% DYNARE CODE
%
%

%----------------------------------------------------------------
% 0. Let's start
%----------------------------------------------------------------

close all;

%----------------------------------------------------------------
% 1. Defining variables
%----------------------------------------------------------------

var rpn c ph ct cn ch cf cstar cnstar ctstar chstar cfstar pt pf pn phstar pfstar ptstar p pstar pnstar yn ynstar yt ytstar tot rer news3 news2 news1 vt;
varexo et;

parameters at ah ahstar phi ytbar rho rhon rhot sigmat rho_v;

%----------------------------------------------------------------
% 2. Calibration
%----------------------------------------------------------------

at = 0.2;
ah  = 0.85;
ahstar = 0.15;

phi = 1.5;
ytbar = 0.5;

rho = 1.5;

rhon = 0.99;
rhot = 0.99;
sigma = 1.5;
sigmat = 0.0075;

rho_v=0.9;

%----------------------------------------------------------------
% 3. Model
%----------------------------------------------------------------

model; 

ch = ah*(ph/pt)^(1-phi)*yt;
cf = (1-ah)*(pf/pt)^(-phi)*(ph/pt)*yt;
cn = (1-at)*(pn/p)^(-rho)*(yt(-1)*ph+yn*pn)/p;
chstar = ahstar*(phstar/ptstar)^(-phi)*(ytstar*pf/ptstar);
cfstar = (1-ahstar)*(pfstar/ptstar)^(1-phi)*ytstar;
cnstar = (1-at)*(pnstar/pstar)^(-rho)*(ytstar*pf+ynstar*pnstar)/pstar;

c=(at^(1/rho)*ct^(1-1/rho)+(1-at)^(1/rho)*cn^(1-1/rho))^(rho/(rho-1));
ct=(ah^(1/phi)*ch^(1-1/phi)+(1-ah)^(1/phi)*cf^(1-1/phi))^(phi/(phi-1));
cstar=(at^(1/rho)*ctstar^(1-1/rho)+(1-at)^(1/rho)*cnstar^(1-1/rho))^(rho/(rho-1));
ctstar=((1-ah)^(1/phi)*chstar^(1-1/phi)+(ah)^(1/phi)*cfstar^(1-1/phi))^(phi/(phi-1));

p = (at*pt^(1-rho)+(1-at)*pn^(1-rho))^(1/(1-rho));
pt = (ah*ph^(1-phi)+(1-ah)*pf^(1-phi))^(1/(1-phi));
pstar = (at*ptstar^(1-rho)+(1-at)*pnstar^(1-rho))^(1/(1-rho));
ptstar = (ah*phstar^(1-phi)+(1-ah)*pfstar^(1-phi))^(1/(1-phi));

ph =1;

yn= cn;
ynstar= cnstar;

ph = phstar;
pf = pfstar;

tot = pf/phstar;
rer = pstar/p;
rpn = pn/pt;

% Complete markets condition

rer=(c/cstar)^(rho);

% Endowments 

yt = (1-rhot)*ytbar + rhot*yt(-1) +news3(-1);
ytstar = (1-rhot)*ytbar + rhot*ytstar(-1);
yn = (1-rhon)*ytbar + rhon*yn(-1);
ynstar = (1-rhon)*ytbar + rhon*ynstar(-1);

% News shocks

news3 = news2(-1);
news2 = news1(-1);
news1 = vt(-1);
vt = rho_v*vt(-1)+et;


end;

%----------------------------------------------------------------
% 4. Initial Values 
%----------------------------------------------------------------

initval;

rpn =1;
c =0.783045; 
ct = 0.335853;
cn = 0.5; 
ch =0.285475;
cf = 0.0503779;
cstar =0.783045;
cnstar =0.5;
ctstar = 0.335853; 
chstar =0.0503779;
cfstar =0.285475;
p =1.66299;
ph =1;
pf =1;
pt = 1;
pn = 1.93268;
pstar =1.66299;
phstar =1; 
pfstar =1;
ptstar = 1;
pnstar = 1.93268;

yt =0.5;
ytstar =0.5;
yn = 0.5;
ynstar = 0.5;

tot =1;
rer =1;
vt=0;
%epsilont = 0;
%epsilontstar = 0;
%epsilonn = 0;
%epsilonnstar = 0;

end;

%----------------------------------------------------------------
% 5. Shocks 
%----------------------------------------------------------------

shocks;
%var epsilont = 1;
%var epsilontstar = 1;
%var epsilonn = 1;
%var epsilonnstar = 1;
var et = 1;
end;

%----------------------------------------------------------------
% 6. Solve the model 
%----------------------------------------------------------------

steady;

stoch_simul(hp_filter = 1600, order = 1);

%----------------------------------------------------------------
% 7. Some Results
%----------------------------------------------------------------

%statistic1 = 100*sqrt(diag(oo_.var(1:6,1:6)))./oo_.mean(1:6);
%dyntable('Relative standard deviations in %',strvcat('VARIABLE','REL. S.D.'),M_.endo_names(1:6,:),statistic1,10,8,4);