# Homotopy Perturbation Method for Fractional Black-Scholes European Option Pricing Equations Using Sumudu Transform

Authors: Asma Ali Elbeleze; Adem Kılıçman; Bachok M. Taib
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

The homotopy perturbation method, Sumudu transform, and He’s polynomials are combined to obtain the solution of fractional Black-Scholes equation. The fractional derivative is considered in Caputo sense. Further, the same equation is solved by homotopy Laplace transform perturbation method. The results obtained by the two methods are in agreement. The approximate analytical solution of Black-Scholes is calculated in the form of a convergence power series
with easily computable components. Some illustrative examples are presented to explain the efficiency and simplicity of the proposed method.


---

## Body

## 1. Introduction

Fractional differential equations have attracted much attention, recently, see, for instance [1–5]. This is mostly due to the fact that fractional calculus provides an efficient and excellent instrument for the description of many practical dynamical phenomena arising in engineering and scientific disciplines such as, physics, chemistry, biology, economy, viscoelasticity, electrochemistry, electromagnetic, control, porous media, and many more, see, for example, [6–9].Many partial differential equations of fractional order have been studied and solved. For example many researchers studied the existence of solutions of the Black-Scholes model using many methods, see [10–14].The homotopy perturbation method was first introduced and applied by He [15–17]. This method has been applied by many authors in many fields, for example, it is applied to nonlinear oscillator [18], nonlinear wave equation [19], nonlinear partial differential equations [20], integro-differential equation of fractional order [21], fuzzy differential equation [22], and other fields [23, 24]. Further homotopy perturbation methods are combined with Laplace transform to solve many problems such as one dimensional nonhomogeneous partial differential equations with a variable coefficient [25], Black-Scholes of fractional order [26], and parabolic partial differential equations [27].  The homotopy perturbation method coupled with Sumudu transform basically illustrates how Sumudu transform can be used to approximate the solutions of the linear and nonlinear differential equations by manipulating the homotopy perturbation method. In [28] Singh et al. studied the solution of linear and nonlinear partial differential equations by using the homotopy perturbation method coupled with Sumudu transform. Further, in [29] the authors proposed the homotopy perturbation method coupled with Sumudu transform to solve nonlinear fractional gas dynamics equation.The Black-Scholes equation is one of the most significant mathematical models for a financial market. It is a second-order parabolic partial differential equation that governs the value of financial derivatives. This Black-Scholes model for the value of an option is described by the following equation:(1)∂v∂t+σx22∂2v∂x2+r(t)x∂v∂x-r(t)v=0,00000(x,t)∈R+×(0,T),0<α≤1,
where v(x,t) is the European call option price at asset price x and at time t,T is the maturity,r(t) is the risk free interest rate, and σ(x,t) represents the volatility function of underlying asset. The payoff functions are
(2)vc(x,t)=max(x-E,0);vp(x,t)=max(E-x,0),
where vc(x,t) and vp(x,t) are the value of the European call and put options, respectively, E denotes the expiration price for the option, and the function max(x,0) gives the large value between x and 0.In this paper, we consider the following fractional Black-Scholes of the form(3)∂αv∂tα+σx22∂2v∂x2+r(t)x∂v∂x-r(t)v=0,00000(x,t)∈R+×(0,T),0<α≤1.In [29] Singh et al. used homotopy perturbation method coupled with Sumudu transform to solve fractional gas dynamics equation. The aim of this paper is to applied the homotopy perturbation method for fractional Black-Scholes equation by using He's polynomials and Sumudu transform.


## 2. Sumudu Transform

The Sumudu transform was first introduced and applied by Watugala [30] in (1998). For further details and properties of Sumudu transform see [31–34]. The Sumudu transform is defined over the set of functions:
(4)A={f(t):∃M,τ1,τ2>0,|f(t)|<Met/τj,=00ift∈(-1)j×[0,∞)}
by the following formula
(5)f-(u)=S[f(t);u]=:∫0∞‍f(ut)e-tdt,u∈(-τ,τ).
The existence and uniqueness was discussed in [35]. For further properties of Sumudu transform and its derivatives, see [36]. Some fundamental further established properties of Sumudu transform can be found in [31].Similarly, this new transform was applied to one-dimensional neutron transport equation [37]. In [34] Kılıçman et al. show that there is a strong relationship between Sumudu and other integral transforms. Further in [33] the Sumudu transform was extended to the distributions, and some of their properties were also studied in [38]. Recently Kılıçman et al. applied this transform to solve system of differentials equations, for more details see [34, 35, 37–39].


## 3. Basic Definitions of Fractional Calculus

In this section, we give some basic definitions and properties of fractional calculus theory which will be used in this paper.Definition 1.
The Riemann-Liouville fractional integral operator of orderα≥0 of a function f∈Cμ, μ≥-1 is defined as follows:
(6)Jαf(x)=1Γ(α)∫0x(x-t)α-1f(t)dt,000000000000000000α>0,t>0
in particular J0f(x)=f(x).For Riemann-Liouville fractional integral, one has(7)Jαxγ=Γ(γ+1)Γ(α+γ+1)xα+γ.Definition 2.
The Caputo fractional derivative off∈C-1m, m∈N is defined as follows:
(8)Dαf(x)=1Γ(m-α)∫0x(x-t)m-α-1fm(t)dt,0000000000000000000000000m-1<α≤m.Lemma 3.
Ifm-1<α≤m, m∈N, f∈Cμm, μ>-1 then the following two properties hold: (1)
Dα[Jαf(x)]=f(x),(2)
Jα[Dαf(x)]=f(x)-∑k=1m-1‍fk(0)(xk/k!).Definition 4.
The Mittag-Leffler functionEα(z) with α>0 is defined by the following series representation, valid in the whole complex plane:
(9)Eα(z)=∑0∞‍znΓ(αn+1).Definition 5.
The Sumudu transform of the Caputo fractional derivative is defined as follows [40]:
(10)S[Dtαf(t)]=u-αS[f(t)]-∑k=0m-1‍u-α+kf(k)(0+),0000000000000000000000000(m-1<α≤m).


## 4. Homotopy Perturbation Method

To illustrate the basic idea of this method, we consider the following nonlinear differential equation:(11)A(u)-f(r)=0,r∈Ω
with boundary conditions
(12)B(u,∂u∂n)=0,r∈Γ,
where A is a general differential operator, B is a boundary operator, f(r) is a known analytic function, and Γ is the boundary of the domain Ω.In general, the operatorA can be divided into two parts L and N, where L is linear, while N is nonlinear. Equation (11) therefor can be rewritten as follows:
(13)L(u)+N(u)-f(r)=0.
By the homotopy technique [41, 42] we construct a homotopy v(r,p):Ω×[0,1]→R which satisfies
(14)H(v,p)=(1-p)[L(v)-L(u0)]+p[A(v)-f(r)]=0000000000000000000000000000000000p∈[0,1],r∈Ω
or
(15)H(v,p)=L(v)-L(u0)+pL(u0)+p[N(v)-f(r)]=0,
where p∈[0,1] is an embedding parameter, and u0 is an initial approximation of (11) which satisfies the boundary conditions.From (14) and (15) we have
(16)H(v,0)=L(v)-L(u0)=0,H(v,1)=A(v)-f(r)=0.
The changing in the process of p from zero to unity is just that of v(r,p) from u0(r) to u(r). In topology this is called deformation, and L(v)-L(u0) and A(v)-f(r) are called homotopic.Now, assume that the solution of (14), (15) can be expressed as
(17)v=v0+pv1+p2v2+⋯.
Setting p=1 results in the approximate solution of (11):
(18)u=limp→1v=v0+v1+v2+⋯.


## 5. Homotopy Perturbation Method Coupled with Sumudu Transform

To illustrate the basic idea of this method, we consider the following nonlinear fractional differential equation:(19)Dtαu(x,t)+L[x]u(x,t)+N[x]u(x,t)0000=q(x,t),t>0,m-1<α≤m,
where Dtα=∂α/∂tα is the fractional Caputo derivative of the function u(x,t), L is the linear differential operator, N is the nonlinear differential operator, and q(x,t) is the source term.Now, applying the Sumudu transform on both sides of (19), we have
(20)S[Dtαu(x,t)]+S[L[x]u(x,t)+N[x]u(x,t)]0000=S[q(x,t)].
Using the differential property of Sumudu transform, we have
(21)S[u(x,t)]000=f(x)-uαS[L[x]u(x,t)+N[x]u(x,t)]00000+uαS[q(x,t)].
Operating with Sumudu inverse on both sides of (21)
(22)u(x,t)=Q(x,t)-S-1[uαS(L[x]u(x,t)+N[x]u(x,t))],
where Q(x,t) represents the term arising from the source term and the prescribed initial conditions.Now, applying the classical homotopy perturbation technique, the solution can be expressed as a power series inp as given below:
(23)u(x,t)=∑n=0∞‍pnun(x,t),
where the homotopy parameter p is considered as a small parameter (p∈[0,1]).We can decompose the nonlinear term as(24)Nu(x,t)=∑n=0∞‍pnHn(u),
where Hn are He’s polynomials of u0,u1,u2,…,un [43–45], and it can be calculated by the following formula:
(25)H(u0,u1,u2,…,un)=1n!∂n∂pn[N(∑n=0∞‍piui)]p=0,n=0,1,2,….
By substituting (23) and (24) and using HPM [15] we get
(26)∑n=1∞‍pnun(x,t)=Q(x,t)-p(S-1[uαS[L∑n=0∞‍pnun(x,t)+∑n=0∞‍pnHn(u)]]).
This is coupling of Sumudu transform and homotopy perturbation method using He’s polynomials. By equating the coefficient of corresponding power of p on both sides, the following approximations are obtained as
(27)p0:u0(x,t)=Q(x,t),p1:u1(x,t)=-S-1(uαS[L[x]u0(x,t)+H0(u)]),p2:u2(x,t)=-S-1(uαS[L[x]u1(x,t)+H1(u)]),p3:u3(x,t)=-S-1(uαS[L[x]u2(x,t)+H2(u)]),p3⋮
Proceeding in the same manner, the rest of the components un(x,t) can be completely obtained, and the series solution is thus entirely determined. Finlay we approximate the solution u(x,t) by truncated series
(28)u(x,t)=limN→∞∑n=0N‍un(x,t).
This series solutions generally converge very rapidly.


## 6. Examples

In this section, we discuss the implementation of the proposed method.Example 6.
We consider the following fractional Black-Scholes option pricing equation as follows:(29)∂αv∂tα=∂2v∂x2+(k-1)∂v∂x-kv,0<α≤1
subject to initial condition
(30)v(x,0)=max(ex-1,0).
Applying Sumudu transform on both sides of (29) subject to initial condition (30), we get
(31)S[v(x,t)]=max(ex-1,0)+uαS[vxx+(k-1)vx-kv].
Operating the inverse Sumudu transform on both sides in (31), we have
(32)v(x,t)=max(ex-1,0)-S-1[uαS(vxx+(k-1)vx-kv)].
Now, applying homotopy perturbation method
(33)∑n=0∞‍pnvn(x,t)=max(ex-1,0)-p(S-1[uαS[∑n=0∞‍pnHn(v)]]),
where
(34)Hn=vnxx+(k-1)vnx+kvn,n∈N.
Equating the corresponding power ofp on both sides in (38), we have
(35)p0:v0(x,t)=max(ex-1,0),p1:v1(x,t)=S-1(uαS[H0(v)])p1:v1(x,t)=-max(ex,0)(-ktα)Γ(α+1)p1:v1(x,t)=+max(ex-1,0)(-ktα)Γ(α+1),p2:v2(x,t)=S-1(uαS[H1(v)])p1:v1(x,t)=max(ex,0)(-ktα)2Γ(2α+1)p1:v1(x,t)++max(ex-1,0)(-ktα)2Γ(2α+1),p1:v1(x,t)⋮pn:vn(x,t)=S-1(uαS[Hn(v)])p1:v1(x,t)=max(ex,0)(-ktα)nΓ(nα+1)p1:v1(x,t)++max(ex-1,0)(-ktα)nΓ(nα+1).
So that the solution v(x,t) of the problem is given by
(36)v(x,t)=limp→1∑i=0∞‍piui(x,t)=max(ex-1,0)Eα(-ktα)+max(ex,0)(1-Eα(-ktα)),
where Eα(z) is Mittag-Leffler function in one parameter. For special case α=1, we get
(37)v(x,t)=max(ex-1,0)e-kt+max(ex,0)(1-e-kt),
which is an exact solution of the given Black-Scholes equation (29) for α=1.
The behaviour ofv(x,t) with respect to x and t when α=1 is given in Figure 1.Figure 1
The surface shows thev(x,t) for (29) with respect to x and t for α=1.Example 7.
We consider the following fractional Black-Scholes option pricing equation as follows:(38)∂αv∂tα+0.08(2+sinx)2x2∂2v∂x2+0.06∂v∂x-0.06v=0,00000000000000000000000000000000000<α≤1
subject to initial condition
(39)v(x,0)=max(x-25e-0.06,0).
Firstly, applying Sumudu transform on both sides of (38) subject to initial condition (39), we get
(40)S[v(x,t)]=max(x-25e-0.06,0)-uαS×[0.08(2+sinx)2x2vxx+0.06vx-0.06v].
Operating the inverse Sumudu transform on both sides in (40), we have
(41)v(x,t)=max(x-25e-0.06,0)-S-1×[uαS(0.08(2+sinx)2x2vxx+0.06vx-0.06v)].
Now, applying the homotopy perturbation method we have
(42)∑n=0∞‍pnvn(x,t)=max(x-25e-0.06-1,0)-p(S-1[uαS[∑n=0∞‍pnHn(v)]]),
where
(43)Hn=0.08(2+sinx)2x2vnxx+0.06vnx-0.06vn,n∈N.
Equating the corresponding power of p on both sides in (42), we have
(44)p0:v0(x,t)=max(x-25e-0.06,0),p1:v1(x,t)p1=S-1(uαS[H0(v)])p1=-x(-0.06tα)Γ(α+1)+max(x-25e-0.06,0)(-0.06tα)Γ(α+1),p2:v2(x,t)p1=S-1(uαS[H1(v)])p1=-x(-0.06tα)2Γ(2α+1)+max(x-25e0.06,0)(-0.06tα)2Γ(2α+1),p1⋮pn:vn(x,t)=S-1(uαS[Hn(v)])p1=-x(-0.06tα)nΓ(nα+1)0000+max(x-25ex,0)(-0.06-tα)nΓ(nα+1).
So that the solution v(x,t) of the problem is given by
(45)v(x,t)=limp→1∑i=0∞‍piui(x,t)=x(1-Eα(-0.06tα))+max(x-25e-0.06,0)Eα(-0.06tα).
This is the exact solution of the given option pricing equation (38). The solution of (38) at the special case α=1 is
(46)v(x,t)=x(1-e0.06t-1,0)+max(x-25e-0.06,0)e-0.06t.
The behaviour of v(x,t) with respect to x and t when α=1 is given in Figure 2.Figure 2
The surface shows thev(x,t) for (38) with respect to x and t for α=1.


## 7. Conclusion

In this paper, the homotopy perturbation Sumudu transform method (HPSTM) is successfully applied for getting the analytical solution of the fractional Black-Scholes option pricing equation. Two examples from the literature [26] are presented. The results of the illustrated examples are in agreement with the results of the method presented in [26]. In conclusion, HPSTM is a very powerful and efficient method to find approximate solutions as well as numerical solutions.



---

*Source: 524852-2013-05-08.xml*