1) For the differential equation y\' = y, y(0) = 1, the value of the solution fu
ID: 3209694 • Letter: 1
Question
1) For the differential equation y' = y, y(0) = 1, the value of the solution function at x=1 is e. Run each of the three numerical algorithms - Euler, Improved Euler, and Runge-Kutta - with N = 50, 100, and 200 points to estimate the value of e to 9 decimal places.
2) The differential equation y' = 9.8 - .00016 (100y + 10y^2 + y^3), y(0) = 0, gives a skydiver's velocity x seconds after jumping out of a plane. Run all three numerical algorithms with 20 points to approximate her velocity after 20 seconds. Keep 5 digits after the decimal since that is where the approximations differ from each other. You do not have to give me all 20 points, just the last one (although it might be interesting to compare the three sets of 20 points).
Explanation / Answer
%%%% Matlab code %%%%
1)
clc;
clear all;
close all;
format long
format compact;
f=@(x,y) y;
N=[50 100 200];
ye(1)=1;
yi(1)=1;
yr(1)=1;
fprintf(' N Eular Method Improved Eular Method RK-4 Method ');
for n=1:length(N)
t=linspace(0,1,N(n));
h=t(3)-t(2);
for k=1:length(t)-1
%%% Eular
ye(k+1)=ye(k)+h*feval(f,t(k),ye(k));
%%%Modified Eular
yp=yi(k)+h*feval(f,t(k),yi(k));
yi(k+1)= yi(k)+h/2*(feval(f,t(k),yi(k))+feval(f,t(k+1),yp));
%%% Rk-4 Method
k1=h*feval(f,t(k),yr(k));
k2=h*feval(f,t(k)+h/2,yr(k)+k1/2);
k3=h*feval(f,t(k)+h/2,yr(k)+k2/2);
k4=h*feval(f,t(k)+h,yr(k)+k3);
yr(k+1)=yr(k)+(k1+2*k2+2*k3+k4)/6;
end
fprintf(' %d %1.9f %1.9f %1.9f ',N(n),ye(k+1),yi(k+1),yr(k+1));
end
OUTPUT:
N Eular Method Improved Eular Method RK-4 Method
50 2.691053247 2.718096008 2.718281825
100 2.704679036 2.718235953 2.718281828
200 2.711483285 2.718270431 2.718281828
2)
%%%Matlab code %%%
clc;
clear all;
close all;
format long
format compact;
f=@(x,y) 9.8-0.00016*(100*y + 10*y^2 + y^3);
N=[20];
ye(1)=0;
yi(1)=0;
yr(1)=0;
fprintf(' time Eular Method Improved Eular Method RK-4 Method ');
for n=1:length(N)
t=linspace(0,20,N(n)+1);
h=t(3)-t(2);
for k=1:length(t)-1
%%% Eular
ye(k+1)=ye(k)+h*feval(f,t(k),ye(k));
%%%Modified Eular
yp=yi(k)+h*feval(f,t(k),yi(k));
yi(k+1)= yi(k)+h/2*(feval(f,t(k),yi(k))+feval(f,t(k+1),yp));
%%% Rk-4 Method
k1=h*feval(f,t(k),yr(k));
k2=h*feval(f,t(k)+h/2,yr(k)+k1/2);
k3=h*feval(f,t(k)+h/2,yr(k)+k2/2);
k4=h*feval(f,t(k)+h,yr(k)+k3);
yr(k+1)=yr(k)+(k1+2*k2+2*k3+k4)/6;
fprintf(' %f %1.5f %1.5f %1.5f ',t(k+1),ye(k+1),yi(k+1),yr(k+1));
end
end
OUTPUT;
time Eular Method Improved Eular Method RK-4 Method
1.000000 9.80000 9.56947 9.63575
2.000000 19.13895 18.16880 18.38507
3.000000 26.92495 24.88411 25.29579
4.000000 32.01113 29.40570 29.94286
5.000000 34.41106 32.14238 32.66947
6.000000 35.24639 33.70213 34.12876
7.000000 35.48886 34.56401 34.86934
8.000000 35.55442 35.03278 35.23475
9.000000 35.57177 35.28565 35.41251
10.000000 35.57634 35.42146 35.49838
11.000000 35.57754 35.49424 35.53973
12.000000 35.57785 35.53320 35.55960
13.000000 35.57794 35.55403 35.56915
14.000000 35.57796 35.56517 35.57373
15.000000 35.57796 35.57113 35.57593
16.000000 35.57796 35.57431 35.57699
17.000000 35.57796 35.57601 35.57750
18.000000 35.57796 35.57692 35.57774
19.000000 35.57796 35.57741 35.57786
20.000000 35.57796 35.57767 35.57791
Related Questions
Navigate
Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.