<html><head><meta name="color-scheme" content="light dark"></head><body><pre style="word-wrap: break-word; white-space: pre-wrap;">%~MoulinÃ© par MaN_auto v.0.33.4 2025-02-11 15:14:10
\documentclass[CRMATH,Unicode,XML,biblatex]{cedram}

\TopicFR{Analyse numÃ©rique}
\TopicEN{Numerical analysis}

\addbibresource{CRMATH_Jamelot_20230366.bib}

\usepackage{amssymb}
\usepackage{mathtools}
\usepackage{tikz}
\usetikzlibrary{calc,intersections}
\usepackage{url}

\DeclareMathOperator{\curls}{curl}
\DeclareMathOperator{\grad}{\mathbf{grad}}
\DeclareMathOperator{\curlv}{\mathbf{curl}}
\DeclareMathOperator{\Grad}{\mathbf{Grad}}
\DeclareMathOperator{\Deltavec}{\mathbf{\Delta}}
\DeclareMathOperator{\divergence}{div}
\let\div=\divergence
\DeclareMathOperator{\card}{card}
\DeclareMathOperator{\vect}{vect}

\usetikzlibrary{calc,intersections}

\newcommand{\eps}{\varepsilon}
\newcommand{\Fcal}{\mathcal{F}}
\newcommand{\Ical}{\mathcal{I}}
\newcommand{\Jcal}{\mathcal{J}}
\newcommand{\Lcal}{\mathcal{L}}
\newcommand{\Pcal}{\mathcal{P}}
\newcommand{\Scal}{\mathcal{S}}
\newcommand{\Tcal}{\mathcal{T}}
\newcommand{\Xcal}{\mathcal{X}}
\newcommand{\bbA}{\mathbb{A}}
\newcommand{\bbB}{\mathbb{B}}
\newcommand{\R}{\mathbb{R}}
\newcommand{\bbL}{\mathbb{L}}
\newcommand{\N}{\mathbb{N}}
\newcommand{\bL}{\mathbf{L}}
\newcommand{\bH}{\mathbf{H}}
\newcommand{\bP}{\mathbf{P}}
\newcommand{\bV}{\mathbf{V}}
\newcommand{\bX}{\mathbf{X}}
\newcommand{\avec}{\mathbf{a}}
\newcommand{\bvec}{\mathbf{b}}
\newcommand{\dvec}{\mathbf{d}}
\newcommand{\evec}{\mathbf{e}}
\newcommand{\fvec}{\mathbf{f}}
\newcommand{\nvec}{\mathbf{n}}
\newcommand{\svec}{\mathbf{s}}
\newcommand{\uvec}{\mathbf{u}}
\newcommand{\vvec}{\mathbf{v}}
\newcommand{\wvec}{\mathbf{w}}
\newcommand{\xvec}{\mathbf{x}}
\newcommand{\hatF}{\hat{F}}
\newcommand{\hatK}{\hat{K}}
\newcommand{\Om}{\Omega}
\newcommand{\pa}{\partial}
\newcommand{\vertiii}[1]{{\left\vert\kern-0.25ex\left\vert\kern-0.25ex\left\vert #1
\right\vert\kern-0.25ex\right\vert\kern-0.25ex\right\vert}}
\newcommand{\ul}{\underline}

\newcommand*{\zerovec}{\mathbf{0}}
\newcommand*{\conf}{\mathrm{c}}
\newcommand*{\CR}{\mathrm{CR}}
\newcommand*{\disc}{\mathrm{disc}}
\newcommand*{\FS}{\mathrm{FS}}
\newcommand*{\LG}{\mathrm{LG}}
\newcommand*{\nonconf}{\mathrm{nc}}
\newcommand*{\PS}{\mathrm{PS}}
\newcommand*{\RT}{\mathrm{RT}}
\newcommand*{\zmv}{\mathrm{zmv}}

\renewcommand{\lesssim}{\mathrel{\raisebox{-2.5pt}{$\stackrel{\raisebox{-1pt}{\scalebox{1}[.8]{&lt;}}}{\sim}$}}}

%%%---------------------------------------------------------------------------

%%% Changer la largeur des adresses
% Jouer avec \hsize :
\makeatletter
\def\@setafterauthor{%
  \vglue3mm%
%  \hspace*{0pt}%
\begingroup\hsize=11cm\advance\hsize\abstractmarginL\raggedright
\noindent
%\hspace*{\abstractmarginL}\begin{minipage}[t]{10cm}
   \leftskip\abstractmarginL
  \normalfont\Small
  \@afterauthor\par
\endgroup
\vskip2pt plus 3pt minus 1pt
}
\makeatother

%%%---------------------------------------------------------------------------

%%% \to macro
\renewcommand*{\to}{\mathchoice{\longrightarrow}{\rightarrow}{\rightarrow}{\rightarrow}}

%%% \mapsto macro
\let\oldmapsto\mapsto
\renewcommand*{\mapsto}{\mathchoice{\longmapsto}{\oldmapsto}{\oldmapsto}{\oldmapsto}}

%%% \longhookrightarrow macro
\newcommand{\longhookrightarrow}{\!\lhook\joinrel\longrightarrow}

%%%---------------------------------------------------------------------------

%\theoremstyle{remark}
%\newtheorem{step}{Step}
\newcounter{stepA}
\newenvironment{stepA}[1][]{\refstepcounter{stepA}\begin{proof}[{Step~\thestepA\if#1\empty\else: #1\fi}]}{\let\qed\relax\end{proof}}

\newcounter{stepB}
\newenvironment{stepB}[1][]{\refstepcounter{stepB}\begin{proof}[{Step~\thestepB\if#1\empty\else: #1\fi}]}{\let\qed\relax\end{proof}}

%%%---------------------------------------------------------------------------

%%% ''Evaluated at'' macro
\newcommand\restr[2]{{% we make the whole thing an ordinary symbol
 \left.\kern-\nulldelimiterspace % automatically resize the bar with \right
 #1 % the function
 \vphantom{\big|} % pretend it's a little taller at normal size
 \right|_{#2} % this is the delimiter
 }}

%%% \abs-like macros
%\usepackage{mathtools}
\DeclarePairedDelimiter{\abs}{\lvert}{\rvert}
\DeclarePairedDelimiter{\norm}{\lVert}{\rVert}
\newcommand{\normsub}{\norm*}					% Use with subscript.
% Star the operator for a scaling version.
% /!\ Attention: use case-sensitive Find&amp;Replace.

%%% ''Such that'' macro
\newcommand*{\tq}{\,\vert\,}
\newcommand*{\btq}{\,\big\vert\,}
% Scaling version (with \left and \right delimiters): \,\middle\vert\,

%%% Divides
% \vert

%%% Separator (non-scaling)
% \mid

%%%---------------------------------------------------------------------------

\graphicspath{{./figures/}}

\newcommand*{\mk}{\mkern -1mu}
\newcommand*{\Mk}{\mkern -2mu}
\newcommand*{\mK}{\mkern 1mu}
\newcommand*{\MK}{\mkern 2mu}

\hypersetup{urlcolor=purple, linkcolor=blue, citecolor=red}

\let\oldtilde\tilde
\renewcommand*{\tilde}[1]{\mathchoice{\widetilde{#1}}{\widetilde{#1}}{\oldtilde{#1}}{\oldtilde{#1}}}
\let\oldhat\hat
\renewcommand*{\hat}[1]{\mathchoice{\widehat{#1}}{\widehat{#1}}{\oldhat{#1}}{\oldhat{#1}}}
\let\oldexists\exists
\renewcommand*{\exists}{\mathrel{\oldexists}}
\let\oldforall\forall
\renewcommand*{\forall}{\mathrel{\oldforall}}

\title{Stability estimates for solving Stokes problem with nonconforming finite elements}
\alttitle{Estimations de stabilitÃ© pour rÃ©soudre le problÃ¨me de Stokes avec des Ã©lÃ©ments finis non conformes}

\author{\firstname{Erell} \lastname{Jamelot}}
\address{UniversitÃ© Paris-Saclay, CEA, Service de Thermo-hydraulique et de MÃ©canique des Fluides, 91191 Gif-sur-Yvette, France}
\thanks{CEA SIMU/SITHY project.}
\email{erell.jamelot@cea.fr}

\subjclass{65N30, 35J57, 76D07}

\begin{abstract}
We propose to analyse the discretization of the Stokes problem with nonconforming finite elements in light of the T-coercivity. First we exhibit a family of operators to prove T-coercivity and we show that the stability constant is equal to the classical one up to a constant which depends on the Babu\v{s}ka--Aziz constant. Then we explicit the stability constants with respect to the shape regularity parameter for order 1 in 2 or 3 dimensions, and order 2 in 2 dimensions. In this last case, we improve the result of the original Fortin--Soulie paper. Second, we illustrate the importance of using a divergence-free velocity reconstruction on some numerical experiments.
\end{abstract}

\begin{altabstract}
Nous proposons d'analyser la discrÃ©tisation du problÃ¨me de Stokes avec des Ã©lÃ©ments finis non conformes Ã&nbsp; la lumiÃ¨re de la T-coercivitÃ©. Tout d'abord, pour prouver la T-coercivitÃ©, nous exhibons une famille d'opÃ©rateurs et nous montrons que la constante de stabilitÃ© est Ã©gale Ã&nbsp; la constante de stabilitÃ© classique, Ã&nbsp; une constante prÃ¨s qui dÃ©pend de la constante de Babu\v{s}ka--Aziz. Par la suite, nous explicitons les constantes de stabilitÃ© par rapport au paramÃ¨tre de rÃ©gularitÃ© de forme pour l'ordre $1$ en dimension~$2$ ou $3$, et l'ordre $2$ en dimension~$2$. Dans ce dernier cas, nous amÃ©liorons le rÃ©sultat de l'article original de Fortin--Soulie. Ensuite nous illustrons l'importance d'utiliser une mÃ©thode de projection conforme dans $\mathbf{H}(\mathrm{div})$ pour certaines expÃ©riences numÃ©riques.
\end{altabstract}


\begin{document}

\maketitle

\section{Introduction}

The Stokes problem describes the steady state of incompressible Newtonian flows. It follows from the Navier--Stokes equations~\cite{GiRa86}. With regard to numerical analysis, the study of Stokes problem helps to build an appropriate approximation of the Navier--Stokes equations. We consider here a discretization with nonconforming finite elements~\cite{CrRa73,FoSo83}. We propose to state the discrete inf-sup condition in light of the T-coercivity (cf.~\cite{Ciar12} for Helmholtz-like problems, see~\cite{JaCi13,CiJK17,Gire18} for the neutron diffusion equation), which allows to estimate the discrete error constant. In Section~\ref{sec:exact_discrete}, we recall the T-coercivity theory~\cite{Ciar12}, which is known to be an equivalent reformulation of the Banach--Ne\v{c}as--Babu\v{s}ka Theorem and we apply it to the continuous Stokes problem. We give details on the triangulation, and we apply the T-coercivity to the discretization of Stokes problem with nonconforming mixed finite elements. For the Stokes problem, in the discrete case, this amounts to finding a Fortin operator. In Section~\ref{sec:NC-disc}, we precise the proof of the well-posedness in the case of order $1$ and $2$ nonconforming mixed finite elements. In Section~\ref{sec:numerical_results}, we illustrate the importance of using a divergence-free velocity on some numerical experiments.

\section{Exact and discrete T-coercivity for Stokes problem}\label{sec:exact_discrete}

\subsection{T-coercivity and application to Stokes problem}\label{sec:Tcoer}

We recall here the T-coercivity theory as written in~\cite{Ciar12}. Consider first the variational problem, where $V$ and $W$ are two Hilbert spaces and $f\in V'$:
\begin{equation}\label{eq:var-pb}
\text{Find $u\in V$ such that $\forall v\in W$, $a(u,v)=\langle f,v\rangle_V$.}
\end{equation}
Classically, we know that Problem~\eqref{eq:var-pb} is well-posed if $a({\,\cdot\,,\cdot\,})$ satisfies the stability and the solvability conditions of the so-called Banach--Ne\v{c}as--Babu\v{s}ka (BNB) Theorem (see e.g.~\cite[Theorem~25.9]{ErGu21-II}). For some models, one can also prove the well-posedness using the T-coercivity theory (cf.~\cite{Ciar12} for Helmholtz-like problems, see~\cite{JaCi13,CiJK17,Gire18} for the neutron diffusion equation).

\begin{defi}\label{def:Tcoer}
Let $V$ and $W$ be two Hilbert spaces and $a({\,\cdot\,,\cdot\,})$ be a continuous and bilinear form over $V\times W$. It is $T$-coercive if
\begin{equation}\label{eq:Tcoer}
\exists T\in\Lcal(V,W), \; \text{ bijective, } \; \exists\alpha_T&gt;0 \; \forall v\in V, \;\bigl\lvert a(v,Tv)\bigr\rvert\geq\alpha_T\norm{v}_V^2.
\end{equation}
\end{defi}

It is proved in~\cite{Ciar12,ChCi13} that the T-coercivity condition is equivalent to the stability and solvability conditions of the BNB Theorem. Whereas the BNB Theorem relies on an abstract inf-sup condition, T-coercivity uses explicit inf-sup operators, both at the continuous and discrete levels. Notice that if the pair $(T,\alpha_T)$ satisfies~\eqref{eq:Tcoer}, then for all $\lambda&gt;0$, the pair $(T_\lambda,\alpha_{T_\lambda})\coloneqq (\lambda T,\lambda \alpha_T)$ also satisfies~\eqref{eq:Tcoer}. Thus, there exists an infinity of pairs $(T,\alpha_T)$ and $(T_\lambda,\alpha_{T_\lambda})_{\lambda&gt;0}$ satisfying~\eqref{eq:Tcoer}.

\begin{theo}[well-posedness]\label{thm:WellPosed}
Let $a({\,\cdot\,,\cdot\,})$ be a continuous bilinear form. Suppose that the form $a({\,\cdot\,,\cdot\,})$ is $T$-coercive. Then Problem~\eqref{eq:var-pb} is well-posed.
\end{theo}

Let $\Om$ be a connected bounded domain of $\R^d$, $d=2, 3$, with a polygonal $(d=2)$ or Lipschitz polyhedral $(d=3)$ boundary $\pa\Om$. We consider Stokes problem:
\begin{equation}\label{eq:Stokes}
\text{Find $(\uvec,p)$ such that }
\left\{
\begin{aligned}
-\nu\Deltavec\uvec+\grad p &amp;=\fvec,\\
\div\uvec&amp;=0,
\end{aligned}
\right.
\end{equation}
with Dirichlet boundary conditions for the velocity $\uvec$ and a normalization condition for the pressure $p$:
\[
\uvec=0\;\text{ on }\;\pa\Om,\;\int_\Om p=0.
\]
The vector field $\uvec$ represents the velocity of the fluid and the scalar field $p$ represents its pressure divided by the fluid density which is supposed to be constant. The first equation of~\eqref{eq:Stokes} corresponds to the momentum balance equation and the second one corresponds to the conservation of the mass. The constant parameter $\nu&gt;0$ is the kinematic viscosity of the fluid. The vector field $\fvec\in\bH^{-1}(\Om)$ represents a body forces divided by the fluid density.

Before stating the variational formulation of Problem~\eqref{eq:Stokes}, we provide some definitions and reminders. Let us set $\bL^2(\Om)=\bigl(L^2(\Om)\bigr)^d$, $\bH^1_0(\Om)=\bigl(H^1_0(\Om)\bigr)^d$, $\bH^{-1}(\Omega)=\bigl(H^{-1}(\Om)\bigr)^d$ its dual space and $L^2_{\zmv}(\Om)=\bigl\{q\in L^2(\Om)\btq\int_{\Om}q=0\bigr\}$. We recall that $\bH(\div; \Om)=\bigl\{\vvec\in\bL^2(\Om)\btq\div\vvec\in\ L^2(\Om)\bigr\}$. Let $h_\Om$ be the diameter of $\Om$. We recall the PoincarÃ©--Steklov inequality~\cite[Lemma~3.24]{ErGu21-I}:
\begin{equation}\label{eq:Poincare}
\exists C_\PS&gt;0 \; \text{ such that } \; \forall v\in H^1_0(\Om),\; \norm{v}_{L^2(\Om)}\leq C_\PS h_\Om \norm{\grad v}_{\bL^2(\Om)}.
\end{equation}
Thanks to this result, in $H^1_0(\Om)$, the semi-norm is equivalent to the natural norm, so that the scalar product reads $(v,w)_{H^1_0(\Om)}=(\grad v,\grad w)_{\bL^2(\Om)}$ and the norm is $\norm{v}_{H^1_0(\Om)}=\norm{\grad v}_{\bL^2(\Om)}$. Let $\vvec, \wvec\in\bH^1_0(\Om)$. We denote by $(v_i)_{i=1}^d$ (resp.\ $(w_i)_{i=1}^d$) the components of $\vvec$ (resp.\ $\wvec$), and we set $\Grad\vvec=\bigl(\pa_j v_i\bigr)_{i,j=1}^d\in\bbL^2(\Om)$, where $\bbL^2(\Om)=\bigl(L^2(\Om)\bigr)^{d\times d}$. We have:
\[
\bigl(\Grad\vvec,\Grad\wvec\bigr)_{\bbL^2(\Om)}=(\vvec,\wvec)_{\bH^1_0(\Om)}=\sum_{i=1}^d(v_i,w_i)_{H^1_0(\Om)}
\]
and
\[
\norm{\vvec}_{\bH^1_0(\Om)}=\left(\sum_{j=1}^d\norm{v_j}_{H^1_0(\Om)}^2\right)^{\frac12}=\norm{\Grad\vvec}_{\bbL^2(\Om)}.
\]
Let us set $\bV=\bigl\{\vvec\in\bH^1_0(\Om)\btq\div\vvec=0\bigr\}$. The vector space $\bV$ is a closed subset of $\bH^1_0(\Om)$. We denote by $\bV^\perp$ the orthogonal of $\bV$ in $\bH^1_0(\Om)$. We recall the following result.

\begin{prop}[{\cite[Corollary~I.2.4]{GiRa86}}]\label{prop:diviso}
The operator $\div\colon\bH^1_0(\Om)\to L^2(\Om)$ is an isomorphism of\/ $\bV^{\perp}$ onto $L^2_{\zmv}(\Om)$. We call $C_{\div}$ the constant such that:
\begin{equation}\label{eq:BAI}
\forall p\in L^2_{\zmv}(\Om), \mathrel{\exists!}\vvec\in\bV^{\perp}\;\text{ such that }\;\div\vvec=p\;\text{ and }\;
\norm{\vvec}_{\bH^1_0(\Om)}\leq C_{\div}\norm{p}_{L^2(\Om)}.
\end{equation}
\end{prop}

The constant $C_{\div}$ depends only on the domain $\Omega$. Recall that we have
\[
\norm{\vvec}_{\bH^1_0(\Om)}^2=\norm{\curlv\vvec}_{\bL^2(\Om)}^2+\norm{\div\vvec}_{L^2(\Om)}^2\geq\norm{p}_{L^2(\Om)}^2,
\]
hence $C_{\div}\geq 1$. Using Proposition~\ref{prop:diviso}, we can define a lifting in $\mathbf{H}^1(\Omega)$ of a function $v\in L^2(\Omega)$:
\begin{coro}\label{coro:svec}
For all $v\in L^2(\Om)$, there exists $\svec\in\bH^1(\Om)$ such that:
\begin{equation}\label{eq:svec}
\div\svec=v\quad\text{ and }\quad\norm{\svec}_{\bL^2(\Om)}+h_\Om\norm{\Grad\svec}_{\bbL^2(\Om)}\leq C_\Om h_\Om\norm{v}_{L^2(\Om)},
\end{equation}
where the dimensionless constant $C_\Om$ depends on $C_{\div}$ and $C_\PS$.
\end{coro}
\begin{proof}
Let $\ul{v}=\int_\Om v/\abs{\Om}$, and $\svec_0\in\bH^1_0(\Om)$ be such that $\div\svec_0=v-\ul{v}$ and $\norm{\Grad\svec_0}_{\bbL^2(\Om)}\leq C_{\div}\norm{v}_{L^2(\Om)}$ (cf.\ Proposition~\ref{prop:diviso}). Let $d'\in\bigl\{1,\dotsc,d\bigr\}$. We consider $\svec\coloneqq \svec_0+\ul{v}(x_{d'}-\ul{x}_{d'})\evec_{d'}\in\bH^1(\Om)$, where $\ul{x}_{d'}=\int_\Om x_{d'}/\abs{\Om}$. We have: $\div\svec=v$ and $\norm{\Grad\svec}_{\bbL^2(\Om)}\leq\widetilde{C}_{\div}\norm{v}_{L^2(\Om)}$, where $\widetilde{C}_{\div}=C_{\div}+1$. Using inequality~\eqref{eq:Poincare}, one can prove that $\norm{\svec}_{\bL^2(\Om)}\leq h_\Om \widetilde{C}_\Om\norm{v}_{L^2(\Om)}$, where $\widetilde{C}_\Omega=C_{\PS}\,C_{\div}+1$. Setting $C_\Om=\widetilde{C}_\Om+\widetilde{C}_{\div}$, we obtain~\eqref{eq:svec}.
\end{proof}

Actually, the constant $C_{\div}$ is such that $C_{\div}=1/\beta(\Om)$ where $\beta(\Om)$, known as the Babu\v{s}ka--Aziz constant, is the inf-sup condition (or Ladyzhenskaya--Babu\v{s}ka--Brezzi condition):
\begin{equation}\label{eq:LBB-CIS}
\beta(\Om)=\inf_{q\in L^2_{\zmv}(\Om)\backslash \{0\}}\sup_{\vvec\in\bH^1_0(\Om)\backslash \{0\}}\frac{(q,\div\vvec)_{L^2(\Om)}}{\norm{q}_{L^2(\Om)} \norm{\vvec}_{\bH^1_0(\Om)}}.
\end{equation}
Generally, the value of $\beta(\Om)$ is not known explicitly. In~\cite{BCDG16}, Bernardi et al.\ established results on the discrete approximation of $\beta(\Om)$ using conforming finite elements. Recently, Gallistl proposed in~\cite{Gall19} a numerical scheme with adaptive meshes for computing approximations to $\beta(\Om)$. In the case of $d=2$, Costabel and Dauge~\cite{CoDa15} established the following bound:
\begin{theo}[{\cite[Theorem~2.3]{CoDa15}}]
Let $\Om\subset\R^2$ be a domain contained in a ball of radius $R$, star-shaped with respect to a concentric ball of radius $\rho$. Then
\begin{equation}\label{eq:betaOm}
\beta(\Om)\geq\frac{\rho}{\sqrt{2} R}\left(1+\sqrt{1-\frac{\rho^2}{R^2}}\right)^{-\frac12}\geq\frac{\rho}{2 R}.
\end{equation}
\end{theo}

Let us detail the bound for some remarkable domains. If $\Om$ is a ball, $\beta(\Om)\geq\frac{1}{2}$ and if $\Om$ is a square, $\beta(\Om)\geq\frac{1}{2 \sqrt{2}}$. Suppose now that $\Om$ is stretched in some direction by a factor $k$, then $\beta(\Om)\geq\frac{1}{2 k}$. Finally, if $\Om$ is L-shaped (resp.\ cross-shaped) such that $L=k l$, where $L$ is the largest length and $l$ is the smallest length of an edge, then $\beta(\Om)\geq\frac{1}{2 \sqrt{2} k}$ (resp.\ $\beta(\Om)\geq\frac{1}{4 k}$).

The variational formulation of Problem~\eqref{eq:Stokes} reads:
\begin{multline}\label{eq:Stokes-FV}
\text{Find $(\uvec,p)\in\bH^1_0(\Om)\times L^2_{\zmv}(\Om)$ such that}\\
\left\{
\begin{aligned}
\nu(\uvec,\vvec)_{\bH^1_0(\Om)}-(p,\div\vvec)_{L^2(\Om)}&amp;=\langle\fvec,\vvec\rangle_{\bH^1_0(\Om)}									&amp;&amp;\forall\vvec\in\bH^1_0(\Om),\\
(q,\div\uvec)_{L^2(\Om)}&amp;=0	&amp;&amp;\forall q\in L^2_{\zmv}(\Om).
\end{aligned}
\right.
\end{multline}
Classically, one proves that Problem~\eqref{eq:Stokes-FV} is well-posed using PoincarÃ©--Steklov inequality~\eqref{eq:Poincare} and Proposition~\ref{prop:diviso}. Check for instance the proof of~\cite[Theorem~I.5.1]{GiRa86}.

Let us set $\Xcal=\bH^1_0(\Om)\times L^2_{\zmv}(\Om)$, which is a Hilbert space which we endow with the following norm:
\begin{equation}\label{eq:normX}
\norm*{(\vvec,q)}_{\Xcal}=\biggl(\norm{\vvec}_{\bH^1_0(\Om)}^2+\nu^{-2} \norm{q}_{L^2(\Om)}^2\biggr)^{\frac12}.
\end{equation}
We consider now the following bilinear symmetric and continuous form:
\begin{equation}\label{eq:aform}
a_S\colon
\left\{
\begin{aligned}
\Xcal\times\Xcal			&amp;\to \R\\
(\uvec',p')\times(\vvec,q)	&amp;\mapsto \nu(\uvec',\vvec)_{\bH^1_0(\Om)}-(p',\div\vvec)_{L^2(\Om)}-(q,\div\uvec')_{L^2(\Om)}.
\end{aligned}
\right.
\end{equation}
We can write Problem~\eqref{eq:Stokes} in an equivalent way as follows:
\begin{equation}\label{eq:Stokes-FVa}
\text{Find $(\uvec,p)\in\Xcal$ such that $a_S \bigl((\uvec,p),(\vvec,q)\bigr)=\langle\fvec,\vvec\rangle_{\bH^1_0(\Om)}\;\forall(\vvec,q)\in\Xcal$.}
\end{equation}
Let us prove that Problem~\eqref{eq:Stokes-FVa} is well-posed using the T-coercivity theory.

\begin{prop}\label{pro:aSTcoer}
The bilinear form $a_S({\,\cdot\,,\cdot\,})$ is $T$-coercive:
\begin{equation}\label{eq:aSTcoer}
\begin{gathered}
\exists T\in\Lcal(\Xcal),\;\text{ bijective, }\;\exists\alpha_T&gt;0, \forall(\uvec',p')\in\Xcal,\\
a_S\Bigl((\uvec',p'), T\bigl((\uvec',p')\bigr)\Bigr)\geq \alpha_T\norm*{(\uvec',p')}_{\Xcal}^2.
\end{gathered}
\end{equation}
\end{prop}
\begin{proof}
We follow here the proof given in~\cite{BoCi23,BaCi24,CiJa24}. Let us consider $(\uvec',p')\in\Xcal$ and let us build $(\vvec^\star,q^\star)=T(\uvec',p')\in\Xcal$ satisfying~\eqref{eq:Tcoer} (with $V=\Xcal$). We need three main steps.
\begin{stepA}
According to Proposition~\ref{prop:diviso}, there exists $\tilde{\vvec}_{p'}\in\bH^1_0(\Om)$ such that: $\div\tilde{\vvec}_{p'}=p'$ in $\Om$ and $\norm{\tilde{\vvec}_{p'}}_{\bH^1_0(\Om)}\leq C_{\div}\norm{p'}_{L^2(\Om)}$. Let us set $\vvec_{p'}=\nu^{-1}\tilde{\vvec}_{p'}$ so that $\div\vvec_{p'}=\nu^{-1} p'$ and
\begin{equation}\label{eq:GiRa}
\norm{\vvec_{p'}}_{\bH^1_0(\Om)}\leq \nu^{-1} C_{\div} \norm{p'}_{L^2(\Om)}.
\end{equation}
Let us set $(\vvec^\star,q^\star)\coloneqq \bigl(\gamma \uvec'-\vvec_{p'},-\gamma p'\bigr)$, with $\gamma&gt;0$. We obtain:
\begin{equation}\label{eq:coer1}
a_S \bigl((\uvec',p'), (\vvec^\star,q^\star)\bigr)=
\nu \gamma \norm{\uvec'}_{\bH^1_0(\Om)}^2+\nu^{-1} \norm{p'}_{L^2(\Om)}^2-\nu (\uvec',\vvec_{p'})_{\bH^1_0(\Omega)}.
\end{equation}
\end{stepA}
\begin{stepA}
In order to bound the last term of~\eqref{eq:coer1}, we use Young inequality and then inequality~\eqref{eq:GiRa}, so that for all $\eta&gt;0$:
\begin{equation}\label{eq:GvGvp}
(\uvec',\vvec_{p'})_{\bH^1_0(\Omega)}\leq
\frac{\eta}{2}\norm{\uvec'}_{\bH^1_0(\Omega)}^2+\frac{\eta^{-1}}{2} \Biggl(\frac{C_{\div}}{\nu}\Biggr)^2 \norm{p'}_{L^2(\Om)}^2.
\end{equation}
\end{stepA}
\begin{stepA}
Using the bound~\eqref{eq:GvGvp} in~\eqref{eq:coer1} and choosing $\eta=\gamma$, we get:
\[
a_S \bigl((\uvec',p'), (\vvec^\star,q^\star)\bigr)\geq\nu\Biggl(\frac{\gamma}{2} \norm{\uvec'}_{\bH^1_0(\Omega)}^2+\nu^{-2} \Biggl(1-\frac{\gamma^{-1}}{2}(C_{\div})^2\Biggr)\norm{p'}_{L^2(\Om)}^2\Biggr).
\]
Consider now $\gamma=(C_{\div})^2$. We obtain:
\[
a_S \bigl((\uvec',p'), (\vvec^\star,q^\star)\bigr)\geq\frac{1}{2}\nu\Bigl((C_{\div})^2 \norm{\uvec'}_{\bH^1_0(\Omega)}^2+\nu^{-2} \norm{p'}_{L^2(\Om)}^2\Bigr).
\]
Reminding that $C_{\div}\geq 1$, it comes:
\begin{equation}\label{eq:Cmin}
a_S \bigl((\uvec',p'), (\vvec^\star,q^\star)\bigr)\geq \frac{\nu}{2} \norm*{(\uvec',p')}_{\Xcal}^2.
\end{equation}
\end{stepA}
We obtain~\eqref{eq:aSTcoer} with $\alpha_T=\frac{\nu}{2}$.

The operator $T$ such that $T\bigl((\uvec',p')\bigr)=(\vvec^\star,q^\star)$ is linear and continuous. We have indeed:
\[
\begin{aligned}
\norm*{T\bigl((\uvec',p')\bigr)}_{\Xcal}^2&amp;\coloneqq \norm{\vvec^\star}_{\bH^1_0(\Om)}^2+\nu^{-2} \norm{q^\star}_{L^2(\Omega)}^2\\
&amp;\:\leq 2 \gamma^2 \norm{\uvec'}_{\bH^1_0(\Om)}^2+2 \norm{\vvec_{p'}}_{\bH^1_0(\Om)}^2+\gamma^2 \nu^{-2} \norm{p'}_{L^2(\Omega)}^2\\
&amp;\:\leq 2 \gamma^2 \norm{\uvec'}_{\bH^1_0(\Omega)}^2+ ( 2 (C_{\div})^2+\gamma^2) \nu^{-2} \norm{p'}_{L^2(\Omega)}^2.
\end{aligned}
\]
We deduce that:
\begin{equation}\label{eq:Tbound}
\norm*{T\bigl((\uvec',p')\bigr)}_{\Xcal}\leq C_{\max} \norm*{(\uvec',p')}_{\Xcal}\;\text{ where }\;C_{\max}=C_{\div} \Bigl(\max\bigl(2+(C_{\div})^2,2 (C_{\div})^2\bigr)\Bigr)^{\frac12}.
\end{equation}
Remark that, given $(\vvec^\star,q^\star)\in\Xcal$, choosing $(\uvec',p')=\bigl(\gamma^{-1}\vvec^\star-\gamma^{-2}\vvec_{q^\star},-\gamma^{-1}q^\star\bigr)$ yields $T\bigl((\uvec',p')\bigr)=(\vvec^\star,q^\star)$. Hence, the operator $T\in\Lcal(\Xcal)$ is bijective.
\end{proof}

We can now prove the following result.
\begin{theo}\label{thm:StokesWellPosed}
Problem~\eqref{eq:Stokes-FVa} is well-posed. It admits one and only one solution such that:
\begin{equation}\label{eq:WellPosed}
\forall\fvec\in\bH^{-1}(\Omega),\;
\left\{
\begin{aligned}
\norm{\uvec}_{\bH^1_0(\Om)}	&amp;\leq\nu^{-1} \norm{\fvec}_{\bH^{-1}(\Omega)},\\
\norm{p}_{L^2(\Om)}			&amp;\leq C_{\div} \norm{\fvec}_{\bH^{-1}(\Omega)}.\\
\end{aligned}
\right.
\end{equation}
\end{theo}
\begin{proof}
According to Proposition~\ref{pro:aSTcoer}, the continuous bilinear form $a_S({\,\cdot\,,\cdot\,})$ is $T$-coercive. Hence, according to Theorem~\ref{thm:WellPosed}, Problem~\eqref{eq:Stokes-FVa} is well-posed. Let us prove~\eqref{eq:WellPosed}. Consider $(\uvec,p)$ the unique solution of Problem~\eqref{eq:Stokes-FVa}. Choosing $\vvec=0$, we obtain that $\forall q\in L^2_{\zmv}(\Om)$, $(q,\div\uvec)_{L^2(\Om)}=0$, so that $\uvec\in\bV$. Now, choosing $\vvec=\uvec$ and using Cauchy--Schwarz inequality, we have: $\nu \norm{\uvec}_{\bH^1_0(\Om)}^2=\langle\fvec,\uvec\rangle_{\bH^1_0(\Om)}\leq\norm{\fvec}_{\bH^{-1}(\Om)} \norm{\uvec}_{\bH^1_0(\Om)}$, so that: $\norm{\uvec}_{\bH^1_0(\Om)}\leq\nu^{-1} \norm{\fvec}_{\bH^{-1}(\Om)}$. Next, we choose in~\eqref{eq:Stokes-FVa} $\vvec=\tilde{\vvec}_p\in\bV^{\perp}$, where $\div\tilde{\vvec}_p=-p$ (see Proposition~\ref{prop:diviso}). Since $\uvec\in\bV$ and $\tilde{\vvec}_p\in\bV^\perp$, we have $(\uvec,\tilde{\vvec}_p)_{\bH^1_0(\Om)}=0$. This gives:
\[
-(p,\div\tilde{\vvec}_p)_{L^2(\Om)}=\norm{p}_{L^2(\Om)}^2=\langle\fvec,\tilde{\vvec}_p\rangle_{\bH^1_0(\Om)}\leq\norm{\fvec}_{\bH^{-1}(\Om)} \norm{\tilde{\vvec}_p}_{\bH^1_0(\Om)}\leq C_{\div} \norm{\fvec}_{\bH^{-1}(\Om)} \norm{p}_{L^2(\Om)},
\]
so that: $\norm{p}_{L^2(\Om)}\leq C_{\div}\norm{\fvec}_{\bH^{-1}(\Om)}$.
\end{proof}

\subsection{Comments on the stability constant}

Using~\eqref{eq:Tbound} in~\eqref{eq:Cmin}, we have:
\[
a_S \bigl((\uvec',p'), (\vvec^\star,q^\star)\bigr)\geq \alpha_T (C_{\max})^{-1}\norm*{(\uvec',p')}_{\Xcal} \norm*{(\vvec^\star,q^\star)}_{\Xcal}
\]
Let us set $\vertiii{T}\coloneqq \sup_{(\uvec',p')\in\Xcal\backslash(\zerovec,0)}\frac{\norm*{T((\uvec',p'))}_{\mathcal{X}}}{\norm*{(\uvec',p')}_{\Xcal}}$. According to~\eqref{eq:Tbound}, we have the bound: $\vertiii{T}\leq C_{\max}$.

Hence, we have:
\[
a_S \bigl((\uvec',p'), (\vvec^\star,q^\star)\bigr)\geq \frac{\alpha_T}{\vertiii{T}}\norm*{(\uvec',p')}_{\Xcal} \norm*{(\vvec^\star,q^\star)}_{\Xcal}
\]
We recover the first Banach--Ne\v{c}as--Babu\v{s}ka condition~\cite[Theorem~25.9, (BNB1)]{ErGu21-II}.

Thus, the T-coercivity approach gives an overestimate of the stability constant $\alpha$ given below:
\[
\frac{\alpha_T}{\vertiii{T}}\geq \alpha\coloneqq \inf_{(\vvec,q)\in\Xcal\backslash(\zerovec,0)} \sup_{(\uvec',p')\in\Xcal\backslash(\zerovec,0)}\frac{a_S \bigl((\uvec',p'), (\vvec,q)\bigr)}{\norm*{(\uvec',p')}_{\Xcal} \norm*{(\vvec,q)}_{\Xcal}}.
\]
It suggests that the stability constant $\alpha$ is proportional to the parameter $\nu$ and depends on the constant $C_{\div}$, therefore on the shape of the domain. More precisely, our estimate gives:
\[
\alpha\leq\frac{\nu}{2}\times
\begin{dcases}
\left(C_{\div} \sqrt{2+(C_{\div})^2}\right)^{-1}&amp;\text{ if }\;1\leq C_{\div}\leq\sqrt{2},\\
\biggl(\sqrt{2} (C_{\div})^2\biggr)^{-1}&amp;\text{ if }\;C_{\div}\geq\sqrt{2}.
\end{dcases}
\]
In our computations, $\alpha_T$ depends on the choice of the parameters $\eta$ and $\gamma$, so that it could be further optimized to minimize $\frac{\alpha_T}{\vertiii{T}}$. Studying the bilinear form $a_S\Bigl((\uvec',p'),T\bigl((\vvec,q)\bigr)\Bigr)$ leads to an alternative variational formulation of Stokes problem, as proposed in~\cite{CiJa24}. It does not depend on the parameters $\eta$ and $\gamma$ because it is coercive. However, the new variational formulation requires a specific treatment of the right-hand side.

\subsection{Conforming discretization and discrete well-posedness}\label{sec:cdisc}

If we were using a conforming discretization to solve Problem~\eqref{eq:Stokes-FVa} (e.g.\ Taylor--Hood finite elements~\cite{TaHo73}), we would use the bilinear form $a_S({\,\cdot\,,\cdot\,})$ to state the discrete variational formulation. Let us call the discrete spaces $\bX_{\conf,h}\subset\bH^1_0(\Om)$ and $Q_{\conf,h}\subset L^2_{\zmv}(\Om)$. Then to prove the discrete T-coercivity, we would need to state the discrete counterpart to Proposition~\ref{prop:diviso}. To do so, we can build a linear operator $\Pi_{\conf}\colon\bX\to\bX_{\conf,h}$, known as Fortin operator, such that (see e.g.~\cite[\S 8.4.1]{BoBF13}):
\begin{align}\label{eq:WellPosed-c-1}
\exists C_\conf\;\text{ such that }\;
\forall\vvec\in\bH^1(\Om),\;\norm{\Grad\Pi_\conf\vvec}_{\bbL^2(\Om)}&amp;\leq C_\conf\norm{\Grad\vvec}_{\bbL^2(\Om)},\\\label{eq:WellPosed-c-2}
\forall\vvec\in\bH^1(\Om),\;(\div\Pi_\conf\vvec,q_h)_{L^2(\Om)}&amp;=(\div\vvec,q_h)_{L^2(\Om)},\;\forall q_h\in Q_{\conf,h}.
\end{align}
Using a nonconforming discretization, we will not use the bilinear form $a_S({\,\cdot\,,\cdot\,})$ to exhibit the discrete variational formulation, but we will need a similar operator to~\eqref{eq:WellPosed-c-1}-\eqref{eq:WellPosed-c-2} to prove the discrete T-coercivity, which is stated in Theorem~\ref{thm:WellPosed-nc}.

\subsection{Discretization notations}\label{sec:disc}

We call $\bigl(O,(x_{d'})_{d'=1}^d\bigr)$ the Cartesian coordinates system, of orthonormal basis $(\evec_{d'})_{d'=1}^d$. Consider $(\Tcal_h)_h$ a simplicial triangulation sequence of $\Om$, where $h$ denotes the mesh size. The triangulations are regular in the sense of Ciarlet. For a triangulation $\Tcal_h$, we use the following index sets:
\begin{itemize}
\item $\Ical_K$ denotes the index set of the elements, such that $\Tcal_h\coloneqq \bigcup_{\ell\in\Ical_K}K_\ell$ is the set of elements.
\item $\Ical_F$ denotes the index set of the facets\footnote{The term facet stands for face (resp.\ edge) when $d=3$ (resp.\ $d=2$).}, such that $\Fcal_h\coloneqq \bigcup_{f\in\Ical_F}F_f$ is the set of facets.
\item[]Let $\Ical_F=\Ical_F^i\cup\Ical_F^b$, where $\forall f\in\Ical_F^i$, $F_f\subset\Om$ and $\forall f\in\Ical_F^b$, $F_f\in\pa\Om$.
\item $\Ical_S$ denotes the index set of the vertices, such that $(S_j)_{j\in\Ical_S}$ is the set of vertices.
\item[]Let $\Ical_S=\Ical_S^i\cup\Ical_S^b$, where $\forall j\in\Ical_S^i$, $S_j\in\Om$ and $\forall j\in\Ical_S^b$, $S_j\in\pa\Om$.
\end{itemize}
We also define the following index subsets:
\begin{itemize}
\item $\forall\ell\in\Ical_K$, $\Ical_{F,\ell}^{(i,b)}=\bigl\{f\in\Ical_F^{(i,b)} \btq F_f\in K_\ell\bigr\},\;\Ical_{S,\ell}=\bigl\{j\in\Ical_S\btq S_j\in K_\ell\bigr\}$.
\item $\forall j\in\Ical_S$, $\Ical_{K,j}=\bigl\{\ell\in\Ical_K\btq S_j\in K_\ell\bigr\}$, $\; N_j\coloneqq \card(\Ical_{K,j})$.
\end{itemize}
For all $\ell\in\Ical_K$, we call $h_\ell$ the diameter of $K_\ell$ and $\rho_\ell$ the diameter of the sphere inscribed in $K_\ell$, and we let: $\sigma_\ell=\frac{h_\ell}{\rho_\ell}$, $h=\max_{\ell\in\Ical_K}h_\ell$. When $(\Tcal_h)_h$ is a shape-regular triangulation sequence (see e.g.~\cite[Definition~11.2]{ErGu21-I}), there exists a constant $\sigma&gt;1$, called the shape regularity parameter, such that for all $h$, for all $\ell\in\Ical_K$, $\sigma_\ell\leq\sigma$. For all $f\in\Ical_F$, $M_f$ denotes the barycenter of $F_f$, and by $\nvec_f$ its unit normal (outward oriented if $F_f\in\pa\Omega$). For all $j\in\Ical_S$, for all $\ell\in\Ical_{K,j}$, $\lambda_{j,\ell}$ denotes the barycentric coordinate of $S_j$ in $K_\ell$, $F_{j,\ell}$ denotes the face opposite to vertex $S_j$ in element $K_\ell$, and $\xvec_{j,\ell}$ denotes its barycenter. We call $\Scal_{j,\ell}$ the outward normal vector of $F_{j,\ell}$ and of norm $\abs{\Scal_{j,\ell}}=\abs{F_{j,\ell}}$.

Let us introduce spaces of piecewise regular elements.

We set $\Pcal_h H^1=\bigl\{v\in L^2(\Om) \btq \forall \ell\in\Ical_K,\;\restr{v}{K_\ell}\in H^1(K_\ell)\bigr\}$, endowed with the scalar product:
\[
(v,w)_h\coloneqq \sum_{\ell\in\Ical_K}(\grad v,\grad w)_{\bL^2(K_\ell)},\quad
\norm{v}_h^2=\sum_{\ell\in\Ical_K}\norm{\grad v}^2_{\bL^2(K_\ell)}.
\]

We set $\Pcal_h\bH^1=\bigl(\Pcal_hH^1\bigr)^d$, endowed with the scalar product:
\[
(\vvec,\wvec)_h\coloneqq \sum_{\ell\in\Ical_K}(\Grad\vvec,\Grad\wvec)_{\bbL^2(K_\ell)},\quad
\norm{\vvec}_h^2=\sum_{\ell\in\Ical_K}\norm{\Grad \vvec}^2_{\bbL^2(K_\ell)}.
\]
Let $f\in\Ical_F^i$ such that $F_f=\pa K_L\cap\pa K_R$ and $\nvec_f$ is outward $K_L$ oriented. The jump of $v\in \Pcal_h H^1$ across the facet $F_f$ is defined as follows: $[v]_{F_f}\coloneqq \restr{v}{K_L}-\restr{v}{K_R}$. For $f\in\Ical_F^b$, we set: $[v]_{F_f}\coloneqq \restr{v}{F_f}$.

We set $\Pcal_h\bH(\div)=\bigl\{\vvec\in \bL^2(\Omega) \btq \forall\ell\in\Ical_K, \restr{\vvec}{K_\ell}\in\bH(\div; K_\ell)\bigr\}$, and we define the operator $\div_h$ such that for all $\vvec\in\Pcal_h\bH(\div)$, $\div_h\vvec\in L^2(\Om)$ is such that:
\[
\forall q\in L^2(\Om),\; \bigl(\div_h\vvec,q\bigr)_{L^2(\Om)}=\sum_{\ell\in\Ical_K}\bigl(\div\vvec,q\bigr)_{L^2(K_\ell)}.
\]

We use the notation $A\lesssim B$ for $A\lesssim CB$ where $A$ and $B$ are scalar quantities and $C$ is a generic positive constant which is independent of the sequence $(\Tcal_h)_h$ and the quantities of interest.

We recall classical finite elements estimates~\cite{ErGu21-I}.
Let $\hatK$ be the reference simplex. For $\ell\in\Ical_K$, we denote by $T_\ell\colon\hatK\to K_\ell$ an affine invertible mapping such that $T_\ell(\hatK)=K_\ell$, $T_\ell(\pa\hatK)=\pa K_\ell$. We set $T_\ell(\hat{\xvec})=\bbB_\ell\hat{\xvec}+\bvec_\ell$, where $\bbB_\ell\in\R^{d\times d}$ and $\bvec_\ell\in\R^d$. Let $J_\ell=\det(\bbB_\ell)$. There holds:
\begin{equation}\label{eq:ch-vb}
\abs{J_\ell}=d!\abs{K_\ell},\quad\norm{\bbB_\ell}=\frac{h_\ell}{\rho_{\hatK}},\quad\norm{{\bbB_\ell}^{-1}}=\frac{h_{\hatK}}{\rho_\ell}.
\end{equation}
Let $f\in\Ical_{F,\ell}$. According to~\cite[Equation~(2.17)]{CrRa73}, we have:
\begin{equation}\label{eq:CrRa73-3.17}
\abs{F_f} \abs{K_\ell}^{-1}\lesssim(\rho_\ell)^{-1}.
\end{equation}
For $v\in L^2(K_\ell)$, we set $\hat{v}_\ell=v\circ T_\ell$. Let $v\in \Pcal_h H^1$. By changing the variable, $\grad \restr{v}{K_\ell}=({\bbB_\ell}^{-1})^T\grad_{\hat{\xvec}}\hat{v}_\ell$, and it holds:
%\begin{align}
%\norm{\grad v}_{\bL^2(K_\ell)}^2	&amp;\leq\norm{{\bbB_\ell}^{-1}}^2 \abs{K_\ell} \norm{\grad_{\hat{\xvec}}\hat{v}_\ell}_{\bL^2(\hat{K})}^2,\label{eq:gradvK}\\
%\norm{\grad_{\hat{\xvec}}\hat{v}_\ell}_{\bL^2(\hat{K})}^2	&amp;\leq\norm{\bbB_\ell}^2 \abs{K_\ell}^{-1} \norm{\grad v}_{\bL^2(K_\ell)}^2.
%\end{align}
%\begin{equation}\label{eq:gradvK}
%\norm{\grad v}_{\bL^2(K_\ell)}^2	\lesssim\norm{{\bbB_\ell}^{-1}}^2 \abs{K_\ell} \norm{\grad_{\hat{\xvec}}\hat{v}_\ell}_{\bL^2(\hat{K})}^2.
%\end{equation}
%\begin{eqnarray}
%\label{eq:gradvK-i}
%\|\grad v\|_{\mathbf{L}^2(K_\ell)}^2
%&amp;\lesssim&amp;
%\|{\mathbb{B}_\ell}^{-1}\|^2\,|K_\ell|\,\|\grad_{\hat{\mathbf{x}}}\hat{v}_\ell\|_{\mathbf{L}^2(\hat{K})}^2,\\
%\label{eq:gradvK-ii}
%\|\grad_{\hat{\mathbf{x}}}\hat{v}_\ell\|_{\bL^2(\hat{K})}^2
%&amp;\lesssim&amp;
%\|\mathbb{B}_\ell\|^2\,|K_\ell|^{-1}\,\|\grad v\|_{\mathbf{L}^2(K_\ell)}^2.
%\end{eqnarray}
\begin{align}
\norm{\grad v}_{\bL^2(K_\ell)}^2	&amp;\lesssim\norm{{\bbB_\ell}^{-1}}^2 \abs{K_\ell} \norm{\grad_{\hat{\xvec}}\hat{v}_\ell}_{\bL^2(\hat{K})}^2,\label{eq:gradvK-i}\\
\norm{\grad_{\hat{\mathbf{x}}}\hat{v}_\ell}_{\bL^2(\hat{K})}^2
&amp;\lesssim
\norm{\mathbb{B}_\ell}^2\,\abs{K_\ell}^{-1}\,\norm{\grad v}_{\mathbf{L}^2(K_\ell)}^2.\label{eq:gradvK-ii}
\end{align}
We will use the following notations:
\begin{gather}\label{eq:vl}
\forall\ell\in\Ical_K, \; \forall v\in L^2(K_\ell), \; \ul{v}_\ell=\int_{K_\ell}v/\abs{K_\ell},\\\label{eq:vf}
\forall f\in\Ical_F, \; \forall v\in L^2(F_f),\; \ul{v}_f=\int_{F_f}v/\abs{F_f}.
\end{gather}
We recall the PoincarÃ©--Steklov inequality in cells:

\begin{prop}[{\cite[Lemma~12.11]{ErGu21-I}}]\label{pro:PSIinC}
For all $\ell\in\Ical_K$ ($K_\ell$ is a convex set), for all $v\in H^1(K_\ell)$:
\begin{equation}\label{eq:PSinC}
\norm{v-\ul{v}_\ell}_{L^2(K_\ell)}\leq\pi^{-1}h_\ell\norm{\grad v}_{\bL^2(K_\ell)}.
\end{equation}
\end{prop}

We will need the following PoincarÃ©--Steklov inequality on faces:

\begin{prop}[{\cite[Lemma~36.8]{ErGu21-II}}]\label{pro:PSIonFaces}
For all $\ell\in\Ical_K$, for all $v\in H^1(K_\ell)$ and for all $f\in\Ical_{F,\ell}$, we have:
\begin{equation}\label{eq:MTI}
\norm{v-\ul{v}_f}_{L^2(F_f)}\lesssim
\Biggl(\frac{\abs{F_f}}{\abs{K_\ell}}\Biggr)^{\frac12} h_\ell \norm{\grad v}_{\bL^2(K_\ell)}
\lesssim(\sigma_\ell)^{\frac{1}{2}}  (h_\ell)^{\frac{1}{2}} \norm{\grad v}_{\bL^2(K_\ell)}.
\end{equation}
\end{prop}
\begin{proof}
We have: $v-\ul{v}_f=(v-\ul{v}_\ell)-\abs{F}^{-1} \int_{F_f}(v-\ul{v}_\ell)$. Hence: $\norm{v-\ul{v}_f}_{L^2(F_f)}\leq 2\norm{v-\ul{v}_\ell}_{L^2(F_f)}$. Changing the variable, using the continuity of the trace operator, we have: $\norm{v-\ul{v}_\ell}_{L^2(F_f)}\lesssim\abs{F_f}^{\frac12}\norm{\hat{v}_\ell-\hat{\ul{v}}_\ell}_{H^1(\hat{K})}$. %Using~\eqref{eq:PSinC} in $\hat{K}$, changing the variable again, we get~\eqref{eq:MTI}, applying~\eqref{eq:CrRa73-3.17}.
Using~\eqref{eq:PSinC} in $\hat{K}$, we obtain $\norm{v-\underline{v}_\ell}_{L^2(F_f)}\lesssim\abs{F_f}^{\frac12}\norm{\grad_{\hat{\mathbf{x}}}\hat{v}_\ell}_{\mathbf{L}^2(\hat{K})}$. Using~\eqref{eq:gradvK-ii} and applying~\eqref{eq:CrRa73-3.17}, we get~\eqref{eq:MTI}.
\end{proof}

For all $D\subset\R^d$, and $k\in\N$, we call $P^k(D)$ the set of order $k$ polynomials on $D$, $\bP^k(D)=\bigl(P^k(D)\bigr)^d$, and we consider the broken polynomial space:
\[
P^k_\disc(\Tcal_h)=\bigl\{q\in L^2(\Om) \btq \forall\ell\in\Ical_K,\;\restr{q}{K_\ell}\in P^k(K_\ell)\bigr\},\quad\bP^k_\disc(\Tcal_h)\coloneqq \bigl(P^k_\disc(\Tcal_h)\bigr)^d.
\]

\subsection{Nonconforming discretization and discrete well-posedness}\label{sec:ncdisc}

The nonconforming finite element method was introduced by Crouzeix and Raviart in~\cite{CrRa73} to solve Stokes problem~\eqref{eq:Stokes}. We approximate the vector space $\bH^1(\Om)$ component by component by piecewise polynomials of order $k\in\N^\star$. Let us consider $X_h$ (resp.\ $X_{0,h}$), the space of nonconforming approximation of $H^1(\Om)$ (resp.\ $H^1_0(\Om)$) of order $k$:
\begin{equation}
\label{eq:Xh}
\begin{gathered}
X_h=\Biggl\{v_h\in P^k_\disc(\Tcal_h) \,\Bigg\vert\, \forall f\in\Ical_F^i, \; \forall q_h\in P^{k-1}(F_f), \; \int_{F_f}[v_h]_{F_f} q_h=0\Biggr\},\\
X_{0,h}=\Biggl\{v_h\in X_h \,\Bigg\vert\, \forall f\in\Ical_F^b, \; \forall q_h\in P^{k-1}(F_f), \; \int_{F_f} v_h q_h=0\Biggr\}.
\end{gathered}
\end{equation}

The condition on the jumps of $v_h$ on the inner facets is often called the patch-test condition. It allows to prove a discrete PoincarÃ©--Steklov inequality, using Corollary~\ref{coro:svec}, and the proof of~\cite[Theorem~D.1]{Saut22}. The proof of~\cite[Lemma~36.6]{ErGu21-II} is similar, but the vector $\svec$ defined in Corollary~\ref{coro:svec} is constructed in~\cite[Lemma~36.6]{ErGu21-II} as the gradient of a scalar function, so that it gives a lower estimate when $\Om$ is nonconvex. Alternative proofs are given in~\cite{Bren03,Vohr05}.

\begin{prop}\label{pro:cps-constant}
The following discrete PoincarÃ©--Steklov inequality holds:
\begin{equation}\label{eq:cps-constant}
\forall v_h\in X_{0,h},\;\norm{v_h}_{L^2(\Omega)}\lesssim \sigma C_\Om h_\Om \norm{v_h}_h.
\end{equation}
\end{prop}
\begin{proof}
Let $v_h\in X_{0,h}$. According to Corollary~\ref{coro:svec}, there exists $\svec\in \bH^1(\Om)$ such that:
\begin{equation}\label{eq:svecvh}
\div\svec=v_h\qquad\text{ and }\qquad\norm{\svec}_{\bL^2(\Om)}+h_\Om\norm{\Grad\svec}_{\bbL^2(\Om)}\leq C_\Om h_\Om\norm{v_h}_{L^2(\Om)}
\end{equation}
We have, by integration by parts:
\begin{equation}\label{eq:cps-constant-0}
\norm{v_h}_{L^2(\Om)}^2=\bigl(v_h,\div\svec\bigr)_{L^2(\Om)}=-\sum_{\ell\in\Ical_K}\bigl(\grad v_h,\svec\bigr)_{\bL^2(K_\ell)}+\sum_{\ell\in\Ical_K}\sum_{f\in\Ical_{F,\ell}}\bigl(v_h,\svec\cdot\nvec_{f,\ell}\bigr)_{L^2(F_f)}.
\end{equation}
The first term can be bounded as follows:
\begin{equation}\label{eq:cps-constant-1}
\bigl(\grad v_h,\svec\bigr)_{\bL^2(K_\ell)}\leq\norm{\grad v_h}_{\bL^2(K_\ell)} \norm{\svec}_{\bL^2(K_\ell)}.
\end{equation}
Due to the patch-test, the second term reads:
\begin{equation}\label{eq:cps-constant-2}
\begin{aligned}
\sum_{\ell\in\Ical_K}\sum_{f\in\Ical_{F,\ell}}\bigl(v_h,\svec\cdot\nvec_{f,\ell}\bigr)_{L^2(F_f)} &amp;=
\sum_{\ell\in\Ical_K}\sum_{f\in\Ical_{F,\ell}}\bigl(v_h-\ul{v}_{h,f},(\svec-\ul{\svec}_f)\cdot\nvec_{f,\ell}\bigr)_{L^2(F_f)}\\
&amp;\leq\sum_{\ell\in\Ical_K}\sum_{f\in\Ical_{F,\ell}}\bigl\lVert v_h-\ul{v}_{h,f}\bigr\rVert_{L^2(F_f)} \bigl\lVert(\svec-\ul{\svec}_f)\cdot\nvec_{f,\ell}\bigr\rVert_{L^2(F_f)}.
\end{aligned}
\end{equation}
Using inequality~\eqref{eq:MTI}, we have:
\begin{equation}\label{eq:cps-constant-3-rg}
\bigl\lVert v_h-\ul{v}_{h,f}\bigr\rVert_{L^2(F_f)} \bigl\lVert(\svec-\ul{\svec}_f)\cdot\nvec_{f,\ell}\bigr\rVert_{L^2(F_f)}\lesssim \sigma_\ell h_{\ell} \norm{\grad v_h}_{\bL^2(K_\ell)} \norm{\Grad\svec}_{\bbL^2(K_\ell)}.
\end{equation}
Using~\eqref{eq:cps-constant-3-rg} in~\eqref{eq:cps-constant-2}, combining the result with~\eqref{eq:cps-constant-1}, inequality~\eqref{eq:cps-constant-0} now reads:
\[
\begin{aligned}
\norm{v_h}_{L^2(\Om)}^2&amp;\lesssim \sum_{\ell\in\Ical_K}\norm{\grad v_h}_{\bL^2(K_\ell)} \Bigl(\norm{\svec}_{\bL^2(K_\ell)}+\sigma_\ell h_{\ell} \norm{\Grad\svec}_{\bbL^2(K_\ell)}\Bigr)\\
&amp;\lesssim \sigma\sum_{\ell\in\Ical_K}\norm{\grad v_h}_{\bL^2(K_\ell)} \Bigl(\norm{\svec}_{\bL^2(K_\ell)}+h_\Om \norm{\Grad\svec}_{\bbL^2(K_\ell)}\Bigr).
\end{aligned}
\]
We obtain~\eqref{eq:cps-constant} using the discrete Cauchy--Schwarz inequality and~\eqref{eq:svecvh}.
\end{proof}

As a consequence of Proposition~\ref{pro:cps-constant}, we have the following result.

\begin{prop}\label{pro:broknorm}
The broken norm $v_h\to\norm{v_h}_h$ is a norm over $X_{0,h}$.
\end{prop}

The space of nonconforming approximation of $\bH^1(\Om)$ (resp.\ $\bH^1_0(\Om)$) of order $k$ is $\bX_h=(X_h)^d$ (resp.\ $\bX_{0,h}=(X_{0,h})^d$). We set $\Xcal_h\coloneqq \bX_{0,h}\times Q_h$ where $Q_h=P^{k-1}_\disc(\Tcal_h)\cap L^2_{\zmv}(\Om)$. We deduce from Proposition~\ref{pro:broknorm} the following result.
\begin{prop}
The broken norm defined below is a norm on $\Xcal_h$:
\begin{equation*}
\bigl\lVert({\,\cdot\,,\cdot\,})\bigr\rVert_{\Xcal_h}\colon
\left\{
\begin{aligned}
\Xcal_h			&amp;\to\R\\
(\vvec_h,q_h)	&amp;\mapsto \left(\norm{\vvec_h}_h^2+\nu^{-2}\norm{q_h}_{L^2(\Om)}^2\right)^{\frac12}
\end{aligned}
\right.
\end{equation*}
\end{prop}

Thus, the product space $\Xcal_h$ endowed with the broken norm $\norm{({\,\cdot\,,\cdot\,})}_{\Xcal_h}$ is a Hilbert space. We consider the discrete continuous bilinear form $a_{S,h}({\,\cdot\,,\cdot\,})$ such that:
\[
a_{S,h}\colon
\left\{
\begin{aligned}
\Xcal_h\times\Xcal_h					&amp;\to\R\\
(\uvec'_h,p'_h)\times(\vvec_h,q_h)	&amp;\mapsto\nu(\uvec'_h,\vvec_h)_h-\bigl(\div_h\vvec_h,p'_h\bigr)_{L^2(\Om)}-\bigl(\div_h\uvec'_h,q_h\bigr)_{L^2(\Om)}
\end{aligned}
\right.
\]
Let us set $\bV_h$ the discrete space of discrete divergence-free velocities :
\begin{equation}\label{eq:Vh}
\bV_h\coloneqq \bigl\{\vvec_h\in\bX_{0,h}\btq \forall q_h\in Q_h, \bigl(\div_h\vvec_h,q_h\bigr)_{L^2(\Om)}=0\bigr\}.
\end{equation}
We recall that the velocities in $\bV_h$ are piecewise divergence-free:
\begin{prop}[{\cite[Lemma~3.1]{BuHa05}}]\label{pro:divvh}
For all\/ $\vvec_h\in\bV_h$, for all $\ell\in\Ical_K$, $\restr{\div\vvec_h}{K_\ell}=0$.
\end{prop}
\begin{proof}
Let $\vvec_h\in\bV_h$. Integrating by parts and using the patch-test, we have:
\[
(\div_h\vvec_h,1)_{L^2(\Om)}= \sum_{\ell\in\Ical_K}\int_{K_\ell}\div\vvec_h=\sum_{f\in\Ical_F}\int_{F_f}[\vvec_h]_{F_f}\cdot\nvec_f=0.
\]
Let $q_h\in P^{k-1}_\disc(\Tcal_h)$ and $\ul{q}_h=\int_\Om q_h/\abs{\Om}$. Then $q_h-\ul{q}_h\in Q_h$ so that: $(\div_h \vvec_h,q_h-\ul{q}_h)_{L^2(\Om)}=0$. Hence, we have: $(\div_h \vvec_h,q_h)_{L^2(\Om)}=(\div_h\vvec_h,\ul{q}_h)_{L^2(\Om)}=0$. Let $\ell\in\Ical_K$. Let $q_h\in P_\disc^{k-1}(\Tcal_h)$ such that $\restr{q_h}{K_\ell}=\restr{\div\vvec_h}{K_\ell}$ and for all $\ell'\in\Ical_K$, $\ell'\neq\ell$, $\restr{q_h}{K_{\ell'}}=0$. We have: $(\div_h \vvec_h,q_h)_{L^2(\Om)}=0$ and $(\div_h \vvec_h,q_h)_{L^2(\Om)}=(\div\vvec_h,q_h)_{L^2(K_\ell)}=\norm{\div\vvec_h}_{L^2(K_\ell)}^2$. Hence $\norm{\div\vvec_h}_{L^2(K_\ell)}^2=0$.
\end{proof}

Let $\Jcal_h\colon X_{0,h}\to Y_{0,h}$, with $Y_{0,h}=\bigl\{v_h\in H^1_0(\Om) \btq \forall\ell\in\Ical_K, \restr{v_h}{K_\ell}\in P^k(K_\ell)\bigr\}$ be the averaging operator described in~\cite[\S 22.4.1]{ErGu21-I}. There exists a constant $C^{\nonconf}_{\Jcal_h}\approx\sigma$ and independent of $h$ such that:
\begin{equation}\label{eq:c-constant}
\forall v_h\in X_{0,h},\;\norm{\Jcal_h v_h}_{H^1_0(\Om)}\leq C^{\nonconf}_{\Jcal_h} \norm{v_h}_h.
\end{equation}
Let $\ell_\fvec\in\Lcal(\Xcal_h,\R)$ be such that for all $(\vvec_h,q_h)\in\Xcal_h$:
\[
\ell_\fvec \bigl((\vvec_h,q_h)\bigr)=
\begin{cases}
(\fvec,\vvec_h)_{\bL^2(\Om)}	&amp;\text{ if }\;\fvec\in\bL^2(\Om),\\
\langle\fvec,\Jcal_h(\vvec_h)\rangle_{\bH^1_0(\Om)} &amp;\text{ if }\;\fvec\in\bH^{-1}(\Om).
\end{cases}
\]
The nonconforming discretization of Problem~\eqref{eq:Stokes-FVa} reads:
\begin{equation}\label{eq:Stokes-FVah}
\text{Find $(\uvec_h,p_h)\in\Xcal_h$ such that $a_{S,h} \bigl((\uvec_h,p_h),(\vvec_h,q_h)\bigr)=\ell_\fvec \bigl((\vvec_h,q_h)\bigr)\;\forall(\vvec_h,q_h)\in\Xcal_h$.}
\end{equation}
To prove that Problem~\eqref{eq:Stokes-FVah} is well-posed, we will also use the T-coercivity theory. We do not need the well-posedness of the continuous problem, i.e.\ Proposition~\ref{pro:aSTcoer}, but we will follow its proof, using a Fortin operator. This operator will be explained later, using the discrete basis functions. We will see that the discrete stability constant depends on this operator (hence polynomial order~$k$).

\begin{prop}\label{pro:aSTcoerh}
Suppose that there exists a Fortin operator $\Pi_{\nonconf}\colon\bH^1(\Om)\to\bX_h$ such that
\begin{align}\label{eq:WellPosed-nc-1}
\exists C_{\nonconf}\;\text{ such that }\; \forall\vvec\in\bH^1(\Om),\;\norm{\Pi_{\nonconf}\vvec}_h&amp;\leq C_{\nonconf}\norm{\Grad\vvec}_{\bbL^2(\Om)},\\\label{eq:WellPosed-nc-2}
\forall\vvec\in\bH^1(\Om),\;(\div_h\Pi_{\nonconf}\vvec,q_h)_{L^2(\Om)}&amp;=(\div\vvec,q_h)_{L^2(\Om)},\;\forall q\in Q_h,
\end{align}
where the constant $C_{\nonconf}$ does not depend on $h$. Then, the bilinear form $a_{S,h}({\,\cdot\,,\cdot\,})$ is $T$-coercive:
\begin{equation}\label{eq:aSTcoerh}
\begin{gathered}
\exists T_h\in\Lcal(\Xcal_h),\;\text{ bijective, }\;\exists\alpha_{T_h}&gt;0, \forall(\uvec_h',p_h')\in\Xcal_h,\\
a_{S,h}\Bigl((\uvec_h',p_h'),T_h\bigl((\uvec_h',p_h')\bigr)\Bigr)\geq \alpha_{T_h} \norm*{(\uvec_h',p_h')}_{\Xcal_h}^2.
\end{gathered}
\end{equation}
\end{prop}
We will exhibit $C_{\nonconf}$ for $k=1$, $d=2, 3$ in Section~\ref{sec:CRNC}, obtaining $C_{\nonconf}=1$; and then for $k=2$, $d=2$ in Section~\ref{sec:FoSo}, obtaining $C_{\nonconf}=\sigma^2+1$. Let us prove Proposition~\ref{pro:aSTcoerh}.
\begin{proof}
We follow the proof of Proposition~\ref{pro:aSTcoer}. Let us consider $(\uvec'_h,p'_h)\in\Xcal_h$ and let us build $(\vvec_h^\star,q_h^\star)\in\Xcal_h$ satisfying~\eqref{eq:Tcoer} (with $V=\Xcal_h$). We need three main steps.
\begin{stepB}\label{fifi}
According to Proposition~\ref{prop:diviso}, there exists $\tilde{\vvec}_{p'_h}\in\bV^{\perp}$ such that $\div\tilde{\vvec}_{p'_h}=p'_h$ in $\Om$ and $\norm{\tilde{\vvec}_{p'_h}}_{\bH^1_0(\Om)}\leq C_{\div}\norm{p'_h}_{L^2(\Om)}$. Let us set $\vvec_{p'_h}=\nu^{-1}\tilde{\vvec}_{p'_h}$. Consider $\vvec_{h,p'_h}=\Pi_{\nonconf}\vvec_{p'_h}$, for all $q_h\in Q_h$ we have:
\begin{equation}\label{eq:GiRa-nc}
\bigl(\div_h\vvec_{h,p'_h},q_h\bigr)_{L^2(\Om)}=\nu^{-1} (p'_h,q_h)_{L^2(\Om)} \quad\text{ and }\quad \norm{\vvec_{h,p'_h}}_h\leq\nu^{-1} C_{\div}^{\nonconf} \norm{p'_h}_{L^2(\Om)}
\end{equation}
where $C_{\div}^{\nonconf}=C_{\nonconf} C_{\div}$. Let us set $(\vvec_h^\star,q_h^\star)\coloneqq \bigl(\gamma_{\nonconf}\uvec'_h-\vvec_{h,p'_h},-\gamma_{\nonconf} p'_h\bigr)$, with $\gamma_{\nonconf}&gt;0$. We obtain:
\begin{equation}\label{eq:coer1-nc}
a_{S,h} \bigl((\uvec'_h,p'_h), (\vvec_h^\star,q_h^\star)\bigr)=
\nu \gamma_{\nonconf}\norm{\uvec'_h}_h^2+\nu^{-1}\norm{p'_h}_{L^2(\Om)}^2-\nu(\uvec'_h,\vvec_{h,p'_h})_h.
\end{equation}
\end{stepB}
\begin{stepB}
In order to bound the last term of~\eqref{eq:coer1-nc}, we use Young inequality and then inequality~\eqref{eq:GiRa-nc} so that for all $\eta_{\nonconf}&gt;0$:
\begin{equation}\label{eq:GvGvp-nc}
\bigl(\uvec'_h,\vvec_{h,p'_h}\bigr)_h\leq
\frac{\eta_{\nonconf}}{2}\norm{\uvec'_h}_h^2+\frac{\eta_{\nonconf}^{-1}}{2} \Biggl(\frac{C_{\div}^{\nonconf}}{\nu}\Biggr)^2 \norm{p'_h}_{L^2(\Om)}^2.
\end{equation}
\end{stepB}
\begin{stepB}
Using the bound~\eqref{eq:GvGvp-nc} in~\eqref{eq:coer1-nc} and choosing $\eta_{\nonconf}=\gamma_{\nonconf}$, we get:
\[
a_{S,h} \bigl((\uvec'_h,p'_h), (\vvec_h^\star,q_h^\star)\bigr)\geq\nu\Biggl(\frac{\gamma_{\nonconf}}{2} \norm{\uvec'_h}_h^2+\nu^{-2}\Biggl(1-\frac{(\gamma_{\nonconf})^{-1}}{2} (C_{\div}^{\nonconf})^2 \Biggr)\norm{p'_h}_{L^2(\Om)}^2\Biggr).
\]
Consider now $\gamma_{\nonconf}=(C_{\div}^{\nonconf})^2$. We obtain:
\[
a_{S,h} \bigl((\uvec'_h,p'_h), (\vvec_h^\star,q_h^\star)\bigr)\geq \alpha_TC_{\min}^{\nonconf} \norm*{(\uvec'_h,p'_h)}_{\Xcal_h}^2,\;\text{ where }\;C_{\min}^{\nonconf}=\min\bigl((C_{\div}^{\nonconf})^2,1\bigr).
\]
We obtain~\eqref{eq:aSTcoerh} with $\alpha_{T_h}=\alpha_T C_{\min}^{\nonconf}$. Suppose that $C_{\nonconf}\geq 1$. Then $\alpha_{T_h}=\alpha_T=\frac{\nu}{2}$.
\end{stepB}
The operator $T_h$ such that $T_h(\uvec_h',p'_h)=(\vvec^\star_h,p^\star_h)$ is linear and continuous. We have indeed:
\[
\norm*{T_h(\uvec_h',p'_h)}_{\Xcal_h}^2=\norm{\vvec_h^\star}_h^2+\nu^{-2} \norm{q_h^\star}_{L^2(\Omega)}^2\leq (C_{\max}^{\nonconf})^2 \norm*{(\uvec_h',p_h')}_{\Xcal_h}^2
\]
where $C_{\max}^{\nonconf}=C_{\div}^{\nonconf}\bigl(\max\bigl(2+(C_{\div}^{\nonconf})^2,2(C_{\div}^{\nonconf})^2\bigr)\bigr)^{\frac12}$.

Remark that, given $(\vvec_h^\star,q_h^\star)\in\Xcal_h$, choosing $(\uvec_h',p_h')=\bigl(\gamma_{\nonconf}^{-1}\vvec_h^\star-\gamma_{\nonconf}^{-2}\vvec_{h,q_h^\star},-\gamma_{\nonconf}^{-1}q_h^\star\bigr)$ yields $T_h\bigl((\uvec_h',p_h')\bigr)=(\vvec_h^\star,q_h^\star)$. Hence, the operator $T_h\in\Lcal(\Xcal_h)$ is bijective.
\end{proof}

\begin{rema}
We recover the first Banach--Ne\v{c}as--Babu\v{s}ka condition~\cite[Theorem~25.9, (BNB1)]{ErGu21-II}:
\[
a_{S,h} \bigl((\uvec_h',p_h'), (\vvec_h^\star,q_h^\star)\bigr)\geq \alpha_{T_h} (C_{\max}^{\nonconf})^{-1} \norm*{(\uvec_h',p_h')}_{\Xcal_h} \norm*{(\vvec_h^\star,q_h^\star)}_{\Xcal_h}.
\]
\end{rema}

We can now prove the discrete counterpart of Theorem~\ref{thm:StokesWellPosed}.
\begin{theo}\label{thm:WellPosed-nc}
Suppose that there exists a Fortin operator $\Pi_{\nonconf}\colon\bH^1(\Om)\to\bX_h$ satifying~\eqref{eq:WellPosed-nc-1}-\eqref{eq:WellPosed-nc-2}. Then Problem~\eqref{eq:Stokes-FVah} is well-posed. It admits one and only one solution $(\uvec_h,p_h)$ such that:
\begin{equation}\label{eq:stabiNC}
\begin{gathered}
\text{if }\;\fvec\in\bL^2(\Om):
\norm{\uvec_h}_h\lesssim \nu^{-1} C_0^{\nonconf} \norm{\fvec}_{\bL^2(\Om)} ,\;\norm{p_h}_{L^2(\Om)}\lesssim 2 C_0^{\nonconf} C_{\div}^{\nonconf} \norm{\fvec}_{\bL^2(\Om)},\\
\text{if }\;\fvec\in\bH^{-1}(\Om):
\norm{\uvec_h}_h\lesssim \nu^{-1} C_{\Jcal_h}^{\nonconf} \norm{\fvec}_{\bH^{-1}(\Om)},\;
\norm{p_h}_{L^2(\Om)}\lesssim 2 C_{\Jcal_h}^{\nonconf} C_{\div}^{\nonconf} \norm{\fvec}_{\bH^{-1}(\Om)},
\end{gathered}
\end{equation}
where $C_0^{\nonconf}=\sigma C_\Om h_\Om$.
\end{theo}
\begin{proof}
Consider $(\uvec_h,p_h)$ the unique solution of Problem~\eqref{eq:Stokes-FVah}. Choosing $\vvec_h=0$, we obtain that $\div_h\uvec_h=0$. Let $\fvec\in\bL^2(\Om)$. Now, choosing $\vvec_h=\uvec_h$ in~\eqref{eq:Stokes-FVah}, using Cauchy--Schwarz inequality, we get that: $\norm{\uvec_h}_h\leq \nu^{-1} \sigma C_\Om h_\Om \norm{\fvec}_{\bL^2(\Om)}$ using inequality~\eqref{eq:cps-constant}. Consider $(\vvec_h,q_h)=(\vvec_{h,p_h},0)$ in~\eqref{eq:Stokes-FVah}, where $\vvec_{h,p_h}=\Pi_{\nonconf}\vvec_{p_h}$ is built as $\vvec_{h,p'_h}$ in Step~\ref{fifi} of the proof of Proposition~\ref{pro:aSTcoerh}, setting $p'_h=p_h$. Suppose that $\fvec\in\bL^2(\Om)$. Notice that $\nu^{-1} \norm{p_h}_{L^2(\Om)}^2=\nu (\uvec_h,\vvec_{h,p_h})_h-(\fvec,\vvec_{h,p_h})_{\bL^2(\Om)}$. Using Cauchy--Schwarz inequality, we have: $\nu^{-1} \norm{p_h}_{L^2(\Om)}^2\leq\nu \norm{\uvec_h}_h \norm{\vvec_{h,p_h}}_h+\norm{\fvec}_{\bL^2(\Om)} \norm{\vvec_{h,p_h}}_{\bL^2(\Om)}$. Using PoincarÃ©--Steklov inequality~\eqref{eq:cps-constant}, hypothesis~\eqref{eq:WellPosed-nc-1}, and the previous estimate on $\norm{\uvec_h}_h$, we have:
\[
\norm{p_h}_{L^2(\Om)}^2\lesssim 2 \sigma C_\Om h_\Om \norm{\fvec}_{\bL^2(\Om)} \norm{\vvec_{h,p_h}}_h\lesssim 2 \sigma C_\Om h_\Om C_{\div}^{\nonconf} \norm{\fvec}_{\bL^2(\Om)} \norm{p_h}_{L^2(\Om)}.
\]
Let $\fvec\in\bH^{-1}(\Om)$. We apply the same reasoning, using inequality~\eqref{eq:c-constant}.
\end{proof}

As a corollary of Theorem~\ref{thm:WellPosed-nc}, the following a priori error estimates follow:
\begin{coro}[{\cite[Theorems~3, 4, 6]{CrRa73}, \cite[Equation~(47)]{FoSo83}}]\label{cor:WellPosed-nc}
Under the assumption of Theorem~\ref{thm:WellPosed-nc}, suppose that $(\uvec,p)\in\bigl(\bH^{1+k}(\Om)\cap\bH^1_0(\Om)\bigr)\times\bigl(H^k(\Om)\cap L^2_{\zmv}(\Om)\bigr)$, we have the estimates:
\begin{align}\label{eq:aprioriEEuh1}
\norm{\uvec-\uvec_h}_h\lesssim\sigma^l h^k 
\bigl(\abs{\uvec}_{\bH^{k+1}(\Om)}+\nu^{-1} \abs{p}_{H^k(\Om)}\bigr),\\\label{eq:aprioriEEph}
\nu^{-1} \norm{p-p_h}_{L^2(\Om)}\lesssim\sigma^l h^k 
\bigl(\abs{\uvec}_{\bH^{k+1}(\Om)}+\nu^{-1} \abs{p}_{H^k(\Om)}\bigr).
\end{align}
Suppose moreover that the domain $\Om$ is convex. Then we have:
\begin{equation}\label{eq:aprioriEEuh0}
\norm{\uvec-\uvec_h}_{\bL^2(\Om)}\lesssim\sigma^{2 l} h^{k+1} 
\bigl(\abs{\uvec}_{\bH^{k+1}(\Om)}+\nu^{-1} \abs{p}_{H^k(\Om)}\bigr).
\end{equation}
The hidden constants depend on $k$ but they don't depend on the mesh. The parameter $\sigma$ is the shape regularity parameter and the exponent $l\in\N^\star$ depends on $k$. When $k=1$, $d=2, 3$, we have $l=1$, and when $k=2$, $d=2$, we have $l=2$.
\end{coro}

When $\Om$ is not convex, the exponent on $h$ in Equation~\eqref{eq:aprioriEEuh0} is equal to $k+s$ where $s\in]0,1[$ depends on $\Om$ (cf.~\cite[Theorem~31.33]{ErGu21-II}).

The main issue with nonconforming mixed finite elements is the construction of the basis functions. In a recent paper, Sauter explains such a construction in two dimensions~\cite[Theorem~1.3]{Saut23}, and gives a bound to the discrete counterpart $\beta_\Tcal(\Om)$ of $\beta(\Om)$ defined in~\eqref{eq:LBB-CIS}:
\begin{equation}\label{eq:LBB-CISh}
\beta_\Tcal(\Om)=\inf_{q_h\in Q_h\backslash \{0\}}\sup_{\vvec_h\in\bX_{0,h}\backslash \{0\}}\frac{(\div_h\vvec_h,q_h)_{L^2(\Om)}}{\norm{q_h}_{L^2(\Om)} \norm{\vvec_h}_h}\geq c_\Tcal \bigl(\log(k+1)\bigr)^{-\alpha},
\end{equation}
where the parameter $\alpha$ is explicit and depends on $k$ and on the mesh topology; and the constant $c_\Tcal$ depends only on the shape-regularity of the mesh.

\section{Examples of nonconforming discretization for Stokes problem}\label{sec:NC-disc}

\subsection{Nonconforming Crouzeix--Raviart mixed finite elements for \texorpdfstring{$k=1$}{k=1}}\label{sec:CRNC}

We study the lowest order nonconforming Crouzeix--Raviart mixed finite elements~\cite{CrRa73}. Let us consider $X_{\CR}$ (resp.\ $X_{0,\CR}$), the space of nonconforming approximation of $H^1(\Om)$ (resp.\ $ H^1_0(\Om)$) of order $1$:
\begin{equation}\label{eq:XCR}
\begin{gathered}
X_{\CR}=\Biggl\{v_h\in P^1_\disc(\Tcal_h) \,\Bigg\vert\,
\forall f\in\Ical_F^i, \; \int_{F_f}[v_h]_{F_f}=0\Biggr\},\\
X_{0,\CR}=\Biggl\{v_h\in X_{\CR} \,\Bigg\vert\, \forall f\in\Ical_F^b, \; \int_{F_f} v_h=0\Biggr\}.
\end{gathered}
\end{equation}
The space of nonconforming approximation of $\bH^1(\Om)$ (resp.\ $\bH^1_0(\Om)$) of order $1$ is $\bX_{\CR}=(X_{\CR})^d$ (resp.\ $\bX_{0,\CR}=(X_{0,\CR})^d$). We set $\Xcal_{\CR}\coloneqq \bX_{0,\CR}\times Q_{\CR}$ where $Q_{\CR}=P^0_\disc(\Tcal_h)\cap L^2_{\zmv}(\Om)$.

We can endow $X_{\CR}$ with the basis $(\psi_f)_{f\in\Ical_F}$ such that:
\[
\forall \ell\in\Ical_K, \; \restr{\psi_f}{K_\ell}=
\begin{cases}
1-d\lambda_{i,\ell}	&amp;\text{ if }\;f\in\Ical_{F,\ell},\\
0					&amp;\text{ otherwise},
\end{cases}
\]
where $i\in\mathcal{I}_S$ is such that $S_i$ is the vertex opposite to $F_f$ in $K_\ell$. We then have $\restr{\psi_f}{F_f}=1$, so that $[\psi_f]_{F_f}=0$ if $f\in\Ical_F^i$, and for all $f, f'\in\Ical_F$, $f'\neq f$, $\int_{F_{f'}}\psi_f=0$.

We have: $X_{\CR}=\vect\bigl((\psi_f)_{f\in\Ical_F}\bigr)$ and $X_{0,\CR}=\vect\bigl((\psi_f)_{f\in\Ical_F^i}\bigr)$.

The Crouzeix--Raviart interpolation operator $\pi_{\CR}$ for scalar functions is defined by:
\[
\pi_{\CR}\colon
\left\{
\begin{aligned}
H^1(\Om)	&amp;\to X_{\CR}\\
v		&amp;\mapsto\sum_{f\in\Ical_F}\pi_fv \psi_{f}
\end{aligned}
\right.
\qquad\text{ where }\;\pi_fv=\frac{1}{\abs{F_f}}\int_{F_f} v.
\]
Notice that $\forall f\in\Ical_F$, $\int_{F_f} \pi_{\CR}v=\int_{F_f} v$. Moreover, the Crouzeix--Raviart interpolation operator preserves the constants, so that $\pi_{\CR}\ul{v}_\Om=\ul{v}_\Om$ where $\ul{v}_\Om=\int_\Om v/\abs{\Om}$. We recall that for $k=1$, the coefficient $C_{\nonconf}$ in~\eqref{eq:WellPosed-nc-1} is equal to $1$:
\begin{lemm}[{\cite[Lemma~2]{ApNS01}}]\label{lem:piCR}
The Crouzeix--Raviart interpolation operator $\pi_{\CR}$ is such that:
\begin{equation}\label{eq:gradvhleqgradv-glo}
\forall v\in H^1(\Om),\;
\norm{\pi_{\CR}v}_h\leq \norm{\grad v}_{\bL^2(\Om)}.
\end{equation}
\end{lemm}
\begin{proof}
We have, integrating by parts twice and using Cauchy--Schwarz inequality:
\[
\begin{aligned}
\grad\pi_{\CR}\restr{v}{K_\ell}&amp;=
\abs{K_\ell}^{-1}\int_{K_\ell}\grad\pi_{\CR}v\\
&amp;=\abs{K_\ell}^{-1}\sum_{f\in\Ical_{F,\ell}}\int_{F_f}\pi_{\CR}v \nvec_f\\
&amp;=\abs{K_\ell}^{-1}\sum_{f\in\Ical_{F,\ell}}\int_{F_f}v \nvec_f\\
&amp;=\abs{K_\ell}^{-1}\int_{K_\ell}\grad v,\\
\bigl\lvert\grad\pi_{\CR}\restr{v}{K_\ell}\bigr\rvert &amp;\leq\abs{K_\ell}^{-\frac12} \norm*{\grad v}_{\bL^2(K_\ell)}\\
\Rightarrow\quad\norm*{\grad \pi_{\CR}v}_{\bL^2(K_\ell)}^2&amp;\leq\norm*{\grad v}_{\bL^2(K_\ell)}^2.
\end{aligned}
\]
Summing these local estimates over $\ell\in\Ical_K$, we obtain~\eqref{eq:gradvhleqgradv-glo}.
\end{proof}

For a vector $\vvec\in\bH^1(\Om)$ of components $(v_{d'})_{d'=1}^d$, the Crouzeix--Raviart interpolation operator is such that: $\Pi_{\CR}\vvec=\bigl(\pi_{\CR}v_{d'}\bigr)_{d'=1}^d$. Let us set $\Pi_f\vvec=\bigl(\pi_fv_{d'}\bigr)_{d'=1}^d$.

\begin{lemm}\label{lem:PiCR}
The Crouzeix--Raviart interpolation operator $\Pi_{\CR}$ can play the role of the Fortin operator:
\begin{align}\label{eq:WellPosed-CR-1}
\forall\vvec\in\bH^1(\Om),\;\norm{\Pi_{\CR}\vvec}_h&amp;\leq\norm{\Grad\vvec}_{\bbL^2(\Om)},\\\label{eq:WellPosed-CR-2}
\forall\vvec\in\bH^1(\Om),\;\bigl(\div_h\Pi_{\CR}\vvec,q_h\bigr)_{L^2(\Om)}&amp;=\bigl(\div\vvec,q_h\bigr)_{L^2(\Om)},\;\forall q\in Q_h.
\end{align}
Moreover, for all $\vvec\in\bP^1(\Om)$, $\Pi_{\CR}\vvec=\vvec$.
\end{lemm}
\begin{proof}
We obtain~\eqref{eq:WellPosed-CR-1} applying Lemma~\ref{lem:piCR} component by component. By integrating by parts, we have $\forall\vvec\in\bH^1(\Om)$, $\forall\ell\in\Ical_K$:
\[
\int_{K_\ell}\div\Pi_{\CR}\vvec =\sum_{f\in\Ical_{F,\ell}}\int_{F_f}\Pi_{\CR}\vvec\cdot\nvec_f =\sum_{f\in\Ical_{F,\ell}}\int_{F_f}\Pi_f\vvec\cdot\nvec_f=\sum_{f\in\Ical_{F,\ell}}\int_{F_f} \vvec\cdot\nvec_f =\int_{K_\ell}\div\vvec,
\]
so that~\eqref{eq:WellPosed-CR-2} is satisfied.
\end{proof}

We can apply the T-coercivity theory to show the following result:
\begin{theo}
Let $\Xcal_h=\Xcal_{\CR}$. Then the continuous bilinear form $a_{S,h}({\,\cdot\,,\cdot\,})$ is $T_h$-coercive and Problem~\eqref{eq:Stokes-FVah} is well-posed.
\end{theo}

\begin{proof}
Using estimates~\eqref{eq:WellPosed-CR-1} and~\eqref{eq:cps-constant}, we apply the proof of Theorem~\ref{thm:WellPosed-nc}.
\end{proof}

Since the constant of the interpolation operator $\Pi_{\CR}$ is equal to $1$, we have $C_{\min}^{\nonconf}=\min\bigl((C_{\div})^2,1\bigr)=1$ and $C_{\max}^{\nonconf}=C_{\max}$: the stability constant of the nonconforming Crouzeix--Raviart mixed finite elements is independent of the mesh. This is not the case for higher order (see~\cite[Theorem~2.2]{CaSa22}).

\subsection{Comments on higher-order methods}

For higher order, we cannot build the interpolation operator component by component, since higher-order divergence moments must be preserved. Thus, for $k&gt;1$, we must build $\Pi_{\nonconf}$ so that for all $\vvec\in\bH^1(\Om)$, for all $\ell\in\Ical_K$, for all $q\in P^{k-1}(K_\ell)$:
\[
\int_{K_\ell}q \div\Pi_{\nonconf}\vvec=\int_{K_\ell}q \div\vvec.
\]
We recall that by integration by parts, we have:
\begin{equation}\label{eq:Green1stId}
\int_{K_\ell}q \div\Pi_{\nonconf}\vvec+\int_{K_\ell}\grad q \cdot\Pi_{\nonconf}\vvec=\int_{\pa K_\ell}q \Pi_{\nonconf}\vvec\cdot\restr{\nvec}{\pa K_\ell}.
\end{equation}
Hence, to obtain a local estimate of $\bigl\lVert\Grad\Pi_{\nonconf}\vvec\bigr\rVert_{\bbL^2(K_\ell)}$, we will need the following Lemma:

\begin{lemm}\label{lem:IPPnormcomp}
Let\/ $\vvec\in \bH^1(K_\ell)$ and $q\in P^{k-1}(K_\ell)\cap L^2_{\zmv}(K_\ell)$. We have:
\begin{equation}\label{eq:IPPnormcomp}
\Biggl\lvert\int_{\pa K_\ell}q (\vvec-\ul{\vvec}_\ell)\cdot\restr{\nvec}{\pa K_\ell}\Biggr\rvert\leq \bigl(\sqrt{d}+1\bigr) \pi^{-1}  h_\ell \norm*{\grad q}_{\bL^2(K_\ell)} \bigl\lVert\Grad\vvec\bigr\rVert_{\bbL^2(K_\ell)}
\end{equation}
\end{lemm}
\begin{proof}
We have by integration by parts, and then using Cauchy--Schwarz inequality:
\begin{align*}
\Biggl\lvert\int_{\pa K_\ell}q (\vvec-\ul{\vvec}_\ell)\cdot\restr{\nvec}{\pa K_\ell}\Biggr\rvert
&amp;\leq\Biggl\lvert\int_{K_\ell}q \div(\vvec-\ul{\vvec}_\ell)\Biggr\rvert+\Biggl\lvert\int_{K_\ell}\grad q \cdot(\vvec-\ul{\vvec}_\ell)\Biggr\rvert\\
&amp;\leq \sqrt{d} \norm{q}_{L^2(K_\ell)} \norm*{\Grad(\vvec-\ul{\vvec}_\ell)
}_{\bbL^2(K_\ell)}+\norm*{\grad q}_{\bL^2(K_\ell)} \norm*{(\vvec-\ul{\vvec}_\ell)}_{\bL^2(K_\ell)}\\
&amp;\leq \bigl(\sqrt{d}+1\bigr) \pi^{-1} h_\ell \norm*{\grad q}_{\bL^2(K_\ell)} \bigl\lVert\Grad\vvec\bigr\rVert_{\bbL^2(K_\ell)},\;\text{ using~\eqref{eq:PSinC} twice.}\qedhere
\end{align*}
\end{proof}

In the next section, we will see that for $k=2$, $d=2$, we will need Lemma~\ref{lem:IPPnormcomp}. For $k\geq 3$, it could be necessary to bound the tangential components of $\vvec-\ul{\vvec}_\ell$. To do so, we would need to preserve curl integrals on $K_\ell$. Indeed, by integration by parts, we have:
\[
\begin{aligned}
\text{For }\;d=2,\; \vvec\in\bH^1(\Om),\; q\in P^{k-1}(K_\ell):&amp;
\int_{K_\ell}q \bigl(\curlv q\cdot\vvec-\curls\vvec q\bigr)=\int_{\pa K_\ell}q \vvec\times\restr{\nvec}{\pa K_\ell}.\\
\text{For }d=3, \vvec\in\bH^1(\Om), \wvec\in \bP^{k-1}(K_\ell):&amp;
\int_{K_\ell}\bigl(\wvec\cdot\curlv\vvec-\curlv \wvec\cdot\vvec\bigr)\!=\!\int_{\pa K_\ell}\bigl(\restr{\nvec}{\pa K_\ell}\times\vvec\times\restr{\nvec}{\pa K_\ell}\bigr)\cdot\bigl(\wvec\times\restr{\nvec}{\pa K_\ell}\bigr).
\end{aligned}
\]

\subsection{Fortin--Soulie mixed finite elements}\label{sec:FoSo}

We consider here the case $k=2$, $d=2$ and we study the so-called Fortin--Soulie mixed finite elements~\cite{FoSo83}. We consider a shape-regular triangulation sequence $(\Tcal_h)_h$.

Let us consider $X_{\FS}$ (resp.\ $X_{0,\FS}$), the space of nonconforming approximation of $H^1(\Om)$ (resp.\ $H^1_0(\Om)$) of order $2$:
\begin{equation}\label{eq:XFS}
\begin{gathered}
X_{\FS}=\Biggl\{v_h\in P^2_\disc(\Tcal_h) \,\Bigg\vert\,
\forall f\in\Ical_F^i, \; \forall q_h\in P^1(F_f), \; \int_{F_f}[v_h]_{F_f} q_h=0\Biggr\},\\
X_{0,\FS}=\Biggl\{v_h\in X_{\FS} \,\Bigg\vert\, \forall f\in\Ical_F^b, \; \forall q_h\in P^1(F_f), \; \int_{F_f}v_h q_h=0\Biggr\}.
\end{gathered}
\end{equation}
The space of nonconforming approximation of $\bH^1(\Om)$ (resp.\ $\bH^1_0(\Om)$) of order $2$ is $\bX_{\FS}=(X_{\FS})^2$ (resp.\ $\bX_{0,{\FS}}=(X_{0,{\FS}})^2$). We set $\Xcal_{\FS}=\bX_{0,\FS}\times Q_{\FS}$ where $Q_{\FS}\coloneqq P^1_\disc(\Tcal_h)\cap L^2_{\zmv}(\Om)$.

The building of a basis for $X_{0,\FS}$ is more involved than for $X_{0,\CR}$ since we cannot use two points per facet as degrees of freedom. Indeed, for all $\ell\in\Ical_K$, there exists a polynomial of order $2$ vanishing on the Gauss--Legendre points of the facets of the boundary $\pa K_\ell$. Let $f\in\Ical_F$. The barycentric coordinates of the two Gauss--Legendre points $(p_{+,f},p_{-,f})$ on $F_f$ are such that:
\[
p_{+,f}=(c_+,c_-),\;  p_{-,f}=(c_-,c_+),\;\text{ where }\;c_{\pm}=\bigl(1\pm 1/\sqrt{3}\bigr)/2.
\]
These points can be used to integrate exactly order three polynomials:
\[
\forall g\in P^3(F_f),\; \int_{F_f} g=\frac{\abs{F_f}}{2}\bigl(g(p_{+,f})+g(p_{-,f})\bigr).
\]
For all $\ell\in\Ical_K$, we define the quadratic function $\phi_{K_\ell}$ that vanishes on the six Gauss--Legendre points of the facets of $K_\ell$ (see Figure~\ref{fig:GLpoints}):
\begin{equation}\label{eq:phiKl}
\phi_{K_\ell}\coloneqq 2-3\sum_{i\in\Ical_{S,\ell}}\lambda_{i,\ell}^2\;\text{ such that }\; \forall f\in\Ical_{F,\ell}, \forall q\in P^1(F_f),\;\int_{F_f}\phi_{K_\ell} q=0.
\end{equation}
\begin{figure}[h!]
\centering
\begin{tikzpicture}
\draw[name path=line 1] (0,0)--(4,0)--(2,3)--(0,0);
\draw[dashed,name path=line 2] (2,1) ellipse (1.6330 and 1.41);
\fill[red,name intersections={of=line 1 and line 2,total=\t}]
   \foreach \s in {1,...,\t}{(intersection-\s) circle (2pt) node [below right] {$p_{\s}$}};
\end{tikzpicture}
\caption{The six Gauss--Legendre points of an element $K_\ell$ and the elliptic function $\phi_{K_\ell}$.}\label{fig:GLpoints}
\end{figure}

We consider the set of the elliptic functions $\phi_{K_\ell}$:
\begin{equation}\label{eq:PhiH}
\Phi_h\coloneqq \bigl\{\phi_h\in L^2(\Om) \btq \forall\ell\in\Ical_K, \; \restr{\phi_h}{K_\ell}=v_{K_\ell}\phi_{K_\ell}, \; v_{K_\ell}\in\R\bigr\}.
\end{equation}
We also define the spaces of $P^2$-Lagrange functions:
\[
X_{\LG}\coloneqq \bigl\{v_h\in H^1(\Om) \btq \forall \ell\in\Ical_K, \; \restr{v_h}{K_\ell}\in P^2(K_\ell)\bigr\},\quad X_{0,\LG}\coloneqq \bigl\{v_h\in X_{\LG} \btq \restr{v_h}{\pa\Om}=0\bigr\}.
\]
The proposition below allows to build a basis for $X_{0,\FS}$:
\begin{prop}[{\cite[Proposition~1]{FoSo83}}]
We have the following decomposition: $X_{\FS}=X_{\LG}+\Phi_h$ with $\dim(X_{\LG}\cap \Phi_h)=1$. Any function of $X_{\FS}$ can be written as the sum of a function of $X_{\LG}$ and a function of $\Phi_h$. This representation can be made unique by specifying one degree of freedom.
\end{prop}

Notice that $\Phi_h\cap X_{\LG}=\vect{(v_\Phi)}$, where for all $\ell\in\Ical_K$, $\restr{v_\Phi}{K_\ell}=\phi_{K_\ell}$. Then, counting the degrees of freedom, one can show that $\dim(X_{\FS})=\dim(X_{\LG})+\dim(\Phi_h)+1$. For problems involving Dirichlet boundary conditions we can thus prove that for $X_{0,\FS}$ the representation is unique and $X_{0,\FS}=X_{0,\LG}\oplus\Phi_h$. We have $X_{\LG}=\vect\bigl((\phi_{S_i})_{i\in\Ical_S},(\phi_{F_f})_{f\in\Ical_F}\bigr)$ where the basis functions are such that:
\[
\forall i,j\in\Ical_S, \forall f,f'\in\Ical_F:\quad\phi_{S_i}(S_j)=\delta_{i,j}, \phi_{S_i}(M_f)=0, \quad\phi_{M_f}(M_{f'})=\delta_{f,f'}, \phi_{M_f}(S_i)=0.
\]
For all $\ell\in\Ical_K$, we will denote by $\bigl(\phi_{\ell,j}\bigr)_{j=1}^6$ the local nodal basis such that:
\[
\bigl(\phi_{\ell,j}\bigr)_{j=1}^3=\bigl(\restr{\phi_{S_i}}{K_\ell}\bigr)_{i\in\Ical_{S,\ell}}\qquad\text{ and } \qquad \bigl(\phi_{\ell,j}\bigr)_{j=4}^6=\bigl(\phi_{F_f}\bigr\rvert_{K_\ell}\bigr)_{f\in\Ical_{F,\ell}}.
\]
The spaces $X_{\FS}$ and $X_{0,\FS}$ are such that:
\begin{equation}\label{eq:XFSbasis}
\begin{gathered}
X_{\FS} =\vect\bigl((\phi_{S_i})_{i\in\Ical_S},(\phi_{F_f})_{f\in\Ical_F},(\phi_{K_\ell})_{\ell\in\Ical_K}\bigr),\\
X_{0,\FS}=\vect\bigl((\phi_{S_i})_{i\in\Ical_S^i},(\phi_{F_f})_{f\in\Ical_F^i},(\phi_{K_\ell})_{\ell\in\Ical_K}\bigr).
\end{gathered}
\end{equation}
We propose here an alternative definition of the Fortin interpolation operator proposed in~\cite{FoSo83}. Let us first recall the Scott--Zhang interpolation operator~\cite{ScZh90,Ciar13}. For all $i\in\Ical_S$, we choose some $\ell_i\in\Ical_{K,i}$, and we build the $L^2(K_{\ell_i})$-dual basis $\bigl(\tilde{\phi}_{\ell_i,j}\bigr)_{j=1}^6$ of the local nodal basis such that:
\[
\forall j, j'\in\bigl\{1,\dotsc,6\bigr\},\;\int_{K_{\ell_i}}\tilde{\phi}_{\ell_i,j} \phi_{\ell_i,j'}=\delta_{j,j'}.
\]
%Let us define the Fortin--Soulie interpolation operator for scalar functions by:
%\begin{equation}\label{eq:piFS}
%\pi_{\FS}\colon
%\left\{
%\begin{aligned}
%H^1(\Om)&amp;\to X_{\FS}\\
%v&amp;\mapsto \sum_{i\in\Ical_S}v_{S_i}\phi_{S_i}+\sum_{f\in\Ical_F}v_{F_f}\phi_{F_f}+\sum_{\ell\in\Ical_K}v_{K_\ell}\phi_{K_\ell}.
%\end{aligned}
%\right.
%\end{equation}
Let us define the Fortin--Soulie interpolation operator for scalar functions by:
\begin{equation}\label{eq:piFS}
\pi_{\FS}\colon
\left\{
\begin{aligned}
H^1(\Om)	&amp;\to X_{\FS}\\
v		&amp;\mapsto\tilde{\pi}v
+\sum_{\ell\in\Ical_K}v_{K_\ell}\phi_{K_\ell}
\end{aligned}\right.
\quad\text{ with }
\quad\tilde{\pi}v=\sum_{i\in\Ical_S}v_{S_i}\phi_{S_i}+\sum_{f\in\Ical_F}v_{F_f}\phi_{F_f}.
\end{equation}
\begin{itemize}
\item The coefficients $(v_{S_i})_{i\in\Ical_S}$ are fixed so that: $\forall i\in\Ical_S$, $ v_{S_i}=\int_{K_{\ell,i}}v \tilde{\phi}_{\ell_i,j_i}$, where $j_i$ is the index such that $\int_{K_{\ell_i}}\tilde{\phi}_{\ell_i,j_i} \restr{\phi_{S_i}}{K_{\ell_i}}=1$.
\item The coefficients $(v_{F_f})_{f\in\Ical_F}$ are fixed so that: $\forall f\in\Ical_F$, $\int_{F_f}\tilde{\pi}v=\int_{F_f}v$.
\item The coefficients $v_{K_\ell}$ are fixed so that: $\int_{K_\ell}\pi_{\FS}v=\int_{K_\ell}v$.
\end{itemize}
The definition~\eqref{eq:piFS} is more general than the one given in~\cite{FoSo83}, which holds for $v\in H^2(\Om)$.

We set $\vvec_{S_i}\coloneqq \bigl(\tilde{\pi}v_1(S_i),\tilde{\pi}v_2(S_i)\bigr)^T$ and $\vvec_{F_f}\coloneqq \bigl(\tilde{\pi}v_1(F_f),\tilde{\pi}v_2(F_f)\bigr)^T$.

We can define two different Fortin--Soulie interpolation operators for vector functions. First, let
\[
\tilde{\Pi}_{\FS}\colon
\left\{
\begin{aligned}
\bH^1(\Om)&amp;\to\bX_{\FS}\\
\vvec&amp;\mapsto\sum_{i\in\Ical_S}\vvec_{S_i}\phi_{S_i}+\sum_{f\in\Ical_F}\vvec_{F_f}\phi_{F_f}+\sum_{\ell\in\Ical_K}\tilde{\vvec}_{K_\ell}\phi_{K_\ell},
\end{aligned}
\right.
\]
where the coefficients $(\tilde{\vvec}_{K_\ell})_{\ell\in\Ical_K}$ are such that:
\begin{equation}\label{eq:Piave}
\forall\ell\in\Ical_K,\;\int_{K_\ell}\tilde{\Pi}_{\FS}\vvec=\int_{K_\ell}\vvec.
\end{equation}
The interpolation operator $\tilde{\Pi}_{\FS}$ preserves the local averages, but it doesn't preserve the divergence. We then define a second interpolation operator which preserves the divergence in a weak sense:
\[
\Pi_{\FS}\colon
\left\{
\begin{aligned}
\bH^1(\Om)&amp;\to\bX_{\FS}\\
\vvec&amp;\mapsto\sum_{i\in\Ical_S}\vvec_{S_i}\phi_{S_i}+\sum_{f\in\Ical_F}\vvec_{F_f}\phi_{F_f}+\sum_{\ell\in\Ical_K}\vvec_{K_\ell}\phi_{K_\ell}.
\end{aligned}
\right.
\]
For all $\ell\in\Ical_K$, the vector coefficient $\vvec_{K_\ell}\in\R^2$ is now fixed so that condition~\eqref{eq:WellPosed-nc-2} is satisfied. We can impose for example that the projection $\Pi_{\FS}\vvec$ satisfies:
\begin{equation}\label{eq:linsyst_vKl}
\int_{K_\ell} T_\ell^{-1}(\xvec) \div\Pi_{\FS}\vvec=\int_{K_\ell} T_\ell^{-1}(\xvec) \div\vvec.
\end{equation}
Notice that due to~\eqref{eq:phiKl}, the patch-test condition is still satisfied.

\begin{prop}\label{pro:rh-FS}
The Fortin--Soulie interpolation operator $\Pi_{\FS}$ is such for all $0\leq s\leq 1$, for all $\vvec\in\bH^{1+s}(\Om)$, we have:
\begin{align}\label{eq:rh-FS-0}
\forall\ell\in\Ical_K,\quad
\bigl\lVert\Grad\bigl(\Pi_{\FS}\vvec-\vvec\bigr)\bigr\rVert_{\bbL^2(K_\ell)}&amp;\lesssim(\sigma_\ell)^2 (h_\ell)^s \abs{\vvec}_{\bH^{1+s}(K_\ell)},
\\\label{eq:rh-FS}
\bigl\lvert\Pi_{\FS}\vvec-\vvec\bigr\rvert_h&amp;\lesssim\sigma^2 h^s \abs{\vvec}_{\bH^{1+s}(\Om)}.
\end{align}
\end{prop}
\begin{rema}
Albeit we are inspired by the proof of~\cite[Lemma~4]{CrRa73}, we changed the transition from Equation~(4.27) to~(4.29) there by using only the properties related to the normal component of the velocity, cf.~\eqref{eq:IPPnormcomp}. In the original proof, one needs a stronger assumption on the regularity of $\vvec$ (namely, $\vvec\in\bigcap_{0&lt;s&lt;s_\Om}\bH^{1+s}(\Om)$ with $s_\Om&gt;\frac12$). Finally, because we do not split the integral over the boundaries of elements into the sum of $d+1$ integrals over the facets, we obtain purely local estimates, which appear to be new for the Fortin--Soulie element in the case of low-regularity fields $\vvec$.
\end{rema}
\begin{proof}[Proof of Proposition~\ref{pro:rh-FS}]
Let $\vvec\in\bH^1(\Om)$. We have:
\begin{equation}\label{eq:gradPiFSv}
\bigl\lVert\Grad\bigl(\Pi_{\FS}\vvec-\vvec\bigr)\bigr\rVert_{\bbL^2(K_\ell)}
\leq
\bigl\lVert\Grad \bigl(\Pi_{\FS}\vvec-\tilde{\Pi}_{\FS}\vvec\bigr)\bigr\rVert_{\bbL^2(K_\ell)}+\norm*{\Grad \bigl(\tilde{\Pi}_{\FS}\vvec-\vvec\bigr)}_{\bbL^2(K_\ell)}.
\end{equation}
Notice that for all $\ell\in\Ical_K$, $\restr{\bigl(\Pi_{\FS}\vvec-\tilde{\Pi}_{\FS}\vvec\bigr)}{K_\ell}=(\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell}) \phi_{K_\ell}$. Using~\eqref{eq:gradvK-i}, we obtain that:
\begin{equation}\label{eq:gradPtvmPv}
\bigl\lVert\Grad \bigl(\Pi_{\FS}\vvec-\tilde{\Pi}_{\FS}\vvec\bigr)\bigr\rVert_{\bbL^2(K_\ell)}
\lesssim \abs{\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell}} \norm*{\grad \phi_{K_\ell}}_{\bL^2(K_\ell)}
\lesssim\norm{{\bbB_\ell}^{-1}} \abs{K_\ell}^{\frac12} \abs{\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell}}.
\end{equation}
Let us estimate $\abs{\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell}}$. On the one hand, we have\footnote{Let $\wvec=\Pi_{\FS}\vvec-\vvec=(w_1,w_2)^T$ and $\xvec\coloneqq (x_1,x_2)^T$. By integration by parts, we have for $d'=1,2$: $\int_{\pa K_\ell}x_{d'} \wvec\cdot\restr{\nvec}{\pa K_\ell}=\int_{K_\ell}x_{d'} \div\wvec+\int_{K_\ell}\wvec\cdot\evec_{d'}$. Due to~\eqref{eq:linsyst_vKl}, $\int_{K_\ell}x_{d'} \div\wvec=0$, so that for $d'=1,2$: $\int_{K_\ell}w_{d'}=\int_{\pa K_\ell}x_{d'} \wvec\cdot\restr{\nvec}{\pa K_\ell}$.}:
\[
\begin{aligned}
\int_{K_\ell}\bigl(\Pi_{\FS}\vvec-\tilde{\Pi}_{\FS}\vvec\bigr)&amp;=\int_{K_\ell}\bigl(\Pi_{\FS}\vvec-\vvec\bigr)\qquad\qquad&amp;&amp;\text{ from~\eqref{eq:Piave},}\\
&amp;=\int_{\pa K_\ell}\xvec \bigl(\Pi_{\FS}\vvec-\vvec\bigr)\cdot\restr{\nvec}{\pa K_\ell}&amp;&amp;\text{ by IBP and using~\eqref{eq:linsyst_vKl}},\\
&amp;=\int_{\pa K_\ell}(\xvec-\ul{\xvec}) \bigl(\Pi_{\FS}\vvec-\vvec\bigr)\cdot\restr{\nvec}{\pa K_\ell}&amp;&amp;\text{ since }\int_{\pa K_\ell}\bigl(\Pi_{\FS}\vvec-\vvec\bigr)\cdot\restr{\nvec}{\pa K_\ell}=0,\\
&amp;=\int_{\pa K_\ell}(\xvec-\ul{\xvec}) \bigl(\tilde{\Pi}_{\FS}\vvec-\vvec\bigr)\cdot\restr{\nvec}{\pa K_\ell}&amp;&amp;\text{ from~\eqref{eq:phiKl}}.
\end{aligned}
\]
On the other hand, it holds: $\int_{K_\ell}\bigl(\Pi_{\FS}\vvec-\tilde{\Pi}_{\FS}\vvec\bigr)=(\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell})\int_{K_\ell}\phi_{K_\ell}=\frac{\abs{K_\ell}}{4} (\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell})$.

Hence we have:
\begin{equation}\label{eq:tobound}
\abs{\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell}}\leq 4 \abs{K_\ell}^{-1} \abs*{\int_{\pa K_\ell}(\xvec-\ul{\xvec}) \bigl(\tilde{\Pi}_{\FS}\vvec-\vvec\bigr)\cdot\restr{\nvec}{\pa K_\ell}}.
\end{equation}
In order to bound the right-hand side of~\eqref{eq:tobound} component by component, we can use Lemma~\ref{lem:IPPnormcomp}, with $q=x_{d'}-\int_{K_\ell}x_{d'}/\abs{K_\ell}$ $(d'=1,2)$, so that $\norm{\grad q}_{\bL^2(K_\ell)}=\abs{K_\ell}^{\frac12}$. We obtain:
\begin{equation}\label{eq:bound-vkl-tvkl-sing}
\abs{\vvec_{K_\ell}-\tilde{\vvec}_{K_\ell}}\leq4 d\times\bigl(\sqrt{d}+1\bigr) \abs{K_\ell}^{-\frac12} \pi^{-1}  h_\ell \norm*{\Grad\bigl(\tilde{\Pi}_{\FS}\vvec-\vvec\bigr)}_{\bbL^2(K_\ell)}.
\end{equation}
Combining~\eqref{eq:gradPtvmPv} and~\eqref{eq:bound-vkl-tvkl-sing}, we have:
\[
\bigl\lVert\Grad\bigl(\Pi_{\FS}\vvec-\tilde{\Pi}_{\FS}\vvec\bigr)\bigr\rVert_{\bbL^2(K_\ell)}
\lesssim\norm{{\bbB_\ell}^{-1}} h_\ell \norm*{\Grad\bigl(\tilde{\Pi}_{\FS}\vvec-\vvec\bigr)}_{\bbL^2(K_\ell)}\lesssim\sigma_\ell \norm*{\Grad\bigl(\tilde{\Pi}_{\FS}\vvec-\vvec\bigr)}_{\bbL^2(K_\ell)}.
\]
For all $\vvec\in \bP^2(K_\ell)$ we have $\tilde{\Pi}_{\FS}(\vvec)=\vvec$ and $\hat{\oldtilde{\Pi}}_{\FS}\hat{\vvec}_\ell=\hat{\vvec}_\ell=\hat{\vvec}_\ell$. Hence, using Bramble--Hilbert/Deny--Lions Lemma~\cite[Lemma~11.9]{ErGu21-I}, we have for $m=0,1$, for all $\vvec\in\bH^{1+m}(\Om)$:
\[
\forall\ell\in\Ical_K,\quad\norm*{\Grad\bigl(\tilde{\Pi}_{\FS}\vvec-\vvec\bigr)}_{\bbL^2(K_\ell)}\lesssim\sigma_\ell (h_\ell)^m \abs{\vvec}_{\bH^{1+m}(K_\ell)}.
\]
We deduce that for $m=0,1$, for all $\vvec\in\bH^{1+m}(\Om)$:
\begin{equation}
\begin{aligned}
\forall\ell\in\Ical_K,\quad\bigl\lVert\Grad\bigl(\Pi_{\FS}\vvec-\vvec\bigr)\bigr\rVert_{\bbL^2(K_\ell)}&amp;\lesssim (\sigma_\ell)^2 (h_\ell)^m \abs{\vvec}_{\bH^{1+m}(K_\ell)},\\
\bigl\lVert\Grad\bigl(\Pi_{\FS}\vvec-\vvec\bigr)\bigr\rVert_{\bbL^2(\Om)}&amp;\lesssim \sigma^2 h^m \abs{\vvec}_{\bH^{1+m}(\Om)}\;\text{ by summation.}
\end{aligned}
\end{equation}
Using interpolation property~\cite[Lemma~22.2]{Tart07}, we obtain~\eqref{eq:rh-FS-0} and~\eqref{eq:rh-FS}.
\end{proof}
Hence, using the triangular inequality, we have:
\[
\bigl\lVert\Grad\Pi_{\FS}\vvec\bigr\rVert_{\bbL^2(K_\ell)}\leq\bigl\lVert\Grad\bigl(\Pi_{\FS}\vvec-\vvec\bigr)\bigr\rVert_{\bbL^2(K_\ell)}+\bigl\lVert\Grad\vvec\bigr\rVert_{\bbL^2(K_\ell)}\lesssim \bigl((\sigma_\ell)^2+1\bigr)\norm{\vvec}_{\bH^1(K_\ell)}.
\]
By summation over $\ell$, we deduce that the coefficient $C_{\nonconf}$ in~\eqref{eq:WellPosed-nc-1} is here equal to $\sigma^2+1$. We recall that the discrete PoincarÃ©--Steklov inequality~\eqref{eq:cps-constant} holds.
\begin{theo}
Let $\Xcal_h=\Xcal_{\FS}$. Then the continuous bilinear form $a_{S,h}({\,\cdot\,,\cdot\,})$ is $T_h$-coercive and Problem~\eqref{eq:Stokes-FVah} is well-posed.
\end{theo}
\begin{proof}
According to Proposition~\ref{pro:rh-FS}, the Fortin--Soulie interpolation operator $\Pi_{\FS}$ satisfies~\eqref{eq:WellPosed-nc-1}-\eqref{eq:WellPosed-nc-2}, so that we can apply the proof of Theorem~\ref{thm:WellPosed-nc}.
\end{proof}

Notice that in the recent paper~\cite{SaTo23}, the inf-sup condition of the mixed Fortin--Soulie finite element is proved directly on a triangle and then using the macro-element technique~\cite{Sten90}, but it seems difficult to use this technique to build a Fortin operator, which is needed to compute error estimates.

The study can be extended to higher orders for $d=2$ using the following papers:~\cite{BaSt07} for $k\geq 4$, $k$ even, \cite{CaSa22a} for $k=3$ and~\cite{CaSa22} for $k\geq 5$, $k$ odd. In~\cite{DiST22}, the authors propose a local Fortin operator for the lowest order Taylor--Hood finite element~\cite{TaHo73} for $d=3$.

\section{Numerical results improving consistency}\label{sec:numerical_results}

\subsection{\texorpdfstring{$\bH(\div)$}{bH(div)}-conforming velocity reconstruction}\label{sec:reco}

Consider Problem~\eqref{eq:Stokes} with data $\fvec=-\grad\phi$, where $\phi\in H^1(\Om)\cap L^2_{\zmv}(\Om)$. The unique solution is then $(\uvec,p)\coloneqq (0,\phi)$. By integrating by parts, the source term in~\eqref{eq:Stokes-FV} reads:
\begin{equation}\label{eq:contIPP}
\forall\vvec\in\bH^1_0(\Om),\;\int_\Om\fvec\cdot\vvec=\int_\Om\phi \div\vvec.
\end{equation}
Recall that the nonconforming space $\bX_h$ defined in~\eqref{eq:Xh} is a subset of $\Pcal_h\bH^1$: using a nonconforming finite element method, the integration by parts must be done on each element of the triangulation, and we have:
\begin{equation}\label{eq:discIPP}
\forall\vvec\in\Pcal_h\bH^1,\;\int_\Om\fvec\cdot\vvec=\bigl(\div_h\vvec,\phi\bigr)_{L^2(\Om)}+\sum_{f\in\Ical_F}\int_{F_f}[\vvec]_{F_f}\cdot\nvec_f \phi.
\end{equation}
Using Lemma~\ref{lem:IPPnormcomp}, we have: $\sum_{f\in\Ical_F}\int_{F_f}[\vvec]_{F_f}\cdot\nvec_f \phi\lesssim h \norm{\vvec_h}_h \norm{\grad\phi}_{\bL^2(\Om)}$. Applying~\eqref{eq:discIPP} to the right-hand side of~\eqref{eq:Stokes-FVah} and choosing $\vvec_h=\uvec_h$, it holds: $\nu \norm{\uvec_h}_h\lesssim h\norm{\grad\phi}_{\bL^2(\Om)}$ (as expected by~\eqref{eq:aprioriEEuh1}). Hence, the term with the jumps acts as a numerical source for the discrete velocity, whose numerical influence is proportional to $h/\nu$. Thus, we cannot obtain exactly $\uvec_h=0$. Linke proposed in~\cite{Link14} to project the test function $\vvec_h\in\bX_h$ on a discrete subspace of $\bH(\div; \Om)$, like Raviart--Thomas or Brezzi--Douglas--Marini finite elements (see~\cite{RaTh77,BrDM85}, or the monograph~\cite{BoBF13}). Let $\Pi_{\div}\colon\bX_{0,h}\to \mathbf{P}^k_{\disc}(\Tcal_h)\cap\bH_0(\div;\Om)$ be some interpolation operator built so that for all $\vvec_h\in\bX_{0,h}$, for all $\ell\in\Ical_K$, $\restr{\bigl(\div\Pi_{\div}\vvec_h\bigr)}{K_\ell}=\restr{\div\vvec_h}{K_\ell}$. Integrating by parts, we have for all $\vvec_h\in\bX_{0,h}$:
\[
\int_\Om\fvec\cdot\Pi_{\div}\vvec_h=\int_\Om\phi \div\Pi_{\div}\vvec_h=\sum_{\ell\in K_\ell}\int_{K_\ell}\phi \div\Pi_{\div}\vvec_h=\sum_{\ell\in K_\ell}\int_{K_\ell}\phi \div\vvec_h=\bigl(\div_h\vvec_h,\phi\bigr)_{L^2(\Om)}.
\]
The projection $\Pi_{\div}$ allows to eliminate the terms of the integrals of the jumps in~\eqref{eq:discIPP}.

Let us write Problem~\eqref{eq:Stokes-FVah} as:
\begin{equation}\label{eq:Stokes-FVah-RT}
\text{Find $(\uvec_h,p_h)\in\Xcal_h$ such that $a_{S,h}\bigl((\uvec_h,p_h),(\vvec_h,q_h)\bigr)=\ell_\fvec\bigl((\Pi_{\div}\vvec_h,q_h)\bigr)\;\forall(\vvec_h,q_h)\in\Xcal_h$.}
\end{equation}
In the case of $\Xcal_h=\Xcal_{\CR}$ and a projection on Brezzi--Douglas--Marini finite elements, the following error estimate holds if $(\uvec,p)\in\bH^2(\Om)\times H^1(\Om)$:
\begin{equation}\label{eq:aprioriEE-PiBDM}
\norm{\uvec-\uvec_h}_{\bL^2(\Om)}\leq \widetilde{C} h^2 \abs{\uvec}_{\bH^2(\Om)},
\end{equation}
where the constant $\widetilde{C}$ is independent of $h$. The proof is detailed in~\cite{BLMS15} for shape-regular meshes and~\cite{AKLM22} for anisotropic meshes. We remark that the error doesn't depend on the norm of the pressure nor on the $\nu$ parameter. We will provide some numerical results to illustrate the effectiveness of this formulation, even with a projection on the Raviart--Thomas finite elements, which, for a fixed polynomial order, are less precise than the Brezzi--Douglas--Marini finite elements.

%For all $\ell\in\Ical_K$, we let $P^k_H(K_\ell)$ be the set of homogeneous polynomials of order $k$ on $K_\ell$.

For $k\in\N^\star$, the space of Raviart--Thomas finite elements can be defined as:
%\[
%\bX_{\RT_k}\coloneqq \Bigl\{\vvec\in\bH(\div; \Om) \,\Big\vert\, \forall \ell\in\Ical_k, \; \restr{\vvec}{K_\ell}=\avec_\ell+b_\ell\xvec, \; (\avec_\ell,b_\ell)\in P^k(K_\ell)^d\times P^k_H(K_\ell)\Bigr\}.
%\]
\[
\bX_{\RT_k}\coloneqq \Bigl\{\vvec\in\bH(\div; \Om) \,\Big\vert\, \forall \ell\in\Ical_k, \; \restr{\vvec}{K_\ell}=\avec_\ell+b_\ell\xvec, \; (\avec_\ell,b_\ell)\in P^k(K_\ell)^d\times P^k(K_\ell)\Bigr\}.
\]
Let $k\leq 1$. The Raviart--Thomas interpolation operator $\Pi_{\RT_k}\colon\bH^1(\Om)\cup\bX_h\to\bX_{\RT_k}$ is defined by:
\begin{equation}\label{eq:PiRTdef}
\forall\vvec\in\bH^1(\Om)\cup\bX_h, \;
\left\{
\begin{aligned}
\forall f\in\Ical_F,\;&amp;\int_{F_f}\Pi_{\RT_k}\vvec\cdot\nvec_f q=\int_{F_f}\vvec\cdot\nvec_f q,\;\forall q\in P^k(F_f),\\
\text{for }\;k=1, \;\forall\ell\in\Ical_K,&amp;\int_{K_\ell}\Pi_{\RT_1}\vvec=\int_{K_\ell}\vvec.
\end{aligned}
\right.
\end{equation}
Note that the Raviart--Thomas interpolation operator preserves the constants. Let $\vvec_h\in\bX_h$. In order to compute the left-hand side of~\eqref{eq:discIPP}, we must evaluate $\restr{(\Pi_{\RT_k}\vvec_h)}{K_\ell}$ for all $\ell\in\Ical_K$. Calculations are performed using the proposition below:
\begin{prop}[{\cite[Lemma~3.11]{Gati14}}]\label{prop:PiRTref}
Let $k\leq 1$. Let $\hat{\Pi}_{\RT_k}\colon\bH^1(\hatK)\to\bP^k(\hatK)$ be the Raviart--Thomas interpolation operator restricted to the reference element, so that:
\begin{equation}\label{eq:PiRTref}
\forall\hat{\vvec}\in\bH^1(\hatK), \;
\left\{
\begin{aligned}
&amp;\forall \hatF\in\pa\hatK,\;
\int_{\hatF}\hat{\Pi}_{\RT_k}\hat{\vvec}\cdot\nvec_{\hatF} \hat{q}=\int_{\hatF}\hat{\vvec}\cdot\nvec_{\hatF} \hat{q},\;\forall \hat{q}\in P^k(\hatF),\\
&amp;\text{for }\;k=1,\;\int_{\hatK}\hat{\Pi}_{\RT_k}\hat{\vvec}=\int_{\hatK}\hat{\vvec}.
\end{aligned}
\right.
\end{equation}
We then have: $\forall\ell\in\Ical_K$,
\begin{equation}\label{eq:PiRTrefPiRT}
\restr{\bigl(\Pi_{\RT_k}\vvec\bigr)}{K_\ell}(\xvec)=
\bbB_\ell \bigl(\hat{\Pi}_{\RT_k}\bbB_\ell{}^{-1}\hat{\vvec}_\ell\bigr)\circ T_\ell{}^{-1}(\xvec)\;\text{ where }\;\hat{\vvec}_\ell=\vvec\circ T_\ell(\hat{\xvec}).
\end{equation}
\end{prop}

The proof is based on the equality of the $\hatF$ and $\hatK$-moments of $\restr{\bigl(\Pi_{\RT_k}\vvec\bigr)}{K_\ell}\circ T_\ell(\hat{\xvec})$ and $\bbB_\ell\bigl( \hat{\Pi}_{\RT_k}\bbB_\ell{}^{-1}\hat{\vvec}_\ell \bigr)(\hat{\xvec})$. For $k=0$, setting for $d'\in\bigl\{1,\dotsc,d\bigr\}$: $\boldsymbol{\psi}_{f,{d'}}\coloneqq \psi_f \evec_{d'}$, we obtain that:
\begin{equation}\label{eq:PiRTCR}
\forall\ell\in\Ical_K , \forall f\in\Ical_{F,\ell} , \; \restr{\bigl(\Pi_{\RT_0}\boldsymbol{\psi}_{f,{d'}}\bigr)}{K_\ell}=\bigl(d \abs{K_\ell}\bigr)^{-1} \bigl(\xvec-\vec{OS}_{f,\ell}\bigr)\Scal_{f,\ell}\cdot\evec_{d'} ,
\end{equation}
where $S_{f,\ell}$ is the vertex opposite to $F_f$ in $K_\ell$.

For $k=1$, the vector $\restr{\bigl(\Pi_{\RT_1}\vvec_h\bigr)}{K_\ell}$ is described by eight unknowns:
\[
\restr{\bigl(\Pi_{\RT_1}\vvec_h\bigr)}{K_\ell}=\bbA_\ell \xvec+(\bvec_\ell\cdot\xvec) \xvec+\dvec_\ell,\;\text{ where }\;\bbA_\ell\in\R^{2\times2}, \bvec_\ell\in\R^2, \dvec_\ell\in\R^2.
\]
We compute only once the inverse of the matrix of the linear system~\eqref{eq:PiRTref}, in $\R^{8\times 8}$.

\subsection{Application with manufactured solutions}\label{sec:resu}

In the Tables~\ref{tab:grad}, \ref{tab:cossin-3} and~\ref{tab:cossin-4}, we call $\eps_0^\nu(\uvec)=\norm{\uvec-\uvec_h}_{\bL^2(\Om)}/\norm{(\uvec,p)}_{\Xcal}$ the velocity error in $\bL^2(\Om)$-norm, where $\uvec_h$ is the solution to Problem~\eqref{eq:Stokes-FVah} (columns $\bX_{\CR}$ and $\bX_{\FS}$) or~\eqref{eq:Stokes-FVah-RT} (columns $\bX_{\CR}+\Pi_{\RT_0}$ and $\bX_{\FS}+\Pi_{\RT_1}$) and $h$ is the mesh size.

We first consider Stokes problem~\eqref{eq:Stokes} in $\Omega=(0,1)^2$ with $\uvec=0$, $p=(x_1)^3+(x_2)^3-0.5$, $\fvec=\grad p=3 \bigl((x_1)^2,(x_2)^2\bigr)^T$. %
%We report in Table~\ref{tab:grad} $\eps_0^\nu(\uvec)$ for $h=5.00\,e-2$ and for different values of $\nu$.
We report in Table~\ref{tab:grad} the velocity error $\varepsilon_0^\nu(\mathbf{u})\coloneqq\nu\norm{\mathbf{u}_h}_{\mathbf{L}^2(\Omega)}/\norm{p}_{L^2(\Omega)}$, for $h=5.00\times 10^{-2}$ and for different values of $\nu$.
%\begin{table}[h]
%\caption{\label{tab:grad}
%Values of $\eps_0^\nu(\uvec)$ for $h=5.00\,e-2$.}
%\begin{tabular}{ccccc}
%\hline $\nu$&amp;$\bX_{\CR}$&amp;$\bX_{\CR}+\Pi_{\RT_0}$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$\\
%\hline $1.00\,e-4$ &amp; $7.96\,e-4$ &amp; $4.59\,e-17$ &amp; $8.81\,e-7$ &amp; $1.54\,e-16$\\
%$1.00\,e-5$ &amp; $7.96\,e-4$ &amp; $4.59\,e-17$ &amp; $8.81\,e-7$ &amp; $1.54\,e-16$\\
%$1.00\,e-6$ &amp; $7.96\,e-4$ &amp; $4.59\,e-17$ &amp; $8.81\,e-7$ &amp; $1.54\,e-16$\\
%\hline
%\end{tabular}
%\end{table}
\begin{table}[h]
\caption{\label{tab:grad}
Values of $\eps_0^\nu(\uvec)$ for $h=5.00\times 10^{-2}$.}
\begin{tabular}{ccccc}
\hline $\nu$&amp;$\bX_{\CR}$&amp;$\bX_{\CR}+\Pi_{\RT_0}$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$\\
\hline &amp;&amp;&amp;&amp;\\
[-1.1em] $1.00\times 10^{-4}$ &amp; $7.96\times 10^{-4}$ &amp; $4.59\times 10^{-17}$ &amp; $8.81\times 10^{-7}$ &amp; $1.54\times 10^{-16}$\\
$1.00\times 10^{-5}$ &amp; $7.96\times 10^{-4}$ &amp; $4.59\times 10^{-17}$ &amp; $8.81\times 10^{-7}$ &amp; $1.54\times 10^{-16}$\\
$1.00\times 10^{-6}$ &amp; $7.96\times 10^{-4}$ &amp; $4.59\times 10^{-17}$ &amp; $8.81\times 10^{-7}$ &amp; $1.54\times 10^{-16}$\\
\hline
\end{tabular}
\end{table}

%Here we have $\norm{(\uvec,p)}_{\Xcal}=\nu\norm{p}_{L^2(\Om)}$. Hence, the $\bL^2(\Om)$-norm of the discrete velocity $\norm{\uvec_h}_{\bL^2(\Om)}$ is proportional to $\nu^{-1}$. Using the projection, we obtain $\eps_0^\nu(\uvec)=0$ close to machine precision.
The $\bL^2(\Om)$-norm of the discrete velocity $\norm{\uvec_h}_{\bL^2(\Om)}$ is proportional to $\nu^{-1}$. Using the projection, we obtain $\eps_0^\nu(\uvec)=0$ close to machine precision.

We now consider Stokes problem~\eqref{eq:Stokes} in $\Omega=(0,1)^2$ with:
\[
\uvec=\Biggl(
\begin{matrix}
\bigl(1-\cos(2\,\pi\,x_1)\bigr)\sin(2\,\pi\,x_2)\\
\bigl(\cos(2\,\pi\,x_2)-1\bigr)\sin(2\,\pi\,x_1)
\end{matrix}\Biggr),\;
p=\sin(2\,\pi\,x_1)\;\text{ and }\;\fvec=-\nu\Delta\uvec+\grad p.
\]
We report in Table~\ref{tab:cossin-3} (resp.~\ref{tab:cossin-4}) the values of $\eps_0^\nu(\uvec)$ in the case $\nu=1.00\times 10^{-3}$ (resp.\ $\nu=1.00\times 10^{-4}$) for mesh sizes. We observe that when there is no projection, $\eps_0^\nu(\uvec)$ is independent of $\nu$, whereas using the projection, $\eps_0^\nu(\uvec)$ is proportional to $\nu$.
\begin{table}[h]
\caption{\label{tab:cossin-3}
Values of $\eps_0^\nu(\uvec)$ for $\nu=1.00\times 10^{-3}$.}
\begin{tabular}{ccccc}
\hline $h$&amp;$\bX_{\CR}$&amp;$\bX_{\CR}+\Pi_{\RT_0}$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$\\
\hline &amp;&amp;&amp;&amp;\\
[-1.1em] $5.00\times 10^{-2}$ &amp; $1.32\times 10^{-3}$ &amp; $2.74\times 10^{-5}$ &amp; $4.73\times 10^{-6}$ &amp; $5.05\times 10^{-7}$ \\
$2.50\times 10^{-2}$ &amp; $3.30\times 10^{-4}$ &amp; $6.93\times 10^{-6}$ &amp; $5.06\times 10^{-7}$ &amp; $6.42\times 10^{-8}$ \\
$1.25\times 10^{-2}$ &amp; $8.25\times 10^{-5}$ &amp; $1.74\times 10^{-6}$ &amp; $6.31\times 10^{-8}$ &amp; $8.10\times 10^{-9}$ \\
$6.25\times 10^{-3}$ &amp; $2.04\times 10^{-5}$ &amp; $4.35\times 10^{-7}$ &amp; $7.44\times 10^{-9}$ &amp; $1.03\times 10^{-9}$ \\
Rate &amp; $h^{2.00}$ &amp; $h^{1.99}$ &amp; $h^{3.08}$ &amp; $h^{2.97}$ \\
\hline
\end{tabular}
%\end{table}
%\begin{table}[h]
%\caption{\label{tab:cossin-4}
%Values of $\eps_0^\nu(\uvec)$ for $\nu=1.00\,e-4$.}
%\begin{tabular}{ccccc}
%\hline $h$&amp;$\bX_{\CR}$&amp;$\bX_{\CR}+\Pi_{\RT_0}$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$\\
%\hline $5.00\,e-2$ &amp; $1.32\,e-3$ &amp; $2.74\,e-6$ &amp; $4.70\,e-6$ &amp; $5.05\,e-8$\\
%$2.50\,e-2$ &amp; $3.30\,e-4$ &amp; $6.93\,e-7$ &amp; $5.10\,e-7$ &amp; $6.43\,e-9$\\
%$1.25\,e-2$ &amp; $8.25\,e-5$ &amp; $1.74\,e-7$ &amp; $6.37\,e-8$ &amp; $8.11\,e-10$\\
%$6.25\,e-3$ &amp; $2.04\,e-5$ &amp; $4.36\,e-8$ &amp; $7.51\,e-9$ &amp; $9.77\,e-11$\\
%Rate &amp; $h^{2.00}$ &amp; $h^{1.99}$ &amp; $h^{3.08}$ &amp; $h^{2.99}$ \\
%\hline
%\end{tabular}
%\end{table}
\end{table}
\begin{table}[h]
\caption{\label{tab:cossin-4}
Values of $\eps_0^\nu(\uvec)$ for $\nu=1.00\times 10^{-4}$.}
\begin{tabular}{ccccc}
\hline $h$&amp;$\bX_{\CR}$&amp;$\bX_{\CR}+\Pi_{\RT_0}$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$\\
\hline &amp;&amp;&amp;&amp;\\
[-1.1em] $5.00\times 10^{-2}$ &amp; $1.32\times 10^{-3}$ &amp; $2.74\times 10^{-6}$ &amp; $4.70\times 10^{-6}$ &amp; $5.05\times 10^{-8}$\\
$2.50\times 10^{-2}$ &amp; $3.30\times 10^{-4}$ &amp; $6.93\times 10^{-7}$ &amp; $5.10\times 10^{-7}$ &amp; $6.43\times 10^{-9}$\\
$1.25\times 10^{-2}$ &amp; $8.25\times 10^{-5}$ &amp; $1.74\times 10^{-7}$ &amp; $6.37\times 10^{-8}$ &amp; \hspace{4pt}$8.11\times 10^{-10}$\\
$6.25\times 10^{-3}$ &amp; $2.04\times 10^{-5}$ &amp; $4.36\times 10^{-8}$ &amp; $7.51\times 10^{-9}$ &amp; \hspace{4pt}$9.77\times 10^{-11}$\\
Rate &amp; $h^{2.00}$ &amp; $h^{1.99}$ &amp; $h^{3.08}$ &amp; $h^{2.99}$ \\
\hline
\end{tabular}
\end{table}

Let us consider Stokes problem~\eqref{eq:Stokes} with a low-regular velocity. Let $\Omega=(0,1)^2$, $S_0=(0.5,0.5)$, and $(r,\theta)$ be the polar coordinates centered on $S_0$. We set:
\[
\uvec=r^\alpha\evec_\theta,\; p=r-\abs{\Omega}^{-1}\int_\Omega r\; \text{ so that }\; \fvec\coloneqq -\nu\Deltavec\uvec+\grad p=\nu(1-\alpha^2)r^{\alpha-2}\evec_\theta+\evec_r.
\]
We report in Table~\ref{tab:vortex} the values of $\eps_0^\nu(\uvec)$ for $\nu=1.00\times 10^{-4}$, and for different mesh sizes, with $\alpha=1$ and $\alpha=0.45$. For $\alpha=1$, $\uvec=(-y,x)^T\in\bH^2(\Om)$. For $\alpha=0.45$, $\uvec\in \bigcap_{0&lt;s&lt;\alpha}\bH^{1+s}(\Om)$, hence $\uvec\notin\bH^2(\Om)$. It seems that the Raviart--Thomas projection is less efficient in that last case.

%\begin{table}[h]
%\caption{\label{tab:vortex}
%Values of $\eps_0^\nu(\uvec)$, regular and low-regular velocity, $\nu=1.00\,e-4$.}
%\begin{tabular}{ccccc}
%\hline &amp;\multicolumn{2}{c}{$\alpha=1$}&amp;\multicolumn{2}{c}{$\alpha=0.45$}\\
%\hline $h$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$\\
%\hline $1.00\,e-1$ &amp; $3.03\,e-5$ &amp; $2.81\,e-6$ &amp; $3.05\,e-5$ &amp; $3.94\,e-6$\\
%$5.00\,e-2$ &amp; $4.34\,e-6$ &amp; $1.54\,e-6$ &amp; $4.57\,e-6$ &amp; $2.15\,e-6$\\
%$2.50\,e-2$ &amp; $4.72\,e-7$ &amp; $2.41\,e-8$ &amp; $9.70\,e-7$ &amp; $8.52\,e-7$\\
%Rate &amp; $h^{3.00}$ &amp; $h^{3.43}$ &amp; $h^{2.48}$ &amp; $h^{1.11}$\\
%\hline
%\end{tabular}
%\end{table}
\begin{table}[h]
\caption{\label{tab:vortex}
Values of $\eps_0^\nu(\uvec)$, regular and low-regular velocity, $\nu=1.00\times 10^{-4}$.}
\begin{tabular}{ccccc}
\hline &amp;\multicolumn{2}{c}{$\alpha=1$}&amp;\multicolumn{2}{c}{$\alpha=0.45$}\\
\hline $h$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$&amp;$\bX_{\FS}$&amp;$\bX_{\FS}+\Pi_{\RT_1}$\\
\hline &amp;&amp;&amp;&amp;\\
[-1.1em] $1.00\times 10^{-1}$ &amp; $3.03\times 10^{-5}$ &amp; $2.81\times 10^{-6}$ &amp; $3.05\times 10^{-5}$ &amp; $3.94\times 10^{-6}$\\
$5.00\times 10^{-2}$ &amp; $4.34\times 10^{-6}$ &amp; $1.54\times 10^{-6}$ &amp; $4.57\times 10^{-6}$ &amp; $2.15\times 10^{-6}$\\
$2.50\times 10^{-2}$ &amp; $4.72\times 10^{-7}$ &amp; $2.41\times 10^{-8}$ &amp; $9.70\times 10^{-7}$ &amp; $8.52\times 10^{-7}$\\
Rate &amp; $h^{3.00}$ &amp; $h^{3.43}$ &amp; $h^{2.48}$ &amp; $h^{1.11}$\\
\hline
\end{tabular}
\end{table}

In order to enhance the numerical results, one could also use a posteriori error estimators to adapt the mesh near point $S_0$ (see~\cite{DaDP95,DoAi05} for $k=1$ and~\cite{AABR12} for $k=2$).

Alternatively, using the nonconforming Crouzeix--Raviart mixed finite element method, one can build a divergence-free basis, as described in~\cite{Hech81} for $k=1$. When $k=1$, following the initial work of~\cite{BeHe00}, one can also add $P^1$-Lagrange basis functions to the space of the discrete pressures as explained in~\cite{JaCS24}. The consistency of the discrete velocity is then improved. Notice that using conforming finite elements, the Scott--Vogelius finite elements~\cite{ScVo85,Zhan05,FGNZ22} produce velocity approximations that are exactly divergence free.

The code used to get the numerical results can be downloaded on GitHub~\cite{Jame22}.

\section{Conclusion}
We analysed the discretization of Stokes problem with nonconforming finite elements in light of the T-coercivity theory. Furthermore, we obtained local stability estimates for order 1 in 2 or 3 dimensions without mesh regularity assumption; and for order 2 in 2 dimensions in the case of a shape-regular triangulation sequence. This local approach, splitting the normal and the tangential components could help to generalize our results to order $k\geq 3$ (using maybe also other internal moment conservation). This is ongoing work. We then provided numerical results to illustrate the importance of using $\bH(\div)$-conforming projection. Further, we intend to extend the study to other mixed finite element methods.

\section*{Acknowledgements}

The author thanks Mahran Rihani and AlbÃ©ric Lefort, and she is grateful to the reviewers for their careful assessment.

\section*{Declaration of interests}

The author does not work for, advise, own shares in, or receive funds from any organization that could benefit from this article, and has declared no affiliation other than her research organization.

\printbibliography
\end{document}</pre></body></html>