A New Approach for the Black–Scholes Model with Linear and Nonlinear Volatilities

Gulen, Seda; Popescu, Catalin; Sari, Murat

doi:10.3390/math7080760

Open AccessArticle

A New Approach for the Black–Scholes Model with Linear and Nonlinear Volatilities

by

Seda Gulen

¹,

Catalin Popescu

²

and

Murat Sari

^3,*

¹

Department of Mathematics, Tekirdag Namik Kemal University, Degirmenalti, Tekirdag 59030, Turkey

²

Department of Business Administration, Petroleum-Gas University, Blvd. Bucuresti, no.39, 100680 Ploiesti, Romania

³

Department of Mathematics, Yildiz Technical University, Istanbul 34220, Turkey

^*

Author to whom correspondence should be addressed.

Mathematics 2019, 7(8), 760; https://doi.org/10.3390/math7080760

Submission received: 24 July 2019 / Revised: 14 August 2019 / Accepted: 16 August 2019 / Published: 19 August 2019

(This article belongs to the Special Issue Financial Mathematics)

Download

Browse Figures

Versions Notes

Abstract

:

Since financial engineering problems are of great importance in the academic community, effective methods are still needed to analyze these models. Therefore, this article focuses mainly on capturing the discrete behavior of linear and nonlinear Black–Scholes European option pricing models. To achieve this, this article presents a combined method; a sixth order finite difference (FD6) scheme in space and a third–order strong stability preserving Runge–Kutta (SSPRK3) over time. The computed results are compared with available literature and the exact solution. The computed results revealed that the current method seems to be quite strong both quantitatively and qualitatively with minimal computational effort. Therefore, this method appears to be a very reliable alternative and flexible to implement in solving the problem while preserving the physical properties of such realistic processes.

Keywords:

Black–Scholes equation; option pricing modelling; European option; volatility; high-order finite difference

1. Introduction

In the last few decades, the corporations looked for important tools in terms of financial securities. As part of the financial securities, options are mainly used to assure assets in order to cover the risks generated in the stock prices changes [1]. To properly understand what alternatives we have for these options, Lesmana and Wang [2] stated that there are two main types of options: European options can only be exercised on a given date (expiry date) and American options can be exercised on or before the expiry date.

An economic issue, which is important for both the traders and the investors, is to have a proper method for price options to determine the appropriate price or theoretical value for a call or a put option. During the time, price options were not so commonly explained by the traders as a sufficient financial instrument and valuation of an option has always been a difficult and challenging task. According to Rad et al. [3], options are priced using mathematical models that are often challenging to solve. Of course, in the last 50 years, there have been some attempts to introduce different methods for price options. For instance, a major step was taken by Black and Scholes [4,5,6] in 1973 to propose a mathematical model for calculating a fair value for an option. According to them [4], the Black–Scholes model is a pure log-normal diffusion model and it leads to a parabolic partial differential equation under Ito’s calculus.

Another important contribution was made by Merton [7], who extended the model equation proposed by Black and Scholes. All three of them demonstrated that their formulae leading to partial differential equations could help to determine a fair value for a call or a put option. One of the main obstacles identified by Lesmana and Wang [2] is that the Black–Scholes option pricing methodology is no longer valid in the presence of transaction costs on trading in the risk-free security or stock. To overcome this drawback, different models were proposed by Leland [8], Boyle and Vorst [9], Kusuoka [10] or Barles and Soner [11]. In addition, Cox and Ross [12] published an article aimed at finding ways to overcome difficulties with payouts and bankruptcies. Another important contribution is mainly related to the Crank–Nicolson scheme used in numerical pricing within the Black–Scholes framework [13]. A good example of this scheme is the work of McCartin and Labadie [14] in the case of pricing vanilla options. In addition, Ankudinova and Ehrhardt [15] used a Crank–Nicolson method and devised a numerical scheme for the linearized Black–Scholes equation using the frozen values of nonlinear volatility. From the point of view of nonlinear volatility, an approach based on the method of lines and the backward Euler scheme was proposed by Company et al. [16]. Other than these, in recent years, Jeong et al. [17], Koleva et al. [18], Manisha [19] and Sari and Bălăcescu [20], Mashayekhi and Hugger [21] and Hout and Valkov [22] used various schemes based on finite differences to analyze the Black–Scholes model. Jeong et al. [17] proposed a finite difference scheme not using the far-field boundary condition. Koleva et al. [18] construct a fourth order compact finite difference scheme in space for European options modeling markets to liquidity shocks. A high-order difference scheme was considered by Manisha [19] for the Black–Scholes equation, which governs four option styles of European-type variable parameters. To analyze the option pricing model, Sari and Balacescu [20] proposed a fourth-order difference and MacCormack schemes. One of the most comprehensive models, the Barles and Soner [11] model was solved with various difference schemes to compare their precision and order of convergence [21]. Hout and Valkov [22] solved the European two-asset options by a finite difference based numerical method considering non-uniform grids. In addition, Hendrick et al. [23] combined a high-order finite difference scheme with the Alternating Direction Implicit scheme for the parabolic partial differential equation representing the European basket option in a sparse grid setting.

Another work by Markolefas [24] demonstrated the ability of a finite element method to formulate discrete approximation models resulting from the Black–Scholes pricing model. Specifically, his/her work demonstrates a version of Galerkin finite element methods to capture accurate and efficient behaviour: the parabolic partial differential equations, with a complex initial, boundary and/or internal conditions, resulting from various option pricing theories. However, Lin and Chen [25] presented a European option pricing model by applying the Model-Order-Reduction method, which is able to reduce the order of the original finite difference method systems.

Likewise, in the case of a nonlinear parabolic partial differential equation that governs the European option pricing in transaction costs, it is worth pointing out the contribution to Lesmana and Wang [2] who proposed an upwind difference scheme for spatial discretization. In addition, there are more than 30 years since the two experts of banking and financial markets developed the first commercially available pricing formula for options linked to the average price of crude oil; they called this option as Asian option [26]. The Asian option, also known as average options, is an option where the return on return is dependent on the arithmetic average of the price of the underlying asset over its lifetime [27,28,29,30]. Kumar et al. [27] described the implementation of the local radial basis function based on a grid-free method for the numerical solution of the Asian option. In addition, Rad et al. [28] described the valuation of the Asian option with a radial basis function approach based on the finite difference method.

In another representative work created by Sin [31], he/she used some numerical methods to compute barrier options considering both Black–Scholes and Heston models. In addition, he/she computed barrier options using the standard Monte Carlo method. The method was improved by Moon [32] in dealing with the problem. The key idea behind his/her approach, the improved one, is to use the exit probability and uniformly distributed random numbers to predict the first instance of hitting the barriers. In addition, Zeng et al. [33] focused on a comparative study of the Monte Carlo simulation and finite difference method for the European option pricing under a regime-switching framework.

In spite of the great efforts made to deal with the problem in the literature, the researchers chose to use relatively low order methods [34,35,36,37,38,39] for the problem.

In recent years, the option market has become increasingly important. This increases the need for higher computational capacity as the complexity of financial models increases in an operating area where real-time market calculations are required. The well-known Black–Scholes model provides a closed form for the values of certain options such as European style. However, for the solutions of more complex models, either there is no closed form solution or, in the absence of tractable closed form, computationally intensive approximate techniques are required. In this case, the finite difference based methods are commonly employed in price options. Therefore, in this article, a sixth order finite difference scheme (FD6) in space and a third-order strong stability preserving Runge–Kutta (SSP-RK3) in time have been combined to obtain effective numerical solutions of the European put option problem that has an exact closed-form solution. This choice provides a direct and accurate estimation of approximation error. High-order approximations reduce the computational effort to achieve the required accuracy. Therefore, they are highly popular in solving real world problems and for all financial markets increased calculation speed is a strategic advantage. Moreover, the choice of time integration is important in many respects such as accuracy, computational effort and stability. The SSP-RK3 scheme, which is a class of high order SSP-time integration techniques, was developed by Gottlieb et al. [40] to investigate hyperbolic conservation laws with a spatial discretization guaranteeing the stability properties expected of the forward Euler method [40]. Therefore, it was considered that the combination of the FD6 and SSP-RK3 could solve the corresponding model with high accuracy. In addition, this method is quite straightforward to write codes in any programming language. Our computations show that the results produced by the current method approximates the exact solution and available solution in the literature [41] very well. To the best knowledge of the authors, this method has not been implemented for the problem represented by the European put option pricing model.

The rest of the paper is organized as follows: in Section 2, the mathematical formulation of the European put option problem is introduced. In Section 3, numerical methods for space and time discretizations are summarized. The results and illustrations are presented in Section 4. Final remarks are reported in Section 5.

2. Mathematical Formulation

The Black–Scholes model is of great importance in option pricing theory. The model was developed by Black and Scholes [4] and previously by Merton [7]. The value of an option

V (S, t)

is a function of stock price S and time t can be determined by the following Black–Scholes partial differential equation,

\frac{\partial V}{\partial t} = L V, t \in (0, T],

(1)

where

L = - \frac{1}{2} σ^{2} S^{2} \frac{\partial^{2}}{\partial S^{2}} - r S \frac{\partial}{\partial S} + r, S \in (0, + \infty),

and

r \geq 0

and

σ

are positive constant parameters, referred to risk-free interest rate and volatility that are both known functions over the life of the option. Here, the stock price S is modeled by a Geometric Brownian motion, i.e., S satisfies the following stochastic differential equation

d S = μ S d t + σ S d W,

(2)

where

μ > 0

is the drift rate of the stock and W is a standard Brownian motion. The Black–Scholes equation is given by Equation (1), which assumes that the market is frictionless, is a trading enviroment and there are no transaction costs and restraints on transactions. In addition, Equation (1) does not include arbitrage opportunities.

The Black–Scholes Equation (1) is very effective for modeling pricing options in a complete market without transaction costs. However, in the presence of transaction cost on trading in the riskless security or stock, the Black–Scholes option pricing methodology is no longer valid, since perfect hedging is impossible [2]. Due to transaction cost where the volatility

σ

can depend on the time t, the stock price S or the derivatives of the option price V itself, the Black–Scholes model (1) is transformed to the following nonlinear equation:

\frac{\partial V}{\partial t} = - \frac{1}{2} \tilde{σ} (t, S, V_{s}, V_{S S}) S^{2} \frac{\partial^{2} V}{\partial S^{2}} - r S \frac{\partial V}{\partial S} + r V, S \in (0, + \infty), t \in (0, T],

(3)

where

\tilde{σ}

is the modified volatility as a function of

t, S, V_{S}, V_{S S}

. The problem represented by the model equation is very complex and depends on too many parameters, which makes the model more realistic, and, in this case, careless use of those parameters can lead to unrealistic results. Among these parameters, volatility is one of the most considerable factors in the pricing model. Volatility must be considered carefully, as the option value depends on the future stock price. The option pricing model with transaction cost has been studied in the literature by some researchers, such as Boyle and Vorst [10] proposed an option price with the volatility of the form

\tilde{σ} = σ_{0} {(1 + c A)}^{1 / 2}, A = \frac{μ}{σ_{0} \sqrt{Δ T}}, c = 1,

(4)

where

μ

is the proportional transaction cost,

Δ t

is the transaction period and

σ_{0}

is the original volatility constant.

Leland [8] obtained that the option price is the solution to Equation (2) with the volatility

{\tilde{σ}}^{2} = σ_{0}^{2} (1 + L e \times s i g n (V_{S S})),

(5)

where

L e

is the Leland number given by

L e = \sqrt{\frac{2}{π}} (\frac{κ}{σ_{0} \sqrt{δ t}}),

where

δ t

and

κ

represent the transaction frequency and transaction cost measure, respectively.

A more complex model was proposed by Barles and Soner [11]. In their model, the nonlinear volatility is given by

{\tilde{σ}}^{2} = σ_{0}^{2} (1 + Ψ [e^{r (T - t)}h a^{2}h S^{2}h V_{S S}]),

(6)

where

a = κ \sqrt{γ N}

with

κ

being the transaction cost parameter,

γ

being a risk aversion factor and N being the number of options to be sold. The function

Ψ

is the solution to the following nonlinear initial value problem:

Ψ^{'} (z) = \frac{Ψ (z) + 1}{2 \sqrt{z Ψ (z)} - z} for z \neq 0 and Ψ (0) = 0 .

(7)

Equation (7) implies that

lim_{x \to \infty} \frac{Ψ (z)}{z} = 1 and lim_{x \to - \infty} Ψ (z) = - 1 .

(8)

This property accepts the treatment of the function

Ψ (.)

as the identity for large arguments and therefore the volatility becomes

{\tilde{σ}}^{2} = σ^{2} (1 + e^{r (T - t)} a^{2} S^{2} V_{S S}) .

(9)

For solving Equations (1) and (3) uniquely, one final condition and two boundary conditions are required. These conditions depend on both types of options (put or call) and whether or not the stock pays cash dividends. The European put option, which gives the holder the right to sell for the strike price K at maturity time T, is studied in this work, and the boundary and final conditions are as follows:

V (S, T) = m a x {K - S, 0}, 0 \leq S < \infty,

(10)

V (S, t) \sim 0, S \to \infty,

(11)

V (0, t) = K e^{r (T - t)}, 0 \leq t \leq T .

(12)

The existence and uniqueness of the solution of the problem were shown in the context of stochastic optimal theory in the study of Barles and Soner [11].

The exact solution of the linear European put option problem is determined by

V (S, t) = K e^{- r (T - t)} N (- γ_{2}) - S e^{- δ (T - t)} N (- γ_{1}),

(13)

with the parameters

γ_{1} = \frac{l n S - l n K + (r - δ + \frac{1}{2} σ^{2}) (T - t)}{σ \sqrt{T - t}},

γ_{2} = γ_{1} - σ \sqrt{T - t},

N (y) = \frac{1}{\sqrt{2 π}} \int_{- \infty}^{y} e^{- \frac{1}{2} x^{2}} d x,

where N is the standard normal cumulative probability distribution function and

δ

is continuous dividend yield [42].

3. The Solution Method

This section is dedicated to a numerical solution of the linear and nonlinear European put option problems. For the sake of brevity, the nonlinear European put option model is discussed with the nonlinear volatility (9) proposed by Barles and Soner [11]. In spatial discretization, a sixth-order finite difference method (FD6) is applied and a third-order strong stability preserving a Runge–Kutta (SSP-RK3) method is considered in temporal discretization. The FD6 scheme based on high order differences achieves the required accuracy. In addition, the SSP-RK3 scheme, which is a class of the SSP-time integration techniques with the spatial discretization, guarantees the stability [40]. Hence, it is expected that the proposed method will solve the problem effectively.

3.1. Spatial Discretization

Spatial derivatives are computed by the FD6 scheme based on the Taylor series expansion [43]. In order to use finite difference approximation, it is started by defining a uniform grid consisting of N points satisfying

0 = x_{0} < x_{1} < \dots < x_{N} = x_{m a x}

. The step size

h = x_{i + 1} - x_{i}

,

i = 0, 1, \dots, N - 1

is equal to each other at any point i. The first derivative

V_{i}^{^{'}}

at point i can be approximated by

V_{i}^{^{'}} = \frac{1}{h} \sum_{j = - L}^{R} a_{j + L} V_{i + j}, 1 \leq i \leq N,

(14)

where R and L represent the number of points on the right-hand side and left-hand side for taken nodes, respectively. Hence, Equation (14) involves

(R + L + 1)

constants,

a_{0}, a_{1}, \dots, a_{R + L}

. R is equal to L for the considered nodes at internal points, but this is not the case for the boundary points. The coefficients

a_{j}

are calculated by using the Taylor expansion near point i. The FD6 scheme uses seven points at interior and boundary points as follows:

V_{i}^{'} (t) = \frac{- 147 V_{i} + 360 V_{i + 1} - 450 V_{i + 2} + 400 V_{i + 3} - 225 V_{i + 4} + 72 V_{i + 5} - 10 V_{i + 6}}{60 h}, i = 1,

V_{i}^{'} (t) = \frac{- 10 V_{i - 1} - 77 V_{i} + 150 V_{i + 1} - 100 V_{i + 2} + 50 V_{i + 3} - 15 V_{i + 4} + 2 V_{i + 5}}{60 h}, i = 2,

V_{i}^{'} (t) = \frac{2 V_{i - 2} - 24 V_{i - 1} - 35 V_{i} + 80 V_{i + 1} - 30 V_{i + 2} + 8 V_{i + 3} - V_{i + 4}}{60 h}, i = 3,

V_{i}^{'} (t) = \frac{- V_{i - 3} + 9 V_{i - 2} - 45 V_{i - 1} + 45 V_{i + 1} - 9 V_{i + 2} + V_{i + 3}}{60 h}, 4 \leq i \leq N - 3,

(15)

V_{i}^{'} (t) = \frac{- 2 V_{i + 2} + 24 V_{i + 1} + 35 V_{i} - 80 V_{i - 1} + 30 V_{i - 2} - 8 V_{i - 3} + V_{i - 4}}{60 h}, i = N - 2,

V_{i}^{'} (t) = \frac{10 V_{i + 1} + 77 V_{i} - 150 V_{i - 1} + 100 V_{i - 2} - 50 V_{i - 3} + 15 V_{i - 4} - 2 V_{i - 5}}{60 h}, i = N - 1,

V_{i}^{'} (t) = \frac{147 V_{i} - 360 V_{i - 1} + 450 V_{i - 2} - 400 V_{i - 3} + 225 V_{i - 4} - 72 V_{i - 5} + 10 V_{i - 6}}{60 h}, i = N .

The sixth order scheme (15) can be written in a more compact form as follows:

V^{'} = \frac{1}{h} A V .

(16)

The second order derivative terms are obtained by applying the first operator twice:

V^{″} = \frac{1}{h} A V^{'},

(17)

where

V = {(V_{1}, V_{2}, . ., V_{N})}^{T}

and

A = [\begin{array}{l} - \frac{147}{60} & \frac{360}{60} & - \frac{450}{60} & \frac{400}{60} & - \frac{225}{60} & \frac{72}{60} & - \frac{10}{60} \\ \frac{10}{60} & - \frac{77}{60} & \frac{150}{60} & - \frac{100}{60} & \frac{50}{60} & - \frac{15}{60} & \frac{2}{60} \\ \frac{2}{60} & - \frac{24}{60} & - \frac{35}{60} & - \frac{80}{60} & - \frac{30}{60} & \frac{8}{60} & - \frac{8}{60} \\ - \frac{1}{60} & \frac{9}{60} & - \frac{45}{60} & 0 & \frac{45}{60} & - \frac{9}{60} & \frac{1}{60} \\ - \frac{1}{60} & \frac{9}{60} & - \frac{45}{60} & 0 & \frac{45}{60} & - \frac{9}{60} & \frac{1}{60} \\ ⋱ & ⋱ & ⋱ & ⋱ & ⋱ & ⋱ \\ - \frac{1}{60} & \frac{9}{60} & - \frac{45}{60} & 0 & \frac{45}{60} & - \frac{9}{60} & \frac{1}{60} \\ \frac{1}{60} & - \frac{8}{60} & \frac{30}{60} & - \frac{80}{60} & \frac{35}{60} & \frac{24}{60} & \frac{- 2}{60} \\ - \frac{2}{60} & \frac{15}{60} & - \frac{50}{60} & \frac{100}{60} & - \frac{150}{60} & \frac{77}{60} & \frac{10}{60} \\ \frac{10}{60} & - \frac{72}{60} & \frac{225}{60} & - \frac{400}{60} & \frac{450}{60} & - \frac{360}{60} & \frac{147}{60} \end{array}] .

3.2. Temporal Discretization

The SSP-RK3 scheme for the discretization of Equation (1) is presented in time. A class of high-order SSP time discretization techniques was introduced by Gottlieb et al. [40] for solving hyperbolic conservation laws with stable spatial discretizations. The SSP methods guarantee the stability properties expected of the forward Euler method [40].

The computational domain for time consists of M points satisfying

0 = t^{0} < t^{1} < \dots < t^{M} = T

. The uniform time step

d t = t^{n + 1} - t^{n}, n = 0, 1, 2, \dots, M - 1

is equidistant at any point n. After applying the FD6 method, Equation (1) can be reduced into a set of ordinary differential equations in time. Then, the equation system can be expressed by

\frac{d V_{i}}{d t} = L V_{i} .

(18)

The SSP-RK3 scheme integrates the semi-discrete Equation (18) from time

t_{0}

to

t_{0} + d t

through the operations

V_{i}^{(1)} = V_{i}^{n} + d t L V_{i}^{n},

V_{i}^{(2)} = \frac{3}{4} V_{i}^{n} + \frac{1}{4} V_{i}^{(1)} + \frac{1}{4} d t L V^{(1)},

(19)

V_{i}^{n + 1} = \frac{1}{3} V_{i}^{n} + \frac{2}{3} V^{(2)} + \frac{2}{3} d t L V^{(2)},

where L is the discretization form of

L

.

For the SSP-RK3, the total variation (TV) of the numerical solution

T V (V) = \sum_{i} ∣ V_{i + 1} - V_{i} ∣

(20)

does not increase in time [40];

T V (V^{k + 1}) \leq T V (V^{k}) .

(21)

4. Numerical Results and Analysis

To show efficiency and accuracy of the proposed methods, numerical experiments on the European put option model have been performed. For the computations through the current schemes, computer codes have been produced in MATLAB^® 2019. The validity of the scheme has been verified through the produced results.

Semi-infinite domain

[0, + \infty)

in Equation (1) is replaced by the bounded interval

[0, S_{m a x})

where

S_{m a x}

is an artifical limit will be chosen larger than three to four times the exercise price. In order to use the difference approximation, it is started by defining a uniform grid consisting of N points satisfying

0 = S_{0} < S_{1} \dots < S_{N} = S_{m a x}

. The step size

d S = S_{i + 1} - S_{i}

,

i = 0, 1,..., N - 1

is equal to each other at any point i. Similarly, the interval

(0, T)

is divided into M subintervals with mesh nodes satisfying

0 = t^{0} < t^{1} <...< t^{M} = T

. The step size

d t = t^{n + 1} - t^{n}

,

n = 0, 1, 2,..., M

is equal to each other at any point n. In the nonlinear model, the square of volatility approximation

{\tilde{σ}}^{2}

given by Equation (9) has been discretized at the node

S_{i}

and time t, being

{\tilde{σ}}_{i}^{2} = σ^{2} (1 + e^{r (T - t)} a^{2} S_{i}^{2} Δ_{i} (V (t))),

(22)

where

Δ_{i} (V (t)) = \frac{V (S_{i + 1}, t) - 2 V (S_{i}, t) + V (S_{i - 1}, t)}{d S^{2}}, i = 0, 1 and i = N, N + 1

(23)

and

Δ_{i} (V (t)) = \frac{1}{180 d S^{2}} (2 V (S_{i - 3}, t) - 27 V (S_{i - 2}, t) + 270 V (S_{i - 1}, t) - 490 V (S_{i}, t) + 270 V (S_{i + 1}, t)

- 27 V (S_{i + 2}, t) + 2 V (S_{i + 3}, t)), i = 3, N - 1

(24)

and has been used in the discretized form of Equation (3).

The parameters

T = 0.25

,

K = 10.0

,

r = 0.1

,

σ = 0.4

,

M = 2000

,

N = 200

[41] and

a = 0.02

are used in all calculations. For the purpose of being realistic, the results obtained by the proposed method are compared with the literature. Complience with the literature and with a realistic problem was taken into consideration in the selection of parameters. In terms of the realistically considered parameters, the computed results are believed to be helpful in making a decision for a viable trading strategy.

The numerical results and exact solution for the linear and numerical results for the nonlinear European put option model are presented in Table 1. The computed results revealed that the current method approximates the exact solution very well and is applicable, effective and easy to use. The qualitative behaviour of the option value is plotted in Figure 1, Figure 2 and Figure 3. Furthermore, the values of the option for the three different values of a are represented in Figure 4. From the figure, it can be seen that qualitative behavior is in agreement with the literature [2].

After applying the methods to the linear Black–Scholes equation and initial and boundary conditions, the resulting system has

M^{N}

unknowns. Since the present methods are explicit, the stability condition is given by

d t \leq \frac{{(d S)}^{2}}{2 a},

(25)

where

a = \frac{1}{2} σ^{2} S^{2}

is the coefficient in front of the second derivative term in the Black–Scholes equation. In the current computation, this condition is provided at each time and spatial step. In the nonlinear model, the nonlinear coefficient

\tilde{σ_{i}}

is frozen at each discretized time level, therefore, the condition (25) is provided.

To show accuracy of the method for the linear model, convergence rates are calculated by

R a t e (∣ ∣ . ∣ ∣_{\infty}) = \frac{|| V_{d s}^{d t} - V_{e x a c t} (S_{i}, t^{j}) {||}_{\infty}}{|| V_{d s / 2}^{d t / 2} - V_{e x a c t} (S_{i}, t^{j}) {||}_{\infty}},

(26)

R a t e (|| . {||}_{2}) = \frac{|| V_{d s}^{d t} - V_{e x a c t} (S_{i}, t^{j}) {||}_{2}}{|| V_{d s / 2}^{d t / 2} - V_{e x a c t} (S_{i}, t^{j}) {||}_{2}},

(27)

where

V_{d s}^{d t}

represents the solution with the spatial mesh size

d s

and the time mesh size

d t

,

V_{e x a c t}

represents the exact solution for the linear model whilst

∣ ∣.∣ ∣_{\infty}

and

∣ ∣.∣ ∣_{2}

are maximum norm and

L^{2}

-norm are given, respectively, by

∣ ∣ V_{d s}^{d t} - V_{e x a c t} ∣ ∣_{\infty} : = max_{1 ⩽ j ⩽ M; 1 ⩽ i ⩽ N} ∣ V_{i}^{j} - V_{e x a c t} (S_{i}, t^{j}) ∣,

(28)

∣ ∣ V_{d s}^{d t} - V_{e x a c t} ∣ ∣_{\infty} : = {(\sum_{1 \leq j \leq M} \sum_{1 \leq i \leq N} ∣ V_{i}^{j} - V_{e x a c t} (S_{i}, t^{j}) ∣^{2} d s d t)}^{1 / 2} .

(29)

To determine these rates, a sequence of meshes generated by halving the mesh sizes of the previous ones by starting from a given coarse mesh is accepted. As seen in Table 2, the current method is seen to be effective to enhance the accuracy of the numerical solution. When M and N are halved, the error reduced by a factor that is equal to ten. The order of convergence of the combination of FD6 and SSPRK3 is about

6.6

in

∣ ∣.∣ ∣_{\infty}

and

4.4

in

∣ ∣.∣ ∣_{2}

.

For the nonlinear model, since the exact solution is unknown, taking fixed value of

d t

, the difference

V (S, t, d s = k) - V (S, t, d s = k / 2)

is plotted in Figure 5. In addition, Table 3 presents the errors of the accuracy of the linear and nonlinear models that are calculated by the

ℓ^{2}

-error given by

E r r o r (t = 0) = {(\frac{1}{d s} \sum_{i = 1}^{N + 1} ∣ v (x_{i}, 0) - V_{i}^{0} ∣^{2})}^{1 / 2} .

(30)

Reference solution is denoted by

v (x_{i}, 0),

which has been obtained by the proposed method on a fine grid

d t

= 6.25 × 10

^{- 5}

and

d s = 0.1

. The numerical solution

V_{i}^{0}

has been calculated with the previously indicated parameters.

To indicate the influence of transaction cost modeled by volatility (9), the difference

V_{n o n l i n e a r} (S, t) - V_{l i n e a r} (S, t)

between the price of the European put option with transaction costs and the price of the European put without transaction costs is plotted in Figure 6. The numerical results present an economically significant price deviation between the linear model and the nonlinear model. It is seen that the difference is not symmetric but decreases closer to the expiry date and is maximal close to the stock price

S = 10

, where the nonlinear price is higher than the linear price.

In Table 4, the required times are compared when the proposed method is applied for the linear and nonlinear option problems. For both situations, the required CPU times increase with respect to the number of spatial and temporal discretizations. Although part of the computational time pertains to the calculation of the nonlinear volatility, the CPU times of the nonlinear model for every N and M discretization point are close to the linear model. The present work provides effective estimates in terms of easy programming and relatively low cost.

5. Conclusions and Recommendations

The complexity and stochastic properties of option pricing problems make it difficult to determine the value of the option. Therefore, accurate approximation methods are needed to understand the behavior of these problems and those are of great importance for scientific developments in financial markets. In this paper, a combination of a sixth-order finite difference scheme in space and a third-order strong stability preserving Runge–Kutta in time has been implemented to obtain effective numerical solutions of the linear and nonlinear European put option models represented by the Black–Scholes equation. The convergence of the solution has been measured by some error norms and it is confirmed that the present method is asymptotically convergent. In addition, the produced results are in good agreement with the literature and the exact solution. Therefore, the proposed method provides a better perspective to describe behavior of the option pricing model represented by the Black–Scholes equation and can be preferred due to reliability and accuracy with minimal computational effort. This study is believed to help researchers who want to uncover challenging financial and stochastic behaviours in modeling. In terms of the realistically considered parameters, it is believed that the computed results are helpful in making decisions for viable trading strategies. The valuation of the European option pricing problem by the proposed method can be compared with realistic option values.

Author Contributions

All authors contributed equally to this paper and they also read and finalized the manuscript.

Funding

This research received no external funding.

Conflicts of Interest

The authors declare no conflicts of interest.

References

Zhao, J.; Davidson, M.; Corless, R.M. Compact finite difference method for American option pricing. J. Comput. Appl. Math. 2007, 206, 306–321. [Google Scholar] [CrossRef] [Green Version]
Lesmana, D.C.; Wang, S. An upwind finite difference method for a nonlinear Black–Scholes equation governing European option valuation under transaction cost. Appl. Math. Comput. 2013, 219, 8811–8828. [Google Scholar] [CrossRef]
Rad, J.A.; Parand, K.; Abbasbandy, S. Local week form meshles techniques based on the radial point interpolation (RPI) method and local boundary integral equation (LBIE) method to evaluate European and American options. Commun. Nonlinear Sci. Numer. Simul. 2015, 22, 1178–1200. [Google Scholar] [CrossRef]
Black, F.; Scholes, M. The pricing of options and other corporate liabilities. J. Polit. Econ. 1973, 81, 637–654. [Google Scholar] [CrossRef]
Yildirim, T. Various Finite Difference Solutions of Option Pricing Models. Master’s Thesis, Yildiz Tecnical University, Istanbul, Turkey, 2016. [Google Scholar]
Widdicks, M.; Duck, P.W.; Andricopoulos, A.D.; Newton, D.P. The Black-Scholes equation revisited: Asymptotic expansions and singular perturbations. Math. Financ. 2005, 15, 373–391. [Google Scholar] [CrossRef]
Merton, R.C. Theory of optional option pricing. Bell J. Econ. Manag. Sci. 1973, 1, 141–183. [Google Scholar] [CrossRef]
Leland, H.E. Option pricing and replication with transaction cost. J. Financ. 1985, 40, 1283–1301. [Google Scholar] [CrossRef]
Boyle, P.P.; Vorst, T. Option replication in dicrete time with transaction cost. J. Financ. 1992, 47, 271–293. [Google Scholar] [CrossRef]
Kusuoka, S. Limit theorem on option replication with transaction costs. Ann. Appl. Probab. 1995, 5, 198–221. [Google Scholar] [CrossRef]
Barles, G.; Soner, H.M. Option pricing with transaction costs and nonlinear Black–Scholes equation. Financ. Stoch. 1998, 2, 799–812. [Google Scholar] [CrossRef]
Cox, C.; Ross, S.A. The valuation of options for alternative stochastic process. J. Financ. Econ. 1976, 3, 145–146. [Google Scholar] [CrossRef]
Tangman, D.Y.; Gopaul, A.; Bhuruth, M. Numerical pricing of options using high order compact finite difference schemes. J. Comput. Appl. Math. 2008, 218, 270–280. [Google Scholar] [CrossRef]
McCartin, B.J.; Labadie, S.M. Accurate and efficient pricing of vanilla stock options via the Crandall-Douglas scheme. Appl. Math. Comput. 2003, 143, 39–60. [Google Scholar] [CrossRef]
Ankudinova, J.; Ehrdardt, M. On the numerical solution of nonlinear Black–Scholes equations. Comput. Math. Appl. 2008, 56, 799–812. [Google Scholar] [CrossRef]
Company, R.; Jodar, L. A numerical method for European option pricing with transaction costs nonlinear equation. Math. Comput. Model. 2009, 50, 910–920. [Google Scholar] [CrossRef]
Jeong, D.; Yoo, M.; Kim, J. Finite difference method for the Black–Scholes equation without boundary conditions. Comput. Econ. 2018, 51, 961–972. [Google Scholar] [CrossRef]
Koleva, M.N.; Mudzimbabwe, W.; Vulkov, L.G. Fourth-order conmpact finite schemes for a parabolic ordinary system of European option pricing liquidity shock model. Numer. Algorithms 2017, 74, 59–75. [Google Scholar] [CrossRef]
Rao, S.C. Numerical solution of generalized Black–Scholes model. Appl. Math. Comput. 2018, 321, 401–421. [Google Scholar]
Sari, M.; Balacescu, A. Discrete algorithms for Black Scholes option pricing economic model. In Business and Applied Economics Book of Abstracts; Institute of Economic Sciences: Belgrade, Serbia, 2018; pp. 152–155. [Google Scholar]
Mashayekhi, S.; Hugger, J. Finite difference schemes for a nonlinear Black–Scholes Model with transaction cost and volatility risk. Acta Math. Univ. Comen. 2015, 84, 255–266. [Google Scholar]
Hout, K.; Valkov, R. Numerical solution of a two-asset option valuation PDE by ADI finite difference discretization. AIP Conf. Proc. 2015, 1648, 1648. [Google Scholar]
Hendrick, C.; Heuer, C.; Ehrhardt, M.; Günther, M. High-order ADI finite difference schemes for parabolic equations in the combination technique with application in finance. J. Comput. Appl. Math. 2017, 316, 175–194. [Google Scholar] [CrossRef]
Markolefas, S. Standard Galerkin formulation with high order Langrange finite elements for option markets pricing. Appl. Math. Comput. 2008, 195, 707–720. [Google Scholar]
Lin, S.B.; Chen, C.D. Applying the Model Order Reduction method to a European option pricing model. Econ. Model. 2013, 33, 533–536. [Google Scholar] [CrossRef]
Falloon, W.; Turner, D. The Evaluation of a Market, Managing Energy Price Risk; Risk Books: London, UK, 1999. [Google Scholar]
Kumar, A.; Tripathi, L.P.; Kadalbajoob, M.K. A numerical study of Asian option with radial basis functions based finite differences method. Eng. Anal. Bound. Elem. 2015, 50, 1–7. [Google Scholar] [CrossRef]
Rad, J.A.; Parand, K.; Ballestra, L.C. Pricing European and American options by radial basis point interpolation. Appl. Math. Comput. 2015, 251, 363–377. [Google Scholar] [CrossRef]
Vecer, J. Black–Scholes representation for Asian option. Math. Financ. 2014, 24, 598–626. [Google Scholar] [CrossRef]
Bayraktar, E.; Xing, H. Pricing Asian options for jump diffusion. Math. Financ. 2011, 21, 117–143. [Google Scholar] [CrossRef]
Sin, K. Numerical Methods for Derivative Pricing with Applications to Barrier Options. Master’s Thesis, University of Waterloo, Waterloo, ON, Canada, 2010. [Google Scholar]
Moon, K.S. Efficient Monte-Carlo algorithm for pricing barrier options. Commun. Korean Math. Soc. 2008, 2, 285–294. [Google Scholar] [CrossRef]
Zeng, X.C.; Guo, I.; Zhu, S.P. Pricing European Options Regime-Switching Assets: A comparative study of Monte Carlo and finite difference approaches. ANZIAM J. 2017, 2, 183–199. [Google Scholar] [CrossRef]
Koleva, M.N.; Vulkov, L.G. On splitting-based numerical methods for nonlinear models of European options. Int. J. Comput. Math. 2016, 93, 781–796. [Google Scholar] [CrossRef]
Hozman, J.; Tichy, T. On the impact of various formulations of the boundary condition within numerical option valuation by DG method. Filomat 2016, 30, 4253–4263. [Google Scholar] [CrossRef] [Green Version]
Anwar, M.N.; Andallah, L.S. A study on numerical solution of Black–Scholes model. J Math. Financ. 2018, 8, 372–381. [Google Scholar] [CrossRef]
Appadu, A.R. Comparison of some finite difference methods for the Black–Scholes equation. AIP Conf. Proc. 2018, 1978, 470098. [Google Scholar]
Jeong, D.; Yoo, M.; Yoo, C.; Kim, J. A hybrid Monte Carlo and finite difference method for option pricing. Comput. Econ. 2019, 53, 111–124. [Google Scholar] [CrossRef]
Cen, Z.; Chen, W. A HODIE finite difference scheme for pricing American options. Adv. Differ. Equ. 2019, 2019, 67. [Google Scholar] [CrossRef]
Gottlieb, S.; Shu, C.W.; Tadmor, E. Strong stability-preserving high-order time discretization methods. SIAM Rev. 2001, 43, 89–112. [Google Scholar] [CrossRef]
Dura, G.; Moşneagu, A.M. Numerical approximation of Black–Scholes equation. Analele Ştiintifice Ale Universitătii Al.I. Cuza Din Iai (S.N) Matematică 2010, 5, 39–64. [Google Scholar] [CrossRef]
Leentvaar, C.C. Numerical Solution of the Black–Scholes Equation with a Small Number of Grid Points. Master’s Thesis, Delft University of Technology, Delft, The Netherlands, 2003. [Google Scholar]
Sari, M.; Gurarslan, G.; Zeytinoglu, A. High-order finite difference schemes for the solution of the generalized Burgers-Fisher equation. Int. J. Numer. Methods Biomed. Eng. 2011, 27, 1296–1308. [Google Scholar] [CrossRef]

Figure 1. Prices of the linear European put option at

t = 0

.

Figure 1. Prices of the linear European put option at

t = 0

.

Figure 2. Prices of the linear European put option for different time values.

Figure 3. Price of the linear European put option.

Figure 4. Valuation of the European put option in both linear and nonlinear cases.

Figure 5. Difference

V_{d s} - V_{d s / 2}

.

Figure 5. Difference

V_{d s} - V_{d s / 2}

.

Figure 6. Barles and Soner model

(a = 0.02)

vs. the linear model.

Figure 6. Barles and Soner model

(a = 0.02)

vs. the linear model.

Table 1. Comparison between exact and numerical solutions of the linear and nonlinear European put option model for various stock price values at

t = 0

.

Table 1. Comparison between exact and numerical solutions of the linear and nonlinear European put option model for various stock price values at

t = 0

.

S	SSP-RK3 + FD6			Dura and Moşneagu [41]	Exact
	Linear Model	Nonlinear Model			Solution
		a = 0.02	a = 0.05		(Linear Model)
4.0	5.753096	5.753096	5.753096	5.753102	5.753100
8.0	1.902102	1.904440	1.915559	1.902102	1.902434
10.0	0.668906	0.673788	0.696650	0.668360	0.669390
16.0	0.005324	0.005498	0.006353	0.005419	0.005386
20.0	1.091776 × 10 $^{- 4}$	1.139578 × 10 $^{- 4}$	1.375594 × 10 $^{- 4}$	1.170806 × 10 $^{- 4}$	1.129336 × 10 $^{- 4}$

Table 2. Convergency results for the linear model.

N	M	${∥.∥}_{ds, \infty}$	Ratio ${∥.∥}_{ds, \infty}$	${∥.∥}_{ds, 2}$	Ratio ${∥.∥}_{ds, 2}$
51	101	0.018822		0.007153
101	1001	0.018035	1.043695	0.002887	2.477658
201	2001	0.002668	6.759370	6.547562 × 10 $^{- 4}$	4.409275
401	4001	4.036584 × 10 $^{- 4}$	6.609549	1.392947 × 10 $^{- 4}$	4.700510

Table 3.

ℓ^{2}

-error for the linear and nonlinear models.

Table 3.

ℓ^{2}

-error for the linear and nonlinear models.

N	M	Error (Linear Model)	Error (Nonlinear Model)
51	101	0.013100	0.012081
101	1001	0.004132	0.003932
201	2001	8.054754 × 10 $^{- 4}$	9.081600 × 10 $^{- 4}$

Table 4. CPU time (seconds).

N	M	CPU (Linear Model)	CPU (Nonlinear Model)	CPU (Nonlinear Model)
			a = 0.02	a = 0.05
51	101	0.136891	0.166548	0.171630
101	2001	3.223077	3.807224	3.830958
201	4001	12.112536	14.122904	14.242410
401	8001	63.284679	72.833058	72.978609

© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Gulen, S.; Popescu, C.; Sari, M. A New Approach for the Black–Scholes Model with Linear and Nonlinear Volatilities. Mathematics 2019, 7, 760. https://doi.org/10.3390/math7080760

AMA Style

Gulen S, Popescu C, Sari M. A New Approach for the Black–Scholes Model with Linear and Nonlinear Volatilities. Mathematics. 2019; 7(8):760. https://doi.org/10.3390/math7080760

Chicago/Turabian Style

Gulen, Seda, Catalin Popescu, and Murat Sari. 2019. "A New Approach for the Black–Scholes Model with Linear and Nonlinear Volatilities" Mathematics 7, no. 8: 760. https://doi.org/10.3390/math7080760

APA Style

Gulen, S., Popescu, C., & Sari, M. (2019). A New Approach for the Black–Scholes Model with Linear and Nonlinear Volatilities. Mathematics, 7(8), 760. https://doi.org/10.3390/math7080760

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

A New Approach for the Black–Scholes Model with Linear and Nonlinear Volatilities

Abstract

1. Introduction

2. Mathematical Formulation

3. The Solution Method

3.1. Spatial Discretization

3.2. Temporal Discretization

4. Numerical Results and Analysis

5. Conclusions and Recommendations

Author Contributions

Funding

Conflicts of Interest

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI