Enveloppe convexe

1. Introduction

On s'intéresse au calcul de l'enveloppe convexe d'un ensemble (noté E) de points du plan, c'est-à-dire le plus petit domaine convexe du plan qui contient tous les points de l'ensemble. Il s'agit en fait d'un polygone convexe dont les sommets sont des points de l'ensemble, comme l'illustre la figure suivante.

Figure pleine page

Soient P_i les sommets du polygone constituant l'enveloppe convexe, ordonnés de manière que le parcours se fasse dans le sens horaire. Pour tout côté P_iP_i+1 du polygone, parcouru dans le sens $P_{i} \to P_{i + 1}$ , tous les autres points de E sont situés à droite de la droite orienté (P_iP_i+1), ou sur cette droite. Autrement dit, pour tout point M de l'ensemble E, on doit avoir l'inégalité :

(\vec{P_{i} M} \land \vec{P_{i} P_{i + 1}}) \cdot \vec{u_{z}} \geq 0 (1)

Nous allons voir deux algorithmes de calcul de l'enveloppe convexe qui utilisent cette propriété.

2. Algorithme direct

La première idée consiste à rechercher directement les paires de points dans l'ensemble E qui vérifient la propriété $(1)$ . On explore donc toutes les paires de l'ensemble. Pour chacune, si tous les autres points de l'ensemble vérifie la condition $(1)$ alors la paire est ajoutée à une liste de côtés.

Il faut pour finir générer le polygone à partir de la liste des côtés.

On voit aisément que le temps de calcul de cet algorithme est O(N³), ce qui le rend pratiquement inutilisable pour les grands ensembles.

3. Algorithme incrémental

Un algorithme incrémental consiste à construire le polygone convexe point par point, en le suivant dans le sens horaire. On commence par ordonner les points de l'ensemble E par abscisse croissante. Notons M₀, M₁,... M_N-1 les points de l'ensemble E ainsi ordonnés. Si deux points ont la même abscisse, on doit les classer par ordonnée croissante. On commence à générer le polygone convexe en partant du point le plus à gauche, c'est-à-dire M₀. Le second point du polygone est M₁ (provisoirement). On commence donc par le polygone P=(M₀,M₁). Dans un premier temps, on cherche à construire la partie supérieure du polygone (car le parcours est dans le sens horaire). Lorsque le dernier point ajouté est M_n, on tente d'ajouter le point M_n+1. On essaye donc le polygone à 3 points P=(M₀,M₁,M₂).

Figure pleine page

Dans le cas montré sur la figure, ce polygone n'est pas correct puisque le point M₂ n'est pas à droite du segment [M₀M₁]. On doit alors enlever l'avant dernier point, ce qui conduit à P=(M₀,M₂). Le point suivant à ajouter est M₃, comme le montre la figure suivante :

Figure pleine page

Là encore, le test $(1)$ échoue, ce qui nous oblige à retirer l'avant dernier point M₂, puis à ajouter M₄. Le polygone P=(M₀,M₃,M₄) n'est toujours pas correct. On enlève donc M₃ et on ajoute M₅. Les trois derniers points du polygone P=(M₀,M₄,M₅) vérifient bien la condition $(1)$ ; il est donc accepté. Le point suivant à ajouter est M₆ :

Figure pleine page

On répète la procédure précédente, consistant à retirer de P l'avant dernier point si la condition $(1)$ n'est pas respectée, puis à ajouter le point situé juste à droite du dernier point ajouté. Dans le cas présent, cela nous conduit à (M₀,M₄,M₆,M₉). La procédure est répétée tant que la liste contient au moins deux points. Les points M₆ et M₈ sont successivement éliminés et on arrive à P=(M₀,M₄,M₉) qui est accepté. Le polygone suivant P=(M₀,M₄,M₉,M₁₀) est accepté. Le dernier point ajouté est le dernier point de l'ensemble E; on a donc obtenu la partie supérieure de l'enveloppe.

Figure pleine page

Pour continuer la génération de l'enveloppe, il faut générer sa partie inférieure en partant du point M₁₀. On procède comme précédemment, mais en ajoutant à chaque étape le point situé juste à gauche du dernier point ajouté.

4. Implémentation

4.a. Génération aléatoire des points

Générer un ensemble E de N points, disposés aléatoirement dans le domaine $[0, 1] \times [0, 1]$ . Utiliser pour cela la fonction random.random(), qui renvoie un nombre à virgule flottante compris dans l'intervalle [0,1[, généré aléatoirement selon une loi uniforme. Les points seront stockés sous la forme [x,y] dans une liste nommée points, que l'on pourra générer avec la fonction append.

Tracer les points avec la fonction matplotlib.pyplot.plot.

4.b. Algorithme direct

[1] Écrire une fonction point_segment(P1,P2,M) qui prend en argument trois points P₁,P₂,M et renvoie True lorsque le point M est à gauche du segment [P₁P₂], c'est-à-dire lorsque la condition $(1)$ n'est pas satisfaite.

[2] Écrire une fonction enveloppe_direct(points) qui détermine l'enveloppe avec la méthode directe, pour une liste de points fournie en argument, et qui la renvoie sous forme d'une liste de côtés. Chaque côté est de la forme [P1,P2] où P1 et P2 sont deux points. Dans un premier temps, on ne cherchera pas à générer le polygone à partir de ses côtés.

[3] Tester la fonction avec la liste de points aléatoires et tracer sur une même figure les points et les côtés de l'enveloppe.

[4] À faire en fin de TP : concevoir et implémenter un algorithme permettant de construire la liste des sommets du polygone dans le sens horaire, à partir de la liste de ses côtés.

4.c. Algorithme incrémental

L'enveloppe sera générée dans une pile. L'étape élémentaire de l'algorithme nécessite de travailler sur les trois derniers points de la pile, qu'il faudra donc au préalable retirer avec la fonction pop.

Pour le tri des points, utiliser la fonction sort, qui permet de trier une liste d'objets de la manière suivante :

def cle(objet):
    return ...
liste.sort(key=cle)

La fonction cle renvoie la valeur de la clé utilisée pour le tri. Dans le cas présent, il s'agit de l'abscisse du point. Pour trier la liste des points par ordonnée croissante puis par abscisse croissante, on écrira donc :

def cle(P):
    return P[1]
points.sort(key=cle)
def cle(P):
    return P[0]
points.sort(key=cle)

[5] Écrire une fonction enveloppe_inc(points) qui renvoie le polygone enveloppe convexe sous la forme d'une liste de points.

[6] Tester la fonction avec la liste de points aléatoires et tracer sur la même figure les points et le polygone.

5. Solution

from matplotlib.pyplot import *
import random

Génération aléatoire des points :

points = []
N = 10
for i in range(N):
    points.append([random.random(),random.random()])

Fonction pour déterminer de quel côté un point se situe par rapport à un segment orienté :

def point_segment(P1,P2,P):
    return (P2[0]-P1[0])*(P[1]-P1[1])-(P2[1]-P1[1])*(P[0]-P1[0])>0

Détermination des côtés de l'enveloppe convexe par la méthode directe :

def enveloppe_direct(points):
    enveloppe = []
    N = len(points)
    for i in range(N):
        for j in range(N):
            if i!=j:
                bord = True
                for k in range(N):
                    if k!=i and k!=j:
                        if point_segment(points[i],points[j],points[k]):
                            bord = False
                            break
                if bord:
                    enveloppe.append([points[i],points[j]])
    return enveloppe

Tracé des points et de l'enveloppe convexe :

figure()
axis([0,1,0,1])
axes().set_aspect('equal')
for i in range(len(points)):
    P = points[i]
    plot([P[0]],[P[1]],"ko")
env = enveloppe_direct(points)
for i in range(len(env)):
    segment = env[i]
    P1 = segment[0]
    P2 = segment[1]
    plot([P1[0],P2[0]],[P1[1],P2[1]],"k-")

Figure pleine page

Détermination de l'enveloppe convexe (pile de points) par la méthode incrémentale :

def enveloppe_inc(points):
    def cle(P):
        return P[1]
    points.sort(key=cle)
    def cle(P):
        return P[0]
    points.sort(key=cle)
    N=len(points)
    enveloppe = [points[0],points[1]]
    for i in range(2,N):
        enveloppe.append(points[i])
        valide = False
        while not(valide) and len(enveloppe)>2:
            P3 = enveloppe.pop()
            P2 = enveloppe.pop()
            P1 = enveloppe.pop()
            if point_segment(P1,P2,P3):
                enveloppe.append(P1)
                enveloppe.append(P3)
            else:
                enveloppe.append(P1)
                enveloppe.append(P2)
                enveloppe.append(P3)
                valide = True
    enveloppe.append(points[N-2])
    for i in range(N-3,-1,-1):
        enveloppe.append(points[i])
        valide = False
        while not(valide) and len(enveloppe)>2:
            P3 = enveloppe.pop()
            P2 = enveloppe.pop()
            P1 = enveloppe.pop()
            if point_segment(P1,P2,P3):
                enveloppe.append(P1)
                enveloppe.append(P3)
            else:
                enveloppe.append(P1)
                enveloppe.append(P2)
                enveloppe.append(P3)
                valide = True
    return enveloppe

Tracé des points et de l'enveloppe convexe :

env=enveloppe_inc(points)
figure()
axis([0,1,0,1])
axes().set_aspect('equal')
for i in range(len(points)):
    P = points[i]
    plot([P[0]],[P[1]],"ko")
env = enveloppe_inc(points)
for i in range(1,len(env)):
    P2 = env[i]
    P1 = env[i-1]
    plot([P1[0],P2[0]],[P1[1],P2[1]],"k-")

Figure pleine page