Skip to content

mjb13028/06_initial_value_ode

master
Switch branches/tags

Name already in use

A tag already exists with the provided branch name. Many Git commands accept both tag and branch names, so creating this branch may cause unexpected behavior. Are you sure you want to create this branch?
Code

Latest commit

 

Git stats

Files

Permalink
Failed to load latest commit information.
Type
Name
Latest commit message
Commit time
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 
 

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.

Analytical vs. Euler

Part c.

Analytical vs. Heun

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.

Analytical vs. Euler

Part c.

Analytical vs. Heun

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.

Analytical vs. Euler

Part c.

Analytical vs. Heun

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');

Height of Plane vs. Distance

About

No description, website, or topics provided.

Resources

Stars

Watchers

Forks

Releases

No releases published

Languages