We consider the incompressible Navier-Stokes equations in a two-dimensional channel. The tangential and normal velocities are assumed to be periodic in the streamwise (horizontal) direction. Moreover, we consider no-slip boundary conditions on the tangential velocity at the top and bottom walls of the channel, and normal velocity actuation at the top and bottom walls. For an arbitrarily large Reynolds number, we design the boundary control inputs to achieve global exponential stabilization, in the $L^2$ sense, of a chosen parabolic Poiseuille profile. Moreover, we design the control inputs such that they have zero mean, but non-zero cubic mean. The zero-mean property is to ensure that the conservation of mass constraint is verified. The non-zero cubic mean property is the key to exploiting the stabilizing effect of nonlinear convection and achieving global stabilization independently of the size of the Reynolds number. This paper is not only the first work where a closed-form feedback law is proposed for global stabilization of parabolic Poiseuille profiles for arbitrary Reynolds number but is also the first generalization of the Cardano-Lyapunov formula, designed initially to stabilize scalar-valued convective PDEs, to a vector-valued convective PDE with a divergence-free constraint on the state.

§ INTRODUCTION

Turbulence is an undesired phenomenon in many fluid engineering applications, including fluid flow through channels. It is characterized by a chaotic movement of the fluid particles that leads, in the case of channel flows, to large friction at the boundary walls. Turbulence also leads, via the generation of vortexes, to the vibration of the structures through which the flow passes, leading to deterioration due to fatigue. Turbulence also plays a crucial role in the prediction of climate change; see e.g. [1, 2], where turbulence is involved in the problem of estimation of basal sliding of the Antarctic ice sheet.

One of the main objectives of fluid flow control is to act on the fluid system, by some control action, to avoid turbulence and impose a desired regular profile on the flow.

Incompressible flows are described by the incompressible Navier-Stokes equations. These equations involve a positive coefficient, which we refer to as the Reynolds number, which is a measure of how turbulent the flow can be in the absence of an adequate control action. Roughly speaking, the larger the Reynolds number is, the more turbulent the flow may be.  Moreover, the Navier-Stokes equations involve a nonlinear term, which we refer to as convection. This nonlinearity is the source of turbulence. When the Reynolds number is high, convection dominates over diffusion and the fluid becomes turbulent. When the Reynolds number is small, the dominant term in the Navier-Stokes equations is diffusion, and the fluid, rather than being turbulent, is laminar. Many attempts have been made in the fluid flow control community to derive boundary feedback laws for the regulation of Navier-Stokes equations towards regular equilibrium profiles, namely, parabolic Poiseuille profiles. The literature on this topic is vast and it is not possible to survey all the obtained results here. To cite just a few references, regulation of the Navier-Stokes equations towards parabolic equilibrium profiles, but for low Reynolds numbers, has been tackled e.g. in [3, 4]. These results aim at maintaining stability or at stability enhancement, as for low Reynolds numbers, parabolic equilibrium profiles are stable in open-loop, i.e. with zero control action. On the other hand, regulation of the linearized Navier-Stokes equations, i.e. without the nonlinear convection, towards parabolic equilibrium profiles, but for arbitrary Reynolds numbers, is considered e.g. in [5, 6, 7, 8, 9]. Other works are based on reduced-order models of the linearized Navier-Stokes equations [10, 11]. The disadvantage of controllers designed based on reduced-order models is that the stability guarantees obtained on the reduced-order models do not necessarily extend to the original PDEs. Optimal control theory has also been a leading approach for the control of the Navier-Stokes equations [12, 13]. However, with optimal control, the feedback is usually not obtained in closed form, and one needs either to linearize the equations or to approach numerically the optimal controller via computationally expensive approaches. We refer e.g. to [14, 15, 16] for reviews and further details.

Obtaining a closed-form feedback law for the nonlinear Navier-Stokes equations that leads to global stability of parabolic equilibrium profiles independently of the size of the Reynolds number is one of the major open problems in fluid flow control. In this paper, we reveal that the primary source of turbulence, nonlinear convection, is the key to tackling this problem. More precisely, we consider a two-dimensional channel of an incompressible fluid, with periodic boundary conditions in the streamwise direction, no-slip boundary conditions, i.e. zero tangential velocity at the top and bottom walls of the channel, and finally normal velocity actuation at the top and bottom walls. We show that, due to nonlinear convection, the derivative of the $L^2$ energy of the regulation error has a cubic dependence in the boundary control input. Then, we exploit this cubic dependence in the control design to enforce global exponential stability of the considered parabolic Poiseuille equilibrium profile in the $L^2$ sense. Moreover, we show that the mean of the control inputs is zero, implying that the net fluid flux through the top and bottom walls of the channel is null. The cubic mean of the control inputs is however not null. As we shall see later on, this is the key property to exploit convection for stabilization.

In addition to the importance of this result in fluid mechanics, we provide here the first generalization of universal control theory to a vector-valued PDE: the Navier-Stokes equations. To understand this last statement, one first needs to recall the seminal works by Artstein [17] and Sontag [18]. In [17], Artstein considers the problem of smooth stabilization of finite-dimensional control systems. One of Artstein's results is that if the control system admits a smooth control Lyapunov function (CLF), then there exists an asymptotically stabilizing feedback law, that is smooth except eventually at the origin. If the system verifies the small control property, then the feedback law is continuous at the origin. Artstein's result is not constructive. The feedback law is constructed by Sontag in [18] for the specific class of control-affine systems, and is known nowadays as Sontag's universal formula. The formula is said to be universal because its structure is independent of the vector fields defining the system. These vector fields are just arguments of the formula. For a long time, this theory of universal control, with explicit feedback law, was specific to control-affine systems. The generalization of Sontag's formula to some classes of control systems admitting a CLF whose derivative is not affine in the control input has been made in [19]. In the latter reference, we have shown that the class of scalar parabolic PDEs with quadratic convection admits, because of convection, the $L^2$ energy as a CLF whose derivative has a cubic dependence in the Dirichlet control input. We constructed a feedback law, the Cardano-Lyapunov formula, that exploits this cubic dependence to perform global stabilization. As the name suggests, scalar parabolic PDEs with quadratic convection can be stabilized by solving cubic equations, via the Cardano root formula. In [20], we have shown that the Cardano-Lyapunov formula is a basis for the construction of a class of inverse optimal stabilizing controllers, that are the unique asymptotically stabilizing minimizers of cost functionals that, rather than being as usual quadratic in the control input, are polynomial. This was the first generalization of inverse optimal control theory to systems admitting a CLF whose derivative is not affine in the control input. The fundamental idea that convection may be helpful for stabilization is not new in the control community. It first appears in [21], for the problem of regulation of constant equilibrium of the scalar-valued Burgers' equation. It is observed in the latter reference that, due to convection, the derivative of the $L^2$ energy of the difference between the state and the equilibrium has a cubic dependence in the Dirichlet control input. This cubic dependence is then exploited to achieve global stabilization. However, since this first paper, very few works have been done where nonlinear convection is considered a helpful term. Some of these works are [22, 23, 24], where quadratic convection is exploited to achieve $L^2$ stabilization via Dirichlet boundary control of the origin of the scalar-valued Kuramoto-Sivashinsky equation. When the PDE is not convective, one needs to rely on diffusion combined with a weighted $L^2$ energy as a CLF candidate [25]. In this paper, we show that the convection in the Navier-Stokes equations is the key to performing global stabilization. The generalization of the Cardano-Lyapunov formula to this vector-valued PDE is however not straightforward. The main difficulty is to handle the divergence-free condition on the state, which accounts for the conservation of mass in the system. The divergence-free condition has the advantage of decoupling the equation of the tangential velocity from the equation of the normal velocity in the $L^2$ energy estimates. However, the price of this decoupling is a strong constraint that the Dirichlet control inputs need to verify. Namely, at each time instant, the mean value of the difference between the control input at the top wall and the control input at the bottom wall needs to be zero. To handle this constraint, we select the control at the top wall and the control at the bottom wall such that both of them have zero mean values. The Cardano-Lyapunov formula developed in [19] has the property of being of constant sign, i.e. either nonpositive or nonnegative, depending on the location of the Dirichlet control input. The only constant-sign function of zero mean is the zero function. As a consequence, the control at the top wall and the control at the bottom wall will not be given by the original Cardano-Lyapunov formula, but by a modification of this formula, that has zero mean value, but non-zero cubic mean value. The fact that the control inputs have a non-zero cubic mean value is, as we reveal in this paper, primordial to exploit the stabilizing effect of convection.

Notation. The divergence of a vector-valued function $\textbf{W}(x,y,t) = (U(x,y,t),V(x,y,t))^{\top}$ with respect to $(x,y)$ is defined by $\nabla \cdot \textbf{W} := U_x + V_y$. We say that $\textbf{W}$ is divergence-free if $\nabla \cdot \textbf{W} = 0$. The Laplacian of $\textbf{W}$ is $\Delta \textbf{W} := U_{xx}+V_{yy}$. The gradient of a scalar field $P(x,y,t)$ with respect to $(x,y)$ is the vector $\nabla P := (P_x,P_y)^{\top}$. Finally, we write $(\textbf{W}\cdot \nabla)\textbf{W} := (UU_x+VU_y,UV_x+VV_y)^{\top}$.

§ PROBLEM STATEMENT AND MAIN RESULT

§.§ Problem statement

We consider hydrodynamic channel flows described by the incompressible Navier-Stokes equations

\begin{equation}\label{NSE}
\left\lbrace
\begin{aligned}
\textbf{W}_t  - \frac{1}{R} \Delta \textbf{W}+(\textbf{W}\cdot \nabla)\textbf{W}+\nabla P &= 0, \\
\nabla \cdot \textbf{W} &= 0,
\end{aligned}
\right.
\end{equation}

where $\textbf{W}(x,y,t) = (U(x,y,t),V(x,y,t))^\top \in \mathbb{R}^2$ is the fluid's velocity field, $P(x,y,t)\in \mathbb{R}$ is the fluid's pressure field, $R>0$ is the Reynolds number, and $(x,y,t) \in [0,L_x]\times [0,L_y]\times (0,+\infty)$, with $L_x, L_y>0$. We refer to $U$ as the tangential velocity of the fluid, and to $V$ as the normal velocity. We assume that (<ref>) is subject to the no-slip boundary conditions

\begin{align}
U(x,0,t) = U(x,L_y,t) = 0.\label{bc1_NSE}
\end{align}

Moreover, we consider periodic boundary conditions in the $x$ direction, namely,

\begin{align}
U(0,y,t) &= U(L_x,y,t), \label{perio1}\\
V(0,y,t) &= V(L_x,y,t), \label{perio2} \\
V_x(0,y,t) &= V_x(L_x,y,t). \label{perio3}
\end{align}

The pressure is assumed to satisfy the boundary condition

\begin{align}
P(0,y,t) = P(L_x,y,t)+aL_x, \label{pressure_b}
\end{align}

where $a\geq 0$ is constant. Finally, we consider normal velocity actuation

\begin{align}
V(x,0,t) &= -V(x,L_y,t) =F(x)\psi (t), \label{control1}
\end{align}

where $(F,\psi)$ are functions to be designed.

The considered boundary conditions are standard in fluid mechanics. The no-slip boundary conditions (<ref>) model the fact that fluid particles, close to a solid boundary, do not move along with the flow since the adhesion force between the fluid particles and the solid boundary particles is stronger than the cohesion force between the fluid particles. The no-slip boundary conditions are used in the vast majority of fluid engineering applications to describe the behavior of the fluid at solid boundaries. The boundary conditions (<ref>)-(<ref>) are frequently employed to approximate large channels by periodic structures. A large channel is then seen as a concatenation of the smaller channels of length $L_x$. The boundary condition (<ref>) models the case where pressure is periodic up to a constant, here $aL_x$. There is a pressure difference, and therefore, a flow through the channel. Finally, the boundary conditions in (<ref>) are the actuation at the top and bottom walls of the channel.
This type of actuation, which is distributed along the
channel's walls, can be implemented using
micro-electro-mechanical systems (MEMS) [16].

When $\psi :=0$, the system (<ref>)-(<ref>) admits as an equilibrium the following parabolic Poiseuille profile

\begin{align}
\bar{U}(y) &:= \frac{R}{2}ay(L_y-y), \label{equ}\\
\bar{V} &:=0, \label{eqv}\\
\bar{P}(x) &:=-ax+b, \label{eqp}
\end{align}

where $b\geq0$. The steady-state velocity field $(\bar{U},\bar{V})$ is schematically represented in Figure <ref>. The objective is to design $(F,\psi)$ to achieve $L^2$ global stabilization of $(\bar{U},\bar{V})$, for any value of $R$. Strictly speaking, for a desired decay rate $\alpha >0$, we aim at achieving

\begin{align}
E(t) \leq E(0)e^{-\alpha t} \quad t\geq 0, \label{exp_stab}
\end{align}

\begin{align}
E &:= \frac{1}{2}\int_{0}^{L_x} \int_{0}^{L_y}\big[u(x,y)^2+v(x,y)^2 \big]\ dxdy, \label{lyapu_2} \\
u&:= U-\bar{U}, \label{u} \\
v &:= V-\bar{V}. \label{v}
\end{align}

[gray!20] (0,0) rectangle (,);

(0,0) rectangle (,);

[thick, black, domain=0:, samples=100] plot (-)*/+ 0.5*, );

in 0.2,0.4,0.6,...,
[black,->] (0.5*,) – (-)*/+ , );

[->] (0,-0.5) – (0,0.5) node[left] $y$;
[->] (-0.5,0) – (0.5,0) node[below] $x$;

[thick, dashed, black] (1.99, 0) – (1.99, );

at (,-0.4) $L_x$;
at (-0.4,) $L_y$;
[black] at (-0.8,/2) $(\bar{U}(y),\bar{V})$;

The velocity equilibrium profile (<ref>)-(<ref>) to be regulated. The actuation is the normal velocity at the bottom and top walls.

Before we state our main result, we shall recall a constraint we have when designing $(F,\psi)$. This constraint comes from the conservation of mass, i.e., the divergence-free condition $\nabla \cdot \textbf{W} = 0$ in (<ref>). The condition $\nabla \cdot \textbf{W} = 0$ means that

\begin{align}
U_x + V_y = 0. \label{div-free}
\end{align}

By integrating both sides of (<ref>) with respect to $x$ and $y$, we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y} \big[U_x+V_y\big] \ dxdy = 0.\label{div-free2}
\end{align}

Using integration by parts, we can rewrite (<ref>) as

\begin{align}
\int_{0}^{L_x} &\big[V(x,L_y)-V(x,0)\big]\ dx \nonumber \\
&~= - \int_{0}^{L_y} \big[U(L_x,y)-U(0,y)\big]\ dy. \label{div-free3}
\end{align}

Since $U(L_x,y)=U(0,y)$, and $V(x,0,t)=-V(x,L_y,t)=F(x)\psi (t)$, then (<ref>) can be rewritten as

\begin{align}
\psi (t)\int_{0}^{L_x}F(x) \ dx = -\psi (t)\int_{0}^{L_x} F(x)\ dx \quad \forall t>0.\label{constraint}
\end{align}

To satisfy (<ref>), we will design $F$ such that it has zero mean. However, as we shall see later on, $F$ needs also to have a non-zero cubic mean, as this will be fundamental for stabilization.

§.§ Main result

We prove in this paper the following result.

Consider the system (<ref>) subject to the boundary conditions (<ref>)-(<ref>). Let $(\bar{U},\bar{V},\bar{P})$ be the parabolic Poiseuille equilibrium profile defined in (<ref>)-(<ref>), and $E$ be the $L^2$ energy of the regulation error defined in (<ref>), with $(u,v)$ the deviation variables in (<ref>)-(<ref>). The equilibrium $(\bar{U},\bar{V})$ is globally exponentially stable in the $L^2$ sense at a prescribed decay rate $\alpha >0$, i.e. for any initial condition $(U(x,y,0),V(x,y,0),P(x,y,0))$ we have $E(t)\leq E(0)e^{-\alpha t}$ for all $t\geq 0$, if

* $F$ is any function such that $F(0)=F(L_x)$, $F'(0)=F'(L_x)$,

\begin{align}
\int_{0}^{L_x}F(x)\ dx = 0, \ \ \text{and} \ \ \int_{0}^{L_x}F(x)^3\ dx \neq 0.
\end{align}

* $\psi$ is given by

\begin{align}
\psi := -\frac{\Psi}{\sqrt[3]{\int_{0}^{L_x}F(x)^3\ dx}}, \label{pPsi}
\end{align}

\begin{align}
\Psi := \sqrt[3]{-\frac{q}{2} + \sqrt{\frac{q^2}{4}+\frac{\beta^3}{27}}}+\sqrt[3]{-\frac{q}{2} - \sqrt{\frac{q^2}{4}+\frac{\beta^3}{27}}}, \label{cardano_thm}
\end{align}

\begin{align}
\beta :=&~ \frac{\int_{0}^{L_x}F(x)\big[p(x,0)+p(x,L_y)\big]\ dx}{\sqrt[3]{\int_{0}^{L_x}F(x)^3\ dx}}, \label{beta_thm} \\
p :=&~ P - \bar{P}, \label{p} \\
q :=&~ \Gamma + \frac{2\sqrt{3}}{9}\big|\beta\big|^{\frac{3}{2}}, \label{q_thm}\\
\Gamma :=&~ \alpha E + \bigg|\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}'v\ dxdy\bigg|. \label{g_thm}
\end{align}

To implement the control law in Theorem <ref>, one needs to select a periodic function $F$ with zero mean but non-zero cubic mean. An example of such a function is provided in the next proposition.

Let $L_x\in \mathbb{R}_{>o}$, $\epsilon \in (0,L_x)$ with $\epsilon \neq L_x/2$ and $\epsilon \neq L_x/4$, and $\theta \in \mathbb{R}^*$. Let $F:[0,L_x]\to \mathbb{R}$ be defined as follows

\begin{equation}
F(x) :=
\left\lbrace
\begin{aligned}
&\frac{\theta (2\epsilon - L_x)}{2\epsilon} \quad &&x\in [0,\epsilon)\cup [L_x-\epsilon,L_x], \\
&\ \ \ \ \theta \quad &&x\in [\epsilon,L_x-\epsilon).
\end{aligned}
\right. \label{f_def}
\end{equation}

Then $F(0)=F(L_x)$, $F'(0)=F'(L_x)$, $\int_{0}^{L_x}F(x)\ dx=0$, and

\begin{align}
\int_{0}^{L_x}F(x)^3\ dx = \theta^3 \bigg[(L_x-2\epsilon) + \frac{(2\epsilon - L_x)^3}{4\epsilon^2}\bigg] \neq 0.  \label{cubic_mean_value}
\end{align}

A schematic representation of the boundary control inputs designed in Theorem <ref> when $F$ is constructed according to Proposition <ref>.

The function $\Psi$ defined by (<ref>) is nonpositive. As a consequence, according to (<ref>), when $\psi \neq 0$, $\psi$ has the same sign as the cubic mean of $F$. A schematic representation of the control inputs when $F$ is designed according to Proposition <ref> is in Figure <ref>.

It is interesting to compare our control strategy, which is based on a Lyapunov-type analysis of the control system (<ref>)-(<ref>), with a popular control strategy in the fluid flow control community, that consists of periodically blowing and sucking fluid at the top and bottom walls of the channel [26, 27]. The latter strategy, which does not require feedback, has been shown numerically to be successful for turbulence attenuation in some cases. The control takes the form of a sinusoidal function $V(x,0,t)=-V(x,L_y,t) := -2A \cos(\omega (x -ct))$, where $A$, $\omega$ and $c$ are tuned experimentally. There are two fundamental differences between our controller and the sinusoidal control input previously mentioned. First, our controller incorporates feedback. It is because of this feedback that we can guarantee global stability results for a chosen parabolic Poiseuille equilibrium profile and any Reynolds number. Second, although our controller has a zero mean, its cubic mean is non-null. The sinusoidal function has both a zero mean and a zero cubic mean. According to the $L^2$ energy estimates in Section <ref>, the convection has a stabilizing effect on the system when the boundary inputs have non-zero cubic mean.

§ PROOF OF THE MAIN RESULT

The first step of the proof consists of writing the system of PDEs describing the evolution of the deviation variables $(u,v,p)$. It is given by

\begin{equation}
\begin{aligned}
&u_t-\frac{1}{R}\Delta u + uu_x + vu_y + \bar{U}u_x + \bar{U}'v+p_x = 0, \\
&v_t-\frac{1}{R}\Delta v + uv_x + vv_y + \bar{U}v_x +p_y = 0, \\
&u_x+v_y = 0. \label{Error_NS}
\end{aligned}
\end{equation}

Moreover, it is subject to the following set of boundary conditions

\begin{align}
u(x,0,t) &= u(x,L_y,t) = 0, \label{err_bo} \\
u(0,y,t) &= u(L_x,y,t), \label{err_b1} \\
v(0,y,t) &= v(L_x,y,t), \label{err_b2} \\
v_x(0,y,t) &= v_x(L_x,y,t), \label{err_b3} \\
p(0,y,t) &= p(L_x,y,t), \label{err_b4} \\
v(x,0,t) &= -v(x,L_y,t) = F(x)\psi(t). \label{err_b6}
\end{align}

Next, we differentiate $E$ with respect to time, and we derive an upper bound on $\dot{E}$. By differentiating $E$ along the strong solutions to (<ref>)-(<ref>), we find

\begin{align}
\dot{E} =&~ \int_{0}^{L_x}\int_{0}^{L_y}\big[uu_t + vv_t\big]\ dxdy \nonumber \\
=&~\int_{0}^{L_x}\int_{0}^{L_y}u\bigg[\frac{1}{R}\Delta u - uu_x - vu_y - \bar{U}u_x- \bar{U}'v \nonumber \\
&~-p_x\bigg]\ dxdy+\int_{0}^{L_x}\int_{0}^{L_y}v\bigg[\frac{1}{R}\Delta v - uv_x \nonumber \\
&~- vv_y - \bar{U}v_x - p_y\bigg]\ dxdy. \label{foo}
\end{align}

Using integration by parts, we obtain

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y} u\Delta u dxdy =&~  \int_{0}^{L_x}\int_{0}^{L_y} u\big[u_{xx}+u_{yy}\big] dxdy \nonumber \\
=&~- \int_{0}^{L_x}\int_{0}^{L_y} \big[u_x^2+u_y^2\big] dxdy \nonumber \\
&~ +\int_{0}^{L_y}\big[uu_x]_{x=0}^{x=L_x}\ dy\nonumber \\
&~+ \int_{0}^{L_x}\big[uu_y]_{y=0}^{y=L_y}\ dx. \label{laplace1}
\end{align}

According to the no-slip boundary conditions $u(x,0)=u(x,L_y)=0$, we have

\begin{align}
\int_{0}^{L_x} \big[uu_y\big]_{y=0}^{y=L_y}\ dx = 0. \label{to_combine1}
\end{align}

Moreover, since $(u,v)$ is divergence-free, then

\begin{align}
u_x(0,y) &= -v_y(0,y), \label{proof_n1}\\
u_x(L_x,y) &= -v_y(L_x,y). \label{proof_n2}
\end{align}

By differentiating both sides of the periodic boundary condition (<ref>) with respect to $y$, we find $v_y(0,y)=v_y(L_x,y)$, which implies, according to (<ref>) and (<ref>), that $u_x(0,y)=u_x(L_x,y)$. Using the periodicity of both $u_x$ and $u$ in the $x$ direction, we conclude that

\begin{align}
\int_{0}^{L_y}\big[uu_x\big]_{x=0}^{x=L_x}\ dy = 0. \label{to_combine2}
\end{align}

Using (<ref>) and (<ref>), we can rewrite (<ref>) as

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y} u\Delta u \ dxdy = -\int_{0}^{L_x}\int_{0}^{L_y} \big[u_x^2+u_y^2\big]\ dxdy. \label{important1}
\end{align}

Similarly, using integration by parts, we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y} v\Delta v dxdy =&~  \int_{0}^{L_x}\int_{0}^{L_y} v\big[v_{xx}+v_{yy}\big] dxdy \nonumber \\
=&~- \int_{0}^{L_x}\int_{0}^{L_y} \big[v_x^2+v_y^2\big] dxdy \nonumber \\
&~ +\int_{0}^{L_y}\big[vv_x]_{x=0}^{x=L_x}\ dy\nonumber \\
&~+ \int_{0}^{L_x}\big[vv_y]_{y=0}^{y=L_y}\ dx. \label{laplace2}
\end{align}

Since $v$ and $v_x$ are periodic in the $x$ direction, then

\begin{align}
\int_{0}^{L_y}\big[vv_x\big]_{x=0}^{x=L_x}\ dy = 0. \label{to_combine3}
\end{align}

Moreover, using the divergence-free condition on $(u,v)$, we conclude that

\begin{align}
v_y(x,0) &= - u_x(x,0), \label{perio_new1}\\
v_y(x,L_y) &= -u_x(x,L_y). \label{perio_new2}
\end{align}

Since $u(x,0)=u(x,L_y)=0$, then (<ref>) and (<ref>) imply that $v_y(x,0)=v_y(x,L_y)=0$. As a consequence, we have

\begin{align}
\int_{0}^{L_x}\big[vv_y\big]_{y=0}^{y=L_y}\ dx = 0. \label{to_combine4}
\end{align}

Using (<ref>) and (<ref>), we can rewrite (<ref>) as

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}v\Delta v\ dxdy = -\int_{0}^{L_x}\int_{0}^{L_y}\big[v_x^2+v_y^2\big]\ dxdy. \label{important2}
\end{align}

Next, using integration by parts, note that we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}up_x\ dxdy =&~ \int_{0}^{L_y}\big[up\big]_{x=0}^{x=L_x}\ dy\nonumber \\
&~- \int_{0}^{L_x}\int_{0}^{L_y}u_xp\ dxdy.
\end{align}

Since $u$ and $p$ are periodic in the $x$ direction, then

\begin{align}
\int_{0}^{L_y}\big[up\big]_{x=0}^{x=L_x}\ dy = 0.
\end{align}

As a consequence, we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}up_x\ dxdy = -\int_{0}^{L_x}\int_{0}^{L_y}u_xp \ dxdy. \label{pressure_add1}
\end{align}

On the other hand, using integration by parts, we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}vp_y\ dxdy =&~ \int_{0}^{L_x}\big[vp\big]_{y=0}^{y=L_y}\ dx \nonumber \\
&~-\int_{0}^{L_x}\int_{0}^{L_y}v_yp\ dxdy. \label{deltap}
\end{align}

Since $v(x,0,t)=-v(x,L_y,t)=F(x)\psi(t)$, then we can rewrite (<ref>) as

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}vp_y\ dxdy =&~ -\psi\int_{0}^{L_x}F(x)p(x,L_y)\ dx \nonumber \\
&~-\psi \int_{0}^{L_x}F(x)p(x,0)\ dx \nonumber \\
&~-\int_{0}^{L_x}\int_{0}^{L_y}v_vp\ dxdy. \label{pressure_add2}
\end{align}

By adding (<ref>) and (<ref>), and using the fact that $u_x+v_y=0$, we conclude that

\begin{align}
\int_{0}^{L_x}&\int_{0}^{L_y}\big[up_x+vp_y\big]dxdy \nonumber \\
=&~-\psi \int_{0}^{L_x}F(x)\big[p(x,L_y)+p(x,0)\big]dx. \label{important3}
\end{align}

Now, we analyze the effect of the terms that come from the convection $(\textbf{W}\cdot\nabla)\textbf{W}$. Using integration by parts, note that we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}u^2u_x\ dxdy = \int_{0}^{L_y}\bigg[ \frac{u^3}{3}\bigg]_{x=0}^{x=L_x}\ dy. \label{eq8}
\end{align}

Since $u$ is periodic in the $x$ direction, then (<ref>) becomes

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}u^2u_x\ dxdy = 0. \label{eq9}
\end{align}

As a consequence, we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}u\big[uu_x+vu_y\big]dxdy = \int_{0}^{L_x}\int_{0}^{L_y}uvu_ydxdy. \label{eq10}
\end{align}

On the other hand, since $(u,v)$ is divergence-free, then

\begin{align}
u\big[(u^2)_x+(uv)_y\big] =&~ 2u^2u_x + u^2 v_y + uvu_y \nonumber \\
=&~u^2\big[u_x+v_y\big] + u\big[uu_x+vu_y\big] \nonumber \\
=&~u\big[uu_x+vu_y\big]. \label{eq11}
\end{align}

Using the identity (<ref>), we obtain

\begin{align}
&\int_{0}^{L_x}\int_{0}^{L_y}u[uu_x+vu_y] \ dxdy = \int_{0}^{L_x}\int_{0}^{L_y}u\big[(u^2)_x \nonumber \\
&~+(uv)_y\big]dxdy= \int_{0}^{L_y}\bigg[\frac{2u^3}{3}\bigg]_{x=0}^{x=L_x} dy\nonumber \\
&~+\int_{0}^{L_x}[u^2v]_{y=0}^{y=L_y} dx - \int_{0}^{L_x}\int_{0}^{L_y} uvu_y\ dxdy.\label{eq12}
\end{align}

Using the periodicity of $u$ in the $x$ direction, and the no-slip boundary conditions $u(x,0)=u(x,L_y)=0$, we can rewrite (<ref>) as

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}u\big[uu_x+vu_y\big]dxdy = - \int_{0}^{L_x}\int_{0}^{L_y}uvu_ydxdy. \label{convect2}
\end{align}

Combining (<ref>) and (<ref>), note that we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}&u\big[uu_x+vu_y\big]\ dxdy= \nonumber \\
&~-\int_{0}^{L_x}\int_{0}^{L_y}u\big[uu_x+vu_y\big]\ dxdy,
\end{align}

which implies that

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}u\big[uu_x+vu_y\big]\ dxdy = 0. \label{important4}
\end{align}

Next, using integration by parts, we obtain

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}v^2v_y \ dxdy =&~ \int_{0}^{L_x} \bigg[\frac{v^3}{3}\bigg]_{y=0}^{y=L_y} \ dx \nonumber \\
=&~ -\frac{2}{3}\psi^3\int_{0}^{L_x}F(x)^3\ dx.\label{eq17}
\end{align}

As a consequence, note that we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y} v\big[uv_x+&vv_y\big]dxdy = -\frac{2}{3}\psi^3\int_{0}^{L_x}F(x)^3\ dx \nonumber \\
&~+\int_{0}^{L_x}\int_{0}^{L_y}vuv_x dxdy. \label{eq19}
\end{align}

On the other hand, since $(u,v)$ is divergence-free, then

\begin{align}
v\big[(v^2)_y + (uv)_x \big]=&~ 2v^2v_y + v^2 u_x + vuv_x \nonumber \\
=&~v^2\big[u_x+v_y\big] + v\big[vv_y+uv_x\big] \nonumber \\
=&~v\big[vv_y+uv_x\big]. \label{eq20}
\end{align}

Using the identity (<ref>), we obtain

\begin{align}
&\int_{0}^{L_x}\int_{0}^{L_y} v\big[uv_x+vv_y\big]dxdy = \int_{0}^{L_x}\int_{0}^{L_y} v\big[(v^2)_y \nonumber \\
&~+ (uv)_x\big]dxdy = \int_{0}^{L_x}\left[\frac{2v^3}{3}\right]_{y=0}^{y=L_y}dx \nonumber \\
&~+ \int_{0}^{L_y}\big[v^2u\big]_{x=0}^{x=L_x}dy-\int_{0}^{L_x}\int_{0}^{L_y}v_xuv\ dxdy. \label{eq21}
\end{align}

Since $u$ and $v$ are periodic in the $x$ direction, then we can rewrite (<ref>) as

\begin{align}
\int_{0}^{L_x}&\int_{0}^{L_y} v\big[uv_x+vv_y\big]dxdy \nonumber \\
=&~ -\frac{4}{3}\psi^3\int_{0}^{L_x}F(x)^3\ dx-\int_{0}^{L_x}\int_{0}^{L_y}v_xuvdxdy. \label{eq19bis}
\end{align}

Combining (<ref>) and (<ref>), note that we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}v_xuv\ dxdy = -\frac{1}{3}\psi^3\int_{0}^{L_x}F(x)^3\ dx.
\end{align}

As a consequence, we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}v\big[uv_x+vv_y\big]\ dxdy = -\psi^3 \int_{0}^{L_x}F(x)^3\ dx. \label{important5}
\end{align}

Finally, using integration by parts, and the fact that $\bar{U}$ is independent of $x$, note that we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}u_x\ dxdy = \int_{0}^{L_y}\bar{U}\bigg[\frac{u^2}{2}\bigg]_{x=0}^{x=L_x}\ dy. \label{uubarux}
\end{align}

Since $u$ is periodic in the $x$ direction, then (<ref>) becomes

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}u_x\ dxdy = 0. \label{important6}
\end{align}

Similarly, using integration by parts, the fact that $\bar{U}$ is independent of $x$, and that $v$ is periodic in the $x$ direction, we have

\begin{align}
\int_{0}^{L_x}\int_{0}^{L_y} v\bar{U}v_x\ dxdy = 0. \label{important7}
\end{align}

Using (<ref>), (<ref>), (<ref>), (<ref>), (<ref>), (<ref>), and (<ref>), note that we have

\begin{align}
\dot{E} \leq&~ -\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}'v\ dxdy \nonumber \\
&~-\bigg[\int_{0}^{L_x}F(x)\big[p(x,0)+p(x,L_y)\big]\ dx\bigg] \psi \nonumber \\
&~-\bigg[\int_{0}^{L_x}F(x)^3\ dx\bigg]\psi^3. \label{edot_ineq2}
\end{align}

Next, we perform the change of variable

\begin{align}
\Psi (t) := -\left[\sqrt[3]{\int_{0}^{L_x}F(x)^3\ dx}\right]\psi(t).
\end{align}

We can rewrite (<ref>) as

\begin{align}
\dot{E}\leq -\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}'v\ dxdy + \beta(p)\Psi + \Psi^3,
\end{align}

where $\beta(q)$ is defined in (<ref>).

Consider now the cubic equation

\begin{align}
\Psi^3+\beta(p)\Psi + q(u,v,p) = 0,\label{cubic_eq_proof}
\end{align}

where $q(u,v,p)$ is defined in (<ref>)-(<ref>). If $\Psi$ is a real root of (<ref>), then

\begin{align}
\dot{E}\leq&~ -\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}'v\ dxdy - q(u,v,p) \nonumber \\
\leq &~ - \int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}'vdxdy - \bigg|\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}'vdxdy\bigg| \nonumber \\
&~ - \frac{2\sqrt{3}}{9}\big|\beta(p)\big|^{\frac{3}{2}} - \alpha E \nonumber \\
\leq &~ - \alpha E,
\end{align}

which would allow us to conclude on global exponential stability at the decay rate $\alpha$. The discriminant of (<ref>) is

\begin{align}
\frac{q^2}{4}+\frac{\beta^3}{27} = \frac{1}{4}\bigg[\Gamma^2+\frac{4\sqrt{3}}{9}\big|\beta\big|^{\frac{3}{2}}\Gamma \bigg] + \frac{1}{27}\bigg[\big|\beta|^3-\beta^3\bigg],
\end{align}

where $\Gamma$ is defined in (<ref>). Note that $\Gamma \geq E$. As a consequence, $q^2/4+\beta^3/27 >0$ as long as $E\neq 0$, i.e. as long as $(u,v)\neq 0$. Therefore, when $(u,v)\neq 0$, the cubic equation (<ref>) admits a unique real root which is given by (<ref>). Now, if $(u,v)=0$, then $q=\beta = 0$. Indeed, when $(u,v)=0$, we have, according to (<ref>), $p_x=0$ and $p_y=0$. It implies that $p$ is constant in $(x,y)$. As a result, note that we have

\begin{align}
\int_{0}^{L_x}F(x)&\big[p(x,0)+p(x,L_y)\big]\ dx = \nonumber \\
&~\big[p(x,0)+p(x,L_y)\big]\int_{0}^{L_x}F(x)\ dx. \label{beta_zero}
\end{align}

Since $F$ has zero mean, then we conclude from (<ref>) that

\begin{align}
\int_{0}^{L_x}F(x)&\big[p(x,0)+p(x,L_y)\big]\ dx = 0,
\end{align}

which implies that $\beta = 0$. On the other hand, $E=0$, and $\int_{0}^{L_x}\int_{0}^{L_y}u\bar{U}'v\ dxdy = 0$, which implies that $q=0$. Therefore, when $(u,v)=0$, the cubic equation (<ref>) admits the unique root $\Psi=0$, which is given by (<ref>).

§ CONCLUSION

We constructed in this paper the first boundary feedback controller that achieves global exponential regulation of the two-dimensional Navier-Stokes equations toward parabolic Poiseuille equilibrium profiles for any Reynolds number. We have shown, through $L^2$ energy estimates, that the source of turbulence, namely nonlinear convection, is the key to stabilizing turbulent flows when the boundary control inputs have zero mean but a non-zero cubic mean. Future research includes not only numerical tests of the proposed controller, but also the generalization of our approach to the Navier-Stokes equations in three-dimensional channels, and to other turbulent flows such as the Hartman flow in magnetohydrodynamic channels [9].

Loria Nešić Loria Nešić

[1]
F. Mourad, E. Witrant, and F. Pattyn, “The Initialization Of Basal Sliding Coefficients For Antarctica A Lyapunov Based Approach,” in Proceedings of the American Control Conference, pp. 3014–3019, 2018.

[2]
F. Mourad, Distributed parameter inversion of basal sliding and diffusion of the Antarctic ice sheet.
PhD thesis, Université Grenoble Alpes, 2020.

[3]
A. Balogh, W.-J. Liu, and M. Krstić, “Stability enhancement by boundary control in 2-D channel flow,” IEEE Transactions on Automatic Control, vol. 46, no. 11, pp. 1696–1711, 2001.

[4]
O. M. Aamo, M. Krstić, and T. R. Bewley, “Control of mixing by boundary feedback in 2D channel flow,” Automatica, vol. 39, no. 9, pp. 1597–1606, 2003.

[5]
R. Vázquez Valenzuela, E. Trélat, and J. M. Coron, “Control for fast and stable laminar-to-high-Reynolds-numbers transfer in a 2D Navier-Stokes channel flow,” Discrete and Continuous Dynamical Systems. Series B, 10 (4), 925-956., 2008.

[6]
R. Vázquez and M. Krstić, “A Closed-Form Feedback Controller for Stabilization of the Linearized 2-D Navier–Stokes Poiseuille System,” IEEE Transactions on Automatic Control, vol. 52, no. 12, pp. 2298–2312, 2007.

[7]
J. Cochran, R. Vázquez, and M. Krstić, “Backstepping boundary control of Navier-Stokes channel flow: a 3D extension,” in Proceedings of the American Control Conference, pp. 6 pp.–, 2006.

[8]
J. Cochran, R. Vázquez, and M. Krstić, “Backstepping Boundary Control of Navier-Stokes Channel Flow: Explicit Gain Formulae in 3D,” in Proceedings of the Conference on Decision and Control, pp. 5329–5334, 2006.

[9]
R. Vázquez and M. Krstić, Control of turbulent and magnetohydrodynamic channel flows: boundary stabilization and state estimation, vol. 2.
Springer Science & Business Media, 2008.

[10]
S. S. Joshi, J. L. Speyer, and J. Kim, “A systems theory approach to the feedback stabilization of infinitesimal and finite-amplitude disturbances in plane Poiseuille flow,” Journal of Fluid Mechanics, vol. 332, pp. 157–184, 1997.

[11]
S. S. Joshi, J. L. Speyer, and J. Kim, “Finite dimensional optimal control of poiseuille flow,” Journal of Guidance, Control, and Dynamics, vol. 22, no. 2, pp. 340–348, 1999.

[12]
T. R. Bewley and S. Liu, “Optimal and robust control and estimation of linear paths to transition,” Journal of Fluid Mechanics, vol. 365, pp. 305–349, 1998.

[13]
T. R. Bewley, P. Moin, and R. Temam, “Dns-based predictive control of turbulence: an optimal benchmark for feedback algorithms,” Journal of Fluid Mechanics, vol. 447, pp. 179–225, 2001.

[14]
T. R. Bewley, “Flow control: new challenges for a new renaissance,” Progress in Aerospace sciences, vol. 37, no. 1, pp. 21–58, 2001.

[15]
M. Gad-el Hak, “Modern developments in flow control,” Appl. Mech. Rev., vol. 47, no. 7, pp. 365–379, 1996.

[16]
O. M. Aamo and T. I. Fossen, “Tutorial on Feedback Control of Flows, Part I: Stabilization of Fluid Flows in Channels and Pipes,” Modeling, Identification and Control, vol. 23, no. 3, pp. 161–226, 2002.

[17]
Z. Artstein, “Stabilization with relaxed controls,” Nonlinear Analysis: Theory, Methods & Applications, vol. 7, no. 11, pp. 1163–1173, 1983.

[18]
E. D. Sontag, “A ‘universal’construction of artstein's theorem on nonlinear stabilization,” Systems & control letters, vol. 13, no. 2, pp. 117–123, 1989.

[19]
M. C. Belhadjoudja, M. Krstić, M. Maghenem, and E. Witrant, “From Sontag’s to Cardano-Lyapunov Formula for Systems Not Affine in the Control: Convection-Enabled PDE Stabilization,” To appear in the Proceedings of the 2024 American Control Conference, Link to the file: https://drive.google.com/file/d/1mR7-QpWwn0DH3a75WIinkhK6VKw8VLnK/view?usp=sharing.

[20]
M. C. Belhadjoudja, M. Krstić, M. Maghenem, and E. Witrant, “Inverse Optimal Cardano-Lyapunov Feedback for PDEs with Convection,” To appear in the Proceedings of the 2024 European Control Conference, Link to the file: https://drive.google.com/file/d/1eH10Q6JqletvBGwqhJZTIjfHRS H4DzgI/view?usp=sharing.

[21]
M. Krstić, “On global stabilization of Burgers’ equation by boundary control,” Systems & Control Letters, vol. 37, pp. 123–141, 1999.

[22]
W.-J. Liu and M. Kristić, “Stability enhancement by boundary control in the Kuramoto–Sivashinsky equation,” Nonlinear Analysis: Theory, Methods & Applications, vol. 43, no. 4, pp. 485–507, 2001.

[23]
M. Maghenem, C. Prieur, and E. Witrant, “Boundary Control of the Kuramoto-Sivashinsky Equation Under Intermittent Data Availability,” in Proceedings of the American Control Conference, 2022.

[24]
M. C. Belhadjoudja, M. Maghenem, E. Witrant, and C. Prieur, “Adaptive stabilization of the Kuramoto-Sivashinsky equation subject to intermittent sensing,” in Proceedings of the American Control Conference, 2023.

[25]
S. Zekraoui, N. Espitia, and W. Perruquetti, “Lyapunov-based nonlinear boundary control design with predefined convergence for a class of 1d linear reaction-diffusion equations,” European Journal of Control, p. 100845, 2023.

[26]
B. K. Lieu, R. Moarref, and M. R. Jovanović, “Controlling the onset of turbulence by streamwise travelling waves. part 2. direct numerical simulation,” Journal of Fluid Mechanics, vol. 663, pp. 100–119, 2010.

[27]
R. Moarref and M. R. Jovanović, “Controlling the onset of turbulence by streamwise travelling waves. part 1. receptivity analysis,” Journal of fluid mechanics, vol. 663, pp. 70–99, 2010.