'fsolve' producing error message

Hi the community,

I've been struggling with the following code for a few days now, not sure what is wrong with it. I would appreciate if you have a look at the code. It is very simple code, a simulation exercise. I am getting the following message

• Error using trustnleqn (line 28)
• Objective function is returning undefined values at initial point. FSOLVE cannot continue.
• Error in fsolve (line 366)
• [x,FVAL,JACOB,EXITFLAG,OUTPUT,msgData]=...
• Error in paper2_v61 (line 76)
• out1=fsolve(@productivity,x0,options);

Thank you!

The code: clc clear all

global A tauc tauw a b sigma alpha S phic phiw rc rw
format short
T=200;

alpha=0.75;   
tauw=1;        
tauc=2;        
rc=3;       
rw=6;      

sigma=3;        
phic=0.5;        
phiw=1;          

a=1;            
b=0.81;         

%Initial values
Ac0=1;           
Aw0=5;

Sw0=7000;
Sc0=35000;
Scb=10^6;
Theta=zeros(T,1);
Ec=zeros(T,1); 
Ew=zeros(T,1);
X=zeros(T,1);
Yc=zeros(T,1); 
Yw=zeros(T,1);
Ac=zeros(T+1,1); 
Aw=zeros(T+1,1);
Ac(1,1)=Ac0;
Aw(1,1)=Aw0;
Pc=zeros(T,1); 
Pw=zeros(T,1);
pc=zeros(T,1); 
pw=zeros(T,1);
cc=zeros(T,1); 
cw=zeros(T,1);
Nc=zeros(T,1); 
Nw=zeros(T,1);
Sc=zeros(T+1,1); 
Sw=zeros(T+1,1);
Sc(1,1)=Sc0;     
Sw(1,1)=Sw0;

x0=[0.5,0.5,100];
for t=1:T

      A=[Ac(t,1),Aw(t,1)];
      S=[Sc(t,1),Sw(t,1)];

      options=optimset('TolX',1e-9,'TolFun',1e-9);
      out1=fsolve(@productivity,x0,options);
      Nc(t,1)=out1(2);
      Nw(t,1)=1-out1(2);
      Theta(t,1)=out1(1);
      Pc(t,1)=out1(3);
      Ac(t+1,1)=(1+tauc*Theta(t,1))*Ac(t,1);
      Aw(t+1,1)=(1+tauw*Theta(t,1))*Aw(t,1);
      X(t,1)=Sc(t,1)*(1/Ac(t+1,1));
      Sc(t,1)=Sc(t,1)+X(t,1);

      Pw(t,1)=(1-Pc(t,1)^(1-sigma))^(1/(1-sigma));

      cc(t,1)=Pc(t,1)^(1/(1-alpha))*Ac(t+1,1)*(1-alpha);
      cw(t,1)=Pw(t,1)^(1/(1-alpha))*Aw(t+1,1)*(1-alpha);

      pc(t,1)=cc(t,1);
      pw(t,1)=cw(t,1);

      Ec(t,1)=(1-(cc(t,1)/phic)^1/rc)*Sc(t,1);
      Ew(t,1)=(1-(cw(t,1)/phiw)^1/rw)*Sw(t,1);

      Sc(t+1)=Sc(t,1)-Ec(t,1);
      Sw(t+1)=Sw(t,1)-Ew(t,1);
  end

Function

   function F=productivity(x)      

      global a b tauc tauw A alpha S sigma phiw phic rc rw

      Ac=(1+tauc*x(1))*A(1);
      Aw=(1+tauc*x(1))*A(2);
      Sc=S(1);
      Sw=S(2);

      cc=x(3)^(1/(1-alpha))*Ac*(1-alpha);
      cw=((1-x(3)^(1-sigma))^(1/(1-sigma)))^(1/(1-alpha))*Aw*(1-alpha);
      E1=((1-(cc/phic)^1/rc)*Sc)/((1-(cw/phiw)^1/rw)*Sw);
      E2=(x(3)/((1-x(3)^(1-sigma))^(1/(1-sigma))))^(alpha/(alpha-1))*(Ac/Aw)^(sigma*(1-alpha)-1)*(cc/cw)^-sigma;

   F=[ -x(1)+a*(x(2)/((1+tauc)*A(1)))^b;
       -x(1)+a*((1-x(2))/((1+tauw)*A(2)))^b;
       E1-E2;     
       ];

Answer 1

[Alan Weiss]

You didn't include your call to fsolve so it is hard to know exactly what you are doing. We don't know your options, your initial point x0, or your fsolve call.

Nevertheless, the error is clear. You have an initial point for fsolve, I suppose it is x0. When evaluating

 productivity(x0)

fsolve found that there were NaN in the result. Either choose an initial point that does not produce this behavior, or debug your program so that it does not produce NaN results.

Alan Weiss

MATLAB mathematical toolbox documentation

[mrrox]

Hi Alan, Many thanks, in fact there is a call for fsolve and x0

x0=[0.5,0.5,100];
and 

out1=fsolve(@productivity,x0,options);

please have a look. Thanks

Answer 2

[Matt J]

When the loop reaches t=46, I find the following

      K>> productivity(x0)

      ans =

         -0.0917
         -0.4998
             Inf

I assume the problem is now clear. productivity() is not returning finite real values at the initial point x0. It is because Sc has reached a super-huge number

        K>> Sc

          Sc =

            3.2735e+301

making E1 numerically infinite

        K>> E1

        E1 =

           Inf

[mrrox]

Hello Matt, thanks again. yes it worked after changing the parameter values and initial values.