3
$\begingroup$

I'm solving the classical Black & Scholes (BS) PDE for a European option using finite difference and the implicit scheme. In other words, I'm trying to solve

$\displaystyle\frac{\partial V}{\partial t} + \frac 12\sigma^2\frac{\partial^2V}{\partial x^2} + \mu\frac{\partial V}{\partial x} - rV = 0$

with boundary condition $V(T,x) = \Phi(e^x) = \max\{e^x - K, 0\}$. The drift parameter is as usual $\mu = r - y - \frac 12 \sigma^2$.

The problem is when my $x$ grid is too narrow, my option prices diverge from the BS. Specifically, I'm solving the ATM case with S0 = K = 2775 and if I use [np.log(2675), np.log(2875)] as the boundary points, the solution has a huge error (~400). However, when I use a wider grid [np.log(1775), np.log(3775)], it works like a charm.

I thought the implicit scheme would be numerically stable enough to not exhibit this kind of behavior, so should I assume there is something wrong with my Python code?

My other parameters are r, y, T, sigma = (0.0278, 0.0189, 1, 0.15).

def buildDs(nx):
    D1 = np.diag(np.ones(nx), 1) - np.diag(np.ones(nx), -1)
    D2 = np.diag(np.ones(nx), 1) + np.diag(np.ones(nx), -1) - 2 * np.eye(nx + 1)

    D1[0, 2] = -1
    D1[0, 1] = 4
    D1[0, 0] = -3

    D1[-1, -1] = 3
    D1[-1, -2] = -4
    D1[-1, -3] = 1


    D2[0, 2] = 1
    D2[0, 1] = -2
    D2[0, 0] = 1

    D2[-1, -1] = 1
    D2[-1, -2] = -2
    D2[-1, -3] = 1

    return csc_matrix(D1), csc_matrix(D2)

def pxFD(Phi, r, mu, sigma,
         xs, ts,
         nx = 2000, nt = 10000):

    dx = (xs[1] - xs[0]) / nx
    dt = (ts[1] - ts[0]) / nt

    xs = np.linspace(xs[0], xs[1], num = nx + 1, endpoint = True)
    ts = np.linspace(ts[0], ts[1], num = nt + 1, endpoint = True)

    V = np.zeros((nt + 1, nx + 1))
    I = identity(nx + 1)

    D1, D2 = buildDs(nx)
    L = 1 / 2 * (sigma / dx) ** 2 * D2 + mu / (2 * dx) * D1 - r * I
    P = I - dt * L

    V[-1] = Phi(xs)
    for j in reversed(range(nt)):
        V[j] = spsolve(P, V[j + 1])    
    return V, xs, ts

def pxFDGBM(Phi, r, y, sigma,
            Ss, ts,
            nx = 2000, nt = 10000):

    mu = r - y - 1 / 2 * sigma ** 2
    xs = np.log(Ss)
    Philog = lambda x: Phi(np.exp(x))

    Vs, xs, ts = pxFD(Phi = Philog,
                      r = r,
                      mu = mu,
                      sigma = sigma,
                      xs = xs,
                      ts = ts,
                      nx = nx,
                      nt = nt)
    return Vs, np.exp(xs), ts
$\endgroup$
5
$\begingroup$

The PDE is defined for $x \in ]-\infty, +\infty[$ but the finite difference scheme requires a truncated domain $[x_{\min}, x_{\max}]$, and the choice of $x_{\min}$ and $x_{\max}$ will affect the quality of the result, regardless of the scheme being explicit, implicit, or mixed.

A good rule of thumb is to choose the truncation $[x_{\min}, x_{\max}]$ such that the probability of $X_t < x_{\min}$ or $X_t > x_{\max}$ is very low under the original SDE, because these are states that won't have much impact on the option value.

In the BS case a truncation at $\pm 5$ standard deviation is usually a good choice. In you example that would be $x_{\min} = \ln(S_0) - 5 \sigma \sqrt{T} \approx \ln(1311)$ and $x_{\max} = \ln(S_0) + 5 \sigma \sqrt{T} \approx \ln(5875)$.

Another thing that can improve accuracy is the choice of the boundary condition at $x_{\min}$ and $x_{\max}$.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.