Time error estimators

for the Chorin-Temam scheme

Nizar Kharrat* — Zoubida Mghazli**

* Département de Mathématiques
(ENIT-LAMSIN), Faculté des Sciences de Bizerte
Tunis
Tunisie
nizar.kharrat@enit.rnu.tn

** Universit Ibn Tofail, Faculté des Sciences (LIRNE-EIMA)
B.P. 133 Kenitra
Maroc
mghazli_zoubida@yahoo.com

**ABSTRACT.** The time-dependent Stokes equations are discretized by the original Chorin's projection method [5] and Temam[15]. According to an idea of [1], we derive time error estimators for velocity and pressure. In particular, the velocity estimator is implemented for adaptation on the time step.

**RÉSUMÉ.** Les équations de Stokes instationnaires sont discrétisées par la méthode de projection classique de Chorin [5] et Temam[15]. En se basant sur une idée de [1], nous construisons des estimateurs sur l'erreur de discrétisation en temps pour la vitesse et la pression. En particulier, l'estimateur associée à la vitesse est mis en œuvre pour l'adaptation sur le pas de temps.

**KEYWORDS** : time-dependent Stokes equations, Chorin's projection method, time error estimators.

**MOTS-CLÉS** : équations de Stokes instationnaires, méthode de projection de Chorin, estimateurs en temps.
# 1. Introduction

For the last forty years, a large number of works have been devoted to the analysis of the numerical difficulty related to the incompressibility constraint in the time-dependent Navier-Stokes equations<fnref n="1" />. The Chorin [5] and Temam [15] scheme, started in the late 1960's, was introduced to overcome such a difficulty. However, this scheme, commonly refered as the projection scheme, suffers from the lack of accuracy on the pressure. The conjecture given by Rannacher [14] and the numerical experiments realized in [13], show in fact, the existence of boundary layers on the domain boundary.

In this paper, we are concerned with the a posteriori analysis of this scheme and we particularly focus on its time error analysis. Our aim is to provide tools to control the time step size. In this approach, we study the discretization error of the time-dependent Stokes problem in two or three dimensionnal bounded domain and for a finite-time interval.

Based on the idea of [1], the methodology that we follow here, is rather similar to [2] for the backward Euler scheme. Its main drawback consists in uncoupling as far as possible the time and the space discretization errors. To this end, we use the space variational formulation of the continuous Stokes problem and that of the two semi-discrete problems defining the projection scheme, namely the prediction step and the projection step. Furthermore, we consider the time continuous Galerkin method to approximate the velocity. For the present case, we use two continuous affine velocities respectively associated to the prediction step and to the projection step. For the pressure, we consider a piecewise constant approach. It should be noted that the splitting feature of the projection method leads to distinguish two residual error estimators. They are respectively associated to the velocity and the pressure. Furthermore, they depend on the fully discrete solution and both of them are, local in time and global with respect to the space variables. For the spatial discretization, we simply use the conforming finite element method.

Then, we prove, up to some terms involving the data and the spatial discretization error, a global upper bound of the error by the Hilbertian sum of these estimators and also a local lower bound of the error by each of them. In particular, we notice that the upper estimate on the pressure estimator is derived independently from the one obtained for the velocity estimator. Moreover, these estimates are obtained with constants independent of any discretization parameter.

Finally, considering a strategy different from the one presented in [1], we implement a simple procedure justifying the efficiency of the estimators for time adaptivity. In particular, we will be concerned only with the velocity estimator.

# 2. The continuous and semi-discrete problems

Let $\Omega$ be a bounded connected domain of $\mathbb{R}^d$ ($d = 2, 3$), with a Lipschitz continuous boundary $\Gamma$. For a given positive real $T$, we consider the time-dependent Stokes problem in the primitive variables

$$ \left\{ \begin{array}{lll} \partial_t \mathbf{u} - \nu \Delta \mathbf{u} + \nabla p = \mathbf{f} & \text{in} & \Omega \times ]0, T], \\ \nabla \cdot \mathbf{u} = 0 & \text{in} & \Omega \times ]0, T], \\ \mathbf{u} = \mathbf{0} & \text{on} & \Gamma \times ]0, T], \\ \mathbf{u}|_{t=0} = \mathbf{u}_0 & \text{in} & \Omega. \end{array} \right. \qquad (1) $$

The unknowns are the velocity $\mathbf{u} = \mathbf{u}(x, t)$ and the pressure $p = p(x, t)$ of the fluid; the data are the density of body forces $\mathbf{f} = \mathbf{f}(x, t)$ and the initial velocity $\mathbf{u}_0 = \mathbf{u}_0(x)$. The kinematic

<footnote n="1">1. for an interesting overview, we refer to [8]</footnote>
viscosity $\nu$ is assumed to be a positive constant. For simplicity, we consider a homogeneous Dirichlet condition.

**Preliminaries:** Let us denote by $\mathfrak{B}$ any separable Banach space and $[a, b]$ any time interval included in $[0, T]$. In the sequel, we denote by $L^2(a, b; \mathfrak{B})$, the space of measurable functions

$v$ from $]a, b[$ in $\mathfrak{B}$ such that $\|v\|_{L^2(a,b;\mathfrak{B})} = \left( \int_a^b \|v(\cdot, s)\|_{\mathfrak{B}}^2 ds \right)^{\frac{1}{2}} < +\infty$, and $C^0(a, b; \mathfrak{B})$

the space of continuous functions from $[a, b]$ in $\mathfrak{B}$. For $k = 0, 1, 2$, we use the Sobolev spaces $H^k(\Omega)$, equipped with their standard norms $\| \cdot \|_k$ and semi-norms $| \cdot |_k$ (the same notation is used for the vector valued functions). As usual, $H^0(\Omega) \equiv L^2(\Omega)$, $L_0^2(\Omega)$ stands for the space of functions in $L^2(\Omega)$ with zero mean value on $\Omega$ and, $H_0^1(\Omega)$ the subspace of $H^1(\Omega)$ with vanishing traces on $\Gamma$. We make use of the subspaces: $H_0(\mathrm{div}, \Omega) = \{v \in L^2(\Omega)^d; \nabla \cdot v \in L^2(\Omega) \text{ with } v|_{\Gamma} \cdot n = 0\}$, the kernels, $V = \{v \in H_0^1(\Omega)^d; \nabla \cdot v = 0\}$ and $H = \{v \in L^2(\Omega)^d; \nabla \cdot v = 0 \text{ with } v|_{\Gamma} \cdot n = 0\}$, where $n$ is the unit normal on $\Gamma$ pointing out of $\Omega$. We also consider the operator $P_H$ as the orthogonal projection from $L^2(\Omega)^d$ onto $H$. This operator plays a key role in the a priori analysis of projection schemes (see for instance [9] and the references therein).

In particular, if $\mathbf{A}$ denotes the Stokes operator [6, Chap. XIX], we have $\mathbf{A}v = -P_H\Delta v, \forall v \in V \cap H^2(\Omega)^d$. To simplify the presentation, we set:

$$X = H_0^1(\Omega)^d, \quad X' = H^{-1}(\Omega)^d, \quad Y = L^2(\Omega)^d \quad \text{and} \quad M = H^1(\Omega) \cap L_0^2(\Omega).$$

The problem (1) admits the variational formulation: find $\mathbf{u}$ in $L^2(0, T; X) \cap C^0(0, T; Y)$ and $p$ in $L^2(0, T; L_0^2(\Omega))$, such that

$$\mathbf{u}(\cdot, 0) = \mathbf{u}_0 \quad \text{a.e. in } \Omega, \qquad (2)$$

and that, for a.e. $t \in ]0, T]$ and for all $(\mathbf{v}, q) \in X \times L_0^2(\Omega)$,

$$\left\{ \begin{array}{l} \langle \partial_t \mathbf{u}, \mathbf{v} \rangle + \nu(\nabla \mathbf{u}, \nabla \mathbf{v}) - (p, \nabla \cdot \mathbf{v}) = \langle \mathbf{f}, \mathbf{v} \rangle, \\ -(\nabla \cdot \mathbf{u}, q) = 0. \end{array} \right. \qquad (3)$$

Moreover, for any $t$ in $]0, T]$ and for all $\mathbf{v}$ in $L^2(0, t; X) \cap C^0(0, t; Y)$, it is useful to define the energy norm $[\mathbf{v}](t) = \left( \|\mathbf{v}(\cdot, t)\|_0^2 + \nu \int_0^t |\mathbf{v}(\cdot, s)|_1^2 ds \right)^{1/2}$. Then, we recall the following stability result [2, Prop. 2.1]:

**Proposition 1.** For any data $(\mathbf{u}_0, \mathbf{f})$ in $\mathbf{H} \times L^2(0, T; X')$, problem (2)-(3) has a unique solution $(\mathbf{u}, p)$, which satisfies for all $t$ in $]0, T]$,

$$[\mathbf{u}](t) \leq \left( \nu^{-1} \|\mathbf{f}\|_{L^2(0,t;X')}^2 + \|\mathbf{u}_0\|_0^2 \right)^{\frac{1}{2}}.$$

Moreover, this solution is such that $\partial_t \mathbf{u} + \nabla p$ belongs to $L^2(0, T; X')$ and satisfies for all $t$ in $]0, T]$,

$$\|\partial_t \mathbf{u} + \nabla p\|_{L^2(0,t;X')} \leq 2 \left( \|\mathbf{f}\|_{L^2(0,t;X')}^2 + \frac{\nu}{2} \|\mathbf{u}_0\|_0^2 \right)^{\frac{1}{2}}.$$

For the sequel, we assume the data $(\mathbf{u}_0, \mathbf{f})$ belong to $\mathbf{H} \times C^0(0, T; X')$. Finally, for all $(\mathbf{v}, \mathbf{w})$ in $Y^2$ and $(a, b)$ in $\mathbb{R}^2$, we make use of the following properties:

$$2(\mathbf{v}, \mathbf{v} - \mathbf{w}) = \|\mathbf{v}\|_0^2 - \|\mathbf{w}\|_0^2 + \|\mathbf{v} - \mathbf{w}\|_0^2, \qquad (4)$$

$$2ab \leq \alpha a^2 + b^2/\alpha \quad (\forall \alpha > 0). \qquad (5)$$
**The projection scheme** Let $N$ be a given integer and $0 = t_0 < t_1 < \dots < t_N = T$ a partition of the interval $[0, T]$ with step sizes $\tau_n = t_n - t_{n-1}$. We denote by $\tau$, the $N$-uple $(\tau_1, \dots, \tau_N)$ and we set $|\tau| = \max_{1 \le n \le N} \tau_n$. Moreover, we assume that partition is regular in the sense of [3, Def. 1.2, Chap. VIII], i.e, there exists a constant $\sigma > 1$, such that

$$ \max_{2 \le n \le N} \frac{\tau_n}{\tau_{n-1}} \le \sigma. \qquad (6) $$

Then, for each $n$, $0 \le n \le N$, such that $\mathbf{u}^0 = \mathbf{u}(\cdot, 0)$, the projection scheme uncouples every iteration step $n$, in two substeps.

At the first step, given $\mathbf{u}^{n-1}$, we seek for a provisional velocity $\tilde{\mathbf{u}}^n$, such that :

$$ \left\{ \begin{array}{ll} \frac{\tilde{\mathbf{u}}^n - \mathbf{u}^{n-1}}{\tau_n} - \nu \Delta \tilde{\mathbf{u}}^n = \mathbf{f}^n & \text{in } \Omega, \\ \tilde{\mathbf{u}}^n = \mathbf{0} & \text{on } \Gamma, \end{array} \right. \qquad (7) $$

where $\mathbf{f}^n$ denotes an approximation to the distribution $\mathbf{f}(\cdot, t_n)$. In the second step, we look for a divergence-free velocity $\mathbf{u}^n$ and a pressure $\Phi^n$, satisfying the equations :

$$ \left\{ \begin{array}{ll} \frac{\mathbf{u}^n - \tilde{\mathbf{u}}^n}{\tau_n} + \nabla \Phi^n = \mathbf{0} & \text{in } \Omega, \\ \nabla \cdot \mathbf{u}^n = 0 & \text{in } \Omega, \\ \mathbf{u}^n \cdot \mathbf{n} = 0 & \text{on } \Gamma. \end{array} \right. \qquad (8) $$

The step (7) is nothing more than an approximation of the viscous part of the Stokes equations, while the step (8) is associated to the incompressibility constraint. The algorithm (7)-(8) admits the following variational formulation :

Find $(\tilde{\mathbf{u}}^n)_{1 \le n \le N}$ in $X^N$ and $(\mathbf{u}^0, (\mathbf{u}^n, \Phi^n)_{1 \le n \le N})$ in $\mathbf{H} \times \mathbf{H}_0(\mathrm{div}, \Omega)^N \times L_0^2(\Omega)^N$, such that

$$ \mathbf{u}^0 = \mathbf{u}(0) \quad \text{a.e. in } \Omega, \qquad (9) $$

and that, for every $n$, $1 \le n \le N$, and for all $(\mathbf{w}, \mathbf{v}, q)$ in $X \times \mathbf{H}_0(\mathrm{div}, \Omega) \times L_0^2(\Omega)$,

$$ (\tilde{\mathbf{u}}^n, \mathbf{w}) + \nu \tau_n (\nabla \tilde{\mathbf{u}}^n, \nabla \mathbf{w}) = (\mathbf{u}^{n-1}, \mathbf{w}) + \tau_n \langle \mathbf{f}^n, \mathbf{w} \rangle, \qquad (10) $$

$$ \left\{ \begin{array}{l} (\mathbf{u}^n, \mathbf{v}) - \tau_n(\nabla \cdot \mathbf{v}, \Phi^n) = (\tilde{\mathbf{u}}^n, \mathbf{v}), \\ -(q, \nabla \cdot \mathbf{u}^n) = 0. \end{array} \right. \qquad (11) $$

The latter step has also the following mixed formulation :

Find $(\mathbf{u}^0, (\mathbf{u}^n, \Phi^n)_{1 \le n \le N})$ in $\mathbf{H} \times \mathbf{H}^N \times M^N$, such that for every $n$, $1 \le n \le N$, and for all $(\mathbf{v}, q)$ in $Y \times M$,

$$ \left\{ \begin{array}{l} (\mathbf{u}^n, \mathbf{v}) + \tau_n(\mathbf{v}, \nabla \Phi^n) = (\tilde{\mathbf{u}}^n, \mathbf{v}), \\ (\mathbf{u}^n, \nabla q) = 0. \end{array} \right. \qquad (12) $$

In particular, (11) is equivalent to a Poisson-Neumann problem with a unknown pressure given by :

$$ \Delta \Phi^n = \frac{1}{\tau_n} \nabla \cdot \tilde{\mathbf{u}}^n \quad \text{in } \Omega, \quad \nabla \Phi^n \cdot \mathbf{n} = 0 \quad \text{on } \Gamma, \qquad (13) $$

and $\mathbf{u}^n$ is given by

$$ \mathbf{u}^n = \tilde{\mathbf{u}}^n - \tau_n \nabla \Phi^n. \qquad (14) $$
In this case, $\Phi^n$ must satisfy a homogeneous Neumann condition on $\Gamma$, which is not necessarily satisfied by the exact pressure. In practice, such a condition often generates boundary layers on $\Gamma$ that exponentially decay in the interior of $\Omega$. Indeed, it is conjectured in [13, 14], that the $L^2$-pressure error is first order in time, for all subdomains strictly included in $\Omega$. With standard arguments, we also prove the following stability result:

**Proposition 2.** There exists a unique solution $(\tilde{\mathbf{u}}^n, \mathbf{u}^n, \Phi^n)_{1\le n\le N}$ of problems (10) and (11) or (12), such that for all $m, 1 \le m \le N$, the following estimates hold:

$$ \begin{aligned} & \|\tilde{\mathbf{u}}^m\|_0^2 + \nu \sum_{n=1}^m \tau_n |\tilde{\mathbf{u}}^n|_1^2 + \sum_{n=1}^m \|\tilde{\mathbf{u}}^n - \mathbf{u}^{n-1}\|_0^2 + \sum_{n=1}^m \|\tilde{\mathbf{u}}^{n-1} - \mathbf{u}^{n-1}\|_0^2 \\ & \le \|\mathbf{u}_0\|_0^2 + \sum_{n=1}^m \frac{\tau_n}{\nu} \|\mathbf{f}^n\|_{X'}^2. \end{aligned} \quad (15) $$

$$ \left( \sum_{n=1}^{m} \tau_n \left\| \frac{\mathbf{u}^n - \mathbf{u}^{n-1}}{\tau_n} + \nabla \Phi^n \right\|_{X'}^2 \right)^{\frac{1}{2}} \le 2 \left( \sum_{n=1}^{m} \tau_n \|\mathbf{f}^n\|_{X'}^2 + \frac{\nu}{2} \|\mathbf{u}_0\|_0^2 \right)^{\frac{1}{2}}. \quad (16) $$

### 3. The upper and the lower bounds of the error

In this section, we introduce two distinct families of time error estimators, defined at each time step as a function of the predicted velocity and of the pressure respectively. We denote by $h > 0$, the mesh size associated to the finite element method and we denote by $(\tilde{\mathbf{u}}_h^n, \mathbf{u}_h^n, \Phi_h^n)_{1\le n\le N}$, the numerical solution obtained by the spatial discretization of the projection algorithm. In addition, with any sequence $(\psi^n)_{0\le n\le N}$ of a given Hilbert space, we associate the function $\psi_\tau$ on $[0, T]$ which is piecewise affine and continuous on the time intervals $[t_{n-1}, t_n]$, $1 \le n \le N$, defined by: $\psi_\tau(t) = \psi^n - \frac{t_n-t}{\tau_n}(\psi^n - \psi^{n-1})$, for $t \in [t_{n-1}, t_n]$. Finally, we introduce the operator $\pi_\tau$ such that, for any function $\psi$ continuous from $[0, T]$ into any Banach space, $\pi_\tau \psi$ denotes the step function which is constant and equal to $\psi(t_n)$ on each interval $[t_{n-1}, t_n]$, $1 \le n \le N$. In order to derive the a posteriori error estimates, we need some regularity on $\Omega$ and the data $(\mathbf{u}_0, \mathbf{f})$.

(R) For any data $g$ in $Y$, the Stokes problem: $\mathbf{A} \mathbf{w} = \mathbf{g}$ on $\Omega$, $\mathbf{w} = 0$ on $\Gamma$, has a unique solution $\mathbf{w}$ in $\mathbf{V} \cap H^2(\Omega)^d$, which satisfies $\|\mathbf{w}\|_2 \le C(\Omega) \|\mathbf{A} \mathbf{w}\|_0$, where $C(\Omega)$ is a constant which depends on $\Omega$.

(D) $\mathbf{u}_0 \in \mathbf{V}$ and $\mathbf{f} \in L^2(0, T; Y) \cap C^0(0, T; X')$.

In this case (cf. [16]), the problem (1) has a unique solution, such that:

$$ \mathbf{u} \in L^2(0, T; H^2(\Omega)^d) \cap C^0(0, T; \mathbf{V}), \quad \partial_t \mathbf{u} \in L^2(0, T; \mathbf{H}) \text{ and } p \in L^2(0, T; M). \quad (17) $$

Besides, we notice from (9), (10) and (12) that:

$$ \mathbf{u}_{\tau}(0) = \mathbf{u}_0 \quad \text{a.e. in } \Omega, \quad (18) $$

and that, for $1 \le n \le N$, for all $t \in [t_{n-1}, t_n]$, for all $(\mathbf{w}, \mathbf{v}, q) \in X \times Y \times M$,

$$ (\partial_t \mathbf{u}_\tau, \mathbf{w}) + \nu(\nabla \tilde{\mathbf{u}}_\tau, \nabla \mathbf{w}) - (\nabla \cdot \mathbf{w}, \Phi^n) = \langle \mathbf{f}^n, \mathbf{w} \rangle - \nu(\nabla(\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_\tau), \nabla \mathbf{w}), \quad (19) $$

and from (11),

$$ \left\{ \begin{array}{l} (\mathbf{u}_\tau - \tilde{\mathbf{u}}_\tau, \mathbf{v}) + \tau_n(\mathbf{v}, \nabla \Phi_\tau^*) = 0, \\ (\mathbf{u}_\tau, \nabla q) = 0, \end{array} \right. \quad (20) $$
where, $\Phi_{\tau}^{*}$ denotes the function :

$$ \Phi_{\tau}^{*}(t) = \frac{\tau_{n-1}}{\tau_{n}} \Phi^{n-1} + \frac{t - t_{n-1}}{\tau_{n}} \left( \Phi^{n} - \frac{\tau_{n-1}}{\tau_{n}} \Phi^{n-1} \right). $$

We notice that, $\Phi_{\tau}^{*}$ defines an affine function associated with $(\Phi^{n})_{1\le n\le N}$, but which is discontinuous at $t_{n-1}$, for $1 \le n \le N$, since for all $t \in ]t_{n-1}, t_{n}]$,

$$ (\Phi_{\tau} - \Phi_{\tau}^{*})(t) = \frac{t_n - t}{\tau_n} \left(1 - \frac{\tau_{n-1}}{\tau_n}\right) \Phi^{n-1}. $$

We also derive from (12), for all $t \in ]t_{n-1}, t_{n}]$ and for all $\mathbf{v}$ in $Y$, the important relation :

$$ (\partial_t \mathbf{u}_\tau, \mathbf{v}) + (\mathbf{v}, \nabla \Phi^n) = (\partial_t \tilde{\mathbf{u}}_\tau, \mathbf{v}) + (\mathbf{v}, \frac{\tau_{n-1}}{\tau_n} \nabla \Phi^{n-1}). \quad (21) $$

Then, by combining (2)-(3) with (18)-(19)-(20) and (21), we deduce that

$$ (\mathbf{u} - \tilde{\mathbf{u}}_{\tau}, \mathbf{u} - \mathbf{u}_{\tau}, p - \pi_{\tau} \Phi_{\tau}) \text{ satisfies } (\mathbf{u} - \mathbf{u}_{\tau})(0) = \mathbf{0} \quad \text{a.e. in } \Omega, $$

and that, for $1 \le n \le N$, for a.e. $t \in ]t_{n-1}, t_{n}]$ and for all $(\mathbf{w}, \mathbf{v}, q) \in X \times Y \times M$,

$$ \langle \partial_t (\mathbf{u} - \mathbf{u}_\tau), \mathbf{w} \rangle + \nu (\nabla (\mathbf{u} - \tilde{\mathbf{u}}_\tau), \nabla \mathbf{w}) - (p - \pi_\tau \Phi_\tau, \nabla \cdot \mathbf{w}) \\ = \langle \mathbf{f} - \pi_\tau \mathbf{f}, \mathbf{w} \rangle + \nu (\nabla (\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_\tau), \nabla \mathbf{w}), \quad (22) $$

$$ \left\{ \begin{array}{l} (\mathbf{u} - \tilde{\mathbf{u}}_{\tau}, \mathbf{v}) + \tau_n (\mathbf{v}, \nabla \Phi_{\tau}^{*}) = (\mathbf{u} - \mathbf{u}_{\tau}, \mathbf{v}), \\ (\mathbf{u} - \mathbf{u}_{\tau}, \nabla q) = 0, \end{array} \right. \quad (23) $$

and moreover,

$$ \langle \partial_t (\mathbf{u} - \mathbf{u}_\tau), \mathbf{v} \rangle - (\nabla \cdot \mathbf{v}), p - \pi_\tau \Phi_\tau = \langle \partial_t (\mathbf{u} - \tilde{\mathbf{u}}_\tau), \mathbf{v} \rangle - (\nabla \cdot \mathbf{v}), p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}. \quad (24) $$

For each $n$, $1 \le n \le N$, we define the error estimators by

$$ \tilde{\zeta}_n = \left( \nu \frac{\tau_n}{3} \right)^{\frac{1}{2}} \left| \tilde{\mathbf{u}}_h^n - \tilde{\mathbf{u}}_h^{n-1} \right|_1, \quad \zeta_n = \frac{1}{\sqrt{3}} \left| \tau_n \Phi_h^n - \tau_{n-1} \Phi_h^{n-1} \right|_1 \quad (25) $$

and their Hilbertian sum by,

$$ \tilde{\zeta}_{n\tau} = \left\{ \sum_{m=1}^{n} \tilde{\zeta}_{m}^{2} \right\}^{\frac{1}{2}}, \quad \zeta_{n\tau} = \left\{ \sum_{m=1}^{n} \zeta_{m}^{2} \right\}^{\frac{1}{2}}. \quad (26) $$

For convenience, we consider $\tilde{e}^n = \tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_h^n$, $e^n = \mathbf{u}^n - \mathbf{u}_h^n$ and $\varepsilon^n = \Phi^n - \Phi_h^n$. By following a similar strategy to that of [2] (see also [1] for the heat equation), we prove in the result below, a global upper bound of the error by the quantities (26) and, a local upper bound by the estimators (25).

**Theorem 1.** Assume the conditions (R) and (D) are satisfied. Then, for each $m \in [1, N]$, there exists $\delta_m \in [0, \frac{1}{4}[ \text{ such that, for any } (\mu_1, \mu_2) \in ]0, 1]^2, \text{ we have the following a posteriori estimate of the error between the solution } \mathbf{u} \text{ of the problem (1) and the functions } (\tilde{\mathbf{u}}_\tau, \mathbf{u}_\tau, \pi_\tau \Phi_\tau)$
associated with the solution $(\mathbf{u}^0, (\tilde{\mathbf{u}}^n, \mathbf{u}^n, \Phi^n)_{1\le n\le N})$ of the projection scheme (7)-(8), given by:

$$ 
\begin{aligned}
& \left( \mathbf{E}_{u\tau}^2(t_m) + \frac{3}{4} \sum_{n=1}^{m} \mathbf{E}_p^2(t_n, \mu_1, \mu_2) \right)^{\frac{1}{2}} \\
& + (2\nu)^{-\frac{1}{2}} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{L^2(0,t_m;X')} \\
& \le 2\sqrt{6} \left( \tilde{\zeta}_{m\tau}^2 + \gamma \zeta_{m\tau}^2 \right)^{\frac{1}{2}} \\
& \quad + 2\sqrt{\frac{3}{\nu}} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{L^2(0,t_m;X')} + (\sqrt{2} + \sqrt{3}) |\tau|^{\frac{1}{2}} C_m(\mathbf{f}, \mathbf{u}_0) \\
& \quad + 4 \left( \sum_{n=1}^{m} \left\{ \nu \tau_n |\tilde{\mathbf{e}}^n - \tilde{\mathbf{e}}^{n-1}|_1^2 + \frac{\gamma}{2} |\tau_n \varepsilon^n - \tau_{n-1} \varepsilon^{n-1}|_1^2 \right\} \right)^{\frac{1}{2}},
\end{aligned} 
\quad (27) $$

where, we denoted by:

$$ \mathbf{E}_{u\tau}(t_m) = \left\{ \|\mathbf{u}(t_m) - \mathbf{u}^m\|_0^2 + \nu \int_0^{t_m} |\mathbf{u} - \tilde{\mathbf{u}}_\tau|^2 dt + \frac{1}{2} [\mathbf{u} - \tilde{\mathbf{u}}_\tau]^2(t_m) \right\}^{\frac{1}{2}}, $$

$$ \mathbf{E}_p(t_n, \mu_1, \mu_2) = \left\{ \tau_n \int_{t_{n-1}}^{t_n} \left( (1-\mu_1)|p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}|_1^2 + (1-\mu_2)|p - \pi_\tau \Phi_\tau|_1^2 \right) dt \right\}^{\frac{1}{2}}, \quad C_m^2(\mathbf{f}, \mathbf{u}_0) = (1+C(\Omega))^2 \left( \sum_{n=1}^{m} \frac{\tau_n}{|\tau|} \|\mathbf{f}\|_{L^2(t_{n-1},t_n;Y)}^2 + 2(\|\mathbf{f}\|_{L^2(0,t_m;Y)}^2 + \nu|\mathbf{u}_0|^2_1) \right), \quad \text{and} \quad \gamma = \frac{1}{4} \left( \frac{\mu_1^{-1} + \mu_2^{-1}}{2} - 1 + \frac{5}{3}\delta_m \right). $$

Then, for each $n$, $1 \le n \le N$, $\tilde{\zeta}_n$ and $\zeta_n$, satisfy the following a posteriori estimates

$$ 
\begin{aligned}
\tilde{\zeta}_n \le \nu^{-\frac{1}{2}} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{L^2(t_{n-1}, t_n; X')} \\
+ \nu^{\frac{1}{2}} \|\mathbf{u} - \tilde{\mathbf{u}}_\tau\|_{L^2(t_{n-1}, t_n; X)} + \nu^{-\frac{1}{2}} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{L^2(t_{n-1}, t_n; X')} \\
+ (\nu \frac{\tau_n}{3})^{\frac{1}{2}} |\tilde{\mathbf{e}}^n - \tilde{\mathbf{e}}^{n-1}|_1.
\end{aligned} 
\quad (28) $$

$$ \frac{\sqrt{3}}{2} \zeta_n \le \frac{1}{\sqrt{2}} |\tau_n \varepsilon^n - \tau_{n-1} \varepsilon^{n-1}|_1 + \mathbf{E}_p(t_n, 0, 0). \quad (29) $$

**Proof.** By taking $\mathbf{w}$ and $\mathbf{v}$ equal to $\mathbf{u} - \tilde{\mathbf{u}}_\tau$ respectively in (22) and (24), and combining them, we get:

$$ 
\begin{aligned}
\frac{1}{2} \frac{d}{dt} \|\mathbf{u} - \tilde{\mathbf{u}}_\tau\|_0^2 + \nu |\mathbf{u} - \tilde{\mathbf{u}}_\tau|_1^2 - (\nabla \cdot (\mathbf{u} - \tilde{\mathbf{u}}_\tau), p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}) = \\
\qquad \langle \mathbf{f} - \pi_\tau \mathbf{f}, \mathbf{u} - \tilde{\mathbf{u}}_\tau \rangle + \nu (\nabla(\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_\tau), \nabla(\mathbf{u} - \tilde{\mathbf{u}}_\tau)).
\end{aligned} 
$$

Next, from (23) for $\mathbf{q}$ in $M$ and $\mathbf{v} = \nabla \mathbf{q}$, we have

$$ (\nabla \cdot (\mathbf{u} - \tilde{\mathbf{u}}_\tau), \mathbf{q}) = -(\mathbf{u} - \tilde{\mathbf{u}}_\tau, \nabla \mathbf{q}) = -(\mathbf{u} - \mathbf{u}_\tau, \nabla \mathbf{q}) - \tau_n (\nabla \Phi_\tau^*, \nabla \mathbf{q}) = \tau_n (\nabla \Phi_\tau^*, \nabla \mathbf{q}), $$

whence

$$ 
\begin{aligned}
\frac{1}{2} \frac{d}{dt} \|\mathbf{u} - \tilde{\mathbf{u}}_\tau\|_0^2 + \nu |\mathbf{u} - \tilde{\mathbf{u}}_\tau|_1^2 - \tau_n (\nabla \Phi_\tau^*, \nabla (p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1})) = \\
\qquad \langle \mathbf{f} - \pi_\tau \mathbf{f}, \mathbf{u} - \tilde{\mathbf{u}}_\tau \rangle + \nu (\nabla (\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_\tau), \nabla (\mathbf{u} - \tilde{\mathbf{u}}_\tau)).
\end{aligned} 
$$

Also, by inserting $\pi_\tau \Phi_\tau$, we obtain

$$ 
\begin{aligned}
\frac{1}{2} \frac{d}{dt} \|\mathbf{u} - \tilde{\mathbf{u}}_\tau\|_0^2 + \nu |\mathbf{u} - \tilde{\mathbf{u}}_\tau|_1^2 - \tau_n (\nabla \Phi_\tau^*, \nabla (p - \pi_\tau \Phi_\tau)) = \\
\qquad \langle \mathbf{f} - \pi_\tau \mathbf{f}, \mathbf{u} - \tilde{\mathbf{u}}_\tau \rangle + \nu (\nabla (\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_\tau), \nabla (\mathbf{u} - \tilde{\mathbf{u}}_\tau)) + (\nabla \Phi_\tau^*, \nabla (\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1})).
\end{aligned} 
$$
By applying (4) and observing moreover that

$$ (\nabla \Phi_{\tau}^{*}, \nabla(\tau_{n} \Phi^{n} - \tau_{n-1} \Phi^{n-1})) = \frac{t - t_{n}}{\tau_{n}^{2}} |\tau_{n} \Phi^{n} - \tau_{n-1} \Phi^{n-1}|_{1}^{2} \\ + \frac{\tau_{n}}{2} (|\Phi^{n}|_{1}^{2} - |\frac{\tau_{n-1}}{\tau_{n}} \Phi^{n-1}|_{1}^{2} + |\Phi^{n} - \frac{\tau_{n-1}}{\tau_{n}} \Phi^{n-1}|_{1}^{2}), $$

we successively deduce

$$ \frac{1}{2} \frac{d}{dt} \|\mathbf{u} - \tilde{\mathbf{u}}_{\tau}\|_{0}^{2} + \nu \|\mathbf{u} - \tilde{\mathbf{u}}_{\tau}\|_{1}^{2} \\ + \frac{\tau_n}{2} (|\frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}|_1^2 - |p|_1^2 + |p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}|_1^2) = \\ \langle \mathbf{f} - \pi_{\tau} \mathbf{f}, \mathbf{u} - \tilde{\mathbf{u}}_{\tau} \rangle + \nu (\nabla(\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_{\tau}), \nabla(\mathbf{u} - \tilde{\mathbf{u}}_{\tau})) \\ + \frac{t - t_{n-1}}{\tau_n} (\nabla(\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}), \nabla(p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1})), \quad (30) $$

and

$$ \frac{1}{2} \frac{d}{dt} \|\mathbf{u} - \tilde{\mathbf{u}}_{\tau}\|_{0}^{2} + \nu \|\mathbf{u} - \tilde{\mathbf{u}}_{\tau}\|_{1}^{2} + \frac{\tau_n}{2} (|\frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}|_{1}^{2} - |p|_{1}^{2} + |p - \pi_{\tau} \Phi_{\tau}|_{1}^{2}) = \\ \langle \mathbf{f} - \pi_{\tau} \mathbf{f}, \mathbf{u} - \tilde{\mathbf{u}}_{\tau} \rangle + \nu (\nabla(\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_{\tau}), \nabla(\mathbf{u} - \tilde{\mathbf{u}}_{\tau})) \\ + \frac{t - t_n}{\tau_n} (\nabla(\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}), \nabla(p - \pi_{\tau} \Phi_{\tau})) \\ + (\frac{t - t_n}{\tau_n^2} + \frac{1}{2\tau_n}) |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2. \quad (31) $$

Besides, we note also, for all $n, 1 \le n \le N$ and for all $t \in [t_{n-1}, t_n]$, that

$$ \nu \int_{t_{n-1}}^{t_n} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_{\tau}|_1^2 dt = \nu \frac{\tau_n}{3} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1^2. \quad (32) $$

By integrating (30) between $t_{n-1}$ and $t_n$, then from (32), the Cauchy-Schwarz inequality and (5), for a given real $\mu_1$ in $[0, 1]$, we derive

$$ \|\mathbf{u}(t_n) - \tilde{\mathbf{u}}^n\|_0^2 - \|\mathbf{u}(t_{n-1}) - \tilde{\mathbf{u}}^{n-1}\|_0^2 + \nu \int_{t_{n-1}}^{t_n} |\mathbf{u} - \tilde{\mathbf{u}}_{\tau}|_1^2 dt \\ + |\tau_{n-1} \Phi^{n-1}|_1^2 + (1 - \mu_1) \tau_n \int_{t_{n-1}}^{t_n} |p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}|_1^2 dt \\ \le \frac{2}{\nu} \int_{t_{n-1}}^{t_n} \|\mathbf{f} - \pi_{\tau} \mathbf{f}\|_{X'}^2 dt + 2\nu \frac{\tau_n}{3} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1^2 \\ + \frac{\mu_1}{3} |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2 + \tau_n \int_{t_{n-1}}^{t_n} |p|_1^2 dt. $$

Similarly, considering also (31), we obtain identical estimate, for a given real $\mu_2$ in $[0, 1]$, while replacing the error $p - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}$ by $p - \pi_{\tau} \Phi_{\tau}$. Consequently, by adding up these estimates, we get

$$ \|\mathbf{u}(t_n) - \tilde{\mathbf{u}}^n\|_0^2 - \|\mathbf{u}(t_{n-1}) - \tilde{\mathbf{u}}^{n-1}\|_0^2 + |\tau_{n-1} \Phi^{n-1}|_1^2 \\ + \nu \int_{t_{n-1}}^{t_n} |\mathbf{u} - \tilde{\mathbf{u}}_{\tau}|_1^2 dt + \frac{1}{2} \mathbf{E}_p^2(t_n, \mu_1, \mu_2) \\ \le 2\nu \frac{\tau_n}{3} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1^2 + \frac{1}{2} \frac{\mu_1^{-1} + \mu_2^{-1}}{3} |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2 \\ + \frac{2}{\nu} \int_{t_{n-1}}^{t_n} \|\mathbf{f} - \pi_{\tau} \mathbf{f}\|_{X'}^2 dt + \tau_n \int_{t_{n-1}}^{t_n} |p|_1^2 dt. \quad (33) $$
From other part, it can be observed from (12) that for all $n \in [1, N]$, the following holds :

$$ \|\tau_n \nabla \Phi^n\|_0^2 = \tau_n (\tilde{u}^n - u^n, \nabla \Phi^n) = \|\tilde{u}^n - u^n\|_0^2 $$

and also that

$$ \|\tau_n \nabla \Phi^n\|_0^2 = \|\mathbf{u}(t_n) - \tilde{\mathbf{u}}^n\|_0^2 - \|\mathbf{u}(t_n) - \mathbf{u}^n\|_0^2. \quad (34) $$

In this case, (33) becomes

$$ \begin{aligned} & \|\mathbf{u}(t_n) - \mathbf{u}^n\|_0^2 - \|\mathbf{u}(t_{n-1}) - \mathbf{u}^{n-1}\|_0^2 + |\tau_n \Phi^n|_1^2 \\ & + \nu \int_{t_{n-1}}^{t_n} |\mathbf{u} - \tilde{\mathbf{u}}_\tau|_1^2 dt + \frac{1}{2} \mathbf{E}_p^2(t_n, \mu_1, \mu_2) \\ & \le 2\nu \frac{\tau_n}{3} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1^2 + \frac{1}{2} \frac{\mu_1^{-1} + \mu_2^{-1}}{3} |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2 \\ & \quad + \frac{2}{\nu} \int_{t_{n-1}}^{t_n} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{X'}^2 dt + \tau_n \int_{t_{n-1}}^{t_n} |p|_1^2 dt. \end{aligned} \quad (35) $$

In particular, we observe that

$$ 3|\tau_n \Phi^n|_1^2 - |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2 = \{|\tau_n \Phi^n|_1^2 - |\tau_{n-1} \Phi^{n-1}|_1^2\} \\ + \{|\tau_n \Phi^n + \tau_{n-1} \Phi^{n-1}|_1^2 - |\tau_{n-1} \Phi^{n-1}|_1^2\}. $$

Thus, setting $R_m = \sum_{n=1}^m \{|\tau_n \Phi^n + \tau_{n-1} \Phi^{n-1}|_1^2 - |\tau_{n-1} \Phi^{n-1}|_1^2\}$, we obtain, after summing with respect to $n$, the inequality

$$ \begin{aligned} & \|\mathbf{u}(t_m) - \mathbf{u}^m\|_0^2 + \nu \int_0^{t_m} |\mathbf{u} - \tilde{\mathbf{u}}_\tau|_1^2 dt + \frac{1}{3}|\tau_m \Phi^m|_1^2 + \frac{1}{3}R_m + \frac{1}{2}\sum_{n=1}^m \mathbf{E}_p^2(t_n, \mu_1, \mu_2) \\ & \le 2\nu \sum_{n=1}^m \frac{\tau_n}{3} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1^2 + \frac{1}{3}(\frac{\mu_1^{-1} + \mu_2^{-1}}{2} - 1) \sum_{n=1}^m |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2 \quad \text{which by} \\ & \quad + \frac{2}{\nu} \int_0^{t_m} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{X'}^2 dt + \sum_{n=1}^m \tau_n \int_{t_{n-1}}^{t_m} |p|_1^2 dt. \end{aligned} $$

taking into account (34), yields

$$ \begin{aligned} & \frac{2}{3} \|\mathbf{u}(t_m) - \mathbf{u}^m\|_0^2 + \frac{1}{3} \|\mathbf{u}(t_m) - \tilde{\mathbf{u}}^m\|_0^2 + \nu \int_0^{t_m} |\mathbf{u} - \tilde{\mathbf{u}}_\tau|_1^2 dt \\ & + \frac{1}{3} R_m + \frac{1}{2} \sum_{n=1}^m \mathbf{E}_p^2(t_n, \mu_1, \mu_2) \\ & \le 2\nu \sum_{n=1}^m \frac{\tau_n}{3} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1^2 + \frac{1}{3} (\frac{\mu_1^{-1} + \mu_2^{-1}}{2} - 1) \sum_{n=1}^m |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2 \quad \text{In ad-} \\ & \quad + \frac{2}{\nu} \int_0^{t_m} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{X'}^2 dt + \sum_{n=1}^m \tau_n \int_{t_{n-1}}^{t_m} |p|_1^2 dt. \end{aligned} $$

dition, we first notice that

$$ |\tau_n \Phi^n + \tau_{n-1} \Phi^{n-1}|_1^2 - |\tau_{n-1} \Phi^{n-1}|_1^2 = Q_n - \frac{1}{4} |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2, $$

where we denote by $Q_n = 2|\tau_n \Phi^n|_1^2 + |\tau_{n-1} \Phi^{n-1}|_1^2 - \frac{3}{4}|\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2$. In particular, we also remark that $Q_n = \frac{1}{2}(|\tau_n \Phi^n|_1^2 - |\tau_{n-1} \Phi^{n-1}|_1^2) + \frac{3}{4}|\tau_n \Phi^n + \tau_{n-1} \Phi^{n-1}|_1^2$. Consequently, if $R_m < 0$ then there exists $\delta_m$ in $]0, \frac{1}{4}]$ such that

$$ \sum_{n=1}^{m} \{2|\tau_n \Phi^n|_1^2 + |\tau_{n-1} \Phi^{n-1}|_1^2 - (1 - \delta_m)|\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2\} = 0. $$
Hence, $R_m = -\delta_m \sum_{n=1}^{m} |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1^2$, by setting $\beta_m = (\frac{\mu_1^{-1}+\mu_2^{-1}}{2} - 1 + \delta_m)$

and using the triangular inequality while inserting (25) then noting (26), we derive

$$ 
\begin{aligned}
& \frac{2}{3} \left( \|\mathbf{u}(t_m) - \mathbf{u}^m\|_0^2 + \nu \int_0^{t_m} |\mathbf{u} - \tilde{\mathbf{u}}_\tau|^2 dt \right) \\
& \qquad + \frac{1}{3} [\mathbf{u} - \tilde{\mathbf{u}}_\tau]^2(t_m) + \frac{1}{2} \sum_{n=1}^m \mathbf{E}_p^2(t_n, \mu_1, \mu_2) \\
& \leq +\frac{2}{3} \sum_{n=1}^m \{2\nu\tau_n |\tilde{\mathbf{e}}^n - \tilde{\mathbf{e}}^{n-1}|_1^2 + \beta_m |\tau_n\varepsilon^n - \tau_{n-1}\varepsilon^{n-1}|_1^2\} \\
& \qquad + \frac{2}{\nu} \int_0^{t_m} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{X'}^2 dt + \sum_{n=1}^m \tau_n \int_{t_{n-1}}^{t_m} |p|_1^2 dt.
\end{aligned}
\quad (36) $$

Furthermore, when $R_m \ge 0$, this inequality is also satisfied with $\delta_m = 0$. Next, in order to estimate the last term in (36), we deduce from the first equation of (1), and (17), that for almost any $t \in [0, T]$, we have $|p|_1 \le \|\mathbf{f}\|_0 + \|\partial_t \mathbf{u}\|_0 + \nu \|\mathbf{u}\|_2$.

Then, by applying the operator $P_H$ to the same equation, we obtain $\partial_t \mathbf{u} + \nu \mathbf{A} \mathbf{u} = P_H \mathbf{f}$.

Consequently, from the condition (R) and the continuity of the operator $P_H$, we then have $\nu \|\mathbf{u}\|_2 \le C(\Omega)\{\|\partial_t \mathbf{u}\|_0 + \|\mathbf{f}\|_0\}$. Whence

$$ \sum_{n=1}^{m} \tau_n \int_{t_{n-1}}^{t_n} |p|_1^2 dt \le 2(1 + C(\Omega))^2 \sum_{n=1}^{m} \tau_n (\|\mathbf{f}\|_{L^2(t_{n-1}, t_n; Y)}^2 + \|\partial_t \mathbf{u}\|_{L^2(t_{n-1}, t_n; Y)}^2). $$

Now, using the classical Faedo Galerkin method (we refer to [6, Chap.XIX, Prop.2] for a similar proof), we derive at each time $t_m$, for any $m, 1 \le m \le N$, the a priori estimate $\|\partial_t \mathbf{u}\|_{L^2(0, t_m; Y)} \le \sqrt{\nu} |\mathbf{u}_0|_1 + \|\mathbf{f}\|_{L^2(0, t_m; Y)}$ and therefore

$$ \sum_{n=1}^{m} \tau_n \int_{t_{n-1}}^{t_n} |p|_1^2 dt \le 2(1 + C(\Omega))^2 \left\{ \sum_{n=1}^{m} \tau_n \|\mathbf{f}\|_{L^2(t_{n-1}, t_n; Y)}^2 + 2|\tau|(\|\mathbf{f}\|_{L^2(0, t_m; Y)}^2 + \nu|\mathbf{u}_0|_1^2) \right\}. $$

Now, we prove a similar upper bound for the error function $\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)$ in the norm of $X'$. In fact, if we observe that

$$ \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{X'} = \sup_{\mathbf{w} \in X} \frac{(\partial_t(\mathbf{u} - \mathbf{u}_\tau), \mathbf{w}) - (p - \pi_\tau \Phi_\tau, \nabla \cdot \mathbf{w})}{|\mathbf{w}|_1}, $$

and using equation (22), we get for any $t \in [t_{n-1}, t_n]$,

$$ \nu^{-\frac{1}{2}} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{X'} \le \nu^{-\frac{1}{2}} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{X'} + \nu^{\frac{1}{2}} |\mathbf{u} - \tilde{\mathbf{u}}_\tau|_1 + \nu^{\frac{1}{2}} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_\tau|_1. $$

Taking the square of this inequality, then integrating between $t_{n-1}$ and $t_n$, from (32) and using the same arguments to establish (36), we derive

$$ (2\nu)^{-\frac{1}{2}} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{L^2(0, t_m; X')} \le \left( \nu \|\mathbf{u} - \tilde{\mathbf{u}}_\tau\|_{L^2(0, t_m; X)}^2 + \frac{1}{\nu} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{L^2(0, t_m; X')}^2 \right. \\
\left. + 2 \sum_{n=1}^m \tilde{\zeta}_n^2 + \frac{2}{3} \sum_{n=1}^m \nu \tau_n |\tilde{\mathbf{e}}^n - \tilde{\mathbf{e}}^{n-1}|_1^2 \right)^{\frac{1}{2}}. $$
For $\mu_1 = \mu_2 = 1$, the second term in the right-hand side of the previous inequality is bounded in (36), and this leads to

$$ (2\nu)^{-\frac{1}{2}} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{L^2(0,t_m;X')} \\ \le (6\tilde{\zeta}_{m\tau}^2 + 2\delta_m \zeta_{m\tau}^2)^{\frac{1}{2}} + \sqrt{\frac{3}{\nu}} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{L^2(0,t_m;X')} + |2\tau|^{\frac{1}{2}} C_m(\mathbf{f}, \mathbf{u}_0) \\ + \sqrt{\frac{2}{3}} \left( \sum_{n=1}^{m} \{3\nu\tau_n |\tilde{\mathbf{e}}^n - \tilde{\mathbf{e}}^{n-1}|_1^2 + \delta_m |\tau_n\varepsilon^n - \tau_{n-1}\varepsilon^{n-1}|_1^2 \} \right)^{\frac{1}{2}}. \quad (37) $$

Finally, we conclude the proof of (27) by combining (36) with (37). Conversely, thanks to triangular inequality, we have

$$ \tilde{\zeta}_n \le \left( \nu \frac{\tau_n}{3} \right)^{\frac{1}{2}} \left( |\tilde{\mathbf{e}}^n - \tilde{\mathbf{e}}^{n-1}|_1 + |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1 \right). \quad (38) $$

To bound the last term, we take first w equal to $\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}_\tau$ in (22), next we integrate between $t_{n-1}$ and $t_n$, then from (32) we deduce

$$ \left(\nu \frac{\tau_n}{3}\right)^{\frac{1}{2}} |\tilde{\mathbf{u}}^n - \tilde{\mathbf{u}}^{n-1}|_1 \le \nu^{-\frac{1}{2}} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{L^2(t_{n-1}, t_n; X')} \\ + \nu^{\frac{1}{2}} \|\mathbf{u} - \tilde{\mathbf{u}}_\tau\|_{L^2(t_{n-1}, t_n; X)} + \nu^{-\frac{1}{2}} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{L^2(t_{n-1}, t_n; X)}. $$

So, inserting this estimate in (38), we get (28). Similarly, we note that

$$ \zeta_n \le \frac{1}{\sqrt{3}} (|\tau_n \varepsilon^n - \tau_{n-1} \varepsilon^{n-1}|_1 + |\tau_n \Phi^n - \tau_{n-1} \Phi^{n-1}|_1). $$

Taking the square of this inequality, we deduce

$$ \zeta_n^2 \le \frac{2}{3} \left( |\tau_n \varepsilon^n - \tau_{n-1} \varepsilon^{n-1}|_1^2 + \tau_n \int_{t_{n-1}}^{t_n} |\Phi^n - \frac{\tau_{n-1}}{\tau_n} \Phi^{n-1}|_1^2 dt \right). \quad (39) $$

This concludes, the proof of the Theorem. $\square$

<figure><img src="image_1.png" /></figure>

The local lower bound of the error (28), is similar to [2, Proposition 3.3]. In particular, (29) concerns the incompressibility part of the Stokes equations while (28) rely on its evolution or diffusive part and is derived independently thereof.

Next, taking the square of (28), multiplying then (39) by $\gamma$ and summing on the $n$, and noting that $\nu \int_0^{t_m} |\mathbf{u} - \tilde{\mathbf{u}}_\tau|^2 dt \le \frac{2}{3} \mathbf{E}_{u\tau}^2(t_m)$, we get

$$ \frac{1}{4}(\tilde{\zeta}_{m\tau}^2 + \gamma\zeta_{m\tau}^2) \le 2 \left( \mathbf{E}_{u\tau}^2(t_m) + (2\nu)^{-1} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{L^2(0,t_m;X')}^2 \right. \\ \left. + \frac{\gamma}{6} \sum_{n=1}^m \mathbf{E}_p^2(t_n, 0, 0) \right) + \nu^{-1} \|\mathbf{f} - \pi_\tau \mathbf{f}\|_{L^2(0,t_m;X')}^2 \\ + \frac{1}{3} \sum_{n=1}^m \left\{ \nu\tau_n |\tilde{\mathbf{e}}^n - \mathbf{e}^{n-1}|_1^2 + \frac{\gamma}{2} |\tau_n\varepsilon^n - \tau_{n-1}\varepsilon^{n-1}|_1^2 \right\}. $$

Consequently, provided that the regularity parameter $\sigma$ is bounded independently of $\tau$, then for $\mu_1$ and $\mu_2$ satisfying $\frac{3}{4}(1 - \mu_1) = \frac{\gamma}{6} = \frac{3}{4}(1 - \mu_2)$, the full error

$$ \left\{ \mathbf{E}_{u\tau}^2(t_m) + \frac{3}{4} \sum_{n=1}^{m} \mathbf{E}_p^2(t_n, \mu_1, \mu_2) \right\}^{\frac{1}{2}} + (2\nu)^{-\frac{1}{2}} \|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau \Phi_\tau)\|_{L^2(0,t_m;X')} \text{ is} $$
equivalent to the Hilbertian sum $(\zeta_{m\tau}^2 + \gamma\zeta_{m\tau}^2)^{\frac{1}{2}}$ up to some terms involving the data $(\mathbf{u}_0, \mathbf{f})$ and the spatial error $\left(\sum_{n=1}^{m} \{\nu\tau_n|\tilde{\mathbf{e}}^n - \tilde{\mathbf{e}}^{n-1}|_1^2 + \frac{\gamma}{2}|\tau_n\varepsilon^n - \tau_{n-1}\varepsilon^{n-1}|_1^2\}\right)^{\frac{1}{2}}$.

In this case, we observe that $\mu_1(=\mu_2)$ is solution of $18\mu^2 - (19 - \frac{5}{3}\delta_m)\mu + 1 = 0$ which admits two solutions in $[0, 1]$. So, selecting $\mu_1 = \frac{1}{36}(19 - \frac{5}{3}\delta_m + \sqrt{(19 - \frac{5}{3}\delta_m)^2 - 72})$, we notice that for $\delta_m \in [0, \frac{1}{4}]$, then $\mu_1(\delta_m) \in [0.97545, 1]$. In particular, we observe that for $\mu_1(0) = 1$ the full error $\mathbf{E}_{u\tau}(t_m) + (2\nu)^{-\frac{1}{2}}\|\partial_t(\mathbf{u} - \mathbf{u}_\tau) + \nabla(p - \pi_\tau\Phi_\tau)\|_{L^2(0,t_m;X')}$ is equivalent to the quantity $\tilde{\zeta}_{m\tau}$ up to some terms involving the data and the spatial error. Which situation is similar to the case of the Euler scheme (see [2, §3]).

### 4. Adaptive algorithm : Time step size control

In the following, we introduce a simple procedure allowing to control the time step size, based on the local lower bound (28)$^2$. For this, for each $n$, $1 \le n \le N$, we introduce the local norm

$$A_n = \left( \nu \int_{t_{n-1}}^{t_n} |\tilde{\mathbf{u}}_{h\tau}|_1^2 dt \right)^{\frac{1}{2}} + \nu^{-\frac{1}{2}} \left( \int_{t_{n-1}}^{t_n} \|\partial_t \mathbf{u}_{h\tau} + \tau_n \nabla \Phi_h^n\|_{X'}^2 dt \right)^{\frac{1}{2}}.$$

The negative norm in the last term of $A_n$ is approximated by $\|\mathbf{v}\|_1$ where $\mathbf{v}$ is the solution of the Laplace equation $-\Delta\mathbf{v} = \partial_t \mathbf{u}_{h\tau} + \tau_n \nabla\Phi_h^n$ with homogeneous Dirichlet boundary conditions at each iteration$^3$. Next, we denote by **tol** the prescribed tolerance in order to bound the error. Then, for a given parameter $\theta$ in $[0, 1]$, and for each $n$, $1 \le n \le N$, we assume the condition $\theta tol \le \frac{\tilde{\zeta}_n}{A_n} \le tol$ is satisfied. So, for a fixed $h$, an initial guess $\tau_1$ and for a fixed value of $\sigma$ (6), we adopt the following procedure :

while $t_n \le T$

$\tau_n = \tau_{n-1}; t_n = t_{n-1} + \tau_n$; solve (7)-(8); compute $\tilde{\zeta}_n$ and $A_n$;

if $tol < \frac{\tilde{\zeta}_n}{A_n}$, $\tau_n = tol \frac{A_n}{\zeta_n} \tau_n$ end if;

if $\frac{\tilde{\zeta}_n}{A_n} < \theta tol$, $\tau_n = \min(\theta tol \frac{A_n}{\zeta_n}, \sigma)\tau_n$ end if;

end while

**A simple test case : the Poiseuille flow.** We assume that, $\Omega = [0, L[\times] - H, +H[$ is a rectangular tube, triangulated by a uniform mesh, made of isosceles rectangular triangles with diameter $h$. We consider the Taylor-Hood element, which uses the polynomial functions of a degree 2 and 1 to approximate respectively the velocity and the pressure. The flow is generated by prescribing a parabolic and horizontal profile of the velocity on $\{0\}\times[-H, +H[$. We assume that, the exact solution is given by :

$$(\mathbf{u}; p) = \alpha(t)((H - y)(H + y), 0; -2\nu(x - L)),$$

with $\alpha(t) = (1+\delta_\epsilon(t)+\cos(r\pi t))^{-1} - (2+\delta_\epsilon(t))^{-1}$; for a fixed $\epsilon$ in $[0, 1]$, $r$ in $\mathbb{N}^*$, and for $T = \frac{6}{r}$, $\delta_\epsilon = \epsilon \times \{\chi_{[0, \frac{2}{r}]} + 0.6 \times \chi_{[\frac{2}{r}, \frac{4}{r}]} + 0.8 \times \chi_{[\frac{4}{r}, T]}\}$. Here, $\chi_{[a,b]}$ defines the characteristic function

2. we also refer to [12] for similar idea and different context

3. we refer to [4, 17] for similar situations
associated with any time interval $[a, b]$. In this case, the exact solution reaches its maximal values at the times $\frac{k}{r}$ in $[0, T]$ with $k$ odd integer. For the numerical test, we take : $L = 0.15$, $H = 0.015$, $\nu = 10^{-3}$, $\epsilon=0.25$, $r=16$, $h = \sqrt{5} \frac{H}{6}$, **tol**=7.5x10⁻⁵, $\tau_1$=1.2x10⁻³, $\sigma$=1.5 and $\theta$=0.5.

Now, for each $n$, $1 \le n \le N$, we reported the time steps $\tau_n$ (Fig.1-(a)) and the global error $[\mathbf{u} - \tilde{\mathbf{u}}_{h\tau}](t_n)$ (Fig.1-(b)), obtained with an adapted and a constant time step. In particular, the constant step corresponds to the average of the variable steps and equal to $9.4 \times 10^{-4}$. Besides, we noticed identical CPU time for both cases. Clearly, we first observe that, the sequence of time

<figure><img src="image_2.png" /><figcaption>**Figure 1.** *Behaviour of (a) $\tau_n$ and (b) $[\mathbf{u} - \tilde{\mathbf{u}}_{h\tau}](t_n)$.*</figcaption></figure>

steps $\tau_n$ reaches its minimum at the peaks level, where the error $[\mathbf{u} - \tilde{\mathbf{u}}_{h\tau}](t_n)$ is maximal. In addition, when time is advancing, the error obtained with the adaptive steps becomes increasingly reduced. The computations are carried out on PC Toshiba (A100-088) intel core-duo (1 GB of RAM) and by using the finite element code Freefem++, see [10].

**5. Conclusion and perspectives**

In this work, we presented a *posteriori* analysis of the original Chorin-Temam scheme by a residual approach. We assumed, the time and the spatial discretization errors are independent. In order to control the time step size, we set up a simple algorithm proving the efficiency of the time error estimators. In particular, only the estimator associated with the velocity seeked at the prediction step is used for the time adaptation. Moreover, we assume that the spatial discretiza- tion error and some other terms involving the data are not considered in the present algorithm. Besides, the present analysis can be extended to other higher order projection schemes<fnref n="4" />. It can be easily developed in the framework of any spatial discretization method viz., the finite volumes or spectral elements. In addition, using a residual approach with a conforming finite-element method, the first analysis given in [11], introduces another family of estimators associated with the spatial discretization error. In practice, they are well adapted for mesh adaptivity. In the same way, they can be also combined with time adaptivity in one single procedure. This will be addressed in a forthcoming paper.

**6. References**

[1] A. BERGAM, C. BERNARDI, Z. MGHAZLI, "A posteriori analysis of the finite element discretization of some parabolic equations", *Math. Comput.*, vol. 74, num. 251, 2004.

<footnote n="4">4. see [8] for a recent overview</footnote>
[2] C. BERNARDI, R. VERFÜRTH, "A posteriori error analysis of the fully discretized time-dependent Stokes equations", M2AN, vol. 38 num. 3, 2004.

[3] C. BERNARDI, Y. MADAY, F. RAPETTI, "Discrétisations variationnelles de problèmes aux limites elliptiques", *Mathématiques et Applications*, SMAI, Springer-Verlag. 2004.

[4] S. BERRONE, "Robust *a posteriori* error estimates for finite element discretizations of the Heat equations with discontinuous coefficients", M2AN, vol. 40 num. 6, 2007.

[5] A. J. CHORIN, "On the convergence of discrete approximations to Navier-Stokes equations", *Mathematics of Computation*, vol. 23, 1969.

[6] R. DAUTRAY, J.-L. LIONS, "Analyse mathématique et calcul numérique pour les sciences et les techniques", *Collection CEA*, Tome 3, Masson, 1985.

[7] V. GIRAULT, P.A. RAVIART, "Finite Element Methods for Navier-Stokes equations, theory and algorithms", Springer-Verlag Berlin Heidelberg, 1986.

[8] J.-L. GUERMOND, P. MINEV, J. SHEN, "An overview of projection methods for incompressible flows", *Comp. Methods Appl. Mech. Engrg.* vol. 195, 2006.

[9] J.-L. GUERMOND, "Some implementation of projection methods for Navier-Stokes equations", M2AN vol. 30, num. 5, 1996.

[10] F. HECHT, O. PIRONNEAU, A. LE HYARIC, K. OHTSUKA, "FreeFem++", www.freefem.org.

[11] N. KHARRAT, Z. MGHAZLI, "Residual error estimators for the time-dependent Stokes equations", C. R. Acad. Sci. Paris, vol. Sér. I, num. 340, 2005.

[12] M. PICASSO, "Adaptive finite elements for a linear parabolic problem", *Comput. Methods Appl. Mech. Engrg.*, vol. 167, 1998.

[13] A. PROHL, "Projection and Quasi-Compressibility Methods for Solving the incompressible Navier-Stokes Equations", Teubner-Verlag, Stuttgart, 1997.

[14] R. RANNACHER, "On Chorin's projection method for the incompressible Navier-Stokes equations", *Lecture Notes in Mathematics. eds : J. G. Heywood, K. Masuda, R. Rautmann, S.A. Solonnikov in : The Navier-Stokes equations II-Theory and Numerical methods. Proc. Oberwolfach*, vol. 1530, 1991.

[15] R. TEMAM, "Sur l'approximation de la solution des équations de Navier-Stokes par la méthode des pas fractionnaires ii", *Arch. Ration. Mech. Anal.*, vol. 33, 1969.

[16] R. TEMAM, "Theory and Numerical analysis of the Navier-Stokes equations", North-Holland, 1977.

[17] R. VERFÜRTH, "Robust *a posteriori* error estimates for non-stationary convection-diffusion equations", *SIAM J. on Numerical Analysis*, vol. 43, num. 4, 2005.