Ionization Potential Depression in Degenerate Plasmas and Pauli Blocking of Multi-Electron Ions

1. Introduction

The properties of atoms and ions immersed in a dense plasma are altered by medium effects. An important example is Debye screening, which modifies the energies of charged particles. Consequently, the ionisation potential changes in a dense plasma [1]. Simple approaches such as the Debye shift, which are valid at low densities, can be extended to higher densities, where other approaches such as the ion sphere model are more suitable. Interpolation formulas [2,3] have been proposed and compared with recent experiments on high-density plasmas [4,5,6,7,8,9]. A more fundamental approach [10], based on a quantum-statistical many-particle theory [11,12], is required for a systematic treatment of these plasma effects. In particular, the dynamic screening by free electrons and the inclusion of the ion structure factor were investigated in [13] and further references given there. More recently, Monte Carlo simulations using path integrals [14,15,16] have been employed.

New experiments with warm, dense matter (WDM) achieve high densities at which the plasma is degenerate, see Refs. [17,18,19,20]. For a plasma with temperature T and free electron density $n_e$, the electron degeneracy parameter is defined as

$$\Theta ={\displaystyle \frac{T}{T_{\mathrm{Fermi}}}}={\displaystyle \frac{2m_e{k}_{B}T}{{\hslash }^2}}{\left(3{\pi }^2{n}_e\right)}^{-2/3}$$

(1)

and describes the region $\Theta \le 1$ in the temperature-density plane in which a classical description is no longer valid. Instead, quantum effects, in particular the Pauli principle as a consequence of the antisymmetry of the many-electron wave function, lead to so-called exchange terms. The associated condition $n_e{\Lambda }_e^3/2\ge 1$ (spin factor 2) with ${\Lambda }_e={(2\pi {\hslash }^2/m_e{k}_{B}T)}^{1/2}$ is known as the condition at which the classical Boltzmann gas approximation is no longer applicable and the quantum description based on the Fermi distribution function must be applied. For example: For $T=100$ eV, the value $\Theta =1$ occurs at a density of $n_e=4.54\times {10}^{24}$ c${\mathrm{m}}^{-3}$. At the National Ignition Facility (NIF) and other facilities, new experiments are being conducted and planned to investigate WDM at free electron densities of up to ${10}^{25}$ c${\mathrm{m}}^{-3}$ and temperatures of up to $k_{B}T={10}^2$ eV, so that the effects of degeneracy are relevant.

A significant issue is the equation of state under such extreme conditions, particularly the composition. The degree of ionisation increases with rising T, but at high densities, pressure ionisation also occurs due to screening effects. Tables that take these effects into account are available, for example OPAL [21,22,23,24]. They are based on a many-particle treatment of screening [11,12,25].

In degenerate systems, new physical phenomena become increasingly significant. Whilst at lower densities, in the classical regime, dynamic screening is the most important medium effect, exchange interactions become increasingly important at these extremely high densities. These Pauli blocking effects are being investigated for Hydrogen plasmas; see [26]. They have been discussed for nuclear matter [27] to describe the dissolution of binding states in nuclear systems at high densities. A recent discussion of Pauli blocking effects in connection with the K-edge shift can be found in [28,29]. In particular, the ionisation potential depression (IPD) and Pauli blocking in degenerate Carbon plasmas at extreme densities were investigated in the article [30]. Due to Pauli blocking, the degree of ionisation at high densities is greater than predicted by standard approaches such as OPAL. These calculations are presented in explicit form in App. Appendix A, where all details are provided in that section so that the calculations can be reproduced immediately without additional sources. Our aim is to investigate the approximations used in these calculations. We mention the discrepancy with density functional theory–molecular dynamics (DFT-MD) simulations [31], which are also discussed in this work.

The influence of Pauli blocking on single-electron ions, i.e. H-like bound states, has already been discussed in several papers; see [26,30,32,33,34] and further references cited therein. In this work, we focus on multi-electron ions, so that other materials and ions can also be addressed. In particular, we discuss the problem of He and He-like bound states. Our aim is to contribute to the study of hot and dense matter consisting of multi-electron components in the region where the electrons are degenerate. We provide a brief overview of the in-medium Schrödinger equation (Sect. Section 2), the composition of the partially ionised plasma (Sect. Section 3) and the treatment of single-electron ions (Sect. Section 4). A new topic is the treatment of two-electron ions (He-like ions) in Sect. Section 5. A discussion of multi-electron ions, as well as a comparison with other approaches, can be found in Sect. Section 6, before we conclude with some final remarks (Sect. Section 7).

In this work, we use Rydberg units, where $\hslash =2m_e=e^2/8\pi {ϵ}_0=1$, so that energies are expressed in Ry, where 1 Ha = 2 Ry = 27.2114 eV = 219474.63 c${\mathrm{m}}^{-1}$. The Bohr radius is $a_{\mathrm{B}}=5.291772\times {10}^{-9}$ cm. We use atomic units $T_{\mathrm{Ry}}$ for the temperature T and $n_{\mathrm{Bohr}}=n{a}_{\mathrm{B}}^3$ for the particle density, so that

$$T_{\mathrm{Ry}}=6.33358\times {10}^{-6}{\displaystyle \frac{T}{\mathrm{K}}},n_{\mathrm{Bohr}}=1.4818471\times {10}^{-25}{n}_{\mathrm{e}}{\mathrm{cm}}^3.$$

(2)

The mass density includes the mass of the nucleus; for ${12}^{\mathrm{C}}$, ${\rho }_{\mathrm{C}}=134.481\times n_{\mathrm{nucl}}{a}_B^3$ g/c${\mathrm{m}}^3$.

2. Multi-Electron Schrödinger Equation

2.1. Free Multi-Electron Schrödinger Equation

We consider warm and dense matter composed of elements with proton number Z, for example $Z=1$ (Hydrogen), $Z=2$ (Helium), or $Z=6$ (Carbon). For neutral matter, the total density of nuclei $n_{\mathrm{nucl}}$ and the total number density of electrons $n_e^{\mathrm{tot}}$ are related as follows: $n_e^{\mathrm{tot}}=Z{n}_{\mathrm{nucl}}$. The interaction is described by Coulomb’s law $Q_1{Q}_2/\left(4\pi {ϵ}_{0}r\right)$, in the momentum representation

$$V_{Q_1,Q_2}(k_1,k_2;k_1^{\prime },k_2^{\prime })={\displaystyle \frac{Q_1{Q}_2}{{\Omega }_0{ϵ}_0{|{\mathbf{k}}_1-{\mathbf{k}}_1^{\prime }|}^2}}{\delta }_{{\mathbf{k}}_1+{\mathbf{k}}_2,{\mathbf{k}}_1^{\prime }+{\mathbf{k}}_2^{\prime }}{\delta }_{{\sigma }_1,{\sigma }_1^{\prime }}{\delta }_{{\sigma }_2,{\sigma }_2^{\prime }}$$

(3)

where $Q_i=Ze$ applies to nuclei with mass $m_Z$ and $Q_i=-e$ applies to electrons with mass $m_e$; ${\Omega }_0$ denotes the system volume. The quantum numbers $k=\mathbf{k},\sigma$ represent the wave number vector and the spin orientation of the interacting particles.

This interaction can give rise to bound states, which are determined by solving the Schrödinger equation. For example, for the isolated electron-proton system ($Z=0$), the energy levels are $E_n^0=-1/n^2$ Ry, where $n=1,2,\dots$. The degeneracy factor of each energy level is $2n^2$. In addition to the bound states, the solution to the Schrödinger equation also includes scattering states corresponding to the energy of relative motion $E\ge 0$.

For $Z>0$, in addition to the neutral atom in which Z electrons are bound to the nucleus, ions $Z_i$ can be formed, in which $Z-i$ electrons are bound to the nucleus, and the charge of the ion is $ie$, where $i=Z_i=0,1,\dots ,Z$. First, we revisit the case of the zero-density limit, where no density effects occur (index “0”). At zero density, the Schrödinger equation for the ion $Z_i$ describes the motion of the center of mass, which approximates the motion of the nucleus with charge $Ze$ and can be separated. The intrinsic motion corresponds approximately to the motion of the $Z_i$ electrons in the Coulomb field of the nucleus at a fixed position,

$$\begin{array}{cc}& \left({\displaystyle \frac{{\hslash }^2{k}_1^2}{2m_e}}+\dots +{\displaystyle \frac{{\hslash }^2{k}_{Z-i}^2}{2m_e}}-E_n^{0,i+}\right){\psi }_n^{0,i+}(k_1,\dots ,k_{Z-i})\hfill \\ & -\sum _{\mathbf{q}}{\displaystyle \frac{Z{e}^2}{{ϵ}_0{q}^2}}\sum _l^{Z-i}{\psi }_n^{0,i+}(k_1,\dots ,k_l+\mathbf{q},\dots ,k_{Z-i})\hfill \\ & +\sum _{\mathbf{q}}{\displaystyle \frac{e^2}{{ϵ}_0{q}^2}}\sum _{l<m}^{Z-i}{\psi }_n^{0,i+}(k_1,\dots ,k_l+\mathbf{q},\dots ,k_m-\mathbf{q},\dots ,k_{Z-i})=0.\hfill \end{array}$$

(4)

When there is more than one electron bound to a nucleus with proton number Z, the solution to the Schrödinger equation cannot be expressed analytically. Approximation methods are used to determine the solution, such as the shell model, in which the interaction between the electrons is replaced by a mean field. For the wave function ${\psi }_n^{0,i+}(k_1,\dots ,k_{Z-i})$ of the ion with charge $+ie$, the antisymmetrized product ansatz (Slater determinant) of single-particle wave functions can be applied, whereby the mean field is described in a self-consistent manner. More generally, density functional theory (DFT) was developed to obtain optimal single-particle wave functions. Numerical methods such as path integral Monte Carlo (PIMC) can also be applied to find solutions that take correlations between the electrons into account. We do not present here the success of atomic physics in solving Eq. (4). Instead, we use the empirical values of the energy levels $E_n^{0,i+}$ of the free ions, charge $Z_{i}e=ie$, excitation state n, and degeneracy factors $g_n^{i+}$, which can be found, for example, in the NIST tables [35].

2.2. In-Medium Single-Particle States

We are interested in the influence of the surrounding plasma on the energy eigenvalues $E_n^{i+}(T,n_{\mathrm{nuc}})$, which depend on the matter density $n_{\mathrm{nuc}}$ and the temperature T of the plasma. The quantum statistical solution to this problem is found using the Green’s function method, which introduces concepts such as the self-energy and the polarization function [11,25]. The perturbation expansion with respect to the Coulomb interaction can be represented by Feynman diagrams.

The single-particle propagator is determined by the self-energy, which depends on the wave number k and the energy. The energy of the free particle c in the quantum state k, $E_k^{c,0}={\hslash }^2{k}^2/\left(2m_c\right)$, is shifted by the self-energy, resulting in the quasiparticle energy. Since the self-energy is generally energy-dependent, the shift ${\Delta }_k^{c,\mathrm{qu}}$ is the value of the self-energy at the quasiparticle energy $E_k^{c,\mathrm{qu}}=E_k^{c,0}+{\Delta }_k^{c,\mathrm{qu}}$. In perturbation theory, the self-energy is expressed as the sum of terms represented by Feynman diagrams.

We consider an approximation of the self-energy given by the lowest-order Feynman diagrams. For the quasi-particle energy, we obtain the approximation

$$E_k^{c,\mathrm{qu}}=E_k^{c,0}+{\Delta }_k^{c,\mathrm{HF}}+{\Delta }_k^{c,\mathrm{scr}}.$$

(5)

The lowest-order contribution to the self-energy is the Hartree-Fock approximation. The Hartree term is zero due to charge neutrality. The Fock term is an exchange term that is relevant only for electrons that can be degenerate. It leads to a quasiparticle shift of the single-particle energies at wave number $\mathbf{k}$

$${\Delta }_k^{e,\mathrm{HF}}=-\sum _{k^{\prime }}{V}_{e,e}(k,k^{\prime };k^{\prime },k)f_e\left(k^{\prime }\right)$$

(6)

with the Fermi distribution function of electrons

$$f_e\left(k\right)={\displaystyle \frac{1}{exp(E_k^{e,0}/k_{B}T-{\mu }_e^0/k_{B}T)+1}}.$$

(7)

The electron chemical potential is given for the ideal (noninteracting) Fermi gas by

$$n_e={\displaystyle \frac{1}{{\pi }^2}}{\int }_0^{\infty }dk{\displaystyle \frac{k^2}{exp(E_k^{e,0}/k_{B}T-{\mu }_e^0/k_{B}T)+1}}.$$

(8)

The Fock term is a consequence of the Pauli principle, around an electron the electrons with the same spin orientation form the Fermi hole. This effect can be neglected in classical plasmas but becomes important in degenerate plasmas, such as the electron gas in metals.

In the next order of perturbation theory, Debye screening emerges, which is already the dominant process in classical plasmas. The screening in Coulomb plasmas is described by the polarization function $\Pi (q,\omega )$, which depends on the wave vector $\mathbf{q}$ and the frequency $\omega$. In the simplest approximation, we obtain the Random Phase Approximation (RPA)

$${\Pi }^{\mathrm{RPA}}(q,\omega )=\underset{ϵ\to 0}{lim}\sum _{c,k}{\displaystyle \frac{f_c\left(k\right)-f_c(\mathbf{k}+\mathbf{q})}{E_k^{c,0}-E_{\mathbf{k}+\mathbf{q}}^{c,0}-\omega +iϵ}}$$

(9)

which, in the static limit, leads to Debye screening of the Coulomb potential. The statically screened interaction $V_{Q_1,Q_2}^{\mathrm{scr}}\left(q\right)$ is given by the Debye interaction in the low-density range,

$$V_{Q_1,Q_2}^{\mathrm{Deb}}\left(\mathbf{q}\right)={\displaystyle \frac{Q_1{Q}_2}{{\Omega }_0{ϵ}_0(q^2+{\kappa }^2)}}$$

(10)

with the Debye screening parameter ($n_i$ is the density of ions, charge $+ie$)

$${\kappa }^2=\sum _i{Z}_i^2{e}^2{n}_i/\left({ϵ}_0{k}_{B}T\right)+{\kappa }_e^2,$$

(11)

$${\kappa }_e^2={\displaystyle \frac{8\pi }{k_{B}T}}{\left({\displaystyle \frac{2\pi {\hslash }^2}{m_e{k}_{B}T}}\right)}^{-3/2}{\displaystyle \frac{e^2}{4\pi {ϵ}_0}}{\displaystyle \frac{1}{\sqrt{\pi }}}{\int }_0^{\infty }dt{\displaystyle \frac{t^{-1/2}}{e^{t-{\mu }_e^0/k_{B}T}+1}}.$$

(12)

A charge $Z_{i}e=+ie$ immersed in the plasma generates a screening cloud. At low densities, the linearized Poisson-Boltzmann equation yields the Debye interaction. The expression for the Debye screening parameter $\kappa$ also includes the contribution of free electrons (${\kappa }_e$), which may be degenerate. The electron contribution is then given by a Fermi integral, which, in the strongly degenerate limit, yields the Thomas-Fermi screening length instead of the Debye screening length; see the references [11,12,21]. The energy shift of a particle with charge $Z_{i}e$ (including electrons, $Z_e=-1$) is

$${\Delta }_k^{i,\mathrm{scr}}\approx {\Delta }_k^{i,\mathrm{Debye}}=-{\displaystyle \frac{\kappa }{2}}{\displaystyle \frac{{\left(Z_{i}e\right)}^2}{4\pi {ϵ}_0}}.$$

(13)

The Debye shifts are merely the lowest-order terms in the perturbation expansion of the polarization function, which is represented by Feynman diagrams. The self-energy shifts have been improved; see [13], where the ion structure factor is included. Alternatively, semi-empirical interpolation expressions [2,3] are introduced to extend the Debye approach to higher densities. In the following, we consider the Stewart-Pyatt expression [3], which is frequently used in practical calculations; see the subsequent Section 2.3 and Appendix A.

2.3. Bound States, Density Effects by Screening of the Interaction

Two-particle correlations are described by the two-particle Green’s function, for which a Bethe–Salpeter equation can be derived. If the interaction with the medium is neglected, the two-particle Green’s function in the low-density limit is obtained by summing the ladder diagrams [25]

$$G_2^0(k_1,k_2;k_1^{\prime },k_2^{\prime };z)=\sum _{nP}\langle k_1,k_2|nP\rangle {\displaystyle \frac{1}{E_{nP}^0-z}}\langle nP:k_1^{\prime },k_2^{\prime }\rangle$$

(14)

where $E_{nP}^0$ and ${\psi }_{nP}^0(k_1,k_2)=\langle k_1,k_2|nP\rangle$ are the eigenvalues and eigenstates of the Schrödinger equation. P denotes the total momentum $\mathbf{P}$ and the total angular momentum, and n is the intrinsic quantum number describing bound and scattering states. The analytical continuation from the Matsubara frequencies into the complex z-plane is considered [11,25]. A two-particle Schrödinger equation describing the interaction of an electron with an ion is given by Eq. (4) for $i=Z-1$.

Similar Schrödinger equations can be formulated for multi-electron ions by considering the corresponding few-particle Green’s functions. In the low-density limit where in-medium effects can be neglected, the solution of the multi-electron Schrödinger equation (4) provides the expression for the Green’s function that generalizes equation (14).

In the case of Carbon, for example, we have binding states for ions consisting of the nucleus and $6-i$ electrons, where the charge of the ions takes the values $i=0,1,\dots ,5$. Instead of solving the corresponding multi-electron Schrödinger equation, we can use the empirical values of the observed energy levels, which can be found, for example, in the NIST tables [35]. The ionization energies $I_0^{0,i}$ are the minimum energies required to remove an electron from the ground state of the Carbon ion $C^{i+}$ (zero density): $I_0^{0,0+}=11.2603$ eV, $I_0^{0,1+}=24.38332$ eV, $I_0^{0,2+}=47.8878$ eV, $I_0^{0,3+}=64.5$ eV, $I_0^{0,4+}=392.087$ eV, $I_0^{0,5+}=489.805$ eV. In addition, these tables also list the excitation energies of the various ions.

We consider the modification of the multi-electron Schrödinger equation (4) in a dense medium. This quantum statistical problem is solved using the Green’s function method [11,25]. We provide a simplified expression; a more detailed discussion can be found in [11,12,36,37]. We replace the dynamic screening with a statically screened interaction. For an ion i with charge $Z_{i}e$ and an electron, which form a cluster with charge $(i-1)e$, we obtain

$$\begin{array}{cc}& \left[E_{k_1}^{e,\mathrm{qu}}+E_{k_2}^{i,\mathrm{qu}}\right]{\psi }_n^{i-1}({\mathbf{k}}_{\mathbf{1}},{\mathbf{k}}_{\mathbf{2}})\hfill \\ & +[1-f_e\left(k_1\right)]\sum _{\mathbf{q}}{V}_{e,i}^{\mathrm{scr}}\left(\mathbf{q}\right){\psi }_n^{i-1}({\mathbf{k}}_{\mathbf{1}}+\mathbf{q},{\mathbf{k}}_{\mathbf{2}}-\mathbf{q})\hfill \\ & =E_{n,k_2}^{i-1,\mathrm{qu}}{\psi }_n^{i-1}({\mathbf{k}}_{\mathbf{1}},{\mathbf{k}}_{\mathbf{2}})\hfill \end{array}$$

(15)

Here, $E_k^{c,\mathrm{qu}}=E_0^{c,\mathrm{qu}}+{\hslash }^2{k}^2/\left(2m_c\right)$ if the shifts are not depending on k. Due to the large mass ratio $m_i/m_e\gg 1$, we set $k_2$ equal to the center of mass momentum, which can be separated. The intrinsic motion follows from the Schrödinger equation in the medium

$$\begin{array}{cc}& \left[E_k^{e,\mathrm{qu}}+E_0^{i,\mathrm{qu}}-E_n^{i-1,\mathrm{qu}}\right]{\psi }_n^{i-1}\left(\mathbf{k}\right)\hfill \\ & +[1-f_e\left(k\right)]\sum _{\mathbf{q}}{V}_{e,i}^{\mathrm{scr}}\left(\mathbf{q}\right){\psi }_n^{i-1}(\mathbf{k}+\mathbf{q})=0.\hfill \end{array}$$

(16)

The statically screened interaction $V_{e,i}^{\mathrm{scr}}\left(q\right)$ is described in the low-density regime by the Debye interaction $V_{e,i}^{\mathrm{Deb}}\left(q\right)$ according to Equ. (10).

The polarization of the plasma by the charge Q leads to a shift in the quasi-particle energies. This part of the self-energy, which is due to correlations, is given at low densities by the Debye term $-\kappa Q^2/\left(8\pi {ϵ}_0\right)$, which does not depend on k. In the classical case, where $f_e\left(k\right)\ll 1$ so that degeneracy effects can be neglected, we expect the contribution to the ionization potential $I_0^{i-1}=E_0^{e,\mathrm{qu}}+E_0^{i,\mathrm{qu}}-E_0^{i-1,\mathrm{qu}}=I_0^{0,i-1}+\Delta I_0^{i-1}$ from the Debye shifts as

$$\Delta I_0^{i-1,\mathrm{Debye}}={\displaystyle \frac{\kappa }{2}}{\displaystyle \frac{e^2}{4\pi {ϵ}_0}}[{(i-1)}^2-i^2-1]=-{\displaystyle \frac{\kappa e^2}{4\pi {ϵ}_0}}i.$$

(17)

This result agrees with the solution of the in-medium Schrödinger equation (16), since, after expanding $e^{-\kappa r}\approx 1-\kappa r$, the Debye potential yields the contribution $-{\displaystyle \frac{\kappa e^2}{4\pi {ϵ}_0}}i$ to the energy $E_n^{i-1,\mathrm{qu}}$, in addition to ${\Delta }^{e,\mathrm{Debye}}+{\Delta }^{i,\mathrm{Debye}}$ according to the self-energy corrections. This reduction in the ionization potential by the term $\Delta I_0^{i-1,\mathrm{Debye}}$ is the ionization potential depression (IPD) in the Debye approximation.

In the case of Hydrogen, it was shown in [11,36,37] that for the bound state, the neutral hydrogen atom, the bound state energy shift caused by the screening of the interaction is exactly compensated by the shift due to the self-energy of the constituents in lowest order of density. For the reaction considered here, ${\mathrm{C}}^{5+}\rightleftharpoons {\mathrm{C}}^{6+}+e$, the Debye shift of the initial bound state is $\Delta E_0^{5+,\mathrm{Debye}}=-\kappa Z_5^2{e}^2/\left(8\pi {ϵ}_0\right)$ while for the final states, ${\Delta }_0^{6+,\mathrm{Debye}}=-\kappa Z_6^2{e}^2/\left(8\pi {ϵ}_0\right)$ for the Carbon nucleus and ${\Delta }^{e,\mathrm{Debye}}=-\kappa e^2/\left(8\pi {ϵ}_0\right)$ for the electron results. The difference yields the reduction in ionization energy as $\Delta I_0^{5+,\mathrm{Debye}}=\kappa (5^2-6^2-1)e^2/\left(8\pi {ϵ}_0\right)=-6\kappa e^2/\left(4\pi {ϵ}_0\right)$.

The Debye shifts are merely the lowest-order terms in the perturbation expansion of the polarization function, which is represented by Feynman diagrams. The self-energy shifts have been refined; see [13]. Alternatively, semi-empirical interpolation expressions are introduced to extend the Debye approximation to higher densities. In practical calculations, the Stewart-Pyatt expression is frequently used, which interpolates between the Debye shift at low densities and the ion sphere model at high densities, where the ions are strongly correlated. The ionization potential depression (IPD) results from the energy shifts of the free-state continuum and the bound-state shift, such that the ionization energy is altered by the medium,

$$\Delta I_0^{i-1,\mathrm{SP}}={\displaystyle \frac{3}{2}}{\displaystyle \frac{(Z_i+1)e^2}{4\pi {ϵ}_0{r}_{\mathrm{IS}}}}\left[{(1+s^3)}^{2/3}-s^2\right]$$

(18)

with the ion sphere radius $r_{\mathrm{IS}}={(3Z_i/\left(4\pi n_e\right))}^{1/3}$ and $s={\left(\kappa r_{\mathrm{IS}}\right)}^{-1}$.

In the parameter range considered here, ${\Delta }_0^{i,\mathrm{SP}}$ depends only weakly on T and i. We use this expression for our calculations, App. Appendix A. Within the framework of a more systematic quantum statistical approach, the Stewart-Pyatt expression [3] for the IPD can be improved using the dynamical structure factor [13]. There, a slightly stronger IPD is observed.

2.4. Degeneracy Effects

Compared to classical plasmas, the Pauli principle, which applies to fermions, gives rise to new effects when the plasmas are degenerate. We use the relation (1), where $n_e$ denotes the free electron density.

For single-particle states, the first-order contribution is the Fock shift; see Eq. (6). The Fock shift of an electron with momentum $\hslash \mathbf{k}$ is given by

$${\Delta }_k^{e,\mathrm{Fock}}=-\sum _q{\displaystyle \frac{e^2}{{\Omega }_0{ϵ}_0{q}^2}}{f}_e(\mathbf{k}+\mathbf{q}).$$

(19)

In the limiting case of strong degeneracy, $T\ll T_{\mathrm{F}}$, we approximate the Fermi distribution function as a step function, $f_e\left(k\right)=\Theta (k_{\mathrm{F}}-k)$. The Fermi wave number is given by $k_{\mathrm{F}}={\left(3{\pi }^2{n}_e\right)}^{1/3}$. The calculation of the Fock shift yields

$${\Delta }_k^{e,\mathrm{Fock}}=-{\displaystyle \frac{e^2}{4\pi {ϵ}_0}}{\displaystyle \frac{1}{\pi k}}\mathrm{Re}\left[k{k}_{\mathrm{F}}+(k_{\mathrm{F}}^2-k^2)arctan\left({\displaystyle \frac{k}{k_{\mathrm{F}}}}\right)\right].$$

(20)

If $k=0$, the following holds for the shift of the continuum edge: ${\Delta }_0^{e,\mathrm{Fock}}=-{\displaystyle \frac{2}{\pi }}{\displaystyle \frac{e^2}{4\pi {ϵ}_0}}{k}_F$. Whilst screening makes the largest contribution to the single-particle shifts in classical plasmas, the Fock shift dominates in the limiting case of strong degeneracy.

In degenerate plasmas, electrons in bound states are also subjected to the Pauli principle. Here we consider the simplest case: an electron bound to a nucleus of charge Z. The two-particle Schrödinger equation without medium corrections yields the well-known hydrogen-like solution for the ground state

$$E_0=-Z^2{\displaystyle \frac{e^4}{{\left(4\pi {ϵ}_0\right)}^2}}{\displaystyle \frac{m_e}{2{\hslash }^2}}=-Z^2\mathrm{Ry}$$

(21)

and

$${\varphi }_0\left(k\right)=8\sqrt{\pi a_Z^3}{\displaystyle \frac{1}{{(1+a_Z^2{k}^2)}^2}},{\psi }_0\left(r\right)={\displaystyle \frac{1}{\sqrt{\pi a_Z^3}}}{e}^{-r/a_Z}$$

(22)

with $a_Z={\displaystyle \frac{4\pi {ϵ}_0}{Z{e}^2}}{\displaystyle \frac{{\hslash }^2}{m_e}}$ (for instance, $a_C=1/6$ in Rydberg units).

At low densities, where the effects of the medium are weak, we can apply perturbation theory to calculate the in-medium shifts in the energies of all quasiparticles, i.e., of the free components and of the bound states. The shift in the bound state energy due to the Fock self-energy is then calculated using the unperturbed solution ${\varphi }_0\left(k\right)$,

$$\Delta E_0^{\mathrm{Fock}}=-\sum _{k,q}{\varphi }_0^2\left(k\right)V_{ee}\left(q\right)f_e(\mathbf{k}+\mathbf{q})=-{\displaystyle \frac{32}{\pi }}{\int }_0^{\infty }dk{\displaystyle \frac{k^2{a}_Z^3}{{(1+a_Z^2{k}^2)}^4}}{\Delta }_e^{\mathrm{Fock}}\left(k\right).$$

(23)

The Pauli blocking shift at $T=0$ is given by (see Ref. [30])

$$\Delta E_0^{\mathrm{Pauli}}=-\sum _{k,q}{\varphi }_0\left(k\right)f_e\left(k\right)V_{ei}\left(q\right){\varphi }_0(\mathbf{k}+\mathbf{q})={\displaystyle \frac{Z{e}^2}{4\pi {ϵ}_0}}{\displaystyle \frac{2}{\pi a_Z}}\left[{\displaystyle \frac{a_Z{k}_F(a_Z^2{k}_F^2-1)}{{(a_Z^2{k}_F^2+1)}^2}}+arctan\left(a_Z{k}_{\mathrm{F}}\right)\right].$$

(24)

For any T, we find

$$\begin{array}{ccc}\hfill \Delta E_0^{\mathrm{Pauli}}& =& {\displaystyle \frac{Z}{2{\pi }^3}}{\int }_0^{\infty }dkk{\int }_0^{\infty }d{k}^{\prime }{k}^{\prime }{\varphi }_0\left(k\right){\varphi }_0\left(k^{\prime }\right)\hfill \\ & & \times {\displaystyle \frac{1}{2}}[f_e\left(k\right)+f_e\left(k^{\prime }\right)]ln\left[{\displaystyle \frac{{(k+k^{\prime })}^2+{\kappa }^2}{{(k-k^{\prime })}^2+{\kappa }^2}}\right].\hfill \end{array}$$

(25)

In the special case $\kappa \to 0$, we have

$$\begin{array}{ccc}\hfill \Delta E_0^{\mathrm{Pauli}}& =& {\displaystyle \frac{1}{2{\pi }^2}}{\int }_0^{\infty }dk{k}^2{\varphi }_0^2\left(k\right)(k^2+Z^2)f_e\left(k\right).\hfill \end{array}$$

(26)

Similar expressions can be derived for the shifts of the excited states. Further details can be found in Sect. Section 4.

Bound states are also regarded as quasi-particles. This is referred to as the chemical picture. In a quantum statistical approach, a corresponding propagator is introduced by equation (14). As a result of the interaction of the bound states with the plasma, the energy eigenvalues are modified, $E_n=E_n^0+\Delta E_n$, where $\Delta E_n=\Delta E_n^{\mathrm{Fock}}+\Delta E_n^{\mathrm{Pauli}}$. The wave functions ${\varphi }_n\left(\mathbf{k}\right)$ are also modified. As the density increases, the influence of the plasma grows, and the bound state energies can merge with the continuum, causing bound states to disappear. This dissolution of the bound states is known as the Mott effect and has important consequences for the macroscopic properties of the plasma; further details can be found in [36,37].

3. Composition

3.1. The Composition of the Partially Ionized Plasma

The total nuclear density $n_{\mathrm{nucl}}^{\mathrm{tot}}$ can be decomposed as follows (the ions are non-degenerate):

$$n_{\mathrm{nucl}}^{\mathrm{tot}}={\displaystyle \frac{1}{{\Omega }_0}}\sum _{i,n,p}{e}^{-\beta E_n^i\left(p\right)+\beta {\mu }_i}=\sum _{i=0}^Z{n}_i.$$

(27)

In this cluster decomposition of the total density, the various channels are denoted by i which is related to the number of particles; the center of mass momentum is $\hslash \mathbf{p}$; and n denotes the intrinsic state of the cluster, such as excited bound states. For example, $i=0,1,\dots ,Z$ represents the ionization state and describes the number of particles (the nucleus and $Z-i$ electrons) for which correlations are taken into account. The sum over the spin orientation and the total angular momentum yields a degeneracy factor. The intrinsic quantum state n encompasses not only the excited bound states, but also the continuum states.

For each component of the plasma, a corresponding chemical potential ${\mu }_i$ can be introduced. In chemical equilibrium, ${\mu }_i$ does not depend on n, and in equation (27) we have

$${\mu }_i={\mu }_{\mathrm{nucl}}+(Z-i){\mu }_e.$$

(28)

The chemical potentials ${\mu }_{\mathrm{nucl}},{\mu }_e$ of the elementary constituents are in turn related to each other by the condition of neutrality. It should be noted that larger clusters $i>Z$ can also be included, such as the formation of ${\mathrm{H}}^{-}$ or molecules; however, this would go beyond the scope of this work.

Since the ions can be regarded as classical particles, the sum/integral over $\mathbf{p}$ is simply evaluated, yielding the following result (where ${\Omega }_0$ is the volume of the system):

$${\displaystyle \frac{1}{{\Omega }_0}}\sum _p{e}^{{\hslash }^2{p}^2/\left(2m_Z{k}_{B}T\right)}=\int {\displaystyle \frac{d^{3}p}{{\left(2\pi \right)}^3}}{e}^{{\hslash }^2{p}^2/\left(2m_Z{k}_{B}T\right)}={\left({\displaystyle \frac{m_Z{k}_{B}T}{2\pi {\hslash }^2}}\right)}^{3/2}={\Lambda }_{\mathrm{nucl}}^{-3},$$

(29)

the mass $m_Z$ of the various ions is nearly identical. The relation (27) can be rewritten as

$$n_{\mathrm{nucl}}^{\mathrm{tot}}={\Lambda }_{\mathrm{nucl}}^{-3}\sum _{i=0}^Z{e}^{\beta {\mu }_i}{\sigma }^i,n_i={\Lambda }_{\mathrm{nucl}}^{-3}{e}^{\beta {\mu }_i}{\sigma }^i$$

(30)

with the intrinsic partition function

$${\sigma }^{0,i}=\sum _n{e}^{-\beta E_n^i}=e^{-\beta E_0^i}\sum _n{e}^{-\beta E_n^{i,*}}$$

(31)

where $E_n^{i,*}=E_n^i-E_0^i$ denotes the excitation energies. For the free ions, these excitation energies are found in tables, e.g. [35], together with the degeneracy factors according to the quantum numbers of angular momentum (including spin).

A further rewriting is used in the literature, introducing the ionization potentials $I_0^i$. The ionization potential is the minimum energy required to extract an electron from the ion i in the ground state, energy $E_0^i(p=0)$. The final state after ionization has the energy $E_0^{i+1}(p=0)+E_k^e$ where $\hslash \mathbf{k}$ gives the momentum of the extracted electron. If $k=0$ has the lowest possible energy of the emitted electron, the edge of the continuum gives the ionization potential

$$I_0^i=E_0^{i+1}\left(0\right)+E_0^e-E_0^i\left(0\right).$$

(32)

Values for the free components are found in the tables, e.g. [35]. In a dense medium, we have to use the quasiparticle energies. For degenerate electrons at zero temperature, only final states above the Fermi energy are possible.

We can replace in Eq. (30)

$$E_0^i=E_{\mathrm{nucl}}+(Z-i)E_0^e-\sum _{i^{\prime }=0}^{Z-i}{I}_0^{Z-i^{\prime }+1},$$

(33)

where the lowest energies of the free components, $E_{\mathrm{nucl}}$ and $E_0^e$, are taken as zero in our calculations. To eliminate both ${\mu }_{\mathrm{nucl}}$ and ${\mu }_e$, a second equation for $n_e^{\mathrm{tot}}$ is required. Similar to eq. (27), we have

$$n_e^{\mathrm{tot}}=\sum _{i=0}^Z(Z-i)n_i+n_e.$$

(34)

Part of electrons are localized within the ions, the other part moves freely. Since the electron mass is small relatively to the ion mass, the free electrons are described by quantum statistics which gives at arbitrary degeneracy the Fermi distribution,

$$n_e={\displaystyle \frac{1}{{\Omega }_0}}\sum _k{f}_e^{\mathrm{qu}}\left(k\right)=2\int {\displaystyle \frac{d^{3}k}{{\left(2\pi \right)}^3}}{\displaystyle \frac{1}{e^{\beta (E_k^{e,\mathrm{qu}}-{\mu }_e)}+1}}$$

(35)

(spin factor 2). Coming back to Eq. (8), which relates the density to the chemical potential in the ideal Fermi gas approximation, the relation (35) takes interaction into account. This expression for the free electron density determines also screening, $n_e^{\mathrm{free}}/n_e^{\mathrm{total}}$ is the ionization degree. The densities of the different components of the plasma are related by the neutrality condition ${\sum }_i{Z}_i{n}_i={\sum }_{i}i{n}_i=n_e$.

3.2. The Intrinsic Partition Function ${\sigma }^i(n,T)$

The intrinsic partition function (31) is obtained if the interaction between the components of the plasma, the free and bound states, is neglected. However, an ideal, noninteracting plasma with occasional reactions to establish chemical equilibrium is not appropriate for a dense plasma where interactions have to be taken into account. This is consistent with the account of excited states which include also scattering states. We are confronted with two problems: (i) the treatment of in-medium effects which leads to the quasiparticle concept, the in-medium shift of energy levels and medium-modified scattering phase shifts, and (ii) the subdivision into free and correlated density of the electrons.

Solving these problems requires a quantum statistical approach where well-defined quantities such as the spectral function are introduced. We give here only the main ideas. To discuss the long-standing problem of the subdivision of the electron density into a bound part and a free part, we consider the second virial coefficient for a short-range potential. An exact expression is given by the Beth-Uhlenbeck formula, for a given channel $\alpha$ (spin, angular momentum, etc.) we have

$$\begin{array}{cc}& {\sigma }_{\alpha }^i\left(T\right)=\sum _n^{\mathrm{bound}}{e}^{-E_{\alpha ,n}/k_{B}T}+{\int }_0^{\infty }{\displaystyle \frac{dE}{\pi }}{e}^{-E/k_{B}T}{\displaystyle \frac{d}{dE}}{\delta }_{\alpha }\left(E\right)\hfill \end{array}$$

(36)

where ${\delta }_{\alpha }\left(E\right)$ is the scattering phase shift in the channel $\alpha$. Both, $E_{\alpha ,n}$ and ${\delta }_{\alpha }\left(E\right)$, can be determined from experiments. After integration by parts and using the Levinson theorem, as well as introducing the quasiparticle energies, we obtain

$$\begin{array}{cc}& {\sigma }_{\alpha }^i\left(T\right)=\sum _n^{\mathrm{bound}}(e^{-E_{\alpha ,n}^{\mathrm{qu}}/k_{B}T}-1)+{\int }_0^{\infty }{\displaystyle \frac{dE}{\pi k_{B}T}}{e}^{-E/k_{B}T}\left\{{\delta }_{\alpha }\left(E\right)-{\displaystyle \frac{1}{2}}sin\left[2{\delta }_{\alpha }\left(E\right)\right]\right\}.\hfill \end{array}$$

(37)

In this second version,

(i) we can assume that the contributions of the continuum remain small and can be neglected. We use this relation to define ${\sigma }^i\left(T\right)={\sum }_n(e^{-E_n/k_{B}T}-1)$ as the contribution of correlations (bound states).

(ii) An additional sin-term appears in the contribution of scattering states. This term is related to the change in the single-particle contribution (35) where quasiparticles are considered. The quasiparticle shifts such as the Hartree-Fock shift contain already part of the interaction, this must be extracted from the second virial coefficient to avoid double counting.

(iii) Also for the two-particle properties, bound state energies and scattering phase shifts, the in-medium corrections have to be considered, in particular the quasiparticle energy shifts of the bound states.

In the case of Coulomb interaction, scattering phase shifts can not be defined in the standard way because of the long-range character of the interaction. This problem has been investigated in plasma physics since a long time, see Refs. [11,12,21]. Already in the approximation where the interaction between the components of the partially ionized plasma is neglected, the definition of the intrinsic partition function ${\sigma }_i^0\left(T\right)$, Eq. (31), contains problems as described by Planck, Larkin, and others [38,39,40]. The sum over all bound states for the Coulomb system is divergent, and truncation must be performed to obtain convergent results.

For the Hydrogen atom, a convergent expression was proposed by Planck and Larkin [38,39,40]

$${\sigma }^i\left(T\right)=e^{-\beta (E_0^{i+1}+E_0^e)}\sum _n^{\mathrm{bound}}\left[e^{\beta I_n^i}-1-\beta I_n^i\right]$$

(38)

where n covers the intrinsic quantum numbers (including spin) of all bound states, $I_n^i=E_0^{i+1}+E_0^e-E_n^i$. For free electrons, the state with lowest energy is $k=0$ so that $E_0^e=0$. For H, $i=0$, is ${\mathrm{H}}^{1+}$ the free nucleus which has also the minimum energy $E_0^{1+}=0$ so that $I_n^0=1/n^2$. For He, $i=0$, is H${\mathrm{e}}^{1+}$ the He ion with one electron, and the minimum energy is $E_0^{1+}=-4$ Ry. The additional last term in the partition function (38) is related to the Debye shift of the quasiparticle energy which is necessary for Coulomb systems, it changes the analytical form of the virial expansion [41,42], and terms of the form $n_e^{1/2}$ appear. The case of Hydrogen is widely discussed in the literature, see Ref. [34] and further references given there.

The remaining part of the continuum of scattering states may be neglected. Note that the subdivision of the intrinsic partition function into a bound state part and a scattering part is not free of ambiguity so that the definition of the ionization degree is model dependent. One possibility for linking the concept of free electron density to physical properties is based on the screening of the Coulomb interaction potential in plasmas [41,42]. The polarisation function $\Pi (q,\omega )$ can provide us with a free electron density defined by its behavior at $\omega =0$ and $q\to 0$.

4. The Single Bound Electron Problem

We investigate the motion of an electron in the screened potential of an ion with charge $Ze$. This model describes the final stage of ionisation in a plasma at increasing density; it defines the transition to a fully ionised plasma. Furthermore, it also describes the mean-field approximation, in which the potential of the nucleus is replaced by an effective potential of the ion. This model is suitable, for example, for a valence electron on top of filled shells. The problem of two electrons moving in the mean-field potential is discussed in Sect. Section 5.

As explained in Sect. Section 2, within the framework of a quantum statistical approach, the in-medium Schrödinger equation is given by (15). As usual, after separation of the momentum of the center of mass motion (which practically coincides with the momentum of the ion) the relative motion remains, which describes the intrinsic structure, i.e., the bound and scattering states, of the electron-ion system. Energy levels and wave functions result from the solution of the corresponding in-medium Schrödinger equation for the relative motion (16),

$$\left[E_k^{e,\mathrm{qu}}+E_0^{i+1,\mathrm{qu}}-E_n^{i,\mathrm{qu}}\right]{\psi }_n^i\left(\mathbf{k}\right)-\left[1-f_e\left(k\right)\right]\int {\displaystyle \frac{d^{3}q}{{\left(2\pi \right)}^3}}{\displaystyle \frac{8\pi Z}{q^2+{\kappa }^2}}{\psi }_n^i(\mathbf{k}+\mathbf{q})=0.$$

(39)

The quasi-particle energy $E_k^{e,\mathrm{qu}}$ (5) includes the Fock shift ${\Delta }_k^{e,\mathrm{HF}}$, Eq. (19), and the shift ${\Delta }_k^{e,\mathrm{scr}}$ due to screening. Dynamic screening has already been considered in the literature; see [36,37]; for simplicity, we assume static screening here, which leads to the Debye term in Eq. (13). In this section, we use Ry units.

The edge of the continuum states $E_{\mathrm{cont}}^{i,\mathrm{qu}}$ from Equ. (39) is given by $k=0$: $E^i,{\mathrm{qu}}_{\mathrm{cont}}=E_0^i+1,\mathrm{qu}+{\Delta }_0^{e,\mathrm{HF}}+{\Delta }_0^{e,\mathrm{scr}}$. For a given T, the edge of the continuum states shifts downward as the density increases. Expressions for the Debye shift of the ion, charge $(Z_i+1)e$, and the electron are given in Eq. (13). The Fock term follows from Eq. (19) for $k\to 0$.

In the non-degenerate case, we consider the screening terms. Since the Debye shifts do not depend on k, they can be combined, and the difference

$$E_0^{e,\mathrm{qu}}+E_0^{i+1,\mathrm{qu}}-E_n^{i,\mathrm{qu}}=I_n^i$$

(40)

is the ionization potential of the state ${\psi }_n^i\left(\mathbf{k}\right)$. In this classical case, the ionization potential depression is given by the Debye shifts, which can be refined using the Stewart-Pyatt expression (18).

In the degenerate case, the Fock self-energy term ${\Delta }_0^{e,\mathrm{HF}}$ appears in the shift of the continuum edge, but ${\Delta }_k^{e,\mathrm{HF}}$ also appears in the calculation of the bound states. If we assume that the k-dependence can be neglected, we can rewrite equation (39) in the form

$$\left[k^2+I_n^i\right]{\psi }_n^i\left(\mathbf{k}\right)-\left[1-f_e\left(k\right)\right]\int {\displaystyle \frac{d^{3}q}{{\left(2\pi \right)}^3}}{\displaystyle \frac{8\pi Z}{q^2+{\kappa }^2}}{\psi }_n^i(\mathbf{k}+\mathbf{q})=0.$$

(41)

As explained above, the appearance of $\kappa$ in the interaction without the Pauli blocking factor $\left[1-f_e\left(k\right)\right]$ already describes the classical part of the IPD, which is enhanced by the Pauli blocking terms in the degenerate case. Neglecting the k-dependence of the Fock shift is justified near the Mott point, where the wavefunction approximates that of a free electron. Although it is possible to account for the k-dependence, it has only a minor effect on the IPD; see Ref. [30].

The integral equation (41) is not Hermitean but can be symmetrized introducing ${\psi }_n^i\left(\mathbf{k}\right)/{(1-f_e\left(k\right))}^{1/2}={\varphi }_n^i\left(\mathbf{k}\right)$. We can assume that the ground state is spherical symmetric (s wave, c.m. momentum 0) so that

$$\left[k^2+I_n^i\right]{\varphi }_n^i\left(k\right)-{(1-f_e\left(k\right))}^{1/2}\int {\displaystyle \frac{d^3{k}^{\prime }}{{\left(2\pi \right)}^3}}{\displaystyle \frac{8\pi Z}{{(\mathbf{k}-{\mathbf{k}}^{\prime })}^2+{\kappa }^2}}{(1-f_e\left(p^{\prime }\right))}^{1/2}{\varphi }_n^i\left(k^{\prime }\right)=0.$$

(42)

Performing the angular integral, we obtain

$$\left[k^2+I_n^i\right]{\varphi }_n^i\left(k\right)-{(1-f_e\left(k\right))}^{1/2}{\displaystyle \frac{Z}{\pi }}{\int }_0^{\infty }d{k}^{\prime }{\displaystyle \frac{k^{\prime }}{k}}ln\left[{\displaystyle \frac{{(k+k^{\prime })}^2+{\kappa }^2}{{(k-k^{\prime })}^2+{\kappa }^2}}\right]{(1-f_e\left(k^{\prime }\right))}^{1/2}{\varphi }_n^i\left(k^{\prime }\right)=0.$$

(43)

or, with ${\varphi }_n^i\left(k\right)=u_n^i\left(k\right)/k$, the symmetric form

$$\left[k^2+I_n^i\right]u_n^i\left(k\right)-{(1-f_e\left(k\right))}^{1/2}{\displaystyle \frac{Z}{\pi }}{\int }_0^{\infty }d{k}^{\prime }ln\left[{\displaystyle \frac{{(k+k^{\prime })}^2+{\kappa }^2}{{(k-k^{\prime })}^2+{\kappa }^2}}\right]{(1-f_e\left(k^{\prime }\right))}^{1/2}{u}_n^i\left(k^{\prime }\right)=0.$$

(44)

We present results for $Z=6$, carbon, and calculate the energy shifts for the final ionization process, the reaction ${\mathrm{C}}^{5+}\leftrightharpoons {\mathrm{C}}^{6+}+e$. The disappearance of the ionization potential determines the emergence of the fully ionized carbon plasma. We begin with the question: How good is perturbation theory? In particular, is there an exact density value (the Mott density) at which the bound states disappear?

To illustrate this, let us consider the case of $C^{5+}$. The Debye shift of the continuum ($C^{6+}+e$) is ${\Delta }^{\mathrm{cont}}=-\kappa Z_{\mathrm{C}}{e}^2/\left(4\pi {ϵ}_0\right)$, see Eq. (17). This result can be refined using the Stewart-Pyatt expression or other results. If we consider the Debye interaction $-Z_{\mathrm{C}}{e}^2/\left(4\pi {ϵ}_{0}r\right)exp(-\kappa r)\approx -Z_{\mathrm{C}}{e}^2/\left(4\pi {ϵ}_{0}r\right)(1-\kappa r)$, perturbation theory yields a shift in the bound state energy of $+Z_{\mathrm{C}}{e}^2/\left(4\pi {ϵ}_0\right)\kappa$. This compensation effect was found in [36,37]. In the case of H, the bound state energy of the neutral H atom remains unchanged, while the bound state energy of ${\mathrm{C}}^{5+}$ and the continuum exhibit the additional Debye shift $-(Z_{\mathrm{C}}-1)\kappa e^2/\left(4\pi {ϵ}_0\right)$. Therefore, the IPD of ${\mathrm{C}}^{5+}$ is mainly determined by the medium-induced modification of the interaction term, i.e., by dynamic screening and Pauli blocking. The situation is illustrated in Fig. 1 of Ref. [30].

We are concerned with the solution of Eq. (44) or, after discretization,

$$\sum _{\beta }\left[k_{\alpha }^2{\delta }_{\alpha \beta }-{(1-f_e\left(k_{\alpha }\right))}^{1/2}{\displaystyle \frac{Zd}{\pi }}ln\left[{\displaystyle \frac{{(k_{\alpha }+k_{\beta })}^2+{\kappa }^2}{{(k_{\alpha }-k_{\beta })}^2+{\kappa }^2}}\right]{(1-f_e\left(k_{\beta }\right))}^{1/2}\right]u_n^i\left(k_{\beta }\right)=-I_n^i{u}_n^i\left(k_{]}\alpha \right),$$

(45)

with $k_{\alpha }=d(\alpha -1/2),\alpha =1,2,\dots N_{\alpha }$, and $\alpha ⟷\beta$. The limit $N_{\alpha }\to \infty ,d\to 0$ must be calculated until convergence is achieved.

It is interesting to compare the exact solution of the in-medium Schrödinger equation with perturbation theory. The perturbation approach based on the undisturbed wave function is not valid near the Mott transition; the wave function is strongly distorted due to medium effects. The Mott transition requires an exact description, regardless of whether the transition is sharp or an adiabatic approximation to the continuum states is present.

The calculations are shown in Fig. Figure 1. For the Debye potential $-2Z/rexp(-\kappa r)$ (Ry units), the lowest-order perturbation theory yields a shift of $2Z\kappa$. The exact solution yields a Mott point at ${\kappa }_{\mathrm{Mott}}=7.1$ instead of the value 3 from perturbation theory. However, there is a sharp transition at which the bound state disappears and a resonance in the continuum appears.

Numerical calculations using Eq. (45) were performed for $N_{\alpha }\le 30$ and various values of d, with the optimal value lying in the range $d=1.5$. The results can be improved by considering larger values of $N_{\alpha }$. To assess the accuracy, we first present results for various $\kappa$ values for ${\mathrm{C}}^{5+}$. In this work, we have not investigated the convergence of our calculations and therefore cannot provide error bars for the calculations. We can only highlight the general trends.

We treat $\kappa$ as a parameter and solve the corresponding Schrödinger equation without the Pauli blocking terms. Perturbation theory yields the shift $\kappa Z{e}^2/\left(4\pi {ϵ}_0\right)$ for the Debye potential, while the solution to the Schrödinger equation exhibits a deviation that is also approximated by the numerical solution. Deviations of the numerical solution from the solution of the differential equation occur for small $\kappa$, since the logarithm in Eq. (45) becomes very sharp and must be resolved.

Furthermore, Fig. Figure 1 illustrates the effect of Pauli blocking. We solve the Schrödinger equation with the Pauli-blocking terms for a given free electron density and $T=100$ eV, which determines the Fermi functions, treating $\kappa$ as a variable. The perturbation theory based on the free bound state does not yield useful shifts, since the ground state is already modified by Pauli blocking.

We present consistent calculations for the IPD $I_0^{5+}$ of C-plasmas, in which, for given values of T and $n_e$ of free electrons, the Fermi functions according to (8) (${\mu }_e$) and the screening parameter $\kappa$ according to (12) are calculated. The results are shown in Figure 2, where various temperatures were considered. The abrupt disappearance of the bound states is evident. The dependence on T is only weak.

The screening of the Coulomb interaction is a dynamic process; see equation (9) for the polarisation function. The use of a statically screened Debye potential in the in-medium Schrödinger equation (39) to describe the contribution of the ions is not justified at high densities. In particular, the ions cannot follow the rapid motion of the electrons. The ions are strongly correlated, so they can make only a small contribution to the screening. Therefore, in our calculations we consider only the electron component of the screening. To include the ion contributions, an improved treatment based on dynamic screening is required.

It should be noted that the solution to Eq. (45) yields also both excited states and continuum states. The contribution of resonances can be determined from the density of states. According to the generalized Beth-Uhlenbeck formula, we can separate the contribution of the single quasiparticles from the correlated part, which contains both bound states and correlations in the continuum. The contribution of continuum correlations can become significant at low temperatures, where the Fermi function can be replaced by the step function. This also solves the problem of the increase in ionization energy at $T=0$, since all states below the Fermi energy are occupied.

Calculations were also performed for other materials with $Z\le 6$. The ionization potentials for the hydrogen-like bound states are shown in Fig. Figure 3 as a function of the free electron density at $k_{B}T=100$ eV. The bond states dissolve due to screening and Pauli blocking at Mott densities; see Table 1.

5. Two-Electron Ions

The multi-electron Schrödinger equation (4) for an ion embedded in a plasma includes medium corrections, the self-energy shifts of the single-particle energies, and the Pauli blocking terms of the interaction. Standard solutions for this few-particle problem introduce an mean potential in which the Coulomb potential of the nucleus is screened by the average density of the other electrons. The total wave function of the electrons is approximated by the antisymmetrised product of single-electron states. We discuss the case of two electrons bound to a nucleus, i.e. He-like ions. To estimate the effect of Pauli blocking, we consider the wave function in momentum space, where the occupation of phase space by the free electrons in the plasma becomes apparent.

We begin with the isolated neutral helium atom, i.e., with vanishing plasma density. Two electrons are bound to a nucleus with $Z=2$; the ground-state energy is -5.80677 Ry [35]. Within the framework of a variational approach, we consider the product ansatz (spin singlet) for the wave function

$${\psi }_0(r_1,r_2)={\displaystyle \frac{1}{\pi a^3}}{e}^{-r_1/a}{e}^{-r_2/a}$$

(46)

With the Hamilton operator, see Eq. (4), the total energy is calculated as

$$E_0={\displaystyle \frac{{\hslash }^2}{m_e}}{\displaystyle \frac{1}{a^2}}-{\displaystyle \frac{e^2}{4\pi {ϵ}_0}}{\displaystyle \frac{27}{8a}}.$$

(47)

It reaches its minimum value of $-3^6/2^7\mathrm{Ry}=-5.6953$ Ry for $a=16/27a_B=0.5926a_B$. This estimate of the binding state energy is slightly higher than the empirical bound state energy of -5.80677 Ry. The difference is due to the simple form of the wave function and the neglect of correlation effects; see Hylleraas [43].

The effective potential contains the Coulomb potential of the nucleus, charge Z, which is screened by the other electron,

$$V_1^{\mathrm{eff}}\left(r_1\right)=-{\displaystyle \frac{Z{e}^2}{4\pi {ϵ}_0{r}_1}}+{\displaystyle \frac{1}{2}}{\displaystyle \frac{e^2}{4\pi {ϵ}_0}}\int d^3{r}_2{\psi }_2^{*}\left(r_2\right){\displaystyle \frac{1}{|{\mathbf{r}}_1-{\mathbf{r}}_2|}}{\psi }_2\left(r_2\right)\approx -{\displaystyle \frac{Z_1^{\mathrm{eff}}{e}^2}{4\pi {ϵ}_0{r}_1}}.$$

(48)

The use on an effective charge $Z_1^{\mathrm{eff}}=27/16=1.6875$ allows to use Hydrogen like orbitals,. It reproduces the wave function determined by the minimum value of the bound state energy, but it also reproduces the value of the bound state energy.

To investigate the influence of the plasma, we consider the case of strong degeneracy, in which the leading terms are the Fock self-energy shift and the Pauli blocking. Since the self-energy shift of a weakly bound electron and a free electron is almost equal, the dissolution of the bound state is mainly due to Pauli blocking. To estimate the shift of the ground state energy, we use the mean-field approximation, in which both electrons move independently on hydrogen-like orbitals. Pauli blocking means that the part of the orbitals, which are already occupied by the free electron states of the plasma, cannot be accessed. The corresponding shift in the bound state energy is given by Eq. (26). In the strongly degenerate case, where the Fermi function can be replaced by the step function, we obtain the expression (24).

A calculation for $k_{B}T=100$ eV is shown in Fig. Figure 4. The bound state energy of H${\mathrm{e}}^0$ at $n_e=0$ is -5.6953 Ry, which is slightly higher than the empirical value -5.80677 Ry. For comparison, the wave function in which only one electron is bound whilst the other is free has the bound state energy -4 Ry (H${\mathrm{e}}^{1+}$); the difference is the ionization potential $I_0^2$ for He.

As the density of free electrons increases, screening and Pauli blocking lead to a shift of the bound state energy, and the two-electron bound state merges into the continuum at $n_e=0.56a_B^{-3}$. We compare this with another test wavefunction in which one electron is bound with the wavefunction $a_1$ ($a_1=a_B/2$ at zero density), whilst the other electron is free ($a_2=\infty$). The bound electron experiences the charge $Z=2$ without screening by the second electron. This state transitions into the continuum at a higher value $n_e=0.71a_B^{-3}$. As can be seen from Fig. Figure 3, the Mott density increases with increasing Z.

The symmetric solution, in which both electrons occupy the same orbit, is no longer the lowest-energy state as the density of free electrons increases. Due to the in-medium shifts, above $n_e^{\mathrm{free}}=0.3a_B^{-3}$, the asymmetric wave function, in which one electron is found at the edge of the continuum, becomes more favorable. If Pauli blocking is taken into account, the Pauli shift becomes smaller, as the wave function is more concentrated in coordinate space and more extended in momentum space; there is less overlap with the Fermi sphere.

A preliminary discussion of the Pauli shift for the neutral He atom in a degenerate plasma was recently presented in [34], in which only the symmetric solution was considered. The approach outlined in this paper describes the dissolution of the neutral He atom in the degenerate plasma at increasing density as a stepwise ionisation process via the formation of the H${\mathrm{e}}^{1+}$ ion. We conclude that the disappearance of the bound state is not the simultaneous ionisation of both electrons, but a step-by-step, sequential ionisation as density increases.

A similar approach can be used to describe the shift of He-like ions in other materials with a charge number Z of the nuclei in the medium. Assuming the symmetric product approach (46), the energy as a function of the parameter a is given by the following (in Ry units):

$$E\left(a\right)={\displaystyle \frac{2}{a^2}}-{\displaystyle \frac{4Z}{a}}+{\displaystyle \frac{5}{4a}},Z^{\mathrm{eff}}={\displaystyle \frac{1}{a^{\mathrm{eff}}}}=Z-{\displaystyle \frac{5}{16}}$$

(49)

is obtained as the minimum of $E\left(a\right)$. For carbon, for example, the energy minimum yields $Z^{\mathrm{eff}}=5.6875$. As in the case of He, the pressure ionisation of ${\mathrm{C}}^{4+}$ in the degenerate plasma occurs step by step, one after the other, via the formation of ${\mathrm{C}}^{5+}$, due to Pauli exclusion principle.

A similar approach can be used to describe the in-medium shift of He-like ions of other materials, with charge number Z of the nuclei. Assuming the symmetric product ansatz (46), the energy as function of the parameter a follows as (Ry units)

$$E\left(a\right)={\displaystyle \frac{2}{a^2}}-{\displaystyle \frac{4Z}{a}}+{\displaystyle \frac{5}{4a}},Z^{\mathrm{eff}}={\displaystyle \frac{1}{a^{\mathrm{eff}}}}=Z-{\displaystyle \frac{5}{16}}$$

(50)

follows as the minimum of $E\left(a\right)$. For instance, for Carbon we find from the minimum of energy. $Z^{\mathrm{eff}}=5.6875$. As in the case of He, the pressure ionisation of ${\mathrm{C}}^{4+}$ in the degenerate plasma due to Pauli blocking occurs stepwise, subsequentally, via the formation of ${\mathrm{C}}^{5+}$.

Ions containing more than two electrons can be treated in a similar way. Using a product ansatz for the many-electron wavefunction, we can introduce an effective mean-field potential. The resulting single-electron orbits describe the wavefunctions of the quasi-particle states, which are partially blocked due to the Pauli principle as the density of free electrons increases. This approach to the many-electron wavefunction, which neglects the formation of correlations, reduces the calculation of energy shifts to the single-electron problem, which evolves in an effective self-consistent potential. A simple case is that of Li-like ions, in which a third electron moves n top of the core of the $1s$ electrons, for example ${\mathrm{C}}^{3+}$. The most important prerequisite for calculating Pauli blocking is knowledge of the many-electron wavefunction of the ion; however, in this work, we cannot enter into this broad field of research, which is dominated by density functional theory (DFT).

6. Comparison: DFT, Averaged Atom Model

We show that theories for classical plasmas lose their validity as soon as the plasma becomes degenerate. Simulation programs such as OPAL [21] do not account for Pauli blocking, but this is an important effect that destroys bound states at high densities. Experiments at high densities suggest higher ionization degrees, as predicted by standard approaches that do not account for the Pauli blocking effect [17,18,19].

DFT offers an approach based on a rigorous quantum mechanical framework that takes the Pauli principle into account. Calculations were carried out [31], but complete ionisation could not be achieved. It remains an open question whether DFT calculations can describe the composition of a partially ionised plasma. DFT is a quasi-particle approach in which correlations between electrons are not adequately taken into account. As is known from the treatment of the ${\mathrm{H}}_2$ molecule, in an approach where electrons move in the effective mean field of the other particles, polar states can arise that are suppressed by Coulomb repulsion. The occurrence of bound states means that an electron is located at a specific ion and is not distributed throughout the entire system. Such correlations may not be correctly described within the DFT formalism.

One approach that models correlation effects in a semi-empirical manner is the average-atom model; see [44]. The introduction of the Wigner–Seitz cell to solve the Schrödinger equation reflects the Pauli principle empirically to a certain extent, as the presence of other electrons reduces the probability that the bound electron will be found outside the Wigner–Seitz cell. Although this model yields plausible results, it does not represent a fundamental approach. The Pauli principle is not formulated in position space, but in the space of quantum states, such as momentum states. Our quantum statistical theory, which uses the Green’s function method, is suitable for a more fundamental approach.

7. Conclusions

We present a detailed study of Pauli blocking in degenerate plasmas that goes beyond perturbation theory. The numerical solution of the in-medium Schrödinger equation yields, for a given temperature, a Mott density that is greater than the value derived from perturbation theory. The disappearance of the bound state remains abruptly. Values for the Mott density of the H-like bound state are given for various plasmas.

We have extended our investigation of Pauli blocking in degenerate plasmas to the case of multi-electron ions. In particular, we have considered He and He-like ions. We have employed a shell model approach, in which a product state is assumed to consist of individual electron orbits moving within an effective mean-field potential. The Planck–Larkin partition function is generalized to multi-electron ions.

We show that the Mott transition of the He-like bound state is not associated with the simultaneous dissociation of the two equivalent electrons as the density increases. It is a stepwise process leading from the two-electron ion via the one-electron ion to the fully ionised ion. The symmetric solution, in which both electrons bound to the ion are equivalent, breaks down into a solution in which one electron is free. The remaining electron is more strongly bound, as there is no screening by the second electron. An abrupt transition from the symmetric to the asymmetric solution is demonstrated.