var yt nt Gammat ct qt vt wt mt ut st Et Ut e_at e_ht e_etat e_bt

zt zrt 

e_ft Jt St EWt EFt;


varexo innovation_at innovation_ht innovation_etat innovation_bt

innovation_ft;


parameters rho beta b eta a sigma chi n q v c y w u E U s h etaz omegab  

omegaf omegah Gzr gzr z zr wprimezr J EW EF S

F rhoa rhoh rhoeta rhob rhof sigmaa sigmah sigmaeta sigmab sigmaf;



beta = 0.98;

sigma = 0.6;

n = 0.9;

q = 0.9;

s = 0.25;

etaz = 1;

Gzr = 0.02;

omegab = 0.9;

omegaf = 0.2;

omegah = 0.2;

a = 1;

rhoa = 0.5;  
rhoh = 0.5;
rhoeta = 0.5;
rhob = 0.5;
rhof = 0.5;
sigmaa = 0.01;  
sigmah = 0.01;
sigmaeta = 0.01;
sigmab = 0.01;
sigmaf = 0.01;






% MODEL EQUATIONS

%model;
model (linear);


%yt=e_at*nt;
%yt=e_at+nt;
yt=e_at+nt+zt;

zt=etaz*zrt;

%chi/qt+e_ht=e_at-wt+rho*beta*Gammat(+1)*(chi/qt(+1)+e_ht(+1));
%-chi/q*qt+h*e_ht=e_at-w*wt+rho*beta*chi/q*(Gammat(+1)-qt(+1))+rho*beta*h*(Gammat(+1)+e_ht(+1));
-chi/q*qt+h*e_ht=(z-w-gzr*zr-gzr*zr*(1+F)-rho*beta*gzr*zr*h)*zrt+(1-Gzr)*z*(e_at+zt)-(1-Gzr)*w*wt-Gzr*F*e_ft+(1-Gzr)*rho*beta*chi/q*(Gammat(+1)-qt(+1))+(1-Gzr)*rho*beta*h*(Gammat(+1)+e_ht(+1));

zr*(e_at+zrt)-wprimezr*zr*zrt+rho*beta*h*(Gammat(+1)+e_ht(+1))+rho*beta*c/q*(Gammat(+1)-qt(+1))+F*e_ft=0;
%zrt=0; 

%ut=1-rho*nt(-1);
ut=-(n/u)*rho*nt(-1);

%nt=rho*nt(-1)+qt*vt;
%n*nt=rho*n*nt(-1)+q*v*(qt+vt);
n*nt=-gzr*zr*(rho*n+q*v)*zrt+(1-Gzr)*(rho*n*nt(-1)+q*v*(qt+vt));

%qt=mt/vt;
qt=mt-vt;

%st=mt/ut;
st=mt-ut;

%mt=sigmam*ut^sigma*vt^(1-sigma);
mt=sigma*ut+(1-sigma)*vt;

%yt=ct+chi*vt+e_ht*vt*qt;
%y*yt=c*ct+chi*v*vt+h*q*v*(qt+vt+e_ht);
y*yt=c*ct+chi*v*vt+h*q*v*(qt+vt+e_ht)+n*F*Gzr/(1-Gzr)*(nt+e_ft)+n*F*gzr*zr/(1-Gzr)^2*zrt;


%Et=wt+beta*Gammat(+1)*(rho*Et(+1)+(1-rho)*Ut(+1));
%E*Et=w*wt+rho*beta*E*(Gammat(+1)+Et(+1))+(1-rho)*beta*U*(Gammat(+1)+Ut(+1));
E*Et=w*wt+rho*beta*( (1-Gzr)*E*(Gammat(+1)+Et(+1))+Gzr*U*(Ut(+1)+Gammat(+1))-zr*gzr*(E-U)*zrt(+1))+(1-rho)*beta*U*(Gammat(+1)+Ut(+1));

%Ut=e_bt+beta*Gammat(+1)*(st(+1)*Et(+1)+(1-st(+1))*Ut(+1));
%U*Ut=b*e_bt+s*beta*E*(Gammat(+1)+Et(+1)+st(+1))+beta*(1-s)*U*(Gammat(11)+Ut(+1))-s*beta*U*st(+1);
U*Ut=b*e_bt+beta*(s*E*(1-Gzr)*(Gammat(+1)+st(+1)+Et(+1))+s*Gzr*U*(Gammat(+1)+st(+1)+Ut(+1))+(1-s)*U*(Gammat(+1)+Ut(+1))-s*U*st(+1)-s*gzr*zr*(E-U)*zrt(+1));

J*Jt=-S/eta*e_etat+(1-eta)/eta*S*St-F*e_ft;

S*St=E*Et-U*Ut;

EF*EFt=-gzr*zr*J*zrt-gzr*zr*F*zrt+(1-Gzr)*J*Jt-Gzr*F*e_ft;

EW*EWt=-gzr*zr*S*zrt-gzr*zr*U*zrt+(1-Gzr)*S*St-Gzr*U*Ut;

%wt=e_etat*(e_at+rho*beta*Gammat(+1)*(chi/qt(+1)+e_ht(+1)))+(1-e_etat)*(e_bt-beta*Gammat(+1)*(rho-st(+1))*(Et(+1)-Ut(+1)));
%w*wt=eta*(1+rho*beta*(chi/q+h))*e_etat + eta*e_at+eta*rho*beta*chi/q*(Gammat(+1)-qt(+1))+eta*rho*beta*h*(Gammat(+1)+e_ht(+1))-eta*(b-(rho-s)*(E-U))*e_etat+(1-eta)*b*e_bt-(1-eta)*beta*(rho-s)*E*(Gammat(+1)+Et(+1))+(1-eta)*beta*(rho-s)*U*(Gammat(+1)+Ut(+1))+(1-eta)*beta*s*(E-U)*st(+1);
w*wt=eta*z*(e_etat+e_at+zt)+eta*F*(e_etat+e_ft)-eta*b*e_etat+(1-eta)*b*e_bt+rho*beta*EF*(Gammat(+1)+EFt(+1))+s*beta*EW*st(+1)+(s-rho)*beta*EW*(Gammat(+1)+EWt(+1));

%Gammat=1;
Gammat=0; 

%e_at=(1-rhoa)*a+rhoa*e_at(-1)+innovation_at;
e_at=rhoa*e_at(-1)+innovation_at;

%e_ht=(1-rhoh)*h+rhoh*e_ht(-1)+innovation_ht;
e_ht=rhoh*e_ht(-1)+innovation_ht;

%e_etat=(1-rhoeta)*eta+rhoeta*e_etat(-1)+innovation_etat;
e_etat=rhoeta*e_etat(-1)+innovation_etat;

%e_bt=(1-rhob)*b+rhob*e_bt(-1)+innovation_bt;
e_bt=rhob*e_bt(-1)+innovation_bt;

e_ft=rhof*e_ft(-1)+innovation_ft;


end;


shocks;
var innovation_at  = sigmaa^2 ;
var innovation_ht = sigmah^2;
var innovation_etat = sigmaeta^2;
var innovation_bt = sigmab^2;
end;

%initval;
%yt=y; 
%nt=n; 
%Gammat=1; 
%ct=c; 
%qt=q; 
%vt=v; 
%wt=w; 
%mt=m; 
%ut=u; 
%st=s; 
%Et=E; 
%Ut=U; 
%e_at=a; 
%e_ht=h; 
%e_etat=eta; 
%e_bt=b;
%end;


model_diagnostics(M_,options_,oo_)


check;
resid(1);
steady;


stoch_simul(irf=40,periods=500, order=1, nograph)yt nt wt ct;