%% Runge Kutta Method Diff. Eq. Example
% Heat transfer problem II.5

%% Runge Kutta 4th order

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

%% RK4

for n=1:N 
    display(['tn=',num2str((n-1)*dt,'%1.1e'),' Tn=',num2str(T,'%1.16e')]);
    f=-alpha*(T^4-Ta^4);
    dT1=dt*f;
    f=-alpha*((T+dT1/2)^4-Ta^4);
    dT2=dt*f;    
    f=-alpha*((T+dT2/2)^4-Ta^4);
    dT3=dt*f;
    f=-alpha*((T+dT3)^4-Ta^4);
    dT4=dt*f;
    T=T+(1/6)*(dT1+2*dT2+2*dT3+dT4);
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)]);

  