var c y w r z g k l I;
varexo e n;
parameters beta phi alpha delta rho theta psi;

/* la calibrazione č presa da Cooley e Prescott (1995) economic growth and business cycle */
alpha = 0.4;
beta = 0.987;
delta = 0.012;
psi = 1;
theta = 0.95;
phi = 0.5;
rho = 0.95;


model;
/* equazione di eulero per il consumo */
(1/c+phi*g) = beta*((1/c(+1)+phi*g(+1))*(1+r(+1)-delta));

/* offerta di lavoro*/
psi*c/(1-l) = w;

/* funzione di produzione */
y = exp(z)*(k(-1)^alpha)*(l^(1-alpha));

/*produttivitā marginale del lavoro */
w = (1-alpha)*exp(z)*(k(-1)^alpha)*(l^(-alpha));

/*produttivitā marginale del capitale */
r = alpha*exp(z)*(k(-1)^(alpha-1))*(l^(1-alpha)) - delta;

/* legge di accumulazione del capitale */
I = k-(1-delta)*k(-1);

/*market clearing condition */
y=c+I+g;

/* shocks */
z = rho*z(-1)+e;
g = theta*g(-1)+n;
end;

write_latex_dynamic_model;
write_latex_static_model;

initval;
k=21.396;
c=1.73154;
l=0.407922;
w=2.92451;
y=1.98829;

r=0;
z=0;
g=0;
end;

% DELETE the previous _steadystate file
delete([M_.fname '_steadystate.m'])
%options_ = rmfield(options_,'steadystate_flag'); 
steady(solve_algo=0);


M_.params(strmatch('K_SS' ,M_.param_names,'exact')) = oo_.steady_state(strmatch('k', M_.endo_names,'exact'));
M_.params(strmatch('C_SS' ,M_.param_names,'exact')) = oo_.steady_state(strmatch('c', M_.endo_names,'exact'));
M_.params(strmatch('L_SS',M_.param_names,'exact')) = oo_.steady_state(strmatch('l',M_.endo_names,'exact'));
M_.params(strmatch('W_SS',M_.param_names,'exact')) = oo_.steady_state(strmatch('w',M_.endo_names,'exact'));
M_.params(strmatch('Y_SS',M_.param_names,'exact')) = oo_.steady_state(strmatch('y',M_.endo_names,'exact'));

fprintf('\n...solve the SS a second time to update the ss ...\n')

steady(solve_algo=0);

if ~exist([M_.fname '_steadystate.m'])
    write_ss_on_dynarefile;
         fprintf('\n...solve the SS a second time to check the _steadystate file...')
         steady(solve_algo=0);
end

resid(1);

%%%% get the SS values %%%%
ss_vector = oo_.steady_state;
n_vars    = length(ss_vector);
ss_cell   = mat2cell(ss_vector,ones(1,n_vars));
tmp       = [M_.endo_names,ones(n_vars,1)*44]';
varnames  = ['[',tmp(:)',']'];
eval([varnames,' = deal(ss_cell{:});'])


disp(' ');disp('%%% Display some SS results : %%%%')
SSexpenditures  = c + I + g;
SSincome        = y;
disp(['expenditures in ss: ',num2str(SSexpenditures)]);
disp(['total income in ss: ',num2str(SSincome)]);
disp(['C/Y: ',num2str(c/y)]);
disp(['I/Y: ',num2str(I/y)]);
disp(['K/Y: ',num2str(k/y)]);
disp(['r (%annual): ',num2str(r*4)]);


shocks;
var e = 0.07;
var n = 0.07;

end;



stoch_simul(periods=150, order=1, irf=100,aim_solver);
forecast;
rplot y;