Abstract

We present enhanced physics-based finite element schemes for two families of turbulence models, the NS-𝜔 models and the Stolz-Adams approximate deconvolution models. These schemes are delicate extensions of a method created for the Navier-Stokes equations in Rebholz (2007), that achieve high physical fidelity by admitting balances of both energy and helicity that match the true physics. The schemes' development requires carefully chosen discrete curl, discrete Laplacian, and discrete filtering operators, in order to permit the necessary differential operator commutations.

1. Introduction

Energy and helicity play fundamental physical roles in the development of turbulent flow. They are conserved in inviscid flow, are delicately balanced when viscous and body forces are present [1, 2], are believed critical in flow structure development [3], and are cascaded jointly through the inertial range until they are removed by viscous forces [4]. Furthermore, their dissipation is closely tied to turbulence decay, and helicity has the topological interpretation that it is nonzero if and only if the velocity field is not reflectionally symmetric, that is, turbulent [3]. Therefore, the accurate treatment of these quantities is paramount for a computation to accurately predict fluid motion.

Unfortunately, often one or both of these quantities is not accurately treated in most models (in their continuous form), and thus cannot be correctly handled in their discretizations. Only a very few models correctly balance both energy and helicity, the most prominent being the Navier-Stokes equations (NSE, the true physical model), NS-𝛼 and related models [5, 6], the Stolz-Adams family of approximate deconvolution models (ADM) [711], and the recently introduced NS-𝜔 family of regularization models [12, 13]. Although these models yield, in a sense, the correct balances of energy and helicity, standard discretizations for these models do not preserve this important feature. For example, the usual Crank-Nicolson method for any of these schemes nonphysically creates helicity via the nonlinearity, thus altering the balance and removing physical meaning from solutions. Correctly accounting for both of these fundamental quantities is difficult; however it is very desirable since incorrect treatment of these quantities can lead to their nonphysical accumulations or degradations, and thus computed solutions that lack physical relevance.

A solution to this problem was first introduced in 2004 by Liu and Wang in [14], where they proposed and tested the EHPS, a finite difference scheme for axisymmetric flow (derived from the vorticity-stream formulation of the NSE) that more accurately treats both energy and helicity than standard schemes, and found excellent results. The idea was later developed for the full 3d NSE in [15], where the main ideas were to use the rotational form with a projected vorticity in a finite element scheme. Since direct computation of the NSE is limited to lower Reynolds numbers, this scheme was extended for use with NS-𝛼 type models [16], where again excellent numerical results were reported. Herein we show how, through carefully chosen definitions of the discrete filter, discrete Laplacian and discrete curl operators, the enhanced-physics based scheme of [15], can be extended for use on the NS-𝜔 and ADM families of models. The NS-𝜔 family of models, in addition to having a sound mathematical foundation [13], provides several advantages over the other classes of energy and helicity preserving models, including a more physical treatment of mass conservation and allowing for more efficient computations [12]. Therefore the development of an enhanced-physics based scheme for this attractive family of models is warranted. For completeness we also develop an enhanced-physics scheme for the ADM, so that each class of models that accurately balances energy and helicity has such a scheme developed for it.

2. The NS-𝜔 Family of Models

Given a domain Ω and endtime 𝑇, the NS-𝜔 model is given by

𝑣𝑡𝑣××𝑣]+𝑞𝜈Δ𝑣=𝑓inΩ×(0,𝑇,(2.1)𝑣=0in[]Ω×0,𝑇,(2.2)𝑣(0)=𝑣0inΩ,(2.3) where 𝑣 denotes velocity, 𝑣0 the initial velocity, 𝑞 pressure, 𝑓 a forcing, 𝜈 the kinematic viscosity, and overbar the Helmholtz filter with filtering radius 𝛼: 𝜙=𝐹𝜙=𝛼2Δ+𝐼1𝜙.(2.4) Extensions of this model to include van Cittert approximate deconvolution, to achieve higher-order consistency error to the NSE, are discussed in [12]. By defining the deconvolution operator 𝐷𝑁=𝑁𝑛=0(𝐼𝐹)𝑛(2.5) and replacing the filtered term 𝑣 in (2.1) by 𝐷𝑁𝑣, where 𝑁 is a positive integer typically chose 1𝑁7, we have the NS-𝜔 deconvolution models: 𝑣𝑡𝑣××𝐷𝑁𝑣],+𝑞𝜈Δ𝑣=𝑓inΩ×(0,𝑇𝑣=0in[],Ω×0,𝑇𝑣(0)=𝑣0inΩ.(2.6) Since the deconvolution operator acts as an approximate inverse to the filter [17], its use creates a more formally accurate class of models. Note also that the 𝑁=0 model is the NS-𝜔 model (2.1)-(2.2) since 𝐷0=𝐼.

NS-𝜔 must be coupled to initial and boundary conditions in order to be a predictor of turbulent flow. As our interest is the accurate balance of energy and helicity, we consider for simplicity a domain Ω to be a 3d rectangular box with periodic boundary conditions. These boundary conditions are not likely to be physically realized, however their use here isolates the fundamental numerical difficulties from the additional complications that arise from more physically relevant boundary conditions. The work herein can be extended to most other boundary conditions, but would need to be done on an individual basis since boundary condition themselves can (and usually do) contribute to the energy and helicity balances.

Recall the balances of energy and helicity for NS-𝜔 from [13]. Denote the 𝐿2 norm and inner product by and (,), respectively.

Lemma 2.1. Given initial condition 𝑣0𝐿2(Ω), endtime 𝑇, and forcing 𝑓𝐿2(0,𝑇;𝐿2(Ω)), the NS-𝜔 model admits the following exact balances for energy and helicity 1Energy2𝑣(𝑇)2+𝜈𝑇01×𝑣(𝑡)𝑑𝑡=2𝑣02+𝑇0(𝑓(𝑡),𝑣(𝑡))𝑑𝑡,(2.7)Helicity𝑣(𝑇),×𝑣(𝑇)+𝛼2×𝑣(𝑇),(×)2𝑣(𝑇)+2𝜈𝑇0×𝑣(𝑡),(×)2𝑣+(𝑡)(×)2𝑣(𝑡),(×)3𝑣=(𝑡)𝑑𝑡𝑣0,×𝑣0+𝛼2×𝑣0,(×)2𝑣0+2𝑇0𝑓(𝑡),×𝑣(𝑡)𝑑𝑡.(2.8)

Remark 2.2. While the energy balance of NS-𝜔 is the same as for the NSE, the helicity balance is the same as for the spatially filtered NSE. It is proven in [13] that up to a filtering radius dependent length scale, energy and helicity are cascaded jointly through the inertial range, just as in true fluid flow [4].

Remark 2.3. Balances for NS-𝜔 with deconvolution can be easily derived by multiplying NS-𝜔-deconvolution by 𝑣 and ×𝐷𝑁𝑣, for energy and helicity respectively, then integrating by parts. The energy balance will be identical to the NS-𝜔 case, while an analogous result for the helicity balance will be obtained where the 𝑣 terms in (2.8) are replaced with 𝐷𝑁1/2𝑣.

3. An Enhanced-Physics Based Finite Scheme for NS-𝜔

The main obstacles in extending the scheme of [15] to NS-𝜔 is the need for the discrete filter and curl operators to commute, and that a form of the vorticity is discretely-divergence free. These difficulties were not present in the extension to NS-𝛼 in [16], and hence this development is fundamentally more difficult. The extension is also quite delicate, but was achieved through the use of a new discrete implementation of the Helmholtz filter that is defined in terms of a discrete curl, which is discretely-divergence free via a suitably chosen projection in its definition. We present the scheme and show its energy and helicity balances after some necessary notation.

Denote the zero mean periodic subspaces of 𝐻1(Ω) and 𝐿2(Ω), respectively, by 𝐻1#(Ω) and 𝐿2#(Ω), let (𝑋,𝑄)(𝐻1#,𝐿2#) be inf-sup stable finite element spaces (e.g., Taylor-Hood elements [1820]) on a regular, conforming mesh of Ω. Let 𝑉 denote the discretely-divergence free subspace of 𝑋: 𝑉=𝑣𝑋,𝑣,𝑞=0𝑞𝑋.(3.1) We now define the discrete curl operator, which will be enforced to be discretely-divergence free by its construction.

Definition 3.1. Define the discrete curl operator to be the 𝐿2 projection of the curl operation into 𝑉: Given 𝜙𝐻1(Ω), curlcurl𝜙 is the unique element in 𝑉 satisfying curl𝜙,𝑣=×𝜙,𝑣𝑣𝑉.(3.2)
Note it is obvious that the linear operator curl is well defined, and also that it is self-adjoint since the usual curl operator is in this setting.
Since the Laplacian operator applied to a divergence free function 𝜒 can be written as
Δ𝜒=×(×𝜒)+(𝜒)=×(×𝜒),(3.3) it is natural to define the discrete Laplacian in terms of the discrete curl operator, as follows.

Definition 3.2. Define the discrete Laplacian Δ𝐻1(Ω)𝑉 by Δ=curlcurl.(3.4)
We define the discrete filter 𝐹 in terms of the discrete Laplacian:
𝐹𝜙=𝜙=𝛼2Δ+𝐼1𝜙.(3.5) Note that, by construction, the filter commutes with the discretely-divergence free projected curl operator and thus also the discrete Laplacian; recall these were the obstacles needed overcome for the development of an enhanced-physics based scheme for NS-𝜔, which we present now.

Algorithm 1 (an enhanced-physics based scheme for NS-𝜔). Given an endtime 𝑇, timestep Δ𝑡>0, filtering radius 𝛼>0, kinematic viscosity 𝜈>0, initial condition 𝑣0𝑉, forcing 𝑓𝐿2(0,𝑇;𝐿2(Ω)), set 𝑀=𝑇/Δ𝑡 and for 𝑛=1,2,,𝑀 find (𝑣𝑛,𝑝𝑛)(𝑋,𝑄) satisfying forall(𝜒,𝑞)(𝑋,𝑄), 1𝑣Δ𝑡𝑛+1𝑣𝑛,𝜒𝑣𝑛+1/2×curl𝑣𝑛+1/2,𝜒𝑝𝑛+1/2,𝜒Δ𝜈𝑣𝑛+1/2,𝜒=𝑓𝑡𝑛+1/2,𝜒,(3.6)𝑣𝑛+1,𝑞=0.(3.7)

Remark 3.3. That solutions exist for Algorithm 1 follows analogous to the case for NS-𝛼 given in [16]. It is straight forward, following similar to [16, 19] to prove optimal convergence to a smooth, weak NSE solution for Taylor-Hood elements.

We now prove the balances of energy and helicity by the scheme. Note the energy balance is, in a sense, a stability estimate.

Theorem 3.4. The scheme of Algorithm 1 for NS-𝜔 admits the following exact energy and helicity balances 1Energy2𝑣𝑀2+𝜈Δ𝑡𝑀1𝑛=0curl𝑣𝑛+1/22=12𝑣02+Δ𝑡𝑀1𝑛=0𝑓𝑡𝑛+1/2,𝑣𝑛+1/2,(3.8)Helicity𝑣𝑀,curl𝑣𝑀+𝛼2curl𝑣𝑀,curl2𝑣𝑀+2𝜈Δ𝑡𝑀1𝑛=0curl𝑣𝑛+1/2,curl2𝑣𝑛+1/2+𝛼2curl2𝑣𝑛+12,curl3𝑣𝑛+1/2=𝑣0,curl𝑣0+𝛼2curl𝑣0,curl2𝑣0+2Δ𝑡𝑀1𝑛=0𝑓𝑡𝑛+1/2,curl𝑣𝑛+1/2.(3.9)

Remark 3.5. The energy and helicity balances of the theorem are discrete analogs to the continuous NS-𝜔 model's balances from Lemma 2.1, and thus are also discrete analogs of the energy balance of the NSE and helicity balance of the spatially filtered NSE respectively.

Proof. The energy balance follows from choosing 𝜒=𝑣𝑛+1/2 in (3.6), 𝑞=𝑝𝑛+1/2 in (3.7), adding the equations then summing over timesteps.
For the helicity balance, begin by choosing 𝜒=curl𝑣𝑛+1/2. This will immediately vanish the nonlinear term since the cross of two vectors is perpendicular to each of them, and will eliminate the pressure term since the discrete curl is discretely-divergence free. This gives us the equation
1𝑣Δ𝑡𝑛+1𝑣𝑛,curl𝑣𝑛+1/2Δ𝜈𝑣𝑛+1/2,curl𝑣𝑛+1/2=𝑓𝑡𝑛+1/2,curl𝑣𝑛+1/2.(3.10) Since the discrete curl and discrete filter commute and are self-adjoint, the first term can be expanded and reduced via 1𝑣Δ𝑡𝑛+1𝑣𝑛,curl𝑣𝑛+1/2=12Δ𝑡𝑣𝑛+1,curl𝑣𝑛+1𝑣𝑛,curl𝑣𝑛+1+𝑣𝑛+1,curl𝑣𝑛𝑣𝑛,curl𝑣𝑛=12Δ𝑡𝑣𝑛+1,curl𝑣𝑛+1𝑣𝑛,curl𝑣𝑛=12Δ𝑡𝑣𝑛+1,curl𝑣𝑛+1+𝛼2curl𝑣𝑛+1,curl2𝑣𝑛+1𝑣𝑛,curl𝑣𝑛𝛼2curl𝑣𝑛,curl2𝑣𝑛.(3.11) Similarly for the viscous term, we have that Δ𝜈𝑣𝑛+1/2,curl𝑣𝑛+1/2=𝜈curl2𝑣𝑛+1/2,curl𝑣𝑛+1/2=𝜈curl2𝑣𝑛+1/2,curl𝑣𝑛+1/2𝛼2𝜈curl2Δ𝑣𝑛+1/2,curl𝑣𝑛+1/2=𝜈curl2𝑣𝑛+1/2,curl𝑣𝑛+1/2+𝛼2𝜈curl3𝑣𝑛+1/2,curl2𝑣𝑛+1/2.(3.12) Summing over timesteps complete the proof.

Remark 3.6. From the proof, the need for the intricate design of the enhanced-physics scheme becomes more clear. Without the use of the discrete curl in the scheme, a test function could not have been chosen to vanish the nonlinearity in the helicity balance, which in turn would cause the nonphysical creation and dissipation of helicity instead of cascading it from the large to small scales. The derivation of differences in helicities at successive timesteps in (3.11) was made possible because the discrete curl and discrete filter commute. Commutation of these operators with the discrete Laplacian was necessary in the analysis of the viscous term.

3.1. Extension to the NS-𝜔 Family of Deconvolution Models

The development of the discrete curl, Laplacian and filter operators allows for an easy extension to the deconvolution family of models (2.4)–(2.6) provided the natural definition for the discrete deconvolution operator is employed: 𝐷𝑁=𝑁𝑛=0𝐼𝐹𝑛.(3.13) Thus the discrete deconvolution can be seen to be positive and a polynomial in the discrete filter, and as such that will be self-adjoint and commute with the discrete curl, discrete Laplacian, and discrete filter. This leads us to the following enhanced-physics based algorithm for NS-𝜔 family of deconvolution models.

Algorithm 2 (an enhanced-physics based scheme for NS-𝜔-deconvolution). Given an endtime 𝑇, timestep Δ𝑡>0, filtering radius 𝛼>0, deconvolution order 𝑁0, kinematic viscosity 𝜈>0, initial condition 𝑣0𝑉, forcing 𝑓𝐿2(0,𝑇;𝐿2(Ω)), set 𝑀=𝑇/Δ𝑡 and for 𝑛=1,2,,𝑀 find (𝑣𝑛,𝑝𝑛)(𝑋,𝑄) satisfying forall(𝜒,𝑞)(𝑋,𝑄), 1𝑣Δ𝑡𝑛+1𝑣𝑛,𝜒𝑣𝑛+1/2×curl𝐷𝑁𝑣𝑛+1/2,𝜒𝑝𝑛+1/2,𝜒Δ𝜈𝑣𝑛+1/2,𝜒=𝑓𝑡𝑛+1/2,𝜒,(3.14)𝑣𝑛+1,𝑞=0.(3.15)

Remark 3.7. As in the case of no deconvolution, proofs of solution existence and convergence follow similarly to work in [16, 21].

We now provide the energy and helicity balances for this algorithm.

Theorem 3.8. The scheme of Algorithm 2 for NS-𝜔-deconvolution admits the following exact energy and helicity balances 1Energy2𝑣𝑀2+𝜈Δ𝑡𝑀1𝑛=0curl𝑣𝑛+1/22=12𝑣02+Δ𝑡𝑀1𝑛=0𝑓𝑡𝑛+1/2,𝑣𝑛+1/2.(3.16)𝐷Helicity𝑁1/2𝑣𝑀,curl𝐷𝑁1/2𝑣𝑀+𝛼2curl𝐷𝑁1/2𝑣𝑀,curl2𝐷𝑁1/2𝑣𝑀+2𝜈Δ𝑡𝑀1𝑛=0curl𝐷𝑁1/2𝑣𝑛+1/2,curl2𝐷𝑁1/2𝑣𝑛+1/2+2𝜈𝛼2Δ𝑡𝑀1𝑛=0curl2𝐷𝑁1/2𝑣𝑛+1/2,curl3𝐷𝑁1/2𝑣𝑛+1/2=𝐷𝑁1/2𝑣0,curl𝐷𝑁1/2𝑣0+𝛼2curl𝐷𝑁1/2𝑣0,curl2𝐷𝑁1/2𝑣0+2Δ𝑡𝑀1𝑛=0𝑓𝑡𝑛+1/2,curl𝐷𝑁1/2𝑣𝑛+1/2.(3.17)

Proof. The energy balance follows exactly as in the case without deconvolution. The helicity balance follows similar to the case of NS-𝜔; choose 𝜒=𝑐𝑢𝑟𝑙𝐷𝑁 to vanish the nonlinearity and pressure, then use commutativity of the differential operators and that the deconvolution operator is positive and self-adjoint to obtain the result.

4. Extension to the Stolz-Adams ADM

The Stolz-Adams ADM family take a similar mathematical form to the NS-𝜔 models, and the tools used to derive an enhanced-physics based scheme are sufficient to devise such a scheme for the ADM. For simplicity, we consider the ADM secondary regularization, often called time relaxation [8, 22]. This term does not correctly balance energy or helicity, however it is negligable away from the smallest resolved scales [22], and so it has little total effect on the overall balances.The ADM can be written in rotational form as

𝑣𝑡𝐷𝑁𝑣𝐷𝑁𝑣+𝑞𝜈Δ𝑣=𝑓,(4.1)𝐷𝑁𝑣=0,(4.2)𝑣(0)=𝑣0,(4.3) where filtering and deconvolution are defined exactly as for NS-𝜔. Its balances for energy and helicity are derived in [10], and follow from multiplying (4.1) by (𝛼2Δ+𝐼)𝐷𝑁𝑣 and (𝛼2Δ+𝐼)(×𝐷𝑁𝑣), respectively, then integrating over the domain, integrating by parts, and reducing. These balances are as follows

1Energy2𝐷𝑁1/2𝑣(𝑇)2+𝛼2×𝐷𝑁1/2𝑣(𝑇)2+𝜈𝑇0×𝐷𝑁1/2𝑣(𝑡)2+𝛼2Δ𝐷𝑁1/2𝑣(𝑡)2=1𝑑𝑡2𝐷𝑁1/2𝑣02+𝛼2×𝐷𝑁1/2𝑣02+𝑇0𝑓(𝑡),𝐷𝑁𝑣(𝑡)+𝛼2×𝑓(𝑡),×𝐷𝑁𝐷𝑣(𝑡)𝑑𝑡,(4.4)Helicity𝑁1/2𝑣(𝑇),×𝐷𝑁1/2𝑣(𝑇)+𝛼2×𝐷𝑁1/2𝑣(𝑇),(×)2𝐷𝑁1/2𝑣(𝑇)+2𝜈𝑇0×𝐷𝑁1/2𝑣(𝑡),(×)2𝐷𝑁1/2𝑣(𝑡)+𝛼2(×)2𝐷𝑁1/2𝑣(𝑡),(×)3𝐷𝑁1/2=𝐷𝑣(𝑡)𝑑𝑡𝑁1/2𝑣0,×𝐷𝑁1/2𝑣0+𝛼2×𝐷𝑁1/2𝑣0,(×)2𝐷𝑁1/2𝑣0+2𝑇0𝑓,×𝐷𝑁𝑣(𝑡)+𝛼2×𝑓(𝑡),(×)2𝐷𝑁𝑣(𝑡)𝑑𝑡.(4.5)

Remark 4.1. Technical arguments in [11] discuss how these balances are equivalent, in a sense, to the energy and helicity balances of the spatially filtered NSE. This is not a coincidence, since the ADM can be derived as an approximation of the spatially filtered NSE [17].

Though perhaps not obvious, the discrete operators which permitted the development of an enhanced-physics scheme for the NS-𝜔 family above are sufficient to devise such a scheme for the ADM. We now present this numerical scheme, which provides discrete analogs for the balances (4.4) and (4.5).

Algorithm 3 (an enhanced-physics based scheme for the ADM). Given an endtime 𝑇, timestep Δ𝑡>0, filtering radius 𝛼>0, deconvolution order 𝑁0, kinematic viscosity 𝜈>0, initial condition 𝑣0𝑉, forcing 𝑓𝐿2(0,𝑇;𝐿2(Ω)), set 𝑀=𝑇/Δ𝑡 and for 𝑛=1,2,,𝑀 find (𝑣𝑛,𝑞𝑛)(𝑋,𝑄) satisfying forall(𝜒,𝑟)(𝑋,𝑄), 1𝑣Δ𝑡𝑛+1𝑣𝑛,𝜒𝐷𝑁𝑣𝑛+1/2×curl𝐷𝑁𝑣𝑛+1/2,𝜒𝑞𝑛+1/2,𝜒Δ𝜈𝑣𝑛+1/2,𝜒=𝑓𝑡𝑛+1/2,𝜒,(4.6)𝑣𝑛+1,𝑟=0.(4.7)

Remark 4.2. Since by definition the discrete curl is discretely-divergence free, the discrete Laplacian must be also. Hence filtered elements in 𝑉 remain in 𝑉, as do deconvolved ones. Thus, (4.7) is sufficient as a discrete enforcement of (4.2).

Remark 4.3. Similar to the NS-𝜔 algorithms, it is straight forward and follows similar to work in [16, 23] that this algorithm admits stable solutions that convergence to weak NSE solutions or spatially filtered NSE solutions.

Theorem 4.4. The scheme of Algorithm 3 for the ADM admits the following exact energy and helicity balances, which are discrete analogs the energy and helicity balances of the continuous ADM 1Energy2𝐷𝑁1/2𝑣𝑀2+𝛼2curl𝐷𝑁1/2𝑣𝑀2+Δ𝑡𝜈𝑀1𝑛=0curl𝐷𝑁1/2𝑣𝑛+1/22+𝛼2Δ𝐷𝑁1/2𝑣𝑛+1/22=12𝐷𝑁1/2𝑣02+𝛼2curl𝐷𝑁1/2𝑣02+Δ𝑡𝑀1𝑛=0𝑓𝑡𝑛+1/2,𝐷𝑁𝑣𝑛+1/2+𝛼2curl𝑓𝑡𝑛+1/2,curl𝐷𝑁𝑣𝑛+1/2,𝐷(4.8)Helicity𝑁1/2𝑣𝑀,curl𝐷𝑁1/2𝑣𝑀+𝛼2curl𝐷𝑁1/2𝑣𝑀,curl2𝐷𝑁1/2𝑣𝑀+2𝜈Δ𝑡𝑀1𝑛=0curl𝐷𝑁1/2𝑣𝑛+1/2,curl2𝐷𝑁1/2𝑣𝑛+1/2+2𝜈𝛼2𝑀1𝑛=0curl2𝐷𝑁1/2𝑣𝑛+1/2,curl3𝐷𝑁1/2𝑣𝑛+1/2=𝐷𝑁1/2𝑣0,curl𝐷𝑁1/2𝑣0+𝛼2curl𝐷𝑁1/2𝑣0,curl2𝐷𝑁1/2𝑣0+2Δ𝑡𝑀1𝑛=0𝑓𝑡𝑛+1/2,curl𝐷𝑁𝑣𝑛+1/2+𝛼2curl𝑓𝑡𝑛+1/2,curl2𝐷𝑁𝑣𝑛+1/2.(4.9)

Proof. For the energy balance, choose 𝜒=𝐹1𝐷𝑁𝑣𝑛+1/2. This vanishes the nonlinear term, since after the 𝐹 and 𝐹1 cancellation of, we have that the cross of two vectors is perpendicular to each of them. The pressure term is also zero, since 𝐷𝑁𝑣𝑛+1/2 is discretely-divergence free. This leaves 1𝑣Δ𝑡𝑛+1𝑣𝑛,𝐹1𝐷𝑁𝑣𝑛+1/2Δ𝜈𝑣𝑛+1/2,𝐹1𝐷𝑁𝑣𝑛+1/2=𝑓𝑡𝑛+1/2,𝐹1𝐷𝑁𝑣𝑛+1/2.(4.10) Now using the definition of 𝐹1 we get 1𝑣Δ𝑡𝑛+1𝑣𝑛,𝐷𝑁𝑣𝑛+1/2𝛼2𝑣Δ𝑡𝑛+1𝑣𝑛,Δ𝐷𝑁𝑣𝑛+1/2Δ𝜈𝑣𝑛+1/2,𝐷𝑁𝑣𝑛+1/2+𝛼2𝜈Δ𝑣𝑛+1/2,Δ𝐷𝑁𝑣𝑛+1/2=𝑓𝑡𝑛+1/2,𝐷𝑁𝑣𝑛+1/2𝛼2𝑓𝑡𝑛+1/2,Δ𝐷𝑁𝑣𝑛+1/2,(4.11) and with the definition of Δ coupled with commutativity of operators, 1𝐷2Δ𝑡𝑁1/2𝑣𝑛+12+𝛼2curl𝐷𝑁1/2𝑣𝑛+12𝐷𝑁1/2𝑣𝑛2𝛼2curl𝐷𝑁1/2𝑣𝑛2+𝜈curl𝐷𝑁1/2𝑣𝑛+1/22+𝛼2𝜈Δ𝐷𝑁1/2𝑣𝑛+1/22=𝑓𝑡𝑛+1/2,𝐷𝑁𝑣𝑛+1/2+𝛼2curl𝑓𝑡𝑛+1/2,curl𝐷𝑁𝑣𝑛+1/2.(4.12) Multiplying through by Δ𝑡 and summing over timesteps now gives the stated energy balance.
For the helicity balance, choose 𝜒=𝐹1curl𝐷𝑁𝑣𝑛+1/2. This vanishes the nonlinear term since after the filter and its inverse cancellation of, the cross of two vectors is perpendicular to each of them. The pressure term also vanishes, since the curl is discretely-divergence free. This leaves
1𝑣Δ𝑡𝑛+1𝑣𝑛,𝐹1curl𝐷𝑁𝑣𝑛+1/2Δ𝜈𝑣𝑛+1/2,𝐹1curl𝐷𝑁𝑣𝑛+1/2=𝑓𝑡𝑛+1/2,𝐹1curl𝐷𝑁𝑣𝑛+1/2.(4.13) Next we decompose 𝐹1, similar to the case of the energy balance. After simplifying, we have 1𝐷2Δ𝑡𝑁1/2𝑣𝑛+1,curl𝐷𝑁1/2𝑣𝑛+1+𝛼2curl𝐷𝑁1/2𝑣𝑛+1,curl2𝐷𝑁1/2𝑣𝑛+11𝐷2Δ𝑡𝑁1/2𝑣𝑛,curl𝐷𝑁1/2𝑣𝑛+𝛼2curl𝐷𝑁1/2𝑣𝑛,curl2𝐷𝑁1/2𝑣𝑛+𝜈curl𝐷𝑁1/2𝑣𝑛+1/2,curl2𝐷𝑁1/2𝑣𝑛+1/2+𝛼2𝜈curl2𝐷𝑁1/2𝑣𝑛+1/2,curl3𝐷𝑁1/2𝑣𝑛+1/2=𝑓𝑡𝑛+1/2,curl𝐷𝑁𝑣𝑛+1/2+𝛼2curl𝑓𝑡𝑛+1/2,curl2𝐷𝑁𝑣𝑛+1/2.(4.14) We complete the proof by multiplying through by 2Δ𝑡 and summing over timesteps.

5. Conclusions

We have developed enhanced-physics based finite element schemes for the NS-𝜔 and Stolz-Adams ADM families of turbulence models that provides physically accurate treatment of discrete energy and helicity. Development of a discretely-divergence-free discrete curl operator, a discrete Laplacian defined in terms of the discrete curl, a discrete filter defined in terms of the discrete Laplacian, and a discrete deconvolution operator defined in terms of the discrete filter, all of which commute, are the key technical tools which allow the discrete schemes to accurately balance two fundamental physical quantities. Future directions of this work will be extensions to commonly used boundary conditions, and efficient implementations of these schemes.