Geophysical Journal International
Geophys. J. Int. (2012) 191, 508–516

doi: 10.1111/j.1365-246X.2012.05611.x

GJI Geomagnetism, rock magnetism and palaeomagnetism

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

On the reﬂection of Alfve´n waves and its implication for Earth’s core modelling
Nathanae¨l Schaeffer,1 Dominique Jault,1,2 Philippe Cardin1 and Marie Drouard1
1ISTerre, Universite´ de Grenoble 1, CNRS, F-38041 Grenoble, France. E-mail: nathanael.schaeffer@ujf-grenoble.fr 2Earth and Planetary Magnetism Group, Institut fu¨r Geophysik, Sonnegstrasse 5, ETH Zu¨rich, CH-8092, Switzerland
Accepted 2012 July 11. Received 2012 July 11; in original form 2011 December 19
SUMMARY Alfve´n waves propagate in electrically conducting ﬂuids in the presence of a magnetic ﬁeld. Their reﬂection properties depend on the ratio between the kinematic viscosity and the magnetic diffusivity of the ﬂuid, also known as the magnetic Prandtl number Pm. In the special case, Pm = 1, there is no reﬂection on an insulating, no-slip boundary, and the incoming wave energy is entirely dissipated in the boundary layer.
We investigate the consequences of this remarkable behaviour for the numerical modelling of torsional Alfve´n waves (also known as torsional oscillations), which represent a special class of Alfve´n waves, in rapidly rotating spherical shells. They consist of geostrophic motions and are thought to exist in the ﬂuid cores of planets with internal magnetic ﬁeld. In the geophysical limit Pm 1, these waves are reﬂected at the core equator, but they are entirely absorbed for Pm = 1. Our numerical calculations show that the reﬂection coefﬁcient at the equator of these waves remains below 0.2 for Pm ≥ 0.3, which is the range of values for which geodynamo numerical models operate. As a result, geodynamo models with no-slip boundary conditions cannot exhibit torsional oscillation normal modes.
Key words: Numerical solutions; Dynamo: theories and simulations; Rapid time variations; Core, outer-core and inner-core; Planetary interiors.

1 INTRODUCTION
Hannes Alfve´n ﬁrst showed the theoretical existence, in an inviscid ﬂuid of inﬁnite electrical conductivity, of hydromagnetic waves that couple ﬂuid motion and magnetic ﬁeld (Alfve´n 1942). The propagation of torsional Alfve´n waves in the Earth’s ﬂuid core was, thereafter, predicted by Braginsky (1970). Such waves arise in rapidly rotating spheres or spherical shells in the presence of a magnetic ﬁeld. In torsional Alfve´n waves, the motions are geostrophic and consist in the rotation ωg(s) of nested cylinders centred on the rotation axis. They, thus, depend only on the distance s to the rotation axis. The period of the fundamental modes of torsional Alve´n waves in the Earth’s ﬂuid core was ﬁrst estimated to be about 60 yr. This timescale was inferred from the analysis of the decadal length of day changes since the ﬁrst half of the 19th century (Jordi et al. 1994) and of the geomagnetic secular variation after 1900 (Braginsky 1984). With hindsight, these time-series were not long enough to show convincingly variations with 60 yr periodicity. Torsional waves with much shorter periods have now been extracted from time-series of core surface ﬂows for the time interval 1955–1985 (Gillet et al. 2010). If this discovery is conﬁrmed, the period of the fundamental modes is of the order of 6 yr and, as such, is much shorter than initially calculated.
Several authors have searched for torsional Alfve´n waves in geodynamo simulations. Using stress-free boundary conditions,
508

Dumberry & Bloxham (2003) and Busse & Simitev (2005) illustrated some parts of the torsional wave mechanism. Dumberry & Bloxham (2003) found that the whole length of the geostrophic cylinders accelerates azimuthally as if they were rigid. The inertial forces, in their simulation, are however, so inﬂuential that they dominate the Lorentz forces. Torsional Alfve´n waves (TAW) have ﬁnally been detected in a set of numerical simulations of the geodynamo with no-slip boundary conditions, for 0.5 ≤ Pm ≤ 10, by Wicht & Christensen (2010) (the magnetic Prandtl number Pm is the ratio of kinematic viscosity over magnetic diffusivity). In both the geophysical (Pm ∼ 10−5) and the numerical studies, there seems to be no reﬂection of the TAW upon their arrival at the equator. However, experimental studies in liquid metals have shown resonance effects on Alfve´n normal modes (Jameson 1964) as well as reﬂection of wave packets (Alboussie`re et al. 2011).
In this paper, we elaborate on the remark that reﬂection of Alfve´n waves is controlled not only by the boundary condition, but also by the magnetic Prandtl number of the ﬂuid in which they propagate (see Jameson 1961, p. 23,24). In the next section, we discuss the governing equations for 1-D Alfve´n waves and the associated boundary conditions for a solid and electrically insulating wall. We remark that for Pm = 1 all the energy of the incident Alfve´n wave is dissipated in a boundary layer, resulting in no reﬂected wave. In the following section, we change geometry to further emphasize our point and brieﬂy present a direct numerical simulation of propagation
C 2012 The Authors Geophysical Journal International C 2012 RAS

and reﬂection of Alfve´n wave in a non-rotating spherical shell. That introduces the section devoted to the geophysical application, where we investigate TAW in the Earth’s core, modelled as a rapidly rotating spherical shell, calculating the energy loss on reﬂection at the Equator as a function of Pm. Finally, we discuss the implications concerning the ability of geodynamo simulations to produce torsional eigenmodes and waves which are expected in the Earth’s core.

2 REFLECTION OF ONE-DIMENSIONAL A L F V E´ N WAV E S

We introduce the problem through the example of Alfve´n waves, transverse to a uniform magnetic ﬁeld in an homogeneous and electrically conducting ﬂuid, hitting a solid wall perpendicular to the imposed magnetic ﬁeld (Roberts 1967). The imposed uniform magnetic ﬁeld B0 is along the x-axis, whereas the induced magnetic ﬁeld b(x, t) and the velocity ﬁeld u(x, t) are transverse to this ﬁeld, along y. Assuming invariance along y- and z-axes, the problem reduce to a 1-D problem, u and b depending only on x. Projecting the Navier–Stokes equation and the induction equation on the y direction (on which the pressure gradient and the non-linear terms do not contribute), one obtains the following equations:

∂t u

=

B0 μ0ρ

∂x

b

+

ν∂xx u,

(1)

∂t b

=

B0∂x u

+

1 μ0σ

∂xx b,

(2)

where μ0 is the magnetic permeability, ρ is the ﬂuid density, ν the kinematic viscosity and σ the electrical conductivity.

2.1 Elsasser variables Introducing the two Elsasser variables h± = u ± b/√μ0ρ, the equation of momentum (1) and the equation of magnetic induction
(2) can be combined into

η+ν

ν−η

∂t h± ∓ VA∂x h± − 2 ∂xx h± = 2 ∂xx h∓,

(3)

where VA = B0/√μ0ρ is the Alfve´n wave speed, and η = (μ0σ )−1

is the magnetic diffusivity. It is already apparent that when ν = η,

the right-hand side of the previous equation vanishes, in which case

h+ and h− are fully decoupled. One can also show that h− travels in

the direction of the imposed magnetic ﬁeld, whereas h+ travels in the opposite direction.

Introducing a length scale L and the timescale L/V A, the previous

equations take the following non-dimensional form:

1

1 Pm − 1

∂t h± ∓ ∂x h± − S ∂xx h± = S Pm + 1 ∂xx h∓,

(4)

where the Lundquist number S and the magnetic Prandtl number Pm are deﬁned as

S = 2VA L η+ν

Pm

=

ν η

.

The propagation of Alfve´n waves requires that the dissipation is
small enough, which is ensured by S 1. The fact that (Pm − 1)/(Pm + 1) = −(Pm−1 − 1)/(Pm−1 +
1) establishes a fundamental symmetry of these equations: when changing Pm into Pm−1, only the sign of the coupling term (right-
hand side of eq. 4) changes.

C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

Reﬂection of Alfve´n waves 509

2.2 Physical boundary conditions and reﬂection of Alfve´n waves

These equations must be completed by boundary conditions. We
assume that the wall is electrically insulating, and that the ﬂuid ve-
locity vanishes at the solid boundary (no-slip boundary condition), which translate to b = 0 and u = 0, leading to h± = 0.
For Pm = 1 the equations for h+ and h− are fully decoupled, regardless of the value of S

1

∂t h± = ±∂x h± + S ∂xx h±.

(5)

In addition, for an insulating solid wall, the boundary condition h± = 0 does not couple h+ and h− either. As a result, reﬂection is not allowed at an insulating boundary when Pm = 1, because reﬂection requires change of travelling direction, and thus transformation of
h+ into h− and vice versa. The energy carried by the wave has to be dissipated in the boundary layer.
For Pm = 1 the equations are coupled: for very small diffusivities (that is large Lundquist number S), the coupling will be effective
only in a thin boundary layer. In addition the coupling will be
more efﬁcient as Pm is further from 1. This gives a mechanism
for reﬂection of Alfve´n waves on an insulating boundary when Pm = 1. Before giving a numerical illustration, it is instructive to consider the boundary conditions in the two limits Pm = 0 and Pm = ∞, with S 1 (dissipationless interior).
In the limit Pm = 0, there is no viscous term and the boundary condition, at the wall x = x0, reduces to

b(x0, t) = 0 ⇒ h+(x0, t) = h−(x0, t).

(6)

There is perfect reﬂection. The incident (+) and reﬂected (−) waves have equal velocities and opposite magnetic ﬁelds. This also corresponds to a stress-free boundary condition for the velocity ﬁeld in combination with an insulating wall (inﬁnitely small vorticity sheet at the wall), leading to perfect reﬂection regardless of the value of Pm used in eq. (4). In this case the boundary condition for the velocity ﬁeld is ∂xu = 0, which translates into ∂x(h+ + h−) = 0 and h+ − h− = 0, effectively coupling h+ and h−.
In the limit Pm = ∞, the boundary condition, at the wall x = x0, reduces instead to

u(x0, t) = 0 ⇒ h+(x0, t) = −h−(x0, t).

(7)

The incident and reﬂected waves have opposite velocities and equal magnetic ﬁelds. This also corresponds to a no-slip boundary condition for the velocity ﬁeld in combination with a perfectly conducting wall (inﬁnitely small current sheet at the wall), leading to perfect reﬂection regardless of the value of Pm used in eq. (4). In this case, the boundary condition for the magnetic ﬁeld is ∂xb = 0, which couples h+ and h−.
Another combination of boundary conditions inhibits reﬂection for Pm = 1: for a stress-free (∂xu = 0) and perfectly conducting wall (∂xb = 0), which translates into ∂xh+ = 0 and ∂xh− = 0, the ﬁelds h+ and h− are decoupled, as for a no-slip insulating wall. Note ﬁnally that a wall with ﬁnite conductivity will allow some weak reﬂection, as illustrated by Fig. 4(h).

2.3 Numerical simulations
We have performed a numerical simulation in a channel 0 ≤ x ≤ x0 with a 1-D ﬁnite difference scheme. The Lundquist number is chosen large enough so that dissipation can be neglected in the interior. The boundary conditions were set to be electrically insulating

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

510 N. Schaeffer et al.

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

Figure 1. Reﬂection coefﬁcient for a 1-D Alfve´n wave packet hitting an insulating boundary with√normal incid√ence, as a function of Pm and for different magnetic Lundquist numbers Lu = V aL/η. The theoretical value for plane waves R(Pm) = (1 − Pm)/(1 + Pm) ﬁts the numerical simulation results perfectly.

and no-slip. The grid is reﬁned next to the boundaries, to have at ltehaicskt n4epsosiδnt=s in√eνaηc/hVbAou(snedearAyplpaeynerd,iwx hAic).h are Hartmann layers of
From the simulation of the travelling wave, we compute the transmission coefﬁcient as the ratio of the velocity amplitude of the reﬂected and incident waves for different values of Pm and S. The results are reported on Fig. 1.
As expected, there is full dissipation for Pm = 1 and energy conservation for Pm 1 or Pm 1. Furthermore, the reﬂection coefﬁcient R is independent of S, and exhibits the expected symmetry R(Pm−1) = −R(Pm). The measured values of R√match perfec√tly the theoretical reﬂection coefﬁcient R(Pm) = (1 − Pm)/(1 + Pm) derived for plane waves, because R depends neither on the pulsation ω, nor on the wavenumber k (see Appendix A).
3 REFLECTION OF A LOCALIZED A L F V E´ N WAV E PA C K E T O N A S P H E R I C A L B O U N DA RY
The peculiar case where no reﬂection occurs is not speciﬁc to the planar, 1-D ideal experiment. Here, we run an axisymmetric simulation in a spherical shell permeated by a non-uniform magnetic ﬁeld, without global rotation. The imposed magnetic ﬁeld is the same as in Jault (2008), and is represented by the dashed ﬁeld lines of Fig. 2. Contrary to the simplest case of the previous section, it is a non-uniform magnetic ﬁeld, which is not perpendicular to the

boundaries. The observed behaviour of Alfve´n wave packets hitting the curved boundaries should therefore apply to many systems.
The numerical pseudo-spectral code is the one used in Gillet et al. (2011), but restrained to axisymmetry. It uses the SHTns library (Schaeffer 2012) for spherical harmonic expansion (Legendre polynomials) in the latitudinal direction, and second-order ﬁnite differences in radius with many points concentrated near the boundaries. It time steps both induction and momentum equation in the sphercial shell using a semi-implicit Crank–Nicholson scheme for the diffusive terms, whereas the coupling and (negligible) nonlinear terms are handled by an Adams–Bashforth scheme (second order in time). The number of radial gridpoints is set to 500 and the maximum degree of Legendre polynomials to 120.
The Alfve´n wave packets are generated mechanically by spinning the conducting inner core for a very short duration (compared to the Alfve´n propagation time). Since the imposed magnetic ﬁeld strength is not uniform, the wave front deforms as it propagates along the ﬁeld lines. When the wave packet hits the outer insulating spherical shell, it does reﬂect and propagates back towards the inner shell for Pm = 0.1 and Pm = 10 but there is no reﬂection for Pm = 1. This is illustrated by the snapshots of Fig. 2.
4 REFLECTION OF TORSIONAL A L F V E´ N WAV E S
Finding evidence of propagation of TAW in the Earth’s ﬂuid core may open a window on the core interior. Properties of TAW in the

C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

Reﬂection of Alfve´n waves 511

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

Figure 2. Snapshot of the azimuthal velocity component of Alfve´n waves propagating in a non-rotating spherical shell. The dashed-lines are the imposed magnetic ﬁeld lines. From left-hand panel to right-hand panel: Panel (a) the incoming waves travelling from the inner shell to the outer shell along magnetic ﬁeld lines; Panel (b) case Pm = 0.1, S = 1800 showing reﬂection with the same sign; Panel (c) case Pm = 1, S = 1000 with total absorption at the wall; Panel (d) case Pm = 10, S = 1800 showing reﬂection with opposite sign.

Earth’s core have, thus, been thoroughly investigated after the initial study of Braginsky (1970). They have been recently reviewed by Jault (2003) and Roberts & Aurnou (2011).

immediately synchronized with the rotation of the solid outer sphere in the course of a spin-up experiment. This is equivalent to a no-slip boundary, as for the 1-D wave case with Pm → ∞.

4.1 Model of torsional Alfve´n waves

To model TAW, magnetic diffusion and viscous dissipation are ne-
glected in the interior of the ﬂuid. The Earth’s ﬂuid core is modelled
as a spherical shell of inner radius ri, outer radius ro and rotation rate . Rapid rotation introduces an asymmetry between the veloc-
ity and magnetic ﬁelds and makes the velocity geostrophic, provided that λ ≡ V A/ ro 1 (Jault 2008). Note that the Lehnert number λ is about 10−4 in the Earth’s core. Geostrophic velocity in a spherical shell consists of the rotation ωg(s) of nested cylinders centred on the rotation axis. It, thus, depends only on the distance s from the
rotation axis (in ro units). A 1-D wave equation for the geostrophic velocity sωg(s) is obtained after elimination of the magnetic ﬁeld b

L ∂2ωg(s) = ∂

∂t2

∂s

L

V˜A2

∂

ωg (s ∂s

)

,

(8)

√ with L = s3H(s) and H (s) = 1 − s2 the half-height of the

geostrophic cylinders, and V˜A2 involves only the z-average of the

squared s-component of the imposed magnetic ﬁeld. Braginsky

(1970) derived (8) rigorously in the geophysical case for which

the viscous Ekman layer is thin compared to the magnetic diffusion

layer located at the top and bottom rims of the geostrophic cylinders.

This condition amounts to Pmλ 1. Then, the velocity remains

geostrophic in the magnetic diffusion layer. We have written the eq.

(8) in its simplest form, when the imposed magnetic ﬁeld is ax-

isymmetric, the mantle is insulating and Ekman friction at the rims

of the geostrophic cylinders is neglected. The eq. (8) needs to be

completed by two boundary conditions, which can be derived when

either Pm 1 or Pm 1.

Interestingly, the eq. (8) may be valid in the limit Pm 1 but also

in the limit Pm 1 (provided Pmλ 1). In the speciﬁc case Pm

1, the appropriate boundary condition on the geostrophic velocity at

the equator (on the inner edge of the Hartmann boundary layer) can

be inferred from the boundary condition on the magnetic ﬁeld. For

an insulating outer sphere, it yields ∂sωg|s=1 = 0 which corresponds to a stress-free boundary, as in the 1-D wave case with Pm → 0. In

the case Pm 1, the appropriate boundary condition is ωg|s=1 = 0 as the angular velocity of the outermost geostrophic cylinder is

C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

4.2 Normal modes

Assuming that ωg varies with time as eict, the eq. (8) can be transformed into a normal mode equation

−c2ωg(s) =

1 L

∂ ∂s

L

V˜A2

∂

ωg (s ∂s

)

.

(9)

Transmission and reﬂection of TAW on the geostrophic cylinder

tangent to the inner core set a special problem that we do not address

here. As an intermediate step, we simply illustrate our discussion

with results for the full sphere case, imposing ∂sωg|s=ε = 0, with ε 1 (we have checked the convergence of the numerical results

as ε → 0). It is of interest to write the solution of this equation in

the case c = 0 and V˜ A uniform

√

ωg (s )

=

1 2

α1

−

1 − s2 − log s2

1 − s2 + 1 + log(s) + α2. (10)

A non-zero solution (uniform rotation ωg(s) = α2) exists for the boundary condition ∂sωg|s=1 = 0 but not for the condition ωg|s=1 = 0 that applies when Pm 1. We are interested in this latter case,
despite its lack of geophysical realism, as contrasting the two boundary conditions sheds light on the nature of the constraint ∂sωg|s=1 = 0 that has always been used in TAW studies.
In the general case (c = 0, non-uniform V˜A), it remains easy to calculate numerically a solution of (9) for 0 < s < 1. We have
successfully checked our numerical results against the eigenvalues
listed in the table C1 of Roberts & Aurnou (2011), that have been obtained analytically for ∂sωg|s=1 = 0 and V˜A = 1. Then, the ﬁrst eigenvalues are (0, 5.28, 8.63, 11.87, 15.07, ..), whereas in the case V˜A = 1 and ωg|s=1 = 0 they are (2.94, 6.35, 9.58, 12.78, 15.95, ..). In the latter case, we recover our previous observation that 0 is not
an eigenvalue.
In contrast with an often-made statement (Buffett 1998; Jault
2003; Roberts & Aurnou 2011), the study of the case Pm 1 shows that it is not required to have ∂sωg|s=1 = 0 to obtain solutions with bounded values of ωg for s ≤ 1. On the other hand, the singularity of ∂sL at s = 1 implies a singularity of ∂sωg (which is O((1 − s)−1/2)

512 N. Schaeffer et al.
as s → 1) . That points to signiﬁcant viscous dissipation once the viscous term is reintroduced.
When Pm is neither very small nor very large, it is not possible to separate the interior region (where (9) applies) and the Hartmann boundary layer.
We can conclude the discussion of normal modes by noting that the solutions for the two cases Pm 1 and Pm 1 differ in a signiﬁcant way at the equator. In both cases, solutions are obtained which satisfy the appropriate boundary conditions and with bounded values of ωg for s in the interval [0, 1]. However, reintroducing dissipation modiﬁes the eigensolutions in the vicinity of the equator and the eigenvalues in the second case only.

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

4.3 Numerical experiments

To determine the reﬂection coefﬁcient of TAW at the equator of the

outer shell, we use a set-up that resembles the Earth’s core. The

code is the same as the one described in Section 3, but this time

with imposed global rotation. The total number of radial points is

typically 1200 and the maximum degree of Legendre polynomials

is set to 360.

For reﬂection to occur, there must be a non-zero imposed mag-

netic ﬁeld Bs at the equator. Hence, we set the simplest poten-

tial quadrupolar ﬁeld (generated from outside the sphere): Bs =

B0s, VA(s

)B=z =Bs2(Bs0)/z√anμd0

Bφ = ρ that

0. is

This large

ensures near the

a local travelling speed reﬂection point (s = 1).

The Lehnert number is small and always set to λ = V A/( ro) =

5 × 10−4, so that λPm is also small.

The initial velocity ﬁeld is along the azimutal direction φ and

depends only on the cylindrical radius s: uφ(s) = sωg(s) = u0s exp ( − (s − s0)2/ 2) with s0 = 0.675. We used two different width
= 0.02 and = 0.063. This initial velocity ﬁeld splits into a

TAW packet propagating inwards that we do not consider here, and

another travelling outwards that we carefully follow and we focus

on the reﬂection of this wave packet at the equator of the outer

shell (s = 1). The Lundquist number S based on the size of the

spherical shell ranges from 6 × 102 to 8 × 104 and the Ekman

number E = ν/ ro2 and magnetic Ekman number Em = η/ ro2 are both always very low and range from 5 × 10−10 to 5 × 10−7 over a wide range of magnetic Prandtl number: from Pm = 10−3 to

Pm = 102.

We measure the extremum of the velocity ﬁeld in the wave packet

before and after the reﬂection, ai and ar respectively, at a ﬁxed radius (s = 0.925 for = 0.02 and s = 0.75 for = 0.063), from

which we compute the corresponding reﬂection coefﬁcient R =

ar/ai, reported in Fig. 3 for an insulating outer shell. We found no signiﬁcant dependence with the Lundquist number S or the width

of the initial pulse (R varies by less than 0.03).

As expected from the discussion of Alfve´n waves equations, the

combination Pm = 1, no-slip boundary condition and insulating

wall corresponds to a special case whereby no reﬂection at all

occurs at the equator (see also Fig. 4g).

However, there are differences with the planar case. First, the

reﬂection coefﬁcient is not symmetric with respect to Pm = 1, as

expected from our discussion of torsional eigenmodes in spherical

geometry in the previous section. For large Pm there is high dissi-

pation and very little reﬂection compared to low Pm. Second, the

reﬂection coefﬁcient is not as large.

Space–time diagrams of the reﬂection of the wave at the equator

are presented in Fig. 4 for a few representative cases. The highest

reﬂection coefﬁcient occurs for the stress-free insulating case at

Figure 3. Reﬂection coefﬁcient for a TAW for insulating and no-slip boundary conditions, as a function of Pm. The Lundquist number is always large (S > 5000 for Pm ≥ 0.01 and S > 600 otherwise). For √reference√, the black curve is the planar Afve´n wave reﬂection coefﬁcient (1− Pm)/( Pm+1), and the red line marks the reﬂection coefﬁcient for a stress-free boundary (obtained with Pm = 1 but which is theoretically independent of Pm and corresponds to a no-slip boundary with Pm → 0).
Pm = 1: from R = 0.86 at S = 1000 to R = 0.88 at S = 1.5 × 104. In this case (Figs 4a and b) one can also see the ampliﬁcation of the velocity ﬁeld very near the boundary, as the magnetic ﬁeld must vanish, doing so by producing the reﬂected wave, just as in the planar case. This is not a boundary layer, but simply the superposition of the incident and reﬂected wave (see also Appendix A). The Hartmann boundary layer is too small to be seen on these plots, but we checked that its size and relative amplitude for velocity and magnetic ﬁelds do match the analytic theory developed in Appendix A.
For Pm = 0.1, the reﬂected wave carries only 16 per cent of the energy, the remaining being dissipated in the boundary layer. The magnetic ﬁeld changes sign at the reﬂection, whereas the velocity keeps the same sign (Figs 4c and d). For Pm = 10, the reﬂected energy drops to 3 per cent and the small reﬂected velocity ﬁeld has opposite sign, whereas the magnetic ﬁeld (barely visible on Fig. 4) keeps the same sign (Figs 4e and f). During its propagation, the incoming wave is also much more damped than for Pm = 0.1, even in the case where S or E have comparable values. This is due to strong dissipation at the top and bottom boundaries, which increases as the wave propagates toward the equator (visible in Fig. 4e) for Pm > 1. This may not be unrelated to the previously discussed singularity for normal modes in the case Pm > 1. A consequence of this large dissipation, is the difﬁculty to clearly identify the reﬂected wave, and to properly deﬁne a reﬂection coefﬁcient. The values reported in Fig. 3 are, thus, not very precise for Pm > 1.
It may also be worth noting that changing the magnetic boundary from insulating to a thin conducting shell allows weak reﬂection for Pm = 1 and no-slip velocity (Fig. 4h), in agreement with the analysis of the governing equations (Section 2.2).
4.4 Energy dissipation and normal modes
We want to emphasize that when no reﬂection occurs, the energy of the wave is dissipated very quickly. However, for liquid metals (Pm 1), only a small amount of the wave energy is absorbed in the event of a reﬂection, but many successive reﬂections can lead to signiﬁcant dissipation. Using the theoretical reﬂection coefﬁcient, we can estimate the timescale of dissipation of an Alfve´n wave due
C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

Reﬂection of Alfve´n waves 513

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

Figure 4. Space–time diagrams of the reﬂection of a TAW for S 104 and = 0.02 recorded in the equatorial plane, near the equator. Top row: stress-free
boundary with Pm = 1 (R = 0.88), (a) the azimuthal angular velocity uφ /s and (b) the azimuthal magnetic ﬁeld bφ (changing sign). Second row: No-slip boundary with Pm = 0.1 (R = 0.40), (c) the azimuthal angular velocity uφ /s and (d) the azimuthal magnetic ﬁeld bφ (changing sign). Third row: No-slip boundary with Pm = 10 (R = −0.17), (e) the azimuthal angular velocity uφ /s (changing sign) and (f) the azimuthal magnetic ﬁeld bφ . Bottom row: (g) azimuthal angular velocity uφ /s for no-slip boundary with Pm = 1 showing no reﬂection (R = 0) for insulating boundary, (h) and little reﬂection when the insulator is
replaced by a solid conductive layer.

to its reﬂections at the boundaries. In the case of an Alfve´n wave turbulence (many wave packets) in a spherical shell of radius L with homogeneous mean energy e, permeated by a magnetic ﬁeld of rms intensity B0, any wave packet will reach the outer insulating boundary once (on average) in the time interval L/V A. When it reﬂects on the boundary, it loses the fraction 1 − R2(Pm) of its energy, where R(Pm) is the reﬂection coefﬁcient (in amplitude). We can then estimate the dissipation rate of energy e due to this process

∂t e ∼ R2(Pm) − 1 L√Bμ0 0ρ e.

(11)

Hence, the timescale of dissipation at the boundaries

τs

=

L VA

1−

1

,

R 2 ( P m)

(12)

which is inversely proportional to the strength of the magnetic ﬁeld,

and depends on the diffusivities only through Pm.

We can compare this to the dissipation of Alfve´n waves of length scale in the bulk of the ﬂuid: τ v = 2 2/(η + ν). It appears that the length scale where surface and bulk dissipation are comparable is

such that

√

L/ = S 1 − R2.

(13)

C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

514 N. Schaeffer et al.

Hence, for the Earth’s core with S ∼ 104, and R ∼ 0.9 (the stress-free value which gives a good approximation of the low Pm value), the dissipation of Alfve´n waves is dominated by the partial absorption at the boundaries for length scales larger than L/45. For numerical simulations of the geodynamo with S ∼ 103 ans R ∼ 0.2, we have L/ ∼ 30.
These timescales are also relevant for torsional normal modes. In 1-D, normal modes are a superposition of waves propagating in opposite directions. Hence, if the dissipation of waves is dominated by their reﬂection, so will it be for the normal modes. From the previous estimation of L/ in the Earth’s core, we expect the dissipation of large wavelength TAW (the ones that can be observed) to be dominated by the effect of reﬂection. Furthermore, to detect a normal mode, its dissipation time must be much larger than its period T = 2π L(cV A)−1. The pulsation c of the ﬁrst torsional normal modes are given in section 4.2 in Alfve´n frequency units, and their dissipation time can be estimated by τ s for the large-scale normal modes. We deﬁne a quality factor for torsional normal modes by

Q

=

τs T

=

c 2π

1

1 −

R

2

.

(14)

Presence of normal modes requires Q 1. Assuming R = 0.9 (stress-free value) in the Earth’s core, we ﬁnd QE 0.8 c and for no-slip numerical simulations of the geodynamo we ﬁnd Qsim < 0.16 c. Considering the largest modes (with c 5 to 15), torsional oscillations could therefore persist in the Earth’s core for a few Alfve´n times, but are completely absent even from the best current geodynamo simulations.

observe eigenmodes, numerical models that use stress-free boundaries (e.g. Kuang & Bloxham 1999; Dumberry & Bloxham 2003; Busse & Simitev 2006; Sreenivasan & Jones 2011) are intrinsically much more suited for the study of torsional normal modes. Quasi-geostrophic dynamo models that can compute dynamo models at very low magnetic Prandtl numbers (Pm < 10−2 in Schaeffer & Cardin 2006), could also provide an interesting tool to study torsional oscillations.
In the case of the Earth’s core, a recent study (Gillet et al. 2010) found no clear evidence for reﬂection at the equator, although this has yet to be conﬁrmed. One might want to invoke turbulent viscosity (see the contrasted views of Deleplace & Cardin (2006) and Buffett & Christensen (2007) in a different context) to explain this fact, leading to an effective Pm close to 1 and inhibiting reﬂection of TAW. This would make numerical models more relevant, but is rather speculative. A solid conductive layer at the top of the core can also have a damping effect on the propagation and reﬂection of torsional waves, and we plan to investigate these matters in a forthcoming study.
ACKNOWLEDGMENTS
The numerical simulations were run at the Service Commun de Calcul Intensif de l’Observatoire de Grenoble (SCCI). We want to thank Mathieu Dumberry and an anonymous reviewer for their help in improving this paper, and Henri-Claude Nataf for useful comments.

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

5 DISCUSSION: IMPLICATION FOR NUMERICAL GEODYNAMO MODELS A N D T H E E A RT H - C O R E
We showed that numerical simulations conducted for Pm ∼ 1 cannot adequately reproduce the boundary conditions for TAW in the Earth’s core (where Pm 1). The small reﬂection coefﬁcient observed for TAW (Fig. 3) means that it is hard to observe TAW reﬂection at the equator in numerical simulations of the geodynamo which currently operate with 0.1 < Pm < 10 (e.g. Takahashi et al. 2008; Sakuraba & Roberts 2009), where the waves are moreover mixed with thermal convection.
As for possible torsional eigenmodes, it is almost impossible to observe them with such low reﬂection coefﬁcients. Unfortunately, that severely limits the ability of geodynamo simulations to exhibit torsional oscillation normal modes, because normal modes require a large reﬂection coefﬁcient to be observable: their period (of order L/V A) must be much larger than the energy dissipation time τ s (see expression 12). A few studies have tried to pin down torsional eigenmodes (Dumberry & Bloxham 2003; Sakuraba & Roberts 2008; Wicht & Christensen 2010) but even though they report waves propagating with the appropriate speed, they report neither reﬂection of these waves, nor eigenmodes.
Another issue for geodynamo models with very low diffusivities, is that the part of the energy carried by Alfve´n waves (regular or torsional) is dissipated very quickly (on an Alfve´n wave timescale), so that an Alfve´n wave turbulence would be damped much faster, and the turbulent state may be far from what we would expect in the Earth’s core.
Changing the boundary condition to stress-free simulates the case Pm = 0 with a high reﬂection coefﬁcient (R = 0.88), but still lower than the planar case. Even though this may still be problematic to

REFERENCES
Alboussie`re, T., Cardin, P., Debray, F., La Rizza, P., Masson, J.P., Plunian, F., Ribeiro, A. & Schmitt, D., 2011. Experimental evidence of Alfve´n wave propagation in a gallium alloy, Phys. Fluids, 23(9), 096 601-1–096601-10, doi:10.1063/1.3633090.
Alfve´n, H., 1942. Existence of electromagnetic-hydrodynamic waves, Nature, 150, 405–406.
Braginsky, S., 1984. Short-period geomagnetic secular variation, Geophys. Astrophys. Fluid Dyn., 30(1–2), 1–78.
Braginsky, S.I., 1970. Torsional magnetohydrodynamic vibrations in the Earth’s core and variations in day length, Geomagn. Aeron., 10, 1–8.
Buffett, B.A., 1998. Free oscillations in the length of day: inferences on physical properties near the core mantle boundary, in The Core–Mantle Boundary Region, Geodynamics Series, Vol. 28, pp. 153–165, eds Gurnis, M., Wysession, M.E., Knittle, E. & Buffett, B.A., AGU, Washington, D.C.
Buffett, B.A. & Christensen, U.R., 2007. Magnetic and viscous coupling at the core-mantle boundary: inferences from observations of the Earth’s nutations, Geophys. J. Int., 171, 145–152.
Busse, F.H. & Simitev, R.D., 2005. Convection in rotating spherical ﬂuid shells and its dynamo states, in Fluid Dynamics and Dynamos in Astrophysics and Geophysics, pp. 359–392, eds Jones, C.A., Soward, A.M., Hughes, D.W. & Weiss, N.O., CRC Press, Boca Raton.
Busse, F.H. & Simitev, R.D., 2006. Parameter dependences of convectiondriven dynamos in rotating spherical ﬂuid shells, Geophys. Astrophys. Fluid Dyn., 100(4), 341–361.
Deleplace, B. & Cardin, P., 2006. Viscomagnetic torque at the core mantle boundary, Geophys. J. Int., 167, 557–566.
Dumberry, M. & Bloxham, J., 2003. Torque balance, Taylor’s constraint and torsional oscillations in a numerical model of the geodynamo, Phys. Earth planet. Inter., 140(1–3), 29–51.
Gillet, N., Jault, D., Canet, E. & Fournier, A., 2010. Fast torsional waves and strong magnetic ﬁeld within the Earth’s core, Nature, 465(7294), 74–77.
C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

Gillet, N., Schaeffer, N. & Jault, D., 2011. Rationale and geophysical evidence for quasi-geostrophic rapid dynamics within the Earth’s outer core, Phys. Earth planet. Inter., 187(3–4), 380–390.
Jameson, A., 1961. Magnetohydrodynamic waves, PhD thesis, University of Cambridge, Cambridge.
Jameson, A., 1964. A demonstration of Alfve´n waves part 1. generation of standing waves, J. Fluid Mech., 19(04), 513–527.
Jault, D., 2003. Electromagnetic and topographic coupling, and LOD variations., in Earth’s Core and Lower Mantle, pp. 56–76, eds Jones, C.A., Soward, A.M. & Zhang, K., Taylor & Francis, Oxford.
Jault, D., 2008. Axial invariance of rapidly varying diffusionless motions in the Earth’s core interior, Phys. Earth planet. Inter., 166, 67–76.
Jordi, C., Morrison, L., Rosen, R., Salstein, D. & Rossello, G., 1994. Fluctuations in the earth’s rotation since 1830 from high-resolution astronomical data, Geophys. J. Int., 117(3), 811–818.
Kuang, W. & Bloxham, J., 1999. Numerical modeling of magnetohydrodynamic convection in a rapidly rotating spherical shell: weak and strong ﬁeld dynamo action, J. Comput. Phys., 153(1), 51–81.
Roberts, P.H., 1967. An Introduction to Magnetohydrodynamics, Elsevier, New York, NY.
Roberts, P.H. & Aurnou, J.M., 2011. On the theory of core-mantle coupling, Geophys. Astrophys. Fluid Dyn., 106, 1–74.
Sakuraba, A. & Roberts, P., 2008. Torsional oscillations in a numerical geodynamo operating in a regime of low ekman and magnetic prandtl numbers, in AGU Fall Meeting Abstracts, pp. A1773+, AGU, Washington, D.C.
Sakuraba, A. & Roberts, P.H., 2009. Generation of a strong magnetic ﬁeld using uniform heat ﬂux at the surface of the core, Nat. Geosci., 2(11), 802–805.
Schaeffer, N., 2012. Efﬁcient Spherical Harmonic Transforms aimed at pseudo spectral numerical simulations, preprint, arXiv:1202.6522.
Schaeffer, N. & Cardin, P., 2006. Quasi-geostrophic kinematic dynamos at low magnetic Prandtl number, Earth planet. Sci. Lett., 245(3–4), 595– 604.
Sreenivasan, B. & Jones, C.A., 2011. Helicity generation and subcritical behaviour in rapidly rotating dynamos, J. Fluid Mech., 688, 5–30.
Takahashi, F., Matsushima, M. & Honkura, Y., 2008. Scale variability in convection-driven MHD dynamos at low ekman number, Phys. Earth planet. Inter., 167(3–4), 168–178.
Wicht, J. & Christensen, U.R., 2010. Torsional oscillations in dynamo simulations, Geophys. J. Int., 181(3), 1367–1380.

A P P E N D I X A : A N A LY T I C A L F V E´ N WAV E S O LU T I O N S I N O N E D I M E N S I O N

A1 Plane wave solutions

Following Jameson (1961, p. 15–18), we look for plane wave

sbo=lut√ioμns0ρoBf

eqs (1) ei (ωt +k x )

and

(2),

substituting

u

=

Uei(ωt +kx)

and

iω + νk2 U = VAik B,

(A1)

iω + ηk2 B = VAikU,

(A2)

which we can combine into

νη k4 + VA2 + i ω(η + ν) k2 − ω2 = 0,

(A3)

for which the exact solutions are

k2

=

−

VA2 2νη

(1

+

2i

)

1±

1

+

4ω2νη VA4(1 + 2i

)2

,

(A4)

where is the reciprocal Lundquist number based on the frequency

= ω(η + ν) . 2Va2

(A5)

C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

Reﬂection of Alfve´n waves 515

In the and also

rωeg√imνηe/wVhA2ere

Alfve´n waves do propagate, we 1 so we can approximate the

have square

1 root

by its ﬁrst-order Taylor expansion, which leads to two solutions

k12 and k22

k12 =

ω2 VA2 (1 + 2i

)−1

k22

=

−

VA2 νη

(1

+

2i

).

(A6)

The solutions k = ±k1 = ±ω/V A(1 − i ), correspond to the prop-

agation in both directions of an Alfve´n wave at the speed V A and

with attenuation on a length scale V A/( ω). The solutions k = ±k2

δ

≡±√i/δν

correspond η/ VA.

to

a

Hartmann

boundary

layer

of

thickness

Finally, from eqs (A1) and (A2) we know that U and B are related

for each k by

B U

=

i kVa iω + ηk2

=

iω + νk2 i kVa

≡ αk ,

(A7)

and for the solutions k = ±k1 and k = ±k2, it reduces to

α±k1 ±1

ν√

α±k2

±

=± η

P m.

(A8)

This means that for the travelling wave solution, U and B have always the same amplitude and the same phase when propagating in the direction opposite to the imposed magnetic ﬁeld, or opposite phase when propagating in the same direction. For the boundary layers, in the limit Pm 1 they involve the velocity ﬁeld alone, whereas for Pm 1 they involve only the magnetic ﬁeld.

A2 Reﬂection coefﬁcient at an insulating wall

To derive the reﬂection coefﬁcient, we consider an insulating wall at x = 0 with an incoming Afve´n wave from the x > 0 region (k = +k1), giving rise to a reﬂected wave (k = −k1). The boundary conditions are matched by a boundary layer (k = +k2) localized near x = 0 (the solution k = −k2 is growing exponentially for x > 0 and has to be rejected for this problem). The solution to this problem
reads

u = eiωt eik1x + Re−ik1x + βeik2x ,

(A9)

b = eiωt αk1 eik1x − Re−ik1x + αk2 βeik2x √μ0ρ,

(A10)

where we have taken into account the fact that α−k1 = −αk1 (see eq. A8).
The boundary conditions u = 0 and b = 0 at x = 0 lead to

1+R+β =0

αk1 (1 − R) + αk2 β = 0,

from which we ﬁnd the amplitude β of the velocity boundary layer contribution, and the reﬂection coefﬁcient R of the amplitude of the velocity component

−2 β=
1 + αk2 /αk1

R = 1 − αk2 /αk1 . 1 + αk2 /αk1

We are √left to evaluate αk2 /αk1 using eqs (A8), which gives

αk2 /αk1 = ν/η at leading order in , and thus

√

1 − Pm

R= √ ,

(A11)

1 + Pm

which is independent of ω and V A. In the case Pm = 1, we then have R = 0 and β = −1 which means that no reﬂection occurs and that the amplitude of the incoming wave is canceled by the boundary layer alone.

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

516 N. Schaeffer et al.
It may be worth emphasizing that, although the boundary layer has the same thickness δ in the velocity and the magnetic ﬁeld components, in the limit Pm → 0, we have β → −2 and αk2 β → 0, so that the boundary layer is apparent only in the velocity ﬁeld component (eq. A9), whereas in the limit Pm → ∞, we have β →

0 and αk2 β → −2, so that the boundary layer is apparent only in the magnetic ﬁeld component (eq. A10).
Finally, we remark that if one sets ν = 0 or η = 0 from the beginning in eqs (A1) and (A2), the solution corresponding to the
boundary layer does not exist anymore.

Downloaded from https://academic.oup.com/gji/article/191/2/508/641474 by guest on 16 January 2024

C 2012 The Authors, GJI, 191, 508–516 Geophysical Journal International C 2012 RAS

