Skip to main content

Main menu

  • Home
  • Articles
    • Current
    • Special Feature Articles - Most Recent
    • Special Features
    • Colloquia
    • Collected Articles
    • PNAS Classics
    • List of Issues
  • Front Matter
    • Front Matter Portal
    • Journal Club
  • News
    • For the Press
    • This Week In PNAS
    • PNAS in the News
  • Podcasts
  • Authors
    • Information for Authors
    • Editorial and Journal Policies
    • Submission Procedures
    • Fees and Licenses
  • Submit
  • Submit
  • About
    • Editorial Board
    • PNAS Staff
    • FAQ
    • Accessibility Statement
    • Rights and Permissions
    • Site Map
  • Contact
  • Journal Club
  • Subscribe
    • Subscription Rates
    • Subscriptions FAQ
    • Open Access
    • Recommend PNAS to Your Librarian

User menu

  • Log in
  • My Cart

Search

  • Advanced search
Home
Home
  • Log in
  • My Cart

Advanced Search

  • Home
  • Articles
    • Current
    • Special Feature Articles - Most Recent
    • Special Features
    • Colloquia
    • Collected Articles
    • PNAS Classics
    • List of Issues
  • Front Matter
    • Front Matter Portal
    • Journal Club
  • News
    • For the Press
    • This Week In PNAS
    • PNAS in the News
  • Podcasts
  • Authors
    • Information for Authors
    • Editorial and Journal Policies
    • Submission Procedures
    • Fees and Licenses
  • Submit
Research Article

Critical vaporization of MgSiO3

Bing Xiao and Lars Stixrude
  1. aDepartment of Earth Sciences, University College London, London WC1E 6BT, United Kingdom

See allHide authors and affiliations

PNAS May 22, 2018 115 (21) 5371-5376; first published May 7, 2018; https://doi.org/10.1073/pnas.1719134115
Bing Xiao
aDepartment of Earth Sciences, University College London, London WC1E 6BT, United Kingdom
  • Find this author on Google Scholar
  • Find this author on PubMed
  • Search for this author on this site
Lars Stixrude
aDepartment of Earth Sciences, University College London, London WC1E 6BT, United Kingdom
  • Find this author on Google Scholar
  • Find this author on PubMed
  • Search for this author on this site
  • For correspondence: l.stixrude@ucl.ac.uk
  1. Edited by Henry J. Melosh, Purdue University, West Lafayette, IN, and approved March 27, 2018 (received for review November 2, 2017)

  • Article
  • Figures & SI
  • Info & Metrics
  • PDF
Loading

Significance

During accretion, giant impacts may partially vaporize planets, with profound consequences for energy and mass transfer between the growing planet and its satellites. We find that giant impacts produce much more supercritical fluid than previously thought, with a much lower liquid–vapor critical temperature. The vapor is depleted in Mg with respect to coexisting liquid and has very different chemical bonding, being dominated by covalent molecules (O2 and SiO) up to the critical point. These molecules also appear in significant fractions within the near-critical liquid phase, making the liquid structure fundamentally different from that near the melting point.

Abstract

Inhomogeneous ab initio molecular dynamics simulations show that vaporization of MgSiO3 is incongruent and that the vapor phase is dominated by SiO and O2 molecules. The vapor is strongly depleted in Mg at low temperature and approaches the composition of the liquid near the critical point. We find that the liquid–vapor critical temperature (6,600±150 K) is much lower than assumed in hydrodynamic simulations, pointing to much more extensive supercritical fluid after the Moon-forming impact than previously thought. The structure of the near-critical liquid is very different from what has been studied previously and includes a significant fraction (10%) of molecular species SiO and O2.

  • ab initio molecular dynamics simulation
  • MgSiO3
  • vaporization
  • giant impact
  • density functional theory

Planets present fundamental challenges to our understanding of materials because of the immense range of pressure and temperature relevant to their evolution and formation. The highest temperatures, exceeding even those at the planetary center, are produced in giant impacts, which appear to be a generic part of planetary accretion. A substantial portion of the planet is vaporized in such events (1, 2); in the case of Earth, the Moon may have condensed from an extended fluid cloud (3, 4). Knowing how much vapor is produced and its composition and speciation is thus essential for understanding the initial stages of planetary evolution. Silicate vaporization is important in other planetary contexts as well, for example, in accreting and steadily evaporating rocky exoplanets (5, 6) and in the passage of meteorites through planetary atmospheres (7).

The conditions of planetary vaporization lie in the realm of warm dense matter, defined by conditions of density and temperature such that interactions among molecules are essential and the electronic structure is not yet fully ionized (e.g., 10−3–101 g⋅cm−3 and 0.5–20 eV) (8). Experimental interest has expanded considerably with the development of methods to produce and accurately characterize matter in this regime, for example on adiabatic release from either shocked or isochorically heated states and because of applications to understanding astrophysical objects, plasma production and characterization, and inertial confinement fusion (8, 9). In silicates most dynamic compression studies have focused on achieving the high pressures characteristic of planetary interiors (10), whereas less attention has focused on expanded states and vaporization (11⇓–13). The behavior of silicates under expansion may be particularly rich as experiments near the triple point suggest that a major change in bonding occurs: from dominantly ionically bonded solid and liquid phases to dominantly covalently bonded vapor, containing, for example, O2 and SiO molecules (14).

Despite the importance of silicate vaporization, our knowledge is still limited because of the very high temperature involved; the liquid–vapor coexistence curve has never been experimentally measured for any silicate composition. The conditions of silicate vaporization extend from the solid–liquid–vapor triple point up to the critical point, which is unknown. In SiO2 previous estimates of the critical point range from 5,000 K to 13,000 K and from 0.07 g⋅cm−3 to 1.1 g⋅cm−3, illustrating the magnitude of uncertainty that remains (15). Predictions based on the extrapolation of experimental data (16⇓–18) are accurate near the triple point, but their range of validity at higher temperatures is unknown. Silicate vaporization in hydrodynamic simulations of giant impact events (4, 19) is based on a semianalytic model, which assumes that the liquid–vapor critical point occurs at 8,800 K and that vaporization is congruent, but the accuracy of these assumptions is untested experimentally near the critical point.

Silicate vaporization presents at least two important challenges to simulation methodology. First, the change in bonding on vaporization means that methods based on classical potentials that are typically tuned to one type of bonding will fail to capture the essential physics. Second is the incongruent nature of silicate vaporization: The chemical composition of the liquid and vapor coexisting in equilibrium may differ (20, 21). This means that methods of determining liquid–vapor coexistence based on the Maxwell equal area construction are incomplete as they assume that liquid and vapor have the same chemical composition. There has been one previous molecular dynamics study of silicate vaporization (SiO2 system) (22), which was based on classical pair potentials and which assumed congruent vaporization based on the Maxwell construction.

Here, we use inhomogeneous ab initio molecular dynamics to simulate silicate vaporization. We have chosen the MgSiO3 system as a prototypical binary silicate that has been extensively studied via first-principles molecular dynamics in the liquid phase (23⇓–25) and that is the dominant chemical component of telluric planets, making up, for example, more than 80% of Earth’s silicate fraction (26).

Our method simulates vaporization directly, mimicking a realizable experimental process (Fig. 1). A slab of liquid is placed in a vacuum. During the dynamical trajectory in the canonical ensemble, we see spontaneous evaporation, leading eventually to a steady-state chemical equilibrium between coexisting liquid and vapor. This approach does not assume congruent vaporization as did the previous molecular dynamics study of silicate vaporization (22). Indeed, we find that the chemical compositions of the vapor and the liquid differ significantly. Our simulations allow us to determine the orthobaric densities of coexisting liquid and vapor, their chemical compositions, the speciation of the vapor phase, and the properties of the interface via direct analysis of the trajectories. Previous two-phase molecular dynamics simulations of liquid–vapor coexistence have been based on classical potentials and have been limited to congruently vaporizing systems, including simple fluids and water (27⇓–29). Further details of our simulations are given in Materials and Methods.

Fig. 1.
  • Download figure
  • Open in new tab
  • Download powerpoint
Fig. 1.

Snapshot from a simulation at 5,500 K (270 atoms, Ω= 6) showing Mg (gold), Si (blue), and O (red) atoms. The box outlines the periodically repeated simulation cell. Molecular species are apparent in the vapor phase, including O2, SiO, SiO2, and atomic species, Mg and O.

Fig. 2 shows the results of typical simulations. The density shows the expected spatial variation: from large density within the liquid slab to much lower density in the vapor region. The density of the liquid is reduced at higher temperature, whereas the density of the vapor is augmented. The width of the Gibbs dividing surface grows with increasing temperature. As expected for near critical fluids, density fluctuations in space and time are apparent in both phases and the magnitude of these fluctuations increases with increasing temperature. However, a steady state has been achieved: There are no secular variations in the mean density of either the liquid or the vapor phase.

Fig. 2.
  • Download figure
  • Open in new tab
  • Download powerpoint
Fig. 2.

Simulated mass density profiles at 4,000 K (Top, Ω=4) and 5,500 K (Bottom, Ω=6). We represent temporal and spatial variations in the density with density profiles from snapshots taken from the simulation at times arranged by rainbow color from 2 ps (red) to 5 ps (blue). The density is coarse grained with a Gaussian of width 2.5 Å and the black line is the best fit of Eq. 6 to the time-averaged density profile. Both simulations have N = 270 atoms.

The densities of coexisting liquid and vapor approach each other with increasing temperature (Fig. 3) and the phase diagram can be represented by the Wegner expansion (30) which yields the location of the critical point (technically the maxcondentherm) at ρc=0.48±0.05 g⋅cm−3 and Tc=6,600±150 K. The density of the liquid decreases monotonically with increasing temperature. At 5,500 K the density (1.29±0.11 g⋅cm−3) is half that at the ambient melting point of this system (31) (2.58 g⋅cm−3 at 1,830 K). The density of the coexisting vapor at 5,500 K is 0.14 ± 0.05 g⋅cm−3. The density of the vapor phase increases with increasing temperature. We find no systematic variation of liquid or vapor densities with computational parameters N or Ω, as expected for orthobaric two-phase coexistence (28).

Fig. 3.
  • Download figure
  • Open in new tab
  • Download powerpoint
Fig. 3.

Liquid–vapor phase diagram from 270-atom (solid symbols) and 135-atom (open symbols) simulations, with Ω = 4 (circles) and 6 (squares). Red, liquid; blue, vapor. The green solid circle is the critical point of the liquid–vapor boundary assumed in recent simulations of the Moon-forming impact (4, 19). The black solid line is the best fit to the Wegner expansion: ρl−ρv=C11−T/Tcβ+C21−T/Tcβ+Δ, 12(ρl+ρv)=ρc+C3(1−T/Tc); with values of the exponents fixed to those expected from renormalization group theory, β=0.325, Δ=12 (30); and the other parameters C1=1.32, C2= 2.36, C3= 1.74. Inset shows the pressure of the two-phase system from our simulations (red symbols), the best fit to the simulation results (red curve), the critical pressure (blue symbol), and the extrapolation of low pressure–temperature experimental data described in the text (black dashed curve).

The pressure along the coexistence line increases monotonically with increasing temperature (Fig. 3). We determine the vapor pressure as the component of the pressure tensor normal to the interface (Pzz) (29). The pressure may be represented by P=exp(A−B/T) with A=11.8±2.0 and B=45,000±14,000 with P in megapascals. This fit yields a critical pressure Pc=140 MPa at Tc=6,600 K. Our pressure–temperature curve agrees remarkably well with an extrapolation from experimental observations at lower temperature: Fitting to the triple point from ref. 21 (0.2 Pa, 1,773 K) and the 1-bar boiling point from ref. 32 (0.1 MPa, 3,350 K) yields A = 12.45 and B = 49,420 and a curve that differs from the best fit to the simulation results by no more than a few megapascals. Incongruent vaporization requires that the vapor pressure line is actually a two-phase coexistence region of finite width. However, a thermochemical modeling study of the silica system (15) indicates that the width of the two-phase region is on the order of 10% of the critical pressure and therefore not resolvable in our simulations within our present uncertainties.

The width of the Gibbs dividing surface between liquid and vapor phases increases systematically with increasing temperature (Fig. 4). We have found that our results can be fitted to the expression (33) w=w01−T/Tc−ν which yields w0=1.4 ±0.4 Å, Tc=6,700±300 K, and ν=0.3±0.2. The value of Tc determined from the fit agrees well with the critical temperature determined from the phase diagram (Fig. 3).

Fig. 4.
  • Download figure
  • Open in new tab
  • Download powerpoint
Fig. 4.

The width of the interface (blue) and the compressibility (red, right-hand axis) as a function of temperature.

As temperature increases, spatial fluctuations in the density of the liquid phase increase in magnitude. We quantify density fluctuations and show that their magnitude is thermodynamically sensible (Fig. 4). Density fluctuations of the liquid may be related to the isothermal compressibility and to the structure factor in the limit of small q as⟨N2⟩−⟨N⟩2⟨N⟩=ρNkBTχT=limq→0S(q),[1]where N is the particle number, angle brackets denote time averages, S(q) is the structure factor, χT is the isothermal compressibility, ρN is the number density, T is temperature, and kB is the Boltzmann constant. We compute the structure factor from the definition (34) S(q)=⟨ρq→ρ−q→⟩/N with ρq→=∑i⁡exp(iq→⋅r→i) given by a sum over the positions r→i of the atoms in the liquid phase, i.e., those for which Eq. 7 is satisfied. To focus on the properties of the liquid phase, we select q→ vectors from the plane parallel to the interfaceq→=2πL(h,k,0)[2]with h and k integers and L the cell dimension in the plane of the interface. To find S(q→0) we fit the structure factor to a quadratic in q in the range 0<q<3Å−1. We have found that this procedure is robust for the relatively large compressibilities of our system: The extrapolated value of the structure factor does not differ significantly from that at the smallest q value, S(2π/L).

The compressibility shows a steep increase with increasing temperature from 0.1 GPa−1 at 4,000 K to 0.5 GPa−1 at 6,500 K (Fig. 4). The compressilbity at 4,000 K is two times larger than that measured experimentally near the triple point (0.05 GPa−1) (35).

The structure of the liquid undergoes important changes with increasing temperature along the vapor–liquid coexistence line (Fig. 5). At 4,000 K, the Si-O, Mg-O, and O-O radial distribution functions are very similar to those seen in previous studies of the homogeneous liquid at similar conditions of temperature and density (23, 36). As temperature increases a new peak emerges in the O-O radial distribution function at a very small distance (1.3 Å). This peak corresponds to the appearance of O2 molecules in the liquid phase.

Fig. 5.
  • Download figure
  • Open in new tab
  • Download powerpoint
Fig. 5.

Cation–oxygen coordination numbers (Top) and the abundance of molecular species in the liquid (Bottom) shown as the proportion of O atoms that are members of O2 molecules and the proportion of Si atoms that are members of an SiO molecule. Bottom Inset shows the partial radial distribution functions of the liquid phase at 6,000 K: Si-O (blue), O-O (red), and Mg-O (gold).

We quantify liquid structure via the coordination numberZαβ=4πρNxβ∫0rminr2gαβ(r)dr,[3]where Zαβ is the mean number of atoms of type β surrounding an atom of type α, gαβ is the corresponding partial radial distribution function, xβ is the concentration of atoms of type β, and rmin is the location of the first minimum in gαβ. We also track ζαβ(Z), the fraction of α atoms surrounded by Z β atoms. All liquid-phase structural quantities are computed from positions of atoms that satisfy Eq. 7.

We find that the Si-O (4) and Mg-O (5) coordination numbers at 4,000 K are very similar to those found in previous experimental and theoretical studies near the triple-point density (23, 36) (Fig. 5). With increasing temperature, ZSiO and ZMgO decrease to values near 3 close to the critical point. The proportion of O atoms that are members of an O2 molecule (ζOO(1)) increases steadily with increasing temperature, as does the number of SiO molecules (ζSiO(1)), both approaching 10% near the critical point (Fig. 5).

The vapor is significantly depleted in Mg compared with the liquid (Fig. 6): At 5,000 K, only 7% of the vapor is Mg, compared with 20% for bulk MgSiO3 composition. The vapor is enriched in O compared with the bulk composition. With increasing temperature, the O content in the vapor decreases, and the Mg content increases, so that the vapor composition approaches the bulk composition at the critical point.

Fig. 6.
  • Download figure
  • Open in new tab
  • Download powerpoint
Fig. 6.

Chemical composition of the vapor phase: atomic fractions of O (red), Si (blue), and Mg (gold). Also shown are the predictions of the MAGMA code (dashed lines). The arrows along the right-hand axis indicate the atomic fractions in the bulk (MgSiO3) composition.

The vapor phase in our simulations consists of several clearly defined and long-lived molecular species, the most important of these being SiO and O2. Because molecules interact frequently, particularly near the critical point, through either collisions or exchanges of atoms, molecules cannot be uniquely defined. Definitions based on bond-length criteria are unsuitable as vibrations are large in amplitude and highly anharmonic so that molecular proportions are sensitive to the choice of bond-length cutoffs. Instead, we use a definition of a molecule that does not rely on arbitrary bond-length distance cutoffs. We define a two-atom cluster as a molecule when the atoms are mutual nearest neighbors. We also permit three-atom molecules if (i) an Si atom i and an O atom j are mutual nearest neighbors and (ii) another O atom k has i as its nearest neighbor and k is the second nearest neighbor of i. The concept of mutual nearest neighbors has been used previously to identify dimers in dense fluid hydrogen (37) and generalized to a consideration of Nth mutual neighbors in the context of simple liquids (38).

We find that the vapor is dominated by SiO (33% at 4,500 K) and O2 (33%) (Fig. 7). The next most abundant species are O (12%), Mg (10%), and SiO2 (9%). With increasing temperature SiO and O2 remain the most abundant species, although their relative proportions decrease somewhat, while the fractions of atomic O and Mg increase. The proportion of SiO2 molecules depends nonmonotonically on temperature, first increasing with temperature, up to 5,500 K, and then decreasing on further heating.

Fig. 7.
  • Download figure
  • Open in new tab
  • Download powerpoint
Fig. 7.

Speciation in the vapor phase. In Top, Middle, and Bottom, monatomic species are indicated by open symbols and molecular species by solid symbols. Top Inset shows the lifetimes of the molecular species compared with their vibrational periods indicated by the tags on the right-hand axis and the vertical blue band joining the periods of the stretching and bending mode of SiO2 (39, 40).

Our results highlight the limitations of the concept of a molecule: At the temperatures under consideration molecular lifetimes may be less than one vibrational period (Fig. 7). While O2 and SiO molecules are long-lived, the lifetime of SiO2 is much less than the period of the bending mode and the lifetime of MgO is much less than the period of the stretching mode. The dynamics and geometry of the SiO2 and MgO “molecules,” and thus their energetics, may differ significantly from what would be expected on the basis of extrapolation from low-temperature experimental data.

As the liquid–vapor coexistence curve has never been determined before via experiment or first-principles simulation, there are few data with which to compare our results. The one previous simulation study of silicate vaporization (on SiO2) (22) used classical potentials and a homogeneous simulation approach, thereby assuming congruent vaporization. The assumption of congruency and the inability of classical potentials to capture the stabilization of covalent species in the vapor phase, in addition to the difference in system composition, likely account for the much higher critical temperature found in the previous study (12,000 K vs. Tc=6,600 K in our study).

Experimental data on vapor phase composition in the MgSiO3 system do not exist. However, our lowest-temperature results are in reasonable accord with available experimental data in the vicinity of the triple point. The silica component appears to be more volatile than the magnesia component in MgSiO3 (21) as well as in a wide variety of multicomponent silicates (41), so that evaporation tends to produce vapor that is Mg depleted, as we find.

We may also compare our results with extrapolations of thermochemical data. For this purpose, we have used the MAGMA code (16), computing the composition of the vapor phase in equilibrium with MgSiO3 liquid. The MAGMA code also predicts vapor depleted in Mg, as we find, and vapor speciation dominated by SiO and O2 with O and SiO2 being next most abundant. We find the composition of the vapor changes more rapidly with increasing temperature than predicted by the MAGMA code, which shows more nearly constant vapor-phase composition (Fig. 6). These differences are due in part to the ideal gas law underlying the MAGMA code, which breaks down as the critical point is approached; our results could be used to expand the scope of such modeling efforts. The MAGMA code also makes a series of assumptions regarding the component species in the liquid phase and does not include SiO or O2 as liquid-phase species. The results of our study may be used to improve the liquid-phase model, for example, by providing information on these molecular components. A recent study of the silica system (15) argued that it is important to include SiO and O2 as liquid-phase species, particularly near the critical point.

The temperature of the liquid–vapor critical point that we find (6,600 K) is much less than what has been assumed in hydrodynamic simulations of planetary accretion (8,800 K) (Fig. 3). The much higher critical temperature assumed in hydrodynamic simulations reflects the properties of the underlying semianalytic thermodynamic model, which is based on the MANEOS approach (42), and the lack of experimental data on near-critical vaporization. Our results point toward the production of much greater amounts of supercritical fluid over a wider range of impact scenarios than previously assumed. The proto-Earth and its surrounding postimpact disk may have been supercritical out to greater distances, possibly encompassing the radius at which the Moon formed. The presence of such an extended one-phase structure may have important implications for understanding the similarity between the Earth and the Moon in O and other isotopes (43) and for isotopic fractionations between the two bodies such as those recently observed in potassium (3).

The incongruent nature of vaporization that we find highlights another shortcoming of the phase diagram used in hydrodynamic simulations: Silicate vaporization is assumed to be congruent. This assumption is made for simplicity, but also based on the experimental observation that forsterite Mg2SiO4 vaporizes congruently near the triple point. However, the Mg/Si ratio of forsterite is very different from that of the bulk silicate Earth, which is much closer to that of our system (Mg/Si = 1.3). Moreover, thermodynamic modeling of the Mg2SiO4–SiO2 system, based on experimental observations, suggests that Earth-like compositions vaporize incongruently to Mg-poor vapor at low temperature (21) as we find. We find that the degree of Mg depletion in the vapor depends strongly on temperature and that vapor and liquid have very similar compositions near the critical point.

Our results open avenues for studying silicates in the lower-density portion of the warm dense-matter regime. Our methods are in principle applicable to any silicate composition, including those that have been studied more extensively by experiment, such as SiO2, and multicomponent systems that capture even more of the richness of planetary vaporization. By the same token, current experimental technology is in principle capable of directly testing our current predictions: Shock and release studies have already measured fluid structure via XANES (X-ray absorption near edge structure) down to 0.9 g⋅cm−3 in SiO2 (13), vapor composition in CaMgSi2O6 (11), and the thermodynamic conditions of one point on the liquid–vapor coexistence line in SiO2 (12).

Materials and Methods

Our molecular dynamics simulations are based on density functional theory using the projector augmented wave (PAW) method (44) as implemented in VASP (45), with valence configurations Si (3s23p2), Mg (3s23p0), and O (2s22p4). We use the PBEsol approximation to the exchange-correlation potential (46) as we have previously found it to yield excellent agreement with experimental measurements of condensed-phase physical properties of oxides (47, 48). We sample the Brillouin zone at the Gamma point and use a plane-wave cutoff of 500 eV, which we find yields pressures and energies that are converged to within 0.01 GPa and 2 meV per atom, respectively. We assume thermal equilibrium between ions and electrons via the Mermin functional (49).

The initial condition is a liquid slab embedded in a vacuum. The simulations are performed in the canonical ensemble using the Nosé–Hoover thermostat (50) with a time step of 1 fs for durations of 4–15 ps. We explored a range of system sizes: 135 atoms or 270 atoms and initial vacuum-to-liquid volume ratios Ω = 4 or 6. Quantities are computed as time averages over the stationary portion of the trajectory, taken to the be the last 50%, and uncertainties are computed with the proper non-Gaussian statistics via the blocking method (51). Initial configurations are chosen from equilibrated homogeneous simulations near zero pressure from our previous work (25). This choice of initial condition helps to minimize behavior that we have observed which complicates analysis: a tendency for the liquid slab to split in two, particularly as the critical point is approached. The homogeneous simulations were performed with N = 135. To initiate our N = 270 simulations, we doubled the liquid slab in the direction normal to the interface. Because of the large cell sizes, these are very demanding simulations: The number of plane waves required is as much as 100 times greater than for homogeneous liquid-phase simulations.

To determine the physical properties and composition of the two coexisting phases, we locate the liquid–vapor interface using the method of ref. 52. Briefly, the interface is defined at every time step by the surface s such that ρ¯(s,t)=c where c is a constant, and ρ¯(r,t) is the coarse-grained density field obtained by convolving atomic positions with a Gaussian. The proximity of atom i to the surface isai=[s(t)−ri(t)]⋅n(t)|s(t)=si*(t),[4]where ri is the position of atom i, si*(t) is the point on the surface nearest to atom i, and n is the surface normal in the direction of the density gradient at that point. We have found that this method yields robust location of the dividing interface even in those cases where the surface is time dependent or when multiple liquid–vapor surfaces appear as a result of splitting of the liquid slab. We also locate the Gibbs dividing surface and determine its properties by replacing the instantaneous interface s(t) with the time-averaged surface ⟨s⟩ in the computation of ai (52).

We compute the mean mass density profileρ(z)=1L2∑imiδ(ai−z),[5]where mi is the mass of atom i and L is the dimension of the cell parallel to the interface. We find that the density profile follows the expected form (53)ρ(z)=12(ρl+ρv)+12(ρl−ρv)tanhzw[6]from which we extract our best-fitting values of ρl and ρv, the density of the liquid and vapor phases, respectively, and w, the width of the dividing surface.

Computation of the ai allows us to identify the phase to which each atom belongs at each time step in the simulation. We assume that atoms belonging to the liquid and vapor phases are defined byai>2w[7]andai<−2w,[8]respectively. This identification allows us to determine several properties of the liquid and vapor phases, including composition, speciation, and structure.

Acknowledgments

We thank Prof. Bruce Fegley (Washington University) for kindly providing the MAGMA code, Zhigang Zhang (Chinese Academy of Sciences) for helpful discussions, and K. Kurosawa and J. Connolly for insightful comments on the manuscript. This project is supported by the European Research Council under Advanced Grant 291432 MoltenEarth (FP7). Calculations were carried out using the IRIDIS cluster partly owned by University College London, ARCHER of the United Kingdom national high-performance computing service, and MARENOSTRUM at the Barcelona Supercomputing Center, Spain of the PRACE Consortium.

Footnotes

  • ↵1To whom correspondence should be addressed. Email: l.stixrude{at}ucl.ac.uk.
  • Author contributions: B.X. and L.S. designed research, performed research, analyzed data, and wrote the paper.

  • The authors declare no conflict of interest.

  • This article is a PNAS Direct Submission.

  • Copyright © 2018 the Author(s). Published by PNAS.

This open access article is distributed under Creative Commons Attribution-NonCommercial-NoDerivatives License 4.0 (CC BY-NC-ND).

References

  1. ↵
    1. Pahlevan K,
    2. Stevenson DJ
    (2007) Equilibration in the aftermath of the lunar-forming giant impact. Earth Planet Sci Lett 262:438–449.
    OpenUrl
  2. ↵
    1. Nakajima M,
    2. Stevenson DJ
    (2015) Melting and mixing states of the Earth’s mantle after the Moon-forming impact. Earth Planet Sci Lett 427:286–295.
    OpenUrl
  3. ↵
    1. Wang K,
    2. Jacobsen SB
    (2016) Potassium isotopic evidence for a high-energy giant impact origin of the Moon. Nature 538:487–490.
    OpenUrl
  4. ↵
    1. Lock SJ,
    2. Stewart ST
    (2017) The structure of terrestrial bodies: Impact heating, corotation limits, and synestias. J Geophys Res Planets 122:950–982.
    OpenUrl
  5. ↵
    1. Visscher C,
    2. Fegley B
    (2013) Chemistry of impact-generated silicate melt-vapor debris disks. Astrophysical J Lett 767:L12.
    OpenUrl
  6. ↵
    1. van Lieshout R, et al.
    (2016) Dusty tails of evaporating exoplanets II. Physical modelling of the KIC 12557548b light curve. Astron Astrophys 596:A32.
    OpenUrl
  7. ↵
    1. Bland PA,
    2. Artemieva NA
    (2003) Efficient disruption of small asteroids by Earth’s atmosphere. Nature 424:288–291.
    OpenUrlCrossRefPubMed
  8. ↵
    1. Lee RW, et al.
    (2003) Finite temperature dense matter studies on next-generation light sources. J Opt Soc Am B 20:770–778.
    OpenUrl
  9. ↵
    1. Clerouin J, et al.
    (2012) A database for equations of state and resistivities measurements in the warm dense matter regime. Phys Plasmas 19:082702.
    OpenUrl
  10. ↵
    1. Spaulding DK, et al.
    (2012) Evidence for a phase transition in silicate melt at extreme pressure and temperature conditions. Phys Rev Lett 108:065701.
    OpenUrlCrossRefPubMed
  11. ↵
    1. Kurosawa K, et al.
    (2010) In-situ spectroscopic observations of silicate vaporization due to >10 km/s impacts using laser driven projectiles. Geophys Res Lett 37:L23203.
    OpenUrl
  12. ↵
    1. Kraus RG, et al.
    (2012) Shock vaporization of silica and the thermodynamics of planetary impact events. J Geophys Res Planets 117:E09009.
    OpenUrl
  13. ↵
    1. Denoeud A, et al.
    (2014) Metallization of warm dense SiO2 studied by XANES spectroscopy. Phys Rev Lett 113:116404.
    OpenUrl
  14. ↵
    1. Schick HL
    (1960) A thermodynamic analysis of the high-temperature vaporization properties of silica. Chem Rev 60:331–362.
    OpenUrlCrossRef
  15. ↵
    1. Connolly JAD
    (2016) Liquid-vapor phase relations in the Si-O system: A calorically constrained van der Waals-type model. J Geophys Res Planets 121:1641–1666.
    OpenUrl
  16. ↵
    1. Fegley B,
    2. Cameron AGW
    (1987) A vaporization model for iron silicate fractionation in the Mercury protoplanet. Earth Planet Sci Lett 82:207–222.
    OpenUrl
  17. ↵
    1. Roine A
    (2002) Outokumpu HSC Chemistry for Windows (Outokumpu Research Oy, Pori, Finland), version 5.
  18. ↵
    1. Schaefer L,
    2. Fegley B
    (2004) A thermodynamic model of high temperature lava vaporization on Io. Icarus 169:216–241.
    OpenUrl
  19. ↵
    1. Nakajima M,
    2. Stevenson DJ
    (2014) Investigation of the initial state of the Moon-forming disk: Bridging SPH simulations and hydrostatic models. Icarus 233:259–267.
    OpenUrlCrossRef
  20. ↵
    1. Nakagawa H,
    2. Asano M,
    3. Kubo K
    (1981) Mass-spectrometric study of the vaporization of lithium metasilicate. J Nucl Mater 102:292–297.
    OpenUrl
  21. ↵
    1. Mysen BO,
    2. Kushiro I
    (1988) Condensation, evaporation, melting, and crystallization in the primitive solar nebula–experimental-data in the system MgO-SiO2-H2 to 1.0 × 10−9 bar and 1870 °C with variable oxygen fugacity. Am Mineral 73:1–19.
    OpenUrlAbstract
  22. ↵
    1. Guissani Y,
    2. Guillot B
    (1996) A numerical investigation of the liquid-vapor coexistence curve of silica. J Chem Phys 104:7633–7644.
    OpenUrlCrossRef
  23. ↵
    1. Stixrude L,
    2. Karki B
    (2005) Structure and freezing of MgSiO3 liquid in Earth’s lower mantle. Science 310:297–299.
    OpenUrlAbstract/FREE Full Text
  24. ↵
    1. de Koker N,
    2. Stixrude L
    (2009) Self-consistent thermodynamic description of silicate liquids, with application to shock melting of MgO periclase and MgSiO3 perovskite. Geophys J Int 178:162–179.
    OpenUrlCrossRef
  25. ↵
    1. Karki BB,
    2. Zhang J,
    3. Stixrude L
    (2013) First principles viscosity and derived models for MgO-SiO2 melt system at high temperature. Geophys Res Lett 40:94–99.
    OpenUrl
  26. ↵
    1. Workman RK,
    2. Hart SR
    (2005) Major and trace element composition of the depleted MORB mantle (DMM). Earth Planet Sci Lett 231:53–72.
    OpenUrl
  27. ↵
    1. Rao M,
    2. Levesque D
    (1976) Surface-structure of a liquid-film. J Chem Phys 65:3233–3236.
    OpenUrl
  28. ↵
    1. Alejandre J,
    2. Tildesley DJ,
    3. Chapela GA
    (1995) Fluid-phase equilibria using molecular-dynamics–The surface-tension of chlorine and hexane. Mol Phys 85:651–663.
    OpenUrl
  29. ↵
    1. Bauer BA,
    2. Patel S
    (2009) Properties of water along the liquid-vapor coexistence curve via molecular dynamics simulations using the polarizable TIP4P-QDP-LJ water model. J Chem Phys 131:084709.
    OpenUrlPubMed
  30. ↵
    1. Vega L,
    2. DeMiguel E,
    3. Rull LF,
    4. Jackson G,
    5. McLure IA
    (1992) Phase-equilibria and critical-behavior of square-well fluids of variable width by Gibbs ensemble Monte-Carlo simulation. J Chem Phys 96:2296–2305.
    OpenUrlCrossRef
  31. ↵
    1. Lange RA,
    2. Carmichael ISE
    (1987) Densities of Na2O-K2O-CaO-MgO-FeO-Fe2O3-Al2O3-TiO2-SiO2 liquids–New measurements and derived partial molar properties. Geochim Cosmochim Acta 51:2931–2946.
    OpenUrlCrossRef
  32. ↵
    1. Ahrens TJ,
    2. O’Keefe JD
    (1972) Shock melting and vaporization of lunar rocks and minerals. Earth Moon Planets 4:214–249.
    OpenUrl
  33. ↵
    1. Beysens D,
    2. Robert M
    (1987) Thickness of fluid interfaces near the critical-point from optical reflectivity measurements. J Chem Phys 87:3056–3061.
    OpenUrl
  34. ↵
    1. Hansen JP
    (1973) Statistical-mechanics of dense ionized matter. 1. Equilibrium properties of classical one-component plasma. Phys Rev A 8:3096–3109.
    OpenUrl
  35. ↵
    1. Rivers ML,
    2. Carmichael ISE
    (1987) Ultrasonic studies of silicate melts. J Geophys Res Solid Earth Planets 92:9247–9270.
    OpenUrl
  36. ↵
    1. Sun N,
    2. Stixrude L,
    3. de Koker N,
    4. Karki BB
    (2011) First principles molecular dynamics simulations of diopside liquid at high pressure. Geochim Cosmochim Acta 75:3792–3802.
    OpenUrlCrossRef
  37. ↵
    1. Tamblyn I,
    2. Bonev SA
    (2010) Structure and phase boundaries of compressed liquid hydrogen. Phys Rev Lett 104:065702.
    OpenUrlCrossRefPubMed
  38. ↵
    1. Larsen RE,
    2. Stratt RM
    (1998) Mutual-nearest-neighbor pairs in fluids. Chem Phys Lett 297:211–216.
    OpenUrl
  39. ↵
    1. Irikura KK
    (2007) Experimental vibrational zero-point energies: Diatomic molecules. J Phys Chem Ref Data 36:389–397.
    OpenUrl
  40. ↵
    1. Pacansky J,
    2. Hermann K
    (1978) Ab initio SCF calculations on molecular silicon dioxide. J Chem Phys 69:963–967.
    OpenUrl
  41. ↵
    1. Hashimoto A
    (1983) Evaporation metamorphism in the early solar nebula–Evaporation experiments on the melt FeO-MgO-SiO2-CaO-Al2O3 and chemical fractionations of primitive materials. Geochem J 17:111–145.
    OpenUrlCrossRef
  42. ↵
    1. Melosh HJ
    (2007) A hydrocode equation of state for SiO2. Meteoritics Planet Sci 42:2079–2098.
    OpenUrlCrossRef
  43. ↵
    1. Young ED, et al.
    (2016) Oxygen isotopic evidence for vigorous mixing during the Moon-forming giant impact. Science 351:493–496.
    OpenUrlAbstract/FREE Full Text
  44. ↵
    1. Kresse G,
    2. Joubert D
    (1999) From ultrasoft pseudopotentials to the projector augmented-wave method. Phys Rev B 59:1758–1775.
    OpenUrlCrossRefPubMed
  45. ↵
    1. Kresse G,
    2. Furthmuller J
    (1996) Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys Rev B 54:11169–11186.
    OpenUrl
  46. ↵
    1. Perdew JP, et al.
    (2008) Restoring the density-gradient expansion for exchange in solids and surfaces. Phys Rev Lett 100:136406.
    OpenUrlCrossRefPubMed
  47. ↵
    1. Holmstrom E,
    2. Stixrude L
    (2015) Spin crossover in ferropericlase from first-principles molecular dynamics. Phys Rev Lett 114:117202.
    OpenUrlPubMed
  48. ↵
    1. Scipioni R,
    2. Stixrude L,
    3. Desjarlais MP
    (2017) Electrical conductivity of SiO2 at extreme conditions and planetary dynamos. Proc Natl Acad Sci USA 114:9009–9013.
    OpenUrlAbstract/FREE Full Text
  49. ↵
    1. Mermin ND
    (1965) Thermal properties of inhomogeneous electron gas. Phys Rev 137:A1441–A1443.
    OpenUrlCrossRef
  50. ↵
    1. Hoover WG
    (1985) Canonical dynamics–Equilibrium phase-space distributions. Phys Rev A 31:1695–1697.
    OpenUrlCrossRefPubMed
  51. ↵
    1. Flyvbjerg H,
    2. Petersen HG
    (1989) Error-estimates on averages of correlated data. J Chem Phys 91:461–466.
    OpenUrlCrossRef
  52. ↵
    1. Willard AP,
    2. Chandler D
    (2010) Instantaneous liquid interfaces. J Phys Chem B 114:1954–1958.
    OpenUrlPubMed
  53. ↵
    1. Widom B
    (1982) Potential-distribution theory and the statistical-mechanics of fluids. J Phys Chem 86:869–872.
    OpenUrlCrossRef
PreviousNext
Back to top
Article Alerts
Email Article

Thank you for your interest in spreading the word on PNAS.

NOTE: We only request your email address so that the person you are recommending the page to knows that you wanted them to see it, and that it is not junk mail. We do not capture any email address.

Enter multiple addresses on separate lines or separate them with commas.
Critical vaporization of MgSiO3
(Your Name) has sent you a message from PNAS
(Your Name) thought you would like to see the PNAS web site.
CAPTCHA
This question is for testing whether or not you are a human visitor and to prevent automated spam submissions.
Citation Tools
Critical vaporization of MgSiO3
Bing Xiao, Lars Stixrude
Proceedings of the National Academy of Sciences May 2018, 115 (21) 5371-5376; DOI: 10.1073/pnas.1719134115

Citation Manager Formats

  • BibTeX
  • Bookends
  • EasyBib
  • EndNote (tagged)
  • EndNote 8 (xml)
  • Medlars
  • Mendeley
  • Papers
  • RefWorks Tagged
  • Ref Manager
  • RIS
  • Zotero
Request Permissions
Share
Critical vaporization of MgSiO3
Bing Xiao, Lars Stixrude
Proceedings of the National Academy of Sciences May 2018, 115 (21) 5371-5376; DOI: 10.1073/pnas.1719134115
del.icio.us logo Digg logo Reddit logo Twitter logo CiteULike logo Facebook logo Google logo Mendeley logo
  • Tweet Widget
  • Facebook Like
  • Mendeley logo Mendeley

Article Classifications

  • Physical Sciences
  • Earth, Atmospheric, and Planetary Sciences
Proceedings of the National Academy of Sciences: 115 (21)
Table of Contents

Submit

Sign up for Article Alerts

Jump to section

  • Article
    • Abstract
    • Materials and Methods
    • Acknowledgments
    • Footnotes
    • References
  • Figures & SI
  • Info & Metrics
  • PDF

You May Also be Interested in

Water from a faucet fills a glass.
News Feature: How “forever chemicals” might impair the immune system
Researchers are exploring whether these ubiquitous fluorinated molecules might worsen infections or hamper vaccine effectiveness.
Image credit: Shutterstock/Dmitry Naumov.
Reflection of clouds in the still waters of Mono Lake in California.
Inner Workings: Making headway with the mysteries of life’s origins
Recent experiments and simulations are starting to answer some fundamental questions about how life came to be.
Image credit: Shutterstock/Radoslaw Lecyk.
Cave in coastal Kenya with tree growing in the middle.
Journal Club: Small, sharp blades mark shift from Middle to Later Stone Age in coastal Kenya
Archaeologists have long tried to define the transition between the two time periods.
Image credit: Ceri Shipton.
Illustration of groups of people chatting
Exploring the length of human conversations
Adam Mastroianni and Daniel Gilbert explore why conversations almost never end when people want them to.
Listen
Past PodcastsSubscribe
Panda bear hanging in a tree
How horse manure helps giant pandas tolerate cold
A study finds that giant pandas roll in horse manure to increase their cold tolerance.
Image credit: Fuwen Wei.

Similar Articles

Site Logo
Powered by HighWire
  • Submit Manuscript
  • Twitter
  • Facebook
  • RSS Feeds
  • Email Alerts

Articles

  • Current Issue
  • Special Feature Articles – Most Recent
  • List of Issues

PNAS Portals

  • Anthropology
  • Chemistry
  • Classics
  • Front Matter
  • Physics
  • Sustainability Science
  • Teaching Resources

Information

  • Authors
  • Editorial Board
  • Reviewers
  • Subscribers
  • Librarians
  • Press
  • Cozzarelli Prize
  • Site Map
  • PNAS Updates
  • FAQs
  • Accessibility Statement
  • Rights & Permissions
  • About
  • Contact

Feedback    Privacy/Legal

Copyright © 2021 National Academy of Sciences. Online ISSN 1091-6490