Convergence rates for the numerical approximation of the 2D stochastic Navier-Stokes equations
Dominic Breit Alan Dodgson


Convergence rates for stochastic Navier–Stokes equations]Convergence rates for the numerical approximation of the 2D stochastic Navier–Stokes equations




Department of Mathematics, Heriot-Watt University, Riccarton Edinburgh EH14 4AS, UK
d.breit@hw.ac.uk, ad335@hw.ac.uk










We study stochastic Navier–Stokes equations in two dimensions with respect to periodic boundary conditions. The equations are perturbed by a nonlinear multiplicative stochastic forcing with linear growth (in the velocity) driven by a cylindrical Wiener process.
 We establish convergence rates for a finite-element based space-time approximation with respect to convergence in probability (where the error is measure in the L^∞_tL^2_x∩ L^2_tW^1,2_x-norm). Our main result provides linear convergence in space
and convergence of order (almost) 1/2 in time. This improves earlier results from [E. Carelli, A. Prohl: Rates of convergence for discretizations of the stochastic incompressible Navier–Stokes equations. SIAM J. Numer. Anal. 50(5), 2467–2496. (2012)] where the convergence rate in time is only (almost) 1/4.
Our approach is based on a careful analysis of the pressure function using a stochastic pressure decomposition.





[2010]65M15, 65C30, 60H15, 60H35




[
    Dominic Breit & Alan Dodgson
    December 30, 2023
================================






































§ INTRODUCTION


In this paper we are concerned with the stochastic Navier–Stokes equations

    {[ =μΔ-(∇)-∇π+Φ() W                 ;               =0                 ;           (0)=_0                  ].

on a filtered probability space (Ω,,(_t)_t≥0,). The equations are perturbed by an (_t)-Wiener process (possibly infinite dimensional) and Φ grows linearly in  (see Section 2.1 for the precise assumptions). The quantity μ>0 is the viscosity of the fluid and _0 is a given (random initial datum).
Here the unknowns are the velocity field :Ω× Q→^N and the pressure π:Ω× Q  →, where Q=(0,T)×𝒪 and 𝒪⊂^N with N=2,3.

The stochastic perturbation in the balance of momentum (<ref>)_1 can take into account for physical, empirical or numerical uncertainties and thermodynamical fluctuations. In addition to that, a main reason why the stochastic Navier–Stokes equations (<ref>) became so popular in fluid mechanical research is their application to turbulence theory (see, for instance, <cit.> and <cit.>). Its mathematical investigation started in the 70's with the pioneering paper of Bensoussan and Temam <cit.>. They provide a semi-deterministic approach based on the flow-transformation. A first fully stochastic theory has been developed by Flandoli and Gatarek <cit.> by showing the existence of a martingale solution. These solutions are weak in the stochastic sense meaning that
the underlying probability space as well the Wiener process W are not a priori known but become an integral part of the solution. In two dimensions, when uniqueness is known, a stochastically strong solution exists (it is defined on a given probability space with a given Wiener process), see <cit.>.
Nowadays there is a huge amount of literature concerning the analysis of (<ref>) and most of the results from the deterministic theory found their stochastic counterpart. For an overview we refer to the recent survey article
<cit.>.

The situation about the numerical approximation of (<ref>) is totally different and only very few results are available.
In <cit.> a fully practical space-time approximation for the three-dimensional stochastic Navier–Stokes equations in a bounded domain is studied. It is shown that the sequence of approximate solutions converges in law (up to a subsequence) to a martingale solution if both discretization parameters tend to zero. This is the best result one can hope for without using some unproven hypotheses about the space-regularity of solutions
(or to be content with local-in-time results). In two dimensions the situation is much better, at least if periodic boundary conditions are considered, that is

    𝒪 = 𝕋^2 = ( (-π, π)|_{ - π, π})^2.

The space-regularity of the unique strong solution is well-known (see for instance <cit.>). Based on this the convergence rates for a finite-element based space-time approximation
is analysed in <cit.>. The result is linear convergence in space
and convergence of order (almost) 1/4 in time. The precise estimate comparing
the solution  and its space-time approximation _h,m reads as

    [1_Ω_Δ t,h(max_1≤ m≤ M(t_m)-_h,m_L^2_x^2   +∑_m=1^M Δ t ∇(t_m)-∇_h,m_L^2_x^2)]
       ≤ c (h^2+(Δ t)^2α)

for any α<1/4. Here we have Ω_Δ t,h⊂Ω with (Ω∖Ω_Δ t,h)→0 as Δ t,h→0. Consequently, one can infer
from (<ref>) convergence in probability (as defined by Printems <cit.>) with asymptotic rates
(almost) 1/4 and 1 respectively. The aim of the present paper is to improve
the convergence rates in time from (almost) 1/4 to (almost) 1/2, see Theorem <ref> for the precise statement.
This is certainly the optimal convergence rate in time in view of the stochastic forcing and should also be considered as the natural result because of the space-regularity of the solution. 


The reason for the low convergence rate in time in <cit.>
is the low time-regularity of the pressure function π in (<ref>). Its appearance can only be avoided when working with finite element functions which are exactly divergence-free. Unfortunately, their construction is quite complicated such that the preferred descretizations (such as the Taylor–Hood, the Crouzeix–Raviart, and the MINI element, see <cit.>) are only 
asymptotically divergence-free. The low regularity of the pressure gradient
arises from the stochastic forcing (in fact, ∇π behaves as W) and it does not seem possible to improve unless the noise is divergence-free (this is quite restrictive as it only allows certain additive or linear multiplicative noise). In the general case (non-divergence-free finite elements and nonlinear multiplicative noise) a more subtle analysis of the pressure function
is required. Our main idea is a decomposition of the pressure into a deterministic and a stochastic component
(such a decomposition first appeared in <cit.>). The deterministic pressure part
behaves as the convective term ⊗. The latter one can be estimated along the lines of <cit.>
(following classical deterministic arguments combined with a discrete stopping time).
In addition a second stochastic integral appears which behaves similarly to the original stochastic integral. Although the time-regularity of this part has not improved, we benefit from the averaging properties of the Itô-integral. Combining these ideas finally leads to the optimal convergence rate in Theorem <ref>.

The paper is organized as follows. In Section 2 we present the mathematical framework, that is the probability setup, the concept of solutions and their qualitative properties. In particular, we give improved (compared to <cit.>) results on the time-regularity of ∇, see Corollary <ref> b). This is needed in Section 3 in order to estimate the error between the continuous solution and the time-discrete solution. The heart of the paper is Section 4 were we estimate the error between the time-discrete solution and the space-time discretization.
Crucial tools are the space-regularity of the time-discrete solution from <cit.>
and the decomposition of the corresponding pressure function.






§ MATHEMATICAL FRAMEWORK






 §.§ Probability setup




Let (Ω,,(_t)_t≥0,) be a stochastic basis with a complete, right-continuous filtration. The process W is a cylindrical Wiener process, that is, W(t)=∑_k≥1β_k(t) e_k with (β_k)_k≥1 being mutually independent real-valued standard Wiener processes relative to (_t) and (e_k)_k≥1 is a complete orthonormal system in a sepa­rable Hilbert space 𝔘.
To give the precise definition of the diffusion coefficient , consider ∈ L^2() and let ():𝔘→ L^2() be defined by Φ()e_k=_k(·,(·)). In particular, we suppose
that _k∈ C^1(×^2) and that the following conditions hold

    ∑_k≥1|_k(x,)|^2 ≤ c(1+||^2)   ,    ∑_k≥1|∇__k(x,)|^2≤ c,
    ∑_k≥1|∇_x _k(x,)|^2    ≤ c(1+||^2),   x∈, ∈^2.

If we are interested in higher regularity some further assumptions are in place and we require additionally _k∈ C^2(×^2) together with

    ∑_k≥1|∇_x^2_k(x,)|^2 ≤ c(1+||^2)   ,    ∑_k≥1|∇^2__k(x,)|^2≤c/1+||^2,
    ∑_k≥1|∇_x∇__k(x,)|^2    ≤ c,   x∈, ∈^2.

We remark that the first inequality of (<ref>) implies

    Φ()_L_2(𝔘;L^2_x)≤ c(1+_L^2_x)  ∀∈ L^2(),

all inequalities from (<ref>) imply

    Φ()_L_2(𝔘;W^1,2_x)≤ c(1+_W^1,2_x)  ∀∈ W^1,2(),

and
(<ref>) finally yields

    Φ()_L_2(𝔘;W^2,2_x)≤ c(1+_W^2,2_x)  ∀∈ W^2,2().

In fact, all our results apply if we replace (<ref>) and (<ref>) by the corresponding norm estimates above.

Furthermore, the conditions imposed on , particularly the first assumption from (<ref>), allow us to define stochastic integrals


Given an (𝔉_t)-progressively measurable process ∈ L^2(Ω;L^2(0,T;L^2())), the stochastic integral 
    t↦∫_0^t()  W

is a well defined process taking values in L^2() (see <cit.> for the detailed construction). Moreover, we can multiply by test functions to obtain
 
    ⟨∫_0^t ()  W,⟩=∑_k≥ 1∫_0^t⟨_k(),⟩ β_k,  ∈ L^2().

Similarly, we can define stochastic integrals with values in W^1,2() and W^2,2() respectively if  belongs to the corresponding class.

The following lemma is a helpful tool to analysis the time-regularity
of stochastic integrals (see, e.g., <cit.> or <cit.>).


 Let ψ∈ L^r(Ω;L^r(0,T;L_2(𝔘,L^2()))), r>2, by an (𝔉_t)-progressively measureable process and W a cylindrical (_t)-Wiener process on 𝔘. Then the paths of the process Z_t:=∫_0^t ψ  W are -a.s. Hölder continuous with exponent α∈(1/r,1/2) and it holds

    [Z_C^α([0,T];L^2())^r]≤ c_α [sup_0≤ t≤ Tψ_L_2(𝔘,L^2())^2]^r/2.









 §.§ The concept of solutions


In dimension two, pathwise uniqueness for weak solutions is known under (<ref>), we refer the reader for instance to Capiński–Cutland <cit.> and Capiński <cit.>. Consequently, we may work with the definition of a weak pathwise solution.



Let (Ω,,(_t)_t≥0,) be a given stochastic basis with a complete right-continuous filtration and an (_t)-cylindrical Wiener process W. Let _0 be an _0-measurable random variable. Then  is called a weak pathwise solution incompressible Navier–Stokes system!weak pathwise solution to (<ref>) with the initial condition _0 provided


  * the velocity field  is (_t)-adapted and

    ∈ C_w([0,T];L^2_())∩ L^2(0,T; W^1,2_())  -a.s.,









  * the momentum equation

    ∫_(t)·-∫__0·
       =∫_0^t∫_⊗:∇  t-μ∫_0^t∫_∇:∇  s
           +∫_0^t∫_Φ()·  W.

holds -a.s. for all ∈ C^∞_() and all t∈[0,T].




Let N=2 and assume that Φ satisfies (<ref>). Let (Ω,,(_t)_t≥0,) be a given stochastic basis with a complete right-continuous filtration and an (_t)-cylindrical Wiener process W. Let _0 be an _0-measurable random variable such that _0∈ L^r(Ω;L^2_div(𝕋^2)) for some r>2. Then there exists a unique weak pathwise solution to (<ref>) in the sense of Definition <ref> with the initial condition _0.


 Now, for ∈ C^∞() we can insert -∇Δ^-1
and obtain

    ∫_(t)·   +∫_0^t∫_μ∇:∇-∫_0^t∫_⊗:∇
       =∫_(0)·
    +∫_0^t∫_π_det 
       +∫_∫_0^tΦ()  W·+∫_∫_0^tΦ^π  W·,

where 

    π_det   =-Δ^-1(⊗),
    Φ^π   =-∇Δ^-1Φ().

This corresponds to the stochastic pressure decomposition introduced in <cit.> (see also <cit.> for a slightly different presentation). However, the situation with periodic boundary conditions we are considering here is much easier as the harmonic component of the pressure disappears. 











 §.§ Regularity of solutions


Let the assumptions of Theorem <ref> be satisfied. 
 

  (a) We have

    [sup_0≤ t≤ T∫_||^2+∫_0^T∫_|∇|^2]^r/2≤ c_r [1+_0_L^2_x^2]^r/2.


  (b) Assume that _0∈ L^r(Ω,W^1,2_div()) for some r≥2. Then we have

    [sup_0≤ t≤ T∫_|∇|^2+∫_0^T∫_|∇^2|^2]^r/2≤ c_r [1+_0_W^1,2_x^2]^r/2.


  (c) Assume that _0∈ L^r(Ω,W^2,2_div())∩ L^5r(Ω,W^1,2_div()) for some r≥2 and that (<ref>) holds. Then we have

    [sup_0≤ t≤ T∫_|∇^2|^2+∫_0^T∫_|∇^3|^2]^r/2≤ c_r [1+_0_W^2,2_x^2+_0_W^1,2_x^10]^r/2.




Part (a) is the standard a priori estimate which follows from applying Itô's formula to the functional f()=1/2_L^2_x^2 (and using Burkholder-Davis-Gundy inequality, assumption (<ref>) and Gronwall's lemma). Note that this is legit in two dimensions since we have ⊗∈ L^2(Q) -a.s. by Ladyshenskaya's inequality.

The proof of (b) and (c) is quite similar to <cit.>. However, we are working with a different setup. So, we decided to give a formal proof although it is certainly known to experts. The proof can be made rigorous by working with a Galerkin-type approximation and show that the following estimates are uniform with respect to the dimension of the ansatz space. Such a procedure is quite standard, so we leave the details to the reader.

In order to show (b) we apply Itô's formula to the function f_γ():=12∂_γ_L^2_x^2 (with γ∈{1,2}) and obtain

    1/2∂_γ(t)_L^2_x^2   =1/2∂_γ_0_L^2_x^2+∫_0^t f'_γ() +1/2∫_0^t f_γ”() ⟨⟨⟩⟩
       =1/2∂_γ_0_L^2_x^2+∫_∫_0^t ∂_γ·∂_γ
       +1/2∫_∫_0^t⟨⟨∫_0^·∂_γ(Φ()  W)⟩⟩
    =:(I)+(II)+(III).

We take the supremum in time, the r/2th power and apply expectations.
Summing over γ, we find  

    (II)   =-(II)_1-(II)_2+(II)_3,
    
    (II)_1   :=μ∫_0^t∫_ |∇^2|^2,
    
    (II)_2   :=∫_0^t∫_∂_γ·∂_γ(Φ()  W),
    
    (II)_3   :=∫_0^t∫_(∇)·Δ.

In two dimensions we have (II)_3=0 by elementary calculations. So we are left with estimating (II)_2 and obtain

    (II)_2   =∑_k∫_∫_0^t∂_γ·∂_γ(Φ()e_k β_k)
       =∑_k∫_∫_0^t∂_γ·∂_γ(_k(·,) β_k)
       =∑_k∫_∫_0^t ∇__k(·,)
    (∂_γ,∂_γ) β_k
       +∑_k∫_∫_0^t∂_γ·∂_γ_k(·,) β_k(σ)
    
    
    
       =:(II)_2^1+(II)_2^2.






On account of assumption (<ref>), Burkholder-Davis-Gundy inequality and Young's inequality we obtain for arbitrary δ>0

    [sup_0≤ t≤ T|(II)^1_2|]^r/2   ≤[sup_0≤ t ≤ T|∫_0^t∑_k∫_∇__k(·,)
    (∂_γ,∂_γ) β_k|]^r/2
       ≤ c [∑_k∫_0^T(∫_∇__k(·,)
    (∂_γ,∂_γ))^2]^r/4
       ≤ c [(∫_0^T(∫_
    |∂_γ|^2)^2]^r/4
       ≤δ [sup_0≤ t≤ T∫_
    |∇|^2]^r/2+ c(δ) [∫_0^T∫_
    |∇|^2]^r/2
       ≤δ [sup_0≤ t≤ T∫_
    |∇|^2]^r/2+c(δ) [1+_0_L^2_x^2]^r/2

using (<ref>) in the last step.
By similar arguments we gain

    [sup_0≤ t≤ T|(II)^2_2|]^r/2   ≤ c [∫_0^T(∫_∂_γ_k(·,)
    ·∂_γ)^2]^r/4
       ≤ c [(∫_0^T(∫_
    |∂_γ|||)^2]^r/4
       ≤ c [sup_(0,T)∫_
    ||^2+∫_0^T∫_
    |∇|^2]^r/2
       ≤ c [1+_0_L^2_x^2]^r/2.

Finally, we have by (<ref>)

    (III)   =1/2∫_∫_0^t  ⟨⟨∫_0^·∂_γ(Φ()  W)⟩⟩
    
       =1/2∑_k∫_∫_0^t  ⟨⟨∫_0^·∂_γ(Φ()e_k) β_k⟩⟩
       ≤1/2∑_k∫_0^t∫_| ∇__k(·,)∂_γ|^2
       +1/2∑_k∫_0^t∫_| ∂_γ_k(·,)|^2
       ≤ c ∫_0^t∫_|∇|^2+c ∫_0^t∫_||^2.

Hence we obtain by 
(<ref>) that

    [sup_0≤ t ≤ T|(III)|]^r/2   ≤ c [1+_0_L^2_x^2]^r/2.

Plugging all together and choosing δ small enough we have shown (<ref>).

The proof of (c) is similar: we simply differentiate once more.
We apply Itô's formula to the function f^β():=12∂^β_L^2_x^2 where β∈ℕ_0^2 is a multi-index of length 2. We obtain

    1/2∂^β(t)_L^2_x^2
       =1/2∂^β_0_L^2_x^2+∫_∫_0^t ∂^β·∂^β
       +1/2∫_∫_0^t⟨⟨∫_0^·∂^β(Φ()  W)⟩⟩
    =:(I)+(II)+(III),

where

    (II)   =-(II)_1+(II)_2-(II)_3,
    
    (II)_1   :=μ∫_0^t∫_ |∂^β∇|^2,
    
    (II)_2   :=∫_0^t∫_∂^β·∂^β(Φ()  W),
    
    (II)_3   :=∫_0^t∫_∂^β((∇))·∂^β.

The main difference is that (II)_3 does not vanish. By <cit.> (with m=2) and Young's inequality we have

    |∫_∂^β(⊗)·∂^β|   ≤ c_W^3,2_x^7/4_W^1,2_x^3/4_L^2_x^1/2
       ≤δ_W^3,2_x^2+c(δ)_W^1,2_x^6_L^2_x^4
       ≤δ_W^3,2_x^2+c(δ)(_W^1,2_x^10+_L^2_x^10),

where δ>0 is arbitrary. This implies

    [sup_0≤ t≤ T|(II)_3|]^r/2   ≤ δ[∫_0^T_W^3,2_x^2]^r/2+c(δ)[sup_0≤ t≤ T_W^1,2_x^2+sup_0≤ t≤ T_L^2_x^2]^5r/2
       ≤ δ[∫_0^T_W^3,2_x^2]^r/2+c(δ)[_0_W^1,2_x]^5r/2

due to (<ref>).
By arguments similar to the proof of (b), using (<ref>), we gain

    [sup_0≤ t≤ T|(II)_2|]^r/2
       ≤δ [sup_0≤ t≤ T^2_W_x^2,2]^r/2+ c(δ) [ sup_0≤ t≤ T_L_x^2^2+sup_0≤ t≤ T^2_W_x^1,2+∫_0^T_W_x^2,2^2]^r/2
       ≤δ [sup_0≤ t≤ T_W_x^2,2^2]^r/2+ c [1+_0_W^1,2_x^2]^r/2

using again (<ref>). Finally, we get
(<ref>) that

    [sup_0≤ t ≤ T|(III)|]^r/2   ≤ c [1+_0_W^1,2_x^2]^r/2

arguing again similarly to (b) and using (<ref>).
Again the claim follows by choosing δ small enough.

Under the assumptions of Lemma <ref> (b) equation (<ref>) is satisfied strongly in the analytical sense. That is
we have

    (t)   =(0)+∫_0^t[μΔ-(∇)
    -∇π_det]+
    ∫_0^t[Φ()+Φ^π]  W

-a.s. for all t∈[0,T], recall equation (<ref>). In the following we analyze the regularity of the pressure components π_det and Φ^π.
In the following the subscript w^* denotes Bochner-measurability with respect to the weak^*-topology.



  * Under the assumptions of Lemma <ref> we have

    π_det∈ L^r/2(Ω,L^2(0,T;L^2()),
    Φ^π∈ L^r(Ω;L^∞_w^*(0,T;L_2(𝔘;L^2()))).
 

  * Under the assumptions of Lemma <ref> (b) we have

    π_det∈ L^r/2(Ω,L^2(0,T;W^1,2()),
    Φ^π∈ L^r(Ω;L^∞_w^*(0,T;L_2(𝔘;W^1,2()))).
 

  * Under the assumptions of Lemma <ref> (c) we have 

    π_det∈ L^r/2(Ω,L^2(0,T;W^2,2()),
    Φ^π∈ L^r(Ω;L^∞_w^*(0,T;L_2(𝔘;W^2,2()))).
 



The key tool is the continuity of Δ^-1 (from W^k-2,p()→ W^k,p() for all k∈ℕ_0 and 1<p<∞). So, the regularity 
transfers from ⊗ to .
We obtain by Ladyshenskaya's inequality

    ∫_0^T∫_|π_det|^2   =∫_0^T∫_|Δ^-1(⊗)|^2≤ c ∫_0^T∫_||^4
       ≤ c ^2_L^∞_tL^2_x∇^2_L^2_tL^2_x≤ c (^4_L^∞_tL^2_x+∇^4_L^2_tL^2_x)

such that

    [∫_0^T∫_|π_det|^2]^r/4   ≤ c [^2_L^∞_tL^2_x+∇^2_L^2_tL^2_x]^r/2<∞

using Lemma <ref> (a). The estimates for (b) and (c) are very similar. We have to compare ∇π with |||∇|, where we have

    |||∇|≤||^2+|∇|^2∈ L^r/2(Ω;L^2(Q))

under the assumptions of Lemma <ref> (b).
Similarly, ∇^2π behaves as |||∇^2|+|∇|^2
which belongs to the same class since Lemma <ref> (c) applies.

Now, we investigate the regularity of Φ^π. We use continuity Δ^-1
to obtain

    [sup_0≤ t≤ TΦ^π_L_2(𝔘;L^2_x)^2]^r/2   =[sup_0≤ t≤ T∇Δ^-1Φ()_L_2(𝔘;L^2_x)^2]^r/2
       =[sup_0≤ t≤ T[∑_k≥1∇Δ^-1_k(·,)_L^2_x^2]^r/2
       ≤ c [sup_0≤ t≤ T∑_k≥1∫_0^T_k(·,)_L^2_x^2]^r/2
       ≤ c [1+sup_0≤ t≤ T_L^2_x^2]^r/2<∞

using (<ref>) and Lemma <ref> (a). Similarly, we have

    [sup_0≤ t≤ TΦ^π_L_2(𝔘;W^1,2_x)^2]^r/2   ≤ c [sup_0≤ t≤ T∑_k≥1_k(·,)_W^1,2_x^2]^r/2
       ≤ c [1+sup_0≤ t≤ T_W^1,2_x^2]^r/2<∞

using (<ref>) and Lemma <ref> (b), as well as

    [sup_0≤ t≤ TΦ^π_L_2(𝔘;W^2,2_x)^2]^r/2   ≤ c [sup_0≤ t≤ T∑_k≥1_k(·,)_W^2,2_x^2]^r/2
       ≤ c [1+sup_0≤ t≤ T_W^2,2_x^2]^r/2<∞

by (<ref>) and the Lemma <ref> (c).

Finally, we investigate the time-regularity of the velocity field.




  * Let the assumptions of Lemma <ref> (b) be satisfied for some r>2. Then we have

    [_C^α([0,T];L^2_x)]^r/2<∞

for all α<1/2.

  * Let the assumptions of Lemma <ref> (c) be satisfied for some r>2. Then we have

    [_C^α([0,T];W^1,2_x)]^r/2<∞

for all α<1/2.



We start with (a) and analyse each term in equation (<ref>) separately. Lemma <ref> (b) implies

    ∫_0^·Δ∈ L^r(Ω;W^1,2(0,T;L^2())).

As seen in the proof of Corollary <ref>  and ⊗ have the same regularity. In particular, Corollary <ref> (b) yields

    ∫_0^·((⊗)+∇)∈ L^r/2(Ω;W^1,2(0,T;L^2())).

Finally, we have

    ∫_0^·Φ()  W∈ L^r(Ω;C^α([0,T];L^2())).

by combing Lemma <ref> with (<ref>). The same conclusion holds for Φ^π using Lemma <ref> (a).
Plugging all together and noting the embedding W^1,2(0,T;X)↪ C^α([0,T];X) for any separably Banach space X the claim follows.

The proof of (b) follows along the same lines using the higher regularity from
Lemma <ref> (c), Corollary <ref> (c) and (<ref>).





 §.§ Discretization in space







We work with a standard finite element set-up for incompressible fluid mechanics, see e.g. <cit.> and <cit.>.
We denote by 𝒯_h a quasi-uniform subdivision of  into triangles of maximal diameter h>0.   
For 𝒮⊂ ^2 and ℓ∈ _0 we denote by
𝒫_ℓ(𝒮) the polynomials on 𝒮 of degree less than or equal
to ℓ. Moreover, we set 𝒫_-1(𝒮):=0. Let us
characterize the finite element spaces V^h() and P^h() as

    V^h()   := _h ∈ W^1,2() : _h|_𝒮∈𝒫_i(𝒮)  ∀𝒮∈𝒯_h,
    
    P^h()   :=π_h ∈ L^2() : π_h|_𝒮∈𝒫_j(𝒮)  ∀𝒮∈𝒯_h.

We will assume that i and j are both natural to get (<ref>) below (this is different from <cit.>, where also j=0 is allowed).
In order to guarantee stability of our approximations we relate V^h() and P^h() by the inf-sup condition, that is we assume that

    sup__h∈ V^h()∫__h π_h/∇_h_L^2_x≥ C π_h_L^2_x   ∀π_h∈ P^h(),

where C>0 does not depend on h. This gives a relation between i and j (for instance the choice (i,j)=(1,0) is excluded whereas (i,j)=(2,0) is allowed).
Finally, we define the space of discretely solenoidal finite element functions by 

    V^h_()   := {_h∈ V^h():  ∫__h  π_h=0 ∀π_h∈ P^h()}.



























































Let Π_h:L^2()→ V_^h() be the L^2()-orthogonal projection onto V_^h(). The following results concerning the approximability of Π_h are well-known (see, for instance <cit.>). There is c>0 independent of h such that we have
  
    ∫_|-Π_h /h|^2+ ∫_∇-∇Π_h ^2   ≤ c ∫_∇^2

for all ∈ W^1,2_() and
  
    ∫_|-Π_h /h|^2+ ∫_∇-∇Π_h ^2   ≤ c h^2 ∫_∇^2
          ^2

for all ∈ W^2,2_(). Similarly, if Π_h^π:L^2()→ P^h() denotes the L^2()-orthogonal projection onto P^h(), we have

    ∫_|p-Π_h^π p/h|^2   ≤ c ∫_∇
          p^2

for all p∈ W^1,2() and

    ∫_|p-Π_h^π p/h|^2   ≤ c h^2 ∫_∇^2
          p^2

for all p∈ W^2,2(). Note that (<ref>) requires the assumption j≥1 in the definition of P^h(), whereas (<ref>) also holds for j=0.




§ TIME-DISCRETIZATION

We consider an equidistant partition of [0,T] with mesh size Δ t=T/M and set t_m=mΔ t. 
Let _0 be an 𝔉_0-mesureable random variable with values in W^1,2_(). We aim at constructing iteratively a sequence of 𝔉_t_m-measurable random variables _m with values in W^1,2_() such that
for every ∈ W^1,2_() it holds true -a.s.

    ∫_   _m· +Δ t(∫_(∇_m)_m·+μ∫_∇_m:∇)
           =∫__m-1·+∫_Φ(_m-1) Δ_mW·,

where Δ_m W=W(t_m)-W(t_m-1). The existence of a unique _m (given _m and Δ_m W) solving (<ref>) is straightforward as it is a stationary Navier–Stokes system. The following result follows from Lemma 3.1 in <cit.>.
 
Assume that _0∈ L^2^q(Ω,W^1,2_div()) for some 1≤ q<∞. Suppose that Φ satisfies (<ref>). Then the iterates (_m)_m=1^M given by (<ref>)
are 𝔉_t_m-measurable. Moreover, the following estimates hold uniformly in M:

    [max_1≤ m≤ M_m^2^q_W^1,2_x+Δ t∑_k=1^M_m_W^1,2_x^2^q-2∇^2_m^2_L^2_x]   ≤ c(q,T,_0),
    [∑_k=1^M_m-_m-1^2_W^1,2_x∇_m^2_L^2_x]   ≤ c(T,_0),
    [(∑_k=1^M_m-_m-1^2_W^1,2_x)^4+(Δ t∑_k=1^M∇_m^2_L^2_x)^4]   ≤ c(T,_0).



Now, for ∈ W^1,2() we can insert -∇Δ^-1∈ W^1,2_() in (<ref>)
and obtain

    ∫__m·   +Δ t(∫_(∇_m)_m·+μ∫_∇_m:∇)
       =∫__m-1·
    +Δ t∫_π_m^det 
       +∫_Φ(_m-1) Δ_m W·+∫_∫_0^tΦ^π_m-1 Δ_m W·,

where 

    π_m^det   =-Δ^-1(_m⊗_m),
    Φ_m-1^π   =-∇Δ^-1Φ(_m-1).



Assume that _0∈ L^8(Ω,W^1,2_div()) and that Φ satisfies (<ref>).
For all m∈{1,...,M} the random variable π_m^det is 𝔉_t_m-measureable, has values in W^1,2() and we have

    [Δ t ∑_m=1^M∇π_m^det_L^2_x^2]≤ C

uniformly in Δ t. 


A corresponding result is shown in <cit.> for the full pressure provided the noise is divergence-free (in this case Φ_m^π vanishes). This is quite restrictive as it means the Φ is linear in .


The 𝔉_t_m-measurability of π_m^det follows directly from the one of _m stated in Lemma <ref>.
By continuity of the operator ∇Δ^-1 on L^2() we have

    ∇π_m^det_L^2_x^2   ≤ c (_m⊗_m)^2_L^2_x≤ c _m_L^4_x^2∇_m_L^4_x^2
       ≤ c _m_L^2_x∇_m_L^2_x^2∇^2_m_L^2_x≤ c _m_L^2_x^2+c ∇_m_L^2_x^4∇^2_m_L^2_x^2

-a.s. using also Ladyshenskaya's inequality v_L^4_x^2≤ c v_L^2_x∇ v_L^2_x which holds in two-dimensions. Now, summing with respect to m, applying expectations and using Lemma <ref> yields the claim.


For all m∈{1,...,M} the random variable Φ^π_m is 𝔉_t_m-measureable, has values in L_2(𝔘;W^1,2()) and we have

    [Δ t ∑_m=1^MΦ^π_m_L_2(𝔘;W^1,2_x)^2]≤ C

uniformly in Δ t. 


As for Lemma <ref> the proof mainly relies on the continuity of ∇Δ^-1 on L^2(). Here, we have by (<ref>)

    Φ^π_m^2_L_2(𝔘;W^1,2_x)   =∑_k≥1∇Δ^-1(Φ(_m)e_k)^2_W^1,2_x
       ≤ c ∑_k≥1Φ(_m)e_k^2_W^1,2_x=c Φ(_m)^2_L_2(𝔘;W^1,2_x)≤ c (1+∇_m^2_W^1,2_x).

Summing over m, applying expectations and using Lemma <ref>
finishes the proof.

Following <cit.> we set for ε>0

    Ω^ε_Δ t={ω∈Ω|max_1≤ m≤ M∇_m^2_L^2_x≤ -εlog (Δ t)}

such that

    (Ω_Δ t^ε)≥ 1-[max_1≤ m≤ M∇_m_L^2_x^2]/-εlog (Δ t)≥1+C/εlogΔ t

using Lemma <ref>. We obtain the following result.

Assume that (<ref>) holds and that _0∈ L^8(Ω,W^1,2_div()) is an 𝔉_0-measureable random variable. Let  be the unique strong solution
to (<ref>) in the sense of Definition <ref>. Assume that we have

    [_C^α([0,T];L^4_x)^4]<∞,  [_C^α([0,T];W^1,2_x)^2]<∞,

for some α∈(0,12); recall Corollary <ref>. Let (_m)_m=1^M be the solution to (<ref>). Then we have the error estimate

    [1_Ω_Δ t^ε(max_1≤ m≤ M(t_m)-_m^2_L^2_x   +Δ t∑_m=1^M  ∇(t_m)-∇_m^2_L^2_x)]≤ c_ε (Δ t)^2α-ε

for any ε>0.




  * Estimate (<ref>) improves the result from <cit.>, where the convergence rate is only α-ε.

  * In the paper <cit.> the time-discretization of the stochastic Navier–Stokes equations is analysed. The corresponding error does not contain an indicator function such as 1_Ω_Δ t^ε. The convergence
 rate is, however, only of logarithmic order.
 


Subtracting (<ref>) and (<ref>) we obtain the equation for the error _m=(t_m)-_m
which reads as

    ∫_   _m· +∫_t_m-1^t_m(∫_((∇)-(∇_m)_m)·+μ∫_∇_m:∇)
           =∫__m-1·+μ∫_t_m-1^t_m∫_(∇(t_m)-∇(t)):∇
           +∫_(∫_t_m-1^t_mΦ()  W-Φ(_m-1) Δ_mW)·

for all ∈ W^1,2_(). Choosing =_m implies

    1/2∫_   |_m|^2 +1/2∫_|_m-_m-1|^2  +μΔ t∫_|∇_m|^2
       =1/2∫_|_m-1|^2 
    -μ∫_t_m-1^t_m∫_(∇(t_m)-∇(t)):∇_m
       +∫_t_m-1^t_m∫_((∇)-(∇_m)_m)·_m
       +∫_(∫_t_m-1^t_mΦ()  W-Φ(_m-1) Δ_mW)·_m.

Iterating this equality and following the arguments in <cit.>
yields

    [   1_Ω_Δ t^ε(max_1≤ m≤ M(t_m)-_m^2_L^2_x+Δ t∑_m=1^M  ∇(t_m)-∇_m^2_L^2_x)]
       ≤ c (Δ t)^2α-ε
    +c [∑_m=1^M∫_t_m-1^t_m∫_|∇(t_m)-∇(t)|^2].

Note that only the first bound from (<ref>) has been used for (<ref>). The second bound from (<ref>), not used in <cit.>,
allows us to estimates the remaining integral by

    [∑_m=1^M∫_t_m-1^t_m∫_|∇(t_m)-∇(t)|^2]≤ c (Δ t)^2α

which finishes the proof.







§ FINITE ELEMENT BASED SPACE-TIME APPROXIMATION


Now we consider a fully practical scheme combining the implicite Euler scheme in time (as in the last section) with a finite element approximation in space. 
For a given h>0 let _h,0 be an 𝔉_0-mesureable random variable with values in V^h_() (for instance Π_h_0). We aim at constructing iteratively a sequence of random variables _h,m with values in V^h_() such that
for every ∈ V^h_() it holds true -a.s.

    ∫_   _h,m· +Δ t∫_((∇_h,m)_h,m-1+(_h,m-1)_h,m)·
       +μ Δ t∫_∇_m:∇=∫__h,m-1·+∫_Φ(_h,m-1) Δ_mW·,

where Δ_m W=W(t_m)-W(t_m-1). We quote the following result concerning the existence of solutions _h,m to (<ref>) from <cit.>.


Let 1≤ q<∞. Assume that _h,0∈ L^2^q(Ω,V_^h()) is an 𝔉_0-mesureable random variable and that [_h,0^2^q_L^2_x|]≤ K uniformly in h for some K>0. Suppose that Φ satisfies (<ref>). Then the iterates (_h,m)_m=1^M given by (<ref>)
are 𝔉_t_m-measurable. Moreover, the following estimate holds uniformly in M and h:

    [max_1≤ m≤ M_h,m^2^q_L^2_x+Δ t∑_k=1^M_h,m^2^q-2_L^2_x∇_h,m^2_L^2_x]   ≤ c(q,T,K).








 §.§ Error analysis



In this subsection we establish convergence with rates of the above defined algorithm. We introduce for ε>0 the sample set

    Ω^ε_h={ω∈Ω|max_0≤ m≤ M(∇_m^4_L^2_x+_h,m_L^2_x^2)≤ -εlog h}

which can be controlled by

    (Ω^ε_h)≥ 1-[max_1≤ m≤ M(∇_m^4_L^2_x+_h,m_L^2_x^2)]/-εlog h≥1+C/εlog h

using Lemma <ref> and <ref>. We also recall the definition of Ω^ε_Δ t in (<ref>). We are now ready to state our main result.


Let _0∈ L^2(Ω,W^1,2_())∩ L^8(Ω;L^2_()) be _0-measurable and assume that Φ satisfies (<ref>). Let  be the unique strong solution to (<ref>) in the sense of Definition <ref>.
Suppose further that

    [_C^α([0,T];L^4_x)^4]<∞,  [_C^α([0,T];W^1,2_x)^2]<∞,

for α∈(0,12); recall Corollary <ref>. Assume that LΔ t≤ (-εlog h)^-1 for some L>0.




Then we have

    [1_Ω^ε_Δ t∩Ω^ε_h(max_1≤ m≤ M(t_m)-_h,m_L^2_x^2   +∑_m=1^M Δ t ∇(t_m)-∇_h,m_L^2_x^2)]
       ≤ c (h^2+(Δ t)^2α-ε),

where (_h,m)_m=1^M is the solution to (<ref>) with _h,0=Π_h_0.
The constant c is (<ref>) is independent of M and h.



We do not expect that it is possible to avoid an indicator function in general (see <cit.> for the numerical approximation of stochastic PDEs with non-Lipschitz nonlinearities).
But it is not clear if our choice of the sample subset is optimal.










The rest of the paper is devoted to the proof of Theorem <ref>.
In fact Theorem <ref> will follow from combining Theorem <ref>
with the following theorem which estimates the error between the time-discrete solution _m to (<ref>)
and the solution _h,m to (<ref>).

Let _0∈ L^2(Ω;W^1,2_())∩ L^8(Ω;L^2_()) be _0-measurable and assume that Φ satisfies (<ref>). Let (_m)_m=1^M be the solution to (<ref>). Assume that LΔ t≤ (-εlog h)^-1 for some L>0.










Then we have

    [1_Ω^ε_h(max_1≤ m≤ M_m-_h,m_L^2_x^2   +Δ t ∑_m=1^M ∇_m-∇_h,m_L^2_x^2)]
       ≤ c (h^2+Δ t),

where (_h,m)_m=1^M is the solution to (<ref>) with _h,0=Π_h_0.
The constant c is (<ref>) is independent of M and h.


In the estimate of <cit.> the error is estimated by h^-3ε(h^2+(Δ t)+h^2Δ t). By our refined pressure analysis we are able to get rid of the term h^2Δ t which leads to a restrictive assumption between space and time-discretization. Additionally, we can remove the factor h^-3ε arising from the convective term, see (<ref>) below. Note that the error term h^2Δ t also appears in <cit.>, where  the finite-element based space-time discretization of the stochastic Stokes equations (that is the linearized version of (<ref>) without convective term) is studied. 




Define the error _h,m=_m-_h,m. Subtracting (<ref>) and (<ref>) we obtain

    ∫_   _h,m· +Δ t∫_μ(∇_m-∇_h,m):∇
       =∫__h,m-1·-Δ t∫_((∇_m)_m-((∇_h,m)_h,m-1+(_h,m-1)_h,m))·
       +∫_(Φ(_m)-Φ(_h,m-1)) Δ_mW·
       -∫_∇Δ^-1Φ(_m-1) Δ_mW·+Δ t∫_π_m^det

for every ∈ V_^h().
Setting =Π_h_h,m and applying the identity ·(-)=1/2(||^2-||^2+|-|^2) (which holds for any ,∈ℝ^n) we gain





















    ∫_   1/2(|Π_h_h,m|^2-|Π_h_h,m-1|^2+|Π_h_h,m-Π_h_h,m-1|^2) +Δ t∫_μ|∇_h,m|^2
       =Δ t∫_μ∇_h,m:∇(_m-Π_h_m)
       -Δ t∫_((∇_m)_m-((∇_h,m)_h,m-1+(_h,m-1)_h,m))·Π_h_h,m
       +Δ t∫_π_m^det Π_h_h,m
       +∫_(Φ(_m)-Φ(_h,m-1)) Δ_mW·Π_h_h,m
       -∫_∇Δ^-1Φ(_m-1) Δ_mW·Π_h_h,m
       =I_1(m)+… +I_5(m).

Young's inequality yields

    I_1(m)   ≤ κ Δ t∫_|∇_h,m|^2+c_κΔ t∫_|∇_m-∇Π_h∇_m|^2
       ≤ κ Δ t∫_|∇_h,m|^2+c_κ h^2 Δ t∫_|∇^2 _m|^2.
 
for every κ>0 using also (<ref>).
The convective term I_2(m) can be decomposed as

    I_2(m)   =I_2^1(m)+I_2^2(m)+I_2^3(m),
    
    I_2^1(m)   =-Δ t∫_(∇_h,m)_m-1·(_m-Π_h_m),
    
    I_2^2(m)   =Δ t∫_(∇_h,m)_h,m-1·(_m-Π_h_m)
       +Δ t∫_(_h,m)_h,m·(_m-Π_h_m),
    
    I_2^3(m)   =-Δ t∫_(∇Π_h_h,m)_h,m-1·_m
       -Δ t∫_(_h,m-1)π_h_h,m·_m,

As in <cit.>
we introduce the sample set

    Ω_h,m^ε={ω∈Ω|max_0≤ n≤ m(∇_n^4_L^2_x+_h,n_L^2_x^2)≤ -εlog h}

and obtain

    1_Ω^ε_h,m-1I_2^2(m)
       ≤ κ1_Ω^ε_h,m-1(∇_h,m-1^2_L^2_x+∇_h,m^2_L^2_x)
       +c_κ Δ t  h^4(max_1≤ n≤ m_h,n_L^2_x^2)∇_m^2_L^2_x∇^2_m_L^2_x^2,
    1_Ω^ε_h,m-1I_2^3(m)
       ≤ κ1_Ω^ε_h,m-1(∇_h,m-1^2_L^2_x+∇_h,m^2_L^2_x)
       +c_κ Δ tlog(h^-ε)(max_1≤ n≤ m1_Ω^ε_h,n-1_h,n_L^2_x)
       +c_κ Δ t∇(_m-_m-1)^2_L^2_x(max_1≤ n≤ m∇_n_L^2_x^2)(max_0≤ n≤ m_h,n_L^2_x^2),

where κ>0 is arbitrary. For I_2^1(m), however, we get a slightly better estimate than in <cit.> since our definition of I_2^1(m) makes use of ^m=0. We get

    1_Ω^ε_h,m-1I_2^1(m)   ≤ κ1_Ω^ε_h,m-1∇_h,m^2_L^2_x+c_κ h^3log(h^-ε)∇^2_m_L^2_x^2

due to (<ref>).

The crucial point in this proof (making the essential difference to <cit.>) are the estimates for the pressure terms
I_3(m) and I_5(m). We will estimate I_3(m) first whereas I_5(m) will only be bounded after iterating (<ref>), see the estimates for ℳ_m,2 below. By (<ref>) we have

    I_3(m)   =Δ t∫_(π_m^det-Π_h^ππ_m^det)Π_h_h,m
       ≤ c_κΔ t  ∫_|π_m^det-Π_h^ππ_m^det|^2+ κΔ t ∫_|∇Π_h_h,m|^2
       ≤  c_κΔ t h^2 ∫_|∇π_m^det|^2+ κΔ t ∫_|∇_h,m|^2,

where κ>0 is arbitrary.
Plugging all together and choosing κ small enough (to absorb the corresponding terms to the left-hand side) we have shown

    1_Ω^ε_h,m-1(∫_|Π_h_h,m|^2  +∫_|Π_h_h,m-Π_h_h,m-1|^2 +Δ t∫_|∇_h,m|^2)
       ≤ 1_Ω^ε_h,m-1∫_|_h,m-1|^2 
    +c Δ tlog(h^-ε)(max_1≤ n≤ m1_Ω^ε_h,n-1_h,n_L^2_x)
       + cΔ t h^2(∫_|∇π_m^det|^2+∫_(1+|∇_m|^2)∫_|∇^2_m|^2)
       +c Δ t  h^4(max_1≤ n≤ m_h,n_L^2_x^2)∫_|∇_m|^2∫_|∇^2_m|^2
       +c Δ t∫_|∇(_m-_m-1)|^2(max_1≤ n≤ m∫_|∇_n|^2)(max_1≤ n≤ m∫_|_h,n|^2)
       +c 1_Ω^ε_h,m-1∫_(∫_t_m-1^t_mΦ(_m-1)-Φ(_h,m-1))  W)·Π_h_h,m
       -c 1_Ω^ε_h,m-1∫_(∫_t_m-1^t_m∇Δ^-1Φ(_m-1)  W)·Π_h_h,m.










Iterating this inequality yields

    1_Ω^ε_h,m-1∫_|Π_h_h,m|^2  +∑_n=1^m1_Ω^ε_h,n-1(∫_|Π_h_h,n-Π_h_h,n-1|^2 +Δ t∫_|∇_h,n|^2)
       ≤ ∫_|_h,0|^2 +c Δ tlog(h^-ε)∑_n=1^m(max_1≤ℓ≤ n1_Ω^ε_h,ℓ-1_h,ℓ^2_L^2_x)
       + c h^2Δ t∑_n=1^m(∫_|∇π_n^det|^2+∫_(1+|∇_n|^2)∫_|∇^2_n|^2)
       +c Δ t  h^4∑_n=1^m(max_1≤ℓ≤ n_h,ℓ_L^2_x^2)∫_|∇_n|^2∫_|∇^2_n|^2
       +c Δ t∑_n=1^m∫_|∇(_n-_n-1)|^2(max_1≤ℓ≤ n∫_|∇_ℓ|^2)(max_1≤ℓ≤ n∫_|_h,ℓ|^2)
       +c ∑_n=1^m1_Ω^ε_h,n-1∫_(∫_t_n-1^t_nΦ(_n-1)-Φ(_h,n-1))  W)·Π_h_h,n
       -c ∑_n=1^m1_Ω^ε_h,n-1∫_(∫_t_n-1^t_n∇Δ^-1Φ(_n-1)  W)·Π_h_h,n.

Now, we explain how to bound line by line in expectation. The error in the initial datum can be bounded by h^2 by (<ref>) and the assumption _0∈ L^2(Ω;W^1,2_()). The second term on the right-hand side can be handled by the discrete Gronwall lemma using the assumption LΔ t≤ (-εlog h)^-1. The expectation of the second line is bounded by h^2 using Lemmas <ref> and <ref> (the estimate for ∇π_m^det is the first main ingredient in our proof). This estimate is better than the corresponding one in <cit.> as only the deterministic part of the pressure appears here.
The expectation of the third line is bounded by

    c h^4([max_1≤ℓ≤ M_h,ℓ_L^2_x^2]^2)^1/2([Δ t ∑_n=1^M∫_|∇_n|^2∫_|∇^2_n|^2]^2)^1/2≤ ch^4

as a consequence of Lemmas <ref> and <ref>.
The fourth line is controlled by Δ t due to the estimate

    [∑_n=1^m∫_|∇(_n-_n-1)|^2(max_1≤ℓ≤ n∫_|∇_ℓ|^2)(max_1≤ℓ≤ n∫_|_h,ℓ|^2)]
       ≤([∑_n=1^M∫_|∇(_n-_n-1)|^2]^2)^1/2([max_1≤ℓ≤ M∫_|∇_ℓ|^2]^4)^1/4×
           ×([max_1≤ℓ≤ M∫_|_h,ℓ|^2]^4)^1/4

and the uniform bounds from Lemmas <ref> and <ref> (with q=3).
It remains to estimate the two stochastic terms

    ℳ_m,1   =
    ∑_n=1^m1_Ω^ε_h,n-1∫_∫_t_n-1^t_n(Φ(_n-1)-Φ(_h,n-1))  W·Π_h_h,n,
    ℳ_m,2   =
    ∑_n=1^m1_Ω^ε_h,n-1∫_∫_t_n-1^t_n(Id-Π_h^π)Δ^-1Φ(_n-1)  W Π_h_h,n.

Note that we used integration by parts and Π_h_h,n∈ V^h_() together with the definition of Π_h^π to rewrite ℳ_m,2 into the form above. Finally, we write

    ℳ_m,1   =
    ∑_n=1^m1_Ω^ε_h,n-1∫_∫_t_n-1^t_n(Φ(_n-1)-Φ(_h,n-1))  W·Π_h_h,n-1
       +∑_n=1^m1_Ω^ε_h,n-1∫_∫_t_n-1^t_n(Φ(_n-1)-Φ(_h,n-1))  W· (Π_h_h,n-Π_h_h,n-1)
       =:ℳ_m,1^1+ℳ_m,1^2

as well as

    ℳ_m,2=
    ∑_n=1^m1_Ω^ε_h,n-1∫_∫_t_n-1^t_n(Id-Π_h^π)Δ^-1Φ(_n-1)  W Π_h_h,n-1
       +∑_n=1^m1_Ω^ε_h,n-1∫_∫_t_n-1^t_n(Id-Π_h^π)Δ^-1Φ(_n-1)  W  (Π_h_h,n-Π_h_h,n-1)
       =:ℳ_m,2^1+ℳ_m,2^2.

These representations have the advantage that ℳ_m,1^1 and ℳ_m,2^2 are martingales (note the index n-1 in the indicator functions).
Consequently, we can apply the Burkholder-Davis-Gundy inequality to estimate them. As far as ℳ_m,1 is concerned we have

    [max_1≤ m≤ M|ℳ_m,1^1|]
       ≤ c [∑_n=1^M1_Ω^ε_h,n-1∫_t_n-1^t_nΦ(_n-1)-Φ(_h,n-1)^2_L_2(𝔘,L^2_x)Π_h_h,n-1^2_L^2_x]^1/2
       ≤ c [max_0≤ n≤ M1_Ω^ε_h,nΠ_h_h,n_L^2_x(∑_n=1^M1_Ω^ε_h,n-1∫_t_n-1^t_nΦ(_n-1)-Φ(_h,n-1)^2_L_2(𝔘,L^2_x))^1/2]
       ≤ κ [max_0≤ n≤ M1_Ω^ε_h,nΠ_h_h,n^2_L^2_x]+ c_κ [∑_n=1^M1_Ω^ε_h,n-1∫_t_n-1^t_n_n-1-_h,n-1_L^2_x^2]
       ≤ κ [max_0≤ n≤ M1_Ω^ε_h,nΠ_h_h,n^2_L^2_x]+ c_κ [Δ t∑_n=1^M1_Ω^ε_h,n-1Π_h_h,n-1_L^2_x^2]
       +c_κ [Δ t∑_n=1^M1_Ω^ε_h,n-1_n-1-Π_h_n-1_L^2_x^2]

Here, we also used (<ref>), (<ref>) as well as Young's inequality for κ>0 arbitrary. Now, the first term can be absorbed for κ small enough.

The second term can be handled by the discrete Gronwall lemma
(note that Ω^ε_h,n⊂Ω^ε_h,n-1 such that 1_Ω^ε_h,n≤1_Ω^ε_h,n-1 a.s.)
Using (<ref>) the last term can be estimated by

    c_κ h^2 [Δ t∑_n=1^M1_Ω^ε_h,n-1∇_n-1_L^2_x^2]

which is bounded by ch^2 using Lemma <ref> (a) (recall that _0∈ L^2(Ω;W^1,2())).
For the term ℳ_m,1^2 we obtain by Cauchy-Schwartz inequality, Young's inequality, Itô-isometry and (<ref>)

    [max_1≤ m≤ M|ℳ_m,1^2|]   ≤ κ [ ∑_n=1^M1_Ω^ε_h,n-1Π_h_h,n-Π_h_h,n-1_L^2_x^2]
       +c_κ [∑_n=1^M1_Ω^ε_h,n-1∫_t_n-1^t_n(Φ(_n-1)-Φ(_h,n-1))  W  _L^2_x^2]
       ≤ κ [ ∑_n=1^M 1_Ω^ε_h,n-1Π_h_h,n-Π_h_h,n-1_L^2_x^2 ]
        + c_κ [∑_n=1^M1_Ω^ε_h,n-1∫_t_n-1^t_n_n-1-_h,n-1_L^2_x^2].

The first term can be absorbed for κ small enough whereas the second one can be estimated as for ℳ_m,1^1.

















Now, we come to the second main ingredient which is the estimate for ℳ_m,2.
We obtain using (<ref>), (<ref>) and continuity of ∇^2Δ^-1







    [max_1≤ m≤ M|ℳ_m,2^1|]
       ≤ c [∑_n=1^M1_Ω_h,n-1∫_t_n-1^t_n(Id-Π_h^π)Δ^-1Φ(_n-1)^2_L_2(𝔘,L^2_x)∇Π_h_h,n-1^2_L^2_x]^1/2
       ≤ c h^2 [max_1≤ n≤ M∇^2Δ^-1Φ(_n-1)_L_2(𝔘,L^2_x)(∑_n=1^M1_Ω_h,n-1∫_t_n-1^t_n∇Π_h_h,n^2_L^2_x)^1/2]
       ≤ c_κ h^4 [max_1≤ n≤ M∇_n-1^2_L^2_x]+ κ [Δ t∑_n=1^M1_Ω_h,n-1∇Π_h_h,n_L^2_x^2].


The first term is bounded by
[max_m|ℳ_m,2^1|]≤ ch^4 using Lemmas <ref> (recall that _0∈ L^2(Ω;W^1,2())). The second term
can be absorbed if κ≪1.







As a consequence of Young's inequality, inverse estimates on V^h(), Itô-ismotry and (<ref>) we infer

    [max_1≤ m≤ M|ℳ_m,2^2|]   ≤ κ h^2 [ ∑_n=1^M 1_Ω_h,n-1∇(Π_h_h,n-Π_h_h,n-1)_L^2_x^2]
       +c_κ h^-2 [∑_n=1^M1_Ω_h,n-1∫_t_n-1^t_n(Id-Π_h^π)Δ^-1Φ(_n-1)  W  _L^2_x^2]
       ≤ c κ [ ∑_n=1^M1_Ω_h,n-1Π_h_h,n-Π_h_h,n-1_L^2_x^2]
       +c_κ h^-2 [∑_n=1^M1_Ω_h,n-1∫_t_n-1^t_n(Id-Π_h^π)Δ^-1Φ(_n-1)_L_2(𝔘;L^2_x)^2]
       ≤ c κ [ ∑_n=1^M1_Ω_h,n-1Π_h_h,n-Π_h_h,n-1_L^2_x^2]
       +c_κ h^2 [∑_n=1^M∫_t_n-1^t_n∇^2Δ^-1Φ(_n-1)_L_2(𝔘;L^2_x)^2].

The first term can be absorbed for κ small enough. 





 Arguing as for ℳ_m,2^1
the last term is bounded by ch^2 [max_n∇_n-1^2_L^2_x].














Plugging all together and noting that Ω^ε_h⊂⋃_n=1^M Ω^ε_h,n shows

    [1_Ω^ε_h(max_1≤ m≤ MΠ_h_h,m_L^2_x^2   +Δ t ∑_m=1^M ∇_h,m_L^2_x^2)]≤ c (h^2+Δ t).


Recalling that _h,m=_m-Π_h_m+Π_h_h,m and using (<ref>) as well as Lemma <ref> (a) gives the claim.



[M]
BeTe A. Bensoussan and R. Temam.
Équations stochastiques du type Navier–Stokes.
J. Funct. Anal. 13, 195–222, 1973.




	
	












Bi B. Birnir: The Kolmogorov–Obukhov Statistical Theory
of Turbulence. J. Nonlinear Sci. 23, 657–688. (2013)






BeBrMi H. Bessaih, Z. Brzeźniak, A. Millet (2014): Splitting up method for the 2D stochastic Navier-Stokes
equations. Stochastic PDE: Anal. Comput. 2, 433–470.
BeMi H. Bessaih, A. Millet (2018): Strong L^2 convergence of time numerical schemes for the stochastic two-dimensional Navier–Stokes equations IMA J. Num. Anal. doi:10.1093/imanum/dry058
Br D. Breit: Existence theory for stochastic power law fluids. J. Math. Fluid Mech. 17, 295–326. (2015)
Br2 D. Breit: Existence theory for generalized Newtonian fluids.
Mathematics in Science and Engineering. Elsevier/Academic Press, London. (2017)
BF F. Brezzi and M. Fortin, Mixed and Hybrid Finite Element Methods, Springer Ser. Comput.
Math. 15, Springer-Verlag, New York, 1991.
BCP
Z. Brzeźniak, E. Carelli, J. A. Prohl (2013): Finite-element-based discretizations of the incompressible Navier–Stokes
equations with multiplicative random forcing, IMA J. Num. Anal. 33, 771–824.



CP
E. Carelli, A. Prohl (2012): Rates of convergence for discretizations of the stochastic incompressible Navier-Stokes equations. SIAM J. Numer. Anal. 50(5), 2467–2496.
CHP E. Carelli, E. Hausenblas, A. Prohl (2012): Time-splitting methods to solve the stochastic incompressible Stokes equation. SIAM J. Numer. Anal. 50, 2917–2939. 
Ca M. Capiński, A note on uniqueness of stochastic Navier-Stokes equations, Univ. Iagell. Acta Math. 30 (1993), 219–228.
CC M.Capiński, N. J. Cutland, Stochastic Navier-Stokes equations, Acta Appl. Math. 25 (1991), 59–85.
   PrZa G. Da Prato, J. Zabczyk (1992): Stochastic Equations in Infinite Dimensions, Encyclopedia
Math. Appl., vol. 44, Cambridge University Press, Cambridge. 





















FlaGat F. Flandoli and D. Gatarek
Martingale and stationary solutions for stochastic Navier-Stokes equations.
Probab. Theory Relat. Fields, 102, 367–391, 1995.






	

GR Girault,V. & Raviart, P. A. (1981) Finite Element Method for Navier–Stokes Equations: Theory and Algorithms. Berlin/Heidelberg/New York: Springer.
GL V. Girault and J.L. Lions, Two-grid finite-element schemes for the steady Navier-Stokes
problem in polyhedra, Port. Math. (N.S.), 58 (2001), pp. 25–57.
GS V. Girault and L.R. Scott, A quasi-local interpolation operator preserving the discrete divergence,
Calcolo, 40 (2003), pp. 1–19.




HR J.G. Heywood and R. Rannacher, Finite element approximation of the nonstationary
Navier–Stokes problem. I. Regularity of solutions and second-order error estimates for
spatial discretization. SIAM J. Numer. Anal. 19, 275–311, 1982.
	Ho1 M. Hofmanová (2013): Degenerate Parabolic Stochastic Partial Differential Equations. Stoch. Pr. Ap. 123 (12), 4294–4336.






KukShi
S. Kuksin and A. Shirikyan.
Mathematics of two-dimensional turbulence, volume 194 of Cambridge Tracts in Mathematics.
Cambridge University Press, Cambridge, 2012.










MikRoz R. Mikulevicius and B. L. Rozovskii. Stochastic Navier-Stokes equations for turbulent flows. SIAM J. Math. Anal., 35(5):1250–1310, 2004.



 


  Pa E. Pardoux (1975): Equations aux dérivées Partielles stochastiques non linéaires monotones. Etude de solutions fortes de
type Itô, Ph.D. thesis, Université Paris Sud.

PrRo  C. Prévôt, M. Röckner (2007): A concise course on stochastic partial differential equations. Lecture Notes in Mathematics, 1905. Springer, Berlin.

Pr J. Printems (2001): On the discretization in time of parabolic stochastic partial differential equations.
Math. Mod. Numer. Anal. 35, 1055–1078.


Ro 
M. Romito.
Some probabilistic topics in the Navier–Stokes equations. 
Recent progress in the theory of the Euler and Navier–Stokes equations.
London Math. Soc. Lecture Note Ser. 430, 175–232,  Cambridge Univ. Press, Cambridge, 2016.


ScoZha90
L. R. Scott and S. Zhang, Finite element interpolation of nonsmooth
  functions satisfying boundary conditions, Math. Comp. 54 (1990),
  no. 190, 483–493.










Ya1 Y. Yan (2004): Semidiscrete Galerkin Approximation for a Linear Stochastic Parabolic Partial Differential Equation Driven by an Additive Noise. Num. Math. 44, 829–847.

Ya2 Y. Yan (2005): Galerkin finite element methods for stochastic parabolic partial differential equations,
SIAM J. Numer. Anal., 43, pp. 1363–1384.

Yo N. Yoshida (2012): Stochastic Shear thickening fluids: strong convergence of the Galerkin approximation and the energy inequality. The Annals of Applied Probability 2012, Vol. 22, No. 3, 1215–1242




