Magnetite-apatite (MtAp) deposits have attracted considerable attention due to their complex genesis and economic importance. These deposits are rich in magnetite ore and can bear significant rare earth elements, but their exact formation mechanisms remain controversial. This study aims to understand the formation processes of MtAp deposits by investigating the role of iron-rich magmatic liquids. Focusing on the El Laco deposit, northern Chile, we follow the hypothesis that iron-rich liquids separate from silicate magma through liquid immiscibility. Building on previous research, this study employs a three-phase 1-D mechanical model to simulate the separation and accumulation of immiscible iron-rich melts within increasingly crystalline magma. The model reproduces the previously suggested transition from isolated droplet settling to an interconnected drainage network and quantifies the relative efficiency of both modes of phase separation. Using scaling analysis, we define porous, mush and suspension flow regimes and construct a regime diagram for three-phase flow. The results indicate that the separation efficiency of immiscible iron-rich melts is maximized in the mush flow regime at intermediate crystallinity. The model-derived accumulation rate of iron-rich melts can be used to estimate the time required to form magnetite deposits of a given scale. Our findings support the physical viability of the liquid immiscibility hypothesis for the genesis of MtAp deposits, offering new insights into the mechanical efficiency of melt separation and contributing to a broader understanding of the formation mechanisms of other valuable deposits that have been linked to immiscible melts.
Magma bodies play a critical role in Earth's geological evolution, influencing volcanic activity, crustal differentiation, and planetary-scale processes. Understanding their thermo-chemical and mechanical evolution requires models that integrate fluid dynamics, phase changes, and chemical transport. This study presents a new numerical model that couples these processes using a multi-phase, multi-component formulation. The model simulates convection, phase segregation, and thermo-chemical evolution across a wide range of scales, from crustal magma chambers to planetary magma oceans. To ensure numerical stability and physical realism, adaptive regularisation schemes are implemented, including eddy diffusivity for higher-dimensional turbulent flows and convective mixing diffusivity for one-dimensional column models. Benchmark tests confirm the accuracy of the numerical scheme, and use cases demonstrate its applicability to scenarios such as fractional crystallisation, wall-rock assimilation, and magma recharge on crustal scales, and magma ocean solidification on planetary scales. By providing an open-source implementation, this work aims to advance our understanding of dynamic magmatic systems and their role in planetary evolution.
The Duluth Complex is a large mafic intrusive system located in northeastern Minnesota emplaced as part of the 1.1-Ga Midcontinent Rift. Several Fe-Ti oxide-bearing ultramafic intrusions are hosted along the Western Margin of the Duluth Complex, and are discordant bodies present in a variety of geometries, hosted in multiple rock types, and dominated by peridotite, pyroxenite, and semi-massive to massive Fe-Ti oxide rock types. Their origin has been debated, and here we present geochemical evidence and modeling that supports a purely magmatic origin for the Titac and Longnose Fe-Ti oxide-bearing ultramafic intrusions. Ilmenite and titanomagnetite textures indicate a protracted cooling process, and delta S-34 values of sulfides reveal little assimilation of the footwall Virginia Formation, a fine-grained pelitic unit that contains sulfide-rich bands. We model the crystallization of a hypothetical parental magma composition to the host intrusion of Longnose using Rhyolite-MELTS and demonstrate that the accumulation of Fe-Ti oxides in the discordant intrusions cannot be explained by density-driven segregation of crystallized Fe-Ti oxides. Instead, we show that the development of silicate liquid immiscibility, occurring by the unmixing of the silicate melt into conjugate Si- and Fe-rich melts, can result in the effective segregation and transportation of the Fe-rich melt. The Fe-rich melt is similar to 2 orders of magnitude less viscous than the Si-rich melt, allowing the Fe-rich melt to be more effectively segregated and transported in the mush regime (crystallinities >50%). This suggests that viscosity, in addition to density, plays a significant role in forming the discordant Fe-Ti oxide-bearing ultramafic intrusions. We propose a genetic model that could also be responsible for the Fe-Ti oxide-rich layers or bands that are hosted within the igneous stratigraphy of mafic intrusions of the Duluth Complex.
Magmatic systems in the Earth's mantle and crust can range from melt-poor partially molten rock to trans-crustal magma mushes with ephemeral lenses of melt-rich suspensions. Most process-based models of magmatic systems, however, are limited to two-phase porous flow at low melt fractions (<20%) or suspension flow at high melt fractions (>60%). A lack of formal extensions to intermediate phase fractions has long hindered investigations into the dynamics of mush flows. To address this knowledge gap and unify two-phase magma flow models, we present a two-dimensional system-scale numerical model of the fluid mechanics of an n-phase system valid at all phase fractions. The numerical implementation uses a finite-difference staggered-grid approach with a dampened pseudo-transient iterative algorithm and is verified using the Method of Manufactured Solutions. Numerical experiments replicate known limits of two-phase flow including rank-ordered porosity wave trains in 1D and porosity wave breakup in 2D in the porous flow regime, as well as particle concentration waves in 1D and mixture convection in 2D in the suspension flow regime. In the mush regime, numerical experiments show strong liquid localisation into pockets and stress-aligned bands. A tentative application to a three-phase, solid-liquid-vapour system demonstrates the broad utility of the n-phase general model and its numerical implementation. The model code is available open source at github.com/kellertobs/pantarhei.
<p>Collisions of planetesimals with pre-differentiated cores aided in the formation of the cores in the terrestrial planets in our Solar System. The differentiation of the interior of planetesimals from a primitive chondritic composition to differentiated metallic core and silicate mantle is therefore an important step in the development of the early Solar System, setting the initial conditions for collisional planetary growth. However, meteoritic evidence of the differentiation stage of planetesimals are rare and challenging to interpret. Therefore, mathematical models are necessary to constrain the timescales and elucidate the physics behind metal-silicate segregation in planetesimals.&#160;We present progress towards a new numerical model which models the thermo-chemical and fluid-mechanical evolution of silicate-metal planetesimals. The model is comprised of four material phases: solid and liquid silicates, and solid and liquid Fe-FeS metal. Hence, the model quantifies the different segregation rates of metal alloys at various stages of planetary melting. The thermo-chemical evolution quantifies the melting and chemical fractionation of the materials using two separate phase diagrams for the silicate and iron systems respectively and calibrated based on enstatite chondrite meteorites as the starting primitive material. The fluid mechanics represents conditions where the silicate liquid phase dominates and captures the settling of solid particles and immiscible liquid metal droplets by hindered Stokes settling. Our model will simulate core formation under a range of initial planetesimal compositions, nebular compositions, and planetesimal sizes. We also intend to implement compatible element partitioning into our model, in particular Hf-W isotope partitioning, to compare our model results with the meteoritic record. Results from our model will allow more robust estimations of the timescales for planetesimal core formation.</p> <div> <div> <div>&#160;</div> </div> <div> <div>&#160;</div> </div> </div>
Despite the first-order importance of crystallisation–differentiation for arc magma evolution, several other processes contribute to their compositional diversity. Among them is the remelting of partly crystallised magmas, also known as cumulate melting or ‘petrological cannibalism’. The impact of this process on the plutonic record is poorly constrained. We investigate a nepheline-normative dyke suite close to the Blumone gabbros, a large amphibole-gabbro unit of the Tertiary Southern Alpine Adamello igneous complex. The compositions of the studied dykes are characterised by low SiO 2 (43–46 wt. %), MgO (5.0–7.2 wt. %), Ni (18–40 μg/g), and high Al 2 O 3 (20.2–22.0 wt. %) contents. Phenocrystic plagioclase in these dykes exhibits major, trace, and Sr isotope compositions similar to Blumone cumulate plagioclase, suggesting a genetic link between the nepheline-normative dykes and the amphibole-gabbro cumulates. We tested this hypothesis by performing saturation experiments on a nepheline-normative dyke composition in an externally heated pressure vessel at 200 MPa between 975 and 1100 °C at fO 2 conditions close to the Ni–NiO buffer. Plagioclase and spinel are near-liquidus phases at and above 1050 °C, contrasting with the typical near-liquidus olivine ± spinel assemblage in hydrous calc-alkaline basalts. The alkaline nature of the dykes results from the abundance of amphibole in the protolith, consistent with melting of amphibole-gabbro cumulates. We modelled the heat budget from the repeated injection of basaltic andesite into a partly crystallised amphibole-gabbro cumulate. The results of this model show that no more than 7% of the cumulate pile reaches temperatures high enough to produce nepheline-normative melts. We propose that such nepheline-normative dykes are a hallmark of hydrous cumulate melting in subvolcanic plumbing systems. Therefore, ne-normative dykes in arc batholiths may indicate periods with high magma fluxes.
Understanding the dynamics of magma ocean crystallisation during planetary cooling can elucidate the initial mantle structure and subsequent evolution of early planetary bodies. However, most studies on magma ocean crystallisation focus on either the thermo-chemistry (e.g., Johnson et al. 2021) or the fluid dynamics of a cooling magma ocean (e.g., Maurice et al. 2017). This precludes investigations into coupled thermo-mechanical processes, such as the effect of convection and phase segregation on chemical differentiation. However, coupled models are challenging to implement due to their numerical complexity and limited experimental constraints on magma ocean crystallisation for model calibration.We develop a two-phase, 6-component model in a 2D rectangular domain based on a multi-phase, multi-component reactive transport model framework (Keller & Suckale, 2019). Magma ocean convection is modelled using Stokes equations while crystal settling is calculated using a form of hindered Stokes law. The fluid mechanics model is coupled with a thermo-chemical model of evolving temperature, phase proportions, and phase compositions to form a reactive transport model, following Keller & Katz (2016). We apply this model to the lunar magma ocean (LMO) by describing the melt and crystal compositions with 6 pseudo-components (approximating forsterite-fayalite, orthopyroxene-clinopyroxene and anorthite-albite mineral systems). To calibrate the melting temperature and composition of each component, we fit data from fractional crystallisation experiments for a Taylor Whole Moon composition (Schmidt & Krättli 2022) using a transitional Markov Chain Monte Carlo method.The 6-component melting model calibrated to experimental data is successfully implemented in the reactive transport model. First results indicate the importance of crystal settling speed and magma convection speed on convective mixing, magma ocean stratification, and crystal cumulate formation. The small size of the Moon and its relatively well-constrained magma ocean history, make the LMO an excellent case study to apply the model. However, with the aid of new experimental data for larger and chemically different planets, such as Mars, this model can provide more general insight into the early evolution of terrestrial bodies.REFERENCES: Maurice et al. (2017) doi:10.1002/2016JE005250, Johnson et al. (2021) doi: 10.1016/j.epsl.2020.116721, Keller & Suckale (2019) doi:10.1093/gji/ggz287, Keller & Katz (2016) doi: 10.1093/petrology/egw030, Schmidt & Krättli (2022) doi:10.1029/2022JE007187
Magmatic systems in the Earth's mantle and crust contain multiple phases including solid crystals, liquid melt and low viscosity fluids. Depending on depth, tectonic setting and chemical composition, magmatic systems can range from partially molten rock at low melt fraction to magma mushes at intermediate melt fraction to magmatic suspensions at high melt fraction. However, the theories underpinning most process-based models of magmatic systems describe magma as a single-phase fluid, or as a two-phase mixture either in the porous flow regime at low melt fractions or the suspension flow regime at high melt fractions. Connections between the two-phase endmember theories are poorly established and hinder investigations into the dynamics of mush flows at intermediate phase fractions, leaving a significant gap in bridging trans-crustal magma processing from source to surface. To address this knowledge gap and unify two-phase magma flow models, we develop a two-dimensional system-scale numerical model of the fluid mechanics of an n-phase system at all phase proportions, based on a recent theoretical model for multi-phase flows in igneous systems. We apply the model to two-phase, solid-liquid mixtures by calibrating transport coefficients to theory and experiments on mixtures with olivine-rich rock and basaltic melt using a Bayesian parameter estimation approach. We verify the model using the Method of Manufactured Solutions and test the scalability for high resolution modelling. We then demonstrate 1D and 2D numerical experiments across the porous, mush and suspension flow regimes. The experiments replicate known phenomena from endmember regimes, including rank-ordered porosity wave trains in 1D and porosity wave breakup in 2D in the porous flow regime, as well as particle concentration waves in 1D and mixture convection in 2D in the suspension flow regime. By extending self-consistently into the mush regime, the numerical experiments show that the weakening solid matrix facilitates liquid localisation into liquid-rich shear bands with their orientation controlled by the solid stress distribution. Although the present model can already be used to investigate three-phase mixtures using conceptually-derived transport coefficients, more rigorous calibration to experiments and endmember theories is needed to ensure accurate time scales and mechanics. With a self-consistent way to examine multi-phase mixtures at any phase proportion, this new model transcends theoretical limitations of existing multi-phase numerical models to enable new investigations into two-phase or higher magma mush dynamics.
<p>The discovery of <sup>18</sup>O-depleted igneous rocks at Krafla, Iceland, suggests that the system interacted with crustal rocks that experienced high-temperature hydrothermal alteration by a meteoric fluid to deviate from the expected mantle signature (&#948;<sup>18</sup>O = 5.5 &#8240;). Such assimilation is documented in low-&#948;<sup>18</sup>O settings worldwide, however, the mechanisms of this dynamic process remain poorly understood. &#160;Due to intense drilling activity and exploration at Krafla, both hydrothermally altered crustal rocks and parental magma are comparably well characterized, making Krafla a great case study for the application of a numerical model that can further advance the understanding of the formation process of low-&#948;<sup>18</sup>O magmas. In this study, we use a new three-phase two-component thermo-chemical-mechanical model to simulate the effect of variable crustal compositions on the assimilation process and the magma chamber dynamics. &#160;We define the simplified square-shaped magma chamber (10 x 10 m) of magma with initially basaltic composition (1250 &#176;C) that assimilates the crustal rock (500 &#176;C) at the top and bottom. Our results indicate that convective behaviour and the formation of cumulate layers can significantly hinder the assimilation process. While the crystal settling Stokes speed scale is the dominant driver for the formation of this boundary layer, depending on the assimilation timescales, the mushy chamber margins are able to grow to sufficient thickness to prohibit additional assimilation of low-&#948;<sup>18</sup>O crustal material. Density and buoyancy contrasts produce three types of convection: chamber convection, layered convection and plume driven convection. Final magma compositions in our preliminary model outputs range from mafic to intermediate but are not able to reach the felsic compositions encountered at Krafla. This suggests that evolution towards the erupted low-&#948;<sup>18</sup>O rhyolitic products involved multiple stages or included additional factors not yet accounted for in our model. Further refining of this and similar thermo-chemical-mechanical model setups may provide important new insights into the assimilation dynamics in the Krafla volcanic field and other low-&#948;<sup>18</sup>O settings worldwide.</p><p>&#160;</p>
Melt extraction from the partially molten mantle is among the fundamental processes shaping the solid Earth today and over geological time. A diversity of properties and mechanisms contribute to the physics of melt extraction. We review progress of the past ∼25 years of research in this area, with a focus on understanding the speed and style of buoyancy-driven melt extraction. Observations of U-series disequilibria in young lavas and the surge of deglacial volcanism in Iceland suggest this speed is rapid compared to that predicted by the null hypothesis of diffuse porous flow. The discrepancy indicates that the style of extraction is channelized. We discuss how channelization is sensitive to mechanical and thermochemical properties and feedbacks, and to asthenospheric heterogeneity. We review the grain-scale physics that underpins these properties and hence determines the physical behavior at much larger scales. We then discuss how the speed of melt extraction is crucial to predicting the magmatic response to glacial and sea-level variations. Finally, we assess the frontier of current research and identify areas where significant advances are expected over the next 25 years. In particular, we highlight the coupling of melt extraction with more realistic models of mantle thermochemistry and rheological properties. This coupling will be crucial in understanding complex settings such as subduction zones. ▪ Mantle melt extraction shapes Earth today and over geological time. ▪ Observations, lab experiments, and theory indicate that melt ascends through the mantle at speeds ∼30 m/year by reactively channelized porous flow. ▪ Variations in sea level and glacial ice loading can cause significant changes in melt supply to submarine and subaerial volcanoes. ▪ Fluid-driven fracture is important in the lithosphere and, perhaps, in the mantle wedge of subduction zones, but remains a challenge to model.
Magnetite-apatite deposits are important sources of iron and other metals. A prominent exam- ple are the magnetite lavas at the El Laco volcano, Northern Chile. Their formation processes remain debated. Here, we test the genetic hypothesis that an Fe-rich melt separated from silicate magma and ascended along collapse-related fractures. We complement recent analy- ses with thermodynamic modelling to corroborate Fe-Si liquid immiscibility evident in melt inclusions at El Laco and present viscometry of Fe- and Si-rich melts to assess the time and length scales of immiscible liquid separation. Using a rock deformation model, we demonstrate that volcano collapse can form failure zones extending towards the edifice flanks along which the ore liquid ascends towards extrusion driven by vapour exsolution despite its high density. Our results support the proposed magmatic genesis for the El Laco deposits. Geochemical and textural similarities indicate magnetite-apatite deposits elsewhere form by similar processes.
The Martian nakhlite meteorites, which represent multiple events that belong to a single magma source region represent a key opportunity to study the evolution of Martian petrogenesis. Here 16 of the 26 identified nakhlite specimens are studied using coupled electron backscatter diffraction (EBSD) and emplacement end‐member calculations. EBSD was used to determine shape preferred orientation of contained augite (high Ca‐clinopyroxene) phenocrysts by considering their crystallographic preferred orientation (CPO). Parameters derived from EBSD, and energy dispersive X‐ray spectroscopy spectra were used in basic emplacement models to assess their dominant mechanism against three end‐member scenarios: thermal diffusion, crystal settling, and crystal convection. Results from CPO analyses indicate low intensity weak‐moderate CPO. In all samples, a consistent foliation within the <001> axes of augite are observed typically coupled with a weaker lineation CPO in one of the other crystallographic axes. These CPO results agree best with crystal settling being the dominant emplacement mechanism for the nakhlites. Modeled crystal settling results identify two distinguishable groups outside of the model's resolution indicating the presence of secondary emplacement mechanisms. Comparison of the two identified groups against CPO, geochemical, and age parameters indicate random variability between individual meteorites. Therefore, coupled CPO and emplacement modeling results identify an overarching characteristic of a dominant crystal settling emplacement mechanism for the nakhlite source volcano despite exhibiting random variation with each discharge through time.
Magnetite-apatite deposits are important sources of iron and other metals. A prominent example are the magnetite lavas at the El Laco volcano, Northern Chile. Their formation processes remain debated. Here, we test the genetic hypothesis that an Fe-rich melt separated from silicate magma and ascended along collapse-related fractures. We complement recent analyses with thermodynamic modelling to corroborate Fe-Si liquid immiscibility evident in melt inclusions at El Laco and present viscometry of Fe- and Si-rich melts to assess the time and length scales of immiscible liquid separation. Using a rock deformation model, we demonstrate that volcano collapse can form failure zones extending towards the edifice flanks along which the ore liquid ascends towards extrusion driven by vapour exsolution despite its high density. Our results support the proposed magmatic genesis for the El Laco deposits. Geochemical and textural similarities indicate magnetite-apatite deposits elsewhere form by similar processes.
Crystals retain an imprint of the dynamic changes within a magma reservoir and hence contain invaluable information about the evolving conditions inside volcanic plumbing systems. However, instead of telling a single, simple story, they comprise overprinted evidence of numerous processes relating to temperature, pressure and composition that drive crystal precipitation and dissolution in magmatic systems. To decipher these different elements in the story that crystals tell, we attempt to identify the observational signatures of a simple, yet ubiquitous process: crystal precipitation and dissolution during magma cooling. To isolate this process in a complex magmatic system with intricate dynamic feedbacks, we assume that synthetic crystals precipitate and dissolve rapidly in response to deviations from thermodynamic equilibrium. In our crystalline‐scale simulations, synthetic crystals drag along the cooler‐than‐ambient melt in which they precipitated and can drive a temperature‐dependent, crystal‐driven convection. We analyze the non‐dimensional conditions for this coupled convection and record the heterogeneous thermal histories that synthetic crystals in this flow regime experience. We show that many synthetic crystals dissolve, loosing their thermal record of the convection. Based on our findings, we suggest that heterogeneity in the thermal history of crystals is more indicative of local, crystal‐scale processes than the overall, system‐wide cooling trend.
Concealed deep beneath the oceans is a carbon conveyor belt, propelled by plate tectonics. Our understanding of its modern functioning is underpinned by direct observations, but its variability through time has been poorly quantified. Here we reconstruct oceanic plate carbon reservoirs and track the fate of subducted carbon using thermodynamic modelling. In the Mesozoic era, 250 to 66 million years ago, plate tectonic processes had a pivotal role in driving climate change. Triassic-Jurassic period cooling correlates with a reduction in solid Earth outgassing, whereas Cretaceous period greenhouse conditions can be linked to a doubling in outgassing, driven by high-speed plate tectonics. The associated 'carbon subduction superflux' into the subcontinental mantle may have sparked North American diamond formation. In the Cenozoic era, continental collisions slowed seafloor spreading, reducing tectonically driven outgassing, while deep-sea carbonate sediments emerged as the Earth's largest carbon sink. Subduction and devolatilization of this reservoir beneath volcanic arcs led to a Cenozoic increase in carbon outgassing, surpassing mid-ocean ridges as the dominant source of carbon emissions 20 million years ago. An increase in solid Earth carbon emissions during Cenozoic cooling requires an increase in continental silicate weathering flux to draw down atmospheric carbon dioxide, challenging previous views and providing boundary conditions for future carbon cycle models.
Fractional crystallization is an essential process proposed to explain worldwide compositional abundances of igneous rocks. It requires crystals to precipitate from the melt and segregate from its residual melt, or experience crystal fractionation. The compositional abundances of volcanic systems show a bell curve distribution suggesting that the process has variable efficiencies. We test crystal fractionation efficiency in convective flow in low to intermediate crystallinity regime. We simulate the physical segregation of crystals from their residual melt at the scale of individual crystals, using a direct numerical method. We find that at low particle Reynolds numbers, crystals sink in clusters. The relatively rapid motion of clusters strips away residual melt. Our results show cluster settling can imprint observational signatures at the crystalline scale. The collective crystal behavior results in a crystal convection that governs the efficiency of crystal fractionation, providing a possible explanation for the bell curve distribution in volcanic systems.
Magma matters. From magmatism facilitating the differentiation of terrestrial planets into core, mantle and crust, to the magmatic activity that modulates plate tectonics and deep volatile cycles to maintain a habitable Earth, to volcanism that causes terrible hazards but also provides rich energy and mineral resources – igneous processes are integral to Earth and other planets. Our understanding of volcanoes and their deep magmatic roots derives from a range of disciplines including field geology, petrology and geochemistry, and geophysical imaging. Observational and experimental studies, however, are hampered by incomplete access to processes that play out across scales ranging from sub-micron size to thousands of kilometres, and from seconds to billions of years. Computational modelling provides tools for investigating igneous processes across these scales. Over the past decade, my research has been focused on advancing the theoretical description and numerical application of multi-phase reaction–transport processes at the volcano to planetary scale. Mixture theory provides a framework to represent the spatially averaged behaviour of a large sample of microscopic phase constituents such as mineral grains, melt films, and vapour bubbles. This approach has been used successfully to model both porous flow of melt percolating through compacting partially molten rock, as well as suspension flow of crystals settling in convecting magma bodies. My recent work has introduced a new modelling framework to bridge the porous and suspension flow limits, and to extend beyound solid-liquid systems to multi-phase systems including several solid, liquid, and vapour phases. These advances provide new insights into the dynamics of crustal mush bodies, the outgassing and eruption of shallow magma reservoirs, and the generation of mineral resources by exsolution of exotic magmatic liquids.