The modified homotopy perturbation method and its application to the dynamics of price evolution in Caputo-fractional order Black Scholes model

Following a financial loss in trades due to lack of risk management in previous models from market practitioners, Fisher Black and Myron Scholes visited the academic setting and were able to mathematically develop an option pricing equation named the Black–Scholes model. In this study, we address the solution of a Caputo fractional-order Black–Scholes model using an analytic method named the modified initial guess homotopy perturbation method. Foremost, the classical Black Scholes model relaxed for European option style is generalized to be of Caputo derivative. The introduced method is established by coupling a power series function of arbitrary order with the renown He’s homotopy perturbation method. The convergence of the method is demonstrated using the fixed point theorem, and its methodology is illustrated by solving a generalized theoretical form of the fractional order Black Scholes model. The applicability of the method is proven by solving three different fractional order Black–Scholes equations derived from different market scenarios and its effectiveness is confirmed as feasible series of arbitrary orders that accelerate fast to the exact solution at an integer order were obtained. The computation of these results was carried out using Mathematica 12 software. Subsequently, the obtained outcomes were utilized in Maple 18 software to conduct a series of numerical simulations. These simulations aimed to analyze the influence of the fractional order on the dynamics of payoff functions regarding the share value as the option approached its expiration date under varying market constraints. In all three scenarios, the results showed that option values decrease as the expiration date approaches the integer order. Furthermore, the comparative outcomes reveal that Caputo fractional order derivatives control the flexibility of the classical Black–Scholes model because its payoff curve exhibits more sensitivity to changes associated with market characteristic parameters, such as volatility and interest rates. We propose that the results of this work should be examined and implemented by mathematicians and economists to better comprehend the influence of Caputo-fractional order derivatives in understanding the dynamics of option price evolution of financial assets.


Background
Options constitute one of the most frequently traded financial assets.The Black-Scholes model is an option pricing model that was proposed by Black and Scholes in [1,2] and Merton in [3,4].This model showed the significance of mathematics in finance.It is a second-order parabolic partial differential equation that is often used to compute the theoretical values of financial options based on market features such as stock price, risk rate, strike price, and volatility rate.
Many researchers have contributed to the understanding and implementation of the Black-Scholes model in mathematical finance.For example, the authors of [5] investigated the numerical solution of the Black-Scholes model.They explored the comprehensive history of Black-Scholes model's, discussed the parameters of the model and their functions, and provided a numerical approach capable of determining the model's outcome.Their findings, which showed convergence with the precise value, were presented and discussed.A numerical method for pricing double-barrier options in a timefractional Black-Scholes model was provided in [6] and an approach suitable for computing the exact solution of Black Scholes model was studied in [7].The stability and convergence of the proposed numerical scheme were shown using Fourier analysis, and the results verified the method's usefulness in solving the Black-Scholes model.The dynamics of option pricing were explored in [4] using a modified Black-Scholes-Merton model for option pricing.In their work, they constructed a conformable model that allows additional flexibility for markets, and an empirical application is carried out.
Fractional calculus has recently piqued the interest of contemporary academics, especially in the modeling of phenomena like diseases and financial theories.This is primarily because it introduces greater flexibility, imparts a long-term memory effect, and lends a more realistic dimension to these phenomena [8][9][10][11].Among the several fractional operators in use, notable ones include the Riemann-Liouville [12,13], Caputo [14,15], Caputo-Fabrizo [16], and Atangana-Baleanu [17] operators.These operators serve as versatile tools for exploring fractal and chaotic phenomena characterized by non-local kernels [18].An in-depth discussion of these properties, advantages, and disadvantages of these operators can be found in [19], which provides several significant results of studies concerning a fractional optimal control problem in systems with timedelay arguments.Similarly, the results provided in [20] on fractional optimal control for variable-order differential systems give supportive insights on the efficacy of mathematical modeling using these derivatives.The impact of Caputo-fractional operators on financial asset option pricing is of special interest in this study.The original Black-Scholes model, which was established using the classical derivative, may be insufficient for examining the behavior of financial asset option prices over time.It is imperative that, for better prediction of future asset price movements, attention be paid to the repeated patterns and trends observed in financial data.In [21], there is empirical evidence backing the repeated patterns and trends in financial markets as associated with fractional order derivatives.As a result, to generalize the derivative of the traditional Black-Scholes model, researchers usually employ a non-localized derivative, such as the Caputo fractional operator proposed by [13] and implemented in [22].
Mathematicians frequently research alternate techniques for obtaining the fractional-order Black-Scholes model's solutions because it is often challenging to obtain the exact solution of the model analytically.Numerous researchers have utilized a variety of numerical approaches and strategies to obtain an approximate solution to the fractional-order Black-Scholes model.The homotopy perturbation method proposed in [23] is a prime example of this approach.This method has been widely employed by different researchers ever since its emergence.For example, [24] modified and applied the homotopy perturbation method to obtain the approximate solution of fractional order Korteweg-De Vries equation.The method proved effective and their results showed rapid convergence to the exact solution.In [4], this same homotopy perturbation method was applied on a fractional order Black-Scholes model described in Caputo sense.It yields an analytic solution in the form of a convergent series with readily computed components.Also, a numerical computation of fractional Black Scholes equation arising in financial market was presented in [25].In their research, they applied the homotopy perturbation method and homotopy analysis method to solve fractional order Black Scholes Model and both methods proved to be highly effective as they both produce convergent series results.This homotopy perturbation approach was also applied in [26] to achieve an approximate analytical solution to a fractional-order Integro-differential equation.Their research indicated that the homotopy perturbation approach is a convergent and easily computable method for solving linear and non-linear fractional order models originating from application areas.A conformable fractional modified homotopy perturbation approach was employed in the study [27] to solve a novel European call option model.Their results showed that the proposed method is an efficient and powerful technique for finding approximate solutions to the fractional Black-Scholes models, which are considered conformable.
Although several studies on techniques for solving the fractional order Black-Scholes model have been published, the majority of the methods are not stable and are limited to different constraints.This motivates the necessity of providing an unconditionally stable and reliable method of solving the fractional-order Black-Scholes model in this study.Hence, we proposed a modified initial guess homotopy perturbation method [24,[28][29][30], which is a simple and effective analytical method to solve and investigate the behavior of the components associated with modern option pricing theory while utilizing less computational effort.
Thus consider the following fractional-order Black-Scholes differential equation.
For an European option, the solution of (1) as explained in [31,32] is defined to be f (S, t) = max f * , 0 .This rep- resents the worth (contingent claim) of the asset where f * = ±ς(S − E) and ς = ±1 respectively stands for the European call and put option prices.The stock price is denoted with S, E = Ke −r(T −t) embed parameters K, r, T, t representing the strike price, the riskless rate, the maturity date and the strike time of the option respectively.Expression (f * , 0) + signifies the maximum value between f * and zero.And the standard deviation (volatility rate) denoted by σ is a normally distributed parameter.
An alternative form of the classical order Black-Scholes model presented in [33,34] was derived using the following transformation: Such that Which yields Modifying (3) by generalizing the classical integer order operator using the Caputo fractional order operator we obtain: (1) where m = 2r σ 2 represents the equilibrium of the free interest rate and stock market volatility [35].

Basic definitions
We discuss some essential features of fractional calculus applied in this study here.
Then for t ∈ [0, T ] the following prop- erties hold: And Is the so called right fractional Caputo derivative which represents the future state of the function ψ(t).For more details on definitions and properties we refer to [19,20,36].
Theorem 1 Suppose there exist a contracting nonlinear mapping γ : δ → ζ specify on two Banach spaces, δ, ζ for all µ, ν ∈ δ , then �γ (µ for some µ 0 ∈ δ which converges to a unique fixed point γ [37]. The proof is defined by using mathematical induction on the contractive property of (C) such that Clearly, the continuity of γ is ensured by the contraction (C).Thus � = lim k→∞ µ k+1 = µ r .

Applied algorithm
In this section, we present the algorithm used to approximate solutions for the Black-Scholes model.

Define the Black-Scholes Model: • Define the input parameters: current stock price (S), strike price (K), time to expiration (T), risk-free interest rate (r), volatility (σ), and option type (call or put).
• Determine whether dividends are involved in the model.

Create an Initial Guess Functional for the Option
Price: • Construct an initial guess functional for the option price given by u(x, t)

If the solution converges:
-Output the approximate option price: -End.
This algorithm presents a systematic process for estimating option prices by employing the Black-Scholes model and the modified He's homotopy perturbation method, and it concludes upon achieving a converged result.Figure 1 illustrates the algorithm's sequential flow, as depicted below.

Homotopy perturbation method
The following non-linear differential equation may be used to explain the basic methodology of He's homotopy perturbation method: Subject to the boundary condition: Δ represents the general differential operator, B stands for the boundary operator, β(μ) is an analytic function and Γ is the boundary operator in the domain Ω.We can separate operator Δ into two parts: where ℓ and η denote the linear and nonlinear operator respectively.( 7) into ( 5) yields ( 5) We can construct a homotopy for (5): simplifying (9) yields; where p is an embedded parameter which undergo deformation process of changing from 0 to 1.
At p = 0, ( 8) We assume a power series solution as follows: Such that the approximate solution of (3) is: The convergence of this series has been proved in [23].Substituting ( 13) into (1212) and equating coefficients of equal powers of p,

Define the Black Scholes Model
The approximate result of each iteration is obtained by solving (13).

Modified initial guess homotopy perturbation method
In this section, we illustrate the methodology and implementation of the proposed modified initial guess homotopy perturbation technique to achieve the approximate solution of a generalized fractional order black Scholes model (FOBSM) in (1) given by: Subject to initial condition: We can construct a correctional functional using the modified initial guess homotopy perturbation method: (15) , and so on.As an assumed solution of (1), Eq. ( 18) must satisfy (1) with unique values of 1 , 2 , 3 , . . .n .Thus we evalu- ate these values, by obtaining the following derivatives: Evaluating ( 16) using (19),

Convergence of solution
We investigated the convergence of the modified initial guess homotopy perturbation technique on the solution obtained for (16).This convergence is only based on contraction of S such that: (25) kr n e −rT Ŵ(nα + 1) .
Clearly, ∞ n=0 r n t n n! = e rt , therefore �µ n − µ e � ∞ ≤ ke −rT ke rt − ke rt ∞ = 0 .Thus, the contraction of µ n (x, t) to μ e as n → ∞ generated by the modified initial guess homotopy perturbation method assures that the method converges unconditionally to the exact solution.

Applications
This section discusses the applications of the modified initial guess homotopy perturbation method in the context of the Black-Scholes equation under various market constraints.
Problem 1 Compute the European call option of a financial asset with drift r = 5  8 and volatility rate σ = 1 2 modelled with the following Caputo fractional order Black Scholes equation:

Subject to initial call option
We can compute the European call option value of the model by applying the modified initial guess homotopy perturbation method.Thus a correctional formula can be constructed such that Evaluating (31) using (30) yields: (29) Evaluating (29) with (33), at t = 0; Which yields:

Such that
To obtain the subsequent approximations, we can construct a homotopy for (29): simplifying (30) yields Substituting ( 23) into (37) and comparing the coefficients of p n , n ≥ 2 , we obtain: And so on.Evaluating (39) using ( 36) we obtain Applying integral operator I α on (42) yields the second approximation obtained as: (33) p 3 : .
Similar to second iteration we can evaluate (40) with (43) to obtain the third approximation: Applying integral operator I α on (44) yields the second approximation obtained as: Subsequent approximations are wisely computed using Mathematica 12 software package such that, Such that the approximate solution

Such that
The integer solution can be obtained when α = 1: Now, modified initial guess agrees that: Thus, at an integer order the exact payoff value for the European call option in ( 29) is: Problem 2 Consider a financial asset with a free interest and stock market volatility rate equilibrium specified (43) f 2 (s, t) = 4s 3 t 2α Ŵ(2α + 1) . (44) . (45) .
by m = 5.Modeled using the Caputo fractional order black Scholes model, Subject to initial condition By modified initial guess, we construct the correctional functional: Substituting (52) into (53), we obtain Equation ( 54) is applied in evaluating (51) such that we obtain: Evaluating (55) at t = 0, and simplifying, we obtain: Such that Constructing a homotopy for (51) to obtain subsequent approximations, Evaluating (58) with (23) and comparing the coefficients of p n , n ≥ 2 , we obtain: And so on.The second approximation can be obtained by evaluating (59) with (57): (51) (55) . (57) p 3 : Applying the integral operator I α on both sides of (62) yields the required approximate value: The third approximation can similarly be computed by evaluating (60) using ( 63) such that, Operator I α can be applied on both sides of (64) to obtain the third approximation such that Repeating the iterative procedure, the following results are achieved as successive approximations: Such that the approximate solution v(s, t) Taking the sum, we obtain The integer solution can be obtained at α = 1: The modified initial guess approach agrees that: Thus at an integer order, the exact solution of the Fractional order Black Scholes model in (51) is: (62) This indeed is the solution of (4) presented in [31] at a value of m = 5 provided that in their obtained result given as: (70) Thus at t = 0, (772 becomes: Solving for λ 1 yields Such that the first approximation yields; (74) We can construct a homotopy for (70) to obtain proceeding approximations such that: Substituting ( 23) into (77) and comparing the coefficients of p n , n ≥ 2 , we obtain: (77) p 3 : And so on.The second approximation result is obtained by evaluating (78) by (76): Applying operator I α on both sides of (81), we obtain Repeating the same process for (79), to obtain u 3 (x, t), (81) Also applying I α on both sides of (83) yields the third approximation obtained as Subsequent approximations are computed such that Such that the approximate solution u(s, t) The sum can be taken such that the payoff function of the put option is: The integer solution can be obtained at α = 1: The modified initial guess approach agrees that: Thus, the exact

Results
The numerical table of results obtained by solving the fractional order Black Scholes differential equation using the proposed approach is presented in this section.In addition, to develop an insight of the price evolution of options, we run numerical simulations of the solved (83) examples using the computational software program Maple 18, and the results are graphically presented.Table 1 shows the numerical results of the MIGHPM, exact solution and the absolute errors.It could be observed that the results produced by the MIGHPM are in good agreement with the exact solution.Although slight deviations of the approximate results from the exact solution were observed as the value of t increases from 0 to 1.For example, the absolute error observed when t = 0, are zero 0 but at t = 0.2 the absolute errors are in order of × 10 −12 to × 10 −10 ; and at t = 1 the absolute error increases to an order which exist in the range of × 10 −4 to × 10 −3 .
Table 2 shows that the absolute error drastically increases as t increases.We observed that, the absolute error is zero when t = 0 and the error ranges from an order of × 10 −13 to × 10 −2 as t progresses from 0 to 1.
Table 3 shows the effectiveness of the MIGHPM in solving fractional order Black Scholes Model.It could be observed that the maximum absolute error recorded exist in the range of × 10 −18 and the absolute error is zero at t = 0 and t = 1.

Discussion
The primary objective of this research is to ascertain the impact of Caputo fractional order derivatives on the price dynamics of financial assets.To achieve this, we examined the behavior of payoff functions concerning fluctuations in stock prices within the range of 0 < α < 1.The results of our investigation are presented as payoff curves.
Figure 2 displays the maturity graphs of option payoffs for problem 1 at various Caputo fractional order levels.Notably, for the assets under consideration, we observed a significant decline in option payout values as the maturity date approached.Specifically, in Fig. 2A, when α = 0.25, the payoff value was 40, but this sharply dropped to 3 as α increased to 0.5, and subsequently decreased further to 1 when α = 0.75, reaching parity with the exact option price.It is worth highlighting that the call value of the financial asset was found to be at its lowest when α took an integer value, signifying a market phenomenon when all constraints were fully enforced.Figure 2E provides a graphical comparison of the call option payoffs, clearly illustrating the diminishing trend as α increases.
Figure 2F, on the other hand, compares the accuracy of our proposed method with the exact solution.Our findings indicate that the proposed method effectively addressed the problem, yielding an approximation result that unconditionally converges to the exact solution.Notably, our outcomes align with the results presented by Esekon et al. [38] in a prior study, where an analytical method was employed to compute the results of the classical Black-Scholes model, underscoring the efficacy of MIGHPM.Additionally, Fig. 3 depict the dynamic behavior of price evolution for the financial asset modeled in problem 2 using the Black-Scholes equation.In Fig. 3A, the payoff value is 6000 when α = 0.25, but this reduces to 2500 as α increases to 0.5 in Fig. 3B, further declining to 300 as α reaches 0.75 in Fig. 3C, ultimately converging to 4 in the exact solution presented in Fig. 3D.These results underscore the importance of exercising caution when dealing with option prices as they approach their expiration dates.Figure 3E offers a comparative plot of these factors, while Fig. 3F reaffirms the effectiveness of our method in solving the fractional order Black-Scholes equation.These findings also align with those of Edeki et al. [33], which introduced a Black-Scholes equation with constant elasticity variance for stock option evaluation.Moving on to problem 3, we present the put option's payout values in Fig. 4. Interestingly, we observed that the asset's put option similarly relies on the derivative of Caputo's fractional order.Figure 4A-C depict a substantial reduction in the asset's payoff value as α approaches the exact solution, in line with the results from [22] that employed the Lagrange multiplier method.Notably, Fig. 4C, D record negative put option payoffs as the option nears its expiration time, indicating investor losses, and the necessity to pay the absolute value of the payoff.Finally, we present graphs comparing the precise and approximate solutions of our proposed method with those generated by the built-in RK4 of Maple 18 in Fig. 4F.Despite some variations in the intervals at which RK4 Maple 18 provides option call values, we found that all three outcomes generally agree.This consistency an implicit finite difference approach for solving the resultant BS-PDE.Theoretical analysis of their research essentially confirms that Caputo fractional order derivatives are the most effective way to simulate option premiums, which deviate significantly from what is predicted by classical differential operator theory.Another significant aspect to keep in mind is that the conventional BS models are known to generate option premium curves with identical curves, which may not accurately reflect the market dynamics.However, the option premium curves produced by the model developed using Caputo fractional derivatives are extremely sensitive to changes in related market characteristics, such as volatility and rates, etc.

Recommendation
Future research could be carried out to investigate the influence of other non-local fractional order operator such as Caputo-Fabrizo [14] and Atangana-Baleanu [15] derivatives on Black Scholes model derived from real market situations using the method applied in this study to better comprehend and forecast the dynamics of option price evolution.

Check for convergence of solution Evaluate the Black Scholes Model using the Initial guess function to obtain the first approximation Construct a He's homotopy perturbation correctional functional Evaluate the correctional functional with the first approximation to obtain the next approximation Compute subsequent approximations using the previous approximation Converges Output the approximate result No Yes Not Convergent End Fig
. 1 Flow chat of the algorithm

Table 1
Numerical Results of problem one at an exact value α = 1

Table 2
Numerical Results of problem Two at an exact value α = 1