-- Simulation numérique en mécanique des fluides --
TP n°9
Méthode de Newton
éléments de correction
simon.marie@lecnam.net
import numpy as np
import matplotlib.pyplot as plt
1 - Méthode de Newton scalaire
La méthode de Newton permet de trouver numériquement la solution de l'équation $f(x)=0$ en partant d'un point de départ $x_0$ et en itérant la suite:
$$ x_{n+1}=x_n-\dfrac{f(x_n)}{f'(x_n)} $$
jusqu'à atteindre un critère de convergence.
Ecrire une fonction permettant de calculer la solution de $f(x)=0$ en connaissant sa dérivée que l'on notera $df$ ainsi que le point de départ que l'on notera $x_0$:
def newton(f,df,x0):
...
return x0
On souhaite utiliser cette méthode pour calculer la solution de l'équation x^2-a=0 Quels doivent être les fonctions $f(x)$ et $f'(x)$ pour que la méthode de Newton donne la solution de cette équation.
def f(x,a=2):
return ...
def df(x,a=2):
return ...
Donner une valeur approchée des nombre suivant: $\sqrt{2}$$\sqrt{3}$ $\sqrt{5}$$\sqrt{7}$
Comarer la convergence de la méthode en fonction de la valeur de $a$
Commentaires
2 - Méthode de Newton vectorielle¶
La méthode de Newton permet de trouver numériquement la solution de l'équation $F(x,y)=0$ en partant d'un point de départ $X_0=(x_0,y_0)$ et en itérant la suite:
$$ \mathbf{J(X_n)}d\mathbf{X_{n+1}}=-\mathbf{F(X_n)} $$
avec:
$$ d\mathbf{X_{n+1}}=\mathbf{X_{n+1}-X_n} $$
jusqu'à atteindre un critère de convergence.
1 - Système non linéaire
Considérons le système suivant:
$$ \left\{\begin{array}{lcr} 9x^2+36y^2+4z^2&=&36\\ x^2-2y^2-20z&=&0\\ x^2-y^2+z^2&=&0 \end{array}\right. $$
Quelle est la Jacobienne du système ?
def F(W):
F1=...
F2=...
F3=...
return np.array([F1,F2,F3])
def J(W):
J=...
return J
def NEWTON(F,dF,W0):
norm=10.
n=0
while (...):
n+=1
...
...
...
print(n,norm)
return W0
En partant du triplet $(1,1,-1)$ donner une valeur approchée de la solution.
W0=...
W=NEWTON(F,J,W0)
print(W)
En évaluant la fonction sur le point $W$ trouvé, assurez-vous que cette valeur (appellé résidu) est très proche de 0.
F(W)
Commentaires....
2 - EDP non linéaire
Considérons maintenant l'EDP non linéaire suivante:
$$ u\dfrac{\partial u}{\partial x}=\nu\dfrac{\partial^2 u}{\partial x^2} $$
Avec les conditions $u(0)=3$ et $u(N)=-1$
Ecrire cette EDP sous forme conservative
Ecrire l'équation sous forme discrète en utilisant des schémas centrés d'ordre 2
Ou souhaite utiliser la méthode de Newton pour résoudre cette EDP non linéaire.
On a ici $F(u)=\dfrac{1}{2}\dfrac{\partial u^2}{\partial x}-\nu\dfrac{\partial^2 u}{\partial x^2}$
avec:
$\dfrac{d F}{d u}=\dfrac{\partial .}{\partial x}-\nu\dfrac{\partial^2 .}{\partial x^2}$
Compléter la cellule suivante pour pouvoir utiliser la méthode de Newton
def F(W,dx,nu):
nx=W.shape[0]
ff=np.zeros(nx)
ff[1:-1]=...
ff[0]=...
ff[-1]=...
return ff
def J(W,dx,nu):
nx=W.shape[0]
J1=...
J2=...
J0=J1+J2
J0[0,0]=...
J0[-1,-1]=...
return J0
def NEWTON(F,dF,W0,dx=1.,nu=0.1,nmax=10000):
norm=10.
n=0
while (norm>1e-3 and n<nmax):
n+=1
DX=...
W0=W0+DX
norm=np.sum(DX**2)**0.5
print(n,norm)
return W0
Calculer la solution de l'EDP ci-dessus en utilisant la méthode de Newton.
W0=...
W=NEWTON(F,J,W0)
Tracer le résidu F(W) (il doit être le plus uniformément nul).
plt.plot(F(W,dx,nu))
Représenter la solution trouvée u(x)
plt.plot(W)
Commentaires...
Conclusion¶
Présentez une synthèse de vos résultats et des étapes de calcul de votre TP.
from IPython.core.display import HTML
style=open('notebooks.css', "r").read()
HTML(style)