Éléments d’optimisation

Pierre Bernhard

1er septembre 2001

2
Table des matières

1 Prologue : rappels d’analyse 5

1.1 La droite réelle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.1.1 Propriétés fondamentales . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.1.2 Intervalles . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
1.1.3 Ouverts et fermés de R, ou “topologie” de R . . . . . . . . . . . . . . . . . 8
1.1.4 Inf, min, sup, max . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
1.2 Espace Rn . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.2.1 Propriétés algébriques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.2.2 Matrices positives . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
1.2.3 Propriétés topologiques . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
1.3 Dérivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
1.3.1 Continuité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14
1.3.2 Dérivées et dérivées partielles . . . . . . . . . . . . . . . . . . . . . . . . . 15
1.3.3 Dérivation en chaı̂ne et dérivées directionnelles . . . . . . . . . . . . . . . . 16
1.4 Existence, unicité, CNS, et toutes ces sortes de choses . . . . . . . . . . . . . . . . . 18
1.4.1 Existence et conditions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
1.4.2 Multiplicateurs, dualité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
1.4.3 Convexité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22

2 Recherche unidimensionnelle 29
2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
2.1.1 Objectif . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
2.1.2 Pente et dérivée numérique . . . . . . . . . . . . . . . . . . . . . . . . . . . 29
2.2 Méthodes directes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.2.1 Dichotomie . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.2.2 Suites de Fibonacci . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
2.2.3 Section dorée . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
2.3 Méthodes indirectes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.3.1 “Backtracking” . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.3.2 Méthode de Newton . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.3.3 Approximation polynômiale . . . . . . . . . . . . . . . . . . . . . . . . . . 35

3 Optimisation dans Rn 39
3.1 Bonnes fonctions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
3.2 Optimisation non contrainte . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39
3.2.1 Relaxation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39

3
4 TABLE DES MATIÈRES

3.2.2 Gradient à pas optimal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41

3.2.3 Méthode de Newton . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
3.3 Optimisation sous contraintes inégalité . . . . . . . . . . . . . . . . . . . . . . . . . 45
3.3.1 Position du problème . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
3.3.2 Gradient projeté . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46
3.3.3 Algorithme d’Uzawa . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48
3.3.4 Pénalisation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50
3.3.5 Méthode du chemin central . . . . . . . . . . . . . . . . . . . . . . . . . . . 52
3.4 Optimisation sous contraintes égalité . . . . . . . . . . . . . . . . . . . . . . . . . . 54
3.4.1 Contraintes affines . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
3.4.2 Contraintes nonlinéaires . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56

4 Programmation linéaire et programmation dynamique 59

4.1 Programmation linéaire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.1.1 Position du problème . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
4.1.2 Étude du polyèdre . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
4.1.3 L’algorithme du simplexe . . . . . . . . . . . . . . . . . . . . . . . . . . . 63
4.1.4 Rudiments de dualité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65
4.2 Programmation dynamique . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
4.2.1 Plus court chemin dans un graphe orienté . . . . . . . . . . . . . . . . . . . 67
4.2.2 Système dynamique et programmation dynamique . . . . . . . . . . . . . . 69
Chapitre 1

Prologue : rappels d’analyse

Les rappels qui suivent sont sans doute inutiles pour la plupart des élèves. Nous les donnons à fin
de référence. Ils se limitent aux concepts d’analyse qui seront utiles à ce cours. Nous conseillons au
lecteur de commencer sa lecture là où commence ce cours, soit au dernier paragraphe “Dérivation”
de ce chapitre, (car il contient des conventions de notations utilisées ensuite) et de ne se reporter aux
paragraphes antérieurs qu’en cas de besoin au cours de l’étude du cours.

1.1 La droite réelle

1.1.1 Propriétés fondamentales
On note par N l’ensemble des entiers. On appelle “rationnel” le rapport de deux entiers. On note
Q l’ensemble des rationnels.
L’ensemble des réels, appelé aussi “droite réelle”, noté R, contient les rationnels, lesquels sont
“denses dans R” (c’est à dire que tout réel peut être approché arbitrairement près par des rationnels).
Sa propriété fondamentale, que ce soit par construction (par “complétion” des rationnels) ou par
axiome, est que les suites de Cauchy convergent. (On dit que R est complet.)
On suppose connu le concept de limite : an → a si quelque soit ε > 0, il existe un entier N tel
que pour tout n supérieur à N , |a − an | < ε. Toutefois, il est plus adroit d’écrire cette définition de
façon à peine différente :

Définition 1.1 (Limite) On dit que la suite {an }n∈N tends vers a, et on note an → a, si
∀ε > 0, ∃N ∈ N : ∀n > N, an ∈ (a − ε, a + ε) .

(Voir ci-dessous notre notation pour les intervalles, mais ici le fait de prendre un intervalle ouvert ou
fermé n’a aucune importance. L’élégance pousse à choisir un ouvert.)
On rappelle à ce propos qu’une suite de Cauchy est une suite {an }n∈N telle que “à condition
d’attendre assez longtemps, les éléments en sont tous arbitrairement proches les uns des autres”.
Mathématiquement :

Définition 1.2 (Suites de Cauchy) Une suite réelle {an }n∈N est appelée suite de Cauchy si
∀ε > 0, ∃N : ∀m, n > N, |an − am | < ε .

Une conséquence fondamentale de la convergence des suites de Cauchy est la célèbre propriété
suivante :

5
6 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Proposition 1.1 (Suites monotones) Les suites monotones bornées convergent.

Démonstration : exercice.
Rappelons enfin une définition :

Définition 1.3 (Point d’accumulation) Le point c est appelé point d’accumulation de la suite {an }
si ∀ε > 0, l’intervalle (c − ε, c + ε) contient un nombre infini de points de la suite.

Bien sûr, une limite est un point d’accumulation, mais une suite peut avoir aucun ou plusieurs points
d’accumulation et (donc) ne pas converger. Au demeurant, on a le résultat facile suivant :

Proposition 1.2 Une suite converge si et seulement si elle a un unique point d’accumulation.

Un concept important est celui de sous-suite. Si {nk }k∈N est une suite strictement croissante
d’entiers tendant vers l’infini, (i.e. croissant au-delà de tout nombre donné), la suite {ank }k∈N est
appelée sous-suite de la suite {an }n∈N .

Proposition 1.3 Les sous-suites d’une suite convergente convergent vers la même limite que la suite.

Les sous-suites sont un outil utile surtout quand on ne sait pas si la suite converge. Elles servent à
démontrer des propriétés des points d’accumulation grâce au fait simple suivant :

Proposition 1.4 Si a est un point d’accumulation de la suite {an }n∈N , il existe une sous-suite qui
converge vers a.

En effet, il suffit de choisir un point xnk de la suite dans chaque intervalle (voir ci-dessous) de lon-
gueur 2/k entourant le point d’accumulation. (Dans Rn , on prendra les boules de rayon 1/k centrées
sur ce point.)

1.1.2 Intervalles
Définition 1.4 (Intervalles fermés et ouverts) On appelle intervalle fermé ou segment [a, b] l’en-
semble des nombres compris entre a et b bornes comprises :

[a, b] = {t ∈ R | a ≤ t ≤ b} .

On appelle intervalle ouvert (a, b) 1 l’ensemble des nombres compris entre a et b bornes non com-
prises :
(a, b) = {t ∈ R | a < t < b}

Nous faisons remarquer les propriétés essentielles suivantes.

Propriété (F) Soit I un intervalle fermé et tn une suite d’éléments de I. Si cette suite converge, soit
t sa limite, alors t appartient à I.
1. Nous avons choisi de privilégier la notation anglosaxonne (a, b) à la notation française ]a, b[ pour deux raisons :
d’une part elle est beaucoup plus répandue dans la littérature, d’autre part, la notation française est sujette à des expressions
peu claires. Le principal inconvénient de notre notation est qu’elle coı̈ncide avec celle du produit scalaire. La confusion est
cependant impossible dans la mesure où un intervalle a des bornes qui sont des nombres réels, un produit scalaire est entre
vecteurs.
1.1. LA DROITE RÉELLE 7

(Les intervalles fermés contiennent les limites de leurs suites convergentes.) Remarquons bien sûr
qu’une suite contenue dans un fermé peut ne pas converger, par exemple la suite définie par t2k = a,
t2k+1 = b. Cette affirmation ne porte que sur les suites convergentes.
Propriété (O) Soit I un intervalle ouvert, et t un élément de I. Il existe un réel ε positif (non nul) tel
que (t − ε, t + ε) ⊂ I.
(Les intervalles ouverts contiennent un intervalle de longueur 2ε non nulle centré en chacun de
leurs points.) Naturellement, ε dépend de t.
Remarque Les demi-droites (−∞, a] et [b, +∞) (c’est à dire respectivement l’ensemble des réels
inférieurs ou égaux à a et l’ensemble des réels supérieurs ou égaux à b) satisfont encore la propriété
(F). De même, les demi-droites ouvertes (−∞, a) et (b, +∞) satisfont la propriété (O). Nous ne les
appellerons pas “intervalles”, réservant ce nom à des intervalles bornés. On peut aussi faire remarquer
que R tout entier satisfait trivialement à la fois les propriétés (F) et (O). Cela n’en fait pas un intervalle
non plus !

Exercice 1.1 (Intervalles fermés emboı̂tés) Soit {In }n∈N une suite d’intervalles fermés emboı̂tés,
i.e. des intervalles [an , bn ] tels que si m > n, alors [am , bm ] ⊂ [an , bn ]. Montrer que les In “ten-
dent”, dans un sens que l’on précisera, vers un intervalle non vide [a, b]. (La propriété des intervalles
fermés emboı̂tés
\
In 6= ∅
n∈N

est parfois prise comme définition axiomatique de R à la place de la convergence des suites de Cauchy.
Elle lui est évidemment équivalente.)
Montrer en outre que si |bn − an | → 0, alors il existe un réel c tel que In → c au sens où pour
toute suite {tn ∈ In } d’éléments des In , tn → c.

Nous énonçons maintenant un théorème fondamental, bien que très facile :

Théorème 1.5 (Bolzano-Weierstrass) Toute suite infinie d’un intervalle fermé (borné, donc) admet
au moins un point d’accumulation. (On dit que les intervalles fermés (bornés) sont compacts.)

Démonstration On procédera par dichotomie et en utilisant le résultat de l’exercice 1.1 ci-dessus.

En effet, si on coupe l’intervalle fermé en deux, l’un au moins des deux demi-intervalles (pris fermé)
contient un nombre infini de points de la suite, et on recommence.

Exercice 1.2 Montrer que si une suite {In }n∈N d’intervalles fermés est telle que toute intersection
finie est non vide, c’est à dire que
\
∀N ∈ N, In 6= ∅,
n≤N

alors l’intersection infinie en est non vide (il existe des points de R qui appartiennent à tous les In ) :
\
In 6= ∅ .
n∈N

Faire un contre-exemple avec des intervalles ouverts, et un avec des demi-droites fermées.
8 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

1.1.3 Ouverts et fermés de R, ou “topologie” de R

Définition 1.5 (Sous-ensembles fermés et ouverts) Tout sous-ensemble de R qui possède la pro-
priété (F) est dit fermé. Tout sous-ensemble qui possède la propriété (O) est dit ouvert.

Proposition 1.6 Le complémentaire d’un fermé est un ouvert, et le complémentaire d’un ouvert est
un fermé.

Proposition 1.7 Les intersections de fermés sont des fermés, les unions finies de fermés sont des
fermés.
Les unions d’ouverts sont des ouverts, les intersections finies d’ouverts sont des ouverts.

Démonstration : exercice.
On donne sans démonstration une caractérisation des ouverts de R :

Proposition 1.8 Les ensembles ouverts de R sont des unions, finies ou dénombrables, d’intervalles
ouverts.

On se méfiera de ce que si les ouverts sont comparativement des objets “simples”, si on en croit
la caractérisation ci-dessus, leurs complémentaires les fermés peuvent être épouvantablement com-
pliqués. En particulier, il serait faux de prétendre que les fermés de R sont des unions finies ou
dénombrables d’intervalles fermés. Un contre-exemple célèbre est donné par l’ensemble de Cantor,
obtenu comme suit.
Exemple : le fermé de Cantor On fabrique une suite de fermés de la façon suivante : C0 = [0, 1],
C1 = [0, 1/3] ∪ [2/3, 1], puis on continue ainsi, chaque Cn étant une union d’intervalles fermés de
longueur 1/3n on enlève le tiers central ouvert de chacun d’entre-eux pour passer à Cn+1 . L’ensemble
C est la limite infinie, ou, si on veut, l’intersection de tous les Cn puisque ces ensembles sont emboı̂tés.
On se convainc facilement que la somme des longueurs des segments retirés vaut 1 (exercice).
On pourrait croire que ne reste que l’ensemble dénombrable des points de division. Qu’il n’en est
rien est montré par le fait qu’on vérifie aussi que 1/4 est intérieur à tous les Cn . En fait, si on adopte
une numération triadique des nombres (i.e. à base 3), on a gardé tous les nombres qui s’écrivent sans
chiffre 1. (On note 0, 111 . . . le nombre 0,2 et donc on le retire, et de même pour tous les nombres a
développement triadique fini se finissant par un 2.)
Les propriétés de cet ensemble sont assez curieuses pour qu’il ait intrigué les mathématiciens.
Nous ne nous y étendrons pas d’avantage, qu’il nous suffise pour nous souvenir que si les ouverts de
R sont “simples”, les fermés peuvent en être compliqués.
Faisons remarquer la propriété essentielle (pour nous) des fermés :

Théorème 1.9 Les fermés bornés de R sont compacts.

Démonstration Soit F un fermé borné de R. Puisqu’il est borné, on peut l’inclure dans un intervalle
fermé borné. Donc toute suite infinie de F comporte au moins un point d’accumulation. Soit donc une
sous-suite convergeant vers ce point d’accumulation. Comme elle est dans F qui est fermé, sa limite
—le point d’accumulation choisi—, est dans F .
Nous avons enfin besoin de deux ou trois éléments de vocabulaire :
1.1. LA DROITE RÉELLE 9

Définition 1.6 (Intérieur, adhérence, frontière) On appelle intérieur d’un ensemble E, et on note
◦
E, l’union de tous les ouverts qu’il contient.
On appelle adhérence ou fermeture de E, et on note Ē, l’intersection de tous les fermés qui le
contiennent.
On appelle frontière de E, et on note ∂E, l’adhérence de E privée de l’intérieur de E.

Sans autre commentaire, ces définitions peuvent sembler arbitraires. En fait, l’intérieur est le
◦
plus grand ouvert contenu dans E. On appelle point intérieur de E tout point de E, ou, de manière
◦
équivalente tout point qui est le centre d’un petit intervalle ouvert contenu dans E. Ainsi, E est l’en-
semble des points intérieurs de E.
De même, Ē est le plus petit fermé contenant E. On appelle point adhérant de E, ou valeur
d’adhérence, tout point qui peut être approché arbitrairement prés par une suite de points de E.
L’adhérence de E est l’ensemble de ses points adhérants.
Un point frontière de E enfin est caractérisé par le fait que tout intervalle centré en ce point
contient des points de E et des points de son complémentaire.
On parle parfois d’extérieur de E pour désigner l’intérieur de son complémentaire, ou de manière
équivalente le complémentaire de son adhérence.

1.1.4 Inf, min, sup, max

Ce cours concerne la minimisation d’une fonction, ou la maximisation, c’est la même chose en
remplaçant la fonction par son opposée. Le minimum sera souvent approché par une suite d’approxi-
mations successives. Il faut donc avoir un peu de vocabulaire concernant cette question.

Définition 1.7 (Inf) Étant donné un ensemble E de réels, on appelle son inf, et on note inf t∈E {t}
son plus grand minorant s’il existe, et on pose inf t∈E {t} = −∞ sinon.

On aura plus souvent la situation où les nombres dont on veut chercher l’inf sont des fonctions d’un
entier (cas d’une suite) ou d’un réel (cas d’une fonction à minimiser sur un sous-ensemble de R) ce
qui mènera à considérer des expressions comme inf n∈N {an }, ou inf t∈F {u(t)}. Ici, {an }n∈N est une
suite réelle, et t ∈ F ⊂ R 7→ u(t) une fonction réelle, parfois notée u(·).
Notons que, suivant l’usage des bonnes imprimeries, Latex écrit pour nous l’inf (ou le min, cf
ci-dessous) en caractères droits, (pourvu qu’on utilise la commande \ inf) même dans les formules
mathématiques, dont les lettres sont plutôt en italique, et dispose l’indice différemment dans les for-
mules détachées du texte 2 :
inf {an } , inf {u(t)} .
n∈N t∈F

Proposition 1.10 Il existe toujours une (des) suite(s) minimisante(s), c’est à dire une sous-suite
{ank }k∈N telle que ank → inf n∈N {an } dans le premier cas, et une suite de réels {tk } de F telle
que u(tk ) → inf t∈F {u(t)} dans le second.

Démonstration Soit ū = inf t∈F {u(t)}. Supposons qu’il ne puisse pas être approché arbitrairement
près par des éléments de {u(t) | t ∈ F }. Il existerait ε > 0 tel qu’aucun élément de cet ensemble ne
2. On remarquera à cette occasion que, quand il ne l’oublie pas, l’auteur de ce texte ponctue même les équations
détachées, ce qui est aussi conforme au bon usage typographique, avec l’accentuation des majuscules, et quelques petits
autres traits telles les ligatures qui distinguent les belles éditions.
10 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

soit dans (ū, ū + ε). Dans ces conditions, ū + ε est lui-même un minorant des {u(t)}, contredisant
la définition de ū comme le plus grand minorant. Ainsi, pour tout n ∈ N, il existe tn tel que u(tn ) ∈
(ū, ū + 1/n), et les u(tn ) forment la suite recherchée.
La démonstration pour le cas d’une suite est identique.
Donc, l’inf peut être approché arbitrairement près par des éléments de l’ensemble considéré. Par
contre, il peut, suivant les cas appartenir ou ne pas appartenir à cet ensemble. Ainsi, par exemple,
inf t∈R {t2 } = 0, et cet inf est atteint en t = 0. On dit alors que c’est un min, et on écrit mint∈R {t2 } =
0. Par contre, inf t>0 {1/t} = 0, mais ici il n’y a pas de min. (C’est l’occasion de rappeler que +∞
n’est pas un nombre réel.) Les suites minimisantes sont toutes le suites tendant vers +∞.

Définition 1.8 (Minimum, minimum strict) S’il existe t∗ ∈ F tel que ∀t ∈ F , u(t) ≥ u(t∗ ), on dit
que c’est un minimum. Si u(t) > u(t∗ ) pour tout t 6= t∗ , le minimum est dit strict.

On a les mêmes définitions et mêmes remarques mutatis mutandis pour le sup, qui est un max s’il
est atteint.

Exercice 1.3 Trouver les infs suivants et dire s’il y a un min (l’un d’entre-eux est difficile) :

ln(n)
inf {ln(n)}, inf { }, sup{sin(n)}, sup{sin(r)}, sup{sin(t)} .
n∈N∗ n∈N∗ n n∈N r∈Q t∈R

Donnons enfin encore une définition utile

Définition 1.9 (Minimum local) Soit u(·) une fonction de R dans R. On dit que t∗ est un minimum
local sur F ⊂ Rn s’il existe un nombre ε > 0 tel que ∀t ∈ F ∩ (t∗ − ε, t∗ + ε), on a u(t) ≥ u(t∗ ).
On dit que ce minimum local est strict si c’est un minimum strict sur F ∩ (t∗ − ε, t∗ + ε)

1.2 Espace Rn
On appelle Rn l’ensemble des n-uplets de nombres (qu’on notera d’habitude en colonne) qu’on
appellera vecteurs :  
x1
 x2 
x= . .
 
 .. 
xn
On notera xt le transposé de x :

xt = ( x1 x2 · · · xn ) .

et de même, la transposée d’une matrice M sera notée M t .

1.2.1 Propriétés algébriques

Algèbre linéaire
Rn est de manière naturelle un espace vectoriel (addition et produit par un scalaire définis élément
à élément). L’étude des propriétés algébriques de la structure d’espace vectoriel nous emmènerait trop
1.2. ESPACE RN 11

loin. Nous y renonçons. Toutefois, la maı̂trise des rudiments de l’algèbre linéaire et du calcul matriciel
est nécessaire pour suivre ce cours.
Les concepts nécessaires comportent la notion de combinaison linéaire, de dépendance et d’in-
dépendance linéaire, d’opérateur linéaire, de matrice, de produit de matrices, de rang d’une matrice
et ses liens avec l’espace image et le noyau de l’application linéaire, de valeur propre et de vecteur
propre.
Le cours évitera la notion de valeur singulière d’une matrice et de décomposition associée. Mais
la connaissance des problèmes de conditionnement est indispensable pour faire du calcul numérique.
On rappelle ici une seule définition :

Définition 1.10 (Rayon spectral) On appelle rayon spectral d’une matrice le module de la valeur
propre de plus grand module,

et le théorème qui rend ce concept important :

Théorème 1.11 La suite des puissances d’une matrice carrée tend vers zéro si et seulement si son
rayon spectral est inférieur à un. Si le rayon spectral est supérieur à un, cette suite diverge.

Structures euclidienne et métrique

On supposera aussi connue la structure euclidienne de Rn , ce qui est un bien grand mot pour
désigner la notion de produit scalaire
n
X
(x, y) = x i yi ,
i=1

de norme
n
!1
2
1 X
kxk = (x, x) =2 x2i ,
i=1

l’inégalité triangulaire
kx + yk ≤ kxk + kyk
et la très utile inégalité de Cauchy-Schwarz

|(x, y)| ≤ kxkkyk ,

où l’inégalité est stricte si les vecteurs x et y ne sont pas colinéaires. (i.e. à moins qu’il n’existe un
nombre a tel que y = ax.)

1.2.2 Matrices positives

Ce paragraphe concerne les formes quadratiques de Rn . Il relève encore de la structure métrique
dans la mesure où la forme quadratique (x, Ax) introduite ci-dessous définit une norme, parfois notée
kxk2A dès que A est positive définie. Il est important pour un cours d’optimisation de connaı̂tre la
notion de matrice positive définie et semi-définie, et d’en connaı̂tre les principales propriétés. C’est
pourquoi nous isolons ces quelques résultats dans un paragraphe à part.
12 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Soient donc ci-dessous A, B,... des matrices carrées symétriques. On considèrera des formes qua-
dratiques de la forme
X n
xt Ax = (x, Ax) = aij xi xj
i,j=1

On fait remarquer au passage que comme xi xj = xj xi , si A n’était pas symétrique, la forme quadra-
tique ci-dessus ne dépendrait que de sa partie symétrique 21 (A + At )

Définition 1.11 (Matrice positive) La matrice A est dite positive définie, ce qu’on note A > 0, si
∀x 6= 0, xt Ax > 0. Elle est dite positive semi-définie, ce qu’on note A ≥ 0, si ∀x, xt Ax ≥ 0.

Remarquons que de telles matrices existent, par exemple l’identité est positive définie. Plus gé-
néralement, de manière évidente, si L est une matrice (rectangulaire) quelconque, Lt L est positive
semi-définie, et définie si L est injective (rang égal à son nombre de colonnes), et LLt est aussi
positive semi-définie, et définie si L est surjective (rang égal à son nombre de lignes.)
De façon très naturelle, on écrira que A > B si A − B > 0, et que A ≥ B si A − B ≥ 0.
On sait que les matrices symétriques admettent des valeurs propres réelles, et sont diagonalisables
(sur une base orthogonale). Soit
A = M t ΛM (1.1)
la forme diagonalisée de A, où Λ est la matrice diagonale des valeurs propres, et M −1 = M t parce
que M est orthogonale. On montre facilement la propriété suivante :

Proposition 1.12 Une matrice symétrique est positive semi-définie si et seulement si ses valeurs
propres sont toutes positives ou nulles. Elle est positive définie si et seulement si ses valeurs propres
sont toutes positives.

Une conséquence est qu’une matrice positive définie est inversible, et son inverse est encore posi-
tive définie. C’est, par contre, un exercice a priori non banal de démontrer le résultat suivant :

Proposition 1.13 Si A > B > 0, alors B −1 > A−1 > 0.

La forme (1.1) permet de montrer divers propriétés intéressantes.

Le changement de base pour passer sous la forme diagonale respectant les normes, on voit immé-
diatement que, si la plus grande valeur propre de A est µ, on a toujours

kAxk ≤ µkxk ,

(µ est la norme d’opérateur de A notée kAk) et donc que (x, Ax) ≤ µkxk2 , et aussi que si la plus
petite valeur propre est ν, on a
(x, Ax) ≥ νkxk2 .
On appelle ν la constante de coercivité de A. Cette dernière inégalité s’écrit aussi A ≥ νI. On a de
même A−1 ≤ 1/νI, de même que kA−1 k = 1/ν, car la forme diagonale montre que la plus grande
valeur propre de A−1 est 1/ν, et sa plus petite valeur propre 1/µ.
Remarquons aussi la propriété suivante :

Proposition 1.14 Si A et B sont des matrices carrées symétriques et 0 ≤ A ≤ B, alors kAk ≤ kBk.
1.2. ESPACE RN 13

En effet, soit x le vecteur propre associé à la plus grande valeur propre de A, de sorte que Ax =
kAkx. On a donc (x, Ax) = kAkkxk2 . Mais par hypothèse, (x, Bx) ≥ (x, Ax), ce qui en majorant
le premier terme par Cauchy Schwarz mène à kBkkxk2 ≥ kAkkxk2 , d’où le résultat annoncé.
En particulier, de A ≤ βI (où I est la matrice identité), on pourra conclure kAk ≤ β.
Appelons naturellement Λ1/2 la matrice diagonale des racines carrées des valeurs propres de A,
et posons
1
L = Λ2 M .
On a immédiatement
A = Lt L . (1.2)
Donc toute matrice positive semi-définie se met sous cette forme. On peut d’ailleurs remarquer que si
certaines valeurs propres sont nulles, la ligne correspondante de L est nulle, et peut être retirée, ce qui
fait de L une matrice rectangulaire, sans perdre la relation (1.2).
Une conséquence utile est la suivante :

Proposition 1.15 Si A ≥ 0 et xt Ax = 0, alors Ax = 0.

Définition 1.12 (Racine carrée) On notera

1 1
A 2 := M t Λ 2 M

Ainsi, A1/2 est symétrique, positive semi-définie, et satisfait (A1/2 )2 = A. Si A > 0, A1/2 > 0.

1.2.3 Propriétés topologiques

On dote naturellement Rn de la norme et de la distance euclidiennes :
n
!1
X 2

kxk = x2k , d(x, y) = ky − xk .

1

On va étendre à Rn les définitions et propriétés topologiques de R, simplement en remplaçant

partout les intervalles ouverts par des boules ouvertes :

Définition 1.13 (Boule ouverte) Étant donnés un vecteur x de Rn et un réel positif r, on appelle
boule ouverte de centre x et de rayon r, et on note B(x, r), l’ensemble des vecteurs y de Rn qui sont
à une distance inférieure à r de x :

B(x, r) = {y ∈ Rn | d(x, y) < r} .

Les propriétés (O) et (F) du paragraphe 1.1.2 servent de définition d’un ouvert et d’un fermé de
Rn . Ce qu’on formalise dans la définition ci-dessous :

Définition 1.14 (Ouverts et fermés de Rn ) La définition 1.5 est étendue à Rn en remplaçant dans la
propriété (O) l’intervalle ouvert de demi-longueur ε par une boule ouverte de rayon ε.

On garde aussi les mêmes définitions pour l’intérieur, l’adhérence, la frontière et l’extérieur d’un
ensemble. De même, les définitions de min et de min local s’étendent immédiatement.
On perd la caractérisation simple des ouverts. Mais on garde le théorème essentiel de Bolzano-
Weierstrass, i.e. la compacité des fermés bornés. Aussi simple qu’il soit, nous indiquons sa démons-
tration en un lemme et un théorème. On passe par les pavés :
14 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Définition 1.15 (Pavés) On appelle pavé fermé un sous-ensemble de Rn de la forme

P = {y ∈ Rn | yk ∈ [xk − rk , xk + rk ] , k = 1, . . . , n} ,

où x est un vecteur de Rn (appelé centre du pavé) et les rk sont des nombres positifs.

À lévidence, un pavé peut aussi être caractérisé par ses deux sommets extrêmes, deux vecteurs a
et b de Rn , et
P = {x ∈ Rn | ak ≤ xk ≤ bk , k = 1, . . . , n} .

Lemme 1.16 Les pavés fermés sont compacts.

Démonstration Étant donnée une suite {x(k)}k∈N , extraire une sous-suite telle que la première co-
ordonnée converge, de cette sous-suite extraire à nouveau une sous-suite telle que la deuxième coor-
donnée converge (dans cette nouvelle sous-suite, la première coordonnée garde la même limite que
dans la première sous-suite) et ainsi de suite n fois.

Théorème 1.17 (Compacts de Rn ) Les fermés bornés de Rn sont compacts.

Démonstration Comme le théorème 1.9, en remplaçant l’intervalle fermé par un pavé fermé.

1.3 Dérivation
1.3.1 Continuité
On rappelle pour mémoire la définition :

Définition 1.16 (Continuité) Une fonction u(·) est dite continue sur un ouvert Ω si pour tout x dans
Ω, et pour toute suite {xn } tendant vers x, u(xn ) tend vers u(x).

“Moralement”, une fonction continue est une fonction dont le graphe peut être dessiné “sans lever
le crayon”. Malheureusement, il peut y avoir des fonctions continues bien trop irrégulières pour qu’on
puisse en dessiner le graphe, puisque, par exemple, on connait des fonctions réelles continues partout
sur un intervalle mais dérivables nulle part. (C’est à dire que le “graphe” n’en admet nulle part de
tangente...)
Les propositions ci-dessous sont peut-être utiles à rappeler aussi :

Proposition 1.18 L’image inverse d’un ouvert par une fonction continue est un ouvert, l’image in-
verse d’un fermé par une fonction continue est un fermé.

On se rend facilement compte que la première propriété ci- dessus est équivalente à la définition
de la continuité. Elle est souvent prise comme définition. En fait, la deuxième propriété lui est trivia-
lement équivalente et pourrait donc aussi bien être retenue comme définition.

Proposition 1.19 L’image (directe) d’un compact par une fonction continue est un compact.

Il serait faux de croire que l’image directe d’un fermé par une fonction continue soit nécessairement
un fermé. Un contre-exemple est fourni par la fonction 1/x qui envoie [1, ∞) ⊂ R (qui est fermé) sur
(0, 1], qui n’est ni ouvert ni fermé.
1.3. DÉRIVATION 15

1.3.2 Dérivées et dérivées partielles

Soit d’abord u(·) : I ⊂ R → R une fonction réelle (scalaire) d’une variable réelle définie sur un
intervalle ouvert I de R. Le lecteur est supposé connaı̂tre la notion de dérivée. Nous choisissons de
réécrire sa définition sous la forme suivante :

Définition 1.17 (Dérivée) La dérivée de u(·) en t est, s’il existe, le nombre (unique : exercice) u0 (t)
tel que l’on ait
u(t + τ ) = u(t) + u0 (t)τ + o(τ ) , (1.3)
où o(·) est une fonction qui tend vers zéro plus vite que son argument. (C’est à dire que o(τ )/τ tend
vers zéro avec τ .)

Si u admet en tout point d’un intervalle ouvert I une dérivée u0 (t) qui est elle même une fonction
continue de t, on dira que “u est (de classe) C 1 ”, ou “u ∈ C 1 ”.
Passons à une fonction de plusieurs variables, u(·) : Ω ⊂ Rn → R d’un ouvert Ω de Rn dans R.
Le lecteur connaı̂t aussi la notion de dérivée partielle. Formellement, on peut écrire que ∂u/∂xi est la
dérivée de l’application partielle

xi 7→ u(x1 , x2 , . . . , xi , . . . , xn ).

∂u
La notation , que nous utiliserons, présente un gros inconvénient qu’il faut souligner ici. C’est
∂xi
l’usage du nom donné à l’argument (ici xi ) dans le nom de la fonction dérivée partielle. Ainsi, com-
ment doit-on écrire l’expression (1.3) pour évaluer l’accroissement de u entre les points
   
y1 y1 + z
 y2   y2 
 ..  et  ?
   
..
 .   . 
yn yn

Nous la noterons (en appelant e1 le vecteur (1, 0, . . . , 0)t )

∂u
u(y + ze1 ) = u(y) + (y)z + o(z) ,
∂x1
et surtout pas
∂u
u(y + ze1 ) = u(y) + (y)z + o(z) .
∂y1
Car on ne saurait changer le nom d’une fonction (ici la dérivée partielle par rapport à la première
variable) à chaque fois qu’on change le nom de son argument. Si on faisait ainsi, que deviendrait cette
dérivée partielle au point (1, 1, . . . , 1)t ?
Pour cette raison, Dieudonné propose de noter les dérivées partielles Di u(·) pour la dérivée par
rapport à la variable de rang i. En cas de besoin, on invite le lecteur à avoir recours à cette notation
qui évite les ambigüités.
Nous placerons d’habitude les vecteurs en colonne, et les dérivées partielles par rapport aux coor-
données d’un vecteur en ligne, réservant la notation u0 (x) à cette présentation :
 
u0 (x) = ∂x ∂u
1
∂u
∂x2
· · · ∂u
∂x n
.
16 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Ainsi nous aurons, utilisant un produit matriciel ordinaire, la formule fondamentale (1.3) préservée :

u(x + h) = u(x) + u0 (x)h + o(khk) . (1.4)

Nous utiliserons la notation ∇u pour désigner le vecteur colonne des dérivées partielles, donc le
transposé de u0 . Ainsi (1.4) s’écrit aussi

u(x + h) = u(x) + (∇u(x), h) + o(khk) ,

où (y, z) désigne le produit scalaire des vecteurs y et z.

Si la fonction u est elle même vectorielle :
 
u1 (x)
 u2 (x) 
u(x) =  ,
 
..
 . 
um (x)

on fera des u0i (x) les lignes de la matrice de type m × n :

 ∂u1 ∂u1 ∂u1
∂x2 · · ·

∂x1 ∂xn
 ∂u2 ∂u2 · · · ∂u2 
 ∂x1 ∂x2 ∂xn
u0 =  ,

..
 . 
∂um ∂um ∂um
∂x1 ∂x2 ··· ∂xn

et la formule fondamentale (1.4) reste correcte.

Dans la notation de Dieudonné,
 
D1 u1 D2 u1 · · · Dn u1
 D1 u2 D2 u2 · · · Dn u2 
u0 =  .
 
..
 . 
D1 um D2 um · · · Dn um

1.3.3 Dérivation en chaı̂ne et dérivées directionnelles

Dérivées en chaı̂ne Soit u(·) et v(·) deux fonctions réelles, et posons w(t) := u(v(t)). On note
d’habitude la fonction w comme u ◦ v. Nous supposons en outre que u et v sont de classe C 1 . Alors
w l’est aussi, et la formule (1.3) donne une façon facile de démontrer (exercice) la formule connue du
lecteur :
w0 (t) = u0 (v(t))v 0 (t) .
Remarquons qu’avec les conventions ci-dessus pour les fonctions vectorielles de variables vecto-
rielles, cette formule demeure. Si v : Rp → Rn et u : Rn → Rm , alors u0 est de type m × n et
v 0 de type n × p, de sorte que le produit matriciel peut bien être fait, et w0 est de type m × p. Tout ceci
découle nécessairement de la formule (1.4).
Dérivée directionnelle Appliquons cela à la dérivée de la “coupe” d’une fonction le long d’une droite.
Soit ξ(t) = x + th une droite de Rn . Ici, x et h sont fixés dans Rn , et t varie dans R. Posons
U (t) = u(ξ(t)). C’est donc la restriction de la fonction u à la droite de direction h passant par x. On
a en appliquant les règles ci-dessus :

U 0 (t) = u0 (ξ(t))h ,
1.3. DÉRIVATION 17

et nous l’utiliserons surtout en t = 0 :

U 0 (0) = (∇u(x), h) . (1.5)
Exercice 1.4 (important) : En déduire que
– la ligne de plus grande pente est parallèle à ∇u (et opposée pour descendre),
– cette direction est orthogonale à la courbe de niveau qui est la courbe {y | u(y) = u(x)}.
.

Différentielles Les règles de calcul précédentes justifient qu’il soit commode d’utiliser la notation
du = u0 (x)dx, puis si x lui-même est fonction de, disons, a, dx = x0 (a)da, d’où en reportant
du = u0 (x)x0 (a)da. Nous ne chercherons pas à justifier plus formellement l’usage de la notation
différentielle qui est très commode. Rappelons que la “rigueur mathématique” ne consiste pas à res-
pecter scrupuleusement des règles formelles de notation, elle consiste à ne pas se tromper. Il suffit
donc de n’utiliser ce type de calcul qu’à bon escient, et de revenir aux théorèmes prouvés quand on
n’est pas sûr de soi.
Intégration par parties Pour une fonction réelle u(·) dérivable en tout t ∈ I où I est un intervalle
ouvert, on a pour tout t1 et t2 de I la célèbre formule d’intégration
Z t2
u(t2 ) − u(t1 ) = u0 (t) dt . (1.6)
t1

Cette formule est moins naı̈ve qu’un long usage ne le fait penser. Elle suppose suffisamment de
régularité de la part de u.
Nous donnons ci-dessous une application importante de cette formule.
Développement au second ordre Soit u(·) une fonction réelle deux fois dérivable continument. (“De
classe C 2 ”.) Nous pouvons utiliser la formule (1.6) en y évaluant u0 (t) à l’aide de la formule (1.3)
appliquée à u0 . Il vient (exercice)
1
u(t + τ ) = u(t) + u0 (t)τ + u00 (t)τ 2 + o(τ 2 ) . (1.7)
2
En fait, on utilisera aussi une forme un peu différente des développements limités, avec le reste de
Lagrange, que nous rappelons ici :
1
u(t + τ ) = u(t) + u0 (t)τ + u00 (t + θτ )τ 2 , (1.8)
2
où θ est un nombre compris entre 0 et 1..
Exercice 1.5 Étendre ces formules à un ordre n quelconque.
Ces formules se généralisent à une fonction scalaire d’une variable vectorielle. Soit donc u(·) :
Ω ⊂ Rn → R. Soit D2 u la matrice symétrique de ses dérivées secondes. On a l’importante formule :
1
u(x + h) = u(x) + u0 (x)h + ht D2 u(x)h + o(khk2 ) , (1.9)
2
ou encore
1
u(x + h) = u(x) + u0 (x)h + ht D2 u(x + θh)h (1.10)
2
Comme toujours, ces formules se déduisent de celles du cas scalaire en regardant la fonction U (t) =
u(x + th).
Rien n’interdit d’étendre ce type de formule à des fonctions vectorielles (facile) et à des ordres
plus élevés. Il y faut des notations plus complexes...ou beaucoup de courage.
18 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

1.4 Existence, unicité, CNS, et toutes ces sortes de choses

Ici seulement commence un cours d’optimisation, fût-il élémentaire.

1.4.1 Existence et conditions

Nous rappelons le théorème fondamental suivant :

Théorème 1.20 (Weierstrass) Une fonction réelle continue sur un compact y atteint son minimum et
son maximum.

Démonstration Soit K un compact de Rn et u(·) une fonction continue de K dans R. Soit encore
u∗ = inf x∈K {u(t)}. Soit {xn } une suite minimisante, i.e. telle que u(xn ) → u∗ quand n → ∞.
Comme K est compact, il existe un point d’accumulation x∗ dans K et une sous-suite {xn0 } qui
converge vers x∗ . Comme {xn0 } est une sous-suite d’une suite minimisante, u(xn0 ) tend, comme
u(xn ), vers u∗ , mais comme u(·) est continue et que {xn0 } → x∗ , u(xn0 ) → u(x∗ ). La limite étant
unique, on en déduit que u(x∗ ) = u∗ . Donc la fonction atteint son minimum en x∗ . La preuve pour le
maximum se fait de la même façon. (Ou en appliquant le résultat juste démontré à −u(·).)
Ce théorème est bien “optimal” au sens où on pourra construire des contre-exemples (exercice)
en renonçant soit au caractère fermé de K, soit à son caractère borné, soit au caractère continu de
la fonction u. Par contre, on a défini des propriétés plus faibles que la continuité, comme la semi-
continuité inférieure ou supérieure, qui suffisent pour avoir l’existence d’un min ou d’un max respecti-
vement. Le contre-exemple construit en renonçant à la continuité ne sera simplement pas semi-continu
inférieurement.

Cas sans contrainte

Le théorème ci-dessus parle d’un compact, donc d’un ensemble fermé. Pourtant, la majeure partie
de l’analyse traditionnelle s’intéresse aux conditions (nécessaires ou suffisantes) prévalant en un mi-
nimum atteint dans un ouvert, ce qu’on appellera le cas sans contrainte. (On verra ci-dessous que la
considération de contraintes fait naturellement sortir de ce cadre.)
Les dérivées ont été inventées pour le théorème suivant :

Théorème 1.21 Soit u(·) une fonction dérivable d’un ouvert Ω dans R. Une condition nécessaire
pour qu’elle atteigne son minimum en x∗ ∈ Ω est que u0 (x∗ ) = 0.

Démonstration On utilise la formule (1.4) pour évaluer u au voisinage de x∗ . Si u0 (x∗ ) 6= 0, il

existe h tel que u0 (x∗ )h < 0 (par exemple h = −εu0 (x∗ )) et de module suffisamment petit pour que
|u0 (x∗ )h| > o(khk) (dans notre exemple, ε suffisamment petit positif), ceci par définition de o(·).
Ainsi, u(x∗ + h) < u(x∗ ), contredisant le fait que x∗ soit un minimum.
Chacun sait que si u0 = 0, on va voir la dérivée seconde pour décider si on a un minimum, un
maximum (local) ou un col. Ceci étant, il n’y a aucune nécessité à ce que la dérivée seconde soit
positive pour avoir un minimum, même strict. La fonction de R dans R u(t) = t4 à l’origine est un
contre-exemple évident. Contre-exemple plus intéressant : la fonction

exp (− t12 ) si t 6= 0 ,

u(t) =
0 si t = 0 ,
1.4. EXISTENCE, UNICITÉ, CNS, ET TOUTES CES SORTES DE CHOSES 19

est continue et infiniment dérivable en 0, où elle a un minimum strict. Pourtant toutes ses dérivées
sont nulles en 0.
Le théorème qu’on peut affirmer cependant est le suivant :
Théorème 1.22 Une condition nécessaire pour qu’une fonction u(·) deux fois continument différen-
tiable atteigne un minimum relatif en x∗ est que la matrice des dérivées secondes D2 u(x∗ ) soit posi-
tive semi définie.

Démonstration Si D2 u(x∗ ) n’est pas positive semi-définie, elle a au moins une valeur propre négative.
On peut donc trouver un vecteur h tel que ht D2 u(x∗ )h < 0. Il suffit alors d’exploiter la formule (1.9)
comme la formule (1.4) dans la condition nécéssaire du premier ordre avec ce vecteur d’accroisse-
ment.
On a enfin comme condition suffisante :
Théorème 1.23 Si u0 (x∗ ) = 0 et u00 (x∗ ) > 0 (dans le cas où x ∈ Rn , comprendre que la matrice
carrée symétrique D2 u(x∗ ) est positive définie), alors x∗ est un minimum local.

Démonstration Il suffit, à nouveau, d’utiliser la formule (1.9).

On étendra le théorème ci-dessus pour les fonctions de R à “première dérivée non nulle d’ordre
pair”. (En se souvenant que ce n’est qu’une condition suffisante de minimum local).

Contraintes
Terminons en examinant ce qu’on peut dire si la variable x, au lieu d’être libre de parcourir tout
Rn , est contrainte à rester dans un ensemble fermé. Commençons par examiner le cas d’une fonction
u(·) d’une seule variable réelle, considérée sur un intervalle fermé I = [a, b]. On peut affirmer de
façon évidente la proposition suivante :

Proposition 1.24 Si u atteint son minimum sur I en t∗ , alors, soit t∗ appartient à l’intérieur de I et
u0 (t∗ ) = 0, soit t∗ = a et u0 (t∗ ) ≥ 0, soit t∗ = b et u0 (t∗ ) ≤ 0.

On laisse le lecteur prouver ce résultat en détail.

Il s’agit de trouver une façon d’exprimer cela qui soit plus élégante, et surtout se généralise à
un fermé de Rn . Dans un premier temps, nous nous limitons à la construction ci-dessous. On verra
comment la simplifier dans le paragraphe sur la convexité.
Nous considérons donc le problème
min u(x)
x∈C
où C est un ensemble fermé donné.
En tout point x de C, on introduit les directions admissibles en x, définies de la façon suivante :
h ∈ Rn est une direction admissible en x si il existe ε > 0 tel que, pour tout t ≤ ε, x + th ∈ C. On
laisse le lecteur démontrer que si x est dans l’intérieur de C, toute direction (tout vecteur de Rn ) est
admissible en x. Par contre, il n’en va pas de même si x est sur la frontière de C.
On a le résultat suivant :
Théorème 1.25 : Si la fonction u(·) est dérivable et atteint son minimum sur C en x∗ ∈ C, nécessai-
rement, on a
(∇u(x∗ ), h) ≥ 0 (1.11)
pour toute direction h admissible en x∗ .
20 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Démonstration Il suffit de considérer la fonction d’une variable scalaire U (t) = u(x∗ + th). Pour
t positif suffisament petit, x∗ + th ∈ C, donc elle doit atteindre son minimum en zéro. On applique
alors la proposition ci-dessus, et la règle de dérivation en chaı̂ne.
◦
Remarquons qu’en fait, si x∗ est intérieur à C, alors u atteint son minimum sur C en x∗ , et il suffit
d’appliquer le théorème 1.21. On a gagné quelquechose si x∗ est un point frontière de C. Mais on n’a
pas à distinguer ces cas pour énoncer le théorème.

1.4.2 Multiplicateurs, dualité

Contraintes inégalité
Nous évoquons trop rapidement cette utilisation essentielle du théorème précédent. Nous n’en
tirerons d’algorithme numérique que dans le cas “convexe” (cf. paragraphe suivant).
Considérons le cas où l’ensemble C des x admissibles est donné par une contrainte scalaire
f (x) ≤ 0. Nous supposerons f dérivable, donc continue. Ainsi, C = {x | f (x) ≤ 0} est l’image
inverse du fermé (−∞, 0] par une fonction continue, donc fermé. Montrons le lemme suivant :

Lemme 1.26 Si x ∈ ∂C, toute direction h telle que (∇f (x), h) < 0, est admissible.

Démonstration Par hypothèse, x ∈ ∂C, soit f (x) = 0 (exercice). Il suffit alors d’utiliser (1.4).
On en déduit le résultat suivant, que nous étendrons ensuite au cas à plusieurs contraintes :

Proposition 1.27 Si x∗ fournit le minimum de u sous la contrainte f (x) ≤ 0, et si au minimum x∗ ,

∇f (x∗ ) 6= 0, nécessairement ∇u(x∗ ) = −p∇f (x∗ ) où p est un nombre positif ou nul, et pf (x∗ ) = 0.

Démonstration En effet, considérons d’abord le cas où on aurait f (x∗ ) < 0. Alors, f étant continue,
x∗ serait intérieur à C, et nécessairement, ∇u(x∗ ) = 0, soit la relation annoncée avec p = 0. Si
maintenant f (x∗ ) = 0, bien sûr cela est encore vérifié si ∇u(x∗ ) = 0. Si non, on déduit du théorème
général et du lemme ci-dessus qu’une condition nécessaire est que, pour tout h tel que (∇f (x∗ ), h) <
0, on ait (∇u(x∗ ), h) ≥ 0. Ceci est bien vérifié si ∇u(x∗ ) = −p∇f (x∗ ) avec p ≥ 0. Que cette
dernière relation soit aussi nécessaire découle du petit calcul suivant. Supposons que cette dernière
condition ne soit pas satisfaite. Alors choisissons

∇u(x∗ ) ∇f (x∗ )
h=− −
k∇u(x∗ )k k∇f (x∗ )k

qui est non nul d’après l’hypothèse même. D’après l’inégalité de Cauchy-Schwarz stricte, (∇f (x∗ ), h)
et (∇u(x∗ ), h) sont tous les deux strictement négatifs, violant donc la condition nécessaire. Enfin, il
découle de ce développement qu’on a toujours soit p = 0 soit f (x∗ ) = 0, et donc toujours pf (x∗ ) = 0.
Le résultat ci-dessus s’étend à m contraintes. Nous n’en donnerons pas la démonstration complète.
Le problème considéré est alors de trouver x∗ ∈ C tel que

min u(x) = u(x∗ )

x∈C

où
C = {x ∈ Rn | fi (x) ≤ 0 , i = 1, . . . , m}
1.4. EXISTENCE, UNICITÉ, CNS, ET TOUTES CES SORTES DE CHOSES 21

De même que nous avons du distinguer ci-dessus les cas f (x∗ ) < 0 et f (x∗ ) = 0, il nous faut
le faire ici contrainte par contrainte. C’est ce que nous faisons en introduisant la terminologie de
contrainte active en x∗ .
Soit donc I(x∗ ) l’ensemble des indices des contraintes actives en x∗ , c’est à dire I = {i ∈
[1 · · · m] | fi (x∗ ) = 0}. Les autres contraintes n’induisent aucune restriction sur les directions admis-
sibles, puisque les fj correspondantes restent négatives dans un voisinage de x∗ par continuité.

Théorème 1.28 Soient u(·) et fi (·), i = 1, . . . , m m + 1 fonctions (scalaires) dérivables. Si x∗

minimise u(·) sous les m contraintes fi (x) ≤ 0, et s’il existe au moins un vecteur h de Rn tel que
∀i ∈ I(x∗ ), (∇fi (x∗ ), h) < 0, alors nécessairement il existe m nombres pi positifs ou nuls tels que
m
X m
X
∇u(x∗ ) + pi ∇fi (x∗ ) = 0 , et pi fi (x∗ ) = 0 .
i=1 i=1

Démonstration Toute direction h de Rn qui satisfait (∇fi (x∗ ), h) < 0 pour toutes les contraintes
actives est admissible (par une extension banale du lemme 1.26). L’affirmation qui nous manque est
donc la suivante, qui généralise le petit calcul effectué dans la proposition :

Lemme 1.29 (Farkas) Soient g0 et gi , i = 1 . . . r des vecteurs de Rn . Toute direction h satisfaisant

(gi , h) < 0, i = 1 . . . r satisfait aussi (g0 , h) ≥ 0 si et seulement si il existe r nombres positifs ou nuls
pi tels que
Xr
g0 = − pi gi .
i=1

Que cela soit suffisant est évident. Le lemme affirme le caractère nécessaire.
Il ne reste plus qu’à appliquer 1.11, et les deux précédents lemmes, (le dernier avec g0 = ∇u(x∗ )
et gi = ∇fi (x∗ ) pour i ∈ I) et à prendre pj = 0 pour j ∈ / I pour obtenir le résultat annoncé. Notons
aussi qu’on a bien par construction, pour tout i, soit fi (x∗ ) = 0 (cas où i ∈ I) soit pi = 0 (cas
contraire).
Remarquons enfin que comme nécessairement, pour tout i, fi (x∗ ) ≤ 0 (x∗ est admissible) et pi ≥
0, chaque produit pi fi (x∗ ) est négatif ou nul. Affirmer que leur somme est nulle est donc équivalent à
affirmer que chacun est nul, ce qui est connu sous le nom de propriété des écarts complémentaires.

Contraintes égalité

On a un résultat analogue pour les contraintes égalité. Nous le citons car il est très important. Sa
démonstration dépend du théorème des fonctions implicites, nous ne l’évoquons pas.

Théorème 1.30 (Multiplicateurs de Lagrange) Si x∗ minimise la fonction dérivable u(·) sous les
contraintes fi (x) = 0, i = 1, . . . , m, et si les vecteurs ∇fi (x∗ ) sont linéairement indépendants (le
jacobien f 0 de f est sujectif en x∗ ), il existe m nombres λi , i = 1, . . . , m, tels que
m
X
∇u(x∗ ) + λi ∇fi (x∗ ) = 0 . (1.12)
i=1
22 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Ce théorème ressemble au précédent, aux notations près : nous avons noté les multiplicateurs λ
au lieu de p. La grande différence est qu’on ne sait rien sur le signe des multiplicateurs λi . C’est une
constante en théorie de la dualité : les multiplicateurs (ou variables duales, ou variables adjointes) sont
signés pour des contraintes inégalité, non signés pour des contraintes égalité. On invite le lecteur à
deviner le théorème qui doit être vrai quand on a des contraintes égalité et des contraintes inégalité
dans le même problème.
La façon classique d’utiliser ce théorème est la suivante. Les n équations (1.12) (une pour chaque
dérivée partielle) permettent, quand tout va bien, de calculer x en fonction des m inconnues λ. En
reportant dans les m équations f (x(λ)) = 0, encore si tout va bien, on peut espérer trouver les λ qui
conviennent.
Ceci n’est qu’une façon de souligner qu’il y a n + m inconnues (x, λ) pour n + m équations. On
écrit souvent ces conditions à l’aide du lagrangien

L(x, λ) := u(x) + (λ, f (x))

sous la forme

∂L
= 0,
∂x
∂L
= 0.
∂λ

Une façon algorithmique d’utiliser cette remarque est de résoudre ce système par la méthode de New-
ton. Cela conduit à l’algorithme de programmation quadratique séquentielle du paragraphe 3.4.2.
Avant de conclure, on donne une forme alternative du théorème de Lagrange :

Théorème 1.30 (Multiplicateurs de Lagrange, 2me forme) Si x∗ minimise la fonction dérivable u(·)
sous les contraintes fi (x) = 0, i = 1, . . . , m, il existe m + 1 nombres λi , i = 0, 1, . . . , m non tous
nuls tels que
Xm
∗
λ0 ∇u(x ) + λi ∇fi (x∗ ) = 0 .
i=1

1.4.3 Convexité
Nous donnons ici une présentation étriquée de quelques résultats de la théorie de la convexité. La
théorie (plus ou moins) complète remplit des livres entiers, et concerne l’analyse des fonctions non
différentiables. C’est dire combien notre définition en termes de dérivée seconde est restrictive.

Fonction convexe

Définition 1.18 (Fonction convexe) Une fonction u(·) deux fois différentiable sera dite convexe si
pour tout x dans le domaine considéré, u00 (x) ≥ 0 pour une fonction de R, ou D2 u(x) ≥ 0 pour une
fonction de Rn . Elle sera dite strictement convexe si les inégalités ci-dessus sont strictes.

Théorème 1.31 Si u(·) est une fonction convexe, alors u0 (x∗ ) = 0 implique que x∗ est un minimum
(la condition nécessaire devient aussi suffisante). Si u est strictement convexe, ce minimum est strict.
1.4. EXISTENCE, UNICITÉ, CNS, ET TOUTES CES SORTES DE CHOSES 23

Démonstration Prenons d’abord le cas d’une fonction de R dans R, et comme de coutume notons t
sa variable, t∗ le point de dérivée nulle : u0 (t∗ ) = 0. Utilisons la formule (1.6) pour évaluer u0 (t). Du
fait que u0 (t∗ ) = 0, on a
Z t
0
u (t) = u00 (s) ds .
t∗
Comme u00
est positive ou nulle pour tout s, (ou positive) il en résulte que u0 (t) est négative ou nulle
(ou négative) pour t ≤ t∗ , et positive ou nulle (ou positive) pour t ≥ t∗ . Utilisons à nouveau la formule
(1.6), mais pour évaluer u cette fois :
Z τ
∗
u(τ ) = u(t ) + u0 (t) dt .
t∗

À nouveau, le signe de u0 (t) nous montre que l’intégrale du membre de gauche ci-dessus est toujours
positive ou nulle, ou strictement positive si u est strictement convexe. Donc u(τ ) ≥ u(t∗ ) —ou
u(τ ) > u(t∗ ) si u est strictement convexe—, ce qu’il fallait démontrer.
Prenons maintenant le cas d’une fonction de Rn . Admettons que le domaine où nous la considérons
est tel que x∗ peut être joint à x par un segment de droite. Pour x fixé, posons donc

ξ(t) = x∗ + t(x − x∗ ) .

On pose U (t) = u(ξ(t)). On se souvient que U 0 (t) = u0 (ξ(t))(x − x∗ ) (et donc que U 0 (0) = 0), et
que
U 00 (t) = (x − x∗)t D2 u(ξ(t))(x − x∗ ) .
Comme D2 u est une matrice positive semi-définie (ou définie), on en déduit que U est convexe (ou
strictement convexe) et on peut lui appliquer la démonstration précédente.

Théorème 1.32 Si u(·) est une fonction convexe, pour tout x et y de Rn on a :

u(y) ≥ u(x) + u0 (x)(y − x) (1.13)

Démonstration Il suffit de remarquer que la fonction

y 7→ u(y) − u(x) − u0 (x)(y − x)

à la même dérivée seconde que y 7→ u(y). Elle est donc aussi convexe. Sa dérivée en y = x est nulle,
elle y est donc minimum par application du théorème précédent, mais elle y est nulle. D’où le résultat
annoncé.
Ce résultat dit simplement que le graphe de la fonction est “au dessus” du plan tangent en un point
quelconque du graphe.

Définition 1.19 (Fonction α-convexe) Une fonction de Rn dans R deux fois continument dérivable
est dite α-convexe (ou coercive de coefficient de coercivité α) où α est un nombre positif, si sa dérivée
seconde est partout supérieure ou égale à α. (Comprendre, dans le cas n > 1, D2 u(x) − αI ≥ 0.)

Théorème 1.33 Pour une fonction u(·) α-convexe, on a pour tout x et y de Rn

[u0 (y) − u0 (x)](y − x) ≥ αky − xk2 , (1.14)

et aussi
α
u(y) ≥ u(x) + u0 (x)(y − x) + ky − xk2 . (1.15)
2
24 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Démonstration Il suffit à nouveau de considérer la fonction U (t) = u(x + t(y − x)), dont la dérivée
seconde est U 00 (t) = ((y − x), u00 (x + t(y − x))(y − x)) ≥ αky − xk2 puis d’utiliser cette minoration
pour minorer U 0 (t) :
Z t
0 0
U (t) = U (0) + U 00 (s) ds ≥ U 0 (0) + tαky − xk2 .
0

L’inéquation, prise avec t = 1 et en utilisant (1.5), est directement (1.14), et en reportant cette mino-
ration dans Z 1
α
U (1) = U (0) + U 0 (t) dt ≥ U (0) + U 0 (0) + ky − xk2
0 2
on obtient (1.15).
Ceci dit que le graphe est non seulement au-dessus du plan tangent (qui est une fonction de dérivée
seconde nulle) mais aussi au-dessus du “paraboloı̈de” tangent au graphe de dérivée seconde α. Une
conséquence importante de cette propriété est la suivante :

Corollaire 1.34 Une fonction u(x) (deux fois différentiable) α-convexe tend vers l’infini quand x tend
vers l’infini, et atteint son minimum sur Rn .

Démonstration exercice (Pour la deuxième assertion, on se ramènera au théorème de Weierstrass.)

Et les inégalités (1.14) et (1.15) ont les importantes conséquences suivantes :

Corollaire 1.35 Si la fonction u(·) est α-convexe et atteint son minimum sur Rn en x∗ , on a pour tout
x ∈ Rn :
k∇u(x)k ≥ αkx − x∗ k (1.16)
et
α
u(x) ≥ u(x∗ ) + kx − x∗ k2 (1.17)
2

Démonstration On utilise (1.14) entre x et x∗ en utilisant le fait que ∇u(x∗ ) = 0 :

(∇u(x), x − x∗ ) ≥ αkx − x∗ k2 ,

et on majore le membre de gauche à l’aide de l’inégalité de Cauchy-Schwarz :

k∇u(x)kkx − x∗ k ≥ αkx − x∗ k2 .

Si kx − x∗ k =
6 0, on simplifie par ce nombre pour obtenir (1.16), qui reste a fortiori vrai si kx − x∗ k =
0. Quant à (1.17), c’est simplement (1.15) avec ∇u(x∗ ) = 0.

Ensemble convexe
(C’est un barbarisme de traiter les ensembles convexes après les fonctions convexes)

Définition 1.20 (Ensemble convexe) Un ensemble est dit convexe si chaque fois qu’il contient deux
points x1 et x2 , il contient la corde qui les joint, i.e. les points {ξ(t) = x1 + t(x2 − x1 ) , t ∈ [0, 1]}.

Nous ne retiendrons, sans démonstration, que deux propriétés des ensembles convexes :
1.4. EXISTENCE, UNICITÉ, CNS, ET TOUTES CES SORTES DE CHOSES 25

Proposition 1.36 Soit C un sous-ensemble convexe d’intérieur non vide de Rn . En tout point x̄ de sa
frontière ∂C, il existe au moins une normale extérieure, c’est à dire un vecteur non nul ν ∈ Rn tel
que,
∀x ∈ C̄, (ν, x − x̄) ≤ 0 .
C’est à dire que tous les points de C sont “de l’autre côté” du plan orthogonal à ν passant par x̄.
Proposition 1.37 Tout point x a une projection x̂ (notée PC (x)) sur C̄ qui est le point de C̄ le plus
proche de x. En outre, si x est extérieur à C, sa projection x̂ appartient à ∂C. Dans ce cas x − x̂ est
une normale extérieure à C en x̂, et cette dernière propriété caractérise la projection.
Notons qu’à l’évidence, si x ∈ C̄, d’après la définition même de la projection sur C, il est sa
propre projection.
Et enfin un théorème pour montrer comment peuvent être utilisés ces concepts de nature géomé-
trique. Pour simplifier l’écriture, nous supposons que C est fermé, ce qui nous dispense d’écrire à
chaque fois que ce serait justifié C̄.
Théorème 1.38 Soit C un ensemble convexe fermé de Rn . La projection sur C est une contraction
au sens large, c’est à dire que quelques soient x1 et x2 , on a kPC (x1 ) − PC (x2 )k ≤ kx1 − x2 k.

Démonstration Notons d’abord qu’en combinant la proposition 1.37 et la définition de la normale,

on a pour tout x extérieur à C, et tout y dans C
(x − PC (x), y − PC (x)) ≤ 0 ,
et cette relation reste vraie si x ∈ C̄, car alors x − PC (x) = 0.
On écrit donc successivement cette relation en prenant l’un des xi pour x et la projection x̂j :=
PC (xj ) de l’autre pour y :
(x1 − x̂1 , x̂2 − x̂1 ) ≤ 0 ,
(x2 − x̂2 , x̂1 − x̂2 ) ≤ 0 .
On change le sens de l’inégalité en changeant le signe d’un des facteurs du produit scalaire dans
chacune des deux lignes ci- dessus, en s’arrangeant pour faire apparaı̂tre le même facteur les deux
fois :
(x1 − x̂1 , x̂1 − x̂2 ) ≥ 0
(x̂2 − x2 , x̂1 − x̂2 ) ≥ 0
et on additionne, pour obtenir
(x1 − x2 − (x̂1 − x̂2 ), x̂1 − x̂2 ) ≥ 0 ,
soit encore
(x1 − x2 , x̂1 − x̂2 ) ≥ kx̂1 − x̂2 k2 .
On majore alors le premier produit scalaire à l’aide de l’inégalité de Cauchy-Schwarz pour obtenir
kx1 − x2 kkx̂1 − x̂2 k ≥ kx̂1 − x̂2 k2 ,
et en simplifiant une fois par kx̂1 − x̂2 k, qu’on peut supposer non nul sans quoi la propriété est
trivialement établie, le résultat annoncé.
Il est utile pour l’intuition de remarquer que c’est de la géométrie euclidienne que nous avons fait
là. Ce résultat servira dans la preuve de convergence d’un des principaux algorithmes d’optimisation
à venir.
26 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE

Inégalité d’Euler, théorème de Kuhn-Tucker

On va maintenant simplifier l’inégalité (1.11) en la spécialisant au cas où l’ensemble des contraintes
C est convexe.

Théorème 1.39 (Inégalité d’Euler) Soit C un convexe fermé de Rn et u(·) une fonction dérivable de
C dans R. Si u atteint son minimum sur C en x∗ , alors,

∀x ∈ C, (u0 (x∗ ), x − x∗ ) ≥ 0 . (1.18)

Démonstration Comme C est supposé convexe, et que x∗ et x appartiennent tous les deux à C, il
en va de même de tous les points du segment [x∗ , x], qui est donc une direction admissible. Il suffit
d’appliquer (1.11)
Nous en déduisons l’important corollaire suivant :

Corollaire 1.40 Pour une fonction convexe, la condition (1.18) est nécessaire et suffisante pour que
x∗ soit son minimum sur le convexe C. En outre, si elle est α-convexe, ce minimum est unique.

Démonstration Il suffit d’uiliser la relation (1.13) pour la première assertion, et (1.15) pour la deuxième.

Théorème 1.41 (Kuhn et Tucker) Soit u(·) une fonction convexe d’un convexe fermé C de Rn dans
R. Soient fi , i = 1 . . . m m fonctions convexes de C dans R. On suppose (hypothèse de Slater) qu’il
existe x0 ∈ C tel que fi (x0 ) < 0, i = 1, . . . m. Alors si u → ∞ quand kxk → ∞ dans l’ensemble
C ∩ {x | fi (x) ≤ 0, i = 1, . . . m} (ou si cet ensemble est borné), u y admet ∗
∗
Pmun minimum x , et il
existe m nombres positifs ou nuls pi , i = 1, . . . m tels que (on note (p, f ) = i=1 pi fi )

∀p ≥ 0, ∀x ∈ C ,
(1.19)
u(x∗ ) + (p, f (x∗ )) ≤ u(x∗ ) + (p∗ , f (x∗ )) ≤ u(x) + (p∗ , f (x))

Réciproquement, s’il existe x∗ ∈ C et p∗ ≥ 0 satisfaisant (1.19), alors fi (x∗ ) ≤ 0 pour i = 1, . . . , m

et x∗ est un minimum de u sur le domaine considéré.

Commentaires Avant de démontrer ce théorème, deux commentaires s’imposent.

1. Les inéquations (1.19) peuvent aussi s’écrire en termes du Lagrangien L(x, p) = u(x) + (p, f )
sous la forme suivante, où on note P + le cône des éléments à coordonnées positives de Rn :

∀(x, p) ∈ C × P + , L(x∗ , p) ≤ L(x∗ , p∗ ) ≤ L(x, p∗ ) .

On dit que (x∗ , p∗ ) forme un point selle sur C × P + .

2. Comme précédemment, le produit scalaire (p, f (x∗ )) est toujours négatif ou nul, puisque p est
à éléments positifs, et les fi (x∗ ) sont négatifs ou nuls. On va très vite voir que l’inégalité de
gauche n’est que la condition des écarts complémentaires (p∗ , f (x∗ )) = 0, qui dit que chacun
des p∗i fi (x∗ ) = 0, les contraintes non actives en x∗ n’entrent pas dans le lagrangien (y ont un
coefficient p∗i nul).
Démonstration Nous ne démontrerons ce théorème que dans le cas simplifié où C est tout Rn . C’est
aussi dans ce cas que l’algorithme que nous en déduirons, l’algorithme d’UZAWA, est praticable.
Le domaine considéré, {x | fi (x) ≤ 0 , i = 1, . . . m}, est convexe fermé. Du fait qu’il est fermé,
si les suites minimisantes y sont bornées, elles ont un point d’accumulation qui est un minimum de u.
1.4. EXISTENCE, UNICITÉ, CNS, ET TOUTES CES SORTES DE CHOSES 27

Du fait qu’il est convexe, la condition de Slater implique que x0 − x∗ est une direction admissible en
x∗ . Nous disposons donc déjà du théorème (1.28). Notons en outre que u étant convexe, les fi aussi, et
les p∗i étant positifs ou nuls, le Lagrangien L est convexe en x. Donc, par le théorème (1.31), x∗ est un
minimum de L(x, p∗ ). Ceci donne l’inégalité de droite du point selle. L’inégalité de gauche découle
immédiatement de la condition des écarts complémentaires qui fait aussi partie du théorème (1.28).
Réciproquement, supposons les inégalités du point selle vérifiées. L’inégalité de gauche donne

∀p ∈ P + (p, f (x∗ )) ≤ (p∗ , f (x∗ )) .

Si les fi (x∗ ) n’étaient pas tous négatifs ou nuls, disons que fj (x∗ ) > 0, en faisant tendre p∗j seul
vers l’infini, on violerait sûrement cette inégalité. Donc aussi (p, f (x∗ )) ≤ 0 pour tout p dans P + .
Si (p∗ , f (x∗ )) était non nul, donc négatif, en remplaçant p∗ par p∗ /2 on violerait à nouveau cette
inégalité. Donc l’inégalité de gauche implique à la fois le caractère admissible de x∗ et la condition
des écarts complémentaires.
Du coup, l’inégalité de droite s’écrit

∀x u(x∗ ) ≤ u(x) + (p∗ , f (x)) .

Considérons un x admissible, c’est à dire tel que fi (x) ≤ 0 pour tout i. Alors, (p∗ , f (x)) ≤ 0, et donc
l’inégalité ci-dessus implique a fortiori u(x∗ ) ≤ u(x), ce qu’il fallait démontrer.
L’interprétation économique de ce théorème est célèbre. Si les contraintes fi (x) ≤ 0 représentent
des ressources disponibles en quantités limitées (une par indice), le lagrangien est ce qu’il faut faire
payer à l’utilisateur, les pi représentant le “prix” des ressources rares. Ce que dit le théorème est que
pour un jeu de prix adéquate, p∗ , la minimisation, sans tenir compte des contraintes de ressources, de
sa fonction de coût ainsi augmentée, mène l’utilisateur à une décision optimale en tenant compte des
contraintes.
28 CHAPITRE 1. PROLOGUE : RAPPELS D’ANALYSE
Chapitre 2

Recherche unidimensionnelle

2.1 Introduction
2.1.1 Objectif
Ce bref chapitre examine une question qui pourrait sembler bien naı̈ve. Elle mérite qu’on s’y
attarde un peu parce qu’elle intervient comme technique intermédiaire dans de nombreux algorithmes
que nous découvrirons ensuite, c’est la “boucle intérieure” de boucles imbriquées, donc celle qu’il
faut soigner le plus.
Soit donc u(·) une fonction d’une seule variable réelle (nous disons “de R dans R”), dont nous
recherchons le minimum sur un intervalle [a, b]. Nous supposons
– que le minimum t∗ recherché est à l’intérieur du segment
– que la fonction u est unimodale sur le segment, c’est à dire d’abord décroissante jusqu’à t∗ ,
puis croissante.
Dans la pratique, choisir a et b pourra poser des problèmes, que nous n’examinons pas ici.
Le problème est de déterminer des algorithmes permettant de trouver t∗ efficacement, c’est à dire
avec une bonne précision mais “pas trop” de calculs.

2.1.2 Pente et dérivée numérique

Suivant les algorithmes proposés, on peut ou non avoir besoin de la “pente” de la fonction, ou
plus précisément du signe ou de la valeur de sa dérivée. Dans bien des cas, la dérivée est aussi facile à
calculer que la fonction elle-même, et ceci ne pose pas de problème. Mais il y a aussi des problèmes
pour lesquels le calcul de la dérivée demande significativement plus de calculs que celui de la fonction.
Remarquons que pour une fonction de Rn dans R, la dérivée consiste en n nombres, soit n fois plus
que la fonction.
Il y a aussi des situations pour lesquelles la dérivée en tant que telle n’est pas connue. Typique-
ment, la fonction u peut n’être donnée que comme un gros programme informatique auquel on fournit
t et qui rend u(t). Il est utile de savoir que sont en train d’apparaı̂tre des outils informatiques —
tels Odyssée— qui prennent en entrée le source d’un programme (FORTRAN pour Odyssée) et
produisent en sortie le source (dans le même langage) d’un programme qui calcule les dérivées par-
tielles de la fonction que définit le programme donné. Mais ces outils sont encore du domaine de la
recherche, ne marchent que sous certaines conditions sur la façon dont est écrit le programme donné,
etc. Supposons, par exemple, que la fonction u soit calculée à l’aide de fonctions intermédiaires ta-
bulées et non disponibles sous forme de combinaison d’opérations et de fonctions “élémentaires”. Il

29
30 CHAPITRE 2. RECHERCHE UNIDIMENSIONNELLE

n’y a aucune chance qu’on puisse en calculer formellement la dérivée.

Dans ces cas, on peut avoir recours à la “dérivation numérique”, un grand mot pour quelque chose
de bien élémentaire. Il s’agit tout simplement d’approximer u0 (t) par (u(t+δ)−u(t))/δ. La difficulté
est de choisir δ, assez petit pour que ceci soit une approximation “raisonnable” de la dérivée, mais
assez grand pour que la différence u(t + δ) − u(t) soit calculée de façon significative. On voit que le
bon choix de δ dépend de la précision avec laquelle sont faits les calculs. Il faut parfois tatonner pour
choisir ce paramètre.
Chaque fois qu’on veut seulement le signe de la dérivée (la réponse à la question “la fonction est-
elle croissante ou décroissante en t ?”), pour savoir de quel côté du minimum on se trouve, le “risque”
est que ce minimum soit entre t et t + δ , et qu’on ne le remarque pas. Notons à ce propos qu’il
est de toutes façons illusoire de chercher le minimum t∗ avec une précision meilleure que celle avec
laquelle on sait distinguer deux valeurs de t par la valeur qu’elles donnent à u(t). (Sauf si la dérivée
est calculable, et avec une meilleure précision que la fonction, cas tout à fait inhabituel.) Mais ceci
impose de choisir δ suffisamment petit.

2.2 Méthodes directes

On appelle “méthodes directes” celles qui ne font pas appel au calcul de la dérivée. Nous étendrons
cela à celles qui ne demandent que le signe de la dérivée.

2.2.1 Dichotomie
La méthode la plus simple, et pas forcément la plus mauvaise, consiste moralement à résoudre
u0= 0 par dichotomie. On arrive ainsi à l’algorithme suivant :
Algorithme Dichotomie
1. a0 := a, b0 := b
2. n = 1
3. m := (an−1 + bn−1 )/2
4. Évaluer u0 (m)
5. si u0 (m) < 0 , an := m , bn := bn−1 ,
si u0 (m) > 0 , an := an−1 , bn := m,
6. Incrémenter n de 1 et retourner au pas 3
On a omis dans l’algorithme ci-dessus le test d’arrêt, qui est un ingrédient nécessaire de tout
algorithme. C’est qu’ici, on sait exactement la précision obtenue au pas n : on peut affirmer que
t∗ ∈ [an , bn ], donc on a une précision de (b − a)/2n . On peut donc déterminer a priori le nombre
de pas à effectuer en fonction de la précision souhaitée. On comparera donc n à ce nombre avant de
l’incrémenter.
On a choisi, dans la description de l’algorithme ci-dessus, une version compréhensible du point de
vue mathématique, et les mathématiques n’aiment pas donner le même nom à des variables différentes.
Un informaticien aurait écrit l’algorithme de la façon plus économe suivante :
Algorithme Faire N fois
1. m := (a + b)/2
2. Si u0 (m) < 0 a := m,
si u0 (m) > 0 b := m
2.2. MÉTHODES DIRECTES 31

3. Recommencer au début
Où N est choisi en fonction de la précision souhaitée. Si cette précision est ε, N croit comme le
logarithme (à base 2) de 1/ε (exercice : calculer N ). Chaque pas demande un calcul de u0 , ou, si on
doit utiliser des dérivées numériques, deux calculs de u. Le nombre de calculs de u à effectuer croit
donc comme 2 log2 (1/ε) = 2 ln(1/ε)/ ln 2.

2.2.2 Suites de Fibonacci

Par une méthode directe, on peut faire mieux que la dichotomie, soit une croissance moins rapide
que 2 log2 (1/ε). Le principe de la méthode, dite des suites de Fibonacci (on verra pourquoi), est le
suivant.
On part du segment [a0 , b0 ], et on suppose qu’on a calculé u(a0 ) et u(b0 ). On choisit deux point
intérieurs c0 < d0 . On calcule encore u(c0 ) et u(d0 ). La proposition de base est la suivante :

Proposition 2.1 Si u(c0 ) < u(d0 ), alors t∗ ∈ [a0 , d0 ], si au contraire u(c0 ) > u(d0 ), alors t∗ ∈
[c0 , b0 ].

En effet, en parcourant le segment [a0 , b0 ], dans le premier cas, la fonction u a nécessairement

commencé à croitre avant d0 , et donc t∗ < d0 , et au contraire dans le deuxième cas, elle a continué à
décroitre audelà de c0 , donc t∗ > c0 .
L’algorithme s’en déduit dans son principe : prendre pour [a1 , b1 ] le nouveau segment contenant
sûrement t∗ , [a0 , d0 ] ou [c0 , b0 ] suivant le cas, et recommencer. Mais au pas suivant, on connait déjà u
en un point intérieur, c0 dans le premier cas et d0 dans le deuxième. Donc on n’aura plus qu’à choisir
un seul autre point intérieur, et calculer u une seule fois, pour itérer.
La difficulté qui subsiste est dans le choix judicieux des points intérieurs à chaque pas. Le premier
choix, naturel, est d’imposer qu’ils soient symétriques par rapport au milieu du segment, ainsi la
longueur du segment restant après l’évaluation de u et application de la proposition sera indépendante
du résultat du test. En outre, le calcul du deuxième nouveau point intérieur à chaque pas est alors
trivial, puisqu’il suffit de prendre le symétrique de celui déjà connu.
Mais cette politique fait dépendre toute la suite des points utilisés du choix des deux premiers (en
fait d’un des deux premiers, l’autre étant fixé par symétrie), et cette dépendance peut mener à des
blocages. Par exemple, le point intérieur connu au pas n (soit cn−1 ou dn−1 suivant le cas) pourrait se
retrouver au milieu, ou tout près du milieu, du segment, une situation qui empêche d’appliquer notre
méthode.
Pour analyser cette question, il faut rentrer un peu dans le détail de cette suite de points.
On supposera qu’on s’est arrangé pour qu’à chaque pas, cn soit plus grand que le milieu (an +
dn )/2 de [an , dn ], et dn soit plus petit que le milieu (cn + bn )/2 de [cn , bn ]. Ainsi, l’algorithme peut
être précisé ainsi :
Algorithme Fibonacci (pas n + 1)
1. Évaluer u au point intérieur manquant (symétrique de celui déjà connu),
2. si u(cn ) < u(dn ),
faire an+1 := an , bn+1 := dn , dn+1 := cn , cn+1 := an+1 + bn+1 − dn+1 ,
si u(cn ) > u(dn ),
faire an+1 := cn , bn+1 := bn , cn+1 := dn , dn+1 := an+1 + bn+1 − cn+1 .
Appelons l0 := b0 − a0 la longueur du segment initial, et de même ln la longueur du segment
[an , bn ]. Du fait de la symétrie, nous avons déjà remarqué que la longueur des segments successifs
32 CHAPITRE 2. RECHERCHE UNIDIMENSIONNELLE

ne dépend pas du résultat du test. Examinons donc deux pas consécutifs en supposant que u(cn−1 ) <
u(dn−1 ) et u(cn ) < u(dn ). Ainsi, [an , bn ] = [an−1 , dn−1 ], et dn = cn−1 , puis [an+1 , bn+1 ] =
[an , dn ] = [an−1 , cn−1 ]. Ainsi, ln+1 = cn−1 − an−1 , tandis que ln = dn−1 − an−1 qui, par symétrie,
est aussi ln = bn−1 − cn−1 . On est donc arrivé à la relation de récurence fondamentale de cette étude :

ln−1 = ln + ln+1 . (2.1)

En faisant tourner cette récurence à l’envers, c’est à dire en posant fk = lN −k pour un nombre N
suffisamment grand, on voit qu’on arrive à la récurrence

fk+1 = fk + fk−1 , (2.2)

qui initialisée en f0 = 0, f1 = 1 donne la suite des nombres de Fibonacci. (exercice : faire un peu de
bibliographie pour retrouver qui était Fibonacci, et pourquoi il s’est intéressé à cette suite.)
On peut facilement calculer les premiers termes de la suite de Fibonacci :

f0 = 0
f1 = 1
f2 = 1
f3 = 2
f4 = 3
f5 = 5
..
.
f16 = 987
f21 = 10946
f26 = 121393

À l’évidence, la récurrence de Fibonacci génère une suite de nombres entiers positifs qui croit
√ très vite
et tends vers l’infini. On peut montrer, et c’est important, qu’elle croit comme ρk+ où ρ+ = ( 5+1)/2
est connu comme le nombre d’or, et est la racine positive de l’équation caractéristique de la récurrence
(2.2) :
ρ2 = ρ + 1 .
√
On aura aussi besoin de r+ = 1/ρ+ = ( 5 − 1)/2, qui satisfait, lui

r2 = 1 − r ,

√ ln+1 = ln−1 − ln , et des

dont on voit qu’elle est l’équation caractéristique de√la récurence (2.1) écrite
racines négatives des mêmes équations ρ− = (1 − 5)/2 et r− = (−1 − 5)/2.
La politique “optimale” consiste à avoir au dernier pas une distance δ entre les deux points
intérieurs, pour être aussi près que possible de diviser le dernier intervalle par deux. Et le résultat
de ce test doit laisser une longueur ε. (Où δ est celui évoqué au titre de la dérivation numérique, et ε
la précision souhaitée.)
On en déduit
lN = ε = lN −1 /2 + δ/2 ,
2.2. MÉTHODES DIRECTES 33

soit lN −1 = 2ε − δ. À partir de là, on peut remonter la suite en utilisant la récurrence (2.1) :

lN = ε = f2 ε
lN −1 = 2ε − δ = f3 ε − f1 δ
lN −2 = 3ε − δ = f4 ε − f2 δ
lN −3 = 5ε − 2δ = f5 ε − f3 δ
..
.
l0 = fN +2 ε − fN δ

On doit donc
Algorithme Fibonacci (complet)
1. Déterminer le plus petit entier N tel que fN +2 ε − fN δ > l0 = b0 − a0 ,
2. calculer ε0 = (l0 + fN δ)/fN +2 ,
3. choisir d0 = a0 + fN +1 ε0 − fN −1 δ et c0 = b0 − fN +1 ε0 + fN −1 δ,
4. dérouler l’algorithme de Fibonacci ci-dessus, de n = 0 à N − 1. (Soit N pas)
Dans cet algorithme, on a évalué u en N + 3 points, pour réduire le segment contenant t∗ dans un
rapport de l’ordre de 1/(ρ+ )N . Il est prudent de suivre la suite des nombres de Fibonacci pour placer
les points intermédiaires à chaque pas plutôt que de se fier au “symétrique”, ce qui évite de laisser
s’accumuler de petites erreurs. En effet, cette procédure est très sensible à de petites erreurs sur la
position des points intérieurs, comme la suite de l’analyse va nous le montrer.

2.2.3 Section dorée

La suite de Fibonacci “à l’envers” engendrée par (2.1) est de la forme
n n
ln = α+ r+ + α− r−

où α+ et α− dépendent des deux termes initiaux. Comme r+ est de module inférieur à 1, le terme
en r+n tend rapidement vers 0. Par contre, r est de module supérieur à 1, et même plus précisément
−
r− < −1, de sorte que le terme en r− n diverge rapidement avec des signes alternés. C’est ce qui

explique la grande sensibilité de l’algorithme “de Fibonacci” au choix de c0 et d0 .

La seule façon de pouvoir poursuivre l’algorithme un nombre arbitraire de pas est de s’assurer que
2 l , . . ., l = r n l .
α− = 0, c’est à dire de choisir l1 = r+ l0 . Ainsi, l2 = (1 − r+ )l0 = r+ 0 n + 0
On a ainsi un algorithme bien plus facile à mettre en œuvre, qui consiste à appliquer l’algorithme
de Fibonacci ci-dessus, en plaçant à chaque pas dn à une distance r+ ln de an , et cn à une distance
(1 − r+ )ln = r+ 2 l , où on rappelle que
n
√
5−1
r+ = ' 0, 618 ,
2√
2 3− 5
1 − r+ = r+ = ' 0, 382 .
2
Cet algorithme est à très peu de chose près aussi efficace que le précédent et beaucoup plus simple.
On n’a pas besoin de déterminer à l’avance le nombre de pas qu’on va effectuer, il suffit de mettre
un test d’arrêt en comparant la longueur du segment [an , bn ] restant après le pas n à la précision ε
désirée. Aussi le récapitulons-nous à fin de référence.
Algorithme Section dorée
34 CHAPITRE 2. RECHERCHE UNIDIMENSIONNELLE

1. Faire [a0 , b0 ] = [a, b] , l0 = b − a , 2l ,

c0 = a + r+ d0 = a + r+ l0 .
0

2. Évaluer u en a0 , b0 , c0 , d0 .
3. Faire n := 0.
4. Calculer ln+1 = r+ (bn − an ).
5. Si u(cn ) < u(dn ),
faire an+1 := an , bn+1 := dn , dn+1 := cn , cn+1 := an+1 + r+ 2l
n+1 ,
si u(cn ) > u(dn ),
faire an+1 := cn , bn+1 := bn , cn+1 := dn , dn+1 := an+1 + r+ ln+1 .
√
6. Si ln+1 < ε ou ( 5 − 2)ln+1 < δ, donner an+1 < t∗ < bn+1 et arrêter,
si non, évaluer u en celui des points cn+1 ou dn+1 où elle n’est pas encore connue,
incrémenter n de 1 et retourner en 4.
exercice : Pourquoi le deuxième critère dans le test d’arrêt ? (Qui limite la précision qu’il est possible
d’obtenir à un peu plus de 4δ. On peut, si vraiment nécessaire, ajouter trois pas de dichotomie.)

2.3 Méthodes indirectes

Comme nous l’avons dit, nous regroupons sous ce titre douteux les méthodes qui reposent de
façon plus essentielle sur le calcul de dérivées. Rappelons que dans bien des cas, la recherche uni-
dimensionnelle est effectuée dans un algorithme de minimisation dans Rn pour trouver le meilleur
point dans une direction de recherche h choisie par ailleurs. Dans ces cas, notre fonction u(t) de ce
chapitre est en fait de la forme u(x + th), et ce qui joue le rôle de u0 (t) est la dérivée directionnelle
(∇u(x + th), h).

2.3.1 “Backtracking”
Sous ce vocable anglosaxon, nous visons une méthode simplissime dont l’objectif n’est pas vrai-
ment de trouver le minimum en t de u(t), mais un point “suffisamment bon”. Cette méthode sera
recommandée dans l’algorithme de Newton “protégé” de l’optimisation multivariable.
On est au voisinage d’un t donné, et on a calculé u0 (t). On cherche un nouveau t0 = t + θ. L’idée
est de partir avec une estimée trop lointaine t + θ0 , avec θ0 u0 (t) < 0, et de reculer jusqu’à ce que la
pente (u(t + θ) − u(t))/θ soit suffisamment grande en valeur absolue, comparée à u0 (t). On va donc
choisir deux nombres positifs r < 0, 5 et ρ < 1 (en pratique on prend ρ vers 0, 8) et de faire
Algorithme
1. choisir θ0 “siffisamment grand”
2. faire n := 0,
3. itérer jusquà ce que u(t + θn ) ≤ u(t) + ru0 (t)θn : θn+1 = ρθn .

2.3.2 Méthode de Newton

La méthode de dichotomie présentée comme méthode “directe” consiste en fait à résoudre l’équa-
tion u0 (x) = 0 par une méthode de dichotomie. On peut bien sûr résoudre cette même équation par
la méthode de Newton. On rappelle que cette méthode itérative consiste à prendre comme prochaine
estimée de la solution d’une équation le point qui serait la solution si la fonction à annuler coı̈ncidait
avec son approximation au premier ordre.
2.3. MÉTHODES INDIRECTES 35

On verra la méthode de Newton dans un cas un peu plus général au chapitre suivant, nous nous
contentons ici d’écrire l’algorithme auquel elle mène pour l’application présente :
Algorithme Newton unidimensionnel
1. Choisir t0 suffisamment proche de t∗ ,
2. faire n := 0,
3. itérer jusqu’à ce que |u0 (tn )| ≤ ε

tn+1 = tn − u00 (tn )−1 u0 (tn ) .

(Ajouter une clause limitant le nombre total d’itérations permis serait une prudence élémentaire.)
On sait que l’algorithme de Newton converge quadratiquement (on le reverra au chapitre suivant),
ce qui est très rapide. Sa grande faiblesse est qu’il est très sensible au choix de la condition initiale, et
peut facilement être mis en défaut si elle est trop éloignée de la solution recherchée. Aussi, on suggère
fréquemment de ne l’utiliser que pour affiner une solution approchée obtenue avec une méthode plus
rustique.
On voit aussi sur la formule que cet algorithme aura des problèmes si u00 (t) est trop proche de 0
au voisinage de t∗ . Il est toujours plus difficile de trouver un minimum très “plat” (avec une très petite
dérivée seconde) qu’un minimum bien marqué. Mais cet algorithme, utilisant explicitement la dérivée
seconde, peut y être plus sensible qu’un autre. Il peut être prudent de tester le module de cette dérivée
seconde (qui devrait rester positive en tout état de cause).

2.3.3 Approximation polynômiale

Approximation parabolique
Dans beaucoup d’applications, notamment celles liées au gradient “à pas optimal” ou “conjugué”,
on connait u et sa dérivée en a. Si dans le cas du gradient conjugué il est important de trouver le
minimum avec une certaine précision, il n’en va pas de même pour le gradient à pas optimal, et
pour quelques autres algorithmes où une approximation moyenne de l’optimum suffit à chaque pas.
(Relaxation,...) Dans ces cas, on peut utiliser la méthode suivante.
Supposons que nous connaissions u et sa dérivée en a. On choisit b tel que [a, b] ait une bonne
chance d’encadrer le minimum t∗ (mais ce n’est pas critique dans cette méthode), on évalue u(b), on
approxime u(t) par la parabole qui aurait même valeur en a et b et même dérivée en a. Et on prend
comme estimée de t∗ le point t̂ où cette parabole atteint son minimum.
Ceci mêne à l’estimée suivante :
(t − a)2
u(t) ' u(a) + (t − a)u0 (a) + α ,
2
soit
u0 (a)
t∗ ' t̂ = a − ,
α
où α est estimée par
(b − a)2
u(b) = u(a) + (b − a)u0 (a) + α
2
ce qui conduit finalement au choix
u0 (a)(b − a)2
t̂ = a − .
2[u(b) − u(a) − (b − a)u0 (a)]
36 CHAPITRE 2. RECHERCHE UNIDIMENSIONNELLE

On ne traite pas cette méthode comme un algorithme au sens propre du terme, car l’idée n’est pas
de l’appliquer itérativement, mais une seule fois dans un algorithme englobant qui requiert d’aller au
minimum dans une direction donnée à chacun de ses pas. On évite ainsi des itérations emboitées, au
prix d’une approximation parfois médiocre de ce minimum.
Au fur et à mesure que l’algorithme englobant se rapproche du minimum recherché, cette méthode
approxime de mieux en mieux le t∗ du pas en cours, comme l’indique le résultat ci-dessous.

Théorème 2.2 Si la fonction u est trois fois continument différentiable, α-convexe sur [a, b], et si ce
segment contient t∗ et a une longueur b − a qui tend vers zéro comme h, |t̂ − t∗ | tend vers zéro comme
h2 . (Donc l’erreur relative |t̂ − t∗ |/h tend vers zéro comme h.)

Démonstration Nous donnons une preuve directe de ce résultat, mais il peut aussi se déduire de la
preuve plus simple que nous donnons du théorème suivant.
Développons u au voisinage de t∗ , en notant u∗ = u(t∗ ) et u00 ∗ = u00 (t∗ ) > α :
1 ∗
u(t) = u∗ + u00 (t − t∗ )2 + ε(h3 )
2
où ε(z) désigne une quantité qui tend vers zéro comme z. On a aussi
∗
u0 (t) = u00 (t − t∗ ) + ε(h2 ) .

Reportons ces expressions dans la formule pour t̂. On trouve facilement

∗
u0 (a)(b − a)2 = u00 (a − t∗ )(b − a)2 (1 + ε(h))

et
1 ∗
u(b) − u(a) − u0 (a)(b − a) = u00 (b − a)2 (1 + ε(h)) ,
2
soit
u00 ∗ (a − t∗ )(b − a)2 (1 + ε(h))
t̂ = a − 1 00 ∗ 2 (1 + ε(h))
= a − (a − t∗ )(1 + ε(h)) = t∗ + ε(h2 ) .
2 u (b − a)
Donc, si cette procédure est utilisée, par exemple, dans un algorithme de gradient, au début, t̂ est
une approximation grossière du t∗ du pas en cours, mais on sait que cela suffit, et au fur et à mesure
que l’algorithme progresse, l’erreur relative tend vers zéro, permettant une bonne convergence de
l’algorithme global.

Approximations cubiques
En fait, la méthode la plus fréquemment utilisée est celle de l’approximation cubique, où la fonc-
tion u est approximée par une cubique, donc. En effet, sa dérivée est alors un polynôme de degré deux,
et on a encore une formule explicite pour son minimum t̂. Ainsi, au prix de calculs guère plus lourds,
on a une estimée meilleure d’un ordre (une erreur relative d’ordre deux), ce qui est excellent.
On prendra donc pour approximation de u

u(t) ' αt3 + βt2 + γt + δ ,

et donc
u0 (t) ' 3αt2 + 2βt + γ ,
2.3. MÉTHODES INDIRECTES 37

ce qui mène à l’approximation t∗ ' t̂ avec

p
β 2 − 3αγ − β
t̂ = .
3α
(On a supposé que γ = u0 (0) < 0 et β = u00 (0) > 0.)
Il reste à calculer α, β et γ. Plusieurs méthodes sont possibles (d’où le pluriel dans le titre de ce
sous-paragraphe).
Soit, en calculant u(b), on a facilement aussi la dérivée u0 (b), et cela fait assez d’information, soit
on considère que calculer une dérivée est plus difficile que la fonction elle même, et on peut alors
préférer calculer u(c) en un point intermédiaire c ∈ [a, b].
Dans le premier cas, on peut évaluer les constantes par les formules suivantes :

2(u(b) − u(a)) − (b − a)(u0 (b) + u0 (a))

α= ,
(b − a)3

u0 (b) − u0 (a)
2β = − 3α(b + a) ,
b−a
γ = u0 (a) − 3αa2 − 2βa ,
voire, pour préserver la symétrie (ce qui peut être numériquement utile) la formule symétrisée pour γ
en faisant 1/2 de cette expression plus la même en b.
Dans le deuxième cas, nous introduisons les quantités

u(a) − u(b) − u0 (a)(a − b)

∆(a, b) =
(a − b)2

et de même pour ∆(a, c), et on peut montrer les formules suivantes :

∆(a, c) − ∆(a, b)
α= ,
b−c

c∆(a, b) − b∆(a, c)
β= − 2aα ,
b−c
γ = u0 (a) − 3αa2 − 2βa .
On a le résultat logique :

Théorème 2.3 Si la fonction u est quatre fois continument différentiable, α-convexe sur [a, b], et si
ce segment contient t∗ et a une longueur b − a qui tend vers zéro comme h, |t̂ − t∗ | tend vers zéro
comme h3 . (Donc l’erreur relative |t̂ − t∗ |/h tend vers zéro comme h2 .)

Démonstration Appelons û(t) notre approximation polynômiale de u, et ε(t) = u(t) − û(t) l’erreur
d’approximation. Le fait essentiel est le suivant :
Proposition 2.4

∀t ∈ [a, b] , ε0 (t) → 0 comme h3 .

38 CHAPITRE 2. RECHERCHE UNIDIMENSIONNELLE

En effet, le développement de Taylor à l’ordre 3 en a avec reste exact nous apprend que u peut sécrire

(t − a)4 (4) 0
u(t) = ū(t) + u (t ) , t0 ∈ [a, t]
24
où ū est un polynôme de degré 3. Si, dans les formules servant à calculer û on remplace u(b) et u(c)
par ū(b) et ū(c), on obtient exactement ū à la place de û. Mais on peut vérifier que û(t) sécrit aussi
(exercice)
(t − a)2 (t − c) (t − a)2 (t − b)
û(t) = ua (t) + u(b) + u(c)
(b − a)2 (b − c) (c − a)2 (c − b)
où ua (t) est un polynôme de degré 3 qui ne dépend que de u(a) et u0 (a), mais pas de u(b) et u(c). (ua
et sa dérivée coincident avec u(a) et u0 (a) en a, et ua s’annulle en b et en c, ce qui au vu des formules
ci-dessus le définit complètement.) Les quantités u(b) et u(c) interviennent (linéairement) dans la
formule ci-dessus avec des coefficients uniforméments bornés quand h → 0. Donc les coefficients du
polynôme û approchent ceux de ū comme ū approche u en b et en c, c’est à dire comme h4 . Donc leur
différence et sa dérivée approchent zéro comme h4 . Ainsi û0 approche ū0 en h4 uniformément en t.
Mais le développement de Taylor de u0 en a nous apprend que ū0 approche u0 comme h3 uniformément
en t. La proposition est démontrée.
Par définition de t̂ on a û0 (t̂) = 0, et donc u0 (t̂) = ε0 (t̂). En outre, u0 (t∗ ) = 0, et donc, comme u
est supposée α-convexe, par (1.14) u0 (t̂) ≥ α|t̂ − t∗ |. Soit bien
1 0
|t̂ − t∗ | ≤ ε (t̂)
α
d’où, avec la proposition, le résultat annoncé.
On remarque que cette méthode de preuve, élémentaire si non élégante, s’étend à un ordre quel-
conque.
Il reste une question ouverte : celui du meilleur choix de c dans la dernière méthode. Y a-t-il
un choix de c qui annule le terme d’ordre 3 dans l’erreur t̂ − t∗ ? C’est peu probable, parce que ce
terme est un polynôme de degré trois en t∗ , et fait donc intervenir quatre coefficients. Mais à notre
connaissance, cette question n’a jamais été regardée. Il faut dire que les calculs demandent à être faits
à la machine, car ils sont gros ! (Il faut développer u au moins à l’ordre cinq.)
Chapitre 3

Optimisation dans Rn

3.1 Bonnes fonctions

Nous présentons ci-dessous divers algorithmes de recherche de minimum dont certains, tel le
gradient à pas optimal, sont très robustes, i.e. convergent dans bien des situations délicates. Cependant,
tant parce qu’on ne peut pas toujours faire beaucoup mieux que par souci de simplicité mathématique,
nous ne démontrerons jamais la convergence que pour une classe de fonctions u assez restreinte, que
nous appellerons les bonnes fonctions (appellation totalement indigène).

Définition 3.1 (Bonnes fonctions) Nous appellerons bonnes fonctions les fonctions u(·) de Rn dans
R α convexes et de dérivée première β-Lipshitz-continue :

∀x, y ∈ Rn , ku0 (y) − u0 (x)k ≤ βky − xk . (3.1)

Si nous définissons la convexité via la positivité de la dérivée seconde, —mais la définition ci-
dessus est plus générale— c’est dire que u(·) doit être deux fois continument dérivable, et qu’il doit
exister deux réels positifs α et β tels que

∀x ∈ Rn , αI ≤ D2 u(x) ≤ βI ,

ou encore que les valeurs propres de D2 u sont comprises pour tout x entre α et β.
Rappelons que u(·) étant α-convexe, elle satisfait outre (3.1) les inégalités (1.14), (1.15), (1.16) et
(1.17), et son minimum est atteint sur Rn .

3.2 Optimisation non contrainte

3.2.1 Relaxation
Nous commençons par un algorithme “direct”, c’est à dire qui ne nécessite pas le calcul des
dérivées de u. Il n’est à recommander que si ce calcul est vraiment difficile, et (si possible) seulement
en petite dimension.
L’algorithme est excessivement simple (voire simpliste), et consiste à faire des minimisations
unidimensionnelles en chacune des variables xi successivement. Un test d’arrêt raisonnable porte sur
la décroissance de u aprés qu’on ait fait cette minimisation sur chacune des coordonnées.
Soit donc xk une estimée de la solution. Nous introduisons les “pas fractionnaires” xk+i/n , que
nous noterons xk,i pour simplifier, de la façon suivante : xk,1 ne diffère de xk que par sa première

39
40 CHAPITRE 3. OPTIMISATION DANS RN

coordonnée obtenue en minimisant u(·) par rapport à cette première coordonnée toutes les autres
étant figées à leur valeur dans xk . On passe ensuite à xk,2 en figeant toutes les coordonnées à leur
valeur dans xk,1 , sauf la deuxième par rapport à laquelle on minimise u, et ainsi de suite. Et on posera
xk+1 = xk,n . Nous décrivons formellement cet algorithme. Nous notons ei le vecteur de base numéro
i de Rn .
Algorithme Relaxation
1. Choisir une estimée initiale x0 et faire k := 0
2. faire xk,0 := xk .
3. pour i = 1 . . . n, calculer xk,i par

u(xk,i ) = min u(xk,i−1 + θei ) .

θ

4. faire xk+1 := xk,n

5. si u(xk ) − u(xk+1 ) < ε, stop, sinon incrémenter k de 1 et retourner en 2.

Théorème 3.1 (Convergence de l’algorithme de relaxation) Si u(·) est une

bonne fonction (au sens de la définition ci-dessus), l’algorithme de relaxation converge vers l’argu-
ment x∗ du minimum.

Démonstration On remarquera d’abord que par l’α-convexité, et à cause de (1.17) par exemple, la
suite {xk } est bornée.
Par construction, on a u(xk,i ) < u(xk,i−1 ), et donc aussi u(xk+1 ) < u(xk ). Comme par α-
convexité, u est bornée inférieurement, ku(xk ) − u(xk+1 )k → 0 quand k → ∞. Plus précisément, la
définition de xk,i implique que u0 (xk,i )ei = 0, de sorte que l’inégalité (1.15) donne
α k,i−1
u(xk,i−1 ) − u(xk,i ) ≥ kx − xk,i k2 .
2
En sommant ces inégalités de i = 1 à n, il vient
α k
u(xk ) − u(xk+1 ) ≥ kx − xk+1 k2 .
2
Donc, en particulier, kxk − xk+1 k tend vers zéro, et donc aussi chacune de ses composantes kxk,i−1 −
xk,i k.
La deuxième inégalité d’α-convexité (1.14) donne, en utilisant ∇u(x∗ ) = 0 :
n
X
αkxk − x∗ k2 ≤ (∇u(xk ), xk − x∗ ) = u0i (xk )(xki − x∗i ) ,
i=1

où u0i = u0 (x)ei désigne bien sur la dérivée partielle en xi . Mais à nouveau, par la définition de xk,i ,
u0i (xk,i ) = 0. De sorte que l’inégalité (3.1) donne

ku0i (xk )k ≤ βkxk − xk,i k .

Nous utilisons cette évaluation dans l’inégalité précédente, pour obtenir

n
X
∗ 2
k
αkx − x k ≤ β kxk − xk,i k kxki − x∗i k .
i=1
3.2. OPTIMISATION NON CONTRAINTE 41

Nous savons que les kxk −xk,i k tendent vers zéro et que les kxki −x∗i k sont bornés. Donc kxk −x∗ k →
0, ce qu’il fallait démontrer.
Cette preuve ne dit pas que la méthode converge “bien”. Elle converge même souvent très mal, et il
n’est guère réaliste d’en espérer une bonne approximation de x∗ . Tout au plus permet-elle d’améliorer
parfois significativement la performance u(x) d’une estimée initiale à moindre frais. Les méthodes
qui suivent sont presque toujours préférables.

3.2.2 Gradient à pas optimal

L’algorithme

L’algorithme de gradient à pas optimal consiste à se déplacer “dans la direction de plus grande
pente”, c’est à dire dans la direction opposée au gradient, jusqu’au point “le plus bas” dans cette
direction. On a ainsi la description formelle suivante :
Algorithme Gradient à pas optimal
1. Choisir une estimée initiale x0 , faire k := 0.
2. Calculer ∇u(xk ). Si k∇u(xk )k < ε stop. Si non,
3. Calculer xk+1 := xk − θk ∇u(xk ) où le pas θk > 0 est déterminé par

u(xk+1 ) = min u(xk − θ∇u(xk ))

θ

4. incrémenter k de 1 et retourner en 2

Théorème 3.2 (Convergence de l’algorithme du gradient à pas optimal)

Soit u(·) une bonne fonction. L’algorithme du gradient à pas optimal converge vers l’argument x∗ du
minimum de u.

Demonstration Nous décomposons cette preuve pour souligner ce que nous apporte chaque hy-
pothèse.

1. La fonction u décroit à chaque pas, mais comme elle est α-convexe, elle est bornée inférieure-
ment. Donc ku(xk ) − u(xk+1 )k → 0.
2. Comme u00 est bornée par βI, et plus précisément grâce à l’inégalité (3.1), la décroissance au pas
k est d’au moins k∇u(xk )k2 /2β, comme l’indique le lemme que nous démontrons séparément
ci-dessous, utilisé avec ĥ = −g/kgk, et donc γ = 1. En conséquence, en tenant compte du (1)
ci-dessus, ∇u(xk ) → 0.
3. Par l’α-convexité, et plus précisément par l’inégalité (1.16), on en déduit que xk → x∗ , ce qu’il
fallait démontrer.

Remarque 3.1 On remarque qu’on n’a pas vraiment besoin de l’α-convexité de u. Il suffit (exercice)
de supposer que u est strictement convexe et tend vers l’infini à l’infini.

Il reste à démontrer le lemme. Nous le démontrons dans un cadre un peu plus large, qui nous
servira par la suite.
42 CHAPITRE 3. OPTIMISATION DANS RN

Lemme 3.3 Soit g = ∇u(x), et ĥ un vecteur de Rn de norme unité, satisfaisant l’inégalité

(g, ĥ) ≤ −γkgk (3.2)
avec 0 < γ ≤ 1.
Soit t+ ∈ R+ déterminé par
u(x + t+ ĥ) = min u(x + tĥ) ,
t∈R+

et soit x+ = x + t+ ĥ.
Sous l’hypothèse (3.1), on a
γ2
u(x) − u(x+ ) ≥ kgk2 (3.3)
2β

Démonstration Introduisons la fonction U (t) = u(x + tĥ). Notons que

U 0 (t) = (∇u(x + tĥ), ĥ)
et donc que U 0 (0) ≤ −γkgk.
Par (3.1), nous avons
k∇u(x + tĥ) − gk ≤ βt ,
soit, avec l’inégalité de Cauchy-Schwarz,
(∇u(x + tĥ) − g, ĥ) ≤ βt
soit encore
U 0 (t) = (∇u(x + tĥ), ĥ) ≤ (g, ĥ) + βt ≤ −γkgk + βt .
On utilise cette minoration pour évaluer U (τ ) :
Z τ
τ2
U (τ ) ≤ U (0) + (βt − γkgk) dt = u(x) + β − γkgkτ .
0 2
Enfin, on utilise cette minoration en τ = βγ kgk pour obtenir

γ γ2
u(x+ ) = min U (τ ) ≤ U ( kgk) ≤ u(x) − kgk2 ,
τ β 2β
ce qu’il fallait démontrer.
Cet algorithme est très robuste, et beaucoup plus efficace que l’algorithme de relaxation. On
démontre sa convergence pour les “bonnes fonctions”, mais c’est un “bon algorithme” au sens où
“une bonne théorie est une théorie qui continue à donner de bons résultats quand on l’utilise en de-
hors des hypothèses sous lesquelles elle a été établie” 1 . Naturellement, en l’absence de convexité, il
converge vers le minimum local dans le “bassin d’attraction” duquel on est parti. Il convient donc
éventuellement de refaire fonctionner l’algorithme avec divers conditions initiales.
Par contre, on ne doit pas attendre une bonne convergence près de l’optimum. C’est à dire que
l’algorithme est efficace pour faire décroitre rapidement la fonction, mais pas pour avoir une bonne
approximation de x∗ . Au point que Claude Lemaréchal a pu écrire dans un autre cours “on devrait
interdire cette méthode”. C’est un peu... totalitaire !. Mais pour approximer finement x∗ , il vaut mieux
finir avec une méthode du second ordre comme celle de l’algorithme que nous présenterons ensuite.
1. Holt Ashley
3.2. OPTIMISATION NON CONTRAINTE 43

Préconditionnement

On remarque que la preuve de convergence demeure si la direction de descente choisie “fait un

angle aigu” avec −∇u, comme le lemme le montre. Cette remarque a de nombreuses applications,
qui tournent autour de ce qu’on appelle le “préconditionnement”.
Supposons qu’on fasse un changement de variable x = Aξ, avec une matrice de changement de
variable A non singulière. On est conduit à considérer la fonction v(ξ) = u(Aξ). On peut montrer
(exercice) que c’est une “bonne fonction”, et qu’on peut donc lui appliquer l’algorithme du gradient.
Cela donne-t-il la même suite de points que l’algorithme initial ? La réponse est non comme nous
allons le voir.
Nous avons v 0 (ξ) = u0 (Aξ)A, soit en transposant ∇v(ξ) = At ∇u(Aξ) et donc, v(ξ − θ∇v) =
u(Aξ − θAAt ∇u). On a donc comme direction de descente −AAt g. Si A est régulière, AAt est posi-
tive définie, et il existe donc γ > 0 tel que (g, AAt g) ≥ γkgk2 . (γ est le carré de la plus petite valeur
singulière de A.) Donc, de la façon dont a été faite la preuve ci-dessus, on déduit immédiatement que
l’algorithme converge encore. Mieux ?, moins bien ? Cela dépend évidemment du choix de A, ou de
la matrice symétrique AAt , qui est appelée matrice de préconditionnement.
Les algorithmes de gradient conjugué, par exemple, peuvent être vus come des algorithmes de
gradient à préconditionnement soigné, et la méthode de Newton protégée ci-dessous comme un pré-
conditionnement “optimal”.
Plus simplement, on recommande en général d’utiliser au moins un préconditionnement diagonal
de la forme ξi = xi /x̄i où les x̄i jouent le rôle de facteur d’échelle, rendant en quelque sorte les ξi
sans dimension, et doivent être choisis de façon que la sensibilité de v(ξ) aux différents ξi soit à peu
près la même. Ce qui est obtenu en prenant des x̄i inversement proportionnels aux (ordre de grandeur
des) ∇i u.

3.2.3 Méthode de Newton

La méthode “pure”

Nous présentons maintenant une méthode “du second ordre”, en ce qu’elle utilise les dérivées
secondes de u, mais aussi en ce qu’elle converge (plus vite que) quadratiquement. C’est dans cette
famille de méthodes qu’il faut chercher si on veut approximer finement l’argument du minimum x∗ .
Le principe est simple, il consiste à résoudre l’équation u0 (x) = 0 par la méthode de Newton.
Nous la rappelons d’abord pour la solution d’une équation g(x) = 0 où g(·) : Rn → Rn .
La méthode de Newton consiste à approximer g au premier ordre autour du dernier point connu,
et à prendre comme prochaine estimée de la solution x∗ de g(x) = 0, le point où cette approximation
linéaire s’annule. Ceci mène à

g(x) ' g(xk ) + g 0 (xk )(x − xk ) ,

soit
xk+1 = xk − [g 0 (xk )]−1 g(xk ) (3.4)

Cette méthode converge quadratiquement comme le montre le résultat suivant :

Théorème 3.4 (Convergence de la méthode de Newton) Si la fonction g est deux fois continument
différentiable au voisinage de x∗ , avec une dérivée première inversible en x∗ , il existe un voisinage de
x∗ dans lequel la méthode de Newton converge quadratiquement.
44 CHAPITRE 3. OPTIMISATION DANS RN

Démonstration Par continuité de g 0 , il existe un voisinage de V de x∗ dans lequel la matrice g 0 (x)

est inversible, et plus précisément a une plus petite valeur singulière bornée inférieurement par un
nombre positif α, de sorte que [g 0 (x)]−1 existe et a une norme bornée supérieurement par 1/α. De
même, dans ce voisinage, g 00 (x) existe et a une norme bornée par un nombre positif γ.
Écrivons l’itération de Newton comme

xk+1 − x∗ = [g 0 (xk )]−1 [g 0 (xk )(xk − x∗ ) − g(xk )]

Par le développement limité (1.10), et en se souvenant que par définition g(x∗ ) = 0, le dernier crochet
ci-dessus peut se ré-écrire comme ci-dessous, coordonnée par coordonnée :
1 k 
gi0 (xk )(xk − x∗ ) − gi (xk ) = (x − x∗ ), D2 gi (xk + θ(x∗ − xk ))(xk − x∗ ) ,
2
de sorte que ce crochet est borné en norme par
γ k
kg 0 (xk )(xk − x∗ ) − g(xk )k ≤ kx − x∗ k2 .
2
Donc, en utilisant la borne sur k(g 0 )−1 k,
γ
kxk+1 − x∗ k ≤ kxk − x∗ k2 .
2α
Il reste à s’assurer que si on part suffisamment prés de x∗ , la suite engendrée ne sort pas de V, ce qui
se fait en bornant la somme des kxk − x∗ k pour kx0 − x∗ k suffisamment petit. Ainsi, le théorème est
démontré.
On voit que si g(x) = ∇u(x), nous avons donné un algorithme de recherche de minimum dans
un ouvert, à savoir h i −1
xk+1 = xk − D2 u(xk ) ∇u(xk ) , (3.5)
et montré sa convergence quadratique si u est trois fois continument différentiable. La grande faiblesse
de cet algorithme très rapide est sa très grande sensibilité aux conditions intiales. C’est pourquoi on
conseille souvent de commencer la recherche du minimum avec une méthode plus robuste comme
celle du gradient, et de finir avec la méthode de Newton.

Remarque 3.2 Une remarque importante est que dans cette version “pure”, on n’a pas besoin d’in-
verser D2 u(xk ) en dépit de ce que semble indiquer la formule (3.5). En effet, il suffit de résoudre le
système linéaire
D2 u(xk )(xk+1 − xk ) = −∇u(xk ) , (3.6)
ce qui peut être significativement moins long. Cela n’évite pas le calcul de la matrice des dérivées
secondes D2 u(xk ) à chaque pas.

La méthode de Newton “protégée”

Une méthode recommandée consiste à considérer h = −[g 0 (xk )]−1 g(xk ) comme une direction
de descente, et, comme précédemment effectuer une recherche unidimensionnelle de minimum dans
cette direction, sachant que l’optimum devrait se trouver près de 1. Le caractère positif défini de
g 0 (x) = D2 u(x), si u est α-convexe, garantit la convergence de cet algorithme en vertu de la preuve
de convergence de l’algorithme de gradient. En pratique, la direction de recherche est si bonne qu’il
suffit de faire du “backtracking” (cf. le paragraphe 2.3.1) depuis le “pas de Newton” 1.
3.3. OPTIMISATION SOUS CONTRAINTES INÉGALITÉ 45

Il n’est pas très sûr que le poids du calcul de [D2 u(xk )]−1 en vaille la peine si on est vraiment trop
loin de l’optimum. Divers aménagement de la méthode de Newton visent à alléger cette étape. Avant
de les évoquer, notons que, toujours si la fonction u est convexe, la matrice ∇g = D2 u à inverser est
positive définie, de sorte qu’on dispose pour cette inversion de la très efficace méthode de Cholesky.

Les méthodes de Newton modifiées, quasi-Newton

Le coût principal en calcul dans la méthode de Newton est dans l’inversion de D2 u(xk ) à chaque
pas. Remarquons que la preuve de convergence demeure inchangée si on remplace D2 u(xk ) par
D2 u(x∗ ), qui lui, est constant, et ne demanderait donc qu’une inversion unique. Bien sûr, on ne connait
pas x∗ , et donc pas non plus D2 u(x∗ ) (sauf si D2 u(x) est constante, mais alors u est une simple
fonction quadratique, et l’algorithme de Newton converge en un pas !) Par contre, la preuve demeure
encore si on remplace D2 u(xk ) par une matrice Hk telle que

kHk − D2 u(x∗ )k ≤ ηkxk − x∗ k .

On aura en effet,
xk+1 − x∗ = Hk−1 [Hk (xk − x∗ ) − g(xk )] ,
et le crochet s’écrit maintenant

Hk (xk − x∗ ) − g(xk ) = D2 u(x∗ )(xk − x∗ ) − g(xk ) + (Hk − D2 u(x∗ ))(xk − x∗ )

dont le module est encore majoré par un terme quadratique en kxk − x∗ k.

On va donc chercher à avoir, à moindre prix, une suite Hk tendant vers D2 u(x∗ ) avec xk . Une idée
simplissime, mais assez efficace si elle est utilisée avec doigté, consiste à ne remettre à jour D2 u(x),
et donc D2 u(x)−1 , que de temps en temps. Typiquement tout les deux à cinq pas.
Des méthodes plus sophistiquées existent, sous le nom de méthodes de “quasi Newton”, qui s’ap-
parentent au gradient conjugué, et cherchent une formule itérative pour approximer Hk−1 . Nous n’en
parlerons pas plus ici.

3.3 Optimisation sous contraintes inégalité

3.3.1 Position du problème
La plupart des problèmes d’optimisation, notamment en théorie de la décision, (contrôle, économie
théorique, recherche opérationnelle, etc.) se présentent avec des contraintes sur les variables de décision
x. D’une manière abstraite, ces contraintes se traduisent par l’existence d’un ensemble C ⊂ Rn de
variables admissibles. On cherche donc x∗ ∈ C tel que

∀x ∈ C u(x) ≥ u(x∗ ) .

L’algorithme du gradient projeté ci-dessous considère le problème sous cette forme, et utilisera la
projection sur C en le supposant convexe fermé. Ceci suppose que cette projection soit facile à faire,
ce qui est rarement le cas.
En pratique, le cas le plus courant est celui où l’ensemble des variables admissibles est défini
indirectement par des contraintes de la forme

C = {x ∈ Rn | fi (x) ≤ 0 , i = 1, 2, . . . , p} ,
46 CHAPITRE 3. OPTIMISATION DANS RN

qu’on écrira aussi f (x) ≤ 0. Les méthodes intéressantes sont alors celles qui sont explicites en fonc-
tion de f .
La théorie de la dualité (multiplicateurs de Lagrange, de Kuhn et Tucker, etc) est à la base de
plusieurs algorithmes visant à répondre à ce type de problème. Nous ne présenterons, dans cette
famille, que l’algorithme d’Uzawa, réputé le plus robuste.
Enfin, par souci de présenter les méthodes les plus employées, nous donnerons des méthodes de
pénalisation, qui échangent une plus grande simplicité théorique contre une efficacité douteuse.

3.3.2 Gradient projeté

Projection sur un convexe simple
On a rappelé dans le chapitre premier l’existence de la projection PC (x) d’un vecteur x sur un
convexe fermé C. Cette projection est une opération simple dans un petit nombre de cas. Typiquement
trois cas :
– C est un pavé, de la forme ai ≤ xi ≤ bi , où les ai et bi sont donnés. Alors la projection se fait
coordonnée par coordonnée et consiste simplement à “ramener xi dans [ai , bi ]”. On peut écrire
cela comme
(PC (x))i = max{ai , min{xi , bi }} .
– C est une boule, de la forme kx − ξk ≤ ρ où le vecteur ξ de Rn et le réel positif ρ sont donnés.
Alors la projection consiste à ramener x dans la boule le long du rayon, soit
ρ
PC (x) = ξ + min{1, }(x − ξ) .
kx − ξk

– C est un demi-espace, de la forme (p, x) ≤ a, où le vecteur p de Rn et le réel a sont donnés. La

projection consiste à ramener x dans le demi espace parallèlement à p :
 
(p, x) − a
PC (x) = x − max 0, p.
(p, p)

Exercice 3.1 Vérifier les formules ci-dessus.

L’algorithme
L’algorithme du gradient projeté décrit ci-dessous n’existe que dans cette version à pas fixe. A
notre connaissance, il n’y a pas de version “à pas optimal”.
Naturellement, si le convexe des contraintes sur lequel on projette est tout Rn , la projection est
l’identité, et donc la preuve ci-dessous montre la convergence de l’algorithme de gradient à pas fixe
sans projection, pour le problème sans contraintes. (Pourvu, toujours, qu’on ait choisi le pas convena-
blement.) En fait, dans le cas sans contrainte, l’algorithme à pas fixe n’est vraiment pas à recomman-
der.
L’algorithme de gradient projeté est fondé sur la remarque suivante. L’inégalité d’Euler pour l’op-
timisation dans un convexe (1.18) nous dit que si x∗ fournit le minimum de u sur C, alors

∀t > 0 PC (x∗ − t∇u(x∗ )) = x∗ . (3.7)

En effet, soit x∗ est intérieur à C, et alors ∇u(x∗ ) = 0, ce qu’implique (3.7) car pour un point
intérieur, et si ∇u(x∗ ) 6= 0, il existe t suffisamment petit pour que x∗ − t∇u(x∗ ) ∈ C, de sorte qu’il
3.3. OPTIMISATION SOUS CONTRAINTES INÉGALITÉ 47

serait sa propre projection. Soit x∗ est un point frontière, et (3.7) est équivalent à l’affirmation que
−∇u(x∗ ) est une normale extérieure à C, ce que dit (1.18).
L’algorithme s’écrit comme une méthode de recherche du point fixe dans (3.7) :
Algorithme Gradient projeté
1. Choisir une éstimée initiale x0 , un pas t > 0, faire k := 0,
2. calculer xk+1 = PC (xk − θ∇u(xk )),
3. si kxk+1 − xk k ≤ ε stop, si non, incrémenter k de 1 et retourner en (2)
On a le résultat de convergence suivant :

Théorème 3.5 (Convergence de l’algorithme du gradient projeté) Si u(·)

est une bonne fonction, et si 0 < θ < 2/β, l’algorithme de gradient projeté converge vers le minimum
x∗ de u sur C.

Remarque 3.3 Le test d’arrêt ne peut pas porter sur le module du gradient de u, dont on ne sait pas a
priori s’il sera grand ou petit. Le test proposé ici vérifie donc qu’il soit “bien orthogonal” au bord de
C si xk est sur ce bord, et petit si-non. En fait, pour l’utilisation avec une fonction u dont la convexité
est douteuse, il vaut mieux faire porter ce test sur la différence kxk+1 − xk−4 k par exemple, c’est à
dire prendre en considération l’évolution de l’algorithme sur plusieurs pas.

Démonstration Puisque nous avons reconnu dans l’algorithme une itération de point fixe (aussi dite
de Picard), montrons que la fonction x 7→ PC (x − θ∇u(x)) est une contraction pour θ choisi comme
dans le théorème. On sait (théorème 1.38) que la projection sur C est une contraction au sens large.
Il suffit donc de montrer que ϕθ (x) := x − θ∇u(x) est une contraction stricte. On a ∇ϕθ (x) =
I − θD2 u(x). Cette matrice est symétrique. Sa norme est donc son rayon spectral, le module de sa
valeur propre de plus grand module. Or ses valeurs propres sont de la forme 1 − θλ(D2 u(x)), où les
λ(D2 u(x)) sont les valeurs propres de D2 u(x). Par hypothèses, celles-ci sont comprises entre α et
β. Il suffit donc de choisir θ de façon que |1 − θα| < 1 et |1 − θβ| < 1, soit θ > 0 et θ < 2/β
respectivement. (Cette dernière condition assurant aussi θ < 1/α.) Ce qui établit la convergence de
l’algorithme sous la condition annoncée.
On peut se demander quel est le θ “optimal”, ou au moins celui pour lequel le module de Lipshitz
de ϕθ est minimal. On voit assez facilement qu’il est tel que 1 − θα = θβ − 1, soit θ = 2/(α + β).
Et alors, le module de Lipshitz de ϕθ est 1 − 2α/(α + β).
On voit que ce nombre est d’autant plus proche de 1 que α est petit. C’est une constante des
problèmes d’optimisation que trouver le minimum est d’autant plus difficile que le module d’alpha
convexité est petit.
En fait, ce théorème souligne surtout la principale difficulté liée à l’utilisation de cet algorithme.
C’est celle du choix du pas θ. En effet, en général, α et β ne sont pas connus —bien heureux s’ils
existent—, et il faut donc des heuristiques d’adaptation du pas θ.
Remarquons d’abord qu’à en croire le calcul précédent, un pas trop petit peut ralentir considéra-
blement la convergence, pas l’empécher comme peut le faire un pas trop grand. Il faut quand même
trouver un moyen d’apprécier le pas qu’on peut se permettre. Une règle qu’on peut proposer est la
suivante. Comparer la décroissance de u obtenue, u(xk ) − u(xk+1 ), à son estimation au premier ordre
u0 (xk )(xk − xk+1 ) = θk∇u(xk )k2 (qui doit être plus grande). Si ces deux quantités coincident “très
bien”, disons à 10% près, on peut sans doute doubler le pas. Si elles coincident mal, disons plus mal
que dans un rapport deux, il faut sans doute diviser le pas par deux.
48 CHAPITRE 3. OPTIMISATION DANS RN

Cet algorithme n’est sans doute pas très performant par rapport à ceux que nous avons vus jus-
qu’ici. Il faut bien comprendre que la minimisation sous contrainte est un problème plus difficile que
la minimisation sans contrainte, et qu’on ne peut espérer trouver des algorithmes aussi efficaces.

3.3.3 Algorithme d’Uzawa

L’algorithme d’Uzawa exploite le point selle du théorème de Kuhn et Tucker. Il va donc chercher à
minimiser le lagrangien par rapport à x et à le maximiser par rapport à la variable duale, disons p. Plus
précisément, il consiste à remettre à jour l’estimée courante de la variable duale par un pas de gradient
—et le gradient du lagrangien par rapport à p est particulièrement simple—, puis à p fixé, à aller
jusqu’au minimum en x. Ce sera donc un “méta algorithme”, puisque nous ne dirons pas comment
effectuer la minimisation en x. Observons pourtant, —et c’est tout ce que la dualité fait pour nous—,
que dans cette minimisation, les contraintes f (x) ≤ 0 ont disparu.
Comme dans le théorème 1.41, nous permettons ici à certaines contraintes de rester “abstraites”
sous la forme x ∈ C tandis que d’autres sont rendues “concretes” sous la forme f (x) ≤ 0. La mini-
misation du lagrangien est alors à effectuer dans C. Si C est tout Rn , on a alors affaire à un problème
de minimisation sans contrainte. Donc un algorithme de gradient à pas optimal, par exemple, est pos-
sible. La dualité a ainsi ramené un problème de minimisation sous contraintes à une suite de problèmes
de minimisation sans contrainte.
On note P+ la projection sur le cône positif de Rn , qui consiste juste à ramener à zéro toute
composante négative du vecteur à projeter.
Algorithme Uzawa
1. Choisir une estimée initiale p0 ≥ 0 Faire k := 0.
2. Calculer xk par
u(xk ) + (pk , f (xk )) = min[u(x) + (pk , f (x))] .
x∈C

3. faire
pk+1 = P+ [pk + ρf (xk )] .
Si kpk+1 − pk k ≤ ε, stop. Sinon, retourner en 2
Cet algorithme est en fait un algorithme de gradient (en p) projeté (sur le cône positif). On
démontre sa convergence d’une façon analogue à la preuve de convergence de cet algorithme. (Cf
théorème 3.5)

Théorème 3.6 (Convergence de l’algorithme d’Uzawa) Si u(·) est une bonne fonction, et f (·) est
Lipshitz-continue de coefficient γ, pour ρ < 2α/γ 2 , la suite {xk } engendrée par l’algorithme d’Uzawa
converge vers l’optimum x∗ .

Démonstration Remarquons d’abord que la condition des écarts complémentaires entraine que pour
tout ρ > 0,
p∗ = P+ [p∗ + ρf (x∗ )] .
Ainsi, par le théorème 1.38, on a

kpk+1 − p∗ k ≤ kpk − p∗ + ρ(f (xk ) − f (x∗ ))k .

En élevant au carré, il vient

kpk+1 − p∗ k2 ≤ kpk − p∗ k2 + 2ρ(pk − p∗ , f (xk ) − f (x∗ )) + ρ2 kf (xk ) − f (x∗ )k2 (3.8)

3.3. OPTIMISATION SOUS CONTRAINTES INÉGALITÉ 49

Nous allons majorer le double produit (le terme central du deuxième membre). Remarquons que, grâce
au fait que les pi sont positifs, u + (p, f ) est encore convexe, et même α-convexe, en x. Par (1.18), la
minoration (1.17) reste valide pour la minimisation dans un convexe fermé. On a donc
α k
u(xk ) + (pk , f (xk )) ≤ u(x∗ ) + (pk , f (x∗ )) − kx − x∗ k2 ,
2
α
u(x∗ ) + (p∗ , f (x∗ )) ≤ u(xk ) + (p∗ , f (xk )) − kxk − x∗ k2 .
2
Sommons membre à membre et faisons passer ce qu’il faut à gauche, pour obtenir

(pk − p∗ , f (xk ) − f (x∗ )) ≤ −αkxk − x∗ k2 .

Ensuite, l’hypothèse que f est Lipshitz de coefficient γ donne

kf (xk ) − f (x∗ )k ≤ γkxk − x∗ k .

En reportant ces deux majorations dans (3.8), il vient

kpk+1 − p∗ k2 ≤ kpk − p∗ k2 + (ρ2 γ 2 − 2ρα)kxk − x∗ k2 .

Si on a choisi ρ < 2α/γ 2 , il existe δ > 0 tel que ρ2 γ 2 − 2ρα < −δ, d’où

kpk+1 − p∗ k2 ≤ kpk − p∗ k2 − δkxk − x∗ k2

qu’on peut ré-écrire

δkxk − x∗ k2 ≤ kpk − p∗ k2 − kpk+1 − p∗ k2 .
Donc, la suite kpk −p∗ k2 est décroissante. Comme elle est composée d’éléments positifs, elle converge,
donc la différence du deuxième membre ci-dessus tend vers zéro. Donc kxk − x∗ k tends vers zéro, ce
qu’il fallait démontrer.
Terminons en évoquant l’interprétation “économique” du théorème de Kuhn et Tucker, et donc de
l’algorithme d’Uzawa. Si les containtes fi (x) ≤ 0 représentent des ressources “rares” dont on ne peut
utiliser que la quantité disponible, et les p∗i associés en sont les prix unitaires à l’équilibre, on voit que
le lagrangien est un coût économique total prenant en compte le “prix” des denrées rares utilisées, et
que l’algorithme consiste tout simplement à diminuer le prix des ressources qu’on n’utilise pas jusqu’à
saturation, et à augmenter le prix de celles pour lesquelles la demande dépasse la disponibilité. Rien
que de très naturel.
Dualisation partielle Une remarque importante doit être faite à ce stade. On a énoncé le théorème de
Kuhn et Tucker pour un problème de la forme

∀x ∈ K u(x) ≥ u(x∗ ) .

avec
K = C ∩ {x ∈ Rn | f (x) ≤ 0} ,
et on a “dualisé” —c’est à dire introduit dans le lagrangien— les seules contraintes “concretes”
f (x) ≤ 0. Les autres sont restées “abstraites” et prises en compte par la condition x ∈ C dans
les inéquations du point selle (1.19) et l’algorithme d’Uzawa.
Si C n’est en effet pas tout Rn , i.e. une partie des contraintes est restée non dualisée, l’étape de
calcul de xk dans l’algorithme est un problème de minimisation sous ces contraintes. Cette minimisa-
tion peut être faite à l’aide d’un autre algorithme qu’Uzawa et la dualité, menant à une “hybridation
50 CHAPITRE 3. OPTIMISATION DANS RN

d’algorithmes”. Si C est un convexe simple, cette minimisation peut être effectuée par un algorithme
de gradient projeté par exemple.
Bien sûr, le choix des contraintes à exprimer d’une façon ou de l’autre (à dualiser ou à ne pas
dualiser) est laissé à l’utilisateur, et il n’est efficace d’utiliser une dualisation partielle qu’en ne laissant
non dualisées que des contraintes simples. Typiquement, si on a un nombre important de contraintes
de borne ai ≤ xi ≤ bi , qui impliqueraient deux fois plus de multiplicateurs. On peut alors préférer
renoncer à les dualiser et les prendre en compte directement dans la minimsation du lagrangien par
projection.

3.3.4 Pénalisation
Les méthodes présentées ici essayent de ramener le problème contraint à un problème non contraint
de façon plus “naı̈ve”, mais qui peut marcher. En particulier, on les utilise de façon non itérative,
contrairement à Uzawa. On ne résoudra donc qu’un, ou quelques, problème(s) de minimisation sans
contrainte. L’idée est la suivante : on va “faire payer” au critère le fait pour x de sortir de C. Et si ce
“prix” est très élevé, l’optimum se trouvera dans C.

Pénalisation extérieure quadratique

Supposons donc que l’ensemble des x admissibles est donné par

fi (x) ≤ 0, i = 1, . . . , p

que nous écrivons aussi f (x) ≤ 0, où f : Rn → Rp . Introduisons en outre la fonction f + (x) appelée
partie positive de f , définie par

fi (x) = max{0 , fi (x)} , i = 1, . . . , p .

Notons la proposition :

Proposition 3.7 Si f est dérivable, kf + k2 l’est aussi, et on a

dkf + (x)k2
= 2(f + (x))t f 0 (x)
dx

Démonstration Notons d’abord que la proposition n’est pas (entièrement) évidente, car f + , elle,
n’est en général pas dérivable aux points où f (x) = 0, donc notamment en x∗ . Prenons le cas d’une
fonction f scalaire. Il suffira ensuite d’appliquer le résultat à chaque composante de f .
Remarquons d’abord que si f (x) > 0 ou f (x) < 0, par continuité l’inégalité reste vraie dans
un voisinage de x. Dans le premier cas f + = f et dans le deuxième f + = 0 dans ce voisinage. La
formule [(f + )2 ]0 = 2(f + )f 0 est alors évidemment vérifiée.
Soit donc x un point où f (x) = 0 (c’est à dire un point frontière de C). Soit ei un vecteur de base
de Rn . On a donc dans le “quotient différentiel”
1 + 1
[(f (x + tei ))2 − (f + (x))2 ] = (f + (x + tei ))2
t t
et
1 + 1
0≤ (f (x + tei ))2 ≤ (f (x + tei ))2 ,
t t
3.3. OPTIMISATION SOUS CONTRAINTES INÉGALITÉ 51

la dernière inégalité parce que dans tous les cas, |f + (x)| ≤ |f (x)|. Dans les inégalités ci-dessus, f
étant dérivable, le terme de droite tends vers 2|f (x)f 0 (x)| = 0, et donc le terme central a une limite,
qui est zéro. On a donc démontré que (f + )2 a en x une dérivée partielle en xi , qui est nulle. Ceci
achève de démontrer la proposition.
Nous considérerons le critère augmenté
1
uε (x) = u(x) + kf + (x)k2 .
ε
La méthode de pénalisation consiste à utiliser la solution xε du problème sans contrainte :

uε (xε ) = minn uε (x)

x∈R

comme approximée de la solution x∗ du problème contraint. Cette méthode est justifiée par le résultat
suivant :

Théorème 3.8 Si u(·) est une bonne fonction, et f (·) est continue, xε → x∗ quand ε → 0.

Démonstration On a manifestement

uε (xε ) ≤ uε (x∗ ) = u(x∗ ) ,

la dernière égalité parce que par définition, f + (x∗ ) = 0. Si u est α-convexe, elle est bornée inférieu-
rement. Donc l’inégalité
1
u(xε ) + f + (xε )2 ≤ u(x∗ )
ε
+
implique que f (xε ) → 0 quand ε → 0. A fortiori, elle implique aussi que u(xε ) reste bornée, donc,
toujours avec l’α-convexité, que xε reste borné. Ainsi, pour toute suite décroissante de εk tendant vers
zéro, les xεk ont des points d’accumulation. Soit x̄ un tel point et ε0 une suite telle que les xε0 → x̄.
Par continuité de f + , f + (x̄) = 0, et x̄ est donc admissible.
On a aussi
u(xε ) ≤ uε (xε ) ≤ u(x∗ ) .
Donc par passage à la limite, (u est continue), u(x̄) ≤ u(x∗ ). Comme x̄ est admissible, il en découle
que u(x̄) = u(x∗ ) et x̄ est optimal. Par l’α-convexité de u, l’optimum est unique. Donc c’est toute la
suite qui converge vers x∗ .
On démontre aussi que, si la matrice f 0 (x∗ ) est injective, ce qui est une hypothèse naturelle, le
produit (1/ε)f + (xε ) tend vers un vecteur λ de Rp , de sorte qu’à l’optimum on a u0 (x∗ ) + λt f 0 (x∗ ) =
0. Ce λ est le multiplicateur de Lagrange (ou de Kuhn et Tucker) associé au problème d’optimisation
sous contraintes.
Cette méthode reste délicate d’emploi pour plusieurs raisons. D’abord, on ne va pas résoudre
beaucoup de fois le problème pénalisé : c’est un travail potentiellement important. Donc quel ε choisir
est non trivial. De plus, si u est trop “plat” au voisinage de C, son rôle dans uε va être masqué par le
terme de pénalisation, nuisant à la précision de l’estimation de x∗ . (Certes, l’hypothèse d’α-convexité
limite ce risque en bornant inférieurement la courbure de u. Mais ceci souligne la dépendance de la
méthode à cette hypothèse qui est d’habitude difficile à tester.)
Autres pénalisations On a pénalisé le critère avec la fonction de pénalisation kf + k2 . Ceci pour avoir
un critère uε dérivable. Le prix à payer est que, génériquement, si le minimum recherché est sur la
frontière de C, on l’approchera par des points extérieurs. Les xε sont non admissibles.
52 CHAPITRE 3. OPTIMISATION DANS RN

On aurait pu prendre f + (x) tout simplement. Alors le critère n’est plus dérivable, mais il est
possible qu’un ε fini permette déjà d’obtenir le minimum exact. Le caractère non différentiable du
critère en x∗ est néanmoins une grave difficulté. En fait, il vaut mieux se référer alors à la théorie de
la dualité.

Pénalisation intérieure

Une autre approche possible est d’utiliser comme fonction de pénalisation une fonction “barrière”,
ϕ(x), qui tend
P vers l’infini quand x tend vers la frontière de C par l’intérieur. On pourrait penser à
ϕ(x) = − Pi 1/ϕi (x) par exemple, mais nous verrons à la section suivante qu’un meilleur choix est
ϕ(x) = − i ln(−fi (x)). Puis on va considérer uε (x) = u(x) + εϕ(x), avec ε assez petit pour que
ceci ne change guère le comportement du critère tant que les fi sont tous loins d’être nuls. On parle
alors d’algorithmes de “points intérieurs”, parce qu’on aborde l’optimum recherché par des points
intérieurs à C.
Soit donc d’une manière un peu plus générale, un critère perturbé uε (x) = u(x) + εϕ(x) ou ϕ(·)
est définie pour les x intérieurs au domaine des contraintes C (c’est à dire les x tels que fi (x) < 0
◦
pour i = 1, . . . p), convexe (continue) positive dans C, et tend vers l’infini quand x → ∂C.
◦
Alors, pour tout ε positif, uε atteint son minimum dans C, en un unique point xε dès lors que u
est strictement convexe.
On a alors le résultat attendu :

Théorème 3.9 Sous les hypothèses ci-dessus, et si u est une bonne fonction, xε → x∗ quand ε → 0.

Demonstration Il faut faire attention qu’on ne peut pas manipuler uε (x∗ ) parce que x∗ peut être (est
généralement) sur la frontière de C et donc ϕ peut n’y être pas définie. Soit donc δ un nombre positif
(arbitrairement petit) et xδ un point itérieur à C tel que u(xδ ) ≤ u(x∗ ) + δ. La suite d’inégalités
ci-dessous est facile à établir :

u(xε ) ≤ uε (xε ) ≤ uε (xδ ) = u(xδ ) + εϕ(xδ ) ≤ u(x∗ ) + δ + εϕ(xδ ) .

En prenant ε ≤ δ/ϕ(xδ ) on en déduit u(xε ) ≤ u(x∗ ) + 2δ. Et comme δ était arbitraire, on en déduit
que u(xε ) tend vers u(x∗ ) quand ε tend vers zéro. Si u est α-convexe, on en déduit par utilisation de
(1.17) que xε tend vers x∗ .
Cette idée est à la base d’une méthode qui est aujourd’hui la meilleure connue si le problème est
réellement convexe et si les dérivées secondes sont faciles à calculer. Nous la présentons ci-dessous.

3.3.5 Méthode du chemin central

Avant d’exposer cette utilisation de la pénalisation intérieure, montrons un résultat qui constitue
la partie facile de la théorie de la dualité.
P Le problème considéré est toujours le même. Posons comme
précédemment L(x, p) = u(x) + i pi fi (x).

Lemme 3.10 Soit x∗ la solution du problème d’optimisation sous contraintes inégalité. Soit p ∈ Rm .
On a
∀p ≥ 0 , min L(x, p) ≤ u(x∗ ) . (3.9)
x
3.3. OPTIMISATION SOUS CONTRAINTES INÉGALITÉ 53

Démonstration
P On a la suite suivante d’inégalités faciles (la seconde vient de ce que pour x ∈ C,
i pi fi (x) ≤ 0) :
min L(x, p) ≤ min L(x, p) ≤ min u(x) = u(x∗ ) .
x x∈C x∈C

Nous utilisons la méthode de pénalisation intérieure avec les fonctions barrières − ln(−fi (x)).
Posons donc X
uε = u − ε ln(−fi (x)) .
i

C’est une fonction fortement convexe. Notons xε son unique minimum sur Rn . Nous affirmons le
théorème facile, mais surprenant :

Théorème 3.11 On a
u(xε ) − pε ≤ u(x∗ ) ≤ u(xε ) .
(Où p est ici le nombre de contraintes scalaires : i = 1, . . . , p.)

Demonstration Le point xε est caractérisé par

X 1
∇uε (xε ) = ∇u(xε ) − ε ∇fi (xε ) = 0 .
fi (xε )
i

Posons alors pi = −ε/fi (xε ). Ce sont des nombres positifs. L’égalité ci-dessus s’écrit
X
∇u(xε ) + pi ∇fi (xε ) = ∇x L(x, p) = 0 .
i

Maintenant, ce L(x, p) est une fonction convexe de x. Donc la condition ci-dessus implique qu’elle
atteint son minimum en xε . Utilisons alors le petit lemme précédent, en remarquant en outre que
pi fi (xε ) = ε. Le théorème en découle immédiatement.
Ainsi non seulement xε approche x∗ quand ε → 0, mais en outre on sait avec quelle précision
∗
u(x ) est atteint.
Voici quel est l’usage qu’on fait de ce résultat dans la méthode du chemin central. On commence
typiquement avec ε = 1 (bien qu’une autre valeur soit bien sûr possible, et peut-être préférable dans-
certains cas.) On résoud le problème min uε sans contrainte par la méthode de Newton. La première
application de cette méthode peut présenter ses difficultés habituelles et demander un peu de soin.
Ensuite, on fait décroı̂tre ε, typiquement d’un ordre de grandeur à chaque fois, et on résoud à nou-
veau le problème sans contrainte en prenant la solution à l’étape précédente comme initialisation de
la méthode de Newton. Cette fois on a une convergence très rapide de cette méthode. Et on continue
à faire décroı̂tre ε jusqu’à ce qu’on ait la précision désirée. (On peut même améliorer le choix du x
initial dans la méthode de Newton en extrapolant la courbe des xε antérieurs.)
On sait que le problème d’optimisation sans contrainte résolu à chaque itération par la méthode
de Newton est de plus en plus mal conditionné au fur et à mesure que ε décroı̂t. Mais pour une raison
pas totalement expliquée, il semble que le choix d’initialisation de la méthode de Newton suggéré
naturellement immunise l’algorithme contre les effets négatifs de ce mauvais conditionnement.
Cette méthode, due à Nesterov et Nemirovsky sur une idée de Karmarkar, puis améliorée par
S. Boyd, a donné, entre les mains de ce dernier, des résultats spectaculaires sur de très nombreux
problèmes. Elle est la méthode à appliquer . . .quand elle s’applique. (La convexité de u et des fi est
essentielle, et il faut que le calcul des dérivées secondes soit raisonnablement simple.)
54 CHAPITRE 3. OPTIMISATION DANS RN

3.4 Optimisation sous contraintes égalité

On considère à présent le cas où les contraintes sont de la forme

C = {x ∈ Rn | fi (x) = 0 , i = 1, . . . , p} . (3.10)

Il n’y a plus lieu de supposer une quelconque convexité pour les fi ni u parce qu’on sort de toutes
façons du cadre convexe. (C’est pour rappeler cela que nous n’avons pas utilisé la notation C —mais
C— pour l’ensemble des x admissibles.)
A quelques indications près, on laissera le soin au lecteur d’imaginer comment hybrider les algo-
rithmes que nous allons discuter avec ceux du cas inégalité si une partie des contraintes est d’un type
et une partie d’un autre.

3.4.1 Contraintes affines

On considère donc maintenant le cas où les fi dans (3.10) sont affines, ou, de manière équivalente,
il nous est donné une matrice F de type p × n et un vecteur f de dimension p, qui définissent les x
admissibles comme devant vérifier
Fx = f (3.11)
En outre, on supposera toujours que F est de rang p, donc surjective, ce qui impose en particulier que
p < n. En effet, si ce n’est pas le cas, soit f est dans l’image de F , mais alors une partie des contraintes
est redondante : on peut supprimer des lignes qui sont linéairement dépendantes des autres, soit f n’est
pas dans l’image de F , et alors il n’y a pas de x admissible.

Algorithmes de gradient
La variété affine admissible est convexe. Donc l’algorithme de gradient projeté peut être conservé
à l’identique. Il reste juste à faire remarquer que la projection est facile à calculer, au moins s’il n’y a
pas trop de contraintes, et s’exprime par

PC (x) = (I − F † F )x + F † f (3.12)

où
F † = F t (F F t )−1
est une inverse à droite (l’inverse de Penrose) de la matrice surjective F .
En fait on peut faire mieux. En effet, la projection sur Ker F du gradient de u est alors le gradient
de la restriction de u à C. On peut donc appliquer un algorithme de gradient à pas optimal à cette
restriction :
Algorithme Gradient projeté à pas optimal
1. Choisir une estimée initiale x0 , faire k := 0.
2. Calculer ∇u(xk ) et hk := (I − F † F )∇u(xk )
3. Si khk k < ε stop. Si non,
4. Calculer xk+1 := xk − θk hk où le pas θk > 0 est déterminé par

u(xk+1 ) = min u(xk − θhk )

θ
3.4. OPTIMISATION SOUS CONTRAINTES ÉGALITÉ 55

5. incrémenter k de 1 et retourner en 2

Remarque 3.4 Il est prudent d’ajouter une correction pour s’assurer que les xk calculés restent bien
dans la variété admissible, ce qui peut être perdu autrement en raison des erreurs d’arrondi. On
peut intercaller à une fréquence à choisir un pas de projection sur C entre deux pas de gradient, soit
le faire systématiquement à tous les pas, remplaçant la formule xk+1 := xk − θk hk par xk+1 :=
(I − F † F )(xk − θk hk ) + F † f .

Les propriétés de convergence de cet algorithme se déduisent de son interprétation comme algo-
rithme de gradient à pas optimal sur la restriction de u à C.

Algorithme d’Uzawa
On a indiqué plus haut que le théorème de Kuhn et Tucker demeure pour des contraintes égalité
affines, à ceci près que le signe des multiplicateurs n’est plus fixé.
En conséquence, l’algorithme d’Uzawa demeure, en supprimant l’opération de projection sur le
cone positif. La preuve de convergence est inchangée.

Programmation quadratique
Un problème classique, dont on verra qu’il joue un rôle par la suite, est le problème d’optimiser
une forme quadratique sous des contraintes affines. Il s’agit donc de minimiser
1
u(x) = xt Qx + q t x
2
où Q est une matrice symétrique positive définie et q un vecteur de Rn , sous les contraintes (3.11).
On laisse le lecteur démontrer (exercice) que la solution s’obtient conjointement avec le multipli-
cateur de Lagrange p en résolvant le système linéaire
Qx + F t p = −q ,
(3.13)
Fx = f

Théorème 3.12 Le systèmes d’équations (3.13) a une solution unique quelque soit F si et seulement
si F est surjective et la restriction de Q au noyau de F est définie (positive ou négative).

Preuve Il suffit de vérifier si le système homogène, c’est à dire obtenu en remplaçant q et f par 0, a
zéro pour seule solution. Soit donc (x, p) satisfaisant le système homogène. Multiplions la première
équation à gauche par xt et tenons compte de la deuxième pour constater qu’on doit avoir xt Qx = 0.
Mais ce x appartient nécessairement au noyau de F . Donc xt Qx = 0 implique x = 0 si et seulement
si la restriction de Q à ce noyau est définie. Mais dans ce cas, on doit aussi avoir F t p = 0, qui implique
p = 0 si et seulement si F est surjective.
Si Q est inversible, le système (3.13) admet la solution
x∗ = Q−1 F t (F Q−1 F t )−1 (F Q−1 q + f ) + Q−1 q .
Cependant, demander l’inversibilité de Q est trop, puisque seule doit être positive définie sa restriction
au noyau de F . On peut donner une formule explicite qui n’utilise que cette hypothèse, en fonction de
la décomposition en valeurs singulières de F :
 
V1
F = U ΣV = U [Σ1 0]
V2
56 CHAPITRE 3. OPTIMISATION DANS RN

(où U et V sont des matrices orthogonales de type p × p et n × n respectivement et Σ1 est la matrice

diagonale des valeurs singulières, V2t engendre Ker F ) comme

x∗ = (I − V2t (V2 QV2t )−1 V2 Q)V1t Σ−1 t t t −1

1 U f − V2 (V2 QV2 ) V2 q .

En pratique, on a construit des algorithmes d’élimination très spécialisés, plus rapides que le calcul
d’une décomposition en valeurs singulières suivi de l’inversion de V2 QV2t .

3.4.2 Contraintes nonlinéaires

Nous sommes maintenant dans les notations de (3.10).

Algorithme à la Uzawa
Le théorème de Kuhn et Tucker n’est plus valide si les contraintes ne sont pas affines. On peut
tenter d’appliquer encore l’algorithme d’Uzawa. On n’est assuré ni de sa convergence, ni du fait que
s’il converge ce soit vers l’optimum. Voici une brève analyse de cette question.
Remarquons que la fonction étendue

ϕ(x) = sup(λ, f (x))

λ

vaut 0 si x est admissible, et +∞ si non. Ainsi le problème posé, de minimiser u(x) sous les
contraintes f (x) = 0, est équivalent au problème de minimiser u(x) + ϕ(x), soit encore de cher-
cher
minn sup(u(x) + (λ, f (x))) .
x∈R λ

On a encore un minsup, et ceci justifie qu’on fasse un algorithme de type gradient ascendant en λ.
Mais ce que calcule l’algorithme d’Uzawa, c’est le

max minn (u(x) + (λ, f (x))) .

λ x∈R

En général, ces deux quantités sont différentes, la deuxième inférieure à la première, ce qu’on appelle
le “saut de dualité”.
Cette méthode ne doit donc être tentée qu’avec circonspection, et si elle converge il convient de
vérifier si f s’annulle au point trouvé. (Ce qui n’est garanti dans le cas convexe que par le fait que
minx supp = maxp minx .)

Programmation quadratique séquentielle

Si les conditions requises sont satisfaites, c’est à dire si les ∇fi (x∗ ) sont linéairement indépendants,
le point recherché est caractérisé par le théorème de Lagrange. C’est à dire qu’il est, avec le multipli-
cateur (vectoriel) de Lagrange λ, solution du système

∇u(x) + ∇f (x)λ = 0 ,
(3.14)
f (x) = 0.

Ici, ∇f (x) désigne la matrice jacobienne de f transposée (f 0 (x))t , de même que ∇u(x) désigne le
transposé de u0 (x).
3.4. OPTIMISATION SOUS CONTRAINTES ÉGALITÉ 57

Une idée naturelle, et fructueuse, consiste à essyer de résoudre ce système d’équations non linéaires
par la méthode de Newton.
Nous introduisons quelques notations à cet effet. On notera Q := Dx2 (u(x)+(λ, f (x))) la matrice
symétrique des dérivées secondes en x du lagrangien, et F := f 0 (x) la matrice jacobienne de f . En
outre, pour alléger encore les notations, on notera f k := f (xk ) et de même pour toutes les fonctions
de x et λ.
Avec ces notations, l’algorithme de Newton s’écrit

Qk xk+1 + (F k )t λk+1 = Qk xk − ∇uk ,

F k xk+1 = F k xk − f k .

On remarque que ce système, où les membres de droite sont connus à l’étape k, et les inconnues sont
xk+1 et λk+1 , a exactement la même forme que le système (3.13). Il peut donc être résolu à l’aide un
algorithme de programmation quadratique —d’où le nom de cette méthode.
Ajoutons que la théorie de la seconde variation montre que les conditions énoncées au théorème
3.12 sont excatement ici la condition de qualification des contraintes d’une part, et la condition suffi-
sante du deuxième ordre pour un minimu local sous contraintes d’autre part.
On peut donc énoncer le théorème :

Théorème 3.13 Si l’optimum x∗ recherché existe, et si la condition de qualification des contraintes

y est satisfaite ainsi que la condition suffisante locale du deuxième ordre, il existe un voisinage de x∗
dans lequel l’algorithme de programmation quadratique séquentielle converge.

Bien sûr, avec ses qualités —très grande vitesse de convergence—, cet algorithme partage les
défauts de la méthode de Newton : faible robustesse, et nécessité de ce fait de partir avec une assez
bonne estimée de x∗ .
58 CHAPITRE 3. OPTIMISATION DANS RN
Chapitre 4

Programmation linéaire et
programmation dynamique

Les deux sujets qu’effleure ce chapitre, à titre d’introduction, ont en commun de se situer à la
frontière de l’optimisation continue et de l’optimisation combinatoire.
La programmation linéaire parle de variables continues sous des contraintes continues, mais le
premier pas dans l’étude de ce problème est de montrer qu’un nombre fini de points sont candidats
à être optimaux, et l’algorithme du simplexe peut être vu comme une façon habile de parcourir cet
ensemble fini.
De son côté, la programmation dynamique sera d’abord présentée comme un problème de re-
cherche du plus court chemin dans un graphe, donc fondamentalement combinatoire. Mais on verra
ensuite qu’un procédé classique ramène un problème en variable d’espace continue (et variable de
temps discrète) à un tel graphe.

4.1 Programmation linéaire

4.1.1 Position du problème
Bien des modèles en ingéniérie et en économie mènent à considérer des problèmes où critère
et contraintes sont linéaires (donc pas α-convexe) et les variables positives. Ces problèmes ont reçu
le (vieux) nom de “programmes linéaires”. Ils ont été étudiés dès les origines de ce qu’on devait
appeler la “recherche opérationnelle”, notamment par G.B. Dantzig dès le début des années 1950 (le
plus ancien article cité remonte à 1951), et le développement des méthodes numériques afférentes est
contemporain de l’apparition des ordinateurs.
Comme modèles simplifiés très naturels de situations concrètes (coûts et ressources consommées
proportionnels aux quantités), ces problèmes ont justifié un effort algorithmique considérable, et ce
jusque dans les années plus récentes, comme le bruit fait par les Bell labs autour de la “méthode de
Karmarkar” l’a montré.
Nous nous contenterons ici de donner un aperçu des résultats de base et des idées sous-jacentes
à l’algorithme du simplexe, comme introduction à l’utilisation des programmes existants. En effet,
écrire un nouveau programme de programmation linéaire, que ce soit par l’algorithme du simplexe
ou une autre méthode, est un exercice strictement réservé aux professionnels de la chose, tant les
“packages” qui existent sont (nombreux et) perfectionnés.

59
60 CHAPITRE 4. PROGRAMMATION LINÉAIRE ET PROGRAMMATION DYNAMIQUE

On utilisera les inégalités entre vecteurs

x≥y ⇔ xi ≥ yi , i = 1, . . . ,

x>y ⇔ x ≥ y et x 6= y ,

x >> y ⇔ xi > yi , i = 1, . . . .

Le problème peut en général être formulé de la façon suivante.

Le critère à minimiser est déterminé par un vecteur c de Rn , et est donné par

n
X
u(x) = (c, x) = ci xi
i=1

Les contraintes sont d’une part

x ≥ 0,

et d’autre part deux jeux de contraintes définies par deux matrices A et B, de dimensions respectives
p × n et q × n, et deux vecteurs a et b de dimension p et q définissant les contraintes égalité et inégalité
par
Ax = a , Bx ≤ b .

Dans les discussions qui suivent, il faut avoir à l’esprit que n, p et q peuvent être de l’ordre de
plusieurs milliers, voire dizaines de milliers.
La première remarque est qu’au moins au plan théorique, on peut remplacer les contraintes iné-
galité par des contraintes égalité et inversement par les artifices suivants. En introduisant q variables
supplémentaires y ∈ Rq , on peut remplacer Bx ≤ b par

Bx + y = b , y ≥ 0.

Au contraire, si on veut privilégier les contraintes inégalité, on peut remplacer Ax = a par

Ax ≤ a et Ax ≥ a .

Bien sur, la première de ces opérations augmente de q le nombre de variables, tandis que la deuxième
augmente de p le nombre de contraintes. Deux opérations qui ne sont pas souhaitables au vu des
dimensions que nous évoquions. Ce sont par contre des artifices utiles pour étudier les propriétés
théoriques du problème.
Nous utiliserons dans l’étude théorique la forme standard égalité caractérisée par m contraintes
égalité signées :
Ax = b , b ≥ 0,

ce qui est toujours possible, quitte à multiplier certaines contraintes par −1, et toujours les contraintes
de positivité
x ≥ 0.
4.1. PROGRAMMATION LINÉAIRE 61

4.1.2 Étude du polyèdre

Pour comprendre le problème de programmation linéaire, il faut étudier l’ensemble des points
défini par les contraintes linéaires et de positivité, appelé polyèdre, que nous prenons sous la forme
standard égalité Ax = b. Nous noterons P cet ensemble.
Nous appelons encore n la dimension de x, sachant qu’elle a pu être augmentée pour donner cette
forme aux contraintes. Nous appelons m le nombre de contraintes, et supposons que m < n. Nous
supposons même plus que cela, à savoir que

rangA = m . (4.1)

En effet, si-non A a des lignes linéairement dépendantes d’autres, et soit la même dépendance linéaire
se retrouve dans les coordonnées de b, et cette ligne est surnuméraire, et peut être omise sans rien
changer au problème, soit les coordonnées de b n’exhibent pas la même dépendance, et le polyèdre
est vide.
Nous introduisons à cette fin la terminologie suivante :

Définition 4.1 (Direction de Rn ) Nous appelons direction l’ensemble des vecteurs portés par une
même demi-droite, c’est à dire multiples positifs d’un d’entre-eux.

Ainsi, soit h ∈ Rn , il définit une direction {θh | θ ∈ R+ }, et tout vecteur de cette direction définit
la même direction.
Une direction sera réputée “admissible” si elle est composée de vecteurs à coordonnées positives
(ce qui est le cas dès qu’un de ses vecteurs l’est) et si elle appartient au noyau de A. Ainsi, si x ∈ P,
x + w ∈ P pour tout w dans cette direction. Ce sont des directions dans lesquelles P est non borné.
(Très précisément, ces directions définissent des “points à l’infini” du polyèdre au sens de la géométrie
projective.)
Le résultat essentiel que nous visons est le suivant :

Théorème 4.1 Les points d’un polyèdre peuvent tous être obtenus comme combinaison convexe d’un
nombre fini de ses points (appelés sommets) plus une somme d’éléments d’un nombre fini de di-
rections (appelées directions admissibles extrêmales ou sommets à l’infini). Réciproquement, toute
combinaison de cette forme appartient au polyèdre.

Ainsi tout polyèdre est caractérisé par un nombre fini de sommets, à distance finie ou à l’infini, et
est constitué par l’ensemble de leurs combinaisons convexes. En particulier, s’il est borné, un polyèdre
se réduit au polytôpe de toutes les combinaisons convexes de ses sommets (en nombre fini).
Pour démontrer ce résultat, nous introduisons la définition suivante :

Définition 4.2 (Points extrémaux) Étant donné un ensemble convexe C, on appelle points extrémaux
de C les points de C qui ne peuvent être représentés comme une combinaison convexe propre d’autres
points de C.

Ainsi, si x̂ est un point extrémal de C, et si x1 et x2 sont deux points de C tels que x̂ = λx1 + (1 −
λ)x2 , alors nécessairement λ = 0 ou λ = 1, et x̂ coincide avec x1 ou x2 .
La même définition s’appliquera aux directions admissibles, sachant que deux vecteurs colinéaires
et de même sens (proportionnels dans un rapport positif) représentent la même direction.
62 CHAPITRE 4. PROGRAMMATION LINÉAIRE ET PROGRAMMATION DYNAMIQUE

Lemme 4.2 Un point admissible (i.e. tel que x ≥ 0 et Ax = b) [resp une direction admissible h, i.e.
telle que h ≥ 0 et Ah = 0] est extrémal[e] si et seulement si les colonnes de A correspondant aux
coordonnées non nulles de x [resp h] sont linéairement indépendantes [resp. ont un défaut de rang
égal à 1].

Les coordonnées non nulles de x sont dites “de base”. On notera que la propriété ci-dessus im-
plique notamment qu’un point extrémal a au plus m coordonnées de base.
Démonstration du lemme Pour simplifier l’écriture, nous supposons que ce sont les p premières
colonnes de A qui sont concernées, et donc les n − p dernières coordonnées de x qui sont nulles.
Nous partitionnons A et x en
 
x̄
A = [Ā Ã] ,
x̃

Soit un point de P qui satisfait la condition énoncée, c’est à dire que x = (x̄ 0), et donc Āx̄ = b.
Si ce point est combinaison convexe propre de deux autres points de P alors, il est intérieur au segment
joignant ces deux points, ce qui veut dire qu’il existe un vecteur w 6= 0 tel que x + w et x − w
appartiennent à C. En particulier, ceci implique que les composantes de w hors base soient nulles (si
non, soit x + w soit x − w aurait des composantes négatives), et donc aussi que Āw̄ = 0, où w̄ désigne
bien sûr les m premières coordonnées de w. Mais par hypothèse, Ā est injective, donc w = 0, ce qui
contredit l’hypothèse.
Réciproquement, soit x un point de P, x̄ l’ensemble de ses coordonnées non nulles, que nous
regroupons au début de la numérotation, et Ā la matrice des colonnes corrsepondantes. Si les colonnes
de Ā ne sont pas linéairement indépendantes (A n’est pas injective), il existe un vecteur w̄ non nul
tel que Āw̄ = 0. Comme les coordonnées de x̄ sont toutes strictement positives, il existe ε assez petit
pour que, en choisissant kw̄k ≤ ε ce qui est toujours loisible, les coordonnées de x̄ + w̄ et de x̄ − w̄
soient encore toutes positives. Alors, en complètant w = (w̄ 0), on voit que Aw = 0, et que x + w
et x − w appartiennent tous les deux à C. Donc x = 1/2(x + w) + 1/2(x − w) n’est pas un point
extrémal de P.
On laisse le lecteur répéter la même preuve pour les directions admissibles extrémales, en se sou-
venant que la différence entre le nombre des vecteurs considéré et le rang du système qu’ils forment,
ou défaut de rang, est la dimension du noyau de la matrice dont ils sont les colonnes.
Cette caractérisation montre qu’il ne saurait y avoir qu’un nombre fini de points extrémaux à
P. Il suffit de tester tous les ensembles de m colonnes de A, et pour ceux qui sont linéairement
indépendants, de verifier si A−1 b est positif. On fera de même avec les ensembles de m + 1 colonnes
de rang m pour chercher les directions admissibles optimales.
Démonstration du théorème Soit x un point de C, et supposons que ce ne soit pas un point extrémal.
Comme précédement, il existe un vecteur w du noyau de A ayant les mêmes coordonnées nulles que
x, et tel que x − w et x + w appartiennent à C. Regardons x + tw, t > 0. Si w n’est pas une direction
admissible de C (il a des coordonnées négatives), pour un certain t+ une des coordonnées de ce vecteur
s’annule, les autres étant encore positives. On peut faire de même avec x − tw. (Si w est une direction
admissible de C, −w ne l’est pas, et réciproquement.) Donc, soit t+ et t− sont tous les deux définis,
et x peut être représenté comme une combinaison convexe de deux vecteurs qui ont une coordonnée
de plus nulle chacun : x = t− /(t+ + t− )x+ + t+ /(t+ + t− )x− , soit on a une représentation comme
une somme x = x− + t− w où w est une direction admissible, et x− a une coordonnée nulle de plus
que x. En répétant ce processus pour les éléments de cette représentation récursivement, on obtient le
théorème. (On ne fait qu’un nombre fini de fois cette opération, puisque le nombre de coordonnées
4.1. PROGRAMMATION LINÉAIRE 63

non nulles diminue chaque fois, et la construction s’arrête quand A n’a plus de noyau —ou un noyau
de dimension 1 pour les directions admissibles—.)
On déduit de ce théorème le corollaire fondamental suivant :
Corollaire 4.3 Si le polyèdre des contraintes n’est pas vide, soit le critère n’a pas d’infimum fini, soit
un des points extrémaux (ou sommets) du polyèdre est solution.

Démonstration Soit il existe une direction admissible w telle que (c, w) < 0. Alors, comme on peut
ajouter tw, t positif arbitraire, à tout point du polyèdre sans en sortir, clairement, le critère (c, x) peut
être rendu arbitrairement grand négatif. Soit, pour toute direction admissible w, (c, w) ≥ 0. Alors, si
un point du polyèdre a w dans sa décomposition comme au théorème ci-dessus, on peut retirer cette
composante. On reste dans le polyèdre, et on améliore le critère. Nous ne considérons donc que des
points combinaison convexe des sommets x1 à xN :
N
X N
X
x= λ i xi , λi ≥ 0, λi = 1 .
i=1 i−1

Alors,
N
X
(c, x) = λi (c, xi ) .
i=1
Si le critère (c, x) est minimum sur C, tous les (c, xi ) sont supérieurs ou égaux à (c, x). Donc seuls
peuvent être non nuls les λi pour lesquels on a l’égalité, et ce sont là des sommets solution. Ou bien —
et c’est le cas générique— un seul sommet est solution, et il ne peut être représenté comme ci-dessus
de manière non dégénérée.

4.1.3 L’algorithme du simplexe

Nous ne donnons ici qu’un bref aperçu du principe de l’algorithme le plus célèbre pour résoudre
numériquement le programme linéaire, l’algorithme du simplexe. Il en existe d’autres aujourd’hui,
plus rapides ... dans la plupart des configurations.
Remarquons d’abord qu’en principe, puisque le nombre de sommets est fini et qu’on sait les
déterminer tous, il suffit de comparer la valeur du critère à tous ces sommets et prendre le meilleur.
Ce qui s’oppose à cette approche naı̈ve est le nombre potentiellement très grand de sommets du
polyèdre. On considère couramment des problèmes où n et m sont tous les deux en milliers. Or il
y a n!/m!(n − m)! combinaisons de m colonnes à tester. C’est à dire, si n = 2000 et m = 1000, un
nombre de l’ordre de 10600 combinaisons.
Il faut faire autrechose. L’algorithme du simplexe va explorer des sommets d’une façon moins
naı̈ve, et d’habitude plus efficace. On sait malheureusement exhiber des problèmes pour lesquels cet
algorithme explore tous les sommets. Mais on sait que ce n’est génériquement pas le cas, et que le
simplexe est génériquement polynomial en m + n.
Le principe de l’algorithme est donc le suivant. Étant donné le choix d’une “base”, c’est à dire de
m colonnes indépendantes de A formant Ā, de sorte que A = [Ā Ã], et un découpage correspondant
des coordonnées de tout x en coordonnées en base x̄ et hors base x̃, on a nécessairement

Āx̄ + Ãx̃ = b ,

soit encore
x̄ = −Ā−1 Ãx̃ + Ā−1 b . (4.2)
64 CHAPITRE 4. PROGRAMMATION LINÉAIRE ET PROGRAMMATION DYNAMIQUE

Donc, en paramétrisant x via x̃ :

(c, x) = (c̃t − c̄t Ā−1 Ã)x̃ + c̄t Ā−1 b ,

que nous réecrivons avec une notation évidente

(c, x) = (w̃, x̃) + c̄t Ā−1 b .

Si l’itérée xk de l’algorithme est un sommet correspondant à cette base, soit x̃k = 0, on va chercher
à améliorer le critère en repérant la coordonnée de w̃ la plus négative, disons w̃M et en donnant une
valeur positive à la coordonnée x̃M correspondante de x. Soit donc eM le vecteur de base numéro M
de Rn , et essayons des x de la forme
x = xk + teM .

On est sûr de faire décroitre le critère ce faisant. Le t maximum permis est atteint la première fois
qu’une autre coordonnée de x̄, tel que calculé par (4.2) passe par zéro. On s’arrête à cette valeur de t,
on fait ainsi rentrer la coordonnée M dans la base et sortir celle qui s’est annulée. Et on itère.
L’algorithme s’arrête quand toutes les coordonnées de w̃ sont positives : on ne peut plus améliorer
le critère.
Il reste à dire comment initialiser l’algorithme. En effet, nous avons indiqué comment passer d’un
sommet du polyèdre à un autre, mais comment trouver un premier sommet ? Nous allons indiquer
comment utiliser le même algorithme pour trouver un sommet initial, et en même temps déterminer
s’il existe un tel sommet, c’est à dire si l’ensemble des états admissibles est non vide. Il peut en effet se
faire que les contraintes Ax = b, x ≥ 0 soient incompatibles, mais ceci même est difficile à découvrir.
La méthode va consister à examiner un autre problème linéaire qui présente la particularité d’avoir
un sommet évident, et de chercher, s’il existe, un sommet du problème d’origine.
Supposons qu’on s’est ramené à b ≥ 0, ce qu’on peut toujours faire en changeant le signe des
lignes de A si nécessaire. Considérons le problème de programmation linéaire portant sur les variables
(positives) (x, w), x ∈ Rn , w ∈ Rm , dont les contraintes sont

Ax + w = b

et le critère
m
X
u(x, w) = wi .
i=1

Comme promis, les contraintes admettent une solution évidente, qui est un sommet : x = 0, w = b.
Et comme ce critère est toujours positif ou nul, son optimum sera zéro si et seulement s’il existe
une solution des contraintes avec w = 0, soit un point admissible du problème d’origine. En outre,
évoluant de point extrémal en point extrémal, le simplexe nous donnera un point extrémal, qui pourra
à son tour servir d’initialisation pour le problème d’origine.
Il y a beaucoup à dire à partir de là (on a écrit des livres entiers). Que faire dans le simplexe si
la nouvelle base n’est pas indépendante ? si deux coordonnées s’annulent en même temps ? etc. On
laisse le lecteur imaginer des parades simples à ces questions.
Plus intéressant : comment simplifier le calcul de Ā−1 pour la nouvelle base en utilisant le fait
qu’on connaissait cette inverse pour une matrice ne différant que par une colonne ?
Pour toutes ces questions, nous renvoyons le lecteur aux livres spécialisés.
4.1. PROGRAMMATION LINÉAIRE 65

4.1.4 Rudiments de dualité

On ne peut pas parler de programmation linéaire sans évoquer la dualité, qui constitue, comme on
va le voir, un outil très utile.
Dans ce numéro, et pour l’élégance des formules obtenues, nous allons supposer qu’on cherche à
maximiser un critère linéaire u = (c, x). Cela revient évidemment à changer c en −c, mais comme
nous n’avons jamais fait d’hypothèse sur le signe des éléments de c, cela est sans conséquence.
Partons donc d’un problème standard, que nous écrivons

max ct x , x ≥ 0, Ax = b ,
x

en notant de façon explicite ct x le produit scalaire (c, x). Supposons que nous connaissions un vecteur
y de Rp (p est le nombre de contraintes) tel que

y t A ≥ ct .

Alors, pour tout x ≥ 0, on a ct x ≤ y t Ax. Mais par ailleurs, pour tout x admissible, Ax = b, de sorte
que l’inégalité ci-dessus donne
ct x ≤ y t b . (4.3)
Avant d’aller plus loin, supposons que la contrainte Ax = b provient en fait d’une contrainte
inégalité, et qu’on a augmenté l’état de “variables d’écart” pour en faire des contraintes égalité, c’est
à dire, avec un abus de notations évident, que le vecteur x ci-dessus est en fait de la forme
 
x
x= ,
ξ

que A est donc de la forme

A = [ A I ],
de sorte que la contrainte est

Ax + ξ = b, ξ ≥ 0, x ≥ 0,

donc équivalente à
Ax ≥ b , x ≥ 0. (4.4)
La “contrainte” sur y se lit alors
y t [ A I ] ≥ [ ct 0 ]
soit
y t A ≥ ct , y≥0 (4.5)
Ainsi, pour tout x admissible au sens de (4.4) et tout y admissible au sens du problème dual dont
l’admissibilité est définie par (4.5), on a (4.3). Le problème

min(b, y)
y

soumis aux contraintes (4.5) est appelé problème dual de celui d’origine (rappelons que nous en avons
ici fait un problème de maximisation sous contraintes inégalité).
On remarque la parfaite symétrie entre problème primal et dual, de sorte que toute affirmation
concernant leur intéraction peut être énoncée dans un sens ou dans l’autre.
66 CHAPITRE 4. PROGRAMMATION LINÉAIRE ET PROGRAMMATION DYNAMIQUE

On voit que si on trouve x et y admissibles pour les problèmes primal et dual tels que (c, x) =
(b, y), alors ils sont nécessairement solution de ces problèmes. Cette remarque est la base de méthodes
efficaces pour résoudre le problème de programmation linéaire, visant à minimiser la différence, tou-
jours positive ou nulle, (b, y) − (c, x) parmi les x et y admissibles.
On voit aussi aisément que si le problème dual a un ensemble d’états admissibles non vide, le
problème primal a son critère borné supérieurement (et mutatis mutandis pour le problème dual si le
problème primal admet des états admissibles).
Ces constatations élémentaires ont des réciproques, que nous ne démontrerons pas :

Théorème 4.4 Chacun des deux problèmes, primal et dual, a un ensemble d’états admissibles non
vide et un suprémum fini (donc une solution) si et seulement si il en va de même de l’autre. Si un des
deux problèmes a un extremum infini, l’autre n’a pas d’état admissible.

Donc les deux problèmes ont une solution ou n’en ont pas simultanément, l’absence de solution
pouvant provenir de l’absence d’états admissibles, ou du fait que le critère n’est pas borné.

Théorème 4.5 Si les problèmes de programmation linéaire primal et dual ont une solution, il existe
x∗ et y ∗ tels que (c, x∗ ) = (b, y ∗ ). Réciproquement, toute paire admissible (x∗ , y ∗ ) satisfaisant cette
égalité est optimale.

Enfin, faisons remarquer que la connaissance de y ∗ , par exemple, permet de trouver x∗ . En effet,
remarquons qu’on a donc toujours, pour tout x et y admissibles

ct x ≤ y t Ax ≤ y t b ,

de sorte que si ct x = y t b, non seulement x et y coincident avec les optimums, mais aussi, toutes les
inégalités ci-dessus sont des égalités. Or, l’égalité

(y ∗ , Ax∗ − b) = 0

alors que y ∗ ≥ 0 et Ax−b ≤ 0 montre que pour toute coordonnée de y ∗ non nulle, la contrainte corres-
pondante en x est “saturée” en x∗ (satisfaite avec l’égalité). De même, l’égalité (ct − (y ∗ )t A)x∗ = 0,
que nous pouvons réécrire
(At y ∗ − c, x∗ ) = 0
alors que At y ∗ − c ≥ 0 et x∗ ≥ 0 implique que pour toute contrainte duale non saturée, la coordonnée
corrsepondante de x∗ est nulle. Comme nous l’avons vu plus tôt, connaı̂tre les contraintes saturées
(les ξi nuls) et les xi nuls suffit à déterminer x∗ par l’équation linéaire qu’on en déduit.
L’utilité de cette théorie est double. D’une part, comme nous l’avons indiqué, elle est le fondement
de méthodes numériques efficaces, ou d’améliorations de l’algorithme du simplexe. D’autre part,
elle permet toujours de choisir si l’on préfère résoudre le problème primal ou dual. L’un a plus de
contraintes (inégalité) que de variables, et l’autre plus de variables que de contraintes. Suivant les
packages de résolution disponibles, l’un peut être préférable à l’autre.

4.2 Programmation dynamique

Nous allons partir d’une vision combinatoire de la programmation dynamique, pour aboutir à une
utilisation dans des problèmes authentiquements “dynamiques” en variables continues, donc voisins
des préoccupations du reste de ce cours.
4.2. PROGRAMMATION DYNAMIQUE 67


( XX 3
 (( XXXX 
((
3 ((((
(( (
A X
 
XXX
XX4X
A
Z XX 2
A J
B Z
2 J1
1
B Z HH
A 1
B Z 2
H
HH
A J
B4
Z
A
  HH J
Z H J
Z
A
2 B
B
Z 4 1

Z B  Z  



Z3 B Z
 
Z ZZ 
Z
B  4 2 
2
Z B  Z 
Z  BB   Z
 
Z ZZ 
Z  1
 

F IGURE 4.1 – Graphe et “poids” des arcs

4.2.1 Plus court chemin dans un graphe orienté

Le problème le plus simple

Nous considérons un graphe orienté, ici de gauche à droite, comme celui de la figure 1, dont
chaque arc est muni d’une “longueur” ou “poids”.
Le problème posé est de trouver le chemin de “longueur” ou “poids” minimum du nœud initial, le
plus à gauche, au nœud terminal, le plus à droite. C’est manifestement un problème “fini” : le graphe
ne comporte qu’une vingtaine de chemins. On pourrait donc tous les lister et choisir le plus court.
Mais on voit bien que le nombre de chemins croit combinatoirement avec la taille du graphe, et une
procédure aussi rustique ne s’étendra pas à des situations beaucoup plus complexes. (Le seul travail
de répertorier tous les chemins devient exorbitant.)
Nous allons indiquer une procédure qui ne croit que linéairement avec le nombre de nœuds, ou
plus précisément comme le produit du nombre de nœuds par le nombre moyen d’arcs par nœud.
Le principe de la méthode est de marquer chaque nœud avec la longueur du chemin minimal de
ce nœud jusqu’à la fin. Cette procédure sera rapide en raison de la remarque banale mais essentielle
qui suit :

Proposition 4.6 (Principe de Bellman) Le chemin optimal a la propriété (dite “principe d’optima-
lité” de Bellman) 1 qu’entre tout nœud N par où il passe et la fin du chemin, il est optimal pour le
problème d’aller de ce nœud N jusqu’à la fin. (Ce que nous appellerons le sous-problème initialisé
en N .)

Démonstration Comme la longueur totale du chemin est la somme de la longueur parcourue du nœud
initial au nœud N plus la longueur de N jusqu’à la fin, si on pouvait trouver un chemin plus court
pour cette dernière longueur, —le sous problème initialisé en N — on pourrait, en le concaténant avec

1. Il s’agit de Richard Bellman, dont on peut contester qu’il ait inventé la programmation dynamique, pas qu’il en ait
compris le premier toute la puissance et toute la généralité.
68 CHAPITRE 4. PROGRAMMATION LINÉAIRE ET PROGRAMMATION DYNAMIQUE

le chemin optimal entre le début et N , trouver un chemin global plus court que le chemin optimal, ce
qui est une contradiction.
Cette démonstration semble être une tautologie tant la propriété est évidente. Pourtant, nous allons
nous servir de ce “principe de Bellman” en le reformulant un peu.
Pour tout nœud, si le chemin optimal passe par ce nœud, de ce nœud jusqu’à la fin, il utilise le
chemin optimal pour ce sous problème. En particulier, depuis un nœud N , si la longueur du che-
min optimal jusquà la fin —c’est à dire la valeur optimale des sous problèmes si non leur solution
complète— est connue pour tous les nœuds immédiatement aval (c’est à dire séparés par un seul arc),
alors résoudre le sous problème initialisé en N est immédiat. En effet, si le chemin optimal (depuis
N ) passe par un certain N 0 aval, la longueur en est la somme de la longueur de l’arc séparant N de
N 0 ajoutée à la valeur optimale du sous-problème initialisé en N 0 . Ainsi, depuis N il suffit de com-
parer ces valeurs, et de retenir la meilleure (la plus petite). On aura ainsi à peu de frais la valeur du
sous-problème initialisé en N .
On obtient ainsi l’algorithme suivant :
Algorithme Programmation dynamique simple
1. Marquer le nœud terminal avec la valeur 0.
2. En tout nœud dont tous les nœuds immédiatement aval sont déjà marqués, faire :
– pour chaque nœud immédiatement aval, calculer la somme de la longueur de l’arc
vers ce nœud et de la valeur de ce nœud.
– Prendre la plus petite de ces sommes pour valeur du nœud courant, et le marquer.
– Marquer le, ou les, arc(s) donnant la valeur retenue.
3. Retourner en (2) jusquà ce que tous les nœuds soient marqués
4. depuis le nœud initial, (comme depuis tout nœud du graphe) tout chemin n’empruntant
que des arcs marqués est optimal, et a une longueur égale à la valeur marquée à ce
nœud.
À titre d’exemple, nous avons dans la figure 2 marqué à chaque nœud sa valeur, et renforcé les
arcs optimaux. On voit que le problème posé n’avait pas une solution unique, mais cette procédure
n’en est nullement affectée.

Extensions du problème du plus court chemin

On peut étendre de nombreuses façons cet algorithme. La plus élémentaire est la suivante. On
peut considérer plusieurs nœuds terminaux et plusieurs nœuds initiaux possibles. De plus, on peut
supposer qu’un “poids” est attaché à chacun des nœuds en plus des arcs.
Dans ce dernier cas, on pourrait aussi bien attacher ce poids du nœud à tout arc qui le rejoint, se
ramenant de manière triviale au problème où seuls les arcs ont un poids. On préfère, pour des raisons
qui apparaı̂tront plus loin, considérer d’une part le poids attaché à chaque nœud terminal, et considérer
que c’est la valeur de ce nœud pour l’algorithme de programmation dynamique, et associer les poids
des autres nœuds aux arcs qui les quittent.
On aboutit ainsi à l’algorithme suivant.
Algorithme Programmation dynamique
1. Marquer les nœuds terminaux avec leur valeur donnée.
2. En tout nœud dont tous les nœuds immédiatement aval sont déjà marqués, faire :
4.2. PROGRAMMATION DYNAMIQUE 69


( 3 XX 3
 (( XXXX 
((
3 ((((
(( (
A X
5 XX 1
 X XX4X
A 
Z XX 2
A J
A2 J1
B Z
1
B Z HH
1
B Z 2
H
HH
A J
B4
Z
A
  HH J
Z H J
Z
A
2 B Z 4 1
6 B 5 1 0
Z B  Z   


Z3 B Z 
 
Z B Z Z 4

Z  2
Z B 2  Z 
Z 
Z
BB  Z 
Z Z  1
Z 
3 2
 

F IGURE 4.2 – Le graphe de la figure 1 avec les valeurs et les arcs optimaux

– pour chaque nœud immédiatement aval, calculer la somme de la longueur de l’arc

vers ce nœud et de la valeur de ce nœud.
– Prendre la plus petite de ces sommes pour valeur du nœud courant, et le marquer.
– Marquer le, ou les, arc(s) donnant la valeur retenue.
3. Retourner en (2) jusquà ce que tous les nœuds soient marqués
4. Tout chemin partant d’un nœud initial de valeur minimum et n’empruntant que des arcs
marqués est optimal, et a une longueur égale à la valeur marquée à ce nœud.
Dans l’exemple de la figure 3, qui a deux nœuds initiaux possibles et trois nœuds terminaux, on
a seulement attaché des poids aux nœuds terminaux, puisque des poids sur les nœuds intermédiaires
auraient seulement augmenté d’autant le poids de chaque arc les quittant, sans augmenter la généralité
de l’exemple.
On donne directement le graphe marqué avec les valeurs et les chemins optimaux. On conseille
au lecteur de refaire l’algorithme lui-même pour constater combien il est simple et rapide. Pourtant ce
graphe a plus de 110 chemins possibles.
De très nombreux problèmes combinatoires peuvent se ramener à un problème de recherche de
chemin de poids minimal dans un graphe. La caractéristique essentielle pour mettre à jour une telle
structure, et l’exploiter, est le sens de parcours unique. En général il provient de ce que le problème
peut être organisé en “étapes” dont l’ordre est imposé par la nature du problème, ou immatériel quant
à la solution cherhée de sorte qu’il peut être fixé arbitrairement. Cet aspect d’étapes successives, ou
dynamique, va être examiné maintenant.

4.2.2 Système dynamique et programmation dynamique

Système dynamique
Nous examinons maintenant un cas particulier extrêmement important. Supposons que les nœuds
du graphe, appelés ici états, peuvent être repérés par un numéro d’étape k ∈ 0, 1, . . . , N , et pour
chaque étape k soit Xk l’ensemble des états possibles à cette étape. L’hypothèse ici est que les arcs
70 CHAPITRE 4. PROGRAMMATION LINÉAIRE ET PROGRAMMATION DYNAMIQUE


( 3 XX 3
 (( XXXX 
((
3 ((((
(( (
A X
5 XX 1 H
 X XX4X
A  H 0 
Z XX 2
A J H
A2 J 1 
1
B Z
1
B Z HH
1
2 H A
2 
B Z
HH J 
B4  J
Z
A 
 2 Z 

H
H 
A   J
2 BB
ZZ 4 H  1
6 5 1 H 0

Z B  
Z   
3 HH


Z3 B 
 Z  
 H
 H
Z  B 
 Z HH

3 Z 
Z 4 3 2 
5 XX Z B 2  Z  1


 XXXZ  BB 
Z
Z
Z  3
Z 
2 X XZ  1 

3 2
 

F IGURE 4.3 – Un graphe à plusieurs nœuds initiaux et terminaux

relient toujours un état d’une étape à un de l’étape suivante, et que ceci correspond au sens de parcours
imposé du graphe. Indiçons les arcs issus d’un nœud (k, x(k)) —où x(k) ∈ Xk est un nœud de l’étape
k— par un indice u ∈ U(k, x(k)) appelé commande. U(k, x(k)) est simplement un ensemble dont
le cardinal est égal au nombre d’arcs quittant le nœud (k, x(k)). On voit que le graphe définit une
équation dynamique de la forme

x(k + 1) = f (k, x(k), u(k)) (4.6)

puisque pour chaque nœud (k, x(k)) et chaque commande u ∈ U(k, x(k)), il définit à quel nœud ou
état conduit cet arc, état toujours situé à l’étape k + 1.
Un chemin dans le graphe est une suite d’états {x(k), k = 0, . . . , N }, appelée trajectoire. Une
trajectoire peut aussi être caractérisée par l’état initial x(0) et une suite de commandes {u(k), k =
0, . . . , N − 1}, qui, via l’équation (4.6), engendre une trajectoire unique.
Notons encore L(k, x, u) le “poids” ou la longueur —nous dirons ici le coût— de l’arc issu du
nœud (k, x) indicé par u, et pour tout nœud terminal x ∈ XN , K(x) le coût attaché à ce nœud. Ainsi,
le coût d’une trajectoire, qu’il s’agit de minimiser, est donné par
N
X −1
J = K(x(N )) + L(k, x(k), u(k)) . (4.7)
k=0

Il y a équivalence complète entre un graphe structuré comme on l’a dit et une équation de la
forme (4.6), et donc entre un problème de plus court chemin dans un tel graphe et le problème de
minimisation de J donné par (4.7) avec la dynamique (4.6). Bien des problèmes seront formulés sous
la forme (4.6),(4.7). Un bon moyen de les résoudre est de faire l’analogie avec le graphe, et de faire
l’algorithme de programmation dynamique sur ce graphe.
Le terme même de “programmation dynamique” vient de là. L’équation (4.6) définit ce qu’on
appelle un système dynamique. L’indice k est généralement interprété comme représentant le temps.
L’équation (4.7) définit une fonctionnelle additive de la trajectoire. On recherche la commande, et
éventuellement l’état initial, qui minimise ce critère ou coût.
4.2. PROGRAMMATION DYNAMIQUE 71

Une forme équationnelle de la programmation dynamique

Ayant décrit le graphe et le critère par des équations, on peut décrire dans ce langage l’algorithme
de la programmation dynamique. Notons V (k, x) le “marquage” associé au nœud (k, x). On l’appelle
plutot ici la fonction “performance”, ou “de Bellman”. L’algorithme que nous avons décrit s’écrit
alors :

∀k ∈ {0, . . . , N − 1} , ∀x ∈ Xk , V (k, x) = min [L(k, x, u) + V (k + 1, f (k, x, u))] , (4.8)

u∈U(k,x)

∀x ∈ XN , V (N, x) = K(x) . (4.9)

L’algorithme de la programmation dynamique consiste donc à appliquer la formule (4.8) pour
calculer V de proche en proche, en commençant par initialiser V avec (4.9), puis en reculant en k. Il
faut avoir deux tableaux des valeurs de V (k, ·) en mémoire, celui qu’on est en train de remplir et celui
qui est utilisé dans le deuxième membre de (4.8). In fine, le tableau des V (0, x) donne pour tout état
initial possible le coût minimum possible.
Il faut aussi en même temps remplir un grand tableau des valeurs de u qui donnent le minimum
à chaque (k, x). Ce tableau donne la stratégie (ou commande en boucle fermée) optimale, en ce que
pour chaque état (k, x) il donne la commande optimale si on se trouve en cet état. Cet aspect est
particulièrement utile si l’équation (4.6) constituait une approximation d’un phénomène physique, de
sorte qu’on est susceptible de constater au cours de la mise en œuvre qu’on est dans un état (k, x(k))
différent de ce que à quoi on s’attendait. L’algorithme ci-dessus a donné une commande conseillée
pour tout état dans le graphe. Certes, l’écart entre le phénomène réel et le modèle fait que cette com-
mande n’est plus tout à fait optimale, mais si cet écart est faible, elle a toutes les chances de rester une
“bonne” commande.

Système à état et continu

Tout le développement de la programmation dynamique a été fait en termes de graphe, avec donc
l’hypothèse constante que dans l’équation (4.6), x et u prennent leurs valeurs dans des ensembles Xk
et U(k, x) finis. Cependant, ces équations ainsi que (4.7) gardent un sens si ces ensembles sont des
sous ensembles de Rn et Rm , disons, respectivement. C’est à dire qu’alors l’état est constitué de n
nombres réels et la commande de m nombres réels, les uns et les autres éventuellement bornés.
Les équations de la programmation dynamique (4.8)(4.9) gardent également un sens. Et on va
démontrer le résultat suivant :

Théorème 4.7 S’il existe une fonction rélle V (·, ·) satisfaisant les équations (4.8) et (4.9), en désignant
par ϕ(k, x) un argument du minimum dans (4.8), si la commande u(k) = ϕ(k, x(k)) est admissible
pour x0 (au sens où elle engendre une trajectoire qui respecte les contraintes x(k) ∈ Xk , ce dont
on peut s’assurer par un choix convenable des U(k, x)), alors cette commande est optimale pour le
problème défini par (4.6)(4.7) initialisé en x(0) = x0 .

Démonstration Soit {u(0), u(1), . . . , u(N − 1)} une commande admissible, engendrant une trajec-
toire {x0 , x(1), . . . , x(N )}. En tout point de cette trajectoire, d’après (4.8), on a

V (k, x(k)) ≤ L(k, x(k), u(k)) + V (k + 1, f (k, x(k), u(k))) ,

ou encore
V (k, x(k)) − V (k + 1, x(k + 1)) ≤ L(k, x(k), u(k)) .
72 CHAPITRE 4. PROGRAMMATION LINÉAIRE ET PROGRAMMATION DYNAMIQUE

Sommons cette inégalité de k = 0 à N − 1, il vient

N
X −1
V (0, x0 ) − V (N, x(N )) ≤ L(k, x(k), u(k)) .
k=0

Utilisons alors (4.9) pour exprimer V (N, x(N )), que nous faisons repasser à droite, il reste
N
X −1
V (0, x0 ) ≤ K(x(N )) + L(k, x(k), u(k)) ,
k=0

soit
V (0, x0 ) ≤ J(x0 , {u}) . (4.10)
Maintenant, si la suite des {u(k)} coincide pour tout k avec ϕ(k, x(k)), les inégalités dans les calculs
ci-dessus sont toutes remplacées par des égalités, et on conclu que

V (0, x0 ) = J(x0 , ϕ) . (4.11)

La comparaison des relations (4.10) et (4.11) établit le théorème.

Système en temps continu

Dans la pratique, le système (4.6) est souvent issu de la discrétisation en temps d’un système en
temps continu, ou système différentiel, de la forme

ẋ = F (t, x, u) .

Si on discrétise ce système avec un pas de temps h, en notant xk = x(kh) et uk = u(kh), on a au

premier ordre
xk+1 = xk + hF (kh, xk , uk )
qui est bien une équation de la forme (4.6), avec

f (k, x, u) = x + hF (kh, x, u) .

De même, le système différentiel peut être muni d’un critère intégral à minimiser, de la forme
Z T
J = K(x(T )) + l(t, x(t), u(t)) dt
0

qui peut être approximé au premier ordre par une somme finie (où T = N h)
N
X −1
J = K(xN ) + hl(kh, xk , uk )
k=0

qui est bien de la forme (4.7).

Ainsi, la programmation dynamique apparait comme un moyen d’aborder l’optimisation d’un
critère intégral pour un système différentiel, un problème connu sous le nom de “commande optima-
le”, ou en Franglais de “contrôle”.
4.2. PROGRAMMATION DYNAMIQUE 73

L’approximation au premier ordre proposée ci-dessus n’est convenable que si on choisit un pas
de temps “assez petit”. De même, l’application pratique demandera souvent qu’on discrétise aussi x
et u, se ramenant en fait à un problème fini. Et encore, cette discrétisation elle-même demande à être
faite avec soin. Enfin, ces discrétisations ne mèneront à un problème faisable en pratique que si les
dimensions des espaces d’état et de commandes sont assez petites pour que le nombre de “nœuds” du
graphe soit raisonnable.
En fait, ce que nous obtenons ici est une discrétisation de l’équation de Hamilton Jacobi Bellman,
l’équation aux dérivées partielles équivalente pour ce problème continu à l’équation de Bellman pour
le problème à temps discret. Une analyse plus approfondie des schémas numériques nécessiterait
beaucoup plus de mathématiques. Nous nous limitons donc à cette aproche naı̈ve.
Mais quelles que soient les limitations de cette méthode, elle reste extrêmement utile dans des
cas où elle s’applique. En particulier par le fait qu’elle se prête à prendre en compte toutes sortes
de contraintes sur les états admissibles, et des données sans bonnes propriétés mathématiques. (Par
exemple, les données peuvent dépendre de fonctions tabulées, etc.)