4th order runge-kutta method to solve two 1st order differential equations
25 views (last 30 days)
Show older comments
%%4th order runge-kutta for Project 2
fy=@(t,y,z) z;
fz=@(t,y,z) (36000/(1500-40*t))-9.81-((0.1962*z^2)/(1500-40*t));
%boundary conditions
t(1) =0;
z(1) =0;
y(1)=0;
j=0;
h=1; %%step size
tfinal=10;
N=tfinal/h;
%Update loop
for j=1:N
t(j+1)=t(j)+h;
k1y=fy(t(j),y(j),z(j));
k1z=fz(t(j),y(j),z(j));
k2y=fy(t(j)+0.5*h , y(j)+0.5*h*k1y , z(j)+0.5*h*k1z);
k2z=fz(t(j)+0.5*h,y(j)+0.5*h*k1y,z(j)+0.5*h*k1z);
k3y=fy(t(j)+0.5*h,y(j)+0.5*h*k2y,z(j)+0.5*h*k2z);
k3z=fz(t(j)+h/2,y(j)+0.5*h*k2y,z(j)+0.5*h*k2z);
k4y=fy(t(j)+h,y(j)+h*k3y,z(j)+h*k3z);
k4z=fz(t(j)+h,y(j)+h*k3y,z(j)+h*k3z);
y(j+1)=y(j)+h/6*(k1y+2*k2y+2*k3y+k4y);
z(j+1)=z(j)+h/6*(k1z+2*k2z+2*k3z+k4z);
end
fprintf('t = %6.4f, y = %6.4f\n', t, y);
so i need to display y values from t=0 to t=10 but the t output is not as expected. and y range shouldnt be a zero in it. please help
[SL changed the last line from a comment in the code to text in a text region.]
0 Comments
Accepted Answer
Alan Stevens
on 13 Jul 2022
Like this?
%%4th order runge-kutta for Project 2
fy=@(t,y,z) z;
fz=@(t,y,z) (36000/(1500-40*t))-9.81-((0.1962*z^2)/(1500-40*t));
%boundary conditions
z(1) =0;
y(1)=0;
h=1; %%step size
tfinal=10;
N=tfinal/h;
%Update loop
t = 0:h:N*h;
for j=1:N
k1y=fy(t(j),y(j),z(j));
k1z=fz(t(j),y(j),z(j));
k2y=fy(t(j)+0.5*h, y(j)+0.5*h*k1y, z(j)+0.5*h*k1z);
k2z=fz(t(j)+0.5*h, y(j)+0.5*h*k1y, z(j)+0.5*h*k1z);
k3y=fy(t(j)+0.5*h, y(j)+0.5*h*k2y, z(j)+0.5*h*k2z);
k3z=fz(t(j)+0.5*h, y(j)+0.5*h*k2y, z(j)+0.5*h*k2z);
k4y=fy(t(j)+h, y(j)+h*k3y, z(j)+h*k3z);
k4z=fz(t(j)+h, y(j)+h*k3y, z(j)+h*k3z);
y(j+1)=y(j)+h/6*(k1y+2*k2y+2*k3y+k4y);
z(j+1)=z(j)+h/6*(k1z+2*k2z+2*k3z+k4z);
end
disp([' t y'])
disp([t' y'])
plot(t,y,'o-',t,z,'*-'),grid
xlabel('t'), ylabel('y & z')
legend('y','z')
More Answers (0)
See Also
Categories
Find more on Numerical Integration and Differential Equations in Help Center and File Exchange
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!