% script file Newton3.m   
   %  User first creates inline functions or mfiles for f(x,y,z), g(x,y,z) and
   % h(x,y,z).  Script will ask user to enter the starting value 
   % as a column vector.  Does 5 iterations. 

     p = input('enter the starting value as a col. vector   ')
     
     % symbolic expressions for f g and h.

     syms x y z;
     ff = sym(f(x,y,z));
     gg = sym(g(x,y,z));
     hh = sym(h(x,y,z));

     % symbolic differentiation 
     ffx = diff(ff,x);
     ffy = diff(ff,y);
     ffz = diff(ff,z);

     ggx = diff(gg,x);
     ggy = diff(gg,y);
     ggz = diff(gg,z);

     hhx = diff(hh,x);
     hhy = diff(hh,y);
     hhz = diff(hh,z);

     fx = inline(vectorize(ffx),'x', 'y', 'z')
     fy = inline(vectorize(ffy), 'x', 'y', 'z')
     fz = inline(vectorize(ffz), 'x', 'y', 'z')
     gx = inline(vectorize(ggx), 'x', 'y', 'z')
     gy = inline(vectorize(ggy), 'x', 'y', 'z')
     gz = inline(vectorize(ggz), 'x', 'y', 'z')
     hx = inline(vectorize(hhx), 'x', 'y', 'z')
     hy = inline(vectorize(hhy), 'x', 'y', 'z')
     hz = inline(vectorize(hhz), 'x', 'y', 'z')

     disp(' ')
     disp('  iterate         p             error (Newton step)  ')
     format compact
     format long
     for n = 1:5
        J = [fx(p(1),p(2), p(3)) fy(p(1), p(2), p(3)) fz(p(1),p(2), p(3)) ; ...
             gx(p(1), p(2), p(3)) gy(p(1), p(2), p(3)) gz(p(1), p(2), p(3)); ...
             hx(p(1), p(2),p(3)) hy(p(1), p(2), p(3)) hz(p(1),p(2), p(3)) ] ;
        F = [f(p(1),p(2),p(3)); g(p(1),p(2),p(3));h(p(1),p(2),p(3))];
        s = -J\F;
        p = p+s;
        n
        [p,s]
        disp(' ')
     end
     format loose
     format short