tol = 1 * 10^-20;

disp(sprintf('*** 5.9c ***'));
xk = [1; 2];
i = 0;
f_xk = [xk(1) + xk(2) - 2*xk(1)*xk(2); xk(1)^2 + xk(2)^2 - 2*xk(1) + 2*xk(2) + 1];
J_xk = [1-2*xk(2) 1-2*xk(1); 2*xk(1)-2 2*xk(2)+2];
sk = J_xk\(f_xk*-1);
xk1 = xk + sk;
disp(sprintf(' '));
disp(sprintf('i   x_k       x_k+1'));
while abs(xk1(1) - xk(1)) > tol
  disp(sprintf('%i   %f  %f', i, xk(1), xk1(1)));
  disp(sprintf('    %f  %f', xk(2), xk1(2)));
  i=i+1;
  xk = xk1;
  f_xk = [xk(1) + xk(2) - 2*xk(1)*xk(2); xk(1)^2 + xk(2)^2 - 2*xk(1) + 2*xk(2) + 1];
  J_xk = [1-2*xk(2) 1-2*xk(1); 2*xk(1)-2 2*xk(2)+2];
  sk = J_xk\(f_xk*-1);
  xk1 = xk + sk;
end
[xk1(1) + xk1(2) - 2*xk1(1)*xk1(2); xk1(1)^2 + xk1(2)^2 - 2*xk1(1) + 2*xk1(2) + 1] % solution should be [0,0]
disp(sprintf(' '));


disp(sprintf('*** 5.9d ***'));
xk = [3;-2];
i=0;
f_xk = [xk(1)^3 - xk(2)^2; xk(1) + xk(1)^2*xk(2) - 2];
J_xk = [3*xk(1)^2  -2*xk(2);  1+2*xk(1)*xk(2)  xk(1)^2];
sk = J_xk\(f_xk*-1);
xk1 = xk+sk;
disp(sprintf(' '));
disp(sprintf('i   x_k       x_k+1'));
while abs(xk1(1) - xk(1)) > tol
  disp(sprintf('%i   %f  %f', i, xk(1), xk1(1)));
  disp(sprintf('    %f  %f', xk(2), xk1(2)));
  i=i+1;
  xk = xk1;
  f_xk = [xk(1)^3 - xk(2)^2; xk(1) + xk(1)^2*xk(2) - 2];
  J_xk = [3*xk(1)^2  -2*xk(2);  1+2*xk(1)*xk(2)  xk(1)^2];
  sk = J_xk\(f_xk*-1);
  xk1 = xk+sk;
end
[xk1(1)^3 - xk1(2)^2; xk1(1) + xk1(1)^2*xk1(2) - 2] % solution should be [0,0]
disp(sprintf(' '));


disp(sprintf('*** 5.9e ***'));
xk = [2*pi;-1*pi]
i=0;
f_xk = [2*sin(xk(1)) + cos(xk(2)) - 5*xk(1); 4*cos(xk(1)) + 2*sin(xk(2)) - 5*xk(2)];
J_xk = [2*cos(xk(1))-5 -sin(xk(2)); -4*sin(xk(1)) 2*cos(xk(2))-5];
sk = J_xk\(f_xk*-1);
xk1 = xk+sk;
disp(sprintf(' '));
disp(sprintf('i   x_k       x_k+1'));
while abs(xk1(1) - xk(1)) > tol
  disp(sprintf('%i   %f  %f', i, xk(1), xk1(1)));
  disp(sprintf('    %f  %f', xk(2), xk1(2)));
  i=i+1;
  xk = xk1;
  f_xk = [xk(1)^3 - xk(2)^2; xk(1) + xk(1)^2*xk(2) - 2];
  J_xk = [3*xk(1)^2  -2*xk(2);  1+2*xk(1)*xk(2)  xk(1)^2];
  sk = J_xk\(f_xk*-1);
  xk1 = xk+sk;
end
[xk1(1)^3 - xk1(2)^2; xk1(1) + xk1(1)^2*xk1(2) - 2] % solution should be [0,0]