3.1 Nonlinear stability of an inviscid columnar rotating flow

We first summarize the results by Wang (Reference Wang2009), for the axisymmetric inviscid nonlinear stability of a general columnar rotating flow between two coaxial cylinders with axial periodic boundary conditions. A general base columnar rotating flow given by (3.1) was considered. Let $\mathscr{E}$ be the total kinetic energy integral of the flow in the domain $V$ , which is given by

(3.3) $$\begin{eqnarray}\displaystyle \mathscr{E}=\frac{1}{2}{\it\rho}\int _{V}(u_{r}^{2}+u_{x}^{2}+u_{{\it\theta}}^{2})\,\text{d}V, & & \displaystyle\end{eqnarray}$$

where $\text{d}V=r\,\text{d}r\,\text{d}x\,\text{d}{\it\theta}=\text{d}y\,\text{d}x\,\text{d}{\it\theta}$ . It was shown that, under the given boundary conditions and the application of integration by parts, $\mathscr{E}$ can also be written as

(3.4) $$\begin{eqnarray}\displaystyle \mathscr{E}=\mathscr{E}({\it\chi},K)=\frac{1}{2}{\it\rho}\int _{V}\left({\it\psi}{\it\chi}+\frac{K^{2}}{2y}\right)\text{d}V. & & \displaystyle\end{eqnarray}$$

Note that the boundary term vanishes in the integration because of the preservation of the volume flow rate. In addition, let $\mathscr{C}_{S}$ be the circulation-dependent integral (Casimir’s integral) given by

(3.5) $$\begin{eqnarray}\displaystyle \mathscr{C}_{S}=\mathscr{C}_{S}(K)={\it\rho}\int _{V}j(K)\,\text{d}V, & & \displaystyle\end{eqnarray}$$

where $j(K)$ is any sufficiently smooth function of $K$ . In the inviscid flow case, $\mathscr{E}({\it\chi},K)$ and $\mathscr{C}_{S}(K)$ are both conserved quantities in time during the perturbation’s evolution, i.e. $\text{d}\mathscr{E}/\text{d}t=0$ and $\text{d}\mathscr{C}_{S}/\text{d}t=0$ for all time $t$ .

Let the reduced Arnol’d function $\mathscr{A}_{rd}$ introduced in Wang (Reference Wang2009) be defined by

(3.6) $$\begin{eqnarray}\displaystyle \mathscr{A}_{rd}({\it\chi},K) & = & \displaystyle \mathscr{E}({\it\chi},K)+\mathscr{C}_{S}(K)\nonumber\\ \displaystyle & = & \displaystyle \frac{1}{2}{\it\rho}\int _{V}\left({\it\psi}{\it\chi}+\frac{K^{2}}{2y}\right)\text{d}V+{\it\rho}\int _{V}j(K)\,\text{d}V.\end{eqnarray}$$

The reduced Arnol’d function $\mathscr{A}_{rd}({\it\chi},K)$ is also conserved in time during the perturbation’s inviscid evolution, i.e. $\text{d}\mathscr{A}_{rd}/\text{d}t=0$ for all time $t$ .

In the formulation of $\mathscr{C}_{S}(K)$ , Wang (Reference Wang2009) chose the function $j(K)$ such that the first derivative of $j(K)$ with respect to $K$ (when evaluated at $K^{C}$ ) satisfies the condition

(3.7) $$\begin{eqnarray}\displaystyle \left.\frac{\text{d}}{\text{d}K}\right|_{K=K^{C}}j(K)=-\frac{K^{C}}{2Y(K^{C})}, & & \displaystyle\end{eqnarray}$$

where $K=K^{C}(y)$ is the circulation function of a base columnar rotating flow state. For a base flow with a monotonically increasing circulation function $K^{C}(y)$ , the function $Y(K^{C})$ is the inverse function of $K^{C}(y)$ , such that $Y(K^{C}(y))=y$ . Note that, for any given base flow, (3.7) is used to formulate the function $j(K)$ , and hence the formulation of $\mathscr{C}_{S}(K)$ .

The difference of the reduced Arnol’d function ${\rm\Delta}\mathscr{A}_{rd}({\it\chi},K):=\mathscr{A}_{rd}({\it\chi},K)-\mathscr{A}_{rd}({\it\chi}^{C},K^{C})$ is also conserved during the perturbation’s inviscid evolution, i.e. $\text{d}{\rm\Delta}\mathscr{A}_{rd}/\text{d}t=0$ for all time $t$ . Through the special choice of $j(K)$ as defined in (3.7), it was found that the first variation of the reduced Arnol’d function vanishes for any base columnar rotating flow. The difference ${\rm\Delta}\mathscr{A}_{rd}({\it\chi},K)$ is then determined by the application of the mean value theorem in the second-order Taylor’s expansion (see details in Wang Reference Wang2009), which is given by

(3.8) $$\begin{eqnarray}\displaystyle {\rm\Delta}\mathscr{A}_{rd}({\it\chi},K)=\tilde{\mathscr{E}}\rfloor _{{\it\xi}}=\frac{1}{2}{\it\rho}\int _{V}\left(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+\left[1-2y\left.\frac{\text{d}}{\text{d}K}\right|_{K=K\rfloor _{{\it\xi}}}\left(\frac{K}{2Y(K)}\right)\right]\tilde{u} _{{\it\theta}}^{2}\right)\text{d}V,\quad & & \displaystyle\end{eqnarray}$$

where $K\rfloor _{{\it\xi}}$ is an intermediate instantaneous function between the circulation function $K^{C}(y)$ of the base columnar rotating flow $\boldsymbol{u}^{C}(y)$ and the circulation function $K(x,y,t)$ of the instantaneous perturbed/total flow $\boldsymbol{u}(x,y,t)$ , that is, $K\rfloor _{{\it\xi}}=K^{C}(y)+{\it\xi}(t)[K(x,y,t)-K^{C}(y)]$ , for some $0<{\it\xi}(t)<1$ . The relationship (3.8) is significant because it involves only the velocity components of the perturbation. It implies that the weighted kinetic energy integral of the perturbation $\tilde{\mathscr{E}}\rfloor _{{\it\xi}}$ (on the right-hand side of (3.8)) is also conserved at all times during the flow evolution. This result leads to the following axisymmetric inviscid nonlinear stability criterion (proven in Wang Reference Wang2009):

Criterion 3.1. If a steady base columnar rotating flow satisfies Rayleigh’s (Reference Rayleigh1916) inviscid linear stability criterion in a strict sense, i.e. $\text{d}(K^{C}(y)^{2})/\text{d}y>0$ for all $y_{1}<y<y_{2}$ , then any finite-amplitude axisymmetric perturbation on the flow must obey the following energy estimate for all times $t\geqslant 0$ :

(3.9) $$\begin{eqnarray}\displaystyle \frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+{\it\lambda}^{-}\tilde{u} _{{\it\theta}}^{2})\,\text{d}V\leqslant \frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r0}^{2}+\tilde{u} _{x0}^{2}+{\it\lambda}^{+}\tilde{u} _{{\it\theta}0}^{2})\,\text{d}V, & & \displaystyle\end{eqnarray}$$

where

(3.10) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle {\it\lambda}^{-}=\min _{\substack{ y_{1}\leqslant y\leqslant y_{2} \\ K_{min}^{C}\leqslant K\rfloor _{{\it\xi}}\leqslant K_{max}^{C}}}\left[1-2y\left.\frac{\text{d}}{\text{d}K}\right|_{K=K\rfloor _{{\it\xi}}}\left(\frac{K}{2Y(K)}\right)\right],\\ \displaystyle {\it\lambda}^{+}=\max _{\substack{ y_{1}\leqslant y\leqslant y_{2} \\ K_{min}^{C}\leqslant K\rfloor _{{\it\xi}}\leqslant K_{max}^{C}}}\left[1-2y\left.\frac{\text{d}}{\text{d}K}\right|_{K=K\rfloor _{{\it\xi}}}\left(\frac{K}{2Y(K)}\right)\right],\end{array}\right\} & & \displaystyle\end{eqnarray}$$

and $\tilde{\boldsymbol{u}}_{0}=(\tilde{u} _{r0},\tilde{u} _{x0},\tilde{u} _{{\it\theta}0})$ is the initial velocity perturbation vector field. If, in addition, ${\it\lambda}^{-}>0$ , then the perturbation $\tilde{\boldsymbol{u}}=(\tilde{u} _{r},\tilde{u} _{x},\tilde{u} _{{\it\theta}})$ (in terms of the norm as defined by the left-hand side of (3.9)) is bounded from above by the initial perturbation $\tilde{\boldsymbol{u}}_{0}=(\tilde{u} _{r0},\tilde{u} _{x0},\tilde{u} _{{\it\theta}0})$ (in terms of the norm as defined by the right-hand side of (3.9)) for all times $t\geqslant 0$ , which means that the base columnar rotating flow $\boldsymbol{u}^{C}(y)$ is nonlinearly stable.

The inequality (3.9) shows the possible gain of the perturbation’s kinetic energy through the exchange of weighted kinetic energy between the perturbation’s velocity components. In the case of a solid-body rotation flow between the two coaxial cylinders with axial periodic boundary conditions, where ${\it\omega}_{1}={\it\omega}_{2}$ and then $B=0$ , we have $K^{C}=2Ay$ and then ${\it\lambda}^{-}={\it\lambda}^{+}=1$ . Then from (3.9), the perturbation’s kinetic energy at a later time $t$ equals the initial perturbation’s kinetic energy for all time $t$ . This means that, in the case of the solid-body rotation flow between two coaxial rotating cylinders with axial periodic boundary conditions, there is no change in the perturbation’s kinetic energy and no possible growth of any perturbation, as was originally found by the linear stability study of Kelvin (Reference Kelvin1880). However, in all other flows that are stable according to Rayleigh’s (Reference Rayleigh1916) inviscid linear stability criterion, including the circular Couette flow and the Lamb–Oseen vortex, there exists a range of flow parameters where $0<{\it\lambda}^{-}<{\it\lambda}^{+}$ . For such flows, the perturbation’s kinetic energy may grow because the perturbation’s azimuthal velocity may convert into radial and axial velocity perturbations. However, this growth is bounded according to (3.9). The inequality (3.9) provides both an estimate and an upper bound for the nonlinear growth of perturbations. From this perspective, (3.9) is also a generalized, axisymmetric, inviscid nonlinear stability criterion for columnar rotating flows.

3.2 Supporting identities for a viscous axisymmetric columnar rotating flow

In this section, we extend the inviscid flow nonlinear stability result by Wang (Reference Wang2009) to a perturbed columnar viscous rotating flow with an axisymmetric finite-amplitude perturbation between two coaxial rotating cylinders, together with axial periodic boundary conditions. We consider the reduced Arnol’d function defined by (3.6) and use the difference given by ${\rm\Delta}\mathscr{A}_{rd}({\it\chi},K):=\mathscr{A}_{rd}({\it\chi},K)-\mathscr{A}_{rd}({\it\chi}^{C},K^{C})$ as a Lyapunov function. We look to determine the operational conditions for which $\text{d}{\rm\Delta}\mathscr{A}_{rd}/\text{d}t\leqslant 0$ for all times $t\geqslant 0$ , and from this to extend the inviscid flow kinetic energy estimates of the perturbation given by the inequality (3.9) to the viscous flow case.

We develop two supporting identities for an axisymmetric viscous columnar rotating flow. In both identities, we consider a base steady, axisymmetric, incompressible, viscous columnar rotating flow in the domain $V$ , where $u_{r}^{C}=u_{x}^{C}=0$ and $u_{{\it\theta}}=u_{{\it\theta}}^{C}(y)$ . The base flow is perturbed according to (3.2).

Identity 3.1. The temporal rate of change of the total kinetic energy integral $\mathscr{E}$ (defined by (3.3)) is given by

(3.11) $$\begin{eqnarray}\displaystyle \frac{\text{d}\mathscr{E}}{\text{d}t}={\it\mu}(I_{1}+I_{2}), & & \displaystyle\end{eqnarray}$$

where $I_{1}$ is the boundary integral containing the first-order perturbation term,

(3.12) $$\begin{eqnarray}\displaystyle I_{1}=4{\rm\pi}\int _{0}^{L}\left[yu_{{\it\theta}}^{C}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\right]_{y=y_{1}}^{y=y_{2}}\,\text{d}x, & & \displaystyle\end{eqnarray}$$

and $I_{2}$ is the volume integral containing the sum of second-order perturbation terms,

(3.13) $$\begin{eqnarray}\displaystyle I_{2}=-\int _{V}\left(\frac{\tilde{u} _{r}^{2}}{2y}+\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{r}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{x}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{{\it\theta}}|^{2}\right)\text{d}V. & & \displaystyle\end{eqnarray}$$

Here, the gradient operator $\boldsymbol{{\rm\nabla}}$ written in axisymmetric cylindrical coordinates is given by

(3.14) $$\begin{eqnarray}\displaystyle \boldsymbol{{\rm\nabla}}:=\frac{\partial }{\partial x}\boldsymbol{e}_{x}+\sqrt{2y}\frac{\partial }{\partial y}\boldsymbol{e}_{r}, & & \displaystyle\end{eqnarray}$$

with $|\boldsymbol{{\rm\nabla}}\tilde{u} _{r}|^{2}:=(\boldsymbol{{\rm\nabla}}\tilde{u} _{r})\boldsymbol{\cdot }(\boldsymbol{{\rm\nabla}}\tilde{u} _{r})=(\partial _{x}\tilde{u} _{r})^{2}+2y(\partial _{y}\tilde{u} _{r})^{2}$ , and similar formulas apply for $\tilde{u} _{x}$ and $\tilde{u} _{{\it\theta}}$ .

Identity 3.2. For any choice of the sufficiently smooth function $j(K)$ , the temporal rate of change of the circulation-dependent integral (Casimir’s integral) $\mathscr{C}_{S}$ (defined by (3.5)) is given by

(3.15) $$\begin{eqnarray}\displaystyle \frac{\text{d}\mathscr{C}_{S}}{\text{d}t}={\it\mu}(S_{1}+S_{2}), & & \displaystyle\end{eqnarray}$$

where $S_{1}$ is the sum of integrals containing the first-order terms of the circulation perturbation $\tilde{K}=K-K^{C}(y)$ ,

(3.16) $$\begin{eqnarray}\displaystyle S_{1} & = & \displaystyle 4{\rm\pi}\int _{0}^{L}\left[yj^{\prime }(K^{C})\frac{\partial \tilde{K}}{\partial y}\right]_{y=y_{1}}^{y=y_{2}}\,\text{d}x\nonumber\\ \displaystyle & & \displaystyle -\,2\int _{V}\left[yj^{\prime \prime \prime }(K^{C})\left(\frac{\text{d}K^{C}}{\text{d}y}\right)^{2}\tilde{K}+2yj^{\prime \prime }(K^{C})\frac{\text{d}K^{C}}{\text{d}y}\frac{\partial \tilde{K}}{\partial y}\right]\text{d}V,\end{eqnarray}$$

and $S_{2}$ is the sum of integrals containing the higher-order perturbation terms of $\tilde{K}$ ,

(3.17) $$\begin{eqnarray}\displaystyle & & \displaystyle \!\!\!\!S_{2}=-\int _{V}j^{\prime \prime }(K)\left(\frac{\partial \tilde{K}}{\partial x}\right)^{2}\,\text{d}V\nonumber\\ \displaystyle & & \displaystyle \quad -\int _{V}y\left[j^{\prime \prime \prime \prime }(K\rfloor _{{\it\zeta}})\!\left(\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\right)^{2}\!\tilde{K}^{2}+4j^{\prime \prime \prime }(K\rfloor _{{\it\zeta}})\!\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\tilde{K}\frac{\partial \tilde{K}}{\partial y}+2j^{\prime \prime }(K\rfloor _{{\it\zeta}})\!\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\right]\text{d}V.\nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

Here, the expressions for the intermediate instantaneous functions $K\rfloor _{{\it\zeta}}$ and $\partial K/\partial y\rfloor _{{\it\zeta}}$ are given by

(3.18) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle K\rfloor _{{\it\zeta}}=K^{C}(y)+{\it\zeta}(t)[K(x,y,t)-K^{C}(y)],\\ \displaystyle \left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}=\frac{\text{d}K^{C}(y)}{\text{d}y}+{\it\zeta}(t)\left[\frac{\partial K(x,y,t)}{\partial y}-\frac{\text{d}K^{C}(y)}{\text{d}y}\right],\quad \text{with }0<{\it\zeta}(t)<1.\end{array}\right\} & & \displaystyle\end{eqnarray}$$

Remark.

For a general choice of the sufficiently smooth function $j(K)$ , the derivation of Identity 3.2 requires the use of the second-order Taylor’s expansion with a remainder. This derivation uses essential calculus techniques and is shown in appendix B. This identity is new and is applicable for any choice of the sufficiently smooth function $j(K)$ . For example, for the special choices such as $j(K)=K$ and $j(K)=K^{2}$ , we find that

(3.19) $$\begin{eqnarray}\displaystyle \frac{\text{d}}{\text{d}t}\left[{\it\rho}\int _{V}K\,\text{d}V\right]=4{\rm\pi}{\it\mu}\int _{0}^{L}\left[y\frac{\partial \tilde{K}}{\partial y}\right]_{y=y_{1}}^{y=y_{2}}\,\text{d}x, & & \displaystyle\end{eqnarray}$$

and

(3.20) $$\begin{eqnarray}\displaystyle \frac{\text{d}}{\text{d}t}\left[{\it\rho}\int _{V}K^{2}\,\text{d}V\right] & = & \displaystyle 8{\rm\pi}{\it\mu}\int _{0}^{L}\left[yK^{C}\frac{\partial \tilde{K}}{\partial y}\right]_{y=y_{1}}^{y=y_{2}}\,\text{d}x-8{\it\mu}\int _{V}y\frac{\text{d}K^{C}}{\text{d}y}\frac{\partial \tilde{K}}{\partial y}\,\text{d}V\nonumber\\ \displaystyle & & \displaystyle -\,2{\it\mu}\int _{V}\left(\frac{\partial \tilde{K}}{\partial x}\right)^{2}\,\text{d}V-4{\it\mu}\int _{V}y\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\,\text{d}V,\end{eqnarray}$$

respectively. The result (3.19) reflects that, in the inviscid flow case, the integrated $K$ in the bulk is a conserved quantity in time; while in the viscous flow case, this quantity changes only due to viscous forces along the walls of the coaxial rotating cylinders. On the other hand, the result (3.20) reflects that, in the inviscid flow case, the integrated $K^{2}$ in the bulk is also a conserved quantity in time; while in the viscous flow case, the changes in this quantity are due to the combination of both viscous forces along the walls of the coaxial rotating cylinders as well as the internal viscous stresses. This second situation reflects a more complicated nonlinear effect of the viscous stresses on flow dynamics.

3.3 Basic formula for Lyapunov stability of the circular Couette flow

Using the identities (3.11) and (3.15), the temporal rate of change of ${\rm\Delta}\mathscr{A}_{rd}$ is given by

(3.21) $$\begin{eqnarray}\displaystyle \frac{\text{d}{\rm\Delta}\mathscr{A}_{rd}}{\text{d}t}=\frac{\text{d}\mathscr{A}_{rd}}{\text{d}t}=\frac{\text{d}\mathscr{E}}{\text{d}t}+\frac{\text{d}\mathscr{C}_{S}}{\text{d}t}={\it\mu}(I_{1}+S_{1}+I_{2}+S_{2}). & & \displaystyle\end{eqnarray}$$

We first compute the term $I_{1}+S_{1}$ in (3.21). From (3.12) and (3.16) we find that

(3.22) $$\begin{eqnarray}\displaystyle I_{1}+S_{1} & = & \displaystyle 4{\rm\pi}\int _{0}^{L}\left[yu_{{\it\theta}}^{C}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\right]_{y=y_{1}}^{y=y_{2}}\,\text{d}x+4{\rm\pi}\int _{0}^{L}\left[yj^{\prime }(K^{C})\frac{\partial \tilde{K}}{\partial y}\right]_{y=y_{1}}^{y=y_{2}}\,\text{d}x\nonumber\\ \displaystyle & & \displaystyle -\,2\int _{V}\left[yj^{\prime \prime \prime }(K^{C})\left(\frac{\text{d}K^{C}}{\text{d}y}\right)^{2}\tilde{K}+2yj^{\prime \prime }(K^{C})\frac{\text{d}K^{C}}{\text{d}y}\frac{\partial \tilde{K}}{\partial y}\right]\text{d}V.\end{eqnarray}$$

We now use the special choice for the function $j(K)$ according to the inviscid approach, given by (3.7), i.e. $j^{\prime }(K^{C})=-K^{C}/2Y(K^{C})$ . The circulation function for the base steady circular Couette flow is $K^{C}(y)=2Ay+B$ ; therefore $Y(K)=(K-B)/2A$ and the various derivatives of $j(K)$ are given by

(3.23a-d ) $$\begin{eqnarray}\displaystyle j^{\prime }(K)=-\frac{AK}{K-B},\quad j^{\prime \prime }(K)=\frac{AB}{(K-B)^{2}},\quad j^{\prime \prime \prime }(K)=\frac{-2AB}{(K-B)^{3}},\quad j^{\prime \prime \prime \prime }(K)=\frac{6AB}{(K-B)^{4}}. & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

Inserting (3.23) into (3.22), and noting that the two boundary terms in (3.22) cancel each other identically under this current choice for the function $j(K)$ , we obtain

(3.24) $$\begin{eqnarray}\displaystyle I_{1}+S_{1} & = & \displaystyle -2\int _{V}\left[\frac{-8A^{3}By}{(K^{C}(y)-B)^{3}}\tilde{K}+\frac{4A^{2}By}{(K^{C}(y)-B)^{2}}\frac{\partial \tilde{K}}{\partial y}\right]\text{d}V\nonumber\\ \displaystyle & = & \displaystyle -2B\int _{V}\left[\frac{-\tilde{K}}{y^{2}}+\frac{1}{y}\frac{\partial \tilde{K}}{\partial y}\right]\text{d}V=0.\end{eqnarray}$$

The last step in (3.24) is obtained through the application of integration by parts with respect to $y$ . It is remarkable that, under this very choice of the function $j(K)$ as defined by (3.7), the sum of the first-order terms in the temporal rate of change of ${\rm\Delta}\mathscr{A}_{rd}$ (given by (3.21)) vanishes trivially in the case of the circular Couette flow.

Substituting the various derivatives of the function $j(K)$ given by (3.23) into the equation for $S_{2}$ according to (3.17), results in

(3.25) $$\begin{eqnarray}\displaystyle \frac{\text{d}{\rm\Delta}\mathscr{A}_{rd}}{\text{d}t} & = & \displaystyle -{\it\mu}\int _{V}\left(\frac{\tilde{u} _{r}^{2}}{2y}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{r}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{x}|^{2}\right)\text{d}V\nonumber\\ \displaystyle & & \displaystyle -\,{\it\mu}\int _{V}\left[\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial x}\right)^{2}+2y\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\right)^{2}\right]\text{d}V-{\it\mu}\int _{V}\frac{AB}{(K-B)^{2}}\left(\frac{\partial \tilde{K}}{\partial x}\right)^{2}\text{d}V\nonumber\\ \displaystyle & & \displaystyle -\,{\it\mu}\int _{V}\left[\frac{6ABy}{(K\rfloor _{{\it\zeta}}-B)^{4}}\left(\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\right)^{2}\tilde{K}^{2}-\frac{8ABy}{(K\rfloor _{{\it\zeta}}-B)^{3}}\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\tilde{K}\frac{\partial \tilde{K}}{\partial y}\right.\nonumber\\ \displaystyle & & \displaystyle +\left.\frac{2ABy}{(K\rfloor _{{\it\zeta}}-B)^{2}}\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\right]\text{d}V.\end{eqnarray}$$

We find from (3.25) that the temporal rate of change of ${\rm\Delta}\mathscr{A}_{rd}$ is not time-invariant, since all the time-invariant quantities that were used in the inviscid flow case (§ 3.1) cease to exist.

To this end, we attempt to use ${\rm\Delta}\mathscr{A}_{rd}$ given by (3.8) as a Lyapunov function. Then, for nonlinear stability in response to axisymmetric perturbations, we require that $\text{d}{\rm\Delta}\mathscr{A}_{rd}/\text{d}t\leqslant 0$ . This implies that ${\rm\Delta}\mathscr{A}_{rd}\leqslant {\rm\Delta}\mathscr{A}_{rd}|_{t=0}$ for all times $t\geqslant 0$ . We clarify here that the equivalence of norms no longer holds since the underlying system is infinite-dimensional. Therefore, a positive second variation of $\mathscr{A}_{rd}$ (when it is evaluated at the circular Couette flow, given its first variation also vanishes identically at this base flow) does not guarantee the local convexity of ${\rm\Delta}\mathscr{A}_{rd}$ around the base flow. However, this is a necessary condition for local convexity (and hence for nonlinear stability), which requires a strict satisfaction of the classical inviscid linear stability criterion of Rayleigh (Reference Rayleigh1916) (when ${\rm\Delta}\mathscr{A}_{rd}$ is used as a Lyapunov function). In the inviscid case, this condition is also known as formal stability, and is also required in the viscous case, since the same Lyapunov function is used. Furthermore, when the mean value theorem is applied in the formulation of ${\rm\Delta}\mathscr{A}_{rd}$ , nonlinear stability also explicitly requires the existence of positive finite value convexity estimates (bounding constants), $0<{\it\lambda}^{-}<{\it\lambda}^{+}<\infty$ , such that the inequalities

(3.26) $$\begin{eqnarray}\displaystyle \Vert \tilde{\boldsymbol{u}}(t)\Vert _{{\it\lambda}^{-}}^{2}\leqslant {\rm\Delta}\mathscr{A}_{rd}|_{t}\leqslant {\rm\Delta}\mathscr{A}_{rd}|_{t=0}\leqslant \Vert \tilde{\boldsymbol{u}}(t=0)\Vert _{{\it\lambda}^{+}}^{2}\leqslant \frac{{\it\lambda}^{+}}{{\it\lambda}^{-}}\Vert \tilde{\boldsymbol{u}}(t=0)\Vert _{{\it\lambda}^{-}}^{2} & & \displaystyle\end{eqnarray}$$

hold for all times $t\geqslant 0$ . More specifically, the norms of the perturbation are given by

(3.27) $$\begin{eqnarray}\displaystyle \Vert \tilde{\boldsymbol{u}}(t)\Vert _{{\it\lambda}^{\mp }}^{2}:=\frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+{\it\lambda}^{\mp }\tilde{u} _{{\it\theta}}^{2})\,\text{d}V, & & \displaystyle\end{eqnarray}$$

where ${\it\lambda}^{\mp }$ represents the minimum/maximum of the coefficient of the square of the azimuthal velocity given by (3.10).

Note. The existence of positive finite value convexity estimates (bounding constants) is sufficient to ensure the local convexity of ${\rm\Delta}\mathscr{A}_{rd}$ around the circular Couette flow. This is with respect to the stability norm $\Vert \cdot \Vert _{{\it\lambda}^{-}}^{2}$ (or any other norm that is equivalent to $\Vert \cdot \Vert _{{\it\lambda}^{-}}^{2}$ ), where the stability norm $\Vert \cdot \Vert _{{\it\lambda}^{-}}^{2}$ is induced directly from the second variation of $\mathscr{A}_{rd}$ (when it is evaluated at the circular Couette flow).

The first term on the right-hand side of (3.25) is always non-positive and supports the nonlinear stability requirement. Furthermore, we regroup the middle term in the second integral on the right-hand side of (3.25) with the third integral and find that

(3.28) $$\begin{eqnarray}\displaystyle -\int _{V}\left[\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial x}\right)^{2}+\frac{AB}{(K-B)^{2}}\left(\frac{\partial \tilde{K}}{\partial x}\right)^{2}\right]\text{d}V=-\int _{V}\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial x}\right)^{2}\left(1+\frac{2ABy}{(K-B)^{2}}\right)\text{d}V. & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

For a non-positive value of these two combined terms, we impose the first sufficient condition for the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ (in the case of the circular Couette flow), i.e.

(3.29) $$\begin{eqnarray}\displaystyle 1+\frac{2ABy}{(K-B)^{2}}\geqslant 0 & & \displaystyle\end{eqnarray}$$

for $y_{1}<y<y_{2}$ . In addition, we also impose the second sufficient condition for the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ (in the case of the circular Couette flow), i.e.

(3.30) $$\begin{eqnarray}\displaystyle & & \displaystyle -\int _{V}\left[\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+2y\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\right)^{2}\right]\text{d}V\nonumber\\ \displaystyle & & \displaystyle \quad -\int _{V}\left[\frac{6ABy}{(K\rfloor _{{\it\zeta}}-B)^{4}}\left(\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\right)^{2}\tilde{K}^{2}-\frac{8ABy}{(K\rfloor _{{\it\zeta}}-B)^{3}}\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\tilde{K}\frac{\partial \tilde{K}}{\partial y}\right.\nonumber\\ \displaystyle & & \displaystyle \quad +\left.\frac{2ABy}{(K\rfloor _{{\it\zeta}}-B)^{2}}\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\right]\text{d}V\leqslant 0.\end{eqnarray}$$

In § 3.4, we analyse the two imposed conditions (3.29) and (3.30) and show that for infinitesimally small perturbations they result in Rayleigh’s (Reference Rayleigh1916) inviscid linear stability criterion for the circular Couette flow. In § 3.5, we analyse the conditions (3.29) and (3.30), and derive a nonlinear stability criterion for the viscous circular Couette flow in response to any finite-amplitude axisymmetric perturbation.

3.4 The linearized stability criterion

We first study the linearized problem where infinitesimally small perturbations are considered, i.e. $K=K^{C}+{\it\epsilon}\tilde{K}$ , where $K^{C}=2Ay+B$ and $0<{\it\epsilon}\ll 1$ . Then, the leading-order $O(1)$ term of the first sufficient condition for the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ (in the case of the circular Couette flow) given by (3.29) is required to satisfy

(3.31) $$\begin{eqnarray}\displaystyle 1+\frac{B}{2Ay}\geqslant 0 & & \displaystyle\end{eqnarray}$$

for $y_{1}<y<y_{2}$ . The condition (3.31) is exactly the classical inviscid flow linear stability criterion of Rayleigh (Reference Rayleigh1916) when applied to the circular Couette flow.

In the leading order $O(1)$ , we also have $K\rfloor _{{\it\zeta}}=2Ay+B$ . By inserting this into the second sufficient condition for the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ (in the case of the circular Couette flow) given by (3.30), we find the following leading-order $O(1)$ condition:

(3.32) $$\begin{eqnarray}\displaystyle -\int _{V}\left[\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+2y\left(\!\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\!\right)^{2}\right]\text{d}V-\frac{B}{A}\int _{V}\left[\frac{3}{2y^{3}}\tilde{K}^{2}-\frac{2}{y^{2}}\tilde{K}\frac{\partial \tilde{K}}{\partial y}+\frac{1}{2y}\left(\!\frac{\partial \tilde{K}}{\partial y}\!\right)^{2}\right]\text{d}V\leqslant 0.\qquad & & \displaystyle\end{eqnarray}$$

Using the relationship $\tilde{K}=\sqrt{2y}\tilde{u} _{{\it\theta}}$ , and the application of integration by parts with the no-slip boundary conditions $\tilde{u} _{{\it\theta}}(x,y_{1},t)=\tilde{u} _{{\it\theta}}(x,y_{2},t)=0$ (or equivalently $\tilde{K}(x,y_{1},t)=\tilde{K}(x,y_{2},t)=0$ ), then the various terms in (3.32) can also be rewritten as

(3.33) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle \int _{V}\left[\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+2y\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\right)^{2}\right]\text{d}V=\int _{V}\left(\sqrt{2y}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}-\frac{\tilde{u} _{{\it\theta}}}{\sqrt{2y}}\right)^{2}\text{d}V,\\ \displaystyle \int _{V}\frac{1}{y^{2}}\tilde{K}\frac{\partial \tilde{K}}{\partial y}\,\text{d}V=\int _{V}\frac{1}{y^{3}}\tilde{K}^{2}\,\text{d}V,\\ \displaystyle \int _{V}\frac{1}{y}\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\text{d}V=\int _{V}\frac{1}{y}\left(\sqrt{2y}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}-\frac{\tilde{u} _{{\it\theta}}}{\sqrt{2y}}\right)^{2}\text{d}V+\int _{V}\frac{1}{y^{3}}\tilde{K}^{2}\,\text{d}V.\end{array}\right\} & & \displaystyle\end{eqnarray}$$

Combining the above simplifications, the leading-order $O(1)$ condition given by (3.32) becomes

(3.34) $$\begin{eqnarray}\displaystyle -\int _{V}\left(1+\frac{B}{2Ay}\right)\left(\sqrt{2y}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}-\frac{\tilde{u} _{{\it\theta}}}{\sqrt{2y}}\right)^{2}\text{d}V\leqslant 0. & & \displaystyle\end{eqnarray}$$

As a result, the leading-order $O(1)$ condition for (3.29) and (3.30) yields the same requirement as in (3.31). In particular, this requirement is the same as the classical inviscid flow linear stability criterion of Rayleigh (Reference Rayleigh1916).

The substitution of $K\rfloor _{{\it\xi}}=2Ay+B$ (the leading-order term for $K\rfloor _{{\it\xi}}$ ) into (3.8) gives

(3.35) $$\begin{eqnarray}\displaystyle \tilde{\mathscr{E}}_{w}=\frac{1}{2}{\it\rho}\int _{V}\left(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+\left(1+\frac{B}{2Ay}\right)\tilde{u} _{{\it\theta}}^{2}\right)\text{d}V, & & \displaystyle\end{eqnarray}$$

where $\tilde{\mathscr{E}}_{w}$ is the linearized form of ${\rm\Delta}\mathscr{A}_{rd}({\it\chi},K)$ in (3.8) (a linearized weighted kinetic energy integral of the perturbation). The time derivative of $\tilde{\mathscr{E}}_{w}$ is the leading-order $O(1)$ time derivative of ${\rm\Delta}\mathscr{A}_{rd}({\it\chi},K)$ , which follows directly from the identity (3.25). It is given by

(3.36) $$\begin{eqnarray}\displaystyle \frac{\text{d}\tilde{\mathscr{E}}_{w}}{\text{d}t} & = & \displaystyle -{\it\mu}\int _{V}\left(\frac{\tilde{u} _{r}^{2}}{2y}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{r}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{x}|^{2}\right)\text{d}V\nonumber\\ \displaystyle & & \displaystyle -\,{\it\mu}\int _{V}\left(1+\frac{B}{2Ay}\right)\left[\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial x}\right)^{2}+\left(\sqrt{2y}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}-\frac{\tilde{u} _{{\it\theta}}}{\sqrt{2y}}\right)^{2}\right]\text{d}V,\end{eqnarray}$$

where all the above perturbations are considered to be infinitesimal, i.e. $u_{r}={\it\epsilon}\tilde{u} _{r}$ , $u_{{\it\theta}}=u_{{\it\theta}}^{C}+{\it\epsilon}\tilde{u} _{{\it\theta}}$ and $u_{x}={\it\epsilon}\tilde{u} _{x}$ with $0<{\it\epsilon}\ll 1$ . In the case when the condition (3.31) is satisfied for $y_{1}<y<y_{2}$ , then this also ensures (i) that $\tilde{\mathscr{E}}_{w}$ is positive definite, and hence becomes a norm for all the axisymmetric infinitesimally small perturbations around the viscous circular Couette base flow, and (ii) the monotonic temporal decay of $\tilde{\mathscr{E}}_{w}$ for all times $t\geqslant 0$ . Hence, $\tilde{\mathscr{E}}_{w}$ is also a Lyapunov function. Thus, under the linearized dynamics, any small axisymmetric perturbation eventually decays monotonically to zero as time increases.

This linearized stability result is identical to that obtained in Kloosterziel (Reference Kloosterziel2010), where a Lyapunov function has been constructed from the set of linearized governing equations. The Lyapunov function used by Kloosterziel (Reference Kloosterziel2010) has the same exact form as (3.35). However, since it is derived directly from the linearized equations, it cannot be extended to study the evolution of finite-amplitude perturbations. The present approach, however, is essentially targeted to the study of the nonlinear stability problem with finite-amplitude perturbations.

3.5 Nonlinear stability criterion

3.5.1 The basic nonlinear stability criterion

We derive the following nonlinear stability criterion:

Criterion 3.2. The viscous circular Couette flow with parameters $A>0$ and $B\geqslant 0$ , such that

(3.37) $$\begin{eqnarray}\displaystyle 0\leqslant \frac{BR_{2}^{2}}{AR_{1}^{4}}\leqslant 3, & & \displaystyle\end{eqnarray}$$

is nonlinearly stable to any finite-amplitude axisymmetric perturbation whose initial total circulation function $K(x,y,t=0)$ satisfies $K_{min}^{C}\leqslant K(x,y,t=0)\leqslant K_{max}^{C}$ (with no additional restrictions on the amplitudes of the initial radial and axial velocity perturbations, respectively). Such finite-amplitude axisymmetric perturbations must also obey the following energy estimate for all times $t\geqslant 0$ :

(3.38) $$\begin{eqnarray}\displaystyle \frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+{\it\lambda}^{-}\tilde{u} _{{\it\theta}}^{2})\,\text{d}V\leqslant \frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r0}^{2}+\tilde{u} _{x0}^{2}+{\it\lambda}^{+}\tilde{u} _{{\it\theta}0}^{2})\,\text{d}V, & & \displaystyle\end{eqnarray}$$

which is independent of fluid viscosity. The minimum and maximum values of the circulation function $K^{C}(y)$ for the base circular Couette flow are $K_{min}^{C}=AR_{1}^{2}+B$ and $K_{max}^{C}=AR_{2}^{2}+B$ , respectively. Convexity estimates (bounding constants) ${\it\lambda}^{-}>0$ and ${\it\lambda}^{+}({>}{\it\lambda}^{-})$ are given by

(3.39a,b ) $$\begin{eqnarray}\displaystyle {\it\lambda}^{-}=1+\frac{(1-{\it\beta}){\it\alpha}^{4}}{{\it\beta}-{\it\alpha}^{2}},\quad {\it\lambda}^{+}=1+\frac{1-{\it\beta}}{({\it\beta}-{\it\alpha}^{2}){\it\alpha}^{2}}, & & \displaystyle\end{eqnarray}$$

where the non-dimensional variables are ${\it\alpha}=R_{1}/R_{2}$ and ${\it\beta}={\it\omega}_{2}/{\it\omega}_{1}$ . The initial velocity perturbation vector field is $\tilde{\boldsymbol{u}}_{0}=(\tilde{u} _{r0},\tilde{u} _{x0},\tilde{u} _{{\it\theta}0})$ .

Proof. Satisfaction of the sufficient conditions (3.29) and (3.30) are both required for the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ (in the case of the circular Couette flow). The condition (3.29) is automatically satisfied given the conditions $A>0$ and $B\geqslant 0$ .

We seek a condition to ensure that (3.30) is also satisfied. By reorganizing the second term in (3.30) and completing the square, this term becomes

(3.40) $$\begin{eqnarray}\displaystyle & & \displaystyle \!\!\!\!\int _{V}\frac{ABy}{(K\rfloor _{{\it\zeta}}-B)^{2}}\left[\frac{6}{(K\rfloor _{{\it\zeta}}-B)^{2}}\left(\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\right)^{2}\tilde{K}^{2}-\frac{8}{(K\rfloor _{{\it\zeta}}-B)}\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\tilde{K}\frac{\partial \tilde{K}}{\partial y}+2\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\right]\text{d}V\nonumber\\ \displaystyle & & \displaystyle \quad =\int _{V}\frac{ABy}{(K\rfloor _{{\it\zeta}}-B)^{2}}\left[\left(\frac{\sqrt{6}}{(K\rfloor _{{\it\zeta}}-B)}\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\tilde{K}-\frac{2\sqrt{2}}{\sqrt{3}}\frac{\partial \tilde{K}}{\partial y}\right)^{2}-\frac{2}{3}\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\right]\text{d}V.\end{eqnarray}$$

Using the relationship $\tilde{K}=\sqrt{2y}\tilde{u} _{{\it\theta}}$ , we also find that

(3.41) $$\begin{eqnarray}\displaystyle \int _{V}\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\text{d}V=\int _{V}\left[\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+2\tilde{u} _{{\it\theta}}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}+2y\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\right)^{2}\right]\text{d}V=\int _{V}\left[\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+2y\left(\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\right)^{2}\right]\text{d}V, & & \displaystyle \nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

since

(3.42) $$\begin{eqnarray}\displaystyle \int _{V}2\tilde{u} _{{\it\theta}}\frac{\partial \tilde{u} _{{\it\theta}}}{\partial y}\,\text{d}V=\int _{V}\frac{\partial (\tilde{u} _{{\it\theta}})^{2}}{\partial y}\,\text{d}V=0, & & \displaystyle\end{eqnarray}$$

because of the no-slip boundary conditions imposed along the walls on the coaxial rotating cylinders.

Inserting (3.40) and (3.41) into (3.30) gives

(3.43) $$\begin{eqnarray}\displaystyle & & \displaystyle -\int _{V}\left(1-\frac{2}{3}\frac{ABy}{(K\rfloor _{{\it\zeta}}-B)^{2}}\right)\left(\frac{\partial \tilde{K}}{\partial y}\right)^{2}\text{d}V\nonumber\\ \displaystyle & & \displaystyle \quad -\int _{V}\frac{ABy}{(K\rfloor _{{\it\zeta}}-B)^{2}}\left(\frac{\sqrt{6}}{(K\rfloor _{{\it\zeta}}-B)}\left.\frac{\partial K}{\partial y}\right\rfloor _{{\it\zeta}}\tilde{K}-\frac{2\sqrt{2}}{\sqrt{3}}\frac{\partial \tilde{K}}{\partial y}\right)^{2}\text{d}V\leqslant 0.\end{eqnarray}$$

As a result (given that $A>0$ and $B\geqslant 0$ ), the condition for $\text{d}{\rm\Delta}\mathscr{A}_{rd}/\text{d}t\leqslant 0$ requires

(3.44) $$\begin{eqnarray}\displaystyle 1-\frac{2ABy}{3(K\rfloor _{{\it\zeta}}-B)^{2}}\geqslant 0 & & \displaystyle\end{eqnarray}$$

for all $(x,y)$ in the domain $D$ and for all times $t\geqslant 0$ . Let $K_{min}=\min K(x,y,t=0)$ and $K_{max}=\max K(x,y,t=0)$ be the minimum and maximum of the initial total circulation function, respectively. By using a maximum principle, we prove in appendix C that

(3.45) $$\begin{eqnarray}\displaystyle K_{min}\leqslant K(x,y,t)\leqslant K_{max} & & \displaystyle\end{eqnarray}$$

for all times $t\geqslant 0$ . This has a clear physical meaning: in the axisymmetric viscous flow case, $K(x,y,t)$ is only transported by the flow motion and is diffused by the viscous effect, but it cannot surpass $K_{min}$ and $K_{max}$ for all times $t\geqslant 0$ . With (3.45), (3.44) becomes

(3.46) $$\begin{eqnarray}\displaystyle 1-\frac{ABR_{2}^{2}}{3(K_{min}-B)^{2}}\geqslant 0. & & \displaystyle\end{eqnarray}$$

By considering $K_{min}=K_{min}^{C}(y_{1})=AR_{1}^{2}+B$ , we obtain the nonlinear stability condition (3.37) for the circular Couette flow.

In addition, the formulas (3.39) are derived directly from (3.10) and are independent of the viscosity. In fact, using $j^{\prime \prime }(K)=AB/(K-B)^{2}$ , $K^{C}=2Ay+B$ , and also by requiring $K_{max}=K_{max}^{C}(y_{2})=AR_{2}^{2}+B$ , (3.10) becomes

(3.47) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle {\it\lambda}^{-}=\min _{y_{1}\leqslant y,~Y\leqslant y_{2}}\left(1+\frac{By}{2AY^{2}}\right)=1+\frac{By_{1}}{2Ay_{2}^{2}},\\ \displaystyle {\it\lambda}^{+}=\max _{y_{1}\leqslant y,~Y\leqslant y_{2}}\left(1+\frac{By}{2AY^{2}}\right)=1+\frac{By_{2}}{2Ay_{1}^{2}},\end{array}\right\} & & \displaystyle\end{eqnarray}$$

where $y_{1}=R_{1}^{2}/2$ and $y_{2}=R_{2}^{2}/2$ . The results in (3.39) follow from (3.47), the definitions of $A$ and $B$ in (2.7), and also the definitions of the non-dimensional parameters ${\it\alpha}$ and ${\it\beta}$ . ☐

Remark.

The maximum principle plays a critical part in deriving Criterion 3.2, as it enables the explicit determination of the conditions for (i)  $0<{\it\lambda}^{-}<{\it\lambda}^{+}<\infty$ , that is, the positivity of finite value convexity estimates (bounding constants), and also (ii) the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ . In particular, since the total circulation function $K(x,y,t)$ is bounded below and above by the respective minimum and maximum circulation values of the circular Couette base flow for all times $t\geqslant 0$ , this then allows the nonlinear stability conditions to be rewritten in terms of only the non-dimensional parameters ${\it\alpha}$ and ${\it\beta}$ . The current convexity estimates are chosen to be the same as those used by Wang (Reference Wang2009) in the study of inviscid flow nonlinear stability; hence the stability region given by ${\it\lambda}^{-}>0$ remains the same as in the inviscid nonlinear stability case, which also coincides with Rayleigh’s (Reference Rayleigh1916) inviscid linear stability region for $B\geqslant 0$ . Under the current formulation, it turns out that the stability region for the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ forms a subset of the stability region given by ${\it\lambda}^{-}>0$ . As a result, the nonlinear stability of the viscous circular Couette flow in response to finite-amplitude axisymmetric perturbations is simply to require the monotonic temporal decay of ${\rm\Delta}\mathscr{A}_{rd}$ (given by (3.37)), as the condition for the positivity of convexity estimates is implicitly satisfied.

Using the non-dimensional parameters ${\it\alpha}=R_{1}/R_{2}$ and ${\it\beta}={\it\omega}_{2}/{\it\omega}_{1}$ , together with the definitions of $A$ and $B$ in (2.7), gives a direct evaluation of condition (3.37) in terms of ${\it\alpha}$ and ${\it\beta}$ ,

(3.48) $$\begin{eqnarray}\displaystyle 0<p({\it\alpha})=\frac{1+3{\it\alpha}^{4}}{1+3{\it\alpha}^{2}}\leqslant {\it\beta}\leqslant 1. & & \displaystyle\end{eqnarray}$$

We note that the inviscid flow linear stability criterion is $K^{C}(y_{2})>K^{C}(y_{1})\geqslant 0$ and is given by the region of the parameters ${\it\alpha}^{2}<{\it\beta}$ . We find that the viscous flow nonlinear stability region given by (3.48) is a subregion of the inviscid linear stability region (see the ${\it\beta}$ versus ${\it\alpha}$ diagram in figure 1). The viscous flow nonlinear stability of the circular Couette flow is a global stability result. Inside this stability region, any finite-amplitude axisymmetric perturbation whose initial total circulation function $K(x,y,t=0)=K^{C}(y)+\tilde{K}(x,y,t=0)$ satisfies $K_{min}^{C}\leqslant K(x,y,t=0)\leqslant K_{max}^{C}$ decays monotonically in time to zero, and the eventual reappearance of the steady circular Couette base flow determined uniquely by the motion of the inner and outer cylinders’ walls is assured (provided that the perturbed flow remains axisymmetric for all times $t\geqslant 0$ ).

Figure 1.

The viscous flow nonlinear stability region for the circular Couette flow: $p({\it\alpha})\leqslant {\it\beta}\leqslant 1$ . Also shown for reference is the inviscid flow linear stability region for the circular Couette flow by Rayleigh’s (Reference Rayleigh1916) criterion: ${\it\alpha}^{2}<{\it\beta}$ .

3.5.2 An extension of the basic nonlinear stability criterion

In this section, we show that the restriction on the initial total circulation function $K(x,y,t=0)$ such that $K_{min}^{C}\leqslant K(x,y,t=0)\leqslant K_{max}^{C}$ is not necessary. There are two cases to consider.

(i) Case 1: for $B=0$ . This particular circular Couette flow is the solid-body rotation flow within the region $V$ . In this case, the difference of the reduced Arnol’d function ${\rm\Delta}\mathscr{A}_{rd}$ (given by (3.8)) simplifies to the kinetic energy integral of the perturbation. Hence, we have

(3.49) $$\begin{eqnarray}\displaystyle {\rm\Delta}\mathscr{A}_{rd}=\tilde{\mathscr{E}}=\frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+\tilde{u} _{{\it\theta}}^{2})\,\text{d}V. & & \displaystyle\end{eqnarray}$$

Furthermore, the condition (3.37) is trivially satisfied, where the convexity estimates are given by ${\it\lambda}^{-}={\it\lambda}^{+}=1>0$ . As a result, in the case when $B=0$ , the circular Couette flow is nonlinearly stable to any finite-amplitude axisymmetric perturbation with no additional restrictions.

(ii) Case 2: for $B>0$ . In this case, we note that $j(K)$ is an analytic function for $K>B$ and becomes singular only at $K=B$ . If the initial total circulation function satisfies $K(x,y,t=0)>B$ , then this would also ensure that the total circulation function $K(x,y,t)$ is greater than $B$ and remains within the analytic domain of $j(K)$ for all times $t\geqslant 0$ (due to a maximum principle in appendix C). Thus, the basic formula (3.25) holds. Consequently (under the current restriction on the initial total circulation function), the derivation of Criterion 3.2 is still valid up to (3.46) (prior to the restriction on $K_{min}$ to be $K_{min}^{C}$ ). Let $\bar{K}_{min}$ be the minimum value of the initial total circulation function $K(x,y,t=0)$ such that (3.46) holds. Then $\bar{K}_{min}$ must satisfy

(3.50) $$\begin{eqnarray}\displaystyle {\textstyle \frac{1}{3}}ABR_{2}^{2}\leqslant (\bar{K}_{min}-B)^{2}. & & \displaystyle\end{eqnarray}$$

We also consider the case when the restriction on the minimum total circulation value is allowed to be extended below the minimum circulation value of the base flow $K_{min}^{C}$ , i.e. $\bar{K}_{min}=K_{min}^{C}-{\rm\Delta}K$ with ${\rm\Delta}K\geqslant 0$ . Then, (3.50) yields an upper bound for ${\rm\Delta}K$ given by

(3.51) $$\begin{eqnarray}\displaystyle 0\leqslant {\rm\Delta}K\leqslant AR_{1}^{2}-\sqrt{{\textstyle \frac{1}{3}}ABR_{2}^{2}}<K_{min}^{C}. & & \displaystyle\end{eqnarray}$$

In terms of the non-dimensional parameters ${\it\alpha}=R_{1}/R_{2}$ and ${\it\beta}={\it\omega}_{2}/{\it\omega}_{1}$ , this upper bound for ${\rm\Delta}K$ can be written as

(3.52) $$\begin{eqnarray}\displaystyle AR_{1}^{2}-\sqrt{{\textstyle \frac{1}{3}}ABR_{2}^{2}}={\it\omega}_{1}R_{1}^{2}q({\it\alpha},{\it\beta}), & & \displaystyle\end{eqnarray}$$

where $q({\it\alpha},{\it\beta})$ is given by

(3.53) $$\begin{eqnarray}\displaystyle q({\it\alpha},{\it\beta})=\left(\frac{{\it\beta}-{\it\alpha}^{2}}{1-{\it\alpha}^{2}}\right)-\frac{1}{{\it\alpha}}\sqrt{\frac{1}{3}\left(\frac{{\it\beta}-{\it\alpha}^{2}}{1-{\it\alpha}^{2}}\right)\left(\frac{1-{\it\beta}}{1-{\it\alpha}^{2}}\right)}. & & \displaystyle\end{eqnarray}$$

Rescaling ${\rm\Delta}K$ with the minimum circulation value of the base flow $K_{min}^{C}={\it\omega}_{1}R_{1}^{2}$ gives

(3.54) $$\begin{eqnarray}\displaystyle 0\leqslant \frac{{\rm\Delta}K}{{\it\omega}_{1}R_{1}^{2}}=\frac{{\rm\Delta}K}{K_{min}^{C}}\leqslant q({\it\alpha},{\it\beta}). & & \displaystyle\end{eqnarray}$$

For a given value of the ratio ${\rm\Delta}K/K_{min}^{C}$ , we can determine the modified nonlinear stability borderline ${\it\beta}_{{\rm\Delta}K/K_{min}^{C}}({\it\alpha})$ for which $q({\it\alpha},{\it\beta})={\rm\Delta}K/K_{min}^{C}$ . For various values of the ratio ${\rm\Delta}K/K_{min}^{C}$ from 0 to 100 %, the corresponding modified stability regions ${\it\beta}_{{\rm\Delta}K/K_{min}^{C}}({\it\alpha})\leqslant {\it\beta}<1$ on the ${\it\alpha}$ versus ${\it\beta}$ plane are displayed in figure 2. It can be seen that the nonlinear stability region decreases as this ratio increases and is completely eliminated when ${\rm\Delta}K/K_{min}^{C}=1$ (which can never occur unless $B=0$ ).

We also determine the convexity estimates ${\it\lambda}^{-}$ and ${\it\lambda}^{+}$ used in the energy estimate (3.38) for the extended initial perturbation domain on the total circulation function given by $\bar{K}_{min}\leqslant K(x,y,t=0)\leqslant \bar{K}_{max}$ , where $\bar{K}_{min}\in [B+\sqrt{ABR_{2}^{2}/3},K_{min}^{C}]$ (determined by rearranging (3.50)) and $\bar{K}_{max}\in [K_{max}^{C},\infty )$ . The quantity $\bar{K}_{max}$ represents the maximum value of the initial total circulation function $K(x,y,t=0)$ such that (3.46) holds, and since $K_{max}$ does not appear in (3.46), hence there is no upper limit on $\bar{K}_{max}$ .

Once the parameters of the circular Couette flow are given, then the intervals for $\bar{K}_{min}$ and $\bar{K}_{max}$ are fixed. By specifying the values of $\bar{K}_{min}$ and $\bar{K}_{max}$ from their respective intervals, then the modified stability region is given by (3.50). The corresponding convexity estimates ${\it\lambda}^{-}$ and ${\it\lambda}^{+}$ are

(3.55) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle {\it\lambda}^{-}=\min _{\substack{ y_{1}\leqslant y\leqslant y_{2} \\ \bar{K}_{min}\leqslant K\rfloor _{{\it\xi}}\leqslant \bar{K}_{max}}}1+2y\frac{AB}{(K\rfloor _{{\it\xi}}-B)^{2}}=1+\frac{ABR_{1}^{2}}{(\bar{K}_{max}-B)^{2}},\\ \displaystyle {\it\lambda}^{+}=\max _{\substack{ y_{1}\leqslant y\leqslant y_{2} \\ \bar{K}_{min}\leqslant K\rfloor _{{\it\xi}}\leqslant \bar{K}_{max}}}1+2y\frac{AB}{(K\rfloor _{{\it\xi}}-B)^{2}}=1+\frac{ABR_{2}^{2}}{(\bar{K}_{min}-B)^{2}},\end{array}\right\} & & \displaystyle\end{eqnarray}$$

which is an extension of (3.10) with $j^{\prime \prime }(K)=AB/(K-B)^{2}$ .

The condition $0<{\it\lambda}^{-}<{\it\lambda}^{+}<\infty$ is satisfied without any further reductions in size of the modified stability region. Furthermore, the initial radial and axial velocity perturbations can be as general as possible with no additional restrictions. In the case when $\bar{K}_{min}=K_{min}^{C}$ and $\bar{K}_{max}=K_{max}^{C}$ , this then recovers the same result as shown in Criterion 3.2.

Figure 2.

The nonlinear stability regions (when $B>0$ ) with respect to any finite-amplitude axisymmetric perturbation whose initial total circulation function satisfies $K_{min}^{C}-{\rm\Delta}K\leqslant K(x,y,t=0)<\infty$ for various levels in terms of ${\rm\Delta}K/K_{min}^{C}$ , where ${\rm\Delta}K\in [0,AR_{1}^{2}-\sqrt{ABR_{2}^{2}/3}]$ .

5.1 Solid-body rotation flow between two coaxial rotating cylinders

For a solid-body rotation flow with ${\it\omega}_{1}={\it\omega}_{2}={\it\omega}$ (i.e. ${\it\beta}=1$ ) between two coaxial rotating cylinders with radii $R_{1}$ and $R_{2}$ ( $0<{\it\alpha}=R_{1}/R_{2}<1$ ), we find from (3.8) that in this special case of ${\it\beta}=1$ we have ${\rm\Delta}\mathscr{A}_{rd}=\tilde{\mathscr{E}}$ . From (3.25), we have the identity (since $B=0$ )

(5.1) $$\begin{eqnarray}\displaystyle \frac{\text{d}\tilde{\mathscr{E}}}{\text{d}t}=-{\it\mu}\int _{V}\left(\frac{\tilde{u} _{r}^{2}}{2y}+\frac{\tilde{u} _{{\it\theta}}^{2}}{2y}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{r}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{x}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{{\it\theta}}|^{2}\right)\text{d}V. & & \displaystyle\end{eqnarray}$$

The identity (5.1) shows that any initial axisymmetric finite-amplitude perturbation around the solid-body rotation flow decays in time. By using the Sobolev inequality (Adams & Fournier Reference Adams and Fournier2003), and the boundary conditions that all the velocity perturbations vanish at the wall, we have

(5.2) $$\begin{eqnarray}\displaystyle \int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+\tilde{u} _{{\it\theta}}^{2})\,\text{d}V\leqslant M_{2}\int _{V}(|\boldsymbol{{\rm\nabla}}\tilde{u} _{r}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{x}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{{\it\theta}}|^{2})\,\text{d}V, & & \displaystyle\end{eqnarray}$$

where $M_{2}>0$ is a constant determined by the size of the cylinders. The constant $M_{2}$ is estimated from the minimization process

(5.3) $$\begin{eqnarray}\displaystyle M_{2}=\left[\min {\displaystyle \frac{\displaystyle \int _{V}|\boldsymbol{{\rm\nabla}}\tilde{u} |^{2}\,\text{d}V}{\displaystyle \int _{V}\tilde{u} ^{2}\,\text{d}V}}\right]^{-1}. & & \displaystyle\end{eqnarray}$$

In this process of minimization, $\tilde{u}$ is searched among all the admissible initial perturbation fields to minimize $M_{2}$ . As a result, we can bound the right-hand side of (5.1) and find that

(5.4) $$\begin{eqnarray}\displaystyle \frac{\text{d}\tilde{\mathscr{E}}}{\text{d}t}\leqslant -\frac{2{\it\nu}}{M_{2}}\tilde{\mathscr{E}}. & & \displaystyle\end{eqnarray}$$

This yields the following analytical upper-bound estimate for the decay of the perturbation’s kinetic energy integral $\tilde{\mathscr{E}}$ :

(5.5) $$\begin{eqnarray}\displaystyle \tilde{\mathscr{E}}(t)\leqslant \tilde{\mathscr{E}}_{0}\text{e}^{-2{\it\nu}t/M_{2}}, & & \displaystyle\end{eqnarray}$$

for all times $t\geqslant 0$ .

For two circular coaxial rotating cylinders with inner and outer wall radii $R_{1}$ and $R_{2}$ , respectively, $M_{2}$ can be found from $M_{2}=R_{2}^{2}/{\it\lambda}_{1}$ , where ${\it\lambda}_{1}$ is the first eigenvalue of the problem

(5.6) $$\begin{eqnarray}\displaystyle {\rm\Delta}w(r^{\dagger })+{\it\lambda}w(r^{\dagger })=0\quad \left({\rm\Delta}w=w_{r^{\dagger }r^{\dagger }}+\frac{w_{r^{\dagger }}}{r^{\dagger }}\right) & & \displaystyle\end{eqnarray}$$

on the annulus ${\it\alpha}=R_{1}/R_{2}<r^{\dagger }=r/R_{2}<1$ with Dirichlet conditions $w({\it\alpha})=0$ and $w(1)=0$ . The general solution of this eigenvalue problem is a linear combination of Bessel functions,

(5.7) $$\begin{eqnarray}\displaystyle w(r^{\dagger })=C_{1}\text{J}_{0}(\sqrt{{\it\lambda}}r^{\dagger })+C_{2}\text{Y}_{0}(\sqrt{{\it\lambda}}r^{\dagger }), & & \displaystyle\end{eqnarray}$$

where $\text{J}_{0}$ and $\text{Y}_{0}$ are the Bessel functions of the first and second kind of order $m=0$ . Applying the boundary conditions yields the following characteristic equation:

(5.8) $$\begin{eqnarray}\displaystyle \text{J}_{0}(\sqrt{{\it\lambda}})\text{Y}_{0}(\sqrt{{\it\lambda}}{\it\alpha})-\text{Y}_{0}(\sqrt{{\it\lambda}})\text{J}_{0}(\sqrt{{\it\lambda}}{\it\alpha})=0. & & \displaystyle\end{eqnarray}$$

The solution of the first eigenvalue ${\it\lambda}_{1}$ of (5.6) (or the first zero of (5.8)) for various values of ${\it\alpha}$ in the range between $0.1$ and $0.8$ is shown in figure 5. It can be seen that ${\it\lambda}_{1}$ increases with ${\it\alpha}$ . The analytical upper bound for the perturbation’s decay rate is then estimated by ${\it\sigma}_{E}^{\ast }=-2{\it\lambda}_{1}{\it\nu}/R_{2}^{2}$ . This upper bound (for the decay rate) is independent of the angular speed of rotation ${\it\omega}$ . Rescaling the time such that $t^{\ast }={\it\omega}t$ and by using $R_{o}={\it\omega}R_{2}^{2}/{\it\nu}$ , we can also rewrite it as ${\it\sigma}_{E}^{\ast }=-2{\it\lambda}_{1}/R_{o}$ . We find that, as the gap between the cylinders becomes smaller ( ${\it\alpha}$ increases towards $1$ ), the flow perturbation decays faster in time, while, as the Reynolds number of the outer cylinder $R_{o}$ increases, the perturbation decays slower in time. In addition, for a fixed rate of rotation of the cylinders (including the case of no rotation) or for cylinders with an accelerating rotational speed from zero up to a given value and then rotation is stopped, the perturbation’s decay rate is approximately the same as this analytical upper-bound estimate as long as the flow stays axisymmetric.

Figure 5.

The dependence of ${\it\lambda}_{1}$ on ${\it\alpha}$ .

To further elaborate on the perturbation’s dynamics, we focus on the example where ${\it\alpha}=0.5$ and then ${\it\lambda}_{1}=39.0138$ . This leads to

(5.9) $$\begin{eqnarray}\displaystyle \tilde{\mathscr{E}}(t)\leqslant \tilde{\mathscr{E}}_{0}\text{e}^{-2{\it\nu}t/M_{2}}\approx \tilde{\mathscr{E}}_{0}\text{e}^{-78.0276t^{\ast }/R_{o}}. & & \displaystyle\end{eqnarray}$$

In this case, the upper bound for the decay rate of the perturbation’s kinetic energy integral $\tilde{\mathscr{E}}$ is estimated by ${\it\sigma}_{E}^{\ast }=-78.0276/R_{o}$ .

We also conducted numerical simulations of the flow dynamics based on (2.3), (2.4a ) and (2.4b ) and the assumed boundary conditions. These simulations are based on a modified code of Rusak et al. (Reference Rusak, Wang, Xu and Taylor2012) and Rusak, Granata & Wang (Reference Rusak, Granata and Wang2015), which includes the viscous terms and axial periodic boundary conditions. For the case of a base flow with ${\it\beta}=1$ and $0<{\it\alpha}<1$ , we used various initial conditions given by

(5.10) $$\begin{eqnarray}\displaystyle \left.\begin{array}{@{}c@{}}\displaystyle {\it\psi}(x,y,0)={\it\epsilon}_{{\it\psi}}\cos \left(\frac{2{\rm\pi}nx}{L}\right)f(y),\\ \displaystyle {\it\chi}(x,y,0)=-{\it\epsilon}_{{\it\chi}}\cos \left(\frac{2{\rm\pi}nx}{L}\right)\left[f_{yy}(y)-\left(\frac{2{\rm\pi}n}{L}\right)^{2}\frac{f(y)}{2y}\right],\\ \displaystyle K(x,y,0)=2{\it\omega}y-{\it\epsilon}_{K}\cos \left(\frac{2{\rm\pi}nx}{L}\right)f(y),\end{array}\right\} & & \displaystyle\end{eqnarray}$$

for $0\leqslant x<L$ and $y_{1}\leqslant y\leqslant y_{2}$ , where ${\it\epsilon}_{{\it\psi}},{\it\epsilon}_{{\it\chi}}$ and ${\it\epsilon}_{K}$ are the perturbation’s amplitudes. The constant $n$ is an integer, and $f(y)=(y-y_{1})^{2}(y_{2}-y)^{2}/c_{m}$ ( $c_{m}$ is a scaling factor such that $\max f(y)=1$ ).

The numerical simulations of the evolution for the non-dimensional perturbation’s kinetic energy integral $\tilde{\mathscr{E}}^{non}=\tilde{\mathscr{E}}/({\it\rho}{\it\omega}_{2}^{2}R_{2}^{2})$ with $R_{o}=100$ , ${\it\alpha}=0.5$ , ${\it\beta}=1$ , ${\it\omega}=1$ , $L=2$ , $R_{2}=1$ and various initial conditions according to (5.10) (with $c_{m}=0.001236$ ) are presented in figure 6. The initial non-dimensional perturbation’s kinetic energy integral $\tilde{\mathscr{E}}^{non}|_{t=0}=0.000242$ is the same for all cases shown. Also shown for reference is the analytic upper bound for the decay of the perturbation’s kinetic energy integral $\tilde{\mathscr{E}}$ according to the right-hand side of (5.9) (the heavy solid line). It can be seen that the case with ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.01236$ and $n=0$ shows a decay of $\tilde{\mathscr{E}}$ (the curve using open circles) that is only slightly below this predicted analytical upper bound (the heavy solid line). In all other cases with various ${\it\epsilon}_{{\it\psi}}$ , ${\it\epsilon}_{{\it\chi}}$ , ${\it\epsilon}_{K}$ and $n$ , $\tilde{\mathscr{E}}$ decays faster than this analytic upper-bound prediction. The decay process for each of these perturbations involves a stage of nonlinear interaction between the flow perturbations that enhances the flow mixing and increases the decay of the perturbations. The numerical simulations demonstrate that this analytic upper-bound decay line provides a simple and relevant criterion for all the processes.

Figure 6.

The computed time history of the non-dimensional perturbation’s kinetic energy integral $\tilde{\mathscr{E}}^{non}=\tilde{\mathscr{E}}/({\it\rho}{\it\omega}_{2}^{2}R_{2}^{2})$ for $R_{o}=100$ , ${\it\alpha}=0.5$ , ${\it\beta}=1$ , ${\it\omega}=1$ , $L=2$ , $R_{2}=1$ for various initial conditions. All perturbations start with the same initial value of $\tilde{\mathscr{E}}^{non}|_{t=0}=0.000242$ . The heavy solid line is the analytic upper bound for the decay of the perturbation’s kinetic energy integral $\tilde{\mathscr{E}}$ according to (5.9). The other curves are for the following cases: open circles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.01236$ and $n=0$ ; open squares, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.01748$ and $n=1$ ; open diamonds, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.01748$ and $n=2$ ; open downward-pointing triangles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00058$ , ${\it\epsilon}_{K}=0.01164$ and $n=0$ ; open upward-pointing triangles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00082$ , ${\it\epsilon}_{K}=0.01626$ and $n=1$ ; open left-pointing triangles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00223$ , ${\it\epsilon}_{K}=0$ and $n=1$ ; dotted line, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00172$ , ${\it\epsilon}_{K}=0$ and $n=0$ .

5.2 Circular Couette flow

We analyse the circular Couette flow for various parameters $0<{\it\alpha}<1$ and $0<{\it\beta}<1$ . The analytical upper-bound estimate for the decay rate of ${\rm\Delta}\mathscr{A}_{rd}$ is derived directly from (3.25). To this end, we first evaluate the term on the left-hand side of (3.46) when $K_{min}=K_{min}^{C}$ , and then rewrite it in terms of the non-dimensional parameters ${\it\alpha}$ and ${\it\beta}$ . This yields

(5.11) $$\begin{eqnarray}\displaystyle d({\it\alpha},{\it\beta})=1-\frac{ABR_{2}^{2}}{3(K_{min}-B)^{2}}=1-\frac{(1-{\it\beta})}{3({\it\beta}-{\it\alpha}^{2}){\it\alpha}^{2}}. & & \displaystyle\end{eqnarray}$$

Thus from (3.25), we have the following upper bound for the temporal rate of change of ${\rm\Delta}\mathscr{A}_{rd}$ :

(5.12) $$\begin{eqnarray}\displaystyle \frac{\text{d}{\rm\Delta}\mathscr{A}_{rd}}{\text{d}t}\leqslant -{\it\mu}\int _{V}(|\boldsymbol{{\rm\nabla}}\tilde{u} _{r}|^{2}+|\boldsymbol{{\rm\nabla}}\tilde{u} _{x}|^{2}+d({\it\alpha},{\it\beta})|\boldsymbol{{\rm\nabla}}\tilde{u} _{{\it\theta}}|^{2})\,\text{d}V, & & \displaystyle\end{eqnarray}$$

which is valid for any initial finite-amplitude axisymmetric perturbation such that $K_{min}^{C}\leqslant K(x,y,t=0)\leqslant K_{max}^{C}$ .

Using the Sobolev inequality (Adams & Fournier Reference Adams and Fournier2003), this leads to

(5.13) $$\begin{eqnarray}\displaystyle \frac{\text{d}{\rm\Delta}\mathscr{A}_{rd}}{\text{d}t}\leqslant -\frac{{\it\mu}}{M_{2}}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+d({\it\alpha},{\it\beta})\tilde{u} _{{\it\theta}}^{2})\,\text{d}V. & & \displaystyle\end{eqnarray}$$

The constant $M_{2}$ can be found again from $M_{2}=R_{2}^{2}/{\it\lambda}_{1}$ , where ${\it\lambda}_{1}$ is the same exact first eigenvalue found in § 5.1 and is a function of ${\it\alpha}$ alone (see figure 5). On the other hand, we also have from (3.8) and (3.10) that

(5.14) $$\begin{eqnarray}\displaystyle 0<{\rm\Delta}\mathscr{A}_{rd}\leqslant \frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+{\it\lambda}^{+}({\it\alpha},{\it\beta})\tilde{u} _{{\it\theta}}^{2})\,\text{d}V. & & \displaystyle\end{eqnarray}$$

Combining (5.13) and (5.14) gives

(5.15) $$\begin{eqnarray}\displaystyle \frac{\text{d}{\rm\Delta}\mathscr{A}_{rd}}{\text{d}t}\leqslant -\frac{d({\it\alpha},{\it\beta})}{{\it\lambda}^{+}({\it\alpha},{\it\beta})}\frac{2{\it\nu}}{M_{2}}{\rm\Delta}\mathscr{A}_{rd}. & & \displaystyle\end{eqnarray}$$

Thus, we have

(5.16) $$\begin{eqnarray}\displaystyle {\rm\Delta}\mathscr{A}_{rd}(t)\leqslant {\rm\Delta}\mathscr{A}_{rd0}\exp \left(-\frac{d({\it\alpha},{\it\beta})}{{\it\lambda}^{+}({\it\alpha},{\it\beta})}\frac{2{\it\nu}}{M_{2}}t\right) & & \displaystyle\end{eqnarray}$$

and

(5.17) $$\begin{eqnarray}\displaystyle \Vert \tilde{\boldsymbol{u}}(t)\Vert _{{\it\lambda}^{-}}^{2} & := & \displaystyle \frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+{\it\lambda}^{-}({\it\alpha},{\it\beta})\tilde{u} _{{\it\theta}}^{2})\,\text{d}V\nonumber\\ \displaystyle & {\leqslant} & \displaystyle \exp \left(-\frac{d({\it\alpha},{\it\beta})}{{\it\lambda}^{+}({\it\alpha},{\it\beta})}\frac{2{\it\nu}}{M_{2}}t\right)\frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r0}^{2}+\tilde{u} _{x0}^{2}+{\it\lambda}^{+}({\it\alpha},{\it\beta})\tilde{u} _{{\it\theta}0}^{2})\,\text{d}V\nonumber\\ \displaystyle & = & \displaystyle \exp \left(-\frac{d({\it\alpha},{\it\beta})}{{\it\lambda}^{+}({\it\alpha},{\it\beta})}\frac{2{\it\nu}}{M_{2}}t\right)\Vert \tilde{\boldsymbol{u}}(t=0)\Vert _{{\it\lambda}^{+}}^{2}.\end{eqnarray}$$

We find that the perturbation $\tilde{\boldsymbol{u}}=(\tilde{u} _{r},\tilde{u} _{x},\tilde{u} _{{\it\theta}})$ at a later time in terms of the stability norm $\Vert \cdot \Vert _{{\it\lambda}^{-}}^{2}$ given by the perturbation’s weighted kinetic energy integral (as defined by the left-hand side of (5.17)) is bounded from above by the initial perturbation $\tilde{\boldsymbol{u}}_{0}=(\tilde{u} _{r0},\tilde{u} _{x0},\tilde{u} _{{\it\theta}0})$ in terms of the stability norm $\Vert \cdot \Vert _{{\it\lambda}^{+}}^{2}$ given by the perturbation’s weighted kinetic energy integral (as defined by the right-hand side of (5.17)). The upper bound for the decay rate of the perturbation’s weighted kinetic energy integral is estimated analytically by ${\it\sigma}_{WE}^{\ast }=-2{\it\lambda}_{1}({\it\alpha})(d({\it\alpha},{\it\beta})/{\it\lambda}^{+}({\it\alpha},{\it\beta}))({\it\nu}/R_{2}^{2})$ . Rescaling the time $t$ such that $t^{\ast }={\it\omega}_{2}t$ and using the Reynolds number based on the outer cylinder $R_{o}={\it\omega}_{2}R_{2}^{2}/{\it\nu}$ , we can write ${\it\sigma}_{WE}^{\ast }=-2{\it\lambda}_{1}({\it\alpha})(d({\it\alpha},{\it\beta})/{\it\lambda}^{+}({\it\alpha},{\it\beta}))(1/R_{o})$ . This upper bound depends on ${\it\alpha}$ (the geometry of the cylinders), ${\it\beta}$ (the ratio of the rates of rotation of the cylinders) and $R_{o}$ (the Reynolds number of the outer cylinder).

Figure 7.

The computed time history of the non-dimensional perturbation’s weighted kinetic energy integral $\tilde{\mathscr{E}}_{{\it\lambda}^{-}}^{non}=\Vert \tilde{\boldsymbol{u}}(t)\Vert _{{\it\lambda}^{-}}^{2}/({\it\rho}{\it\omega}_{2}^{2}R_{2}^{2})$ for $R_{o}=100$ , ${\it\alpha}=0.8$ , ${\it\beta}=0.96$ , ${\it\omega}_{2}=1$ , $L=2$ , $R_{2}=1$ and various initial conditions. All perturbations start with the same initial value, which is given by the initial weighted kinetic energy integral $\tilde{\mathscr{E}}_{{\it\lambda}^{+}}^{non}|_{t=0}=0.00021$ (according to the right-hand side of (5.18)). The heavy solid line is the analytic upper bound for the decay of perturbation according to (5.18). The other curves are for the following cases: open circles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.01755$ and $n=0$ ; open squares, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.02482$ and $n=1$ ; open diamonds, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.02482$ and $n=2$ ; open downward-pointing triangles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00067$ , ${\it\epsilon}_{K}=0.01410$ and $n=0$ ; open upward-pointing triangles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00092$ , ${\it\epsilon}_{K}=0.02000$ and $n=1$ ; open left-pointing triangles, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00155$ , ${\it\epsilon}_{K}=0$ and $n=1$ ; dotted line, ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0.00111$ , ${\it\epsilon}_{K}=0$ and $n=0$ .

For the case with ${\it\alpha}=R_{1}/R_{2}=0.8$ and ${\it\beta}={\it\omega}_{2}/{\it\omega}_{1}=0.96$ , the computation gives that ${\it\lambda}_{1}=246.4298$ , ${\it\lambda}^{-}=1.0512$ , ${\it\lambda}^{+}=1.1953$ and $d=0.9349$ . Then, (5.17) becomes

(5.18) $$\begin{eqnarray}\displaystyle \frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r}^{2}+\tilde{u} _{x}^{2}+1.0512\tilde{u} _{{\it\theta}}^{2})\,\text{d}V\leqslant \text{e}^{-385.49t^{\ast }/R_{o}}\frac{1}{2}{\it\rho}\int _{V}(\tilde{u} _{r0}^{2}+\tilde{u} _{x0}^{2}+1.1953\tilde{u} _{{\it\theta}0}^{2})\,\text{d}V.\quad & & \displaystyle\end{eqnarray}$$

The upper bound for the decay rate of the perturbation’s weighted kinetic energy integral is estimated analytically by ${\it\sigma}_{WE}^{\ast }=-385.49/R_{o}$ .

Numerical simulation results of the evolution of the non-dimensional perturbation’s weighted kinetic energy integral $\tilde{\mathscr{E}}_{{\it\lambda}^{-}}^{non}=\Vert \tilde{\boldsymbol{u}}(t)\Vert _{{\it\lambda}^{-}}^{2}/({\it\rho}{\it\omega}_{2}^{2}R_{2}^{2})$ according to the left-hand side of (5.18) for $R_{o}=100$ , ${\it\alpha}=0.8$ , ${\it\beta}=0.96$ , ${\it\omega}_{2}=1$ , $L=2$ , $R_{2}=1$ and various initial conditions similar to those in (5.10) (with $c_{m}=0.000066$ ) are presented in figure 7. The initial non-dimensional perturbation’s weighted kinetic energy integral $\tilde{\mathscr{E}}_{{\it\lambda}^{+}}^{non}|_{t=0}$ according to the right-hand side of (5.18) is the same for all cases shown, which equals $0.00021$ . Also shown for reference is the analytic upper-bound line for the decay of the perturbation (the heavy solid line) according to the right-hand side of (5.18). It can be seen that the case with ${\it\epsilon}_{{\it\psi}}={\it\epsilon}_{{\it\chi}}=0$ , ${\it\epsilon}_{K}=0.01755$ and $n=0$ shows a decay of the non-dimensional perturbation’s weighted kinetic energy integral (the curve with open circles) that is below this predicted analytic upper bound (the heavy solid line). In all other cases studied, the non-dimensional perturbation’s weighted kinetic energy integral decays faster than this upper-bound prediction. Again, each of these decay processes involves a stage of nonlinear interaction between the flow perturbations that enhances the flow mixing and increases the decay of the perturbations. The numerical simulations demonstrate that this analytic upper-bound decay line provides a simple and relevant criterion for all the processes.