Abstract

We investigate the influence of different thermal and velocity boundary conditions on numerical geodynamo models. We concentrate on the implications for magnetic field morphology, heat transport scaling laws, force balances and generation mechanisms. The field morphology most strongly depends on the local Rossby number, but there is some variation in the dipolarity of the field with boundary condition. Scaling laws also depend on the boundary conditions, but a diffusivity-free scaling is a good first order approximation for all our dipolar models. Our multipolar models, however, obey different scaling laws from dipolar models implying a different force balance in these models. We find that our dipolar models have a stronger degree of Lorentz–Coriolis balance compared to our multipolar models which have a stronger degree of Lorentz-inertial balance.The models with a stronger Lorentz–Coriolis dominance can be generated by either αω, α2ω or α2 mechanisms whereas the models with a stronger Lorentz-inertial balance are all α2 dynamos. These results imply that some caution is necessary when extrapolating results from dynamo models to Earth-like parameters since the choice of boundary conditions can have important effects.

1 INTRODUCTION

Earth's magnetic field is generated in the planet's fluid outer core through dynamo action. In this process, convection and differential rotation of an electrically conducting fluid maintain the magnetic field against its ohmic decay. In the last two decades, 3-D numerical geodynamo models driven by convection in rotating spherical shells have been successful in generating some of the basic properties of the observed geomagnetic field, including the field's dipole dominance, intensity and reversal statistics. This is despite the fact that, due to computational constraints, models operate in a parameter regime that is not representative of planetary interiors (Glatzmaier & Roberts 1996; Kuang & Bloxham 1997; Takahashi et al. 2005; Miyagoshi et al. 2010). A question arises as to how these chosen model parameters, as well as the boundary conditions, affect the numerical solutions. Specifically, since we are interested in determining core properties such as the core–mantle boundary (CMB) heat flux and ohmic dissipation to understand the Earth's heat budget, are dynamo models getting these right? In order to address this question, scaling studies of the influence of model parameters and boundary conditions on dynamo models are required.

1.1 Dynamo scaling laws

As computational capabilities improve, dynamo studies in previously unexplored parameter regimes become feasible. For example, Sakuraba & Roberts (2009) and Miyagoshi et al. (2011) have studied dynamo models at an Ekman number (ratio of viscous to Coriolis force) as low as 10−7. However, the Earth's Ekman number is about 10−15 and hence this example illustrates the many orders of magnitude separating model and planetary parameters. The most practical way to proceed, currently, is to run suites of models using numerically reasonable ranges of parameter values to develop scaling laws between the input and output parameters. The scaling laws must then be assumed to hold over the many orders of magnitude separating model and planetary parameters, and used to infer properties of the geodynamo.

Various studies have developed scaling laws from dynamo simulations for characteristics such as magnetic field strength, heat transport and dipolarity (e.g. Christensen & Tilgner 2004; Aubert 2005; Christensen & Aubert 2006; Olson & Christensen 2006; Takahashi et al. 2008; Aubert et al. 2009; Christensen 2010; Stelzer & Jackson 2013; Yadav et al. 2013). Some important results from these scaling law studies are that (1) the morphology of the magnetic field depends on the local Rossby number (a proxy for the importance of rotation on the convective motions) and (2) diffusivity-free scalings work well for predicting characteristics such as field strength (although a recent analysis calls this into question Stelzer & Jackson 2013). This has lead to the hypothesis that the resulting magnetic field strength depends only on the convective power available to drive the dynamo rather than the balance between Lorentz and Coriolis forces (the so called ‘magnetostrophic balance’ popular in previous literature). The scaling laws developed from these numerical simulations also appear to work well in explaining the field strength for some planets and stars (Christensen et al. 2009).

Most of the models in these simulations use fixed temperature and no-slip boundary conditions although a few of the studies do investigate the effects of some other mechanical and thermal boundary conditions. Several of the studies also incorporate an insulating inner core rather than a conducting inner core, an approximation that can affect force balances in the models (Dharmaraj & Stanley 2012).

1.2 The influence of boundary conditions

Various studies have demonstrated that the choice of thermal and mechanical boundary conditions in numerical models can affect the resulting flows and fields as well as the force balances in the models.

1.2.1 Thermal boundary conditions

Kutzner & Christensen (20002002) studied the magnetic field strength and reversal frequencies of dynamo models for a range of Rayleigh numbers with no-slip velocity boundary conditions and different thermal boundary conditions, with and without internal heating. They found that the choice of thermal boundary conditions made little difference. However, this study utilized a solid inner core that was electrically insulating as opposed to the more physically realistic choice of finitely electrically conducting, and recent work (Dharmaraj & Stanley 2012) has demonstrated that this can affect the results.

Another study on the influence of different thermal boundary conditions showed that at lower Ekman numbers, fixed heat flux models produce Earth-like, strong dipolar magnetic fields, unlike fixed temperature models that produce weak magnetic fields (Sakuraba & Roberts 2009). The difference in the solutions is attributed to large scale circulation, that is strong differential rotation is promoted by fixed heat flux models, which in turn coherently organizes the magnetic field and enhances the dipole component. In the fixed flux models, more heat is transported by convection in the equatorial region than in the polar regions, whereas the fixed temperature models redistribute excess heat more homogeneously at the boundaries.

Hori et al. (2010) studied the influence of thermocompositional boundary conditions on dynamos driven by internal heating. They showed that the impact of the thermal boundary conditions depended on the convective driving source. For example, changing the outer boundary condition only affected the solutions when the driving force was dominated by internal buoyancy sources like secular cooling. Similarly, changing the inner boundary condition only affected the solutions when the dominant driving force was at the inner core boundary, such as the release of light elements or latent heat from a growing inner core. This study was also performed using electrically insulating inner cores, and as mentioned earlier, this may affect the results.

1.2.2 Velocity boundary conditions

Studies provide a variety of results for the influence of different velocity boundary conditions on the dynamo. For example, Kuang & Bloxham (1997) and Kuang (1999) studied the influence of different velocity boundary conditions for a dynamo model with fixed heat flux thermal boundary conditions. They demonstrated that the velocity boundary condition significantly affects force balances in the models. In contrast, Christensen et al. (1999), using fixed temperature boundary conditions, demonstrated that although the choice of velocity boundary conditions influenced the zonal flows that develop in the core, the magnetic field morphology is broadly similar in both cases. Miyagoshi et al. (2010), also using fixed temperature boundaries, demonstrated that at lower Ekman numbers (near the limit of what is numerically achievable currently), the zonal flows are more similar for different mechanical boundary conditions.

Most of the variety in these results is due to the fact that other parameters in these studies are different, for example, different thermal boundary conditions, or control parameter values, or the treatment of the inner core. It would therefore be advantageous for a single study to systematical vary all the boundary conditions to determine the implications.

1.3 Dynamo generation mechanisms

Dynamo models have also been used to investigate the magnetic field generation mechanism in Earth's core (Christensen & Olson 1998; Christensen et al. 1999; Kuang 1999; Olson et al. 1999). The process by which differential rotation can produce toroidal magnetic field is known as the ‘ω effect’. The process by which helical convective upwelling produces both toroidal and poloidal magnetic fields is known as the ‘α effect’. If the toroidal magnetic field is dominantly produced by the ω effect then the dynamo is known as the ‘αω dynamo’, if it is dominantly produced by α effect then the dynamo is known as the ‘α2 dynamo’ and if α and ω effects are equally important, the dynamo is known as the ‘α2ω dynamo’. For a schematic illustration of the generation of poloidal from toroidal magnetic field and vice versa we refer the reader to Olson et al. (1999).

The generation mechanism of the dynamo depends on the force balance. Non-magnetic convection models demonstrate that the flow in a rapidly rotating spherical shell is organized by the Coriolis force into axial columns and the main force balance is between the Coriolis force and the pressure gradient force. In a dynamo model, the addition of a strong magnetic field can influence the fluid motion, breaking the axial column and the competition between the Coriolis, viscous and inertial forces has to accommodate the Lorentz force. It is often assumed that the influence of magnetic field will be important when the Elsasser number (ratio of Lorentz to Coriolis force) > O(1). However, recent studies have shown that the traditional definition of the Elsasser number is not a good representation of the Lorentz to Coriolis force balance. Dharmaraj & Stanley (2012) demonstrated that a comparison of the axisymmetric ϕ component of the Lorentz and Coriolis forces gives a better estimate of the Lorentz–Coriolis balance whereas Soderlund et al. (2012) showed that calculating volume integrals of the forces can do so as well.

In this paper, we use numerical dynamo models to investigate the influence of different thermal and velocity boundary conditions on scaling laws and dynamo generation mechanisms. In Section 2, we describe our numerical model and in Section 3 we present our results. The conclusions are then given in Section 4.

2 NUMERICAL MODEL

We use the Kuang and Bloxham numerical dynamo model (Kuang & Bloxham 1999) that solves the governing equations for a dynamo in a 3-D spherical, rotating, Boussinesq, electrically conducting fluid shell surrounding a solid inner core. Non-dimensionalization of the governing equations is done using the radius of the outer core ro as a length scale, the magnetic diffusion time |$\tau = r_{\rm o}^2/\eta$| as a time scale [where η = 1/(σμ0) is the magnetic diffusivity, σ is the electrical conductivity and μ0 is the magnetic permeability of free space], the magnetostrophic balance intensity |$B = \sqrt{(2 \Omega \mu _0 \rho \eta )}$| as a magnetic field scale (where Ω is the rotation rate and ρ is the fluid density), and hTro or ▵T/[rio(1 − rio)] as a temperature scale for the fixed heat flux or fixed temperature thermal boundary conditions, respectively (where hT is the temperature gradient at the inner core boundary, ▵T is the temperature difference from the inner core to outer core boundary and rio is the inner to outer core radius ratio which is held constant at rio = 12/35 for all models). The non-dimensional equations are:
\begin{eqnarray} {\rm Ro}_M \left( \frac{\partial }{\partial t} + \vec{v} \cdot \nabla \right) \vec{v} + \hat{z} \times \vec{v} &=& -\nabla p + \skew3\vec{J} \times \skew3\vec{B} + E \nabla ^2 \vec{v} + {\rm Ra} \Theta \vec{r} \nonumber\\ \end{eqnarray}
(1)
\begin{eqnarray} \left( \frac{\partial }{\partial t} - \nabla ^2 \right) \skew3\vec{B} &=& \nabla \times (\vec{v} \times \skew3\vec{B}) \end{eqnarray}
(2)
\begin{eqnarray} \left( \frac{\partial }{\partial t} - q_k \nabla ^2 \right) \Theta &=& -\vec{v} \cdot \nabla [T_0(r) + \Theta], \end{eqnarray}
(3)
where |$\vec{v}$| is the velocity field, |$\skew3\vec{B}$| is the magnetic field, Θ is the temperature perturbation field, T0(r) is the background state temperature solution and |$\hat{z}$| is a unit vector in the direction of the rotation axis. As discussed in Kuang & Bloxham (1999), this model only includes the axisymmetric inertial term in eq. (1). We demonstrate in the Supporting Information that this does not significantly affect our results.
The definition of the non-dimensional numbers used in the above equations are as follows:
\begin{eqnarray} \text {Magnetic Rossby number: } {\rm Ro}_{M} &=& \frac{\eta }{2 \Omega r_{\rm o}^2} \end{eqnarray}
(4)
\begin{eqnarray} \text {Ekman number: } E &=& \frac{\nu }{2 \Omega r_{\rm o}^2} \end{eqnarray}
(5)
\begin{eqnarray} \text {Modified Rayleigh number: } {\rm Ra}_{F} &=& \frac{\alpha _T g_{\rm o} h_T r_{\rm o}^2}{2 \Omega \eta } \end{eqnarray}
(6)
\begin{eqnarray} \text {(for fixed heat flux models)} \nonumber \\ \text {Modified Rayleigh number: } {\rm Ra}_{T} &=& \frac{\alpha _T g_{\rm o} \triangle T r_{\rm o}}{2 \Omega \eta r_{{\rm io}}(1-r_{{\rm io}})} \end{eqnarray}
(7)
\begin{eqnarray} \text {(for fixed temperature models)} \nonumber \\ \text {Roberts number: } q_k &=& \frac{\kappa }{\eta }, \end{eqnarray}
(8)
where αT is the coefficient of thermal expansion, go is the gravitational acceleration at ro, ν is the kinematic viscosity and κ is the thermal diffusivity.
All models have spherical harmonic degree Lmax = 58 and order mmax = 53. The number of radial grid points in the inner core is |$N_{r_{\rm i}}=36$| and in the outer core is either |$N_{r_{\rm o}}=64$| or 93 (see Tables S1 and S2 for specific values for each model run). We employ hyperdiffusivities to work at moderately supercritical Rayleigh numbers. The hyperdiffusivity is applied in such a way that the thermal and magnetic diffusivities and the kinematic viscosity depend on the spherical harmonic degree l as:
\begin{eqnarray} \kappa (l) &=& \kappa _0 [1+0.06(l-l_0)^2] \end{eqnarray}
(9)
\begin{eqnarray} \eta (l) &=& \eta _0 [1+0.06(l-l_0)^2] \end{eqnarray}
(10)
\begin{eqnarray} \nu (l) &=& \nu _0 [1+0.05(l-l_0)^2], \end{eqnarray}
(11)
where κ0, η0 and ν0 are constant values. l0 is a constant value that lies in the range 1 to Lmax and determines the minimum degree at which hyperdiffusivity onsets. We tried different l0's in the range of 0–65 and found that as long as l0 > 20 there was little effect on the solutions. For all the models presented in this paper we use l0 ≥ 30. As demonstrated in the Supporting Information, our use of hyperdiffusivity does not significantly affect our results.

The control parameters, boundary conditions and results are summarized in Table S1 for dynamo models and Table S2 for pure thermal convection models with no magnetic field. The boundary conditions given in the supporting tables are applied at both the inner and outer core boundaries. In all the dynamo models presented in this paper, the inner and outer cores have the same electrical conductivities. In our models, the inner core is free to oscillate in response to magnetic and viscous torques. Since we compare models with different thermal and velocity boundary conditions, for convenience, we will use the abbreviations: FFSF for fixed heat flux and stress-free models, FFNS for fixed heat flux and no-slip models, FTSF for fixed temperature and stress-free models and FTNS for fixed temperature and no-slip models. To differentiate between dynamo models and pure thermal convection models with no magnetic field, we use a suffix ‘C’ for ‘convection’ in the boundary condition abbreviation mentioned above. The simulations were run for many magnetic dipole diffusion times and the results presented were analysed after the initial transients decayed.

3 RESULTS

We have generated 149 dynamo models and 96 non-magnetic convection models with different combinations of boundary conditions and control parameters as presented in Tables S1 and S2. qk = 1 for all our models. We consider three choices for (E, Rom) pairs: (1) E = Rom = 2 × 10− 5, (2) E = 6 × 10− 5, Rom = 2 × 10− 5 and (3) E = 2 × 10− 4, Rom = 1 × 10− 4 and then consider a range of Ra and boundary conditions for each (E, Rom) pair. If one is interested in the corresponding values for the Prandtl (Pr) number and magnetic Prandtl (Pm) number, they can be found using:
\begin{eqnarray} {\rm Pm} & = & \frac{\nu }{\eta } = \frac{E}{{\rm Ro}_m} \end{eqnarray}
(12)
\begin{eqnarray} {\rm Pr} & = & \frac{\nu }{\kappa } = \frac{{\rm P}{\rm m}}{q_k}. \end{eqnarray}
(13)
We begin by examining dynamo models to characterize the class of magnetic fields that they produce. We examine the heat transfer of dipole dominated dynamos by obtaining scaling laws, and compare them to non-magnetic convection models while examining if there is a dependence on boundary conditions. We investigate if there are any differences between dipolar dynamo models with different boundary conditions. Finally to understand the difference between dipolar and multipolar dynamos we examine the force balance and generation mechanisms of characteristic models.

3.1 Dynamo regimes and scaling laws

3.1.1 Dynamo regimes

Inspired by the work of Christensen & Aubert (2006) and Olson & Christensen (2006), we investigate whether the magnetic field morphology depends on flow scale for our dynamo models. We considered the subset of our models that matched the same criteria as used in Christensen & Aubert (2006), namely, models with Ekman number less than or equal to 6.5 × 10−5 (using our definition of the Ekman number) and Nusselt number (Nu) greater than 2. In Fig. 1, we plot the relative dipole strength defined as:
\begin{equation} f_{{\rm dip}} = \left[ \frac{\displaystyle \int _S \skew3\vec{B}_{l=1}(r=r_{\rm o})\cdot \skew3\vec{B}_{l=1}(r=r_{\rm o}) \, {\rm d}S}{\displaystyle \int _S \left(\sum \limits _{l=1}^{12} \left[\skew3\vec{B}_l(r=r_{\rm o})\cdot \skew3\vec{B}_l(r=r_{\rm o})\right]\right) {\rm d}S} \right]^{1/2} \end{equation}
(14)
(where S is the CMB surface) versus the local Rossby number:
\begin{equation} {\rm Ro}_l = {\rm Ro} \frac{l_v}{\pi }. \end{equation}
(15)
Ro = UD is the standard Rossby number based on the spherical shell thickness D and a characteristic speed U. We use the Christensen & Aubert (2006) choice for U as the rms of the velocity in the shell. The local Rossby number incorporates the characteristic length scale of the flow by dividing the Rossby number by the characteristic half-wavelength of the flow (π/lv). lv is the mean spherical harmonic degree of the kinetic energy spectrum and is given by:
\begin{eqnarray} l_v &=& \frac{\sum l \langle v_l \cdot v_l \rangle }{E_k}, \end{eqnarray}
(16)
where vl is the component of the flow at degree l, 〈〉 denotes the rms, the sum is taken over all degrees and Ek is the non-dimensional kinetic energy:
\begin{eqnarray} E_k &=& \int _V |\vec{v}|^2 \rm{d}V. \end{eqnarray}
(17)
We find a similar dependence of dynamo dipolarity on Rol as previous studies. Irrespective of the control parameters and boundary conditions, most models follow the transition from dipolar (fdip > 0.35) to non-dipolar (fdip < 0.35) as Rol increases, with a critical value for the transition around Rol = 0.12. Some of our models are bistable, that is for the same set of control parameters and boundary conditions, the field morphology strongly depends on the initial magnetic field. Bistability is especially prominent in fixed flux and stress-free models as also observed by Yadav et al. (2013), Simitev & Busse (2009), Simitev & Busse (2012), Sasaki et al. (2011), Schrinner et al. (2012) and Gastine et al. (2012). Ignoring the weak bistable solutions, our results imply that the transition from dipolar to multipolar dynamos is governed by the local Rossby number and the transition point is independent of the boundary conditions chosen in the models. However it should be mentioned that for FFNS boundary conditions, it is harder to draw such a conclusion since we have no models at Rol > 0.12 that meet the Nu > 2, E < 6.5 × 10− 5 criteria.
Figure 1.

Relative dipole strength and magnetic energy versus local Rossby number: (a) relative dipole strength versus local Rossby number, and (b) magnetic energy versus local Rossby number. The bistable models are connected with a black dashed line. Panels (c) and (d) are (a) and (b) with the bistable models removed, Nu > 2 and E < 6.5 × 10−5. Red denotes fixed flux and stress-free (FFSF) models, magenta denotes fixed flux and no-slip (FFNS) models, blue denotes fixed temperature and stress-free (FTSF) models, cyan denotes fixed temperature and no-slip (FTNS) models, and a green plus sign inside the symbols denotes bistable models.

For a given Rol in the dipolar regime, the dipolarity fdip does weakly depend on the boundary conditions. As can be observed in Fig. 1, fixed temperature models are generally more dipolar than fixed flux models. For the fixed flux models, stress-free models are more dipolar than no-slip models (ignoring the weak bistable solutions).

We also plot the time averaged non-dimensional magnetic energy:
\begin{eqnarray} E_B &=& \int _V |\skew3\vec{B}|^2 {\rm d}V \end{eqnarray}
(18)
versus the local Rossby number in Figs 1(b) and (d). Note that based on our choice of non-dimensionalization, EB/V is the Elsasser number. Also, in order to compare the magnitudes of the dimensional energies, the ratio of the dimensional energies is related to the ratio of the non-dimensional energies as:
\begin{eqnarray} \frac{E^d_k}{E^d_B} &=& {\rm Ro}_m \frac{E_k}{E_m} , \end{eqnarray}
(19)
where the superscript d refers to dimensional quantities. Although there is significant scatter in Figs 1(b) and (d), and the energy depends on properties other than just Rol, it appears that for low Rol the magnetic energy increases rapidly but then stalls around the constant value of EB = 10 for models approaching the transition zone. Once past the transition zone, the energy begins increasing again. This suggests that there may be some new dynamics or force balances occurring for models near the transition region between dipolar and multipolar solutions. This is important since the geodynamo may have Rol in this region (Olson & Christensen 2006).

3.1.2 Heat transport

In order to examine the influence of boundary conditions on heat transport, we begin by searching for power laws. First, we plot Nu, which is the ratio of total heat flux to conductive heat flux at the CMB versus the Rayleigh number for all our models in Figs 2(a: dynamo models) and (b: non-magnetic convection models) that meet the Christensen & Aubert (2006) criteria. Since the Rayleigh number definitions for fixed heat flux and fixed temperature boundary conditions are different in our models (i.e. eqs 6 and 7), we calculate the RaT value for the fixed flux models using the relationship:
\begin{eqnarray} {\rm Ra}_{T} = \frac{{\rm Ra}_{F}}{{\rm Nu}} \end{eqnarray}
(20)
and use RaT for our horizontal axis. In both plots the Nusselt number increases as the Rayleigh number increases, as expected, but there is a lot of scatter and the results do not fall on a single straight line. Dynamo models seem to diverge at high Rayleigh numbers, whereas non-magnetic convection models seem to converge. The solid red and blue lines represent the best fitting lines for fixed flux and fixed temperature models, respectively, with fdip > 0.35 [i.e. dipolar dynamo models for comparison with Christensen & Aubert (2006)]. The multipolar models at higher Ra number do not obey the same scaling laws as the dipolar models.
Figure 2.

Heat transfer power laws: dynamo models are shown in the left-hand column and convection models with no magnetic field are shown in the right-hand column. The first row plots Nusselt number versus Rayleigh number. The second row plots Nusselt number versus Rayleigh number normalized by it's critical value. The third row plots the modified Nusselt number versus modified flux based Rayleigh number. Symbols as in Fig. 1 without the non-dipolar bistable models. The brown plus signs inside symbols denote multipolar models which have not been included when determining the best-fitting lines since they appear to diverge from the trend. The blue and red lines denote power law fits for fixed temperature (blue and cyan symbols) and fixed heat flux (red and magenta symbols) models, respectively.

The effect of different boundary conditions can be seen more clearly by plotting Nu versus the supercriticality (Ra/RaC) rather than versus Ra. In Figs 2(c) and (d), we see a clear difference between the fixed flux and fixed temperature models, even though there is a significant amount of scatter, especially in the fixed temperature models. Also, the power-law exponents (i.e. the line slope) of all the dipolar dynamo models are similar to the non-magnetic convection models, and the dipolar dynamo models obey similar scaling laws irrespective of the velocity boundary conditions.

Next, for comparison to the diffusivity-free power laws of Christensen & Aubert (2006), we plot the modified Nusselt number (that is independent of thermal diffusivity) versus the modified flux based Rayleigh number (that is independent of diffusivities and based on the advected heat flux) in Figs 2(e: dynamo models) and (f: non-magnetic convection models), where:
\begin{eqnarray} \rm {Modified\; Nusselt \#: } {\rm Nu}^* &=& ({\rm Nu}-1) \ q_k \ {\rm Ro}_M \ \frac{2}{(1-r_{{\rm io}})^2} \end{eqnarray}
(21)
\begin{eqnarray} \rm {Modified\; Flux-based\; Rayleigh \#: } {\rm Ra}_Q^* &=& 4r_{{\rm io}} \ {\rm Ra} \ {\rm Ro}_M \ {\rm Nu}^*. \end{eqnarray}
(22)
We use the Rayleigh number relationship given above in eq. (20) to plot all the models with different thermal boundary conditions on the same horizontal axis. Compared to Figs 2(a)–(d), the scatter is significantly reduced and we obtain power laws similar to Christensen & Aubert (2006). The power law dependence of dipolar dynamo models is shown in Table 2. The power law exponent of our fixed temperature dipolar dynamo models with Nu > 2 is 0.47, similar to Stelzer & Jackson (2013) [a more rigorous analysis of Christensen & Aubert (2006) data set]. Again, it should be noted that the multipolar models clearly do not follow the trend.

3.2 Field morphology

In order to understand the influence of different boundary conditions on the models, we select specific ‘characteristic’ models with the following identical parameters: E = 6 × 10−5, RoM = 2 × 10− 5 and qk=1 but varying boundary conditions and dipolarity. Since the Rayleigh number is defined differently for fixed flux and fixed temperature models, we choose our characteristic models at similar Nusselt numbers. The characteristic dipolar dynamo models at Nu ≈ 5 are: FFSF-21, FFNS-15, FTSF-15 and FTNS-24, and the characteristic multipolar dynamo models at Nu ≈ 11.5 are: FTSF-20 and FTNS-33. The corresponding non-magnetic convection models at similar Nusselt number are: FFSFC-16, FFNSC-14, FTSFC-09, FTNSC-08, FTSFC-13 and FTNSC-15. The control parameters and boundary conditions of the models are given in Table 1.

Table 1.

Characteristic model results. The column headings are: model ID, Rayleigh number (Ra), supercriticality of Rayleigh number (⁠|$\frac{{\rm Ra}}{{\rm Ra}_C}$|⁠), magnetic energy (EB), kinetic energy (Ek), ratio of poloidal to toroidal magnetic energy (⁠|$\frac{E_{B_P}}{E_{B_T}}$|⁠), relative dipole strength (fDip), Rossby number (Ro), local Rossby number (Rol), Nusselt number (Nu), Match Coefficient of the Lorentz to Coriolis force (MLC) and the Match Coefficient of α and ω effects (Mα and Mω), respectively. All models use hyperdiffusivity parameter l0 = 30. All models were run for at least three magnetic dipole diffusion times and results are based on analysis after the system equilibrated.

Dynamo models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EBEk|$\frac{E_{B_P}}{E_{B_T}}$|fDipRoRolNuMLCMαMω
FFSF-212.50e+0487.418.1e+003.2e+058.8e−018.3e–011.1e–024.6e–025.360.450.390.23
FFNS-152.00e+0482.994.5e+001.6e+051.0e+007.5e–018.1e–033.0e–025.210.420.520.17
FTSF-155.00e+0314.494.6e+003.1e+051.3e+009.3e–011.1e–025.5e–025.060.360.410.27
FTNS-245.00e+0318.056.6e+001.7e+051.3e+008.1e–018.3e–033.4e–025.050.470.390.32
FTSF-202.00e+0457.971.2e+014.8e+068.1e-011.7e–014.4e–021.5e–0111.740.090.600.19
FTNS-333.00e+04108.301.3e+013.3e+069.1e−019.7e–023.6e–021.4e–0111.750.020.670.12
Convection models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EkRoRolNu
FFSFC-163.50e+04122.381.5e+077.6e–028.6e–024.83
FFNSC-142.50e+04103.735.9e+051.5e–023.7e–025.20
FTSFC-097.50e+0321.746.5e+065.1e–025.4e–024.87
FTNSC-085.00e+0318.054.2e+051.3e–024.1e–024.65
FFSFC-132.00e+0469.932.6e+063.2e–024.3e–023.96
FTNSC-153.00e+04108.301.0e+076.3e–021.8e–0111.63
Dynamo models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EBEk|$\frac{E_{B_P}}{E_{B_T}}$|fDipRoRolNuMLCMαMω
FFSF-212.50e+0487.418.1e+003.2e+058.8e−018.3e–011.1e–024.6e–025.360.450.390.23
FFNS-152.00e+0482.994.5e+001.6e+051.0e+007.5e–018.1e–033.0e–025.210.420.520.17
FTSF-155.00e+0314.494.6e+003.1e+051.3e+009.3e–011.1e–025.5e–025.060.360.410.27
FTNS-245.00e+0318.056.6e+001.7e+051.3e+008.1e–018.3e–033.4e–025.050.470.390.32
FTSF-202.00e+0457.971.2e+014.8e+068.1e-011.7e–014.4e–021.5e–0111.740.090.600.19
FTNS-333.00e+04108.301.3e+013.3e+069.1e−019.7e–023.6e–021.4e–0111.750.020.670.12
Convection models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EkRoRolNu
FFSFC-163.50e+04122.381.5e+077.6e–028.6e–024.83
FFNSC-142.50e+04103.735.9e+051.5e–023.7e–025.20
FTSFC-097.50e+0321.746.5e+065.1e–025.4e–024.87
FTNSC-085.00e+0318.054.2e+051.3e–024.1e–024.65
FFSFC-132.00e+0469.932.6e+063.2e–024.3e–023.96
FTNSC-153.00e+04108.301.0e+076.3e–021.8e–0111.63
Table 1.

Characteristic model results. The column headings are: model ID, Rayleigh number (Ra), supercriticality of Rayleigh number (⁠|$\frac{{\rm Ra}}{{\rm Ra}_C}$|⁠), magnetic energy (EB), kinetic energy (Ek), ratio of poloidal to toroidal magnetic energy (⁠|$\frac{E_{B_P}}{E_{B_T}}$|⁠), relative dipole strength (fDip), Rossby number (Ro), local Rossby number (Rol), Nusselt number (Nu), Match Coefficient of the Lorentz to Coriolis force (MLC) and the Match Coefficient of α and ω effects (Mα and Mω), respectively. All models use hyperdiffusivity parameter l0 = 30. All models were run for at least three magnetic dipole diffusion times and results are based on analysis after the system equilibrated.

Dynamo models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EBEk|$\frac{E_{B_P}}{E_{B_T}}$|fDipRoRolNuMLCMαMω
FFSF-212.50e+0487.418.1e+003.2e+058.8e−018.3e–011.1e–024.6e–025.360.450.390.23
FFNS-152.00e+0482.994.5e+001.6e+051.0e+007.5e–018.1e–033.0e–025.210.420.520.17
FTSF-155.00e+0314.494.6e+003.1e+051.3e+009.3e–011.1e–025.5e–025.060.360.410.27
FTNS-245.00e+0318.056.6e+001.7e+051.3e+008.1e–018.3e–033.4e–025.050.470.390.32
FTSF-202.00e+0457.971.2e+014.8e+068.1e-011.7e–014.4e–021.5e–0111.740.090.600.19
FTNS-333.00e+04108.301.3e+013.3e+069.1e−019.7e–023.6e–021.4e–0111.750.020.670.12
Convection models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EkRoRolNu
FFSFC-163.50e+04122.381.5e+077.6e–028.6e–024.83
FFNSC-142.50e+04103.735.9e+051.5e–023.7e–025.20
FTSFC-097.50e+0321.746.5e+065.1e–025.4e–024.87
FTNSC-085.00e+0318.054.2e+051.3e–024.1e–024.65
FFSFC-132.00e+0469.932.6e+063.2e–024.3e–023.96
FTNSC-153.00e+04108.301.0e+076.3e–021.8e–0111.63
Dynamo models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EBEk|$\frac{E_{B_P}}{E_{B_T}}$|fDipRoRolNuMLCMαMω
FFSF-212.50e+0487.418.1e+003.2e+058.8e−018.3e–011.1e–024.6e–025.360.450.390.23
FFNS-152.00e+0482.994.5e+001.6e+051.0e+007.5e–018.1e–033.0e–025.210.420.520.17
FTSF-155.00e+0314.494.6e+003.1e+051.3e+009.3e–011.1e–025.5e–025.060.360.410.27
FTNS-245.00e+0318.056.6e+001.7e+051.3e+008.1e–018.3e–033.4e–025.050.470.390.32
FTSF-202.00e+0457.971.2e+014.8e+068.1e-011.7e–014.4e–021.5e–0111.740.090.600.19
FTNS-333.00e+04108.301.3e+013.3e+069.1e−019.7e–023.6e–021.4e–0111.750.020.670.12
Convection models: E = 6.00e–05 , RoM = 2.00e–05, qk = 1, |$N_{r_0}=64$|
ModelRa|$\frac{{\rm Ra}}{{\rm Ra}_C}$|EkRoRolNu
FFSFC-163.50e+04122.381.5e+077.6e–028.6e–024.83
FFNSC-142.50e+04103.735.9e+051.5e–023.7e–025.20
FTSFC-097.50e+0321.746.5e+065.1e–025.4e–024.87
FTNSC-085.00e+0318.054.2e+051.3e–024.1e–024.65
FFSFC-132.00e+0469.932.6e+063.2e–024.3e–023.96
FTNSC-153.00e+04108.301.0e+076.3e–021.8e–0111.63
Table 2.

Scaling laws for dynamo and non-magnetic convection models. The subscript C denotes convection models. ▵a and ▵b are standard errors in the offset or Y-intercept (a) and slope (b), respectively. These power laws are obtained for models with Nu > 2 and E < 3 × 10−4. Dipolar dynamo models have fDip > 0.35 and multipolar dynamo models have fDip < 0.35. Ignore the polarity for convection models.

a ± ▵ab ± ▵baC ± ▵aCbC ± ▵bC
|${\rm Nu} = (a \pm \triangle a) \ {\rm Ra}_{FT}^{(b \pm \triangle b)}$|
All0.034 ± 0.0040.59 ± 0.010.017 ± 0.0050.65 ± 0.03
Dipolar fixed temperature0.02 ± 0.0030.65 ± 0.020.019 ± 0.0070.64 ± 0.04
Dipolar fixed flux0.013 ± 0.0030.71 ± 0.030.016 ± 0.0090.66 ± 0.07
Multipolar fixed temperature0.38 ± 0.10.34 ± 0.03
|${\rm Nu} = (a \pm \triangle a) \ (\frac{{\rm Ra}}{{\rm Ra}_C})^{(b \pm \triangle b)}$|
All1.6 ± 0.20.33 ± 0.031.6 ± 0.40.3 ± 0.06
Dipolar fixed temperature1.1 ± 0.060.57 ± 0.030.94 ± 0.10.58 ± 0.04
Dipolar fixed flux0.82 ± 0.020.41 ± 0.0070.66 ± 0.040.43 ± 0.01
Multipolar fixed temperature3.7 ± 0.50.25 ± 0.03
|${\rm Nu}^* = (a \pm \triangle a) \ {\rm Ra}_Q^{*(b \pm \triangle b)}$|
All0.027 ± 0.0020.43 ± 0.0070.033 ± 0.0040.47 ± 0.01
Dipolar fixed temperature0.04 ± 0.0040.47 ± 0.0080.028 ± 0.0040.45 ± 0.02
Dipolar fixed flux0.053 ± 0.0070.49 ± 0.010.052 ± 0.010.51 ± 0.02
Multipolar fixed temperature0.0074 ± 0.0010.28 ± 0.02
a ± ▵ab ± ▵baC ± ▵aCbC ± ▵bC
|${\rm Nu} = (a \pm \triangle a) \ {\rm Ra}_{FT}^{(b \pm \triangle b)}$|
All0.034 ± 0.0040.59 ± 0.010.017 ± 0.0050.65 ± 0.03
Dipolar fixed temperature0.02 ± 0.0030.65 ± 0.020.019 ± 0.0070.64 ± 0.04
Dipolar fixed flux0.013 ± 0.0030.71 ± 0.030.016 ± 0.0090.66 ± 0.07
Multipolar fixed temperature0.38 ± 0.10.34 ± 0.03
|${\rm Nu} = (a \pm \triangle a) \ (\frac{{\rm Ra}}{{\rm Ra}_C})^{(b \pm \triangle b)}$|
All1.6 ± 0.20.33 ± 0.031.6 ± 0.40.3 ± 0.06
Dipolar fixed temperature1.1 ± 0.060.57 ± 0.030.94 ± 0.10.58 ± 0.04
Dipolar fixed flux0.82 ± 0.020.41 ± 0.0070.66 ± 0.040.43 ± 0.01
Multipolar fixed temperature3.7 ± 0.50.25 ± 0.03
|${\rm Nu}^* = (a \pm \triangle a) \ {\rm Ra}_Q^{*(b \pm \triangle b)}$|
All0.027 ± 0.0020.43 ± 0.0070.033 ± 0.0040.47 ± 0.01
Dipolar fixed temperature0.04 ± 0.0040.47 ± 0.0080.028 ± 0.0040.45 ± 0.02
Dipolar fixed flux0.053 ± 0.0070.49 ± 0.010.052 ± 0.010.51 ± 0.02
Multipolar fixed temperature0.0074 ± 0.0010.28 ± 0.02
Table 2.

Scaling laws for dynamo and non-magnetic convection models. The subscript C denotes convection models. ▵a and ▵b are standard errors in the offset or Y-intercept (a) and slope (b), respectively. These power laws are obtained for models with Nu > 2 and E < 3 × 10−4. Dipolar dynamo models have fDip > 0.35 and multipolar dynamo models have fDip < 0.35. Ignore the polarity for convection models.

a ± ▵ab ± ▵baC ± ▵aCbC ± ▵bC
|${\rm Nu} = (a \pm \triangle a) \ {\rm Ra}_{FT}^{(b \pm \triangle b)}$|
All0.034 ± 0.0040.59 ± 0.010.017 ± 0.0050.65 ± 0.03
Dipolar fixed temperature0.02 ± 0.0030.65 ± 0.020.019 ± 0.0070.64 ± 0.04
Dipolar fixed flux0.013 ± 0.0030.71 ± 0.030.016 ± 0.0090.66 ± 0.07
Multipolar fixed temperature0.38 ± 0.10.34 ± 0.03
|${\rm Nu} = (a \pm \triangle a) \ (\frac{{\rm Ra}}{{\rm Ra}_C})^{(b \pm \triangle b)}$|
All1.6 ± 0.20.33 ± 0.031.6 ± 0.40.3 ± 0.06
Dipolar fixed temperature1.1 ± 0.060.57 ± 0.030.94 ± 0.10.58 ± 0.04
Dipolar fixed flux0.82 ± 0.020.41 ± 0.0070.66 ± 0.040.43 ± 0.01
Multipolar fixed temperature3.7 ± 0.50.25 ± 0.03
|${\rm Nu}^* = (a \pm \triangle a) \ {\rm Ra}_Q^{*(b \pm \triangle b)}$|
All0.027 ± 0.0020.43 ± 0.0070.033 ± 0.0040.47 ± 0.01
Dipolar fixed temperature0.04 ± 0.0040.47 ± 0.0080.028 ± 0.0040.45 ± 0.02
Dipolar fixed flux0.053 ± 0.0070.49 ± 0.010.052 ± 0.010.51 ± 0.02
Multipolar fixed temperature0.0074 ± 0.0010.28 ± 0.02
a ± ▵ab ± ▵baC ± ▵aCbC ± ▵bC
|${\rm Nu} = (a \pm \triangle a) \ {\rm Ra}_{FT}^{(b \pm \triangle b)}$|
All0.034 ± 0.0040.59 ± 0.010.017 ± 0.0050.65 ± 0.03
Dipolar fixed temperature0.02 ± 0.0030.65 ± 0.020.019 ± 0.0070.64 ± 0.04
Dipolar fixed flux0.013 ± 0.0030.71 ± 0.030.016 ± 0.0090.66 ± 0.07
Multipolar fixed temperature0.38 ± 0.10.34 ± 0.03
|${\rm Nu} = (a \pm \triangle a) \ (\frac{{\rm Ra}}{{\rm Ra}_C})^{(b \pm \triangle b)}$|
All1.6 ± 0.20.33 ± 0.031.6 ± 0.40.3 ± 0.06
Dipolar fixed temperature1.1 ± 0.060.57 ± 0.030.94 ± 0.10.58 ± 0.04
Dipolar fixed flux0.82 ± 0.020.41 ± 0.0070.66 ± 0.040.43 ± 0.01
Multipolar fixed temperature3.7 ± 0.50.25 ± 0.03
|${\rm Nu}^* = (a \pm \triangle a) \ {\rm Ra}_Q^{*(b \pm \triangle b)}$|
All0.027 ± 0.0020.43 ± 0.0070.033 ± 0.0040.47 ± 0.01
Dipolar fixed temperature0.04 ± 0.0040.47 ± 0.0080.028 ± 0.0040.45 ± 0.02
Dipolar fixed flux0.053 ± 0.0070.49 ± 0.010.052 ± 0.010.51 ± 0.02
Multipolar fixed temperature0.0074 ± 0.0010.28 ± 0.02

3.2.1 Surface magnetic field

The difference in the observed magnetic field can be seen using surface magnetic power spectra. In a current free region (e.g. outside the Earth's core), the magnetic field can be represented as the gradient of the scalar potential |$(\skew3\vec{B}=-\nabla \Phi )$|⁠. The scalar potential can be expanded in spherical harmonics as follows:
\begin{eqnarray} \Phi &=& a \sum _{l=1}^{\infty } \sum _{m=0}^{l} \!\left( \frac{a}{r} \right)^{l+1} P_l^m (\cos \theta ) \left[ g_l^m \cos (m \phi ) \!+\! h_l^m \sin (m \phi )\right], \end{eqnarray}
(23)
where a is the radius of planet, (r, θ, ϕ) are spherical coordinates, |$P_l^m$| are Schmidt-normalized associated Legendre polynomials and |$g_l^m$| and |$h_l^m$| are Gauss coefficients for each spherical harmonic degree l and order m. The power p in each degree and order can be defined in terms of the mean square field intensity as follows:
\begin{eqnarray} p(l,m,r) &=& (l+1)\left (\frac{a}{r}\right )^{2l+4} \left[ (g_l^m)^2 +(h_l^m)^2\right]. \end{eqnarray}
(24)
In Fig. 3, we show the resulting average surface magnetic power spectra over a magnetic dipole diffusion time for our characteristic dipolar and multipolar dynamo models. As expected, the dipolar models result in axially dipolar dominated fields while the multipolar models result in non-axisymmetric, non-dipolar magnetic fields. Notably, the spectrum does not strongly depend on the thermal or velocity boundary conditions.
Figure 3.

Surface magnetic power spectra: Left-hand plot shows power versus degree normalized to the dipole (l=1) power and the right-hand plot shows power versus order normalized to axisymmetric (m = 0) power. The m = 0 power in the right-hand figure is artificially plotted at m = 0.4 in order to use a log axis.

3.2.2 Magnetic field

In this section, we examine the characteristics of the dynamo models inside the dynamo generation region to understand the difference in the potential fields observed outside as described in the previous section. The magnetic energies of the characteristic multipolar models (that are found at higher Rol) are approximately twice as strong as their respective characteristic dipolar models with similar boundary conditions. The toroidal magnetic field intensities are approximately equal in magnitude to the poloidal magnetic field intensity in all the characteristic dipolar and multipolar models. The multipolar dynamo models are dominated by their non-axisymmetric toroidal and poloidal magnetic fields, whereas the dipolar dynamo models have approximately equal contribution from their axisymmetric and non-axisymmetric components.

Next we examine the magnetic field morphology of dipolar and multipolar dynamos with different boundary conditions at the CMB and inside the dynamo generation region. Although Figs 4 and 5 present snapshots in time of the resulting magnetic fields of the characteristic dipolar and multipolar models, their features are similar at other times. Fig. 4 presents the radial magnetic field at the CMB. In the dipolar regime, all the models produce strong magnetic field that is mainly dipolar with some non-axisymmetric features, whereas the models in the multipolar regime produce intense small-scale non-dipolar and non-axisymmetric magnetic fields at the CMB as expected. Even though the fixed flux models produce dipolar dynamos, they produce intense patches at the poles, whereas the fixed temperature models are more homogeneous.

Figure 4.

Radial magnetic field at the core mantle boundary: for (a) FFSF-21, (b) FFNS-15, (c) FTSF-15, (d) FTNS-24, (e) FTSF-20 and (f) FTNS-33 models. Red and blue colours represent different field directions. The units are non-dimensional.

Figure 5.

Axisymmetric magnetic field: for (a) FFSF-21, (b) FFNS-15, (c) FTSF-15, (d) FTNS-24, (e) FTSF-20 and (f) FTNS-33 models. The left-hand (right) plot shows the meridional slice of the axisymmetric toroidal (streamlines of poloidal) magnetic field. In the toroidal plots, red shades represent prograde direction and blue shades represent retrograde direction. In the poloidal plots, red is clockwise and blue is anticlockwise. The units are non-dimensional.

Fig. 5 presents axisymmetric meridional slices of the toroidal and poloidal magnetic fields. The magnetic field morphologies of the models in the dipolar regime are different from the models in the multipolar regime. The dipolar models produce large scale toroidal flux patches that are antisymmetric at the equator while the multipolar models produce small scale flux patches that are scattered in the bulk of the outer core. The magnetic field morphology does not strongly depend on the choice of thermal or velocity boundary conditions.

3.2.3 Velocity field

In this section, we examine the characteristics of the flow velocities inside the dynamo generation region and compare them to non-magnetic convection models at similar parameters and boundary conditions. The kinetic energies of the characteristic multipolar dynamo models are an order of magnitude bigger than the dipolar dynamo models. All the characteristic dynamo models are dominated by their toroidal kinetic energies, but the dominance is more prominent in the stress-free models. The characteristic dipolar dynamo models are dominated by the non-axisymmetric component of the toroidal kinetic energies, whereas the characteristic multipolar dynamo models have equal contribution from the axisymmetric and non-axisymmetric components.

The kinetic energies of the characteristic non-magnetic convection models are an order of magnitude bigger than their corresponding characteristic dynamo models. All the characteristic non-magnetic convection models are also dominated by their toroidal kinetic energies. Again, the dominance is more prominent in the stress-free models. Unlike the dynamo models, the toroidal energy is dominated by its axisymmetric component for all non-magnetic convection models.

Next, we examine the velocity field morphology of characteristic dynamo models and compare them to non-magnetic convection models. Figs 6 and 7 present typical snapshots of the axisymmetric velocity fields for dynamo models and non-magnetic convection models, respectively. The velocity fields of the models in the dipolar regime are different from both the models in the multipolar regime and the convection models with no magnetic fields. The zonal flows of the multipolar models are stronger than the dipolar models. The toroidal velocity fields of the convection models with no magnetic field are approximately constant on cylinders parallel to the rotation axis, a sign of strong Taylor–Proudman influence. The velocity field morphologies of the multipolar models are similar to the convection models (also seen by Aubert 2005) even though they produce strong magnetic fields. In this case, the strong magnetic fields of multipolar models do not influence the large scale velocity fields.

Figure 6.

Axisymmetric velocity field: for (a) FFSF-21, (b) FFNS-15, (c) FTSF-15, (d) FTNS-24, (e) FTSF-20 and (f) FTNS-33 dynamo models. The toroidal velocity fields are shown on the left-hand side where red (blue) denotes prograde (retrograde) circulation and streamlines of poloidal velocity fields are shown on the right-hand side where red (blue) denotes clockwise (anticlockwise) circulation.

Figure 7.

Axisymmetric velocity field: for (a) FFSFC-16, (b) FFNSC-14, (c) FTSFC-09, (d) FTNSC-08, (e) FTSFC-13 and (f) FTNSC-15 convection models with no magnetic field. The toroidal velocity fields are shown on the left-hand side where red (blue) denotes prograde (retrograde) circulation and streamlines of poloidal velocity fields are shown on the right-hand side where red (blue) denotes clockwise (anticlockwise) circulation.

On the other hand, the dipolar models have strong differential rotation inside the tangent cylinder and a retrograde equatorial jet outside the tangent cylinder close to the outer core boundary; a sign that Lorentz forces are strong enough to break the Taylor–Proudman constraint. The meridional circulation of the dipolar dynamo models are predominantly equatorially antisymmetric throughout the core. The intensity of the zonal flows of the dynamo models are weaker than the convection models without magnetic field due to magnetic braking. Notably, the velocity field morphologies do not strongly depend on the choice of thermal or velocity boundary conditions.

3.3 Force balances

In this section, we examine the force balances of the characteristic dipolar and multipolar dynamo models. As demonstrated in our previous work (Dharmaraj & Stanley 2012), using the traditional definition of the Elsasser number to determine the force balance is inadequate and we therefore directly plot the axisymmetric ϕ component of the Lorentz and Coriolis forces in Fig. 8. In order to quantify the force balance we define the ‘Match Coefficient’ (MLC) as follows:
\begin{eqnarray} M_{{\rm LC}} &=& \int \frac{F_{\rm L} \ * \ F_{\rm C}}{\left(F_{\rm L}^2 \ + \ F_{\rm C}^2\right)/2} \ \rm{d}V, \end{eqnarray}
(25)
where |$F_{\rm L} = (\skew3\vec{J} \times \skew3\vec{B})_{\phi }$| is the ϕ component of the Lorentz force and |$F_{\rm C} = (\hat{z} \times \vec{v})_{\phi }$| is the ϕ component of the Coriolis force. The Match Coefficient ranges between −1 (forces are equal in magnitude and opposite in direction) and +1 (forces are equal in magnitude and direction). Note that if FL and FC are both zero then the expression for MLC gives an indefinite number. If this occurs at a given point, we assign the integrand a value of 1 since this corresponds to a perfect match. The more standard normalized correlation measure is defined as:
\begin{eqnarray} M_{{\rm LC}}^{{\rm stand}} &=& \frac{\int F_{\rm L} \ * \ F_{\rm C} \,\rm{d}V}{\sqrt{\int F_{\rm L}^2 \,\rm{d}V * \int F_{\rm C}^2 \,\rm{d}V} }. \end{eqnarray}
(26)
Figure 8.

Force balance: for (a) FFSF-21, (b) FFNS-15, (c) FTSF-15, (d) FTNS-24, (e) FTSF-20 and (f) FTNS-33 models. The left- and right-hand plots show the time averaged axisymmetric ϕ components of the Lorentz and Coriolis forces, respectively. The Match Coefficient for the Lorentz and Coriolis forces is also given. For comparison, the more standard measure of normalized correlation coefficients (eq. 26) for the models are (a) 0.7, (b) 0.52, (c) 0.72, (d) 0.44, (e) 0.15 and (f) −0.09. Note that the reported values are calculated by computing the instantaneous match coefficients and correlations and then averaging those values over time.

We also report these values in the Fig. 8 caption. We prefer the Match Coefficient defined in eq. (25) since it computes point-by-point normalization rather than an overall normalization.

As seen in Fig. 8, the models in the dipolar regime have a much higher degree of magnetostrophic balance than the models in the multipolar regime. This figure only shows the degree of correlation of the ϕ-averaged forces, that is those that accelerate a zonal flow, whereas our Match Coefficient measures a more detailed balance including the non-axisymmetric components. Since the multipolar models are dominated by their non-axisymmetric components we examined the force balance due only to the non-axisymmetric components (not shown here) and found that they were also not in magnetostrophic balance. It is evident that the strong, large scale fields of the dipolar models are a result of a stronger degree of magnetostrophic force balance in these models. In contrast, the multipolar models also produce intense magnetic fields due to the presence of strong Lorentz forces, but the Lorentz force is not balanced by the Coriolis force (and hence not in significant magnetostrophic balance). In the next section, we examine the power budget to investigate the force balances in these models.

3.4 Axisymmetric power budget

Similar to the procedure adopted by Aubert (2005), we multiply the axisymmetric ϕ component of the momentum equation by the ϕ component of the velocity field to get
\begin{eqnarray} \frac{1}{2} \frac{\partial }{\partial t} \bar{v_{\phi }}^2 \!&=&\! \overline{[\skew3\vec{J} \times \skew3\vec{B}]_{\phi }}\bar{v_{\phi }} -\bar{v_s}\bar{v_{\phi }} - {\rm Ro}_M \overline{[(\vec{v} \cdot \nabla ) \vec{v}]_{\phi }} \bar{v_{\phi }} + E \overline{[\nabla ^2 \vec{v}]_{\phi }} \bar{v_{\phi }} ,\!\!\!\nonumber\\ \end{eqnarray}
(27)
where an overbar denotes the axisymmetric component. Fig. 9 plots the time average of the four terms on the right-hand side of eq. (27), namely, the Lorentz, Coriolis, inertial and viscous terms, respectively. For the models in the dipolar regime, the main power source is the Coriolis term whereas the Lorentz and viscous terms act as the main sinks. In contrast, for the models in the multipolar regime, the main power source is the inertial term (Reynolds stresses) and the Lorentz and viscous terms act as sinks. In this case the Coriolis term acts as both a source and a sink in different regions of the outer core. The source is strongest at the location of the equatorial jet and inside the tangent cylinder in the dipolar models. In the multipolar models, the source is strongest outside the tangent cylinder in columns corresponding to where strong zonal flows are seen in Fig. 6. It is evident that the Lorentz force is balanced by the Coriolis in the dipolar dynamo models whereas it is balanced by the inertial force in the multipolar dynamo models. There is no strong dependence on the boundary conditions. For comparison, we plot the Coriolis, inertial and viscous terms of the corresponding non-magnetic convection models in Fig. 10. It is evident that the inertial and viscous terms in the convection models are more similar to the magnetic multipolar models demonstrating that the magnetic fields are not significantly altering the zonal flows in these cases.
Figure 9.

Axisymmetric power budget: for FFSF-21, FFNS-15, FTSF-15, FTNS-24, FTSF-20 and FTNS-33 dynamo models from the first to last row, respectively. The first, second, third and fourth columns show the time averaged axisymmetric ϕ components of the Lorentz, Coriolis, inertial and viscous terms, respectively. In these plots, red (blue) denotes sources (sinks) of power.

Figure 10.

Axisymmetric power budget: for the FFSFC-16, FFNSC-14, FTSFC-09, FTNSC-08, FTSFC-13 and FTNSC-15 convection models with no magnetic field from the first to last row, respectively. The first, second, and third columns show the time averaged axisymmetric ϕ components of the Coriolis, inertial and viscous terms, respectively. In these plots, red (blue) denotes sources (sinks) of power.

The force balance can also be quantified by comparing the total forces over a volume integral in a spherical shell as employed by Soderlund et al. (2012). Fig. 11 shows the force integrals and their respective ratios as a function of Nu for the different boundary conditions employed. In these plots, we fix the parameters to E = 6 × 10− 5, RoM = 2 × 10− 5, qk = 1. In this parameter space the Coriolis force dominates all the models irrespective of the thermal or velocity boundary conditions. At Nu ≥ 8, the Lorentz force is balanced by the inertial force, that is the ratio of the Lorentz to inertial force is approximately 1. At Nu < 8, the Lorentz force does not seem to be balanced by either the Coriolis force or the inertial force. These results are in contrast to our Matching Coefficient study. This could be due to the fact that the Matching Coefficient only considers the ϕ component of the forces, or because the total force magnitude does not determine the force balance, and instead, a point-by-point comparison needs to be made (as is done with the Match Coefficient).

Figure 11.

The rms force integrals: for models with E = 6 × 10−5, RoM = 2 × 10−5 and qk = 1. (a) FFSF, (b) FFNS, (c) FTSF and (d) FTNS. The top plots are the time averaged volume integrals of the rms Lorentz, Coriolis, inertial and viscous forces. The bottom plots are the time averaged ratios of the Lorentz force to the Coriolis and inertial forces, and the ratio of inertial force to Coriolis and viscous forces.

In Fig. 12, we only consider the ϕ component in the volume integrals. The degree of Lorentz–Coriolis balance seen in Fig. 8 is still not evident here, suggesting that total force magnitudes do not fully represent the force balances. Using Match Coefficients can better predict the force balance and the dipolarity of the field.

Figure 12.

The rms force integrals of the ϕ component: for models with E = 6 × 10−5, RoM = 2 × 10−5 and qk = 1. (a) FFSF, (b) FFNS, (c) FTSF and (d) FTNS. The top plots are the time averaged volume integrals of the ϕ component rms Lorentz, Coriolis, inertial and viscous forces. The bottom plots are the time averaged ratios of the ϕ component of the Lorentz force to the Coriolis and inertial forces, and the ratio of inertial force to Coriolis and viscous forces.

3.5 Generation mechanisms

Since the dipolar and multipolar models have different force balances, we investigate if they have different dynamo generation mechanisms. Poloidal magnetic field is generated by the α effect, through shearing and twisting of toroidal magnetic field lines in a helical flow. Toroidal magnetic field can be generated by either the α effect or by the ω effect, through shearing of field lines by large scale differential rotation. These different mechanisms can be seen mathematically by decomposing the generation term of the magnetic induction equation (eq. 2) as follows:
\begin{eqnarray} \nabla \times (\vec{v} \times \skew3\vec{B}) &=& \underbrace{(\skew3\vec{B} \cdot \nabla )\vec{v}}_{{\rm Creation}} - \underbrace{(\vec{v} \cdot \nabla )\skew3\vec{B}}_{{\rm Advection}} \end{eqnarray}
(28)
\begin{eqnarray} &&{(\skew3\vec{B} \cdot \nabla ) \vec{v} =}\nonumber\\ && \underbrace{\overbrace{(\skew3\vec{B}_T \cdot \nabla ) \vec{v}_P}^{{\rm Creates} \ \skew3\vec{B}_P} + \overbrace{(\skew3\vec{B}_P \cdot \nabla ) \vec{v}_P + (\skew3\vec{B}_P \cdot \nabla ) \vec{v}_T^{{\rm Naxi}} + (\skew3\vec{B}_T \cdot \nabla ) \vec{v}_T^{Naxi}}^{{\rm Creates} \ \skew3\vec{B}_T}}_{\alpha \ {\rm effect}} \nonumber \\ &&\quad +\, \underbrace{\overbrace{(\skew3\vec{B}_P \cdot \nabla ) \vec{v}_T^{{\rm Axi}} + (\skew3\vec{B}_T \cdot \nabla ) \vec{v}_T^{{\rm Axi}}}^{{\rm Creates} \ \skew3\vec{B}_T}}_{\omega \ {\rm effect}} , \end{eqnarray}
(29)
where Naxi and Axi refer to the non-axisymmetric and axisymmetric components, respectively. From eq. (29) we see that only the α effect creates |$\skew3\vec{B}_P$|⁠, hence, if we want to categorize the dynamo type, it is sufficient to compare the |$\skew3\vec{B}_T$| generation components.
In Olson et al. (1999) and Christensen & Olson (1998), the authors determined that their dynamos were α2 type, that is helical motions are responsible for both the poloidal and toroidal field generation. Using the same method as these studies with identical model parameters and boundary conditions, we obtain the same results. Here we determine the dynamo type for all our models by calculating the Match Coefficient of both α and ω effects with the total generation term for the toroidal field:
\begin{eqnarray} M_{\alpha } &=& \int \frac{F_{\alpha } \ * \ F_{\skew3\vec{B}_T}}{\left(F_{\alpha }^2 \ + \ F_{\skew3\vec{B}_T}^2\right)/2} \,\rm{d}V \end{eqnarray}
(30)
\begin{eqnarray} M_{\omega } &=& \int \frac{F_{\omega } \ * \ F_{\skew3\vec{B}_T}}{\left(F_{\omega }^2 \ + \ F_{\skew3\vec{B}_T}^2\right)/2} \,\rm{d}V, \end{eqnarray}
(31)
where Fα is the α effect that creates the axisymmetric ϕ component of |$\skew3\vec{B}_T$| (i.e. the second, third and fourth terms on the right-hand side of eq. 29), Fω is the ω effect that creates the axisymmetric ϕ component of |$\skew3\vec{B}_T$| (i.e. the fifth and sixth terms on the right-hand side of eq. 29), and |$F_{\skew3\vec{B}_T}$| is the sum of α and ω effects that create the axisymmetric ϕ component of |$\skew3\vec{B}_T$| (i.e. the sum of terms two through six on the right-hand side of eq. 29).

Fig. 13 presents these Match Coefficients as a function of Rol for different boundary conditions at fixed control parameters of E = 6 × 10− 5, RoM = 2 × 10− 5, qk = 1. At larger values of Rol, all of the models produce α2 dynamos. Intuitively, this is expected since a larger Rol implies a larger Ra and hence a stronger convective forcing. As Rol is decreased, the models with no-slip boundary conditions show a clear trend of increasing importance of ω effect and decreasing importance of α effect in generating the toroidal field. So the dynamo type in the no-slip models varies smoothly from αω to α2ω to α2 as Rol increases. The results for stress-free models are not as simple. At high Rol, the behaviour is similar to the no-slip models, but after Rol decreases to the point where the Match Coefficients of both effects are similar, the monotonic trends do not continue. There is large scatter in the Match Coefficient trends at lower Rol and it does not appear that the ω effect ever significantly dominates over the α effect. The stress-free dynamos, therefore, have either α2 or α2ω mechanisms.

Figure 13.

Match coefficient of α and ω effects versus local Rossby number: for models with E = 6 × 10− 5, RoM = 2 × 10− 5 and qk = 1. (a) FFSF, (b) FFNS, (c) FTSF and (d) FTNS. The green dotted line marks the transition of α effect dominance over ω effect.

Fig. 14 shows details of the field generation mechanisms. The time averaged toroidal magnetic field, creation of |$\skew3\vec{B}_T$|⁠, and the creation of |$\skew3\vec{B}_T$| due to α and ω effects (terms in eq. 29) for models in order of increasing Rol, namely, FTNS-18, FTNS-22 and FTNS-32 are presented. It is evident that model FTNS-18 is an αω type dynamo, model FTNS-22 is an α2ω type dynamo and model FTNS-32 is an α2 type dynamo.

Figure 14.

Dynamo creation mechanism: for models FTNS-18, FTNS-22 and FTNS-32 from first to last row, respectively. The first column shows the time averaged axisymmetric toroidal magnetic field. The second, third and fourth columns show the time averaged axisymmetric ϕ component of the total creation term of the toroidal magnetic field, toroidal magnetic field created by α effect and toroidal magnetic field created by the ω effect, respectively.

These results imply that axially dipolar dominated dynamos can be generated by α2, α2ω and αω mechanisms. In our models, all of the multipolar dynamos are generated by an α2 mechanism. This appears counterintuitive since the zonal flows in the high Rol models are stronger (see Fig. 6), however, as demonstrated in Fig. 9, the source of these strong zonal flows is the inertial force rather than the Coriolis force. It therefore appears that zonal flows must result from dynamos whose dominant force balance is between the Lorentz and Coriolis forces (rather than the Lorentz and inertial forces) to generate an axially dipolar dominated field in our models. These results also imply that the dominant Lorentz–Coriolis balance can occur with either generation mechanism for the toroidal field.

4 CONCLUSIONS

In this paper, we investigated the influence of different thermal and velocity boundary conditions on numerical dynamo models. We found that the resulting magnetic field morphologies most strongly depend on the local Rossby number (similar to previous work), but there are some secondary influences on the dipolarity of the field by the boundary conditions. Namely, at a given Rol, fixed temperature models are more dipolar than fixed flux models and stress-free models are more dipolar than no-slip models.

A single diffusivity-free scaling law for heat transport does seem feasible for all the dipolar models studied, although there is some variation among the different boundary conditions that may be important when extrapolating results of scaling laws to Earth-like parameters. The influence of different boundary conditions is particularly evident on Nu versus Ra/Rac scalings. Our multipolar models appear to obey different scaling laws from the dipolar models. Further investigation of scaling laws for multipolar models would benefit from a larger range in Ra, but this is numerically difficult presently.

The different scalings for multipolar and dipolar models is directly related to the different force balances in these models. Even though the magnetic energies of the dipolar and multipolar models were similar, dipolar models had a stronger degree of Lorentz–Coriolis balance whereas the multipolar dynamos possessed a dominant Lorentz-inertia balance. We demonstrated that plotting the axisymmetric Lorentz and Coriolis forces and estimating the Match Coefficient is a better way to establish the force balance than solely calculating the Elsasser number or comparing volume integrals of the forces.

We also demonstrated that a specific force balance does not imply a specific generation mechanism. A dominant Lorentz–Coriolis balance can be generated by either an α2, α2ω or an αω mechanism and result in dipolar dominated fields. In contrast, our multipolar fields, with dominant Lorentz-inertia balance, are all α2 type.

The authors thank U. Christensen and an anonymous reviewer for useful comments that improved this manuscript. S. Stanley acknowledges funding for this project by the Natural Sciences and Engineering Research Council of Canada (NSERC). Computations were performed on the gpc supercomputer at the SciNet HPC Consortium. SciNet is funded by the Canada Foundation for Innovation under the auspices of Compute Canada, the Government of Ontario, Ontario Research Fund – Research Excellence, and the University of Toronto.

REFERENCES

Aubert
J.
Steady zonal flows in spherical shell dynamos
J. Fluid Mech.
2005
, vol. 
542
 (pg. 
53
-
67
)
Aubert
J.
Labrosse
S.
Poitou
C.
Modelling the palaeo-evolution of the geodynamo
Geophys. J. Int.
2009
, vol. 
179
 (pg. 
1414
-
1428
)
Christensen
U.R.
Dynamo scaling laws and applications to the planets
Space Sci. Rev.
2010
, vol. 
152
 
1–4
(pg. 
565
-
590
)
Christensen
U.R.
Aubert
J.
Scaling properties of convection-driven dynamos in rotating spherical shells and application to planetary magnetic fields
Geophys. J. Int.
2006
, vol. 
166
 (pg. 
97
-
114
)
Christensen
U.R.
Olson
P.
A dynamo model interpretation of geomagnetic field structures
Geophys. Res. Lett.
1998
, vol. 
25
 
10
(pg. 
1565
-
1568
)
Christensen
U.R.
Tilgner
A.
Power requirement of the geodynamo from ohmic losses in numerical and laboratory dynamos
Nature
2004
, vol. 
429
 
6988
(pg. 
169
-
171
)
Christensen
U.R.
Olson
P.
Glatzmaier
G.A.
Numerical modelling of the geodynamo: a systematic parameter study
Geophys. J. Int.
1999
, vol. 
138
 (pg. 
393
-
409
)
Christensen
U.R.
Holzwarth
V.
Reiners
A.
Energy flux determines magnetic field strength of planets and stars
Nature
2009
, vol. 
457
 (pg. 
167
-
169
)
Dharmaraj
G.
Stanley
S.
Effect of inner core conductivity on planetary dynamo models
Phys. Earth planet. Inter.
2012
, vol. 
212–213
  
doi:10.1016/j.pepi.2012.09.003
Gastine
T.
Duarte
L.
Wicht
J.
Dipolar versus multipolar dynamos: the influence of the background density stratification
A&A
2012
, vol. 
546
 (pg. 
A19:1
-
13
)
Glatzmaier
G.A.
Roberts
P.H.
An anelastic evolutionary geodynamo simulation driven by compositional and thermal convection
Physica D
1996
, vol. 
97
 (pg. 
81
-
94
)
Hori
K.
Wicht
J.
Christensen
U.R.
The effect of thermal boundary conditions on dynamos driven by internal heating
Phys. Earth planet. Inter.
2010
, vol. 
182
 (pg. 
85
-
97
)
Kuang
W.
Force balances and convective state in the earth's core
Phys. Earth planet. Inter.
1999
, vol. 
116
 (pg. 
65
-
79
)
Kuang
W.
Bloxham
J.
An Earth-like numerical dynamo model
Nature
1997
, vol. 
389
 (pg. 
371
-
374
)
Kuang
W.
Bloxham
J.
Numerical modeling of magnetohydrodynamic convection in a rapidly rotating spherical shell: weak and strong field dynamo action
J. Comp. Phys.
1999
, vol. 
153
 (pg. 
51
-
81
)
Kutzner
C.
Christensen
U.R.
Effects of driving mechanisms in geodynamo models
Geophys. Res. Lett.
2000
, vol. 
27
 (pg. 
29
-
32
)
Kutzner
C.
Christensen
U.R.
From stable dipolar towards reversing numerical dynamos
Phys. Earth planet. Inter.
2002
, vol. 
131
 (pg. 
29
-
45
)
Miyagoshi
T.
Kageyama
A.
Sato
T.
Zonal flow formation in the Earth's core
Nature
2010
, vol. 
463
 (pg. 
793
-
796
)
Miyagoshi
T.
Kageyama
A.
Sato
T.
Formation of sheet plumes, current coils, and helical magnetic fields in a spherical magnetohydrodynamic dynamo
Phys. Plasmas
2011
, vol. 
18
 pg. 
072901
  
doi:10.1063/1.3603822
Olson
P.
Christensen
U.
Dipole moment scaling for convection-driven planetary dynamos
Earth planet. Sci. Lett.
2006
, vol. 
250
 (pg. 
561
-
571
)
Olson
P.
Christensen
U.
Glatzmaier
G.A.
Numerical modeling of the geodynamo: mechanisms of field generation and equilibration
J. geophys. Res.
1999
, vol. 
104
 (pg. 
10 383
-
10 404
)
Sakuraba
A.
Roberts
P.H.
Generation of a strong magnetic field using uniform heat flux at the surface of the core
Nat. Geosci.
2009
, vol. 
2
 (pg. 
802
-
805
)
Sasaki
Y.
Takehiro
S.
Kuramoto
K.
Hayashi
Y.
Weak-field dynamo emerging in a rotating spherical shell with stress-free top and no-slip bottom boundaries
Phys. Earth planet. Inter.
2011
, vol. 
188
 
3–4
(pg. 
203
-
213
)
Schrinner
M.
Petitdemange
L.
Dormy
E.
Dipole collapse and dynamo waves in global direct numerical simulations
Astrophys. J.
2012
, vol. 
752
 pg. 
121
  
doi:10.1088/0004-637X/752/2/121
Simitev
R.D.
Busse
F.H.
Bistability and hysteresis of dipolar dynamos generated by turbulent convection in rotating spherical shells
EPL
2009
, vol. 
85
  
doi:10.1209/0295-5075/85/19001
Simitev
R.D.
Busse
F.H.
Bistable attractors in a model of convection-driven spherical dynamos
Phys. Scr.
2012
, vol. 
86
  
doi:10.1088/0031-8949/86/01/018409
Soderlund
K.M.
King
E.M.
Aurnou
J.M.
The influence of magnetic fields in planetary dynamo models
Earth planet. Sci. Lett.
2012
, vol. 
333–334
 (pg. 
9
-
20
)
Stelzer
Z.
Jackson
A.
Extracting scaling laws from numerical dynamo models
Geophys. J. Int.
2013
, vol. 
2013
 
Takahashi
F.
Matsushima
M.
Honkura
Y.
Simulations of a quasi-Taylor state geomagnetic field including polarity reversals on the Earth simulator
Science
2005
, vol. 
309
 (pg. 
459
-
461
)
Takahashi
F.
Matsushima
M.
Honkura
Y.
Scale variability in convection-driven MHD dynamos at low Ekman number
Phys. Earth planet. Inter.
2008
, vol. 
167
 (pg. 
168
-
178
)
Yadav
R.K.
Gastine
T.
Christensen
U.R.
Scaling laws in spherical shell dynamos with free-slip boundaries
Icarus
2013
, vol. 
225
 (pg. 
185
-
193
)

SUPPORTING INFORMATION

Additional Supporting Information may be found in the online version of this article:

Table S1. Dynamo model results.

Table S2. Results of convection models.

Table S3. Models without hyperdiffusion (Supplementary Data).

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the article.

Supplementary data