An improved virial expansion for the low-density limit of the electrical conductivity $ \sigma ( T , n )$ of hydrogen as the simplest ionic plasma is presented. Quantum statistical methods provide exact values for the lowest virial coefficients, which serve as a benchmark for analytical approaches to electrical conductivity as well as for numerical results from density functional theory-based molecular dynamics simulations (DFT-MD) or path-integral Monte Carlo simulations. The correction factor introduced by Reinholz *et al.* [Phys. Rev. E **91**, 043105 (2015)] is applied to describe the inclusion of electron–electron collisions in DFT-based calculations of transport coefficients. As a benchmark, the first virial coefficient is correctly described with this approach. The value of the second virial coefficient is discussed, and questions about its value according to DFT-MD simulations are addressed.

## I. INTRODUCTION

*T*and particle number density $ n e = n i = n$ (we consider a charge-neutral hydrogen plasma consisting of electrons

*e*, with mass

*m*, and protons

_{e}*i*, with mass

*m*), the electrons can be bound or degenerate so that a quantum statistical approach is required. Other ionic plasmas with atomic nuclei, charge number

_{i}*Z*, can be treated in a similar way but are not considered here. Instead of the variables

*T*and

*n*, which characterize the state of the plasma, the dimensionless plasma parameter

^{1}using kinetic theory. The Boltzmann equation valid in the low-density limit and the corresponding Fokker–Planck equation were solved, with the result

^{2},

*V*is the volume, $P$ the total momentum of electrons, and the

*t*dependence is according to the Heisenberg picture

*ω*= 0 so that $ \sigma dc ( T , n ) = \sigma ( 0 ; T , n )$. This case is considered in this work; for the general

*ω*-dependent conductivity, see Ref. 3.

*η*is a parameter of absolute convergence, and the limit value $ \eta \u2192 + 0$ must be carried out.

Equilibrium correlation functions can be evaluated using the methods of quantum statistics, such as the Green's function method.^{4,5} This is a perturbative approach in which an expansion with respect to the Coulomb interaction is performed.^{6–8} By applying partial summations describing, e.g., quasiparticles, screening, and the formation of bound states, useful results are obtained for a wide range of plasma parameters. The application to electrical conductivity was also shown in Ref. 9. For the relation to kinetic theory, see Ref. 3. However, exact results are only obtained in special limiting cases. In the following, we discuss the virial expansion of the inverse conductivity of the hydrogen plasma, which contains these exact results. They can be employed as benchmarks for any approach.

A possibility to avoid the perturbation-theoretical treatment of the electron–ion interaction is provided by the density functional theory (DFT), which is already successfully used in condensed matter physics. The application to the calculation of electrical conductivity was made within the framework of the DFT-MD (molecular dynamic) simulations.^{10} The ion dynamics are treated classically using molecular dynamics (MD) simulations, while the electrons are treated as quantum particles, with the wave equations being solved using the DFT formalism. The evaluation of the correlation function (4) is performed numerically using the Kubo–Greenwood formula. However, the treatment of the electron–electron (*e*–*e*) interaction is approximated by a suitable choice of the exchange-correlation functional *E _{xc}*.

^{11}This problem is also a subject of the present work.

The rigorous treatment of the *e*–*e* interaction is possible with the help of path-integral Monte Carlo (PIMC) simulations (see Refs. 12–14 and references therein). Note that the equilibrium correlation function in Eq. (4) is determined by the corresponding spectral function, which can be calculated using the imaginary time (Matsubara) technique.^{4,5} However, the calculations are very complex, and the current shortcomings of this approach include the relatively small number of particles and the sign problem for fermions.

In this work, we focus on two aspects: the virial expansion of the resistivity of hydrogen plasmas in the low-density region and how DFT-MD simulations can be improved to describe the contribution of electron–electron collisions to plasma conductivity. In Sec. II, we discuss the virial expansion of the resistivity. The DFT-MD simulations of plasma conductivity are presented in Sec. III, and ways to account for *e*–*e* collisions are discussed in Sec. IV. Finally, we draw conclusions in Sec. V.

## II. BENCHMARKS FROM THE VIRIAL EXPANSION OF THE RESISTIVITY

### A. The standard virial expansion

For the sake of simplicity, we consider hydrogen plasmas. The results can be transferred to other plasmas with ionic charge *Ze*, since the low-density limit considered here is dominated by the Coulomb force and short-range interactions contribute to higher orders of the virial expansion of the resistivity.

^{3})

^{1}who calculated the value $ \gamma E = 0.5816$, the first virial coefficient

*T*. For the second virial coefficient, the high-temperature limit follows from the quantum Lenard–Balescu (QLB) approach

^{7,15–17}

^{16,18}

*T*are very difficult (see Refs. 16 and 18). In principle, values for $ \rho 2 ( T )$ can be obtained from PIMC simulations as soon as calculations for the conductivity of hydrogen plasmas are possible. Thus, the calculation of the temperature dependence of $ \rho 2 ( T )$ is a challenge for future PIMC simulations.

*x*= 0 is

*ρ*

_{1}and the slope is

*ρ*

_{2}. To perform the extrapolation to

*x*= 0, the variable

*x*must be sufficiently small for linear behavior to be observed.

This method was applied to several approaches for determining the electrical conductivity of hydrogen plasmas, in particular to the data presented in Ref. 20. As shown in Ref. 19, several semi-empirical expressions for plasma conductivity fail to fulfill this benchmark.

^{21}An approximate expression for this third virial term was given in Ref. 6. This higher order virial coefficient can be extracted by analyzing the effective second virial coefficient

*n*is

### B. Generalized virial expansions

^{7}The reference density is $ n ref ( T ) = n / \Gamma 3$. Expressions for the second virial coefficient $ \rho 2 ( T )$ are derived for the classical case in Refs. 16 and 18,

*T*remains open.

The choice of the reference density influences the higher virial coefficients. In principle, $ n \u0302 ref ( T )$ can be chosen such that the next virial coefficient disappears, $ \rho \u0302 2 ( T ) = 0$, as shown at the end of this section [see Eq. (28)]. The logarithmic term in the virial expansion gives us the freedom to choose $ n ref ( T )$ in such a way that the contribution of the higher order virial terms is minimized.

*b*(

*T*,

*n*), Eq. (23), is related to the collision integral evaluated for a screened potential, see Refs. 5, 6, and 15 to regularize the divergent Coulomb potential. The screening parameter $ \kappa 2 = \nu n e 2 / ( \u03f5 0 k B T )$ is determined by the density parameter

*νn*. Various approximations are discussed, static screening by electrons only results in

*ν*= 1, screening by electrons and ions results in

*ν*= 2, for dynamic screening of the electron–ion interaction the value $ \nu = 1.4727$ was derived in Ref. 7. We can scale

*b*(

*T*,

*n*) by a factor

*y*without changing the analytical structure of the virial expansion. It is possible to set the second virial coefficient to zero at $ T \u2192 \u221e$ if the parameter

*b*is replaced by $ b \u0302 = b / y$ with

*y*= 1.5158. Then, we have

*b*by $ b \u0302$. We have $ lim T \u2192 \u221e \rho \u0302 2 ( T ) = 0$, and in the high-temperature limit

*T*-dependent

*y*(

*T*) can be introduced to compensate for the second virial coefficient for all

*T*. These modifications of the virial expansion lead to changes in the higher virial coefficients. Reducing the contributions of the higher order virial coefficients would improve the convergence of the virial expansion and could extend the range of application. Therefore, the modification of the virial expansion is an interesting topic that still needs further investigation.

As example, the values for the conductivity for specific values of *T*, *n* are shown in Table I. The values for low densities and high *T* are considered where the parameter $ x , x \u0303 = 1 / \Lambda $ or $ x \u0302 = 1 / \Lambda \u0302$ is small so that the virial expansion with error bars below few percent is valid. The largest discrepancy occurs for the relatively high density 10 g/cm^{3} and the low temperature 100 eV. The method of virial plots can also be used to extract the virial coefficients $ \rho 1 ( T ) , \rho 2 ( T ) , \rho \u0303 2 ( T )$, or $ \rho \u0302 2 ( T )$ from experimental data or numerical simulations. An example is discussed in Sec. III, see also Ref. 19.

ρ (g cm^{−3})
. | T (eV)
. | Γ . | Θ . | $ 1 / \u2009 ln \u2009 ( \Theta / \Gamma )$ . | σ ( $ \Omega \u2212 1$ cm^{−1})
. | $ \sigma \Lambda $ ( $ \Omega \u2212 1$ cm^{−1})
. | $ 1 / \Lambda $ . | $ \sigma \Lambda \u0302$ ( $ \Omega \u2212 1$ cm^{−1})
. | $ 1 / \Lambda \u0302$ . |
---|---|---|---|---|---|---|---|---|---|

1 | 2 000 | 0.009 777 3 | 77.284 | 0.1114 | 3.5849 $ \xd7 10 6$ | 3.5849 $ \xd7 10 6$ | 0.1003 | 3.5849 $ \xd7 10 6$ | 0.1046 |

1 | 1 000 | 0.019 554 6 | 38.642 | 0.1318 | 1.4825 $ \xd7 10 6$ | 1.4825 $ \xd7 10 6$ | 0.1164 | 1.4825 $ \xd7 10 6$ | 0.1224 |

1 | 700 | 0.027 935 1 | 27.049 | 0.1454 | 951 335 | 951 229 | 0.127 | 951 281 | 0.1341 |

1 | 400 | 0.048 886 5 | 15.457 | 0.1737 | 483 512 | 483 446 | 0.1481 | 483 412 | 0.1578 |

1 | 200 | 0.097 773 | 7.7284 | 0.2288 | 218 811 | 218 659 | 0.1862 | 218 582 | 0.2018 |

1 | 100 | 0.195 546 | 3.8642 | 0.3352 | 107 445 | 107 035 | 0.2504 | 106 826 | 0.2789 |

1 | 70 | 0.279 351 | 2.7049 | 0.4405 | 78 668.2 | 77 910.1 | 0.3035 | 77 531 | 0.3456 |

1.67 | 20 000 | 0.001 16 | 549.05 | 0.076 53 | 7.9375 $ \xd7 10 7$ | 7.9375 $ \xd7 10 7$ | 0.0711 | 7.9375 $ \xd7 10 7$ | 0.073 27 |

1.67 | 10 000 | 0.002 32 | 274.53 | 0.085 61 | 3.1236 $ \xd7 10 7$ | 3.1236 $ \xd7 10 7$ | 0.078 88 | 3.1236 $ \xd7 10 7$ | 0.081 55 |

1.67 | 1 000 | 0.0232 | 27.453 | 0.1413 | 1.5818 $ \xd7 10 6$ | 1.5817 $ \xd7 10 6$ | 0.1239 | 1.5817 $ \xd7 10 6$ | 0.1306 |

1.67 | 100 | 0.232 | 2.7453 | 0.4047 | 125 498 | 124 570 | 0.2865 | 124 103 | 0.324 |

10 | 1 000 | 0.042 129 1 | 8.3252 | 0.1892 | 2.0643 $ \xd7 10 6$ | 2.0638 $ \xd7 10 6$ | 0.1591 | 2.0636 $ \xd7 10 6$ | 0.1704 |

10 | 100 | 0.421 291 | 0.832 52 | 1.468 | 303 434 | 274 960 | 0.5529 | 263 134 | 0.687 |

100 | 1 000 | 0.090 764 4 | 1.7936 | 0.3352 | 3.3977 $ \xd7 10 6$ | 3.3847 $ \xd7 10 6$ | 0.2504 | 3.3781 $ \xd7 10 6$ | 0.2789 |

ρ (g cm^{−3})
. | T (eV)
. | Γ . | Θ . | $ 1 / \u2009 ln \u2009 ( \Theta / \Gamma )$ . | σ ( $ \Omega \u2212 1$ cm^{−1})
. | $ \sigma \Lambda $ ( $ \Omega \u2212 1$ cm^{−1})
. | $ 1 / \Lambda $ . | $ \sigma \Lambda \u0302$ ( $ \Omega \u2212 1$ cm^{−1})
. | $ 1 / \Lambda \u0302$ . |
---|---|---|---|---|---|---|---|---|---|

1 | 2 000 | 0.009 777 3 | 77.284 | 0.1114 | 3.5849 $ \xd7 10 6$ | 3.5849 $ \xd7 10 6$ | 0.1003 | 3.5849 $ \xd7 10 6$ | 0.1046 |

1 | 1 000 | 0.019 554 6 | 38.642 | 0.1318 | 1.4825 $ \xd7 10 6$ | 1.4825 $ \xd7 10 6$ | 0.1164 | 1.4825 $ \xd7 10 6$ | 0.1224 |

1 | 700 | 0.027 935 1 | 27.049 | 0.1454 | 951 335 | 951 229 | 0.127 | 951 281 | 0.1341 |

1 | 400 | 0.048 886 5 | 15.457 | 0.1737 | 483 512 | 483 446 | 0.1481 | 483 412 | 0.1578 |

1 | 200 | 0.097 773 | 7.7284 | 0.2288 | 218 811 | 218 659 | 0.1862 | 218 582 | 0.2018 |

1 | 100 | 0.195 546 | 3.8642 | 0.3352 | 107 445 | 107 035 | 0.2504 | 106 826 | 0.2789 |

1 | 70 | 0.279 351 | 2.7049 | 0.4405 | 78 668.2 | 77 910.1 | 0.3035 | 77 531 | 0.3456 |

1.67 | 20 000 | 0.001 16 | 549.05 | 0.076 53 | 7.9375 $ \xd7 10 7$ | 7.9375 $ \xd7 10 7$ | 0.0711 | 7.9375 $ \xd7 10 7$ | 0.073 27 |

1.67 | 10 000 | 0.002 32 | 274.53 | 0.085 61 | 3.1236 $ \xd7 10 7$ | 3.1236 $ \xd7 10 7$ | 0.078 88 | 3.1236 $ \xd7 10 7$ | 0.081 55 |

1.67 | 1 000 | 0.0232 | 27.453 | 0.1413 | 1.5818 $ \xd7 10 6$ | 1.5817 $ \xd7 10 6$ | 0.1239 | 1.5817 $ \xd7 10 6$ | 0.1306 |

1.67 | 100 | 0.232 | 2.7453 | 0.4047 | 125 498 | 124 570 | 0.2865 | 124 103 | 0.324 |

10 | 1 000 | 0.042 129 1 | 8.3252 | 0.1892 | 2.0643 $ \xd7 10 6$ | 2.0638 $ \xd7 10 6$ | 0.1591 | 2.0636 $ \xd7 10 6$ | 0.1704 |

10 | 100 | 0.421 291 | 0.832 52 | 1.468 | 303 434 | 274 960 | 0.5529 | 263 134 | 0.687 |

100 | 1 000 | 0.090 764 4 | 1.7936 | 0.3352 | 3.3977 $ \xd7 10 6$ | 3.3847 $ \xd7 10 6$ | 0.2504 | 3.3781 $ \xd7 10 6$ | 0.2789 |

The method of virial plots has also been applied to thermodynamic variables describing the equation of state of plasmas (see Ref. 23). There, higher order virial coefficients and generalized virial expansions are discussed. Higher order virial coefficients are of interest to describe degeneracy as known from the Lee–More model (see, e.g., Ref. 22).

## III. THE KUBO–GREENWOOD FORMULA AND DFT-MD SIMULATIONS

*T*and

*n*. In the Born–Oppenheimer approximation, the ions are treated as classical particles whose time-dependent configuration is described by molecular dynamics simulations. They act like an external potential on the electron subsystem. Together with the kinetic energy and Coulomb's

*e*–

*e*interaction, a many-electron Schrödinger equation is formulated, which determines the many-body wave functions. In density functional theory, the correlated many-electron wave function is replaced by the antisymmetrized product of optimized single-electron states, the Kohn–Sham orbitals. For the corresponding Kohn–Sham equations, the Hamiltonian contains the kinetic energy of a non-interacting reference system, the mean-field Coulomb interaction energy, and an exchange-correlation energy

*E*that accounts to some approximation for the effects of

_{xc}*e*–

*e*correlations and antisymmetrization. A simple approximation uses expressions for

*E*that are derived for the uniform electron gas. For frequency-dependent electrical conductivity, see (4) for $ z = \omega + i \eta $, the Kubo–Greenwood formula

_{xc}^{2,24}

^{19,25–29}Kohn–Sham wave functions $ \Psi i , k$ from density functional theory calculations are used to calculate the transition matrix elements of the electron momentum operator $ p \u0302 \alpha , \u2009 \alpha = { x , y , z}$. The Fermi–Dirac distribution $ f ( \u03f5 )$ accounts for the average occupation at energy

*ϵ*, and the summation over momentum space

*k*contains the

*k*-point weights

*w*.

_{k}Due to the finite size of the simulation box, the *δ*-function in Eq. (29) must be approximated by a Gaussian with finite width, which also prevents the direct calculation of the dc conductivity at *ω* = 0. Therefore, the dynamic conductivity is extrapolated to the limit $ \omega \u2192 0$ using a Drude fit. The extrapolation procedure to *ω* = 0 can be improved by using a frequency-dependent collision frequency.^{8}

One of the main shortcomings of the DFT-MD approach is that the many-particle interaction is replaced by a mean-field potential. When using product wave functions for the many-electron system, correlations are excluded. The exchange-correlation energy density functional reflects the Coulomb interaction in a certain approximation, e.g., as it exists in the homogeneous electron gas. However, it becomes problematic in the low-density limit, where correlations, especially *e*–*e* collisions, are important.

DFT-MD simulations have been successfully used to calculate the transport properties in the warm–dense matter region for various materials. The transition to condensed matter, the liquid or solid state, is adequately described, in particular for the degenerate electron system. The question whether the contribution of *e*–*e* collisions to the conductivity is correctly described also in the low-density range has long been the subject of controversial debate. The *e*–*e* Coulomb interaction is included in the exchange-correlation energy of the DFT functional so that it was argued that the *ab initio* approach also describes the Coulomb interaction of the electrons.^{10,30} However, it was pointed out^{8} that DFT-MD is not able to reproduce the Spitzer limiting value of plasma conductivity, what leads to significant deviations in the low-density range. The correct behavior of the plasma conductivity at low density was determined using generalized linear response theory (gLRT).

^{19,31}see Fig. 1. Instead of the correct virial expansion (8) for the resistivity of the hydrogen plasma,

*k*,

*l*run over all

*N*particles. The pseudopotential is introduced in the framework of the DFT formalism, a simple version of this Kohn–Sham pseudopotential is the statically screened Debye–Hückel potential $ V c d ( r ) = e c e d e \u2212 \kappa r / ( 4 \pi \u03f5 0 r )$, which is often used for analytical calculations (see Sec. IV C). The Lorentz model plasma has the first virial coefficient

^{1,6,15}

*E*is also included. However, dynamical screening is not described by

_{xc}*E*, in contrast to the quantum Lenard–Balescu equation.

_{xc}We would like to mention that in the case of thermal conductivity, it has been shown that the contribution of *e*–*e* collisions is not taken into account in DFT-MD simulations^{10,22,30} and results in an additional term. Other approaches such as generalized linear response theory can be used to solve this problem.

The DFT-MD simulations have been used very successfully to calculate the transport properties of warm and dense matter (WDM) in the region of degeneracy. However, electron–electron collisions are not treated adequately. For degenerate systems, *e*–*e* collisions are suppressed by Pauli blocking so that their contribution can be neglected and the results of DFT-MD are applicable. The systematic inclusion of *e*–*e* collisions is indispensable for the description of plasmas at $ \Theta > 1$.

## IV. THE ACCOUNT OF ELECTRON–ELECTRON COLLISIONS IN NUMERICAL SIMULATIONS OF PLASMA CONDUCTIVITY

### A. The correction factor for Lorentz plasmas

We refer to a Lorentz plasma as a model system of electrons and ions, in which the interaction of electrons with ions is considered [see Eq. (32)]. The electron–electron interaction is only considered in the mean-field approximation so that the electron–ion interaction is screened, but the effect of electron–electron collisions on the transport coefficients is neglected. For this Lorentz model plasma, the solution of the linearized Boltzmann equation is found using the relaxation time approach. The relaxation time ansatz is possible for elastic collisions, which is the case in the limit $ m i / m e \u226b 1$, but not for *e*–*e* collisions.

For the Lorentz model plasma, the electrical conductivity $ \sigma Lorentz ( T , n )$ can be calculated beyond perturbation theory. The effective interaction of electrons with the ion subsystem is treated, for example, using the density functional theory. Similar approaches are given by the average atom model or the Lee–More model. These approaches are applicable in the high-density, degenerate case ( $ \Theta \u2264 1$), in which electron–electron collisions are suppressed due to Pauli blocking. They fail to describe real plasmas such as hydrogen plasmas in the low-density range, as *e*–*e* collisions are not taken into account.

*n*by Θ) approaches 1 for $ \Theta \u2264 1$ and the value $ 0.591 \pi 3 / 2 / 2 5 / 2 = 0.582$ for $ \Theta \u226b 1$.

Exact expressions for the conductivity $ \sigma ( T , n )$ of plasmas are obtained from the generalized linear response theory (gLRT) in the form of correlation functions. The kinetic theory (KT) follows as a special case.^{3} However, the analytical evaluation of the correlation functions is only possible in some approximations using Green function diagram techniques. For example, the dynamic screening of the electrons is treated in the random phase approximation (RPA), the ionic structure factor is also approximated, and the Born approximation can be improved by a T-matrix approach.

*Z*, a smooth interpolation formula was given, which reads for hydrogen (

*Z*= 1) in the low-density, classical limit [see Eq. (C.20) of Ref. 8]

*x*.

Within the framework of a refined calculation, we must investigate the effects of dynamical screening, ion–ion structure factor, and strong collisions (see Refs. 7, 15, 33, and 34). These effects are of minor importance for the correction due to *e*–*e* collisions, both in the high-density limit where $ lim \Theta \u2192 0 R e e ( \Theta ) = 1$ and in the low-density limit, where they only occur in higher orders of the virial expansion. Although these corrections are small, we analyze them in Sec. IV B.

### B. Example: DFT-MD simulations for hydrogen plasmas

As shown in Fig. 1, the DFT-MD simulations for the effective resistivity $ \rho 1 eff , DFT \u2212 MD ( T , n )$ are on a straight line pointing to the Lorentz limit (see Ref. 19). The corrected values $ \rho 1 eff , ee \u2212 corr ( T , n ) = \rho 1 eff , DFT \u2212 MD / R e e$ are shown in Table II. The correction factor *R _{ee}*, Eq. (35), leads to a considerable increase (up to a factor $ \u2248 2$ in the low-density range) of the resistivity and solves the drawback of the correct

*ρ*

_{1}benchmark by construction. However, the slope determining the second virial coefficient seems to remain too large compared to the high-temperature limit, Eq. (11), obtained from QLB or gLRT approaches.

^{15,17}While the slope of the DFT-MD simulations is 0.9965, the slope of the

*e*–

*e*corrected virial line is 0.732. The slope of the virial line in Fig. 1 determines the second virial coefficient $ \rho 2 ( T )$. The

*e*–

*e*corrected value is closer to the high-temperature benchmark 0.4917, but remains different.

n (g/cm^{3})
. | k (eV)
. _{B}T | Γ . | Θ . | $ 1 / \u2009 ln \u2009 ( \Theta / \Gamma )$ . | $ \sigma DFT \u2212 MD$ ( $ 10 6 / \Omega $ m) . | $ \rho 1 eff , DFT \u2212 MD$ . | $ \sigma ee \u2212 corr$ ( $ 10 6 / \Omega $ m) . | $ \rho 1 eff , ee \u2212 corr$ . |
---|---|---|---|---|---|---|---|---|

2 | 75 | 0.3285 | 1.8257 | 0.583 02 | 11.44 | 1.073 | 10.034 | 1.2234 |

2 | 100 | 0.246 37 | 2.4343 | 0.436 57 | 15.26 | 0.9269 | 12.434 | 1.1375 |

3 | 100 | 0.282 03 | 1.8577 | 0.530 47 | 16.85 | 1.020 | 14.412 | 1.1925 |

3 | 150 | 0.188 02 | 2.7866 | 0.370 92 | 25.67 | 0.8603 | 20.157 | 1.0956 |

4 | 150 | 0.206 94 | 2.3003 | 0.415 22 | 27.39 | 0.9026 | 22.058 | 1.1208 |

n (g/cm^{3})
. | k (eV)
. _{B}T | Γ . | Θ . | $ 1 / \u2009 ln \u2009 ( \Theta / \Gamma )$ . | $ \sigma DFT \u2212 MD$ ( $ 10 6 / \Omega $ m) . | $ \rho 1 eff , DFT \u2212 MD$ . | $ \sigma ee \u2212 corr$ ( $ 10 6 / \Omega $ m) . | $ \rho 1 eff , ee \u2212 corr$ . |
---|---|---|---|---|---|---|---|---|

2 | 75 | 0.3285 | 1.8257 | 0.583 02 | 11.44 | 1.073 | 10.034 | 1.2234 |

2 | 100 | 0.246 37 | 2.4343 | 0.436 57 | 15.26 | 0.9269 | 12.434 | 1.1375 |

3 | 100 | 0.282 03 | 1.8577 | 0.530 47 | 16.85 | 1.020 | 14.412 | 1.1925 |

3 | 150 | 0.188 02 | 2.7866 | 0.370 92 | 25.67 | 0.8603 | 20.157 | 1.0956 |

4 | 150 | 0.206 94 | 2.3003 | 0.415 22 | 27.39 | 0.9026 | 22.058 | 1.1208 |

The *e*–*e* correction factor, Eq. (35), was obtained from a two-moment approximation. The inclusion of higher moments of the distribution function would change the *e*–*e* correction factor and the corresponding values of $ \rho 1 , \rho 2 ( T )$. These corrections are expected to be small, see the good convergence with increasing number of moments reported in Ref. 15.

The virial coefficient $ \rho 2 ( T )$ is determined by the screening of the Coulomb potential. The Born approximation can be applied in the high temperature range $ T Ha \u226b 1$ so that strong collisions leading to a T-matrix approach are not relevant. A reason for the *ρ*_{2} discrepancy can also be found in the treatment of screening, see the analysis of the slope parameter in Ref. 31. The electrons should be treated as dynamic screening, in lowest approximation by the RPA expression for the polarization function, to obtain the QLB benchmark Eq. (11). It is not clear to what extent the dynamic screening by electrons is already implemented in the DFT-MD simulations by the choice of the exchange-correlation part *E _{xc}* of the energy-density functional.

If $ \sigma Lorentz ( T , n )$ which is required for the calculation of *R _{ee}*, Eq. (35), treats the dynamical screening of the electron–ion interaction differently than the DFT-MD approach $ \sigma DFT \u2212 MD ( T , n )$, the use of the correction factor is not consistent since the slope parameter $ \rho 2 ( T )$ is not corrected. With respect to the contribution of the ions to the screening, the ion structure factor describes the static screening, which is well approximated in the DFT-MD simulations.

In conclusion, the correction factor (35) given in Ref. 8 reproduces the general behavior for $ \Theta \u2264 1$ on $ R e e \u2248 1$ and is a reasonable approximation to the correct Spitzer result in the low-density limit. Only a two-moment approximation of the electron distribution function is considered to calculate *R _{ee}*, Eq. (35). It can be improved by considering higher moments of the electron distribution function so that the exact value of the Spitzer benchmark is reproduced. With respect to the second virial coefficient, the specific choice of the

*e*–

*i*interaction, which defines the Lorentz model plasma as a reference system, is discussed in Sec. IV C.

### C. Systematic approaches

While electron–ion collisions are well understood to calculate the conductivity of plasmas, the inclusion of electron–electron collisions is notoriously difficult. One reason for this is that a relaxation time approach is only valid for elastic collisions where the single-particle energy is conserved, what holds to a good approximation for *e*–*i* collisions, but not for *e*–*e* collisions. Several approaches have been worked out to implement *e*–*e* collisions. At low densities, this problem is solved within the framework of kinetic theory. There, a linearized Boltzmann equation can be derived for a weak electric field.^{1,5} The perturbation of the single electron distribution function results from the Kohler variational principle. In dense plasmas, gLRT is used to express the conductivity in terms of equilibrium correlation functions.^{3,6}

To calculate the conductivity of degenerate electrons in condensed matter, the Ziman formula is known, which considers the force–force autocorrelation function. Calculations of the collision frequency based on the solution of the electron–ion scattering problem are improved by a renormalization factor $ r ( \omega )$ to obtain the correct dynamic collision frequency of the plasma.^{8,34,35} While the Ziman formula is frequently used to describe degenerate Coulomb plasmas such as electrons in liquid metals or solids, the renormalization factor allows to treat also non-degenerate plasmas, especially in the low-density range. The method of including *e*–*e* correlations by means of a renormalization factor $ r ( \omega )$ is equivalent to using the correction factor *R _{ee}*, Eq. (34), and requires analogous approximations.

An improved approach to treat dense plasmas is to start from the Lorentz model plasma, in which the *e*–*i* interaction is solved, and the electron distribution function is rigorously obtained from the relaxation time ansatz for any degree of degeneracy. Starting from the linearized Boltzmann equation, thermoelectric transport coefficients for Lorentz model plasmas are calculated, see the Lee–More model.^{36} For a dense system, the averaged atom model was worked out, which describes the average influence of the surrounding plasma on the atom. The interaction of the electrons with the ion system is better described by the DFT approach, and the transition to condensed matter is well established. The electron–electron interaction is only considered in mean-field approximation to introduce an effective potential for the single-electron quasiparticle states.

The success of relaxation-time models for transport coefficients to treat the electron–ion interaction has led to develop a generalization of these models which includes *e*–*e* collisions. For example, as in Refs. 30 and 37, the relaxation-time model by Starrett^{38} includes Pauli blocking and accounts for correlations using a mean-force scattering cross section for electron–ion collisions, but electron–electron collisions are accounted for only through the correction formula for *R _{ee}* proposed in Ref. 8.

A possible solution for this particular problem of including *e*–*e* collisions in the determination of plasma conductivity can be the use of path-integral Monte Carlo (PIMC) simulations to calculate $ \sigma ( T , n )$. In PIMC simulations, *e*–*e* correlations are treated in an exact way. This enables the immediate calculation of $ R e e ( T , \Theta )$ according to (34) and the construction of consistent interpolation formulas. In particular, the determination of virial coefficients would be of interest. Such PIMC simulations are not yet possible. Highly accurate results are only available for the uniform electron gas.^{23} In addition, there is the unsolved sign problem, and PIMC simulations are currently performed for too small numbers of particles (a few tens). It is expected that significant results for the conductivity of hydrogen plasmas from PIMC simulations will be available in the near future. This will overcome the drawback of DFT-MD simulations where *e*–*e* correlations are not treated precisely.

We propose to use the virial expansion (8) for the determination of the conductivity of the hydrogen plasma in order to analyze the consistency of different theoretical approaches using the virial coefficients as a benchmark. Having $ \rho 1 , \rho 2 ( T ) ,$ etc. available, we can find improved versions of the correction factor $ R e e ( T , n )$. Since these correction factors relate the plasma conductivity $ \sigma ( T , n )$ to the conductivity of a Lorentz model plasma $ \sigma Lorentz ( T , n )$, see Eq. (34), these correction factors depend on the definition of the Lorentz model Hamiltonian. They depend, in particular, on the definition of the effective electron–ion potential, which includes dynamic or static screening by electrons, or on the particular choice of the exchange-correlation energy *E _{xc}* in DFT-MD calculations. If the determination of the conductivity and virial expansion for the Lorentz model plasma is known, a corresponding correction factor can be defined.

*Z*, so that it can also be applied to other substances.

## V. CONCLUSIONS

We propose an exact virial expansion (8) for the plasma conductivity to analyze the consistency of theoretical approaches. The virial coefficients serve as a benchmark for different approaches and can be better visualized using virial plots. For example, various analytical calculations of the dc conductivity $ \sigma ( T , n )$, which were presented in Ref. 20, do not fulfill this exact requirement and do not provide accurate results in the range of small densities. The results of DFT-MD simulations, which are currently considered to be the most reliable, are checked by virial plots, and shortcomings are found in the region of small densities. By benchmarking with the virial expansion (8) for $ x \u2192 0$, future PIMC simulations can also be tested. These *ab initio* simulations become a computational challenge in the low-density region, but the virial expansion enables extrapolation into this range. For the construction of interpolation formulas, see Ref. 18, knowledge of the virial coefficients is an important ingredient.

A particular problem concerns the treatment of the *e*–*e* interaction in DFT-MD simulations. To a certain extent, the *e*–*e* interaction is included in the exchange-correlation energy *E _{xc}* in order to obtain an effective mean field for optimal single-particle states. However, with the help of virial expansion, it was shown that

*e*–

*e*collisions are not included in the DFT-MD calculations. By introducing a correction factor

*R*, which is obtained from a model plasma with statically screened Coulomb interaction, the first virial coefficient

_{ee}*ρ*

_{1}can be reproduced. In this work, the problems with the second virial coefficient $ \rho 2 ( T )$, which is determined by the treatment of dynamic screening, are emphasized. Since the dynamic screening is not considered in the DFT-MD simulations, the second virial coefficient $ \rho 2 QLB$, Eq. (11), is also not obtained by the correction factor

*R*. A numerical calculation can be considered in the context of PIMC simulations to determine the second virial coefficient $ \rho 2 ( T )$.

_{ee}A correction factor $ R e e DFT \u2212 MD ( T , n )$ was introduced, which takes into account the specific choice of exchange-correlation energy in the DFT approach. The high-temperature limit of the virial expansion was given in Eq. (38), and the full *T* dependence would be of interest. The DFT-MD simulations are very accurate at higher densities since *e*–*e* correlations are strongly suppressed due to Pauli blocking in the degenerate range ( $ \Theta \u2264 1$). With the knowledge of $ R e e DFT \u2212 MD ( T , n )$, results for the conductivity of hydrogen plasmas can also be obtained from DFT-MD calculations for the non-degenerate range. In addition to the analytical evaluation of correlation functions for higher moments of the distribution function in the framework of gLRT, the numerical calculation in the context of PIMC simulations can be considered to determine the second virial coefficient $ \rho 2 ( T )$.

The approach described here can also be applied to other transport properties such as thermal conductivity, thermopower, viscosity, and diffusion coefficients.^{22} It is also interesting to extend the virial expansion to substances other than hydrogen, where different ions can be formed and the electron–ion interaction is no longer purely Coulombic.

## ACKNOWLEDGMENTS

The author would like to thank M. Bethkenhagen, M. French, R. Redmer, H. Reinholz, and M. Schörner for valuable and fruitful discussions.

## AUTHOR DECLARATIONS

### Conflict of Interest

The author has no conflicts to disclose.

### Author Contributions

**G. Roepke:** Conceptualization (lead); Investigation (lead); Writing – original draft (lead).

## DATA AVAILABILITY

The data that support the findings of this study are available within the article.

## REFERENCES

*Quantum Theory of Many-Particle Systems*

*Nonequilibrium Statistical Physics*