function x = newton_solver(f,m,x0,dfdx) %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % % function x = newton_solver(f,m,x0,dfdx) % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % % % Solves the equation f(x) = m using x0 as % % start guess where f is a string with the name% % of the file in which f is defined. % % The jacobian can be specified % % in the file 'dfdx'. If omitted, a numerical % % approximation is used % % % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% tol=0.01; x=x0; fx=feval(f,x0); err=1; errmax=100; if nargin == 4 % Use exact jacobian while (err>tol) fx = feval(f,x); J = feval(dfdx,x); dx = -J\(fx-m); x = x+dx; fx = feval(f,x); err=norm(fx-m,'inf'); if(err>errmax) error('Divergent Newton solver (file: newton_solver.m)') end end else % Use approximate jacobian while (err>tol) fx = feval(f,x); J = appr_dfdx(f,x); dx = -J\(fx-m); while (norm(feval(f,x+dx)-m)>norm(fx-m)) % added 000211 dx=dx./2 end % end added x = x+dx; fx = feval(f,x); err=norm(fx-m,'inf'); if(err>errmax) error('Divergent Newton solver (file: newton_solver.m)') end end end