Ab initio Gibbs ensemble Monte Carlo simulations of the liquid-vapor equilibrium and the critical point of sodium.

The ab initio (ai) Gibbs ensemble (GE) Monte Carlo (MC) method coupled with Kohn-Sham density functional theory is successful in predicting the liquid-vapour equilibrium of insulating systems. Here we show that the aiGEMC method can be used to study also metallic systems, where the excited electronic states play an important role and cannot be neglected. For this we include the electronic free energy in the formulation of the effective energy of the system to be used in the acceptance criteria for the MC moves. The application of this aiGEMC method to sodium yields a good agreement with available experimental data on the liquid-vapour equilibrium densities. We predict a critical point for sodium at 2338 ± 108 K and 0.24 ± 0.03 g cm-3. The liquid structure stemming from aiGEMC simulations is very similar to the one from ab initio molecular dynamics. Since this method can determine phase transition without computing the Gibbs free energy, it may offer a new possibility to study other materials with a reasonable computational cost.

Anharmonicity and Phase Diagram of Magnesium Oxide in the Megabar Regime.

With density functional molecular dynamics simulations, we computed the phase diagram of MgO from 50 to 2000 GPa up to 20 000 K. Via thermodynamic integration (TDI), we derive the Gibbs free energies of the B1, B2, and liquid phases and determine their phase boundaries. With TDI and a pseudo-quasi-harmonic approach, we show that anharmonic effects are important and stabilize the B1 phase in particular. As a result, the B1-B2 transition boundary in the pressure-temperature plane exhibits a steep slope. We predict the B1-B2-liquid triple point to occur at approximately T=10000 K and P=370 GPa, which is higher in pressure than was inferred with quasiharmonic methods alone. We predict the principal shock Hugoniot curve to enter the B2 phase stability domain but only over a very small range of parameters. This may render it difficult to observe this phase with shock experiments because of kinetic effects.

Magnesium oxide at extreme temperatures and pressures studied with first-principles simulations.

We combine two first-principles computer simulation techniques, path integral Monte Carlo and density functional theory molecular dynamics, to determine the equation of state of magnesium oxide in the regime of warm dense matter, with densities ranging from 0.35 to 71 g cm-3 and temperatures ranging from 10 000 K to 5 × 108 K. These conditions are relevant for the interiors of giant planets and stars as well as for shock wave compression measurements and inertial confinement fusion experiments. We study the electronic structure of MgO and the ionization mechanisms as a function of density and temperature. We show that the L-shell orbitals of magnesium and oxygen hybridize at high density. This results in a gradual ionization of the L-shell with increasing density and temperature. In this regard, MgO behaves differently from pure oxygen, which is reflected in the shape of the MgO principal shock Hugoniot curve. The curve of oxygen shows two compression maxima, while that of MgO shows only one. We predict a maximum compression ratio of 4.66 to occur for a temperature of 6.73 × 107 K. Finally, we study how multiple shocks and ramp waves can be used to cover a large range of densities and temperatures.

Path integral Monte Carlo simulations of dense carbon-hydrogen plasmas.

Carbon-hydrogen plasmas and hydrocarbon materials are of broad interest to laser shock experimentalists, high energy density physicists, and astrophysicists. Accurate equations of state (EOSs) of hydrocarbons are valuable for various studies from inertial confinement fusion to planetary science. By combining path integral Monte Carlo (PIMC) results at high temperatures and density functional theory molecular dynamics results at lower temperatures, we compute the EOSs for hydrocarbons from simulations performed at 1473 separate (ρ, T)-points distributed over a range of compositions. These methods accurately treat electronic excitation effects with neither adjustable parameter nor experimental input. PIMC is also an accurate simulation method that is capable of treating many-body interaction and nuclear quantum effects at finite temperatures. These methods therefore provide a benchmark-quality EOS that surpasses that of semi-empirical and Thomas-Fermi-based methods in the warm dense matter regime. By comparing our first-principles EOS to the LEOS 5112 model for CH, we validate the specific heat assumptions in this model but suggest that the Grüneisen parameter is too large at low temperatures. Based on our first-principles EOSs, we predict the principal Hugoniot curve of polystyrene to be 2%-5% softer at maximum shock compression than that predicted by orbital-free density functional theory and SESAME 7593. By investigating the atomic structure and chemical bonding of hydrocarbons, we show a drastic decrease in the lifetime of chemical bonds in the pressure interval from 0.4 to 4 megabar. We find the assumption of linear mixing to be valid for describing the EOS and the shock Hugoniot curve of hydrocarbons in the regime of partially ionized atomic liquids. We make predictions of the shock compression of glow-discharge polymers and investigate the effects of oxygen content and C:H ratio on its Hugoniot curve. Our full suite of first-principles simulation results may be used to benchmark future theoretical investigations pertaining to hydrocarbon EOSs and should be helpful in guiding the design of future experiments on hydrocarbons in the gigabar regime.

Equation of state and shock compression of warm dense sodium-A first-principles study.

As one of the simple alkali metals, sodium has been of fundamental interest for shock physics experiments, but knowledge of its equation of state (EOS) in hot, dense regimes is not well known. By combining path integral Monte Carlo (PIMC) results for partially ionized states [B. Militzer and K. P. Driver, Phys. Rev. Lett. 115, 176403 (2015)] at high temperatures and density functional theory molecular dynamics (DFT-MD) results at lower temperatures, we have constructed a coherent equation of state for sodium over a wide density-temperature range of 1.93-11.60 g/cm3 and 103-1.29×108 K. We find that a localized, Hartree-Fock nodal structure in PIMC yields pressures and internal energies that are consistent with DFT-MD at intermediate temperatures of 2×106 K. Since PIMC and DFT-MD provide a first-principles treatment of electron shell and excitation effects, we are able to identify two compression maxima in the shock Hugoniot curve corresponding to K-shell and L-shell ionization. Our Hugoniot curves provide a benchmark for widely used EOS models: SESAME, LEOS, and Purgatorio. Due to the low ambient density, sodium has an unusually high first compression maximum along the shock Hugoniot curve. At beyond 107 K, we show that the radiation effect leads to very high compression along the Hugoniot curve, surpassing relativistic corrections, and observe an increasing deviation of the shock and particle velocities from a linear relation. We also compute the temperature-density dependence of thermal and pressure ionization processes.

The phase diagrams of beryllium and magnesium oxide at megabar pressures.

We performab initiosimulations of beryllium (Be) and magnesium oxide (MgO) at megabar pressures and compare their structural and thermodynamic properties. We make a detailed comparison of our two recently derived phase diagrams of Be (Wuet al2021Phys. Rev.B104014103) and MgO (Soubiran and Militzer 2020Phys. Rev. Lett.125175701) using the thermodynamic integration technique, as they exhibit striking similarities regarding their shape. We explore whether the Lindemann criterion can explain the melting temperatures of these materials through the calculation of the Debye temperature at high pressure. From our free energy calculations, we find that the melting line of both materials is well represented by the Simon-Glazel fitTm(P) =T0(1 +P/a)1/c, whereT0= 1564 K,a= 15.8037 GPa andc= 2.4154 for Be, whileT0= 3010 K,a= 10.5797 GPa andc= 2.8683 for the MgO in the B1. For the B2 phase, we use the valuesa= 26.1163 GPa andc= 2.2426. Both materials exhibit negative Clapeyron slopes on the boundaries between the two solid phases that are strongly affected by anharmonic effects, which also influence the location of the solid-solid-liquid triple point. We find that the quasi-harmonic approximation underestimates the stability range of the low-pressure phases, namely hcp for Be and B1 for MgO. We also compute the phonon dispersion relations at low and high pressure for each of the phases of these materials, and also explore how the phonon density of states is modified by temperature. Finally, we derive secondary shock Hugoniot curves in addition to the principal Hugoniot curve for both materials, and study their offsets in pressure between solid and liquid branches.

Equivalence between pressure- and structure-defined ionization in hot dense carbon.

The determination of the ionization of a system in the hot dense regime is a long standing issue. Recent studies have shown inconsistencies between standard predictions using average atom models and evaluations deduced from electronic transport properties computed with quantum molecular dynamics simulations [Bethkenhagen et al., Phys. Rev. Res. 2, 023260 (2020)]2643-156410.1103/PhysRevResearch.2.023260. Here, we propose a definition of the ionization based on its effect on the plasma structure as given by the pair distribution function (PDF), and on the concept of effective one-component plasma (eOCP). We also introduce a definition based on the total pressure and on a modelization of the electronic pressure. We show the equivalence of these definitions on two studies of carbon along the 100 eV isotherm and the 10 g/cm^{3} isochor. Simulations along the 100 eV isotherm are obtained with the newly implemented Ext. First principles molecular dynamics (Fpmd) method in Abinit for densities ranging from 1 to 500 g/cm^{3}and along the 10 g/cm^{3} isochor with the recently published Spectral quadrature DFT (Sqdft) simulations, between 8 and 860 eV. The resulting ionizations are compared to the predictions of the average-atom code Qaam which is based on the muffin-tin approximation. A disagreement between the eOCP and the actual PDFs (non-OCP behavior) is interpreted as the onset of bonding in the system.

Analyzing Melts and Fluids from Ab Initio Molecular Dynamics Simulations with the UMD Package.

We have developed a Python-based open-source package to analyze the results stemming from ab initio molecular-dynamics simulations of fluids. The package is best suited for applications on natural systems, like silicate and oxide melts, water-based fluids, and various supercritical fluids. The package is a collection of Python scripts that include two major libraries dealing with file formats and with crystallography. All the scripts are run at the command line. We propose a simplified format to store the atomic trajectories and relevant thermodynamic information of the simulations, which is saved in UMD files, standing for Universal Molecular Dynamics. The UMD package allows the computation of a series of structural, transport and thermodynamic properties. Starting with the pair-distribution function it defines bond lengths, builds an interatomic connectivity matrix, and eventually determines the chemical speciation. Determining the lifetime of the chemical species allows running a full statistical analysis. Then dedicated scripts compute the mean-square displacements for the atoms as well as for the chemical species. The implemented self-correlation analysis of the atomic velocities yields the diffusion coefficients and the vibrational spectrum. The same analysis applied on the stresses yields the viscosity. The package is available via the GitHub website and via its own dedicated page of the ERC IMPACT project as open-access package.

First-principles equation of state database for warm dense matter computation.

We put together a first-principles equation of state (FPEOS) database for matter at extreme conditions by combining results from path integral Monte Carlo and density functional molecular dynamics simulations of the elements H, He, B, C, N, O, Ne, Na, Mg, Al, and Si as well as the compounds LiF, B_{4}C, BN, CH_{4}, CH_{2}, C_{2}H_{3}, CH, C_{2}H, MgO, and MgSiO_{3}. For all these materials, we provide the pressure and internal energy over a density-temperature range from ∼0.5 to 50 g cm^{-3} and from ∼10^{4} to 10^{9} K, which are based on ∼5000 different first-principles simulations. We compute isobars, adiabats, and shock Hugoniot curves in the regime of L- and K-shell ionization. Invoking the linear mixing approximation, we study the properties of mixtures at high density and temperature. We derive the Hugoniot curves for water and alumina as well as for carbon-oxygen, helium-neon, and CH-silicon mixtures. We predict the maximal shock compression ratios of H_{2}O, H_{2}O_{2}, Al_{2}O_{3}, CO, and CO_{2} to be 4.61, 4.64, 4.64, 4.89, and 4.83, respectively. Finally we use the FPEOS database to determine the points of maximum shock compression for all available binary mixtures. We identify mixtures that reach higher shock compression ratios than their end members. We discuss trends common to all mixtures in pressure-temperature and particle-shock velocity spaces. In the Supplemental Material, we provide all FPEOS tables as well as computer codes for interpolation, Hugoniot calculations, and plots of various thermodynamic functions.

Electrical conductivity and magnetic dynamos in magma oceans of Super-Earths.

Super-Earths are extremely common among the numerous exoplanets that have been discovered. The high pressures and temperatures in their interiors are likely to lead to long-lived magma oceans. If their electrical conductivity is sufficiently high, the mantles of Super-Earth would generate their own magnetic fields. With ab initio simulations, we show that upon melting, the behavior of typical mantle silicates changes from semi-conducting to semi-metallic. The electrical conductivity increases and the optical properties are substantially modified. Melting could thus be detected with high-precision reflectivity measurements during the short time scales of shock experiments. We estimate the electrical conductivity of mantle silicates to be of the order of 100 Ω-1 cm-1, which implies that a magnetic dynamo process would develop in the magma oceans of Super-Earths if their convective velocities have typical values of 1 mm/s or higher. We predict exoplanets with rotation periods longer than 2 days to have multipolar magnetic fields.

First-principles equation of state and shock compression predictions of warm dense hydrocarbons.

We use path integral Monte Carlo and density functional molecular dynamics to construct a coherent set of equations of state (EOS) for a series of hydrocarbon materials with various C:H ratios (2:1, 1:1, 2:3, 1:2, and 1:4) over the range of 0.07-22.4gcm^{-3} and 6.7×10^{3}-1.29×10^{8}K. The shock Hugoniot curve derived for each material displays a single compression maximum corresponding to K-shell ionization. For C:H = 1:1, the compression maximum occurs at 4.7-fold of the initial density and we show radiation effects significantly increase the shock compression ratio above 2 Gbar, surpassing relativistic effects. The single-peaked structure of the Hugoniot curves contrasts with previous work on higher-Z plasmas, which exhibit a two-peak structure corresponding to both K- and L-shell ionization. Analysis of the electronic density of states reveals that the change in Hugoniot structure is due to merging of the L-shell eigenstates in carbon, while they remain distinct for higher-Z elements. Finally, we show that the isobaric-isothermal linear mixing rule for carbon and hydrogen EOS is a reasonable approximation with errors better than 1% for stellar-core conditions.

Degassing cascades in a shear-thinning viscoelastic fluid.

We report the experimental study of the degassing dynamics through a thin layer of shear-thinning viscoelastic fluid when a constant air flow is imposed at its bottom. The fluid is an aqueous solution of cetyltrimethylammonium bromide (CTAB) and sodium salicylate (NaSal). Over a large range of parameters, the air is periodically released through a series of successive bubbles, hereafter named cascades. Each cascade is followed by a continuous degassing, lasting for several seconds, corresponding to an open channel crossing the fluid layer. The periodicity between two cascades does not depend on the injected flow rate. Inside one cascade, the properties of the overpressure signal associated with the successive bubbles vary continuously. The pressure threshold above which the fluid starts flowing, fluid deformation and pressure drop due to degassing through the thin fluid layer can be simply described by a Maxwell model.