Résolution d'équations avec CoffeeScript

Dans cet article on va résoudre l'équation x3=50 par dichotomie et par l'algorithme de Raphson, en mode "normal" et en mode "multiprécision". Pour ce dernier cas, on utilisera la bibliothèque Big.js. Mais d'abord on va voir comment la notion de dichotomie peut être introduite en dehors du contexte des résolutions d'équations, avec un jeu qui peut être vu en seconde comme activité préliminaire:

Le jeu de la valoche

Présentation du jeu

Chaque matin, un animateur radio annonce un montant aléatoire dont on sait uniquement qu'il est inférieur ou égal à 100 €, et donné par un nombre entier. Puis il téléphone à un auditeur au hasard et lui demande s'il connaît le montant de la "valoche". Mais l'auditeur n'était à ce moment précis, pas en train d'écouter la radio. L'animateur lui propose de deviner en quelques essais, et lui dit si ses tentatives tombent en dessous ou au dessus du montant réel. En combien de tentatives maximum peut-on être sûr de deviner le montant ?

L'algorithme est intéressant parce que

Et si le temps le permet, on peut demander aux élèves de programmer un robot pour qu'il gagne au jeu s'il est contacté par l'animateur (on peut aussi demander d'enrichir l'algorithme ci-dessus pour qu'il affiche le nombre de tentatives):

Inversion des rôles

Cette version également, est intéressante, puisqu'en plus des aspects précédents, elle incite assez naturellement à chercher quels sont les nombres qui nécessitent beaucoup d'essais avec la dichotomie. Et permet d'explorer les nombres maximaux autres que 100.

Résolution d'équation par dichotomie

Avec Coffeescript

Pour résoudre l'équation x3=50, on va la réécrire x3-50=0, et déterminer deux nombres a et b tels que a3-50<0 et b3-50>0. La dichotomie permet de rapprocher a et b jusqu'à ce qu'ils soient tous deux proches de la solution de l'équation. En renommant ƒ la fonction x3-50, on économise du code et on peut facilement modifier le script pour résoudre d'autres équations (attention toutefois au cas où ƒ est décroissante, il y a plus de choses à changer).

Attention à ne pas aller au-delà de 10-15, la précision de la machine ne permettant pas de distinguer a et b dans ce cas. Pour faire mieux il faut utiliser une bibliothèque de précision comme Big.js:

Résolution de l'équation en grande précision

On voit que ça prend un certain temps pour trouver que la racine cubique de 50 vaut environ 3,6840314986403866057798228335798072219172581174381826692564614945772359755087770427815054848656019769...

Algorithme de Newton-Raphson

Le principe est que la résolution de l'équation ƒ(x)=0 revient à résoudre x-ƒ(x)/ƒ'(x)=x puisque ces équations sont équivalentes. Or pour résoudre une équation du type g(x)=x, on peut simplement itérer g: Le théorème des accroissements finis permet de prouver que, plus la dérivée de g est petite, et plus la convergence est rapide. Or, en la solution de l'équation g(x)=x, la dérivée de g est nulle, on peut difficilement faire plus petit...

Le problème avec la méthode de Raphson est qu'elle nécessite la dérivée de ƒ: En calculant cette dérivée par un algorithme, on court le double risque de rallonger le temps de calcul, et de perdre en précision. On a donc tout intérêt à calculer la dérivée par calcul formel avant de rédiger l'algorithme. Dans le cas présent, ƒ'(x)=3x² donc g(x)=x-(x3-50)/(3x²)=2x/3+50/(3x²). En fait il est possible de décomposer la suite en deux suites adjacentes. On va utiliser la suite des itérés de x-ƒ(x)/ƒ'(x) pour calculer rapidement la racine cubique de 50: On choisit une valeur initiale quelconque (par exemple u0=1) et on pose un+1=2un/3+50/(3un2):

Avec CoffeeScript

La convergence est assez rapide pour permettre d'explorer empiriquement (déterminer par exemple qu'il faut boucler 11 fois pour avoir 15 décimales correctes). L'avantage des suites adjacentes est qu'on a une condition de sortie de boucle (comme dans la dichotomie ci-dessus).

En grande précision

On recommence mais en version Big.js:

On boucle 11 fois, on affiche la valeur de u trouvée, puis on calcule g(u)-u, et on affiche le résultat: On constate qu'il n'est pas tout-à-fait nul. Alors on essaye avec une autre valeur que 11, et (rapidement si on opère par dichotomie comme dans la première partie !) on trouve que 12 tours de boucle suffisent pour avoir la racine cubique à 10-99 près. En commençant par Big.DP = 1000, il faut moins d'une minute de dichotomie sur le nombre de boucles pour déterminer que 16 tours de boucle donnent la racine cubique de 50 à 1000 décimales près:

3,6840314986403866057798228335798072219172581174381826692564614945772359755087770427815054848656019772860747983678083780522156700088223538177216810146750284113828284571623221900563081617469139352929220190907374863012262111864041050751715759460820960781632078837083253346646354147835146855588917409805541998473269280043103799096042474408909658828725533094320987865961933203374632268723094722834639987079803810816015336357513187484518197229575763388142873072828323977028712073361152104104452595910283657000624189150839340881551495850374065352139642214191968398590792406127926965628277910390858895146284586904456039658612496852464581800173508811409229032019467813550329572419934668327725869638735344743140753075599138089049989752830287261897227622172445223435264854105519837917257265114987324240550792456823332877430417318356003669713158512113084960418844448070440191051639513346565793687092991581756975404821394174627036698035099783373544660522364050690770510456247295512737504458235303880227798220929657