I am trying to solve my system with 5 equations/ and 5 unknowns but it seems like my code is giving me a wrong answer an

[answerhappygod]

I am trying to solve my system with 5 equations/ and 5 unknowns
but it seems like my code is giving me a wrong answer and it takes
all iterations it does not stop. please anyone can help to fix. I
am using newton Raphson method c = ci - inv(g)*f
here is my code
clear all
close
clc
bK=1*10^-3; %bulk concentration
nj= 5; %number of node points
x= 0.5; % the thickness
h=x/nj;
a=1*10^-5;
b=0.0;
k=10^-5;
v=0.01;
c1(1) = 0.001;
c2(1) = 0.001;
c3(1) = 0.001;
c4(1) = 0.001;
c5(1) = 0.001;
c0=[c1(1);c2(1);c3(1);c4(1);c5(1)];
tol=[1e-6;1e-6;1e-6;1e-6;1e-6];
c_old = c0;
for i=1:1000
ci=[0.01;
0.01;
0.01;
0.01;
0.01];
f = zeros (5,1)
f(1,1)=-a+b*c1*((-c3+4*c2-3*c1)/(h)^2)-k*c1;
f(2,1)=v*((c3-c1)/2*h)-b*((c3-c1)/2*h)^2-a+b*c2*(c3-2*c2+c1)/(h)^2;

f(3,1)=v*((c4-c2)/2*h)-b*((c4-c2)/2*h)^2-a+b*c3*(c4-2*c3+c2)/(h)^2;

f(4,1)=v*((c5-c3)/2*h)-b*((c5-c3)/2*h)^2-a+b*c4*(c5-2*c4+c3)/(h)^2;

f(5,1)=c5-bK;
g = zeros (5,5)
g(1,1)=(3*(a + b*c1))/(2*h) - k + (b*(3*c1 - 4*c2 +
c3))/(2*h);
g(1,2)=-(2*(a + b*c1))/h;
g(1,3)=(a + b*c1)/(2*h);
g(1,4)=0;
g(1,5)=0;
g(2,1)=-(a + b*c2)/h^2 - v/(2*h) - (b*(2*c1 - 2*c3))/(4*h^2);
g(2,2)=(2*(a + b*c2))/h^2 - (b*(c1 - 2*c2 + c3))/h^2;
g(2,3)=v/(2*h) - (a + b*c2)/h^2 + (b*(2*c1 - 2*c3))/(4*h^2);
g(2,4)=0;
g(2,5)=0;
g(3,1)=0;
g(3,2)=-(a + b*c3)/h^2 - v/(2*h) - (b*(2*c2 - 2*c4))/(4*h^2);
g(3,3)= (2*(a + b*c3))/h^2 - (b*(c2 - 2*c3 + c4))/h^2;
g(3,4)= v/(2*h) - (a + b*c3)/h^2 + (b*(2*c2 - 2*c4))/(4*h^2);
g(3,5)=0;
g(4,1)=0;
g(4,2)=0;
g(4,3)=-(a + b*c4)/h^2 - v/(2*h) - (b*(2*c3 - 2*c5))/(4*h^2);
g(4,4)=(2*(a + b*c4))/h^2 - (b*(c3 - 2*c4 + c5))/h^2;
g(4,5)=v/(2*h) - (a + b*c4)/h^2 + (b*(2*c3 - 2*c5))/(4*h^2);
g(5,1)=0;
g(5,2)=0;
g(5,3)=0;
g(5,4)=0;
g(5,5)=1;
cn = ci - (g\f);

if abs(cn - ci) > tol

ci = cn;
i = i+1;
if abs(cn - ci) < tol

ci = cn;
break;
end
end



end
plot(cn)
cn
i