dT = 0.5;
T = 0:dT:8;
Y1(1) = 1;
for i = 1:length(T)-1
k1 = Y1(i)*(sin(T(i)))^3;
Y1(i+1) = Y1(i) + k1*dT;
k2 = Y1(i)*(sin(T(i+1)))^3;
Y1(i+1) = Y1(i) + 0.5*( k1 + k2 )*dT;
end
plot(T, Y1, '-s')
hold on
Y2(1) = 1;
for i = 1:length(T)-1
k1 = Y2(i)*(sin(T(i)))^3;
y_star = Y2(i) + 0.5*k1*dT;
k2 = y_star*(sin((T(i)+0.5*dT)))^3;
y_star = Y2(i) + 0.5*k2*dT;
k3 = y_star*(sin((T(i)+0.5*dT)))^3 ;
y_star = Y2(i) + k3*dT;
k4 = y_star*sin((T(i+1)))^3 ;
Y2(i+1) = Y2(i) + (1/6)*(k1 + 2*k2 + 2*k3 + k4)*dT;
end
plot(T, Y2, '-d')
Y3(1) = 1;
Y3(2) = Y1(2);
for i = 1:length(T)-2
k1 = Y3(i)*(sin(T(i)))^3;
k2 = Y3(i)*(sin(T(i+1)))^3 ;
Y3(i+2) = Y3(i+1) + 0.5*dT*(3*k2 - k1);
end
plot(T, Y3, '-p')
xlabel('t'), ylabel('y')
legend('Improved Euler', '4th order Runge-Kutta', '2nd order Adams-Bashforth')
Best Answer