# Accuracy, Robustness, and Efficiency of the Linear Boundary Condition for the Black-Scholes Equations

Authors: Darae Jeong; Seungsuk Seo; Hyeongseok Hwang; Dongsun Lee; Yongho Choi; Junseok Kim
Journal: Discrete Dynamics in Nature and Society
 (2015)

---

## Abstract

We briefly review and investigate the performance of various boundary conditions such as Dirichlet, Neumann, linear, and partial differential equation boundary conditions for the numerical solutions of the Black-Scholes partial differential equation. We use a finite difference method to numerically solve the equation. To show the efficiency of the given boundary condition, several numerical examples are presented. In numerical test, we investigate the effect of the domain sizes and compare the effect of various boundary conditions with pointwise error and root mean square error. Numerical results show that linear boundary condition is accurate and efficient among the other boundary conditions.


---

## Body

## 1. Introduction

We briefly review and perform a comparison study for Dirichlet, Neumann, linear, and partial differential equation (PDE) boundary conditions (BCs) for the Black-Scholes (BS) partial differential equations. LetS and t denote the price of the underlying asset and time, respectively. Then, the value of the option u
(
S
,
t
) is governed by the BS equation [1]:
(1)
∂
u
S
,
t
∂
t
+
1
2
σ
S
2
∂
2
u
S
,
t
∂
S
2
+
r
S
∂
u
S
,
t
∂
S
=
r
u
(
S
,
t
)
,
for
S
,
t
∈
R
+
×
0
,
T
,
where σ is a constant volatility of the asset and r
>
0 is a constant riskless interest rate. The final condition is the payoff function Λ
(
S
) at expiry T:
(2)
u
(
S
,
T
)
=
Λ
(
S
)
.
It is difficult to find the analytic solutions of (1) and (2) for exotic options. Therefore we need to use a numerical approximation. To obtain an approximation of the option value, one can compute a solution of the BS equations (1) and (2) using a finite difference method (FDM) [2–8], finite element method [9–11], finite volume method [12–14], a fast Fourier transform [15–17], and also their optimal BC [18]. Here, we consider FDM approach to solve the BS PDE. See [1, 19–24] for more details for FDM applied to option pricing. To solve the BS PDE using the FDM, we need to truncate the infinite domain to a finite domain and have to use an artificial BC. To reduce large errors in the numerical solution due to this approximation of the BCs, the truncated domain must be large enough. Thus the purpose of this paper is to investigate the effects of several BCs on the numerical solutions for the BS PDE (1) and (2).The outline of the paper is as follows. In Section2, we formulate the generalized version of the BS PDE. And we explain five cases which are Dirichlet I, Dirichlet II, Neumann, linear, and PDE BCs. Section 3 presents the results of the numerical experiments and the conclusions are drawn in Section 4.


## 2. Numerical Solution

Letx
=
S be the value of the underlying asset price and let τ
=
T
-
t be the time to expiry; then (1) becomes a more natural initial value problem
(3)
∂
u
(
x
,
τ
)
∂
τ
=
1
2
σ
x
2
∂
2
u
x
,
τ
∂
x
2
+
r
x
∂
u
x
,
τ
∂
x
-
r
u
(
x
,
τ
)
,
for
x
,
τ
∈
Ω
×
0
,
T
with an initial condition u
(
x
,
0
)
=
Λ
(
x
) for x
∈
Ω
=
(
0
,
L
). Here we truncate the infinite domain into a finite domain since the infinite domain cannot be represented in the computer [23].


### 2.1. Discretization with Finite Difference Method

We apply the FDM for solving (3) numerically. Let us first discretize the given computational domain Ω
=
[
0
,
L
] as a uniform grid with a spatial step size h
=
L
/
N
x and a temporal step size Δ
τ
=
T
/
N
t, where N
x is the number of subintervals (see Figure 1) and N
t is the number of time steps.Figure 1
Uniform grid with a spatial step sizeh.Let us denote the numerical approximation of the solution as(4)
u
i
n
≡
u
(
x
i
,
τ
n
)
=
u
(
i
h
,
n
Δ
τ
)
,
where i
=
0,1
,
…
,
N
x and n
=
0,1
,
…
,
N
t. By applying the fully implicit-in-time and space-centered difference scheme to (3), we have
(5)
u
i
n
+
1
-
u
i
n
Δ
τ
=
σ
2
x
i
2
2
u
i
-
1
n
+
1
-
2
u
i
n
+
1
+
u
i
+
1
n
+
1
h
2
+
r
x
i
u
i
+
1
n
+
1
-
u
i
-
1
n
+
1
2
h
-
r
u
i
n
+
1
.
We can rewrite (5) by
(6)
α
i
u
i
-
1
n
+
1
+
β
i
u
i
n
+
1
+
γ
i
u
i
+
1
n
+
1
=
b
i
,
where α
i
=
r
x
i
/
2
h
-
σ
2
x
i
2
/
2
h
2, β
i
=
1
/
Δ
τ
+
σ
2
x
i
2
/
h
2
+
r, γ
i
=
-
r
x
i
/
2
h
-
σ
2
x
i
2
/
2
h
2, and b
i
=
u
i
n
/
Δ
τ.In order to solve the linear system (6), we need to know u
0
n and u
N
x
n for all n
=
0
,
…
,
N
t. At x
=
0, we simply set u
0
n
=
0. Next, we present five different boundary conditions for specifying the values of u
N
x
n.


### 2.2. Boundary Conditions

In solving the BS equation numerically, there are many BCs. In this section, we introduce five BCs, which are used to solve the BS equation. We focus on two options: European call option and cash-or-nothing option. The payoff functions are given by(7)
u
x
,
0
=
max
⁡
x
-
K
,
0
,
(8)
u
(
x
,
0
)
=
C
if
x
>
K
0
otherwise
,
for the European call option and the cash-or-nothing option, respectively. Here K is the strike price and C denotes the return value at expiration if the option is in-the-money. We have closed-form solutions for these options. For the European call option, the closed-form solution of the BS equation is
(9)
u
x
,
τ
=
x
N
d
1
-
K
e
-
r
τ
N
d
2
,
∀
x
∈
0
,
L
,
∀
τ
∈
0
,
T
,
(10)
d
1
=
ln
⁡
x
/
K
+
r
+
(
1
/
2
)
σ
2
τ
σ
τ
,
d
2
=
d
1
-
σ
τ
,
where N
(
d
)
=
(
1
/
2
π
)
∫
-
∞
d
‍
exp
⁡
-
x
2
/
2
d
x is the cumulative distribution function for the standard normal distribution [25]. For the cash-or-nothing option, the closed form solution is
(11)
u
(
x
,
τ
)
=
C
e
-
r
τ
N
(
d
2
)
.


#### 2.2.1. Dirichlet I Boundary Condition

From (9), we can observe that N
(
d
1
) and N
(
d
2
) are close to one when x is large enough. Therefore, Dirichlet BC for the European call option is defined by
(12)
u
(
L
,
τ
)
=
L
-
K
e
-
r
τ
for a sufficiently large L. By (12), we use u
N
x
n
+
1
=
x
N
x
-
K
e
-
r
(
n
+
1
)
Δ
τ. Using this, we can rewrite (6) as
(13)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
γ
N
x
-
1
(
L
-
K
e
-
r
(
n
+
1
)
Δ
τ
)
.
Likewise, we have u
N
x
n
+
1
=
C
e
-
r
(
n
+
1
)
Δ
τ for the cash-or-nothing option.


#### 2.2.2. Dirichlet II Boundary Condition

The other Dirichlet boundary condition is simply setting the boundary values to be fixed all the time with the payoff value. For the max call option, we setu
N
x
n
+
1
=
L
-
K:
(14)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
γ
N
x
-
1
(
L
-
K
)
.
For the cash-or-nothing option, we set u
N
x
n
+
1
=
C. This boundary condition was used in [26].


#### 2.2.3. Neumann Boundary Condition

To specify values for the derivative of the solution at the boundary of the spatial domain, we have the following equation from (9):
(15)
∂
u
(
x
,
τ
)
∂
x
=
N
(
d
1
)
.
Therefore, for a sufficiently large L, we assume
(16)
∂
u
(
L
,
τ
)
∂
x
=
1
,
which we call by Neumann BC for European call option. Equation (16) can be discretized as (
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
)
/
h
=
1; that is, u
N
x
n
+
1
=
u
N
x
-
1
n
+
1
+
h. Therefore, we have
(17)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
+
γ
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
h
γ
N
x
-
1
.
Similarly, applying BC in terms of the first derivative of x for the cash-or-nothing option, we can get
(18)
∂
u
(
x
,
τ
)
∂
x
=
C
e
-
(
r
τ
+
d
2
2
/
2
)
x
σ
2
π
τ
from (11). Therefore, we can assume
(19)
∂
u
(
L
,
τ
)
∂
x
=
0
for a sufficiently large L. This represents Neumann BC for cash-or-nothing option. We now replace (19) with the one-sided derivative (
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
)
/
h
=
0; that is, u
N
x
n
+
1
=
u
N
x
-
1
n
+
1.


#### 2.2.4. Linear Boundary Condition

Linear boundary condition assumes that the second derivative of the option value with respect to the underlying asset pricex vanishes to zero for the large value of the asset price. To demonstrate this, we consider (9) at the right end of the domain. The second derivative value is given by
(20)
∂
2
u
(
x
,
τ
)
∂
x
2
=
e
-
d
1
2
/
2
σ
2
π
τ
x
.
If the asset price approaches to the large value L, we assume
(21)
∂
2
u
(
L
,
τ
)
∂
x
2
=
0
,
which we call by linear BC for European call option. Equation (21) can be discretized as (
u
N
x
-
1
n
+
1
-
2
u
N
x
n
+
1
+
u
N
x
+
1
n
+
1
)
/
h
2
=
0. By this relation, we obtain the boundary value as u
N
x
+
1
n
+
1
=
2
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1. By substituting this in (3), we get
(22)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
1
β
N
x
-
1
γ
N
x
-
1
0
⋯
0
α
N
x
-
γ
N
x
β
N
x
+
2
γ
N
x
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
1
n
+
1
u
N
x
n
+
1
=
b
1
b
2
⋮
b
N
x
-
1
b
N
x
.
Similarly, for the cash-or-nothing option, we have
(23)
∂
2
u
(
x
,
τ
)
∂
x
2
=
-
C
d
1
e
-
(
r
τ
+
d
2
2
/
2
)
σ
2
τ
2
π
x
2
from (11). Therefore, for a sufficiently large L, we have
(24)
∂
2
u
(
L
,
τ
)
∂
x
2
=
0
.
In this case, we can obtain the same BC as European call option.


#### 2.2.5. PDE Boundary Condition

Next we consider the PDE BC [23, 27]. Using the BS equation itself, we derive the BC. We use one-sided discretization:
(25)
u
N
x
n
+
1
-
u
N
x
n
Δ
τ
=
σ
2
x
N
x
2
2
u
N
x
n
+
1
-
2
u
N
x
-
1
n
+
1
+
u
N
x
-
2
n
+
1
h
2
+
r
x
N
x
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
h
-
r
u
N
x
n
+
1
,
where the drift and volatility terms are discretized by using one-sided derivatives. The linear system of N
x equations can be written in the following matrix form:
(26)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
1
β
N
x
-
1
γ
N
x
-
1
0
⋯
α
N
x
β
N
x
γ
N
x
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
1
n
+
1
u
N
x
n
+
1
=
b
1
b
2
⋮
b
N
x
-
1
b
N
x
,
where α
N
x
=
-
σ
2
x
N
x
2
/
2
h
2, β
N
x
=
σ
2
x
N
x
2
/
h
2
+
r
x
N
x
/
h, γ
N
x
=
1
/
Δ
τ
-
σ
2
x
N
x
2
/
2
h
2
-
r
x
N
x
/
h
+
r, and b
N
x
=
u
N
x
n
/
Δ
τ.


## 2.1. Discretization with Finite Difference Method

We apply the FDM for solving (3) numerically. Let us first discretize the given computational domain Ω
=
[
0
,
L
] as a uniform grid with a spatial step size h
=
L
/
N
x and a temporal step size Δ
τ
=
T
/
N
t, where N
x is the number of subintervals (see Figure 1) and N
t is the number of time steps.Figure 1
Uniform grid with a spatial step sizeh.Let us denote the numerical approximation of the solution as(4)
u
i
n
≡
u
(
x
i
,
τ
n
)
=
u
(
i
h
,
n
Δ
τ
)
,
where i
=
0,1
,
…
,
N
x and n
=
0,1
,
…
,
N
t. By applying the fully implicit-in-time and space-centered difference scheme to (3), we have
(5)
u
i
n
+
1
-
u
i
n
Δ
τ
=
σ
2
x
i
2
2
u
i
-
1
n
+
1
-
2
u
i
n
+
1
+
u
i
+
1
n
+
1
h
2
+
r
x
i
u
i
+
1
n
+
1
-
u
i
-
1
n
+
1
2
h
-
r
u
i
n
+
1
.
We can rewrite (5) by
(6)
α
i
u
i
-
1
n
+
1
+
β
i
u
i
n
+
1
+
γ
i
u
i
+
1
n
+
1
=
b
i
,
where α
i
=
r
x
i
/
2
h
-
σ
2
x
i
2
/
2
h
2, β
i
=
1
/
Δ
τ
+
σ
2
x
i
2
/
h
2
+
r, γ
i
=
-
r
x
i
/
2
h
-
σ
2
x
i
2
/
2
h
2, and b
i
=
u
i
n
/
Δ
τ.In order to solve the linear system (6), we need to know u
0
n and u
N
x
n for all n
=
0
,
…
,
N
t. At x
=
0, we simply set u
0
n
=
0. Next, we present five different boundary conditions for specifying the values of u
N
x
n.


## 2.2. Boundary Conditions

In solving the BS equation numerically, there are many BCs. In this section, we introduce five BCs, which are used to solve the BS equation. We focus on two options: European call option and cash-or-nothing option. The payoff functions are given by(7)
u
x
,
0
=
max
⁡
x
-
K
,
0
,
(8)
u
(
x
,
0
)
=
C
if
x
>
K
0
otherwise
,
for the European call option and the cash-or-nothing option, respectively. Here K is the strike price and C denotes the return value at expiration if the option is in-the-money. We have closed-form solutions for these options. For the European call option, the closed-form solution of the BS equation is
(9)
u
x
,
τ
=
x
N
d
1
-
K
e
-
r
τ
N
d
2
,
∀
x
∈
0
,
L
,
∀
τ
∈
0
,
T
,
(10)
d
1
=
ln
⁡
x
/
K
+
r
+
(
1
/
2
)
σ
2
τ
σ
τ
,
d
2
=
d
1
-
σ
τ
,
where N
(
d
)
=
(
1
/
2
π
)
∫
-
∞
d
‍
exp
⁡
-
x
2
/
2
d
x is the cumulative distribution function for the standard normal distribution [25]. For the cash-or-nothing option, the closed form solution is
(11)
u
(
x
,
τ
)
=
C
e
-
r
τ
N
(
d
2
)
.


### 2.2.1. Dirichlet I Boundary Condition

From (9), we can observe that N
(
d
1
) and N
(
d
2
) are close to one when x is large enough. Therefore, Dirichlet BC for the European call option is defined by
(12)
u
(
L
,
τ
)
=
L
-
K
e
-
r
τ
for a sufficiently large L. By (12), we use u
N
x
n
+
1
=
x
N
x
-
K
e
-
r
(
n
+
1
)
Δ
τ. Using this, we can rewrite (6) as
(13)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
γ
N
x
-
1
(
L
-
K
e
-
r
(
n
+
1
)
Δ
τ
)
.
Likewise, we have u
N
x
n
+
1
=
C
e
-
r
(
n
+
1
)
Δ
τ for the cash-or-nothing option.


### 2.2.2. Dirichlet II Boundary Condition

The other Dirichlet boundary condition is simply setting the boundary values to be fixed all the time with the payoff value. For the max call option, we setu
N
x
n
+
1
=
L
-
K:
(14)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
γ
N
x
-
1
(
L
-
K
)
.
For the cash-or-nothing option, we set u
N
x
n
+
1
=
C. This boundary condition was used in [26].


### 2.2.3. Neumann Boundary Condition

To specify values for the derivative of the solution at the boundary of the spatial domain, we have the following equation from (9):
(15)
∂
u
(
x
,
τ
)
∂
x
=
N
(
d
1
)
.
Therefore, for a sufficiently large L, we assume
(16)
∂
u
(
L
,
τ
)
∂
x
=
1
,
which we call by Neumann BC for European call option. Equation (16) can be discretized as (
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
)
/
h
=
1; that is, u
N
x
n
+
1
=
u
N
x
-
1
n
+
1
+
h. Therefore, we have
(17)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
+
γ
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
h
γ
N
x
-
1
.
Similarly, applying BC in terms of the first derivative of x for the cash-or-nothing option, we can get
(18)
∂
u
(
x
,
τ
)
∂
x
=
C
e
-
(
r
τ
+
d
2
2
/
2
)
x
σ
2
π
τ
from (11). Therefore, we can assume
(19)
∂
u
(
L
,
τ
)
∂
x
=
0
for a sufficiently large L. This represents Neumann BC for cash-or-nothing option. We now replace (19) with the one-sided derivative (
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
)
/
h
=
0; that is, u
N
x
n
+
1
=
u
N
x
-
1
n
+
1.


### 2.2.4. Linear Boundary Condition

Linear boundary condition assumes that the second derivative of the option value with respect to the underlying asset pricex vanishes to zero for the large value of the asset price. To demonstrate this, we consider (9) at the right end of the domain. The second derivative value is given by
(20)
∂
2
u
(
x
,
τ
)
∂
x
2
=
e
-
d
1
2
/
2
σ
2
π
τ
x
.
If the asset price approaches to the large value L, we assume
(21)
∂
2
u
(
L
,
τ
)
∂
x
2
=
0
,
which we call by linear BC for European call option. Equation (21) can be discretized as (
u
N
x
-
1
n
+
1
-
2
u
N
x
n
+
1
+
u
N
x
+
1
n
+
1
)
/
h
2
=
0. By this relation, we obtain the boundary value as u
N
x
+
1
n
+
1
=
2
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1. By substituting this in (3), we get
(22)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
1
β
N
x
-
1
γ
N
x
-
1
0
⋯
0
α
N
x
-
γ
N
x
β
N
x
+
2
γ
N
x
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
1
n
+
1
u
N
x
n
+
1
=
b
1
b
2
⋮
b
N
x
-
1
b
N
x
.
Similarly, for the cash-or-nothing option, we have
(23)
∂
2
u
(
x
,
τ
)
∂
x
2
=
-
C
d
1
e
-
(
r
τ
+
d
2
2
/
2
)
σ
2
τ
2
π
x
2
from (11). Therefore, for a sufficiently large L, we have
(24)
∂
2
u
(
L
,
τ
)
∂
x
2
=
0
.
In this case, we can obtain the same BC as European call option.


### 2.2.5. PDE Boundary Condition

Next we consider the PDE BC [23, 27]. Using the BS equation itself, we derive the BC. We use one-sided discretization:
(25)
u
N
x
n
+
1
-
u
N
x
n
Δ
τ
=
σ
2
x
N
x
2
2
u
N
x
n
+
1
-
2
u
N
x
-
1
n
+
1
+
u
N
x
-
2
n
+
1
h
2
+
r
x
N
x
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
h
-
r
u
N
x
n
+
1
,
where the drift and volatility terms are discretized by using one-sided derivatives. The linear system of N
x equations can be written in the following matrix form:
(26)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
1
β
N
x
-
1
γ
N
x
-
1
0
⋯
α
N
x
β
N
x
γ
N
x
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
1
n
+
1
u
N
x
n
+
1
=
b
1
b
2
⋮
b
N
x
-
1
b
N
x
,
where α
N
x
=
-
σ
2
x
N
x
2
/
2
h
2, β
N
x
=
σ
2
x
N
x
2
/
h
2
+
r
x
N
x
/
h, γ
N
x
=
1
/
Δ
τ
-
σ
2
x
N
x
2
/
2
h
2
-
r
x
N
x
/
h
+
r, and b
N
x
=
u
N
x
n
/
Δ
τ.


## 2.2.1. Dirichlet I Boundary Condition

From (9), we can observe that N
(
d
1
) and N
(
d
2
) are close to one when x is large enough. Therefore, Dirichlet BC for the European call option is defined by
(12)
u
(
L
,
τ
)
=
L
-
K
e
-
r
τ
for a sufficiently large L. By (12), we use u
N
x
n
+
1
=
x
N
x
-
K
e
-
r
(
n
+
1
)
Δ
τ. Using this, we can rewrite (6) as
(13)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
γ
N
x
-
1
(
L
-
K
e
-
r
(
n
+
1
)
Δ
τ
)
.
Likewise, we have u
N
x
n
+
1
=
C
e
-
r
(
n
+
1
)
Δ
τ for the cash-or-nothing option.


## 2.2.2. Dirichlet II Boundary Condition

The other Dirichlet boundary condition is simply setting the boundary values to be fixed all the time with the payoff value. For the max call option, we setu
N
x
n
+
1
=
L
-
K:
(14)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
γ
N
x
-
1
(
L
-
K
)
.
For the cash-or-nothing option, we set u
N
x
n
+
1
=
C. This boundary condition was used in [26].


## 2.2.3. Neumann Boundary Condition

To specify values for the derivative of the solution at the boundary of the spatial domain, we have the following equation from (9):
(15)
∂
u
(
x
,
τ
)
∂
x
=
N
(
d
1
)
.
Therefore, for a sufficiently large L, we assume
(16)
∂
u
(
L
,
τ
)
∂
x
=
1
,
which we call by Neumann BC for European call option. Equation (16) can be discretized as (
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
)
/
h
=
1; that is, u
N
x
n
+
1
=
u
N
x
-
1
n
+
1
+
h. Therefore, we have
(17)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
2
β
N
x
-
2
γ
N
x
-
2
0
⋯
0
α
N
x
-
1
β
N
x
-
1
+
γ
N
x
-
1
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
2
n
+
1
u
N
x
-
1
n
+
1
=
b
1
b
2
⋮
b
N
x
-
2
b
N
x
-
1
-
h
γ
N
x
-
1
.
Similarly, applying BC in terms of the first derivative of x for the cash-or-nothing option, we can get
(18)
∂
u
(
x
,
τ
)
∂
x
=
C
e
-
(
r
τ
+
d
2
2
/
2
)
x
σ
2
π
τ
from (11). Therefore, we can assume
(19)
∂
u
(
L
,
τ
)
∂
x
=
0
for a sufficiently large L. This represents Neumann BC for cash-or-nothing option. We now replace (19) with the one-sided derivative (
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
)
/
h
=
0; that is, u
N
x
n
+
1
=
u
N
x
-
1
n
+
1.


## 2.2.4. Linear Boundary Condition

Linear boundary condition assumes that the second derivative of the option value with respect to the underlying asset pricex vanishes to zero for the large value of the asset price. To demonstrate this, we consider (9) at the right end of the domain. The second derivative value is given by
(20)
∂
2
u
(
x
,
τ
)
∂
x
2
=
e
-
d
1
2
/
2
σ
2
π
τ
x
.
If the asset price approaches to the large value L, we assume
(21)
∂
2
u
(
L
,
τ
)
∂
x
2
=
0
,
which we call by linear BC for European call option. Equation (21) can be discretized as (
u
N
x
-
1
n
+
1
-
2
u
N
x
n
+
1
+
u
N
x
+
1
n
+
1
)
/
h
2
=
0. By this relation, we obtain the boundary value as u
N
x
+
1
n
+
1
=
2
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1. By substituting this in (3), we get
(22)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
1
β
N
x
-
1
γ
N
x
-
1
0
⋯
0
α
N
x
-
γ
N
x
β
N
x
+
2
γ
N
x
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
1
n
+
1
u
N
x
n
+
1
=
b
1
b
2
⋮
b
N
x
-
1
b
N
x
.
Similarly, for the cash-or-nothing option, we have
(23)
∂
2
u
(
x
,
τ
)
∂
x
2
=
-
C
d
1
e
-
(
r
τ
+
d
2
2
/
2
)
σ
2
τ
2
π
x
2
from (11). Therefore, for a sufficiently large L, we have
(24)
∂
2
u
(
L
,
τ
)
∂
x
2
=
0
.
In this case, we can obtain the same BC as European call option.


## 2.2.5. PDE Boundary Condition

Next we consider the PDE BC [23, 27]. Using the BS equation itself, we derive the BC. We use one-sided discretization:
(25)
u
N
x
n
+
1
-
u
N
x
n
Δ
τ
=
σ
2
x
N
x
2
2
u
N
x
n
+
1
-
2
u
N
x
-
1
n
+
1
+
u
N
x
-
2
n
+
1
h
2
+
r
x
N
x
u
N
x
n
+
1
-
u
N
x
-
1
n
+
1
h
-
r
u
N
x
n
+
1
,
where the drift and volatility terms are discretized by using one-sided derivatives. The linear system of N
x equations can be written in the following matrix form:
(26)
β
1
γ
1
0
⋯
0
α
2
β
2
γ
2
⋯
0
⋮
⋱
⋱
⋱
⋮
0
⋯
α
N
x
-
1
β
N
x
-
1
γ
N
x
-
1
0
⋯
α
N
x
β
N
x
γ
N
x
u
1
n
+
1
u
2
n
+
1
⋮
u
N
x
-
1
n
+
1
u
N
x
n
+
1
=
b
1
b
2
⋮
b
N
x
-
1
b
N
x
,
where α
N
x
=
-
σ
2
x
N
x
2
/
2
h
2, β
N
x
=
σ
2
x
N
x
2
/
h
2
+
r
x
N
x
/
h, γ
N
x
=
1
/
Δ
τ
-
σ
2
x
N
x
2
/
2
h
2
-
r
x
N
x
/
h
+
r, and b
N
x
=
u
N
x
n
/
Δ
τ.


## 3. Numerical Experiments

To compare the above five BCs, we perform the numerical tests with European vanilla call option and cash-or-nothing option. Figures2(a) and 2(b) show European call option and cash-or-nothing option prices at τ
=
0 and τ
=
1 on Ω
=
[
0,300
], respectively. Here, σ
=
0.35 and r
=
0.05 are used. In the Appendix, we provide MATLAB codes for the closed-form solutions of these options.Valuesu
(
x
,
τ
) of (a) European call option and (b) cash-or-nothing option at τ
=
0 and τ
=
1, respectively.
(a)
(b)In the following sections, unless otherwise specified, we use strike priceK
=
100, cash C
=
1, the risk-free interest rate r
=
0.05, and volatility σ
=
0.35.


### 3.1. Convergence Test

First, we investigate the convergence of the numerical solutions with different BCs. To show this, we calculate the root mean square error with respect to spatial and temporal step sizes. Root mean square error (RMSE) is defined by(27)
RMSE
=
1
N
∑
x
j
∈
[
0.8
K
,
1.2
K
]
‍
u
j
n
-
v
j
n
2
,
where N is the number of grid points on [
0.8
K
,
1.2
K
]. Here, u is exact solution and v is numerical solution.Table1 represents the RMSE of European vanilla call option with five different BCs for varying h and Δ
τ at T
=
1. As shown in Table 1, we can observe that the RMSE converges to zero with decreasing space and time step sizes. From now on, we will use the Δ
τ
=
0.025 and h
=
0.5 for the following numerical tests.Table 1
RMSE for European call option with varyingh and Δ
τ with five different BCs.
h
Δ
τ
0.2
0.1
0.05
0.025
Dirichlet I
2.0
1.4525
e
-
3
1.4168
e
-
3
1.3986
e
-
3
1.3894
e
-
3
1.0
2.5706
e
-
4
2.5075
e
-
4
2.4754
e
-
4
2.4592
e
-
4
0.5
4.5474
e
-
5
4.4359
e
-
5
4.3790
e
-
5
4.3504
e
-
5
Dirichlet II
2.0
1.4560
e
-
3
1.4197
e
-
3
1.4010
e
-
3
1.3917
e
-
3
1.0
2.5768
e
-
4
2.5126
e
-
4
2.4797
e
-
4
2.4631
e
-
4
0.5
4.5583
e
-
5
4.4445
e
-
5
4.3867
e
-
5
4.3574
e
-
5
Neumann
2.0
1.4533
e
-
3
1.4173
e
-
3
1.3989
e
-
3
1.3896
e
-
3
1.0
2.5719
e
-
4
2.5085
e
-
4
2.4759
e
-
4
2.4595
e
-
4
0.5
4.5490
e
-
5
4.4375
e
-
5
4.3799
e
-
5
4.3509
e
-
5
Linear
2.0
1.4535
e
-
3
1.4174
e
-
3
1.3989
e
-
3
1.3896
e
-
3
1.0
2.5724
e
-
4
2.5086
e
-
4
2.4756
e
-
4
2.4595
e
-
4
0.5
4.5505
e
-
5
4.4377
e
-
5
4.3800
e
-
5
4.3509
e
-
5
PDE
2.0
1.4532
e
-
3
1.4173
e
-
3
1.3989
e
-
3
1.3896
e
-
3
1.0
2.5715
e
-
4
2.5084
e
-
4
2.4759
e
-
4
2.4595
e
-
4
0.5
4.5476
e
-
5
4.4373
e
-
5
4.3799
e
-
5
4.3509
e
-
5


### 3.2. Pointwise Error on Different Domain SizeL

To study the effects of the various BCs, we perform the numerical tests on different domain sizeL. For comparison, we compute the pointwise error as the absolute difference between the numerical solution u
num and the closed-form solution u
ex; that is, error
=
|
u
ex
-
u
num
|.First, we take the payoff of a European vanilla call option (7). In Figure 3, left and right columns represent the option values and the pointwise errors versus x, respectively. In Figure 3, we compare the numerical results with three different domain sizes L
=
150, 200, and 300. As domain size is large, there is no difference between the numerical results with five different boundary conditions and analytic solution at x
=
L. Also, pointwise error in [
0.8
K
,
1.2
K
] decreases when the domain size increases.Comparison of European call option values (left column) and pointwise error (right column) with five BCs on different domain: (a)L
=
150, (b) L
=
200, and (c) L
=
300, respectively.
(a)
L
=
150
(b)
L
=
200
(c)
L
=
300As second example, we test cash-or-nothing option. Similar to the previous case, we have same results as shown in Figure4. When L is large, the pointwise error decreases. However, when the domain size is not sufficiently large, we have gap between the analytic and numerical solutions.Comparison of cash-or-nothing option values (left column) and pointwise error (right column) with five BCs on different domain: (a)L
=
150, (b) L
=
200, and (c) L
=
300, respectively.
(a)
L
=
150
(b)
L
=
200
(c)
L
=
300Therefore, it is important to choose the domain sizeL. Also, to reduce the numerical error, we need to choose the proper boundary condition.


### 3.3. Root Mean Square Error

The following numerical tests illustrate RMSE of various BCs with differentr and σ. First we consider European call option with h
=
0.5, Δ
t
=
0.025, T
=
1, and L
=
300. Figures 5(a) and 5(b) represent the RMSE against the interest rate r and the volatility σ, respectively. Here, we use σ
=
0.35 and r
=
0.05 in Figures 5(a) and 5(b), respectively. As we can observe from Figure 5(a), the RMSE increases as r is large. Also, the RMSE versus σ has similar behavior in Figure 5(a). Through this test, we can see that Dirichlet I and Neumann BCs are more efficient than other BCs in case of r and σ.RMSE of various BCs versus (a)r and (b) σ with European call option.
(a)
(b)Now, we consider cash-or-nothing option withh
=
0.5, Δ
t
=
0.025, T
=
1, and L
=
300. Figure 6 shows the RMSE of various BCs versus r and σ. Here, we use σ
=
0.35 and r
=
0.05 in Figures 6(a) and 6(b), respectively. In this test, we can observe that Dirichlet II and Neumann BCs have lower RMSE than the other BCs for various r and σ.RMSE of various BCs versus (a)r and (b) σ with cash-or-nothing option.
(a)
(b)Figure7 illustrates RMSE of various BCs versus time T with European call option and cash-or-nothing option, respectively. As shown in Figure 7(a), we can see that the choice of BCs is important when τ is large. In addition, we can confirm that Neumann BC is efficient to reduce numerical error. Similarly, we observe that Neumann BC is efficient when L is small or T is large in Figure 7(b).RMSE of various BCs versusT with (a) European call option and (b) cash-or-nothing option.
(a)
(b)


## 3.1. Convergence Test

First, we investigate the convergence of the numerical solutions with different BCs. To show this, we calculate the root mean square error with respect to spatial and temporal step sizes. Root mean square error (RMSE) is defined by(27)
RMSE
=
1
N
∑
x
j
∈
[
0.8
K
,
1.2
K
]
‍
u
j
n
-
v
j
n
2
,
where N is the number of grid points on [
0.8
K
,
1.2
K
]. Here, u is exact solution and v is numerical solution.Table1 represents the RMSE of European vanilla call option with five different BCs for varying h and Δ
τ at T
=
1. As shown in Table 1, we can observe that the RMSE converges to zero with decreasing space and time step sizes. From now on, we will use the Δ
τ
=
0.025 and h
=
0.5 for the following numerical tests.Table 1
RMSE for European call option with varyingh and Δ
τ with five different BCs.
h
Δ
τ
0.2
0.1
0.05
0.025
Dirichlet I
2.0
1.4525
e
-
3
1.4168
e
-
3
1.3986
e
-
3
1.3894
e
-
3
1.0
2.5706
e
-
4
2.5075
e
-
4
2.4754
e
-
4
2.4592
e
-
4
0.5
4.5474
e
-
5
4.4359
e
-
5
4.3790
e
-
5
4.3504
e
-
5
Dirichlet II
2.0
1.4560
e
-
3
1.4197
e
-
3
1.4010
e
-
3
1.3917
e
-
3
1.0
2.5768
e
-
4
2.5126
e
-
4
2.4797
e
-
4
2.4631
e
-
4
0.5
4.5583
e
-
5
4.4445
e
-
5
4.3867
e
-
5
4.3574
e
-
5
Neumann
2.0
1.4533
e
-
3
1.4173
e
-
3
1.3989
e
-
3
1.3896
e
-
3
1.0
2.5719
e
-
4
2.5085
e
-
4
2.4759
e
-
4
2.4595
e
-
4
0.5
4.5490
e
-
5
4.4375
e
-
5
4.3799
e
-
5
4.3509
e
-
5
Linear
2.0
1.4535
e
-
3
1.4174
e
-
3
1.3989
e
-
3
1.3896
e
-
3
1.0
2.5724
e
-
4
2.5086
e
-
4
2.4756
e
-
4
2.4595
e
-
4
0.5
4.5505
e
-
5
4.4377
e
-
5
4.3800
e
-
5
4.3509
e
-
5
PDE
2.0
1.4532
e
-
3
1.4173
e
-
3
1.3989
e
-
3
1.3896
e
-
3
1.0
2.5715
e
-
4
2.5084
e
-
4
2.4759
e
-
4
2.4595
e
-
4
0.5
4.5476
e
-
5
4.4373
e
-
5
4.3799
e
-
5
4.3509
e
-
5


## 3.2. Pointwise Error on Different Domain SizeL

To study the effects of the various BCs, we perform the numerical tests on different domain sizeL. For comparison, we compute the pointwise error as the absolute difference between the numerical solution u
num and the closed-form solution u
ex; that is, error
=
|
u
ex
-
u
num
|.First, we take the payoff of a European vanilla call option (7). In Figure 3, left and right columns represent the option values and the pointwise errors versus x, respectively. In Figure 3, we compare the numerical results with three different domain sizes L
=
150, 200, and 300. As domain size is large, there is no difference between the numerical results with five different boundary conditions and analytic solution at x
=
L. Also, pointwise error in [
0.8
K
,
1.2
K
] decreases when the domain size increases.Comparison of European call option values (left column) and pointwise error (right column) with five BCs on different domain: (a)L
=
150, (b) L
=
200, and (c) L
=
300, respectively.
(a)
L
=
150
(b)
L
=
200
(c)
L
=
300As second example, we test cash-or-nothing option. Similar to the previous case, we have same results as shown in Figure4. When L is large, the pointwise error decreases. However, when the domain size is not sufficiently large, we have gap between the analytic and numerical solutions.Comparison of cash-or-nothing option values (left column) and pointwise error (right column) with five BCs on different domain: (a)L
=
150, (b) L
=
200, and (c) L
=
300, respectively.
(a)
L
=
150
(b)
L
=
200
(c)
L
=
300Therefore, it is important to choose the domain sizeL. Also, to reduce the numerical error, we need to choose the proper boundary condition.


## 3.3. Root Mean Square Error

The following numerical tests illustrate RMSE of various BCs with differentr and σ. First we consider European call option with h
=
0.5, Δ
t
=
0.025, T
=
1, and L
=
300. Figures 5(a) and 5(b) represent the RMSE against the interest rate r and the volatility σ, respectively. Here, we use σ
=
0.35 and r
=
0.05 in Figures 5(a) and 5(b), respectively. As we can observe from Figure 5(a), the RMSE increases as r is large. Also, the RMSE versus σ has similar behavior in Figure 5(a). Through this test, we can see that Dirichlet I and Neumann BCs are more efficient than other BCs in case of r and σ.RMSE of various BCs versus (a)r and (b) σ with European call option.
(a)
(b)Now, we consider cash-or-nothing option withh
=
0.5, Δ
t
=
0.025, T
=
1, and L
=
300. Figure 6 shows the RMSE of various BCs versus r and σ. Here, we use σ
=
0.35 and r
=
0.05 in Figures 6(a) and 6(b), respectively. In this test, we can observe that Dirichlet II and Neumann BCs have lower RMSE than the other BCs for various r and σ.RMSE of various BCs versus (a)r and (b) σ with cash-or-nothing option.
(a)
(b)Figure7 illustrates RMSE of various BCs versus time T with European call option and cash-or-nothing option, respectively. As shown in Figure 7(a), we can see that the choice of BCs is important when τ is large. In addition, we can confirm that Neumann BC is efficient to reduce numerical error. Similarly, we observe that Neumann BC is efficient when L is small or T is large in Figure 7(b).RMSE of various BCs versusT with (a) European call option and (b) cash-or-nothing option.
(a)
(b)


## 4. Conclusion

In this paper, we reviewed and studied the performance of the five different boundary conditions such as Dirichlet, Neumann, linear, and partial differential equation boundary conditions for the numerical solutions of the BS partial differential equation. We used a finite difference method to numerically solve the BS equation with the five different BCs. To show the efficiency of the given boundary condition, several numerical examples such as a convergence test, domain size effect, and parameter effect are presented. Numerical results suggested that the linear boundary condition is accurate and efficient among the other boundary conditions. As a future research, we will investigate the BCs on the multidimensional BS equations.



---

*Source: 359028-2015-03-02.xml*