une boucle ? une suite !
Qu'est-ce qu'une boucle?
Une héroine d'un conte pour enfants ?
Ou bien ceci ?
r ← r0 pour k de 1 jusque n faire r ← f(r) fpour retourne r
Mais non, bien sûr, C'est une suite:
$$r_{k+1}=f(r_k)=f^{k+1}(r_0) \ {\rm avec } \ r_0=r0$$
L'algorithme peut d'ailleurs être écrit récursivement:
boucle(r : type 1, n : entier) : type 2 si n <= 0 alors retourne r sinon retourne boucle(f(r), n - 1) fsi
| Orbite d'un système dynamique |
|---|
| L'orbite issue de $r_0$ pour la fonction $f$ est la suite $(r_n)$ définie par $r_n=f^n(r_0)$ |
En programmation fonctionnelle, on traduit une boucle par un pliage. Par exemple avec Haskell.
Prelude> let boucle invariantBoucle ini ni nf pas = foldl (\ r k -> invariantBoucle r) ini [ni, ni + pas .. nf] Prelude> boucle (\ x -> 2*x) 1 1 10 1 1024
On peut avoir l'orbite avec la fonction iterate qui renvoie une liste infinie : c'est bien l'orbite au sens mathématique du terme.
let orbite = iterate
Pour n'avoir que les n premiers éléments de l'orbite, on utilise la fonction take:
Prelude> take 11 ( orbite (\ x -> 2*x) 1 ) [1,2,4,8,16,32,64,128,256,512,1024]
On peut avoir presque la même chose sur XCAS.
On utilise l'opérateur de composition itérée: f@@n est alors l'itérée n fois de la composée de f par elle-même.
0>> orbite := (f,r0,n) -> (f@@n)(r0) 1>> f := x -> 2*x 2>> orbite(f,1,10) 1024
Le programmeur voudrait savoir si sa boucle va effectivement renvoyer ce qui est
prévu dans sa spécification et surtout en combien de temps: s'il s'agit de
calculer la vitesse d'une voiture en fonction de la distance qui la sépare de la
voiture qui est en face, il faut faire vite.
exemple d'étude: l'algorithme de Babylone
Il y a presque 4000 ans, les petits Babyloniens calculaient la racine carrée de
2 avec 6 bonnes décimales en utilisant un algorithme élémentaire.
On peut lire $1 + \frac{24}{60} + \frac{51}{60^2}+ \frac{10}{60^3}\approx 1,414 213$
r ← 1.0 pour k de 1 jusque n faire r ← (r + 2/r) * 0.5 fpour retourne r
Vous pouvez le traduire instantanément dans n'importe quelle langage avec une boucle pour. Évitez autant que faire se peut les
tant que pour ne pas risquer de lancer des boucles infinies.
Avec XCAS, c'est immédiat: programmer en XCAS c'est écrire en pseudo-langage...mais pour de vrai !
boucle(n):={ r := 1; pour k de 1 jusque n faire r := (r + 2/r) * 0.5 fpour; retourne r }:;
En récursif:
boucle r n si n <= 0 alors retourne r sinon retourne boucle((r + 2 / r) * 0.5, n - 1) fsi
La traduction sur XCAS est carément du mot à mot:
boucle(r,n):={ si n <= 0 alors retourne r sinon retourne boucle((r + 2 / r) * 0.5, n - 1); fsi }:;
On peut aussi utiliser notre fonction boucle créée avec Haskell:
Prelude> boucle (\x -> (x + 2 / x) * 0.5) 1 1 6 1 1.414213562373095
Premier problème : est-ce que la suite des valeurs calculées par la boucle
converge vers $\sqrt{2}$?###
Nous aurons besoin de deux théorèmes que nous admettrons et dont nous ne
donnerons que des versions simplifiées.
|Théorème de la limite monotone|
|---------------------------------------------|
|Toute suite croissante majorée (ou décroissante minorée) converge|
Un théorème fondamental est celui du point fixe. Il faudrait plutôt parler des
théorèmes du point fixe car il en existe de très nombreux avatars qui portent
les noms de mathématiciens renommés: Banach, Borel,
Brouwer, Kakutani, Kleene,…Celui qui nous
intéresserait le plus en informatique est celui de Knaster-Tarski. On
peut en trouver une présentation claire dans
Les démonstrations et les algorithmes.
Ils donnent tous des conditions d'existence de points fixes (des solutions de
$f(x)=x$ pour une certaine fonction).
Nous nous contenterons de la version light vue au lycée.
|Un théorème light du point fixe|
|------------------------------------------|
|Soit I un intervalle fermé de $\mathbb R$, soit $f$ une fonction de Ι vers I et soit
$(r_n)$ une suite d'éléments de I telle que $r_{n+1}=f(r_n)$ pour tout entier
naturel $n$. |
|SI $(r_n)$ est convergente ALORS sa limite est UΝ point fixe de $f$
appartenant à I.|
Cela va nous aider à étudier notre suite définie par $r_{n+1}= f(r_n)=\frac{r_n +
\frac{2}{r_n}}{2}$ et $r_0=1$.
On peut alors démontrer par récurrence que
- pour tout entier naturel non nul $n$, on a $r_n \geqslant \sqrt{2}$
- la suite est décroissante
- $\sqrt{2}$ est l'unique point fixe positif de $f$
On peut alors conclure et être assuré que notre algorithme va nous permettre d'approcher $\sqrt{2}$.
Deuxième problème: quelle est la vitesse de convergence?
Ici, cela se
traduit par « combien de décimales obtient-t-on en plus à chaque
itération?».
Introduisons la notion d'ordre d'une suite:
| Définition: Ordre d'une suite - Constante asymptotique d'erreur |
|---|
| Soit $(r_n)$ une suite convergeant vers $\ell$. S'il existe un entier $k>0$ tel que: |
$$
\displaystyle\lim_{n\to +\infty}\frac{|r_{n+1}-\ell|}{|r_n-\ell|^k}=C
$$
avec $C\neq 0$ et $C\neq +\infty$ alors on dit que $(r_n)$ est d'ordre $k$ et
que C est la constante asymptotique d'erreur.|
Déterminons l'ordre et la constante asymptotique d'erreur de la suite
de Babylone donnant une approximation de $\sqrt{2}$.
On démontre que
$$
\left|r_{n+1} - \sqrt{2}\right| = \left| \frac{(r_n-\sqrt{2})^2}{2r_n}\right|
$$
Ainsi la suite est d'ordre 2 et sa constante asymptotique est $\frac{1}{2\sqrt{2}}$.
Ici, on peut même avoir une idée de la vitesse sans étude asymptotique. En effet, on a $r_n \geqslant 1$ donc
$$
\left|r_{n+1} - \sqrt{2}\right| \leqslant \left| (r_n-\sqrt{2})^2\right|
$$
Les élèves de terminale étant fort habitué(e)s à manipuler des pH et des concentrations de $H_3O^+$, on peut leur demander de démontrer
que le nombre de décimales d'un nombre positif plus petit que 1 est $-\log_{10}x$.
En posant $d_n=-\log_{10}|r_n-\sqrt{2}|$ on obtient donc que $d_{n+1} \geqslant 2d_n$: à chaque tour de boucle, on double au minimum le
nombre de bonnes décimales.
La précision de la machine étant de 16 décimales si les nombres sont codés en binaires sur 64 bits, il faut au maximum 6 itérations pour
obtenir la précision maximale.
Ou, si vous préférez avec XCAS:
racine(a,xk,n):={ si n == 0 alors xk sinon racine(a,(xk + a / xk) * 0.5, n - 1); fsi }:;
ou toujours avec XCAS/giac mais en ligne de commande:
0>> DIGITS := 16:; "Done" // Time 0 1>> racine(a,xk,n) := si n <= 0 alors xk; sinon racine(a,(xk + a / xk) * 0.5, n - 1); fsi; // Interprète racine // Succès lors de la compilation racine racine: définition récursive "Done" // Time 0 2>> racine(2,1,5) 1.414213562373095 // Time 0
ou encore plus simplement avec l'opérateur de composition itérée @@.
Cela correspond donc au pliage de la programmation fonctionnelle : on peut tout faire avec XCAS...
Ainsi:
0>> DIGITS := 16 1>> f := x -> (x + 2/x) * 0.5 2>> racine := n -> (f@@n)(1) 3>> racine(6) 1.414213562373095
Avec Python:
def racine(a,xk,n): """ Calcule une approximation de a à partir d'une première valeur par défaut xk On effectue n itérations de l'algorithme de Babylone (ou de Héron) """ if n == 0: return xk else: return racine(a, (xk + a / xk) * 0.5 , n - 1)
Une application : on peut ainsi calculer rapidement sur un petit système
embarqué (voire un smart phone) la vitesse de la
voiture en fonction de la distance la séparant de la voiture qui la précède.
Voici par exemple une application disponible sur iPhone:
En multiprécision
Il se peut que vous ayez besoin d'aller au-delà de la précision de la machine (environs 16 décimales pour 64 bits).
Vous pouvez vous tourner vers des logiciels travaillant en précision infinie naturellement comme par exemple XCAS/giac. Il suffit
alors de fixer DIGITS à une valeur voulue:
0>> DIGITS := 100:; 1>> racine(a,xk,n):= si n <= 0 alors xk; sinon racine(a,(xk + a / xk) * 0.5, n - 1); fsi; 2>> racine(2,1,8) // avec le doublement de la précision, 8 devrait suffire : 1.414213562373095048801688724209698078569671875376948073176679737990732478462107038850387534327641573 // on compare avec la valeur donnée par l'approximation native de XCAS: 3>> evalf(sqrt(2),100) 1.414213562373095048801688724209698078569671875376948073176679737990732478462107038850387534327641573
En fait, XCAS/giac utilise la bibliothèque C MPFR de calcul multi-précision mise au point par une équipe française et qui est devenue la référence mondiale.
On peut l'utiliser via Sage dont la syntaxe de programmation est celle de Python. On précise la taille de la mantisse
utilisée dans la représentation des nombres: 53 en 64 bits ou 113 en 128 bits.
On travaille sur 128 bits (on donne la précision = taille de la mantisse) RN = RealField(113) # utilise MPFR de C def racine(a,xk,n): """ Calcule une approximation de a à partir d'une première valeur par défaut xk On effectue n itérations de l'algorithme de Babylone (ou de Héron) On travaille sur une mantisse de prec bits fixée globalement via la bibliothèque MPFR """ if n == 0: return RN(xk) else: return racine(RN(a), RN(xk + a / xk) * RN(0.5) , n - 1)
7 itérations suffisent:
sage: racine(2,1,7) 1.41421356237309504880168872420970 sage: sqrt(2).n(113) 1.41421356237309504880168872420970
On pourra se référer à ce poly de cours pour plus de précisions sur le traitement des nombres à virgule flottante.