Robust error estimates for stabilized finite element approximations of the two dimensional Navier-Stokes equations with application to implicit large eddy simulation
aa r X i v : . [ m a t h . NA ] A p r ROBUST ERROR ESTIMATES FOR STABILIZED FINITE ELEMENTAPPROXIMATIONS OF THE TWO DIMENSIONALNAVIER-STOKES EQUATIONS WITH APPLICATION TO IMPLICITLARGE EDDY SIMULATION
ERIK BURMAN ∗ Abstract.
We consider error estimates in weak parametrised norms for stabilized finite elementapproximations of the two-dimensional Navier-Stokes’ equations. These weak norms can be relatedto the norms of certain filtered quantities, where the parameter of the norm, relates to the filterwidth. Under the assumption of the existence of a certain decomposition of the solution, into largeeddies and fine scale fluctuations, the constants of the estimates are proven to be independent of boththe Reynolds number and the Sobolev norm of the exact solution. Instead they exhibit exponentialgrowth with a coefficient proportional to the maximum gradient of the large eddies. The errorestimates are on a posteriori form, but using Sobolev injections valid on finite element spaces andthe properties of the stabilization operators the residuals may be upper bounded uniformly, leadingto robust a priori error estimates.
Key words.
Navier-Stokes’ equations, stability, error estimates, large eddy simulation, finiteelement methods, stabilization
AMS subject classifications.
1. Introduction.
In this paper we will be interested in stabilized finite elementmethods in the context of so called implicit large eddy simulation (ILES), see [2]. Thisis a numerical approach to the computation of turbulent flow where no modelling ofthe Reynolds stresses is performed on the continuous level. Instead the Navier-Stokes’equations are approximated numerically using a method that dissipates sufficient en-ergy on the scale of the mesh size. This eliminates the buildup of energy that createsspurious oscillations in any energy conservative approximation method. It has beenargued that the truncation error of such methods by itself may act as a subgrid model[1, 18] and there exists numerical evidence that ILES methods work for the simulationof two dimensional turbulence, provided back scatter effects are not strong [16]. Thereis also numerical evidence of the potential for adaptive LES/DNS driven by adaptive,stabilized finite element simulations, see [12, 13] and [21].Our objective in this paper is to provide a numerical analysis for stabilized finiteelement methods under minimal regularity assumptions and to provide sufficient con-ditions on the exact solution for the derivation of rigorous error estimates that areindependent of both the Reynolds number and Sobolev norms of the exact solution. Itis well known that provided the exact solution is sufficiently smooth the approximatesolution u h of the Navier-Stokes’ equations on velocity-pressure form can be provedto satisfy estimates of the type k u − u h k L (Ω) . e k∇ u k ∞ h | u | L ( I ; H (Ω)) , (1.1)if a consistent stabilized finite element method with piecewise affine approximation isused. Here u denotes the flow velocity and h the mesh size. See [11, 7] for examples ofanalyses of Navier-Stokes’ equations on velocity-pressure form and [15, 17] for analyseson velocity-vorticity form. We also give a proof of (1.1) for one the methods proposedherein in appendix. Here we use the notation a . b for a ≤ Cb with C a constant ∗ Department of Mathematics, University College London, Gower Street, London, UK–WC1E6BT, United Kingdom; ( [email protected] ) 1
E. BURMAN independent of the physical parameters of the problem, unless they can be expectedto have O (1) contribution, it can also include some dependence on initial data, thatmay be assumed to be O (1). We will also use a ∼ b for a . b and b . a . Note thatthere is no explicit dependence on the viscosity in the estimate (1.1). For this estimateto be useful the included Sobolev norms must be small, which rarely is the case inthe high Reynolds number regime and hence the dependence of the viscosity entersin an implicit manner. The purpose of the present paper is to propose an alternativeapproach, where the estimate is indeed independent of the Reynolds number, both inthe sense that the estimate is free from inverse power of the viscosity in the upperbound, but also that the dependence on unknown Sobolev norms of the exact solutionis strongly reduced. It does not seem possible to eliminate this dependence completely,due to the possible presence of backscatter. Observe that the presence of k∇ u k ∞ inthe exponential of (1.1) reflects the effect of diverging characteristics in the transportequation and is present already in the linear convection–diffusion equation at highP´eclet numbers. We may define a timescale for the flow separation, τ := k∇ u k − ∞ .The reason this time scale becomes so small is that it will be the smallest timescaleof the flow and hence equal to the micro time, because of fine scale fluctuations ofthe velocity. From the physical point of view it is argued that LES will be successfulfor flows where both the quantities of interest and the rate-controlling processes aredetermined by the resolved large scales , see Pope [20] for a discussion. We will use thisidea as a starting point for our assumptions on the flow.To derive error estimates for a numerical method we need the following:– continuous dependence on data, independent of the exact solution;– some smooth quantity that we can apply approximation estimates to.At a first glance both these prerequisites appear to fail for the two-dimensional Navier-Stokes’ equation. The first fails because of the presence of the exponential factor andthe second fails because Sobolev norms of the exact solution can be huge for smallviscosities. The following three points allow us to break this deadlock:1. the use of a parametrized weak norm, corresponding to measuring the errorin filtered quantities of the solution;2. introduction of an assumption on the structure of the exact solution that issufficient for an implicit large eddy simulation to be robust;3. a stabilized finite element method, giving enhanced a priori control of residualquantities in the high Reynolds regime.The idea of measuring error in filtered quantities was considered in [9, 10], but theestimates were not robust in the Reynolds number and the constant included highorder Sobolev norms of the exact solution. In [4] weak norm estimates were usedin order to derive robust estimates for the Burgers’ equation, where the constant inthe right hand side only depended on initial data. The second point, which was notnecessary in the case of the Burgers’ equation, reflects the difficulty to characterisethe solution structure in higher dimension and the ensuing need of some structuralassumptions in order to rule out strong backscatter effects. The third point allows usto control residual quantities independent of the viscosity.As a first approximation it is reasonable to assume that for a solution to beamenable to large eddy simulation, there are relatively smooth eddies, with largeassociated Reynolds number, containing the bulk of the energy and small scale fluctu-ations that may vary rapidly in space, but carries a negligible part of the energy. Tomake this precise, we assume that there exists a decomposition of the exact solutionin the spatially slowly varying large scales and an arbitrarily rough fine scale, with obust error estimates for the two dimensional Navier-Stokes equations u = ¯ u + u ′ , ¯ u, u ′ ∈ W , ∞ (Ω) . We then assume that the Reynolds number associated to the large scales, Re may belarge, but k ¯ u k W , ∞ ( Q ) ∼
1, whereas for the fine scales k u ′ k W , ∞ ( Q ) may be large, butthe energy small. To give a precise meaning to small here, we introduce a global timescale for the flow, defined using the large scales τ F := k ¯ u k − W , ∞ ( Q ) ∼ . This is in agreement with the statement that rate controlling processes are deter-mined by the large scale. Using the viscosity coefficient we may make the followingassumption on the energy content of the small scales k u ′ k ∞ ∼ ν/τ F . The length scale based on | u ′ | and τ F writes l ′ := | u ′ | τ F ∼ | u ′ | ν | u ′ | and it follows that the small scale Reynolds number is Re ′ := u ′ l ′ ν ∼ | u ′ | ν | u ′ | ν ∼ . (1.2)Alternatively one may assume that the fine scale Reynolds number is one and thatthe large scale characteristic time, is the globally relevant time scale and then derivethe bound on the energy. We will refer to the above as the large eddy assumption .Under this assumption we prove the following bound on the approximate velocitiessup t ∈ (0 ,T ) k u − u h k L (Ω) . h . (1.3)The hidden constant in the above estimate only depends on initial data (maximuminitial vorticity) and the mesh geometry, but is independent of the Reynolds numberand Sobolev norms of the exact solution. The discussion is limited to two spacedimensions and hence we do not properly speaking address the question of turbulentflows.Let us end this introductory discussion by emphasising that what we compute is anapproximation to the solution of the Navier-Stokes’ equations. For this approximatesolution we can prove that provided τ F is not too small, corresponding to slowlyvarying large scale velocity field, the filtered part of the vorticity is stable underperturbations resulting in robust error estimates in weak norms for vorticity. Usingthese estimates we may then control the L -norm of the velocity error as shown above.Herein our main concern will be the high mesh Reynolds number case Re h := U hν > , where U := k u h ( · , k L ∞ (Ω) ∼ E. BURMAN
Reynolds number is low, other approaches than those presented herein might be moreappropriate. Let us also point out that another feature of our estimates is that theyprovide the first error estimates with an order in h for nonlinear stabilization schemes,satisfying a discrete maximum principle, in two space dimensions.We will consider the Navier-Stokes’ equations written on vorticity-velocity form.Let Ω be the unit square and assume that the boundary conditions are periodic inboth cartesian directions. The L -scalar product over some space-time domain willbe denoted ( · , · ) X with associated norm k · k X where the subscript may be dropped for X = Ω. Define the time interval I := (0 , T ) and the space-time domain Q := Ω × I .The equations then writes, ∂ t ω + ∇· ( uω ) − ν ∆ ω = 0 , in Q, − ∆Ψ = ω in Q, (1.4) u = rot Ψ in Q,ω ( x,
0) = ω , with ω ∈ L ∞ (Ω). The associated weak formulation takes the form, for t >
0, find( ω, Ψ) ∈ H (Ω) × H (Ω) ∩ L ∗ (Ω) such that( ∂ t ω, v ) Ω + a ( u ; ω, v ) = 0 , (1.5)( ∇ Ψ , ∇ Φ) Ω = ( ω, Φ) Ω , (1.6) u = rot Ψ , for all ( v, Φ) ∈ H (Ω) × H (Ω) ∩ L ∗ (Ω), where the semi-linear form a ( · ; · , · ) is definedby a ( u ; ω, v ) := ( ∇· ( uω ) , v ) Ω + ( ν ∇ ω, ∇ v ) Ω . This problem is known to be well-posed, but a priori error estimates on the solutionare in general strongly dependent on the viscosity coefficient reflecting the possiblepoor stability of the equations in the high Reynolds number regime.
2. Finite element discretization.
Let {T h } h> be a family of affine, simplicialDelaunay meshes of Ω. We assume that the meshes are kept fixed in time and that thefamily {T h } h> is quasi-uniform. Mesh faces are collected in the set F . For a smoothenough function v that is possibly double-valued at F ∈ F with F = ∂T − ∩ ∂T + ,we define its jump at F as [[ v ]] := v | T − − v | T + , and we fix the unit normal vectorto F , denoted by n F , as pointing from T − to T + . The arbitrariness in the sign of[[ v ]] is irrelevant in what follows. Define V h to be the standard space of piecewisepolynomial, continuous periodic functions, V kh := { v h ∈ H (Ω) : v h | K ∈ P k ( K ); ∀ K ∈ T h ; v h periodic in x and y } . The set of gradients of functions in V kh will be denoted by W k − h := { w h = ∇ v h ; v h ∈ V kh } . Let L := L (Ω) and set L ∗ := { q ∈ L ; R Ω q = 0 } . Let V ∗ := V kh ∩ L ∗ . We let π L denote the L -projection on V kh and π V the H -projection( ∇ π V u, ∇ v h ) Ω = ( ∇ u, ∇ v h ) Ω ∀ v h ∈ V h and Z Ω ( π V u − u ) d x = 0 . obust error estimates for the two dimensional Navier-Stokes equations π L and π V , k π L u − u h k + h k∇ ( π L u − u ) k ≤ c h s | u | s , with 1 ≤ s ≤ k + 1 (2.1)and k π V u − u h k + h k∇ ( π V u − u ) k ≤ c h s | u | s , with 1 ≤ s ≤ k + 1 . (2.2)We consider continuous finite elements with k = 1 to discretize the vorticity ω in spaceand k = 1 , u h | K := rot Ψ h := ( ∂ y Ψ h , − ∂ x Ψ h ). Note that using this definition ∇ · u h = 0 in Ω,i.e. the discrete velocity is globally divergence free. We discretise in space using astabilized finite element method. For t > ω h , Ψ h ) ∈ V h × V l ∗ , with l = 1 , ∂ t ω h , v h ) M + a ( u h ; ω h , v h ) + s ( u h ; ω h , v h ) = 0 , (2.3)( ∇ Ψ h , ∇ Φ h ) Ω − ( ω h , Φ h ) Ω = 0 , (2.4) u h − rot Ψ h = 0 , for all ( v h , Φ h ) ∈ V h × V ∗ and with initial data w := π L ω ( · , s ( · ; · , · ) denotes astabilization operator that is linear in its last argument and ( ∂ t ω h , v h ) M denotes thebilinear form defining the mass matrix, this operator either coincides with ( · , · ) Ω or isdefined as ( · , · ) Ω approximated using nodal quadrature, i.e. so called mass lumping.We will assume the stabilization term satisfiesinf v h ∈ V h k h ( u h · ∇ ω h − v h ) k . s ( u h , ω h ; ω h ) . h ( U + k u h k L ∞ (Ω) ) k∇ ω h k (2.5)and s ( u h , ω h ; v h ) . h ( U + k u h k L ∞ (Ω) ) s ( u h , ω h ; ω h ) k∇ v h k . (2.6)The formulation (2.3)-(2.4) satisfies the following stability estimates Lemma 2.1. sup t ∈ I k ω h ( · , t ) k M + 2 k ν ∇ ω h k Q + 2 Z I s ( u h ; ω h , ω h ) d t ≤ k ω h ( · , k M , (2.7) and if exact integration is used for ( · , · ) M , k u h ( · , T ) k M + 2 k ν ω h k Q = k u h ( · , k M − Z T s ( u h ; ω h , Ψ h ) d t, (2.8) k u h ( · , t ) k L ∞ (Ω) ≤ c q k ω h ( · , t ) k L q (Ω) , q > and for l = 1 , Z T k∇ ∂ t ω h k d t . Z T ( h − ( U + k u h k L ∞ (Ω) ) s ( u h , ω h , ω h ) + νh − k∇ ω h k ) d t. (2.10) Proof . Inequality (2.7) is immediate by taking v h = ω h in (2.3). Inequality (2.8)is obtained by taking v h = Ψ h in the equation (2.3) and deriving the equation (2.4) E. BURMAN in time and taking Φ h = ω h . For the inequality (2.9), consider the auxiliary problem, − ∆ ˜Ψ = ω h in Ω and note that by [22] there holds k u h ( · , t ) k L ∞ (Ω) . k ˜Ψ( · , t ) k W , ∞ (Ω) and adapting the analysis of [19] we have for the (simpler) case or periodic boundaryconditions, k ˜Ψ( · , t ) k W , ∞ (Ω) ≤ c q k ω h ( · , t ) k L q (Ω) , q > . To prove (2.10) finally we introduce a function ξ h ∈ V h such that( ξ h , v h ) M = ( ∇ ∂ t ω h , ∇ v h ) Ω , ∀ v h ∈ V h , it follows by taking v h = ξ h and using the Cauchy-Schwarz inequality followed by aninverse inequality that ( ξ h , ξ h ) M ∼ k ξ h k . h − k ∂ t ∇ ω h k . (2.11)Observe that by norm equivalence on discrete spaces the L -norm defined using nodalquadrature is equivalent to the consistent L -norm. Taking v h = ξ h in (2.3) yields k ∂ t ∇ ω h k = − ( u h · ∇ ω h , ξ h ) Ω − ( ν ∇ ω h , ∇ ξ h ) Ω − s ( u h , ω h , ξ h ) . We may then apply the Cauchy-Schwarz inequality in the second term of the righthand side and (2.6) in the last term, followed by inverse inequalities on k∇ ξ h k andthe estimate (2.11). For the first term we write, using the properties of ξ h and thebound | ( v h , ξ h ) M − ( v h , ξ h ) Ω | . ( h |∇ v h | , |∇ ξ h | ) Ω , | ( u h · ∇ ω h , ξ h ) Ω , | . | ( u h · ∇ ω h − v h , ξ h ) Ω | + | ( ∂ t ∇ ω h , ∇ v h ) Ω | + ( h |∇ v h | , |∇ ξ h | ) Ω Since both u h and ∇ ω h are constant per element ∇ v h | K = ∇ ( v h − u h · ∇ ω h ) | K . Usinginverse inequalities and the bound (2.11) on ξ h we have | ( u h ·∇ ω h − v h , ξ h ) | + | ( ∂ t ∇ ω h , ∇ v h ) | +( h |∇ v h | , |∇ ξ h | ) . h − k ∂ t ∇ ω h kk u h ·∇ ω h − v h k . The claim follows by the inequality (2.5), (2.11) and finally by integrating in time.It follows from (2.8) that the method is energy consistent if s ( u h ; ω h , Ψ h ) = 0. Takingthe difference of the formulations (1.5) - (1.6) (with v = v h ) and (2.3) - (2.4) andsetting e ω = ω − ω h and e Ψ = Ψ − Ψ h , the following consistency relation holds( ∂ t e ω + u ·∇ e ω + rot e Ψ ·∇ ω h , v h ) Ω + ( ν ∇ e ω , ∇ v h ) Ω = ( ∂ t ω h , v h ) M − ( ∂ t ω h , v h ) Ω + s ( u h , ω h ; v h ) in Q, (2.12)( ∇ e Ψ , ∇ Φ h ) Ω − ( e ω , Φ h ) Ω = 0 in Q. As mentioned in the introduction, if the solution ( u, ω ) is smooth one may prove anerror estimate that is robust with respect to ν using standard linear theory and per-turbation arguments. For the methods we consider herein, this result is an extensionof the works in [17] and [7] and we state it here only with the dominant terms present. obust error estimates for the two dimensional Navier-Stokes equations Proposition 2.2.
Let ( u, ω ) be a smooth solution of (1.5) - (1.6) and ( u h , ω h ) be the solution of (2.3) - (2.4) , where the stabilization operator satisfies the additionalweak consistency property s ( u h ; π L ω, π L ω ) ≤ c ( u, ω ) h then for l = 1 , k ( u − u h )( · , T ) k + k ( ω − ω h )( · , T ) k . c ω ( h | ω | L ( I ; H (Ω)) + h l | Ψ | L ∞ ( I ; H l +1 (Ω)) ) where c ω := e k∇ ω k L ∞ ( Q ) T . In addition there holds for the stabilization operator s ( u h ; ω h , ω h ) . c ω ( h | ω | L ( I ; H (Ω)) + h l | Ψ | L ∞ ( I ; H l +1 (Ω)) ) . Observe that the exponential factor here depends on k∇ ω k L ∞ (Ω) , compared to k∇ u k L ∞ (Ω) in (1.1). This is the prize we pay for estimating the L -error of thevorticity. As we shall see below, the use of weaker norms for the estimation of ω h allows us to revert back to the exponential factor of (1.1) and under the large eddyassumption, the exponential growth is moderate.
3. Dual problem.
From the consistency relation (2.12) we deduce the following(homogeneous) perturbation formulation for the evolution of ( e ω , e Ψ )( ∂ t e ω + u ·∇ e ω + rot e Ψ ·∇ ω h , ϕ ) Ω + ( ν ∇ e ω , ∇ ϕ ) Ω = 0 in Q, (3.1)( ∇ e Ψ , ∇ ϕ ) Ω − ( e ω , ϕ ) Ω = 0 in Q, (3.2)where ϕ , ϕ are the solutions to a dual adjoint perturbation equation related to thecontinuous equation (1.5)-(1.6) and the discretization (2.3)-(2.4). Since the jump ofthe tangential derivative of ω h is zero, we may integrate by parts in (3.1), to arriveat the dual adjoint problem − ∂ t ϕ − u ·∇ ϕ − ϕ − ν ∆ ϕ = 0 in Q, (3.3) − ∆ ϕ − ∇ ω h · rot ϕ = 0 in Q, (3.4) ϕ ( x, T ) = ξ ( x ) in Ω , (3.5)where ξ ( x ) is some initial data to be fixed later, the choice of ξ determines thequantity of interest.A key result for the present analysis is the following stability estimate for the dualadjoint solution Proposition 3.1.
The following stability estimate holds for the solution ( ϕ , ϕ ) of (3.3) - (3.5) , sup t ∈ I k∇ ϕ ( · , t ) k + k ν D ϕ k Q . e TτF k∇ ξ k (3.6) Z I k∇ ϕ ( · , t ) k d t ≤ e TτF Z I k ω h k L ∞ (Ω) d t k∇ ξ k (3.7) where τ F is defined in the proof. If the large eddy assumption holds τ F ∼ . E. BURMAN
Proof . First multiply (3.3) by − ∆ ϕ and (3.4) by ϕ and integrate over Q ∗ :=Ω × ( t ∗ , T ), where t ∗ is a time to be chosen. By summing the two relations we obtain( ∂ t ϕ , ∆ ϕ ) Q ∗ | {z } I + ( u ·∇ ϕ , ∆ ϕ ) Q ∗ | {z } I + ( ∇ ω h · rot ϕ , ϕ ) Q ∗ | {z } I + k ν ∆ ϕ k Q ∗ | {z } I = 0 . We will now treat the terms I - I term by term. First note that by integration byparts first in space and then integration in time we have I = − Z Tt ∗ ddt k∇ ϕ ( · , t ) k d t = 12 k∇ ϕ ( · , t ∗ ) k − k∇ ξ k . The second term is handled using the decomposition of u in the large scale and finescale component and then an integration by parts only in the large scale part. Here ∇ S u denotes the symmetric part of the gradient of the vector u . I = − (( ∇ S ¯ u −
12 ( ∇· ¯ u ) I × ) ∇ ϕ , ∇ ϕ ) Q ∗ − ( u ′ ·∇ ϕ , ∆ ϕ ) Q ∗ ≤ Z Q (Λ(¯ u, u ′ , ν ) ∇ ϕ ) T ·∇ ϕ d x d t + 12 k ν ∆ ϕ k Q ∗ , where Λ(¯ u, u ′ , ν ) is a two by two, symmetric matrix defined by,Λ(¯ u, u ′ , ν ) = −∇ S ¯ u + 12 ∇· ¯ u I × + 12 ν u ′ T u ′ . We now define the global timescale τ F of the flow by( τ F ) − := inf ¯ u ∈ L ∞ ( I ; W , ∞ (Ω)) k σ + p (Λ(¯ u, u ′ , ν )) k L ∞ ( Q ) . Here σ + p denotes the largest positive eigenvalue of the matrix. This results in anontrivial minimization problem in L ∞ . We leave the precise study of this problemfor further work and here simply observe that by computing the eigenvalues of thesymmetric part of the gradient tensor we may write( τ F ) − ≤ inf ¯ u J (¯ u, u ′ )where J (¯ u, u ′ ) := sup t ∈ I (cid:16)q k ( ∂ x ¯ u − ∂ x ¯ u ) + ( ∂ x ¯ u + ∂ x ¯ u ) k L ∞ (Ω) + ν − k u ′ k L ∞ (Ω) (cid:17) . We observe that the global stability does not depend on the divergence componentor the rotational of ¯ u , only on the other two components of the velocity gradientmatrix. Since u ′ = u − ¯ u , it follows that we can minimize over all large scale vectorfields ¯ u ∈ [ W , ∞ (Ω)] and the infimum value obtained is the optimal timescale of theflow. Under the assumptions made in the introduction, that k ¯ u k W , ∞ (Ω) ∼ ν − k u ′ k L ∞ (Ω) ∼
1, for all t , we immediately deduce that τ F ∼ ∇· rot ϕ = 0 and ∇ ϕ · rot ϕ =0 we have I = − ( ω h ∇· rot ϕ , ϕ ) Q ∗ − ( ω h rot ϕ , ∇ ϕ ) Q ∗ = 0 . obust error estimates for the two dimensional Navier-Stokes equations I − I we have k∇ ϕ ( · , t ∗ ) k L + k ν ∆ ϕ k Q ∗ ≤ τ − F k∇ ϕ k Q ∗ + k∇ ξ k L . The inequality for ϕ follows after a Gronwall’s inequality and by taking the supremumover t ∗ , resulting in sup t ∈ I k∇ ϕ ( · , t ) k + k ν D ϕ k Q . e TτF k∇ ξ k . Elliptic regularity has been used for the second term.For the bound on ϕ multiply equation (3.4) by ϕ and integrate over Ω, k∇ ϕ ( · , t ) k = − ( ω h rot ϕ , ∇ ϕ ) Ω ≤ k ω h ( · , t ) k L ∞ (Ω) k∇ ϕ ( · , t ) kk∇ ϕ ( · , t ) k . Then divide by k∇ ϕ ( · , t ) k , integrate in time and use that Z I k ω h ( · , t ) k L ∞ (Ω) k∇ ϕ ( · , t ) k d t ≤ Z I k ω h ( · , t ) k L ∞ (Ω) d t sup t ∈ I k∇ ϕ ( · , t ) k . Finally use equation (3.6) to bound the term in k∇ ϕ ( · , t ) k .Note the dependence of ω h in the bound (3.7). This appearance of a finite elementfunction in the stability estimate shows that the global stability depends on the mono-tonicity of the approximation scheme. However as we shall see, strict monotonicity isnot necessary, only L ∞ -control of the vorticity.
4. A posteriori and a priori error estimates for the abstract method.
Let e ω = ω − ω h and let the filtered error ˜ e ω be defined as the solution to the problem − δ ∆˜ e ω + ˜ e ω = e ω . (4.1)We introduce a norm on ˜ e ω such that |k ˜ e ω k| δ := k δ ∇ ˜ e ω k + k ˜ e ω k = ( e ω , ˜ e ω ) Ω . Thisnorm coincides with the L -norm for δ = 0 and is related to the H − -norm for δ = 1.By choosing δ := δ ( h ), i.e. by reducing the filter width with the mesh size, we obtaina family of norms that become stronger as the mesh size is reduced.Using the above norm and the relations (2.12), (3.3)-(3.4) as well as the stabilityresult of Proposition 3.1 we may derive a posteriori estimates for the filtered quantity˜ e ω . We here derive the result for the abstract finite element element method (2.3)-(2.4)and then show how these estimates can be transformed into a priori error estimates,depending on the properties of the stabilization operator s ( u h , ω h ; v h ). The use ofweak norms and stabilized finite element methods in the following estimates draws onideas from [14] and [3, 4]. Theorem 4.1. (A posteriori error estimates) |k ˜ ω − ˜ ω h k| δ . e TτF (cid:18) hδ (cid:19) X i =0 R i , (4.2) with R := k ( ω − ω h )( · , k , R := Z T inf v h ∈ V h k h ( u h ·∇ ω h − v h ) k d t, E. BURMAN R := min( h, ν T ) k ν [[ n F ·∇ ω h ]] k F× I , R := Z T k ω h ( · , t ) k L ∞ (Ω) d t min( c sup t ∈ I k Ψ h ( · , t ) k ∆ , , c h sup t ∈ I k ω h ( · , t ) k ) where k Ψ h ( · , t ) k ∆ ,s := k h s [[ n F · ∇ Ψ h ( · , t )]] k F + inf v h ∈ V lh X K ∈T h k h + s (∆Ψ h ( · , t ) − v h ) k K ! , R := h Z T k ∂ t ∇ ω h k d t and R := ( U + k u h k L ∞ ( Q ) ) Z T s ( u h ; ω h , ω h ) d t. The term R is omitted if the consistent mass matrix is used. For the velocities wehave the estimate, for all t ∈ I , k ( u − u h )( · , t ) k . (cid:16) k Ψ h ( · , t ) k ∆ , + |k (˜ ω − ˜ ω h )( · , t ) k| (cid:17) (4.3) where the second term in the right hand side may be a posteriori bounded by taking δ = 1 in (4.2) .Proof . By the definition of ˜ e ω we have, taking ξ = ˜ e ω in (3.5), |k ˜ e ω k| δ = ( e ω ( T ) , ϕ ( T )) Ω + ( e ω , − ∂ t ϕ − u ·∇ ϕ − ν ∆ ϕ ) Q + ( e Ψ , − ∆ ϕ − ∇ ω h · rot ϕ ) Q = ( e ω (0) , ϕ (0)) Ω + ( ∂ t e ω + u ·∇ e ω + rot e Ψ ·∇ ω h , ϕ ) Ω + ( ν ∇ e ω , ∇ ϕ ) Ω + ( ∇ e Ψ , ∇ ϕ ) Ω − ( e ω , ϕ ) Ω . Using now the consistency relation (2.12) we obtain |k ˜ e ω k| δ = ( e ω (0) , ( ϕ − π L ϕ )( · , Ω + ( ∂ t e + u ·∇ e + rot e Ψ ·∇ ω h , ϕ − π L ϕ ) Q + ( ν ∇ e, ∇ ( ϕ − π L ϕ )) Q + ( ∇ e Ψ , ∇ ( ϕ − Π ϕ )) Q − ( e, ϕ − Π ϕ ) Q − ( ∂ t ω h , π L ϕ ) M,Q + ( ∂ t ω h , π L ϕ ) Q − s ( u h , ω h ; π L ϕ ) , where Π : H (Ω) V lh will be taken as either π L or π V . Using the equations (1.5)-(1.6) and the definitions of the projections π L and π V we deduce for Π := π V , |k ˜ e ω k| δ = ( e ω (0) , ( ϕ − π L ϕ )( · , Ω − ( u h ·∇ ω h − v h , ϕ − π L ϕ ) Q − ( ν ∇ ω h , ∇ ( ϕ − π L ϕ ) Q + ( ω h , ϕ − π V ϕ ) Q − ( ∂ t ω h , π L ϕ ) M,Q + ( ∂ t ω h , π L ϕ ) Q − Z T s ( u h , ω h ; π L ϕ ) d t, and similarly for Π := π L , |k ˜ e ω k| δ = ( e ω ( · , , ( ϕ − π L ϕ )( · , Ω − ( u h ·∇ ω h − v h ) , ϕ − π L ϕ ) Q − ( ν ∇ ω h , ∇ ( ϕ − π L ϕ ) Q − ( ∇ Ψ h , ∇ ( ϕ − π L ϕ )) Q − ( ∂ t ω h , π L ϕ ) M,Q + ( ∂ t ω h , π L ϕ ) Q − Z T s ( u h , ω h ; π L ϕ ) d t. obust error estimates for the two dimensional Navier-Stokes equations π V , |k ˜ e ω k| δ . (cid:18) hδ (cid:19) (cid:16) k e ω ( · , k + Z T inf v h ∈ V h k h ( u h ·∇ ω h − v h ) k d t + min( h,ν T ) k ν [[ n F ·∇ ω h ]] k F× I + c h sup t ∈ I k ω h ( · , t ) k Z T k ω h ( · , t ) k L ∞ (Ω) d t + h Z T k ∂ t ∇ ω h k d t + ( U + k u h k L ∞ ( Q ) ) Z T s ( u h ; ω h , ω h ) d t (cid:17) × (sup t ∈ I k δ ∇ ϕ ( · , t ) k + k δν D ϕ k Q ) . If Π := π L the fourth term on the right hand side is replaced using( ∇ Ψ h , ∇ ( ϕ − π L ϕ )) Q . (cid:18) hδ (cid:19) c sup t ∈ I k Ψ h ( t ) k ∆ , Z T k δ ∇ ϕ ( · , t ) k d t, followed by the bound (3.7) on ϕ . The estimate (4.2) now follows by taking theminimum of the two expressions and noting that by (3.6)sup t ∈ I k δ ∇ ϕ ( · , t ) k + k δν D ϕ k Q . e TτF |k ˜ e ω k| δ . The velocity estimate (4.3) is obtained by noting that, with e Ψ := Ψ − Ψ h , k u − u h k := k∇ e Ψ k = ( ∇ e Ψ , ∇ ( e Ψ − π L e Ψ )) + ( e ω , π L e Ψ ) . Using the equation (1.6) we have k∇ e Ψ k = ( ∇ Ψ h , ∇ ( e Ψ − π L e Ψ )) + ( ω, e Ψ ) − ( ω h , π L e Ψ )= ( ∇ Ψ h , ∇ ( e Ψ − π L e Ψ )) + ( ω − ω h , e Ψ )Let ˜ e be the solution of (4.1) with δ = 1. Then k u − u h k = ( ∇ Ψ h , ∇ ( e Ψ − π L e Ψ )) Ω + ( ∇ ˜ e ω , ∇ e Ψ ) Ω + (˜ e ω , e Ψ ) Ω . By an integration by parts in the first term, followed by a Cauchy-Schwarz inequalityand the Poincar´e-Friedrichs inequality in the last term we may write k u − u h k . k h F [[ ∇ Ψ h ]] k F k h − ( e Ψ − π L e Ψ ) k F + X K ∈T h k h (∆Ψ h − v h ) k K ! k h − ( e Ψ − π L e Ψ ) k + |k (˜ ω − ˜ ω h ) k| k ( u − u h ) k . By elementwise trace inequalities and the approximation property (2.1) we have k h − ( e Ψ − π L e Ψ ) k F + k h − ( e Ψ − π L e Ψ ) k . k u − u h k E. BURMAN by which we conclude.If the stability properties of the stabilized method are sufficient, these a posteriorierror estimates translate into a priori error estimates. We propose two strategies forthis. One using stability concepts based on Sobolev injections for discrete spaces andone based on monotonicity, applicable to monotone stabilized finite element methodsand monotone implicit large eddy methods. The advantage of the former is that itallows the derivation of a priori error estimates for quasi linear terms s ( u h ; ω h , v h ) andthe use of the consistent mass matrix. The latter technique on the other hand allowsfor the derivation of a priori error estimates with precise control of the constants in theestimates. We will use the notion of the discrete maximum principle (DMP) and theassociated, DMP-property of the forms defining a finite element method introducedin [6]. Proposition 4.2.
Assume that the mass ( · , · ) M is evaluated exactly and that inaddition to (2.5) and (2.6) the following stability estimate holds for all t > , k ω h k L ∞ (Ω) . c ( h )( k ω h k + s ( u h , ω h , ω h ) ) . (4.4) Then there holds for all ǫ > , |k (˜ ω − ˜ ω h )( T ) k| δ . e TτF (cid:18) hδ (cid:19) ( c ǫ h − ǫ + c ( h ) h ) and k ( u − u h )( · , T ) k . inf v h ∈ W l − h k ( u − v h )( · , T ) k + e TτF h ( c ǫ h − ǫ + c ( h ) h ) . Proof . First we recall that k ω h ( · , k ≤ k ω ( · , k . Then by (2.5) and (2.7) Z T inf v h ∈ V h k h ( u h ·∇ ω h − v h ) k d t . T (cid:16)Z T s ( u h ; ω h , ω h ) d t (cid:17) . T k ω ( · , k . Using an elementwise trace inequality and (2.7) we also havemin( h, ν T ) k ν [[ n F ·∇ ω h ]] k F× I ≤ h k ν ∇ ω h k Q . h k ω ( · , k . For R we use the discrete Sobolev injection (4.4) to deduce h sup t ∈ I k ω h ( · , t ) k Z T k ω h ( · , t ) k L ∞ (Ω) d t . h k ω h ( · , k c ( h ) Z T ( k ω h k + s ( u h , ω h , ω h ) ) d t . h c ( h ) k ω h ( · , k . The only remaining term is the stabilization term, which is not innocent since we donot have an a priori bound on the factor k u h k L ∞ ( Q ) . Here we use (2.9) to deduce, forall t > q > k u h k L ∞ (Ω) Z T s ( u h ; ω h , ω h ) d t ≤ c q k ω h k L q (Ω) T (cid:16)Z T s ( u h ; ω h , ω h ) d t (cid:17) obust error estimates for the two dimensional Navier-Stokes equations k u h k L ∞ ( Q ) Z T s ( u h ; ω h , ω h ) d t . c q h − qq T sup t ∈ I k ω h ( · , t ) kk ω h ( · , k and the estimate follows taking ǫ = ( q − /q . Note that the constant c q explodes as q → L -norm for the velocities follows as before from the vorticityestimate using, with C h denoting the Cl´ement interpolant, k u − u h k = k∇ e Ψ k = ( ∇ e Ψ , ∇ (Ψ − C h Ψ)) − ( e ω , Ψ h − C h Ψ) ≤ k u − u h kk∇ (Ψ − C h Ψ) k − ( − ∆˜ e ω + ˜ e ω , ( C h Ψ − Ψ h )) ≤ k u − u h kk∇ (Ψ − C h Ψ) k + |k ˜ ω − ˜ ω h k| k C h Ψ − Ψ h k H (Ω) . We conclude by using the H -stability of the Cl´ement interpolant, a Poincar´e inequal-ity and finally by dividing both sides with k u − u h k . Proposition 4.3. (A priori error estimate using monotonicity) Assume that Re h > , that the mass ( · , · ) M is evaluated using nodal quadrature, that the form a ( u h ; ω h , v h ) + s ( u h ; ω h , v h ) has the DMP property as defined in [6] and that (2.5) - (2.6) are satisfied as well as the assumptions of Lemma 2.1 Then there holds |k (˜ ω − ˜ ω h )( T ) k| δ . e TτF (cid:18) hδ (cid:19) and k ( u − u h )( · , T ) k . inf v h ∈ W l − h k ( u − v h )( · , T ) k + e TτF h . Proof . The terms R − R are bounded as in the proof of Proposition 4.2. Sinceby assumption the spatial discretization of (2.3) has the DMP property and the mass-matrix is evaluated using nodal quadrature, we know from [5, 6] that k ω h k L ∞ ( Q ) = k ω h ( · , k L ∞ (Ω) . Hence by (2.9) k u h k L ∞ ( Q ) ≤ c ∞ k ω h ( · , k L ∞ (Ω) . We may then use these L ∞ -boundstogether with the stabilities of Lemma 2.1 to upper bound the remaining residualquantities of (4.2). Using (2.9) and (2.7) we immediately have h sup t ∈ I k ω h ( · , t ) k Z T k ω h ( · , t ) k L ∞ (Ω) d t . h T k ω h ( · , kk ω h ( · , k L ∞ (Ω) For the residual term resulting from the mass lumping we have using the stability(2.10) h Z T k ∂ t ∇ ω h k d t . T ( U + k u h k L ∞ ( Q ) ) (cid:16)Z T ( s ( u h ; ω h , ω h ) + k ν ∇ ω h k ) d t (cid:17) ≤ T U k ω ( · , k . E. BURMAN
The remaining contribution from the stabilization is bounded as before using themaximum principle and (2.7). The proof of the L -norm estimate on the velocities isidentical to that of Proposition 4.2Note that only the proof of Proposition 4.3 uses the assumption Re h >
5. Stabilized finite element methods.
The estimates of Theorem 4.1 holdsfor any finite element method on the form (2.3)-(2.4). Indeed by taking s ( · ; · , · ) ≡ k ω h k L ∞ ( Q ) . When theconsistent mass matrix is used one may proved that (4.4) holds giving onceagain a priori estimates, at the price of a logarithmic factor.2. high order stabilization, we propose to stabilize the jump of the streamlinederivative. This scheme does not yield a maximum principle, so the residualscan not be completely a priori bounded. The scheme has some interestingconservation properties for two-dimensional Navier-Stokes’ computations thatwe will point out. If a nonlinear stabilization term is added and mass-lumpingis used the solution may be made monotone and the a priori error estimateof Proposition 4.3 holds, this time with the possibility of higher order conver-gence in the smooth portion of the flow. Finally if the consistent mass matrixis used and stabilization is added also in the crosswind direction, an estimateof the type (4.4) can be shown to hold leading to a priori error bounds usingProposition (4.2). We consider first the stabi-lization method obtained by penalizing the jumps of the streamline derivative overelement faces. We use the exact mass matrix in (2.3) and the stabilizing operator s sd ( u h , ω h , v h ) := γ X F ∈F U − ( h F [[ u h ·∇ ω h ]] , [[ u h ·∇ v h ]]) F . (5.1)For this formulation the following stability estimates hold Lemma 5.1. sup t ∈ I k ω h ( · , t ) k + 2 k ν ∇ ω h k Q + 2 γU − k h F [[ u h ·∇ ω h ]] k F ≤ k ω h ( · , k (5.2) and if the consistent mass matrix is used, k u h ( · , T ) k + 2 k ν ω h k Q = k u h ( · , k (5.3) obust error estimates for the two dimensional Navier-Stokes equations Proof . the proof of (5.2) is an immediate consequence of (2.7) and the definition(5.1). The inequality (5.3) follows by observing that s sd ( u h , ω h , Ψ h ) = γ X F ∈F ( U − h F [[ u h ·∇ ω h ]] , [[rot Ψ h ·∇ Ψ h ]]) F = 0 . Observe that the method dissipates enstrophy but conserves energy exactly as thephysics of the problem suggests. Using known results on interpolation between dis-crete spaces it is also straightforward to show (see [8]),inf v h ∈ V h k h ( u h · ∇ ω h − v h ) k . s sd ( u h , ω h , ω h ) . Unfortunately this stabilization operator can not be shown to satisfy (4.4). For thiswe need the stabilization to act also in the crosswind direction. We therefore proposethe following two stabilization operators, the first is the standard artificial viscositymethod s av ( u h ; ω h , ω h ) := ( γh ( U + | u h | ) U − ∇ ω h , ∇ v h ) (5.4)and the second is a modification of (5.1) where also the crosswind gradient is penalizeddefined by s cd ( u h , ω h , v h ) := s sd ( u h , ω h , v h ) + γ X K ∈T h U h µK Z ∂K [[ n F · ∇ ω h ]][[ n F · ∇ v h ]] d s (5.5)Observe that the first part of s cd ensures the satisfaction of (2.5) and as we shall seethe second part is necessary for (4.4) to hold. Proposition 5.2.
Both stabilization operators (5.4) and (5.5) satisfy (2.5) and (2.6) . The stabilization operator s av ( · ; · , · ) satisfy (4.4) with c ( h ) ∼ h − (1 + | log ( h ) | ) and s cd ( · ; · , · ) satisfy (4.4) with c ( h ) ∼ h − µ (1 + | log ( h ) | ) , µ > .Proof . The proofs of (2.5) - (2.6) are consequences of the Cauchy-Schwarz in-equality and in the case of s cd trace inequalities. To prove (4.4) we note that in twospace dimensions there holds (see [23]), k ω h k L ∞ (Ω) . (1 + | log( h ) | ) k ω h k H (Ω) . This allows us to conclude for s av . For s cd we use that k∇ ω h k ≤ (cid:16) X F ∈F Z F | [[ ∇ ω h · n F ]] || ω h | d s (cid:17) . A Cauchy-Schwarz inequality followed by a trace inequality in the right hand sideleads to k∇ ω h k . (cid:16) X K ∈T h h − µ k ω h k K k h µ [[ n f ·∇ ω h ]] k ∂K (cid:17) . h − µ ( k ω h k + s cd ( u h ; ω h , ω h ) ) . Since the assumptions of Proposition 4.2 are satisfied, we may conclude that the6
E. BURMAN method (2.3)-(2.4) using the stabilization (5.4) statisfy the a priori error bounds for ǫ > |k (˜ ω − ˜ ω h )( T ) k| δ . e TτF (cid:18) hδ (cid:19) ( c ǫ h − ǫ + 1 + | log( h ) | )and k ( u − u h )( · , T ) k . inf v h ∈ W l − h k ( u − v h )( · , T ) k + e TτF (cid:18) hδ (cid:19) ( c ǫ h − ǫ + 1 + | log( h ) | ) . Similarly we have the following estimates if the stabilization (5.5) is used. |k (˜ ω − ˜ ω h )( T ) k| δ . e TτF (cid:18) hδ (cid:19) ( c ǫ h − ǫ + (1 + | log( h ) | ) h − µ − )and k ( u − u h )( · , T ) k . inf v h ∈ W l − h k ( u − v h )( · , T ) k + e TτF h ( c ǫ h − ǫ + (1 + | log( h ) | ) h − µ − ) . We see that if we take µ = 1 in (5.5) we get the same order for the two methods,however if we want the method to have optimal convergence for smooth solutions wechoose µ = 2 and l = 2, resulting in an a priori convergence order of O ( h ) in thenon-smooth case. Since the consistent mass matrix is non-monotonewe herein only consider methods using lumped mass. Monotone methods can also bedesigned using a nonlinear switch that changes the local quadrature as a function ofthe solution ω h so that the consistent mass is used away from local extrema to reducethe dispersion error known to haunt mass-lumping schemes, such methods are beyondthe scope of the present paper. A monotone method using linear artificialviscosity is obtained by taking (see [6]) s ( u h , ω h , v h ) := γ X K (max( U , k u h k L ∞ ( K ) ) h K X F ∈ ∂K ( ∇ ω h × n F , ∇ v h × n F ) F . (5.6)Then the estimates (2.7)-(2.8) hold and we observe that there exists positive constants c , c such that c k| u h | h ∇ ω h k ≤ s ( u h , ω h , ω h ) ≤ c k| u h | h ∇ ω h k . Let the mass matrix be evaluated using nodal quadrature so that the matrix cor-responding to ( · , · ) M is diagonal. We may use the theory of [5, 6] to prove thatthe operator a ( ω h , v h ) + s ( u h , ω h , v h ) has the DMP-property and hence the followingdiscrete maximum principle holds k ω h k L ∞ ( Q ) = k ω h ( · , k L ∞ (Ω) . This requires the parameter γ to be chosen large enough, however it does not requireany additional acute condition on the mesh, since the discretization of the Laplace op-erator results in an M-matrix on Delaunay meshes. Since by the maximum principle, k u h k L ∞ ( Q ) . k ω h ( · , k L ∞ (Ω) we have k| u h | h ∇ ω h k Q . k u h k L ∞ ( Q ) Z T s ( u h , ω h , ω h ) d t . k ω h ( · , k M (5.7) obust error estimates for the two dimensional Navier-Stokes equations v h = 0. It is straightforward to prove also (2.6). Comparingwith Proposition 4.3 we conclude that the assumptions are satisfied and hence thatthe Proposition holds for (2.3)-(2.4) with stabilization given by (5.6) and the massmatrix evaluated using nodal quadrature. Here we assume that l = 1 so that both ω h and Ψ h are discretized using piecewise affine elements. We propose a stabilizationterm consisting of one linear part and one nonlinear part. The role of the nonlinearpart is to ensure that the form a ( · ; · , · ) + s ( · ; · , · ) has the DMP property. The linearpart is necessary to ensure that the inequality (2.5) holds. We define s ( u h ; ω h , v h ) := s sd ( u h ; ω h , v h ) (5.8)+ γ X K h X F ∈ ∂K R F ( u h , ω h )(sign( ∇ ω h × n F ) , ∇ v h × n F ) F (5.9)where R F ( u h , ω h ) := k u h k L ∞ (∆ F ) (1 + U − k u h k L ∞ (∆ F ) ) m F ([[ n F · ∇ ω h ]])with ∆ F := ∪ K ∈T h ; K ∩ F = ∅ K and m F ([[ n F · ∇ ω h ]]) = max F ′ ∈F F ′ ∈ ∂K ′ ; K ′ ∩ F = F k [[ n F · ∇ ω h ]] k F ′ . It is shown in [6] that with this definition a ( u h ; ω h , v h ) + s ( u h ; ω h , v h ) has the DMP-property for γ large enough. Since the bounds (2.5)-(2.10) also hold, the assumptionsof Proposition 4.3 are satisfied and its estimates hold. We conclude that for themethods defined by mass lumping and the stabilization operators (5.6) or (5.8) thefollowing estimates hold |k (˜ ω − ˜ ω h )( T ) k| δ . e TτF (cid:18) hδ (cid:19) and k ( u − u h )( · , T ) k . inf v h ∈ W l − h k ( u − v h )( · , T ) k + e TτF (cid:18) hδ (cid:19) .
6. Conclusion.
We have shown that under a certain structural assumption onthe solution of the two dimensional Navier-Stokes’ equation one may derive robusterror estimates with an order in h , independent of both the Reynolds number andhigh order Sobolev norms of the exact solution. Robustness is obtained for a classof stabilized finite element methods. The estimates are both on a posteriori form,and on a priori form, providing an upper bound on the error. Due to the strongassumptions on the mesh the present a posteriori error estimates are not immediatelysuitable for use in adaptive algorithms, but a more detailed analysis may allow themesh assumptions to be relaxed. If the solution is smooth we also prove that optimalconvergence may be obtained, provided the stabilization operator is weakly consistentto the right order.Observe that it is natural that the LES estimate has much poorer convergenceorder, since we may assume no smoothness of the exact solution. Even the large scalesare assumed to have moderate gradients only.8 E. BURMAN
We show how several stabilized methods enter the framework, both first andsecond order accurate ones. The interest of the first order artificial viscosity methodis primarily its close relationship to the vertex centered finite volume method. Notealso that the estimates with an order proposed herein for nonlinear monotone schemesto the best of our knowledge are the first of their kind in the literature.It appears that for implicit large eddy simulations both the estimate (1.1) forsmooth solutions and the estimate (1.3) for rough solutions derived herein are desirableproperties for the theoretical justification of a method.Future work will focus on numerical investigations both in two and three spacedimensions. Of particular interest is to study the stability of the incompressible Eulerequations to see if the limit estimate with no allowed small scales is sharp.
Acknowledgment.
Partial funding for this research was provided by EPSRC(Award number EP/J002313/1).
Appendix A. Proof of Proposition 2.2.
We introduce the discrete errors,with I h denoting the Lagrange interpolant e h,ψ := Ψ h − I h Ψ and e h,ω := ω h − π L ω. First consider the second equation (2.4) and use Galerkin orthogonality k∇ e h, Ψ k = ( ∇ (Ψ − I h Ψ) , ∇ e h, Ψ ) − ( ω − ω h , e h, Ψ ) . Applying Poincar´es inequality followed by Cauchy Schwarz inequality we obtain thefollowing bound for Ψ in terms of the error in the vorticity k∇ e h, Ψ k . k∇ (Ψ − I h Ψ) k + k ω − π L ω h k + k e h,ω k . Consider now the equation (2.3) taking v h = e h,ω and observing that there holds12 ddt k e h,ω k + s ( u h ; e h,ω , e h,ω )= ( ∂ t ( ω − π L ω ) , e h,ω ) + ( ω, u · ∇ e h,ω ) − ( π L ω, u h · ∇ e h,ω ) − s ( u h ; π L ω, e h,ω ) . By integration by parts in time we see that the first term on the right hand side iszero, by the orthogonality of the L -projection. We then add and subtract u h in thesecond term on the right hand side to obtain12 ddt k e h,ω k + s ( u h ; e h,ω , e h,ω ) = ( ω, ( u − u h ) · ∇ e h,ω )+ ( ω − π L ω, u h · ∇ e h,ω ) − s ( u h ; π L ω, e h,ω ) = I + II + III.
In the first term on the right hand side we now reintegrate by parts and use Cauchy-Schwarz inequality, I ≤ k ω k W , ∞ k u − u h kk e h,ω k . k ω k W , ∞ ( k∇ (Ψ − I h Ψ) k + k ω − ω h k ) k e h,ω k . k ω k W , ∞ ( k∇ (Ψ − I h Ψ) k + k ω − π L ω k + k e h,ω k ) . obust error estimates for the two dimensional Navier-Stokes equations L -projection to retract somefunction v h and then apply (2.5), II = ( ω − π L ω, u h · ∇ e h,ω − v h ) ≤ c k h − ( ω − π L ω ) k s ( u h ; e h,ω , e h,ω ) ≤ ch s − k ω k H s + 14 s ( u h ; e h,ω , e h,ω ) . For the stabilization term finally we apply the Cauchy-Schwarz inequality and anarithmetic-geometric inequality to obtain
III = s ( u h ; π L ω, e h,ω ) ≤ s ( u h ; π L ω, π L ω ) + 14 s ( u h ; e h,ω , e h,ω ) . (A.1)Then we observe that by adding and subtracting I h rot Ψ we may write s ( u h ; π L ω, π L ω ) . s ( u h − I h rot Ψ; π L ω, π L ω ) + s ( I h rot Ψ; π L ω, π L ω )and using the definition (5.1) and the stability of the L -projection on quasi uniformmeshes, we have, s ( u h − I h rot Ψ; π L ω, π L ω ) . X F ∈F Z F h | u h − I h rot Ψ | |∇ π L ω h | d s . k∇ ω k L ∞ (Ω) h k u h − I h rot Ψ k . k∇ ω k L ∞ (Ω) h ( k∇ Ψ − I h ∇ Ψ) k + k∇ (Ψ − I h Ψ) k + k ω − π L ω k + k e h,ω k )and then s ( I h rot Ψ; π L ω, π L ω ) ≤ k I h rot Ψ k L ∞ ( Q ) X K h K (cid:0) k∇ ( ω − π L ω ) k K + h K k∇ ( ω − π L ω ) k (cid:1) ≤ k u k L ∞ ( Q ) Ch s − K k ω k L ( I ; H s (Ω)) . We conclude by collecting the upper bounds for the terms I − III , applying approx-imability and Gronwall’s lemma that k e h,ω ( · , T ) k + Z T s ( u h ; e h,ω , e h,ω ) . exp( cT k ω k W , ∞ ( Q ) ) k ω k W , ∞ ( Q ) × ( h l k Ψ k L ( I ; H l +1 (Ω) + h k +1 k ω k L ( I ; H (Ω)) ) . Here we assumed h k ω k W , ∞ ( Q ) . k u k L ∞ ( Q ) (thatis upper bounded by k ω ( · , k L ∞ (Ω) ). It follows that for l = 2 and sufficiently smoothsolutions we have k ( ω − ω h )( · , T ) k + k ( u − u h )( · , T ) k . h . REFERENCES[1] A. Aspden, N. Nikiforakis, S. Dalziel, and J. B. Bell. Analysis of implicit LES methods.
Commun. Appl. Math. Comput. Sci. , 3:103–126, 2008. E. BURMAN[2] J. P. Boris. On large eddy simulation using subgrid turbulence models comment 1. In J. L.Lumley, editor,
Whither Turbulence? Turbulence at the Crossroads , Lecture Notes inPhysics, page 344353. Berlin Springer Verlag, 1990.[3] E. Burman. A posteriori error estimation for interior penalty finite element approximations ofthe advection-reaction equation.
SIAM J. Numer. Anal. , 47(5):3584–3607, 2009.[4] E. Burman. Computability of filtered quantities for the Burgers’ equation. Technical ReportarXiv:1111.1182, 2012.[5] E. Burman and A. Ern. The discrete maximum principle for stabilized finite element methods.In Franco Brezzi, Annalisa Buffa, Stefania Corsaro, and Almerico Murli, editors,
NumericalMathematics and Advanced Applications , pages 557–566. Springer Milan, 2003.[6] E. Burman and A. Ern. Stabilized Galerkin approximation of convection-diffusion-reactionequations: discrete maximum principle and convergence.
Math. Comp. , 74(252):1637–1652(electronic), 2005.[7] E. Burman and M. A. Fern´andez. Continuous interior penalty finite element method for thetime-dependent Navier-Stokes equations: space discretization and convergence.
Numer.Math. , 107(1):39–77, 2007.[8] E. Burman and P. Hansbo. Edge stabilization for Galerkin approximations of convection-diffusion-reaction problems.
Comput. Methods Appl. Mech. Engrg. , 193(15-16):1437–1453,2004.[9] A. Dunca and V. John. Finite element error analysis of space averaged flow fields defined by adifferential filter.
Math. Models Methods Appl. Sci. , 14(4):603–618, 2004.[10] A. Dunca, V. John, and W. Layton. Approximating local averages of fluid velocities: Theequilibrium Navier-Stokes equations.
Appl. Numer. Math. , 49(2):187–205, 2004.[11] P. Hansbo and A. Szepessy. A velocity-pressure streamline diffusion finite element methodfor the incompressible Navier-Stokes equations.
Comput. Methods Appl. Mech. Engrg. ,84(2):175–192, 1990.[12] J. Hoffman. Computation of mean drag for bluff body problems using adaptive DNS/LES.
SIAM J. Sci. Comput. , 27(1):184–207 (electronic), 2005.[13] J. Hoffman and C. Johnson. Stability of the dual Navier-Stokes equations and efficient com-putation of mean output in turbulent flow using adaptive DNS/LES.
Comput. MethodsAppl. Mech. Engrg. , 195(13-16):1709–1721, 2006.[14] P. Houston, J. A. Mackenzie, E. S¨uli, and G. Warnecke. A posteriori error analysis for numericalapproximations of Friedrichs systems.
Numer. Math. , 82(3):433–470, 1999.[15] C. Johnson and J. Saranen. Streamline diffusion methods for the incompressible Euler andNavier-Stokes equations.
Math. Comp. , 47(175):1–18, 1986.[16] J. Kent, J. Thuburn, and N. Wood. Assessing implicit large eddy simulation for two-dimensionalflow.
Quarterly Journal of the Royal Meteorological Society , 138(663):365–376, 2012.[17] J.-G. Liu and C.-W. Shu. A high-order discontinuous Galerkin method for 2D incompressibleflows.
J. Comput. Phys. , 160(2):577–596, 2000.[18] L. G. Margolin, W. J. Rider, and F. F. Grinstein. Modeling turbulent flow with implicit LES.
J. Turbul. , 7:Paper 15, 27 pp. (electronic), 2006.[19] V. Maz’ya. On the boundedness of first derivatives for solutions to the Neumann-Laplaceproblem in a convex domain.
J. Math. Sci. (N. Y.) , 159(1):104–112, 2009. Problems inmathematical analysis. No. 40.[20] S. B. Pope. Ten questions concerning the large-eddy simulation of turbulent flows.
New Journalof Physics , 6(1):35, 2004.[21] J. Principe, R. Codina, and F. Henke. The dissipative structure of variational multiscalemethods for incompressible flows.
Comput. Methods Appl. Mech. Engrg. , 199(13-16):791–801, 2010.[22] R. Rannacher and R. Scott. Some optimal error estimates for piecewise linear finite elementapproximations.
Math. Comp. , 38(158):437–445, 1982.[23] R. Scott. Optimal L ∞ estimates for the finite element method on irregular meshes. Math.Comp. , 30(136):681–697, 1976.[24] V. Thom´ee.
Galerkin finite element methods for parabolic problems , volume 25 of