%% This is the MATLAB program to find the steady state in BZ
%% It takes the parameters from the DYNARE BZfirst.mod and gives the steady state
%%
%%
function [ys,params,check] = tt_steadystate(ys,exo,M_,options_)

NumberOfParameters = M_.param_nbr;
for ii = 1:NumberOfParameters
  paramname = M_.param_names{ii};
  eval([ paramname ' = M_.params(' int2str(ii) ');']);
end
% initialize indicator
check = 0;


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


R            = 1/betta;
Lambda       = betta;

% Targets
p_mom        =  0.03;     %Steady state default. 3% default quarterly. 
Q=1;
etaq=1;
p            = p_mom;
rho          = 1 + 0.02/4; %0.02/4; %    %steady state premium
phie_mom     = 2; 
phie         = phie_mom;



npar=22;
x0=ones(npar,1);
x1=[0.98,0.3, 0.49, 0.008, 0.98, 0.30, 0.2, 0.49,1.02];
x0=[x0;x1'];
fun = @(theta) DataSS(theta,Q,etaq,delt,z,alphha,ksie,xpri,phie, cA, Databar, betta, gammma, chil,epsl,p_mom,rho,phie_mom,R);
options=optimset('MaxFunEvals',100000000000000000000);
theta     = fsolve(fun,x0,options);
F= DataSS(theta,Q,etaq,delt,z,alphha,ksie,xpri,phie, cA, Databar, betta, gammma, chil,epsl,p_mom,rho,phie_mom,R);

As=theta(1);
Roc=theta(2);
K=theta(3);
Data=theta(4); 
MPK=theta(5);
Lab=theta(6); 
A=theta(7); 
DR=theta(8);
C=theta(9); 
Y=theta(10); 
Le=theta(11);
Ne=theta(12); 
I=theta(13); 
Atil=theta(14); 
L_over_K=theta(15); 
Ce_over_K=theta(16); 
C_over_K=theta(17); 
Y_over_K=theta(18); 
Y_over_C=theta(19); 
Ne_over_K=theta(20); 
Z=theta(21);
Ce=theta(22);
sigmae=theta(23);
sigma_omega=theta(24);
fnGam=theta(25);
fnG=theta(26);
DGam=theta(27);
DG=theta(28);
mon=theta(29); 
omega  = exp(theta(30));
Ra=theta(31); 
sigma_omega_ss=sigma_omega;
%Rest variables
W=(1-alphha)*Y/Lab;
Uc=(C*(1-gammma))^(-1) -betta*gammma*(C*(1-gammma))^(-1);
Abar=1;
EfLab=Lab*(1-xpri);
rhoRex=rho*(R+(1-mon)*fnG*DR);
psi=1;
YY=1;   Loansloans=1; A=1; YY=1; Loansloans=1; KK=1; CC=1; LabLab=1; II=1; RR=1; RkRk=1; RlRl=1;
WW=1; ZZ=1;  NeNe=1; UcUc=1; QQ=1;phiephie=1;  TT=1; GG=1; rhorho=1; 
pp=1; rhoRexrhoRex=1; spreadRkRk=RkRk-RR;
spreadRlR=RlRl-RR; premium =Ra-R;  S=0; Sdash=0; premiumpremium=Ra-R;
sigma_pref=1;
zeta=1;


spread = Ra-R;

sigma_omega= sigma_omega_ss;
bankrupt = p;
X=1;



YY=1;
KK=1;
II=1;
CC=1;
WW=1;
LabLab=1;
QQ=1;
RR=1;
RkRk=1;
RnRn=1;
PIPPIP=1;
spreadspread=0;
NeNe=1;
rhorho=1;
psipsi=1;
phiephie=1;
OUTGAPOUTGAP=1;
LeLe =1;
%% END OF THE MODEL SPECIFIC BLOCK.

%% DO NOT CHANGE THIS PART.
%%
%% Store model parameters
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

end                                                  % End of the loop.
%%