A B-Spline Quasi Interpolation Crank–Nicolson Scheme for Solving the Coupled Burgers Equations with the Caputo–Fabrizio Derivative

In this paper, a Crank–Nicolson finite difference scheme based on cubic B-spline quasi-interpolation has been derived for the solution of the coupled Burgers equations with the Caputo–Fabrizio derivative. The first- and second-order spatial derivatives have been approximated by first and second derivatives of the cubic B-spline quasi-interpolation. The discrete scheme obtained in this way constitutes a system of algebraic equations associated with a bi-pentadiagonal matrix. We show that the proposed scheme is unconditionally stable. Numerical examples are provided to verify the efficiency of the method.


Introduction
e coupled Burgers equations are coupled partial differential equations which are capable of describing realitic polydispersive supensions. e coupled Burgers equation perdicts an interesting phenomenon, which is called phase shifts [1]. is equation is one of the fundamental models in fluid mechanics and arises in gas dynamics, chromatography, and flood waves in rivers [2]. e coupled viscous Burges equation is given by with the initial conditions and the boundary conditions u(a, t) � f 1 (a, t), where η is a real constant and α and β are arbitrary constants depending on the system parameters such as Peclet number, Stokes velocity of particles due to gravity, and Brownian diffusivity [3]. Spline is a special function defined piecewise by polynomials. e spline approximation first appeared in a paper by Schoenberg [4]. Spline interpolation is a form of interpolation where the interpolant is a special type of piecewise polynomial called a spline. Applications of spline function in fractional partial differential equations can be found in [5][6][7][8][9][10][11][12][13][14][15].
x ∈ [a, b], t ∈ [0, T], 0 < c < 1, with initial and the boundary conditions where c is order of time fractional derivative. Also, η, α 1 , and α 2 are those ones we said before.
For the description of memory and some physical properties of various materials and processes, modeling with fractional derivatives is very appropriate. is is the main benefit of fractional derivatives in comparison with classical integer order models, in which such effects are missed. In recent years, the coupled system of Burgers equations with fractional derivatives has been the focus of attention. For example, in [27], the Adomian decomposition method is directly extended to study the coupled Burgers equations with time and space fractional derivatives. Khan et al. [28] proposed the generalized differential transform method (GDTM) and homotopy perturbation method (HPM) for time fractional Burgers and coupled Burgers equations. e fractional variational iteration method (FVIM) to solve a time and space fractional coupled burgers equations is given by Prakash et al. [29]. In [30], a q-homotopy analysis transform method (q-HATM) for time and space fractional coupled Burgers equations is introduced. Aminikhah and Malekzadeh [31] introduced a new homotopy perturbation method for system of variable coefficient coupled Burgers equations with time fractional derivative. In [32], the Laplace-Adomian decomposition method (LADM), the Laplace-variational iteration method (LVIM), and the reduced differential transform method (RDTM) are proposed to solve the one-and two-dimensional fractional coupled Burgers equations. Albuohimad and Adibi derived a hybrid spectral exponential Chebyshev method (HSECM) for time fractional coupled Burgers equations [33]. In [34], authors investigate the fractional coupled viscous Burgers equation involving Mittag-Leffler kernel. In [35], the generalized twodimensional differential transform method (DTM) was applied to solve the coupled Burgers equations with space and time fractional derivatives. Ozdemir et al. used the Gegenbauer wavelets-based computational methods to find the approximate solutions of the coupled system of Burgers equations with time fractional derivative [36].
Our aim is to propose a Crank-Nicolson finite difference scheme using cubic B-spline quasi-interpolation to solve time fractional coupled viscous Burgers equations. e firstand second-order spatial derivatives have been approximated by first and second derivatives of the cubic B-spline quasi-interpolation.
is approximations have not been used for the fractional coupled Burgers equations before. e paper is organized as follows. In Section 2, we present some basic definitions and concepts of quasiinterpolants. In Section 3, using the quasi-interpolant and Crank-Nicolson finite difference method, we obtain a numerical scheme. e stability of this method is studied in Section 4. In Section 5, some numerical examples are proposed. Finally, conclusions are given in Section 6.

Univariate Spline Quasi-Interpolants
In this section, we introduce the basic concepts about B-spline and univariate B-spline quasi-interpolants that we will use in Section 3.
According to [37], let P 1 d ≔ space of univariate polynomials of degree at most d, (9) and Ω � [a, b] be an interval that has been partitioned into subintervals via a set of points We define the space of univariate polynomial splines of smoothness r and degree d with knots Δ as where 0 ≤ r < d are given integers. We have For a formal proof of this fact, see eorem 4.4 of [38]. Given 0 ≤ r < d and Δ � x i k+1 i�0 , the associated extended partition Δ e is defined to be y i n+d+1 i�0 , where n is the dimension of S r d (Δ) given in (11): 2 Mathematical Problems in Engineering a � y 1 � · · · � y d+1 , y n+1 � · · · � y n+d+1 � b, Given an extended partition Δ e , let for i � 1, . . . , n + d, and let We call these the normalized B-splines of order m (or degree m − 1) associated with the extended partition Δ e .
In [37], univariate B-spline quasi-interpolants can be defined as a formula of the form where N i n+d i�1 are the B-splines forming a basis of S r d (Δ). Quasi-interpolants have been heavily studied in the literature. Some basic ideas and sources for further information can be found in [38]. For a good approximations, we need to make sure it reproduces polynomials, i.e., Qp � p for all p ∈ P 1 d . For each i � 1, . . . , n, we assume that the coefficient λ i is a linear functional defined on C[a, b] that can be computed from samples of f at some set of points According to [39], the error of a quasi-interpolation satisfies where D d y � [y d+1 , y n+1 ], D x is the union of the supports of all B-splines N i , i ∼ x and ‖f (d+1) ‖ ∞,D x denote the maximum norm of f (d+1) on D x and h(x) � max y∈D x |y − x| that ∼ is used to indicate proportionality. If the local mesh ratio is bounded, i.e., if the quotients of the lengths of adjacent knot intervals are ≤r y , then the error of the derivatives on the knot intervals (y l , y l+1 ) can be estimated by en, (17) defines a linear operator mapping C[a, b] into S 2 3 (Δ) with Qp � p for all cubic polynomials p. For approximate derivatives of f by derivatives of Q 3 f up to the order h 3 , we can evaluate the value of f ′ and f 〞 at . . , n. By solution of the linear systems, we obtain where D 1 , D 2 ∈ R (n+1)×(n+1) and is obtained as follows:

Numerical Scheme
We consider a grid . e values of the function u at the grid points are denoted u k j � u(x j , t k ) and U k j is the approximate solution at the point (x j , t k ).
A discrete approximation to the CF 0 D c t u(x, t) at (x j , t k ) can be obtained by the following approximation [40]: where Proof. Using the Taylor series expansion with integral remainder, we have Since, and e − σΔt(k− l) < 1 for l � 1, 2, . . . , k and t k � kΔt; hence, the result will be achieved. Now, using eorem 1, we obtain We introduce some lemmas which will be used in numerical scheme and stability analysis. □ Lemma 1 (see [41]). Suppose u(t) ∈ C 3 [0, t k+1 ]; then, we have Lemma 2 (see [42]). For the definition M j , we have M j > 0 and M j+1 < Mj, ∀j ≤ k.

Mathematical Problems in Engineering
erefore, in each time step we solve the following bi-

Stability Analysis
To study the stability analysis of the proposed scheme, we use the Fourier method. In applying the Fourier stability method, the nonlinear terms are temporarily frozen, since the stability analysis is strictly only applicable to linear equations. us, we have linearized the nonlinear terms uu x and (uv) x in equation (4) by freezing u and v as a local constants β 1 and β 2 , respectively. We have Substituting approximations (22) and (32) yield the following difference equation: Let U k j and V k j be the approximate solutions of (43), and define So, we have 6 Mathematical Problems in Engineering where ] � (ηβ 1 + α 1 β 2 ). Now, we define the grid functions as follows: We expand the ζ k (x) and ξ k (x) into the following Fourier series expansions: where (48) Applying the Parseval equality, we have (51) Now, we suppose that where σ x � (2lπ/L). Substituting the above relations into (45) leads to (53)

Mathematical Problems in Engineering
Set We have (56) (55), then there are positive constants C k such that

Theorem 2. If A k is the solution of equation
Proof. We use the mathematical induction for proof. For k � 1, we have Using the convergence of the series on the right-hand side of equation (49), we know that there exists a positive constant P 2 such that So, Now, suppose that By equation (56), we have Now, assume that So, □ Remark 1. From (5), with similar way, there are positive constants Q k such that Theorem 3. e finite difference schemes (39) and (40) are unconditionally stable for c ∈ (0, 1).
Proof. According to eorem 2 and Remark 1, using (51), we obtain so that which shows that schemes (39) and (40) are unconditionally stable.

Numerical Results
In this section, we provide two examples to illustrate efficiency of schemes (39) and (40) where e n j � u n j − U n j and e n j � v n j − V n j . We evaluate the convergence order with the following formula: We solve this problem with the method developed in this article with several values of T and c. Graphs of numerical solution and exact solution at different times have been demonstrated in Figure 4. Figure 4 shows that the proposed method is efficient. Table 5 gives the approximation errors for t � 0.5, 1, 1.5, 2, 2.5, 3, 4 with different α 1 , α 2 and CPU times. We choose h � (1/160) for space step size to obtain the numerical results. Figure 5 shows the comparison of numerical solution and exact solution for c � 0.1 at t � 1. Also, Figure 6 shows the pointwise errors and contour plot of numerical solution for u(x, t).

Conclusion
In this article, we constructed a Crank-Nicolson finite difference scheme based on cubic B-spline quasi-interpolation to solve time fractional coupled Burgers equations. By the Fourier series method, we proved that this scheme is unconditionally stable. Numerical examples have been carried out to show the convergence orders and applicability of the scheme and error norms are calculated with respect to different space step sizes. From the error tables and graphs of exact and numerical solution, we can say that our method has a good accuracy. For the numerical computations, we have used Matlab.

Data Availability
All results have been obtained by conducting the numerical procedure and the ideas can be shared for the researchers.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.