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

**Fig. 4** | Partial melting in the early mantle. The current adiabatic profile of the mantle corresponds to a surface $T_p$ between 1,600 K and up to 1,750 K (Supplementary Fig. 1). At the Archaean-to-Proterozoic transition, $T_p$ was significantly higher, at 1,800–1,850 K (ref. ³¹). For a $T_p$ even higher than 1,800–1,850 K, a large part of the upper mantle should be partially molten. The transition between the younger, relatively hot and largely molten mantle and the older, colder and less molten mantle could have played a major role in the change of mantle dynamics at the transition. Ol, olivine; Wd, wadsleyite; Rw, ringwoodite; Bg, bridgmanite; Fp, ferropericlase.

## Mantle changes at the Archaean-to-Proterozoic transition

The Archaean-to-Proterozoic transition is associated with major changes in the mantle dynamics with a progressive establishment of modern plate tectonics³²,³⁵. An increasing number of arc signatures in greenstone belts younger than 3.5 Ga suggests that slab subduction became a more-common process after this transition. This global change could be correlated to a transition in the upper-mantle mechanical properties. Indeed, it was shown that the coupling strength between tectonic plates and the underlying mantle is inversely proportional to the viscosity at the lithosphere–asthenosphere interface⁴¹⁻⁴³. A hotter mantle during the Archaean would favour the formation of a partially molten layer and the mechanical weakening at a mid-upper-mantle depth. In addition, shear deformation can produce laminated lithologies and aligned melt accumulations, with further weakening⁴⁴. Therefore, the presence of a partially molten layer could have favoured the early geodynamic regime dominated by a stagnant lid.

With secular mantle cooling, the degree of partial melting in the upper mantle eventually decreased. At one point, the state of the upper mantle became comparable to the present-day situation, in which partial melting still takes place, but at relatively low F values²². The decrease of F below a certain threshold value should be associated with a mechanical hardening, which favours mechanical coupling at the lithosphere–asthenosphere interface. This transition could have induced the change of geodynamic regime at the Archaean-to-Proterozoic transition.

## Methods

Methods, including statements of data availability and any associated accession codes and references, are available at https://doi.org/10.1038/s41561-017-0053-9.

Received: 26 June 2017; Accepted: 13 December 2017;
Published online: 22 January 2018

## References

1. Solomatov, V. S. in *Treatise on Geophysics* Vol. 9 (ed. G. Schubert) 81–104 (Elsevier, Amsterdam, 2015).

2. Romanowicz, B. A global tomographic model of the shear attenuation in the upper mantle. *J. Geophys. Res. Solid Earth* **100**, 12375–12394 (1995).

3. Tauzin, B., Debayle, E. & Wittlinger, G. Seismic evidence for a global low-velocity layer within the Earth's upper mantle. *Nat. Geosci.* **3**, 718–721 (2010).

4. Lay, T., Garnero, E. J. & Williams, Q. Partial melting in a thermo-chemical boundary layer at the base of the mantle. *Phys. Earth Planet. Inter.* **146**, 441–467 (2004).

5. Ferot, A. & Bolfan-Casanova, N. Water storage capacity in olivine and pyroxene to 14 GPa: implications for the water content of the Earth's upper mantle and nature of seismic discontinuities. *Earth Planet. Sci. Lett.* **349**, 218–230 (2012).

6. Aubaud, C., Hauri, E. H. & Hirschmann, M. M. Hydrogen partition coefficients between nominally anhydrous minerals and basaltic melts. *Geophys. Res. Lett.* **31**, L20611 (2004).

7. Herzberg, C., Raterron, P. & Zhang, J. New experimental observations on the anhydrous solidus for peridotite KLB-1. *Geochem. Geophys. Geosyst.* **1**, 1051 (2000).

8. Iwamori, H., McKenzie, D. & Takahashi, E. Melt generation by isentropic mantle upwelling. *Earth Planet. Sci. Lett.* **134**, 253–266 (1995).

9. Tronnes, R. G. & Frost, D. J. Peridotite melting and mineral–melt partitioning of major and minor elements at 22–24.5 GPa. *Earth Planet. Sci. Lett.* **197**, 117–131 (2002).

10. Walter, M. J. Melting of garnet peridotite and the origin of komatiite and depleted lithosphere. *J. Petrol.* **39**, 29–60 (1998).

11. Agee, C. B., Li, J., Shannon, M. C. & Circone, S. Pressure–temperature phase-diagram for the Allende meteorite. *J. Geophys. Res. Solid Earth* **100**, 17725–17740 (1995).

12. Litasov, K. & Ohtani, E. Phase relations and melt compositions in CMAS-pyrolite–H₂O system up to 25 GPa. *Phys. Earth Planet. Inter.* **134**, 105–127 (2002).

13. Katsura, T., Yoneda, A., Yamazaki, D., Yoshino, T. & Ito, E. Adiabatic temperature profile in the mantle. *Phys. Earth Planet. Inter.* **183**, 212–218 (2010).

14. Laporte, D., Toplis, M. J., Seyler, M. & Devidal, J. L. A new experimental technique for extracting liquids from peridotite at very low degrees of melting: application to partial melting of depleted peridotite. *Contrib. Mineral. Petrol.* **146**, 463–484 (2004).

15. Sifre, D. et al. Electrical conductivity during incipient melting in the oceanic low-velocity zone. *Nature* **509**, 81–85 (2014).

16. Sato, H. & Ida, Y. Low frequency electrical impedance of partially molten gabbro: the effect of melt geometry on electrical properties. *Tectonophysics* **107**, 105–134 (1984).

17. Partzsch, G. M., Schilling, F. R. & Arndt, J. The influence of partial melting on the electrical behavior of crustal rocks: laboratory examinations, model calculations and geological interpretations. *Tectonophysics* **317**, 189–203 (2000).

18. Maumus, J., Bagdassarov, N. & Schmeling, H. Electrical conductivity and partial melting of mafic rocks under pressure. *Geochim. Cosmochim. Acta* **69**, 4703–4718 (2005).

19. Andrault, D. et al. Melting of subducted basalt at the core–mantle boundary. *Science* **344**, 892–895 (2014).

20. Davis, F. A., Hirschmann, M. M. & Humayun, M. The composition of the incipient partial melt of garnet peridotite at 3 GPa and the origin of OIB. *Earth Planet. Sci. Lett.* **308**, 380–390 (2011).

21. Yoshino, T., Laumonier, M., McIsaac, E. & Katsura, T. Electrical conductivity of basaltic and carbonatite melt-bearing peridotites at high pressures: implications for melt distribution and melt fraction in the upper mantle. *Earth Planet. Sci. Lett.* **295**, 593–602 (2010).

22. Chantel, J. et al. Experimental evidence supports mantle partial melting in the asthenosphere. *Sci. Adv.* **2**, e1600246 (2016).

23. Hirschmann, M. M. Mantle solidus: experimental constraints and the effects of peridotite composition. *Geochem. Geophys. Geosyst.* **1**, C000070 (2000).

24. Simon, F. & Glatzel, G. Fusion–pressure curve. *Z. Anorg. Allg. Chem.* **178**, 309–316 (1929).

25. Zhang, J. & Herzberg, C. Melting experiments on anhydrous peridotite KLB-1 from 5.0 to 22.5 GPa. *J. Geophys. Res.* **99**, 17729–17742 (1994).

26. Hirose, K. & Kushiro, I. Partial melting of dry peridotite at high-pressures—determination of compositions of melts segregated from peridotite using aggregates of diamonds. *Earth Planet. Sci. Lett.* **114**, 477–489 (1993).
27. Béjina, F., Jaoul, O. & Liebermann, R. C. Diffusion in minerals at high pressure: a review. *Phys. Earth Planet. Inter.* **139**, 3–20 (2003).

28. Hirschmann, M. M., Asimow, P. D., Ghiorso, M. S. & Stolper, E. M. Calculation of peridotite partial melting from thermodynamic models of minerals and melts. III. Controls on isobaric melt production and the effect of water on melt production. *J. Petrol.* **40**, 831–851 (1999).

29. Dalton, C. A., Langmuir, C. H. & Gale, A. Geophysical and geochemical evidence for deep temperature variations beneath mid-ocean ridges. *Science* **344**, 80–83 (2014).

30. Litasov, K. D., Shatskiy, A. & Ohtani, E. Melting and subsolidus phase relations in peridotite and eclogite systems with reduced C–O–H fluid at 3–16 GPa. *Earth Planet. Sci. Lett.* **391**, 87–99 (2014).

31. Herzberg, C., Condie, K. & Korenaga, J. Thermal history of the Earth and its petrological expression. *Earth Planet. Sci. Lett.* **292**, 79–88 (2010).

32. Kamber, B. S. The evolving nature of terrestrial crust from the Hadean, through the Archaean, into the Proterozoic. *Precambrian Res.* **258**, 48–82 (2015).

33. Arndt, J., Barnes, S. J. & Lesher, C. M. *Komatiite* (Cambridge Univ. Press, Cambridge, 2008).

34. Korenaga, J. Initiation and evolution of plate tectonics on Earth: theories and observations. *Annu. Rev. Earth Planet. Sci.* **41**, 117–151 (2013).

35. van Hunen, J. & Moyen, J.-F. in *Annual Review of Earth and Planetary Sciences* Vol. 40 (ed. R. Jeanloz) 195–219 (Annual Reviews, Palo Alto, CA, 2012).

36. Moore, W. B. & Webb, A. A. G. Heat-pipe Earth. *Nature* **501**, 501–505 (2013).

37. Michaut, C. & Jaupart, C. Secular cooling and thermal structure of continental lithosphere. *Earth Planet. Sci. Lett.* **257**, 83–96 (2007).

38. Monteux, J., Andrault, D. & Samuel, H. On the cooling of a deep terrestrial magma ocean. *Earth Planet. Sci. Lett.* **448**, 140–149 (2016).

39. Matsukage, K. N., Jing, Z. C. & Karato, S. Density of hydrous silicate melt at the conditions of Earth’s deep upper mantle. *Nature* **438**, 488–491 (2005).

40. Sakamaki, T., Suzuki, A. & Ohtani, E. Stability of hydrous melt at the base of the Earth’s upper mantle. *Nature* **439**, 192–194 (2006).

41. Sizova, E., Gerya, T., Brown, M. & Perchuk, L. L. Subduction styles in the Precambrian: insight from numerical experiments. *Lithos* **116**, 209–229 (2010).

42. Foley, B. J., Bercovici, D. & Elkins-Tanton, L. T. Initiation of plate tectonics from post-magma ocean thermochemical convection. *J. Geophys. Res. Solid Earth* **119**, 8538–8561 (2014).

43. Doglioni, C., Ismail-Zadeh, A., Panza, G. & Riguzzi, F. Lithosphere–asthenosphere viscosity contrast and decoupling. *Phys. Earth Planet. Inter.* **189**, 1–8 (2011).

44. Jin, Z. M., Green, H. W. & Zhou, Y. Melt topology in partially molten mantle peridotite during ductile deformation. *Nature* **372**, 164–167 (1994).

### Acknowledgements

We thank A. Bouhifd, N. Cayzer, M. Guitreau, T. Kawamoto, T. Komabayashi, D. Laporte, M. Laumonier, H. Martin, S. Parman, B. Reynard and F. Schiavi for help and fruitful discussions. This work was supported by the French National Research Agency (ANR) contract 'OxyDeep'. This research was financed by the French Government Laboratory of Excellence initiative no. ANR-10-LABX-0006, the Région Auvergne and the European Regional Development Fund. This is Laboratory of Excellence ClerVolc contribution number 279.

### Author Contributions

L.H. and G.P. synthesized the glass starting material. G.P. and G.M. performed the EC measurements. G.P., D.A., G.M., N.B.C., D.N., N.G. and A.K. performed the X-ray diffraction measurements, which were subsequently treated by D.A. G.P. and N.B.C. determined the water content in samples using infrared spectroscopy. D.A. wrote the manuscript with help from G.P., G.M. and J.M. All of the authors discussed and commented on various versions of the manuscript.

### Additional information

*Supplementary information* is available for this paper at https://doi.org/10.1038/s41561-017-0053-9.

*Reprints and permissions information* is available at www.nature.com/reprints.

*Correspondence and requests for materials* should be addressed to D.A.

*Publisher’s note:* Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
## Methods

**Preparation of starting material.** Our sample consisted of a mixture of oxides with the composition of the chondritic mantle after core segregation (Supplementary Table 1). We selected this composition because it could represent well both the early Earth's mantle and the deep mantle today. Our composition includes the minor elements Ti, Cr and Na, in addition to the major Si, Mg, Al, Fe and Ca elements. This is an important difference with some previous works based on simplified composition (for example, the CMAS composition¹²). No trace elements were added, as their effect on the solidus is negligible as long as their solubility limit in solid phases is not exceeded¹⁵.

Under upper-mantle P-T conditions, we expect the sample mineralogy dominated, on the one hand, by olivine and its high-pressure polymorphs (wadsleyite (Wd) and ringwoodite (Rw)) and, on the other hand, a variable amount of pyroxenes and majoritic garnet. In our experiments, these lithologies were confirmed by in situ X-ray diffraction (see below). Final-phase fractions may be slightly different than those present in the mantle, first because the true mantle composition remains relatively open and second because of potential chemical heterogeneities. The chemical activities of the major and minor elements should not be largely different if the model composition changes between peridotite, pyrolite or chondritic-type materials. Therefore, these compositions should present similar solidus temperatures.

The mixture of oxides was ground in an agate mortar for approximately 1 h to ensure a homogeneous powder. Then, glass spheres of ~2.5 mm in diameter were prepared using the high-temperature aerodynamic levitation technique associated with CO₂ laser heating at CNRS¹⁶. This containerless method uses a conical nozzle that enables the diffusion of a regulated gas flow below the sample, which eliminates completely any contact between the sample and the container. In this process, the powder is first cold-pressed under a load of 6 tons. Weighted pieces of the recovered agglomerate are then melted at a temperature above 2,273 K, which gives a spherical shape to the sample. After a couple of minutes, samples were quenched to room temperature with a typical cooling rate of 500 K s⁻¹. The gas used for the sample levitation was pure Ar, with less than 0.05 ppm of oxygen; however, the vessel remained open to air during the synthesis at high temperature. This prevents the formation of metallic Fe at very high temperatures when using a reduced atmosphere.

For each experiment, an initial sphere of glass 2–3 mm in diameter was reshaped in the form of a cylinder and directly inserted into the capsule of the high-pressure assembly. Such a starting material provides consistent advantages when performing studies of the melting behaviour. First, the amount of moisture adsorbed at the sample surface is drastically reduced compared with that adsorbed by a powder (Supplementary Fig. 2). Second, the glass presents an ideal chemical homogeneity (as confirmed by electron microprobe analyses) compared with the intrinsic heterogeneities typical of mixtures of oxide powders.

**Determination of the water content in our samples.** Infrared spectroscopy was performed to determine the H₂O content in our samples before and after the melting experiment. Infrared spectra were measured in air on double-polished sections with thicknesses of 600 µm (for the glass starting material) or 900 µm (for the recovered samples). Spectra were collected between 700 and 7,000 cm⁻¹ with 2 cm⁻¹ steps and 300 scans for the sample and 100 scans for the background. In the starting glass, we measured less than 1 ppm wt H₂O, which indicates a starting material nominally dry. For samples recovered from high-pressure–high-temperature conditions, the grain size was significantly smaller (Fig. 2) than the infrared beam. We performed several measurements within the same sample, which appeared to be consistent with each other. From the broad asymmetric band centred at 3,500 cm⁻¹ and the OH peaks related to point defects in olivine⁷ at 3,572, 3,527 and 3,357 cm⁻¹, we calculate a H₂O content of ~90 ppm (Supplementary Fig. 2). These observations confirm that the H₂O contamination remained very limited during the melting experiments. The sample H₂O content on melting is critical because the presence of volatile elements can induce a significant decrease in the melting temperature. Still, due to the presence of 90 ppm water in our samples, we corrected all experimental solidus temperatures by +35 to +55 K, based on the previously reported cryoscopic relation¹⁷ (Supplementary Table 2).

**EC measurements.** We measured the EC of our samples up to their melting point in a 1,500-ton 6-8 split cylinder Kawai-type multi-anvil apparatus at the Laboratory of Magmas and Volcanoes in Clermont-Ferrand. Quasi-hydrostatic pressures of 5 or 20 GPa were obtained using eight tungsten-carbide cubic anvils with edges truncated to 11 or 4 mm, respectively (Supplementary Figs. 3 and 4). We used an octahedral pressure medium with an edge length of 18 or 10 mm made of semisintered MgO with 5% Cr₂O₃. The pressure calibration for these assemblages has an uncertainty of ~0.5 GPa (ref. ⁴⁸). The heater was made of a 50 µm thick Re foil folded into a cylindrical shape and inserted in the ZrO₂ sleeve, which served as a thermal insulator. We applied an alternating current with a frequency of 100 Hz to generate temperatures up to more than 2,500 K. A cylindrical sleeve of MgO was inserted between the heater and the sample to prevent chemical reactions with the Re heater, ensure electrical insulation of the sample and maintain the sample shape on melting.

Temperatures were measured using type C thermocouple (WRe₅%–WRe₂₆%). We added a third WRe₅% wire to the opposite side of the cylindrical sample to allow measurements of the sample impedance. Two Re discs were placed on the top and bottom of the cylindrical sample. The discs served as electrodes for EC measurements. As the determination of temperature gradients is essential in this work, we modelled thermal gradients in the cell assembly based on a finite-element method¹⁹ for thermocouple temperatures between 1,473 and 2,273 K (Supplementary Fig. 3b). As expected, the hottest part is found at the centre of the sample (extreme top-left position in Supplementary Fig. 3b). Each temperature contour (thin black lines) represents a temperature interval of 50 K in Supplementary Fig. 3b. In all the simulations, the temperature difference (ΔT<sub>simul</sub>) between the sample centre and the thermocouple does not exceed 100 K. We therefore apply a temperature correction of ΔT<sub>simul</sub> to all our temperature measurements (Supplementary Table 2).

The sample EC was determined using the impedance spectroscopy method in the frequency range of 10⁶ to 10¹¹ Hz (ref. ⁵⁰) (Supplementary Fig. 4). The insulation resistance (detection limit) of the assembly was determined at similar pressure–temperature conditions prior to the actual experiments. Polycrystalline samples are characterized by a combination of resistor and constant-phase element (R-C/ CPE) circuits and the resistance can be obtained by fitting the impedance spectra to appropriate equivalent circuits. Once the sample resistance has been determined, conductivity can be calculated using the sample diameter and length measured after each experiment, assuming the sample geometry remains unchanged during the experiment. The activation enthalpy (ΔH) of each conduction mechanism can be obtained by fitting data to the Arrhenius equation, where σ is the EC (S m⁻¹), T is the absolute temperature, σ₀ is the pre-exponential factor (S m⁻¹) and k is the Boltzmann constant (J K⁻¹):

$$ \sigma = \sigma_0 e^{(-\Delta H/kT)} $$

For each experiment, we performed a series of heating and cooling cycles, similar to those in previous works⁵¹ (Supplementary Fig. 5a). In a first heating cycle, below 573 K, a low activation enthalpy is likely to result from an extrinsic mechanism associated with the presence of free protons (adsorbed moisture, mostly around the sample). At the desired pressure, a sample was kept at 573 K for more than 12 h. While maintaining 573 K, the electrical resistance of the sample was measured along regular intervals and the next heating cycle started once the sample resistance reached a steady value, which is often 1–2 orders of magnitude higher than the resistance measured at the beginning of the heating cycle. This crucial step ensures the removal of the adsorbed moisture from the sample and the surrounding assembly materials. The adsorbed moisture could otherwise interfere with measurements at higher temperatures.

The crystallization of our samples was observed between 1,370 and 1,570 K, thus at temperatures at least 300 K below the solidus temperature (Supplementary Table 2). The temperature was kept constant above the crystallization temperature for 2–3 h to allow full crystallization of the initial glass samples. The EC of the crystalline samples was measured in several heating and cooling cycles until the heating and cooling paths were reproducible.

The crystallized samples show an activation energy of 1.02 eV (Fig. 1, green line). This is typical of a conduction mechanism with electron hole hopping between Fe²⁺ and Fe³⁺ (small polaron conduction)⁵⁰. At a temperature of ~1,700 K for a nominal pressure of 5 GPa, we observe a first-order change of slope in the conductivity profile (Fig. 1, blue line). The change of slope can be interpreted as the first occurrence of an interconnected partial melt. The experimental uncertainty for the determination of the solidus temperature is about ±25 degree. The errors associated with EC measurements were less than 5%, and were mainly associated with temperature estimations, data-fitting errors and the estimation of sample dimensions. A good reproducibility is observed for EC profiles obtained for our different samples (Supplementary Fig. 5b). Our EC results compare well with ECs of both solid minerals and melts reported in previous studies.

**X-ray diffraction measurements.** In situ X-ray diffraction experiments were performed using the 1,200-ton DIA-type multi-anvil apparatus installed at the PSICHE beamline of the SOLEIL synchrotron in Gif-sur-Yvette. Pressures up to 28 GPa (Supplementary Table 2) were generated using a 7/3 multi-anvil assembly configuration. We used 14 mm WC cubes as secondary anvils (Supplementary Fig. 6). Similar to EC measurements, we used semisintered MgO with a 5% Cr₂O₃ octahedron pressure media and a cylindrical ZrO₂ sleeve in which a 50 µm thick Re-foil heater was inserted. The cylindrical sample of chondritic glass was inserted into an MgO capsule. The Re furnace was oriented horizontally to allow an X-ray beam path parallel to the Re tube. To apply electrical power, two Re electrodes were placed in contact with the furnace. Temperature was measured using a W/Re thermocouple that touched the Re furnace. The thermocouple signal was filtered to remove the common electrical mode at 100 Hz. We used alternatively gaskets made of pyrophilite or a mixture of boron and epoxy. The sample pressure was determined based on the pressure–volume–temperature equation of state of MgO of its capsule, with an accuracy of ~0.5 GPa.

At PSICHE, an under-vacuum wiggler produces a white X-ray beam with energy range from 20 to more than 80 keV. The X-ray beam is focused vertically
between the anvil gap of the multi-anvil apparatus using a horizontal mirror located at 18.5 m from the source and 4 m from the sample. The resulting vertical full-width at half-maximum of the beam at the sample position is approximately 20 µm, which corresponds to almost one order of magnitude flux gain compared with the unfocused beam. The beam is then collimated horizontally to a typical width of 50 µm. Two different modes were used: (1) samples were first aligned in the beam using X-ray radiography, for which the press was scanned vertically in the X-ray beam and the image was then reconstructed numerically (Supplementary Fig. 6c), and (2) diffraction was performed in an energy-dispersive mode using a Caesar-type diffractometer³. A best resolution in interplanar distance *d*<sub>hkl</sub> is achieved when the diffracted X-ray beam is collimated horizontally to ~50 µm by two slit systems located after the sample. A typical diffraction 2θ angle of 8.0° is convenient to register the major diffraction lines of our silicate samples. The resolution in diffraction angle is about 10⁻³ degrees.

On heating, the sample temperature was kept at 700 K for several minutes to remove the adsorbed moisture in the transmitting media. This procedure also mimics the experimental procedure followed for EC measurements. We then increased the temperature until sample crystallization was observed, based on a clear evolution of X-ray diffraction (Supplementary Table 2). At this point, we wait several minutes to finalize the grain growth until no significant change with time was observed on the diffraction pattern. At temperatures below the solidus, the absence of melt prevents the atomic diffusion required for additional grain growth in such a multiphasic system (that is, the Ostwald ripening effect² (Fig. 3)).

For each phase, relative intensities between different diffraction peaks follow quite well those theoretically expected based on the atomic packing in the unit-cell lattice (Supplementary Fig. 6). Such an achievement is typical of very small grains in a sample free of stresses. Still, the use of the energy-dispersive X-ray diffraction technique (compared with the angle-dispersive one) prevents the refinement of phase fractions for the coexisting olivine (or its high-pressure polymorphs), garnet and pyroxenes phases present in our samples. In addition, the limited energy resolution of the Ge detector produces peaks that overlap and the identification of minor phases sometimes becomes difficult. In contrast, this method enables the acquisition of a decent X-ray-diffraction pattern every few seconds, which offers new criteria for the determination of the solidus temperature.

**Analysis of recovered samples.** Recovered samples were cut in two using a diamond-wire saw. A sample half was embedded in epoxy, polished to a mirror-like surface and coated with a thin layer of carbon using a low-vacuum sputter-coater machine. Textural and semiquantitative chemical analyses were performed using a JEOL JSM-5910 LV SEM equipped with a PGT microanalysis system at Laboratory of Magmas and Volcanoes and a Carl Zeiss SIGMA HD VP Field

Emission SEM with an Oxford AZtec ED X-ray analysis operated at the University of Edinburgh. Textural analyses with a field-emission gun SEM were performed at a high accelerating voltage (15 kV) and a small working distance (6–10 mm). These conditions are ideal to detect small features in samples with a low degree of partial melting. For the semiquantitative analysis based on X-ray fluorescence, we instead used a working distance of 20 mm. Energy was also lowered to 5 kV to minimize the penetration of electrons in samples and, consequently, to optimize the spatial resolution.

**Data availability.** The data that support the findings of this study are available from the corresponding author on request.

## References

45. Corgne, A., Liebske, C., Wood, B. J., Rubie, D. C. & Frost, D. J. Silicate perovskite-melt partitioning of trace elements and geochemical signature of a deep perovskitic reservoir. *Geochim. Cosmochim. Acta* **69**, 485–496 (2005).

46. Hennet, L. et al. Aerodynamic levitation and laser heating: applications at synchrotron and neutron sources. *Eur. Phys. J. Spec. Top.* **196**, 151–165 (2011).

47. Tenner, T. J. & Hirschmann, M. M. & Humayun, M. The effect of H₂O on partial melting of garnet peridotite at 3.5 GPa. *Geochem. Geophys. Geosyst.* **13**, Q03016 (2012).

48. Boujibar, A. et al. Metal-silicate partitioning of sulphur, new experimental and thermodynamic constraints on planetary accretion. *Earth Planet. Sci. Lett.* **391**, 42–54 (2014).

49. Hernlund, J., Leinenweber, K., Locke, D. & Tyburczy, J. A. A numerical model for steady-state temperature distributions in solid-medium high-pressure cell assemblies. *Am. Miner.* **91**, 295–305 (2006).

50. Yoshino, T. Laboratory electrical conductivity measurement of mantle minerals. *Surv. Geophys.* **31**, 163–206 (2010).

51. Manthilake, M. et al. Electrical conductivity of wadsleyite as a function of temperature and water content. *Phys. Earth Planet. Inter* **174**, 10–18 (2009).

52. Wang, Y. B. et al. A new technique for angle-dispersive powder diffraction using an energy-dispersive setup and synchrotron radiation. *J. Appl. Crystallogr.* **37**, 947–956 (2004).

53. Yamazaki, D., Kato, T., Ohtani, E. & Toriumi, M. Grain growth rates of MgSiO₃ perovskite and periclase under lower mantle conditions. *Science* **274**, 2052–2054 (1996).
# Evidence for anorthositic crust formed on an inner solar system planetesimal

P. Frossard<fnref n="1" />*, M. Boyet<fnref n="1" />, A. Bouvier<fnref n="2" />,<fnref n="3" />, T. Hammouda<fnref n="1" />, J. Monteux<fnref n="1" />

## Abstract

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

During the first million years of solar system history, planetesimals experienced extensive melting powered by the radioactive decay of $^{26}$Al (Lee et al., 1977). To date, the only known anorthositic crust on a solar system body is that of the Moon, formed by plagioclase flotation on top of the magma ocean (Wood et al., 1970). Here we show evidence from the ungrouped achondrite meteorite Northwest Africa (NWA) 8486 that an anorthositic crust formed on a planetesimal very early in solar system history (<1.7 Ma). NWA 8486 displays the highest anomalies in Eu and Sr found in achondrites so far and, for the first time, this characteristic is also identified in clinopyroxene. Elemental modeling, together with calculated timescales for crystal settling, show that only

the melting of an anorthosite can produce NWA 8486 within the first 5 million years of solar system history. Our results indicate that such a differentiation scenario was achievable over short timescales within the inner solar system, and must have contributed to the making and elemental budget of the terrestrial planets.

OPEN ACCESS

doi: 10.7185/geochemlet.1921

Received 18 March 2019 | Accepted 14 August 2019 | Published 7 October 2019

## Introduction

Over 80 ungrouped achondrites have been found in the past 20 years, enriching our collections with new types of meteoritic samples. Northwest Africa (NWA) 8486 is one of several paired stones, along with the ungrouped achondrites NWA 7325 and 8014, that have been found in the Sahara desert (Ruzicka et al., 2017). They are plagioclase-rich cumulative gabbros that experienced remelting and fast cooling (Yang et al., 2019), with a very peculiar calcic and magnesian mineralogy. NWA 8486 formed under reduced conditions, with an oxygen fugacity ($fO_2$) of 3.2 log units below the iron-wüstite (IW) buffer (Sutton et al., 2017), and therefore likely originates from the inner solar system. Based on Cr, Ti and O isotopic compositions, NWA 7325/8486 have affinities with both acapulcoite-lodranite and ureilite groups (Barrat et al., 2015; Weber et al., 2016; Goodrich et al., 2017). The Pb-Pb isochron age of 4563.4 ± 2.6 Ma, Al-Mg age of 4563.09 ± 0.26 Ma and initial Mg isotopic composition for NWA 7325 indicate that the parent material of this meteorite had to be formed within 1-2 Myr after Solar System formation (Koefoed et al., 2016).

## Eu and Sr Anomalies in NWA 7325/8486

We investigated major and trace element composition of NWA 8486 through *in situ* and in solution analyses. The mineral compositions of NWA 8486 are consistent with those of NWA

7325 (Barrat et al., 2015; Weber et al., 2016; Goodrich et al., 2017), with $An_{88.7\pm3.0}Ab_{11.2\pm3.0}$ plagioclase, $Wo_{45.4\pm0.5}En_{53.4\pm0.5}Fs_{1.2\pm0.1}$ clinopyroxene and $Fo_{97.1\pm0.3}$ olivine. The modal compositions reported for different fragments of NWA 7325 vary (Barrat et al., 2015; Weber et al., 2016; Goodrich et al., 2017), illustrating the heterogeneous distribution of the minerals at the centimetre scale. The fragment of NWA 8486 studied here contains the highest pyroxene modal content at 52 %, the lowest plagioclase content at 44 %, with the remaining 4 % consisting of olivine, sulphides and metal. The range reported for NWA 7325 is 30-44 % for pyroxene, 54-60 % for plagioclase and 2-15 % for olivine (Barrat et al., 2015; Weber et al., 2016; Goodrich et al., 2017).

We report trace element abundances in mineral (plagioclase, pyroxene and olivine) and a whole rock powder of NWA 8486. Minerals were analysed either *in situ* or in solution after mechanical separation (Supplementary Information). The whole rock composition slightly differs from that of NWA 7325 from Barrat et al. (2015) as a consequence of different modal compositions (Fig. 1). Incompatible and moderately volatile elements are depleted in this meteorite, below 0.5 × CI chondrites for the whole rock with the exceptions of Eu and Sr. Positive Eu and Sr anomalies (Eu/Eu* and Sr/Sr*) are present in each mineral phase with different magnitudes, from 1.8 to 6.5 in pyroxene and 450 to 1039 in plagioclase (Supplementary Information). Their amplitude in the whole rock is much higher than that measured in lunar anorthosites (Fig. 1).

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

<footnote n="1">1. Université Clermont Auvergne, CNRS, IRD, OPGC, Laboratoire Magmas et Volcans, F-63000 Clermont-Ferrand, France</footnote>

<footnote n="2">2. Department of Earth Sciences, Centre for Planetary Science and Exploration, University of Western Ontario, Ontario, Canada</footnote>

<footnote n="3">3. Bayerisches Geoinstitut, Universität Bayreuth, Germany</footnote>

<footnote>* Corresponding author (email: paul.frossard@uca.fr)</footnote>
<figure><img src="image_340.png" /><figcaption>**Figure 1** Trace element composition of NWA 8486 normalised to CI chondrites (Anders and Grevesse, 1989). All fractions exhibit Eu and Sr positive anomalies. NWA 7325 has a lower content for most incompatible elements compared to NWA 8486 owing to its modal enrichment in plagioclase. Lunar anorthosites are reported for comparison (data from Haskin et al., 1973 and Norman et al., 2003).</figcaption></figure>

Europium changes valence in a continuum from a 3+ state to a 2+ state around the IW buffer (where Eu³⁺/ΣEu = 0.5) and in these conditions plagioclase preferentially incorporates Eu²⁺ in its lattice, resulting in an enrichment in Eu over the other trivalent rare earth elements (REE³⁺). The similar plagioclase/melt partition coefficients of Sr²⁺ and Eu²⁺ arise from very close ionic radii. NWA 7325/8486 formed within an fO₂ of IW-3.2 ± 0.2, which is more reduced than for lunar rocks that formed in a range of IW-2 to IW (Wadhwa, 2008; Sutton et al., 2017). Despite their higher fO₂, lunar anorthosite plagioclases show similar Eu anomalies to reduced meteorites (Fig. S-4). Thus, even considering variations due to parental melt composition, there should be little difference between NWA 8486 and lunar anorthosite plagioclase compositions if they were formed under similar conditions. Although positive Eu and Sr anomalies are common in plagioclase, the present data are the first report of positive Eu and Sr anomalies in clinopyroxene in an achondritic meteorite. Europium partitioning in clinopyroxene changes depending on the composition of the system with Eu²⁺ being either similarly or less compatible than Eu³⁺ (Karner et al., 2010). Therefore, the peculiar trace element inventory of NWA 8486 cannot be solely attributed to redox conditions.

## Anorthosites as Source Rocks of NWA 7325/8486

NWA 8486 stands out as a unique meteorite when compared to the composition of other achondrites and lunar anorthosites. Pyroxene, plagioclase and whole rock data from all known types of achondrites define separate fields in an EuN/SmN versus SrN/NdN diagram (with N indicating concentrations normalised to corresponding CI chondrite abundances) (Fig. 2). The composition for NWA 8486 is shifted towards both higher EuN/SmN and SrN/NdN ratios compared to other planetary compositions. The rare felsic (Si-rich, evolved) achondrites GRA 06128/9 and Almahata Sitta (ALM-A clast) do not appear any different from other achondrites in terms of Eu and Sr

enrichments. Early occurrences of such felsic rocks have been related to partial melting of chondritic material instead of the multi-stage processes involved to form felsic crust on Earth (Day et al., 2009; Bischoff et al., 2014). Therefore, the particular composition of NWA 7325/8486 must be related to the melting of a specific source which was already enriched in Eu (relative to other REEs) and Sr. Barrat et al. (2015) suggested that impact melting of a gabbroic crust might produce such anomalies. In view of a homogeneous distribution of ²⁶Al in the solar system, Barrat et al. (2015) and Koefoed et al. (2016) calculated that the differentiation event for NWA 7325's parent body was ~1-2 Myr after formation of the solar system, while its crystallisation (from Pb-Pb chronometry) occurred within about 3 Myr after that. Considering the diopside-anorthite binary system at 1 atm (Osborn, 1942), the source of NWA 7325/8486 parent melt needs to be already enriched in plagioclase to reach the high plagioclase modal compositions for the meteorite, at least 48 % plagioclase to match NWA 8486. The modal composition of NWA 8486 lies on the eutectic point of the diagram. Thus, basaltic sources with low plagioclase content can produce eutectic compositions, but higher anorthite modes cannot be reached. We performed melting models of varied compositions and compared the compositions of these liquids with those in equilibrium with NWA 8486 minerals (Supplementary Information). The enrichment in Eu and the general depletion of incompatible elements are not reproduced with cumulate eucrite (Moore County) compositions (Fig. 3). The only possibility found to produce a parental magma enriched in both Eu and Sr is melting of a plagioclase-, Mg- and Ca-rich rock. Low degree melting (below 5 %) of an anorthosite (98 % anorthite) produces liquids corresponding to NWA 7325/8486 modes, but does not yield both high enough Eu anomalies and REE depletion. However, a higher pyroxene content can produce a larger amount of melt of eutectic composition. NWA 8486 features can be reproduced by 20-50 % melting of a pyroxene-rich anorthosite (Apollo noritic anorthosite 62236 with 20 % pyroxene) (Fig. 3). Anorthosites are therefore the most likely source for the parental magma of the NWA 7325/8486 suite.

<figure><img src="image_341.png" /></figure>
<figure><img src="image_342.png" /><figcaption>**Figure 2** Eu/Sm and Sr/Nd ratios of achondrites and NWA 8486 normalised to CI chondrites (Anders and Grevesse, 1989). Overlap between the different fields is due to whole rocks mainly composed of pyroxene or plagioclase (e.g., aubrites, ureilites, lunar anorthosites). For each field group, NWA 8486 exhibits higher Eu<sub>N</sub>/Sm<sub>N</sub> and Sr<sub>N</sub>/Nd<sub>N</sub>. Only lunar anorthosites are similar to NWA 8486 whole rock, but they contain much more plagioclase than NWA 8486. See Supplementary Information for data sources.</figcaption></figure>

<figure><img src="image_343.png" /><figcaption>**Figure 3** REE composition of liquids in equilibrium with NWA 8486 minerals (in grey dashed lines) compared to liquids modelled with compositions of the cumulate eucrite Moore County and the pyroxene-rich lunar anorthosite Apollo 62236, normalised to CI chondrites (Anders and Grevesse, 1989). A non-modal melting is considered, in agreement with petrological constraints, of eutectic proportions of 42 % plagioclase and 58 % pyroxene (Osborn, 1942). The range of composition of the liquids from 1 % to 50 % degree of melting is represented for each source composition. NWA 8486 whole rock composition is shown (black line) for comparison. See Supplementary Information for details on the model.</figcaption></figure>

<figure><img src="image_344.png" /></figure>
<figure><img src="image_345.png" /><figcaption>**Figure 4** Time for plagioclase crystal ascent during the late stage of a magma ocean. This time is calculated for a range of viscosities (0.1 to 10 Pa.s), crystal diameters (a; 100 µm to 1 cm) and density contrasts ($\Delta\rho$; 50 to 500 kg/m$^3$) between the crystal and the liquid phase. The time for a crystal to reach the surface of the magma ocean is represented with a colour scale ranging between 10 kyr (blue) and 2 Myr (red). Both panels show that the time of ascent of the crystal is a few tens of thousands of years in magma ocean conditions, except for small crystals around 100 µm and high viscosities of >10 Pa.s for which time of ascent is longer.</figcaption></figure>

## Timescales of Anorthosite Formation on Early Solar System Planetesimals

Anorthosites are associated with large magma systems that experience fractional crystallisation, in which plagioclase floats at the surface. Lunar anorthosites are derived from a magma ocean (MO) that formed subsequently to the Moon-forming impact event, but their petrogenesis model is still debated. Specific conditions are required for anorthositic crust formation (Albarède and Blichert-Toft, 2007). The body needs to be rather dry, as plagioclase formation is delayed in presence of water and would appear late in the crystallisation sequence (Elkins-Tanton, 2008). The size of the planetary body is also important as plagioclase is stable up to 1 GPa, which corresponds to a depth of 75 km for a Mars-like body or 200 km for a Moon-like body (Albarède and Blichert-Toft, 2007). Studies on Mercury's MO have shown that in very reduced conditions anorthite is denser than the melt and then may not float (Brown and Elkins-Tanton, 2009; Vander Kaaden and McCubbin, 2015). Crystal settling that enables anorthositic crust formation may be hindered by highly turbulent convection on small parent bodies with low gravity. Nonetheless, some authors have suggested that crystal settling can occur in such environments (Taylor et al., 1993; Elkins-Tanton, 2012).

Geochemical modelling indicates that the NWA 7325/8486 source was most likely an anorthosite. Although Mercury may not be an analogue to NWA 7325/8486 parent body, it may represent an end member in terms of density and composition of the MO. Calculations have been carried out to constrain the timescales of formation of an anorthositic crust in a MO. The plagioclase formation occurs generally when the MO is crystallised by 70-80 vol. % (Supplementary Information). Assuming that the remaining molten material (≈20 vol. %) forms a buoyant layer over the solidified mantle and that 15 vol. % of the inner body is made of a metallic core, we can estimate the depth of the molten layer relative to the radius of the body (Supplementary Information). Considering that the rising velocity for plagioclase crystals in a MO

scales with a Stokes' Law settling velocity (Martin and Nokes, 1989), the time for a plagioclase crystal to reach the surface in a shallow MO can be assessed. The calculation yields that for the MO viscosity range (10⁻¹-10¹ Pa.s, Dygert et al., 2017) and density contrasts between crystal and magma (50-500 kg/m³, Brown and Elkins-Tanton, 2009; Dygert et al., 2017), a crystal of at least 100 µm in diameter will reach the surface of the MO in a few tens of thousands of years (Fig. 4). The time of ascent of plagioclase is much smaller than the age of the first differentiation on NWA 7325/8486 estimated at 1.7 Ma using the ²⁶Al-²⁶Mg systematics (Koefoed et al., 2016). Therefore, we believe that NWA 7325/8486 is the first evidence of the formation and re-processing of an anorthositic crust on a planetesimal very early in the solar system.

If anorthosite could form quickly on planetesimals and magma oceans were common in the first few million years of the solar system, where are the anorthosites in the meteorite record? The reason for the absence of anorthositic crusts on achondrite parent bodies probably lies in their thermal and chemical evolution. The conditions necessary for magma oceans to form with anorthositic crusts may have been limited by the time and place of accretion (Greenwood et al., 2012). Inner solar system planetesimals that accreted early were enriched in ²⁶Al and were thus more likely to have experienced magma ocean conditions (Grimm and McSween Jr, 1993). Iron meteorites are evidence of these processes. Their parent bodies were rather small (20 to 200 km; Chabot and Haack, 2006) and possibly formed within the terrestrial accretion zone below 1 AU before being scattered into the main asteroid belt. Silicate layers from these bodies have a low probability of survival in the chaotic early inner solar system (Bottke et al., 2006, O'Neill and Palme, 2008). It is likely, then, that anorthositic crusts were formed in the inner solar system but rarely preserved. Material from these planetesimal silicate layers could have been subsequently added to the accreting terrestrial planets, providing a non-chondritic source for refractory lithophile elements.

<figure><img src="image_346.png" /></figure>
Acknowledgements

We thank Dr. L. Garvie (Center for Meteorite Studies, Arizona State University) for the loan of NWA 8486, J.-L. Devidal for assistance with electron microprobe and LA-ICP-MS analyses, and J.-M. Hénot for assistance with SEM analyses. We also thank T. Withers for comments and J.-A. Barrat for enriching discussions on geochemical modelling. Helen Williams is acknowledged for editorial handling and Addi Bischoff and Stephen Elardo for constructive and thorough reviews that greatly improved the manuscript. This project has received funding from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (Grant Agreement No 682778 – ISOREE to MB), and NSERC Discovery Grant, Canada Research Chair, and the Canada Foundation for Innovation JELF programs (to AB). This is Laboratory of Excellence ClerVolc contribution number 369.

Editor: Helen Williams

Additional Information

Supplementary Information accompanies this letter at http://
www.geochemicalperspectivesletters.org/article1921.

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

This work is distributed under the
Creative Commons Attribution
Non-Commercial No-Derivatives
4.0 License, which permits unre-

stricted distribution provided the original author and source
are credited. The material may not be adapted (remixed, trans-
formed or built upon) or used for commercial purposes without
written permission from the author. Additional information is
available at http://www.geochemicalperspectivesletters.org/
copyright-and-permissions.

Cite this letter as: Frossard, P., Boyet, M., Bouvier, A., Hammouda, T., Monteux, J. (2019) Evidence for anorthositic crust formed on an inner solar system planetesimal. *Geochem. Persp. Let.* 11, 28–32.

References

ALBARÈDE, F., BLICHERT-TOFT, J. (2007) The split fate of the early Earth, Mars, Venus, and Moon. *Comptes Rendus Geosciences* 339, 917-927.

ANDERS, E., GREVESSE, N. (1989) Abundances of the elements: Meteoritic and solar. *Geochimica et Cosmochimica Acta* 53, 197-214.

BARRAT, J.-A., GREENWOOD, R.C., VERCHOVSKY, A.B., GILLET, P., BOLLINGER, C., LANGLADE, J.A., LIORZOU, C., FRANCHI, I.A. (2015) Crustal differentiation in the early solar system: Clues from the unique achondrite Northwest Africa 7325 (NWA 7325). *Geochimica et Cosmochimica Acta* 168, 280-292.

BISCHOFF, A., HORSTMANN, M., BARRAT, J.-A., CHAUSSIDON, M., PACK, A., HERWARTZ, D., WARD, D., VOLLMER, C., DECKER, S. (2014) Trachyandesitic volcanism in the early solar system. *Proceedings of the National Academy of Sciences* 111, 12689-12692.

BOTTKE, W.F., NESVORNÝ, D., GRIMM, R.E., MORBIDELLI, A., O'BRIEN, D.P. (2006) Iron meteorites as remnants of planetesimals formed in the terrestrial planet region. *Nature* 439, 821-824.

BROWN, S.M., ELKINS-TANTON, L.T. (2009) Compositions of Mercury's earliest crust from magma ocean models. *Earth and Planetary Science Letters* 286, 446-455.

CHABOT, N., HAACK, H. (2006) Evolution of Asteroidal Cores. In: Lauretta, D.S., McSween Jr., H.Y. (Eds.) *Meteoritics and the Early Solar System II*. The University of Arizona Press, Tucson, 747-771.

DAY, J.M., ASH, R.D., LIU, Y., BELLUCCI, J.J., RUMBLE III, D., McDONOUGH, W.F., WALKER, R.J., TAYLOR, L.A. (2009) Early formation of evolved asteroidal crust. *Nature* 457, 179-182.

DYGERT, N., LIN, J., MARSHALL, E.W., KONO, Y., GARDNER, J.E. (2017) A low viscosity lunar magma ocean forms a stratified anorthitic flotation crust with mafic poor and rich units. *Geophysical Research Letters* 44, 11,282-11,291.

ELKINS-TANTON, L.T. (2008) Linked magma ocean solidification and atmospheric growth for Earth and Mars. *Earth and Planetary Science Letters* 271, 181-191.

ELKINS-TANTON, L.T. (2012) Magma oceans in the inner solar system. *Annual Review of Earth and Planetary Science* 40, 113-139.

GOODRICH, C.A., KITA, N.T., YIN, Q.-Z., SANBORN, M.E., WILLIAMS, C.D., NAKASHIMA, D., LANE, M.D., BOYLE, S. (2017) Petrogenesis and Provenance of Ungrouped Achondrite Northwest Africa 7325 from Petrology, Trace Elements, Oxygen, Chromium and Titanium Isotopes, and Mid-IR Spectroscopy. *Geochimica et Cosmochimica Acta* 203, 381-403.

GREENWOOD, R., FRANCHI, I., GIBSON, J., BENEDIX, G. (2012) Oxygen isotope variation in primitive achondrites: The influence of primordial, asteroidal and terrestrial processes. *Geochimica et Cosmochimica Acta* 94, 146-163.

GRIMM, R.E., MCSween Jr, H.Y. (1993) Heliocentric zoning of the asteroid belt by aluminum-26 heating. *Science* 653-655.

HASKIN, L. A., HELMKE P. A., BLANCHARD D. P., JACOBS J. W., TELANDER K. (1973) Major and trace elements abundances in samples from lunar highlands. *Proceedings of the Lunar Science Conference* 4, 1275-1296.

KARNER, J.M., PAPIKE, J.J., SUTTON, S.R., BURGER, P.V., SHEARER, C.K., LE, L., NEWVILLE, M., CHOI, Y. (2010) Partitioning of Eu between augite and a highly spiked martian basalt composition as a function of oxygen fugacity (IW-1 to QFM): Determination of Eu²⁺/Eu³⁺ ratios by XANES. *American Mineralogist* 95, 410-413.

KOEFOED, P., AMELIN, Y., YIN, Q.-Z., WIMPENNY, J., SANBORN, M.E., IIZUKA, T., IRVING, A.J. (2016) U-Pb and Al-Mg systematics of the ungrouped achondrite Northwest Africa 7325. *Geochimica et Cosmochimica Acta* 183, 31-45.

LEE, T., PAPANASTASSIOU, D., WASSERBURG, G. (1977) Aluminum-26 in the early solar system-Fossil or fuel. *Astrophysical Journal* 211, L107-L110.

MARTIN, D., NOKES, R. (1989) A fluid-dynamical study of crystal settling in convecting magmas. *Journal of Petrology* 30, 1471-1500.

NORMAN, M.D., BORG, L.E., NYQUIST, L.E., BOGARD, D.D (2003) Chronology, geochemistry, and petrology of a ferroan noritic anorthosite clast from Descartes breccia 67215: Clues to the age, origin, structure, and impact history of the lunar crust. *Meteoritics and Planetary Science* 38, 645-661.

O'NEILL, H.S.C., PALME, H. (2008) Collisional erosion and the non-chondritic composition of the terrestrial planets. *Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences* 366, 4205-4238.

OSBORN, E.F. (1942) The system CaSiO₃-diopside-anorthite. *American Journal of Science* 240, 751-788.

RUZICKA, A., GROSSMAN, J., BOUVIER, A., AGEE, C.B. (2017) The Meteorit-
ical Bulletin, No. 103. *Meteoritics and Planetary Science* 52, 1014.

SUTTON, S.R., GOODRICH, C.A., WIRICK, S. (2017) Titanium, Vanadium and Chromium Valences in Silicates of Ungrouped Achondrite (NWA) 7325 and Ureilite Y-791538 Record Highly-Reduced Origins. *Geochimica et Cosmochimica Acta* 204, 313-330.

TAYLOR, G.J., KEIL, K., MCCOY, T., HAACK, H., SCOTT, E.R. (1993) Asteroid differentiation: Pyroclastic volcanism to magma oceans. *Meteoritics* 28, 34-52.

VANDER KAADEN, K.E., MCCUBBIN, F.M. (2015) Exotic crust formation on Mercury: Consequences of a shallow, FeO-poor mantle. *Journal of Geophysical Research: Planets* 120, 195-209.

WADHWA, M. (2008) Redox conditions on small bodies, the Moon and Mars. *Reviews in Mineralogy and Geochemistry* 68, 493-510.

WEBER, I., MORLOK, A., BISCHOFF, A., HIESINGER, H., WARD, D., JOY, K., CROWTHER, S., JASTRZEBSKI, N., GILMOUR, J., CLAY, P. (2016) Cosmochemical and spectroscopic properties of Northwest Africa 7325—A consortium study. *Meteoritics and Planetary Science* 51, 3-30.

WOOD, J.A., DICKEY JR, J., MARVIN, U.B., POWELL, B. (1970) Lunar anortho-
sites and a geophysical model of the moon. *Geochimica et Cosmochimica
Acta Supplements* 1, 965-988.

YANG, J., ZHANG, C., MIYAHARA, M., TANG, X., GU, L., LIN, Y. (2019) Evidence for early impact on hot differentiated plantesimal from Al-rich micro-inclusions in ungrouped achondrite Northwest Africa 7325. *Geochimica et Cosmochimica Acta* 258, 310-335.

<figure><img src="image_348.png" /></figure>
# Deep mantle origin of large igneous provinces and komatiites

Rémy Pierru<fnref n="1" />*, Denis Andrault<fnref n="1" />, Geeth Manthilake<fnref n="1" />, Julien Monteux<fnref n="1" />, Jean Luc Devidal<fnref n="1" />, Nicolas Guignot², Andrew King², Laura Henry²

Large igneous provinces (LIPs) resulted from intraplate magmatic events mobilizing volumes of magma up to several million cubic kilometers. LIPs and lavas with deep mantle sources have compositions ranging from komatiites found in Archean greenstone belts to basalts and picrites in Phanerozoic flood basalt and recent oceanic islands. In this study, we identify the mantle conditions appropriate to each type of lava based on an experimental study of the melting of pyrolite. The depth of the mantle source decreases from 600 to 700 km for the oldest komatiites to 100 to 300 km for picrites and basalts, and the extent of mantle melting ranges from 10 to 50%. We develop a geodynamical model that explains the origin of the hot mantle plumes capable of generating these melting P-T conditions. Within a superadiabatic temperature gradient persisting in the deep mantle, the ascent of hot mantle plumes creates excess temperatures up to 250 to 300 K by adiabatic decompression.

Copyright © 2022
The Authors, some
rights reserved;
exclusive licensee
American Association
for the Advancement
of Science. No claim to
original U.S. Government
Works. Distributed
under a Creative
Commons Attribution
NonCommercial
License 4.0 (CC BY-NC).

## INTRODUCTION

### The mantle potential temperature

The mantle temperature (T) is generally defined by the potential surface temperature ($T_p$), which corresponds to the extrapolation to the surface of the adiabatic T gradient. The mantle $T_p$ can be refined from the analysis of primary magmas formed along accretionary ridges (1) and is calculated after a calibration of, e.g., the melt MgO content against several parameters including T and pressure (P) (2, 3). Currently, $T_p$ of 1350° ± 50°C is suggested by primary basalts with MgO and FeO contents at 10 to 13 weight % (wt%) and 6.5 to 8.0 wt%, respectively (4). Still, there are “hot” primary magmas (basaltic-picritic) requiring $T_p$ as high as ~1600°C (5). These are “hot spots,” such as Hawaii, and they could be related to instabilities in the thermal boundary layer in the lowermost mantle (6).

The Archean mantle $T_p$ has important implications for understanding the geodynamics of the early Earth and the onset of plate tectonics. Recent geodynamic models have recognized the importance of the reevaluation of Archean mantle $T_p$ for crustal formation mechanisms in stagnant lava environments (7). Average $T_p$ estimates for large igneous provinces (LIPs) and Hawaii is $\Delta T = +200$ K or higher than the average mantle $T_p$. In Portugal, where komatiites occur, $\Delta T$ can exceed +300 K (8).

Regarding the accuracy of the $T_p$ determination, note that the agreement between the compositions of experimental melts and primary magmas on the field are in very good agreement only if the source is not very deep (P lower than 6 to 7 GPa). At higher P, the accuracy of $T_p$ is more questionable, especially for Archean sources of ultramafic rocks. This is, in part, because experimental conditions are less well controlled at high P and the composition of the melt at equilibrium is also more difficult to measure. For example, estimates of the $T_p$ required to produce a primary magma with MgO >20 wt% diverge by ~100 K with the PRIMELT calculation (3). In a more general manner, any uncertainty on the melting diagram is directly transferred into uncertainties on the knowledge of the melting condition in the mantle.

In addition to this, anomalously high water contents in the mantle could produce comparable effects as high T, as water lowers the solidus. Where $T_p$ estimates from dry material show rather consistent results, melts produced from a wet mantle do not seem to provide usable computational solutions (9). More generally, the likely presence of chemical heterogeneities increases the uncertainty on the $T_p$ determination locally (10).

### A broad range of primitive lavas

Mafic-ultramafic tholeiitic lavas such as basalts, picrites, and komatiites (11) represent the deepest and hottest magmas that erupted on Earth. Three main families of komatiites are distinguished depending on their chemical composition; nicknames have been given to them according to locality of the main outcrops (12): Al-depleted komatiites (ADK; Barberton type), Al-undepleted komatiites (AUK; Munro type), and Al-enriched komatiites/picrites (AEK; Gorgona type). These rocks show a large diversity of chemical composition resulting from different geological contexts including (i) mantle source, composition, depth, extent of melting (Φ), and possible presence of fluids and (ii) fractional crystallization, contamination, or assimilation possibly occurring during the ascent of melts to the lithosphere (13). There is a clear evolution of the composition of LIPs and komatiites through time over the past ~4 billion years (Ga; Fig. 1 and fig. S1). Unlike basalts and picrites that can be found at all ages, ADK seems to be only formed during the Archean from -4 to -2.8 Ga, while AUK erupted until about -1.9 Ga in Paleo-Proterozoic. Association of basalts and picrites with komatiites can be found, for example, in the form of dykes into Archean greenstone belts (14). Last, AEK formed from the Archean to the Phanerozoic. The most recent komatiites have been dated at ~90 Ma (15).

Picrites are rare intrusive igneous lavas enriched in magnesium (>12% MgO) (16) compared to the basalts erupting today on Earth. These lavas seem to be products of mantle melting at P above ~3 GPa and extent of melting (Φ) close to 20 to 30% (17). However, several other possible origins were proposed such as high rate of melting of lherzolitic mantle, partial melting of mantle with high Mg content, or olivine enrichment through processes of magmatic differentiation.

<footnote n="1">¹Université Clermont Auvergne, CNRS, IRD, OPGC, LMV, F-63000 Clermont-Ferrand, France. ²Synchrotron SOLEIL, Gif-sur-Yvette, France.
*Corresponding author. Email: r pierru30@gmail.com</footnote>
<figure><img src="image_349.png" /><figcaption>Fig. 1. Age of the different lavas found in LIPs: References are the same as those used for fig. S1.</figcaption></figure>

Komatiites contain up to 34 wt% MgO, present very high liqui-
dus *T* at atmospheric *P* (18) and have a close chemical relationship
with mantle peridotites. For these reasons, komatiites must have formed
at relatively high extents of melting in a mantle warmer than today
(+200 to 300 K) (1). Alternatively, komatiites could be hydrated
magmas formed at *T* not much higher than the surrounding mantle
(5); the source of volatiles would be non-degassed primordial man-
tle or subducted oceanic crust.

1) ADKs present the very high MgO contents (≥30 wt%) and CaO/Al₂O₃ ratios (>~1.5), low Al₂O₃/TiO₂ ratios (<10 to 15), and high Gd/Yb ratios, compared to chondrites. They are also enriched in trace elements and light rare earth elements (REEs) but depleted in heavy REE (HREE). However, their composition could not be reproduced experimentally from melting of typical mantle materials. Still, we note that CaO/Al₂O₃ ratios greater than 1.5 in melts were successfully reproduced by melting either KLB-1 peridotite at extent lower than 1% (19), peridotite with addition of up to 10% water (20), or ad hoc dry or hydrated CaO-MgO-Al₂O₃-SiO₂ (CMAS) compositions (21). The composition of ADK suggests melting with a residual Al-bearing phase (majoritic garnet and bridgmanite), which requires *P* over 15 to 17 GPa (22–25). Zirconium and hafnium, which are high in ADK, act as incompatible elements for melting processes at moderate Φ only (26), suggesting ADK magmas produced at lower Φ (~30%) than the other Archean komatiites.

2) AUKs present low CaO/Al₂O₃ ratios (~1 or less), high Al₂O₃/TiO₂ ratios (~15 to 20), low MgO content (<30 wt%) compared to ADK, and an REE pattern close to chondrites. A robust characteristic is moderate to strong depletion of HREE and other immobile incompatible trace elements. The degree of depletion approaches, but rarely meets or exceeds, that of modern mid-ocean ridge N-type basalts (MORB). AUKs are less rich in zirconium than ADKs (12). According to batch melting models based on the available partition coefficients [e.g., (18)], AUK could result from high extent of melting (~50%). The composition of AUKs from the late Archean are difficult to explain by a single melting event, as they present conflicting nonchondritic Al₂O₃/TiO₂ and CaO/Al₂O₃ ratios implying the presence of residual garnet and quasi-chondritic Gd/Yb ratio implying no residual garnet (23). Alternatively, a deep mantle source undergoing high extent of melting could produce a melt with high

MgO content, however depleted in Al₂O₃, before the slow uplift of
such melt would induce progressive enrichment in Al₂O₃ and de-
pletion in MgO [e.g., (2, 12)]. All arguments suggest AUK mantle
sources considerably shallower than for ADKs. This difference, to-
gether with a lower MgO content (<30 wt%), indicates lower forma-
tion *T* for AUKs than ADKs.

3) AEKs have similar CaO/Al₂O₃ ratios and REE patterns to AUKs (12). However, they present high Al₂O₃/TiO₂ ratios ranging from the chondritic value of ~20 to as high as 50. AEKs are richer in Al₂O₃ than other komatiites and contain less MgO (>20 wt%). Accordingly, they formed by mantle melting at relatively low *P* (3 to 4 GPa), under conditions where garnet is absent from the residue. If the source were to be a mantle fertile peridotite, then the average *Φ* should be 30 to 40% (26). In terms of chemistry, AEK makes a transition from komatiites to picrites, as both compositions are similar.

Possible origin of the primitive lavas

The correlation between komatiite compositions, timing of the erup-
tions, and the depth of the primary sources indicates that over time,
there is less magma originating from the deeper mantle regions
(19, 23, 27). This evolution could be explained by the secular mantle
cooling; however, this hypothesis remains partially speculative. The
debate is still vivid about the magmatic processes capable of gener-
ating LIPs with large volume of mafic to ultramafic magmas in short
pulses (typically less than a million years). In addition, komatiitic
ultramafic lavas require large extent of melting at high *T* (28, 29).
The greenstone belts, for example, may represent a small fragment
only of a much larger LIP (11). LIPs and komatiites could be associ-
ated with the rise of hot mantle plumes from the deep mantle (30),
potentially induced by the formation of a subducted-slab graveyard
(31). This plume could generate a positive anomaly in *T* and/or
volatile element (H₂O, CO₂) concentration and favor major mantle
melting. Other models propose shallower mechanisms such as mantle
warming beneath supercontinents (32), plate-related processes (29),
or lithospheric delamination (33).

The isotopic signature of LIPs and komatiites inform on the na-
ture of their source. Similar to modern hotspots [e.g., Galapagos,
Iceland, Hawaii, and Samoa (34)], the LIP signature overlaps with
many mantle components such as depleted MORB mantle (DMM),
ocean island basalts (OIB), Focal zone (Fozo), etc. (35). FOZO is interpreted as being located in the lower mantle (36). Together with high ³He/⁴He anomalies measured in some LIPs, it provides evidence for a source originating from deep mantle regions containing primordial ³He (28). In contrast, the DMM component represents the upper mantle depleted by ~3 Ga of MORB extraction. In addition, reservoirs enriched mantle 1 (EM1), enriched mantle 2 (EM2), and hight-µ (HIMU) are likely located in the sublithospheric mantle. The overlaps found in LIPs suggest the involvement of deep sublithospheric or asthenospheric mantle components, as well as crustal or lithospheric contamination. Concerning komatiites, isotopic studies of Hf, Sr, and Nd show the overlap with the OIB component [e.g., (37)]. Together with high levels of He isotopes (38), it indicates the presence in the komatiite source of material from the deep mantle. With the high *T* required to form komatiitic magmas, this information strongly supports that komatiites formed by melting of mantle plumes coming from the very deep mantle.

Up to now, the variation of magma chemistry with *P*, *T*, and *Φ* is only well constrained experimentally for the first 100- to 200-km mantle depths (23, 39–44). Higher *P* have also been investigated, however, for compositions different from those expected for the sources of LIPs and komatiites [e.g., CaO-MgO-Al₂O₃-SiO₂ system (21)] and limited ranges of *Φ* [e.g., around 1% (19)]. Previous experiments proposed global tendencies for the change of melt composition with *P* and *T*; however, the

ranges of *P*-*T*-*Φ* covered by experimental studies remain insufficient to accurately model mantle melting at depths greater than 200 km.

## RESULTS

### Melting diagram of pyrolite

#### Solidus temperatures

To bring new constraints on the origin of LIPs and komatiites, we undertook an experimental determination of the melting diagram and of the melt composition for the first 700-km depth (i.e., ~26 GPa) for a pyrolitic mantle composition (table S1). Our experimental strategy was to synthesize samples at controlled *Φ* values using in situ techniques: electrical conductivity and x-ray diffraction (figs. S2 to S5 and table S2). Pyrolite solidus is found similar, within ~100 K, to the chondritic-type mantle measured with similar methods (45). This is substantially lower than previous reports on peridotite solidus (fig. S6) (19, 23, 27, 42, 46). As developed in (45), reasons for this discrepancy with previous works cannot be a difference in composition alone. It could result from several factors: (i) Use of simplified sample compositions, for example, CMAS (21). It is well accepted that additional elements induce lower solidus *T* for a typical mantle mineralogy. (ii) Use of coarse grain powders, instead of a glass sample in our study. Slow atomic diffusion may impede the formation of a solidus melt composition at grain interfaces. (iii) Lack of in situ

<figure><img src="image_350.png" /><figcaption>**Fig. 2.** Progressive increase of the extent of melting associated with disappearance of residual minerals. Left frames show experiments realised at 4 to 5 GPa, central frames 15 to 17 GPa, and right frames 22 to 25 GPa. (A and B) MA105, (C) T2, (D) MA119, (E and F) MA116, (G) MA149, (H) MA143, and (I) MA144. At *P* around 4 to 5 GPa, several mineral phases [clinopyroxene (Cpx) and garnet (Gt)] disappear at small extent of melting. With increasing *T*, orthopyroxene (Opx) is enriched in compatible elements (i.e., MgO) and presents grain size exceeding that of olivine. At higher *T*, melting of Opx is associated with a major growth of olivine grain, which is the liquidus phase. At *P* around 15 to 17 GPa, wadsleyite (Wad) and garnet continue to coexist up to very high melting extent. Garnet appears to be the liquidus phase of pyrolite at *P* > 15 GPa. For *P* > 22 GPa, ringwoodite (Rw) and garnet disappear at a relatively low concentration (<25%). After that, bridgmanite (Bg) starts to melt and ferropericlase (Mw) becomes the liquidus phase.</figcaption></figure>
melting criterion and characterization method insufficiently sensitive to detect thin melt films in the recovered samples. Small pockets of melt can be destroyed during the polishing steps. In addition, melt films can be difficult to identify using the scanning electron microscope (SEM) (47). (iv) Last, important thermal gradient within the sample charge can lead to *T* and *Φ* errors.

Apart from (45) that is in good agreement with our melting curves, none of the studies cited above used in situ techniques to detect the onset of melting (such as electrical conductivity and x-rays), and none probes the evolution of $Φ$ at higher $T$ (such as a falling sphere). In previous works, melting was generally detected on the basis of “trials and errors” and the solidus and liquidus $T$ were sometimes identified using the $T$ gradient in a single sample. However, major thermal gradients induce melt migration; the presence of chemical gradients within the capsule disables the determination of the melting behavior of the starting material. For example, a solid residue identified in an unmelted zone of the sample can result from the compaction of refractory minerals, which melting $T$ can be far above the solidus and liquidus of the bulk composition.

### Extent of melting

To probe the melting properties at higher $Φ$, we developed an original falling sphere technique to detect in situ the occurrence of the rheological transition upon heating (fig. S5). When $Φ$ becomes 30 to 40%, a Re sphere falls at an irregular speed in the melt located between the solid grains. It indicates the change of rheological behavior from solid type to liquid type (48). The samples were quenched at the rheological transition or at a slightly higher $T$ on the basis of thermocouple reading. The extent of melting was determined for each sample by numerical analyses (Fig. 2 and fig. S7). To achieve a maximum precision, we correlated backscattering electron images and chemical maps obtained by SEM; melt regions present a maximum concentration of incompatible elements. A rapid increase of $Φ$ up to 30 to 40% melt within the first 100 K above the solidus was observed (Figs. 3 and 4). Then, $Φ$ increases less steeply at higher $T$. At high extent of melting, the grain size of the residual solid phases becomes

10 to 30 µm or more, with well-developed polygonal textures, especially for olivine, orthopyroxene, and bridgmanite.

### Liquidus temperatures

Liquidus $T$ can primarily be constrained on the basis of the extrapolation to 100% of the relationship between $Φ$ and $T$ (Fig. 3). Additional constraints come from $T$ measured for (i) very fast and linear trajectory of the falling sphere indicating the absence of solid grains [such as in (49)] and (ii) complete disappearance of x-ray diffraction peaks from the sample. (iii) We also monitored the reappearance of crystallites upon cooling and (iv) microstructure of a number of “trial and error” syntheses at variable $T$, indicating whether the sample has been fully molten or not. All evidences converge on a pyrolitic-mantle liquidus 275 ± 50 K higher than the solidus at all upper mantle depths. Our experiments are compatible with previous works showing that for the pyrolitic/peridotitic composition, majorite replaces olivine at the liquidus between 13 and 16 GPa, and ferropericlase replaces majorite at ~22 GPa (table S3) (19, 23, 46).

### Comparison between experimental melts and natural lavas

Chemical compositions of the melts were determined using the electron microprobe (Supplementary Materials, Figs. 5 and 6, fig. S8, and table S4). From the correlation between Al₂O₃ content and CaO/Al₂O₃ ratio in experimental melts (19, 23), three major trends arise: (i) High extent of melting yields low CaO/Al₂O₃ ratio and low Al₂O₃ content. The pyrolitic melt (i.e., Φ=100%) presents CaO/Al₂O₃ ~ 0.45 and Al₂O₃ ~ 4.5 wt%. (ii) Low P favors high Al₂O₃ content and (iii) high P favors high CaO/Al₂O₃ ratio. In this figure, we superimpose a compilation of the composition of numerous basalts, picrites, and komatiites of different ages and locations. The domains where the different LIPs and komatiites overlap with experimental melts provide strong constraints on their mantle origin (see the Supplementary Information for a discussion of the role of alteration on the CaO/Al₂O₃ ratio in LIPs and komatiites).

The comparison suggests a minimum depth of 550 to 600 km (20 to 23 GPa) to produce the old ADK (Fig. 6). Their compositions

<figure><img src="image_351.png" /><figcaption>**Fig. 3.** Temperature dependence of extent of melting. The volumetric extent of melting ($Φ$) is determined at various $P$ and $T$ on the basis of the joint analyses of electron microphotographs and chemical maps of the recovered samples. High melting extents are also constrained from mass balance calculations. Extrapolation to $Φ$=100% of the quasi-linear relation found at high $Φ$ yields a good estimate of the liquidus $T$ as a function of $P$ (distinguished by colors).</figcaption></figure>

<figure><img src="image_352.png" /><figcaption>**Fig. 4.** Melting diagram of the pyrolitic mantle. The solidus is determined up to 25 GPa using in situ electrical conductivity (red circles) and x-ray diffraction (pink circle). The rheological transition at 30 to 40% of melt (gray dashed line) and the liquidus are determined using various methods, including a new falling sphere technique (Supplementary Materials). Both solidus and liquidus present a change of slope at ~15 GPa, when the high-$P$ polymorphs of olivine become stable.</figcaption></figure>
<figure><img src="image_353.png" /><figcaption>**Fig. 5.** Variation diagrams showing oxide abundances of melt versus extent of melting (wt%). (Left) Chemistry of our melts at various pressures and degrees of mantle melting. Averages of several electron microprobe analyses in the melt regions of our samples are plotted as colored symbols, for *P* between 4 and 25 GPa. One SD is <5% for each melt composition. (Right) Chemistry of melt's literature at various pressures and degrees of mantle melting. Each color represents a different study for the pyrolitic/peridotitic upper mantle: (19, 21, 23, 39–44).</figcaption></figure>

overlap with melts at Φ values ranging between 10 and 25% and 20 and 50% for highest and lowest *P* (thus with decreasing CaO/Al₂O₃ ratio), respectively. The presence of residual Al-bearing garnet at moderate *P* can induce low Al content in the melt. However, the CaO/

Al₂O₃ ratio cannot be very high in this case because garnet also incorporates CaO easily. At lower mantle pressures, bridgmanite and davemaoite (Dvm; CaSiO₃-perovskite) are rich in Al and Ca, respectively. Melting of Dvm above solidus yields primary melts enriched
<figure><img src="image_354.png" /><figcaption>**Fig. 6.** Correlation between Al₂O₃ contents and CaO/Al₂O₃ ratio of experimental melts and natural lavas. Melt compositions are regrouped by experimental *P* ranges (colored symbols). We also report melt composition from previous experimental reports: Light or dark gray circles indicate studies reproducing the CaO/Al₂O₃ ratio from melting of peridotite/pyrolite under mantle relevant (23, 27, 39–43, 46, 79) or more specific (19–21) experimental conditions, respectively. Colored areas correspond to natural lavas of Barberton-type (ADK, light green), Munro-type (AUK, dark green), Gorgona-type (AEK), and ancient (blue) and recent (red) picrites and basalts [see (54)]. On the basis of the degree of partial melting (Φ) measured in our samples (Fig. 3), we draw contour lines corresponding to Φ ~ 20, 40, and 60% (dashed gray lines). Comparing the composition of LIPs/komatiites of different ages (fig. S1) and experimental melts at varying *P* and Φ indicates the mantle source for each type of lava.</figcaption></figure>

in CaO and poor in Al₂O₃. In these mantle conditions, melts with high MgO and SiO₂ contents and high CaO/Al₂O₃ ratios can be produced at relatively low extent of melting. These results are in relatively good agreement with previous works (22–25), including (50) that suggests that komatiite genesis involves a high extent of melting (>30%) at low *P* < 5 to 6 GPa, after adiabatic decompression of mantle diapirs that begin to melt in the transition zone (*P* > 15 GPa) or at even higher depth (*P* > 20 GPa). In all cases, *T*ₚ at the time of ADK formation should be appreciably higher than that prevailing today (1).

The intermediate position of AUK suggests mantle sources of 300- to 550-km depth (i.e., 10 to 20 GPa), as proposed previously [e.g., (19, 23, 27)]. To produce Al-undepleted melts, the presence of a residual garnet is not mandatory. There is a smooth transition between the ADK and AUK between 20 and 25 GPa, which explains the coexistence of AUK and ADK in some Munro-type outcrops. Because of similar sources, the two melts present similar major element contents, but the *P* difference can still induce substantial differences in minor and trace elements. Formation of the AUK occurred primarily during the Late Archean and could therefore correspond to lower *T* than the Barberton komatiites (11, 12). On the basis of our experiments at AUK, lavas are compatible with Φ from 15% to more than 60%, which makes AUK compositions quite close to bulk pyrolite. AUK can be produced in a wide range of *P* from ~5 to 10 GPa to ~25 to 30 GPa, when Φ is very high (>50%).

AEK, picrites, and basalts align along near-vertical trends at CaO/Al₂O₃ ratio between 0.5 and 1 and with Al₂O₃ contents varying from ~8 to ~17 wt%. AEK can be produced by melting at 200- to 400-km depth [i.e., 7 to 13 GPa, in agreement with (2, 19, 23)] and

picrites and basalts at 100 to 200 km (i.e., 4 to 5 GPa). Our results suggest moderate extent of melting (Φ < 50%) for recent Gorgona-type komatiites and picrites; melts enriched in Al₂O₃ suggest mantle melting without garnet in the solid residue. Last, picritic and basaltic lavas, which represent LIPS and recent hotspots such as Hawaii or Reunion Island, correspond to melting extents lower than 30% (2, 17).

The melt contents in other major elements confirm the trends presented above. We draw correlation diagrams between them, which show the effect of *P* on the melt composition (Fig. 5). The bulk sample composition, pyrolite, is a singular point in these diagrams because the melt composition evolves toward it with increasing Φ, irrespective of the experimental *P*. On these diagrams, we superimpose the compositions previously reported for basalts, picrites, and komatiites (fig. S8): (i) CaO versus SiO₂ diagram. The high SiO₂ content obtained at the highest experimental *P* and moderate Φ explains well the composition of ADKs. (ii) Al₂O₃ versus CaO and Al₂O₃ versus MgO diagrams. At relatively low Φ (i.e., at high CaO or low MgO content), the sequence of natural lavas (ADK, AUK, AEK, and picrites/basalts) is well compatible with the *P*-induced decrease of Al₂O₃ content in melts (Fig. 5). (iii) CaO versus FeO diagram. The increasing FeO content with increasing *P* at relatively low Φ is also well compatible with the sequence of natural lavas.

## DISCUSSION

### Mantle conditions for the origin of LIPs and komatiites

We now combine the evolution of Φ between the solidus and liquidus of pyrolite at various mantle depths (from Fig. 3) with the chemical overlap between experimental melts and natural lavas (from Fig. 6). The result (Fig. 7) highlights the likely source region for primary magmas at the origin of LIPs and komatiites. In addition, we plot, the current mantle geotherm (51) and its range of uncertainties (gray area), as well as mantle geotherms with *T*ₚ 100, 200, and 300 K higher than today (gray lines) (1). We acknowledge that these *T* profiles do not formally correspond to mantle *T* profiles; as they go through regions between the solidus and liquidus, the effect of latent heat of mantle fusion should be included. Still, the mantle depth at which an adiabatic profile crosses the solidus informs truthfully on the onset of mantle melting along a decompression path. It appears clear that *T*ₚ > 200 to 300 K higher than today induces melting in the *P* domain corresponding to ADK. Secular mantle cooling can explain the progressive disappearance of ADK, with preservation of AUK lavas when *T*ₚ was >150 to 200 K above that prevailing today (Fig. 1). Last, a mantle plume >50 to 150 K hotter than the average mantle would yield the production of picrites and deep basalts, which are currently found at oceanic islands (e.g., La Reunion) and mid-oceanic ridges, as well as AEK. These hot mantle plumes could still be uprising today (52), producing supervolcanoes like Toba (Indonesia) or Yellowstone (the United States).

Previous works were based on the modeling of mantle cooling using the Urey ratio (53) or thermodynamical modeling using the PRIMELT method (1, 3). Other less extreme scenarios exist on the basis of similar calculation from the Urey ratio or petrological modeling technique to derive primary magma compositions using large datasets of non-arc basaltic lavas (9, 54, 55). Our estimates of Archean temperatures are in better agreement with the most recent studies. These discrepancies can be explained by several critical points in the thermodynamical calculations. Different petrological models give different parameterizations for the MgO content in the melt
<figure><img src="image_355.png" /><figcaption>**Fig. 7.** Mantle P-T regions are sources of the different LIPs and komatiites. Solidus (blue line) and liquidus (red line) of pyrolite are reported together with the current mantle geotherm (gray zone), and *T* profiles that are 100, 200, and 300 K higher than the geotherm (light gray lines; these *T* profiles are not corrected for the effect of latent heat of fusion between the solidus and liquidus). Φ values between solidus and liquidus are constrained from our experiments (Fig. 3). The likely source region of the different types of lavas is determined on the basis of their composition (Figs. 5 and 6). The low-*P* basalts and picrites plotted in the figure correspond to flood basalts and OIBs found in the Phanerozoic. MgO-rich picrites (i.e., corresponding to higher *P* and melting extents) found on the islands of Gorgona and Curaçao present a chemistry close to komatiites and are therefore classified in the Gorgona-type group. The thicker blue line serves as a guide to the eyes for the ascent of a hot mantle plume that could be the origin of the oldest komatiites (see Fig. 8).</figcaption></figure>

produced by decompression melting, as a function of *P* and *T* [e.g., (56)]. The thermodynamical models rely heavily on the available experimental data, and before our work, these data were limited to *P* below 7 to 8 GPa concerning the formation of high MgO magmas (2, 23). We show that Archean komatiites can be produced at higher *P* than previously thought. A major reason of this discrepancy is our lower solidus profile; it yields differences of 100 to 200 K lower than previous works for the *T*<sub>p</sub> value of the primitive mantle. In addition, we provide an accurate determination of the rate of melting as a function of *P* and *T* for the anhydrous mantle. This yields a more accurate dataset on the melt chemistry over a wide range of mantle P-T. Still, we observe relatively good agreement with the literature, indicating that the key point is the accuracy of the *T* measurement within the sample, and the relationship between melting rate and *T* between solidus and liquidus.

**Which mantle thermal state produces primitive ultramafic magmas?**

The production of a large volume of mafic to ultramafic magma is generally explained by the rising of a hot mantle plume (29, 30). However, if the mantle *T* profile follows an adiabatic gradient [corresponding to *T*<sub>p</sub> of 1600 K or higher (57)], then substantial thermal anomalies can only be generated if the ascending material comes from the hot thermal boundary layer in contact with the core. Such

mechanism can explain the present-day OIB (28, 58); the hot material can eventually form a dome and rise toward the surface. Still, the thickness of the hot-layer atop the core-mantle boundary (CMB) could only be several hundred kilometers, yielding mantle plume with a reduced volume. These smaller plumes are unlikely to generate great LIPs.

*T*<sub>p</sub> higher in the past than that prevailing today would yield ubiquitous and persistent mantle melting in a wide range of mantle depths (Fig. 7). This thermal state is not consistent with the sporadic character of the LIPs. In addition, if the extent of melting should become larger than ~40%, as required to explain the composition of some komatiites, then this would indicate a layer of very low viscosity within the upper mantle for a very long time, which is a new concept with implications to be investigated. Last, warmer adiabatic *T* profiles with a mushy upper mantle are not found in line with the cooling history of Earth, when the starting point is the MO solidification after the Moon-forming impact [see e.g., (59)]. Alternatively, major thermal anomalies within the deep mantle could more simply explain the properties of LIPs and komatiites.

The simplest mechanism to produce large mantle plumes with positive *T* anomalies is a superadiabatic *T* profile in the deep solid mantle (Fig. 8). This favors mantle instabilities and overturns; large volumes of deep mantle can rise to the surface and become hotter than the surrounding mantle because they cool less during the adiabatic decompression than the superadiabatic gradient prevailing in the mantle. We recently demonstrated that super adiabatic gradients should prevail in the deep mantle for 1 Ga or more after the major Moon-forming impact (59). This thermal state is a necessary step between (i) the solidification of the global magma ocean, which yields a *T* profile close to the mantle solidus and (ii) the adiabatic profile of the present-day mantle with *T*<sub>p</sub> of 1600 to 1650 K (if an adiabatic *T* gradient is really achieved today). Along this cooling stage, *T* should decrease by 500 to 1000 K in the lowermost mantle. It implies the extraction of a huge amount of heat, especially because of other sources of energy available, such as radioactive disintegration.

**Superadiabatic mantle gradients could be at the origin of LIPs**

Following previous work (59), we recalculate the evolution of the mantle *T* profile over time after the crystallization of a pyrolitic magma ocean. A most appreciable difference in the parameters is the solidus of pyrolite, which presents a progressive increase with increasing depths in the upper mantle (fig. S6). Because of this curvature, the last mantle region to solidify upon cooling is at mid-upper mantle depths (fig. S9), instead of just below the surface for the chondritic-type mantle composition used in the previous work. It results that the cooling rate of the pyrolitic mantle is significantly slower, with a *T* profile remaining largely superadiabatic in the deep mantle for 2 to 3 Ga or more (Fig. 8 and fig. S9). Together with the fact that melting happens solely in the first 1000 km of the mantle (fig. S12), such long time scale provides a primary justification for using the properties of the pyrolite mantle for our discussions; while the bulk primitive mantle could very well be of chondritic-type composition after the core segregation, the mean mantle involved in the convection processes 2 to 3 Ga later in the shallow mantle is believed to be similar to pyrolite composition.

We also calculate the P-T path that uprising material would follow along adiabatic decompression paths, taking into account the latent heat of fusion [e.g., (50)] when the plume enters the P-T domain between the solidus and the liquidus. The adiabatic uprising
<figure><img src="image_356.png" /><figcaption>**Fig. 8. Thermal state in the ancient mantle.** Panel (A) reports the solidus and liquidus of the pyrolitic mantle (green lines), the present-day reference adiabatic T gradient for Tp = 1600 K (light blue), ancient adiabatic T gradients going through a mushy zone in the upper mantle (blue lines), and ancient T profiles with ages from 4.0 to 1.0 Ga before present (red lines). Panel (B) is a zoom of the figure A between 0 and 1000 km and 1500 to 2500 K. The T profiles result from modeling the cooling of a pyrolitic mantle composition from a global magma ocean, following previous work (59). The thicker blue line is a guide for eyes for the ascent of a hot mantle plume. With an initial temperature typical of an ancient T profile, for example, 2850 K at 2300-km depth in a ~3-Ga-old mantle [blue dot in (A)], a rapid decompression path along its corresponding adiabatic T gradient (bold blue line) yields a mantle plume with a major T excess in the shallow mantle, compared to the average T profile. This plume would cross the solidus at ~800-km depth and melt partially at levels of 25 to 30%, ~50%, and 30 to 40%, at mantle depths of 660 km, 300 to 500 km, and 100 to 200 km, respectively. This information is given by the position of the thick blue line between the solidus and liquidus in (B), following the expected extent of melting as a function of P and T (Fig. 3). These conditions are compatible with the formation of Barberton-type komatiites (Fig. 7). Our model of liquid komatiite petrogenesis by adiabatic decompression is similar to that of picrite or MORB petrogenesis, except that it involves a superadiabatic T gradient to generate hot mantle plume coming from great mantle depths.</figcaption></figure>

of a mantle blob from 2300-km depth and 2850 K (which corresponds to mantle conditions about 3 Ga ago; blue dot in Fig. 8) would cross the solidus at ~800-km depth and melt partially at levels of 25 to 30%, ~50%, and 30 to 40% at depths of 660 km, 300 to 500 km, and 100 to 200 km, respectively (fig. S8). These melting conditions are very compatible with those needed to produce ADK, AUK, and early AEK-picrite lavas, respectively (Fig. 6). Therefore, the prevalence of superadiabatic gradients in Archean and early Proterozoic mantle

offers new explanations for the sources of LIPs and komatiites. Within the first ~1000-km depth, the magnitude of the T excess of the mantle plume, compared to the average mantle, depends on the level of mantle superadiabaticity of the T profile in the deep mantle and the depth from which the mantle plume originates. With secular cooling, the level of superadiabaticity decreases, inducing smaller T anomalies in the shallow mantle.

Another major aspect that still needs to be addressed is the volume of uprising material, which needs to be several million cubic kilometers to produce major LIPs. While the presence of superadiabatic mantle gradients could provide an explanation for the presence of hot plumes in the upper mantle, the dynamics of such hot mantle could be different from today. A steep T profile in the past implies high Rayleigh number value and turbulent mantle convection. Unsteady convection cells could have developed with characteristic length scales much smaller than the mantle thickness (60). At early time, uprising of hot material from the deep Earth may have occurred sporadically along heat pipes (61). Mantle cooling over the Hadean and the Archean should have progressively increased the length scale of mantle convection, due to the decrease of the Rayleigh number, before global mantle convection is progressively established (62). Under these circumstances, the production of LIPs and komatiites and slab subduction would be precursors of the establishment of global mantle convection prevailing in modern Earth.

## MATERIALS AND METHODS

### Sample preparation

As starting material, we used a typical composition of pyrolitic mantle. This model is constructed by adding some MORB to mantle peridotite. It is dedicated to represent a primitive undepleted upper mantle (63). Our pyrolitic composition (table S1) is not far from KLB-1 (27) and fertile oceanic peridotites (64). A mixture of synthetic powders was homogenized in an agate mortar, before glass spheres of diameter ~2.5 mm were produced using an aerodynamic levitation system (65). The use of a regulated flow of pure Ar yields a convenient $Fe^{3+}/(Fe^{2+}+Fe^{3+})$ ratio of 7 to 8% for such type of sample composition (66).

We determine the water content in the pyrolitic glass using infrared spectroscopy, based on the Beer-Lambert law and absorption coefficients given for basaltic (67) and basanitic (68) glasses. These coefficients are valid when using the height of a water band located at, e.g., 3450 cm⁻¹. With water contents ranging between 18 and 22 parts per million, we consider the glass starting material nominally dry.

Glass spheres are reshaped by polishing in the form of a cylinder to be directly inserted in the high-*P* assembly of the multi-anvil press (MAP). The cylinder of glass offers several advantages including a rigid shape during compression and a very small surface/volume ratio for moisture absorption, compared to powders. We keep the sample length less than 1 mm to reduce *T* gradients during high *P* and *T* experiments.

### Electrical conductivity measurements

Electrical conductivity measurements are performed up to 20 GPa at the Laboratoire Magmas et Volcans using the same experimental procedure as in our previous work (45). We use tungsten carbide cubes with edge truncation of 11, 6, or 4 mm, octahedron *P* medium made of $Cr_2O_3$-doped MgO with edge length of 18, 14, or 10 mm,
respectively, and cylindrical furnace made of a Re foil of 25-µm thickness. The heater is surrounded by a zirconia sleeve, which serves as thermal insulator. The sample is inserted in an MgO capsule; MgO is softer than alumina parts, which minimizes the deformation of the molten sample. Before preparation of the assembly, all ceramic parts are fired at 1000°C for 1 hour and stored in a vacuum furnace at 150°C together with the glass sample. The sample length is determined before and after the experiment using a high-precision digital gauge and the SEM. We apply a correction based on elastic parameters to account for the change of sample length due to compression and thermal expansion. *P* is determined on the basis of a press-load versus sample-*P* calibration using the well-known phase transition in mantle minerals. The uncertainty evolves from 0.5 to 1 GPa with increasing *P* (69).

Rhenium disks electrodes with diameter of 1 mm are positioned at the top and bottom of the sample. They maintain the electrical connection between the sample and the electrodes and favor the spreading of the current at the sample surface. On one sample side, we insert the two wires of a C-type thermocouple (W<sub>0.95</sub>Re<sub>0.05</sub>-W<sub>0.74</sub>Re<sub>0.26</sub>) through the MgO cylinder and in contact with the Re disk. A third W<sub>0.95</sub>Re<sub>0.05</sub> wire is positioned at the other sample side. On the basis of this configuration, the sample *T* can be measured on one sample side and impedance spectra are collected between the two W<sub>0.95</sub>Re<sub>0.05</sub> wires. Because of the proximity with the Re furnace, we estimate our sample fO<sub>2</sub> close to the Re-ReO<sub>x</sub> buffer.

Absence of a metallic capsule around the silicate sample can induce the escape of some melt from the sample core, especially when achieving high extent of melting. This melt could react with the MgO sleeve. Fortunately, this effect is reduced when a stronger cohesion of the cell assembly is achieved at high *P* (70). In addition, we perform fast heating after the solidus *T* is detected (see below). Using the trial-and-error method, we optimized the experimental duration to obtain a good textural and chemical equilibrium together with minimum melt migration.

Electrical impedance of samples was recorded using the ModuLab MTS Impedance/Gain-Phase Analyzer in a range of frequencies between 1 MHz and 1 Hz [see details in (45)]. After compression to the target *P*, we performed different cycles of heating and cooling to remove adsorbed moisture from the assembly components, until a reproducible sample conductivity is achieved. Then, *T* is raised up to ~1300 K, which is above the glass transition of our pyrolite glass. The sample resistivity increases again as a sign of final dehydration and crystallization. The grain size is up to 2 to 10 µm, based on the microstructure of quenched samples. After this step, the sample is cooled to 600 to 800 K before the final heating up to the melting *T*. Our major criterion for the detection of solidus is the change of slope in the *T* dependence of sample conductivity (fig. S4).

## X-ray measurements

We performed x-ray diffraction and contrast imaging in situ in the MAP operating at the PSICHE beamline of the SOLEIL synchrotron using the same experimental procedure as in our previous work [see details in (45)]. *P* up to 27 GPa is generated using WC cubes with 3- or 4-mm truncation edges mounted with 7- or 10-mm-length Cr-doped MgO octahedra. High *T* is provided by a Re tubular furnace surrounded by a ZrO<sub>2</sub> sleeve for thermal insulation. Re electrodes with a donut shape are used to connect each furnace extremity to one cube truncation. Heating is obtained by applying an alternating current at 100 Hz.

*T* is monitored by a C-type thermocouple (W/Re) in contact with the furnace at the center of the assembly. All pieces are positioned horizontally in parallel with the x-ray beam.

The beamline configuration is similar to that for a previous study dedicated to measuring the silicate melt viscosity (49). At the sample position, we use either (i) a relatively large beam of about 2 mm by 2 mm to record the sample image as a function of time using a fast charge-coupled device camera (>100 frames per second) located behind the MAP or (ii) a beam collimated to 50 µm by 200 µm to characterize the sample mineralogy and monitor *P*. X-ray diffraction patterns are recorded in energy-dispersive mode using a Ge solid-state detector coupled to the Caesar-type diffractometer. A 2θ diffraction angle of 8.0° (±10⁻³) is convenient to register the major diffraction lines of our silicate samples. On the basis of uncertainties of the determination of the MgO cell parameter and on the MgO equation of state itself, the uncertainty in the *P* determination is ~0.3 GPa (71).

## The new falling sphere technique

We investigate the rheological transition occurring in the sample at an extent of melting (Φ) of 30 to 40% using x-ray imaging and an original “falling sphere” method (fig. S5). For this, using laser ablation, we drill a small hole on the lateral surface of the cylindrical sample at equal distance from the two ends. One or two rhenium spheres with diameters of 50 to 80 µm are inserted in the hole. We verify that Re spheres are located on the upper side of sample when inserting the assembly in the MAP. The sample *P* is determined on the basis of the x-ray diffraction signal of the MgO disks located on each side of the sample.

Upon heating above the solidus, we first observe some instabilities in the position of the sphere likely produced by an increasing Φ in the sample. With further *T* increase, the sphere starts to fall within the sample. Its movement appear to be chaotic, indicating the sphere bouncing on the surrounding mineral grains. The sphere does not fall at a constant speed and takes several seconds to go across the sample. Its complete fall takes orders of magnitude longer than a sphere fall through a fully molten sample at similar *P* (49). The slow falling of the sphere evidences a rheological transition between solid-like to liquid-like behavior. This transition occurs when the amount of melt is sufficient to disrupt the rigid cohesion between grains forming the solid matrix. Still, the high fraction of solid grains prevents a rapid fall of the sphere. We have heated some of our samples to even higher *T*, and in this case, the sphere falls through the entire sample in milliseconds (~200 ms). Using this technique, we have monitored the rheological transition at eight different *P* from 5 to 24 GPa and at *T* ranging from ~1855 to ~2250 K (table S2) The extent of melting is determined a posteriori on the recovered samples.

## Temperature gradients and uncertainties

Estimation of thermal gradients is crucial for an accurate determination of the pyrolite melting *T*. They can yield to heterogeneities in partially molten samples, although the sample volume is less than 1 mm³. Our recovered samples do not present significant zoning when using the 18/11 cell assembly, neither in Φ nor in the chemical composition of the different phases. The situation is not as ideal when using the smaller 14/6, 10/4, or 7/3 assemblies. We limit artifacts possibly generated by chemical segregations by duplicating experiments several times and crosschecking the petrological and chemical analyses of different samples. For a few samples, the texture and
sample mineralogy were slightly different at the two capsule ends because of non-negligible *T* gradients. If the overall chemical homogeneity is well preserved, then these samples can still be used to determine the melt composition, based on chemical analyses of the different sample areas.

We model the *T* field for our different cell assemblies for thermocouple readings between 1773 and 2273 K (fig. S10) (72). We obtain a *T* difference of less than 100 K between the sample and the thermocouple located in contact with (i) an electrode at the top of the sample for electrical conductivity experiments and (ii) the tubular furnace around the sample for x-ray experiments. The difference can only be higher if the sample is significantly off the center of the assembly (and especially for 10/4 and 7/3 assemblies), which can be checked on the recovered samples. The coldest zones are the top and bottom of the cylindrical sample, while the central plane of the sample records the highest temperatures. We applied a correction to all *T* measurements based on the modeled *T* gradients.

For conductivity measurements, the sample volume is ~0.8 mm³ and maximum *T* difference of 60 or 70 to 80 K is calculated using 18/11 or 14/6 and 10/4 assemblies, respectively. For x-ray measurements, maximum *T* difference is 80, or 90 K, for 10/4 assembly with a sample volume of 0.4 to 0.5 mm³, or 7/3 assembly with a sample of 0.2 to 0.3 mm³, respectively.

## Microstructure of recovered samples

Texture and mineralogy of recovered samples is characterized using the SEM. Samples are mounted in an epoxy resin and polished to mirror-like surface. Longitudinal section (along the furnace axis) is the most convenient to check the sample homogeneity and estimate the extent of melting and abundance of the different minerals. Observations are carried out using both JEOL (JSM-5910 LV at LMV) and field-emission gun ZEISS (supra 55VP at 2MAtech, Aubière, France) SEMs with accelerating voltage of 15 kV and working distance (WD) from 6 to 11.6 mm. The use of small WD enables the detection of small veins of melts typical of low Φ.

Observations confirm that the initial shape of the samples is well preserved, even at high melting extents. Images show a good textural equilibrium for all our samples, based on well-crystallized polygonal grains with diameters of several micrometers and grain boundary angles close to 120° (Fig. 2 and figs. S3 and S7). The particles size does not change significantly with *P* but rather with increasing *T*. The presence of cracks in some samples could result from tensile stresses during quench and decompression.

Several chemical maps (Si, Fe, Mg, Ca, Al, and K) have been recorded using energy-dispersive spectroscopy to image the elemental distribution between the melt and the solid residue. Mg, Al, and Si provide a rapid identification of major minerals; olivine and its high-*P* polymorphs, pyroxenes, and garnet. Fe, Ca, and K provide clues for the identification of the network of melt pockets.

## Determination of the extent of melting

We use high-resolution electron microphotographs to determine Φ in our samples. At the solidus *T*, samples show nanometric traces of melts at the grain triple junctions and along grain boundaries. In this case, the size of melt pockets or films is generally too small to enable a precise measurement of the melt composition, due to chemical interferences with the adjacent minerals. At *T* of 20 to 50 K or more above the solidus, the interconnected melt presents a very low dihedral angle (fig. S3). To prevent artifacts in the determination of

Φ, we used the MATLAB software to correlate the backscattered electron images with chemical maps of the incompatible elements, which present a maximum concentration in the melt regions. We perform several calculations with adaptable pixel sizes to check the robustness of our method; uncertainties on calculated Φ are <2 and <5% for low and high extents of melting, respectively (Fig. 2 and table S3). These Φ values are in very good agreement with other Φ values calculated on the basis of the chemical composition of all the phases coexisting in a same sample. Still, the mass balance technique provides larger uncertainties compared to image analyses. When a thermal gradient is detectable inside a sample, the extent of melting is studied in the different regions and the *T* in each zone is estimated from modeled *T* gradients (fig. S10).

## Chemical analyses

Chemical compositions of melts are determined by electron probe microanalyses. We used alternatively Cameca SX100 and SXFiveTactis electron probes with a same accelerating voltage of 15 kV and a beam current of 20 nA. Quantitative analyses are obtained for eight elements (Si, K, Fe, Na, Ti, Al, Ca, and Mg) with counting times of 20 s on peak maxima and 2 × 10 s on the background. As standards, we use natural minerals and synthetic oxides (Si and Ca, wollastonite; Mg, forsterite; Al, Al₂O₃; Fe, fayalite; Na, albite; K, orthoclase; Ti and Mn, MnTiO₃; Cr, Cr₂O₃; and Ni, NiO) and pure metal (Fe, Mn, Mg, Ni, and Cr) for metal alloys. Mineral phases and vein of melts are analyzed with a focused beam (2 to 5 µm), while large melt pools (with texture of quench) are analyzed with a beam defocused to diameter of ~10 µm. In this case, standards are also analyzed with the defocused beam.

Regions of melt could be adequately localized on the basis of chemical maps obtained using SEM. In this work, we do not attempt to analyze the incipient melt composition (the smallest Φ). Instead, we cover the most extensive P-T-Φ domain typical of current and primitive upper mantle. As logically expected, the melt composition evolves toward pyrolite, the sample bulk composition, with increasing Φ (Fig. 5). Chemical trends, however, vary significantly with experimental P. Noticeable trends are the following:

1) CaO and Al₂O₃ contents in melt decrease significantly with increasing *T* because of their known incompatible character. CaO presents a major decrease from 10 to 15 to ~3 wt% at all *P*, except at 23 to 25 GPa where it becomes slightly less incompatible. Al₂O₃ content also decreases largely with *T* from ~17 to ~5% at 4 to 5 GPa.

2) The incompatible character of Al disappears with increasing *P*, as expected from the presence of residual garnet above the solidus for *P* higher than ~7 GPa. At low Φ, the Al content in the melt increases from ~2 to ~18 wt% with *P* from 4 to 5 GPa to 23 to 25 GPa, respectively.

3) Our melt FeO contents are low, indicating a moderately incompatible behavior. They present a relative stagnation up to Φ ~ 60%, independently of the *P*, except at 23 to 25 GPa. It suggests solid-melt partition coefficients close to unity with the solid phases that disappear before Φ ~ 60% and an incompatible behavior only relative to the mineral phase that remain stable up the liquidus. We note that our FeO results diverge from previous works performed at very low melting extents, as significantly higher Fe content is reported in the incipient melt [e.g., (73, 74)].

4) MgO presents a major increase of its content from 15 to 25 wt% to ~40 wt% with increasing Φ. Still, its compatible character decreases with increasing *P*; its content in the melt at low Φ evolves
from ~15 to ~25 wt% with increasing P. At high Φ, melts produced at P higher than ~7 GPa contain 25 to 30% MgO, similar to komatiitic liquids.

5) SiO₂ contents in the melt stagnate at 43 to 48 wt% at all Φ and below 23 to 25 GPa. Therefore, upper mantle melts should not be particularly depleted or enriched in SiO₂ compared to pyrolite. The highest melt SiO₂ content is 50 to 52% at 23 to 25 GPa and Φ lower than 20%. This effect can be related to ferropericlase being the liquidus phase in this *P* range [as reported in (75) for peridotite].

6) Simultaneous low MgO and stagnant or high SiO₂ contents yield high SiO₂/MgO ratio in the first melts to form [e.g., (76)].

7) High Ti, Na, and K concentrations are measured in the first drops of melts (table S4), and these contents decrease fast with increasing Φ, as expected for strongly incompatible elements.

## Geodynamical modeling

Following our previous work (59), we calculate the time evolution of the pyrolitic-mantle geotherm during the first ~3 Ga of Earth's history. First, we set the primordial mantle *T* gradient parallel and slightly above the mantle solidus. It corresponds to the thermal state of the mushy mantle when the magma ocean has solidified to Φ of ~40% at all depths, several thousand years after the major Moon-forming impact (77). We note that this model could rise discussions; still, such highly superadiabatic *T* gradient is very likely to happen after the global magma ocean solidification. Then, we use a one-dimensional numerical model in spherical geometry accounting for turbulent convective heat transfer to calculate the extent of cooling at all mantle depths. This is the same calculations as recently reported for a chondritic-type mantle [see details in (59)].

We adopt the physical parameters typical of the pyrolitic mantle model. The composition establishes a number of thermodynamical parameters; e.g., equation of state, heat capacity, and latent heat of fusion, which are all calculated using a mixing model based on the molar fractions of enstatite (37.63%), forsterite (51.64%), fayalite (5.39%), anorthite (3.13%), and diopside (2.22%). Concerning the pyrolite solidus and liquidus, we used (i) the solidus and liquidus measured below ~27 GPa in this study (Fig. 4), (ii) the solidus measured recently in LMV using the laser-heated diamond anvil cell, which is very similar to the solidus of the chondritic-type mantle (78), and (iii) the same liquidus as reported for the chondritic-type mantle (78) but shifted by ~200 K to match the liquidus *T* of pyrolite found at ~25 GPa and ~2500 K in this study.

We observe that the liquidus *T* profile in the lower mantle does not have a substantial impact on the results of our calculations. When fitted using Simon-Glatzel equations $T_m(P) = T_0(P/a + 1)^{1/c}$, where $T_m(P)$ and *P* are melting *T*(*K*) and *P*(*GPa*), respectively, it yields (i) ($T_0 = 1622.7$; *a* = 237.79; *c* = 0.33615) and ($T_0 = 349.82$; *a* = 0.10937; *c* = 2.9646) for the pyrolite solidus at *P* below and above 24 GPa, respectively, and (ii) ($T_0 = 1931.2$; *a* = 222.27; *c* = 0.48148) and ($T_0 = 1498.9$; *a* = 10.760; *c* = 2.3668) for the pyrolite liquidus at *P* below and above 18 GPa, respectively (Fig. 8).

For mantle viscosity, we consider equation 19B of (59) with $A_{\text{diff}} = 10^{-15} \text{ Pa}^{-1} \text{ s}^{-1}$, $V_{\text{diff}} = 10^{-6} \text{ m}^3 \text{ mol}^{-1}$, $E_{\text{diff}} = 200 \text{ kJ mol}^{-1}$, and no dichotomy between upper and lower mantle viscosity. We also considered radiogenic heating as in our previous work.

In addition to the calculated pyrolitic-mantle geotherm as a function of time, we also calculate several mantle adiabatic *T* profiles for the deep mantle. They are primarily constrained by the pyrolite equation of state; however, complications arise at *T* between the

solidus and the liquidus, where the degree of mantle partial melting can change with depth. This is taken into account by considering the latent heat of fusion for the variation of Φ with depth, as already described in previous works [e.g., (50)]. It yields adiabatic *T* profiles steeper for the mushy mantle between the solidus and liquidus compared to subsolidus mantle conditions (Fig. 8).

## SUPPLEMENTARY MATERIALS

Supplementary material for this article is available at https://science.org/doi/10.1126/sciadv.abo1036

## REFERENCES AND NOTES

1. C. Herzberg, K. Condie, J. Korenaga, Thermal history of the Earth and its petrological expression. *Earth Planet. Sci. Lett.* **292**, 79–88 (2010).

2. C. Herzberg, M. J. O'Hara, Plume-associated ultramafic magmas of phanerozoic age. *J. Petrol.* **43**, 1857–1883 (2002).

3. C. Herzberg, P. D. Asimow, Petrology of some oceanic island basalts: PRIMELT2.XLS software for primary magma calculation. *Geochem. Geophys. Geosystems* **9**, Q09001 (2008).

4. C.-T. Lee, P. Luffi, T. Plank, H. Dalton, W. Leeman, Constraints on the depths and temperatures of basaltic magma generation on Earth and other terrestrial planets using new thermobarometers for mafic magmas. *Earth Planet. Sci. Lett.* **279**, 20–33 (2009).

5. T. Grove, S. Parman, Thermal evolution of the Earth as recorded by komatiites. *Earth Planet. Sci. Lett.* **219**, 173–187 (2004).

6. I. H. Campbell, R. W. Griffiths, The changing nature of mantle hotspots through time: Implications for the chemical evolution of the mantle. *J. Geol.* **100**, 497–523 (1992).

7. A. Piccolo, B. Kaus, R. White, R. Palin, Generation of earth's early continents from a relatively cool archean mantle. *Geochem. Geophys. Geosystems* **20**, 1679–1697 (2019).

8. S. Matthews, K. Wong, O. Shorttle, M. Edmonds, J. MacLennan, Do olivine crystallization temperatures faithfully record mantle temperature variability? *Geochem. Geophys. Geosystems* **22**, e2020GC009157 (2021).

9. J. Ganne, X. Feng, Primary magmas and mantle temperatures through time. *Geochem. Geophys. Geosystems* **18**, 872–888 (2017).

10. M. Yoshida, A new conceptual model for whole mantle convection and the origin of hotspot plumes. *J. Geodyn.* **78**, 32–41 (2014).

11. R. E. Ernst, *Large Igneous Provinces* (Cambridge University Press, 2014).

12. J. Arndt, S. J. Barnes, C. M. Lesher, *Komatiites* (Cambridge University Press, 2008).

13. J. M. D. Day, Hotspot volcanism and highly siderophile elements. *Chem. Geol.* **341**, 50–74 (2013).

14. S. A. Gibson, Major element heterogeneity in archean to recent mantle plume starting-heads. *Earth Planet. Sci. Lett.* **195**, 59–74 (2002).

15. B. Jahn, G. Gruau, A. Glikson, Komatiites of the Onverwacht Group, S. Africa: REE geochemistry, Sm/Nd age and mantle evolution. *Contrib. to Mineral. Petrol.* **80**, 25–40 (1982).

16. M. Ile Bas, IUGS reclassification of the high-mg and picritic volcanic rocks. *J. Petrol.* **41**, 1467–1470 (2000).

17. C. Herzberg, M. J. O'Hara, Phase equilibrium constraints on the origin of basalts, picrites, and komatiites. *Earth Sci. Rev.* **44**, 39–79 (1998).

18. D. Green, Genesis of archean peridotitic magmas and constraints on archean geothermal gradients and tectonics. *Geology* **3**, 15–18 (1975).

19. C. Herzberg, J. Zhang, Melting experiments on anhydrous peridotite KLB-1: Compositions of magmas in the upper mantle and transition zone. *J. Geophys. Res.* **101**, 8271–8295 (1996).

20. R. Nomura, K. Hirose, K. Uesugi, Y. Ohishi, A. Tsuchiyama, A. Miyake, Y. Ueno, Low core-mantle boundary temperature inferred from the solidus of pyrolite. *Science* **343**, 522–525 (2014).

21. K. Litasov, E. Ohtani, Phase relations and melt compositions in CMAS-pyrolite-H₂O system up to 25 GPa. *Phys. Earth Planet. Inter.* **134**, 105–127 (2002).

22. T. Inoue, R. P. Rapp, J. Z. Zhang, T. Gasparik, D. J. Weidner, T. Irfune, Garnet fractionation in a hydrous magma ocean and the origin of Al-depleted komatiites: Melting experiments of hydrous pyrolite with REEs at high pressure. *Earth Planet. Sci. Lett.* **177**, 81–87 (2000).

23. M. J. Walter, Melting of garnet peridotite and the origin of komatiite and depleted lithosphere. *J. Petrol.* **39**, 29–60 (1998).

24. M. Boyet, M. Garçon, N. Arndt, R. W. Carlson, Z. Konc, Residual liquid from deep magma ocean crystallization in the source of komatiites from the ICDP drill core in the Barberton Greenstone Belt. *Geochim. Cosmochim. Acta* **304**, 141–159 (2021).

25. C. Herzberg, Generation of plume magmas through time: An experimental perspective. *Chem. Geol.* **126**, 1–16 (1995).
26. A. C. Kerr, G. F. Marriner, N. T. Arndt, J. Tarney, A. Nivia, A. D. Saunders, R. A. Duncan, The petrogenesis of Gorgona komatiites, picrites and basalts: New field, petrographic and geochemical constraints. *Lithos* **37**, 245–260 (1996).

27. E. Takahashi, Melting of a dry peridotite KLB-1 up to 14 GPa—Implications on the origin of the peridotitic upper mantle. *J. Geophys. Res. Solid Earth Planets* **91**, 9367–9382 (1986).

28. V. Courtillot, A. Davaille, J. Besse, J. Stock, Three distinct types of hotspots in the Earth's mantle. *Earth Planet. Sci. Lett.* **205**, 295–308 (2003).

29. G. Foulger, J. Natland, D. Presnall, D. Anderson, *Plates, Plumes, and Paradigms* (Geological Society of America, 2005).

30. C. A. Dalton, C. H. Langmuir, A. Gale, Geophysical and geochemical evidence for deep temperature variations beneath mid-ocean ridges. *Science* **344**, 80–83 (2014).

31. S. Maruyama, M. Santosh, D. Zhao, Superplume, supercontinent, and post-perovskite: Mantle dynamics and anti-plate tectonics on the core–mantle boundary. *Gondw. Res.* **11**, 7–37 (2007).

32. N. Coltice, H. Bertrand, P. Rey, F. Jourdan, B. R. Phillips, Y. Ricard, Global warming of the mantle beneath continents back to the Archaean. *Gondw. Res.* **15**, 254–266 (2009).

33. L. T. Elkins-Tanton, Continental magmatism, volatile recycling, and a heterogeneous mantle caused by lithospheric gravitational instabilities. *J. Geophys. Res. Solid Earth* **112**, B03405 (2007).

34. C. Hemond, N. T. Arndt, U. Lichtenstein, A. W. Hofmann, N. Oskarsson, S. Steinthorsson, The heterogeneous Iceland plume: Nd-Sr-O isotopes and trace element constraints. *J. Geophys. Res.* **98**, 15833–15850 (1993).

35. F. Jourdan, A. Marzoli, H. Bertrand, S. Cirilli, L. Tanner, D. Kontak, G. McHone, P. Renne, G. Bellieni, $^{40}$Ar/$^{39}$Ar ages of CAMP in North America: Implications for the Triassic–Jurassic boundary and the $^{40}$K decay constant bias. *Lithos* **110**, 167–180 (2009).

36. M. G. Jackson, M. D. Kurz, S. R. Hart, R. K. Workman, New Samoan lavas from Ofu Island reveal a hemispherically heterogeneous high 3He/4He mantle. *Earth Planet. Sci. Lett.* **264**, 360–374 (2007).

37. J. Blichert-Toft, N. T. Arndt, Hf isotope compositions of komatiites. *Earth Planet. Sci. Lett.* **171**, 439–451 (1999).

38. T. Matsumoto, A. Seta, J. Matsuda, M. Takebe, Y. Chen, S. Arai, Helium in the Archaean komatiites revisited: Significantly high 3He/4He ratios revealed by fractional crushing gas extraction. *Earth Planet. Sci. Lett.* **196**, 213–225 (2002).

39. F. A. Davis, M. M. Hirschmann, M. Humayun, The composition of the incipient partial melt of garnet peridotite at 3 GPa and the origin of OIB. *Earth Planet. Sci. Lett.* **308**, 380–390 (2011).

40. M. B. Baker, E. M. Stolper, Determining the composition of high-pressure mantle melts using diamond aggregates. *Geochim. Cosmochim. Acta* **58**, 2811–2827 (1994).

41. R. J. Kinzler, Melting of mantle peridotite at pressures approaching the spinel to garnet transition: Application to mid-ocean ridge basalt petrogenesis. *J. Geophys. Res. B Solid Earth* **102**, 853–874 (1997).

42. E. Takahashi, T. Shimazaki, Y. Tsuzaki, H. Yoshida, Melting study of a peridotite KLB-1 to 6.5 GPa, and the origin of basaltic magmas. *Philos. Trans. Royal Soc. Phys. Sci. Eng.* **342**, 105–120 (1993).

43. K. Hirose, I. Kushiro, Partial melting of dry peridotite at high-pressures—Determination of compositions of melts segregated from peridotite using aggregates of diamonds. *Earth Planet. Sci. Lett.* **114**, 477–489 (1993).

44. D. Novella, D. J. Frost, The composition of hydrous partial melts of garnet peridotite at 6 GPa: Implications for the origin of group II kimberlites. *J. Petrol.* **55**, 2097–2123 (2014).

45. D. Andrault, G. Pesce, G. Manthilake, J. Monteux, N. Bolfan-Casanova, J. Chantel, D. Novella, N. Guignot, A. King, J.-P. Itie, L. Hennet, Deep and persistent melt layer in the Archaean mantle. *Nat. Geosci.* **11**, 139–143 (2018).

46. R. G. Trønnes, D. J. Frost, Peridotite melting and mineral–Melt partitioning of major and minor elements at 22–24.5 GPa. *Earth Planet. Sci. Lett.* **197**, 117–131 (2002).

47. C. Herzberg, P. Raterson, J. Zhang, New experimental observations on the anhydrous solidus for peridotite KLB-1. *Geochim. Cosmochim. Acta* **1**, 1051–104 (2000).

48. T. Scott, D. L. Kohlstedt, The effect of large melt fraction on the deformation behavior of peridotite. *Geochim. Cosmochim. Acta* **246**, 177–187 (2006).

49. L. Xie, A. Yoneda, D. Yamazaki, G. Manthilake, Y. Higo, Y. Tange, N. Guignot, A. King, M. Scheel, D. Andrault, Formation of bridgmanite-enriched layer at the top lower-mantle during magma ocean solidification. *Nat. Commun.* **11**, 548 (2020).

50. G. H. Miller, E. M. Stolper, T. J. Ahrens, The equation of state of a molten komatiite 2, Application to komatiite petrogenesis and the Hadean mantle. *J. Geophys. Res.* **96**, 11849–11864 (1991).

51. T. Katsura, A. Yoneda, D. Yamazaki, T. Yoshino, E. Ito, Adiabatic temperature profile in the mantle. *Phys. Earth Planet. Inter.* **183**, 212–218 (2010).

52. A. Davaille, A simultaneous generation of hotspots and superswells by convection in a heterogeneous planetary mantle. *Nature* **402**, 756–760 (1999).

53. J. Korenaga, Urey ratio and the structure and evolution of Earth's mantle. *Rev. Geophys.* **46**, RG2007 (2008).

54. K. C. Condie, R. C. Aster, J. van Hunen, A great thermal divergence in the mantle beginning 2.5 Ga: Geochemical constraints from greenstone basalts and komatiites. *Geosci. Front.* **7**, 543–553 (2016).

55. G. F. Davies, Effect of plate bending on the Urey ratio and the thermal evolution of the mantle. *Earth Planet. Sci. Lett.* **287**, 513–518 (2009).

56. H. Iwamori, D. McKenzie, E. Takahashi, Melt generation by isentropic mantle upwelling. *Earth Planet. Sci. Lett.* **134**, 253–266 (1995).

57. C. Herzberg, P. D. Asimow, N. Arndt, Y. Niu, C. M. Lesher, J. G. Fitton, M. J. Cheadle, A. D. Saunders, Temperatures in ambient mantle and plumes: Constraints from basalts, picrites, and komatiites. *Geochem. Geophys. Geosystems* **8**, 2006GC001390 (2007).

58. C. G. Farnetani, A. W. Hofmann, Dynamics and internal structure of the Hawaiian plume. *Earth Planet. Sci. Lett.* **295**, 231–240 (2010).

59. J. Monteux, D. Andrault, M. Guitreau, H. Samuel, S. Demouchy, A mushy Earth's mantle for more than 500 Myr after the magma ocean solidification. *Geophys. J. Int.* **221**, 1165–1181 (2020).

60. N. Arndt, A. Davaille, Episodic earth evolution. *Tectonophysics* **609**, 661–674 (2013).

61. W. B. Moore, A. A. G. Webb, Heat-pipe earth. *Nature* **501**, 501–505 (2013).

62. M. Maurice, N. Tosi, H. Samuel, A. C. Plesa, C. Huttig, D. Breuer, Onset of solid-state mantle convection and mixing during magma ocean solidification. *J. Geophys. Res. Planets* **122**, 577–598 (2017).

63. A. E. Ringwood, in *International Series in the Earth's and Planetary Sciences*, (McGraw-Hill, 1975), pp. 189–194.

64. W. F. McDonough, S. S. Sun, The composition of the earth. *Chem. Geol.* **120**, 223–253 (1995).

65. L. Hennet, V. Cristiglio, J. Kozaily, I. Pozdnyakova, H. E. Fischer, A. Bytchkov, J. W. E. Drewitt, M. Leydier, D. Thiaudiere, S. Gruner, S. Brassamin, D. Zanghi, G. J. Cuello, M. Koza, S. Magazu, G. N. Greaves, D. L. Price, Aerodynamic levitation and laser heating: Applications at synchrotron and neutron sources. *Eur. Phys. J. Spec. Top.* **196**, 151–165 (2011).

66. D. Andrault, M. Muñoz, G. Pesce, V. Cerantola, A. Chumakov, I. Kantor, S. Pascarelli, R. Rüffer, L. Hennet, Large oxygen excess in the primitive mantle could be the source of the great oxygenation event. *Geochem. Perspect. Lett.* **6**, 5–10 (2018).

67. M. Mercier, A. Di Muro, N. Métrich, D. Giordano, O. Belhadj, C. W. Mandeville, Spectroscopic analysis (FTIR, Raman) of water in mafic and intermediate glasses and glass inclusions. *Geochim. Cosmochim. Acta* **74**, 5641–5656 (2010).

68. T. A. Shishkina, R. E. Botcharnikov, F. Holtz, R. R. Almeev, A. M. Jazwa, A. A. Jakubiak, Compositional and pressure effects on the solubility of H₂O and CO₂ in mafic melts. *Chem. Geol.* **388**, 112–129 (2014).

69. A. Boujibar, D. Andrault, M. A. Bouhifid, N. Bolfan-Casanova, J.-L. Devidal, N. Trcera, Metal-silicate partitioning of sulphur, new experimental and thermodynamic constraints on planetary accretion. *Earth Planet. Sci. Lett.* **391**, 42–54 (2014).

70. A. Pommier, K. Leinenweber, D. L. Kohlstedt, C. Qi, E. J. Garnero, S. J. Mackwell, J. A. Tyburczy, Experimental constraints on the electrical anisotropy of the lithosphere-asthenosphere system. *Nature* **522**, 202–206 (2015).

71. Y. Tange, E. Takahashi, Y. Nishihara, K. I. Funakoshi, N. Sata, Phase relations in the system MgO-FeO-SiO₂ to 50 GPa and 2000 degrees C: An application of experimental techniques using multianvil apparatus with sintered diamond anvils. *J. Geophys. Res. Solid Earth.* **114**, B02214 (2009).

72. J. Hernlund, K. Leinenweber, D. Locke, J. A. Tyburczy, A numerical model for steady-state temperature distributions in solid-medium high-pressure cell assemblies. *Am. Mineral.* **91**, 295–305 (2006).

73. K. N. Matsukage, Z. Jing, S. I. Karato, Density of hydrous silicate melt at the conditions of Earth's deep upper mantle. *Nature* **438**, 488–491 (2005).

74. P. Condamine, D. J. Frost, in *Annual Report* (Bayerisches Geoinstitute, 2018).

75. E. Ito, A. Kubo, T. Katsura, M. J. Walter, Melting experiments of mantle materials under lower mantle conditions with implications for magma ocean differentiation. *Phys. Earth Planet. Inter.* **143–144**, 397–406 (2004).

76. C. Liebske, D. J. Frost, Melting phase relations in the MgO-MgSiO₃ system between 16 and 26 GPa: Implications for melting in Earth's deep interior. *Earth Planet. Sci. Lett.* **345**, 159–170 (2012).

77. J. Monteux, D. Andrault, H. Samuel, On the cooling of a deep terrestrial magma ocean. *Earth Planet. Sci. Lett.* **448**, 140–149 (2016).

78. D. Andrault, N. Bolfan-Casanova, G. Lo Nigro, M. A. Bouhifid, G. Garbarino, M. Mezouar, Melting curve of the deep mantle applied to properties of early magma ocean and actual core-mantle boundary. *Earth Planet. Sci. Lett.* **304**, 251–259 (2011).

79. W. Wang, E. Takahashi, Subsolidus and melting experiments of K-doped peridotite KLB-1 to 27 GPa: Its geophysical and geochemical implications. *Most* **105**, 2855–2868 (2000).

80. H. S. Smith, A. J. Erlank, Geochemistry and petrogenesis of komatiites from the Barberton greenstone belt, South Africa. *Komatiites*, 347–397 (1982).

81. C. Herzberg, Depth and degree of melting of komatiites. *J. Geophys. Res.* **97**, 4521–4540 (1992).

82. Y. Lahaye, N. Arndt, Alteration of a komatiite flow from Alexo, Ontario, Canada. *J. Petrol.* **37**, 1261–1284 (1996).

83. M. J. de Wit, R. A. Hart, R. J. Hart, The Jamestown Ophiolite Complex, Barberton mountain belt: A section through 3.5 Ga oceanic crust. *J. African Earth Sci* (1983) **6**, 681–730 (1987).
84. D. R. Pyke, A. J. Naldrett, O. R. Eckstrand, Archean ultramafic flows in Munro Township, Ontario. *GSA Bull.* **84**, 955–978 (1973).

85. B. Aitken, L. Echeverría, Petrology and geochemistry of komatiites and tholeiites from Gorgona Island, Colombia. *Contrib. to Mineral. Petrol.* **86**, 94–105 (1984).

86. E. Nisbet, N. Arndt, M. Bickle, W. Cameron, C. Chauvel, M. Cheadle, E. Hegner, A. Martin, R. Renner, E. Roedder, Uniquely fresh 2.7 Ga Komatiites from the Belingwe Greenstone-Belt, Zimbabwe. *Geology* **15**, 1147–1150 (1987).

87. M. Bickle, A. Martin, E. G. Nisbet, Basaltic and peridotitic komatiites and stromatolites above a basal unconformity in the Belingwe greenstone belt, Rhodesia. *Earth Planet. Sci. Lett.* **27**, 155–162 (1975).

88. Q. Xie, R. Kerrich, J. Fan, HFSE/REE fractionations recorded in three komatiite-basalt sequences, Archean Abitibi greenstone belt: Implications for multiple plume sources and depths. *Geochim. Cosmochim. Acta* **57**, 4111–4118 (1993).

89. E. Nisbet, M. Bickle, A. Martin, The mafic and ultramafic lavas of the belingwe greenstone belt, Rhodesia. *J. Petrol.* **18**, 521–566 (1977).

90. J. G. Arth, N. T. Arndt, A. J. Naldrett, Genesis of archean komatiites from Munro Township, Ontario: Trace-element evidence. *Geology* **5**, 590–594 (1977).

91. J. A. Hallberg, D. A. C. Williams, Archaean mafic and ultramafic rock associations in the Eastern Goldfields region, Western Australia. *Earth Planet. Sci. Lett.* **15**, 191–200 (1972).

92. S. S. Sun, R. W. Nesbitt, Petrogenesis of Archaean ultrabasic and basic volcanics: Evidence from rare earth elements. *Contrib. to Mineral. Petrol.* **65**, 301–325 (1978).

93. R. Nesbitt, S. S. Sun, A. C. Purvis, Komatiites: Geochemistry and genesis. *Can. Mineral.* **17**, 165–186 (1979).

94. B.-M. Jahn, B. Auvray, S. Blais, R. Capdevila, J. Cornichet, F. Vidal, J. Hameurt, Trace element geochemistry and petrogenesis of finnish greenstone belts. *J. Petrol.* **21**, 201–244 (1980).

95. R. Nesbitt, S.-S. Sun, Geochemistry of Archaean spinifex-textured peridotites and magnesian and low-magnesian tholeiites. *Earth Planet. Sci. Lett.* **31**, 433–453 (1976).

96. D. A. C. Williams, R. G. Furnell, Reassessment of part of the Barberton type area, South Africa. *Precambrian Res.* **9**, 325–347 (1979).

97. G. Gruau, C. Chauvel, N. T. Arndt, J. Cornichet, Aluminum depletion in komatiites and garnet fractionation in the early Archean mantle: Hafnium isotopic constraints. *Geochim. Cosmochim. Acta* **54**, 3095–3101 (1990).

98. T. Chaudhuri, M. Satish-Kumar, R. Mazumder, S. Biswas, Geochemistry and Sm-Nd isotopic characteristics of the Paleoarchean Komatiites from Singhbhum Craton, Eastern India and their implications. *Precambrian Res.* **298**, 385–402 (2017).

99. N. Arndt, Komatiites, kimberlites, and boninites. *J. Geophys. Res.* **108**, 2293 (2003).

100. N. Arndt, A. Fowler, Textures in komatiites and variolitic basalts, in *The Precambrian Earth: Tempos and Events* (Elsevier, 2004).

101. P. A. Sossi, S. M. Eggins, R. W. Nesbitt, O. Nebel, J. M. Hergt, I. H. Campbell, H. S. C. O'Neill, M. Van Kranendonk, D. R. Davies, Petrogenesis and geochemistry of archean komatiites. *J. Petrol.* **57**, 147–184 (2016).

102. Y. Lahaye, N. Arndt, G. Byerly, C. Chauvel, S. Fourcade, G. Gruau, The influence of alteration on the trace-element and Nd isotopic compositions of komatiites. *Chem. Geol.* **126**, 43–64 (1995).

103. S. Sun, (1984).

104. F. F. Guo, S. Svetov, W. D. Maier, E. Hanski, S. H. Yang, Z. Rybnikova, Geochemistry of komatiites and basalts in Archean greenstone belts of Russian Karelia with emphasis on platinum-group elements. *Miner. Depos.* **55**, 971–990 (2020).

105. M. J. Bickle, N. T. Arndt, E. G. Nisbet, J. L. Orpen, A. Martin, R. R. Keays, R. Renner, Geochemistry of the igneous rocks of the Belingwe greenstone belt: Alteration, contamination and petrogenesis, in *The Geology of the Belingwe Greenstone Belt, Zimbabwe* (CRC Press, 2020) pp. 175–213.

106. I. Vlastélic, T. Staudacher, M. Semet, Rapid change of lava composition from 1998 to 2002 at Piton de la Fournaise (Réunion) inferred from Pb isotopes and trace elements: Evidence for variable crustal contamination. *J. Petrol.* **46**, 79–107 (2005).

107. P. M. Holm, R. C. O. Gill, A. K. Pedersen, J. G. Larsen, N. Hald, T. F. D. Nielsen, M. F. Thirlwall, The tertiary picrites of West Greenland: Contributions from 'Icelandic' and other sources. *Earth Planet. Sci. Lett.* **115**, 227–244 (1993).

108. J. F. G. Wilkinson, H. D. Hensel, The petrology of some picrites from Mauna Loa and Kilauea volcanoes, Hawaii. *Contrib. to Mineral. Petrol.* **98**, 326–345 (1988).

109. S. Révillon, N. T. Arndt, E. Hallot, A. C. Kerr, J. Tarney, Petrogenesis of picrites from the Caribbean Plateau and the North Atlantic magmatic province. *Lithos* **49**, 1–21 (1999).

110. A. Kerr, J. Mahoney, Oceanic plateaus: Problematic plumes, potential paradigms. *Chem. Geol.* **241**, 332–353 (2007).

111. A. C. Kerr, *Oceanic Plateaus* (Elsevier, ed. 2, 2014), vol. 4.

112. Z. Zhang, J. J. Mahoney, J. Mao, F. Wang, Geochemistry of picritic and associated basalt flows of the western emeishan flood basalt province, China. *J. Petrol.* **47**, 1997–2019 (2006).

113. D. Wang, S. Karato, Electrical conductivity of talc aggregates at 0.5 GPa: Influence of dehydration. *Phys. Chem. Miner.* **40**, 11–17 (2013).

114. H. Brandstätter, I. Hanzu, M. Wilkening, Myth and reality about the origin of inductive loops in impedance spectra of lithium-ion electrodes—A critical experimental approach. *Electrochim. Acta* **207**, 218–223 (2016).

115. D. Klotz, Negative capacitance or inductive loop?—A general assessment of a common low frequency impedance feature. *Electrochem. Commun.* **98**, 58–62 (2019).

116. T. Yoshino, A. Shimojuku, S. Shan, X. Guo, D. Yamazaki, E. Ito, Y. Higo, K. I. Funakoshi, Effect of temperature, pressure and iron content on the electrical conductivity of olivine and its high-pressure polymorphs. *J. Geophys. Res. Solid Earth* **117**, 1–15 (2012).

117. T. Yoshino, Laboratory electrical conductivity measurement of mantle minerals. *Surv. Geophys.* **31**, 163–206 (2010).

118. J. Zhang, C. Herzberg, Melting experiments on anhydrous peridotite KLB-1 from 5.0 to 22.5 GPa. *J. Geophys. Res.* **99**, 17729–17742 (1994).

**Acknowledgments:** We thank two anonymous reviewers for their thoughtful comments. We thank L. Hennet for the preparation of the pyrolitic glasses and L. Xie, Q. Dumont, and M. Garçon for help and fruitful discussions. **Funding:** This research was financed by the French Government Laboratory of Excellence initiative no. ANR-10-LABX-0006 and the Université Clermont Auvergne. This is the laboratory of excellence Clervolc contribution no. 561.

**Author contributions:** R.P. prepared all experiments. R.P. and G.M. performed the conductivity measurements. R.P., D.A., G.M., N.G., A.K., and L.H. performed the experiments at SOLEIL. R.P. and J.L.D. performed the sample characterization and chemical analyses of melts. J.M. and D.A. elaborated the geodynamic modeling. R.P. and D.A. wrote the manuscript, which was revised by G.M. and J.M. and validated by all coauthors. **Competing interests:** The authors declare that they have no competing interests. **Data and materials availability:** All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials.

Submitted 13 January 2022

Accepted 17 August 2022

Published 2 November 2022

10.1126/sciadv.abo1036
Deep mantle origin of large igneous provinces and komatiites

Rmy Pierru, Denis Andrault, Geeth Manthilake, Julien Monteux, Jean Luc Devidal, Nicolas Guignot, Andrew King, and
Laura Henry

Sci. Adv., 8 (44), eabo1036.

DOI: 10.1126/sciadv.abo1036

View the article online

https://www.science.org/doi/10.1126/sciadv.abo1036

Permissions

https://www.science.org/help/reprints-and-permissions

Use of this article is subject to the [Terms of service](#).
# 8 Cas des satellites de glace

La plupart des processus physiques évoqués lors des premières étapes de la formation des planètes se produisent également dans la partie extérieure du Système Solaire, lors de la formation des satellites de glace. Dans la majorité de ces corps, le manteau est principalement composé de glace d'eau, tandis que le noyau est essentiellement constitué de matériaux rocheux. Cependant, des différences de composition et de structure interne existent entre les principaux satellites glacés de Jupiter et de Saturne, suggérant des histoires d'accrétion et de différenciation distinctes (par exemple Tobie et al. (2013) et Monteux et al. (2018)). Ces satellites se seraient formés dans des disques circum-planétaires à partir d'agrégats de glace et de roche (Canup and Ward, 2002; Estrada et al., 2008). Bien que la variété des structures internes des lunes glacées soit encore débattue, le rôle essentiel de l'accrétion est largement accepté pour comprendre les différences entre les lunes glacées.

## 8.1 Impact géant et déformation du noyau d'Encelade

Malgré sa petite taille (R=252 km), Encelade, une des lunes de Saturne, est l'un des corps géologiquement les plus actifs du système solaire. Son activité endogène surprenante est caractérisée par une province très active au pôle sud, d'où des éruptions de vapeur d'eau et de grains de glace émanant de crêtes tectoniques chaudes ont été observées par la sonde Cassini (Porco et al., 2006; Hansen et al., 2006; Waite et al., 2006; Spencer et al., 2006). Cette activité est associée à une énorme puissance thermique estimée entre 5 et 15 GW (Spencer and Nimmo, 2013), ce qui implique un intérieur chaud, compatible avec une couche d'eau liquide sous une couche de glace et un intérieur différencié (Nimmo and Kleine, 2007; Schubert et al., 2007). Les modèles de dissipation des marées peuvent expliquer pourquoi l'activité est concentrée aux pôles, où la dissipation devrait être maximale (Tobie et al., 2008; Běhounková et al., 2010). Cependant, il n'y a toujours pas d'explication satisfaisante sur la raison pour laquelle cette activité est localisée uniquement au sud, et non au nord.

Il a été proposé que la dichotomie entre les hémisphères nord et sud soit le résultat d'une asymétrie dans la forme du noyau rocheux McKinnon (2013). En raison de la faible pression et de la température modérée attendues dans le noyau d'Encelade, de grandes anomalies de topographie peuvent en effet être conservées sur de très longues périodes et peuvent expliquer pourquoi les activités de convection dans la couche de glace sont confinées uniquement au pôle sud (Showman et al., 2013). Outre la dépression polaire sud, les anomalies de topographie du noyau pourraient expliquer, au moins en partie, l'existence d'autres grandes dépressions observées à des latitudes modérées et non corrélées à des limites géologiques (Schenk and McKinnon, 2009). McKinnon (2013) a proposé trois hypothèses pour expliquer l'irrégularité possible du noyau rocheux d'Encelade : la fusion lors de l'accrétion de la région externe d'Encelade associée à une instabilité de degré un, l'accrétion de proto-lunes glacées autour de morceaux de roche irréguliers et l'assemblage par collision de deux proto-lunes précédemment différenciées.
<figure><img src="image_357.png" /></figure>

FIGURE 30 – Morphologie du noyau rocheux en fonction de la taille de l'impacteur et de la vitesse d'impact ($v_{esc}$ = 240 m/s). Pour chaque morphologie, le cercle rouge représente la forme sphérique avant impact du noyau impacté. Au-dessus de la courbe noire en tirets, la topographie induite par l'impact est négative. En dessous de cette ligne théorique critique, la topographie induite par l'impact est positive. Au-dessus de la courbe noire en pointillés, des noyaux très fortement déformés se forment et la fluidisation acoustique peut contribuer à leur forme finale. Cependant la déformation est trop importante et probablement non compatible avec la morphologie d'Encelade. Nous limitons le temps de simulation post-impact à une heure, ce qui signifie que pour les grandes vitesses d'impact ($\ge$ 6 km/s) et les grands rayons d'impact ($\ge$ 75 km), le matériau rocheux excavé du noyau d'Encelade et en orbite autour de la lune se déplace toujours avec une vitesse significative à la fin de la simulation. (D'après Monteux et al. (2016b)).

Dans une étude développée au LPG Nantes, nous avons étudié les conséquences des collisions avec de grands impacteurs sur la forme du noyau. Nous avons effectué des simulations d'impact à l'aide du code iSALE2D (Wünnemann et al., 2006) en considérant de grands impacteurs différenciés ayant un rayon compris entre 25 et 100 km et des vitesses d'impact comprises entre 0.24 et 2.4 km/s. Nos simulations ont montré que les principaux paramètres contrôlant la forme du noyau rocheux d'Encelade après l'impact sont le rayon et la vitesse de l'impacteur et, dans une moindre mesure, la présence d'un océan d'eau interne ainsi que la porosité
et la résistance du noyau rocheux. Pour les impacts de faible énergie, les impacteurs ne tra-
versent pas complètement le manteau glacé. L'enfoncement et l'étalement ultérieurs du noyau
rocheux de l'impacteur entraînent une anomalie topographique positive du noyau. Dans le cas
d'impacts modérément énergétiques, les impacteurs pénètrent complètement dans le manteau
glacé, ce qui induit une topographie négative du noyau entourée d'une anomalie positive de
moindre amplitude. La profondeur et l'étendue latérale de la zone excavée sont principalement
déterminées par le rayon et la vitesse de l'impacteur. Dans le cas d'impacts très énergétiques,
le noyau rocheux est fortement déformé et le corps entier est susceptible d'être perturbé. L'ex-
plication de la forme irrégulière du noyau d'Encelade sur de grandes longueurs d'onde par des
impacts impliquerait de multiples collisions à faible vitesse ($<$ 2.4 km/s) avec des impacteurs
différenciés déca-kilométriques, ce qui n'est possible qu'après la période LHB.

8.2 Accrétion et évolution thermique des satellites de glace

Des différences de composition et de structure interne existent entre les principaux satel-
lites glacés de Jupiter et de Saturne, suggérant des histoires d'accrétion et de différenciation
distinctes (Kirk and Stevenson, 1987; Mueller and McKinnon, 1988; Mosqueira and Estrada,
2003; Barr and Canup, 2008). Le facteur de moment d'inertie élevé déduit des mesures gra-
vimétriques de Galileo (C/MR² = 0.346) (Anderson et al., 2001) suggère que la séparation
glace-roche pourrait être incomplète à l'intérieur de Callisto, une des lunes de Jupiter. En re-
vanche, Ganymède a un moment d'inertie beaucoup plus faible (C/MR² = 0.31) (Anderson
et al., 2001) et montre des signes d'activité endogène passée (Pappalardo et al., 2004). Une
séparation complète de la glace et de la roche est suggérée pour Ganymède, ainsi que la for-
mation d'un noyau métallique, à l'origine d'un champ magnétique intrinsèque relativement
intense (Kivelson et al., 1998).

Avec une taille et une masse similaires, Titan, lune de Saturne, pourrait être un cas inter-
médiaire entre Callisto et Ganymède. Son moment d'inertie estimé à partir des mesures de
gravité de Cassini (Iess et al., 2010, 2012) suggère que l'intérieur de Titan est plus différencié
que Callisto mais probablement beaucoup moins que Ganymède. Comme Callisto, Titan pour-
rait encore posséder une couche de mélange glace-roche entre un noyau rocheux et un manteau
extérieur riche en glace, à moins que le noyau rocheux ne soit principalement composé de mi-
néraux fortement hydratés (Sohl et al., 2010; Castillo-Rogez and Lunine, 2010). Le fait que
l'intérieur de Callisto et peut-être de Titan puissent encore contenir une couche de mélange
glace-roche suggère que le satellite a pu éviter une fusion importante lors de l'accrétion et de
l'évolution ultérieure.

Il est possible d'éviter la fusion si l'accumulation de l'énergie d'accrétion est inefficace,
c'est-à-dire si l'énergie est rayonnée à un taux comparable au taux d'accrétion (Schubert et al.,
1981; Squyres et al., 1988; Kossacki and Leliwa-Kopystyński, 1993; Coradini et al., 1995;
Grasset and Sotin, 1996; Barr and Canup, 2008, 2010). L'échelle de temps d'accrétion devrait
être supérieure à 1 Myr pour éviter une fusion importante et donc une différenciation de Cal-
listo, tandis qu'une échelle de temps d'accrétion aussi courte que 1000-10 0000 ans pourrait être possible pour Ganymède. Cependant, ces échelles de temps dépendent de la façon dont le dépôt de chaleur et le refroidissement sont traités.

En utilisant une approche unidimensionnelle comme initialement développée pour l'accrétion des planètes terrestres (Safronov, 1978; Kaula, 1979), la croissance d'un corps peut être modélisée comme un dépôt de couches successives de matière. Dans ce cas, l'effet thermique d'un impact n'est pas considéré individuellement, mais est moyenné sur toute la surface et intégré. Cette approche est valable tant que les impacteurs restent petits (> 1 km) et sont distribués de manière aléatoire à la surface. Cela pourrait être le cas au tout début du processus d'accrétion, mais les impacteurs de plus de 1 km sont probablement devenus de plus en plus abondants à la fin de la phase d'accrétion (Estrada et al., 2009). On peut également s'attendre à des impacteurs de plus de 100 km (Sekine and Genda, 2012; Dwyer et al., 2013). Pour des impacts aussi volumineux, il est nécessaire de disposer d'une description détaillée de chaque impact, y compris du dépôt et du transfert d'énergie.

Afin de mieux comprendre l'évolution thermique des satellites de glace en croissance, nous avons développé un modèle numérique tri-dimensionnel à partir du code numérique Oedipus développé par Choblet et al. (2007). Ce modèle numérique retrace, à partir d'une population réaliste d'impacteurs, l'évolution thermique d'un satellite de glace pendant son accrétion. Mes résultats montrent que la durée de l'accrétion ou les tailles caractéristiques des impacteurs influent peu sur l'état thermique du corps pendant son accrétion (Monteux et al., 2014). Cependant, l'efficacité avec laquelle l'énergie cinétique de l'impacteur est transmise au corps impacté sous forme de chaleur est un paramètre clé qui gouverne le taux de fusion final du satellite et ainsi sa structure interne. Par ailleurs, mes modèles montrent qu'il est probable que Callisto soit actuellement différenciée.

## 8.3 Le rôle de l'accrétion dans l'évolution de l'atmosphère de Titan

Titan, le plus gros satellite de Saturne, est aussi le seul satellite du système solaire à posséder une atmosphère massive. La masse et la composition de l'atmosphère massive de Titan, qui est actuellement dominée par $N_2$ et $CH_4$, ont probablement varié tout au long de son histoire en raison d'une combinaison de processus exogènes et endogènes. Les processus de formation de cette atmosphère riche en azote et méthane sont encore assez mal contraints (Tobie et al., 2013). Des études en laboratoire ont permis de quantifier la quantité de $N_2$ dégazé lors d'un impact sur un objet riche en $NH_3$ et de mettre en évidence l'importance des paramètres d'accrétion sur l'évolution de l'épaisseur de cette atmosphère (Sekine et al., 2011).

J'ai développé avec N. Marounina un modèle numérique pour déterminer la quantité de gaz qui peut être dégagée lors d'impacts géants se produisant au cours de l'accrétion de Titan et pendant le premier milliard d'années d'évolution (Marounina et al., 2015). Dans cette étude, nous étudions l'évolution de l'atmosphère de Titan pendant le late heavy bombardment (LHB) en modélisant la perte et l'apport compétitifs de volatils par les impacts cométaires
<figure><img src="image_358.png" /></figure>

FIGURE 31 – Coupes équatoriales du champ de température (colonnes de gauche et du milieu) et représentations topographiques 3D (à droite) du satellite de glace en croissance en fonction du temps (de haut en bas). La colonne de gauche représente l'évolution d'un cas avec "accrétion froide" où, jusqu'à la fin du régime intermédiaire, $\gamma_{li} = \gamma_{lay} = 0.1$, $x_{m,li} = 10\%$ tandis que la colonne du milieu représente d'un cas avec "accrétion chaude" où $\gamma_{li} = \gamma_{lay} = 0.3$ et $x_{m,li} = 33\%$. L'échelle de couleur de la température est saturée en blanc pour la température au point de fusion de la glace (> 273 K). Entre chaque régime (précoce, intermédiaire, tardif), le champ de température est interpolé sur une grille plus large. Dans le régime tardif, $\gamma_{li} = 0.1$, $\gamma_{lay} = 0.3$ ; $x_{m,li} = 33\%$, $r_{min} = 10$ km et $r_{max} = 100$ km pour les colonnes de gauche et du milieu (D'après Monteux et al. (2014)).

et leurs conséquences sur l'équilibre atmosphérique. Nous modélisons le dégazage et le re-
<figure><img src="image_359.png" /><figcaption>FIGURE 32 – (a) Distribution du NH₃ à la surface de Titan après un épisode de LHB (b) distribution des tailles d'impact cumulées pendant ce LHB et (c) taux de dégazage à la surface de Titan après le LHB pour la distribution illustrée dans le panel (a) (D'après Marounina et al. (2015).</figcaption></figure>

froidissement de la croûte suite à des impacts géants et tout en tenant compte de l'érosion de l'atmosphère lors des impacts (Svetsov, 2007; Shuvalov, 2009). Pour des albédos de surface compris entre 0.1 et 0.7, nous avons examiné l'émergence d'une atmosphère pendant le LHB ainsi que l'évolution d'une atmosphère primitive de différentes masses et compositions avant cet événement, en tenant compte de la conversion NH3-N2 de la croûte induite par l'impact et du dégazage subséquent ainsi que de l'érosion atmosphérique induite par l'impact. En considérant une population d'impacteurs caractéristique du LHB, nous avons montré que la génération d'une atmosphère riche en N2 avec une masse équivalente à celle d'aujourd'hui nécessite une fraction de masse d'ammoniac de 2-5%, dépendant des albédos de surface, dans une couche glacée d'au moins 50 km sous la surface, impliquant un intérieur indifférencié au moment du LHB. A l'exception des albédos de surface élevés (≥ 0.7) où la plupart du N2 libéré reste gelé à la surface, nos calculs indiquent que les impacts à grande vitesse ont conduit à une forte érosion atmosphérique. Pour un Titan différencié avec une fine croûte enrichie en ammoniac (65 km) et un albédo inférieur à 0.6, toute atmosphère préexistante au LHB devrait être plus de 5 fois plus massive qu'actuellement, afin de maintenir une atmosphère équivalente à celle d'aujourd'hui. Ceci implique qu'une atmosphère massive s'est formée sur Titan lors de son accrétion ou que l'atmosphère riche en azote a été générée après le LHB.
Référence personnelle

Tobie, G., Lunine, J., **Monteux, J.**, Mouzis, O., Nimmo, F., 2013, *The origin and evolution of Titan. In: Mueller-Wodarg, I., Griffith, C., Cravens, T., Lellouch, E. (Eds.), *Titan: Surface, Atmosphere and Magnetosphere. LIEN HAL*.

**Monteux, J.**, G., Tobie, G., Choblet, M., Le Feuvre, 2014b, *Can large icy Moons accrete undifferentiated, Icarus, Volume 237, Pages 377-387*. **LIEN HAL**.

Marounina, N., Tobie, G., Carpy, S., **Monteux, J.**, Charnay, B., Grasset, O., 2015, *Evolution of Titan's Atmosphere during the Late Heavy Bombardment, Icarus*, 257, 324-335. **LIEN HAL**.

**Monteux, J.**, Collins G.S., Tobie G., Choblet G., 2016b, *Consequences of large impacts on Enceladus' core shape. Icarus vol.264, p.300-310*. **LIEN HAL**.

**Monteux, J.**, Golabek G.J., Rubie D.C., Tobie G., Young E.D., 2018, *Water and the interior structure of terrestrial planets and icy bodies. Space Science Reviews*. **LIEN HAL**.
# The Origin and Evolution of Titan

G. TOBIE, J.I. LUNINE, J. MONTEUX, O. MOUSIS AND F. NIMMO

## Abstract

Saturn's moon Titan is the second largest natural satellite in the solar system, and the only one that possesses a substantial atmosphere. The origin of the massive nitrogen atmosphere and the source of the present-day methane, which is continuously destroyed by UV-driven photochemistry, has long been debated. Data collected by the *Cassini-Huygens* mission since its arrival to Saturn in 2004 now provide key constraints on Titan's atmosphere composition, surface morphology and interior structure, which restrict the possible scenarios of formation and evolution. In the present chapter, after reviewing our present knowledge about the interior structure and composition based on the *Cassini-Huygens* data, we present the different physical and chemical processes that have potentially affected the origin and evolution of Titan. In particular, we argue that all along Titan's evolution, from accretion to present, interactions between the icy shell, the internal water ocean and the deep rock-dominated interior may have affected the evolution of the atmosphere. In spite of considerable progresses in our understanding of Titan thanks to *Cassini-Huygens*, many questions still remain, and we conclude on how they may be addressed by future exploration missions.

## 2.1 INTRODUCTION

Although Titan is similar in terms of mass and size to Jupiter's moons, Ganymede and Callisto, it is the only one harboring a massive atmosphere. Moreover, unlike the Jovian system populated with four large moons, Titan is the only large moon around Saturn. The other Saturnian moons are much smaller and have an average density at least 25% less than Titan's uncompressed density and much below the density expected for a Solar composition (Johnson and Lunine, 2005), although with a large variation from satellite to satellite. Both Jupiter's and Saturn's moon systems are thought to have formed in a disk around the growing giant planet. However, the difference in architecture between the two systems probably reflects different disk characteristics and evolution (e.g. Sasaki et al., 2010), and in the case of Saturn, possibly the catastrophic loss of one or more Titan-sized moons (Canup, 2010). Moreover, the presence of a massive atmosphere on Titan as well as the emission of gases from Enceladus' active south polar region (Waite

et al., 2009) suggest that the primordial building blocks that comprise the Saturnian system were probably more volatile-rich than Jupiter's.

The composition of the present-day atmosphere, dominated by nitrogen, with a few percent methane and lesser amounts of other species, probably does not directly reflect the composition of the primordial building blocks and is rather the result of complex evolutionary processes involving internal chemistry and outgassing, impact cratering, photochemistry, escape, crustal storage and recycling, and other processes. The low $^{36}$Ar/$N_2$ ratio measured by the GCMS on *Huygens* (Niemann et al., 2005, 2010) suggest that nitrogen was not brought to Titan in the form of $N_2$, but rather in the form of $NH_3$ (Owen, 1982). The argument goes as follows: Ar and $N_2$ have similar volatility and affinity with water ice. Thus, if the primary carrier of Titan's $N_2$ were molecular nitrogen itself, the Ar/$N_2$ ratio on Titan should be within an order of magnitude of the solar composition ratio of about 0.1 (Lunine and Stevenson, 1985), which is about 500,000 times larger than the observed ratio. This indicates that nitrogen has been brought in the form of easily condensable compounds, most likely ammonia, which then has been converted *in situ* by impact-driven chemistry (McKay et al., 1988; Sekine et al., 2011; Ishimaru et al., 2011) or by photochemistry (Atreya et al., 1978). As will be further discussed in the chapter, both conversion processes required restrictive conditions to occur, and therefore specific evolutionary scenarios for Titan.

Photochemical processes are also known to lead to an irreversible destruction of methane, implying that a source of replenishment must exist to explain its few percent abundance (vertically averaged) in the atmosphere. In the absence of such a replenishment, the atmospheric methane should disappear over a timescale of several tens of millions of years. Before the *Cassini-Huygens* mission, the atmospheric methane was proposed to be in equilibrium with a global surface liquid reservoir (e.g. Lunine et al., 1983). Even though liquid reservoirs have been identified on Titan's surface by *Cassini* in the form of lakes and seas at high latitudes (Stofan et al., 2007), their estimated total volume is too small to sustain methane in the atmosphere on geologic timescales (Lorenz et al., 2008a). This suggests that a subsurface reservoir of methane exists and replenishes the atmospheric reservoir episodically or continuously. Different sources have been proposed, but despite a number
of new constraints from the *Cassini-Huygens* mission, we lack a strong discriminant among them. The ¹²C/¹³C ratio measured in methane by the GCMS on *Huygens* (Niemann et al., 2010) is similar to the ¹²C/¹³C ratio measured in several solar system objects, and close to the solar value, possibly reflecting an atmospheric methane inventory that has not been affected by atmospheric escape. The analysis of Mandt et al. (2009) indicates that if the current atmospheric methane is the remnant of methane injected into the atmosphere more than about two hundred million years ago the ¹²C/¹³C should be decreased by at least 10% compared to the initial value. The current atmospheric methane may then have been recently delivered to the surface, from whatever reservoir. The detection of a significant amount of ⁴⁰Ar (the decay product of ⁴⁰K initially contained in silicate phases) further indicates that exchanges with the rocky interior have occurred (Niemann et al., 2010).

In this chapter we argue that exchange between the interior and atmosphere of Titan is likely to have occurred throughout Titan's evolution, at least during accretion and in the following billion years, and possibly until the present. A better understanding of these exchange processes requires good knowledge of the interior structure. Based on its mean density of 1881 kg.m⁻³, Titan has long been known to be composed of a roughly 60%-40% mixture of rock and ice by mass. Gravimetric and altimetric data collected by the *Cassini* spacecraft permit us to constrain the distribution of these components in the interior. However, as is further discussed in Section 2, the solution remains non-unique and the different possible structures of the interior imply a range of formation and evolution scenarios. Other data from the *Cassini-Huygens* mission suggest the presence of a water ocean, possibly doped with ammonia or other agents depressing the ocean crystallization point, beneath the solid icy shell. The existence of such an internal water ocean at present has important consequences for the evolution of Titan as well as for its astrobiological potential.

In Section 2.2, we review the present-day constraints on the structure and dynamics of Titan's interior and present different possible interior models. The present-day composition and structure of the atmosphere and surface are already discussed in other chapters and will not be repeated here. A good knowledge of the present-day state of Titan is a necessary preliminary to understanding its origin and evolution. Section 2.3 presents some theoretical constraints on the formation of the Saturnian system and describes the accretion process(es) of Titan. Coupled evolution of the interior, surface and atmosphere, including internal differentiation and generation and recycling of the atmosphere, are then discussed in Section 2.4. Finally, we conclude in Section 2.5 by listing a series of remaining questions and how they may be addressed by future exploration missions.

## 2.2 PRESENT-DAY INTERIOR STRUCTURE AND DYNAMICS

Models of Titan's formation (Section 2.3) and subsequent evolution (Section 2.4) are constrained by its present-day interior structure and dynamics. In this Section, we focus on three questions which are of particular relevance to Titan's history.

1. *To what extent is Titan differentiated?* The degree of differentiation depends mainly on the thermal conditions during and just after accretion. Thus, Titan's degree of differentiation constrains the conditions under which it formed and evolved during the first billion years.

2. *Is the ice shell convecting?* The long-term thermal evolution of Titan's interior is determined by the rate at which heat can be transferred across its icy outer layer(s). Thus, whether or not convection is taking place controls how rapidly Titan cools. If an ocean is present, the lifetime of the ocean (and thus its habitability) is controlled by the behaviour of the overlying ice shell. The ice shell also controls transport between the interior and the atmosphere, as well as potential tectonic and cryovolcanic activities.

3. *Does Titan possess a sub-surface ocean?* If an ocean is present, Titan's astrobiological potential is increased. Moreover, an ocean can have profound geophysical consequences; for instance, it will likely increase the rate of tidal dissipation in the overlying ice shell.

Prior to the *Cassini-Huygens* mission, there were essentially no observational constraints on any of these questions. Below we discuss the data sets relevant to the internal structure and dynamics of Titan.

### 2.2.1 Interior structure inferred from the gravity field data

Measurements of Titan's gravity field are made by precise radio tracking of the *Cassini* spacecraft as it passes close to Titan. Unfortunately, because of the lack of a scan platform, gravity tracking precludes the acquisition of most other observations, therefore the number of flybys dedicated to these measurements is limited. Furthermore, Titan's extended atmosphere can result in non-gravitational forces on the spacecraft, which complicates the analysis and requires that the closest approach distance is farther from Titan than what is desirable for an optimal resolution of the gravity field. Despite these obstacles, Titan's degree-two gravity coefficients have been determined with relatively high confidence, as well as the degree-three coefficients but with much less confidence (Iess et al., 2010). The degree-two spherical harmonic coefficients, which are the dominant coefficients for a body subjected to tidal and rotational forces, $J_2 = (-C_{20})$ and $C_{22}$ are tabulated in Table 2.1. Because of the number of flybys, the coefficients were determined in-
dependently, without requiring a value of 10/3 for the ratio of $J_2$ to $C_{22}$. This situation differs from bodies where fewer flybys have been carried out, such as Callisto (Anderson et al., 2001) and Rhea (Mackenzie et al., 2008).

The degree-two coefficients provide very important constraints on the radial mass distribution within the body. It can be shown mathematically that they are directly related to the principal moments of inertia (MoI), A, B and C:

$$C_{20}^{\text{obs}} = \frac{(A+B)/2 - C}{M_T \bar{a}^2} \quad ; \quad C_{22}^{\text{obs}} = \frac{B - A}{4M_T \bar{a}^2} \quad (2.1)$$

with $M_T$ Titan's mass and $\bar{a}$ the mean equatorial radius. Thus, measuring the $C_{20}$ and $C_{22}$ coefficients only gives the difference between the principal moments of inertia. In order to determine their absolute value, a third condition is needed. For Mars and the Earth, measurement of the precession rate (which gives $(C - \bar{A})/C$, where C and $\bar{A}$ are the polar and averaged equatorial moments of inertia) allows retrieval of the exact values of the three principal moments of inertia. In the case of Titan, such a measurement is not available. To infer the moment of inertia, one approach is to use the Radau-Darwin approximation, which assumes that the body responds to tidal and rotational forces as a fluid body (i.e. it has no long-term strength, and therefore it is in hydrostatic equilibrium). This approximation assumes that the body surface as well as each interface between internal layers correspond to an iso-potential surface. The potential at each interface is the sum of the tidal potential, centrifugal potential and gravitational potential. If the body is in hydrostatic equilibrium, the degree-two gravitational potential is proportional to the tidal and centrifugal potential. The response of the body to the tidal and centrifugal degree-two potential is classically described with the fluid Love number $k_2^f$. In hydrostatic equilibrium,

$$C_{20}^{\text{th}} = -\left(\frac{\alpha}{2} + \frac{q}{3}\right) k_2^f \quad ; \quad C_{22}^{\text{th}} = \frac{\alpha}{4} k_2^f \quad (2.2)$$

where $\alpha = M_S/M_T \times (R_T/D_{ST})^3$, is the tidal parameter, and $q = \omega^2 R_T^3/GM_T$ is the rotational parameter (with $M_S$ Saturn's mass, $R_T$ Titan's mean radius, $D_{ST}$ Saturn-Titan mean distance, $\omega$ Titan's spin rate). In the case of a synchronously rotating body in a circular orbit, $\alpha = q$ and therefore $C_{20}/C_{22} = -10/3$ if the body is hydrostatic. For Titan, $\alpha = q = 3.9555 \times 10^{-5}$.

In the Radau-Darwin approximation, the normalized moment of inertia about the spin axis, $C/M_T R_T^2$ can then be calculated from the fluid Love number, inferred from the gravity coefficients assuming that they result from a fluid response of the body and in the limit of small density variations:

$$\frac{C}{M_T R_T^2} = \frac{2}{3} \left[ 1 - \frac{2}{5} \left( \frac{5}{k_2^f + 1} - 1 \right)^{1/2} \right] \quad (2.3)$$

In the case of Titan, the ratio between the observed $C_{20}^{\text{obs}}$ and $C_{22}^{\text{obs}}$ is close to -10/3, suggesting that Titan is relatively close to hydrostatic equilibrium. We can therefore use the

<table><caption>Table 2.1. Observed ($C_{lm}$) un-normalized degree-two gravity coefficients. Observed values are from Iess et al. (2010), where SOL1 and SOL2 refer to the multi-arc and global solutions, respectively, and errors are given as $\pm 2\sigma$. Un-normalized surface topography coefficients $C_{lm}^h$ are derived from data detailed in Zebker et al. (2009b) and quoted from Nimmo and Bills (2010). Normalized moment of inertia $C/M_T R_T^2$ was calculated assuming hydrostatic equilibrium using equation 2.3.</caption><thead><tr><th></th><th>C<sub>lm</sub> × 10<sup>-6</sup><br>SOL1</th><th>C<sub>lm</sub> × 10<sup>-6</sup><br>SOL2</th><th>C<sub>lm</sub><sup>h</sup><br>m</th></tr></thead><tbody><tr><td>l = 2, m = 0</td><td>-31.80±0.80</td><td>-33.46±1.26</td><td>-358</td></tr><tr><td>l = 2, m = 2</td><td>9.98±0.08</td><td>10.02±0.14</td><td>63</td></tr><tr><td>ratio</td><td>-3.19±0.08</td><td>-3.34±0.13</td><td>-5.68</td></tr><tr><td>k<sub>2</sub><sup>f</sup>(C<sub>20</sub>)</td><td>0.965</td><td>1.015</td><td>-</td></tr><tr><td>k<sub>2</sub><sup>f</sup>(C<sub>22</sub>)</td><td>1.010</td><td>1.013</td><td>-</td></tr><tr><td>C/M<sub>T</sub> R<sub>T</sub><sup>2</sup>(C<sub>20</sub>)</td><td>0.335</td><td>0.342</td><td>-</td></tr><tr><td>C/M<sub>T</sub> R<sub>T</sub><sup>2</sup>(C<sub>22</sub>)</td><td>0.341</td><td>0.342</td><td>-</td></tr></tbody></table>

degree two gravity coefficients to first order to estimate the fluid Love number and the corresponding moment of inertia. Depending on whether we use the $C_{20}$ or $C_{22}$ coefficients as a reference to determine the fluid Love number, the normalized moment of inertia varies between 0.335 and 0.342. This may be compared with the value of 0.4 for a uniform sphere, 0.394 for the Moon, 0.33 for the Earth, 0.31 for Ganymede and 0.355 for Callisto (assuming hydrostatic equilibrium (Anderson et al., 1996, 2001)). The moment of inertia factor of Titan is therefore intermediate between Ganymede and Callisto, suggesting a modest increase of density toward the center.

As already mentioned, this interpretation however assumes that the $C_{20}^{\text{obs}}$ or $C_{22}^{\text{obs}}$ coefficient is entirely determined by the dynamical distortion of the satellite under the action of tidal and rotational forces. In reality, part of the gravity signal may be attributed to non-hydrostatic contributions such as uncompensated relief, ice shell thickness variations or lateral density variations. The fact that the $C_{20}^{\text{obs}}/C_{22}^{\text{obs}}$ ratio is not exactly -10/3 and that the degree three gravity coefficients are not zero clearly indicates that there are some non-hydrostatic contributions to the observed signals and that the gravity field of Titan is not only attributed to the dynamical distortion of a radially-distributed mass interior. The degree-three coefficients indicate that non negligible lateral variations in the mass distribution exist on Titan. These should affect the estimation of the moment of inertia. If we assume that the non-hydrostatic contributions to degree two are comparable to the inferred degree 3 terms, a few to 10% of the observed degree two signals may be attributed to the non-hydrostatic (Iess et al., 2010, Supplementary Information). According to the Radau-Darwin equation, a 5%-10% overestimate of the hydrostatic parts would shift $C/M_T R_T^2$ from 0.341 to 0.334-0.327.
In summary, the interpretation of the gravity data indicates that the normalized moment of inertia of Titan is of the order of 0.33-0.34. This MoI value combined with the mean density of Titan indicates that the density of the core is low, which may be explained by incomplete separation of rock from ice (less et al., 2010) or the presence of highly hydrated silicate minerals (Castillo-Rogez and Lunine, 2010; Fortes, 2011). Based on the present data, the presence of a small iron core (< 500 km) cannot be ruled out (Castillo-Rogez and Lunine, 2010; Fortes, 2011), but the low internal temperatures (< 900 - 1000K) required to maintain either an ice-rock mixture or highly hydrated silicate minerals are incompatible with the segregation of iron and the formation of an iron core. The most likely present-day structures is either a pure rock core, mostly anhydrous, surrounded by a mixture of ice and rock, surrounded by pure water ice (Fig. 2.5, top and middle), or a low density rock core consisting mainly, but not entirely, of hydrated silicate minerals surrounded by H₂O (Fig. 2.5, bottom). In both cases, a liquid water layer may be present between a high-pressure ice mantle and an outer ice shell. We will discuss the possible evolution scenarios leading to these interior structures in more details in Section 2.4.

## 2.2.2 Global shape: Constraints on the structure and thermal state of the ice shell

Just as with gravity, the degree-two shape of a synchronous body can be used to determine its moment of inertia - as long as the body is in hydrostatic equilibrium. As with gravity, the degree-two shape coefficients of such a body, $C_{20}^h$ and $C_{22}^h$, will be in the ratio -10/3. Thus, with sufficiently good topography we can check whether or not the hydrostatic condition is satisfied.

Titan's shape has been measured (Zebker et al., 2009b) by a combination of traditional radar altimetry (Zebker et al., 2009a), and a new technique which uses overlaps in the SAR-imaging swaths to determine topography (Stiles et al., 2009). In contrast to the results obtained from gravity, $C_{20}^h$ and $C_{22}^h$ were not in the expected -10/3 ratio (see Table 2.1).

This result is apparently paradoxical: how can Titan have nearly hydrostatic gravity but non-hydrostatic topography? The hydrostatic gravity rules out scenarios in which Titan's shape is derived from an earlier epoch (a "fossil bulge" similar to that of the Moon). One possibility is that Titan's outer ice shell - presuming the existence of an ocean - is of varying thickness or density and is (Airy or Pratt) isostatically compensated. An isostatic shell would give rise to only very small gravity anomalies (see below), while the shell thickness variations would result in a long-wavelength non-hydrostatic shape.

Nimmo and Bills (2010) explored a version of the Airy scenario, in which shell thickness variations were caused by spatial variations in tidal heating in the ice shell (cf. Ojakangas and Stevenson, 1989). They found that the observed topography could be fitted with a mean ice shell

of thickness $D \approx 100$ km, heated from below at a rate of 4-5 mW m⁻² (appropriate to a chondritic Titan) and from within at a rate about one-tenth as large (due to tidal dissipation). The tidal Love number $h_2$, which describes the response of the body to the eccentricity tides, was found to be $\approx 1.2$, consistent with theoretical predictions (Sohl et al., 2003). Given the relatively modest rate of tidal dissipation, the e-folding timescale to damp Titan's eccentricity is quite long (2-3 Gyr). Alternatively, Choukroun and Sotin (2012) proposed that the flattened shape of Titan may result from the accumulation of dense hydrocarbon clathrates at high latitudes owing to enhanced ethane precipitation. Their Pratt isostatic compensation model may explain the ~300 m difference between expected and measured polar radii if the ethane-rich clathrate crust is about 3-km thick.

Although the variable-thickness or variable-density shell is assumed isostatic, it will still have a small effect on the degree-2 gravity coefficients. Correcting for this effect before using equation (2.3) changes the resulting moment of inertia by about 2% (the MoI factor is still between 0.33 and 0.34). An isostatic shell is probably a good assumption: a shell with an elastic thickness of 50 km would still be about 95% compensated at degree-2 (Turcotte et al., 1981), assuming a Young's modulus of 9 GPa.

Perhaps the most surprising result of these studies is that, if it is correct, Titan's ice shell must be conductive, or at most weakly convecting, in order to support thickness or density variations on the long term. In the thickness-variation model, surface topography generates pressure gradients at depth; if the ice shell were convecting, the ice viscosity would be low enough to flow in response to these pressure gradients and progressively remove the surface topography (Stevenson, 2000). Whether or not an ice shell convects depends on whether its Rayleigh number $Ra$ of the whole ice layer exceeds some critical value $Ra_{cr}$ (e.g. Mitri and Showman, 2008):

$$Ra = \frac{\rho g \alpha \Delta T D^3}{\kappa \eta_b} \approx 3 \times 10^8 \left( \frac{D}{100 \text{km}} \right)^3 \left( \frac{10^{14} \text{ Pa s}}{\eta_b} \right) \geq Ra_{cr} \quad (2.4)$$

Here $\rho, \alpha, \kappa, \eta_b$ are the density, thermal expansivity, thermal diffusivity and basal viscosity of the ice shell, respectively, $\Delta T$ is the temperature contrast across the shell ($\Delta T \approx 170 - 180$ K between the surface and the ocean) and numerical values assumed are from Nimmo and Bills (2010). The value of $Ra_{cr}$ depends on how sensitive viscosity is to temperature changes, and the geometry of the shell (Solomatov, 1995; Barr and McKinnon, 2007), and for Titan is $\approx 4 \times 10^7$ in the stagnant lid regime. Equation (2.4) shows that whether or not the shell convects is a sensitive function of the viscosity of the base of the ice shell. For pure ice, the viscosity in turn depends mainly on the temperature of the ocean beneath, and secondarily on the ice grain size. The ocean temperature is controlled by how much ammonia is present; a sufficiently ammonia-rich ocean and/or large ice grain sizes would prevent convection because the ice would
be too viscous (e.g. Tobie et al., 2005b).

Nevertheless, it is possible that the ice shell was convective in the past and that convection ceased as the ocean crystallized and the ammonia concentration increased owing to reduced ocean volume (Tobie et al., 2005b; Mitri and Showman, 2008). The elevated orbital eccentricity of Titan indicates however that this convective period should have lasted less than about 0.5-1 billion of years (Tobie et al., 2005b). In contrast to the Galilean moons where the Laplace resonance excites the eccentricity of Europa and Io, there is no major resonance in the Saturnian system to force Titan's eccentricity. Unless some unknown process has recently forced the eccentricity, the latter should have monotonically decayed as a function of time since the formation of the system owing to tidal dissipation at the surface and in the interior (Sears, 1995; Sohl et al., 1995; Tobie et al., 2005b). The dissipation rate is reduced either if the interior is totally rigid (with no internal ocean) (Sohl et al., 1995) or if the outer ice shell above the ocean is thin and/or cold (i.e. conductive) (Tobie et al., 2005b). Indeed, the volume wherein tidal energy is dissipated is much bigger in a convective ice shell (where about half of the shell has a temperature larger than 240 K) than in a conductive ice shell where the warm temperature is confined at the bottom of the ice shell.

If Titan has an internal ocean at present, it implies that the ice shell remained conductive and relatively thin (< 30 – 40 km) during most of Titan's evolution, in order to reduce the eccentricity damping (Tobie et al., 2005b). This is possible if the viscosity of the icy shell is large enough due to large grain sizes or to low oceanic temperatures (owing to ammonia Tobie et al. (2005b) or methanol Deschamps et al. (2010)) and/or if the icy shell contains large fraction of low conductivity materials such as methane clathrate (Tobie et al., 2006). Following this requirement, the icy shell has never been convecting or has been convecting only during a limited period of time. A possible scenario may be that the ice shell started convecting during the last roughly billion years and then stopped later because of ocean crystallization as the concentration in ammonia or methanol in the ocean became too high and the temperature too low to allow convection in the overlying viscous ice shell (e.g. Mitri and Showman, 2008). Indeed, even for initially moderate ammonia fraction (about 2-3% relative to the total mass of H₂O), the concentration significantly increases as the ocean crystallizes and the ocean temperature gets below 250 K when the ammonia concentration in the crystallizing ocean reaches 10 wt% (Tobie et al., 2005b). At such a temperature, the viscosity at the base of the ice shell exceeds 5.10¹⁵ – 10¹⁶ Pa.s, which is too large to promote convection even for a 100 km-thick ice shell.

### 2.2.3 Geophysical evidence for an internal water ocean

For the Galilean moons, the most convincing evidence for subsurface oceans came from analysis of the magnetic data acquired by the Galileo mission during close flybys. The tilt of Jupiter's magnetic field relative to its spin pole means that the Galilean satellites experience a time-varying field and resulting induction currents. This effect has been used to probe the conductivity structure of the satellites, and deduce the presence of oceans for the three outer satellites (Khurana et al., 1998; Kivelson et al., 2002), and a partially-molten mantle for Io (Khurana et al., 2011). Unfortunately, Saturn's spin and magnetic poles are almost exactly coplanar, so the induction technique is very difficult to use at Titan (Saur et al., 2010). A recent study of Titan's magnetic field did not see any induction effects, and did not detect any permanent internal magnetic field (Wei et al. 2010). Although it is possible to have a molten core without generating a magnetic field (e.g. Venus, Io), this result is perhaps consistent with Titan's interior being cold and not fully differentiated.

On Titan, another technique based on the measurements of the electric field during the *Huygens* descent by the Permittivity, Wave and Altimetry (PWA) instrument (Béghin et al., 2007) provides however some information on the subsurface electrical conductivity. The measurements of low frequency waves and atmospheric conductivity by the PWA instrument revealed the existence of a Schumann-like resonance trapped within Titan's atmospheric cavity. The observed signal may be triggered and sustained by strong electric currents induced in the ionosphere by Saturn's magnetospheric plasma flow (Simões et al., 2007; Béghin et al., 2009, 2010). The characteristics of this trapped-resonant mode imply the presence of a conductive layer at about 30-60 km below the surface (Béghin et al., 2010), which may correspond to the ice/ocean interface provided that the ocean is sufficiently conductive, possibly doped with small amounts of ammonia and/or salts (Béghin et al., 2010).

As explained in Section 2.2.2, the observed topography may be explained by variations in ice shell thickness, resulting from inhomogeneous crystallization of an underlying ocean. If this analysis is correct, it implies that an internal ocean is still present in Titan's interior and that it is currently slowly crystallizing. This constitutes additional indirect evidence for an internal ocean inside Titan.

Another piece of evidence is provided by Titan's spin state. A non-spherical body, such as a synchronous satellite, precesses about its rotation axis at a rate determined by the ratio, $(C - \bar{A})/C$. Since the gravity coefficient $J_2 = (C - \bar{A})/M_T R_T^2$, if $J_2$ and the precession rate can be measured, the moment of inertia $C$ may be determined directly without assuming hydrostatic equilibrium. Although measuring a satellite's precession rate directly is difficult, measuring the obliquity (the angle between spin pole and
orbit pole) is relatively easy. If a satellite is dissipative, tidal torques drive the obliquity to a state in which the spin pole and orbit pole remain coplanar as the latter precesses around a fixed (invariable) pole (e.g. Peale, 1969). If this state, known as a damped Cassini state, exists, then measuring the obliquity is sufficient to determine the satellite moment of inertia (e.g. Bills and Nimmo, 2011). In practice, other small effects (such as Saturn's spin pole precession) have to be taken into account, but they do not alter the underlying physical picture.

Titan's spin pole location has been established through analysis of SAR images. It is very nearly coplanar with the invariable and orbit poles (suggesting occupation of a Cassini state), and corresponds to an obliquity of 0.32° (Stiles et al., 2008, 2010). If Titan occupies a damped Cassini state, the implied value of $C/(M_T R_T^2)$ is $0.45 \pm 0.02$ (Bills and Nimmo, 2011). This result is clearly unreasonable, in that it suggests a body which is under-dense in the interior, and also disagrees with the results derived from gravity. The solution is probably a combination of two factors. First, the fact that Titan's spin pole is not quite coplanar with the orbit and invariable poles suggests some degree of present-day excitation of the obliquity. Two possible candidate sources of excitation are the atmosphere (Tokano et al., 2011) or fluid motion in a subsurface ocean, both of which can transfer angular momentum to the ice shell. The second factor is that Titan is probably not acting as a rigid body (which the Cassini state analysis assumes), and that instead the ice shell is partially decoupled from the interior by an ocean. If the shell were totally decoupled, its normalized moment of inertia would be 0.67. The derived value of 0.45 suggests some degree of coupling, perhaps due to pressure or gravitational torques (Baland et al., 2011). Fortunately, continued monitoring of the position of Titan's spin pole should be able to disentangle these two effects, but the current results are strongly indicative of a subsurface ocean.

Another potential source of evidence for a subsurface ocean would be the existence of non-synchronous rotation (NSR) of the ice shell (Ojakangas and Stevenson, 1989). It was initially thought that Titan's shell was exhibiting such behaviour (Lorenz et al., 2008b), but this observation turned out to be erroneous and there is currently no evidence for NSR (Stiles et al., 2010).

Probably, the most promising technique to confirm the existence of a subsurface ocean would be to measure the tidal fluctuation from accurate gravity or topography measurements. The Radio Science Experiment on *Cassini* has in theory the capability to measure the tide-generated gravity signals (Rappaport et al., 1997, 2008). However, the atmospheric drag due to the extended atmosphere as well as the existence of relatively large signals in harmonics of degree three (Iess et al., 2010), suggesting significant signal at degrees 4 and 5 (Sotin et al., 2010), make the analysis more difficult than initially anticipated. In spite of these difficulties, a recent analysis of the *Cassini* gravity data ob-

tained during six close flybys indicates that Titan responds to the variable tidal field exerted by Saturn at a detectable level, providing further evidence for a global ocean at depth on Titan (Iess et al., 2012). A future mission with dedicated multiple flybys or even better in orbit around Titan will be required to confirm this detection and to get more precise informations on the ocean and ice shell thicknesses.

### 2.2.4 Summary

At the start of this Section we posed three questions; we will conclude by providing some provisional answers.

1. **To what extent is Titan differentiated?** Based on the gravity results (Section 2.2.1), the normalized moment of inertia of Titan is estimated between 0.33 and 0.34. This clearly indicates that Titan is not fully differentiated and does not have an iron core like Ganymede. Earlier we suggested two possible internal structures – one in which rock and ice are incompletely separated, the other involving a rocky core that is partially hydrated – that satisfy the rather high moment of inertia. We discuss this further in Sections 2.3.4 and 2.4.1.1.

2. **Is the ice shell convecting?**

A conductive ice shell is favoured by the model for Titan's long-wavelength shape (Section 2.2.2) and is consistent with the requirement that its eccentricity not damp too rapidly. To avoid convection, the ice shell must have a high viscosity, which can be achieved if the ocean beneath is ammonia-rich and therefore cold. Convection may have occurred in the past before the ammonia concentration increased, but this convective period must have lasted less than 0.5-1 billion years in order to explain Titan's present-day elevated orbital eccentricity.

3. **Does Titan possess a sub-surface ocean?** Electric signals measured in Titan's atmosphere by *Huygens* suggest the presence of a conductive layer at depth, possibly a water ocean doped with small amounts of salts or ammonia to increase the electrical conductivity. Titan's spin state and tidal gravity response are also suggestive of a decoupled shell (Section 2.2.3) while its long-wavelength shape can be explained if a subsurface ocean is assumed (Section 2.2.2). Given the existence of subsurface oceans on the comparably-sized Galilean satellites, Callisto and Ganymede, it would be surprising if Titan did not possess an ocean. However, one important difference is that Titan's ocean may be colder and more ammonia-rich (see below).
## 2.3 ORIGIN OF THE SATURNIAN SYSTEM AND ACCRETION OF TITAN

Titan was formed as a regular satellite in a disk that was the outgrowth of the formation of Saturn itself. The conditions under which Saturn formed therefore had a direct impact on the formation of Titan. They determined its composition, the architecture of the whole system as well as the final assemblage that led to the formation of Titan and of the other regular moons. The present section is structured around three main questions which are of particular relevance to Titan's formation.

1. What was the composition of the primordial building blocks that formed Titan and the Saturnian system? Titan and the other Saturnian moons formed from a mixture of rock and ice initially present in the vicinity of Saturn. We will present in this section different theoretical and observational constraints that can be used to estimate the composition of these ice-rock blocks, notably their volatile content.

2. What makes the Jupiter and Saturn systems so different? Contrary to the Jovian system, there is no gradual variation of density and hence rock/ice ratio in the Saturnian moons as a function of planet distance. Titan is the unique big moon, the second largest moon Rhea having a mass less than 2% of Titan's mass. Another particular characteristic of Saturn is its extended ring system. In this Section, we will describe the main phases of satellite formation in a disk around giant planets and try to identify which processes differed between Jupiter and Saturn.

3. Can Titan accrete undifferentiated? The Cassini gravity data suggests that Titan is moderately differentiated at present (section 2.2.1), suggesting that the differentiation process has been limited. One possibility is that Titan accreted relatively cold, which would have delayed the differentiation. To address that question, we will review the physics of accretion and determine under which conditions Titan might have remained cold during accretion.

### 2.3.1 Formation of ices in Saturn's environment

The primitive building blocks that formed Titan and the rest of the satellite system were presumably composed of a complex assemblage of various ices, hydrates, silicate minerals and organics. An important issue associated with the formation of giant planets and their satellite systems concerns the condensation of water ice and other ice compounds as the solar nebula cools down. The condensation sequence depends on the initial gas phase, the cooling rate of the disk and the stability of the gas molecules in solid phases.

The composition of the initial gas phase of the disk can be defined by assuming that the abundances of all elements, including oxygen, are protosolar (Asplund et al., 2009) and

that O, C, and N exist only in the form of H₂O, CO, CO₂, CH₃OH, CH₄, N₂, and NH₃. The abundances of CO, CO₂, CH₃OH, CH₄, N₂ and NH₃ are then determined from the adopted CO/CO₂/CH₃OH/CH₄ and N₂/NH₃ gas phase molecular ratios. Once the abundances of these molecules are fixed, the remaining O gives the abundance of H₂O. Concerning the distribution of elements in the main volatile molecules, the ratio of CO/CO₂/CH₃OH/CH₄ can be set to 70/10/2/1 in the gas phase of the disk, values that are consistent with the ISM measurements considering the contributions of both gas and solid phases in the lines of sight (Mousis et al., 2009c). In addition, S is assumed to exist in the form of H₂S, with H₂S/H₂ = 0.5 × (S/H₂)⊙, and other refractory sulfide components (Pasek et al., 2005). N₂/NH₃ is assumed equal to 10/1 in the nebula gas-phase, a value compatible with thermochemical models of the solar nebula (Lewis and Prinn, 1980) and with observations of cometary comae (e.g. Hersant et al., 2008).

To first order, the process by which volatiles are trapped in icy planetesimals can be calculated assuming thermo-dynamic equilibrium. The sequence of clathration<fnref n="1" /> and condensation can be determined by using the equilibrium curves of stoichiometric hydrates, clathrates and pure condensates, and the thermodynamic path (hereafter cooling curve) detailing the evolution of temperature and pressure between roughly 5 and 20 AU, a distance range largely encompassing the migration path followed by proto-Jupiter and proto-Saturn during their formation in the solar nebula (Alibert et al., 2005b). We refer the reader to the works of Papaloizou and Terquem (1999) and Alibert et al. (2005a) for a full description of the turbulent model of accretion disk used here. For the pressure range expected in the solar nebula, most of the gas compounds, except CO₂, are enclathrated at temperatures higher than the temperatures they condense in the form of pure condensates. This has fundamental consequences for the trapping of volatiles in primordial building blocks.

As illustrated in Figure 2.1, the cooling curve intercepts the equilibrium curves of the different ices at particular temperatures and pressures. For each ice considered, the domain of stability is the region located below its corresponding equilibrium curve. The clathration process stops when no more crystalline water ice is available to trap the volatile species. As CO₂ is the only species that crystallizes at a higher temperature than its associated clathrate, solid CO₂ is assumed to be the only existing condensed form of CO₂ in this environment. In addition, CH₃OH is assumed to be in the form of pure ice, since, to the best of our knowledge, no experimental data concerning the equilibrium curve of its associated clathrate have been reported in the literature.

<footnote n="1">¹ Clathrate hydrates are crystalline water-based structures that form cages in which gas molecules can be trapped (Lunine and Stevenson, 1985; Choukroun et al., 2011), and for this reason they play a crucial role in the volatile trapping. They also play a key role in the evolution of the volatile inventory on Titan as it will be further discussed in Sections 2.4.1.3 and 2.4.1.4.</footnote>
<figure><img src="image_360.png" /></figure>

**Figure 2.1** Equilibrium curves of NH₃-H₂O hydrate, H₂S, Xe, CH₄ and CO clathrates (solid lines), CH₃OH, CO₂, Kr, CO, Ar and N₂ pure condensates (dotted lines), and thermodynamic path followed by the primordial nebula between 5 and 20 AU as a function of time, respectively, assuming a full efficiency of clathration. Species remain in the gas phase above the equilibrium curves. Below, they are trapped as clathrates or simply condense. The equilibrium curves of hydrates and clathrates derive from Lunine and Stevenson (1985)'s compilation of published experimental work, in which data are available at relatively low temperatures and pressures. The equilibrium curves of pure condensates used here derive from the compilation of laboratory data provided by Lide (2004).

For simplicity, the clathration efficiency is assumed to be total, implying that guest molecules had the time to diffuse through porous water-ice solids before their growth into planetesimals and their accretion by proto-planets or proto-satellites. This statement remains plausible only if collisions between planetesimals have exposed essentially all the ice to the gas over time scales shorter or equal to planetesimal lifetimes in the nebula (Lunine and Stevenson, 1985). In the present case, which is based on the protosolar abundances of Asplund et al. (2009), NH₃, H₂S, Xe, CH₄ and ~31.5% of CO form NH₃-H₂O hydrate and H₂S, Xe, CH₄ and CO clathrates with the available water in the outer nebula. The remaining CO, as well as N₂, Kr, Ar and Ne (not shown on the figure), whose clathration normally occurs at lower temperatures, remain in the gas phase until the nebula cools enough to allow the formation of pure condensates. As the gas phase composition of the disk is assumed not to vary with heliocentric distance, the calculated clathration conditions remain the same in the 5–20 AU range of the nebula (see Marboeuf et al. (2008) for details). Once crystallized, these ices will agglomerate and form planetesimals large enough to decouple from the nebular gas and will be accreted by the forming giant planets and satellites.

Many works published in the last decade and detailing the formation conditions of ices in the solar nebula (Gautier et al., 2001; Hersant et al., 2008; Mousis and Alibert, 2006; Alibert and Mousis, 2007, e.g.) have used the equilibrium

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

**Figure 2.2** Mole fraction of volatiles engaged in H₂S (a), Xe (b), CH₄ (c) and CO (d) clathrates. Grey and dark bars correspond to structure I and structure II clathrates, respectively, which are two possible clathrate structures with slightly different cage sizes (e.g. Sloan, 1998).

curves established experimentally for single guest molecules (similar to the curves shown on Fig. 2.2). This simplified approach permits an understanding to first order of the condensation sequence in the nebula, but in reality, clathrates forming in a gas mixture always simultaneously trap several guest molecules (e.g. Lunine and Stevenson, 1985; Sloan, 1998). Here we consider the effect of simultaneous trapping of guest species on the predicted composition of volatiles trapped in the water ice.

The composition of the clathrate phases can be constrained by using a statistical thermodynamic model derived from the approach of van der Waals and Platteeuw (1959) (see Lunine and Stevenson, 1985; Mousis et al., 2010, for details). This thermodynamic approach is based on the use of intermolecular potentials, which themselves depend on parameters describing the interaction between the guest molecule and the cage made of water molecules, called “Kihara parameters”. The afore-mentioned formation sequence of the different ices occurring during the cooling of the solar nebula is considered (see Fig. 2.1). As a consequence, any volatile already trapped or condensed at a higher temperature than the formation temperature of the clathrate under consideration is excluded from the coexisting gas phase composition. This implies that CO₂, Xe, CH₄, CO, Kr, Ar and N₂ are considered as possible guests in the case of H₂S clathrate. On the other hand, only N₂, Ar, Kr can become guests in CO clathrate. Figure 2.2 represents the mole fraction *f* of volatiles engaged in structure I and structure II clathrates *a priori* dominated by H₂S, Xe, CH₄ and CO and formed in the primordial nebula. This figure shows that H₂S, Xe, CH₄ and CO are the dominating guest species if structure I clathrates form in the nebula. However, in the case of structure II clathrate formation, CO becomes a minor compound in the clathrate that is expected to be dominated by this molecule. Indeed, Kr becomes more abundant than CO in this clathrate, irrespective of their initial abundances in the gas phase