% The model to be solved is a version of Hansen and Sargent (1987).
%          
%          max E sum{ beta^t (log(c(t)) -a1*(n1(t)-n2(t)) -a2*n2(t)
%          s.t.
%              c(t) + i(t) = xi(t) k^theta {h1*(n1(t)-d/2*(n1(t)-n1(t-1))
%                            ^2)^(1-theta)+h2*n2(t)^(1-theta)} 
%              k(t+1)=(1-delta)/lambda k(t) + (1/lambda) i(t)
%              xi(t+1)=rho0+rho1*xi(t)+rho2*xi(t-1)+epsilon(t+1)
%
%  States:       [1 xi(t) k(t) n1(t-1)]'
%  Controls:     [n1(t) n2(t) i(t)]'
%  Vars in u(.): [xi(t) k(t) n1(t-1) n1(t) n2(t) i(t)]'
%
% This file is called by 'gaussmle.m'.
%

d=15; 
h1=.46; 
h2=.13; 
theta=.36;
aa=2; 

delta=.023; 
lambda=1.005186;
rho0=.05; 
rho1=.95; 
rho2=0;

sigma=.0075;
beta=.99; 



P=[..
     3.142873885834261e+00
     2.429905208221457e-01
     2.774547493434945e-01
     4.699367578603676e+00
     7.577733793429229e-01
     1.150064812098385e+00
    -1.528591531851431e-01
     2.266116310096529e-02
     2.347237224297064e+00
     9.917699037252606e-01];
P=[..
    1.913445505006352e+00
     2.311074502102888e-07
     5.386396789976379e-01
     7.328212542279084e+00
     9.067369548678742e-01
     1.115650126949900e+00
    -1.181068856623376e-01
     2.266067700557153e-02
     3.090020152495124e+00
     9.917723796283989e-01];
P=[..
     1.480055615731623e+01
     4.972739146866751e-01
     3.746246720263898e-10
     3.023992649783557e+00
     2.130502534675324e+00
     1.266854875457478e+00
    -2.774770724345758e-01
     2.057701550613851e-02
     8.941300077988720e+00
     9.906373665804667e-01];
P=[..
     7.433662912667610e+00
     4.986093881297829e-01
     .1
     3.036942324822586e+00
     8.882764300182777e+01
    -2.131586778838403e-01
     7.943914260054218e-01
     2.241288294026891e-02
     2.965557069670502e+00
     9.915112701265081e-01];

    % 4.126601479755413e-05

d=P(1); h1=P(2); h2=P(3); aa=P(4); rho0=P(5); rho1=P(6); rho2=P(7);
delta=P(8); sigma=P(9); beta=P(10);

param1=[d,h1,h2,theta]';
param11=param1;
param2=aa;
param=[param1',param2 lambda rho0 rho1 rho2 delta sigma beta]';

a1=-aa*log(1-h1);
a2=-aa*log(1-h1-h2);
%
A=[1     0    0   0      0;
   rho0 rho1 rho2 0      0;
   0     1    0   0      0;
   0     0    0 (1-delta)/lambda  0;
   0     0    0   0      0];
B=[0 0 0;0 0 0;0 0 0;0 0 1/lambda;1 0 0];
D=[0 1 0 0 0;
   0 0 0 1 0;
   0 0 0 0 1;
   0 0 0 0 0;
   0 0 0 0 0;
   0 0 0 0 0];
E=[0 0 0;
   0 0 0;
   0 0 0;
   1 0 0;
   0 1 0;
   0 0 1];
Sigma=zeros(5);
Sigma(2,2)=sigma^2;
g=eye(5);
%

N=1;
IC=0;
nov1=6;
ct=1;
