Download
<< Back
OutputCode
%% Q2 Solution
%% a) Linear Oscillator
%
% We solve linear oscillator with $a = 0$, $b = 9$, $c = 10$, $\omega = 1$
% and $x_0=(2,1)^T$ using Milne-Simpson 2, Nystrom 2, Adams-Moulton 2, and ODE23
figure
axis([0 10 -4 4]);
hold on;
xlabel('x_1');
ylabel('x_2');
[t,x]=ms2(@oscillator,0,10, [2;1], .1);
plot(t,x(:,1),'ro','DisplayName','Milne-Simpson 2');
[t,x]=ny2(@oscillator,0,10, [2;1], .1);
plot(t,x(:,1),'mx','DisplayName','Nystrom 2');
[t,x]=am2(@oscillator,0,10, [2;1], .1);
plot(t,x(:,1),'b+','DisplayName','Adams-Moulton 2');
[t,x]=ode23(@oscillator, [0,10], [2;1]);
plot(t,x(:,1),'ks','DisplayName','ode23');
legend('Location', 'NorthEast');
%% b) Linear System - Adams-Moulton 2-step
%
% We solve the stiff linear system with $x_0=(2,1)^T$
% using Milne-Simpson 2, Nystrom 2, Adams-Moulton 2, and ODE23
figure
axis([0 10 -4 4]);
hold on;
xlabel('x_1');
ylabel('x_2');
[t,x]=ms2(@linsystem,0,0.05, [2;1], .001);
plot(t,x(:,1),'ro','DisplayName','Milne-Simpson 2');
[t,x]=ny2(@linsystem,0,0.05, [2;1], .001);
plot(t,x(:,1),'mx','DisplayName','Nystrom 2');
[t,x]=am2(@linsystem,0,0.05, [2;1], .001);
plot(t,x(:,1),'b+','DisplayName','Adams-Moulton 2');
[t,x]=ode23(@linsystem, [0,0.05], [2;1]);
plot(t,x(:,1),'ks','DisplayName','ode23');
ylim([-5 15]);
xlim([-.1 .1]);
legend('Location', 'NorthEast');
%% Convergence Analysis for Milne-Simpson and Nystrom
fprintf('Milne-Simpson 2-step = ');
conv_analysis(@ms2);
fprintf('Nystrom 2-step = ');
conv_analysis(@ny2);