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

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



%%%% Model equations to be entered here

% Estimated Parameters Initialization
beta1=1.5;
h=0.5;
gamma_l=2;
kappa=4;
theta_w=0.6;
chi_w=0.5;
alpha=0.45;
theta_p=0.6;
chi=0.5;
gamma_y=2;
gamma_pi=2;
Lambda_mu1=0.62;
yta=10;
gamma_w=0.5;
rho_d=0.5;
rho_phi=0.5;

eps=1.5;

c=0;
lam=0;
ud=0;
uphi=0;
z=0;
R=0;
m=0;
pi=0;
pi_s=0;
v=0;
r=0;
q=0;
i=0;
w=0;
w_s=0;
f=0;
ld=0;
omega=0;
g_1=0;
g_2=0;
yd=0;
bigpsi=0;
mu=0;
k=0;
Delta_p=0;
Delta_w=0;
A=0;
l=0;

laml=0;
pil=0;
ql=0;
rl=0;
w_sl=0;
fl=0;
g_1l=0;
g_2l=0;
pi_sl=0;
xi_lam=0;
xi_pi=0;
xi_q=0;
xi_r=0;
xi_ws=0;
xi_f=0;
xi_g1=0;
xi_g2=0;
xi_pis=0;
xi_i=0;

eps_d=0;
eps_phi=0;
eps_mu=0;
eps_a=0;
eps_w=0;

laml=0;
pil=0;
rl=0;
ql=0;
w_sl=0;
fl=0;
g_1l=0;
g_2l=0;
pi_sl=0;
il=0;

dy=0;
di=0;
dpi=0;
dR=0;
domega=0;

% Fixed Prior
beta = 1.0/(beta1/100.0+1.0);
Lambda_mu=Lambda_mu1 / 100;
muss=exp(Lambda_mu);
piss=1.072;
lss=0.25;
zss=exp(0.0238);

% Steady State

delta = 1.0-(1.0-0.0678)*zss*muss; %depreciation rate
r = zss*muss/beta-(1.0-delta); %gamma1
R = piss*zss/beta; %s.s. nominal interest rate
Ass = exp(log(zss)*(1-alpha)-log(muss)*alpha); %s.s. growth rate of neutral technology shock
gamma_1=r;

pi_w = ((1.0-theta_w*piss^((1.0-chi_w)*(yta-1.0))*zss^(yta-1.0))/(1.0-theta_w))^(1.0/(1.0-yta)); %s.s. ratio of optimal wage to aggregate wage

pi_s = ((1.0-theta_p*piss^((1.0-chi)*(eps-1.0)))/(1.0-theta_p))^(1.0/(1.0-eps)); %s.s. ratio of optimal price to aggregate price
Delta_p = (1.0-theta_p)*pi_s^(-eps)/(1.0-theta_p*piss^((1.0-chi)*eps)); %s.s. relative price dispersion
bigpsi = ((eps-1.0)/eps)*pi_s*(1.0-beta*theta_p*piss^((1.0-chi)*eps))/(1.0-beta*theta_p*piss^((1.0-chi)*(eps-1.0))); %s.s. real marginal cost
w = (1.0-alpha)*(bigpsi*(alpha/r)^alpha)^(1.0/(1.0-alpha)); %s.s. real aggregate wage
Delta_w = (1.0-theta_w)*pi_w^(-yta)/(1.0-theta_w*zss^yta*piss^((1.0-chi_w)*yta)); %s.s. relative wage dispersion

ld = lss/Delta_w; %s.s. effective labor service
k = ld*(alpha/(1.0-alpha))*(w/r)*zss*muss; %s.s. capital-labor ratio
%F = Delta_p/5.873-(Ass/zss)*(k/ld)^(alpha-1.0);

eps0=1.5;
x = fsolve(@(x)eps_solver(x,beta,alpha,theta_p,chi,Lambda_mu1),eps0,optimset('Display','off'));
eps=x

yd=k/5.873;
w_s = w*pi_w; %s.s. real optimal wage
%//Omega = k_Ld; %auxiliary variable for s.s. capital-labor ratio
%//k_tilde = Omega*Ld; %s.s. capital stock
%//i_tilde = i_k_tilde*k_tilde; %s.s. investment
%//y_d_tilde = k_tilde/k_y_tilde; %s.s. output
c=yd-i;
%//c_tilde = y_d_tilde*(1.0-i_tilde/y_d_tilde); %s.s. consumption
lam = (1.0-h/zss)^(-1.0)*c^(-1.0); %s.s. shadow price
nu_m = lam*(1.0-beta/(zss*piss)); %s.s. ratio related to money demand 
%//LHS = (1.0-beta*theta_w*zss^(yta*(1.0+gamma_l))*piss^(yta*(1.0-chi_w)*(1.0+gamma_l)))/(1.0-beta*theta_w*zss^(yta-1.0)*piss^((1.0-chi_w)*(yta-1.0)));
%//psi = (LHS*((yta-1.0)/yta)*wss_s*lamss)/(piss_w^(-yta*gamma_l)*ldss^gamma_l); %weight on labor service
i=k*0.0678;

omega=zss*piss;

rss=r
Rss=R
piss_s=pi_s
Deltass_p=Delta_p
bigpsiss=bigpsi
wss=w
kss=k
piss_w=pi_w
Deltass_w=Delta_w
ldss=ld
wss_s=w_s
css=c
lamss=lam
nu_mss=nu_m
omegass=zss*piss


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);
for i = 1:NumberOfEndogenousVariables
  varname = deblank(M_.endo_names(i,:));
  eval(['ys(' int2str(i) ') = ' varname ';']);
end
else
ys=zeros(length(lgy_),1);
for i = 1:length(lgy_)
    ys(i) = eval(lgy_(i,:));
end
check = 0;
end