Finite Difference/Fourier Spectral for a Time Fractional Black–Scholes Model with Option Pricing

We study the fractional Black–Scholes model (FBSM) of option pricing in the fractal transmission system. In this work, we develop a full-discrete numerical scheme to investigate the dynamic behavior of FBSM. )e proposed scheme implements a known L1 formula for the α-order fractional derivative and Fourier-spectral method for the discretization of spatial direction. Energy analysis indicates that the constructed discrete method is unconditionally stable. Error estimate indicates that the 2 − α-order formula in time and the spectral approximation in space is convergent with order O(Δt2− α + N1− ), where m is the regularity of u and Δt and N are step size of time and degree, respectively. Several numerical results are proposed to confirm the accuracy and stability of the numerical scheme. At last, the present method is used to investigate the dynamic behavior of FBSM as well as the impact of different parameters.


Introduction
e classical option pricing model is proposed by Black and Scholes [1], which is based on the assumption that stocks and options are in an "ideal state" in the market and Samuelson's model [2]: where S be the the stock value, B(t) is the Brownian motion with the unit variance, and μ and σ are two constants. But in many cases, fractional Brownian motion is more accurate than integer order [3]. On the other hand, more and more diffusion processes were found to be non-Fickian [4,5], and the fractional order stochastic differential equation is considered as an extension of the stochastic differential equation. One view is that fractional order option trading equation is regarded as nonrandom growth process caused by Brownian motion. erefore, Jumarie [6] and Liang et al. [7] considered fractional Brownian motion in Samuelson's model equation: dS � μSdt + σSdB(t, α), for double barrier options are often biased when price changes are considered as fractal transmission systems. Chen et al. [9] revealed that it would be better to use fractional order Black-Scholes equation to explain the pricing in fractal transmission systems. ere is a lot of work in the modeling and calculation of fractional equations. Yang et al. [10,11] developed a new definition of fractional derivative. e advantage of this definition is that it does not contain singular kernel. Inc et al. [12,13] studied the isolated solutions of a class of fractional equations with Kerr law nonlinearity by Riccati-Bernoulli method. e exact dark optical and periodic singular soliton solution is obtained. Singhet al. [14][15][16] solved a series of fractional equations by using homotopy analysis technique and Laplace transform algorithm.
As we all know, as an effective method, L1 formula [17] has been widely used in the calculation of fractional differential equations. Langlands and Henry [18], Sun and Wu [19], and Lin and Xu [20] discussed the error estimate for L1 scheme. De Staelena and Hendybc [21] constructed a numerical method of fourth-order finite difference in space and 2 − α in time. Stability, uniqueness, and error estimates are analyzed. Zhang et al. [22] presented a discrete implicit finite scheme to solve time fractional Black-Scholes model. ey discussed the stability and error estimation of numerical schemes by Fourier analysis. Unfortunately, their analysis methods are local approach. Due to the importance of FBSM, it is necessary to reconstruct an efficient numerical method and analyze global stability and error estimates.
In this work, we will develop an efficient full-discrete scheme to approximate the Black-Scholes model with α-order fractional derivative. We apply L1 method to discretize the direction of time and Fourier-spectral method to discretize the direction of space. Using the energy analysis method, we discuss the stability and error estimate of the fully discrete numerical method. e detailed analysis shows that the scheme is unconditionally energy stable, and the error estimates indicate that our full-discrete scheme can achieve 2 − α-order accuracy in time and exponential accuracy in space direction. Finally, some numerical examples are conducted to support the theoretical claims. At the same time, the dynamic behavior of FBSM is studied by the proposed method.
We organize the rest of the paper as follows. Section 2 will briefly introduce the FBSM. In Section 3, we develop a time-discrete method for FBSM and then present its discrete energy law. In Section 4, we will study the error estimate of the full-discrete scheme. In Section 5, we present accuracy/ stability tests and numerous numerical examples to demonstrate the validity of the full-discrete method. In addition, we will discuss the properties of the solution of the FBSM. Some concluding remarks are given in Section 6.

Black-Scholes Model
In this work, we will consider the following time fractional Black-Scholes model: where V is the price of the option, S the price of the underlying asset, r the interest rate, and σ the volatility of the stock price, 0 < α ≤ 1; the time fractional derivative z α V(S, t)/zt α is defined by is is a linear parabolic partial differential equation which has been studied extensively.
We transform the problem to an initial value problem by using the time to mature t: � T − τ, and we then set x � ln S, V(S, τ) � u(e x , T − τ); we can rewrite (4) as where ξ � 1/2σ 2 , ω � r − ξ, and with the following boundary (barrier) and initial conditions, In order to solve the above model by numerical method, it is necessary to truncate the original unbounded region into a finite interval. erefore, we will consider problem (8) in bounded interval (0, 2π). en, we will study the following problem: Remark 1. In fact, one can choose homogeneous or inhomogeneous boundary conditions. It all depends on the actual option price. We have tested it, and it does not make any difference in actual numerical examples.

2 − αOrder Numerical Method
Here, we will develop the time-discrete method for equation (10). First, given a positive integer K, set Δt � T/K be the time step size, and denote τ n � nΔt, (0 ≤ n ≤ K) as the mesh point. en, we introduce an L1 method to discrete the Caputo fractional derivative of order α: Lemma 1 (see [20,23,24]). e coefficients of formula (13) satisfy the following properties: en, we can obtain the following time-discrete scheme: where a 0 � Δt α Γ(2 − α). It should be noted that if n � 0, we can rewrite the above equation as First of all, we have the following energy stability results for time-discrete (15).

Theorem 1.
e time-discrete scheme (15) is unconditionally stable. It satisfies the following energy dissipation law: Proof. When n � 0, computing the L 2 inner product of (16) with 2u 1 , we obtain It is easy to verify that the following formula is correct: us, Giving up some positive terms, we have Assume the following inequality holds: Next, we will show ‖u n+1 ‖ 2 ≤ ‖u 0 ‖ 2 is still valid. If j � n + 1, taking the L 2 inner product of (15) with 2u n+1 , we derive Note the fact that us, we get is yields (17).

Error Estimate for Full Discretization
In this part, we will study the Fourier-spectral method for the time-discrete method (15). First, we define S N as the polynomial space. Define π N : L 2 (Ω) ⟶ S N be the L 2 -projection operator which satisfies We have the following estimate [25]: en, we can develop the following full-discrete scheme: We now present the stability results of the fully discrete scheme (28).
Next, we begin to analyze the error estimates of the fulldiscrete scheme (28). Define the following error function: From [20,19,24,23], we know that R n+1 satisfies We also define the following error functions: e n u � π N u t n − u n N , e n u � u t n − π N u t n , e n u � e n u + e n u � u t n − u n N . (32) Lemma 2. For a 0 and b n , we have the following results: Proof.
Note that erefore, we obtain □ Remark 2. In [23,24], readers will also find similar parameter estimation. is result is very useful for us to analyze error estimate.

Theorem 3.
For the constructed numerical scheme (28), we have the following error estimate: Proof. For n � 0, we can write equation (28) as Subtracting (38) from (10) at τ 1 , we note that en, we have (40) (41) Dropping some positive terms, we find (43) Next, we will prove that it holds also for k � n + 1. Subtracting (28) from a reformulation of (10) at t n+1 , we find (45) Note that b n− j ≥ b n , thus Note that is ends the proof.

Numerical Examples
In this section, several numerical examples will be present to confirm the accuracy and applicability of the full-discrete scheme (28). We consider a rectangular computed domain of [0, 2π]. In order to better simulate the periodic boundary conditions, the Fourier-spectral method will be used to discretize space direction.

Verification of Convergence of Numerical Method.
First, in order to conduct a time accuracy test, an exact solution will be constructed to evaluate the convergence of the full-discrete scheme (28).
Example 1. We consider the following FBSM with α-order Caputo derivative: where It is easy to verify that the exact solution will be u � τ 2 sin x.
We set N � 128. e default values for the parameters are set as T � 1, ξ � 1, ω � 0, and r � 1. In Table 1, we show the temporal convergence orders of various time steps. As can be seen from Table 1, our full-discrete scheme is close to 2 − α-order accuracy in time, which is confirm with the result in eorem 3.
Fix N � 128, in Figure 1, we give L 2 error for different α. It is obvious that our numerical scheme has good convergence in time direction. Let T � 0.1, Δt � 10 − 6 , ξ � r � 0.5, ω � 0, and u 0 � cos 8x. Figure 2 shows that the full-discrete     Mathematical Problems in Engineering 7 scheme (28) has excellent convergence behavior in space direction.

Effect of Various Parameters.
is section is devoted to investigate the dynamic behavior of FBSM equation with different α. In the following numerical experiments, we fix ξ � 0.5, w � 0.1, r � 0.6, N � 64, Δt � 0.01, and u 0 � |sin x|. From Figures 3-5, we know that α has certain influence on the solution, as α increases, and the solution becomes smoother. In order to test the influence of ξ on the option price, we set r � 0.6, u 0 � sin(x/2) and let ξ change at the same time. Figure 6 shows that the parameter ξ has a significant effect on the price of options, and there will be an inflection point around x � 2.8. Finally, we investigated the influence of r on the option price, and the results are shown in Figure 7, it can be seen that when r increases, the option price also increases.

Conclusion
In this paper, a new full-discrete numerical method is developed to solve the FBSM. An efficient 2 − α-order and unconditionally energy stable method is constructed by combining the L1 approach in time and Fourier method in space direction. It is proved that the full-discrete converges to the order O(Δt 2− α + N − s + N − m ) globally. Numerical examples demonstrate the robustness and accuracy of the developed full-discrete method, numerically. Finally, we also study the properties of the solution of the FBSM.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.

Authors' Contributions
Juan He carried out an efficient numerical approach to time fractional Black-Scholes model. Aiqing Zhang helped to draft the manuscript. All authors read and approved the final manuscript.