%% Euler's and Midpoint Method Diff. Eq. Example
% Heat transfer problem II.5

%% Explicit Euler Solution

T=2500;
Ta=250;
alpha=4e-12;
dt=1;
N=10;

for n=1:N 
    display(['tn=',num2str((n-1)*dt,'%1.1e'),' Tn=',num2str(T,'%1.16e')]); %,' T=',num2str(Ta,'%1.16e'),' err=',num2str(err,'%1.16e')]);
    f=-alpha*(T^4-Ta^4);
    T=T+dt*f;
end;
    Tr=1758.26337470;
    err=T-Tr;
    display(['tn=',num2str((n)*dt,'%1.1e'),' Tn=',num2str(T,'%1.16e'),' T=',num2str(Tr,'%1.16e'),' err=',num2str(err,'%1.16e'),char(10)]);

%% Midpoint Method Solution

T=2500;
Ta=250;
alpha=4e-12;
dt=2;
N=5;

for n=1:N 
    display(['tn=',num2str((n-1)*dt,'%1.1e'),' Tn=',num2str(T,'%1.16e')]); %,' T=',num2str(Ta,'%1.16e'),' err=',num2str(err,'%1.16e')]);
    f=-alpha*(T^4-Ta^4);
    Tp=T+0.5*dt*f;
    fp=-alpha*(Tp^4-Ta^4);
    T=T+dt*fp;
end;
   Ta=1758.26337470;
   err=T-Ta;
   display(['tn=',num2str((n)*dt,'%1.1e'),' Tn=',num2str(T,'%1.16e'),' T=',num2str(Ta,'%1.16e'),' err=',num2str(err,'%1.16e')]);
