Mathematical analysis and symmetric fractional-order reduction method for diffusion-wave equations

1. Problem Statement

The paper addresses the numerical solution of time-fractional diffusion-wave equations (also known as fractional wave equations) of order $\alpha \in (1, 2)$:
$$\partial_t^\alpha u - u_{xx} = g(x, t), \quad (x, t) \in \Omega \times (0, T)$$
with initial conditions $u(x,0)=a_0(x)$, $\partial_t u(x,0)=a_1(x)$, and homogeneous Dirichlet boundary conditions.

Key Challenges:

• Low Regularity: The solutions often exhibit weak singularities near $t=0$, particularly when the initial velocity $a_1(x)$ lacks smoothness (e.g., $a_1(x) \notin C^2(\Omega)$) or when the source term $g(x,t)$ is singular.
• Limitations of Standard Reduction: Traditional order reduction methods (introducing $v = u_t$) often require $a_1(x)$ to be sufficiently smooth (e.g., $C^2$) to define terms like $\Delta a_1(x)$, which restricts their applicability to problems with discontinuous or low-regularity initial data.
• Uniform Meshes: Standard numerical schemes on uniform temporal grids fail to achieve optimal convergence rates for solutions with singularities at $t=0$.

2. Methodology

The authors propose a Symmetric Fractional-Order Reduction (SFOR) method combined with nonuniform temporal meshes (graded meshes) and L-type discretization schemes.

A. Symmetric Fractional-Order Reduction (SFOR)

Instead of the standard reduction $v = \partial_t u$, the authors utilize a decomposition based on Lemma 1.1:
$$\partial_t^\alpha u(t) = \partial_t^{\alpha/2} (\partial_t^{\alpha/2} u(t)) - u'(0)\omega_{2-\alpha}(t)$$
where $\omega_\beta(t) = t^{\beta-1}/\Gamma(\beta)$. This leads to two equivalent reformulated systems:

1. Model (1.4): Introduces an auxiliary variable $v = \partial_t^{\alpha/2} u$.
   $$\partial_t^{\beta} v - u_{xx} = a_1(x)\omega_{2-\alpha}(t) + g(x, t)$$
   Note: This model retains the singular source term $a_1(x)\omega_{2-\alpha}(t)$ but avoids requiring $\Delta a_1(x)$.

2. Model (1.5): Introduces a modified auxiliary variable $z = \partial_t^{\beta} u - [a_1(x) + \Gamma(2-\alpha)f(x)]\omega_{2-\alpha/2}(t)$.
   $$\partial_t^{\beta} z - u_{xx} = 0$$
   Note: This model absorbs the singularity into the definition of $z$, resulting in a homogeneous equation, provided $\nabla a_1(x)$ exists.

These formulations allow the method to handle cases where $a_1(x)$ has weaker regularity (e.g., $a_1 \in L^2$ or $H^1$) compared to previous methods requiring $a_1 \in C^2$.

B. Numerical Discretization

The authors apply two classical schemes on nonuniform graded meshes defined by $t_k = T(k/N)^r$ ($r \ge 1$):

• L1 Scheme: First-order approximation ($\theta = 0$).
• L2-1$\sigma$ Scheme: Second-order approximation ($\theta = \beta/2 = \alpha/4$).

The discrete Caputo derivative is approximated using convolution sums with discrete kernels $A^{(n)}_{n-k}$. The paper rigorously analyzes the positivity and monotonicity of these kernels on nonuniform grids to ensure stability.

3. Key Contributions

1. Extension of SFOR to Low Regularity: The paper extends the Symmetric Fractional-Order Reduction method to cases where the initial velocity $a_1(x)$ is not in $C^2(\Omega)$. It demonstrates that Models (1.4) and (1.5) are numerically preferable to the standard reduction (Model 1.6) when $a_1(x)$ is discontinuous or has limited smoothness, as they do not require the computation of $\Delta a_1(x)$.
2. Optimal Mesh Grading Parameters: The authors derive optimal mesh grading parameters ($r$) tailored to the specific regularity of the solution and the chosen scheme:
   • For L1 scheme on Model (1.4): $r = \frac{4-\alpha}{2-\alpha}$.
   • For L1 scheme on Model (1.5): $r = \max\{\frac{4-\alpha}{\alpha}, 2\}$.
   • For L2-1$\sigma$ scheme on Model (1.4): $r = \frac{4}{2-\alpha}$.
   • For L2-1$\sigma$ scheme on Model (1.5): $r = \max\{\frac{4}{\alpha}, \frac{8}{4-\alpha}\}$.
3. Rigorous Regularity and Stability Analysis:
   • Proved stability estimates for Models (1.4), (1.5), and (1.6) under different regularity assumptions for $a_1(x)$.
   • Established detailed regularity theory for the solution $u$ and auxiliary variables $v, z$ in Sobolev spaces $D((-\Delta)^\gamma)$, showing how the singularity strength depends on the smoothness of initial data.
4. Unconditional Convergence: Provided unconditional convergence proofs for both L1 and L2-1$\sigma$ schemes in $H^1(\Omega)$ and $L^2(\Omega)$ norms, even for discontinuous initial data.

4. Results

• Theoretical Convergence Rates:

  • The L1 scheme achieves an optimal convergence rate of $O(N^{-\min\{r(1-\beta), 2-\beta\}})$ in the $H^1$ norm. With the optimal $r$, this recovers the expected order (e.g., $2-\beta$).
  • The L2-1$\sigma$ scheme achieves an optimal convergence rate of $O(N^{-\min\{r(1-\beta), 2\}})$ (effectively second-order) in the $H^1$ norm when the mesh parameter $r$ is chosen correctly.
  • Similar optimal rates were derived for the $L^2$ norm when $a_1(x)$ is discontinuous.

• Numerical Experiments:

  • Example 1: Tested with smooth and non-smooth initial data ($a_0, a_1$).
    • Uniform meshes ($r=1$) resulted in convergence orders $< 1$.
    • Graded meshes with optimal $r$ restored the theoretical convergence orders (e.g., $\approx 1.375$ for L1 with $\alpha=1.25$, and $\approx 2$ for L2-1$\sigma$).
  • Example 2: Tested with discontinuous initial velocity $a_1(x) = \chi_{(0, \pi/2]}(x)$.
    • The L1 scheme on Model (1.4) successfully achieved the predicted convergence rates, validating the method's robustness for non-smooth data where standard reduction methods might fail or require regularization.
  • Tables 1-7: Comprehensive data confirms that the proposed schemes match the theoretical error bounds across various $\alpha \in (1, 2)$ and mesh parameters.

5. Significance

• Relaxed Regularity Assumptions: This work removes the restrictive requirement that initial velocity must be $C^2$ smooth, making the numerical treatment of fractional wave equations applicable to a broader class of physical problems involving shocks or discontinuous sources.
• Systematic Framework: It provides a systematic protocol for selecting the appropriate SFOR model (1.4 vs. 1.5) and mesh grading parameter $r$ based on the regularity of the input data.
• Computational Efficiency: By identifying the exact optimal mesh grading parameters, the method ensures that high-order accuracy is achieved without unnecessary computational cost (over-refining the mesh).
• Theoretical Foundation: The paper fills a gap in the literature by verifying the stability and convergence of SFOR methods on nonuniform grids for low-regularity solutions, addressing previous analytical limitations regarding the positivity of discrete coefficients.

In summary, the paper presents a robust, theoretically grounded numerical framework for solving fractional diffusion-wave equations with singular solutions, significantly advancing the state-of-the-art in handling low-regularity initial data through symmetric reduction and adaptive meshing.