# Numerical dynamos with outer boundary heat flux inferred from probabilistic tomography-consequences for latitudinal distribution of magnetic flux

Hagay Amit, Frédéric Deschamps, Gael Choblet

► To cite this version:

Hagay Amit, Frédéric Deschamps, Gael Choblet. Numerical dynamos with outer boundary heat flux inferred from probabilistic tomography-consequences for latitudinal distribution of magnetic flux. Geophysical Journal International, 2015, 203 (2), pp.840-855. 10.1093/gji/ggv332. hal-02319531

HAL Id: hal-02319531

https://hal.science/hal-02319531

Submitted on 19 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Numerical dynamos with outer boundary heat flux inferred from probabilistic tomography—consequences for latitudinal distribution of magnetic flux

Hagay Amit,<fnref n="1" /> Frédéric Deschamps<fnref n="2" /> and Gaël Choblet<fnref n="1" />

Accepted 2015 August 10. Received 2015 August 7; in original form 2014 December 1

## SUMMARY

Mantle control on the geodynamo is often simulated using numerical dynamos with imposed outer boundary heat flux inferred from lower mantle tomography, assuming that seismic and thermal anomalies in the lowermost mantle are highly correlated. However, non-thermal effects might perturb this idealized linear seismic–thermal mapping. Here we use a probabilistic tomography model to isolate the thermal part of the seismic anomaly in order to impose a more realistic core–mantle boundary heat flux pattern on the outer boundary of numerical dynamo simulations. We demonstrate that on time average these dynamo models have more low-latitude convective and magnetic activity than corresponding models with conventional tomographic heat flux. In addition, the low-latitude magnetic flux and kinetic energy contributions are more time-dependent in the dynamo models with a probabilistic tomography heat flux, and thus may recover the observed latitudinal distribution of geomagnetic flux on the core–mantle boundary, which we propose as a morphological criterion for Earth-like dynamo models.

**Key words:** Dynamo: theories and simulations; Magnetic field; Heat flow; Mantle processes; Core, outer core and inner core; Seismic tomography.

## 1 INTRODUCTION

Various observed features of the geodynamo have been recovered by imposing heterogeneous heat flux patterns on the outer boundary of numerical dynamo simulations (e.g. Olson & Christensen 2002; Christensen & Olson 2003; Gubbins et al. 2007; Aubert et al. 2007, 2008; Davies et al. 2008). In these simulations, the prescribed heat flux patterns are proportional to seismic tomography models of the lowermost mantle, under the assumption that seismic and thermal anomalies in this region are highly correlated. However, observational and numerical studies of the lower mantle suggest that the structures and dynamics there are far more complex, with other (non-thermal) contributions perturbing the idealized linear seismic–thermal mapping (e.g. Trampert et al. 2004; Lay et al. 2008; Nakagawa & Tackley 2008; Bull et al. 2009; Tackley 2011; Deschamps et al. 2012). Clearly an appropriate core–mantle boundary (CMB) heat flux model is required in order to model adequately core–mantle thermal interactions and the resulting persistent geodynamo features.

The primary observational constraint for any dynamo model, with an homogeneous or heterogeneous outer boundary pattern, is the morphology of the geomagnetic field on the CMB. Fig. 1 shows two models of the radial geomagnetic field on the CMB inverted from surface observatories and dedicated satellites (Jackson et al. 2000; Olsen et al. 2010). The field morphology is characterized

by regions of concentrated geomagnetic flux (Gubbins 2003; Christensen et al. 2010). Two intense normal polarity flux patches appear at high latitudes of each hemisphere. Because a spherical harmonic degree and order 2 pattern dominates lower-mantle seismic anomalies (e.g. Masters et al. 2000), it has been suggested that these two pairs of intense high-latitude geomagnetic flux patches are maintained by the lower-mantle heterogeneity (Bloxham 2002; Olson & Christensen 2002; Gubbins et al. 2007; Aubert et al. 2008). The underlying concept is that while the tangent cylinder determines the latitudes of these patches (Olson et al. 1999), the mantle heterogeneity dictates the statistically preferred longitudes (Gubbins 2003), although on short timescales the patches may oscillate between one preferred location to another (Amit et al. 2010). Jackson (2003) noted the presence of intense flux patches at the equatorial region as well, in particular below Africa and the Atlantic Ocean (Fig. 1). The existence of these low-latitude features is less well understood, especially because numerical dynamos are often characterized by zonal equatorial upwelling (Aubert 2005; Amit & Olson 2006) that disperses field lines at low latitudes.

To assess quantitatively what dynamo model can be considered as Earth-like, Christensen et al. (2010) defined four morphological criteria based on observed geomagnetic field models: Relative axial dipole power, equatorial symmetry, zonality and flux concentration. While all these criteria are global, objective and reflect robust geomagnetic field morphological features, the latitudinal distribution

<footnote n="1">¹CNRS UMR 6112, Université de Nantes, Laboratoire de Planétologie et de Géodynamique, 2 rue de la Houssinière, Nantes, F-44000, France.</footnote>

<footnote>E-mail: Hagay.Amit@univ-nantes.fr</footnote>

<footnote n="2">²Institute of Earth Sciences, Academia Sinica, 128 Academia Road Sec. 2, Nangang, Taipei 11529, Taiwan</footnote>
<figure><img src="image_2.png" /><figcaption>**Figure 1.** Models of the radial geomagnetic field on the CMB expanded until spherical harmonic degree 14. (a) From *CHAOS-4* for the year 2010 (Olsen *et al.* 2010); (b) From *gufm1* for the year 1980 (Jackson *et al.* 2000); (c,d) Corresponding zonal $B_r$ (solid) and rms of non-zonal $B_r$ (dashed) versus latitude in mT. For *CHAOS-4* 2010 the low/high-latitude magnetic flux ratios are $B_r^{rms}$ = 0.805 and $B_r^{max}$ = 0.913; For *gufm1* 1980 $B_r^{rms}$ = 0.735 and $B_r^{max}$ = 0.787. When *CHAOS-4* 2010 is truncated at $ℓ_{max}$ = 10, the ratios are $B_r^{rms}$ = 0.766 and $B_r^{max}$ = 1.037.</figcaption></figure>

of flux is not represented, in particular the two pairs of high-latitude intense flux patches observed near the interception of the tangent cylinder with the CMB (Fig. 1) that are supported by theory and laboratory experiments (Aurnou *et al.* 2003). In addition, the low-latitude intense flux patches (again Fig. 1) are also robust. Their mobility accounts for the most striking secular variation signatures (Finlay & Jackson 2003; Aubert *et al.* 2013; Livermore *et al.* 2013). The specific latitudes where these geomagnetic flux patches are found, that is, near the tangent cylinder and at low latitudes, are not constrained by the criteria of Christensen *et al.* (2010). Therefore, a dynamo model might score high without recovering this particular and important observed latitudinal distribution of geomagnetic flux.

Another observational constraint on dynamo models may come from the geomagnetic secular variation. Aubert *et al.* (2013) recovered the westward drift of the modern field with their dynamo model. On long timescales, Davies & Constable (2014) proposed that the axial dipole frequency spectrum of a dynamo model should obey the same power law fit as that of the paleomagnetic dipole over the past 2 Myr (Ziegler *et al.* 2011). Davies & Constable (2014) applied their temporal variability criterion in combination with the criteria of Christensen *et al.* (2010) on dynamo models with either homogeneous or tomographic CMB heat flux to assess how Earth-like is a dynamo model.

The dynamo models investigated by Christensen *et al.* (2010) have an homogeneous outer boundary condition. An important ingredient that may affect geomagnetic field morphology is boundary control. In particular, mantle control on the geodynamo may help explain the presence of intense geomagnetic flux patches at low latitudes. The much explored tomographic pattern has little potential in that respect because it is dominated by a $Y_2$ term with no significant latitudinal preference. In order to incorporate possible compositional and small-scale effects in the lower mantle (e.g. Bull *et al.* 2009), Amit & Choblet (2012) constructed a CMB heat flux pattern with hot thermal ridges at the edges of large-scale hot and

cold lowermost mantle structures. They found that such hot ridges may yield local upwelling and a flow barrier that result in localized concentration of magnetic flux, most effectively at low latitudes. However, the deviation of the CMB heat flux pattern proposed by Amit & Choblet (2012) from that of conventional tomography is mild, and its effects are moderate as well. A stronger impact on the resulting field morphology could be obtained by a more detailed description of lowermost mantle thermo-chemical heterogeneities and consequently a stronger departure from conventional tomography pattern.

Low shear wave velocity provinces (LLSVPs), which are the prominent features in conventional mantle tomography models and give rise to the dominant $Y_2$ distribution, are unlikely purely thermal structures. They are better explained if a compositional component is also present (Trampert *et al.* 2004; Deschamps *et al.* 2012; Mosca *et al.* 2012), despite the detailed nature of compositional anomalies being still debated. One hypothesis is that LLSVPs are enriched in iron, for example, due to early partial differentiation of the mantle, in which case these would be reservoirs of dense, undegassed material. Another hypothesis is that LLSVPs consist of piles of recycled material. However, equation of state modelling of the mantle aggregates indicates that unless they are very hot, LLSVPs are unlikely to consist exclusively of recycled oceanic crust (Deschamps *et al.* 2012). Hybrid scenarios have recently been proposed in which LLSVPs of primitive origin may be periodically re-fed with small amounts of recycled oceanic crust reaching the bottom of the mantle (Tackley 2011; Li *et al.* 2014). Post-perovskite, if present, may also play a key role in the structure of the lowermost mantle. Its large Clapeyron slope implies that it is not stable in hot regions. The transition from perovskite to post-perovskite is a good candidate for the D'' discontinuity, which is observed locally in the deep mantle. Cobden & Thomas (2013) showed that it would explain D'' reflections in regions where P and S waveforms have different polarities, but that other explanations (including chemical changes
and anisotropy) may also be possible in regions where both *P* and *S* waveforms have positive polarities.

By taking into account seismic normal modes, Ishii & Tromp (1999) provided an independent constraint on mantle density distributions, which can be used to resolve different contributions to seismic anomalies. In probabilistic tomography studies (Trampert et al. 2004; Mosca et al. 2012), anomalies in density, shear wave velocity and bulk-sound (or compressional) velocity are first derived from a seismic data set that includes normal mode data using a Monte-Carlo (neighbourhood algorithm) inversion. A key advantage of the Monte-Carlo inversion is that it provides full probability density functions (PDFs) for anomalies in density and seismic velocities at each node of the model. In a second step, these PDFs of seismic anomalies are combined with appropriate equation of state modelling of the mantle aggregate and mineral physics database (including uncertainties in these parameters), giving direct access to PDFs of thermochemical anomalies, including the thermal anomalies. It turns out that the patterns of thermal anomalies in probabilistic tomography models are remarkably different than their seismic counterparts. In the study by Trampert et al. (2004), the thermal heterogeneity resembles a $Y_2^1$ term, in contrast to the $Y_2^2$ dominance of the seismic anomaly. Mosca et al. (2012) used updated mineral physics constraints and confined their vertical resolution to a thinner layer above the CMB. The most striking feature of their thermal anomaly pattern is the concentration of intense structures at low latitudes. Note, however, that at the bottom of the mantle the standard deviation in temperature anomalies is at most locations larger than the mean temperature anomaly.

Applying probabilistic tomography is problematic because it contains only even degrees, a consequence of normal modes being used as an independent constraint on the density structure. To resolve this issue we design a simple procedure to incorporate odd degrees by mapping the ratio between thermal and seismic (in our case, shear wave velocity) anomalies from probabilistic tomography images. This mapping is then combined with a conventional seismic tomography pattern, which contains odd degrees, to obtain a re-scaled thermal anomaly that includes a complete spectrum up to degree 6.

Here we use the probabilistic tomography model of Mosca et al. (2012) in conjunction with the conventional tomography model of Masters et al. (2000) to obtain a more realistic CMB heat flux pattern, which we impose on the outer boundary of numerical dynamo simulations. In Section 2, we present the formalism to construct the CMB heat flux model from probabilistic tomography, and we attempt to assess its validity using some synthetic tests. The numerical dynamo models and the statistical measures of the latitudinal magnetic and kinetic distributions are described in Section 3. In Section 4, we demonstrate that these probabilistic tomography dynamo models have more low-latitude convective and magnetic activity than conventional tomography dynamo models and thus may improve the recovery of the observed latitudinal distribution of geomagnetic flux on the CMB. We discuss and summarize our main findings in Section 5.

# 2 A CMB HEAT FLUX MODEL FROM PROBABILISTIC TOMOGRAPHY

## 2.1 Numerical procedure

Probabilistic tomography provides distributions of shear wave velocity ($\delta\ln V_s$) and temperature ($\delta T$) anomalies, which we use to construct the CMB heat flux pattern. Interestingly, the PDFs provided

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Probability density functions (PDFs) for (a) seismic velocity anomalies ($\delta\ln V_s^e$) and (b) temperature anomalies ($\delta T^e$) in the lowermost mantle (sampling the layer 2750–2891 km) at two different locations (from Mosca et al. 2012). These PDFs approximate well Gaussian distributions, here defined with a mean and a standard deviation. (c) PDFs of the ratio $\delta T^e/\delta\ln V_s^e$.</figcaption></figure>

by probabilistic tomography approximate well Gaussian distributions (Mosca et al. 2012), and can thus be defined by a mean and a standard deviation. Figs 2(a) and (b) represent distributions of $\delta\ln V_s$ and $\delta T$ in the lowermost mantle (sampling the layer 2750–2891 km) at two different locations. A limitation of the probabilistic tomography models published so far is that only even spherical harmonic degrees up to 6 are included (Trampert et al. 2004; Mosca et al. 2012). While the truncation at degree 6 is unlikely to modify the large-scale structure, the missing odd degrees may introduce a substantial bias in the maps of thermal anomalies. To overcome this problem we map the ratio between shear wave velocity (even degrees) and thermal anomalies (even degrees) using probabilistic tomography. We then apply this mapping to a conventional lower mantle seismic tomography model (even and odd degrees) to obtain
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Construction of a CMB heat flux model based on probabilistic tomography. (a) Probabilistic tomography model of shear wave velocity anomalies $\delta \ln V_s^e$ (maximum likelihood) in the lowermost mantle (2750–2891 km) until spherical harmonic degree 6, but with even degrees only (Mosca et al. 2012). (b) Probabilistic tomography model of temperature anomalies $\delta T^e$ (maximum likelihood) at the lowermost mantle until spherical harmonic degree 6, again with even degrees only (Mosca et al. 2012). (c) Conventional tomography model of seismic anomalies $\delta \ln V_s^C$ at the lowermost mantle until spherical harmonic degree 6 (Masters et al. 2000). (d) Re-scaled model of temperature anomalies $\delta T^P$ at the lowermost mantle until spherical harmonic degree 6 with both even and odd degrees, constructed from (a), (b) and (c) (see the text for details).</figcaption></figure>

a thermal anomaly distribution with both even and odd degrees. Mathematically we write

$$ \delta T^P = \frac{\delta T^e}{\delta \ln V_s^e} \delta \ln V_s^C, \qquad (1) $$

where the superscript ‘C’ denotes a conventional tomography model. Here we use the model of Masters et al. (2000) truncated at spherical harmonic degree 6 (Fig. 3c). The superscript ‘e’ denotes a probabilistic tomography model with only even degrees content. Physically eq. (1) implicitly assumes that the even degree’s thermal contribution to the even seismic velocity available from probabilistic tomography is comparable to the unknown odd degree’s thermal contribution to the odd seismic part (see synthetic tests of this approximation in Section 2.2). We use the model of Mosca et al. (2012) including degrees 2, 4 and 6 for the seismic (Fig. 3a) and thermal (Fig. 3b) anomalies. The resulting rescaled thermal anomaly including all degrees until $\ell_{max} = 6$ is denoted by superscript ‘P’ (Fig. 3d).

The scaling procedure defined in eq. (1) may trigger numerical instabilities at regions where $\delta \ln V_s^e$ from probabilistic tomography is low. To avoid this difficulty, at each node of the probabilistic tomography model we calculated the ratios of the full PDFs of $\delta T^e$ and $\delta \ln V_s^e$ (rather than the ratios of their mean values), resulting in a full PDF of the ratio $\delta T^e / \delta \ln V_s^e$. Because PDFs of $\delta \ln V_s^e$ and $\delta T^e$

are Gaussian distributions to a good approximation, the PDF of the ratio $\delta T^e / \delta \ln V_s^e$ obeys the following analytical integration:

$$ g_{\delta T^e / \delta \ln V_s^e}(y) = \int_{-\infty}^{\infty} |x| g_{\delta T^e}(yx) g_{\delta \ln V_s^e}(x) dx, \qquad (2) $$

where $g_{\delta T^e}$ and $g_{\delta \ln V_s^e}$ are the PDFs of $\delta T^e$ and $\delta \ln V_s^e$ respectively, $x$ is a dummy variable representing the input distributions (in our case $\delta T^e$ and $\delta \ln V_s^e$) and $y$ is a dummy variable representing the output distribution (in our case $\delta T^e / \delta \ln V_s^e$). Note that because we assume Gaussian distributions for both $\delta T^e$ and $\delta \ln V_s^e$, the output ratio $\delta T^e / \delta \ln V_s^e$ follows a Cauchy distribution. In practice, we calculated this integral within bounds $x_1$ to $x_2$ defined by

$$ x_1 = \min(\delta T^e - 5\sigma_{\delta T^e}, \delta \ln V_s^e - 5\sigma_{\delta \ln V_s^e}) \qquad (3) $$

and

$$ x_2 = \max(\delta T^e + 5\sigma_{\delta T^e}, \delta \ln V_s^e + 5\sigma_{\delta \ln V_s^e}), \qquad (4) $$

where $\sigma_{\delta T^e}$ and $\sigma_{\delta \ln V_s^e}$ are the corresponding standard deviations.

Fig. 2(c) shows resulting distributions of $\delta T^e / \delta \ln V_s^e$ at two locations, one where $\delta \ln V_s^e$ is close to zero (blue curve; 60°N, 255°E; $\delta \ln V_s^e = 0.02$ per cent), and another where $\delta \ln V_s^e$ is large (red curve; 15°N, 135°E; $\delta \ln V_s^e = 1.4$ per cent). To calculate the rescaled temperature anomalies, we then simply multiply the velocity anomalies from Masters et al. (2000) (for which uncertainties
<figure><img src="image_5.png" /><figcaption>**Figure 4.** Correlation between $\delta T^p$ and $\delta T$ as a function of the amplitude of the standard deviation. $\delta T^p$ was calculated based on our numerical procedure with the following synthetic setup: $\delta T \propto \delta \ln V_s^C$, $\delta \ln V_s^e$ is the even part of $\delta \ln V_s^C$ and $\delta T^e$ is the even part of $\delta T$. The dashed horizontal line denotes the reference correlation between $\delta T^e$ and $\delta T$.</figcaption></figure>

are not available) by the values of $\delta T^e / \delta \ln V_s^e$ at the maximum of the PDFs (i.e. maximum likelihood).

Finally, inferring the pattern of CMB heat flux anomaly $\delta q$ from the lower-mantle thermal anomaly $\delta T$ is rather straightforward. Note that because temperature differences on the mantle side of the CMB are by far larger than on the core side, $\delta q$ is antiproportional to $\delta T$, that is, $\delta q^C \propto -\delta T^C$ and $\delta q^P \propto -\delta T^P$, so that cold (warm) lower-mantle regions correspond to large (low) CMB heat flux respectively. The amplitude of the CMB heat flux anomaly is determined by the control parameter $q^*$ (see Section 3.1).

## 2.2 Synthetic tests

It is important to assess the validity of the procedure for obtaining the CMB heat flux model. As a first idealized test we consider the limit of perfectly correlated seismic and thermal anomalies, that is, $\delta T \propto \delta \ln V_s^C$ (where $\delta T$ is the true thermal anomaly). We then take the even part of $\delta \ln V_s^C$ to get $\delta \ln V_s^e$, and likewise $\delta T^e$ is assumed to be the even part of $\delta T$. If the standard deviation is zero, eq. (1) is solved analytically, in this case giving $\delta T^p = \delta T$, that is, perfect recovery. Next we apply eq. (1) to our procedure with variable standard deviation amplitudes to calculate $\delta T^p$. The correlation between our calculated $\delta T^p$ and the true $\delta T$ is shown in Fig. 4 as a function of the standard deviation amplitude. For comparison, an horizontal line denotes the reference correlation between $\delta T^e$ and $\delta T$. Only for the largest standard deviation amplitude $\delta T^v$ is better correlated with $\delta T$ than the model produced by our procedure; for all other standard deviation amplitudes, our procedure provides larger correlations between $\delta T^p$ and $\delta T$.

The above idealized test has two unrealistic features. First, the full thermal and seismic anomalies are not expected to be perfectly correlated (which is precisely the basis of this study). Second, although the probabilistic tomography (purely even) seismic anomaly and the even part of the conventional tomography model have a significant spatial correlation of 0.5, these two patterns are not identical. The next test accounts for these more realistic conditions. Here we use $\delta \ln V_s^e$ and $\delta T^e$ from the probabilistic tomography model, and thus $\delta T^p$ is also the same as that used in this paper. We then construct a synthetic true thermal anomaly $\delta T$ by summing $\delta T^v$ plus a random odd part. Fig. 5(a) shows the correlation between $\delta T^p$ and $\delta T$ (red) compared to the correlation between $\delta T^e$ and $\delta T$ (green).

Overall, the two correlations are comparable with slight advantage to $\delta T^v$ (Fig. 5b) and larger dispersion for $\delta T^p$. When the relative odd part is smaller, both correlations are higher, as expected. The vertical lines in Fig. 5 denote the relative odd power in $\delta \ln V_s^C$ and $\delta T^p$. If these vertical lines represent the relative odd power in $\delta T$ of Earth's lower mantle, Fig. 5(a) suggest maximal correlations of ~0.5 between $\delta T^p$ and $\delta T$.

Two points are important to emphasize concerning the second test. First, although the correlation between $\delta T^p$ and $\delta T$ is comparable to that of $\delta T^e$ to $\delta T$, our approach is conceptually better because it includes all degrees (odd and even). Second, a random odd part of $\delta T$ is a worst case scenario; eq. (1) is expected to predict an odd part which is closer to the real temperature than just a random distribution. This motivates additional synthetic tests with possibly some geophysical merit. In Fig. 5(c), the synthetic $\delta T$ is constructed from the sum of $\delta T^e$ plus a non-random odd part pattern with a random amplitude. When the pattern of the odd part of $\delta \ln V_s^C$ is added (solid), corresponding to thermochemical even part and purely thermal odd part, $\delta T^e$ always exhibits larger correlations, possibly due to the inconsistent physical origins of the even and odd parts in this synthetic case. However, the difference between the two correlations diminishes with increasing odd power (Fig. 5d). When the pattern of the odd part of $\delta T^p$ is added (dashed), the correlation of $\delta T^p$ becomes superior for relative odd power larger than 0.15. The peak correlation between $\delta T^p$ and $\delta T$ of 0.9 is reached in the plausible range of relative odd power of 0.15-0.45 delimited by the vertical lines.

The latter tests (Fig. 5c) yield significantly larger correlations than the test with random odd part (Fig. 5a). Moreover, in Fig. 5(c) the correlations degrade with increasing odd power at a slower rate compared to in Fig. 5(a). This demonstrates that introducing plausibility in the synthetic $\delta T$ provides superior correlations, that is, our procedure is radically advantageous in computing the missing odd part compared to purely random addition.

# 3 NUMERICAL DYNAMO MODELS

## 3.1 Setup

We use 3-D self-consistent numerical dynamos in a spherical shell. We implement the simulation MAGIC (Wicht 2002), originally coded by Gary Glatzmaier. The control parameters are the heat flux based (modified) Rayleigh number $Ra = \alpha g_0 q_0 D^4 / \kappa k \nu$, the Ekman number $E = \nu / \Omega D^2$, the Prandtl number $Pr = \nu / \kappa$, the magnetic Prandtl number $Pm = \nu / \lambda$ and a heat flux anomaly amplitude $q^* = (q_{max} - q_{min}) / 2q_0$, where $\alpha$ is thermal expansion coefficient, $g_0$ is gravity on the outer boundary, $q_0$ is mean heat flux across the outer boundary, $D$ is shell thickness, $k$ is thermal conductivity, $\kappa$ is thermal diffusivity, $\nu$ is kinematic viscosity, $\Omega$ is rotation rate and $\lambda$ is magnetic diffusivity. We use $D$, $D^2 / \nu$, $Dq_0 / k$ and $\sqrt{\rho \mu_0 \lambda \Omega}$ to scale length, time, co-density (Braginsky & Roberts 1995) and magnetic field, respectively, where $\rho$ is the fluid density and $\mu_0$ is the permeability of free space. The inner to outer radii ratio is 0.35. The boundary conditions are rigid for the velocity and insulating for the magnetic field. The inner boundary has fixed co-density, and the outer boundary has a prescribed heat flux pattern. The co-density source/sink term in the energy equation is set to zero corresponding to generic thermochemical convection (Aubert et al. 2009). The grid includes 41 radial, 192 longitudinal and 96 latitudinal points, and in spectral space the fields were expanded up to degree and order 64. The dynamo model's control parameters are given in Table 1.
<figure><img src="image_6.png" /><figcaption>**Figure 5.** Left: correlation between $\delta T^P$ and $\delta T$ (red) and between $\delta T^e$ and $\delta T$ (green) for synthetic $\delta T$ as a function of the power of odd harmonics (degrees 1, 3 and 5) normalized by the global power (degrees 1–6). $\delta T^P$ and $\delta T^e$ are the ones used in this study. Right: difference between the two correlations. Top: 60 000 cases with $\delta T$ constructed from $\delta T^e$ plus a random odd part. Bottom: $\delta T$ constructed from $\delta T^e$ plus the pattern of the odd part of $\delta T^P$ with a random amplitude (solid), and $\delta T$ constructed from $\delta T^e$ plus the pattern of the odd part of $\delta T^P$ with a random amplitude (dashed). The blue and purple vertical lines denote the relative power of odd harmonics in $\delta T^P$ and $\delta T^e$, respectively.</figcaption></figure>

For more details on the governing equations and numerical method see Christensen & Wicht (2007).

All models are in a parametric regime in which dipole polarity reversals are absent, and overall the magnetic field morphology in a snapshot is dominated by an axial dipole component. Two patterns are imposed as outer boundary conditions. The first is based on a conventional lower mantle tomography model (Masters *et al.* 2000) and is termed ‘C’ in Table 1. The second is derived from probabilistic tomography (see above) and is termed ‘P’.

## 3.2 Statistical measures of latitudinal distribution

The latitudinal distributions of two quantities are calculated. Primarily the magnetic flux on the outer boundary is examined for comparison with observed geomagnetic field models on the CMB. In addition, the latitudinal distribution of kinetic energy is quantified both at the top of the free stream and throughout the volume of the shell. Because the core flow generates and maintains the magnetic field, it may be expected that similar latitudinal distributions

will be found for the magnetic flux and the kinetic energy. However, because the particular field-flow interactions might be different at low- and high latitudes (Finlay & Amit 2011), it is worth-while examining both quantities.

To quantify the ratio of low- to high-latitude distributions of magnetic field and flow, we divide the spherical surface and spherical shell into two regions. High latitudes are defined above a critical latitude in both hemispheres, and low latitudes are defined below this critical latitude. We arbitrarily choose as a critical latitude 30° because it conveniently separates equal spherical areas of low- and high latitudes. In addition, since high-latitude intense flux patches have branches that often extend to mid-latitudes, the choice of a relatively low-latitude of 30° as a separation avoids interpretation of such branches as intense low-latitude flux. Overall, the main objective is not to interpret the absolute value of these ratios, but rather to compare ratios of pairs of dynamo models with different outer boundary heat flux patterns and otherwise identical dynamo internal parameters, so in that respect the choice of a separation latitude is probably not critical. We nevertheless examine the sensitivity of our results to the choice of the critical latitude.
<table><caption>**Table 1.** Dynamo model setup and results. In cases denoted by 'C' the imposed outer boundary heat flux pattern is from conventional tomography (Masters *et al. 2000*), in the 'P' cases the pattern is derived from probabilistic tomography. The magnitude of the heat flux heterogeneity is $q^*$. The internal dynamo control parameters *Ra* and *Pm* are given. In all cases $E = 3 \times 10^{-4}$ and $Pr = 1$. Also given is the magnetic Reynolds number *Rm* computed from the time-average kinetic energy in the volume of the shell. The low/high-latitude magnetic flux rms/maxima ratios (time-average and its standard deviation) are $B_{r_{lh}}^{rms}$ and $B_{r_{lh}}^{max}$, respectively. The corresponding rms ratios for exclusively normal and reversed flux are $B_{N_{lh}}^{rms}$ and $B_{R_{lh}}^{rms}$, respectively. The low/high-latitude kinetic energy ratios at the volume of the shell and just below the Ekman boundary layer (time-average and its standard deviation) are $KE_{lh}^{V}$ and $KE_{lh}^{S}$, respectively. The low/high-latitude magnetic flux ratios are also given for the geomagnetic field models *gufm1* (Jackson *et al. 2000*) in the period 1840–1990, *CHAOS-4* in the period 2000–2010 (Olsen *et al. 2010*) as well as the archeomagnetic field model *CALS3k.3* over the past three millennia (Korte *et al. 2009*). Note that the $B_{R_{lh}}^{rms}$ value for the archeomagnetic field model is non-applicable because in some snapshots high latitudes are completely deprived of reversed flux and the ratio becomes infinity.</caption><thead><tr><th>Case</th><th>q*</th><th>Ra</th><th>Pm</th><th>Rm</th><th>B<sub>r</sub><sup>rms</sup><sub>lh</sub></th><th>B<sub>r</sub><sup>max</sup><sub>lh</sub></th><th>B<sub>N</sub><sup>rms</sup><sub>lh</sub></th><th>B<sub>R</sub><sup>rms</sup><sub>lh</sub></th><th>KE<sub>lh</sub><sup>V</sup></th><th>KE<sub>lh</sub><sup>S</sup></th></tr></thead><tbody><tr><td>C1</td><td>0.7</td><td>1 × 10<sup>6</sup></td><td>3</td><td>128</td><td>0.439 ± 0.05</td><td>0.305 ± 0.08</td><td>0.352 ± 0.04</td><td>0.722 ± 0.11</td><td>0.335 ± 0.04</td><td>0.290 ± 0.04</td></tr><tr><td>P1</td><td>0.7</td><td>1 × 10<sup>6</sup></td><td>3</td><td>122</td><td>0.524 ± 0.08</td><td>0.408 ± 0.13</td><td>0.422 ± 0.06</td><td>0.909 ± 0.20</td><td>0.338 ± 0.06</td><td>0.314 ± 0.06</td></tr><tr><td>C2</td><td>1.34</td><td>1 × 10<sup>6</sup></td><td>3</td><td>132</td><td>0.425 ± 0.05</td><td>0.283 ± 0.07</td><td>0.342 ± 0.04</td><td>0.682 ± 0.12</td><td>0.354 ± 0.07</td><td>0.287 ± 0.06</td></tr><tr><td>P2</td><td>1.34</td><td>1 × 10<sup>6</sup></td><td>3</td><td>123</td><td>0.563 ± 0.09</td><td>0.461 ± 0.16</td><td>0.461 ± 0.07</td><td>0.968 ± 0.21</td><td>0.368 ± 0.07</td><td>0.361 ± 0.08</td></tr><tr><td>C3</td><td>0.7</td><td>3 × 10<sup>6</sup></td><td>3</td><td>236</td><td>0.641 ± 0.11</td><td>0.539 ± 0.18</td><td>0.504 ± 0.09</td><td>0.991 ± 0.20</td><td>0.244 ± 0.05</td><td>0.206 ± 0.04</td></tr><tr><td>P3</td><td>0.7</td><td>3 × 10<sup>6</sup></td><td>3</td><td>236</td><td>0.702 ± 0.11</td><td>0.641 ± 0.19</td><td>0.546 ± 0.09</td><td>1.147 ± 0.23</td><td>0.270 ± 0.04</td><td>0.228 ± 0.03</td></tr><tr><td>C4</td><td>0.7</td><td>1 × 10<sup>6</sup></td><td>10</td><td>382</td><td>0.499 ± 0.07</td><td>0.405 ± 0.13</td><td>0.401 ± 0.05</td><td>0.741 ± 0.13</td><td>0.417 ± 0.06</td><td>0.420 ± 0.06</td></tr><tr><td>P4</td><td>0.7</td><td>1 × 10<sup>6</sup></td><td>10</td><td>376</td><td>0.575 ± 0.09</td><td>0.533 ± 0.17</td><td>0.467 ± 0.07</td><td>0.844 ± 0.16</td><td>0.463 ± 0.06</td><td>0.505 ± 0.09</td></tr><tr><td>gufm1</td><td>na</td><td>na</td><td>na</td><td>na</td><td>0.675 ± 0.04</td><td>0.725 ± 0.05</td><td>0.658 ± 0.04</td><td>2.130 ± 0.60</td><td>na</td><td>na</td></tr><tr><td>CHAOS-4</td><td>na</td><td>na</td><td>na</td><td>na</td><td>0.797 ± 0.01</td><td>0.882 ± 0.02</td><td>0.772 ± 0.01</td><td>1.539 ± 0.03</td><td>na</td><td>na</td></tr><tr><td>CALS3k.3</td><td>na</td><td>na</td><td>na</td><td>na</td><td>0.470 ± 0.07</td><td>0.648 ± 0.17</td><td>0.428 ± 0.07</td><td>na</td><td>na</td><td>na</td></tr></tbody></table>

Two statistical measures are defined to assess the ratio of low- to high-latitude magnetic flux on the outer boundary at any snapshot. The first is an rms ratio:

$$B_{r_{lh}}^{\text{rms}} = \sqrt{\frac{S_h}{S_l} \frac{\int_0^{2\pi} \int_{\theta_c}^{\pi-\theta_c} B_r^2 dS}{\int_0^{2\pi} \int_0^{\theta_c} B_r^2 dS + \int_0^{2\pi} \int_{\pi-\theta_c}^{\pi} B_r^2 dS}}, \quad (5)$$

where $B_r$ is the radial magnetic field on the outer boundary, $\theta_c$ is the critical northern co-latitude that separates low- and high latitudes, $S_l$ and $S_h$ are the low- and high-latitude spherical surface areas, respectively (for a critical latitude of 30°, $S_l = S_h$) and $dS = r^2 \sin \theta d\phi d\theta$ is a spherical surface increment on the outer boundary. For the magnetic flux, we use $r = r_o$ where $r_o$ is the outer radius of the spherical shell. The rms ratio (5) ranges between zero (no flux at low latitudes) to infinity (all flux concentrated at low latitudes). In eq. (5), the first integral is over longitude and the second over co-latitude. The second statistical measure for the magnetic flux latitudinal distribution is a maxima ratio of the absolute radial field:

$$B_{r_{lh}}^{\text{max}} = \frac{\max(|B_r|^{\pi-\theta_c})}{\max(|B_r|^{\theta_c}, \ |B_r|^{\pi-\theta_c})}. \quad (6)$$

In eq. (6) the ranges represent co-latitudes for all longitudes.

The ratio of low- to high-latitude unsigned flux (5) reflects the latitudinal distribution of magnetic flux, but it is insensitive to the latitudinal distribution of polarity, which is important for defining an Earth-like field morphology. We therefore examine the ratios of low- to high-latitude magnetic flux applied exclusively to normal and reversed flux:

$$B_{N_{lh}}^{\text{rms}} = \sqrt{\frac{S_h}{S_l} \frac{\int_0^{2\pi} \int_{\theta_c}^{\pi-\theta_c} B_N^2 dS}{\int_0^{2\pi} \int_0^{\theta_c} B_N^2 dS + \int_0^{2\pi} \int_{\pi-\theta_c}^{\pi} B_N^2 dS}} \quad (7)$$

$$B_{R_{lh}}^{\text{rms}} = \sqrt{\frac{S_h}{S_l} \frac{\int_0^{2\pi} \int_{\theta_c}^{\pi-\theta_c} B_R^2 dS}{\int_0^{2\pi} \int_0^{\theta_c} B_R^2 dS + \int_0^{2\pi} \int_{\pi-\theta_c}^{\pi} B_R^2 dS}} \quad (8)$$

where the normal flux is defined by

$$B_N = \begin{cases} B_r & \text{if } \theta < \pi/2 \text{ and } B_r < 0 \text{ or } \theta > \pi/2 \text{ and } B_r > 0 \\ 0 & \text{otherwise} \end{cases} \quad (9)$$

and the reversed flux is defined by

$$B_R = \begin{cases} B_r & \text{if } \theta < \pi/2 \text{ and } B_r > 0 \text{ or } \theta > \pi/2 \text{ and } B_r < 0 \\ 0 & \text{otherwise.} \end{cases} \quad (10)$$

To quantify the ratio of low- to high-latitude of the flow activity we integrate the kinetic energy KE. Because the effects of rapid rotation produce flows that are quasi-invariant in the direction of the rotation axis (e.g. Busse 1970), we separate here low- and high latitudes in the shell volume using cylindrical coordinates. The low and high-latitude kinetic energies are integrated over all grid points with $s > s_c$ and $s < s_c$, respectively, where $s$ is the radial distance from the rotation axis in cylindrical coordinates and $s_c = r_o \sin \theta_c$ is the critical distance of a co-centric cylinder that intersects the CMB at co-latitude $\theta_c$. We calculate the kinetic energy ratio in the volume of the spherical shell as follows:

$$KE_{lh}^{V} = \frac{V_h}{V_l} \frac{\int_{0}^{2\pi} \int_{\theta_c}^{\pi-\theta_c} \int_{r_o}^{r_c} KE \, dV}{\int_{0}^{2\pi} \int_{\theta_c}^{\pi-\theta_c} \int_{r_l}^{r_c} KE \, dV + \int_{0}^{2\pi} \int_{\theta_c}^{\pi-\theta_c} \int_{r_l}^{r_c} KE \, dV}, \quad (11)$$

where $V_h$ and $V_l$ are the volumes of the high- and low latitudes, respectively, $r_c = s_c/\sin \theta$ and $dV = r^2 \sin \theta d\phi d\theta dr$ is a spherical volume increment. In eq. (11), the integrals are performed over longitude, co-latitude and the radial direction, respectively. Finally, as in eq. (5), we also calculate the ratio of low- to high-latitude of the kinetic energy at the top of the free stream just below the upper Ekman boundary layer:

$$KE_{lh}^{S} = \frac{S_h}{S_l} \frac{\int_{0}^{2\pi} \int_{\theta_c}^{\pi-\theta_c} KE \, dS}{\int_{0}^{2\pi} \int_{\theta_c}^{\pi-\theta_c} KE \, dS + \int_{0}^{2\pi} \int_{\pi-\theta_c}^{\pi} KE \, dS} \quad (12)$$

In eq. (12), the integrals are performed over longitude and co-latitude, respectively.
<figure><img src="image_7.png" /><figcaption>**Figure 6.** Non-dimensional CMB heat flux anomaly patterns imposed on the numerical dynamo models. (a) Based on conventional tomography (Masters *et al.* 2000). (b) Inferred from probabilistic tomography. Red/blue corresponds to positive/negative anomalies, respectively.</figcaption></figure>

## 4 RESULTS

Fig. 6 shows the two CMB heat flux patterns used and compared in this study. The well-explored conventional tomography pattern (Fig. 6a and cases ‘C’ in Table 1) is characterized by a $Y_2^2$ dominance with two large CMB heat flux structures below the Americas and east Asia/Oceania. These features are very large-scale, practically extending across all latitudes. No significant latitudinal preference is observed in the conventional tomography pattern. The pattern derived from probabilistic tomography (Fig. 6b and cases ‘P’ in Table 1) is characterized by smaller scales. Three positive anomalies appear below Central America, Indonesia and east Australia, that is, most of the large CMB heat flux is concentrated at low latitudes.

Figs 7 and 8 compare the convective activities in two snapshots from cases C2 and P2. These two cases have identical control parameters and differ only by their CMB heat flux patterns (see Table 1). The specific snapshots chosen here have close to average latitudinal distribution of kinetic energy, reflecting the typical convective state in these models. In both cases, the zonal reduced co-density $T(\theta, r) - \bar{T}(r)$ (where $\bar{T}(r)$ is the mean over a spherical surface at radial distance $r$) is characterized by downwellings at the edge of the tangent cylinder and equatorial rising plumes, but in case C2 this plume is intense (Fig. 7a), whereas in case P2 it is weaker (Fig. 8a). The meridional flow in case C2 is strongly concentrated at high latitudes (see streamline concentration inside the tangent cylinder in Fig. 7b), whereas in case P2 it is much more balanced with comparable convective activity at mid- and even low latitudes (Fig. 8b). The radial vorticity maps at the top of the free stream show in both cases north-south elongated structures, evidence for the dominance

of rotational effects in both dynamo models. However, in case C2 significant radial vorticity structures are concentrated at the polar regions (Fig. 7c), whereas again in case P2 the radial vorticity latitudinal distribution is much more balanced with much less activity inside the tangent cylinder (Fig. 8c). Overall, in these two dynamo models the low- to high-latitude ratio of kinetic energy in the spherical shell volume is comparable. However, at the top of the free stream the ratio is significantly larger, by 26 per cent, when a CMB heat flux pattern derived from probabilistic tomography is imposed (Table 1).

Next we compare the resulting magnetic fields on the outer boundary (Figs 9 and 10). For comparison with observations, we apply the statistical measures to the radial fields low-pass filtered around spherical harmonic degree 10 (Figs 9b and 10b). Again we show snapshots that reflect the typical latitudinal distribution, in this case of magnetic flux. In both cases, high-latitude flux patches prevail, but in case P2 there is an additional intense magnetic flux patch below the equatorial Indian Ocean (Figs. 10a and b), where anomalously positive heat flux is present (Fig. 6b). The more balanced latitudinal distribution of magnetic flux in case P2 is also apparent in the profiles of zonal radial field and rms of non-zonal radial field (compare Figs 9c and 10c). In case C2, larger values of the zonal radial field characterize high latitudes, while low latitudes are roughly comparable in cases C2 and P2. Overall, the low- to high-latitude ratios in case P2 are larger than in case C2 by 32 per cent and 63 per cent based on rms and maxima measures respectively (Table 1).

The statistical measures (5) and (6) consider the unsigned magnetic flux. The geomagnetic field contains prominent reversed flux
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Images of a snapshot from case C2. (a) Non-dimensional zonal reduced co-density $T(\theta, r) - \bar{T}(r)$. (b) Non-dimensional zonal flow (colours) and meridional streamlines (contours, solid/dashed denote anticlockwise/clockwise circulation, respectively). (c) Non-dimensional radial vorticity at the top of the free stream (red/blue denote anticlockwise/clockwise circulation, respectively). This snapshot exhibits low/high-latitude kinetic energy ratios similar to the time-average for this case: $KE_{\theta\theta}^V = 0.333$ and $KE_{\theta\theta}^S = 0.267$.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 8.** As in Fig. 7 for case P2 with identical control parameters to those of C2 but a non-dimensional CMB heat flux pattern inferred from probabilistic tomography. This snapshot exhibits low/high-latitude kinetic energy ratios similar to the time-average for this case: $KE_{\theta\theta}^V = 0.405$ and $KE_{\theta\theta}^S = 0.416$.</figcaption></figure>

patches (Fig. 1). In order to assess whether a dynamo model is Earth-like or not, it is important to compare the latitudinal distribution of signed flux as well. In Table 1, we calculated the low- to high-latitude ratios of exclusively normal and reversed flux (eqs (7) and (8), respectively). The dynamo models differ in their distribution of signed flux in a similar manner as their difference in terms of unsigned flux ratios, with systematically larger relative low-latitude contributions in the P cases, both normal and reversed, for each C and P pair. While the unsigned flux ratios $B_{\theta\theta}^{\text{rms}}$ in the large Ra models reach that of gufm1, the normal flux ratio $B_{\theta\theta}^{\text{rms}}$ of the dynamo models is somewhat lower than the observed.

We tested the sensitivity of our results to the choice of separation latitude. The results in Table 1 were obtained with critical latitude of 30°. Using critical latitude of 35°, we find very small changes to the time-average ratios. For example, with critical latitude of 35° we obtain in case C2 $B_{\theta\theta}^{\text{rms}} = 0.414 \pm 0.05$ and $B_{\theta\theta}^{\text{max}} = 0.297 \pm 0.08$, in case P2 $B_{\theta\theta}^{\text{rms}} = 0.559 \pm 0.09$ and $B_{\theta\theta}^{\text{max}} = 0.499 \pm 0.18$, and for gufm1 $B_{\theta\theta}^{\text{rms}} = 0.670 \pm 0.04$ and $B_{\theta\theta}^{\text{max}} = 0.743 \pm 0.03$. Compari-

son with the ratios obtained with critical latitude of 30° (Table 1) reveals variations of ~1–3 per cent and ~2–8 per cent for $B_{\theta\theta}^{\text{rms}}$ and $B_{\theta\theta}^{\text{max}}$, respectively. We conclude that the sensitivity of our results to the choice of $\theta_c$ is very weak in both the dynamo models and the geomagnetic field models, in particular concerning the rms ratio.

Next we examine the time-series of the rms and maxima ratios of low- to high-latitude magnetic flux ratios (eqs (5) and (6), respectively) for cases C2 and P2 (Figs. 11a and b, respectively). As expected, the pointwise maxima ratio is a noisier quantity characterized by larger time-dependence (see standard deviations in Table 1). In case P2, the time-average values are larger than in case C2, that is, low-latitude relative magnetic flux is larger in case P2. Moreover, the temporal fluctuations are by far larger in case P2. Fig. 12 shows the radial magnetic field in a snapshot of case P2 in which the relative intensity of low-latitude flux patches is high. In this snapshot, negative/positive flux patches still dominate high latitudes of the northern/southern hemisphere, respectively, so that the axial dipole is still prominent (see the solid line in Fig. 12b). However, the flux
<figure><img src="image_10.png" /><figcaption>**Figure 9.** A snapshot of the non-dimensional radial magnetic field on the outer boundary from case C2, non-filtered (a) and filtered (b). (c) Non-dimensional zonal $B_r$ (solid) and rms of non-dimensional non-zonal $B_r$ (dashed) versus latitude. The snapshot in (b) exhibits low/high-latitude magnetic flux ratios similar to the time-average for this case: $B_{rh}^{\text{rms}} = 0.470$ and $B_{rh}^{\text{max}} = 0.281$.</figcaption></figure>

patches at low latitudes are comparably intense (see the dashed line in Fig. 12b), especially below the Atlantic hemisphere (Fig. 12a).

Table 1 gives the kinetic energy and magnetic flux latitudinal distribution ratios for dynamo models with different control parameters. Increasing $q^*$ with a conventional tomography heat flux pattern (Fig. 6a) results in slight decrease in relative low-latitude magnetic flux ratios and comparable kinetic energy ratios (compare cases C1 and C2), whereas with the pattern derived from probabilistic tomography (Fig. 6b) increasing $q^*$ gives significant increase in relative low-latitude magnetic flux and kinetic energy (compare cases P1 and P2). Increasing $Ra$ increases low-latitude magnetic flux and decreases low-latitude kinetic energy in both cases (compare cases C1 and C3, and compare cases P1 and P3). Increasing $Pm$ which results in larger magnetic Reynolds number $Rm$ increases low-latitude

dynamics in both cases (compare cases C1 and C4, and compare cases P1 and P4). Most importantly, in all four pairs the P cases exhibit more balanced latitudinal distribution of flow and magnetic field, that is, relatively stronger magnetic flux and convective activity at low latitudes.

It is worth comparing the low- to high-latitude magnetic flux ratios in our dynamo models with those of the observed geomagnetic field. Fig. 13 shows the magnetic rms and maxima ratios over the historical era based on the field model *gufm1* of Jackson *et al.* (2000). The increase of the ratios with time might reflect an increase in the resolution of the models due to more available geomagnetic measurements. Furthermore, over the short period of 150 yr the ratios drift rather than oscillate (as in Fig. 11), so their standard deviations are probably not very meaningful for uncertainty estimates.
<figure><img src="image_11.png" /><figcaption>**Figure 10.** As in Fig. 9 for the parallel case P2 with a non-dimensional CMB heat flux pattern inferred from probabilistic tomography. The snapshot in (b) exhibits low/high-latitude magnetic flux ratios similar to the time-average for this case: $B_r^{\text{rms}} = 0.517$ and $B_{\ell h}^{\text{max}} = 0.480$.</figcaption></figure>

The observed ratios of ~0.6–0.7 are comparable to the time-average values of case P3 (Table 1). Accounting for temporal fluctuations, any of the P cases may occasionally recover and even exceed these observed values (e.g. Figs 11b and 12).

Table 1 compares the latitudinal distributions of geomagnetic flux among the historical field model *gufm1* (Jackson *et al.* 2000) in the period 1840–1990, the modern field model *CHAOS-4* in the period 2000–2010 (Olsen *et al.* 2010) and the archeomagnetic field model *CALS3k.3* over the past three millennia (Korte *et al.* 2009). The low- to high-latitude unsigned flux rms and maxima ratios as well as the normal flux ratio are all largest in *CHAOS-4* and lowest in *CALS3k.3*, suggesting that increasing data coverage at more recent models improves the resolution of low-latitude features in a more significant manner than at high latitudes.

Note that the reversed flux ratio $B_{R\ell h}^{\text{rms}}$ in *gufm1* and *CHAOS-4* is strongly biased by the intense positive flux patch below equatorial Atlantic which is considered reversed in the northern hemisphere. This problem can be circumvented by a topological algorithm that maps the magnetic equator (Terra-Nova *et al.* 2015), but such an algorithm has not yet been implemented for dynamo models. In addition, for the archeomagnetic field model the reversed flux ratio is non-applicable because in some snapshots high latitudes are completely deprived of reversed flux and the ratio thus becomes infinity. Overall, the results for *CALS3k.3* should be treated with caution due to the model’s low spatial and temporal sampling. Nevertheless, these results provide an insight into millennial timescale variability and the effect of reduced spatial resolution.
<figure><img src="image_12.png" /><figcaption>**Figure 11.** Time-series of low/high-latitude magnetic flux ratios for cases C2 (a) and P2 (b). Black/red curves are rms/maxima ratios, respectively.</figcaption></figure>

## 5 DISCUSSION

To the best of our knowledge, this is the first time that a CMB heat flux pattern inferred from probabilistic tomography is imposed on the outer boundary of numerical dynamo simulations. The advantage of this approach is clearly the decomposition of the seismic anomaly to different origins, for our purposes the isolation of the lowermost mantle thermal anomaly. However, caution should be taken when constraining geodynamo models with such patterns due to imperfections of our procedure to incorporate odd degrees and due to the uncertainties associated with probabilistic tomography.

Our synthetic tests demonstrate that if a random odd thermal part is added to the probabilistic even thermal anomaly, our model performs comparably well while we benefit from a complete spectrum (Figs 5a and b); If a thermal part with a non-random (and geophysically plausible) pattern but random amplitude is added, our model may significantly improve the recovery of the true thermal anomaly (Figs 5c and d). Overall, these tests demonstrate that our procedure is far better than a purely random determination of the missing odd part. In addition, the correlation between our model and the real mantle temperature is expected to be high if the odd power is relatively low, as is expected for the Earth.

A more adequate test with both seismic and thermal synthetic patterns is more difficult to design because there is a physical relation

between these two quantities (as well as their standard deviations) prescribed by the mantle dynamics (Nakagawa & Tackley 2008). Unfortunately this relation is not trivial, thus considering a random synthetic thermal pattern in conjunction with a specific seismic pattern is dynamically inconsistent. Despite these uncertainties, our model seems geodynamically sensible because the high spatial correlation between the purely even thermal pattern of probabilistic tomography and our model is consistent with the similarity between the purely even seismic pattern of probabilistic tomography and the dominantly even pattern of global seismic tomography.

More rigorous tests may be performed based on the assessment of the relations among seismic, thermal and chemical anomalies, including the partitioning to even and odd parts, in a large set of ‘Earth-like’ mantle convection simulations. Such a study will improve the assessment of the validity of our procedure and may possibly lead to an improved procedure. However, even such an approach cannot be considered as a proof due to uncertainties in mantle convection modelling (e.g. Tackley 2012).

A measure of the uncertainties in the seismic shear wave velocity anomalies $\delta\ln V_s^e$ is given by the standard deviation in its (Gaussian) PDF derived from the neighbourhood algorithm and normal mode data. Uncertainties in the temperature anomalies $\delta T^e$ further account for uncertainties in the mineral physics data set and reference radial models of mantle temperature and composition. In
<figure><img src="image_13.png" /><figcaption>**Figure 12.** A snapshot of the non-dimensional radial magnetic field on the outer boundary from case P2, non-filtered (a). (b) Non-dimensional zonal $B_r$ (solid) and rms of non-dimensional non-zonal $B_r$ (dashed) versus latitude. This snapshot exhibits large low/high-latitude magnetic flux ratios of $B_r^{\text{rms}}_{1h} = 0.905$ and $B_r^{\text{max}}_{1h} = 1.010$.</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>**Figure 13** As in Fig. 11 for the geomagnetic field model *gufm1* for the period 1840–1990.</figcaption></figure>

Mosca *et al.* (2012), the rms in standard deviations of $\delta\ln V_s^e$ and $\delta T^e$ are 1.0 per cent and 340 K, respectively. These uncertainties add up when calculating the ratio $\delta\ln V_s^e/\delta T^e$. Overall, the rms in standard deviation of the rescaled temperature anomalies $\delta T^p$ is about 540 K.

Ideally, these uncertainties may be used to define a set of CMB heat flux patterns and associated numerical dynamo models. In practice, this would require an exaggerated computational time.

Apart from the above formal uncertainties, several theoretical and observational developments may alter the maps of thermal anomalies used in this study. Improved normal mode catalogues (Deuss *et al.* 2013), potentially including odd degrees, may lead to refined maps of density and seismic velocity anomalies, leading in turn to new maps of thermal and compositional anomalies. Small fractions of recycled oceanic crust may be incorporated in LLSVPs (Tackley 2012; Li *et al.* 2014), an effect which is not taken into account when inferring thermal anomalies from probabilistic tomography. Electrical conductivity, which is sensitive to both temperature and composition, is a promising tool to probe the thermochemical structure of the mantle. Combined with conventional or probabilistic tomography models, it may provide finer constraints on the temperature. Synthetic reconstructions of electrical conductivity from probabilistic tomography and available mineral physics suggest that a belt of high conductivity peaking in LLSVPs may be present in the lowermost mantle (Deschamps 2015). However, due to the lack of high-frequency magnetic field variation data over long periods (1 yr and higher), current models of electrical conductivity are limited to 1500 km depth (Kelbert *et al.* 2009; Semenov & Kuvshinov 2012). In addition, mineral physics measurements of electrical conductivity of mantle minerals at lower-mantle conditions remain sparse and should be completed. Finally, small-scale structures, in particular the ultra low velocity zones (ULVZs) observed by seismologists (Garnero & Helmberger 1995) may locally modify the temperature predicted by global tomography. The nature of ULVZs is still
debated. Hypotheses include pockets of partial melts (Williams &
Garnero 1996) and lenses of iron-rich post-perovskite (Mao et al.
2006), which would have different implications in terms of the ther-
mal distribution.

In addition, lateral variations in the lower mantle thermal con-
ductivity may affect the CMB heat flux. Recent ab initio calcu-
lations suggest that thermal conductivity of perovskite decreases
with both increasing temperature and iron content, but that these
corrections are small (Tang et al. 2014). If confirmed, this would
slightly enhance the heat flux heterogeneity amplitude derived from
the mantle thermo-chemical structure. Ammann et al. (2014) found
that the thermal conductivity of post-perovskite at lower mantle
conditions is about 50 per cent higher than that of perovskite, which
would further increase the amplitude of lateral variations in thermal
conductivity. Overall, a complete CMB heat flux modelling should
incorporate the lateral variability of the thermal conductivity of the
lowermost mantle in addition to that of the temperature.

While the seismic pattern of probabilistic tomography (Fig. 3a)
is rather similar to the conventional seismic tomography pattern
(Fig. 3c), the thermal anomaly is distinctive (Fig. 3d). The conven-
tional tomography pattern in Fig. 3(a) is characterized by two belts
of low seismic shear wave velocity, one along the Americas and the
other along east Asia, both extending through all latitudes. In con-
trast, the large heat flux structures of the probabilistic tomography
pattern in Fig. 3(d) are concentrated at low latitudes. If the latter
pattern indeed reflects the CMB thermal heterogeneity, previous
inferences based on conventional tomography should be revised, in
particular concerning the possible impact of mantle heterogeneity
on core dynamics at low latitudes.

We demonstrated that the probabilistic tomography dynamo mod-
els (P cases) indeed contain relatively stronger convective and mag-
netic activity at low latitudes than the conventional tomography
dynamo models (C cases). For a given set of internal dynamo pa-
rameters, the ratio of kinetic energy at low- to high-latitude is sys-
tematically larger in the P cases than in the C cases, both in the shell
volume and at the top of the shell. Likewise, the ratio of low- to
high-latitude magnetic flux (unsigned as well as exclusively normal
or reversed) on the outer boundary is systematically larger in the
P cases than in the C cases. In addition, the standard deviations of
these ratios are either comparable or larger in the P cases (Table 1),
facilitating occasional bursts of intense magnetic flux patches at low
latitudes that surpass the corresponding ratios in the observed geo-
magnetic field models (Fig. 13). These results persist when varying
Ra and Pm. When increasing the CMB heat flux heterogeneity am-
plitude q*, the magnetic flux ratio decreases for the C cases but sig-
nificantly increases for the P cases (Table 1), further demonstrating
the efficiency of the probabilistic pattern in producing low-latitude
intense flux patches as opposed to the tendency of the conventional
pattern to intensify the high-latitude flux patches. In addition, if
thermal conductivity decreases with increasing temperature (Tang
et al. 2014) and increases at cold mantle downwellings where post-
perovskite is expected to be present (Ammann et al. 2014), a larger
q* is expected for a given lower-mantle thermal distribution. In that
case, a lower q* C case should be compared with a larger q* P case,
thus enhancing the contrast between the two CMB heat flux models
in terms of low-latitude contributions to magnetic and convective
activities.

Dynamo models with homogeneous CMB heat flux exhibit downwellings at the edge of the tangent cylinder (Olson et al. 1999) and zonal equatorial upwelling (Aubert 2005; Amit & Olson 2006). In the P cases, the large heat flux at the equatorial region produces mantle-driven downwellings so the homogeneous dynamo equato-

rial upwelling is suppressed (Fig. 8), resulting in more convective activity and magnetic flux concentration at low latitudes. Such dynamics may explain the different latitudinal distribution of magnetic flux and convective activity for different dynamo internal control parameters (Table 1). For a given CMB heat flux (either conventional or probabilistic), increasing Ra strengthens high-latitude cyclones at the expense of the equatorial upwelling, thus weakens low-latitude convective activity. The weaker equatorial upwelling results in less efficient magnetic dispersion and stronger low-latitude magnetic flux. Increasing Pm also yields stronger low-latitude magnetic flux, suggesting that diffusion is more important at low latitudes of the dynamo models.

We emphasize that we do not claim that the dynamo models presented here are more Earth-like in terms of their field morphology than other models previously published (e.g. Christensen et al. 2010). To find such optimal Earth-like model, a much more extensive parametric study is required. Here we merely examine the possible impact of a CMB heat flux inferred from probabilistic tomography on the resulting magnetic field morphology, and point to its tendency to reproduce one observed geomagnetic feature, that is, low-latitude intense flux patches (Jackson 2003). The effect of the probabilistic tomography pattern on the magnetic field is compared with that of conventional tomography CMB heat flux, rather than searching for an optimal combination of internal control parameters. Overall, the observed latitudinal distribution of magnetic flux on the CMB (Fig. 13) may be combined with the criteria of Christensen et al. (2010) through a $\chi^2$ evaluation as a new criterion for Earth-like dynamo models.

We maintained moderate $q^*$ values in our dynamo models, thus in these models rotational effects still dominate in the form of nearly axially invariant flows, and the magnetic field is still dominated by an axial dipole and high-latitude intense flux patches. The boundary heterogeneity mildly modifies the background homogeneous dynamo, most notably by strengthening low-latitude magnetic and convective activities. However, some encouraging results using the conventional tomographic pattern might be lost using our alternative probabilistic tomography pattern. In particular, the P cases do not recover the preferential longitudes of high-latitude flux patches (Gubbins et al. 2007; Aubert et al. 2008), but the existence of such non-axisymmetric field features in the time-average paleomagnetic field (Johnson & Constable 1995; Kelly & Gubbins 1997) is controversial (Carlut & Courtillot 1997) while at intermediate timescales these patches were found mobile in tomographic dynamo models (Bloxham 2002; Amit et al. 2010) and in the archeomagnetic field (Korte & Holme 2010; Amit et al. 2011).

Efforts to isolate the thermal contribution to the lowermost man-
tle seismic anomalies are essential to better understand structures
and dynamical processes on both sides of the CMB. On the man-
tle side, it is crucial to understand the interplay between thermal,
compositional and mineralogical heterogeneities. On the core side,
it is essential to impose an optimal thermal outer boundary con-
dition on geodynamo simulations. Overall, progress in elucidating
lower-mantle dynamics may help improve models of CMB heat flux,
provide more realistic constraint on numerical dynamo models and
allow for better recovery of observed geodynamo related properties.

ACKNOWLEDGEMENTS

We thank two anonymous reviewers for constructive sugges-
tions that significantly improved the paper. FD was supported by
Academia Sinica grant 102-CDA-M02 for this research.
## REFERENCES

Amit, H. & Choblet, G., 2012. Mantle-driven geodynamo features - effects of compositional and narrow D'' anomalies, *Phys. Earth planet. Inter.*, **190-191**, 34-43.

Amit, H. & Olson, P., 2006. Time-average and time-dependent parts of core flow, *Phys. Earth planet. Inter.*, **155**, 120-139.

Amit, H., Aubert, J. & Hulot, G., 2010. Stationary, oscillating or drifting mantle-driven geomagnetic flux patches?, *J. geophys. Res.*, **115**, B07108, doi:10.1029/2009JB006542.

Amit, H., Korte, M., Aubert, J., Constable, C. & Hulot, G., 2011. The time-dependence of intense archeomagnetic flux patches, *J. geophys. Res.*, **116**, B12106, doi:10.1029/2011JB008538.

Ammann, M.W., Walker, A.M., Stackhouse, S., Wookey, J., Forte, A.M., Brodholt, J.P. & Dobson, D.P., 2014. Variation of thermal conductivity and heat flux at the Earth's core mantle boundary, *Earth planet. Sci. Lett.*, **390**, 175-185.

Aubert, J., 2005. Steady zonal flows in spherical shell fluid dynamos, *J. Fluid Mech.*, **542**, 53-67.

Aubert, J., Amit, H. & Hulot, G., 2007. Detecting thermal boundary control in surface flows from numerical dynamos, *Phys. Earth planet. Inter.*, **160**, 143-156.

Aubert, J., Amit, H., Hulot, G. & Olson, P., 2008. Thermo-chemical wind flows couple Earth's inner core growth to mantle heterogeneity, *Nature*, **454**, 758-761.

Aubert, J., Finlay, C.C. & Fournier, A., 2013. Bottom-up control of geomagnetic secular variation by the Earth's inner core, *Nature*, **502**, 219-223.

Aubert, J., Labrosse, S. & Poitou, C., 2009. Modelling the paleo-evolution of the geodynamo, *Geophys. J. Int.*, **179**, 1414-1428.

Aurnou, J., Andreadis, S., Zhu, L. & Olson, P., 2003. Experiments on convection in Earth's core tangent cylinder, *Earth Planet. Sci. Lett.*, **212**, 119-134.

Bloxham, J., 2002. Time-independent and time-dependent behaviour of high-latitude flux bundles at the core-mantle boundary, *Geophys. Res. Lett.*, **29**, doi:10.1029/2001gl014543.

Braginsky, S.I. & Roberts, P.H., 1995. Equations governing convection in Earth's core and the geodynamo, *Geophys. Astrophys. Fluid Dyn.*, **79**, 1-97.

Bull, A., McNamara, A. & Ritsema, J., 2009. Synthetic tomography of plume clusters and thermochemical piles, *Earth Planet. Sci. Lett.*, **278**, 152-162.

Busse, F.H., 1970. Thermal instabilities in rapidly rotating systems. *J. Fluid Mech.*, **44**, 441-460.

Carlut, J. & Courtillot, V., 1997. How complex is the time-averaged geomagnetic field over the past 5 myr?, *Geophys. J. Int.*, **134**, 527-544.

Christensen, U. & Olson, P., 2003. Secular variation in numerical geodynamo models with lateral variations of boundary heat flow, *Phys. Earth planet. Inter.*, **138**, 39-54.

Christensen, U. & Wicht, J., 2007. Numerical dynamo simulations, in *Treatise on Geophysics*, Vol. 8, ed. Olson, P., Elsevier Science.

Christensen, U., Aubert, J. & Hulot, G., 2010. Conditions for Earth-like geodynamo models, *Earth planet. Sci. Lett.*, **296**, 487-496.

Cobden, L. & Thomas, C., 2013. The origin of d''reflections a systematic study of seismic array data sets, *Geophys. J. Int.*, **194**, 1091-1118.

Davies, C.J. & Constable, C.G., 2014. Insights from geodynamo simulations into long-term geomagnetic field behaviour, *Earth planet. Sci. Lett.*, **404**, 238-249.

Davies, C.J., Gubbins, D., Willis, A.P. & Jimack, P.K., 2008. Time-averaged paleomagnetic field and secular variation: predictions from dynamo solutions based on lower mantle seismic tomography, *Phys. Earth planet. Inter.*, **169**, 194-203.

Deschamps, F., 2015. Lower mantle electrical conductivity inferred from probabilistic tomography, *Terr. Atmos. Ocean. Sci.*, **26**, 27-40.

Deschamps, F., Cobden, L. & Tackley, P.J., 2012. The primitive nature of large low shear-wave velocity provinces, *Earth planet. Sci. Lett.*, **349-350**, 198-208.

Deuss, A., Ritsema, J. & van Heijst, H., 2013. A new catalogue of normal-mode splitting function measurements up to 10 mHz, *Geophys. J. Int.*, **193**, 920-937.

Finlay, C.C. & Amit, H., 2011. On flow magnitude and field-flow alignment at Earth's core surface, *Geophys. J. Int.*, **186**, 175-192.

Finlay, C.C. & Jackson, A., 2003. Equatorially dominated magnetic field change at the surface of Earth's core, *Science*, **300**, 2084-2086.

Garnero, E.J. & Helmberger, D.V., 1995. A very slow basal layer underlying large-scale low-velocity anomalies in the lower mantle beneath the Pacific: evidence from core phases, *Earth planet. Sci. Lett.*, **91**, 161-176.

Gubbins, D., 2003. Thermal core-mantle interactions: theory and observations, in *Earth's Core: Dynamics, Structure and Rotation*, eds Dehant, V., Creager, K., Karato, S. & Zatman, S., AGU Geodynamics Series American Geophysical Union.

Gubbins, D., Willis, P.W. & Sreenivasan, B., 2007. Correlation of Earth's magnetic field with lower mantle thermal and seismic structure, *Phys. Earth Planet. Inter.*, **162**, 256-260.

Ishii, M. & Tromp, J., 1999. Normal-mode and free-air gravity constraints on lateral variations in velocity and density of Earth's mantle, *Science*, **285**, 1231-1236.

Jackson, A., 2003. Intense equatorial flux spots on the surface of the Earth's core, *Nature*, **424**, 760-763.

Jackson, A., Jonkers, A.R.T. & Walker, M.R., 2000. Four centuries of geomagnetic secular variation from historical records, *Phil. Trans. R. Soc. Lond.*, A **358**, 957-990.

Johnson, C.L. & Constable, C.G., 1995. The time averaged geomagnetic field as recorded by lava flows over the past 5 myr, *Geophys. J. Int.*, **122**, 489-519.

Kelbert, A., Schultz, A. & Egbert, G., 2009. Global electromagnetic induction constraints on transition-zone water content variations, *Nature*, **460**, 1003-1007.

Kelly, P. & Gubbins, D., 1997. The geomagnetic field over the past 5 million years, *Geophys. J. Int.*, **128**, 315-330.

Korte, M. & Holme, R., 2010. On the persistence of geomagnetic flux lobes in global field models, *Phys. Earth planet. Inter.*, **182**, 179-186.

Korte, M., Donadini, F. & Constable, C., 2009. The geomagnetic field for 0-3 ka: 2. a new series of time-varying global models, *J. geophys. Res.*, **10**, Q06008, doi:10.1029/2008GC002297.

Lay, T., Herlund, J. & Buffett, B., 2008. Core-mantle boundary heat flow, *Nature Geosci.*, **1**, 25-32.

Li, M., McNamara, A.K. & Garnero, E.J., 2014. Chemical complexity of hotspots caused by cycling oceanic crust through mantle reservoirs, *Nature Geosci.*, **7**, 366-370.

Livermore, P.W., Hollerbach, R. & Jackson, A., 2013. Electromagnetically driven westward drift and inner-core superrotation in Earth's core, *Proc. Nat. Acad. Sci.*, **110**(40), 15914-15918.

Mao, W.L., Mao, H.-K., Sturhahn, W., Zhao, J., Prakapenka, V.B., Meng, Y., Shu, J., Fei, Y. & Hemley, R.J., 2006. Iron-rich post-perovskite and the origin of ultralow-velocity zones, *Science*, **312**, 564-565.

Masters, G., Laske, G., Bolton, H. & Dziewonski, A., 2000. The relative behavior of shear velocity, bulk sound velocity, and compressional velocity in the mantle: implications for chemical and thermal structure, in *Earths Deep Interior*, Vol. 117, eds Karato, S., Forte, A., Liebermann, R., Masters, G. & Stixrude, L., AGU Monograph.

Mosca, I., Cobden, L., Deuss, A. & Ritsema, J. and Trampert J., 2012. Seismic and mineralogical structures of the lower mantle from probabilistic tomography, *J. geophys. Res.*, **117**, doi:10.1029/2011JB008851.

Nakagawa, T. & Tackley, P.J., 2008. Lateral variations in cmb heat flux and deep mantle seismic velocity caused by a thermal-chemical-phase boundary layer in 3D spherical convection, *Earth planet. Sci. Lett.*, **271**, 348-358.

Olson, P. & Christensen, U., 2002. The time averaged magnetic field in numerical dynamos with nonuniform boundary heat flow, *Geophys. J. Int.*, **151**, 809-823.

Olson, P., Christensen, U.R. & Glatzmaier, G.A., 1999. Numerical modeling of the geodynamo: mechanisms of field generation and equilibration, *J. Geophys. Res.*, **104**, 10383-10404.

Olsen, N., Luehr, H., Sabaka, T.J., Michaelis, I., Rauberg, J. & Tøffner-Clausen, L., 2010. CHAOS-4 – a high-resolution geomagnetic field model derived from low-altitude CHAMP data. In AGU Fall Meeting, Abstract GP21A-0992.
Semenov, A. & Kuvshinov, A., 2012. Global 3-D imaging of mantle con-
ductivity based on inversion of observatory C-responses–II. Data analysis
and results, *Geophys. J. Int.*, **191**, 965–992.

Tackley, P.J., 2011. Living dead slabs in 3-D: the dynamics of compositionally-stratified slabs entering a ‘slab graveyard’ above the core-mantle boundary, *Phys. Earth planet. Inter.*, **188**, 150–162.

Tackley, P.J., 2012. Dynamics and evolution of the deep mantle resulting
from thermal, chemical, phase and melting effects, *Earth Sci. Rev.*, **110**,
1–25.

Tang, X., Ntam, M.C., Dong, J., Rainey, E.S.G. & Kavner, A., 2014. The thermal conductivity of Earth's lower mantle, *Geophys. Res. Lett.*, **41**, 2746–2752.

Terra-Nova, F., Amit, H., Hartmann, G.A. & Trinidade, R.I.F., 2015. The time dependence of reversed archeomagnetic flux patches, *J. geophys. Res.*, **120**(2), 691–704.

Trampert, J., Deschamps, F., Resovsky, J. & Yuen, D., 2004. Probabilistic tomography maps chemical heterogeneities throughout the lower mantle, *Science*, **306**, 853–856.

Wicht, J., 2002. Inner-core conductivity in numerical dynamo simulations,
*Phys. Earth planet. Inter.*, **132**, 281–302.

Williams, Q. & Garnero, E.J., 1996. Seismic evidence for partial melt at the base of the Earth mantle, *Science*, **273**, 1528–1530.

Ziegler, L.B., Constable, C.G., Johnson, C.L. & Tauxe, L., 2011. PADM2M: a penalized maximum likelihood model of the 0–2 Ma palaeomagnetic axial dipole moment, *Geophys. J. Int.*, **184**(3), 1069–1089.