Numerical Methods II (APM3711), UNISA, Complete Exam Pack – assignment solutions, worked examples, and exam preparation material.

1 APM3711/203


2



SEMESTER 1
ASSIGNMENT 02 03
DUE DATE: 20
10JUNE
APRIL202
UNIQUE NR: 229231


QUESTION 1

Solve the boundary–value problem

y 00 + x2 y 0 4xy = 2x3 + 6x2 2; y(0) = 0; y(1) = 2

by using the shooting method. Use the modi…ed Euler method (with only one correction at each step),
and take h = 0:2: Start with an initial slope of y 0 (0) = 1:9 as a …rst attempt and y 0 (0) = 2:1 as a second
attempt. Then interpolate.
Compare the result with the analytical solution y = x4 x2 + 2x:


SOLUTION
To solve the following boundary–value problem by using the shooting method:

y 00 + x2 y 0 4xy = 2x3 + 6x2 2 ((1))

y (0) = 0; y (1) = 2;

we use the modi…ed Euler method with h = 0:2 and compare the result with the analytic solution

y = x4 x2 + 2x:


Set
dy
z= = y0
dx
then
dz
z0 =
= y 00 :
dx
The second–order di¤erential equation (1) can thus be written as a system of two coupled …rst–order
di¤erential equations:

y0 = z ((2))
0 2 3 2
z = 4xy x z + 2x + 6x 2 = g (x; y; z) ((3))

with boundary conditions
y (0) = 0; y (1) = 2:

To solve (2) –(3) we must estimate z (0) = y 0 (0). For a chosen estimate z0 the algorithm is

Predictor Corrector
zp = zi + hzi0 zi+1 = zi + h [zi0 + z 0 (xi + h; yp ; zp )] =2
yp = yi + hyi0 yi+1 = yi + h [yi0 + y 0 (xi + h; yp ; zp )] =2

, 2
APM3711/203

3 APM3711/202


Sequence of calculations

x0 = 0
y0 = 1
z0 = estimate
zi0 = g (xi ; yi ; zi )
zp = zi + hzi0
yi0 = zi
yp = yi + hyi0
zp0 = g (xi + h; yp ; zp )
zi+1 = zi + h zi0 + zp0 =2
yp0 = zp
yi+1 = yi + h yi0 + yp0 =2
xi+1 = xi + h

Note that one could improve the algorithm slightly by using the corrected value zi+1 instead of zp to
calculate yp0 :

The program below uses the estimates

z0 = 1:9 (…rst attempt)
z0 = 2:1 (second attempt)

The third estimate is calculated by using the extrapolation formula:
G2 G1
z0 = G 1 + (D R1 )
R2 R1
where

G1 = 1:9; G2 = 2:1; D = 2:0
R1 = y (1) calculated when z0 = 1:9
R2 = y (1) calculated when z0 = 2:1:

This procedure is repeated, each time using the previous two estimates, until the di¤erence of 2:0 and the
value of y (1) calculated for the last extrapolated estimate z0 is less than a chosen tolerance.

Since (1) is a linear boundary–value problem we expect the problem to be solved after only one extrapola-
tion. This is con…rmed by the results. Note that although the exact values for y (0) and y (1) are obtained
the intermediate values are not very accurate. This is due to the inaccuracy of the modi…ed Euler method
for the large step size h = 0:2.

PROGRAM A3_1(output);
uses printer;

, 3
APM3711/203
4


CONST
xinitial = 0.0;
xfinal = 1.0;
yinitial = 0.0;
yfinal = 2.0;
zinit1 = 1.9;
zinit2 = 2.1;
h = 0.2;
tolerance = 1e-7;
jmax = 10;
VAR
x, y, z : real;
i, imax : integer;
j : 1..jmax;
zi, yf : array[1..jmax] of real;
fst : text;


FUNCTION fy(x,y,z : real) : real;
(* calculates y’ *)
BEGIN
fy := z;
END; {fy}


FUNCTION fz(x,y,z : real) : real;
(* calculates z’ *)
BEGIN
fz := 4*x*y + sqr(x)*(2*x + 6 - z) - 2;
END; {fz}


PROCEDURE Calculate(zinitial : real;
VAR yend : real);
(* solves the differential equation
using modified Euler, using the estimate
zinitial for the initial value of z *)
VAR
x, y, z, yp, zp, fy0, fz0 : real;
i : integer;


BEGIN
x := xinitial;
y := yinitial;
z := zinitial;
writeln(fst);

, 4
5 APM3711/203
APM3711/202


writeln(fst,’ zinitial = ’,zinitial:10:6);
writeln(fst,’ x y z’);
writeln(fst,x:12:6, y:12:6, z:12:6);
FOR i := 1 to imax DO
BEGIN
fy0 := fy(x,y,z);
yp := y + h*fy0;
fz0 := fz(x,y,z);
zp := z + h*fz0;
x := x + h;
y := y + h*(fy0 + fy(x,yp,zp))/2;
z := z + h*(fz0 + fz(x,yp,zp))/2;
writeln(fst,x:12:6, y:12:6, z:12:6);
END; {for i}
yend := y;
writeln(fst,’ error in final y = ’,
(yfinal - yend):10:6);
END; {Calculate}


PROCEDURE PrintExact;
VAR i:integer;
BEGIN
writeln(fst);
writeln(fst,’ Exact solution’);
writeln(fst,’ x y z’);
FOR i := 0 to imax DO
BEGIN
x := xinitial + i*h;
y := sqr(x)*(sqr(x) - 1) + 2*x;
z := x*(4*sqr(x) - 2) + 2;
writeln(fst,x:12:6, y:12:6, z:12:6);
END; {for i}
END; {PrintExact}


BEGIN {program}
assign(fst,’c:napm311nas00_3_1.dat’);
rewrite(fst);
writeln(fst); writeln(fst);
writeln(fst,’ ***** ASSIGNMENT 3, QUESTION 1’,
’ *****’);
imax := round((xfinal - xinitial)/h);
zi[1] := zinit1;