MTHM003 ANALYSIS AND COMPUTATION FOR FINANCE

Overview and Submission information

General guidelines. Please write down your answers clearly quoting numerical values to 4 significant

figures, unless otherwise stated. Full marks achieved in this assignment will contribute 15% of

the final module mark. Usual university penalties apply for late submission unless mitigation is

approved. Feedback will be written on your work, and general comments will also appear on

ELE.

Submission format. Please submit a single PDF file summarising your answers to the questions

(e.g., computed quantities, graphs, plots, etc.) and including all supporting MATLAB code via

eBART submission. Any queries, please ask at the front desk for help. Analytical calculations

and discussions (when asked for) can be done by hand-writing.

Allocation of marks. When solving any particular problem, marks will be awarded for clearly

commented (MATLAB) codes, and carefully reasoned mathematical solutions. Any diagrams,

graphs, etc., should be clear to understand, with good comment or clear labelling. For questions

requiring numerical simulation/output, marks will primarily be awarded for demonstrating use of

MATLAB: little credit will be given for other solution schemes.

Academic Misconduct. All work submitted including functions and scripts must be your

own and not copied from other students or other sources (e.g., web pages). Any evidence

of possible plagiarism will be investigated following university procedures. You are welcome to

discuss ideas and possible methods, of course, with each other and with lecturers, etc., but the

work you hand in must be your own.

Questions

1. Chaotic dynamics in a simple system of ODEs.

(a) (10 marks) Write a MATLAB function RK4.m to solve a system of ODEs using the 4th

order Runge–Kutta method. Recall that to obtain the solution at time tj+1 from the state

yj at time tj, this iterative scheme requires you to calculate four variables

k1 = f(yj, tj), k2 = f

(

yj +

h

2

k1, tj +

h

2

)

,

k3 = f

(

yj +

h

2

k2, tj +

h

2

)

,

k4 = f(yj + hk3, tj + h),

and time step using

yj+1 = yj +

h

6

(k1 + 2k2 + 2k3 + k4).

Your implementation should use the same inputs and outputs as the explicit Euler and

Euler–Heun methods given in the lecture notes.

该文档是极速PDF编辑器生成，

如果想去掉该提示,请访问并下载：

http://www.jisupdfeditor.com/

(b) (5 marks) The Ro¨ssler equations are a system of ODEs known to produce chaotic dynamics.

The equations are related to the Lorenz equations describing atmospheric convection. For

a column state vector X = (x, y, z)T, time t and parameters a, b and c the equations are

defined as follows:

x˙ = −y − z,

y˙ = x+ ay,

z˙ = b+ z(x− c).

For c2 > 4ab the equations have two fixed points X± = (x±, y±, z±)

T given by

x± = −ay±,

y± = −z±,

z± =

c±√c2 − 4ab

2a

.

(i) Fix the parameters to a = −0.1, b = 0.2 and c = 5.7. Starting from the initial

condition Xˆ = (0, 0, 0.25)T numerically solve the Ro¨ssler equations using your RK4

implementation for a time interval t ∈ [0, 200] with a step size of h = 0.005 (use the

same step size for all questions). Plot each of the state variables against time.

(ii) You should find that the trajectory converges to one of the fixed points, give the value

of this fixed point to 4 decimal places as calculated in MATLAB through simulation;

check your computed value against the above formulae.

(c) (10 marks)

(i) For the parameter values a = 0.1, b = 0.2 and c = 5.7, the same initial condition as in

(b), and a time interval t ∈ [0, 800], numerically solve the Ro¨ssler equations. Discard

any transient behaviour (say t < 300) and plot the trajectory in (x, y, z) coordinates.

You should find a closed curve (a periodic orbit).

(ii) Taking the Poincare´ section P = {x, z : y = 0, x > 0} (the positive-x part of the plane

y = 0), find all points at which the trajectory passes through P . Plot the trajectory in

(x, z) coordinates, marking any points along the trajectory that cross through P . At

a = 0.1 these points should all coincide exactly. Give the value of the periodic orbit’s

minimal period T to 4 decimal places.

(iii) Plot each of the state variables against time, marking the points where the trajectory

crosses through P .

(d) (5 marks) Repeat part (c) for a = 0.13, a = 0.15, a = 0.1525 and a = 0.2. Comment on

how the dynamics change as a increases.

Hints:

1. The Ro¨ssler equations can be defined as an anonymous function in MATLAB to work with

your function RK4.m as follows:

[email protected](t,x)[...] % function defining right-hand side

2. You can check the equations are defined correctly by calling ross(0,Xpl) or ross(0,Xmi)

with Xpl and Xmi being the fixed points given above, which should return the column vector

(0, 0, 0)T.

3. In part (c) the points should coincide exactly. You can improve the accuracy by interpolating

from points along the trajectory either side of P to find more accurate values of x and z

when crossing through at y = 0.

[30]

2. Fokker–Planck equation.

A simple model for asset prices St for t ≥ 0 is given by the following stochastic differential equation

(SDE):

dSt = µStdt+ σStdWt, S0 = s0,

where Wt is a standard Wiener process, µ is the drift coefficient or expected return and σ is

the volatility. You are provided with a MATLAB script sde example cw2.m that generates 500

realisations of this SDE and plots the distribution of the trajectories at final time T = 2 as a

histogram for the parameter values µ = 0.2, σ = 0.25 and initial conditions s0 drawn from a

normal distribution N(ν, λ2) centered at ν = 0.5 with standard deviation λ = 0.05.

The probability density p(t, s) of trajectories of the above SDE over the spatial variable asset price

s ∈ (0,∞) satisfies the Fokker–Planck equation (FPE)

∂tp(t, s) = −∂s[µsp(t, s)] + 1

2

∂ss[σ

2s2p(t, s)],

a partial differential equation. The FPE can be solved numerically using a forward in time, central

in space (FTCS) finite difference method as described below:

Discretise the time interval t ∈ [0, 2] as tj = t0 + j∆t for j = 0, 1, . . . , nt, where the final

time point is T = t0 + nt∆t = 2.

For s ∈ [0, 3] let sk = (k − 1)∆s with k = 1, 2, . . . , ns where ns = 3/∆s + 1. While s can

take any value in (0,∞), here we truncate the solution assuming p(t, s) = 0 for s ≥ 3 for

the given time interval.

Notation: pjk ≈ p(tj, sk) for j = 0, 1, . . . , nt and k = 1, 2, . . . , ns

Initial condition (IC): p0k = p(0, sk) = (1/

√

2piλ) exp[−(sk − ν)2/(2λ2)] for k = 1, 2, . . . , ns,

that is, the normal distribution N(ν, λ2)

Dirichlet boundary conditions (BCs): pj1 = p(tj, s1) = 0 and pjns = p(tj, sns) = 0 for

j = 0, 1, . . . , nt

Recursively define (loop over j and k):

p(j+1)k = pjk − µ∆t

2∆s

[sk+1pj(k+1) − sk−1pj(k−1)] + σ

2∆t

2∆2s

[s2k+1pj(k+1) − 2s2kpjk + s2k−1pj(k−1)]

(a) (20 marks) Write MATLAB code for the above FTCS method and simulate the FPE with

∆t = 0.0001, ∆s = 0.01, parameter values, initial condition and boundary conditions as

given above.

(b) (10 marks) Demonstrate that your implementation works and is accurate by producing the

following plots:

该文档是极速PDF编辑器生成，

如果想去掉该提示,请访问并下载：

http://www.jisupdfeditor.com/

– Plot your solution p(t, s) as a surface or as a 2D colour map (use, for example, surf or

imagesc).

– Plot the distributions p(t, s) at the initial and final time points, and compare these

with the distribution of St at the final time T = 2 as produced by the given script

sde example cw2.m.

Hints:

Your solution to the FPE will be a matrix p of size (nt+1)×ns. Before using the recursive equation

given above in a for-loop you can pre-allocate some values in p using IC for the first row, and

using BCs for the first and final columns. [30]

3. The diffusion equation and the Black–Scholes equation.

(a) (15 marks) Use the MATLAB PDE solver pdepe to solve the diffusion equation

∂u

∂t

=

∂2u

∂x2

on the interval x ∈ [0, 1] with the initial condition

u(0, x) = u0(x) = − sin 5pix

2

subject to the boundary conditions

u(t, 0) = 0

and

∂u

∂x

(t, 1) = 0.

Use 41 equally spaced mesh points in the x-direction and integrate between t = 0 and t = 0.1

using 21 equally spaced time steps, and follow the method given in the example in section 12

of the lecture notes (edit the functions pdefun, pdeic and pdebc for use with the MATLAB

function pdepe). For output, plot u(t, x) at the times t = 0, t = 0.025, t = 0.05, t = 0.075

and t = 0.1 on one graph.

(b) (25 marks)

(i) Solve the transformed Black–Scholes equation in the European put case

∂u

∂τ

=

∂2u

∂x2

with initial condition

u(0, x) = u0(x) = max

{

exp

[

1

2

(k − 1)x

]

− exp

[

1

2

(k + 1)x

]

, 0

}

and boundary conditions

u(τ, x) = exp

[

1

2

(k − 1)x+ 1

4

(k − 1)2τ

]

as x→ −∞

and

u(τ, x) = 0 as x→ +∞

using the MATLAB solver pdepe. The constant k is given below. We cannot apply

the boundary conditions at x = ±∞, however, doing so at x = ±xb works well with

the choice xb = 3, for example. Plot your results (using, for example, the MATLAB

functions surf or imagesc).

To recover the financial variables from the transformed variables, the following trans-

formations are used:

S = Eex, t = T − 2τ

σ2

, P = Eu exp

[

−1

2

(k − 1)x− 1

4

(k + 1)2τ

]

.

(ii) Solve the transformed equations for the case E = 10, r = 0.1 and σ = 0.45 with

T = 1/3, where k = 2r/σ2. Explain why the scaled time τ runs from 0 to 0.03375.

Your output should be P at integer values of S from S = 0 to S = 16 (noting that

P0 = Ee

−kτ can be determined analytically; explain why this is the case). Some of

these values are given in the table below for comparison.

S P

0.00 9.6722

2.00 7.6722

4.00 5.6723

6.00 3.6977

(iii) Experiment by increasing the number of mesh points in the x-direction, and by in-

creasing the value of xb, which gives the points where the boundary conditions are

applied.

It is important to annotate your code with comments, explaining what each of the function

m-files does and how these function m-files refer to the general equation solved by pdepe.

Annotate also the driving m-file script, explaining in particular the role of xb and the use of

pdeval in extracting P in terms of S for the output.

[40]

欢迎咨询51作业君