%This script simulates the path of capital stock from the NR method and
%compares it to the known solution. Writen by Philip Shaw, 2022;
clear all
global alpha beta frac T
T=100; %set the time frame for the simulation
alpha=.33; %production function parameter
beta=.98; %preference parameter
frac=.05; %start the economy at a fraction of steady state capital
KT=(alpha*beta)^(1/(1-alpha)); %steady state capital stock
k0=frac*KT; %initial capital stock
k01=k0; %store initial capital stock for later use
x0=ones(1,T-2)*.16; %initialize guess for capital stock sequence
[xstar]=NR(@(x)eulereq(x),x0); %final solution from N-R approximation
 %simulate capital stock over time with starting and terminal values
 %starting with known solution
 for i=1:T-1
     kt(i,1)=alpha*beta*k0^alpha; %known solution
     k0=kt(i); %update capital stock for iterations
 end
kapprox=[k01 xstar' KT]; %the approximate solution
plot(1:1:T,[k01 xstar' KT],1:1:T,[k01 kt']) %plot approximate and true solution
legend('Approximate Solution','Known Solution')
xlabel('Time')
ylabel('Capital Stock')
title('Evolution of Capital Stock Over Time')