06_initial_value_ode
Problem 2
Script:
% Part a.
y_analytical = @(t ) exp (-t );
% Part b.
dt = 0.01 ; % step size
a = figure (1 );
set (0 , ' defaultAxesFontsize' , 16 )
set (0 , ' defaultTextFontsize' , 16 )
set (0 , ' defaultLineLineWidth' , 2 )
t = [0 : dt : 3 ]' ;
y_n = zeros (size (t ));
y_n (1 ) = 1 ; % initial condition
for i = 2 : length (t )
dy = ode2 (t , y_n (i - 1 ));
y_n (i ) = y_n (i - 1 ) + dt * dy ;
end
plot (t , y_analytical (t ), ' k-' , t , y_n , ' o-' );
title (' Analytical vs. Euler' )
legend (' Analytical' , ' Euler' , ' Location' , ' Northeast' )
saveas (a , ' euler.png' );
% Part c.
b = figure (2 );
y_heun = heun_sol_order1 (@(t ,y ) ode2 (t ,y ), dt , 1 , [0 3 ]);
plot (t , y_analytical (t ), ' k-' , t , y_heun , ' o-' );
title (' Analytical vs. Heun' )
legend (' Analytical' , ' Heun' , ' Location' , ' Northeast' )
saveas (b , ' heun.png' );
Part b.
Part c.
Problem 3
Script:
% Part a.
y_analytical = @(t ) cos (3 * t );
% Part b.
dt = 0.001 ; % step size
a = figure (1 );
set (0 , ' defaultAxesFontsize' , 16 )
set (0 , ' defaultTextFontsize' , 16 )
set (0 , ' defaultLineLineWidth' , 2 )
t = [0 : dt : 3 ]' ;
y_n = zeros (length (t ), 2 );
y_n (1 , : ) = [1 0 ]; % initial condition
for i = 2 : length (t )
dy = ode3 (t , y_n (i - 1 , : ));
y_n (i , : ) = y_n (i - 1 , : ) + dt * dy ;
end
plot (t , y_analytical (t ), ' k-' , t , y_n (: ,1 ), ' o-' );
title (' Analytical vs. Euler' )
legend (' Analytical' , ' Euler' , ' Location' , ' Northeast' )
saveas (a , ' euler3.png' );
% Part c.
b = figure (2 );
dt = 0.001 ;
t = [0 : dt : 3 ];
y_heun = heun_sol_order2 (@(t ,y ) ode3 (t ,y ), dt , [1 0 ], [0 3 ]);
plot (t , y_analytical (t ), ' k-' , t , y_heun (: ,1 ), ' o-' );
title (' Analytical vs. Heun' )
legend (' Analytical' , ' Heun' , ' Location' , ' Northeast' )
saveas (b , ' heun3.png' );
Part b.
Part c.
Problem 4
Script:
% part a.
g = 9.81 ; % m/s^2
cd = 0.25 ; % kg/m
m = 60 ; % kg
v_term = sqrt (m * g / cd );
x_analytical = @(t ) ((v_term ^ 2 )/g )*log (abs (cosh ((g * t )/v_term ))) + 100 ;
% part b.
dt = 0.01 ; % step size
a = figure (1 );
set (0 , ' defaultAxesFontsize' , 16 )
set (0 , ' defaultTextFontsize' , 16 )
set (0 , ' defaultLineLineWidth' , 2 )
t = [0 : dt : 12 ]' ;
y_n = zeros (length (t ), 2 );
y_n (1 , : ) = [100 0 ]; % initial condition
for i = 2 : length (t )
dy = ode4 (t , y_n (i - 1 , : ));
y_n (i , : ) = y_n (i - 1 , : ) + dt * dy ;
end
plot (t , x_analytical (t ), ' k-' , t , y_n (: ,1 ), ' o-' );
title (' Analytical vs. Euler' )
legend (' Analytical' , ' Euler' , ' Location' , ' Northeast' )
saveas (a , ' euler4.png' );
% Part c.
b = figure (2 );
dt = 0.01 ;
t = [0 : dt : 12 ];
y_heun = heun_sol_order2 (@(t ,y ) ode4 (t ,y ), dt , [100 0 ], [0 12 ]);
plot (t , x_analytical (t ), ' k-' , t , y_heun (: ,1 ), ' o-' );
title (' Analytical vs. Heun' )
legend (' Analytical' , ' Heun' , ' Location' , ' Northeast' )
saveas (b , ' heun4.png' );
Part b.
Part c.
Problem 5
Part a.
function dy = phugoid_ode (t ,y )
% glider equations describing phugoid path
% y = [v, theta, x, y]
g = 9.81 ; % m/s^2
vt = 5.5 ; % m/s
cl = 5.2 ;
cd = 1 ;
dy = zeros (size (y ));
dy (1 ) = - g * sin (y (2 )) - cd / cl * g / vt ^ 2 * y (1 )^2 ;
dy (2 ) = -(g / y (1 ))*cos (y (2 )) + g / vt ^ 2 * y (1 );
dy (3 ) = y (1 )*cos (y (2 ));
dy (4 ) = y (1 )*sin (y (2 ));
end
Part b.
% Part b.
dt = 0.1 ;
t = [0 : dt : 20 ]' ;
y_n1 = zeros (length (t ), 4 );
y_n1 (1 , : ) = [10 0 0 2 ]; % initial condition
for i = 2 : length (t )
dy = phugoid_ode (t (i - 1 ), y_n1 (i - 1 , : ));
y_n1 (i , : ) = y_n1 (i - 1 , : ) + dt * dy ;
end
dt = 0.01 ;
t = [0 : dt : 20 ]' ;
y_n2 = zeros (length (t ), 4 );
y_n2 (1 , : ) = [10 0 0 2 ]; % initial condition
for i = 2 : length (t )
dy = phugoid_ode (t (i - 1 ), y_n2 (i - 1 , : ));
y_n2 (i , : ) = y_n2 (i - 1 , : ) + dt * dy ;
end
a = figure (1 );
set (0 , ' defaultAxesFontsize' , 16 )
set (0 , ' defaultTextFontsize' , 16 )
set (0 , ' defaultLineLineWidth' , 2 )
[t23 , y23 ] = ode23 (@(t ,y ) phugoid_ode (t ,y ), [0 ,20 ], [10 0 0 2 ]);
plot (y_n1 (: ,3 ), y_n1 (: , 4 ), ' .' , y_n2 (: ,3 ), y_n2 (: , 4 ), ' o' , y23 (: ,3 ), y23 (: ,4 ), ' -' );
title (' Height of Plane vs. Distance' )
xlabel (' Distance (m)' )
ylabel (' Height (m)' )
legend (' dt = 0.1' , ' dt = 0.01' , ' analytical' , ' Location' , ' Northeast' )
saveas (a , ' problem5.png' );