Masanori Matsui

The MD simulation of the equation of state of MgO: Application as a pressure calibration standard at high temperature and high pressure

Abstract Molecular dynamics (MD) simulation is used to calculate the elastic constants and their temperature and pressure derivatives, and the T-P-V equation of state of MgO. The interionic potential is taken to be the sum of pairwise additive Coulomb, van der Waals attraction, and repulsive interactions. In addition, to account for the observed large Cauchy violation of the elastic constants of MgO, the breathing shell model (BSM) is introduced in MD simulation, in which the repulsive radii of O ions are allowed to deform isotropically under the effects of other ions in the crystal. Quantum correction to the MD pressure is made using the Wigner-Kirkwood expansion of the free energy. Required energy parameters, including oxygen breathing parameters, were derived empirically to reproduce the observed molar volume and elastic constants of MgO, and their measured temperature and pressure derivatives as accurately as possible. The MD simulation with BSM is found to be very successful in reproducing accurately the measured molar volumes and individual elastic constants of MgO over a wide temperature and pressure range. The errors in the simulated molar volumes are within 0.3% over the temperature range between 300 and 3000 K at 0 GPa, and within 0.1% over the pressure range from 0 up to 50 GPa at 300 K. The simulated bulk modulus is found to be correct to within 0.7% between 300 and 1800 K at 0 GPa. Here we present the MD simulated T-P-V equation of state of MgO as an accurate internal pressure calibration standard at high temperatures and high pressures

Molecular dynamics study of the structures and bulk moduli of crystals in the system CaO-MgO-Al2O3-SiO2

Molecular dynamics (MD) simulations have been used to calculate the structures and bulk moduli of crystals in the system CaO-MgO-Al2O3-SiO2 (CMAS) using an interatomic potential model (CMAS94), which is composed of pairwise additive Coulomb, van der Waals, and repulsive interactions. The crystals studied, total of 27, include oxides, Mg meta- and ortho-silicates, Al garnets, and various Ca or Al bearing silicates, with the coordination number of cations ranging 6 to 12 for Ca, 4 to 12 for Mg, 4 to 6 for Al, and 4 and 6 for Si. In spite of the simplicity of the CMAS94 potential and the diversity of the structural types treated, MD simulations are quite satisfactory in reproducing well the observed structural data, including the crystal symmetries, lattice parameters, and average and individual nearest neighbour Ca-O, Mg-O, Al-O, and Si-O distances. In addition MD simulated bulk moduli of crystals in the CMAS system compare well with the observed values.

Molecular dynamics simulation of structures, bulk moduli, and volume thermal expansivities of silicate liquids in the system CaO-MgO-Al2O3-SiO2

Materials in the system CaO-MgO-Al2O3-SiO2 are important constituents of the Earths lower crust and mantle. Silicate liquids in this geophysically important system have been studied using molecular dynamics(MD) simulation with an empirical interatomic potential (CMAS94). MD simulations are quite satisfactory in reproducing well the observed structure and pressure-volume-temperature equation-of-state parameters of molten enstatite (MgSiO3), wollastonite (CaSiO3), diopside (CaMgSi2O6), and anorthite (CaAl2Si2O8) at 1900 K and 0 GPa. However, the MD simulated bulk modulus of molten forsterite (Mg2SiO4) at 2300 K and 0 GPa, K0=18.0(6) GPa, is found to be much smaller than the value, K0=∼60 GPa at similar temperature and pressure conditions, estimated previously based on melting curve analyses of forsterite. In an attempt to investigate the possible occurrence of the density inversion between magmatic liquids and residual crystals in the upper mantle conditions, as proposed by Stolper et al. [1981], we have further applied the MD technique with the CMAS94 potential to the diopside system at 1900 K as an example, and have found that the density inversion between crystal and liquid in this system actually occurs at approximately 11 GPa.

Comparison between the Au and MgO pressure calibration standards at high temperature

[1]xa0We compare the MgO pressure scales previously developed using molecular dynamics simulations based on empirical and non-empirical interatomic potentials, to check the accuracy of the reported pressure scales at high temperatures. We find excellent agreement between the two independent pressure scales over wide temperature and pressure ranges, with simulated pressure discrepancies between the two for specified relative volumes being within 0.5 GPa in the temperature and pressure ranges up to 2000 K and 30 GPa. Based on recent simultaneous pressure measurements using both the MgO and Au scales with a multi-anvil apparatus at 1873 K, we then compare the resulting pressures based on the MgO scale with those based on the Andersons Au scale, which is recently widely used to estimate pressures in high temperature and high pressure experiments. We find the Andrersons Au scale underestimates pressures by 1.4 ± 0.3 GPa relative to the MgO scale at 1873 K.

Computer simulation of the MgSiO3 polymorphs

Six polymorphs of MgSiO3 have been studied using molecular dynamic (MD) simulation techniques, based on the empirical potential (MAMOK), which is composed of terms to describe pairwise additive Coulomb, van der Waals attraction, and repulsive interactions. Crystal structures, bulk moduli, volume thermal expansivities, and enthalpies were simulated for the known MgSiO3 polymorphs; orthoenstatite, clinoenstatite, protoenstatite, garnet, ilmenite, and perovskite. The simulated values compare very well with the available experimental data, and the results are quite satisfactory in view of the diversity of the crystal structures of the six polymorphs, the wide range of simulated properties, and the simplicity of the MAMOK potential. MD simulation was further successfully used to study the possibile existence of a post-protoenstatite phase at high temperature, and a C2/c phase at high pressure, both phases being suggested or inferred previously from experimental works.

Breathing shell model in molecular dynamics simulation: Application to MgO and CaO

Molecular dynamics (MD) simulation is used to calculate the elastic constants of both MgO and CaO at zero pressure, and their temperature dependences, as well as the temperature–pressure–volume equation of states of the two oxides. The interionic potential is taken to be the sum of pairwise additive Coulomb, van der Waals, and repulsive interactions. In order to account for the observed large departures from the Cauchy relation of the elastic constants of the two oxides, the breathing shell model (BSM) is introduced in MD simulation, in which the repulsive radii of O ions are allowed to deform isotropically under the effects of other ions in the crystal, with each core and breathing shell being linked by a harmonic spring with force constant k. Required energy parameters, including k, were derived empirically to reproduce the observed molar volumes and elastic constants of the two oxides at ambient conditions, and their temperature dependences as accurately as possible. The MD simulation with BSM is very ...

The case for a body-centered cubic phase (α′) for iron at inner core conditions

Abstract The molecular dynamics (MD) method is used to simulate the structures and elastic constants of liquid and solid iron at Earths core conditions. The configurational energy of the iron system is approximated as a sum of Morse-type pair interactions between Fe ions. The contribution to pressure from the electronic thermal energy is taken from band structure calculations reported for solid iron. The three energy parameters for the Morse potential are determined empirically to reproduce (1) the volume compression data of hexagonal close-packed (hcp) iron up to 304 GPa at 300 K, (2) the density ϱ, bulk modulus K T , and volume thermal expansivities α of liquid iron at 1900 K, (3) the three elastic constants of face-centered cubic (fcc) iron at 1428 K, and (4) the α values of both fcc and hep iron at high pressures. Then we apply the MD method with the optimized potential to simulate ϱ, K S , μ (rigidity), and α of liquid, fcc, hcp, and bcc (body-centered cubic) iron at pressures up to 400 GPa, and temperatures up to 8000 K. At these P-T conditions, the fcc, hcp and bcc iron phases have very similar simulated values for both ϱ and K S , and the fcc and hcp phases also have essentially the same μ values. However, the predicted μ values are found to be much smaller for the bcc phase than for the fcc or hcp phase. We further apply the MD technique to simulate the 200 and 240 GPa phase transitions at high temperatures previously found from shock compression experiments by Brown and McQueen [Brown J.M., McQueen, R.G., 1986. Phase transitions, Gruneisen parameter, and elasticity for shocked iron between 77 GPa and 400 GPa. J. Geophys. Res. 91, 7485–7494]. When we interpret the 200 GPa and the subsequent 240 GPa transitions as the changes from hcp to bcc, and from bcc to liquid, in which the bcc phase and the other phases are considered to be nonmagnetic at these very high-pressure and high-temperature conditions, the MD simulated values of the longitudinal sound velocity drops at the two transitions compare well with the corresponding shock compression data. This prompts us to suggest that the bcc phase, unlike the hcp and fcc phases, duplicates the longitudinal sound velocity in the 200 to 240 GPa region as found by Brown and McQueen. This, in a consequence, suggests possible existence of bcc iron in the inner core conditions.

Computer simulation of the Mg2SiO4 phases with application to the 410 km seismic discontinuity

Abstract Molecular dynamics (MD) simulation is used to predict the structure and elasticity of the olivine (α), modified-spinel (β) and spinel (γ) forms of Mg2SiO4 at high temperatures and high pressures. The interionic potential is taken to be the sum of pairwise additive Coulomb, van der Waals attraction, and repulsive interactions. In order to take account of non-central forces in crystals, the breathing shell model (BSM) is used for simulation, in which the repulsive radii of O ions are allowed to deform isotropically under the effects of other ions in the crystal. The same potential model is used for the three Mg2SiO4 phases. Required energy parameters, including the breathing parameters of O ions, were obtained empirically using the measured structural and elastic properties of the three phases. The MD simulation with BSM is found to be very successful in reproducing accurately the observed values for the three phases, including the structural parameters and individual elastic constants at ambient conditions, the temperature and pressure derivatives of bulk and shear moduli, and the volume thermal expansivity and volume compression over wide temperature and pressure ranges. We further apply MD simulation to predict the density and seismic velocity contrasts between α- and β-Mg2SiO4 at high temperature and high pressure conditions corresponding to the 410 km mantle discontinuity, and compare the simulated results with seismologically observed data. The simulated velocity contrasts for both P- and S-waves support a previous estimate of maximum 50% by volume for the olivine content at the 410 km discontinuity. In contrast, the simulated density difference between the two phases, when compared with the seismologically reported density jump (PREM) at the 410 km discontinuity, requires implausibly too high olivine content of about 90% by volume at this discontinuity. This inconsistency between the simulated density and seismic velocity data may indicate that the bulk composition at the 410 km depth mantle is different from that usually considered such as the pyrolite model, or that the seismic density data at the 410 km mantle contains significant overestimation.

Comparison between the lattice dynamics and molecular dynamics methods: Calculation results for MgSiO3 perovskite

The lattice dynamics (LD) and molecular dynamics (MD) methods have been used to calculate the structure, bulk modulus, and volume thermal expansivity of MgSiO3 perovskite, in order to investigate the reliability of the two simulation techniques over a wide range of temperature and pressure conditions. At an intermediate temperature of 500 K and zero pressure, the LD and MD values are in exellent agreement for both the structure and bulk modulus of MgSiO3 perovskite. At high temperatures and zero pressure, however, the LD method, which is based on the quasi-harmonic approximation, increasingly overestimates the molar volume of MgSiO3 perovskite because of the neglect of higher-order anharmonic terms. At the high temperatures and high pressures prevailing in the lower mantle, the errors in the LD values for both the molar volume and bulk modulus, relative to the MD values, are generally small or negligible. However, since anharmonicity decreases substantially with pressure but increases rapidly with temperature, the error in the LD simulated volume thermal expansivity is serious, especially in the lower pressure region.

Molecular dynamics simulation of MgSiO3 perovskite and the 660-km seismic discontinuity

Molecular dynamics (MD) simulation is applied to MgSiO3 perovskite (Pv), in which breathing shell model (BSM) is developed for oxygen ions to take account of non-central forces in crystals. MD simulation with BSM is found to be very successful in accurately reproducing the observed structural and elastic properties of MgSiO3 Pv over wide temperature and pressure ranges where experimental data are available. Based on newly developed BSM for MgSiO3 Pv together with previously obtained BSM for each MgO periclase (Pc) and Mg2SiO4 spinel (Sp), we use MD simulation to predict the density contrast between (Mg,Fe)2SiO4 Sp and (Mg,Fe)SiO3 Pv plus (Mg,Fe)O magnesiowustite (Mw) at high-temperature and high-pressure conditions corresponding to the 660-km seismic discontinuity, and compare the resulting density contrast with seismologically measured data. Comparison of the MD density contrast (8.4%) with a seismic value (5%) recently obtained by analysing reflected P and S waves from the 660-km discontinuity requires an (Mg,Fe)2SiO4 content at this discontinuity to be 60% by volume, which is in harmony with about 60 vol.% for the (Mg,Fe)2SiO4 content in a pyrolite composition; while the density contrast from PREM (9.3%) demands an unreal (Mg,Fe)2SiO4 content of more than 100% at this discontinuity.