Finite Difference methods

Finite Difference Methods (Forward, Backward and Crank-Nicolson)
The Black-Scholes PDE can be transformed into the heat equation: $$u_t- u_{xx}=0,$$ with known boundary conditions and initial condition $u(x,0)=f(x)$. The heat equation can be approximated using finite difference methods.

Let

• $x\in [XS,XE]$, and $t\in [0,\alpha]$, where $\alpha=\frac{1}{2} r ~T~ \sigma^2$, $XS << 1$ and $XE>>1$.
• i.e. $XS:=-10$ and $XS:=2$, assuming the share price is bounded by the interval $S_t=e^{-10} ~E$ and $S_t=e^2 ~E$, during the lifetime of the option.
• Discretize $x$ over $[XS,XE]$ with $x_i=XS+i~h, ~\forall i \in[0,J]$, with an increment $h=\frac{XE-XS}{J}$, where $J$ is the number of points required over the $x$-space.
• Discretize $t$ over $[0,\alpha]$ with $t_n=n ~\Delta t, ~\forall n \in[0,N]$, with an increment $\Delta t=\frac{\alpha}{N}$, where $N$ is the number of points required over the $t$-space.
• For a put option, we have the initial condition of $$u(x_i,0)=\max \left(e^{0.5(w-1)X_i}-e^{0.5(w+1)X_i},0 \right),$$ and boundary conditions $u(XE,t_n)=0$ and $V(x_1,t_n)=e^{0.5(w-1)x_1+(w+1)^2~t_n/4}$, where $\displaystyle {w=\frac{2r}{\sigma^2}}$.

The equation of the forward method:
$$\frac{u(x_i,t_{n+1})-u(x_i,t_n)}{\Delta t}-~ D^2_x \left[u(x_i,t_n) \right]=0,$$ where $$D^2_x \left[u(x_i,t_n) \right]=\frac{u(x_{i+1},t_{n})-2 u(x_i,t_n)+u(x_{i-1},t_n)}{h^2},$$ $h=x_{i+1}-x_i$ and $\Delta t=t_{n+1}-t_n$. This gives $$u_i^{n+1}=\frac{\Delta t}{h^2} \left (u_{i+1}^n + u_{i-1}^n \right)+ \left(1-2 \frac{\Delta t}{h^2} \right) u_i^n.$$ It can be shown that the forward numerical scheme is only stable when $$(1-2 \lambda)>0 \Longleftrightarrow N>\frac{r T \sigma^2 J^2}{(XE-XS)^2},$$ where $\lambda=\frac{\Delta t}{h^2}$.


The equation of the backward method:
$$\frac{u(x_i,t_{n+1})-u(x_i,t_n)}{\Delta t}-~ D^2_x[u(x_i,t_{n+1})]=0,$$ where $$D^2_x[u(x_i,t_{n+1})]=\frac{u(x_{i+1},t_{n+1})-2 u(x_i,t_{n+1})+u(x_{i-1},t_{n+1})}{h^2},$$ $h=x_{i+1}-x_i$ and $\Delta t=t_{n+1}-t_n$. This gives $$-\frac{\Delta t}{h^2} u_{i+1}^{n+1}+\left(1+2 \frac{\Delta t}{h^2} \right) u_i^{n+1}-\frac{\Delta t}{h^2} u_{i-1}^{n+1}=u_{i}^n.$$ It can be shown that the backward numerical scheme is always stable. Thus, we can find the solution by solving the system of linear equation
$$M_{(J-1)\times (J-1)} ~U^{n+1}_{(J-1)\times 1} = V^{n}_{(J-1) \times 1},$$ where $$\begin{align*} U^{n+1}_{(J-1)\times 1}= \left( \begin{array}{c} u_1^{n+1}\\ u_2^{n+1}\\ u_3^{n+1}\\ \vdots \\ u_{J-2}^{n+1} \\ u_{J-1}^{n+1} \end{array} \right),~ V^{n}_{(J-1)\times 1}= \left( \begin{array}{c} u_1^{n}+\lambda u_0^{n+1}\\ u_2^{n}\\ u_3^{n}\\ \vdots \\ u_{J-2}^{n}\\ u_{J-1}^{n}+\lambda u_{J}^{n+1} \end{array} \right) \end{align*}$$

$$\begin{align*} M_{(J-1)\times (J-1)} = \left( \begin{array}{ccccccc} 1+2\lambda & -\lambda& 0 &0&0&\cdots&0 \\ -\lambda&1+2\lambda & -\lambda& 0 &0&\cdots&0 \\ 0& -\lambda&1+2\lambda & -\lambda& 0 &\cdots&0 \\ 0 &\cdots&\cdots&\cdots&\cdots&\cdots&0 \\ 0 &\cdots&\cdots&\cdots&\cdots&\cdots&0 \\ 0 &\cdots&\cdots&\cdots& -\lambda&1+2\lambda & -\lambda \\ 0 &\cdots&\cdots&\cdots&0& -\lambda&1+2\lambda \\ \end{array} \right).~ \end{align*}$$

Gauss elimination can then be used to find $u_i^{n+1}$ from $u_i^{n}$, $\forall n \in [1,N-1]$. Note that solving for $U^{n+1}_{(J-1)\times 1} = M^{-1}_{(J-1)\times (J-1)} V^{n}_{(J-1) \times 1}$ is computationally expensive when calculating the inverse of a tridiagonal sparse matrix.


The equation of Crank-Nickolson method:
$$\frac{u(x_i,t_{n+1})-u(x_i,t_n)}{\Delta t}-~ D^2_x[\frac{u(x_i,t_{n})+u(x_i,t_{n+1})}{2}]=0,$$ where $$D^2_x[\frac{u(x_i,t_{n})+u(x_i,t_{n+1})}{2}]=\frac{1}{2} \left(D^2_x[u(x_i,t_{n})]+D^2_x[u(x_i,t_{n+1})] \right).$$ It can be shown that Crank-Nickolson numerical scheme is always stable. Thus, we can find the solution by solving the system of linear equation
$$M_{(J-1)\times (J-1)} ~U^{n+1}_{(J-1)\times 1} = V^{n}_{(J-1) \times 1},$$ where $$\begin{align*} U^{n+1}_{(J-1)\times 1}= \left( \begin{array}{c} u_1^{n+1}\\ u_2^{n+1}\\ u_3^{n+1}\\ \vdots \\ u_{J-2}^{n+1} \\ u_{J-1}^{n+1} \end{array} \right),~ V^{n}_{(J-1)\times 1}= \left( \begin{array}{c} f_1^{n}+\lambda u_0^{n+1}\\ f_2^{n}\\ f_3^{n}\\ \vdots \\ f_{J-2}^{n}\\ f_{J-1}^{n}+\lambda u_{J}^{n+1} \end{array} \right), \end{align*}$$

$$\begin{align*} M_{(J-1)\times (J-1)} = \left( \begin{array}{ccccccc} 1+2\lambda & -\lambda& 0 &0&0&\cdots&0 \\ -\lambda&1+2\lambda & -\lambda& 0 &0&\cdots&0 \\ 0& -\lambda&1+2\lambda & -\lambda& 0 &\cdots&0 \\ 0 &\cdots&\cdots&\cdots&\cdots&\cdots&0 \\ 0 &\cdots&\cdots&\cdots&\cdots&\cdots&0 \\ 0 &\cdots&\cdots&\cdots& -\lambda&1+2\lambda & -\lambda \\ 0 &\cdots&\cdots&\cdots&0& -\lambda&1+2\lambda \\ \end{array} \right),~ \end{align*}$$
where $f_i^n:=\lambda u_{i+1}^n+(1-2\lambda) u_i^n+\lambda u_{i-1}^n$. Gauss elimination can then be used to find $u_i^{n+1}$ from $u_i^{n}$, $\forall n \in [1,N-1]$.


Ameican Options.
For American options, we have to include the property of the early exercise of the option. For a Put option, the early exercise has a value of $$g_i^{n+1}=e^{\frac{t_{n+1}}{4}(w+1)^2} \mathrm{max} \left(e^{\frac{w-1}{2}x_i}-e^{\frac{w+1}{2}x_i},0 \right)$$ .
Taking the forward method as an example, we get $$u_i^{n+1}= \mathrm{max} \left( g_i^{n+1}, \frac{\Delta t}{h^2} \left (u_{i+1}^n + u_{i-1}^n \right)+ \left(1-2 \frac{\Delta t}{h^2} \right) u_i^n \right).$$

Forward Finite Differences

Backward Finite Differences

Crank-Nicolson Finite Differences

---

References

1. Burden, R. L., & Faires, J. D. (2005). Numerical analysis. Thomson Brooks/Cole.
2. Haidar, H. (2008). Pricing options for special prices movement,” MSc dissertation, University of Sussex.