clear all;
close all;

% Define the structural parameters of the model,
% i.e. policy invariant preference and technology parameters
% alpha : capital's share of output
% beta  : time discount factor
% delta : depreciation rate
% sigma : risk-aversion parameter, also intertemp. subst. param.
% gamma : unconditional expectation of the technology parameter

alpha = .35;
beta  = .98;
delta = .025;
sigma = 2;
gamma = 5;

% Find the steady-state level of capital as a function of 
%   the structural parameters
%   kstar = ((1/beta - 1 + delta)/(alpha*gamma))^(1/(alpha-1));

% Define the number of discrete values k can take
g = 3;

k = [2.85
     3.00
     3.15];

% (a)
% Compute a (3 × 3) dimensional consumption matrix c
%   for all the (3 × 3) values of k_t and k_t+1.

% vypocet v cyklu
for i = 1 : g
   for j = 1 : g
      c(i,j) = gamma*k(i)^alpha + (1-delta)*k(i) - k(j);
      % i is the counter for the state variable k_t
      % j is the counter for tue control variable k_t+1
   end
end

% efektivni vypocet



% vypocet v cyklu
% % Compute a (3 × 3) dimensional utility matrix u
% %   for all the (3 × 3) values of k_t and k_t+1.
% for i = 1 : g
%    for j = 1 : g
%       u(i,j) = (c(i,j)^(1-sigma) - 1)/(1-sigma);    
%       % i is the counter for the state variable k_t
%       % j is the counter for tue control variable k_t+1
%    end
% end

% efektivni vypocet



% (b)

v = [167.6
     168.1
     168.6];
 
% Transform the v vector into a 3x3 matrix where each row is identical
% First define a coloumn vector of ones
enere = [];
         

% Then premultiply the transformed v vector with this vector
tmp1 = [];

% alternativne pomoci fce repmat
tmp0 = [];

% Display the result and delete the variable
disp(''); % blank line
disp('v transformed:')
disp(tmp1);
clear tmp1

% Compute the sum {u(k_t,k_t+1) + beta*v(k_t+1)}
tmp2 = [];

% Display the result and delete the variable
disp(''); % blank line
disp('The matrix within the curly brakcets:')
disp(tmp2);
clear tmp2

% Find the largest element on each row,
%  i.e. for each k_t find the utility maximizing k_t+1

Tv = [];


disp(''); % blank line
disp('The resulting value from the maximizing choice for each k_t:')
disp(Tv);


% Extra 

% Find what element of k_t+1 which gave the optimal value

disp(''); % blank line
disp('Element of the k_t+1 vector chose:')
disp(i);

% Find what element of k_t+1 which gave the optimal value
kdecrule = [];
disp(''); % blank line
disp('Decision rule for k_t+1:')
disp(kdecrule);