i noticed that there is an error in the sign of third term of x update equation (for higher order Newton iteration): the line x := (x -(f(x)/D(f)(x)) -((((D@@2)(f))*((f(x))^2))/(2*(D(f)(x))^3))) should be x := x -(f(x)/D(f)(x)) + ((((D@@2)(f))*((f(x))^2))/(2*(D(f)(x))^3)); then following modified code finds the solution NEWTON_RAPH:= proc(f,z,M) local x,A; x := z: for A to M do: x := x -(f(x)/D(f)(x)) + ((((D@@2)(f))*((f(x))^2))/(2*(D(f)(x))^3)); end do: evalf(x): proc: f:= x-> x^2-5: NEWTON_RAPH(f,1,10); evalf(%);