**5. Dynamics of counter-rotating vortices**

As the final example, we consider the test problem of interacting counter-rotating vortices that involves both their nonlinear dynamics and, as a by-product, sound generation. For small viscosity, the direct simulation of vortex dynamics and acoustics by solving the compressible Navier-Stokes equations on a Eulerian grid is a challenging problem because of the thin vorticity filaments that are generated as the process evolves in time. These are difficult to capture because of numerical dissipation-dispersion problems mentioned in the introduction. In the literature, examples of flow simulations have Reynolds number, as defined based on the velocity circulation, in the range of 1000-4000 (e.g., Inoue, 2002). Eldridge (2007) manages to accurately compute the problem of dynamics and acoustics of counter-rotating vortex pairs at a high Reynolds number, Re=10000 with the use of a Lagrangian vortex particle method. In the latter, the governing fluid flow equations are solved in a non-conservative form and the advected vortex solution is regularly reinitialised on a Eulerian grid to reduce the complexity of thin vorticity filaments and stabilise the solution. In the present subsection, the problem of counter-rotating vortices is solved on a fixed Eulerian grid for the range of Reynolds numbers, Re=5000-10000 with the conservative Navier-Stokes CABARET method.

Fig. 8 shows the problem setup. Four viscous-core counter-rotating vortices are initiated in an open domain. Each of the vortices has a constant velocity circulation at infinity , 0 and a Gaussian distribution of the vorticity with the core radius 0*r*

$$\varphi = \frac{1.25\Gamma}{\pi r\_0^2} e^{-1.25\left(\frac{r}{r\_0}\right)^2} \tag{5}$$

where 0.24 *a* that corresponds to the vortex core Mach number 0 *M* 0.3 and 10 *x y* 3 0 *r* is half-distance between the adjacent vortex centres. In non-dimensional variables, the flow parameters at infinity are taken to be *p* 1, 1 and 0.2 . The initial location of the centre of mass of the system corresponds to **x**=0.

Direct Numerical Simulations of Compressible Vortex Flow Problems 33

resulting velocity field is substituted to the momentum equations in conservation variables that are then integrated numerically with the use of the isentropic flow relation between

(a) (b) Fig. 9. Computational grid: (a) full domain view and (b) zoom in the beginning of the grid

Fig.10 shows the time evolution of vorticty field of the system as the vortex pairs sleep through each other. The grid density is 6 grid cells per vortex radius and the Reynolds number of the CABARET simulation is 9400. This particular Reynolds number is chosen as the best match for the reference Lagrangian particle solution of Eldridge (2007) that corresponded Re=10000. The 5% difference in the Reynolds numbers may be attributed to the differences in numerical approximation of viscous terms in the momentum and energy equations in the conservative CABARET and the non-conservative vortex particle method. For Re=9400, the vortex pairs merge after 4 sleep-through events. The sleep-though events correspond to Fig10(c),(d),(e),(f). For Re=5000, the coalescence happens earlier in comparison with the Re=9400 case: for the lower Re-number the vortices coalesce already after 3 sleepthroughs. This change can be compared with results for the fully inviscid follow case calculation that was conducted with the same CABARET Euler method. In the latter case the vortices manage to undergo 6 sleep-through events before coalescence. The capability of the CABARET model to capture the qualitative differences between the Re=9400 and the fully

inviscid solution is noted as a good indication of how low-dissipative the method is.

To zoom into the flow details, Fig.11 compares instantaneous vorticity contours obtained from the CABARET solution at two grid resolutions, 6 and 12 cells per core radius, with the reference Lagrangian particle solution from Eldridge (2007) at one time moment. This time moment corresponds to the last vortex sleep-through before the coalescence. At this time the vortices are very close and strongly interact with each other through fine vorticity filaments. For all three solutions, the same contour levels are plotted that show a good agreement

It is also interesting to compare the 2D solution computed with the results of the vortical structure visualization obtained experimentally by Bricteux et al 2011 for a 3D high-Reynolds number jet**.** In this experimental work a moving window technique is used in the framework of Particle Image Velocimetry (PIV) method to visualise vortex paring in the jet shear layer. Fig.12 shows the results of the visualisation in the jet symmetry plane that

For a quantitative comparison, the centre mass velocity of the vortex system and the values of first few vortex slip-through periods are compared with the reference analytical solution

appear qualitatively very similar to the results of the 2D simulation (cf. fig.10a,b,d,g).

pressure and density.

stretching zone.

down to a small detail.

for point vortices in inviscid flow.

Fig. 8. Problem configuration for two counter-rotating vortex pairs.

The vortex system undergoes a jittering motion with one vortex pair sleeping through the other and taking turns. The centre of mass travels at a positive subsonic speed, which corresponds to the right horizontal direction in Fig.8. After each sleep-through, the distance between the vortices decreases until they finally coalesce and continue the movement as a single core. By using a point-vortex approximation with neglecting viscous effects, the dynamics of the vortex system before coalescence can be described by the classical analytical solution of Hicks (1922). This solution gives the following expressions for the centre mass velocity *U* and the slip through period *Tp* of the vortex system:

$$\overline{\mathcal{U}} = \frac{\delta\_y \cdot \Gamma}{4\pi \delta\_x^2} \left[ \frac{k^2 E(k)}{E(k) - (1 - k^2) K(k)} \right], \\ T\_p = \frac{32\pi \delta\_x^2}{\Gamma k^2 (1 + k)} \left( E(k) - (1 - k^2) K(k) \right) \tag{6}$$

where *k*=1/[1+(*<sup>y</sup> / x*)2]=0.5 and *Kk Ek* ( ), ( ) are the complete elliptic integrals of the first and second kind. For the specified parameters of the model, the centre mass velocity and slipthrough period are

$$T\_p = \frac{54.46\delta^2}{\Gamma} = 12.21, \ \overline{\mathcal{U}} = \frac{0.1437 \cdot \Gamma}{\delta} = 0.1282\ . \tag{7}$$

The velocity of the centre mass corresponds to the Mach number *M* 0.1083 .

Because of the problem symmetry with regard to axis y=0, one-half of the computational domain is considered with the symmetry boundary condition. The problem is solved in the inertial frame of reference which velocity with regard to the absolute frame equals the velocity of the centre of mass,*U* . The latter is available from the analytical point-vortex solution. Because of this choice the centre of mass is approximately stationary in the reference coordinate system. The latter is helpful for minimising the size of the computational domain.

The computational domain of size 440 x 220 (axially times vertically) is covered by a Cartesian grid that has a uniform grid spacing in the central block (60 x *30*). Exponential grid stretching is applied near the outer boundaries to reduce numerical reflections (Fig.9). Three grid resolutions are considered: 6, 9 and 12 grid cells per vortex core radius.

Initial conditions in conservation variables are computed in the following way. By combining the initial vorticity distribution (5) with the solenoidal velocity field condition, the velocity field is computed from solving the Laplace equation for velocity potential. The 32 Advanced Fluid Dynamics

2

The vortex system undergoes a jittering motion with one vortex pair sleeping through the other and taking turns. The centre of mass travels at a positive subsonic speed, which corresponds to the right horizontal direction in Fig.8. After each sleep-through, the distance between the vortices decreases until they finally coalesce and continue the movement as a single core. By using a point-vortex approximation with neglecting viscous effects, the dynamics of the vortex system before coalescence can be described by the classical analytical solution of Hicks (1922). This solution gives the following expressions for the centre mass

Fig. 8. Problem configuration for two counter-rotating vortex pairs.

velocity *U* and the slip through period *Tp* of the vortex system:

*x*

through period are

computational domain.

The computational domain of size 440

.

The velocity of the centre mass corresponds to the Mach number *M* 0.1083 .

 x 220

Three grid resolutions are considered: 6, 9 and 12 grid cells per vortex core radius.

Cartesian grid that has a uniform grid spacing in the central block (60

22 2

*<sup>y</sup> <sup>x</sup> <sup>p</sup>*

2 2

*kEk U T <sup>E</sup> <sup>k</sup> <sup>k</sup> <sup>K</sup> <sup>k</sup> Ek k Kk k k*

where *k*=1/[1+(*<sup>y</sup> / x*)2]=0.5 and *Kk Ek* ( ), ( ) are the complete elliptic integrals of the first and second kind. For the specified parameters of the model, the centre mass velocity and slip-

<sup>2</sup> 54.46 0.1437 12.21, 0.1282 *T U <sup>p</sup>*

Because of the problem symmetry with regard to axis y=0, one-half of the computational domain is considered with the symmetry boundary condition. The problem is solved in the inertial frame of reference which velocity with regard to the absolute frame equals the velocity of the centre of mass,*U* . The latter is available from the analytical point-vortex solution. Because of this choice the centre of mass is approximately stationary in the reference coordinate system. The latter is helpful for minimising the size of the

grid stretching is applied near the outer boundaries to reduce numerical reflections (Fig.9).

Initial conditions in conservation variables are computed in the following way. By combining the initial vorticity distribution (5) with the solenoidal velocity field condition, the velocity field is computed from solving the Laplace equation for velocity potential. The

( ) <sup>32</sup> , ( ) (1 ) ( ) 4 ( ) (1 ) ( ) (1 )

 

(6)

2

. (7)

2

(axially times vertically) is covered by a

 x *30*

). Exponential

resulting velocity field is substituted to the momentum equations in conservation variables that are then integrated numerically with the use of the isentropic flow relation between pressure and density.

Fig. 9. Computational grid: (a) full domain view and (b) zoom in the beginning of the grid stretching zone.

Fig.10 shows the time evolution of vorticty field of the system as the vortex pairs sleep through each other. The grid density is 6 grid cells per vortex radius and the Reynolds number of the CABARET simulation is 9400. This particular Reynolds number is chosen as the best match for the reference Lagrangian particle solution of Eldridge (2007) that corresponded Re=10000. The 5% difference in the Reynolds numbers may be attributed to the differences in numerical approximation of viscous terms in the momentum and energy equations in the conservative CABARET and the non-conservative vortex particle method.

For Re=9400, the vortex pairs merge after 4 sleep-through events. The sleep-though events correspond to Fig10(c),(d),(e),(f). For Re=5000, the coalescence happens earlier in comparison with the Re=9400 case: for the lower Re-number the vortices coalesce already after 3 sleepthroughs. This change can be compared with results for the fully inviscid follow case calculation that was conducted with the same CABARET Euler method. In the latter case the vortices manage to undergo 6 sleep-through events before coalescence. The capability of the CABARET model to capture the qualitative differences between the Re=9400 and the fully inviscid solution is noted as a good indication of how low-dissipative the method is.

To zoom into the flow details, Fig.11 compares instantaneous vorticity contours obtained from the CABARET solution at two grid resolutions, 6 and 12 cells per core radius, with the reference Lagrangian particle solution from Eldridge (2007) at one time moment. This time moment corresponds to the last vortex sleep-through before the coalescence. At this time the vortices are very close and strongly interact with each other through fine vorticity filaments. For all three solutions, the same contour levels are plotted that show a good agreement down to a small detail.

It is also interesting to compare the 2D solution computed with the results of the vortical structure visualization obtained experimentally by Bricteux et al 2011 for a 3D high-Reynolds number jet**.** In this experimental work a moving window technique is used in the framework of Particle Image Velocimetry (PIV) method to visualise vortex paring in the jet shear layer. Fig.12 shows the results of the visualisation in the jet symmetry plane that appear qualitatively very similar to the results of the 2D simulation (cf. fig.10a,b,d,g).

For a quantitative comparison, the centre mass velocity of the vortex system and the values of first few vortex slip-through periods are compared with the reference analytical solution for point vortices in inviscid flow.

Direct Numerical Simulations of Compressible Vortex Flow Problems 35

<sup>1</sup>

<sup>2</sup>

within the order of accuracy the point-vortex model, (r0/

Re=5000 and 9400 on the grids of different resolution.

adjacent vortex centres.

to *a t* / ~ 110 

*a t* / ~ 90 

control surface at distance of 20

**Grid spacing, h** *r*<sup>0</sup> /6 *r*<sup>0</sup> /9 *r*<sup>0</sup> /12

*U* 0.1271 0.1277 0.1279

*Tp* 10.85 10.9 10.9

*Tp* 7.9 8.1 8.2

from the vortex centre of vortices. The control surface is

)2=0.16 and thus characterises how

Table 1. Integral characteristics of the vortex system as obtained from the numerical solution Here (T*p)1* is the time period between the first and the second vortex sleep-through and (T*p)2*  is the time period before the second and the third sleep-through. The agreement for the meanflow velocity between the point-vortex theory value (0.1282) and the numerical values on 3 different grids is within 1%. For the sleep-through period, the values obtained on different grids are converged within 0.5%. The numerically predicted time period, however, is 10-15% shorter in comparison with the point-vortex theory (12.21). This discrepancy is

non-compact the viscous vortex core is in comparison with the distance between the

In addition to the near-filed, the far-field pressure field has been computed on a circle

located in the same reference system as the centre of mass that moves at a small subsonic speed, *M* 0.1083 with respect to the absolute frame. Fig.13 shows the pressure signals obtained at the control points corresponding to 300 and 900 angle to the flow direction for

The pressure fluctuations are defined with the reference to the pressure field value at infinity, *p* 1 . The peaks of the pressure signatures correspond to the vortex sleepthrough events and the number of the peaks corresponds to the total number of vortex sleep-throughs, respectively. The phase of intense vortex interaction during the vortex pairing is followed by a "calming" period that corresponds to the vortex roll-up after the coalescence. In comparison with the pre-coalescence time history that is dominated by large-

For the higher Re-number case, the amplitude of the last acoustic "burst" that corresponds

 has some 20% higher amplitude in comparison with other peaks. This loud acoustic event corresponds to the last vortex sleep-through, which takes place at

 and which is well-captured on the grids of different resolution. After the vortex coalescence, the increase of Reynolds number from 5000 to 9400 also leads to a notable prolongation of small-scale acoustic fluctuations in the post-coalescence phase. These effects may be associated with the small spatial structures that are generated shortly before the

For the pre-coalescence period of vortex evolution, the numerical solutions that correspond to the grids of different resolution are converged within 1-2% for both Reynolds numbers. For the post-coalescence time history, the grid convergence for the high Reynolds-number case, Re=9400, slows down in comparison with the Re=5000 case. For both Re-number cases, however, the CABARET solution on the grid resolution 12 cells per vortex radius appears

time scales the post-coalescence signal is dominated by small-time-scale events.

vortices coalesce (e.g., fig.10f,g) and which are more sensitive to viscous dissipation.

adequate to capture the fine pressure field fluctuations well.

Fig. 10. Vorticity distribution of leapfrogging vortices at several consecutive time moments: *a t* /0 (a), 18 (b), 36 (c), 54 (d), 72 (e), 90 (f), 108 (g), 120 (h).

Fig. 11. Vorticity distribution of the leapfrogging vortex pairs at *a t* / 90 for (a) CABARET solution with the grid density of 6 cells per vortex core radius, (b) CABARET solution on the grid with 12 cells per vortex core radius, and (c) the reference vortex particle method solution from Eldridge (2007).

Fig. 12. PIV of vorticity distribution of leapfrogging vortex rings in the symmetry plane of a high Reynolds number jet obtained with a moving window technique by Bricteux et al 2011.

34 Advanced Fluid Dynamics

(а) (b) (c)

(d) (e) (f)

(g) (h)

Fig. 10. Vorticity distribution of leapfrogging vortices at several consecutive time moments:

(a) (b) (c)

CABARET solution with the grid density of 6 cells per vortex core radius, (b) CABARET solution on the grid with 12 cells per vortex core radius, and (c) the reference vortex particle

Fig. 12. PIV of vorticity distribution of leapfrogging vortex rings in the symmetry plane of a high Reynolds number jet obtained with a moving window technique by Bricteux et al 2011.

for (a)

Fig. 11. Vorticity distribution of the leapfrogging vortex pairs at *a t* / 90

method solution from Eldridge (2007).

0 (a), 18 (b), 36 (c), 54 (d), 72 (e), 90 (f), 108 (g), 120 (h).

*a t* /


Table 1. Integral characteristics of the vortex system as obtained from the numerical solution

Here (T*p)1* is the time period between the first and the second vortex sleep-through and (T*p)2*  is the time period before the second and the third sleep-through. The agreement for the meanflow velocity between the point-vortex theory value (0.1282) and the numerical values on 3 different grids is within 1%. For the sleep-through period, the values obtained on different grids are converged within 0.5%. The numerically predicted time period, however, is 10-15% shorter in comparison with the point-vortex theory (12.21). This discrepancy is within the order of accuracy the point-vortex model, (r0/)2=0.16 and thus characterises how non-compact the viscous vortex core is in comparison with the distance between the adjacent vortex centres.

In addition to the near-filed, the far-field pressure field has been computed on a circle control surface at distance of 20 from the vortex centre of vortices. The control surface is located in the same reference system as the centre of mass that moves at a small subsonic speed, *M* 0.1083 with respect to the absolute frame. Fig.13 shows the pressure signals obtained at the control points corresponding to 300 and 900 angle to the flow direction for Re=5000 and 9400 on the grids of different resolution.

The pressure fluctuations are defined with the reference to the pressure field value at infinity, *p* 1 . The peaks of the pressure signatures correspond to the vortex sleepthrough events and the number of the peaks corresponds to the total number of vortex sleep-throughs, respectively. The phase of intense vortex interaction during the vortex pairing is followed by a "calming" period that corresponds to the vortex roll-up after the coalescence. In comparison with the pre-coalescence time history that is dominated by largetime scales the post-coalescence signal is dominated by small-time-scale events.

For the higher Re-number case, the amplitude of the last acoustic "burst" that corresponds to *a t* / ~ 110 has some 20% higher amplitude in comparison with other peaks. This loud acoustic event corresponds to the last vortex sleep-through, which takes place at *a t* / ~ 90 and which is well-captured on the grids of different resolution. After the vortex coalescence, the increase of Reynolds number from 5000 to 9400 also leads to a notable prolongation of small-scale acoustic fluctuations in the post-coalescence phase. These effects may be associated with the small spatial structures that are generated shortly before the vortices coalesce (e.g., fig.10f,g) and which are more sensitive to viscous dissipation.

For the pre-coalescence period of vortex evolution, the numerical solutions that correspond to the grids of different resolution are converged within 1-2% for both Reynolds numbers. For the post-coalescence time history, the grid convergence for the high Reynolds-number case, Re=9400, slows down in comparison with the Re=5000 case. For both Re-number cases, however, the CABARET solution on the grid resolution 12 cells per vortex radius appears adequate to capture the fine pressure field fluctuations well.

Direct Numerical Simulations of Compressible Vortex Flow Problems 37

Bogey, C. and Bailly, C., "A family of low dispersive and low dissipative explicit schemes for flow and noise computations", J. Comput. Physics, 194 (2004), pp. 194-214. Bogey C. and Bailly C., "Influence of nozzle-exit boundary-layer conditions on the flow and acoustic fields of initially laminar jets", J. Fluid Mech., Vol.25, 2010, pp507-540. Boris, J.P., Book, D.L., and Hain, K., "Flux-corrected transport: Generalization of the

Bricteux, L., Schram C., Duponcheel M., Winckelmans, G., "Jet flow aeroacoustics at

Cockburn B and Shu CW, "Runge–Kutta Discontinuous Galerkin Methods for Convection-Dominated Problems", Journal of Scientific Computing, 16(3): 173-261, Sept 2001. Colonius T., Lele S.K., and Moin P., "The scattering of sound waves by a vortex: numerical simulations and analytical solutions", J Fluid Mech (1994), 260, pp 271-298. Colonius T and Lele SK., "Computational aeroacoustics: progress on nonlinear problems of sound generation." Progress in Aerospace sciences, 2004, 40, pp. 345-416. Dritschel, D.G., Polvani, L.M. and Mohebalhojeh, A.R.: The contour-advective semi-

Eldridge, J.D., "The dynamics and acoustics of viscous two-dimensional leapfrogging

Fritsch G. and Giles M., "Second-order effects of unsteadiness on the performance of

Georges, T.M. "Acoustic ray paths through a model vortex with a viscous core",

Goloviznin V.M. and Samarskii, A.A. "Difference approximation of convective transport

Goloviznin V.M. and Samarskii, A.A., "Some properties of the CABARET scheme",

Goloviznin, V.M. "Balanced characteristic method for systems of hyperbolic conservation

Harten, A., Engqist, B., Osher, S., and Chakravarthy, S. "Uniformly High Order Accurate Essentially Non-Oscillatiry Schemes III", J. Comput. Phys, 71, (1987), pp. 231-303. Hicks, W.M. "On the mutual threading of vortex rings", Proceedings of the Royal Society of

Hirsh, C., "Numerical computation of internal and external flows", vol. 2, John Wiley &

Howe, M.S., "The generation of sound by aerodynamic sources in an inhomogeneous steady

Inoue, O "Sound generation by the leapfrogging between two coaxial vortex rings", Physics

Iserles, A. "Generalized Leapfrog Methods", IMA Journal of Numerical Analysis, 6 (1986), 3,

Karabasov, S.A., Berloff, P.S. and Goloviznin, V.M. "CABARET in the Ocean Gyres", J.

with spatial splitting of time derivative", Mathematical Modelling, Vol. 10, No 1,

J.Acoust.Soc. of America, Vol. 51, No. 1 (Part 2) pp. 206-209 (1972).

Mathematical Modelling, Vol. 10, No 1, 1998, pp. 101–116.

flow", J. Fluid Mech., Vol 67, No. 3, 1975, pp. 597-610.

laws", Doklady. Mathematics, 2005, vol. 72, no1, pp. 619-62313.

Aeroacoustics Conference), 6-8 June, Portland, Oregon, 2011

Re=93000: comparison between experimental results and numerical predictions", AIAA-2011-2792, 17th AIAA/CEAS Aeroacoustics Conference (32nd AIAA

Lagrangian algorithm for the shallow water equations. Mon. Wea. Rev. 127(7), pp.

method", J. Comput. Phys, 31,(1975), 335-350.

vortices", J. Sound Vib., 301 (2007) 74–92.

turbomachines", ASME Paper 92-GT-389; 1992.

**7. References** 

1551–1565 (1999).

pp. 86-100.

Sons, 1998.

381-392.

London A 10 (1922) 111–131.

of Fluids 14 (9) (2002) 3361–3364.

Ocean Model., 30 (2009), рр. 155–168.

Fig. 13. Acoustic pressure signals at different observer angles to the flow: 300 (a),(c) and 900 (b),(d) for Re=5000 (a),(b) and Re=9400 (c),(d).
