%
% ex9 is a matlab file that simulates
% the solutions of the rbc model
%

% ----- save output to file ----- %

diary ex9.mout
diary off
delete ex9.mout
diary ex9.mout

% ----- display date and time of computation ----- %

%format long
date
time0 = clock;

% ----- define the parameters ----- %

psi      = .7/.3;
alpha    = 0.40;
delta    = 0.025;
beta     = 0.98;
thetabar = 1;
phi      = 0.95;

% ----- calibrated parameter ----- %

C = 0.0270;

% ----- define steady states ----- %

cbar = 0.8248;
nbar = 0.2480;
kbar = 9.3177;

% ----- enter matrix solutions ----- %

PP = [0.9467];

QQ = [1.3745];

RR = [ 0.6108; ...
        -0.2889 ];

SS = [-6.4210; ...
       10.1685 ];


% ----- generate 8000 simulated exogenous variables ----- %

% ---- initial log(theta) ---- %

ltheta(:,1) = 0;

% ---- for loop and random number generator ---- %

for i = 1:8000
 ltheta(:,i+1) = phi*ltheta(:,i) +  C*normrnd(0,1) ;
end

% ----- generate log(k/kbar) from states ----- %

% ---- initial log(k/kbar) ---- %

lk(:,1) = 0;

% ---- for loop ---- %

for i = 1:8000
 lk(:,i+1) = PP*lk(:,i) + QQ*ltheta(:,i) ;
end

% ----- generate log(c/cbar) and log(n/nbar) from states ----- %

% ---- for loop ---- %

for i = 1:8000
 y(:,i) = RR*lk(:,i) + SS*ltheta(:,i) ;
end

lc = y(1,:);
ln = y(2,:);

% ----- generate levels from steady states ----- %

k     = kbar*exp(lk);
c     = cbar*exp(lc);
n     = nbar*exp(ln);
theta = thetabar*exp(ltheta);

% ---- output ---- %

for i = 1:8000
 out(:,i) = theta(:,i)*k(:,i)^alpha*n(:,i)^(1-alpha) ;
end

% ---- investment ---- %

for i = 1:8000
 inv(:,i) = k(:,i+1) - (1-delta)*k(:,i) ;
end

% ----- put solutions into matrix ----- %

X = [out(:,1000:8000)',inv(:,1000:8000)',c(:,1000:8000)',n(:,1000:8000)'] ;

% ----- compute std deviation matrix ----- %

std(X)

% ----- compute correlation matrix ----- %

corrcoef(X)

% ----- exit ----- %

comptime = etime(clock, time0)
diary off