Discussion :

[crippella]

bonjour,
j'ai un mini projet en matlab sur la methode de newton raphson sur les équations non liniare j'ai fait un essai mais j'ai pas obtenu la bonne solution
j'ai pris comme équation:f(x)=2x+exp(x/x+1)
voila le programme :
code:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
x=-10:0.1:10;
n=50;
y=2*x+exp(x/(x+1));
   z=diff(y)./diff(x);
y(1)=1;
for i=1:n;
   y(i+1)=y(i)-y(i)/z(i);
end
plot(x,y)

aidez moi s'il vous plait
merci d'avance

[kmaniche]

Attention à l'écriture du code :

* => .*
/ => ./

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
 
x=-10:0.1:10;
n=50;
y=2.*x+exp(x./(x+1));
 
z=diff(y)./diff(x);
 
y(1)=1;
 
 for i=1:n;
    y(i+1)=y(i)-y(i)./z(i);
 end
plot(x,y)

[Jean Dumoncel]

@ crippella : un conseil, étudie bien la méthode de newton-raphson avant de passez au code. Tu en trouveras une explication claire et détaillée sur wikipedia. Déroules la méthode au brouillon avec ton exemple.

A chaque itération, tu dois calculer :

C'est à l'intérieur de ta boucle que tu dois calculer f(xn) et f'(xn).

[tubaas]

@Magelan : je n epense pas.
il faut d'abord calculer f' avant de la prendre en x(n)
par contre il ne me semble pas que

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

    y(i+1)=y(i)-y(i)./z(i);

colle avec la formule puisque mélange de xn et f(xn)

edit: ayant zyeuté vite fait wikipédia je pense qu'il faut écrire

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
 
vzero=zeros(1,n+1);
vzero(1)=1
for i =1:n
vzero(i+1)=vzero(i)-y(i)./z(i);
end

[FR119492]

Salut!
Tu as besoin de l'expression analytique de f'(x). La fonction diff ne convient pas.
Jean-Marc Blanc

[Jean Dumoncel]

@Jean-Marc Blanc : oui, c'est juste. Par contre on peut prendre une approximation de la dérivée en un point, mais là on s'éloigne un peu de la méthode de Newton-Raphson.

@tubaas : je t'assure qu'il ne faut pas faire ce calcul avant. Regardes l'exemple de la page de wikipedia pour t'en persuader. Si f est la fonction étudiée, dans ton code, tu utilises f(i) et f'(i) à chaque itération, ce qui n'a pas de sens. Il faut calculer f(x(n)) et f'(x(n)), et donc on a besoin de x(n) pour les calculer!

[tubaas]

je ne pense pas qu'on aie besoin de l'expression analytique.
Matlab sert à faire du calcul numérique donc il y a forcément des approximations et des erreurs à minimaliser. le but est de trouver les bons paramètres pour concilier une convergence rapide vers la solution avec une erreur acceptable

@magelan: c'est exact, autant pour moi.
cela dit, on a besoin de xn mais ce n'est pas une erreur de calculer f' une fois pour toute et ensuite on la prend en xn, j'insiste (c'est vrai qu'en fait ça ne change rien on pourrait calculer dans la boucle la dérivée mais mon insistance est due au fait que dans le premier code, y était modifié dans la boucle alors ce serait une catastrophe de calculer la dérivée sur un y modifié).

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
 
dx=0.01;
x=-10:dx:10;
n=50;
y=2*x+exp(x./(x+1));
 
z=diff(y)./diff(x);
z(length(x))=z(length(x)-1);%pour avoir une expression de la dérivée en tous les points
 
vzero=zeros(1,n+1);
vzero(1)=1;
for i =1:n
pind=find(abs(vzero(i)-x)<=dx/2,1,'first');%on trouve l'indice correspondant à la meilleure approximation possible de xn
vzero(i+1)=vzero(i)-y(pind)/z(pind);
%ou sinon
%vzero(i+1)=vzero(i)-y(pind)/((y(pind+1)-y(pind-1))/(x(pind+1)-(x(pind-1)));
%tu as en fait plusieurs choix pour approximer ta dérivée en ce point si tu le fais dans la boucle :
%(y(pind+1)-y(pind))/(x(pind+1)-(x(pind))
%(y(pind)-y(pind-1))/(x(pind)-(x(pind-1))
%(y(pind+1)-y(pind-1))/(x(pind+1)-(x(pind-1))
%attention cependant à gérer les limites de ton domaine (pind=1 ou pind=length(x))
%de manière particulière ce que je ne fais pas là
%c'est vrai que tu as ainsi l'avantage de pouvoir comparer différentes méthodes de calcul de dérivée
end

ceci semble converger vers -0.3154 en améliorant dx (ie : diminuant !, icic avec dx=0.0001).
dans cet algo tu dois jouer sur dx et n et sur vzero(1) !.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
pind=find(abs(-0.3154-x)<=dx/2,1,'first');
>> y(pind)
 
ans =
 
    3.749694754406896e-005
 
>> y(pind+1)
 
ans =
 
    3.720912274263188e-004
 
>> y(pind-1)
 
ans =
 
   -2.971079354017192e-004

je t'ai mis un zoom de ta fonction y qui me semblait intéressant après l'avoir ploté

[phryte]

Bonjour.
<un exemple de l'algorithme :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
clear
syms x
fonc = 2*x+exp(x/(x+1));
foncprime=diff(fonc);
erreur = 1e-6;
n = 0;x = -10;
hold on
while abs(eval(fonc)) > erreur
 n = n+1; 
if eval(foncprime) ~= 0
 xn = x - eval(fonc)/eval(foncprime);%Formule
else  disp('Dérivée nulle'); end
plot([n-1 n],[x xn],'-k*','MarkerSize',6,'LineWidth',1)
x = xn; end;grid

[tubaas]

un exemple de l'algorithme :

attention ! : si on dispose de la symbolic math toolbox (j'ai voulu tester )

[Jean Dumoncel]

@ tubaas : décidément tu es tétu...

Quelque inconvénient de ton algo :
tu limites la recherche entre -10 ET 10
Tu trouves une solution correcte pour dx=0.0001, donc dans x tu places 200000 valeurs et si n=50, au maximum tu n'en utilises que 50...
et plus tu augmente la taille de ton x, plus la recherche prendra du temps.

Fais-le en calculant uniquement la dérivée à chaque itération, (ça prends 7 lignes de code), tu verrras les différences...

[tubaas]

je suis d'accord avec vos remarques et que le code est optimisable.
cependant je présente juste des solutions, je ne fais pas tout le travail, dont le plus important, l''optimisation' (entendu ici au sens de calcul numérique : choix des paramètres et de méthodes de calcul de dérivée qui dépendent du problème).
ensuite à moins d'utiliser tic toc je ne verrai pas la différence car j'ai la chance que ce soit quasi instantané (je sais c'est pas un argument ça !)
de plus ce qui est long c'est la boucle sur n
je peux la réduire car la convergence est beaucoup plus rapide.
je suis obligé de remplir x au max selon la précision sur la valeur que je recherche.
si elle est de 10^(-3) je DOIS faire un pas de 0.001 car la limitation vient de là (enfin il me semble, étant donné que mes vzéros devront finalement être approximé par une valeur du vecteur x)
qu'importe le nombre de valeurs que j' utilise : il faut qu'il y ait une bonne approximation de la valeur que je recherche, et si pour cela je prends 100000000000000000000000 de valeurs pour n'en utiliser que 50, c'est que j'y suis contraint, sinon il n'y aura pas la valeur que je recherche.

edit: image de la convergence de vzero, dx=0.0001, n=10, vzero(1)=-1

edit2: je suis d'accord qu'il est mieux de dériver dans la boucle car mon vecteur z et non x sera trop consistant, peut-être est-ce ce que vous vouliez-dire, cependant ce n'est pas faux de dériver avant car si l'on applique le même principe pour le calcul de la dérivée (y(pind+1)-y(pind))/(x(pind+1)-x(pind)), le résultat sera identique je pense

edit3: et oui je suis têtu désolé

edit4 : pour n grand, il vaut mieux (pour le temps d'exécution mesuré avec tic toc) calculer f' (z) au préalable. le point d'équilibre se trouve aux alentours de n=280. mais comme on a vu que l'on convergeait très vite...

EDIT5 : après explication par mp, je dois toutes mes excuses à Magelan dont la méthode qui consiste à construire uniquement les termes de la suite convergente est bien plus efficace et propre. je n'avais tout simplement pas compris où il voulait en venir.

[FR119492]

Bonjour à tous!
C'est vrai qu'on peut utiliser une approximation de la dérivée, mais ce n'est plus la méthode de Newton-Raphson, et c'est cette méthode que l'énoncé de son problème impose à Crippella.
Jean-Marc Blanc

[membreComplexe12]

C'est vrai qu'on peut utiliser une approximation de la dérivée, mais ce n'est plus la méthode de Newton-Raphson,

tout à fait d'accord, au faite matlab ne fait pas de calcul formel? il ne pourrais pas donner l'expression analytique? (je suis debutant je ne connais pas toutes les onctions de matlab)

[phryte]

Bonjour.

matlab ne fait pas de calcul formel?

Regarde un peu plus haut.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
clear
syms x
fonc = 2*x+exp(x/(x+1));
foncprime=diff(fonc);

[membreComplexe12]

donc la reponse est non? c'est pour cela que vous utilisez diff qui est une approximation de la derivée en un point.

[phryte]

La réponse est oui, diff effectue le calcul analytique :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
fonc =
 
2*x+exp(x/(x+1))
 
 
foncprime =
 
2+(1/(x+1)-x/(x+1)^2)*exp(x/(x+1))

[Jean Dumoncel]

Comme cela a déjà été précisé dans ce topic, MATLAB peut faire du calcul symbolique si l'on possède la toolbox math symbolique.


Comme beaucoup de fonction, diff est une fonction surchargée :
il existe entre autre une version de référence qui fait le calcul numérique de la dérivée
(Differences and approximate derivatives)

et une autre dans la toolbox math symbolique (Differentiate symbolic expression).

[membreComplexe12]

La réponse est oui, diff effectue le calcul analytique :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
fonc =
 
2*x+exp(x/(x+1))
 
 
foncprime =
 
2+(1/(x+1)-x/(x+1)^2)*exp(x/(x+1))

a ok, en regardant rapidement l'aide (trop rapidement apparemment) j'avais compris le contraire!

merci de ta reponse
donc si on peut avoir l'expression de la derivée analytiquement on effectue bien une methode de Newton Raphson et non une methode de quasi newton..


edit=> ok merci aussi magellan pour tes precisions....


a+