
% Provides steady state values to dynamic problem solver
 
function [ys,params,check] = dyn_main_steadystate(ys,exo,M_,options_)
%function [ys,params,check] =NK_baseline_steadystate(ys,exo,M_,options_) is an example in https://git.dynare.org/Dynare/dynare/-/wikis/BreakingFeaturesIn4.6
  
NumberOfParameters = M_.param_nbr;
for ii = 1:NumberOfParameters
  paramname = M_.param_names{ii};
  eval([ paramname ' = M_.params(' int2str(ii) ');']);
end
% initialize indicator
check = 0;


  % Initial guess
   %     X0 = [0.34; 0.34; 0.02; 0.02; 1.8; 1.8];
 X0 = [0.34; 0.34; 1; 1; 2.5; 2.5]; %last used by Diego

  
  options = optimset('fsolve');
  options = optimset(options,'Display','off');
  X  = fsolve(@(X) dyn_det(X),X0,options);
  
  nH = X(1);
  nF = X(2);
  xH = X(3);
  xF = X(4);
  wH = X(5);
  wF = X(6);

  [cH, YH, WH, qH, pH, dH, mH, yH, lHy, lHq, IH, UcH, kH. RH, ...
   cF, YF, WF, qF, pF, dF, mF, yF, lFy, lFq, IF, UcF, kF, RF] = sseqs(X);

    bH  = 0; %eq 43 in dyn_main
    bF  = 0;  

    r   = bet; %eq 16 in dyn_main
    
    zH = 1; %eq 19 in dyn_main
    zF = 1; 
    
  % Terms of trade %eq 25 in dyn_main
    totH = 1;
    totF = 1;
    
  % Price index (P) % eq 26 in dyn_main
    PH = 1; %(alfH*IH*pH^nu + (1-alfH)*IF*tau*pF^nu)^(1/nu);
    PF = 1; %(alfF*IF*pF^nu + (1-alfF)*IH*tau*pH^nu)^(1/nu);
  
  % Real exchange rate %eq 27 in dyn_main
    rerH = PF./PH;
    rerF = PH./PF;
    
  % Exports %eq 21 in dyn_main
    exH = tau*IH.*mF.*pH;
    exF = tau*IF.*mH.*pF;
    
  % GDP %eq20 in dyn_main
    gdpH = YH + exH - exF; %gdpH=IH*pH*(dH+mF);
    gdpF = YF + exF - exH; %gdpF=IF*pF*(dF+mH);
    
  % Net exports %eq28 in dyn_main
    nxH = (exH-exF)./gdpH;
    nxF = (exF-exH)./gdpF;
    
  % Measured aggregates (coincide in ss)
    mYH = YH; %eq22 dyn_main
    mYF = YF;
    mcH = cH; %eq24 dyn_main
    mcF = cF;
    mxH = xH; %eq 23 dyn_main
    mxF = xF;
    
  % Relative output (RY) %eq29 in dyn_main
    RYH = gdpH./gdpF;
    RYF = gdpF./gdpH;
    
  % Relative consumption (RC) %eq 30 in dyn_main
    RCH = mcH./mcF;
    RCF = mcF./mcH;
    
  % Quality adjusted variables %eq 31-42 in dyn_main
    qadjH = 1;
    qYH   = YH;
    qcH   = cH;
    qxH   = xH;
    qtotH = totH;
    qPH   = PH;
    qrerH = rerH;
    qexH  = exH;
    qgdpH = gdpH;
    qnxH  = nxH;
    qRYH  = RYH;
    qRCH  = RCH;
    qadjF = 1;
    qYF   = YF;
    qcF   = cF;
    qxF   = xF;
    qtotF = totF;
    qPF   = PF;
    qrerF = rerF;
    qexF  = exF;
    qgdpF = gdpF;
    qnxF  = nxF;
    qRYF  = RYF;
    qRCF  = RCF;
   

params=NaN(NumberOfParameters,1);
for iter = 1:length(M_.params) %update parameters set in the file
  eval([ 'params(' num2str(iter) ') = ' M_.param_names{iter} ';' ])
end

NumberOfEndogenousVariables = M_.orig_endo_nbr; %auxiliary variables are set automatically
for ii = 1:NumberOfEndogenousVariables
  varname = M_.endo_names{ii};
  eval(['ys(' int2str(ii) ') = ' varname ';']);
end