Trying to get my runge-kutta to plot two graphs one of the original differential and one exact solution

This is what I have so far but it just sits there and doesnt graph anything. Is there something else I need to do or is there something I am missing?

clear,clc

function [x,y] = rk4(ti,tf,npts,yi,xi,F_tx,F_ty)

h=0.1;

ti = 0;

tf = 24;

t=ti:h:tf;

xi= 50;

yi= 50;

F_tx=@(x,y)(75-75*exp^(-.25*t));

F_ty=@(x,y)(0.25*(75-yi));

x=zeros(1,length(t));

x(1)=xi;

y=zeros(1,length(t));

y(1)=yi;

%t= linspace(ti,tf,length(t));

for i=1:(length(t)-1)

kx1 = F_tx(x(i),y(i));

ky1 = F_ty(x(i),y(i));

kx2 = F_tx(x(i)+0.5*h*kx1,y(i)+0.5*h*ky1);

ky2 = F_ty(x(i)+0.5*h*kx1,y(i)+0.5*h*ky1);

kx3 = F_tx((x(i)+0.5*h*kx2),y(i)+0.5*h*ky2);

ky3 = F_ty((x(i)+0.5*h*kx2),(y(i)+0.5*h*ky2));

kx4 = F_tx((x(i)+kx3*h),y(i)+ky3*h);

ky4 = F_ty((x(i)+kx3*h),(y(i)+ky3*h));

x(i+1) = x(i) + h*(kx1+2*kx2+2*kx3+kx4)/6;

y(i+1) = y(i) + h*(ky1+2*ky2+2*ky3+ky4)/6;

end

figure(1)

plot(t,y)

hold on

plot(t,x)

end

Answer 1

[Abraham Boayue]

So you meant to solve dT/dt = 0.25(75-T) with initial condition of T(0) = 0? If so, then the exaction solution is T(t) = 75-75e^(-0.25t). Let me know if this is what you wanted.

function [y, t,N] = rk4(f,to,tf,ya,h)
N= ceil((tf-to)/h);
halfh = h / 2;
y(1,:) = ya;

t(1) = to;
h6 = h/6;
for i = 1 : N
    t(i+1) = t(i) + h;
    th2 = t(i) + halfh;
    s1 = f(t(i), y(i,:));
    s2 = f(th2, y(i,:) + halfh * s1);
    s3 = f(th2, y(i,:) + halfh * s2);
    s4 = f(t(i+1), y(i,:) + h * s3);
    y(i+1,:) = y(i,:) + (s1 + s2+s2 + s3+s3 + s4) * h6;
end
%% This is your mfile. Run this alone after you have saved the above function in the same 
% folder as the mfile.
yo = 0;
h  = 0.1;
to = 0;
tf = 24;
Exact =@(t)(75-75*exp(-.25*t));
f = @(t,y)(0.25*(75-y));
[y, t,N] = rk4(f,to,tf,yo,h);
plot(t,Exact(t),'linewidth',2.5)
hold on
plot(t,y,'--','linewidth',2.5)
b = xlabel('t');
set(b,'Fontsize',15);
b = ylabel('y');
set(b,'Fontsize',15);
a= title('Exact Solution and Numerical solution');
set(b,'Fontsize',15);
legend('Exact Soln.','Numerical Soln.','location','best')
grid

[Zachary David]

so the original differential was dT/dt=0.25(75-T)

I think i obtained correctly the exact solution by using the integrating factor method which is where i got:

T(t)=75-75e^(-0.25t)

obviously by the look of the graph it doesnt really look right if i made the step size smaller it might look different but even then the Runge-Kutta is suppose to be very close as it is.

[Zachary David]

OK i definitely had the integration wrong sorry.

[Abraham Boayue]

The the integration was done correctly. Your asnswer for the exact solution is right.