这个作业是用Matlab完成数值计算相关的问题

EMAT20920: Numerical Methods in MATLAB

Question 1: Root-finding

(a) For each of the equations (i)–(iii) below, state

• the number of solutions of the equation in the given range,

• a bracket for each of the solutions,

• the values of the solutions, as accurately as possible (using a built-in MATLAB root-finding

routine).

You should report not just the solutions, but also the code used to find them, and any supporting

material necessary to justify your solutions.

(i) x

2

e

x−1 = tanh(2x − 1) for x ∈ R,

(ii) tan

2

5/2x

3

!

= tan

2x

3

for −2π < x < 2π,

(iii) 50000x

5 + 1000x

4 + 35460x

3 + 286290x

2 − 332910x − 39340

300x

4 − 299

= 499 for x ∈ R,

(10 marks)

(b) The bisection method — which finds a root of a continuous function f(x) in an interval [a, b] on

which f changes sign — can be modified by using a point other than the midpoint m of the interval

as the estimate of the root. One alternative is to instead use the point c where the straight line

joining the points (a, f(a)) and (b, f(b)) crosses the x-axis; that is

c =

af(b) − bf(a)

f(b) − f(a)

. (1)

EMAT20920 2019-20 1 Coursework Assessment

Use this modified bisection method (which uses c in equation (1), instead of the midpoint) to find

the root of f(x) = cos(x) − x in the interval [0, 1], accurate to 13 decimal places. Use the results

to plot a graph that shows the order of convergence of the method, and use the graph to find the

order.

(5 marks)

(c) The iteration

xn =

xn−2f(xn−1) − xn−1f(xn−2)

f(xn−1) − f(xn−2)

(2)

can be used, as n → ∞, to find a root x = limn→∞

xn of the function f, given initial values x0 and x1.

Use equation (2) to find the root of f(x) = cos(x) − x, given initial values x0 = 0 and x1 = 1,

accurate to 13 decimal places. Use the results to plot a graph that shows the order of convergence

of the method, and use the graph to find the order.

(5 marks)

Question 2: Numerical calculus

(a) The expression

b − a

18 ”

5f

a + b

2

−

r

3

5

b − a

2

!

+ 8f

a + b

2

+ 5f

a + b

2

+

r

3

5

b − a

2

!# (3)

can be used to approximate the integral R b

a

f(x) dx.

Use this definition to create a MATLAB function that calculates a numerical approximation to

an integral R β

α

f(x) dx — given an integrand f, region [α, β], and absolute error tolerance —

using a composite version of equation (3); i.e., where the interval [α, β] is divided into equal-width

subintervals, on each of which we use equation (3). Use your function to find the degree of accuracy

and order of the composite rule, explaining your reasons in each case.

(7 marks)

(b) The continuous random variable X has probability density function

fX(x) =

C(k)x

2

log(1 + kx) e−x/k x > 0

0 otherwise

where k > 0 is a constant, and C(k) = 1/

R ∞

0

x

2

log(1 + kx) e−x/k dx. Use MATLAB to plot a

(single) graph of the mean, median, and mode of X, as a function of k.

(6 marks)

EMAT20920 2019-20 2 Coursework Assessment

(c) The expressions

−f(a + 2h) + 8f(a + h) − 8f(a − h) + f(a − 2h)

12h

, (4)

and

−11f(a) + 18f(a + h) − 9f(a + 2h) + 2f(a + 3h)

6h

(5)

can both be used, as h → 0, to approximate the derivative of a function f(x) at x = a.

Use equations (4) and (5) to approximate the derivative of f(x) = x sin(1/x) at x = 3/π, and plot

a graph of the error between the approximate and true value of the derivative in both cases, for a

range of values of h sufficient to illustrate the different sources of error in the formulae. Describe

the results that you find, and explain which is the better formula and why.

(7 marks)

Question 3: Ordinary differential equations

The following set of ODEs is a simple model of infectious disease dynamics in a population:

dS

dt

= −βSI, (6)

dI

dt

= βSI − µI, (7)

dR

dt

= µI. (8)

The variables represent the proportions of the population at time t who are 1. healthy but susceptible

to the disease, S(t), 2. infected with the disease, I(t), and 3. recovered from the disease, R(t). Recovery

from the disease confers immunity, and there is no death (either from the disease, or from natural causes).

The parameter β is the infection rate, and µ is the recovery rate.

(a) Mathematically (without using MATLAB) show that the total population S + I + R is a conserved

quantity, i.e. that it remains fixed for all time. Explain why that means we can set S + I + R = 1,

and why we only need to integrate equations (6) and (7) to solve the SIR system.

(2 marks)

(b) Suppose that at the intial time t = 0, 1% of the population is infected with the disease and everyone

else is susceptible, so S(0) = 0.99 and I(0) = 0.01. Set β = 2 and µ = 1. Use the forward Euler

method with stepsize h = 0.1 to solve the SIR model up to time t = 15. Plot a graph of S(t), I(t),

R(t) against t (in a single graph), with nice-looking labels, distinguishable lines, etc., and explain

what you see.

(5 marks)

EMAT20920 2019-20 3 Coursework Assessment

(c) Write MATLAB code that uses both the forward Euler and the fourth-order Runge-Kutta methods

to solve the SIR model, with the same parameter values and initial conditions as in part (b) above,

for a given stepsize h. Use it to plot a single graph of the final R(t) value against stepsize h, for

different values of h between 100 and 10−3

. Explain what you see.

(5 marks)

(d) Use an in-built MATLAB ODE solver to solve the SIR model, with the same parameter values and

initial conditions as in part (b) above, and find the time for which I(t) is a maximum.

(3 marks)

(e) Fix the recovery rate µ = 1. Solve the SIR model, using an in-built MATLAB ODE solver, until

the dynamics reaches a steady state, for various values of β, and plot a graph of the final R(t) value

against β. Explain the graph you have obtained.

## 发表评论