# Revisiting the linear instabilities of plane channel flow between compliant walls

Smail Lebbal, Frédéric Alizard, Benoît Pier

► To cite this version:

Smail Lebbal, Frédéric Alizard, Benoît Pier. Revisiting the linear instabilities of plane channel flow between compliant walls. Physical Review Fluids, 2022, 7 (2), pp.023903. 10.1103/PhysRevFluids.7.023903. hal-03337441v3

HAL Id: hal-03337441

https://hal.science/hal-03337441v3

Submitted on 11 Feb 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

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

Distributed under a Creative Commons Attribution 4.0 International License
Revisiting the linear instabilities of plane channel
flow between compliant walls

Smail Lebbal,¹ Frédéric Alizard,²,\* and Benoît Pier¹

¹Univ Lyon, Ecole Centrale de Lyon, CNRS, Univ Claude Bernard Lyon 1, INSA Lyon, LMFA, UMR5509,
69130 Ecully, France

²Univ Lyon, Univ Claude Bernard Lyon 1, CNRS, Ecole Centrale de Lyon, INSA Lyon, LMFA, UMR5509,
69622 Villeurbanne France

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

(Received 4 September 2021; accepted 19 January 2022; published 11 February 2022)

The present investigation revisits the linear stability of Poiseuille channel flow interacting with compliant walls. The results obtained include the dynamics of Tollmien-Schlichting (TS) modes as well as flow-induced surface-instability (FSI) modes, in the form of both traveling-wave flutter (TWF) and divergence modes. The compliant wall model consists of a spring-backed plate with a viscous substrate deformable in the vertical direction [C. Davies and P. W. Carpenter, *J. Fluid Mech.* **352**, 205 (1997)]. At the interface between the fluid and the walls, the continuity of velocities and stresses, including both viscous and pressure contributions, are taken into account. The FSI modes (both varicose and sinuous) and TS modes are then reinterpreted in the light of the two principal nondimensional control parameters: the Reynolds number (Re), which characterizes the base flow, and the reduced velocity ($V_R$), which measures the response of the flexible wall to hydrodynamic loading [E. De Langre, La Houille Blanche, 3, 14 (2000)]. The characteristics of TS and FSI modes are systematically investigated over a large control-parameter space, including wall dissipation, spring stiffness, and flexural rigidity. We observe that TWF modes are primarily governed by $V_R$ and largely independent of the Reynolds number. It is found that the instability is generally dominated by the TWF mode of varicose symmetry. Divergence and TS modes are both affected by $V_R$ and Re, confirming that these modes belong to a different class. The onset of the divergence mode is mainly observed for the sinuous motion, when increasing the dissipation. To provide physical insight into the mechanisms driving these instabilities, the perturbative energy equations for both FSI and TS modes are analyzed for a wide range of wall parameters and wavenumbers.

DOI: 10.1103/PhysRevFluids.7.023903

I. INTRODUCTION

The constant scientific interest to extend the laminar regime for industrial applications has led to the development of compliant walls since the beginning of the 20th century. In particular, researchers focused on finding optimum wall properties, aiming to delay the laminar-turbulent transition.

In the biological context, interactions between fluid and elastic forces associated with a deformable channel or tube lead to a variety of physiologically significant phenomena. In particular, deformability plays a prominent role in blood flow as well as peristaltic transport, for example, through the intestines and the urogenital tract (see [1,2] for a review).

Such interest arose from Gray's paradox [3]. Gray showed that, to overcome the friction drag of a swimming dolphin subjected to a turbulent flow around its body, the muscles have to be capable

<footnote>*frederic.alizard@univ-lyon1.fr</footnote>
of generating a power at least seven times greater than that of other kinds of mammalian muscle. Hence, Gray has suggested that the dolphin skin is able to delay the laminar-turbulent transition. In the laminar regime, indeed, the power developed by muscles would still conform to that of other types of mammalian muscle.

Kramer [4] in the 1950s developed a compliant coating trying to mimic the dolphin's skin. The author claimed that he was able to reduce the drag of a torpedolike model by as much as 60%. Later, the conclusions of both Gray's hypotheses and Kramer's experiments were in part questioned. Russian and American experiments since the 1980s failed to reproduce results provided by Kramer. In addition, scientists observed that the lower turbulence level around the dolphin's swimming body could also be attributed to local pressure gradients. It was also suggested that a reduction in the friction drag may result from the fact that the dolphin leaps out of the water for breathing (see [5–7] for reviews).

However, Gaster's experiments in 1988 [8] gave new hope in using compliant walls to delay laminar-turbulent transition. He showed that the growth rates of artificially generated Tollmien-Schlichting (TS) waves are inhibited when using appropriate coatings. These experiments have given a strong impulse to theoretical developments aiming to tackle this problem.

A major difficulty arises from the design of compliant wall models that are able to couple fluid and solid dynamics. These models may be separated into two categories: surface based and volume based (see [7,9] for a review). The first class of models is less computationally demanding and considers an infinitely thin wall interacting with the fluid through an interface condition. In this case, the wall is defined as a thin plate mounted on springs and dampers. The wall parameters are classically the spring stiffness, the tension, the bending stiffness, its mass, and the damping coefficient.

For the second family of models (i.e., volume based), the wall material is fully described to include single- or multilayer coatings (see [10] and [9] for recent reviews). However, Duncan [11] has shown that a surface model can be used to quantitatively describe many aspects of instabilities and wave propagation on the surface of an elastic and incompressible coating.

Here, the surface-based approach is adopted. Within this approach and following the lead of Benjamin [12–14], which relies on the theory developed by Miles for water waves [15], Carpenter and Garrad [16,17] focused on the stability of boundary-layer flows over Kramer-type compliant walls. They provided some confirmation of the transition-delaying potential of compliant coatings. According to Carpenter and Garrad [16,17], instability modes can be classified into two categories: fluid based (TS mode) and solid based (flow-induced surface instabilities, or fluid-structure instabilities, referenced as FSI hereafter). The last category includes both the traveling-wave flutter (TWF) modes and the (almost static) divergence modes. For the divergence mode, scientists are still arguing about its precise nature. It is either interpreted as an absolute instability [18] or it may also result from a modal instability with a nearly vanishing phase velocity when increasing the wall dissipation [19].

While the surface-based model has been extensively used for external flows such as boundary layers, very few studies have been conducted on internal flows (see, for instance, the recent review of Kumaran [9]). Among them, Davies and Carpenter [20,21] investigated linear instability waves that emerge when a viscous incompressible flow interacts with infinite or finite spring-backed plates in a plane channel.

Considering only sinuous perturbations, Davies and Carpenter [20] derived a theory for the motion of the walls and obtained neutral stability curves for both FSI and TS modes and some values of flexural rigidity of the plate, spring stiffness, and damping. The same year, these authors investigated the dynamics of TS waves that propagate over a compliant panel of finite length, using linearized numerical simulations [21]. For low external disturbances, Davies and Carpenter [21] then suggested that the use of multiple-panel compliant walls could maintain the laminar flow regime at all Reynolds numbers in a plane channel flow. This seminal work was also limited to the sinuous symmetry. These prior analyses have been extended to the varicose motion of the walls for the same case with infinite compliant walls by Nagata and Cole [22]. Despite the fact that this
investigation focused only on a narrow region of the parameter space to make definite statements, the authors found that the varicose TWF mode can be amplified prior to its sinuous counterpart. A similar analysis was recently conducted to address the effect of pulsatile Poiseuille flow through a compliant channel [23] for modal and nonmodal instabilities. Although the latter extended results discussed above to pulsatile flows and short-time perturbations, it has not brought new elements for the asymptotic regime of the steady flow case. Finally, Gajjar and Sibanda [24] investigated a viscous flow in a channel in which only one wall is compliant. The same spring-plate model is used by the authors and the nonlinear evolution of disturbances have been obtained in the limit of large Reynolds numbers within the nonlinear critical layer theory.

The aforementioned analyses address the amplification of FSI and TS modes for compliant channel flows but it should be noted that the models used above have some limitations. First, these analyses use the Reynolds number as a control parameter for both FSI and TS modes. This approach is not entirely satisfactory because a change in the Reynolds number also leads to a change in the compliant wall characteristics. Therefore, the reduced velocity, defined as the ratio between the characteristic time of the solid and that of the fluid, has been advocated by de Langre [25] to measure the strength of the fluid-structure interactions independently of the Reynolds number. Second, a rigorous treatment of the pressure at the wall is key for a complete analysis of fluid-structure interaction problems. On one hand, a consistent formulation was indeed derived by Davies and Carpenter [20] under the form of an Orr-Sommerfeld-like equation accounting for compliant walls. However, this formulation only applies to perturbations of sinuous symmetry. On the other hand, Nagata and Cole [22] and Tsigklifis and Lucey [23] used a primitive formulation combined with a spectral collocation method without any specific treatments for the pressure and the additional boundary conditions arising from the wall equation. In this case, spurious eigenvalues are likely to occur [26]. Third, as suggested by Nagata and Cole [22], the varicose symmetry is as important as the sinuous one for FSI modes. Nevertheless, as far as we know, existing studies of varicose modes only explore a narrow parameter range. As a consequence, the potential gain in using compliant walls to delay laminar-turbulent transition in a plane channel flow could be challenged significantly.

The purpose of the present paper is therefore to revisit the instability modes that emerge when a viscous flow interacts with two parallel spring-backed plates and hope to address the inherent limitations mentioned above. Thus, we reconsider this problem by using the reduced velocity to describe FSI modes. A general formulation is implemented, free of spurious modes, suitable for both sinuous and varicose symmetries, and taking into account the exact hydrodynamic forces acting on the walls. A wide range of wall parameters has been explored to highlight their influence on both FSI and TS modes. Only two-dimensional perturbations are considered since the Squire theorem holds for compliant walls [27].

The paper is organized as follows. Section II presents the model and governing equations used for the fluid-structure interaction problem. In particular, the dimensionless parameters and linearized equations will be introduced. After the numerical methods are presented in Sec. III, Sec. IV provides some physical insight into the influence of wall parameters onto FSI and TS modes for both the sinuous and varicose symmetries. For that purpose, kinetic energy budgets are also computed. Conclusions and prospects are given in the last section.

## II. PROBLEM FORMULATION

### A. Fluid-structure interaction model and interface conditions

In the entire paper, we restrict our analysis to the two-dimensional problem. Using a Cartesian coordinate system (x, y) with unit vectors (e<sub>x</sub>, e<sub>y</sub>), we consider an incompressible Newtonian fluid with dynamic viscosity μ and density ρ between two spring-backed deformable plates located at y = ξ<sup>±</sup>(x, t), which are allowed to move only in the y direction. A schematic diagram of the configuration is shown in Fig. 1.
<figure><img src="image_3.png" /><figcaption>FIG. 1. Channel flow with infinite spring-backed flexible walls.</figcaption></figure>

The flow between the walls is governed by the Navier-Stokes equations:

$$ \frac{\partial \mathbf{u}}{\partial t} + \mathbf{u} \cdot \nabla \mathbf{u} = -\frac{1}{\rho} \nabla p + \nu \Delta \mathbf{u}, \\ \nabla \cdot \mathbf{u} = 0, \tag{1} $$

where $\mathbf{u} = (u, v)$ and $p$ represent the velocity and pressure fields, $u$ ($v$) denotes the streamwise (wall-normal) velocity component, and $\nu \equiv \mu/\rho$ is the kinematic viscosity of the fluid.

The movement of the flexible plates obeys the following equations, derived through Newton's second law:

$$ m \frac{\partial^2 \zeta^\pm}{\partial t^2} + d \frac{\partial \zeta^\pm}{\partial t} + \left( B \frac{\partial^4}{\partial x^4} - T \frac{\partial^2}{\partial x^2} + K \right) \zeta^\pm = f^\pm. \tag{2} $$

Here, $m$ denotes the mass per unit area of the plates, $d$ their damping coefficient, $B$ the flexural rigidity, $T$ the wall tension, $K$ the spring stiffness, and $f^\pm$ represents the $y$ component of the forces acting on the plates. These forces are obtained as

$$ f^{\pm} = e_y \cdot f^{\pm} \quad \text{with} \quad f^{\pm} = (\bar{\tau}^{\pm} - \delta p^{\pm} \mathbf{I}) \cdot \mathbf{n}^{\pm}. \tag{3} $$

Here $\bar{\tau}^{\pm}$ denotes the viscous stress tensor at the walls and $\delta p^{\pm}$ the transmural surface pressure, i.e., the difference between the surface pressure inside and outside of the channel, and $\mathbf{n}^{\pm}$ is the unit vector normal to the walls pointing towards the fluid.

The viscous stress tensor at the walls has the following expression:

$$ \bar{\tau}^{\pm} = \begin{pmatrix} 2\mu \frac{\partial u}{\partial x} & \mu \left( \frac{\partial u}{\partial y} + \frac{\partial v}{\partial x} \right) \\ \mu \left( \frac{\partial u}{\partial y} + \frac{\partial v}{\partial x} \right) & 2\mu \frac{\partial v}{\partial y} \end{pmatrix} \Big|_{y=\zeta^{\pm}}, \tag{4} $$
and the normal vectors to the interface $\mathbf{n}^{\pm} = (n_x^{\pm}, n_y^{\pm})$ are obtained as

$$n_x^{\pm} = \pm \frac{\partial \zeta^{\pm}}{\partial x} \frac{1}{\sqrt{1 + \left(\frac{\partial \zeta^{\pm}}{\partial x}\right)^2}} \quad \text{and} \quad n_y^{\pm} = \mp \frac{1}{\sqrt{1 + \left(\frac{\partial \zeta^{\pm}}{\partial y}\right)^2}} \quad (5)$$

Combining Eqs. (3)–(5) yields the y component of the hydrodynamic forces acting on the compliant walls,

$$f^{\pm} = \left[ \pm \mu \left( \left. \frac{\partial u}{\partial y} \right|_{y=\zeta^{\pm}} + \left. \frac{\partial v}{\partial x} \right|_{y=\zeta^{\pm}} \right) \zeta_x^{\pm} \mp 2\mu \left. \frac{\partial v}{\partial y} \right|_{y=\zeta^{\pm}} \pm \delta p^{\pm} \right] / \sqrt{1 + \left( \left. \frac{\partial \zeta^{\pm}}{\partial x} \right) \right)^2, \quad (6)$$

which governs the wall dynamics (2) since the wall movement is constrained to occur only in the y direction. Note that in this approach there are no hydrodynamic forces acting on the plate from the outside, except for a pressure.

Finally, the no-slip conditions on both walls lead to the kinematic conditions prevailing at the moving boundaries:

$$u = 0 \quad \text{and} \quad v = \frac{\partial \zeta^{\pm}}{\partial t} \quad \text{for} \quad y = \zeta^{\pm}. \quad (7)$$

Thus, the fluid-structure interaction problem is completely determined by the coupling of the fluid equations (1), the wall equations (2), and boundary conditions (7).

## B. Dimensionless control parameters

The present compliant-channel flow configuration is characterized by nine dimensional parameters: the volumetric flow rate $Q$ (m$^3$ s$^{-1}$); the half height $h$ (m) of the channel; the fluid density $\rho$ (kg m$^{-3}$); the kinematic viscosity $\nu$ (m$^2$ s$^{-1}$); the mass of the plate per unit area, $m$ (kg m$^{-2}$); the damping coefficient of the wall, $d$ (kg m$^{-2}$ s$^{-1}$); the bending stiffness of the plate, $B$ (kg m$^2$ s$^{-2}$); the wall tension $T$ (kg s$^{-2}$); and the spring stiffness $K$ (kg m$^{-2}$ s$^{-2}$). Hence, this system may be described by six dimensionless parameters. A useful parameter to characterize fluid-structure interaction phenomena is the reduced flow velocity $V_R$, defined [25] as the ratio of a characteristic time of the structure to a characteristic time of the flow. Using time scales based on spring stiffness,

$$\tau_K = \sqrt{\frac{m}{K}},$$

and flow advection,

$$\tau_Q = \frac{4h^2}{Q},$$

the reduced velocity is obtained as $V_R = \tau_K / \tau_Q$. Other choices based on different characteristic times would be possible. For $V_R \ll 1$, the influence of the wall compliance is negligible, while $V_R \gg 1$ corresponds to very soft walls. Hence, the resulting six nondimensional control parameters are

$$\left\{ \begin{array}{lll} \text{Re} = \frac{Q}{\nu}, & V_R = \frac{Q}{4h^2} \sqrt{\frac{m}{K}}, & \Gamma = \frac{m}{\rho h}, \\ d_* = \frac{d}{\sqrt{mK}}, & B_* = \frac{B}{Kh^4}, & T_* = \frac{T}{h^2 K}. \end{array} \right. \quad (8)$$

Here Re is the Reynolds number based on channel diameter and average flow velocity, and $\Gamma$ is the mass ratio between the compliant walls and the fluid. The three nondimensional wall parameters $d_*$, $B_*$, and $T_*$ are all relative to the spring stiffness $K$, which serves as reference for the reduced velocity $V_R$. One may notice that in several previous studies [20,22,28] fluid quantities are used to
build nondimensional parameters for the wall. As underlined by Domaradzki and Metcalfe [29], this may correspond to nonphysical situations where at each Reynolds number a different compliant wall and a different fluid are considered. In the present work, the Reynolds number may be modified using $v$ and/or $Q$ without changing the wall properties. Hereafter, the mass ratio is kept constant at $\Gamma = 2$ and we only consider walls without tension, $T = 0$. We fix the three dimensional parameters at $\rho = 1$, $h = 1$, and $\tau_Q = 1$.

### C. Formulation of the linearized model

This entire study considers the dynamics of small-amplitude perturbations, obeying the linearized version of the governing equations around a steady base state. The unperturbed base configuration thus consists of a parabolic Poiseuille flow $U(y) = (U_b(y), 0)$ between parallel walls located at $y = \pm h$ [see Fig. 1(a)]. This flow is driven by a pressure $P_b(x) = P_0 - Gx$ of constant streamwise gradient, and we assume a pressure outside the channel walls always equal to $P_b(x)$, so as to equilibrate the forces acting on both sides of the walls for unperturbed conditions.

The total flow fields are then decomposed as $\mathbf{u}(x, y, t) = \mathbf{U}(y) + \mathbf{u}'(x, y, t)$ and $p(x, y, t) = P_b(x) + p'(x, y, t)$. The wall displacement is written as $\zeta^{\pm}(x, t) = \pm h + \eta^{\pm}(x, t)$.

Considering that the perturbation components $\mathbf{u}', p'$, and $\eta^{\pm}$ are of small amplitude, the governing equations may be linearized about the base state.

Since the base state is steady and homogeneous in the streamwise direction, the perturbation to the velocity fields, pressure fields, and normal displacements are expressed in normal-mode form as

$$ \mathbf{u}'(x, y, t) = \hat{\mathbf{u}}(y)e^{i(\alpha x - \omega t)}, \quad p'(x, y, t) = \hat{p}(y)e^{i(\alpha x - \omega t)}, \quad \text{and} \quad \eta^{\pm}(x, t) = \hat{\eta}^{\pm}e^{i(\alpha x - \omega t)}, \quad (9) $$

with $\alpha$ the streamwise wavenumber and $\omega$ the frequency. Hereafter, we adopt a temporal viewpoint where $\alpha \in \mathbb{R}$ and $\omega = \omega_r + i\omega_i \in \mathbb{C}$, with $\omega_i$ the temporal amplification rate of the mode and $\omega_r$ its circular frequency. Substitution of this decomposition into the Navier-Stokes equations and linearization about the base flow leads to

$$ \begin{cases} -i\omega\rho\hat{\mathbf{u}} + \rho i\alpha\hat{\mathbf{u}}\mathbf{U}_b + \rho\hat{v}\frac{d\mathbf{U}_b}{dy} = -i\alpha\hat{p} + \mu\left(\frac{d^2}{dy^2} - \alpha^2\right)\hat{u}, \\ -i\omega\rho\hat{v} + \rho i\alpha\hat{v}\mathbf{U}_b = -\frac{d\hat{p}}{dy} + \mu\left(\frac{d^2}{dy^2} - \alpha^2\right)\hat{v}, \end{cases} \quad (10) $$

together with the divergence-free condition

$$ i\alpha\hat{\mathbf{u}} + \frac{d\hat{v}}{dy} = 0. \quad (11) $$

Linearization of the wall equations (2) and (6) yields

$$ \begin{cases} -\omega^2 m \hat{\eta}^+ - i\omega d \hat{\eta}^+ + (B\alpha^4 + T\alpha^2 + K) \hat{\eta}^+ = +\hat{p}(h) + \mu \left( i\alpha \left. \frac{dU_b}{dy} \right|_h \hat{\eta}^+ \right) - 2\mu \left. \frac{d\hat{v}}{dy} \right|_{+h}, \\ -\omega^2 m \hat{\eta}^- - i\omega d \hat{\eta}^- + (B\alpha^4 + T\alpha^2 + K) \hat{\eta}^- = -\hat{p}(-h) - \mu \left( i\alpha \left. \frac{dU_b}{dy} \right|_{-h} \hat{\eta}^- \right) + 2\mu \left. \frac{d\hat{v}}{dy} \right|_{-h}. \end{cases} \quad (12) $$

Following Shankar and Kumaran [30], the boundary conditions at the perturbed interface are implemented using Taylor series about the base state at $y = \pm h$. At linear order, the flow velocity at the walls reads

$$ \mathbf{u}(x, y = \zeta^{\pm}, t) = \mathbf{u}'(x, y = \pm h, t) + \eta^{\pm} \frac{dU_b}{dy} \Big|_{\pm h} \mathbf{e}_x. \quad (13) $$
Thus, the kinematic boundary conditions (7) become

$$ \hat{u}(\pm h) + \hat{\eta}^{\pm} \frac{dU_b}{dy} \Big|_{\pm h} = 0 \quad \text{and} \quad \hat{v}(\pm h) = -i\omega\hat{\eta}^{\pm}. \qquad (14) $$

Then, by using Eqs. (14) and the divergence-free condition, the right-hand side of Eqs. (12) can be further simplified, leading to

$$ \left\{ \begin{array}{l} -\omega^2 m \hat{\eta}^{+} - i \omega d \hat{\eta}^{+} + (B \alpha^4 + T \alpha^2 + K) \hat{\eta}^{+} = + \hat{p}(h) - \mu \frac{d\hat{v}}{dy} \Big|_{h}, \\ -\omega^2 m \hat{\eta}^{-} - i \omega d \hat{\eta}^{-} + (B \alpha^4 + T \alpha^2 + K) \hat{\eta}^{-} = - \hat{p}(-h) + \mu \frac{d\hat{v}}{dy} \Big|_{-h}. \end{array} \right. \qquad (15) $$

Equations (10), (11), and (15) completely govern the dynamics of small-amplitude perturbations and take into account the linearized fluid-structure coupling as derived from the exact hydrodynamic forces.

### III. NUMERICAL METHODS

In this work, we follow the general framework described by Manning et al. [26] for avoiding spurious eigenvalues. First of all, we rewrite Eqs. (15) using velocity components at the boundaries. For illustration purposes, only the upper wall is here considered.

Using the condition $-i\omega\hat{\eta}^{+} = \hat{v}(h)$, we obtain

$$ -i\omega\hat{v}(h) = W_1\hat{v}(h) + W_2\hat{u}(h) + \frac{1}{m}p(h), \qquad (16) $$

with

$$ W_1 = -\frac{d}{m} - 2\frac{\mu}{m}\frac{d}{dy} \quad \text{and} \quad W_2 = \frac{(B\alpha^4 + T\alpha^2 + K)}{m\left.\frac{dU_b}{dy}\right|_h} - \frac{\mu}{m}i\alpha. $$

The kinematic condition $\hat{u}(h) + \hat{\eta}^{+} \frac{dU_b}{dy}|_h = 0$ is recast as

$$ -i\omega\hat{u}(h) + \hat{v}(h) \frac{dU_b}{dy} \Big|_h = 0. \qquad (17) $$

The velocity components and pressure are discretized in the y-direction using a Chebyshev collocation method. To avoid spurious pressure modes, we consider the so-called $\mathbb{P}_N - \mathbb{P}_{N-2}$ approximation in which the pressure is approximated with a polynomial of degree $N-2$ while the velocity is discretized with a polynomial of degree $N$ [31]. From this point, we note the vectors con-

taining the unknowns: $\mathbf{V} = (\hat{u}_0, \hat{u}_N, \hat{v}_0, \hat{v}_N, \hat{u}_1, \dots, \hat{u}_{N-1}, \hat{v}_1, \dots, \hat{v}_{N-1})$ and $\mathbf{P}_I = (\hat{p}_1, \dots, \hat{p}_{N-1})$, where we separate the boundary values ($\mathbf{V}_{BC}$) from the interior points ($\mathbf{V}_I$). Hence, the discretized counterpart of the previous continuous model [Eqs. (10), (15), and (14)] reads

$$ \underbrace{\mathbf{A}_1}_{4} \quad \underbrace{\mathbf{A}_2}_{2(N-2)} \quad \underbrace{\mathbf{A}_3}_{N-2} \quad \mathbf{V}_{BC} \quad \mathbf{V}_I \quad \mathbf{P}_I $$
$$ = i\omega \begin{bmatrix} \mathbf{I} & \mathbf{0} & \mathbf{0} \\ \mathbf{0} & \mathbf{I} & \mathbf{0} \\ \mathbf{0} & \mathbf{0} & \mathbf{0} \end{bmatrix} \begin{bmatrix} \mathbf{V}_{\text{BC}} \\ \mathbf{V}_{\text{I}} \\ \mathbf{P}_{\text{I}} \end{bmatrix} \begin{bmatrix} 4 \\ 2(N-2) \\ N-2 \end{bmatrix} \quad (18) $$

where the divergence-free condition is imposed on the interior points and **I** denotes the identity matrix. The derivative matrices based on Chebyshev polynomials are either expressed on the interior points only (for the pressure) or all the nodes (for the velocity components) [31]. Note that the boundary equations involve the pressure at $y = \pm h$: these values are readily obtained by polynomial interpolation with spectral accuracy, corresponding to matrix $A_3$.

The discrete counterpart of the divergence-free condition reads $C_2V_I = 0$. Hence, from

$$ \mathbf{0} = C_2(i\omega\mathbf{V}_I) = C_2\mathbf{B}_1\mathbf{V}_{\text{BC}} + C_2\mathbf{B}_2\mathbf{V}_I + C_2\mathbf{B}_3\mathbf{P}_I, $$

the vector $\mathbf{P}_I$ can be expressed as a function of $\mathbf{V}_I$ and $\mathbf{V}_{\text{BC}}$:

$$ \mathbf{P}_I = \underbrace{-(C_2\mathbf{B}_3)^{-1}C_2\mathbf{B}_1}_{\mathbf{M}_{\text{BC}}} \mathbf{V}_{\text{BC}} - \underbrace{(C_2\mathbf{B}_3)^{-1}C_2\mathbf{B}_2}_{\mathbf{M}_I} \mathbf{V}_I. $$

Thus eliminating the pressure, the system (18) is recast as

$$ \begin{bmatrix} \mathbf{A}_1 + \mathbf{A}_3\mathbf{M}_{\text{BC}} \\ \mathbf{B}_1 + \mathbf{B}_3\mathbf{M}_{\text{BC}} \end{bmatrix} \begin{bmatrix} \mathbf{A}_2 + \mathbf{A}_3\mathbf{M}_I \\ \mathbf{B}_2 + \mathbf{B}_3\mathbf{M}_I \end{bmatrix} \begin{bmatrix} \mathbf{V}_{\text{BC}} \\ \mathbf{V}_I \end{bmatrix} = i\omega \begin{bmatrix} \mathbf{I} & \mathbf{0} \\ \mathbf{0} & \mathbf{I} \end{bmatrix} \begin{bmatrix} \mathbf{V}_{\text{BC}} \\ \mathbf{V}_I \end{bmatrix}. \quad (19) $$

System (19) still contains $N-2$ null eigenvalues due to the divergence-free constraint. We can further reduce Eq. (19) by eliminating the streamwise velocity components at the interior points

[26]. Indeed, using $\mathbf{V}_I = (\underbrace{\hat{u}_1, \dots, u_{\hat{N}-1}}_{\mathbf{U}}, \hat{v}_1, \dots, v_{\hat{N}-1})$, the divergence-free condition $C_2V_I = 0$ becomes $i\alpha\mathbf{U} + C_2\mathbf{V}\mathbf{V} = 0$. Thus, for $\alpha \neq 0$, the streamwise velocity $\mathbf{U}$ is obtained as a function of $\mathbf{V}$. This then leads to a discrete version of the Orr-Sommerfeld equation for the fluid-structure interaction problem of the form

$$ \begin{bmatrix} \bullet & \bullet \\ \bullet & \bullet \end{bmatrix} \begin{bmatrix} \mathbf{V}_{\text{BC}} \\ \mathbf{V} \end{bmatrix} = i\omega \begin{bmatrix} \mathbf{V}_{\text{BC}} \\ \mathbf{V} \end{bmatrix}. \quad (20) $$

This system may be further reduced [32], by considering perturbations of either sinuous or varicose symmetry and using only half of the channel together with derivative operators appropriate for the symmetry of each component of the different flow fields.

Apart from the fact that the above algebraic transformations remove spurious eigenvalues, they also drastically reduce the computational effort. The system is either solved using the QZ algorithm from the LAPACK library or an Arnoldi technique provided by the ARPACK software. The numerical procedure is validated and discussed in the Appendix. The number of collocation points is varied from 100 to 300 as the Reynolds number is increased.

## IV. LINEAR STABILITY RESULTS

After the formulation of the linear fluid-structure interaction problem and the presentation of the numerical methods, we are now in a position to analyze its dynamics. The different classes of modes and their dependence on the control parameters are investigated in detail. Specific attention is devoted to provide physical insight through total energy budget analyses.

### A. Spectra and classes of modes

A typical spectrum is shown in Fig. 2, corresponding to a base configuration at $Re = 7000$, $V_R = 1$, $B_* = 4$, $d_* = 0.2$, and perturbations with wavenumber $\alpha = 0.6$. Since the base state is symmetric in $y$, the entire spectrum consists of the same number of modes of either varicose (red symbols) or sinuous (blue) symmetry. The Orr-Sommerfeld modes ($\times$ and $+$) are essentially due to the base
<figure><img src="image_4.png" /><figcaption>FIG. 2. Eigenvalue spectrum for perturbations with $\alpha = 0.6$ at $Re = 7000$, $V_R = 1$, $B_\star = 4$, and $d_\star = 0.2$. Orr-Sommerfeld modes (+, sinuous; ×, varicose) are located on three main branches (A, P, and S) and dominated by TS mode (circle); four modes ($\langle \downarrow \rangle$ and $\rangle$) are due to fluid-structure interactions. Modes of varicose and sinuous symmetry are shown in red and blue, respectively.</figcaption></figure>

Poiseuille flow and organized in three branches (classically labeled A, P, and S [33]), as for rigid channels. This part of the spectrum is dominated by the TS mode (indicated by the circle in Fig. 2) of sinuous symmetry. In the present configuration, the coupling between the fluid and wall equations leads to four additional eigenvalues (indicated by $\langle \downarrow \rangle$ and $\rangle$) and referenced hereafter as FSI modes. Two of these FSI modes travel upstream ($\omega_r < 0$), while the other pair of FSI modes propagates along the flow direction ($\omega_r > 0$); each of these pairs consists of a sinuous and a varicose mode.

To gain a better understanding of these eigenmodes, we monitor changes in the spectrum resulting from the variation of some control parameters. A few typical scenarios are shown in Fig. 3.

FSI modes strongly depend on the reduced velocity $V_R$, as shown in Fig. 3(a). When $V_R \to 0$, which corresponds to approaching the rigid-walls case, the growth rates of the FSI modes reach neutrality ($\omega_i \to 0$) while their phase velocities tend to infinity ($\omega_r \to \pm\infty$). For the range of base state configurations shown here, the upstream propagating FSI modes are always stable, albeit with a weak decay rate, and their (negative) phase velocities reach very small values as $V_R$ is increased: this behavior is characteristic of divergence modes, as observed for boundary-layer flows along highly damped walls [19]. On the other hand, the downstream propagating FSI modes are strongly destabilized as $V_R$ is increased; these modes are identified as TWF modes. In contrast with the FSI mode dynamics, the TS mode only weakly depends on $V_R$ since it is mainly driven by the shear flow [see inset in Fig. 3(a)]. Consistently with the definition of $V_R$, the eigenvalue corresponding to the TS mode matches the one found for the rigid case when $V_R \to 0$. The other Orr-Sommerfeld modes on the A branch also display only a weak dependence on $V_R$, while those on the P and S branches appear to be mostly unaffected.

The influence of the wall-damping parameter $d_\star$ is shown in Figs. 3(b) and 3(c) for $V_R = 1$ and $V_R = 2$, respectively. The growth rate $\omega_i$ of the downstream propagating TWF modes is seen to significantly decrease with wall dissipation $d_\star$. Thus wall damping has a strongly stabilizing effect on both sinuous and varicose TWF modes. For $V_R = 1$, Fig. 3(b) shows that wall damping has a similar stabilizing effect on the upstream propagating FSI modes. However, at the larger value $V_R = 2$ of the reduced velocity [Fig. 3(c)], these FSI modes are nearly stationary divergence modes. In that regime, an increase in wall damping $d_\star$ results in an increase of their negative growth rate, and the phase velocity is seen to vanish in the limit of large wall damping $d_\star$. Here destabilization of
<figure><img src="image_5.png" /><figcaption>(a)</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>(b)</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>(c)</figcaption></figure>

FIG. 3. Influence of $V_R$ and $d_*$ on eigenspectra for $\alpha = 0.6$ at Re = 7000 and $B_* = 4$. (a) $V_R = 0.01, 0.5, 1, 1.5,$ and 2 (in blue, black, red, green, and brown, respectively) at $d_* = 0$. (b) $d_* = 0, 0.2, 0.5, 1,$ and 2 (in black, blue, green, red, and brown, respectively) at $V_R = 1$. (c) $d_* = 0, 0.4, 0.8, 1, 2, 4, 6, 10,$ and 20 (in black, blue, red, green, orange, light blue, grey, pink, and brown, respectively) at $V_R = 2$. Some mode trajectories are also represented.

the divergence modes occurs for large values of $d_*$. One may recall that while the sinuous TWF and divergence modes were investigated numerically by Davies and Carpenter [20], the varicose TWF and divergence modes have not been explored for the plane channel flow.

Thus, the influence of the various control parameters may be summarized as follows: the TS mode is temporally damped by an increase of $V_R$, but its temporal amplification rate may be amplified with an increase in wall damping for a certain range of streamwise wavenumbers. A different behavior is found when the TWF modes are considered: their growth rate increases with $V_R$, while it decreases with $d_*$. This is consistent with the mode classification given by Benjamin [14] (i.e., class A for TS modes and class B for TWF modes). Moreover, the divergence mode is seen to be amplified for high values of wall-damping parameter. Finally, the effect of the wall compliance is seen to be negligible for both P and S branches.

To further identify the different types of modes with respect to Benjamin's classification, the eigenfunctions of a few selected modes are shown in Fig. 4.
<figure><img src="image_8.png" /><figcaption>FIG. 4. Eigenfunctions for Re = 7000 and $B_\star = 4$. (a) TS (sinuous) mode for $\alpha = 1$, $d_\star = 4$, and $V_R = 0.5$. (b) Sinuous TWF mode for $\alpha = 0.6$, $d_\star = 0$, and $V_R = 1$. (c) Varicose TWF mode for $\alpha = 0.6$, $d_\star = 0$, and $V_R = 1$. (d) Sinuous divergence mode for $\alpha = 0.6$, $d_\star = 20$, and $V_R = 2$. (e) Varicose divergence mode for $\alpha = 0.6$, $d_\star = 20$, and $V_R = 2$. The phase angles of eigenfunctions are denoted as $\phi_u$ and $\phi_v$ for $u'$ and $v'$, respectively. The eigenfunctions are normalized to the unit kinetic energy norm.</figcaption></figure>

One may recall that the mechanism whereby both TWF modes and TS modes grow involves the action of the streamwise velocity base flow gradient along the wall-normal position working against the Reynolds stresses. In the absence of an inflection point, both modes involve a phase shift of the disturbance velocity at some distance from the wall [34]. However, mechanisms are quite different
for each of these modes. Within a large-Reynolds-number asymptotic theory, viscous effects are
only present in the vicinity of the viscous wall layer or the critical layer, located at the wall-normal
position $y_c$ where $U_b(y_c) = \omega_r/\alpha$ [34]. The inviscid approximation is therefore accurate in the other
regions of the flow. For the TWF mode, the instability mechanism is essentially driven by the wall.
Indeed, the instability is amplified if the work done by the pressure disturbance on the wall is positive
when averaged over one period [12]. As first proven by Miles [15] for water waves, Benjamin [13]
shows for the boundary-layer-flow case that it results from a phase shift between disturbance wall
pressure and wall displacement. In particular, Benjamin [13] proves that the pressure at the wall is
associated with the integrated effect of the velocity perturbations along the wall direction and is a
consequence of a phase shift for velocity components near the critical layer. We note hereafter $\phi_u$
and $\phi_v$ as the phases of the respective velocity components $\hat{u}$ and $\hat{v}$ of the associated eigenfunctions.
For TWF modes [see Figs. 4(b) and 4(c)], the essential phase shift (i.e., $\phi_u - \phi_v \neq \pi/2$) occurs
near the critical layer in the limit of large Reynolds numbers (class B modes). For the TS mode
[Fig. 4(a)], the phase shift is rather associated with the viscous wall layer. They belong to class
A modes and are stabilized through a transfer of energy to the wall. Hence it is essential, when
addressing the different classes of modes, to monitor the phases of the perturbation components.
Figure 4(b) shows that the sinuous TWF mode exhibits a clear phase shift near the critical layer
in agreement with class B modes. It is consistent with the theoretical investigation of Davies and
Carpenter [20] for the same flow case. The linear behavior along with the wall-normal position of
$\hat{u}$ is associated with a displacement of the Poiseuille solution when the walls are shifted with $\eta$. It
is easily verified that the small deviation from the Poiseuille solution due to sinuous motion of the
walls is proportional to $-2\eta\dot{y}$ (see also [35]). Interestingly, the varicose case [Fig. 4(c)], not studied
by Davies and Carpenter [20], also exhibits a phase shift.

Finally, the divergence mode [Figs. 4(d) and 4(e)] exhibits a phase shift in the viscous wall layer.
Nevertheless, due to the low velocity phase, the viscous wall layer and the critical layer are not
well separated. As a consequence, the theoretical framework derived by Davies and Carpenter [20]
cannot be applied for these modes.

For all FSI modes (both TWF and divergence), velocity and pressure fluctuations are concentrated near the wall. For the TS mode, as for the rigid case, the streamwise velocity component peaks at the critical layer.

B. Temporal growth and instability onset

The stability features of compliant channel flow configurations depend on a large number of
parameters. In the previous section we have shown that the reduced velocity $V_R$ is the main control
parameter governing fluid-structure interactions, but growth rates of the different classes of modes
may also depend significantly on Reynolds number Re, wall dissipation $d_*$, and flexural rigidity $B_*$.
In the present section, we will map out the stability characteristics by monitoring the dispersion
relation in the ($V_R, \alpha$) plane for selected values of the other relevant control parameters. By
taking advantage of the base flow symmetry, sinuous and varicose perturbations may be efficiently
computed and their properties are here discussed separately. We first consider perturbations of
varicose symmetry, which are generally the modes most amplified by fluid-structure interactions.

Figure 5(a) shows isolines of the temporal growth rate $\omega_i$ of the leading varicose eigenmode in the ($V_R, \alpha$) plane at Re = 5000, $d_* = 10$, and $B_* = 4$. The neutral curve ($\omega_i = 0$) exhibits two distinct minima at $V_R \approx 1.6$ and $V_R \approx 3.9$. Beyond onset, a finite range of wavenumbers $\alpha$ display positive temporal growth rates. The associated phase speeds $\omega_r/\alpha$ are given in Fig. 5(b). It is found that near the first minimum $V_R \approx 1.6$, modes travel with vanishingly small phase speeds. In contrast, near the second minimum $V_R \approx 3.9$, modes travel with phase speeds of the order of the mean base flow velocity. This behavior is characteristic of divergence and TWF modes, which are thus each found to dominate the perturbation dynamics in distinct regions of the ($V_R, \alpha$) plane for these parameter settings. Figure 6 illustrates the influence of wall dissipation $d_*$, for the same values of Re and $B_*$. It is observed that energy dissipation in the wall only weakly influences the temporal growth
<figure><img src="image_9.png" /><figcaption>FIG. 5. Dispersion relation for leading varicose mode in ($V_R$, $\alpha$) plane for $Re = 5000$, $d_\star = 10$, and $B_\star = 4$: contours of (a) growth rate $\omega_i$ and (b) phase velocity $\omega_i/\alpha$.</figcaption></figure>

rate for $0 \le d_\star \le 1$, while stronger stabilization occurs for $d_\star > 1$. At these high dissipation rates, the growth rate and the range of unstable wavenumbers are greatly reduced; however, the critical value of $V_R$ for onset of instability (denoted $V_R^c$ hereafter) remains of the same order of magnitude. Monitoring the neutral curves more precisely reveals that $V_R^c$ increases from $V_R^c \simeq 0.45$ for $d_\star = 0$ to reach a maximal value of about $V_R^c \simeq 1.6$ for $d_\star = 10$ and decreases again for larger values of $d_\star$ ($V_R^c \simeq 1.0$ for $d_\star = 100$). This nonmonotonous effect of wall dissipation on instability onset is due to a change of the nature of the leading eigenmode: while the instability is dominated by the TWF mode at low values of $d_\star$, the unstable dynamics is governed by the divergence mode for strong dissipation in the compliant wall as already suggested by Fig. 5.

<figure><img src="image_10.png" /><figcaption>FIG. 6. Temporal growth rates $\omega_i$ for the varicose instability at $Re = 5000$ and $B_\star = 4$ for $d_\star = 0, 0.6, 1, 10, 60$, and $100$. The contour levels are 0, 0.1, 0.2, and 0.3.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>FIG. 7. Marginal curves for onset of varicose instability at $B_{\star} = 4$ for Re = 10000, 40000, 80000, 160000, and 240000 (in black, blue, red, green, and brown respectively) and (a) $d_{\star} = 0$, (b) $d_{\star} = 10$, and (c) $d_{\star} = 100$.</figcaption></figure>

This observation is further illustrated in Fig. 7, where neutral curves are shown for various Reynolds numbers and wall dissipations. In particular, for Re = 10000 and $d_{\star} = 10$ [Fig. 7(b)], we observe the coexistence of the divergence and TWF modes. The critical reduced velocities $V_R^c$ for the divergence and TWF modes are $\approx 2$ and $\approx 4$, respectively. Neutral curves associated with the TWF mode [Fig. 7(a)] are seen to be almost independent of the Reynolds number. For the divergence mode, we observe a destabilizing effect of viscosity [Fig. 7(c)]: the critical value $V_R^c$ for onset of instability is seen to increase with Reynolds number, which is consistent with the fact that divergence modes are intimately connected to the viscous wall layer. Nevertheless, in the regime dominated by divergence modes, the temporal amplification rates reach much lower values than those prevailing for TWF modes at low values of $d_{\star}$ (see Fig. 6).

After the previous extensive discussions of results for varicose perturbations we now focus on the sinuous symmetry. The linear dynamics of sinuous eigenmodes is very similar to that of their varicose counterpart, except that the (sinuous) Tollmien-Schlichting modes may also display positive growth rates.

The growth-rate isolines in the ($V_R$, $\alpha$) plane of Fig. 8(a) are obtained for Re = 8000, $B_{\star} = 4$, and $d_{\star} = 10$. This Reynolds number (based on channel diameter and mean fluid velocity) is slightly in excess of the critical value Re_c $\simeq 7696$ for Tollmien-Schlichting instability developing in rigid channel flow. As for the varicose case, the neutral curve exhibits two minima. The first appears near $V_R \approx 0$ while the second minimum is close to $V_R \approx 2$. The associated phase speeds $\omega_r/\alpha$ given in Fig. 8(b) indicate that divergence, TWF, and TS modes can be involved in this regime.

Figure 9 shows the evolution of sinuous temporal instability characteristics with the wall dissipation parameter $d_{\star}$. For $d_{\star} = 0$ and low values of $V_R$, corresponding to near-rigid compliant walls, there exists a narrow band of unstable wavenumbers near $\alpha = 1$ associated with unstable
<figure><img src="image_12.png" /><figcaption>FIG. 8. Dispersion relation for leading sinuous mode in (V<sub>R</sub>, α) plane for Re = 8000, d<sub>*</sub> = 10, and B<sub>*</sub> = 4: contours of (a) growth rate ω<sub>i</sub> and (b) phase velocity ω<sub>r</sub>/α.</figcaption></figure>

Tollmien-Schlichting modes. At larger values of V<sub>R</sub>, stronger fluid-structure coupling leads to the destabilization of (sinuous) TWF modes: the temporal growth rates ω<sub>i</sub> and the range of unstable wavenumbers rapidly increase with V<sub>R</sub>. Note that there exists a narrow region near V<sub>R</sub> = 1 where both TS and TWF modes are stable. As d<sub>*</sub> is increased from 0 to 0.6, it is observed that the regions corresponding to unstable TS and TWF modes merge, giving rise to the so-called transition mode.

Except for the merging of TS and TWF instabilities, the instability features of sinuous perturba-
tions shown in Fig. 9 are very similar to those observed for varicose perturbations. For d\* > 1,

<figure><img src="image_13.png" /><figcaption>FIG. 9. Temporal growth rates ωᵢ for the sinuous instability at Re = 8000 and B* = 4 for d* = 0, 0.6, 1, 10, 60, and 100. The contour levels are 0, 0.1, 0.2, and 0.3.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>FIG. 10. Marginal curves for onset of sinuous instability at $B_{\star} = 4$ for Re = 10000, 40000, 80000, 160000, and 240000 (in black, blue, red, green, and brown respectively) and (a) $d_{\star} = 0$, (b) $d_{\star} = 10$, and (c) $d_{\star} = 100$.</figcaption></figure>

wall dissipation $d_{\star}$ significantly reduces the temporal growth rates and the range of unstable wavenumbers. For $d_{\star} \ge 10$, Fig. 9 shows that divergence and transition modes coexist. In addition, the critical reduced velocity $V_R^c$ for the onset of the divergence mode is only weakly influenced by the wall dissipation ($V_R^c \approx 2$). The critical curves shown in Fig. 10 correspond to marginal ($\omega_i = 0$) sinuous instability for a range of Reynolds numbers and $d_{\star} = 0, 10$, and $100$. Due to the fact that, for Re > 7696, TS instability prevails at low values of $V_R$, down to $V_R = 0$, a critical value of reduced velocity $V_R^c$ cannot be defined for onset of sinuous instability. For the sinuous instability, it may be hard to distinguish between TS or TWF modes since branch switching occurs as some parameters are continuously varied. To better illustrate this phenomenon, Fig. 11 shows the dispersion relation for both branches in the range $0.55 < \alpha < 1.05$ for $d_{\star} = 0.10, 0.13$, and $0.16$, at Re = 10000, $B_{\star} = 1$, and $V_R = 1$. For $\alpha < 0.8$, the branches with largest temporal growth rate [upper branches in Fig. 11(a)] are of TS type while the other branches are always stable in that wavenumber range and can be identified as TWF modes. For $\alpha > 0.8$, the unstable branch displays a growth rate $\omega_i$ rapidly increasing with $\alpha$, and is found to correspond to a mode of the TWF type, while the other branch is strongly stabilized at these wavenumbers. Due to the branch switching that occurs near $\alpha = 0.8$ and $d_{\star} = 0.13$, the unstable TWF branch prevailing for $\alpha > 0.8$ is continuously connected to the TS branch when $d_{\star} > 0.13$ while it is continued as a stable TWF mode for $\alpha < 0.8$ when $d_{\star} < 0.13$. This behavior is further illustrated in Fig. 12. For a small amount of wall dissipation, the TS mode is seen to be damped as $V_R$ is increasing [Fig. 12(a)], for all the Reynolds numbers considered. In particular, for $d_{\star} = 0$, the critical reduced velocity of TS mode suppression is varying from $V_R \approx 0.8$ for Re = 10000 to $V_R \approx 4$ for Re = 40000. However, as $d_{\star}$ is increased beyond 0.13, a transition mode emerges and it is no longer possible to distinguish between TWF and TS modes [Figs. 12(b)–(d)].
<figure><img src="image_15.png" /><figcaption>FIG. 11. Dispersion relation for the two leading sinuous modes at Re = 10000, B<sub>*</sub> = 1, and V<sub>R</sub> = 1:
(a) temporal growth rate ω<sub>i</sub> and (b) frequency ω<sub>r</sub>. Branch switching occurs near α = 0.8 and d<sub>*</sub> = 0.13.</figcaption></figure>

Finally, we present some results for the destabilization of the divergence mode and its dependence on wall dissipation $d_*$ and Reynolds number Re. We compute the critical value $V_R^c$ for onset of instability in the ($V_R$, $\alpha$) plane at fixed values of the other control parameters. A Newton-Raphson search algorithm with an adaptive step has been implemented to automatically identify the start of the neutral curve in the ($V_R$, $\alpha$) plane for different values of $d_*$ and Re. The critical value $V_R^c$ (and associated wavenumber $\alpha^c$) for onset of instability is then obtained when $\frac{dV_R}{d\alpha}$ vanishes along the neutral curve.

Critical curves as functions of $d_*$ for different Re values are plotted in Figs. 13(a) and 13(b), respectively, for the divergence modes of varicose and sinuous symmetry. Flexural rigidity is kept constant at $B_* = 1$. Both sets of curves clearly indicate the stabilizing effect of the Reynolds number on the divergence mode, for both the varicose and sinuous cases. For the varicose case, unstable divergence modes only occur at relatively high wall dissipation, $d_* > 6$ for Re = 5000 and $d_* > 9.5$ for Re = 80000. For the sinuous case, divergence mode instability already starts for $d_*$ in the range 2–4, with a weaker Reynolds number dependence. In both varicose and sinuous cases, the critical $V_R^c$ appears to asymptote towards a finite limit for large values of $d_*$.

As discussed above, the combined effects of the different wall parameters are conveniently summarized by monitoring either the variations of the critical reduced velocity $V_R^c$ or the critical Reynolds number $Re^c$ depending on whether the FSI or TS mode is considered. First, we focus on FSI instability modes of varicose symmetry. To that purpose, we consider a high value of the Reynolds number in order to investigate only the influence of the wall properties.

Results are reported in Fig. 14 for $0 \le d_* \le 60$ and $0.2 \le B_* \le 10$. The inset in the figure shows that for small values of $d_*$, flexural rigidity has a moderately stabilizing effect on TWF modes: $V_R^c$ increases as $B_*$ is increased for fixed $d_*$. In the range $5 < d_* < 15$, onset of instability is seen to display an almost universal behavior with a linear relationship between $V_R^c$ and $d_*$, almost independent of $B_*$. In this regime, instability always occurs by a TWF mode. For larger wall dissipation rates, $d_* > 20$, the dynamics is dominated by the divergence modes, associated with vanishing phase velocities. In this latter regime, the critical values $V_R^c$ weakly depend on flexural rigidity $B_*$ and decrease with increasing wall dissipation $d_*$. For large values of $d_*$, a limit value of $V_R^c \approx 2.85$ is asymptotically reached, independently of $B_*$. Note that the crossover from
<figure><img src="image_16.png" /><figcaption>FIG. 12. Neutral curves for the TS mode for Re = 10000, 20000, 30000, and 40000 (in black, red, green, and blue, respectively) and $B_{\star} = 1$ in the plane $(\alpha, V_{R})$. (a) $d_{\star} = 0$, (b) $d_{\star} = 0.5$, (c) $d_{\star} = 1$, and (d) $d_{\star} = 10$.</figcaption></figure>

the TWF-dominated instabilities (low $d_{\star}$) to the divergence-dominated instabilities (high $d_{\star}$) also depends on $B_{\star}$.

Finally, we address the influence of the wall properties on the stability of TS modes by monitoring the critical Re$^c$. We restrict our analysis to the case without wall damping. Results are summarized in Fig. 15, where the critical Reynolds number for the TS mode is displayed as a function of $V_{R}$ for $B_{\star} = 1, 2$, and 4. The figure shows that the TS mode is significantly stabilized with increasing wall compliance (i.e., increasing $V_{R}$), and this stabilizing effect is enhanced at low flexural rigidity $B_{\star}$. In

<figure><img src="image_17.png" /><figcaption>FIG. 13. Divergence mode. Critical reduced velocity distribution $V_{R}^{c}$ with the wall-damping parameter and for Re = 5000, 10000, 20000, 40000, and 80000 with $B_{\star} = 1$: (a) the varicose case and (b) the sinuous case.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>FIG. 14. Evolution of critical reduced velocity $V_R^e$ with wall-damping parameter $d_*$, for a range of flexural rigidity values $B_*$ (from 0.2 to 10). Perturbations of varicose symmetry are considered.</figcaption></figure>

the figure, onset of TWF instability is indicated by the nearly vertical lines; recall that these modes are almost independent of the Reynolds number. Dotted lines correspond to sinuous TWF modes while dashed lines have been computed for varicose TWF modes. It is immediately apparent from the figure that the potential effect of a compliant wall regarding transition delay is clearly limited by the destabilization of TWF modes. Indeed, the region of the ($V_R$, $Re$) plane corresponding to stable configurations is located below the $Re^c$ curves for onset of TS instability and delimited at large values of $V_R$ by onset of TWF instability. The results of the present investigation reveal that the extent of the stable region is drastically reduced by taking into account TWF perturbations of varicose symmetry, which display growth rates in excess of their sinuous counterparts. For instance,

<figure><img src="image_19.png" /><figcaption>FIG. 15. The critical Reynolds number associated with TS instability is plotted as a function of $V_R$ for $B_* = 1$, 2, and 4 in solid lines. Nearly vertical lines indicate onset of TWF instability for perturbations of sinuous (dotted lines) and varicose (dashed lines) symmetry. The colored regions correspond to stable configurations at $B_* = 2$, allowing only sinuous perturbations (light red) or all types of perturbations (solid red).</figcaption></figure>
when only sinuous perturbations are taken into account, the upper right corner of the shaded region in Fig. 15 indicates that stable configurations would be possible up to Re = 8370, for $B_\star = 2$. However, when perturbations of all symmetries are allowed, the stable domain is reduced to the solid red region, with maximal Reynolds number of Re = 7970, only marginally larger than the critical Reynolds number of Re = 7696 corresponding to the rigid case. Similar conclusions hold for different values of $B_\star$. Moreover, since it has been shown that wall dissipation has a destabilizing effect on TS waves, it is not expected that using $d_\star > 0$ could increase the extent of the stable regions in parameter space. Thus it seems impossible to significantly delay instability onset by wall compliance.

### C. Comparison with asymptotic theories for $d_\star = 0$

Davies and Carpenter [20] derived an analytical expression for the wall pressure for $d_\star = 0$ in the limit of small $\alpha$ and high Reynolds numbers, for modes of sinuous symmetry. This pressure, denoted $\hat{p}(\alpha, c, U_{\text{ref}})$, is obtained as a function of the wavenumber $\alpha$, the phase velocity $c = \omega_r/\alpha$, and a reference value for the fluid velocity, $U_{\text{ref}} = Q/2h$, and includes the effects of both the critical and viscous layers. The wall pressure is obtained as an expansion in $\alpha$ up to $\alpha^2$: $\hat{p} = p_0 + \alpha^2 p_1$, where only the term $p_1$ includes the effect of viscous and critical layers. For the viscous layer, the approximation is carried out up to $O((\alpha \text{Re})^{-1/2})$. Neglecting the viscous stress at the wall, they obtain the dispersion relation:

$$m(c^2 - c_0^2) + \hat{p}(\alpha, c, U_{\text{ref}}) + i\left(\frac{c}{\alpha}\right)d = 0, \quad (21)$$

with $c_0 = \sqrt{\frac{1}{m}(B\alpha^2 + T + \frac{K}{\alpha^2})}$ the free wave speed for the wall. Davies and Carpenter [20] express the onset of instability with the Reynolds number. Here, we suggest that it is more appropriate to use the reduced velocity $V_R$. The resulting neutral curves are shown in Fig. 16 for $B_\star = 4$. For the sinuous case, Fig. 16(a) shows a very good agreement between the analytical model and complete numerical resolution of the full system up to $\alpha \approx 3$. In order to remove the Reynolds number effect, a numerical solution for Re = $1 \times 10^6$ has been carried out. For the latter case, the grid mesh is increased up to $N = 300$ in order to correctly capture both viscous and critical layers. Figure 16(a) shows that the departure from the theoretical model is due to the expansion in terms of streamwise wavenumber up to $\alpha^2$; hence this approximation is no longer valid for $\alpha$ greater than 3.5 (not shown in [16]). In addition, one observes that neglecting the viscous stress at the wall and in the analytical expression of the pressure yields an almost perfect approximation of the exact dispersion relation. In Fig. 16(b), comparisons with varicose cases are shown. While the theoretical model is derived only for the sinuous symmetry, it is interesting to notice that in the limit of high Reynolds numbers, the model associated with merely the critical layer gives a quite accurate description of the varicose symmetry for $\alpha$ varying from 0.8 to 3. For streamwise wavenumbers greater than 4, the varicose and sinuous neutral curves fall in one single curve for all Reynolds numbers. For $\alpha < 0.8$, the varicose case exhibits a more complex Reynolds number dependence. Since the varicose mode always dominates over the sinuous mode, the critical value $V_R^c$ for onset of instability is associated with a varicose perturbation for all configurations that have been considered in the present study. The latter observation is also in agreement with results provided by Nagata and Cole [22]. In addition, it is found that the Reynolds number has a slight stabilizing effect. In Fig. 16, we provide a comparison with the theoretical model derived by Huang [36] for the varicose symmetry, only based on the critical layer. The figure shows a good agreement for moderate values of $\alpha$ between the model and the numerical simulation.

### D. Energy budgets

This final section addresses the energy transfer mechanisms between the different components of the compliant channel flow configurations, in order to shed further light on the dynamics and on
<figure><img src="image_20.png" /><figcaption>FIG. 16. Neutral curves for $d_* = 0$ and $B_* = 4$. For the full system of equations, the Reynolds numbers are fixed to $Re = 4 \times 10^4$, $8 \times 10^4$, $2.4 \times 10^5$, and $1 \times 10^6$. For the analytical dispersion relation $Re$ is fixed to $1 \times 10^6$. (a) Sinuous case and (b) sinuous and varicose cases comparison.</figcaption></figure>

the fundamental mechanisms promoting instability. This section is greatly influenced by the work of Domaradski and Metcalfe [29] and Carpenter and Morris [37].

For the compliant channel flow configuration under consideration, the total energy of the system is the sum of three components:

$$E_{tot} = E_{FK} + E_{WK} + E_{WP}, \qquad (22)$$

where $E_{FK}$ represents the fluid kinetic energy, while the wall energy consists of both kinetic and potential contributions, $E_{WK}$ and $E_{WP}$, respectively.

The global fluid kinetic energy is obtained by integration over the channel diameter:

$$E_{FK} = \int_{-h}^{+h} e(y)dy, \qquad (23)$$

where

$$e(y) = \rho \hat{u}(y) \cdot \hat{u}(y)^* \equiv \rho[\hat{u}(y)\hat{u}(y)^* + \hat{v}(y)\hat{v}(y)^*] \qquad (24)$$

denotes the local kinetic energy of the flow, averaged over $x$ at a given wall-normal position $y$, using the notations introduced in Sec. II C. The temporal variation of the local kinetic energy then follows
from the governing equation (10) as

$$
\begin{equation}
\begin{aligned}
2\omega_i e(y) = & -\rho[\hat{u}(y)\hat{v}(y)^* + \hat{u}(y)^*\hat{v}(y)] \underbrace{\frac{dU_b(y)}{dy}}_{\mathcal{P}: \text{Reynolds stress work against the mean shear}} - \underbrace{\frac{d}{dy}[\hat{p}(y)\hat{v}(y)^* + \hat{p}(y)^*\hat{v}(y)]}_{\mathcal{P}: \text{pressure diffusion}} \\
& -2\mu \underbrace{\left[ \frac{d\hat{u}(y)}{dy} \frac{d\hat{u}(y)^*}{dy} + \frac{d\hat{v}(y)}{dy} \frac{d\hat{v}(y)^*}{dy} + \alpha^2 (\hat{u}(y)\hat{u}(y)^* + \hat{v}(y)\hat{v}(y)^*) \right]}_{\mathcal{E}: \text{viscous dissipation}} \\
& +\mu \underbrace{\left[ \hat{u}(y) \frac{d\hat{u}(y)^*}{dy} + \hat{u}(y)^* \frac{d\hat{u}(y)}{dy} + \hat{v}(y) \frac{d\hat{v}(y)^*}{dy} + \hat{v}(y)^* \frac{d\hat{v}(y)}{dy} \right]}_{\mathcal{D}: \text{viscous diffusion}}
\end{aligned}
\tag{25}
\end{equation}
$$

and is the result of four distinct mechanisms as indicated in the above equation. Integration of this expression over the channel diameter leads to the equivalent equation governing the evolution of the total fluid kinetic energy:

$$
\begin{equation}
\begin{aligned}
2\omega_i E_{FK} = & -\int_{-h}^{+h} \rho[\hat{u}(y)\hat{v}(y)^* + \hat{u}(y)^*\hat{v}(y)] \frac{dU_b(y)}{dy} dy - [\hat{p}(y)\hat{v}(y)^* + \hat{p}(y)^*\hat{v}(y)]_{-h}^{+h} \\
& -2\mu \int_{-h}^{+h} \left[ \frac{d\hat{u}(y)}{dy} \frac{d\hat{u}(y)^*}{dy} + \frac{d\hat{v}(y)}{dy} \frac{d\hat{v}(y)^*}{dy} + \alpha^2 (\hat{u}(y)\hat{u}(y)^* + \hat{v}(y)\hat{v}(y)^*) \right] dy \\
& +\mu \left[ \hat{u}(y) \frac{d\hat{u}(y)^*}{dy} + \hat{u}(y)^* \frac{d\hat{u}(y)}{dy} + \hat{v}(y) \frac{d\hat{v}(y)^*}{dy} + \hat{v}(y)^* \frac{d\hat{v}(y)}{dy} \right]_{-h}^{+h}.
\end{aligned}
\tag{26}
\end{equation}
$$

While interaction with the base shear flow and viscous dissipation prevails throughout the channel cross section, pressure and viscous diffusion only contribute at the boundaries and transfer energy between the fluid and the compliant walls.

The kinetic and potential energies associated with the walls are obtained as

$$
E_{WK} = m|\omega|^2(|\hat{\eta}^+|^2 + |\hat{\eta}^-|^2) \quad \text{and} \quad E_{WP} = (B\alpha^4 + T\alpha^2 + K)(|\hat{\eta}^+|^2 + |\hat{\eta}^-|^2), \quad (27)
$$

respectively. Using wall equations (15), together with boundary conditions (14), yields the temporal
variation of the wall energy as

$$
\begin{equation}
\begin{aligned}
2\omega_i (E_{WK} + E_{WP}) = & -\underbrace{2d|\omega|^2(|\hat{\eta}^+|^2 + |\hat{\eta}^-|^2)}_{E_0} + \underbrace{[\hat{p}(y)\hat{v}(y)^* + \hat{p}(y)^*\hat{v}(y)]}_{E_1}^{+h} \\
& -\mu \underbrace{\left[\hat{v}(y) \frac{d\hat{v}(y)^*}{dy} + \hat{v}(y)^* \frac{d\hat{v}(y)}{dy}\right]}_{E_2}^{+h}.
\end{aligned}
\tag{28}
\end{equation}
$$

Thus, changes in total wall energy are seen to be the result of either dissipation within the wall (E₀) or energy exchange at the interface between the fluid and the compliant walls: work done by the pressure force (E₁) or the normal viscous stress (E₂). Both terms E₁ and E₂ also appear in Eq. (26) but with opposite sign; these fluid-structure interaction terms only account for an exchange of energy between the fluid and the walls but do not modify the total energy of the system.

The temporal variation of the total energy (22) is then obtained by adding (26) and (28), which
leads to the following integrated total energy budget:

$$
2\omega_i E_{tot} = - \underbrace{\int_{-h}^{+h} \rho [\hat{u}(y)\hat{v}(y)^* + \hat{u}(y)^*\hat{v}(y)] \frac{dU_b(y)}{dy} dy}_{C_1: \text{energy exchange with the base flow}}
$$
$$
\begin{equation}
\begin{aligned}
& \underbrace{-\mu \left[ \left( \frac{d\hat{u}(y)}{dy} \hat{\eta}^* + \frac{d\hat{u}(y)^*}{dy} \hat{\eta} \right) \frac{dU_b(y)}{dy} \right]_{-h}^{+h}}_{C_2: \text{energy exchange with the base flow at the walls}} \\
& -2\mu \underbrace{\int_{-h}^{+h} \left[ \frac{d\hat{u}(y)}{dy} \frac{d\hat{u}(y)^*}{dy} + \frac{d\hat{v}(y)}{dy} \frac{d\hat{v}(y)^*}{dy} + \alpha^2 (\hat{u}(y)\hat{u}(y)^* + \hat{v}(y)\hat{v}(y)^*) \right] dy}_{C_3: \text{viscous dissipation}} \\
& -\underbrace{2d|\omega|^2 (|\hat{\eta}^+|^2 + |\hat{\eta}^-|^2)}_{C_4: \text{wall damping}}.
\end{aligned}
\tag{29}
\end{equation}
$$

Hence, the only mechanisms that contribute to variations of the total energy are interactions with the base flow and dissipation (see also [19,35]). Energy transfer from or to the base flow occurs in the bulk ($C_1$) as well as at the boundaries ($C_2$), and energy dissipation takes place both in the fluid ($C_3$) and in the compliant walls ($C_4$). Note that the kinematic boundary conditions (14) have been used to bring to the fore the role of the base flow shear in the exchange term $C_2$. As underlined by Carpenter and Morris [37], the contribution $C_2$ arises from the interaction of the displaced mean flow and shear stress. In the literature, the terms $C_1$ and $C_2$ are often labeled as irreversible energy transfer from the base flow to the perturbation; depending on the signs and phases of the different components in $C_1$ and $C_2$, they may have a destabilizing or a stabilizing influence.

The energy budget (29) may be used to recover the temporal growth rate as

$$ \omega_i = \underbrace{\hat{C}_1 + \hat{C}_2 + \hat{C}_3 + \hat{C}_4}_{\Sigma}, \qquad (30) $$

where the different contributions have been renormalized by the total energy, $\hat{C}_i \equiv C_i/(2E_{tot})$.

Now that we have identified the different components that contribute to the variation of the perturbation energy, we proceed to analyze their role in the dynamics of the different classes of modes that prevail in the present configuration. In sequence we will address TS, divergence, and TWF modes and discuss the corresponding total energy budget as well as the spatial structure of the different contributions.

First we consider the stabilization mechanism of the TS mode as the reduced velocity $V_R$ is increased. To that purpose, we investigate configurations with $Re = 10000$, $B_* = 1$, and $d_* = 0$. For these typical control parameter values, the TS mode is stable for $V_R$ in excess of approximately 0.85, as shown by the black curve in Fig. 12(a). To study the influence of $V_R$ on the energy transfer mechanisms, the most unstable TS mode is considered as $V_R$ is varied; i.e., the streamwise wavenumber $\alpha$ is chosen to maximize the temporal growth rate $\omega_i$ for each value of $V_R$.

Figure 17(a) plots the components of the total energy budget as $V_R$ is increased. The associated growth rate $\omega_i$ is also reported in the same figure. The excellent agreement between the curves of $\omega_i$, derived from the eigenvalue problem, and of $\Sigma$, right-hand side of Eq. (30), gives confidence that the computation of the different energy terms is correctly implemented. The curves in Fig. 17(a) show that for the range of $V_R$ considered here, the major destabilizing contribution is due to the action of the basic velocity gradient working against the Reynolds stress ($\hat{C}_1$). However, as $V_R$ is increased, this production term $\hat{C}_1$ is observed to decrease and to be partially balanced by $\hat{C}_2$. This suggests two stabilizing mechanisms associated with the compliant wall: one reducing the bulk production term $\hat{C}_1$ and due to a modification of the perturbation velocity profiles, and another one directly connected to the wall term $\hat{C}_2$ and due to the displaced mean flow that acts as a dissipative term here. The viscous dissipation ($\hat{C}_3$) is seen to weaken as $V_R$ is increased, but the overall stabilizing influence dominates for increasing $V_R$.

The relative importance of the different components of the total energy (22) are displayed in Fig. 17(b) as $V_R$ is varied, using the notation $\hat{E}_{FK} = E_{FK}/E_{tot}$, $\hat{E}_{WK} = E_{WK}/E_{tot}$, $\hat{E}_{WP} = E_{WP}/E_{tot}$, and $\hat{E}_W = \hat{E}_{WK} + \hat{E}_{WP}$. This plot shows that, as $V_R$ is increased, a small part of the fluid kinetic
<figure><img src="image_21.png" /><figcaption>FIG. 17. Energy of most unstable TS mode as $V_R$ is increased at $Re = 10000$, $B_\star = 1$, and $d_\star = 0$. (a) Evolution of the integrated total energy budget and comparison with temporal growth rate $\omega_i$. (b) Different contributions to the total energy for $d_\star = 0$.</figcaption></figure>

energy is indeed transferred to the wall, mainly as potential energy. This is consistent with a class A mode.

We now analyze the spatial structure of the fluid kinetic energy budget (25). Figure 18 shows the wall-normal profiles of the different contributions for the most unstable TS mode computed at different values of $V_R$ and $d_\star$. Note that these profiles have been normalized to unit fluid kinetic energy $E_{FK}$. Figures 18(a) and 18(b), corresponding to $V_R = 0$ and $V_R = 0.03$ at $d_\star = 0$, show that the production term $\mathcal{P}$ is significantly modified by an increase of $V_R$. In particular, the amplitude of $\mathcal{P}$ decreases with $V_R$ and it exhibits a small region of negative production above the critical layer. In this region, the energy is transferred from the wave to the mean flow leading to a decrease of the total energy associated with the fluctuation. A similar observation is made by Metcalfe and Domaradski [29] for the case of a laminar boundary layer stabilized by a compliant membrane. Comparison of Figs. 18(a) and 18(b) also illustrates the importance of the viscous diffusion term $\mathcal{D}$ in redistributing energy produced by the Reynolds stress as $V_R$ is increased, whereas the pressure diffusion term $\Pi$ has a minor influence. This shows that under the action of the viscous diffusion term, the energy produced by the work of the Reynolds stress is transferred towards the wall where it is dissipated by viscosity ($\varepsilon$). Figure 18(a) also reveals that the production term $\mathcal{P}$ increases near the wall with the emergence of a second peak as $V_R$ is increased. Due to the pressure diffusion and viscous diffusion terms, which are negatively correlated with $\mathcal{P}$, this additional production does not result in a destabilizing effect.

The influence of the dissipation within the compliant wall is illustrated in Figs. 18(b)–(d) for $V_R = 1$. When $d_\star$ is increased from 0 [Fig. 18(b)] to 0.14 [Fig. 18(c)], the TS mode is destabilized (not shown here for the sake of conciseness). In Fig. 18(c), the pressure-diffusion term $\Pi$ is positive across the entire channel, which leads to work of the pressure force at the walls; the viscous diffusion profile $\mathcal{D}$ is mostly unchanged, but the peak near the wall of the production term $\mathcal{P}$ becomes the dominant feature. The increase of the dissipation leads to an increase in the production term in the viscous layer which is consistent with the TS mode (i.e., the phase shift mainly occurs in this region).

Finally, for a stronger wall damping $d_\star = 0.2$ [Fig. 18(d)], the pressure diffusion term $\Pi$ is almost identical to the production term $\mathcal{P}$, in contrast with what is observed for $d_\star = 0.14$ and $d_\star = 0$. In particular, $\Pi$ exhibits positive values near the wall which are associated with work of the pressure at the wall. This indicates that the mode is then in strong interaction with the wall and could probably be classified as a TWF mode. In particular, at $d_\star = 0.15$, a branch switching occurs and a collapse between TS and TWF modes is observed. This behavior is illustrated in Fig. 19. In the figure, $\alpha$ is chosen to maximize $\omega_i$ for each value of $d_\star$. On one hand, for $d_\star < 0.15$ where the TS mode exists, the production term $\hat{C}_1$ is increasing with the dissipation. It has as a consequence to increase the
<figure><img src="image_22.png" /><figcaption>FIG. 18. Fluid kinetic energy budget profiles for most unstable TS mode at Re = 10000, $B_\star = 1$, and (a) $d_\star = 0$, $V_R = 0.03$, (b) $d_\star = 0$, $V_R = 1$, (c) $d_\star = 0.14$, $V_R = 1$, and (d) $d_\star = 0.2$, $V_R = 1$. The wall-normal position of the critical layer is shown as a dashed line.</figcaption></figure>

temporal amplification rate of the TS mode in agreement with class A modes. On the other hand, for $d_\star > 0.15$, the temporal amplification rate is decreasing with the wall dissipation. This further indicates that beyond $d_\star = 0.15$, the mode is mainly associated with a TWF instability.

After the discussion of the TS modes, we now address the total energy budgets prevailing for divergence modes, of both sinuous and varicose symmetry. The evolution of the energy transfer mechanisms is monitored as $d_\star$ is varied, since this is the main control parameter influencing the dynamics of divergence modes. Figure 20 shows data computed over the range $0 < d_\star < 10$ at Re = 10000, $B_\star = 1$, and $V_R = 2$. The evolution with $d_\star$ of the different terms of the energy budget for sinuous [Fig. 20(a)] and varicose [Fig. 20(c)] modes reveals that the term $\hat{C}_2$ dominates for both symmetries. This term accounts for the energy exchange with the base flow due to wall displacement, which is therefore identified as the main mechanism promoting instability of the
<figure><img src="image_23.png" /><figcaption>FIG. 19. Energy of most unstable sinuous modes. Evolution of the integrated total energy budget with $d_*$ for $V_R = 1$. TS mode is shown as solid lines. Sinuous TWF mode is shown in dashed lines. At $d_* \approx 0.15$ modes collapse.</figcaption></figure>

divergence modes. It also confirms the influence of the Reynolds number onto the divergence mode as observed in the previous section. The importance of viscous effects for divergence modes has also been observed by Carpenter and Morris [37] for the boundary-layer flow case. The destabilizing effect of $\hat{C}_2$ is partially balanced by the work of the Reynolds stress against the basic shear $\hat{C}_1$, the viscous dissipation $\hat{C}_3$, and the wall dissipation $\hat{C}_4$. Hence, it illustrates the dual nature of viscous effects for the divergence mode. On one hand, it promotes the instability by propagating the

<figure><img src="image_24.png" /><figcaption>FIG. 20. Energy of most unstable divergence mode as $d_*$ is increased at $Re = 10000$, $B_* = 1$, and $V_R = 2$. (a, c) Evolution of integrated total energy budget and (b, d) breakdown of total energy into its components for (a, b) sinuous and (c, d) varicose modes.</figcaption></figure>
<figure><img src="image_25.png" /><figcaption>FIG. 21. Fluid kinetic energy budget profiles for most unstable divergence mode at $Re = 10000$, $B_* = 1$, and $V_R = 2$ and (a) $d_* = 5$ sinuous mode, (b) $d_* = 5$ varicose mode, and (c) $d_* = 9$ varicose mode.</figcaption></figure>

energy production from the wall displacement into the flow domain, but on the other hand, viscosity also plays its usual dissipative role. Inspection of the different contributions to the total energy [Figs. 20(b) and 20(d)] shows that while the fluid energy remains the main factor in both situations, the varicose modes involve significantly more wall energy than their sinuous counterparts. Due to the slow dynamics of these modes, the wall energy is almost entirely made up of potential energy for both symmetries. It should also be mentioned that due to the long wavelengths of the divergence modes, the main contribution to the wall potential energy $E_{WP}$ [Eq. (27)] is here due to the stiffness while the flexural rigidity only plays a marginal role.

The cross-channel profiles of the fluid kinetic energy budgets for the divergence modes are shown in Fig. 21. The figure demonstrates that the production, diffusion, and dissipation contributions are all localized in the near-wall region. Comparison of Figs. 21(a) and 21(b) reveals a notable difference between sinuous and varicose modes at $d_* = 5$: while the production term $\mathcal{P}$ is negative throughout the channel cross section for the varicose mode, the sinuous mode exhibits a small region near the wall with positive values of the production $\mathcal{P}$. Interestingly, when $d_*$ is further increased up to $d_* = 9$, the plots in Fig. 21(c) show that then the varicose production term $\mathcal{P}$ also exhibits a weakly positive region near the wall. This change of sign of the production term appears to approximately coincide with the onset of divergence instability, i.e., change of sign of the growth rate $\omega_i$ plotted in Fig. 20(c). Hence, while the energy budget is dominated by the term $\hat{C}_2$, it seems that the instability is also significantly influenced by $\mathcal{P}$.

Finally we consider TWF modes. For this final class of modes, we restrict our analysis to the varicose symmetry which has been observed to always dominate over the sinuous symmetry. The evolution of the total energy budgets with $V_R$ is shown in Fig. 22(a) for $Re = 10000$, $B_* = 1$, and $d_* = 0$. These plots indicate that the budget is driven by the production term $\hat{C}_1$, in agreement with previous analyses [20,29]. However, it is interesting to notice that energy transfer $\hat{C}_2$ from the base flow to the perturbation via the boundaries also promotes the instability. Both contributions $\hat{C}_1$ and $\hat{C}_2$ have a destabilizing influence. The distribution of the total energy among its different components plotted in Fig. 22(b) shows that for small values of $V_R$ the wall is the most energetic component of this fluid-structure system, but for $V_R > 1.6$ the kinetic energy of the fluid overcomes the wall contribution. Thus, the present scenario appears to be exactly opposite of the situation prevailing for TS modes: as $V_R$ is increased, the ratio $E_{FK}/E_W$ decreases for TS modes while it increases for TWF modes. This is in accordance with the classification of TS modes as class A and TWF modes as class B modes. Finally, monitoring the potential and kinetic components of the wall energy ($E_{WP}$ and $E_{WK}$, respectively) shows that the wall energy is essentially due to the kinematic contribution for TWF modes, which is in contrast to the situation prevailing for divergence modes.

In Fig. 22(c), the fluid kinetic energy profile is shown across the channel diameter for $V_R = 3$. The role of viscosity is found to be mainly concentrated in two distinct areas. The viscous diffusion term $\mathcal{D}$ exhibits a first shallow peak around the critical layer (indicated by the dashed line), while
<figure><img src="image_26.png" /><figcaption>FIG. 22. Energy of most unstable varicose TWF mode as $V_R$ is increased at $Re = 10000$, $B_* = 1$, and $d_* = 0$. (a) Evolution of integrated total energy budget. (b) Distribution of total energy as fluid and wall contributions. (c) Cross-channel profile of the fluid kinetic energy budget for $V_R = 3$. The wall-normal position of the critical layer is shown as a dashed line.</figcaption></figure>

a second stronger peak emerges near the wall due to the viscous layer. This is consistent with the theoretical model developed by Davies and Carpenter [20]. However, the dominant production term $\mathcal{P}$ is associated with the work of the Reynolds stress against the mean shear and does not exhibit a maximum near the critical layer. Its profile presents rather an inflection point near the critical layer, which is in contrast with observations by Metcalfe and Domaradski [29] for the flat plate boundary layer. We also observe that both pressure diffusion $\Pi$ and viscous diffusion $\mathcal{D}$ display large positive values near the wall. Thus the energy production due to the basic shear flow is transferred to the wall through the action of both diffusion processes.

To conclude this section on energy transfer mechanisms, the influence of wall dissipation $d_*$ is reported in Fig. 23 for the most unstable varicose TWF modes at $Re = 10000$, $B_* = 1$, and $V_R = 3$. The plots in Fig. 23(a) show that for small values of $d_*$ the stabilization of the TWF mode is essentially due to the increasing energy dissipation $\hat{C}_4$ in the compliant walls. For larger values of $d_*$, both the work of the Reynolds stress against the basic shear $\hat{C}_1$ and the wall dissipation $\hat{C}_4$ evolve so as to stabilize the TWF mode. This illustrates that two distinct mechanisms are at play to diminish the growth rate for the TWF mode. The curves in Fig. 23(b) illustrate the redistribution of the total energy due to wall dissipation: as $d_*$ is increased, the wall contribution (essentially kinetic energy) to the total energy decreases almost linearly in favor of the fluid contribution. The effect of the wall damping on the wall-normal profiles of the fluid kinetic energy is shown in Fig. 23(c). While these profiles at $d_* = 4$ are similar to those prevailing without wall dissipation [see Fig. 22(c)], it is observed that the production $\mathcal{P}$ exhibits a lower amplitude than the pressure diffusion term $\Pi$. As $d_*$ increases, the energy transfer from the fluid towards the wall is no longer sufficient to balance the higher wall dissipation and, as a result, the proportion of fluid energy increases in the total energy of the perturbation.
<figure><img src="image_27.png" /><figcaption>FIG. 23. Energy of most unstable varicose TWF mode as $d_*$ is increased at $Re = 10000$, $B_* = 1$, and $V_R = 3$. (a) Evolution of integrated total energy budget. (b) Distribution of total energy as fluid and wall contributions. (c) Cross-channel profile of the fluid kinetic energy budget for $d_* = 4$. The wall-normal position of the critical layer is shown as a dashed-line.</figcaption></figure>

V. CONCLUDING REMARKS

In this paper, we have revisited the linear stability problem of a fluid interacting with a compliant channel. The walls are modeled as spring-backed deformable plates including a damping mechanism. A general numerical method free of spurious modes is derived to tackle this problem. By taking advantage of the base flow symmetries, varicose and sinuous eigenmodes are computed separately and efficiently. A dimensional analysis has been carried out to identify the most physically relevant control parameters. Traveling wave flutter (TWF), divergence, and Tollmien-Schlichting (TS) modes are then recovered. The main conclusions of the current study are as follows:

(1) TWF waves are mainly driven by the reduced velocity $V_R$ which measures the strength of the coupling between the fluid and the compliant walls. The Reynolds number effect is not significant for these modes.

(2) For TWF modes, the perturbations of varicose symmetry are observed to be destabilized first when increasing $V_R$. While a decrease in the flexural rigidity is accompanied by a stabilization of the TS mode, it enhances the amplification of TWF modes. This is especially true for the varicose TWF modes. In particular, we show that the compliant wall potential for laminar-turbulent transition delay is almost negligible when the varicose symmetry is taken into account. This finding extends the previous linear stability results obtained for the same flow case but either limited to the sinuous modes [20,21] or focusing on a quite narrow region of parameter space [22].

(3) For the varicose symmetry, a linear relationship is found between the critical reduced velocity and the wall dissipation, independently of the flexural rigidity. For large values of $V_R$, the dynamics is eventually dominated by divergence modes for a certain amount of wall dissipation and for both symmetries. While a stabilizing effect is observed for the divergence mode when increasing the Reynolds number, sinuous perturbations are amplified at a lower wall dissipation
than their varicose counterpart. Such a distinction between varicose and sinuous symmetries was not identified in previous studies.

We have also uncovered physical mechanisms responsible for the amplification of the different categories of modes. This was accomplished by evaluating total energy budgets. First, for $d_{\star} = 0$, kinetic energy budgets indicate that the TS mode is stabilized as $V_R$ is increasing under the action of two mechanisms. On one hand, the work of the Reynolds stress against the mean flow exhibits a region above the critical layer which is negative. Here, energy is transferred from the wave to the base flow. The second mechanism is associated with an irreversible energy transfer from the fluctuation to the mean flow due to the mean flow displacement near the walls. On the other hand, the vertical distribution of the production term associated with the TWF mode is observed to peak near the walls and exhibits an inflection point close to the critical layer for both symmetries. Hence, it shows some differences with results provided by Domaradzki and Metcalfe [29] for the boundary layer where a peak is observed near the critical layer position. Furthermore, the energy transfer from the wall to the fluid is mainly attributed to a pressure diffusion term.

Second, the dominant term for divergence modes appears to be associated with an irreversible energy exchange due to the interaction of the displaced mean flow and shear stress.

Extension of the present study to nonmodal stability analyses can be considered in a future work. Hœpfner et al. [35] have investigated the same problem for three-dimensional perturbations but numerical oscillations were observed by the authors when increasing the number of eigenmodes used to compute transient growth. The formulation adopted here is free of spurious modes and we believe that the problem can be solved using the numerical procedure developed in the present analysis. In addition, the role of the reduced velocity has not been investigated by Hœpfner et al. [35]. Finally, it should also be interesting to extend our study to pulsatile flow [38,39] and the pipe geometry, which cover more biologically significant phenomena.

## ACKNOWLEDGMENTS

Fédération lyonnaise de modélisation et sciences numériques (FLMSN) is gratefully acknowledged for providing access to high-performance computing resources.

## APPENDIX: VALIDATION AND INFLUENCE OF BOUNDARY CONDITIONS

For validation purposes, we adopt choices made by Davies and Carpenter [20] for dimensionless wall equations. Here, the Reynolds number is based on $h$, $v$, and the centerline velocity. The dimensionless wall parameters are

$$ \Gamma = \frac{m}{\rho h}, \quad B_{\star} = \frac{B}{h \rho v^2}, \quad K_{\star} = \frac{Kh^3}{\rho v^2}, \quad d_{\star} = \frac{dh}{\rho v}. $$

Davies and Carpenter [20] used an Orr-Sommerfeld equation for solving the fluid-structure interaction problem. The pressure at the lower wall is recovered through

$$ p(-1) = \frac{1}{2} \left( -i\omega \int_{-1}^{1} \hat{v} \, dy + i\alpha \int_{-1}^{1} U_b \hat{v} \, dy \right), $$

where the viscous terms have been neglected. For consistency, the authors have also neglected the effect of the viscous stress at the wall for the coupling between the fluid and the walls. We recall that this formulation is correct only for the sinuous configuration. Guaus and Bottaro [28] used a primitive formulation for the same problem and considered only the pressure at the wall for the force associated with the fluid acting at the walls. However, the pressure is computed without neglecting the viscous stress along the wall-normal direction. In the present study, we have not neglected the effect of the viscous stress on both the pressure and the force acting on the walls. In Fig. 24, we show the neutral curve for the sinuous TWF mode for $K_{\star} = 10^7$, $B_{\star} = 4K$, $d_{\star} = 0$, and $\Gamma = 2$. The figure shows an almost perfect agreement between our approach and the one given by Davies
<figure><img src="image_28.png" /></figure>

FIG. 24. Neutral curves. Sinuous configuration. Effects of the boundary conditions. $K_{\star} = 10^{7}$, $B_{\star} = 4K_{\star}$, $d_{\star} = 0$, and $\Gamma = 2$. DC, results extracted from [20]; BC1, the formulation adopted in the present paper; BC2, the pressure at the wall is expressed as in [20] and we neglect the viscous stress at the walls; BC3, formulation adopted in [28].

and Carpenter [20]. It supports the hypotheses made by the previous authors. Nevertheless, our approach is more general because it also allows to deal with the varicose symmetry. In addition, the system also reduces to an Orr-Sommerfeld problem with only one velocity component. Figure 24 also shows that the formulation made by Guaus and Bottaro [28] is not consistent and leads to discrepancies near the critical Reynolds number.

[1] J. B. Grotberg and O. E. Jensen, Biofluid mechanics in flexible tubes, Annu. Rev. Fluid Mech. 36, 121 (2004).

[2] P. G. Larose and J. B. Grotberg, Flutter and long-wave instabilities in compliant channels conveying developing flows, J. Fluid Mech. 331, 37 (1997).

[3] J. Gray, Studies in animal locomotion. VI. The propulsive powers of the dolphin, J. Exp. Biol. 13, 192 (1936).

[4] M. Kramer, Boundary layer stabilization by distributed damping, J. Aeronaut. Sci. 24, 459 (1957).

[5] D. M. Bushnell and K. J. Moore, Drag reduction in nature, Annu. Rev. Fluid Mech. 23, 65 (1991).

[6] F. E. Fish and G. V. Lauder, Passive and active flow control by swimming fishes and mammals, Annu. Rev. Fluid Mech. 38, 193 (2006).

[7] M. Gad-el-Hak, Boundary layer interactions with compliant coatings: An overview, Appl. Mech. Rev. 39, 511 (1986).

[8] M. Gaster, Is the dolphin a red herring? in Turbulence Management and Relaminarisation (Springer-Verlag, Berlin, 1988), pp. 285–304.

[9] V. Kumaran, Stability and the transition to turbulence in the flow through conduits with compliant walls, J. Fluid Mech. 924, P1 (2021).

[10] R. Patne and V. Shankar, Stability of flow through deformable channels and tubes; implications of consistent formulation, J. Fluid Mech. 860, 837 (2019).

[11] J. H. Duncan, A comparison of wave propagation on the surfaces of simple membrane walls and elastic coatings bounded by a fluid flow, J. Sound Vib. 119, 565 (1987).

[12] B. T. Benjamin, Shearing flow over a wavy boundary, J. Fluid Mech. 6, 161 (1959).

[13] B. T. Benjamin, Effects of a flexible boundary on hydrodynamics stability, J. Fluid Mech. 9, 513 (1960).

[14] B. T. Benjamin, The threefold classification of unstable disturbances in flexible surfaces bounding inviscid flows, J. Fluid Mech. 16, 436 (1963).
[15] J. W. Miles, On the generation of surface waves by shear flows, J. Fluid Mech. 3, 185 (1957).

[16] P. W. Carpenter and A. D. Garrad, The hydrodynamic stability of flow over Kramer-type compliant surfaces. Part 1. Tollmien-Schlichting instabilities, J. Fluid Mech. 155, 465 (1985).

[17] P. W. Carpenter and A. D. Garrad, The hydrodynamic stability of flow over Kramer-type compliant surfaces. Part 2. Flow-induced surface instabilities, J. Fluid Mech. 170, 199 (1986).

[18] O. Wipler and U. Ehrenstein, Numerical simulation of linear and nonlinear disturbance evolution in a boundary layer with compliant walls, J. Fluids Struct. 14, 157 (2000).

[19] M. Malik, M. Skote, and R. Bouffanais, Growth mechanisms of perturbations in boundary layers over a compliant wall, Phys. Rev. Fluids 3, 013903 (2018).

[20] C. Davies and P. W. Carpenter, Instabilities in a plane channel flow between compliant walls, J. Fluid Mech. 352, 205 (1997).

[21] C. Davies and P. W. Carpenter, Numerical simulation of the evolution of Tollmien-Schlichting waves over finite compliant panels, J. Fluid Mech. 335, 361 (1997).

[22] M. Nagata and T. R. Cole, On the stability of plane Poiseuille flow between compliant boundaries, in *Proceedings of the International Conference on Computational Methods and Experimental Measurements IX*, edited by G. M. Carlomagno, University of Naples di Napoli, Italy and C.A. Brebbia, Wessex Institute of Technology, United Kingdom (WIT Press, 1999), pp. 231–240.

[23] K. Tsigklifis and A. D. Lucey, Asymptotic stability and transient growth in pulsatile Poiseuille flow through a compliant channel, J. Fluid Mech. 820, 370 (2017).

[24] J. S. B. Gajjar and P. Sibanda, The hydrodynamic stability of channel flow with compliant boundaries, Theor. Comput. Fluid Dyn. 8, 105 (1996).

[25] E. de Langre, *Fluides et Solides* (École polytechnique, 91120 Palaiseau).

[26] M. L. Manning, B. Bamieh, and J. M. Carlson, Descriptor approach for eliminating spurious eigenvalues in hydrodynamic equations, arXiv:0705.1542.

[27] J. M. Rotenberry and P. G. Saffman, Effect of compliant boundaries on weakly nonlinear shear waves in channel flow, SIAM J. Appl. Math. 50, 361 (1990).

[28] A. Guaus and A. Bottaro, Instabilities of the flow in a curved channel with compliant walls, Proc. R. Soc. A 463, 2201 (2007).

[29] J. A. Domaradzki and R. W. Metcalfe, Stabilization of laminar boundary layers by compliant membranes, Phys. Fluids 30, 695 (1987).

[30] V. Shankar and V. Kumaran, Stability of wall modes in fluid flow past a flexible surface, Phys. Fluids 14, 2324 (2002).

[31] R. Peyret, *Spectral Methods for Incompressible Viscous Flow* (Springer, Berlin, 2002).

[32] L. Trefethen, *Spectral Methods in MATLAB* (Society for Industrial and Applied Mathematics, Philadelphia, 2000).

[33] P. J. Schmid and D. S. Henningson, *Stability and Transition in Shear Flows*, Applied Mathematical Sciences Vol. 142 (Springer, Berlin, 2001).

[34] P. G. Drazin and W. H. Reid, *Hydrodynamic Stability* (Cambridge University Press, Cambridge, U.K., 1981).

[35] J. Hœpffner, A. Bottaro, and J. Favier, Mechanisms of non-modal energy amplification in channel flow between compliant walls, J. Fluid Mech. 642, 489 (2010).

[36] L. Huang, Reversal of the Bernoulli effect and channel flutter, J. Fluids Struct. 12, 131 (1998).

[37] P. W. Carpenter and P. J. Morris, The effect of anisotropic wall compliance on boundary-layer stability and transition, J. Fluid Mech. 218, 171 (1990).

[38] B. Pier and P. J. Schmid, Linear and nonlinear dynamics of pulsatile channel flow, J. Fluid Mech. 815, 435 (2017).

[39] B. Pier and P. J. Schmid, Optimal energy growth in pulsatile channel and pipe flows, J. Fluid Mech. 926, A11 (2021).