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

options= optimset('Display','Final', 'TolX',1e-10, 'TolFun', 1e-10 ,'MaxFunEvals', 200000,'MaxIter', 200000);

Pi 	= 1;
Pihf 	= 1;
Pihfst 	= 1;
Pii 	= 1;
Piist 	= 1;
Pif 	= 1; 
E 	= 1;
Dh 	= 1;
Di 	= 1;
Q 	= betaa;
R 	= 1/betaa;
Af 	= 0;
Ai 	= 0;

Ast 	= 0;
rst =(1/betast) ;
Pist = 1; 
Piffst = 1; 
Pint_st = 1;
Piff_st = 1;

Pintst = 1;
Dff = 1;
Dnt = 1;
Yst = ((1-alphast)*(1/mst))^(1/(1+psii));
Ynt = (alphast*(1/mst))^(1/(1+psii));
Cst = Ynt^(alphast)*Yst^(1-alphast);
F1st = Yst^(psii+2)/((1-alphast)*(1 - betast*thetaa));
F2st = Yst/(1- betast*thetaa);
V1st = Ynt^(psii+2)/(alphast*(1 - betast*thetaa));
V2st = Ynt/(1- betast*thetaa);
Gpst = (Ynt/(alphast*Cst))^(1/(alphast-1));

thetstar=Gpst;

x0=[1.1615;0.8547;0.3709;0.7385;0.5137;0.6514;0.9036;0.1329];

[lam] = fsolve(@steadystate, x0, options,Cstarr,thetstar,psii,m,mi,alphaa,a,deltaa,chii,rhoo,Rw);

C=lam(1);
S=lam(2);
Thet=lam(3);
Yh=lam(4);
Yi=lam(5);
Ni=lam(6);
Nf=lam(7);
G=lam(8);

MCh	= 1/m;
F1	= Yh*MCh/(1-thetaf*betaa);
F2	= Yh/(1-thetaf*betaa);
MCi	= 1/mi;
V1	= Yi*MCi/(1-thetai*betaa);
V2	= Yi/(1-thetai*betaa);

U = log(C)-((Dh*Yh/(G^chii)) + (Di*Yi/((deltaa*G)^chii)))^(1+psii)/(1+psii);
W = U/(1-betaa);

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