2 Mathematical formulation

Consider a Boussinesq fluid sphere of radius $r_{o}$ with constant thermal diffusivity $\unicode[STIX]{x1D705}$ , thermal expansion coefficient $\unicode[STIX]{x1D6FC}$ and kinematic viscosity $\unicode[STIX]{x1D708}$ . The fluid sphere rotates uniformly with constant angular velocity $\hat{\boldsymbol{z}}\unicode[STIX]{x1D6FA}$ in the presence of its own gravitational field $-\unicode[STIX]{x1D6FE}\boldsymbol{r}$ , where $\unicode[STIX]{x1D6FE}$ is a positive constant and $\boldsymbol{r}$ is the position vector. The whole sphere is heated by a uniform distribution of heat sources (see, for example, Chandrasekhar Reference Chandrasekhar1961; Roberts Reference Roberts1968; Busse Reference Busse1970; Soward Reference Soward1977; Jones et al. Reference Jones, Soward and Mussa2000; Sanchez et al. Reference Sanchez, Garcia and Net2016), producing the unstable conducting temperature gradient $-\unicode[STIX]{x1D6FD}\boldsymbol{r}$ , $\unicode[STIX]{x1D6FD}$ being a positive constant. The problem of convection is governed by the dimensionless equations:

where $\unicode[STIX]{x1D6E9}$ represents the deviation of the temperature from its static distribution, $p$ is the total pressure, and $\boldsymbol{u}$ is the velocity of convection in spherical polar coordinates $(r,\unicode[STIX]{x1D703},\unicode[STIX]{x1D719})$ with unit vectors ( $\hat{\boldsymbol{r}},\hat{\unicode[STIX]{x1D73D}},\hat{\unicode[STIX]{x1D753}}$ ). In deriving (2.1)–(2.3), we have employed $r_{o}$ as the length scale, $1/\unicode[STIX]{x1D6FA}$ as the unit of time and $\unicode[STIX]{x1D6FD}r_{o}^{4}\unicode[STIX]{x1D6FA}/\unicode[STIX]{x1D705}$ as the unit of temperature fluctuation $\unicode[STIX]{x1D6E9}$ . The three non-dimensional parameters, the Rayleigh number $R$ , the Prandtl number $Pr$ and the Ekman number $E$ , are defined as

The amplitude of convection in this study is assumed to be sufficiently small in order to neglect the nonlinear terms in (2.1) and (2.2).

Two different sets of boundary condition are considered in the asymptotic analysis. The first is perfectly conducting (isothermal), impenetrable and stress-free conditions given by

and the second is no-slip and isothermal conditions described by

We first derive an asymptotic solution of (2.1)–(2.3) at $0<E\ll 1$ describing torsional convection for stress-free conditions (2.5), and then extend the asymptotic analysis to that for no-slip conditions (2.6).

3 Asymptotic solution for the stress-free condition

Our asymptotic analysis is based on an assumption that torsional convection in rapidly rotating spheres is dominated by spherical inertial oscillation while buoyancy forces appear only at the next order to drive it against viscous damping. This assumption leads to an asymptotic expansion for $0<E\ll 1$ in the form

where $\text{i}=\sqrt{-1}$ , $\unicode[STIX]{x1D70E}$ denotes the half-frequency (the frequency $\unicode[STIX]{x1D714}=2\unicode[STIX]{x1D70E}$ ) of torsional convection with $0<|\unicode[STIX]{x1D70E}|<1$ , $\widehat{\boldsymbol{u}}$ , $\widehat{p}$ and $\unicode[STIX]{x1D70E}_{1}$ represent small perturbations, caused by viscous effects, to the leading-order solution $\boldsymbol{u}_{0},p_{0}$ and $\unicode[STIX]{x1D70E}_{0}$ , and satisfy $|\widehat{\boldsymbol{u}}|\ll |\boldsymbol{u}_{0}|$ and $|\unicode[STIX]{x1D70E}_{1}|\ll |\unicode[STIX]{x1D70E}_{0}|$ , and $\widetilde{\boldsymbol{u}}$ denotes a weak viscous boundary flow on the bounding spherical surface due to stress-free condition (2.5).

Substitution of the asymptotic expansion (3.1)–(3.4) into (2.1)–(2.3) yields the leading-order problem that describes non-dissipative thermal-inertial oscillation and is governed by the equations

subject to the boundary conditions

It can be shown, after some manipulations, that solutions of (3.5)–(3.7) are expressible as

where $j_{l}$ denotes the spherical Bessel function of the first kind, $P_{l}$ is the Legendre function, $\unicode[STIX]{x1D6FD}_{lq}$ with $q=1,2,3,\ldots$ are solutions of $j_{l}(\unicode[STIX]{x1D6FD}_{lq})=0$ and ordered such that $0<\unicode[STIX]{x1D6FD}_{l1}<\unicode[STIX]{x1D6FD}_{l2}<\unicode[STIX]{x1D6FD}_{l3}<\cdots \,$ , $k=1,2,3,\ldots$ can be regarded as a parameter of the solution which is to be determined by the next-order problem, ${\mathcal{C}}_{kij}$ is defined as

and $\unicode[STIX]{x1D70E}_{0}$ is a solution of

which has $2k$ real distinct solutions within $0<|\unicode[STIX]{x1D70E}_{0}|<1$ and is also to be determined by the next-order problem. The method for deriving ${\mathcal{C}}_{kij}$ and $\unicode[STIX]{x1D70E}_{0}$ can be found in Zhang et al. (Reference Zhang, Earnshaw, Liao and Busse2001).

For a given $E$ and $Pr$ , evidently, the solution $(p_{0},\boldsymbol{u}_{0},\unicode[STIX]{x1D6E9}_{0})$ associated with the unknown values of $k$ and $\unicode[STIX]{x1D70E}_{0}$ is not determined by the leading-order problem. Its determination is through the solvability condition of the next-order problem described by

subject to the stress-free boundary condition

where we have neglected higher-order terms such as $E\unicode[STIX]{x1D6FB}^{2}\widehat{\boldsymbol{u}}$ since $0<E\ll 1$ and $|\widehat{\boldsymbol{u}}|\ll |\boldsymbol{u}_{0}|$ , but the boundary-layer flow $\widetilde{\boldsymbol{u}}$ , albeit weak, must be retained such that $(\boldsymbol{u}_{0}+\widetilde{\boldsymbol{u}})$ satisfies the stress-free boundary condition. At this order, thermal effects are coupled with the non-dissipative thermal-inertial oscillation, driving torsional convection against viscous dissipation. Since (3.16), whose right-hand side is related to the leading-order solution $(\boldsymbol{u}_{0},\unicode[STIX]{x1D6E9}_{0})$ , is inhomogeneous, it requires a solvability condition whose real part, after making use of the stress-free condition, yields an expression for the Rayleigh number $R$ :

where $f^{\ast }$ denotes the complex conjugate of $f$ and

The physically preferred solution $(p_{0},\boldsymbol{u}_{0},\unicode[STIX]{x1D6E9}_{0})$ corresponds to the smallest Rayleigh number $R$ , denoted by the critical Rayleigh number $R_{c}$ , at which convective instabilities can occur, which is determined by minimizing $R$ in (3.19) over a manifold of $\boldsymbol{u}_{0}$ and $\unicode[STIX]{x1D6E9}_{0}$ with different values of $k$ and $\unicode[STIX]{x1D70E}_{0}$ .

Figure 1. The Rayleigh number $R$ plotted as a function of the Prandtl number $Pr$ at $E=10^{-4}$ : (a) for the stress-free condition and (b) for the no-slip condition. The curve labelled with A represents solutions for retrogradely travelling waves with the azimuthal wavenumber $m=1$ ; the curve B for axisymmetric torsional oscillation; the curve C for progradely travelling waves with $m=1$ ; and the curve D for retrogradely travelling waves with $m=2$ . The mode D is not preferred with the stress-free condition and, hence, is not shown in (a).

After determining the critical Rayleigh number $R_{c}$ together with the values of $k$ and $\unicode[STIX]{x1D70E}_{0}$ using (3.19), the imaginary part of the solvability condition is then used to calculate the critical half-frequency $\unicode[STIX]{x1D70E}_{c}$ via the expression

where

which is always negative, suggesting that viscous effects always reduce the size of $\unicode[STIX]{x1D70E}_{0}$ . Evidently, the asymptotic solution (3.19)–(3.22) is valid only for $0<E\ll 1$ with sufficiently small $Pr$ .

Figure 2. (a,b) Contours of the radial flow $\hat{\boldsymbol{r}}\boldsymbol{\cdot }\boldsymbol{u}$ and (c,d) contours of the azimuthal flow $\hat{\unicode[STIX]{x1D753}}\boldsymbol{\cdot }\boldsymbol{u}$ in a meridional plane at two different oscillation states of torsional convection in rapidly rotating spheres. Solid contours indicate $\hat{\boldsymbol{r}}\boldsymbol{\cdot }\boldsymbol{u}>0$ (or $\hat{\unicode[STIX]{x1D753}}\boldsymbol{\cdot }\boldsymbol{u}>0$ ) while dashed contours correspond to $\hat{\boldsymbol{r}}\boldsymbol{\cdot }\boldsymbol{u}<0$ (or $\hat{\unicode[STIX]{x1D753}}\boldsymbol{\cdot }\boldsymbol{u}<0$ ).

The result of the asymptotic solution obtained from (3.19)–(3.22) for torsional convection is fully consistent with that of the recent numerical simulation (Sanchez et al. Reference Sanchez, Garcia and Net2016). Figure 1(a) presents some typical asymptotic solutions obtained at $E=10^{-4}$ for different values of $Pr$ . Of three curves in figure 1(a), the curve labelled with A represents an asymptotic solution for non-axisymmetric, equatorially symmetric, retrogradely travelling waves with the azimuthal wavenumber $m=1$ , the curve C for progradely travelling waves with $m=1$ , and both the curves A and C are computed according to the asymptotic formula given by Zhang (Reference Zhang1994) for non-axisymmetric travelling waves. The curve B, calculated using (3.19) with $k=1$ and $\unicode[STIX]{x1D70E}_{0}=1/\sqrt{5}$ , represents the new branch of torsional convection that is physically preferred in the range $7.0\times 10^{-4}\leqslant Pr\leqslant 1.35\times 10^{-3}$ . The structure of torsional convection is displayed in figure 2 at two different states of oscillation. Asymptotic relations (3.19)–(3.21) at different values of $E$ also show largely similar behaviours, reaffirming the numerical finding (Sanchez et al. Reference Sanchez, Garcia and Net2016) that the most unstable mode of convection for $Pr/E\approx 10$ and $0<E\ll 1$ is in the form of torsional oscillation. Moreover, the leading-order asymptotic solution for the pressure $p$ and the velocity $\boldsymbol{u}$ of torsional convection is given by

which perhaps represents the simplest analytical solution of thermal convection in rapidly rotating fluid spheres with $0<E\ll 1$ . It is remarkable that a physically realizable solution of rotating convection in a spherical geometry can be analytically so simple.

It should be noticed that the asymptotic solution given by (3.19)–(3.22) is also reasonably accurate. As a typical example for $E=10^{-4}$ and $Pr=10^{-3}$ , asymptotic relations (3.19)–(3.22) give rise to $R_{c}=0.782$ and $\unicode[STIX]{x1D70E}_{c}=0.447$ while our fully numerical simulation, using a similar spectral method discussed by Zhang (Reference Zhang1992), produces $R_{c}=0.760$ and $\unicode[STIX]{x1D70E}_{c}=0.447$ . Even for a moderately small Ekman number $E=10^{-3}$ with $Pr=0.01$ , (3.19)–(3.22) still yield the most unstable mode of torsional convection with $R_{c}=7.82$ and $\unicode[STIX]{x1D70E}_{c}=0.444$ while our numerical analysis gives $R_{c}=7.33$ and $\unicode[STIX]{x1D70E}_{c}=0.444$ . They demonstrate that a good quantitative agreement is achieved between the numerical simulation and the analytical solution. In short, expressions (3.19)–(3.22) together with (3.23)–(3.26) constitute the complete asymptotic solution of torsional convection in rapidly rotating fluid spheres with the stress-free condition for $0<E\ll 1$ .

4 Asymptotic analysis for the no-slip condition

Whether torsional convection is also physically preferred in rapidly rotating spheres with the no-slip boundary condition can be answered by performing a similar, but more complicated, asymptotic analysis. The asymptotic expansion at $0<E\ll 1$ for the no-slip boundary condition is still given by (3.1)–(3.4) but needs to be interpreted in a different way. While $\widehat{\boldsymbol{u}}$ and $\widehat{p}$ represent small interior perturbations to the leading-order solution $\boldsymbol{u}_{0}$ and $p_{0}$ , $\widetilde{\boldsymbol{u}}$ denotes a strong viscous boundary flow with $|\widetilde{\boldsymbol{u}}|=O(|\boldsymbol{u}_{0}|)$ and an explicit solution for both the viscous boundary layer $\widetilde{\boldsymbol{u}}$ and its influx is required in the asymptotic analysis of the no-slip problem.

Substitution of asymptotic expansion (3.1)–(3.4) into (2.1)–(2.3) produces the leading-order interior problem describing non-dissipative thermal-inertial oscillation, which remains the same as that for the stress-free condition. At the next order, however, the governing equations for the secondary flow $\widehat{\boldsymbol{u}}$ in the interior of the sphere become

where higher-order terms, such as $E\unicode[STIX]{x1D6FB}^{2}\widehat{\boldsymbol{u}}$ , are neglected and the interior secondary flow $\widehat{\boldsymbol{u}}$ is subject to the boundary condition

It can be shown, through a standard manipulation (Greenspan Reference Greenspan1968), that the boundary-layer solution $\widetilde{\boldsymbol{u}}$ on the bounding surface of the sphere takes the form of

where $\unicode[STIX]{x1D709}=(1-r)/\sqrt{E}$ with $0<E\ll 1$ and

The boundary-layer solution (4.4) breaks down at the critical latitude, $\unicode[STIX]{x1D703}_{c}=\cos ^{-1}(\unicode[STIX]{x1D70E}_{0})$ , where the thickness of the layer changes from $O(E^{1/2})$ to $O(E^{2/5})$ (Stewartson & Roberts Reference Stewartson and Roberts1963). But it is anticipated that the existence of the critical latitude does not significantly affect the primary properties of convection at $0<E\ll 1$ because the influx from this small singular region is insignificant compared with the dominant influx from the rest of the spherical boundary layer. With the boundary-layer flow $\widetilde{\boldsymbol{u}}$ given by (4.4), an analytical expression for its influx can be written as

where

are real and functions of $\unicode[STIX]{x1D703}$ . Expression (4.7) is required in deriving the solvability condition required for inhomogeneous equations (4.1)–(4.2), which can be written in the form

Inserting (4.7) into (4.9) and taking its real part, we obtain an expression for the Rayleigh number $R$ :

where $P_{0}(\unicode[STIX]{x1D703})=\left[p_{0}(r,\unicode[STIX]{x1D703})\right]_{r=1}$ is a real function of $\unicode[STIX]{x1D703}$ , and the integral involving $\unicode[STIX]{x1D6E9}_{0}$ is given by (3.20). The smallest Rayleigh number can be determined by minimizing $R$ in (4.10) over different $p_{0},\boldsymbol{u}_{0}$ and $\unicode[STIX]{x1D70E}_{0}$ in connection with different values of $k$ . After determining the smallest Rayleigh number using (4.10), we then use the imaginary part of the solvability condition to determine the correction $\unicode[STIX]{x1D70E}_{1}$ , which is

where the integral involving $\unicode[STIX]{x1D6E9}_{0}$ is given by (3.22).

Our extensive calculation, based on asymptotic expression (4.10) over a large manifold of $p_{0},\boldsymbol{u}_{0}$ and $\unicode[STIX]{x1D6E9}_{0}$ with different $k$ and $\unicode[STIX]{x1D70E}_{0}$ , suggests that convection in the form of torsional oscillation cannot be physically preferred in any parameter regime of the problem; and this suggestion is also in agreement with our extensive numerical analysis. For the purpose of easy comparison, four different asymptotic solutions obtained at $E=10^{-4}$ for the no-slip condition are presented in figure 1(b) for different values of the Prandtl number $Pr$ : the curve labelled with A represents the asymptotic solution for non-axisymmetric, equatorially symmetric, retrogradely travelling waves with the azimuthal wavenumber $m=1$ , the curve C for progradely travelling waves with $m=1$ and the curve D for retrogradely travelling waves with $m=2$ , all of which are computed according to the asymptotic formula given by Zhang (Reference Zhang1995) for non-axisymmetric travelling waves; and the curve B represents the axisymmetric torsional mode calculated from (4.10). It is clear from figure 1(b) that the axisymmetric oscillation mode with the no-slip condition, in contrast to the curve B in figure 1(a) for the stress-free condition, is not physically preferred at any values of $Pr$ . Our calculations using different values of $E$ show largely similar behaviours, suggesting that the type of velocity boundary condition plays a critical role in determining the primary properties of convection for small $Pr$ in rapidly rotating spheres with $0<E\ll 1$ .

5 Summary and remarks

We have derived two asymptotic solutions for convection-driven torsional oscillation in rapidly rotating fluid spheres with the stress-free and no-slip boundary conditions; we have shown that the asymptotic solution described by (3.19)–(3.22) and (3.23)–(3.26) is in satisfactory quantitative agreement with the corresponding numerical solution; we believe that the convective flow given by (3.23)–(3.26) represents the simplest analytical solution in closed form in rotating spherical convective systems with $0<E\ll 1$ ; and our asymptotic analysis confirms, as suggested by the recent numerical simulation (Sanchez et al. Reference Sanchez, Garcia and Net2016), that convection in the form of axially symmetric and equatorially antisymmetric torsional oscillation is physically preferred for $0<E\ll 1$ in a special range of small Prandtl number with the stress-free boundary condition. However, our asymptotic analysis indicates that torsional convection cannot be physically preferred in rapidly rotating spheres with the no-slip boundary condition, revealing an essential role played by the velocity boundary condition. Of course, this does not rule out the possible relevance of the equatorially antisymmetric torsional mode in the nonlinear regime as indicated by the result of direct numerical simulation (Landeau & Aubert Reference Landeau and Aubert2000).

A satisfactory quantitative agreement achieved between the asymptotic and numerical results in turn validates the hypotheses/approaches adopted in our asymptotic analysis: the leading-order velocity and pressure of torsional convection for $0<E\ll 1$ are dominated by spherical inertial oscillation; the temperature is passively driven by the inertial oscillation and, hence, the leading-order problem describes non-dissipative thermal-inertial oscillation; there exist no simple general asymptotic scalings for $0<E\ll 1$ and, consequently, we have to tackle a set of partial differential equations that govern the problem; and different types of velocity boundary condition require different asymptotic approaches.

An interesting but intriguing property is that the convection-driven oscillatory flow $\boldsymbol{u}_{0}$ cannot generate, through its nonlinear interaction, any mean zonal flow in the vicinity of convective instabilities in rapidly rotating fluid spheres. This is because the Reynolds stress in connection with $\boldsymbol{u}_{0}$ is

where $0<|\unicode[STIX]{x1D70E}_{0}|<1$ . It is therefore important to understand the weakly and strongly nonlinear properties of torsional convection through either perturbation analysis or direct nonlinear numerical simulation in the parameter regime suggested by the asymptotic analysis in this paper. Finally, it is worth mentioning that the convective flow described by (3.24)–(3.26) – which is simple, time-dependent, fully analytical and in closed form – makes it highly desirable to study if this flow is capable of sustaining dynamo action in rapidly rotating spheres. Since the solution (3.24)–(3.26) represents a physically realizable convective flow in a rotating spherical geometry for $0<E\ll 1$ , its corresponding dynamo solution offers valuable insight into the basic generation mechanism of planetary fields.