We present an equation of state for the solid and liquid phases of lithium fluoride that covers a wide range of conditions from ambient pressure and temperature to the high pressures and temperatures exhibited in shock- and ramp-compression studies. The particular solid phase we have focused on in this work is the B1 phase. We have followed an approach where the pressure and heat-capacity functions of both phases are fit to experimental data and our own quantum molecular dynamics simulations and are then integrated in a thermodynamically consistent way to obtain the corresponding free-energy functions. This approach yields a two-phase equation of state that provides better overall agreement with experimental data than other equations of state for lithium fluoride, such as SESAME 7271v3, LEOS 2240, and the model presented by Smirnov. The last of these is a three-phase equation of state that predicts a B1–B2 transition along the shock Hugoniot at a pressure of about 140 GPa. This solid–solid transition has been a topic of speculation and debate in the literature for over 50 years, culminating in the work of Smirnov, who has developed the only potentially viable equation of state that allows for this transition. We explain why the proposed B1–B2 transition at 140 GPa is not consistent with recent velocimetry data.

Lithium fluoride (LiF) is an alkali halide with a wide range of applications. LiF is an important component of salt mixtures that serve as coolants and carriers for molten salt nuclear reactors because it is relatively non-volatile, chemically inert, and has favorable transport properties.1 Fluoride compounds have been studied in planetary science as surrogates of the oxides and silicates that are of more immediate interest to understanding the behavior of the Earth’s mantle.2 LiF is also commonly used in static-compression experiments [e.g., in diamond anvil cells (DACs)] as both a pressure standard (calibrant) and a pressure-transmitting medium.3–5 In these experiments, a material is hydrostatically compressed at a fixed temperature to reveal information about how its volume varies as a function of pressure along a particular isotherm. The role of the pressure-transmitting medium is to help produce the desired hydrostatic-pressure profile by minimizing pressure gradients to ensure isotropic, uniform loading. Finally, LiF is often selected to be the window material in dynamic-compression experiments (e.g., shock compression described by a Hugoniot curve6 or ramp compression along a quasi-isentrope7). The purpose of the window in these experiments is to help maintain the desired high-pressure state in the sample of interest for a longer duration of time. An important goal of these experiments is to examine the material response of the sample with interferometric velocimetry techniques in which light traverses through the window and reflects off of the sample. LiF is a popular choice of window material in large part because it remains transparent over a wide range of pressures. For this reason, its optical properties are of great interest within high-energy-density science (HEDS).8–15 

It is essential to have a reliable equation of state (EOS) for the applications described above. We present such an EOS for the solid and liquid phases of LiF that covers a wide range of conditions from ambient pressure and temperature to the high pressures (say, a few hundred GPa) and temperatures exhibited in dynamic-compression studies. We develop this EOS by following the same approach we have used to develop a two-phase EOS for water.16 In this approach, we assume particular functional forms for the pressure and the heat capacity, adjust the parameters in these forms to fit experimental data and predictions from atomistic simulations, and integrate the pressure and heat-capacity functions in a thermodynamically consistent way to obtain the corresponding free-energy functions. The solid phase of interest in our work is the B1 (sodium chloride-type) phase, which consists of a face-centered-cubic lattice of lithium cations with an interpenetrating face-centered-cubic lattice of fluoride anions. There has been discussion in the literature10,14,15,17–19 about the possibility of a B1–B2 transition at sufficiently high pressures and temperatures. The B2 (cesium chloride-type) phase has two interpenetrating body-centered-cubic lattices. This solid–solid transition has not been confirmed to occur in any of the published experiments, and in fact, we explain later why it is not consistent with recent velocimetry measurements. For this reason, we have decided to focus only on the B1 and liquid phases in this work.

We have organized this paper in the following manner. Section II outlines our EOS-development methodology and introduces some of the functional forms that we have elected to use. Section III shows our results for the solid (B1) phase, including comparisons with other theoretical models and data from shock-, ramp-, and shock-ramp-compression experiments. Our results for the liquid phase are presented in Sec. IV. This section includes a discussion on the two-phase segment of the principal Hugoniot, which is the portion that overlaps the melt curve and represents states in which the shocked material exists as a two-phase mixture. We show that the widely used SESAME 7271v313,20 and LEOS 224021–23 equations of state do not agree well with experimental determinations of the melt curve, nor do they properly model the two-phase segment of the Hugoniot. To further compare our equation of state with these other available models, we present hydrodynamic simulations using different equations of state where we have shock compressed LiF to pressures as high as 300 GPa. Section IV also analyzes the issue of the B1–B2 transition. Section V concludes this paper with a summary of our main points.

The goal of this study is to develop a model for the Helmholtz energy F = F(V, T) of LiF as a function of the molar volume V and temperature T. We aim to develop two separate free energy functions F(V, T), one for the liquid phase and the other for the solid (B1) phase. We follow the methodology described in our previous work16 on a particular high-pressure solid phase of water known as ice VII. In this approach, we do not prescribe a model for F(V, T) and find other thermodynamic properties such as the pressure P = −(∂F/∂V)T and isochoric heat capacity CV=T(2F/T2)V by taking derivatives of F. Rather, we work with these derivatives directly and integrate them with respect to V or T to infer what the corresponding F(V, T) is. This requires that we assign particular functional forms to P = P(V, T) and CV = CV(V, T), being mindful of the fact that these quantities need to obey thermodynamic consistency relations such as

(1)

and fit the model parameters to available experimental or atomistic simulation data.

We decompose the pressure into two terms:

(2)

where Pisotherm(V) = P(V, Tref) is the pressure along an isotherm at an arbitrary reference temperature Tref and Pthermal(V, T) is the thermal contribution to the pressure. In order to satisfy Eq. (1) and the Maxwell relation (∂P/∂T)V = (∂S/∂V)T, the heat capacity CV = T(∂S/∂T)V and the entropy S = −(∂F/∂T)V must obey16 

(3)
(4)

where Vref is an arbitrary reference volume and CV(Vref, T) is an arbitrary function that depends on temperature only. Furthermore, it can be shown that the Helmholtz energy F must be given by

(5)

in which

(6)
(7)

Here, Fref = F(Vref, Tref) is the reference energy.

Let us examine the particular case where the thermal contribution to pressure assumes the following empirical form in which it varies linearly with temperature:

(8)

where η and Vthermal are constants. As we demonstrate later, our quantum molecular dynamics (QMD) simulations suggest that this particular form is applicable to both phases of LiF for the conditions of interest in this study. Since P has a linear dependence on temperature so that (2P/T2)V=0, we see from (3) that CV in this case can depend only on temperature,

(9)

This simplifying assumption of modeling the thermal pressure Pthermal with a functional form that depends linearly on temperature turns out to be very helpful toward developing the EOS because it decouples the pressure from the heat capacity. This means that the pressure and heat capacity—two important derivatives of the free energy—can be determined independently of each other. Despite its simplicity, we have found that the set of equations above is capable of producing an accurate two-phase EOS for water.16 As a result, we expect that it will also be sufficient for LiF, which is chemically more simple and does not exhibit much of the anomalous behavior famously displayed by water. We expect that Eq. (8) will be valid for most non-reactive, single-component materials over a relatively limited range of conditions. It will not be valid in very high-pressure and high-temperature regimes, such as those examined by Driver and Militzer,24 where the free energy is dominated by contributions from electronic excitations. The restriction to single-component systems also excludes porous materials since the air that fills the pores may be treated as a second component that complicates the pressure–volume–temperature behavior.

In summary, Eqs. (2)–(7) involve four quantities that may be determined independently of each other: Fref, Pisotherm, Pthermal, and CV = CV(T). Once these four quantities have been set, the Helmholtz energy F may be obtained by evaluating Eqs. (5)–(7). It is straightforward to then derive all other thermodynamic properties from F, such as the internal energy E = F + TS or the Gibbs energy G = F + PV. Thus the EOS is completely specified once the four quantities have been determined. The rest of this study explains how they may be determined from available data to develop models for the liquid and solid phases of LiF. Each phase will be described by a different set of Fref, Pisotherm, Pthermal, and CV. Since much more data are available for the solid, we focus on this phase first, before extending the EOS to describe the liquid. We present the following itemized summary of the procedure that we have followed as a reference guide for the reader:

  1. Fit P(V, T) of the solid phase to PV isotherms from static-compression experiments and PT isochores from quantum molecular dynamics simulations conducted as part of this study.

  2. Represent CV(T) of this phase with a Debye-like form and fit this form to experimental heat capacity data.

  3. Set the reference energy Fref to an arbitrarily chosen value. Integrate P(V, T) and CV(T) with respect to volume and temperature according to Eqs. (5)–(7) to obtain the Helmholtz energy F(V, T) of the solid.

  4. Fit P(V, T) of the liquid phase to PT isochores from our quantum molecular dynamics simulations. These atomistic simulations represent the few pieces of information available on the volumetric behavior of this phase; no such static-compression data have been published, for instance.

  5. Find what Fref and CV(T) of the liquid have to be in order to reproduce the experimentally determined LiF melt curve.

  6. Combine P(V, T) and CV(T) together to produce F(V, T) of this phase.

One source of data to which we may fit the pressure P(V, T) = Pisotherm(V) + Pthermal(V, T) is isotherms from static-compression experiments, such as those performed in diamond anvil cells (DACs). Since there is a relative abundance of PV isotherm data at or near room temperature (which we approximate to be 300 K), we choose Tref = 300 K and represent the pressure Pisotherm at this temperature with the third-order Birch–Murnaghan equation,25 

(10)

where Bref is the isothermal bulk modulus at Vref and Bref is the pressure derivative of Bref at Vref. This equation may be rearranged to

(11)

in which the coefficients a5/3, a7/3, and a3 are constants, given by

(12)
(13)
(14)

Equation (11) is more useful than the traditional Birch–Murnaghan form in (10) for practical purposes because it can more readily be differentiated and integrated with respect to volume. We adjust Vref, Bref, and Bref to fit data along the 300 K isotherm (Fig. 1).

FIG. 1.

The pressure–volume–temperature behavior of solid LiF: (a) comparison of static-compression data collected along two isotherms with our EOS fit to these isotherms; (b) a magnified version of (a) at lower pressures where more data are available; (c) comparison of isochores obtained from our quantum molecular dynamics (QMD) simulations with our EOS fit to these results. The five isotherms (300 K, 560–573 K, 873 K, 985 K, 1073 K) depicted in (b) cover the represented range of conditions.3,5,26–28 The QMD simulations are discussed in more detail in the text.

FIG. 1.

The pressure–volume–temperature behavior of solid LiF: (a) comparison of static-compression data collected along two isotherms with our EOS fit to these isotherms; (b) a magnified version of (a) at lower pressures where more data are available; (c) comparison of isochores obtained from our quantum molecular dynamics (QMD) simulations with our EOS fit to these results. The five isotherms (300 K, 560–573 K, 873 K, 985 K, 1073 K) depicted in (b) cover the represented range of conditions.3,5,26–28 The QMD simulations are discussed in more detail in the text.

Close modal

We use higher-temperature isotherms from static-compression experiments as well as results from our own quantum molecular dynamics (QMD) simulations to determine the parameters η and Vthermal in the thermal pressure of Eq. (8). The QMD simulations are a valuable complement since the highest experimental pressure achieved so far for LiF at elevated temperatures is 37 GPa [see Liu et al.3 in Fig. 1(a)], which is rather limited considering our stated interest in having the EOS be applicable to pressures of at least a few hundred GPa. Our QMD simulations involve ab initio density functional theory (DFT) simulations carried out with a revised version of the Perdew–Burke–Ernzerhof (PBE) exchange–correlation functional29 referred to as PBEsol.30 This functional is especially designed for densely packed solids such as LiF. Our system is a supercell composed of 4 × 4 × 4 = 64 atoms. To sample the Brillouin zone, we use the (1/4, 1/4, 1/4) k-point. We have verified that this supercell and k-point sampling lead to converged results for the energy and pressure by comparing the QMD results on the 1000 K isotherm to the results from the Quasi-harmonic model31 at the same temperature using an 8 × 8 × 8 Monkhorst–Pack mesh32 for the k-point integration. In our QMD computations, the electronic-structure calculations from DFT provide the potential energy surface through which the positions of the ions are evolved through NVT molecular dynamics simulations performed with time steps of 0.25 fs and cutoff radii of 1.7 and 1.1 bohr for the lithium and fluoride ions, respectively. Table I lists values for all the aforementioned parameters used in calculating the pressure of the solid phase with our EOS. We note that the QMD results justify our assumption in Eq. (8) that the thermal pressure varies linearly with temperature.

TABLE I.

Parameters for computing the pressure with our solid LiF equation of state.

TrefVrefBrefηVthermal
(K)(cm3/mole)(GPa)Bref(GPa/K)(cm3/mole)
300 9.8 69.0 4.27 −0.005 48.0 
TrefVrefBrefηVthermal
(K)(cm3/mole)(GPa)Bref(GPa/K)(cm3/mole)
300 9.8 69.0 4.27 −0.005 48.0 

Some discussion on the uncertainty in the parameters in the table is warranted. Naturally, this is a reflection of the uncertainty in the data to which the EOS has been fit. For a given pressure along the 300 K isotherm, the reported volume differs by as much as about 5% among the studies depicted in Fig. 1. Therefore it is reasonable to place error bars of ±5% on our Bref and Bref. (There is far less uncertainty in Vref because this parameter determines the density at ambient conditions, which is a well-established quantity.) The uncertainty in the static-compression data is a result of several factors, including the limited precision of the experimental tools, the degree to which hydrostatic conditions are maintained, and the choice of pressure standard and pressure medium. For example, Cynn et al.28 used copper as a pressure standard and did not use a medium for their 300 K isotherm. For their 560 K isotherm, they used platinum for their standard and argon as the medium, and the heating was achieved through a 600 W heating coil. By contrast, the studies by Liu et al.3 and Yagi26 involve gold and sodium chloride pressure standards, respectively. Questions have been raised about Yagi’s choice of pressure standard,3,27 which might explain why their 1073 K isotherm exhibits a large curvature and overlaps with the 985 K isotherm from Liu et al. at higher pressures. The results in Fig. 1 represent our best attempt to fit as much of the data as we can, but it is not possible to fit all of the points due to the inconsistencies among the different data sets. Given this reality, we expect the uncertainty on our η and Vthermal to be at least ±5% as well.

Now that Pisotherm and Pthermal have been established, the next step is to determine the isochoric heat capacity CV = CV(Vref, T). One way to represent this quantity, which depends only on temperature in our EOS, is through the Debye integral

(15)

where R is the gas constant and θD is the Debye temperature. The Debye integral can be approximated fairly well by the much simpler expression

(16)

This functional form is similar to the integral in (15) in that they both ensure that CV does not exceed the value of 6R dictated by the Dulong–Petit law. However, the Dulong–Petit law is based on the presumption that each atom in the LiF lattice may be treated as a three-dimensional harmonic oscillator that vibrates independently of all the other oscillating atoms in the lattice. In reality, there are other effects, such as anharmonicity and electronic excitations (which both become more prominent with temperature), that allow CV to exceed the 6R Dulong–Petit limit. For the temperatures of interest in this study, we expect these effects to be small but not negligible. Therefore we employ a slight variation of (16) for the heat capacity,

(17)

in which a and θD are adjustable constants where the former may be larger than but not very different from 6. We set the values of these constants (a = 6.55 and θD = 234 K) to match isobaric heat capacity CP measurements reported by Douglas and Dever33 and Andersson and Bäckström34 (Fig. 2). The former study has measured CP as a function of temperature at ambient pressure, while the latter has determined CP as a function of pressure up to 1 GPa at 296 K. While it is true that the trend in the data reported by Douglas and Dever cannot be fully captured by a Debye representation of the heat capacity, the maximum deviation between our EOS and their results is less than 4%. Moreover, CV is fairly insensitive to volume in the range of conditions we have examined in our QMD simulations [Fig. 2(c)]. This gives some justification for our choice of functional forms in Eqs. (8) and (17), in which CV depends only on temperature. The maximum deviation compared to the QMD simulations is less than 5%.

FIG. 2.

Heat capacity of solid LiF nondimensionalized by the gas constant R, including comparisons between our EOS and experimental data:33,34 (a) the isobaric heat capacity CP vs. temperature at ambient pressure; (b) CP vs. pressure at 296 K; (c) comparison between CV predicted by our QMD simulations and that obtained from our EOS; and (d) the Debye-like exponential in (17) plotted alongside the polynomial in (18). The two dashed lines in (b) represent the lower- and upper-bound results from Andersson and Bäckström.34 

FIG. 2.

Heat capacity of solid LiF nondimensionalized by the gas constant R, including comparisons between our EOS and experimental data:33,34 (a) the isobaric heat capacity CP vs. temperature at ambient pressure; (b) CP vs. pressure at 296 K; (c) comparison between CV predicted by our QMD simulations and that obtained from our EOS; and (d) the Debye-like exponential in (17) plotted alongside the polynomial in (18). The two dashed lines in (b) represent the lower- and upper-bound results from Andersson and Bäckström.34 

Close modal

Purely for practical purposes [i.e., to avoid having to compute the exponential integral Ei(x) when evaluating the free energy F and entropy S, which are integrals of CV], we approximate the exponential function (17) with a high-order polynomial,

(18)

where the coefficients are given in Table IV of the  Appendix. Figure 2(d) verifies that this polynomial gives a close representation of Eq. (17) over a wide temperature range.

Now that we have established the pressure P(V, T) and the heat capacity CV(T), our EOS for the solid phase will be completely determined once we specify the value of the reference free energy Fref = F(Vref, Tref). Since this reference sets the zero of the energy scale, which is an arbitrary quantity, we choose Fref = 0 for simplicity. If we apply Eqs. (5)–(7), the resulting Helmholtz energy F(V, T) = Fisotherm(V) + Fthermal(V, T) is given by

(19)
(20)

Before proceeding to Sec. III C, we take this opportunity to clarify why we do not set θD in Eq. (17) to the ambient-pressure Debye temperature of 740 K reported in the literature.35 It would be sensible to use this value if we were to entirely represent the free energy with a Debye model. [Furthermore, in order to accommodate high pressures, we should have the Debye temperature be dependent on volume, encapsulate this volume dependence through a Grüneisen parameter Γ(V), and set 740 K as the value of the Debye temperature at V = Vref.] But as evident in Eqs. (19) and (20), our free-energy models include non-Debye contributions, namely, the integrals of Pisotherm and Pthermal with respect to volume. Thus the value of 234 K that we have set for θD in Eq. (17) does not represent the true Debye temperature of solid LiF. The θD that appears in that equation is merely a parameter that we have adjusted to obtain good agreement with QMD results on CV and experimental data on CP. In fact, if we were to set θD in (17) to be 740 K, the CP produced by our EOS would be far too small due to the presence of the non-Debye terms.

The EOS for the solid phase has now been determined since the Helmholtz energy F(V, T) has been established. This EOS cannot fully reproduce the fine-scale behavior within a particular narrow range of conditions (i.e., it will generally not be as accurate as an EOS that is specifically targeted toward those conditions), but we have demonstrated that it gives a good representation of the average behavior over a fairly wide range of conditions. We can apply it to compute various thermodynamic quantities. Of great interest to high-energy-density science (HEDS) are the properties along the principal Hugoniot, which is a curve through V, T space that represents the loci of states behind the shock wave when a material is shock compressed from ambient conditions. The Hugoniot satisfies the so-called Rankine–Hugoniot equation,6 which is an energy-balance equation that relates the states behind and ahead of the shock,

(21)

where V0, P0, and E0 are the volume, pressure, and internal energy of the material ahead of the shock, respectively, which is at ambient conditions. In addition, the shock speed us and particle speed up may obtained by solving

(22)
(23)

in which M is the molecular mass. These two equations represent mass and momentum balances across the shock interface, respectively. Another quantity of interest along the Hugoniot is the bulk sound speed cs, which is given by

(24)

Here, the isentropic bulk modulus BS is calculated from

(25)

Setting P0 = 10−4 GPa = 1 bar and T0 = 300 K so that V0 = V(P0, T0) and E0 = E(P0, T0), we obtain the results shown in Fig. 3.

FIG. 3.

The principal Hugoniot of solid LiF from our EOS and from other studies8,11,17,19,36–39 plotted as a function of different variables: (a) density vs. pressure, (b) temperature vs. pressure, (c) sound speed cs vs. pressure, and (d) shock speed us vs. particle speed up. For clarity, we have labeled whether the results presented in (b) come from experiments or theory. The range for up in (d) corresponds to the range in pressure in the other three figures.

FIG. 3.

The principal Hugoniot of solid LiF from our EOS and from other studies8,11,17,19,36–39 plotted as a function of different variables: (a) density vs. pressure, (b) temperature vs. pressure, (c) sound speed cs vs. pressure, and (d) shock speed us vs. particle speed up. For clarity, we have labeled whether the results presented in (b) come from experiments or theory. The range for up in (d) corresponds to the range in pressure in the other three figures.

Close modal

There is good agreement among all the studies (both theoretical and experimental) portrayed in Fig. 3 for all properties except the temperature. This discrepancy regarding the Hugoniot temperature is an issue not only with LiF but with many other materials as well. It arises due to the fact that: (1) there is virtually no experimental information about the heat capacity at high pressures and temperatures for any material and (2) obtaining reliable temperature measurements in shock-compression studies remains a daunting challenge because it requires a determination of the equilibrium state of the bulk material—and not the surface—on fast (sub-microsecond) time scales. In other words, there is quite a bit of leeway in choosing how to model CV when developing an EOS because there is not much information on CV itself. Furthermore, one usually does not have much information about the temperature, which is the quantity that is most sensitive to CV and can therefore help the most to discriminate among different CV models. To further illustrate this point, Fig. 3 includes results from two other equations of state: SESAME 7271v3, which is the most current version13,20 of SESAME 7271, and LEOS 2240, which is based on the quotidian EOS (QEOS) formalism.21–23 These are perhaps the two most commonly used equations of state for LiF within the HEDS community. Section IV provides a more detailed comparison between our EOS and these other equations of state. For now, we note that there is fairly good agreement among all of these models for the density. Although we have not shown it, there is also good agreement for the sound, shock, and particle speeds presented in Figs. 3(c) and 3(d). But this is not the case with the temperature. For example, the difference in the predicted temperature between our EOS and SESAME 7271v3 at 140 GPa is about 3% (roughly 100 K) even though there is virtually perfect agreement in terms of the density.

Kormer36 had inferred the temperature through radiance measurements conducted via optical pyrometry. The temperatures that this study reports, which represent the only currently available experimental data on the solid-phase Hugoniot temperature of LiF, are significantly higher than predictions from the theoretical studies in Fig. 3(b). In fact, the EOS developed by Kirillov, Kormer, and Sinitsyn37 are far more in line with the other theoretical studies and significantly differ from Kormer’s own experimental results. Interestingly, in a later study40 that we discuss in more detail in Sec. IV, Kormer et al. observed the opposite trend above 10 000 K; the measured temperatures in this later study are much lower than the temperatures predicted by their EOS.37 This EOS takes into account two contributions to the Helmholtz energy: an ionic term that is represented by a Mie–Grüneisen EOS with a constant CV of 6R, and an electronic term that involves a constant bandgap of 11.5 eV. Due to the large value of the bandgap, the electronic contribution is found to be negligible for the temperature range in Fig. 3(b), which is why we have labeled the EOS of Kirillov et al.37 in that figure as simply being a Mie–Grüneisen EOS.

In order to resolve the temperature discrepancy, Kormer worked together with Zel’dovich41 to develop a kinetic model for the effect of shock compression on electron transfer across the bandgap to the conductivity band. This model predicts that at high temperatures (>10 000 K), the experimental temperatures inferred from pyrometry tend to be much lower than EOS predictions because the radiation from the bulk, where electrons are equilibrated with the ionic lattice, is screened by “cold” electrons. The model also predicts that at low temperatures (<4000 K), the experimental temperatures obtained from pyrometry are much higher than EOS predictions because there is some mechanism (perhaps plastic deformation) that results in the production of hot electrons into the conduction band. These electrons are not equilibrated with the lattice due to the fast time scales of the shock compression, and so the temperature inferred from measuring the radiation emitted by these electrons will be significantly higher than the true, equilibrium temperature.

There are other reasons to view Kormer’s results with some caution. Because LiF is a poor light emitter, it is unclear whether the radiation he measured was even emitted by LiF at all; it may have come from elsewhere in his setup. A similar pyrometry study conducted by Radousky and Ross42 on xenon, which also has a large bandgap like LiF, found no evidence of non-equilibrium electronic effects of the kind proposed by Zel’dovich and Kormer. Radousky and Ross report good agreement between their measurements and temperatures predicted by an earlier theoretical model.43 Kormer was also limited in the technology available at the time that his study was conducted: high-resolution, multi-channel pyrometers were not yet available, for instance. In summary, the only reported experimental data on the solid-phase Hugoniot temperature of LiF likely do not represent true readings of the equilibrium temperature. We encourage future experimental studies on this topic in order to provide the crucial information needed to discriminate among the available equations of state.

Another common source of EOS data in HEDS is ramp-compression experiments. The thermodynamic path in ramp compression is nearly isentropic, and for this reason, it is often referred to as a quasi-isentrope. It can be shown that for a single-phase material, the differential change in temperature dT along an arbitrary path is related to its change in entropy dS and pressure dP along that path according to

(26)

Since T and CP are both positive, isentropic compression (dS = 0) results in a smaller temperature increase than does shock compression (dS > 0), as can be seen from comparing Figs. 4(b) with 3(b). Figure 4 shows results for ramp compression of solid LiF initiated from ambient conditions,13 where the resulting path is well-approximated by the principal isentrope. In addition, the figure presents results from shock-ramp studies,44,45 in which the material is shocked to some state on the principal Hugoniot and then ramp compressed to some (usually) much higher pressure. This shock-plus-ramp loading allows one to access high-pressure states that are much colder than the same pressure point on the Hugoniot. Such a loading path is well-suited for the purpose of exploring solidification kinetics,7,46 for example. The maximum deviation between our EOS predictions and that of Davis et al. in Fig. 4(a) is less than 0.5%. A similar level of agreement is obtained with the work of Seagle et al.44 in Fig. 4(c), if we take into account the uncertainty in their results. The maximum deviation in density between our EOS and that of Ali et al.45 in Fig. 4(d) is less than 1%. Our results for the sound speed lie within their reported uncertainties for pressures lower than about 700 GPa. We note that this close agreement has been achieved without fitting any of our EOS parameters to the data in Fig. 4.

FIG. 4.

Comparison of predictions from our solid LiF EOS (solid curves) with data from ramp- and shock-ramp-compression experiments: [(a) and (b)] the principal isentrope in terms of density vs. pressure and temperature vs. pressure, respectively; [(c) and (d)] results for two different shock-ramp studies. We have plotted in (a) the empirical correlation that Davis et al.13 have used to fit their isentrope data. Seagle et al.44 have shock compressed their LiF sample to about 100 GPa and then ramp compressed it using the Z machine at the Sandia National Laboratories47 to a final pressure of nearly 200 GPa. The dotted curves in (c) represent the lower and upper bounds on their uncertainty. Ali et al.45 have shock compressed their sample to roughly 20 GPa and ramp compressed it with a laser drive9 at the National Ignition Facility to a final pressure of over 900 GPa.

FIG. 4.

Comparison of predictions from our solid LiF EOS (solid curves) with data from ramp- and shock-ramp-compression experiments: [(a) and (b)] the principal isentrope in terms of density vs. pressure and temperature vs. pressure, respectively; [(c) and (d)] results for two different shock-ramp studies. We have plotted in (a) the empirical correlation that Davis et al.13 have used to fit their isentrope data. Seagle et al.44 have shock compressed their LiF sample to about 100 GPa and then ramp compressed it using the Z machine at the Sandia National Laboratories47 to a final pressure of nearly 200 GPa. The dotted curves in (c) represent the lower and upper bounds on their uncertainty. Ali et al.45 have shock compressed their sample to roughly 20 GPa and ramp compressed it with a laser drive9 at the National Ignition Facility to a final pressure of over 900 GPa.

Close modal

As in the case of the solid, our first step toward determining the EOS of the liquid phase is to fit the pressure P(V, T) = Pisotherm(V) + Pthermal(V, T) to the functional forms in Eqs. (11) and (8) for Pisotherm(V) and Pthermal(V, T), respectively. We have stated in Sec. II that static-compression data have not been published for liquid LiF. As a result, we rely on QMD simulations to determine P(V, T) of this phase. By fitting the EOS parameters to these simulations (see Table II), we obtain the results shown in Fig. 5. The QMD simulations suggest that the simplifying assumption in Eq. (8) of Pthermal(V, T) depending linearly on temperature is valid for the liquid phase as well.

TABLE II.

Parameters for computing the pressure with our liquid LiF equation of state.

TrefVrefBrefηVthermal
(K)(cm3/mole)(GPa)Bref(GPa/K)(cm3/mole)
1118.0 13.5 26.0 4.33 −0.0074 18.5 
TrefVrefBrefηVthermal
(K)(cm3/mole)(GPa)Bref(GPa/K)(cm3/mole)
1118.0 13.5 26.0 4.33 −0.0074 18.5 
FIG. 5.

The pressure–volume–temperature behavior of liquid LiF: (a) isochores obtained from our QMD simulations (the details of which are discussed in Sec. III A), along with our EOS fit to these results; (b) comparison between our EOS and QMD simulations from Clérouin et al.8 for five different points. The temperatures of these points (which represent the liquid-phase segment of the principal Hugoniot in their study) are 3920 K, 5250 K, 8580 K, 9355 K, and 15 800 K.

FIG. 5.

The pressure–volume–temperature behavior of liquid LiF: (a) isochores obtained from our QMD simulations (the details of which are discussed in Sec. III A), along with our EOS fit to these results; (b) comparison between our EOS and QMD simulations from Clérouin et al.8 for five different points. The temperatures of these points (which represent the liquid-phase segment of the principal Hugoniot in their study) are 3920 K, 5250 K, 8580 K, 9355 K, and 15 800 K.

Close modal

One could choose the reference temperature Tref of the EOS to again be 300 K, but there would be no practical way to validate the parameters at this temperature (even with atomistic simulations) because the liquid is extremely unstable at 300 K and would therefore very quickly transform to solid if somehow undercooled to this temperature. We instead choose Tref to be 1118 K, which is roughly the melting temperature of LiF at ambient pressure. The value of Vref at this temperature can be inferred by applying the Clausius–Clapeyron relation to measurements of the enthalpy of fusion ΔH presented by Douglas and Dever33 together with the slope of the melt curve (which we discuss in Sec. IV B). The slope of the isochores in Fig. 5(a) is determined by adjusting η and Vthermal, while the relative spacing of the isochores is controlled by Bref and Bref. We expect the uncertainty in these four parameters to be at least ±5%, given the sparsity of available information on the pressure–volume–temperature behavior of the liquid phase.

The final step in the development of the liquid EOS is to find CV(T) = CV(Vref, T) of this phase. We do this by following the procedure described in our earlier work.16 Essentially, it involves setting Fref = 0 and adjusting CV(T) so that Fthermal(V, T) in Eq. (7), which is the only contribution to the Helmholtz energy that depends on CV, acquires the behavior needed for the Gibbs energy G of the solid and liquid phases to be equal along the desired melt curve. For this purpose, we have found that it is sufficient to represent the CV-dependent term in Fthermal(V, T) with a third-order polynomial in T,

(27)

where the polynomial coefficients are listed in Table III and Tref = 1118 K, as mentioned above. Once this polynomial has been determined, the EOS of the liquid becomes completely specified by the free-energy functions in Eqs. (5)–(7).

TABLE III.

Polynomial coefficients for computing the heat-capacity dependent contribution to the Helmholtz energy of liquid LiF with Eq. (27).

F0F1F2F3
(kJ/mole)[(kJ/mole)/K][(kJ/mole/)K2][(kJ/mole/)K3]
−32.200 291 −0.092 093 −1.863 249 × 10−5 7.396 478 × 10−10 
F0F1F2F3
(kJ/mole)[(kJ/mole)/K][(kJ/mole/)K2][(kJ/mole/)K3]
−32.200 291 −0.092 093 −1.863 249 × 10−5 7.396 478 × 10−10 

The resulting two-phase EOS yields the results presented in Fig. 6. We obtain good agreement with experimental determinations of the melt curve from Boehler et al.48 and Jackson.2 Both of these experimental studies have obtained their results through static-compression techniques. By contrast, the SESAME 7271v3 and LEOS 2240 equations of state13,20–23 introduced in Sec. III C do not match the experiments well except in certain narrow ranges in pressure where their associated melt curve happens to intersect the data. Perhaps more importantly, SESAME 7271v3 and LEOS 2240 do not correctly represent the two-phase segment of the principal Hugoniot. This portion of the Hugoniot denotes the states where the shocked material exists as a two-phase mixture. Since an EOS represents the equilibrium behavior (i.e., we neglect kinetics7,46), the two-phase segment is manifested when the solid-phase Hugoniot intersects the melt curve, which occurs at around 141 GPa in our EOS. This pressure is referred to as the solidus. As explained by Kormer36 and as illustrated in Fig. 6(b), the two-phase segment of the Hugoniot overlaps the melt curve. More and more liquid is formed at the expense of the solid with increasing pressure along this segment. At some sufficiently high pressure, all of the solid will have transformed to a single-phase liquid. This pressure, which is about 254 GPa in our EOS, is referred to as the liquidus. The Hugoniot again deviates significantly from the melt curve at pressures above the liquidus.

FIG. 6.

Melt curve and principal Hugoniot of LiF: (a) melt curve and solid-phase Hugoniot from our EOS compared with melt-curve data from two static-compression experiments2,48 and melt-curve predictions from the SESAME 7271v3 and LEOS 2240 equations of state, (b) melt-curve data from Boehler et al.48 plotted alongside the melt curve and Hugoniot from three different equations of state (our EOS, SESAME 7271v3, and LEOS 2240), and [(c) and (d)] our EOS predictions for the two-phase segment of the Hugoniot. This segment overlaps the melt curve and represents the loci of shocked states where the material exists as a two-phase mixture. We show the mole fraction of the liquid in (c), while (d) illustrates the density and sound speed of the two-phase mixture, plus the density and sound speed of the individual phases. The curve labeled the “equilibrium Hugoniot” in (b) is comprised of two single-phase segments (solid for P < 141 GPa and liquid for P > 254 GPa) and the two-phase segment (141 GPa < P < 254 GPa). The metastable extension of the solid-phase Hugoniot is included in (b) for comparison purposes.

FIG. 6.

Melt curve and principal Hugoniot of LiF: (a) melt curve and solid-phase Hugoniot from our EOS compared with melt-curve data from two static-compression experiments2,48 and melt-curve predictions from the SESAME 7271v3 and LEOS 2240 equations of state, (b) melt-curve data from Boehler et al.48 plotted alongside the melt curve and Hugoniot from three different equations of state (our EOS, SESAME 7271v3, and LEOS 2240), and [(c) and (d)] our EOS predictions for the two-phase segment of the Hugoniot. This segment overlaps the melt curve and represents the loci of shocked states where the material exists as a two-phase mixture. We show the mole fraction of the liquid in (c), while (d) illustrates the density and sound speed of the two-phase mixture, plus the density and sound speed of the individual phases. The curve labeled the “equilibrium Hugoniot” in (b) is comprised of two single-phase segments (solid for P < 141 GPa and liquid for P > 254 GPa) and the two-phase segment (141 GPa < P < 254 GPa). The metastable extension of the solid-phase Hugoniot is included in (b) for comparison purposes.

Close modal

It is clear from Fig. 6(b) that neither SESAME 7271v3 nor LEOS 2240 properly models the two-phase segment of the Hugoniot because in these equations of state, there is no overlap with the corresponding melt curve. This is entirely expected in the case of LEOS 2240. The QEOS formalism21 on which LEOS 2240 is based is not a true multiphase EOS methodology in the sense described in our study here. For instance, QEOS does not allow for distinct liquid and solid phases because the enthalpy change ΔH and volume change ΔV of melting are both zero. (One can also examine the isotherms produced by this EOS to confirm the lack of phase transitions; there are no density discontinuities along these isotherms.) The melt curve in LEOS 2240, which is modeled with the Lindemann melt criterion,21,49 is therefore a rather artificial construct since it is not a true phase boundary that delineates distinct phases. SESAME 7271v3 also has a Lindemann melt curve, but unlike LEOS 2240, this EOS is a true multiphase equation of state. Nevertheless, the two-phase segment of the Hugoniot in SESAME 7271v3 does not overlap with the melt curve from that EOS.

The Rankine–Hugoniot equations in (21)–(23) must be obeyed along each of the three segments of the Hugoniot that we have described above (solid, liquid, and two-phase). In all cases, E0, V0, and P0 refer to the solid at ambient conditions, and it is clear what E, V, and P represent in the case of the single-phase segments. Along the two-phase segment, E and V denote the internal energy and molar volume, respectively, of the mixture. If we label the mole fraction of phase 1 as x, then

(28)
(29)

where the subscripts indicate the individual phases. For a given pressure P along the two-phase segment (i.e., along the melt curve), we substitute (28) and (29) into (21) and solve (21) to determine the value of x. That is, there is a unique proportion of solid to liquid that defines these two-phase Hugoniot states. We then obtain the shock speed us and particle speed up of the mixture by solving (22) and (23). It turns out that the sound speed cs of a two-phase mixture is not a simple linear combination of the sound speeds of the individual phases. One can rigorously show50 that

(30)

Here, for compactness of notation, we have defined the variables x1 = x, x2 = 1 − x, ΔV = V1V2, and ΔS = S1S2. All the quantities that appear in this equation can be directly evaluated from the EOS with exception of (Sj/P)S. This derivative reflects the change in the entropy of phase j with respect to pressure while holding the total entropy S = xS1 + (1 − x)S2 of the mixture fixed. This derivative must be evaluated iteratively. Figures 6(c) and 6(d) show that as expected, solving the Rankine–Hugoniot equations leads to a result where the mole fraction of the liquid evolves from 0 to 1 with increasing pressure along the two-phase segment. This evolution from single-phase solid to single-phase liquid is reflected in the behavior of mixture properties such as the density and the sound speed.

To further compare our EOS with SESAME 7271v3 and LEOS 2240, we have performed one-dimensional shock-and-release simulations using the single-fluid, multi-material Ares code.51 The simulations model the planar impact of a projectile consisting of a 5 mm lexan sabot (LEOS 5060) and a 1 mm tantalum flyer plate (SESAME 90210a) onto a target consisting of a 0.5 mm aluminum baseplate (LEOS 130) and a 5 mm LiF window. Typical densities are used for all materials. The simulations are run fully Lagrangian using a uniform mesh resolution of 1 μm per zone. No strength, heat-conduction, or phase-transition-kinetics7,46 models are enabled for these example calculations. All simulations are terminated at 500 ns, and a typical time step is 20 ps. These simulations use the monotonic scalar artificial-viscosity model in Ares.52 The most common way in which equations of state are implemented into the simulations is to have Ares read in tabular versions. Toward this purpose, we have generated a tabular version of our EOS using the Multiphase Equation of State (MEOS) table-generating code developed by the Lawrence Livermore National Laboratory. The table indicates the Helmholtz energy of the most stable phase at a particular density and the temperature. All other thermodynamic properties are computed from the free energy. We have used MEOS to create a tabular version of our EOS with a two-dimensional grid of 600 points linearly spaced apart in density from 2.4 to 6.5 g/cm3, and 1200 points linearly spaced apart in temperature from 273 to 10 000 K.

The projectile is initialized with a large velocity directed toward the target to drive a strong shock at the start of the simulation. For each of the three LiF equations of state, we have performed simulations with impact velocities of 8 km/s and 10 km/s. These impact velocities result in shock pressures of approximately 200 and 300 GPa in the LiF window, respectively. Figure 7 shows a comparison of the pressure, density, temperature, and particle velocity recorded by a Lagrangian tracer placed 2 mm into the LiF from the LiF–aluminum interface. When the shock arrives at the tracer in the LiF, a large jump is seen in the quantities plotted in Fig. 7. These values are maintained until the arrival of a release wave originating from the sabot–flyer-plate interface. As seen in Fig. 7, the pressure, density, and particle velocity histories are quite similar for the three LiF equations of state. This is not the case, however, with the temperature. At 200 GPa along the Hugoniot, LiF exists as a two-phase mixture according to our EOS [Fig. 6(b)]. By contrast, SESAME 7271v3 predicts that LiF remains as a single-phase solid, and although LEOS 2240 predicts that it is a single-phase liquid at 200 GPa, this last EOS treats the latent heat ΔH as being zero. As a result, the Hugoniot temperature predicted by our EOS, which has a ΔH whose value at ambient pressure has been adjusted to agree with experimental data,33 is much lower than that predicted by the other two. This temperature disparity continues to persist at 300 GPa, where the pressure is high enough that LiF exists as a single-phase liquid according to all three equations of state.

FIG. 7.

Results from the Ares52 simulations of LiF shock compression described in the text: (a) density vs. time, (b) pressure vs. time, (c) particle velocity vs. time, and (d) temperature vs. time. For each of the three LiF equations of state depicted here, we have simulated shock compression at two different impact velocities: 8 km/s (solid curves) and 10 km/s (dashed curves). These impact velocities result in shock pressures of about 200 and 300 GPa, respectively. The shock reaches the location of the LiF tracer at about 165 ns (190 ns) for the faster (slower) impact scenario.

FIG. 7.

Results from the Ares52 simulations of LiF shock compression described in the text: (a) density vs. time, (b) pressure vs. time, (c) particle velocity vs. time, and (d) temperature vs. time. For each of the three LiF equations of state depicted here, we have simulated shock compression at two different impact velocities: 8 km/s (solid curves) and 10 km/s (dashed curves). These impact velocities result in shock pressures of about 200 and 300 GPa, respectively. The shock reaches the location of the LiF tracer at about 165 ns (190 ns) for the faster (slower) impact scenario.

Close modal

In this section, we address the controversy alluded to in the Introduction, which is the question of whether the B2 solid phase can be thermodynamically stable. First of all, we note that every ab initio DFT study that has investigated this topic10,14,15,19 predicts that the cold free energy (the Gibbs energy at 0 K) of the B2 phase is greater than that of the B1 phase by more than 0.3 eV/atom (about 60 kJ/mole of LiF). A couple of these studies have examined pressures as high as 1000 GPa, and thus there is no question in the literature that the B1 phase is more thermodynamically stable at 0 K over all relevant pressures. The controversy is whether thermal effects can stabilize the B2 phase relative to B1 at some sufficiently high temperature. We review the evidence that has been presented in the literature regarding the possibility of a B1–B2 transition and present our arguments for why we believe that this evidence is not conclusive.

The earliest study we are aware of that has suggested the possibility of this solid–solid transition is a publication by Kormer et al.17 from 1965, where a few different alkali halides are examined along their principal Hugoniot curves to pressures as high as 500 GPa. They report a slight discontinuity in the slope of their linear usup relation over the course of this range. This should be reflected in a curvature change in other quantities as well, such as the density and sound speed.53–55 Their study does not involve X-ray diffraction or any other type of structural characterization, but they suggest the B1–B2 transition as being a possible cause of this discontinuity. (Kormer et al. pointed out, however, that the discontinuity they observe for LiF is much smaller than that exhibited by the other alkali halides in their study, especially sodium chloride. It is sufficiently small that they concede that it would be justifiable to represent the entire range of their data with a single usup relation.) Figure 8 indicates that the slope discontinuity observed by Kormer et al. could be explained by taking into account the melting that occurs along the Hugoniot. At pressures above about 180 GPa, the equilibrium Hugoniot from our EOS has a curvature that is somewhat different from that of the metastable extension of our solid-phase (B1) Hugoniot.

FIG. 8.

The principal Hugoniot of LiF plotted as a function of different variables: (a) density vs. pressure, (b) temperature vs. pressure, (c) bulk sound speed cs vs. pressure, and (d) shock speed us vs. particle speed up. The metastable extension of our B1 solid-phase Hugoniot is included for comparison purposes. We contend that the B2 solid phase is not thermodynamically stable at the conditions depicted in these figures for reasons explained in the text. The results from Clérouin et al.8 represent predictions from their QMD simulations for both the solid and the liquid. The phase boundaries from Smirnov19 represent three different phase transitions: B1–liquid (for pressures <120 GPa), B1–B2 (the lower curve at pressures >120 GPa), and B2–liquid (the upper curve at pressures >120 GPa). The 1969 study by Kormer et al.40 presents two data points for a given pressure [see the filled squares in (b)]: these correspond to temperatures inferred from optical pyrometry in the blue (440 nm wavelength) and red (650 nm) regions of the electromagnetic spectrum. The range for up in (d) corresponds to the range in pressure in the other three figures.

FIG. 8.

The principal Hugoniot of LiF plotted as a function of different variables: (a) density vs. pressure, (b) temperature vs. pressure, (c) bulk sound speed cs vs. pressure, and (d) shock speed us vs. particle speed up. The metastable extension of our B1 solid-phase Hugoniot is included for comparison purposes. We contend that the B2 solid phase is not thermodynamically stable at the conditions depicted in these figures for reasons explained in the text. The results from Clérouin et al.8 represent predictions from their QMD simulations for both the solid and the liquid. The phase boundaries from Smirnov19 represent three different phase transitions: B1–liquid (for pressures <120 GPa), B1–B2 (the lower curve at pressures >120 GPa), and B2–liquid (the upper curve at pressures >120 GPa). The 1969 study by Kormer et al.40 presents two data points for a given pressure [see the filled squares in (b)]: these correspond to temperatures inferred from optical pyrometry in the blue (440 nm wavelength) and red (650 nm) regions of the electromagnetic spectrum. The range for up in (d) corresponds to the range in pressure in the other three figures.

Close modal

In a later publication, Belonoshko et al.18 performed classical molecular dynamics simulations on the melting behavior of LiF. Their simulations predict that at pressures above 100 GPa, the B2 phase is more stable than B1 and that the B2 melt curve is higher in temperature than the B1 melt curve. Based on this result, Belonoshko et al. speculated that allowing for the B1–B2 transition could resolve an apparent discrepancy between the static-compression data from Boehler et al. to which we have fit our EOS melt curve and the shock-compression data from Kormer36 regarding the melt curve. This discrepancy is as follows. In his 1968 study, Kormer reported that the liquidus of the Hugoniot occurs at a pressure and temperature of 280 GPa and 6000 K, respectively. This is much higher than the roughly 4000 K temperature that would be expected from extrapolating the melt curve of Boehler et al. to 280 GPa. (We may consider our melt curve at this pressure as a good approximation to this extrapolation of Boehler et al.) We note that the claims made by Belonoshko et al. are questionable because their simulations do not agree well with the experimental data on the melt curve. For instance, at least three different studies2,33,48 have determined that the melting temperature at ambient pressure is approximately 1118 K, yet Belonoshko et al. predicted that melting at this pressure occurs at around 900 K. The disagreement could be due to the fact that their system contains only 4096 atoms (modern, two-phase simulations involve millions of atoms56) and/or due to inaccuracies with their interatomic potential. This discrepancy grows with pressure; at 100 GPa, their melting temperature is about 1000 K lower than that from Boehler et al.48 As a result of difficulties with measuring melt curves in DACs (e.g., is the sample uniformly heated?), it is certainly possible that the uncertainty in the results of Boehler et al. may be larger than that indicated in Fig. 6(a). However, Boehler et al. mentioned that they determine the melt temperature through two different methods (visual observation and detection of a discontinuous change in the absorption of laser radiation) and that both methods yield identical results. Thus it is unlikely that their results would be in error by more than 30%, which is the magnitude of the uncertainty needed for the results of Belonoshko et al. to fall within the experimental error bars.

Nevertheless, in order to bring this idea of a B1–B2 transition floated by Belonoshko et al. to fruition, Smirnov19 performed a follow-up study in which he models LiF using DFT for the electrons and a Debye model for the lattice. He adjusted the Debye temperatures of the different phases to obtain good agreement with the Hugoniot density–pressure data in Fig. 8(a), with the assumption that the data from the 1965 study by Kormer et al. correspond to the B2 phase. By doing this, he obtained a model which predicts that the B2 phase is more thermodynamically stable than B1 at pressures and temperatures above 120 GPa and 2500 K, respectively. Moreover, the B1 principal Hugoniot in his study intersects his B1–B2 transition curve at a pressure of about 140 GPa, as shown in Fig. 8(b). Smirnov’s models predict that the solid–solid transition that occurs at this pressure involves a decrease in the longitudinal sound speed of more than 10% and an even larger (about 50%) decrease in the tranverse sound speed, leading to a more than 20% decrease in the bulk, thermodynamic sound speed.

We do not believe that the aforementioned B1–B2 transition proposed by the theoretical studies of Smirnov and Belonoshko et al. is justified because it contradicts some recent experimental data. The strongest evidence against it comes from the sound speed measurements by Liu et al.39 They observed the onset of a phase transition at pressures of around 140 GPa along the Hugoniot. According to their study, this pressure range represents a set of two-phase states in which LiF undergoes a phase transition where the longitudinal sound speed approaches that of the bulk, thermodynamic sound speed. The fact that these two velocities become equal to each other (due to a vanishing shear modulus) indicates that the phase transition in question is a shock-induced melt transition and not a solid–solid transition. Furthermore, Liu et al. did not observe the large decrease in the bulk sound speed predicted by Smirnov around this pressure range. Such a noticeable decrease in the sound speed should also be reflected in a sizeable decrease in the shock speed, but this is not observed either in the experiments. For example, the highest up (about 5.8 km/s) data point by Rigg et al.11 [see Fig. 8(d)] corresponds to a pressure of about 200 GPa, yet even at this pressure, there is no indication of the B1–B2 transition. Given the velocimetry data that contradict the existence of this phase transition, we conclude that the B1 solid phase is likely to be more stable than B2 for the conditions shown in Fig. 8.

There is good agreement among the studies in Fig. 8 for all properties except temperature. We have explained in Sec. III C, however, that the temperatures reported by Kormer et al.—which are the only published experimental temperature results on the LiF Hugoniot—are likely to not be reliable. Before concluding this section, we would like to provide some clarification on the temperature reported by Holland and Ahrens,57 which has sometimes been incorrectly interpreted in the literature as being experimental data. The primary goal of Holland and Ahrens is to ascertain the Hugoniot temperature of iron, using LiF as a window, from information about the temperature of the iron/LiF interface (which they infer through radiance measurements) and the temperature of the LiF bulk. Although the details are not clear, their study suggests that they determine the bulk temperature of LiF along its Hugoniot through a Mie–Grüneisen EOS where the Grüneisen parameter is fitted to the well-established trajectory of the LiF Hugoniot in density–pressure space.

We conclude this section by pointing out that the equilibrium Hugoniot produces better agreement with the experimental data than if we were to neglect the phase transition and follow the metastable extension of the solid-phase Hugoniot instead. For instance, taking into account the melting that occurs for pressures above the solidus (≈140 GPa) leads to a small decrease in the slope of our usup relation, allowing for the resulting Hugoniot to lie within the error bars of the highest up data point from Kormer et al. The density of the equilibrium Hugoniot is also higher than that of the metastable solid Hugoniot, yielding better agreement with the 1965 and 1969 studies by Kormer et al. It may seem contradictory that allowing for melting would lead to a higher density. This apparent contradiction can be resolved by realizing that for a given pressure, the temperature of the equilibrium Hugoniot is much lower than that of the metastable solid-phase Hugoniot, as is evident in Fig. 8(b). The temperature disparity arises because the equilibrium Hugoniot accounts for the fact that some of the energy from shock compression goes toward providing the latent heat required to transform the solid to the liquid.

We have developed an equation of state for both the solid (B1) and liquid phases of LiF. Our models for both phases have been fit to isochores from our own quantum molecular dynamics simulations, and in addition, the solid-phase EOS has been fit to isotherms from static-compression experiments (Figs. 1 and 5). Although the data to which we have fit our EOS is limited to a few hundred GPa, our EOS also yields good agreement with experiments in which LiF is ramp compressed along a quasi-isentrope to pressures as high as about 900 GPa (Fig. 4). The heat capacity of the solid is modeled with a Debye-like function that we have calibrated with experimental and simulation data on that property (Fig. 2), while the liquid-phase heat capacity is modeled in a way so as to reproduce experimental data on the melt curve [Fig. 6(a)]. By contrast, the melt curves from the widely used SESAME 7271v313,20 and LEOS 224021–23 equations of state do not match the experiments well except in certain narrow ranges in pressure where their associated melt curve happens to intersect the data.

An important part of this work is to compare our predictions for the principal Hugoniot with results from the literature (Figs. 3 and 6–8). We have examined a number of properties along the Hugoniot: the density, temperature, sound speed, shock speed, and particle speed. Our EOS agrees well with the experimental data on all of these properties except for the temperature. In the solid phase, the various available theoretical models (including our EOS) disagree by as much as 400 K, depending on the pressure [Fig. 3(b)]. But this level of disagreement is much smaller than the discrepancy that these studies have with the work of Kormer et al., who report experimentally determined temperatures that are more than 1000 K higher than what any of the theoretical models predict for a given pressure. In the liquid phase, there is significant disagreement among the various theoretical studies regarding the Hugoniot temperature, in addition to the discrepancy between these studies and the work of Kormer et al.,36,40 as shown in Figs. 6(b) and 8(b). The work of Kormer et al. still remains the only set of experimental studies on the LiF Hugoniot temperature despite the fact that it was conducted in the 1960s. We have explained in Sec. III C why we believe that their work should be viewed with some caution. In fact, they themselves have expressed skepticism about whether their radiance measurements actually reflect the equilibrium, bulk Hugoniot temperature. A theoretical kinetics model was developed by Zel’dovich, Kormer, and Urlin to give plausible explanations for what the measurements may have detected if it is not the radiation from the bulk material.41 

Developing experimental methods to accurately and precisely measure the temperature of materials under dynamic compression remains an active area of research.58,59 This issue is closely related to the fact that there is virtually no experimental information about the heat capacity at high pressures and temperatures for any material. Given this lack of information, the various theoretical models have chosen to model the heat capacity in different ways, which has lead to the large spread described above in the predictions for the Hugoniot temperature. It is clear that reliable measurements are needed to further constrain these models. We encourage such future experiments in order to obtain this crucial information.

Nevertheless, the information that is currently available does provide some basis to choose among the different models. We have elected to fit our EOS to experimental data on the melt curve, including an extrapolation of these data beyond the 100-GPa limit of the currently available measurements.48 (It would certainly be useful in the future to have data beyond 100 GPa to validate this extrapolation.) By doing so, our resulting EOS happens to agree with recent results from Liu et al.,39 who report the onset of melting along the Hugoniot at a pressure of about 140 GPa. We can see from Fig. 6(b) that the SESAME 7271v3 and LEOS 2240 equations of state do not agree with that of Liu et al. since their Hugoniot curves do not intersect their melt curves at the appropriate pressure. Moreover, these equations of state do not properly model the Hugoniot because the two-phase segment of this curve in both of these models does not overlap the corresponding melt curve.

We have also addressed the issue of the B1–B2 transition that has been speculated in the literature. This solid–solid transition was first suggested as a possibility by Kormer et al.,17 who notice a slight discontinuity in the slope of their Hugoniot over the range from 0 to 500 GPa. We have demonstrated that the melting which occurs along the Hugoniot does result in a slight change in curvature (Fig. 8). Thus one does not need to postulate a solid–solid transition to explain the observations from Kormer et al. In fact, this Hugoniot that takes into account melting (we have referred to this curve as the “equilibrium Hugoniot”) yields better agreement with the experimental data than if we were to neglect the melting and follow the metastable extension of the solid-phase Hugoniot instead. We have also explained in Sec. IV D that the particular B1–B2 transition suggested by Smirnov,19 which was inspired by earlier work from Belonoshko et al.,18 is likely to not be valid since it does not agree with velocimetry data from a few different studies11,17,39 along the Hugoniot.

This work was performed under the auspices of the U.S. Department of Energy by the Lawrence Livermore National Laboratory under Contract No. DE-AC52-07NA27344. We acknowledge informative discussions with L. X. Benedict, R. Chau, N. C. Holmes, R. A. Managan, J. M. Montgomery, and H. B. Radousky. We thank E. Chen, J. Leek, C. J. Prisbrey, and C. J. Wu for help with the Multiphase Equation of State (MEOS) table-generating code. We also thank J. L. Belof and P. A. Sterne for help with obtaining the SESAME 7271/7271v3 tables and the accompanying documentation. We are grateful to S. J. Ali, D. E. Fratanduono, L. E. Kirsch, and R. F. Smith for shock-ramp-compression data obtained from the National Ignition Facility.

TABLE IV.

Coefficients for computing the isochoric heat capacity CV of solid LiF with the polynomial in (18). The units are such that the resulting polynomial evaluates to CV/R, where R is the gas constant.

CoefficientValue
C0 3.283 789 698 998 
C1 7.107 689 614 097 × 10−3 
C2 −8.163 106 441 508 × 10−6 
C3 5.591 980 980 961 × 10−9 
C4 −2.472 769 836 397 × 10−12 
C5 7.446 338 187 500 × 10−16 
C6 −1.583 860 723 719 × 10−19 
C7 2.438 088 233 136 × 10−23 
C8 −2.755 804 802 825 × 10−27 
C9 2.300 698 341 833 × 10−31 
C10 −1.414 461 220 788 × 10−35 
C11 6.317 748 058 643 × 10−40 
C12 −1.992 438 441 636 × 10−44 
C13 4.203 112 948 497 × 10−49 
C14 −5.319 483 121 230 × 10−54 
C15 3.052 977 702 298 × 10−59 
CoefficientValue
C0 3.283 789 698 998 
C1 7.107 689 614 097 × 10−3 
C2 −8.163 106 441 508 × 10−6 
C3 5.591 980 980 961 × 10−9 
C4 −2.472 769 836 397 × 10−12 
C5 7.446 338 187 500 × 10−16 
C6 −1.583 860 723 719 × 10−19 
C7 2.438 088 233 136 × 10−23 
C8 −2.755 804 802 825 × 10−27 
C9 2.300 698 341 833 × 10−31 
C10 −1.414 461 220 788 × 10−35 
C11 6.317 748 058 643 × 10−40 
C12 −1.992 438 441 636 × 10−44 
C13 4.203 112 948 497 × 10−49 
C14 −5.319 483 121 230 × 10−54 
C15 3.052 977 702 298 × 10−59 
1.
H.
Luo
,
S.
Xiao
,
S.
Wang
,
P.
Huai
,
H.
Deng
, and
W.
Hu
, “
Molecular dynamics simulation of diffusion and viscosity of liquid lithium fluoride
,”
Comput. Mater. Sci.
111
,
203
208
(
2016
).
2.
I.
Jackson
, “
Phase relations in the system LiF–MgF2 at elevated pressures: Implications for the proposed mixed-oxide zone of the earth’s mantle
,”
Phys. Earth Planet. Inter.
14
,
86
94
(
1977
).
3.
J.
Liu
,
L.
Dubrovinsky
,
T.
Boffa Ballaran
, and
W.
Crichton
, “
Equation of state and thermal expansivity of LiF and NaF
,”
High Pressure Res.
27
,
483
489
(
2007
).
4.
I.
Uts
,
K.
Glazyrin
, and
K. K. M.
Lee
, “
Effect of laser annealing of pressure gradients in a diamond-anvil cell using common solid pressure media
,”
Rev. Sci. Instrum.
84
,
103904
(
2013
).
5.
H.
Dong
,
S. M.
Dorfman
,
C. H.
Holl
,
Y.
Meng
,
V. B.
Prakapenka
,
D.
He
, and
T. S.
Duffy
, “
Compression of lithium fluoride to 92 GPa
,”
High Pressure Res.
34
,
39
48
(
2014
).
6.
Y. B.
Zel’dovich
and
Y. P.
Raizer
,
Physics of Shock Waves and High-Temperature Hydrodynamic Phenomena
, Annotated Edition (
Dover Publications
,
Mineola
,
2012
), p.
946
.
7.
P. C.
Myint
and
J. L.
Belof
, “
Rapid freezing of water under dynamic compression
,”
J. Phys. Condens. Matter
30
,
233002
(
2018
).
8.
J.
Clérouin
,
Y.
Laudernet
,
V.
Recoules
, and
S.
Mazevet
, “
Ab initio study of the optical properties of shocked LiF
,”
Phys. Rev. B
72
,
155122
(
2005
).
9.
D. E.
Fratanduono
,
T. R.
Boehly
,
M. A.
Barrios
,
D. D.
Meyerhofer
,
J. H.
Eggert
,
R. F.
Smith
,
D. G.
Hicks
,
P. M.
Celliers
,
D. G.
Braun
, and
G. W.
Collins
, “
Refractive index of lithium fluoride ramp compressed to 800 GPa
,”
J. Appl. Phys.
109
,
123521
(
2011
).
10.
Z.-H.
Sun
,
J.
Dong
, and
Y.-W.
Xia
, “
First-principles calculations of the structural, electronic, and optical properties of LiF up to 300 GPa
,”
Physica B
406
,
3660
3665
(
2011
).
11.
P. A.
Rigg
,
M. D.
Knudson
,
R. J.
Scharff
, and
R. S.
Hixson
, “
Determining the refractive index of shocked [100] lithium fluoride to the limit of transmissibility
,”
J. Appl. Phys.
116
,
033515
(
2014
).
12.
C. D.
Spataru
,
L.
Shulenburger
, and
L. X.
Benedict
, “
Ab initio many-body Green’s function calculations of optical properties of LiF at high pressures
,”
Phys. Rev. B
92
,
245117
(
2015
).
13.
J.-P.
Davis
,
M. D.
Knudson
,
L.
Shulenburger
, and
S. D.
Crockett
, “
Mechanical and optical response of [100] lithium fluoride to multi-megabar dynamic pressures
,”
J. Appl. Phys.
120
,
165901
(
2016
).
14.
R. E.
Jones
and
D. K.
Ward
, “
Estimates of crystalline LiF thermal conductivity at high temperature and pressure by a Green–Kubo method
,”
Phys. Rev. B
94
,
014309
(
2016
).
15.
X.-W.
Sun
,
Z.-J.
Liu
,
W.-L.
Quan
,
T.
Song
,
R.
Khenata
, and
S.
Bin-Omran
, “
High-pressure and high-temperature physical properties of LiF studied by density functional theory calculations and molecular dynamics simulations
,”
J. Phys. Chem. Solids
116
,
209
215
(
2018
).
16.
P. C.
Myint
,
L. X.
Benedict
, and
J. L.
Belof
, “
Free energy models for ice VII and liquid water derived from pressure, entropy, and heat capacity relations
,”
J. Chem. Phys.
147
,
084505
(
2017
).
17.
S. B.
Kormer
,
M. V.
Sinitsyn
,
A. I.
Funtikov
,
V. D.
Urlin
, and
A. V.
Blinov
, “
Investigation of the compressibility of five ionic compounds at pressures up to 5 Mbar
,”
Sov. Phys. JETP
20
,
811
819
(
1965
), see http://www.jetp.ac.ru/cgi-bin/e/index/e/20/4/p811?a=list.
18.
A. B.
Belonoshko
,
R.
Ahuja
, and
B.
Johansson
, “
Molecular dynamics of LiF melting
,”
Phys. Rev. B
61
,
11928
(
2000
).
19.
N. A.
Smirnov
, “
Ab initio calculations of the thermodynamic properties of LiF crystal
,”
Phys. Rev. B
83
,
014109
(
2011
).
20.
S. D.
Crockett
and
S.
Rudin
, “
Lithium fluoride equations of state (SESAME 7271)
,” Technical Report No. LA-UR-06-8401,
Los Alamos National Laboratory
,
2006
, see http://permalink.lanl.gov/object/view?what=info:lanl-repo/lareport/LA-UR-06-8401.
21.
R. M.
More
,
K. H.
Warren
,
D. A.
Young
, and
G. B.
Zimmerman
, “
A new quotidian equation of state (QEOS) for hot dense matter
,”
Phys. Fluids
31
,
3059
3078
(
1988
).
22.
D. A.
Young
,
Phase Diagrams of the Elements
(
University of California Press
,
Berkeley and Los Angeles
,
1991
).
23.
D. A.
Young
and
E. M.
Corey
, “
A new global equation of state model for hot, dense matter
,”
J. Appl. Phys.
78
,
3748
3755
(
1995
).
24.
K. P.
Driver
and
B.
Militzer
, “
First-principles simulations of warm dense lithium fluoride
,”
Phys. Rev. E
95
,
043205
(
2017
).
25.
A.
Patiño Douce
,
Thermodynamics of the Earth and Planets
(
Cambridge University Press
,
Cambridge
,
2011
), p.
722
.
26.
T.
Yagi
, “
Experimental determination of thermal expansivity of several alkali halides at high pressures
,”
J. Phys. Chem. Solids
39
,
563
571
(
1978
).
27.
R.
Boehler
and
G. C.
Kennedy
, “
Thermal expansion of LiF at high pressures
,”
J. Phys. Chem. Solids
41
,
1019
1022
(
1980
).
28.
H.
Cynn
,
Z.
Jenei
,
M. J.
Lipp
, and
W. J.
Evans
, “
Pressure-volume isotherms of lithium fluoride at room temperature and at 560 K
” (unpublished).
29.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
, “
Generalized gradient approximation made simple
,”
Phys. Rev. Lett.
77
,
3865
3868
(
1996
).
30.
J. P.
Perdew
,
A.
Ruzsinszky
,
G. I.
Csonka
,
O. A.
Vydrov
,
G. E.
Scuseria
,
L. A.
Constantin
,
X.
Zhou
, and
K.
Burke
, “
Restoring the density-gradient expansion for exchange in solids and surfaces
,”
Phys. Rev. Lett.
100
,
136406
(
2008
).
31.
M. T.
Dove
,
Introduction to Lattice Dynamics
, 1st ed. (
Cambridge University Press
,
Cambridge
,
2005
), p.
280
.
32.
H. J.
Monkhorst
and
J. D.
Pack
, “
Special points for Brillouin-zone integrations
,”
Phys. Rev. B
13
,
5188
5192
(
1976
).
33.
T. B.
Douglas
and
J. L.
Dever
, “
Lithium fluoride: Heat content from 0 to 900°, the melting point and heat of fusion
,”
J. Am. Chem. Soc.
76
,
4824
4826
(
1954
).
34.
S.
Andersson
and
G.
Bäckström
, “
Thermal conductivity and heat capacity of single-crystal LiF and CaF2 under hydrostatic pressure
,”
J. Phys. C.: Solid State Phys.
20
,
5951
5962
(
1987
).
35.
L. E. A.
Jones
, “
High-temperature behaviour of the elastic moduli of LiF and NaF: Comparison with MgO and CaO
,”
Phys. Earth Planet. Inter.
13
,
105
118
(
1976
).
36.
S. B.
Kormer
, “
Optical study of the characteristics of shock-compressed condensed dielectrics
,”
Sov. Phys. Usp.
11
,
229
254
(
1968
).
37.
G. A.
Kirillov
,
S. B.
Kormer
, and
M. V.
Sinitsyn
, “
Equilibrium radiation of shock-compressed ionic crystals
,”
JETP Lett.
7
,
290
291
(
1968
), see http://www.jetpletters.ac.ru/ps/1696/article_25789.shtml.
38.
39.
Q.
Liu
,
X.
Zhou
,
X.
Zeng
, and
S. N.
Luo
, “
Sound velocity, equation of state, temperature and melting of LiF single crystals under shock compression
,”
J. Appl. Phys.
117
,
045901
(
2015
).
40.
S. B.
Kormer
,
M. V.
Sinitsyn
, and
A. I.
Kuryapin
, “
Nonequilibrium radiation from shock compressed ionic crystals at temperatures above 1 eV. I
,”
Sov. Phys. JETP
28
,
852
854
(
1969
), see http://jetp.ac.ru/cgi-bin/dn/e_028_05_0852.pdf.
41.
Y. B.
Zel’dovich
,
S. B.
Kormer
, and
V. D.
Urlin
, “
Nonequilibrium radiation from shock compressed ionic crystals at temperatures above 1 eV. II
,”
Sov. Phys. JETP
28
,
855
859
(
1969
), see http://www.jetp.ac.ru/cgi-bin/dn/e_028_05_0855.pdf.
42.
H. B.
Radousky
and
M.
Ross
, “
Shock temperature measurements in high density fluid xenon
,”
Phys. Lett. A
129
,
43
46
(
1988
).
43.
M.
Ross
and
A. K.
McMahan
, “
Condensed xenon at high pressure
,”
Phys. Rev. B
21
,
1658
1664
(
1980
).
44.
C. T.
Seagle
,
J.-P.
Davis
, and
M. D.
Knudson
, “
Mechanical response of lithium fluoride under off-principal dynamic shock-ramp loading
,”
J. Appl. Phys.
120
,
165902
(
2016
).
45.
S. J.
Ali
,
L. E.
Kirsch
,
D. G.
Braun
,
R. F.
Smith
,
D. E.
Fratanduono
, et al., “
Shock-ramp compression of lithium fluoride to 900 GPa
” (unpublished).
46.
P. C.
Myint
,
A. A.
Chernov
,
B.
Sadigh
,
L. X.
Benedict
,
B. M.
Hall
,
S.
Hamel
, and
J. L.
Belof
, “
Nanosecond freezing of water at high pressures: Nucleation and growth near the metastability limit
,”
Phys. Rev. Lett.
121
,
155701
(
2018
).
47.
C. A.
Hall
,
J. R.
Asay
,
M. D.
Knudson
,
W. A.
Stygar
,
R. B.
Spielman
,
T. D.
Pointon
,
D. B.
Reisman
,
A.
Toor
, and
R. C.
Cauble
, “
Experimental configuration for isentropic compression of solids using pulsed magnetic loading
,”
Rev. Sci. Instrum.
72
,
3587
3595
(
2001
).
48.
R.
Boehler
,
M.
Ross
, and
D. B.
Boercker
, “
Melting of LiF and NaCl to 1 Mbar: Systematics of ionic solids at extreme conditions
,”
Phys. Rev. Lett.
78
,
4589
4592
(
1997
).
49.
F. A.
Lindemann
, “
The calculation of molecular vibration frequencies
,”
Physik. Z
11
,
609
612
(
1910
).
50.
P. C.
Myint
,
L. X.
Benedict
,
C. J.
Wu
, and
J. L.
Belof
, “
The sound speed in two-phase mixtures
” (unpublished).
51.
R. M.
Darlington
,
T. L.
McAbee
, and
G.
Rodrigue
, “
A study of ALE simulations of Rayleigh–Taylor instability
,”
Comput. Phys. Commun.
135
,
58
73
(
2001
).
52.
R. B.
Christensen
, “
Godunov methods on a staggered mesh—An improved artificial viscosity
,” Technical Report No. UCRL-JC-105269,
Lawrence Livermore National Laboratory
,
1990
, see https://e-reports-ext.llnl.gov/pdf/219547.pdf.
53.
J. M.
Brown
and
R. G.
McQueen
, “
Melting of iron under core conditions
,”
Geophys. Res. Lett.
7
,
533
536
, https://doi.org/10.1029/gl007i007p00533 (
1980
).
54.
J. H.
Nguyen
and
N. C.
Holmes
, “
Melting of iron at the physical conditions of the Earth’s core
,”
Nature
427
,
339
342
(
2004
).
55.
J. H.
Nguyen
,
M. C.
Akin
,
R.
Chau
,
D. E.
Fratanduono
,
W. P.
Ambrose
,
O. V.
Fat’yanov
,
P. D.
Asimow
, and
N. C.
Holmes
, “
Molybdenum sound velocity and shear modulus softening under shock compression
,”
Phys. Rev. B
89
,
174109
(
2014
).
56.
L. A.
Zepeda-Ruiz
,
B.
Sadigh
,
A. A.
Chernov
,
T.
Haxhimali
,
A.
Samanta
,
T.
Oppelstrup
,
S.
Hamel
,
L. X.
Benedict
, and
J. L.
Belof
, “
Extraction of effective solid–liquid interfacial free energies for full 3D solid crystallites from equilibrium MD simulations
,”
J. Chem. Phys.
147
,
194704
(
2017
).
57.
K. G.
Holland
and
T. J.
Ahrens
, “
Properties of LiF and Al2O3 to 240 GPa for metal shock temperature measurements
,” in
Properties of Earth and Planetary Materials at High Pressure and Temperature
, Geophysical Monograph (
American Geophysical Union
,
1998
), Chap. 34, pp.
335
343
.
58.
Y.
Ping
,
F.
Coppari
,
D. G.
Hicks
,
B.
Yaakobi
,
D. E.
Fratanduono
,
S.
Hamel
,
J. H.
Eggert
,
J. R.
Rygg
,
R. F.
Smith
,
D. C.
Swift
,
D. G.
Braun
,
T. R.
Boehly
, and
G. W.
Collins
, “
Solid iron compressed up to 560 GPa
,”
Phys. Rev. Lett.
111
,
065501
(
2013
).
59.
M. C.
Akin
,
P. C.
Myint
,
E. L.
Shi
, and
R.
Chau
, “
The temperature of Fe at 3 Mbar
,” in
APS March Meeting
,
Los Angeles, CA, USA
,
2018
, see https://meetings.aps.org/Meeting/MAR18/Session/K38.8.