Benchmark solution for the stability of plane Couette flow with net throughflow

This paper investigates the stability of an incompressible viscous fluid flow between relatively moving horizontal parallel plates in the presence of a uniform vertical throughflow. A linear stability analysis has been performed by employing the method of normal modes and the resulting stability equation is solved numerically using the Chebyshev collocation method. Contrary to the stability of plane Couette flow (PCF) to small disturbances for all values of the Reynolds number in the absence of vertical throughflow, it is found that PCF becomes unstable owing to the change in the sign of growth rate depending on the magnitude of throughflow. The critical Reynolds number triggering the instability is computed for different values of throughflow dependent Reynolds number and it is shown that throughflow instills both stabilizing and destabilizing effect on the base flow. It is seen that the direction of throughflow has no influence on the stability of fluid flow. A comparative study between plane Poiseuille flow and PCF has also been carried out and the similarities and differences are highlighted.

stabilizes the basic channel flow. Later, Fransson and Alfredsson 21 dealt with the same problem with a slightly different approach, wherein they considered the stability of a shear flow allowing it to be both a function of the mean velocity distribution and the Reynolds number separately. In particular, they used the maximal channel velocity as their velocity scale, and thus separated the effects of base velocity magnitude from those of the base velocity distribution. Interestingly, they demonstrated that the effect of throughflow could be either stabilizing or destabilizing depending on the values of governing parameters. Guha and Frigaard 22 presented a linear stability of the PCPF in the presence of throughflow. They observed that the destabilization of the basic flow as the throughflow velocity is increased and then stabilization at large values of throughflow velocity and also, the energy analysis was executed to decipher the mechanism of instability. The stability of PPF with net throughflow and an adverse temperature gradient was studied by Bajaj 23 and showed that the effect of throughflow is to stabilize or destabilize the base flow depending upon its magnitude as well as direction and the Prandtl number. Recently, Shankar and Shivakumara 24 investigated the effect of a uniform vertical throughflow on the stability of porous-Poiseuille flow and found that the throughflow dependent Reynolds number instills both stabilizing and destabilizing effects on the base flow.
To the best of our knowledge, the effect of a uniform vertical throughflow on the hydrodynamic stability of PCF against small disturbances has not been investigated so far, and the objective of this investigation is to furnish the missing information. It is well known that the PCF is always stable for all the values of Reynolds number so far as small disturbances are concerned with no throughflow. At first glance, it is natural to assume that the PCF in the presence of throughflow will also be stable against small disturbances for all the values of Reynolds number. However, the present investigation has revealed that this conjecture is not always true because the basic velocity profile is distorted by the effect of throughflow and permits the unstable mode of disturbance. The stability of fluid flow is analyzed by solving the stability eigenvalue problem numerically using the Chebyshev collocation method.

Formulation of the problem
We are considering an incompressible viscous fluid bounded by two horizontal rigid, however permeable, plates which are moving in the opposite directions with a uniform speed U 0 and in the presence of a uniform vertical throughflow of strength W 0 . The flow is considered to be only due to the motion of the plates as well as throughflow and not due to the applied pressure gradient. The sketch in Fig. 1 is three-dimensional wherein the vertical coordinate z spans the limited interval [−h, h] , while the ( x , y)-coordinates are unbounded and varies over all possible real values. A dimensionless form of the governing equations, given by the local mass balance and the local momentum balance, can be formulated as (Drazin and Reid 3 ) where the coordinates (x, y, z) are non-dimensionalized in the usual way by the half-width of the channel h , the fluid velocity � q = ( u, v, w ) by U 0 , the pressure p by ρU 2 0 and the time t by h/U 0 , where U 0 is the velocity at the center line of the channel and ρ is the fluid density. Here, R = U 0 h/ν is the Reynolds number, while ν is the kinematic viscosity.
The boundaries are rigid but permeable and are moving with a uniform speed in the opposite directions. Hence, the formulation of the problem sketched in Fig. 1 leads to the dimensionless boundary conditions

Base flow
The base flow is fully developed, unidirectional, steady and laminar with a constant vertical throughflow of mag- The plots of u b (z) as a function of z for different values of R T are displayed in Fig. 2. In the absence of throughflow, the basic velocity reduces to the well-known linear profile as expected i.e. lim R T →0 u b = z (Drazin and Reid 3 ). The presence of vertical throughflow, however, changes the basic flow drastically, which has important implications on the stability of fluid flow. It is seen from the figure that an increase in the value of |R T | is to compress the basic velocity profiles toward the boundaries to which the flow is directed. In particular, a very thin boundary layer structure develops close to z = ±1 2 when |R T | is very large. Figure 2 shows an evident symmetry between positive and negative values of R T , when the absolute value of this parameter is fixed. In fact, there is symmetry of solution (5) yielding an invariance under the transformation

Linear stability analysis
To discuss the linear stability of fluid flow, the field variables are split into a basic state and an infinitesimal disturbance as where the primes denote the perturbed quantities. Substituting Eq. (7) into Eqs. (1) and (2) and neglecting the nonlinear terms, we obtain www.nature.com/scientificreports/ The stability problem involves the growth of arbitrary infinitesimal perturbations, but it has been generally assumed that such perturbations can be resolved into independent wave-like components. Each component is supposed to satisfy the linearized equations of motion and boundary conditions separately. So we consider a typical wave component with where α and β are the wave numbers in the x and y-directions, respectively and c = c r + ic i is the complex wave speed. It is to be understood that the real part has to be taken to get physical quantities, this being permissible for the linear problem i.e. the flow is stable, neutrally stable or unstable for c i < 0, c i = 0 or c i > 0 , respectively. Equation (10) serves to separate the variables and reduce the linearized governing Eqs. (8) and (9) from partial to ordinary differential equations and the pressure is eliminated from the momentum equation so that we are led to the following stability equation where D = d/dz . The term with coefficient R T has arisen due to the vertical transport of perturbation velocity. The appropriate boundary conditions are Equations (11) and (12) forms an eigenvalue problem which is symmetric, on account of Eq. (5) and it is left invariant by the transformation

Computation of the eigenvalue
The resulting eigenvalue problem can be solved with high numerical accuracy by employing the Chebyshev collocation method. This mathematical technique is well documented in the literature (Canuto et al. 25 ). In the collocation method, the kth order Chebyshev polynomial is and the Chebyshev collocation points are where j = 0 and j = N correspond to the upper and lower wall boundaries, respectively. The field variable w is approximated in terms of Chebyshev polynomials as where w j is a constant. Using (16), the stability Eqs. (11) and (12) are discretized in terms of Chebyshev polynomials to get z j = cos(πj/N), j = 0, 1, 2, . . . , N,

Numerical results and discussion
The  Table 1 and note that the numerical accuracy improves with increasing N in Eq. (16). It can also be pointed out from the table that the numerical scheme needs a good number of collocation points to achieve the convergence with increasing R T . Growth rate. The sign of c i yields significant information about the stability characteristics of the fluid flow.
To accomplish this, the growth rate c i versus the wave number α for different values of R and |R T | is illustrated in Fig. 3a-c. In the absence of throughflow (i.e. R T = 0 ), it is seen that c i is negative for all values of the Reynolds number considered indicating that no linear perturbation mode grows with time so that the basic flow is always linearly stable (Fig. 3a) confirming the result documented in Drazin and Reid 3 . For R T = 10 , however, it is observed that the sign of c i changes from negative to positive signaling the flow switching over from stable to unstable depending on the value of R and the same is evident from Fig. 3b. A closer inspection of this figure Table 2. Comparison of (R c , α c , c c ) between the present study and Fransson and Alfredsson 21 for different values of R T .   www.nature.com/scientificreports/ enclosed by the corresponding neutral stability curve. The region below and outside the curve is that defining the parametric conditions of linear stability (c i < 0) . As a consequence, the absolute minimum of the neutral stability curve yields the so-called critical value of α, R and c which are denoted by α c , R c and c c , respectively. The physical meaning of the critical conditions is that no linear instability is possible for R < R c . From the figures, it is intriguing to note that the two-dimensional disturbances ( β = 0 ) are not the first to go unstable for all values of α . A closer look at the figures, it is evident that three-dimensional disturbances ( β = 0 ) become unstable earlier to that of β = 0 at lower range of values of α (Fig. 4b). Nonetheless, the observed trend reverses with increasing α and thereby the critical Reynolds number at which the instability triggers first turns out to be the least for only two-dimensional disturbances (i.e. validity of Squire's theorem). The reason for this conflicting behavior may be attributed to the fact that we are not plotting the neutral curves against the overall wave number of the roll and probably one should consider the orientations of the roll as well to overcome this inconsistency. Figure 5 shows the neutral stability curves for different values of |R T | with β = 0 and it is evident that the vertical throughflow ensures both stabilizing and destabilizing effect. This means that the neutral stability curve moves both downward and upward when |R T | increases so as the pushing of instability region. It is also observed that the normal modes activating the instability correspond to increasing wave numbers as R T increases, which physically means very small wavelengths. Moreover, the sign of critical wave speed c c remains to be positive (negative) for all positive (negative) values of R T indicating the perturbation modes are travelling.

Stability of PCF.
The critical triplets (R c , α c , c c ) for PCF are plotted in Fig. 6 as a function of R T and the results obtained for PPF are also shown for the sake of comparison. The direction of throughflow does not alter R c (Fig. 6a, Table 3) and hence it has no impact on the stability of fluid flow. From Fig. 6a, it is evident that there exists a threshold value of |R T |= 3.353 prior to which the flow is always stable because the perturbations exhibit a negative growth rate as the throughflow effect is not strong enough to induce instability. Beyond this threshold value, however, a self-excited mode of disturbance exists for PCF; a result of contrast observed in the absence of throughflow (Drazin and Reid 3 ). This may be due to the occurrence of Tollmien-Schlichting wave with increasing |R T | . Moreover, the value of R c decreases rapidly and assumes a minimum value of 167,087.609 at |R T | = 4.899 . That is, there exists a range of |R T |(3.353 ≤ |R T | ≤ 4.899 ) in which R c decreases (i.e. flow is destabilizing) but R c starts increasing (i.e. flow is stabilizing) again with further increase in |R T |(> 4.899) and attains an asymptotic value 27,189.31 |R T | for |R T | > 25 . This may be due to the transport of momentum that opposes the stabilizing effect of viscous diffusion, and this is large when |R T | increases, and hence, higher values of R c are needed to initiate instability.
The corresponding values of critical wave number α c for PCF are presented in Fig. 6b. The critical wave number is independent of the direction of throughflow and an increasing function of |R T | . In particular, α c increases from |R T |= 3.353 and shows a steady rise in their values, eventually attaining a limiting value α c = 0.155469|R T | for sufficiently large values of |R T | . We have shown that the principal effect of the vertical throughflow is to modify the steady velocity distribution, which tends to exponential boundary-layer type with thickness of the order of α |R T | with increasing |R T | . Therefore, we may expect that the critical wavelength will be of the same order of  www.nature.com/scientificreports/ magnitude as the boundary-layer thickness. Thus the critical wave number will increase with increasing |R T | indicating the effect of throughflow is to reduce the size of convection cells. The curves of critical wave speed c c start from the value of |R T | beyond which the flow becomes unstable (Fig. 6c). Initially, the curve of c c exhibits a decreasing trend and thereafter remains invariant with further increase in |R T | . The sign of c c becomes negative (positive) when |R T | is negative (positive) but its modulus does not change. Thus the cells move up (down) the layer when the solutions correspond to c c < 0 (c c > 0).

Comparative study between PPF and PCF.
In the case of PPF both the plates are at rest and an external pressure gradient drives the flow in the horizontal direction whereas in PCF both the plates are considered to be moving in the opposite directions with a uniform speed U 0 , and the relative motion drives the shearing action in the fluid between the plates. As a result, the basic velocity for these two flows is found to be different due to asymmetric boundary conditions. However, the resulting stability equation is the same for both PPF and PCF even in the presence of throughflow. In the case of no throughflow ( R T = 0), the basic state becomes unstable when the Reynolds number is larger than 5772.22 (Orszag 1 ) for PPF whereas in the case of PCF, the basic flow is always stable (Drazin and Reid 3 ). It is observed from Fig. 6a that the critical Reynolds number of PPF lies well below the critical Reynolds number of PCF for all values of |R T | , indicating PCF has a more stabilizing effect than PPF. The vertical throughflow, irrespective of its direction, imparts only stabilizing effect on PPF, in general. In contrast, the throughflow indicts both destabilizing and stabilizing effect on PCF. The critical wave numbers for PPF are higher than PCF, indicating that the size of the cell width of PPF is smaller than that of PCF throughout the domain of |R T | (Fig. 6b). The critical wave speed values of PPF are lower than those of PCF for upward throughflow, while an opposite trend is observed for downward throughflow (Fig. 6c). Although the sign of c c remains positive for all values of |R T | for PPF, its sign changes for PCF and it depends on the direction of throughflow as well. Figure 6b,c also show that there is an abrupt variation in their values in the case of PPF and this may be due to the occurrence of bi-modal structure within a limited range of |R T | in the neutral stability curve. In the beginning, the left minimum (lower wave number) determines the critical mode; however, a slight increase in |R T | from 1.29 to 1.3 makes the right minimum (higher wave number) to be lower than the left one to dominate the instability. The switching over of instability mode is the cause for this sudden change in these critical values. While in the case of PCF uni-modal structure is observed and hence there is no change in the mode of instability.

Conclusions
The change in the hydrodynamic stability of plane Couette flow (PCF) due to the presence of a uniform vertical throughflow has been investigated. A modal stability analysis of small-amplitude disturbances has been performed and the resulting stability eigenvalue problem is solved numerically using the Chebyshev collocation method. Some interesting and unexpected results are unveiled on the stability characteristics of the fluid flow. The numerical approach has shown that that three-dimensional disturbances are always more stable than twodimensional ones. Suitable symmetries in the stability eigenvalue problem is found allowing one to gather information on the regime |R T | < 0 from the results of |R T | > 0 . It is found that PCF is always stable in the absence of throughflow and also for values of |R T | ≤ 3.352 and beyond this value it becomes unstable as the growth rate turns out to be positive from negative. In the parameter space 3.353 ≤ |R T | ≤ 4.899 , the flow gets destabilized manifesting itself as a minimum in the (R T , R c )-plane and stabilized reversely for values of |R T | > 4.899 and attains an asymptotic value 27,189.31 |R T | for sufficiently large values of |R T | . Besides, the instability is found to set in always via travelling-wave mode and PCF has a more stabilizing effect than PPF in the presence of vertical throughflow. The size of convection cells decrease with increasing |R T | and the cells may move up ( c c < 0 ) or down (c c > 0) the layer depending on the direction of throughflow.