12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19-22, 2022

SELF-EXCITED OSCILLATIONS IN VARIABLE DENSITYCOUNTER-CURRENT ROUND JETS

Karol Wawrzak
Department of Thermal Machinery
Czestochowa University of Technology
Al. Armii Krajowej 21, 42-201 Czestochowa, Poland
karol.wawrzak@pcz.pl

Andrzej Boguslawski
Department of Thermal Machinery
Czestochowa University of Technology
Al. Armii Krajowej 21, 42-201 Czestochowa, Poland
andrzej.boguslawski@pcz.pl

ABSTRACTThe paper focuses on a global instability phenomenon invariable density counter-current round jets. The analysis isperformed using a large eddy simulation and the computations
are carried out employing a high-order numerical code. The
configuration of two co-axial nozzles, where suction appliedin an annular nozzle is a driving force for the counterflow,
is investigated for different velocity and density ratios. The
research have shown the emergence of two global modes depending of the velocity and density ratios.

INTRODUCTION
Round jets have been studied theoretically, experimentally as well as numerically for decades. Their widespreadpresence in many industrial devices has inspired scientists
to find a way to control them to obtain desired properties.
One of the most interesting phenomena emerging in round
jets, which considerably influences their behavior, is the selfexcited global instability phenomenon triggered by absolutely
unstable local flow regions (Huerre & Monkewitz, 1990). The
theoretical predictions of Monkewitz & Sohn (1988) and Jendoubi & Strykowski (1994) showed that the absolute instability can be induced in variable density and counter-current jet
configurations. Jendoubi & Strykowski (1994) identified two
distinct absolutely unstable modes called Mode I and ModeII. The main features of the modes are: (i) Mode I has maximum fluctuation in the shear layer region, while Mode II has
maximum on the jet centreline; (ii) both modes become moreunstable as the velocity ratio I = U2 /U1 (U1 - velocity of thejet, U2 - velocity of the counterflow) is increased, the density
ratio S = ρ jet /ρambient is decreased or the Mach number is decreased; (iii) Mode I is absolutely unstable only provided thatI > 0 whereas Mode II does not require any counterflow below
critical density ratio Scr ≈ 0.7.
The results of linear stability analysis were confirmed
in experimental (Strykowski & Niccum, 1991; Strykowski& Wilcoxon, 1993; Kyle & Sreenivasan, 1993; Hallberg &Strykowski, 2006) and numerical works (Lesshafft et al., 2006,
2007; Foysi et al., 2010; Boguslawski et al., 2016; Wawrzak

Figure 1.

Schematic view on computational configuration.

et al., 2021). Despite the fact that the research cited above
univocally showed the occurrence of the global instability inthe counter-current jets, there is a lack of an analysis demonstrating a clear distinction between Mode I and Mode II in
the available literature. Moreover, there are no investigations
devoted to hot round jets with the counter-flow to the best authors’ knowledge. The present paper is aimed at extending theknowledge concerning the global instability in variable density

1

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19-22, 2022

1.2

1.2

S1.0-I0.000
S1.0-I0.050
S1.0-I0.100
S1.0-I0.150
S1.0-I0.200
S1.0-I0.250
S1.0-I0.275
S1.0-I0.300

jet
1.0

1

0.8
0.8

Ux/U1

Ux/U1

0.6
0.6

0.4
0.4

0.2

co - flow
0.2

suction

0.0
-0.2
0.0

0

0.5

r/D

1.0

0

5

1.5

x/D

10

15

(a)
0.4

(a)

S1.0-I0.000
S1.0-I0.050
S1.0-I0.100
S1.0-I0.150
S1.0-I0.200
S1.0-I0.250
S1.0-I0.275
S1.0-I0.300

0.35

1.2
0.3

0.25

ρ/ρambientux’/U1

1.0

0.2

0.15

0.8

0.1

0.6

0.05

0

0.4

0.0

0

5

x/D

10

15

(b)
0.5

1.0

Figure 3. Time-averaged mean (a) and fluctuating (b) velocity profiles along the jet axis for S = 1.

1.5

r/D
(b)

Figure 2.

scales resulting from the LES filtering are modelled using the
sub-grid model of Vreman (Vreman, 2004).
The calculations are performed for the configuration
shown schematically in Fig. 1, which corresponds to the experimental configuration of Strykowski & Niccum (1991). In
the present simulations, the internal geometry of the nozzles
is not included in the computational domain. Instead, the jet
is modeled through velocity and density profiles at the inlet
plane. The computational domain is a simple rectangular box
with dimensions Lx = 15D, Ly = Lz = 6D, where D is the diameter of the central main jet and x is the axial direction.
The inlet boundary conditions are specified in terms of themean density and instantaneous velocity profiles. The inlet velocity profile is computed by adding fluctuating velocity components to the mean values. The velocity fluctuations are computed applying the digital filtering method proposed by Klein
et al. (2003) whereas the inlet mean velocity is described by
the hyperbolic-tangent profile. The inlet turbulence intensity is
assumed at the level of 0.001U1 . A low-velocity co-flow equal
to 0.05U1 is added outside the jet region to mimic natural entrainment from the surroundings. The inlet density profile isobtained from the Busemann-Crocco relation (Schlichting &Gersten, 1979). The mean velocity and density profiles as-

Mean velocity (a) and density (b) profiles.

jets by adding a counterflow.

NUMERICAL ALGORITHM AND COMPUTATIONAL PARAMETERS
We consider an incompressible variable density flow forwhich the filtered continuity and Navier-Stokes equations are
solved using large eddy simulation (LES). The LES calculations are performed using an academic high-order solver
SAILOR based on the low Mach number approximation. The
code employed in this study was used previously in various
computations including isothermal, non-isothermal, excitedand reacting jets (Boguslawski et al., 2017). A solution algorithm is based on a projection method for the pressurevelocity coupling. The time integration is performed bya predictor-corrector (Adams-Bashforth/Adams-Moulton) approach. The model equations are discretized in space using
the 6th order compact difference method for half-staggered
meshes (Tyliszczak, 2014; Tyliszczak A., 2016).The sub-grid2

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19-22, 2022

Amplitude/U1

10

0

10

-1

10

-2

x/D= 1.0

1.2

S1.0-I0.000
S1.0-I0.050
S1.0-I0.100
S1.0-I0.150

10-3

S0.7-I0.000
S0.7-I0.050
S0.7-I0.100
S0.7-I0.150
S0.7-I0.200
S0.7-I0.250
S0.7-I0.300

1

10-4
0.8

-5

0.0

0.5

Amplitude/U1100

1.0 St

D

1.5

2.0

x/D= 1.0

10-1

-3

10

-4

10

-5

0.0

3.0

S1.0-I0.200
S1.0-I0.250
S1.0-I0.275
S1.0-I0.300

10-2
10

2.5

Ux/U1

10

0.6

0.4

0.2

0.5

1.0 St

D

1.5

2.0

2.5

3.0
0

Figure 4. Evolution of the spectra of the axial velocity fluctuations registered at the jet axis for S = 1.

0

x/D

5

10

15

(a)
0.4

sumed at the inlet plane are presented in Fig. 2.
The Reynolds number defined as ReD = U1 D/ν, where νis the kinematic viscosity, is equal to ReD = 5000. The analysis is performed for a jet with the shear layer thickness characterized by the parameter D/θ = 40, various velocity ratios
I = 0 − 0.3, and density ratios equal to S = 0.5, 0.7, 1. It should
be noted that a relatively thick shear layer has been chosen to
avoid interactions with the self-sustained convective oscillations reported by Boguslawski et al. (2013, 2019) that appear
for D/θ > 50.
The computations were carried out on the mesh consisting of 192 × 128 × 128 points in x, y, z directions, respectively.
Preliminary computations showed that such a mesh density is
sufficient to obtain virtually grid independent results.

S0.7-I0.000
S0.7-I0.050
S0.7-I0.100
S0.7-I0.150
S0.7-I0.200
S0.7-I0.250
S0.7-I0.300

0.35

0.3

ux’/U1

0.25

0.2

0.15

0.1

0.05

0

0

x/D

5

RESULTS

10

15

(b)

As mentioned above the analysis was performed for threedensity ratios. The first density ratio S = 1 corresponds to thehomogeneous jet for which a counterflow (I > 0) is required
to trigger both global modes. The next density ratio S = 0.7 isslightly above the critical density ratio Scr below which global
Mode II do not require any counterflow (I = 0) for a jet withshear layer characterized by D/θ = 40. The last density ratio
S = 0.5 is below critical Scr for global Mode II without counterflow.
Figure 3 presents the profiles of the time-averaged axial
mean and fluctuating velocity obtained for S = 1 and velocity ratio varying from I = 0 to 0.3. Note that the particular cases are hereafter denoted as SS − II , e.g. S1.0 − I0.10
refers to S = 1.0 and I = 0.1. It can be seen that the jets with
I = 0.0 − 0.15 are characterised by the mean velocity profiles
weakly dependent on I parameter. The profiles are nearly identical and characteristic for a circular free jet with the potential
core length approximately equal to 5 − 6D. Despite a more
rapid decay of the mean velocity for I = 0.2 and I = 0.25, all
the profiles are similar. A stronger impact of velocity ratio I
on the mean profiles is observed for I = 0.275 and I = 0.3.
Local minima at the axial distances x/D ≈ 2 appear for boththe cases. Changes of the mean velocity profiles are connectedwith alteration of the profiles of the velocity fluctuations. For
I = 0 − 0.25 their shapes are similar and typical for round free
jets. In these cases, the fluctuations start growing at x/D ≈ 3
and achieve a level of 18% of U1 . A significant change is observed for the two last I = 0.275, 0.3, for which the growth of

Figure 5. Time-averaged mean (a) and fluctuating (b) velocity profiles along the jet axis for S = 0.7.

Amplitude/U1

100

x/D= 1.0

10-1
10

-2

10

-3

10

-4

10-5
0.0

Amplitude/U1

10

0.5

1.0 StD

1.5

S0.7-I0.000
S0.7-I0.050
S0.7-I0.100
S0.7-I0.150

2.0

2.5

3.0

0

10

-1

10

-2

x/D= 1.0

S0.7-I0.200
S0.7-I0.250
S0.7-I0.300

10-3
10-4
10-5
0.0

0.5

1.0 St

D

1.5

2.0

2.5

3.0

Figure 6. Evolution of the spectra of the axial velocity fluctuations registered at the jet axis for S = 0.7.

the fluctuations starts right from the nozzle exit. Additionally,
for I = 0.3 the profile of fluctuations is characterised by the
local maximum that correspond to the local minimum of the3

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19-22, 2022
10

Amplitude/U1

1.2
S0.5-I0.000
S0.5-I0.050
S0.5-I0.100
S0.5-I0.150
S0.5-I0.200
S0.5-I0.250
S0.5-I0.300

1

10
10

-2

x/D= 1.0

10-4
10

-5

0.0

0.5

100

Amplitude/U1

10-1

0.4

0.2

0

5

x/D

10

10

-2

10

-3

10

-4

10

-5

0.0

1.0 St

D

1.5

2.0

x/D= 1.0

0.6

0S0.5-I0.000
S0.5-I0.050
S0.5-I0.100
S0.5-I0.150

10-3

0.8

Ux/U1

0

-1

0.5

1.0 St

D

1.5

2.5

3.0

S0.5-I0.200
S0.5-I0.250
S0.5-I0.300

2.0

2.5

3.0

Figure 8. Evolution of the spectra of the axial velocity fluctuations registered at the jet axis for S = 0.5.

15

(a)
0.4
S0.5-I0.000
S0.5-I0.050
S0.5-I0.100
S0.5-I0.150
S0.5-I0.200
S0.5-I0.250
S0.5-I0.300

0.35

0.3

profile is characteristic for an emergence of global Mode II in
the flow. For higher velocity ratios I = 0.2, 0.25 and I = 0.3
one can observe the damping of the fluctuations at the jet centerline which indicates a transition from global Mode II to
global Mode I. It should be noted that S = 0.7 is only slightly
above the critical density ratio below which global Mode II
does not require any counterflow for the jet with D/θ = 40.
For S = 0.7 the velocity ratio I = 0.15 is a critical velocity
ratio for global Mode II while the I = 0.2 is a critical velocity
ratio for global Mode I. The obtained results are in good agreement with linear stability analysis that showed that an increasein velocity ratio I decreases the critical density ratio for bothglobal modes.
Looking at the spectra presented in Fig. 6 one can observe that a decrease in the density ratio slightly lowers thenon-dimensional frequency StD in the case of global Mode I.
Global Mode I is characterised by StD = 0.75, 0.81 and 0.9
for S0.7 − I0.2 , S0.7 − I0.25 and S0.7 − I0.3 . Such an observation
is in good agreement with theoretical predictions of Jendoubi& Strykowski (1994). Regarding global Mode II (S0.7 − I0.15 )the frequency is lower than for Mode I and characterised by
StD = 0.6.
Figure 7 presents the profiles of the time-averaged axial
mean and fluctuating velocity obtained for the last density ratio considered, i.e. S = 0.5. In this case, as expected, the global
instability is observed for all the velocity ratios considered. As
can be deduced from the velocity profiles both global modes
appeared in the flow. For the cases S0.5 − I0.0 , S0.5 − I0.05 andS0.5 − I0.10 global Mode II is observed while starting fromI = 0.15 there is a transition from Mode II to Mode I. Thechange of the S to the lower value again decreased the criticalvelocity ratio for both modes. Another important dependence
is well seen in the evolution of spectra presented in Fig. 8,
namely the decrease in density ratio lower the frequency of themode. However, this dependence is rather weak. Taking into
consideration the dependence of the frequency on the velocity
ratio one can observe that for all the density ratio S = 0.5, 0.7.1
the increase in I parameter causes an increase in StD value forboth modes. A similar impact of density and velocity ratios
on the frequency of the global modes was predicted by linear
stability analysis (Jendoubi & Strykowski, 1994).
Figure 9 shows iso-surfaces of the Q-parameter and theaxial velocity contours for both global modes. The Qparameter is commonly used to indicate organised vorticalstructures and is defined as Q = 1/2(|Ω|2 − |S|2 ), where Ωux’/U1

0.25

0.2

0.15

0.1

0.05

0

0

5

x/D

10

15

(b)

Figure 7. Time-averaged mean (a) and fluctuating (b) velocity profiles along the jet axis for S = 0.5.

mean velocity profile. Such profiles are characteristic for occurrence of the global Mode I in the jet flow as was presentedby Strykowski & Niccum (1991) and Wawrzak et al. (2021).
Therefore I = 0.275 can be considered as a critical velocity
ratio Icr for global Mode I for S = 1.
The emergence of the global oscillations is confirmed by
Fig. 4 presenting the amplitude spectra of the axial velocity
fluctuations registered at the axial location x/D = 1.0. The
results obtained for the cases S1.0 − I0.275 and S1.0 − I0.30 reveal the existence of well-defined distinct peaks that indicate
the formation of strong vortical structures. The peaks occur at
StD = 0.85 and 0.96 for the cases S1.0 − I0.275 and S1.0 − I0.30 ,respectively.
Figure 5 presents the profiles of the time-averaged axial
mean and fluctuating velocity obtained for S = 0.7 and velocity
ratio varying from I = 0 to 0.3. In this case, the sudden change
in the time-averaged profiles is observed for a lower velocity
ratio than in the case of S = 1, i.e. for I = 0.15. Hence, theI = 0.15 can be considered as a critical velocity ratio. However, the fluctuating velocity profile for the case S0.7 − I0.15
qualitatively differ from those obtained for S = 1. The perturbation growth is observed close to the nozzle exit and the fluctuations reaches 30%U1 at the axial distance x/D ≈ 4. Such a4

12th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19-22, 2022
is identified. In this case, the mixing processes are drasticallyintensified already from the inlet plane. This may be desirablefrom a practical application point of view, e.g. in combustion
chamber, where the fuel/oxidiser mixing is important.

CONCLUSION
The paper presents the LES study on counter-current variable density jets characterized by three density ratios: (i) S = 1
much above; (ii) S = 0.7 very close to; (iii) S = 0.5 below thecritical value for which the global instability is triggered in
the jet without the counterflow. Wide range of the velocity
ratio I controlling the suction were considered that allow us
to observe both globally unstable modes. Consistently to theliterature, the critical velocity ratio for both modes decreases
along with a decreae of the density ratio. Nevertheless, thevariable density counter-current round jets still require significanlty deeper analysis.

ACKNOWLEDGEMENTS
The research was supported by the Polish National Science Centre, Project No. 2016/21/B/ST8/00414. The computations were carried out using PL-Grid Infrastructure.The
authors are grateful for stimulating discussions to Professor
Luc Vervisch and Professor Laurent Gicquel concerning possible applications of the global modes in counter-current jets for
combustion enhancement performed within the project fundedby the National Agency for Academic Exchange (NAWA)
within the International Academic Partnerships Programme
No. PPI/APM/2019/1/00062.

(a)

REFERENCES
Boguslawski, A., Tyliszczak, A., Drobniak, S. & Asendrych,
D. 2013 Self-sustained oscillations in a homogeneousdensity round jet. J. Turbul. 14 (4), 25–52.
Boguslawski, A., Tyliszczak, A., Wawrzak, A. & Wawrzak, K.
2017 Numerical simulation of free jets. Int . J. Numer. Meth.
Heat Fluid 27 (5), 1056–1063.
Boguslawski, A., Tyliszczak, A. & Wawrzak, K. 2016 Largeeddy simulation predictions of absolutely unstable roundhot jet. Phys. Fluids 28 (2), 025108.
Boguslawski, A., Wawrzak, K. & Tyliszczak, A. 2019 A new
insight into understanding the crow and champagne preferred mode: a numerical study. J. Fluid Mech. 869, 385–
416.
Foysi, H., Mellado, J.P. & Sarkar, S. 2010 Large-eddy simulation of variable-density round and plane jets. Int. J. Heat
Fluid Flow 31 (3), 307–314.
Hallberg, M.P. & Strykowski, P.J. 2006 On the universalityof global modes in low-density axisymmetric jets. J. Fluid
Mech. 569, 493.
Huerre, P. & Monkewitz, P.A. 1990 Local and global instabilities in spatially developing flows. Annu. Rev. Fluid Mech.
22 (1), 473–537.
Jendoubi, S. & Strykowski, P.J. 1994 Absolute and convective
instability of axisymmetric jets with external flow. Phys.
Fluids 6 (9), 3000–3009.
Klein, M., Sadiki, A. & Janicka, J. 2003 A digital filter based
generation of inflow data for spatially developing direct numerical or large eddy simulations. J. Comput. Phys. 186 (2),652–665.
Kyle, D.M. & Sreenivasan, K.R. 1993 The instability and(b)

Figure 9. Iso-surfaces of instantaneous value of the Qparameter with the axial velocity contours at the main crosssection for the cases: S0.5 − I0.00 (a) and S0.5 − I0.20 (b).

is the vorticity tensor and S is the rate of strain tensor. The
case S0.5 − I0.00 , for which global Mode II is observed, is characterized by a formation of strong vortices close to the nozzle.
They pair up and consequently break up further downstream.
The formation of vortices even closer to the nozzle and a significantly faster break-up of the structures is observed for thesecond case presented in Fig. 9(b) for which the global Mode I

512th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19-22, 2022
due to global oscillations in jets with annular counterflow.
AIAA J. 31 (3), 564–570.
Tyliszczak, A. 2014 A high-order compact difference algorithm for half-staggered grids for laminar and turbulent incompressible flows. J. Comput. Phys. 276, 438–467.
Tyliszczak A. 2016 High-order compact difference algorithm
on half-staggered meshes for low Mach number flows.
Comput. Fluids 127, 131–145.
Vreman, A.W. 2004 An eddy-viscosity subgrid-scale modelfor turbulent shear flow: Algebraic theory and applications.
Phys. Fluids 16 (10), 3670–3681.
Wawrzak, K., Boguslawski, A. & Tyliszczak, A. 2021 A numerical study of the global instability in counter-current homogeneous density incompressible round jets. Flow Turbul.
Combust. pp. 1–35.

breakdown of a round variable-density jet. J. Fluid Mech.
249, 619–664.
Lesshafft, L., Huerre, P. & Sagaut, P. 2007 Frequency selection
in globally unstable round jets. Phys. Fluids 19 (5), 054108.
Lesshafft, L., Huerre, P., Sagaut, P. & Terracol, M. 2006 Nonlinear global modes in hot jets. J. Fluid Mech. 554, 393–
409.
Monkewitz, P.A. & Sohn, K. 1988 Absolute instability in hot
jets. AIAA J. 26 (8), 911–916.
Schlichting, H & Gersten, K 1979 Boundary-layer theory newyork: Macgraw hill.
Strykowski, P.J. & Niccum, D.L. 1991 The stability of countercurrent mixing layers in circular jets. J. Fluid Mech. 227,309–343.
Strykowski, P.J. & Wilcoxon, R.K. 1993 Mixing enhancement

6