Structure of Periodic Flows through a Channel with a Suddenly Expanded and Contracted Part ()
1. Introduction
Flow paths with abrupt changes in cross-sectional area can be found in a variety of objects such as fluid machinery and chemical plants. Examples of laminar flow paths include injection molds and etched fabrication of printed circuit boards. Such channels are prone to flow separation and self-sustained oscillations. In recent years, transport phenomena in laminar flow paths have attracted attention for applications in medical engineering, such as compact analytical devices and artificial organs. One such application is the use of autonomous oscillations of fluids to promote heat and mass mixing [1] [2] .
In a 2-D channel with only a suddenly expanding or suddenly contracting section, the flow is steady at a Reynolds number on the order of 100 and transitions abruptly from a steady flow to an irregular oscillatory flow [3] - [7] . In contrast, in channels with cavities consisting of both suddenly expanding and suddenly contracting sections, the transition from steady flow to periodic oscillatory flow is known by 2-D numerical analysis and experiments [8] [9] . This oscillatory flow is a Tollmien-Schlichting (T-S) wave excited by the Kelvin-Helmholtz (K-H) instability in the cavity section [10] . The critical Reynolds number for the change from steady to oscillatory flow depends on the boundary conditions as well as the expansion ratio of the flow channel and the aspect ratio of the cavity section. There are two types of boundary condition settings depending on the treatment of the inlet and outlet surfaces: cyclic boundary condition and inlet/outlet boundary condition. In the latter case, the velocity distribution is specified at the inflow surface, and the outflow surface is used as a free boundary condition. Assuming that the cavity is single and symmetrical about the channel center axis, with an expansion ratio Ex = 3 and aspect ratio As = 7/3, the critical Reynolds number Rec = 843 for the inlet/outlet boundary condition [11] [12] , while Rec = 48.7 for the periodic boundary condition [13] . The definitions of the parameters are described below.
Just as linear stability analysis is used for a plane Poiseuille flow, linear stability analysis can be applied to a two-dimensional flow in a suddenly expanding and contracting channel [14] . In a plane Poiseuille flow, the wavenumber component of the disturbance in the main flow direction is included in the exponential part. In contrast, in a flow in a suddenly expanding and contracting channel, the angular frequency and time-decay rate can be obtained, but not the wavenumber directly, because the waveform distribution is in two variables, one in the main flow direction and the other in the vertical direction. Nevertheless, under periodic boundary conditions, the wave number is defined by the channel length. In fact, the time series distribution of velocity at an arbitrary point revealed that the disturbance is a monochromatic wave [15] . However, in the inflow/outflow boundary condition, there are no restrictions on the direction of the main flow, which can lead to multivalence in the wavenumber of the excited disturbance. Although vibration characteristics have been investigated in the downstream entrance region, only the time evolution of velocity at a specific point has been obtained [11] .
In this study, two-dimensional numerical simulations were performed for the flow in a suddenly expanding and suddenly contracting channel under inflow and outflow boundary conditions. By analyzing the structure of the obtained periodic oscillatory flow, the various characteristics of the oscillatory flow, especially the specific structure of the wave number, are clarified.
2. Methods
2.1. Governing Equations
In this study, a two-dimensional incompressible viscous fluid is considered. The governing equations are the continuity Equation (1) and the Navier-Stokes Equation (2), expressed in terms of pressure p and velocity u = (u, v) with time t and space x = (x, y) as parameters.
(1)
(2)
where the density ρ and kinematic viscosity μ/ρ are constants.
Figure 1 shows the channel geometry, related parameters, and coordinate system. The x-axis is in the main flow direction, and the y-axis is perpendicular to the x-axis. The origin O of the coordinate system is set on the center axis of the flow path at the center of the cavity section. The channel geometry was set to an expansion ratio of Ex = 3h/h = 3 and an aspect ratio of As = L0/3h = 7/3, in accordance with previous studies. The lengths of the entrance regions were set to be upstream L1/h = 13 and downstream L2/h = 33, which are sufficiently long so that boundary conditions do not affect the flow.
The Reynolds number is defined as Re = ρUmax h/2μ. Here, the reference velocity Umax is the maximum value of the theoretical solution for a plane Poiseuille flow, and the characteristic length h/2 is half the channel height in the auxiliary section before and after the cavity section. By fixing each value at Umax = 1 [m/s] and h/2 = 1 [m], the physical quantities were made virtually dimensionless, and Re was determined by varying the kinematic viscosity coefficient μ/ρ [m2/s].
2.2. Numerical Solutions
In this study, OpenFOAM was used as the simulation tool [16] , and unsteady numerical simulations were performed. Many numerical simulations using this software have been reported so far [17] [18] [19] . The application to flow in a
![]()
Figure 1. A geometry with parameters and coordinate system.
3-D suddenly expanding duct in a transition zone has been reported [20] [21] [22] .
OpenFOAM discretizes the governing equations using the finite volume method. The numerical scheme used is the Crank–Nicolson method for time, and central difference and Gaussian integration for space. Although the Quadratic Upstream Interpolation for Convective Kinematics (QUICK) method could have been chosen for the convection term, it was found to be unsuitable due to significant numerical oscillations in the calculations. The algorithm was the Pressure-Implicit with Splitting of Operators (PISO) method. The solvers are the Geometric-Algebraic Multi-Grid (GAMG) method for the pressure and Gauss-Seidel for the velocity. Iterations were performed until the residuals reached 10−10.
The boundary conditions are described below. The inflow condition is a theoretical solution for a plane Poiseuille flow with a maximum velocity of Umax = 1 [m/s]. The outflow boundary condition was set to zero for the mean pressure and the Sommerfeld radial boundary condition for the velocity [23] . The velocity boundary condition was compared between the free boundary condition and the Sommerfeld radial boundary condition. The free boundary condition increased the root mean square (RMS) of the velocity component v by up to about 20% compared to the Sommerfeld radial boundary condition. The no-slip boundary condition was specified for the wall. The other boundary conditions were set to zero gradient in the normal direction.
The computational domain is meshed with equispaced-structured mesh elements at a spacing of Δx = Δy. The mesh spacing is expressed as a percentage of the channel height h = 2 [m]. In this study, a mesh with Δx/h = 1/40 was mainly used. The time interval was set to Δt = 0.005 [s]. The validity of these settings is discussed in detail in Section 3.4.
3. Results
3.1. Flow Patterns
The flow in a suddenly expanding and contracting channel transitions from a steady flow to an oscillatory flow above a critical Reynolds number. In this study, the critical Reynolds number was determined to be Rec = 1027. The validity of this value is discussed below in Section 3.5.
Figure 2(a) and Figure 2(b) show examples of streamlines based on the velocity u at a certain instant. This was represented by finding the contour lines of the flow function. At Re = 1030, which is near the critical Reynolds number, there are only a few oscillations that can be read from the stream plots, while at Re = 1200, the meandering of the main flow can be clearly seen.
The streamlines of time-averaged velocity um shown in Figure 2(c) and Figure 2(d) resemble a stable symmetric steady flow, with the mainstream slightly widening as it moves downstream in the cavity section and abruptly shrinking before the steeply contracting section.
![]()
Figure 2. Streamlines (a) and (b) for velocity u, (c) and (d) for time mean velocity um, and (e) and (f) for disturbance velocity u' at an instant. Although cellular vortices are confirmed in (e) and (f), their shapes are distorted as wavenumber of disturbance varies depending on a position. The disturbance shows a tendency to decrease its wavenumber at a center of the channel.
Figure 2(e) and Figure 2(f) show streamlines of the time-varying velocity u'. At Re = 1030, cellular vortices like those of a plane Poiseuille flow were observed in the downstream entrance region, the vortices gradually became distorted, and then vortex pairs appeared as if they were filling the gaps. At Re = 1200, the distortion of the vortices becomes more and more pronounced, and the vortex pattern cannot be called cellular after a while downstream from the steeply decreasing section.
In a plane Poiseuille flow, a cellular vortex pattern appears when streamlines are drawn in a fluctuating velocity field. If the periodic boundary condition is imposed in the main flow direction, the wavenumber is an integer multiple of the inverse of the channel length, and the boundary between the vortices is perpendicular to the wall surface. If periodic boundary conditions are imposed in the main flow direction, even for the flows in a suddenly expanding and contracting channel, the possible wave number is defined by the length of the channel. A non-periodic boundary condition with Dirichlet boundary conditions at the inlet and outlet results in an arbitrariness in wavenumber. Therefore, the wavenumber changes in a complex manner as one moves downstream. A similar case is that in a parallel-plate flow with multiple cavities, two types of waves with different phase velocities are generated at about Re = 58, resulting in an irregular arrangement of cellular vortices [24] .
3.2. Time-Averaged Velocity Distributions
Figure 3(a) plots the time-averaged velocity component um. At x = 7, which is
![]()
Figure 3. Scatter diagrams of time mean velocity um along x-axis at Re = 1200. In the entrance region on the downstream side, a position y where vm is maximum is not constant. In conjunction with Figure 7 described later, the phase of disturbance velocity is influenced by the change in time mean velocity.
the boundary between the cavity and the entrance region, the velocity reaches a minimum at 0 ≤ y ≤ 0.5, which is close to the center axis of the channel, while the velocity increases at 0.6 ≤ y ≤ 0.9. Moreover, a bulge is observed around x = 7, and the closer to y = 1, the larger the rise is. The mean velocity distribution approaches a uniform distribution at 5 ≤ x ≤ 7, and then the distribution becomes parabolic as it proceeds downstream through a steeply decreasing area.
Figure 3(b) and Figure 3(c) show the distribution of time-averaged velocity component vm as in Figure 3(a). The vm has a more regular distribution than um. The mean velocity component is vm = 0 at y = 0 along the channel center axis, reflecting the fact that the disturbances are all symmetric. In the cavity section, the position of the vm maxima is consistent regardless of y. The farther away from the channel center, the larger the maximum value of vm is. The location of the vm maxima in the cavity section is related to the vortex configuration. For example, the maxima are particularly large at x ≈ 5 and x ≈ 6.5, which correspond to the periphery of the most downstream vortex, where the velocity component v is particularly large. In the downstream entrance region, vm reaches a maximum around y = 0.5, though the position is not stable. The smooth curves begin to intersect after x = 6, and gradually decay while intertwining with each other after passing the steeply expanding section.
3.3. Frequency Spectrum
To characterize the oscillatory flow, the spectrum of velocity component v was analyzed by discrete Fourier transform (DFT). To obtain a dimensionless frequency interval of Δf = 0.001 with N = 2000 samples, the data output interval was set to t = 0.5 in dimensionless time, requiring a data for 1000 dimensionless times. In the actual computation, a data for 3000 dimensionless times was set aside to check the degree of convergence, and a spectral distribution was created based on a data for the latter 1000 dimensionless times. The time required was approximately three days with one processor for Δx = 1/40, and approximately one week with four processors for Δx = 1/80. To further refine the frequency of disturbance, the number of samples was set to 100 (N = 1901 - 2000), and DFT were performed on all samples. Therefore, the frequency with the largest amplitude of the fundamental wave was selected.
Figure 4 shows the spectral distribution of the amplitude A with respect to the frequency f obtained by the DFT of the velocity component v at (x, y) = (10, 0.5). Depending on the number of samples, though the distribution may be continuous at the root of the peak, a discrete spectral distribution can be obtained by selecting an appropriate number of samples. It was found that the disturbance includes not only the fundamental wave but also harmonics. The frequency of the fundamental wave tends to increase slightly in proportion to Re. As Re increases, the amplitude of the disturbance also increases.
In the following descriptions, the velocity components of the nth harmonic, which constitute the velocity components u and v, are denoted as Un and Vn.
![]()
Figure 4. Spectrum of amplitude A of velocity components v with respect to frequency f at a point (x, y) = (10, 0.5). Frequency f(V1) and amplitude A(V1) of a fundamental V1 are shown in these figures, which have a discrete spectrum and contains harmonics as well as a fundamental.
The frequency f, amplitude A, and phase θ of the harmonics are defined as dependent variables of Un and Vn.
3.4. Discretization Error Evaluation
Spatial discretization errors are caused by the structure and fineness of the mesh and have a significant impact on the characteristics of the computation results. In this section, we evaluate the spatial discretization error due to mesh size, referring to the work of Ferziger and Perić [25] . If the physical quantity of interest converges monotonically when the mesh size is systematically made finer, then the lattice dependence is sufficiently small. In the case of equally spaced meshes, it is best to reduce the mesh size by half. Although discretization errors are inevitable in simulation results, if the physical quantities of interest tend to converge monotonically, the reproducibility of qualitative features is at least guaranteed. Based on the results obtained from the three models with systematically finer mesh sizes, the scheme order o can be obtained. When the scheme order o is a positive real number, the results converge monotonically to the true value.
Table 1 shows the mesh dependence evaluated by making the mesh size progressively finer at Re = 1200. The meshes are equally spaced with Δx = Δy, and three mesh sizes are available: Δx = 1/20, 1/40, and 1/80. The time step was set to
![]()
Table 1. Comparison of discrete precision of three kinds of mesh size Δx with disturbance velocity v' at Δt = 0.005, (x, y) = (10, 0.5) and Re = 1200. Amplitude of disturbance velocity A(v'), its fundamental harmonic A(V1), second one A(V2), and third one A(V3) are compared.
Δt = 0.005. In addition to the amplitude A(v') of the fluctuation velocity component v', the amplitudes A(V1), A(V2), and A(V3) of the fundamental, second harmonic, and third harmonic waves were selected as the physical quantity of interest, and the data at (x, y) = (10, 0.5) were used. Based on these, the order ο of the scheme was calculated according to Equation (3) [26] [27] :
(3)
Here, r represents the rate of increase of the mesh, r = 2. The φ is an arbitrary physical quantity, and the amplitude A(v') of the fluctuating velocity component or the amplitude A(Vn) of the nth harmonic was assigned to φ.
The order of the fluctuation velocity component v' was 1.34. For each harmonic, the order o of the fundamental V1 was 5.11, which is sufficiently large, while the accuracy dropped sharply to 1.32 for the second harmonic V2 and −2.11 for the third harmonic V3. Because the higher harmonics are more strongly affected by the spatial discretization accuracy, the higher frequency disturbances are at the same time higher wavenumbers.
Table 2 examines the effect of time step Δt. Using the method for evaluating discrete errors in space, the convergence was evaluated for decreasing time increments of Δt = 0.01, 0.005, and 0.0025. The mesh spacing was Δx = 1/40, and the other conditions are the same as in Table 1. In obtaining the order o of the scheme, l in Equation (3) was replaced by the smallest time increment Δt = 0.0025. The amplitude A(v') of the fluctuating velocity component converges to a constant value, and thus can be said to adequately reproduce the characteristics of the target. The amplitude A(V1) of the fundamental wave decreased once and then began to increase as the time step Δt decreased. The order o of the scheme could not be determined, because the change was not monotonic. Although the order o of the amplitude A(V3) of the third harmonic exceeds that of the amplitude A(v') of the fluctuation velocity component, this is considered to be a fake order, considering that the spatial discretization accuracy is not sufficiently ensured as shown in Table 1. It is judged that the variability velocity component after the third harmonic lacks reproducibility qualitatively.
![]()
Table 2. Comparison of discrete precision of three kinds of time increment Δt with disturbance velocity v' at Δx = 1/40, (x, y) = (10, 0.5) and Re = 1200. This is same as Table 1 about the compared objects.
3.5. Bifurcation Diagram
Figure 5 is a bifurcation diagram based on the DFT result of the velocity component v at position (x, y) = (10, 0.5). The mesh spacing was set to Δx = 1/40. The initial value is of an oscillatory flow at Re larger than the target Re, and the computation is continued until the oscillations are sufficiently damped. The condition for convergence was that the peak value of the velocity component v be within 1% of the amplitude in the range of 1000 dimensionless times.
The critical Reynolds number was obtained by drawing an approximate curve using the seventh-order least-squares method based on the plotted points of the amplitude A(V1) of the fundamental wave, and its value was found to be Rec = 1027. In the article by Mizushima et al. [12] , the critical Reynolds number is given as Rec = 843. Comparing the two numerical solution methods, the previous study discretized the time term with first-order accuracy and the other terms with second-order accuracy, while the present study adopted a second-order accuracy scheme for all terms. Although the time interval of the present study is Δt = 0.005 compared to Δt = 0.0005 in the previous study, the steady flow condition was maintained even when Δt = 0.0005 was used at Re = 1000, which refers that the order of the scheme was dominant. Furthermore, as shown in Table 2, the time increments have little influence on the results. Therefore, the dependence on the value of the critical Reynolds number is not due to a lack of discretization accuracy, but rather to the characteristics of the numerical solution method.
The approximate curves for the second and third harmonics are drawn by means of a type of cubic spline interpolation, the Akima spline interpolation, assuming that the critical Reynolds number obtained from the fundamental is common to the other harmonics. From this approximate curve, it can be read that the fraction of harmonics in the disturbance increases with the Reynolds number. From the viewpoint of bifurcation theory, when a stationary flow such as a plane Poiseuille flow becomes unstable and transitions to an oscillatory flow, it is the Hop bifurcation, in which a branch of the oscillatory flow vertically branches off from a branch of the stationary flow at a critical point. A similar bifurcation structure should exist for suddenly expanding and suddenly contracting channel flows. In fact, the bifurcation diagram based on the amplitude A(V1) of the fundamental suggests that it is the Hop bifurcation. However, for the harmonics, the bifurcation is an incomplete Hop bifurcation, with a branch
![]()
Figure 5. A bifurcation diagram based on amplitude of the fundamental harmonic A(V1), the second one A(V2) and the third one A(V3) of velocity components v with respect to Reynolds numbers Re. The critical Reynolds number is of Rec = 1027. An approximate curve of the fundamental wave was drawn by least-squares method with seven parameters. In proportion to Reynolds numbers, not only the amplitude of velocity components v but also a proportion of harmonics in disturbance increases. When a steady flow becomes unstable and transitions to an oscillatory flow, it is a Hopf bifurcation, but with respect to the harmonics it becomes an incomplete Hopf bifurcation.
of the oscillatory flow extending tangential to the branch of the stationary flow. If this result is taken as is, the cascade down of disturbances is activated in proportion to the Reynolds number, with a very little activity near the critical Reynolds number. Nevertheless, high-frequency disturbances have high wavenumber at the same time, and there is a concern about the effects of numerical errors due to mesh coarseness.
3.6. Amplitude Distributions
Figure 6 shows the amplitude distribution of each harmonic with respect to the velocity components u and v. By applying the DFT to the velocity distribution over the entire flow field, the characteristics of the disturbance were evaluated.
The amplitudes A(Un) and A(Vn) of the harmonics are even or odd functions with respect to the channel center axis at the symmetry of the oscillatory flow. For the fundamental wave, A(U1) is an odd function, and A(V1) is an even function. Such disturbances with parity are called antisymmetric disturbances and are expected to be initially unstable with respect to the flow field. These properties are consistent with those of a plane Poiseuille flow. For the nth harmonic, A(Un) is even-functional for even-order harmonics, and A(Vn) is even-functional for odd-order harmonics.
To evaluate the contribution of each harmonic to the overall disturbance, the ratio of the maximum amplitude of each harmonic in the flow field to the representative velocity was compared. At Re = 1030, the maximum values of the amplitude distribution of the fundamental wave are A(U1) = 0.0377 and A(V1) = 0.0415, while those of the second harmonic are A(U2) = 0.00208 and A(V2) =
![]()
Figure 6. Distributions of amplitude for some frequency components of disturbance. There are points where the amplitude temporarily becomes very weak at a downstream side of the channel, which are like a clause generated by superposition of waves. The disturbance is a synthesized wave which consists of waves having different wavenumbers.
0.00152, less than 10% of those of the fundamental wave. At Re = 1200, the maximum values of the variation velocity component u' are 0.217, 0.071, 0.030, 0.019, 0.011, and 0.005 from the fundamental to the sixth harmonic, in that order, while those of the fluctuation velocity component v' are 0.197, 0.042, 0.027, 0.013, 0.009, and 0.006. The fraction of fundamental waves in the disturbance is high. Furthermore, the amplitude distribution of the fundamental is qualitatively in good agreement with the RMS distribution of the fluctuation velocity components u' and v'. Therefore, the behavior of the disturbance is represented by the fundamental waves.
The maximum amplitude for the fluctuation component of u' is located near x ≈ 10, slightly downstream of the suddenly contracting part. The amplitude peaks extend from the corners of the suddenly contracting section. This amplitude distribution is assumed to be caused by the compression of the vibration in the y-direction during the development of a boundary layer from the wall surface due to viscous action after the vibration flow excited in the cavity section enters the downstream parallel plate space. In contrast, the amplitude of v' peaks at x = 7, when the flow path is suddenly contracting. This is thought to be because the vibration in the y-direction is suppressed between the downstream parallel plates, resulting in a gradual damping.
3.7. Phase Distributions
As the wavenumber k is the reciprocal of the wavelength λ, the wavenumber k can be obtained from the relation Δθ/Δx = 2π/λ = k based on the phase θ shift between the two points. The phase is obtained by the DFT. Based on the phase distribution, the wavenumber distribution of each harmonic can be estimated by reading the rate of change of the phase Δθ/Δx.
Figure 7 shows the distribution of phase θ obtained by applying the DFT to the velocity components u and v as in Figure 6. However, if the phase θ is expressed as it is in a contour plot, a discontinuity surface is created between π and −π, making the plot difficult to understand. Therefore, the phase is treated as sin θ using a trigonometric function so that continuous contours can be drawn.
As the period of rotation of the vortex in the cavity corresponds only to the period of the fundamental wave, the second and subsequent harmonics are judged to have been generated by cascading down. Moreover, as the frequency increases, the spacing of the phase distribution becomes narrower, while the disturbance of the high-frequency component is simultaneously higher in wavenumber. This is supported by the process of cascading down the disturbances defined by the trigonometric functions of the convective term in the Navier-Stokes equations [28] .
For the velocity component u, the phase velocity of the disturbance increases in the center of the channel, whereas the phase velocity is almost constant and does not change near the wall. Therefore, a boundary parallel to the channel can be observed in the phase distribution, which is clearly visible at Re = 1200 and
![]()
Figure 7. Distributions of phase for some frequency components of disturbance. Since a fault line is confirmed in the phase distributions of the fundamental wave, the fundamental wave is composed of at least two modes: a maxima of an amplitude distribution along the y-axis is near a wall and at a center of the channel.
even at Re = 1030, where the oscillations are very small. The fundamental wave has a clear boundary at y ≈ 0.7, and multiple boundaries can be seen at higher harmonics.
For the velocity component v, the phase distribution becomes non-uniform with respect to the Reynolds number. At Re = 1030, the contour lines of the phase distribution are approximately perpendicular to the channel wall. At Re = 1200, the contours of the phase distribution run obliquely, and their slopes vary irregularly. The phase of the disturbance changes randomly as the flow progresses. The phase distribution is almost constant near the wall, while it changes more irregularly as the flow approaches the center of the channel. The fact that the phase distribution varies from place to place clearly indicates that the wavenumber is not constant.
4. Discussion
Based on the above results, the wavenumber of disturbance in the oscillatory flow is discussed. As the channel under consideration here is not under a periodic boundary condition, the wavenumber can take any value. In fact, the phase distribution in Figure 7 shows that the wavenumber is not uniform; varies from place to place. This chapter discusses the mechanism of this phenomenon.
The fact that the phase varies from place to place suggests that the disturbance in the oscillatory flow is not a monochromatic wave. Furthermore, as shown in Figure 6, the amplitude distribution of the fundamental wave of the velocity component v has a point where the amplitude temporarily becomes very small at the center of the downstream channel. This is like a node generated by the superposition of stationary waves. Since a single wavenumber does not produce a node, the disturbance is a composite wave composed of multiple wavenumbers. Even in the cavity section, the distribution of wavenumbers is not uniform and includes several combinations of wavenumbers that are not integer multiples at this point. As shown in the stream diagram in Figure 2, the cavity generates multiple vortex pairs that are responsible for the generation of disturbances. The effect of their different sizes is thought to generate multiple disturbances with different wave numbers.
5. Conclusions
In this study, two-dimensional numerical simulations of the flow in a suddenly expanding and suddenly contracting channel were performed to reproduce periodic flow in the range of 1020 < Re ≤ 1200. The following conclusions were obtained regarding the oscillatory characteristics of the periodic flow.
The frequency of the periodic flow is a discrete spectrum. The fundamental wave accounts for most of it, though harmonics are also included. The proportion of harmonics in the frequency spectrum increases with Reynolds number. In the cavity section, only the fundamental oscillations are excited by the Kelvin-Helmholtz instability, while the harmonics are generated by cascading down.
Disturbances consist of multiple wavenumbers with arbitrary frequency components. The amplitude and phase of the disturbance are extracted by Fourier transform in the form of a composite wave at a specific frequency. The amplitude distribution has nodes due to the superposition of waves with different wavenumbers. The phase distribution is non-uniform and tends to be lower in the center of the channel. Therefore, the streamlines of the velocity field of the disturbance are distorted so that the cellular vortex protrudes at the center of the channel.
As the phase distribution of the fundamental wave shows boundaries, the disturbance of the fundamental wave consists of at least two types: one with a small attenuation rate and the maximum value of the y-directional velocity distribution near the wall, and the other with the maximum value at the center of the flow path. These correspond to the wall mode and the center mode in the linear stability analysis of a plane Poiseuille flow.
The magnitude of the wavenumber is linked to the change in the time-averaged velocity. The time-averaged velocity component um promotes low wavenumber in the center of the channel as it develops into a parabolic distribution in the entrance region. The time-averaged velocity component vm does not have a fixed maximum position in the entrance region, which makes the phase distribution of the time-varying velocity component v' irregular.
In summary, the periodic flow excited in the cavity section consists of multiple disturbances with different wavenumbers and velocity distributions for all frequency components, and the wavenumbers change irregularly in the downstream region due to differences in attenuation rates and changes in mean velocity.