%----------------------------------------------------------------
% 0. Close all graphic windows
%----------------------------------------------------------------

close all;

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

var  cf ci w nf ni k r t a b p q yf yi g gstar gsmall;

varexo ea eb; 

%----------------------------------------------------------------
% 2. Defining parameters and calibration
%----------------------------------------------------------------                                   
                                             
parameters  tao taostar phi kai alpha delta miu v gamma cbeta rhoa rhob
cf_ss ci_ss w_ss nf_ss ni_ss k_ss t_ss a_ss q_ss b_ss p_ss yf_ss yi_ss g_ss gstar_ss gsmall_ss r_ss;

tao = 0.21;
taostar = 0.16;
phi = 2;
kai = 0.1;
alpha = 0.32;
delta = 0.1;
miu = 0.3;
v = 0.3;
gamma = 2;
cbeta = 0.99;
rhoa = 0.9;
rhob = 0.8;
sigmaea = 0.03;
sigmaeb =0.01;

%----------------------------------------------------------------
% 3. Characteristic equations of the steady state
%----------------------------------------------------------------

q_ss = 1;
r_ss = (1/cbeta -1 + delta)/(1-taostar);
k_ss = nf_ss*a_ss*((1/cbeta -1 + delta)/(alpha*(1-taostar)))^(1/alpha-1);
w_ss = (1-alpha)*((1/cbeta -1 + delta)/(alpha*(1-taostar)))^(alpha/alpha-1);
yf_ss = (k_ss^(alpha))*(nf_ss)^((1-alpha));
yi_ss = (ni_ss)^(v);
nf_ss = (1 + ((w_ss*(1-tao))/(v*(b_ss^(v))*(1-(p_ss*phi*tao)))))^(1/(v-1));
ni_ss = (1-nf_ss)^(v);
cf_ss = yf_ss + k_ss*(r_ss - 1);
ci_ss = yi_ss*phi*tao*(1-p_ss);
p_ss = 1/(phi*tao);
gsmall_ss = yi_ss;
g_ss = (tao*w_ss*nf_ss)+(taostar*r_ss*k_ss)+(yi_ss);
gstar_ss = (1-kai)*g_ss;
t_ss = (1-kai-miu)+(kai*miu)*g_ss;
a_ss = 1;
b_ss = 0.75;


%----------------------------------------------------------------
% 4. Model
%----------------------------------------------------------------

model;

(cf+ci) = (cbeta*(cf(+1)+ci(+1))^(-gamma)*((1-delta)+((1-taostar)*r(+1))))^(1/(-gamma));      %1. Euler Equation
w = (-(v*q*b*(b*(1-nf))^(v-1))*(1-(p*phi*tao)))/(1-tao);      %2. Wages HH
cf = yf+(r*k(-1))-k;      %3. Market clearing condition cf
ci = (1-p)*phi*tao*yi;      %4. Market clearing condition ci
yi+yf = ci+cf+k+gsmall-(r*k(-1));      %5. Total production
1 = ni+nf;      %6. Labour endowment
yf = (k(-1)^(alpha))*((nf*a)^(1-alpha));      %7. Formal production function
yi = (b*ni)^(v);      %8. Informal production function
r = alpha*(k(-1)^(alpha-1))*((nf*a)^(1-alpha));      %9. Markup equation rt
w = (1-alpha)*(k(-1)^(alpha))*((nf*a)^(-alpha));     %10. Markup equation wt
log(a) = rhoa*log(a(-1))+ea;      %11. Formal technology
log(b) = rhob*log(b(-1))+eb;      %12. Informal technology
gsmall = p*phi*tao*yi;      %13. Government revenue from punishment
g = (tao*w*nf)+(taostar*r*k(-1))+gsmall;      %14. Total gov. revenue
gstar = (1-kai)*g;      %15. Gov. rev. after corruption
p = ni*gstar*(1+miu);      %16. Probability of auditing an informal
t = (1-miu)*gstar;      %17. Final transfer to household



end;


%----------------------------------------------------------------
% 5. Computation
%----------------------------------------------------------------

initval;

cf = cf_ss;
ci = ci_ss;
w = w_ss;
nf = nf_ss;
ni = ni_ss;
k = k_ss;
r = r_ss;
t = t_ss;
a = a_ss;
b = b_ss;
p = p_ss;
q = q_ss;
yf = yf_ss;
yi = yi_ss;
g = g_ss;
gstar = gstar_ss;
gsmall = gsmall_ss;
ea = 0;
eb = 0;
 

end;

check;

steady;

%----------------------------------------------------------------
% 6. Shocks
%----------------------------------------------------------------

shocks;
var ea = sigmaea^2;
var eb = sigmaeb^2;
var ea,eb = 0*sigmaea*sigmaeb;
end;


stoch_simul (periods=1000, order = 1, irf=10) cf ci nf ni yf yi t;