1 Introduction
The seismic F-layer is defined by a decrease in the compression wave (P wave) velocity gradient in an approximately 150–200 km-thick layer in Earth's outer core located just above the inner core boundary (ICB). The anomalous gradient in this region is large enough that it has recently been incorporated into a number of global seismic models, as shown in Fig. 1. Most observations indicate that the F-layer is global, that is, it surrounds the entire inner core (Cormier, 2009; Cormier et al., 2011; Souriau and Poupinet, 1991; Zou et al., 2008), although lateral variations in its properties remain an open possibility (Yu et al., 2005).
Interpretations of the F-layer attribute its anomalous seismic velocity gradient to an approximately monotonic increase in the heavy element concentration (Fe and Ni) with depth, or equivalently, a corresponding decrease in the light element concentration (O, Si, S…) with depth at the base of the liquid outer core (Gubbins et al., 2008). An increase in iron content relative to light elements that accounts for the reduced P wave velocity gradient there has the opposite effect on outer core density (Badro et al., 2007), producing a density increase with depth, i.e., a stable compositional stratification.
The possibility of stable density stratification at the base of the outer core raises important questions about energy transfer and dynamics in the core. According to the standard model of core energetics (Labrosse, 2003), as the core cools, solidification at the ICB partitions heavy elements into the solid and lighter elements into the liquid, providing the primary source of buoyancy for driving convection in the liquid outer core. Under these conditions, the region above the ICB would be expected to have neutral or slightly unstable density stratification due to its elevated light element concentration. Adding a layer with stable stratification complicates this picture in several ways. First, it begs the question of how the F-layer formed. Second, it appears to limit the downward transport of light elements from the outer core to the ICB, thereby inhibiting compositional convection in the outer core.
Several mechanisms have been offered to explain the formation of the F-layer, each with far-reaching implications for the dynamics in the core. One possibility is that the F-layer is a relic of the core formation process. This mechanism is based on the idea that the light element abundances in core-forming metals evolved with time as the Earth accreted in such a way that the core was built with an initial radial stratification consisting of progressively less dense alloys (Hernlund et al., 2013). According to this scenario, the F-layer consists of the remnants of the initial stratification. Another proposal is that the F-layer is maintained by iron solidifying at the top of the layer then re-melting as it precipitates through the layer (Gubbins et al., 2008).
By far the most provocative mechanism, and the one we focus on in this study, assumes that the F-layer is maintained through the interaction of separated melting and solidifying regions distributed over the ICB (Alboussière et al., 2010). Because the ICB is a phase change boundary, substantial lateral variations in temperature must be present there for melting and freezing to occur simultaneously. The core-mantle boundary (CMB) is one possible source of these lateral temperature variations. Experiments (Sumita and Olson, 1999, 2002) and numerical simulations (Aubert et al., 2008) have shown that temperature anomalies generated by strongly heterogeneous CMB heat flux can be transmitted from the CMB to the ICB by outer core convection. Gubbins et al. (2011) demonstrated that, under proper conditions, these CMB-generated temperature anomalies can produce a pattern of simultaneous melting and freezing on the ICB.
The other possibility is convection in the solid inner core. The simplest form of inner core convection consists of melting and solidification in separate hemispheres, leading to a lateral translation of the solid material in the inner core from the freezing hemisphere to the melting one, the so-called inner core translation mode (Alboussière et al., 2010; Monnereau et al., 2010). In addition to the substantial interest in its energetics and dynamics, inner core translation offers a plausible explanation for the observed east–west variations in seismic anisotropy in the inner core (Bergman et al., 2010; Deuss et al., 2010; Niu and Wen, 2001; Sun and Song, 2008; Tanaka and Hamaguchi, 1997), thereby linking its seismic structure to its growth.
The onset of subsolidus convective instabilities in the inner core, including the translational mode, has been examined using several approaches (Buffett, 2009; Cottaar and Buffett, 2012; Deguen and Cardin, 2011; Deguen et al., 2013; Jeanloz and Wenk, 1988; Mizzon and Monnereau, 2013; Weber and Machetel, 1992). The main requirement for convective instability in the inner core is an adverse radial density gradient. Inner core translation is the preferred mode of instability at high viscosity, as it involves little or no solid-state deformation. Cellular (i.e., higher mode) convection favored at viscosities less than about 31018 Pa.s can also produce localized melting and solidification, but less efficiently than the translation mode (Deguen et al., 2013; Mizzon and Monnereau, 2013).
A systematic investigation of the influence of inner core translation on the outer core by Davies et al. (2013) examined thermal convection driven by a spherical harmonic degree and order one pattern of heat flux applied at the ICB. They found that the flow transitions from the usual columnar-style convection when the ICB is homogeneous (Sumita and Olson, 2000) to a pattern dominated by larger-scale mostly prograde (eastward) spiraling jets when the ICB heterogeneity is strong. In cases where the ICB heterogeneity was large enough to simulate melting (corresponding to negative ICB heat flux in their model) the spiral jets found by Davies et al. (2013) resulted in large hemispherical differences in azimuthal velocity everywhere in the outer core, including below the CMB.
Several previous studies have examined dynamo action with spherical harmonic degree one ICB heterogeneity. Results of these studies include dipole eccentricity (Olson and Deguen, 2012) as well as east-west asymmetry in the secular variation of the magnetic field (Aubert, 2013; Aubert et al., 2013). Aubert et al. (2013) have shown that the differences in the geomagnetic secular variation observed between Atlantic and Pacific hemispheres can be produced by a relatively small amount of hemispherically asymmetric inner core buoyancy flux, provided it is properly oriented. In particular, they found that the observed asymmetry of the geomagnetic secular variation is best explained if the buoyancy flux is maximum in the eastern hemisphere, implying westward inner core translation, which is the opposite direction from the original interpretations of the inner core anisotropy (Bergman et al., 2010; Geballe et al., 2013; Monnereau et al., 2010), though in agreement with other interpretations of the observed seismic properties of the inner core (Cormier and Attanayake, 2013; Cormier et al., 2011).
In this study we consider the effects of asymmetric inner core growth including translation in the context of thermochemical convection and dynamo action, using hemispheric variations of codensity at the ICB to represent inner core melting and solidification. First we use numerical simulations of rotating non-magnetic convection to examine the changes in flow and stratification in the outer core that accompany increasing amounts of this type of codensity heterogeneity. We then consider these effects on the geomagnetic field using numerical dynamos that include hemispherical inner core codensity variations. Lastly, we derive an analytical gravity current model for interpreting our numerical simulations and relating them to the observed F-layer structure.
2 Thermochemical convection with hemispherical inner core growth
To model the effects of hemispherically asymmetric inner core growth including melting and solidification on the ICB, we define codensity C in the outer core as the sum of temperature and light element concentration,
(1) |
(2) |
At the CMB we assume no-slip velocity conditions and a uniform positive (destabilizing) codensity flux F0, so that
(3) |
with denoting the dimensionless version of (3). At the ICB, we assume no-slip velocity conditions and we specify a distribution of codensity modeling the combination of uniform plus heterogeneous inner core growth. Specifically, we set
(4) |
where and C′(θ, ϕ) are the mean and laterally varying codensity at the ICB, θ and ϕ denoting colatitude and longitude, respectively. To represent hemispherical ICB melting and solidification, we set
(5) |
We can now define additional dimensionless parameters to characterize the ICB heterogeneity. One new control parameter is a lateral Rayleigh number based on the amplitude of the ICB codensity heterogeneity:
(6) |
A parameter that measures the system response to ICB heterogeneity is the buoyancy flux ratio, . Using the codensity flux at the ICB
(7) |
(8) |
where and denote the negative and positive Fi-values and the integrals are evaluated over the ICB area. The ratios Ra*/Ra and parameterize the extent of melting regions versus solidification regions on the ICB. For large values of these ratios, melting and solidification each cover approximately one-half of the ICB. For intermediate values, the melting region is finite but covers less than one-half of the ICB. For small values of Ra*/Ra or for , (4) is positive everywhere on the ICB and there is no melting. Cases with small values of Ra*/Ra with positive codensity everywhere on the ICB are meant to simulate uneven growth of the inner core that is subcritical to the translation instability. In the other extreme, cases with very large Ra*/Ra are meant to simulate the effects of the translation instability.
3 Non-magnetic thermochemical convection with ICB heterogeneity
We have carried out two sets of non-magnetic thermochemical convection calculations using the MAGIC code (Wicht, 2002) with Y1,1-type ICB conditions, varying the parameters Ra and Ra*. Fig. 2 shows the resulting structure of the convection in the equatorial plane at E = 3 × 10−5 and Pr = 1 from calculations with numerical resolution (nr, nθ, nϕ) = (129, 192, 384). From top to bottom, the images show the effects of increasing amounts of ICB heterogeneity, represented by increasing values of the ratio Ra*/Ra. The first three columns in Fig. 2 are snapshots of codensity, axial (z-component) vorticity, and azimuthal (ϕ-component) of velocity; the last column is the time average of the streamfunction of the flow.
The top row of images in Fig. 2 show a moderate amount of hemispherical asymmetry in all variables at Ra*/Ra = 0.125. Since Fi is positive over the entire ICB in this calculation, there is no melting implied and nothing that would correspond to the F-layer. Nevertheless, convection is substantially more vigorous above the ICB in the hemisphere where C′ > 0. The azimuthal velocities are primarily retrograde (westward), are highest near the solidification pole and near the CMB there is a broad velocity minimum (CMB quiet zone) about 90° east of the solidification pole.
The two middle rows in Fig. 2 show the effects of increasing ICB heterogeneity including regions with Fi < 0, analogous to melting. At Ra*/Ra = 0.25 (2nd row in Fig. 2) the hemispherical differences are magnified compared to the Ra*/Ra = 0.125 case, and a thin spherical cap with positive codensity gradient (negative codensity flux) is clearly evident above the ICB in the hemisphere where C′ < 0. Convection above this stably stratified cap is noticeably attenuated, and the vorticity pattern in the cap indicates the presence of a thin gravity current just above the ICB. In addition, the region with the most intense convection and the CMB quiet zone are both shifted farther to the west. These effects are amplified in the Ra*/Ra = 0.5 case shown in the third row of Fig. 2, where the high density gravity current covers more than half of the ICB and the CMB quiet zone is located almost directly above the solidification pole. In the extreme case with Ra* finite and shown in the bottom row of Fig. 2, the convection is maintained by ICB heterogeneity alone, the gravity current covers most of the ICB, convection is restricted to a localized region, and the CMB quiet zone has been shifted west of the solidification pole.
The winding of the spiral in the azimuthal flow with increasing Ra*/Ra is most evident in the time average streamfunction in Fig. 2, as the longitude of the CMB quiet zone shifts progressively westward. The tendency for spiral tightening with increasing Ra* and decreasing E was first reported by Davies et al. (2013). Interestingly, however, the direction of the dominant azimuthal flow in the equatorial plane of their calculations is prograde (eastward), possibly because they model thermal convection, whereas our convection is primarily compositional. Our azimuthal flows also differ somewhat from the azimuthal flows reported by Aubert (2013) who used heterogeneous codensity flux conditions on the ICB (rather than the heterogeneous codensity ICB conditions we use) and found a prograde and a retrograde spiral jet of nearly equal strength. In light of these differences, it is clear that the large-scale azimuthal flows are sensitive to the choice of boundary conditions and the type of convective forcing assumed, not to mention external factors such as angular momentum exchange with the mantle.
Nevertheless, there are qualitative similarities between our azimuthal flow patterns and the generally westward flow in the equatorial plane of the outer core inferred by Pais and Jault (2008) and Gillet et al. (2009) from the geomagnetic secular variation. Based on the location of the CMB quiet zone, the Ra*/Ra = 0.125 case in Fig. 2 is closest to the model preferred by Aubert et al. (2013) for explaining asymmetry in the geomagnetic secular variation, which would imply no ICB melting.
We have made a second set of non-magnetic calculations in order to determine the minimum amount of melting needed to produce the F-layer. Fig. 3 shows the correlation of the time average versus Ra* from 13 cases at E = 10−4 and various Ra*-values. The inset shows the criterion used to identify stable stratification, namely, a sign reversal in the average codensity gradient above the ICB. Filled and open symbols denote the presence or absence of an F-layer according to this prescription. Although there is considerable scatter, the results in Fig. 3 are consistent with the interpretation that a certain minimum amount of melting, corresponding to , is needed for a stable layer to form.
4 Dynamo action with inner core melting
The previous studies of dynamo action with hemispheric inner core growth imposed relatively modest ICB heterogeneity, generally too small in amplitude to produce negative buoyancy flux and therefore too small to simulate inner core melting. Accordingly, we consider here the influence of large amplitude ICB heterogeneity, as would produce an extensive melting region on the ICB.
Fig. 4 shows the structure of a numerical dynamo with Ra = 4 × 106, Ra* = 2.3 ×107, , E = 10−4, ɛ = –1 along with Pr = Pm = 1, where Pm is the magnetic Prandtl number. The upper row of equatorial plane images are, from left to right, a snapshot of the codensity, and the corresponding snapshot of the streamfunction with z-vorticity contours, and the time average of the streamfunction with z-vorticity contours. As with previous equatorial sections, the melting pole is on the left, the solidification pole is on the right, the solid line marks the melting equator, and the gray shading indicates inner core age increasing from light to dark. The lower row shows the time average of the radial magnetic field on the CMB and a polar view of the radial magnetic field on the CMB.
For this dynamo, Ra*/Ra is large, corresponding to rapid translation, melting is extensive, and the influence of ICB heterogeneity is pervasive. The codensity pattern in the outer core in Fig. 4 reflects the hemispherical ICB codensity pattern, but in this case the contrast between the stable cap on the melting side and the fine scale convection near the solidification pole is particularly evident, and highlights the contrasting dynamics in stable versus unstable regions.
Streamlines in the snapshot and time average images show predominantly westward azimuthal flow throughout most of the outer core, with the CMB quiet zone lying slightly west of the solidification pole, as was found for large Ra* non-magnetic convection. Note the presence of a stagnation point in the snapshot and time average streamline patterns just above the melting pole, marking the center of divergent flow in the dense layer above the ICB, with westward flow to its west and eastward flow to its east. The location of this streamline divergence marks the center of the dense gravity current that in this case extends past the melting equator, covering nearly 80% of the ICB.
The magnetic field structure shows more radial coherence than the azimuthal velocity. Fig. 4 shows the high intensity field is shifted westward by about 90° at the CMB with respect to the solidification pole. More significantly, the dipole axis is offset from the rotation axis. The best-fitting dipole axis, determined using the method of James and Winch (1967), has virtually zero tilt in this case, a consequence of the north-south mirror symmetry in the large-scale magnetic field. However, the best-fitting dipole axis is offset from the center by a distance equal to 0.2ro in the general direction of the inner core melting equator, to the west of the solidification pole and to the east of the melting pole. In terms of offset direction, this is qualitatively consistent with the findings of Olson and Deguen (2012), although the westward shifts of the dipole axes they found were generally smaller than shown here, differences can be attributed to differences in dynamo model parameters and boundary conditions.
5 Gravity currents in the F-layer
Figs. 2 and 4 indicate that formation of the F-layer by a heterogeneous distribution of melting and solidification on the ICB can be conceptualized in terms of a thin gravity current spreading slowly over a spherical surface, subject to a heterogeneous influx of mass representing inner core melting and solidification.
To model these dynamics, we assume the density anomaly in the F-layer decreases linearly with distance above the ICB according to
(9) |
(10) |
The particular form of the friction coefficient depends on the mechanism of resistance to flow in the layer. Resistance to flow from the geomagnetic field is likely to be most important in Earth's electrically conducting core, whereas viscosity resists the flow in our non-magnetic convection calculations. In general, the friction coefficient can be written as a sum of these two effects:
(11) |
Resistance to flow in the layer from the Lorentz force due to a uniform radial magnetic field with intensity B0 yields
(12) |
(13) |
Assuming incompressible flow, conservation of mass for the layer can be written as
(14) |
(15) |
(16) |
and
(17) |
(18) |
The inner core translation hypothesis generally assumes that the melting and solidification poles lie in the equatorial plane, in the simplest case, 180° apart in longitude. Suppose the melting pole is located at (θm = π/2, ϕm), the solidification pole is located at (θs = π/2, ϕs = ϕm + π), and the melting/solidification pattern is symmetric about the the axis defined by these two poles. If λ denotes the great circle angular distance from the melting pole, then according to (17) and (18),
(19) |
(20) |
The seismic observations argue for a global F-layer, formed by a dense stably stratified fluid enveloping the entire inner core, including the solidifying regions. In the context of our gravity current model, the transition from partial coverage to full inner core coverage is defined by the condition h2 = 0 at λ = π. According to (20), the thickness of a global layer at its melting pole λ = 0 must therefore satisfy
(21) |
Assuming σ = 106 S/m, B0 = 10−3 T, and ρ0 = 1.2 × 104 kg/m3 implies in the F-layer. Using g′ = 0.03 m2/s and assuming W = 3 × 10−10 m/s (equivalent to a 300-Myr inner core translation time), the critical value of h0 based on the r.h.s. of (21) is of the order of ten meters. Even this might be an overestimate, since we have ignored the possibility of magnetic flux expulsion from the gravity current. Indeed, the seismically observed F-layer thickness greatly exceeds the minimum thickness needed for full coverage of the ICB for any plausible combination of core properties, which is consistent with its inferred global structure.
When the flow resistance is due to viscosity, as in our non-magnetic convection calculations, the friction coefficient kv defined by (13) depends on the layer thickness h and the structure of the gravity current is slightly different than derived above assuming magnetic field resistance. For a viscous gravity current, the analog equation to (15) for the evolution of the layer thickness is
(22) |
(23) |
For purely hemispherical melting/solidification, the layer thickness varies over the ICB according to
(24) |
(25) |
We can use (25) to calculate the relationship between the angular coverage by the gravity current and its dimensionless thickness at the melting pole. For partial coverage, these parameters are related by
(26) |
Fig. 5 shows the relationship between the angular coverage λmax and melting pole thickness from the non-magnetic convection calculations shown in Fig. 2 compared to the gravity current model predictions using (26). The angular coverage in the convection calculations is defined as the average great circle distance at which the stable stratification vanishes, and the melting pole thickness is estimated using a straight line fit to the codensity profile there.
Considering the various sources of error in these estimates, along with the highly idealized gravity current model, the rather small discrepancy in Fig. 5 is reassuring, although perhaps fortuitous. In particular, the stable layer thickness in the convection calculations is expected to differ from theoretical predictions because diffusion and entrainment effects, both of which tend to disperse the dense layer, are present in the convection calculations but are ignored in the gravity current model. Calculations with greater F-layer coverage have more extensive stably stratified regions where entrainment is presumably weak, whereas cases with less F-layer coverage have weaker, less extensive stratified regions and therefore are more prone to entrainment. Accordingly, cases with limited F-layers are more eroded relative to the non-entraining gravity current model. Finally, it needs to be emphasized that although the relative thicknesses h* of the stable layer in our convection calculations are comparable to the observed thickness of the seismic F-layer relative to the outer core depth ro – ri, this is probably another fortuitous coincidence that results from the choice of the diffusivity in our convection models, which is far larger than appropriate for Earth's core.
6 Conclusions
Our calculations show that asymmetric release of buoyancy at the inner core boundary related to asymmetric inner core growth produces local and global perturbations in the core and the geomagnetic field. Moderately asymmetric inner core growth without melting produces moderate asymmetry in the core flow. In contrast, strongly asymmetric inner core growth with extensive melting, as would accompany rapid inner core translation, produces global-scale asymmetry in the core flow, and a highly eccentric magnetic field structure, perturbations that would be readily observable in the geomagnetic field, its secular variation, and probably also in the paleomagnetic field. The asymmetric structure of our numerical dynamo best conforms to the present-day dipole offset direction and the relatively quiet Pacific geomagnetic secular variation if the central longitude in Fig. 4 approximately corresponds to 180* in the core, which would imply inner core translation from the western to the eastern hemisphere. However, this interpretation rests on very shaky ground, because the observed offset of the geomagnetic dipole is known to be time variable (Olson and Deguen, 2012) and because the longitude shifts in our models depend sensitively on poorly constrained core properties.
Inner core melting, represented in this study by negative buoyancy flux at the ICB, leads to formation of a stably stratified region above portions of the ICB, analogous to the seismic F-layer. Unlike the observed F-layer, the stable layers in our numerical calculations never entirely envelop the inner core, although they cover much of it in our most extreme cases. Our gravity current model explains this discrepancy in terms of the unrealistically high flow resistance in the numerical calculations, compared to the flow resistance expected in the outer core.
Lastly, we speculate that full coverage of the inner core by the F-layer as inferred seismically and predicted theoretically on the basis of its thickness presents conceptual problems for the long-term maintenance of the geodynamo via inner core growth. Full F-layer coverage is expected to inhibit the inward transport of light elements from the outer core to the ICB, meanwhile permitting their outward escape. Consequently, it is possible that over time the F-layer will so restrict the resupply of light elements to the ICB that the compositional buoyancy force in the outer core diminishes. Resupply of light elements could take place if portions of the ICB remain uncovered, or if the stable layer somehow remains permeable to their inward transport. Otherwise, the F-layer could eventually starve the outer core of its main source of power for convection.
Disclosure of interest
The authors declare that they have no conflicts of interest concerning this article.
Acknowledgments
This research was supported by Frontiers in Earth System Dynamics grant EAR-1135382 from the U.S. National Science Foundation.