Newton Raphson method multiple varibles

I'm trying to solve a system of nonlinear equations for x,y and z. But my jocobian matrix and function vectors contain Nan values and im not sure why. It doesnt seem to reconginze the values in either matrix and I'm not sure why.

% Problem 2c
% System of Nonlinear Equations , Newton-Raphson Method
clear;clc

P=2; % Maybe need to convert this
K_p=-2.88^10;

% To solve sequence of nonlinear equations using Newton Raphson
% Initial conditions
X0=[1;1;1/2]; %initial condition value of x1;x2;x3

maxIter=1e3; % maximum iteration...SHOULD BE HIGHER
tolX=1e-3; % Tolerance for error...SHOULD BE HIGHER

%% Computation using Newton Raphson
X=X0;
Xold=X0;

% f=zeros([1,3])
% j=zeros([3,3])
% f=[];
% j=[]
for i=1:maxIter

%     [f,j]=myfunction_to_solve_systemofNonlinearEquations(X);
x=X(1);
y=X(2);
z=X(3);

%Define f(x)
fval(1,1)=(x.*y.^(1/2))./z.*2.^(1/2)-(-2.88^10); % note p=2atm and K_p=-2.88^10
fval(2,1)=(z+x)./(x+1/2.*y+z)-1/2;
fval(3,1)=x+y+z-(1-1e-5);  % Wont converge if set to 1 so close to 1

%Define jacobian
jac=[(y^(1/2).*2^(1/2))./z, (2^(1/2).*x)./(2.*z.*y^(1/2)),(-2^(1/2).*x.*y^(1/2))./z^2;
    (2.*y)/(2.*z+2.*x+y)^2,-(z+x)/(2.*(z+x+y/2)^2),(2.*y)/(2.*z+2.*x+y)^2;
    1, 1, 1];

    X=X-inv(jac)*fval;
    err(:,i)=abs(X-Xold);
    Xold=X;
    if (err(:,i)<tolX)
        break; % end for loop
        
    end


end
display(X)
fprintf('Number of iterations: %d \n' ,i)

I keep getting a Warning: Matrix is singular, close to singular or badly scaled. Results may be inaccurate.

RCOND = NaN.

and a anser of :

X =

NaN

NaN

NaN

Number of iterations: 1000

[Torsten]

If you output jac, you will notice that the first and third column approach each other quite rapidly.

I didn't examine why this happens, but it makes the Jacobian singular and Newton's method fail.

You might want to examine your (rewritten) problem with MATLAB's "fsolve" first:

format long
fun = @(x,y,z)[x*sqrt(2*y)+2.88^10*z;(z+x)-1/2*(x+0.5*y+z);x+y+z-1];
x0 = [1;1;1/2];
options = optimset('MaxFunEvals',10000,'MaxIter',10000);
sol = fsolve(@(x)fun(x(1),x(2),x(3)),x0,options)

Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.

sol = 3×1

0.333343138071722 0.666666666666667 -0.000009804738388

fun(sol(1),sol(2),sol(3))

ans = 3×1

1.0e-14 * -0.716093850883226 -0.005551115123126 -0.011102230246252