Abstract

We propose some new mixed finite element methods for the time-dependent stochastic Stokes equations with multiplicative noise, which use the Helmholtz decomposition of the driving multiplicative noise. It is known (Langa, J. A., Real, J. & Simon, J. (2003) Existence and regularity of the pressure for the stochastic Navier--Stokes equations. Appl. Math. Optim., 48, 195--210) that the pressure solution has low regularity, which manifests in suboptimal convergence rates for well-known inf-sup stable mixed finite element methods in numerical simulations; see Feng X., & Qiu, H. (Analysis of fully discrete mixed finite element methods for time-dependent stochastic Stokes equations with multiplicative noise. arXiv:1905.03289v2 [math.NA]). We show that eliminating this gradient part from the noise in the numerical scheme leads to optimally convergent mixed finite element methods and that this conceptual idea may be used to retool numerical methods that are well known in the deterministic setting, including pressure stabilization methods, so that their optimal convergence properties can still be maintained in the stochastic setting. Computational experiments are also provided to validate the theoretical results and to illustrate the conceptual usefulness of the proposed numerical approach.

1. Introduction

This paper is concerned with fully discrete mixed finite element approximations of the following time-dependent stochastic Stokes equations with multiplicative noise for viscous incompressible fluids covering the domain |$D=(0,L)^d$| for |$d=2,3$|⁠:
(1.1a)
 
(1.1b)
 
(1.1c)
where |$ \textbf {u}$| and |$p$|⁠, respectively, denote the velocity field and the pressure of the fluid, which are spatially periodic with period |$L>0$| in each coordinate direction; |$\mathbf {u}_0$| and |$\mathbf {f}$| denote, respectively, the prescribed initial velocity and body force, which are spatially periodic (see Section 2 for the details). For the sake of simplicity and ease of presentation, we assume |$\{W(t); t\geq 0\}$| to be an |${\mathbb R}$|-valued Wiener process; see Section 2 for further details.

When |$ \textbf {B} \equiv \textbf {0}$|⁠, (1.1) is the well-known (deterministic) Stokes system; one motivation for studying (1.1a)–(1.1b) with ‘random force’ |$\mathbf {f}+ \textbf {B}( \textbf {u})\frac {\textrm{d}W}{\textrm{d}t}$| is to develop mathematical models of this type for turbulent fluids (Bensoussan, 1995; Hairer & Mattingly, 2006). In addition to their importance in applied sciences and engineering, the Stokes equations are a well-known partial differential equation (PDE) model with saddle point structure, which requires special numerical discretizations to construct optimally convergent methods. It should be noted that although the involved deterministic Stokes operator is linear, system (1.1a)–(1.1b) is nonlinear due to the nonlinear function |$ \textbf {B}$|⁠.

The numerical analysis of the deterministic Stokes problem is well established in the literature; see Brezzi & Fortin (1991), Girault & Raviart (1986), Heywood & Rannacher (1982). Well-known numerical methods include: exactly divergence-free methods, which approximate the velocity in exactly divergence-free finite element spaces; mixed finite element methods, where the (discrete) inf-sup condition is the key criterion that distinguishes stable pairings of finite element ansatz spaces for the velocity (with more degrees of freedom) and the pressure (with fewer degrees of freedom); and mixed methods, which allow a more flexible, broader application if compared to exactly divergence-free methods, thus putting them in the center of research on numerical methods for saddle point problems in recent decades. Another class of related numerical methods are stabilization methods, which were initiated in Hughes et al. (1986), where the incompressibility constraint (1.1b) is relaxed into
(1.2)
This relaxation allows for stable pairings of equal order (nodal-based) finite element ansatz spaces for both velocity and pressure (putting |$\varepsilon = {\mathcal O}(h^2)$|⁠, where |$h>0$| is the spatial mesh size). We remark that optimal-order error estimates had been obtained for all three classes of finite element methods in the deterministic setting (cf. Girault & Raviart, 1986; Brezzi & Fortin, 1991), where
This work contributes to the numerical analysis of the stochastic Stokes problem (1.1) (i.e., |$ \textbf {B} \neq \textbf {0}$|⁠). By Langa et al. (2003), even for smooth datum functions |$\mathbf {u}_0$| and |$\mathbf {f}$|⁠, the (temporal) regularity of the pressure |$p\in L^1\bigl (\varOmega ; W^{-1,\infty }(0,T;H^1(D)/\mathbb {R} )\bigr )$| is limited in general due to the driving noise. In order to motivate its impact on the pressure, we here discuss the related question of |$k$|-independent stability estimates for the pair of random variables |$( \textbf {u}^{n+1}, p^{n+1})$| of the following time-implicit discretization of (1.1) on a uniform mesh of |$[0,T]$| with mesh size |$k>0$|⁠:
(1.3a)
 
(1.3b)
where |$\varDelta _{n+1} W\!:=\!W(t_{n+1})-W(t_{n})\thicksim {\mathcal N}(0,k)=\sqrt {k}{\mathcal N}(0, 1)$| and |$\mathbf {f}^{n+1} = \mathbf {f}(t_{n+1},\cdot ) \in L^2(\varOmega , L^2_\textrm{per}(D;\mathbb {R}^d))$|⁠. A crucial observation for the motivation of this paper is that the pressure gradient on the left-hand side is scaled by |$k$|⁠, while the noise term is |$\mathcal{O}(\sqrt {k})$|⁠. Let us assume that estimate (3.5) in Lemma 3.2 for |$\{\varDelta \textbf {u}^n\}_{n}$| taking values in |$L^2(D; {\mathbb R}^d)$| has already been shown, and we now look for a uniform bound for |$\{\nabla p^n\}_n$| taking values in |$L^2(D; {\mathbb R}^d)$|⁠. The strategy for deriving such a stability estimate is to fix one |$\omega \in \varOmega $| and to multiply (1.3) by |$\nabla p^{n+1}(\omega )$|⁠: all the terms that involve the velocity vanish due to the incompressibility property and the periodic boundary condition, and we end up with
(1.4)
Note that the term on the right-hand side does not vanish since |$ \textrm {div}\, \textbf {B}( \textbf {u}^n)\neq 0$| for a general (Lipschitz) nonlinear mapping |$ \textbf {B}$|⁠. We now take expectations |${\mathbb E}[\cdot ]$| on both sides, sum over all time steps and use (3.5), the facts that |$ \textbf {B}( \textbf {u}^n)$| and |$\varDelta _{n+1} W$| are independent and |${\mathbb E}\bigl [\vert \varDelta _{n+1} W\vert ^2\bigr ] \leq Ck$| and Young’s inequality (with |$\alpha>0$|⁠) to obtain the estimate
Taking |$\alpha = \frac {1}{k}$| allows the last term on the right-hand side to be absorbed into the one on the left, but the remaining term is |$\sum _{n=0}^{N-1} {\mathbb E}\bigl [ \Vert \textbf {B} ( \textbf {u}^n )\Vert ^2\bigr ] \propto {\mathcal O}(k^{-1})$|⁠; therefore, we end up with the following |$k$|-dependent estimate:
(1.5)

The above consideration crucially affects the error analysis of a space-time discretization of (1.1a)–(1.1b):

  • Exactly divergence-free methods require restricted settings of data, including the dimension, topology and regularity of the spatial domain |$D$|⁠. However, an optimal-order error estimate can be proved for the velocity approximation; see Carelli & Prohl (2012), which uses the fact that no pressure is involved in the analysis.

  • The error estimate for the velocity approximation of inf-sup stable mixed finite element methods in Feng & Qiu (2019) was obtained based on a stability bound of type (1.5) to bound the related best-approximation error for the pressure that appears in (an auxiliary temporal discretization of) (1.1), thus leading to a sub-optimal error estimate for the velocity of order |${\mathcal O}(k^{\frac 12} + h {k}^{-\frac 12})$|⁠. The computational studies in Feng & Qiu (2019) suggest that this error bound is sharp.

The first goal of the paper is to construct optimally convergent inf-sup stable mixed finite element methods, with ‘minimum’ extra effort. Our main idea, which is partly borrowed from Carelli et al. (2012), is to perform the Helmholtz decomposition for the noise term at each time step first and then to determine the new velocity and pressure iterates simultaneously via the mixed finite element method. Below we shall use the semidiscrete time-stepping scheme (1.3) to motivate our strategy. Introducing the Helmholtz decomposition of |$\mathbf {B}$| as:
(1.6)
and setting |$r^{n+1}:= p^{n+1} - k^{-1} \varDelta _{n+1} W \xi ^{n}$|⁠, then (1.3) can be rewritten as
(1.7a)
 
(1.7b)
In contrast to estimate (1.5) for |$p^{n+1}$|⁠, it can be shown that the new pressure |$r^{n+1}$| satisfies the following improved stability estimate (see Lemma 3.2):
(1.8)
which is a consequence of the divergence-free property of the modified noise term (i.e., the last term on the right-hand side of (1.7a)). Conceptually, this improved stability for the new pressure |$r^{n+1}$| is obtained by removing the stochastic pressure  |$\xi ^n$| from the driving noise in (1.3a). As will be detailed in Section 4, any inf-sup stable mixed finite element discretization of (1.7) then gives optimally convergent velocity approximations (see Theorem 4.5, whose proof essentially relies on (1.8)). We also present optimal-error estimates for (temporal averages of) the pressure approximations in |$L^2$|⁠, which improve corresponding suboptimal estimates in Feng & Qiu (2019).
We therefore conclude by saying that it is essential to identify the proper role of the semidiscrete pressures, namely |$\{p^n\}_n$| in (1.3) vs. |$\{r^n\}$| in (1.7), for inf-sup stable mixed finite element methods for (1.1) in order to construct optimally convergent mixed methods. Moreover, this insight also suggests how to construct optimally convergent stabilization methods for (1.1) that circumvent the inf-sup stability criterion for mixed element methods, and hence allow a more efficient discretization such as
(1.9a)
 
(1.9b)
for which |$\varepsilon = {\mathcal O}(h^2)$| will be shown to be the optimal choice in Section 5. The error analysis in Section 5 verifies optimal-order convergence for a standard finite element discretization of (1.9) which employs the same finite element space for approximating both |$ \textbf {u}_{\varepsilon }^{n+1}$| and |$r^{n+1}_\varepsilon $|⁠; see Theorem 5.3. Corresponding computational studies in Section 6 support the conclusion that the choice of pressure in the stabilization is crucial for achieving an optimally convergent stabilization method for (1.1).

The remainder of this paper is organized as follows. In Section 2 we give exact assumptions on the data in (1.1) and recall the definition and known properties of the (strong) variational solution for problem (1.1). In Sections 3 and 4 we analyze the Helmholtz decomposition enhanced Euler–Maruyama time-stepping scheme (1.6)–(1.7) and its mixed finite element approximations and establish optimal convergence for both. Section 5 establishes optimal convergence for the stabilized scheme (1.9) and its equal-order finite element approximations. Two-dimensional numerical experiments and computational studies are given in Section 6 to validate the theoretical error bounds and to computationally evidence that a proper selection of the pressure for the construction of optimally convergent mixed methods is indeed necessary.

2. Preliminaries

2.1 Notation

Standard function and space notation will be adopted in this paper. For example, |$H^\ell _\textrm{per}(D, {\mathbb R}^d)\, (\ell \geq 0)$| denotes the subspace of the Sobolev space |$H^\ell (D, {\mathbb R}^d)$| consisting of |${\mathbb R}^d$|-valued periodic functions with period |$L$| in each spatial coordinate direction and |$(\cdot ,\cdot ):=(\cdot ,\cdot )_D$| denotes the standard |$L^2$|-inner product, with induced norm |$\Vert \!\cdot\! \Vert $|⁠. Let |$(\varOmega ,{\mathcal {F}}, \{{\mathcal {F}}_t\},{\mathbb {P}})$| be a filtered probability space with the probability measure |${\mathbb {P}}$|⁠, the |$\sigma $|-algebra |${\mathcal {F}}$| and the continuous filtration |$\{{\mathcal {F}}_t\} \subset {\mathcal {F}}$|⁠. For a random variable |$v$| defined on |$(\varOmega ,{\mathcal {F}}, \{{\mathcal {F}}_t\},{\mathbb {P}})$|⁠, let |${\mathbb E}[v]$| denote the expected value of |$v$|⁠. For a vector space |$X$| with norm |$\|\cdot \|_{X}$|⁠, and |$1 \leq p < \infty $|⁠, we define the Bochner space |$\bigl (L^p(\varOmega ,X); \|v\|_{L^p(\varOmega ,X)} \bigr )$|⁠, where |$\|v\|_{L^p(\varOmega ,X)}:=\bigl ({\mathbb E} [ \Vert v \Vert _X^p]\bigr )^{\frac 1p}$|⁠. Throughout this paper, unless stated otherwise we shall use |$C$| to denote a generic positive constant that may depend on |$T$|⁠, the datum functions |$\mathbf {u}_0$| and |$\mathbf {f}$| and the domain |$D$| but is independent of |$k$| and the mesh parameter |$h$|⁠.

We also define
We recall from Girault & Raviart (1986) that the (orthogonal) Helmholtz projection |$ \textbf {P}_{{\mathbb H}}: L^2_\textrm{per}(D; {\mathbb R}^d)$|  |$\rightarrow {\mathbb H} $| is defined by |$ \textbf {P}_{{\mathbb H}} \textbf {v} = \pmb {\eta }$| for every |$ \textbf {v} \in L^2_\textrm{per}(D; {\mathbb R}^d)$|⁠, where |$(\pmb {\eta }, \xi ) \in {\mathbb H} \times H^1_\textrm{per}(D)/\mathbb {R}$| is a unique tuple such that
and |$\xi \in H^1_\textrm{per}(D)/\mathbb {R}$| solves the following Poisson problem (cf. Babutzka & Kunstmann, 2018):
(2.1)
In this paper we denote the Stokes operator by |$ \textbf {A}:= \textbf {P}_{{\mathbb H}} \varDelta : H^2(D; {\mathbb R}^d) \rightarrow {\mathbb H}$|⁠.
We assume that |$ \textbf {B}: L^2(\varOmega ;H^1_\textrm{per}(D; {\mathbb R}^d)) \rightarrow L^2(\varOmega ; H^1_\textrm{per}(D; {\mathbb R}^d))$| is Lipschitz continuous and has linear growth, i.e., there exists a constant |$C> 0$| such that for all |$ \textbf {v}, \textbf {w} \in L^2_\textrm{per}(D; {\mathbb R}^d)$|⁠,
(2.2a)
 
(2.2b)
 
(2.2c)
where |$\mathcal {D}\mathbf {B}$| denotes the Gateaux derivative of |$\mathbf {B}$| and |$\|\cdot \|_*$| is its operator norm.

2.2 Variational formulation of the stochastic Stokes equations

We first recall the solution concept for (1.1) and refer to Chow (2007), Da Prato & Zabczyk (1992) for its existence and uniqueness.

 

Definition 2.1
Given |$(\varOmega ,{\mathcal {F}}, \{{\mathcal {F}}_t\},{\mathbb {P}})$|⁠, let |$W$| be an |${\mathbb R}$|-valued Wiener process on it. Suppose |$ \textbf {u}_0\in L^2(\varOmega , {\mathbb V})$| and |$\mathbf {f} \in L^2(\varOmega ;L^2((0,T);L^2_\textrm{per}(D;\mathbb {R}^d)))$|⁠. An |$\{{\mathcal {F}}_t\}$|-adapted stochastic process |$\{ \textbf {u}(t) ; 0\leq t\leq T\}$| is called a variational solution of (1.1) if |$ \textbf {u} \in L^2\bigl (\varOmega ; C([0,T]; {\mathbb V})) \cap L^2\bigl (0,T;H^2_\textrm{per}(D;{\mathbb R}^d)\bigr )$| and satisfies |${\mathbb {P}}$|-a.s. for all |$t\in (0,T]$|⁠,
(2.3)

The following estimates from Carelli et al. (2012), Feng & Qiu (2019) establish the Hölder continuity in time of the variational solution in various spatial norms.

 

Theorem 2.2
Additionally suppose |$ \textbf {u}_0 \!\in\! L^2\bigl (\varOmega ; {\mathbb V} \cap H^2_\textrm{per}(D; {\mathbb R}^d)\bigr )$| and |$\mathbf {f} \!\in\! L^2(\varOmega ,C^{\frac 12}([0,T]);H^1_\textrm{per}(D;\mathbb {R}))$|⁠. There exists a constant |$C>0$| such that the variational solution to problem (1.1) satisfies for |$s,t \in [0,T]$|⁠,
(2.4a)
 
(2.4b)

 

Remark 2.3

To avoid the technicality of tracking the required ‘minimum’ assumptions on |$\mathbf {u}_0$| and |$\mathbf {f}$| for each stability and/or error estimate, unless it is stated otherwise we shall implicitly make the ‘maximum’ assumption |$ \textbf {u}_0 \in L^2\bigl (\varOmega ; {\mathbb V} \cap H^2_\textrm{per}(D; {\mathbb R}^d)\bigr )$| and |$\mathbf {f} \in L^2(\varOmega ,C^{\frac 12}([0,T]);H^1_\textrm{per}(D;\mathbb {R}))$| in the rest of the paper.

We also note that (2.4b) has only been proved for the periodic boundary condition case in the literature, which is the main reason for restricting our attention to the periodic boundary condition case in this paper as well.

2.3 Definition and role of the pressure

Definition 2.1 only addresses the velocity |$\mathbf {u}$| in the stochastic PDE (1.1); a corresponding pressure that satisfies a proper formulation (see Theorem 2.4 below) may be constructed after the existence of a velocity field |$ \textbf {u}$| has been established. We therefore consider processes
Evidently, |$\mathbf {U} \in L^2\bigl (\varOmega , L^2(0,T; H^2_\textrm{per}(D,\mathbb {R}^d))\bigr )$| and |$\mathbf {G}\in L^2\bigl (\varOmega , L^2(0,T; L^2_\textrm{per}(D,\mathbb {R}^d))\bigr )$|⁠, and (2.3) therefore implies
(2.5)
By the Helmholtz decomposition (Langa et al., 2003, Theorem 4.1 and Remark 4.3), there exists a unique |$P \in {L^2\bigl (\varOmega , L^2(0,T; H^1_\textrm{per}(D))/\mathbb {R}\bigr )}$| such that
(2.6)
in the distributional sense. It is shown in Langa et al. (2003, Section 5) that its distributional time derivative |$p:= \partial _t P \in L^1\bigl (\varOmega ; W^{-1,\infty }(0,T; H^1_\textrm{per}(D)/\mathbb {R})\bigr )$|⁠. As a consequence, we have the following result.

 

Theorem 2.4
Let |$\{ \textbf {u}(t) ; 0\leq t\leq T\}$| be the variational solution of (1.1). There exists a unique adapted process |$P\in {L^2\bigl (\varOmega , L^2(0,T; H^1_\textrm{per}(D)/\mathbb {R})\bigr)}$| such that |$(\mathbf {u}, P)$| satisfies |${\mathbb {P}}$|-a.s. for all |$t\in (0,T]$|⁠,
(2.7a)
 
(2.7b)

System (2.7) can be regarded as a mixed formulation for the stochastic Stokes system (1.1), where the (time-averaged) pressure |$P$| is defined. Below, we also define another time-averaged ‘pressure’
where we use the Helmholtz decomposition |$ \textbf {B}(\mathbf {u}(t)) = \pmb {\eta }(t) + \nabla \xi (t)$|⁠, where |$\xi \in H^1_\textrm{per}(D)/\mathbb {R}$|  |${\mathbb P}\mbox {-a.s.}$| such that
(2.8)
Then (2.6) can be rewritten as
(2.9)
The time averaged ‘pressure’ |$\{R(t); 0\leq t\leq T\}$| will also be a target process to be approximated in our numerical methods.

3. Semidiscretization in time

In this section we study the stability and convergence properties of a Helmholtz-decomposition-enhanced Euler–Maruyama time discretization scheme that is based on (1.7), where the stochastic pressure is removed from the noise term via the Helmholtz decomposition, but its |${\mathbb V}$|-valued velocity approximation |$\{ \textbf {u}^{n+1}\}_n$| still solves the original Euler–Maruyama scheme (1.3).

3.1 Formulation of the time-stepping scheme

In the following, let |$N$| be a positive integer, |$k = \frac {T}{N}$| and |$t_n= nk$| for |$n = 0, 1, \ldots , N$| be a uniform mesh that covers |$[0,T]$|⁠.

Algorithm 1

Let |$ \textbf {u}^0= \textbf {u}_0$|⁠. For |$n=0,1,\ldots , N-1$| do the following steps:

Step 1: Find |$\xi ^{n} \in L^2\bigl (\varOmega ,H^1_\textrm{per}(D)/\mathbb {R}\bigr )$| by solving
(3.1)
Step 2: Set |${\pmb {\eta }}^n:= \textbf {B}( \textbf {u}^n)-\nabla \xi ^{n}$| and find |$( \textbf {u}^{n+1},r^{n+1}) \in L^2\bigl (\varOmega , {\mathbb V} \times L^2_\textrm{per}(D)/\mathbb {R}\bigr )$| by solving
(3.2a)
 
(3.2b)

Step 3: Define |$p^{n+1}:= r^{n+1} + k^{-1} \xi ^n \varDelta _{n+1} W $|⁠.

 

Remark 3.1
By the elliptic regularity theory (see Girault & Raviart,1986, p. 13), the solution of (3.1) is in |$\xi ^{n} \in L^2\bigl (\varOmega ,H^2_\textrm{per}(D)/\mathbb {R}\bigr )$| and satisfies Lebesgue-a.e.
(3.3a)
Moreover, there exists a constant |$C>0$| such that
(3.4)

The solvability of Algorithm 1 is clear because a linear coercive elliptic PDE problem is solved at each step. Step 1 in Algorithm 1 requires a Poisson problem (3.1) to be solved, which only slightly increases the computational cost if a fast solver is used to solve them. The iterates |$\{(\mathbf {u}_n, r_n)\}_n$| and |$\{p_n\}_n$| defined in Steps 2 and 3 aim to approximate |$\{ (\mathbf {u}(t),r(t)); 0\leq t\leq T\}$| and |$\{ p(t); 0\leq t\leq T\}$|⁠, respectively. See Section 3.4 for details.

3.2 Stability estimates

In this subsection we present some stability estimates for the time-stepping scheme given in Algorithm 1. All these estimates, in particular the estimate for |$\{ \nabla r^{n+1}\}_n$|⁠, will play an important role in establishing optimal-order error estimates for the fully mixed finite element discretization to be given in the next section.

 

Lemma 3.2
Let |$\{ ( \textbf {u}^{n+1}, r^{n+1})\}_n$| be generated by Algorithm 1. There exists a constant |$C>0 $|⁠, such that
(3.5)
 
(3.6)

 

Proof.
Since a proof of estimate (3.5) can be found in Carelli & Prohl (2012, Lemma 3.1) and in Feng & Qiu (2019, Lemma 2.1), we only give a sketch of the proof for estimate (3.6) here. To prove (3.6), multiplying (1.7a) with |$\nabla r^{n+1}$| and integrating over |$D$| and using (1.7b) and |$ \textrm {div}\, \eta ^n = 0$| yields
hence,
Taking the expectation and applying the summation operator |$\sum _{n=0}^{N-1}$| on both sides gives the desired estimate. The proof is complete.

3.3 Error estimate for the velocity approximation

Since the velocity approximation |$\{ \textbf {u}^{n+1}\}_n$| generated by Algorithm 1 also solves the original Euler–Maruyama time-stepping scheme (1.3), the following optimal-order error estimate for |$\{ \textbf {u}^{n+1}\}_n$| was established in Carelli & Prohl (2012), Feng & Qiu (2019).

 

Theorem 3.3
Let |$\{ ( \textbf {u}^{n+1}, r^{n+1})\}_n$| be generated by Algorithm 1. There exists a constant |$C>0$|⁠, such that
(3.7)

We note that the proof of the above error estimate crucially uses the fact that |$\mathbf {u}^n$| is exactly divergence-free for each |$0\leq n\leq N$|⁠.

3.4 Error estimates for the pressure approximations

An optimal-order error estimate was obtained in Feng & Qiu (2019) for |$\{P(t_n)\}_n$| via the Euler–Maruyama time-stepping scheme (1.3). For the reader’s convenience, we give its proof here.

 

Theorem 3.4
Let |$\{p^n; 1\leq n\leq N\}$| be the pressure in (1.3) and |$\{P(t); 0\leq t\leq T\}$| be defined in Theorem 2.4. There exists a constant |$C>0$|⁠, such that
(3.8)

 

Proof.
Consider (1.3a) and take the sum over steps |$0 \leq n \leq m-1$|⁠. We denote |$ \textbf {U}^{m}:= k \sum _{n=0}^{m-1} \textbf {u}^{n+1}$| and |$P^m:= k \sum _{n=0}^{m-1} p^{n+1}$| and therefore obtain
(3.9)
We subtract this equation from (2.6) at time |$t = t_m$| and denote |$ \textbf {E}_{\textbf {U}}^m:= \textbf {U}(t_m) - \textbf {U}^m\in L^2(\varOmega ; H^1_\textrm{per}(D;\mathbb {R}^d))$| and |$E_{P}^m:= P(t_m) - P^m \in L^2(\varOmega ; L^2_\textrm{per}(D))$|⁠. By the stability of the divergence operator, there exists |$\beta>0$|⁠, such that
(3.10)
Taking squares on both sides and then applying expectations, Theorem 3.3 in combination with Hölder’s inequality leads to
By Ito’s isometry and (2.2a), as well as (2.4a), and Theorem 3.3, we find the bounds
and by using the Cauchy–Schwarz inequality we get
which leads to the desired estimate (3.8).

We now consider the pressure |$R^m:= k \sum _{n=0}^{m-1} r^{\,n+1}$|⁠, where |$\{r^n\}_{n}$| is defined by Algorithm 1. Using the new notation, (3.9) can be written as
(3.11)
We again subtract this equation from (2.9) at time |$t = t_m$| and adapt the error notation in (3.10),
where V is the same as above, hence, |${\mathbb E}\bigl [ \Vert {\texttt {V}}\Vert ^2\bigr ] \leq Ck$|⁠, and
By a stability result for the Poisson problems (2.8), (3.1) and property (2.2a), we easily obtain, thanks to Theorem 3.3,
Similarly, we get |${\mathbb E}\bigl [ \Vert {\texttt {VI}}_{\texttt {2}}\Vert ^2\bigr ] \leq Ck$|⁠. Therefore, we have the following corollary.

 

Corollary 3.5
Let |$\{r^n; 1\leq n\leq N\}$| be the discrete process from Algorithm 1. There exists a constant |$C>0$| such that
(3.12)

4. Fully discrete, inf-sup stable mixed finite element method

In this section we discretize Algorithm 1 in space via an inf-sup stable mixed finite element method. We choose the prototypical Taylor--Hood mixed finite element (see, e.g., Girault & Raviart, 1986; Brezzi & Fortin, 1991) as an example and give a detailed error analysis for the resulting fully discrete method, but we remark that the convergence analysis below also applies to general inf-sup stable mixed finite elements.

4.1 Preliminaries

Let |${\mathcal {T}}_h$| be a quasi-uniform triangular or rectangular mesh of |$D\subset \mathbf {R}^d$| with mesh size |$0<h \ll 1$|⁠. We define the following finite element spaces:
where |$P_\ell (K; {\mathbb R}^d)$| (⁠|$\ell =1,2$|⁠) denotes the set of |${\mathbb R}^d$|-valued polynomials of degree less than or equal to |$\ell $| over the element |$K\in {\mathcal {T}}_h$|⁠. Note that |$S_{h}=W_{h}$|⁠, which could be relaxed to |$S_{h}\subseteq W_{h}$| when higher-order Taylor--Hood elements are used.
We recall that the pair |$(\mathbb {X}_{h}, {W}_{h})$| satisfies the (discrete) inf-sup condition: there exists an |$h$|-independent constant |$\gamma>0$| such that
(4.1)
Next let |$\rho _h: L^2_\textrm{per}(D)\to W_h$| resp. |$\mathcal {R}_h: H^1_\textrm{per}(D)/\mathbb {R}\to S_h$| denote the |$L^2$|-resp. the Ritz-projection operators which are defined by
(4.2)
 
(4.3)
Then the following approximation properties are well known (cf. Girault & Raviart, 1986; Ern & Guermond, 2004; Falk, 2008):
(4.4)
 
(4.5)
for |$s=1,2$|⁠. Here, |$C$| is a positive constant independent of |$h$|⁠.
We also consider the space |${\mathbb V}_h\subset {\mathbb X}_h$| of discretely divergence-free functions,
and define the |$L^2_\textrm{per}(D; {\mathbb R}^d)$|-projection operator |$\mathbf {P}_h: L^2_\textrm{per}(D;\mathbb {R}^d) \rightarrow {\mathbb V}_h$| by
The following approximation properties are well known (cf. Heywood & Rannacher, 1982):
(4.6)
for |$s=1,2$|⁠. Here, |$C$| is again a positive constant independent of |$h$|⁠.

4.2 Formulation of the fully discrete mixed finite element method

The fully discrete, inf-sup stable finite element below is a spatial discretization of Algorithm 1. We note that since |${\mathbb V}_h \not \subset {\mathbb V}$|⁠, in general, the mixed finite element discretization requires improved stability estimates for the semidiscrete pressure |$\{ r^{n+1}\}_n$| as given in Lemma 3.2 in order to ensure optimal convergence properties.

Algorithm 2

Let |$ \textbf {u}_h^0\in L^2(\varOmega ; {\mathbb X}_h)$|⁠. For |$n=0,1,\ldots , N-1$|⁠, we do the following steps:

Step 1: Determine |${\xi ^{n}_h} \in L^2(\varOmega ; S_h)$| by solving
(4.7)
Step 2: Set |${\pmb {\eta }}^n_h:= \textbf {B}( \textbf {u}^n_h)-\nabla \xi ^{n}_h$|⁠. Find |$( \textbf {u}^{n+1}_h,r^{n+1}_h) \in L^2\bigl (\varOmega , {\mathbb V}_h \times W_h\bigr )$| by solving
(4.8a)
 
(4.8b)

Step 3: Define the |$W_h$|-valued random variable |$p^{n+1}_h = r^{n+1}_h + k^{-1} \xi ^{n}_h \varDelta _{n+1} W $|⁠.

 

Remark 4.1

Because of (4.7), we have |$({\pmb {\eta }}^n_h, \nabla \phi _h) = 0$| for all |$\phi _h \in S_h$|⁠, |${\mathbb P}$|-a.s. We also note that each of Step 1 and Step 2 solves a linear problem that is clearly well posed; in particular, the well-posedness of (4.8) is ensured by the inf-sup property (4.1) of the mixed finite element spaces |${\mathbb X}_h$| and |$W_h$|⁠.

4.3 Error estimate for the velocity approximation

The main result of this section is to prove the following optimal estimate for the velocity error |$ \textbf {u}^n- \textbf {u}^n_h$|⁠.

 

Theorem 4.2
Suppose that
Let |$\{ ( \textbf {u}^n, r^n); 1\leq n\leq N\}$| and |$\{ ( \textbf {u}^n_h, r_h^n); 1\leq n\leq N \}$| be, respectively, the solutions of Algorithms 1 and 2. Then there exists a constant |$C>0$| such that
(4.9)

 

Proof.
Define |$ \textbf {e}^n_{\textbf {u}} = \textbf {u}^n - \textbf {u}^n_h$| and |$e^n_r = r^n - r^n_h$|⁠. It is easy to check that |$\{( \textbf {e}^n_{\textbf {u}},e^n_r)\}_n$| satisfies the following error equations |${\mathbb P}$|-a.s. for all tuples |$( \textbf {v}_h, q_h) \in {\mathbb X}_h \times W_h$|⁠:
(4.10)
 
(4.11)
Now for any fixed |$\omega \in \varOmega $|⁠, setting |$ \textbf {v}_h = \mathbf {P}_h \textbf {e}^{n+1}_{\textbf {u}}(\omega ) \in {\mathbb V}_h$| in (4.10) yields
(4.12)
We now estimate each term on the left-hand side of (4.12) from below. First, by the definition of |$\mathbf {P}_h$| we get
(4.13)
Next, using again the fact that |$\mathbf {P}_h \textbf {u}^{n+1}_h = \textbf {u}^{n+1}_h$| and the Schwarz inequality, we obtain
(4.14)
where we have used (4.6) to get the last inequality.
For the next term in (4.12), using the fact that |$\mathbf {P}_h \textbf {e}^{n+1}_{\textbf {u}}$| takes values in |${\mathbb V}_h$|⁠, and estimates (4.4), (4.5) and (4.6), we get
(4.15)
Finally, we bound the only term on the right-hand side of (4.12) from above. By the independence of the increments |$\{ \varDelta _{n+1} W\}_n$|⁠, and its distribution, we get
(4.16)
and because of (2.2) and (4.6), there holds
(4.17)
To control |$\| \nabla (\xi ^n -\xi ^n_h)\|$|⁠, we recall the definitions of |$\xi ^n$| and |$\xi ^n_h$| to get
Setting |$\phi _h = \mathcal {R}_h [\xi ^n - \xi ^n_h]= (\xi ^n -\xi ^n_h) - (\xi ^n-\mathcal {R}_h \xi ^n )$|⁠, properties (2.2a) and (4.5) yield
Hence, by (3.4) in Remark 3.1, (4.6) and (2.2c) we get
(4.18)
Therefore,
(4.19)
We insert estimates (4.13)–(4.19) into (4.12), take the expectation, and apply the summation operator |$\sum _{n=0}^{m}$| for any |$0\leq m \leq N-1$| to conclude
(4.20)
where we have used estimates (3.5) and (3.6) to obtain the last inequality.
Applying the discrete Gronwall inequality to (4.20) then leads to
(4.21)

Finally, the desired estimate (4.9) follows from an application of the triangle inequality on |$ \textbf {e}^{m+1}_{\textbf {u}} = ( \textbf {u}^{m+1} -\mathbf {P}_h \textbf {u}^{m+1}) + \mathbf {P}_h \textbf {e}^{m+1}_{\textbf {u}}$| and using (4.21) and (4.6). The proof is complete.

4.4 Error estimates for the pressure approximations

In this subsection, we derive some error estimates for both |$r^n-r^n_h$| and |$p^n-p^n_h$|⁠. The argumentation parallels that in Section 3.4, and uses the inf-sup condition (4.1) in particular.

 

Theorem 4.3
Suppose that
Let |$\{ ( \textbf {u}^n, r^n); 1\leq n\leq N\}$| and |$\{ ( \textbf {u}^n_h, r_h^n); 1\leq n\leq N \}$| be, respectively, the solutions of Algorithms 1 and 2. There exists a constant |$C>0$|⁠, such that

 

Proof.
Summing (4.10) (after lowering the index by 1) over |$1\leq n\leq m \leq N $| leads to
By (4.1) we conclude (compare with (3.10))
Taking expectations after squaring both sides, using estimates (4.17), (4.18), (4.19) and Theorem 4.2 yield the desired result.

The following result now is a simply corollary of Theorem 4.3.

 

Corollary 4.4
Let |$\{p^n_h\}_n$| be the solution in Algorithm 2. Then there exists a constant |$C>0$|⁠, such that

4.5 Space-time error estimates for Algorithm 2

Theorems 3.3, 3.4, 4.2 and 4.3 and Corollaries 3.5 and 4.4 now provide the following global error estimates.

 

Theorem 4.5
Let |$( \textbf {u}, P)$| solve (1.1), and |$\{ ( \textbf {u}^n_h, r_h^n, p_h^n); 1\leq n\leq N \}$| solve Algorithm 2. There exists a constant |$C>0$| such that
for all |$1 \leq m \leq N$|⁠.

5. Stabilization methods for (1.1)

The scheme in Section 4 requires inf-sup stable pairings |$({\mathbb X}_h, W_h)$| of which the Taylor–Hood mixed finite element is one example. By recalling its definition in Section 4.1, we observe that the dimension of |${\mathbb X}_h$| exceeds that of |$W_h$|⁠. The motivation for the stabilization methods in Hughes et al. (1986) is to relax the inf-sup stability criterion for pairings of ansatz spaces in order to allow for equal-order ansatz spaces for both velocity and pressure approximates; see Brezzi & Fortin (1991), Ern & Guermond (2004), Girault & Raviart (1986), Hughes et al. (1986) for further details.

Below we replace |${\mathbb X}_h$| defined in Section 4.1 by
to which we associate the |$L^2_\textrm{per}$|-projection operator |$ \textbf {Q}_h: L^2_\textrm{per}(D; {\mathbb R}^d) \rightarrow {\mathbb Y}_{h}$| by
which satisfies the following approximation property (cf. Ern & Guermond, 2004):
(5.1)
for |$s=1,2$|⁠. Here, |$C>0$| is a constant independent of |$h$|⁠. Moreover, let |$W_h$| be the same as in Section 4, and let |$\widetilde {\mathcal R}_h$| denote the Ritz projection from |$H^1_\textrm{per}(D)/\mathbb {R}$| to |$W_{h}$|⁠. Again, we take |$S_h=W_h$| in this section.
In this section we consider the equal-order pairing |$({\mathbb Y}_{h}, W_{h})$| to discretize (1.1) based on (1.9a)–(1.9b), which violates the inf-sup condition. In fact, the following estimate is known to hold (cf. Hughes et al., 1986): there exists |$\delta>0$| independent of |$h>0$|⁠, such that
(5.2)

Equation (5.2) can be regarded as the reason why this pairing still performs optimally when applied to the Stokes problem, where |$\varepsilon = {\mathcal O}(h^2)$|⁠. Below we show that such a strategy can again be successful for the stochastic Stokes problem (1.1) if the proper pressure is chosen for the perturbation and that using the Helmholtz projection of the noise provides such an approach.

To prepare for the analysis, we start with a modification of Algorithm 1 that perturbs the incompressibility constraint.

Algorithm 3

Let |$0<\varepsilon \ll 1$| and |$ \textbf {u}_{\varepsilon }^0 = \textbf {u}_0$|⁠. For |$n=0,1,\ldots , N-1$|⁠, do the following steps:

Step 1: Find |$\xi ^{n}_{\varepsilon } \in L^2\bigl (\varOmega ,H^1_\textrm{per}(D)/\mathbb {R}\bigr )$| by solving
(5.3)
Step 2: Set |${\pmb {\eta }}^n_{\varepsilon }:= \textbf {B}( \textbf {u}^n_\varepsilon )-\nabla \xi ^{n}_\varepsilon $|⁠, and find |$( \textbf {u}^{n+1}_\varepsilon ,r^{n+1}_\varepsilon ) \in L^2\bigl (\varOmega , H^1_\textrm{per}(D; {\mathbb R}^d) \times H^1_\textrm{per}(D)/\mathbb {R}\bigr )$| by solving
(5.4a)
 
(5.4b)

Step 3: Define |$p^{n+1}_\varepsilon := r^{n+1}_\varepsilon +k^{-1} \xi ^{n}_\varepsilon \varDelta _{n+1} W$|⁠.

Because each step involves a coercive linear problem, Algorithm 3 has a unique solution. The lemma below establishes an energy estimate for the solution.

 

Lemma 5.1
Let |$\{ (u^n_\varepsilon ,r^n_\varepsilon ) \}_{n\geq 0}$| be the solution of Algorithm 3. Then there exists a constant |$C>0$|⁠, such that
(5.5)

 

Proof.
By fixing one |$\omega \in \varOmega $| and choosing |$( \textbf {v}, q ) = \bigl ( \textbf {u}^{n+1}_\varepsilon (\omega ), r^{n+1}_\varepsilon (\omega )\bigr )$| in (5.4a), we then obtain the identity
(5.6)
Taking expectations, applying the summation operator |$\sum _{n=0}^{m}$| for any |$0\leq m \leq N-1$| and using the independence of the increments |$\{\varDelta _{n+1}W\}_n$| yields
(5.7)
By (5.3) and (2.2b) we have
therefore the first term on the right-hand side of (5.7) can be bounded as
We insert these auxiliary estimates into (5.6), take expectations, sum over all time steps and use the discrete Gronwall inequality to get the desired estimate.

Note that the estimate for |$\{ \nabla r^{n+1}\}$| is scaled by |$\varepsilon>0$|⁠, which is too weak in the following to verify optimal error estimates for a spatial discretization of Algorithm 3. The following stability result therefore sharpens the estimate (5.5); its proof crucially exploits again the fact that each |${\pmb {\eta }}^n_\varepsilon $| is an |${\mathbb H}$|-valued random variable.

 

Lemma 5.2
Let |$\{ (u^n_\varepsilon ,r^n_\varepsilon ) \}_n$| be the solution of Algorithm 3. Then there exists a constant |$C>0$|⁠, such that
(5.8)

 

Proof.
Step 1: We adapt the argumentation from Carelli et al. (2012, Theorem 3.1) and interpret problem (1.9)—with |${\pmb {\eta }}^n$| being replaced by |${\pmb {\eta }}^n_\varepsilon $|—as perturbation of problem (1.7). Subtracting the corresponding equations of both systems and denoting |$ \textbf {e}^{n+1}_{\textbf {u}}:= \textbf {u}^{n+1} - \textbf {u}^{n+1}_\varepsilon $| resp. |$e^{n+1}_r:= r^{n+1} - r^{n+1}_{\varepsilon }$| yields
(5.9a)
 
(5.9b)
Now fix one |$\omega \in \varOmega $|⁠, test (5.9a) with |$ \textbf {v} = \textbf {e}_{\textbf {u}}^{n+1}(\omega )$|⁠, and (5.9b) with |$q = e^{n+1}_r(\omega )$|⁠, and afterwards sum both equations; we then conclude
(5.10)
Using Young’s inequality, hiding one part of the last term in the corresponding term on the left-hand side, using the independence of |$ \varDelta _{n+1}W$| and |$ \textbf {e}^n_{\textbf {u}}$|⁠, |$ \varDelta _{n+1}W$| as well as of |$[{\pmb {\eta }}^n - {\pmb {\eta }}^n_{\varepsilon }]$|⁠, and utilizing (2.2), we obtain
Subtracting (5.3) from (3.1) and using (2.2) we get
Hence, we then conclude from (5.10) with the help of the discrete Gronwall inequality, and (3.6) that
Consequently, by (3.6) we conclude that
Step 2: Fix one |$\omega \in \varOmega $| in (1.9a), multiply the equation by |$-\varDelta \textbf {u}^{n+1}_\varepsilon (\omega )$|⁠, integrate, perform an integration by parts on the last term and use the periodicity of |${\pmb {\eta }}^n_{\varepsilon }$| and |$ \textbf {u}^{n+1}_\varepsilon $|⁠; we get
(5.11)
To bound the last term above, we use the Schwarz inequality, the fact that |${\pmb {\eta }}^{n}_{\varepsilon } = \mathbf {B}(\mathbf {u}^{n}_{\varepsilon }) - \nabla \xi ^n_{\varepsilon }$|⁠, (2.2) and (3.4) to get
(5.12)
Substituting (5.12) into (5.11), summing over all time steps, and using (5.5) and the result of Step 1 we get the desired estimate (5.8). The proof is complete.

From Step 1 of the above proof we also obtain the following result.

 

Theorem 5.3
Let |$\{ ( \textbf {u}^{n+1}, r^{n+1})\}_n$| and |$\{ ( \textbf {u}^n_\varepsilon ,r^n_\varepsilon ) \}_n$| be the solutions of Algorithms 1 and 3, respectively. Then there exists a constant |$C>0$|⁠, such that
(5.13)

We are ready to bound the error between the pressures |$\{r^n\}_n$| and |$\{r^n_\varepsilon \}_n$|⁠; the proof of it uses (5.9) after summation in time, and follows along the lines of the proof of Theorem 3.4, using the stability of the divergence operator (cf. estimate (3.10)), and Theorem 5.3.

 

Theorem 5.4
Let |$\{r^n; 1\leq n\leq N\}$| be generated by Algorithm 1 and |$\{r^n_\varepsilon ; 1\leq n\leq N\}$| by Algorithm 3. There exists a constant |$C>0$|⁠, such that for |$m=1,2,\ldots , N$|⁠,

The above estimate and the algebraic relation in Step 3 of Algorithm 3 immediately imply the following estimate for |$p^n -p^n_\varepsilon $|⁠.

 

Corollary 5.5
Let |$\{p^n; 1\leq n\leq N\}$| be generated by Algorithm 1 and |$\{p^n_\varepsilon ; 1\leq n\leq N\}$| by Algorithm 3. There exists a constant |$C>0$|⁠, such that for |$m=1,2,\ldots , N$|⁠,

Next we present the following modification of Algorithm 2.

Algorithm 4

Let |$0<\varepsilon \ll 1$| and |$ \textbf {u}_{\varepsilon ,h}^0\in L^2(\varOmega ; {\mathbb Y}_{h})$|⁠. For |$n=0,1,\ldots , N-1$|⁠, do the following steps:

Step 1: Determine |${\xi ^{n}_{\varepsilon ,h}} \in L^2(\varOmega ; S_{h})$| from
(5.14)
Step 2: Set |${\pmb {\eta }}^n_{\varepsilon ,h}:= \textbf {B}( \textbf {u}^n_{\varepsilon ,h})-\nabla \xi ^{n}_{\varepsilon ,h}$|⁠. Find |$( \textbf {u}^{n+1}_{\varepsilon ,h},r^{n+1}_{\varepsilon ,h}) \in L^2\bigl (\varOmega , {\mathbb Y}_{h} \times W_{h} \bigr )$| by solving
(5.15a)
 
(5.15b)

Step 3: Define the |$W_{h}$|-valued random variable |$p^{n+1}_{\varepsilon ,h} = r^{n+1}_{\varepsilon ,h} + k^{-1}\xi ^{n}_{\varepsilon ,h} \varDelta _{n+1} W $|⁠.

The main result of this section is the following estimate for the velocity error |$ \textbf {u}^n_{\varepsilon } - \textbf {u}^n_{\varepsilon ,h}$|⁠.

 

Theorem 5.6
Suppose
Let |$\{ ( \textbf {u}^n_\varepsilon , r^n_\varepsilon ); 1\leq n\leq N\}$| and |$\{ ( \textbf {u}^n_{\varepsilon ,h}, r_{\varepsilon ,h}^n); 1\leq n\leq N \}$| be the solutions of Algorithms 3 and 4, respectively. Then there exists a constant |$C>0$|⁠, such that
(5.16)
This estimate suggests that |$\varepsilon = {\mathcal O}(h^2)$| is the optimal choice of |$\varepsilon $|⁠.

 

Proof.
Let |$ \textbf {e}^{n+1}_{\textbf {u}}:= \textbf {u}^{n+1}_\varepsilon - \textbf {u}^{n+1}_{\varepsilon ,h}$| and |$e^{n+1}_r:= r^{n+1}_\varepsilon - r^{n+1}_{\varepsilon ,h}$|⁠. Then |$ \{( \textbf {e}^n_{\textbf {u}},e^n_r)\}_n$| satisfies the following error equations |${\mathbb P}$|-a.s. for all tuples |$( \textbf {v}_h, q_h) \in {\mathbb Y}_{h} \times W_{h}$|⁠:
(5.17)
 
(5.18)
Now consider (5.17)–(5.18) for |$\omega \in \varOmega $| fixed, and choose
we then deduce
(5.19)
We can adopt the corresponding arguments in (4.13) and (4.14), and use Lemma 5.2 to treat the first two terms in (5.19), and also the argument around (4.16) may easily be adopted to the present setting. But a different treatment is required to deal with the last term on the left-hand side of (5.19) because it involves the error in the pressure. We rewrite this term as follows:
We estimate |${\texttt {II}}_2$| with the help of (4.5) and using Lemma 5.2,
Integrating by parts in |${\texttt {II}}_3$|⁠, using (5.1) and again Lemma 5.2 yields
Because of (5.18) we have
Putting the above auxiliary estimates together, we obtain that there exists some |$h$|- and |$\varepsilon $|-independent constant |$C>0$| such that
(5.20)
for every |$0 \leq m \leq N$|⁠. The desired estimate (5.16) then follows from an application of the discrete Gronwall inequality. The proof is complete.

The last result gives an estimate for the pressure approximation error. Using (5.2), equation (5.17) after taking a summation in |$n$|⁠, and Lemma 5.2, we obtain the following theorem.

 

Theorem 5.7
Let |$\{r^n_\varepsilon ; 1\leq n\leq N\}$| be the pressure in Algorithm 3 and |$\{r^n_{\varepsilon ,h}; 1\leq n\leq N\}$| be the pressure in Algorithm 4. There exists a constant |$C>0$|⁠, such that for all |$1 \leq m \leq N$|⁠,

The above estimate and the algebraic relation in Step 3 of Algorithm 4 immediately imply the following estimate for |$p^n_\varepsilon -p^n_{\varepsilon ,h}$|⁠.

 

Corollary 5.8
Let |$\{p^n_\varepsilon ; 1\leq n\leq N\}$| be the pressure in Algorithm 3 and |$\{p^n_{\varepsilon ,h}; 1\leq n\leq N\}$| be the pressure in Algorithm 4. There exists a constant |$C>0$|⁠, such that for all |$1 \leq m \leq N$|⁠,

To sum up the results in this section, we have shown the following error estimates for Algorithm 4.

 

Theorem 5.9
Let |$( \textbf {u}, P)$| be the solution of (1.1) and |$\{ ( \textbf {u}^n_{\varepsilon ,h}, r_{\varepsilon ,h}^n, p_{\varepsilon ,h}^n); 1\leq n\leq N \}$| be the solution of Algorithm 4. There exists a constant |$C>0$|⁠, such that

6. Computational experiments

We present computational results to validate the theoretical error estimates in Theorems 4.5 and 5.9, and evidence how crucial the numerical treatment of the pressure part in the noise is to obtain an optimally convergent mixed method for (1.1). Our computations are done using the software packages FreeFem++ (Hecht et al. 2008) and Matlab, and the physical domain of all experiments is taken to be |$D=(0,1)^2$|⁠, i.e., |$L=1$|⁠.

Specifically, we use Algorithm 2 to compute the solution of the following initial-(Dirichlet) boundary value problem:
(6.1a)
 
(6.1b)
 
(6.1c)
 
(6.1d)
and use Algorithm 4 to compute the solution of the pressure stabilization of the above system which is obtained by replacing (6.1b) by the following two conditions:
(6.2a)
 
(6.2b)
where |$\partial _{\textbf {n}} p$| stands for the normal derivative of |$p$|⁠.
Test 1. Let |$ \textbf {u}_0=(0,0), \mathbf {f} = (1,1)^\textrm{T} $| and |$ \textbf {B}(u_1,u_2) = \bigl ((u_1^2+1)^{\frac 12}, (u_2^2+1)^{\frac 12} \bigr )$|⁠, which is nonsolenoidal. We choose |$W$| in (1.1) to be an |${\mathbb R}^J$|-valued Wiener process, with increment
(6.3)
where |$ \textbf {x}=(x_1,x_2) \in D, \, \xi _{j_1j_2}^n \sim \mathcal{N}(0,1)$|⁠, |$\lambda _{j_1j_2} = \frac {1}{j_1^2 + j_2^2}$| and
(6.4)
We use the following parameters: |$J = 4$| and |$T = 1$|⁠, and take |$N_p = 501$| to be the number of realizations in this test.
Let |$k_0$| and |$k$| denote the fine and regular time-step sizes that are used to generate the numerical true solution and a computed solution, clearly |$k_0\ll k$|⁠. Moreover, |$( \textbf {u}_h^n(\tau ), p_h^n(\tau ))$| denotes the numerical solution at time step |$t_n$| using the time-step size |$\tau $|⁠; below, |$\tau =k_0$| or |$k$|⁠. For any |$1\leq n \leq N$| we use the following numerical integration formulas:
and
The definitions of |${\mathcal {E}}_{r,\textrm{av}}^N$| and |${\mathcal {E}}_{r,0}^N$| are similar.

We then implement Algorithm 2 and verify the convergence orders of the time and spatial discretizations proved in Theorem 4.5.

To generate a numerical exact solution for computing the orders of convergence, we use |$k_0=\frac {1}{600}$| and |$h_0=\frac {1}{100}$| as fine mesh sizes to compute such a solution. Then to compute the convergence order of the time discretization for the velocity, we fix |$h=\frac {1}{100}$| and then compute the numerical solution with following time mesh sizes: |$k=\frac 15,\frac {1}{10}, \frac {1}{20}, \frac {1}{40}$|⁠. The errors in the |$L^2$|-norm (⁠|${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|⁠) and |$H^1$|-norm (⁠|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|⁠) are shown in Table 1. The numerical results verify the convergence order |$\mathcal{O}(k^{\frac 12})$| that is stated in Theorem 4.5.

Table 1

Algorithm 2: time discretization errors for the velocity |$\{ \textbf {u}^n_h\}_n$|

|$k$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.162530.25558
|$1/10$|0.115210.4960.180500.5018
|$1/20$|0.081450.50020.125800.5209
|$1/40$|0.057300.50730.087580.5225
|$k$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.162530.25558
|$1/10$|0.115210.4960.180500.5018
|$1/20$|0.081450.50020.125800.5209
|$1/40$|0.057300.50730.087580.5225
Table 1

Algorithm 2: time discretization errors for the velocity |$\{ \textbf {u}^n_h\}_n$|

|$k$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.162530.25558
|$1/10$|0.115210.4960.180500.5018
|$1/20$|0.081450.50020.125800.5209
|$1/40$|0.057300.50730.087580.5225
|$k$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.162530.25558
|$1/10$|0.115210.4960.180500.5018
|$1/20$|0.081450.50020.125800.5209
|$1/40$|0.057300.50730.087580.5225

Tables 2 and 3 display, respectively, the |$L^2$|-norm errors |$({\mathcal {E}}_{\alpha ,\textrm{av}}^N)$| and |$({\mathcal {E}}_{\alpha ,0}^N)$| (⁠|$\alpha =r$| and |$p$|⁠) of the time-averaged pressure approximations using time mesh sizes |$k=\frac 15,\frac {1}{10}, \frac {1}{20}, \frac {1}{40}$|⁠. The numerical results indicate the convergence rate |$\mathcal{O}(k^{\frac 12})$| that was predicted in Theorem 4.5. We also present the standard |$L^2$|-norm errors |${\mathcal {E}}_{r,0}^N $| and |${\mathcal {E}}_{p,0}^N$| in Tables 2 and 3, respectively, for comparison purposes, for which we observe a significantly slower rate. It should be noted that our convergence theory does not conclude such convergence behavior.

Table 2

Algorithm 2: time discretization errors for the pressure |$\{r^n_h\}_n$|

|$k$||${\mathcal {E}}_{r,\textrm{av}}^N$|Order|${\mathcal {E}}_{r,0}^N $|Order
|$1/5$|0.063520.08013
|$1/10$|0.044860.50190.062310.3629
|$1/20$|0.031610.50490.048420.3639
|$1/40$|0.022190.51020.037340.3745
|$k$||${\mathcal {E}}_{r,\textrm{av}}^N$|Order|${\mathcal {E}}_{r,0}^N $|Order
|$1/5$|0.063520.08013
|$1/10$|0.044860.50190.062310.3629
|$1/20$|0.031610.50490.048420.3639
|$1/40$|0.022190.51020.037340.3745
Table 2

Algorithm 2: time discretization errors for the pressure |$\{r^n_h\}_n$|

|$k$||${\mathcal {E}}_{r,\textrm{av}}^N$|Order|${\mathcal {E}}_{r,0}^N $|Order
|$1/5$|0.063520.08013
|$1/10$|0.044860.50190.062310.3629
|$1/20$|0.031610.50490.048420.3639
|$1/40$|0.022190.51020.037340.3745
|$k$||${\mathcal {E}}_{r,\textrm{av}}^N$|Order|${\mathcal {E}}_{r,0}^N $|Order
|$1/5$|0.063520.08013
|$1/10$|0.044860.50190.062310.3629
|$1/20$|0.031610.50490.048420.3639
|$1/40$|0.022190.51020.037340.3745
Table 3

Algorithm 2: time discretization errors for the pressure approximation |$\{p^n_h\}_n$|

|$k$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.002170.0967
|$1/10$|0.001540.49470.07220.3211
|$1/20$|0.001090.49860.05790.3184
|$1/40$|0.000770.50140.04610.3288
|$k$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.002170.0967
|$1/10$|0.001540.49470.07220.3211
|$1/20$|0.001090.49860.05790.3184
|$1/40$|0.000770.50140.04610.3288
Table 3

Algorithm 2: time discretization errors for the pressure approximation |$\{p^n_h\}_n$|

|$k$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.002170.0967
|$1/10$|0.001540.49470.07220.3211
|$1/20$|0.001090.49860.05790.3184
|$1/40$|0.000770.50140.04610.3288
|$k$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.002170.0967
|$1/10$|0.001540.49470.07220.3211
|$1/20$|0.001090.49860.05790.3184
|$1/40$|0.000770.50140.04610.3288
Table 4

Algorithm 2: spatial discretization errors for the velocity approximation |$\{ \textbf {u}^n_h\}_n$|

|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.079810.50832
|$1/10$|0.040340.98440.253151.0057
|$1/20$|0.020161.00070.126620.9995
|$1/40$|0.010071.00140.063221.0021
|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.079810.50832
|$1/10$|0.040340.98440.253151.0057
|$1/20$|0.020161.00070.126620.9995
|$1/40$|0.010071.00140.063221.0021
Table 4

Algorithm 2: spatial discretization errors for the velocity approximation |$\{ \textbf {u}^n_h\}_n$|

|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.079810.50832
|$1/10$|0.040340.98440.253151.0057
|$1/20$|0.020161.00070.126620.9995
|$1/40$|0.010071.00140.063221.0021
|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$|Order|${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|Order
|$1/5$|0.079810.50832
|$1/10$|0.040340.98440.253151.0057
|$1/20$|0.020161.00070.126620.9995
|$1/40$|0.010071.00140.063221.0021

To verify the convergence rate |$\mathcal{O}(h)$| for the velocity approximation, we fix |$k = \frac {1}{200}$| and use different spatial mesh sizes |$h = \frac 15,\frac {1}{10}, \frac {1}{20}, \frac {1}{40}$| to compute the errors |${\pmb {\mathcal {E}}}_{\mathbf {u},0}^n$| and |${\pmb {\mathcal {E}}}_{\mathbf {u},1}^n$|⁠. Table 4 contains the computational results, which verify a first-order convergence rate for both, as stated in Theorem 4.5.

To verify the convergence rate for the pressure approximation, we fix |$k=\frac {1}{200}$| and use different spatial mesh sizes |$h = \frac 15,\frac {1}{10}, \frac {1}{20}, \frac {1}{40}$|⁠. Tables 5 and 6 display the error |${\mathcal {E}}_{p,\textrm{av}}^N$| of the pressure approximation. It is evident that |${\mathcal {E}}_{p,\textrm{av}}^N$| converges linearly in |$h$| as stated in Theorem 4.5. For comparison purposes, we also compute the error |${\mathcal {E}}_{p,0}^N $| and include it in Tables 5 and 6. The numerical results suggest that the error |${\mathcal {E}}_{p,0}^N $| converges with a slower rate.

Table 5

Algorithm 2: spatial discretization errors for the pressure approximation |$\{r^n_h\}_n$|

|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.042890.30972
|$1/10$|0.021450.99970.235720.3939
|$1/20$|0.010711.00220.179770.3901
|$1/40$|0.005341.00380.136200.3972
|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.042890.30972
|$1/10$|0.021450.99970.235720.3939
|$1/20$|0.010711.00220.179770.3901
|$1/40$|0.005341.00380.136200.3972
Table 5

Algorithm 2: spatial discretization errors for the pressure approximation |$\{r^n_h\}_n$|

|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.042890.30972
|$1/10$|0.021450.99970.235720.3939
|$1/20$|0.010711.00220.179770.3901
|$1/40$|0.005341.00380.136200.3972
|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.042890.30972
|$1/10$|0.021450.99970.235720.3939
|$1/20$|0.010711.00220.179770.3901
|$1/40$|0.005341.00380.136200.3972
Table 6

Algorithm 2: spatial discretization errors for the pressure |$\{p^n_h\}_n$|

|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.1276200.44524
|$1/10$|0.0681610.90480.365040.2865
|$1/20$|0.0360680.91820.297070.2973
|$1/40$|0.0192620.90490.241890.2965
|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.1276200.44524
|$1/10$|0.0681610.90480.365040.2865
|$1/20$|0.0360680.91820.297070.2973
|$1/40$|0.0192620.90490.241890.2965
Table 6

Algorithm 2: spatial discretization errors for the pressure |$\{p^n_h\}_n$|

|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.1276200.44524
|$1/10$|0.0681610.90480.365040.2865
|$1/20$|0.0360680.91820.297070.2973
|$1/40$|0.0192620.90490.241890.2965
|$h$||${\mathcal {E}}_{p,\textrm{av}}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.1276200.44524
|$1/10$|0.0681610.90480.365040.2865
|$1/20$|0.0360680.91820.297070.2973
|$1/40$|0.0192620.90490.241890.2965
Test 2. In this test, we use Algorithm 2 to solve the driven cavity problem with stochastic forcing, which is described by system (1.1a)–(1.1b) with the following nonhomogeneous boundary condition:

Let |$W$| be the Wiener process as in (6.3), and |$ \textbf {B} \equiv (1,1)^\textrm{T} $|⁠, i.e., the noise is additive. We use the following parameters in the test: |$T =1$|⁠, |$h=\frac {1}{20}$|⁠, |$k=0.01$| and the number of realizations is |$N_p = 1001$|⁠.

Figure 1 displays the expected values of the computed stochastic velocity |$ \textbf {u}^N_h$| and pressure |$p^N_h$|⁠; as expected, they behave similarly to their deterministic counterparts. On the other hand, individual realizations of the computed stochastic velocity |$ \textbf {u}^N_h$| and pressure |$p^N_h$| given in Figs 24 show quite different behaviors from their deterministic counterparts.

Test 3. In this test we study the stabilization method in Section 5. Specifically, we implement Algorithm 4 with the same function |$ \textbf {B}$| as in Test 1, and |$\{W(t); 0\leq t\leq T\}$| is chosen as an |$\mathbb {R}$|-valued Wiener process. We also add a constant forcing term |$ \textbf {f} \equiv (1,1)^\textrm{T} $| to (1.1a) in order to construct an exact solution to system (1.1). We also take |$ \textbf {u}_0 = (0,0)$|⁠, |$T = 1$|⁠, the number of realizations |$N_p = 800$| and the minimum time step |$k_0 = \frac {1}{4096}$|⁠. The computations are done on a uniform mesh of |$D$| with mesh size |$h = \frac {1}{100}$|⁠.

In order to verify the optimal convergence rate |$\mathcal{O}(h)$| of Theorem 5.6, we fix |$k = \frac {1}{256}$| and |$\varepsilon =h^2$|⁠, and then compute the numerical solutions for different values of |$h$|⁠. The standard |$L^2$|-errors |${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$| and |${\mathcal {E}}_{p,0}^N $| for the velocity and pressure approximations are presented in Table 7. The numerical results verify the first-order convergence rate for the spatial approximation of the velocity as stated in Theorem 5.6.

For comparison purposes, we also implement the ‘standard’ stabilization method, which is based on (1.2) instead of (1.9a)–(1.9b), with the same noise and parameters as above. Table 8 displays the |$L^2$|-errors |${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$| and |${\mathcal {E}}_{p,0}^N $| of the velocity and pressure approximations. The numerical results indicate that the velocity approximation is also convergent but at a slower rate. This confirms the advantages of the proposed Helmholtz decomposition enhanced stabilization method (Algorithm 4) over the ‘standard’ stabilization method.

Test 2: (a) The expected value of $\{ \textbf {u}^N_h\}_n$. (b) Level-lines of the expected value of $\{p^N_h\}_n$. (c) The streamlines of the expected
value of $\{ \textbf {u}^N_h\}_n$.
Fig. 1.

Test 2: (a) The expected value of |$\{ \textbf {u}^N_h\}_n$|⁠. (b) Level-lines of the expected value of |$\{p^N_h\}_n$|⁠. (c) The streamlines of the expected value of |$\{ \textbf {u}^N_h\}_n$|⁠.

First realization of (a) the velocity $\{ \textbf {u}^N_h\}_n$, (b) the pressure $\{p^N_h\}_n$ and (c) the streamline of $\{ \textbf {u}^N_h\}_n$.
Fig. 2.

First realization of (a) the velocity |$\{ \textbf {u}^N_h\}_n$|⁠, (b) the pressure |$\{p^N_h\}_n$| and (c) the streamline of |$\{ \textbf {u}^N_h\}_n$|⁠.

Table 7

Algorithm 4: spatial discretization errors for the velocity |$\{ \textbf {u}^n_{\varepsilon ,h}\}_n$| and pressure |$\{p^n_{\varepsilon ,h}\}_n$|

|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0183920.147406
|$1/10$|0.0090831.01780.0929130.6658
|$1/20$|0.0040951.14930.0526110.8205
|$1/40$|0.0022790.84540.0447230.2344
|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0183920.147406
|$1/10$|0.0090831.01780.0929130.6658
|$1/20$|0.0040951.14930.0526110.8205
|$1/40$|0.0022790.84540.0447230.2344
Table 7

Algorithm 4: spatial discretization errors for the velocity |$\{ \textbf {u}^n_{\varepsilon ,h}\}_n$| and pressure |$\{p^n_{\varepsilon ,h}\}_n$|

|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0183920.147406
|$1/10$|0.0090831.01780.0929130.6658
|$1/20$|0.0040951.14930.0526110.8205
|$1/40$|0.0022790.84540.0447230.2344
|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0183920.147406
|$1/10$|0.0090831.01780.0929130.6658
|$1/20$|0.0040951.14930.0526110.8205
|$1/40$|0.0022790.84540.0447230.2344
Table 8

Standard stabilization method: spatial discretization errors for the velocity |$\{ \textbf {u}^n_{\varepsilon ,h}\}_n$| and pressure |$\{p^n_{\varepsilon ,h}\}_n$|

|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0376580.735843
|$1/10$|0.0255860.55760.888352-0.2717
|$1/20$|0.0193420.40360.5798180.6155
|$1/40$|0.0114120.76110.4426910.3893
|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0376580.735843
|$1/10$|0.0255860.55760.888352-0.2717
|$1/20$|0.0193420.40360.5798180.6155
|$1/40$|0.0114120.76110.4426910.3893
Table 8

Standard stabilization method: spatial discretization errors for the velocity |$\{ \textbf {u}^n_{\varepsilon ,h}\}_n$| and pressure |$\{p^n_{\varepsilon ,h}\}_n$|

|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0376580.735843
|$1/10$|0.0255860.55760.888352-0.2717
|$1/20$|0.0193420.40360.5798180.6155
|$1/40$|0.0114120.76110.4426910.3893
|$h$||${\pmb {\mathcal {E}}}_{\mathbf {u},0}^N$|Order|${\mathcal {E}}_{p,0}^N $|Order
|$1/5$|0.0376580.735843
|$1/10$|0.0255860.55760.888352-0.2717
|$1/20$|0.0193420.40360.5798180.6155
|$1/40$|0.0114120.76110.4426910.3893
Second realization of (a) the velocity $\{ \textbf {u}^N_h\}_n$, (b) the pressure $\{p^N_h\}_n$ and (c) the streamline of $\{ \textbf {u}^N_h\}_n$.
Fig. 3.

Second realization of (a) the velocity |$\{ \textbf {u}^N_h\}_n$|⁠, (b) the pressure |$\{p^N_h\}_n$| and (c) the streamline of |$\{ \textbf {u}^N_h\}_n$|⁠.

Third realization of (a) the velocity $\{ \textbf {u}^N_h\}_n$, (b) the pressure $\{p^N_h\}_n$ and (c) the streamline of $\{ \textbf {u}^N_h\}_n$.
Fig. 4.

Third realization of (a) the velocity |$\{ \textbf {u}^N_h\}_n$|⁠, (b) the pressure |$\{p^N_h\}_n$| and (c) the streamline of |$\{ \textbf {u}^N_h\}_n$|⁠.

Acknowledgements

After this paper was finished, reference Breit & Dodgson (2019) was brought to our attention by Prof. D. Breit. The finite element method proposed in Breit & Dodgson (2019) is essentially equivalent to Algorithm 2 of this paper although they are different algorithmically.

Funding

National Science Foundation (DMS-1620168 and DMS-2012414 to X.F.; DMS-1620168 and DMS-2012414 to L.V.).

Dedicated to the memory of John W. Barrett

References

Babutzka
,
J.
&
Kunstmann
,
P. C.
(
2018
)
|${L}^p$|-Helmholtz decomposition on periodic domains and applications to Navier-Stokes equations
.
J. Math. Fluid Mech.
,
20
,
1093–1121
.

Breit
,
D.
&
Dodgson
,
A.
(
2019
)
Convergence rates for the numerical approximation of the 2D stochastic Navier-Stokes equations
.
arXiv:1906.11778v2 [math.NA]
.

Bensoussan
,
A.
(
1995
)
Stochastic Navier–Stokes equations
.
Acta Appl. Math.
,
38
,
267
304
.

Brezzi
,
F.
&
Fortin
,
M.
(
1991
)
Mixed and Hybrid Finite Element Methods
.
New York
:
Springer
.

Carelli
 
E.
,
Hausenblas
 
E.
&
Prohl
 
A.
(
2012
)
Time-splitting methods to solve the stochastic incompressible Stokes equations
.
SIAM J. Numer. Anal.
,
50
,
2917
2939
.

Carelli
 
E.
&
Prohl
 
A.
(
2012
)
Rates of convergence for discretizations of the stochastic incompressible Navier–Stokes equations
.
SIAM J. Numer. Anal.
,
50
,
2467
2496
.

Chow
,
P.-L.
(
2007
)
Stochastic Partial Differential Equations
.
2nd edn
.
Boca Raton, Florida
:
Chapman and Hall/CRC Press
.

Da Prato
,
G.
&
Zabczyk
,
J.
(
1992
)
Stochastic Equations in Infinite Dimensions
.
Cambridge, UK
:
Cambridge University Press
.

Ern
,
A.
&
Guermond
,
J.-L.
(
2004
)
Theory and Practice of Finite Elements
.
New York
:
Springer
.

Falk
,
R.
(
2008
)
A Fortin operator for two-dimensional Taylor-Hood elements
.
ESAIM: Math. Model. Num. Anal.
,
42
,
411
424
.

Feng
 
X.
, &
Qiu
 
H.
(
2019
)
Analysis of fully discrete mixed finite element methods for time-dependent stochastic Stokes equations with multiplicative noise
.
arXiv:1905.03289v2 [math.NA]
.

Girault
,
V.
&
Raviart
,
P.-A.
(
1986
)
Finite Element Methods for Navier–Stokes Equations
.
New York
:
Springer
.

Hairer
,
M.
&
Mattingly
,
J. C.
(
2006
)
Ergodicity of the 2D Navier–Stokes equations with degenerate stochastic forcing
.
Ann. of Math.
,
164
,
993
1032
.

Hecht
 
F.
,
LeHyaric
 
A.
, &
Pironneau
 
O.
 
2008
 
Freefem++ version 2.24–1
, available at http://www.freefem.org/ff++.

Heywood
,
J. G.
&
Rannacher
,
R.
(
1982
)
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
.

Hughes
,
T. J. R.
,
Franca
,
L. P.
&
Balestra
,
M.
(
1986
)
A new finite element formulation for computational fluid mechanics: V. Circumventing the Babuska–Brezzi condition: A stable Petrov–Galerkin formulation of the Stokes problem accomodating equal order interpolation
.
Comp. Meth. Appl. Mech. Eng.
,
59
,
85
99
.

Langa
,
J. A.
,
Real
,
J.
&
Simon
,
J.
(
2003
)
Existence and regularity of the pressure for the stochastic Navier–Stokes equations
.
Appl. Math. Optim.
,
48
,
195
210
.

Prohl
,
A.
(
1997
)
Projection and Quasi-compressibility Methods for Solving the Incompressible Navier–Stokes Equations
.
Stuttgart
:
B.G. Teubner
.

This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://dbpia.nl.go.kr/journals/pages/open_access/funder_policies/chorus/standard_publication_model)