Loading [MathJax]/jax/element/mml/optable/BasicLatin.js
Research article

Fat tails arise endogenously from supply/demand, with or without jump processes

  • We show that the quotient of Levy processes of jump-diffusion type has a fat-tailed distribution. An application is to price theory in economics, with the result that fat tails arise endogenously from modeling of price change based on an excess demand analysis resulting in a quotient of arbitrarily correlated demand and supply whether or not jump discontinuities are present. The assumption is that supply and demand are described by drift terms, Brownian (i.e., Gaussian) and compound Poisson jump processes. If P1dP/dt (the relative price change in an interval dt) is given by a suitable function of relative excess demand, (DS)/S (where D and S are demand and supply), then the distribution has tail behavior F(x)xζ for a power ζ that depends on the function G in P1dP/dt=G(D/S). For G(x)|x|1/q one has ζ=q. The empirical data for assets typically yields a value, ζ˜=3, or ζ[3,5] for many financial markets.

    Many theoretical explanations have been offered for the disparity between the tail behavior of the standard asset price equation and empirical data. This issue never arises if one models price dynamics using basic economics methodology, i.e., generalized Walrasian adjustment, rather than the usual starting point for classical finance which assumes a normal distribution of price changes. The function G is deterministic, and can be calibrated with a smaller data set. The results establish a simple link between the decay exponent of the density function and the price adjustment function, a feature that can improve methodology for risk assessment.

    The mathematical results can be applied to other problems involving the relative difference or quotient of Levy processes of jump-diffusion type.

    Citation: Gunduz Caginalp. Fat tails arise endogenously from supply/demand, with or without jump processes[J]. AIMS Mathematics, 2021, 6(5): 4811-4846. doi: 10.3934/math.2021283

    Related Papers:

    [1] Yinlin Ye, Hongtao Fan, Yajing Li, Ao Huang, Weiheng He . An artificial neural network approach for a class of time-fractional diffusion and diffusion-wave equations. Networks and Heterogeneous Media, 2023, 18(3): 1083-1104. doi: 10.3934/nhm.2023047
    [2] Kexin Li, Hu Chen, Shusen Xie . Error estimate of L1-ADI scheme for two-dimensional multi-term time fractional diffusion equation. Networks and Heterogeneous Media, 2023, 18(4): 1454-1470. doi: 10.3934/nhm.2023064
    [3] Narcisa Apreutesei, Vitaly Volpert . Reaction-diffusion waves with nonlinear boundary conditions. Networks and Heterogeneous Media, 2013, 8(1): 23-35. doi: 10.3934/nhm.2013.8.23
    [4] Yaxin Hou, Cao Wen, Yang Liu, Hong Li . A two-grid ADI finite element approximation for a nonlinear distributed-order fractional sub-diffusion equation. Networks and Heterogeneous Media, 2023, 18(2): 855-876. doi: 10.3934/nhm.2023037
    [5] Zhe Pu, Maohua Ran, Hong Luo . Effective difference methods for solving the variable coefficient fourth-order fractional sub-diffusion equations. Networks and Heterogeneous Media, 2023, 18(1): 291-309. doi: 10.3934/nhm.2023011
    [6] Feiyang Peng, Yanbin Tang . Inverse problem of determining diffusion matrix between different structures for time fractional diffusion equation. Networks and Heterogeneous Media, 2024, 19(1): 291-304. doi: 10.3934/nhm.2024013
    [7] Panpan Xu, Yongbin Ge, Lin Zhang . High-order finite difference approximation of the Keller-Segel model with additional self- and cross-diffusion terms and a logistic source. Networks and Heterogeneous Media, 2023, 18(4): 1471-1492. doi: 10.3934/nhm.2023065
    [8] Aslam Khan, Abdul Ghafoor, Emel Khan, Kamal Shah, Thabet Abdeljawad . Solving scalar reaction diffusion equations with cubic non-linearity having time-dependent coefficients by the wavelet method of lines. Networks and Heterogeneous Media, 2024, 19(2): 634-654. doi: 10.3934/nhm.2024028
    [9] Caihong Gu, Yanbin Tang . Global solution to the Cauchy problem of fractional drift diffusion system with power-law nonlinearity. Networks and Heterogeneous Media, 2023, 18(1): 109-139. doi: 10.3934/nhm.2023005
    [10] Verónica Anaya, Mostafa Bendahmane, David Mora, Ricardo Ruiz Baier . On a vorticity-based formulation for reaction-diffusion-Brinkman systems. Networks and Heterogeneous Media, 2018, 13(1): 69-94. doi: 10.3934/nhm.2018004
  • We show that the quotient of Levy processes of jump-diffusion type has a fat-tailed distribution. An application is to price theory in economics, with the result that fat tails arise endogenously from modeling of price change based on an excess demand analysis resulting in a quotient of arbitrarily correlated demand and supply whether or not jump discontinuities are present. The assumption is that supply and demand are described by drift terms, Brownian (i.e., Gaussian) and compound Poisson jump processes. If P1dP/dt (the relative price change in an interval dt) is given by a suitable function of relative excess demand, (DS)/S (where D and S are demand and supply), then the distribution has tail behavior F(x)xζ for a power ζ that depends on the function G in P1dP/dt=G(D/S). For G(x)|x|1/q one has ζ=q. The empirical data for assets typically yields a value, ζ˜=3, or ζ[3,5] for many financial markets.

    Many theoretical explanations have been offered for the disparity between the tail behavior of the standard asset price equation and empirical data. This issue never arises if one models price dynamics using basic economics methodology, i.e., generalized Walrasian adjustment, rather than the usual starting point for classical finance which assumes a normal distribution of price changes. The function G is deterministic, and can be calibrated with a smaller data set. The results establish a simple link between the decay exponent of the density function and the price adjustment function, a feature that can improve methodology for risk assessment.

    The mathematical results can be applied to other problems involving the relative difference or quotient of Levy processes of jump-diffusion type.



    Time-fractional partial differential equations (TFPDEs) have attracted considerable attention due to their ability to model memory and nonlocal properties. It has been established that TFPDEs are crucial mathematical and physical models for describing a wide range of anomalous phenomena and complex systems in the fields of natural science and engineering [1,2,3,4,5]. Successful applications of TFPDEs include signal processing [6], Powell-Eyring fluid [7], fluid mechanics [8,9], and robotics [10]. The multi-term TFPDEs have been found useful in describing many complex natural phenomena, such as magnetic resonance imaging (MRI) [11], viscoelastic mechanical models [12], and oxygen delivery through capillaries to a tissue (see [13] and the references therein). The multi-term time-fractional mixed sub-diffusion and diffusion wave equations belong to an important class of multi-term TFPDEs that have three major subclasses: the time-fractional sub-diffusion equations (TFSDEs), whose fractional order belongs to (0,1); the time fractional diffusion-wave equations (TFDWEs), whose fractional order belongs to (1,2); and TFDEs whose fractional order belongs to (0,2). Some researchers studied the analytical solutions of these equations, such as the authors of [14] who studied the solution of multi-term TFDWEs using the method of separation of variables. In [15], the authors used the Laguerre polynomials for obtaining the approximate solution of 1D TFDWE. The authors of [16] obtained the analytic solution of the multi-term space-time fractional advection-diffusion equation using a method based on spectral representation of fractional Laplacian operator. In most cases, the analytical solution of the these equations is provided in terms of special functions such as the multivariate Mittag-Leffler function and the Fox H-function [17,18], which are very complex and challenging to evaluate. Therefore, to solve these equations, numerical methods would be preferable, particularly in situations where analytical solutions are not present. Various numerical schemes have been developed to solve TFPDEs. Liu et al.[19] developed a first-order finite difference approach with stability constraints for both time and spatial derivatives. The multi-term time-fractional mixed sub-diffusion with a variable coefficient was introduced by Zhao [20] using the finite element method. The solution of the wave-diffusion equation based on the Legendre spectral method and finite difference method for discretization of space and time derivatives was established in [21]. Shen and Gu [22] developed two novel finite difference schemes for 2D time fractional mixed diffusion and wave-diffusion equations. Agrawal [23] utilized the Laplace transform and finite sine transform for obtaining the solution of the wave-diffusion equation. In [24], the authors coupled spatial extrapolation method with the Crank-Nicholson method to obtain the solution to the fractional diffusion problem. In [25], a fully discrete spectral method was employed for the solution of a novel multi-term time-fractional mixed diffusion and diffusion-wave equation. The mesh-less method based on radial basis function was employed by [26] for the approximate solution of multi-term time-fractional mixed sub-diffusion. In [27], the authors coupled the dual reciprocity method and the improved singular boundary method for multi-term fractional wave-diffusion equations. Ye et al. [28] solved the 2D and 3D multi-term time-space fractional diffusion equations via a series expansion. Li et al.[29] analyzed the temporal asymptotic behavior with well-posedness for the solution of the multi-term time-fractional diffusion equation. The authors of [30] employed the spectral method for the numerical solution of the two-dimensional multi-term time-fractional mixed sub-diffusion. Sun [31] examined the simultaneous inversion of the potential term and the fractional orders in a multi-term time-fractional diffusion equation. Spectral methods are among the most popular techniques for discretizing spatial variables in partial differential equations, which have the reputation of producing extremely accurate approximations of sufficiently smooth solutions [32,33]. Over time, spectral methods have been widely used in many fields, such as quantum mechanics [34], fluid dynamics [35], weather forecasting [36], and heat conduction [37], due to their high-order accuracy. These techniques haven't, however, been applied to many problems where the finite-difference and finite-element approaches are still frequently used because of some disadvantages. One disadvantage is that using spectral methods to discretize partial differential equations results in the solution of large systems of linear or nonlinear equations that require full matrices. On the other hand, finite-difference and finite-element approaches produce sparse matrices, which can be handled by suitable techniques to significantly reduce the computational complexity. Another drawback of spectral methods is that they face difficulties for problems defined on complex geometries. Although, some authors have attempted to use spectral methods in complex geometries[38,39]. In this work, our aim is to obtain the approximate solution of the multi-term time-fractional mixed sub-diffusion using the pseudospectral method in space coupled with the Laplace transform in time. The pseudospectral method is an easy method to implement and offers low computational cost with high accuracy. Also, it has been proven that for time discretization, the Laplace transform is suitable for solving various initial and boundary value problems and can be utilized as an alternative to the classical time stepping techniques [40]. In finite difference time stepping techniques, the accuracy can be obtained for extremely short time steps, which results in high computational cost and crucial stability constraints. A large number of valuable works are available in the literature, which couples Laplace transform in time with another method in space (see [41,42] and references therein). In this article, we consider a two-dimensional multi-term time-fractional mixed sub-diffusion and diffusion-wave equation of the following form [43]:

    pk=1ϱ1,kαku(ˉξ,t)tαk+ϱ2Dtu(ˉξ,t)+qr=1ϱ3,rβru(ˉξ,t)tβr+ϱ4u(ˉξ,t)=ϱ5Δu(ˉξ,t)+ϱ6γΔu(ˉξ,t)tγ+G(ˉξ,t),ˉξΓ,t[0,1] (1.1)

    with boundary conditions

    LBu(ˉξ,t)=h(ˉξ,t),ˉξΓ,t[0,1] (1.2)

    and initial condition

    u(ˉξ,0)=g1(ˉξ),ut(ˉξ,0)=g2(ˉξ),ˉξΓ (1.3)

    where G(ˉξ,t),g1(ˉξ),g2(ˉξ),h(ˉξ,t) are given continuous functions, ϱ1,k>0,ϱ3,r>0,ϱl>0, l=2,4,5,6, the coefficients are not simultaneously equal to zero, 1<α1<α2<...<αp<2,0<β1<β2<...<βq<1,0<γ<1, ˉξ=(ξ,ζ),Dt=t,Δ=2ξ2+2ζ2,LB is the boundary operator, ΓR2 is the domain, Γ is its boundary, and αu(ˉξ,t)tα, βu(ˉξ,t)tβ, and γu(ˉξ,t)tγ are fractional derivatives in Caputo's sense, defined as in [44].

    αu(ˉξ,t)tα=1Γ(nα)t0nu(ˉξ,ν)νn(tν)αn+1dν,n1<αn

    In literature, there are many physical processes whose characteristics cannot be described by utilizing the time fractional sub-diffusion equations or the time fractional wave diffusion equation independently. In order to avoid such limitations, the multi-term time-fractional mixed sub-diffusion can be utilized as an effective mathematical model [43]. These models are particularly beneficial at representing anomalous diffusion processes and capturing the properties of media, including power-law frequency dependence [30]. Additionally, they are proficient in the modeling of various types of viscoelastic damping, modeling the unsteady flow of a fractional Maxwell fluid [45], and describing the behavior of an Oldroyd-B fluid, which has been used to simulate the response of many dilute polymeric liquids [22].

    In this section, the Laplace transform is employed on the multi-term time fractional mixed sub-diffusion and diffusion-wave equation in order to reduce the problem into an equivalent time-independent problem. The Laplace transform of u(ˉξ,t) is denoted and defined as in [44]

    ˆu(ˉξ,s)=0estu(ˉξ,t)dt

    The Laplace transform of Caputo's derivative αu(ˉξ,t)tα of fractional order α(m1,m] is defined as in [44]

    L{αu(ˉξ,t)tα}=sαˆu(ˉξ,s)m1j=0sαj1u(j)(ˉξ,0)

    Applying Laplace transform on Eqs (1.1)–(1.3), we have

    pk=1ϱ1,k(sαkˆu(ˉξ,s)sαk1u(ˉξ,0)sαk2ut(ˉξ,0))+ϱ2(sˆu(ˉξ,s)u(ˉξ,0))+qr=1ϱ3,r(sβrˆu(ˉξ,s)sβr1u(ˉξ,0))+ϱ4ˆu(ˉξ,s)=ϱ5Lˆu(ˉξ,s)+ϱ6(sγLˆu(ˉξ,s)sγ1Lu(ˉξ,0))+ˆG(ˉξ,s)

    and

    LBˆu(ˉξ,s)=ˆh(ˉξ,s)

    which implies

    (pk=1ϱ1,ksαkI+ϱ2sI+qr=1ϱ3,rsβrI+ϱ4Iϱ5Lϱ6sγL)ˆu(ˉξ,s)=ˆH(ˉξ,s) (1.4)
    LBˆu(ˉξ,s)=ˆh(ˉξ,s) (1.5)

    where ˆH(ˉξ,s)=pk=1ϱ1,ksαk1g1(ˉξ)+pk=1ϱ1,ksαk2g2(ˉξ)+ϱ2g1(ˉξ)+qr=1ϱ3,rsβr1g1(ˉξ)ϱ6sγ1Lg1(ˉξ)+ˆG(ˉξ,s), L=Δ, and I is the identity operator. In the proposed approach, first we discretize the linear differential operator L via the pseudospectral method, then the full-discrete system (Eqs (1.4) and (1.5)) is solved in Laplace transform space. Finally the solution of Eqs (1.1)–(1.3) is obtained via inversion of the Laplace transform. Our solution method is highly suitable for parallel computations; it belongs to the class of parallel in time methods for TFPDEs (see [46]). The pseudospectral method is discussed in next section.

    The pseudospectral method is an accurate and extremely precise approach for the numerical solution of TFPDEs. The pseudospectral method uses the Chebyshev points for collocation that are nonzero over the entire domain, whereas finite element methods use basis functions that are nonzero only on small subintervals [35]. In the pseudospectral method, the expansion of solutions is characterized by some global basis functions, such as Lagrange's interpolation polynomials. The concept of interpolation and differentiation matrices [47] is very important in the pseudospectral method. To provide a comprehensive description of the pseudospectral method based on Lagrange's interpolation polynomial basis, the main aspects of the differentiation matrices are reviewed in the next section.

    In the pseudospectral method, the solution is considered over [1,1] and the data {(ξj,ˆu(ξj))}mj=0 is interpolated by the Lagrange polynomial (LP) ηj(ξ) of degree m [33,48]

    Im(ξ)=mj=0ηj(ξ)ˆuj

    where ηj(ξ) is the LP at the point ξj(j=0,1,...,m), which is

    ηj(ξ)=(ξξ0)...(ξξj1)(ξξj+1)...(ξξm)(ξjξ0)...(ξjξj1)(ξjξj+1)...(ξjξm) (1.6)

    where ˆuj=ˆu(ξj). The Chebyshev nodes are defined as

    ξj=cos(jπm),j=0,1,...,m (1.7)

    and are used to discretize the domain [1,1]. The 1st derivative ˆu(ξ)ξ is approximated as

    ˆu(ξ)ξΦmˆu

    and the matrix Φm has elements in the following form

    [Φm]i,j=ηj(ξi),i,j=1,2,...,m

    The non-diagonal elements of [Φm]i,j have the following form

    [Φm]i,j=αjαi(ξiξj),ij

    where α1j=mij(ξiξj), and the elements on diagonal entries are expressed as

    [Φm]i,j=mj=0,ji[Φm]i,j,i=0,1,2,..,m

    The elements of the μth-order derivative of Φμm are analytically obtained as

    [Φμm]i,j=ημj(ξi),i,j=1,2,...,m

    In [47,49], more accurate evaluation of the differentiation matrices is elaborated. For the differentiation matrices, Welfert [47] deduced a handy recursion relation as follows

    [Φμm]kj=μξkξj(αjαk[Φ(μ1)m]kk[Φ(μ1)m]kj),kj

    Chebyshev points are remarkable for their ability to enable spectral convergence when appropriate analytic assumptions are made.

    Let us consider the square domain Γ=[1,1]2. Generally, the point in Γ is denoted by ˉξ and is expressed as

    ˉξij=(cos(πim),cos(πjm)),i,j=0,1,2,...,m

    The LPs associated to ˉξij can be written as

    ηij(ˉξ)=ηi(ξ)ηj(ζ)

    where ηij(ˉξij)=δij. The 2nd-order derivatives of the LPs Eq (1.6) are given as

    2ηij(ˉξrs)ξ2=ηi(ξr)ηj(ζs)=[Φ2m]riδjs,
    2ηij(ˉξrs)ζ2=ηi(ξr)ηj(ζs)=δri[Φ2m]sj

    where Φ2m represents the second order differentiation matrix based on Chebyshev points. Applying the operator L on the LPs Eq (1.6) based on the points {ˉξrs} is

    L(ηij(ˉξrs))=([Φ2m]riδjs+δri[Φ2m]sj)

    Consequently, the approximation L via the pseudospectral method is obtained as

    LM=(ImΦ2m+Φ2mIm)

    The conditions in Eq (1.5) are incorporating by considering the interpolation matrix LM and considering all points ˉξ. Furthermore, the rows of LM in correspondence with boundary nodes are replaced with unit vectors that have a one in accordance with the diagonal elements of LM. Hence, the boundary conditions LBˆu(ˉξ,s)=ˆh(ˉξ,s) in Eq (1.5) will be implemented directly [33]. By rearranging the columns and rows of the matrix LM, the following block matrix is obtained

    LΓ=[WR0I]

    where W and I are the nonzero square blocks of order (mmB)×(mmB) and (mB×mB), respectively. However, mB represents the boundary points. The solution of the system Eqs (1.4) and (1.5) can be attained by solving the linear block system

    LΓˆu(ˉξ,s)=[ˆH(ˉξ,s)ˆh(ˉξ,s)]

    where ˆH and ˆh collect the values at interior and boundary collocation points, respectively. In the final step, the inverse Laplace transform is used to obtain the solution of Eqs (1.1)–(1.3).

    The interpolation operator is based on Chebyshev collocation points Eq (1.7) and Lagrange polynomials Eq (1.6) as follows:

    Im:C(Γ)Pm,Im(u)=mj=0u(ξj)ηj(ξ).

    The steps proposed in [50] will be followed for the construction of the interpolation polynomial error bound, for a constant Mm, which satisfies the following estimate

    Im(u)Mmu,uC[1,1] (1.8)

    Additionally,

    Im(u)=u,uPm

    It is possible to show that for Chebyshev interpolation,

    Mm=ln(1+m)π2+1(m+1)

    Based on m, the constant of stability increases very slowly. For uCm+1[1,1], we have the following bound [50]

    uIm(u)12mΓ(m+2)u(m+1) (1.9)

    Theorem 1. [50] If Eqs (1.8) and (1.9) hold, and uCm+1[1,1], we have

    u()Im(u)()2(m)(M()m+1)Γ(m+2)u(m+1),=0,1,2,...,m (1.10)

    with the stability constant

    M()m=1Γ(+1)[MmΓ(m+1)Γ(m+1)]

    Proof: By using Eqs (1.9) and (1.10) for Eq (1.1) in the one-dimension case, the governing operator is written as L=2ξ2; thus, the error estimate is given as

    E=(pk=1ϱ1,kDαktu+ϱ2Dtu+qr=1ϱ3,rDβrtu+ϱ4uϱ5Luϱ6Dγt(Lu))(pk=1ϱ1,kDαktIm(u)+ϱ2Im(Dtu)+qr=1ϱ3,rDβrtIm(u)+ϱ4Im(u)ϱ5LIm(u)ϱ6Dγt(LIm(u)))=pk=1ϱ1,kDαkt(uIm(u))+ϱ2(Dt(uIm(u)))+qr=1ϱ3,rDβrt(uIm(u))+ϱ4(uIm(u))ϱ5L(uIm(u))ϱ6Dγt(L(uIm(u)))pk=1ϱ1,kDαkt(uIm(u))+ϱ2Dt(uIm(u))+qr=1ϱ3,rDβrt(uIm(u))+ϱ4(uIm(u))+ϱ5L(uIm(u))+ϱ6Dγt(L(uIm(u)))pk=1|ϱ1,k|Dαkt(uIm(u))+|ϱ2|Dt(uIm(u))+qr=1|ϱ3,r|Dβrt(uIm(u))+|ϱ4|uIm(u)+|ϱ5|uξξIm(u)ξξ+|ϱ6|Dγt(uξξIm(u)ξξ)Epk=1|ϱ1,k|Dαkt(uIm(u))+|ϱ2|Dt(uIm(u))+qr=1|ϱ3,r|Dβrt(uIm(u))+|ϱ4|12mΓ(m+2)u(m+1)+|ϱ5|2(M(2)m+1)Γ(m)(12)m1u(m+1)+|ϱ6|Dγt(uξξIm(u)ξξ)

    since the time derivatives are accurately computed. Therefore, the error bound of Dαt(uIm(u)),Dt(uIm(u)),Dβt(uIm(u)) and uIm(u) are of the same order and the error bound of Dγt(uξξIm(u)ξξ) and uξξIm(u)ξξ are of the same order. Thus, the following stability bound is obtained:

    ENu(m+1)

    where the constant N is determined by computing the coefficients of u(m+1). A similar stability estimate can be obtained for two-dimensional problems by using the tensor product interpolation operators.

    In this section, we utilize numerical inverse Laplace methods to transform the pseudospectral method solution ˆu(ˉξ,s) from the Laplace domain to the time domain u(ˉξ,t) as follows:

    u(ˉξ,t)=12πiρ+iρiesτˆu(ˉξ,s)ds,ρ>ρ0 (1.11)

    Here, ˆu(ˉξ,s) needs to be inverted and the converging abscissa is ρ0 and ρ>ρ0. Thus, the open half plane Re(s)<ρ contains all the singularities of ˆu(ˉξ,s). Our goal is to approximate the integral defined in Eq (1.11). In most cases, the integral defined in Eq (1.11) is quite difficult to be evaluated analytically; therefore, a numerical method needs to be employed. The integral defined in Eq (1.11) may be evaluated by several numerical algorithms available in literature. Every approach has a distinct application and is appropriate for a certain type of function. In this article, we use the two popular inversion algorithms, the improved Talbot's method [51] and the Stehfest's method [52], which are presented as follows.

    Here, we utilize the Talbot's method to approximate u(ˉξ,t)

    u(ˉξ,t)=12πiρ+iρiestˆu(ˉξ,s)ds=12πiΛestˆu(ˉξ,s)ds,Re(s)>0 (1.12)

    where Λ is a suitably chosen contour. In the Talbot's method, the numerical quadrature is applied to the integral in Eq (1.12). The trapezoidal and midpoint rule are the two effective rules used in conjunction with the deformation of contour [53]. The purpose of contour deformation is to handle the exponential factor. Particularly, the path of integration can be deformed to a Hankel contour, i.e., a contour whose real part starts at in the 3rd quadrant and winds around all the singularities of the transform function going again to in the 2nd quadrant. The exponential component decays rapidly on such contour, making the Bromwich integral appropriate for approximation and using the trapezoidal and midpoint rule. Cauchy's theorem justifies such deformation, provided that the contour remains in the region where the transformed function ˆu(ˉξ,s) is analytic. Furthermore, some mild restrictions are required in the left complex plane on the decay of the transformed function ˆu(ˉξ,s)[54,55]. We consider the contour in the following parametric form as[51]

    Λ:s=s(ω),πωπ,Res(±π)=

    We have

    s(ω)=mTtϖ(ω),ϖ(ω)=η1+η2ωcot(η3ω)+η4iω (1.13)

    where η1,η2,η3, and η4 are to be selected by the user. Plugging Eq (1.13) in Eq (1.12), we have

    u(ˉξ,t)=12πiππes(ω)tˆu(ˉξ,s)(s(ω))s(ω)dω (1.14)

    Midpoint rule with uniform step h=2πmT is utilized to approximate Eq (1.14) as

    uApp(ˉξ,t)1mTimTk=1es(ωk)tˆu(s(ωk))s(ωk),ωk=π+(k12)h (1.15)

    The error analysis of the improved Talbot's approach is based on the following theorem.

    Theorem 1.1. [51] Let ωk be defined as in Eq (1.15). Let g:ΞC be an analytic function in the set

    Ξ={ωC:π<Reω<π,ands<Imω<r}

    when r,s>0, then

    ππg(ω)dω2πmTmTj=1g(ωk)=Ψ(ϕ)+Ψ+(ψ)

    where

    Ψ+(ϕ)=12(π+iϕπ+π+iϕπ+iϕ+ππ+iϕ)(1+itan(mTω2))g(ω)dω

    and

    Ψ(ψ)=12(πiψπ+πiψπiψ+ππiψ)(1itan(mTω2))g(ω)dω

    0<ϕ<r,0<ψ<s,andmTeven if mT is an odd number; we can replace tan(mTω2) with cot(mTω2) if g(ω) is real valued, that is, g(ˉω)=¯g(ω) and if r and s can be taken to be equal, then

    Ψ(ψ)=Ψ+(ψ)+Ψ(ψ)=Reπ+iψπ+iψ(1+itan(mTω2))g(ω)dω

    By examining the complex tangent function's behavior, this may be bounded as

    |Ψ(ψ)|4πCexp(rmT)1

    The above process has been done for even mT and C and r are some positive constants. A similar approach can be used for an odd mT.

    The optimal values of the parameters included in Eq (1.13) can be utilized to find the optimal contour of integration, which is necessary to obtain the best results. The authors of [51] have obtained the optimal values of the parameters as follows

    η1=0.61220,η2=0.50170,η3=0.64070,andη4=0.26450

    The corresponding error estimate is given as

    Errest=| uApp(ˉξ,t)u(ˉξ,t)|=O(e(1.358)mT)

    One of the most effective and straightforward techniques for Laplace transform inversion is the Gaver-Stehfest approach. The latter part of the 1960s saw its design. Because of its simplicity and efficacy, it has become more and more popular in a variety of fields, including computational physics, finance, economics, and chemistry. The basis of the Gaver-Stehfest approach is the series of Gaver approximants, as found by Gaver [56]. Since the convergence of the Gaver approximants was essentially logarithmic, acceleration was necessary. A linear acceleration approach was provided by Stehfest [52] using the Salzer acceleration method. Using a series of functions, the Gaver-Stehfest method approximates u(ˉξ,t) as

    uApp(ˉξ,t)=ln2tmSk=1θkˆu(ˉξ,ln2tk) (1.16)

    where θk are given as

    θk=(1)mS2+kmin(k,mS2)n=k+12nmS2(2n)!(mS2n)!n!(n1)!(kn)!(2nk)! (1.17)

    Solving the system Eqs (1.4) and (1.5) for the corresponding Laplace parameters s=ln2tk,k=1,2,3,...,mS, the approximate solution uApp(ˉξ,t) of the problem in Eq (1.1) can be obtained via Eq (1.16). The Gaver-Stehfest algorithm has a few noteworthy qualities, such as: (ⅰ) u(ˉξ,t) are linear in the context of values of ˆu(ˉξ,s); (ⅱ) the values of ˆu(ˉξ,s) are required merely for real value of s; (ⅲ) the procedure of determining the coefficients is fairly effortless; (ⅳ) for constant functions, this approach results in significantly precise approximations, i.e., if uc,thenuAppc for all mS1. In literature, this methodology has been employed by many researchers in [57,58], in which it is revealed that this strategy converges promptly to uApp(ˉξ,t) (given u(ˉξ,t) is non-oscillatory).

    Convergence

    The convergence of uApp(ˉξ,t) has been derived by the author in [57]. The results are based on the following theorem.

    Theorem 1.2. Let u:(0,)R be a locally integrable function. Let s>0, define the Laplace transform ˆu(ˉξ,s), and let uApp(ˉξ,t) be the numerical solution as given by Eq (1.16).

    1. uApp(ˉξ,t) converges given u(ˉξ,t) near t.

    2. Let for some real number ϕ and 0<ψ<0.25,

    ψ0|u(ˉξ,tlog2(1/2+η))+u(ˉξ,tlog2(1/2η))2ϕ|η1dη<

    then uApp(ˉξ,t)ϕ as mS+.

    3. Let u(ˉξ,t) be of bounded variation near t, then

    uApp(ˉξ,t)u(ˉξ,t+0)+u(ˉξ,t0)2,asmS+

    Corollary 1.3. Using the assumptions of the above theorem, if

    u(ˉξ,t+η)u(ˉξ,t)=O(|η|ϑ)

    η, and some ϑ, then uApp(ˉξ,t)u(ˉξ,t) as mS+.

    Moreover, the authors in [59] conducted a number of experiments to determine how parameters affected the numerical scheme's correctness. Their conclusions are as follows: "For ν1 significant digits, select a mS positive integer 2.2ν1: After setting the system precision to ν2=1.1mS, compute θi,1imS for a given mS using Eq (1.17). Next, compute the uApp(ˉξ,t) in Eq (1.16) for the provided transformed function ˆu(ˉξ,s) and the argument t." These conclusions indicate that the error is 10(ν1+1)uApp(ˉξ,t)u(ˉξ,t)u(ˉξ,t)10ν1, where mS=2.2ν1 [40].

    In this section, three test problems are considered to assess the effectiveness and accuracy of the Laplace transformed pseudospectral method. We compute the maximum absolute error (Err), the relative error (Err2), and the root mean square error (Errrms) between the numerical and the exact solutions, which are defined as follows:

    Err=max
    {Err}_{2} = \sqrt{\frac{\sum^{m}_{k = 1}(u(\mathit{\boldsymbol{\bar{\xi}}}_k,t)-u_{App}(\mathit{\boldsymbol{\bar{\xi}}}_k,t))^2}{\sum^{m}_{k = 1}(u(\mathit{\boldsymbol{\bar{\xi}}}_k,t))^2}},
    {Err}_{rms} = \sqrt{\frac{\sum^{m}_{k = 1}(u(\mathit{\boldsymbol{\bar{\xi}}}_k,t)-u_{App}(\mathit{\boldsymbol{\bar{\xi}}}_k,t))^2}{m}}

    where u_{App}(\mathit{\boldsymbol{\bar{\xi}}}_k, t) is the approximate solution and u(\mathit{\boldsymbol{\bar{\xi}}}_k, t) is the exact solution.

    Problem 1

    We consider Eq (1.1) with exact solution

    \begin{equation*} u(\xi, \zeta ,t) = (t^{3}+1) \sin (\pi \xi) \sin(\pi \zeta) \end{equation*}

    and the source term is

    \begin{split} \mathcal{G}(\xi, \zeta ,t)= &\sin(\pi \xi)\sin(\pi \zeta) \bigg[ \varrho_{1,1} \frac{6}{\Gamma(4-\alpha)}t^{3-\alpha}+3\varrho_{2} t^{2}+ \varrho_{3,1} \frac{6}{\Gamma(4-\beta)}t^{3-\beta}\\ &+(\varrho_{4}+2\varrho_{5} \pi^{2})(t^{3}+1)+2\varrho_{6}\pi^{2} \frac{6}{\Gamma(4-\gamma)}t^{3-\gamma} \bigg] \end{split}

    where \varrho_{1, 1} = \varrho_{2} = \varrho_{3, 1} = \varrho_{4} = \varrho_{5} = \varrho_{6} = 1 and (\xi, \zeta) \in [-1, 1]^{2}. The initial-boundary conditions are obtained from the analytical solution. The Err_{2}, Err_{\infty}, and Err_{rms} errors of the proposed method for problem 1 obtained via the Stehfest's and the Talbot's methods by using various values of m, \; m_{T}, \; m_{S}, \; \alpha, \; \beta, \; \gamma, and t = 1 are presented in Tables 1 and 2, respectively. The numerical solution to problem 1 is shown in Figure 1. Figure 2 presents a comparison of error norms Err_{2}, Err_{\infty}, and Err_{rms} computed using the Stehfest's and the Talbot's methods for various values of m_{S} and m_{T} with \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m = 30, and t = 1 . The comparison of Err_{2}, Err_{\infty}, and Err_{rms} computed by the Stehfest's and the Talbot's methods for different values of t with \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m_{S} = 14, \; m_{T} = 30, and m = 30 is shown in Figure 3. The comparison of Err_{2}, \; Err_{\infty}, and Err_{rms} computed using the Stehfest's and the Talbot's methods versus \alpha with \beta = 0.7, \; \gamma = 0.3, \; m = 30, \; m_{S} = 14, \; m_{T} = 30, and t = 1 is presented in Figure 4. The comparison of Err_{2}, Err_{\infty}, and Err_{rms} computed via the Stehfest's and the Talbot's methods versus \beta with \alpha = 1.5, \gamma = 0.7, m = 30, m_{S} = 14, m_{T} = 30, and t = 1 is presented in Figure 5. The comparison of Err_{2}, Err_{\infty}, and Err_{rms} computing by the Stehfest's and the Talbot's methods versus \gamma with \alpha = 1.5, \beta = 0.7, m = 30, m_{S} = 14, m_{T} = 30, and t = 1 is presented in Figure 6. The plot of Err_{2} obtained via the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and \beta \in [0, 1] with \gamma = 0.3, \; t = 1, \; m_{S} = 16, m_{T} = 30, and m = 30 is shown in Figure 7. In Figure 8, the plots of Err_{\infty} obtained using the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and \gamma \in [0, 1] with \beta = 0.7, \; t = 1, \; m_{S} = 16, \; m_{T} = 30, and m = 30 are shown. Similarly, Figure 9 presents the plots of Err_{rms} obtained using the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and t \in [0, 1] with \beta = 0.7, \gamma = 0.3, m_{S} = 16, \; m_{T} = 30, and m = 30. It is observed that the accuracy of the Stehfest's method steadily decreases for m_{S}\geq 14 . From the results presented in tables and figures, we conclude that the proposed method is accurate, stable, and efficient. The numerical results demonstrate that the Talbot's method is more accurate than the Stehfest's method.

    Table 1.  The Err_{2}, Err_{\infty}, Err_{rms} computed by the Stehfest's method for problem 1.
    (\alpha, \beta, \gamma) m m_{S} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    (1.3, 0.2, 0.9) 22 14 4.7216 \times 10^{-04} 9.6894 \times 10^{-05} 2.1462 \times 10^{-05} 0.199250
    25 5.3662 \times 10^{-04} 1.0751 \times 10^{-04} 2.1465 \times 10^{-05} 0.326975
    28 7.5644 \times 10^{-04} 1.0128 \times 10^{-04} 2.7016 \times 10^{-05} 0.556057
    30 10 6.9931 \times 10^{-02} 5.9791 \times 10^{-03} 2.3310 \times 10^{-03} 0.575029
    12 1.5653 \times 10^{-03} 1.3408 \times 10^{-04} 5.2178 \times 10^{-05} 0.663883
    14 7.8051 \times 10^{-04} 1.3790 \times 10^{-04} 2.6017 \times 10^{-05} 0.755536
    (1.6, 0.5, 0.7) 22 14 4.0428 \times 10^{-04} 4.4012 \times 10^{-05} 1.8376 \times 10^{-05} 0.180055
    25 4.8390 \times 10^{-04} 4.4246 \times 10^{-05} 1.9356 \times 10^{-05} 0.324530
    28 6.4225 \times 10^{-04} 8.5592 \times 10^{-05} 2.2938 \times 10^{-05} 0.576450
    30 10 6.9931 \times 10^{-02} 5.9791 \times 10^{-03} 2.3310 \times 10^{-03} 0.581491
    12 1.5680 \times 10^{-03} 1.3409 \times 10^{-04} 5.2266 \times 10^{-05} 0.688463
    14 7.4051 \times 10^{-04} 1.3304 \times 10^{-04} 2.4684 \times 10^{-05} 0.719391
    (1.9, 0.8, 0.5) 22 14 3.8530 \times 10^{-04} 4.2171 \times 10^{-05} 1.7514 \times 10^{-05} 0.184399
    25 4.5504 \times 10^{-04} 4.3708 \times 10^{-05} 1.8202 \times 10^{-05} 0.329158
    28 5.2746 \times 10^{-04} 6.0001 \times 10^{-05} 1.8838 \times 10^{-05} 0.556485
    30 10 6.9931 \times 10^{-02} 5.9791 \times 10^{-03} 2.3310 \times 10^{-03} 0.514627
    12 1.5660 \times 10^{-03} 1.3402 \times 10^{-04} 5.2201 \times 10^{-05} 0.630183
    14 6.5840 \times 10^{-04} 1.3505 \times 10^{-04} 2.1947 \times 10^{-05} 0.760430

     | Show Table
    DownLoad: CSV
    Table 2.  The Err_{2}, Err_{\infty}, Err_{rms} computed by the Talbot's method for problem 1.
    (\alpha, \beta, \gamma) m m_{T} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    (1.3, 0.2, 0.9) 22 30 1.4423 \times 10^{-10} 3.8758 \times 10^{-11} 6.5561 \times 10^{-12} 1.431363
    25 2.2738 \times 10^{-10} 6.1642 \times 10^{-11} 9.0953 \times 10^{-12} 2.823492
    28 3.0634 \times 10^{-10} 5.6744 \times 10^{-11} 1.0941 \times 10^{-11} 4.853607
    30 20 5.1049 \times 10^{-06} 4.3647 \times 10^{-07} 1.7016 \times 10^{-07} 4.409938
    25 1.1297 \times 10^{-08} 9.6706 \times 10^{-10} 3.7658 \times 10^{-10} 5.576429
    30 3.7908 \times 10^{-10} 8.4039 \times 10^{-11} 1.2636 \times 10^{-11} 6.435613
    (1.6, 0.5, 0.7) 22 30 8.8800 \times 10^{-11} 1.6472 \times 10^{-11} 4.0364 \times 10^{-12} 1.432723
    25 1.3312 \times 10^{-10} 3.5030 \times 10^{-11} 5.3248 \times 10^{-12} 2.644643
    28 1.9449 \times 10^{-10} 3.9329 \times 10^{-11} 6.9462 \times 10^{-12} 4.498548
    30 20 5.1049 \times 10^{-06} 4.3647 \times 10^{-07} 1.7016 \times 10^{-07} 4.512210
    25 1.1277 \times 10^{-08} 9.6519 \times 10^{-10} 3.7589 \times 10^{-10} 5.250997
    30 2.5877 \times 10^{-10} 8.3971 \times 10^{-11} 8.6256 \times 10^{-12} 6.315079
    (1.9, 0.8, 0.5) 22 30 6.3703 \times 10^{-11} 1.7695 \times 10^{-11} 2.8956 \times 10^{-12} 1.465253
    25 9.4756 \times 10^{-11} 2.4410 \times 10^{-11} 3.7902 \times 10^{-12} 2.624655
    28 1.5070 \times 10^{-10} 5.0718 \times 10^{-11} 5.3821 \times 10^{-12} 4.675835
    30 20 5.1049 \times 10^{-06} 4.3647 \times 10^{-07} 1.7016 \times 10^{-07} 4.469602
    25 1.1281 \times 10^{-08} 9.6487 \times 10^{-10} 3.7602 \times 10^{-10} 5.300393
    30 2.0811 \times 10^{-10} 4.6520 \times 10^{-11} 6.9370 \times 10^{-12} 6.311030

     | Show Table
    DownLoad: CSV
    Figure 1.  Numerical solution of problem 1 obtained by proposed scheme.
    Figure 2.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus m_{S} and for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus m_{T} for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, t = 1, and m = 30 .
    Figure 3.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus t for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, m_{S} = 14, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus t for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, m_{T} = 30, and m = 30 .
    Figure 4.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus \alpha for \beta = 0.7, \; \gamma = 0.3, \; m_{S} = 14, \; t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus \alpha for \beta = 0.7, \; \gamma = 0.3, \; m_{T} = 30, \; t = 1, and m = 30 .
    Figure 5.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus \beta for \alpha = 1.5, \gamma = 0.3, m_{S} = 14, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus \beta for \alpha = 1.5, \gamma = 0.3, m_{T} = 30, t = 1, and m = 30 .
    Figure 6.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Stehfest's method for \alpha = 1.5, \; \beta = 0.7, \; m_{S} = 14, \; t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Talbot's method for \alpha = 1.5, \; \beta = 0.7, \; m_{T} = 30, \; t = 1, and m = 30 .
    Figure 7.  (a) Plot of Err_{2} versus \alpha and \beta obtained by Stehfest's method for \gamma = 0.3, t = 1, m_{S} = 16, and m = 30. (b) Plot of Err_{2} versus \alpha and \beta obtained by Talbot's method for \gamma = 0.3, t = 1, m_{T} = 24, and m = 30 .
    Figure 8.  (a) Plot of Err_{\infty} versus \alpha and \gamma obtained by Stehfest's method for \beta = 0.7, \; t = 1, \; m_{S} = 16, and m = 30. (b) Plot of Err_{\infty} versus \alpha and \gamma obtained by Talbot's method for \beta = 0.7, \; t = 1, \; m_{T} = 24, and m = 30 .
    Figure 9.  (a) Plot of Err_{rms} versus \alpha and t obtained by Stehfest's method for \beta = 0.7, \gamma = 0.3, m_{S} = 16, and m = 30. (b) Plot of Err_{rms} versus \alpha and t obtained by Talbot's method for \beta = 0.7, \gamma = 0.3, m_{T} = 24, and m = 30 }.

    Problem 2

    We consider Eq (1.1) with exact solution

    \begin{equation*} u(\xi, \zeta ,t) = (t^{3}+1)(1- \xi^{2}- \zeta^{2}) \end{equation*}

    and the source term is

    \begin{split} \mathcal{G}(\xi, \zeta ,t)= &(1- \xi^{2}- \zeta^{2}) \bigg[ \varrho_{1,1} \frac{6}{\Gamma(4-\alpha)}t^{3-\alpha}+3\varrho_{2} t^{2}+ \varrho_{3,1} \frac{6}{\Gamma(4-\beta)}t^{3-\beta}\\ &+\varrho_{4} (t^{3}+1)\bigg] +4\varrho_{5} (t^{3}+1)+ \varrho_{6} \frac{6}{\Gamma(4-\gamma)}t^{3-\gamma} \end{split}

    where \varrho_{1, 1} = \varrho_{2} = \varrho_{3, 1} = \varrho_{4} = \varrho_{5} = \varrho_{6} = 1, (\xi, \zeta) \in [-1, 1]^{2} and the initial-boundary conditions are obtained from the analytical solution.

    The Err_{2}, Err_{\infty}, and Err_{rms} errors of the proposed method for problem 2 obtained via the Stehfest's and the Talbot's methods by using various values of m, \; m_{T}, \; m_{S}, \; \alpha, \; \beta, \; \gamma, and t = 1 are presented in Tables 3 and 4, respectively.

    Table 3.  The Err_{2}, Err_{\infty}, Err_{rms} computed by Stehfest's method for problem 2.
    (\alpha, \beta, \gamma) m m_{S} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    (1.3, 0.2, 0.9) 22 14 4.4852 \times 10^{-03} 1.6695 \times 10^{-03} 2.0387 \times 10^{-04} 0.188688
    25 6.8895 \times 10^{-03} 1.5726 \times 10^{-03} 2.7558 \times 10^{-04} 0.394322
    28 1.7303 \times 10^{-02} 6.2685 \times 10^{-03} 6.1797 \times 10^{-04} 0.524445
    30 10 9.4263 \times 10^{-02} 5.9793 \times 10^{-03} 3.1421 \times 10^{-03} 0.568557
    12 2.1084 \times 10^{-03} 3.2664 \times 10^{-04} 7.0279 \times 10^{-05} 0.658355
    14 2.4313 \times 10^{-02} 9.5089 \times 10^{-03} 8.1042 \times 10^{-04} 0.818657
    (1.6, 0.5, 0.7) 22 14 2.2303 \times 10^{-03} 1.0201 \times 10^{-03} 1.0138 \times 10^{-04} 0.200221
    25 6.1739 \times 10^{-03} 2.3026 \times 10^{-03} 2.4695 \times 10^{-04} 0.316503
    28 9.0251 \times 10^{-03} 2.3816 \times 10^{-03} 3.2232 \times 10^{-04} 0.510553
    30 10 9.4257 \times 10^{-02} 5.9791 \times 10^{-03} 3.1419 \times 10^{-03} 0.550572
    12 2.0884 \times 10^{-03} 2.1157 \times 10^{-04} 6.9612 \times 10^{-05} 0.663605
    14 1.3078 \times 10^{-02} 3.9249 \times 10^{-03} 4.3594 \times 10^{-04} 0.735535
    (1.9, 0.8, 0.5) 22 14 1.7605 \times 10^{-03} 4.8804 \times 10^{-04} 8.0021 \times 10^{-05} 0.181867
    25 3.2429 \times 10^{-03} 1.5817 \times 10^{-03} 1.2971 \times 10^{-04} 0.318878
    28 7.1727 \times 10^{-03} 2.6411 \times 10^{-03} 2.5617 \times 10^{-04} 0.512775
    30 10 9.4250 \times 10^{-02} 5.9791 \times 10^{-03} 3.1417 \times 10^{-03} 0.585836
    12 2.1968 \times 10^{-03} 2.7998 \times 10^{-04} 7.3227 \times 10^{-05} 0.644230
    14 1.1617 \times 10^{-02} 4.9286 \times 10^{-03} 3.8722 \times 10^{-04} 0.708296

     | Show Table
    DownLoad: CSV
    Table 4.  The Err_{2}, Err_{\infty}, Err_{rms} computed by Talbot's method for problem 2.
    (\alpha, \beta, \gamma) m m_{T} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    (1.3, 0.2, 0.9) 22 30 2.2778 \times 10^{-09} 8.4692 \times 10^{-10} 1.0354 \times 10^{-10} 1.387058
    25 5.2077 \times 10^{-09} 2.4683 \times 10^{-09} 2.0831 \times 10^{-10} 2.893077
    28 5.9240 \times 10^{-09} 1.8102 \times 10^{-09} 2.1157 \times 10^{-10} 4.540104
    30 20 6.8804 \times 10^{-06} 4.3647 \times 10^{-07} 2.2935 \times 10^{-07} 4.290310
    25 1.5504 \times 10^{-08} 1.7167 \times 10^{-09} 5.1681 \times 10^{-10} 5.406541
    30 9.2636 \times 10^{-09} 2.3673 \times 10^{-09} 3.0879 \times 10^{-10} 6.253001
    (1.6, 0.5, 0.7) 22 30 1.1655 \times 10^{-09} 4.5895 \times 10^{-10} 5.2979 \times 10^{-11} 1.413894
    25 2.4327 \times 10^{-09} 8.7358 \times 10^{-10} 9.7309 \times 10^{-11} 2.551745
    28 5.9798 \times 10^{-09} 1.5282 \times 10^{-09} 2.1357 \times 10^{-10} 4.802169
    30 20 6.8807 \times 10^{-06} 4.3647 \times 10^{-07} 2.2936 \times 10^{-07} 4.319773
    25 1.5690 \times 10^{-08} 2.3208 \times 10^{-09} 5.2301 \times 10^{-10} 5.494691
    30 9.9438 \times 10^{-09} 3.3832 \times 10^{-09} 3.3146 \times 10^{-10} 6.327990
    (1.9, 0.8, 0.5) 22 30 9.5101 \times 10^{-10} 3.0484 \times 10^{-10} 4.3228 \times 10^{-11} 1.500904
    25 1.6298 \times 10^{-09} 5.9318 \times 10^{-10} 6.5192 \times 10^{-11} 2.605735
    28 3.3339 \times 10^{-09} 1.4952 \times 10^{-10} 1.1907 \times 10^{-10} 4.379358
    30 20 6.8805 \times 10^{-06} 4.3647 \times 10^{-07} 2.2935 \times 10^{-07} 4.246011
    25 1.5398 \times 10^{-08} 1.2926 \times 10^{-09} 5.1325 \times 10^{-10} 5.247654
    30 6.2433 \times 10^{-09} 3.0606 \times 10^{-09} 2.0811 \times 10^{-10} 6.329043

     | Show Table
    DownLoad: CSV

    The numerical solution of problem 2 is shown in Figure 10. Figure 11 presents a comparison of error norms Err_{2}, \; Err_{\infty}, and Err_{rms} computed via the Stehfest's and the Talbot's methods versus m_{S} and m_{T} with \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m = 30, and t = 1 . The comparison of Err_{2}, Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods versus t with \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m_{S} = 14, \; m_{T} = 30, and m = 30 is shown in Figure 12.

    Figure 10.  Numerical solution of problem 2 obtained by proposed scheme.
    Figure 11.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus m_{S} and for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus m_{T} for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, t = 1, and m = 30 .
    Figure 12.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus t for \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m_{S} = 14, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus t for \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m_{T} = 30, and m = 30 .

    The comparison of Err_{2}, \; Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods versus \alpha with \beta = 0.7, \; \gamma = 0.3, \; m = 30, \; m_{S} = 14, \; m_{T} = 30, and t = 1 is presented in Figure 13. The comparison of Err_{2}, Err_{\infty}, and Err_{rms} computed using the Stehfest's and the Talbot's methods versus \beta with \alpha = 1.5, \gamma = 0.7, m = 30, m_{S} = 14, m_{T} = 30, and t = 1 is shown in Figure 14.

    Figure 13.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus \alpha for \beta = 0.7, \gamma = 0.3, m_{S} = 14, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus \alpha for \beta = 0.7, \gamma = 0.3, m_{T} = 30, t = 1, and m = 30 .
    Figure 14.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus \beta for \alpha = 1.5, \; \gamma = 0.3, \; m_{S} = 14, \; t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus \beta for \alpha = 1.5, \; \gamma = 0.3, \; m_{T} = 30, \; t = 1, and m = 30 .

    The comparison of Err_{2}, Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods versus \gamma with \alpha = 1.5, \beta = 0.7, m = 30, m_{S} = 14, m_{T} = 30, and t = 1 is presented in Figure 15. The plots of Err_{2} obtained using the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and \beta \in [0, 1] with \gamma = 0.3, \; t = 1, \; m_{S} = 16, m_{T} = 30, and m = 30 are shown in Figure 16.

    Figure 15.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Stehfest's method for \alpha = 1.5, \beta = 0.7, m_{S} = 14, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Talbot's method for \alpha = 1.5, \beta = 0.7, m_{T} = 30, t = 1, and m = 30 .
    Figure 16.  (a) Plot of Err_{2} versus \alpha and \beta obtained by Stehfest's method for \gamma = 0.5, \; t = 1, \; m_{S} = 16, and m = 30. (b) Plot of Err_{2} versus \alpha and \beta obtained by Talbot's method for \gamma = 0.5, \; t = 1, \; m_{T} = 24, and m = 30 .

    In Figure 17, the plots of Err_{\infty} obtained using the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and \gamma \in [0, 1] with \beta = 0.7, \; t = 1, \; m_{S} = 16, \; m_{T} = 30, and m = 30 are shown. Similarly, Figure 18 shows the plot of Err_{rms} obtained via the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and t \in [0, 1] with \beta = 0.7, \gamma = 0.3, m_{S} = 16, \; m_{T} = 30, and m = 30. Comparable performances like the one we witnessed in Problem 2 are observed.

    Figure 17.  (a) Plot of Err_{\infty} versus \alpha and \gamma obtained by Stehfest's method for \beta = 0.8, t = 1, m_{S} = 16, and m = 30. (b) Plot of Err_{\infty} versus \alpha and \gamma obtained by Talbot's method for \beta = 0.8, t = 1, m_{T} = 24, and m = 30 .
    Figure 18.  (a) Plot of Err_{rms} versus \alpha and t obtained by Stehfest's method for \beta = 0.8, \gamma = 0.5, m_{S} = 16, and m = 30. (b) Plot of Err_{rms} versus \alpha and t obtained by Talbot's method for \beta = 0.8, \gamma = 0.5, m_{T} = 24, and m = 30 .

    Problem 3

    We consider Eq (1.1) with analytical solution as

    \begin{equation*} u(\xi, \zeta ,t) = t^{2} \sin(1- \xi) (e^{\xi}-1) \sin(1- \zeta)(e^{\zeta}-1) \end{equation*}

    and the source term is

    \begin{equation*} \begin{split} \mathcal{G}(\xi, \zeta ,t) & = \sin(1- \xi) (e^{\xi}-1) \sin(1- \zeta)(e^{\zeta}-1) \bigg[ \varrho_{1,1} \frac{2}{\Gamma(3-\alpha)}t^{2-\alpha}+2\varrho_{2} t+ \varrho_{3,1} \frac{2}{\Gamma(3-\beta)}t^{2-\beta}+\varrho_{4} t^{2}\bigg]\\ & -\bigg(\varrho_{5} t^{2}+ \varrho_{6} \frac{2}{\Gamma(3-\gamma)}t^{2-\gamma} \bigg) \bigg[ sin(\zeta-1) (e^{\zeta}-1)(2cos(\xi-1)e^{\xi}+sin(\xi-1))\\ & +sin(\xi-1) (e^{\xi}-1)(2cos(\zeta-1)e^{\zeta}+sin(\zeta-1))\bigg] \end{split} \end{equation*}

    where \varrho_{1, 1} = \varrho_{2} = \varrho_{3, 1} = \varrho_{4} = \varrho_{5} = \varrho_{6} = 1 and \mathit{\boldsymbol{\bar{\xi}}} = (\xi, \zeta) \in [-1, 1]^{2}. The initial-boundary conditions are obtained from the analytical solution. The Err_{2}, Err_{\infty}, and Err_{rms} errors obtained using the proposed method for problem 3 obtained via the Stehfest's and the Talbot's methods by using various values of m, \; m_{T}, \; m_{S}, \; \alpha, \; \beta, \; \gamma, and t = 1 are presented in Tables 5 and 6, respectively.

    Table 5.  The Err_{2}, Err_{\infty}, Err_{rms} computed by Stehfest's method for problem 3.
    (\alpha, \beta, \gamma) m m_{S} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    (1.3, 0.2, 0.9) 22 14 3.3480 \times 10^{-05} 8.6337 \times 10^{-06} 1.5218 \times 10^{-06} 0.254196
    25 9.7449 \times 10^{-05} 3.7684 \times 10^{-05} 3.8980 \times 10^{-06} 0.363822
    28 2.3165 \times 10^{-04} 1.0850 \times 10^{-04} 8.2733 \times 10^{-06} 0.585754
    30 10 2.2116 \times 10^{-04} 1.8488 \times 10^{-05} 7.3719 \times 10^{-06} 0.542104
    12 1.7804 \times 10^{-04} 1.4958 \times 10^{-05} 5.9348 \times 10^{-06} 0.687321
    14 3.1175 \times 10^{-04} 7.7178 \times 10^{-05} 1.0392 \times 10^{-06} 0.802970
    (1.6, 0.5, 0.7) 22 14 3.0799 \times 10^{-05} 1.0421 \times 10^{-05} 1.4000 \times 10^{-06} 0.191420
    25 8.6136 \times 10^{-05} 3.5275 \times 10^{-05} 3.4454 \times 10^{-06} 0.339131
    28 1.4053 \times 10^{-04} 4.7890 \times 10^{-05} 5.0190 \times 10^{-06} 0.559366
    30 10 2.2101 \times 10^{-04} 1.8488 \times 10^{-05} 7.3670 \times 10^{-06} 0.583434
    12 1.8068 \times 10^{-04} 1.4953 \times 10^{-05} 6.0227 \times 10^{-06} 0.673867
    14 1.2604 \times 10^{-04} 4.4833 \times 10^{-05} 4.2014 \times 10^{-06} 0.733930
    (1.9, 0.8, 0.5) 22 14 2.0609 \times 10^{-05} 5.6164 \times 10^{-06} 9.3675 \times 10^{-07} 0.185813
    25 5.6168 \times 10^{-05} 1.6845 \times 10^{-05} 2.2467 \times 10^{-06} 0.337030
    28 4.7054 \times 10^{-05} 2.1159 \times 10^{-05} 1.6805 \times 10^{-06} 0.513103
    30 10 2.2103 \times 10^{-04} 1.8488 \times 10^{-05} 7.3676 \times 10^{-06} 0.648050
    12 1.7979 \times 10^{-04} 1.4953 \times 10^{-05} 5.9930 \times 10^{-06} 0.627418
    14 1.3026 \times 10^{-04} 4.7709 \times 10^{-05} 4.3420 \times 10^{-06} 0.718273

     | Show Table
    DownLoad: CSV
    Table 6.  The Err_{2}, Err_{\infty}, Err_{rms} computed by Talbot's method for problem 3.
    (\alpha, \beta, \gamma) m m_{T} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    (1.3, 0.2, 0.9) 22 30 2.7733 \times 10^{-11} 1.2506 \times 10^{-11} 1.2606 \times 10^{-12} 1.406847
    25 6.7413 \times 10^{-11} 2.5993 \times 10^{-11} 2.6965 \times 10^{-12} 2.649332
    28 1.5151 \times 10^{-11} 4.4017 \times 10^{-11} 5.4112 \times 10^{-12} 4.415065
    30 20 6.8451 \times 10^{-08} 5.7242 \times 10^{-09} 2.2817 \times 10^{-09} 4.279803
    25 1.6659 \times 10^{-10} 3.2206 \times 10^{-11} 5.5531 \times 10^{-12} 5.292333
    30 1.5071 \times 10^{-10} 4.2789 \times 10^{-11} 5.0236 \times 10^{-12} 6.272503
    (1.6, 0.5, 0.7) 22 30 2.4943 \times 10^{-11} 1.2897 \times 10^{-11} 1.1338 \times 10^{-12} 1.503324
    25 5.4494 \times 10^{-11} 1.8470 \times 10^{-11} 2.1798 \times 10^{-12} 2.687497
    28 1.1024 \times 10^{-11} 3.6371 \times 10^{-11} 3.9370 \times 10^{-12} 4.498761
    30 20 6.8453 \times 10^{-08} 5.7242 \times 10^{-09} 2.2818 \times 10^{-09} 4.177472
    25 1.4326 \times 10^{-10} 3.6404 \times 10^{-11} 4.7752 \times 10^{-12} 5.148521
    30 1.3951 \times 10^{-10} 6.3866 \times 10^{-11} 4.6503 \times 10^{-12} 6.443956
    (1.9, 0.8, 0.5) 22 30 7.6747 \times 10^{-12} 2.4080 \times 10^{-12} 3.4885 \times 10^{-13} 1.476408
    25 3.0024 \times 10^{-11} 1.6923 \times 10^{-11} 1.2010 \times 10^{-12} 2.663378
    28 4.6899 \times 10^{-11} 1.8591 \times 10^{-11} 1.6750 \times 10^{-12} 4.540820
    30 20 6.8448 \times 10^{-08} 5.7242 \times 10^{-09} 2.2816 \times 10^{-09} 4.224953
    25 1.3981 \times 10^{-10} 3.3329 \times 10^{-11} 4.6602 \times 10^{-12} 5.341039
    30 5.1002 \times 10^{-11} 2.4527 \times 10^{-11} 1.7001 \times 10^{-12} 6.258475

     | Show Table
    DownLoad: CSV

    The numerical method solution to problem 3 is shown in Figure 19. Figure 20 represents a comparison of error norms Err_{2}, Err_{\infty}, and Err_{rms} computed by the Stehfest's and the Talbot's methods versus m_{S} and m_{T} with \alpha = 1.5, \beta = 0.7, \; \gamma = 0.3, \; m = 30 , and t = 1 . The comparison of Err_{2}, Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods versus t with \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m_{S} = 14, \; m_{T} = 30, and m = 30 is shown in Figure 21.

    Figure 19.  Numerical solution of problem 3 obtained by proposed scheme.
    Figure 20.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus m_{S} and for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus m_{T} for \alpha = 1.5, \beta = 0.7, \gamma = 0.3, t = 1, and m = 30 .
    Figure 21.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus t for \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m_{S} = 14, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus t for \alpha = 1.5, \; \beta = 0.7, \; \gamma = 0.3, \; m_{T} = 30, and m = 30 .

    The comparison of Err_{2}, \; Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods for different values of \alpha with \beta = 0.7, \; \gamma = 0.3, \; m = 30, m_{S} = 14, \; m_{T} = 30, and t = 1 is shown in Figure 22. The comparison of Err_{2}, Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods for different values of \beta with \alpha = 1.5, \gamma = 0.7, m = 30, m_{S} = 14, m_{T} = 30, and t = 1 is presented in Figure 23.

    Figure 22.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus \alpha for \beta = 0.7, \gamma = 0.3, m_{S} = 14, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus \alpha for \beta = 0.7, \gamma = 0.3, m_{T} = 30, t = 1, and m = 30 .
    Figure 23.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus \beta for \alpha = 1.5, \; \gamma = 0.3, \; m_{S} = 14, \; t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus \beta for \alpha = 1.5, \; \gamma = 0.3, \; m_{T} = 30, \; t = 1, and m = 30 .

    The comparison of Err_{2}, Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods versus \gamma with \alpha = 1.5, \beta = 0.7, m = 30, m_{S} = 14, m_{T} = 30, and t = 1 is presented in Figure 24. The plot of Err_{2} using the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and \beta \in [0, 1] with \gamma = 0.3, \; t = 1, \; m_{S} = 16, m_{T} = 30, and m = 30 is shown in Figure 25.

    Figure 24.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Stehfest's method for \alpha = 1.5, \beta = 0.7, m_{S} = 14, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Talbot's method for \alpha = 1.5, \beta = 0.7, m_{T} = 30, t = 1, and m = 30 .
    Figure 25.  (a) Plot of Err_{2} versus \alpha and \beta obtained by Stehfest's method for \gamma = 0.5, \; t = 1, \; m_{S} = 16, and m = 30. (b) Plot of Err_{2} versus \alpha and \beta obtained by Talbot's method for \gamma = 0.5, \; t = 1, \; m_{T} = 24, and m = 30 .

    In Figure 26, the plots of Err_{\infty} obtained using the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and \gamma \in [0, 1] with \beta = 0.7, \; t = 1, \; m_{S} = 16, \; m_{T} = 30, and m = 30 are shown. Similarly, Figure 27 shows the plots of Err_{rms} using the Stehfest's and the Talbot's methods for \alpha \in [1, 2] and t \in [0, 1] with \beta = 0.7, \gamma = 0.3, m_{S} = 16, \; m_{T} = 30, and m = 30.

    Figure 26.  (a) Plot of Err_{\infty} versus \alpha and \gamma obtained by Stehfest's method for \beta = 0.8, t = 1, m_{S} = 16, and m = 30. (b) Plot of Err_{\infty} versus \alpha and \gamma obtained by Talbot's method for \beta = 0.8, t = 1, m_{T} = 24, and m = 30 .
    Figure 27.  (a) Plot of Err_{rms} versus \alpha and t obtained by Stehfest's method for \beta = 0.8, \; \gamma = 0.5, \; m_{S} = 16, and m = 30. (b) Plot of Err_{rms} versus \alpha and t obtained by Talbot's method for \beta = 0.8, \; \gamma = 0.5, \; m_{T} = 24, and m = 30 .

    It is observed that the accuracy of the Stehfest's method steadily decreases for m_{S}\geq 14 . From the results presented in tables and figures, we conclude that the proposed method is accurate, stable, and efficient. The numerical results demonstrates that the Talbot's method is more accurate than the Stehfest's method.

    Problem 4

    We consider a limiting case of Eq (1.1) with analytical solution as

    \begin{equation*} u(\xi, \zeta ,t) = E_{\gamma}\bigg(- \frac{\pi^{2}}{2} t^{\gamma}\bigg) \cos \big( \frac{\pi}{2} \xi\big) \cos \big( \frac{\pi}{2} \zeta \big) \end{equation*}

    and the source term is

    \begin{equation*} \begin{split} \mathcal{G}(\xi, \zeta ,t) = E_{\gamma}\bigg(- \frac{\pi^{2}}{2} t^{\gamma}\bigg) \cos \big( \frac{\pi}{2} \xi\big) \cos \big( \frac{\pi}{2} \zeta \big) \end{split} \end{equation*}

    where \varrho_{1, 1} = \varrho_{3, 1} = \varrho_{2} = 0, \varrho_{6} = \varrho_{4} = \varrho_{5} = 1, (\xi, \zeta) \in [-1, 1]^{2}, and the initial-boundary conditions are obtained from the analytical solution. E_{\gamma}(t) = \sum_{n = 0}^{\infty} \frac{t^{n}}{\Gamma(n \gamma+1)}, \; (\gamma \in \mathbb{C}, Re({\gamma}) > 0) is the Mittag Leffler function. The Err_{2}, Err_{\infty}, and Err_{rms} errors of the proposed method for problem 4 by using various values of m, \; m_{T}, \; m_{S}, \; \gamma, and t = 1 are presented in Tables 7 and 8, respectively. The numerical solution of problem 4 is shown in Figure 28. Figure 29 presents a comparison of error norms Err_{2}, Err_{\infty}, and Err_{rms} obtained via the Stehfest's and the Talbot's methods versus m_{S} and m_{T} with \gamma = 0.7, \; m = 30, and t = 1 . The comparison of Err_{2}, Err_{\infty}, and Err_{rms} obtained using the Stehfest's and the Talbot's methods versus t with \gamma = 0.7, \; m_{S} = 14, \; m_{T} = 25, and m = 30 is shown in Figure 30. The comparison of Err_{2}, Err_{\infty}, and Err_{rms} using the Stehfest's and the Talbot's methods versus \gamma with m = 30, m_{S} = 14, m_{T} = 25, and t = 1 is presented in Figure 31. The plot of Err_{2} using the Stehfest's and the Talbot's methods for \gamma \in [0.5, 1] and t \in [0, 1] with m_{S} = 14, m_{T} = 25, and m = 30 is shown in Figure 32. In Figure 33, the plots of Err_{\infty} computed using the Stehfest's and the Talbot's methods with \gamma \in [0.5, 1], \; t \in [0, 1], \; m_{S} = 14, \; m_{T} = 25, and m = 30 are shown. Similarly, Figure 34 shows the plot of Err_{rms} obtained using the Stehfest's and the Talbot's methods with \gamma \in [0.5, 1], \; t \in [0, 1], \; m_{S} = 14, \; m_{T} = 25, and m = 30. The computational results undeniably demonstrate that Talbot's method is greater in precision than Stehfest's method.

    Table 7.  The Err_{2}, Err_{\infty}, Err_{rms} computed by Stehfest's method for problem 4.
    m m_{S} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    \gamma=0.5 21 14 7.5454 \times 10^{-05} 1.0208 \times 10^{-05} 3.5930 \times 10^{-06} 0.155965
    24 8.7095 \times 10^{-05} 1.0340 \times 10^{-05} 3.6290 \times 10^{-06} 0.263790
    27 9.6493 \times 10^{-05} 1.0212 \times 10^{-05} 3.5738 \times 10^{-06} 0.475612
    30 08 3.3892 \times 10^{-05} 3.2475 \times 10^{-06} 1.1297 \times 10^{-06} 0.433870
    10 1.1589 \times 10^{-04} 1.1104 \times 10^{-05} 3.8630 \times 10^{-06} 0.562005
    12 1.1026 \times 10^{-04} 1.0576 \times 10^{-05} 3.6754 \times 10^{-06} 0.658822
    \gamma=0.7 21 14 8.2377 \times 10^{-06} 1.7617 \times 10^{-06} 3.9227 \times 10^{-07} 0.156614
    24 7.5928 \times 10^{-06} 1.5342 \times 10^{-06} 3.1637 \times 10^{-07} 0.263522
    27 1.4674 \times 10^{-05} 3.3973 \times 10^{-06} 5.4349 \times 10^{-07} 0.459460
    30 08 2.2064 \times 10^{-04} 2.1141 \times 10^{-05} 7.3547 \times 10^{-06} 0.428611
    10 2.3103 \times 10^{-04} 2.2138 \times 10^{-05} 7.7011 \times 10^{-06} 0.567407
    12 3.1044 \times 10^{-05} 2.9706 \times 10^{-06} 1.0348 \times 10^{-06} 0.623242
    \gamma=0.9 21 14 9.4376 \times 10^{-05} 1.2762 \times 10^{-05} 4.4941 \times 10^{-06} 0.160795
    24 1.0828 \times 10^{-04} 1.2942 \times 10^{-05} 4.5116 \times 10^{-06} 0.288674
    27 1.2360 \times 10^{-04} 1.2863 \times 10^{-05} 4.5779 \times 10^{-06} 0.424751
    30 08 1.2649 \times 10^{-03} 1.2120 \times 10^{-04} 4.2162 \times 10^{-05} 0.462303
    10 1.8471 \times 10^{-03} 1.7699 \times 10^{-04} 6.1570 \times 10^{-05} 0.521154
    12 6.5184 \times 10^{-04} 6.2457 \times 10^{-05} 2.1728 \times 10^{-05} 0.645028

     | Show Table
    DownLoad: CSV
    Table 8.  The Err_{2}, Err_{\infty}, Err_{rms} computed by Talbot's method for problem 4.
    m m_{T} Err_{2} Err_{\infty} Err_{rms} CPU (sec.)
    \gamma=0.5 21 30 7.5633 \times 10^{-05} 1.0211 \times 10^{-05} 3.6016 \times 10^{-06} 0.565917
    25 9.0040 \times 10^{-05} 1.0253 \times 10^{-05} 3.6016 \times 10^{-06} 1.337206
    29 1.0445 \times 10^{-04} 1.0278 \times 10^{-05} 3.6016 \times 10^{-06} 2.600535
    30 10 1.0893 \times 10^{-04} 1.0438 \times 10^{-05} 3.6311 \times 10^{-06} 1.103470
    12 1.0800 \times 10^{-04} 1.0348 \times 10^{-05} 3.6000 \times 10^{-06} 1.386152
    14 1.0804 \times 10^{-04} 1.0353 \times 10^{-05} 3.6015 \times 10^{-06} 1.538342
    \gamma=0.7 21 30 3.2260 \times 10^{-11} 4.3311 \times 10^{-12} 1.5362 \times 10^{-12} 0.598598
    25 3.8574 \times 10^{-11} 4.3535 \times 10^{-12} 1.5430 \times 10^{-12} 1.348554
    29 4.3999 \times 10^{-11} 4.2745 \times 10^{-12} 1.5172 \times 10^{-12} 2.666525
    30 12 7.2355 \times 10^{-08} 6.9329 \times 10^{-09} 2.4118 \times 10^{-09} 1.266167
    16 2.8863 \times 10^{-10} 2.7649 \times 10^{-11} 9.6211 \times 10^{-12} 1.678309
    20 4.5249 \times 10^{-11} 4.3492 \times 10^{-12} 1.5083 \times 10^{-12} 2.078595
    \gamma=0.9 21 30 5.4576 \times 10^{-12} 1.7432 \times 10^{-12} 2.5989 \times 10^{-13} 0.579733
    25 6.9966 \times 10^{-12} 1.5988 \times 10^{-12} 2.7987 \times 10^{-13} 1.395903
    29 1.3880 \times 10^{-11} 4.5894 \times 10^{-12} 4.7862 \times 10^{-13} 2.534978
    30 12 8.2775 \times 10^{-08} 7.9314 \times 10^{-09} 2.7592 \times 10^{-09} 1.266089
    16 3.3689 \times 10^{-10} 3.2289 \times 10^{-11} 1.1230 \times 10^{-11} 1.622688
    20 3.5520 \times 10^{-12} 8.5037 \times 10^{-13} 1.1840 \times 10^{-13} 1.994124

     | Show Table
    DownLoad: CSV
    Figure 28.  Numerical solution of problem 4 obtained by proposed scheme.
    Figure 29.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus m_{S} and for \gamma = 0.7, t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus m_{T} for \gamma = 0.7, t = 1, and m = 30. It can be seen that the observed errors are in good agreement with theoretical error.
    Figure 30.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Stehfest's method versus t for \gamma = 0.7, m_{S} = 14, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} obtained by Talbot's method versus t for \gamma = 0.7, m_{T} = 25, and m = 30 .
    Figure 31.  (a) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Stehfest's method for m_{S} = 14, \; t = 1, and m = 30. (b) Plots of Err_{2}, Err_{\infty}, and Err_{rms} versus \gamma obtained by Talbot's method for m_{T} = 25, \; t = 1, and m = 30 .
    Figure 32.  (a) Plot of Err_{2} versus \gamma and t obtained by Stehfest's method for m_{S} = 14, and m = 30. (b) Plot of Err_{2} versus \gamma and t obtained by Talbot's method for m_{T} = 25, and m = 30 .
    Figure 33.  (a) Plot of Err_{\infty} versus \gamma and t obtained by Stehfest's method for m_{S} = 14, and m = 30. (b) Plot of Err_{\infty} versus \gamma and t obtained by Talbot's method for m_{T} = 25, and m = 30 .
    Figure 34.  (a) Plot of Err_{rms} versus \gamma and t obtained by Stehfest's method for m_{S} = 14, and m = 30. (b) Plot of Err_{rms} versus \gamma and t obtained by Talbot's method for m_{T} = 25, and m = 30 .

    In the current work, we have developed an efficient and stable numerical method, which combines the numerical Laplace transform method in time with the pseudospectral method in space for the numerical solution of the two-dimensional time fractional multi-term mixed sub-diffusion and diffusion-wave equation. The proposed method offers an excellent approach to the solution process of the considered equation. In our technique, first, the Laplace transform was employed, which reduced the problem into an elliptic problem in the Laplace transform domain, then the pseudospectral method was utilized to obtain the approximate solution to the transformed problem. Finally, the improved Talbot's method and the Stehfest's method were used to convert the obtained solution in Laplace transform domain back into time domain. The improved Talbot's method and the Stehfest's method provide a numerical inversion process that is accurate, stable, easy to implement, and does not suffer from stability issues which occur with finite difference time stepping methods. From the results presented in tables and figures, it can be observed that the Stehfest's method provides optimal results for m_{S} < 14, and for m_{S} \geq 14, its accuracy decreases. Furthermore, the obtained numerical results clearly demonstrate that the Talbot's method is considerably more precise than the Stehfest's method.

    The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

    The authors D. Santina and N. Mlaiki would like to thank the Prince Sultan University for paying the publication fees for this work through TAS LAB. This research is supported by Prince Sultan University, Saudi Arabia.

    The authors declare that they have no conflicts of interest.



    [1] M. Aoki, H. Yoshikawa, Reconstructing Macroeconomics: A Perspective from Statistical Physics and Combinatorial Stochastic Processes, Cambridge University Press, 2006.
    [2] L. Bachelier, Théorie de la spéculation, Ann. Sci. École Norm. Sup., 17 (1900), 21–86.
    [3] W. A. Brock, C. H. Hommes, Heterogeneous beliefs and routes to chaos in a simple asset pricing model, Journal of Economic Dynamics and Control, 22 (1998), 1235–1274. doi: 10.1016/S0165-1889(98)00011-6
    [4] J. P. Bouchaud, M. Potters, Theory of financial risks. From Statistical Physics to Risk Management, Cambridge University Press, 2000.
    [5] F. Black, M. Scholes, The pricing of options and corporate liabilities, J. Political Economy, 81 (1973), 637–654. doi: 10.1086/260062
    [6] G. Caginalp, D. Balevonich, Asset flow and momentum: Deterministic and stochastic equations, Phil. Trans. Royal Soc. A., 357 (1999), 2119–2113. doi: 10.1098/rsta.1999.0421
    [7] C. Caginalp, G. Caginalp, The quotient of normal random variables and application to asset price fat tails, Physica A, 499 (2018), 457–471. doi: 10.1016/j.physa.2018.02.077
    [8] C. Caginalp, G. Caginalp, Price equations with symmetric supply/demand; implications for fat tails, Economics Letters, 176 (2019), 79–82. doi: 10.1016/j.econlet.2018.12.037
    [9] G. Caginalp, G. B. Ermentrout, A kinetic thermodynamics approach to the psychology of fluctuations in financial markets, Appl. Math. Lett., 3 (1990), 17–19.
    [10] C. Carroll, The method of endogenous gridpoints for solving dynamic stochastic optimization problems, Economics Letters, 91 (2006), 312–320. doi: 10.1016/j.econlet.2005.09.013
    [11] N. Champagnat, M. Deaconu, A. Lejay, N. Navet, S. Boukherouaa, An empirical analysis of heavy-tails behavior of financial data: The case for power laws, HAL archives-ouvertes, 2013.
    [12] R. Cont, P. Tankov, Financial modeling with jump processes, Chapman & Hall/CRC, London, 2004.
    [13] J. Daníelsson, B. N. Jorgensen, G. Samorodnitsky, M. Sarma, C. G. de Vries, Fat tails, VaR and subadditivity, J. Econometrics, 172 (2013), 283–291. doi: 10.1016/j.jeconom.2012.08.011
    [14] E. Díaz-Francés, F. J. Rubio, On the existence of a normal approximation to the distribution of the ratio of two independent normal random variables, Stat. Papers, 54 (2013), 309–323. doi: 10.1007/s00362-012-0429-2
    [15] M. DeSantis, D. Swigon, Slow-fast analysis of a multi-group asset flow model with implications for the dynamics of wealth, PloS one, 13 (2018), e0207764. doi: 10.1371/journal.pone.0207764
    [16] C. Eom, T. Kaizoji, S. H. Kang, L. Pichl, Bitcoin and investor sentiment: statistical characteristics and predictability, Physica A, 514 (2019), 511–521. doi: 10.1016/j.physa.2018.09.063
    [17] C. Eom, T. Kaizoji, E. Scalas, Fat tails in financial return distributions revisited: Evidence from the Korean stock market, Physica A, 526 (2019), 121055. doi: 10.1016/j.physa.2019.121055
    [18] A. Erdelyi, Asymptotic Expansions, Dover Publications, 1956.
    [19] L. C. Evans, Partial Differential Equations, 2nd Ed., Am. Math. Soc., 2010.
    [20] G. W. Evans, S. Honkapohja, Stochastic gradient learning in the cobweb model, Economics Letters, 61 (1998), 333–337. doi: 10.1016/S0165-1765(98)00176-1
    [21] E. Fama, The Behavior of Stock-Market Prices, J. Business, 38 (1965), 34–105. doi: 10.1086/294743
    [22] J. P. Fouque, G. Papanicolaou, R. Sircar, K. Sølna, Multiscale stochastic volatility for equity, interest rate, and credit derivatives, Cambridge University Press, 2011.
    [23] X. Gabaix, Power laws in economics: An introduction. Journal of Economic Perspectives, 30 (2016), 185–206.
    [24] X. Gabaix, P. Gopikrishnan, V. Plerou, H. E. Stanley, Institutional investors and stock market volatility, The Quarterly Journal of Economics, 121 (2006), 461–504. doi: 10.1162/qjec.2006.121.2.461
    [25] S. Gjerstad, J. Dickhaut, Price formation in double auctions, Games and economic behavior, 22 (1998), 1–29. doi: 10.1006/game.1997.0576
    [26] S. Gjerstad, The competitive market paradox, Journal of Economic Dynamics and Control, 31 (2007), 1753–1780. doi: 10.1016/j.jedc.2006.07.001
    [27] S. Gjerstad, Price dynamics in an exchange economy, Economic Theory, 52 (2013), 461–500. doi: 10.1007/s00199-011-0651-5
    [28] R. Geary, The frequency distribution of the quotient of two normal variates, J. Royal Stat. Soc., 93 (1930), 442–446. doi: 10.2307/2342070
    [29] P. Gopikrishnan, M. Meyer, L. Amaral, H. E. Stanley, Inverse cubic law for the distribution of stock price variations, Eur. Phys. J. B, 3 (1998), 139–140.
    [30] J. Henderson, R. Quandt, Microeconomic theory - A mathematical approach, McGraw-Hill, New York, 1980.
    [31] J. Hirshleifer, A. Glazer, D. Hirshleifer, Price theory and applications: decisions, markets, and information, Cambridge Univ. Press, 2005.
    [32] J. C. Hull, Options, Futures and Other Derivatives, 9th Ed., Pearson, 2015.
    [33] T. D. Jeitschko, Equilibrium price paths in sequential auctions with stochastic supply, Economics Letters, 64 (1999), 67–72. doi: 10.1016/S0165-1765(99)00066-X
    [34] D. W. Jansen, C. G. De Vries, On the frequency of large stock returns: Putting booms and busts into perspective, The review of economics and statistics, (1991), 18–24.
    [35] K. Kawamura, The structure of bivariate Poisson distribution, Kodai Math Sem. Rep., 25 (1973), 246–256. doi: 10.2996/kmj/1138846776
    [36] M. Kemp, Extreme events - robust portfolio construction in the presence of fat tails, Wiley Finance, Hoboken, NJ, 2011.
    [37] P. J. Klenow, B. A. Malin, Microeconomic evidence on price-setting, Handbook of monetary economics, 3 (2010), 231–284. doi: 10.1016/B978-0-444-53238-1.00006-5
    [38] M. Kirchler, J. Huber, Fat tails and volatility clustering in experimental asset markets, Journal of Economic Dynamics and Control, 31 (2007), 1844–1874. doi: 10.1016/j.jedc.2007.01.009
    [39] T. Lux, The stable Paretian hypothesis and the frequency of large returns: an examination of major German stocks, Applied financial economics, 6 (1996), 463–475. doi: 10.1080/096031096333917
    [40] B. LeBaron, R. Samanta, Extreme value theory and fat tails in equity markets, 2005. Available at SSRN 873656.
    [41] B. Mandelbrot, Sur certain prix speculatifs: faits empirique et modele base sur les processes stables additifs de Paul Levy, Comptes Rendus de l'Academie de Sciences, 254 (1962), 3968–3970.
    [42] A. M. Mood, F. A. Graybill, D. C. Boes, Introduction to the theory of statistics, 3rd Ed., McGraw-Hill, 1973.
    [43] B. Mandelbrot, R. Hudson, The misbehavior of markets: A fractal view of financial turbulence, Basic Books, New York, 2007.
    [44] G. Marsaglia, Ratios of normal variables, J. Stat. Softw., 16 (2006), 1–10.
    [45] H. Merdan, M. Alisen, A mathematical model for asset pricing, App. Math. Comp., 218 (2011), 1449–1456. doi: 10.1016/j.amc.2011.06.028
    [46] P. Milgrom, Discovering Prices: Auction Design in Markets with Complex Constraints, Columbia University Press, 2017.
    [47] R. N. Mantegna, H. E. Stanley, Introduction to econophysics: correlations and complexity in finance, Cambridge University press, 1999.
    [48] T. Pham-Gia, N. Turkkan, E. Marchand, Density of the ratio of two normal random variables and applications, Communications in Statistics-Theory and Methods, 35 (2006), 1569–1591. doi: 10.1080/03610920600683689
    [49] C. Plott, K. Pogorelskiy, Call market experiments: efficiency and price discovery through multiple calls and emergent Newton adjustments, American Economic Journal: Microeconomics, 9 (2017), 1–41.
    [50] D. Porter, S. Rassenti, A. Roopnarine, V. Smith, Combinatorial auction design, Proceedings of the National Academy of Sciences, 100 (2003), 11153–11157. doi: 10.1073/pnas.1633736100
    [51] S. T. Rachev, Handbook of Heavy Tailed Distributions in Finance, Handbooks in Finance, Book 1, Elsevier, 2003.
    [52] M. Raberto, S. Cincotti, S. M. Focardi, M. Marchesi, Agent-based simulation of a financial market, Physica A, 299 (2001), 319–327. doi: 10.1016/S0378-4371(01)00312-0
    [53] S. T. Rachev, C. Menn, F. J. Fabozzi, Fat-tailed and skewed asset return distributions: implications for risk management, portfolio selection, and option pricing, John Wiley & Sons, 2005.
    [54] Z. Rachev, B. Rocheva-Iotovo, S. Stoyanov, Capturing fat tails, Risk, 23 (2010), 72.
    [55] A. W. Rathgeber, J. Stadler, S. Stöckl, Financial modelling applying multivariate Lévy processes: New insights into estimation and simulation, Physica A, 532 (2019), 121386. doi: 10.1016/j.physa.2019.121386
    [56] Y. Si, S. Nadarajah, X. Song, On the distribution of quotient of random variables conditioned to the positive quadrant, Communications in Statistics-Theory and Methods, 49 (2020), 2514–2528. doi: 10.1080/03610926.2019.1576893
    [57] N. Taleb, G. Daniel, The problem is beyond psychology: The real world is more random than regression analyses, International Journal of Forecasting, Forthcoming, 2011.
    [58] Y. Tong, The multivariate normal distribution, Springer-Verlag, New York, 1990.
    [59] D. Watson, M. Getz, Price theory and its uses, 5th Ed., University Press of America, 1981.
    [60] E. R. Weintraub, Microfoundations, Cambridge University Press, 1979.
  • This article has been cited by:

    1. Ujala Gul, Zareen A. Khan, Salma Haque, Nabil Mlaiki, A Local Radial Basis Function Method for Numerical Approximation of Multidimensional Multi-Term Time-Fractional Mixed Wave-Diffusion and Subdiffusion Equation Arising in Fluid Mechanics, 2024, 8, 2504-3110, 639, 10.3390/fractalfract8110639
    2. Yubing Jiang, Hu Chen, Error estimate of GL‐ADI scheme for 2D multiterm nonlinear time‐fractional subdiffusion equation, 2024, 47, 0170-4214, 14588, 10.1002/mma.10290
    3. Salma Aljawi, Ahmad Aloqaily, Nabil Mlaiki, A local meshless method for the numerical solution of multi-term time-fractional generalized Oldroyd-B fluid model, 2024, 10, 24058440, e34061, 10.1016/j.heliyon.2024.e34061
    4. Jia Tan, Tian-jun Wang, Spectral collocation method coupled with domain decomposition for exterior problems of the Fisher equation, 2024, 175, 08981221, 289, 10.1016/j.camwa.2024.10.003
    5. Zeeshan Ali, Kamal Shah, Thabet Abdeljawad, Amjad ali, Analyzing wave dynamics of Burger-Poisson fractional partial differential equation, 2025, 26668181, 101153, 10.1016/j.padiff.2025.101153
    6. E. Aourir, H. Laeli Dastjerdi, Kamal Shah, Thabet Abdeljawad, Computational scheme for the numerical solution of fractional order pantograph delay-integro-differential equations via the Bernstein approach, 2025, 26668181, 101194, 10.1016/j.padiff.2025.101194
    7. Razia Begum, Sajjad Ali, Kamal Shah, Thabet Abdeljawad, Fractional modeling of a communicable through air droplets based disease in Mittag–Leffler kernel sense and its numerical simulation, 2025, 11, 2363-6203, 10.1007/s40808-025-02461-5
  • Reader Comments
  • © 2021 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(3140) PDF downloads(100) Cited by(0)

Figures and Tables

Figures(1)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog