
function [ys,check]=dsge_steadystate(ys,exo)

global M_ 

% read out parameters to access them with their name
NumberOfParameters = M_.param_nbr;
for ii = 1:NumberOfParameters
  paramname = deblank(M_.param_names(ii,:));
  eval([ paramname ' = M_.params(' int2str(ii) ');']);
end
% initialize indicator
check = 0;

%delta1=1/betta-(1-delta0);
A=1;
Z=1;
v=1;
u=1;
S=1;
pssi=psis;
tauc=taucs;
taun=tauns;
tauk=tauks;
pi=pi_target;
delta1 = ((1+pi_target)/betta-pi_target)/betta + delta0 -1;
%i = i_target;
%i=1/betta*(1+pi)-1;
i = rho3/(1-rho2)*cm;
%i = i_target;
q = 1/((1+pi)*(1+i-pi));
%i_bc = 1/betta*(1+pi)-1;
i_bc = i + 0.5;
%R = q * (1-betta*(1-delta0))/(betta*(1-tauk));
%R = q/((1-tauk)) * ((1+i-pi)/betta-1+delta0);
R = q*delta1/(1-tauk);
pi_star=(((1+pi)^(1-epsilonp)-fip*(1+pi)^(etap*(1-epsilonp)))/(1-fip))^(1/(1-epsilonp))-1;
v_p=(1-fip)*((1+pi_star)/(1+pi))^(-epsilonp)/(1-fip*(1+pi)^(epsilonp*(1-etap)));
mc=(epsilonp-1)/epsilonp*(1+pi_star)/(1+pi)*(1-fip*betta*(1+pi)^(epsilonp*(1-etap)))/(1-fip*betta*(1+pi)^((1-epsilonp)*(etap-1)));
w=(1-alfa)*mc*(alfa*mc/R)^(alfa/(1-alfa));
w_star=w*((1-fiw*(1+pi)^((epsilonw-1)*(1-etaw)))/(1-fiw))^(1/(1-epsilonw));
N_d=((epsilonw-1)/epsilonw*w_star/(psis*v)*(1-taun)/(1+tauc)/(1-cb)*(1-betta*cb)*((1-omegga)*(alfa*mc/R)^(alfa/(1-alfa))-delta0*(alfa*mc/R)^(1/(1-alfa)))^(-1)*(w/w_star)^(-epsilonw*ka)*((1-fiw*betta*(1+pi)^(epsilonw*(1+ka)*(1-etaw)))/(1-fiw*betta*(1+pi)^((epsilonw-1)*(1-etaw)))))^(1/(1+ka));
% C =  N_d*((1-omegga)*A*(alfa*mc/R)^(alfa/(1-alfa))-delta0*(alfa*mc/R)^(1/(1-alfa)));
K_hat=(alfa*mc/R)^(1/(1-alfa))*N_d;
K=K_hat;
I=delta0*K;
Y=A*K_hat^alfa*N_d^(1-alfa)/v_p;
M = cm * Y;
B = c2*Y - M;
G= omegga*Y;
BC=I-B;
C = Y - I - G;
lamb = 1/C*(1-betta*cb)/(1-cb)/(1+tauc);
i_m = 1/lamb - (v*cgamma_m*(M^chi_m)/lamb)-1;
x1=lamb*mc*Y/(1-fip*betta*(1+pi)^((1-etap)*epsilonp));
x2=lamb*Y/(1-fip*betta*(1+pi)^((etap-1)*(1-epsilonp)));
h1=psis*(w/w_star)^(epsilonw*(1+ka))*N_d^(1+ka)/(1-fiw*betta*(1+pi)^(epsilonw*(1+ka)*(1-etaw)));
h2=(1-taun)*lamb*(w/w_star)^epsilonw*N_d/(1-fiw*betta*(1+pi)^((epsilonw-1)*(1-etaw)));

check=0;


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

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