2. Numerical method

2.1. Problem statement

The computational domains are presented in figure 1. A fluid column of stroke length $L(t)$ is pushed through a two-dimensional nozzle with a height $H$ at the exit plane. The origin of coordinates is located at the centre of the nozzle exit, and the symmetry axis is aligned with the $x$ direction. The distances from the nozzle exit to the upstream nozzle inlet and downstream exit boundary are both fixed at $40H$ . The outer boundary is placed at a lateral distance of $6H$ from the symmetry axis. A squared-off lip at the orifice is used to avoid numerical singularities.

Figure 1. Computational domains (not to scale) for (a) straight nozzle and (b) orifice nozzle. Sample meshes (every sixth mesh point is shown in each direction) for (c) straight nozzle and (d) orifice nozzle.

No-slip condition is enforced on the nozzle wall and the outer orifice plate. The pressure outlet with zero gauge pressure is employed on the outer and downstream exits. The piston motion is modelled as a time-dependent velocity inlet boundary (Rosenfeld, Rambod & Gharib Reference Rosenfeld, Rambod and Gharib1998; Zhu et al. Reference Zhu, Zhang, Gao and Yu2023a,b ). A uniform velocity profile is specified at the nozzle inlet to control piston velocity $U_p(t)$ . An impulsive velocity programme at the nozzle exit is defined as

(2.1) \begin{equation} U_{0}(t) \begin{cases} U & \mbox {for }\quad t \leqslant \tau , \\[5pt] 0 & \mbox {for }\quad t \gt \tau , \\ \end{cases} \end{equation}

where $U$ is the constant velocity at the nozzle exit and $\tau$ is the jet discharge time. The maximum stroke length and time-averaged exit velocity are $L_m=\int _{0}^{\tau } U_{0} \textrm{d}t=U\tau$ and $\overline {U}_{0}=(\tau )^{-1}\int _{0}^{\tau } U_{0} \textrm{d}t=U$ . The non-dimensional formation time is defined as $t^*=\overline {U}_{0}t/H$ .

As listed in table 1, a total of 15 cases are simulated. Water is used as the working fluid, with density $\rho =998.2\,\mathrm{kg}\,\mathrm{m}^{-3}$ and dynamic viscosity $\mu =1.003\times 10^{-3}\,\mathrm{kg}\,\mathrm{m}^{-1},\mathrm{s}^{-1}$ . The same exit height $H=2\,\mathrm {cm}$ and jet velocity $U=5\,\mathrm{cm}\,\mathrm{s}^{-1}$ for both the straight nozzle and orifice nozzle provide the same momentum flux at the nozzle exit. To avoid the possible occurrence of Kelvin–Helmholtz instability during jet discharge, a Reynolds number ( $Re_H=\overline {U}_{0}H/\nu$ ) of 1000 is chosen for all cases, where $\nu$ is the kinematic viscosity. The onset of vortex pinch-off for the non-parallel flow cases would occur earlier than for the parallel flow cases (Limbourg & Nedić Reference Limbourg and Nedić2021c ). The chosen stroke ratios for the orifice cases should be large enough to study the formation number in the present study. Therefore, the chosen stroke ratios for orifice cases are slightly smaller than those for cases with the straight nozzle.

Table 1. Cases considered.

2.3. Verification of numerical method

Cases for respective independent tests are listed in table 2. The tests are evaluated for the streamwise vortex trajectory ( $x^*_v=x_v/H$ ) for cases SN250 and ON175. For the domain-independent tests shown in figure 2(a), the convergence of numerical results can be observed as the domain size increases. The trajectory for the chosen domains deviates less than 2 % from that for the larger domains $40H\times 8H$ . For the grid-independent tests in figure 2(b) and time-step-independent tests in figure 2(c), the numerical results are not sensitive to the mesh and time-step choices. Comparing with the denser mesh $6001\times 901$ or smaller time step $0.0005\,\mathrm {s}$ , the errors for the developed vortex pair are generally less than 4 %. In short, the error and convergence for the present domain size, mesh and time step suggest that the numerical set-up is sufficient to describe the vortex motion adequately. To validate the numerical results, the evolution of the transverse trajectory ( $y^*_v=y_v/H$ ) starting at the nozzle edge is plotted in figure 2(d). The numerical result for the straight nozzle satisfies the $2/3$ power law (red dashed line) suggested by the roll-up theory of a two-dimensional vortex sheet developed at the edge of a semi-infinite flat plate (Saffman Reference Saffman1978).

Table 2. Grid sizes used for independent tests.

Figure 2. Numerical tests on vortex trajectories for (a) domain sizes, (b) grid sizes and (c) time steps. (d) Validation with Saffman theory. Abbreviations: SDS/LDS, CM/FM and STS/LTS represent small/large domain size, coarser/finer mesh and small/large time step.

2.4. Vortex invariants and vortex boundary

The main dynamical properties of the two-dimensional vortex motion are the circulation, the hydrodynamic impulse and the kinetic energy. The quantities can be obtained by the area integration inside the vortex core on the upper semi-plane (Batchelor Reference Batchelor1967; Saffman Reference Saffman1992), such as

(2.2) \begin{equation} \varGamma _v=\int \int \omega \,\textrm{d} x\,\textrm{d} y , \end{equation}

(2.3) \begin{equation} I_v=\rho \int \int y\omega \,\textrm{d}x\,\textrm{d}y, \end{equation}

(2.4) \begin{equation} E_v=\frac {1}{2}\rho \int \int \omega \psi \,\textrm{d}x\,\textrm{d}y, \end{equation}

where $\omega$ is the azimuthal vorticity and $\psi$ is the stream function. It is noted that the integrals of the vortex pair are different from those of the vortex ring (Gharib et al. Reference Gharib, Rambod and Shariff1998). The cutoff level to determine the vortex core boundary is defined as 5 % of the local maximum vorticity (Rosenfeld et al. Reference Rosenfeld, Rambod and Gharib1998; Zhao et al. Reference Zhao, Frankel and Mongeau2000; Sau & Mahesh Reference Sau and Mahesh2007). Sensitivity analysis had been made by the cutoff level of 2 %. The differences of the vortex invariants are within 3 % for both nozzle configurations. Therefore, the vortex invariants are found to be insensitive to the choice of the cutoff level. In addition, the pressure-based method is used to further identify the rear vortex boundary (Lawson & Dawson Reference Lawson and Dawson2013). The features of the vortex boundary in the pressure field are found to be accurate by comparing with the Lagrangian vortex boundary (Baskaran & Mulleners Reference Baskaran and Mulleners2022). Based on the maximum in the centreline pressure, the trailing pressure maximum at the vortex rear boundary is determined along the symmetry axis (e.g. the blue dashed lines in figure 3 a, b). Moreover, the vortex centre (red dashed lines) and front boundary (yellow dashed lines) can also be identified by this pressure-based method. This method can locate the rear boundary before a clear vorticity separation appears between the leading vortex and the trailing jet, as shown in figure 3(c). Based on the cutoff of vorticity level together with the pressure-based method, the vortex boundary can be identified over a wide range of non-dimensional time for the calculation of vortex invariants throughout the entire studies.

Figure 3. Pressure-based vortex boundary (case SN250). (a,b) The pressure normalized by the maximum and (c,d) the vorticity contours normalized by local maximum (from 0.05 to 1, with 10 levels). Results at (a,c) $t^*=20$ and (b,d) $t^*=50$ . Blue, red and yellow dashed lines represent the rear, the centre and the front of the vortex.

3. Slug model for total invariants

The invariants of the motion were introduced in unbounded inviscid incompressible flows, such as circulation, impulse and energy. The classical slug model was originally developed to predict total circulation shedding from the sharp edge of an axisymmetric straight nozzle (Didden Reference Didden1979). By definition, the rate of circulation change is equal to the vorticity flux at the nozzle exit plane during jet discharge, i.e.

(3.1) \begin{equation} \frac {\textrm{d}\varGamma _{slug}}{\textrm{d}t}=\int _{0}^{\infty }u_{x}\omega \,\textrm{d}y =\int _{0}^{\infty }u_{x} \left (\frac {\partial u_{y}}{\partial x}-\frac {\partial u_{x}}{\partial y} \right ) \textrm{d}y=\int _{0}^{\infty }u_{x}\frac {\partial u_{y}}{\partial x} \textrm{d}y + \frac {1}{2}U_{outlet}^2, \end{equation}

where $U_{outlet}$ is the centreline velocity at the nozzle exit. Considering the uniform assumption of streamwise velocity and the parallel-flow assumption for a straight nozzle, the transverse velocity component $u_{y}$ would be zero, leading to an approximated form:

(3.2) \begin{equation} \textrm{d}\varGamma _{slug}=\frac {1}{2}U^2\textrm{d}t. \end{equation}

This slug model for the total circulation was extended to estimate the hydrodynamic impulse and kinetic energy for a parallel starting jet (Gharib et al. Reference Gharib, Rambod and Shariff1998), such as

(3.3) \begin{equation} \textrm{d}I_{slug}=\frac {1}{4} \pi \rho U^2 {D}^2 \textrm{d}t, \end{equation}

(3.4) \begin{equation} \textrm{d}E_{slug}=\frac {1}{8} \pi \rho U^3 {D}^2 \textrm{d}t, \end{equation}

where $D$ is the diameter of a circular nozzle. However, this model deviates from the real situation with a time-dependent streamwise velocity profile due to the growth of the boundary layer inside the nozzle (Rosenfeld et al. Reference Rosenfeld, Rambod and Gharib1998). Thus, this model performs worse at a later time for parallel starting jets. In addition, this model was less successful in the cases of non-parallel flows generated by a circular orifice (Krieg & Mohseni Reference Krieg and Mohseni2013). The reason was attributed to the significant transverse velocity associated with the over-pressure effect at the nozzle exit (Krueger Reference Krueger2005).

To overcome the disadvantage of the slug model in orifice-generated flows, a modification for axisymmetric flows was proposed based on the contraction coefficient defined at the vena contracta originally for two-dimensional jets (Limbourg & Nedić Reference Limbourg and Nedić2021a ). The contraction coefficient is defined as the ratio of the area of the vena contracta $A_{\star }$ to the area of the nozzle exit $A$ , such as

(3.5) \begin{equation} {C_c} = \frac {A_{\star }}{A} = \frac {{D}_{\star }^2}{D^2} = \frac {U}{U_{\star }} = \frac {L}{L_{\star }}. \end{equation}

According to mass conservation, the quantities at the nozzle exit can be transformed into those at the downstream vena contracta (subscript $\star$ ). The contraction-based slug model for impulsively axisymmetric starting jets can thus be expressed as

(3.6) \begin{equation} \varGamma _{\star } = \frac {1}{2}L_{\star }U_{\star } = \frac {1}{2}LU \times 1/C^2_c, \end{equation}

(3.7) \begin{equation} I_{\star } = \frac {1}{4} \pi \rho L_{\star }U_{\star }{D}_{\star }^2 = \frac {1}{4} \pi \rho LU{D}^2 \times 1/C_c, \end{equation}

(3.8) \begin{equation} E_{\star } = \frac {1}{8} \pi \rho L_{\star }U_{\star }^2{D}_{\star }^2 = \frac {1}{8} \pi \rho LU^2{D}^2 \times 1/C^2_c. \end{equation}

Determination of the contraction coefficient $C_c$ for orifice cases was related to the orifice-to-tube diameter ratio (Limbourg & Nedić Reference Limbourg and Nedić2021a ). It is noted that this contraction-based model with suggested value $C_c=0.9$ was available for parallel flows with considerable boundary layer (Limbourg & Nedić Reference Limbourg and Nedić2021b ). Thus, this contraction-based model may be used to predict total invariants for both an axisymmetric straight nozzle and orifice nozzle.

For impulsively two-dimensional flows (Afanasyev Reference Afanasyev2006), the slug-flow model can be written as

(3.9) \begin{equation} \varGamma _{slug2} = \frac {1}{2}LU, \end{equation}

(3.10) \begin{equation} I_{slug2} = \rho LHU, \end{equation}

(3.11) \begin{equation} E_{slug2} = \frac {1}{2}\rho LHU^2. \end{equation}

With the consideration of contracted flows generated by the two-dimensional orifice nozzle, the two-dimensional slug model can be rewritten as

(3.12) \begin{equation} \varGamma _{\star 2} = \frac {1}{2}L_{\star }U_{\star }=\frac {1}{2}LU\times 1/C^2_c, \end{equation}

(3.13) \begin{equation} I_{\star 2} = \rho L_{\star }H_{\star }U_{\star } = \rho LHU\times 1/C_c, \end{equation}

(3.14) \begin{equation} E_{\star 2} = \frac {1}{2}\rho L_{\star }H_{\star }U_{\star }^2 = \frac {1}{2}\rho LHU^2\times 1/C^2_c, \end{equation}

where the contraction coefficient follows the relations

(3.15) \begin{equation} {C_c} = \frac {H_{\star }}{H} = \frac {U}{U_{\star }} = \frac {L}{L_{\star }}. \end{equation}

The slug model (3.9)–(3.11) and the contraction-based slug model (3.12)–(3.14) for two-dimensional flows are evaluated in § 5.

4. Two-dimensional formation number

According to the definition of critical stroke ratio for the pinch-off process (Gharib et al. Reference Gharib, Rambod and Shariff1998), the existence of a trailing jet after jet discharge is crucial to describe the maximum energy state of the leading vortex. To examine the effect of stroke ratio on the existence of a trailing jet, comparisons of the vorticity fields at the instant slightly larger than the maximum stroke ratio ( $t^*=L_m/H+2.5$ ) and at a much later moment ( $t^*=L_m/H+12.5$ ) are presented in figures 4 and 5 for the straight nozzle and orifice nozzle, respectively. Vorticity fields are normalized by the local maximum. As shown in figures 4(a,c,e) and 5(a,c,e), the starting jet includes a leading vortex ring and a trailing jet without any separations between them at around the maximum stroke ratio. The trailing jet is eventually absent when the stroke ratio is small (see figures 4 b and 5 b). On the other hand, shear layers separating from the leading vortex pair can be observed in figures 4(d) and 5(d) when the stroke ratio is large. This result implies a critical formation number at around 10 for the two-dimensional starting jets generated by either the straight nozzle or orifice nozzle. With the increase in stroke ratio, the vortex pair grows in size and becomes more difficult to disconnect completely from the trailing jet in the vorticity field (figures 4 f and 5 f).

Figure 4. Normalized vorticity fields for the straight nozzle at stroke ratios of (a,b) 7.5, (c,d) 10 and (e,f) 15. The vorticity field has 10 levels from 0.05 to 1.

Figure 5. Normalized vorticity fields for the orifice nozzle at stroke ratio of (a,b) 7.5, (c,d) 10 and (e,f) 15. The vorticity field has 10 levels from 0.05 to 1.

The absence of a clear pinch-off in a previous simulation (Pedrizzetti Reference Pedrizzetti2010) might be due to the limited stroke ratio (up to 10). As stated in experimental results (Gharib et al. Reference Gharib, Rambod and Shariff1998), a complete pinch-off might require non-dimensional time up to two times the formation number for an axisymmetric jet. A clear pinch-off for vortex pairs before the termination of jet discharge was also not observed in the present studies. Therefore, cases with larger stroke ratio could be meaningful for better understanding of the whole pinch-off process in starting jets. Additionally, a clear pinch-off for vortex dipoles was also not observed by the visualization in the experiment for stroke ratio up to 51 (Afanasyev Reference Afanasyev2006). This might be due to the low Reynolds numbers (of the order of hundreds). It is also noted that a clear pinch-off also could not be found in axisymmetric jets at low Reynolds numbers even though the vortex ring has attained the maximum circulation level (Bi & Zhu Reference Bi and Zhu2020).

The appearance of a trailing jet in the present simulations is different from the initial absence of a trailing jet in the experiments of non-parallel planar starting jets generated by a rectangular nozzle with an aspect ratio of 40 (Steinfurth & Weiss Reference Steinfurth and Weiss2020). In their experiments, the leading vortex ring could absorb all vorticity from the nozzle exit up to a formation time of 12. Vortex pinch-off was not observed. This may be due to the initial absence of a trailing jet. The result was explained by the effects of high normalized over-pressure (greater than 1) at the nozzle exit. The present orifice cases are with low normalized over-pressure $p^*_{over} = (p_{exit}-p_{\infty }) /\rho U^2 \approx 0.5$ , where $p_{exit}$ and $p_{\infty }$ are the pressure at the nozzle exit centre and the pressure at the ambient fluid. It is noted that the absence of a trailing jet was also not observed in experiments of non-parallel axisymmetric jets with low over-pressure of about 0.4 by Krieg & Mohseni (Reference Krieg and Mohseni2013). Therefore, the appearance of a trailing jet and absence of significant vortex growth in the present studies might be due to the low magnitude of over-pressure.

Figure 6. Development of normalized vortex circulation (a,c) and maximum vortex circulation versus stroke ratio (b,d) for straight nozzle (a,b) and orifice nozzle (c,d).

To quantify the growth of the vortex pair, the evolution of vortex circulation is presented in figure 6. Before the trailing pressure maximum appears, the vortex circulation can be approximated by the total circulation. The vortex circulation is normalized as $\varGamma ^*_{v} = \varGamma _{v}/UH$ . As shown in figure 6(a,c), the leading vortex pair continuously absorbs vorticity from the trailing jet even after the termination of piston motion. Due to the self-induced velocity, the vortex pair eventually disconnects from the trailing jet and starts its decaying stage. This process also can be found in the formation of the leading vortex ring (Rosenfeld et al. Reference Rosenfeld, Rambod and Gharib1998). For the orifice cases, the shear layer is more unstable than that for the straight nozzle. Although the Kelvin–Helmholtz instability is not observed during the jet discharge, secondary vortices in the trailing jet appear for cases ON150 and ON175 due to the perturbation caused by the termination of jet discharge (see figure 5 f). Thus, the vortex circulation can increase rapidly by vortex pairing.

Figure 7. Development of normalized vortex energy (a,c) and vortex energy at $t^*=60$ versus stroke ratio (b,d) for straight nozzle (a,b) and orifice nozzle (c,d).

The maximum vortex circulation for each case is plotted in figure 6(b,d). For straight nozzle cases, the growth rate of vortex circulation becomes smaller with the increase of stroke ratio. This convergent tendency suggests a limiting situation in vortex pair formation, which is analogous to tube-generated vortex ring formation (Gharib et al. Reference Gharib, Rambod and Shariff1998). For the orifice cases, the increase of slope can be explained by the vortex pairing with the trailing vortex. After the pairing, the leading vortex moves faster due to the increased circulation strength. This process leads to a significant decrease in the growth of vortex circulation after $L_m/H=15$ . This result therefore agrees with the trend of the limiting situation in vortex formation analogous to orifice-generated vortex ring formation (Limbourg & Nedić Reference Limbourg and Nedić2021c ).

To further examine the state of the vortex pair, the evolution of vortex energy is presented in figure 7. The vortex energy is normalized as $E^{**}_v=E_v/\rho \varGamma ^2_v$ . Large and small values of the normalized energy are represented by thin and thick vortex pairs (Pierrehumbert Reference Pierrehumbert1980). As shown in figure 7(a,c), the vortex core develops into a thick core with decreasing non-dimensional energy. The thickest state of a vortex pair is generally enhanced and postponed by the increasing stroke ratio. The vortex pairing process of large-stroke-ratio orifice cases (cases ON150 and ON175) greatly thickens the leading vortex core within a short time. The trailing jet for case ON175 is more energetic on account of the longer jet discharging time, leading to an earlier vortex pairing. During the post-formation stage where the leading vortex becomes stable and decays steadily after a complete pinch-off, vortex pairs generated by non-parallel jets settle to a constant state of energy faster than those of parallel jets. The convergence of vortex energy with respect to a nearly constant vortex impulse after a complete pinch-off from the trailing jet satisfies the suggestion of Kelvin’s variational principle. For the energy at $t^*=60$ , the vortex energy for straight nozzle cases converges to a universal value of about 0.068 (figure 7 b). For the orifice cases, the universal energy is at about 0.076 (figure 7 d). The universal energy of steady vortex rings was also observed in experiments (Gharib et al. Reference Gharib, Rambod and Shariff1998). A thicker straight-nozzle-generated vortex core (i.e. smaller non-dimensional energy) than the orifice-generated vortex core also could be found in starting axisymmetric jets (Rosenfeld, Katija & Dabiri Reference Rosenfeld, Katija and Dabiri2009; Krieg & Mohseni Reference Krieg and Mohseni2021). This might be explained by the thicker shear layer of straight nozzle cases and slower vortex motion due to the absence of the over-pressure effect.

On account of the aforementioned similarities between vortex ring formation and vortex pair formation, the formation number (i.e. the onset of pinch-off) will be further studied in two-dimensional flows. Five criteria are presented here, including three local analyses (total circulation, vortex velocity and induced velocity criteria) and two global analyses (entrainment and circulation ratio criteria). The formation numbers will be identified by the circulation criterion first and subsequently they will be further analysed by the other four criteria.

4.1. Circulation criterion

The formation number $F$ may be determined by the circulation criterion (Gharib et al. Reference Gharib, Rambod and Shariff1998) when the total circulation $\varGamma ^*_t=\varGamma _t/UH$ becomes equal to the maximum vortex circulation (Rosenfeld et al. Reference Rosenfeld, Rambod and Gharib1998). In terms of locally maximum vortex circulation $\varGamma ^*_{vmax}$ for each case, the difference between cases SN225 and SN250 is 5 % while that for cases ON150 and ON175 is about 4 %. The differences are small enough to approximate maximum vortex circulation with errors of about 5 % by averaging these two cases (i.e. the two highest values in vortex circulation) for both nozzles. As shown in figure 8, the formation numbers (red arrows) are found by the intersection of total circulation (black solid line) and maximum vortex circulation (blue dashed line). The formation numbers are about 13.6 and 9.3 for the straight nozzle and orifice nozzle, respectively. These are much larger than the 4 and 2 for axisymmetric jets (Gharib et al. Reference Gharib, Rambod and Shariff1998; Limbourg & Nedić Reference Limbourg and Nedić2021c ).

For the non-parallel jet, the formation number of 7 found in Gao & Yu (Reference Gao and Yu2016a ) is smaller than that found in the present study. The smaller value may be attributed to the trailing jet instability at higher Reynolds number (Zhao et al. Reference Zhao, Frankel and Mongeau2000; Gao & Yu Reference Gao and Yu2012).

Figure 8. Circulation criterion to determine formation number for (a) straight nozzle (total circulation obtained from case SN250 and vortex circulation obtained from the average of cases SN225 and SN250) and (b) orifice nozzle (total circulation obtained from case ON175 and vortex circulation obtained from the average of cases ON150 and ON175).

4.2. Entrainment criterion

To examine the entrainment property of starting jets, the area of scalar-carrying fluid is shown in figure 9. The area $A_{sc}^*=A_{sc}/H^2$ is computed inside the jet boundary which is defined above a threshold of 1 % for scalar concentration (Sau & Mahesh Reference Sau and Mahesh2007). The jet area increases linearly for the impulsive piston motion. After the termination of piston motion, the slopes become nearly constant (figure 9 a,c). The increase can be attributed to diffusion and mixing of ambient fluids. A constant slope for the rate of area change for $L_m/H \leqslant F$ and a decreased slope for $L_m/H \gt F$ (Sau & Mahesh Reference Sau and Mahesh2007) are found. The different growth rates of $\textrm{d}A^*_{sc}/\textrm{d}t^*$ in terms of stroke ratio can be explained by the fact of the different mixing abilities between the leading vortex and the trailing jet. More and more vorticity would be accumulated in the trailing jet with the increase of stroke ratio for $L_m/H \gt F$ . As a consequence, the relative contribution of the leading vortex to the total entrainment decreases as the stroke ratio increases. As shown in figure 9(b,d), the growth rates of $\textrm{d}A^*_{sc}/\textrm{d}t^*$ decrease significantly after the critical formation time. This situation in two-dimensional flows is similar to the situation in axisymmetric flows (Sau & Mahesh Reference Sau and Mahesh2007). Based on the intersection of the fitted lines (blue dashed lines), the formation numbers (13.4 for the straight nozzle and 9.6 for the orifice nozzle) can be determined with a prior knowledge of the formation number. The differences between the circulation criterion and the entrainment criterion are about 1 % and 3 % for the straight nozzle and the orifice nozzle, respectively.

Figure 9. Development of normalized area of scalar-carrying fluid (a,c) and normalized growth rate of the area versus stroke ratio (b,d) for straight nozzle (a,b) and orifice nozzle (c,d).

4.3. Kinematic criterion

A kinematic criterion in terms of the translational velocity was also proposed in previous studies for axisymmetric jets (Mohseni & Gharib Reference Mohseni and Gharib1998; Shusser & Gharib Reference Shusser and Gharib2000). Pinch-off starts when the vortex ring velocity equals a specific velocity (e.g. the jet velocity near the leading vortex). A specific value $N$ is defined as the ratio of this particular velocity to the space-averaged velocity at the nozzle exit. This specific value $N$ in parallel jets is 0.5 from a Hamiltonian theory for vortex rings (Mohseni & Gharib Reference Mohseni and Gharib1998) or about 0.6 from the analytical model with a closure from experimental results (Shusser & Gharib Reference Shusser and Gharib2000). By the fourth-order central-difference approximation, the vortex velocity can be determined from the vortex trajectory. Moreover, the normalized vortex translational velocity is defined as $U^*_{v}=U_{v}/U$ . Based on the formation numbers found from the circulation criterion, the value $N$ can be identified to be about 0.49 (figure 10 a) and 0.58 (figure 10 b) for the straight nozzle and orifice nozzle, respectively. Thus, our results for vortex pair formation can generally agree with the suggestions (Mohseni & Gharib Reference Mohseni and Gharib1998; Shusser & Gharib Reference Shusser and Gharib2000) and measurements (Limbourg & Nedić Reference Limbourg and Nedić2021d ) for axisymmetric flows.

Figure 10. Development of vortex translational velocity for (a) straight nozzle (case SN250) and (b) orifice nozzle (case ON175).

Recently, another kinematic criterion was proposed for axisymmetric jets (Krieg & Mohseni Reference Krieg and Mohseni2021). Pinch-off starts when the characteristic vortex velocity reaches a characteristic feeding velocity of $2U$ . The characteristic vortex velocity is also commonly referred to as the induced velocity $U^*_{cmax}$ ( $=U_{cmax}/U$ ). Velocity $U^*_{cmax}$ is the maximum velocity along the centreline. The streamwise location of the induced velocity is determined by the streamwise trajectory of the vortex pair. Based on the formation numbers obtained from the circulation criterion, the induced velocity criteria are found to be about $1.75U$ and $2.12U$ for the straight nozzle (figure 11 a) and orifice nozzle (figure 11 b), respectively.

Figure 11. Development of induced velocity for (a) straight nozzle (case SN250) and (b) orifice nozzle (case ON175).

4.4. Maximum circulation criterion

The ratio of the maximum vortex circulation to maximum total circulation $\varGamma ^*_{vmax}/\varGamma ^*_{tmax}$ is defined to examine the efficiency of starting jets on vortex formation. The ratio decreases as the stroke ratio increases due largely to the existence of the trailing jet from small formation time (such as $t^*=1$ ). Based on the formation number obtained from the circulation criterion, the circulation ratio criterion is defined when the ratio reaches a critical value. A universal ratio of about 72 % is found to be appropriate for both nozzle configurations in two-dimensional flows (see figure 12). Therefore, this criterion is not similar to the kinematic criterion due to its independence of the nozzle configuration. The formation number represents a critical circulation ratio for vortex formation in starting jets. This criterion agrees well with the experiment for non-parallel two-dimensional starting jets with a critical ratio of about 71 % (Gao & Yu Reference Gao and Yu2016a ).

Considering the contraction effect, the transformed formation number has been used to unify the classical formation number for the straight nozzle and orifice nozzle (Limbourg & Nedić Reference Limbourg and Nedić2021b ), such as

(4.1) \begin{equation} F_{\star } = \frac {U_{\star }t_{critical}}{H_{\star }} = \frac {Ut_{critical}}{H}\times 1/C^2_c = F\times 1/C^2_c, \end{equation}

where the contraction coefficient $C_c$ is 0.9 suggested for the straight nozzle (Limbourg & Nedić Reference Limbourg and Nedić2021b ) or 0.76 for the orifice nozzle. The contraction coefficient $C_c$ is obtained from the steady simulation and identified from the lowest transversal position along the streamline passing the inner corner of the orifice plate. Quantity $t_{critical}$ is the physical time corresponding to the critical formation time. Based on (4.1), the present formation numbers of 13.6 and 9.3 ( $F$ ) will become 16.8 and 16.1 ( $F_{\star }$ ) for the straight nozzle and orifice nozzle, respectively. Therefore, the transformed formation number ( $F_{\star }$ ) of 16.5 with deviation $\pm 0.4$ can effectively unify the critical formation time for the straight nozzle case and orifice case. In additional, with the use of the contraction coefficient, two kinematics criteria also can be universal for both nozzle configurations. For the first kinematic criterion (§ 4.3), the critical translational velocity of about $0.44U_{\star }$ is universal for the straight nozzle and orifice nozzle. This conclusion agrees well with the experiments for axisymmetric jets (Limbourg & Nedić Reference Limbourg and Nedić2021d ). For the second kinematic criterion (§ 4.3), the critical induced velocity of about $1.60U_{\star }$ is also found to be nearly independent of the nozzle configuration. The deviations for both straight nozzle and orifice nozzle are within $0.03U_{\star }$ .

Figure 12. Circulation ratio versus stroke ratio for (a) straight nozzle and (b) orifice nozzle.

4.5. Effects of Reynolds number

Some studies indicated that the development of vortex circulation (i.e. the formation number and the scaling law for axial and radial trajectories) is nearly independent of the Reynolds number for laminar vortex rings (Didden Reference Didden1979; Gharib et al. Reference Gharib, Rambod and Shariff1998). The insignificant influence on the formation number of laminar vortex rings has been examined by a series of simulations over a range of Reynolds numbers from 500 to 5000. It is noted that the Reynolds number would play an important role to the velocity profile at the nozzle exit plane (Rosenfeld et al. Reference Rosenfeld, Rambod and Gharib1998, Reference Rosenfeld, Katija and Dabiri2009). To evaluate the dependence of formation number on the Reynolds number, several more cases are simulated as listed in table 3. To avoid any possible developments of transitional or turbulent flows, the Reynolds numbers are chosen up to 2000. Due to more shear-layer instabilities for higher Reynolds numbers (e.g. 1600 and 2000), finer meshes with an average grid spacing of $H/180$ in streamwise and transverse directions are used to provide more precise results for the straight nozzle and the orifice nozzle. The evolutions of total circulation and vortex circulation are presented in figure 13. According to the circulation criterion, the formation numbers are identified for different Reynolds numbers, as summarized in table 3.

Table 3. Cases studied for the dependence of the formation number on the Reynolds number.

Figure 13. Effects of the Reynolds number on the evolution of total circulation and vortex circulation for (a) straight nozzle and (b) orifice nozzle.

As can be seen in figure 13(a), the total circulation generated by the straight nozzle would be larger for smaller Reynolds number (Rosenfeld et al. Reference Rosenfeld, Rambod and Gharib1998). This can be attributed to the viscous effect which produces a thicker boundary layer within the nozzle with higher vorticity. Similar reason contributes to higher vortex circulation initially. However, the vortex circulation for case SN250-Re500 gradually becomes smaller than that for other cases with the straight nozzle due to the viscous diffusion. The greater total circulation and smaller maximum vortex circulation for case SN250-Re500 bring down the formation number slightly to 12. As indicated by table 3, the maximum vortex circulation is identical for cases SN250 and SN250-Re1600. Therefore, the influence of Reynolds number on the local maximum of vortex circulation is limited. However, the formation number for case SN250 is smaller due to larger total circulation increase rate. The final increase of vortex circulation generated by the vortex pairing with the secondary vortex appearing after the jet discharge is excluded to identify the formation number for cases SN250-Re1600 and SN250-Re2000. This is because the local maximum of vortex circulation has already been achieved. This phenomenon is different from that of the non-parallel flow cases which still have a continuous increase of vortex circulation before the vortex pairing. The difference may be explained by the stronger shear-layer instability for orifice cases. Stronger secondary vortex for case SN250-Re2000 contributes to a faster pinch-off and smaller vortex circulation in comparison with case SN250-Re1600. However, due to the smaller total circulation increase rate, the formation number for case SN250-Re2000 is maintained. It is noted that the formation number could be increased slightly with the Reynolds number for the straight nozzle due to a decrease in total circulation. In comparison with the formation number of 13.6 for the straight nozzle, the differences produced by the Reynolds numbers from 500 to 2000 are within 12 %. This range agrees well with the comparable difference of about 10 % in experiments and simulations for axisymmetric jets (Gharib et al. Reference Gharib, Rambod and Shariff1998; Rosenfeld et al. Reference Rosenfeld, Rambod and Gharib1998).

Different from the situation for the straight nozzle, the total circulation for the orifice cases could be enhanced with the increasing Reynolds number due to greater transverse velocity gradient at higher Reynolds numbers (see figure 13 b). The increase of vortex circulation becomes smaller as the Reynolds number increases. This leads to similar vortex circulation for cases ON175, ON175-Re1600 and ON175-Re2000 initially and smaller vortex circulation for case ON175-Re500. Therefore, the formation number for case ON175-Re500 is smaller than that at Reynolds number of 1000. Higher Reynolds number induces more shear-layer instabilities (e.g. double vortex pairing for case ON175-Re1600 and triple vortex pairing for case ON175-Re2000). This leads to a faster pinch-off and smaller formation number identified by the local maximum vortex circulation during a constant state after the first vortex pairing. The onset of pinch-off can appear earlier by up to 27 % (in comparison with the formation number of 9.3) due to the appearance of secondary vortices in the trailing jet. This result agrees well with previous findings in simulations for axisymmetric jets (Zhao et al. Reference Zhao, Frankel and Mongeau2000). Additionally, the formation numbers for cases ON175-Re1600 and ON175-Re2000 also agree with the range found by earlier experiments with the orifice nozzle and Reynolds numbers from 1490 to 3380 (Gao & Yu Reference Gao and Yu2016a ).

5. Two-dimensional flows versus axisymmetric flows

After determining the range of formation numbers in the preceding sections, we focus on the difference in the formation number between the axisymmetric and two-dimensional cases. Slower vortex motion in two-dimensional flows has been attributed to be the reason (Pedrizzetti Reference Pedrizzetti2010; Domenichini Reference Domenichini2011). For more quantitative discussions on the differences between axisymmetric and two-dimensional flows, two axisymmetric cases (case ASN100 for the straight nozzle and case AON100 for the orifice nozzle) with stroke ratio $L_m/D=10$ are further simulated for comparisons. As shown in figure 14, the translational velocity of the vortex rings is far greater than the translational velocity of the vortex pairs, as may have been expected. The oscillation of vortex velocity may be attributed to the vorticity redistribution during the pinch-off process and the non-circular geometry of the vortex core at the post-formation process (Danaila & Helie Reference Danaila and Helie2008). To examine the factors influencing the translational velocity, the self-induced velocity of a curved vortex filament can be estimated by the Biot–Savart law (Lim & Nickels Reference Lim and Nickels1995):

(5.1) \begin{equation} \mathbf {u} = \frac {\varGamma \kappa }{4\pi }\ln \left (\frac {R_c}{a} \right ) \mathbf {b}, \end{equation}

where $\kappa$ is the local curvature; $R_c$ is the radius of curvature; $a$ is the radius of the vortex core; and $\mathbf {b}$ is the unit vector along the binormal direction. At the initial formation time ( $t^*\lt 1$ ), the vortex core is small enough to provide a constant term $\ln ({R_c}/{a})$ (Wu, Ma & Zhou Reference Wu, Ma and Zhou2007). The vortex circulation may be evaluated by the total circulation in the initial formation. Based on the numerical results, the vortex pair would have larger circulation but smaller translational velocity during the initial formation in comparison with the vortex ring. This result suggests that the local curvature is more important than the vortex circulation in estimating vortex translational velocity even though they are of the same order in the Biot–Savart law.

Figure 14. Development of (a) normalized streamwise trajectory and (b) normalized translational velocity of the leading vortex.

The velocity formula for the steady vortex can be used to further understand the faster motion of vortex rings. For two-dimensional flows, the translational velocity of a vortex pair considered as two point vortices can be expressed by (Saffman Reference Saffman1992)

(5.2) \begin{equation} U_v = \frac {\varGamma _v}{4 \pi y_v}. \end{equation}

For axisymmetric flows, the vortex ring velocity depends on its core thickness when the curvature of the vortex tube is considered. The translational velocity of a thin-cored vortex ring can be expressed by a second-order formula as (Fraenkel Reference Fraenkel1972)

(5.3) \begin{equation} U_v = \frac {\varGamma _v}{4 \pi y_v} \underbrace {\left \{ \ln \left ( \frac {8}{\varepsilon } \right )-\frac {1}{4}+\varepsilon ^2 \left [\frac {15}{32}-\frac {3}{8}\ln \left ( \frac {8}{\varepsilon } \right ) \right ] \right \}}_{\text {I}}, \end{equation}

where mean core radius $\varepsilon$ is defined as the ratio of the core radius $a$ to ring radius $y_v$ . A universal energy of a steady thin vortex ring was found to be 0.33 (Gharib et al. Reference Gharib, Rambod and Shariff1998; Limbourg & Nedić Reference Limbourg and Nedić2021c ), leading to a specific ratio $\varepsilon$ of about 0.5 in the family of steady vortex rings (Norbury Reference Norbury1973). Thus, the term I in equation (5.3) represents the effect of core size and is equal to approximately 2.4. Based on the results of vortex circulation and transverse trajectory, the velocities at $t^*=12$ and $t^*=30$ are calculated for vortex pairs and vortex rings as listed in table 4. The magnitude of the calculated velocity at both instants approximately agrees with the translational velocity obtained from numerical results (figure 14 b). At $t^*=12$ , the circulation is larger for vortex pairs than for vortex rings while the difference in transverse trajectory is negligible. However, the velocity of the vortex ring is larger. This suggests the important effect of the local curvature.

Table 4. Estimated translational velocity using equations (5.2) and (5.3).

Figure 15. Evolution of normalized total invariants: (a,d) circulation, (b,e) impulse and (c,f) energy for straight nozzle (a–c) and orifice nozzle (d–f).

As expected, total invariants are different between axisymmetric flows and two-dimensional flows. They are reflected in the velocity distribution over the nozzle exit. This influence on total invariants and the accuracy of the contraction-based slug model (equations (3.6)–(3.8) for axisymmetric flows and equations (3.12)–(3.14) for two-dimensional flows) are evaluated in figure 15. The contraction coefficient $C_c$ for the axisymmetric flow is found to be 0.73 from the steady simulation. Therefore, the difference of $C_c$ for the two-dimensional flow and axisymmetric flow is negligible, in accord with previous discussions (Limbourg & Nedić Reference Limbourg and Nedić2021a ,Reference Limbourg and Nedić b ). Here, the results for the axisymmetric flow are approximately predicted by the two-dimensional $C_c$ .

First, an evaluation of circulation prediction is conducted. It is noted that the formulas of circulation prediction are identical for two-dimensional and axisymmetric flows (equations (3.6) and (3.12)). As can be seen in figure 15(a,d), the contraction-based slug model can approximately predict the total circulation of axisymmetric flows at small time (e.g. $t^*\lt 4$ ). A better estimation than that of the classical slug model (corresponding to a smaller slope of invariant growth) can be observed as in a previous study (Limbourg & Nedić Reference Limbourg and Nedić2021a ). However, the development of boundary-layer thickness brings increasing deviation between the numerical result and the contraction-based model. Due to a reduction of boundary-layer thickness, the contraction-based model can appropriately predict total circulation in two-dimensional flows. The relative error for the prediction of the straight nozzle case at $t^*=10$ (which is the end of the piston motion and is around the critical formation time) is $-11\%$ . In addition, the contraction-based slug model can also provide approximate predictions for orifice cases in both two-dimensional and axisymmetric flows. The relative error for the two-dimensional prediction of the orifice case at $t^*=10$ is $-18\%$ .

For the prediction of total impulse $I^*_t=I_t/\rho D^3\it U$ and total energy $E^*_t=E_t/\rho D^3\it U^2$ for axisymmetric flows, the contraction-based slug model (blue dashed lines in figure 15) can generally provide accurate descriptions. Although the initial over-pressure effect on total impulse cannot be included, the error for prediction is insignificant at later time (e.g. $t^*\gt 6$ ). However, due to the two-dimensional nature, the contraction-based slug model for two-dimensional flows (red dashed lines) predicts larger slopes than that of axisymmetric flows for total impulse $I^*_t=I_t/\rho H^2U$ and total energy $E^*_t=E_t/\rho H^2U^2$ . This overestimated behaviour of the slug model is not found in two-dimensional flows. In reality, the total momentum and total energy are lower for two-dimensional flows due to slower vortex motion. In comparison with the difference in total circulation between two-dimensional and axisymmetric flows, larger differences can be observed for total impulse and total energy. This can be attributed to the two-dimensional characteristic for circulation and the three-dimensional characteristic for impulse and energy in axisymmetric flows.

To improve the contraction-based slug model for two-dimensional flows, a modified coefficient $C_m$ is proposed for impulse and energy formulas. Equations (3.13) and (3.14) are rewritten as

(5.4) \begin{equation} I_{\star 2m} = \rho LHU\times 1/C_c\times 1/C_m, \end{equation}

(5.5) \begin{equation} E_{\star 2m} = \frac {1}{2}\rho LHU^2\times 1/C^2_c\times 1/C_m, \end{equation}

where the modified coefficients are 2.2 and 2.0 for the straight nozzle and orifice nozzle, respectively. As shown in figure 15, the modified contraction-based slug model for two-dimensional flows (red dotted lines) can appropriately predict the total impulse and total energy. The relative errors of the total invariants predicted by the modified contraction-based slug model (equations (3.12), (5.4) and (5.5)) can be found in table 5. Errors are obtained by time averaging between $t^*=5$ and $t^*=10$ .

Table 5. Time-averaged relative errors for the modified contraction-based slug model.

6. Analytical estimation of two-dimensional formation number

For axisymmetric flows, an asymptotic matching procedure was proposed to analytically estimate the formation number. The total invariants described by the slug model were matched to a model of steadily isolated vortex rings such as the Norbury–Fraenkel model (Norbury Reference Norbury1973) or Kaplanski–Rudi model (Kaplanski & Rudi Reference Kaplanski and Rudi2005). Ultimately, the critical time scale was obtained by an additional closure assumption, for example, a universal vortex energy in non-dimensional form (Gharib et al. Reference Gharib, Rambod and Shariff1998), a specific vortex velocity normalized by the initial jet velocity (Mohseni & Gharib Reference Mohseni and Gharib1998; Shusser & Gharib Reference Shusser and Gharib2000) or a volume constraint of vortex ring bubble (Linden & Turner Reference Linden and Turner2001). In addition, the contraction-based slug model was used to predict the formation number (Limbourg & Nedić Reference Limbourg and Nedić2021d ).

For two-dimensional flows, a family of steady vortex pairs (Pierrehumbert Reference Pierrehumbert1980) may be used to further predict the formation number. To the best of the authors’ knowledge, the asymptotic matching procedure has not been studied for two-dimensional starting jets. A parameter $\overline {A}_0=A_0/A_1$ was typically used to identify the member in the Pierrehumbert model, where $A_0$ is the minimum distance of the core boundary from the symmetry axis and $A_1$ is the maximum distance from the symmetry axis. However, distance $A_0$ only can be identified when the leading vortex pair moves far away from the trailing jet. This is owing to the influence of the trailing jet during jet discharge from the nozzle. Instead, another geometrical parameter $\gamma =y_v/A_1$ can be used to determine the vortex member in the current study. As shown in figure 16(a), vortex pairs generated by various nozzle configurations and stroke ratios eventually reach a universal value at $\gamma =0.49$ (averaged from all cases). This result generally agrees with the conclusion from figure 7 and the suggestion of a universal vortex ring generated by various starting jets (Gharib et al. Reference Gharib, Rambod and Shariff1998).

Figure 16. (a) Evolution of parameter $y_v/A_1$ and (b) comparison of simulated vortices with those from the Pierrehumbert family.

The Pierrehumbert model provides a one-to-one correspondence between the geometrical parameter $\gamma$ and the non-dimensional vortex energy $E^{**}_v$ (square symbols in figure 16 b). The relation may be approximated by a quartic polynomial (blue dotted line) as

(6.1) \begin{equation} E^{**}_{v} = 10.875\gamma ^4 - 26.434\gamma ^3 + 24.032\gamma ^2 - 8.898\gamma + 1.186. \end{equation}

It is noted that the Pierrehumbert vortex with $\gamma =0.49$ and $E^{**}_v=0.11$ has overestimated the non-dimensional vortex energy of the simulated vortex pairs with differences of 62 % and 45 % for the straight nozzle and orifice nozzle, respectively.

Combining the proposed slug model (equations (3.12), (5.4) and (5.5)) and the definition of vortex energy normalization ( $E^{**}_v=E_v/\rho \varGamma ^2_v$ ), the relation between stroke ratio and vortex energy can be represented as

(6.2) \begin{equation} \frac {L}{H}=\frac {2C^2_c}{C_mE^{**}_v}. \end{equation}

The required stroke ratio to generate a particular vortex pair can be estimated from (6.1) and (6.2). The limiting cases of vortex pairs are a thick axis-touching vortex for $\gamma \rightarrow 0.413$ and a thin vortex for $\gamma \rightarrow 1$ . For the limiting cases at $\gamma \rightarrow 0.413$ , the limiting stroke ratio should be at about 11.5 for the straight nozzle and 9.0 for the orifice nozzle. These results suggest the minimum stroke ratio to generate the theoretically thickest vortex pair.

To estimate the formation number, an additional assumption is required. The pinch-off would start when the vortex pair becomes steady. Two models on the basis of this assumption are proposed in terms of the vortex impulse and vortex velocity.

First, an impulse-based model to predict the onset of pinch-off is proposed based on momentum conservation and the assumption of the steady state of a vortex pair. The vortex impulse (Batchelor Reference Batchelor1967) can be expressed as

(6.3) \begin{equation} I_v = \rho y_v \varGamma _v. \end{equation}

When the vortex pair is assumed to be bounded in a circular vortex bubble moving with velocity $U_v$ , the vortex impulse can be approximately expressed as

(6.4) \begin{equation} I_v = \rho \pi A^2_1 U_v. \end{equation}

Combining equations (5.2), (6.3) and (6.4) yields

(6.5) \begin{equation} \frac {y^2_v}{A^2_1}=\frac {1}{4}. \end{equation}

Therefore, a critical geometrical parameter $\gamma =0.5$ can be obtained from the impulse-based model for a steady vortex pair. This result agrees with the numerical results shown in figure 16(a). Based on equation (6.1), the critical value of vortex energy is at about 0.12. Inserting this result into equation (6.2) produces critical stroke ratios of 6.1 and 4.8 for the straight nozzle and orifice nozzle, respectively. Although the matching process of vortex impulse allows a prediction of the geometrical parameter, the predictions of the formation numbers are far away from the numerical results due to the error of the Pierrehumbert model in terms of the vortex energy.

Second, the velocity-based model proposes the onset of pinch-off when the vortex velocity reaches the jet velocity near the vortex rear boundary (Shusser & Gharib Reference Shusser and Gharib2000). With the assumption of jet half-width equalling the transverse trajectory of the vortex pair, the momentum conservation of jet flows allows

(6.6) \begin{equation} HU=2y_vU_v. \end{equation}

Combining equations (3.12), (5.4) and (6.3) yields

(6.7) \begin{equation} y_v=\frac {2C_cH}{C_m}. \end{equation}

Substituting equation (6.7) into equation (6.6) yields

(6.8) \begin{equation} U_v=\frac {C_mU}{4C_c}. \end{equation}

Combining equations (3.12), (5.2) and (6.7) yields

(6.9) \begin{equation} U_v=\frac {C_mLU}{16\pi C^3_c H}. \end{equation}

Combining equations (6.8) and (6.9) yields

(6.10) \begin{equation} \frac {L}{H} = 4\pi C^2_c. \end{equation}

Therefore, critical stroke ratios of 10.2 and 7.3 are found for the straight nozzle and orifice nozzle, respectively. The relative errors of formation number prediction are presented in table 6. In comparison with the dynamical model in terms of the vortex impulse, the kinematic model would be more appropriate to predict the formation numbers and reduce the errors by half.

Table 6. Relative errors for analytical models.