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

SHOCK-LIKE FOCUSING OF INERTIAL WAVES - THE LOCALIZEDGENERATION OF TURBULENCE

Jie Liu
Ahmed Mohamed
Chair of Fluid Dynamics,
Fluid mechanics and acoustics laboratory
Department of Mechanical Engineering
Ecole Centrale de Lyon
Technical University of Darmstadt
69130 Ecully, France
64287 Darmstadt, Germany
ahmed.mohamed@ec-lyon.fr
liu@fdy.tu-darmstadt.deAlexandre Delache
Fabien Godeferd
Fluid mechanics and acoustics laboratory Fluid mechanics and acoustics laboratory
Ecole Centrale de Lyon
Ecole Centrale de Lyon
69130 Ecully, France
69130 Ecully, France
alexandre.delache@ec-lyon.fr
fabien.godeferd@ec-lyon.fr

Yongqi Wang
Martin Oberlack
Chair of Fluid Dynamics,
Chair of Fluid Dynamics,
Department of Mechanical Engineering Department of Mechanical Engineering
Technical University of Darmstadt
Technical University of Darmstadt
64287 Darmstadt, Germany
64287 Darmstadt, Germany
wang@fdy.tu-darmstadt.deoberlack@fdy.tu-darmstadt.de

ABSTRACT

turbulence. The optimal angle θ p is determined here by using
numerical simulation and theoretical analysis.

In this work, the focusing of inertial waves emitted by avertically vibrating torus in a uniformly rotating fluid has been
investigated by numerical simulation and theoretical analysis.
The primary objective of this study is to explore the localized
generation of turbulence around the focal point.

GOVERNING EQUATIONS
We consider a uniformly rotating fluid with angular velocity Ω , density ρ and kinematic viscosity ν. For low-amplitude
disturbances, the departures of pressure p and velocity u fromthe undisturbed rotating flow are governed by the linearized
equations of motion in the rotating system (Greenspan, 1968)

CONTEXTIn the ocean, due to stably stratified fluid, the energy injected by tidal forces at the surface generates some internal
gravity waves, which interact with the topography. The transfer of this energy into local turbulence occurs due to a focus of internal gravity waves through two mechanisms (Klymak et al., 2012): wave-topography interaction (remote turbulence) and wave-wave interaction away from the topography(local turbulence). An idealized configuration such as an oscillating torus can reproduce this focus of waves at the apexof the propagation cone with an angle θ (see figure 1). Several studies using internal gravity waves or inertial waves in arotating fluid have analyzed this torus configuration by using
numerical simulations of inertial waves (Duran-Matute et al.,2013), experiments of internal gravity waves (Shmakova &Flór, 2019), or the theory of internal gravity waves (Voisin
et al., 2011). Based on the theoretical work on an oscillating
sphere by Voisin et al. (2011), Shmakova & Flór (2019) assumed a preferential angle θ p at which the amplification of theamplitude in the focal region is maximal for internal gravity
waves. At this angle θ p , the maximum energy of the waves istransferred to the focus zone and could be converted into local∂u
1Ω × u = − ∇p + ν∇2 u + f ,
+ 2Ω
∂t
ρ
∇ · u = 0.

(1)
(2)

Here p represents the reduced pressure, which includes the hydrodynamic pressure and the centrifugal force. It is known thatinertial wave can be generated in the rotating flow as a resultof the restoring action of the Coriolis force. The propagation
direction of the wave energy is determined only by the ratio
between the forcing frequency ω f and the rotation rate Ω ofthe fluid as following (Sagaut & Cambon, 2008)
θIW = cos−1 (ω f /2Ω).

(3)

In this work, by using a vertically oscillating ring of Diractype to simulate the oscillation of a torus, we force the inertial waves to focus into the apex of the propagation cone (see112th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19–22, 2022
and the rotation rate Ω of the fluid. This equation can be analytically solved using the Fourier transform. Substituting thesolution for p into the momentum equation (1) yields the solution for velocity u . The focusing phenomenon caused by thethe forcing f is included in the particular solutions, which is
the key topic of the present analysis.
Considering the forcing condition (5), the wave-like equation
for the pressure (4) in a viscous fluid is given by2∂∂2p− Ek∇2 ∇2 p + 4 2
∂t
∂z


2 
2
=Λ0 4 − σ − iEk∇
δ (r − 1)δ ′ (z)e−iσt .
Figure 1. Schematic drawing of the oscillating torus and thefocusing effect of the resulting inertial waves, where a is theradius of the tube and b the major radius of the torus witha ≪ b, Ω is the constant angular velocity of the rotating flow.
c g denotes the group velocity, at which the wave energy propagates, and θ is the propagation angle.

Using the Fourier transform in cylindrical coordinate system,
we obtain the inhomogeneous solution for the pressure in
Fourier spaceh2 i
4 − σ + iEkκ 2
kz J0 (kr )p̂(kr , kz ,t) = i2πΛ0 e−iσt. (7)2
σ + iEkκ 2 κ 2 − 4kz2

figure 1).
Several non-dimensional numbers have been used to characterize the oscillation of the torus in a similar way to what has been
done in stratified fluid (Ermanyuk et al., 2017): the KeuleganCarpenter number Ke = Ab , which represents the ratio of theforcing amplitude to the minor radius of the torus, the Stokes

Additionally, the transformation of the momentum equation(1) in z-direction into Fourier space results in the following
scalar equation

ω f b2
ωf A
ν , the Rossby number Ro = 2bΩ (Davidson,ω f A·2b2013), and the Reynolds number Re = ν . We adapt this

number St =number in our case of a ring with Dirac-type forcing by changing b into a.
Next, the theoretical analysis and numerical simulation ofwave focusing will be presented in detail.

(8)

Substituting (7) into (8), the inhomogeneous solution ofuz in Fourier space follows as−iσtûz = i2πΛ0 eLet b, Ω−1 , U characterize the typical length, time andrelative fluid velocity to the rotating flow. With this, the linearized equation system can be reduced to a dimensionless∂2
form. Furthermore, applying the operator ( ∂t
2 ∇·) to (1) andtaking the continuity equation (2) into account lead to a single
wave-like equation for the pressure p in the formJ0 (kr )kr2 σ + iEkκ 2.
2
σ + iEkκ 2 κ 2 − 4kz2

(9)

Applying the Fourier transform in the cylindrical coordinate
system to the wave-like equation (4) and then inserting the solution of the pressure p into the momentum equation (1), we
obtain the analytical solution of the vertical velocity
Λ0 ∞ 3
k J0 (kr )J0 (rkr ) · ei|z|α f1 e−|z|α f2
(10)
2 0 r
2
2
2
2
[σ (4 − σ ) f1 + 4Ekkr f2 ] + i[4Ekkr f1 − σ (4 − σ ) f2 ]×dkr
σ (4 − σ 2 )3/2 ( f12 + f22 )uz = e−iσt2∂− Ek∇2 ∇2 p + (2nn · ∇)2 p = 4nn · (nn · ∇) f +
∂t



2
∂
∂
2nn ·
− Ek∇2 (∇ × f ) +− Ek∇2 (∇ · f ), (4)
∂t
∂t

Z

where α = √ σ4−σ 2

where n is a unit vector indicating the direction of the rotation
axis. For the following analysis, we define n as the z-axis. Theν
Ekman number Ek = Ωb
2 is a gross measure of the viscous
force compared to the Coriolis force.
In our case, see figure 1, the forcing f has only one component
in z-direction, which can be described using the Dirac functionin the cylindrical coordinate system as
fz (r, ϕ, z,t) = Aδ (r − 1)δ (z)e−iσt .

∂+ Ekκ 2 ûz = −ikz p̂ + fˆz ,∂twhere fˆz = 2πΛ0 J0 (kr )e−iσt .

METHODS
Theoretical analysis



(6)

, and

vs
u 
2


u
kr t
64Ek2 kr4
1024Ek2 kr4
64Ek2 kr4
f1,2 = √1− 2+ 2± 1− 2
.
2 )3
2 )4
2 )3
σ
(4
−
σ
σ
(4
−
σ
σ
(4
−
σ
2

Numerical simulations
Besides the analysis in the previous section, the governing equations (1) and (2) are solved using the standard pseudo
spectral method in a periodic domain with 2π length cube andspatial resolution of 5123 points.
We perform 25 linear direct numerical simulations DNS
by varying the angle of inertial waves propagation to know if
there is a propagation angle maximizes the energy in the focus(5)Here A is a constant line density of the force along the torus,
σ is the ratio between the oscillation frequency ω f of the torus

212th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19–22, 2022
becomes increasingly blurred, the amplitude of the waves isgreatly reduced and the position of the maximum amplitude ofuz moves in the negative direction of the z-axis.
In order to be able to analyze the viscous effect described in(10) more intuitively, we further simplify it by introducing theassumption of small Ekman number, i.e. Ek ≪ 1. Using Taylor expansion and neglecting all terms with O(Ek2 ), f1 and f2
can be simplified to the following approximations
f1 ≈ kr ,(11a)

16Ekkr3f2 ≈.
σ (4 − σ 2 )2

(11b)

Correspondingly, the solution of uz in (10) can be simply expressed as
Figure 2. Inertial-wave cone excited by the oscillating diractorus as depicted by isosurface of uz according to the analyticalsolution (10) for inviscid fluid, i.e. Ek= 0.

uz ≈ e−iσt

Z k̄r
0

3

Aei(|z|αkr −ψ) e−|z|αβ kr dkr

(12)

with
We keep the viscosity constant ν = 10−4 , the ring ra-zone.
dius a = 1, the amplitude of forcing A = 0.01 and frequency ofω A·aforcing ω f = 5, in that case Re = fν = 500, the KeuleganCarpenter number Ke = 0.01 ≪ 1 where the forcing amplitude
is assumed to be smaller than the radius b, and the Stokes number St = 50000 ≫ 1 as a result of the choice of small viscosity.
Nevertheless we vary the angle of cone in the range of θ = 10◦
to 70◦ by varying the rotating rate 2Ωz = 5.077 to 14.61 andωf A
at the same time Ro = aΩvaries from 0.019 to 0.0068 (table 1). Note that, inertial waves appear in low Reynolds andlow Rossby regime, where in the regime of high Reynolds andlow Rossby, waves-turbulence interactions occur according toGodeferd & Moisy (2015).
In the second set of simulations, we perform 5 linear simulations using a solid body with a three-dimensional torus
shape. By using a penalisation method as done by DuranMatute et al. (2013), we force the torus body to oscillate vertically and we vary the propagation angle.
The parameters of this set is chosen to be, the viscosity ν = 10−5 , the amplitude of forcing A = 0.0157 and thefrequency of forcing ω f = 1.256. With using torus geometry as following, minor radius of b = 0.0942 and major radiusa = 0.848, by varying the rotating rate 2Ωz = 1.996 to 2.44,
we vary the angle θ between 51◦ and 59◦ , see table 1. The
resulting non-dimensional numbers, i.e. the Reynolds number Re = 372, the Rossby number which varies between 0.085
and 0.105, the Stokes number St = 1115, and the KeuleganCarpenter number Ke = 0.167.

Λ0 α 2
k J0 (kr )J0 (rkr ),2σ r
4σ Ekkr2
ψ = arctan,(4 − σ 2 )2
16Ek
β=
.
σ (4 − σ 2 )2
A=

(13a)
(13b)
(13c)

Here A denotes the amplitude as in the inviscid fluid, ψ represents the phase shift and β a factor of viscous effect. It is to be
noted that both ψ and β are affected not only by the viscosity
but also by the forcing frequency. The larger ψ and β are, thewider is the wave beam and weaker is the focusing effect.
At the focal point, where r = 0 and |z| = 1/α, the effectof viscous attenuation is actually dominated by the factor
3e−β kr , which decays rapidly as kr increases. Fig. 4 shows,3
the exponential function e−β kr converges to zero at a smaller
wavenumber kr with increasing value of Ek as well as σ → 0
and σ → 2. For this we postulate the existence of a cutoff wavenumber k̄r so that viscous attenuation eliminates all
waves with the wavenumbers that are significantly greater thank̄r . In numerical evaluation, an attenuation factor e−c can bepredefined, where the constant c is used to control the precision of the integration result. As an example, the choice of
c = 10 corresponds to an attenuation factor e−c ≈ 4.5 × 10−5 ,which may be sufficiently small for applications. The value of3
cutoff wavenumber k̄r of the function e−β kr corresponding to−c
the chosen attenuation factor e can be reasonably estimatedby

RESULTS
Viscous spreading on the focus zone 1/3 
1/3ccσ (4 − σ 2 )2
=
.
k̄r =
β
16Ek

We evaluated the integral with respect to kr in (10) numerically and plotted the solution. Figure 2 shows the isosurfaceof uz for the inviscid case Ek= 0. Figure 3 shows the distribution of the vertical velocity amplitude according to (10)in the vertical plane (r, z) for σ = 1 and the Ekman number
Ek= 0, 10−5 , 10−4 and 10−3 , respectively. Due to symmetry,
we only present the propagation of the wave beam in one quadrant. The result shows that as the value of the Ekman number
increases, the effect of viscous attenuation becomes stronger,
which leads to a rapid dissipation of energy during the propagation. As a result, with increasing viscosity the focal zone

(14)

Equation (14) shows that for the same non-dimensional forcing frequency σ , the cutoff wavenumber k̄r decreases with ascaling proportional to Ek1/3 . Furthermore, for the same valueof Ek, the cutoff wavenumber k̄r also tends to zero as the σ → 0
and σ → 2, which corresponds to disappearance of focusing
phenomenon.
Fig. 5 shows the change of the amplitude of the vertical velocity at the focal point with the Ekman numbers Ek and the3

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

ω A

DNS

2Ωz

θ

f
Ro = (a)Ω

1

5.077

10◦

0.019

5.32

20◦

0.018

3

5.77

30◦

0.017

4

6.52

40◦

0.0153

5
..
.

6.62
..
.

41◦
..
.

0.0151
..
.

24

10

60◦

0.01

25

14.61

70◦

0.0068

2

ω A

DNS

2Ωz

θ

f
Ro = (b)Ω

26

1.9968

51◦

0.105

27

2.088

53◦

0.1

28

2.1908

55◦

0.095

29

2.307

57◦

0.09

2.44

59◦

0.085

30

Table 1. Simulation parameters for Linear DNS for ring of Dirac forcing on the left from number 1 to number 25, and solid body
torus forcing on the right from number 26 to number 30.

45

4.5

40

4

35

3.5

30

3

25

2.5

20

2

15

1.5

10

1

5

0.5

3
2.5
2

0.7
0.6
0.5
0.4

1.5

0.3
1

0.2
0.5

0.1

Figure 3. Amplitude of the vertical velocity in the vertical plane according to equation (10) for σ = 1 in a viscous fluid with the
Ekman number (a) Ek= 0, (b) Ek= 0.10−5 , (c) Ek= 10−4 and (d) Ek= 10−3 , respectively.

3

Figure 4. Plot of the function e−β kr defined in equation (13c) with the change of σ and kr for the Ekman number (a) Ek = 10−5 , (b)Ek = 10−4 and (c) Ek = 10−3 , respectively.

z direction. Furthermore, σ → 2 corresponds also to the casewith an infinitely large viscous effect as described in (13c),so that the energy is rapidly dissipated during propagation andcannot be focused effectively. Additionally, For forcing frequency close to the angular velocity of rotating fluid, i.e. σ
is close to 1, the amplitude of the vertical velocity at the focal point reaches its maximum, which corresponds to a wave

dimensionless forcing frequency σ . It is easy to understandthat the focusing phenomenon disappears when σ tends to 0,
because in this case the waves propagate in vertical direction
along the cylindrical surface on which the oscillating torus islocated according to Eq. (3) and thus do not focus on the z-axis.
However, the disappearance of wave focusing at σ → 2 is limited to the case where the external force is present only in the412th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19–22, 2022

Figure 5. Amplitude of the vertical velocity at the focal point
in a viscous fluid with a slight lowering compared its positionin an inviscid fluid (r = 0, z = 1/α) as a function of the dimensionless frequency σ from 0 to 2 and the Ekman number
Ek from 10−5 to 10−3 .

Figure 6. 2D cut of the vertical velocity uz contour from DNS
using the parameters number 19 in table 1.

propagation angle slightly less than 60◦ due to viscous effects.
Since the theoretical solution is very complicated, this preferential angle will be further investigated by numerical simulation.

Transfer to vertical velocity
To compare with the theoretical results, we first simulate the vertical velocity distribution. Figure 6 shows a twodimensional cut of the vertical velocity uz contour in the vertical plane, where the rays of inertial waves form a double conewith two focal points. This phenomenon is consistent with thetheoretical results. Furthermore, for each simulation forced
with a vertically oscillating ring of Dirac forcing, we extract
the vertical velocity uz of the focus zone as a function of time.
The theoretical position of the focus point is Z f = a · tan(θ ),so for each simulation we extract a small cube of the velocitiesof the theoretical focus point and its neighborhood, the cubehas 4 grid cells in each direction. For each point, we obtain thesignal of vertical velocity. Extracting the amplitudes and thenaveraging them, we obtain the amplitude reported on figure 7
for each angle θ . We observe that there is a vertical velocity
focusing maximum where the angle θ = 55◦ maximizes thevertical velocity at the focus zone. For the case of torus bodyforcing simulations, the averaged vertical velocity amplitude
for each simulation is obtained, and a maximum of the focusing of the velocity is observed in the simulation at θ = 57◦ .
In the stratified fluid the dispersion relation is given as
θIGW = ±N sin−1 (ω/N), where N is the Brunt–Väisälä frequency, while in the rotating fluid the dispersion relation isθIW = ± cos−1 (ω/2Ω) as defined in equation (3). In the workof Shmakova & Flor (2019), it was demonstrated that the angle at which the maximum energy is extracted into the focus
zone in the case of a stratified fluid is θIGW = 35◦ for internal
gravity wave, so it is coherent to have θIW = θIGW − 90◦ .

Figure 7. The normalized amplitude of the vertical velocity
component uz in the focus zone as a function in the waves
propagation angle θ . Comparison of the analytical solution in(10) for Ek= 10−4 with DNS of the parameters from number
1 to number 25 listed in table 1.
major radius. The geostrophic mode appears also in the ocean,
where the Coriolis force is balanced by the pressure gradient(Davidson, 2013).

CONCLUSION AND OUTLOOKBoth our analytical and numerical result show that therays of inertial wave emitted by a vertically vibrating torus
form a double cone symmetric about the torus with two focal points, where the wave amplitude has a maximum due towave focusing. This phenomenon is in good agreement withthe experimental and numerical study by Duran-Matute et al.
(2013). Furthermore, we have analyzed the the effect of viscous spreading on the focus zone. The result shows that asthe Ekman number increases, the wave beam becomes wider,
the focusing effect becomes weaker, and focus zone moves towards the negative direction of the z- axis compared with thecase of an inviscid fluid. Moreover, our numerical simulationshows that the inertial waves are able to extract the maximum
energy of focusing at an angle of θIW = 55◦ . Based on the
results obtained in this work, we are presently analyzing the
focusing of inertial waves in the fully nonlinear regime andexplore the localized generation of turbulence around the fo-

Geostrophic mode
The solution of the vertical velocity using the torus bodyforcing shows that a stationary geostrophic mode develops in
the flow, vertically elongated above and below the torus body,see figure 8. This geostrophic mode is dominating in the mean
flow. In the mean flow, we observe two clockwise and counterclockwise rotating cylinders with radius similar to the torus

512th International Symposium on Turbulence and Shear Flow Phenomena (TSFP12)
Osaka, Japan (Online), July 19–22, 2022
Press.
Duran-Matute, M., Flór, J.-B., Godeferd, F. S. & Jause-Labert,
C. 2013 Turbulence and columnar vortex formation throughinertial-wave focusing. Phys. Rew. E 87, 041001(R).
Ermanyuk, Evgeny V, Shmakova, ND & Flór, J-B 2017 Internal wave focusing by a horizontally oscillating torus. Journal of Fluid Mechanics 813, 695–715.
Godeferd, Fabien S & Moisy, Frédéric 2015 Structure and dynamics of rotating turbulence: a review of recent experimental and numerical results. Applied Mechanics Reviews
67 (3).
Greenspan, H. P. 1968 The theory of rotating fluids. Cambridge
University Press.
Klymak, J. M., Legg, S., Alford, M. H., Buijsman, M., Pinkel,
R. & Nash, J. D. 2012 The direct breaking of internal waves
at steep topography. Oceanography society 25, 041001(R).
Sagaut, P. & Cambon, C. 2008 Homogeneous turbulence dynamics. Cambridge University Press.
Shmakova, N. D. & Flór, J.-B. 2019 Nonlinear aspects of focusing internal waves. J. Fluid Mech. 862, R4.
Voisin, B., Ermanyuk, E. & Flór, J. 2011 Internal wave generation by oscillation of a sphere, with application to internal
tides. J. Fluid Mech. 666, 308–357.

Figure 8. Vertical averaged velocity on the vertical plane andhorizontal averaged velocity on the horizontal plane from DNS
using the parameters number 29 in table 1.

cal point.

REFERENCES
Davidson, Peter Alan 2013 Turbulence in rotating, stratified
and electrically conducting fluids. Cambridge University

6