function [xn,iter]=newtonm(x0,f,x)
format long;
epsilon=1e-10;
maxval=10000;
J=jacobian(f,x);
d=1;
iter=0;
while (d>epsilon)
JJ=subs(J,x,x0);
F=-subs(f,x,x0);
JJ=double(JJ);
F=double(F);
if abs(det(JJ))<epsilon
error('Newton is singular,choose a new initial value')
disp(['iterations = ' num2str(iter)]);
break;
end
if abs(F)>maxval
error('Solution diverge');
disp(['iterations = ' num2str(iter)]);
break;
end
xn=x0+inv(JJ)*F;
d=norm(xn-x0,inf);
x0=xn;
iter=iter+1;
end
end
