Differential Equations. Computer Practical 1

Differential Equations. Computer Practical 1

(Week 2. Explicit and Implicit Euler methods for finding numerical solutions to the initial value problems; first-order ODEs)

February 6, 2014

Aim: To learn how to solve numerically elementary initial value problems using Explicit and Implicit Euler method

Outcomes: written report (must be uploaded to the Blackboard) Due date: Friday, February 21st, by 23.59.

1 Explicit Euler method

Consider the following initial value problem
dy =f(x,y), y(x0)=y0,

dx

where the function f : R2 → R is continuous. Let us suppose that we are interested in finding an approximate solution to this equation.

Let x0,x1,…,xk,…,xn be the set of points from R. For simplicity we consider the case in which points x0,x1,…,xk,…,xn are evenly spaced. Let xk+1 − xk = ∆ be the distance between the points. Implicit Euler method is formally described as follows:

yk+1 =yk +∆·f(xk,yk) (1) The values of yk, k = 1, . . . , n − 1 are approximations of the values of y(xk) of the actual solution at xk.

(N.B. those who have DiPrima’s book are advised to look at Chapter 8.1).

1.1 Implicit Euler method

Implicit Euler method is very similar to the explicit one. Yet, there is one difference: instead of the function f(xk, yk) in (1) the function f(xk+1, yk+1) is used:

yk+1 =yk +∆f(xk+1,yk+1). (2) One needs to solve equation (2) for yk+1 at every k-th iteration of this process.

2 Problems

Consider the following initial value problem
dy = 5y(1−y), y(0) = 0.1, x ∈ [0,2]. (3)

dx

2.1 Explicit Euler

C1. (4 marks) Please implement explicit Euler integration scheme with ∆ = 0.03 (equation (1)) for (4).

C2. (3 marks) Solve (4) analytically.

C3. (3 marks) Let y(xk) be the value of system’s (4) solution at x = xk. Calculate errors of integration ek = yk − y(xk) between the values of yk (numerically estimated solutions) and y(xk) (derived analytically) for ∆ = 0.03, ∆ = 0.02, and ∆ = 0.01. Compute maximal values of integration error for every ∆ (0.03, 0.02, and 0.01) as

emax =max|y(xk)−yk|. k

Please observe and report how the maximal values of integration error change with ∆.

C1. (4 marks) C2. (3 marks) C3. (3 marks)

dx
Please implement implicit Euler integration scheme with ∆ = 0.05 for (4).

Solve (4) analytically

Please calculate (maximal) errors of integration ek = yk − y(xk) between the values of yk (numer- ically estimated solutions) and y(xk) (derived analytically) for ∆ = 0.05, ∆ = 0.1, and ∆ = 0.01. Compute maximal integration errors for every ∆ (0.05, 0.1, and 0.01).

2.2 Implicit Euler

Consider the following initial value problem
dy =−y+(cos(x)+2)y2, y(0)=0.4, x∈[0,4] (4)

3 Examples
3.1 Explicit Euler integration

Implement in MATLAB explicit Euler integration method for the following initial value problem: dy =−y, y(0)=1, x∈[0,π/4]

dx

A possible way to follow is to use basic “while” loop available in MATLAB.

%this is the left boundary of the integration interval [0,pi/4]
x0=0;

%this is the right boundary of the integration interval [0,pi/4]
x1=pi/4;

%definition of Delta (please see formula (1))
Delta=0.01

%here we set up a grid of points in [0,pi/4]
x=x0:Delta:x1;

% y is an array in which the output of the integration is stored;
% y(1) is the first point in this array, and hence it must satisfy
% the initial condition: y(1)=y0=1;

i=1; y(i)=1;

% while loop starts here
while (x(i)<x1)

    y(i+1)=y(i)+Delta*(-y(i));
    i=i+1;

end;
% Vuala, we finished! Let’s plot the estimates y_k of the solution y(x_k) now
plot(x,y);

3.2 Implicit Euler integration

Implement implicit Euler integration method (in MATLAB) for finding numerical solutions of the follow- ing initial value problem:

dy =−y+tan(x), y(0)=1, x∈[0,π/4]. dx

A possible way to follow is to use basic while loop available in MATLAB (you are welcome to exploit other possibilities).

x1=pi/4;
%
x0=0;
%
Delta=0.01
x=x0:Delta:x1;
i=1;

y(i)=1;
while (x(i)<x1)
y(i+1)=(1+Delta)^(-1)*(y(i)+Delta*tan(x(i+1)));
i=i+1;
end;
plot(x,y);

Table 1: Supplementary marking chart