Skip to content
No description, website, or topics provided.
MATLAB
Branch: master
Clone or download
This branch is even with mjb13028:master.

Latest commit

Fetching latest commit…
Cannot retrieve the latest commit at this time.

Files

Permalink
Type Name Latest commit message Commit time
Failed to load latest commit information.
HW6Problem2.m
HW6Problem3.m
HW6Problem3.m~
HW6Problem4.m
HW6Problem4.m~
HW6Problem5.m
HW6Problem5.m~
README.md
README.md~
euler2.png
euler3.png
euler4.png
heun2.png
heun3.png
heun4.png
heun_sol.m
heun_sol.m~
heun_sol_order1.m
heun_sol_order2.m
ode2.m
ode3.m
ode4.m
problem5.png

README.md

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

You can’t perform that action at this time.