232
6 R´ esolution des ´ equations et des syst` emes non lin´ eaires
Programme 49 - newthorn : M´ ethode de Newton-Horner avec raffinement
function [xn,iter,root,itrefin]=newthorn(A,n,tol,x0,nmax,iref)
%NEWTHORN M´ ethode de Newton-Horner avec raffinement.
% [XN,ITER,ROOT,ITREFIN]=NEWTHORN(A,N,TOL,X0,NMAX,IREF) tente
% de calculer toutes les racines d’un polynˆ ome de degr´ e N et de
% coefficients A(1),...,A(N). TOL est la tol´ erance de la m´ ethode.
% X0 est la donn´ ee initiale. NMAX est le nombre maximum d’it´ erations.
% Si IREF vaut 1, la proc´ edure de raffinement est activ´ ee.
apoly=A;
for i=1:n, it=1; xn(it,i)=x0+sqrt(-1)*x0; err=tol+1; Ndeg=n-i+1;
if Ndeg == 1
it=it+1; xn(it,i)=-A(2)/A(1);
else
while ittol
[px,B]=horner(A,Ndeg,xn(it,i)); [pdx,C]=horner(B,Ndeg-1,xn(it,i));
it=it+1;
if pdx ˜=0
xn(it,i)=xn(it-1,i)-px/pdx;
err=max(abs(xn(it,i)-xn(it-1,i)),abs(px));
else
fprintf(’ Arrˆ et dˆ u `
a une annulation de p’’ ’);
err=0; xn(it,i)=xn(it-1,i);
end
end
end
A=B;
if iref==1
alfa=xn(it,i); itr=1; err=tol+1;
while err>tol*1e-3 & itr
[px,B]=horner(apoly,n,alfa); [pdx,C]=horner(B,n-1,alfa);
itr=itr+1;
if pdx˜=0
alfa2=alfa-px/pdx;
err=max(abs(alfa2-alfa),abs(px));
alfa=alfa2;
else
fprintf(’ Arrˆ et dˆ u `
a une annulation de p’’ ’);
err=0;
end
end
itrefin(i)=itr-1; xn(it,i)=alfa;
end
iter(i)=it-1; root(i)=xn(it,i); x0=root(i);
end
return
6 R´ esolution des ´ equations et des syst` emes non lin´ eaires
Programme 49 - newthorn : M´ ethode de Newton-Horner avec raffinement
function [xn,iter,root,itrefin]=newthorn(A,n,tol,x0,nmax,iref)
%NEWTHORN M´ ethode de Newton-Horner avec raffinement.
% [XN,ITER,ROOT,ITREFIN]=NEWTHORN(A,N,TOL,X0,NMAX,IREF) tente
% de calculer toutes les racines d’un polynˆ ome de degr´ e N et de
% coefficients A(1),...,A(N). TOL est la tol´ erance de la m´ ethode.
% X0 est la donn´ ee initiale. NMAX est le nombre maximum d’it´ erations.
% Si IREF vaut 1, la proc´ edure de raffinement est activ´ ee.
apoly=A;
for i=1:n, it=1; xn(it,i)=x0+sqrt(-1)*x0; err=tol+1; Ndeg=n-i+1;
if Ndeg == 1
it=it+1; xn(it,i)=-A(2)/A(1);
else
while it
[px,B]=horner(A,Ndeg,xn(it,i)); [pdx,C]=horner(B,Ndeg-1,xn(it,i));
it=it+1;
if pdx ˜=0
xn(it,i)=xn(it-1,i)-px/pdx;
err=max(abs(xn(it,i)-xn(it-1,i)),abs(px));
else
fprintf(’ Arrˆ et dˆ u `
a une annulation de p’’ ’);
err=0; xn(it,i)=xn(it-1,i);
end
end
end
A=B;
if iref==1
alfa=xn(it,i); itr=1; err=tol+1;
while err>tol*1e-3 & itr
itr=itr+1;
if pdx˜=0
alfa2=alfa-px/pdx;
err=max(abs(alfa2-alfa),abs(px));
alfa=alfa2;
else
fprintf(’ Arrˆ et dˆ u `
a une annulation de p’’ ’);
err=0;
end
end
itrefin(i)=itr-1; xn(it,i)=alfa;
end
iter(i)=it-1; root(i)=xn(it,i); x0=root(i);
end
return
