Jean-Baptiste Wacheul

Michael Le Bars

Experiments on fragmentation and thermo-chemical exchanges during planetary core formation

Keywords: Planet collision, Two phase flow equilibration, Laboratory experiments

The initial thermo-chemical state of telluric planets was largely controlled by mixing following the collision of differentiated proto-planets. Up to now, most models of planet formation simply assume that the iron core of the impactors immediately broke up to form an "iron rain" within a large-scale magma ocean, leading to the rapid equilibration of the whole metal with the whole mantle. Only recent studies have focused on resolving the fluid mechanics of the problem, with the aim to define more relevant diffusion-advection models of thermal and chemical exchanges within and between the two fluids.

Furthermore, the influence of the viscosity ratio on this dynamical process is generally neglected, whilst it is known to play a role in the breakup of the initial iron diapirs and in the shape of the resulting droplets. Here we report the results of analog laboratory experiments matching the dynamical regime of the geophysical configuration. High speed video recording allows us to describe and characterize the fluid dynamics of the system, and temperature measurements allow us to quantify the diffusive exchanges integrated during the fall of the liquid metal. We find that the early representation of this flow as an iron rain is far from the experimental results. The equilibration coefficient at a given depth depends both on the initial size of the metal diapir

Introduction

The aggregation time scale of moon-size planets in a protoplanetary disk is on the order of 10 to 100 Myr [START_REF] Chambers | Planetary accretion in the inner Solar System[END_REF][START_REF] O'brien | Terrestrial planet formation with strong dynamical friction[END_REF]. It is comparable to the typical time scale of differentiation of terrestrial bodies inferred from meteorite analysis [START_REF] Lee | Hf-W Isotopic Evidence for Rapid Accretion and Differentiation in the Early Solar System[END_REF][START_REF] Kleine | The W isotope evolution of the bulk silicate Earth: constraints on the timing and mechanisms of core formation and accretion[END_REF] and numerical simulations [START_REF] Short Neumann | Differentiation and core formation in accreting planetesimals[END_REF]. This implies that Earth and other terrestrial bodies most likely formed from the collision of already differentiated proto-planets [START_REF] Yoshino | Core formation in planetesimals triggered by permeable flow[END_REF]. These events were by themselves energetic enough to cause a local melting of the mantle [START_REF] Safronov | The heating of the Earth during its formation[END_REF][START_REF] Kaula | Thermal evolution of Earth and Moon growing by planetesimal impacts[END_REF][START_REF] Reese | Fluid dynamics of local martian magma oceans[END_REF][START_REF] Monteux | Thermomechanical adjustment after impacts during planetary growth[END_REF]. In addition, there was enough 26 Al for the decay heat to elevate the mantle's temperature above solidus across depths of several hundreds of kilometers [START_REF] Merk | Numerical Modeling of 26 Al-Induced Radioactive Melting of Asteroids Considering Accretion[END_REF], and thus produce a deep global magma ocean [START_REF] Tonks | Core formation by giant impacts[END_REF].

After the impact, the liquid iron from the impactors' core sank through the molten silicate mantle allowing efficient equilibration, the details of which being a matter of fluid dynamics. Once formed, the core and mantle of terrestrial planets are weakly coupled regarding heat and element exchanges, because the surface of exchange is relatively small and because transfers can only take place through diffusion and very slow advection in the solid mantle.

Consequently, the internal dynamics of these reservoirs strongly depend on their initial thermochemical states [START_REF] Samuel | Heat partitioning in terrestrial planets during core formation by negative diapirism[END_REF][START_REF] King | Heat partitioning in metal-silicate plumes during Earth differentiation[END_REF]. Also, our understanding of the event chronology taking place during planetary accretion relies on radio-nucleid: loose constraints on metal/silicate partitioning at the very end of accretion don't allow to discriminate between a wide range of possible scenarii [START_REF] Kleine | The W isotope evolution of the bulk silicate Earth: constraints on the timing and mechanisms of core formation and accretion[END_REF]Allègre et al., 2008).

Past studies provide us with a partial story on the flow occurring after the impact. Provided that the initial mass of liquid metal can be seen as a unique source referred to as a "diapir", it falls as a turbulent thermal (i.e. an isolated buoyant mass of fluid, in which gravitational potential energy is converted to turbulent motion, causing mixing: see [START_REF] Deguen | Experiments on turbulent metal-silicate mixing in a magma ocean[END_REF], breaks-up within a depth of a few diapir's initial radii to capillary-sized droplets [START_REF] Ichikawa | Direct numerical simulation of an iron rain in the magma ocean[END_REF][START_REF] Samuel | A re-evaluation of metal diapir breakup and equilibration in terrestrial magma oceans[END_REF][START_REF] Wacheul | Laboratory experiments on the breakup of liquid metal diapirs[END_REF], except for the largest diapirs [START_REF] Dahl | Turbulent mixing of metal and silicate during planet accretion and interpretation of the Hf-W chronometer[END_REF] or if planetesimals cores undergo massive mixing during their impact [START_REF] Kendall | Differentiated planetesimal impacts into a terrestrial magma ocean: Fate of the iron core[END_REF]. Turbulent thermals are dominated by collective effects between droplets, until they transition to an "iron rain" regime when diluted enough [START_REF] Bush | Particle clouds in homogeneous and stratified environments[END_REF][START_REF] Deguen | Experiments on turbulent metal-silicate mixing in a magma ocean[END_REF]. The size and speed of droplets in this last regime allow droplets to fully equilibrate with the ambient silicate well before reaching the bottom of the magma ocean [START_REF] Stevenson | Fluid dynamics of core formation[END_REF][START_REF] Rubie | Mechanisms of metal-silicate equilibration in the terrestrial magma ocean[END_REF][START_REF] Ichikawa | Direct numerical simulation of an iron rain in the magma ocean[END_REF][START_REF] Ulvrová | Compositional and thermal equilibration of particles, drops, and diapirs in geophysical flows[END_REF][START_REF] Samuel | A re-evaluation of metal diapir breakup and equilibration in terrestrial magma oceans[END_REF]. However, the typical length scale required for significant equilibration is very important because it sets the depth, i.e. the pressure of equilibration between the metal and silicate.

Recent studies [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF] suggest that the smallest scales of the induced turbulence, even before fragmentation, could drastically reduce this equilibration length scale. It is then crucial to precisely understand how the fluid mechanics of this flow governs diffusion.

We report in the present article the results of laboratory experiments on a fluid system analog to the liquid silicate/ liquid iron system, where we performed measurements of both the fluid dynamics and the thermal diffusive exchanges over the depth of our experimental domain. The set-up and relevant parameters are described in sections 2 and 3 respectively. The fluid dynamics of sedimentation and fragmentation are addressed in section 4. Equilibration results are analysed in section 5 in comparison with various models from the literature. Finally, conclusions and future works are shortly reviewed in section 6.

Experimental set-up

As shown in figure 1, our set-up consists in a 45x45x100cm polycarbonate tank filled with a mixture of water and UCON T M oil of density 1050 kg.m -3 , that is used here to mimic the molten silicate. Its viscosity can be adjusted over 5 orders of magnitude depending on the mass fraction of water, in order to reproduce the range of viscosity ratios of the liquid iron/silicate magma system at low pressure [START_REF] Rubie | Mechanisms of metal-silicate equilibration in the terrestrial magma ocean[END_REF]. To mimic the impacting iron core, we use a balloon filled with galinstan, a gallium alloy that is liquid at room temperature. Its surface energy is 0.718 J.m -2 , and its density 6440 kg.m -3 . Before each experimental run, the mixture water/UCON T M oil is vigorously mixed in order to homogenize the temperature. A balloon is filled with a given mass of galinstan and precisely weighted. A thermocouple is inserted in the balloon before closing it. The balloon is then placed in a hot container in order to heat the galinstan. The viscosity of the mixture water/UCON T M decreases with temperature, therefore the temperature of the galinstan is elevated only by about ten degrees Celsius to avoid lubrication due to higher temperature in the thermal boundary layer. Eventually, the balloon is fixed at the top of the tank where it is immersed in the mixture water/UCON T M oil and popped immediately. The balloon takes only a few thousandths of seconds to retract and the galinstan then starts falling through the viscous fluid. At the bottom of the tank (70 cm below the surface) sits a square funnel made of polystyrene collecting the liquid metal and channeling it towards an adiabatic chamber. When the galinstan is collected, its residual speed mixes completely the temperature and the whole mass of metal reaches an almost constant mean temperature in a few seconds. A thermocouple is placed a few millimeters above the bottom of the adiabatic chamber in order to measure this final temperature. In addition, a third thermocouple is fixed in the middle of the tank one centimeter away from the wall in order to measure the temperature in the bulk of the ambient fluid.

For each experimental run, two events are clearly identifiable in the time series of the thermocouples, as seen in figure 2: the release of the galinstan and the moment it reaches the adiabatic chamber at the bottom of the tank.

From theses signals we compute the temperature drop of the galinstan with respect to ambient fluid during its fall. Each run of the experiment is also recorded by a high speed camera PHOTRON Fastcam at 1000 frames per second with a resolution of 800x1240 pixels.

Parameters and non dimensional numbers

It is still unclear how the initial size of impacting planetesimal converts into the typical size of liquid iron structures at the surface of the magma ocean just after the impact [START_REF] Canup | Simulations of a late lunar-forming impact[END_REF]. However, it is very likely that this typical size depends on the angle of the impact, since stretching and dilution of the core become very important when the impact is far from the vertical direction [START_REF] Kendall | Differentiated planetesimal impacts into a terrestrial magma ocean: Fate of the iron core[END_REF]. In our approach, we consider this typical size as a free parameter, and correspondingly we vary the initial radius of our balloon of galinstan. The other free parameter is the viscosity of the ambient fluid, because the magma viscosity is very dependent on temperature, composition and pressure [START_REF] Karki | Viscosity of MgSiO3 Liquid at Earth's Mantle Conditions: Implications for an Early Magma Ocean[END_REF]. It is adjusted in our set-up by changing the amount of water in the ambient fluid.

Symbol

Definition

Value / range in our experiment (measured experimentally) The two-phase buoyancy-driven flow, theoretically described by the Navier-Stokes and Laplace pressure equations, is then fully characterized by five dimensionless numbers, derived from the parameters listed in table 1:

Re = ρ a .u.R 0 µ a ; W e = ρ a .u 2 .R 0 σ ; Bo = (ρ m -ρ a ).g.R 2 0 σ ; (1)
R µ = µ a µ m ; R ρ = (ρ m -ρ a ) ρ a . (2)
From left to right and up to down, the Reynolds number is the ratio of the flow's inertia over the viscous forces; the Weber number is the ratio of the flow's inertia over the surface tension; the Bond number is the ratio of the buoyancy forces over the surface tension; and the last two numbers are the viscosity ratio and the density anomaly.

We suppose (and we will prove later using experimental results) that the inertial forces are in balance with the buoyancy forces over most part of the flow. This balance sets the order of magnitude of the flow speed, known as the Newtonian terminal velocity and corresponding to the inviscid free fall

speed u ∼ (ρ m -ρ a ) ρ a .g.R 0 1 2 . ( 3
)
This scaling implies that the Bond number is proportional to the Weber number, their ratio being equal to the drag coefficient

Cd = (ρ m -ρ a ).g.R 0 ρ a .u 2 . ( 4
)
Cd will be determined experimentally in the following.

The ranges of dimensionless parameters explored in our set-up are given in is 10 -2 -10 5 , corresponding to regimes which are satisfyingly matched by the experiments. Note that the normalized density anomaly of our fluid system is 5, larger than the corresponding planetary value of about 1.3. In both cases however, the density contrast is globally of order 1, meaning that ambient and metal densities are of the same order of magnitude, with clear non-Boussinesq effects: we thus expect both cases to be in the same dynamical regime.

Note finally that although the proper way to describe the fall and fragmentation of the liquid metal in a viscous fluid involves all dimensionless numbers described above, the experimental results below will also be described in terms of the dimensional parameters R 0 and µ a , which are the two adjustable parameters in our set-up.

Dynamics of the fall

Series of snapshots from raw videos of the experiments are presented in figure 3 for different viscosity ratios. These pictures show several features of importance when considering the diffusion advection process over the whole fall. First, we observe that the flow is always turbulent after falling over a distance equal to 2 radii of the initial diapir, which also corresponds to the length needed to develop strong corrugation of the surface. Our study focuses on the following dynamics. It is also clear from figure 3 that while the ambient viscosity changes over 3 orders of magnitude, the overall falling velocity of the metal remains about the same. This means in our view that our experiments are always in the Newtonian regime, i.e. that their falling dynamics is dominated by inertia and largely independent of viscous dissipation. However, it is also clear that the fragmentation process strongly depends on the ambient viscosity: in our view, this is a consequence of the stress continuity at the metal / ambient fluid interface. Hence, and even if strictly speaking, changing the ambient viscosity leads in our set-up to changes both in the Reynolds number and in the viscosity ratio, observed results are interpreted in the following as a function of the viscosity ratio only, tacitly assuming that the Reynolds number is always large (as for planetary applications).

Dilution of the liquid metal

The dilution of the galinstan mass during its fall seems consistent with the representation of [START_REF] Deguen | Experiments on turbulent metal-silicate mixing in a magma ocean[END_REF], at least for the three lower viscosity ratios explored here (i.e. R µ ≤ 72.5 at least). For these viscosity ratios, the diapir starts deforming on the sides, presenting oscillations of its interface with the ambient fluid because of Kelvin-Helmholtz instabilities (KH). Vortices grow as an evolution of these instabilities and incorporate fluid into the region occupied by the liquid metal. This region grows like a spheroid with a radius apparently increasing linearly with depth, even after the galinstan has broken up in droplets. This is consistent with the classical hypothesis of [START_REF] Morton | Turbulent gravitational convection from maintained and instantaneous sources[END_REF] concerning turbulent entrainment into a buoyant vol-ume of fluid, which states that the rate of entrainment is proportional to the mean fall speed. It follows from this assumption that the radius of the region occupied by the buoyant fluid grows linearly with the depth of this region, regardless of the speed at which it falls (see [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF] for details on the non-Boussinesq turbulent thermal model). We tested this hypothesis using the images extracted from the videos and binarized according to an intensity threshold. We measured the barycenter of the black pixels (i.e.

containing metal) and we assumed that its depth is the same as the actual center of mass of the liquid metal. We also measured the standard deviation of the distance between black pixels and the barycenter of the black pixels, and assumed that it is proportional to the spherical equivalent radius of the region occupied by the galinstan. This radius and its derivative with respect to the depth of the center of mass are shown in figure 4. About 40 experimental videos containing an average of 700 frames were used for each ambient fluid's viscosity to calculate this derivative. For the two lowest viscosity ratios, this derivative is almost constant with depth and reaches a mean value of α = dr dz 0.19, which is close to the values usually found for the coefficient of entrainment for turbulent thermals of particle suspension, ranging from 0.1 to 0.35 [START_REF] Bush | Particle clouds in homogeneous and stratified environments[END_REF]. In summary, for these viscosity ratios the apparent expansion of the liquid metal before and after the breakup is compatible with the scenario of an evolution as turbulent thermals.

Figure 4 also shows a different behavior for the larger viscosity ratio:

there, the initial z-derivative of the radius is smaller than for the other cases, and strongly increases during the fall. The falling diapir clearly does not evolve as a turbulent thermal. This confirms the direct observations on figure 3 : for this large viscosity ratio, the KH instabilities are rapidly saturated and the ambient fluid is not incorporated efficiently in the galinstan. Instead, the galinstan and the mixture water/UCON T M are stretched and form scales that slide on each other, accommodating the shear while minimizing viscous damping (see a close-up in figure 5). The scales are eventually drained in such a manner that very few galinstan separates from the bulk mass. Massive and unique breakup event is not the general breakup mode, in contrast to the 3 other viscosity ratios.

The evolution of the radius for R µ = 72.5 appears to be intermediate between the two extreme trends described above. Nevertheless, the breakup mode is much more similar to the lowest viscosity case. A closer look at the videos of the experiments at R µ = 72.5 gives an explanation for why the apparent evolution of the radius is leaning both ways, while the majority of the liquid metal breaks up in a burst: indeed this case shares with the one at large viscosity ratio long stretched sheets that develop on the side, are peeled where numerous ligaments are stretched and breakup (jellyfish regime) and a regime where fragmentation is mainly due to the growth of Rayleigh-Taylor instabilities (RT piercing regime). For our lower radii, the fragmentation indeed proceeds in a very similar way to the jellyfish regime. For our higher radii with the higher viscosity ratio, the fragmentation proceeds with RTinstabilities and erosion by shear, which can be seen as a transition between the jellyfish and the RT piercing regimes. However, for our higher radii with the 3 lower viscosity ratios, the growth of RT instabilities is generally very limited and they are often advected at the back of the liquid metal mass. Their full development is more occasional than in the experiments of [START_REF] Landeau | Experiments on the fragmentation of a buoyant liquid volume in another liquid[END_REF]. This may be due to the higher density ratio in our experiments, inducing a much higher initial acceleration that could reduce the growth rate of the instabilities and enhance the advection of these perturbations. This effect clearly deserves more study, but is beyond the scope of the present paper.

Fall velocity

From the same series of movies, we computed the velocity of the metal diapir (or of the droplets cloud after breakup) by assuming it equates the speed of the barycenter of the black pixels. Some examples of computed barycenter depth as a fonction of time are shown in figure 6. The velocity is almost constant over the entire fall of each diapir, including after its breakup. We expect that over a longer distance, droplets may finally behave independently and not as a cloud anymore, in which case the global barycenter velocity would be meaningless. Note that this effect is quantified by the so-called Rouse number, equal to the ratio of a single particle terminal velocity in a quiet infinite medium over the velocity of the flow it is caught in. The Rouse number determines if a cloud of droplet behave collectively (Rouse 1) or if the droplets fall independently as they would in an infinite medium (Rouse 1). For our experiments, the Rouse number is between 1/3 and 1/2.

We show on figure 7 our measures of the drag coefficient Cd, i.e. the square of the ratio between the inviscid free fall speed scaling over the measured mean velocity, as defined in equation (4). Cd is almost constant over the explored range of Reynolds number, including experiments with the high- This corresponds to Re = 14, with W e = 31; this is an exemple of a diapir that does not burst although it is above the critical radius for breakup limited by surface tension [START_REF] Hinze | Fundamentals of the hydrodynamic mechanism of splitting in dispersion processes[END_REF].

est ambient viscosity: inertia always dominates the flow. The mean value of Cd is Cd = 3.7 ± 1. This value is larger than those for rigid objects at corresponding Reynolds numbers, which are typically around or below 1 [START_REF] Clift | Bubbles, drops, and particles, 3rd Edition[END_REF]. Such a high value, corresponding to a lower effective falling velocity, is presumably related to the higher area of the diapir's cross flow surface since it deforms easily.

For comparison, the theoretical speed of a turbulent thermal developped in [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF] is also shown in figure 6. The agreement with our experimental results is reasonably good, at least within the uncertainty of our data. We nevertheless claim that our experiments do not exhibit the deceleration predicted by the turbulent thermal model. Two explanations for this discrepancy could be provided. First, as mentioned above, the Rouse number is between 1/3 and 1/2, which means that drops behave collectively but they may not impart their momentum to the ambient with a maximum efficiency as a cloud of droplets with a very low Rouse number would. Sec- ondly, we also argue that some physics may be missing in the thermal model, especially before break up. There, engulfment of the ambient fluid happens mostly at the rear of the diapir (see figure 3). This flux produces an expansion of the liquid metal sheet that works against surface tension. This corresponds to an over-pressure in the liquid that is entrained in the liquid metal, much similar to a balloon. This pressure could force the previously entrained liquid out of the envelope made by the liquid metal sheet. This could result in a linear expansion of the region occupied by the liquid metal with depth, simultaneously with a non-zero flux of momentum associated with the entrainement. This is the main difference with a turbulent thermal, whose equations of movement are based on the assumption that entrainment adds mass with no momentum, the momentum of the buoyant liquid being globally redistributed, leading to a global speed reduction. This explanation for the constant speed cannot be tested for the moment because it requires precise velocimetry of the ambient fluid surrounding the diapir before it breaks up. Nevertheless, it does not question the validity of the assumptions made for the diffusive-advective process of heat in an unbroken turbulent thermal made by [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF], as they rely on the turbulent nature of the flow and on the apparent radius of the liquid metal mass. We will test these scalings in section 5.

Breakup distance

For the highest viscosity case, fragmentation proceeds with "slices" of liquid metal that are peeled away, and by successive division of the initial diapir in several branches, which makes the breakup progressive. On the contrary, the experiments with R µ ≤ 72.5 present one big burst as a consequence of the intense stretching that the galinstan undergoes when incorporating the ambient fluid. The metal becomes a very thin sheet highly convoluted, which is suddenly pierced with many holes. Those holes grow until forming ligament when two holes collide. Each ligament then breaks up as a result of Rayleigh-type capillary instabilities [START_REF] Lasheras | Break-up and atomization of a round water jet by a high-speed annular air jet[END_REF]. This mode of fragmentation allows us to define a quantitative breakup criterion. We define the breakup as the moment when most of the liquid metal sheets have retracted from the holes. In order to compute a breakup distance, we measured several quantities on the binarized images extracted from the video that we assume are affected by this changing topology. These are (see an example figure 8): the quantity of black pixels, the absolute value of its time derivative, and the quantity of white pixels that are enclosed in a black pixel region (a hole). We then compute the following variable as a function of depth:

K(z) = nb a bl * nb b ho * dnb bl dt c , (5)
where nb bl is the number of black pixels, nb ho is the number of white pixels in the regions enclosed by black pixels, and a, b and c are positive exponents.

We finally define the breakup depth as the first depth of the barycenter of black pixels at which this variable K reaches half its maximum. From our experience, fixing a threshold value of half the maximum of the variable K is necessary to avoid spurious maxima that sometimes appear, for instance when the cloud of droplets reaches the bottom of the set-up. The variation of this measured depth is reported in figure 9, non-dimentionnalised by the initial radius of the diapir. The choice of the breakup criterion and of the weighting exponent a, b and c is arbitrary: it was chosen so as to minimize the noise from the experimental data. But even though it is arbitrary, the results presented in figure 9 are only weakly sensitive to this choice. Thus this simple criterion based on light intensity and topology seems robust enough to give meaningful trends.

The smallest diapirs give much dispersed results, because of the weaker probability to record a hole visible on the video. Increasing the mass gives more chance to record holes, which makes our criterium better suited for large diapirs, with a better statistical convergence. Yet a significant variability is still observed, even for the largest diapirs: it corresponds in our view to the natural variability of such non-linear dynamics. Note also again that for our lowest Reynolds data points (i.e. highest viscosity), the breakup distance does not correspond to a clear transition from a single coherent mass to a cloud of droplets: our criteria nevertheless indicates a characteristic distance for a significant change of topology, that will deserve additional study. Here we focus below on the three largest Reynolds numbers that give a clear trend:

334 the breakup distance L bk is proportional to the initial radius, as expected in the literature [START_REF] Deguen | Experiments on turbulent metal-silicate mixing in a magma ocean[END_REF][START_REF] Samuel | A re-evaluation of metal diapir breakup and equilibration in terrestrial magma oceans[END_REF]. We also see here that it weakly depends on the viscosity ratio. From figure 9, one could use a mean value L bk /R 0 = 13 ± 2.

To explain this proportionality, let us assume that the region occupied by the galinstan stays dense (in the mathematical sense) when incorporating ambient fluid in the region of space it occupies on the video (i.e. the largest domain around the barycenter of black pixels where metal is present). Then, the mean stretching undergone by a galinstan fluid particle has to be a growing function of the volume fraction of the ambient fluid in the region occupied by the galinstan. In that case, the diapir deforms and incorporates more and more ambient fluid, until it reaches a critical volume fraction: there, structures are stretched to such an extent that capillary instabilities grow faster than stretching. Let's assume there is a critical volume fraction Φ bk that does not depend on the initial radius. By conservation of the total volume of galinstan, it can be written with the initial radius R 0 and the radius of the region occupied by the galinstan when the breakup occurs R bk

Φ bk = R 0 R bk 3 ⇒ R bk = R 0 . Φ -1 3 bk (6)
As shown in section 4.1, the depth and the radius are proportional, with α = 0.19 the mean derivative of the radius with respect to depth, determined in figure 4. The breakup distance L bk can the be written has

L bk = α -1 (R bk -R 0 ) = α -1 . R 0 .(Φ -1 3 bk -1), (7)
With our experimental results, we find Φ bk = 0.0239 ± 0.01. Viscosity is expected to slow down the dynamics of thin fluid sheets [START_REF] Villermaux | The viscous Savart 42 sheet[END_REF]: this implies a smaller critical volume fraction for a larger viscosity contrast, hence a larger breakup distance, in agreement with the the small global increase at low Reynolds number that might be noticed in figure 9. But we don't have yet a sufficient resolution to quantify this effect. Obviously, the details of the break-up have to be addressed by further studies of the small scale fluid mechanics happening during this extremely fast event. For now, we suggest the use of the breakup distance shown in figure 9 for simulations of planet accretion that consider a post-impact initial condition with a given length scale. The same analysis of the droplets radius distribution as in our previous study [START_REF] Wacheul | Laboratory experiments on the breakup of liquid metal diapirs[END_REF] has been conducted. A better imaging of the experiment in the present set-up enables us to analyze the three lowest viscosity ratios. For R µ = 691, a converged distribution of sizes for stable drops was unaccessible, as only a small fraction of the liquid metal was contained in stable structures. Distributions of sizes for droplets coming from a 20 mm radius diapir are shown in figure 10: it shows a mean radius around 3.5 ± 1 mm, with no clear or systematic trend with the viscosity. This value agrees well with the length scale at which surface tension equilibrates the dynamic pressure due to turbulence l σ [START_REF] Hinze | Fundamentals of the hydrodynamic mechanism of splitting in dispersion processes[END_REF])

Droplet size

l σ = R 0 .W e -3 5 . ( 8
)
For our experiment, this length scale is about 3 mm. Note that this value is commonly used to estimate the mean radius of fragments resulting from turbulence in a diluted two phase flow [START_REF] Kolmogorov | On the breakage of drops in a turbulent flow[END_REF][START_REF] Hinze | Fundamentals of the hydrodynamic mechanism of splitting in dispersion processes[END_REF][START_REF] Chen | Drop size distribution in agitated liquid-liquid systems[END_REF]. As already suggested in [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF], we support the use of this scaling rather than the classical "laminar" capillary length scaling

l lam σ = R 0 .W e -1 2 4.5mm. (9)

Equilibration

Let us define the experimental degree of equilibration:

T * (z) = T initial -T (z) T initial -T ambient (10)
where T initial is the initial temperature of the galinstan, T ambient is the temperature of the mixture water/UCON T M , and T (z) is the mean temperature of the galinstan at depth z.

The experimental degrees of equilibration are shown in figure 11. For reference, this figure also shows the prediction for the degree of equilibration assuming total thermal equilibration between the mass of metal and a finite volume of the ambient fluid. The calculation was made with two volumes:

the volume of ambient fluid that a turbulent thermal with a coefficient of entrainment α = 0.15 incorporates, corresponding to the most conservative estimate of α from our experiment, and the volume of ambient fluid that a turbulent thermal with a coefficient of entrainment α = 0.19 incorporates, which is the average value of α from our experiments. The fact that the experimental degrees of equilibration are all below the predicted ones proves that no full equilibration can be expected, therefore a good understanding of the fluid mechanics is needed to model the exchanges between the two phases, as already claimed by [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF].

Most models of terrestrial planet formation require a length scale of equilibration as they consider potential temperatures and partition coefficients as a function of depth (see e.g. [START_REF] Boujibar | Metalsilicate partitioning of sulphur, new experimental and thermodynamic constraints on planetary accretion[END_REF]. The length scale of equilibration allows the quantification of the precise conditions in temperature and pressure relevant for the equilibration of the liquid metal with the silicate magma. In order to compute such a length scale from our experiments, we use the diffusion equation, stating that the heat flux lost by a diapir is proportional to the temperature difference between its interior and the ambient fluid, and inversely proportional to the size of the diffusive boundary layer separating the two. Taking into account the large thermal diffusivity of the metal and the expected intense flows inside the low viscosity metal drops, we assume that their interior is well-mixed and at a unique mean temperature T (z). Consequently, temperature gradients are set by the diffusion in the boundary layer of size l dif around each drop. Then heat equation writes :

d (V.ρ m .Cp m .T (z)) dt = V.ρ m .Cp m . u . dT (z) dz = -S . λ a . (T (z) -T ambient ) l dif ( 11
)
where V is the total volume of metal, Cp m its heat capacity, S the total metal / ambient fluid surface of exchange, and λ a the thermal conductivity of the ambient fluid. The time evolution of the mean metal temperature, or equivalently its depth evolution since the sinking velocity u is constant, then depends on the heat flux times the surface of exchange, divided by the (constant) volume of metal. At this point, we don't know the exact evolution of the surface of exchange as a function of time, which would help to provide a better model of equilibration. However, since the other parameters can be considered as constant, we will compute a mean experimental length scale of equilibration, by implicitly assuming a surface of exchanges averaged along the fall. With this approximation, the differential equation for the mean temperature of the liquid metal can be simply written:

dT * dz = (1 -T * )/L eq (12)
whose solution writes

T * = 1 -exp -z/Leq . (13)
Here the equilibration length scale is

L eq = l dif . ρ m .Cp m λ a . u . V S (14)
On the basis of our temperature measurements of galinstan at the top and bottom of the set-up, and using equation ( 13), we compute the experimental length scale of equilibration with the following formula:

L eq = - H ln (1 -T * (H)) ( 15
)
where H is the total height of the galinstan's fall.

The equilibration length scales computed from our experiments are presented in figure 12. For comparison, the equilibration length scales according to three different models are computed from the theoretical formula ( 14). Hypotheses for these three models are given in table 3, regarding the three key parameters: the large scale advection speed u, the thickness of the thermal boundary layer l dif , and the typical ratio of volume over surface area V /S.

The first model corresponds to the classical iron rain. It assumes that equilibration is made between the ambient fluid and drops of laminar capillary size (see equation 9) falling at their terminal velocity (i.e. equation 3 using the laminar capillary size). The size of the thermal boundary layer l dif is determined by equilibrating advection at the size of the drops with diffusion through l dif [START_REF] Samuel | A re-evaluation of metal diapir breakup and equilibration in terrestrial magma oceans[END_REF]. Besides, the volume to surface ratio is also given by the laminar capillary size. The corresponding scaling for the length scale of equilibration is then given in terms of global dimensionless numbers by

L eq = c 1
Cp m . ρ m Cp a . ρ a . P e 1/2 . W e -7/8 . R 0 (16)

where c 1 is a non-dimensional coefficient and P e = u.R 0 /D a is the Peclet number of the diapir, D a = λ a /ρ a Cp a being the thermal diffusivity of the ambient fluid. Interestingly, rewriting L eq in terms of primary parameters, we find

L eq = c 1 Cp m . ρ m Cp a . ρ a . 1 D 1/2 a . ρ m -ρ a ρ a g -5/8 . σ ρ a 7/8 , (17)
hence an equilibration length independent of the initial diapir radius, as well as of the viscosity of the ambient fluid: within this model, the equilibration length in our system should be constant over all experimental runs.

On the contrary, the two other models take into account the large scale evolution of the flow as a turbulent thermal, falling at the terminal velocity in the Newton regime of the initial diapir (see section 4.2). Therefore, they assume that turbulence is fully developed and controls diffusive exchanges.

l dif is determined by equilibrating diffusion with advection at the smallest advective scale of the turbulent flow, i.e. the Kolmogorov length scale l kol (see details in [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF]:

l kol = R 0 . Re -3 4 . ( 18
)
Assuming a Kolmogorov type turbulence cascade, the velocity fluctuations u kol at this scale are

u kol = u l kol R 0 1 3 . ( 19
)
Hence,

l dif = P e -1 2 . Re -1 4 . R 0 . (20)
The difference between the two "turbulent thermal" models considered here comes from the assumption for the interfacial area between the ambient fluid and the liquid metal. The model developped by [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF] deals with the liquid metal in its unbroken form, i.e. in its state of crumpled sheet. It assumes that the interface is advected by an isotropic turbulence in the same way as an iso-surface of concentration would. The theory of turbulence then provides the fractal dimension of the interface D f ractal = 8/3 (see [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF], for a detailed analysis). The exact value of the area of the interface between the two fluids is then a unique function of the cutoff scale l below which the fractal aspect of the interface is lost: S = 4π R 2 0 (l/R 0 ) 2-D f ractal . We take for this cutoff length the turbulent capillary length given by (8), i.e. the scale at which surface tension equilibrates the dynamics pressure due to turbulence. In this framework, the interfacial area scales as

S = 4πR 2 0 . W e 2 5 . (21)
Then, the equilibration length is

L eq = c 2 . Cp m .ρ m Cp a .ρ a . P e 1 2 . W e -2 5 . Re -1 4 . R 0 (22)
where c 2 is a constant, or using the primary parameters

L eq = c 2 Cp m . ρ m Cp a . ρ a . 1 ρ 1/4 a D 1/2 a . ρ m -ρ a ρ a g -11/40 . σ ρ a 2/5 . µ 1/4 a . R 23/40 0 . (23)
The last model considered here describes the flow as a turbulent thermal where the buoyant fluid is a suspension of droplets. The interfacial area is just the cumulative surface of the drops and the volume to surface ratio should be well approximated by the mean radius of the drops, corresponding to l σ defined by equation ( 8). Note that we also computed the Sauter diameters from our size distributions, but it does not show any significant difference with the above simple approximation, which is thus used for simplicity. Then, the length scale of equilibration is

L eq = c 3 . Cp m .ρ m Cp a .ρ a . P e 1 2 . W e -3 5 . Re -1 4 . R 0 (24)
where c 3 is a constant, or using the primary parameters

L eq = c 3 Cp m . ρ m Cp a . ρ a . 1 ρ 1/4 a D 1/2 a . ρ m -ρ a ρ a g -19/40 . σ ρ a 3/5 . µ 1/4 a . R 7/40 0 . (25)
The equilibration length for each model is presented in figure 12 together with the experimental measurements. We focus here on the slope of the experimental data compared with the slope predicted by the models. The Table 3: Models tested in this study and their assumptions constants c 1 , c 2 and c 3 have been adjusted for each viscosity ratio to best fit the experimental data. Although the iron rain model is a widely used approximation [START_REF] Karato | Core formation and chemical equilibrium in the Earth -I. Physical considerations[END_REF][START_REF] Rubie | Mechanisms of metal-silicate equilibration in the terrestrial magma ocean[END_REF][START_REF] Ulvrová | Compositional and thermal equilibration of particles, drops, and diapirs in geophysical flows[END_REF][START_REF] Samuel | A re-evaluation of metal diapir breakup and equilibration in terrestrial magma oceans[END_REF], it is incapable of explaining the observed dependence on both the viscosity ratio and the initial radius. In contrast, the turbulent thermal models include a "turbulent" assumption for the thickness of the boundary layer that gives the relevant trend for the influence of the viscosity ratio. Note that this kind of assumption is also used in chemical engineering when considering the mass transfer of bubbly flows in stirred tanks [START_REF] Petitti | Simulation of coalescence, break-up and mass transfer in a gasliquid stirred tank with CQMOM[END_REF]. Besides, the iron rain model is based on the value of the falling velocity at the size of one drop, while the two other models consider the global advection speed based on the falling velocity of the initial diapir, in agreement with our experimental results. This explains the dependence of the equilibration length on the initial radius. One could argue that the diapir always breaks up, forms a droplet suspension, and dilutes to a point where the drops fall independently [START_REF] Bush | Particle clouds in homogeneous and stratified environments[END_REF]: the iron rain should then be the relevant asymptotic behaviour. However, our experiments suggest that unless the length scale of initial structures is already on the same order of magnitude as the capillary length, the transient before forming a diluted rain takes place over a significant part of the fall (see additional discussion in [START_REF] Deguen | Experiments on turbulent metal-silicate mixing in a magma ocean[END_REF], and that a large part of the equilibration is achieved before the flow is diluted enough to be considered as an iron rain. Therefore it is very unlikely that the iron rain scaling is valid for a post impact flow in planets.

It is important to notice that while close, neither the turbulent thermal

Conclusion

To conclude, our measurements of thermal exchanges between a liquid metal and an ambient fluid integrated over the fall show a dependence on the viscosity ratio between the two fluids and on the initial size of the diapir.

We suggest for advanced models of iron/silicate equilibration the use of the turbulent thermal solution for the dilution of the iron diapir in the magma ocean [START_REF] Morton | Turbulent gravitational convection from maintained and instantaneous sources[END_REF][START_REF] Deguen | Experiments on turbulent metal-silicate mixing in a magma ocean[END_REF], and the use of the corresponding turbulent scalings for the length scale of equilibration: respectively (22) before the breakup (modified from [START_REF] Deguen | Turbulent metal-silicate mixing, fragmentation, and equilibration in magma oceans[END_REF], and (24) after the breakup. When building such models, a special attention should be paid when choosing the size of the initial structures as these are still discussed [START_REF] Kendall | Differentiated planetesimal impacts into a terrestrial magma ocean: Fate of the iron core[END_REF] and seem to have a relatively high influence on equilibration. The size of these remnant structures will also determine the initial speed of the liquid metal, which is in our experiment the typical speed of the flow when most of equilibration happens. On the one hand, if one considers the remnant structures as the whole core of the planetesimal, it would probably inherit its initial speed from the impact velocity, which should be close to the escape velocity [START_REF] Dahl | Turbulent mixing of metal and silicate during planet accretion and interpretation of the Hf-W chronometer[END_REF]. On the other hand, if the structures are very well diluted after the impact [START_REF] Kendall | Differentiated planetesimal impacts into a terrestrial magma ocean: Fate of the iron core[END_REF] or, as recently suggested by [START_REF] Kraus | Impact vaporization of planetesimal cores in the late stages of planet formation[END_REF], the iron is vaporized upon impact in a large proportion, the initial speed of the metal would probably be the terminal speed of the remnant structures, which obviously depends on their size. In addition, our experiment also suggests that breakup distance should be set to approximately 13 times the size of the initial diapirs. Since the physics of the exchanges appear to be changing dramatically after the breakup, the size of the initial diapirs is further more important.

Finally, future experiments will be required to better characterize the surface to volume of the droplets distribution and to quantify the exchanges as a function of depth during the fall, in order to complement the integrated measures made in this study.

Acknowledgments

We

Figure 1 :

 1 Figure 1: Schematic and picture of the experimental set-up.

Figure 2 :

 2 Figure 2: Example of temperature measurements: time series of the temperature difference between the ambient fluid recorded at the initial location of the balloon and in the adiabatic container at the bottom of the tank.

Figure 3 :

 3 Figure 3: Series of snapshots taken 0.13 s apart from videos of experiments with a 31 mm galinstan diapir falling through a viscous fluid with different viscosities: from top to bottom, R µ = 0.383, R µ = 8.37, R µ = 72.5 and R µ = 691, corresponding respectively to Re = 34000, Re = 1700, Re = 195 and Re = 20, with W e = 51 in all cases.

  away, and destabilized by Rayleigh-Taylor instability. These structures are common at R µ = 691 and rare at R µ = 72.5, but they still make the area occupied by the liquid metal larger. Despite the corresponding rather large discrepancy shown in measurements of entrainement rate in figure 4, we conclude that the turbulent entrainment hypothesis stands for R µ = 72.5 as for the two lowest viscosity ratios.Using a similar set-up (but without significant variation of the viscosity ratio),[START_REF] Landeau | Experiments on the fragmentation of a buoyant liquid volume in another liquid[END_REF] mapped a regime diagram of the breakup regime. Following their classification, our experiments lie between a regime

Figure 4 :

 4 Figure 4: Top: Radius of the region occupied by the liquid metal measured from video analysis as a function of the depth z of the center of mass of the liquid metal for all experiments. Bottom: z-derivative of the radius of the region occupied by the liquid metal averaged between all the experiments at the same viscosity ratio. The color scale for both plots corresponds to the different viscosity ratios, as given in the bottom plot.

Figure 5 :

 5 Figure 5: Snapshot of the fall of a 24 mm diapir in the most viscous ambient fluid R µ = 691.

Figure 6 :

 6 Figure 6: Depth of the galinstan center of mass as a function of time measured during one run for each viscosity ratio, starting from a 31 mm initial radius diapir. From the lowest viscosity ratio to the highest, this corresponds to Re = 34000, Re = 1700, Re = 195 and Re = 20, with W e = 51 for each case. For comparison, we also show the results of the turbulent thermal model equation of Deguen et al. (2014) integrated in time, using the parameters found by Deguen et al. (2014) (dashed dotted line). The initial speed of the turbulent thermal is set to the speed calculated from the videos.

Figure 7 :

 7 Figure 7: Top: Drag coefficient measured from video analysis as a function of the initial radius of the diapir for each experiment. Bottom: Drag coefficient as a function of the Reynolds number.

Figure 8 :

 8 Figure 8: Number of black pixels (blue) and number of white pixels enclosed in a region of black pixels (green) in each frame as a function of time for several videos of the fall of a 0.031 m radius diapir in a water/UCON T M oil corresponding to R µ ≤ 72.5. The inset shows a colorized snapshot from one video, imaging those two computed quantities using the same color code.

Figure 9 :

 9 Figure 9: Top: Breakup distance (divided by the initial radius) computed from video analysis as a function of the initial radius of the diapir for each experiment. Bottom: Breakup distance as a function of the Reynolds number.

Figure 10 :

 10 Figure 10: Average histogram of the droplets sizes (bin size 1mm) measured on the videos of a 20 mm radius diapir for the three lowest viscosity ratios. From the lowest viscosity ratio to the highest, this corresponds to Re = 15900, Re = 790, and Re = 91, with in each case W e = 18. Each of these histograms results from 6 runs, analyzing for each run the last 20 frames before the diapir reaches the bottom of the experiment.

Figure 11 :

 11 Figure 11: Distance from total equilibration (i.e. 1 minus the experimental degree of equilibration at the bottom of our experiment) from our measurements (symbols) and from the theoretical prediction assuming total equilibration of the liquid metal with the ambient fluid incorporated in the thermal, for both a conservative and the average value of the entrainment coefficient (lines).

  model nor the turbulent cloud of droplets model actually predicts the correct influence of the diapir initial radius (see figure 12 bottom). We think it is due to the approximation of considering a constant / averaged surface of exchange during the whole fall. Each model is indeed relevant with part of the dynamics only. The crumpled sheet turbulent thermal of Deguen et al. (2014) is well fitted for the flow before the breakup, while the cloud of droplets provides a good description of the flow after the breakup. Hence the measured equilibration length is intermediate between our two estimates. To go any further would require additional experiments, focusing on the time evolution of the surface of exchange, including the initial transient.

  acknowledge support from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement No. 681835-FLUDYCO-ERC-2015-CoG). We thank our reviewers H. Samuel and R. Deguen for useful comments and interesting discussions, that helped to significantly improve the first version of this paper. Allègre, C. J., Manhès, G., Göpel, C., Mar. 2008. The major differentiation of the Earth at 4.45 Ga. Earth and Planetary Science Letters 267 (1-2),

Figure 12 :

 12 Figure 12: Top: experimental equilibration lengths (circles) and equilibration lengths computed according to the "iron rain" (black dashed line) and turbulent thermal models (fuchsia dashed line) plotted as a function of the viscosity ratio. Bottom: experimental equilibration lengths (circles) and equilibration lengths computed according to the "iron rain" model (black dashed line), the turbulent thermal model of Deguen et al. (2014) (dashed dotted lines), and the cloud of droplet turbulent thermal model (solid lines), plotted as a function of the initial diapir radius.

Table 1 :

 1 Definition and experimental range of the variables used in this paper.

table 2

 2

	. As an example for comparison, a planetesimal core of radius 100 km
	falls at a typical velocity of 10 3 m.s -1 , giving a Reynolds number of about
	10 11 and a Weber number of about 10 15 . While it is impossible to reproduce
	exactly the planetary regime in the laboratory, our set-up allows exploring
	the relevant range of dynamical parameters, i.e. Re	1 meaning that the
	flow is highly turbulent and dominated by inertia, and W e	1 meaning
	that the initial diapir is prone to fragmentation. Our experiment can also be
	regarded as a close up to the scales where surface tension and viscous forces
	both matter and affect the advection-diffusion process, see also our previous
	study Wacheul et al. (2014). The expected range of viscosity ratio on Earth

Table 2 :

 2 Values of the dimensionless numbers explored with our set-up.