Opérations sur les polyèdres

4 Opérations sur les polyèdres

Pour les opérations sur les polyèdres SPPoC repose presque entièrement sur la bibliothèque polyédrique connue sous le nom de PolyLib [23]. À l'origine, la PolyLib a été développée à l'IRISA de Rennes par Doran Wilde dans le cadre du projet Alpha . Une seconde version de la PolyLib a vu le jour en collaboration avec Philippe Clauss et Vincent Loechner de l'ICPS à Strasbourg [5]. Cette version permet de manipuler des polyèdres paramétrés, en particulier de trouver leurs sommets et de compter le nombre de points entiers qu'ils contiennent.

L'interface entre SPPoC et la PolyLib utilise le fait que des routines écrites en C peuvent être utilisées en Objective Caml par simple édition de liens. Le plus difficile est d'écrire les fonctions de conversion des structures Caml en structures C et vice-versa.

4.1 Construction des polyèdres

La PolyLib manipule en fait des unions de polyèdres, chaque polyèdre étant défini par un couple formé de ses contraintes et d'une représentation duale à base de sommets, de rayons et de droites. Ces informations sont redondantes, la forme duale pouvant être obtenue à partir des contraintes et réciproquement, mais elles permettent d'effectuer les calculs plus rapidement. Pour la même raison SPPoC manipule des listes de couples.

Il est possible de créer un polyèdre directement à partir des deux informations mais c'est déconseillé car aucun contrôle de cohérence n'est effectué. On peut cependant vouloir utiliser cette possibilité pour des raisons de performance. Il est beaucoup plus pratique de créer un polyèdre à partir d'un simple système de contraintes :

# let p = Polyhedron.of_system <:s<3*i+j>=10, 1<=i<=n, 1<=j<=m>>;;
val p : SPPoC.Polyhedron.t =
  {({j <= m, j >= 1, i <= n, i >= 1, 3*i+j >= 10},
    {ray(m), ray(n), vertex(i+7*j+7*m+n),
     ray(j+m), ray(i+n), vertex(3*i+j+m+3*n)})}

La présence de rayons dans la forme duale (calculée à partir des contraintes) nous apprend qu'il ne s'agit pas d'un polyèdre borné. En effet, le translaté d'un polyèdre suivant la direction d'un de ses rayons est inclus dans le polyèdre original.

Il est aussi possible de créer un polyèdre en spécifiant sa forme duale. Prenons l'exemple d'une bande diagonale dans un espace à deux dimensions dont les bords passent par les points (0,1) et (1,0). Sa représentation duale est constituée de ces deux sommets et de la ligne de vecteur directeur (1,1) :

# let p = Polyhedron.of_rays <:r<vertex(i), vertex(j), line(i+j)>>;;
val p : SPPoC.Polyhedron.t =
  {({i-j <= 1, i-j >= -1}, {line(i+j), vertex(j), vertex(-j)})}

La représentation de p fait bien apparaître les deux inéquations définissant la bande.

De plus SPPoC permet la création de polyèdres ne contenant aucun point ou tous les points d'un espace défini par le nom des variables associées à ses dimensions. Enfin, des fonctions sont disponibles pour créer directement des unions de polyèdres et pour récupérer la liste des contraintes ou la liste des représentations duales de telles unions. Les unions de polyèdres sont représentées par une liste de polyèdres. La calculatrice SPPoC présente toujours des unions préalablement simplifiées par la PolyLib : les polyèdres des listes sont disjoints.

4.2 Opérations ensemblistes

SPPoC reprend les opérations sur les polyèdres offertes par la PolyLib : intersection, union, soustraction, image ou pré-image par une fonction affine, calcul de l'enveloppe convexe, simplification dans le contexte d'une autre union de polyèdres, test d'inclusion et enfin test de vacuité.

Il est important de comprendre que SPPoC construit les polyèdres en ne tenant compte que des variables qui interviennent dans les contraintes ou la représentation duale. Supposons qu'un utilisateur souhaite calculer l'intersection entre une droite horizontale dans le plan et un polyèdre plus complexe de dimension deux. La définition de la droite peut se faire comme suit :

# let d = Polyhedron.of_system <:s< i=4 >>;;
val d : SPPoC.Polyhedron.t = {({i = 4}, {vertex(4*i)})}

Ceci définit un polyèdre de dimension 1 (il est possible de s'en convaincre en demandant la liste des variables ; elle se limite à i). Le convexe de dimension deux peut, lui, se définir par :

# let p = Polyhedron.of_system <:s<2*i+j>1, i-4*j<10, j<=20>>;;
val p : SPPoC.Polyhedron.t =
  {({j <= 20, i-4*j <= 10, 2*i+j >= 1},
    {vertex((14*i-19*j)/9), vertex((-19*i+40*j)/2), vertex(90*i+20*j)})}

Lors du calcul d'intersection de la droite d et du polyèdre p, SPPoC se rend compte que les deux polyèdres n'ont pas la même dimension et cherche à les plonger dans un espace commun. En l'occurrence un espace à deux dimensions (associées aux variables i et j). Le résultat de l'intersection est un segment de droite :

# Polyhedron.inter d p;;
- : SPPoC.Polyhedron.t =
{({i = 4, j <= 20, 2*j >= -3}, {vertex((8*i-3*j)/2), vertex(4*i+20*j)})}

Dans tous les cas, il est possible de trouver un espace commun : SPPoC récupère les ensembles de variables des unions de polyèdres, en calcule l'union V et plonge les deux unions de polyèdres dans un espace ayant autant de dimensions que V a d'éléments. L'identification des dimensions communes se fait par égalité des noms de variables les caractérisant. Il n'y a pas de renommage automatique.

Pour les images et pré-images de polyèdres par une fonction affine, un mécanisme semblable est mis en place. Pour comprendre l'utilité de l'auto-ajustement en présence de fonctions, prenons l'exemple de la projection du polyèdre ci-dessous selon la direction associée à i :

# let q = Polyhedron.of_rays
          <:r< vertex(i+5j), vertex(2i-j),
               vertex(3j), vertex(2i) >>;;
val q : SPPoC.Polyhedron.t =
  {({i <= 2, 5*i+j <= 10, 2*i-j >= -3, 2*i+j >= 3},
    {vertex(i+5*j), vertex(2*i-j), vertex(3*j), vertex(2*i)})}

La fonction à utiliser est <:fct< i <- i >>. ³ Notre fonction de projection se retrouve donc avec un espace de départ et d'arrivée à une dimension. SPPoC va ajuster automatiquement la dimension de cet espace lors du calcul de l'image :

# Polyhedron.image q <:fct< i <- i >;;
- : SPPoC.Polyhedron.t = {({i >= 0, i <= 2}, {vertex(0), vertex(2*i)})}

Si, au contraire, on souhaite plonger q dans un espace avec un plus grand nombre de dimensions (par exemple, en utilisant <:fct< i,j,k <- i,j,k >>), SPPoC va automatiquement ajuster la dimension de q avant d'appliquer la fonction. L'utilisateur n'a donc pas, en principe, à se préoccuper des dimensions de ses fonctions et de ses polyèdres. Cependant SPPoC permet de forcer le plongement d'une fonction ou d'une union de polyèdres dans un espace de dimension supérieure.

Par contre, l'utilisateur doit absolument utiliser les mêmes variables pour référencer les dimensions des espaces dans lesquels sont définis ses polyèdres. Si pour une raison ou une autre, ce n'est pas possible lors de la définition des polyèdres, SPPoC permet, par la suite, de renommer les variables liées aux dimensions.

4.3 Opérations sur les polyèdres paramétriques

Même sans opération spécifique, la PolyLib permet de manipuler des polyèdres paramétriques en ajoutant des dimensions supplémentaires pour les paramètres. Par exemple il est possible de définir un carré paramétré par la taille n de ses cotés :

# let carre = Polyhedron.of_system <:s< 0<=i<=n, 0<=j<=n >>;;
val carre : SPPoC.Polyhedron.t =
  {({j <= n, j >= 0, i <= n, i >= 0},
    {vertex(0), ray(n), ray(i+n), ray(i+j+n), ray(j+n)})}

Il est possible de faire des opérations ensemblistes sur ces polyèdres paramétriques. Par contre, la représentation duale de ces objets est déconcertante ; elle ne donne plus les sommets des polyèdres bornés. En effet, pour << paramétrer >> le polyèdre, des dimensions représentant les paramètres sont ajoutées. La représentation duale tient compte de ces dimensions et, comme les paramètres sont généralement non bornés, on se retrouve avec des lignes ou des rayons suivant les dimensions paramétriques.

Il existe une version étendue de la PolyLib qui est celle à laquelle SPPoC propose une interface et qui permet de trouver les sommets d'un polyèdre paramétrique et de calculer le nombre de points entiers contenus dans un tel polyèdre. Voici un appel à l'extension de la PolyLib qui donne les sommets d'un polyèdre paramétrique sous sa forme habituelle, c'est à dire eux-mêmes paramétrés :

# Polyhedron.vertices
    carre ( Polyhedron.of_system <:s< n>=0 >> );;
- : (SPPoC.Polyhedron.t * SPPoC.Expr.t list list) list =
[{({n >= 0}, {ray(n), vertex(0)})}, [[n; n]; [n; 0]; [0; n]; [0; 0]]]

On reconnait bien ici les sommets du carré : (0, 0), (0, n), (n, 0) et (n, n), ainsi que la définition du domaine du paramètre n : n³0.

Pour compter le nombre de points entiers contenus dans un polyèdre paramétrique borné⁴ l'extension de la PolyLib utilise les polynômes de Ehrhart. Les polynômes de Ehrhart sont des polynômes à coefficients périodiques. En fait les coefficients sont des couples formés d'une liste de nombres et d'un paramètre. La valeur du coefficient est le nombre de la liste de rang le paramètre modulo la taille de la liste. Intuitivement, ces coefficients périodiques sont utiles pour construire des expressions dépendant de paramètres en évitant une définition au cas par cas suivant les valeurs de ces paramètres. Prenons l'exemple du polynôme nm+[1/2,2/3]_n m + [0,1/3]_m n + 1, sa valeur pour n=2 et m=3 est 2× 3+1/2× 3 + 1/3× 2 +1. Deux exemples de calcul du nombre de points entiers sont donnés ci-après. Les résultats sont, bien entendu, donnés en fonction des paramètres.

# Polyhedron.enum
    carre 
    (Polyhedron.of_system <:s< n>=0 >>);;
- : SPPoC.Polyhedron.expr_quast =
[{({n >= 0}, {ray(n), vertex(0)})}, 1+2*n+n**2]
# let rectangle = 
    Polyhedron.of_system <:s< 0<=i<=n, 0<=j<=n/3 >>;;
val p : SPPoC.Polyhedron.t =
  {(3*j <= n, j >= 0, i <= n, i >= 0,
    vertex(0), ray(n), ray(i+n), ray(3*i+j+3*n), ray(j+3*n))}
# Polyhedron.enum
     rectangle 
    (Polyhedron.of_system <:s< n>=0 >>);;
- : SPPoC.Polyhedron.expr_quast =
[{({n >= 0}, {ray(3*n), vertex(0)})},
 1/3*n**2+n*[4/3, 1, 2/3]n+[1, 2/3, 1/3]n]

L'opération de comptage de points entiers est coûteuse en mémoire et en temps. En appelant directement la fonction de comptage de la PolyLib , seuls des polyèdres de petite dimension peuvent être traités. Pour parvenir à traiter des polyèdres plus complexes, SPPoC offre une fonction de plus haut niveau Polyhedron.enumeration. Cette fonction utilise le fait que la plupart des polyèdres manipulés peuvent s'écrire sous la forme d'un produit cartésien. Un résultat bien connu est que le cardinal d'un produit cartésien est le produit des cardinaux des ensembles utilisés dans le produit. Le fait de lancer le calcul du nombre de points sur plusieurs polyèdres de faible dimension plutôt que sur le polyèdre original permet de limiter l'explosion combinatoire. L'exemple le plus frappant est celui d'un hyper-cube de dimension d ; un tel polyèdre est le produit cartésien de d polyèdres de dimension 1.

L'algorithme de décomposition d'un polyèdre en produit cartésien consiste à partitionner l'ensemble des contraintes définissant le polyèdre. On regroupe deux contraintes dans le même sous-ensemble si et seulement si elles ont au moins une variable en commun. À la partition de l'ensemble des contraintes correspond donc une partition de l'ensemble des variables (toutes les variables apparaissent dans les contraintes puisque le polyèdre doit être borné). Le polyèdre original est le produit cartésien des polyèdres définis par les sous-ensembles de contraintes sur les sous-ensembles de variables correspondants.