diff --git a/HW6Problem2.m b/HW6Problem2.m index c186f79..a56c1b2 100644 --- a/HW6Problem2.m +++ b/HW6Problem2.m @@ -1,20 +1,29 @@ +%Part a. +y_analytical = @(t) exp(-t); + %Part b. -dt = [0.001]; -figure() +dt = 0.01; %step size +a = figure(1); +set(0, 'defaultAxesFontsize', 16) +set(0, 'defaultTextFontsize', 16) +set(0, 'defaultLineLineWidth', 2) hold on -for k = 1:length(dt) - 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(k)*dy; - end - - plot(t, exp(-t), 'k-', t, y_n, 'o-'); +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. -figure() -y_heun = heun_sol_order1(@(t,y) ode2(t,y), 1, 1, [0 3]); -plot([0:1:3], y_heun); +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'); diff --git a/HW6Problem2.m~ b/HW6Problem2.m~ new file mode 100644 index 0000000..6f374b4 --- /dev/null +++ b/HW6Problem2.m~ @@ -0,0 +1,27 @@ +%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) +hold on +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), 1, 1, [0 3]); +plot([0:1:3], y_heun); +plot(t, y_analytical(t), 'k-', t, y_heun, 'o-'); diff --git a/HW6Problem4.m b/HW6Problem4.m new file mode 100644 index 0000000..fd9e0c7 --- /dev/null +++ b/HW6Problem4.m @@ -0,0 +1,27 @@ +%part a. +t = [0:.1:12]; +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; + + +plot(t, x_analytical(t)) + +%part b. +%dy = ode4(0, [1;0]); +dt = [0.1 0.01]; +figure() +hold on +for k = 1:length(dt) + t = [0:dt(k):12]'; + y_n = zeros(length(t), 2); + y_n(1, :) = [100 0]; %initial condition + for i = 2:length(t) + dy = ode4(t(i-1), y_n(i-1, :)); + y_n(i, :) = y_n(i-1, :) + dt(k)*dy; + end + %[t23, y23] = ode23(@(t,y) phugoid_ode(t,y), [0,20], [10 0 0 2]); + plot(t, x_analytical(t), '.', t, y_n(:,1), 's'); +end \ No newline at end of file diff --git a/HW6Problem4.m~ b/HW6Problem4.m~ new file mode 100644 index 0000000..1fce921 --- /dev/null +++ b/HW6Problem4.m~ @@ -0,0 +1,27 @@ +%part a. +t = [0:.1:12]; +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; + + +plot(t, x_analytical(t)) + +%part b. +dy = ode4(0, [1;0]); +dt = [0.1 0.01]; +figure() +hold on +for k = 1:length(dt) + t = [0:dt(k):12]'; + y_n = zeros(length(t), 4); + y_n(1, :) = [10 0 0 2]; %initial condition + for i = 2:length(t) + dy = phugoid_ode(t(i-1), y_n(i-1, :)); + y_n(i, :) = y_n(i-1, :) + dt(k)*dy; + end + [t23, y23] = ode23(@(t,y) phugoid_ode(t,y), [0,20], [10 0 0 2]); + plot(y_n(:,3), y_n(:, 4), '.', y23(:,3), y23(:,4), 's'); +end \ No newline at end of file diff --git a/HW6Problem5.m b/HW6Problem5.m new file mode 100644 index 0000000..35bead7 --- /dev/null +++ b/HW6Problem5.m @@ -0,0 +1,16 @@ +%Part b. +dy = ode3(0, [1;0]); +dt = [0.1 0.01]; +figure() +hold on +for k = 1:length(dt) + t = [0:dt(k):20]'; + y_n = zeros(length(t), 4); + y_n(1, :) = [10 0 0 2]; %initial condition + for i = 2:length(t) + dy = phugoid_ode(t(i-1), y_n(i-1, :)); + y_n(i, :) = y_n(i-1, :) + dt(k)*dy; + end + [t23, y23] = ode23(@(t,y) phugoid_ode(t,y), [0,20], [10 0 0 2]); + plot(y_n(:,3), y_n(:, 4), '.', y23(:,3), y23(:,4), 's'); +end \ No newline at end of file diff --git a/HW6Problem5.m~ b/HW6Problem5.m~ new file mode 100644 index 0000000..ceabf14 --- /dev/null +++ b/HW6Problem5.m~ @@ -0,0 +1,19 @@ +%Part b. +dy = ode3(0, [1;0]); +dt = [0.1]; +figure() +hold on +for k = 1:length(dt) + t = [0:dt(k):20]'; + y_n = zeros(length(t), 4); + y_n(1, :) = [10 0 0 2]; %initial condition + for i = 2:length(t) + dy = phugoid_ode(t, y_n(i-1, :)); + y_n(i, :) = y_n(i-1, :) + dt(k)*dy; + end + [t23, y23] = ode23(@(t,y) phugoid_ode(t,y), [0,20], [100 0 0 2]); + plot(y_n(:,3), y_n(:, 4)); + %plot(t, cos(3*t), 'k-', t, y_n(:,1), 'o-'); +end + +%Part c. \ No newline at end of file diff --git a/README.md b/README.md index c3494ff..2a77514 100644 --- a/README.md +++ b/README.md @@ -1 +1,39 @@ -# 06_initial_value_ode \ No newline at end of file +# 06_initial_value_ode + +## Problem 2 + +### Script: +```matlab +%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) +hold on +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 a. +![Analytical vs. Euler](./euler.png) \ No newline at end of file diff --git a/README.md~ b/README.md~ new file mode 100644 index 0000000..c3494ff --- /dev/null +++ b/README.md~ @@ -0,0 +1 @@ +# 06_initial_value_ode \ No newline at end of file diff --git a/euler.png b/euler.png new file mode 100644 index 0000000..d929569 Binary files /dev/null and b/euler.png differ diff --git a/heun.png b/heun.png new file mode 100644 index 0000000..cdfc2f6 Binary files /dev/null and b/heun.png differ diff --git a/ode4.m b/ode4.m new file mode 100644 index 0000000..d395ccd --- /dev/null +++ b/ode4.m @@ -0,0 +1,9 @@ +function dx = ode4(t,x) + %differential equation from hw6 problem 4 + g = 9.81; %m/s^2 + cd = 0.25; %kg/m + m = 60; %kg + dx = zeros(size(x)); + dx(1) = x(2); + dx(2) = g - (cd/m)*x(2)^2; +end \ No newline at end of file