Les différentes erreurs lorsqu'on calcule à l'aide d'un ordinateur
Derivation numérique
Dans ce document, on s’intéresse au différentes erreurs qui interviennent lorsqu’on effectue des calculs à l’aide d’un ordinateur. On prend pour exemple le calcul numérique de la dérivée d’une fonction, qu’elle soit définie analytiquement ou numériquement (par une table de valeurs).
Dérivation d’une fonction analytique — idée simple
Introduction
Le nombre dérivé $f’(x)$ en un point $x \in \mathbb{R}$, d’une fonction1 :
$$
\begin{array}{l|rcl}
f: & \mathbb{R} & \longrightarrow & \mathbb{R} \cr
& x & \longmapsto & f(x)
\end{array}
$$
est défini par :
$$f^{\prime}(x) = \lim\limits_{h \to 0} \dfrac{f(x+h) - f(x)}{h}$$
Un ordinateur ne pouvant manipuler qu’un nombre limité de chiffres, tous les nombres envisageables en mathématique ne le sont pas en informatique. La notion de limite n’a en particulier aucune signification en informatique ; on recherche donc une valeur approchée de la dérivée en $x$.
Le procédé le plus simple consiste à considérer que si $h$ est suffisamment petit il n’est pas nécessaire de passer à la limite. On considère alors que :
$$ f^{\prime}(x) \cong \dfrac{f(x+h) - f(x)}{h} \qquad \qquad \text{(1.1.)}$$
Comment évaluer l’approximation réalisée ?
L’approximation (1.1.) s’appuie sur la valeur du paramètre $h$. Il semble évident que cette valeur doit être la plus petite possible de façon à limiter le plus possible la non prise en compte de la limite dans la formule (1.1.). Cette intuition est cependant fausse, nous allons voir que :
la valeur de $h$ ne doit pas être plus petite que le plus petit nombre représentable dans le langage choisi ;
la valeur de $h$ doit être assez grande pour que $f(x)$ et $f(x+h)$ soient reconnues comme deux valeurs différentes, afin d’être à l’abri d’erreurs d’arrondi ;
la valeur de $h$ doit être assez petite afin que les erreurs de troncature soient négligeables.
En Python, on accède facilement à certaines de ces informations :
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
>>>importsys>>>sys.float_infosys.float_info(max=1.7976931348623157e+308,max_exp=1024,max_10_exp=308,min=2.2250738585072014e-308,min_exp=-1021,min_10_exp=-307,dig=15,mant_dig=53,epsilon=2.220446049250313e-16,radix=2,rounds=1)>>>sys.float_info.min# Plus petit float (positif) accessible2.2250738585072014e-308>>>1./1.2250738585072014e-3088.162773150824861e+307>>>sys.float_info.epsilon# plus petit intervalle entre deux float 2.220446049250313e-16>>>1./1e-309# Nombre non représentableinf# Résultat d'une division par 0
Calcul numérique et erreurs
Les opérations numériques ne sont jamais exactes et toujours entachées d’erreurs. Ces erreurs sont généralement petites pour chaque opération mais, comme un ordinateur réalise un très grand nombre d’opérations, leur accumulation peut conduire à des résultats complètement erronés.
Erreur d’arrondi
Les erreurs d’arrondi sont liées à la représentation en mémoire des nombres : tout nombre réel n’est connu qu’avec une précision donnée de $n$ chiffres significatifs.
Afin d’illustrer la notion d’erreur d’arrondi, on suppose, dans cette partie, que les calculs sont effectués à l’aide d’un calculateur à trois chiffres significatifs.
Si $f(x) = x^2$ on calcule que $f’(7)=14,0$. L’utilisation de l’expression approchée (1.1.) donne :
Une petite perturbation sur la valeur de $f(x)$ (erreur d’arrondi) peut induire une grande variation de $f’(x)$.
Si on note $\delta$ la précision relative de la machine (par exemple : 7 chiffres significatifs donnent $\delta = 10^{-7}$), l’erreur absolue sur l’évaluation d’une fonction $f$ en un point $x$ est de l’ordre de $\delta\ |f(x)|$. On peut donc estimer l’erreur sur le numérateur de l’expression (1.1.) par : $\delta\ |f(x+h)| + \delta\ |f(x)| \approx 2 \delta\ |f(x)|$.
L’erreur absolue sur l’évaluation du nombre dérivé par l’expression (1.1.) est donc :
$$
E_a = 2 \delta\ |\dfrac{f(x)}{h}|
$$
Retour à l’exemple ($\delta = 10^{-3}$ et $f(x)=x^2$) :
# -*- coding: utf-8 -*-"""
Mise en évidence de l'erreur d'arrondi lors de la dérivation
"""importmatplotlib.pyplotaspltfromtypingimportListdefdf(f:"function",a:float,h:float)->float:"""
Retourne le nombre dérivé de la fonction f à l'abscisse
a en appliquant le schéma des différences finies à
droite avec le pas h.
"""return(f(a+h)-f(a))/hdefaffichage(x:List[float],y:List[float],titre:str,x_label:str,y_label:str,save:bool=False)->None:"""
Trace erreur = f(pas).
"""fig=plt.figure()ax=fig.add_subplot(1,1,1)ax.set_title(titre)ax.set_xlabel(x_label)ax.set_ylabel(y_label)ax.plot(x,y)plt.grid()if(save):plt.savefig("Graphique-{}.png".format(compteur))plt.show()defmain():"""
Fonction principale
"""f=lambdax:x**2# fonction à dérivera=2# nombre dérivé est calculé en x = 2f_prime=4# valeur théorique du nombre dérivéhs=[0.05/2**iforiinrange(0,45)]# valeurs des pas hdfs=[df(f,a,h)forhinhs]# valeurs des nombres dérivéserreurs=[abs(df-f_prime)/f_prime*100fordfindfs]# évolution de l'erreur#affichage(dxs, df, "Valeur de la dérivée")affichage(hs,erreurs,"Évolution de l'erreur dans le calcul de la dérivéeen fonction du pas","h","erreur relative (%)",False)foriinrange(len(hs)):print("Pas : {:.5e}, Nombre calculé : {:.9f}, erreur relative (%) : {:.9e}".format(hs[i],dfs[i],erreurs[i]))main()
Soit une fonction $f$ possédant $n+1$ dérivées continues sur l’intervalle $[a,b]$ et soient $x$ et $x_0$ deux points de cet intervalle. On peut alors écrire :
$$f(x) = p_n(x) + R_{n+1}(x)$$
avec
$$
\begin{array}{rcl}
p_n(x) & = & f(x_0) + \dfrac{x - x_0}{1!} f^{\prime}(x_0) + \ldots + \dfrac{(x - x_0)^n}{n!} f^{(n)}(x_0), \cr
R_{n+1}(x) & = & \dfrac{1}{n!} \int_{x_0}^{x} (x-t)^n f^{(n+1)}(t) dt \cr
& = & \dfrac{(x - x_0)^{(n+1)}}{(n+1)!} f^{(n+1)}(\epsilon), \hspace{0.5cm} x_0 \leqslant \epsilon \leqslant x. \cr
\end{array}
$$
En limitant le développement de Taylor à l’ordre 2, on peut écrire :
$$f(x+h) = f(x) + f^{\prime}(x)h + f^{\prime\prime}(\epsilon)\dfrac{h^2}{2} \hspace{0.5cm} |x - \epsilon| < |h|$$
Les erreurs de troncature sont liées à la précision de l’algorithme utilisé. Si, par exemple, une fonction est approchée par son développement de Taylor, l’erreur de troncature sera obtenue par une évaluation du reste du développement. Son contrôle sera obtenu par une majoration de ce reste.
Erreur de méthode
Les erreurs de méthode interviennent lorsqu’une expression est mal équilibrée et/ou mélange des valeurs très différentes. C’est un problème qui engendre des erreurs d’arrondi. Dans la plupart des cas, une erreur de méthode impose de modifier l’algorithme.
Exemple 1. Détermination de la position du milieu d’un segment.
Le calcul effectué spontanément : $$m = \dfrac{a+b}{2}$$ peut conduire à des erreurs si la somme $a+b$ dépasse les capacités de représentation des nombres du calculateur, alors que la valeur de $m$ est représentable. Le calcul à préférer est : $$m = a + \dfrac{b-a}{2}$$
1
2
3
4
5
6
7
8
>>>a=1e308>>>b=1.7e308>>>(a+b)/2# a + b > 1.7976931348623157e+308 (valeur maximale) inf# on obtient une erreur>>>a+(b-a)/2# (b -a ) < 1.7976931348623157e+308 (valeur maximale) 1.35e+308# on obtient bien la moyenne
Exemple 2. Détermination des racines du polynôme du second degré.
Par exemple, la recherche des racines du polynôme $10^{-8} x^2 - 0,8 x + 10^{-8}$, dont les deux racines sont $x_1 \cong 0.8\cdot 10^{8}$ et $x_2 \cong 1.25\cdot 10^{-8}$, nécessite le calcul du discriminant $\Delta = 0,64 - 4\cdot 10^{-16}$ dont le résultat n’est pas toujours correct car $4\cdot 10^{-16}$ peut être considéré comme négligeable par certains calculateurs, comparé à 0,64.
Pour obtenir les racines, il est donc nécessaire de modifier l’algorithme en considérant, par exemple, la somme et le produit des racines.
\end{description}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
>>>a=1e-8>>>b=-0.8>>>c=1E-8>>>delta=b**2-4*a*c>>>frommathimportsqrt>>>x1=(-b+sqrt(delta))/(2*a)# Première racine >>>x179999999.99999999>>>x2=(-b-sqrt(delta))/(2*a)# Deuxième racine >>>x21.1102230246251565e-08>>>somme_racines=somme_racines=-b/a# Ici on utilise les propriétés des racines, pas les calculs précédents>>>produit_racines=c/a# Idem>>>x1+x2,somme_racines# Semble identiques (80000000.0,80000000.0)>>>x1*x2,produit_racines(0.8881784197001251,1.0)# Problème !
Erreur totale
L’erreur totale, dans l’exemple étudié dans ce chapitre, est définie par :
$$
E = E_a + E_t
$$
Dans le cadre de cette étude3, l’erreur totale a donc pour expression :
$$ E = 2 \delta |\dfrac{f(x)}{h}| + \dfrac{|h|}{2} |f^{\prime\prime}(\epsilon)|$$
Rechercher la valeur approchée la plus précise possible de l’expression (1.1.) nécessite donc de minimiser $E$. $$E = 2 \delta |\dfrac{f(x)}{h}| + \dfrac{|h|}{2} |f^{\prime\prime}(\epsilon)|$$ est minimale pour $$|h| = 2 \sqrt{\delta |\dfrac{f(x)}{f^{\prime\prime}(\epsilon)}|}$$
Démonstration :
Soit $g(s) = \dfrac{a}{s} + bs$ avec $a = 2 \delta |f(x)|$ et $b = \dfrac{1}{2} |f^{\prime\prime}(\epsilon)|$. Puisque $a$ et $b$ sont positifs, $g^{\prime}(s) = 0 \Longleftrightarrow -\dfrac{a}{s^2} +b =0 \Longleftrightarrow s = \sqrt{\dfrac{a}{b}}$
Retour à l’exemple : $f(x) = x^2$ et $\delta = 10^{-3}$
$f^{\prime\prime}(x) = 2, h = \sqrt{2 \delta f(x)} = \sqrt{2\delta} |x|$. Pour $x=7$, $h$ a pour valeur : $h = 10\sqrt{\delta} \approx 0,3$.
Convergence et stabilité
On dit qu’un algorithme est stable si le résultat des opérations est peu sensible à l’accumulation des erreurs.
L’étude de stabilité se fait par des méthodes expérimentales et parfois théoriques.
Dérivation d’une fonction analytique, les différents schémas
L’expression :
$$f^{\prime}(x) = \dfrac{f(x+h) - f(x)}{h}$$
est appelée différence finie à droite.
On peut aussi choisir une différence finie à gauche :
$$f^{\prime}(x) = \dfrac{f(x) - f(x-h)}{h}$$
et une différence finie centrée autour de $x$ :
$$f^{\prime}(x) = \dfrac{f(x+h) - f(x-h)}{2h}$$
Dans toutes ces expressions l’erreur d’arrondi reste la même. Par contre, l’erreur de troncature est modifiée dans la dernière expression.
Démonstration :
À partir de la formule de Taylor,
$$
\begin{array}{rcl}
f(x+h) & = & f(x) + hf^{\prime}(x) + \frac{h^2}{2}f^{\prime\prime}(x) + \frac{h^3}{6}f^{\prime\prime\prime}(\epsilon) \hspace{1cm} \text{où } \epsilon \in [x,x+h] \cr
f(x-h) & = & f(x) - hf^{\prime}(x) + \frac{h^2}{2}f^{\prime\prime}(x) - \frac{h^3}{6}f^{\prime\prime\prime}(\eta) \hspace{1cm} \text{où } \eta \in [x-h,x] \cr
\end{array}
$$
Donc $$f^{\prime}(x) - \dfrac{f(x+h) - f(x-h)}{2h} = -\dfrac{h^2}{12}[f^{\prime\prime\prime}(\epsilon) + f^{\prime\prime\prime}(\eta)]$$
Finalement, l’erreur de troncature est diminuée lorsqu’on utilise la différence finie centrée autour de $x$.
"""
Différences finies centrées, ordre 02 en h.
Mise en évidence des erreurs de troncature et d'arrondis.
"""frommathimportsin,cosdefdf(x):returncos(x)deferr_centree(x,f,df,h):returnabs(df(x)-(f(x+h/2)-f(x-h/2))/h)deferr_progressive(x,f,df,h):returnabs(df(x)-(f(x+h)-f(x))/h)deferr_retrograde(x,f,df,h):returnabs(df(x)-(f(x)-f(x-h))/h)defmain():f=lambdax:sin(x)# Fonction étudiéedf=lambdax:cos(x)# Fonction dérivée de la fonction étudiéex0=1# Point où le nombre dérivé est calculéhs=[1/(10**i)foriinrange(1,15)]# Liste des pas hprint("h \t\t Diff. finies Prog. Diff. finies Retro. Diff. finies centrées")forhinhs:print("{:.1e}\t{:.4e}\t{:.4e}\t{:.4e}".format(h,err_progressive(x0,f,df,h),err_retrograde(x0,f,df,h),err_centree(x0,f,df,h)))main()
On peut être amené à calculer la dérivée d’une fonction définie par une table de valeurs (résultat de mesures expérimentales, par exemple). Une bonne compréhension des résultats des formules précédentes (en particulier celle de la différence finie centrée autour de $x$) nécessite de prendre en compte :
que le pas $h$ est maintenant imposé (et pas forcément constant) ;
que les inévitables erreurs de mesure vont être « amplifiées » par la division par $h$ (qui est « petit »).
Des méthodes existent pour l’importance du second point, ci-dessus, mais elles ne seront pas abordées cette année.
Dérivée seconde d’une fonction
Pour déterminer une expression de la dérivée seconde, on utilise à nouveau la formule de Taylor, en « poussant » à l’ordre 4 le développement :
$$
\begin{array}{rcl}
f(x+h) & = & f(x) + hf’(x) + \frac{h^2}{2}f^{\prime\prime}(x) + \frac{h^3}{6}f^{\prime\prime\prime}(x) + \frac{h^4}{24}f^{\prime\prime\prime\prime}(\epsilon) \text{ où } \epsilon \in [x,x+h]\cr
f(x-h) & = & f(x) - hf’(x) + \frac{h^2}{2}f^{\prime\prime}(x) - \frac{h^3}{6}f^{\prime\prime\prime}(x) + \frac{h^4}{24}f^{\prime\prime\prime\prime}(\eta) \text{ où } \eta \in [x-h,x] \cr
\end{array}
$$