How can I get values for n=500, 1000 and 2000 ? only getting values for n=4000!
1 view (last 30 days)
Show older comments
m Euler_1 Euler_2 RK_4
____ ________ __________ __________
500 0 0 0
1000 0 0 0
2000 0 0 0
4000 0.093706 2.0728e-05 8.9537e-13
clc
clear all
close all
%% ----------1-STAGE-EULER----------
n = [500,1000,2000,4000];
T = table();
for p = 1:4
steps = n(p);
T.m(p) = steps;
end
init_cond = [1 1i];
G = 1;
M = 1;
t_f = 2*pi/norm(1i);
dt = t_f/steps;
t = 0:dt:t_f;
z = init_cond;
f = @(t,z) [z(2),-G*M*z(1)/norm(z(1))^3];
for i = 1:steps
k1 = f(t(i),z(i,:));
z(i+1,:) = z(i,:) + dt.*k1;
end
T.Euler_1(p) = abs(init_cond(1) - z(end,1));
clearvars -except T steps n p Euler_1
%% ----------2-STAGE-EULER----------
init_cond = [1 1i];
G = 1;
M = 1;
t_f = 2*pi/norm(1i);
dt = t_f/steps;
t = 0:dt:t_f;
z = init_cond;
f = @(t,z) [z(2),-G*M*z(1)/norm(z(1))^3];
for i = 1:steps
k1 = f(t(i),z(i,:));
k2 = f(t(i)+dt,z(i,:)+dt*k1);
z(i+1,:) = z(i,:) + dt/2.*(k1 + k2);
end
T.Euler_2(p) = abs(init_cond(1) - z(end,1));
clearvars -except T steps n p Euler_2
%% ----------4-STAGE-CLASSICAL-RK----------
init_cond = [1 1i];
G = 1;
M = 1;
t_f = 2*pi/norm(1i);
dt = t_f/steps;
t = 0:dt:t_f;
z = init_cond;
f = @(t,z) [z(2),-G*M*z(1)/norm(z(1))^3];
for i = 1:steps
k1 = f(t(i),z(i,:));
k2 = f(t(i) + 0.5*dt,z(i,:) + 0.5*dt.*k1);
k3 = f((t(i) + 0.5*dt),(z(i,:) + 0.5*dt.*k2));
k4 = f((t(i) + dt),(z(i,:) + k3.*dt));
z(i+1,:) = z(i,:) + (1/6)*(k1 + 2*k2 + 2*k3 + k4).*dt;
end
T.RK_4(p) = abs(init_cond(1) - z(end,1));
clearvars -except T steps n p
figure
loglog(T.m,T.Euler_1)
hold on
loglog(T.m,T.Euler_2)
loglog(T.m,T.RK_4)
xlabel('n'), ylabel('Error'), grid on
title('Error Mag. vs. Number of Points (n)')
legend('1 Step Euler','2 Step Euler','Runge-Kutta 4')
T
0 Comments
Accepted Answer
Alan Stevens
on 16 Oct 2020
Your assignments
T.Euler_1(p) = abs(init_cond(1) - z(end,1));
etc. are all outside of the for p = ... loop, hence they just retain the last value of p, namely 4, and only store values in T.Euler_1(4) etc. overwriting anything previously stored there, hence ending with only the n = 4000 values.
More Answers (0)
See Also
Categories
Find more on Entering Commands in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!