On certain aspects of the THERMOS toolkit for modeling experiments

TL;DR

This paper reviews the THERMOS toolkit's capabilities in modeling plasma radiative properties, highlights recent upgrades for dense plasma effects, and demonstrates its application through experimental spectrum analysis.

Contribution

It introduces recent enhancements to the THERMOS toolkit, including ionization potential lowering effects, and showcases its use in analyzing experimental plasma spectra.

Findings

Successful modeling of chlorine spectra

Effective analysis of silicon plasma transmission

Validation against experimental data

Abstract

The THERMOS toolkit has been developed to calculate radiative properties of plasmas. This article contains a brief survey of some of its key features used by calculation of opacities and emissivities and by analysis of specific experiments. The code has recently been upgraded to account for the effect of ionization potential lowering in dense plasmas. The functionality of the code is illustrated for several cases from the 10th NLTE Code Comparison Workshop, in particular, for the experimental spectra of chlorine [1] and for the measured transmission of a silicon plasma [2].

Figures

Equations

$$\Delta I^{SP}(Z)=3(Z+1)/2r_0,$$

$$\Delta I^{EK}(Z)=(Z+1)/r_{EK},$$

$$T_{\omega }=\exp (-{\varkappa }_{\omega }\rho L),$$

$$F_{\omega }=\pi \frac{j_{\omega }}{\rho {\varkappa }_{\omega }}\left[1-e^{-\tau }+\tau e^{-\tau }-{\tau }^2{E}_1(\tau )\right],$$

$$I_{\omega }=B_{\omega }=K\frac{15\sigma }{{\pi }^5}\frac{{\omega }^3}{\exp (\omega /T_{rad})-1},$$

$$\mu \frac{d{I}_{\omega }}{dz}=j_{\omega }-\rho {\varkappa }_{\omega }{I}_{\omega },$$

$$\sum _{j^{\prime }{s}^{\prime }}\left(n_{j^{\prime }{s}^{\prime }}{w}_{j^{\prime }{s}^{\prime }\to js}-n_{js}{w}_{js\to j^{\prime }{s}^{\prime }}\right)=0,$$

$$I_{\omega }(z=0,\mu >0)=B_{\omega },$$

$$I_{\omega }(z=L,\mu <0)=0,$$

$$F_{\omega }^{+}=2\pi \underset{0}{\overset{\frac{\pi }{2}}{\int }}{I}_{\omega }(z=L)\cos \theta \sin \theta d\theta .$$

Full text

††journal: Journal of LaTeX Templates

1 Introduction

Originally, the THERMOS code was developed at KIAM RAS (Keldysh Institute of Applied Mathematics of the Russian Academy of Sciences) by the research team of A.F. Nikiforov, V.B. Uvarov and V.G. Novikov to calculate the properties of equilibrium plasmas in a wide range of temperatures and densities. It is based on the non-relativistic modified Hartree-Fock-Slater model in the average atom approximation and is described in the book “Quantum Statistical Models of High-Temperature Plasma” [3]. Spectral opacity tables, calculated with the THERMOS code, have been extensively used to solve various problems of radiation dominated plasmas; see for example Refs. [4, 5, 6, 7].

The development of new powerful facilities, such as NIF and Z, gave rise to a large amount of experimental data for non-equilibrium plasmas that resulted in revision of many existing models. Additional strong motivation came also from the development of extreme ultraviolet light sources, where simulations require high accuracy of line positioning and new methods for computing the spectral tables. These problems gave a major impetus for further development of the THERMOS toolkit.

A number of codes, developed earlier by our research team for calculating the radiative properties of plasmas, have recently been consolidated into the THERMOS toolkit software package. In result, it has become a modern instrument that allows calculation of the emission and absorption spectra in a wide range of temperatures and densities, both in the approximation of local thermodynamic equilibrium (LTE), and for the non-equilibrium (NLTE) plasmas with an arbitrary radiation field. The newly implemented models take into account the density effects by allowing for the ionization potential depression (IPD) according to the Stewart-Pyatt [8] and/or the Ecker-Kröll [9] formulae, and for the ion broadening of spectral lines in the two-level approximation [10]. Also, a new code module has been added where the equations of radiation transfer and level kinetics are solved self-consistently in a steady-state approximation across a multi-layer slab sample.

The THERMOS toolkit is a regular participant at the NLTE Code Comparison Workshops [11], and its calculation results demonstrate a fair agreement with the other well-known codes from all over the world. At the 10th NLTE Workshop (2017, San Diego, CA, USA), the experimentally obtained emission and transmission spectra for steady-state cases of the aluminum, chlorine, and silicon plasmas were proposed as benchmark problems for simulation. Some of the simulation results, obtained with the THERMOS toolkit for this workshop, were recently presented at the RPHDM18 [12] conference and are described in some detail below.

2 THERMOS toolkit overview

The THERMOS toolkit includes codes for numerical simulation of the LTE and non-LTE plasmas, see Fig. 1. Calculations for an LTE plasma can be carried out “ab initio” in the average atom approximation. Currently, three variants of this model are available: the Thomas-Fermi model (THERMOS_TF), the modified Hartree-Fock-Slater (THERMOS_HFS) model, and the quasizone model (THERMOS_QZM) [13]). In addition, for a low density plasma, the Saha-Boltzmann (THERMOS_SB) statistics can be used. For the non-LTE cases, the system of level-kinetics equations is solved in the quasi-stationary approach with a fixed ambient radiation field by using the collisional-radiative equilibrium (THERMOS_CRE) model. The THERMOS_SB and THERMOS_CRE models require an external atomic database.

Calculation of the atomic database for a particular substance is based on the isolated-ion model with automatic selection of non-relativistic configurations of ions according to given criteria. The toolkit includes a module that allows to improve the positions of spectral lines and the oscillator strengths by using external data from detailed atomic codes, such as RCG [14] and FAC [15], as well as the known experimental data. In addition to that, a special technique has been developed for averaging the atomic data on a given photon energy grid — the Radiative Unresolved Spectra Atomic Model or RUSAM [16], which is aimed at reduction of the calculation time with little to none detriment to the accuracy.

3 Solid-density aluminum plasma

The aluminum cases, discussed at the 10th NLTE Workshop, were aimed to explore the collisional-radiative kinetics and the K-shell emission in a solid-density plasma relevant to the XFEL experiments [17]. To compare how the density effects are accounted for in different codes, it was proposed to calculate the emissivity of the aluminum plasma at temperatures $T_{e}=10$, 30, 100, and 300 eV, and the electron number densities $N_{e}=2\times 10^{23}$ and $5\times 10^{23}$ cm*-3*. To compare THERMOS with other codes, we have selected points with the electron density $N_{e}=5\times 10^{23}$ cm*-3* and the temperatures $T_{e}=30$ and 100 eV.

There is some uncertainty with regard to which of the existing IPD models would be the most appropriate for modeling a dense aluminum plasma [17, 18]. In the THERMOS toolkit, we have a choice between two IPD models. In the Stewart-Pyatt (SP) approximation [8], the ionization potential depression $\Delta I$ (in atomic units) for an ion with a net charge $Z$ is given by

[TABLE]

where the dimensionless radius $r_{0}$ of the electrically neutral atomic cell is defined by the relation $(4\pi/3)(r_{0}a_{0})^{3}N_{i}=1$; here $N_{i}=N_{e}/Z$ is the ion number density, $a_{0}=0.529\times 10^{-8}$ cm is the Bohr radius. The Ecker-Kröll (EK) approximation [9] uses the expression

[TABLE]

where $r_{EK}^{3}=r_{0}^{3}/(1+Z_{0})$ is the modified atomic сell radius depending on the mean ion charge $Z_{0}$.

Initially, the atomic database was calculated for isolated aluminum ions, where only the excited states with the principal quantum numbers $n\leq n_{max}=6$ were taken into account. Then, as the IPD correction was applied, all the excited states with the binding energies below $\Delta I$ were discarded, which caused $n_{max}$ to decrease. In the present approximation, the radius $r_{EK}$ defines by the ion charge $Z$.

The calculated dependence of $n_{max}$ on the ion charge $Z$ for the electron density $N_{e}=5\times 10^{23}$ cm*-3* and temperature $T_{e}=100$ eV is shown in Fig. 2. This figure shows also the results by other participants of the 10th NLTE Workshop (gray color). One observes a broad scatter of the results. The THERMOS toolkit yields $n_{max}=2$ for the ions with charges $3\leq Z<9$ within the Stewart-Pyatt approximation, and $n_{max}=2$ for all the ions $Z\geq 3$ when the Ecker-Kröll formula is used. Ions with $Z<3$ are not present under the considered conditions. Both approximations produce very close results, and we cannot give preference to any of them.

Figure  3 shows the emissivity of the aluminum plasma, calculated with the THERMOS toolkit, in comparison with the other results from the 10th NLTE Workshop. Generally, a good agreement with the results by other scientific groups is observed, although the THERMOS spectrum is shifted towards lower photon energies due to the lower mean ion charge (see Fig. 4).

To clarify the cause of the above discrepancy, we examined the case with the electron density $N_{e}=5\times 10^{23}$ cm*-3* and the temperature $T_{e}=30$ eV. Under these conditions the plasma can be treated as LTE, and the THERMOS_HFS code is applicable for performing the reference calculation. As is seen in Fig. 5, whereas the THERMOS_CRE produced a mean charge of $Z_{0}=3.2$, the THERMOS_HFS yielded $Z_{0}=4.3$ — which is closer to the 10th NLTE results. The emissivity of the aluminum plasma, calculated with the THERMOS_CRE and THERMOS_HFS and compared with other results from the 10th NLTE Workshop, is displayed in Fig. 6. This figure indicates that the THERMOS_HFS code does not describe the line positions correctly, although the bound-free cross-sections seem to be in order; also, the broadening of the photoabsorption edges due to the density effects is clearly seen in contrast to the corresponding THERMOS_CRE results. Here we conclude that in solid-density plasmas the simplest approximations [19] for the ionization and recombination rates, used in the present THERMOS_CRE version, are not applicable and need to be reconsidered.

4 ORION experiment with the chlorine plasma

The simulation of the emissivity of the chlorine plasma was related to the recent high-resolution measurements on the OHREX spectrometer fielded at the Orion facility, as is shown in Fig. 5(b) of Ref. [1]. According to the experimental data, it was not possible to obtain accurate values for the temperature and density of the plasma, though the ranges of their variation were evaluated as $T_{e}\approx 400$–600 eV and $N_{e}\approx 10^{21}$–$10^{23}$ cm*-3*. At the Workshop it was proposed to determine the best fit to the experiment. To reach this goal, it was necessary to calculate the emissivity and to investigate the widths and intensities of the K-lines in a dense chlorine plasma within the given ranges of temperature and density.

In the experimental emission spectrum of the chlorine plasma, shown in Fig. 7, two groups of lines can be distinguished. The first one corresponds to the $1s3p\rightarrow 1s^{2}$ transitions in the helium-like Cl*+15* ion. The second group corresponds to the lithium-like Cl*+14* ion transitions $1s2s3p\rightarrow 1s^{2}2s$ and $1s2p3p\rightarrow 1s^{2}2p$. The best fit solution was found by adjusting the ratio of the intensities of the Cl*+14* satellite lines $I(1s2s3p\rightarrow 1s^{2}2s)/I(1s2p3p\rightarrow 1s^{2}2p)$, calculated by solving the system of the level kinetics equations, to its measured value. In result, the best agreement was achieved for the electron density $N_{e}=10^{23}$ cm*-3* and the temperature $T_{e}=600$ eV.

Under these conditions the chlorine plasma is optically thick, so one can assume that the measured radiation originates from a thin surface layer. Figure 8 shows two transmission spectra of the chlorine plasma at $N_{e}=10^{23}$ cm*-3* and $T_{e}=600$ eV, calculated for two values of the layer thickness $L$, namely, for $L=4\times 10^{-5}$ and 10*-4* cm. The transmission of a layer of thickness $L$ at a photon energy $\omega$ was calculated as

[TABLE]

where $\varkappa_{\omega}$ (cm2/g) is the opacity, and $\rho$ is the mass density. Figure 8 clearly indicates that for $L=4\times 10^{-5}$ cm there already is a significant amount of reabsorption in the $1s3p\rightarrow 1s^{2}$ line, while for greater thicknesses $L\geq 10^{-4}$ cm the reabsorption must be taken into account for practically all the observed spectral lines.

The radiation flux, emitted by a uniform slab with a thickness $L$ at a photon energy $\omega$, can be calculated as

[TABLE]

where $E_{k}(\tau)=\int\limits_{\tau}^{\infty}t^{-k}e^{-t}dt$, $\tau=L\rho\varkappa_{\omega}$ is the optical thickness, and $j_{\omega}$ and $\varkappa_{\omega}$ are, respectively, the spectral emissivity and the opacity. The radiation flux, calculated in this way for $L=4\times 10^{-5}$ cm with the $j_{\omega}$ and $\varkappa_{\omega}$ provided by the THERMOS_CRE code for $N_{e}=10^{23}$ cm*-3* and $T_{e}=$600 eV, is shown in Fig. 7 together with the experimental data. It is seen that a fair agreement with the measured emission spectrum can be achieved by using the approximation of a thin uniform layer.

5 Transmission of a silicon plasma slab

At the 10th NLTE Workshop, a special interest was paid to the results of Ref. [2], where the experimental transmission spectrum for soft X-rays through a slab of photoionized silicon plasma was presented — under the conditions that are close to astrophysical environment. In such a plasma, radiative processes dominate over the collisional ones due to the presence of an intense external radiation field, which pumps the autoionization states and strongly enhances the role of the Auger transitions.

For the numerical modeling of this experiment, it was proposed to consider a homogeneous slab of a silicon plasma, whose thickness $L$, temperature $T_{e}$, and the electron density $N_{e}$ are treated as free parameters; for the intensity $I_{\omega}$ of the ambient radiation field, it was proposed to use the diluted black-body spectrum

[TABLE]

where $K$ is the dilution factor, $\omega$ is the photon energy, $T_{rad}$ is the radiation temperature, and $\sigma$ is the Stefan-Boltzmann constant; $I_{\omega}$ and $B_{\omega}$ are defined per unit photon energy interval. The plasma parameters were supposed to be determined by the best fit to the measured transmission spectrum. The analysis, carried out in Ref. [2], allowed to outline the bounds for the search area: the plasma size should be in the range $L\simeq 0.1$–1 cm, the electron density $N_{e}\simeq 10^{19}$ cm*-3*, and the temperature $T_{e}\simeq 30$ eV, with the mean ion charge being around $Z_{0}\approx 10$.

First we consider one of the cases proposed for the preliminary analysis, namely, $L=0.1$ cm, $N_{e}=10^{19}$ cm*-3*, $T_{e}=30$ eV, $K=1$, and $T_{rad}=63$ eV. For further discussion, this case is denoted the “uniform slab”, which means that the electron temperature and density, as well as the radiation energy density $\displaystyle U_{\omega}=\frac{4\pi}{c}B_{\omega}$ are constant over the plasma volume. Calculation using the THERMOS_CRE code for a uniform slab with the radiation field (3) resulted in the mean charge $Z_{0}=11.8$. This value does not match the estimated optimum of $Z_{0}~{}\approx$ 10, and, accordingly, the calculated transmission spectrum (the blue curve in Fig. 11) is rather far from the experimental one. In the course of discussions at the Workshop, it was decided to diminish the role of the ambient radiation field by setting $K=0.1$, which led to the mean ion charge lying closer to 10 and to a better agreement with the experimental transmission spectrum. Here, the dilution factor of 0.1 made it possible to artificially take into account the effects of the reabsorption of radiation in a dense plasma.

At the same time, a more thorough analysis of the opacity of the silicon plasma revealed that the major part of the external radiation flux was absorbed in a thin edge layer, and was unable to penetrate deep into the slab. For the relevant values of the temperature, electron density and radiation field, a layer with a thickness in excess of $10^{-3}$ cm is opaque for the photons in the investigated wavelength range (6.64–6.96 Å). Hence, the spatial distribution of the radiation field must be treated as non-uniform. In order to take into account the effect of reabsorption in the plasma slab, we have modified the problem statement and the numerical model.

In the new formulation, the initial slab with $L=0.1$ cm is divided into $N$ layers ($z_{i-1}<z<z_{i},\,i=1,N$). In each layer, the electron density and temperature are kept the same as in the initial case of a uniform slab, i.e. $N_{e}=10^{19}$ cm*-3* and $T_{e}=30$ eV. Also, the radiation energy density $U_{\omega}=\displaystyle\frac{1}{c}\int_{4\pi}I_{\omega}d\Omega$ is assumed to be constant within each layer, but, in contrast to $N_{e}$ and $T_{e}$, it is allowed to vary from layer to layer. This case is called the “multi-layer slab” (see Fig. 9).

The radiation intensity $I_{\omega}(z,\theta)$ is found as a solution of the transfer equation

[TABLE]

where $\mu=\cos\theta$. The emissivity $j_{\omega}(z,n_{js})$ and the opacity $\varkappa_{\omega}(z,n_{js})$ at given $T_{e}$, $N_{e}$, and $U_{\omega}$ are determined by the populations $n_{js}(z)$ of the ion species with the ionization degree $j$ in a quantum state $s$. The values of $n_{js}(z)$ are found as a solution of the system of rate equations

[TABLE]

where $w_{js\to j^{\prime}s^{\prime}}(U_{\omega})$ is the total probability of transition from an ion state $js$ to a state $j^{\prime}s^{\prime}$ due to the relevant elementary processes [3].

To calculate the plasma radiative properties, one has to solve self-consistently the transfer equation (4) together with the rate equations (5). The boundary condition for Eq. (4) takes into account the incident flux of the external radiation,

[TABLE]

[TABLE]

where $B_{\omega}$ is given by Eq. (3) with $K=1$ and $T_{rad}=$63 eV. The self-consistent solution of equations (4) and (5) for $I_{\omega}(z)$ and $n_{js}(z)$ is obtained by successive iterations. At each iteration, the transfer equation (4) is solved exactly for given piece-wise constant coefficients $j_{\omega}(z)$ and $\rho\varkappa_{\omega}(z)$ — in accordance with the discretization of $z$. At the initial iteration $\nu=0$, the level populations in the first layer $z_{0}<z<z_{1}$ are calculated with the $\displaystyle U_{\omega}=\frac{4\pi}{c}B_{\omega}$, and in the remaining $N-1$ layers at $z_{1}<z<z_{N}$ with the $U_{\omega}=0$. All the subsequent iterations $\nu\geq 1$ use the values $U_{\omega}^{(\nu-1)}$, calculated at the previous step $\nu-1$. The iterations are assumed to have converged when the relative accuracy of $10^{-3}$ is reached in the $L_{1}$ norm. The radiation flux at the right boundary is calculated as

[TABLE]

The reference calculation, compared with the experiment, was performed for a slab divided into $N=5$ layers. All the layer thicknesses were equal to $L/N$. The obtained solution has a non-uniform distribution of the radiation intensity in the slab. Figure 10 compares the two spectral fluxes, emitted from the right boundary $z=L$, as calculated from Eq. (8) in the multi-layer approximation (the red curve), and for the uniform slab from Eq. (2) (the blue curve). The two spectra significantly differ from each other because in the “multi-layer slab” case the mean ion charge decreases from $Z_{0}=11.6$ in the leftmost layer to $Z_{0}=9.7$ in the rightmost layer, whereas for the “uniform slab” it does not depend on $z$ and is everywhere equal to $Z_{0}=11.8$.

Figure 11 shows the transmission spectra for the considered silicon plasma slab: the two theoretical spectra, calculated respectively in the multi-layer approximation (the red curve) and for the uniform slab (the blue curve), are compared with the measured spectrum (the black curve) from Ref. [2]. One clearly sees that the “multi-layer slab” spectrum gives a significantly better agreement with the experiment than the transmission of a uniform slab because the non-uniformity of the radiation field along the $z$-axis is properly accounted for. The deviations in the positions and strengths of individual lines are explained by the use of non-relativistic configurations when solving the rate equations (5).

6 Summary

Some of the key features of the THERMOS toolkit for modeling radiative plasma properties are described. The impact of the plasma density effects is analyzed for the aluminum and chlorine plasma cases. The THERMOS simulation results demonstrate a reasonable agreement with the experimental data, as well as with the other codes that participated in the 10th NLTE Workshop. The model of self-consistent calculation of radiation transfer and level kinetics in a multi-layer plasma slab was used to account for the radiation reabsorption, important for the case of non-equilibrium silicon plasma. In this way, a good agreement with the experimental transmission spectrum was achieved, except for the line strengths of certain ion species.

Acknowledgment

Authors thank the participants and organizers of the NLTE and RPHDM conferences for productive discussions of the problems arising from non-equilibrium plasma simulations. This work was funded by the Russian Science Foundation through grant No. 14-11-00699. Calculations have been performed at HPC MVS-10P (JSCC RAS) and HPC K100 (KIAM RAS).