Baojin Su¹, Yanan Wang¹, Jingwen Qi¹, Yousen Li^2*
1Qihe County Yongfeng Experimental School, Dezhou, China.
²Middle School Affiliated to Shandong University, Jinan, China.
DOI: 10.4236/jamp.2022.1010210   PDF    HTML   XML   114 Downloads   548 Views   Citations

Abstract

In this paper, we introduce high-order finite volume methods for the multi-term time fractional sub-diffusion equation. The time fractional derivatives are described in Caputo’s sense. By using some operators, we obtain the compact finite volume scheme have high order accuracy. We use a compact operator to deal with spatial direction; then we can get the compact finite volume scheme. It is proved that the finite volume scheme is unconditionally stable and convergent in L_∞-norm. The convergence order is O(τ^2-α + h⁴). Finally, two numerical examples are given to confirm the theoretical results. Some tables listed also can explain the stability and convergence of the scheme.

Keywords

Multi-Term Time Fractional Sub-Diffusion Equation, High-Order Compact Finite Volume Scheme, Stable, Convergent

Share and Cite:

1. Introduction

During the past several decades, the study of fractional partial differential equations has attracted many scholars’ attention [1] [2] [3] [4]. Fractional order partial differential equations can provide mathematical tools to describe many phenomena, such as engineering [5] [6], chemistry, physics and so on. Fractional order partial differential equations are different from classical partial differential equations. The fractional integrals and derivatives satisfy the nonlocal properties. Although many important works about theoretical analysis have been carried on, the analytic solutions cannot be obtained exactly in most fractional partial differential equations. Even if their solutions can be found, they are in the form of series, which are difficult to evaluate. We can refer to some recent related works [7] [8] [9]. So the numerical investigation of the fractional partial differential equations based on convergence and stability analyses has been an important topic in recent years [10] - [22]. Liu et al. have discussed the stability and convergence of fractional partial differential equations by using a new energy method [10] [11] [12] [13] [14]. Yuste and Acedo found a finite difference method that can solve the time fractional diffusion equation by using the forward-Euler scheme, and they discussed the stability and convergence of the scheme [15] [16].

For the time-fractional partial differential equations, there also have been lots of works. Henry and Langlands presented an implicit finite difference scheme by using the L1 scheme to approximate the time-fractional derivative, and they analyzed the stability and convergence by Fourier analysis [17]. Sun and Wu have derived the error estimate of the L1 formula which can approximate the Caputo derivative and derived a fully discrete scheme for the diffusion-wave equation [18]. Gao and Sun have derived the L1 approximation for the time-fractional derivative, and they constructed a compact finite difference scheme for the sub-diffusion equation [19]. Zhuang et al. constructed a Crank-Nicolson-type difference scheme for sub-diffusion equations which have a variable time, and they proved that this scheme is stable and convergent with the discrete H1 norm [13] [20]. Also, the maximum norm error estimate has been obtained [21]. In addition, more researchers found other numerical schemes such as the finite element method [22] and others. Tang has discussed convergence and superconvergence of fully discrete finite element for time fractional optimal control problems [23]. Wang et al. have derived the local discontinuous Galerkin method for the Time-Fractional KdV equation [24]. It is noted that only one item time fractional order is included in the study of time fractional sub-diffusion equation. In fractional physics, especially diffusion movement, the concept of brown movement is extended because of the generalization of the Gauss probability function. The scope of nuclear magnetic resonance is expanding by added resolving power, so one item time fractional order cannot explain this kind of problem.

The finite volume method is also called the control volume method. Let’s take a brief look at the idea of this method. First, we mesh the space and there are non-repetitive control volumes near each grid point. Second, we integrate the equations separately on each control volume. And then, we approximate the first-order partial derivatives with the function values of nodes. This method is usually used to solve integer order equations, and it can also be used to solve the fractional equation. The finite volume method has advantages of integral conservation [25].

Consider the following one-dimensional multi-term time fractional sub-diffusion equation with homogeneous source term on the interval $\left[\mathrm{0,}L\right]$

$\{\begin{array}{ll}{}_0{}^{C}D{}_t^{\alpha }u\left(x,t\right)+{}_0{}^{C}D{}_t^{\beta }u\left(x,t\right)=\frac{{\partial }^{2}u}{\partial x^2}+f\left(x,t\right),\hfill & \left(x,t\right)\in \left(0,L\right)\times \left(0,T\right],\hfill \\ u\left(0,t\right)=\mu \left(t\right),u\left(L,t\right)=\nu \left(t\right),\hfill & t\in \left[0,T\right],\hfill \\ u\left(x,0\right)=\varphi \left(x\right),\hfill & x\in \left[0,L\right],\hfill \end{array}$ (1)

where $0<\beta <\alpha <1$ and the operator ${}_0{}^{C}D{}_t^{\alpha }$ denotes the Caputo fractional derivative of order $\alpha$ defined by

${}_0{}^{C}D{}_t^{\alpha }f\left(t\right)=\frac{1}{\Gamma \left(1-\alpha \right)}{\displaystyle {\int }_0^t}\frac{f^{\prime }\left(\xi \right)}{{\left(t-\xi \right)}^{\alpha }}\text{d}\xi \mathrm{.}$

$\Gamma (\cdot )$ is the gamma function, and $\mu \left(t\right)\mathrm{,}\nu \left(t\right)\mathrm{,}\varphi \left(x\right)$ and $f\left(x\mathrm{,}t\right)$ are the known functions, $\mu \left(0\right)=\varphi \left(0\right)$, $\nu \left(0\right)=\varphi \left(L\right)$.

For the numerical approximation, take two positive integers M, N and let $h=\frac{L}{M},\tau =\frac{T}{N}$. Define $x_i=ih\left(0\le i\le M\right)$, $t_n=n\tau \left(0\le n\le N\right)$, ${\Omega }_h=\left\{x_i|0\le i\le M\right\}$, ${\Omega }_{\tau }=\left\{t_n|0\le i\le N\right\}$, then the computational domain $\left[0,L\right]\times \left[0,T\right]$ is covered by ${\Omega }_h\times {\Omega }_{\tau }$. Define primal partition $I_h$, grid element $I_i=\left[x_i\mathrm{,}{x}_{i+1}\right]\left(0\le i\le M-1\right)$. Define dual partition $I_h^{\ast }$, $I_0^{\ast }=\left[0,x_{\frac{1}{2}}\right]$, $I_i^{\ast }=\left[x_{i-\frac{1}{2}},x_{i+\frac{1}{2}}\right]\left(1\le i\le M-1\right)$, $I_M^{\ast }=\left[x_{M-\frac{1}{2}},L\right]$. Suppose $u=\left\{u_i^n|0\le i\le M,0\le n\le N\right\}$ is a grid function on ${\Omega }_{h\tau }={\Omega }_h\times {\Omega }_{\tau }$. For every grid function u, we define the following notations:

${\delta }_x{u}_{i-\frac{1}{2}}^n=\frac{1}{h}\left(u_i^n-u_{i-1}^n\right).$

The main goal of this paper is to construct a high-order compact finite scheme and establish the corresponding error estimate. The remainder of the article is arranged as follows. In Section 2, the compact finite volume scheme is derived. In Section 3, the existence and uniqueness, stability and convergence of the finite volume scheme are proved. In Section 4, two numerical examples are given to demonstrate the theoretical results. Finally, we obtain a brief conclusion.

2. The Derivation of the Compact Finite Volume Scheme

We need to follow some lemmas in the derivation of the compact finite volume scheme.

Lemma 1. If $g\left(x\right)\in C^3\left[\mathrm{0,}L\right]$ , and $x_{i+\frac{1}{2}}=\left(i+\frac{1}{2}\right)h$ , $0\le i\le M-1$ . Then

${\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}g\left(x\right)\text{d}x=\frac{h}{24}\left[g\left(x_{i-1}\right)+22g\left(x_i\right)+g\left(x_{i+1}\right)\right]+C{h}^4,$

where $C={\displaystyle {\int }_{-\frac{1}{2}}^{\frac{1}{2}}}\frac{g^{\left(3\right)}\left(\theta \right)}{6}\xi \left({\xi }^2-1\right)\text{d}\xi$ , $\theta \in \left(x_{i-\frac{1}{2}},x_{i+\frac{1}{2}}\right)$.

Proof. We use Lagrange interpolation to approximate $g\left(x\right)$ at the points $\left(x_{i-1}\mathrm{,}g\left(x_{i-1}\right)\right)\mathrm{,}\left(x_i\mathrm{,}g\left(x_i\right)\right)$ and $\left(x_{i+1}\mathrm{,}g\left(x_{i+1}\right)\right)$. We can obtain that

$g\left(x\right)=L_2\left(x\right)+R_2\left(x\right),$

where

$\begin{array}{c}{L}_2\left(x\right)=\frac{\left(x-x_i\right)\left(x-x_{i+1}\right)}{\left(x_{i-1}-x_i\right)\left(x_{i-1}-x_{i+1}\right)}g\left(x_{i-1}\right)+\frac{\left(x-x_{i-1}\right)\left(x-x_{i+1}\right)}{\left(x_i-x_{i-1}\right)\left(x_i-x_{i+1}\right)}g\left(x_i\right)\\ +\frac{\left(x-x_{i-1}\right)\left(x-x_i\right)}{\left(x_{i+1}-x_{i-1}\right)\left(x_{i+1}-x_i\right)}g\left(x_{i+1}\right)\end{array}$

and

$R_2\left(x\right)=\frac{g^{\left(3\right)}\left(\theta \right)}{6}\left(x-x_{i-1}\right)\left(x-x_i\right)\left(x-x_{i+1}\right),\theta \in \left(x_{i-\frac{1}{2}},x_{i+\frac{1}{2}}\right).$

Integrating $g\left(x\right)$ by parts, we can have

${\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}g\left(x\right)\text{d}x={\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}\left(L_2\left(x\right)+R_2\left(x\right)\right)\text{d}x.$

To obtain the approximation order, we use substitution to simplify integral terms. For example, let $\xi =\frac{x-ih}{h}$, so the integral terms can be simplified into

$L_2\left(\xi \right)=\frac{1}{2}\xi \left(\xi -1\right)g\left(x_{i-1}\right)-\left({\xi }^2-1\right)g\left(x_i\right)+\frac{1}{2}\xi \left(\xi +1\right)g\left(x_{i+1}\right)$

and

$R_2\left(\xi \right)=\frac{g^{\left(3\right)}\left(\theta \right)h^3}{6}\xi \left({\xi }^2-1\right)\mathrm{.}$

Such that

${\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}{L}_2\left(x\right)\text{d}x=h{\displaystyle {\int }_{-\frac{1}{2}}^{\frac{1}{2}}}{L}_2\left(\xi \right)\text{d}\xi =\frac{h}{24}\left(g\left(x_{i-1}\right)+22g\left(x_i\right)+g\left(x_{i+1}\right)\right)$

and

${\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}{R}_2\left(x\right)\text{d}x=h{\displaystyle {\int }_{-\frac{1}{2}}^{\frac{1}{2}}}{R}_2\left(\xi \right)\text{d}\xi =h^4{\displaystyle {\int }_{-\frac{1}{2}}^{\frac{1}{2}}}\frac{g^{\left(3\right)}\left(\theta \right)}{6}\xi \left({\xi }^2-1\right)\text{d}\xi .$ □

Define grid function space $U_h=\left\{g|g=\left(g_0,g_1,\cdots ,g_M\right)\right\}$. For every $g\in U_h$, define the integral operator as follows:

${\left(Sg\right)}_i=\frac{h}{24}\left(g_{i-1}+22g_i+g_{i+1}\right)\mathrm{,}\left(1\le i\le M-1\right)\mathrm{.}$

with, ${\left(Sg\right)}_0=g_0$, ${\left(Sg\right)}_M=g_M$.

Lemma 2. If $g\left(x\right)\in C^5\left[\mathrm{0,}L\right]$ , and $x_{i+\frac{1}{2}}=\left(i+\frac{1}{2}\right)h$ , $0\le i\le M-1$ . Then

$\frac{1}{24}\left[g^{\prime }\left(x_{i-\frac{3}{2}}\right)+22g^{\prime }\left(x_{i-\frac{1}{2}}\right)+g^{\prime }\left(x_{i+\frac{1}{2}}\right)\right]=\frac{1}{h}\left[g\left(x_i\right)-g\left(x_{i-1}\right)\right]+C{h}^4\mathrm{,}$

where

$\begin{array}{c}C={\displaystyle {\int }_0^1}\{\frac{1}{144}\left[g^{\left(5\right)}\left(x_{i-\frac{1}{2}}-\lambda h\right)+g^{\left(5\right)}\left(x_{i-\frac{1}{2}}+\lambda h\right)\right]{\left(1-\lambda \right)}^3\\ -\frac{1}{48}\left[g^{\left(5\right)}\left(x_{i-\frac{1}{2}}-\frac{\lambda }{2}h\right)+g^{\left(5\right)}\left(x_{i-\frac{1}{2}}+\frac{\lambda }{2}h\right)\right]{\left(1-\lambda \right)}^4\}\text{d}\lambda .\end{array}$

Especially,

$3g^{\prime }\left(x_{\frac{1}{2}}\right)-5g^{\prime }\left(x_{\frac{3}{2}}\right)+4g^{\prime }\left(x_{\frac{5}{2}}\right)-g^{\prime }\left(x_{\frac{7}{2}}\right)=\frac{1}{h}\left[g\left(x_1\right)-g\left(x_0\right)\right]+O\left(h^4\right),$

Proof. Based on the Taylor expansion, we can obtain the following equations

$\begin{array}{c}g\left(x_i\right)=g\left(x_{i-\frac{1}{2}}\right)+\frac{h}{2}{g}^{\prime }\left(x_{i-\frac{1}{2}}\right)+\frac{h^2}{8}{g}^{″}\left(x_{i-\frac{1}{2}}\right)+\frac{h^3}{48}{g}^{\left(3\right)}\left(x_{i-\frac{1}{2}}\right)\\ +\frac{h^4}{384}{g}^{\left(4\right)}\left(x_{i-\frac{1}{2}}\right)+\frac{1}{24}{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_i}}{g}^{\left(5\right)}\left(s\right){\left(x_i-s\right)}^4\text{d}s,\end{array}$ (2)

$\begin{array}{c}g\left(x_{i-1}\right)=g\left(x_{i-\frac{1}{2}}\right)-\frac{h}{2}{g}^{\prime }\left(x_{i-\frac{1}{2}}\right)+\frac{h^2}{8}{g}^{″}\left(x_{i-\frac{1}{2}}\right)-\frac{h^3}{48}{g}^{\left(3\right)}\left(x_{i-\frac{1}{2}}\right)\\ +\frac{h^4}{384}{g}^{\left(4\right)}\left(x_{i-\frac{1}{2}}\right)+\frac{1}{24}{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i-1}}}{g}^{\left(5\right)}\left(s\right){\left(x_{i-1}-s\right)}^4\text{d}s,\end{array}$ (3)

$\begin{array}{c}{g}^{\prime }\left(x_{i-\frac{3}{2}}\right)=g^{\prime }\left(x_{i-\frac{1}{2}}\right)-h{g}^{″}\left(x_{i-\frac{1}{2}}\right)+\frac{h^2}{2}{g}^{\left(3\right)}\left(x_{i-\frac{1}{2}}\right)-\frac{h^3}{6}{g}^{\left(4\right)}\left(x_{i-\frac{1}{2}}\right)\\ +\frac{1}{6}{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i-\frac{3}{2}}}}{g}^{\left(5\right)}\left(s\right){\left(x_{i-\frac{3}{2}}-s\right)}^3\text{d}s,\end{array}$ (4)

$\begin{array}{c}{g}^{\prime }\left(x_{i+\frac{1}{2}}\right)=g^{\prime }\left(x_{i-\frac{1}{2}}\right)+h{g}^{″}\left(x_{i-\frac{1}{2}}\right)+\frac{h^2}{2}{g}^{\left(3\right)}\left(x_{i-\frac{1}{2}}\right)+\frac{h^3}{6}{g}^{\left(4\right)}\left(x_{i-\frac{1}{2}}\right)\\ +\frac{1}{6}{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}{g}^{\left(5\right)}\left(s\right){\left(x_{i+\frac{1}{2}}-s\right)}^3\text{d}s.\end{array}$ (5)

To obtain the approximation order, we use substitution to simplify integral re-mainder. For example, let $s=\left(i-\frac{1}{2}\right)h+\frac{\lambda }{2}h$, so the integral remainder of (2) can be simplified into

$\frac{h^5}{48}{\displaystyle {\int }_0^1}{g}^{\left(5\right)}\left(x_{i-\frac{1}{2}}+\frac{\lambda }{2}h\right){\left(1-\lambda \right)}^4\text{d}\lambda \mathrm{.}$

We simplify the integral remainder of (2)-(5) by using the same way, then

$\begin{array}{l}\frac{1}{24}\left[g^{\prime }\left(x_{i-\frac{3}{2}}\right)+22g^{\prime }\left(x_{i-\frac{1}{2}}\right)+g^{\prime }\left(x_{i+\frac{1}{2}}\right)\right]-\frac{1}{h}\left[g\left(x_i\right)-g\left(x_{i-1}\right)\right]\\ =\frac{h^4}{144}{\displaystyle {\int }_0^1}\left[g^{\left(5\right)}\left(x_{i-\frac{1}{2}}-\lambda h\right)+g^{\left(5\right)}\left(x_{i-\frac{1}{2}}+\lambda h\right)\right]{\left(1-\lambda \right)}^3\text{d}\lambda \end{array}$

$\begin{array}{l}-\frac{h^4}{48}{\displaystyle {\int }_0^1}\left[g^{\left(5\right)}\left(x_{i-\frac{1}{2}}-\frac{\lambda }{2}h\right)+g^{\left(5\right)}\left(x_{i-\frac{1}{2}}+\frac{\lambda }{2}h\right)\right]{\left(1-\lambda \right)}^4\text{d}\lambda \\ =h^4{\displaystyle {\int }_0^1}\{\frac{1}{144}\left[g^{\left(5\right)}\left(x_{i-\frac{1}{2}}-\lambda h\right)+g^{\left(5\right)}\left(x_{i-\frac{1}{2}}+\lambda h\right)\right]{\left(1-\lambda \right)}^3\\ -\frac{1}{48}\left[g^{\left(5\right)}\left(x_{i-\frac{1}{2}}-\frac{\lambda }{2}h\right)+g^{\left(5\right)}\left(x_{i-\frac{1}{2}}+\frac{\lambda }{2}h\right)\right]{\left(1-\lambda \right)}^4\}\text{d}\lambda .\end{array}$ □

Define fractional index grid function space ${\bar{U}}_h=\left\{g\mathrm{<mo>\; |\; </mo>}g=\left(g_{\frac{1}{2}}\mathrm{,}{g}_{\frac{3}{2}}\mathrm{,}\cdots \mathrm{,}{g}_{M-\frac{1}{2}}\right)\right\}$. For every $g_{i-\frac{1}{2}}\in {\bar{U}}_h$, define compact operator as follows:

${\left(Ag\right)}_{i-\frac{1}{2}}=\frac{1}{24}\left(g_{i-\frac{3}{2}}+22g_{i-\frac{1}{2}}+g_{i+\frac{1}{2}}\right)\mathrm{,}\left(2\le i\le M-1\right)\mathrm{.}$

with, ${\left(Ag\right)}_0=g_0$, ${\left(Ag\right)}_M=g_M$.

Lemma 3. [26] If $0<\alpha <1$ , $g\in C^2\left[\mathrm{0,}{t}_n\right]$ . Then,

${}_0{}^{C}D{}_t^{\alpha }g\left(t_n\right)=\frac{{\tau }^{-\alpha }}{\Gamma \left(2-\alpha \right)}\left[a_0^{\alpha }g\left(t_n\right)-{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)g\left(t_k\right)-a_{n-1}^{\alpha }g\left(0\right)\right]+R_{t_n}^{\alpha }$

and

$\left|R_{t_n}^{\alpha }\right|\le \frac{1}{\Gamma \left(2-\alpha \right)}\left[\frac{1}{4}+\frac{\alpha }{\left(1-\alpha \right)\left(2-\alpha \right)}\right]\underset{0\le t\le t_n}{\max }\left|g^{″}\left(t\right)\right|{\tau }^{2-\alpha }\mathrm{,}$

where $a_k^{\alpha }={\left(k+1\right)}^{1-\alpha }-k^{1-\alpha }$.

Define the L1 approximation operator as follows

$D_{\tau }^{\alpha }{g}_i^n=\frac{{\tau }^{-\alpha }}{\Gamma \left(2-\alpha \right)}\left[a_0^{\alpha }g\left(x\left(i\right),t_n\right)-{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)g\left(x\left(i\right),t_k\right)-a_{n-1}^{\alpha }g\left(x\left(i\right),0\right)\right],$

$0\le i\le M\mathrm{,}0\le n\le N\mathrm{,}$

where the definition of $a_k^{\alpha }$ is as same as the Lemma 3.

Let us now construct a compact finite volume scheme for problem (1). On ${\Omega }_h\times {\Omega }_{\tau }$, we now define the grid functions

$U_i^n=u\left(x_i,t_n\right),f_i^n=f\left(x_i,t_n\right),0\le i\le M,0\le n\le N.$

Suppose $u\left(x\mathrm{,}t\right)\in C_{x\mathrm{,}t}^{\mathrm{5,2}}\left(\left[\mathrm{0,}L\right]\times \left[\mathrm{0,}T\right]\right)$, which symbols of $u\in C^5\left(x\right)$ and $u\in C^2\left(t\right)$ meanwhile. We consider the Equation (1) at the point $\left(x\mathrm{,}{t}_n\right)$, and we can have

${}_0{}^{C}D{}_t^{\alpha }u\left(x,t_n\right){+}_0^C{D}_t^{\beta }u\left(x,t_n\right)=\frac{{\partial }^{2}u\left(x,t_n\right)}{\partial x^2}+f\left(x,t_n\right),x\in \left(0,L\right),1\le n\le N.$ (6)

Integrating (6) in intervals, we can have

$\begin{array}{l}{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}\left({}_0{}^{C}D{}_t^{\alpha }u\left(x,t_n\right)+{}_0{}^{C}D{}_t^{\beta }u\left(x,t_n\right)\right)\text{d}x\\ ={\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}\left(\frac{{\partial }^{2}u\left(x,t_n\right)}{\partial x^2}+f\left(x,t_n\right)\right)\text{d}x,1\le i\le M-1,1\le n\le N.\end{array}$ (7)

We can obtain

$\begin{array}{l}{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}\left({}_0{}^{C}D{}_t^{\alpha }u\left(x,t_n\right)+{}_0{}^{C}D{}_t^{\beta }u\left(x,t_n\right)\right)\text{d}x\\ =\frac{\partial u\left(x_{i+\frac{1}{2}},t_n\right)}{\partial x}-\frac{\partial u\left(x_{i-\frac{1}{2}},t_n\right)}{\partial x}+{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}f\left(x,t_n\right)\text{d}x,1\le i\le M-1,1\le n\le N.\end{array}$ (8)

By using Lemma 2, we can get

$\begin{array}{l}A{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}\left({}_0{}^{C}D{}_t^{\alpha }u\left(x,t_n\right)+{}_0{}^{C}D{}_t^{\beta }u\left(x,t_n\right)\right)\text{d}x\\ ={\delta }_x{U}_{i+\frac{1}{2}}^n-{\delta }_x{U}_{i-\frac{1}{2}}^n+A{\displaystyle {\int }_{x_{i-\frac{1}{2}}}^{x_{i+\frac{1}{2}}}}f\left(x,t_n\right)\text{d}x+r_i^n,1\le i\le M-1,1\le n\le N,\end{array}$ (9)

where $r_i^n=C{h}^4$ and the definition of operator A is as same as Lemma 2.

Using numerical integration formula to deal with spatial integral and L1 interpolation to discretize the time fractional derivative, then we have from Lemma 1 and Lemma 3.

$\begin{array}{l}\frac{1}{24}S\left(D_{\tau }^{\alpha }{U}_{i-1}^n+D_{\tau }^{\beta }{U}_{i-1}^n\right)+\frac{22}{24}S\left(D_{\tau }^{\alpha }{U}_i^n+D_{\tau }^{\beta }{U}_i^n\right)+\frac{1}{24}S\left(D_{\tau }^{\alpha }{U}_{i+1}^n+D_{\tau }^{\beta }{U}_{i+1}^n\right)\\ ={\delta }_x{U}_{i+\frac{1}{2}}^n-{\delta }_x{U}_{i-\frac{1}{2}}^n+\frac{1}{24}S{f}_{i-1}^n+\frac{22}{24}S{f}_i^n+\frac{1}{24}S{f}_{i+1}^n+R_i^n,1\le n\le N,\end{array}$ (10)

where $R_i^n=C_1\left({\tau }^{2-\alpha }\right)+C_2\left(h^4\right)$ and $C_1\mathrm{,}{C}_2$ are constants which are independent of $h\mathrm{,}\tau$. The definition of operator S is as same as Lemma 3.

Notice that $u\left(0,t\right)=\mu \left(t\right),u\left(L,t\right)=\nu \left(t\right)$ and $u\left(x,0\right)=\varphi \left(x\right)$, so we can have

$U_i^0=\varphi \left(x_i\right),1\le i\le M-1,$ (11)

$U_0^n=\mu \left(t_n\right),U_M^n=\nu \left(t_n\right),0\le n\le N.$ (12)

Therefore, we leave out the infinitesimal, then

$\begin{array}{l}\frac{1}{24}S\left(D_{\tau }^{\alpha }{u}_{i-1}^n+D_{\tau }^{\beta }{u}_{i-1}^n\right)+\frac{22}{24}S\left(D_{\tau }^{\alpha }{u}_i^n+D_{\tau }^{\beta }{u}_i^n\right)+\frac{1}{24}S\left(D_{\tau }^{\alpha }{u}_{i+1}^n+D_{\tau }^{\beta }{u}_{i+1}^n\right)\\ ={\delta }_x{u}_{i+\frac{1}{2}}^n-{\delta }_x{u}_{i-\frac{1}{2}}^n+\frac{1}{24}S{f}_{i-1}^n+\frac{22}{24}S{f}_i^n+\frac{1}{24}S{f}_{i+1}^n,1\le n\le N,\end{array}$ (13)

$u_i^0=\varphi \left(x_i\right),1\le i\le M-1,$ (14)

$u_0^n=\mu \left(t_n\right),u_M^n=\nu \left(t_n\right),0\le n\le N.$ (15)

These above equations are the compact finite volume scheme of question (1).

3. Analysis of the Compact Finite Volume Scheme

In this section, we will prove the existence and uniqueness, stability and convergence.

First, we introduce the norms in the space $\gamma$. Let $g=\left(g_0\mathrm{,}{g}_1\mathrm{,}\cdots \mathrm{,}{g}_M\right)\in \gamma$. Denote

${\Vert g\Vert }_{\infty }=\underset{0<i<M}{\max }\left|g_i\right|.$

Now, we introduce some important lemmas which can be used in the following verifications.

Lemma 4. Suppose $\alpha \in \left(\mathrm{0,1}\right)$ , and $a_l^{\alpha }$ is defined by Lemma 3, $l=0,1,\cdots$ , then

1) $1=a_0^{\alpha }>a_1^{\alpha }>a_2^{\alpha }>\cdots >a_l^{\alpha }$ ; $a_l^{\alpha }\to 0$, when $l\to \infty$ ;

2) $\left(1-\alpha \right)l^{-\alpha }<a_{l-1}^{\alpha }<\left(1-\alpha \right){\left(l-1\right)}^{-\alpha }$, $l\ge 1$.

Lemma 5. Let $p=\left(p_0,p_1,\cdots ,p_M\right)\in \gamma$ , $q=\left(q_0,q_1,\cdots ,q_M\right)\in \gamma$ , then

${\Vert p+q\Vert }_{\infty }\le {\Vert p\Vert }_{\infty }+{\Vert q\Vert }_{\infty }\mathrm{.}$

Theorem 1. The solution of the compact finite volume scheme (13)-(15) is existent and unique.

Proof. Suppose $u^n=\left(u_0^n\mathrm{,}{u}_1^n\mathrm{,}\cdots \mathrm{,}{u}_M^n\right)\in \gamma$. The numerical solution of $u^0$ can be obtained by (14). If the numerical solutions of $u^0\mathrm{,}{u}^1\mathrm{,}\cdots \mathrm{,}{u}^{n-1}$ are existing and unique, we can obtain non-homogeneous linear equations about $u^n$ from (13) and (15). So if we prove the existence and uniqueness of $u^n$, we only need to prove that homogeneous linear equations only have a zero solution. Define $S_{\alpha }={\tau }^{\alpha }\Gamma \left(2-\alpha \right)$, $S_{\beta }={\tau }^{\beta }\Gamma \left(2-\beta \right)$.

$\begin{array}{l}\left(\frac{S}{24S_{\alpha }}+\frac{S}{24S_{\beta }}\right)u_{i-1}^n+\left(\frac{22S}{24S_{\alpha }}+\frac{22S}{24S_{\beta }}\right)u_i^n+\left(\frac{S}{24S_{\alpha }}+\frac{S}{24S_{\beta }}\right)u_{i+1}^n\\ ={\delta }_x{u}_{i+\frac{1}{2}}^n-{\delta }_x{u}_{i-\frac{1}{2}}^n,1\le i\le M-1,\end{array}$ (16)

$u_0^n=u_M^n=0.$ (17)

Now, we prove that (13) and (15) only have a zero solution. Equation (16) can be rewritten into

$\begin{array}{l}\left(\frac{S}{24S_{\alpha }}+\frac{S}{24S_{\beta }}\right)u_{i-1}^n+\left(\frac{22S}{24S_{\alpha }}+\frac{22S}{24S_{\beta }}+\frac{2}{h}\right)u_i^n+\left(\frac{S}{24S_{\alpha }}+\frac{S}{24S_{\beta }}\right)u_{i+1}^n\\ =\frac{1}{h}\left(u_{i-1}^n+u_{i+1}^n\right),1\le i\le M-1.\end{array}$

According to Lemma 5, equation mentioned above also can be transformed into

$\left(\frac{S}{S_{\alpha }}+\frac{S}{S_{\beta }}+\frac{2}{h}\right){\Vert u^n\Vert }_{\infty }=\frac{2}{h}{\Vert u^n\Vert }_{\infty }.$

So ${\Vert u^n\Vert }_{\infty }=0$, and we can obtain $u^n=0$.

According to inductive principle, Equations (13)-(15) have a unique solution. □

Theorem 2. Suppose the finite volume scheme (13)-(15) has a solution. We record it as $\left\{v_i^n\mathrm{<mo>\; |\; </mo>0}\le i\le M\mathrm{,0}\le n\le N\right\}$ , then

${\Vert v^k\Vert }_{\infty }\le {\Vert v^0\Vert }_{\infty }+c\underset{1\le m\le k}{\max }{\Vert f^m\Vert }_{\infty },1\le k\le N,$ (18)

where $c=\frac{1}{2}\max \left\{T^{\alpha }\Gamma \left(1-\alpha \right),T^{\beta }\Gamma \left(1-\beta \right)\right\}$, ${\Vert f^m\Vert }_{\infty }=\underset{1\le i\le M-1}{\max }\left|f_i^m\right|$.

Proof. Equation (13) can be written into

$\begin{array}{l}\frac{S}{24S_{\alpha }}\left[a_0^{\alpha }{v}_{i-1}^n-{\displaystyle \underset{k\mathrm{<mo>\; =\; </mo>1}}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)v_{i-1}^k-a_{n-1}^{\alpha }{v}_{i-1}^0\right]\\ +\frac{S}{24S_{\beta }}\left[a_0^{\beta }{v}_{i-1}^n-{\displaystyle \underset{k\mathrm{<mo>\; =\; </mo>1}}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)v_{i-1}^k-a_{n-1}^{\beta }{v}_{i-1}^0\right]\\ +\frac{22S}{24S_{\alpha }}\left[a_0^{\alpha }{v}_i^n-{\displaystyle \underset{k\mathrm{<mo>\; =\; </mo>1}}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)v_i^k-a_{n-1}^{\alpha }{v}_i^0\right]\\ +\frac{22S}{24S_{\beta }}\left[a_0^{\beta }{v}_i^n-{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)v_i^k-a_{n-1}^{\beta }{v}_i^0\right]\end{array}$

$\begin{array}{l}+\frac{S}{24S_{\alpha }}\left[a_0^{\alpha }{v}_{i+1}^n-{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)v_{i+1}^k-a_{n-1}^{\alpha }{v}_{i+1}^0\right]\\ +\frac{S}{24S_{\beta }}\left[a_0^{\beta }{v}_{i+1}^n-{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)v_{i+1}^k-a_{n-1}^{\beta }{v}_{i+1}^0\right]\\ ={\delta }_x{v}_{i+\frac{1}{2}}^n-{\delta }_x{v}_{i-\frac{1}{2}}^n+\frac{1}{24}S{f}_{i-1}^n+\frac{22}{24}S{f}_i^n+\frac{1}{24}S{f}_{i+1}^n,\\ 1\le i\le M-1,1\le n\le N.\end{array}$

Equation mentioned above also can be written into

$\begin{array}{l}\left(\frac{S{a}_0^{\alpha }}{24S_{\alpha }}+\frac{S{a}_0^{\beta }}{24S_{\beta }}\right)v_{i-1}^n+\left(\frac{22S{a}_0^{\alpha }}{24S_{\alpha }}+\frac{22S{a}_0^{\beta }}{24S_{\beta }}+\frac{2}{h}\right)v_i^n+\left(\frac{S{a}_0^{\alpha }}{24S_{\alpha }}+\frac{S{a}_0^{\beta }}{24S_{\beta }}\right)v_{i+1}^n\\ =\frac{S}{24S_{\alpha }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)v_{i-1}^k+a_{n-1}^{\alpha }{v}_{i-1}^0\right]\\ +\frac{S}{24S_{\beta }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)v_{i-1}^k+a_{n-1}^{\beta }{v}_{i-1}^0\right]\\ +\frac{22S}{24S_{\alpha }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)v_i^k+a_{n-1}^{\alpha }{v}_i^0\right]\end{array}$

$\begin{array}{l}+\frac{22S}{24S_{\beta }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)v_i^k+a_{n-1}^{\beta }{v}_i^0\right]\\ +\frac{S}{24S_{\alpha }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)v_{i+1}^k+a_{n-1}^{\alpha }{v}_{i+1}^0\right]\\ +\frac{S}{24S_{\beta }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)v_{i+1}^k+a_{n-1}^{\beta }{v}_{i+1}^0\right]\\ +\frac{1}{h}\left(v_{i-1}^n+v_{i+1}^n\right)+\frac{1}{24}S{f}_{i-1}^n+\frac{22}{24}S{f}_i^n+\frac{1}{24}S{f}_{i+1}^n,\\ 1\le i\le M-1,1\le n\le N.\end{array}$

And then, we add absolute values to both sides of the equation. By using Lemma 5, we can obtain

$\begin{array}{l}\left(S\frac{a_0^{\alpha }}{S_{\alpha }}+S\frac{a_0^{\beta }}{S_{\beta }}+\frac{2}{h}\right){\Vert v^n\Vert }_{\infty }\\ \le S\frac{1}{S_{\alpha }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right){\Vert v^k\Vert }_{\infty }+a_{n-1}^{\alpha }{\Vert v^0\Vert }_{\infty }\right]\end{array}$

$\begin{array}{l}+S\frac{1}{S_{\beta }}\left[{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right){\Vert v^k\Vert }_{\infty }+a_{n-1}^{\beta }{\Vert v^0\Vert }_{\infty }\right]+\frac{2}{h}{\Vert v^n\Vert }_{\infty }+S{\Vert f^n\Vert }_{\infty }\\ ={\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left[S\frac{1}{S_{\alpha }}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)+S\frac{1}{S_{\beta }}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)\right]{\Vert v^k\Vert }_{\infty }\\ +\left[S\frac{a_{n-1}^{\alpha }}{S_{\alpha }}+S\frac{a_{n-1}^{\beta }}{S_{\beta }}\right]{\Vert v^0\Vert }_{\infty }+\frac{2}{h}{\Vert v^n\Vert }_{\infty }+S\frac{a_{n-1}^{\alpha }}{S_{\alpha }}\cdot \frac{S_{\alpha }}{2a_{n-1}^{\alpha }}{\Vert f^n\Vert }_{\infty }\\ +S\frac{a_{n-1}^{\beta }}{S_{\beta }}\cdot \frac{S_{\beta }}{2a_{n-1}^{\beta }}{\Vert f^n\Vert }_{\infty }.\end{array}$ (19)

By using Lemma 4, we can obtain

$\frac{S_{\alpha }}{2a_{n-1}^{\alpha }}\le \frac{{\tau }^{\alpha }\Gamma \left(2-\alpha \right)}{2h\left(1-\alpha \right)n^{-\alpha }}=\frac{t_n^{\alpha }\Gamma \left(1-\alpha \right)}{2h},$ (20)

$\frac{S_{\beta }}{2a_{n-1}^{\beta }}\le \frac{{\tau }^{\beta }\Gamma \left(2-\beta \right)}{2h\left(1-\beta \right)n^{-\beta }}=\frac{t_n^{\beta }\Gamma \left(1-\beta \right)}{2h}.$ (21)

Bring (20) and (21) into (19), then

$\begin{array}{l}\left(A\frac{a_0^{\alpha }}{S_{\alpha }}+A\frac{a_0^{\beta }}{S_{\beta }}\right){\Vert v^n\Vert }_{\infty }\\ \le {\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left[S\frac{1}{S_{\alpha }}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)+S\frac{1}{S_{\beta }}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)\right]{\Vert v^k\Vert }_{\infty }\\ +\left(S\frac{a_{n-1}^{\alpha }}{S_{\alpha }}+S\frac{a_{n-1}^{\beta }}{S_{\beta }}\right)\left({\Vert v^0\Vert }_{\infty }+c{\Vert f^n\Vert }_{\infty }\right),1\le n\le N,\end{array}$ (22)

where $c=\frac{1}{2}\max \left\{T^{\alpha }\Gamma \left(1-\alpha \right),T^{\beta }\Gamma \left(1-\beta \right)\right\}$.

Next, we use mathematical induction to prove (18).

According to (22), we can obtain when $n=1$,

$\left(S\frac{a_0^{\alpha }}{S_{\alpha }}+S\frac{a_0^{\beta }}{S_{\beta }}\right){\Vert v^1\Vert }_{\infty }\le \left(S\frac{a_0^{\alpha }}{S_{\alpha }}+S\frac{a_0^{\beta }}{S_{\beta }}\right)\left({\Vert v^0\Vert }_{\infty }+c{\Vert f^1\Vert }_{\infty }\right)\mathrm{.}$

It also can be written as

${\Vert v^1\Vert }_{\infty }={\Vert v^0\Vert }_{\infty }+c{\Vert f^1\Vert }_{\infty }.$

Therefore when $k=1$, Equation (18) satisfies conditions.

Suppose Equation (22) also satisfies conditions when $k=1,2,\cdots ,n-1$. According to (22), we can obtain

$\begin{array}{l}\left(S\frac{a_0^{\alpha }}{S_{\alpha }}+S\frac{a_0^{\beta }}{S_{\beta }}\right){\Vert v^n\Vert }_{\infty }\\ \le {\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left[S\frac{1}{S_{\alpha }}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)+S\frac{1}{S_{\beta }}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)\right]\left({\Vert v^0\Vert }_{\infty }+c\underset{1\le m\le k}{\max }{\Vert f^m\Vert }_{\infty }\right)\\ +\left(S\frac{a_{n-1}^{\alpha }}{S_{\alpha }}+S\frac{a_{n-1}^{\beta }}{S_{\beta }}\right)\left({\Vert v^0\Vert }_{\infty }+c{\Vert f^n\Vert }_{\infty }\right)\end{array}$

$\begin{array}{l}\le \{{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left[S\frac{1}{S_{\alpha }}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)+S\frac{1}{S_{\beta }}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)\right]\\ +\left(S\frac{a_{n-1}^{\alpha }}{S_{\alpha }}+S\frac{a_{n-1}^{\beta }}{S_{\beta }}\right)\}\cdot \left({\Vert u^0\Vert }_{\infty }+c\underset{1\le m\le n}{\max }{\Vert f^m\Vert }_{\infty }\right)\\ =\left(S\frac{a_0^{\alpha }}{S_{\alpha }}+S\frac{a_0^{\beta }}{S_{\beta }}\right)\left({\Vert v^0\Vert }_{\infty }+c{\Vert f^n\Vert }_{\infty }\right).\end{array}$

It also can be written as

${\Vert v^n\Vert }_{\infty }\le {\Vert v^0\Vert }_{\infty }+c{\Vert f^n\Vert }_{\infty }\mathrm{.}$

So when $k=n$, Equation (18) satisfies conditions.

According to the inductive principle, theorem 2 holds. □

Theorem 3. The solution of the compact finite volume scheme (13)-(15) is convergent. Suppose $\left\{U_i^n\mathrm{<mo>\; |\; </mo>0}\le i\le M\mathrm{,0}\le n\le N\right\}$ is the solution of equation (1), $\left\{u_i^n\mathrm{<mo>\; |\; </mo>0}\le i\le M\mathrm{,0}\le n\le N\right\}$ is the solution of the compact finite volume scheme (13)-(17). Define

$e_i^n=U_i^n-u_i^n,0\le i\le M,0\le n\le N,$

then

${\Vert e^n\Vert }_{\infty }\le c\left({\tau }^{2-\alpha }+h^4\right)\mathrm{,}1\le n\le N\mathrm{.}$ (23)

Proof. We subtract (18)-(22) from (13)-(15), then we can obtain error equations.

$\begin{array}{l}S\frac{1}{S_{\alpha }}\left[a_0^{\alpha }{e}_i^n-{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\alpha }-a_{n-k}^{\alpha }\right)e_i^k-a_{n-1}^{\alpha }{e}_i^0\right]\\ +S\frac{1}{S_{\beta }}\left[a_0^{\beta }{e}_i^n-{\displaystyle \underset{k=1}{\overset{n-1}{\sum }}}\left(a_{n-k-1}^{\beta }-a_{n-k}^{\beta }\right)e_i^k-a_{n-1}^{\beta }{e}_i^0\right]\\ ={\delta }_x{e}_{i+\frac{1}{2}}^n-{\delta }_x{e}_{i-\frac{1}{2}}^n+r_i^n,1\le i\le M-1,1\le n\le N,\end{array}$ (24)

$e_i^0=0,1\le i\le M-1,$ (25)

$e_0^n=0,e_M^n=0,0\le n\le N.$ (26)

By using Theorem 2, we can obtain

${\Vert e^n\Vert }_{\infty }\le {\Vert e^0\Vert }_{\infty }+c\underset{1\le k\le n}{\max }{\Vert r^k\Vert }_{\infty }\mathrm{,}1\le k\le N\mathrm{,}$

where $c=\frac{1}{2}\max \left\{T^{\alpha }\Gamma \left(1-\alpha \right),T^{\beta }\Gamma \left(1-\beta \right)\right\}$.

According to Lemma 5, it’s easy to obtain

${\Vert r^k\Vert }_{\infty }\le c_1\left({\tau }^{2-\alpha }+h^4\right)\mathrm{.}$

Such that, the result

${\Vert e^n\Vert }_{\infty }\le c\left({\tau }^{2-\alpha }+h^4\right)\mathrm{,}1\le n\le N\mathrm{,}$

is proved.

4. Numerical Experiments

In this section, we report on some numerical results to show the convergence orders and effectiveness of our finite volume scheme.

Example 1. Let $L=1$, $T=1$. In order to obtain the order of convergence of the finite volume scheme, we refer to the exact solution of the problem (1) is

$u\left(x,t\right)=t^2\sin \left(\pi x\right).$

We can obtain the corresponding source term $f\left(x\mathrm{,}t\right)$ and the initial and boundary conditions with $\alpha =0.5$, $\beta =0.2$, which are respectively

$f\left(x,t\right)=\frac{\Gamma \left(3\right)\sin \left(\pi x\right)}{\Gamma \left(3-\alpha \right)}{t}^{2-\alpha }+\frac{\Gamma \left(3\right)\sin \left(\pi x\right)}{\Gamma \left(3-\beta \right)}{t}^{2-\beta }+{\pi }^2{t}^2\sin (\; \pi \; x\; )$

and

$\mu \left(t\right)=0,\nu \left(t\right)=0,\varphi \left(x\right)=0.$

Denote the maximum errors

$E_{\infty }\left(h\mathrm{,}\tau \right)=\underset{1\le i\le M}{\max }\left|U_i^N-u_i^N\right|\mathrm{.}$

Table 1 and Table 2 show the maximum errors and the convergence order of the finite volume scheme. Suppose

$E_{\infty }\left(h,\tau \right)=O\left(h^p+{\tau }^q\right).$

If $\tau$ is small enough, then $E_{\infty }\left(h\mathrm{,}\tau \right)\approx O\left(h^p\right)$. Consequently,

$\frac{E_{\infty }\left(2h,\tau \right)}{E_{\infty }\left(h,\tau \right)}\approx 2^p,\log 2\frac{E_{\infty }\left(2h,\tau \right)}{E_{\infty }\left(h,\tau \right)}\approx p.$

So p is the convergence order with respect to the spatial step-size. Similarly, we can obtain

$q\approx \log 2\frac{E_{\infty }\left(h,2\tau \right)}{E_{\infty }\left(h,\tau \right)}$

for small enough h. Denote

$Order1=\log 2\frac{E_{\infty }\left(h,2\tau \right)}{E_{\infty }\left(h,\tau \right)},Order2=\frac{E_{\infty }\left(2h,\tau \right)}{E_{\infty }\left(h,\tau \right)}.$

In order to verify the integral conservation of the scheme, we define $Mas{s}_h={\displaystyle \underset{i=0}{\overset{M}{\sum }}}{u}_i^{N}h$. Under different meshes, $Mas{s}_h$ keep about the same size, which we can explain the integral conservation.

In Table 1 and Table 2, we compare the exact solution and the numerical solution. In order to test the convergence order of our scheme in temporal direction, we fix sufficiently small spatial step size $h=\frac{1}{1000}$ and vary the temporal step sizes. Table 1 list the numerical results for different temporal step sizes. In order to test the convergence order of our finite volume scheme in spatial direction, we fix sufficiently small temporal step sizes $\tau =\frac{1}{1000}$ and vary the spatial step sizes. Table 2 list the numerical results for different spatial step sizes. Figure 1 show the effect of the numerical solution and exact solution at fixed $h=\frac{1}{1000}$ and $\tau =\frac{1}{64}$. Figure 2 show the effect of the numerical solution and exact solution at

fixed $\tau =\frac{1}{1000}$ and $h=\frac{1}{64}$. In Figure 3, the red line represents temporal error order and the blue line represents spatial error order.

In order to observe error orders more intuitively, we plot a figure about error orders which the slope represents the error order. We can observe that the temporal error order is about $2-\alpha$ and the spatial error order is about $h^4$.

Example 2. Let $L=1$, $T=1$. As before, we refer to the exact solution of the problem (1) is

$u\left(x,t\right)=t^{2}x\left(1-x\right).$

It is also not difficult to obtain the corresponding source term $f\left(x,t\right)$ and the initial and boundary conditions with $\alpha =0.9$, $\beta =0.1$, which are respectively

$f\left(x,t\right)=\frac{\Gamma \left(3\right)x\left(1-x\right)}{\Gamma \left(3-\alpha \right)}{t}^{2-\alpha }+\frac{\Gamma \left(3\right)x\left(1-x\right)}{\Gamma \left(3-\beta \right)}{t}^{2-\beta }+2t^2$

and

$\mu \left(t\right)=0,\nu \left(t\right)=0,\varphi \left(x\right)=0.$

In order to test the order of convergence of the finite volume scheme, we compute the maximum norm errors of the numerical solution, which is defined as same as Example 1. Denote

$Order1=\log 2\frac{E_{\infty }\left(h,2\tau \right)}{E_{\infty }\left(h,\tau \right)},Order2=\frac{E_{\infty }\left(2h,\tau \right)}{E_{\infty }\left(h,\tau \right)}.$

In Table 3, we compute the problem for a longer time by fixing

$N=8,16,32,64,128$, and still choosing $h=\frac{1}{1000}$. In Table 4, we compute the problem for a longer space by fixing $M=8,16,32,64,128$, and still choosing $\tau =\frac{1}{1000}$. Figure 4 shows the effect of the numerical solution and exact solution at fixed $h=\frac{1}{1000}$ and $\tau =\frac{1}{64}$. Figure 5 shows the effect of the numerical solution and exact solution at fixed $\tau =\frac{1}{1000}$ and $h=\frac{1}{64}$. In Figure 6, the red line represents temporal error order and the blue line represents spatial error order.

In Example 1, we plot a figure about error orders which slope represents the error order. We can observe that the temporal error order is about $2-\alpha$ and the spatial error order is about $h^4$.

According to these tables, we can obtain that the compact finite volume scheme is convergent with different spatial step sizes and time step sizes. Under

different space steps, $Mas{s}_h$ keep about the same size, which we can get a conclusion that the scheme is conserved.

5. Conclusion

In this article, we constructed a compact finite volume scheme for the multi-term time fractional sub-diffusion equation. Indeed, we use some important operators to deal with the multi-term time fractional sub-diffusion equation. By using a compact operator, we obtain a high order accuracy scheme. We proved the existence and uniqueness, stability and $L_{\infty }$ convergence of our scheme. Two numerical results show that the scheme is conserved and convergent with the order $O\left({\tau }^{2-\alpha }+h^4\right)$. Some tables and figures also can prove the compact finite volume scheme is stable and convergent.

Funding

This work is supported by the Project of Shandong Education Development Promotion Association and the Project of Dezhou Institute of Educational Sciences.

Authors’ Contributions

Baojin Su carried out the main part of this article. All authors read and approved the final manuscript.

Acknowledgements

The authors would like to express their sincere thanks to the referees for the valuable comments and suggestions which helped to improve the original paper.

Supported

This work is supported by Project of Shandong Education Development Promotion Association, Project of Dezhou Institute of Educational Sciences.