Electrical Resistivity Tomography Investigation of Permafrost Conditions in a Thermokarst Site in Fairbanks, Alaska

Abdallah Basiru; Shishay T Kidanu; Sergei Rybakov; Nicholas Hasson; Moustapha Kebe; Emmanuel Osei Acheampong; Abdallah Basiru; Shishay T Kidanu; Sergei Rybakov; Nicholas Hasson; Moustapha Kebe; Emmanuel Osei Acheampong

doi:10.3934/geosci.2024001

AIMS Geosciences

2024, Volume 10, Issue 1: 1-27. doi: 10.3934/geosci.2024001

Previous Article Next Article

Research article Topical Sections

Electrical Resistivity Tomography Investigation of Permafrost Conditions in a Thermokarst Site in Fairbanks, Alaska

1.
Department of Civil, Geological, and Environmental Engineering, University of Alaska Fairbanks, 1760 Tanana Loop Duckering 245, USA
2.
Geophysical Institute, University of Alaska-Fairbanks, 903 N Koyukuk Drive, Fairbanks, AK, 99709, USA
3.
Department of Mining and Minerals Engineering, University of Alaska Fairbanks, 1760 Tanana Loop Duckering 245, USA
4.
Department of Geological Sciences, Ohio University, Clippinger Laboratories 316 Ohio University Athens, OH 45701, USA

Received: 12 September 2023 Revised: 23 November 2023 Accepted: 01 December 2023 Published: 09 January 2024

The degradation of permafrost poses severe environmental threats to communities in cold regions. As near-surface permafrost warms, extensive topographic variability is prevalent in the Arctic and Sub-Arctic communities. Geologic hazards such as thermokarst are formed due to varying rates of permafrost degradation, resulting in ground subsidence. This gradual subsidence or abrupt collapse of the earth causes a danger to existing infrastructure and the economic activities of communities in cold regions. Understanding the causes of thermokarst development and its dynamics requires imaging its underground morpho-structures and characterizing the surface and subsurface controls. In this study, we conducted a two-dimensional (2D) electrical resistivity tomography (ERT) survey to characterize the permafrost conditions in a thermokarst prone site located in Fairbanks, Alaska. To increase the reliability in the interpretability of the ERT data, borehole data and the depth-of-investigation (DOI) methods were applied. By using the 2D and three-dimensional (3D) ERT methods, we gained valuable information on the spatial variability of transient processes, such as the movement of freezing and thawing fronts. Resistivity imaging across the site exhibited distinct variations in permafrost conditions, with both low and high resistive anomalies observed along the transects. These anomalies, representing taliks and ice wedges, were characterized by resistivity values ranging from 50 Ωm and above 700 Ωm, respectively. The results from this study showed the effectiveness of ERT to characterize permafrost conditions and thermokarst subsurface morpho-structures. The insights gained from this research contribute to a better understanding of the causes and dynamics of thermokarst, which can be instrumental for engineers in developing feasible remedial measures.

Keywords:

Citation: Abdallah Basiru, Shishay T Kidanu, Sergei Rybakov, Nicholas Hasson, Moustapha Kebe, Emmanuel Osei Acheampong. Electrical Resistivity Tomography Investigation of Permafrost Conditions in a Thermokarst Site in Fairbanks, Alaska[J]. AIMS Geosciences, 2024, 10(1): 1-27. doi: 10.3934/geosci.2024001

Related Papers:

[1]	Muhammad Asim Khan, Norma Alias, Umair Ali . A new fourth-order grouping iterative method for the time fractional sub-diffusion equation having a weak singularity at initial time. AIMS Mathematics, 2023, 8(6): 13725-13746. doi: 10.3934/math.2023697
[2]	Mubashara Wali, Sadia Arshad, Sayed M Eldin, Imran Siddique . Numerical approximation of Atangana-Baleanu Caputo derivative for space-time fractional diffusion equations. AIMS Mathematics, 2023, 8(7): 15129-15147. doi: 10.3934/math.2023772
[3]	Junying Cao, Zhongqing Wang, Ziqiang Wang . Stability and convergence analysis for a uniform temporal high accuracy of the time-fractional diffusion equation with 1D and 2D spatial compact finite difference method. AIMS Mathematics, 2024, 9(6): 14697-14730. doi: 10.3934/math.2024715
[4]	Ajmal Ali, Tayyaba Akram, Azhar Iqbal, Poom Kumam, Thana Sutthibutpong . A numerical approach for 2D time-fractional diffusion damped wave model. AIMS Mathematics, 2023, 8(4): 8249-8273. doi: 10.3934/math.2023416
[5]	Zeshan Qiu . Fourth-order high-precision algorithms for one-sided tempered fractional diffusion equations. AIMS Mathematics, 2024, 9(10): 27102-27121. doi: 10.3934/math.20241318
[6]	Zhichao Fang, Ruixia Du, Hong Li, Yang Liu . A two-grid mixed finite volume element method for nonlinear time fractional reaction-diffusion equations. AIMS Mathematics, 2022, 7(2): 1941-1970. doi: 10.3934/math.2022112
[7]	Yanjie Zhou, Xianxiang Leng, Yuejie Li, Qiuxiang Deng, Zhendong Luo . A novel two-grid Crank-Nicolson mixed finite element method for nonlinear fourth-order sin-Gordon equation. AIMS Mathematics, 2024, 9(11): 31470-31494. doi: 10.3934/math.20241515
[8]	Krunal B. Kachhia, Jyotindra C. Prajapati . Generalized iterative method for the solution of linear and nonlinear fractional differential equations with composite fractional derivative operator. AIMS Mathematics, 2020, 5(4): 2888-2898. doi: 10.3934/math.2020186
[9]	Abdul-Majeed Ayebire, Saroj Sahani, Priyanka, Shelly Arora . Numerical study of soliton behavior of generalised Kuramoto-Sivashinsky type equations with Hermite splines. AIMS Mathematics, 2025, 10(2): 2098-2130. doi: 10.3934/math.2025099
[10]	Lei Ren . High order compact difference scheme for solving the time multi-term fractional sub-diffusion equations. AIMS Mathematics, 2022, 7(5): 9172-9188. doi: 10.3934/math.2022508

Abstract

1. Introduction

Fractional partial differential equations (FPDEs) have attracted considerable attention in various fields. Though research shows that many phenomena can be described by FPDEs such as physics ^[1], engineering ^[2], and other sciences ^[3,4]. However, finding the exact solutions of FPDEs by using current analytical methods such as Laplace transform, Green's function, and Fourier-Laplace transform (see ^[5,6] for examples) are often difficult to achieve^[7]. Thus, proposing numerical methods to find approximate solutions of these equations has practical importance. Due to this fact, in recent years a large number of numerical methods have been proposed for solving FPDEs, for instances see ^{[8,9,10,11,12]} and the references therein.

The time fractional diffusion-wave equation is obtained from the classical diffusion-wave equation by replacing the second order time derivative term with a fractional derivative of order $\alpha$ , $1 < \alpha < 2$ , and it can describe the intermediate process between parabolic diffusion equations and hyperbolic wave equations. Many of the universal mechanical, acoustic and electromagnetic responses can be accurately described by the time fractional diffusion-wave equation, see ^[13,14] for examples. The fourth order space derivative arises in the wave propagation in beams and modeling formation of grooves on a flat surface, thus considerable attention has been devoted to fourth order fractional diffusion-wave equation and its applications, see ^[15]. In this paper, the following nonlinear time fractional diffusion-wave equation with fourth order derivative in space and homogeneous initial boundary conditions will be considered

$\begin{align} \frac{\partial^2u(x, t)}{\partial t^2} + {}^{C}_{0}D_{t}^{\alpha}u(x, t) +K_{c}\frac{\partial^{4}u(x, t)}{\partial x^4} = \frac{\partial^2 u(x, t)}{\partial x^2} + g(u)+ f(x, t), \end{align}$

(1.1)

where $1 < \alpha < 2$ , $f(x, t)$ is a known function, $g(u)$ is a nonlinear function of $u$ with $g(0) = 0$ and satisfies the Lipschitz condition, and ${}^{C}_{0}D_{t}^{\alpha}u(x, t)$ denotes the temporal Caputo derivative with order $\alpha$ defined as

$\begin{align*} {}^{C}_{0}D_{t}^{\alpha}u(x, t) = \frac{1}{\Gamma(2-\alpha)}\int_{0}^{t}(t-s)^{1-\alpha}\frac{\partial^2u(x, s)}{\partial s^2}ds. \end{align*}$

Recently, there exist many works on numerical methods for time fractional diffusion-wave equations (TFDWEs), see ^{[16,17,18,19,20,21,22]} and the references therein. Chen et al. ^[17] proposed the method of separation of variables with constructing the implicit difference scheme for fractional diffusion-wave equation with damping. Heydari et al. ^[19] have proposed Legendre wavelets (LWs) for solving TFDWEs where fractional operational matrix of integration for LWs was derived. Bhrawy et al. ^[16] have proposed Jacobi tau spectral procedure combined with the Jacobi operational matrix for solving TFDWEs. Ebadian et al. ^[18] have proposed triangular function (TFs) methods for solving a class of nonlinear TFDWEs where fractional operational matrix of integration for the TFs was derived. Mohammed et al. ^[21] have proposed shifted Legengre collocation scheme and sinc function for solving TFDWEs with variable coefficients. Zhou et al. ^[22] have applied Chebyshev wavelets collocation for solving a class of TFDWEs where fractional integral formula of a single Chebyshev wavelets in the Riemann-Liouville sense was derived. Khalid et al. ^[20] have proposed the third degree modified extended B-spline functions for solving TFDWEs with reaction and damping terms. Some other numerical methods were presented for solving time fractional diffusion equations, one can see ^{[23,24,25,26]} and the references therein.

To the best of our knowledge, there is no existing numerical method which can be used to solve Eq (1.1) neither directly nor by transferring Eq (1.1) into an equivalent integro-differential equation. Thus, the aim of this study is devoted to constructing the high order numerical schemes to solve Eq (1.1), and carrying out the corresponding numerical analysis for the proposed schemes. Herein, we firstly transform Eq (1.1) into the equivalent partial integro-differential equations by using the integral operator. Secondly, the Crank-Nicolson technique is applied to deal with the temporal direction. Then, we use the midpoint formula to discretize the first order derivative, use the weighted and shifted Gr $\ddot{\text{u}}$ nwald difference formula to discretize the Caputo derivative, and apply the second order convolution quadrature formula to approximate the first order integral. The classical central difference formula, the fourth order Stephenson scheme, and the fourth order compact difference formula are applied for spatial approximations.

The rest of this paper is organized as follows. In Section 2, some preparations and useful lemmas are provided and discussed. In Section 3, the finite difference scheme is constructed and analyzed. In Section 4, the compact finite difference scheme is deduced, and the convergence and the unconditional stability are strictly proved. Numerical experiments are provided to support the theoretical results in Section 5. Finally, some concluding remarks are given.

2. Preliminaries

Lemma 2.1. (see Lemma 6.2 in ^[27]) Eq (1.1) is equivalent to the following partial integro-differential equation,

$\begin{align} \frac{\partial u(x, t)}{\partial t} + {}^{C}_{0}D_{t}^{\alpha-1}u(x, t) + K_{c}\cdot{}_{0}J_{t}\frac{\partial^{4}u(x, t)}{\partial x^4} = {}_{0}J_{t} \frac{\partial^2 u(x, t)}{\partial x^2} + {}_{0}J_{t}g(u) + F(x, t), \end{align}$

(2.1)

where $F(x, t) = {}_{0}J_{t}f(x, t)$ and ${}_{0}J_{t}$ is first order integral operator, i.e., ${}_{0}J_{t}u(\cdot, t) = \int_{0}^{t}u(\cdot, s)ds$ .

To discretize Eq (2.1), we introduce the temporal step size $\tau = T/N$ with a positive integer $N$ , $t_n = n\tau$ , and $t_{n+1/2} = (n + 1/2)\tau$ . Similarly, define the spatial step size $h = L/M$ with a positive integer $M$ , and denote $x_i = ih$ . Then, define a grid function space $\Theta_h = \{v_i^n|\text{ }0 \leq n \leq N, 0 \leq i \leq M, v^n_0 = v^n_M = 0\}$ , and introduce the following notations, inner product, and norm, i.e., for $u^n, v^n \in \Theta_h$ , we define

$\begin{align*} \Delta_xu_{i}^n = & \frac{1}{2h}\left(u_{i+1}^{n} - u_{i-1}^n\right), \quad\delta^2_xu_{i}^n = \frac{1}{h^2}\left(u_{i-1}^n - 2u_{i}^n + u^{n}_{i+1}\right), \\ \langle u^n, v^n\rangle = &h\sum\limits_{i = 1}^{M-1}u^n_{i}v^n_{i}, \quad||u^n||^2 = \langle u^n, u^n\rangle, \end{align*}$

$\mathcal{H}u_i^n = \left\{ \begin{array}{ll} \left(1+\frac{h^2}{12}\delta_x^2\right)u_i^n = \frac{1}{12}\left(u_{i-1}^n + 10u_{i}^n + u_{i+1}^n\right), & \mbox{ $1 \leq i\leq M-1$}, \\ u^n_i, & \mbox{ $i = 0$ or $M$}.\\ \end{array} \right.$

Lemma 2.2. (see Lemmas 2.2 and 2.3 in ^[28]) If $u(\cdot, t) \in C^2([0, T])$ and $0 < \gamma < 1$ , then it holds

$\begin{align*} {}_{0}J_{t}u(\cdot, t_{n+1/2}) = \frac{1}{2}\left[{}_{0}J_{t}u(\cdot, t_{n+1})+{}_{0}J_{t}u(\cdot, t_{n})\right] + O(\tau^2). \end{align*}$

Furthermore, if $u(\cdot, t)\in C^3([0, T])$ , then we have

$\begin{align*} u_t(\cdot, t_{n+1/2}) = \frac{u(\cdot, t_{n+1}) - u(\cdot, t_{n})}{\tau}+ O(\tau^2) = \delta_tu(\cdot, t_{n+1/2})+ O(\tau^2) \end{align*}$

and

$\begin{align*} {}^{C}_{0}D^{\gamma}_{t}u(\cdot, t_{n+1/2}) = \frac{1}{2}\left({}^{C}_{0}D^{\gamma}_{t}u(\cdot, t_{n+1})+{}^{C}_{0}D^{\gamma}_{t}u(\cdot, t_{n})\right) + O(\tau^2). \end{align*}$

Lemma 2.3. (see Theorem 4.1 in ^[29]) Let $\{\omega_{k}\}$ be the weights from generating function $\left(3/2 - 2z + z^2/2\right)^{-1}$ , i.e., $\omega_k = 1-3^{-(k+1)}$ . If $u(\cdot, t) \in C^2\left([0, T]\right)$ and $u(\cdot, 0) = u_t(\cdot, 0) = 0$ , then we have

$\begin{equation*} {}_{0}J_{t_{n+1}}u(\cdot, t) -\tau \sum\limits_{k = 0}^{n+1}\omega_{n+1-k}u(\cdot, t_{k}) = O(\tau^2). \end{equation*}$

Lemma 2.4. (see Theorem 2.4 in ^[30]) For $u(\cdot, t) \in L^1(\mathbb{R})$ , ${}_{-\infty}^{RL}D_{t}^{\gamma +2}u(\cdot, t)$ and its Fourier transform belong to $L^1(\mathbb{R})$ , if we use the weighted and shifted Gr $\ddot{\mathit{\text{u}}}$ nwald difference operator to approximate the Riemann-Liouville derivative, then it holds

$\begin{align*} {}_{0}^{RL}D^{\gamma}_{0}u(\cdot, t_{k+1}) = \tau^{-\gamma}\sum\limits_{j = 0}^{k+1}\sigma_{j}^{(\gamma)}u(\cdot, t_{k+1-j}) + O(\tau^2), \quad 0 < \gamma < 1, \end{align*}$

where

$\begin{align*} \sigma_{0}^{(\gamma)} = \frac{2 + \gamma}{2}c_{0}^{(\gamma)}, \; \sigma_j^{(\gamma)} = \frac{2 + \gamma}{2}c_j^{(\gamma)} - \frac{\gamma}{2}c_{j-1}^{(\gamma)}, \mathit{\text{}} j\geq1, \end{align*}$

and $c_j^{(\gamma)} = (-1)^j\left(\begin{array}{c} \gamma \\ j \end{array}\right)$ for $j\geq 0$ .

Lemma 2.5. (see Lemma 1.2 in ^[31]) Suppose $u(x, \cdot)\in C^4([x_{i-1}, x_{i+1}])$ , let $\zeta(s) = u^{(4)}(x_i+sh, \cdot) + u^{(4)}(x_i-sh, \cdot)$ , then

$\begin{align*} \delta_x^2u(x_i, \cdot) = \frac{u(x_{i-1}, \cdot) - 2u(x_{i}, \cdot) + u(x_{i+1}, \cdot)}{h^2} = u_{xx}(x_i, \cdot) + \frac{h^2}{24}\int_{0}^{1}\zeta(s)(1-s)^3ds. \end{align*}$

Lemma 2.6. (see Page 6 of ^[32]) Assume that $u(x, \cdot)\in C^8\left([0, L]\right)$ with $u(0, \cdot) = u(L, \cdot) = u_x(0, \cdot) = u_x(L, \cdot) = 0$ , and define the operator $\delta_x^4$ by

$\begin{align*} \delta_x^4u_i^n = \frac{12}{h^2}\left(\Delta_xv_i^n-\delta_x^2u_i^n\right), \end{align*}$

where $v_i^n$ is a compact approximation of $u_x(x_i, t_n)$ , i.e.,

$\begin{align*} \frac{1}{6}v_{i-1}^n + \frac{2}{3}v_{i}^n + \frac{1}{6}v_{i+1}^n = \Delta_xu_i^n. \end{align*}$

Then, we have the following approximation

$\begin{align*} \delta_x^4 u_i^n = \frac{\partial^4u(x_i, t_n)}{\partial x^4} + O(h^4). \end{align*}$

Furthermore, let $u^n = \left(u_1^n, u_2^n, \cdots, u_{M-1}^n\right)^T$ , then the matrix representation of the operator $\delta_x^4$ is

$\begin{align*} \mathbf{S} u^n = \frac{6}{h^4}\left(3\mathbf{K}\mathbf{P}^{-1}\mathbf{K} + 2\mathbf{D}\right)u^n, \end{align*}$

where

$\begin{align*} \mathbf{K} = \left( \begin{array}{ccccc} 0 &1 & & & \\ -1 & 0 & 1 & & \\ & \ddots & \ddots & \ddots & \\ & & -1 & 0 & 1 \\ & & & -1 & 0 \\ \end{array} \right)_{(M-1)\times(M-1)}\mathit{\text{, }} \mathbf{P} = \left( \begin{array}{ccccc} 4 &1 & & & \\ 1 & 4 & 1 & & \\ & \ddots & \ddots & \ddots & \\ & & 1 & 4 & 1 \\ & & & 1 & 4 \\ \end{array} \right)_{(M-1)\times(M-1)}, \end{align*}$

and $\mathbf{D} = 6\mathbf{I} - \mathbf{P}$ with the identity matrix $\mathbf{I}$ .

Lemma 2.7. (see Lemma 3.3 in ^[32]) The matrix $\mathbf{S}$ defined in Lemma 2.6 is symmetric positive definite.

It follows from Lemma 2.7, there is an invertible matrix $\mathbf{B}$ such that, $\mathbf{S} = \mathbf{B}^T\mathbf{B}$ . Then for $w^n, v^n \in \Theta_h$ , we have

$\begin{align} \langle \mathbf{S}w^n, v^n\rangle = \langle \mathbf{B}^T\mathbf{B}w^n, v^n\rangle = \langle\mathbf{B}w^n, \mathbf{B}v^n\rangle. \end{align}$

(2.2)

The following lemma is required when we use compact operator $\mathcal{H}$ to increase the spatial accuracy.

Lemma 2.8. (see Lemma 1.2 in ^[31]) Suppose $u(x, \cdot) \in C^6\left([x_{i-1}, x_{i+1}]\right)$ , $1 \leq i\leq M-1$ , and $\zeta(s) = 5\left(1-s\right)^3 - 3\left(1-s\right)^5$ . Then it holds that

$\begin{align*} &\frac{1}{12}\left[u_{xx}(x_{i-1}, \cdot) + 10u_{xx}(x_{i}, \cdot) + u_{xx}(x_{i+1}, \cdot)\right] - \frac{1}{h^2}\left[u(x_{i-1}, \cdot) - 2u(x_{i}, \cdot) + u(x_{i+1}, \cdot)\right]\\ & = \frac{h^4}{360}\int_{0}^{1}\left[u^{(6)}(x_i-sh, \cdot) + u^{(6)}(x_i+sh, \cdot)\right]\zeta(s)ds. \end{align*}$

In order to linearize the nonlinear function $g(u)$ , we can easily get the following lemma by Taylor expansions.

Lemma 2.9. Assume that $u(\cdot, t)\in C^1([0, T]) \cap C^2((0, T])$ , then the following approximation holds

$\begin{align*} u(\cdot, t_{n+1}) = 2u(\cdot, t_n) - u(\cdot, t_{n-1}) + O(\tau^2). \end{align*}$

3. Derivation and analysis of the finite difference scheme

3.1. The derivation of the finite difference scheme

In this subsection, a finite difference scheme with the accuracy $O(\tau^2 + h^2)$ for nonlinear Problem (2.1) is constructed.

Assume that $u(x, t) \in C_{x, t}^{8, 3}([0, L]\times[0, T])$ , and $u(\cdot, 0) = u_t(\cdot, 0) = 0$ . Consider Eq (2.1) at the point $u(x_i, t_{n+1/2})$ , we have

$\begin{align*} \left.\frac{\partial u(x_i, t)}{\partial t}\right|_{t = t_{n+1/2}} = & -{}^{C}_{0}D^{\alpha-1}_{t_{n+1/2}}u(x_i, t) - K_c\cdot {}_{0}J_{t_{n+1/2}}\frac{\partial^4 u(x_i, t)}{\partial x^4} + {}_{0}J_{t_{n+1/2}}\frac{\partial^2 u(x_i, t)}{\partial x^2} \\ &+ {}_{0}J_{t_{n+1/2}}g(u(x_i, t)) + F(x_i, t_{n+1/2}). \end{align*}$

The Crank-Nicolson technique and Lemma 2.2 for the above equation yield

$\begin{align} \frac{u(x_i, t_{n+1})-u(x_i, t_{n})}{\tau} = &-\frac{1}{2}\left[{}^{C}_{0}D^{\alpha-1}_{t_{n+1}}u(x_i, t)+{}^{C}_{0}D^{\alpha-1}_{t_{n}}u(x_i, t)\right]\nonumber-\frac{K_c}{2}\left[ {}_{0}J_{t_{n+1}}\frac{\partial^4 u(x_i, t)}{\partial x^4} + {}_{0}J_{t_{n}}\frac{\partial^4 u(x_i, t)}{\partial x^4} \right]\nonumber\\ &+\frac{1}{2}\left[{}_{0}J_{t_{n+1}}\frac{\partial^2 u(x_i, t)}{\partial x^2}+{}_{0}J_{t_{n}}\frac{\partial^2 u(x_i, t)}{\partial x^2}\right]+\frac{1}{2}\left[{}_{0}J_{t_{n+1}}g(x_i, t)+{}_{0}J_{t_{n}}g(x_i, t)\right]\nonumber\\ &+ F(x_i, t_{n+1/2})+O(\tau^2). \end{align}$

(3.1)

Let $u(x_i, t_n) = u_i^n$ . Since the initial values are $0$ , thus the Riemann $-$ liouville derivative is equivalent to Caputo derivative. We apply Lemmas 2.3 and 2.4 to discretize the first order integral operator and Caputo derivative in Eq (3.1) respectively, apply Lemma 2.6 to discretize $\frac{\partial^4u(x_i, t)}{\partial x^4}$ , and Lemma 2.5 to discretize $\frac{\partial^2u(x_i, t)}{\partial x^2}$ , then we get

$\begin{align} \frac{u_i^{n+1}-u_i^{n}}{\tau} = & -\frac{\tau^{1-\alpha}}{2}\left[\sum\limits_{k = 0}^{n+1}\sigma_k^{(\alpha-1)}u_i^{n+1-k}+\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}u_i^{n-k}\right]-\frac{K_c\tau}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^4u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^4u_i^{n-k}\right]\\ &+\frac{\tau}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^2u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^2u_i^{n-k}\right]+\frac{\tau}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_kg(u_i^{n+1-k}) + \sum\limits_{k = 0}^{n}\omega_kg(u_i^{n-k})\right]\\ &+F_i^{n+\frac{1}{2}}+(R_1)_i^{n+1}, \end{align}$

(3.2)

where $(R_1)_i^{n+1} = O(\tau^2 + h^2 + h^4) = O(\tau^2 + h^2).$

It is clear that Eq (3.2) is a nonlinear system with respect to the unknown $u_i^{n+1}$ . To linearly solve Eq (3.2), we use $u_i^1 = u_i^0 + \tau(u_t)_i^0 +O(\tau^2)$ and Lemma 2.9 to linearize Eq (3.2) for $n = 0$ and $1\leq n \leq N-1$ , respectively, and then multiply Eq (3.2) by $\tau$ , i.e.,

$\begin{align} u_i^{1} - u_i^{0} = &-\frac{\tau^{2-\alpha}}{2}\left[\sum\limits_{k = 0}^{1}\sigma_k^{(\alpha-1)}u_i^{1-k}+\sigma_0^{(\alpha-1)}u_i^0\right]-\frac{K_c\tau^2}{2}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^4u_i^{1-k} + \omega_0\delta_x^4u_i^0\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^2u_i^{1-k} + \omega_0\delta_x^2u_i^0\right]+\frac{\tau^2}{2}\left[\omega_0g(u_i^0+\tau(u_t)_i^0) + \omega_1g(u_i^0) + \omega_0g(u_i^0)\right]\\ &+\tau F_i^{n+\frac{1}{2}}+O(\tau^{3}+\tau h^2) \end{align}$

(3.3)

and

$\begin{align} u_i^{n+1} - u_i^{n} = &-\frac{\tau^{2-\alpha}}{2}\left[\sum\limits_{k = 0}^{n+1}\sigma_k^{(\alpha-1)}u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}u_i^{n-k}\right]-\frac{K_c\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^4u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^4u_i^{n-k}\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^2u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^2u_i^{n-k}\right]+\frac{\tau^2}{2}\left[\sum\limits_{k = 1}^{n+1}\omega_kg(u_i^{n+1-k}) + \sum\limits_{k = 0}^{n}\omega_kg(u_i^{n-k}) \right]\\ &+\frac{\tau^2 \omega_0}{2}g(2u_i^n-u_i^{n-1}) + \tau F_i^{n+\frac{1}{2}} + O(\tau^{3} + \tau h^2), \text{ for } 1\leq n \leq N-1. \end{align}$

(3.4)

Noting $(u_t)_i^0 = 0$ , neglecting the truncation error term $O(\tau^{3} + \tau h^2)$ in both above equations, and replacing the $u_i^n$ with its numerical solution $U_i^n$ , we deduce the following finite difference scheme for Problem (2.1)

$\begin{align} U_i^{1} - U_i^{0} = &-\frac{\tau^{2-\alpha}}{2}\left[\sum\limits_{k = 0}^{1}\sigma_k^{(\alpha-1)}U_i^{1-k}+\sigma_0^{(\alpha-1)}U_i^0\right]-\frac{K_c\tau^2}{2}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^4U_i^{1-k} + \omega_0\delta_x^4U_i^0\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^2U_i^{1-k} + \omega_0\delta_x^2U_i^0\right]+\frac{\tau^2}{2}\left[\omega_0g(U_i^0) + \omega_1g(U_i^0) + \omega_0g(U_i^0)\right]\\ &+ \tau F_i^{n+\frac{1}{2}} \end{align}$

(3.5)

and

$\begin{align} U_i^{n+1} - U_i^{n} = &-\frac{\tau^{2-\alpha}}{2}\left[\sum\limits_{k = 0}^{n+1}\sigma_k^{(\alpha-1)}U_i^{n+1-k} + \sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}U_i^{n-k}\right]-\frac{K_c\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^4U_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^4U_i^{n-k}\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^2U_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^2U_i^{n-k}\right]+\frac{\tau^2}{2}\left[\sum\limits_{k = 1}^{n+1}\omega_kg(U_i^{n+1-k}) + \sum\limits_{k = 0}^{n}\omega_kg(U_i^{n-k}) \right]\\ &+\frac{\tau^2 \omega_0}{2}g(2U_i^n-U_i^{n-1}) + \tau F_i^{n+\frac{1}{2}}, \text{ for } 1\leq n \leq N-1. \end{align}$

(3.6)

Remark 3.1. In case of $g(u) = f(x, t) = 0$ , the only solution of the finite difference Scheme (3.5) and (3.6) is zero solution.

3.2. Analysis of the finite difference Scheme (3.5) and (3.6)

In this subsection, the convergence and stability of the finite difference Scheme (3.5) and (3.6) will be discussed. For convenience, let $C$ be a generic constant, whose value is independent of discretization parameters and may be different from one line to another. To begin, we provide two lemmas that will be used in our convergence and stability analysis.

Lemma 3.2. (see Proposition 5.2 in ^[33] and Lemma 3.2 in ^[34]) Let $\{\omega_k\}$ and $\{\sigma_k^{(\alpha-1)}\}$ be the weights defined in Lemmas 2.3 and 2.4, respectively. Then for any positive integer $K$ and real vector $\left(V_1, V_2, \cdots, V_K\right)^{T}$ , the inequalities

$\begin{align*} \sum\limits_{n = 0}^{K-1}\left(\sum\limits_{j = 0}^{n}\omega_jV_{n+1-j}\right)V_{n+1}\geq 0 \end{align*}$

and

$\begin{align*} \sum\limits_{n = 0}^{K-1}\left(\sum\limits_{j = 0}^{n}\sigma_j^{(\alpha -1)}V_{n+1-j}\right)V_{n+1}\geq 0 \end{align*}$

hold.

Lemma 3.3. (see Lemma 4.2.2 in ^[35]) For any grid function $w^n, v^n \in \Theta_h$ , it holds

$\begin{align*} \langle \delta_x^2w^n, v^n \rangle = -\langle \delta_x w^n, \delta_x v^n \rangle. \end{align*}$

Theorem 3.4. Assume $u(x, t)\in C_{x, t}^{8, 3}\left([0, L] \times[0, T]\right)$ and $u(\cdot, 0) = u_t(\cdot, 0) = 0$ , and let $u(x, t)$ be the exact solution of Eq (2.1) and $\{U_i^n|\mathit{\text{}}0\leq i \leq M, 1 \leq n \leq N\}$ be the numerical solution for Scheme (3.7) and (3.8). Then, for $1 \leq n \leq N$ , it holds that

$\begin{align*} \|u^n - U^n\| \leq C(\tau^{2} + h^2). \end{align*}$

Proof. Let us start by analyzing the error of (3.6). Subtracting Eq (3.6) from Eq (3.4), we have

$\begin{align*} e_i^{n+1} - e_i^n = & -\frac{\tau^{2-\alpha}}{2}\left[\sum\limits_{k = 0}^{n+1}\sigma_k^{(\alpha-1)}e_i^{n+1-k} + \sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}e_i^{n-k} \right]\\ &-\frac{K_c\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^4e_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^4e_i^{n-k}\right]+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^2e_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^2e_i^{n-k}\right]\\ &+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}\left(\omega_{k+1} + \omega_k\right)\left[g(u_i^{n-k}) - g(U_i^{n-k})\right]+\frac{\tau^2 \omega_0}{2}\left[g(2u_i^{n} - u_i^{n-1}) - g(2U_i^{n} - U_i^{n-1})\right]\\ &+O(\tau^{3} + \tau h^2), \end{align*}$

where $e_i^n = u_i^n - U_i^n$ . Since $e_i^0 = 0$ , the above equation becomes

$\begin{align*} e_i^{n+1} - e_i^n = & -\frac{\tau^{2-\alpha}}{2}\left[\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}(e_i^{n+1-k} + e_i^{n-k})\right]-\frac{K_c\tau^2}{2}\left[\sum\limits_{k = 0}^{n}\omega_k\delta_x^4\left(e_i^{n+1-k} + e_i^{n-k}\right) \right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n}\omega_k\delta_x^2\left(e_i^{n+1-k} + e_i^{n-k}\right)\right]+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}\left(\omega_{k+1} + \omega_k\right)\left[g(u_i^{n-k}) - g(U_i^{n-k})\right]\\ &+\frac{\tau^2 \omega_0}{2}\left[g(2u_i^{n} - u_i^{n-1}) - g(2U_i^{n} - U_i^{n-1})\right] + O(\tau^{3} + \tau h^2). \end{align*}$

Multiplying the both sides of the above equation by $h(e_i^{n+1} + e_i^{n})$ and summing over $1\leq i \leq M-1$ . Then using Lemmas 3.3, 2.6, and Eq (2.2), we have

$\begin{align*} \|e^{n+1}\|^2 -\|e^n\|^2 = & -\frac{\tau^{2-\alpha}}{2}\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}\langle e^{n+1-k}+e^{n-k}, e^{n+1} + e^n\rangle\\ &- \frac{K_c\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \mathbf{B}(e^{n+1-k} + e^{n-k}), \mathbf{B}(e^{n+1} + e^n)\rangle\\ &-\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \delta_x(e^{n+1-k} + e^{n-k}), \delta_x(e^{n+1} + e^n)\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle g(u^{n-k})-g(U^{n-k}), e^{n+1} + e^n\rangle\\ &+\frac{\tau^2 \omega_0}{2}\langle g(2u^n -u^{n-1})- g(2U^n -U^{n-1}), e^{n+1} + e^n \rangle\\ &+\langle O(\tau^{3} + \tau h^2), e^{n+1} + e^n\rangle. \end{align*}$

Summing the above equation over $n$ from $1$ to $J-1$ leads to

$\begin{align} \|e^{J}\|^2 -\|e^1\|^2 = & -\frac{\tau^{2-\alpha}}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}\langle e^{n+1-k}+e^{n-k}, e^{n+1} + e^n\rangle\\ &- \frac{K_c\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}\omega_k\langle \mathbf{B}(e^{n+1-k} + e^{n-k}), \mathbf{B}(e^{n+1} + e^n)\rangle\\ &-\frac{\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}\omega_k\langle \delta_x(e^{n+1-k} + e^{n-k}), \delta_x(e^{n+1} + e^n)\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle g(u^{n-k})-g(U^{n-k}), e^{n+1} + e^n\rangle\\ &+\frac{\tau^2 \omega_0}{2}\sum\limits_{n = 1}^{J-1}\langle g(2u^n -u^{n-1})- g(2U^n -U^{n-1}), e^{n+1} + e^n \rangle\\ &+\sum\limits_{n = 1}^{J-1}\langle O(\tau^{3} + \tau h^2), e^{n+1} + e^n\rangle. \end{align}$

(3.7)

Now, we turn to analyze $\|e^1\|$ . Subtracting Eq (3.5) from Eq (3.3), and by the similar deductions as above, we can derive that

$\begin{align} \|e^1\|^2 = & -\frac{\tau^{2-\alpha}}{2}\sigma_0^{(\alpha-1)}\langle e^{1}+e^{0}, e^{1} + e^0\rangle - \frac{K_c\tau^2}{2}\omega_0\langle \mathbf{B}(e^{1} + e^{0}), \mathbf{B}(e^{1} + e^0)\rangle\\ &-\frac{\tau^2}{2}\omega_0\langle \delta_x(e^{1} + e^{0}), \delta_x(e^{1} + e^0)\rangle+\tau^2\omega_0\langle g(u^{0})-g(U^{0}), e^{1} + e^0\rangle\\ &+\frac{\tau^2 \omega_1}{2}\langle g(u^0)- g(U^0), e^{1} + e^0 \rangle +\langle O(\tau^{3} + \tau h^2), e^{1} + e^0\rangle. \end{align}$

(3.8)

Sum up Eq (3.7) and Eq (3.8), and apply Lemma 3.2, it deduces that

$\begin{align} \|e^J\|^2 \leq& \frac{\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle g(u^{n-k})-g(U^{n-k}), e^{n+1} + e^n\rangle\\ &+\frac{\tau^2 \omega_0}{2}\sum\limits_{n = 1}^{J-1}\langle g(2u^n -u^{n-1})- g(2U^n -U^{n-1}), e^{n+1} + e^n \rangle\\ &+\tau^2\omega_0\langle g(u^{0})-g(U^{0}), e^{1} + e^0\rangle+\frac{\tau^2 \omega_1}{2}\langle g(u^0)- g(U^0), e^{1} + e^0 \rangle\\ &+C\sum\limits_{n = 1}^{J-1}\langle \tau^{3} + \tau h^2, e^{n+1} + e^n\rangle. \end{align}$

(3.9)

Using the Lipschitz condition of $g$ and exchanging the order of two summations in the above inequality, we have

$\begin{align} \|e^{J}\|^2 \leq& C\tau^{2}\sum\limits_{k = 0}^{J-1}\sum\limits_{n = k}^{J-1}(\omega_{n+1-k} + \omega_{n-k})\|e^k\|\|e^{n+1}+e^n\|+C\tau^2 \sum\limits_{n = 1}^{J-1}\|e^n\|\|e^{n+1} + e^n\|+ C\sum\limits_{n = 1}^{J-1}(\tau^{3} +\tau h^2)\|e^{n+1} + e^n \|. \end{align}$

(3.10)

Assuming ${\|e^P\| = \max\limits_{0 \leq p \leq N} \|e^p\|}$ . Since ${\tau \sum\limits_{n = k}^{N}\left(\omega_{n+1-k} + \omega_{n-k}\right)}$ is bounded (see ^[29]), then the above inequality yields

$\begin{align} \|e^P\| \leq C \tau \sum\limits_{k = 0}^{P-1}\|e^k\| + C(\tau^{2} + h^2). \end{align}$

(3.11)

Once the discrete Gronwall inequality has been applied to Inequality (3.11), we arrive at the estimate

$\begin{align*} \|e^P\| \leq C(\tau^{2} + h^2), \end{align*}$

thus the proof is completed.

Theorem 3.5. Let $\{U_i^n|\mathit{\text{}}0\leq i \leq M, 0\leq n \leq N\}$ be the numerical solution of Scheme (3.5) and (3.6) for Problem (2.1). Then for $1 \leq K \leq N$ , it holds

$\begin{align} \|U^K\| \leq C\left(\max\limits_{0\leq n \leq N} \|g(U^n)\| +\max\limits_{0\leq n \leq N-1} \|F^{n+\frac{1}{2}}\| \right). \end{align}$

(3.12)

Proof. Multiplying (3.6) by $h(U_i^{n+1}+U_i^n)$ and summing up for $i$ from $1$ to $M - 1$ , we have

$\begin{align*} \|U^{n+1}\|^2 -\|U^n\|^2 = & -\frac{\tau^{2-\alpha}}{2}\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}\langle U^{n+1-k}+U^{n-k}, U^{n+1} + U^n\rangle\\ &- \frac{K_c\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \delta_x^4(U^{n+1-k} + U^{n-k}), U^{n+1} + U^n\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \delta^2_x(U^{n+1-k} + U^{n-k}), U^{n+1} + U^n\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle g(U^{n-k}), U^{n+1} + U^n\rangle\\ &+\frac{\tau^2 \omega_0}{2}\langle g(2U^n -U^{n-1}), U^{n+1} + U^n \rangle\\ & - \frac{K_c\tau^2}{2}\omega_{n+1}\langle \delta_x^4U^{0}, U^{n+1} + U^n\rangle -\frac{\tau^{2-\alpha}}{2}\sigma_{n+1}^{(\alpha-1)}\langle U^{0}, U^{n+1} + U^n\rangle \\ &+\frac{\tau^2}{2}\omega_{n+1}\langle \delta^2_xU^{0}, U^{n+1} + U^n\rangle+\tau\langle F^{n+\frac{1}{2}}, U^{n+1} + U^n\rangle. \end{align*}$

Note that Eq (1.1) is equipped with the homogeneous initial conditions, thus it deduces

$\begin{align*} \|U^{n+1}\|^2 -\|U^n\|^2 = & -\frac{\tau^{2-\alpha}}{2}\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}\langle U^{n+1-k}+U^{n-k}, U^{n+1} + U^n\rangle\\ &- \frac{K_c\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \delta_x^4(U^{n+1-k} + U^{n-k}), U^{n+1} + U^n\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \delta^2_x(U^{n+1-k} + U^{n-k}), U^{n+1} + U^n\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle g(U^{n-k}), U^{n+1} + U^n\rangle\\ &+\frac{\tau^2 \omega_0}{2}\langle g(2U^n -U^{n-1}), U^{n+1} + U^n \rangle \\ &+\tau\langle F^{n+\frac{1}{2}}, U^{n+1} + U^n\rangle. \end{align*}$

Applying the similar deductions to get Eq (3.9), it achieves that

$\begin{align} \|U^{J}\|^2 \leq& C \tau \sum\limits_{k = 0}^{J-1}\|g(U^k)\|\left(\|U^{n+1}\|+\|U^{n}\|\right) +\frac{\tau^2}{2} \omega_0\sum\limits_{n = 1}^{J-1}\|g(2U^{n}-U^{n-1})\|\left(\|U^{n+1}\|+\|U^n\|\right)\\ &+C\tau \sum\limits_{n = 1}^{J-1}\|F^{n+\frac{1}{2}}\|\left(\|U^{n+1}\|+\|U^n\|\right). \end{align}$

(3.13)

One can estimate $\| g(2U^n -U^{n-1}) \|$ as the following

$\begin{align} \| g(2U^n -U^{n-1})\| = & \| g(2U^n -U^{n-1})-g(U^n)+g(U^n)\|, \\ \leq&\| g(2U^n -U^{n-1})-g(U^n)\| + \| g(U^n)\|, \\ \leq& C(\| U^n\| + \|U^{n-1}\|) + \| g(U^n)\|. \end{align}$

(3.14)

Substituting Eq (3.14) into Eq (3.13) and using Young's inequality, then we have

$\begin{align} \|U^J\|^2 \leq C\tau \sum\limits_{n = 0}^{J-1}\|U^n\|^2 + C\max\limits_{0\leq n\leq N}\|g(U^n)\|^2 + C\max\limits_{0\leq n\leq N-1}\|F^{n+\frac{1}{2}}\|^2. \end{align}$

(3.15)

By applying the Gronwall inequality to (3.15), it becomes

$\begin{align*} \|U^J\|^2 \leq C\left(\max\limits_{0\leq n\leq N}\|g(U^n)\|^2 + \max\limits_{0\leq n\leq N-1}\|F^{n+\frac{1}{2}}\|^2\right), \end{align*}$

and this completes the proof.

4. Derivation and analysis of the compact finite difference scheme

4.1. The derivation of the compact finite difference scheme

In this subsection, a compact finite difference scheme with accuracy $O(\tau^2 + h^4)$ for nonlinear Problem (2.1) is presented.

Now let us act on both sides of Eq (3.1) with the compact operator $\mathcal{H}$ . Then, by using Lemma 2.8, we obtain

$\begin{align} \mathcal{H}\left[\frac{u(x_i, t_{n+1})-u(x_i, t_{n})}{\tau}\right] = &-\frac{1}{2}\mathcal{H}\left[{}^{C}_{0}D^{\alpha-1}_{t_{n+1}}u(x_i, t)+{}^{C}_{0}D^{\alpha-1}_{t_{n}}u(x_i, t)\right]\\ &- \frac{K_c}{2}\mathcal{H}\left[ {}_{0}J_{t_{n+1}}\frac{\partial^4 u(x_i, t)}{\partial x^4} + {}_{0}J_{t_{n}}\frac{\partial^4 u(x_i, t)}{\partial x^4} \right]\\ &+ \frac{1}{2}\left[{}_{0}J_{t_{n+1}}\delta_x^2 u(x_i, t)+{}_{0}J_{t_{n}}\delta_x^2 u(x_i, t)\right] \\ &+ \frac{1}{2}\mathcal{H}\left[{}_{0}J_{t_{n+1}}g(x_i, t)+{}_{0}J_{t_{n}}g(x_i, t)\right]\\ &+ \mathcal{H}F_i^{n+\frac{1}{2}}+O(\tau^2+h^4) . \end{align}$

(4.1)

Apply the similar deductions to get Eqs (3.3) and (3.4), it achieves

$\begin{align} \mathcal{H}\left[u_i^{1} - u_i^{0}\right] = &-\frac{\tau^{2-\alpha}}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{1}\sigma_k^{(\alpha-1)}u_i^{1-k}+\sigma_0^{(\alpha-1)}u_i^0\right]-\frac{K_c\tau^2}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^4u_i^{1-k} + \omega_0\delta_x^4u_i^0\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^2u_i^{1-k} + \omega_0\delta_x^2u_i^0\right] + \frac{\tau^2}{2}\mathcal{H}\left[\omega_0g(u_i^0) + \omega_1g(u_i^0) + \omega_0g(u_i^0)\right] \\ &+ \tau \mathcal{H}F_i^{n+\frac{1}{2}} +O(\tau^3 +\tau h^4) \end{align}$

(4.2)

and

$\begin{align} \mathcal{H}\left[u_i^{n+1} - u_i^{n}\right] = &-\frac{\tau^{2-\alpha}}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{n+1}\sigma_k^{(\alpha-1)}u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}u_i^{n-k}\right] - \frac{K_c\tau^2}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^4u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^4u_i^{n-k}\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^2u_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^2u_i^{n-k}\right]+\frac{\tau^2}{2}\mathcal{H}\left[\sum\limits_{k = 1}^{n+1}\omega_kg(u_i^{n+1-k}) + \sum\limits_{k = 0}^{n}\omega_kg(u_i^{n-k}) \right]\\ &+\frac{\tau^2 \omega_0}{2}\mathcal{H}g(2u_i^n-u_i^{n-1}) +\tau \mathcal{H}F_i^{n+\frac{1}{2}} + O(\tau^3 +\tau h^4), \text{ for } 1\leq n \leq N-1. \end{align}$

(4.3)

Neglecting the truncation error term $O(\tau^{3} + \tau h^4)$ in both above equations, and replacing the $u_i^n$ with its numerical solution $U_i^n$ , we deduce the following compact finite difference scheme for Problem (2.1)

$\begin{align} \mathcal{H}\left[U_i^{1} - U_i^{0}\right] = &-\frac{\tau^{2-\alpha}}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{1}\sigma_k^{(\alpha-1)}U_i^{1-k}+\sigma_0^{(\alpha-1)}U_i^0\right]-\frac{K_c\tau^2}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^4U_i^{1-k} + \omega_0\delta_x^4U_i^0\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{1}\omega_k\delta_x^2U_i^{1-k} + \omega_0\delta_x^2U_i^0\right]+\frac{\tau^2}{2}\mathcal{H}\left[\omega_0g(U_i^0) + \omega_1g(U_i^0) + \omega_0g(U_i^0)\right]\\ & + \tau \mathcal{H}F_i^{n+\frac{1}{2}} \end{align}$

(4.4)

and

$\begin{align} \mathcal{H}\left[U_i^{n+1} - U_i^{n}\right] = &-\frac{\tau^{2-\alpha}}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{n+1}\sigma_k^{(\alpha-1)}U_i^{n+1-k} + \sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}U_i^{n-k}\right]\\ &-\frac{K_c\tau^2}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^4U_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^4U_i^{n-k}\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^2U_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^2U_i^{n-k}\right]\\ &+\frac{\tau^2}{2}\mathcal{H}\left[\sum\limits_{k = 1}^{n+1}\omega_kg(U_i^{n+1-k}) + \sum\limits_{k = 0}^{n}\omega_kg(U_i^{n-k}) \right]\\ &+\frac{\tau^2 \omega_0}{2}\mathcal{H}g(2U_i^n-U_i^{n-1}) +\tau \mathcal{H}F_i^{n+\frac{1}{2}}, \text{ for } 1\leq n \leq N-1. \end{align}$

(4.5)

Remark 4.1. In case of $g(u) = f(x, t) = 0$ , the only solution of the compact finite difference Scheme (4.4) and (4.5) is zero solution.

4.2. Analysis of the compact finite difference Scheme (4.4) and (4.5)

In this subsection, we turn to analyze the convergence and stability of the compact finite difference Scheme (4.4) and (4.5). Firstly, we provide the following lemmas, which will be used in our convergence and stability analysis.

Lemma 4.2. (see Lemma 5 in ^[36]) Let $\{\sigma_k^{(\alpha-1)}\}$ be the weighted coefficients defined in Lemma 2.4, then for any positive integer $n$ and $w^n \in \Theta_h$ , it holds that

$\begin{align*} \sum\limits_{m = 0}^{n}\sum\limits_{k = 0}^{m}\sigma_{k}^{(\alpha-1)}\langle \mathcal{H}w^{m-k}, w^{m}\rangle \geq 0. \end{align*}$

Lemma 4.3. (see Lemma 4.2 in ^[37]) For any grid function $w^n \in \Theta_h$ , we have

$\begin{align*} \frac{2}{3}\|w^n\|^2 \leq \langle \mathcal{H}w^n, w^n\rangle \leq \|w^n\|^2. \end{align*}$

Theorem 4.4. Assume $u(x, t)\in C_{x, t}^{8, 3}\left([0, L] \times[0, T]\right)$ and $u(\cdot, 0) = u_t(\cdot, 0) = 0$ , and let $u(x, t)$ be the exact solution of Eq (2.1) and $\{U_i^n|\mathit{\text{}}0\leq i \leq M, 1 \leq n \leq N\}$ be the numerical solution for Scheme (4.4) and (4.5). Then, for $1 \leq n \leq N$ , it holds that

$\begin{align*} \|u^n - U^n\| \leq C(\tau^{2} + h^4). \end{align*}$

Proof. Let us start by analyzing the error of (4.5). Subtracting Eq (3.5) from Eq (4.3), we have

$\begin{align*} \mathcal{H}\left[e_i^{n+1} - e_i^n\right] = & -\frac{\tau^{2-\alpha}}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{n+1}\sigma_k^{(\alpha-1)}e_i^{n+1-k} + \sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}e_i^{n-k} \right]-\frac{K_c\tau^2}{2}\mathcal{H}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^4e_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^4e_i^{n-k}\right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n+1}\omega_k\delta_x^2e_i^{n+1-k} + \sum\limits_{k = 0}^{n}\omega_k\delta_x^2e_i^{n-k}\right]+\frac{\tau^2}{2}\mathcal{H}\sum\limits_{k = 0}^{n}\left(\omega_{k+1} + \omega_k\right)\left[g(u_i^{n-k}) - g(U_i^{n-k})\right]\\ &+\frac{\tau^2 \omega_0}{2}\mathcal{H}\left[g(2u_i^{n} - u_i^{n-1}) - g(2U_i^{n} - U_i^{n-1})\right]+O(\tau^{3} + \tau h^4), \end{align*}$

where $e_i^n = u_i^n - U_i^n$ . Since $e_i^0 = 0$ , the above equation becomes

$\begin{align*} \mathcal{H}\left[e_i^{n+1} - e_i^n\right] = & -\frac{\tau^{2-\alpha}}{2}\left[\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}\mathcal{H}(e_i^{n+1-k} + e_i^{n-k})\right]-\frac{K_c\tau^2}{2}\left[\sum\limits_{k = 0}^{n}\omega_k\mathcal{H}\delta_x^4\left(e_i^{n+1-k} + e_i^{n-k}\right) \right]\\ &+\frac{\tau^2}{2}\left[\sum\limits_{k = 0}^{n}\omega_k\delta_x^2\left(e_i^{n+1-k} + e_i^{n-k}\right)\right] +\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}\left(\omega_{k+1} + \omega_k\right)\mathcal{H}\left[g(u_i^{n-k}) - g(U_i^{n-k})\right]\\ &+\frac{\tau^2 \omega_0}{2}\mathcal{H}\left[g(2u_i^{n} - u_i^{n-1}) - g(2U_i^{n} - U_i^{n-1})\right] + O(\tau^{3} + \tau h^4). \end{align*}$

Multiplying the both sides of the above equation by $h(e_i^{n+1} + e_i^{n})$ and summing over $1\leq i \leq M-1$ . Then using Lemmas 2.6, 3.2, 4.2, and Eq (2.2), we have

$\begin{align*} \|e^{n+1}\|^2 -\|e^n\|^2 \leq& -\frac{\tau^{2-\alpha}}{2}\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}\langle \mathcal{H}(e^{n+1-k}+e^{n-k}), e^{n+1} + e^n\rangle\\ &- \frac{K_c\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \mathcal{H}\mathbf{B}(e^{n+1-k} + e^{n-k}), \mathbf{B}(e^{n+1} + e^n)\rangle\\ &-\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}\omega_k\langle \delta_x(e^{n+1-k} + e^{n-k}), \delta_x(e^{n+1} + e^n)\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle \mathcal{H}\left(g(u^{n-k})-g(U^{n-k})\right), e^{n+1} + e^n\rangle\\ &+\frac{\tau^2 \omega_0}{2}\langle \mathcal{H} \left(g(2u^n -u^{n-1})- g(2U^n -U^{n-1})\right), e^{n+1} + e^n \rangle\\ &+C\langle \tau^{3} + \tau h^4, e^{n+1} + e^n\rangle. \end{align*}$

Summing the above inequality over $n$ from $1$ to $J-1$ leads to

$\begin{align} \|e^{J}\|^2 -\|e^1\|^2 \leq& -\frac{\tau^{2-\alpha}}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}\sigma_k^{(\alpha-1)}\langle \mathcal{H}(e^{n+1-k}+e^{n-k}), e^{n+1} + e^n\rangle\\ &- \frac{K_c\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}\omega_k\langle \mathcal{H}B(e^{n+1-k} + e^{n-k}), B(e^{n+1} + e^n)\rangle\\ &-\frac{\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}\omega_k\langle \delta_x(e^{n+1-k} + e^{n-k}), \delta_x(e^{n+1} + e^n)\rangle\\ &+\frac{\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle \mathcal{H}\left(g(u^{n-k})-g(U^{n-k})\right), e^{n+1} + e^n\rangle\\ &+\frac{\tau^2 \omega_0}{2}\sum\limits_{n = 1}^{J-1}\langle \mathcal{H} \left(g(2u^n -u^{n-1})-g(2U^n -U^{n-1})\right), e^{n+1} + e^n \rangle\\ &+C\sum\limits_{n = 1}^{J-1}\langle \tau^{3} + \tau h^4, e^{n+1} + e^n\rangle. \end{align}$

(4.6)

Now, we turn to analyze $\|e^1\|$ . From Eqs (4.4), (4.2), and by the similar deductions as above, we can derive that

$\begin{align} \|e^1\|^2 \leq& -\frac{\tau^{2-\alpha}}{2}\sigma_0^{(\alpha-1)}\langle \mathcal{H}(e^{1}+e^{0}), e^{1} + e^0\rangle - \frac{K_c\tau^2}{2}\omega_0\langle \mathcal{H}\mathbf{B}(e^{1} + e^{0}), \mathbf{B}(e^{1} + e^0)\rangle\\ &-\frac{\tau^2}{2}\omega_0\langle \delta_x(e^{1} + e^{0}), \delta_x(e^{1} + e^0)\rangle +\frac{\tau^2\omega_{1}}{2}\langle \mathcal{H}\left(g(u^{0})-g(U^{0})\right), e^{1} + e^0\rangle\\ &+\tau^2 \omega_0\langle \mathcal{H} \left(g(u^0)- g(U^0)\right), e^{1} + e^0 \rangle + C\langle \tau^{3} + \tau h^4, e^{1} + e^0\rangle. \end{align}$

(4.7)

Sum up Eqs (4.6) and (4.7), and apply Lemmas 3.2 and 4.2, it deduces that

$\begin{align*} \|e^J\|^2 \leq& \frac{\tau^2}{2}\sum\limits_{n = 1}^{J-1}\sum\limits_{k = 0}^{n}(\omega_{k+1} + \omega_k)\langle \mathcal{H}\left(g(u^{n-k})-g(U^{n-k})\right), e^{n+1} + e^n\rangle\nonumber\\ &+\frac{\tau^2 \omega_0}{2}\sum\limits_{n = 1}^{J-1}\langle \mathcal{H} \left(g(2u^n -u^{n-1})- g(2U^n -U^{n-1})\right), e^{n+1} + e^n \rangle\nonumber\\ &+\frac{\tau^2\omega_{1}}{2}\langle \mathcal{H}\left(g(u^{0})-g(U^{0})\right), e^{n+1} + e^n\rangle \nonumber\\ &+\tau^2 \omega_0\langle \mathcal{H} \left(g(u^0)- g(U^0)\right), e^{n+1} + e^n \rangle +C \sum\limits_{n = 1}^{J-1}\langle \tau^{3} + \tau h^4, e^{n+1} + e^n\rangle. \end{align*}$

According to the same technique as for dealing with (3.9), we can achieve

$\begin{align*} \|e^P\| \leq C(\tau^{2} + h^4), \end{align*}$

thus completes the proof.

Theorem 4.5. Let $\{U_i^n|\mathit{\text{}}0\leq i \leq M, 0\leq n \leq N\}$ be the numerical solution of Scheme (4.4) and (4.5) for Problem (2.1). Then for $1 \leq K \leq N$ , it holds

$\begin{align*} \|U^K\| \leq C\left(\max\limits_{0\leq n \leq N} \|g(U^n)\| +\max\limits_{0\leq n \leq N-1} \|F^{n+\frac{1}{2}}\| \right). \end{align*}$

5. Numerical experiments

In this section, we carry out numerical experiments to verify the theoretical results and demonstrate the performance of our new schemes. All of the computations are performed by using a MATLAB on a computer with Intel(R) Core(TM) i5-8265U CPU 1.60GHz 1.80GHz and 8G RAM.

Example 5.1. Consider the following problem with exact solution $u(x, t) = t^{2+\alpha} \sin^2(\pi x)$

$\begin{align*} \frac{\partial^2u(x, t)}{\partial t^2} + {}^{C}_{0}D_{t}^{\alpha}u(x, t) +\frac{\partial^{4}u(x, t)}{\partial x^4} = \frac{\partial^2 u(x, t)}{\partial x^2}+ f(x, t) + g(u), \end{align*}$

where $T = 1$ , $0 < x < 1$ , $0 < t \leq T$ , and $1 < \alpha < 2$ . The nonlinear function $g(u) = u^2$ and $f(x, t)$ is

$\begin{align*} f(x, t) = (2&+\alpha)(1+\alpha) t^{\alpha}\sin^2(\pi x) + \frac{\Gamma(3+\alpha)}{2}t^{2}\sin^2(\pi x) - 8\pi^4 t^{2+\alpha}\cos(2\pi x)\\ -& 2\pi^2 t^{2+\alpha}\cos(2\pi x) - t^{2(2+\alpha)}\sin^{4}(\pi x). \end{align*}$

It is clear that $u(x, t)$ satisfies all smoothness conditions required by Theorems 3.4 and 4.4, so that both of our schemes can be applied in this example. In Figures 1 and 2, we compare the exact solution with the numerical solution of finite difference Scheme (3.5) and (3.6) and compact finite difference Scheme (4.4) and (4.5). We easily see that the exact solution can be well approximated by the numerical solutions of our schemes.

Figure 1. The comparison of numerical solution of Scheme (3.5) and (3.6) with the exact solution for

$\tau = h = 0.01$ and

$\alpha = 1.6$ .

DownLoad: Full-Size Img PowerPoint

Figure 2. The comparison of numerical solution of the compact finite difference Scheme (4.4) and (4.5) with the exact solution for

$\tau = h = 0.01$ and

$\alpha = 1.6$ .

DownLoad: Full-Size Img PowerPoint

First, we in , and show that the errors, time and space convergence order $\thickapprox 2$ and CPU times (second) of the finite difference Scheme (3.5) and (3.6) for $\alpha = 1.25, 1.5, 1.75$ . The average CPU time, expressed as the mean time (mean) for $\alpha = 1.25, 1.5, 1.75$ . Specifically, tests the case that when $\tau = h$ . In , we set $h = 0.001$ , a value small enough such that the spatial discretization errors are negligible as compared with the temporal errors, and choose different time step size. In , we set $\tau = 0.001$ , a value small enough such that the temporal discretization errors are negligible as compared with the spatial errors, and choose different space step size. From all scenarios above, we conclude that the temporal and spatial convergence order is 2. It verifies Theorem 3.4.

Table 1. The errors, CPU times (second) for different

$\alpha$ , and numerical convergence orders of Scheme (3.5) and (3.6) for different

$\tau = h$ .

$\tau=h$	$\alpha=1.25$		$\alpha=1.5$		$\alpha = 1.75$		CPU time
$\tau=h$	error	order	error	order	error	order	mean
$1/5$	$6.6627\times 10^{-2}$		$7.8031\times 10^{-2}$		$8.9815\times 10^{-2}$		$0.0896$
$1/10$	$1.8412\times 10^{-2}$	$1.8555$	$2.1839\times 10^{-2}$	$1.8371$	$2.5456\times 10^{-2}$	$1.8190$	$0.0973$
$1/20$	$4.8132\times 10^{-3}$	$1.9355$	$5.7273\times 10^{-3}$	$1.9310$	$6.6917\times 10^{-3}$	$1.9275$	$0.0994$
$1/40$	$1.2137\times 10^{-3}$	$1.9876$	$1.4621\times 10^{-3}$	$1.9698$	$1.7210\times 10^{-3}$	$1.9591$	$0.1359$

| Show Table

DownLoad: CSV

Table 2. The errors, CPU times (second) for different

$\alpha$ , and temporal numerical convergence orders of Scheme (3.5) and (3.6) for

$h = 0.001$ and different

$\tau$ .

$\tau$	$\alpha=1.25$		$\alpha=1.5$		$\alpha = 1.75$		CPU time
$\tau$	error	order	error	order	error	order	mean
$1/5$	$7.0844\times 10^{-2}$		$8.2130\times 10^{-2}$		$9.3783\times 10^{-2}$		$0.5852$
$1/10$	$1.9012\times 10^{-2}$	$1.8977$	$2.2432\times 10^{-2}$	$1.8724$	$2.6040\times 10^{-2}$	$1.8486$	$1.0501$
$1/20$	$4.9405\times 10^{-3}$	$1.9442$	$5.8537\times 10^{-3}$	$1.9381$	$6.8169\times 10^{-3}$	$1.9335$	$2.4071$
$1/40$	$1.2435\times 10^{-3}$	$1.9903$	$1.4917\times 10^{-3}$	$1.9724$	$1.7504\times 10^{-3}$	$1.9615$	$6.7799$

| Show Table

DownLoad: CSV

Table 3. The errors, CPU times (second) for different

$\alpha$ , and spatial numerical convergence orders of Scheme (3.5) and (3.6) for

$\tau = 0.001$ and different

$h$ .

$h$	$\alpha=1.25$		$\alpha=1.5$		$\alpha = 1.75$		CPU time
$h$	error	order	error	order	error	order	mean
$1/5$	$4.7813\times 10^{-3}$		$4.7510\times 10^{-3}$		$4.7111\times 10^{-3}$		$2.2382$
$1/10$	$6.1943\times 10^{-4}$	$2.9484$	$6.1518\times 10^{-4}$	$2.9492$	$6.0963\times 10^{-4}$	$2.9501$	$2.2952$
$1/20$	$1.2773\times 10^{-4}$	$2.2778$	$1.2654\times 10^{-4}$	$2.2815$	$1.2503\times 10^{-4}$	$2.2857$	$2.5410$
$1/40$	$2.8950\times 10^{-5}$	$2.1415$	$2.8363\times 10^{-5}$	$2.1575$	$2.7665\times 10^{-5}$	$2.1761$	$3.4293$

| Show Table

DownLoad: CSV

On the other hand, we check the numerical convergence orders and CPU times (second) in time and space of the compact finite difference Scheme (4.4) and (4.5) for $\alpha = 1.25, 1.5, 1.75$ in and , respectively. The average CPU time, expressed as the mean time (mean) for $\alpha = 1.25, 1.5, 1.75$ . As expected, the numerical results reflect that the compact finite difference has a convergence order of $2$ and $4$ in time and space, respectively, which verifies our Theorem 4.4.

Table 4. The errors, CPU times (second) for different

$\alpha$ , and temporal convergence orders of Scheme (4.4) and (4.5) for

$h = 0.001$ and different

$\tau$ .

$\tau$	$\alpha=1.25$		$\alpha=1.5$		$\alpha = 1.75$		CPU time
$\tau$	error	order	error	order	error	order	mean
$1/5$	$7.0844\times 10^{-2}$		$8.2129\times 10^{-2}$		$9.3783\times 10^{-2}$		$0.9501$
$1/10$	$1.9012\times 10^{-2}$	$1.8978$	$2.2432\times 10^{-2}$	$1.8724$	$2.6040\times 10^{-2}$	$1.8486$	$2.3622$
$1/20$	$4.9407\times 10^{-3}$	$1.9441$	$5.8538\times 10^{-3}$	$1.9381$	$6.8169\times 10^{-3}$	$1.9335$	$7.6793$
$1/40$	$1.2436\times 10^{-3}$	$1.9901$	$1.4919\times 10^{-3}$	$1.9723$	$1.7506\times 10^{-3}$	$1.9612$	$28.9326$

| Show Table

DownLoad: CSV

Table 5. The errors, CPU times (second) for different

$\alpha$ , and spatial convergence orders of Scheme (4.4) and (4.5) for

$\tau = 0.0005$ and different

$h$ .

$h$	$\alpha=1.25$		$\alpha=1.5$		$\alpha = 1.75$		CPU time
$h$	error	order	error	order	error	order	mean
$1/5$	$3.8110\times 10^{-3}$		$3.7871\times 10^{-3}$		$3.7555\times 10^{-3}$		$12.5566$
$1/10$	$2.5308\times 10^{-4}$	$3.9125$	$2.5141\times 10^{-4}$	$3.9130$	$2.4922\times 10^{-4}$	$3.9135$	$14.0490$
$1/20$	$2.2087\times 10^{-5}$	$3.5183$	$2.1851\times 10^{-5}$	$3.5243$	$2.1557\times 10^{-5}$	$3.5312$	$18.1726$
$1/40$	$1.8261\times 10^{-6}$	$3.5964$	$1.7163\times 10^{-6}$	$3.6703$	$1.5904\times 10^{-6}$	$3.7607$	$43.9104$

| Show Table

DownLoad: CSV

6. Conclusions

We in this paper constructed two linearized finite difference schemes for time fractional nonlinear diffusion-wave equations with the space fourth-order derivative. The equations were transformed into equivalent partial integro-differential equations. Then, the Crank-Nicolson technique, the midpoint formula, the weighted and shifted Gr $\ddot{\text{u}}$ nwald difference formula, the second order convolution formula, the classical central difference formula, the fourth-order approximation and the compact difference technique were applied to construct the two proposed schemes. The finite difference Scheme (3.5) and (3.6) has the accuracy $O(\tau^2 + h^2)$ . The compact finite difference Scheme (4.4) and (4.5) has the accuracy $O(\tau^2 + h^4)$ . It should be mentioned that our schemes require the exact solution $u(\cdot, t) \in C^3([0, T])$ , while it requires $u(\cdot, t)\in C^4([0, T])$ if one discretizes Eq (1.1) directly to get the second order accuracy in time. Theoretically, the convergence and the unconditional stability of the two proposed schemes are proved and discussed. All of the numerical experiments can support our theoretical results.

Acknowledgments

This research is supported by Natural Science Foundation of Jiangsu Province of China (Grant No. BK20201427), and by National Natural Science Foundation of China (Grant Nos. 11701502 and 11871065).

Conflict of interest

The authors declare that they have no competing interests.

References

[1]	Jorgenson MT (2013) Thermokarst Terrains, In: Shroder J, Giardino R, and Harbor J, Eds., Treatise on Geomorphology, San Diego, Academic Press, 313–324. https://doi.org/10.1016/B978-0-12-374739-6.00215-3
[2]	Harris C, Arenson LU, Christiansen HH, et al. (2009) Permafrost and climate in Europe: Monitoring and modelling thermal, geomorphological and geotechnical responses. Earth Sci Rev 92: 117–171. https://doi.org/10.1016/j.earscirev.2008.12.002 doi: 10.1016/j.earscirev.2008.12.002
[3]	Jorgenson MT, Romanovsky V, Harden J, et al. (2010) Resilience and vulnerability of Permafrost to climate change. Can J For Res 40: 1219–1236. https://doi.org/10.1139/X10-060 doi: 10.1139/X10-060
[4]	Osterkamp TE, Romanovsky VE (1999) Evidence for warming and thawing of discontinuous permafrost in Alaska. Permafrost Periglac 10: 17–37. https://doi.org/10.1002/(SICI)1099-1530(199901/03)10:1<17:AID-PPP303>3.0.CO;2-4 doi: 10.1002/(SICI)1099-1530(199901/03)10:1<17:AID-PPP303>3.0.CO;2-4
[5]	Romanovsky VE, Smith SL, Christiansen HH (2010) Permafrost thermal state in the polar northern hemisphere during the international polar year 2007–2009: A synthesis. Permafrost Periglac 21: 106–116. https://doi.org/10.1002/ppp.689 doi: 10.1002/ppp.689
[6]	Romanovsky VE, Osterkamp TE (2000) Effects of unfrozen water on heat and mass transport processes in the active layer and permafrost. Permafrost Periglac 11: 219–239. https://doi.org/10.1002/1099-1530(200007/09)11:3<219:AID-PPP352>3.0.CO;2-7 doi: 10.1002/1099-1530(200007/09)11:3<219:AID-PPP352>3.0.CO;2-7
[7]	Shur YL, Jorgenson MT (2007) Patterns of permafrost formation and degradation in relation to climate and ecosystems. Permafrost Periglac 18: 7–19. https://doi.org/10.1002/ppp.582 doi: 10.1002/ppp.582
[8]	Osterkamp TE, Jorgenson MT, Schuur EAG, et al. (2009) Physical and ecological changes associated with warming permafrost and thermokarst in Interior Alaska. Permafrost Periglac 20: 235–256. https://doi.org/10.1002/ppp.656 doi: 10.1002/ppp.656
[9]	Olenchenko VV, Gagarin LA, Khristoforov II, et al. (2017) The structure of a site with thermo-suffosion processes within bestyakh terrace of the lena river, according to geophysical data. Kriosfera Zemli 21: 16–26.
[10]	Briggs MA, Campbell S, Nolan J, et al. (2017) Surface Geophysical Methods for Characterising Frozen Ground in Transitional Permafrost Landscapes. Permafrost Periglac 28: 52–65. https://doi.org/10.1002/ppp.1893 doi: 10.1002/ppp.1893
[11]	Shiklomanov NI, Nelson FE (1999) Analytic representation of the active layer thickness field, Kuparuk River Basin, Alaska. Ecol Modell 123: 105–125. https://doi.org/10.1016/S0304-3800(99)00127-1 doi: 10.1016/S0304-3800(99)00127-1
[12]	Shur Y, Hinkel KM, Nelson FE (2005) The transient layer: Implications for geocryology and climate-change science. Permafrost Periglac 16: 5–17. https://doi.org/10.1002/ppp.518 doi: 10.1002/ppp.518
[13]	Kneisel C, Hauck C, Fortier R, et al. (2008) Advances in geophysical methods for permafrost investigations. Permafrost Periglac 19: 157–178. https://doi.org/10.1002/ppp.616 doi: 10.1002/ppp.616
[14]	Conaway CH, Johnson CD, Lorenson TD, et al. (2020) Permafrost Mapping with Electrical Resistivity Tomography: A Case Study in Two Wetland Systems in Interior Alaska. J Environ Eng Geophys 25: 199–209. https://doi.org/10.2113/JEEG19-091 doi: 10.2113/JEEG19-091
[15]	Douglas TA, Jorgenson MT, Kanevskiy MZ, et al. (2008) Permafrost Dynamics at the Fairbanks Permafrost Experimental Station Near Fairbanks, Alaska. University of Alaska, Fairbanks, 373–377.
[16]	Fortier R, LeBlanc AM, Allard M, et al. (2008) Internal structure and conditions of permafrost mounds at Umiujaq in Nunawik, Canada, inferred from field investigation and electrical resistivity tomography. Can J Earth Sci 45: 367–387. https://doi.org/10.1139/E08-004 doi: 10.1139/E08-004
[17]	Swarzenski PW, Johnson CD, Lorenson TD, et al. (2016) Seasonal Electrical Resistivity Surveys of a Coastal Bluff, Barter Island, North Slope Alaska. J Environ Eng Geophys 21: 37–42. https://doi.org/10.2113/JEEG21.1.37 doi: 10.2113/JEEG21.1.37
[18]	Voytek EB, Rushlow CR, Godsey SE, et al. (2016) Identifying hydrologic flowpaths on arctic hillslopes using electrical resistivity and self potential. Geophysics 81: WA225–WA232. https://doi.org/10.1190/GEO2015-0172.1 doi: 10.1190/GEO2015-0172.1
[19]	Parshley L (2022) Climate and Science. As permafrost thaws, the ground beneath Alaska is collapsing. Available from: https://grist.org/science/alaska-permafrost-thawing-ice-climate-change/.
[20]	Péwé TL (1975) Quaternary Geology of Alaska, Washington DC: US Government Printing Office.
[21]	Péwé TL (1954) Effect of Permafrost on Cultivated Fields, Fairbanks Area, Alaska, US Government Printing Office, 315–351.
[22]	Prindle LM, Smith PS, Katz FJ (1913) A Geologic Reconnaissance of The Fairbanks Quadrangle, Alaska Detailed Description of the Fairbanks District Account of Lode Mining Near Fairbanks, Washington: Government Printing Office.
[23]	Kanevskiy, M, Ping, C, Shur, Y, et al. (2011). Permafrost of Northern Alaska. The Proceedings of the Twenty-First (2011). Int Offshore Polar Eng Conf, 1179–1186.
[24]	Alaska Department of Natural Resources. "DNR Well Log—WELTS, " Well Logs at UAF. Available from: https://dnr.alaska.gov/welts.
[25]	Douglas TA, Jorgenson MT, Brown DRN, et al. (2016) Degrading permafrost mapped with electrical resistivity tomography, airborne imagery and LiDAR, and seasonal thaw measurements. Geophysics 81: WA71–WA85. https://doi.org/10.1190/geo2015-0149.1 doi: 10.1190/geo2015-0149.1
[26]	Jorgenson MT, Racine CH, Walters JC, et al. (2001) Permafrost Degradation and Ecological Changes Associated with a Warming Climate in Central Alaska. Clim Change 48: 551–579.
[27]	Frolov AD, Zykov DY, Snegirev AM, et al. (1998) Principal problems, progress, and directions of Geophysical investigations in permafrost regions. 7th International Conference on Permafrost. Université Laval, Canada, 305–311.
[28]	Asare A, Appiah-Adjei EK, Owusu-Nimo F, et al. (2022) Lateral and vertical mapping of salinity along the coast of Ghana using Electrical Resistivity Tomography: The case of Central Region. Results Geophys Sci 12: 100048. https://doi.org/10.1016/j.ringps.2022.100048 doi: 10.1016/j.ringps.2022.100048
[29]	Pazzi V, Morelli S, Fanti R (2019) A Review of the Advantages and Limitations of Geophysical Investigations in Landslide Studies. Int J Geophys 2019: 1–27. https://doi.org/10.1155/2019/2983087 doi: 10.1155/2019/2983087
[30]	Hoekstra P, Sellmann PV, Delaney A (1975) Ground and Airborne Resistivity Surveys of Permafrost Near Fairbanks, Alaska. Geophysics 40: 641–656. https://doi.org/10.1190/1.1440555 doi: 10.1190/1.1440555
[31]	Hauck C, Kneisel C (2008) Applied geophysics in periglacial environments. Cambridge: Cambridge University Press.
[32]	Fortier R, Allard M, Seguin MK (1994) Effect of physical properties of frozen ground on electrical resistivity logging. Cold Reg Sci Technol 22: 361–384. https://doi.org/10.1016/0165-232X(94)90021-3 doi: 10.1016/0165-232X(94)90021-3
[33]	Mathys T, Hilbich C, Arenson LU, et al. (2022) Towards accurate quantification of ice content in permafrost of the Central Andes—Part 2: An upscaling strategy of geophysical measurements to the catchment scale at two study sites. Cryosphere 16: 2595–2615. https://doi.org/10.5194/tc-16-2595-2022 doi: 10.5194/tc-16-2595-2022
[34]	Loke MH (2000) Electrical imaging surveys for environmental and engineering studies. A practical guide to 2. Available from: https://www.abem.se.
[35]	Gharibi M, Bentley LR (2005) Resolution of 3-D electrical resistivity images from inversions of 2-D orthogonal lines. J Environ Eng Geophys 10: 339–349. https://doi.org/10.2113/JEEG10.4.339 doi: 10.2113/JEEG10.4.339
[36]	Kidanu S, Varnavina A, Anderson N, et al. (2020) Pseudo-3D-electrical resistivity tomography imaging of subsurface structure of a sinkhole—A case study in Greene County, Missouri. AIMS Geosci 6: 54–70. https://doi.org/10.3934/geosci.2020005 doi: 10.3934/geosci.2020005
[37]	Golden Software (2019) Surfer ® Powerful contouring, gridding & surface mapping system Full User's Guide. Colorado. Available from: https://www.GoldenSoftware.com.
[38]	Loke MH, Chambers JE, Rucker DF, et al. (2013) Recent developments in the direct-current geoelectrical imaging method. J Appl Geophys 95: 135–156. https://doi.org/10.1016/j.jappgeo.2013.02.017 doi: 10.1016/j.jappgeo.2013.02.017
[39]	Degroot-Hedlin C, Constable S (1990) Occam' s inversion to generate smooth, two-dimensional models from magnetotelluric data. Geophysics 55: 1613–1624. https://doi.org/10.1190/1.1442813 doi: 10.1190/1.1442813
[40]	Loke MH, Barker RD (1996) Rapid least-squares inversion of apparent resistivity pseudosections by a quasi-Newton method^l. Geophys Prospect 44: 131–152. https://doi.org/10.1111/j.1365-2478.1996.tb00142.x doi: 10.1111/j.1365-2478.1996.tb00142.x
[41]	Loke MH (2004) Tutorial: 2-D and 3-D electrical imaging surveys. Available from: https://www.geotomosoft.com.
[42]	Hauck C, Mühll DV (2003) Inversion and interpretation of two-dimensional geoelectrical measurements for detecting permafrost in mountainous regions. Permafrost Periglac 14: 305–318. https://doi.org/10.1002/ppp.462 doi: 10.1002/ppp.462
[43]	Loke MH (2021) Tutorial: 2-D and 3-D electrical imaging surveys. Geotomo Software, 1–232. Available from: https://www.geotomosoft.com.
[44]	Oldenburg DW, Li Y (1999) Estimating depth of investigation in dc resistivity and IP surveys. Geophysics 64: 403–41. https://doi.org/10.1190/1.1444545 doi: 10.1190/1.1444545
[45]	Roy A, Apparao A (1971) Depth of Investigation in Direct Current Methods. Geophysics 36: 943–959. https://doi.org/10.1190/1.1440226 doi: 10.1190/1.1440226
[46]	Edwards LS (1977) A Modified Pseudosection for Resistivity and IP. Geophysics 42: 1020–1036. https://doi.org/10.1190/1.1440762 doi: 10.1190/1.1440762
[47]	Szalai S, Novák A, Szarka L (2009) Depth of investigation and vertical resolution of surface geoelectric arrays. J Environ Eng Geophys 14: 15–23. https://doi.org/10.2113/JEEG14.1.15 doi: 10.2113/JEEG14.1.15
[48]	Marescot L, Loke MH, Chapellier D, et al. (2003) Assessing reliability of 2D resistivity imaging in mountain permafrost studies using the depth of investigation index method. Near Surf Geophys 1: 57–67. https://doi.org/10.3997/1873-0604.2002007 doi: 10.3997/1873-0604.2002007
[49]	Hilbich C, Fuss C, Hauck C (2011) Automated time-lapse ERT for improved process analysis and monitoring of frozen ground. Permafrost Periglac 22: 306–319. https://doi.org/10.1002/ppp.732 doi: 10.1002/ppp.732
[50]	Loke MH, Dahlin T (2002) A comparison of the Gauss-Newton and quasi-Newton methods in resistivity imaging inversion. J Appl Geophys 49: 149–162. https://doi.org/10.1016/S0926-9851(01)00106-9 doi: 10.1016/S0926-9851(01)00106-9
[51]	Liljedahl AK, Boike J, Daanen R, et al. (2016) Pan-Arctic ice-wedge degradation in warming permafrost and its influence on tundra hydrology. Nature Geosci 9: 312–318. https://doi.org/10.1038/ngeo2674 doi: 10.1038/ngeo2674
[52]	Farquharson L, Anthony KW, Bigelow N, et al. (2016) Facies analysis of yedoma thermokarst lakes on the northern Seward Peninsula, Alaska. Sediment Geol 340: 25–37. https://doi.org/10.1016/j.sedgeo.2016.01.002 doi: 10.1016/j.sedgeo.2016.01.002
[53]	Douglas TA, Jones MC, Heimstra AC (2014) Sources and Sinks of Carbon in Boreal Ecosystems of Interior Alaska: A review. Elementa 2: 000032. https://doi.org/10.12952/journal.elementa.000032 doi: 10.12952/journal.elementa.000032

This article has been cited by:

1.	Emadidin Gahalla Mohmed Elmahdi, Jianfei Huang, A linearized finite difference scheme for time–space fractional nonlinear diffusion-wave equations with initial singularity, 2022, 0, 1565-1339, 10.1515/ijnsns-2021-0388
2.	Emadidin Gahalla Mohmed Elmahdi, Jianfei Huang, EFFICIENT NUMERICAL SOLUTION OF TWO-DIMENSIONAL TIME-SPACE FRACTIONAL NONLINEAR DIFFUSION-WAVE EQUATIONS WITH INITIAL SINGULARITY, 2022, 12, 2156-907X, 831, 10.11948/20210444
3.	Chaeyoung Lee, Seokjun Ham, Youngjin Hwang, Soobin Kwak, Junseok Kim, An explicit fourth-order accurate compact method for the Allen-Cahn equation, 2024, 9, 2473-6988, 735, 10.3934/math.2024038
4.	Emadidin Gahalla Mohmed Elmahdi, Yang Yi, Jianfei Huang, Two linearized difference schemes on graded meshes for the time-space fractional nonlinear diffusion-wave equation with an initial singularity, 2025, 100, 0031-8949, 015215, 10.1088/1402-4896/ad95c4

Reader Comments

Your name:*

Email:*
© 2024 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)