function [ys,check1]=quest1_steadystate(junk, exs)

global M_ options_

for j=1:size(M_.param_names,1)
  eval([deblank(M_.param_names(j,:)),' = M_.params(j);']) 
  assignin('base',deblank(M_.param_names(j,:)),M_.params(j)); 
end

check1=0;
if isempty(DDYN), DDYN=1; end,
if isempty(SIMFLAG), SIMFLAG=0; end,

 

E_EX_RW = E_EX_INOMW - GPW0;
E_EX_R = 1/BETAE-1;
if SIMFLAG==0,
  E_EX_RW = E_EX_R;
  E_EX_INOMW = E_EX_RW + GPW0;
end
  
GTFP0      =   (ALPHAE+ALPHAGE-1)/ALPHAE*GY0-(2-ALPHAE-ALPHAGE)/ALPHAE*GPCPI0; 
GYW0 = GY0;
% DELTAKE    =   (DELTAE+GPOP0)/(1+GPOP0);
% DELTAKGE    =   (DELTAGE+GPOP0)/(1+GPOP0);
DELTAKE    =   (DELTAE+GPOP0);
DELTAKGE    =   (DELTAGE+GPOP0);
M_.params(strmatch('E_EX_RW',M_.param_names,'exact')) = E_EX_RW;
M_.params(strmatch('GYW0',M_.param_names,'exact')) = GYW0;
M_.params(strmatch('GTFP0',M_.param_names,'exact')) = GTFP0;
% end;
%SE=0.5*(SESUM+SEDIFF);
SWE=SE;
E_ETA =(1-TAUE);    
% E_DEPRAT = DEPRAT0;
% E_DEPRAT1 = DEPRAT0;
% E_DEPRAT2 = DEPRAT0;
E_GC = GY0;        
E_GCL = GY0+GPOP0;        
E_GCLC = GY0;        
E_GCNLC = GY0;          
E_GE = GP0 - GPW0;       
E_GER =0;       
E_GEX=GY0;        
E_GEXL = GY0+GPOP0+DGEX;        
E_GG =GY0;        
E_GGL =GY0+GPOP0;        
E_GI =GY0+GPCPI0;   
E_GIG = E_GI;
E_GIL =GY0+GPCPI0+GPOP0;        
E_GIM=GY0;        
E_GIML = GY0+GPOP0+DGIM;        
E_GK =E_GI;        
E_GKG =E_GI;        
E_GL  =0;      
E_GPOPA = GPOP0;
E_GTFP = (ALPHAE+ALPHAGE-1)/ALPHAE*GY0-(2-ALPHAE-ALPHAGE)/ALPHAE*GPCPI0;      
E_GTFPUCAP = ALPHAE*E_GTFP;
E_GUC = 0;       
E_GUCAP = 0;     
E_GWRY = 0;      
E_GY=GY0;       
E_GYL=GY0+GPOP0;       
E_GTAX = GY0+GP0;
E_GTR = GY0;
E_GYPOT=GY0;      
E_GYW=GYW0;        
E_INOM  = E_EX_R+GP0;
E_INOMW  =E_EX_INOMW;
E_LOL = 0;       
% E_NPART = NPART0;       
E_BWRY = (E_EX_RW - E_EX_R)/RPREME;
TBYN = (-E_INOM+GPOP0+GP0+GY0)*E_BWRY;
E_TBYN=TBYN;
E_LYWY = LYWY0;
M_.params(strmatch('LYWY0',M_.param_names,'exact')) = LYWY0;
YWY = exp(LYWY0);

% E_LIK = log(DELTAKE+GY0+GPCPI0)-(GY0+GPCPI0);   
% E_LIKG = log(DELTAKGE+GY0+GPCPI0)-(GY0+GPCPI0);   
% IK = (DELTAKE+GY0+GPCPI0)/exp(GY0+GPCPI0);   
E_LIK = log(DELTAKE+GY0+GPCPI0);   
E_LIKG = log(DELTAKGE+GY0+GPCPI0);   
IK = (DELTAKE+GY0+GPCPI0);   
if DDYN==0,
    E_LYKPPI  = E_LIK-log(ISN);     
    YKN = IK/ISN; 
else    
    YKN = (E_EX_R+RPREMK+DELTAE)/(1-TP)/(1-TAUE)/(1-ALPHAE);
    E_LYKPPI=log(YKN);
    ISN=IK/YKN;
    M_.params(strmatch('ISN',M_.param_names,'exact')) = ISN;
    assignin('base','ISN',ISN);
end

if SIMFLAG == 1,
  UCAP0 = (1-TAUE)*(1-ALPHAE)*exp(E_LYKPPI)/A1E;
  set_param_value('UCAP0',UCAP0);
  % assignin('base','UCAP0',UCAP0);
  % M_.params(strmatch('UCAP0',M_.param_names,'exact')) = UCAP0;
  set_param_value('SIMFLAG',-1);
  if abs(E_EX_RW - E_EX_R)<1.e-8,
    pname={'L0'};
    eqind=10;
%     L0 = fzero(@(x) quest1_steadystate_util(x, pname, eqind),L0);
    L0 = fzero(@(x) quest1_steadystate_util(x, pname),L0);
    LER0=0;
  else
    while abs(quest1_steadystate_util)>options_.dynatol,
    pname={'L0'};
    eqind=10;
    L0 = fzero(@(x) quest1_steadystate_util(x, pname, eqind),L0);
    pname={'LER0'};
    eqind=30;
    LER0 = fzero(@(x) quest1_steadystate_util(x, pname, eqind),LER0);
    end
%     if quest1_steadystate_util>options_.dynatol,
%     opt=optimset('fminsearch');
% %     opt.Display='iter';
%     opt.TolX=1.e-4;
%     opt.TolFun=1.e-4;
%     opt.MaxFunEvals=1000;
%     pname={'L0','LER0'};
%     xx = fminsearch(@(x) quest1_steadystate_util(x, pname),[L0; LER0],opt);
% %     xx = dynare_solve('quest1_steadystate_util',[L0 LER0],0,{'L0','LER0'},[10 30]);
% 
%     L0=xx(1);
%     LER0=xx(2);
%     end
    set_param_value('LER0',LER0);
  end
  set_param_value('L0',L0);
  set_param_value('SIMFLAG',1);


end
E_L =L0;
E_LL =log(L0);
E_LL0 =log(L0);    

if SIMFLAG==0 | abs(E_EX_RW - E_EX_R)<1.e-8,
  LER0=0;
end
E_ER=exp(LER0);
E_LER = LER0;

PMP = (E_ER^ALPHAX);    
E_LPMP = log(PMP);    
%PCP = PMP^(1-SE);
PCP = (SE+(1-SE)*PMP^(1-SIGIME))^(1/(1-SIGIME));
E_LPCP= log(PCP);
%E_LPCP = (1-SE)*E_LPMP - SE*(1-SE)*E_LPMP^2*(SIGIME-1)/2;
%PCP=exp(E_LPCP);
E_LPCPM0 = E_LPCP-E_LPMP;
PXP = 1; %PCP^(S0*(1-SXDE));
E_LPXP = 0; %(1-SXDE)*S0*E_LPCP; 
% E_LIMYN = log(1-SE);
% E_LEXYN = log(1-SWE);   

IMYN = (1-SE)*exp(E_LPCP-E_LPMP)^(SIGIME-1)*(1-TBYN);
E_LIMYN = log(IMYN);
EXYN = (1-SE)*exp(E_LER*(ALPHAX*SE))^SIGEXE*exp(E_LYWY)^ALPHAX;       
E_LEXYN = log(EXYN);   


E_LIGSN =  log(IGSN);
E_ZEPS_IG = 0;

E_LGSN = log(GSN);       
E_LISN = log(ISN);      
CSN = 1-(IGSN+GSN+ISN+TBYN);
M_.params(strmatch('CSN',M_.param_names,'exact')) = CSN;
E_LCSN = log(CSN);      
E_LBGYN=log(BGTAR);
YWR = 1/((1-TAUE)*ALPHAE/L0*(1+E_LOL))  ;
E_LYWR = log(YWR) ;
E_LWPTU = 0;
E_TW = 0.2;
%E_TW=TW0*(E_WPTU)^TW1;     

% TW0=E_TW/(E_WPTU)^TW1;    

M_.params(strmatch('TW0',M_.param_names,'exact')) = TW0;
E_TRW = TRSN;
E_WS = exp(E_LL-E_LYWR);
E_TAXYN  = (E_INOM-GP0-GY0-GPOP0)*BGTAR+IGSN+GSN+TRSN*E_WS ...
     -(E_TW+SSC)*E_L/YWR -TP*(1-E_L/YWR) -TVAT*CSN;

CLCSN = ((1-E_TW-SSC)*E_L/YWR + TRSN*E_WS - E_TAXYN)/(1+TVAT);      
E_LCLCSN = log(CLCSN);      
CNLCSN = (CSN-CLCSN*SLC*SLCFLAG)/(1-SLC*SLCFLAG);
E_LCNLCSN = log(CNLCSN);
UCTERM = ((1-SLCFLAG*SLC)*(CNLCSN*(1-HABE))^(-SIGC)+SLCFLAG*SLC*CLCSN^(-SIGC))/((1-SLCFLAG*SLC)*(CNLCSN*(1-HABE))^(1-SIGC)+SLCFLAG*SLC*CLCSN^(1-SIGC));
A = (THETAE-1)/THETAE*(1-E_TW-SSC)/(1+TVAT)*E_WS*UCTERM/L0;
if SIMFLAG==0,
  OMEGE = A/(KAPPAE*(L0*(1-HABLE))^(KAPPAE-1)+A*(L0*(1-HABLE))^KAPPAE);
end
UCYN = (CNLCSN*(1-HABE))^(-SIGC)*(1-OMEGE*(L0*(1-HABLE))^KAPPAE)^(1-SIGC);     
E_LUCYN = log(UCYN);
UCLCYN = (CLCSN)^(-SIGC)*(1-OMEGE*(L0*(1-HABLE))^KAPPAE)^(1-SIGC);     
E_LUCLCYN = log(UCLCYN);
E_VL = (CNLCSN*(1-HABE))^(1-SIGC)*(1-OMEGE*(L0*(1-HABLE))^KAPPAE)^(-SIGC)*KAPPAE*OMEGE*(L0*(1-HABLE))^(KAPPAE-1);
E_VLLC = (CLCSN)^(1-SIGC)*(1-OMEGE*(L0*(1-HABLE))^KAPPAE)^(-SIGC)*KAPPAE*OMEGE*(L0*(1-HABLE))^(KAPPAE-1);

E_LYGAP = 0 ;   
E_LYGAP1 = 0 ;   
E_LCY = E_LCSN-E_LPCP;
E_LGY = E_LGSN-E_LPCP;
E_LWS = log(E_WS);

E_CLCSN = exp(E_LCLCSN);
E_TRYN = E_TRW*exp(E_LL-E_LYWR) ;
E_LTRYN = log(E_TRYN);
E_TRTAXYN = E_TRW*exp(E_LL-E_LYWR) - E_TAXYN;
E_BGYN = exp(E_LBGYN);
E_GSN = exp(E_LGSN);
E_WSW = (1-E_TW-SSC)*E_WS;
E_MRY = (1+E_INOM)^(-ZETE);
E_PHI=GP0;        
E_PHIC=GP0;       
E_PHIPI=GPCPI0;       
E_PHIM=GP0;       
E_PHIML=GP0+DGPM;       
E_PHIW=GPW0;       
E_PHIX=GP0;       
E_PHIXL=GP0+DGPX;       
E_Q = 1 ;        
E_R = E_EX_R;
E_TI = 0;
E_UCAP  =UCAP0;     
E_UCAP0  =UCAP0;     
E_WPHI  = GP0+GY0;
E_WRPHI  = GY0;
E_ZEPS_C = 0;     
E_ZEPS_CLC = 0;     
E_ZEPS_EQ = 0;     
E_ZEPS_ETA = 0;   
E_ZEPS_ETAM = 0;    
E_ZEPS_ETAX = 0;   
E_ZEPS_EX = 0;    
E_ZEPS_G = 0;     
E_ZEPS_I = 0;     
E_ZEPS_IM = 0;     
E_ZEPS_L = 0;     
E_ZEPS_M = 0;     
E_ZEPS_PC = 0;    
E_ZEPS_PPI = 0;
E_ZEPS_POP = 0;    
E_ZEPS_PX = 0;    
E_ZEPS_RPREME = 0;
E_ZEPS_RPREMK = 0;
E_ZEPS_SLC = 0;   
E_ZEPS_TFP = 0;   
E_ZEPS_TR = 0;   
E_ZEPS_W = 0;     
E_ZEPS_Y = 0;     
E_ZEPS_YW = 0;     
E_DBGYN = 0;     
E_ZEPS_YGAP = 0;     
E_ZPHIT = GP0;     

ys=[
E_BGYN       
E_BWRY       
E_CLCSN      
E_DBGYN      
E_LER        
E_ETA        
E_GC         
E_GCL        
E_GCLC       
E_GCNLC      
E_GE         
E_GEX        
E_GEXL       
E_GG         
E_GGL        
E_GI         
E_GIG        
E_GIL        
E_GIM        
E_GIML       
E_GK         
E_GKG        
E_GL         
E_GSN        
E_GTAX       
E_GTFP       
E_GTFPUCAP   
E_GTR        
E_GUC        
E_GUCAP      
E_GWRY       
E_GY         
E_GYL        
E_GYPOT      
E_GYW        
E_INOM       
E_INOMW      
E_LL         
E_LL0        
E_LBGYN      
E_LCSN       
E_LCLCSN     
E_LCNLCSN    
E_LEXYN      
E_LGSN       
E_LIGSN      
E_LIMYN      
E_LIK        
E_LIKG       
E_LISN       
E_LOL
E_LPCP       
E_LPMP       
E_LPXP       
E_LTRYN      
E_LUCYN      
E_LUCLCYN    
E_LYGAP      
E_LYKPPI     
E_LYWR       
E_LYWY       
E_MRY        
E_PHI        
E_PHIC       
E_PHIPI      
E_PHIM       
E_PHIML      
E_PHIW       
E_PHIX       
E_PHIXL      
E_Q          
E_R          
E_TAXYN      
E_TBYN       
E_TRTAXYN    
E_TRW        
E_TRYN       
E_TW         
E_UCAP       
E_UCAP0      
E_VL         
E_VLLC       
E_WPHI       
E_WRPHI      
E_WS         
E_WSW        
E_ZEPS_C     
E_ZEPS_ETA   
E_ZEPS_ETAM  
E_ZEPS_ETAX  
E_ZEPS_EX    
E_ZEPS_G     
E_ZEPS_IG    
E_ZEPS_L     
E_ZEPS_M     
E_ZEPS_PPI   
E_ZEPS_RPREME
E_ZEPS_RPREMK
E_ZEPS_TR    
E_ZEPS_W     
E_ZPHIT      
E_LCY
E_LGY
E_LWS
];

if SIMFLAG==0,

A1E = (1-TAUE)*(1-ALPHAE)*exp(E_LYKPPI)/UCAP0;
assignin('base','A1E',A1E);
M_.params(strmatch('A1E',M_.param_names,'exact')) = A1E;

M_.params(strmatch('OMEGE',M_.param_names,'exact')) = OMEGE;
assignin('base','OMEGE',OMEGE);
if DDYN==0,
    RPREMK = (1-TAUE)*(1-ALPHAE)*(1-TP)*exp(E_LYKPPI) - DELTAE-E_EX_R;
    M_.params(strmatch('RPREMK',M_.param_names,'exact')) = RPREMK;
    assignin('base','RPREMK',RPREMK);
end
end

%PAR1 = 1/((1-TP)*(1/exp(E_LISN)-(1/exp(E_LISN))*(1/exp(E_LYWR))- SHCRED*E_INOM*exp(E_LIK)));
%PAR1 = exp(E_LIK)/((1-TP)*(exp(E_LYKPPI) -exp(E_LYKPPI)/exp(E_LYWR)-SHCRED*(E_INOM+RPREMK)));
%M_.params(strmatch('PAR1',M_.param_names,'exact')) = PAR1;

% exs = M_.exo_nbr;
% ex_ = zeros(1,exs);
% 
% [fmax, iff] = max(abs(quest1_static(ys, ex_, M_.params)));
% 
% if fmax>1.e-8; %options_.dynatol;
%     disp(['The analytic steady state is not correct! ',num2str(fmax)])
%     disp(['Equation ',int2str(iff),', discrepancy ',num2str(fmax)])
%     check1=1;
% end

% if SIMFLAG==1,
% exs = M_.exo_nbr;
% ex_ = zeros(1,exs);
% disp(num2str(max(abs( quest1_static(ys, ex_, M_.params))))),
% end