AN EMBEDDED DISCONTINUOUS GALERKIN METHOD FOR THE OSEEN EQUATIONS

In this paper, the a priori error estimates of an embedded discontinuous Galerkin method for the Oseen equations are presented. It is proved that the velocity error in the L(Ω) norm, has an optimal error bound with convergence order k +1, where the constants are dependent on the Reynolds number (or ν−1), in the diffusion-dominated regime, and in the convection-dominated regime, it has a Reynolds-robust error bound with quasi-optimal convergence order k + 1/2. Here, k is the polynomial order of the velocity space. In addition, we also prove an optimal error estimate for the pressure. Finally, we carry out some numerical experiments to corroborate our analytical results. Mathematics Subject Classification. 65N12, 65N22, 65N30, 76D07. Received April 21, 2021. Accepted September 20, 2021.

popular in simulating the Navier-Stokes equations [10,14,23,30]. Unfortunately, DG methods have a large number of degrees of freedom compared to the continuous finite element methods. To address this problem, the hybridized discontinuous Galerkin (HDG) method has been developed by introducing new degrees of freedom on the element interfaces. The globally coupled unknowns after static condensation are defined only on element interfaces rather than the interior of the cells, thereby obviating the common criticism of DG methods.
In recent years, many H(div)-conforming HDG methods have been developed [11,21,22,26]. The HDG method introduced in [26] uses discontinuous facet function spaces for the trace velocity and pressure approximations. A simple modification is possible to use continuous facet function spaces for the trace velocity, which was wellknown as an embedded-hybridized discontinuous Galerkin (E-HDG) method [28]. The E-HDG method has less globally coupled degrees of freedom than the HDG method due to the use of continuous facet function spaces on a given mesh. The HDG and E-HDG methods provide an exactly divergence-free and H(div)-conforming velocity field, in which the velocity error bound is pressure-robust [28]. The number of global coupled degrees of freedom can be reduced even further by using continuous facet function spaces for the trace velocity and pressure approximations, which was well-known as an embedded discontinuous Galerkin (EDG) method. Using continuous facet function spaces, the number of global degrees of freedom is the same as that of a continuous Galerkin method on a given mesh. These properties, usually associated with DG methods, can be obtained with the same number of global degrees of freedom as a continuous Galerkin method on the same mesh, which has been presented for the incompressible Navier-Stokes equations [19]. The EDG method combines attractive properties of continuous and discontinuous Galerkin methods. The EDG method for the Stokes problem results in a pointwise divergence-free approximate velocity on cells. Unfortunately, the EDG method is not pressurerobust, because the approximate velocity is not H(div)-conforming [28]. In [28], it is proved that using continuous facet function spaces can significantly reduce the CPU time and the number of iterations to convergence than using discontinuous facet function spaces, by using the preconditioner in [27]. In view of the computational efficiency of the EDG method, a pressure-robust EDG method have been developed for the Stokes equations by introducing a local reconstruction operator [20]. Finally, we can see that these types of methods mentioned above has been applied to the incompressible Navier-Stokes equations [19,26]. In recent years, the space-time formulations of these types of methods have been developed for the incompressible Navier-Stokes equations on moving domains [16,17].
In this paper, we consider the embedded discontinuous Galerkin method for the Oseen equations. It is proved that a Reynolds-dependent error bound for the 2 error of the velocity, has an optimal convergence order + 1 in the diffusion-dominated regime. In the convection-dominated regime, it has a Reynolds-robust error bound with quasi-optimal convergence order + 1/2 for the 2 error of the velocity. The analysis here also covers the pressure-robust HDG and E-HDG methods. This work was motivated by the EDG method, which combines attractive properties of continuous and discontinuous Galerkin methods.
The structure of the paper is as follows: In Section 2, we introduce the EDG, HDG and E-HDG methods for the Oseen equations. Some preliminaries are presented in Section 3. In Section 4, the error estimates for the velocity are given and for the pressure in Section 5. Numerical studies are presented in Section 6. We end in Section 7 with some conclusions and future directions.

Notation
Let { ℎ } 0<ℎ≤1 be a family of triangulations of the domain Ω, which are shape-regular and quasi-uniform. For each triangulation ℎ , define mesh size ℎ = ∈ ℎ ℎ , with ℎ the diameter of each element ∈ ℎ . Denote the set of all facets and the mesh skeleton by ℱ ℎ and 0 , respectively. Let ℱ ℎ = ℱ ∪ ℱ , in which ℱ and ℱ are the subset of interior facets and boundary facets, respectively. We denote the boundary of a cell by , and the outward unit normal vector on by . The set of faces of the element is denoted by ℱ . ( ) ( ≥ 0) is the space of polynomials of degree at most on a domain . Introduce the trace operator : (Ω) → −1/2 (ℱ ℎ ) ( ≥ 1), by restricting functions in (Ω) to ℱ ℎ . Next, introduce the following discontinuous finite element spaces on Ω: and the following discontinuous facet finite element spaces on 0 : with ≥ 1.

Weak formulation
Next, based on the embedded, hybridized, and embedded-hybridized discontinuous Galerkin methods for the Stokes equations in [25,28] and for Navier-Stokes equations in [19,26], we present a straightforward extension to the Oseen equations. Firstly, define the following multilinear forms: Here, > 0 is the velocity penalty parameter, and it needs to be chosen sufficiently large to ensure stability [25]. The unified formulation of the embedded, hybridized, and embedded-hybridized discontinuous Galerkin methods reads as follows: given for ∀( ℎ , ℎ ) ∈ ⋆ ℎ , and ⋆ ℎ is given by We denote function pairs in ⋆ ℎ and ⋆ ℎ by boldface, for example, ℎ = ( ℎ ,¯ℎ) ∈ ⋆ ℎ and ℎ = ( ℎ ,¯ℎ) ∈ ⋆ ℎ . The HDG method uses discontinuous velocity and pressure facet function spaces. For the E-HDG method, the velocity facet function space is continuous and the pressure facet function space is discontinuous, and for the EDG method, both velocity and pressure facet functions are continuous. These methods yield the pointwise divergence-free velocity fields on cells. For the HDG and E-HDG methods, ℎ ∈ (div, Ω), and for the EDG method, the normal component of the velocity is only weakly continuous across cell facets, see [28]. In addition, these methods are local mass and momentum conservative. The local momentum conservation is in terms of the numerical flux ℎ = ℎ ⊗ + (¯ℎ − ℎ ) ⊗ +ˆ, ℎ , where is a function that takes on a value of either one or zero as is defined in [19], andˆ, ℎ =¯ℎ − ∇ ℎ − ℎ (¯ℎ − ℎ ) ⊗ , see [19].

Preliminaries
In this section, we give some preliminaries. First, introduce the following composite function spaces , the jump · and average { · } operators across the cell boundary , ∀ ∈ ℎ , are defined by Define the broken Sobolev space We define the following norms on ⋆ (ℎ) and ⋆ (ℎ), respectively where ||| · ||| ′ and ||| · ||| are equivalent on ⋆ ℎ ([25], (28)), and We introduce the following seminorm Introduce the following discrete Poincaré inequality where is a positive constant independent of ℎ, see Proposition A.2 of [6]. Next, we give the following energy stability estimate. By testing (2.4) with ( ℎ , ℎ ), using the discrete coercivity in Lemma 3.1 and the stability of ℎ (3.10) on the left-hand side, and Cauchy-Schwarz inequality and Young's inequality on the right-hand side, we can obtain the estimate Let 0 and 1 , ∞, we have the local inverse inequality Lemma 1.138 of [9] The following continuous and discrete trace inequalities will be used We denote the Lagrange interpolant of order of a continuous function by ℐ ℎ , then the following bound holds Theorem 4.4.4 of [4]: where > / when 1 < ≤ ∞ and ≥ when = 1. Given a domain ⊂ R ( = 2, 3) and, for an integer ≥ 0, we define the 2 -orthogonal projector : Now, we present the stability and boundedness of the multilinear forms which ensures that the method is well-posed, and the consistency of the method. The following lemmas hold true for the EDG, E-HDG and HDG formulations, see [18,25,28] for more details. > 0 and > 0, independent of ℎ and , such that for all ℎ ∈ ⋆ ℎ and ∈ ⋆ (ℎ), Proof. The proof of this lemma was provided in the proof of Lemma 4.8 in [25].
Proof. See the proof of Lemma 4.7 of [25] for the HDG formulation. It is also true for the EDG and E-HDG formulations by following the proof.

Error estimates for the velocity
In this section, we present the velocity error analysis for the EDG, E-HDG and HDG methods of the Oseen problem in a unified setting.
Proof. Firstly, by subtracting (2.4) from (3.11), we can obtain By using (4.3) and taking ( ℎ , ℎ ) = ( , ), then we make the error splitting to obtain where we use the equation ℎ ( , ) = ℎ ( , ) = 0 from (4.1). By applying the discrete coercivity of ℎ in Lemma 3.1 and the stability of ℎ (3.10) on the left-hand side of (4.5), we obtain For the second term of the right-hand side of (4.6), for the HDG and E-HDG methods, ℎ ( , ) = 0 with , which is pointwise divergence-free and H(div)-conforming, and for the EDG method, we apply Young's inequality to obtain Denote Ψ = ℎ ( ; , ), then we have For the term Ψ 1 , applying H¨lder's inequality and Young's inequality, we have (4.8) For the term Ψ 2 , applying Young's inequality and the triangle inequality, we can obtain Collecting the above estimates, we can obtain Inserting (4.7) and (4.9) into (4.6), we can finish the proof.
Next, by combining Lemma 4.1, the following Theorem gives an optimal error estimate with convergence order + 1 for the velocity 2 error, in which the constants are dependent on the negative power of the viscosity.
Proof. First, note the condition that when < 1, ∀ ∈ , ∀ ∈ ℎ , namely, = | · |ℎ , ∀ ∈ , ∀ ∈ ℎ . Using the condition, Lemma 4.1, the interpolation estimates (4.2), (3.5) and (3.6), and the triangle inequality, we can conclude the proof. For the continuous interior penalty method in [5], an optimal 2 -error estimate for the velocity is proved when the local Reynolds number is low, by giving an adjoint problem with a regularity hypothesis.
The Raviart-Thomas interpolation operator satisfies the following approximation properties.
Proof. By using (4.12) and taking ( ℎ , ℎ ) = ( , ), we make the error splitting to obtain (4.13) Note that ℎ ( , ) = 0 with , which is pointwise divergence-free and H(div)-conforming. Then, we apply the discrete coercivity of ℎ in Lemma 3.1 and the stability of ℎ (3.10) on the left-hand side of (4.13). On the right-hand side of (4.13), we use Cauchy-Schwarz inequality to get (4.14) For the second term of the right-hand side of (4.14), by applying boundedness of ℎ in Lemma 3.1 and Young's inequality, we have Then, for the third term of the right-hand side of (4.14), as proceeded in (4.7), we have Next, we give a bound for the term Φ = ℎ ( ; , ). Firstly, Then, by using (4.10), we have where we recall the definition of 0 , and For the term Φ 1 , using (4.17) and applying H¨lder's inequality, inverse inequality and Young's inequality, we have (4.18) For Φ 2 , we apply Young's inequality to obtain Collecting the above estimates, we can obtain Thus, inserting (4.15), (4.16) and (4.19) into (4.14), we can finish the proof.
Next, the following Theorem gives the main results of the paper, which states that the velocity errors are Renolds-robust with the quasi-optimal convergence order + 1/2 at small viscosity.
Then, there exists a constant 2 > 0, independent of ℎ and −1 , such that when Here, the constant 2 is dependent on the physical parameters and .

Error estimates for the pressure
In this section, we give an error estimate for the pressure, which is based on the velocity error estimate in previous section.
Here, the constant is dependent on the physical parameters and .
Proof. Next, we use (4.12) to split the errors for the velocity and the pressure, respectively, then we have By the discrete inf-sup stability (3.9), (5.1), and noting that ℎ ( ℎ , ) = 0, we have By boundedness of ℎ in Lemma 3.1 and boundedness of ℎ (3.7), Using Cauchy-Schwarz inequality and the Poincaré inequality (3.1), we have Note that we will not track the dependency of the error estimates on unless it's necessary. For the term Φ ′ 1 , we have For the term Φ ′ 2 , applying Cauchy-Schwarz inequality, the triangle inequality and the trace inequality, we obtain Collecting the above estimates, and using the triangle inequality, the trace inequality and inverse inequality, we can obtain By the triangle inequality and (5.7), we can conclude the proof.

Remark 5.2.
Here, we will comment on the regularity for the exact pressure for the HDG and E-HDG methods. The related pressure doesn't appear in the velocity error estimates for the two methods, thus, the higher regularity ∈ +1 (Ω) is not needed in Theorems 4.2 and 4.6. In fact, we only need the regularity ∈ (Ω) to prove the optimal error estimates of the pressure in Theorem 5.1, where the interpolant of the pressure can be replaced byĪ = ℐ | ℱ ℎ ∈¯ℎ. Here, ℐ is the continuous Scott-Zhang interpolant [4].

Numerical studies
In this section, the numerical experiments are devised to verify our theoretical estimates. In particular, we recover the convergence rates for the velocity and the pressure. For the ample numerical performance of these types of methods, we can refer to the references [19,[26][27][28].
In our implementation, the velocity penalty parameter is chosen to be 6 2 . We choose the polynomial order = 2. We use the meshes, in which a sequence of the regular triangulations with diagonals (from bottom right to top left), with the same number of subdivisions on each coordinate direction was generated, then each of these triangulations was barycentrically refined. Here, 'EOC' represents the average estimated orders of convergence and 'ndof' represents the number of global degrees of freedom on the element interfaces. All numerical experiments are implemented in this NGSolve software [29]. Table 1. Velocity errors in the 2 -norm of the HDG, E-HDG and EDG methods with varying .    Table 2. Pressure errors in the 2 -norm of the EDG method with varying . Since this analytical solution is independent of the viscosity , it allows us to easily assess the approximation properties for a wide range of the Reynolds numbers by changing the value of . We use the above-mentioned meshes with = 6, 12, 24, 48 subdivisions in each coordinate direction. We test the convergence orders of the velocity errors in the 2 -norm for the HDG, E-HDG and EDG methods with varying viscosity, respectively. From Table 1, we can observe the same convergence behavior for the HDG, E-HDG and EDG methods. The velocity error has an optimal convergence rate for large values of , and for small values of , it has the quasioptimal convergence rate, as we predicted in Theorems 4.2 and 4.6. In addition, we can fix the mesh size to observe that as the viscosity decreases, the velocity errors become larger and larger, and when the viscosity is small enough, the velocity errors are almost unchanged, namely, the velocity errors are independent of the small viscosity. We can also note that the HDG method with more degrees of freedom has no obvious superiority in terms of accuracy than the E-HDG and EDG methods on a given mesh. In addition, we also test the convergence orders of the pressure errors in the 2 -norm of the EDG method with varying . From Table 2, we can observe that the pressure errors have a super-convergence order and better accuracy at small viscosity. The HDG and E-HDG methods have the same convergence behavior as that of the EDG method (for brevity not shown here).
Next, we test the HDG, E-HDG and EDG methods with different pressure by varying . We set = 1 − 3. The above-mentioned mesh with = 50 subdivisions in each coordinate direction, is used. Here, | ℎ | nj = Table 3. Velocity errors in the 2 -norm of the HDG, E-HDG and EDG methods with varying the pressure.  Table 3 that the velocity errors for the HDG and E-HDG methods are indeed independent of the pressure. For the EDG method, the velocity error isn't pressure-robust, and the velocity isn't H(div)-conforming.

Conclusions and future directions
In this paper, we analyzed the EDG method for the Oseen equations. The a priori error estimates are given for the velocity and the pressure. The velocity error in the 2 (Ω) norm is Reynolds-robust with convergence order of + 1/2 in the convection-dominated regime. In the diffusion-dominated regime, the 2 error of the velocity has an optimal convergence order of + 1. In addition, the convergence rate of the pressure is also proved. The analysis here also covers the pressure-robust HDG and E-HDG methods.
In view of the attractive properties of the EDG method, a further work is to recover the pressure-robustness of the EDG method for the Oseen problem. As we see, in [20], a simple reconstruction operator is introduced to recover the pressure-robustness of the EDG method for the Stokes equations. Next, we will consider whether the simple reconstruction operator is applicable to the EDG method for the Oseen problem, or whether a new reconstruction operator can be introduced, in which the 2 (Ω) error of the velocity is proved to be Reynoldsrobust and pressure-robust with convergence order of + 1/2 in the convection-dominated regime.