An extension of Herstein's theorem on Banach algebra

Abu Zaid Ansari; Suad Alrehaili; Faiza Shujat; Abu Zaid Ansari; Suad Alrehaili; Faiza Shujat

doi:10.3934/math.2024201

AIMS Mathematics

2024, Volume 9, Issue 2: 4109-4117. doi: 10.3934/math.2024201

Previous Article Next Article

Research article Special Issues

An extension of Herstein's theorem on Banach algebra

1.
Department of Mathematics, Faculty of Science, Islamic University of Madinah, Saudi Arabia
2.
Department of Mathematics, Faculty of Science, Taibah University, Madinah, Saudi Arabia

Received: 04 October 2023 Revised: 24 November 2023 Accepted: 02 January 2024 Published: 12 January 2024
MSC : 16N60, 16B99, 16W25

Let $\mathcal{A}$ be a $(p+q)!$ -torsion free semiprime ring. We proved that if $\mathcal{H}, \mathcal{D} : \mathcal{A}\to \mathcal{A}$ are two additive mappings fulfilling the algebraic identity $2\mathcal{H}(a^{p+q}) = \mathcal{H}(a^p) a^q+ a^p \mathcal{D}(a^q)+\mathcal{H}(a^q) a^p+ a^q \mathcal{D}(a^p)$ for all $a\in \mathcal{A}$ , then $\mathcal{H}$ is a generalized derivation with $\mathcal{D}$ as an associated derivation on $\mathcal{A}$ . In addition to that, it is also proved in this article that $\mathcal{H}_1$ is a generalized left derivation associated with a left derivation $\delta$ on $\mathcal{A}$ if they fulfilled the algebraic identity $2\mathcal{H}_1(a^{p+q}) = a^p \mathcal{H}_1(a^q)+ a^q \delta(a^p)+a^q \mathcal{H}_1(a^p)+ a^p \delta(a^q)$ for all $a \in \mathcal{A}$ . Further, the legitimacy of these hypotheses is eventually demonstrated by examples and at last, an application of Banach algebra is presented.

Keywords:

Citation: Abu Zaid Ansari, Suad Alrehaili, Faiza Shujat. An extension of Herstein's theorem on Banach algebra[J]. AIMS Mathematics, 2024, 9(2): 4109-4117. doi: 10.3934/math.2024201

Related Papers:

[1]	Yunzhang Zhang, Xinghui Yong . Analysis of the linearly extrapolated BDF2 fully discrete Modular Grad-div stabilization method for the micropolar Navier-Stokes equations. AIMS Mathematics, 2024, 9(6): 15724-15747. doi: 10.3934/math.2024759
[2]	Eunjung Lee, Dojin Kim . Stability analysis of the implicit finite difference schemes for nonlinear Schrödinger equation. AIMS Mathematics, 2022, 7(9): 16349-16365. doi: 10.3934/math.2022893
[3]	Zuliang Lu, Ruixiang Xu, Chunjuan Hou, Lu Xing . A priori error estimates of finite volume element method for bilinear parabolic optimal control problem. AIMS Mathematics, 2023, 8(8): 19374-19390. doi: 10.3934/math.2023988
[4]	Danuruj Songsanga, Parinya Sa Ngiamsunthorn . Single-step and multi-step methods for Caputo fractional-order differential equations with arbitrary kernels. AIMS Mathematics, 2022, 7(8): 15002-15028. doi: 10.3934/math.2022822
[5]	Rubayyi T. Alqahtani, Jean C. Ntonga, Eric Ngondiep . Stability analysis and convergence rate of a two-step predictor-corrector approach for shallow water equations with source terms. AIMS Mathematics, 2023, 8(4): 9265-9289. doi: 10.3934/math.2023465
[6]	Philsu Kim, Seongook Heo, Dojin Kim . A high-order convergence analysis for semi-Lagrangian scheme of the Burgers' equation. AIMS Mathematics, 2023, 8(5): 11270-11296. doi: 10.3934/math.2023571
[7]	Danxia Wang, Ni Miao, Jing Liu . A second-order numerical scheme for the Ericksen-Leslie equation. AIMS Mathematics, 2022, 7(9): 15834-15853. doi: 10.3934/math.2022867
[8]	Narcisse Batangouna . A robust family of exponential attractors for a time semi-discretization of the Ginzburg-Landau equation. AIMS Mathematics, 2022, 7(1): 1399-1415. doi: 10.3934/math.2022082
[9]	Ruby, Vembu Shanthi, Higinio Ramos . A numerical approach to approximate the solution of a quasilinear singularly perturbed parabolic convection diffusion problem having a non-smooth source term. AIMS Mathematics, 2025, 10(3): 6827-6852. doi: 10.3934/math.2025313
[10]	Sara S. Alzaid, Pawan Kumar Shaw, Sunil Kumar . A study of Ralston's cubic convergence with the application of population growth model. AIMS Mathematics, 2022, 7(6): 11320-11344. doi: 10.3934/math.2022632

Abstract

1. Introduction

The grad-div stabilization method was first introduced for the finite element approximation of the incompressible Navier-Stokes equations by adding the grad-div stabilization term to the momentum equation in ^[12]. For the most of Lagrange finite element pairs, such as the Taylor-Hood element and MINI element, the grad-div stabilization term is nonzero, and there is no point-wise mass conservation. Thus, this stabilization term is used to improve the mass conservation and to reduce the velocity error caused by the pressure error. The influence of the grad-div stabilization term on the accuracy of the numerical solution was studied in ^[40] for the steady Stokes problem.

It is well known that the standard Galerkin finite element methods are not suitable for the numerical simulation of the incompressible flow problems in the case of the high Reynolds number or the small viscosity. Different stabilized methods have been developed to overcome the instabilities of numerical simulation, such as Galerkin least square methods in ^[13,25], the residual-free bubbles methods in ^[14,15], the large eddy simulation methods in ^[45], the sub-grid scale methods in ^[23,31], the variational multiscale methods in ^{[26,27,34,49]} and the streamline upwind Petrov-Galerkin (SUPG) method ^[5]. We note that the grad-div stabilization method is also a powerful tool for the high Reynolds flow problems. There have an amount of works for this issue. For example, the SUPG and grad-div stabilization methods were studied for the generalized Oseen problem in ^[37], where it was concluded that the SUPG method is less important for the inf-sup stable pair of velocity and pressure due to the choice of stabilization parameter. Moreover, if one considers only the grad-div stabilization method, it was acknowledged that the grad-div method can lead to stabilized results. The grad-div method for the rotation form of the Navier-Stokes equations can be found in ^[32], where numerical results showed that the difference between the skew-symmetric and the rotation form of the nonlinearity is from the increased error in the Bernoulli pressure, which results in the increasing of velocity error. The use of the grad-div stabilization term ameliorates the effect and reduces the velocity error. For the time-dependent Oseen problem, the semi-discrete and fully discrete schemes based on the grad-div stabilization method were analyzed in ^[17]. Under the assumption of sufficiently smooth solutions and based on the specific Stokes projection, the authors proved the optimal error estimates for velocity and pressure with constants that are independent of the viscosity. Subsequently, this observation was extended to the time-dependent Navier-Stokes equations in ^[18], where the constants in the optimal error bounds do not depend on the inverse power of viscosity. Thus, the results in ^[17] and ^[18] gave a theoretical confirmation about stable simulations of the high Reynolds number flows by using the grad-div stabilization method. Recently, the analysis technique in ^[18] was extended to the second-order backward differentiation formula (BDF2) scheme with variable time-step size of the Navier-Stokes equations in ^[21]. A review study about error analysis of the grad-div stabilization methods can be found in ^[20].

Since there are some advantages of the grad-div stabilization method in numerical simulations of the incompressible flows, it has been applied to some related and coupled systems, such as the magnetohydrodynamic (MHD) equations ^[11] and the time-dependent dual-porosity-Stokes equations ^[33]. In this paper, we will consider the following time-dependent penetrative convection equations in the nondimensional form:

$\begin{align} \frac {\partial {{\textbf{u}}}}{\partial t}-\nu\Delta{{\textbf{u}}}+({{\textbf{u}}}\cdot\nabla){{\textbf{u}}}+\nabla p-(\gamma_{1}\theta+\gamma_{2}\theta^{2}){\textbf{i}}_{3} = {{\textbf{f}}},\quad &{\mbox{for}}\quad ({\textbf{x}},t)\in\ \Omega\times(0,T], \end{align}$

(1.1)

$\begin{align} \nabla\cdot{{\textbf{u}}} = 0,\quad &{\mbox{for}}\quad ({\textbf{x}},t)\in\ \Omega\times(0,T] , \end{align}$

(1.2)

$\begin{align} \frac {\partial \theta}{\partial t}-\kappa\Delta\theta+({{\textbf{u}}}\cdot\nabla)\theta = g,\quad &{\mbox{for}}\quad({\textbf{x}},t)\in \ \Omega\times(0,T], \end{align}$

(1.3)

$\begin{align} {{\textbf{u}}} = {\textbf{0}},\quad\theta = 0,\quad &{\mbox{for}}\quad({\textbf{x}},t)\in \ \partial\Omega\times(0,T], \end{align}$

(1.4)

$\begin{align} {{\textbf{u}}}({\textbf{x}}, 0) = {{\textbf{u}}}_{0}({\textbf{x}}),\quad\theta({\textbf{x}},0) = \theta_{0}({\textbf{x}}),\quad&{\mbox{for}}\quad{\textbf{x}}\in\Omega, \end{align}$

(1.5)

where $\Omega$ is a bounded and convex polyhedral domain in $\mathbb R^3$ with boundary $\partial\Omega$ , and $(0, T]$ is a finite time interval. In the above system, the unknown $({{\textbf{u}}}, p, \theta)$ denotes the velocity of the fluid, the pressure and temperature, respectively. $\nu > 0$ and $\kappa > 0$ represent the viscosity coefficient and thermal conductivity parameter. ${{\textbf{f}}}$ and $g$ are two given functions. ${\textbf{i}}_{3}$ is the unit basis vector given by ${\textbf{i}}_{3} = (0, 0, 1)^{\mathrm{T}}$ , and $(\gamma_{1}\theta+\gamma_{2}\theta^{2}){\textbf{i}}_{3}$ in (1.1) denotes the buoyancy term. The constraint $\nabla\cdot{{\textbf{u}}} = 0$ represents the incompressibility of the fluid.

The penetrative convection equations (1.1)–(1.3) used to describe the motion in which convection in a thermally unstable region extends into an adjacent stable region. It has a wide range of applications in the fields of atmospheric dynamics, atmospheric fronts, katabatic winds, dense gas dispersion, natural ventilation, cooling of electronic equipment ^{[29,36,39,42]}. The stabilities of penetrative convection equations have been studied in ^{[7,8,16,41,46]}, and numerical simulations have been reported in ^[3,38]. Convergence and error estimates of finite element fully discrete schemes were studied in ^[44] and ^[6], where the BDF2 and Crank-Nicolson schemes were used, respectively. However, the high Reynolds number was not considered in ^[6,44]. By considering a linear Boussinesq approximation, i.e., $\gamma_2 = 0$ in (1.1), a grad-div stabilization finite element method was studied in ^[19]. In addition, a rigorous error analysis was given and error estimates were derived in ^[19], where the constants in error bounds are independent of inverse powers of the viscosity and thermal conductivity coefficients.

In this paper, based on the grad-div stabilization method, we consider and study two finite element schemes by using the first-order backward Euler and BDF2 formulas to discrete the time derivatives, respectively. The proposed schemes are both linearized semi-implicit schemes where we use the implicit-explicit method and the standard extrapolation formula to deal with nonlinear terms. Moreover, these schemes are unconditionally stable without any condition of the time step size and mesh size. Following the analysis techniques developed in ^[18,19], we derive error bounds of the velocity and temperature in which the constants are independent of negative powers of the viscosity and thermal conductivity coefficients.

The rest of this paper is organized as follows. In next section, we state some notations and recall some known results in finite element theory. In sections three and four, the first-order Euler and second-order BDF2 schemes with the grad-div stabilization are proposed, respectively. Unconditional stability of numerical schemes and error analysis are also presented in these sections. In section five, we give numerical results to support the theory analysis and numerical stability of the grad-div stabilization method for the high Reynolds flows.

2. Materials and methods

For $k\in\mathbb N^+$ and $1\leq p\leq +\infty$ , let $L^{p}(\Omega)$ and $W^{k, p}(\Omega)$ denote the standard Lebesgue space and Sobolev space, respectively. When $p = 2$ , $W^{k, 2}(\Omega)$ is the Hilber space $H^k(\Omega)$ . The norms in $L^{p}(\Omega)$ , $W^{k, p}(\Omega)$ and $H^k(\Omega)$ are defined as the classical senses (cf.^[1]) and are denoted by $\|\cdot\|_{L^{p}}$ , $\|\cdot\|_{W^{k, p}}$ and $\|\cdot\|_{H^{k}}$ . We define $H_0^{k}(\Omega)$ to be the subspace of $H^{k}(\Omega)$ of functions with zero trace on $\partial\Omega$ . The dual space of $H_0^1(\Omega)$ is denoted by $H^{-1}(\Omega)$ . The boldface notations ${{\textbf{L}}}^p(\Omega)$ , ${{\textbf{W}}}^{k, p}(\Omega)$ and ${{\textbf{H}}}^k(\Omega)$ are used to denote the vector-value spaces $L^p(\Omega)^3$ , $W^{k, p}(\Omega)^3$ and $H^k(\Omega)^3$ , respectively. The corresponding norms are denoted by $\|\cdot\|_{{{\textbf{L}}}^{p}}$ , $\|\cdot\|_{{{\textbf{W}}}^{k, p}}$ and $\|\cdot\|_{{{\textbf{H}}}^{k}}$ .

For simplicity, we denote ${{\textbf{V}}} = {{\textbf{H}}}_{0}^{1}(\Omega)$ , $Y = H_{0}^{1}(\Omega)$ ,

$\begin{align} {{\textbf{V}}}_0 = \{{{\textbf{u}}}\in{{\textbf{V}}},\quad \nabla\cdot {{\textbf{u}}} = 0 \ {\mbox{in}} \ \Omega\},\quad Q = L_{0}^{2}(\Omega) = \{q\in\ L^{2}(\Omega),\quad \int_{\Omega}q{\mbox d}{\textbf{x}} = 0\}. \end{align}$

Introduce the trilinear forms $c_1(\cdot, \cdot, \cdot)$ and $c_2(\cdot, \cdot, \cdot)$ , which are given by

$\begin{align*} &c_1({{\textbf{u}}}_{1},{{\textbf{u}}}_{2},{{\textbf{u}}}_{3}) = \int_{\Omega} \left( ({{\textbf{u}}}_{1}\cdot \nabla){{\textbf{u}}}_{2}\cdot{{\textbf{u}}}_{3}+\frac{1}{2} (\nabla\cdot{{\textbf{u}}}_{1}){{\textbf{u}}}_{2}\cdot{{\textbf{u}}}_{3} \right) d{\textbf{x}}\\ &c_2({{\textbf{u}}}_{1},\theta_{1},\theta_{2}) = \int_{\Omega} \left( ({{\textbf{u}}}_{1}\cdot\nabla\theta_{1} )\theta_{2}+\frac{1}{2}(\nabla\cdot{{\textbf{u}}}_{1})\theta_{1}\theta_{2}\right) d{\textbf{x}} \end{align*}$

for ${{\textbf{u}}}_{i}\in{{\textbf{V}}}$ with $i = 1, 2, 3$ and $\theta_{j}\in Y$ with $j = 1, 2$ . By the integration by parts, it is easy to check that the trilinear forms $c_1(\cdot, \cdot, \cdot)$ and $c_2(\cdot, \cdot, \cdot)$ satisfy the skew-symmetric properties, i.e.,

$\begin{align*} c_1({{\textbf{u}}}_{1},{{\textbf{u}}}_{2},{{\textbf{u}}}_{3}) = -c_1({{\textbf{u}}}_{1},{{\textbf{u}}}_{3},{{\textbf{u}}}_{2})\quad{\mbox{and}} \quad c_2({{\textbf{u}}}_{1},\theta_{1},\theta_{2}) = - c_2({{\textbf{u}}}_{1},\theta_{2},\theta_{1}), \end{align*}$

which give

$\begin{align} c_1({{\textbf{u}}}_{1},{{\textbf{u}}}_{2},{{\textbf{u}}}_{2}) = 0 \quad{\mbox{and}} \quad c_2({{\textbf{u}}}_{1},\theta_{1},\theta_{1}) = 0. \end{align}$

(2.1)

For given ${{\textbf{f}}}\in {{\textbf{L}}}^2(\Omega), g\in L^2(\Omega)$ , ${{\textbf{u}}}_0\in{{\textbf{L}}}^2(\Omega)$ and $\theta_0\in L^2(\Omega)$ , in terms of ther skew-symmetric properties (2.1), there {holds} the following energy inequalities:

$\begin{align} \|\theta(t)\|_{L^2}^2+\kappa{ \int}_0^t \|\nabla\theta(s)\|_{{{\textbf{L}}}^2}^2 ds \leq \|\theta_0\|_{L^2}^2 + C { \int}_0^t\|g(s)\|_{L^2}^2 ds, \end{align}$

(2.2)

and

$\begin{align} \begin{split} & \|{{\textbf{u}}}(t)\|_{{{\textbf{L}}}^2}^2+\nu{ \int}_0^t \|\nabla{{\textbf{u}}}(s)\|_{{{\textbf{L}}}^2}^2 ds \\ \leq & \|{{\textbf{u}}}_0\|_{{{\textbf{L}}}^2}^2 + C { \int}_0^t ( \|{{\textbf{f}}}(s)\|_{{{\textbf{L}}}^2}^2 + \|\theta(s)\|_{L^2}^2 + \|\theta(s)\|_{L^2}^2 \|\theta(s)\|_{L^3}^2 ) ds \\ \leq & \|{{\textbf{u}}}_0\|_{{{\textbf{L}}}^2}^2 + C { \int}_0^t ( \|{{\textbf{f}}}(s)\|_{{{\textbf{L}}}^2}^2 + \|g(s)\|_{L^2}^2 ) ds + C \left( { \int}_0^t\|g(s)\|_{L^2}^2 ds \right)^2 \end{split} \end{align}$

(2.3)

for $0\leq t\leq T$ .

Let $0 = t_0 < t_1 < \cdots < t_N = T$ be a uniform partition of the time interval $[0, T]$ with time step $\tau = T/N$ and $t_n = n\tau$ with $0\leq n\leq N$ . Denote ${{\textbf{u}}}^n = {{\textbf{u}}}(t_n)$ , $p^n = p(t_n)$ , $\theta^n = \theta(t_n)$ , ${{\textbf{f}}}^{n} = {{\textbf{f}}}(t_{n})$ and $g^{n} = g(t_{n})$ . For any sequence of functions $\{f^n\}_{n = 0}^N$ , denote the backward Euler discretization and BDF2 discretization by

$D_{1,\tau}f^{n+1} = (f^{n+1}-f^{n})/\tau \quad {\mbox{for}} \ \ 0\leq n \leq N-1,$

and

$D_{2,\tau}f^{n+1} = ({3f^{n+1}-4f^{n}+f^{n-1}})/2\tau, \quad \widehat {f}^{n} = {2f^{n}-f^{n-1}}\quad {\mbox{for}} \ \ 1\leq n \leq N-1,$

where $\widehat {f}^{n}$ is the standard extrapolation formula. For the discrete time derivative $D_{2, \tau}$ , there holds the telescope formula ^[35]:

$\begin{align} \begin{split} (D_{2,\tau}f^{n+1}, f^{n+1}) = & \frac{1}{4\tau} (\|f^{n+1}\|_{L^2}^{2}-\|f^{n}\|_{L^2}^{2}+\|\widehat {f}^{n+1}\|_{L^2}^{2}-\|\widehat {f}^{n}\|_{L^2}^{2}) \\ & + \frac{1}{4\tau} \|f^{n+1}-2f^{n}+f^{n-1}\|_{L^2}^{2} . \end{split} \end{align}$

(2.4)

We define the finite element spaces as follows. Let $\mathcal {T}_h = \{K_j\}_{j = 1}^L$ be a quasi-uniform tetrahedron partition of $\Omega$ with the mesh size $h = \max \, \{diam \, K_j, \ j = 1, \cdots, L\}$ . In finite element discretization of (1.1)–(1.3), we use the Taylor-Hood ( $P_2-P_1$ ) finite element for the velocity ${{\textbf{u}}}$ and pressure $p$ , and the piece-wise quadratic finite element for the temperature $\theta$ . The corresponding finite element subspaces of ${{\textbf{V}}}, Q$ and $Y$ are denoted by ${{\textbf{V}}}_h,$ $Q_h$ and $Y_h$ , respectively, i.e.,

$\begin{align*} {{\textbf{V}}}_h & = \{{{\textbf{v}}}_h\in C(\overline\Omega)\cap {{\textbf{V}}}, \ {{\textbf{v}}}_h|_K \in (P_2(K))^3 , \ \forall \ K\in \mathcal T_h \} , \\ Q_h & = \{q_h\in C(\overline\Omega)\cap H^1(\Omega), \ q_h|_K\in P_1(K), \ \forall \ K\in\mathcal T_h , \ { \int}_\Omega q_h dx = 0\} , \\ Y_h & = \{r_h\in C(\overline\Omega)\cap Y, \ r_h|_K\in P_2(K), \ \forall \ K\in\mathcal T_h \} . \end{align*}$

It is {well known} that the discrete inf-sup condition holds for the inf-sup stable element such that

$\begin{align} \beta_0 \|q_h\|_{L^2} \leq \sup\limits_{{{\textbf{v}}}_h\in {{\textbf{V}}}_h} \frac {(\nabla\cdot{{\textbf{v}}}_h, q_h)}{\|\nabla{{\textbf{v}}}_h\|_{{{\textbf{L}}}^2}}, \quad \forall \ q_h\in Q_h, \end{align}$

(2.5)

where $\beta_{0} > 0$ is independent of the mesh size $h$ .

Denote the space of discrete divergence-free functions by

${{\textbf{V}}}_{0h} = \{{{\textbf{v}}}_{h}\in{{\textbf{V}}}_{h},\quad (\nabla\cdot{{\textbf{v}}}_{h},q_{h}) = 0 \quad\forall \ q_{h}\in Q_{h}\}.$

In the sequel, let $\pi_h$ and $R_h$ be the $L^2$ orthogonal and Ritz orthogonal projections onto $Q_h$ and $Y_h$ , respectively. Then the following error bounds and stability hold ^[4,43]:

$\begin{align} \|p-\pi_{h}p\|_{L^2} & \leq Ch^2 \|p\|_{H^2} \quad \forall \ p\in H^2(\Omega), \end{align}$

(2.6)

$\begin{align} \|\psi-R_{h}\psi\|_{L^2}+h \|\nabla(\psi-R_{h}\psi)\|_{{{\textbf{L}}}^2} & \leq Ch^3\|\psi\|_{H^3} \quad\forall \ \psi\in H^{3}(\Omega), \end{align}$

(2.7)

$\begin{align} \|R_h\psi\|_{W^{1,\infty}} & \leq C\|\psi\|_{H^3} \quad\forall \ \psi\in H^{3}(\Omega). \end{align}$

(2.8)

Let $I_{h}$ be the Lagrange interpolation operator onto ${{\textbf{V}}}_h$ . The following bound can be found in ^[4]:

$\begin{align} \|{{\textbf{u}}}-I_{h}{{\textbf{u}}}\|_{{{\textbf{W}}}^{m,p}}\leq Ch^{n-m}\|{{\textbf{u}}}\|_{{{\textbf{W}}}^{n,p}} \quad 0\leq m\leq n\leq 3 , \end{align}$

(2.9)

where $n > 3/p$ when $1 < p\leq\infty$ and $n\geq 3$ when $p = 1$ .

Next, we recall some known results about finite element approximations of Stokes problem in ^[17,18]. Consider the Stokes problem:

$\begin{align} \begin{split} -\nu\Delta {{\textbf{u}}} + \nabla p = {{\textbf{F}}} & \quad {\mbox{in}} \ \Omega, \\ \nabla\cdot{{\textbf{u}}} = 0 &\quad {\mbox{in}} \ \Omega, \\ {{\textbf{u}}} = 0 & \quad {\mbox{on}} \ \partial \Omega. \end{split} \end{align}$

(2.10)

Let us denote $({{\textbf{u}}}_h, p_h)\in{{\textbf{V}}}_h\times Q_h$ , the mixed finite element approximation solution to the Stokes problem (2.10). Then one has error estimates ^[22,30]

$\begin{align} \|\nabla({{\textbf{u}}}-{{\textbf{u}}}_h)\|_{{{\textbf{L}}}^2} & \leq C \left( \inf\limits_{{{\textbf{v}}}_h\in{{\textbf{V}}}_h} \|{{\textbf{u}}}-{{\textbf{v}}}_h\|_{{{\textbf{H}}}^1} + \nu^{-1} \inf\limits_{q_h\in Q_h}\|p-q_h\|_{L^2} \right), \end{align}$

(2.11)

$\begin{align} \|p-p_h\|_{L^2} & \leq C \left( \nu \inf\limits_{{{\textbf{v}}}_h\in{{\textbf{V}}}_h} \|{{\textbf{u}}}-{{\textbf{v}}}_h\|_{{{\textbf{H}}}^1} + \inf\limits_{q_h\in Q_h}\|p-q_h\|_{L^2} \right), \end{align}$

(2.12)

$\begin{align} \|{{\textbf{u}}}-{{\textbf{u}}}_h\|_{{{\textbf{L}}}^2} & \leq C h \left( \inf\limits_{{{\textbf{v}}}_h\in{{\textbf{V}}}_h} \|{{\textbf{u}}}-{{\textbf{v}}}_h\|_{{{\textbf{H}}}^1} + \nu^{-1} \inf\limits_{q_h\in Q_h}\|p-q_h\|_{L^2} \right), \end{align}$

(2.13)

where $C > 0$ is independent of $h, \nu$ and $\kappa$ . It is clear that the error estimates of the velocity depend on the negative power of $\nu$ . To prove error bounds with constants being independent of the negative power of $\nu$ , we recall the specific projection of $({{\textbf{u}}}, p)$ introduced in ^[17] that we will denote by ${{\textbf{s}}}_h$ defined by

$(\nabla{{\textbf{s}}}_h,\nabla{{\textbf{v}}}_h) = (\nabla{{\textbf{u}}},\nabla{{\textbf{v}}}_h), \quad {{\textbf{v}}}_h\in {{\textbf{V}}}_{0h}.$

Let $({{\textbf{u}}}, p, \theta)$ be the solution to the system (1.1)–(1.5) with ${{\textbf{u}}}\in {{\textbf{L}}}^\infty((0, T];{{\textbf{V}}}\cap{{\textbf{H}}}^3(\Omega))$ , $p\in L^\infty((0, T];Q\cap H^2(\Omega))$ , $\theta\in L^\infty((0, T];Y\cap H^3(\Omega))$ and ${{\textbf{u}}}_t\in {{\textbf{L}}}^\infty((0, T];{{\textbf{H}}}^1(\Omega))$ . We note that $({{\textbf{u}}}, 0)$ is the solution to the Stokes problem (2.10) with the right-hand side ${{\textbf{F}}}$ given by

$\begin{align} {{\textbf{F}}} = {{\textbf{f}}}-{{\textbf{u}}}_t-({{\textbf{u}}}\cdot\nabla){{\textbf{u}}} -\nabla p +(\gamma_{1}\theta+\gamma_{2}\theta^{2}){\textbf{i}}_{3}. \end{align}$

(2.14)

Let $({{\textbf{s}}}_h, l_h)\in {{\textbf{V}}}_h\times Q_h$ be the corresponding mixed finite element approximation of $({{\textbf{u}}}, 0)$ . Then, from (2.11)–(2.13), we have

$\begin{align} \|{{\textbf{u}}}-{{\textbf{s}}}_{h}\|_{{{\textbf{L}}}^2}+h\|{{\textbf{u}}}-{{\textbf{s}}}_{h}\|_{{{\textbf{H}}}^{1}}&\leq Ch^3\|{{\textbf{u}}}\|_{{{\textbf{H}}}^{3}}, \end{align}$

(2.15)

$\begin{align} \|l_{h}\|_{L^2}&\leq C\nu h^2 \|{{\textbf{u}}}\|_{{{\textbf{H}}}^3}, \end{align}$

(2.16)

where $C > 0$ is independent of $h$ and $\nu$ . In terms of the inverse inequality in the theory of the finite element method and $\|I_{h}{{\textbf{u}}}\|_{{{\textbf{L}}}^\infty}\leq C\|{{\textbf{u}}}\|_{{{\textbf{L}}}^\infty}$ , we have

$\begin{align} \begin{split} \|{{\textbf{s}}}_{h}\|_{{{\textbf{L}}}^\infty}&\leq\|{{\textbf{s}}}_{h}-I_{h}{{\textbf{u}}}\|_{{{\textbf{L}}}^\infty}+\|I_{h}{{\textbf{u}}}\|_{{{\textbf{L}}}^\infty} \\ &\leq C(h^{-3/2}\|{{\textbf{s}}}_{h}-I_{h}{{\textbf{u}}}\|_{{{\textbf{L}}}^2}+\|{{\textbf{u}}}\|_{{{\textbf{L}}}^\infty})\\ &\leq C ( h^{-3/2}\|{{\textbf{s}}}_{h}-{{\textbf{u}}}\|_{{{\textbf{L}}}^2}+h^{-3/2}\|{{\textbf{u}}}-I_{h}{{\textbf{u}}}\|_{{{\textbf{L}}}^2}+\|{{\textbf{u}}}\|_{{{\textbf{H}}}^3})\\ &\leq C\|{{\textbf{u}}}\|_{{{\textbf{H}}}^3}, \end{split} \end{align}$

(2.17)

where $C > 0$ is independent of $h, \nu$ and $\kappa$ . Furthermore, there holds

$\begin{align} \|\nabla({{\textbf{u}}}-{{\textbf{s}}}_{h})\|_{{{\textbf{L}}}^\infty}\leq C\|\nabla{{\textbf{u}}}\|_{{{\textbf{L}}}^\infty}\leq C\|{{\textbf{u}}}\|_{{{\textbf{H}}}^3} \quad {\mbox{and}}\quad \|\nabla {{\textbf{s}}}_{h}\|_{{{\textbf{L}}}^\infty}\leq C\|\nabla{{\textbf{u}}}\|_{{{\textbf{L}}}^\infty}\leq C\|{{\textbf{u}}}\|_{{{\textbf{H}}}^3}, \end{align}$

(2.18)

where $C > 0$ is independent of $h, \nu$ and $\kappa$ .

Finally, we recall the discrete Gronwall inequality established in ^[24].

Lemma 2.1. Let $a_k, b_k, c_k$ and $\gamma_k$ , for integers $k\geq 0$ , be the nonnegative numbers such that

$\begin{align} a_n+\tau \sum\limits_{k = 0}^n b_k\leq \tau \sum\limits_{k = 0}^n \gamma_ka_k+ \tau \sum\limits_{k = 0}^n c_k+B \qquad {\mbox{for}} \ n\geq 0. \end{align}$

(2.19)

Suppose that $\tau\gamma_k < 1$ for all $k$ and set $\sigma_k = (1-\tau \gamma_k)^{-1}$ , then

$\begin{align} a_n+\tau \sum\limits_{k = 0}^n b_k\leq \exp(\tau \sum\limits_{k = 0}^n \gamma_k \sigma_k)(\tau \sum\limits_{k = 0}^n c_k+B ) \qquad {\mbox{for}} \ n\geq 0. \end{align}$

(2.20)

Remark 2.1. If the first sum on the right in (2.19) extends only up to $n-1$ , then the estimate (2.20) holds for all $\tau > 0$ with $\sigma_k = 1$ .

3. The Euler grad-div stabilization finite element approximation

In this section, we consider the first-order Euler fully discrete scheme for the problems (1.1)–(1.5). This Euler scheme is a semi-implicit scheme, i.e., one only solves two linear systems at each time discrete level. Based on the grad-div stabilization method, we propose the following first-order Euler finite element scheme for $0\leq n\leq N-1$ .

Step Ⅰ: Find $\theta_h^{n+1}\in Y_h$ by

$\begin{align} (D_{1,\tau}\theta_{h}^{n+1},\psi_{h})+\kappa(\nabla \theta_{h}^{n+1},\nabla\psi_{h})+c_2({{\textbf{u}}}_{h}^{n},\theta_{h}^{n+1},\psi_{h}) = (g^{n+1},\psi_{h}) \end{align}$

(3.1)

with the initial iteration values $\theta_{h}^{0} = R_h\theta_{0}$ , ${{\textbf{u}}}_{h}^{0} = I_h{{\textbf{u}}}_{0}$ and for any $\psi_h\in Y_h$ .

Step Ⅱ: Find $({{\textbf{u}}}_{h}^{n+1}, p_{h}^{n+1})\in {{\textbf{V}}}_{h}\times Q_{h}$ by

$\begin{eqnarray} (D_{1,\tau}{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})+ \nu(\nabla {{\textbf{u}}}_{h}^{n+1},\nabla {{\textbf{v}}}_{h})+c_1({{\textbf{u}}}_{h}^{n},{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})-(\nabla\cdot{{\textbf{v}}}_{h}, p_{h}^{n+1}) +(\nabla\cdot{{\textbf{u}}}_{h}^{n+1}, q_{h}) \\ +\beta(\nabla\cdot{{\textbf{u}}}_{h}^{n+1},\nabla\cdot{{\textbf{v}}}_{h}) -\gamma_{1}(\theta_{h}^{n}\textbf i_{3},{{\textbf{v}}}_{h})-\gamma_{2}((\theta_{h}^{n})^{2}\textbf i_{3},{{\textbf{v}}}_{h}) = ({{\textbf{f}}}^{n+1},{{\textbf{v}}}_{h}) \end{eqnarray}$

(3.2)

for any $({{\textbf{v}}}_{h}, q_{h})\in {{\textbf{V}}}_{h}\times Q_{h}$ , where $\beta > 0$ is the stabilization parameter. The detailed discussion about the choice of $\beta$ can be found in ^[28].

Next lemma gives the stability of the scheme (3.1)–(3.2), where the discrete energy inequalities can be viewed as the discrete version of the energy estimates (2.2)–(2.3). Moreover, since (3.1)–(3.2) are both linearized problems, the discrete energy inequalities imply the existence and uniqueness of numerical solution $({{\textbf{u}}}_h^{n+1}, p_h^{n+1}, \theta_h^{n+1})$ to the Euler scheme (3.1)–(3.2).

Lemma 3.1. For $0\leq n\leq N-1$ and all $\tau > 0$ , $h > 0$ , the finite element scheme (3.1) and (3.2) has a unique solution $\theta_{h}^{n+1}\in Y_{h}$ and $({{\textbf{u}}}_{h}^{n+1}, p_{h}^{n+1})\in {{\textbf{V}}}_{h}\times Q_{h}$ . Moreover, the discrete energy inequalities hold:

$\begin{align} \|\theta_{h}^{m+1}\|_{L^{2}}^{2}+\kappa\tau\sum\limits_{n = 0}^{m}\|\nabla{\theta}_{h}^{n+1}\|_{{{\textbf{L}}}^{2}}^{2} \leq C\tau\sum\limits_{n = 0}^{N-1}\|g^{n+1}\|_{L^{2}}^{2}+\|\theta_{h}^{0}\|_{L^{2}}^{2}, \end{align}$

(3.3)

and

$\begin{align} \begin{split} & \|{{\textbf{u}}}_{h}^{m+1}\|_{{{\textbf{L}}}^{2}}^{2}+\tau\nu \sum\limits_{n = 0}^{m}\|\nabla {{{\textbf{u}}}}_{h}^{m+1}\|_{{{\textbf{L}}}^{2}}^{2} +2\tau\beta \sum\limits_{n = 0}^{m}\|\nabla\cdot {{{\textbf{u}}}}_{h}^{m+1}\|_{L^{2}}^{2} \\ \leq & C\tau\sum\limits_{n = 0}^{N-1}(\|{{\textbf{f}}}^{n+1}\|_{{{\textbf{L}}}^{2}}^{2}+\|g^{n+1}\|_{L^{2}}^{2} ) + C (\|{{\textbf{u}}}_{h}^{0}\|_{{{\textbf{L}}}^{2}}^{2} + \|\theta_{h}^{0}\|_{L^{2}}^{2} ) +C \left(\tau\sum\limits_{n = 0}^{N-1}\|g^{n+1}\|_{L^{2}}^{2}+\|\theta_{h}^{0}\|_{L^{2}}^{2} \right)^2 \end{split} \end{align}$

(3.4)

for all $0\leq m\leq N-1$ , where $C > 0$ is independent of $\tau$ and $h$ .

Proof. Setting $\psi_{h} = 2\tau\theta_{h}^{n+1}$ in (3.1), using the skew-symmetric property (2.1), the Hölder inequality and Young inequality, we have

$\begin{align*} \|\theta_{h}^{n+1}\|_{L^2}^{2}-\|\theta_{h}^{n}\|_{L^2}^{2}+ \|\theta_{h}^{n+1}-\theta_{h}^{n}\|_{L^2}^{2}+2\tau\kappa\|\nabla\theta_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} &\leq2\tau\|g^{n+1}\|_{L^2}\|\theta_{h}^{n+1}\|_{L^2}\\ &\leq C\tau\|g^{n+1}\|_{L^2}^{2}+\tau\kappa\|\nabla\theta_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}. \end{align*}$

Summing up the above inequality from $n = 0$ to $n = m$ gives (3.3).

Similarly, taking $({{\textbf{v}}}_{h}, q_h) = 2\tau({{\textbf{u}}}_{h}^{n+1}, p_h^{n+1})$ in (3.2), we have

$\begin{align*} &\|{{\textbf{u}}}_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}-\|{{\textbf{u}}}_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}+\|{{\textbf{u}}}_{h}^{n+1}-{{\textbf{u}}}_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}+2\tau\nu\|\nabla {{\textbf{u}}}_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}+2\tau\beta\|\nabla\cdot{{\textbf{u}}}_{h}^{n+1}\|_{L^2}^{2}\\ &\leq C\tau\|{{\textbf{f}}}^{n+1}\|_{{{\textbf{L}}}^2}\|{{\textbf{u}}}_{h}^{n+1}\|_{{{\textbf{L}}}^2}+C\tau\|\theta_{h}^{n}\|_{L^2}\|{{\textbf{u}}}_{h}^{n+1}\|_{{{\textbf{L}}}^2} +C\tau \|\theta_{h}^{n}\|_{L^2} \|\theta_{h}^{n}\|_{L^3}\|{{\textbf{u}}}_{h}^{n+1}\|_{{{\textbf{L}}}^6}\\ &\leq \tau\nu\|\nabla{{\textbf{u}}}_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + C\tau(\|{{\textbf{f}}}^{n+1}\|_{L^2}^{2}+\|\theta_{h}^{n}\|_{L^2}^{2}+\|\theta_{h}^{n}\|_{L^2}^{2}\|\theta_{h}^{n}\|_{L^3}^{2}). \end{align*}$

Summing up the above inequality from $n = 0$ to $n = m$ and noticing (3.3), we get (3.4).□

Denote

$e_{h}^{n} = {{\textbf{s}}}_{h}^{n}-{{\textbf{u}}}_{h}^{n}, \quad \theta^{n}-\theta_{h}^{n} = (\theta^{n}-R_h\theta^{n})+(R_h\theta^{n}-\theta_{h}^{n}) : = \eta_{\theta}^{n}+e_{\theta}^{n}{.}$

The ${{\textbf{s}}}_h(t)$ is the finite element solution to the Stokes problem (2.10) with the right-hand side (2.14) and ${{\textbf{s}}}_{h}^{n} = {{\textbf{s}}}_{h}(t_n)$ .

Thanks to the discrete inf-sup condition (2.5), the finite element scheme (3.2) is equivalent to: Find ${{\textbf{u}}}_{h}^{n+1}\in {{\textbf{V}}}_{0h}$ by

$\begin{eqnarray} (D_{1,\tau}{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})+ \nu(\nabla {{\textbf{u}}}_{h}^{n+1},\nabla {{\textbf{v}}}_{h})+c_1({{\textbf{u}}}_{h}^{n},{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h}) +\beta(\nabla\cdot{{\textbf{u}}}_{h}^{n+1},\nabla\cdot{{\textbf{v}}}_{h}) \\ -\gamma_{1}(\theta_{h}^{n}\textbf i_{3},{{\textbf{v}}}_{h})-\gamma_{2}((\theta_{h}^{n})^{2}\textbf i_{3},{{\textbf{v}}}_{h}) = ({{\textbf{f}}}^{n+1},{{\textbf{v}}}_{h}) \quad\quad \forall \ {{\textbf{v}}}_h\in {{\textbf{V}}}_{0h}. \end{eqnarray}$

(3.5)

For simplicity, we take ${{\textbf{s}}}_h^{0} = {{\textbf{u}}}_{h}^{0}$ . For $0\leq n\leq N-1$ , ${{\textbf{s}}}_h^{n+1}$ and $\theta^{n+1}$ satisfy

$\begin{align} \begin{split} &(D_{1,\tau}{{\textbf{s}}}_h^{n+1},{{\textbf{v}}}_{h})+\nu(\nabla {{\textbf{s}}}_h^{n+1},\nabla {{\textbf{v}}}_{h}) +c_1({{\textbf{s}}}_h^{n},{{\textbf{s}}}_h^{n+1},{{\textbf{v}}}_{h})+ \beta(\nabla\cdot {{\textbf{s}}}_h^{n+1},\nabla\cdot{{\textbf{v}}}_{h}) \\ & -\gamma_{1}(\theta^{n}\textbf i_{3},{{\textbf{v}}}_{h}) -\gamma_{2}((\theta^{n})^{2}\textbf i_{3},{{\textbf{v}}}_{h})\\ & = ({{\textbf{f}}}^{n+1},{{\textbf{v}}}_{h})+(\nabla\cdot{{\textbf{v}}}_{h}, p^{n+1}-\pi_{h}p^{n+1}) \\ & +\beta(\nabla\cdot({{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1}),\nabla\cdot{{\textbf{v}}}_{h}) +(D_{1,\tau}{{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}_{t}(t_{n+1}),{{\textbf{v}}}_{h})\\ & +c_1({{\textbf{s}}}_h^{n},{{\textbf{s}}}_h^{n+1},{{\textbf{v}}}_{h})-c_1({{\textbf{u}}}^{n+1},{{\textbf{u}}}^{n+1},{{\textbf{v}}}_{h})\\ &+\gamma_{1}(\theta^{n+1}-\theta^{n},\textbf i_{3}\cdot {{\textbf{v}}}_{h}) + \gamma_{2} ((\theta^{n+1})^{2}-(\theta^{n})^{2}, \textbf i_{3}\cdot {{\textbf{v}}}_{h}) \quad \forall \ {{\textbf{v}}}_h\in{{\textbf{V}}}_{0h}, \end{split} \end{align}$

(3.6)

and

$\begin{eqnarray} (D_{1,\tau}\theta^{n+1},\psi_{h})+\kappa(\nabla \theta^{n+1},\nabla\psi_{h})+c_2({{\textbf{u}}}^{n+1},\theta^{n+1},\psi_{h}) = (D_{1,\tau}\theta^{n+1}-\theta_t (t_{n+1}),\psi_{h}) \\ +(g^{n+1},\psi_{h}) \quad \forall \ \psi_h\in Y_h. \end{eqnarray}$

(3.7)

In (3.6), we use $\nabla\cdot{{\textbf{u}}}^{n+1} = 0$ in $\Omega$ and $(\nabla\cdot{{\textbf{v}}}_h, \pi_h p^{n+1}) = 0$ due to ${{\textbf{v}}}_h\in {{\textbf{V}}}_{0h}$ .

Next, we give error equations. Subtracting (3.5) from (3.6) leads to

$\begin{align} (D_{1,\tau}e_{h}^{n+1},{{\textbf{v}}}_{h})+\nu(\nabla e_{h}^{n+1},\nabla{{\textbf{v}}}_{h}) +\beta(\nabla\cdot e_{h}^{n+1},\nabla\cdot{{\textbf{v}}}_{h}) = \sum_{j = 1}^6 (J_{j}^{n+1},{{\textbf{v}}}_{h}) \quad \forall \ {{\textbf{v}}}_h\in{{\textbf{V}}}_{0h}, \end{align}$

(3.8)

where

$\begin{eqnarray*} (J_{1}^{n+1},{{\textbf{v}}}_{h})& = &( D_{1,\tau}{{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}_{t}(t_{n+1}),{{\textbf{v}}}_{h}) , \\ (J_{2}^{n+1},{{\textbf{v}}}_{h})& = &(\nabla\cdot{{\textbf{v}}}_{h},(p^{n+1}-\pi_h p^{n+1})+\beta\nabla\cdot({{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1})), \\ (J_{3}^{n+1},{{\textbf{v}}}_{h})& = &c_1({{\textbf{s}}}_h^{n},{{\textbf{s}}}_h^{n+1},{{\textbf{v}}}_{h})-c_1({{\textbf{u}}}^{n+1},{{\textbf{u}}}^{n+1},{{\textbf{v}}}_{h}), \\ (J_{4}^{n+1},{{\textbf{v}}}_{h})& = &c_1({{\textbf{u}}}_{h}^{n},{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h}) -c_1({{\textbf{s}}}_h^{n},{{\textbf{s}}}_h^{n+1},{{\textbf{v}}}_{h}),\\ (J_{5}^{n+1},{{\textbf{v}}}_{h})& = &\gamma_{1}(\theta^{n+1}-\theta^{n},\textbf i_{3}\cdot {{\textbf{v}}}_{h}) + \gamma_{2} ((\theta^{n+1})^{2}-(\theta^{n})^{2}, \textbf i_{3}\cdot {{\textbf{v}}}_{h}),\\ (J_{6}^{n+1},{{\textbf{v}}}_{h})& = &\gamma_{1}((\theta^{n}-\theta_{h}^{n}),\textbf i_{3}\cdot{{\textbf{v}}}_{h})+ \gamma_{2}(((\theta^{n})^{2}-(\theta_{h}^{n})^{2}),\textbf i_{3}\cdot{{\textbf{v}}}_{h}). \end{eqnarray*}$

Subtracting (3.1) from (3.7), we have

$\begin{align} (D_{1,\tau}e_{\theta}^{n+1},\psi_{h})+\kappa(\nabla e_{\theta}^{n+1},\nabla\psi_{h})+c_2({{\textbf{u}}}_h^n,e_{\theta}^{n+1},\psi_h) = \sum_{j = 1}^{3}(\lambda_{j}^{n+1},\psi_{h}) \quad \forall \ \psi_h\in Y_{h}, \end{align}$

(3.9)

where

$\begin{eqnarray*} (\lambda_{1}^{n+1},\psi_{h})& = &( D_{1,\tau}\theta^{n+1}-\theta_t(t_{n+1}),\psi_{h}),\\ (\lambda_{2}^{n+1},\psi_{h})& = &-(D_{1,\tau}\eta_{\theta}^{n+1},\psi_h),\\ (\lambda_{3}^{n+1},\psi_{h})& = &c_2({{\textbf{u}}}_{h}^{n},R_h\theta^{n+1},\psi_{h}) -c_2({{\textbf{u}}}^{n+1},\theta^{n+1},\psi_{h}). \end{eqnarray*}$

The main result in this section is presented in the following theorem.

Theorem 3.1. Let $({{\textbf{u}}}, p, \theta)$ be the unique solution to the time-dependent penetrative convection problems (1.1)–(1.5) and satisfy the following regularity assumptions:

$\begin{align} \begin{cases} {{\textbf{u}}}_0\in {{\textbf{H}}}^3(\Omega), \ \theta_0\in H^3(\Omega), \ {{\textbf{f}}}\in L^{\infty}((0,T];{{\textbf{L}}}^{2}(\Omega)),\quad g\in L^{\infty}((0,T];L^{2}(\Omega)) , \\ {{\textbf{u}}}\in {{\textbf{L}}}^{\infty}((0,T];{{\textbf{H}}}^{3}(\Omega)),\quad p\in L^{\infty}((0,T];H^{2}(\Omega)), \quad \theta\in L^{\infty}((0,T];H^{3}(\Omega)) , \\ {{\textbf{u}}}_t\in {{\textbf{L}}}^{\infty}((0,T];{{\textbf{H}}}^{1}(\Omega))\cap {{\textbf{L}}}^{2}((0,T];{{\textbf{H}}}^{2}(\Omega)),\\ \theta_t\in L^{\infty}((0,T];H^{1}(\Omega))\cap L^{2}((0,T];H^{2}(\Omega)), \\ {{\textbf{u}}}_{tt}\in {{\textbf{L}}}^{2}((0,T];{{\textbf{L}}}^{2}(\Omega)),\quad \theta_{tt} \in L^{2}((0,T];L^{2}(\Omega)). \end{cases} \end{align}$

(3.10)

Then, for the first-order Euler grad-div scheme (3.1)–(3.2) under the time step condition $\tau\leq Ch^2$ , when $h$ and $\tau$ are sufficiently small, we have the following error estimate:

$\begin{align} \|e_h^{m+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_\theta^{m+1}\|_{L^2}^{2} +\tau\sum\limits_{n = 0}^{m} \left ( \nu\|\nabla e_h^{m+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{m+1}\|_{{{\textbf{L}}}^2}^{2} \right ) \leq C_0^2(h^{4}+\tau^2) \end{align}$

(3.11)

with $0\leq m\leq N-1$ , where $C_0 > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Proof. Taking ${{\textbf{v}}}_h = 2\tau e_h^{n+1}$ in (3.8) and summing up the resulting equation from $n = 0$ to $n = m$ , we have

$\begin{align} \|e_{h}^{m+1}\|_{{{\textbf{L}}}^2}^{2}+2\tau\sum\limits_{n = 0}^{m} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{n+1}\|_{L^2}^{2} \right ) = 2\tau\sum\limits_{n = 0}^{m} \sum_{j = 1}^6 (J_{j}^{n+1},e_{h}^{n+1}) . \end{align}$

(3.12)

We estimate the right-hand side of (3.12) as follows. For $J_1^{n+1}$ , we split it as

$\begin{align*} J_1^{n+1} = ( D_{1,\tau}{{\textbf{s}}}_h^{n+1}-D_{1,\tau}{{\textbf{u}}}^{n+1})+ (D_{1,\tau}{{\textbf{u}}}^{n+1} -{{\textbf{u}}}_{t}(t_{n+1})) . \end{align*}$

According to the Taylor's formula and the regularity assumption (3.10), one has

$\begin{align*} \tau\sum\limits_{n = 0}^{m} \|D_{1,\tau}{{\textbf{u}}}^{n+1} -{{\textbf{u}}}_{t}(t_{n+1})\|_{{{\textbf{L}}}^2}^2 \leq C\tau^2. \end{align*}$

Notice that

$\begin{align*} D_{1,\tau}{{\textbf{s}}}_h^{n+1}-D_{1,\tau}{{\textbf{u}}}^{n+1} = \frac 1\tau { \int}_{t_n}^{t_{n+1}} \partial_t({{\textbf{s}}}_h(t)-{{\textbf{u}}}(t)) dt . \end{align*}$

From the Hölder inequality, we have

$\begin{align*} \|D_{1,\tau}{{\textbf{s}}}_h^{n+1}-D_{1,\tau}{{\textbf{u}}}^{n+1}\|_{{{\textbf{L}}}^2} & \leq \frac 1{\tau} { \int}_{t_n}^{t_{n+1}} \| \partial_t({{\textbf{s}}}_h(t)-{{\textbf{u}}}(t))\|_{{{\textbf{L}}}^2} dt \\ & \leq \frac {Ch^2}{\tau} { \int}_{t_n}^{t_{n+1}} \| \partial_t{{\textbf{u}}}(t)\|_{{{\textbf{H}}}^2} dt \\ & \leq \frac {Ch^2}{\tau^{1/2}} \left( { \int}_{t_n}^{t_{n+1}} \| \partial_t{{\textbf{u}}}(t)\|_{{{\textbf{H}}}^2}^2 dt \right)^{1/2} . \end{align*}$

Thus,

$\begin{align} 2\tau\sum\limits_{n = 0}^{m} (J_{1}^{n+1},e_{h}^{n+1}) \leq C \tau \sum_{n = 0}^{m} \| e_h^{n+1}\|_{{{\textbf{L}}}^2}^{2} + C(\tau^2+h^{4}), \end{align}$

(3.13)

where $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ .

For $J_2^{n+1}$ , it follows from the Hölder inequality and (2.6) that

$\begin{align*} (J_{2}^{n+1},e_{h}^{n+1}) \leq & \left( \|p^{n+1}-\pi_h p^{n+1}\|_{L^2} + \beta \|{{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1}\|_{{{\textbf{H}}}^1} \right) \|\nabla\cdot e_h^{n+1}\|_{L^2} \\ \leq & C h^2 \|\nabla\cdot e_h^{n+1}\|_{L^2}{,} \end{align*}$

then

$\begin{align} 2\tau \sum\limits_{n = 0}^{m} (J_{2}^{n+1},e_{h}^{n+1}) \leq Ch^{4} + \epsilon_1\beta \tau \sum_{n = 0}^{m} \|\nabla\cdot e_h^{n+1}\|_{L^2}^{2}, \end{align}$

(3.14)

where $C > 0$ is independent of $h, \tau, \nu, \kappa$ and $\epsilon_1 > 0$ is some small constant determined later.

For $J_3^{n+1}$ , in terms of the skew-symmetric property (2.1), (2.18), the Hölder inequality and the regularity assumption (3.10), there holds

$\begin{align*} (J_{3}^{n+1},e_{h}^{n+1}) & = c_1({{\textbf{s}}}_h^{n}-{{\textbf{u}}}^{n},{{\textbf{s}}}_h^{n+1},e_{h}^{n+1})+c_1({{\textbf{u}}}^{n}-{{\textbf{u}}}^{n+1},{{\textbf{s}}}_h^{n+1},e_{h}^{n+1}) \\ &\quad +c_1({{\textbf{u}}}^{n+1},{{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1},e_{h}^{n+1}) \\ &\leq (\|{{\textbf{s}}}_h^{n}-{{\textbf{u}}}^{n}\|_{{{\textbf{H}}}^1}+\|{{\textbf{u}}}^{n+1}-{{\textbf{u}}}^{n}\|_{{{\textbf{H}}}^1})\|{{\textbf{s}}}_h^{n+1}\|_{{{\textbf{W}}}^{1,\infty}}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}\\ &\quad + \|{{\textbf{u}}}^{n+1}\|_{{{\textbf{L}}}^{\infty}}\|{{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1}\|_{{{\textbf{H}}}^1}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}\\ & \leq C (h^2+\tau ) \|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}, \end{align*}$

where $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ . Thus, the third term can be estimated by

$\begin{align} 2\tau \sum\limits_{n = 0}^{m} (J_{3}^{n+1},e_{h}^{n+1}) \leq C (h^{4}+\tau^2) + C \tau \sum_{n = 0}^{m} \| e_h^{n+1}\|_{{{\textbf{L}}}^2}^2, \end{align}$

(3.15)

where $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ .

A similar method to bound $J_4^{n+1}$ leads to

$\begin{align*} (J_{4}^{n+1},e_{h}^{n+1}) & = c_1 (e_{h}^{n},{{\textbf{s}}}_h^{n+1},e_{h}^{n+1}) \\ & \leq \|{e}_{h}^{n}\|_{{{\textbf{L}}}^2}\|\nabla{{\textbf{s}}}_h^{n+1}\|_{{{\textbf{L}}}^\infty}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2} +\|{{\textbf{s}}}_h^{n+1}\|_{{{\textbf{L}}}^\infty}\|\nabla\cdot e_{h}^{n}\|_{L^2}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}, \end{align*}$

and

$\begin{align} 2\tau \sum\limits_{n = 0}^{m} (J_{4}^{n+1},e_{h}^{n+1}) \leq \epsilon_1\beta \tau \sum_{n = 0}^{m} \|\nabla\cdot e_h^{n+1}\|_{L^2}^{2} + C\tau \sum_{n = 0}^{m}(\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}), \end{align}$

(3.16)

where $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ .

For $J_5^{n+1}$ , using the Hölder inequality, Young inequality and the regularity assumption (3.10), we have

$\begin{align*} |(J_{5}^{n+1},e_{h}^{n+1})|&\leq C(1+\|\theta^{n}+\theta^{n+1}\|_{L^\infty}) \int_{t_{n}}^{t_{n+1}}\|\theta_{t}(t)\|_{L^2}dt \|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}\\ &\leq C\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}+C\tau^2{,} \end{align*}$

then

$\begin{align} 2\tau \sum\limits_{n = 0}^{m} (J_{5}^{n+1},e_{h}^{n+1}) \leq C \tau \sum_{n = 0}^{m} \| e_h^{n+1}\|_{{{\textbf{L}}}^2}^2 + C\tau^2, \end{align}$

(3.17)

where $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ .

By a similar method, we estimate the last term by

$\begin{align*} (J_{6}^{n+1},e_{h}^{n+1}) &\leq \|\eta_{\theta}^{n}+e_{\theta}^{n}\|_{L^2}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}+ \|(\eta_{\theta}^{n}+e_{\theta}^{n})\|_{L^{2}}\|\theta^{n}+\theta_{h}^{n}\|_{L^{\infty}}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}\\ &\leq C (1+\|\theta_{h}^{n}\|_{L^{\infty}}^2) (h^{4}+\|e_{\theta}^{n}\|_{L^2}^2)+ C (\|e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}) , \end{align*}$

where we use the Sobolev imbedding inequality and $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ . Taking the sum gives

$\begin{align} \begin{split} 2\tau \sum\limits_{n = 0}^{m} (J_{6}^{n+1},e_{h}^{n+1}) \leq & C h^{4} \tau \sum_{n = 0}^{m} (1+\|\theta_{h}^{n}\|_{L^{\infty}}^2) + C \tau \sum_{n = 0}^{m} (\|e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}) \\ & + C\tau \sum_{n = 0}^{m} (1+\|\theta_{h}^{n}\|_{L^{\infty}}^2) \| e_{\theta}^{n}\|_{L^2}^2 , \end{split} \end{align}$

(3.18)

where $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ .

Substituting the estimates (3.13)–(3.18) into (3.12), we obtain

$\begin{align} \begin{split} & \|e_{h}^{m+1}\|_{{{\textbf{L}}}^2}^{2}+2\tau\sum\limits_{n = 0}^{m} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{n+1}\|_{L^2}^{2} \right ) \\ \leq & C(h^{4}+\tau^2) + C h^{4} \tau \sum_{n = 0}^{m} (1+\|\theta_{h}^{n}\|_{L^{\infty}}^2) + C \tau \sum_{n = 0}^{m} (\|e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}) \\ & +2\epsilon_1\beta \tau \sum_{n = 0}^{m} \|\nabla\cdot e_h^{n+1}\|_{L^2} + C\tau \sum_{n = 0}^{m} (1+\|\theta_{h}^{n}\|_{L^{\infty}}^2) \|e_{\theta}^{n}\|_{L^2}^2 , \end{split} \end{align}$

(3.19)

where $C > 0$ is independent of $h, \tau, \nu$ and $\kappa$ .

Next, we estimate $e_\theta^{n+1}$ . Taking $\psi_h = 2\tau e_\theta^{n+1}$ in (3.9), summing up the resulting equation from $n = 0$ to $n = m$ and noticing the skew-symmetric property (2.1), we have

$\begin{align} \|e_{\theta}^{m+1}\|_{L^2}^{2}+2\kappa \tau\sum\limits_{n = 0}^{m} \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} = 2\tau\sum\limits_{n = 0}^{m} \sum_{j = 1}^3 (\lambda_{j}^{n+1}, e_\theta^{n+1}) , \end{align}$

(3.20)

where we noted $e_{\theta}^0 = 0$ . By the Taylor's formula and the Hölder inequality, it is easy to prove that

$\begin{align} \begin{split} 2\tau\sum\limits_{n = 0}^{m} (\lambda_{1}^{n+1}, e_\theta^{n+1}) \leq C \tau \sum_{n = 0}^{m} \| e_\theta^{n+1}\|_{L^2}^2 + C\tau^2 , \\ \end{split} \end{align}$

(3.21)

$\begin{align} \begin{split} 2\tau\sum\limits_{n = 0}^{m} (\lambda_{2}^{n+1}, e_\theta^{n+1}) &\leq C \tau \sum_{n = 0}^{m} \| e_{\theta}^{n+1}\|_{L^2}^2 + Ch^{4} \tau \sum_{n = 0}^{m} \|D_{1,\tau}\theta^{n+1}\|_{H^1}^2 \\ &\leq C \tau \sum_{n = 0}^{m} \| e_{\theta}^{n+1}\|_{L^2}^2 + Ch^{4} \end{split} \end{align}$

(3.22)

by noticing

$\begin{align*} \tau \sum_{n = 0}^{m} \|D_{1,\tau}\theta^{n+1}\|_{H^1}^2 \leq \frac 1 \tau \sum_{n = 0}^{m} \left( { \int}_{t_n}^{t_{n+1}} \|\theta_t(t)\|_{H^1}dt \right)^2 \leq C, \end{align*}$

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

For last term $\lambda_3^{n+1}$ , in terms of the following splitting:

$\begin{align*} (\lambda_{3}^{n+1}, e_{\theta}^{n+1}) & = ({{\textbf{u}}}_h^n\cdot\nabla R_h\theta^{n+1},e_{\theta}^{n+1}) +\frac{1}{2}((\nabla\cdot{{\textbf{u}}}_h^n)R_h\theta^{n+1},e_{\theta}^{n+1}) -({{\textbf{u}}}^{n+1}\cdot\nabla\theta^{n+1},e_{\theta}^{n+1})\\ & = -(e_h^n\cdot\nabla R_h\theta^{n+1},e_{\theta}^{n+1}) -({{\textbf{s}}}_h^n\cdot\nabla\eta_{\theta}^{n+1},e_{\theta}^{n+1}) +(({{\textbf{s}}}_h^n-{{\textbf{u}}}^{n+1})\cdot\nabla\theta^{n+1},e_{\theta}^{n+1})\\ &-\frac{1}{2}((\nabla\cdot e_h^n)R_h\theta^{n+1},e_{\theta}^{n+1}) +\frac{1}{2}(\nabla\cdot({{\textbf{s}}}_{h}^n-{{\textbf{u}}}^n)R_h\theta^{n+1},e_{\theta}^{n+1}){,} \end{align*}$

we estimate it by

$\begin{align*} (\lambda_{3}^{n+1}, e_{\theta}^{n+1}) \leq & ( \|e_h^n\|_{{{\textbf{L}}}^2} \|\nabla R_h\theta^{n+1}\|_{{{\textbf{L}}}^\infty} +\|{{\textbf{s}}}_{h}^{n}\|_{{{\textbf{L}}}^\infty}\|\nabla\eta_{\theta}^{n+1}\|_{{{\textbf{L}}}^2} + \|{{\textbf{s}}}_h^n-{{\textbf{u}}}^{n+1}\|_{{{\textbf{L}}}^3}\|\nabla\theta^{n+1}\|_{{{\textbf{L}}}^6}) \|e_{\theta}^{n+1}\|_{L^2} \\ +&\frac{1}{2}(\|\nabla\cdot e_h^n\|_{L^2}\|R_h\theta^{n+1}\|_{L^{\infty}} +\|{{\textbf{s}}}_h^n-{{\textbf{u}}}^n\|_{{{\textbf{H}}}^1}\|R_h\theta^{n+1}\|_{L^{\infty}})\|e_{\theta}^{n+1}\|_{L^2}\\ \leq & C (h^{4}+\tau^2) + C (\| e_{\theta}^{n+1}\|_{L^2}^2 + \|e_h^n\|_{{{\textbf{L}}}^2}^2) + \frac {\epsilon_1\beta }2 \|\nabla \cdot e_h^{n}\|_{L^2}^2 . \end{align*}$

Furthermore, we get

$\begin{align} \begin{split} 2\tau \sum\limits_{n = 0}^{m} (\lambda_{3}^{n+1}, e_\theta^{n+1}) \leq & C (h^{4}+\tau^2) + C\tau \sum\limits_{n = 0}^{m} (\| e_{\theta}^{n+1}\|_{L^2}^2 + \|e_h^n\|_{{{\textbf{L}}}^2}^2) + \epsilon_1\beta \tau \sum_{n = 0}^{m} \|\nabla \cdot e_h^{n}\|_{L^2}^2, \end{split} \end{align}$

(3.23)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Substituting the estimates (3.21)–(3.23) into (3.20), we get

$\begin{align} \begin{split} & \|e_{\theta}^{m+1}\|_{L^2}^{2}+2\kappa \tau\sum\limits_{n = 0}^{m} \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} \\ \leq & C(h^{4}+\tau^2) + C\tau \sum\limits_{n = 0}^{m} (\| e_{\theta}^{n+1}\|_{L^2}^2 + \|e_h^n\|_{{{\textbf{L}}}^2}^2) + \epsilon_1\beta \tau \sum_{n = 0}^{m} \|\nabla \cdot e_h^{n}\|_{L^2}^2, \end{split} \end{align}$

(3.24)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Taking the sum of (3.19) and (3.24), we have

$\begin{align} \begin{split} & \|e_{h}^{m+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{m+1}\|_{L^2}^{2}+2\tau\sum\limits_{n = 0}^{m} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{n+1}\|_{L^2}^{2} \right ) \\ \leq & C \tau \sum_{n = 0}^{m} (\|e_{\theta}^{n+1}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} +\|e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}) + C(h^{4}+\tau^2) \\ & + C h^{4} \tau \sum_{n = 0}^{m} (1+\|\theta_{h}^{n}\|_{L^\infty}^2) + 3\epsilon_1\beta \tau \sum_{n = 0}^{m} \|\nabla\cdot e_h^{n+1}\|_{L^2} + C\tau \sum_{n = 0}^{m} (1+\|\theta_{h}^{n}\|_{L^\infty}^2) \|e_{\theta}^{n}\|_{L^2}^2 , \end{split} \end{align}$

(3.25)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

To uniformly bound the norm $\|\theta_h^n\|_{L^\infty}$ , we use the method of mathematical induction. Taking $m = 0$ in (3.25) and noting $e_h^0 = 0$ , $e_\theta^0 = 0$ and $\|\theta_h^0\|_{L^\infty} = \|R_h\theta_0\|_{L^\infty}\leq \|\theta_0\|_{H^3}$ , we get

$\begin{align} \begin{split} & \|e_{h}^{1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{1}\|_{L^2}^{2}+2\tau \left ( \nu\|\nabla e_{h}^{1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{1}\|_{L^2}^{2} \right ) \\ \leq & C \tau (\|e_{\theta}^{1}\|_{L^2}^{2}+\|e_{h}^{1}\|_{{{\textbf{L}}}^2}^{2} ) + C(h^4+\tau^2) + 3\epsilon_1\beta \tau \|\nabla\cdot e_h^{1}\|_{L^2}. \end{split} \end{align}$

(3.26)

For sufficiently small $\tau$ with $C\tau < 1/2$ , we select small parameter $\epsilon_1\leq 2/3$ in (3.26). We then conclude that there exists some $C_1 > 0$ , which is independent of $h, \tau$ , $\nu$ and $\kappa$ such that

$\begin{align*} \|e_{h}^{1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{1}\|_{L^2}^{2}+\tau \left ( \nu\|\nabla e_{h}^{1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{1}\|_{{{\textbf{L}}}^2}^{2}\right ) \leq C_1^2 (h^4+\tau^2). \end{align*}$

Thus, the error estimate (3.11) holds for $m = 0$ by taking $C_0 > C_1$ . Now, we assume that the error estimate (3.11) holds for $m\leq k-1$ with $1\leq k\leq N-1$ . Under the time step condition $\tau\leq Ch^2$ , we have

$\begin{align} \|e_{\theta}^{n}\|_{L^2} \leq CC_0h^2 \qquad \forall \ 1\leq n\leq k. \end{align}$

(3.27)

By the inverse inequality, we get

$\begin{align*} \|e_{\theta}^{n}\|_{L^\infty} \leq Ch^{-3/2}\|e_{\theta}^{n}\|_{L^2} \leq CC_0h^{1/2} \leq C \end{align*}$

for sufficiently small $h$ with $C_0h^{1/2} \leq 1$ . Thus,

$\begin{align} \|\theta_h^{n}\|_{L^\infty} \leq \|e_\theta^n\|_{L^\infty} + \|\theta^n\|_{L^\infty}\leq C \quad {\mbox{and}}\quad \tau \sum_{n = 1}^{k} \| \theta_{h}^{n}\|_{L^\infty}^2 \leq C. \end{align}$

(3.28)

As a result, resetting $m$ by $k$ in (3.25), we have

$\begin{align} \begin{split} & \|e_{h}^{k+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{k+1}\|_{L^2}^{2}+2\tau\sum\limits_{n = 0}^{k} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{n+1}\|_{L^2}^{2} \right ) \\ \leq & C(h^4+\tau^2) + C \tau \sum_{n = 0}^{k} (\|e_{\theta}^{n+1}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} +\|e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}) \\ & + 3\epsilon_1\beta \tau \sum_{n = 0}^{k} \|\nabla\cdot e_h^{n+1}\|_{L^2} . \end{split} \end{align}$

(3.29)

We select small parameter $\epsilon_1\leq 2/3$ in (3.29), then

$\begin{align*} & \|e_{h}^{k+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{k+1}\|_{L^2}^{2}+\tau\sum\limits_{n = 0}^{k} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} \right ) \\ \leq & C(h^4+\tau^2) + C \tau \sum_{n = 0}^{k} (\|e_{\theta}^{n+1}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} +\|e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}). \end{align*}$

Applying the discrete Gronwall inequality to the above inequality, we conclude that there exists some $C_2 > 0$ , which is independent of $h, \tau$ , $\nu$ and $\kappa$ such that

$\begin{align*} \|e_{h}^{k+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{k+1}\|_{L^2}^{2}+\tau\sum\limits_{n = 0}^{k} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_{\theta}^{n+1}\|_{{{\textbf{L}}}^2}^{2} \right ) \leq C_2^2(h^4+\tau^2) \end{align*}$

for sufficiently small $\tau$ . Thus, we prove that the error estimate (3.11) also holds for $m = k$ by taking $C_0\geq C_2$ , and we finish the mathematical induction. □

Based on (3.11) in Theorem 3.1, we get the following error estimate for the first-order Euler grad-div finite element scheme (3.1)–(3.2).

Theorem 3.2 Under the assumptions in Theorem 3.1, we have

$\begin{align} \|{{\textbf{u}}}^{m+1}-{{\textbf{u}}}_h^{m+1}\|_{{{\textbf{L}}}^2}^{2}+\|\theta^{m+1}-\theta_h^{m+1}\|_{L^2}^{2} \leq C(h^{4}+\tau^2) \end{align}$

(3.30)

with $0\leq m\leq N-1$ , where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

4. The BDF2 grad-div stabilization finite element approximation

In this section, we consider the second-order BDF2 fully discrete scheme. Based on the extrapolation method and grad-div stabilization method, we propose the following linearized second-order finite element scheme for $1\leq n\leq N-1$ :

Step Ⅰ: Find $\theta_h^{n+1}\in Y_h$ by

$\begin{align} (D_{2,\tau}\theta_{h}^{n+1},\psi_{h})+\kappa(\nabla \theta_{h}^{n+1},\nabla\psi_{h})+c_2(\widehat{{{\textbf{u}}}}_{h}^{n},\theta_{h}^{n+1},\psi_{h}) = (g^{n+1},\psi_{h}) \quad \forall \ \psi_h\in Y_h. \end{align}$

(4.1)

Step Ⅱ: Find $({{\textbf{u}}}_{h}^{n+1}, \theta_{h}^{n+1})\in {{\textbf{V}}}_{h}\times Q_{h}$ by

$\begin{eqnarray} (D_{2,\tau}{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})+\nu(\nabla{{\textbf{u}}}_{h}^{n+1},\nabla {{\textbf{v}}}_{h})+c_1(\widehat {{{\textbf{u}}}}_{h}^{n},{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})-(\nabla\cdot{{\textbf{v}}}_{h},p_{h}^{n+1}) +(\nabla\cdot{{\textbf{u}}}_{h}^{n+1},q_{h})\\ +\beta(\nabla\cdot{{\textbf{u}}}_{h}^{n+1},\nabla\cdot{{\textbf{v}}}_{h})-(\gamma_{1}\widehat {\theta}_{h}^{n} + \gamma_2(\widehat \theta_{h}^{n})^{2} ,{\textbf{i}}_{3} \cdot {{\textbf{v}}}_{h}) = ({{\textbf{f}}}^{n+1},{{\textbf{v}}}_{h}) \quad \forall \ ({{\textbf{v}}}_{h},q_{h}) \in{{\textbf{V}}}_{h} \times Q_{h}. \end{eqnarray}$

(4.2)

Remark 4.1. Numerical solution $({{\textbf{u}}}_h^1, \theta_h^1)\in {{\textbf{V}}}_{0h}\times Y_h$ can be solved from the Euler finite element scheme (3.1)–(3.2) in the above section. Usually, since the one-step error is first higher order, then one has

$\begin{align} \|e_{h}^{1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{1}\|_{L^2}^{2}+\tau \left ( \nu\|\nabla e_{h}^{1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_{\theta}^{1}\|_{{{\textbf{L}}}^2}^{2} + \beta\|\nabla\cdot e_h^1\|_{L^2}^2 \right ) \leq C_3^2(h^{4}+\tau^4), \end{align}$

(4.3)

where $C_3 > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

By a similar proof for Lemma 3.1, we can get the following unconditional stabilities of the numerical scheme (4.1)–(4.2), which imply the existence and uniqueness of numerical solution $({{\textbf{u}}}_h^{n+1}, p_h^{n+1}, \theta_h^{n+1})$ to the BDF2 scheme (4.1)–(4.2).

Lemma 4.1. For $1\leq n\leq N-1$ and all $\tau > 0$ , $h > 0$ , the BDF2 finite element scheme (4.1) and (4.2) has a unique solution $\theta_{h}^{n+1}\in Y_{h}$ and $({{\textbf{u}}}_{h}^{n+1}, p_{h}^{n+1})\in {{\textbf{V}}}_{h}\times Q_{h}$ . Moreover, the discrete energy inequalities hold:

$\begin{align} \|\theta_{h}^{m+1}\|_{L^{2}}^{2}+\|\widehat {\theta}_{h}^{m+1}\|_{L^{2}}^{2} +2\kappa\tau\sum\limits_{n = 1}^{m}\|\nabla\theta_{h}^{n+1}\|_{{{\textbf{L}}}^{2}}^{2}\leq C\tau\sum\limits_{n = 1}^{m}\|g^{n+1}\|_{L^{2}}^{2} +C\left( \|\theta_{h}^{1}\|_{L^{2}}^{2} + \| {\theta}_{h}^{0}\|_{L^{2}}^{2} \right), \end{align}$

(4.4)

and

$\begin{align} \begin{split} &\|{{\textbf{u}}}_{h}^{m+1}\|_{{{\textbf{L}}}^{2}}^{2}+\|\widehat {{\textbf{u}}}_{h}^{m+1}\|_{{{\textbf{L}}}^{2}}^{2}+2\tau\nu\sum\limits_{n = 1}^{m}\|\nabla {{\textbf{u}}}_{h}^{n+1}\|_{{{\textbf{L}}}^{2}}^{2}+4\tau\beta\sum\limits_{n = 1}^{m}\|\nabla\cdot{{\textbf{u}}}_{h}^{n+1}\|^{2}_{{{\textbf{L}}}^2}\\ \leq & C\tau\sum\limits_{n = 1}^{m}(\|{{\textbf{f}}}^{n+1}\|_{{{\textbf{L}}}^{2}}^{2}+\|g^{n+1}\|_{L^{2}}^{2})+C(\|\theta_{h}^{1}\|_{L^{2}}^{2} +\|{\theta}_{h}^{0}\|_{L^{2}}^{2}) \\ & +C\left (\tau\sum\limits_{n = 1}^{m}\|g^{n+1}\|_{L^{2}}^{2}+\|\theta_{h}^{1}\|_{L^{2}}^{2} +\|{\theta}_{h}^{0}\|_{L^{2}}^{2} \right )^{2} \end{split} \end{align}$

(4.5)

for all $1\leq m\leq N-1$ , where $C > 0$ is independent of $h$ and $\tau$ .

In terms of the discrete inf-sup condition (2.5), we rewrite the discrete variational problem (4.2) as: find ${{\textbf{u}}}_{h}^{n+1}\in {{\textbf{V}}}_{0h}$ by

$\begin{eqnarray} (D_{2,\tau}{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})+\nu(\nabla{{\textbf{u}}}_{h}^{n+1},\nabla {{\textbf{v}}}_{h}) +c_1(\widehat {{{\textbf{u}}}}_{h}^{n},{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})+\beta(\nabla\cdot{{\textbf{u}}}_{h}^{n+1},\nabla\cdot{{\textbf{v}}}_{h}) -(\gamma_{1}\widehat {\theta}_{h}^{n}+\gamma_2 (\widehat \theta_{h}^{n})^2, {\textbf{i}}_{3}\cdot{{\textbf{v}}}_{h}) = ({{\textbf{f}}}^{n+1},{{\textbf{v}}}_{h}) \quad \forall \ {{\textbf{v}}}_h\in{{\textbf{V}}}_{0h}. \end{eqnarray}$

(4.6)

On the other hand, let ${{\textbf{s}}}_{h}^{0} = {{\textbf{u}}}_{h}^{0}$ . For $1\leq n\leq N-1$ , ${{\textbf{s}}}_{h}^{n+1}$ and $\theta^{n+1}$ satisfy

$\begin{align} \begin{split} &(D_{2,\tau}{{\textbf{s}}}_{h}^{n+1},{{\textbf{v}}}_{h})+\nu(\nabla {{\textbf{s}}}_{h}^{n+1},\nabla {{\textbf{v}}}_{h})+c_1(\widehat{{{\textbf{s}}}}_{h}^{n},{{\textbf{s}}}_{h}^{n+1},{{\textbf{v}}}_{h}) +\beta(\nabla\cdot {{\textbf{s}}}_{h}^{n+1},\nabla\cdot{{\textbf{v}}}_{h})\\ &-\gamma_{1}(\widehat {\theta}^{n},\textbf i_{3}\cdot{{\textbf{v}}}_{h}) -\gamma_{2}((\widehat \theta^{n})^{2},\textbf i_{3}\cdot{{\textbf{v}}}_{h})\\ & = ({{\textbf{f}}}^{n+1},{{\textbf{v}}}_{h})+(\nabla\cdot{{\textbf{v}}}_{h}, p^{n+1}-\pi_{h}p^{n+1})\\ &+\beta(\nabla\cdot({{\textbf{s}}}_{h}^{n+1}-{{\textbf{u}}}^{n+1}),\nabla\cdot{{\textbf{v}}}_{h}) +(D_{2,\tau} {{\textbf{s}}}_{h}^{n+1}-{{\textbf{u}}}_{t}(t_{n+1}),{{\textbf{v}}}_{h})\\ & +c_1(\widehat{{{\textbf{s}}}}_{h}^{n},{{\textbf{s}}}_{h}^{n+1},{{\textbf{v}}}_{h})-c_1({{\textbf{u}}}^{n+1},{{\textbf{u}}}^{n+1},{{\textbf{v}}}_{h})\\ &+\gamma_{1}(\theta^{n+1}-\widehat {\theta}^{n}, \textbf i_{3}\cdot{{\textbf{v}}}_{h}) +\gamma_{2}((\theta^{n+1})^{2}-(\widehat \theta^{n})^{2},\textbf i_{3}\cdot{{\textbf{v}}}_{h}) \quad \forall \ {{\textbf{v}}}_h\in{{\textbf{V}}}_{0h}, \end{split} \end{align}$

(4.7)

and

$\begin{eqnarray} (D_{2,\tau}\theta^{n+1},\psi_{h})+\kappa(\nabla \theta^{n+1},\nabla\psi_{h})+c_2({{\textbf{u}}}^{n+1},\theta^{n+1},\psi_{h}) = (g^{n+1},\psi_{h}) +(D_{2,\tau}\theta^{n+1}-\theta_{t}^{n+1},\psi_{h}) \quad \forall \ \psi_h\in Y_h. \end{eqnarray}$

(4.8)

In (4.7), we use $\nabla\cdot{{\textbf{u}}}^{n+1} = 0$ in $\Omega$ and $(\nabla\cdot{{\textbf{v}}}_h, \pi_hp^{n+1}) = 0$ .

Next, we give error equations corresponding to the BDF2 scheme. Subtracting (4.6) from (4.7) leads to

$\begin{align} (D_{2,\tau}e_{h}^{n+1},{{\textbf{v}}}_{h})+\nu(\nabla e_{h}^{n+1},\nabla{{\textbf{v}}}_{h}) +\beta(\nabla\cdot e_{h}^{n+1},\nabla\cdot{{\textbf{v}}}_{h}) = \sum_{j = 1}^6 (I_{j}^{n+1},{{\textbf{v}}}_{h}) \quad \forall \ {{\textbf{v}}}_h\in{{\textbf{V}}}_{0h}, \end{align}$

(4.9)

where

$\begin{eqnarray*} (I_{1}^{n+1},{{\textbf{v}}}_{h})& = &( D_{2,\tau}{{\textbf{s}}}_{h}^{n+1}-{{\textbf{u}}}_{t}(t_{n+1}),{{\textbf{v}}}_{h}) , \\ (I_{2}^{n+1},{{\textbf{v}}}_{h})& = &(\nabla\cdot{{\textbf{v}}}_{h},(p^{n+1}-\pi_h p^{n+1})+\beta\nabla\cdot({{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1})), \\ (I_{3}^{n+1},{{\textbf{v}}}_{h})& = &c_1(\widehat {{{\textbf{s}}}}_{h}^{n},{{\textbf{s}}}_{h}^{n+1},{{\textbf{v}}}_{h})-c_1({{\textbf{u}}}^{n+1},{{\textbf{u}}}^{n+1},{{\textbf{v}}}_{h}), \\ (I_{4}^{n+1},{{\textbf{v}}}_{h})& = &c_1(\widehat{{{\textbf{u}}}}_{h}^{n},{{\textbf{u}}}_{h}^{n+1},{{\textbf{v}}}_{h})-c_1(\widehat {{{\textbf{s}}}}_{h}^{n},{{\textbf{s}}}_{h}^{n+1},{{\textbf{v}}}_{h}),\\ (I_{5}^{n+1},{{\textbf{v}}}_{h})& = &\gamma_{1}(\theta^{n+1}-\widehat {\theta}^{n}, \textbf i_{3}\cdot{{\textbf{v}}}_{h}) +\gamma_{2}( (\theta^{n+1})^{2}- (\widehat \theta^{n})^{2}, \textbf i_{3}\cdot{{\textbf{v}}}_{h}),\\ (I_{6}^{n+1},{{\textbf{v}}}_{h})& = &\gamma_{1}((\widehat \theta^{n}-\widehat \theta_{h}^{n}),\textbf i_{3}\cdot{{\textbf{v}}}_{h})+ \gamma_2(((\widehat\theta^{n})^{2}-(\widehat\theta_{h}^{n})^{2}), \textbf i_{3}\cdot{{\textbf{v}}}_{h}). \end{eqnarray*}$

Subtracting (4.1) from (4.8) leads to

$\begin{align} (D_{2,\tau}e_{\theta}^{n+1},\psi_{h})+\kappa(\nabla e_{\theta}^{n+1},\nabla\psi_{h})+c_2(\widehat{{\textbf{u}}}_h^n,e_{\theta}^{n+1},\psi) = \sum_{j = 1}^{3}(X_{j}^{n+1},\psi_{h}) \quad \forall \ \psi_h\in Y_{h}, \end{align}$

(4.10)

where

$\begin{eqnarray*} (X_{1}^{n+1},\psi_{h})& = &( D_{2,\tau}\theta^{n+1}-\theta_t(t_{n+1}),\psi_{h}),\\ (X_{2}^{n+1},\psi_{h})& = &- (D_{2,\tau}\eta_\theta^{n+1},\psi_{h}),\\ (X_{3}^{n+1},\psi_{h})& = &c_2(\widehat {{\textbf{u}}}_{h}^{n},R_h\theta^{n+1},\psi_{h}) -c_2( {{\textbf{u}}}^{n+1},\theta^{n+1},\psi_{h}). \end{eqnarray*}$

The main result in this section is presented in the following theorem.

Theorem 4.1. Under the assumptions in (3.10), we further assume that

$\begin{align} {{\textbf{u}}}_{ttt} \in {{\textbf{L}}}^2((0, T]; {{\textbf{L}}}^2(\Omega)), \quad \theta_{ttt} \in L^2((0, T]; L^2(\Omega)). \end{align}$

(4.11)

Then, for the second-order BDF2 grad-div scheme (4.1)–(4.2), under the time step condition $\tau\leq Ch$ , when $h$ and $\tau$ are sufficiently small, we have the following error estimate:

(4.12)

with $0\leq m\leq N-1$ , where $\widehat C_0 > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Proof. Taking ${{\textbf{v}}}_h = 4\tau e_h^{n+1}$ in (4.9) and summing up the resulting equation from $n = 1$ to $n = m$ , we have

(4.13)

where $C > 0$ is independent of $h, \tau$ , $\nu, \kappa$ and we use

$\begin{align*} \|\widehat e_h^1\|_{{{\textbf{L}}}^2} \leq 2 \|e_h^1\|_{{{\textbf{L}}}^2} +\|e_h^0\|_{{{\textbf{L}}}^2} = 2\|e_h^1\|_{{{\textbf{L}}}^2} . \end{align*}$

We estimate the right-hand side of (4.13) as follows. Since

$\begin{align*} I_1^{n+1} = ( D_{2,\tau}{{\textbf{s}}}_{h}^{n+1}-D_{2,\tau}{{\textbf{u}}}^{n+1})+ (D_{2,\tau}{{\textbf{u}}}^{n+1} -{{\textbf{u}}}_{t}(t_{n+1})) , \end{align*}$

by using the Taylor's formula and the regularity assumption (4.11), one has

$\begin{align*} D_{2,\tau}{{\textbf{u}}}^{n+1} -{{\textbf{u}}}_{t}(t_{n+1}) = \frac{1}{2\tau}\int_{t_{n-1}}^{t^{n+1}} \left( 2(t-t_{n})^{2}-\frac{1}{2}(t-t_{n-1})^{2}\right)\partial_{ttt}{{\textbf{u}}}(t) dt, \end{align*}$

which results in

$\begin{align*} \tau\sum\limits_{n = 1}^{m} \|D_{2,\tau}{{\textbf{u}}}^{n+1} -{{\textbf{u}}}_{t}(t_{n+1})\|_{{{\textbf{L}}}^2}^2 \leq C\tau^4, \end{align*}$

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ . From the Hölder inequality, we have

$\begin{align*} &\quad \|D_{2,\tau}{{\textbf{s}}}_{h}^{n+1}-D_{2,\tau}{{\textbf{u}}}^{n+1}\|_{{{\textbf{L}}}^2} \\ &\leq \frac {3}{2\tau}\int_{t_{n}}^{t^{n+1}} \|\partial_{t}({{\textbf{s}}}_{h}(t)-{{\textbf{u}}}(t))\|_{{{\textbf{L}}}^2}dt+ \frac{1}{2\tau}\int_{t_{n-1}}^{t^{n}} \|\partial_{t}({{\textbf{s}}}_{h}(t)-{{\textbf{u}}}(t))\|_{{{\textbf{L}}}^2}dt\\ &\leq \frac {Ch^2}{\tau} \int_{t_{n-1}}^{t^{n+1}} \|\partial_{t}{{\textbf{u}}}(t)\|_{{{\textbf{H}}}^2}dt\\ &\leq \frac {Ch^2}{\tau^{1/2}} \left( { \int}_{t_n-1}^{t_{n+1}} \| \partial_t{{\textbf{u}}}(t)\|_{{{\textbf{H}}}^2}^2 dt \right)^{1/2}. \end{align*}$

Thus, we can obtain

$\begin{align} 4\tau\sum\limits_{n = 1}^{m} (I_{1}^{n+1},e_{h}^{n+1}) \leq C \tau \sum_{n = 1}^{m} \| e_h^{n+1}\|_{{{\textbf{L}}}^2}^{2} + C(\tau^4+h^{4}), \end{align}$

(4.14)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

From the Hölder inequality and (2.6), we estimate $(I_{2}^{n+1}, e_{h}^{n+1})$ by

$\begin{align*} (I_{2}^{n+1},e_{h}^{n+1}) \leq & \left( \|p^{n+1}-\pi_h p^{n+1}\|_{L^2} + \beta \|{{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1}\|_{{{\textbf{H}}}^1} \right) \|\nabla\cdot e_h^{n+1}\|_{L^2} \\ \leq & C h^2 \|\nabla\cdot e_h^{n+1}\|_{L^2}. \end{align*}$

Then, by the Young inequality, it is easy to see that

$\begin{align} 4\tau \sum\limits_{n = 1}^{m} (I_{2}^{n+1},e_{h}^{n+1}) \leq Ch^{4} + \epsilon_2\beta \tau \sum_{n = 1}^{m} \|\nabla\cdot e_h^{n+1}\|_{L^2}^{2} , \end{align}$

(4.15)

where $C > 0$ is independent of $h, \tau, \nu, \kappa$ and $\epsilon_2 > 0$ is a small constant determined later.

For the term of $I_3^{n+1}$ , according to the Hölder inequality, (2.6), (2.18), the regularity assumption (3.10) and the Taylor formula with integral remainder

$\begin{align} \widehat {{{\textbf{u}}}}^{n}-{{\textbf{u}}}^{n+1} = \int_{t_{n-1}}^{t_{n+1}}\{2(t-t_{n})_+ -(t-t_{n-1})\}\partial_{tt}{{\textbf{u}}}(t) dt, \end{align}$

(4.16)

where $(t-t_{n})_+ = \max\{t-t_n, 0\}$ , there holds

$\begin{align*} (I_{3}^{n+1},e_{h}^{n+1}) & = c_1(\widehat {{{\textbf{s}}}}_{h}^{n}-\widehat {{{\textbf{u}}}}^{n},{{\textbf{s}}}_{h}^{n+1},e_{h}^{n+1}) +c_1(\widehat {{{\textbf{u}}}}^{n}-{{\textbf{u}}}^{n+1},{{\textbf{s}}}_{h}^{n+1},e_{h}^{n+1}) \\ &\quad +c_1({{\textbf{u}}}^{n+1},{{\textbf{s}}}_{h}^{n+1}-{{\textbf{u}}}^{n+1},e_{h}^{n+1})\\ &\leq ( \|\widehat {{{\textbf{s}}}}_{h}^{n}-\widehat {{{\textbf{u}}}}^{n}\|_{{{\textbf{H}}}^1} + \|\widehat {{{\textbf{u}}}}^{n}-{{\textbf{u}}}^{n+1}\|_{{{\textbf{L}}}^2} )\|{{\textbf{s}}}_h^{n+1}\|_{{{\textbf{W}}}^{1,\infty}} \|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}\\ & \quad + C \|{{\textbf{u}}}^{n+1}\|_{{{\textbf{L}}}^{\infty}} \| {{\textbf{s}}}_h^{n+1}-{{\textbf{u}}}^{n+1}\|_{{{\textbf{H}}}^1} \|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}\\ &\leq C(h^2+(\tau^{3}\int_{t_{n-1}}^{t^{n+1}}\|\partial_{tt}{{\textbf{u}}}(t)\|_{{{\textbf{L}}}^2}^{2}dt)^{\frac{1}{2}})\|e_{h}^{n+1}\|_{{{\textbf{L}}}^{2}}. \end{align*}$

Thus, we have

$\begin{align} 4\tau \sum\limits_{n = 1}^{m} (I_{3}^{n+1},e_{h}^{n+1}) \leq C (h^{4}+\tau^4) + C \tau \sum_{n = 1}^{m} \| e_h^{n+1}\|_{{{\textbf{L}}}^2}^2, \end{align}$

(4.17)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

A similar method leads to

$\begin{align*} (I_{4}^{n+1},e_{h}^{n+1}) & = c_1 (\widehat e_{h}^{n},{{\textbf{s}}}_{h}^{n+1},e_{h}^{n+1}) \\ & \leq \|{\widehat e}_{h}^{n}\|_{{{\textbf{L}}}^2}\|\nabla{{\textbf{s}}}_{h}^{n+1}\|_{{{\textbf{L}}}^\infty}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2} +\|{{\textbf{s}}}_{h}^{n+1}\|_{{{\textbf{L}}}^\infty}\|\nabla\cdot \widehat e_{h}^{n}\|_{L^2}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2} \\ & \leq \|{\widehat e}_{h}^{n}\|_{{{\textbf{L}}}^2}\|\nabla{{\textbf{s}}}_{h}^{n+1}\|_{{{\textbf{L}}}^\infty}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2} +\|{{\textbf{s}}}_{h}^{n+1}\|_{{{\textbf{L}}}^\infty}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2} (2\|\nabla\cdot e_{h}^{n}\|_{L^2}+\|\nabla\cdot e_{h}^{n-1}\|_{L^2}) \end{align*}$

and

$\begin{align} \begin{split} 4\tau \sum\limits_{n = 1}^{m} (I_{4}^{n+1},e_{h}^{n+1}) \leq & \epsilon_2\beta \tau \sum_{n = 1}^{m} (\|\nabla\cdot e_h^{n}\|_{L^2}^{2} + \|\nabla\cdot e_h^{n-1}\|_{L^2}^{2} ) \\ & + C\tau \sum_{n = 1}^{m}(\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}+\|\widehat e_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}), \end{split} \end{align}$

(4.18)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

For $I_5^{n+1}$ , using (4.16), the Hölder inequality and the regularity assumption (3.10), we have

$\begin{align*} (I_{5}^{n+1},e_{h}^{n+1}) &\leq C \|e_h^{n+1}\|_{{{\textbf{L}}}^2} (\tau^{3}\int_{t_{n-1}}^{t^{n+1}}\|\partial_{tt}\theta(t)\|_{L^2}^{2}dt)^{\frac{1}{2}}, \end{align*}$

{then}

$\begin{align} 4\tau \sum\limits_{n = 1}^{m} (I_{5}^{n+1},e_{h}^{n+1}) \leq C \tau \sum_{n = 1}^{m} \| e_h^{n+1}\|_{{{\textbf{L}}}^2}^{2} + C\tau^4, \end{align}$

(4.19)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

We estimate the last term $I_6^{n+1}$ by

$\begin{align*} (I_{6}^{n+1},e_{h}^{n+1}) &\leq \|\widehat \eta_{\theta}^{n}+\widehat e_{\theta}^{n}\|_{L^2}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}+ \|(\widehat\eta_{\theta}^{n}+\widehat e_{\theta}^{n})\|_{L^{2}} \|\widehat\theta^{n}+\widehat\theta_{h}^{n}\|_{L^{\infty}}\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}\\ &\leq C (1 +\|\widehat \theta_{h}^{n}\|_{L^\infty}^2) ( h^{4}+\|\widehat e_{\theta}^{n}\|_{L^2}^2) +C ( \|\widehat e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}) , \end{align*}$

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ . Taking the sum of the above inequality gives

$\begin{align} \begin{split} 4\tau \sum\limits_{n = 0}^{m} (I_{6}^{n+1},e_{h}^{n+1}) &\leq C h^{4} \tau \sum_{n = 1}^{m} (1+\|\widehat \theta_{h}^{n}\|_{L^\infty}^2) + C \tau \sum_{n = 1}^{m} ( \|\widehat e_{\theta}^{n}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}) \\ &+C\tau \sum_{n = 1}^{m} (1+\|\widehat \theta_{h}^{n}\|_{L^\infty}^2) ( \| e_{\theta}^{n}\|_{L^2}^2 +\| e_{\theta}^{n-1}\|_{L^2}^2 ), \end{split} \end{align}$

(4.20)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Substituting the estimates (4.14)–(4.20) into (4.13), we obtain

$\begin{align} \begin{split} & \|e_{h}^{m+1}\|_{{{\textbf{L}}}^2}^{2}+\|\widehat e_{h}^{m+1}\|_{{{\textbf{L}}}^2}^{2}+4\tau\sum\limits_{n = 1}^{m} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{n+1}\|_{L^2}^{2} \right ) \\ \leq & C(h^{4}+\tau^4) + C h^{4} \tau \sum_{n = 1}^{m} \|\widehat \theta_{h}^{n}\|_{L^\infty}^2 + C \tau \sum_{n = 1}^{m} ( \|\widehat e_{\theta}^{n}\|_{L^2}^{2}+\|\widehat e_{h}^{n}\|_{{{\textbf{L}}}^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}) \\ & +3\epsilon_2\beta \tau \sum_{n = 1}^{m} \|\nabla\cdot e_h^{n+1}\|_{L^2}^2 + C\tau \sum_{n = 1}^{m} (1+\|\widehat\theta_{h}^{n}\|_{L^\infty}^2) ( \|e_{\theta}^{n}\|_{L^2}^2 +\| e_{\theta}^{n-1}\|_{L^2}^2 ), \end{split} \end{align}$

(4.21)

where $C > 0$ is independent of $h, \tau$ , $\nu, \kappa$ and the estimate (4.3) is used.

Next, we estimate $e_\theta^{n+1}$ . Taking $\psi_h = 4\tau e_\theta^{n+1}$ in (4.10), summing up the resulting equation from $n = 1$ to $n = m$ and noticing the skew-symmetric property (2.1), we have

$\begin{align} \|e_{\theta}^{m+1}\|_{L^2}^{2}+\|\widehat e_{\theta}^{m+1}\|_{L^2}^{2}+4\kappa \tau\sum\limits_{n = 1}^{m} \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} \leq 4\tau\sum\limits_{n = 1}^{m} \sum_{j = 1}^3 (X_{j}^{n+1}, e_\theta^{n+1})+C\|e_{\theta}^1\|_{L^2}^{2} , \end{align}$

(4.22)

where we noted $e_{\theta}^0 = 0$ .

By the regularity assumptions (3.10) and (4.11), it is easy to prove that

$\begin{align*} (X_{1}^{n+1}, e_\theta^{n+1}) & \leq \frac{1}{2\tau}\int_{t_{n-1}}^{t^{n+1}} \left( 2(t-t_{n})^{2}-\frac{1}{2}(t-t_{n-1})^{2}\right) \|\partial_{ttt}\theta(t)\|_{L^2} dt\|e_{\theta}^{n+1}\|_{L^2},\\ \tau \sum_{n = 1}^{m} \|D_{2,\tau}\theta^{n+1}\|_{H^1}^2 & \leq \frac C \tau \sum_{n = 1}^{m} \left( { \int}_{t_{n-1}}^{t_{n+1}} \|\theta_t(t)\|_{H^1}dt \right)^2 \leq C{,} \end{align*}$

then we have

$\begin{align} 4\tau\sum\limits_{n = 1}^{m} (X_{1}^{n+1}, e_\theta^{n+1})&\leq C \tau \sum_{n = 1}^{m} \| e_\theta^{n+1}\|_{L^2} + C\tau^4 , \end{align}$

(4.23)

$\begin{align} \begin{split} 4\tau\sum\limits_{n = 1}^{m} (X_{2}^{n+1}, e_\theta^{n+1}) & \leq C \tau \sum_{n = 1}^{m} \| e_{\theta}^{n+1}\|_{L^2}^2 + Ch^{4} \tau \sum_{n = 1}^{m} \|D_{2,\tau}\theta^{n+1}\|_{H^1}^2 \\ & \leq C \tau \sum_{n = 1}^{m} \| e_{\theta}^{n+1}\|_{L^2}^2 + Ch^{4}, \end{split} \end{align}$

(4.24)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

For last term $X_3^{n+1}$ , we estimate it by

$\begin{align*} (X_{3}^{n+1}, e_{\theta}^{n+1}) & = (\widehat{{{\textbf{u}}}}_h^n\cdot\nabla R_h\theta^{n+1},e_{\theta}^{n+1}) +\frac{1}{2}((\nabla\cdot\widehat{{{\textbf{u}}}}_h^n)R_h\theta^{n+1},e_{\theta}^{n+1}) -({{\textbf{u}}}^{n+1}\cdot\nabla\theta^{n+1},e_{\theta}^{n+1})\\ & = -(\widehat{e}_h^n\cdot\nabla R_h\theta^{n+1},e_{\theta}^{n+1}) -(\widehat{{{\textbf{s}}}}_h^n\cdot\nabla\eta_{\theta}^{n+1},e_{\theta}^{n+1}) +((\widehat{{{\textbf{s}}}}_h^n-\widehat{{{\textbf{u}}}}^{n})\cdot\nabla\theta^{n+1},e_{\theta}^{n+1})\\ &+((\widehat{{{\textbf{u}}}}^n-{{\textbf{u}}}^{n+1})\cdot\nabla\theta^{n+1},e_{\theta}^{n+1})\\ &-\frac{1}{2}((\nabla\cdot \widehat{e}_h^n)R_h\theta^{n+1},e_{\theta}^{n+1}) +\frac{1}{2}(\nabla\cdot(\widehat{{{\textbf{s}}}}_{h}^n-\widehat{{{\textbf{u}}}}^n)R_h\theta^{n+1},e_{\theta}^{n+1}) , \end{align*}$

where the skew-symmetric property (2.1) is used. We have

$\begin{align} \begin{split} 4\tau \sum\limits_{n = 1}^{m} (X_{3}^{n+1}, e_\theta^{n+1}) &\leq (\|\widehat{e}_h^n\|_{{{\textbf{L}}}^2}\|\nabla R_h\theta^{n+1}\|_{{{\textbf{L}}}^\infty} +\|\widehat{{{\textbf{s}}}}_h^n\|_{{{\textbf{L}}}^{\infty}}\|\nabla\eta_{\theta}^{n+1}\|_{{{\textbf{L}}}^2})\|e_{\theta}^{n+1}\|_{L^2}\\ &+(\|\widehat{{{\textbf{s}}}}_h^n-\widehat{{{\textbf{u}}}}^{n}\|_{{{\textbf{L}}}^2}+\|\widehat{{{\textbf{u}}}}^n-{{\textbf{u}}}^{n+1}\|_{{{\textbf{L}}}^2})\|\nabla\theta^{n+1}\|_{{{\textbf{L}}}^\infty}\|e_{\theta}^{n+1}\|_{L^2}\\ &+(\|\nabla\cdot \widehat{e}_h^n\|_{L^2}+\|\widehat{{{\textbf{s}}}}_{h}^n-\widehat{{{\textbf{u}}}}^n\|_{{{\textbf{H}}}^1})\|R_h\theta^{n+1}\|_{L^\infty}\|e_{\theta}^{n+1}\|_{L^2}\\ &\leq C(h^{4}+\tau^4) + \epsilon_2\beta \tau \sum_{n = 1}^{m} \|\nabla \cdot e_h^{n}\|_{L^2}^2\\ & + C\tau \sum\limits_{n = 1}^{m} (\| e_{\theta}^{n+1}\|_{L^2}^2 + \|\widehat e_h^n\|_{{{\textbf{L}}}^2}^2), \end{split} \end{align}$

(4.25)

where the $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Substituting the estimates (4.23)–(4.25) into (4.22) and using (4.3), we get

$\begin{align} \begin{split} & \|e_{\theta}^{m+1}\|_{L^2}^{2}+4\kappa \tau\sum\limits_{n = 1}^{m} \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} \\ \leq & C(h^{4}+\tau^4) + C\tau \sum\limits_{n = 1}^{m} (\| e_{\theta}^{n+1}\|_{L^2}^2 + \|\widehat e_h^n\|_{{{\textbf{L}}}^2}^2) +\epsilon_2\beta \tau \sum_{n = 1}^{m} \|\nabla \cdot e_h^{n}\|_{L^2}^2 , \end{split} \end{align}$

(4.26)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Taking the sum of (4.21) and (4.26), we have

$\begin{align} \begin{split} & \|e_{h}^{m+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{m+1}\|_{L^2}^{2}+4\tau\sum\limits_{n = 1}^{m} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{n+1}\|_{L^2}^{2} \right ) \\ & \leq C(h^{4}+\tau^4) + C \tau \sum_{n = 1}^{m} (\|e_{\theta}^{n+1}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \|\widehat e_h^n\|_{{{\textbf{L}}}^2}^2+\|\widehat e_\theta^n\|_{L^2}^2) + C h^{4} \tau \sum_{n = 1}^{m} \|\widehat \theta_{h}^{n}\|_{L^\infty}^2\\ & +4\epsilon_2\beta \tau \sum_{n = 1}^{m} \|\nabla\cdot e_h^{n+1}\|_{L^2}^2 + C\tau \sum_{n = 1}^{m} (1+\|\widehat \theta_{h}^{n}\|_{L^\infty}^2) ( \| e_{\theta}^{n}\|_{L^2}^2+ \| e_{\theta}^{n-1}\|_{L^2}^2), \end{split} \end{align}$

(4.27)

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Next, we bound uniformly the norm $\|\theta_h^n\|_{L^\infty}$ by the method of mathematical induction as in the proof of Theorem 3.1.

According to (4.3), we can see that (4.12) holds for $m = 0$ if we take $\widehat C_0\geq C_3$ . Now, we assume that (4.12) also holds for $m\leq k-1$ with $1\leq k\leq N-1$ . Under the condition $\tau\leq Ch$ and the inverse inequality, one has

$\begin{align*} \|e_\theta^n\|_{L^\infty}\leq C h^{-3/2}\|e_\theta^n\|_{L^2} \leq C\widehat C_0 h^{1/2} \leq C \quad \forall \ 1\leq n\leq k \end{align*}$

for sufficiently small $h$ with $\widehat C_0 h^{1/2}\leq 1$ , which further gives

$\begin{align} \|\widehat \theta_h^n\|_{L^\infty} \leq \|\widehat e_\theta^n\|_{L^\infty} + \| \theta^n\|_{L^\infty} \leq C \quad {\mbox{and}} \quad \tau \sum_{n = 1}^{k} \|\widehat \theta_{h}^{n}\|_{L^\infty}^2 \leq C. \end{align}$

(4.28)

To finish the mathematical induction, we need to prove that (4.12) also holds for $m\leq k$ . As a result of (4.28), we have

$\begin{align*} & \|e_{h}^{k+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{k+1}\|_{L^2}^{2}+4\tau\sum\limits_{n = 1}^{k} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \beta \|\nabla\cdot e_{h}^{n+1}\|_{L^2}^{2} \right ) \\ & \leq C(h^4+\tau^4) + C \tau \sum_{n = 1}^{k} (\|e_{\theta}^{n+1}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \|\widehat e_h^n\|_{{{\textbf{L}}}^2}^2+\|\widehat e_\theta^n\|_{L^2}^2) \\ & +4\epsilon_3\beta \tau \sum_{n = 1}^{k} \|\nabla\cdot e_h^{n+1}\|_{L^2}^2 \end{align*}$

by setting $m = k$ in (4.27). Select small $\epsilon_2\leq 1$ , then we have

$\begin{align*} & \|e_{h}^{k+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{k+1}\|_{L^2}^{2}+\tau\sum\limits_{n = 1}^{k} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_\theta^{n+1}\|_{{{\textbf{L}}}^2}^{2} \right ) \\ \leq & C(h^4+\tau^4) + C \tau \sum_{n = 1}^{k} (\|e_{\theta}^{n+1}\|_{L^2}^{2}+\|e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2}+ \|\widehat e_h^n\|_{{{\textbf{L}}}^2}^2), \end{align*}$

where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

Applying the discrete Grönwall inequality to the above inequality, we conclude that there exists some $C_4 > 0$ , which is independent of $h, \tau$ , $\nu$ and $\kappa$ , such that

$\begin{align*} \|e_{h}^{k+1}\|_{{{\textbf{L}}}^2}^{2}+\|e_{\theta}^{k+1}\|_{L^2}^{2}+\tau\sum\limits_{n = 1}^{k} \left ( \nu\|\nabla e_{h}^{n+1}\|_{{{\textbf{L}}}^2}^{2} + \kappa \|\nabla e_{\theta}^{n+1}\|_{{{\textbf{L}}}^2}^{2} \right ) \leq C_4^2(h^4+\tau^4) \end{align*}$

for sufficiently small $\tau$ . Thus, we prove that the error estimate (4.12) also holds for $m\leq k$ by taking $\widehat C_0\geq C_4$ , and we finish the mathematical induction.□

Based on (4.12) in Theorem 4.1, we get the following error estimate for the second-order BDF2 grad-div finite element scheme (4.1)–(4.2).

Theorem 4.2 Under the assumptions in Theorem 4.1, we have

$\begin{align} \|{{\textbf{u}}}^{m+1}-{{\textbf{u}}}_h^{m+1}\|_{{{\textbf{L}}}^2}^{2}+\|\theta^{m+1}-\theta_h^{m+1}\|_{L^2}^{2} \leq C(h^{4}+\tau^4) \end{align}$

(4.29)

with $0\leq m\leq N-1$ , where $C > 0$ is independent of $h, \tau$ , $\nu$ and $\kappa$ .

5. Numerical results

In this section, we present numerical results of the second-order BDF2 grad-div stabilization method for the penetrative convection problems (1.1)–(1.5) with small viscosity coefficient $\nu$ . For simplicity, we only consider the 2D problem with the domain $\Omega = [0, 1]^2$ . To check the convergence rates, we take the appropriate ${{\textbf{f}}}$ and $g$ such that the exact solution $({{\textbf{u}}}, p, \theta)$ is given by

${{\textbf{u}}}_1 = 10 x^2(x-1)^2 y(y-1)(2y-1)\exp(-t),\quad {{\textbf{u}}}_2 = -10 x(x-1)(2x-1) y^2 (y-1)^2\exp(-t) ,$

$p = 10(2x-1)(2y-1)\exp(-t), \quad \theta = \sin(\pi x)\sin(\pi y)\exp(-t).$

In numerical computation, we take the physical parameters $\kappa = 10^{-1}, \gamma_{1} = 10^{-1}, \gamma_{2} = 10^{-1}$ , the stabilized parameter $\beta = 0.1$ and the final time $T = 1$ . In addition, we denote

$\begin{eqnarray*} \|{{\textbf{u}}}-{{\textbf{u}}}_h\|_{{{\textbf{L}}}^2} & = & \|{{\textbf{u}}}^N-{{\textbf{u}}}_h^N\|_{{{\textbf{L}}}^2}, \quad \|\theta-\theta_h\|_{L^2} = \|\theta^N-\theta_h^N\|_{L^2} , \\ \|{{\textbf{u}}}-{{\textbf{u}}}_h\|_{l^2(V)} & = & \left( \tau \sum_{n = 1}^N \|\nabla({{\textbf{u}}}^n-{{\textbf{u}}}^n_h)\|_{{{\textbf{L}}}^2}^2 \right)^{1/2}, \\ \|\theta-\theta_h\|_{l^2(Y)} & = & \left( \tau \sum_{n = 1}^N \|\nabla(\theta^n-\theta^n_h)\|_{{{\textbf{L}}}^2}^2 \right)^{1/2}. \end{eqnarray*}$

In terms of error estimates in Theorems 4.1 and 4.2, we have the second-order convergence rates of the velocity and temperature if we select $\tau = \mathcal O(h)$ . We use a uniform mesh on $\Omega$ with the spatial mesh size $h$ in each direction. By taking gradually decreasing mesh sizes $h = 1/4, 1/8, \cdots, 1/128$ , we choose different iteration numbers $N = 1/h$ such that the time step size $\tau = h$ . We present numerical results with small viscosity coefficients $\nu = 10^{-3}$ and $\nu = 10^{-4}$ in Tables 1 and 2, respectively, from which we can see that the second-order convergence rates are reached. Thus, these numerical results are in good agreement with the theoretical analysis.

Table 1. Numerical errors and convergence rates with

$\nu = 10^{-3}$ and

$\tau = h$ .

$h$	$\\|{{\textbf{u}}}-{{\textbf{u}}}_{h}\\|_{{{\textbf{L}}}^2}$	rate	$\\|{{\textbf{u}}}-{{\textbf{u}}}_h\\|_{l^2(V)}$	rate	$\\|\theta-\theta_{h}\\|_{L^2}$	rate	$\\|\theta-\theta_{h}\\|_{l^2(Y)}$	rate
1/4	4.28913E-03	$\quad$	6.57526E-02	$\quad$	2.09701E-03	$\quad$	4.88982E-02	$\quad$
1/8	9.74180E-04	2.14	2.16564E-02	1.60	3.44116E-04	2.61	1.23962E-02	1.98
1/16	2.34556E-04	2.05	4.99515E-03	2.12	6.70734E-05	2.36	3.11457E-03	1.99
1/32	5.84811E-05	2.00	1.03426E-03	2.27	1.51182E-05	2.15	7.79712E-04	2.00
1/64	1.46777E-05	1.99	2.10634E-04	2.30	3.61134E-06	2.07	1.94978E-04	2.00
1/128	3.67926E-06	2.00	4.68535E-05	2.17	8.84215E-07	2.03	4.87451E-05	2.00

| Show Table

DownLoad: CSV

Table 2. Numerical errors and convergence rates with

$\nu = 10^{-4}$ and

$\tau = h$ .

$h$	$\\|{{\textbf{u}}}-{{\textbf{u}}}_{h}\\|_{{{\textbf{L}}}^2}$	rate	$\\|{{\textbf{u}}}-{{\textbf{u}}}_h\\|_{l^2(V)}$	rate	$\\|\theta-\theta_{h}\\|_{L^2}$	rate	$\\|\theta-\theta_{h}\\|_{l^2(Y)}$	rate
1/4	4.75610E-03	$\quad$	7.67595E-02	$\quad$	2.10166E-03	$\quad$	4.89096E-02	$\quad$
1/8	1.20543E-03	1.98	3.84952E-02	1.00	3.45921E-04	2.60	1.23987E-02	1.98
1/16	2.87010E-04	2.07	1.36084E-02	1.50	6.79439E-05	2.35	3.11550E-03	1.99
1/32	6.89986E-05	2.06	3.67556E-03	1.89	1.54451E-05	2.14	7.80036E-04	2.00
1/64	1.73110E-05	1.99	8.06317E-04	2.19	3.70433E-06	2.06	1.95075E-04	2.00
1/128	4.24985E-06	2.03	1.45108E-04	2.47	9.07167E-07	2.03	4.87662E-05	2.00

| Show Table

DownLoad: CSV

On the other hand, the second-order convergence rates in the $L^2$ norm is sub-optimal since we use the $P_2$ element to approximate ${{\textbf{u}}}$ and $\theta$ . Theoretically, the optimal error estimate should be

$\begin{align} \|{{\textbf{u}}}-{{\textbf{u}}}_{h}\|_{{{\textbf{L}}}^2}+\|\theta-\theta_{h}\|_{L^2} \leq C (h^3+\tau^2) \end{align}$

(5.1)

for $1\leq n\leq N$ . To check (5.1), we take $\tau = h^{3/2}$ such that

$\begin{align} \|{{\textbf{u}}}-{{\textbf{u}}}_{h}\|_{{{\textbf{L}}}^2}+\|\theta-\theta_{h}\|_{L^2} \leq C h^3. \end{align}$

(5.2)

By taking different mesh sizes $h = 1/2^2, 1/3^2, \cdots, 1/7^2$ , we give numerical results in and with $\nu = 10^{-3}$ and $\nu = 10^{-4}$ , respectively. We can see that the almost third-order convergence rates in $L^2$ norm are reached. Thus, the $L^2$ error estimates in (4.29) is not optimal. How to prove the optimal error estimates in $L^2$ norm will be reported in future work.

Table 3. Numerical errors and convergence rates with

$\nu = 10^{-3}$ and

$\tau = h^{3/2}$ .

$h$	$\\|{{\textbf{u}}}-{{\textbf{u}}}_{h}\\|_{{{\textbf{L}}}^2}$	rate	$\\|\theta-\theta_{h}\\|_{L^2}$	rate
1/4	3.18961E-03	$\quad$	1.39164E-03	$\quad$
1/9	3.03406E-04	2.90	1.22372E-04	3.00
1/16	5.02327E-05	3.13	2.18932E-05	2.99
1/25	1.15939E-05	3.29	5.75192E-06	3.00
1/36	3.31503E-06	3.43	1.93196E-06	2.99
1/49	1.12367E-06	3.51	7.67713E-07	2.99

| Show Table

DownLoad: CSV

Table 4. Numerical errors and convergence rates with

$\nu = 10^{-4}$ and

$\tau = h^{3/2}$ .

$h$	$\\|{{\textbf{u}}}-{{\textbf{u}}}_{h}\\|_{{{\textbf{L}}}^2}$	rate	$\\|\theta-\theta_{h}\\|_{L^2}$	rate
1/4	3.65795E-03	$\quad$	1.39237E-03	$\quad$
1/9	5.23458E-04	2.40	1.22415E-04	3.00
1/16	1.16969E-04	2.60	2.18917E-05	2.99
1/25	3.32275E-05	2.82	6.22161E-06	2.82
1/36	1.12122E-05	2.98	2.25287E-06	2.79
1/49	4.52035E-06	2.95	9.47731E-07	2.81

| Show Table

DownLoad: CSV

6. Conclusions

In this paper, we studied the first-order Euler and second-order BDF2 finite element schemes for the approximation of the the time-dependent penetrative convection equations. In designing of numerical scheme, we used the grad-div stabilization method to overcome instabilities from the high Reynolds number. Main advantages of the proposed schemes are two folds. One is that they both are unconditionally stable without any condition of the time step and mesh size. Another is that one only needs to solve linearized systems at each time step. In terms of the analysis technique in ^[18,19], uniform error estimates in $L^2$ norm were derived in which the constants are independent of inverse powers of the viscosity coefficient and thermal conductivity coefficient.

In some related models, the nonlinear term $2\nu\mathcal D({{\textbf{u}}}):\mathcal D({{\textbf{u}}})$ , which is used to describe the viscous dissipation, appears in the temperature equation (1.3) ^[2,9,47]. Here, $\mathcal D({{\textbf{u}}})$ is the strain tensor having components $\mathcal D_{ij}({{\textbf{u}}}) = \frac 12 (\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i})$ . However, this strongly nonlinear term will result in much more difficulties in well-posedness analysis and numerical analysis. Thus, we consider the Eq (1.3) without the viscous dissipation in this paper. On the other hand, different modified Boussinesq approximations for some practical problems have been considered, such as the heat and mass transfer models in ^[10,48]. In future works, we can use the grad-div stabilization method for these modified models.

On the other hand, for the reason of simplicity, we use the zero Dirichlet boundary condition for the temperature. We remark that uniform error estimates derived in this paper can be extended to the problem with the mixed boundary condition for the temperature in ^[19], where the mixed boundary condition is composed by zero Dirichlet and Neumann boundary conditions.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

This work was supported by National Natural Science Foundation of China (No. 11771337).

Conflict of interest

On behalf of all authors, the corresponding author states that there is no conflict of interest.

References

[1]	S. Ali, On generalized left derivations in rings and Banach algebras, Aequationes Math., 81 (2011), 209–226. https://doi.org/10.1007/s00010-011-0070-5 doi: 10.1007/s00010-011-0070-5
[2]	M. Ashraf, S. Ali, On generalized Jordan left derivations in rings, B. Korean Math. Soc., 452 (2008), 253–261. https://doi.org/10.4134/BKMS.2008.45.2.253 doi: 10.4134/BKMS.2008.45.2.253
[3]	M. Bresar, Jordan mappings of semiprime rings, J. Algebra, 127 (1089), 218–228. http://dx.doi.org/10.1016/0021-8693(89)90285-8 doi: 10.1016/0021-8693(89)90285-8
[4]	J. Cusack, Jordan derivations in rings, P. Am. Math. Soc., 53 (1975), 321–324. https://doi.org/10.1090/S0002-9939-1975-0399182-5 doi: 10.1090/S0002-9939-1975-0399182-5
[5]	B. L. M. Ferreira, W. Feng, Mixed $$ -Jordan-type derivations on $$ -algebras, J. Algebra Appl., 22 (2023), 2350100. https://doi.org/10.1142/S0219498823501001 doi: 10.1142/S0219498823501001
[6]	B. L. M. Ferreira, H. Guzzo, R. N. Ferreira, F. Wei, Jordan derivations of alternative rings, Commun. Algebra, 48 (2020), 717–723. https://doi.org/10.1080/00927872.2019.1659285 doi: 10.1080/00927872.2019.1659285
[7]	I. N. Herstein, Derivations in prime rings, P. Am. Math. Soc., 8 (1957), 1104–1110. https://doi.org/10.1090/S0002-9939-1957-0095864-2 doi: 10.1090/S0002-9939-1957-0095864-2
[8]	B. E. Johnson, A. M. Sinclair, Continuity of derivations and a problem of Kaplansky, Am. J. Math., 90 (1968), 1068–1073. https://doi.org/10.2307/2373290 doi: 10.2307/2373290
[9]	E. C. Posner, Derivations in prime rings, P. Am. Math. Soc., 1957, 1093–1100. https://doi.org/10.1090/S0002-9939-1957-0095863-0
[10]	G. Scudo, A. Z. Ansari, Generalized derivations on Lie ideals and power values on prime rings, Math. Slovaca, 65 (2015), 975–980. https://doi.org/10.1515/ms-2015-0066 doi: 10.1515/ms-2015-0066
[11]	I. M. Singer, J. Wermer, Derivations on commutative normed spaces, Math. Ann., 129 (1995), 435–460.
[12]	M. P. Thomos, The image of a derivation is contained in the radical, Ann. Math., 128 (1988), 435–460. https://doi.org/10.2307/1971432 doi: 10.2307/1971432
[13]	J. Vukman, On left Jordan derivations on rings and Banach algebras, Aequationes Math., 75 (2008), 260–266. https://doi.org/10.1007/s00010-007-2872-z doi: 10.1007/s00010-007-2872-z
[14]	S. M. A. Zaidi, M. Ashraf, S. Ali, On Jordan ideals and left ( $\theta, \theta$ )-derivation in prime rings, Int. J. Math. Math. Sci., 37 (2004), 1957–1965. https://doi.org/10.1155/S0161171204309075 doi: 10.1155/S0161171204309075
[15]	B. Zalar, On centralizers of semiprime rings, Comment. Math. Univ. Ca., 32 (1991), 609–614.
[16]	J. Zhu, C. Xiong, Generalized derivations on rings and mappings of $P$ -preserving kernel into range on Von Neumann algebras, Acta. Math. Sin., 41 (1998), 795–800. https://doi.org/10.12386/A1998sxxb0133 doi: 10.12386/A1998sxxb0133

This article has been cited by:

耕耘田, Temporal Error Estimate of First-Order Fractional Step Algorithm for Unsteady Penetrative Convection Model, 2024, 13, 2324-7991, 3039, 10.12677/aam.2024.137289

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)

通讯作者: 陈斌, bchen63@163.com

1.
沈阳化工大学材料科学与工程学院沈阳 110142

AIMS Mathematics

1.8 3.1

Metrics

Article views(1608) PDF downloads(77) Cited by(2)

Preview PDF

Download XML

Export Citation

Article outline

Show full outline

AIMS Mathematics

An extension of Herstein's theorem on Banach algebra

Related Papers:

Abstract

1. Introduction

2. Materials and methods

3. The Euler grad-div stabilization finite element approximation

4. The BDF2 grad-div stabilization finite element approximation

5. Numerical results

6. Conclusions

Use of AI tools declaration

Acknowledgments

Conflict of interest

References

This article has been cited by:

Reader Comments

通讯作者: 陈斌, bchen63@163.com

Metrics

Other Articles By Authors

Catalog

AIMS Mathematics

An extension of Herstein's theorem on Banach algebra

Related Papers:

Abstract

1. Introduction

2. Materials and methods

3. The Euler grad-div stabilization finite element approximation

4. The BDF2 grad-div stabilization finite element approximation

5. Numerical results

6. Conclusions

Use of AI tools declaration

Acknowledgments

Conflict of interest

References

This article has been cited by:

Reader Comments

通讯作者: 陈斌, bchen63@163.com

Metrics

Other Articles By Authors

Related pages

Tools

Export File

Citation

Format

Content

Catalog