# Classical Molecular dynamic codes for hot dense plasmas: The BinGo code suite

Annette Calisti, Sandrine Ferri, Caroline Mossé, Bernard Talin

► To cite this version:

Annette Calisti, Sandrine Ferri, Caroline Mossé, Bernard Talin. Classical Molecular dynamic codes for hot dense plasmas: The BinGo code suite. High Energy Density Physics, 2024, 10.1016/j.hedp.2024.101084. hal-04454313

HAL Id: hal-04454313

https://hal.science/hal-04454313v1

Submitted on 9 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Classical Molecular dynamic codes for hot dense plasmas: the BinGo code suite.

A. Calisti<sup>a</sup>, S. Ferri<sup>a</sup>, C. Mossé<sup>a</sup>, B. Talin<sup>a</sup>

<sup>a</sup>Aix Marseille Université, CNRS, PIIM UMR 7345, 13397, Marseille, France.

## Abstract

The purpose of this paper is to illustrate our contribution to general plasma physics studies obtained since the 90s with multiple versions and adaptations of the classical molecular dynamics (CMD) simulation interactive code called BinGo. After a description of the particulars of the CMD simulation models and the BinGo code suite, some applications are discussed for illustration. These results validate the CMD simulation as a powerful tool of investigation for hot dense plasmas.

**Keywords:** Classical Molecular Dynamics, Simulations, Dense plasmas

## 1. Introduction

Classical molecular dynamics (CMD) was born years ago with computers [1, 2]. Its primary motivation was to provide alternative statistical results (numerical experiments) and at the same time, overcome complexity (N-body problem) in statistical physics calculations. These developments resulted into helpful discussions and comparisons. Nowadays, the motivation remains the same but benefits from the extraordinary growing capabilities of computers giving to CMD the status of a real investigation tool. Due to the interaction length between charged particles, CMD is well adapted to plasmas for a wide domain of composition, temperature and density conditions. By going from the simulation of one-component plasmas to simulations of electron-ion two-component plasmas, CMD is an unavoidable tool with a great versatility which allows one to consider various applications in plasmas.

The BinGo code suite is one of these CMD codes. It has been developed for more than two decades in order to accompany studies of plasma spectroscopic diagnostics. It has been used not only to generate the field histories required for the simulation of spectral line shapes in plasmas but it has

also been very useful for various studies of particle correlation effects on plasma statistical properties. Among these, the effects of particle correlations on the broadening due to the radiator motion (Doppler broadening) in the framework of a detailed analysis of the different processes that contribute to the spectral broadening of the Ni-like Ag XUV laser line, have been studied using BinGo applied to a one-component Yukawa plasma [3]. The study of electron dynamics around an ion impurity provided results for comparisons with theoretical modeling of field correlation and distribution functions measured at an ion [4, 5]. The CMD code enabled us to explore micro field statistical properties in the context of line shape studies for plasma diagnostic by spectroscopy in both cases, an impurity (or a small percentage) of hydrogen-like helium in protons [6] and neutral hydrogen in an infinite system of interacting electrons and protons [7]. The most recent investigations have been performed with BinGo-TCP, a version of BinGo designed for two component (TCP) electron ion plasma simulations. This version led us to explore warm-dense and hot-dense plasma conditions with simulations of dynamics structure factors [8] and the evaluation of the mechanism of ionization potential depression in dense plasmas [9], respectively.

The purpose of this paper is not to make a review article but to illustrate our contribution to general plasma physics studies obtained since the 90s with

Email address: annette.calisti@univ-amu.fr (A. Calisti)

Preprint submitted to High Energy Density Physics

January 25, 2024
multiple versions and adaptations of the interactive
CMD code called BinGo. In the following, after a
section devoted to give the particulars of the CMD
simulation models and the BinGo code suite, we will
come back on some of the previously cited applica-
tions for illustration. These advances and results
validate the CMD simulation as a powerful investi-
gation tool for hot dense plasmas.

## 2. CMD basis

### 2.1. Generalities

Classical molecular dynamics is a technique
shared by many physicists and its detailed descrip-
tion is not needed here (for more details see for
example [10] and references there in). It consists
in the simulation of the movement of interacting
atoms or molecules treated as classical non rela-
tivistic pointlike particles. In its standard form,
CMD deals with a three dimensional infinite sys-
tem in thermodynamical equilibrium composed of
joint identical cubic cells. All the step by step cal-
culation work is done into a unique cubic cell with
the so called minimum image convention [11] im-
plying that the acceleration on a given particle is
calculated considering a box centered on this par-
ticle. The number of interacting particles consid-
ered at a time and moving according to Newton's
laws integrated using for instance the velocity Ver-
let algorithm [12], is constrained by computer ca-
pabilities. Different adaptations of this basic CMD
have been developed. A non exhaustive list of ex-
amples is given here after for illustration: mirror
walls approximation for the boundaries of the sim-
ulated system [13], rectangular periodic boxes, non
periodic finite systems [14, 15], systems with an in-
terface [16, 17], etc.

The standard CMD simulations are performed
for thermodynamic state characterized by a fixed
number of atoms, *N*, a fixed volume, *V*, and a
fixed energy, *E*. These extensive variables (*N*,
*V*, *E*) of the microcanonical ensemble should be
strictly conserved (i.e., time-independent) during
the run of the simulation. The corresponding in-
tensive variables, namely the chemical potential,
*ν*, the pressure, *P*, and the temperature, *T*, are
not conserved. In a simulation at equilibrium, the
corresponding instantaneous observables fluctuate
around well-defined average values and in a non-
equilibrium simulation, these quantities may un-
dergo a systematic drift. For studies dedicated to

equilibrated system, the preparation, i.e., the ther-
modynamical equilibration at a given temperature
of an assembly of interacting particles is a crucial
step to be performed before any further relevant
statistical samplings. In general, tricks like ther-
malization are used to accelerate this preparation
following practical considerations. A number of
methods have been developed to keep the temper-
ature constant while using the microcanonical en-
semble, e.g. velocity rescaling, the Andersen ther-
mostat, the Nosé-Hoover thermostat, Nosé-Hoover
chains, the Berendsen thermostat or Langevin dy-
namics [18]. A system will be considered as equi-
librated, if the intensive variables do not show any
drift and if the total energy of the particles in
the box undergoes negligible numerical fluctuations
that do not affect the static or dynamic measure-
ments being performed when the system evolves.

### 2.2. Interaction potentials

A molecular dynamics simulation requires the
definition of the potential function that governs in-
teractions between the particles inside the simula-
tion box. This issue will be considered differently
according to the studied plasma conditions.

#### 2.2.1. The one component plasma (OCP and Yukawa OCP) systems

The one component plasma (OCP) model, which
consists in dense systems of charged particles of
one species interacting through a Coulomb poten-
tial, embedded in a uniform background of oppo-
site charge ensuring over-all electrical neutrality, is
of great astrophysical interest. It provides an ex-
cellent model for describing superdense, completely
ionized matter typical of white dwarfs or the outer
layers of neutron stars [19]. This model is valid
mainly when the plasma conditions are such as the
electrons can be considered as a highly degener-
ate electron gas and then as a rigid uniform back-
ground. When it is necessary to account for the
electron gas polarization by the ionic charge distri-
bution, the bare Coulomb potential is replaced by
an effective screened potential [20]. This is called
one-component Yukawa systems. Depending of the
degree of electron degeneracy, the screening length
will be chosen differently, e.g. the Debye length for
ions in a pure classical electron gas or the Thomas-
Fermi screening length in the case of partially de-
generate electrons.

In all cases, the interaction length in CMD is lim-
ited, by construction, to the box size *a*, defined as
$a = (N/n)^{1/3}$, with $n$ the density of charges. For pure Coulomb potential or if the interaction length is larger than the simulation box, it is necessary to account for particles outside the simulation box. In these cases, different methods can be used such as Ewald sum, particle mesh Ewald (PME) or particle-particle particle-mesh (P3M) techniques [21].

### 2.2.2. The two component ion and electron plasma systems

A step further towards more realistic plasma models is done simulating ions and electrons together provided that an ion-electron potential finite at short distance prevent electron-ion collapse. Such effective regularized potentials are designed to approximately account for known quantum properties appropriate for the kind of investigated plasma, e.g. quantum interference and diffraction effects at small distances involved in electron-electron or electron-ion collisions for hydrogenic systems [22] or the correct ionization energy of the ions in the ground state when the electrons are situated on top of the ions for nano plasma properties [23]. The computer cost of these simulations is important as the time steps have to be short enough for a good description of electron motion and the total number of time steps has to be large enough to allow ions to move significantly in order to build, self consistently, their spatial structure. However, they enable us to study electron dynamics, which is the key to understanding a number of plasma mechanisms.

As mentioned above, the TCP simulations can be performed provided that an effective regularized ion-electron potential is used. The choice of this potential can affect the measurement of certain properties, in particular those that are strongly dependent on quantum effects. Since the pioneering work of Hansen and McDonald in the 1980s [24] who first used such potentials in classical MD simulations, a great deal of work has been devoted to determining the effective potentials describing classical systems composed of ionic charges and free electrons. The classical approach is guided by an analysis of the quantum diffraction of an electron by a ionic charge allowing to determine an equivalent classical regularized potential. This potential depends on temperature and because the regularization a potential energy minimum becomes effective for an electron at null distance from the ion. Several regularized potentials have been developed by different authors [25–29] in which a common regularization distance used is the DeBroglie wave length.

A widely used regularized electron-ion potential, i.e., finite at short distances is defined as [26]:

$$V_{ie}(r) = -\frac{Z_i e^2}{r} \left(1 - \exp\left(-\frac{r}{\delta_{ie}}\right)\right), \quad (1)$$

where the regularization distance is the DeBroglie wave length of the relative electron motion, $\delta_{ie} = \hbar/(2\pi\mu_{ie}k_B T)^{1/2}$ with $\mu_{ie}$ the ion-electron reduced mass. Here, the spatial distribution of the bound electrons is not considered, their charges are included in $Z_i$. Thus, it is important to make sure that in the MD simulations most of the electrons are free. The simulated electrons occupying the quasi-periodical orbits with $r \ll r_0$, $r_0$ being the inter ionic distance, around individual ions can be considered bound. It has been shown in [13] that the plasma conditions must be selected such that the regularization distance $\delta_{ie}$ is large enough to suppress the formation of the classical bound states of electrons while satisfying $\delta_{ie} \lesssim r_0$ so not to affect the free electron density at $r \sim r_0$. An estimation of the probability to find an electron within the volume $r \lesssim \delta_{ie}$ shows that for $\Gamma_{ie} \gtrsim 0.5$, with $\Gamma_{ie} = \frac{Z_i e^2}{k_B T r_0}$ the ion-electron coupling parameter, the bound electrons represent a significant fraction of all simulated electrons. Thus, the MD simulations clearly do not correspond to the formulation of the initial problem, since the simulated $Z_i$ and the electronic density are effectively reduced. Furthermore, for plasma conditions such that bound electrons exist, the simulation equilibrium state is difficult to reach. The filling of negative energy states by electron trapping, i.e. the adjustment of electrons to ions, results from three-body collisions involving one ion and two electrons, which are rare events. Nevertheless, this preliminary phase is necessary as at start the ions and electrons are put together at random in the simulation box and the ion-electron interaction is switched on inducing a non equilibrium state of the ion-electron plasma. Because the pairwise exchange of kinetic and potential energy, this non equilibrium state, would evolve towards a state far from the density temperature conditions expected regarding the objective of the simulation.

Ion-ion and electron-electron potentials are Coulomb potentials,

$$V_{ii,ee}(r) = \frac{Z_{i,e}^2 e^2}{r}. \quad (2)$$

Accounting in some way for electron degeneracy can be achieved using an electron-electron poten-
tial finite at short distances to account for quantum diffraction effects plus an effective Pauli potential with symmetry effects [30]:

$$V_{ee}(r) = \frac{e^2}{r} \left(1 - \exp\left(-\frac{r}{\delta_{ee}}\right)\right) + k_B T(\ln 2) \exp\left(-\frac{r^2}{\pi \delta_{ee}^2 (\ln 2)}\right). \quad (3)$$

The domain of application of these pair-wise potentials is limited, for a given temperature, to not too high densities as they do not treat many-body effects. For an analysis of this issue and of the techniques to model many-body effects see [31] and references therein. Implementing those potentials in MD leads necessarily to the Maxwell momentum distribution and not to the momentum distribution of an electron gas governed by Fermi statistics. This may be corrected by using a momentum-dependent potential [32–34].

### 2.2.3. Virtual neutral plasma systems

For various purposes, it is necessary to simulate virtual plasmas composed of electrons and ions in different ionization states. Depending on the plasma conditions, e.g. when they are such that the electron-ion interactions are strong or those generated by direct target interaction with intense laser beams, it is required to account, in the simulations, for bound states as well as quantum mechanical processes to create and destroy them. This permits to avoid non-physical Coulomb collapse and to extend the applicability of classical MD simulations, which are ideally suited for fully-ionized, non-degenerate plasma, to a widest range of plasma conditions. Various approaches have been developed and proposed by different groups, to perform molecular dynamics simulations in such systems. A powerful method is the quantum molecular dynamics (QMD), in which the dynamics of the ions is solved through classical MD techniques according to the electronic structure obtained at each time step using density functional methods [35]. This method is widely used to generate benchmark data but it is very computationally expensive and it is limited to low temperature plasma conditions. To extend QMD to higher temperatures, orbital free density functional theory molecular dynamics (OF-DFT-MD) simulations have been developed. In these simulations based on the Thomas-Fermi approximation, the electronic fluid is represented by a free energy entirely determined by the local density [36, 37]. Both these methods use the

Born-Oppenheimer approximation, i.e., they suppose that the electrons adjust instantaneously to ions, and do not solve explicitly the electron dynamics. Moreover, even though the OF-DFT-MD simulation is more efficient than QMD in terms of computational time, it still does not allow to simulate a large number of particles over a long period of time. This has been a strong motivation for developing alternative methods. Among these, it is worthwhile to mention the ddcMD code which has been designed to simulate high energy density plasmas undergoing thermonuclear burn and which is the backbone of the Cimarron project based at Lawrence Livermore National Laboratory [38]. This code is a massively parallel molecular dynamics code in which quantum processes such as emission and absorption of X-ray photons, thermonuclear fusion or collisional ionization and recombination have been incorporated by using the "Small Ball" (SB) method [39]. The latter combines a classical MD technique as defined previously to solve the many-body classical dynamics when interparticle distances are greater than the SB radius and a Monte Carlo technique to solve the two-body quantum mechanics governed by probabilities for interparticle distances lower than the SB radius.

Another example concerns studies on cluster dynamics for which TCP ion-electron MD simulation incorporates some quantum mechanical effects through the interaction pseudopotential between ions and electrons [23]. In these simulations, ions and electrons are treated as classical particles represented by a Gaussian of fixed width leading to a inter-particle potential of the form:

$$V_{ij}(r) = Z_i Z_j e^2 \frac{\text{erf}(r/\sigma_{ij})}{r}, \quad (4)$$

where $\sigma_{ij}$ is given by $\sigma_{ij}^2 = \sigma_i^2 + \sigma_j^2$ with $\sigma_i$ the Gaussian width of the particle $i$. Depending on the choice of $\sigma_i$, core and valence electrons can then be treated at the same level in the simulations. For example, in [15] devoted to a study of laser excited small sodium clusters, $\sigma_i$ is chosen in order to reproduce the ionization energy for sodium in the ground state when the electrons are situated on top of the ions.

In our BinGo-TCP model, we have also chosen to simulate ions and electrons as classical particles and to incorporate the ionization potentials of isolated ions as a minimum quantum information in a regularized potential to model ionization and collisional recombination processes. The BinGo-TCP code is
designed to deal with neutral mixtures composed of ions of the same atom with different charge states and electrons. Inside the simulation box the ion charges can change from one to another according to the density temperature conditions. For that purpose, a regularized electron-ion potential, depending on the ion charge $Z_i$ is defined as:

$$V_{ie}(r) = -Z_i e^2 e^{-r/\lambda}(1 - e^{-r/\delta(Z_i)})/r, \quad (5)$$

where the regularization distance $\delta(Z_i)$ is chosen to reproduce the ionization energy $E_i$ of the unperturbed ion of charge $Z_i$ in the ground state when the electron is located at the ion ($r=0$),

$$\delta(Z_i) = Z_i e^2 / E_i. \qquad (6)$$

$\delta(Z_i)$ also plays the role of radius of the ion of charge $Z_i$. The screening factor present in these potentials $e^{-r/\lambda}$ where $\lambda$ is half the simulation box size, helps to smooth the small fluctuations of forces arising with the periodic boundary conditions. It could affect the mechanisms controlling the particle motion in the simulation box and thus, the measure of some plasma properties such as non-equilibrium properties. Concerning the examples of applications hereafter, it has been checked that the results do not depend on the choice of $\lambda$ providing that the box size is large enough (a few times the natural plasma screening length).

The choice of this potential associated to the knowledge of the position and velocity of individual particles at each time step, allows us to design a collisional ionization - recombination process. Ionization-recombination mechanisms rely on an approximate analysis of collisional events between one ion and 1 or 2 electrons. In dense plasmas, the concept of collisions is not straightforward as the interaction involves all particles within the screening length. The main idea, in our model, is to extract from the data provided by the MD simulation (positions and velocities of particles) a local characterization of the plasma around a given ion "A" in order to infer if the conditions are favorable to a ionization or recombination of this ion. For that purpose, notions of hot shell and cold shell of "A" have been introduced. The first mutual nearest neighbor (MNN) and the next nearest neighbor (NNN) electrons of "A" are first identified. According to the location of the MNN electron - it must be located in between $\delta(Z_i)$ and $\sqrt{2} \times \delta(Z_i)$ - and the sign of the total energy of the MNN and NNN electrons - calculated accounting for all the interactions in the simulation box -, one then evaluates if

locally the plasma, at one step of its evolution, is favorable to the ionization (positive energy) or the recombination (negative energy) of "A". This test results into a pre-ionization, i.e., an increase by 1 of the ion charge and the addition of one electron located at the ion, or a recombination, i.e., a decrease by 1 of the ion charge while the MNN electron is removed. This local discontinuity over one time step is then accounted for by the whole system through a normal evolution. The pre-ionized state, i.e., an ion with a trapped electron can be converted into an ionized state through multiple collisions. In this approach the ionization will be considered as completed when a new hot shell surrounds the ion. In the mean time the ion is considered as excited or multi-excited if there are more than one trapped electrons in the ion potential.

At start, ions of same charge and electrons are put together in the simulation box and the ion-electron interactions are switched on inducing a state of the ion-electron plasma far from equilibrium. The ionization/recombination process allows the evolution of the charge state population towards a stationary state depending on temperature, density and composition of the plasma and it favors the setting up of a population of electrons temporary trapped in the ion wells. During this phase of equilibration, the electron temperature is constrained by velocity rescaling while the ionization recombination process broadens the charge distribution. During this short phase, no significant change of ion temperature and positions can occur. At the end of the preparation phase, the system follows a quasi stable evolution with stationary trapped and free electron populations and is ready for extracting data samplings dedicated to plasma investigations.

The equilibration of the charge states is illustrated in Figure 1, for a Carbon plasma at solid density and 50 eV. The balancing of carbon charge state populations is presented for two simulations with different initial conditions. The solid lines show the results obtained for a simulation with all charges equal to z=6 at the initial time and the dashed lines show the same results but for z=1 at t=0. We can see that thanks to the ionisation-recombination process, whatever the choice of initial state, after a few femto-seconds the concentrations of charge states stabilise around average values. Once this state of equilibrium is reached, each ion continues to undergo ionisation or recombination, but only infrequently, and the charge distribution no longer changes. It has been shown that, in
<figure><img src="image_2.png" /><figcaption>Figure 1: Time evolution of the concentration of charge states in a Carbon plasma at solid density and 50 eV during and after the equilibration phase of the simulation. Full lines: the initial ionic charges are set at z=6, Dashed lines: the initial ionic charges are set at z=1.</figcaption></figure>

general, these simulated charge distributions compare well with those calculated using the FLYCHK code [40] for the same plasma conditions.

Our aim in developing the BinGo-TCP code was to avoid excess of complexity of models in order to preserve a straightforward interpretation of physics from simulation results. One gains the ability to describe the ion-electron coupling accounting for mixtures of ions undergoing changes of their charge states. Note that, the coupling of electrons with radiation is ignored, the notion of discrete energy for the ionic excited states is here replaced by its continuous equivalent and that, the lowest energy of an electron depends on all the charges including the closest ion.

### 3. Example of applications

#### 3.1. One component Yukawa systems

##### 3.1.1. Effects of particle correlations on spectral line shapes.

The intrinsic spectral profile of a given line in a plasma is determined predominantly by spontaneous emission rates, electron collisional rates, Stark broadening and Doppler broadening [41] with possible complications due to ion turbulence [42] and additional ion-ion interactions [42, 43]. In medium with gain, the observed profile is modified by radiative transport effects in being narrowed approximately as the square root of the gain-length product in the small signal regime [44]. As the laser saturates, if the intrinsic profile is dominated by inhomogeneous rather than homogeneous broadening

mechanisms, the line can be re-broaden to its intrinsic width. This points out clearly the importance to have a good representation of the intrinsic profile together with a good understanding of the different mechanisms responsible for broadening.

It is usual to consider that the inhomogeneous broadening mechanisms are caused by the local inhomogeneities of the medium such as Doppler shifts, quasi-static electric micro fields or turbulence and that homogeneous broadenings are mainly due to electronic collision and/or spontaneous emission. Nevertheless, in the same manner as the quasi-static electric field approximation can be inappropriate in the Stark broadening theory due to ion dynamics [45], the collision-free formalism involved in the standard Doppler effect calculation could fail if velocities change over time scales of the same order or shorter than the effective radiative lifetime of the oscillator (i.e. the inverse of homogeneous spectral linewidth). In other words, if the velocities change before the light emission happens, it is no longer possible to consider ions on straight trajectories and a calculation, taking into account collisions, has to be done. In some circumstances, the breakdown of the collision-free approximation can result in the effective "homogenization" of the ordinarily inhomogeneous Doppler profile by collisional redistribution when the mean velocity-changing collision interval $t_c$ is less than the effective radiative lifetime or in the collisional narrowing of the Doppler profile when $t_c$ is shorter than the effective Doppler correlation time (Dicke narrowing [46]).

A detailed analysis of the different broadening mechanisms of the spectral profiles of the $4d-4p$ ($J=0-1$) lasing line in Ni-like Ag ($\lambda=13.9$ nm) has been performed by using the PPP line shape code [47] and the BinGo code in its Yukawa OCP version. The plasma temperatures and densities that are consistent with the ones required for collisional excitation pumping of Ni-like Ag in laser-produced plasmas, give rise to strong coupling plasma regime where correlations between particles could no longer be ignored in the case of transient XUV lasers. The electron temperature can vary between 200 and 700 eV depending on the parameters of the pulse of the laser used to create the plasma and the ionic temperature is in between 20 and 50 eV for transient XUV laser and of order of 200 eV for QSS XUV lasers.

A study of the accuracy of the free-particle Doppler approximation, versus densities and tem-
peratures has been done by using CMD simulations, for the Ni-like Ag laser line and for conditions relevant to transient and quasi-steady-state (QSS) XUV lasers. Significant differences in the behavior of intrinsic line profiles between these two laser regimes have been demonstrated [3].

If one accounts for the emitter motion, the general expression of the line profile reads [45]:

$$I(\omega) = \Re e \frac{1}{\pi} \int_{0}^{\infty} dt e^{i\omega t} \times \\ < e^{i(\vec{k} \cdot \vec{r}(t) - \vec{k} \cdot \vec{r}(0))} \vec{d}(t) \cdot \vec{d}(0) > \quad (7)$$

where <<>> denotes an ensemble average over the emitter plus plasma system, $\vec{d}$ is the radiator dipole operator and $k = 2\pi/\lambda$ with $\lambda$ the wavelength of the considered line. The factors $e^{\pm i \vec{k} \cdot \vec{r}}$ account for the radiator's center-of-mass motion. Broadenings due to the interaction of the emitting ion with surrounding particles and to emitter motion are statistically dependent in the general case. Broadening due to interactions results from a modification of the internal state of the atomic oscillator. Both this internal state and the velocity of translational motion of the emitter can be altered in the same collision.

In this study, interactions with the electronic component of the plasma dominate, giving rise to a phase shift of the atomic oscillator. This is due to electronic collisions which change substantially the phase without altering the velocity of the emitter owing to the great difference of masses. So, it is quite accurate to ignore correlations between the ion translation $\vec{r}(t)$ and the dipole moment $\vec{d}(t)$:

$$I(\omega) = \Re e \frac{1}{\pi} \int_{0}^{\infty} dt e^{i\omega t} \times \\ < e^{i(\vec{k} \cdot \vec{r}(t) - \vec{k} \cdot \vec{r}(0))} > < \vec{d}(t) \cdot \vec{d}(0) > . \quad (8)$$

The line shape appears as the Fourier-transformed of a product of two correlation functions, the radiator dipole operator correlation function, $C(t) = <\vec{d}(t) \cdot \vec{d}(0)>$, and the self-structure factor, $S_s(k, t) = <e^{i(\vec{k} \cdot \vec{r}(t) - \vec{k} \cdot \vec{r}(0))}>$ which can be calculated independently. The line shape is then the convolution of the profile due to emitter interactions with the plasma, $I_{int}(\omega)$ and the profile due to emitter motion, $I_D(\omega)$.

In [3], the PPP line shape code, a multi-electron radiator line broadening code developed to calcu-

late theoretical spectral line profiles for a general emitter in a plasma, is used to calculate $I_{int}(\omega)$.

The self-structure factor, $S_s(k, t)$ is well known in the free-particle limit resulting from the hypothesis that each radiating ion moves at constant velocity $\vec{r}(t) = \vec{v}t$ with a Maxwellian distribution of velocities, and is given by:

$$S_s(k, t) = e^{-k^2 t^2 / 2\beta m}, \quad (9)$$

with $\beta = 1/k_B T$ and $m$ the ion mass. A straightforward way to take into account interactions between ions in the calculation of the line shape, is to use a classical molecular dynamics simulation techniques (MD) to compute $S_s(k, t)$. The plasma model consists of classical point ions interacting together through a coulombic potential screened by electrons and localized in a cubic box of side $L$ with periodic boundary conditions. Newton's equations of particle motion are integrated by using a velocity-Verlet algorithm using a time-step consistent with energy conservation. Due to periodic boundary conditions, $k = 2\pi/\lambda$ must satisfy $k_{x,y,z} = n_{x,y,z}2\pi/L$, $n_{x,y,z}$ being an integer number. The number of particles, $N$ (thus $L$), is chosen to find $k$ as close as possible to that of the considered laser line.

Integrating the Newton's equation gives access to the positions and velocities of the ions as a function of time and thus to the associated static and dynamic statistical properties such as structure factors, velocity correlation functions, diffusion coefficients, ion-ion collision rates etc.

Calculations of spectral $4d-4p$ line Doppler profiles accounting for ion correlations for different electronic densities, in cases of transient and QSS XUV lasers with ($Ti = 20eV, Te = 200eV$) and ($Ti = Te = 200eV$), respectively, have been performed.

In the cases of transient XUV lasers, where the ionic temperature is relatively low and the electronic densities are in between $N_e = 5 \times 10^{19}$ and $N_e = 7 \times 10^{20}$ cm$^{-3}$, the plasma coupling parameter $\Gamma = Z^2e^2/(r_0kT_i)$, $r_0$ being the ion sphere radius, $r_0 = (3/(4\pi N_i))^{1/3}$, takes values from ~6 to ~14. The corresponding plasmas are strongly coupled making the concept of binary collision between ions compromised because the ions are always in interaction. Collective effects (multiple collisions) are expected on the profiles. When accounting for ionic correlations, the line profiles are not only narrowed but a structure appears in the wings more or
less far from the center depending on the density. The profile is described in [43, 46, 48] as a narrow Lorentzian central peak superimposed on a broader Gaussian plateau. We will see in the following that it is not so simple. The frequencies, $\Delta\omega_{osc}$, delimiting the two behaviors correspond to the frequencies of oscillations of the velocity autocorrelation functions, $C_v(t) = \langle \vec{v}(t) \cdot \vec{v}(0) \rangle$ obtained by MD simulations.

The spectral line shapes carry the mark of the plasma oscillations together with the multiple collisions which change velocities. To some extend, it is possible to say that for frequencies greater than $\Delta\omega_{osc}$ (for times shorter than the time necessary for the velocity to change sign), the velocity is constant and the profile is given by the Doppler free particle limit while for frequencies smaller than $\Delta\omega_{osc}$, the time of interest is large enough for having a noticeable changing in the velocity direction resulting in a strong narrowing of the profile. The description of the narrow central peak is rather complicated in the cases of interest here because the concept of collision fails in those strongly coupled plasmas. It has been shown that even though the two limits, free-particle and diffusion limit, seem to be applicable for low and high densities respectively, the profiles are never completely Gaussian nor Lorentzian and can be well fitted by a Voigt profile.

In the case of QSS XUV laser study, the plasma coupling parameter falls in between 0.7 and 1.4 due to a higher ionic temperature. One can expect that the notion of collision is meaningful and that collisional models are an alternative to simulation. A first noticeable difference comes from the $C_v(t)$ function which, in all the cases here, does not oscillate and can be described by an exponential decay function. The effects of narrowing are less important than for the transient laser nevertheless the shape of the line profiles is modified. Due to the fact that Doppler effect leads to much more broader lines, the plasma frequency $\omega_{pi}$, which is a marker of collective effects is always inferior to the Doppler width ($\Delta\omega_D \sim 9.37 \times 10^{-3}$), so no structure appears on the line profile, but here again the central part of the profile is neither lorentzian nor gaussian.

Supposing that the concept of collisions makes sense here, we have deduced from the diffusion constant obtained by MD simulations through the relation $D = \frac{1}{3} \int_{0}^{\infty} C_v(t)dt$, an effective mean free path, $l = 3\Gamma^{1/2}D^*r_0$ with $D^* = D/\omega_{pi}r_0^2$. This permits to define a collision frequency with $1/\tau_{coll} = \bar{v}/l$

<figure><img src="image_3.png" /><figcaption>Figure 2: Comparisons of the 4d-4p laser line profiles obtained by MD simulation (red circles) and the strong collision model (blue line), for $N_e = 7 \times 10^{20} \text{cm}^{-3}$ and $T_i = T_e = 200$ eV. The black dash line represents the free particle limit.</figcaption></figure>

where $\bar{v} = \sqrt{2kT_i/m}$. These rates have been used in a strong collision model according to [43]. In this model, it is assumed that the modifications of the profiles are due to collisions such as the velocity $v$ of a particle after collision is independent of its velocity $v'$ before collision. According to this, the effect of velocity changing collisions is modeled by a stationary Markov process giving rise to the following expression for the line profile:

$$I(\omega) = \frac{1}{\pi} \Re e \frac{\int \frac{W(v)dv}{\nu_{coll} + i(\omega - \bar{k} \cdot \vec{v})}}{1 - \nu_{coll} \int \frac{W(v)dv}{\nu_{coll} + i(\omega - \bar{k} \cdot \vec{v})}}, \quad (10)$$

with $W(v)$ the Maxwellian distribution function of velocities and $\nu_{coll}$ the velocity changing rate. Note that the profile is just a function of the profile given by the Doppler free-particle limit and a unique collision rate. Figure 2 is the comparison of the profile obtained by MD simulation with the result of eq.10 for $N_e = 7 \times 10^{20} \text{cm}^{-3}$ and $T_i = T_e = 200$ eV.

The same excellent agreement have been obtained for all the considered densities.

## 3.2. Electron and ion plasmas

### 3.2.1. Microfield statistical properties.

In plasmas, the Stark broadening of spectral lines results from the interaction of the emitter's internal degrees of freedom with the ionic and electronic micro-fields created by the plasma. Modeling this broadening is a complex problem involving a complex combination of atomic physics, statistical physics and plasma physics, the most difficult problem being to identify the emitter environment completely and correctly. In particular, taking into
account fluctuations in the electric fields produced on the emitter by moving ions and electrons has been the subject of constant interest since the 1960s [45]. With the exception of a few well-known cases, the calculation of line profiles in plasmas requires the use of numerical codes. There are several codes of this type, of varying degrees of complexity, which necessarily differ in terms of scope and accuracy [49]. Among these codes, numerical simulation is playing an increasingly important and unique role. In numerical simulations, the Schrödinger equation describing the time evolution of the emitter wave functions in the presence of a time-dependent electric field is solved and then averaged over a statistically representative set of electric fields to obtain the line profile.

Today, thanks to advances in computer technology, the electrical microfields derived from force calculations on individual charges can be simulated using classical molecular dynamics (CMD) techniques applied to plasmas with two components, ions and electrons, where all the interactions between the charges are taken into account in the motion of the particles [7, 50–53].

CMD provides the stochastic functions of the microfields needed to solve the emitter evolution equation. As an example, Figure 3 shows the electric micro field components, $E_{x,y,z}(t)$ measured at a neutral Hydrogen in a electron-proton plasma at $n_e = 10^{18}$ cm$^{-3}$ and $T_e = T_i = 1$ eV, simulated by TCP CMD. The time history in second of the three components of the total field, $\vec{E}(t) = \vec{E}_i(t) + \vec{E}_e(t)$ have been plotted together with the components of the ionic, $\vec{E}_i(t)$, and the electronic, $\vec{E}_e(t)$, fields. The differences in the dynamics of ions and electrons due to their mass differences, are clearly seen.

The spectral line shapes obtained by using simulated electric fields to solve the emitter evolution equation are used as ideal experimental data. Unfortunately, this technique is very time-consuming and is therefore limited to simple atomic systems. To overcome this problem, approximate models have been developed that can reproduce the simulated profiles, taking into account correlations and charge movements, and that can also be used to study complex emitters.

In line-shape models, the plasma environment of an emitter is taken into account by means of certain statistical properties such as the distribution or correlation functions of electric micro fields, positions, velocities, and so on.

<figure><img src="image_4.png" /><figcaption>Figure 3: Time history (in second) of micro field components in a Hydrogen plasma at $n_e = 10^{18}$ cm$^{-3}$ and $T_e = T_i = 1$ eV. $E_x$: full; $E_y$: dot; $E_z$: dash.</figcaption></figure>

With TCP simulations involving ions and electrons one gets static field distributions of the total field at ionic emitters, i.e., the electron field plus the ion field. Studies of the total field distribution functions have been carried out using HNC and APEX techniques [51]. However, these functions have no practical interest for line shape modeling as in the quasi-static approximation, only the slow varying component of the electric field contributes to the quasi-static Stark splitting of atomic states. The study of the electron dynamics shows that the low frequency does not only concern the ions but also a set of electrons which follow the movement of the ions. This is illustrated in Figure 4, which shows the correlations between directionalities of micro-field components,

$$C_{a,b}(\tau) = \int dt \vec{\Phi}_a(t) \cdot \vec{\Phi}_b(t + \tau), \quad (11)$$

where,

$$\vec{\Phi}_a(t) = \frac{\vec{E}_a(t)}{E_a(t)}, \quad (12)$$

and the indices $a$ and $b$ represent either electrons $e$ or ions $i$.

It can be seen that the correlation function of the electron electric field has two characteristic times, one very short and the other longer, of the same order as the correlation time of the ion electric field. Comparisons between correlation functions obtained accounting for all the interactions between charges, accounting for interactions between charges of same sign only and without interactions (ideal case) show that the low-frequency plateau in
<figure><img src="image_5.png" /><figcaption>Figure 4: Correlation functions of the micro field directionalities in the same conditions as Figure 3. Full line: All interactions between charges; Dash line: Interactions between charges of same sign only; Dot line: Ideal case.</figcaption></figure>

the electron field is due to interactions between ions and electrons (see Figure 4).

To extract its low-frequency component from the total electric field, in ref. [7], the following definition has been proposed:

$$ \begin{aligned} \vec{E}(t) &= \vec{E}_i(t) + \vec{E}_e(t) \\ &= \vec{E}_{Slow}(\Delta t; t) + \vec{E}_{Fast}(\Delta t; t), \end{aligned} \quad (13) $$

with

$$ \vec{E}_{Slow}(\Delta t; t) = \frac{1}{\Delta t} \int_{-\frac{\Delta t}{2}}^{\frac{\Delta t}{2}} \vec{E}(t - t') dt'. \quad (14) $$

The choice of $\Delta t$ must be a compromise between the plasma and atomic time scale constraints. In the case of an Hydrogen plasma, in [7], it has been shown that the choice of $\Delta t$ lower than the time corresponding to the end of the temporal plateau of the field correlation function did not make it possible to recover the well-known low-frequency distribution [54], in which it is assumed that ions interact with each other through an effective potential which includes electron-ion shielding. An increase in $\Delta t$ beyond this time led to an average of the low-frequency field component and this, due to the ion (proton) dynamics. The same simulations have then been performed with frozen ions for removing ion dynamics and as it is shown in Figure 5, an increase of $\Delta t$ leads correctly to the Hooper's low frequency distribution.

In [52], the question of how best to extract the quasistatic (low-frequency) microfield from a classical molecular dynamics simulation is explored in

<figure><img src="image_6.png" /><figcaption>Figure 5: Slow field distribution components at neutrals.</figcaption></figure>

some detail, in the case of a dense carbon plasma. It was shown in this article that, for a given time interval, the time-averaged microfield does not change significantly, suggesting that the separation of the slow and fast components of the microfield has been completed. After defining a "slow gradient of the distribution" to measure changes in the distributions as a function of $\Delta t$, it appears that the minimum gradient is obtained for a $\Delta t$ that corresponds to the end of the temporal plateau of the field correlation function. According to the fact that atomic processes in plasmas depend on the time average of the microfields, atomic time scale characteristics have to be taken into account too.

### 3.2.2. Ionization potential depression.

The radiative properties of an atom or an ion surrounded by a plasma are modified through various mechanisms. Depending on plasma conditions, the electrons supposedly occupying the upper quantum levels of radiators no longer exist. All the charges of the plasma contribute to the lowering of the energy required to free a bound electron. This mechanism is known as ionization potential depression (IPD). The evaluation of the IPD has important implications for dense plasma physics and in particular for the detailed prediction of the dense plasma equation of state and radiative opacity in stellar interiors, inertial confinement fusion research, or planetary interiors.

This evaluation deals with highly complex N-body coupled systems, involving particles with different dynamics and attractive ion-electron forces. There are several theoretical IPD models that apply to various conditions. Two of them, the Stewart-Pyatt (SP) [55] and the Ecker-Kröll (EK) [56] mod-
els which apply across a wide range of densities,
allow a general discussion of experiments.

In experiments performed at LCLS (Stanford) on
Aluminum [57] the observation of the K-α fluores-
cence and the measurement of the position of the
K-edge of ions show that the formula of Ecker and
Kröll is more adequate than the formula of Stewart
and Pyatt. Nevertheless, the agreement is not sat-
isfactory for the highest ion charges. On the other
hand, in an experiment performed at the Orion
laser facility (UK) [58], with a plasma at higher
temperatures, [500-700] eV, and densities in the
range [1-10] g.cm⁻³, the Aluminum K-shell spec-
trum shows a better agreement with calculations
performed with the Stewart–Pyatt IPD rather than
the Ecker–Kröll one. These two main experiments
have renewed interest for this issue and have stimu-
lated many theoretical investigations of IPD (see for
instance references [59–62] and references therein).

Here, the study of IPD is illustrated and dis-
cussed for aluminum and magnesium plasmas in the
conditions of both Ciricosta's and Hoarty's experi-
ments.

To get closer of the conditions of Ciricosta's ex-
periments, one can use CMD to simulate a two com-
ponent plasma of ions at room temperature and
solid density, and electrons in pseudo equilibrium
with the cold ion population. The argument here is
that the electron adjustment to the ions can occur
in a time that does not allow the ion population to
be heated by the electrons.

<figure><img src="image_7.png" /><figcaption>Figure 6: Ionization energy probability density function for various ion charges in an Aluminum plasma in the Ciricosta's experiment conditions.</figcaption></figure>

The BinGo-TCP code has then been used to sim-
ulate an Aluminum plasma with typically the ion

temperature, $T_i$ = 300 K, the electron temperature,
$T_e$ = 50 eV and a ionic density $\rho$ = 2.7 g.cm⁻³.

By taking advantage of the characteristics of the
ionization protocol, it is possible, when the ion is
in a pre-ionization state, to measure the energy re-
quired to ionize an electron from the ground state of
an ion, accounting for all interactions with the sur-
rounding charges. Due to the fluctuating local en-
vironment of the ions, this ionization energy is rep-
resented by a probability density function (PDF).
Figure 6 shows these PDFs for the various ion
charges present in the plasma.

By defining the ionisation potential depression
as the difference between the ionisation energy of
the isolated ion (see [63], for example) and the ion-
isation energy obtained by simulation, we obtain
an IPD distribution function that can be approxi-
mated by a Gaussian. In Figure 7, we have plot-
ted the average values calculated over the simulated
IPD distributions (red stars) and compared them
with the SP and EK models (dash-dotted and full
lines, respectively) and experimental values (blue
stars). The pink shaded area shows the dispersion
of the IPD values linked to the widths of the distri-
butions.

<figure><img src="image_8.png" /><figcaption>Figure 7: Ionization potential depression for an Aluminum plasma at solid density and $T_i$ = 300 K and $T_e$ = 50 eV. Comparisons of simulation results - red stars -, models - SP dash-dotted line, EK full line - and experimental data - blue stars -. The pink shaded region represents the dispersion due to the PDF.</figcaption></figure>

The simulation results fall in between the two
models and are in qualitative agreement with the
experimental data.

Still using the experiments of Ciricosta and col-
leagues as support, we have simulated a Magnesium
plasma at solid density (ρ = 1.74 g.cm⁻³). The re-
sults are shown in Figure 8.
<figure><img src="image_9.png" /><figcaption>Figure 8: Same as Figure 7, for a Magnesium plasma at solid density and $T_i = 300$ K and $T_e = 50$ eV.</figcaption></figure>

In Figure 9 (a) and (b), the Magnesium results are shown with those for Aluminium to highlight the effects of density. As expected, the decrease in density leads to a decrease in IPD. Here again, the simulation results are in good qualitative agreement with the experiments, although the density effects appear to be greater in the simulations than in the experiments. It should be noted that radiation is not taken into account in our simulations, and although the conditions chosen are close to those used in experiments, we are not trying to reproduce the experimental results exactly.

<figure><img src="image_10.png" /><figcaption>Figure 9: Ionization potential depression for Magnesium and Aluminum plasmas at solid density. Density effects (a)- Ciri-costa's experimental data, (b)- simulation results.</figcaption></figure>

The results of simulations carried out under conditions close to those of Hoarty's experiment are shown in Figure 10. Here ionic and electronic temperatures are equal. Aluminum Plasmas such that $\rho = 6 \text{ g.cm}^{-3}$ have been simulated for two temper-

<figure><img src="image_11.png" /><figcaption>Figure 10: Ionization potential depression in an Aluminum plasmas at $\rho = 6 \text{ g.cm}^{-3}$ and $T_i = T_e = 550 \text{ eV}$ (dark blue) and $T_i = T_e = 200 \text{ eV}$ (light blue). In red are plotted the ionization energies for isolated ions (stars) and surrounded ions (points).</figcaption></figure>

atures, $T_i = T_e = 550 \text{ eV}$ and $T_i = T_e = 200 \text{ eV}$. Figure 10 shows the simulation results compared to SP and EK models. The IPD are plotted in blue whereas, in red, the ionisation energies for isolated Aluminum ions and for Aluminum ions immersed in a plasma at 550 eV are compared, in order to visualise the effects of local fluctuations with as little confusion as possible in the graph.

According to the simulations, a decrease in temperature leads to an increase in the IPD. In our model, changes in temperature lead to a change in the average ionic charge and therefore a change in the electron density. Here, we have measured, $\bar{Z} = 10.82$ and $\bar{Z}^2 = 117.4$ for the case $T_i = T_e = 550$ eV and $\bar{Z} = 10.2$ and $\bar{Z}^2 = 104.18$ for the case $T_i = T_e = 200$ eV, giving $n_e = 1.45 \times 10^{24} \text{ cm}^{-3}$ and $n_e = 1.37 \times 10^{24} \text{ cm}^{-3}$ and $\Gamma \simeq 5.3$ and $\Gamma \simeq 6$, respectively (Γ being the plasma coupling parameters). The correlation effects between particles are therefore stronger when the temperature is lower, which may partly explain the trend observed in the simulated results. The calculations performed with the EK and SP models do not show the same trend. The results obtained with the SP model show no significant difference with the change in temperature. Only the results at $T_i = T_e = 550$ eV are shown so as not to overload the graph. The EK model shows a slight increase in the IPD with temperature.

The strength of our simulations lies in their abil-
ity to take into account all the interactions between particles as well as the dynamics of charged particles, ions and electrons. These effects are neglected in most models in the literature. The role and importance of charge fluctuations in IPD have already been the subject of research. In [64], it is suggested that fluctuations, which are neglected in models but are essential to describe the absorption of energy by a system, are important enough to have an impact on the interpretation of experimental results. More recently, the effects of ionisation potential depression and fluctuations on electron impact ionisation have been studied [65].

These investigations will be continued in order to gain a better understanding of the influence of particle correlation and fluctuation effects on atomic properties in plasmas.

## 4. Conclusion

We have presented the BinGo series of codes based on classical molecular dynamics techniques and designed primarily for the study of dense plasmas. CMD simulations appear to be a very versatile and powerful tool for studying these plasmas. A few years ago we developed the BinGo-TCP code, which has the ability to simulate neutral plasmas containing electrons and ions of different charges that interact with each other. Our aim in developing the BinGo-TCP code was to avoid excess of complexity of models in order to preserve a straightforward interpretation of physics from simulation results. The results obtained with our TCP-MD simulation code are very encouraging. They compare well with results obtained by other models or simulations and experimental results when they exist. Our simulations provide data for further discussion on IPD models and more generally on the charge correlation effects on the plasma properties.

## References

[1] B. J. Alder and T. E. Wainwright, J. Chem. Phys. 27 (1957) 1208; J. Chem. Phys. 31 (1959) 459.

[2] A. Rahman, Phys. Rev. 136 (1964) A405.

[3] A. Calisti, S. Ferri, C. Mossé, B. Talin, A. Klisnick, L. Meng, D. Benredjem and O. Guilbaud, High Energy Density Physics 9 (2013) 516-522.

[4] B. Talin, A. Calisti, J. W. Dufty, I. V. Pogorelov, Phys. Rev. E 77 (2008) 036410.

[5] A. Calisti, E. Dufour, B. Talin and J. W. Dufty, AIP Conf. Proc. 645 (2002) 252.

[6] A. Calisti, A. Demura, M. A. Gigosos, D. González-Herrero, C. A. Iglesias, V. S. Lisitsa, E. Stambulchik, Atoms 2 (2014) 259-276.

[7] A. Calisti, S. Ferri, C. Mossé, B. Talin, M. A. Gigosos, M. A. González, High Energy Density Physics 7 (2011) 197-202.

[8] A. Calisti, S. Ferri, M. Marciante, B. Talin, High Energy Density Physics 13 (2014) 1-8.

[9] A. Calisti, S. Ferri and B. Talin, Contrib. Plasma Phys. 55, Issue 5 (2015) 360-365.

[10] D. C. Rapaport, The Art of Molecular Dynamics Simulation, Cambridge University Press, 2004.

[11] N. Metropolis et al., J. Chem. Phys. 21 (1953) 1087.

[12] L. Verlet, Phys. Rev. 159 (1967) 98-103.

[13] D. V. Fisher and Y. Maron, Eur. Phys. J. D 14 (2001) 349-359.

[14] M. Marciante, C. Champenois, A. Calisti, M. Knoop, Appl. Phys. B 107 (2012) 1117-1123.

[15] T. Raitza, H. Reinholz, G. Röpke and I. Morozov, J. Phys. A: Math. Theor. 42 (2009) 214048.

[16] Lin-Lin Wang, Ajith Perera and Hai-Ping Cheng, Phys. Rev. B 68 (2003)115409.

[17] C. Ticknor, S. D. Herring, F. Lambert, L. A. Collins, and J. D. Kress, Phys. Rev. E 89 (2014) 013108.

[18] P. H. Hünenberger Adv. Polym. Sci. 173 (2005) 105-149.

[19] J. P. Hansen, I. R. McDonald, E. L. Pollock, Phys. Rev. A 11, 3 (1975)1225.

[20] S. Galam and J. P. Hansen, Phys. Rev. A 14, 2 (1976) 816.

[21] R. W. Hockney, J. W. Eastwood, Computer Simulation using Particles. McGraw Hill Company, New York, 1981.

[22] H. Minoo, M. M. Gombert, and C. Deutsch, Phys. Rev. A 23 (1981) 924.

[23] M. Belkacem, F. Megi, P.-G. Reinhard, E. Suraud, and G. Zwicknagel, Eur. Phys. J. D 40, (2006) 247-255.

[24] J. P. Hansen and I. H. McDonald, Phys. Rev. Lett. 41 (1978) 1379; Phys. Rev. A 23 (1981) 2041.

[25] T. Dunn and A. A. Broyles, Phys. Rev. 157 (1967) 156.

[26] C. Deutsch, Phys. Lett. 60A (1977) 317.

[27] G. Kelbg, Ann. Physik 13 (1964) 354; 14 (1964) 394.

[28] W. Ebeling, G. E. Norman, A. A. Valuev and I. Valuev, Contrib. Plasma Phys. 39 (1999) 61.

[29] A. Filinov, M. Bonitz, W. Ebeling, J. Phys. A (Math. Gen.) 36 (2003) 5957-5962.

[30] C. Deutsch, M. M. Gombert and H. Minoo, Phys. Lett. A 66 (1978) 381; Phys. Lett. 72 (1979) 481.

[31] C. S. Jones and M. S. Murillo, High Energy Density Phys. 3 (2007) 379-394.

[32] C. Dorso, S. Duarte, J. Randrup, Phys. Lett. B 188 (1987) 287.

[33] H. Wagenknecht, W. Ebeling, A. Förster, Contrib. Plasma Phys. 41, 1 (2001) 15-25.

[34] K. LaGattuta, Opt. Express 8 (2001) 401.

[35] R. Car, M. Parrinello, Phys. Rev. Lett. 55 (1985) 2471.

[36] G. Zérah, J. Clérouin, and E. L. Pollock, Phys. Rev. Lett. 69 (1992) 446.

[37] F. Lambert, J. Clérouin, G. Zérah, Phys. Rev. E 73 (2006) 016403.

[38] F. R. Graziani et al., High Energy Density Phys. 8 (2012) 105-131.

[39] S. P. Hau-Riege, J. Weisheit, J. I. Castor, R. A. London, H. Scott and D. F. Richards, New Journal of Physics 15 (2013) 015011.

[40] H.-K. Chung, M.H. Chen, W.L. Morgan, Yu. Ralchenko
and R.W. Lee, High Energy Density Phys. 1 (2005) 3-12.

[41] H. R. Griem, Spectral line broadening by plasmas, Academic, New York, 1974.

[42] H. R. Griem, Phys. rev. A 33 (1986) 3580.

[43] S. G. Rautian and I. I. Sobel’man, Usp. Fiz. Nauk 90 (1966) 209.

[44] G. J. Pert, JOSA B 11 (1994) 1425.

[45] H. R. Griem, Principles of Plasma Spectroscopy, Cambridge University Press, 1997.

[46] R. H. Dicke, Phys. Rev. 89 (1953) 472.

[47] A. Calisti et al., Phys. Rev A 42 (1990) 5433.

[48] J. A. Koch et al., Phys. Rev. A 50 (1994) 1877.

[49] https://plasma-gate.weizmann.ac.il/slsp/

[50] E. Stambulchik, D. V. Fisher, Y. Maron, H. R. Griem, and S. Alexiou, High Energy Density Phys. 3 (2007) 272.

[51] H. B. Nersisyan, C. Toepffer, and G. Zwicknagel, Contrib. Plasma Phys. 50 (2010) 193.

[52] S. Hau-Riege and J. Weisheit, Phys. Rev. E 91 (2015) 033106.

[53] M. A. Gigosos, D. González-Herrero, N. Lara, R. Florido, A. Calisti, S. Ferri, and B. Talin, Phys. Rev. E 98, 3, (2018) 033307.

[54] C. F., Jr. Hooper, Phys. Rev. 165 (1968) 215.

[55] J. Stewart and K., Jr. Pyatt, Astrophys. J. 144 (1966) 1203.

[56] G. Ecker and W. Kröll, Phys. Fluids 6 (1963) 62.

[57] O. Ciricosta et al., Phys. Rev. Lett. 109 (2012) 065002; O. Ciricosta et al., Nat. Commun. (2016) 7:11713.

[58] D. J. Hoarty et al., Phys. Rev. Lett. 110 (2013) 265003.

[59] S. K. Son, R. Thiele, Z. Jurek, B. Ziaja and R. Santra, Phys. Rev. X 4 (2014) 031004.

[60] S. Vinko, O. Ciricosta and J. Wark, Nat. Commun. 5 (2014) 3533.

[61] C. Lin, G. Röpke, W. D. Kraeft and H. Reinholz, Phys. Rev. E 96 (2017) 013202; C. Lin, Phys. Plasmas 26 (2019) 122707.

[62] J. L. Zeng, Y. Li, C. Gao and J. Yuan, Astron. Astrophys. 634 (2020) A117; J. L. Zeng et al., Results in Physics 40 (2022) 105836.

[63] R. D Cowan, The theory of atomic structure and spectra, University of California Press, Berkeley. Los Angeles. London, 1981.

[64] C. A. Iglesias and P. A. Sterne, High Energy Density Physics 9 (2013) 103 - 107.

[65] D. Benredjem, J. C. Pain, A. Calisti, S. Ferri, Phys. Rev. E 108 (2023) 035207.