In order to show how the Convex package works, we will implement Kushnirenko's bound for the sum of the Milnor numbers of the singularities of a complex polynomial (see A. G. Kouchnirenko, Polyèdres de Newton et nombres de Milnor, Invent. Math. 32 (1976), 1–31). But don't worry - you do not need to know what a singularity is, let alone its Milnor number; just believe that it is something interesting and that it can be expressed in terms of convex geometry. We will explain how to do so along the way.

Loading the package

After invoking Maple, we first load the package:
Convex version 1.2.0, Copyright (C) 1999-2016 Matthias Franz This package is distributed under the GNU General Public License See for more information
[affinehull, ambientdim, arecompatible, boundary, codim, contains, containsrelint, convhull, corank, corners, crosspolytope, cube, cyclicpolytope, delaunay, dim, directsum, distance, domain, dotprod, draw, dual, edges, emptypcomplex, emptypolyhedron, emptypolytope, facefan, faces, facets, fan, flagf, flagh, fullcone, fullpolyhedron, furthestdelaunay, fvector, genhpolynomial, genhvector, hilbertbasis, homology, hplanes, hspacenos, hspaces, hvector, image, incidencematrix, incidentfacets, incidentrays, intersection, isaffine, isbounded, iscomplete, iscontained, isempty, isface, isfulldim, islinear, ispointed, ispolytopal, isquasipolytopal, isregular, issimple, issimplex, issimplicial, issimplicial1, join, lensspace, lineality, linearhull, lines, maximal, maximum, minimal, minimum, minkowskisum, modz, newtonpolytope, normalfan, pcomplex, permutahedron, plotdata, polar, poshull, posorthant, pred, preimage, projspace, proximum, randompolytope, rank, raynos, rays, readpoly, recession, regularpart, regularsubdiv, relint, simplicialsubdiv, skeleton, stdsimplex, stellarsubdiv, succ, support, surface, torsion, transversalfan, traverse, traverse2, vertexnos, vertices, volume, voronoi, wprojspace, writepoly, zerocone, zerofan]

The Newton polytope of a polynomial

Kushnirenko's formula uses the Newton polytope of a polynomial p. It is defined as follows: One looks at all monomials appearing in p (with non-zero coefficients) and considers the exponent of each such monomial as an integral point in Euclidean space. Their convex hull is the Newton polytope of p.

Let us illustrate this with an example. We start with the following polynomial:

p := (x-y*z+x^2*y)^3-(z+2*y)^4+(x*y*z)^3;
p := (x^2*y-y*z+x)^3-(z+2*y)^4+x^3*y^3*z^3
The Convex package comes already with the function newtonpolytope to compute Newton polytopes.
P := newtonpolytope(p, [x, y, z]);
P := POLYTOPE(3,3,6,8)
Let us explain the last line: Convex provides the data type POLYTOPE for polytopes. An object of this type stores complete information about the corresponding polytope (essentially the vertices, the facets ad the incidence matrix). Apart from trivial examples, this is much more information than one could reasonably display. Therefore, only a summary of the polytope's properties is shown, namely (in this order) the ambient dimension, the dimension, the number of vertices and of facets. The same applies to other types defined by Convex. For each type, you can find the exact meaning of the short display on the help page for that type. (A click on the type gets you there.)

To know more about a polytope, one uses some of the functions defined by the Convex package. For instance, let us see what the vertices of P are:

[[0, 4, 0], [0, 0, 4], [6, 3, 0], [3, 3, 3], [0, 3, 3], [3, 0, 0]]
This is indeed what we have expected. The vertex [3, 0, 0] corresponds to the monomial x^3 in p because x is the first variable in the ordering we have chosen.

If we are curious, we may also draw the Newton polytope:

[Maple plot]

Kushnirenko's formula

Kushnirenko's formula involves the convex hull of the Newton polytope and the origin. (In other words, one perturbs the constant coefficient of the polynomial if it happens to be zero.) We therefore replace P by this convex hull. This is done by the function convhull.
P := convhull(P, [0, 0, 0]);
P := POLYTOPE(3,3,7,8)
The function convhull actually lies at the heart of the Convex package since it allows to define polytopes as the convex hull of some points and, more generally, all kinds of polyhedra as the convex hull of points, rays, lines and/or previously defined polyhedra. (The function newtonpolytope is based on it as well.) There is also a dual function, called intersection, which computes polyhedra given as the intersection of halfspaces, hyperplanes and/or other polyhedra. We will use it below.

Kushnirenko's bound on the sum of the Milnor numbers is an alternating sum over the volume of the intersection of P with all coordinate hyperplanes in the ambient space. Here "volume" means the induced Euclidean volume in a, say, i-dimensional coordinate hyperplane, which is 0 if the intersection is not i-dimensional. Each such volume is scaled by (-1)^(n-i)*i!, where n is the ambient dimension of P (which is 3 in our example).

Hence, we need to know the volume of the intersection of P with all coordinate hyperplanes. Since the polytope lies in the positive octant (or "positive orthant" in general), we can equally well intersect it with the faces of the positive orthant, considered as a polyhedron. Let's compute these faces.

PO := posorthant(3);
PO := CONE(3,3,0,3,3)
The function posorthant returns the positive orthant in the space of given dimension, considered as a cone. For various reasons, the Convex package distinguishes between cones (with are of type CONE) on the one hand and polytopes (of type POLYTOPE) and other polyhedra (of type POLYHEDRON) on the other.

The following call converts the CONE PO to a POLYHEDRON, which is then stored in the same variable.

PO := convert(PO, affine);
PO := POLYHEDRON(3,3,0,[1, 3],[3])
Next we compute all faces of the positive orthant:
fl := faces(PO);
fl := Array(-1 .. 3,{-1 = [PFACE(0,4)], 0 = [PFACE(1,3)], 1 = [PFACE(2,2), PFACE(2,2), PFACE(2,2)], 2 = [PFACE(3,1), PFACE(3,1), PFACE(3,1)], 3 = [PFACE(4,0)]},datatype = anything,storage = rectangular,order = Fortran_order)
They are returned as an Array whose i-th entry is the list of all i-dimensional faces. The faces of a polyhedron are of type PFACE, not POLYHEDRON. The difference is that a PFACE knows of which polyhedron it is a face. This for example permits calculations in face lattices. (The internal representation of this type is also different.) Note that in the above result there is a face of dimension -1. This is the empty face, which would not be there if we had taken the CONE posorthant(3). In our example, we will not need this face.

Now we run through all non-empty faces of PO and sum up the volumes of their intersections with P:

V := 0: # the sum of all volumes for i from 0 to 3 do Vi := 0; # the sum of the i-dimensional volumes for f in fl[i] do Q := intersection(P, convert(f, POLYHEDRON)); if dim(Q) = i then Vi := Vi+volume(maximal(Q)) fi od; V := V+(-1)^(3-i)*i!*Vi od: V;
The call convert(f, POLYHEDRON) converts a PFACE to a POLYHEDRON, and the function intersection computes the intersection of its arguments. The function volume computes the Euclidean volume of a polytope. By convention, this is 0 if the polytope is not full-dimensional. But this is not what we want. We use volume(maximal(Q)) instead, which gives the volume in the affine plane generated by Q. This is the i-dimensional volume we are interested in if Q is of dimension i. Summing up all volumes gives the desired result.

The new function

We collect the code we have written and make a new function out of it. Moreover, we make the function independent of the dimension and ambient dimension of the Newton polytope. Here is the result:
kushnirenko_bound := proc(p, indets) local n, P, fl, i, V, Vi, f, Q; n := nops(indets); # = ambient dimension of Newton polytope P := convhull(newtonpolytope(p, indets), [0$n]); fl := faces(convert(posorthant(n), affine)); V := 0; for i from 0 to dim(P) do Vi := 0; for f in fl[i] do Q := intersection(P, convert(f, POLYHEDRON)); if dim(Q) = i then Vi := Vi+volume(maximal(Q)) fi od; V := V+(-1)^(n-i)*i!*Vi od; V end:
Let's try it out with the polynomial p from before:
kushnirenko_bound(p, [x, y, z]);
Now in higher dimensions:
kushnirenko_bound((x+y+z+u+v+w)^3, [x, y, z, u, v, w]);