The dynamics of an inertial confinement fusion capsule, or of a stellar or planet interior, obey a very similar set of equations: magneto-radiative-hydrodynamic equations. The solutions of these equations, however, depend entirely on the transport properties associated with the different materials at play. To properly model the dynamics of these systems, it is necessary to determine with high accuracy the transport coefficients of several materials over a large range of thermodynamic conditions. Experimental capabilities in this respect are still limited due to the nature of the microphysics at play and the extreme conditions involved. Numerical simulations are thus necessary, and in this respect, molecular dynamics simulations based on density functional theory offer exquisite possibilities to constrain the transport properties in the warm to hot dense matter regime. In this paper, we report the methodology used to extract different transport properties based on molecular dynamics performed with the software Abinit. The examples shown are based on the specific cases identified for the purpose of the second charged-particle transport code comparison workshop.

Topics
Density functional theory, Molecular dynamics, Thermal conductivity, Electrical conductivity, Electronic transport, Transport properties, Planetary interiors, Plasma confinement, Hydrodynamics, Classical statistical mechanics

The modeling of high energy density phenomena, such as inertial confinement fusion (ICF) capsule implosion,1 planetary interiors,2 and stellar interiors,3 relies on magneto-radiation-hydrodynamic simulations. However, to be representative of the microscale physics, it is crucial to capture the proper transport properties of the systems from the initial conditions and throughout their evolution. Transport phenomena play a very important role in the aforementioned processes. For instance, heat conductivity in ICF cavities not only determines the efficiency of the energy transfer and deposit within the cavity4 but also impacts the growth of hydrodynamic instabilities at the ablator–fuel interface.5 In planetary interiors, electrical conductivity is key to characterize the existence of a dynamo process,6 and interdiffusion and self-diffusion coefficients are of first importance to determine the cooling efficiency.7–9 In stars, mixing between the core and the stably stratified layers is driven by convective overshoot and the transport coefficients strongly influence the efficiency of the mixing of the energy transport.10 

Density and temperature conditions can vary significantly whether we consider planets, stars, or ICF capsules; the spatial and temporal scales vary significantly as well. The microphysics at play is, however, similar. The use of theoretical models, such as Spitzer–Härm, for thermal conductivities,11 Chapman–Enskog for viscosities,12 or Drude for electrical conductivities,13 provides asymptotic limits. Serious efforts have been made in the recent past to improve our understanding and quantitative characterization of the transport coefficients in the warm dense regime. Average-atom models provided a fair evaluation of electrical conductivities for single species systems.14 To include mixing effects as well as ionic contributions for the viscosity and diffusion coefficients, N-body simulation tools are required.

The advent of molecular dynamics simulations based on density functional theory (DFT), whether using an orbital-free or a Mermin–Kohn–Sham (KS) scheme, offered an exquisite opportunity to explore the N-body effects on transport properties. For instance, orbital-free molecular dynamics simulations of asymmetric mixtures allowed to exhibit very different diffusion mechanisms depending on the concentration of light and heavy elements.15 At much lower temperatures, KS-DFT molecular dynamics shed some light on several processes relevant for planetary interiors. Diffusion coefficients have been determined for the interior of gaseous giants and allowed to constrain the efficiency of the core erosion in Jupiter.16 They also allowed to provide an explanation of the depth of the zonal winds on Uranus and Neptune by invoking ionic electrical conductivities.17 KS-DFT also allowed the computation of precise electrical conductivities for several mixtures,18,19 which can be used in geodynamo simulations.

Until recently, KS-DFT molecular dynamics simulations could only be performed at moderate temperatures because of the orbital wall20 faced by this method. The development of hybrid methods, such as the mixed deterministic and stochastic DFT approach21 or the extended method,22–24 alleviates this issue. It is now possible to simulate warm dense to hot dense matter in a computationally tractable manner using KS-DFT. This provides insight into the role of electronic structure calculations in the computation of transport properties.

In the present paper, we show the current capabilities of KS-DFT molecular dynamics in computing transport coefficients for warm dense matter using the Abinit code.25–27 We explain the strategies developed to design the simulations, to extract, and to post-process the data so as to provide the relevant transport coefficients. We illustrate the discussion using different cases as defined for the second charged-particle transport code comparison workshop.28 

To compute the different transport coefficients, we relied on molecular dynamics (MD) simulations in the isokinetic ensemble. The ions are advanced in time with a Verlet algorithm, and we ensured conservation of temperature by utilizing a velocity rescaling thermostat:29 both of which are implemented in Abinit. The time step was adapted to ensure that the trajectory of each atom was properly sampled with respect to the chemical bond oscillation frequency and to the mean-free path. The simulation cell contained 64–256 atoms. Previous studies showed that such a simulation cell would have limited finite-size effects on bulk properties at least [such as equation of states (EOS)].30 Ionic transport properties are more prone to finite-size effects,31,32 but they should be acceptable within the bounds of the other uncertainties. However, small simulation cells allow for longer simulations which contributes to reducing the uncertainty on the ionic transport coefficients as discussed in the following sections. Our simulations were at least 3000 time steps long (up to 15 000 for the longest).

While the nuclei are all treated as classical particles, the electronic density is self-consistently computed at each time step using the Kohn–Sham density functional theory (KS-DFT)33 with the Mermin finite-temperature functional form.34 Following previous work, we opted for either the generalized gradient approximation (GGA) Perdew, Burke and Ernzerhof exchange-correlation functional35 or for a local density approximation (LDA) (see Table I for further details). These functionals are originally designed for 0 K calculations which constitute one source of uncertainty. Temperature-dependent exchange-correlation functionals are under active development,36–39 and codes like Abinit must be adapted to make use of such exchange-correlation functionals. We used specifically designed projector augmented-wave (PAW) pseudo-potentials.27,40 They were designed using the Atompaw code41 and have been utilized for similar plasma conditions as presented here (see, e.g., Refs. 42–44) We used short cutoff radii to prevent the PAW spheres from overlapping including most of the electrons in the valence (see Table I for details). With such hard pseudo-potentials, the energy cutoff of the plane wave expansion had to be set above 50 Ha. To sample the Brillouin zone, we used the Baldereschi mean-value point.45 We chose the number of Kohn–Sham orbitals N band KS in order to have a last level occupancy of ∼ 10 − 5⁠. For the highest temperatures, we used the extended DFT version of the Abinit code23 to increase the efficiency of the calculations by reducing the number of Kohn–Sham orbitals. The extended DFT model replaces high energy orbitals with pure plane waves, assuming their electrons as free fermions, while keeping a Kohn–Sham description of inner shell electrons with the full plane wave basis set. Much fewer numbers of KS orbitals need to be considered, reducing the computational cost of the calculation. It allowed for simulations containing more particles that were carried out for longer time lengths at multi-eV temperatures.

TABLE I.

KS-DFT molecular dynamics input parameters for ionic transport properties analysis. Mean-value and gamma points are, respectively, referred to as MP and Γ. GGA stands for the Perdew–Burke–Ernzerhof (PBE) exchange-correlation functional.35  N band KS stands for the number of KS orbitals (excluding extended DFT pure plane waves). E cut PW stands for the energy cutoff of the plane wave expansion. r PAW is the radius of the PAW spheres, N electron valence is the number of electrons kept in the valence (while other electrons are frozen in the core), and ⟨ V max ovlp ⟩ stands for the mean voluminal overlap ratio of the PAW spheres.

H1 C1 CH1 Al1 Cu1 HCu1 CH2 Au1 C2
ρ (g/cm3)  1.0  10.0  1.0  2.7  8.97  1.8  0.9  19.32  100.0 
k B T (eV)  2.0  2.0  2.0  1.0  1.0  1.0  7.8  10.0  0.2 
N atom  256 H  64 C  48 H–48 C  108 Al  108 Cu  50 H/50 Cu  48 H/48 C  108 Au  108 C 
k →-point(s)  MP  23  23  MP  23  MP  23  Γ  33 
N band KS  512  512  512  1024  2048  1024  1024  3072  512 
E cut PW (Ha)  62.0  320.0  60.0  50.0  30.0  15.0  60.0  50.0  400.0 
Δ t (atu)  5.0  8.0  5.0  100.0  50.0  50.0  3.0  10.0  5.0 
XC Func.  GGA  GGA  GGA  GGA  LDA  GGA  GGA  LDA  GGA 
r PAW (in Bohrs)  0.7  0.3  0.8–1.1  1.5  2.0  1.0–2.0  0.8–1.1  2.25  0.3 
N electron valence/atom  1  6  1–4  11  19  1–19  1–4  33  6 
N timestep  18 787  7475  22 494  4148  3359  1429  8334  5948  5370 
Total time (ps)  2.272  1.446  2.721  10.034  4.063  1.728  0.605  1.439  0.649 
⟨ V max ovlp ⟩    0.00%  8.76%  0.00%  6.04%  20.17%  15.58%  13.12%  0.00% 
Notes            M H ≡ M Cu  Ext. DFT  Ext. DFT  FCC crystal 
H1 C1 CH1 Al1 Cu1 HCu1 CH2 Au1 C2
ρ (g/cm3)  1.0  10.0  1.0  2.7  8.97  1.8  0.9  19.32  100.0 
k B T (eV)  2.0  2.0  2.0  1.0  1.0  1.0  7.8  10.0  0.2 
N atom  256 H  64 C  48 H–48 C  108 Al  108 Cu  50 H/50 Cu  48 H/48 C  108 Au  108 C 
k →-point(s)  MP  23  23  MP  23  MP  23  Γ  33 
N band KS  512  512  512  1024  2048  1024  1024  3072  512 
E cut PW (Ha)  62.0  320.0  60.0  50.0  30.0  15.0  60.0  50.0  400.0 
Δ t (atu)  5.0  8.0  5.0  100.0  50.0  50.0  3.0  10.0  5.0 
XC Func.  GGA  GGA  GGA  GGA  LDA  GGA  GGA  LDA  GGA 
r PAW (in Bohrs)  0.7  0.3  0.8–1.1  1.5  2.0  1.0–2.0  0.8–1.1  2.25  0.3 
N electron valence/atom  1  6  1–4  11  19  1–19  1–4  33  6 
N timestep  18 787  7475  22 494  4148  3359  1429  8334  5948  5370 
Total time (ps)  2.272  1.446  2.721  10.034  4.063  1.728  0.605  1.439  0.649 
⟨ V max ovlp ⟩    0.00%  8.76%  0.00%  6.04%  20.17%  15.58%  13.12%  0.00% 
Notes            M H ≡ M Cu  Ext. DFT  Ext. DFT  FCC crystal 
View Large

In order to compute the electronic transport properties, we had to select several ionic configurations from the molecular dynamics trajectories and perform a more converged DFT calculation. For instance, we increased the Brillouin zone sampling to a 23 Monkhorst–Pack grid.46 We also needed to increase the number of Kohn–Sham orbitals to include a large number of weakly occupied states. For these calculations we opted out of the extended DFT and performed a full Kohn–Sham calculation. We then applied a Kubo–Greenwood formalism to the eigenstates and eigenvalues to extract the real part of the electrical conductivity.

Molecular dynamics directly provides the time evolution of the energy and of the stresses. To use that as an evaluation of the equation of state, one has to perform a time-average of these quantities. For instance, we evaluate the total pressure as being equal to
(1)
with Tr ( σ ̂ ( t ) ) = σ xx ( t ) + σ yy ( t ) + σ zz ( t ) ,   σ i i ( t ) the ii component of the stress tensor at time t, t0 and N t the initial time step and the number of time steps of the average window, N tot the total number of atoms, k B the Boltzmann constant, T the temperature, and V the total volume. An example of a time evolution of the instantaneous pressure can be seen in Fig. 1.
FIG. 1.
Pressure evolution of case CH1. (Top) Time evolution of the instantaneous pressure in the molecular dynamics. The dashed line indicates the average value over the time range considered for the statistics. (Bottom) Time evolution of the pressure autocorrelation function. An exponential decay fit is also plotted. The first vertical dashed line shows the characteristic decay time. The second dashed line shows the approximated decorrelation time considered.
View largeDownload slide

Pressure evolution of case CH1. (Top) Time evolution of the instantaneous pressure in the molecular dynamics. The dashed line indicates the average value over the time range considered for the statistics. (Bottom) Time evolution of the pressure autocorrelation function. An exponential decay fit is also plotted. The first vertical dashed line shows the characteristic decay time. The second dashed line shows the approximated decorrelation time considered.

FIG. 1.
Pressure evolution of case CH1. (Top) Time evolution of the instantaneous pressure in the molecular dynamics. The dashed line indicates the average value over the time range considered for the statistics. (Bottom) Time evolution of the pressure autocorrelation function. An exponential decay fit is also plotted. The first vertical dashed line shows the characteristic decay time. The second dashed line shows the approximated decorrelation time considered.
View largeDownload slide

Pressure evolution of case CH1. (Top) Time evolution of the instantaneous pressure in the molecular dynamics. The dashed line indicates the average value over the time range considered for the statistics. (Bottom) Time evolution of the pressure autocorrelation function. An exponential decay fit is also plotted. The first vertical dashed line shows the characteristic decay time. The second dashed line shows the approximated decorrelation time considered.

Close modal

To estimate the statistical uncertainty on the pressure, we use the standard deviation and divide it by the number of uncorrelated steps. The latter is estimated based on the auto-correlation function (ACF) of the pressure fluctuations. We fit the ACF with an exponential decay providing a characteristic decay time. We considered that the pressure was uncorrelated after twice the decay time. An example can be seen in Fig. 1.

We used the trajectories computed during the molecular dynamics simulation to extract the average pair distribution function. For each time step, we compute the distance in between all the particles, including their periodic images, and construct a histogram. By assigning each species appropriately we obtain the pair distribution function using the following expression:
(2)
where α and β refer to two particle types, V is the total volume, t0 and N t are the initial time step and the number of time steps of the average window, N α is the number of particles of type α, Δ r is the bin size in radial distance, r → i α ( t ) is the position vector at time t of the i α th atom of type α, and r → i β * ( t ) is the position vector at time t of the closest periodic image of the i β th atom of type β.

In Fig. 2, we plotted an example of a pair distribution function (PDF) for the equimolar CH mixture, with the C–C, C–H, and H–H respective curves. The PDF does not provide direct information on the transport properties. It does, however, provide details on the structure of the system. Figure 2 shows, for instance, that hydrogen atoms are mostly dissociated, which means their diffusion is mostly that of an atom. The first peak of the C–C PDF would, however, suggest that some carbon–carbon bonding, possibly transient, does exist, which will reduce the magnitude of the carbon interdiffusion.

FIG. 2.
Pair distribution function associated with case CH1.
View largeDownload slide

Pair distribution function associated with case CH1.

FIG. 2.
Pair distribution function associated with case CH1.
View largeDownload slide

Pair distribution function associated with case CH1.

Close modal
To determine the self-diffusion and interdiffusion coefficients, we used two different strategies. The first one relies on the mean squared displacement (MSD). From the fluctuation–dissipation theorem, we can determine that D α⁠, the diffusion coefficient of species α, is given, in three dimensions, as
(3)
with ⟨ ⟩ α being the average over all particles of type α. Numerically speaking, we proceeded in two steps. We first compute the time evolution of the mean squared displacement as follows:
(4)
where the position vector is the physical one and not periodically reintroduced in the simulation box. It is also important to keep the center-of-mass fixed to ensure we compute the diffusion coefficient and not the global movement of the system. As can be seen in Fig. 3, the time evolution is relatively noisy. We can, however, improve the statistics by using the multiple origin technique as follows:47 
(5)
where N t 0 is the number of origins we used. In most cases, we incremented the origin by one time step and used a time window, a third of the total length. To estimate the uncertainty of the MSD, we computed its standard deviation as a function of time. The time origins being correlated and because we merely have a time window of the same duration as the diffusion time, we only divided the standard deviation by N α⁠.
FIG. 3.
Mean squared displacement of aluminum nuclei in the Al1 case. In purple, the MSD was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation and the dashed lines depict an affine fit of the evolutions.
View largeDownload slide

Mean squared displacement of aluminum nuclei in the Al1 case. In purple, the MSD was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation and the dashed lines depict an affine fit of the evolutions.

FIG. 3.
Mean squared displacement of aluminum nuclei in the Al1 case. In purple, the MSD was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation and the dashed lines depict an affine fit of the evolutions.
View largeDownload slide

Mean squared displacement of aluminum nuclei in the Al1 case. In purple, the MSD was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation and the dashed lines depict an affine fit of the evolutions.

Close modal

The second step is to extract the diffusion coefficient from the MSD time evolution. For this purpose, we used a linear regression to an affine curve. We let the algorithm find the most optimal start to the fit to minimize the uncertainty on the diffusion coefficient. We used this criterion to restrict the fit to the diffusion part omitting the ballistic section of the curve. In Fig. 3, we depicted the result of the fit by the dashed lines.

The second strategy resides in using the fluctuation–dissipation theorem in velocity space. We can then link the diffusion coefficient D α to the velocity auto-correlation function as follows:
(6)
where v → ( t ) is the velocity vector at instant t of a given particle with respect to the center-of-mass. Similar to the MSD, we compute the velocity auto-correlation function (VACF) γ v , α using the multi-origin technique as follows:
(7)
Similar to the MSD, we estimate the uncertainty on the VACF by computing the standard deviation on γ v , α and divided it by N t 0 ¯⁠, where we estimated the number of uncorrelated origins N t 0 ¯ considering that they were uncorrelated after two decaying times of the VACF. To compute the diffusion coefficient, we simply integrate the VACF over the time windows and average over the last few time steps. The uncertainty on the integral is difficult to estimate as the uncertainties over the VACF are correlated. We, thus, simply integrated the uncertainty on the VACF.

Integrated VACF and MSD give very consistent results as previously shown in Ref. 48. For instance, the cases plotted in Figs. 3 and 4 give a diffusion coefficient of 9.5 × 10 − 8 m2/s for the MSD and 9.9 × 10 − 8 m2/s for the integrated VACF. It is also interesting to note that we computed these quantities within an isokinetic ensemble which is not considered as the most well-suited ensemble for such calculations because it rescales the velocity components at each time step. It is known that the choice of thermostat is essential for the quality of the ionic transport properties estimates.49 Nevertheless, based on the comparison displayed in Ref. 28, it is clear that the isokinetic ensemble with our methodology gives very satisfactory results. It is possible that the cell size, the short chosen time steps, and the multi-origin technique help compensate for the crude thermalization associated with the isokinetic ensemble.

FIG. 4.
Normalized (top) and integrated (bottom) velocity autocorrelation function of aluminum nuclei in the Al1 case. In purple, the VACF was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation.
View largeDownload slide

Normalized (top) and integrated (bottom) velocity autocorrelation function of aluminum nuclei in the Al1 case. In purple, the VACF was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation.

FIG. 4.
Normalized (top) and integrated (bottom) velocity autocorrelation function of aluminum nuclei in the Al1 case. In purple, the VACF was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation.
View largeDownload slide

Normalized (top) and integrated (bottom) velocity autocorrelation function of aluminum nuclei in the Al1 case. In purple, the VACF was computed with a single origin, the blue case is with multiple origins. The shaded areas correspond to the one-σ deviation.

Close modal
Another important quantity for hydrodynamic simulations is the viscosity, which we can also determine from molecular dynamics. The fluctuation–dissipation theorem also links the viscosity with the stress auto-correlation function (SACF) as follows:
(8)
where s i j ( t ) is one of the off diagonal components of the stress tensor (including the kinetic contribution) taken at instant t, V is the total volume, and T is the temperature. Numerically, we compute the SACF in a very similar manner to the VACF using the multi-origin technique and averaging the results of the three off diagonal components. We, thus, compute the SACF γs using the following expression:
(9)
where we denoted δ s i j , t 0 ( t ) as the value of the ij component of the stress tensor at instant t when we remove its average value over the time window starting at t0. Doing so allows to remove any slight residue in the off diagonal stress tensor component which would lead to a diverging integral. The uncertainty on the SACF is computed as the standard deviation divided by 3 N t 0 ¯⁠, where we estimated the number of uncorrelated origins N t 0 ¯ considering that they were uncorrelated after two decaying times of the SACF.
The SACF tends to be much more noisy than the VACF. This is because the SACF is derived from an intensive property, namely, the stress tensor, which does not scale with the number of particles in the system, in contrast with the MSD and VACF. We, thus, opted for a fit of the SACF to offer a less noisy estimate of the viscosity as mentioned in Ref. 50. For this purpose, we approximated the SACF by16 
(10)
where A, B, τ1, τ2, and ω are adjusted on the SACF.
After integration, we obtain the viscosity from the fit parameters as follows:
(11)
We use an error propagation method to determine the uncertainty on the viscosity. An example of SACF and its integral as well as the corresponding fit are plotted in Fig. 5.
FIG. 5.
Normalized (top) and integrated (bottom) stress tensor autocorrelation function of gold nuclei in the Au1 case. The full line and shaded area correspond to the multiple origin calculation of the average SACF. The dashed line is the result of the fit presented in Eq. (10).
View largeDownload slide

Normalized (top) and integrated (bottom) stress tensor autocorrelation function of gold nuclei in the Au1 case. The full line and shaded area correspond to the multiple origin calculation of the average SACF. The dashed line is the result of the fit presented in Eq. (10).

FIG. 5.
Normalized (top) and integrated (bottom) stress tensor autocorrelation function of gold nuclei in the Au1 case. The full line and shaded area correspond to the multiple origin calculation of the average SACF. The dashed line is the result of the fit presented in Eq. (10).
View largeDownload slide

Normalized (top) and integrated (bottom) stress tensor autocorrelation function of gold nuclei in the Au1 case. The full line and shaded area correspond to the multiple origin calculation of the average SACF. The dashed line is the result of the fit presented in Eq. (10).

Close modal

Density functional theory within the Kohn–Sham–Mermin framework gives access to the electronic structure of the complete system. It, however, requires a higher level of convergence regarding the electronic density to offer reliable results. This is why we performed new DFT calculations on a series of snapshots extracted from the initial molecular dynamics. We used a higher number of bands to improve the convergence on the density of states and on the transition matrix elements of use in the Kubo–Greenwood formula (see below). For the MD calculations performed within the extended DFT framework,23 the snapshots were first converged with the extended approach then with the standard Kohn–Sham orbitals. The detailed parameters of the simulations are available in Table II. For all the cases reported in Ref. 28, we restricted ourselves to non spin-polarized calculations.

TABLE II.

KS-DFT snapshots of input parameters and electronic transport property analysis. N snapshot stands for the number of snapshots selected in the molecular dynamics trajectory, at least separated by the decorrelation time (Fig. 1). Δ t is the time between each snapshot, where Δ t is necessarily greater than the decorrelation time.

H1 C1 CH1 Al1 Cu1 CH2 Au1 C2
N snapshot  15  6  6  6  6  6  6  6 
Δ t (ps)  0.121  0.193  0.242  1.209  0.605  0.0435  0.145  0.0726 
N atom  256 H  64 C  48 H/48 C  108 Al  108 Cu  48 H/48 C  108 Au  108 C 
k →-point(s)  23  33  33  MP  43  33  23  83 
N band KS  512  4096  2048  8000  4096  4096  12 288  1024 
E cut PW (Ha)  62.0  320.0  60.0  50.0  40.0  60.0  50.0  400.0 
H1 C1 CH1 Al1 Cu1 CH2 Au1 C2
N snapshot  15  6  6  6  6  6  6  6 
Δ t (ps)  0.121  0.193  0.242  1.209  0.605  0.0435  0.145  0.0726 
N atom  256 H  64 C  48 H/48 C  108 Al  108 Cu  48 H/48 C  108 Au  108 C 
k →-point(s)  23  33  33  MP  43  33  23  83 
N band KS  512  4096  2048  8000  4096  4096  12 288  1024 
E cut PW (Ha)  62.0  320.0  60.0  50.0  40.0  60.0  50.0  400.0 
View Large
The density of states (DOS) was computed for each snapshot and averaged after aligning the electronic chemical potential to zero. To transform the discrete eigenvalues, we used normalized Gaussian functions with a half-height width ς of 0.5 eV typically. The averaged density of states D is then given as
(12)
where N t is the number of snapshots, w k → is the relative weight of the k →-points, N b is the number of Kohn–Sham orbitals, ϵ i , k → ( t ) is the ith eigenstate at k →-point k → and instant t, and μ ( t ) is the electronic chemical potential at instant t. The occupied DOS is determined from the full DOS weighted by a Fermi–Dirac factor. An example of the resulting DOS is plotted in Fig. 6. In this figure, we can clearly see the 2s and 2p peaks while the 3s states are already merged with the continuum.
FIG. 6.
Electronic density of states of aluminum for case Al1 ( ρ = 2.7 g/cm3, k B T = 1 eV). The zero of energy corresponds to the chemical potential. Low energy spikes correspond to 2s and 2p electron core shells, while higher energy shells are merged in the continuum, starting at ε ≈ − 12.5 eV.
View largeDownload slide

Electronic density of states of aluminum for case Al1 (⁠ ρ = 2.7 g/cm3, k B T = 1 eV). The zero of energy corresponds to the chemical potential. Low energy spikes correspond to 2s and 2p electron core shells, while higher energy shells are merged in the continuum, starting at ε ≈ − 12.5 eV.

FIG. 6.
Electronic density of states of aluminum for case Al1 ( ρ = 2.7 g/cm3, k B T = 1 eV). The zero of energy corresponds to the chemical potential. Low energy spikes correspond to 2s and 2p electron core shells, while higher energy shells are merged in the continuum, starting at ε ≈ − 12.5 eV.
View largeDownload slide

Electronic density of states of aluminum for case Al1 (⁠ ρ = 2.7 g/cm3, k B T = 1 eV). The zero of energy corresponds to the chemical potential. Low energy spikes correspond to 2s and 2p electron core shells, while higher energy shells are merged in the continuum, starting at ε ≈ − 12.5 eV.

Close modal
From the detailed Kohn–Sham electronic structure, it is possible to determine the Onsager coefficients via a Kubo–Greenwood formula as implemented in Abinit code.42 We have51,
(13)
where e is the electron charge, m e is the electron mass, ω is the pulsation, f is the Fermi–Dirac function, h e is the electronic enthalpy, and ⟨ Ψ j , k → | ∇ α | Ψ i , k → ⟩ are the transition matrix elements as computed from the KS orbitals. From this relationship, we obtain the real part of the electrical conductivity as follows:51,52
(14)
and the electronic contribution to the thermal conductivity,51 
(15)
These quantities are computed for each snapshot and then averaged over the different snapshots to obtain an average quantity. Some large fluctuations can be observed for some systems, especially for those that are not fully metallic. We estimate the uncertainty based on the standard deviation and the number of snapshots, which we consider as uncorrelated. Examples of electrical and thermal conductivities are plotted in Figs. 7 and 8, respectively.
FIG. 7.
Electrical conductivity as a function of the excitation energy for aluminum (top, case Al1), copper (middle, case Cu1), and CH (bottom, case CH2). The full line and the shaded areas are the conductivity and one-sigma deviation as obtained from DFT. We fitted the results with a Drude fit [Eq. (16)], an exponent of the ω τ term [Eq. (17) labeled “Alpha-fit”], and a Drude–Smith expression [Eq. (18)].
View largeDownload slide

Electrical conductivity as a function of the excitation energy for aluminum (top, case Al1), copper (middle, case Cu1), and CH (bottom, case CH2). The full line and the shaded areas are the conductivity and one-sigma deviation as obtained from DFT. We fitted the results with a Drude fit [Eq. (16)], an exponent of the ω τ term [Eq. (17) labeled “Alpha-fit”], and a Drude–Smith expression [Eq. (18)].

FIG. 7.
Electrical conductivity as a function of the excitation energy for aluminum (top, case Al1), copper (middle, case Cu1), and CH (bottom, case CH2). The full line and the shaded areas are the conductivity and one-sigma deviation as obtained from DFT. We fitted the results with a Drude fit [Eq. (16)], an exponent of the ω τ term [Eq. (17) labeled “Alpha-fit”], and a Drude–Smith expression [Eq. (18)].
View largeDownload slide

Electrical conductivity as a function of the excitation energy for aluminum (top, case Al1), copper (middle, case Cu1), and CH (bottom, case CH2). The full line and the shaded areas are the conductivity and one-sigma deviation as obtained from DFT. We fitted the results with a Drude fit [Eq. (16)], an exponent of the ω τ term [Eq. (17) labeled “Alpha-fit”], and a Drude–Smith expression [Eq. (18)].

Close modal
FIG. 8.
Thermal conductivity of CH (case CH2) at 0.9 g/cm3 and 7.8 eV as a function of the excitation energy. The DC value obtained in this case is within the error bar of a recent experimental value of 465 ± 233 W/m/K.55
View largeDownload slide

Thermal conductivity of CH (case CH2) at 0.9 g/cm3 and 7.8 eV as a function of the excitation energy. The DC value obtained in this case is within the error bar of a recent experimental value of 465 ± 233 W/m/K.55 

FIG. 8.
Thermal conductivity of CH (case CH2) at 0.9 g/cm3 and 7.8 eV as a function of the excitation energy. The DC value obtained in this case is within the error bar of a recent experimental value of 465 ± 233 W/m/K.55
View largeDownload slide

Thermal conductivity of CH (case CH2) at 0.9 g/cm3 and 7.8 eV as a function of the excitation energy. The DC value obtained in this case is within the error bar of a recent experimental value of 465 ± 233 W/m/K.55 

Close modal

For the purpose of magneto-hydrodynamic simulations, one is interested in the low frequency or even the direct current value of the electrical and thermal conductivities. Due to the finite-size of the simulation cell, the energy spectrum of the Kohn–Sham orbitals is discretized and has a finite resolution, which prevents a direct calculation of the DC conductivities. Despite recent efforts made into computing more accurately the intra-band contributions,53 it is necessary to extrapolate the conductivity values toward the lowest frequencies.

It is possible to use physically informed extrapolation methods, such as a Drude model, for conducting systems. In this case, we fit a Drude-like model of the form
(16)
on the low-frequency range of the conductivity, typically from 0.5 to 3 eV, and adjust the σ0 and τ parameters. This expression, however, fails to reproduce the conductivity behavior of transition metals, such as copper, for which the d-electrons create a large spike in conductivity at low energy. For such systems, it is more realistic to use an exponent of the ω τ term as follows:
(17)
For some of the cases explored, especially the plastics, it is also possible to use a Drude–Smith functional form,54 
(18)
Figure 7 presents three examples of conductivity profiles and the fit with different models. It is to be noted that we restrict the fit to the low photon energies since we are interested in obtaining the direct current (DC) conductivities. We adapted the fitted region to the type of functional form of the fit in order to minimize the residuals. We typically used the conductivity values over the range from 0.5 to 3–5 eV.

For the thermal conductivity, it is possible to use a fit similar to Eq. (17). However, because it includes higher order Onsager coefficients, the thermal conductivity tends to be more continuous near the origin. Figure 8 shows an example for the equimolar CH mixture at a high temperature, around 7.8 eV, and a density of 0.9 g/cm3. By simply averaging the DC value obtained for each snapshot, we obtain a thermal conductivity of 650 ± 10 W/m/K. This value compares very well with the recent experimental determination of 465 ± 233 W/m/K,55 which shows the capabilities of the KS-DFT-MD method in computing accurate transport coefficients.

From the Onsager coefficients, it is also possible to extract optical properties. In particular, from the frequency-dependent electrical conductivity, one can compute the absorption coefficient α. A first step is to compute the imaginary part of the conductivity using a Kramers–Kronig relationship and to then determine the optical index and its related quantities (reflectivity and absorption coefficients).52 These quantities were not explicitly requested for the workshop but they are linked to the opacities, which are related to radiative transport. It is also important to note that during shock experiments it is very often an optical quantity that is measured and we can use ab initio results to relate the optical quantities to the DC conductivity.56 An example of the absorption coefficient is plotted in Fig. 9. It corresponds to the Al1 case. As a comparison, we also added the result obtained at k B T = 1 eV and ρ = 0.27 g/cm3. We can see that a decrease in the density by a factor of 10 results in a decrease in absorption by a factor of nearly 10 as well. One could also note the pronounced peak in the absorption near 5 eV, which is mainly due to transitions from the (semi-)bound states near the continuum to the continuum, the 3s states being partly recombined at low density.

FIG. 9.
Frequency-dependent absorption coefficient as a function of the incident photon energy for Al at different densities and k B T = 1 eV. The result was averaged over six snapshots in the two cases.
View largeDownload slide

Frequency-dependent absorption coefficient as a function of the incident photon energy for Al at different densities and k B T = 1 eV. The result was averaged over six snapshots in the two cases.

FIG. 9.
Frequency-dependent absorption coefficient as a function of the incident photon energy for Al at different densities and k B T = 1 eV. The result was averaged over six snapshots in the two cases.
View largeDownload slide

Frequency-dependent absorption coefficient as a function of the incident photon energy for Al at different densities and k B T = 1 eV. The result was averaged over six snapshots in the two cases.

Close modal

To summarize, KS-DFT molecular dynamics offers a unified approach to compute ionic and electronic transport properties as well as standard thermodynamic quantities. However, the final determination of these properties requires some additional assumptions and techniques to circumvent difficulties encountered when using this theoretical scheme. We bypassed the orbital wall encountered at high temperatures using the extended DFT method, allowing for multi-eV molecular dynamics simulations with comparable accuracy as conventional KS-DFT. Additionally, to bypass the lack of statistics when determining diffusion coefficients, and viscosity, we used the multiple origin technique. While thermodynamic quantities are not so sensible to finite size effects, this is not the case for electronic transport properties such as direct current and thermal conductivities. We used physically informed extrapolation methods to determine if the finite size effects are negligible or not. These different techniques coupled with the high accuracy of the KS-DFT molecular dynamics method showed accurate results, which we were able to compare to recent experiments and other codes via the second charged-particle transport code comparison workshop.

We are very thankful to Robin Piron for fruitful discussions and to the three anonymous referees for their insightful comments.

The authors have no conflicts to disclose.

Augustin Blanchet: Conceptualization (equal); Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Validation (equal); Visualization (equal); Writing – original draft (supporting); Writing – review & editing (equal). Vanina Recoules: Conceptualization (equal); Data curation (equal); Formal analysis (equal); Methodology (supporting); Validation (equal); Visualization (equal); Writing – original draft (supporting); Writing – review & editing (equal). François Soubiran: Conceptualization (equal); Data curation (equal); Formal analysis (lead); Investigation (equal); Methodology (lead); Validation (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead). Mikael Tacu: Conceptualization (equal); Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Validation (equal); Visualization (equal); Writing – original draft (supporting); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.

1.
S.
Liberatore
,
P.
Gauthier
,
J.-L.
Willien
,
P.-E.
Masson-Laborde
,
F.
Philippe
,
O.
Poujade
,
E.
Alozy
,
R.
Botrel
,
G.
Boutoux
,
T.
Caillaud
,
C.
Chicanne
,
C.
Chollet
,
A.
Debayle
,
S.
Depierreux
,
W.
Duchastenier
,
M.
Ferri
,
O.
Henry
,
P.
Hoch
,
S.
Laffite
,
O.
Landoas
,
L.
Le-Deroff
,
E.
Lefebvre
,
G.
Legay
,
I.
Marmajou
,
C.
Meyer
,
K.
Molina
,
O.
Morice
,
E.
Peche
,
P.
Prunet
,
R.
Riquier
,
R.
Rosch
,
V.
Tassin
,
X.
Vaisseau
, and
B.
Villette
, “
First indirect drive inertial confinement fusion campaign at laser megajoule
,”
Phys. Plasmas
30
,
122707
(
2023
).
https://doi.org/10.1063/5.0176446
Google Scholar
Crossref
Search ADS
 
2.
N.
Schaeffer
, “
Efficient spherical harmonic transforms aimed at pseudo-spectral numerical simulations
,”
Geochem. Geophys. Geosyst.
14
,
751
(
2013
).
https://doi.org/10.1002/ggge.20071
Google Scholar
Crossref
Search ADS
 
3.
B.
Paxton
,
L.
Bildsten
,
A.
Dotter
,
F.
Herwig
,
P.
Lesaffre
, and
F.
Timmes
, “
Modules for experiments in stellar astrophysics (MESA)
,”
Astrophys. J.
192
,
3
(
2011
).
https://doi.org/10.1088/0067-0049/192/1/3
Google Scholar
Crossref
Search ADS
 
4.
R.
Devriendt
 and
O.
Poujade
, “
Classical molecular dynamic simulations and modeling of inverse bremsstrahlung heating in low Z weakly coupled plasmas
,”
Phys. Plasmas
29
,
073301
(
2022
).
https://doi.org/10.1063/5.0091662
Google Scholar
Crossref
Search ADS
 
5.
B.
Hammel
,
S.
Haan
,
D.
Clark
,
M.
Edwards
,
S.
Langer
,
M.
Marinak
,
M.
Patel
,
J.
Salmonson
, and
H.
Scott
, “
High-mode Rayleigh-Taylor growth in NIF ignition capsules
,”
High Energy Density Phys.
6
,
171
–
178
(
2010
).
https://doi.org/10.1016/j.hedp.2009.12.005
Google Scholar
Crossref
Search ADS
 
6.
F.
Soubiran
 and
B.
Militzer
, “
Electrical conductivity and magnetic dynamos in magma oceans of super-earths
,”
Nat. Commun.
9
,
3883
(
2018
).
https://doi.org/10.1038/s41467-018-06432-6
Google Scholar
Crossref
Search ADS
PubMed
 
7.
J.
Leconte
 and
G.
Chabrier
, “
A new vision of giant planet interiors: Impact of double diffusive convection
,”
Astron. Astrophys.
540
,
A20
(
2012
).
https://doi.org/10.1051/0004-6361/201117595
Google Scholar
Crossref
Search ADS
 
8.
B.
Militzer
,
F.
Soubiran
,
S. M.
Wahl
, and
W. B.
Hubbard
, “
Understanding Jupiter's interior
,”
J. Geophys. Res.
121
,
1552
, https://doi.org/10.1002/2016JE005080 (
2016
).
Google Scholar
Crossref
Search ADS
 
9.
N.
Nettelmann
,
K.
Wang
,
J. J.
Fortney
,
S.
Hamel
,
S.
Yellamilli
,
M.
Bethkenhagen
, and
R.
Redmer
, “
Uranus evolution models with simple thermal boundary layers
,”
Icarus
275
,
107
(
2016
).
https://doi.org/10.1016/j.icarus.2016.04.008
Google Scholar
Crossref
Search ADS
 
10.
B.
Freytag
,
H.-G.
Ludwig
, and
M.
Steffen
, “
Hydrodynamical models of stellar convection. The role of overshoot in DA white dwarfs, a-type stars, and the sun
,”
Astron. Astrophys.
313
,
497
(
1996
).
Google Scholar
 
11.
L.
Spitzer
 and
R.
Härm
, “
Transport phenomena in a completely ionized gas
,”
Phys. Rev.
89
,
977
(
1953
).
https://doi.org/10.1103/PhysRev.89.977
Google Scholar
Crossref
Search ADS
 
12.
S.
Chapman
 and
T. G.
Cowling
,
The Mathematical Theory of Non-Uniform Gases
(
Cambridge University Press
,
1970
).
Google Scholar
 
13.
P.
Drude
, “
Zur elektrontheorie der metalle
,”
Ann. Phys.
306
,
566
(
1900
).
https://doi.org/10.1002/andp.19003060312
Google Scholar
Crossref
Search ADS
 
14.
G.
Faussurier
,
C.
Blancard
,
P.
Renaudin
, and
P. L.
Silvestreli
, “
Electrical conductivity of warm expanded Al
,”
Phys. Rev. B
73
,
075106
(
2006
).
https://doi.org/10.1103/PhysRevB.73.075106
Google Scholar
Crossref
Search ADS
 
15.
C.
Ticknor
,
J. D.
Kress
,
L. A.
Collins
,
J.
Clérouin
,
P.
Arnault
, and
A.
Decoster
, “
Transport properties of an asymmetric mixture in the dense plasma regime
,”
Phys. Rev. E
93
,
063208
(
2016
).
https://doi.org/10.1103/PhysRevE.93.063208
Google Scholar
Crossref
Search ADS
PubMed
 
16.
F.
Soubiran
,
B.
Militzer
,
K. P.
Driver
, and
S.
Zhang
, “
Properties of hydrogen, helium, and silicon dioxide mixtures in giant planet interiors
,”
Phys. Plasmas
24
,
041401
(
2017
).
https://doi.org/10.1063/1.4978618
Google Scholar
Crossref
Search ADS
 
17.
D.
Soyuer
,
F.
Soubiran
, and
R.
Helled
, “
Constraining the depth of the winds on Uranus and Neptune via ohmic dissipation
,”
Mon. Not. R. Astron. Soc.
498
,
621
–
638
(
2020
).
https://doi.org/10.1093/mnras/staa2461
Google Scholar
Crossref
Search ADS
 
18.
M.
French
,
A.
Becker
,
W.
Lorenzen
,
N.
Nettelmann
,
M.
Bethkenhagen
,
J.
Wicht
, and
R.
Redmer
, “
Ab initio simulations for material properties along the Jupiter adiabat
,”
Astrophys. J. Suppl. Ser.
202
,
5
(
2012
).
https://doi.org/10.1088/0067-0049/202/1/5
Google Scholar
Crossref
Search ADS
 
19.
M.
Preising
,
M.
French
,
C.
Mankovich
,
F.
Soubiran
, and
R.
Redmer
, “
Material properties of Saturn's interior from ab initio simulations
,”
Astrophys. J. Suppl. Ser.
269
,
47
(
2023
).
https://doi.org/10.3847/1538-4365/ad0293
Google Scholar
Crossref
Search ADS
 
20.
A.
Blanchet
,
M.
Torrent
, and
J.
Clérouin
, “
Requirements for very high temperature Kohn–Sham DFT simulations and how to bypass them
,”
Phys. Plasmas
27
,
122706
(
2020
).
https://doi.org/10.1063/5.0016538
Google Scholar
Crossref
Search ADS
 
21.
V.
Sharma
,
L. A.
Collins
, and
A. J.
White
, “
Stochastic and mixed density functional theory within the projector augmented wave formalism for simulation of warm dense matter
,”
Phys. Rev. E
108
,
L023201
(
2023
).
https://doi.org/10.1103/PhysRevE.108.L023201
Google Scholar
Crossref
Search ADS
PubMed
 
22.
S.
Zhang
,
H.
Wang
,
W.
Kang
,
P.
Zhang
, and
X. T.
He
, “
Extended application of Kohn–Sham first-principles molecular dynamics method with plane wave approximation at high energy-from cold materials to hot dense plasmas
,”
Phys. Plasmas
23
,
042707
(
2016
).
https://doi.org/10.1063/1.4947212
Google Scholar
Crossref
Search ADS
 
23.
A.
Blanchet
,
J.
Clérouin
,
M.
Torrent
, and
F.
Soubiran
, “
Extended first-principles molecular dynamics model for high temperature simulations in the Abinit code: Application to warm dense aluminum
,”
Comput. Phys. Commun.
271
,
108215
(
2022
).
https://doi.org/10.1016/j.cpc.2021.108215
Google Scholar
Crossref
Search ADS
 
24.
P.
Hollebon
 and
T.
Sjostrom
, “
Hybrid Kohn-Sham + Thomas-Fermi scheme for high-temperature density functional theory
,”
Phys. Rev. B
105
,
235114
(
2022
).
https://doi.org/10.1103/PhysRevB.105.235114
Google Scholar
Crossref
Search ADS
 
25.
X.
Gonze
,
B.
Amadon
,
G.
Antonius
,
F.
Arnardi
,
L.
Baguet
,
J.-M.
Beuken
,
J.
Bieder
,
F.
Bottin
,
J.
Bouchet
,
E.
Bousquet
,
N.
Brouwer
,
F.
Bruneval
,
G.
Brunin
,
T.
Cavignac
,
J.-B.
Charraud
,
W.
Chen
,
M.
Côté
,
S.
Cottenier
,
J.
Denier
,
G.
Geneste
,
P.
Ghosez
,
M.
Giantomassi
,
Y.
Gillet
,
O.
Gingras
,
D. R.
Hamann
,
G.
Hautier
,
X.
He
,
N.
Helbig
,
N.
Holzwarth
,
Y.
Jia
,
F.
Jollet
,
W.
Lafargue-Dit-Hauret
,
K.
Lejaeghere
,
M. A. L.
Marques
,
A.
Martin
,
C.
Martins
,
H. P. C.
Miranda
,
F.
Naccarato
,
K.
Persson
,
G.
Petretto
,
V.
Planes
,
Y.
Pouillon
,
S.
Prokhorenko
,
F.
Ricci
,
G.-M.
Rignanese
,
A. H.
Romero
,
M. M.
Schmitt
,
M.
Torrent
,
M. J.
van Setten
,
B. V.
Troeye
,
M. J.
Verstraete
,
G.
Zérah
, and
J. W.
Zwanziger
, “
The Abinit project: Impact, environment and recent developments
,”
Comput. Phys. Commun.
248
,
107042
(
2020
).
https://doi.org/10.1016/j.cpc.2019.107042
Google Scholar
Crossref
Search ADS
 
26.
F.
Bottin
,
S.
Leroux
,
A.
Knyazev
, and
G.
Zérah
, “
Large-scale ab initio calculations based on three levels of parallelization
,”
Comput. Mater. Sci.
42
,
329
–
336
(
2008
).
https://doi.org/10.1016/j.commatsci.2007.07.019
Google Scholar
Crossref
Search ADS
 
27.
M.
Torrent
,
N.
Holzwarth
,
F.
Jollet
,
D.
Harris
,
N.
Lepley
, and
X.
Xu
, “
Electronic structure packages: Two implementations of the projector augmented wave (PAW) formalism
,”
Comput. Phys. Commun.
181
,
1862
–
1867
(
2010
).
https://doi.org/10.1016/j.cpc.2010.07.036
Google Scholar
Crossref
Search ADS
 
28.
L. J.
Stanek
,
A.
Kononov
,
S. B.
Hansen
,
B. M.
Haines
,
S.
Hu
,
P. F.
Knapp
,
M. S.
Murillo
,
L. G.
Stanton
,
H. D.
Whitley
,
S. D.
Baalrud
 et al, “
Review of the second charged-particle transport coefficient code comparison workshop
,”
Phys. Plasmas
31
,
052104
(
2024
).
https://doi.org/10.1063/5.0198155
Google Scholar
Crossref
Search ADS
 
29.
P.
Minary
,
G. J.
Martyna
, and
M. E.
Tuckerman
, “
Algorithms and novel applications based on the isokinetic ensemble. I. Biophysical and path integral molecular dynamics
,”
J. Chem. Phys.
118
,
2510
–
2526
(
2003
).
https://doi.org/10.1063/1.1534582
Google Scholar
Crossref
Search ADS
 
30.
F. M. C.
Soubiran
 and
B.
Militzer
, “
Anharmonicity and phase diagram of magnesium oxide in the megabar regime
,”
Phys. Rev. Lett.
125
,
175701
(
2020
).
https://doi.org/10.1103/PhysRevLett.125.175701
Google Scholar
Crossref
Search ADS
PubMed
 
31.
L. J.
Stanek
,
R. C.
Clay
,
M.
Dharma-Wardana
,
M. A.
Wood
,
K. R.
Beckwith
, and
M. S.
Murillo
, “
Efficacy of the radial pair potential approximation for molecular dynamics simulations of dense plasmas
,”
Phys. Plasmas
28
,
032706
(
2021
).
https://doi.org/10.1063/5.0040062
Google Scholar
Crossref
Search ADS
 
32.
I.-C.
Yeh
 and
G.
Hummer
, “
System-size dependence of diffusion coefficients and viscosities from molecular dynamics simulations with periodic boundary conditions
,”
J. Phys. Chem. B
108
,
15873
–
15879
(
2004
).
https://doi.org/10.1021/jp0477147
Google Scholar
Crossref
Search ADS
 
33.
W.
Kohn
 and
L. J.
Sham
, “
Self-consistent equations including exchange and correlation effects
,”
Phys. Rev.
140
,
A1133
–
A1138
(
1965
).
https://doi.org/10.1103/PhysRev.140.A1133
Google Scholar
Crossref
Search ADS
 
34.
N. D.
Mermin
, “
Thermal properties of the inhomogeneous electron gas
,”
Phys. Rev.
137
,
A1441
–
A1443
(
1965
).
https://doi.org/10.1103/PhysRev.137.A1441
Google Scholar
Crossref
Search ADS
 
35.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
, “
Generalized gradient approximation made simple
,”
Phys. Rev. Lett.
77
,
3865
–
3868
(
1996
).
https://doi.org/10.1103/PhysRevLett.77.3865
Google Scholar
Crossref
Search ADS
PubMed
 
36.
V. V.
Karasiev
,
T.
Sjostrom
,
J.
Dufty
, and
S.
Trickey
, “
Accurate homogeneous electron gas exchange-correlation free energy for local spin-density calculations
,”
Phys. Rev. Lett.
112
,
076403
(
2014
).
https://doi.org/10.1103/PhysRevLett.112.076403
Google Scholar
Crossref
Search ADS
PubMed
 
37.
S.
Groth
,
T.
Dornheim
,
T.
Sjostrom
,
F. D.
Malone
,
W.
Foulkes
, and
M.
Bonitz
, “
Ab initio exchange-correlation free energy of the uniform electron gas at warm dense matter conditions
,”
Phys. Rev. Lett.
119
,
135001
(
2017
).
https://doi.org/10.1103/PhysRevLett.119.135001
Google Scholar
Crossref
Search ADS
PubMed
 
38.
V. V.
Karasiev
,
J. W.
Dufty
, and
S.
Trickey
, “
Nonempirical semilocal free-energy density functional for matter under extreme conditions
,”
Phys. Rev. Lett.
120
,
076401
(
2018
).
https://doi.org/10.1103/PhysRevLett.120.076401
Google Scholar
Crossref
Search ADS
PubMed
 
39.
V. V.
Karasiev
,
D. I.
Mihaylov
, and
S. X.
Hu
, “
Meta-GGA exchange-correlation free energy density functional to increase the accuracy of warm dense matter simulations
,”
Phys. Rev. B
105
,
L081109
(
2022
).
https://doi.org/10.1103/PhysRevB.105.L081109
Google Scholar
Crossref
Search ADS
 
40.
P. E.
Blöchl
,
C. G.
Van de Walle
, and
S. T.
Pantelides
, “
First-principles calculations of diffusion coefficients: Hydrogen in silicon
,”
Phys. Rev. Lett.
64
,
1401
–
1404
(
1990
).
https://doi.org/10.1103/PhysRevLett.64.1401
Google Scholar
Crossref
Search ADS
PubMed
 
41.
N.
Holzwarth
,
A.
Tackett
, and
G.
Matthews
, “
A projector augmented wave (PAW) code for electronic structure calculations, Part i: Atompaw for generating atom-centered functions
,”
Comput. Phys. Commun.
135
,
329
–
347
(
2001
).
https://doi.org/10.1016/S0010-4655(00)00244-7
Google Scholar
Crossref
Search ADS
 
42.
N.
Brouwer
,
V.
Recoules
,
N.
Holzwarth
, and
M.
Torrent
, “
Calculation of optical properties with spin-orbit coupling for warm dense matter
,”
Comput. Phys. Commun.
266
,
108029
(
2021
).
https://doi.org/10.1016/j.cpc.2021.108029
Google Scholar
Crossref
Search ADS
 
43.
N.
Jourdain
,
V.
Recoules
,
L.
Lecherbourg
,
P.
Renaudin
, and
F.
Dorchies
, “
Understanding xanes spectra of two-temperature warm dense copper using ab initio simulations
,”
Phys. Rev. B
101
,
125127
(
2020
).
https://doi.org/10.1103/PhysRevB.101.125127
Google Scholar
Crossref
Search ADS
 
44.
G.
Huser
,
V.
Recoules
,
N.
Ozaki
,
T.
Sano
,
Y.
Sakawa
,
G.
Salin
,
B.
Albertazzi
,
K.
Miyanishi
, and
R.
Kodama
, “
Experimental and ab initio investigations of microscopic properties of laser-shocked ge-doped ablator
,”
Phys. Rev. E
92
,
063108
(
2015
).
https://doi.org/10.1103/PhysRevE.92.063108
Google Scholar
Crossref
Search ADS
 
45.
A.
Baldereschi
, “
Mean-value point in the Brillouin zone
,”
Phys. Rev. B
7
,
5212
–
5215
(
1973
).
https://doi.org/10.1103/PhysRevB.7.5212
Google Scholar
Crossref
Search ADS
 
46.
H. J.
Monkhorst
 and
J. D.
Pack
, “
Special points for Brillouin-zone integrations
,”
Phys. Rev. B
13
,
5188
–
5192
(
1976
).
https://doi.org/10.1103/PhysRevB.13.5188
Google Scholar
Crossref
Search ADS
 
47.
D.
Frenkel
 and
B.
Smit
,
Understanding Molecular Simulation: From Algorithms to Applications
(
Academic Press
,
2002
).
Google Scholar
 
48.
C.
Kim
,
O.
Borodin
, and
G. E.
Karniadakis
, “
Quantification of sampling uncertainty for molecular dynamics simulation: Time-dependent diffusion coefficient in simple fluids
,”
J. Comput. Phys.
302
,
485
–
508
(
2015
).
https://doi.org/10.1016/j.jcp.2015.09.021
Google Scholar
Crossref
Search ADS
 
49.
J. E.
Basconi
 and
M. R.
Shirts
, “
Effects of temperature control algorithms on transport properties and kinetics in molecular dynamics simulations
,”
J. Chem. Theory Comput.
9
,
2887
–
2899
(
2013
).
https://doi.org/10.1021/ct400109a
Google Scholar
Crossref
Search ADS
PubMed
 
50.
E. R.
Meyer
,
J. D.
Kress
,
L. A.
Collins
, and
C.
Ticknor
, “
Effect of correlation on viscosity and diffusion in molecular-dynamics simulations
,”
Phys. Rev. E
90
,
043101
(
2014
).
https://doi.org/10.1103/PhysRevE.90.043101
Google Scholar
Crossref
Search ADS
 
51.
B.
Holst
,
M.
French
, and
R.
Redmer
, “
Electronic transport coefficients from ab initio simulations and application to dense liquid hydrogen
,”
Phys. Rev. B
83
,
235120
(
2011
).
https://doi.org/10.1103/PhysRevB.83.235120
Google Scholar
Crossref
Search ADS
 
52.
S.
Mazevet
,
M.
Torrent
,
V.
Recoules
, and
F.
Jollet
, “
Calculations of the transport properties within the paw formalism
,”
High Energy Density Phys.
6
,
84
–
88
(
2010
).
https://doi.org/10.1016/j.hedp.2009.06.004
Google Scholar
Crossref
Search ADS
 
53.
L.
Calderín
,
V. V.
Karasiev
, and
S. B.
Trickey
, “
Kubo-greenwood electrical conductivity formulation and implementation for projector augmented wave datasets
,”
Comput. Phys. Commun.
221
,
118
–
142
(
2017
).
https://doi.org/10.1016/j.cpc.2017.08.008
Google Scholar
Crossref
Search ADS
 
54.
T. L.
Cocker
,
D.
Baillie
,
M.
Buruma
,
L. V.
Titova
,
R. D.
Sydora
,
F.
Marsiglio
, and
F. A.
Hegmann
, “
Microscopic origin of the Drude–Smith model
,”
Phys. Rev. B
96
,
205439
(
2017
).
https://doi.org/10.1103/PhysRevB.96.205439
Google Scholar
Crossref
Search ADS
 
55.
S.
Jiang
,
O. L.
Landen
,
H. D.
Whitley
,
S.
Hamel
,
R.
London
,
D. S.
Clark
,
P.
Sterne
,
S. B.
Hansen
,
S. X.
Hu
,
G. W.
Collins
, and
Y.
Ping
, “
Thermal transport in warm dense matter revealed by refraction-enhanced x-ray radiography with a deep-neural-network analysis
,”
Commun. Phys.
6
,
98
(
2023
).
https://doi.org/10.1038/s42005-023-01190-4
Google Scholar
Crossref
Search ADS
 
56.
M.
Guarguaglini
,
F.
Soubiran
,
J.-A.
Hernandez
,
A.
Benuzzi-Mounaix
,
R.
Bolis
,
E.
Brambrink
,
T.
Vinci
, and
A.
Ravasio
, “
Electrical conductivity of warm dense silica from double-shock experiments
,”
Nat. Commun.
12
,
840
(
2021
).
https://doi.org/10.1038/s41467-021-21046-1
Google Scholar
Crossref
Search ADS
PubMed
 
© 2024 Author(s). All article content, except where otherwise noted, is licensed under a Creative Commons Attribution-NonCommercial-NoDerivs 4.0 International (CC BY-NC-ND) license (https://creativecommons.org/licenses/by-nc-nd/4.0/).