function [ys,check] = model_steadystate(ys,exe);
global M_ lgy_ options_

if isfield(M_,'param_nbr') == 1
NumberOfParameters = M_.param_nbr;
for i = 1:NumberOfParameters
  paramname = deblank(M_.param_names(i,:));
  eval([ paramname ' = M_.params(' int2str(i) ');']);
end
check = 0;
end

ETA       = 0.022;
EPSILON   = 0.39; 
ZETAi     = 0.99; 
ZETAs     = 0.985;
BETTA     = 0.97;
OMEGA     = 0.1;
KAPPA_bar = 0.75; 
N         = 0.5;
SIGMA     = 2; 
RHO_a     = 0.9;
RHO_kappa = 0.9; 
SIGMANUi  = 0.005;
SIGMANUk  = 0.011;
Z1 =  1 ;

BETAi_s = BETTA;
BETAs_s = BETTA;
rb_s = 1/BETTA;
ks_s = 1/(1+3*N);
ki_s = 4*ks_s;
q_s = (BETTA*OMEGA*Z1*(ks_s)^(OMEGA-1))/(1-BETTA);
bi_s = KAPPA_bar*q_s*ki_s;
bs_s = -(N*bi_s)/(1-N);
w_s = EPSILON*(ki_s)^(EPSILON-1);
r_s   = (1-EPSILON)*(ki_s)^-EPSILON;
cs_s = w_s + Z1 *(ks_s)^OMEGA + (1-rb_s)*bs_s;
ci_s = w_s + (r_s * ki_s) +(1-rb_s)*bi_s;
mu_s = (BETTA/(1-KAPPA_bar))*(((q_s + r_s)/q_s)-rb_s);
a_s  = 1; 
kappa_t_s = 1;

BETAi = BETTA;
BETAs = BETTA;
ci = log(ci_s);
cs = log(cs_s);
ki = log(ki_s);
ks = log(ks_s);
r = log(r_s);
rb = log(rb_s);
w = log(w_s);
a = log(a_s);
kappa_t = log(kappa_t_s);
bi = bi_s;
bs = bs_s;
%the value of bs_s is negative, therefore no log()
mu = log(mu_s);
q = log(q_s);


for iter = 1:length(M_.params)
  eval([ 'M_.params(' num2str(iter) ') = ' M_.param_names(iter,:) ';' ])
end

if isfield(M_,'param_nbr') == 1

if isfield(M_,'orig_endo_nbr') == 1
NumberOfEndogenousVariables = M_.orig_endo_nbr;
else
NumberOfEndogenousVariables = M_.endo_nbr;
end
ys = zeros(NumberOfEndogenousVariables,1);
if options_.linear
    
else
for i = 1:NumberOfEndogenousVariables
  varname = deblank(M_.endo_names(i,:));
  eval(['ys(' int2str(i) ') = ' varname ';']);
end
end

else
ys=zeros(length(lgy_),1);
for i = 1:length(lgy_)
    ys(i) = eval(lgy_(i,:));
end
check = 0;
end