Malika Boutiba | Selma Baghli-Bendimerad* | Abbès Benaïssa
© 2022 IIETA. This article is published by IIETA and is licensed under the CC BY 4.0 license (http://creativecommons.org/licenses/by/4.0/).
OPEN ACCESS
In this paper, we apply to a class of partial differential equation the finite element method when the problem is involving the Riemann-Liouville fractional derivative for time and space variables on a bounded domain with bounded conditions. The studied equation is obtained from the standard time diffusion equation by replacing the first order time derivative by a for 0<a<1 and for the second standard order space derivative by b for 1<b<2 respectively. The existence of the unique solution is proved by the Lax-Milgram Lemma. We present here three schemes to approximate numerically the time derivative and use the finite element method for the space derivative using the Hadamard finite part integral and the Diethlem's first degree compound quadrature formula, the second approach is based on the link between Riemann-Liouville and Caputo fractional derivative, when the third method was based on the approximation of the Riemann-Liouville by the Grunwald-Letnikov fractional derivative. For the approximation of the space fractional derivative, the finite element method is introduced for all the three approaches. Finally, to check the effectiveness of the three methods, a numerical example was given.
partial differential equations, Riemann-Liouville fractional derivative, finite element method, Lax-Milgram Lemma, Hilbert spaces, numerical analysis
We are interesting here by the existence of the unique solution of the Riemann-Liouville fractional derivative of the partial differential equations of the form.
${ }_0^{R L} D_t^\alpha u(x, t)={ }_0^{R L} D_x^\beta u(x, t)+f(x, t)$$(x, t) \in 0:=[0,1] \times[0, T]$ (1)
under the initial and border conditions:
$u(x, 0)=u_0(x), x \in \Omega:=[0,1]$ (2)
$u(0, t)=u(1, t)=0, t \in I:=[0, T]$ (3)
where, 0<a<1 and 1<b<2, T>0, ${ }_0^{R L} D_t^\alpha u(x, t)$ and ${ }_0^{R L} D_x^\beta u(x, t)$ denote the left Riemann-Liouville fractional derivative of order a and b respectively, $f: \Omega \times I \rightarrow \mathbb{R}$ and $u_0: \Omega \rightarrow \mathbb{R}$ are given functions.
Fractional partial differential equations have many applications in various fields, for example electro-magnetic, viscoelastic mechanics, fractal media, mathematical biology, and chemistry and so on. Analytical solution of these problems has been studies using Fourier and Laplace transforms, Green's functions [1-7].
Some different numerical solution methods are proposed to resolve the space and time fractional partial differential equations [8-23].
Meerschaert and Tadjeran [24] give a finite difference approximation to resolve the space fractional dispersion-advection equation involving Riemann-Liouville fractional derivative.
Zeng et al. [19] proposed two finite difference element methods to approximate the time-fractional sub diffusion equation with the Caputo fractional derivative. Li and Xu [25] proposed a finite difference spectral approximation method for the time-fractional derivative for diffusion equations. Zheng et al. [26] proposed a note on the finite element method for the space-fractional derivative for diffusion equation. Zheng et al. [21] proposed a novel high order for space-time spectral method to the time-fractional derivative for Fokker-Planck equation.
However, the references of the numerical and analysis methods for space-time fractional derivative for partial differential equations are less limited then numerical and analysis methods for the partial differential equations with only fractional derivative. Shen et al. [18] presented implicit and explicit difference approximation for the space-time Riesz-Caputo fractional derivative for diffusion-advection equation.
Hejazi et al. [27] presented a finite volume method to resolve the time-space on two sides’ fractional derivative for dispersion-advection equation.
In this paper, we study the existence of the unique solution of a partial differential equation class involving fractional derivative theoretically and numerically by the finite element method with three different time discretization methods: Hadamard-Diethlem quadrature formula, Riemann-Liouville and Caputo link and Grunwald-Letnikov approximation. A numerical example is given finally to verify the consistence of the theoretical and numerical results as conclusion.
Here are some fractional integrals and fractional derivatives order definitions.
Definition 2.1 For any positive integer n and s∈]$n-1, n[$, the fractional integration is defined by
$\mathrm{I}_{\mathrm{t}}^s \mathrm{v}(\mathrm{t})=\frac{1}{\Gamma(\mathrm{s})} \int_0^{\mathrm{t}}(\mathrm{t}-\tau)^{s-1} \mathrm{v}(\tau) \mathrm{d} \tau$, for all $\mathrm{t} \in[0, \mathrm{~T}]$
where, $\Gamma(\mathrm{z})=\int_0^{+\infty} \mathrm{x}^{\mathrm{z}-1} \mathrm{e}^{-\mathrm{x}}$ is Euler’s gamma function and n is the integer part of s.
Definition 2.2 For any positive integer n and s∈]$n-1, n[$, the left fractional derivative is defined by:
${ }_0^R D_t^s v(t)=\frac{1}{\Gamma(n-s)} \frac{d^n}{d t^n} \int_0^t(t-\tau)^{n-s-1} v(\tau) d \tau$
and the right fractional derivative is defined by
${ }_t^R D_t^s v(t)=\frac{(-1)^n}{\Gamma(n-s)} \frac{d^n}{d t^n} \int_t^T(t-\tau)^{n-s-1} v(\tau) d \tau$
for all t∈[0,T].
Definition 2.3 For any positive integer n and s∈]$n-1, n[$, the Caputo fractional derivative is defined by
${ }^c D_0^s v(t)=\frac{1}{\Gamma(n-s)} \frac{d}{d t} \int_0^t(t-\tau)^{n-\alpha-1} v(\tau) d \tau$
for all t∈[0,T].
Definition 2.4 For any positive integer n and s∈]$n-1, n[$, the Grunwald-Letnikov fractional derivative of a function u is defined by
${ }_0^{G L} D_t^S u(t)=\lim _{\Delta t \rightarrow 0} \frac{1}{\Delta t^s} \sum_{k=0}^n \delta_k^s u\left(t_{n-k}\right)$
where, $\delta_k^s=\frac{(-1)^k \Gamma(s+1)}{\Gamma(s+1-k) \Gamma(k+1)}$.
We introduce now some notations and definitions some functional spaces equipped with their norms, semi-norms and inner products, which are used hereafter.
Let Λ be a domain which may stand for I, Ω and $\mathcal{O}$. $L^2(\Lambda)$ is the space of measurable function whose square is Lebesgue integral in Λ, his inner product and norm are defined by:
$(u, v)_{\Lambda}=\int_{\Lambda} u v d \Lambda,\|u\|_{0, \Lambda}=(u, u)_{\Lambda}^{\frac{1}{2}}$
for all $u, v \in L^2(\Lambda)$.
We define some functional spaces with their norms:
$H^1(\Lambda)=\left\{w \in L^2(\Lambda), \frac{d w}{d x} \in L^2(\Lambda)\right\}$, $H_0^1(\Lambda)=\left\{w \in L^2(\Lambda),\left.w\right|_{\partial \Lambda}\right\}$, $H^m(\Lambda)=\left\{w \in L^2(\Lambda), \frac{d^k w}{d x^k} \in L^2(\Lambda), \forall k \leq m\right\}$,
where, the inner product and the corresponding norm of $H_0^1(\Lambda)$ is respectively define by:
$(u, w)_1=(u, w)+\left(\frac{d u}{d x}, \frac{d w}{d x}\right),\|w\|_1=(w, w)_1^{1 / 2}$.
We use the standard norm of $L^2(\Lambda)$ with the norm
$\|w\|_1=\left(\|w\|^2+d\left\|\frac{d w}{d x}\right\|^2\right)^{1 / 2}$.
Let $H^m(\Lambda)$ and $H_0^m(\Lambda)$ be the usual Sobolev spaces $W^{m, 2}(\Lambda)$ with usual norms denoted by $\|.\|_{m \cdot \Lambda}:$
$\|w\|_m=\left(\sum_{k=0}^m\left\|\frac{d^k w}{d x^k}\right\|_0^2\right)^{1 / 2}$.
Let $C_0^{\infty}(\Lambda)$ be the smooth functions space with compact support in Λ.
We also need to give some definition of some Sobolev spaces: For a non-negative real number s for the Sobolev space X with the norm$\|.\|_X$, let:
$H^s(I, X)=\left\{v ;\|v(., t)\|_X \in H^s(I)\right\}$
equipped with the norm
$\|v\|_{H^s(I, X)}=\|\| v(., t)\left\|_X\right\|_{s, I}$.
In particular, when X is $H^\sigma(\Omega)$ or $H_0^\sigma(\Omega), \sigma \geq 0$ the norm of the space HS (I, X) will be denoted by$\|.\|_{\sigma, s, O}$.
In the rest of this paper, while no confusion would arise about the domain symbols Ω, I or O is omitted from the notation.
Let Q=(a,b) which may stand for I or Ω, we define the spaces:
Definition 2.5 For any real s≥ 0, define the space${ }^l H_0^s(Q)$ as the closure of the $C_0^{\infty}(Q)$ with respect to the norm $\|.\|_{} l_{H_0^s(Q)}$, that is:
${ }^l H^s(Q)=\left\{v ;\|v\|_{l_{H^s}(Q)}<\infty\right\}$
with the norm
$\|v\|_{l_{H^s}(Q)}=\left(\|v\|_{0, Q}^2+|v|^2 l_{H^s(Q)}\right)^{\frac{1}{2}}$
and the semi-norm
$|v|_{l_{H^s}(Q)}=\|\|^{R L} D_Z^S v \|_{0, Q}$.
Definition 2.6 For any real s≥0, define the space r $H_0^s(Q)$ as the closure of the $C_0^{\infty}(Q)$ with respect to the norm $\|.\| r_{H_0^S}(Q)$, that is: ${r}_{H^s}(Q)=\left\{v ;\|v\|_{r_{H}{ }^s(Q)}<\infty\right\}$
with the norm
$\|v\|_{{r}_{H^s}(Q)}=\left(\|v\|_{0, Q}^2+|v|_{{r}_{H^s}(Q)}^2\right)^{\frac{1}{2}}$
and the semi-norm
$|v|_{rH^s(Q)}=\left\|_z^{R L} D^s v\right\|_{0, Q}$.
In the above notation “l” and “r” are used to indicate the left and right fractional spaces and its corresponding norm and semi-norm.
Definition 2.7 In the usual Sobolev space$H_0^S(Q)$, we also define the semi-norm
$|v|_{H_0^S(Q)}^*=\left(\frac{\left({ }^{R L} D_z^s v,{ }_z^{R L} D^S v\right)_Q}{\cos (\pi s)}\right)^{\frac{1}{2}}$
for all $v \in H_0^s(Q)$.
Lemma 2.1 For $s>0, s \neq n-\frac{1}{2}$, the spaces ${ }^l H^s(Q),{ }^r H^s(Q)$ and $H_0^s(Q)$ are equal in the sense their semi-norms are all equivalent to $|v|_{H_0^S(Q)}^*$.
Lemma 2.2 For $0<s<2, s \neq 1, w \in H^{\frac{s}{2}}(Q)$, we have:
${ }^{R L} D_z^s w={ }^{R L} D_z^{\frac{s}{2} R L} D_z^{\frac{s}{2}} w$.
Lemma 2.3 For $0<s<2, s \neq 1, w, v \in H^{\frac{s}{2}}(Q)$, we have:
$\left({ }^{R L} D_z^s w(z), v(z)\right)_Q=\left({ }^{R L} D_z^{\frac{s}{2}} w(z),{ }_z^{R L} D^{\frac{s}{2}} v(z)\right)_Q$,$\left({ }_z^{R L} D_{:}^S w(z), v(z)\right)_Q=\left({ }_z^{R L} D^{\frac{s}{2}} w(z),{ }^{R L} D_z^{\frac{s}{2}} v(z)\right)_Q$.
We define the space
$B^{s, \sigma}(\mathcal{O})=H^s\left(I, L^2(\Omega)\right) \cap L^2\left(I, H_0^\sigma(\Omega)\right)$
with the norm
$\|v\|_{B^{s, \sigma}}=\left(\|v\|_{H^s\left(I ; L^2(\Omega)\right)}^2+\|v\|_{L^2\left(I, H_0^\sigma(\Omega)\right)}^2\right)^{\frac{1}{2}}$
where, $\mathcal{O}=\Omega \times I$.
The weak formulation of our problem (1)-(3) is:
for $f \in L^2(\Omega)$ find $u \in B^{\frac{\alpha}{2} \frac{\beta}{2}}(\mathcal{O})$ such that:
$\mathcal{A}(u, v)=\mathcal{F}(v)$ (4)
for all $v \in B^{\alpha, \beta}(\mathcal{O})$ where $\mathcal{A}(.,$. is the bilinear form given by:
$\begin{aligned} & \mathcal{A}(u, v)=\left({ }_0^{R L} D_t^\alpha u, v\right)_{L^2(\mathcal{O})}+\left({ }_0^{R L} D_x^\beta u, v\right)_{L^2(0)^{\prime}} \\ = & \left({ }_0^{R L} D_t^{\frac{\alpha}{2}} u,{ }_t^{R L} D_T^{\frac{\alpha}{2}} v\right)_{L^2(0)}+\left({ }_0^{R L} D_x^{\frac{\beta}{2}} u,{ }_0^{R L} D_x^{\frac{\beta}{2}} v\right)_{L^2(0)}\end{aligned}$
and the functional $\mathcal{F}($. is given by:
$\mathcal{F}(v)=(f, v)_{L^2(o)}$
Theorem 3.1 For 0<α<1 and $f \in L^2(\mathcal{O})$, the problem (1)-(3) has a unique solution. Furthermore, we have the following stability result:
$\|u\|_B \frac{\alpha}{2}, \frac{\beta}{2}(\mathcal{O}) \leq C\|f\|_{L^2(\mathcal{O})}$. (5)
Proof: The existence of the unique solution is done by the well-known Lax-Milgram Lemma. It consists to prove the cœrcivity and continuity of the bilinear for $\mathcal{A}$ and the continuity of the linear functional $\mathcal{F}$ which is easy to prove.
1. The cœrcivity: Using above Lemmas, $\forall v \in B^{\frac{\alpha}{2}, \frac{\beta}{2}}(\mathcal{O})$ we have:
$\begin{aligned} & \mathcal{A}(v, v)=\left({ }_0^{R L} D_t^{\frac{\alpha}{2}} v,{ }_t^{R L} D_T^{\frac{\alpha}{2}} v\right)_{L^2(\mathcal{O})}+\left({ }_0^{R L} D_x^{\frac{\beta}{2}} v,{ }_0^{R L} D_x^{\frac{\beta}{2}} v\right)_{L^2(\mathcal{O})} \\ &=\cos \left(\frac{\pi \alpha}{2}\right) \\ &+\cos \left(\frac{\pi \beta}{2}\right)\left({ }_0^{R L} D_x^{\frac{\beta}{2}} v,{ }_0^{R L} D_x^{\frac{\beta}{2}} v\right)_{L^2(\mathcal{O})} \\ & \geq C\|v\|_{B^{\frac{\alpha}{2}}, \frac{\beta}{2}(\mathcal{O})} .\end{aligned}$
2. The continuity: $\forall v \in B^{\frac{\alpha}{2}, \frac{\beta}{2}}(\mathcal{O})$ we have:
$|\mathcal{A}(u, v)| \leq\left\|{ }_0^{R L} D_t^{\frac{\alpha}{2}} u\right\|_{L^2(\mathcal{O})}\left\|{ }_t^{R L} D_T^{\frac{\alpha}{2}} v\right\|_{L^2(\mathcal{O})}$$+\left\|{ }_0^{R L} D_x^{\frac{\beta}{2}} u\right\|_{L^2(O)}\left\|{ }_0^{R L} D_x^{\frac{\beta}{2}} v\right\|_{L^2(\mathcal{O})}$$\leq\|u\|_{H^{\frac{\alpha}{2}}\left(I ; L^2(\Omega)\right)}\|v\|_{H^{\frac{\alpha}{2}}\left(I ; L^2(\Omega)\right)}$$+\|u\|_{H^{\frac{\beta}{2}}\left(I ; L^2(\Omega)\right)}\|v\|_{H^{\frac{\beta}{2}}\left(I ; L^2(\Omega)\right)}$$\leq\|u\|_{B^{\frac{\alpha}{2}, \frac{\beta}{2}}(\mathcal{O})}\|v\|_{B^{\frac{\alpha}{2}}, \frac{\beta}{2}(\mathcal{O})}$.
By the Lax-Milgram Lemma, the problem (1)-(3) has a unique solution in $B^{\frac{\alpha}{2}, \frac{\beta}{2}}(\mathcal{O})$.
To prove the stability, we take v=u in (4) and using the cœrcivity result we get:
$\|u\|_{B^{\frac{\alpha}{2}}, \frac{\beta}{2}(O)} \leq C\|f\|_{L^2(\mathcal{O})}$.
The proof is complete.
4.1 Nodal base functions and their fractional derivatives properties
Let Ω=[0,1] be a finite domain. Let Ωh be a uniform partition of Ω, which is given by:
$0=x_0<x_1<\cdots<x_{m-1}<x_m=1$
where, m is a positive integer. $h=\frac{1}{m}=x_i-x_{i-1}$ and each $\Omega_{\mathrm{i}}=\left[x_{i-1}-x_i\right]$ for $i=1,2, \ldots, m$.
We also define the space Sh as the set of piecewise-linear polynomial define on Ωh.
$S_h=\left\{v ;\left.v\right|_{\Omega_{\mathrm{i}}} \in P_1\left(\Omega_{\mathrm{i}}\right),\{v\} \in C(\Omega)\right\}$,
where, $P_1\left(\Omega_i\right)$ is the space of linear polynomial define on Ωi.
The nodal base functions $\phi_i, i=0,1, \ldots, m$ of Sh are given as
$\phi_i= \begin{cases}\frac{x-x_i}{h}, & x \in\left[x_{i-1}, x_i\right], \\ \frac{x_{i+1}-x}{h}, & x \in\left[x_i, x_{i+1}\right], \\ 0, & \text { elsewhere }\end{cases}$
and
$\phi_0=\left\{\begin{array}{cc}\frac{x_1-x}{h}, & x \in\left[x_0, x_1\right], \\ 0, & \text { elsewhere },\end{array}\right.$$\phi_m=\left\{\begin{array}{cc}\frac{x-x_{m-1}}{h}, & x \in\left[x_m, x_{m-1}\right], \\ 0, & \text { elsewhere. }\end{array}\right.$
Lemma 4.1 For i=1, 2,…, m-1, we have:
$\left(\phi_i(x), \phi_j(x)\right)= \begin{cases}1, & |j-1|=1, \\ 4, & j=i, j=0,1,2, \ldots, m, \\ 0, & \text { otherwise. }\end{cases}$
Lemma 4.2 For i=1, 2,…, m-1 and 0<λ<1, we have:
${ }^{R L} D_x^\lambda \phi_i(x)$$=\mu \begin{cases}0, & a \leq x \leq x_{i-1} \\ \left(x-x_{i-1}\right)^{1-\lambda}, & x_{i-1} \leq x \leq x_i \\ \left(x-x_{i-1}\right)^{1-\lambda}-2\left(x-x_i\right)^{1-\lambda}, & x_i \leq x \leq x_{i+1}(x) \\ \left(x-x_{i-1}\right)^{1-\lambda}-2\left(x-x_i\right)^{1-\lambda}+\left(x-x_{i+1}\right)^{1-\lambda} \\ & x_{i+1} \leq x \leq b\end{cases}$
where, $\mu=\frac{1}{h \Gamma(2-\lambda)}$.
Lemma 4.3 For i=1, 2,…, m-1, we have:
$\int_{x_{i-1}}^{x_i}{ }^{R L} D_x^\lambda \phi_j(x) d x=\kappa \cdot p_{i, j}^{(1)}$$\int_{x_j}^{x_{i+1}}{ }_{R L} D_x^\lambda \phi_j(x) d x=\kappa \cdot p_{i, j}^{(2)}$
where,
$p_{i, j}^{(1)}=\kappa \begin{cases}c_{i-j}, & j \leq i-2, \\ 2^{2-\lambda}-3, & j=i-1, \\ 1, & j=i, \\ 0, & j \geq i+1\end{cases}$
and
$p_{i, j}^{(2)}=\kappa \begin{cases}c_{i-j+1}, & j \leq i-1, \\ 2^{2-\lambda}-3, & j=i, \\ 1, & j=i+1 \\ 0, & j \geq i+2\end{cases}$
where, $\kappa=\frac{h^{1-\lambda}}{\Gamma(3-\lambda)}$ and
$c_i=-(i-2)^{2-\lambda}+3(i-1)^{2-\lambda}-3 i^{2-\lambda}+(i+1)^{2-\lambda}$.
4.2 The General Time Discretization Scheme (GTDS)
We can rewrite the system (1)-(3) as:
${ }_0^{R L} D_t^\alpha\left[u(t)-u_0\right]+A u(t)=f(t)$ (6)
where, $A={ }_0 D_x^\beta$..
Consider that the fractional Riemann-Liouville derivative can be written in the form of a finite-part integral:
${ }_0^{R L} D_t^\alpha u(t)=\frac{1}{\Gamma(-\alpha)} \int_0^t \frac{u(\tau)}{(t-\tau)^{1+\alpha}} d \tau$
where, this integral can be evaluated as a Hadamard finite part integral.
Let us consider tj=j⁄n,j=0,…,n a partition of [0,1], then we have:
${ }_0^{R L} D_t^\alpha\left[u-u_0\right]\left(t_j\right)=\frac{t_j^{-\alpha}}{\Gamma(-\alpha)} \int_0^1 g(w) w^{-1-\alpha} d w$
where, $g(w)=u\left(t_j-t_j w\right)-u(0)$..
We use now the Diethelm's first degree compound quadrature formula with equidistant nodes 0,1⁄j,2⁄j,…,1 to obtain the following approximation of the fractional derivative:
$\int_0^1 g(w) w^{-1-\alpha} d w=\frac{t_j^{-\alpha}}{\Gamma(-\alpha)} \sum_{k=0}^j \alpha_{k j} g\left(\frac{k}{j}\right)+R_j(g)$ (7)
where, the remainder term Rj (g) satisfy:
$\left\|R_j(g)\right\| \leq C_i^{\alpha-2} \sup _{0 \leq t \leq T}\left\|g^{\prime \prime}(t)\right\|$
and the weights αkj (for j≥1) are given by:
$\alpha_{k j}=\frac{j^\alpha}{\alpha(1-\alpha)}\left\{\begin{array}{r}-1, for k=0, \\ 2 k^{1-\alpha}-(k-1)^{1-\alpha}-(k+1)^{1-\alpha}, \\ \text { for } k=1,2, \ldots, j-1, \\ (\alpha-1) k^{-\alpha}-(k-1)^{1-\alpha}+k^{1-\alpha} \\ \text { for } k=i\end{array}\right.$
So (7) can be reducing to:
${ }_0^{R L} D_t^\alpha\left[u-u_0\right]\left(t_j\right)$$=\frac{\Delta t^{-\alpha}}{\Gamma(2-\alpha)} \sum_{k=0}^j \beta_{k j}\left[u\left(t_j-t_k\right)\right.$$-u(0)]+\frac{t^{-\alpha}}{\Gamma(-\alpha)} R_j(g)$ (8)
where,
$\beta_{k j}=\left\{\begin{array}{rr}1, & \text { for } k=0, \\ -2 k^{1-\alpha}+(k-1)^{1-\alpha}+(k+1)^{1-\alpha}, \\ \text { for } k=1,2, \ldots, j-1, \\ -(\alpha-1) k^{-\alpha}+(k-1)^{1-\alpha}-k^{1-\alpha}, \\ & \text { for } k=i.\end{array}\right.$
for t=tj, (8) is given for j=0, …, n as:
$\frac{\Delta t^{-\alpha}}{\Gamma(2-\alpha)} \sum_{k=0}^j \beta_{k j}\left[u\left(t_j-t_k\right)-u(0)\right]+A u\left(t_j\right)=f\left(t_j\right)-\frac{t^{-\alpha}}{\Gamma(-\alpha)} R_j(g)$. (9)
Noting uj the approximation of u(tj ) and f(tj )=fj, so we have:
$\sum_{k=0}^j \beta_{k j}\left(u^{j-k}-u^0\right)+p A u^j=p f_j$ (10)
where, $p=\Delta t^\alpha \Gamma(2-\alpha)$.
The weak formulation of the problem (1)-(3) can also given by:
$\left({ }_0^{R L} D_t^\alpha u, v\right)+\left({ }_0^{R L} D_x^\lambda u, \frac{\partial v}{\partial x} v\right)=(f, v), 0<\lambda<1$. (11)
So the semi-discretization is:
$\left(u^j, v\right)+p\left({ }_0^{R L} D_x^\lambda u^j, \frac{\partial v}{\partial x}\right)=-\sum_{k=0}^j \beta_{k j}\left(u^{j-k}, v\right)+p\left(f_j, v\right)+\sum_{k=0}^j \beta_{k j}\left(u^0, v\right)$. (12)
To get the full discretization, we set$u^j=\sum_{s=0}^m u_s^j \phi_s(x)$, and choose every test function v to be $\phi_i, i=1, \ldots, m-1$ in (12), we get:
$\begin{aligned} & \sum_{s=0}^m u_s^j\left(\phi_s, \phi_i\right)+p \sum_{s=0}^m u_s^j\left({ }_0^{R L} D_x^\lambda \phi_s, \frac{\partial \phi_i}{\partial x}\right) \\ &=p\left(f_j, \phi_i\right)-\sum_{k=0}^{j-1} \beta_{k j} \sum_{s=0}^m u_s^{j-k}\left(\phi_s, \phi_i\right) \\ &+\sum_{k=0}^{j-1} \beta_{k j}\left(u^0, \phi_i\right)\end{aligned}$
for i=1,…,m-1.
Applying Lemma 4.1 and Lemma 4.3, we obtain:
$\begin{aligned} \frac{h}{6}\left(u_{i-1}^j+4 u_i^j\right. & \left.+u_{i+1}^j\right)+r\left[\sum_{s=0}^m\left(p_{s i}^{(1)}-p_{s i}^{(2)}\right) u_s^j\right] \\ & =-\frac{h}{6} \sum_{k=0}^{j-1} \beta_{k j}\left(u_{i-1}^j+4 u_i^j+u_{i+1}^j\right) \\ & +p\left(f_j, \phi_i\right)+\sum_{k=0}^{j-1} \beta_{k j}\left(u^0, \phi_i\right)\end{aligned}$ (13)
where, $r=\frac{p \kappa}{h}$.
4.3 Caputo's Method (CM)
We need the following lemma in order to give the second scheme.
Lemma 4.4 For 0<α<1, we have:
${ }_0^{R L} D_t^\alpha u(x, t)={ }^C D_t^\alpha u(x, t)+\frac{u(x, 0) t^{-\alpha}}{\Gamma(1-\alpha)}$.
This lemma gives the link between the Riemann-Liouville fractional derivative and the Caputo fractional derivative.
So the problem (1)-(3) can be written as:
${ }_0^C D_t^\alpha u(x, t)-{ }_0^{R L} D_x^\beta u(x, t)=f(x, t)-\frac{u_0(x) t^{-\alpha}}{\Gamma(1-\alpha)}$. (14)
The time-fractional derivative is estimated [25, 28]:
${ }_0^C D_t^\alpha u\left(x, t_n\right)$$=\frac{\Delta t^{1-\alpha}}{\Gamma(1-\alpha)} \sum_{k=0}^{n-1} b_k \frac{u\left(x, t_{n-k}\right)-u\left(x, t_{n-k-1}\right)}{\Delta t}$ (15)
where, $b_k=(k+1)^{1-\alpha}-k^{1-\alpha}$.
The weak form of problem (1)-(3) is given by:
$\left({ }_0^c D_t^\alpha u, v\right)+\left({ }_0^{R L} D_x^\lambda u, \frac{\partial v}{\partial x^\lambda}\right)=(f, v)-\frac{t^{-\alpha}\left(u_0, v\right)}{\Gamma(1-\alpha)}$ (16)
for all $v \in H_0^1(\Omega)$ where 0<λ<1.
Denoting u^nthe approximation solution of u(x,tn ) and fn (x)=f(x,tn). Using (15), we can write the weak form (16) at tnas:
$\begin{aligned}\left(u^n, v\right) & +p\left({ }_0^{R L} D_x^\lambda u^n, \frac{\partial v}{\partial x}\right) \\ & =\sum_{k=1}^{n-1} w_k\left(u^{n-k}, v\right)\left(b_{n-1}\right. \\ & \left.+(\alpha-1) n^{-\alpha}\right)\left(u_0, v\right)+p\left(f_n, v\right)\end{aligned}$ (17)
where, $p=\Gamma(2-\alpha) \Delta t^\alpha$ and $w_k=b_{k-1}-b_k$.
Let $u^n=\sum_{j=0}^m u_j^n \phi_j(x)$, choosing v to be $\phi_i(x)$ for i=1,…,m-1 and applying Lemma 4.1 and Lemma 4.3, we get the full discretization of the problem:
$\frac{h}{6}\left(u_{i-1}^n+4 u_i^n+u_{i+1}^n\right)+r\left[\sum_{j=0}^m\left(p_{i j}^{(1)}-p_{i j}^{(2)}\right) u_j^n\right]$$=\frac{h}{6} \sum_{k=1}^{j-1} w_k\left(u_{i-1}^n+4 u_i^n+u_{i+1}^n\right)+p\left(f_n, \phi_i\right)$$+\left(b_{n-1}+(\alpha-1) n^{-\alpha}\right)\left(u^0, \phi_i\right)$. (18)
4.4 The Grunwald-Letnikov Method (GLM)
The Grunwald-Letnikov fractional derivative can be approximated the Riemman-Liouville fractional derivative [29, 30], which is giving by the following formula:
$\left.{ }_0^{R L} D_t^\alpha u(t)\right|_{t=t_n} \approx \Delta t^{-\alpha} \sum_{k=0}^n \delta_k^{(\alpha)} u\left(t_{n-k}\right)$
where, $\delta_k^{(\alpha)}=\frac{(-1)^k \Gamma(\alpha+1)}{\Gamma(\alpha+1-k) \Gamma(k+1)}$..
Consider the following semi discretization of the problem (13):
$\Delta t^{-\alpha} \sum_{k=0}^n \delta_k^{(\alpha)}\left(u^{n-k}, v\right)+\left({ }_0^{R L} D_x^\lambda u^n, \frac{\partial v}{\partial x}\right)=(f n, v)$ (19)
Let 〖$\varepsilon_k^{(\alpha)}=-\delta_k^{(\alpha)}$, then we have:
$\left(u^n, v\right)+\Delta t^\alpha\left({ }^{R L} D_x^\lambda u^n, \frac{\partial v}{\partial x}\right)=\sum_{k=1}^{n-1} \varepsilon_k^{(\alpha)}\left(u^{n-k}, v\right)$$+\Delta t^\alpha\left(f_n, v\right)+\varepsilon_k^{(\alpha)}\left(u^0, v\right)$. (20)
To get the full discretization, we set $u_j=\sum_{j=0}^m u_j^m \phi_j(x)$ and we choose v to be ϕi (x) for i=1,…,m-1, then we have
$\begin{aligned} \sum_{j=0}^m u_j^n\left(\phi_j, \phi_i\right) & +\Delta t^\alpha \sum_{j=0}^m u_j^m\left({ }_0^{R L} D_x^\lambda \phi_j, \frac{\partial \phi_i}{\partial x}\right) \\ = & \sum_{k=1}^{n-1} \varepsilon_k^{(\alpha)} \sum_{j=0}^m u_j^{n-k}\left(\phi_j, \phi_i\right) \\ & +\Delta t^\alpha\left(f_n, \phi_i\right)+\varepsilon_n^{(\alpha)}\left(u^0, \phi_i\right)\end{aligned}$. (21)
Applying Lemma 4.1 and Lemma 4.3, we obtain:
$\begin{aligned} \frac{h}{6}\left(u_{i-1}^n\right. & \left.+4 u_i^n+u_{i+1}^n\right)+r_1 \sum_{j=0}^m\left(p_{i j}^{(1)}-p_{i j}^{(2)}\right) u_j^n \\ = & \sum_{k=1}^{n-1} \varepsilon_k^{(\alpha)}\left(u_{i-1}^{n-k}+4 u_i^{n-k}+u_{i+1}^{n-k}\right) \\ & +\Delta t^\alpha\left(f_n, \phi_i\right)+\varepsilon_n^{(\alpha)}\left(u^0, \phi_i\right)\end{aligned}$ (22)
where, $r_1=\frac{\Delta t^\alpha \kappa}{h}$.
In order to verify the efficiency of the three numerical methods, we present the next example. The L2 norm is used to investigate the error. All the numerical result are evaluated at T=1.
${ }_0^{R L} D_t^a u(x, t)-{ }_0^{R L} D_x^{\lambda+1} u(x, t)=f(x, t), \quad(x, t) \in[0,1]^2$
subject to the initial and boundary conditions:
$\begin{array}{cc}u(x, 0)=0, & x \in[0,1] \\ u(0, t)=u(1, t)=0, & t \in[0, T].\end{array}$
where, $f(x, t)=\frac{t^{1-\alpha}}{\Gamma(2-\alpha)} x^3+\frac{6 t}{\Gamma(3-\alpha)} x^{2-\lambda}$.
The exact solution is u(x,t)=tx3.
Table 1 shows exclusively the error between the exact and numerical solution using the GTDS and CM, and the Table 2 for the GLM, for a fixed value of $\alpha=0.5 $ and a fixed step time $\Delta t=1 / 100$ , and different value of $\lambda$. For the Table 3 and Table 4, we fixed $\lambda=0.75$ and $h \lambda=1 / 1000$ and the result are gotten for a different value of $\alpha$.
Table 1. The error of GTDS and CM (×10-6) for α=0.5 and Δt=1/100
|
h |
λ=0.25 |
λ=0.5 |
λ=0.75 |
λ=0.95 |
|
1/8 |
0.08686 |
0.04526 |
0.0196 |
0.03116 |
|
1/16 |
0.01895 |
0.01029 |
0.4755 |
0.7853 |
|
1/32 |
0.4429 |
0.247 |
0.117 |
0.1971 |
|
1/64 |
0.1076 |
6.083 |
0.2904 |
4.9937 |
|
1/128 |
2.659 |
1.513 |
7.241 |
1.235 |
Table 2. The error of GLM (×10-5) for α=0.5 and Δt=1/100
|
h |
λ=0.5 |
λ=0.5 |
λ=0.75 |
λ=0.95 |
|
1/8 |
0.977 |
0.496 |
0.1814 |
0.2822 |
|
1/16 |
0.3067 |
0.1571 |
4.691 |
4.947 |
|
1/32 |
0.1783 |
9.457 |
4.077 |
95.91 |
|
1/64 |
0.153 |
8.373 |
4.418 |
2.417 |
|
1/128 |
0.147 |
8.15 |
4.529 |
2.789 |
Table 3. The error of GTDS and CM (×10-8) for λ=0.75 and h=1/10000
|
Δt |
α=0.25 |
α=0.5 |
α=0.75 |
α=0.95 |
|
1/20 |
1.11759 |
1.18435 |
1.21105 |
1.16581 |
|
1/40 |
1.11757 |
1.18456 |
1.214 |
1.16724 |
|
1/80 |
1.11756 |
1.18462 |
1.21177 |
1.16782 |
|
1/160 |
1.11754 |
1.18468 |
1.21179 |
1.16805 |
|
1/320 |
1.11761 |
1.18467 |
1.1182 |
1.16815 |
Table 4. The error of GLM (×10-5) for λ=0.75 and h=1/10000
|
Δt |
α=0.25 |
α=0.5 |
α=0.75 |
α=0.95 |
|
1/20 |
0.1507 |
0.2257 |
0.1886 |
5.046 |
|
1/40 |
7.57 |
0.1136 |
9.511 |
2.521 |
|
1/80 |
3.794 |
5.703 |
4.79 |
1.264 |
|
1/160 |
1.9 |
2.86 |
2.409 |
63.52 |
|
1/320 |
0.9865 |
0.01404 |
0.01037 |
0.2398 |
In this paper, we studied the finite element method for solving a class of Riemann-Liouville space-time fractional partial differential equations. In the first step we approximated the time fractional derivative using the Hadamard finite part integral and the Diethlem's first degree compound quadrature formula; the second approach was based on the link between Riemann-Liouville and Caputo fractional derivative, when the third method was based on the approximation of the Riemann-Liouville by the Grunwald-Letnikov fractional derivative. For the approximation of the space fractional derivative the finite element method was introduced for all the three approaches. Finally to check the effectiveness of the three methods a numerical example was given.
The authors thank the General Directorate for Scientific Research and Technological Development (DGRSDT). The authors thank also the anonymous reviewers for their helpful remarks.
[1] Chaurasia, V.B.L., Singh, J. (2012). Solution of a time-space fractional diffusion equation by integral transform method. Tamsui Oxford Journal of Information and Mathematical Sciences, 28(2): 153-164.
[2] Povstenko, Y.Z. (2014). Fundamental solutions to time-fractional advection diffusion equation in a case of two space variables. Mathematical Problems in Engineering, 2014: Article ID: 705364. https://doi.org/10.1155/2014/705364
[3] Povstenko, Y., Klekot, J. (2016). The Dirichlet problem for the time-fractional advection-diffusion equation in a line segment. Boundary Value Problems, 2016(89): 1-8. https://doi.org/10.1186/s13661-016-0597-4
[4] Mainardi, F., Pagnini, G. (2002). Space-time fractional diffusion: exact solutions and probability interpretation. Waves and Stability in Continuous Media, 296-301. https://doi.org/10.1142/9789812777331_0037
[5] Salim, T.O., El-Kahlout, A. (2009). Analytical solution of time-fractional advection dispersion equation. Applications and Applied Mathematics, 4(1): 176-188.
[6] Stern, R., Effenberger, F., Fichtner, H., Schäfer, T. (2013). The space-fractional diffusion-advection equation: analytical solutions and critical assessment of numerical solutions. Fractional Calculus and Applied Analysis, 17(1): 171-190. https://doi.org/10.2478/s13540-014-0161-9
[7] Zhuang, P., Liu, F., Anh, V., Turner, I. (2008). New solution and analytical techniques of the implicit numerical method for the anomalous subdiffusion equation. SIAM Journal on Numerical Analysis, 46(2): 1079-1095. https://doi.org/10.1137/060673114
[8] Chen, J., Liu, F., Burrage, K., Shen, S. (2013). Numerical techniques for simulating a fractional mathematical model of epidermal wound healing. Journal of Applied Mathematics and Computing, 41(1): 33-47. https://doi.org/10.1007/s12190-012-0591-7
[9] Diethelm, K. (1997). Generalized compound quadrature formulae for finite-part integrals. IMA Journal of Numerical Analysis, 17(3): 479-493.https://doi.org/10.1093/imanum/17.3.479
[10] Feng, L.B., Zhuang, P., Liu, F., Turner, I., Gu, Y.T. (2016). Finite element method for space-time fractional diffusion equation. Numerical Algorithms, 72(3): 749-767. https://doi.org/10.1007/s11075-015-0065-8
[11] Ford, N.J., Xiao, J., Yan, Y. (2011). A finite element method for time fractional partial differential equations. Fractional Calculus and Applied Analysis, 14(3): 454-474. https://doi.org/10.2478/s13540-011-0028-2
[12] Ford, N., Xiao, J., Yan, Y. (2011). A finite element method for time fractional partial differential equations. Fractional Calculus and Applied Analysis, 14(3): 454-474. https://doi.org/10.2478/s13540-011-0028-2
[13] Hou, D., Hasan, M.T., Xu, C. (2018). Müntz spectral methods for the time-fractional diffusion equation. Computational Methods in Applied Mathematics, 18(1): 43-62. https://doi.org/10.1515/cmam-2017-0027
[14] Li, X., Xu, C. (2009). A space-time spectral method for the time fractional diffusion equation. SIAM Journal on Numerical Analysis, 47(3): 2108-2131. https://doi.org/10.1137/080718942
[15] Li, X., Xu, C. (2010) Existence and uniqueness of the weak solution of the space-time fractional diffusion equation and a spectral method approximation. Communications in Computational Physics, 8(5): 1016-1051. https://doi.org/10.4208/cicp.020709.221209a
[16] Li, C., Chen, A. (2017). Numerical methods for fractional partial differential equations. International Journal of Computer Mathematics, 95(2): 1-60. https://doi.org/10.1080/00207160.2017.1343941
[17] Li, B., Luo, H., Xie, X. (2018). A time-spectral algorithm for fractional wave problems. Journal of Scientific Computing, 77(2): 1164-1184. https://doi.org/10.48550/arXiv.1708.02720
[18] Shen, S., Liu, F., Anh, V. (2011). Numerical approximations and solution techniques for the space-time Riesz-Caputo fractional advection-diffusion equation. Numerical Algorithms, 56(3): 383-403. https://doi.org/10.1007/s11075-010-9393-x
[19] Zeng, F., Li, C., Liu, F., Turner, L. (2013). The use of finite difference/element approaches for solving the time-fractional subdiffusion equation. SIAM Journal on Scientific Computing, 35(6): A2976-A3000. https://doi.org/10.1137/130910865
[20] Zhao, J., Xiao, J., Xu, Y. (2013). Stability and convergence of an effective finite element method for multiterm fractional partial differential equations. Abstract and Applied Analysis, 2013(3):https://doi.org/10.1155/2013/857205
[21] Zheng, M., Liu, F., Turner, I., Anh, V. (2015). A novel high order space-time spectral method for the time fractional Fokker-Planck equation. SIAM Journal on Scientific Computing, 37(2): A701-A724. https://doi.org/10.1137/140980545
[22] Edeki, S.O., Akinlabi, G.O., Adeosun, S.A. (2016). Analytic and numerical solutions of time-fractional linear Schrödinger equation. Communications in Mathematics and Applications, 7(1): 1-10.https://doi.org/10.26713/cma.v7i1.327
[23] Ali Shah, N., Saleem, S., Chung, J.D. (2021). Numerical analysis of time-fractional diffusion equations via a novel approach. Journal of Function Spaces, 2021:Article ID: 9945364. https://doi.org/10.1155/2021/9945364
[24] Meerschaert, M.M., Tadjeran, C. (2004). Finite difference approximations for fractional advection-dispersion flow equations. Journal of Computational and Applied Mathematics, 172: 65-77. https://doi.org/10.1016/j.cam.2004.01.033
[25] Li, Y., Xu, C. (2007). Finite difference/spectral approximations for the time-fractional diffusion equation. Journal of Computational Physics, 225(2): 1533-1552.https://doi.org/10.1016/j.jcp.2007.02.00111
[26] Zheng, Y., Li, C., Zhao, Z. (2010). A note on the finite element method for the space-fractional advection diffusion equation. Computers & Mathematics with Applications, 59(5): 1718-1726. https://doi.org/10.1016/j.camwa.2009.08.071
[27] Hejazi, H., Moroney, T., Liu, F. (2013). A finite volume method for solving the two-sided time-space fractional advection-dispersion equation. Central European Journal of Physics, 11(10): 1275-1283.https://doi.org/10.2478/s11534-013-0317-y
[28] Jiang, Y., Ma, J. (2011). High-order finite element methods for time-fractional differential equations. Journal of Computational and Applied Mathematics, 235(11): 3285-3290. https://doi.org/10.1016/j.cam.2011.01.011
[29] Podlubny, I. (1999). Fractional differential equations. Academic Press, New York.
[30] Li, C., Chen, A. (2018). Numerical methods for fractional partial differential equations. International Journal of Computer Mathematics, 95(6-7): 1048-1099. https://doi.org/10.1080/00207160.2017.1343941