Processing math: 63%
Research article Special Issues

A recursive filter for a class of two-dimensional nonlinear stochastic systems

  • Received: 24 August 2024 Revised: 12 December 2024 Accepted: 30 December 2024 Published: 24 January 2025
  • A recursive filtering problem on minimum variance is investigated for a type of two-dimensional systems incorporating noise and a random parameter matrix in the measurement equation, along with random nonlinearity. It methodically describes random variables using statistical characteristics, placing a strong emphasis on the application of random multivariate analysis and computational techniques. A bidirectional time-sequence recursive filter is designed to achieve unbiasedness and reduce error variance effectively. This involves deriving the gain matrix through a completion of squares method and solving a complex difference equation with two independent variances. To facilitate the online implementation of this filter, various formulations and an algorithm are proposed. A numerical study demonstrates the effectiveness of the design in practical applications.

    Citation: Shulan Kong, Chengbin Wang, Yawen Sun. A recursive filter for a class of two-dimensional nonlinear stochastic systems[J]. AIMS Mathematics, 2025, 10(1): 1741-1756. doi: 10.3934/math.2025079

    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
  • A recursive filtering problem on minimum variance is investigated for a type of two-dimensional systems incorporating noise and a random parameter matrix in the measurement equation, along with random nonlinearity. It methodically describes random variables using statistical characteristics, placing a strong emphasis on the application of random multivariate analysis and computational techniques. A bidirectional time-sequence recursive filter is designed to achieve unbiasedness and reduce error variance effectively. This involves deriving the gain matrix through a completion of squares method and solving a complex difference equation with two independent variances. To facilitate the online implementation of this filter, various formulations and an algorithm are proposed. A numerical study demonstrates the effectiveness of the design in practical applications.



    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:

    utνΔu+(u)u+p(γ1θ+γ2θ2)i3=f,for(x,t) Ω×(0,T], (1.1)
    u=0,for(x,t) Ω×(0,T], (1.2)
    θtκΔθ+(u)θ=g,for(x,t) Ω×(0,T], (1.3)
    u=0,θ=0,for(x,t) Ω×(0,T], (1.4)
    u(x,0)=u0(x),θ(x,0)=θ0(x),forxΩ, (1.5)

    where Ω is a bounded and convex polyhedral domain in R3 with boundary Ω, and (0,T] is a finite time interval. In the above system, the unknown (u,p,θ) denotes the velocity of the fluid, the pressure and temperature, respectively. ν>0 and κ>0 represent the viscosity coefficient and thermal conductivity parameter. f and g are two given functions. i3 is the unit basis vector given by i3=(0,0,1)T, and (γ1θ+γ2θ2)i3 in (1.1) denotes the buoyancy term. The constraint 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., γ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.

    For kN+ and 1p+, let Lp(Ω) and Wk,p(Ω) denote the standard Lebesgue space and Sobolev space, respectively. When p=2, Wk,2(Ω) is the Hilber space Hk(Ω). The norms in Lp(Ω), Wk,p(Ω) and Hk(Ω) are defined as the classical senses (cf.[1]) and are denoted by Lp, Wk,p and Hk. We define Hk0(Ω) to be the subspace of Hk(Ω) of functions with zero trace on Ω. The dual space of H10(Ω) is denoted by H1(Ω). The boldface notations Lp(Ω), Wk,p(Ω) and Hk(Ω) are used to denote the vector-value spaces Lp(Ω)3, Wk,p(Ω)3 and Hk(Ω)3, respectively. The corresponding norms are denoted by Lp, Wk,p and Hk.

    For simplicity, we denote V=H10(Ω), Y=H10(Ω),

    V0={uV,u=0 in Ω},Q=L20(Ω)={q L2(Ω),Ωqdx=0}.

    Introduce the trilinear forms c1(,,) and c2(,,), which are given by

    c1(u1,u2,u3)=Ω((u1)u2u3+12(u1)u2u3)dxc2(u1,θ1,θ2)=Ω((u1θ1)θ2+12(u1)θ1θ2)dx

    for uiV with i=1,2,3 and θjY with j=1,2. By the integration by parts, it is easy to check that the trilinear forms c1(,,) and c2(,,) satisfy the skew-symmetric properties, i.e.,

    c1(u1,u2,u3)=c1(u1,u3,u2)andc2(u1,θ1,θ2)=c2(u1,θ2,θ1),

    which give

    c1(u1,u2,u2)=0andc2(u1,θ1,θ1)=0. (2.1)

    For given fL2(Ω),gL2(Ω), u0L2(Ω) and θ0L2(Ω), in terms of ther skew-symmetric properties (2.1), there {holds} the following energy inequalities:

    θ(t)2L2+κt0θ(s)2L2dsθ02L2+Ct0g(s)2L2ds, (2.2)

    and

    u(t)2L2+νt0u(s)2L2dsu02L2+Ct0(f(s)2L2+θ(s)2L2+θ(s)2L2θ(s)2L3)dsu02L2+Ct0(f(s)2L2+g(s)2L2)ds+C(t0g(s)2L2ds)2 (2.3)

    for 0tT.

    Let 0=t0<t1<<tN=T be a uniform partition of the time interval [0,T] with time step τ=T/N and tn=nτ with 0nN. Denote un=u(tn), pn=p(tn), θn=θ(tn), fn=f(tn) and gn=g(tn). For any sequence of functions {fn}Nn=0, denote the backward Euler discretization and BDF2 discretization by

    D1,τfn+1=(fn+1fn)/τfor  0nN1,

    and

    D2,τfn+1=(3fn+14fn+fn1)/2τ,ˆfn=2fnfn1for  1nN1,

    where ˆfn is the standard extrapolation formula. For the discrete time derivative D2,τ, there holds the telescope formula [35]:

    (D2,τfn+1,fn+1)=14τ(fn+12L2fn2L2+ˆfn+12L2ˆfn2L2)+14τfn+12fn+fn12L2. (2.4)

    We define the finite element spaces as follows. Let Th={Kj}Lj=1 be a quasi-uniform tetrahedron partition of Ω with the mesh size h=max{diamKj, j=1,,L}. In finite element discretization of (1.1)–(1.3), we use the Taylor-Hood (P2P1) finite element for the velocity u and pressure p, and the piece-wise quadratic finite element for the temperature θ. The corresponding finite element subspaces of V,Q and Y are denoted by Vh, Qh and Yh, respectively, i.e.,

    Vh={vhC(¯Ω)V, vh|K(P2(K))3,  KTh},Qh={qhC(¯Ω)H1(Ω), qh|KP1(K),  KTh, Ωqhdx=0},Yh={rhC(¯Ω)Y, rh|KP2(K),  KTh}.

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

    β0qhL2supvhVh(vh,qh)vhL2, qhQh, (2.5)

    where β0>0 is independent of the mesh size h.

    Denote the space of discrete divergence-free functions by

    V0h={vhVh,(vh,qh)=0 qhQh}.

    In the sequel, let πh and Rh be the L2 orthogonal and Ritz orthogonal projections onto Qh and Yh, respectively. Then the following error bounds and stability hold [4,43]:

    pπhpL2Ch2pH2 pH2(Ω), (2.6)
    ψRhψL2+h(ψRhψ)L2Ch3ψH3 ψH3(Ω), (2.7)
    RhψW1,CψH3 ψH3(Ω). (2.8)

    Let Ih be the Lagrange interpolation operator onto Vh. The following bound can be found in [4]:

    uIhuWm,pChnmuWn,p0mn3, (2.9)

    where n>3/p when 1<p and n3 when p=1.

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

    νΔu+p=Fin Ω,u=0in Ω,u=0on Ω. (2.10)

    Let us denote (uh,ph)Vh×Qh, the mixed finite element approximation solution to the Stokes problem (2.10). Then one has error estimates [22,30]

    (uuh)L2C(infvhVhuvhH1+ν1infqhQhpqhL2), (2.11)
    pphL2C(νinfvhVhuvhH1+infqhQhpqhL2), (2.12)
    uuhL2Ch(infvhVhuvhH1+ν1infqhQhpqhL2), (2.13)

    where C>0 is independent of h,ν and κ. It is clear that the error estimates of the velocity depend on the negative power of ν. To prove error bounds with constants being independent of the negative power of ν, we recall the specific projection of (u,p) introduced in [17] that we will denote by sh defined by

    (sh,vh)=(u,vh),vhV0h.

    Let (u,p,θ) be the solution to the system (1.1)–(1.5) with uL((0,T];VH3(Ω)), pL((0,T];QH2(Ω)), θL((0,T];YH3(Ω)) and utL((0,T];H1(Ω)). We note that (u,0) is the solution to the Stokes problem (2.10) with the right-hand side F given by

    F=fut(u)up+(γ1θ+γ2θ2)i3. (2.14)

    Let (sh,lh)Vh×Qh be the corresponding mixed finite element approximation of (u,0). Then, from (2.11)–(2.13), we have

    ushL2+hushH1Ch3uH3, (2.15)
    lhL2Cνh2uH3, (2.16)

    where C>0 is independent of h and ν. In terms of the inverse inequality in the theory of the finite element method and IhuLCuL, we have

    shLshIhuL+IhuLC(h3/2shIhuL2+uL)C(h3/2shuL2+h3/2uIhuL2+uH3)CuH3, (2.17)

    where C>0 is independent of h,ν and κ. Furthermore, there holds

    (ush)LCuLCuH3andshLCuLCuH3, (2.18)

    where C>0 is independent of h,ν and κ.

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

    Lemma 2.1. Let ak,bk,ck and γk, for integers k0, be the nonnegative numbers such that

    an+τnk=0bkτnk=0γkak+τnk=0ck+Bfor n0. (2.19)

    Suppose that τγk<1 for all k and set σk=(1τγk)1, then

    an+τnk=0bkexp(τnk=0γkσk)(τnk=0ck+B)for n0. (2.20)

    Remark 2.1. If the first sum on the right in (2.19) extends only up to n1, then the estimate (2.20) holds for all τ>0 with σk=1.

    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 0nN1.

    Step Ⅰ: Find θn+1hYh by

    (D1,τθn+1h,ψh)+κ(θn+1h,ψh)+c2(unh,θn+1h,ψh)=(gn+1,ψh) (3.1)

    with the initial iteration values θ0h=Rhθ0, u0h=Ihu0 and for any ψhYh.

    Step Ⅱ: Find (un+1h,pn+1h)Vh×Qh by

    (D1,τun+1h,vh)+ν(un+1h,vh)+c1(unh,un+1h,vh)(vh,pn+1h)+(un+1h,qh)+β(un+1h,vh)γ1(θnhi3,vh)γ2((θnh)2i3,vh)=(fn+1,vh) (3.2)

    for any (vh,qh)Vh×Qh, where β>0 is the stabilization parameter. The detailed discussion about the choice of β 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 (un+1h,pn+1h,θn+1h) to the Euler scheme (3.1)–(3.2).

    Lemma 3.1. For 0nN1 and all τ>0, h>0, the finite element scheme (3.1) and (3.2) has a unique solution θn+1hYh and (un+1h,pn+1h)Vh×Qh. Moreover, the discrete energy inequalities hold:

    θm+1h2L2+κτmn=0θn+1h2L2CτN1n=0gn+12L2+θ0h2L2, (3.3)

    and

    um+1h2L2+τνmn=0um+1h2L2+2τβmn=0um+1h2L2CτN1n=0(fn+12L2+gn+12L2)+C(u0h2L2+θ0h2L2)+C(τN1n=0gn+12L2+θ0h2L2)2 (3.4)

    for all 0mN1, where C>0 is independent of τ and h.

    Proof. Setting ψh=2τθn+1h in (3.1), using the skew-symmetric property (2.1), the Hölder inequality and Young inequality, we have

    θn+1h2L2θnh2L2+θn+1hθnh2L2+2τκθn+1h2L22τgn+1L2θn+1hL2Cτgn+12L2+τκθn+1h2L2.

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

    Similarly, taking (vh,qh)=2τ(un+1h,pn+1h) in (3.2), we have

    un+1h2L2unh2L2+un+1hunh2L2+2τνun+1h2L2+2τβun+1h2L2Cτfn+1L2un+1hL2+CτθnhL2un+1hL2+CτθnhL2θnhL3un+1hL6τνun+1h2L2+Cτ(fn+12L2+θnh2L2+θnh2L2θnh2L3).

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

    Denote

    enh=snhunh,θnθnh=(θnRhθn)+(Rhθnθnh):=ηnθ+enθ.

    The sh(t) is the finite element solution to the Stokes problem (2.10) with the right-hand side (2.14) and snh=sh(tn).

    Thanks to the discrete inf-sup condition (2.5), the finite element scheme (3.2) is equivalent to: Find un+1hV0h by

    (D1,τun+1h,vh)+ν(un+1h,vh)+c1(unh,un+1h,vh)+β(un+1h,vh)γ1(θnhi3,vh)γ2((θnh)2i3,vh)=(fn+1,vh) vhV0h. (3.5)

    For simplicity, we take s0h=u0h. For 0nN1, sn+1h and θn+1 satisfy

    (D1,τsn+1h,vh)+ν(sn+1h,vh)+c1(snh,sn+1h,vh)+β(sn+1h,vh)γ1(θni3,vh)γ2((θn)2i3,vh)=(fn+1,vh)+(vh,pn+1πhpn+1)+β((sn+1hun+1),vh)+(D1,τsn+1hut(tn+1),vh)+c1(snh,sn+1h,vh)c1(un+1,un+1,vh)+γ1(θn+1θn,i3vh)+γ2((θn+1)2(θn)2,i3vh) vhV0h, (3.6)

    and

    (D1,τθn+1,ψh)+κ(θn+1,ψh)+c2(un+1,θn+1,ψh)=(D1,τθn+1θt(tn+1),ψh)+(gn+1,ψh) ψhYh. (3.7)

    In (3.6), we use un+1=0 in Ω and (vh,πhpn+1)=0 due to vhV0h.

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

    (D1,τen+1h,vh)+ν(en+1h,vh)+β(en+1h,vh)=6j=1(Jn+1j,vh) vhV0h, (3.8)

    where

    (Jn+11,vh)=(D1,τsn+1hut(tn+1),vh),(Jn+12,vh)=(vh,(pn+1πhpn+1)+β(sn+1hun+1)),(Jn+13,vh)=c1(snh,sn+1h,vh)c1(un+1,un+1,vh),(Jn+14,vh)=c1(unh,un+1h,vh)c1(snh,sn+1h,vh),(Jn+15,vh)=γ1(θn+1θn,i3vh)+γ2((θn+1)2(θn)2,i3vh),(Jn+16,vh)=γ1((θnθnh),i3vh)+γ2(((θn)2(θnh)2),i3vh).

    Subtracting (3.1) from (3.7), we have

    (D1,τen+1θ,ψh)+κ(en+1θ,ψh)+c2(unh,en+1θ,ψh)=3j=1(λn+1j,ψh) ψhYh, (3.9)

    where

    (λn+11,ψh)=(D1,τθn+1θt(tn+1),ψh),(λn+12,ψh)=(D1,τηn+1θ,ψh),(λn+13,ψh)=c2(unh,Rhθn+1,ψh)c2(un+1,θn+1,ψh).

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

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

    {u0H3(Ω), θ0H3(Ω), fL((0,T];L2(Ω)),gL((0,T];L2(Ω)),uL((0,T];H3(Ω)),pL((0,T];H2(Ω)),θL((0,T];H3(Ω)),utL((0,T];H1(Ω))L2((0,T];H2(Ω)),θtL((0,T];H1(Ω))L2((0,T];H2(Ω)),uttL2((0,T];L2(Ω)),θttL2((0,T];L2(Ω)). (3.10)

    Then, for the first-order Euler grad-div scheme (3.1)–(3.2) under the time step condition τCh2, when h and τ are sufficiently small, we have the following error estimate:

    em+1h2L2+em+1θ2L2+τmn=0(νem+1h2L2+κem+1θ2L2)C20(h4+τ2) (3.11)

    with 0mN1, where C0>0 is independent of h,τ, ν and κ.

    Proof. Taking vh=2τen+1h in (3.8) and summing up the resulting equation from n=0 to n=m, we have

    em+1h2L2+2τmn=0(νen+1h2L2+βen+1h2L2)=2τmn=06j=1(Jn+1j,en+1h). (3.12)

    We estimate the right-hand side of (3.12) as follows. For Jn+11, we split it as

    Jn+11=(D1,τsn+1hD1,τun+1)+(D1,τun+1ut(tn+1)).

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

    τmn=0D1,τun+1ut(tn+1)2L2Cτ2.

    Notice that

    D1,τsn+1hD1,τun+1=1τtn+1tnt(sh(t)u(t))dt.

    From the Hölder inequality, we have

    D1,τsn+1hD1,τun+1L21τtn+1tnt(sh(t)u(t))L2dtCh2τtn+1tntu(t)H2dtCh2τ1/2(tn+1tntu(t)2H2dt)1/2.

    Thus,

    2τmn=0(Jn+11,en+1h)Cτmn=0en+1h2L2+C(τ2+h4), (3.13)

    where C>0 is independent of h,τ,ν and κ.

    For Jn+12, it follows from the Hölder inequality and (2.6) that

    (Jn+12,en+1h)(pn+1πhpn+1L2+βsn+1hun+1H1)en+1hL2Ch2en+1hL2,

    then

    2τmn=0(Jn+12,en+1h)Ch4+ϵ1βτmn=0en+1h2L2, (3.14)

    where C>0 is independent of h,τ,ν,κ and ϵ1>0 is some small constant determined later.

    For Jn+13, in terms of the skew-symmetric property (2.1), (2.18), the Hölder inequality and the regularity assumption (3.10), there holds

    (Jn+13,en+1h)=c1(snhun,sn+1h,en+1h)+c1(unun+1,sn+1h,en+1h)+c1(un+1,sn+1hun+1,en+1h)(snhunH1+un+1unH1)sn+1hW1,en+1hL2+un+1Lsn+1hun+1H1en+1hL2C(h2+τ)en+1hL2,

    where C>0 is independent of h,τ,ν and κ. Thus, the third term can be estimated by

    2τmn=0(Jn+13,en+1h)C(h4+τ2)+Cτmn=0en+1h2L2, (3.15)

    where C>0 is independent of h,τ,ν and κ.

    A similar method to bound Jn+14 leads to

    (Jn+14,en+1h)=c1(enh,sn+1h,en+1h)enhL2sn+1hLen+1hL2+sn+1hLenhL2en+1hL2,

    and

    2τmn=0(Jn+14,en+1h)ϵ1βτmn=0en+1h2L2+Cτmn=0(en+1h2L2+enh2L2), (3.16)

    where C>0 is independent of h,τ,ν and κ.

    For Jn+15, using the Hölder inequality, Young inequality and the regularity assumption (3.10), we have

    |(Jn+15,en+1h)|C(1+θn+θn+1L)tn+1tnθt(t)L2dten+1hL2Cen+1h2L2+Cτ2,

    then

    2τmn=0(Jn+15,en+1h)Cτmn=0en+1h2L2+Cτ2, (3.17)

    where C>0 is independent of h,τ,ν and κ.

    By a similar method, we estimate the last term by

    (Jn+16,en+1h)ηnθ+enθL2en+1hL2+(ηnθ+enθ)L2θn+θnhLen+1hL2C(1+θnh2L)(h4+enθ2L2)+C(enθ2L2+en+1h2L2),

    where we use the Sobolev imbedding inequality and C>0 is independent of h,τ,ν and κ. Taking the sum gives

    2τmn=0(Jn+16,en+1h)Ch4τmn=0(1+θnh2L)+Cτmn=0(enθ2L2+en+1h2L2)+Cτmn=0(1+θnh2L)enθ2L2, (3.18)

    where C>0 is independent of h,τ,ν and κ.

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

    em+1h2L2+2τmn=0(νen+1h2L2+βen+1h2L2)C(h4+τ2)+Ch4τmn=0(1+θnh2L)+Cτmn=0(enθ2L2+en+1h2L2+enh2L2)+2ϵ1βτmn=0en+1hL2+Cτmn=0(1+θnh2L)enθ2L2, (3.19)

    where C>0 is independent of h,τ,ν and κ.

    Next, we estimate en+1θ. Taking ψh=2τen+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

    em+1θ2L2+2κτmn=0en+1θ2L2=2τmn=03j=1(λn+1j,en+1θ), (3.20)

    where we noted e0θ=0. By the Taylor's formula and the Hölder inequality, it is easy to prove that

    2τmn=0(λn+11,en+1θ)Cτmn=0en+1θ2L2+Cτ2, (3.21)
    2τmn=0(λn+12,en+1θ)Cτmn=0en+1θ2L2+Ch4τmn=0D1,τθn+12H1Cτmn=0en+1θ2L2+Ch4 (3.22)

    by noticing

    τmn=0D1,τθn+12H11τmn=0(tn+1tnθt(t)H1dt)2C,

    where C>0 is independent of h,τ, ν and κ.

    For last term λn+13, in terms of the following splitting:

    (λn+13,en+1θ)=(unhRhθn+1,en+1θ)+12((unh)Rhθn+1,en+1θ)(un+1θn+1,en+1θ)=(enhRhθn+1,en+1θ)(snhηn+1θ,en+1θ)+((snhun+1)θn+1,en+1θ)12((enh)Rhθn+1,en+1θ)+12((snhun)Rhθn+1,en+1θ),

    we estimate it by

    (λn+13,en+1θ)(enhL2Rhθn+1L+snhLηn+1θL2+snhun+1L3θn+1L6)en+1θL2+12(enhL2Rhθn+1L+snhunH1Rhθn+1L)en+1θL2C(h4+τ2)+C(en+1θ2L2+enh2L2)+ϵ1β2enh2L2.

    Furthermore, we get

    2τmn=0(λn+13,en+1θ)C(h4+τ2)+Cτmn=0(en+1θ2L2+enh2L2)+ϵ1βτmn=0enh2L2, (3.23)

    where C>0 is independent of h,τ, ν and κ.

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

    em+1θ2L2+2κτmn=0en+1θ2L2C(h4+τ2)+Cτmn=0(en+1θ2L2+enh2L2)+ϵ1βτmn=0enh2L2, (3.24)

    where C>0 is independent of h,τ, ν and κ.

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

    em+1h2L2+em+1θ2L2+2τmn=0(νen+1h2L2+κen+1θ2L2+βen+1h2L2)Cτmn=0(en+1θ2L2+en+1h2L2+enθ2L2+enh2L2)+C(h4+τ2)+Ch4τmn=0(1+θnh2L)+3ϵ1βτmn=0en+1hL2+Cτmn=0(1+θnh2L)enθ2L2, (3.25)

    where C>0 is independent of h,τ, ν and κ.

    To uniformly bound the norm θnhL, we use the method of mathematical induction. Taking m=0 in (3.25) and noting e0h=0, e0θ=0 and θ0hL=Rhθ0Lθ0H3, we get

    e1h2L2+e1θ2L2+2τ(νe1h2L2+κe1θ2L2+βe1h2L2)Cτ(e1θ2L2+e1h2L2)+C(h4+τ2)+3ϵ1βτe1hL2. (3.26)

    For sufficiently small τ with Cτ<1/2, we select small parameter ϵ12/3 in (3.26). We then conclude that there exists some C1>0, which is independent of h,τ, ν and κ such that

    e1h2L2+e1θ2L2+τ(νe1h2L2+κe1θ2L2)C21(h4+τ2).

    Thus, the error estimate (3.11) holds for m=0 by taking C0>C1. Now, we assume that the error estimate (3.11) holds for mk1 with 1kN1. Under the time step condition τCh2, we have

    enθL2CC0h2 1nk. (3.27)

    By the inverse inequality, we get

    enθLCh3/2enθL2CC0h1/2C

    for sufficiently small h with C0h1/21. Thus,

    θnhLenθL+θnLCandτkn=1θnh2LC. (3.28)

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

    ek+1h2L2+ek+1θ2L2+2τkn=0(νen+1h2L2+κen+1θ2L2+βen+1h2L2)C(h4+τ2)+Cτkn=0(en+1θ2L2+en+1h2L2+enθ2L2+enh2L2)+3ϵ1βτkn=0en+1hL2. (3.29)

    We select small parameter ϵ12/3 in (3.29), then

    ek+1h2L2+ek+1θ2L2+τkn=0(νen+1h2L2+κen+1θ2L2)C(h4+τ2)+Cτkn=0(en+1θ2L2+en+1h2L2+enθ2L2+enh2L2).

    Applying the discrete Gronwall inequality to the above inequality, we conclude that there exists some C2>0, which is independent of h,τ, ν and κ such that

    ek+1h2L2+ek+1θ2L2+τkn=0(νen+1h2L2+κen+1θ2L2)C22(h4+τ2)

    for sufficiently small τ. Thus, we prove that the error estimate (3.11) also holds for m=k by taking C0C2, 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

    um+1um+1h2L2+θm+1θm+1h2L2C(h4+τ2) (3.30)

    with 0mN1, where C>0 is independent of h,τ, ν and κ.

    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 1nN1:

    Step Ⅰ: Find θn+1hYh by

    (D2,τθn+1h,ψh)+κ(θn+1h,ψh)+c2(ˆunh,θn+1h,ψh)=(gn+1,ψh) ψhYh. (4.1)

    Step Ⅱ: Find (un+1h,θn+1h)Vh×Qh by

    (D2,τun+1h,vh)+ν(un+1h,vh)+c1(ˆunh,un+1h,vh)(vh,pn+1h)+(un+1h,qh)+β(un+1h,vh)(γ1ˆθnh+γ2(ˆθnh)2,i3vh)=(fn+1,vh) (vh,qh)Vh×Qh. (4.2)

    Remark 4.1. Numerical solution (u1h,θ1h)V0h×Yh 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

    e1h2L2+e1θ2L2+τ(νe1h2L2+κe1θ2L2+βe1h2L2)C23(h4+τ4), (4.3)

    where C3>0 is independent of h,τ, ν and κ.

    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 (un+1h,pn+1h,θn+1h) to the BDF2 scheme (4.1)–(4.2).

    Lemma 4.1. For 1nN1 and all τ>0, h>0, the BDF2 finite element scheme (4.1) and (4.2) has a unique solution θn+1hYh and (un+1h,pn+1h)Vh×Qh. Moreover, the discrete energy inequalities hold:

    θm+1h2L2+ˆθm+1h2L2+2κτmn=1θn+1h2L2Cτmn=1gn+12L2+C(θ1h2L2+θ0h2L2), (4.4)

    and

    um+1h2L2+ˆum+1h2L2+2τνmn=1un+1h2L2+4τβmn=1un+1h2L2Cτmn=1(fn+12L2+gn+12L2)+C(θ1h2L2+θ0h2L2)+C(τmn=1gn+12L2+θ1h2L2+θ0h2L2)2 (4.5)

    for all 1mN1, where C>0 is independent of h and τ.

    In terms of the discrete inf-sup condition (2.5), we rewrite the discrete variational problem (4.2) as: find un+1hV0h by

    (D2,τun+1h,vh)+ν(un+1h,vh)+c1(ˆunh,un+1h,vh)+β(un+1h,vh)(γ1ˆθnh+γ2(ˆθnh)2,i3vh)=(fn+1,vh) vhV0h. (4.6)

    On the other hand, let s0h=u0h. For 1nN1, sn+1h and θn+1 satisfy

    (D2,τsn+1h,vh)+ν(sn+1h,vh)+c1(ˆsnh,sn+1h,vh)+β(sn+1h,vh)γ1(ˆθn,i3vh)γ2((ˆθn)2,i3vh)=(fn+1,vh)+(vh,pn+1πhpn+1)+β((sn+1hun+1),vh)+(D2,τsn+1hut(tn+1),vh)+c1(ˆsnh,sn+1h,vh)c1(un+1,un+1,vh)+γ1(θn+1ˆθn,i3vh)+γ2((θn+1)2(ˆθn)2,i3vh) vhV0h, (4.7)

    and

    (D2,τθn+1,ψh)+κ(θn+1,ψh)+c2(un+1,θn+1,ψh)=(gn+1,ψh)+(D2,τθn+1θn+1t,ψh) ψhYh. (4.8)

    In (4.7), we use un+1=0 in Ω and (vh,πhpn+1)=0.

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

    (D2,τen+1h,vh)+ν(en+1h,vh)+β(en+1h,vh)=6j=1(In+1j,vh) vhV0h, (4.9)

    where

    (In+11,vh)=(D2,τsn+1hut(tn+1),vh),(In+12,vh)=(vh,(pn+1πhpn+1)+β(sn+1hun+1)),(In+13,vh)=c1(ˆsnh,sn+1h,vh)c1(un+1,un+1,vh),(In+14,vh)=c1(ˆunh,un+1h,vh)c1(ˆsnh,sn+1h,vh),(In+15,vh)=γ1(θn+1ˆθn,i3vh)+γ2((θn+1)2(ˆθn)2,i3vh),(In+16,vh)=γ1((ˆθnˆθnh),i3vh)+γ2(((ˆθn)2(ˆθnh)2),i3vh).

    Subtracting (4.1) from (4.8) leads to

    (D2,τen+1θ,ψh)+κ(en+1θ,ψh)+c2(ˆunh,en+1θ,ψ)=3j=1(Xn+1j,ψh) ψhYh, (4.10)

    where

    (Xn+11,ψh)=(D2,τθn+1θt(tn+1),ψh),(Xn+12,ψh)=(D2,τηn+1θ,ψh),(Xn+13,ψh)=c2(ˆunh,Rhθn+1,ψh)c2(un+1,θn+1,ψh).

    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

    utttL2((0,T];L2(Ω)),θtttL2((0,T];L2(Ω)). (4.11)

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

    em+1h2L2+em+1θ2L2+τmn=0(νem+1h2L2+κem+1θ2L2)ˆC20(h4+τ4) (4.12)

    with 0mN1, where ˆC0>0 is independent of h,τ, ν and κ.

    Proof. Taking vh=4τen+1h in (4.9) and summing up the resulting equation from n=1 to n=m, we have

    em+1h2L2+ˆem+1h2L2+4τmn=1(νen+1h2L2+βen+1h2L2)4τmn=16j=1(In+1j,en+1h)+Ce1hL2, (4.13)

    where C>0 is independent of h,τ, ν,κ and we use

    ˆe1hL22e1hL2+e0hL2=2e1hL2.

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

    In+11=(D2,τsn+1hD2,τun+1)+(D2,τun+1ut(tn+1)),

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

    D2,τun+1ut(tn+1)=12τtn+1tn1(2(ttn)212(ttn1)2)tttu(t)dt,

    which results in

    τmn=1D2,τun+1ut(tn+1)2L2Cτ4,

    where C>0 is independent of h,τ, ν and κ. From the Hölder inequality, we have

    D2,τsn+1hD2,τun+1L232τtn+1tnt(sh(t)u(t))L2dt+12τtntn1t(sh(t)u(t))L2dtCh2τtn+1tn1tu(t)H2dtCh2τ1/2(tn+1tn1tu(t)2H2dt)1/2.

    Thus, we can obtain

    4τmn=1(In+11,en+1h)Cτmn=1en+1h2L2+C(τ4+h4), (4.14)

    where C>0 is independent of h,τ, ν and κ.

    From the Hölder inequality and (2.6), we estimate (In+12,en+1h) by

    (In+12,en+1h)(pn+1πhpn+1L2+βsn+1hun+1H1)en+1hL2Ch2en+1hL2.

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

    4τmn=1(In+12,en+1h)Ch4+ϵ2βτmn=1en+1h2L2, (4.15)

    where C>0 is independent of h,τ,ν,κ and ϵ2>0 is a small constant determined later.

    For the term of In+13, according to the Hölder inequality, (2.6), (2.18), the regularity assumption (3.10) and the Taylor formula with integral remainder

    ˆunun+1=tn+1tn1{2(ttn)+(ttn1)}ttu(t)dt, (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 .

    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 Tables 3 and 4 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

    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.

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

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

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



    [1] A. Habibi, Two-dimensional Bayesian estimate of images, Proc. IEEE, 60 (1972) 877–883. https://doi.org/10.1109/PROC.1972.8787
    [2] R. P. Roesser, A discrete state-space model for linear image processing, IEEE Trans. Autom. Contr., 20 (1975), 1–10. https://doi.org/10.1109/TAC.1975.1100844 doi: 10.1109/TAC.1975.1100844
    [3] E. Fornasini, G. Marchesini, State-space realization theory of two-dimensional filters, IEEE Trans. Autom. Control, 21 (1976), 484–492. https://doi.org/10.1109/TAC.1976.1101305 doi: 10.1109/TAC.1976.1101305
    [4] J. W. Woods, C. H. Radewan, Kalman filtering in two dimensions, IEEE Trans. Inf. Theory, 23 (1977), 473–482. https://doi.org/10.1109/TIT.1977.1055750 doi: 10.1109/TIT.1977.1055750
    [5] J. W. Woods, V. K. Ingle, Kalman filtering in two dimensions: Further results, IEEE Trans. Acoust., Speech, Signal Process., 29 (1981), 188–197. https://doi.org/10.1109/TASSP.1981.1163533 doi: 10.1109/TASSP.1981.1163533
    [6] T. Katayama, M. Kosana, Recursive filtering algorithm for a two-dimensional system, IEEE Trans. Autom. Control, 24 (1979), 130–132. https://doi.org/10.1109/TAC.1979.1101956 doi: 10.1109/TAC.1979.1101956
    [7] M. \breve{S}ebek, Polynomial solution of 2D Kalman Bucy filtering problem, IEEE Trans. Autom. Control, 37 (1992), 1530–1533. https://doi.org/10.1109/9.256417 doi: 10.1109/9.256417
    [8] A. Concetti, L. Jetto, Two-dimensional recursive filtering algorithm with edge preserving properties and reduced numerical complexity, IEEE Trans. Circuits Syst. II-Analog Digital Signal Process., 44 (1997), 587–591. https://doi.org/10.1109/82.598429 doi: 10.1109/82.598429
    [9] X. Chen, C. Yang, The state estimation of the stochastic 2D FMII models, Acta Autom. Sin., 27 (2001), 131–135.
    [10] Y. Zou, M. Sheng, N. Zhong, S. Xu, A generalized Kalman filter for 2D discrete systems, Circuits Syst. Signal Process., 23 (2004), 351–364. https://doi.org/10.1007/s00034-004-0804-x doi: 10.1007/s00034-004-0804-x
    [11] J. Liang, F. Wang, Z. Wang, Minimum-variance recursive filtering for two-dimensional systems with degraded measurements: Boundedness and Monotonicity, IEEE Trans. Autom. Control, 64 (2019), 4153–4166. https://doi.org/10.1109/TAC.2019.2895245 doi: 10.1109/TAC.2019.2895245
    [12] F. Wang, Z. Wang, J. Liang, X. Liu, Robust finite horizon filtering for 2-D systems with randomly varying sensor delays, IEEE Trans. Syst., Man, Cybern., Syst., 50 (2020), 220–232. https://doi.org/10.1109/TSMC.2017.2788503 doi: 10.1109/TSMC.2017.2788503
    [13] F. Wang, Z. Wang, J. Liang, C. Silvestre, Recursive locally minimum-variance filtering for two-dimensional systems: When dynamic quantization effect meets random sensor failure, Automatica, 148 (2023), 110762. https://doi.org/10.1016/j.automatica.2022.110762 doi: 10.1016/j.automatica.2022.110762
    [14] F. Wang, J. Liang, J. Lam, J. Yang, C. Zhao, Robust filtering for 2-D systems with uncertain-variance noises and weighted try-once-discard protocols, IEEE Trans. Syst., Man, Cybern., Syst., 53 (2023), 2914–2924. https://doi.org/10.1109/TSMC.2022.3219919 doi: 10.1109/TSMC.2022.3219919
    [15] F. Wang, Z. Wang, J. Liang, Q. Ge, S. X. Ding, Recursive filtering for two-dimensional systems with amplify-and-forward relays: Handling degraded measurements and dynamic biases, Inform. Fusion, 108 (2024), 102368. https://doi.org/10.1016/j.inffus.2024.102368 doi: 10.1016/j.inffus.2024.102368
    [16] P. Zhang, C. Zhu, B. Yang, Z. Wang, M. Hao, Event-triggered ultimately bounded filtering for two-dimensional discrete-time systems under hybrid cyber attacks, J. Franklin Inst., 361 (2024), 683–711. https://doi.org/10.1016/j.jfranklin.2023.12.019 doi: 10.1016/j.jfranklin.2023.12.019
    [17] S. Kong, Y. Sun, H. Zhang, Optimal Kalman-like filtering for a class of nonlinear stochastic systems, J. Ocean Eng. Sci., 8 (2023), 500–507. https://doi.org/10.1016/j.joes.2022.03.002 doi: 10.1016/j.joes.2022.03.002
    [18] W. D. Koning, Optimal estimation of linear discrete-time systems with stochastic parameters, Automatica, 20 (1984), 113–115. https://doi.org/10.1016/0005-1098(84)90071-2 doi: 10.1016/0005-1098(84)90071-2
    [19] J. Hu, Z. Wang, H. Gao, Recursive filtering with random parameter matrices, multiple fading measurements and correlated noise, Automatica, 49 (2013), 3440–3448. https://doi.org/10.1016/j.automatica.2013.08.021 doi: 10.1016/j.automatica.2013.08.021
    [20] W. Wang, J. Zhou, Optimal linear filtering design for discrete-time systems with cross-correlated stochastic parameter matrices and noises, IET Control Theory Appl., 11 (2017), 3353–3362. https://doi.org/10.1049/iet-cta.2017.0425 doi: 10.1049/iet-cta.2017.0425
  • This article has been cited by:

    1. 耕耘 田, 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
  • © 2025 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

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(590) PDF downloads(47) Cited by(0)

Figures and Tables

Figures(2)  /  Tables(1)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog