COMP5930M - Scientific Computing
Coursework 1
October 19, 2021
Deadline
09:00, Monday 15th November
Task
The three numbered sections of this document describe problems which can be modelled by
nonlinear equations. In each case a proper formulation of the problem as a nonlinear system
is required with appropriate initial data. You will use one or more of the algorithms we have
covered in the course to produce a numerical solution.
MATLAB scripts referred to in this document can be downloaded from the VLE under Learning
Resources / Coursework. MATLAB implementations of some algorithms have been provided
as part of the course but you can implement your solutions in any other language (e.g. Python)
if you prefer.
You should submit your answers as a single PDF document via the VLE before the stated
deadline. Source code submitted as part of your answers should be included in the document.
All functions should include appropriate help information that describe the purpose and use of
that function.
Standard late penalties apply for work submitted after the deadline.
You should properly cite any sources of information used outside of the course material.
Disclaimer
This is intended as an individual piece of work and, while discussion of the work is
encouraged, plagiarism in any form is strictly prohibited.
1. A system of 2 nonlinear equations
[14 marks total]
The following system of nonlinear equations is considered:
x21 − 10x1 + 25 = 5
x1x2 − exp(x2) = −10
(a) Implement in MATLAB code a function that represents this problem as a system of
nonlinear equations in the form F(x) = 0.
(b) Using the Matlab function newtonSys.m and the numerical Jacobian approximation
fdJacobian.m, solve the problem with Newton’s method starting from the initial
guesses x0 = (−2, 3) and x0 = (4,−1) and tolerance tol = 10−12.
State in both cases the solutions obtained xN , the final number of iterations taken
N , and the value of |F(xn)| at each iteration n = 1, 2, . . . , N . Include the values of
all the input parameters used to call newtonSys.m.
(c) Solve the same problem starting from the initial guess x0 = (1, 1) and tolerance
tol = 10−12. What do you observe happening to Newton’s method? State the solution
obtained xN , and the final number of iterations taken N .
Use the function plotContours.m to plot graphically the iterates in the (x1, x2) plane.
Include the plot as a figure in your report.
(d) Implement the gradient descent method:
xk+1 = xk − αkJ(xk)T f(xk), k = 0, 1, . . . , N
as a MATLAB solver following the same design as the newtonSys.m. Include a line-
search step using the step-halving (Armijo) method to select the step size 0 < αk ≤ 1
in each step that guarantees the algorithm only takes descent steps.
Solve the problem using gradient descent with line search starting from the ini-
tial guess x0 = (1, 1), tolerance tol = 10
−6, and maximum number of iterations
kmax = 100 . Use a maximum of 10 step-halving iterations in the line search step.
What do you observe happening to the gradient descent method?
Use the function plotContours.m to plot graphically the iterates in the (x1, x2)
plane. Zoom in on the solution to analyse further. Include the plot as a figure in
your report.
Figure 1: A series of reactors
2. Chemical engineering
[14 marks total]
Figure 1 shows a series of n continuously-stirred reactor vessels, each of volume V, that
are used for a particular chemical reaction. Materials are fed into reactor 1 at rate G. The
chemical of interest has concentration xi−1 at the inlet and concentration ai at the outlet
to reactor i in the series, which feeds into reactor i+ 1.
Given a reaction rate constant k, we can write a nonlinear equation for the steady-state
concentration in the reactor as
βx2i + xi − xi−1 = 0
where β = kVG .
This leads to a system of n nonlinear equations for the series of n reactors.
We assume that values for V , G, a0 and k are known constants.
(a) For the case n = 5 reactors, write out the nonlinear equation system in the form
F(U) = 0 defining F and U. Implement the system in MATLAB code.
(b) Derive the analytical formulas for all elements of the Jacobian matrix for this problem.
(c) Implement the Jacobian computation in MATLAB code.
(d) Produce a numerical solution for the n = 5 case with the following parameters:
V = 1.0, G = 35.0, k = 0.6, a0 = 6.0. using Newton’s method with the exact
Jacobian and tol = 10−12. Explain which initial guess x0 you used. Tabulate the
values of |F| at each iteration.
(e) The steady-state problem above only models the final concentration of the tanks
assuming that an infinite time is allowed to pass. The unsteady version of the problem
is: find the time-dependent concentrations xi(t) s.t.
dxi
dt
= βx2i + xi − xi−1, for t > 0
given an initial state xi(0) = x
0
i for for i = 1, . . . , 5.
Discretise this unsteady problem in time using the backward Euler (implicit) method
and a time step of ∆t. Write the nonlinear equations to be solved in each time step.
How does the Jacobian matrix need to be modified in the unsteady problem compared
to the steady-state problem?
3. Compressible flow in a pipe
[12 marks total]
Weymouth’s equation is a model relating the pressure drop in a pipe to the flow rate of
the compressible gas it is transporting. The model can be expressed as
Q0 = 433.54
T0
P0
(
P 21 − P 22
LσT
) 1
2
d2.667η (1)
where:
Q0 is the gas flow rate
L is the pipe length
d is the pipe diameter
T0 is the standard temperature
T is the actual gas temperature
P0 is the standard pressure
P1 is the upstream (inflow) pressure
P2 is the downstream (outflow) pressure
σ is the specific gravity of gas
η is the efficiency factor
The following data are provided:
Q0 = 2e6, L = 0.1894, T0 = 520, T = 530, P0 = 14.7, P2 = 21.7, σ = 0.7, η = 0.7.
(a) Define a nonlinear equation, F (P1) = 0, for solving the unknown upstream pressure P1
from Weymouth’s equation assuming that all the other variables are known. Provide
a MATLAB implementation of the function.
(b) Compute the analytical derivative F ′(P1) for the nonlinear equation defined in (a).
(c) Assume we are given a pipe of diameter d = 4.026 and solve for the required upstream
pressure P1 using the Newton’s method. All the other variables should be fixed at
the values above. Explain which initial guess x0 you used, the stopping tolerance
chosen, the final iterate xN obtained, and the number of iterations.
(d) Define a function P1(d, L) for solving the upstream pressure P1 for a given pipe
diameter d and length L from Weymouth’s equation. All the other variables should
be fixed at the values above. Provide a MATLAB implementation of the function.
Solve the required upstream pressure P1 the for a range of diameters d ∈ [2, 5], length
L ∈ [0.1, 0.2] using Newton’s method, and plot the surface of the function P1(d, L) in
the two-dimensional coordinates (d, L). Include the figure in your report. Provide a
MATLAB implementation of the function.
Hint: Use the MATLAB functions meshgrid and surf to plot the solution surface in
3D, where the (x, y)-axes are given by the values (d, L) and the z-axis represents the
corresponding value of P1 = (d, L).
Learning objectives
• Formulation of nonlinear equation systems.
• Implementation of nonlinear systems in MATLAB.
• Application of standard solution algorithms for nonlinear systems.
• Convergence properties of Newton’s method and modifications to improve robustness.
• Numerical solution of simple ordinary differential equations
Marking scheme
There are 40 marks in total. This piece of work is worth 20% of the final module grade.
1. A system of 2 nonlinear equations [14 marks total]
(a) [2 marks] nonlinear system and MATLAB code
(b) [4 marks] convergence and solution from two initial guesses
(c) [3 marks] convergence analysis from third initial guess
(d) [5 marks] convergence of gradient descent method
2. Chemical engineering [14 marks total]
(a) [2 marks] nonlinear system and MATLAB code
(b) [3 marks] analytical Jacobian
(c) [2 marks] Jacobian in MATLAB
(d) [3 marks] numerical solution and convergence
(e) [4 marks] unsteady problem analysis
3. Compressible flow in a pipe [12 marks total]
(a) [2 marks] nonlinear system in MATLAB
(b) [2 marks] analytical derivative
(c) [4 marks] solution of the nonlinear problem and MATLAB function
(d) [4 marks] MATLAB function for plotting the surface and the plot
学霸联盟