Méthode de Newton

La méthode de Newton consiste à chercher une racine de fonction. La seule condition est que la fonction soit dérivable sur tout l'intervalle considéré.

Considérons la fonction $f: x \mapsto 0,5\,x^4 - 2\,x^3 + x^2 + x +1$

Commençons par tracer la fonction. On peut bien sûr utiliser un grapheur comme desmos. Mais on peut aussi tracer la courbe avec Python et matplotlib.pyplot.

import matplotlib.pyplot as plt
import numpy as np

def f(x):
    return 0.5*x**4 - 2*x**3 + x**2 + x + 1

# tracé du graphique de f sur [-1;3.5]
xs = np.linspace(-1, 3.5 , num = 100)
ys = f(xs)  # création du tableau pour y

plt.plot(xs,ys)

Quelques commentaires :

• numpy permet de gérer des tableaux spéciaux.
• Ligne 8, xs est justement un tableau numpy qui contient 100 valeurs allant de -1 à 3.5 régulièrement espacées.
• Ligne 9, ys est encore un tableau numpy – la ligne de commande toute simple ys = f(xs) fait que ys est un tableau constitué des valeurs f(x) pour x contenu dans xs. C'est toute la puissance de numpy que de permettre ce genre de raccourci.
• Ligne 11, avec matplotlib, tracé d'une courbe passant par l'ensemble des cordonnées (x;y) avec les x pris dans xs et les y pris dans ys. La courbe est donc une suite de segments mais comme ils sont courts, on a l'impression d'une courbe arrondie – même système que sur une calculatrice !

Il s'agit en gros de suivre le sens de la pente. La dérivée d'une fonction me dit comment la fonction varie. Par exemple, supposons que $x = 1$, alors $y = f(1) = 1,5$. Si je veux diminuer $y$, je dois diminuer $f(x)$ et pour savoir comment diminuer $f(x)$ j'ai besoin de $f'(x)$.

$f'(x) = 2\,x^3 - 6\,x^2 + 2\,x + 1 \Rightarrow f'(1) = -1 < 0$

On voit donc que $f$ décroît en $x = 1$ et donc pour diminuer $f(x)$ il faut augmenter $x$.

Faut-il l'augmenter beaucoup ? On fait comme si $f$ allait diminuer à une vitesse constante. Puisque $y = f(x) = 1,5$ et $f'(x) = -1$, il faut que $x$ augmente de $-\frac{f(x)}{f'(x)} = 1,5$. On va donc aller voir en $x = 2,5$.

Sur ce graphique, $x_0 = 1$ est le point de départ. En utilisant la dérivée en $x_0$, on peut tracer la tangente en vert et calculer la valeur $x_1 = 2,5$.

$x_1$ n'est pas une racine. Alors on recommence.

ENTRÉES : réel x0, fonction f, réel positif seuil_erreur
SORTIE : estimation d'une racine x de la fonction f
DÉBUT
    soit x = x0
    TANT QUE l'écart entre f(x) et 0 est supérieur à seuil_erreur RÉPÉTER
        x prend la valeur x - f(x) / f'(x)
    FIN
    RENVOYER x
FIN

L'algorithme utilise $f'(x)$. Dans notre exemple, nous connaissions $f'(x)$ car nous pouvions faire nous même le calcul suivant :

$$f(x) = 0,5\,x^4 - 2\,x^3 + x^2 + x +1 \Rightarrow f'(x) = 2\,x^3 - 6\,x^2 + 2\,x + 1$$

Mais alors, faudra-t-il que nous fournissions au programme l'expression de $f$ et celle de $f'$ ? Ne serait-il pas préférable que la machine détermine automatiquement $f'(x)$ ? Mais c'est sans doute beaucoup plus difficile… – C'est possible néanmoins !

Pour contourner cette difficulté, je vous propose de ne pas calculer $f'(x)$ avec une formule exacte mais de nous contenter d'une approximation – Après tout, nous cherchons une racine approximative, alors quitte à approximer…

Rappelez vous qu'en théorie $\displaystyle f'(x) = \lim_{\varepsilon \to 0} \frac{f(x + \varepsilon) - f(x)}{\varepsilon}$. C'est une idée assez abstraite puisque $\varepsilon$ doit devenir aussi proche que l'on veut de $0$ mais sans atteindre $0$.

On peut se contenter d'une approximation en prenant $\varepsilon$ petit. Par exemple :

EPSILON = 1e-10

# exemple de définition de f
def f(x):
    return 0.5*x**4 - 2*x**3 + x**2 + x + 1

# définition de la dérivée approximative
def der_f(fonction, x):
    """
    x: réel, antécédent de la fonction
    fonction: fonction de x
    renvoie une approximation de fonction'(x)
    """
    return (fonction(x + EPSILON) - fonction(x)) / EPSILON

Reprenons ce que nous avons déjà :

EPSILON = 1e-10

# exemple de définition de f
def f(x):
    return 0.5*x**4 - 2*x**3 + x**2 + x + 1

# définition de la dérivée approximative
def der_f(fonction, x):
    """
    x: réel, antécédent de la fonction
    fonction: fonction de x
    renvoie une approximation de fonction'(x)
    """
    return (fonction(x + EPSILON) - fonction(x)) / EPSILON
    
# Votre travail
def newton(x0, fonction, seuil_erreur):
    """
    x0: réel, valeur initiale de la recherche
    fonction: fonction de x
    seuil_erreur: réel positif,
      si écart entre f(x) et 0 < seuil_erreur, on considère x comme une racine
    renvoie une racine approximative de la fonction
    """
    
    # à vous

# Exemple d'utilisation
r = newton(1, f, 1e-5)
print("Une racine approximative de f est {}".format(r))