Interface between MuPAD and Polylib:
Computation of Ehrhart's polynomials

François Thomasset
INRIA - Rocquencourt - France
April 2003

The Polyhedral Library (PolyLib for short) operates on objects made up of unions of polyhedra of any dimension. It was first developed by Doran Wilde at IRISA, in Rennes, France, in connection with the ALPHA project. It is now available from the site of University Louis Pasteur, at Strasbourg: 

	  http://icps.u-strasbg.fr/polylib/

We use it here to compute the Ehrhart polynomial of a polyhedron defined under MuPAD. The so-called "Ehrhart polynomial" gives a closed form to the number of integer points in a polyhedron defined by a set of inequalities; the polyhedron is allowed to be defined in terms of symbolic parameters. More details can be found on the same site at Strasbourg : 

     http://icps.u-strasbg.fr/Ehrhart/

For instance, from the following input, describing a polyhedron parameterized by p: 

    [x-1>=0, -2*x+p>=0], // the polyhedron
    [-p+80>=0],		 // context
    [x],		 // list of iterators
    [p],		 // list of parameters
    []			 // substitution list for the parameters, e.g. [p=10]
the present interface will return to MuPAD the expression below:
              -- table(                                       --
              |             p                                  |
              |    _value = - + _periodic([0, -1/2], p),       |
	      |             2                                  |
              |    _cond = _cond = 0 <= p - 2 and 0 <= 80 - p  |
              -- )                                            --
This means a list of alternatives, each alternative being implemented as a table: for each alternative, when _cond is evaluated to TRUE, _value yields the expected expression, a function of the parameters. In the present case, there is only one alternative, whose condition is implied by the context. When given to the pretty printing function prt_poly, this gives an "Ehrhart polynomial":
	1/2*p + _periodic([0, -1/2], p)
Here _periodic takes 2 parameters: a list l and a parameter p; the result is the k-th element of the list (counting up from 0), where k=p (modulo length of list). E.g. we get (-1/2) when p = 11.
Note that since MuPAD count the list elements from 1, we have to add 1 to this k to find the correct position, as done by the function _periodic provided in test_count.mupad.

Note that whenever a parameter is given a value in the substitution list, it gets evaluated within the interface.

To install and use this interface:

  1. Install MuPad in some directory MuPadDIR (see the web site http://www.mupad.de).
    Make sure that MuPadDIR/share/bin is in your search path.
    Then saying:
    mupad
    in your shell starts a mupad session.
    Please use a version at least as recent as 2.0.0 of MuPad.
    Don't forget to register yourself as a MuPad user at www.sciface.com/register.shtml.
  2. Decide of a directory where the polylib library is installed ; define in your environment a variable named COMMON pointing to this directory. This might be e.g. /usr/local on a machine where you can become root.
    It is expected that the library reside in the following directory:
    $COMMON/polylib
  3. In your environment, setup your library path to the place where polylib libraries are installed:
    setenv LD_LIBRARY_PATH $COMMON/lib/
  4. Make sure that mupad and mmg are in your execution path.
    (mmg is the dynamic module generator of MuPAD; it will be used to compile the interface).
    Go to the Count subdirectory in this distribution; if a file count.mdm exists, remove it ('make clean'); then type 'make'. A file count.mdm should be generated.
  5. Start mupad.
    File test_count.mupad in this distribution gives an example.
    In MuPAD, define the READPATH variable, so that it contains the place where you installed the Count subdirectory from this distribution: 
       >> READPATH := your_directory."/Count":
    The MuPAD procedure constraints2poly takes a friendly description of the polyhedron, and converts it to a description readable by the interface.
    The arguments of constraints2poly are: The result can be forwarded to the interface:
       >> constraints_ := constraints2poly ([x-1>=0, -2*x+p>=0],
				        [-p+80>=0], [x], [p], []) :
   >> module(count):
   >> count::count (op(constraints_,1), op(constraints_,2),[p]) ;

    The result of count can then be given to prt_poly.
Bug reports, comments: Francois.Thomasset@inria.fr