Poly 435
Poly 435
Poly 435
31 janvier 2023
Table des matières
3 ALGORITHMES D’OPTIMISATION 59
3.1 Algorithmes de type gradient (sans contraintes) . . . . . . . . . . . . 60
3.1.1 Algorithme de gradient à pas optimal . . . . . . . . . . . . . . 60
3.1.2 Algorithme de gradient à pas fixe . . . . . . . . . . . . . . . . 63
3.2 Généralisations et autres algorithmes de type gradient . . . . . . . . . 66
3.2.1 Vitesse de convergence . . . . . . . . . . . . . . . . . . . . . . 67
3.2.2 Algorithme de Nesterov . . . . . . . . . . . . . . . . . . . . . 71
3.2.3 Algorithme de la boule pesante . . . . . . . . . . . . . . . . . 78
3.2.4 Algorithme du gradient conjugué . . . . . . . . . . . . . . . . 80
3.2.5 Algorithme de sous-gradient . . . . . . . . . . . . . . . . . . . 83
3.2.6 Gradient stochastique . . . . . . . . . . . . . . . . . . . . . . 86
3.2.7 Algorithme proximal . . . . . . . . . . . . . . . . . . . . . . . 89
3.3 Méthode de Newton . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
3.3.1 Cas de la dimension finie . . . . . . . . . . . . . . . . . . . . . 93
i
ii TABLE DES MATIÈRES
Préface
G. Allaire, A. Ern
Paris, le 1er février 2023
vi Préface
Chapitre 1
INTRODUCTION À
L’OPTIMISATION ET AU
CONTRÔLE
Ce chapitre est une introduction aux deux sujets de ce cours qui, quoique dif-
férents et indépendants, partagent de nombreux points communs. L’objectif est de
présenter les motivations à ces deux sujets et d’illustrer leur portée et leur diversité
à travers de nombreux exemples concrets.
1.1 Motivations
L’optimisation et le contrôle sont des sujets très anciens qui connaissent un nou-
vel essor depuis l’apparition des ordinateurs et dont les méthodes s’appliquent dans
de très nombreux domaines : économie, gestion, planification, logistique, automa-
tique, robotique, conception optimale, sciences de l’ingénieur, traitement du signal,
etc. L’optimisation et le contrôle couvrent ainsi un champ scientifique relativement
vaste, qui touche aussi bien au calcul des variations qu’à la recherche opérationnelle
(en lien avec les processus de gestion ou de décision). Nous ne ferons souvent qu’ef-
fleurer ces sujets car il faudrait un polycopié complet pour chacun d’eux si nous
voulions les traiter à fond.
Avant de considérer l’optimisation ou le contrôle d’un phénomène physique ou
d’un système industriel, il faut déjà passer par une étape de modélisation qui
permet de représenter cette réalité (et éventuellement de la simplifier si elle est trop
complexe) par un modèle mathématique. Dans ce qui suit, nous considérerons des
exemples de problèmes d’optimisation et de contrôle où les modèles peuvent être
de nature très différente. Dans le cas le plus simple des problèmes d’optimisation,
le modèle sera une simple équation algébrique et il s’agira simplement d’optimiser
une fonction définie sur un espace de dimension finie (disons Rn ). Une deuxième
catégorie de problèmes correspond au cas où la fonction à optimiser dépend de la
solution d’une équation différentielle ordinaire (autrement dit, cette fonction est
définie sur un espace de dimension infinie, par exemple l’espce C[0, T ] des fonctions
1
2 CHAPITRE 1. INTRODUCTION
continues sur l’intervalle fermé [0, T ] en temps). On parle alors de contrôle (ou de
commande) optimale, et les applications sont très nombreuses en automatique et
robotique. La troisième et dernière catégorie correspond à l’optimisation de fonctions
de la solution d’une équation aux dérivées partielles. Il s’agit alors de la théorie
du contrôle optimal des systèmes distribués qui a de nombreuses applications, par
exemple en conception optimale ou pour la stabilisation de structures mécaniques.
Il ne nous sera pas possible dans ce cours de niveau introductif d’aborder le cas du
contrôle des systèmes distribués. Remarquons néanmoins que ces catégories ne sont
pas hermétiquement cloisonnées puisqu’après discrétisation spatiale une équation
aux dérivées partielles se ramène à un système d’équations différentielles ordinaires
et, qu’après discrétisation temporelle, une équation différentielle ordinaire se ramène
à un système d’équations algébriques.
On peut aussi séparer l’optimisation en deux grandes branches aux méthodes
fort différentes selon que les variables sont continues ou discrètes. Typiquement, si
l’on minimise une fonction f (x) avec x ∈ Rn , il s’agit d’optimisation en variables
continues, tandis que si x ∈ Zn on a affaire à de l’optimisation combinatoire ou
en variables discrètes. Malgré les apparences, l’optimisation en variables continues
est souvent plus “facile” que l’optimisation en variables discrètes car on peut utiliser
la notion de dérivée qui est fort utile tant du point de vue théorique qu’algorith-
mique. L’optimisation combinatoire est naturelle et essentielle dans de nombreux
problèmes de la recherche opérationnelle. C’est un domaine où, à coté de résultats
théoriques rigoureux, fleurissent de nombreuses “heuristiques” essentielles pour ob-
tenir de bonnes performances algorithmiques. Dans ce cours de niveau introductif
nous traiterons majoritairement d’optimisation continue et nous renvoyons au cours
de troisième année [5] pour l’optimisation combinatoire.
Quant aux problèmes de contrôle, on peut également en distinguer deux branches
principales selon que l’objectif soit d’amener le système en un état fixé (on parle
alors de contrôlabilité ou de commandabilité), ou de minimiser une fonction-
nelle évaluant le coût de l’action sur le système. Ce coût résulte bien souvent d’un
compromis entre la réalisation de certaines performances (comme l’atteinte d’une
cible ou le fait de s’en rapprocher) et le coût afférent à la réalisation de ce contrôle
(et dû par exemple à la consommation énergétique, à la poussée des moteurs, etc.)
On parle alors de problèmes de contrôle optimal.
Pour finir cette brève introduction nous indiquons le plan de la suite du cours. Le
reste de ce chapitre présente à travers des exemples les applications de l’optimisation
et du contrôle. Les Chapitres 2, 3 et 4 sont consacrés aux problèmes d’optimisation.
Le Chapitre 2 va principalement porter sur les aspects théoriques. On y étudiera la
question de l’existence et de l’unicité de solutions à des problèmes d’optimisation,
que ce soit en dimension finie ou infinie. En particulier, nous verrons le rôle crucial
de la convexité pour obtenir des résultats d’existence en dimension infinie. Nous
y verrons aussi les conditions d’optimalité, reposant sur les dérivées des fonctions
optimisées, qui permettent de caractériser les solutions possibles. Le Chapitre 3
développera les algorithmes numériques qui découlent de ces conditions d’optimalité.
Le Chapitre 4 traite de la programmation linéaire dans le cas continu ou discret.
Dans ce dernier cas, cela constitue une très brève, et très biaisée, introduction aux
1.2. EXEMPLES EN OPTIMISATION 3
N
M X
!
X
inf cij vij
(vij )
i=1 j=1
Une variante consiste à autoriser des valeurs réelles positives de aij ∈ R+ . Ce type
de problèmes est appelé problème d’affectation (il intervient dans des contextes
industriels plus sérieux comme l’affectation des équipages et des avions dans une
compagnie aérienne). Bien que ce ne soit pas forcément la meilleure manière de
poser le problème, on peut l’écrire sous une forme voisine de l’Exemple 1.2.1. Les
variables de décision sont notées vij qui vaut 1 s’il y a mariage entre la femme i et
l’homme j et 0 sinon. On veut maximiser
N X N
!
X
sup aij vij
(vij ) i=1 j=1
sous les contraintes (qui assurent que chaque femme trouvera un mari et chaque
homme une épouse)
N
X N
X
vij = 0 ou 1, vij = 1, vij = 1 pour 1 ≤ i, j ≤ N.
j=1 i=1
est plus délicat mais peut encore se résoudre en le voyant comme un problème de
flot sur un graphe (voir la Section 4.3). Nous ne faisons qu’évoquer ces techniques
dans ce cours et nous renvoyons au cours de troisième année [5] pour plus de détails.
•
c’est-à-dire qui maximise l’utilité sous une contrainte de budget maximal (voir la
Sous-section 2.6.2 pour la résolution). •
Voici maintenant un exemple issu du domaine de l’apprentissage machine (ou
machine learning) qui connait un développement spectaculaire ces dernières années.
Comme pour l’Exemple 1.2.4 des moindres carrés, il est d’usage de régulariser le
vecteur des paramètres w et, à l’aide d’un coefficient ` > 0, de considérer une
version augmentée de (1.4)
n
1X `
inf P (hw,τ (xi ), yi ) + kwk22 (1.5)
w∈R ,τ ∈R n
d 2
i=1
où k.k2 est la norme euclidienne. On peut aussi remplacer cette norme par une
norme l1 si l’on souhaite trouver des vecteurs w creux. Il s’agit donc d’un problème
d’apprentissage supervisé (puisqu’on dispose de labels yi fournis par un expert pour
apprendre le modèle à partir des données xi ) dont la solution optimale (w∗ , τ ∗ ) donne
un algorithme de classification, c’est-à-dire
de prédiction
du label pour une nouvelle
donnée x grâce à la formule y = sgn hw∗ ,τ ∗ (x) .
Le problème (1.5) peut aussi s’interpréter comme une séparation des données en
deux classes les plus éloignées possibles au sens de la distance Euclidienne dans Rd :
le rapport 2/kwk2 est la distance entre les deux hyperplans affines w · x − τ = +1 et
w · x − τ = −1 qui séparent les données si celles-ci vérifient yi (w · xi − τ ) ≥ 1. Dans ce
cas, cette technique relève de ce que l’on appelle les machines à vecteurs de support
ou séparateurs à vaste marge (en anglais, support vector machine, ou SVM). •
Exemple 1.2.8 (calcul des variations) Pour fixer les idées, considérons une poutre
uni-dimensionnelle qui au repos est représentée par le segment de droite (0, L) où
L > 0 est la longueur de la poutre. On note x ∈ (0, L) les points de ce segment. Sous
l’action d’une force f (x), le point x de la poutre se déplace perpendiculairement au
segment d’une distance u(x). On peut trouver la position d’équilibre de la poutre
en résolvant un problème de minimisation d’énergie. L’énergie mécanique totale se
décompose en deux termes : d’une part une énergie de déformation (le terme quadra-
tique en la dérivée u0 (x) ci-dessous), d’autre part une énergie potentielle (le terme
proportionnel à la force ci-dessous). Autrement dit, il s’agit du problème d’optimi-
sation
1 L
Z Z L
0 2
inf µ|u (x)| dx − f (x) u(x) dx,
u∈C 1 (0,L) 2 0 0
où µ > 0 est un coefficient de rigidité de la poutre. Si la poutre est encastrée à une de
ses extrémités (ou aux deux) on rajoutera les contraintes u(0) = 0 et/ou u(L) = 0.
Plus généralement, on peut vouloir résoudre des problèmes du type
Z L
0
inf
1
j u (x), u(x) dx (1.6)
u∈C (0,L) 0
où j(λ, u) est une fonction régulière sur R × R. Nous verrons que l’existence d’une
solution au problème (1.6) n’est absolument pas une évidence et qu’il faut des condi-
tions particulières sur l’intégrande j pour s’en assurer. Par ailleurs, la définition de
l’espace sur lequel on minimise est aussi une question très délicate dont nous ne
dirons rien ici et qui conduit à des développements importants en analyse (disons
seulement que l’espace C 1 (0, L) doit être remplacé par des espaces plus généraux,
de type Sobolev, utilisant la théorie des distributions, voir [1], [9], [15]). •
Ce système de contrôle se récrit sous la forme (1.7) avec f (t, x, u) = Ax+Bu, si bien
qu’il y a d = 2 degrés de liberté et k = 1 variable de contrôle. En guise de condition
initiale, on prescrit la position et la vitesse du tram à t = 0. Par exemple, le tram
est initialement à l’arrêt (V (0) = 0) à la gare située en X(0) = 0. En se fixant
un horizon temporel T > 0, le problème de la contrôlabilité du tram consiste à se
demander s’il existe un contrôle u : [0,T ] → R capable d’amener le tram au temps
T en une position X1 avec une vitesse V1 (par exemple, X1 peut être la position
de la gare suivante, et dans ce cas on prescrit V1 = 0). Nous verrons dans ce cours
que la réponse à cette question est toujours positive. On peut également poser cette
question en rajoutant une contrainte sur le contrôle, par exemple que celui-ci soit en
tout temps à valeurs dans un ensemble compact, comme [−1, 1], ce qui décrit le fait
que les moteurs du tram ne peuvent fournir qu’une accélération bornée. On peut
également, tout en conservant cette contrainte sur le contrôle, chercher à atteindre
la cible le plus rapidement possible. Par exemple, s’il s’agit d’amener le tram d’une
gare à la suivante, on s’attend à ce qu’un première phase d’accélération soit suivie
par une phase de décélération jusqu’à l’arrêt complet du tram à la prochaine gare.
•
Remarquons que la fonction x(t) dépend de la variable u à travers (1.8). Comme les
commandes admissibles sont éventuellement limitées (la puissance d’un moteur est
souvent bornée...), on introduit un convexe fermé non vide K de Rk qui représente
l’ensemble des contrôles admissibles. Le problème est donc de résoudre
inf J(u).
u(t)∈K, t∈[0,T ]
Il faudra, bien sûr, préciser dans quels espaces fonctionnels on minimise J(u) et on
définit la solution x de (1.8) (voir le Chapitre 6 pour la résolution). •
c(Y ) u `
(0, 0) X
l’axe des X coïncide avec la berge de départ et l’axe des Y est transverse au canal.
La barque a une vitesse d’amplitude constante notée v, et le courant a une vitesse
c(Y ). On suppose que c(Y ) > v, pour tout Y ∈ [0, `] ; il s’agit d’une hypothèse dite
de courant fort car l’amplitude du courant est toujours supérieure à la vitesse de la
barque. La configuration est illustrée à la figure 1.1. Le contrôle est l’angle u ∈ [0, 2π]
de la vitesse de la barque par rapport à l’axe des X, la vitesse étant considérée dans
le repère du courant. L’état de la barque est décrit par le couple x := (X, Y ) ∈ R2
donnant les coordonnées de la barque dans le repère fixé. La trajectoire de la barque
est régie par la dynamique suivante :
v cos(u(t)) + c(Y (t))
ẋ(t) = f (x(t), u(t)) = , ∀t ∈ [0,T ], (1.9)
v sin(u(t))
et la condition initiale est x(0) = (0, 0). Nous pouvons considérer (entre autres) deux
problèmes de contrôle optimal pour atteindre la berge opposée :
1. Minimiser le déport latéral : le critère à minimiser est
J1 (u) = x(T ),
et on rajoute la contrainte de cible Y (T ) = `. On obtient le problème de
minimisation sous contraintes
inf J1 (u).
u(t)∈[0,2π], t∈[0,T ], Y (T )=`
On notera que dans ces deux problèmes, le temps T n’est pas fixé et qu’il s’agit d’une
inconnue supplémentaire. Ces deux problèmes peuvent être résolus en appliquant le
principe du minimum de Pontryaguine que nous verrons dans ce cours. •
12 CHAPITRE 1. INTRODUCTION
Chapitre 2
ASPECTS THÉORIQUES DE
L’OPTIMISATION
13
14 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
est l’ensemble des éléments admissibles du problème, ou bien que K définit les
contraintes s’exerçant sur le problème considéré. Enfin, le critère, ou la fonction
coût, ou la fonction objectif, à minimiser, noté J, est une fonction définie sur K
à valeurs dans R. Le problème étudié sera donc noté
Lorsque l’on utilise la notation inf pour un problème de minimisation, cela indique
que l’on ne sait pas, a priori, si la valeur du minimum est atteinte, c’est-à-dire s’il
existe v ∈ K tel que
J(v) = inf J(v).
v∈K⊂V
Si l’on veut indiquer que la valeur du minimum est atteinte, on utilise de préférence
la notation
min J(v),
v∈K⊂V
mais il ne s’agit pas d’une convention universelle (quoique fort répandue). Pour
les problèmes de maximisation, les notations sup et max remplacent inf et min,
respectivement. Précisons quelques définitions de base.
Définition 2.1.1 On dit que u est un minimum (ou un point de minimum) local
de J sur K si et seulement si
Lemme 2.1.3 Pour tout problème d’optimisation, si K est non vide, il existe au
moins une suite minimisante.
Démonstration. Notons m ∈ R ∪ {−∞} la valeur de l’infimum de J sur K.
Par définition, pour tout v ∈ K, J(v) ≥ m. Supposons qu’il n’existe aucune suite
(un )n∈N ∈ K telle que limn→+∞ J(un ) = m. Cela revient à dire qu’il existe > 0 tel
que, pour tout v ∈ K, J(v) ≥ m + . Mais cela est une contradiction avec le fait que
m est définie comme la plus grande constante qui minore J sur K. 2
2.2. OPTIMISATION EN DIMENSION FINIE 15
Remarque 2.2.2 Notons que la propriété (2.2), qui assure que toute suite minimi-
sante de J sur K est bornée, est automatiquement vérifiée si K est borné. Lorsque
l’ensemble K n’est pas borné, cette condition exprime que, dans K, J est infinie à
l’infini. •
Exercice 2.2.1 Montrer par des exemples que le fait que K est fermé ou que J est
continue est en général nécessaire pour l’existence d’un minimum. Donner un exemple
de fonction continue et minorée de R dans R n’admettant pas de minimum sur R.
Exercice 2.2.2 Montrer que l’on peut remplacer la propriété “infinie à l’infini” (2.2)
par la condition plus faible
inf J(v) < lim inf J(v) .
v∈K R→+∞ kvk≥R
Exercice 2.2.3 Montrer que l’on peut remplacer la continuité de J par la semi-
continuité inférieure de J définie par
Exercice 2.2.4 Montrer qu’il existe un minimum pour les Exemples 1.2.1 et 1.2.3.
16 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
Exemple 2.3.1 Soit l’espace de Hilbert (de dimension infinie) des suites de carré
sommable dans R
( +∞
)
X
`2 (R) = x = (xi )i≥1 tel que x2i < +∞ ,
i=1
P+∞
muni du produit scalaire hx, yi = i=1 xi yi . On considère la fonction J définie sur
`2 (R) par
+∞ 2
2 X
2 xi
J(x) = kxk − 1 + .
i=1
i
pour lequel nous allons montrer qu’il n’existe pas de point de minimum. Vérifions
tout d’abord que
inf J(x) = 0.
x∈`2 (R)
Introduisons la suite xn dans `2 (R) définie par xni = δin pour tout i ≥ 1, où δin est
le symbole de Kronecker qui vaut 1 si i = n et 0 sinon. On vérifie aisément que
1
J(xn ) = → 0 quand n → +∞.
n
Comme J est positive, on en déduit que xn est une suite minimisante et que la valeur
du minimum est nulle. Cependant, il est évident qu’il n’existe aucun x ∈ `2 (R) tel
que J(x) = 0. Par conséquent, il n’existe pas de point de minimum pour (2.3).
On voit dans cet exemple que la suite minimisante xn “part à l’infini” et n’est pas
compacte dans `2 (R) (bien qu’elle soit bornée). •
Voici maintenant un exemple modèle qui, malgré son caractère simplifié, est
très représentatif de problèmes réalistes et pratiques de minimisation d’énergies de
transition de phases en science des matériaux.
2.3. EXISTENCE D’UN MINIMUM EN DIMENSION INFINIE 17
h/n
0 k/n 1
Figure 2.1 – Suite minimisante un pour l’Exemple 2.3.2.
Remarque 2.3.1 Le lecteur attentif pourrait faire remarquer que l’espace V dans
l’Exemple 2.3.2 n’est pas un espace de Hilbert. Mais ce n’est pas là que se place
la difficulté et le même résultat est vrai si on remplace V par l’espace de Sobolev
H 1 (0, 1) qui, muni du même produit scalaire, est bien un espace de Hilbert de
dimension infinie, voir par exemple [1].
18 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
Ainsi, sous les conditions (2.2) et (2.6), le problème (2.5) admet une solution.
Malheureusement, la condition (2.6) est inutilisable car invérifiable en général !
On peut cependant la vérifier pour une classe particulière de problèmes, très im-
portants en théorie comme en pratique : les problèmes de minimisation convexe.
Comme nous le verrons dans la Sous-section 2.3.3, si V est un espace de Hilbert,
K un convexe fermé de V , et que J est une fonction convexe et continue sur
K, alors (2.6) a lieu et le problème (2.5) admet une solution. Les motivations pour
introduire ces conditions sont, d’une part, que les hypothèses de convexité sont
souvent naturelles dans beaucoup d’applications, et d’autre part, qu’il s’agit d’une
des rares classes de problèmes pour lesquels la théorie est suffisamment générale et
complète. Mais ceci ne signifie pas que ces conditions sont les seules qui assurent
l’existence d’un minimum ! Néanmoins, en dehors du cadre convexe développé dans
les sous-sections suivantes, des difficultés du type de celles rencontrées dans les
contre-exemples précédents peuvent survenir.
Définition 2.3.2 On dit qu’une fonction J définie sur un ensemble convexe non
vide K ∈ V et à valeurs dans R est convexe sur K si et seulement si
J(θu + (1 − θ)v) ≤ θJ(u) + (1 − θ)J(v) ∀ u, v ∈ K , ∀ θ ∈ [0, 1] . (2.7)
De plus, J est dite strictement convexe si l’inégalité (2.7) est stricte lorsque u 6= v
et θ ∈]0, 1[.
2.3. EXISTENCE D’UN MINIMUM EN DIMENSION INFINIE 19
Remarque 2.3.3 Si J est une application définie sur K à valeurs dans R, on appelle
épigraphe de J l’ensemble Epi(J) = {(λ, v) ∈ R × K , λ ≥ J(v)}. Alors J est
convexe si et seulement si Epi(J) est une partie convexe de R × V . •
Exercice 2.3.1 Soient J1 et J2 deux fonctions convexes sur V, λ > 0, et ϕ une fonction
convexe croissante sur un intervalle de R contenant l’ensemble J1 (V ). Montrer que
J1 + J2 , max(J1 , J2 ), λJ1 et ϕ ◦ J1 sont convexes.
Exercice 2.3.2 Soit (Li )i∈I une famille (éventuellement infinie) de fonctions affines
sur V . Montrer que supi∈I Li est convexe sur V . Réciproquement, soit J une fonction
convexe continue sur V . Montrer que J est égale au supLi ≤J Li où les fonctions Li sont
affines.
Pour les fonctions convexes il n’y a pas de différence entre minima locaux et
globaux comme le montre le résultat élémentaire suivant.
Proposition 2.3.4 Si J est une fonction convexe sur un ensemble convexe K, tout
point de minimum local de J sur K est un minimum global et l’ensemble des points
de minimum est un ensemble convexe (éventuellement vide).
Si de plus J est strictement convexe, alors il existe au plus un point de minimum.
Démonstration. Soit u un minimum local de J sur K. D’après la Définition 2.1.1,
nous pouvons écrire
∃ δ > 0 , ∀ w ∈ K , kw − uk < δ =⇒ J(w) ≥ J(u) . (2.8)
Soit v ∈ K. Pour θ ∈]0, 1[ suffisamment petit, wθ = θv +(1−θ)u vérifie kwθ −uk < δ
et wθ ∈ K puisque K est convexe. Donc, J(wθ ) ≥ J(u) d’après (2.8), et la convexité
de J implique que J(u) ≤ J(wθ ) ≤ θJ(v) + (1 − θ)J(u), ce qui montre bien que
J(u) ≤ J(v), c’est-à-dire que u est un minimum global sur K.
D’autre part, si u1 et u2 sont deux minima et si θ ∈ [0, 1], alors w = θu1 + (1 −
θ)u2 est un minimum puisque w ∈ K et que
inf J(v) ≤ J(w) ≤ θJ(u1 ) + (1 − θ)J(u2 ) = inf J(v) .
v∈K v∈K
Le même raisonnement avec θ ∈]0, 1[ montre que, si J est strictement convexe, alors
nécessairement u1 = u2 . 2
Une propriété agréable des fonctions convexes “propres” (c’est-à-dire qui ne
prennent pas la valeur +∞) est qu’elles sont continues.
Lemme 2.3.5 Soit v0 ∈ V et J une fonction convexe majorée sur la boule unité de
centre v0 . Alors J est continue sur cette boule ouverte.
Démonstration. Sans perte de généralité, par translation et addition on peut se
ramener au cas où v0 = 0 et J(0) = 0. Soit v 6= 0, kvk < 1 et M la majoration de J
sur la boule unité. Par convexité de J pour θ = kvk, on a
v v
J(v) = J θ + (1 − θ)0 ≤ θJ( ) + (1 − θ)J(0) ≤ M kvk
kvk kvk
20 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
Définition 2.3.6 On dit qu’une fonction J, définie sur un ensemble convexe K, est
fortement convexe (ou α-convexe) s’il existe α > 0 tel que, pour tout u, v ∈ K et
tout θ ∈ [0, 1],
αθ(1 − θ)
J(θu + (1 − θ)v) ≤ θJ(u) + (1 − θ)J(v) − ku − vk2 . (2.9)
2
Evidemment, les fonctions fortement convexes sont également strictement convexes.
Dans la Définition 2.3.6, il est possible de ne tester la forte convexité de J que pour
la valeur θ = 1/2, comme le montre l’exercice suivant.
Exercice 2.3.3 Soit J une fonction continue sur un ensemble convexe K et α > 0
tels que, pour tout u, v ∈ K,
u+v J(u) + J(v) α
J ≤ − ku − vk2 . (2.10)
2 2 8
Montrer que J est fortement ou α-convexe.
Exercice 2.3.4 Montrer qu’une fonction J, définie sur un ensemble convexe K, est
fortement ou α-convexe si et seulement si la fonction J(v) − α2 kvk2 est convexe sur K.
Exercice 2.3.6 Soit V un espace de Hilbert (avec les notations usuelles pour son
produit scalaire et sa norme). Soit L(v) une forme linéaire continue sur V , et soit
a(u, v) une forme bilinéaire continue et symétrique sur V × V . Soit la fonction J définie
sur V par
1
J(v) = a(v, v) − L(v).
2
Montrer que J est convexe sur V si a est positive, c’est-à-dire a(v, v) ≥ 0 pour tout
v ∈ V . Montrer que J est fortement convexe sur V si a est coercive, c’est-à-dire s’il
existe C > 0 tel que a(v, v) ≥ Ckvk2 pour tout v ∈ V .
2.3. EXISTENCE D’UN MINIMUM EN DIMENSION INFINIE 21
Le résultat suivant sera essentiel dans l’obtention d’un résultat d’existence d’un
minimum en dimension infinie. En particulier, il permet de conclure qu’une fonction
J fortement convexe et continue sur un ensemble K convexe fermé non vide est
“infinie à l’infini” dans K, c’est-à-dire vérifie la propriété (2.2).
Si de plus J est fortement convexe sur K, alors il existe deux constantes γ > 0 et
η ∈ R telles que
J(v) ≥ γkvk2 + η ∀ v ∈ K . (2.12)
On en déduit
α α
J(v) ≥ kvk2 − hv, v0 i + L(v) + C1 ,
4 2
avec C1 = (α/4)kv0 k2 + L(v0 ) − J(v0 ) + 2δ. D’après l’inégalité de Cauchy-Schwarz
appliqué à hv, v0 i et la continuité de L, i.e. |L(v)| ≤ kLkV 0 kvk (voir la Définition
8.1.10), il vient
α 2 αkv0 k α
J(v) ≥ kvk − kLkV 0 + kvk + C1 ≥ kvk2 + η ,
4 2 8
ce qui montre que la suite (un ) est de Cauchy, et donc converge vers une limite u,
qui est nécessairement un minimum de J sur K puisque J est continue et K fermé.
L’unicité du point de minimum a été montrée dans la Proposition 2.3.4. Enfin, si
v ∈ K, (u + v)/2 ∈ K car K est convexe, d’où, toujours grâce à (2.10),
α 2 J(u) J(v) u+v J(v) − J(u)
ku − vk ≤ + −J ≤ ,
8 2 2 2 2
u+v
car J ≥ J(u). 2
2
Il est possible de généraliser en grande partie le Théorème 2.3.8 au cas de
fonctions J qui sont seulement convexes (et non pas fortement convexes). Cependant,
autant la démonstration du Théorème 2.3.8 est élémentaire, autant celle du théorème
suivant est délicate. Elle repose en particulier sur la notion de convergence faible
que l’on peut considérer comme “hors-programme” dans le cadre de ce cours.
Remarque 2.3.11 Le Théorème 2.3.9 reste vrai si l’on suppose simplement que V
est un espace de Banach réflexif (i.e. que le dual de V 0 est V ). •
Nous indiquons brièvement comment on peut démontrer le Théorème 2.3.9 dans
le cas d’un espace de Hilbert séparable (c’est-à-dire qui admet une base hilbertienne
dénombrable). On définit la notion de convergence faible dans V (pour plus de
détails à ce sujet, voir les chapitres 3 et 5 de [9]).
Définition 2.3.12 On dit qu’une suite (un )n∈N de V converge faiblement vers u ∈ V
si
∀ v ∈ V , lim hun , vi = hu, vi .
n→+∞
Soit (ei )i≥1 une base hilbertienne de V . Si on note uni = hun , ei i les composantes
dans cette base d’une suite un , uniformément bornée dans V , il est facile de vérifier
que la Définition 2.3.12 de la convergence faible est équivalente à la convergence
de toutes les suites de composantes (uni )n≥1 pour i ≥ 1.
Comme son nom l’indique la convergence faible est une notion “plus faible”
que la convergence usuelle dans V . En effet, par simple application de l’inégalité de
Cauchy-Schwarz, on voit que si un converge vers u, c’est-à-dire si limn→+∞ kun −
uk = 0, alors un converge faiblement vers u. Réciproquement, en dimension infinie il
existe des suites qui convergent faiblement mais pas au sens usuel (que l’on appelle
parfois “convergence forte” par opposition). Par exemple, la suite un = en converge
faiblement vers zéro, mais pas fortement puisqu’elle est de norme constante égale à
1. L’intérêt de la convergence faible vient du résultat suivant.
Lemme 2.3.13 De toute suite (un )n∈N bornée dans V on peut extraire une sous-
suite qui converge faiblement.
Démonstration. Comme la suite un est bornée, chaque suite d’une composante
uni est bornée dans R. Pour chaque i, il existe donc une sous-suite, notée uni i , qui
converge vers une limite ui . Par un procédé d’extraction diagonale de suites, on
0
obtient alors une sous-suite commune n0 telle que, pour tout i, uni converge vers ui .
0
Ce qui prouve que un converge faiblement vers u (on vérifie que u ∈ V ). 2
Lemme 2.3.14 Une partie convexe fermée K de V est l’intersection des demi-
espaces fermés qui contiennent K.
24 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
Démonstration. Il est clair que K est inclus dans l’intersection des demi-espaces
fermés qui le contiennent. Réciproquement, supposons qu’il existe un point u0 de
cette intersection qui n’appartiennent pas à K. On peut alors appliquer le Théorème
8.1.12 de séparation d’un point et d’un convexe et construire ainsi un demi-espace
fermé qui contient K mais pas u0 . Ceci est une contradiction avec la définition de
u0 , donc u0 ∈ K. 2
Lemme 2.3.15 Soit K un ensemble convexe fermé non vide de V . Alors K est
fermé pour la convergence faible.
De plus, si J est convexe et semi-continue inférieurement sur K (voir l’Exercice
2.2.3 pour cette notion), alors J est aussi semi-continue inférieurement sur K pour
la convergence faible.
Démonstration. Par définition, si un converge faiblement vers u, alors L(un )
converge vers L(u). Par conséquent, un demi-espace fermé de V est fermé pour la
convergence faible. Le Lemme 2.3.14 permet d’obtenir la même conclusion pour K.
D’après les hypothèses sur J, l’ensemble Epi(J) (défini à la Remarque 2.3.3)
est un convexe fermé de R × V , donc il est aussi fermé pour la convergence faible.
On en déduit alors facilement le résultat : si la suite (v n ) tend faiblement vers v
dans K, alors lim inf n→+∞ J(v n ) ≥ J(v). 2
2.4 Différentiabilité
Jusqu’ici nous ne nous sommes intéressés qu’aux questions d’existence de mini-
mum aux problèmes d’optimisation. Mais il importe aussi de caractériser les points
de minimum, autant d’un point de vue théorique que pratique, afin de les calculer.
Pour ce faire, on utilise des conditions d’optimalité, c’est-à-dire des conditions
nécessaires et parfois suffisantes de minimalité. Ces conditions d’optimalité s’écrivent
avec les dérivées de la fonction objectif et des éventuelles contraintes. Nous allons
donc rappeler comment on calcule ces dérivées ou différentielles de fonctions définies
sur des espaces de Hilbert.
Avant d’en venir à ces considérations techniques, nous motivons l’étude de ces
conditions d’optimalité en rappelant le cas, très simple, du calcul des minima d’une
fonction dérivable J(x) définie sur l’intervalle [a, b] ⊂ R à valeurs réelles. Il est bien
connu que si x0 est un point de minimum local de J sur l’intervalle [a, b], alors on a
Remarque 2.4.2 La Définition 2.4.1 est en fait valable si V est seulement un espace
de Banach (on n’utilise pas de produit scalaire dans (2.15)). Cependant, si V est un
espace de Hilbert, on peut préciser la relation (2.15) en identifiant V et son dual V 0
grâce au Théorème de représentation de Riesz 8.1.11. En effet, il existe un unique
p ∈ V tel que hp, wi = L(w), donc (2.15) devient
|o(w)|
J(u + w) = J(u) + hp, wi + o(w) avec lim =0. (2.16)
w→0 kwk
On note aussi parfois p = J 0 (u), ce qui peut prêter à confusion... La formule (2.16)
est souvent plus “naturelle” que (2.15), notamment si V = Rn ou V = L2 (Ω). •
Dans la plupart des applications, il suffit souvent de déterminer la forme linéaire
continue L = J 0 (u) ∈ V 0 car on n’a pas besoin de l’expression explicite de p = J 0 (u) ∈
V lorsque V 0 est identifié à V . En pratique, il est plus facile de trouver l’expression
explicite de L que celle de p, comme le montrent les exercices suivants.
Exercice 2.4.2 Soit a une forme bilinéaire symétrique continue sur V × V . Soit L une
forme linéaire continue sur V . On pose J(u) = 21 a(u, u) − L(u). Montrer que J est
dérivable sur V et que hJ 0 (u), wi = a(u, w) − L(w) pour tout u, w ∈ V .
Remarque 2.4.3 Il existe d’autres notions de différentiabilité, plus faible que celle
au sens de Fréchet. Par exemple, on rencontre souvent la définition suivante. On
dit que la fonction J, définie sur un voisinage de u ∈ V à valeurs dans R, est
différentiable au sens de Gâteaux en u s’il existe L ∈ V 0 tel que
J(u + δw) − J(u)
∀w ∈ V , lim = L(w) . (2.18)
δ→0, δ>0 δ
On parle aussi de différentiabilité directionnelle et w est la direction de dérivation
dans (2.18). L’intérêt de cette notion est que la vérification de (2.18) est plus aisée
que celle de (2.15). Cependant, si une fonction dérivable au sens de Fréchet l’est
aussi au sens de Gâteaux, la réciproque est fausse, même en dimension finie, comme
le montre l’exemple suivant dans R2
x6
J(x, y) = pour (x, y) 6= (0, 0) , J(0, 0) = 0 .
(y − x2 )2 + x8
Convenons que, dans ce qui suit, nous dirons qu’une fonction est dérivable lorsqu’elle
l’est au sens de Fréchet, sauf mention explicite du contraire. •
Examinons maintenant les propriétés de base des fonctions convexes dérivables.
Remarque 2.4.6 Les conditions (2.20) et (2.23) ont une interprétation géométrique
simple : la première signifie que la fonction convexe J(v) est toujours au dessus de
son plan tangent en u (considéré comme une fonction affine de v), tandis que la
deuxième affirme que J(v) est au dessus d’une fonction quadratique en v qui lui est
tangente en u. Les conditions (2.21) et (2.24) sont des propriétés de monotonie ou
de croissance de J 0 . •
dans laquelle on fait tendre θ vers 0 pour obtenir (2.23). Pour obtenir (2.24) il
suffit d’additionner (2.23) avec lui-même en échangeant u et v. Montrons que (2.24)
implique (2.22). Pour u, v ∈ V et t ∈ R, on pose ϕ(t) = J(u + t(v − u)). Alors ϕ est
dérivable et donc continue sur R, et ϕ0 (t) = hJ 0 (u + t(v − u)), v − ui, de sorte que,
d’après (2.24)
ϕ0 (t) − ϕ0 (s) ≥ α(t − s)kv − uk2 si t ≥ s . (2.25)
Soit θ ∈]0, 1[. En intégrant l’inégalité (2.25) de t = θ à t = 1 et de s = 0 à s = θ, on
obtient
αθ(1 − θ)
θϕ(1) + (1 − θ)ϕ(0) − ϕ(θ) ≥ kv − uk2 ,
2
c’est-à-dire (2.22). 2
Exercice 2.4.5 Montrer qu’une fonction J dérivable sur V est strictement convexe si
et seulement si
ou encore
hJ 0 (u) − J 0 (v), u − vi > 0 ∀ u, v ∈ V avec u 6= v .
ko(w)kW
f (u + w) = f (u) + L(w) + o(w) avec lim =0. (2.26)
w→0 kwkV
Définition 2.4.7 Soit J une fonction de V dans R. On dit que J est deux fois
dérivable en u ∈ V si J est dérivable dans un voisinage de u et si sa dérivée J 0 (u)
est dérivable en u. On note J 00 (u) la dérivée seconde de J en u qui vérifie
ko(w)kV 0
J 0 (u + w) = J 0 (u) + J 00 (u)w + o(w) avec lim = 0.
w→0 kwkV
Telle qu’elle est définie la dérivée seconde est difficile à évaluer en pratique car
w → J 00 (u)w est un opérateur linéaire continu de V dans V 0 . Heureusement, en
faisant agir J 00 (u)w sur v ∈ V on obtient une brave forme bilinéaire continue sur
V ×V que l’on notera J 00 (u)(w, v) en lieu et place de (J 00 (u)w) v. Au vu de la formule
de Taylor d’ordre deux ci-dessous, en pratique on calcule plutôt J 00 (u)(w, w).
Lemme 2.4.8 Si J est une fonction deux fois dérivable de V dans R, elle vérifie
1
J(u + w) = J(u) + hJ 0 (u), wi + J 00 (u)(w, w) + o(kwk2 ), (2.27)
2
o(kwk2 )
avec limw→0 kwk2
= 0 où J 00 (u) est identifiée à une forme bilinéaire continue sur
V ×V.
Démonstration. On écrit un développement de Taylor avec reste exact
Z 1
J(u + w) = J(u) + hJ 0 (u + sw), wids
0
Nous laissons au lecteur le soin de finir les détails pour arriver à (2.27). 2
Exercice 2.4.6 Soit a une forme bilinéaire symétrique continue sur V × V . Soit L une
forme linéaire continue sur V . On pose J(u) = 21 a(u, u) − L(u). Montrer que J est deux
fois dérivable sur V et que J 00 (u)(v, w) = a(v, w) pour tout u, v, w ∈ V . Appliquer ce
résultat aux exemples des Exercices 2.4.3, 2.4.4.
Lorsque J est deux fois dérivable on retrouve la condition usuelle de convexité :
si la dérivée seconde est positive, alors la fonction est convexe.
Exercice 2.4.7 Montrer que si J est deux fois dérivable sur V les conditions des
Propositions 2.4.4 et 2.4.5 sont respectivement équivalentes à
Lemme 2.4.10 Soit J une fonction différentiable de V dans R telle que sa dérivée
J 0 est Lipschitzienne de constante L > 0. Alors
L
J(v) ≤ J(u) + hJ 0 (u), v − ui + kv − uk2 pour tout v, u ∈ V.
2
Remarque 2.4.11 Une interprétation du Lemme 2.4.10 est qu’une fonction diffé-
rentiable à dérivée L-Lipschitzienne peut être majorée par une fonction quadratique
qui lui est tangente au point u. C’est en quelque sorte une version opposée (ma-
joration au lieu de minoration) de la Remarque 2.4.6 pour une fonction α-convexe
différentiable (notons néanmoins qu’il n’y a pas d’hypothèse de convexité dans le
Lemme 2.4.10). Par conséquent, les fonctions fortement convexes à dérivée Lipschit-
zienne sont remarquables car elles sont encadrées par deux fonctions quadratiques
qui sont tangentes au point u. Cette propriété très utile sera exploitée dans les
preuves de convergence des algorithmes d’optimisation. •
Démonstration. Ecrivons une formule de Taylor avec reste exact
Z 1
J(v) = J(u) + hJ 0 (u + t(v − u)), v − uidt.
0
Par conséquent,
Z 1
0
J(v) − J(u) − hJ (u), v − ui ≤ kJ 0 (u + t(v − u)) − J 0 (u)kkv − ukdt,
0
Exercice 2.4.8 Soit J une fonction de classe C 2 de V dans R. On suppose qu’il existe
une constante L > 0 telle que J 00 (u)(w, w) ≤ Lkwk2 pour tout u, w ∈ V . Montrer que
J 0 est Lipschitzienne de constante L > 0.
kx − xK k = min kx − yk.
y∈K
Exercice 2.5.2 On reprend l’Exemple 1.2.4 du problème “aux moindres carrés”. Mon-
trer que ce problème admet toujours une solution et écrire l’équation d’Euler correspon-
dante. Montrer aussi que la solution est unique si et seulement si KerA = {0}.
2.5. CONDITIONS D’OPTIMALITÉ 31
Ax − b = B ∗ p avec p ∈ Rm .
1 L
Z Z L
0 2
inf J(u) = µ|u (x)| dx − f (x) u(x) dx,
u∈K 2 0 0
où f (x) est une fonction continue sur [0, L], µ > 0 et K est défini par
Vérifier que K est convexe et calculer la dérivée directionnelle J 0 (u)(w). Montrer que le
point de minimum u∗ ∈ K de J (s’il existe et s’il est de classe C 2 ) vérifie la condition
nécessaire suivante
−µu00∗ (x) = f (x)
si 0 < x < L,
u∗ (0) = u∗ (L) = 0.
En déduire une formule pour ce point de minimum et donc qu’il est unique.
Exercice 2.5.5 Soit K un convexe fermé non vide de V , soit a une forme bilinéaire
symétrique continue coercive sur V , et soit L une forme linéaire continue sur V . Montrer
que J(v) = 21 a(v, v)−L(v) admet un unique point de minimum dans K, noté u. Montrer
que u est aussi l’unique solution du problème (appelé inéquation variationnelle)
Exercice 2.5.6 Soit J1 et J2 deux fonctions convexes continues sur une partie convexe
fermée non vide K ⊂ V . On suppose que J1 seulement est dérivable. Montrer que u ∈ K
est un minimum de J1 + J2 si et seulement si
Les remarques suivantes, qui sont des applications simples du Théorème 2.5.1,
vont nous donner l’intuition de la notion de “multiplicateur de Lagrange” qui sera
développée à la Sous-section suivante.
P = {v ∈ V , hai , vi = 0 pour 1 ≤ i ≤ M } ,
32 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
avec des réels λi (qui seront appelés multiplicateurs de Lagrange dans le Théorème
2.5.6).
Exercice 2.5.8 Soit (a1 , . . . , aM ) une famille libre dans V . Supposons que K est défini
par
K = {v ∈ V , hai , vi ≤ 0 pour 1 ≤ i ≤ M } .
Vérifier que K est un cône convexe fermé, ce qui signifie que K est un ensemble convexe
fermé tel que λv ∈ K pour tout v ∈ K et tout λ ≥ 0. Montrer que la condition
d’optimalité (2.29) implique que
En déduire que si u vérifie (2.29), alors il existe des réels positifs ou nuls λi ≥ 0 tels que
M
X
u∈K et ∃ λ1 , . . . , λM ≥ 0 , J 0 (u) + λ i ai = 0 , (2.33)
i=1
et, de plus, λi = 0 si hai , ui < 0. Nous verrons que ces réels λi sont encore appelés
multiplicateurs de Lagrange au Théorème 2.5.17.
Exercice 2.5.9 Montrer que K(v) = V si v est intérieur à K. Montrer que K(v0 ) =
{0} si K = K0 ∪ {v0 } avec la distance d(v0 , K0 ) > 0.
L’intérêt du cône des directions admissibles réside dans le résultat suivant, qui
donne une condition nécessaire d’optimalité.
hJ 0 (u), wi ≥ 0 ∀ w ∈ K(u) .
Démonstration. Soit une direction admissible w ∈ K(u). Il existe donc une suite
vn − u
lim v n = u , lim εn = 0 , lim =w
n→+∞ n→+∞ n→+∞ εn
Pour n suffisamment grand,
Contraintes d’égalité
Dans ce premier cas on suppose que K est donné par
K = {v ∈ V , F (v) = 0} , (2.36)
où F (v) = (F1 (v), ..., FM (v)) est une application de V dans RM , avec M ≥ 1. La
condition nécessaire d’optimalité prend alors la forme suivante.
Théorème 2.5.6 Soit u ∈ K où K est donné par (2.36). On suppose que J est dé-
rivable en u ∈ K et que les fonctions (Fi )1≤i≤M sont continûment dérivables dans un
voisinage de u. On suppose de plus que les vecteurs Fi0 (u) 1≤i≤M sont linéairement
Figure 2.3 – Variété K (en bleu) et son hyperplan tangent K(u) en un point u (en
rose). Ce dessin correspond à V = R3 et M = 1 contrainte.
Démonstration. Montrons d’abord que le cône des directions admissibles K(u) est
précisément l’espace tangent à la variété K, définie par (2.36), au point u (voir la
Figure 2.3), c’est-à-dire
K(u) = KerF 0 (u), (2.38)
avec, par définition
M n o
⊥
\
0
KerF (u) = [Fi0 (u)] = w ∈ V , hFi0 (u), wi = 0 pour i = 1, . . . , M .
i=1
0 = hF 0 (u), wi + o(1),
c’est-à-dire, à la limite n → +∞, que w appartient à KerF 0 (u). Pour démontrer l’in-
clusion inverse, on utilise le Lemme 2.5.7 ci-dessous avec, pour tout w ∈ KerF 0 (u),
F(ε, v) = F (v + εw). Pour v0 = u, on vérifie sans peine les hypothèses du Lemme
2.5.7, en particulier car la différentielle partielle Dv F(0, u) = F 0 (u) est égale à la
matrice de lignes Fi0 (u) 1≤i≤M qui est surjective de V dans RM puisque ses lignes
sont libres. Pour tout ε > 0, on choisit v = u + εw qui vérifie F (v) = o(ε) puisque
hF 0 (u), wi = 0. On en déduit donc que la suite v ε fournie par le Lemme 2.5.7
est telle que kv − vε kV = o(ε), c’est-à-dire qu’elle est admissible pour la direction
w dans la définition de K(u). Par conséquent, w ∈ K(u), ce qui démontre que
K(u) = KerF 0 (u).
Comme nous venons de montrer que K(u) est un espace vectoriel, on peut
prendre successivement w et −w dans la Proposition 2.5.5, ce qui conduit à
M
⊥
\
0 0
hJ (u), wi = 0 ∀ w ∈ KerF (u) = [Fi0 (u)] ,
i=1
c’est-à-dire que J 0 (u) est engendré par les Fi0 (u) 1≤i≤M (notons que les multiplica-
teurs de Lagrange sont définis de manière unique). Une autre démonstration (plus
géométrique) est proposé dans la preuve de la Proposition 2.5.13. 2
Lemme 2.5.7 Soit F(ε, v) une fonction de R×V dans RM . On suppose qu’il existe
v0 ∈ V tel que F(0, v0 ) = 0 et que F(ε, v) est continûment différentiable dans
un voisinage de (0, v0 ) (c’est-à-dire différentiable de différentielle continue). Si la
différentielle (partielle) V 3 h → hDv F(0, v0 ), hi ∈ RM est surjective, alors il existe
un autre voisinage V de (0, v0 ) et une constante C > 0 tels que pour tout (ε, v) ∈ V
il existe vε ∈ V qui vérifie
Une démonstration du Lemme 2.5.7 peut être trouvée, dans un cadre plus gé-
néral (théorème de l’application surjective), dans [4].
Remarque 2.5.8 Lorsque les vecteurs Fi0 (u) 1≤i≤M sont linéairement indépen-
dants (ou libres), on dit que l’on est dans un cas régulier. Dans le cas contraire,
on parle de cas non régulier et la conclusion du Théorème 2.5.6 est fausse comme
le montre l’exemple suivant.
Prenons V = R, M = 1, F (v) = v 2 , J(v) = v, d’où K = {0}, u = 0, F 0 (u) = 0 :
il s’agit donc d’un cas non régulier. Comme J 0 (u) = 1, (2.37) n’a pas lieu. •
Pour bien comprendre la portée du Théorème 2.5.6, nous l’appliquons sur l’Exemple
1.2.3
1
min J(x) = Ax · x − b · x ,
x∈ KerB 2
36 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
Exercice 2.5.10 Généraliser les résultats ci-dessus pour cette variante de l’Exemple
1.2.3
1
min J(x) = Ax · x − b · x ,
Bx=c 2
où c ∈ Rm est un vecteur donné.
Exercice 2.5.11 Soit A une matrice carrée d’ordre n, symétrique. On veut caractériser
et calculer les solutions de
inf
n
J(x) = Ax · x,
x∈R ,kxk=1
où kxk est la norme euclidienne de x. Appliquer le Théorème 2.5.6 et en déduire que les
points de minimum de J sur la sphère unité sont des vecteurs propres de A associés à
la plus petite valeur propre.
Exercice 2.5.12 Selon la légende rapportée par Virgile dans l’Énéide, la fondation de
la ville de Carthage par la reine Didon conduisit à un problème typique du calcul des
variations. Il s’agissait de trouver la plus grande surface possible (pour la ville) s’appuyant
sur le rivage (une droite) et de frontière terrestre de longueur donnée. La réponse est
intuitivement un demi disque. Mathématiquement, le problème est de trouver la courbe
plane, définie par le graphe d’une fonction y(x) pour x ∈ (0, ξ), de longueur fixée l ≥ 0
qui enclôt avec le segment (0, ξ) reliant ses deux extrémités l’aire maximum. Autrement
dit, on résout Z ξ
sup y(x) dx,
ξ≥0,y∈C 1 [0,ξ] 0
2.5. CONDITIONS D’OPTIMALITÉ 37
sup b · x et sup b · x
Ax·x≤1 Ax·x=1
sont équivalents et qu’ils ont une solution. Utiliser le Théorème 2.5.6 pour cal-
culer cette solution et montrer qu’elle est unique.
2. On introduit un ordre partiel dans l’ensemble des matrices symétriques définies
positives d’ordre n en disant que A ≥ B si et seulement si Ax · x ≥ Bx · x
pour tout x ∈ Rn . Déduire de la question précédente que, si A ≥ B, alors
B −1 ≥ A−1 .
Exercice 2.5.15 En théorie cinétique des gaz les molécules de gaz sont représentées
en tout point de l’espace par une fonction de répartition f (v) dépendant de la vitesse
microscopique v ∈ RN . Les quantités macroscopiques, comme la densité du gaz ρ, sa
vitesse u, et sa température T , se retrouvent grâce aux moments de la fonction f (v)
Z Z Z
1 2 N 1
ρ= f (v) dv , ρu = v f (v) dv , ρu + ρT = |v|2 f (v) dv .
RN R N 2 2 2 RN
(2.39)
Boltzmann a introduit l’entropie cinétique H(f ) définie par
Z
H(f ) = f (v) log f (v) dv .
RN
Montrer que H est strictement convexe sur l’espace des fonctions f (v) > 0 mesurables
telle que H(f ) < +∞. On minimise H sur cet espace sous les contraintes de moment
38 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
(2.39), et on admettra qu’il existe un unique point de minimum M (v). Montrer que ce
point de minimum est une Maxwellienne définie par
|v − u|2
ρ
M (v) = exp − .
(2πT )N/2 2T
Remarque 2.5.11 On pourrait croire que le Lemme 2.5.10 est seulement une as-
tuce ou un jeu d’écriture pour faire disparaitre la contrainte de manière artificielle. Il
n’en est rien et la notion de Lagrangien est extrêmement utile pour les algorithmes
numériques de minimisation. Donnons en un rapide aperçu. Nous verrons que tous
les algorithmes d’optimisation sont itératifs, c’est-à-dire que l’on construit une suite
de solutions approchées un qui convergerait vers une solution u. A chaque itération,
on définit la nouvelle itérée un+1 à partir de la précédente un . En général, il est très
2.5. CONDITIONS D’OPTIMALITÉ 39
difficile, voire impossible, de garantir que les solutions approchées un vérifient exac-
tement les contraintes F (un ) = 0. L’utilisation du Lagrangien permet de contourner
cette difficulté. Imaginons, par exemple dans le cas d’une seule contrainte (M = 1),
que F (un ) > 0. Alors, si on choisit un multiplicateur de Lagrange µn > 0, la minimi-
sation (totale ou partielle) sans contrainte de L(v, µn ) permet d’obtenir une nouvelle
itérée un+1 qui minimise au mieux la fonction objectif et la violation de la contrainte
(la positivité de µn conduit à minimiser F (v), comme J(v)). Inversement, si on avait
F (un ) < 0, on aurait choisi un multiplicateur de Lagrange µn < 0 et la minimisation
sans contrainte de L(v, µn ) conduirait à maximiser F (v) (tout en minimisant tou-
jours J(v)), c’est-à-dire à réduire la violation de la contrainte. Tout ceci sera précisé
dans la Sous-section 3.4.2 lorsque sera présenté l’algorithme d’Uzawa. •
Le résultat suivant permet de donner une interprétation concrète à la notion
de multiplicateur de Lagrange. Plus précisément, les multiplicateurs de Lagrange
λi dans le Théorème 2.5.6 donnent la sensibilité de la valeur minimale de J aux
variations du niveau des contraintes Fi (au signe près). Par exemple, en économie, les
multiplicateurs de Lagrange mesurent des prix ou des coûts marginaux, alors qu’en
mécanique ils fournissent des forces de réaction correspondant à des déplacements
imposés.
Pour ∈ RM , on considère le problème d’optimisation
inf J(v), (2.41)
v∈V, F (v)=
où le niveau des contraintes est paramétré par . On note Jmin () la valeur minimum
dans (2.41).
0
J et F1 , ..., FM sont deux fois continûment dérivables
pose que les fonctions
M
et que
les vecteurs Fi (u) 1≤i≤M sont linéairement indépendants. Soit λ ∈ R le multipli-
cateur de Lagrange défini par le Théorème 2.5.6. Alors tout minimum local u de J
sur K vérifie
M
! M
⊥
X \
00 00
J (u) + λi Fi (u) (w, w) ≥ 0 ∀ w ∈ K(u) = [Fi0 (u)] . (2.42)
i=1 i=1
En dérivant on obtient
et
fi00 (t) = Fi00 u(t) u0 (t), u0 (t) + hFi0 u(t) , u00 (t)i pour 1 ≤ i ≤ M.
Grâce à (2.43) on peut éliminer les dérivées premières et u00 (0) pour obtenir (en
sommant les deux dernières équations)
M
!
X
λi Fi00 (u) + J 00 (u) u0 (0), u0 (0) ≥ 0,
i=1
qui n’est rien d’autre que (2.42) lorsque u0 (0) parcourt K(u).
L’existence de tels chemins admissibles u(t) et le fait que l’ensemble des u0 (0)
décrit la totalité du cône des directions admissibles K(u) est une conséquence du
théorème des fonctions implicites que l’on peut appliquer grâce à l’hypothèse que la
famille Fi0 (u) 1≤i≤M est libre (voir par exemple [4]). 2
Contraintes d’inégalité
Dans ce deuxième cas on suppose que K est donné par
Théorème 2.5.17 On suppose que K est donné par (2.44), que les fonctions J et
F1 , . . . , FM sont dérivables en u et que les contraintes sont qualifiées en u. Alors, si
u est un minimum local de J sur K, il existe λ1 , . . . , λM ≥ 0, appelés multiplicateurs
de Lagrange, tels que
M
X
0
J (u)+ λi Fi0 (u) = 0, λi ≥ 0, λi = 0 si Fi (u) < 0, ∀ i ∈ {1, . . . , M }. (2.46)
i=1
(On peut montrer que K̃(u) n’est autre que le cône K(u) des directions admissibles,
voir [4]). Soit w une direction admissible satisfaisant (2.45), w ∈ K̃(u), et un réel
δ > 0. Nous allons montrer que u + ε(w + δw) ∈ K pour tout réel ε > 0 assez petit.
Il faut examiner trois cas de figure.
1. Si i ∈
/ I(u), on a Fi (u) < 0 et Fi (u + ε(w + δw)) < 0 par continuité de Fi si
ε est assez petit.
2. Si i ∈ I(u) et hFi0 (u), wi < 0, alors
et
n M
X o
K̂ = q ∈ V , ∃ λ1 , . . . , λM ≥ 0 , q = − λi ai .
i=1
hp, wi ≥ 0 ∀ w ∈ K =⇒ p ∈ K̂ .
une suite d’éléments de K̂ (donc avec λni ≥ 0 ∀ i ∀ n), convergeant vers une limite
q ∈ V . Alors il est clair que chaque suite (λni ) converge dans R+ vers une limite
λi ≥ 0 (pour 1 ≤ i ≤ M ) puisque les vecteurs (ai )1≤i≤M forment une base de
l’espace qu’ils engendrent. On a donc q = − M
P
i=1 λi ai ∈ K̂, qui est donc fermé.
Si les vecteurs (ai )1≤i≤M sont linéairement dépendants, nous procédons par ré-
currence sur leur nombre M . La propriété est évidente lorsque M = 1, et nous sup-
posons qu’elle est vraie lorsque le nombre de vecteurs ai est inférieur PM à M . Comme
les vecteurs (ai )1≤i≤M sont liés, il existe une relation de la forme i=1 µi ai = 0, avec
au moins un des coefficients µi qui est strictement positif. Soit alors q = − M
P
i=1 λi ai
XM
un élément de K̂. Pour tout t ≤ 0, on peut aussi écrire q = − (λi + tµi )ai , et on
i=1
peut choisir t ≤ 0 pour que
Remarque 2.5.22 Comme pour le cas des contraintes d’égalité, la notion de La-
grangien est extrêmement utile pour les algorithmes numériques de minimisation.
Si la contrainte F (un ) ≤ 0 est violée pour une solution approchée un , alors, pour
un multiplicateur de Lagrange µn ≥ 0, dont les composantes correspondant à une
contrainte violée sont strictement positives tandis que les composantes correspon-
dant à une contrainte inactive sont nulles, la minimisation sans contrainte de L(v, µn )
permet d’obtenir une nouvelle itérée un+1 qui minimise au mieux la fonction objectif
et la violation de la contrainte. Cela sera précisé dans la Sous-section 3.4.2 lorsque
sera présenté l’algorithme d’Uzawa. •
La condition nécessaire d’optimalité (2.46) du Théorème 2.5.17 est bien la
stationnarité du Lagrangien de la Définition 2.5.20 puisque
∂L
(u, λ) = J 0 (u) + λ · F 0 (u) = 0,
∂v
2.5. CONDITIONS D’OPTIMALITÉ 45
Ax − b + B ∗ p = 0 , p≥0, p · (Bx − c) = 0.
Exercice 2.5.18 Soit f une fonction définie sur un ouvert borné Ω. Pour > 0 on
considère le problème de régularisation suivant
Z
Z inf |∇u|2 dx,
Ω
u ∈ V0 , |u − f |2 dx ≤ 2
Ω
Nous pouvons alors énoncer les conditions nécessaires d’optimalité sur l’en-
semble (2.53).
46 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
des coefficients αi découle de l’inversibilité de la matrice (hFi0 (u), Fj0 (u)i)ij . Il est
clair cependant que (2.45) n’implique pas (2.56). •
On vérifie facilement que (2.57) entraîne (2.54), c’est-à-dire que les contraintes sont
qualifiées. •
2.6. POINT-SELLE, THÉORÈME DE KUHN ET TUCKER, DUALITÉ 47
2.6.1 Point-selle
De manière abstraite, V et Q étant deux espaces de Hilbert réels, un Lagrangien
L est une application de V × Q (ou d’une partie U × P de V × Q) dans R. Dans
le cadre du Théorème 2.5.6 sur les contraintes d’égalité (ou plutôt de la Définition
2.5.9), nous avions U = V , P = Q = RM et L(v, q) = J(v)+q ·F (v). La situation est
un peu différente dans le cadre du Théorème 2.5.17 sur les contraintes d’inégalité, où
pour le même Lagrangien L(v, q) = J(v) + q · F (v) il faut prendre U = V , Q = RM
et P = (R+ )M .
48 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
Le Théorème 2.5.17 a donné une condition nécessaire d’optimalité. Dans cette sous-
section nous allons voir que cette condition est aussi suffisante si les contraintes et
la fonction coût sont convexes. En effet, la Proposition 2.6.2 affirme que, si (u, p)
est un point-selle du Lagrangien, alors u réalise le minimum de J sur K. Pour un
problème de minimisation convexe avec des contraintes d’inégalités convexes, nous
allons établir une réciproque de ce résultat, c’est-à-dire que, si u réalise le minimum
de J sur K, alors il existe p ∈ (R+ )M tel que (u, p) soit point-selle du Lagrangien.
On suppose désormais que J, F1 , . . . , FM sont convexes continues sur V .
(R+ )M tel que (u, p) soit un point-selle du Lagrangien L sur V × (R+ )M ou, de
manière équivalente, tel que
M
X
0
F (u) ≤ 0 , p ≥ 0 , p · F (u) = 0 , J (u) + pi Fi0 (u) = 0 . (2.63)
i=1
Comme dans le cas des contraintes d’égalité (voir le Lemme 2.5.12), on peut
donner une interprétation concrète des multiplicateurs de Lagrange λi dans le Théo-
rème 2.5.17 : ils donnent la sensibilité de la valeur minimale de J aux variations du
niveau des contraintes Fi (au signe près). Pour un niveau des contraintes ∈ RM ,
on considère le problème d’optimisation
Lemme 2.6.6 On se place sous les hypothèses du Théorème 2.6.4 de Kuhn et Tu-
cker et on suppose que, pour petit, le problème (2.64) admet une solution unique,
notée u . On suppose aussi que la valeur minimale Jmin () de (2.64) est différen-
tiable par rapport à en 0. Alors, si on note u la solution pour = 0 et λ ∈ RM le
multiplicateur de Lagrange associé, on a, pour 1 ≤ i ≤ M ,
∂Jmin
λi = − (0).
∂i
Démonstration. D’après le Théorème 2.5.17, la solution u du problème (2.64) non
perturbé, c’est-à-dire pour = 0, vérifie la condition d’optimalité pour un certain
multiplicateur de Lagrange λ ≥ 0
2.6.3 Dualité
Donnons un bref aperçu de la théorie de la dualité pour les problèmes d’optimi-
sation. Nous l’appliquerons au problème de minimisation convexe avec contraintes
d’inégalité de la sous-section précédente. Nous avons associé à ce problème de
minimisation un problème de recherche d’un point-selle (u, p) pour le Lagrangien
L(v, q) = J(v) + q · F (v). Mais nous allons voir que, à l’existence d’un point-selle
(u, p) du Lagrangien, on peut associer inversement non pas un mais deux problèmes
d’optimisation (plus précisément, un problème de minimisation et un problème de
maximisation), qui seront dits duaux l’un de l’autre. Nous expliquerons ensuite sur
deux exemples simples en quoi l’introduction du problème dual peut être utile
pour la résolution du problème d’origine, dit problème primal (par opposition au
dual).
Revenons un instant au cadre général de la Définition 2.6.1.
Pour v ∈ U et q ∈ P , posons
J (v) = sup L(v, q) G(q) = inf L(v, q) . (2.68)
q∈P v∈U
Remarque 2.6.8 Bien sûr, sans hypothèses supplémentaires, il peut arriver que
J (v) = +∞ pour certaines valeurs de v ou que G(q) = −∞ pour certaines valeurs
de q. Mais l’existence supposée du point-selle (u, p) dans la Définition 2.6.7 nous
assure que les domaines de J et G (i.e. les ensembles {v ∈ U , J (v) < +∞} et
{q ∈ P , G(q) > −∞} sur lesquels ces fonctions sont bien définies) ne sont pas
vides, puisque (2.67) montre que J (u) = G(p) = L(u, p). Les problèmes primal et
dual ont donc bien un sens. Le résultat suivant montre que ces deux problèmes sont
étroitement liés au point-selle (u, p). •
Si le sup et l’inf sont atteints dans (2.72) (c’est-à-dire qu’on peut les écrire max
et min, respectivement), on voit alors que (2.72) traduit la possibilité d’échanger
l’ordre du min et du max appliqués au Lagrangien L. Ce fait (qui est faux si L
n’admet pas de point selle) explique le nom de min-max qui est souvent donné à un
point-selle. •
Démonstration. Soit (u, p) un point-selle de L sur U × P . Notons L∗ = L(u, p).
Pour v ∈ U , il est clair d’après (2.68) que J (v) ≥ L(v, p), d’où J (v) ≥ L∗ d’après
(2.67). Comme J (u) = L∗ , ceci montre que J (u) = inf v∈U J (v) = L∗ . On montre
de la même façon que G(p) = supq∈P G(q) = L∗ .
Réciproquement, supposons que (2.71) a lieu et posons L∗ = J (u). La définition
(2.68) de J montre que
L(u, q) ≤ J (u) = L∗ ∀q ∈ P . (2.73)
De même, on a aussi :
L(v, p) ≥ G(p) = L∗ ∀v ∈ U , (2.74)
et on déduit facilement de (2.73)-(2.74) que L(u, p) = L∗ , ce qui montre que (u, p)
est point-selle. 2
2.6. POINT-SELLE, THÉORÈME DE KUHN ET TUCKER, DUALITÉ 53
En effet, pour tout v ∈ U et q ∈ P , L(v, q) ≥ inf v0 ∈U L(v 0 , q), donc supq∈P L(v, q) ≥
supq∈P inf v0 ∈U L(v 0 , q), et puisque ceci est vrai pour tout v ∈ V , inf v∈V supq∈P L(v, q) ≥
supq∈P inf v0 ∈U L(v 0 , q), ce qui donne (2.75). La différence (positive) entre les deux
membres de l’inégalité (2.75) est appelée saut de dualité. •
Exercice 2.6.1 Soit F une fonction bornée non constante de R dans R. On définit sur
R×R le Lagrangien L(v, q) = F (v+q). Vérifier que pour ce Lagrangien l’inégalité (2.75)
est stricte avec ses deux membres finis.
Le Théorème 2.6.9 de dualité affirme que l’existence d’un point selle pour un
Lagrangien L(v, q) est équivalente à l’échange du min en v et du max en q pour son
optimisation. On peut alors se poser la question des conditions sur le Lagrangien
pour l’existence d’un point selle. Le Théorème 2.6.4 de Kuhn et Tucker donne un
premier exemple d’existence d’un point selle pour un Lagrangien du type L(v, q) =
J(v) + q · F (v), défini sur V × (R+ )M , avec J fortement convexe et F1 , . . . , FM
convexes. En effet, dans ce cas il existe un unique minimum global de J(v) sous les
contraintes F (v) ≤ 0 (si l’ensemble K, défini par (2.62) n’est pas vide). Le résultat
suivant donne une réponse plus générale sous l’hypothèse cruciale que le Lagrangien
est convexe en v et concave en q.
L(vθopt
n
, (1 − θn )q ∗ + θn q) ≤ L(v, (1 − θn )q ∗ + θn q).
L(ṽ, q ∗ ) ≤ L(v, q ∗ ).
Application
Nous appliquons ce résultat de dualité au problème précédent de minimisation
convexe avec contraintes d’inégalité convexes
ce qui montre que le problème primal inf v∈V J (v) est exactement le problème d’ori-
gine (2.78) ! D’autre part, la fonction G(q) du problème dual est bien définie par
(2.68), car (2.68) est ici un problème de minimisation convexe. De plus, G(q) est une
fonction concave car elle est l’infimum de fonctions affines (voir l’Exercice 2.3.2).
Par conséquent, le problème dual
sup G(q) ,
q∈(R+ )M
est un problème de maximisation concave plus simple que le problème primal (2.78)
car les contraintes sont linéaires ! Cette particularité est notamment exploitée dans
des algorithmes numériques (cf. l’algorithme d’Uzawa). Une simple combinaison des
Théorèmes de Kuhn et Tucker 2.6.4 et de dualité 2.6.9 nous donne le résultat suivant.
min L(v, p) .
v∈V
Précisons qu’avec les hypothèses faites il n’y a pas a priori d’unicité des solutions
pour tous ces problèmes. Précisons aussi que pour obtenir l’existence du minimum
u dans le Corollaire 2.6.13 il suffit d’ajouter une hypothèse de forte convexité ou de
comportement infini à l’infini sur J.
56 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
a une solution unique puisque v → L(v, q) est une fonction fortement convexe. Cette
solution vérifie ∂L
∂v
(v, q) = Av − b + B ∗ q = 0, soit v = A−1 (b − B ∗ q). On obtient donc
Certes, la fonctionnelle à maximiser dans (2.83) n’a pas une allure particulière-
ment sympathique. Il s’agit encore d’un problème avec fonctionnelle quadratique
et contraintes affines. Cependant, le Corollaire 2.6.13 nous assure qu’il a une so-
lution. On peut voir d’ailleurs que cette solution n’est pas forcément unique (sauf
si la matrice B est de rang M car la matrice BA−1 B ∗ est alors définie positive).
Mais l’avantage important du problème dual (2.83) vient du fait que les contraintes
(q ≥ 0) s’expriment sous une forme particulièrement simple, bien plus simple que
pour le problème primal ; et nous verrons à la Section 3.4 que cet avantage peut
être utilisé pour mettre au point un algorithme de calcul de la solution du problème
primal. •
Remarque 2.6.15 Une autre manière, très simple, d’éliminer des contraintes d’in-
égalité non-linéaires est l’introduction de variables supplémentaires, dites d’écart (ou
"slack variables" en anglais). On remplace les contraintes Fi (v) ≤ 0, pour 1 ≤ i ≤ M ,
par
Fi (v) + zi = 0 et zi ≥ 0 ,
puis on optimise par rapport au couple (v, z) ∈ V × (R+ )M avec ces nouvelles
contraintes d’égalité et d’inégalité (mais très simples pour z). Cette manière de faire
n’est donc pas très économe car on augmente le nombre de variables et de contraintes
(nous reverrons cette idée en programmation linéaire, voir le chapitre 4). •
2.6. POINT-SELLE, THÉORÈME DE KUHN ET TUCKER, DUALITÉ 57
Exercice 2.6.2 On reprend l’Exemple 1.2.8, déjà étudié à l’Exercice 2.5.4. En mé-
canique il est bien connu que la minimisation de l’énergie “en déplacements” J(u) de
l’Exercice 2.5.4 est équivalent à la minimisation d’une autre énergie, dite complémen-
taire dont la signification physique est tout aussi importante que celle de J(u). Cette
énergie complémentaire est définie en terme d’un champ de contraintes (mécaniques)
τ (x) par
1 L 2
Z
G(τ ) = |τ | dx. (2.84)
2 0
Elle s’accompagne d’une contrainte sur τ qui doit être statiquement admissible,
c’est-à-dire vérifier −τ 0 = f dans (0, L). Autrement dit, on considère le problème de
minimisation sous contrainte
1 L 2
Z
inf G(τ ) = |τ | dx . (2.85)
−τ 0 =f dans (0,L) 2 0
Montrer que la fonction duale D(v) correspondante n’est rien d’autre que l’opposée
de l’énergie −J(u) de l’Exercice 2.5.4. En admettant que −J(u) admette un point de
maximum u et que (2.85) admette un point de minimum σ, montrer que (σ, u) est un
point selle du Lagrangien et que σ = u0 .
58 CHAPITRE 2. ASPECTS THÉORIQUES DE L’OPTIMISATION
Chapitre 3
ALGORITHMES
D’OPTIMISATION
59
60 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
Remarque 3.0.1 Nous nous limitons aux seuls algorithmes déterministes et nous
ne disons rien des algorithmes de type stochastique (recuit simulé, algorithmes géné-
tiques, etc.). Outre le fait que leur analyse fait appel à la théorie des probabilités (que
nous n’abordons pas dans ce cours), leur utilisation est très différente. Pour sché-
matiser simplement, disons que les algorithmes déterministes sont les plus efficaces
pour la minimisation de fonctions convexes, tandis que les algorithmes stochastiques
permettent d’approcher des minima globaux (et pas seulement locaux) de fonctions
non convexes (à un prix toutefois assez élevé en pratique). Nous n’évoquons pas, non
plus, les algorithmes qui n’utilisent pas de dérivées car ils sont très peu efficaces en
pratique. •
un+1 = un − µn wn , (3.2)
J 0 (un )
montre que le choix de la direction wn = permet de trouver la plus petite
kJ 0 (un )k
valeur de J(un+1 ) si on néglige le terme de reste (c’est-à-dire en l’absence d’autres
informations comme les dérivées supérieures ou les itérées antérieures).
Ce simple constat conduit à choisir parmi les méthodes du type (3.2), qui sont
appelées "méthodes de descente", la méthode de gradient ou de plus grande pente
(steepest descent en anglais) dans laquelle wn est positivement proportionnelle à
3.1. ALGORITHMES DE TYPE GRADIENT (SANS CONTRAINTES) 61
Notons que l’on n’a pas normalisé le vecteur gradient dans (3.3). Cet algorithme
converge comme l’indique le résultat suivant.
Théorème 3.1.1 On suppose que J est α-convexe différentiable et que J 0 est Lip-
schitzien sur tout borné de V , c’est-à-dire que
∀ M > 0, ∃ LM > 0, kvk + kwk ≤ M ⇒ kJ 0 (v) − J 0 (w)k ≤ LM kv − wk . (3.5)
Alors l’algorithme de gradient à pas optimal converge : quel que soit u0 , la suite (un )
définie par (3.3) et (3.4) converge vers la solution u de (3.1).
Démonstration. La fonction j(µ) = J un − µ J 0 (un ) est fortement convexe (si
un+1 = un − µn P J 0 (un )
conduit à
si µn > 0 est petit et J 0 (un ) 6= 0. Il peut être intéressante d’utiliser un tel précondi-
tionnement du gradient pour améliorer la convergence. Par exemple, si les ordres
de grandeur des différentes composantes de u, et donc de J 0 (u), sont très différents,
on peut utiliser une telle matrice P diagonale pour remettre à la même échelle ces
composantes. Plus généralement, le choix d’une telle matrice revient à changer le
produit scalaire de V avec lequel on identifie le gradient. De ce point de vue, cette
approche se généralise aux espaces de Hilbert V de dimension infinie : on change la
direction de descente si on change le produit scalaire de V . •
est la règle d’Armijo-Goldstein qui, pour des paramètres 0 < p− < p+ < 1 fixés,
demande que le pas µ̃n > 0 vérifie
où on rappelle que j 0 (0) < 0 tant qu’on n’a pas convergé exactement. Il existe
d’autres règles du même type, comme celle de Wolfe (voir [6]). •
3.1. ALGORITHMES DE TYPE GRADIENT (SANS CONTRAINTES) 63
à partir d’une initialisation u0 ∈ V . Cette méthode est donc plus simple que l’al-
gorithme de gradient à pas optimal, puisqu’on fait à chaque étape l’économie de
la résolution de (3.4). Le résultat suivant affirme que, si le pas de descente µ est
suffisamment petit, alors l’algorithme (3.8) converge.
Théorème 3.1.6 On suppose que J est α-convexe différentiable et que J 0 est Lip-
schitzien sur V , c’est-à-dire qu’il existe une constante L > 0 telle que
Alors, si 0 < µ < 2α/L2 , l’algorithme de gradient à pas fixe converge : quel que soit
u0 , la suite (un ) définie par (3.8) converge vers la solution u de (3.1)
p
kun − uk ≤ γ n ku0 − uk avec γ = 1 − 2αµ + µ2 L2 < 1.
64 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
n n 0 n+1
Démonstration. Posons v = u − u. Comme J (u) = 0, on a v = vn −
0 n 0
µ J (u ) − J (u) , d’où il vient
2
kv n+1 k2 = kv n k2 − 2µhJ 0 (un ) − J 0 (u), un − ui + µ2
J 0 (un ) − J 0 (u)
(3.10)
2 2 n 2
≤ 1 − 2αµ + L µ kv k ,
d’après (3.9) et l’α-convexité. Si 0 < µ < 2α/L2 , il est facile de voir que 1 − 2αµ +
L2 µ2 ∈]0, 1[, et la convergence se déduit de (3.10). De manière équivalente, la même
démonstration montre que l’application v 7→ v − µJ 0 (v) est strictement contractante
lorsque 0 < µ < 2α/L2 , donc elle admet un unique point fixe (qui n’est autre que
u) vers lequel converge la suite un . 2
Remarque 3.1.7 L’algorithme de gradient à pas fixe est plus simple que celui à pas
optimal puisqu’il ne nécessite pas de résoudre une optimisation uni-dimensionnelle
(3.4) à chaque itération. Par contre, on ne connait en général pas de bonne estimation
des constantes α et L et on ne sait pas a priori si un choix de pas de descente µ > 0
est suffisamment petit. •
kJ 0 (un )k ≤ C0 γ n ,
Proposition 3.1.9 Soit une fonction J(u) définie sur V = RN , strictement convexe,
différentiable, infinie à l’infini et telle que J 0 est Lipschitzien, c’est-à-dire vérifie
(3.9) pour une constante L > 0. Alors, si 0 < µ < 2/L, l’algorithme de gradient à
pas fixe converge : quel que soit u0 , la suite (un ) définie par (3.8) converge vers la
solution u de (3.1).
Démonstration. Puisque un+1 = un − µJ 0 (un ) on applique le Lemme 3.1.10 avec
v = un et v ∗ = un+1 , ce qui conduit à
On en déduit que, pour 0 < µ < 2/L, la suite J(un ) est décroissante et, comme J
est infinie à l’infini, la suite un est bornée. Comme V = RN il existe une sous-suite
0 0 0 0
convergente, notée un , de la suite un . De la formule un +1 − un = −µJ 0 (un ), on
0
déduit que limn0 →+∞ J 0 (un ) = 0. Comme J 0 est continue (puisque Lipschitzien), la
0
limite de la sous-suite un est un zéro de J 0 . Or, J est strictement convexe donc
admet au plus un point de minimum, c’est-à-dire que J 0 ne s’annule qu’en un seul
point u. Par conséquent, toutes les sous-suites convergentes de un convergent vers
la même limite u. Autrement dit, toute la suite un converge vers u qui est l’unique
point de minimum de J. 2
Démonstration. En vertu du Lemme 2.4.10 la fonction J(v ∗ ) est majorée par une
fonction quadratique de v ∗
L ∗
J(v ∗ ) ≤ J(v) + hJ 0 (v), v ∗ − vi + kv − vk2 .
2
On remplace v ∗ − v par −µJ 0 (v) pour obtenir le résultat. 2
Remarque 3.1.11 Le Lemme 3.1.10 est une garantie de décroissance d’un pas de
l’algorithme du gradient, sans aucune condition de convexité sur la fonction. •
Si on ne suppose même plus que la fonction J est convexe, alors on ne peut pas
garantir que la suite un converge vers un point de minimum.
Exercice 3.1.6 Soit une fonction J(u) différentiable, infinie à l’infini et telle que J 0
est Lipschitzien, c’est-à-dire vérifie (3.9) pour une constante L > 0. On ne suppose pas
que J est convexe. A l’aide des mêmes raisonnements que pour la Proposition 3.1.9,
montrer que l’algorithme du gradient à pas fixe 0 < µ < 2/L vérifie que J(un ) décroit
et J 0 (un ) tend vers 0.
L’exercice suivant donne un exemple de fonction J pour lequel l’algorithme
de gradient converge seulement vers un point d’inflexion et pas vers un point de
minimum.
Exercice 3.1.7 Soit la fonction J(x) = 21 x21 + 14 x42 − 21 x22 définie sur R2 . Vérifier que
J n’est pas convexe et que l’ensemble de ses points critiques (où son gradient s’annule)
sont O = (0, 0), P = (0, −1) et Q = (0, 1). Montrer que P et Q sont des points de
minimum tandis que O est un point selle. Montrer que, si on initialise un algorithme de
gradient au point x0 = (1, 0), alors quelque soit la règle de choix du pas de descente
µn > 0 (fixe, optimal ou autre) qui assure la convergence la limite est le point O qui
n’est pas un minimum, même local, de J.
(on parle alors de méthode multi-pas). Nous en donnons quelques exemples simples
pour illustrer la richesse de l’algorithmique en optimisation, sachant que certains
algorithmes sont spécialisés pour des classe de fonctions coûts bien particulières.
1 − α/L
kun − uk ≤ γ nopt ku0 − uk avec γ opt = < 1.
1 + α/L
Démonstration. On écrit à nouveau
L’hypothèse (3.9) sur J 0 Lipschitzien implique J 00 ≤ L Id, tandis que la forte convexité
de J implique J 00 ≥ α Id. Par conséquent, la matrice symétrique B n vérifie
et dans cette majoration la valeur optimale du pas est µ opt = 2/(L + α) qui conduit
L−α
à γ opt = L+α . 2
68 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
Remarque 3.2.2 Comme l’a montré l’Exercice 3.1.5 sur une fonction quadratique,
on ne peut pas espérer améliorer la vitesse de convergence de l’algorithme du gradient
à pas fixe obtenue à la Proposition 3.2.1, qui est donc bien optimale. De même,
l’Exercice 3.1.5 montre que l’on ne peut pas élargir l’ensemble des pas, 0 < µ < 2/L,
pour lesquels l’algorithme du gradient à pas fixe converge. •
Une propriété intéressante des fonctions J, α-convexes à dérivée L-Lipschitzienne,
est que la décroissance de la norme de la dérivée J 0 (un ) est équivalente à la décrois-
sance de la norme de l’erreur (un − u) où u est l’unique point de minimum de J
et un est la suite des itérées produite par n’importe quel algorithme d’optimisation.
Autrement dit, il suffit d’étudier la convergence vers zéro de la dérivée J 0 (un ) pour
avoir une idée de la convergence vers zéro de l’erreur (un − u), sans même connaitre
le point de minimum. C’est une conséquence du résultat suivant.
Démonstration. Notons que J 0 (u) = 0. La première inégalité est alors une consé-
quence de (2.24) et de l’inégalité de Cauchy-Schwarz, tandis que la deuxième inéga-
lité est simplement l’application de la Définition 2.4.9 sur la L-Lipschitzianité de J 0 .
2
Si la fonction J que l’on minimise n’est pas fortement (ou α-) convexe, alors la
vitesse de convergence est bien moindre, devenant algébrique au lieu de géométrique.
Proposition 3.2.4 Soit une fonction J(u) définie sur V = RN , strictement convexe,
de classe C 2 , infinie à l’infini et telle que J 0 est Lipschitzien, c’est-à-dire vérifie (3.9)
pour une constante L > 0. Alors l’algorithme de gradient à pas fixe converge pour
0 < µ < 2/L, c’est-à-dire que quel que soit u0 , la suite (un ) définie par (3.8) converge
vers la solution u de (3.1), et on a
κ ku0 − uk2
0 ≤ J(un ) − J(u) ≤ avec κ = .
n+1 µ(1 − Lµ/2)
Remarque 3.2.5 Noter que la Proposition 3.2.4 ne dit rien de la vitesse de conver-
gence de kun − uk et que la vitesse de convergence est ici définie comme la vitesse
de convergence de J(un ) vers J(u). Cette dernière est algébrique (comme 1/n),
ce qui est beaucoup moins rapide que la convergence géométrique du cas fortement
convexe (comme γ 2n avec 0 < γ < 1, d’après l’Exercice 3.1.4). On ne peut pas avoir
une estimation de kun − uk, indépendante de J, car les fonctions convexes, mais pas
fortement convexes, peuvent être très plates près de leur minimum et la convergence
de un vers u peut être aussi lente que l’on veut selon la forme de J. •
Démonstration. On a déjà démontré à la Proposition 3.1.9 que un converge vers
l’unique point de minimum u de J. On calcule la norme au carré de la relation
et l’on majore le dernier terme à l’aide du Lemme 3.2.6 ci-dessous. Comme J 0 (u) = 0,
on en déduit
et donc que la suite kun − uk est décroissante pour 0 < µ < 2/L. Par ailleurs, la
convexité de J en un donne
ku0 − uk2
∆n+1 ≤ ∆n − κ−1 (∆n )2 avec κ = . (3.13)
µ(1 − Lµ/2)
κ
Montrons alors par récurrence sur n que ∆n ≤ n+1
, ce qui est le résultat recherché.
Pour n = 0, on reprend (3.11) et on majore
ku0 − uk
kJ 0 (u0 )k = kJ 0 (u0 ) − J 0 (u)k ≤ Lku0 − uk ≤
µ(1 − Lµ/2)
Lemme 3.2.6 Soit une fonction J(u) définie sur V = RN , strictement convexe, de
classe C 2 et telle que J 0 est L-Lipschitzien au sens de la Définition 2.4.9 pour une
constante L > 0. Alors
qui vérifie A−1 ≥ 1/L Id. On déduit alors de la formule de Taylor ci-dessus
1 0
kJ (v) − J 0 (w)k2 ≤ hA−1 (J 0 (v) − J 0 (w)) , J 0 (v) − J 0 (w)i = hJ 0 (v) − J 0 (w), v − wi.
L
En fait, le résultat reste vrai pour une fonction J(u) convexe mais pas strictement
convexe : il suffit de lui ajouter kuk2 avec > 0 puis de passer à la limite → 0
dans le résultat. 2
Proposition 3.2.8 Soit J une fonction convexe différentiable, telle que J 0 est L-
Lipschitzienne, c’est-à-dire vérifie (3.9). Soit un pas de descente µ tel que 0 < µ ≤
1/L. On suppose qu’il existe un point de minimum u de J sur V . Alors la suite un
de l’agorithme de Nesterov, définie par (3.14), vérifie
κ
0 ≤ J(un ) − J(u) ≤ avec κ = 2µ−1 ku0 − uk2 + 4(J(u0 ) − J(u)).
(n + 2)2
Remarque 3.2.9 Avec quelques hypothèses supplémentaires on peut montrer que
la suite un converge vers u (voir l’Exercice 3.2.1 ci-dessous). La vitesse de convergence
de la Proposition 3.2.8 améliore sensiblement (1/n2 au lieu de 1/n) celle obtenue
dans la Proposition 3.2.4. Comme on le verra à la Remarque 3.2.10 cette vitesse
est optimale pour les fonctions convexes. Néanmoins, pour des fonctions fortement
convexes, quitte à changer la définition précise du coefficient d’extrapolation an ,
on peut encore améliorer cette vitesse de convergence et retrouver une convergence
géométrique (voir la Proposition 3.2.11). •
Démonstration. Soit pn = (an − 1)(un−1 − un ) de sorte que v n+1 = un − pn /an+1 .
L’astuce de Nesterov est d’étudier la convergence de la suite un − pn plutôt que celle
de un ou v n . On calcule
un+1 − pn+1 = un+1 − (an+1 − 1)(un − un+1 ) = an+1 un+1 − (an+1 − 1)un
= an+1 (v n+1 − µJ 0 (v n+1 )) − (an+1 − 1)un = un − pn − an+1 µJ 0 (v n+1 ).
On en déduit
kun+1 − pn+1 − uk2 = kun − pn − uk2 + a2n+1 µ2 kJ 0 (v n+1 )k2
+2an+1 µhJ 0 (v n+1 ), u − v n+1 i (3.15)
+2(an+1 − 1)µhJ 0 (v n+1 ), pn i,
72 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
en ayant utilisé an+1 (v n+1 +pn −un ) = (an+1 −1)pn dans la dernière ligne. On majore
la deuxième ligne de (3.15) par convexité de J
kun+1 − pn+1 − uk2 + 2µa2n+1 (J(un+1 ) − J(u)) ≤ kun − pn − uk2 + 2µa2n (J(un ) − J(u))
donc, comme a0 = 1 et p0 = 0,
dont l’unique point de minimum est u = (1, ..., 1, 0, ..., 0), qui a ses p premières
composantes égales à 1 et les suivantes égales à 0, avec J(u) = 0. On initialise
un algorithme du gradient avec u0 = 0. Tous les algorithmes de gradient que nous
étudions dans cette section (pas optimal ou fixe, Nesterov, boule pesante ou gradient
conjugué) vérifient la propriété suivante : à chaque itération n, la solution approchée
un appartient à l’espace affine u0 + En où En est le sous-espace vectoriel engendré
par {J 0 (u0 ), J 0 (u1 ), ..., J 0 (un−1 )}. Par conséquent, une borne inférieure est
J(un ) ≥ min J(v).
v∈u0 +En
Or, dans l’exemple ci-dessus, il est facile de voir, par récurrence sur n, que J 0 (un−1 )
appartient à l’espace engendré par les n premiers vecteurs de la base canonique,
{e1 , ..., en }, et par conséquent que En ⊂ Vect{e1 , ..., en } donc
J(un ) ≥ min J(v) ≥ min J(v).
v∈u0 +En v∈u0 +Vect{e1 ,...,en }
Autrement dit, sur cet exemple (très défavorable !) l’algorithme de gradient (quelle
que soit sa version) ne peut propager l’information que la solution u a ses p premières
composantes égales à 1, à partir d’une initialisation nulle u0 = 0, qu’à la "vitesse"
d’une composante par itération. Un calcul facile permet de calculer le minimum sur
le sous-espace affine u0 + Vect{e1 , ..., en } : ce point de minimum, noté wn , a ses
i
dernières N − n composantes nulles tandis que les n premières sont win = 1 − n+1
1
pour 1 ≤ i ≤ n. Par conséquent (w1n − 1)2 = (win − wi−1 n
)2 = (wnn )2 = (n+1)2 pour
On choisit alors t = an+1 , où an est la suite définie par (3.19). Grâce aux propriétés
du Lemme 3.2.14, on a an+1 ≥ 1 ainsi que
avec
1
∆n = a2n (J(un ) − J(u)) + ku − un − (an − 1)(un − un−1 )k2 ,
2µ
76 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
c’est-à-dire que
n
Y
∆n ≤ ∆0 (1 − ak αµ).
k=1
√
Comme 0 ≤ 1 − ak αµ ≤ 1 et que limn→+∞ an = 1/ q, on en déduit que ∆n
tend vers 0 et l’algorithme converge. De plus, si on choisit le pas µ = 1/L, comme
a2
1 − a1n = (1 − qan ) n−1
a2
, on a
n
n n n
Y Y a2n Y 1 √
(1 − ak αµ) = (1 − qak ) = 2 (1 − ) ≤ a2n (1 − q)n ,
k=1 k=1
a0 k=1 ak
√ n
d’où l’on déduit en particulier J(un ) − J(u) ≤ (1 − q) . Alors la forte convexité
de J et le fait que J 0 (u) = 0 donne
α n √
ku − uk2 ≤ J(un ) − J(u) ≤ (1 − q)n ,
2
ce qui conclut la démonstration. 2
Lemme 3.2.13 Soit J une fonction α-convexe sur V dont la dérivée J 0 est Lip-
schitzienne de constante L. Pour tout v, w ∈ V et pour tout 0 < µ ≤ 1/L, en notant
v ∗ = v − µJ 0 (v), on a
1 − αµ 1 ∗
J(w) + kw − vk2 ≥ J(v ∗ ) + kv − wk2 . (3.21)
2µ 2µ
Démonstration. Pour tout v, w ∈ V la forte convexité de J s’écrit
α
J(w) ≥ J(v) + hJ 0 (v), w − vi + kw − vk2 (3.22)
2
auquel on ajoute un terme kw − vk2 /(2µ) de part et d’autre pour obtenir
1 − αµ 1
J(w) + kw − vk2 ≥ J(v) + hJ 0 (v), w − vi + kw − vk2 . (3.23)
2µ 2µ
On note Q(w) la fonction convexe quadratique à droite de (3.23)
1
Q(w) = J(v) + hJ 0 (v), w − vi + kw − vk2
2µ
qui admet v ∗ = v − µJ 0 (v) comme unique point de minimum. Puisque Q(w) est
quadratique et que Q0 (v ∗ ) = 0, on a
1 ∗
Q(w) = Q(v ∗ ) + kv − wk2 . (3.24)
2µ
En combinant (3.23) et (3.24) on obtient
1 − αµ 1 1
J(w)+ kw −vk2 ≥ J(v)+hJ 0 (v), v ∗ −vi+ kv ∗ −vk2 + kv ∗ −wk2 . (3.25)
2µ 2µ 2µ
3.2. GÉNÉRALISATIONS ET AUTRES ALGORITHMES DE TYPE GRADIENT77
car µ ≤ 1/L. 2
avec la donnée initiale u(0) = u0 (la notation u̇ signifie la dérivée par rapport à la
variable de temps t). En effet, si on assimile le pas de descente µ > 0 à un pas de
temps, alors un peut être interprété comme une approximation au temps tn = nµ de
la solution u(t) (voir [1] pour des détails sur cette notion). L’équation différentielle
ordinaire (3.30) modélise la trajectoire d’un point qui descendrait le long du potentiel
J(u) (on peut vérifier que J(u(t)) est décroissant le long de la trajectoire ; c’est ce
qu’on appelle une fonction de Lyapunov). Si la fonction J est fortement convexe, on
peut aussi montrer que la trajectoire u(t) converge, quand t tend vers l’infini, vers
l’unique point de minimum de J.
Expliquons maintenant la motivation et l’origine de l’algorithme de la boule
pesante dans ce contexte. Si la fonction J n’est pas convexe et présente des minima
locaux, la solution u(t) peut converger vers un tel minimum local, loin d’un minimum
global. Pour améliorer cette situation défavorable, une idée est de changer l’équation
différentielle, pour y rajouter un terme d’inertie, qui permettrait au point matériel
solution de cette nouvelle équation de sortir d’un minimum local en remontant la
pente grâce à l’inertie acquise. Autrement dit, l’équation différentielle de cette boule
pesante serait
mü(t) + u̇(t) = −J 0 (u(t)),
où m > 0 serait la masse de cette boule. Par analogie avec une discrétisation explicite
en temps (qui nécessite 3 points en temps pour approcher la dérivée seconde), on
obtiendrait
un+1 = un − µ J 0 (un ) + ν(un − un−1 ), (3.29)
où µ > 0 et ν > 0 sont deux paramètres positifs. Bien sûr, il faut joindre à cette
récurrence pour n ≥ 1 deux initialisations u1 , u0 ∈ V . Pour simplifier l’analyse on
se place en dimension finie.
c’est-à-dire que
d’où l’on déduit la convergence de z n vers 0 pour δ < (1 − ρmax )/2 (et aussi le fait
que si un est proche de u à un certain rang, cela restera vrai pour les itérations
suivantes).
Or ρ est valeur propre dans (3.31) si Hx = λx avec λ = (1 + ν − ρ − ν/ρ)/µ
où λ est une valeur propre de la Hessienne H. Autrement dit, ρ est une racine du
polynome
ρ2 − (1 + ν − µλ)ρ + ν = 0,
dont le discriminant est ∆ = (1+ν −µλ)2 −4ν. On cherche des conditions suffisantes
pour que les racines vérifient |ρ| < 1. On choisit ν ≥ 0 et comme ν est le produit
des racines du polynome, il faut √ que ν < 1. Si ∆ ≤ 0, alors les deux racines sont
complexes et de module égal à ν < 1. Si ∆ > 0, alors un calcul facile montre que
les deux racines appartiennent à l’intervalle (−1, +1) si et seulement si 0 < µ <
2(1 + ν)/λ pour toute valeur propre λ de la Hessienne H. A cause de l’hypothèse
(3.30) les valeurs propres de H vérifient λ ∈ [α, L]. Par conséquent, les conditions
0 < ν < 1 et 0 < µ < 2(1 + ν)/L assurent la convergence de l’algorithme. Pour
80 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
trouver des valeurs optimales des paramètres, une analyse (un peu longue et pas
détaillée ici) montre
√ qu’il est bon de se placer dans le cas des racines complexes de
module égal à ν, ce qui arrive lorsque ∆ ≤ 0 . On vérifie que ∆ ≤ 0 est équivalent
à √ √
(1 − ν)2 ≤ µλ ≤ (1 + ν)2 ,
qui doit être vrai pour toute valeur propre λ ∈ [α, L] de H. C’est donc vrai en
particulier si √ √
√ p
ν ≥ αµ − 1, ν ≥ 1 − Lµ.
√ √ √ √
On choisit µ = 2/( α+ L) pour que les deux bornes inférieures de √ν coincident
√ √
et on prend la plus petite valeur qui en résulte, c’est-à-dire ν = √α− L
√ .
α+ L
2
On voit donc, comme 0 < α/L < 1, que la méthode de la boule pesante avec ces
paramètres optimaux converge plus vite que les algorithmes précédents. Malheu-
reusement, les valeurs de α et L ne sont souvent pas connues en pratique et il est
difficile d’ajuster les valeurs optimales des paramètres µ et ν. •
n kJ 0 (un+1 )k2
β = .
kJ 0 (un )k2
3.2. GÉNÉRALISATIONS ET AUTRES ALGORITHMES DE TYPE GRADIENT81
Expliquons d’où viennent ces formules en quelques mots. Tout d’abord, la recherche
linéaire pour trouver le pas µn est classique et reminiscente de l’algorithme du gra-
dient à pas optimal. La formule pour β n s’inspire du cas où J(u) est quadratique.
Si J(u) est une fonction non-linéaire quelconque, l’algorithme du gradient conju-
gué fonctionnera néanmoins, approximativement comme dans le cas quadratique,
au voisinage d’un point de minimum non dégénéré (c’est-à-dire où la Hessienne est
définie positive). Nous allons donc formellement supposer que la Hessienne J 00 (u)
est indépendante de u et égale à une matrice constante A. L’idée essentielle est de
demander à ce que les directions de descente pn et pn+1 soient conjuguées par rap-
port à la matrice A, c’est-à-dire que pn+1 · Apn = 0. Comme expliqué ci-dessous (voir
l’Exercice 3.2.2), cette propriété est à la source de la convergence de l’algorithme du
gradient conjugué pour une fonction quadratique et on essaye donc de la reproduire
dans le cas général.
On souhaite donc que la valeur de β n soit telle que
pn+1 · Apn = 0.
Par ailleurs,
d’où l’on déduit la formule de Polak-Ribière en reportant (3.34) et (3.35) dans (3.33).
Un calcul similaire conduit à la formule de Fletcher-Rieves. Remarquons que ces
deux formules ne sont pas égales dans le cas général mais qu’elles coincident pour
une fonction J(u) quadratique.
En fait, la méthode du gradient conjugué a d’abord été inventée pour résoudre
des systèmes linéaires dont la matrice est symétrique, définie positive, c’est-à-dire
pour minimiser sur RN une fonction quadratique
1
J(x) = Ax · x − b · x
2
avec b ∈ RN et A une matrice symétrique définie positive. C’est dans ce cas qua-
dratique que l’algorithme du gradient conjugué est le plus efficace puisqu’on peut
82 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
krn k2 krn+1 k2
µn = et βn = .
Apn · pn krn k2
Exercice 3.2.2 Soit A une matrice symétrique définie positive. Soit x0 ∈ RN . Pour
la récurrence (3.36) démontrer que rn = Axn − b. Montrer que la suite pn vérifie
Apn · pj = 0 pour 0 ≤ j < n (on dit que la suite pn est conjuguée par rapport à A).
Vérifier que les formules de Fletcher-Rieves et Polak-Ribière coincident dans ce cas.
Le fait que la suite des directions de descente pn soit conjuguée par rapport
à A donne son nom à l’algorithme. C’est une propriété cruciale pour l’efficacité de
la méthode car elle assure que la nouvelle direction de descente pn “travaille” de
manière orthogonale (au sens du produit scalaire hx, yiA = Ax · y) par rapport aux
directions précédentes. Le nouveau gain de minimisation ne nuit pas aux gains passés
(pas d’effet de convergence en zig-zag) et c’est la raison du résultat de convergence
exacte en moins de N itérations. Néanmoins, la convergence approchée peut avoir
lieu en beaucoup moins d’itérations comme le montre le résultat suivant que nous
admettrons.
Proposition 3.2.18 Soit A une matrice symétrique réelle définie positive, de va-
leurs propres ordonnées 0 < λ1 ≤ · · · ≤ λN . Soit x la solution exacte du système
Ax = b. Soit (xn )n≥0 la suite de solutions approchées du gradient conjugué, définie
par (3.36). Alors
p !n
p 1− λ1 /λN
kxn − xk ≤ 2 λN /λ1 p kx0 − xk2 .
1 + λ1 /λN
Exercice 3.2.3 Soit J(x) = |x| définie de R dans R. Calculer son sous-gradient ∂J(0).
Une classe importante de fonctions convexes, qui ne sont pas différentiables
mais admettent un sous-gradient que l’on sait calculer, est celle des fonctions qui
sont des maxima de fonctions convexes différentiables. C’est une situation courante
puisque, par exemple, on sait, d’après l’Exercice 2.3.2, que toute fonction convexe
est le supremum des fonctions affines qui la minorent. Plus généralement, pour un
espace (non vide) de paramètres Λ, on définit
J (x) = sup Jλ (x), (3.38)
λ∈Λ
Lemme 3.2.22 Supposons Λ fini, soit J la fonction définie par (3.38) avec des
fonctions Jλ (x) convexes et différentiables, et pour tout x ∈ Rn , posons Γ(x) = {λ ∈
Λ | Jλ (x) = J (x)}. Alors, pour tout λ ∈ Γ(x), Jλ0 (x) est un sous-gradient de J au
point x.
Démonstration. Pour tout λ ∈ Γ(x), et pour tout y ∈ Rn , on a J (y) − J (x) ≥
Jλ (y) − Jλ (x) ≥ Jλ0 (x) · (y − x), par convexité de Jλ , ce qui démontre le résultat. 2
Remarque 3.2.23 Les fonctions non régulières du type (3.38) se retrouvent dans au
moins deux contextes importants. D’une part, elles sont utilisées dans les méthodes
de relaxations Lagrangiennes pour l’optimisation combinatoire ou en variables en-
tières (voir [5]). D’autre part, elles correspondent à l’approche, dite du pire des
cas, en optimisation robuste. Supposons que l’on veuille optimiser un système décrit
par une variable x. Mais la description de ce système est entachée d’incertitudes,
d’erreurs ou de données inconnues, caractérisées par λ. On peut identifier chaque
valeur de λ à un scenario possible de fonctionnement du système. Si l’on veut optimi-
ser pour tous les scenarios possibles, l’approche du pire des cas consiste à minimiser
le maximum par rapport à la variable λ. Notons que c’est une approche assez pes-
simiste (on prévoit le pire !) et qu’elle est parfois remplacée par une approche en
moyenne où on optimise la moyenne de Jλ (x) sur Λ (ou bien une somme pondérée
de la moyenne et de la variance). •
L’algorithme de sous-gradient pour minimiser la fonction convexe J consiste,
pour une initialisation x0 ∈ Rn , à construire la suite
ρk
xk+1 = xk − pk , (3.39)
kpk k
où pk est un sous-gradient quelconque de J au point xk , et où ρk > 0 est une suite
de réels strictement positifs telle que
X X
ρk → 0, ρi = +∞, ρ2i < +∞ . (3.40)
i∈N i∈N
Remarque 3.2.25 Les hypothèses de la Proposition 3.2.24 sont vérifiées pour les
fonctions J définies par (3.38) avec un ensemble Λ fini, chacune des fonctions Jλ (x)
étant convexe différentiable et l’une d’entre elle infinie à l’infini. •
3.2. GÉNÉRALISATIONS ET AUTRES ALGORITHMES DE TYPE GRADIENT85
On en déduit tout d’abord que kxi+1 −x∗ k est bornée à cause de la dernière condition
de (3.40), donc la suite xk est bornée par une constante M indépendante de k. On
applique alors l’hypothèse que J est Lipschitzienne à la définition du sous-gradient
pk
pk · (y − xk ) ≤ J(y) − J(xk ) ≤ LM ky − xk k
d’où l’on déduit en variant y autour de xk que kpk k ≤ LM . Par conséquent (3.42)
implique que
LM kx0 − x∗ k2 + ∞ 2
P
k=0 ρk
min (J (xk ) − J (x∗ )) ≤ Pi ,
0≤k≤i 2 k=0 ρk
ce qui démontre que lim inf k→+∞ J (xk ) = J (x∗ ) puisque la série des ρk diverge. Soit
xk0 une sous-suite de xk telle que limk→+∞ J (xk0 ) = J (x∗ ). Comme elle est bornée
dans Rn , il existe encore une sous-suite, toujours notée xk0 par souci de simplicité,
qui converge vers une limite x0∗ qui, par continuité de J est un point de minimum
de J. Cette limite x0∗ n’est peut-être pas le point de minimum x∗ introduit ci-dessus
mais rien ne nous empêche de remplacer x∗ par x0∗ et tout le raisonnement reste
le même car les itérées xk et pk ne dépendent pas du choix x∗ ou x0∗ . Par souci de
simplicité dans les notations, on réécrit x0∗ = x∗ . Montrons que toute la suite xk
converge en fait vers x∗ .
Pour δ > 0 petit fixé, il existe un indice k 0 , suffisamment grand, tel que
X
kxk0 − x∗ k2 ≤ δ/2 et ρ2i ≤ δ/2.
i≥k0
Proposition 3.2.29 On suppose que F (x) est α-convexe et qu’il existe une constante
C, indépendante de n, telle que pour tout x ∈ Rd
n
1X 0
kfi (x)k2 ≤ C(1 + kx − x∗ k2 ) (3.45)
n i=1
Remarque 3.2.30 Les hypothèses sur F sont vérifiées dans le cas de l’Exemple
1.2.6 (car chacune des dérivées fi0 (x) est bornée) ou bien si chaque fonction fi est
quadratique et fortement convexe, uniformément par rapport à i. •
Démonstration. On réécrit (3.44) sous la forme
kxk+1 − x∗ k2 = kxk − x∗ k2 − 2µk (xk − x∗ ) · fi0k (xk ) + (µk )2 kfi0k (xk )k2 . (3.46)
On prend l’espérance de (3.46) par rapport à la seule variable aléatoire qui gouverne
le choix de l’indice ik , en notant que xk ne dépend pas de cette variable aléatoire,
pour obtenir
n
k+1 ∗ 2
k ∗ 2 k k ∗ 0 k (µk )2 X 0 k 2
E kx − x k = kx − x k − 2µ (x − x ) · F (x ) + kf (x )k
n i=1 i
où, cette fois-ci, E désigne l’espérance pour toutes les variables aléatoires des indices
de i0 à ik , et avec
Yk
k
Π = ρj .
j=0
ce qui prouve la convergence en moyenne (un résultat similaire avec d’autres ex-
posants s’obtient lorsque α ≥ 1/2). On peut démontrer la convergence presque
sûrement mais cela dépasse la cadre de ce cours. 2
Comme pour l’algorithme du sous-gradient la convergence est particulièrement
lente (algébrique au lieu de géométrique pour le gradient à pas fixe). Du coup, il
n’est pas clair que cet algorithme soit préférable à celui du gradient. Mais il ne faut
pas oublier que si n est très grand, le coût d’une itération de gradient stochastique
est approximativement n plus fois faible que le coût d’une itération du gradient
puisqu’on n’évalue qu’une seule dérivée fi0 au lieu de n (pour les deux algorithmes
on calcule aussi régulièrement la fonction objectif F pour vérifier sa décroissance et la
bonne convergence des itérations). Par ailleurs, il existe des stratégies heuristiques du
choix du pas de descente µk qui peuvent améliorer en pratique le choix théorique ci-
dessus. En particulier, l’algorithme du gradient stochastique est souvent plus rapide
lors des premières itérations que l’algorithme du gradient et surtout plus insensible
au "bruit" (qu’il soit numérique ou dans les données).
Plus τ est grand et plus la solution x aura des composantes nulles. Au contraire,
à la limite τ → 0+ , il s’agit d’un problème de moindres carrés standard dont on
sait qu’il admet toujours au moins une solution qui est unique si KerA = {0} (cf.
Exercice 2.5.2).
Remarque 3.2.32 Expliquons sur un cas particulier très simple pourquoi la norme
`1 dans (3.49) conduit à ce que la solution ait de nombreuses composantes nulles.
On choisit p = n et A = Id. Dans ce cas, la minimisation de (3.49) se fait “à la
main”, composante par composante. La fonction de R dans R,
1
xi 7→ |xi − bi |2 + τ |xi |,
2
est fortement convexe, dérivable partout sauf en 0, et admet donc un unique point
de minimum x∗i qui vérifie, soit x∗i = 0, soit la condition d’optimalité
On en déduit facilement l’unique point de minimum qui est donné par l’opérateur
de contraction (shrinkage, en anglais)
∗ 0 si |bi | < τ,
xi = Sτ (bi ) = (3.50)
bi − sgn(bi )τ si |bi | ≥ τ.
L’interprétation est simple si l’on pense que la donnée b est bruitée et que τ corres-
pond au seuil maximal du bruit. Si |bi | < τ , la donnée est indistinguable du bruit et
on obtient x∗i = 0. Sinon, on obtient x∗i ≈ bi à une correction de l’ordre de τ près.
La solution x∗ ∈ Rn est donc d’autant plus parcimonieuse que le seuil τ est élevé. •
Pour résoudre (3.49) on ne peut pas employer une méthode de gradient classique
car la norme `1 , x 7→ kxk1 , n’est pas dérivable en tout point dont une des compo-
santes s’annule. On pourrait utiliser l’algorithme du sous-gradient car la fonction
objectif est convexe mais on a vu que cet algorithme converge lentement. Il y a
donc de la place pour concevoir un algorithme plus efficace. Plus généralement, on
considère la minimisation de la somme de deux fonctions convexes sur Rn
si xk+1 en est le point de minimum. Or, si (3.51) n’a pas de sens pour une fonction
J2 non différentiable, au contraire (3.52) admet bien un unique point de minimum,
noté xk+1 , dès que J2 est convexe, puisque la fonction à minimiser est fortement
convexe. Ainsi, (3.52) permet de définir l’application proximale de τ J2
proxτ J2 (xk ) = xk+1 . (3.53)
Avec une initialisation x0 ∈ Rn , (3.53) définit un algorithme de minimisation pour
J2 , appelé algorithme proximal. On vérifie bien qu’il s’agit d’un algorithme de
descente car J2 (xk+1 ) < J2 (xk ), sauf si xk+1 = xk auquel cas la suite xm est station-
naire pour m ≥ k. Cependant, cet algorithme est toujours de type conceptuel car
pour une fonction convexe quelconque J2 on ne sait pas résoudre le problème de mi-
nimisation (3.52) dont la solution est requise pour définir l’algorithme. En pratique,
cet algorithme est donc limité à des fonctions simples pour lesquelles on sait résoudre
“à la main” (3.52). C’est précisément le cas pour J2 (x) = kxk1 , comme on vient de
le voir à la Remarque 3.2.32. Dans ce cas, on a xk+1 = proxτ k·k1 (xk ) = Sτ (xk ), où
Sτ est l’opérateur de contraction défini par (3.50).
Finalement, on propose un algorithme de type explicite-implicite pour minimiser
J1 + τ J2 : soit une initialisation x0 ∈ Rn , pour k ≥ 0, on calcule
= xk − τ J10 (xk ),
( k+1/2
x
(3.54)
xk+1 = proxτ J2 (xk+1/2 ),
où on a appliqué une itération de l’algorithme du gradient à pas fixe à J1 suivi d’une
itération de l’algorithme proximal à J2 avec le même pas fixe τ pour les deux. Pour
le cas particulier (3.49), cet algorithme proximal est
k+1 k ∗ k
x = Sτ x − τ A (Ax − b) .
Vérifier que IK est une fonction convexe et que l’application proximale proxIK est
simplement l’opérateur de projection orthogonale sur K.
Exercice 3.2.5 Soit J une fonction convexe minorée par une fonction affine sur Rn .
Pour τ > 0, on définit la régularisation de Moreau-Yosida de J, notée Jτ , par
1
Jτ (x) = minn kx − yk2 + J(y).
y∈R 2τ
Montrer que Jτ est convexe, différentiable et que minx∈Rn J(x) = minx∈Rn Jτ (x) (indi-
cation : on pourra utiliser la forte convexité de la fonction dont la minimisation donne
Jτ pour montrer d’abord que l’application proximale est différentiable). Vérifier que
c’est-à-dire qu’une itération de l’algorithme proximal pour J est équivalent à une itération
de l’algorithme du gradient à pas fixe pour Jτ .
Exercice 3.2.6 Soit J1 une fonction fortement convexe différentiable sur Rn . Soit K
un convexe fermé de Rn et J2 = IK , la fonction indicatrice de K définie comme à
l’Exercice 3.2.4. Vérifier que l’algorithme explicite-implicite (3.54) n’est rien d’autre que
l’algorithme de gradient projeté de la Sous-section 3.4.1.
Exercice 3.2.7 Soit A une matrice n×n symétrique réelle, définie positive, et b ∈ Rn .
On définit la fonction quadratique J sur Rn par
1
J(x) = Ax · x − b · x.
2
Déterminer l’opérateur proximal proxJ .
Lemme 3.2.33 Soit J(x) une fonction convexe différentiable de Rn dans R qui
admet un unique point de minimum x∗ . Pour τ > 0 et x0 ∈ Rn on définit l’algorithme
proximal par la suite, indexée par k ≥ 0,
k+1 k 1 k 2
x = proxτ J (x ) = arg minx∈Rn kx − x k + τ J(x) .
2
dont on déduit que J(xk ) converge vers J(x∗ ) et que la suite (xk ) est bornée dans
0
Rn . Par conséquent, pour une sous-suite, xk converge vers une limite x∞ et, comme
J est continue et que son point de minimum est unique, il vient que x∞ = x∗ et que
toute la suite xk converge vers x∗ .
Si J est α-convexe, alors on peut améliorer (3.55) car x 7→ J(x)+kx−xk k2 /(2τ )
est (α + 1/τ )-convexe
1 1 1
J(x) + kx − xk k2 ≥ J(xk+1 ) + kxk+1 − xk k2 + (α + 1/τ )kx − xk+1 k2 .
2τ 2τ 2
On choisit encore x = x∗ et on combine cette inégalité avec la forte convexité (2.23)
de J à son point de minimum x∗ (où sa dérivée s’annule)
α k+1
J(xk+1 ) ≥ J(x∗ ) + kx − x∗ k 2 .
2
En sommant les deux dernières inégalités on obtient
d’où le résultat. 2
On imagine ainsi qu’elle peut être plus efficace qu’un algorithme du premier ordre,
puisqu’elle utilise, en plus du gradient, la Hessienne de la fonction, mais qu’elle
risque d’être plus compliquée à mettre en oeuvre puisqu’il faut justement calculer
ces dérivées secondes. Comme dans les sections précédentes on considère un problème
d’optimisation sans contrainte, sauf dans la dernière sous-section.
c’est-à-dire
−1
u = v − (F 0 (v)) F (v) + O kv − uk2 .
Rappelons que l’on ne calcule pas l’inverse de la matrice F 0 (un ) dans (3.56) mais que
l’on résout un système linéaire. Du point de vue de l’optimisation, la méthode de
Newton s’interprète de la manière suivante. Soit J une fonction de classe C 3 de RN
dans R, et soit u un minimum local de J. En notant F = J 0 , la méthode précédente
permet de résoudre la condition nécessaire d’optimalité J 0 (u) = 0. Plus précisément,
on peut aussi voir la méthode de Newton comme une méthode de minimisation. A
cause du développement de Taylor
1
J(w) = J(v) + J 0 (v) · (w − v) + J 00 (v)(w − v) · (w − v) + O kw − vk3 , (3.57)
2
on peut approcher J(w) au voisinage de v par une fonction quadratique. La méthode
de Newton consiste alors à minimiser cette approximation quadratique et à itérer.
Le minimum de la partie quadratique du terme de droite de (3.57) est donné par
94 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
Supposons que toutes les itérées jusqu’à un sont restées proches de u, au sens où
ku − un k ≤ δ, donc F 0 (un ) est inversible. Comme F (u) = 0, on déduit de (3.56)
−1
un+1 − u = un − u − (F 0 (un )) (F (un ) − F (u))
c’est-à-dire 2n
−1 n
n
ku − uk ≤ C 0
Cku − uk ≤ C −1 (C)2
S n+1 = S n + C n
où C n est une matrice de rang 1 qui dépend de un , un+1 , F (un ), F (un+1 ), choisie
de manière à ce que S n − (F 0 (un ))−1 converge vers 0. Pour plus de détails sur ces
méthodes de quasi-Newton nous renvoyons à [6] et [13]. •
min kAn u − bn k2 ,
u∈RN
c’est-à-dire que un+1 est une solution de l’équation normale, (∇F (un ))∗ ∇F (un )un =
(∇F (un ))∗ bn . Si l’on suppose que Ker∇F (un ) = {0}, alors la matrice (∇F (un ))∗ ∇F (un )
est inversible et on vérifie aisément que un+1 est donnée par
−1
un+1 = un − (∇F (un ))∗ ∇F (un ) (∇F (un ))∗ F (un )
Remarque 3.3.4 Si N = M et Ker∇F (un ) = {0}, alors ∇F (un ) est une matrice
inversible et la formule ci-dessus se simplifie en
−1
un+1 = un − ∇F (un ) F (un )
alors la fonction J n (w), définie par (3.59), est fortement convexe, donc admet un
unique point de minimum un+1 caractérisé par la condition d’optimalité du premier
ordre
(J n )0 (un+1 )(w) = 0 pour tout w ∈ V,
3.3. MÉTHODE DE NEWTON 97
L’équation (3.60) n’est rien d’autre qu’une formulation variationnelle sur V qu’on
peut résoudre à l’aide du lemme de Lax-Milgram grâce à l’hypothèse de coercivité
de J 00 (un ) (voir [1]). C’est typiquement de cette façon que sont résolues les équations
aux dérivées partielles non-linéaires.
Si les vecteurs (G01 (u), ..., G0M (u)) sont linéairement indépendants, la condition né-
cessaire d’optimalité du Théorème 2.5.6 est
M
X
J 0 (u) + λi G0i (u) = 0 , Gi (u) = 0 1 ≤ i ≤ M , (3.62)
i=1
la matrice F 0 (u, λ) est inversible. On remarque que (3.63) est l’inégalité stricte dans
la condition d’optimalité d’ordre 2 de la Proposition 2.5.13. Il est donc naturel de
faire l’hypothèse (3.63) qui permet d’utiliser l’algorithme de Newton. On peut ainsi
démontrer la convergence de cette méthode (voir [6]). Il est intéressant d’interpréter
cet algorithme comme une méthode de minimisation. On introduit le Lagrangien
L(v, µ) = J(v) + µ · G(v), ses dérivées par rapport à v, L0 et L00 , et on vérifie que
l’équation
−1
(un+1 , λn+1 ) = (un , λn ) − (F 0 (un , λn )) F (un , λn )
est la condition d’optimalité pour que un+1 soit un point de minimum du problème
quadratique à contraintes affines
avec
1
n
Q (w) = L(u , λ ) + L (u , λ ) · (w − u ) + L00 (un , λn )(w − un ) · (w − un )
n n 0 n n n
,
2
hu − u − µJ 0 (u) , v − ui ≥ 0 ∀ v ∈ K .
(3.68)
u = PK u − µJ 0 (u) .
(3.69)
Plus précisément, (3.69) est en fait équivalent à (3.67) et caractérise donc la solu-
tion u de (3.66). L’algorithme de gradient à pas fixe avec projection (ou plus
simplement de gradient projeté) est alors défini par l’itération
un+1 = PK un − µJ 0 (un ) ,
(3.70)
où µ > 0 est un pas de descente fixe. Autrement dit, (3.70) est une simple adaptation
de l’algorithme de gradient à pas fixe auquel on a ajouté une étape de projection
pour ne pas quitter l’ensemble admissible K.
Théorème 3.4.1 On suppose que J est α-convexe différentiable et que J 0 est Lip-
schitzien sur V (de constante L, voir (3.9)). Alors, si 0 < µ < 2α/L2 , l’algorithme
de gradient à pas fixe avec projection converge : quel que soit u0 ∈ K, la suite (un )
définie par (3.70) converge vers la solution u de (3.66).
Démonstration. La démonstration reprend celle du Théorème 3.1.6 où l’on a mon-
tré que l’application v 7→ v − µJ 0 (v) est strictement contractante lorsque 0 < µ <
2α/L2 , c’est-à-dire qu’il existe γ ∈]0, 1[ tel que
Puisque la projection
0
PK est faiblement contractante d’après (8.2), l’application
v 7→ PK v − µJ (v) est strictement contractante, ce qui prouve la convergence de
la suite (un ) définie par (3.70) vers la solution u de (3.66). 2
PN
Exercice 3.4.1 Soit V = RN et K = {x ∈ RN tel que i=1 xi = 1}. Expliciter
l’opérateur de projection orthogonale PK et interpréter dans ce cas la formule (3.69) en
terme de multiplicateur de Lagrange.
100 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
Le Théorème 3.4.1 montre que l’algorithme de gradient à pas fixe avec pro-
jection est applicable à une large classe de problèmes d’optimisation convexe avec
contraintes, dès lors que l’on sait calculer l’opérateur de projection orthogonale PK .
Mais malheureusement, il est très difficile, non seulement de connaitre explicitement
l’opérateur de projection PK , mais même de construire un algorithme simple et ef-
ficace pour l’évaluer, dans le cas d’un convexe fermé K quelconque dans V . Malgré
sa simplicité apparente l’algorithme du gradient projeté est donc souvent un leurre
du point de vue pratique par défaut d’une connaissance explicite de cet opérateur
de projection PK .
Donnons néanmoins deux exceptions importantes pour lesquelles on sait calculer
explicitement PK . Par souci de simplicité dans la présentation, on se limite au cas
de la dimension finie V = RN , mais ces deux exemples se généralisent sans difficulté
à des espaces de Hilbert de dimension infinie. On considère tout d’abord le cas d’un
sous-espace affine
K = {x ∈ RN Bx = c},
avec B matrice M × N de rang maximal, rang(B) = M ≤ N , et c ∈ RM . On trouve
alors facilement que
Remarque 3.4.2 Ce qui rend très attractif, mais très difficile à mettre en oeuvre,
l’algorithme du gradient projeté est qu’il s’agit d’un algorithme faisable, c’est-à-
dire qu’à chaque itération la solution approchée un appartient bien à l’ensemble
K qui définit les contraintes. A l’opposé de cette approche il existe une classe
d’algorithmes infaisables où la suite de solutions approchées un n’appartient pas
à K mais, en cas de convergence, leur limite u appartient bien à K. C’est évidem-
ment beaucoup plus facile à mettre en oeuvre mais, en général, avant la convergence
on n’a pas une solution, même de qualité médiocre, qui satisfait les contraintes. La
section suivante donne des exemples de tels algorithmes infaisables. •
3.4. ALGORITHMES DE TYPE GRADIENT (AVEC CONTRAINTES) 101
on déduit que (p − q) · F (u) ≥ 0 pour tout q ∈ (R+ )M , d’où on tire, pour tout réel
µ > 0,
(p − q) · p − p + µF (u) ≤ 0 ∀ q ∈ (R+ )M ,
p = PRM
+
(p + µF (u)) ∀ µ > 0 , (3.76)
PRM+
désignant la projection de RM sur (R+ )M .
Au vu de cette propriété et de la seconde inégalité dans (3.75), nous pouvons
introduire l’algorithme d’Uzawa : à partir d’un élément quelconque p0 ∈ (R+ )M ,
on construit les suites (un ) et (pn ) déterminées par les itérations
Théorème 3.4.3 On suppose que J est α-convexe différentiable, que F est convexe
et Lipschitzienne de V dans RM , c’est-à-dire qu’il existe une constante L telle que
et que l’ensemble K = {v ∈ V tel que F (v) ≤ 0} est non vide. On suppose que
les contraintes sont qualifiées pour la solution u de (3.73). Si 0 < µ < 2α/L2 ,
l’algorithme d’Uzawa converge : quel que soit l’élément initial p0 , la suite (un ) définie
par (3.77) converge vers la solution u du problème (3.73).
Démonstration. Tout d’abord, l’existence et l’unicité d’une solution u de (3.73)
découle du Théorème 2.3.8 puisque l’on minimise une fonction fortement convexe
sur un convexe fermé non vide K. Comme les contraintes sont qualifiées en u, le
Théorème de Kuhn et Tucker 2.6.4 affirme alors qu’il existe un multiplicateur de
Lagrange p ∈ (R+ )M tel que (u, p) est un point-selle du Lagrangien (3.74) sur
V × (R+ )M . De même, pn étant fixé, le problème de minimisation dans (3.77) a bien
une solution unique un . D’après l’Exercice 2.5.6, les inéquations d’Euler satisfaites
par u et un s’écrivent
soit
krn+1 k2 ≤ krn k2 + 2µrn · F (un ) − F (u) + µ2 kF (un ) − F (u)k2 .
Si 0 < µ < 2α/L2 , on peut trouver β > 0 tel que L2 µ2 − 2µα < −β, d’où
Ceci montre alors que la suite krn k2 est décroissante : le membre de droite de (3.82)
tend donc vers 0, ce qui entraîne que un tend vers u. 2
où µn > 0 est un pas de descente positif. Autrement dit, (3.85) recherche un point
selle en alternant un pas de minimisation en v et un pas de maximisation en q.
Théorème 3.4.4 On suppose que J est α-convexe différentiable, que F est convexe
différentiable de dérivée bornée par CF > 0. On suppose aussi que, pour tout p ∈
(R+ )M , la fonction v 7→ J 0 (v) + p · F 0 (v) est L-Lipschitzienne sur V . Alors il existe
un point-selle (u, p) du Lagrangien (3.74) sur V × (R+ )M et u est la solution unique
du problème (3.73). Si le pas de descente vérifie, pour C 2 = max(2L2 + CF2 , 2CF2 ),
αγ kun − uk2
µn = avec 0 < γ < 1,
C 2 kun − uk2 + kpn − pk2
alors, l’algorithme d’Arrow-Hurwicz converge : quel que soit les éléments initiaux
p0 , u0 , la suite (un ) définie par (3.85) converge vers la solution u du problème (3.73).
Remarque 3.4.5 Bien que le Théorème 3.4.4 prescrive une valeur variable de µn ,
en pratique on utilise une valeur fixe 0 < µ < α/C 2 car on ne connait pas le point
selle (u, p) et on ne peut donc pas calculer µn .
Démonstration. Les hypothèses impliquent qu’il existe un point-selle (u, p) du
Lagrangien (3.74) sur V × (R+ )M . Rappelons que d’après (3.76) on a p = PRM +
(p +
n
µF (u)) pour tout réel µ > 0. Dans ce qui suit on écrit µ au lieu de µ par souci de
simplification. Comme la projection PK est faiblement contractante d’après (8.2),
on déduit de la deuxième équation de (3.85) que
∂L n n ∂L
kun+1 − uk2 = kun − uk2 − 2µh (u , p ), (un − u)i + µ2 k (un , pn )k2 . (3.87)
∂v ∂v
Or, comme J 0 (u) + p · F 0 (u) = 0, on a
∂L n n
(u , p ) = J 0 (un ) + p · F 0 (un ) − (J 0 (u) + p · F 0 (u)) + (pn − p) · F 0 (un ),
∂v
3.4. ALGORITHMES DE TYPE GRADIENT (AVEC CONTRAINTES) 105
donc
∂L n n 2
(u , p )k ≤ 2 kJ 0 (un ) + p · F 0 (un ) − (J 0 (u) + p · F 0 (u))k2 + k(pn − p) · F 0 (un )k2
k
∂v
≤ 2 L2 kun − uk2 + CF2 kpn − pk2 ,
car J 0 (v)+p·F 0 (v) est L-Lipschitzienne et F 0 est borné. Par conséquent, en sommant
(3.86) et (3.87) et en appliquant le Lemme 3.4.6 pour v = un et q = pn , on trouve
avec C 2 = max(2L2 + CF2 , 2CF2 ) et en ayant remplacé µ par sa valeur µn donnée dans
l’énoncé. On en déduit que la suite kun − uk2 + kpn − pk2 est décroissante positive
donc convergente, et par ailleurs que le rapport kun − uk4 / (kun − uk2 + kpn − pk2 )
tend vers zéro. Si la suite kun − uk2 + kpn − pk2 tend vers zéro, alors la suite (un , pn )
converge vers le point selle (u, p). Sinon, alors c’est la suite kun − uk4 qui tend vers
zéro, c’est-à-dire que la suite un converge vers u (mais on ne peut rien dire pour la
suite pn ). 2
Lemme 3.4.6 Soit le Lagrangien L(v, q) = J(v) + q · F (v), défini sur V × (R+ )M ,
avec J et F convexes différentiables et J α-convexe. On suppose que (u, p) est un
point selle du Lagrangien. Alors, pour tout (v, q) ∈ V × (R+ )M ,
α ∂L ∂L
ku − vk2 ≤ h (v, q), (v − u)i − (v, q) · (q − p) .
2 ∂v ∂q
Démonstration. Pour q ≥ 0 le Lagrangien est α-convexe
∂L α
L(u, q) ≥ L(v, q) + h (v, q), (u − v)i + ku − vk2 .
∂v 2
Comme (u, p) est un point selle, L(u, q) ≤ L(u, p) et on en déduit
α ∂L
ku − vk2 ≤ h (v, q), (v − u)i + L(u, p) − L(v, q) . (3.88)
2 ∂v
Or L(v, q) = L(v, p) + (q − p) · F (v) et ∂L
∂q
(v, q) = F (v). La deuxième inégalité du
point selle, L(u, p) ≤ L(v, p) permet enfin de conclure à partir de (3.88). 2
stratégie de résolution qui doit, pour sa partie numérique, utiliser l’un des algo-
rithmes de la Sous-section 3.1 pour l’optimisation sans contrainte. Cette résolution
numérique peut d’ailleurs soulever des difficultés, car le problème “pénalisé” est sou-
vent “mal conditionné”. Néanmoins l’avantage de cette méthode de pénalisation est
sa facilité de mise en oeuvre, même si elle est parfois décevante pour la qualité de
ses résultats numériques.
Nous nous plaçons pour simplifier dans le cas où V = RN , et nous considérons
de nouveau le problème de minimisation convexe
inf J(v) , (3.89)
F (v)≤0
est non vide. En notant u l’unique solution de (3.89) et, pour ε > 0, uε l’unique
solution de (3.90), on a alors
lim uε = u .
ε→0
Exercice 3.4.3 En plus des hypothèses de la Proposition 3.4.7, on suppose que les
fonctions J et F1 , . . . , FM sont continûment différentiables. On note de nouveau I(u)
l’ensemble des contraintes actives en u, et on suppose que les contraintes sont qualifiées
en u au sens de la Définition 2.5.15. Enfin, on suppose que les vecteurs Fi0 (u) i∈I(u)
où le deuxième terme est appelée une fonction “barrière” dont le choix n’est pas
unique. De manière générale, une fonction barrière est une fonction définie sur R−
qui tend vers +∞ quand son argument tend vers zéro. Le but de cette barrière est
d’empêcher une suite de solutions approchées de (3.93) de s’approcher de la limite
Fi (v) = 0. Ainsi les contraintes F (v) < 0 ne sont jamais actives et peuvent être
ignorées en pratique. Ainsi, si on utilise un algorithme de gradient pour le problème
sans contrainte (3.93) avec une initialisation qui vérifie strictement les contraintes
et un pas de descente suffisamment petit, alors la suite des itérées va rester dans le
domaine strictement faisable F (v) < 0. Pour que l’effet de la barrière soit de plus en
plus négligeable, on fait progressivement tendre le paramètre ε > 0 vers zéro (notez
la différence fondamentale avec le coefficient 1/ε en pénalisation extérieure). Nous
étudierons un exemple de pénalisation intérieure à la Sous-section 4.2.3. •
108 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
λ∗ = λn + µF (v n ). (3.98)
De (3.98) on tire deux informations. D’une part, si la suite λn tend vers λ∗ , il n’est
pas nécessaire de faire tendre le coefficient de pénalisation µ vers l’infini pour avoir
la satisfaction progressive de la contrainte F (v) = 0. D’autre part, (3.98) suggère
une formule de mise à jour du multiplicateur de Lagrange
λn+1 = λn + µF (v n ).
Cette formule est très semblable à celle de l’algorithme d’Uzawa, voir (3.85), sauf
qu’ici le pas µ n’est pas petit. Au total l’algorithme du Lagrangien augmenté
est le suivant : on choisit µ > 0, on initialise λ0 ∈ RM et on construit deux suites v n
et λn par
Laug (v n , λn , µ) = min Laug (v, λn , µ),
v∈RN
λn+1 = λn + µF (v n ).
Par ailleurs, de temps en temps, et un nombre limité de fois (voir le chapitre 17
de [26]), on augmente la valeur du coefficient de pénalisation µ. Mais, comme le
rappelle le résultat ci-dessous, il n’est pas nécessaire de faire tendre µ vers l’infini
pour converger.
Alors, il existe µ0 > 0 tel que, pour tout µ ≥ µ0 , v ∗ est un point de minimum local
de v 7→ Laug (v, λ∗ , µ) (sans contrainte).
Remarque 3.4.10 L’intérêt du Lemme 3.4.9 est de montrer que, si l’on connait la
valeur du multiplicateur de Lagrange optimal, alors la minimisation, sans contrainte
et avec une pénalisation finie de la contrainte, du Lagrangien augmenté conduit à
une solution du problème d’origine (3.94). Cela justifie, en quelque sorte, le fait qu’il
n’est pas nécessaire de faire tendre vers l’infini le coefficient de pénalisation µ pour
que l’algorithme du Lagrangien augmenté converge, du moins si on a une bonne
estimation du multiplicateur de Lagrange. •
Démonstration. Si v ∗ un point de minimum local de (3.94) et que les contraintes
sont qualifiées, alors la condition d’optimalité du premier ordre du Théorème 2.5.6
donne
J 0 (v ∗ ) + λ∗ · F 0 (v ∗ ) = 0 et F (v ∗ ) = 0,
110 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
ce qui implique que L0aug (v ∗ , λ∗ , µ) = 0 (ici et dans toute la suite de cette preuve le
signe 0 indique une dérivée en v). Calculons la dérivée seconde en v du Lagrangien
augmenté
avec
L00 (v ∗ , λ∗ )(w, w) = J 00 (v ∗ )(w, w) + λ∗ · F 00 (v ∗ )(w, w)
et où il n’y a pas de dérivée seconde dans le terme de pénalisation car F (v ∗ ) = 0.
L’hypothèse (3.99) montre que le premier terme à droite dans (3.100) est une forme
quadratique définie positive sur KerF 0 (v ∗ ), tandis que le second terme est une forme
quadratique définie positive sur le sous-espace de dimension M engendré par les
colonnes de F 0 (v ∗ ) (qui est aussi l’orthogonal de KerF 0 (v ∗ )). Notons que sur ce sous-
espace de dimension finie L00aug (v ∗ , λ∗ , µ)(w, w) peut prendre des valeurs négatives
mais qu’on peut les minorer par −µ0 kwk2 où −µ0 est la plus petite valeur propre
de la matrice M × M qui représente cette forme quadratique sur ce sous-espace. Il
suffit alors de choisir µ > µ0 pour que la somme des deux termes soit strictement
positif pour w 6= 0. Au total, la dérivée seconde du Lagrangien augmenté est définie
positive en v ∗ qui est un point critique. C’est donc un point de minimum local. 2
où J(v) et (F1 (v), ..., FM (v)) = F (v) sont des fonctions régulières de RN dans R.
Les remarques qui suivent s’appliquent de la même manière aux problèmes avec
contraintes d’inégalité, moyennant des modifications évidentes.
Si l’application directe des algorithmes d’optimisation proposés dans ce chapitre
est trop compliquée ou coûteuse pour (3.101), une stratégie courante consiste à rem-
placer ce dernier par une succession de problèmes approchés, obtenus, par exemple,
par développement de Taylor des fonctions J et F autour de la solution approchée
précédente. L’idée sous-jacente est qu’il est plus facile de résoudre le problème appro-
ché que le problème exact. Ces approximations n’ayant en général qu’un caractère
local, il faut itérer cette stratégie en faisant un nouveau développement de Taylor
au point de minimum obtenu sur le précédent problème approché.
Une méthode d’approximation successive peut donc se présenter comme suit.
Etant donnée une initialisation v 0 ∈ V = RN , on construit une suite de solutions
approchées v n ∈ V , n ≥ 1, où v n est le point de minimum du problème
pour que (3.102) admette une solution unique v n , mais aussi pour que le calcul de
v n soit facile.
Expliquons le principe de cette méthode dans un cas très simple, en l’absence
de contraintes d’égalité, c’est-à-dire F ≡ 0. On choisit Fn−1 ≡ 0 et
1
Jn−1 (v) = J(v n−1 ) + J 0 (v n−1 ) · (v − v n−1 ) + kv − v n−1 k2 , (3.103)
2µ
avec un paramètre µ > 0. Notons que les fonctions J et Jn−1 sont tangentes en v n−1
et donc que Jn−1 est bien une approximation de J en ce point.
1 n
J 0 (v n−1 ) + (v − v n−1 ) = 0,
µ
qui n’est rien d’autre qu’un programme linéaire, comme étudié dans la Section 4.2 et
pour lequel on dispose d’algorithmes extrêmement efficaces. Une difficulté immédiate
dans la résolution de (3.104) est que sa valeur minimum peut être −∞ et qu’il n’y a
pas de solution optimale. Notons que, sous une condition de qualification standard,
(F10 (v n−1 ), ..., FM
0
(v n−1 )) famille libre de RN , l’ensemble admissible de (3.104) n’est
112 CHAPITRE 3. ALGORITHMES D’OPTIMISATION
pas vide. C’est pourquoi, en pratique, cette méthode s’accompagne d’une contrainte
supplémentaire, dite de région de confiance, qui prend la forme
kv − v n−1 k ≤ δ, (3.105)
où δ > 0 est un paramètre qui définit la taille du voisinage de v n−1 dans lequel
(3.104) est une bonne approximation de (3.101). La norme Pdans (3.105) peut être
N
soit la norme kvk∞ = max1≤i≤N |vi |, soit la norme kvk1 = i=1 |vi |, ce qui dans les
deux cas préserve le fait que le problème approché est un programme linéaire. Ce
dernier a alors nécessairement au moins une solution optimale puisque l’ensemble
admissible est désormais borné.
le point de minimum de (3.101), ce qui entraine que (3.106) admet au moins une
solution optimale. Par contre, la matrice J 00 (v n−1 ) n’a aucune raison d’être positive
en général. Néanmoins, si v n−1 n’est pas un point de minimum, il peut être nécessaire
de recourir à nouveau à une contrainte de région de confiance du type de (3.105).
Pour plus de détails nous renvoyons à [26].
où J10 (v) est une matrice N2 × N et J20 (v2 ) est un vecteur de RN2 puisque J2 est une
fonction de RN2 dans R. En posant v2 = J1 (v), h2 = J10 (v)h + o(h) et en composant
les deux dérivées on obtient
J2 ◦ J1 (v + h) = J2 ◦ J1 (v) + J20 ◦ J1 (v) · J10 (v)h + o(h) pour tout h ∈ RN .
3.6. MODÉLISATION, STRUCTURES ET ALGORITHMES SPÉCIFIQUES 115
Exercice 3.6.1 P Montrer que pour calculer la dérivée directionnelle J 0 (v) · h la formule
(3.110) nécessite m i=1 Ni Ni+1 multiplications, tandis que pour calculer le gradient com-
plet J 0 (v) la formule (3.111) nécessite exactement le même nombre de multiplications.
est déterminé par les paramètres v. Cette contrainte d’égalité est typiquement un
modèle qui permet de calculer y en fonction de v. Par exemple, y est la solution
d’une équation différentielle où v est un paramètre (un coefficient, un terme source,
etc.). Nous rencontrerons cette situation en théorie du contrôle optimal. En tout
état de cause, il faut comprendre que la résolution du modèle pour trouver y en
fonction de v est une opération non explicite et coûteuse.
Etant donné une fonction objectif J(v, y), différentiable de RN × V dans R, et
une fonction contrainte F (v, y), différentiable de RN × V dans V , on considère le
problème d’optimisation
Bien sûr, il est possible d’avoir d’autres contraintes mais nous nous concentrons ici
sur la contrainte de modèle qui relie y à v. Nous faisons l’hypothèse suivante sur ce
modèle.
Hypothèse (H). Pour tout v ∈ RN , il existe une unique solution y(v) ∈ V de
la contrainte F (v, y) = 0, cette solution v 7→ y(v) est différentiable de RN dans V
et l’opérateur linéarisé z 7→ ∂F
∂y
(v, y(v))(z) est continu, inversible de V dans V et
d’inverse continu de V dans V .
Sous cette hypothèse, le problème (3.112) est équivalent au problème sans
contrainte
inf J(v)
e = J(v, y(v)). (3.113)
v∈RN
Lemme 3.6.3 Sous l’hypothèse (H) la fonction Je est différentiable sur RN et, pour
tout h ∈ RN , on a
0 ∂J ∂J 0
J (v) · h =
e (v, y(v)) · h + (v, y(v)) , y (v) · h ,
∂v ∂y
un simple système linéaire. Cet exemple est tout à fait représentatif car il correspond,
par exemple, à la discrétisation d’une équation différentielle linéaire dont la solution
y dépend d’une variable de contrôle v. La difficulté avec le Lemme 3.6.3 est qu’il
ne donne pas une formule pour le gradient Je0 (v) mais seulement pour la dérivée
directionnelle dans la direction h. Pour obtenir le gradient il faut calculer chacune
de ses composantes en résolvant (3.114) pour h successivement égal à chacun des
vecteurs ei de la base canonique de RN , ce qui revient à résoudre N fois ce système
linéaire, opération trop coûteuse dans de nombreux cas pratiques où M et N sont
grands. •
Le Lagrangien est affine en p et sa dérivée partielle par rapport à p n’est rien d’autre
que la contrainte F (v, y). Cela motive l’étude d’une autre dérivée partielle du La-
grangien, celle par rapport à y, qui vaut, pour tout z ∈ V ,
∂L ∂J ∂F
h (v, y, p), zi = h (v, y), zi + hp, (v, y)(z)i. (3.115)
∂y ∂y ∂y
Lemme 3.6.5 Le problème adjoint est défini par : trouver p ∈ V solution de
∗
∂F ∂J
(v, y(v)) p = − (v, y(v)). (3.116)
∂y ∂y
Sous l’hypothèse (H), ce problème admet une unique solution p ≡ p(v) ∈ V , appelé
état adjoint. Autrement dit, l’état adjoint est l’unique p(v) ∈ V tel que
∂L
h (v, y(v), p(v)), zi = 0 pour tout z ∈ V.
∂y
∂F
Démonstration. Comme ∂y
(v, y(v))
est un opérateur linéaire continu de V dans
∗
∂F
V , il admet un adjoint, noté ∂y (v, y(v)) . D’après l’hypothèse (H) il est de plus
inversible d’inverse continu, et c’est vrai aussi pour son adjoint. Puisque J est diffé-
rentiable, le second membre de (3.116) appartient à V et on en déduit donc l’exis-
tence et l’unicité de l’état adjoint p(v) solution de (3.116). On remplace alors y par
y(v) et p par p(v) dans (3.115) pour obtenir que la dérivée partielle du Lagrangien
par rapport à y s’annule au point (v, y(y), p(v)). 2
∂L
Je0 (v) · h = (v, y(v), p(v)) · h,
∂v
ce qui donne la formule (3.117). 2
3.6.3 Décomposition-coordination
On considère un problème d’optimisation sous contraintes d’égalité
n
X n
X
J(v) = Ji (vi ) et F (v) = Fi (vi ),
i=1 i=1
De nombreux problèmes pratiques peuvent se mettre sous cette forme et, pour fixer
les idées, nous allons considérer le cas de l’optimisation du fonctionnement d’une
3.6. MODÉLISATION, STRUCTURES ET ALGORITHMES SPÉCIFIQUES 119
entreprise composée de n divisions ou unités qui ont chacune leur propre variable
de décision vi et qui doivent minimiser leur coût de fonctionnement Ji (vi ) avec
leur contrainte de production Fi (vi ), indépendants de celui des autres unités j 6= i,
mais en partageant le niveau global de contraintes imposées, comme un budget, des
ressources humaines ou des niveaux de production qui sont agrégées au niveau de
l’entreprise entière.
Si on note p le multiplicateur de Lagrange associé à la contrainte F (v) = F0 ,
pour p ∈ RM et v ∈ RN , on introduit le Lagrangien
n
X
L(v, p) = J(v) + p · (F (v) − F0 ) = (Ji (yi ) + p · Fi (v)) − p · F0 .
i=1
L’algorithme d’Uzawa (3.77) appliqué au problème (3.118) consiste, pour une initia-
lisation p0 ∈ RM , à construire les suites (uk ) et (pk ) déterminées par les itérations,
pour k ≥ 0,
L(uk , pk ) = inf L(v, pk ) ,
v∈RN
(3.120)
k+1 k k
p = p + µ(F (u ) − F0 ) ,
où µ > 0 est un pas positif fixé. L’hypothèse de décomposition des fonctions J et
F rend en fait la minimisation du Lagrangien, à p fixé, particulièrement facile car il
suffit de construire uk = (uk1 , ..., ukn ), où chaque composante uki est la solution d’un
problème de minimisation
Chacun des problèmes (3.121) est de plus petite taille que (3.118) et indépendant
des autres composantes j 6= i du problème. Dans ce contexte, l’algorithme d’Uzawa
est appelé algorithme de décomposition par les prix.
Une interprétation usuelle de cet algorithme est la suivante. Supposons, pour
simplifier, que les fonctions Ji et Fi sont positives : elles correspondent au coût
de production et au niveau de production, respectivement, de l’unité i. A chaque
itération k, un organe central de l’entreprise fixe un prix −pk pour la valeur de
la production (le vecteur pk a toutes ses composantes négatives). Du coup, chaque
unité veut maximiser son gain, c’est-à-dire minimiser son coût moins ses revenus
correspondant au prix −pk de ce qu’elle produit. Autrement dit, chaque unité résout
le problème (3.121) sans se soucier des autres unités. A la fin de l’iteration k, la
direction de l’entreprise regarde si elle a atteint le niveau de production souhaité F0 .
Si c’est le cas, la procédure a convergé. Sinon, on corrige les prix en les augmentant
si le niveau de production est trop bas, F (uk ) < F0 , ou en les diminuant si la
production est trop forte, F (uk ) > F0 (on rappelle que le prix est −pk ).
converge, c’est-à-dire que l’algorithme de décomposition par les prix converge dans
ce cas particulier. On suppose ici que F (v) est affine, ce qui permet de remplacer la
contrainte d’égalité par deux contraintes (opposées) d’inégalité, qui sont toutes les
deux convexes. •
Une autre approche de décomposition applicable au problème (3.119) est la dé-
composition par les quantités. Expliquons-en le principe en reprenant l’exemple
de l’entreprise composée de n unités. Il s’agit toujours d’un algorithme itératif et à
chaque itération k la direction de l’entreprise confie à chaque unité i le soin de pro-
duire une quantité Fik ∈ RM , en ayant pris soin de choisir ses quantités qui vérifient
la contrainte globale
Xn
Fik = F0 .
i=1
Chaque unité doit donc résoudre un problème d’optimisation sous contrainte
inf Ji (vi ), (3.122)
vi ∈Rni , Fi (vi )=Fik
dont on suppose qu’il admet une solution unique uki ∈ Rni . On suppose aussi qu’il
existe un unique multiplicateur de Lagrange pki ∈ RM pour la contrainte dans
(3.122). On note Gi (Fik ) la valeur du minimum dans (3.122), qui dépend du ni-
veau de contraintes Fik ∈ RM . On sait par le Lemme 2.5.12 que, si la fonction Gi
est différentiable, alors on a
pki = −∇Gi (Fik ).
Autrement dit, la sensibilité de la valeur du minimum par rapport au niveau de
contraintes est égale à l’opposé du multiplicateur de Lagrange (qui s’interprète clas-
siquement comme l’opposé d’un prix marginal). Ainsi, si deux unités i 6= j ont
des multiplicateurs de Lagrange différents pki 6= pkj , alors la direction de l’entreprise
peut changer l’allocation de production à ces deux unités, en Fik − µ(pkj − pki ) et
Fjk + µ(pkj − pki ), pour µ > 0 petit, ce qui conduit à une diminution de la somme des
minima dans (3.122) pour i et j qui vaut au premier ordre −µ|pkj − pki |2 < 0. On
en conclut que tant que les multiplicateurs de Lagrange des unités sont différents,
l’entreprise peut améliorer son minimum en changeant l’allocation des quantités à
produire. Concrètement, pour Ppasser à l’itération (k + 1), la direction de l’entreprise
forme le prix moyen pk+1 = ni=1 pki /n et propose la nouvelle allocation
Fik+1 = Fik + µ(pki − pk+1 ).
On vérifie aisément que la contrainte globale ni=1 Fik+1 = F0 est toujours vérifiée
P
et que, pour µ > 0 petit, le gain dans le minimum de (3.119) est égal au premier
ordre à n
X
−µ |pki − pk+1 |2 .
i=1
On améliore donc systématiquement le minimum avec cette méthode de décompo-
sition par les quantités qui, par ailleurs, vérifie toujours la contrainte F (uk ) = F0
au cours des itérations, contrairement à la méthode de décomposition par les prix
(voir [10] pour plus de détails).
Chapitre 4
PROGRAMMATION LINÉAIRE
4.1 Introduction
Ce chapitre est consacré à la programmation linéaire qui permet de résoudre
efficacement les problèmes d’optimisation où les contraintes et le critère s’expriment
linéairement en fonction des variables (voir l’Exemple 1.2.1). Ce type de problème est
extrêmement fréquent dans le domaine de la recherche opérationnelle (ou plus
simplement RO). La RO est l’ensemble des méthodes scientifiques qui permettent de
modéliser et analyser des situations complexes afin de prendre des décisions, sinon
optimales, du moins les plus efficaces possibles. C’est un domaine au carrefour des
mathématiques, de l’informatique et de l’ingénierie (au sens très large) où le mot
“opérationnel” fait référence à la planification des opérations militaires car la RO
est née pendant la seconde guerre mondiale. Depuis elle s’est très largement civilisée
et elle est employée dans de très nombreuses applications pour l’industrie et les
services, en économie et sciences de la gestion
La RO ne se limite pas du tout à la programmation linéaire et, bien au contraire,
utilise des outils très divers, issus de plusieurs champs scientifiques : optimisa-
tion (continue et combinatoire), probabilités (algorithmes stochastiques), théorie
des jeux, mathématiques discrètes, théorie des graphes, théorie de la complexité
(informatique), et programmation par contraintes. Nous n’étudierons pas tous ces
aspects de la RO, pas plus que nous ne toucherons aux questions de modélisation
ou à la conception d’“heuristiques", lorsque des méthodes rigoureuses font défaut,
ce qui représente souvent une part importante de la pratique en RO... Nous nous
contentons ici de donner un éclairage très partiel sur les apports de l’optimisation,
à travers la programmation linéaire, dans ce vaste domaine. Pour une introduction
à la partie mathématisée de la RO nous renvoyons le lecteur vers l’ouvrage [5], issu
d’un cours de troisième année à l’Ecole Polytechnique, ou bien vers [27].
La Section 4.2 est une présentation classique de la programmation linéaire. Au
delà des (nombreux) cas où le modèle d’optimisation (fonction objectif et contraintes)
est linéaire, rappelons que la programmation linéaire est aussi une brique de base
dans de nombreux algorithmes d’optimisation qui linéarisent les problèmes à ré-
soudre. La Section 4.3 est plus spécifique à la RO car elle étudie quelques exemples
où l’on recherche des solutions à valeurs entières, et non pas seulement réelles, de
121
122 CHAPITRE 4. PROGRAMMATION LINÉAIRE
programmes linéaires. Lorsqu’on impose cette contrainte additionnelle (et très forte)
de solutions entières, on parle alors de programmation linéaire en nombres entiers,
voire d’optimisation combinatoire car il serait possible en théorie, mais illusoire en
pratique, d’énumérer tous les candidats possibles pour trouver la solution optimale.
Là aussi nous ne ferons qu’effleurer cette vaste question et nous renvoyons à [5] pour
plus de détails.
inf c · x.
x∈Rn tel que Ax≥b, A0 x=b0
peut se mettre sous la forme standard (4.1) quitte à changer la taille des données.
En effet, remarquons tout d’abord que les contraintes d’égalité A0 x = b0 sont évi-
demment équivalentes aux contraintes d’inégalité A0 x ≤ b0 et A0 x ≥ b0 . On peut
donc se restreindre au cas suivant (qui ne contient que des contraintes d’inégalité)
inf c · x. (4.2)
x∈Rn tel que Ax≥b
inf c · x. (4.3)
(x,λ)∈R(n+m) tel que Ax=b+λ, λ≥0
x3
x2
x1
qui est bien sous forme standard (mais avec plus de variables). Il n’y a donc aucune
perte de généralité à étudier le programme linéaire standard (4.1).
Nous avons déjà donné une motivation concrète de la programmation linéaire au
début du Chapitre 1 (voir l’Exemple 1.2.1). Considérons pour l’instant un exemple
simple qui va nous permettre de comprendre quelques aspects essentiels d’un pro-
gramme linéaire
min x1 + 4x2 + 2x3 . (4.5)
x1 ≥0,x2 ≥0,x3 ≥0
2x1 +x2 +3x3 =6
Sur la Figure 4.1 nous avons tracé l’ensemble des (x1 , x2 , x3 ) qui vérifient les contrain-
tes : c’est un triangle plan T . C’est un fermé compact de R3 , donc la fonction continue
x1 +4x2 +2x3 y atteint son minimum que l’on note M . Pour déterminer ce minimum
on peut considérer la famille de plans parallèles x1 + 4x2 + 2x3 = v paramétrée par
v. En augmentant la valeur de v à partir de −∞, on “balaie” l’espace R3 jusqu’à
atteindre le triangle T , et le minimum M est obtenu lorsque le plan “touche” ce
triangle. Autrement dit, tout point de minimum de (4.5) est sur le bord du triangle
T . Une autre façon de le voir est de dire que la fonction x1 + 4x2 + 2x3 a un gradient
non nul dans T donc ses extréma se trouvent sur le bord de T . Pour l’exemple (4.5)
le point de minimum (unique) est le sommet (3, 0, 0) de T . Nous verrons qu’il s’agit
d’un fait général : un point de minimum (s’il existe) peut toujours se trouver en un
des sommets de l’ensemble géométrique des vecteurs x qui vérifient les contraintes.
Il “suffit” alors d’énumérer tous les sommets afin de trouver le minimum : c’est
précisément ce que fait (de manière intelligente) l’algorithme du simplexe que nous
verrons dans la prochaine sous-section.
Pour établir cette propriété en toute généralité pour le programme linéaire stan-
dard (4.1), nous avons besoin de quelques définitions qui permettent de préciser le
124 CHAPITRE 4. PROGRAMMATION LINÉAIRE
vocabulaire.
Définition 4.2.1 L’ensemble Xad des vecteurs de Rn qui satisfont les contraintes
de (4.1), c’est-à-dire
Xad = {x ∈ Rn tel que Ax = b, x ≥ 0} ,
est appelé ensemble des solutions admissibles. On appelle sommet ou point ex-
trémal de Xad tout point x ∈ Xad qui ne peut pas se décomposer en une combinai-
son convexe (non triviale) de deux autres points de Xad , c’est-à-dire que, s’il existe
y, z ∈ Xad et θ ∈]0, 1[ tels que x = θy + (1 − θ)z, alors y = z = x.
Lemme 4.2.3 Il existe au moins une solution optimale (ou point de minimum) du
programme linéaire standard (4.1) si et seulement si la valeur du minimum est finie
−∞ < inf c · x < +∞.
x∈Rn tel que Ax=b, x≥0
Démonstration. Soit (xk )k≥1 une suite minimisante de (4.1). On introduit la ma-
trice A définie par ∗
c
A= .
A
La suite Axk appartient au cône suivant
( n )
X
C= xi Ai avec xi ≥ 0 ,
i=1
Ax = BxB + N xN .
Lemme 4.2.5 Les sommets du polyèdre Xad sont exactement les solutions basiques.
Démonstration. Si x ∈ Xad est une solution basique, dans une certaine base de
Rn on a xP = (x1 , ..., xm , 0, ..., 0), A = (B, N ) avec B = (b1 , ..., bm ), une base de Rm
telle que m i=1 xi bi = b. Supposons qu’il existe 0 < θ < 1 et y, z ∈ Xad tels que
x = θy + (1 − θ)z. Nécessairement, les n − m dernières P composantes P de y et z sont
nulles et, comme y et z appartiennent à Xad , on a m y
i=1 i ib = b et m
i=1 zi bi = b.
Par unicité de la décomposition dans une base, on en déduit que x = y = z, et donc
x est un sommet de Xad .
Réciproquement, si x est un sommet de Xad , on note k le nombre Pk de ses com-
posantes non nulles, et après un éventuel réarrangement on a b = i=1 xi ai où les
(ai ) sont les colonnes de A. Pour montrer que x est une solution basique il suffit
de prouver que la famille (a1 , ..., ak ) est libre dans Rm (on obtient une base B en
Pk cette famille). Supposons que ce ne soit pas le cas : il existe alors y 6= 0
complétant
tel que i=1 yi ai = 0 et (yk+1 , ..., yn ) = 0. Comme les composantes (x1 , ..., xk ) sont
strictement positives, il existe > 0 (petit) tel que (x + y) ∈ Xad et (x − y) ∈ Xad .
Le fait que x = (x + y)/2 + (x − y)/2 contredit le caractère extrémal de x, donc x
est une solution basique. 2
Le résultat fondamental suivant nous dit qu’il est suffisant de chercher une
solution optimale parmi les sommets du polyèdre Xad .
Proposition 4.2.6 S’il existe une solution optimale du programme linéaire stan-
dard (4.1), alors il existe une solution optimale basique.
k
X
b= x i ai ,
i=1
126 CHAPITRE 4. PROGRAMMATION LINÉAIRE
où les (ai ) sont les colonnes de A. Si la famille (a1 , ..., ak ) est libre dans Rm , alors x
est une solution optimale basique. Si (a1 , ..., ak ) est lié, alors il existe y 6= 0 tel que
k
X
yi ai = 0 et (yk+1 , ..., yn ) = 0.
i=1
Comme les composantes (x1 , ..., xk ) sont strictement positives, il existe > 0 tel que
(x ± y) ∈ Xad . Comme x est un point de minimum, on a nécessairement
c · x ≤ c · (x ± y),
max x1 + 2x2
x1 ≥0,x2 ≥0
Exercice 4.2.2 Montrer que l’on peut choisir la matrice A de taille m×n et le vecteur
b ∈ Rm de telle façon que Xad soit le cube unité [0, 1]n−m dans le sous-espace affine de
dimension n − m défini par Ax = b. En déduire que le nombre de sommets de Xad est
alors 2n−m .
une solution optimale (ce qui est garanti si le programme linéaire admet effective-
ment une solution optimale). L’algorithme du simplexe ne se contente pas d’énumé-
rer tous les sommets, il décroît la valeur de la fonction c · x en passant d’un sommet
au suivant.
On considère le programme linéaire standard (4.1). Rappelons qu’un sommet
(ou solution basique) de l’ensemble des solutions admissibles Xad est caractérisé par
une base B (m colonnes libres de A). Après permutation de ses colonnes, on peut
écrire
A = (B, N ) et x = (xB , xN ),
de sorte qu’on a Ax = BxB + N xN . Toute solution admissible peut s’écrire xB =
B −1 (b − N xN ) ≥ 0 et xN ≥ 0. Le sommet associé à B est défini (s’il existe) par
xN = 0 et xB = B −1 b ≥ 0. Si on décompose aussi c = (cB , cN ) dans cette base, alors
on peut comparer le coût d’une solution admissible quelconque x avec celui de la
solution basique x
c · x − c · x = cB · B −1 (b − N xN ) + cN · xN − cB · B −1 b = (cN − N ∗ (B −1 )∗ cB ) · xN . (4.6)
Proposition 4.2.8 Supposons que la solution basique associée à B est non dégéné-
rée, c’est-à-dire que B −1 b > 0. Une condition nécessaire et suffisante pour que cette
solution basique associée à B soit optimale est que
c̃N = cN − N ∗ (B −1 )∗ cB ≥ 0. (4.7)
c · x − c · x = c̃N · xN ≥ 0 ,
puisque xN ≥ 0. Donc la condition (4.7) est suffisante pour que x soit optimal.
Réciproquement, supposons qu’il existe une composante i de c̃N qui soit strictement
négative, (c̃N · ei ) < 0. Pour > 0 on définit alors un vecteur x() par xN () = ei et
xB () = B −1 (b − N xN ()). Par construction Ax() = b et, comme B −1 b > 0, pour
des valeurs suffisamment petites de on a x() ≥ 0, donc x() ∈ Xad . D’autre part,
x(0) = x et, comme > 0, on a
ce qui montre que x n’est pas optimal. Donc la condition (4.7) est nécessaire. 2
Algorithme du simplexe
— Initialisation (phase I) : on cherche une base initiale B 0 telle que la solution
basique associée x0 soit admissible
(B 0 )−1 b
0
x = ≥ 0.
0
xk () = (xkB (), xkN ()) avec xkN () = ei , xkB () = (B k )−1 (b − ai ).
— Soit on peut choisir > 0 aussi grand que l’on veut avec xk () ∈ Xad .
Dans ce cas, le minimum du programme linéaire est −∞.
— Soit il existe une valeur maximale k ≥ 0 et un indice j tels que la j-ème
composante de xk (k ) s’annule. On obtient ainsi une nouvelle solution
admissible basique
xk+1 = xk (k ),
correspondant à une nouvelle base B k+1 déduite de B k en remplaçant sa
j-ème colonne par la colonne ai . La solution admissible xk+1 a un coût
inférieur ou égal à celui de xk .
Il reste un certain nombres de points pratiques à préciser dans l’algorithme du
simplexe. Nous les passons rapidement en revue.
Dégénérescence et cyclage
On a toujours c·xk+1 ≤ c·xk , mais il peut y avoir égalité si la solution admissible
basique xk est dégénérée, auquel cas on trouve que k = 0 (si xk n’est pas dégénérée,
la démonstration de la Proposition 4.2.8 garantit une inégalité stricte). On a donc
changé de base sans améliorer le coût : c’est le phénomène du cyclage qui peut
empêcher l’algorithme de converger. Il existe des moyens de s’en prémunir, mais en
pratique le cyclage n’apparaît jamais.
En l’absence de cyclage, l’algorithme du simplexe parcourt un sous-ensemble
des sommets de Xad en diminuant de façon stricte le coût. Comme il y a un nombre
4.2. PROGRAMMATION LINÉAIRE 129
Initialisation
Comment trouver une solution admissible basique lors de l’initialisation ? (Rap-
pelons que la condition d’admissibilité xB = B −1 b ≥ 0 n’est pas évidente en géné-
ral.) Soit on en connaît une à cause de la structure du problème. Par exemple,
pour le problème (4.4) qui possède m variables d’écart, − Idm est une base de la
matrice “globale” des contraintes d’égalité de (4.4). Si de plus b ≤ 0, le vecteur
(x0+ , x0− , λ0 ) = (0, 0, −b) est alors une solution admissible basique pour (4.4).
Dans le cas général, on introduit une nouvelle variable y ∈ Rm , un nouveau
vecteur coût k = (1, ..., 1) et un nouveau programme linéaire
Inversion de la base
Tel que nous l’avons décrit l’algorithme du simplexe demande l’inversion à
chaque étape de la base B k , ce qui peut être très coûteux pour les problèmes de
130 CHAPITRE 4. PROGRAMMATION LINÉAIRE
grande taille (avec beaucoup de contraintes puisque l’ordre de B k est égal au nombre
de contraintes). On peut tirer parti du fait que B k+1 ne diffère de B k que par une co-
lonne pour mettre au point une meilleur stratégie. En effet, si c’est la j-ème colonne
qui change, on a
1 l1
... ..
. 0
..
1 .
B k+1 = B k E k avec E k = lj ,
.
.. 1
.. ...
0 .
ln 1
et E k est facile à inverser
−l1
1
.. ..
. . 0
1 −lj−1
k −1 1
(E ) = 1 .
lj
−lj+1 1
.. ..
0 . .
−ln 1
On utilise donc la formule, sous forme factorisée,
(B k )−1 = (E k−1 )−1 (E k−2 )−1 · · · (E 0 )−1 (B 0 )−1 .
Exercice 4.2.3 Résoudre par l’algorithme du simplexe le programme linéaire
min x1 + 2x2
x1 ≥0, x2 ≥0, x3 ≥0, x4 ≥0, x5 ≥0
inf c · x. (4.9)
x∈Rn tel que Ax=b, x≥0
Remarquons qu’en pratique la contrainte x > 0 n’en est pas une car elle n’est jamais
active : quand on minimise (4.11) on ne peut pas “s’approcher” du bord de x > 0
sous peine de faire “exploser” le potentiel π(x) vers +∞. Il s’agit donc d’un exemple
de méthode de pénalisation intérieure.
Le principe de l’algorithme de trajectoire centrale est de minimiser (4.11) par
une méthode de Newton pour des valeurs de plus en plus petites de µ. En effet,
lorsque µ tend vers zéro, le problème pénalisé (4.11) tend vers le programme linéaire
(4.9).
Cette hypothèse est nécessaire car on peut construire des exemples où, Xad n’étant
pas borné, le programme linéaire (4.9) admet une solution optimale x0 mais le pro-
blème pénalisé (4.11) n’admet pas de solution optimale (voir l’Exercice 4.2.6 pour
un tel contre-exemple). •
Démonstration. La fonction µπ(x) + c · x est strictement convexe sur Rn+ donc,
s’il existe un point de minimum à (4.11), il est unique. Par contre, cette fonction
0
n’est pas infinie à l’infini donc on a besoin de l’hypothèse que Xad est borné. Par
0
ailleurs, Xad est convexe mais pas fermé, donc on ne peut pas utiliser directement
0
le Théorème 2.2.1 bien que J soit continue et Xad borné. Soit xn ∈ Xad0
une suite
0
minimisante de (4.11). Comme Xad = Xad est borné, on peut en extraire une sous-
suite, toujours notée xn , qui converge vers une limite xµ dans Xad . En fait, toutes
les composantes de xµ sont strictement positives et xµ ∈ Xad0
car sinon
ce qui serait une contradiction avec le fait que la valeur minimale de (4.11) est
0
bornée supérieurement puisque Xad n’est pas vide. Comme la fonction µπ(x) + c · x
est continue, on a donc obtenu
Autrement dit, xµ est la solution optimale de (4.11). On écrit les conditions d’op-
timalité (comme toujours dans ce chapitre, on suppose que la matrice A, de taille
m × n, est de rang m ≤ n, donc les contraintes sont qualifiées) pour xµ : il existe un
multiplicateur de Lagrange pµ ∈ Rm tel que
1
Axµ = b, c−µ + A∗ pµ = 0, (4.12)
xµ
où 1/x désigne le vecteur de Rn de composantes 1/xi . On veut passer à la limite,
quand µ → 0, dans (4.12). Comme l’ensemble Xad est borné, la suite xµ est bornée
et on peut en extraire une sous-suite qui converge vers une limite x0 ∈ Xad car
Xad est fermé. Par contre, en général la suite pµ n’est pas bornée, quand µ → 0, et
nous allons simplement montrer que la suite A∗ pµ est bornée dans Rn . On déduit de
(4.12) que
A∗ pµ + c ≥ 0. (4.13)
Au vu de (4.13) les composantes de A∗ pµ sont bornées inférieurement. Montrons
0
qu’elles sont aussi bornées supérieurement. Soit y ∈ Xad qui est supposé non vide.
µ
On multiplie la contrainte Ay = b par p pour obtenir
y · A ∗ pµ = b · pµ .
c · xµ − nµ + A∗ pµ · xµ = c · xµ − nµ + pµ · b = 0, (4.14)
4.2. PROGRAMMATION LINÉAIRE 133
donc
y · A∗ pµ = nµ − c · xµ ≤ C,
car la suite xµ est bornée. Comme y > 0, les composantes de A∗ pµ sont bornées
supérieurement. Autrement dit, A∗ pµ est une suite bornée dans Rn . On peut donc
en extraire une sous-suite qui converge. Par ailleurs, elle appartient à ImA∗ qui est
un sous-espace vectoriel fermé, donc il existe p0 ∈ Rm tel que la sous-suite A∗ pµ
converge vers A∗ p0 .
On passe alors à la limite dans la première égalité de (4.12), dans (4.13) et
(4.14)
Ax0 = b, A∗ p0 + c ≥ 0, c · x0 + A∗ p0 · x0 = 0.
Il s’agit des conditions d’optimalité de (4.9), donc en vertu du théorème de Kuhn et
Tucker (voir le Théorème 4.2.13 plus bas) x0 est la solution optimale de (4.9) et p0
est un multiplicateur de Lagrange pour ce problème. Si on suppose que (4.9) admet
une solution unique, alors toute la suite xµ converge vers x0 lorsque µ tend vers zéro.
2
Exercice 4.2.6 Donner un exemple de programme linéaire (4.9) qui admet une solution
optimale x0 mais pour lequel l’ensemble admissible Xad n’est pas borné et le problème
pénalisé (4.11) n’admet pas de solution optimale. Donner un autre exemple d’ensemble
0
admissible Xad non vide mais tel que Xad est vide.
4.2.4 Dualité
La théorie de la dualité (déjà évoquée lors de la Sous-section 2.6.3) est très utile
en programmation linéaire. Considérons à nouveau le programme linéaire standard
que nous appellerons primal (par opposition au dual)
inf c · x, (4.15)
x∈Rn tel que Ax=b, x≥0
sup p · b. (4.18)
p∈Rm tel que A∗ p−c≤0
134 CHAPITRE 4. PROGRAMMATION LINÉAIRE
Les programmes linéaires (4.15) et (4.18) sont dits en dualité. L’intérêt de cette
notion vient du résultat suivant qui est un cas particulier du Théorème de dualité
2.6.13.
Théorème 4.2.13 Si (4.15) ou (4.18) a une valeur optimale finie, alors il existe
x ∈ Xad solution optimale de (4.15) et p ∈ Pad solution optimale de (4.18) qui
vérifient
min c·x =c·x=p·b= max p·b (4.19)
x∈Rn tel que Ax=b, x≥0 p∈Rm tel que A∗ p−c≤0
Ax = b, x ≥ 0, A∗ p − c ≤ 0, x · (c − A∗ p) = 0. (4.20)
Si (4.15) ou (4.18) a une valeur optimale infinie, alors l’ensemble des solutions
admissibles de l’autre problème est vide.
c·x=b·p
alors x est solution optimale de (4.15) et p de (4.18). Ces deux propriétés permettent
de trouver facilement des bornes pour les valeurs optimales de (4.15) et (4.18), et
de tester si un couple (x, p) est optimal. •
Démonstration. Supposons que Xad et Pad sont non vides. Soit x ∈ Xad et p ∈ Pad .
Comme x ≥ 0 et A∗ p ≤ c, on a
c · x ≥ A∗ p · x = p · Ax = p · b,
Supposons maintenant que l’un des deux problèmes primal ou dual admet une
valeur optimale finie. Pour fixer les idées, admettons qu’il s’agisse du problème dual
(un argument symétrique fonctionne pour le problème primal). Alors, le Lemme
4.2.3 affirme qu’il existe une solution optimale p de (4.18). Si Xad n’est pas vide,
on se retrouve dans la situation précédente ce qui finit la démonstration. Montrons
donc que Xad n’est pas vide en utilisant encore le Lemme de Farkas 2.5.19. Pour
p ∈ Rm , on introduit les vecteurs de Rm+1
b p
b̃ = et p̃ = .
−b · p 1
Comme b̃ · p̃ ≤ 0 pour tout p̃ ∈ C, le Lemme de Farkas 2.5.19 nous dit qu’il existe
x̃ ∈ Rn tel que x̃ ≥ 0 et b̃ = Ãx̃, c’est-à-dire que x̃ ∈ Xad qui n’est donc pas vide.
Finalement, supposons que la valeur optimale du problème primal est (4.15)
−∞. Si Pad n’est pas vide, pour tout x ∈ Xad et tout p ∈ Pad , on a c · x ≥ b · p.
En prenant une suite minimisante dans Xad on obtient b · p = −∞, ce qui absurde.
Donc Pad est vide. Un raisonnement similaire montre que, si la valeur optimale de
(4.15) est infinie, alors Xad est vide. 2
Exercice 4.2.7 Utiliser la dualité pour résoudre “à la main” (et sans calculs !) le pro-
gramme linéaire
min 8x1 + 9x2 + 4x3 + 6x4
x1 ≥0, x2 ≥0, x3 ≥0, x4 ≥0
Exercice 4.2.8 Trouver le problème dual de (4.15) lorsqu’on dualise aussi la contrainte
x ≥ 0, c’est-à-dire qu’on introduit le Lagrangien
L(x, p, q) = c · x + p · (b − Ax) − q · x
Exercice 4.2.9 Vérifier que le problème dual de (4.18) est à nouveau (4.15).
Montrer que le problème dual peut se mettre sous la forme suivante, avec q ∈ Rm
sup −b · q . (4.22)
q≥0
A∗ q+c≥0
(c + A∗ q) · v = 0 et (b − Av) · q = 0. (4.23)
Les deux égalités de (4.23) sont appelées conditions des écarts complémentaires
(primales et duales, respectivement).
Lemme 4.3.1 Pour une matrice inversible A ∈ Zn×n les deux propriétés suivantes
sont équivalentes :
1. detA = ±1 ;
2. pour tout b ∈ Zn , on a A−1 b ∈ Zn .
Démonstration. Au vu des formules de Cramer pour la résolution du système
linéaire Ax = b, il est clair que si detA = ±1, alors A−1 b ∈ Zn pour tout second
membre b ∈ Zn . Pour prouver la réciproque montrons tout d’abord que A−1 ∈ Zn×n .
En effet, si on choisit b comme le i-ème vecteur de la base canonique de Rn , on en
déduit que la i-ème colonne de A−1 , qui coïncide avec A−1 b, est à coefficients entiers.
En variant 1 ≤ i ≤ n, on obtient donc A−1 ∈ Zn×n . En particulier, cela implique que
detA−1 ∈ Z et, comme 1 = detA detA−1 , cela montre que detA divise 1, c’est-à-dire
que detA = ±1. 2
Ce lemme motive la définition suivante qui est centrale dans cette section.
Définition 4.3.2 On dit qu’une matrice carrée A ∈ Zn×n est unimodulaire quand
detA = ±1, et qu’une matrice générale A ∈ Zm×n est totalement unimodulaire
quand toute sous-matrice carrée extraite de A est de déterminant ±1 ou 0.
En prenant des sous-matrices 1 × 1, on voit en particulier que les coefficients
d’une matrice totalement unimodulaire valent nécessairement ±1 ou 0. L’introduc-
tion des matrices totalement unimodulaires est motivée par le résultat suivant.
Remarque 4.3.4 La Proposition 4.3.3 est très importante en pratique car elle im-
plique que si l’on résout (4.24) par l’algorithme du simplexe "usuel", alors la solution
que l’on obtiendra (ainsi que toutes les solutions intermédiaires) seront entières car
le simplexe énumère des solutions basiques. Il n’y a rien à changer à l’algorithme
du simplexe pour obtenir des solutions entières, du moins sous les hypothèses de la
Proposition 4.3.3.
La Proposition 4.3.3 ne dit surtout pas que toutes les solutions optimales de
(4.24) sont entières. D’ailleurs, il ne peut en être ainsi, à moins que la solution
optimale ne soit unique, car, par linéarité du coût c · x, tout barycentre de solutions
optimales d’un programme linéaire est aussi solution optimale.
Il est intéressant de remarquer que l’on ne fait aucune hypothèse sur le vecteur
c qui apparait dans le coût. Celui-ci peut ne pas être à valeurs entières mais, dans
tous les cas, les solutions basiques (et donc une solution optimale au moins) seront
entières. •
Démonstration. Rappelons qu’une solution basique correspond à un réarrange-
ment des colonnes de A tel que A = (B, N ), où B est une matrice inversible m×m et
N est une matrice m×(n−m), et que dans cette base on a x = (xB , xN ) avec xN = 0,
xB ≥ 0 et BxB = b. Comme A est totalement unimodulaire, B est unimodulaire et
xB ∈ Zm . Comme xN = 0, on a bien x ∈ Zn . 2
Bien sûr, pour que la Proposition 4.3.3 tienne toutes ses promesses, il faut véri-
fier qu’il existe "suffisamment" de matrices totalement unimodulaires. Il n’existe pas
de caractérisation des matrices totalement unimodulaires mais on connait certaines
conditions suffisantes. On donne ici une telle condition, très utile en pratique, due
à Poincaré.
Exercice 4.3.1 Soit une matrice A ∈ Zm×n telle que les premiers éléments de chacune
de ses colonnes sont des 1 alors que tous les autres éléments sont des 0. Autrement dit,
pour tout 1 ≤ j ≤ n, il existe nj ∈ N tel que aij = 1 pour 1 ≤ i ≤ nj et aij = 0
4.3. VERS LA PROGRAMMATION LINÉAIRE EN NOMBRES ENTIERS 139
Nous faisons une brève présentation des problèmes de flots qui ne servent ici
que de prétexte à illustrer la résolution de programmes linéaires en nombres entiers
et, accessoirement, à donner des exemples de matrices vérifiant les hypothèses du
Lemme 4.3.5 de Poincaré. Pour plus de détails sur les problèmes de flots et notam-
ment pour d’autres algorithmes permettant de les résoudre, nous renvoyons à [5],
[19]. Les problèmes de flots sont des problèmes d’optimisation sur des graphes, qui
permettent de modéliser des réseaux (de transport, de communication, ou sociaux).
Pour illustrer notre propos, imaginons qu’il s’agit d’un réseau de transport.
On introduit donc un graphe orienté G = (N , A) où N est l’ensemble des
nœuds, numérotés de 1 à m, et A ⊂ N × N est l’ensemble des arcs ou arêtes
qui relient deux nœuds distincts. Les nœuds sont par exemple des villes et les arcs
des lignes de transport qui les relient. On note (i, j) l’arc qui relier le nœud i au
nœud j. Le graphe est dit orienté car l’arc (i, j) n’est pas le même que l’arc (j, i)
(voir la Figure 4.2) : cela permet de distinguer le sens du trajet, ce qui est utile
dans de nombreuses applications. Tous les nœuds ne sont pas forcément reliés entre
eux par un arc et on note n le nombre total d’arcs. Chaque arc (i, j) ∈ A possède
une capacité maximale uij ∈ R+ ∪ {+∞} (la quantité maximale de biens ou de
personnes pouvant circuler sur cette ligne de transport) ainsi qu’un coût cij ∈ R
(le prix unitaire du voyage sur cette ligne). Remarquez que, si la capacité maximale
est bien positive, on laisse la possibilité d’un coût négatif en cas de subventions de
certaines lignes... En chaque nœud i du graphe il est aussi possible d’avoir la donnée
d’un flot extérieur bi ∈ R, c’est-à-dire une quantité de biens ou de personnes qui
arrivent en i par d’autres moyens de transport que ceux modélisés par le graphe.
Si bi > 0, ce flot exogène est dit entrant, tandis que si bi < 0, le flot exogène est
sortant.
On appelle flot un vecteur x ∈ Rn , défini sur sur A par (i, j) 7→ xij , vérifiant
140 CHAPITRE 4. PROGRAMMATION LINÉAIRE
Dans la suite, nous supposerons systématiquement que cette condition (4.27) est
vérifiée. Un flot est dit admissible s’il satisfait en plus les contraintes de capacité
Il s’agit d’un programme linéaire qui admet toujours une solution si les coûts sont
positifs, cij ≥ 0, ou bien si les capacités sont finies, uij < +∞. Pour donner une forme
plus explicite de (4.29), sous la forme d’un programme linéaire, nous réécrivons la
loi des nœuds de Kirchoff (4.25) comme Ax = b, où la matrice A ∈ Rm×n , appelée
matrice d’incidence nœuds-arcs de G, est définie par
−1 si i = k,
Ai,(j,k) = 1 si i = j, (4.30)
0 sinon.
La i-ème ligne de A est la loi de Kirchoff au nœud i, tandis que les élément non nuls
de la (j, k)-ème colonne de A indiquent à quels nœuds (entrant ou sortant selon le
signe) est relié l’arc (j, k). Le problème de flot (4.29) est donc équivalent à
X
min cij xij . (4.31)
x∈Rn ,Ax=b, 0≤x≤u
(i,j)∈A
Ce problème (4.31) n’est pas sous la forme standard (4.1) à cause de la contrainte de
capacité maximale x ≤ u. Néanmoins, grâce à l’introduction d’une variable d’écart
4.3. VERS LA PROGRAMMATION LINÉAIRE EN NOMBRES ENTIERS 141
et aux clients j
M
X
0= vij − rj pour 1 ≤ j ≤ N.
i=1
On choisit donc les flots externes comme bi = si pour les entrepôts et bj 0 = −rj 0 pour
les clients (notez le changement de signe de rj pour prendre en compte le fait que
les clients reçoivent alors que les entrepôts expédient). Avec le choix des capacités
maximales infinies, uij = +∞, le problème de transport de l’Exemple 1.2.1 est donc
un
PMproblème PNde flot à coût minimum. Si les stocks sont supérieurs à la demande,
s
i=1 i > j=1 rj , alors il faut rajouter un client imaginaire qui reçoit tout le stock
non consommé par les autres clients avec un coût de transport nul. Le problème est
à nouveau un problème de flot à coût minimum. •
Il arrive souvent en pratique que l’on cherche des solutions entières d’un pro-
blème de flot. Par exemple, pour le problème de transport de l’Exemple 1.2.1, les
marchandises à livrer peuvent être des colis, et livrer une fraction de colis peut ne
pas avoir de sens. Il est donc naturel de se demander si un problème de flot à coût
minimum, dont les données sont entières, a des solutions optimales qui sont aussi
entières. Une réponse positive à cette question est apportée par le résultat suivant.
142 CHAPITRE 4. PROGRAMMATION LINÉAIRE
Théorème 4.3.7 On suppose que les flots extérieurs ou exogènes sont entiers, bi ∈
Z, ainsi que les capacités, uij ∈ N ∪ {+∞}. Alors les solutions basiques de (4.32)
sont entières. En particulier, le problème (4.29) de flot à coût minimal admet une
solution optimale entière x ∈ Nn .
Lemme 4.3.9 La matrice d’incidence nœuds-arcs A d’un graphe, définie par (4.30),
est totalement unimodulaire.
Démonstration. Par construction la matrice A a exactement un seul 1 et un seul
−1 dans chacune de ses colonnes car chaque arc possède un unique nœud de départ
et un unique nœud d’arrivée. On conclut grâce au Lemme 4.3.5 de Poincaré. 2
Terminons cette section en donnant un exemple, très important en pratique, de
problème de flot à coût minimum : le problème du flot maximal, qui s’intéresse
seulement aux capacités et pas aux coûts. Dans le graphe G on distingue deux
nœuds, notés s et p, appelés respectivement source et puits (voir la Figure 4.2), et
on suppose que s n’a pas de prédécesseur, autrement dit {i ∈ N | (i, s) ∈ A} = ∅,
et p n’a pas de successeur, c’est-à-dire {i ∈ N | (p, i) ∈ A} = ∅. Soit v ∈ R+ . On
appelle flot admissible de s à p de valeur v une solution x de Ax = b, 0 ≤ x ≤ u
avec
v si i = s,
bi = −v si i = p, (4.34)
0 sinon.
Exercice 4.3.2 Montrer que le problème du flot maximal est bien un exemple de
problème de flot à coût minimal. Pour cela, on rajoutera un arc reliant p à s au graphe
4.3. VERS LA PROGRAMMATION LINÉAIRE EN NOMBRES ENTIERS 143
Les deux dernières contraintes égalité indiquent que chaque femme trouvera un mari
et chaque homme une épouse (comme on maximise la réussite des mariages, on aurait
pu écrire ces contraintes comme des inégalités aussi). Les matrices v = (vij )1≤i,j≤N
qui vérifient ces contraintes sont précisément les matrices correspondantes aux per-
mutations σ ∈ SN .
Si maintenant on relâche (ou relaxe) la contrainte vij = 0 ou 1, en 0 ≤ vij ≤ 1,
les matrices qui vérifient les contraintes relaxées
N
X N
X
0 ≤ vij ≤ 1, vij = 1, vij = 1 pour 1 ≤ i, j ≤ N (4.36)
j=1 i=1
144 CHAPITRE 4. PROGRAMMATION LINÉAIRE
où on a effectué la même astuce que dans la Remarque 4.3.6 pour que la contrainte
de mariage pour chaque homme j ressemble à la loi des nœuds de Kirchoff. On peut
alors résumer ces contraintes d’égalité par Av = b avec b ∈ R2N le vecteur dont les N
2
premières composantes sont égales à 1 et les N dernières égales à −1 et A ∈ Z2N ×N
est une matrice dont chacune des colonnes a exactement un élément égal à 1, un
autre égal à −1 et tous les autres nuls. Le programme linéaire relaxé est donc
N XN
!
X
max aij vij . (4.37)
v≥0, Av=b
i=1 j=1
Remarque 4.3.10 Il existe une autre preuve du caractère entier des solutions opti-
males de (4.37) (voir [5]). L’ensemble des matrices doublement stochastiques, c’est-
à-dire qui vérifient (4.36), est un ensemble convexe compact dont on peut montrer
(théorème de Birkhoff - von Neumann) que les points extrêmes sont précisément les
matrices de permutations. On peut aussi montrer (théorème de Minkowski) que l’en-
semble des matrices doublement stochastiques est l’enveloppe convexe de ses points
extrêmes. Par linéarité de la fonction coût, on en déduit qu’il existe des solutions qui
sont des matrices de permutation, autrement dit des solutions optimales entières. •
Chapitre 5
CONTRÔLABILITÉ DES
SYSTÈMES DIFFÉRENTIELS
u : [0,T ] → Rk
nous permet d’agir sur le système afin d’en modifier l’état. On dit que u est le
contrôle. Une fois le contrôle u fixé, (5.1) est un problème de Cauchy. Afin
d’expliciter le fait que la trajectoire x, solution de (5.1), dépend du contrôle u, nous
la noterons souvent xu , et nous écrirons (5.1) sous la forme
u ∈ L1 ([0,T ]; Rk ),
145
146 CHAPITRE 5. CONTRÔLABILITÉ DES SYSTÈMES DIFFÉRENTIELS
et nous serons parfois amenés à faire des hypothèses un peu plus fortes sur le contrôle,
comme par exemple que u prend ses valeurs dans un sous-ensemble fermé non-vide
U de Rk , ce que nous noterons u ∈ L1 ([0,T ]; U ) ; nous ferons parfois des hypothèses
d’intégrabilité plus forte en temps, comme par exemple L2 ([0,T ]; U ) ou L∞ ([0,T ]; U ).
Rappelons à toutes fins utiles que l’espace L1 ([0,T ]; Rk ) est équipé de la norme
Z T
kukL1 ([0,T ];Rk ) = |u(s)|Rk ds,
0
où |·|Rk désigne la norme euclidienne sur Rk . (On peut remplacer la norme euclidienne
par toute autre norme sur Rk .) Rappelons que la notation ∗ désigne la transposition
des vecteurs ou des matrices ; on écrit donc x∗ y pour le produit scalaire entre deux
vecteurs et Z ∗ pour la transposée de la matrice Z.
Définition 5.1.1 (Systèmes de contrôle linéaires) On dit que (5.2) est un sys-
tème de contrôle linéaire. On dit que ce système est autonome (ou station-
naire) lorsque les matrices A et B ne dépendent pas du temps. Plus généralement,
on dit que le système de contrôle linéaire est instationnaire lorsqu’il s’écrit sous
la forme
Si une fonction F est absolument continue sur [0,T ], alors elle est continue sur [0,T ]
et elle est dérivable presque partout, de dérivée égale à f .
5.1. CONTRÔLABILITÉ DES SYSTÈMES LINÉAIRES 147
On notera que cette expression a bien un sens pour u ∈ L1 ([0,T ]; Rk ) car la fonction
s 7→ e(t−s)A est bornée sur [0,T ].
Démonstration. (1) Supposons d’abord que rang(C) < d. Par conséquent les lignes
de C sont liées et il existe un vecteur Ψ ∈ Rd , Ψ 6= 0, tel que
Ψ∗ B = Ψ∗ AB = · · · = Ψ∗ Ad−1 B = 0 (∈ Rk ),
où Ψ∗ désigne le transposé de Ψ (Ψ∗ est un vecteur ligne). D’après le théorème de
Cayley-Hamilton, il existe des réels s0 , · · · , sd−1 tels que
Ad = s0 Id + · · · + sd−1 Ad−1 ,
où Id est la matrice identité dans Rd×d . On en déduit par récurrence que Ψ∗ Ak B =
0 pour tout k ∈ N, puis que Ψ∗ etA B = 0 pour tout t ∈ [0,T ]. Par conséquent,
Ψ∗ Φ(u) = 0 pour tout contrôle u, i.e., l’application Φ ne peut être surjective.
(2) Réciproquement, si l’application Φ n’est pas surjective, il existe un vecteur Ψ ∈
Rd , Ψ 6= 0, tel que
Z T
∗
Ψ e(T −s)A Bu(s) ds = 0, ∀u ∈ L∞ ([0,T ]; Rk ).
0
∗
En choisissant le contrôle u(s) = B ∗ e(T −s)A Ψ, qui est bien dans L∞ ([0,T ]; Rk ), on
en déduit que
Ψ∗ etA B = 0 (∈ Rk ), ∀t ∈ [0,T ].
En t = 0, il vient Ψ∗ B = 0, puis en dérivant par rapport à t, il vient Ψ∗ AB = 0 et
ainsi de suite ; d’où
Ψ∗ B = Ψ∗ AB = · · · = Ψ∗ Ad−1 B = 0 (∈ Rk ).
La matrice C ne peut donc être de rang maximal. 2
5.1. CONTRÔLABILITÉ DES SYSTÈMES LINÉAIRES 149
ou encore ẍ(t) = − R
L
1
ẋ(t) − LC x(t) + L1 u(t) On obtient le système de contrôle linéaire
(avec d = 2, k = 1) sous la forme
0 1 0 x(t)
Ẋ(t) = −1 −R X(t) + 1 u(t), X(t) = .
LC L L
ẋ(t)
si bien que Ψ∗ e(T −s)A B = 0 pour tout s ∈ [0,T ]. Par la formule de Duhamel, on
obtient Ψ∗ (xu (T ) − eT A x0 ) = 0, ce qui montre que xu (T ) est dans un hyperplan
affine. Par conséquent, le système n’est pas contrôlable. 2
est inversible.
Démonstration. Identique au cas autonome. 2
La matrice KT n’est donc pas inversible, si bien que le système (5.9) n’est pas
contrôlable. Le problème vient du fait que la matrice R(s)−1 B(s) est indépendante
de s. En revanche, si le vecteur B était constant (et non-nul), le système serait
contrôlable car B et AB seraient alors des vecteurs orthogonaux non-nuls, si bien
que la matrice de Kalman C = (B, AB) serait de rang plein.
A1 A2
d(x1, A2)
x1
d(x2, A1)
x2
où xi est la trajectoire associée au contrôle ui , i ∈ {1, 2}. Posons u(s) = θu1 (s) +
(1 − θ)u2 (s), pour tout s ∈ [0,t]. La fonction u est mesurable et cette fonction est à
valeurs dans U grâce à la convexité du sous-ensemble U . De plus, par linéarité, la
5.1. CONTRÔLABILITÉ DES SYSTÈMES LINÉAIRES 153
Lemme 5.1.19 (Lyapunov) Soit t > 0. Soit une fonction f ∈ L1 ([0,t]; Rn ). Alors,
le sous-ensemble Z
f (s) ds | E ⊂ [0,t] mesurable (5.14)
E
est un sous-ensemble convexe de Rn .
154 CHAPITRE 5. CONTRÔLABILITÉ DES SYSTÈMES DIFFÉRENTIELS
Remarque 5.1.20 On peut montrer que l’ensemble atteignable pour des contrôles
à valeurs dans U est le même que pour des contrôles à valeurs dans conv(U ) (l’en-
veloppe convexe de U ).
x
u ≡ +1
x(t) ∈ A[−1,1](t, x0) = A{−1,1}(t, x0)
x0 t
u ≡ −1
Figure 5.2 – Ensemble atteignable par un point matériel dont on contrôle la vitesse
dans U = [−1, 1].
Remarque 5.2.3 Les hypothèses du Lemme 5.2.2 sont bien vérifiées dans le cas
linéaire. L’ensemble atteignable varie donc continûment en temps dans ce cas.
Démonstration. Soit > 0. On va montrer qu’il existe δ > 0 tel que
∀t1 , t2 ∈ [0,T ], |t1 − t2 | ≤ δ =⇒ d(A1 , A2 ) ≤ ,
où A1 = A(t1 , x0 ) et A2 = A(t2 , x0 ). Supposons pour fixer les idées que t2 > t1 . Soit
x2 ∈ A2 . Il existe donc un contrôle u ∈ L∞ ([0,t2 ]; U ) tel que
Z t2
x2 = x 0 + f (s, x(s), u(s)) ds.
0
156 CHAPITRE 5. CONTRÔLABILITÉ DES SYSTÈMES DIFFÉRENTIELS
Ceci montre que d(x2 , A1 ) ≤ |x2 − x1 |Rd ≤ C|t2 − t1 |. On raisonne de même pour
x1 ∈ A1 , ce qui conclut la preuve. 2
(iv) pour tout (t, x) ∈ [0,T ]×Rd , l’ensemble des vecteurs vitesse K(t, x) := {f (t, x, u) | u ∈
U } est un sous-ensemble convexe de Rd .
Alors, pour tout t ∈ [0,T ], l’ensemble atteignable A(t, x0 ) est un sous-ensemble com-
pact de Rd .
Remarque 5.2.5 Les hypothèses du Lemme 5.2.4 sont bien vérifiées dans le cas
linéaire avec U convexe. L’ensemble atteignable est donc compact dans ce cas.
Démonstration. On se place dans l’espace de Hilbert V = L2 ([0,T ]; Rd ) et on
va montrer la compacité de l’ensemble atteignable A(T, x0 ). La preuve utilise des
notions de topologie faible dans les espaces de Hilbert (voir la Sous-section 2.3.3
pour cette notion).
(1) Soit (yn )n∈N une suite d’éléments de A(T, x0 ) ⊂ Rd . Soit (un )n∈N une suite de
contrôles dans L∞ ([0,T ]; U ) et (xn )n∈N la suite de trajectoires correspondantes dans
AC([0,T ]; Rd ) menant de x0 à yn . Posons gn (s) = f (s, xn (s), un (s)) ds pour tout
n ∈ N et s ∈ [0,T ]. On a
Z t
xn (t) = x0 + gn (s) ds, ∀t ∈ [0,T ] et yn = xn (T ).
0
D’après les hypothèses, la suite (gn )n∈N est bornée dans V . En invoquant le Lemme 2.3.13
sur la compacité faible dans les espaces de Hilbert, on en déduit qu’à une sous-suite
près, la suite (gn )n∈N converge vers une fonction g ∈ V pour la topologie faible. On
définit la trajectoire x ∈ AC([0,T ]; Rd ) en posant
Z t
x(t) = x0 + g(s) ds, ∀t ∈ [0,T ].
0
5.2. CONTRÔLABILITÉ DES SYSTÈMES NON-LINÉAIRES 157
Rt Rt
Par convergence faible, on a 0
gn (s) ds = (gn , 1[0,t] )V → (g, 1[0,t] )V = 0
g(s) ds, i.e.,
En particulier, on a donc
lim yn = x(T ).
n→+∞
Il reste à montrer que la trajectoire x(t) peut bien être engendrée par un contrôle
u ∈ L∞ ([0,T ]; U ).
(2) Posons θn (s) = f (s, x(s), un (s)) et introduisons l’ensemble
de sorte que (θn )n∈N est une suite de Θ. Par hypothèse, K(s, x(s)) est un sous-
ensemble convexe de Rd pour tout s ∈ [0,T ]. On en déduit que Θ est un sous-
ensemble convexe de V . De plus, Θ est fermé dans V car la convergence dans V
implique la convergence p.p. d’une sous-suite, et K(s, x(s)) est fermé dans Rd . Grâce
au Lemme 2.3.15 sur la fermeture faible des convexes dans les espaces de Hilbert,
on en déduit que Θ est faiblement fermé dans V . De plus, comme la suite (θn )n∈N
est bornée dans V , on déduit du Lemme 2.3.13 qu’elle converge faiblement, à une
sous-suite près, vers une fonction θ ∈ Θ. Il existe donc une fonction u : [0,T ] → U
telle que θ(s) = f (s, x(s), u(s)) p.p. dans [0,T ], et la fonction u peut être choisie
mesurable (cf. la Section 8.2 pour plus de précisions sur ce point). Pour tout ϕ ∈ V ,
on a
Z T Z T Z T
gn (s)ϕ(s) ds = θn (s)ϕ(s) ds+ (f (s, xn (s), un (s))−f (s, x(s), un (s))ϕ(s) ds.
0 0 0
(5.16)
Comme |f (s, xn (s), un (s))−f (s, x(s), un (s))|Rd ≤ C|xn (s)−x(s)|Rd et |xn (s)−x(s)|Rd
tend vers zéro p.p. dans [0,T ], le deuxième terme au membre de droite de (5.16) tend
vers zéro (invoquer le théorème de convergence dominée de Lebesgue). En outre, par
RT RT
convergence faible, on a 0 g(s)ϕ(s) ds = 0 θ(s)ϕ(s) ds, i.e., g(s) = θ(s) p.p. dans
[0,T ]. En conclusion, on a bien g(s) = f (s, x(s), u(s)) p.p. sur [0,T ]. 2
δu ∈ L∞ ([0,T ]; Rk ),
Remarque 5.2.9 Dans le Lemme 5.2.8 la notion de différentielle de Fréchet est celle
pour une application définie sur un espace de Banach B et la notation hL, δui désigne
le produit de dualité entre une application linéaire continue L (qui appartient au
dual B 0 ) et un élément δu ∈ B. La Définition 2.4.1 ne considérait que des applications
définies sur un espace de Hilbert mais, comme expliqué dans la Remarque 2.4.2, il
5.2. CONTRÔLABILITÉ DES SYSTÈMES NON-LINÉAIRES 159
ẋu+δu (t) − ẋu (t) = f (t, xu+δu (t), u(t) + δu(t)) − f (t, xu (t), u(t))
= ∂f
∂x
(t, xu (t), u(t))(xu+δu (t) − xu (t)) + ∂f
∂u
(t, xu (t), u(t))δu(t) + o(δu)
= Au (t)(xu+δu (t) − xu (t)) + Bu (t)δu(t) + o(δu),
LE SYSTÈME
LINÉAIRE-QUADRATIQUE
V = L2 ([0,T ]; Rk ). (6.2)
161
162 CHAPITRE 6. LE SYSTÈME LINÉAIRE-QUADRATIQUE
le critère
Z T Z T
1 ∗ 1 1
J(u) = u(t) Ru(t) dt + exu (t)∗ Qexu (t) dt + exu (T )∗ Dexu (T ), (6.3)
2 0 2 0 2
ξ(0)
x(t) x(T )
x0 ξ(T )
ξ(t)
avec
1 T
Z
JR (u) = u(t)∗ Ru(t) dt,
2 0
1 T
Z
1
JQD (u) = exu (t)∗ Qexu (t) dt + exu (T )∗ Dexu (T ).
2 0 2
J 0 (u) = Ru + B ∗ p ,
Remarque 6.2.3 [État adjoint] Attention, il n’y a pas de condition initiale sur
p, mais une condition finale en T . Par ailleurs, dans la littérature, la convention
est parfois de définir l’état adjoint comme un vecteur ligne pb := p∗ . Dans ce cas, le
p(t)A−ex (t)∗ Q, pour tout t ∈ [0,T ],
système différentiel rétrograde s’écrit dtd pb(t) = −b
et pb(T ) = ex (T ) D. Enfin, le contrôle optimal est u(t) = −R−1 B ∗ pb(t)∗ .
∗
Comme les matrices D et Q sont positives et que la matrice R est définie positive,
on en déduit que B ∗ p(t) = 0 sur [0,T ]. Donc, u(t) = 0, ce qui implique que x(t) = 0,
et ce qui implique enfin que p(t) = 0.
qui réalise une pondération au sens des moindres carrés entre l’atteinte de la cible
nulle sur [0,T ] et le fait que le contrôle ne soit pas trop grand dans L2 ([0,T ]; R). Ce
problème rentre dans le cadre du système LQ introduit à la section 6.1 en posant
A = 0, B = 1, R = 1, Q = 1, D = 0, ξ ≡ 0.
u(t) = −p(t),
dp
(t) = −x(t), ∀t ∈ [0,T ], p(T ) = 0.
dt
On a donc
d x(t) 0 −1 x(t) tZ cosh(t) − sinh(t)
= , e = ,
dt p(t) −1 0 p(t) − sinh(t) cosh(t)
| {z }
=Z
si bien que
x0
x0 tanh(T ) x0/ cosh(T )
t
T
−x0 tanh(T )
On notera que l’état adjoint initial est, à ce stade, encore inconnu. Afin de le dé-
terminer, on utilise la condition en t = T sur l’état adjoint, à savoir p(T ) = 0. On
obtient facilement que p(0) = x0 tanh(T ). En conclusion, l’extrémale s’écrit
1
x(t) = x0 cosh(T )
cosh(T − t),
1
p(t) = x0 cosh(T )
sinh(T − t),
1
u(t) = −p(t) = −x0 cosh(T )
sinh(T − t).
Cette extrémale est illustrée à la Figure 6.2.
dx ∂H
(t) = Ax(t) + Bu(t) + f (t) = (t, x(t), p(t), u(t)), (6.9a)
dt ∂p
dp ∂H
(t) = −A∗ p(t) − Q(x(t) − ξ(t)) = − (t, x(t), p(t), u(t)), (6.9b)
dt ∂x
et d’autre part que
∂H
(t, x(t), p(t), u(t)) = 0. (6.10)
∂u
Comme la fonction v 7→ H(t, x, p, v) est fortement convexe en v ∈ Rk pour tout
triplet (t, x, p) fixé dans [0,T ] × Rd × Rd , l’équation (6.10) ne signifie rien d’autre
que
u(t) = arg minv∈Rk H(t, x(t), p(t), v), ∀t ∈ [0,T ].
Il s’agit du principe du minimum de Pontryaguine (PMP) dans le cas parti-
culier du système LQ. Résumons ce résultat sous la forme d’une proposition.
avec
dx ∂H
(t) = (t, x(t), p(t), u(t)) = Ax(t) + Bu(t), x(0) = x0 , (6.11a)
dt ∂p
dp ∂H
(t) = − (t, x(t), p(t), u(t)) = −A∗ p(t) − Qex (t), p(T ) = Dex (T ), (6.11b)
dt ∂x
où ex (t) = x(t) − ξ(t).
b := −H et aboutir
Remarque 6.3.3 [Convention de signe] On aurait pu définir H
à un principe du maximum pour H.
b
∂H
(t, x, p, u) = 0.
∂t
Dans ce cas onn dit que le Hamiltonien est autonome.
Théorème 6.4.1 On suppose que dérive et cible sont nulles. Il existe une unique
matrice P ∈ C 1 ([0,T ]; Rd×d ) solution de l’équation de Riccati
et on a
p(t) = P (t)x(t), ∀t ∈ [0,T ],
si bien que le contrôle optimal s’écrit sous forme de boucle fermée :
et on a X (0) = Id.
(2) Inversibilité de X (t). Nous allons montrer que la matrice X (t) est inversible pour
tout t ∈ [0,T ]. Pour ce faire, on raisonne par l’absurde. Soit s ∈ [0,T ] et 0 6= x0 ∈ Rd
tels que x(s) = X (s)x0 = 0. On a nécessairement s > 0 car X (0) = Id. De plus, on
a vu que
d
p(t)∗ x(t) = −x(t)∗ Qx(t) − (B ∗ p(t))∗ R−1 B ∗ p(t).
dt
En intégrant de s à T , et comme x(s) = 0, il vient
Z T
∗ ∗ ∗ ∗ −1 ∗
0 = (Dx(T )) x(T ) + x(t) Qx(t) + (B p(t)) R B p(t) dt ≥ 0.
s
On a donc dx dt
(t) = Ax(t) et x(s) = 0 ; d’où x(t) = 0 sur [s, T ]. De même, comme on a
dp ∗
dt
(t) = −A p(t) et p(T ) = Dx(T ) = 0, il vient p(t) = 0 sur [s, T ]. On en déduit que
(x, p) vérifie un système différentiel linéaire avec conditions finales x(T ) = p(T ) = 0.
6.4. ÉQUATION DE RICCATI : FEEDBACK 171
Ceci implique que x(t) = p(t) = 0 sur [0, T ] ; en particulier, on obtient x0 = 0, d’où
la contradiction.
(3) Équation de Riccati. On pose
P (t) = P(t)X (t)−1 , ∀t ∈ [0,T ].
Par construction, on a P ∈ C 1 ([0,T ]; Rd×d ). De plus, on constate que
dp dP dx
(t) = (t)x(t) + P (t) (t)
dt dt dt
dP −1 ∗
= (t) + P (t)A − P (t)BR B P (t) x(t),
dt
et par ailleurs, on a également dp
dt
(t) = −A∗ p(t) − Qx(t). On en déduit que
dP ∗ −1 ∗
(t) + P (t)A + A P (t) − P (t)BR B P (t) + Q x(t) = 0,
dt
pour tout t ∈ [0,T ] et pour tout x0 ∈ Rd . Pour tout t ∈ [0,T ] fixé, le vecteur
x(t) décrit Rd lorsque x0 décrit Rd (car X (t) est inversible). Par conséquent, la
fonction t 7→ P (t) est bien solution de l’équation de Riccati pour tout t ∈ [0,T ].
En raisonnant de manière analogue, on constate que p(T ) = Dx(T ) = P (T )x(T ).
Comme x(T ) décrit Rd lorsque x0 décrit Rd , on conclut que P (T ) = D.
(4) Propriétés de P (t). La fonction t 7→ P (t) est solution d’un système différentiel
quadratique. La non-linéarité satisfait donc une condition de Lipschitz locale, ce qui
assure l’unicité de la solution. L’unicité prouve que P (t) est symétrique pour tout
t ∈ [0,T ] car la fonction t 7→ P (t)∗ satisfait la même équation. Afin d’établir la
positivité de P (t) pour tout t ∈ [0,T ], on raisonne comme suit. Soit x ∈ Rd . Posons
x0 = X (t)−1 x de sorte que x = x(t) où x est la trajectoire optimale issue de x0 .
Comme la fonction t 7→ p(t)∗ x(t) est décroissante, il vient
x∗ P (t)x = x(t)∗ P (t)x(t) ≥ x(T )∗ Dx(T ) ≥ 0,
ce qui montre que P (t) est semi-définie positive. Enfin, si la matrice D est définie
positive, cela entraîne x(T ) = 0, d’où x = X (t)X (T )−1 x(T ) = 0, i.e., la matrice
P (t) est alors définie positive.
(5) Valeur optimale du critère. Il vient
1 T
Z
1
J(u) = x(t) Qx(t) + u(t) Ru(t) dt + x(T )∗ Dx(T )
∗ ∗
2 0 2
Z T
1 1
= x(t)∗ Qx(t) − p(t)∗ Bu(t) dt + x(T )∗ Dx(T )
2 0 2
Z T
1 1
= x(t)∗ Qx(t) − p(t)∗ Bu(t) dt + p(T )∗ x(T )
2 0 2
Z T
1 d 1
p(t)∗ x(t) dt + p(T )∗ x(T )
= −
2 0 dt 2
1 1 1
= p(0)∗ x(0) = x(0)∗ P (0)x(0) = x∗0 P (0)x0 ,
2 2 2
ce qui conclut la preuve. 2
172 CHAPITRE 6. LE SYSTÈME LINÉAIRE-QUADRATIQUE
On note R(t) = e(T −t)A la résolvante associée à ce système différentiel (de dimension
2d et telle que R(T ) = Id). On pose
R1 (t) R2 (t)
R(t) = ∈ R(2d)×(2d) ,
R3 (t) R4 (t)
où les quatre blocs sont à valeurs dans Rd×d . On a x(t) = R1 (t)x(T ) + R2 (t)p(T )
et p(t) = R3 (t)x(T ) + R4 (t)p(T ). Or p(T ) = Dx(T ), si bien qu’en posant XT (t) =
R1 (t) + R2 (t)D et PT (t) = R3 (t) + R4 (t)D, il vient x(t) = XT (t)x(T ) et p(t) =
PT (t)x(T ). En conclusion, la matrice P (t) solution de l’équation de Riccati s’obtient
également à partir de la résolvante du système linéaire de taille 2d ci-dessus en posant
Cette expression est intéressante en pratique car elle évite de devoir résoudre un
système différentiel non-linéaire.
A = 0, B = 1, R = 1, Q = 1, D = 0, ξ ≡ 0.
On obtient P (t) = tanh(T − t). Le contrôle optimal se met alors sous forme de
boucle fermée
ce qui permet de retrouver l’expression ci-dessus liant u(t) à x(t). Enfin, la valeur
optimale du critère est J(u) = 21 x20 P (0) = 21 x20 tanh(T ).
Chapitre 7
PRINCIPE DU MINIMUM DE
PONTRYAGUINE
173
174 CHAPITRE 7. PRINCIPE DU MINIMUM DE PONTRYAGUINE
Nous allons formuler quelques hypothèses (en général, raisonnables) sur les dif-
férents ingrédients intervenant dans la formulation du problème de contrôle opti-
mal (7.4), à savoir la fonction f pour la dynamique et les fonctions g et h pour le
critère. Commençons par les hypothèses sur la dynamique. On suppose que
(a) f ∈ C 0 ([0,T ] × Rd × U ; Rd ) et f est de classe C 1 par rapport à x ;
(b) ∃C, |f (t, y, v)|Rd ≤ C(1 + |y|Rd + |v|Rk ), ∀t ∈ [0,T ], ∀y ∈ Rd , ∀v ∈ U ;
(c) Pour tout R > 0, ∃CR , | ∂f ∂x
(t, y, v)|Rd×d ≤ CR (1 + |v|Rd ), ∀t ∈ [0,T ], ∀y ∈
B(0, R), ∀v ∈ U .
Dans ces hypothèses, C et CR désignent des constantes génériques indépendantes
de (t, y, v), CR dépendant du rayon R de la boule fermée B(0, R) ; par la suite, nous
utiliserons les symboles C et CR avec la convention que les valeurs de C et de CR
peuvent changer à chaque utilisation tant qu’elles restent indépendantes du temps,
de l’état du système et de la valeur du contrôle. L’objectif des trois hypothèses ci-
dessus est d’assurer, pour tout contrôle u ∈ U, l’existence et l’unicité de la trajectoire
associée xu ∈ AC([0,T ]; Rd ).
Lemme 7.1.1 (Existence et unicité des trajectoires) Dans le cadre des hypo-
thèses (a), (b), (c) ci-dessus, pour tout contrôle u ∈ U, il existe une unique trajectoire
associée xu ∈ AC([0,T ]; Rd ) solution de (7.1).
Démonstration. Il s’agit d’une conséquence de la version locale du théorème de
Cauchy–Lipschitz avec une dynamique mesurable en temps uniquement (cf. le Théo-
rème 8.3.6). On considère le système dynamique ẋ(t) = F (t, x(t)) avec la fonction
F : [0,T ] × Rd → Rd telle que F (t, x) = f (t, x, u(t)). La fonction F est mesurable en
t, et elle est continue en x. De plus, F est localement lipschitzienne par rapport à x
puisque l’on a, pour tout t ∈ [0,T ] et tout x1 , x2 ∈ B(0, R),
|F (t, x1 ) − F (t, x2 )|Rd ≤ C0 (t)|x1 − x2 |Rd , C0 (t) = sup ∂f
∂x
(t, y, u(t)) d×d .
R
y∈B(0,R)
Théorème 7.2.2 (PMP) Si u ∈ U est un contrôle optimal, i.e., si u est une so-
lution de (7.4), alors en notant x = xu ∈ AC([0,T ]; Rd ) la trajectoire associée au
contrôle u et en définissant l’état adjoint p ∈ AC([0,T ]; Rd ) solution de
dp ∂h
(t) = −A(t)∗ p(t) − b(t), ∀t ∈ [0,T ], p(T ) = (x(T )) ∈ Rd , (7.7)
dt ∂x
où pour tout t ∈ [0,T ],
∂f ∂g
A(t) = (t, x(t), u(t)) ∈ Rd×d , b(t) = (t, x(t), u(t)) ∈ Rd , (7.8)
∂x ∂x
on a, p.p. t ∈ [0,T ],
u(t) ∈ arg min H(t, x(t), p(t), v), (7.9)
v∈U
Remarque 7.2.3 L’état adjoint p est solution d’un système linéaire (à (x, u) fixés)
instationnaire et rétrograde en temps. Ce système, ainsi que la condition finale sur
p(T ), sont bien définis grâce aux hypothèses (a) et (d) ci-dessus. De plus, ce système
admet une unique solution car la fonction b est bien intégrable en temps grâce à
l’hypothèse (f) et la fonction A est dans L1 ([0,T ]; Rd×d ) grâce à l’hypothèse (c).
Alors, on a inf u∈U J(u) = 0 et il n’existe pas de contrôle optimal. Pour le montrer,
on considère pour tout n ∈ N∗ la suite minimisante de contrôles
un(t)
xn(t)
Figure 7.1 – Illustration du Contre-exemple 7.2.6 : contrôle issu d’une suite mini-
misante et trajectoire associée.
Cette inégalité, toujours grâce au Théorème 2.5.1, équivaut au fait que u soit mini-
miseur sur K de la fonctionnelle
Z T
˜
J(u) = ∗
p(t) Bu(t) + g(t, x(t), u(t)) dt.
0
7.3. APPLICATION AU SYSTÈME LQ AVEC CONTRAINTES 179
En posant Φ(t, v) := p(t)∗ Bv + g(t, x(t), v), nous avons donc établi que
Z T Z T
Φ(t, u(t)) dt ≤ Φ(t, u(t)) dt, ∀u ∈ K.
0 0
Supposons par l’absurde que Φ(t, u(t)) > minv∈U Φ(t, v) sur un sous-ensemble de
[0,T ] de mesure strictement positive. Comme Φ(t, u(t)) ≥ minv∈U Φ(t, v) puisque
u(t) ∈ U par hypothèse, ceci implique que
Z T Z T
min Φ(t, v) dt < Φ(t, u(t)) dt.
0 v∈U 0
En posant û(t) := arg minv∈U Φ(t, v) sur [0,T ], on peut montrer (voir le Théorème
8.2.9 de sélection mesurable dont la démonstration sort du cadre de ce cours) que
la fonction û ainsi définie est bien mesurable. Elle est de plus à valeurs dans U par
construction, et comme U est borné par hypothèse, û est bien de carré sommable.
En conclusion, û ∈ K, ce qui fournit la contradiction attendue puisqu’il vient
Z T Z T Z T Z T
Φ(t, û(t)) dt = min Φ(t, v) dt < Φ(t, u(t)) dt ≤ Φ(t, û(t)) dt.
0 0 v∈U 0 0
Ainsi u(t) est minimiseur instantané de v 7→ p(t)∗ Bv + g(t, x(t), v), ce qui n’est rien
d’autre que minimiser le Hamiltonien par rapport à v. 2
K = L2 ([0,T ]; U ). (7.11)
Lemme 7.3.1 Il existe une unique solution au problème (7.12), i.e., la fonctionnelle
J définie par (7.13) admet un unique minimiseur sur le sous-ensemble K défini
par (7.11).
dp
(t) = −A∗ p(t) − Qex (t), ∀t ∈ [0,T ], p(T ) = Dex (T ), (7.14)
dt
1 1
H(t, x, p, u) = p∗ (Ax + Bu) + u∗ Ru + (x − ξ(t))∗ Q(x − ξ(t)). (7.15)
2 2
En appliquant le PMP (cf. le Théorème 7.2.2), on en déduit qu’une condition néces-
saire d’optimalité est que, p.p. t ∈ [0,T ], u(t) est un minimiseur de H(t, x(t), p(t), v)
sur U , i.e.,
u(t) ∈ arg minv∈U H(t, x(t), p(t), v). (7.16)
7.3. APPLICATION AU SYSTÈME LQ AVEC CONTRAINTES 181
Remarque 7.3.3 Le fait que le contrôle optimal u ∈ K soit une fonction lipschit-
zienne du temps montre que pour le système LQ avec contraintes, il n’y a pas de
phénomènes de type bang-bang pour le contrôle optimal.
Démonstration. (1) La fonctionnelle J étant convexe et différentiable sur V , une
condition nécessaire et suffisante d’optimalité pour le problème (7.12) est l’inéqua-
tion d’Euler (cf. le Théorème 2.5.1)
hJ 0 (u), v − uiV ≥ 0, ∀v ∈ K.
où la fonctionnelle
Z T
∗ ∗ 1 ∗
Jp : V → R, Jp (v) = v(t) B p(t) + v(t) Rv(t) dt
0 2
182 CHAPITRE 7. PRINCIPE DU MINIMUM DE PONTRYAGUINE
(2) Montrons que la fonction u] (t) ainsi définie est lipschitzienne en t sur [0,T ]. Soit
t1 , t2 ∈ [0,T ]. On a
|u] (t2 ) − u] (t1 )|Rk = |δu] |Rk ≤ λmin (R)−1 kB ∗ kRk×d |p(t2 ) − p(t1 )|Rd ,
où λmin (R) > 0 désigne la plus petite valeur propre de la matrice R. Comme la
fonction t 7→ p(t) est de classe C 1 en t, cela montre que la fonction t 7→ u] (t) est
lipschitzienne en t.
(3) En conclusion, la fonction u] : [0,T ] → Rk est mesurable (car lipschitzienne), de
carré sommable et à valeurs dans U . On a donc u] ∈ K. De plus, comme u(t) ∈ U
p.p. t ∈ [0,T ], l’inégalité suivante est satisfaite p.p. t ∈ [0,T ] :
1 1
u(t)∗ B ∗ p(t) + u(t)∗ Ru(t) ≥ u] (t)∗ B ∗ p(t) + u] (t)∗ Ru] (t).
2 2
En intégrant cette inégalité de 0 à T , il vient
Jp (u) ≥ Jp (u] ).
u
Soit ∈ U un contrôle optimal, de trajectoire associée (a, r)∗ . Comme ∂f ∂x
(x, u) =
ϕ(u) 0 ∂g
γu 0 et ∂x (x, u) = 0, l’état adjoint p = (pa , pr )∗ : [0,T ] → R2 est tel que
dp
a (t) = −ϕ(u(t))pa (t) − γu(t)pr (t),
dt ∀t ∈ [0,T ], (7.24)
dp
r (t) = 0,
dt
et la condition finale sur l’état adjoint est
On a donc
dpa
(t) = −ϕ(u(t))pa (t) + γu(t), pr (t) ≡ −1, ∀t ∈ [0,T ]. (7.26)
dt
Par ailleurs, le Hamiltonien est autonome (cf. la Définition 7.2.1) et s’écrit sous la
forme
H(x, p, u) = pa ϕ(u)a + γpr ua. (7.27)
184 CHAPITRE 7. PRINCIPE DU MINIMUM DE PONTRYAGUINE
La condition de minimisation (7.9) s’écrit, en utilisant le fait que a(t) 6= 0 pour tout
t ∈ [0,T ],
u(t) ∈ arg minv∈[0,1] ψ(t)v, (7.28)
où la fonction de commutation est donnée par
a
u
r
t∗ t
pa
Figure 7.2 – Trajectoire, état adjoint et contrôle optimal pour le modèle de ruche.
— Cas 2. t∗ > 0 (ce qui correspond au cas d’un horizon temporel T rela-
tivement grand) ; le contrôle optimal est u ≡ 0 sur [0, t∗ [ et u ≡ 1 sur
]t∗ , T ]. En effet, le contrôle u vérifie bien le PMP car dpdta (t) = (β − α)pa (t),
pa (t∗ ) = − αγ , d’où pa (t) = − αγ e(β−α)(t−t∗ ) < − αγ sur [0, t∗ ], si bien que la
fonction de commutation est positive, ce qui correspond bien à u(t) = 0.
L’ensemble {t ∈ [0,T ] | ψ(t) = 0} est réduit au singleton {t∗ } et est donc de
mesure nulle.
Une illustration de la trajectoire, de l’état adjoint et du contrôle optimal est présen-
tée à la Figure 7.2 dans le cas où il y a une commutation.
On rappelle les hypothèses qui avaient été introduites afin de garantir l’existence
et l’unicité d’une trajectoire xu pour un contrôle donné u ∈ U (cf. en particulier le
Lemme 7.1.1) et le fait que la fonctionnelle J(u) est bien définie :
1. f ∈ C 0 ([0,T ] × Rd × U ; Rd ) et f est de classe C 1 par rapport à x ;
2. ∃C, |f (t, y, v)|Rd ≤ C(1 + |y|Rd + |v|Rk ), ∀t ∈ [0,T ], ∀y ∈ Rd , ∀v ∈ U ;
3. Pour tout R > 0, ∃CR , | ∂f
∂x
(t, y, v)|Rd×d ≤ CR (1 + |v|Rd ), ∀t ∈ [0,T ], ∀y ∈
B(0, R), ∀v ∈ U ;
4. g ∈ C 0 ([0,T ] × Rd × U ; R) et g est de classe C 1 par rapport à x ; de plus,
h ∈ C 1 (Rd ; R) ;
5. Pour tout R > 0, ∃CR , |g(t, y, v)| ≤ CR (1 + |v|Rk ), ∀t ∈ [0,T ], ∀y ∈ B(0, R),
∀v ∈ U ;
∂g
6. Pour tout R > 0, ∃CR , | ∂x (t, y, v)|Rd ≤ CR (1 + |v|Rk ), ∀t ∈ [0,T ], ∀y ∈
B(0, R), ∀v ∈ U ;
7. Les fonctions g et h sont minorées respectivement sur [0,T ] × Rd × U et sur
Rd .
Dans ces hypothèses, C et CR désignent des constantes génériques indépendantes
de (t, y, v), CR dépendant du rayon R de la boule fermée B(0, R) ; comme précé-
demment, nous continuons à utiliser les symboles C et CR avec la convention que
les valeurs de C et de CR peuvent changer à chaque utilisation tant qu’ils restent
indépendants du temps, de l’état du système et de la valeur du contrôle.
Rappelons enfin l’énoncé du PMP (cf. le Théorème 7.2.2).
Théorème 7.5.1 (PMP) Si u ∈ U est un contrôle optimal, i.e., si u est une solu-
tion de (7.34), alors, en notant x = xu ∈ AC([0,T ]; Rd ) la trajectoire associée à u,
et en définissant l’état adjoint p ∈ AC([0,T ]; Rd ) solution de
dp ∂h
(t) = −A(t)∗ p(t) − b(t), ∀t ∈ [0,T ], p(T ) = (x(T )), (7.35)
dt ∂x
avec A(t) = ∂f ∂x
(t, x(t), u(t)) ∈ Rd×d et b(t) = ∂g
∂x
(t, x(t), u(t)) ∈ Rd pour tout t ∈
[0,T ], on a, p.p. t ∈ [0,T ],
On rappelle enfin qu’un triplet (x, p, u) satisfaisant les conditions ci-dessus est
appelé une extrémale et que le PMP ne fournit qu’une condition nécessaire
d’optimalité ; en revanche, il ne dit rien sur l’existence d’un contrôle optimal et il
ne fournit pas a priori de condition suffisante.
Démonstration. Nous allons nous contenter de donner une esquisse de la preuve,
en insistant sur les idées principales sans nécessairement fournir tous les détails tech-
niques pour certains résultats intermédiaires. Ce qui compte ici est donc davantage
7.5. PMP : ESQUISSE DE PREUVE 187
si bien que
ẏδ (s) = A(s)yδ (s), ∀s ∈ Iδ+ , yδ (t + δ) = f (t, x(t), v) − f (t, x(t), u(t)),
∂f
où on rappelle que A(s) = ∂x
(s, x(s), u(s)). On en déduit que
xδ (s) − x(s) = δyδ (s) + Φδ (s), ∀s ∈ Iδ+ , Φδ = o(δ) unif. sur Iδ+ .
En effet, on a vu que Φδ (t + δ) = o(δ) et Φ̇δ (s) = Ψδ (s) + A(s)Φδ (s), pour tout
s ∈ Iδ+ , où Ψδ (s) = o(δ) uniformément sur Iδ+ , car
Ψδ (s) = f (s, xδ (s), u(s)) − f (s, x(s), u(s)) − A(s)(xδ (s) − x(s)).
188 CHAPITRE 7. PRINCIPE DU MINIMUM DE PONTRYAGUINE
v xδ (t + δ) xδ
u
x(t) x
x(t + δ)
Iδ Iδ+ T Iδ Iδ+ T
0 ≤ g(t, x(t), v) − g(t, x(t), u(t)) + p(t)∗ (f (t, x(t), v) − f (t, x(t), u(t))),
Définition 8.1.1 Un espace de Hilbert réel est un espace vectoriel sur R, muni
noté hx, yi, qui est complet pour la norme associée à ce produit
d’un produit scalaire, p
scalaire, notée kxk = hx, xi. (On rappelle qu’un espace vectoriel normé est complet
si toute suite de Cauchy est une suite convergente dont la limite appartient à cet
espace.)
Dans tout ce qui suit nous noterons V un espace de Hilbert réel, et hx, yi son
produit scalaire associé.
xK ∈ K, hxK − x, xK − yi ≤ 0 ∀y ∈ K. (8.1)
191
192 ANNEXE : RAPPELS
kPK x − PK yk ≤ kx − yk ∀ x, y ∈ V . (8.2)
xW ∈ W, hxW − x, zi = 0 ∀z ∈ W.
Montrons que y n est une suite de Cauchy. En utilisant la symétrie du produit scalaire,
il vient
1 1 1
kx − (y n + y p )k2 + k (y n − y p )k2 = (d2n + d2p ).
2 2 2
1 n
n p
Or, par convexité de K, (y + y )/2 ∈ K, et kx − 2 (y + y p )k2 ≥ d2 . Par conséquent
ce qui montre que y n est une suite de Cauchy. Comme V est un espace de Hilbert, il
est complet, donc la suite y n est convergente vers une limite xK . Par ailleurs, comme
K est fermé, cette limite xK appartient à K. Par conséquent, on a d = kx − xK k.
Comme toute la suite minimisante est convergente, la limite est forcément unique,
et xK est le seul point de minimum de miny∈K kx − yk.
Soit xK ∈ K ce point de minimum. Pour tout y ∈ K et θ ∈ [0, 1], par convexité
de K, xK + θ(y − xK ) appartient à K et on a
kx − xK k2 ≤ kx − xK k2 + θ2 ky − xK k2 − 2θhx − xK , y − xK i,
0 ≤ −2hx − xK , y − xK i + θky − xK k2 .
Remarque 8.1.7 L’existence d’une base hilbertienne dénombrable n’est pas ga-
rantie pour tous les espaces de Hilbert. Néanmoins, on peut construire des bases
hilbertiennes pour les espaces de Hilbert séparables (i.e. qui contiennent une famille
dénombrable dense). •
On écrit alors X
x= hx, en ien .
n≥1
Démonstration. S’il existe une suite (xn )n≥1 de réels telle que limp→+∞ pn=1 xn en =
P
x, alors par projection sur en (et comme cette suite est par définition indépendante
de p) on a xn = hx, en i, ce qui prouve l’unicité de la suite (xn )n≥1 . Montrons mainte-
nant son existence. Par définition d’une base hilbertienne, pour tout x ∈ V et pour
tout > 0, il existe y, combinaison linéaire finie des (en )n≥1 , tel que kx − yk < .
Grâce au Théorème 8.1.3 on peut définir une application linéaire Sp qui, à tout
point z ∈ V , fait correspondre Sp z = zW , où zW est la projection orthogonale sur
le sous-espace vectoriel W engendré par les p premiers vecteurs (en )1≤n≤p . En vertu
de (8.1), (z − Sp z) est orthogonal à tout élément de W , donc en particulier à Sp z.
On en déduit que
kzk2 = kz − Sp zk2 + kSp zk2 , (8.4)
ce qui implique
kSp zk ≤ kzk∀z ∈ V.
Comme Sp z est engendré par les (en )1≤n≤p , et que (z − Sp z) est orthogonal à chacun
des (en )1≤n≤p , on vérifie facilement que
p
X
Sp z = hz, en ien .
n=1
qui n’est rien d’autre que la formule de sommation (8.3), dite de Parseval. 2
Définition 8.1.9 Soit V et W deux espaces de Hilbert réels. Une application linéaire
A de V dans W est dite continue s’il existe une constante C telle que
kAxkW ≤ CkxkV ∀x ∈ V.
La plus petite constante C qui vérifie cette inégalité est la norme de l’application
linéaire A, autrement dit
kAxkW
kAk = sup .
x∈V,x6=0 kxkV
Définition 8.1.10 Soit V un espace de Hilbert réel. Son dual V 0 est l’ensemble des
formes linéaires continues sur V , c’est-à-dire l’ensemble des applications linéaires
continues de V dans R. Par définition, la norme d’un élément L ∈ V 0 est
|L(x)|
kLkV 0 = sup .
x∈V,x6=0 kxk
d’où le résultat désiré avec y = L(z0 )z0 (l’unicité est évidente). D’autre part, on a
et
|L(x)| hx, z0 i
kLkV 0 = sup = L(z0 ) sup .
x∈V,x6=0 kxk x∈V,x6=0 kxk
Le maximum dans le dernier terme de cette égalité est atteint par x = z0 , ce qui
implique que kLkV 0 = kyk. 2
Théorème 8.1.12 (Séparation d’un point et d’un convexe) Soit K une par-
tie convexe non vide et fermée d’un espace de Hilbert V , et x0 ∈
/ K. Alors il existe
un hyperplan fermé de V qui sépare strictement x0 et K, c’est-à-dire qu’il existe une
forme linéaire L ∈ V 0 et α ∈ R tels que
x ∈ C 1 ([0,T ]; Rd ). (8.12)
Cette solution satisfait donc le système différentiel (8.10) pour tout t ∈ [0,T ].
On voit donc que même si l’application f est régulière en u, le fait que le contrôle
ne dépende pas continûment du temps fait que l’application F ne sera pas néces-
sairement continue en t. Afin de traiter cette situation, on dispose de la variante
suivante du Théorème 8.3.1 (la preuve utilise des arguments analogues à ceux évo-
qués ci-dessus). On renvoie le lecteur à la Définition 5.1.2 pour la notion de fonction
absolument continue.
∃C0 ∈ L1 ([0,T ]; R+ ),
p.p. t ∈ [0,T ], ∀x1 , x2 ∈ Rd , |F (t, x1 ) − F (t, x2 )|Rd ≤ C0 (t)|x1 − x2 |Rd . (8.16)
x ∈ AC([0,T ]; Rd ). (8.17)
Cette solution, qui est dérivable p.p. sur [0,T ], satisfait le système différentiel (8.10)
pour presque tout t ∈ [0,T ] ; elle vérifie également
Z t
x(t) = x0 + F (s, x(s)) ds, ∀t ∈ [0,T ]. (8.18)
0
203
204 BIBLIOGRAPHIE
[15] EKELAND I., TEMAM R., Analyse convexe et problèmes variationnels, Dunod,
Paris (1974).
[16] FLETCHER R., Practical methods of optimization. Vol. 1. John Wiley &
Sons, Ltd., Chichester, 1980. Unconstrained optimization, A Wiley-Interscience
Publication.
[17] GAUTHIER T., Calcul différentiel et analyse complexe, cours de 2ème année à
l’Ecole Polytechnique (2020).
[18] GOLSE F., LASZLO Y., PACARD F., VITERBO C., Analyse réelle, cours de
1ère année à l’Ecole Polytechnique (2019).
[19] GRIVA I., NASH S., SOFER A., Linear and Nonlinear Optimization, Second
Edition, SIAM, Philadelphia (2009).
[20] HERMES H., LASALLE J. P., Functional analysis and time optimal control.
Academic Press, New York-London, 1969. Mathematics in Science and Engi-
neering, Vol. 56.
[21] ISIDORI A., Nonlinear control systems. Communications and Control Engi-
neering Series. Springer-Verlag, Berlin, third edition, 1995.
[22] LEE E. B., MARKUS L., Foundations of optimal control theory. Robert E.
Krieger Publishing Co., Inc., Melbourne, FL, second edition, 1986.
[23] LIONS P.-L., Contrôle de modèles dynamiques. Cours polycopié. École Poly-
technique, 2016.
[24] MASSOT M., SERIES L., BREDEN M., PICHARD T., Introduction à l’analyse
numérique : des fondements mathématiques à l’expérimentation avec Jupyter,
cours de 2ème année à l’Ecole Polytechnique (2020).
[25] NESTEROV Y., Lectures on convex optimization, Springer Optimization and
Its Applications, 137, Springer, Cham, 2018.
[26] NOCEDAL J., WRIGHT S., Numerical optimization, Springer Series in Ope-
rations Research and Financial Engineering, New York (2006).
[27] PADBERG M., Linear optimization and extensions, Springer, Berlin (1999).
[28] ROCKAFELLAR R. T., WETS R. J.-B., Variational analysis, volume 317 of
Grundlehren der Mathematischen Wissenschaften [Fundamental Principles of
Mathematical Sciences]. Springer-Verlag, Berlin, 1998.
[29] SCHRIJVER A., Theory of linear and integer programming, Wiley, New York
(1986).
[30] SONTAG E. D., Mathematical control theory, volume 6 of Texts in Applied
Mathematics. Springer-Verlag, New York, second edition, 1998. Deterministic
finite-dimensional systems.
[31] TRELAT E., Contrôle optimal. Mathématiques Concrètes. Vuibert, Paris, 2005.
Théorie & applications.
[32] VINTER R., Optimal control. Systems & Control : Foundations & Applications.
Birkhäuser Boston, Inc., Boston, MA, 2000.
Index
205
206 INDEX
méthode de Newton, 93
méthode du gradient conjugué, 81
minimum global, 14
minimum local, 14
multiplicateurs de Lagrange, 38, 44
nœud, 139
nombres entiers, 136
opérateur, 194
parcimonie, 5, 89
pénalisation, 105, 131
point-selle, 48
polyèdre, 124
primal, 52, 133
principe du minimum de Pontryaguine,
168, 176
problème dual, 52, 133
problème primal, 52
programmation linéaire séquentielle, 111
programmation quadratique séquentielle,
112
programme linéaire, 122
projection orthogonale, 191
proximal, 90
résolvante, 150
rétroaction, 169
recherche opérationnelle, 3, 121
unimodulaire, 137