Quintic B-spline method for time-fractional superdiffusion fourth-order differential equation

The main objective of this paper is to obtain the approximate solution of superdiffusion fourth-order partial differential equations. Quintic B-spline collocation method is employed for fractional differential equations (FPDEs). The developed scheme for finding the solution of the considered problem is based on finite difference method and collocation method. Caputo fractional derivative is used for time fractional derivative of order $$\alpha$$α, $$1<\alpha <2$$1<α<2. The given problem is discretized in both time and space directions. Central difference formula is used for temporal discretization. Collocation method is used for spatial discretization. The developed scheme is proved to be stable and convergent with respect to time. Approximate solutions are examined to check the precision and effectiveness of the presented method.


Introduction
The idea of differentiation and integration to non integer order has been studied in fractional calculus. The subject fractional calculus is the generalized form of classical calculus. It is as old as classical calculus but getting more attention these days due to its applications in various fields, such as science and engineering. Oldham and Spanier [1], Podlubny [3], and Miller and Ross [2] provide the history and a comprehensive treatment of this subject.
In recent decades, it has been observed by many mathematicians, scientists, and researchers that non-integer operators (differential and integral) play an important role in describing the properties of physical phenomena. Fractional derivatives and integrals efficiently describe the history and inherited properties of different processes and equipments. It has further been observed that fractional models are more efficient and accurate than already developed classical models. The mathematical modeling of real world problems, such as earthquake, traffic flow, fluid flow, signal processing, and viscoelastic problems, results in FPDEs.
A fourth-order linear PDE is defined as The fourth-order problems play a significant role for the development of engineering and modern science. For example, window glasses, floor systems, bridge slabs, etc. are modeled by fourth-order PDEs. It also represents beam problem where v is the transversal displacement of the beam, m is the ratio of flexural rigidity of the beam to its mass per unit length, t is time, and x is space variable, h(x, t) is the dynamic driving force acting per unit mass. The considered equation, is the superdiffusion fourth-order partial differential equation, which has the following form where t 0 ; t f represents initial and final time, respectively. The time-fractional superdiffusion fourth-order problems play an important role in modern science and engineering. For example, airplane wings and transverse vibrations of sustained tensile beam can be modeled as plates with different boundary supports which are successfully governed by superdiffusion fourth-order differential equations.
The given problem is solved with initial conditions Collocation method is widely used to provide the solution of partial differential equations. Depending on the situation, some times, it is useful to find the solution of FPDEs at various locations in the given problem domain, then the spline solutions guarantee to give the information of spline interpolation between mesh points. In many situations, it is harder to calculate the analytical solutions of fractional PDEs using analytical methods. Due to this reason, the mathematicians are motivated to develop numerical methods for the approximate solution of FPDEs. The approximate methods for solving FPDEs have become popular among the researchers in the last 10 years. A variety of approximate techniques were used for finding the solution of FPDEs, such as the FDM, FEM, FVM, and the collocation method etc.
Many authors applied different numerical methods to find the solution of fractional PDEs.
Siddiqi and Arshed [7] developed a numerical scheme for the solution of time-fractional fourth-order partial differential equation with fractional derivative of order a ð0\a\1Þ.
Lin and Xu [6] developed a numerical scheme for timefractional diffusion equation with fractional derivative of order a; ð0\a\1Þ. Yang et al. [11] developed a novel numerical scheme for solving fourth-order partial integrodifferential equation with a weakly singular kernel. Yang et al. [8] developed a quasi-wavelet-based numerical method for solving fourth-order partial integro-differential equations with a weakly singular kernel. Khan et al. [9], proposed a numerical method for the solution of time-fractional fourthorder differential equations with variable coefficients using Adomian decomposition method (ADM) and He's variational iteration method (HVIM). Zhang and Han [10] developed a quasi-wavelet-based numerical method for time-dependent fractional partial differential equation. Zhang et al. [12] proposed quintic B-spline collocation method for the numerical solution of fourth-order partial integro-differential equations with a weakly singular kernel. Khan et al. [4] used sextic spline solution for solving fourthorder parabolic partial differential equation. Dehghan [5] used finite difference techniques for the numerical solution of partial integro-differential equation arising from viscoelasticity. In [15], Bhrawy and Abdelkawy used fully spectral collocation approximation for solving multi-dimensional fractional Schrodinger equations. Bhrawy et al. [16] proposed a new formula for fractional integrals of Chebyshev polynomials and application for solving multiterm fractional differential equations. Bhrawy [17] developed a new spectral algorithm for a time-space fractional partial differential equations with subdiffusion and superdiffusion. In [18], Bhrawy et al. gave a review of operational matrices and spectral techniques for fractional calculus. Bhrawy et al. [19] proposed Chebyshev-Laguerre Gauss-Radau collocation scheme for solving time fractional sub-diffusion equation on a semi-infinite domain. Bhrawy [20] developed a space-time collocation scheme for modified anomalous subdiffusion and nonlinear superdiffusion equations. In [21], Doha et al. developed an efficient Legendre Spectral Tau matrix formulation for solving fractional subdiffusion and reaction sub-diffusion equations. Bhrawy [22] proposed a highly accurate collocation algorithm for 1 ? 1 and 2 ? 1 fractional percolation equations. In [23], Bhrawy and Zaky used fractional-order Jacobi tau method for a class of time-fractional PDEs with variable Coefficients. This paper is designed to determine the approximate solution of superdiffusion fourth order PDEs. The approximate solution is based on central difference formula and B-spline collocation method.
The paper is divided into six sections. The quintic B-spline basis function is given in Sect. 2. The finite difference approximation for time discretization of the given problem is discussed in Sect. 3. The stability and convergence of temporal discretization are established. Space discretization based on B-spline collocation technique is discussed in Sect. 4. Approximate solutions are obtained in Sect. 5 which support the theoretical results. The conclusion is given is Sect. 6.

Quintic B-spline
The problem domain [a, b] has been subdivided into N elements having uniform step size h with knots z j ; To define basis functions for quintic B-spline, it is first needed to introduce ten additional knots, The quintic B-spline basis functions Q j ðzÞ, j ¼ À2; . . .; N þ 2 are defined as The approximate solution V kþ1 ðzÞ at k þ 1 time level to the exact solution v kþ1 ðzÞ is obtained as a linear combination of quintic B-spline basis functions Q j ðzÞ as where c j are unknowns. Using Eq. (2.1) and basis functions, the approximate solution V and its derivatives at the nodes z j , can be calculated as

Time discretization
Caputo fractional derivative is used for temporal discretization. Let is the time step size and t f be the final time. Caputo fractional derivative is discretized using central difference approximation as The time discrete operator is represented by D a and defined as D a vðz; t kþ1 Þ : Then Eq. (3.1) can be written as The second order approximation of Caputo derivative as discussed in [13] is given as where E is a constant.
The above equation can be rewritten as where v kþ1 ðxÞ ¼ vðx; t kþ1 Þ and a 0 ¼ Cð3 À aÞDt a , Using the following initial conditions and the boundary conditions ðI; II; IIIÞ The developed method is a three time level method. To implement the developed method, it is first required to calculate ðv À1 Þ and ðv 0 Þ.
v À1 ðzÞ ¼ vðz; t 0 Þ À Dtv t ðz; t 0 Þ: v À1 ðzÞ ¼ f 0 ðzÞ À Dtf 1 ðzÞ: In particular, for k ¼ 0, the scheme simply leads to r kþ1 is the error term can be defined as ; where L 2 X ð Þ is the space of measurable functions whose square is Lebesgue integrable in X. The inner products of L 2 X ð Þ and H 2 ðXÞ are defined, respectively, by and the corresponding norms by The norm : Instead of using the above standard H 2 -norm, it is preferred to define k:k 2 where a 0 ¼ Cð2 À aÞDt a . For the convergence and stability analysis, the following weak formulation of Eqs. (3.4) and (3.7) is needed, i.e., finding v kþ1 2 H 2 0 ðXÞ, such that for all w 2 H 2 0 ðXÞ, where k:k 2 is defined in (3.10).
Proof For k ¼ 0, let w ¼ v 1 in Eq. (3.11), Equation (3.11) takes the following form Integrating the above equation, we have Due to boundary conditions on w all the boundary contributions disappeared.
Using the inequality kwk 0 kwk 2 and Schwarz inequality Eq. (3.15) takes the following form Suppose that the result hold for w ¼ v l i.e.

ð3:16Þ
Taking w ¼ v kþ1 in Eq. (3.11), it can be written as Integrating the above equation, we get Using the Schwarz inequality and the inequality kwk 0 kwk 2 the above equation can be written as where E is a constant and Proof Let s k ¼ vðz; t k Þ À v k ðzÞ, for k ¼ 1, by combining Eqs. (1.1), (3.11), the error equation becomes using Schwarz inequality and the inequality kwk 0 kwk 2 Let v ¼ s 1 , noting s 0 ¼ 0, the above equation takes the following form using kwk 0 kwk 2 and Schwarz inequality the above equation, for v ¼ s kþ1 , can be written as

Space discretization
Collocation method with quintic B-spline basis functions is employed for space discretization. Putting Eq. (2.1) into Eq. (3.4), we obtain the following equation The above equation can be rewritten as þ a 0 q kþ1 : The above equation can be rewritten as a system of linear equations þ a 0 q kþ1 : This system has N þ 1 equations in N þ 5 unknowns c kþ1 À2 ; c kþ1 À1 ; . . .; c kþ1 Nþ1 ; c kþ1 Nþ2 . Unique solution of the above system is obtained by eliminating the unknowns c kþ1 À2 ; c kþ1 À1 ; c kþ1 Nþ1 , and c kþ1 Nþ2 using boundary conditions.

Numerical results
In The error norms L 1 and L 2 are calculated for N ¼ 100 and different temporal step sizes. The temporal rate of convergence at t f ¼ 1:0 is presented in Table 1. The temporal rate of convergence obtained by the proposed method support the theoretical results.   To indicate the effects of the proposed method for larger time level K, the exact solution and the approximate solution are plotted using N ¼ 50, K ¼ 1000 and Dt ¼ 0:00001 as shown in Fig. 1. It is clear from the Fig. 1 that the numerical solution is highly consistent with the exact solution, which indicates that the proposed method is very effective. In Fig. 2, the exact solution is represented by solid line and the numerical solution is represented by dotted line at K ¼ 1000 time level.
The temporal rate of convergence of L 2 norm error as a function of the time step Dt for a ¼ 1:50 is shown in Fig. 3. The exact solution of the problem is vðz; tÞ ¼ t þ 1 ð Þ p 5 sin pz þ 1 p 5 cos pz À 1 p 5 cos 3pz : The error norms L 1 and L 2 are calculated for N ¼ 80 and different temporal step sizes. The temporal rate of convergence at t f ¼ 1:0 is presented in Table 2. The temporal rate of convergence obtained by the proposed method support theoretical results.
To indicate the effects of the proposed method for larger time level, the exact solution and the approximate solution are plotted using N ¼ 80, K ¼ 1000, and Dt ¼ 0:00001 as shown in Fig. 4. It is clear from the Fig. 4 that the The temporal rate of convergence of L 2 norm error as a function of the time step Dt for a ¼ 1:75 is shown in Fig. 6.

Conclusion
The solution of superdiffusion fourth order equation has been developed. Central difference formula and quintic Bspline are used for discretizing time and space variables, respectively. The proposed method is stable and the approximate results approach the analytical results with order OðDt 2 Þ. The approximate results obtained by the suggested method support the theoretical estimates.
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://crea tivecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.