# Ionization potential depression and ionization balance in dense plasmas

###### Abstract

Theoretical modelling of ionization potential depression and the related ionization equilibrium in dense plasmas, in particular in warm/hot dense matter, represents a significant challenge due to ionic coupling and electronic degeneracy effects. We present a quantum statistical model based on dynamical structure factors for the ionization potential depression, where quantum exchange and dynamical correlation effects in plasma environments are consistently and systematically taken into account in terms of the concept of self-energy. Under the condition of local thermodynamic equilibrium, the charge state distribution (or ionic fraction) characterized by the ionization balance is obtained by solving the coupled Saha equations. Calculations for the ionization potential depression of different chemical elements are performed with the electronic and ionic structure factors. The ionic structure factors are determined by solving the Ornstein-Zernike equation in combination with the hypernetted-chain equation. As a further application of our approach, we present results for the charge state distribution of aluminium plasmas at several temperatures and densities.

## I Introduction

Within the plasma community, one of the well-known many-body effects is ionization potential depression (IPD) or continuum lowering. The IPD significantly alters the ionization balance and the corresponding charge state distribution, which has a strong influence on transport, optical, and thermodynamic properties of the system. This is of essential importance not only for plasma physics but also for astrophysics, planetary science, and solid state physics. However, an accurate description of IPD in an interactive many-body environment is extremely complicated, because strong correlations and quantum degeneracy have to be taken into account consistently. Different semi-empirical models have been proposed for the IPD. In particular, the Ecker-Kröll (EK) EK63 () and the Stewart-Pyatt (SP) model SP66 () have been widely applied in plasma physics for simulations of physical properties and the analysis of experimental observations.

With newly developed experimental facilities, it has become possible to explore warm/hot dense matter and materials in the high-energy density regime, whereby the physical systems are found to be strongly coupled and nearly degenerate. The experimental outcomes can be used to benchmark physical assumptions and theoretical models commonly applied within plasma physics. Over the last few years, new experiments related to the IPD have been performed with high-intense laser beams in LCLS at SLAC National Accelerator Laboratory ciricosta12 (); ciricosta16 (), at ORION in the UK Hoarty13 (), and in NIF at Lawrence Livermore National Laboratory Fletcher14 (); Kraus16 (). These precise experiments have revealed the lack of a consistent picture about IPD, since no self-consistent explanation for these experiments can be drawn from any of the commonly accepted IPD models. Therefore, a more fundamental understanding of the phenomenon IPD in warm/hot dense regime is required. Recently, attempts for a better understanding of the new experimental data have been made using different numerical approaches and simulation methods. Two-step Hartree-Fock calculations STJZS14 (), simulations based on the finite-temperature density functional theory VCW14 (); Hu17 (), classical molecular dynamics simulations Calisti15 (), and Monte-Carlo simulations Stransky16 () have been worked out. Although these state-of-the-art modelling methods are very successful to explain properties of warm dense matter, they are restricted due to the large computational cost. Besides these simulations methods, other analytic improvements have also been proposed, such as fluctuation model IS13 () and atomic-solid-plasma model Rosmej18 ().

In view of the demands from the theoretical as well as experimental aspect, we have, starting from a quantum statistical theory, derived an analytic self-consistent approach to IPD, which is demonstrated to be valid in a wide range of temperatures and densities LRKR17 (). The influence of the surrounding plasma on an embedded ion/atom is described by the self energy (SE) which contains the dielectric function KKER86 (); KSKB05 (). The dielectric function can be expressed in terms of the dynamic structure factor (SF) according to the fluctuation-dissipation theorem. Since the dynamic SF comprises the time and space correlations of particles in the plasma, a microscopic understanding of the IPD in this quantum statistical approach is obvious. The developed quantum statistical model for IPD is based on the assumption of a two-component plasma model with ionic and electronic subsystem. In a realistic physical problem, the investigated system usually consists of different ion species WCK17 (), for example, H-C mixture in laboratory physics Kraus16 (); LR12 (). Because of the large mass and charge asymmetry, the light element moves more liberally between the highly charged, strongly correlated heavy components, which generally causes an additional dynamical screening effect for the calculation of IPD.

In this work we give a detailed description of our quantum statistical approach for IPD, which is based on the model developed in the previous work LRKR17 () and is extended to describe a multicomponent plasma. The present work is organized as follows: in Sec. II, we outline the basis of single-particle SE in the GW approximation, where G denotes the dressed Green’s function and W indicates the dynamically screened interaction potential. In terms of the SE, definition of IPD within the quantum statistical theory is introduced in Sec. II.3. In the subsequent Sec. III, we demonstrate that the IPD is connected to charge-charge dynamical SF, which describes dynamical correlation effects in plasmas. Quantum exchange (statistical) correlations including the Fock contribution and the Pauli blocking effect are also taken into account in this section. Using the effective ionization potential defined within the developed IPD model, the derivation of the coupled Saha equations in the chemical picture is discussed in Sec. IV. In the high-temperature ideal plasma limit, the Debye-Hückel model is exactly reproducible from our approach, as shown in Sec. V.1. As applications of our method, the IPD for different chemical elements are investigated in Sec. V.2, where comparisons with experimental observations ciricosta16 () are performed. Then we apply the developed IPD model to calculate the charge state distribution of Al plasmas corresponding to the thermodynamical conditions in experiments of Hoarty et al. Hoarty13 () in Sec. V.3. Finally, conclusions are drawn in Sec. VI.

## Ii Self energy and ionization potential depression

### ii.1 Plasma parameters

We consider a multicomponent mixture consisting of free electrons with charge and mass as well as ions of different species with charge and mass in a volume . is the elementary charge. The partial particle number density for ion species is and the total particle number density for all ions is . The corresponding number concentration is then given by . According to the charge neutrality, the electron density is with the mean ionization degree of ions . Additionally, we introduce the effective charge number of plasma ions , which effectively describes the plasma as a whole and therefore regards the ionic perturbers in plasma as a single ionic species.

In a many-body environment the motion of particles is correlated with the motion of their nearby particles. The coupling strength of such correlation is represented by the dimensionless plasma parameter

(1) |

which is taken as the ratio of the average unscreened interaction potential between type and type to the thermodynamic kinetic energy characterized by the temperature . In general, ions and electrons in multicomponent charged particle systems can have different temperatures with and , respectively. The electron-ion interaction temperature has been constructed with different ansatzes, see Refs. GRHGR07 (); HFBKRY17 (). The averaged interparticle distance reads

(2) |

Additionally, the quantum degeneracy for the electron subsystem is defined via the ratio of the thermodynamic kinetic energy and the Fermi energy as follows

(3) |

### ii.2 Single-particle self energy in GW approximation

We calculate the single-particle SE (SPSE) in the GW-approximation. In the Lehman representation KKER86 (), the Green’s function for species in the energy-momentum space, i.e. , can be expressed in terms of the spectral function

(4) |

where is the Matsubara frequencies with for fermions and for bosons. Similarly, the screened interaction potential can be written in the spectral representation via the inverse dielectric function

(5) |

with the Coulomb interaction . The inverse dielectric function describes the response of a many-body system to an external perturbation. Then the SPSE in GW approximation

(6) |

can be decomposed into a Hartree-Fock (HF) contribution due to quantum exchange effects and a correlation one because of dynamical interactions

(7) |

where the HF SE is given by

(8) |

and the correlation part of the SPSE reads

(9) |

In the following, the sum over momentum , i.e. , is replaced by the integral . Performing the summation over the Matsubara frequencies yields

(10) | |||

(11) |

where the distribution function reads

(12) |

with the upper sign for fermions denoted as , and the lower sign for bosons with . is the chemical potential of species . The dynamical effects within mediums are described by the function that is given by

(13) | ||||

The HF SE, Eq. (11), has no dependence on the frequency and is a real quantity. The first contribution of the HF SE is denoted as Hartree term which vanishes for a homogeneous system because of charge neutrality KKER86 (). The second contribution is the so-called Fock term, which arises from exchange correlation of identical particles and has no classical counterpart. We only calculate the Fock term in this work, but still refer this contribution as Hartree-Fock

(14) |

The correlation contribution, i.e. Eq. (10), can be split into a real and an imaginary part by means of the analytic continuation (i.e. Wick rotation) with KKER86 (); KSKB05 (); SL13 (). After the analytic continuation the function can be rewritten as

(15) | ||||

with the use of Dirac’s identity ( denotes the principal value). Generally, we have for the non-degenerate ions, in particular, for the ion involved in the ionization process. Obviously, the correlation contribution of the SPSE, i.e. Eq. (10), can be decomposed into a real part (related to the shift of eigenstates) and an imaginary part (connected to the broadening of eigenstates).

### ii.3 Definition of ionization potential depression within the quantum statistical theory

Generally, the IPD in a medium is defined as the change of the ionization potential with respect to the isolated case. It can be extracted from the effective binding energy, which is obtained by solving the Schrödinger equation with an effective interaction potential. The IPD acquired from the solution of the standard Schrödinger equation (i.e. eigenenergies for both scattering and bound eigenstates) is a real quantity. However, the quantum eigenstates of the investigated system in a many-body environment are visibly broadened due to the fast oscillating part of the microfield generated by the surrounding charged particles. In order to account for both the shift and the broadening of quantum eigenstates, the Bethe-Salpeter equation (or the in-medium Schrödinger equation) has to be solved KKER86 (). Therefore, the commonly defined IPD can be extended to a complex quantity, which includes the standard IPD and an additional contribution due to broadening effects.

Alternatively, the problem of generalized IPD (GIPD) can also be tackled within the Green’s function technique for the quantum statistical theory LRKR17 (). In the framework of quantum statistical theory, the modifications of the atomic/ionic properties are described by the SE. In the chemical picture, the investigated system undergoing ionization is treated as different ionic species before (ion ) and after the ionization (ion plus an ionized electron). Apparently, the SE of the investigated system before and after the ionization is significantly changed, while the surrounding environment is assumed to be not changed during the ionization process if the relaxation effect is negligible. Under this assumption, the GIPD can be defined via the difference between the SEs of the corresponding investigated system before and after the ionization. In the following, the letters in the subscripts and in the summation represent both ions and electrons. If only ionic species are involved in equations, the Greek letters will be used.

The SPSE of species is expressed via the dressed propagator (Green’s function) and the screened interaction potential KKER86 (); KSKB05 () (see also Eq. (6))

(16) |

Then the frequency- and momentum-dependent GIPD is given by

(17) |

The commonly defined IPD is described by the real part of the above introduced GIPD

(18) |

and the broadening of the IPD is given as

(19) |

The electronic SE contains two contribution: the SPSE of the ionized electron as a free particle and energy modification of the ionized electron in its parent ion as well as the restriction of phase space occupation in the bound-free ionization process, i.e.

(20) |

Concerning bound-free transition of the ionized electron, the following many-particle effects have to be considered in the SE :(1) energy levels of the bound states are shifted due to dynamical interaction (collisional shift) and quantum exchange effect (Fock shift) between bound electron and free electrons in plasmas; (2) bound electrons can not be ionized to those states that are already occupied by the free electrons (known as Pauli blocking); (3) the energy levels are also broadened because of random collision of bound electron with its surrounding charged particles (known as pressure broadening). In this work, we only consider the Pauli blocking and the Fock shift of bound states in , i.e.

(21) |

Then the GIPD can be defined as

(22) |

where the momentum- and frequency-dependence are suppressed. The basic quantity describing the GIPD is the SPSE , which can be evaluated within different approximations. In the next section, we will discuss the GW approximation for the SPSE , the Pauli blocking of the ionized electron , and the Fock shift of bound states .

Based on the SPSE, we can define the frequency- and momentum-dependent GIPD in this section and its reduced version in the next section. Inserting Eq. (7) in combination with Eq. (8) and Eq. (II.2) into Eq. (22) yields

(23) |

with the Hartree-Fock contribution

(24) |

and correlation contribution from dynamical interactions

(25) |

As mentioned before that all HF SEs are real quantities, hence contributes only to the commonly defined IPD. Here we collect all contributions from statistical correlation in the expression of GIPD

(26) |

Inserting Eq. (10) into Eq. (25) and introducing the following expression

(27) |

to describe atomic properties of the ion involved in the ionization reaction, the dynamical correlation part of the GIPD, i.e. Eq. (25), can be rewritten as

(28) |

With the expression (15) the real and imaginary part of the dynamical correlation contribution (28) take the following expressions

(29) | ||||

(30) |

The essential quantities determining the GIPD are the spectral function and the dielectric function . The spectral function of particle will be discussed in detail in the subsequent subsection. The dielectric function describes dielectric response of a material to an external field. A large amount of physical properties, such as stopping power, optical spectra, and conductivity, are directly connected to the dielectric function KKER86 (); KSKB05 (). The dielectric function is extremely complicate to be determined. A famous approximation for the dielectric function is the random phase approximation AB84 (), which describes the effective response in a non-interacting gas. Improvements based on the random phase approximation can be performed, such as local field correlation UI80 (); FWR10 (), Born-Mermin ansatz Mermin65 (), and the extended Born-Mermin ansatz RSWR99 (); SRWRPZ01 (). Attempts to develop dielectric function for two-component plasma have been proposed by different authors SRWRPZ01 (); Roepke98 (); RW98 (); AAADT14 (). However, all these improvements are inadequate to describe the effective charge response in a multicomponent plasmas under warm dense matter conditions. In this work we express the dielectric function in terms of the charge-charge dynamical SF according to the fluctuation-dissipation theorem. Furthermore, the proposed approach in this work can be improved by taking higher-order cluster contributions, e.g. through considering the multiple interaction in the T-matrix approximation LRKR17 (). Such improvements will not be handled within the context of the present investigation.

## Iii Ionization potential depression: statistical and dynamical correlations in many-body systems

In this section we demonstrate that the IPD in a multi-component Coulomb system can be directly related to the spatial distribution and temporal fluctuation of the plasma environment, i.e. the dynamical SFs. The SPSE incorporated in the GIPD includes the HF term and the correlation contribution. Correspondingly, the GIPD can be decomposed into two contributions, i.e. contribution from the statistical (or exchange) correlation and from the dynamical correlation, respectively. The dynamical correlation contribution describes dynamical interactions between the investigated system (denoted as impurity) and its surrounding charged environment (treated as perturber). The exchange correlation includes the HF contribution of the continuum edge , the Pauli blocking , and Fock shift for bound states stemming from spin statistics of identical particles. In this section we will discuss the statistical and dynamical correlations in many-body systems, where the real part of GIPD, i.e. the commonly defined IPD, is of central relevance. The uncertainty of IPD, characterized by broadening of the continuum edge and energy levels, will be briefly discussed in the present study. More details will be presented in a forthcoming work.

### iii.1 Approximation for spectral function

The spectral function contains all information about the dynamical behavior of a particle in an interacting many-body environment and satisfies the normalization condition . It is related to the SPSE as follows KKER86 (); KSKB05 (); SL13 ()

(31) |

As shown in Eq. (6), the SPSE is further connected to the spectral function, so that the spectral function and the SPSE have to be determined self-consistently. In the case of non-interacting gases or in the case of negligible width () of weakly interacting gases, the Lorentz form of spectral function, i.e. Eq. (31), can be replaced by a simple -shape

(32) |

Such simplified treatment of the spectral function results in the approximation of the SPSE, where denotes the undressed Green’s function for free particles. Note that we actually do not know which approximation is better for the evaluation of SE. It is well know in condense matter physics that the self-consistent GW approximation usually results in a good quasi-energy but overestimates the energy gap, while approximation with screened potential in RPA level gives a better energy gap Holm99 (); KNO16 (). For a consistent approach, vertex corrections have also to be included in the calculation.

In the following, we perform the approximation for the SPSE to calculate the IPD in multicomponent plasmas. Then we have for the function , i.e. Eq. (II.3), in the expression of GIPD (28)

(33) | ||||

As discussed in Ref. LRKR17 () (for details also see Appendix A), we take the momentum after performing classical dispersion relation for the investigated ion , i.e. . In other words, we define the GIPD as . In order to achieve an analytic expression for the IPD, a further simplification in the derivation is to take the classic limit in the propagator in Eq. (29) for the real part of the correlation contribution. Under these approximations we have

(34) |

with which the real part of dynamical interaction contribution in the GIPD, i.e. Eq. (29), is given by

(35) |

Similarly, the broadening of the IPD, i.e. Eq. (30), can be described by

(36) |

Detailed derivations of Eq. (34), Eq. (35) and Eq. (III.1) are given in Appendix A.

### iii.2 IPD due to statistical correlation: Hartree-Fock contribution and Pauli blocking

In this subsection we discuss the statistical correlation contribution of the GIPD

(37) |

As described in Sec. II.3, here the Hartree-Fock contribution to the continuum edge and the energy shift of bound states due to Pauli blocking and Fock shift are taken into account.

#### iii.2.1 Hartree-Fock contribution to the continuum edge

Within the approximation, the HF term of the SPSE contained in the HF contribution in the GIPD, i.e. , is given by

(38) |

Inserting the explicit expression for the distribution function , i.e. (12), into Eq. (38) yields

(39) |

where the upper/lower sign corresponds to the case of fermions/bosons according to spin statistics and is the polylogarithm function. For non-degenerate ions, the Fermi-Dirac or the Bose-Einstein distribution can be approximated by the Maxwell-Boltzmann distribution. This approximation yields the following expression of HF SE

(40) |

Consequently, the HF contribution to the continuum edge, i.e. , is given by

(41) |

with the electronic contribution

(42) |

and the ionic contribution

(43) |

Here the chemical potentials of ideal gases for ions are used. For non-degenerate ions, is applied, where is the statistical weight of the ground state of ionization state . For electronic chemical potential , the expression (102) is utilized in the calculation. The results for a plasma at temperature eV are shown in Fig. 1. It can be seen that the HF SE of the ion and its next ionization stage compensates with each other, so that the whole ionic contribution is negligible in comparison to the electronic contribution. Even in the case of low density for , the ionic contribution amounts to of the electronic one. Therefore, in most case only the electronic contribution, i.e. Eq. (42), has to be taken into account for the determination of IPD in plasmas.

#### iii.2.2 Energy shift of bound states

The impact of the statistical correlation is not only reflected in the reduction of the continuum edge but also in the shift of discrete energy levels . We assume that the optical electron are ionized from the outermost shell of ion . Therefore, in this section the subscript is used to denote the bound electron in the outermost shell of the ground state of charge state . The shift consists of the Fock shift Roepke19 ()

(44) |

and the Pauli blocking shift Roepke19 ()

(45) |

The bound electron is assumed to move in the mean field produced by the nucleus and other bound electrons. The effective charge number experienced by a bound electron in the outermost shell of charge state in the isolated case is described by . It can be calculated within the screened hydrogenic model FBR08 (), where a many-electron atom/ion is approximated by a hydrogen-like system. Assuming that the selected bound electron occupies the -like state of the hydrogen-like system with effective core charge number , then the corresponding wave function reads

(46) |

with . In the momentum space the wave function is given by

(47) |

Inserting this wave function into the Fock shift (44) as well as into the Pauli blocking term (45), we obtain

(48) |

and

(49) |

Gathering both contribution and referring it as Pauli-Fock contribution, we arrive at

(50) |

with and

(51) |

It can be demonstrated that the bound electrons have an important influence on the physical properties, in particular impacted by the Pauli blocking effect in strongly degenerate plasmas Roepke19 (); Lin19 (). A more detailed description demands a systematic investigation of the internal structure and the knowledge of the interaction between the bound electrons in complex many-electron systems, which is not intended in the present work.

### iii.3 IPD due to dynamical correlation: real part of dynamical interaction contribution

In this section we discuss the dynamical interaction contribution of the GIPD, where the charge-charge dynamical SF are introduced according to the fluctuation-dissipation theorem to describe the influence of the plasma environment on the investigated atomic/ionic system. The IPD in plasmas is proved to be directly determined by the spatial arrangement and the temporal fluctuation of surrounding particles.

#### iii.3.1 Fluctuation-dissipation theorem

The response of an interactive plasma system to external perturbations are totally determined by the dielectric function, which contains the complete information on the ions and the free electrons in this interacting system. It is directly connected to the density-density response function between particles of species and via the following relation KSKB05 (); GR09 ()

(52) |

which describes the induced density fluctuations of species owing to the influence of an external field on particles of species . Additionally, the detailed spatial and temporal structure of a many-body system is elaborately described by its density-density dynamical SF. Such dynamical SF determines many transport and optical properties, such as stopping power, the equation of state, the spectral lines, IPD and the corresponding ionization balance. Taking into account the fact that the dynamical SFs are related to their corresponding density-density correlation functions via Fourier transformation, the partial density-density dynamical SF for different components in plasmas can be defined in terms of the density-density response function via the following expression Hoell07 ()

(53) |

Therefore, for a multi-component plasma the fluctuation-dissipation theorem can be described by means of the following relation

(54) |

The free-bound dynamical SF accounting for the correlation between bound and free electrons and the electron-electron dynamical SF are related to the ionic dynamical SF and the dynamical SF of fast-moving free electrons Chihara00 (); WVGG11 ()

(55) | ||||

(56) |

Within the framework of the linear response theory the screening function can be expressed in terms of the electronic dielectric function Chihara00 (); GVWG10 (); Chapman15 ()

(57) |

Taking electron-ion interactions to be Coulomb interaction potentials, the long-wavelength limit of the screening function within the linear response is given by GVWG10 (); Chapman15 ()

(58) |

which is proportional to the charge number of the test particle . In the high density regime, such long-wavelength approximation hidden in the RPA dielectric function might be inapplicable to describe the finite-wavelength screening, so that the full version of RPA dielectric function has to be utilized in the calculation of the screening function Chapman15 ().

Consequently, the fluctuation-dissipation theorem can be recast in terms of the total charge-charge dynamical SF . The effective charge-charge response of a multi-component charged plasma to an immersed impurity is then described by

(59) |

with respect to the total inverse screening length . Detailed explanation of this expression is given in the Appendix B. The inverse screening parameter is exhaustively discussed in the next subsection (also see the Appendix C). The total charge-charge dynamical SF is given by

(60) | ||||

For the screening function , the long-wavelength approximation, i.e. Eq. (58), will be used in the following calculations. Introducing

(61) |

the effective charge-charge dynamical SF for the total many-body system can be then rewritten as

(62) |

with the ionic charge-charge dynamical SF

(63) |

#### iii.3.2 Non-linear screening effect: effective inverse screening length

Assuming that the electrostatic potential takes the following form

(64) |

according to the Boltzmann distribution KKER86 (), the mean particle density of species in the vicinity of the test particle with charge number is described by

(65) |

In a classical plasma, the inverse screening parameter in Eq. (64) is reasonably described by the inverse Debye length with and . In a strongly coupled, non-ideal plasma, the electronic screening is determined by the Thomas-Fermi length. Similarly, the ionic screening in high density plasmas can not be described by the Debye length any more, since non-linear effect relating to strong ionic coupling has to be taken into account. In the following, we determine the effective screening parameter according to the perfect screening rule (charge neutrality) Nordholm84 ()

(66) |

where indicates the charge number of a particular ion after the ionization process has taken place, i.e. is the spectroscopic symbol. The screening cloud is assumed to be spherically symmetric and is defined as

(67) |

Inserting this expression into the Eq. (66), an effective inverse screening parameter can be determined, if the detailed charge state distribution is known. However, such procedure is very complicated to be performed, since integration over a series of transcendental functions has to be worked out. Regarding the ionic mixture effectively as a single ionic species with charge , the screening cloud can be approximated by the following expression

(68) |

The derivation of this approximation is given in detail in Appendix C. Introducing the following dimensionless variables

and the impurity-perturber coupling strength

(69) |

with respect to the ionic radius

(70) |

we obtain the following closed equation for the effective screening parameter in terms of the impurity-perturber coupling strength