% Endogenous variables 
var Cih Cif Csh Csf kihh kihf kiff kifh ksh ksf Bih Bif Bsh Bsf q1 q2 Rh Rf Wh Wf Rkh Rkf Ah Af Zh Zf L Yh Yf Jh Jf BETAi BETAs;

% Exogenous variables
varexo nu1 nu2 nu3 nu4;

parameters Epsilon KAPPA ETA N RHO SIGMA OMEGA DELTA r ZETAi ZETAs; 

%Parameters values
Epsilon = 0.36;     %Capital share in final goods production
ETA = 0.022;        %Discount function
N = 0.5;            %Proportion of investors
RHO = 0.9;         % Productivity shock persistence for investors in home country producing final good
SIGMA = 2;          % coefficient of relative risk aversion
ZETAi = 0.975;      %Investors discount function
ZETAs = 0.9825;     %Savers discountfunction
OMEGA = 0.1;        %Capital share in home production
DELTA = 0.0226;     %capital depreciation rate
r = 0.04; 	        %world interest rate


%Steady state variables 
BETAi_ss = 1/(1+r);
BETAs_ss = 1/(1+r);
Rh_ss = (1/BETAs_ss);
Rf_ss = (1/BETAs_ss);
Rkh_ss =(1/BETAs_ss);
Rkf_ss = (1/BETAs_ss);
kihh_ss = 0.5*(Rh_ss/Epsilon)^(1/(Epsilon-1));
kihf_ss = 0.5*(Rh_ss/Epsilon)^(1/(Epsilon-1));
kiff_ss = 0.5*(Rh_ss/Epsilon)^(1/(Epsilon-1));
kifh_ss = 0.5*(Rh_ss/Epsilon)^(1/(Epsilon-1));
ksh_ss = (1/(1-N))*(1-(Rh_ss/Epsilon)^(1/(Epsilon-1)));
ksf_ss = (1/(1-N))*(1-(Rh_ss/Epsilon)^(1/(Epsilon-1)));
Yh_ss =(Rh_ss/Epsilon)^(Epsilon/(Epsilon-1));
Yf_ss =(Rh_ss/Epsilon)^(Epsilon/(Epsilon-1));
Jh_ss = (ksh_ss )^(OMEGA);
Jf_ss = (ksf_ss )^(OMEGA);
Cih_ss = Yh_ss ;
Cif_ss =Yf_ss; 
Csh_ss = Jh_ss;
Csf_ss = Jf_ss; 
Wh_ss = (1-Epsilon)*(Rh_ss/Epsilon)^(1/(Epsilon-1));
Wf_ss =  (1-Epsilon)*(Rh_ss/Epsilon)^(1/(Epsilon-1));
Bih_ss = 0;
Bif_ss = 0;
Bsh_ss  =0;
Bsf_ss = 0;
q1_ss =1;
q2_ss = 1;
Ah_ss = 1;
Af_ss = 1;
Zh_ss = 1;
Zf_ss = 1;
L_ss = 1;

model;

%  discount factor 
BETAi= ZETAi*(1+ exp(Cih))^(-ETA);
%discount factor 
 BETAs= ZETAs*(1+ exp(Cif))^(-ETA);
%equation 2  /budget constraint investment home country
exp(Cih)+(exp(q1)*exp(kihh))+(exp(q2)*exp(kihf))= exp(Wh)+(exp(q1)+exp(Rkh))*exp(kihh(-1))+(exp(q2)+exp(Rkf))*exp(kifh(-1))+ exp(Bih)-exp(Rh(-1))*exp(Bih(-1));

%equation 2' /budget constraint investment foreign country
exp(Cif)+ (exp(q1)*exp(kiff))+(exp(q2)*exp(kifh))= exp(Wf)+(exp(q1)+exp(Rkf))*exp(kiff(-1))+(exp(q2)+exp(Rkf))*exp(kihf(-1))+exp(Bif)-exp(Rf(-1))*exp(Bif(-1));

%equation 4 /FOC capital 1 home country 
exp(Cih)^(-SIGMA)=BETAi *exp(Cih(+1))^(-SIGMA)*(exp(q1(+1))+exp(Rkh(+1)))/exp(q1);

%equation 4' / FOC capital 1 foreign country
exp(Cif)^(-SIGMA) = BETAi * exp(Cif(+1))^(-SIGMA)*(exp(q1(+1))+exp(Rkf(+1)))/exp(q1);

%equation 5 / FOC capital 2 home country
exp(Cih)^(-SIGMA) = BETAi * exp(Cih(+1))^(-SIGMA)*(exp(q2(+1))+exp(Rkh(+1)))/exp(q2);

%equation 5'  / FOC capital 2 foreign country
exp(Cif)^(-SIGMA) = BETAi * exp(Cif(+1))^(-SIGMA)*(exp(q2(+1))+exp(Rkf(+1)))/exp(q2);

%equation 6    / FOC bond home country
exp(Cih)^(-SIGMA) = BETAi * exp(Cih(+1))^(-SIGMA)*exp(Rh(+1));

%equation 6'   / FOC bond foreign country
exp(Cif)^(-SIGMA) = BETAi * exp(Cif(+1))^(-SIGMA)*exp(Rf(+1));

%equation9  / budget constraint savers home country
exp(Csh)+ (exp(q1)*exp(ksh))= exp(Wh)+ exp(q1)*exp(ksh(-1))+ exp(Zh)*exp(ksh(-1))^(OMEGA)+ exp(Bsh)- (exp(Rh(-1))*exp(Bsh(-1)));

%equation 9'   / budget constraint saves foreign country
exp(Csf)+ exp(q1)*exp(ksf)= exp(Wh)+ exp(q1)*exp(ksh(-1))+exp(Zh)*exp(ksh(-1))^(OMEGA)+ exp(Bsh)- exp(Rh(-1))*exp(Bsh(-1));

%equation 11   / FOC capital home country
exp(Csh)^(-SIGMA) = BETAs*exp(Csh(+1))^(-SIGMA)*(exp(q1(+1))+ OMEGA*exp(Zh(+1))*exp(ksh)^(OMEGA-1))/exp(q1);


%equation 11'  / FOC capital foreign country
exp(Csf)^(-SIGMA) = BETAs*exp(Csf(+1))^(-SIGMA)*(exp(q2(+1))+ OMEGA*exp(Zf(+1))*exp(ksf)^(OMEGA-1))/exp(q2);

% Equation 12   / FOC bond home country (for savers)
exp(Csh)^(-SIGMA) = BETAs*exp(Csh(+1))^(-SIGMA)*exp(Rh);

% Equation 12' / FOC bond foreign country (for savers)
exp(Csf)^(-SIGMA) = BETAs*exp(Csf(+1))^(-SIGMA)* exp(Rf);

%Equation 15  / FOC with respect to L home country
exp(Wh)= exp(Ah)*(1-Epsilon)*((N*(exp(kihh)+exp(kifh)))/exp(L))^(Epsilon);

%Equation 15' / FOC with respect to L foreign country
exp(Wf)= exp(Af)*(1-Epsilon)*((N*(exp(kiff)+exp(kihf)))/exp(L))^(Epsilon);

%Equation 16  / FOC with respect to K home country
exp(Rkh)= exp(Ah)*(1-Epsilon)*((N*(exp(kihh)+exp(kifh)))/exp(L))^(Epsilon);

%Equation 16'   / FOC with respect to K foreign country
exp(Rkf)=exp(Af)*(1-Epsilon)*((N*(exp(kiff)+exp(kihf)))/exp(L))^(Epsilon);

%equation 17    / AR(1) home investors
Ah= RHO*Ah(-1)+ nu1;
  
% equation 18  / AR(1) foreign investors
Af= RHO*Af(-1)+ nu2;

%The stochastic process for home goods productivity in the home and foreign countries
% equation 19   / AR(1) Home savers
Zh= RHO*Zh(-1)+ nu3;

% equation20  / AR(1) foreign savers
Zf= RHO*Zf(-1)+ nu4;

%Equation 21    / bond market home country
N*exp(Bih)+(1-N)*exp(Bsh)=0;

%Equation 22 / fixed asset equilibrium home country
N*exp(kihh)+N*exp(kifh)+(1-N)*exp(ksh)=1;

%Equation 21'   / bond market foreign country
N*exp(Bif)+(1-N)*exp(Bsf)=0;

%Equation 22'   / fixed asset equilibrium foreign country
N*exp(kiff) + N*exp(kihf)+ (1-N)*exp(ksf)=1;

% Production functions for home and foreign investors
% equation 13    / production function home investors
exp(Yh)= exp(Ah)*(N*(exp(kihh(-1))+exp(kifh(-1))))^Epsilon*exp(L)^(1-Epsilon);
%Kh =kihh+kifh
% equation 13'    / production function foreign investors
exp(Yf)= exp(Af)*(N*(exp(kiff(-1))+exp(kihf(-1))))^Epsilon*exp(L)^(1-Epsilon);
%Kf = kiff+kihf
% Production functions for home and foreign savers
%equation 14     / production function home savers 
exp(Jh)= exp(Zh)*(exp(ksh(-1)))^(OMEGA);
% equation 14'    / production function foreign savers
exp(Jf)= exp(Zf)*(exp(ksf(-1)))^(OMEGA);

% equation 23    / world market clearing condition
N*exp(Cih)+ N*exp(Cif)+(1-N)*(exp(Csh)+exp(Csf))= exp(Yh)+exp(Yf)+(1-N)*exp(Jh)+(1-N)*exp(Jf);

end;

initval;
Rh = log(Rh_ss);
Rf = log(Rf_ss);
Rkh = log(Rkh_ss);
Rkf = log(Rkf_ss);
kihh = log(kihh_ss);
kihf = log(kihf_ss);
kiff = log(kiff_ss);
kifh =log(kifh_ss);
ksh = log(ksh_ss);
ksf = log(ksf_ss);
Yh = log(Yh_ss);
Yf = log(Yf_ss);
Jh = log(Jh_ss);
Jf = log(Jf_ss);
Cih = log(Cih_ss) ;
Cif = log(Cif_ss); 
Csh = log(Csh_ss);
Csf =log(Csf_ss);
Wh = log(Wh_ss);
Wf = log(Wf_ss);
L = log(L_ss);
Bih = 0;
Bif = 0;
Bsh  =0;
Bsf = 0;
q1 = log(q1_ss);
q2 = log(q2_ss);
Ah = log(Ah_ss);
Af = log(Af_ss);
Zh = log(Zh_ss);
Zf = log(Zf_ss);
nu1 =0;
nu2 =0;
nu3 =0;
nu4 = 0;
end;
shocks;
var nu1 =1;
var nu2 =1;
var nu3 =1;
var nu4 =1;
end;
steady; 

stoch_simul(hp_filter=1600,order=1,irf=40);