A hybrid mechanism for enhanced core-mantle boundary chemical interaction

,

Direct physical entrainment of core material into rising mantle plumes (Petford et al., 2005;Otsuka & Karato, 2012) may seem to be the most straightforward way to explain the isotopic observations.However, such an exchange may be limited by the high density and low viscosity of the liquid outer core.Furthermore, the absence of a correlated enhancement of siderophile element abundances in lavas bearing low 182 W/ 184 W and high 3 He/ 4 He is inconsistent with direct transport of metal upward into the mantle (Mundl-Petermeier et al., 2017).This latter constraint may only be reconciled if metals and silicates are allowed to undergo chemical interaction in the CMB region, while the metals are left behind as the reacted silicates are subsequently borne upward to the shallow mantle (Mundl-Petermeier et al., 2020).
Although it was originally proposed as a mechanism for producing a high electrical conductivity layer that provides magnetic coupling of core and mantle angular momentum, intrusion of metal into pore spaces inside mantle rocks at CMB dynamic topography lows (Kanda & Stevenson, 2006) may satisfy these constraints.In order for this mechanism to work, liquid metal must "wet" grain boundaries in the rock (Takafuji et al., 2004;Mann et al., 2008) to allow both efficient intrusion and subsequent compaction and expulsion of metals back into the core as material is transported away from the topographic lows where immersion and mixing occurs.Such compaction at the CMB has been shown to be very efficient unless grain sizes are very small, of order 10 µm or less (Buffett et al., 2000).Owing to the small length scales involved, of order the grain size, chemical equilibration inside a metal-silicate "mush" may be expected to occur on time scales much shorter than mantle convection flows (Hernlund & McNamara, 2015).Because CMB dynamic topography of order ∼1 km is expected (Olson et al., 1987), consistent with seismological constraints (Sze & van der Hilst, 2003;Tanaka, 2010), mantle circulation may expose ∼10 21 kg of mantle to silicate-metal interaction every time the CMB is refreshed by mantle convection.While this is small in comparison to the total mass of the Earth's mantle (4×10 24 kg), if the mantle side of the CMB is replaced ∼100 times over Earth's history, then the cumulative amount of exposed mantle material rises to of order ∼1% of the silicate Earth, which may be sufficient to account for occasional observations of core flavors in surface lavas (Hernlund & McNamara, 2015).
In this paper we investigate a scenario in which a metal-silicate "mush" layer is formed by metal intrusion at CMB topographic lows, permitting the mixing, equilibration, and subsequent unmixing of silicates and metals in a Kanda-Stevenson-like mushy layer at the CMB.We additionally consider the potential for weakening and lateral gravitational collapse of the layer, as well as its consequent feedbacks with mantle convection.In particular, we are interested in exploring the degree to which collapse of a mushy layer is able to alter mantle convection circulation in the CMB region and enhance the degree of interaction between core and mantle materials.Using a coupled model of mantle convection and layer collapse, we show that this hybrid "soft CMB" mechanism becomes effective if the viscosity of the metal mush layer is ∼10 5 times smaller than the viscosity of the deep mantle, for which a secondary circulation arises around CMB topographic lows and may begin to exert a strong influence on deep mantle dynamics.

The "Soft CMB" Mechanism
The CMB is depressed into the core in the vicinity of mantle downwelling flows as a consequence of deviatoric stresses derived from buoyancy-driven mantle convection.
The expected dynamic topography at the CMB is of order ∼1 km (Olson et al., 1987).
At CMB pressure-temperature conditions, a liquid iron-alloy is expected to "wet" solid grain boundaries and intrude between the grains to form an inter-connected network (Takafuji et al., 2004;Mann et al., 2008).Combined with the excess fluid pressure head induced in topographic lows, this drives intrusion of metal upward into submerged basal mantle rock (Kanda & Stevenson, 2006).The amount of metal that may be ingested into the mushy region is limited to the disaggregation fraction since solids must maintain a continuous touching network in order to transmit a contrasting pressure gradient relative to liquid metal, and may only penetrate into the mantle by an amount similar to the magnitude of CMB topography (i.e., ∼1 km).
A metal mush mixture formed in CMB topographic lows will be buoyant with respect to the underlying core, and may become rheologically weakened, thus raising the possibility of gravitational collapse.Lateral spreading of metal mush draws more mantle down from above to maintain the dynamic topography dictated by large scale mantle convection (Fig. 1a).By creating a non-linear feedback in the system, such collapse enhances circulation of mantle rock into and through the mushy layer (Fig. 1b).The com- bined effects result in a "softening" of the lower boundary condition for mantle convection in downwelling regions.

Model
Although the dynamics of mantle convection can be highly complex, here we focus on building a basic illustrative model by assuming isoviscous mantle convection of an incompressible Boussinesq fluid in a Cartesian geometry.Normal stresses σ zz exerted by the convective flows v on the CMB raise a dynamic topography h given by, where g is the acceleration of gravity and, where φ is the volume fraction of liquid metal (here it is assumed constant) intruded into the submerged portions of the metal-rock mush (i.e., where h < 0), and ρ m and ρ r are the densities of liquid metal (≈9900 kg m −3 ) and mantle rock (≈5500 kg m −3 ), respectively.The quantity ρ mix is the density of the mushy metal-rock mixture.
The model is started from a quasi-steady convection solution with a downwelling in the middle of the domain and upwellings at the edges.We assume that decompaction and infiltration of metal into submerged rock occurs on time scales much shorter than the residence time of mantle rocks at the CMB.In addition, the reverse process of compaction and expulsion of metal back to the core as the mush moves laterally away from depressions occurs on similarly short time scales.We expect variations in the mushy layer to occur over lateral length scales L that are much larger in comparison to h.In other words, since h/L 1, we apply the "thin-layer approximation" from lubrication theory to describe gravitational collapse of the mushy layer (Reynolds, 1886;Hier-Majumder & Revenaugh, 2010;Hernlund & Bonati, 2019).Gravitational collapse of the mushy layer can be approximated as a diffusion process with where t is time, µ is the (assumed constant) viscosity of the mushy layer, and Since mantle viscous forces are assumed to maintain the equilibrium dynamic topography described by Equation ( 1), keeping h constant for a given buoyancy-driven convection flow, the effect of lateral spreading in the layer is to draw down solid mantle from above.We equate u z+ = (1 − φ)∂h/∂t, where u z+ is the induced draw-down velocity of silicate solids from above at the top of the mushy layer.The factor (1 − φ) accounts for the solid flux into the mushy region that is a mixture of both solids and metals.A secondary collapse-driven flow u thus develops in the mantle that is coupled to gravitational spreading of the mushy layer described by the equation of u z+ at the lower boundary.With the assumption of a linear rheology, the collapse-driven Stokes flow u can be solved separately from buoyancy-driven convection v at each time step, after which they are combined to obtain a total effective velocity v ef f = v + u that is used to advect temperature in the mantle.See the supporting information for more details and the full set of equations.
We neglect the small variations in boundary topography when solving for v, for which we assume free-slip (i.e., tangential stress-free) and impenetrable (i.e., v z (x, y, z = 0) = 0) boundary conditions at the CMB.However, we will need to obtain an estimate of the vertical velocity due to buoyancy-driven flow by itself (independently of collapse-driven flow) at the top of the layer, for which we use, where the same ∂v z /∂z is used to compute h in Eq. ( 1).This will be used to measure the relative contributions of solid flux through the metal-rock mush due to buoyancydriven convection in order to compare it to collapse-driven flux.

Results
We solved for mantle convection flow in 2D Cartesian geometry with a Rayleigh number for Ra = 10 4 −10 6 , where α is the thermal expansivity, ∆T is the super-adiabatic temperature change across the mantle, H is the mantle thickness, η is the reference viscosity of the mantle, and κ is the thermal diffusivity.We vary Ra by changing η while holding other quantities constant.The values used for the parameters are described in Table S1.Two different viscosity contrasts ξ = µ/η are considered here: 10 −5 and 10 −6 .
Larger values (i.e., higher mushy layer viscosities) do not yield any significant collapsedriven flow.These ratios capture the behavior at the point where collapse-driven flux through the mushy layer becomes comparable in magnitude to buoyancy driven-flux due to large-scale convection.
The temperature field, mushy layer thickness, and streamlines for both buoyancy- The pattern of u and v (Fig. 2 c-d) do not change significantly over the parameter ranges considered here.However, their amplitudes are sensitive to the input parameters.This leads to an enhancement of solid flux through the mushy layer that can be quantified as a "gain" G defined as: where F cd and F bd are the mass fluxes due to the collapse-driven and buoyancy-driven flows respectively, and S is the mantle-mushy layer interface.A plot of G as a function of Ra for two viscosity ratios are shown in Figure 3.The gain decreases moderately as Ra increases, while an order of magnitude decrease in ξ leads to an order of magnitude increase in G.

Discussion
The models show that collapse-driven flux reaches parity with buoyancy-driven flux through the mushy layer for ξ ∼10 −5 .As shown in Figure 3, there is a modest decrease in G with increasing Ra, such that this basic conclusion is unlikely to change significantly (at the order of magnitude level) even allowing for broad uncertainties in lowermost mantle properties.G decreases with Ra because mantle viscosity (η) is used as the control variable for convective vigor, thus a reduction in viscosity (increase in Ra) decreases the magnitude of deviatoric stresses acting on the CMB topography more so than flow velocities increase with Ra (v ∝ Ra 2/3 ).This reduction in topography has a strong effect on gravitational collapse due to the non-linear dependence upon h 4 in the diffusion operator of Equation (3).The value of ξ is also an important variable that determines which type of flow dominates the system.We can scale the flux of each flow-type according to the velocities near the CMB as such: u ∼ ∆ρgh 4 /(µL 2 ) and v ∼ δρgHh/η where δρ is the density variations caused by buoyancy forces.Comparing the two velocities gives From Equation ( 1), we obtain a scaling for h according to the densities as follows: h/H ∼ δρ/∆ρ.This is plugged back into Equation ( 7) to eliminate the density ratio and H which finally gives the following scaling for G Equation ( 8) tells us that the gain depends largely on the aspect ratio of the mushy layer and the viscosity contrast between the two domains.A preliminary estimate can be made to determine at which value of ξ the collapse-driven flow becomes dominant (i.e., u/v ≥ 1) by estimating the order of magnitude for each term.Previously h was estimated to be ∼ 10 3 m, while in the numerical models, L ∼ 10 6 m.Combining these values together, we see that flow due to the collapse of the mushy layer becomes dominant when ξ ≤∼ 10 −6 , in good agreement with our results.This implies that once the mushy layer becomes rheologically weak past a certain threshold, the positive feedback on the downwellings begins to dominate flows in the CMB region.
Figure 3 shows a clear negative trend between log 10 G and log 10 Ra that indicates a reduced enhancement of flow into the mushy layer with increasing convective strength of the mantle.In our calculations, the half-width at half the maximum amplitude of the layer was used to approximate the horizontal length scale L. From the numerical models, the mushy layer becomes smaller and narrower with increasing Ra.The following  relations describing the dimensions of the mushy layer with Ra were obtained: h ∼ HRa −0.2325 and L ∼ HRa −0.1750 (see figures S2a and S2c in the supplementary material).Plugging these values into Equation (8) shows that for a constant ξ, G ∼ Ra −0.1150 .This exponent is similar, though slightly larger, than what is obtained in our numerical models (Fig. 3).
The efficacy of the soft CMB mechanism as measured by G dominantly depends on the viscosity ratio ξ between the metal mush and the solid mantle.The viscosity of the mush mixture is expected to decrease as metal fraction increases and drops to values similar to liquid metal above the disaggregation limit (when grains are no longer forming a continuous skeletal touching network).However, the ability for metal to intrude into the pore spaces depends on the existence of a grain-touching network and therefore this limit is never reached under the present assumptions.The key factor is the decrease in mixture viscosity µ corresponding to the maximum infiltration capacity for the mush, at the point where it is no longer able to draw in additional metal.While a ξ of order 10 −5 or smaller is certainly plausible in this scenario, the grain scale dynamics of this process and the effects on mixture viscosity are complex and difficult to constrain, even within several orders of magnitude.
The model presented here is relatively simple and is intended to introduce the basic idea of the soft CMB mechanism.Numerous other complications are expected to influence the efficacy of this mechanism.Variable viscosity, particularly temperature dependence, can have a strong influence on the lower boundary layer for mantle convection and needs to be considered in future studies.Furthermore, chemical reactions between rock and metal following exchange in a mush can change their densities and lead to enhanced convection and/or accumulation of layers on either side of the CMB, depending on whether reactions decrease or increase their densities.Finally, the long time evo-lution with these and other effects also needs to be studied in greater detail, rather than simply considering a snapshot.
In summary, the soft CMB mechanism, whereby chemical interaction in a metalrock mushy layer induced by CMB dynamic topography is enhanced by gravitational collapse, appears to be a viable mechanism for increasing core-mantle chemical exchange.
Further study of this mechanism may generate predictions that can be tested against seismological and other observations.The possibility that hybrid processes like these, which are produced by collaboration of simpler processes occurring across a broad range of scales, additionally serves to illustrate the capacity for nature to find degrees of freedom that often escape our attention.

Introduction
This supporting information contains the description of equations used in the numerical models and gravitational collapse of the mushy layer.

Text S1.
For an isoviscous mantle that is approximated as an incompressible, Boussinesq fluid with an infinite Prandtl number, the governing equations of motion in Cartesian geometry are (e.g., Schubert et al., 2001): where η is the reference viscosity of the mantle,  is the velocity field due to buoyancy variations,  " is the pressure field,  # is the reference density of the mantle,  is the thermal expansion coefficient,  is the temperature,  # is the reference temperature, and  is the acceleration due to gravity.As mentioned in the main text, convective flows exert normal stresses  $$ on the CMB which creates dynamic topography ℎ.In our calculations, we assume that the mantle follows a Newtonian rheology such that This viscous stress deforms the CMB to produce topography variations according to Equation (1) in the main text.
The rate of change of the mushy-layer thickness %& %' in Equation (3) in the main text has contributions from both the silicate and metallic components, and we let the combined vertical velocity of both components at the mantle-mushy layer interface be  ; $ (, ) such that %& %' =  ; $ (, ).Since only the silicate component is exchanged across this interface, the metallic component has to be subtracted from the total velocity such that  ; $ −  ; $ = (1 − ) ; $ where  ; $ is the contribution from the liquid metal.Let  $( be the velocity of the silicate component.Hence,  $( = (1 − ) ; $ .Plugging this result back into Equation (3), the velocity of the silicate component is Gravitational collapse of the mushy layer creates a secondary flow in the overlying mantle  that superimposes on the existing buoyancy-driven flow.This velocity field can be obtained by solving the Stokes' and continuity equations: where  + is the pressure field arising from the flow .These two equations have to be solved with the inhomogeneous boundary condition given by Equation (S4) at each timestep.The two velocity fields  and  are summed together to advect any relevant scalar fields in the domain.In the instance of temperature, the temperature field  will evolve according to the following equation: where  ,--=  +  .
To compute the material flux due to buoyancy-driven flows  ./, we first only approximate  $ (, , ) about the point  = 0 using the Maclaurin series Since the boundaries are impenetrable  $ (, , 0) = 0, the first term on the right side disappears.Ignoring the terms  ! and larger, we can approximate  N $ (, ) =  $ (, , ℎ) as Since  N $ contains contributions from both silicate and metallic components, the latter has to be subtracted from  N $ to obtain the velocity of the silicate component  $( .Hence, This allows us to compute  ./using the following equation where  is the area of the interface between the mushy layer and overlying mantle.Equation (S11) is written with a "rectified" velocity because the flux into the layer is only non-zero when the  $( is less than 0.

Figure 1 .
Figure 1.(a) Schematic illustration of the hybrid mechanism.h represents the thickness of the mushy layer, which is similar to the amplitude of dynamic topography.Black arrows illustrate downwelling mantle flow that induces dynamic topography.(b) Flow chart of the soft CMB mechanism explicitly showing the feedback loop.
driven flow Ψ v and flow due to the gravitational collapse of the mushy layer Ψ u for Ra = 10 4 and ξ = 10 −6 are shown in Figures2a, 2b, 2c, and 2d, respectively.The buoyancydriven flow follows typical convective flow patterns, whereas for the collapse-driven flow, we observe a secondary circulation pattern in the vicinity of the downwelling just above the CMB.The secondary circulation arises from gravitational collapse of the mushy layer and we can see from the streamlines (Fig.2d) that downwelling flows are indeed enhanced, especially close to the CMB.

Figure 2 .
Figure 2. Results for Ra = 10 4 and ξ = 10 −6 at steady state.(a) Temperature field.(b) Mushy layer profile and thickness induced by deviatoric stresses at the CMB.(c) Streamlines of buoyancy-driven flow with black arrows indicating the direction of the flow.(d) Streamlines of collapse-driven flow at the CMB with black arrows indicating the direction of the flow.

Figure 3 .
Figure 3. Plot of G against Ra in log scale.Lines show the least squares linear fit of log 10 G with log 10 Ra.The slope m represents the exponent in the following expression G ∝ Ra m .Black circles and blue circles correspond to ξ = 10 −6 and ξ = 10 −5 respectively.

Figure S1 .
Figure S1.Plot of the initial temperature field during each calculation.In our illustrative models, buoyancy-driven convection is run to steady-state and then collapse-driven flow is initiated.This eliminates transient effects associated with particular choices for initial conditions.

Figure S2 .
Figure S2.In our numerical calculations, the topographic depressions produced by buoyancy forces in the mantle are symmetric about the mid-point (see Figure2bin the main text).We use