Solution of Assignment 2 (MATH246 Spring 2024)
clearvars; clf; format short g
Problem 1(a)
Solve the homogeneous problem 
: Plugging in 
gives
: Plugging in gives
Hence we need 
. The quadratic formula gives 
, yielding the fundamental set of solutions
. The quadratic formula gives , yielding the fundamental set of solutions 
Find a particular solution: For 
we obtain
we obtain
This gives the linear system 
which has the solution
which has the solution[-3 1;-1 -3]\sym([0;1]) % A\b solves linear system A x = b
We obtain the particular solution 
.
.Hence the general solution of 
is 
:
is :
We determine 
from the initial conditions:
from the initial conditions:
gives 
which has the following solution 
:
:[1 1; 1 -2]\sym([1.1; 2.3])
Hence the solution of the initial value problem is
Problem 1(b)
Note that the coefficient functions are not continuous at 
. Since the initial condition is given at 
we will obtain solutions which exist on the interval 
, i.e., for 
.
. Since the initial condition is given at we will obtain solutions which exist on the interval , i.e., for .Solve the homogeneous problem : Plugging in 
gives
: Plugging in gives
Hence we need 
. We get the solutions 
, yielding the fundamental set of solutions
. We get the solutions , yielding the fundamental set of solutions
Find a particular solution: For 
we obtain
we obtain
which gives 
and 
. Hence we get the particular solution 
.
and . Hence we get the particular solution .Hence the general solution of is :
is :
We determine from the initial conditions:
from the initial conditions: gives 
which has the following solution :
:[1 1; 1 -1]\[ 1-1/sym(15) ; 1-4/sym(15) ]
Hence the solution of the initial value problem is
Problem 2(a)
We already solved the homogeneous ODE in 1(a). For the particular solution we use 
and get
and get
This gives 
and 
, hence 
yielding 
.
and , hence yielding .We determine from the initial conditions:
from the initial conditions: gives 
which has the following solution :
:[1 1; 1 -2]\[5/sym(4); 3/sym(2)]
Hence the solution of the initial value problem is
Problem 2(b)
We rewrite the problem as a first order system of ODEs 
:
:
f = @(t,y) [ y(2) ; t-y(2)+2*y(1) ];
with initial condition 
:
:t0 = 0; y0 = [1;1];
Perform 1 step of improved Euler with 
:
:h = 1/2;
s1 = f(t0,y0)
yE = y0 + h*s1
s2 = f(t0+h,yE)
y1 = y0 + h*(s1+s2)/2
We obtain for 
the approximation 
, hence 
is approximated by 1.625.
the approximation , hence is approximated by 1.625.Problem 2(c)
The code for the function IEulermethod is at the end of this document.
We want to see how the error improves when we double n, and print out
yex = @(t) -t/2-1/4+4/3*exp(t)-1/12*exp(-2*t); % solution from (a)
T = 2;
exact = yex(T);
olderr = NaN;
for n=2.^(3:8)
[ts,ys] = IEulerMethod(f,[0,T],y0,n);
err = ys(end,1)-exact;
fprintf('n=%3g: error=%12e, factor=%g\n',n,err,olderr/err)
olderr = err;
end
The Improved Euler Method has convergence order 2: this means that the error at the final time T satisfies
This upper bound gets divided by 4 when we double the number n. We see that the actual errors decrease in the same way, as the computed factors get closer and closer to 4.
Problem 2(d)
format long g % show 15 digits of results
T = 8;
exact = yex(T);
[ts,ys] = ode45(f,[0,T],y0);
Y = ys(end,1)
err = Y - exact
opt = odeset('RelTol',3e-14,'AbsTol',1e-14);
[ts,ys] = ode45(f,[0,T],y0,opt);
Y = ys(end,1)
err = Y - exact
Problem 3(i)(a)
We rewrite the problem as a first order system of ODEs :
:
f = @(t,y) [ y(2) ; -y(2)+2*y(1) ];
vectorfield(f,-2:.25:2,-2:.25:2); hold on
for a=-2:2
for b=-2:2
for T=[8,-8]
[ts,ys] = ode45(f,[0 T],[a;b]); % find solution for t from 0 to T
plot(ys(:,1),ys(:,2),'b'); % plot trajectory in the phas plane
end
end
end
hold off; axis([-2 2 -2 2]); xlabel('y_1'); ylabel('y_2')
Problem 3(i)(b)
We found in 1(a): the general solution of the homogeneous ODE is
, 
This gives for 
(1) For 
we get 
for 
. Hence 
goes from 0 to ∞, so the solution starts at the origin and moves along a half line going in the direction 
(dashed blue line below).
we get for . Hence goes from 0 to ∞, so the solution starts at the origin and moves along a half line going in the direction (dashed blue line below).(2) For 
we get 
: a half line starting at the origin and going in the direction 
(dotted blue line below).
we get : a half line starting at the origin and going in the direction (dotted blue line below).(3) For 
we get 
for . Hence 
goes from ∞ to 0, so the trajectory is a half line going in the direction 
, the solution moves along this half line towards the origin (dashed green line below)
we get for . Hence goes from ∞ to 0, so the trajectory is a half line going in the direction , the solution moves along this half line towards the origin (dashed green line below)(4) For 
we get 
for (dotted green line below)
we get for (dotted green line below)plot([0 2],[0 2],'b--'); hold on
plot([0 -2],[0 -2],'b:');
plot([0 2],[0 -4],'g--');
plot([0 -2],[0 4],'g:');
vectorfield(f,-2:.5:2,-2:.5:2);
legend('(1)','(2)','(3)','(4)','Location','eastoutside')
axis([-2 2 -2 2]); hold off;
title('Four trajectories')
Problem 3(ii)(a)
We rewrite the problem as a first order system of ODEs :
:
f = @(t,y) [ y(2) ; -1/2*y(1)-9/4*y(2) ];
vectorfield(f,-2:.25:2,-2:.25:2); hold on
for a=-2:1:2
for b=-2:1:2
for T=[8,-8]
[ts,ys] = ode45(f,[0 T],[a;b]);
plot(ys(:,1),ys(:,2),'b')
end
end
end
hold off; axis([-2 2 -2 2]); xlabel('y_1'); ylabel('y_2')
Problem 3(ii)(b)
We get the characteristic equation 
. This gives 
, so the general solution of the homogeneous ODE is
. This gives , so the general solution of the homogeneous ODE is
, 
This gives for
(1) For we get for . Hence goes from ∞ to 0, this is a half line going in the direction (dashed blue line below), the solution is moving toward the origin.
we get for . Hence goes from ∞ to 0, this is a half line going in the direction (dashed blue line below), the solution is moving toward the origin.(2) For we get : a half line starting at the origin and going in the direction (dotted blue line below).
we get : a half line starting at the origin and going in the direction (dotted blue line below).(3) For we get 
for . Hence 
goes from ∞ to 0, so the trajectory is a half line going in the direction 
, the solution moves along this half line towards the origin (dashed green line below)
we get for . Hence goes from ∞ to 0, so the trajectory is a half line going in the direction , the solution moves along this half line towards the origin (dashed green line below)(4) For we get 
for (dotted green line below)
we get for (dotted green line below)plot([0 2],[0 -4],'b--'); hold on
plot([0 -2],[0 4],'b:');
plot([0 8],[0 -2],'g--');
plot([0 -8],[0 2],'g:');
vectorfield(f,-2:.25:2,-2:.25:2);
legend('(1)','(2)','(3)','(4)','Location','eastoutside')
axis([-2 2 -2 2]); hold off;
title('Four trajectories')
Function IEulerMethod for problem 2
function [ts,ys] = IEulerMethod(f,interval,y,n)
% [ts,ys] = EulerMethod(f,[t0,T],y0,n)
% initial value problem y'=f(t,y), y(t0) = y0
% use n steps of Euler method to approximate solution in interval [t0,T]
t = interval(1); T = interval(2);
h = (T-t)/n;
m = size(y,1); % size of vector y
ts = zeros(n+1,1); ys = zeros(n+1,m); % initialize vector for t-values, array for y-values
ts(1) = t; ys(1,:) = y.'; % y.' is transpose, gives row vector
for k=1:n
s1 = f(t,y);
yE = y + h*s1;
s2 = f(t+h,yE);
t = t + h;
y = y + h*(s1+s2)/2;
ts(k+1) = t; ys(k+1,:) = y.'; % store t,y in vector ts, array ys
end
end