Name: Ariful IslamTusher
Roll: 2002030
25.1
(a) Analytical solution :
syms t y
ode = y*t^3 - 1.5*y;
y_analytical = dsolve('Dy = y*t^3 - 1.5*y', 'y(0) = 1', 't');
t_values = 0:0.01:2;
y_values_analytical = double(subs(y_analytical, t_values));
figure;
plot(t_values, y_values_analytical, '-o');
title('Analytical Solution');
xlabel('t');
ylabel('y');
grid on;
(b) Euler’s method:
h1 = 0.5;
h2 = 0.25;
t_values_euler1 = 0:h1:2;
t_values_euler2 = 0:h2:2;
y_values_euler1 = zeros(size(t_values_euler1));
y_values_euler2 = zeros(size(t_values_euler2));
y_values_euler1(1) = 1;
y_values_euler2(1) = 1;
for i = 2:length(t_values_euler1)
y_prime = y_values_euler1(i-1)*t_values_euler1(i-1)^3 - 1.5*y_values_euler1(i-1);
y_values_euler1(i) = y_values_euler1(i-1) + h1*y_prime;
end
for i = 2:length(t_values_euler2)
y_prime = y_values_euler2(i-1)*t_values_euler2(i-1)^3 - 1.5*y_values_euler2(i-1);
y_values_euler2(i) = y_values_euler2(i-1) + h2*y_prime;
end
figure;
plot(t_values_euler1, y_values_euler1, '-o', t_values_euler2, y_values_euler2, '-o');
legend('h = 0.5', 'h = 0.25');
title("Euler's Method");
xlabel('t');
ylabel('y');
grid on;
(c) Fourth order RK method:
h_rk = 0.5;
t_values_rk = 0:h_rk:2;
y_values_rk = zeros(size(t_values_rk));
y_values_rk(1) = 1;
for i = 2:length(t_values_rk)
k1 = h_rk * (y_values_rk(i-1)*t_values_rk(i-1)^3 - 1.5*y_values_rk(i-1));
k2 = h_rk * ((y_values_rk(i-1) + k1/2)*(t_values_rk(i-1) + h_rk/2)^3 - 1.5*(y_values_rk(i-1) +
k1/2));
k3 = h_rk * ((y_values_rk(i-1) + k2/2)*(t_values_rk(i-1) + h_rk/2)^3 - 1.5*(y_values_rk(i-1) +
k2/2));
k4 = h_rk * ((y_values_rk(i-1) + k3)*(t_values_rk(i-1) + h_rk)^3 - 1.5*(y_values_rk(i-1) + k3));
y_values_rk(i) = y_values_rk(i-1) + (k1 + 2*k2 + 2*k3 + k4)/6;
end
figure;
plot(t_values_rk, y_values_rk, '-o');
title('Fourth-order Runge-Kutta Method');
xlabel('t');
ylabel('y');
grid on;
25.19:
(a) Euler’s Method:
% Given parameters
g = 9.81; % acceleration due to gravity (m/s^2)
cd = 0.225; % drag coefficient (kg/m)
m = 90; % mass (kg)
initial_height = 1000; % initial height (m)
% Time parameters
t_start = 0;
t_end = 100; % Choose a sufficiently large end time
h = 0.1; % time step
% Initialize arrays
t_values_euler = t_start:h:t_end;
v_values_euler = zeros(size(t_values_euler));
y_values_euler = zeros(size(t_values_euler));
% Initial conditions
v_values_euler(1) = 0; % initial velocity
y_values_euler(1) = initial_height; % initial height
% Euler's method
for i = 2:length(t_values_euler)
v_prime = g - (cd / m) * v_values_euler(i-1)^2;
v_values_euler(i) = v_values_euler(i-1) + h * v_prime;
y_values_euler(i) = y_values_euler(i-1) - h * v_values_euler(i-1);
% Break if the object hits the ground
if y_values_euler(i) <= 0
break;
end
end
% Plot results
figure;
subplot(2, 1, 1);
plot(t_values_euler, v_values_euler, '-o');
title("Euler's Method - Velocity vs Time");
xlabel('Time (s)');
ylabel('Velocity (m/s)');
grid on;
subplot(2, 1, 2);
plot(t_values_euler, y_values_euler, '-o');
title("Euler's Method - Distance vs Time");
xlabel('Time (s)');
ylabel('Distance (m)');
grid on;
(b) Forth-order RK method:
% Initialize arrays
t_values_rk = t_start:h:t_end;
v_values_rk = zeros(size(t_values_rk));
y_values_rk = zeros(size(t_values_rk));
% Initial conditions
v_values_rk(1) = 0; % initial velocity
y_values_rk(1) = initial_height; % initial height
% Fourth-order Runge-Kutta method
for i = 2:length(t_values_rk)
k1 = h * (g - (cd / m) * v_values_rk(i-1)^2);
k2 = h * (g - (cd / m) * (v_values_rk(i-1) + k1/2)^2);
k3 = h * (g - (cd / m) * (v_values_rk(i-1) + k2/2)^2);
k4 = h * (g - (cd / m) * (v_values_rk(i-1) + k3)^2);
v_values_rk(i) = v_values_rk(i-1) + (k1 + 2*k2 + 2*k3 + k4)/6;
y_values_rk(i) = y_values_rk(i-1) - h * v_values_rk(i-1);
% Break if the object hits the ground
if y_values_rk(i) <= 0
break;
end
end
% Plot results
figure;
subplot(2, 1, 1);
plot(t_values_rk, v_values_rk, '-o');
title("Fourth-order Runge-Kutta Method - Velocity vs Time");
xlabel('Time (s)');
ylabel('Velocity (m/s)');
grid on;
subplot(2, 1, 2);
plot(t_values_rk, y_values_rk, '-o');
title("Fourth-order Runge-Kutta Method - Distance vs Time");
xlabel('Time (s)');
ylabel('Distance (m)');
grid on;
