%%%%%%%%%%%%%%%%% DETERMINISTIC CAKE EATING MODEL %%%%%%%%%%%%%%%%%%%%%
%%%  THIS FILE SOLVES THE DETERMINISTIC CAKE EATING MODEL PRESENTED
%%%  IN LECTURE ONE OF RECURSIVE MACROECONOMICS OF MARTIN ELLISON 
%%%  (AVAILABLE ALSO IN CHAPTER 2 OF THE ADDA/COOPER BOOK, PAGE 16)

%--------------------------------------------------------------------------
%                           DEFINE PARAMETERS
%--------------------------------------------------------------------------

clear
dimIter=30;                     % number of iterations
beta=0.9;                       % discount factor
beta1=0.5;                      % a different value of the discount factor
K=0:0.01:1;                     % grid over cake size
[rowk,colK]=size(K);

%--------------------------------------------------------------------------
%                          VALUE FUNCTION ITERATION
%--------------------------------------------------------------------------

V=zeros(colK,dimIter);          % initialize the value function. Rows are cake size, columns the iteration
V1=zeros(colK,dimIter);          % initialize the value function. Rows are cake size, columns the iteration
for iter=1:dimIter              % loop over all current cake sizes
   
    aux=zeros(colK,colK)+NaN;
    aux1=zeros(colK,colK)+NaN;
    
    for ik=1:colK               % loop over next period cake size
        for ik2=1:(ik-1)
            aux(ik,ik2)=log(K(ik)-K(ik2))+beta*V(ik2,iter);
            aux1(ik,ik2)=log(K(ik)-K(ik2))+beta1*V1(ik2,iter);   
        end
    end
    V(:,iter+1)=max(aux')';     % optimizing over size of next period cake
    V1(:,iter+1)=max(aux1')';   % optimizing over size of next period cake
end

%--------------------------------------------------------------------------
%       CALCULATE OPTIMAL CONSUMPTION AS FCT. OF INITIAL CAKE SIZE
%--------------------------------------------------------------------------

[Val,Ind]=max(aux');
optK=K(Ind);
optK=optK+0*Val;
optC=K'-optK';

[Val1,Ind1]=max(aux1');
optK1=K(Ind1);
optK1=optK1+0*Val;
optC1=K'-optK1';

%--------------------------------------------------------------------------
%                              GRAPH RESULTS
%--------------------------------------------------------------------------

figure(1)
plot(K,V);
xlabel('Size of Cake');
ylabel('Value Function');

%figure(3)
%plot(K,V1);
%xlabel('Size of Cake');
%ylabel('Value Function');

figure(2)
plot(K,optC,'LineWidth',2)
hold on
plot(K,K','--r',...
    'LineWidth',2)
%hold on
%plot(K,optC1,'-g',...
 %   'LineWidth',2)

xlabel('Size of Cake');
ylabel('Optimal Consumption');
text(0.4,0.65,'45 degree line','FontSize',18)
text(0.4,0.13,'Optimal Consumption','FontSize',18)