clear all
close all

% set parameters
phi = 0.75;
sigma = 1; 
nu =1;
pssi =1;
epss =10;
theta =1;
raua = 0.95;
raum = 0;
sigmam = 0.01;
sigmaa = 0.01;
betta = 0.95;


options = optimset('Display','Final','TolX',1e-2,'TolFun',1e-2);

x0 = [1,1,0.5,0.5,0.5,1,1,1,0.1,0.2,1,0,0,0];
x = fsolve(@(x)funn(x,betta, pssi, nu, sigma,phi,epss,theta),x0, options)
Cs = x(1);
is = x(2);
Ns = x(3);
ws = x(4);
ms = x(5);
mc_ss = x(6);
Ys = x(7); 
vs = x(8); 
x1s = x(9);
x2s = x(10); 
As =x(11);
pis = x(12);
pihs = x(13);
deltams = x(14);


save param_trial betta phi nu epss  sigma pssi theta raua raum sigmam sigmaa  Cs is Ns ws ms mc_ss Ys vs x1s x2s As pis pihs deltams

dynare trial noclearall nolog

figure
subplot(2,2,1)
plot(100*C_ea,'-k','Linewidth',2)
title('log(C)')

subplot(2,2,2)
plot(400*pi_ea,'-k','Linewidth',2)
title('pi')

subplot(2,2,3)
plot(400*w_ea,'-k','Linewidth',2)
title('w')

subplot(2,2,4)
plot(100*m_ea,'-k','Linewidth',2)
title('log(m)')



function F = funn(x,betta, pssi, nu, sigma,phi,epss,theta)



F(1) = exp(x(1))^(-sigma) -  betta* exp(x(1))^(-1)* (1+x(2)); 
F(2) = pssi*exp(x(3))^(nu) - exp(x(1))^(-sigma)*exp(x(4));
F(3) = exp(x(5)) - theta*(1+x(2))/x(2) *exp(x(1))^(sigma);
F(4) = exp(x(6)) - exp(x(4)) / exp(0);
F(5) = exp(x(1))- exp(x(7));
F(6) = exp(x(7)) - exp(0)*exp(x(3))/x(8);
F(7) = x(8)- 1;
F(8) = 1 - epss/(epss-1) * exp(x(9))/exp(x(10));
F(9) = exp(x(9)) - exp(x(1))^(-1) *exp(x(6)) *exp(x(7)) - phi*betta*exp(x(9)); 
F(10) = exp(x(10)) - exp(x(1))^(-1)  *exp(x(7)) - phi*betta*exp(x(10));
F(11) = exp(x(11)) - 1;
F(12) = x(12) - 0;
F(13) = x(13)-x(12);
F(14) = x(14) - 0; 

end
