% computes the steady state of MSSF.mod. 
% Taken from an example file by Kevin Moran
% May 2016, Abdellah Manadir

function [ys,check] = MSSF_steadystate(ys,exe)
global M_

%% DO NOT CHANGE THIS PART.
%%
%% Here we load the values of the deep parameters in a loop.
%%
NumberOfParameters = M_.param_nbr;                            % Number of deep parameters.
for i = 1:NumberOfParameters                                  % Loop...
  paramname = deblank(M_.param_names(i,:));                   %    Get the name of parameter i. 
  eval([ paramname ' = M_.params(' int2str(i) ');']);         %    Get the value of parameter i.
end                                                           % End of the loop.  
check = 0;
%%
%% END OF THE FIRST MODEL INDEPENDENT BLOCK.


%% THIS BLOCK IS MODEL SPECIFIC.
%%
%% Here the user has to define the steady state.
%%


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% steady-state computations:

lSa_ss = log(gam) ;

lSx_ss = 0 ;

lSup_ss = 0 ;

lSuw_ss = 0 ;

SMp_ss = elas/(1-elas) ;

SMw_ss = elasw/(1-elasw);

SR_ss = gam/bet ;

SRn_ss = SR_ss ;

Spi_ss = 1 ;

SP_ss = 1 ;

SPstar_ss = 1 ;

Sfi_ss = 1/SMp_ss ;

SU_ss = 1 ;

Sq_ss = 1/(1-(xi*0.5)*((gam -1)^2) - gam*(1-bet)*xi*(gam -1)) ;

SRh_ss = ((gam/bet) - (1-delta))*Sq_ss ;

SW_ss = (Sfi_ss/((SRh_ss^alfa)*(alfa^(-alfa))*((1-alfa)^(alfa - 1))))^(1/(1-alfa)) ;

SWstar_ss = (((1- (probaw/(gam^(1-elasw))))/(1 - probaw))^(1/(1-elasw)))*SW_ss ;

SWh_ss = ((1-bet*probaw*(gam^elasw))/(1 - bet*probaw*(gam^(elasw-1))))*(SWstar_ss/SMw_ss) ;

GY_ss= gyshare ;

SYL = (SW_ss)/(Sfi_ss*(1-alfa)) ;

SKL = gam*(SYL^(1/alfa)) ;

SCL = (1- GY_ss)*SYL - ((1 - ((1-delta)/gam))/(1 - (xi*0.5)*((gam-1)^2)))*SKL ;

SlambdaL = (1 - (h/gam))*SCL ;

SL_ss = (SWh_ss/SlambdaL)^(1/(1+phi)) ;

Slambda_ss = SlambdaL*SL_ss ;

SC_ss = SCL*SL_ss ;

SY_ss = SYL*SL_ss;

SK_ss = SKL*SL_ss ;

SI_ss = ((1 - ((1- delta)/gam))/(1 - (xi*0.5)*((gam-1)^2)))*SK_ss ;

A1_ss = ((SW_ss^elasw)*SL_ss)/(Slambda_ss*(1 - bet*probaw*(gam^(elasw -1)))) ;

A2_ss = ((SW_ss^elasw)*SL_ss*SWh_ss*SMw_ss)/(Slambda_ss*(1 - bet*probaw*(gam^elasw ))) ;

A3_ss = (SY_ss/(Slambda_ss*(1-bet*probap))) ;

A4_ss= A3_ss ;



% Now ready to assign these "*_ss values to all variables


Slambda = Slambda_ss

SC     = SC_ss ;

SR     = SR_ss ; 
 
SL    = SL_ss  ;   

SWh   = SWh_ss ;

SP    = SP_ss ;

SI     = SI_ss ;

SRh         = SRh_ss ;

SU          = SU_ss ;

SK          = SK_ss ;

Sq          = Sq_ss ;

A1         = A1_ss ;

A2         = A2_ss ;

A3        = A3_ss ;

A4        = A4_ss ;

SW         = SW_ss ;

SWstar      =SWstar_ss

SMw         = SMw_ss ;


SY          = SY_ss ;

SMp         = SMp_ss ;

Sfi         = Sfi_ss ;

SPstar     = SPstar_ss ;

SRn        = SRn_ss ;

Spi         =Spi_ss ;

lSuw        = lSuw_ss ;

lSup        = lSup_ss ;

lSa         = lSa_ss ;

GY          = GY_ss ;

lSx         = lSx_ss ;       



%%
%% END OF THE MODEL SPECIFIC BLOCK.


%% DO NOT CHANGE THIS PART.
%%
%% Here we define the steady state values of the endogenous variables of
%% the model.
%%
NumberOfEndogenousVariables = M_.endo_nbr;                    % Number of endogenous variables.
ys = zeros(NumberOfEndogenousVariables,1);                    % Initialization of ys (steady state).
for i = 1:NumberOfEndogenousVariables                         % Loop...
  varname = deblank(M_.endo_names(i,:));                      %    Get the name of endogenous variable i.                     
  eval(['ys(' int2str(i) ') = ' varname ';']);                %    Get the steady state value of this variable.
end                                                           % End of the loop.
%%
%% END OF THE SECOND MODEL INDEPENDENT BLOCK.