Solving Nonlinear Equations Error

28 vues (au cours des 30 derniers jours)
Tim
Tim le 4 Juil 2011
I have a problem with solving two nonlinear equations using the newton method. My problem is that the code works with some values but with others I only receive the following errors:
??? Error using ==> erfc Input must be real and full.
Error in ==> normcdf at 94 p(todo) = 0.5 * erfc(-z ./ sqrt(2));
Error in ==> test4 at 14 f(1) = x*normcdf((log(x/18.5)+(0.005 +0.5*y^2))/y) - 18.5 * exp(-0.005)*normcdf((log(x/185)+(0.03 +0.5*y^2))/y)-y-0.5;
My Code is:
% Newton Raphson solution of two nonlinear algebraic equations
% set up the iteration
error1 = 1.e8;
xx(1) = 19.0; % initial guesses
xx(2) = 0.02;
iter=0;
itermax=3000.
% begin iteration
while error1>1.e-12
iter=iter+1;
x = xx(1);
y = xx(2);
% calculate the functions
f(1) = x*normcdf((log(x/18.5)+(0.005 +0.5*y^2))/y) - 18.5 * exp(-0.005)*normcdf((log(x/185)+(0.03 +0.5*y^2))/y)-y-0.5;
f(2) =(x/0.5)*normcdf((log(x/18.5)+(0.005 +0.5*y^2))/y)*y-1.42;
% calculate the Jacobian
J(1,1) = normcdf((log(x/18.5)+(0.005 +0.5*y^2))/y);
J(1,2) = x*1*(exp((-((log(x/18.5)+(0.005 +0.5*y^2))/y)^2)/2))/((sqrt(2)*pi*1)*x*y);
J(2,1) = normcdf((log(x/18.5)+(0.005 +0.5*y^2))/y)*y +(x*1*(exp((-((log(x/18.5)+(0.005 +0.5*y^2))/y)^2)/2))/((sqrt(2)*pi*1)*x*y));
J(2,2) = x*normcdf((log(x/18.5)+(0.005 +0.5*y^2))/y)+x*y*(exp((-((log(x/18.5)+(0.005 +0.5*y^2))/y)^2)/2))/((sqrt(2)*pi)*x*y)*((log(x/(18.5*exp(-0.005)))))/(y^2)+0.5;
% solve the linear equations
y = -J\f';
% move the solution, xx(k+1) - xx(k), to xx(k+1)
xx = xx + y';
% calculate norms
error1=sqrt(y(1)*y(1)+y(2)*y(2));
error(iter)=sqrt(f(1)*f(1)+f(2)*f(2));
ii(iter)=iter;
if (iter > itermax)
error1 = 0.;
s=sprintf('****Did not converge within %3.0f iterations.****',itermax);
disp(s)
end
% check if error1 < 1.e-12
end
x = xx(1);
y = xx(2);
f(1) = x*normcdf((log(x/18.5)+(0.005 +0.5*y^2))/y) - 18.5 * exp(-0.005)*normcdf((log(x/185)+(0.005 +0.5*y^2))/y)-y-0.5;
f(2) =(x/0.5)*normcdf((log(x/18.5)+(0.005 +0.5*y^2))/(y*1))*y-1.42;

Réponses (0)

Catégories

En savoir plus sur Symbolic Math Toolbox dans Help Center et File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by