Free energy and configurational entropy of liquid silica: Fragile-to-strong crossover and polyamorphism

Recent molecular dynamics (MD) simulations of liquid silica, using the “BKS” model [Van Beest, Kramer, and van Santen, Phys. Rev. Lett. 64, 1955 (1990)], have demonstrated that the liquid undergoes a dynamical crossover from super-Arrhenius, or “fragile” behavior, to Arrhenius, or “strong” behavior, as temperature $T$ is decreased. From extensive MD simulations, we show that this fragile-to-strong crossover (FSC) can be connected to changes in the properties of the potential energy landscape, or surface (PES), of the liquid. To achieve this, we use thermodynamic integration to evaluate the absolute free energy of the liquid over a wide range of density and $T$. We use this free energy data, along with the concept of “inherent structures” of the PES, to evaluate the absolute configurational entropy $S_c$ of the liquid. We find that the temperature dependence of the diffusion coefficient and of $S_c$ are consistent with the prediction of Adam and Gibbs, including in the region where we observe the FSC to occur. We find that the FSC is related to a change in the properties of the PES explored by the liquid, specifically an inflection in the $T$ dependence of the average inherent structure energy. In addition, we find that the high $T$ behavior of $S_c$ suggests that the liquid entropy might approach zero at finite $T$, behavior associated with the so-called Kauzmann paradox. However, we find that the change in the PES that underlies the FSC is associated with a change in the $T$ dependence of $S_c$ that elucidates how the Kauzmann paradox is avoided in this system. Finally, we also explore the relation of the observed PES changes to the recently discussed possibility that BKS silica exhibits a liquid-liquid phase transition, a behavior that has been proposed to underlie the observed polyamorphism of amorphous solid silica.

Figure 1

State points simulated in the (a) $V\text{−}D$ and (b) $V\text{−}T$ planes. All points give results for equilibrated liquids.

Figure 2

Test of $T$ dependence of basin shape. IS’s from three $T$, and for three isochores, are rapidly heated in order to confine the sampling to a single basin. Using velocity scaling, $T$ is increased from $0$ to $7000\mathrm{K}$ over $100\text{fs}$. Each curve is an average over $10$ runs. The curves for the same isochore are approximately the same, indicating that the anharmonic contributions to the vibrational energy can be assumed, for the present purposes, to be the same for each basin.

Figure 3

$\Omega$ as a function of $T$ for (a) isochore A, (b) isochore D, and (c) isochore I. (See Table for the definition of the isochore labels.) Also shown is the standard deviation about the mean value based on 100 samples. We note that a difference in $\Omega$ of $0.01$ yields a change in entropy of $0.24\mathrm{J}∕\text{mol}\mathrm{K}.$ For the purpose of this work, we therefore consider $\Omega$ to be constant along isochores. (d) $\Omega \left(T=4000\mathrm{K}\right)$ as a function of $V$; this is the value of $\Omega$ used in our calculations.

Figure 4

$\Omega$ as a function of $e_{\text{IS}}$ for (a) isochore A, (b) isochore D, and (c) isochore I. The error bars represent the standard deviation about the mean value based on 100 samples.

Figure 5

Isotherm of $\mid P_{\text{LJ}}^{\text{ex}}\mid$ for the LJ fluid at $T=4000\mathrm{K}$. The solid line is the fit to Eq. (23) with $M=8$. The dashed line is the curve given by Eq. (20). The data are shown on a log-log plot to simplify the comparison of the data to Eq. (20) at large $V$. The cusp near $V=10{\text{cm}}^3∕\text{mol}$ is due to the fact that $P_{\text{LJ}}^{\text{ex}}$ changes sign.

Figure 6

Plot of the integrand $I=⟨{\mathcal{U}}_{\text{BKS}}-{\mathcal{U}}_{\mathrm{LJ}}⟩$, in Eq. (26). The solid line is a cubic spline fit to the data.

Figure 7

Thermodynamic properties of liquid BKS silica along the $T_0=4000\mathrm{K}$ isotherm: (a) $S$, (b) $P$, and (c) $U$.

Figure 8

$U\left(V,T=0\right)$ for quartz, coesite, and stishovite. Solid lines are common tangent constructions used to find the system ground state energy for specific values of $V$ (as indicated by dotted lines).

Figure 9

(a) Isochores of $D$, the diffusion coefficient of Si atoms. (b) Test of the AG relation along isochores of $D$; the legend indicates $\rho$ in $\mathrm{g}∕{\text{cm}}^3$. If the data fall on a straight line, the AG relation is satisfied. For comparison, the inset shows $D_{\mathrm{O}}$, the diffusion coefficient of O atoms, along the $\rho =3.01$ (triangles) $\rho =2.31{\mathrm{g}∕\text{cm}}^3$ isochores (stars).

Figure 10

(a) Isochores of $D∕T$. (b) Test of the AG relation along isochores of $D∕T$; the legend indicates $\rho$ in ${\mathrm{g}∕\text{cm}}^3$. For comparison, the inset shows $D_{\mathrm{O}}$ along the $\rho =3.01$ (triangles) $\rho =2.31{\mathrm{g}∕\text{cm}}^3$ isochores (stars). The lines in the main panel are fits of a straight line to the data, used to obtained the values of $A$ and ${\mu }_0$ shown in Fig. 15.

Figure 11

(a) Arrhenius plot of $D_{\mathrm{O}}$ (open symbols) and $D$ (filled symbols) along the $\rho =3.01$ (circles) and $\rho =2.31{\mathrm{g}∕\text{cm}}^3$ (diamonds) isochores. (b) Ratio $D_{\mathrm{O}}∕D$ as a function of $1∕T$ along the isochores shown in (a).

Figure 12

$e_{\text{IS}}\left(T\right)$ along isochores. The lines are fits to each isochore of a fifth order polynomial in $T$ with no linear term (i.e., zero slope at $T=0$); these are the polynomial curves we use to estimate the integral term in Eq. (9). The legend indicates $\rho$ in ${\mathrm{g}∕\text{cm}}^3$. The inset shows $e_{\text{IS}}$ versus $1∕T$ for three isochores spanning the density range. The low density isochore shows a marked departure from the relation $e_{\mathrm{IS}}\sim 1∕T$ at low $T$. The symbols used in the inset correspond to the same $\rho$ as in the main panel.

Figure 13

Detail of isochoric $e_{\text{IS}}$ behavior for three $\rho$ spanning the range of our calculations. The thick horizontal lines show the value of the $T=0$ crystal energies obtained from Fig. 8.

Figure 14

(a)–(i) $S_c$ along isochores; each panel is labeled by the density in $\mathrm{g}∕{\text{cm}}^3$. Each curve is obtained using Eq. (9), by integrating the fitted curves for $e_{\mathrm{IS}}$ shown in Fig. 12. Dotted curves are fits to a two-state model [43, 44]. (j) $1∕S_c$ as a function of $T$.

Figure 15

Estimates of the parameters (a) ${\mu }_0$ and (b) $A$, that occur in Eq. (1). These estimates are based on the straight-line fits to the data shown in Fig. 10(b).

Figure 16

Isotherms of the $V$ dependence of $D$ (filled circles), $S_c$ (squares), and $S$ (triangles) at various $T$. For plotting purposes, $S$ has been shifted down so that $S\left(V_0\right)=S_c\left(V_0\right)$. The shifts in $S$ for the various panels are (a) $-66.6042\mathrm{J}∕\text{mol}\mathrm{K}$, (b) $-70.6031\mathrm{J}∕\text{mol}\mathrm{K}$, (c) $-77.8919\mathrm{J}∕\text{mol}\mathrm{K}$, and (d) $-84.3673\mathrm{J}∕\text{mol}\mathrm{K}$.

Figure 17

Comparison of $U$ and the contributions to $C_V$. The top panels (a)–(c) show isochores of $U$ spanning the studied density range; also shown are estimates for $U$ of the crystalline state of the system in the harmonic approximation, derived from the $T=0$ estimate of $U$ in Fig. 8 and extended to higher $T$ using a straight line of slope $3R∕2$. The bottom panels (d)–(f) show the contributions to $C_V$. The solid line is $C_V-\left(3R∕2\right)\left(N-1\right)∕N$; the dot-dashed line is $C_V^{\text{anh}}$; and the dashed line is $C_V^{\text{IS}}$. Each curve is obtained by differentiating the function fitted to the data for the corresponding energy.

Figure 18

Location of the line of $C_V$ maxima (asterisks) in the $V\text{−}T$ plane. Also shown are points on the “temperature of maximum density” (TMD) line (triangles), at which the isobaric expansivity changes sign; and the location of maxima in the contribution of $e_{\text{IS}}$ to $C_V$ (squares). The diamond indicates where evidence for liquid-liquid phase separation was found in Ref. 48. The regions labeled I and II are referred to in the text.