newton.m

来自「共轭梯度法」· M 代码 · 共 29 行

function   s=newton(F,x,d)
x=sym(x);
d=sym(d);
a=sym('a');
E=0.0000001;
p=sym(findsym(F));
n=length(p);
q=x+a*d;
j=length(q);
if n==j
for i=1:1:n
    F=subs(F,q(i),p(i));
end
F=collect(F);
end
F1=diff(F,[a]);
F2=diff(F1,[a]) ; 
a=0.0000001;
f1=subs(F1);
f2=subs(F2);
while abs(f1)>=E
    a=a-f1/f2;
    f1=subs(F1);
    f2=subs(F2);
end    
s=a;