I. IntroductionEn informatique, la programmation dynamique est une méthode algorithmique pour résoudre des problèmes d'optimisation. Le concept a été introduit au début des années 1950 par Richard Bellman.
La programmation dynamique consiste à résoudre un problème en le décomposant en sous-problèmes, puis à résoudre les sous-problèmes, des plus petits aux plus grands en stockant les résultats intermédiaires.
Pour mieux comprendre cette méthode algorithmique on partira d'un exemple simple de calcul du nème terme de la suite de Fibonacci, puis on cherchera à l'aide de la programmation dynamique à optimiser le calcul des différences divisées.
Enfin, on montrera comment évaluer un polynôme d'interpolation de Newton en mémorisant ses différences divisées successives.
Note importante :
Par la suite, pour faciliter leur comptage, les appels récursifs des fonctions seront représentés sur les lignes d'un triangle (partiel ou entier), ces lignes seront numérotées de 0 à n.
La position des appels situés sur une même ligne ne respectera pas forcément leur chronologie.
II. Principe de la méthode - Suite de Fibonacci
Prenons comme exemple le calcul du nème terme de la suite de Fibonacci.
II-A. Définition de la suite de Fibonacci
Elle est définie par la relation de récurrence :
Fi = Fi-1 + Fi-2
Avec comme premiers termes :
F0 = 0
F1 = 1
II-B. Fonction récursive - Méthode intuitive
Partant de la définition précédente, on peut alors écrire une fonction récursive en Python qui évalue le nème terme de cette suite :
| Code Python : | Sélectionner tout |
1 2 3 4 5 | def fib_v1(n): if n <= 1: # si 1er ou 2e terme : n=0 ou n=1 return n # retourne n else: # sinon, appels récursifs basés sur la relation de récurrence return fib_v1(n-1) + fib_v1(n-2) |
Représentons maintenant à l'aide d'un schéma les différents appels récursifs de fib(5) :
Note importante : pour mieux visualiser les choses on a regroupé les appels récursifs identiques sur chaque ligne.
On a donc besoin de réaliser 15 appels (1 appel initial et 14 appels récursifs) pour évaluer fib(5), les calculs étant effectués à la remontée des appels.
L'algorithme est en fait en temps exponentiel.
On constate qu'il y a :
- 3 appels de fib(0) ;
- 5 appels de fib(1) ;
- 3 appels pour fib(2) ;
- 2 appels pour fib(3).
On va donc chercher à mémoriser les valeurs de fib(i) pour éviter de les recalculer par la suite.
II-C. Programmation dynamique - Méthode ascendante
Dans la méthode ascendante, on va "de bas en haut" : il s'agit du sens contraire à celui utilisé par les méthodes récursives. On commence par calculer des solutions pour les sous-problèmes les plus petits, puis, de proche en proche, on calcule les solutions des problèmes en utilisant le principe d'optimalité et on mémorise les résultats dans un tableau.
f[0] = 0
f[1] = 1
f[2] = f[1] + f[0] = 1 + 0 = 1
f[3] = f[2] + f[1] = 1 + 1 = 2
...
f[n] = f[n-1] + f[n-2]
On obtient alors la fonction itérative :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 | def fib_v2(n): if n==0: # si 1er terme : n=0 return n # retourne n et sort # initialisation de la liste des termes avec des zéros : f = [0, 0, ..., 0] f=[0]*(n+1) # initialisation des 2 pemiers termes de la liste : f = [0, 1, ..., 0] # f[0]= 0 inutile car on a déjà : f[0] = 0 f[1] = 1 for i in range(2, n+1): # parcours des indices : 2, ..., n f[i] = f[i-1] + f[i-2] # mémorisation du terme d'indice i : Fi = Fi-1 + Fi-2 return f[n] # renvoie le terme d'indice n : Fn # return f # retourne la séquence de Fibonacci : [0, 1, 1, 2, ..., Fn] |
Note : on peut également retourner la séquence de Fibonacci f au lieu du nombre f[n].
Représentons maintenant l'enchaînement des calculs avec mémorisation des résultats intermédiaires pour fib(5) :
On constate que le calcul de Fib(n) est cette fois en temps linéaire.
Note importante :
Dans notre exemple d'optimisation on utilise la méthode dite « ascendante », mais il existe également une autre méthode dite « descendante » nécessitant d'écrire un algorithme récursif.
III. Évaluation des différences divisées
III-A. Définition
Les différences divisées d'ordre 0, 1 et 2 sont définies comme suit :
On obtient ainsi le schéma récursif suivant :
Note importante :
Elles interviennent notamment dans le polynôme d'interpolation de Newton N(x) associé à k+1 points (x0, y0), …,(xk, yk) et défini par :
N(x) = [y0] + [y0, y1](x - x0) + [y0, y1, y2](x - x0)(x - x1) + ... + [y0, y1, ..., yk](x - x0)(x - x1)...(x - xk-1)
III-B. Fonction récursive - Méthode intuitive
En utilisant la relation de récurrence :
[yi, ..., yi+n] = ([yi+1, ..., yi+n] - [yi, ..., yi+n-1])/(xi+n - xi)
avec :
[yi]= yi
On peut alors écrire la fonction récursive :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 | def difference_divisee_v1(x , y, i, n): # x, y : listes des valeurs # i, n : indice mini et ordre de la différence divisée [yi, yi+1, ..., yi+n] if n==0: # si on est à l'ordre 0 return y[i] # on renvoie yi : [yi] = yi else: # sinon, appels récursifs selon la relation de récurrence [yi, yi+1, ..., yi+n] = ([yi+1, yi+1, ..., yi+n] - [yi, yi+1, ..., yi+n-1]) / (xi+n - xi) return (difference_divisee_v1(x , y , i+1, n-1) - difference_divisee_v1(x , y , i, n-1))/(x[i+n] - x[i]) |
Représentons maintenant à l'aide d'un schéma les différents appels récursifs nécessaires pour évaluer [y0, ..., y3] :
On a donc besoin de réaliser 23+1 - 1 = 15 appels (1 appel initial et 14 appels récursifs) pour évaluer [y0, ..., y3], les calculs étants effectués à la remontée des appels.
Par conséquent, l'évaluation des différences divisées à l'ordre n nécessite 2n+1 - 1 appels, et "seulement" 2n - 1 calculs car les valeurs des différences divisées d'ordre 0 (ligne n) sont déjà connues.
Donc l'algorithme est en temps exponentiel en n.
On constate également qu'il y a 3 appels de [y1] et [y2], et 2 pour [y1, y2].
On va donc chercher à mémoriser les valeurs de [yi, ..., yi+k] pour éviter de les recalculer par la suite.
III-C. Programmation dynamique - Méthode ascendante
On reprend la relation de récurrence :
[yi, ..., yi+n] = ([yi+1, ..., yi+n] - [yi, ..., yi+n-1])/(xi+n - xi)
avec :
[yi]= yi
Si maintenant on part de la base, on peut alors évaluer la différence divisée en remontant le triangle ligne par ligne et en mémorisant à chaque fois les résultats intermédiaires :
On a donc besoin de réaliser 3×4/2 = 6 calculs pour évaluer [y0, ..., y3].
Donc le calcul des différences divisées à l'ordre n nécessite cette fois seulement n(n+1)/2 opérations.
On obtient ainsi la fonction itérative :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 | def difference_divisee_v2(x , y , j, n): # x, y : listes des valeurs # j, n : indice mini et ordre de la différence divisée [yj, yj+1, ..., yj+n] dy = y.copy() # on copie la liste y dans dy : différences divisées à l'ordre 0 for k in range(1,n+1): # on parcourt les ordres suivants des différences divisées : ordre 1, 2, ..., n for i in range(j, j + n-k+1): # on parcourt les indices de x et dy jusqu'à j+n-k h = (x[i+k] - x[i]) # on détermine la différence entre xi+k et xi pour le calcul de la différence divisée à l'ordre k dyi = (dy[i+1] - dy[i])/h # on évalue la différence divisée à l'ordre k : [yi, ..., yi+k]) dy[i] = dyi # on remplace l'élément d'indice i de la liste dy par la valeur de la différence divisée dyi return dy[j] # on renvoie la différence divisée |
Note : on a voulu également optimiser l'espace mémoire en utilisant une liste à une seule dimension pour enregistrer les résultats intermédiaires.
IV. Évaluation du polynôme d'interpolation de Newton
Nous avons vu précédemment que le polynôme d'interpolation de Newton N(x) associé à k+1 points (x0, y0), …,(xk, yk) est défini par :
N(x) = [y0] + [y0, y1](x - x0) + [y0, y1, y2](x - x0)(x - x1) + ... + [y0, y1, ..., yk](x - x0)(x - x1)...(x - xk-1)
Pour éviter de recalculer à chaque fois les différences divisées successives [y0], [y0, y1], [y0, y1, y2], ..., on doit d'abord les mémoriser en remontant le triangle ligne par ligne :
On peut ainsi écrire une fonction Python qui génère la liste des différences divisées à partir des 2 listes de valeurs x et y :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 | def differences_divisees(x,y): n=len(x)-1 # nombre de valeurs de la série - 1, il correspond au degré du polynôme, exemple pour n=2 : a2x^2 + a1x + a0, où a2, a1 et a0 représentent les coefficients dy = y.copy() # on copie la liste y dans dy : différences divisées à l'ordre 0 diff_divisees = [dy[0]] # on initialise la liste permettant d'enregistrer les différences divisées : [y0] for k in range(1,n+1): # on parcourt les ordres suivants des différences divisées : ordre 1, 2, ..., n for i in range(n-k+1): # on parcourt les indices de x et dy jusqu'à n-k h = (x[i+k] - x[i]) # on détermine la différence entre xi+k et xi pour le calcul de la différence divisée à l'ordre k dyi = (dy[i+1] - dy[i])/h # on évalue la différence divisée à l'ordre k : [yi, ..., yi+k]) dy[i] = dyi # on remplace l'élément d'indice i de la liste dy par la valeur de la différence divisée dyi diff_divisees = diff_divisees + [dy[0]] # on ajoute à la liste la différence divisée d'ordre k. return diff_divisees # on renvoie la liste des différences divisées |
Évaluons maintenant le polynôme d'interpolation en se basant sur sa formule générale.
La fonction evaluer_polynome_newton prend en arguments 2 listes x et y, représentant la série de points (xi,yi), et renvoie la valeur du polynôme pour valeur_x :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 | def evaluer_polynome_newton(x, y, valeur_x): diff_divisees = differences_divisees(x,y) # on génère la liste des différences divisées : [[y0], [y0, y1], [y0, y1, y2], ...] k=len(x)-1 # indice maxi de la liste des x valeur_polynome_newton = y[0] # on initialise la variable avec y[0] produit_facteurs=1 # on initialise le produit de facteurs for i in range(k): # parcours des indices : 0, 1, ..., k-1 produit_facteurs = produit_facteurs*(valeur_x - x[i]) # on évalue le produit de facteurs : (x-x0)(x-x1)...(x-xi) valeur_polynome_newton = valeur_polynome_newton + diff_divisees[i+1]*produit_facteurs # on évalue progressivement le polynôme de Newton return valeur_polynome_newton # renvoie la valeur du polynôme de Newton |
V. Module complet de test
Il contient l'ensemble des fonctions décrites dans ce billet et le code permettant de comparer les versions récursives et optimisées :
| Code Python : | Sélectionner tout |
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 | import time def fib_v1(n): if n<=1: # si 1er ou 2e terme : n=0 ou n=1 return n # retourne n else: # sinon, appels récursifs basés sur la relation de récurrence return fib_v1(n-1) + fib_v1(n-2) def fib_v2(n): if n==0: # si 1er terme : n=0 return n # retourne n et sort # initialisation de la liste des termes avec des zéros : f = [0, 0, ..., 0] f=[0]*(n+1) # initialisation des 2 pemiers termes de la liste : f = [0, 1, ..., 0] # f[0]= 0 inutile car on a déjà : f[0] = 0 f[1] = 1 for i in range(2, n+1): # parcours des indices : 2, ..., n f[i] = f[i-1] + f[i-2] # mémorisation du terme d'indice i : Fi = Fi-1 + Fi-2 return f[n] # renvoie le terme d'indice n : Fn # return f # retourne la séquence de Fibonacci : [0, 1, 1, 2, ..., Fn] def difference_divisee_v1(x , y, i, n): # x, y : listes des valeurs # i, n : indice mini et ordre de la différence divisée [yi, yi+1, ..., yi+n] if n==0: # si on est à l'ordre 0 return y[i] # on renvoie yi : [yi] = yi else: # sinon, appels récursifs selon la relation de récurrence [yi, yi+1, ..., yi+n] = ([yi+1, yi+1, ..., yi+n] - [yi, yi+1, ..., yi+n-1]) / (xi+n - xi) return (difference_divisee_v1(x , y , i+1, n-1) - difference_divisee_v1(x , y , i, n-1))/(x[i+n] - x[i]) def difference_divisee_v2(x , y , j, n): # x, y : listes des valeurs # j, n : indice mini et ordre de la différence divisée [yj, yj+1, ..., yj+n] dy = y.copy() # on copie la liste y dans dy : différences divisées à l'ordre 0 for k in range(1,n+1): # on parcourt les ordres suivants des différences divisées : ordre 1, 2, ..., n for i in range(j, j + n-k+1): # on parcourt les indices de x et dy jusqu'à j+n-k h = (x[i+k] - x[i]) # on détermine la différence entre xi+k et xi pour le calcul de la différence divisée à l'ordre k dyi = (dy[i+1] - dy[i])/h # on évalue la différence divisée à l'ordre k : [yi, ..., yi+k]) dy[i] = dyi # on remplace l'élément d'indice i de la liste dy par la valeur de la différence divisée dyi return dy[j] # on renvoie la différence divisée def differences_divisees(x,y): n=len(x)-1 # nombre de valeurs de la série - 1, il correspond au degré du polynôme, exemple pour n=2 : a2x^2 + a1x + a0, où a2, a1 et a0 représentent les coefficients dy = y.copy() # on copie la liste y dans dy : différences divisées à l'ordre 0 diff_divisees = [dy[0]] # on initialise la liste permettant d'enregistrer les différences divisées : [y0] for k in range(1,n+1): # on parcourt les ordres suivants des différences divisées : ordre 1, 2, ..., n for i in range(n-k+1): # on parcourt les indices de x et dy jusqu'à n-k h = (x[i+k] - x[i]) # on détermine la différence entre xi+k et xi pour le calcul de la différence divisée à l'ordre k dyi = (dy[i+1] - dy[i])/h # on évalue la différence divisée à l'ordre k : [yi, ..., yi+k]) dy[i] = dyi # on remplace l'élément d'indice i de la liste dy par la valeur de la différence divisée dyi diff_divisees = diff_divisees + [dy[0]] # on ajoute à la liste la différence divisée d'ordre k. return diff_divisees # on renvoie la liste des différences divisées def evaluer_polynome_newton_v1(x, y, valeur_x): k=len(x)-1 # indice maxi de la liste des x valeur_polynome_newton = y[0] # on initialise la variable avec y[0] produit_facteurs =1 # on initialise le produit de facteurs for i in range(k): # parcours des indices : 0, 1, ..., k-1 produit_facteurs = produit_facteurs*(valeur_x - x[i]) # on évalue le produit de facteurs : (x-x0)(x-x1)...(x-xi) valeur_polynome_newton = valeur_polynome_newton + difference_divisee_v1(x,y,0,i+1)*produit_facteurs # on évalue progressivement le polynôme de Newton return valeur_polynome_newton # renvoie la valeur du polynôme de Newton def evaluer_polynome_newton_v2(x, y, valeur_x): diff_divisees = differences_divisees(x,y) # on génère la liste des différences divisées : [[y0], [y0, y1], [y0, y1, y2], ...] k=len(x)-1 # indice maxi de la liste des x valeur_polynome_newton =y[0] # on initialise la variable avec y[0] produit_facteurs=1 # on initialise le produit de facteurs for i in range(k): # parcours des indices : 0, 1, ..., k-1 produit_facteurs = produit_facteurs*(valeur_x - x[i]) # on évalue le produit de facteurs : (x-x0)(x-x1)...(x-xi) valeur_polynome_newton = valeur_polynome_newton + diff_divisees[i+1]*produit_facteurs # on évalue progressivement le polynôme de Newton return valeur_polynome_newton # renvoie la valeur du polynôme de Newton n=40 print("Comparaison de la durée d'exécution des 2 fonctions fib_v1 et fib_v2") print("") print("Calcul de fib(" + str(n) + ") - version n°1 :") start_time = time.time() fib_n = fib_v1(n) print("fibo_v1(" + str(n) + ") = " + str(fib_n)) print("Durée d'exécution : %s secondes" % (time.time() - start_time)) print("") print("Calcul de fib(" + str(n) + ") - version n°2 :") start_time = time.time() fib_n = fib_v2(n) print("fibo_v2(" + str(n) + ") = " + str(fib_n)) print("Durée d'exécution : %s secondes" % (time.time() - start_time)) |
Le code affiche :
Comparaison de la durée d'exécution des 2 fonctions fib_v1 et fib_v2
Calcul de fib(40) - version n°1 :
fibo_v1(40) = 102334155
Durée d'exécution : 108.11277770996094 secondes
Calcul de fib(40) - version n°2 :
fibo_v2(40) = 102334155
Durée d'exécution : 0.01198887825012207 secondes
VI. Conclusion
En partant d'un exemple simple consistant à évaluer le nème terme de la suite de Fibonacci, nous avons mieux compris le principe de la programmation dynamique.
Ceci nous a ensuite permis d'optimiser le calcul des différences divisées et celui du polynôme d'interpolation de Newton.
Enfin, nous avons pu vérifier sur nos fonctions que cette méthode algorithmique apporte un net gain de temps d'exécution.
Sources :
https://fr.wikipedia.org/wiki/Programmation_dynamique
https://en.wikipedia.org/wiki/Dynamic_programming
https://fr.wikipedia.org/wiki/Suite_de_Fibonacci
https://fr.wikipedia.org/wiki/Diff%C..._divis%C3%A9es
https://fr.wikipedia.org/wiki/Triangle_de_Pascal