Fsolve within a for loop

Hi. I need to get each column of z_ss in each iteration j after fsolve. The problem is I do not know how to adjust the function in the separate file accordingly, such that in the first iteration it is considering z_ss(:,1) , in the second z_ss(:,2) ect...

This is the code:

clear 
clc 

r_annual = 6.5;
rss = r_annual/400;
delta = 0.021; 
sigma = 0.15;
alpha = 0.36; 
rho = 0.95;
sigma_theta = 0.01;
theta_ss = 1;  %(1)
M = 2;

g  = 0:0.1:0.2;

z_ss = ones(M,length(g));
gamma = zeros(length(g),1);
beta = zeros(length(g),1);

for j = 1:length(g)
    gamma(j) = 1+g(j);
    beta(j) = gamma(j)/(1+rss);
    params = [sigma alpha rss beta(j) r_annual delta rho gamma(j) M];

options=optimset('display','iter');%,'TolFun',1e-4,'TolX', 1e-4);
options=optimset(options,'MaxFunEvals',10000,'TolX',1e-8);

z_steady = @(z_ss)steadygrow(params, z_ss);

z_ss(:,j) = fsolve(z_steady, z_ss, options); %(2)

end;

and the function is :

function H = steadygrow(params, z_ss)
%params = [sigma alpha rss beta(j) r_annual delta rho gamma(j) M];

h = zeros(params(9),1);
s1 = zeros(params(9),1);
s2 = zeros(params(9),1);

for i = 1:params(9)
    s1(i) = (1-normcdf(z_ss(i)-params(1)))*(params(4)*(1-params(6)))^(i-1)*(params(8))^(-i);
    s2(i) = (1-normcdf(z_ss(i)))*(params(4)*(1-params(6)))^(i-1);
end;

      S1 = log(sum(s1));
      S2 = log(sum(s2));

      for i =1:params(9)
           h(i) = z_ss(i) - 1/params(1)*( log(1-params(2)) + S1 - S2 + 1/2*params(1)^2 + i * log (params(8)) );
      end;
          H = h;
  end

Any help would be appreciated.

[Torsten]

You mean how to supply z_ss(:,j) as initial guess for z_ss in the j'th call to fsolve ?

Answer 1

[Matt J]

I'm not completely sure what you want for the initial guess, x0, but the following should be a starting point.

options=optimset('display','iter');%,'TolFun',1e-4,'TolX', 1e-4);
options=optimset(options,'MaxFunEvals',10000,'TolX',1e-8); 

 for j = 1:length(g)

      gamma(j) = 1+g(j);
      beta(j) = gamma(j)/(1+rss);
      params = [sigma alpha rss beta(j) r_annual delta rho gamma(j) M];

     z_steady = @(x)steadygrow(params, x);

     x0=ones(M,1);

     z_ss(:,j) = fsolve(z_steady, x0, options); %(2)

end;

[Matt J]

Adelia's response moved here

Thank you. That worked :)