Construction of Interval Shannon Wavelet and Its Application in Solving Nonlinear Black-Scholes Equation

Interval wavelet numerical method for nonlinear PDEs can improve the calculation precision compared with the common wavelet. A new interval Shannon wavelet is constructed with the general variational principle. Compared with the existing interval wavelet, both the gradient and the smoothness near the boundary of the approximated function are taken into account. Using the new interval Shannon wavelet, a multiscale interpolation wavelet operator was constructed in this paper, which can transform the nonlinear partial differential equations intomatrix differential equations; this can be solved by the coupling technique of the wavelet precise integration method (WPIM) and the variational iteration method (VIM). At last, the famous Black-Scholes model is taken as an example to test this new method. The numerical results show that this method can decrease the boundary effect greatly and improve the numerical precision in the whole definition domain compared with Yan’s method.


Introduction
The nonlinear Black-Scholes equations have been increasingly attracting interest over the last two decades, since they provide more accurate values by taking into account more realistic assumptions, such as the transaction costs, risks from an unprotected portfolio, large investor's preferences, or illiquid markets, which may have an impact on the stock price, the volatility, the drift, and the option price itself [1].In contrast with the traditional Black-Scholes equation, it is almost impossible to find the analytical solution of the nonlinear Black-Scholes equations.Most of the numerical schemes focus on the finite difference methods [2,3].As the matter of fact, the wavelet precise integration method (WPIM) is a simple and effective method for linear partial differential equations proposed by Mei et al. [4].For linear steady structural dynamic systems, their numerical results at the integration points are almost equal to those of the exact solution in machine accuracy.Yan proposed a wavelet precise integration method for the nonlinear Black-Scholes model [5] recently by combining the variational iteration method (VIM) and homotopy perturbation method (HPM) [6].
Since the definition domain of wavelet transformation is an infinite interval, the boundary effect would occur when being applied for resolving the engineering problems with bounded interval, for example, ordinary differential equations (ODEs).Consequently, it will decrease the precision and computational efficiency of the solution.Nevertheless, the boundary effect can be eliminated effectively by constructing an interval wavelet using numerical methods.There are several ways available to construct an interval wavelet.Mei et al. [7] proposed a general construction method of interval wavelet based on the restricted variational principle.According to this method, there exist some interval wavelets that have been constructed to solve PDEs in engineering [8][9][10][11].
The main purpose of this paper is to deduce a novel interval wavelet based on the restricted variational principle and construct an interval wavelet precise integration method (WPIM) for the nonlinear Black-Scholes model.According to WPIM, the nonlinear differential equation should be transformed to a system of ordinary differential equations with the multiscales wavelet interpolation operator [12], and then the nonlinear PDEs become a system of nonlinear differential equations.The matrix differential equation (MDE) can be solved almost exactly by the coupling technique of the precise integration method and the variational iteration method

Construction of Interval Shannon Wavelet
Based on Restricted Variational Principle Mei et al. [7] proposed a kind of construction method of the interval interpolation wavelet based on the restricted variational principle.It does not matter with the representation of the wavelet function.In this section, a new interval wavelet function with this method will be constructed.It takes into account the gradient and the smoothness of the approximated function.This is helpful to improve the calculation precision of the wavelet numerical method.The representation of Shannon wavelet [18][19][20][21] is based upon approximating the Dirac delta function as a bandlimited function and is given by The Shannon wavelet possesses many excellent numerical properties such as interpolating, relative sparse, and orthogonal properties.A perceived disadvantage of ( 1) is that it tends to zero quite slowly as || → ∞.A direct consequence of this is that there are a large number of grid points will contribute to the derivatives calculation of approximated function.It is for this reason that Hoffman et al. [22] have suggested using the Shannon-Gabor wavelet as follows: It has been proved that (2) can improve the localized and asymptotic behavior of the Shannon scaling function. is the width parameter (or called window size).However, the presence of the Gaussian window destroys the orthogonal properties possessed by the Shannon wavelet, effectively worsening the approximation to a Dirac delta function.So, the Shannon wavelet representation of the Dirac delta function is adopted in this paper, and it is shown that this representation ensures that the approach is identical to the weighted residual approach.The corresponding discrete formula of the Shannon scaling function is In the following, the construction of the interval Shannon wavelet is based on the general variational principle [7].Consider one-dimensional problem (), whose domain of definition Ω is [, ].Assuming that the number of discrete points is 2  +1,  ∈ , the th discrete point of variable  could be written as If the Shannon wavelet with interpolation property is employed as the trial function   , that is, any element in the matrices K and G  could be expressed, respectively, as where Ñ is the Lagrange basic function.
As the definition domain of wavelet transformation is a double infinite interval, the wavelet coefficients are large near the endpoints of the bounded signal, which increases the computational error.To obtain a high computational precision, the additional condition of scalar fonctionelle in boundary Γ could be changed to Thus where   is a known value.According to the general variational principle, we obtain where Correspondingly, (7) changed to the form The Substituting ( 13) into (10) results in where In the matrices G1 and G2,   and   could be calculated, respectively, as According to G   = 0 and the relational expression ( 13), ( 14) could be condensed in the following way: Therefore, the function () can be represented as where Subsequently, as interval interpolation basic functions, the interval Shannon wavelet can be obtained from (18) as follows: where  is the support domain of the wavelet function; that is, sup  = [−, ].
The interval wavelet constructed in [8] is different from (20), which is represented as follows: where where  is the amount of the external collocation points, the amount of discrete points in the definition domain is 2  + 1 ( ∈ ), and [ min ,  max ] is the definition domain of the approximated function.
It is obvious that ( 21) is just the Lagrange extrapolation representation.The parameter  appeared in ( 21) is relative to the smoothness and gradient of the function to be approximated.In theory, the increase of  can improve the approximate precision.But it was pointed out in [23] that the increase of  can also bring the increase of the condition number of the problem to be solved, and this can decrease the numerical precision greatly.In addition, the larger  can also bring more calculation.
Table 1 shows the maximum absolute errors of the numerical approximation to curves with different gradient near the boundary corresponding to different value of .It illustrates that  is not only related to the smoothness, but also to the gradient of original function near the boundary.For the same , the greater the original function's gradient is, the bigger the error is.To the different gradient original functions in Table 1, the reasonable  should be  ≥ 3,  ≥ 4,  ≥ 4, and  ≥ 8, respectively.In contrast, the interval wavelet proposed in this paper is not sensitive to the gradient of the function to be approximated, and the calculation amount is the same as the interval wavelet in [7] with  = 1.As the matter of fact, the proposed interval wavelet is constructed based on that the derivative function should be continuous, and the smaller  avoids the Gibbs phenomenon appeared in the Lagrange interpolator.This is the reason why it is not sensitive to the gradient value.

Black-Scholes Equation and Its Interval Wavelet Approximation
The famous Black-Scholes equation can be represented as where () denotes the underlying asset,  ∈ (0, ),  denotes the expiry date,  is the volatility (measures the standard deviation of the returns), and  is the riskless interest rate.
In general, the parameter  is taken as a constant since the transaction cost is taken as zero.Obviously, the  is not really a constant, and then we can obtain the nonlinear Black-Scholes equation as follows: where σ denotes a nonconstant volatility.In order to solve the problem, it is necessary to perform a variable transformation as follows: Substituting (25) into (24), the following equation can be obtained: where Initial condition is And the boundary conditions are According to the transformation relation ( 25), it is easy to understand that the point  = 0 is corresponding to the strike price  = .Obviously, the initial solution curve is smooth except around the position  = 0, where a sharp steep wave happened.So, an adaptive numerical method is suited to this problem.
According to the wavelet precision integration method, the nonlinear Black-Scholes equation can be approximately represented as d  (  , ) where  = 0, 1, 2, . . ., 2  .Let ] . (31) Then the system of (30) can be expressed as the matrix format: It is obvious that  0 is a constant matrix according to (32). 1 (  )   is the nonlinear term in (32), which can be linearized subsectionally in terms of first-order Taylor series expansion; that is, the nonlinear terms can be taken as linear ones in the time step (  , +1 ), and then (32) can be rewritten as where Letting  2 = d 1 /d = − diag(d  /d) 1 and combining with (32), the expression of  1 can be obtained as The precise time integration scheme for the linear ordinary differential equations (33) can be expressed as where and then solving (33) is changed into the calculation of matrix . the exponential matrix  can be calculated accurately by PIM as follows: Let Δ = /2  , where  is a positive integer (usually take  = 20, and then Δ = /1048576).Taking into account that  is a small time interval and Δ is a very smaller value, then which is the Taylor series expansion of exp(M 0 Δ).In order to calculate the matrix  more accurately, it is necessary to factorize the matrix  as After doing  times of factorization as above, a more accurate solution of  can be obtained.The calculation of the exponent matrix (ℎ) at different time steps is needed in solving nonlinear equations through iteration based on the precise integration method, and the algorithm of the matrix (ℎ) presented here can obtain all the matrices at different time steps for once.
The numerical solution of the nonlinear Black-Scholes equation is shown in Figure 1, where the riskless rate  = 0.1, the volatility  = 0.2, and strike price  = 100.With the increase of the parameter , the smoothness of the solution function near the strike price  becomes better and better, and the amount of the collocation points is decreasing correspondingly.This illustrates the performance of the adaptive wavelet numerical method.
Compared with the wavelet precision integration method employed in Yan's method [5], the amount of the collocation points in the interval wavelet precision integration method proposed in this paper is decreased greatly, especially around the boundary of the definition domain.This illustrates that the interval wavelet precision integration method can eliminate the boundary effect completely.It is obviously that the decrease of the collocation points can improve the efficiency of the algorithm.The performance of Mei's interval wavelet [4] relates to the external interpolation points (), and the choice scheme of  was not given in [3].It was pointed out in [23] that the bigger  usually introduces bigger condition number which can decrease the precision greatly.This can be illustrated in Figures 1(d) and 1(e).As  = 1, the boundary effect can be eliminated effectively.On the contrary, larger  almost cannot restrict it.
In addition, the interval wavelet can also improve the numerical precision compared with the common wavelet method.For convenience, we take the linear Black-Scholes equation as an example, which has exact analytical solution as follows: where is the call price,  is the underlying asset price,  is the strike price,  is the riskless rate,  is the maturity,  is the volatility, and ( 1 ) express the normal distribution.
All the comparisons in Table 2 are made qualitatively by comparing the calculation precision in the same time step and space mesh grid size.The first measure of error  1 is given by which provides a measure of the accuracy of the solution near the boundary.The second measure of error  2 is given by In Table 2, the parameter  is taken as 1 in Mei's method.The error  1 of Mei's method almost equals that of Yan's method, and the error e 2 of Mei's method is smaller evidently than Yan's method.This indicates that the interval wavelet method   proposed by Mei can eliminate the boundary effect to some extent as  = 1, but it cannot decrease the maximum of the absolute error evidently.In addition, with the increasing of the parameters  and , the errors of both Yan and Mei's method become larger and larger evidently.This can be explained by the condition number theory [23].On the contrary, the interval wavelet method proposed in this paper is not sensitive to the parameters  and , and the precision is the best among all above methods.

Conclusion
The solution of the nonlinear Black-Scholes equation has an evident sharp steep wave near the boundary, and this can introduce evident boundary effect in the wavelet approximation such as Yan's method.The interval wavelet approximation proposed in this paper makes sure that the derivative function near the boundary of the definition domain is continued, which is helpful to eliminate the boundary effect and improve the numerical precision in the whole definition domain.Compared with Mei's interval wavelet method, the interval wavelet proposed in this paper is not sensitive to the amount of the levels and the time step and does not have to choose the external collocation points.So, it is a simple and robust tool for solving the Black-Scholes equation.In fact, we can construct various interval wavelet functions aiming at different problems in engineering based on the general variational principle.

Figure 1 :
Figure 1: Evolution of the call option price with the parameter .

8 Mathematical
Problems in Engineering interpolation properties of Ñ and N imply that a and b are the vectors consisting of the approximate values of  in discrete points   and extended discrete points outside the domain of [, ], respectively.Typically, periodic, symmetric, and zero extension methods could be employed.Assuming that there are  = 2  discrete points in [, ], that is,  0 ,  1 , . . .,  2 , and  extended points in both outsides of the domain [, ], respectively, that is,  −−1 ,  − , . . .,  −1 and  +1 ,  +2 , . . .,  + , the function  can be extended outside of the definition domain using the symmetry method as below:

Table 2 :
The numerical precision comparison.