Pattern preservation during the decay and growth of localized wave packet
in two-dimensional channel flow

Localized structures are found to be the key characteristics of subcritical transitions in wall-bounded shear flows [1, 2, 3], e.g. puffs in pipe flows and turbulent bands in plane Poiseuille flow (PPF), plane Couette flow (PCF), and circular Couette flow. For the two-dimensional PPF, the subcritical transition may occur at Reynolds numbers ($Re$) larger than 89, the threshold for energy instability [4, 5], and less than 5772, the critical value of linear instability [6], and the counterpart of turbulent band is localized Tollmien-Schlichting waves or localized wave packet (LWP) [7, 8, 9, 10]. A possible localization mechanism is the subharmonic instability [11, 12] of upper-branch Tollmien-Schlichting waves (TSW) [13, 14], which do not exist for Reynolds numbers less than 2600 [15, 11, 16]. By applying the continuation of invariant solution in a periodic domain with a periodicity length $L$ about 16.7$\pi$, the modulated Tollmien-Schlichting waves (MTSWs) generated through Hopf bifurcation from TSW is tracked to a critical Reynolds number $Re_{c}$ about 2285 [12]. When the periodicity length is increased to $32\pi$ , the minimal Reynolds number traced by the continuation method for LWP is found to be 2335 [17], above which the lower branch is unstable, corresponding to the edge state [18], and the upper branch undergoes a Hopf bifurcation at a higher Reynolds number. With time-dependent numerical simulations in a long channel ($L=20\pi$), single localized wave packet was obtained with $Re_{c}\sim 2330$ [8], indicating a saddle-node bifurcation. Such kind of bifurcation is also found in PCF [19] and pipe flow [20]. Later simulations of two-dimensional PPF reveal $Re_{c}$ of 2332 with finer resolutions and much longer channel lengths, e.g. $L=1600$ [21] and $L=400$ [22], indicating that the periodic MTSW obtained previously as $Re>2332$ will evolve to LWP when the channel is long enough.

An important issue for subcritical transition is the persistency of localized structures, which can be determined by analyzing their lifetimes $\tau$. For PCF, it is shown experimentally that the lifetime of turbulent patches scales as $\tau\propto(Re_{c}-Re)^{-1}$ with $Re_{c}\approx 323\pm 2$ [23], which agrees with the numerical observation $Re_{c}=324\pm 1$ [24]. The same scaling relation was found for puffs in the pipe flows by experiments [25] and numerical simulations [26], representing that the localized turbulence can sustain when $Re>Re_{c}$. Later experiments with longer observation time showed that the mean lifetimes of puff before decaying and splitting are superexponential functions of $Re$, indicating that puffs are transient phenomena at moderate Reynolds numbers [27, 28]. The lifetime statistics of turbulent bands in tilt and narrow computational domains have been studied as well, and the Reynolds numbers determined by the crossover of the mean decay time and splitting time for PPF and PCF are 965 [29] and 325 [30], respectively. For large channel domains, however, it is illustrated numerically that the transition starts with a sparse turbulence state, where the localized turbulent bands are sustained by a dynamic balance between band extension and band breaking [31, 32]. The lower bound or threshold Reynolds number of the sparse turbulence state is $Re_{c}=660$, which is confirmed experimentally by velocity field measurements [33, 34] and flow visualizations [35]. Above $Re\approx 700$, most turbulent bands are found numerically to have the same orientation and the spanwise symmetry is restored and the two-dimensional directed-percolation is retrieved only above $Re\approx 1000$ [36]. For two-dimensional PPF, however, the scaling law of LWP’s lifetime with $Re$ hasn’t been studied so far.

The spatio-temporal evolution of localized flow structures, e.g. decay and growth, is another important feature of the subcritical transitions. For PCF, the decay of turbulent band was studied in large domains with under-resolution simulations, and it was found that the turbulent band shrinks slowly and laminar holes are formed in the bands [37, 38]. At moderate Reynolds numbers, direct numerical simulations illustrate that an isolated turbulent band decays in a style of longitudinal contraction with a statistically constant velocity, and experiences transient growths before the eventual relaminarization [39]. In PPF, turbulent band breaking occurs and the fragments relax to streaks before vanishing [31]. The growth of transitional characteristic structures has a plenty of scenarios. An isolated turbulent band in PCF may experience oblique extension, transverse split, and longitudinal split at moderate Reynolds numbers [40], and their combined effects lead to a labyrinthine growth [41], which can also be triggered directly by localized perturbations [24]. For the turbulent bands in PPF, the oblique extension happens at moderate Reynolds numbers, and the parallel or longitudinal split and the transverse split occur at relatively high Reynolds numbers [31, 36]. In the two-dimensional PPF, however, it is still unknown whether LWPs evolve in similar ways to the three-dimensional cases during the decay and growth processes.

This paper is organized as follows. In Sec. II the numerical methods are described, and the resolution and domain size independence are verified for LWP. Different decaying periods are revealed at low Reynolds in Sec. III, and the lifetime scaling law is explained theoretically with a pattern preservation approximation and the disturbance kinetic energy equation in Sec. IV. The growth of LWP from the lower branch or edge state to the upper branch is studied numerically and compared with theoretical predictions in Sec. V, and some concluding remarks are provided in Sec. VI.

The incompressible two-dimensional plane Poiseuille flow is numerically simulated with a spectral code [42]. The flow rate is kept constant, and the half channel height $h$ and 1.5 times of the bulk velocity $U_{m}$ are chosen as the characteristic length and velocity scales, respectively. No-slip boundary conditions at the side walls ($y=\pm 1$) and periodic boundary condition in the streamwise direction $x$ are applied. $N_{x}$ Fourier modes and $N_{y}$ Chebyshev modes are used in the $x$ and $y$ directions, respectively. The Reynolds number is defined as $Re=1.5U_{m}h/\nu$, where $\nu$ is the kinematic viscosity of the fluid. For details of the simulation methods, we refer to the previous papers [31, 32].

Following the previous studies [22], the length $l_{p}$, the center coordinate $x_{c}$, and the convection velocity $c_{p}$ of LWP are defined as

where $e_{k}$ is the kinetic energy of disturbance velocity ${\bf{u^{\prime}}}=(u^{\prime},v^{\prime})=(u-U_{0},v^{\prime})$ and $U_{0}(y)=1-y^{2}$ is the basic flow solution.

The computational domain length $L$ is just the space interval between neighbor LWPs due to the periodic boundary condition, and its effects on $l_{p}$ are shown in Fig. 1, where the same initial field shown in Fig. 1(a) is used. With the increase of $L$ and the spatial resolution ($N_{x}\times N_{y}$), $l_{p}$ converges at different Reynolds numbers [see Fig. 1(c)], reflecting the independence of LWP on the computational domain size and resolution. In the previous studies, ($N_{x}\times N_{y}$)=$(200\times 33)$ with $L=20\pi$ [8] and $(256\times 33)$ per 100 units length [21] were used. According to Fig. 1, larger $L$ and higher spatial resolution are required to obtain converged solutions at higher $Re$s. For example, $l_{p}$ is about 33 at $Re=4500$, and $L=100$ is so short that the interaction among the LWP and its neighbors leads to decay, a similar phenomenon is observed for turbulent bands in three dimensional channel flows [32]. Consequently, $L=400$ and $(N_{x}\times N_{y})=(2048\times 65)$ are used hereafter. In a translational coordinate moving with the LWP’s velocity $c_{p}$, the time series of normal velocities $v^{\prime}$ are recorded at different locations shown by the symbols in Fig. 1(b) , and the corresponding power spectra are calculated. These spectra share the same dominant frequency [see Fig. 1(d)], which is identified as a global frequency $\omega_{g}$ by spatio-temporal stability analyses [22] and increases mildly with $Re$ in a similar manner to $c_{p}$ as shown in the inset.

In order to study the decay properties and lifetime characteristics of LWP, the LWPs obtained at $Re=$ 2335, 2400, and $3000$ are used as the initial disturbances of simulations within the range of $Re=2298\sim 2333$. The total disturbance kinetic energy $E_{k}=\int e_{k}dxdy$ is calculated, and the lifetime $\tau$ of each case is defined as the period with $E_{k}\geq 10^{-5}$. Though the initial disturbance kinetic energies are different, it is shown in Fig. 2 that the temporal evolutions of $E_{k}$ illustrate some common features. Firstly, $E_{k}$ obtained at $Re=2333$ reaches a constant value eventually, while other cases at lower $Re$s decay monotonically: the higher $Re$ is, the longer the LWP lifetime is. Secondly, for LWPs with high initial $E_{k}$, e.g. the cases shown in Fig. 2(b), there are three characteristic periods along the $E_{k}$ evolution curves, i.e. the initial steep descent period, the middle plateau period, and the final steep descent period. The last two periods are characterized by the first and the second inflection points labeled by the square and the triangle symbols in Figs. 2(a) and 2(b), respectively. $E_{k1}$ and $E_{k2}$ are the corresponding disturbance kinetic energies. Finally, the lifetime increases faster than the superexponential scaling when $Re$ is close to a critical value $Re_{c}$ indicated by a vertical dash dot line in Fig. 2(c). Note that the lifetimes of transient puffs in pipe flows scale superexponentially with $Re$ [27]. As shown in Fig. 2(d), $\tau^{-2}$ of different initial disturbances decreases linearly and can be extrapolated to a value $Re_{c}=2332.35\pm 0.35$, i.e. $\tau^{-2}\propto(Re_{c}-Re)$ or $\tau\propto(Re_{c}-Re)^{-1/2}$, indicating that the lifetime tends to infinity when $Re$ increases to $Re_{c}$ and LWP becomes a temporally persistent structure as $Re>Re_{c}$. This critical Reynolds number is consistent with the previous values obtained with initial random disturbances [21] and localized disturbances [22], where LWPs are arranged initially with large intervals.

The lifetimes of LWPs change for different Reynolds numbers and initial disturbances, but the common features shown in Fig. 2 suggest that the decay process is completely different from the chaotic behaviors of puffs in pipe flows. It is shown in Figs. 3(a-d) that the flow patterns obtained at different decay moments look similar, i.e. localized wave packets whose streamwise wavelengths do not change much with time but whose perturbation amplitudes and kinetic energy decrease monotonically. In addition, $E_{k}$ at the inflection points on the $E_{k}$ curves are nearly the same for different Reynolds numbers and different initial perturbation energies [see Figs. 2(a) and 2(b)]. Therefore, LWPs seem to follow the same route to decay, preserving the flow pattern but fading the kinetic energy, and hence a pattern preservation approximation is proposed: for an evolving LWP, the domain integral properties of the velocity field do not depend on the Reynolds number and the initial perturbation strength, but are single valued functions of $E_{k}$.

Applying the no-slip boundary conditions on the sidewalls and periodic conditions in the streamwise direction, the growth rate $\sigma$ of $E_{k}$ is defined according to the Reynolds-Orr equation,

where $\mathcal{P}$ and $\mathcal{D}$ are the kinetic energy production term and the enstrophy of the disturbance field, respectively, and $\mathcal{D}/Re$ represents the viscous dissipation.

Along the decay processes, both $\mathcal{P}$ and $\mathcal{D}$ are calculated for different Reynolds numbers and different initial kinetic energies and are shown in Figs. 3(e) and 3(f), respectively. The curves collapse with each other, indicating that both $\mathcal{P}$ and $\mathcal{D}$ are single-valued, monotonically increasing functions of $E_{k}$. Consequently, the pattern preservation approximation is validated. It is noted that the $\mathcal{D}$ and $\mathcal{P}$ curves shown in Fig. 3 have small-amplitude oscillations, which are related to the periodic production of travelling waves at the downstream end of LWP [22] during the decay process. In order to eliminate the interference of these oscillations, the growth rate $\sigma$, $\mathcal{D}/\mathcal{P}$, and $d\mathcal{D}/d\mathcal{P}$ are calculated with smoothed $E_{k}$ curves hereafter, where the oscillation components have been filtered out.

When $Re$ is close to $Re_{c}$, or $Re=Re_{c}(1-\delta)$ with $0<\delta\ll 1$, we have

If a LWP is saturated at $Re_{sat}$, we have $\sigma=0$ and $\mathcal{D}/\mathcal{P}=Re_{sat}$ according to Eq. (2). Considering that $\sigma=0$ at $Re_{c}$ and $Re_{c}$ is the minimum of $Re_{sat}$, where $dRe_{sat}/dE_{k}=d(\mathcal{D}/\mathcal{P})/dE_{k}=0$, or $d\mathcal{D}/d\mathcal{P}=\mathcal{D}/\mathcal{P}$, the subcritical transition threshold $Re_{c}$ can be derived analytically based on the $Re$-independent properties ($\mathcal{D}$ and $\mathcal{P}$) as

As shown in Fig. 3(g), the crossover point of $d\mathcal{D}/d\mathcal{P}$ (thick solid line) and $\mathcal{D}/\mathcal{P}$ (thick dotted line) is consistent with $Re_{c}$ obtained in simulations [see Fig. 2(d)].

When $Re<Re_{c}$, $dE_{k}/dt<0$ and the inflection points in the $E_{k}$ curves require

Considering that $d\mathcal{D}/d\mathcal{P}$ only depends on $E_{k}$, Eq. (5) suggests that the $E_{k}$s at the inflection points ($E_{k1}$ and $E_{k2}$) are nearly constants, a feature confirmed by the numerical data shown in Figs. 2(a) and 2(b). As illustrated by the vertical dashed line in Fig. 3(g), $E_{k1}\simeq E_{kc}$, a result consistent with Eqs. (4) and (5). Consequently, the growth rate at the first inflection point is

a relation confirmed by simulations shown in Fig. 4(a).

Considering the growth rate at the second inflection point $\sigma_{2}=\mathcal{P}(E_{k2})-\mathcal{D}(E_{k2})/Re$, a linear relation between $\sigma_{2}$ and $1/Re$ is expected and is confirmed by the numerical data shown in Fig. 4(a), where the fitted dashed line represents $\sigma_{2}=0.000382-1.288/Re$ and the coefficients agree reasonably with the values of $\mathcal{P}(E_{k2})$ and $-\mathcal{D}(E_{k2})$, 0.000304 and $-1.02$, respectively.

When $Re$ approaches to $Re_{c}$, $\sigma$ in Eq. (3) is approximately a single valued function of $E_{k}$ and can be expanded around the first inflection point as

where $A$ is approximated as a positive constant when $1/Re\simeq 1/Re_{c}$ is used. Note that at the inflection point, $(d\sigma/dE_{k})\mid_{E_{k}=E_{k1}}=[(d\sigma/dt)/(dE_{k}/dt)]\mid_{E_{k}=E_{k1}}=0$. When $\delta\rightarrow 0$, substituting Eq. (6) into Eq. (7) leads to a simplified equation for the reduced kinetic energy $E_{k}-E_{k1}$,

the classical partial differential equation for saddle-node bifurcation, i.e. when $Re<Re_{c}$ there are only decaying solutions for $E_{k}-E_{k1}$, and when $Re>Re_{c}$ there is a stable upper branch solution and an unstable lower-branch one. After rescaling the time and the reduced kinetic energy with LWP’s lifetime scale $\hat{\tau}$ as $t=\hat{\tau}\tilde{t}$ and $E_{k}-E_{k1}=\tilde{E_{k}}/(A\hat{\tau})$, Eq. (8) can be transformed into the unified form $d\tilde{E}_{k}/d\tilde{t}=-1-\tilde{E}_{k}^{2}$ for the decay cases by setting $1/\hat{\tau}^{2}=(Re_{c}-Re)A\mathcal{D}(E_{kc})/Re_{c}^{2}$, i.e. $\hat{\tau}\propto(Re_{c}-Re)^{-1/2}$, the lifetime scaling law found in the simulations shown in Fig. 2(d). Equation (8) suggests that there will be long-lasting transient process (the middle plateau period) and the lifetime will diverge when the system is close to the saddle-node bifurcation, manifesting a phenomenon known as “critical slowing-down” [43]. Similar energy plateau periods have been reported for PCF [44].

Since $E_{k1}\sim 10^{-1}\ll 1$ and $|d^{3}\sigma/dE_{k}^{3}|\ll|d^{2}\sigma/dE_{k}^{2}|$ at $E_{k}=E_{k1}$ due to the smooth feature of $\sigma$ curves shown in Fig. 4(b), $|E_{k}-E_{k1}|$ needs not to be much smaller than $E_{k1}$ in order to ignore the higher order terms in Eq. (7). It is shown in Fig. 4(b) that all the decay rate differences $-(\sigma-\sigma_{1})$ at different Reynolds numbers with different initial perturbations ($Re_{i}$) collapse together and agree well with the quadratic relation Eq. (7) in the range of $0.5E_{k1}<E_{k}<2E_{k1}$, which covers nearly the whole lifetime of LWP as shown in Fig. 4(c). Consequently, the lifetime for a case with the initial disturbance kinetic energy $E_{k}^{i}$ can be evaluated with Eq. (7) as

Note that $E_{k1}/\sqrt{-\sigma_{1}}\simeq E_{k1}Re_{c}/\sqrt{(Re_{c}-Re)\mathcal{D}(E_{k1})}\gg 1$ as $Re$ approaches $Re_{c}$. Specifically,

Considering that $Re\sim 10^{3}$ and $\delta\ll 1$, the condition $E_{k}^{i}-E_{k1}\gg\sqrt{\mathcal{D}(E_{k1})\delta/(ARe)}$ is easy to be satisfied for $E_{k}^{i}\neq E_{k1}$ when $\delta$ is small enough. According to Eq. (10), the lifetime of LWP with high $E_{k}^{i}$ should be about twice as long as that with $E_{k}^{i}=E_{k1}$, i.e. $\tau/\tau_{1}\simeq 2$, a prediction consistent with the simulations shown in Fig. 4(d), indicating that the initial disturbance kinetic energy does affect the lifetime of LWP.

When $Re>Re_{c}$, it is shown in Fig. 5(a) that the weak LWPs, which are obtained numerically at different times during a decay process, grow directly to a saturated state with the same $E_{k,\mathrm{UB}}$, the corresponding upper-branch LWP solution. Especially, when $E_{k}$ of the initial LWP is lower than a specific value $E_{k,\mathrm{LB}}$, the LWP will decay as shown by dotted curves in Fig. 5(a), where the bisection method [18] is used to track the lower-branch or edge state. The lower and upper branches of LWP calculated by the time-dependent simulations are shown by the empty and filled circles in Fig. 5(b), respectively, where the solutions obtained by the continuation method [17] are shown as references.

Both branches of LWP are steady solutions, and according to Eq. (2) the corresponding $\sigma=0$ and $\mathcal{D}/\mathcal{P}=Re$. The horizontal dashed line in Fig. 3(g) is tangent to the $\mathcal{D}/\mathcal{P}$ curve, indicating the minimum $Re$ for steady LWP. Moving this line up will intersect with the $\mathcal{D}/\mathcal{P}$ curve at two points, whose $E_{k}$s should correspond to $E_{k,\mathrm{LB}}$ and $E_{k,\mathrm{UB}}$, respectively. Based on the $\mathcal{D}$ and $\mathcal{P}$ obtained during the decay process ($Re<Re_{c}$), the predicted $E_{k}$s for the upper and the lower branches are shown by the grey curve in Fig. 5(b), agreeing with the simulation data. Especially, the predicted $E_{k,\mathrm{LB}}$ is consistent with the simulation value not only near $Re_{c}$ but also at high Reynolds numbers, e.g. $Re=4000$, indicating that the pattern preservation approximation is applicable for LWPs with low $E_{k}$s in a wide range of $Re$. The flow patterns at different times during a growth process from the lower branch to the upper branch are shown in Figs. 5(c-f) and share the same features, localized short-wavelength waves with a strong downstream end and a long and decaying upstream tail.

The kinetic energy production term $\mathcal{P}$ and the disturbance enstrophy $\mathcal{D}$ during the growth processes at different Reynolds numbers are shown in Figs. 6(a) and 6(b), respectively. The curves for moderate Reynolds numbers agree acceptably with each other, confirming the availability of the pattern preservation approximation. When both the Reynolds number and the disturbance kinetic energy are high, e.g. $E_{k}>0.15$ and $Re>2500$, discrepancies among the $\mathcal{P}$ and $\mathcal{D}$ curves of growing cases can be observed, reflecting the Reynolds number effect on strong LWPs and the limitation of the approximation. When the decay process starts from a LWP saturated at a higher $Re$, there is always an offset adaptation. This offset effect becomes significant as $E_{k}>0.15$, and leads to the underestimation of $\mathcal{P}$ and $\mathcal{D}$ [Fig. 6(a) and 6(b)], the dependence of decay rate on the initial $E_{k}(Re_{i})$ [Fig. 4(b)], and the underprediction of the upper branch with high $E_{k}$ as shown in Fig. 5(b).

When LWP grows from the lower branch to the upper branch, there is one and only one inflection point in the $E_{k}$ curve, corresponding to a disturbance kinetic energy $E_{k3}$ with the largest growth rate $\sigma_{3}$ as shown in Fig. 5(a). This phenomenon can be explained briefly with the pattern preservation approximation as follows. At the inflection point, Eq. (5) is still applicable, i.e. $d\mathcal{D}/d\mathcal{P}=Re$, and according to Fig. 3(g) there is only one intersection point between a horizontal line for $Re>Re_{c}$ and the $d\mathcal{D}/d\mathcal{P}$ curve in the region above the $\mathcal{D}/\mathcal{P}$ curve, which confines the $E_{k}$ range between the lower and the upper branches. Using the similar arguments to those for Eq. (7), the growth rate can be expanded around this inflection point as well and we have $\sigma\simeq\sigma_{3}+B(E_{k}-E_{k3})^{2}$, where $B$ may be further expanded as $B=(1/2)(d^{2}\sigma/dE_{k}^{2})_{E_{k}=E_{k3}}=(1/2)(d^{2}\sigma/dE_{k}^{2})_{E_{k}=E_{kc}}+...$. As shown in Fig. 6(c), the growth rate differences $\sigma-\sigma_{3}$ obtained numerically at different Reynolds numbers coincide with each other and agree with the quadratic relation. The fitted coefficient $B=-0.0496$, which is very close to -0.048, the value of $(1/2)(d^{2}\sigma/dE_{k}^{2})_{E_{k}=E_{kc}}$.

According to the above discussions, just based on the pattern preservation approximation and the $\mathcal{D}$ and $\mathcal{P}$ data obtained from the decaying cases, both $Re$ and $E_{k}$ of the whole lower branch, the turning point, and the upper branch with $E_{k}<0.15$ of this saddle-node bifurcation can be predicted directly by the disturbance kinetic energy equation. The reason why the flow pattern can be preserved during the decay and growth processes lies in the nonlinear self-sustaining mechanism of LWP [22]. In the mean-flow field after deducting the basic flow, there is a finite-amplitude vortex dipole at the downstream end of LWP, which provides an unstable region to produce TS-type waves propagating and decaying in the upstream stable region, forming the long upstream tail. On the other hand, the Reynolds stresses contributed by the TS-type waves strengthen the vortex dipole and hence LWP can be self-sustained. Since the length of vortex dipole is confined by the channel height and the wavelength of the TS-type wave in the long upstream tail are mainly determined by the dispersion relation of the basic flow, i.e. both length scales are intrinsically related to the channel height, the geometric shape of LWP weakly depends on the Reynolds number, and its pattern can be retained.

Near the onset of turbulence in three-dimensional PCF and PPF, the growth of disturbance kinetic energy may be caused by band extension and band split, and the decrease of $E_{k}$ may correspond to band breaking or longitudinal contraction. For two-dimensional PPF, however, the variation of $E_{k}$ does not change much the flow pattern of LWP. By applying the pattern preservation approximation, where the integral kinematic properties are independent of $Re$ and only functions of $E_{k}$, it is shown that the disturbance kinetic energy equation can be transformed into the classical differential equation for saddle-node bifurcation. As a result, LWP’s lifetime can be derived analytically and is shown to be proportional to $(Re_{c}-Re)^{-1/2}$, a scaling law found in numerical simulations. Especially, the Reynolds number and $E_{k}$ of both the upper branch ($E_{k}<0.15$) and the lower branch predicted by the disturbance kinetic energy equation are consistent with the simulation data in a wide range of Reynolds number, indicating that the pattern preservation is an intrinsic character of LWP. Considering that the two-dimensional PPF is a simplified model for the three-dimensional case and retains required ingredients for the subcritical transition, e.g. solitary localized structure at moderate Reynolds numbers, structure split at relative high Reynolds numbers, and small-scale quasi-periodic components coupled with a large-scale mean flow, the present results are expected to be helpful in understanding the fully three-dimensional transitions in channel flows and wait for experimental validations through film or other two-dimensional flow techniques.