Project Topics
Get Complete Project Material File(s) Now! »
Three-dimensional stability analysis
Mathematical model
We consider an incompressible swirling flow of a Newtonian fluid in a rotating, straight, finite-length circular pipe. The axial and radial distances are scaled with the pipe radius R such that the pipe non-dimensional radius is unity and the pipe non-dimensional length is L. The flow domain is given by VC = {(x, r, θ ) : 1 ≤ x ≤ L, 0 ≤ r ≤ 1, 0 ≤ θ < 2π}. The velocity components are scaled with the characteristic axial speed entering the pipe, denoted as U. The gauge pressure from the chamber-base reference pressure is scaled by the characteristic dynamic pressure ρU2 and is denoted as p. The time t is scaled with the ratio of pipe radius to inlet characteristic speed. The flow inside the domain VC is described by the non-dimensional continuity and axial, radial, and circumferential momentum equations (refer to (2.2)).
youThis system is subjected to certain boundary conditions that may reflect a physical setting of a flow in a straight circular rotating pipe that is generated by a vortex generator at a steady, continuous and smooth operation ahead of the pipe inlet. In the present study we focus on a special setting where the inlet profile is a steady solid-body rotating flow with a uniform axial velocity. The non-dimensional axial, radial and azimuthal velocities at the pipe inlet, x = 0, are fixed for all time t ≥ 0, Here ω ≥ 0 is the inlet-flow swirl ratio. These inlet conditions may represent the nature of a swirling flow generated by a rotating honeycomb ahead of the rotating pipe (see Leclaire Sipp [26]). Slightly different inlet conditions that may represent the case of a vortex generator that is built of guiding vanes have been studied in Wang & Rusak ([66], [68]). The flow is treated as being locally inviscid at the pipe wall and at the pipe outlet, see for example, Fletcher [16] for a discussion of boundary conditions. The boundary condition at the pipe outlet, x = L, is described by a zero axial gradient of the pressure for all time t ≥ 0, This condition is supported by the experimental data shown in Umeh et al. [63] where the mean swirling flow has established a parallel flow state after several radii from the pipe inlet. The pipe wall is modelled as an inviscid, slip boundary. At the pipe wall, r = 1, we impose the no-penetration condition for all time t ≥ 0, Beran [4] and Lopez [35] have used a similar wall condition in their numerical study of viscous swirling flows in pipes. Under the boundary conditions specified above, the solid-body rotation flow with a uniform axial velocity and no radial velocity, is a steady-flow solution of the problem for every swirl ratio ω and all Re; this is the base flow for the present linear stability analysis. The base-flow non-dimensional velocity vector and gauge pressure are denoted u0 = (1, 0, ωr) and p0(r) = ω2r2/2.
Substituting (3.5) into the governing equations (2.2) leads, at order O(ε), to the equations of the linearized flow dynamics inside the domain VC for all time t ≥ 0, We consider a suitable general mode analysis of the linearized flow problem in a straight circular pipe that is composed of an azimuthal normal mode in θ (since the problem preserves azimuthal symmetry) and a non-axial Fourier mode in x and r that must obey the non-periodic inlet-outlet conditions. This mode is given in the form:
Here σ may be a complex number that relates to the growth rate of the perturbation, i is the imaginary unit and m is an integer representing the perturbation’s azimuthal wavenumber. The perturbation admits a separation of variables solution in terms of (x, r), θ and t. Substituting the general mode shape (3.7) into (3.6a)-(3.6d) results for each m in The system (3.8) forms the flow linear stability equations of a mode with an azimuthal wavenumber m. Since the base flow is axisymmetric, this approach simplifies the perturbation equations from three-dimensional to two-dimensional in terms of the radial and axial coordinates, (x, r). These equations are subjected to the following boundary conditions for every azimuthal wavenumber m (refer to figure 3.1), In addition, the perturbation mode functions W , U , V and P must satisfy the compatibility conditions at the pipe centreline due to the singular nature of the coordinate system, see Ash & Khorrami [1], In addition, the perturbation mode functions W , U , V and P must satisfy the compatibility conditions at the pipe centreline due to the singular nature of the coordinate system, see Ash & Khorrami [1], The perturbation equations (3.8) with the boundary conditions (3.10a) – (3.10c) and the compatibility conditions (3.11a) – (3.11c) form the flow linear stability problem for the solu-tion of the eigenvalue σ and the mode’s shape functions {W (x, r), U (x, r), V (x, r), P(x, r)} as a function of the swirl ratio ω, the pipe non-dimensional length L and the perturbation’s azimuthal wavenumber m. For each m and every ω, L and Re, this problem has a continuum of solutions of perturbation modes, where some have real values of σ and shape functions and others have complex values of σ and shape functions. For each m, we are seeking the range of ω where the real part σr of σ becomes positive, indicating an unstable perturbation mode. Also, for each m, we denote ωm,k as a critical swirl level where a certain mode k of perturbation changes its stability.
Remark 1: In the inviscid flow case, the perturbation equations (3.8) are also the governing equations for the classical stability problem where a vortex flow in an infinitely-long pipe is considered (Kelvin [24]). This set of equations together with the imposed far-field boundedness condition of the perturbations at ±∞ and the wall conditions (3.10c) form a well-posed eigenvalue problem and thus uniquely determine all the eigenmodes and their eigenvalues. The normal modes used in the classical linear stability theory are but the only possible solutions of this eigenvalue problem. Then, the perturbation mode functions W , U , V and P must be each written as multiples of an axial Fourier mode, ei(αx+mθ )+σt (where α is a real axial wavenumber, m is an integer azimuthal wavenumber, and σ is a complex growth rate) and a radial mode shape that depends only on the radial distance r. This forms a complete separation-of-variables solution of the perturbation in x, r, θ and t, separately. In this way, the solution of the perturbation modes is even further simplified for all m to the solution of a system of ordinary differential equations in r alone, with only the centreline (r = 0) compatibility and wall (r = 1) conditions. These eigenvalue problems determine the radial mode shapes and the growth rates of the normal modes. However, the classical normal mode solutions do not satisfy the present boundary conditions (3.10a)-(3.10c) and cannot be used in the current analysis which forms a different problem. More general mode shapes that depend on both (x, r) (with no separation-of-variables shapes) are necessary in the present analysis.
Remark 2: Applying a complex-conjugate operation on the stability equations (3.8) and conditions (3.10a) – (3.10c) and (3.11a) – (3.11c) shows that the growth rates and eigenmodes for a case with a fixed m < 0 are respectively the complex-conjugate growth rates and eigenmodes for a case with m > 0. Therefore, in the present analysis, for all m, the perturbation modes with m < 0 are identical to their counterparts with m > 0. This is in contrast to the classical normal mode analysis of the stability problem, where the growth rates and eigenmodes for each m < 0 are different from those for m > 0. In the classical stability problem, the perturbations to a solid-body rotation flow are rotating and travelling waves with different rotational and axial speeds that are free to move along the pipe, either downstream, in the direction of the base axial flow or upstream, with a different speed, against the direction of the base axial flow. In the present flow set-up, the perturbations are rotating waves about the centreline in the direction of the base flow rotation, but are subjected to the relatively-fixed inlet and outlet conditions that modify the axial travelling-wave nature of the perturbations. For example, in the axisymmetric (m = 0) case, the perturbations are described by a standing growing wave in a definite range of ω around the first critical swirlω0,1.
Remark 3: In the inviscid flow problem with a similar physical set-up and boundary conditions, the special case of m = 0 was initially introduced and studied by Wang & Rusak
They found that the linear stability nature of a vortex in a finite-length pipe is strikingly different from the predictions of classical linear stability theory based on the normal mode analyses by Kelvin [24], Rayleigh [44], Synge [58] and Howard & Gupta [21]. This reflects the physical fact that a vortex flow in a finite-length pipe is sensitive to the confined conditions at the pipe inlet and outlet. As a result, the set of small-disturbance eigenmodes that govern the flow is completely different from the classical normal modes. At high swirl ratios ω, this leads to a flow instability when ω is above the first critical swirl ω0,1. The m = 0 unstable modes have been extensively studied by Wang & Rusak ([66], [67], [68]), Gallaire & Chomaz [17], Vyazmina et al. [64], Meliga & Gallaire [40], Wang & Rusak [70], Rusak et al. [47], Wang et al. [71], Rusak & Wang [46], due to the fact that they are closely related to the axisymmetric vortex breakdown process.
Remark 4: We note that the m = 0 case is special in its mathematical treatment; it also admits a separation-of-variables solution that is different from the classical normal modes, further simplifies the analysis and leads to an ordinary differential equation that can be solved using an accurate Runge-Kutta integration method. However, the cases with m = 1, 2, …
modes do not admit such a similar approach; therefore a numerical solution of the linear flow stability problem is an appropriate approach. The goal of this article is to extend the m = 0 stability theory to the cases with a general m mode and to complete the stability analysis of both inviscid and viscous solid-body rotation flows in a straight, circular finite-length pipe.
Numerical method and validation
The stability equations (3.8) with the boundary conditions (3.10a) – (3.10c) and the com-patibility conditions (3.11a) – (3.11c) are solved numerically by a finite-difference method. The computational domain 0 ≤ x ≤ L and 0 ≤ r ≤ 1 is discretized by a uniform mesh with M ×N grid points in the x and r directions, respectively. A second-order central-difference scheme is used in discretizing the first- and second-order partial derivatives in x and r. The obtained discretized equations form an eigenvalue-problem for the solution of σ in terms of the base flow swirl ratio ω, Re, L and azimuthal wavenumber m. A standard Matlab eigenvalue problem solver is used to determine the eigenvalues and eigenfunctions.
In this paper we show results of stability computations for modes with m = 0, m = 1 and m = 2 in the range of swirl ratios ω where these modes become unstable as swirl ratio is increased. Similar computations can be conducted for higher-order modes with |m| ≥ 3. We compare the computed results and determine the competition between the various perturbation modes.
In order to assess the accuracy of our numerical scheme, we first compare the present computed results for the benchmark inviscid flow case where m = 0 with results using the well-established method for the stability characteristics of axisymmetric modes developed in Wang & Rusak [66]; see additional results in Wang et al. [71]. We note that in the Wang
Rusak [66] model the inlet conditions considered a fixed vanishing azimuthal vorticity and the outlet condition considered a vanishing radial velocity, whereas in the present study a vanishing inlet radial velocity condition and a vanishing outlet axial pressure gradient are applied. We used the same streamfunction vorticity solver to solve the m = 0 problem. The ordinary differential equation (ODE) eigenvalue problem with the present boundary conditions has been solved by a shooting iterative method to machine numerical accuracy. These accurate growth rates are shown in figure 3.2 (the solid line). Also shown in figure 3.2 are the presently-computed results of σ of the least stable mode with m = 0 in the range of swirl ratios between ω = 1.932 and ω = 1.986 (the open circles). In this range of swirl ratio, of the least-stable mode is real. It is evident that the present computations using a mesh size of 600 ×250 give accurate results of the growth rates for m = 0.
In order to validate the computed results, we also conducted a thorough mesh convergence study by using various sizes of meshes for the modes m = 0, 1, 2 in both the inviscid and viscous flow cases. Tables 3.1 and 3.2 compare for various meshes the computed results of σ of the most-unstable linear mode of the solid-body rotation flow with m = 1 for a representative inviscid flow case with ω = 2.1 and L = 6. From table 3.1 we find that the computed σ change by only ∼ 10−6 when the mesh radial resolution is increased from
N = 100 to N = 300 with a fixed M = 1000. This demonstrates that a mesh with N = 100 in the radial direction is sufficient for an accurate computation of σ. From table 3.2 we find that the computed σ change by ∼ 6 ×10−4 when the mesh axial resolution is increased from M = 200 to M = 1600 with a fixed N = 100 (approximately 2.4% error in computing the real part of σ and less than 0.1% error in computing the imaginary part of σ). This demonstrates that a mesh of M ×N = 800 ×100 is sufficient for an accurate computation of
of the various m = 1 modes in the inviscid flow case. Similarly, we find that a mesh of M ×N = 800 ×100 is also sufficient for an accurate computation of σ of the various m = 0 and m = 2 modes in the inviscid flow case.
Tables 3.3 and 3.4 compare for various meshes the computed results of σ of the most-unstable linear mode with m = 0 for a representative viscous flow case with ω = 2.15, Re = 1000 and L = 6. These tables demonstrate that a mesh of M ×N = 800 ×100 is again sufficient for an accurate computation of σ of the m = 0 modes in the viscous flow case.
Similarly, tables 3.5 and 3.6 compare for various meshes the computed results of the growth rate of the most-unstable linear mode with m = 1 for a representative viscous flow case with ω = 2.1, Re = 1000 and L = 6. These tables demonstrate again that a mesh of M ×N = 800 ×100 is sufficient for an accurate computation (within 1% or less) of σ of the various m = 1 modes in the viscous-flow case. Comparing the results for same ω = 2.1 from the inviscid-flow case in table 3.2 and from the Re = 1000 case in table 3.6, we find that the real part of σ decreases with the addition of the viscous effect while the imaginary part stays nearly the same.
Comparison for various meshes for the computed results of the growth rate of the mostunstable linear mode with m = 2 for a representative viscous-flow case with ω = 3.64, Re = 500 and L = 6 are shown in table 3.7 and table 3.8. From table 3.7 we find that the growth rate is changed by 2 ×10−4 as the radial mesh resolution is increased from N = 100 to N = 300 with a fixed M = 1000, and from table 3.8 we find the computed growth rate is changed by ∼ 4 × 10−4 when the axial mesh resolution is increased from M = 400 to M = 800 with a fixed N = 100. This demonstrates that a mesh with M ×N = 800 ×100 is sufficient for an accurate computation of the growth rates of the m = 2 modes in the viscous flow. Moreover, in later chapters, we will further demonstrate that the viscous effect is rather passive and nearly the same for all ω in the range of studied swirl ratios.
Tables 3.9 and 3.10 compare for various meshes the computed results of σ of the most-unstable linear mode of the solid-body rotation flow with m = 0 for a representative inviscid flow case with ω = 2.35 and L = 2. From table 3.9 we find that the computed σ change by only ∼ 10−5 when the mesh radial resolution is increased from N = 100 to N = 300 with a fixed M = 1000. This demonstrates that a mesh with N = 100 in the radial direction is sufficient for an accurate computation of σ. From table 3.10 we find that the computed change by ∼ 3 × 10−5 when the mesh axial resolution is increased from M = 200 to M = 1600 with a fixed N = 100. This demonstrates that a mesh of M = 800 × N = 100 is sufficient for an accurate computation of σ of the various m = 0 modes in the viscous-flow case.
Table 3.11 and table 3.12 compare for various meshes the computed results of σ of the most-unstable linear mode with m = 1 for a representative viscous flow case with ω = 2.3, Re = 500 and L = 2. The computed σ change by ∼ 1 ×10−4 when the mesh radial resolution is increased from N = 100 to N = 300 with a fixed M = 1000 as shown in table 3.11.
The inviscid-flow perturbation modes on a solid-body rotation flow
The growth rates and shapes of modes of the inviscid solid-body rotation flow are presented here. For each perturbation mode m = 0, m = 1 and m = 2, the growth rate σr as a function of the swirl ratio ω is first given, followed by the velocity vector field and the streaklines of the mode.
The m=0 mode
We first consider the axisymmetric (m = 0) perturbation modes. A similar linear stability problem was first studied analytically by Wang & Rusak [66]. The present stability equations (3.8) result in a system of linear partial differential equations for the solution of the eigenvalue. The problem has a continuum of solutions of perturbation modes that turn in sequence from being asymptotically stable to unstable as ω is increased above the first critical swirl ω0,1, where ω0,1 = ωB2 + π2/(4L2) (here ωB = 1.9159 is Benjamin’s [3] special swirl). The growth rates σr (the real part of σ) as a function of the swirl ratio ω of the various m = 0 perturbation modes on the solid-body rotation flow in a pipe with L = 6 are shown in figure 3.3.
Fig. 3.3 The growth rates σr (the real part of σ) as a function of the swirl ratio ω of the various inviscid axisymmetric (m = 0) perturbation modes on a solid-body rotation flow in a pipe with L = 6.
It is found that that all modes have real values of σ that are negative and real shape functions when 0 < ω < ω0,1 = 1.9337, i.e., all the m = 0 modes are asymptotically stable when the swirl ratio is below the first m = 0 critical swirl. The σ of the least-stable mode is real and changes from negative to positive, and the flow becomes unstable as ω is increased from below the first critical swirl ω0,1 = 1.9337 to above it. The growth rate σr of the k = 1 mode increases in the range ω0,1 = 1.9337 < ω < 1.971, reaches a maximum value of 0.00475 at ω = 1.971 and then decreases in the range 1.971 < ω < ω0,2
The growth rate σr of the second (k = 2) mode changes from negative to zero as ω increases to the second critical swirl ω0,2 = 1.9862. We note that the real growth rates of the two least-stable modes merge to form the first fold point at ω0,2 = 1.9862, where for both modes σr = 0. All the other perturbation modes have real negative σ and real shape functions and are asymptotically stable in the range ω0,1 < ω < ω0,2. The two least-stable modes combine together at ω0,2 = 1.9862 to form a pair of modes with complex-conjugate values of σ and mode shape functions in the range ω > ω0,2. Only the real part σr of these two modes is shown in figure 2; it increases in the range ω0,2 = 1.9852 < ω < 2.045, reaches a maximum of 0.0095 at ω = 2.045, and then deceases in the range 2.045 < ω < ω0,3,q where ω0,3 = ωB2 + 9π2/(4L2) = 2.0706. These two complex-conjugate modes have an oscillatory behaviour of σr as ω increases above ω0,3 = 2.0706.
The σr of the third (k = 3) least-stable mode changes from negative to positive as ω increases from below the third critical swirl ω0,3 = 2.0706 to above it. The growth rate and shown that the critical swirl levels with even values of n are also fold points of the σr plots, and complex-conjugate modes appear at these fold swirl levels. Also, the critical swirl levels with odd values of n are points of exchange of stability of a mode with real σ and shape functions. The dominant unstable m = 0 perturbation at each swirl ratio is the perturbation mode with the greatest value of σr at that swirl level.
Figures 3.4 and 3.5 show a representative example of the solution of the most-unstable m = 0 (axisymmetric) mode shape functions at ω = 2.13 (between ω0,3 and ω0,4) in a pipe with L = 6 for which σ = 0.025 (on solution branch k = 3 of m = 0 case). Figure 3.4 presents the computed planar projected velocity vector, (1 + εW , εU ), of the mode in the flow domain; here W and U are real and ε is the perturbation’s amplitude such that there is a near-stagnation point at the pipe outlet. Figure 3.5 describes a snapshot of the streaklines of the mode consisting of particles released at the pipe inlet (x = 0) on a circle with r = 0.15. It can be seen that the mode describes an axisymmetric axial growing wave with a wavelength of approximately 3.6 and an amplitude that grows in x. Similar figures can be given for all other computed m = 0 modes.
The m=1 mode
In parallel to the axisymmetric (m = 0) modes, the stability problem with m = 1 also has a continuum of solutions of spiral perturbation modes that are studied in this paper for the first time. The growth rates σr (the real part of σ) and frequency σi (the imaginary part of σ) as a function of the swirl ratio ω of the various inviscid m = 1 perturbation modes on the solid-body rotation flow in a pipe with L = 6 are shown in figure 3.6. Each of the different modes is characterized by a continuous solution as a function of ω, with a similar axial wavelength of the shape functions that increases with the mode number k. The modes are numbered in ascending numerical order k = 1, k = 2, and so on. Figure 3.6a describes σr of the modes k = 1 and k = 2 in the range 1.5 ≤ ω ≤ 2.2 and figure 3.6b describes the corresponding σi of these modes. Figure 3.6c describes σr of the modes k = 3 and k = 4 in the range 1.9 ≤ ω ≤ 2.5 and figure 3.6d describes the corresponding σi of these modes.
It is found that each one of these complex modes turns in sequence from being asymptot-ically stable to unstable as ω is increased. Mode k = 1 changes its stability at the swirl level ω1,1 = 1.6134 and reaches a local maximum positive value of σr = 0.000253 at ω = 1.629 (figure 3.6a). Then σr of this mode continues with oscillations as ω is further increased and has maximum points of σr = 0.00335 at ω = 1.771, of σr = 0.00308 at ω = 1.954 and so on. Similarly, mode k = 2 changes its stability at the swirl level ω1,2 = 1.766 and reaches a maximum of σr = 0.00765 at ω = 1.826 (figure 3.6a). The σr of this mode continues with oscillations as ω is further increased with a maximum point of σr = 0.0139 at ω = 1.978. The higher-order modes with k = 3 and k = 4 exhibit a similar behaviour (figure 3.6c). The mode k = 3 becomes unstable at swirl level ω1,3 = 1.976 and the mode k = 4 becomes unstable at the swirl level ω1,4 = 2.199. We note that each mode reaches a maximum of σr that increases with the mode number k. This mode becomes the dominant unstable mode in the range of swirl around its maximum of σr.
Comparing with the m = 0 results, the growth rates σ are all complex in the m = 1 case since the equations contain iω terms which reflect the rotation of the azimuthal wave. The imaginary part σi of the modes describes the rate of rotation of the mode around the pipe centreline. For all modes σi is negative, indicating a rotation of the mode in the direction of the base flow rotation. For this case with m = 1 the values of σi as a function of ω are close to −ω, i.e., the rate of rotation of the mode is close to the negative rate of rotation of the base flow (see reference dash-dot blue line σi = −ω in figures 3.6b and 3.6d). When |σi| > ω, the mode rotates around the centreline with an angular speed that is greater than that of the base flow and, when |σi| < ω, the mode rotates around the centreline with an angular speed that is less than that of the base flow. For example, mode k = 1 has a faster rate of rotation than the base flow when ω < 1.613, and has a slower rate of rotation than the base flow when ω > 1.613. Similarly, mode k = 2 has a faster rate of rotation than the base flow when ω < 1.766, and has a slower rate of rotation than the base flow when ω > 1.766. We find that for each of the modes k = 1 through k = 4 with m = 1, the value of σ + iω = 0 at the swirl levels ω1,k where σr = 0.
Table of contents
List of figures
List of tables
1 Introduction
1.1 Vortex stability theory
1.2 The researches conducted in the thesis
2 Vortex stability theory and vortex breakdown phenomenon
2.1 Vortex breakdown phenomenon
2.2 Classical study of swirling flow stability
2.3 The classical vortex stability theory and vortex breakdown
2.4 A ground-breaking stability analysis of axisymmetric swirling flow in a finite-length pipe
3 Three-dimensional stability analysis
3.1 Mathematical model
3.2 Numerical method and validation
3.3 The inviscid-flow perturbation modes on a solid-body rotation flow
3.4 The viscous-flow perturbation modes on a solid-body rotation flow
3.5 The neutral-stability lines of the perturbation modes
3.6 Kinetic energy transfer mechanism
4 The Lamb-Oseen vortex
4.1 The stability equations for Lamb-Oseen vortex
4.2 Numerical method and validation
4.3 The viscous-flow perturbation modes on a Lamb-Oseen flow
4.4 Kinetic energy transfer mechanism of the perturbation modes
5 Conclusion
References
GET THE COMPLETE PROJECT
Stability analysis of swirling flow in a finite-length pipe
