## Abstract

We investigate numerically the flow of an electrically conducting fluid confined in a spherical shell, with the outer sphere stationary, the inner sphere rotating and a magnetic field imposed parallel to the rotation axis. We compute both the axisymmetric basic states and their non-axisymmetric instabilities. Two distinct instability classes emerge, one connected to previous non-magnetic results, the other to previous strongly magnetic results. Both instabilities arise from the basic state's meridional circulation, but are otherwise very different from one another, and are separated by a region of stability that persists even for large Reynolds numbers. Finally, we compute the fully three-dimensional nonlinear equilibration of both instabilities. The second class exhibits a rich variety of secondary bifurcations, involving mode transitions between different azimuthal wave numbers.

## 1. Introduction

Spherical Couette flow is the flow induced in a fluid-filled spherical shell by fixing the outer sphere and rotating the inner one. The problem is defined by only two parameters, the radius ratio of the two spheres and the Reynolds number measuring the inner sphere's rotation rate. Despite this seeming simplicity, spherical Couette flow can generate an enormous variety of flow patterns and transitions, the study of which is an important part of classical fluid dynamics and hydrodynamic stability theory.

Now suppose the fluid is electrically conducting, and a magnetic field is externally imposed, thereby introducing fundamentally new, magnetohydrodynamic effects. Specifically, we will impose a uniform field parallel to the inner sphere's axis of rotation. For sufficiently strong fields, a shear layer emerges on the so-called tangent cylinder , the cylinder parallel to the field and just touching the inner sphere. The thickness of this layer scales as *M*^{−1/2}, where the Hartmann number *M* is a measure of the field strength (Starchenko 1998).

This problem was first studied numerically by Hollerbach (2000*a*), who obtained results in broad agreement with Starchenko's asymptotic scalings, but also showed that the electromagnetic boundary conditions imposed on the inner and/or outer spheres play a crucial role. Only insulating boundaries yield a shear layer on ; conducting boundaries yield a very strong counter-rotating jet. Hollerbach & Skinner (2001) then investigated the non-axisymmetric instabilities of these shear layers and jets, and showed them to be essentially Kelvin–Helmholtz-type instabilities, in which the initially axisymmetric layers and jets roll up into a series of vortices.

The work presented here is a sequel to Hollerbach and Skinner, hereafter referred to as paper I, but now focusing on the weakly magnetic regime, as opposed to the strongly magnetic regime considered in I. There are two reasons for reconsidering this problem. First, Hollerbach *et al*. (2006), hereafter referred to as paper II, considered non-axisymmetric instabilities in non-magnetic spherical Couette flow, and obtained results completely different from the strongly magnetic results in I. It is therefore of interest to map out the transition from one regime to the other, and indeed we will find some quite unexpected behaviour in the intermediate regime.

The second reason for returning to this problem is the work of Sisan *et al*. (2004), who did liquid metal experiments in exactly this configuration, and obtained a very intriguing series of non-axisymmetric structures. We suggest that their results may perhaps be turbulent analogues of the non-axisymmetric instabilities computed here.

Finally, we note that, instead of a uniform axial field, one could also impose a dipolar magnetic field. This problem was considered numerically by Hollerbach (1994) and Dormy *et al*. (1998), asymptotically by Kleeorin *et al*. (1997) and experimentally by Nataf *et al*. (2006, 2008). Some of the resulting instabilities were also calculated by Hollerbach *et al*. (2007). These papers show that the dipolar problem is different in many ways from the axial problem, so we will not consider dipolar fields here.

## 2. Equations

As derived in I, the non-dimensional equations to be solved are(2.1)(2.2)where ** U** is the fluid flow and

**is the induced magnetic field. The Reynolds number is a measure of the inner sphere's rotation rate**

*b**Ω*, and the Hartmann number

*M*=

*B*

_{0}

*r*

_{i}/(

*μρνη*)

^{1/2}is a measure of the imposed magnetic field . The boundary conditions associated with equation (2.1) are(2.3)where the inner and outer spheres' radii are chosen as

*r*

_{i}=1,

*r*

_{o}=3.

For equation (2.2), we will consider insulating boundary conditions only. As noted above, in the large Hartmann number limit, conducting boundaries yield a counter-rotating jet rather than a shear layer on . The gradual emergence of this counter-rotating jet in the basic state is likely to make the resulting transitions between different instabilities even more complicated than the insulating results presented here.

As in I, we solve these equations in three different variants. We begin by computing the axisymmetric basic states. Next we linearize equations (2.1) and (2.2) about the given basic state, thereby decoupling the different azimuthal wave numbers *m*. Because all the basic states turn out to be equatorially symmetric (perturbations of the opposite symmetry were introduced, but always decayed away), each wave number *m* decouples further into two distinct equatorial symmetry classes, satisfying either(2.4)which we call equatorially symmetric, or(2.5)which we call equatorially antisymmetric. Finally, having solved for the linear onset of the different instability modes, we study their nonlinear equilibration in the supercritical regime. All calculations were done using the code described by Hollerbach (2000*b*), in which the radial structure is expanded in Chebyshev polynomials, the angular structure in spherical harmonics and the time integration is a modified second-order Runge–Kutta method.

The highest resolutions used were 120 modes in *r* and spherical harmonics up to *l*=240 for the two-dimensional basic state and linear onset calculations, and 60 modes in *r* and spherical harmonics up to *l*=120, *m*=30 for the fully three-dimensional calculations. These resolutions were sufficient to compute the linear onset up to *Re*=3000, and the nonlinear equilibration up to *Re*=1200, and were carefully checked to ensure that all solutions were fully resolved.

## 3. The basic state

Figure 1 shows examples of the basic state, at a fixed *Re*=250, and Hartmann number increasing from 1 to 64. For *M*=1, the flow is still much the same as in the non-magnetic problem considered in II. The primary flow is the angular velocity driven directly by the rotation of the inner sphere. As the fluid rotates with the sphere, inertia tends to fling it outward, thereby driving a secondary meridional circulation consisting of a narrow jet right on the equator and a return flow in the rest of the shell. Over most of the inner sphere, the meridional circulation is therefore directed inward. It is this feature that explains the boundary-layer structure observed in the angular velocity, with advection inward just balancing diffusion outward (Howarth 1951). Only at the equator are the angular velocity contours also drawn outward by the equatorial jet.

As the Hartmann number is increased, magnetohydrodynamic effects become more and more important. The angular velocity becomes increasingly concentrated into a narrow shear layer on the tangent cylinder , together with Hartmann layers on the inner and outer boundaries. This alignment in the axial direction is caused by the magnetic tension along the field lines. Fluid inside is coupled to both boundaries, and hence rotates at a rate intermediate between 0 at the outer boundary and 1 at the inner, whereas fluid outside is coupled only to the outer boundary, and hence remains at rest. Note also how the meridional circulation is very strongly suppressed as *M* increases.

## 4. Linear onset of instability

In I, we considered this problem for Hartmann numbers from 32 to 1000. The instability curves from fig. 4*a* of I are shown here as the right half of figure 2. We see how the critical wave number increases from *m*=2 to 5; the corresponding critical Reynolds numbers also increase, roughly as *M*^{0.66} according to I.

The left half of figure 2 shows the new results obtained by considering the weakly magnetic regime *M*=1–32. We first note how the previous instability curves turn dramatically upward again, with an increasingly large range of wave numbers also becoming unstable almost simultaneously.

Next, coming from the left, at *M*=1, *m*=2 and 3 become unstable at almost the same critical Reynolds numbers. These are precisely the non-magnetic instabilities previously considered in II. (That *m*=2 and 3 occur at virtually the same *Re* is an accident of this particular radius ratio; as shown in II, for slightly wider shells *m*=2 is preferred, for narrower shells *m*=3 and higher.) If we now increase *M*, these instability curves also turn sharply upward.

Separating the two different instability classes is a quite extraordinary strip of stability; for example, *M*=20, *Re*=3000 is stable, whereas *M*=20, *Re*=300 is unstable. It would be very interesting to know whether this region of stability extends to indefinitely large Reynolds numbers, and, if so, whether it continues to follow the *Re*∼*M*^{2} scaling suggested by the results here. Unfortunately, numerical limitations make the large *Re* regime increasingly difficult to explore.

One further aspect of these two instability classes worth mentioning concerns their equatorial symmetry: the modes on the left are antisymmetric (equation (2.5)) and the modes on the right are symmetric (equation (2.4)). In the remainder of this paper, we will study these two instability classes in more detail, relating the left modes to the previous non-magnetic results in II, and the left half of the right modes to the previous strongly magnetic results in I.

## 5. Antisymmetric modes

Figure 3 shows the solutions at the four points on the *m*=2 curve indicated in figure 2. As before in figure 1, the top row shows the angular velocity and the bottom row shows the meridional circulation. The superimposed grey shading (the same in both rows) is the azimuthally integrated kinetic energy of the instability, that is, the quantity , as a function of *r* and *θ*. Based on where the instability is concentrated, it is clearly associated with the radially outflowing jet in the meridional circulation. As demonstrated in II, the essence of this instability is that the radial jet adopts a wavy structure, alternately slightly above and below the equatorial plane.

As shown in figure 2, once *M* is increased beyond approximately 5, the critical Reynolds number starts to increase very rapidly. Given how effectively the magnetic field suppresses the meridional circulation, it is hardly surprising that it would also suppress this instability of it. There is no fundamental, qualitative change in the nature of the instability though as one moves along the instability curve into the increasingly strong field regime.

Figure 4 shows the nonlinear equilibration of this instability, which is indeed seen to consist of a slightly wavy structure in the radial jet, with an *m*=2 periodicity in *ϕ*. Comparing this figure with fig. 16 in II, we see that the flow is again very similar to the non-magnetic case.

## 6. Symmetric modes

Figure 5 shows the solutions at the four points on the symmetric curves indicated in figure 2. The instabilities are again associated with the meridional circulation, but not so much with the radial jet as with its outer edge, where it turns away from the equatorial plane to begin its recirculation in the rest of the volume. These instabilities are thus quite different not only from those considered in §5, but also from the Kelvin–Helmholtz-type shear layer instabilities considered in I, to which they are continuously connected though. On the left half of these symmetric curves, one obtains the meridional circulation instabilities here; on the right half, one obtains the shear layer instabilities in I.

To study the nonlinear equilibration of these modes, we fix *Re*=1200 and increase *M*. Figure 6 illustrates the various bifurcations that occur. The first wave number to become unstable is *m*=5, at *M*=17.2. In the supercritical regime, the fully three-dimensional solution then contains only multiples of this basic wave number. This mode remains stable until *M*=20.4, at which point it becomes unstable to perturbations of other wave numbers. Perturbations having the opposite equatorial symmetry were also introduced, but all solutions throughout this entire section always remained equatorially symmetric.

After a brief transient evolution, a new solution emerges that contains only multiples of *m*=4. Increasing *M* further, this mode in turn becomes unstable at *M*=21.8, and switches to a solution containing only multiples of *m*=3. Alternatively, if one decreases *M* at any point in this progression, both the *m*=4 and 3 modes eventually switch back to the *m*=5 mode, as shown in figure 6.

Figure 7 shows these *m*=5, 4 and 3 solutions. Figure 7*a* shows the stream function of the vertically integrated horizontal flow, as in I. Figure 7*b* shows *U*_{z} in the plane *z*=1, the plane just touching the inner sphere's pole. The small dashed circle in the centre of each panel indicates the location of the tangent cylinder. We see that the instabilities are indeed concentrated well outside the tangent cylinder, and consist primarily of a fluctuation in *U*_{z}, corresponding precisely to the grey shading in figure 5.

Returning to figure 6, the *m*=3 mode is stable up to *M*=24.1, at which point it switches to a mode containing only multiples of *m*=2. Figure 8 shows these solutions, at *M*=25, 35 and 50. Note, in particular, how the instability is gradually forced inward, until by *M*=50 it is concentrated virtually on the tangent cylinder. These are now precisely the shear layer instabilities considered in I, indicating that, in the nonlinear regime as well, there is a smooth progression from the meridional circulation instabilities for small *M* to the shear layer instabilities for large *M*.

Finally, if we are on the *m*=2 branch and decrease *M* again, at *M*=23.0 we find a Hopf bifurcation to an oscillatory solution. That is, all of these solutions so far have been oscillatory in the sense that they drift in longitude, but in a suitably co-rotating frame they appear steady. This new solution is genuinely time dependent. Figure 9*a* shows how the kinetic energy in wave numbers *m*=2 and 4 fluctuates periodically in time. Figure 10 shows these solutions, at the six times indicated by the dots in figure 9. The oscillation is seen to consist of a fluctuation in the relative contributions from *m*=2 and 4. Note also how the last two panels are indeed exactly the same as the first two, merely shifted in longitude, indicating the periodic nature of these solutions.

Returning to figure 9, figure 9*b* shows how the energy in *m*=1 and 3 perturbations decays in time, indicating that this solution is stable to odd wave number perturbations. If *M* is further reduced, we find that at *M*=22.4, these perturbations start to grow rather than decay. In this way, the stability boundary of these oscillatory solutions can be determined just as accurately as any of the other bifurcation points. For *M*<22.4, one switches back to the *m*=3 mode again.

## 7. Conclusion

In the process of connecting the previous non-magnetic results in II with the strongly magnetic results in I, we have discovered that the entire sequence of gradually increasing *M* for fixed *Re* involves quite a number of different regimes: first comes the basically non-magnetic radial jet instability of §5, then a narrow region of stability, then the return flow instabilities of §6. These gradually merge into the shear layer instabilities of I, until for even larger Hartmann numbers the basic state would become stable again.

Returning to the experiments of Sisan *et al*. (2004), they performed exactly this sequence of fixing *Re* and increasing *M*, and also obtained a number of different non-axisymmetric solutions and transitions between different azimuthal wave numbers and equatorial symmetries. However, their Reynolds numbers were *O*(10^{6}), and correspondingly even the ‘basic states’ from which their various solutions arose were already fully turbulent. It is unfortunately not possible therefore to make a direct comparison with their experiments and the work presented here.

Nevertheless, one might conjecture that at least some of their results may be turbulent analogues of the solutions presented here, just as in classical Taylor–Couette flow, for example, it is possible to have turbulent Taylor vortices that are otherwise similar to laminar Taylor vortices at much lower Reynolds numbers (e.g. Koschmieder 1979; Takeda 1999). It would certainly be well worth doing further work in magnetohydrodynamic spherical Couette flow, both computational and experimental, to bridge the gap in parameter space that currently exists between the two.

## Footnotes

- Received January 3, 2009.
- Accepted March 2, 2009.

- © 2009 The Royal Society