Corrigendum: Models for the computation of opacity of mixtures

. Equation (4.1) in our paper (2013 New J. Phys. 15 015012) is incorrect. Its corrected version is provided here and justiﬁed. It applies to all the models. Abstract. We compare four models for the partial densities of the components of mixtures. These models yield different opacities as shown on polystyrene, acrylic and polyimide in local thermodynamical equilibrium (LTE). Two of these models, the ‘whole volume partial pressure’ model (M1) and its modiﬁcation (M2) are not thermodynamically consistent (TC). The other two models are TC and minimize free energy. M3, the ‘partial volume equal pressure’ model, uses equality of chemical potential. M4 uses commonality of free electron density. The latter two give essentially identical results in LTE, but M4’s convergence is slower. M4 is easily generalized to non-LTE conditions. Non-LTE effects are shown by the variation of the Planck mean opacity of the mixtures with temperature and density.


Introduction
Equation (4.1) in [1] gives the opacity (in cm 2 g −1 ) of a mixture as a linear combination of the opacities of the constituents. It is written there as κ(ρ, T ) = 1 ρ l q l ρ l κ l (ρ l , T ), (4.1) where q l is the number fraction of the element l in the mixture and ρ l is the model-dependent effective density of the latter.

Correction
The mistake in equation (4.1) stems from the confusion between the real density (hereafter ρ l ) and the effective density (hereafterρ l ). The 'real density' is an expression of the number of absorbing atoms, while the 'effective density' describes the physical state (e.g. the 'size') of the atoms caused by their environment. The latter is model dependent, as shown in [1]. This 1 Author to whom any correspondence should be addressed.
distinction is clear when one writes the linear attenuation coefficient as where N l is the number of atoms of element l in the unit volume V and σ is the cross-section (in cm 2 per atom) for absorption at the effective densityρ l and temperature T. Introducing the atomic weight A l (g per atom), we write where M l is the mass of the l atoms in the volume V and κ is the opacity in cm 2 g −1 . Writing now the real (also known as 'partial') density for element l like in paragraph 3.2 of [1]: where M is the mass of the mixture, m l is the mass fraction of element l and ρ mix is the total density of the mixture, we finally obtain a correction to equation (4.1) valid for all the models described in [1]: κ(ρ mix , T ) = 1 ρ mix l ρ l κ l (ρ l , T ) = l m l κ l (ρ l , T ).

Remark
With this correction, we always have l ρ l = ρ mix . But lρ l can be greater than ρ mix as shown in [1]. Thus, the error in (4.1) could lead to unphysical situations.

Introduction
In order to simulate the experiments of the Centre for Radiative Shocks (CRASH) project [1], opacities of various materials are needed. Among them are several kinds of 'plastics': polystyrene (CH), polyethylene (CH 2 ), acrylic (C 5 H 8 O 2 ) and polyimide (C 22 H 10 O 5 N 2 ). The aim of the project is not only to simulate the experiment, but also to quantify the uncertainty of the simulation by comparing different models and checking what influence these have on the final result. It is with this in mind that we carried out this work. When several elements 2 are present in a plasma at finite temperature and density, they interact and arrive at an equilibrium that may modify their properties as compared with their state when isolated. Various models have been proposed in order to find the modification of these properties. For instance, one can assume that the elements experience equal pressure or that all elements fill up the entire volume or that equilibrium is achieved when they all possess the same chemical potential (for a definition of the latter see section 2.1). All the models are built based on the assumption that all the elements are at the same temperature. The way the computations are performed in practice is by computing opacity databases separately for the various elements on a grid of temperatures and densities. In this work, we do not consider molecular and/or solid-state effects. Consequently, the way the mixture effects are taken into account is by finding partial densities for each element, for each temperature and density of the mixture. Therefore, each model consists of different ways of computing these partial densities. In this work, we consider four models: 1. M1, 'the partial pressure model', in which each element fills the entire volume and has a partial pressure according to Dalton's law of mass action. This simple model has been widely used for low-density astrophysical plasmas. 2. M2, an 'ad hoc' model which is shown as an example of the sensitivity of the assumptions.
It is the same as M1, but arbitrarily modified by requesting equality of chemical potentials instead of the strict application of Dalton's law (see section 2.1). 3. M3, for which elements are assumed to occupy partial volumes at equal pressure, with the requirement of equality of chemical potential (see section 2.1). 4. M4, which is the same as M3, but for which equality of free electron density replaces the requirement of equality of chemical potential. The advantage of the latter is that it can easily be generalized to non-local thermodynamical equilibrium (non-LTE) plasma conditions in the framework of the RADIOM model (see section 5.1 and [13]).
We will show by examples that M3 and M4 are equivalent and give essentially identical results. M1 and M2 give results differing from the latter. They are not thermodynamically consistent (TC) at finite densities.
In section 2, we briefly describe the computation procedure with the Super Transition Array (STA) and MIX (see section 2.2) codes. In section 3, the different models are described in some detail. Examples of the results are given in section 4. Approximate non-LTE effects are shown in section 5. We conclude, in section 6, by stating that M3 for computing opacities for the CRASH project is the model of choice for LTE plasmas, and M4 is equivalent and can give an estimate of non-LTE effects when needed.

The Super Transition Array code
For the opacity of each element, the STA model and code are used. Originally proposed and programmed by Bar-Shalom and Oreg [2,3], the STA code has recently been updated and made more stable [4].
This opacity code is based on the ion sphere model, with a Wigner-Seitz radius, R WS = N −3 i (N i is the ion volume density). The space between the touching, but not interacting, ion spheres is filled with a constant potential (the so-called 'jellium'), corresponding to a bath of free electrons. STA is working in the chemical picture: we consider the plasma as consisting of ions (nucleus + bound electrons) of each element and free electrons shared among all ions. This approach is in line with INFERNO and subsequent models [5,6]. Like in [5,6], the STA code starts by modelling the ions in the average atom framework [7] with a self-consistent distribution of bound and free electrons with all wave functions (bound and free) solutions of the Dirac equation. The bound and continuum electron densities are determined by populating the wave functions according to their statistical weight modulated by the Fermi distribution function f (ε, µ) = [1 + exp(ε − µ)/kT ] −1 with the chemical potential µ varied to ensure charge neutrality within R WS : where P and Q are the large and small components of the Dirac wave functions, and κ is the relativistic orbital quantum number. The free electrons are to be found between the ion spheres, but are also penetrating the ion spheres to satisfy equation (2.1). In the ion sphere model, the interactions between the positive ions are neglected. Consequently, the electron chemical potential obtained from equation (2.1) is equal to the ion chemical potential. When the above system is solved self-consistently, taking into account pressure ionization, one obtains the chemical potential µ-in the 'chemical picture' framework-of the average atom at temperature T and mass density ρ. The average ion charge Z * can be defined as where Z * is also a function of T and ρ.
The STA model then goes further and expands the average atom model on a Saha-Boltzmann distribution of ion charges and super-configurations, each of which has its optimized energy. These expansions and optimizations are intended to yield a description of the spectra in terms of STAs, which is finer than the average atom one. The resulting charge states are then in equilibrium with the same chemical potential µ, and the average charge Z * is conserved. In the following, the expressions 'chemical potential' and 'ion charge' refer to these average quantities µ and Z * .
The STA code is well adapted to the computation of opacities of LTE plasmas and has been verified in many experiments. Figure 1 shows a comparison of an emission spectrum of tungsten embedded in CH plastic, at a temperature T = 80 eV and electronic density N e = 3 × 10 20 e cm −3 . The experiment was performed at the Naval Research Laboratory, Washington, DC [2]. In this and more recent versions of STA, the introduction of configuration interactions (CIs) considerably improves the agreement with experiment.
The STA code generates spectrally resolved opacities (usually 15 000 frequency points) on a logarithmically equidistant grid of temperatures and matter densities. Typically, for the CRASH project, we compute 200 temperatures from 0.03 to 1000 eV, times 200 densities (ρ) from 10 −6 to 100 g cm −3 . For each (T, ρ) point, the code gives the chemical potential µ, the average charge Z * and the average Rosseland and Planck means. These data are used by the MIX code. The present version of the STA code does not compute any equation of state (EOS) data. A recent idea using the virial theorem in the Dirac equation [8] has not yet been implemented.

The MIX code
MIX takes as the input the files created by STA, and an input file describing the mixture and various options. It reads an additional file with data for photon groups (bins) in order to create multi-group average Rosseland and Planck opacities. The MIX extract for each element from the STA records a block of all computed density records for one temperature. It then computes the partial densities for the chosen model in order to retrieve the corresponding records from the opacity database within this temperature block. If needed, it will interpolate between two density records on the spectrally resolved opacities and then compute the group opacities. Then it proceeds to the next temperature.

Generalities
In the following, we consider a volume V of the mixture, of mass M, containing L elements, with N l atoms of element l, their mass being M l . We further define for each element the mass fraction m l = M l /M and the number fraction q l = N l / N l .
The mixture is assumed to be in a plasma state. We do not take into account, in any of the following models, molecular or solid state effects. All the elements are in the atomic or ionic state, as individually computed with the STA code. The four models considered in this work use the same STA databases; they differ only in their approach of equilibrium and consequently in the way they evaluate partial densities. Complete equilibrium consists of the following: • thermal equilibrium: all the elements should be at the same temperature; • mechanical equilibrium: the elements should be at the same pressure; • chemical equilibrium: the elements should have the same chemical potential; • electrical equilibrium: in order to keep the overall plasma neutrality, the free electron density N e should balance the positive ion density Within the average atom framework the chemical potential is set to ensure the neutrality of the Wigner-Seitz ion sphere according to relation (2.1). Thus the chemical equilibrium and the electrical equilibrium are bound together. The basic difference between the models consists in the way they satisfy-or not-these requirements. From these, one deduces the partial densities of the elements constituting the mixture. The corresponding spectrally resolved opacities are then retrieved at these partial densities and added to yield the group average opacities (see section 4.1).

M1: the 'whole volume partial pressure' model
The simplest model consists in supposing that each and every element occupies the entire volume V. This model has been widely used for astrophysical plasmas, which are usually of very low density [9]. In this model, only the thermal equilibrium is enforced. Assuming that the components of the mixture are ideal gases, each component has a partial pressure proportional to their number fraction, following Dalton's law: In our notation, the partial densities ρ l are where m l is the mass fraction of element l and ρ is the density of the mixture. The partial densities are additive: This model assumes that the various elements do not disturb, or interact with, each other. This would be valid at very low density, where the mean free path between atomic collisions is very much larger than the atomic radii. For these plasmas, the free electron density is usually a given, caused by a very dominant element such as hydrogen, and the ion charges of the other elements-often considered as traces-are computed accordingly, so the electrical equilibrium is enforced. At finite temperature, and at the moderate to high densities considered in this work, we cannot assume that each component, or their mixture, behaves like an ideal gas. Consequently, this model is, for plasmas considered in this work, not TC, because the elements are not at the same pressure. In this model, the partial densities are independent of the temperature, are strictly proportional to the total density of the mixture and are just given by equation (3.1). We mention this model here and display some of its results in order to quantify the difference with more realistic models.

M2: the modified partial density mix model
In model M1, the elements do not interact, and only the thermal and electrical equilibria are enforced. In model M2, as in M1, we still suppose that the partial densities are additive, but the constraint of equal chemical potential replaces the mass fraction condition of equation ( This introduces some interaction between atoms. This is an ad hoc modification. It is tentative to recover the chemical equilibrium in addition to the thermal equilibrium common to all models. This model leaves the question of partial pressure unanswered. We do not have an additional theoretical justification for it. As in the case of M1, we show M2 results in order to quantify the differences. We will show that in some cases, the results from M2 are actually worse than those of M1. Thus we show that enforcing arbitrarily incomplete equilibrium conditions may lead to wrong results. Equation (3.4) is solved iteratively by the Brent method [10], which is an accelerated bisection method. Figure 2 shows the chemical potential for the elements H, C, O and N at T = 30 eV as a function of density. For any given µ, the partial densities are easily found or interpolated.

M3: partial volumes, equal pressure and chemical potential
We now assign to each element a partial volume V l , part of the total volume V, and require equal chemical potential to yield equilibrium. We have It was shown [11] that the assumption of equilibrium with partial volumes implies equal pressure even in the case of non-ideal gases. This fact adds mechanical equilibrium to the thermal equilibrium present in this model like in all models. Furthermore, it was proven recently [12] that this model minimizes free energy. Equation (3.5) is one equation with L unknowns. In order to solve it, additional information, or constraint, is required. For model M3, we will require equality of chemical potential. We have then also for this model chemical equilibrium, and according to the remark in section 3.1, also electrical equilibrium. Thus all the equilibrium conditions mentioned in section 3.1 are enforced, making this model TC. We have now, for each temperature T and mixture density ρ, one equation with the unknown µ: The STA databases give tables of µ l (ρ) for each element at any given temperature T, as shown in figure 2. These tables need to be inverted in order to obtain ρ(µ). Equation (3.6) is then solved by the Brent method [10]. This procedure supposes that the inverted tables are single valued, i.e. that STA tables µ l (ρ) are monotonic. It may so happen that because of numerical errors this condition is locally not satisfied. In order to deal with this problem, we introduced in MIX a binomial smoothing process of the order of 4 [10]. An example of smoothed values is plotted in figure 2.

M4: partial volumes, equal free electron densities
In this model, as in M3, we attempt to solve equation (3.6). However, instead of the equality of chemical potential, we invoke the following argument: the plasma consists of ions of different elements embedded in a sea of free electrons. Therefore, the average charge of each element Z * l has to be consistent with the same electronic density N e . In section 3.4 the electrical equilibrium was achieved as a consequence of the chemical equilibrium. In model M4, it is the reverse. Comparing the results of M3 and M4 will thus show the consistency of the STA and MIX algorithms. In addition, with this model, one can estimate the non-LTE effects; see section 5.
For each element we have N e (l) = N l Z * l . The unknown N e is obtained by inverting the tables Z l * (ρ l ) in a similar manner as the µ tables for model M3. It may have similar numerical problems, and therefore the same smoothing procedure is applied. Equation (3.8) is solved iteratively by the Brent method as above.

Opacity of the mixture
Once the partial densities are found for every element of the mixture by one of the models or algorithms mentioned in the previous section, these partial densities are used for locating the spectral opacity record at the given temperature in their respective database. In most cases, the partial density is not exactly the same as that of the densities used for computing the STA database, and linear interpolation on a logarithmic scale of ρ is applied. The resulting opacity κ (in cm −2 g −1 ) of the mixture is then easily shown to be given by where q l is the number fraction of element l in the mixture and κ l is the opacity of element l at partial density ρ l .

Examples in local thermal equilibrium
In this section, a few examples involving polystyrene, acrylic and polyimide will be shown. figure 3 the partial densities of C and H in CH at T = 50 eV and ρ = 10 −4 g cm −3 are displayed. The models M3 and M4 give the same results. Note that the partial densities for these models can be larger than the total density because they are not additive. For M1 and M2, the partial densities are additive; their sum is equal to 10 −4 g cm −3 . The convergence of M3 and M4 for CH in the same conditions is shown in figure 4. The absolute values of the increments of chemical potential (red line) and of partial densities of C (blue) and H (yellow) are displayed as a function of the convergence step. The convergence of M4 for the same case is somehow slower. It needs 30 iterations (steps) to converge to the same values of the partial densities as for M3 with 13 iterations. The convergence of the density is described here as the increment of the function

Polystyrene (CH) at T = 50 eV. In
In figure 5 the Planck mean opacity for CH at 50 eV versus mixture density is shown. M3 and M4 are nearly identical. M2 (black) is even further from M3 and M4 than M1 (green). As mentioned above, M3 and M4 are TC and could be taken as the reference. Consequently, we consider that the modification introduced in M2 is not an improvement.

Acrylic (C 5 H 8 O 2
). In figure 6, we compare the four models for acrylic at T = 50 eV. The Rosseland mean opacity is displayed as a function of acrylic density. M3 (blue) and M4 (red) are nearly identical. M2 (green) is farther than M1 (black).
In figure 7 the rate of convergence of M2 and M4 is compared for acrylic at T = 50 eV and density ρ = 10 −4 g cm −3 . The convergence of M4 is measured with the function f(ne) defined in equation (4.2). It appears that convergence on electronic density is slower than that on the chemical potential (M2 and M3). In both cases, the convergence of C and O is identical.  In figure 9, the mean Planck opacity of polyimide at 31 eV is shown as a function of density. M3 (blue line) and M4 (red line) are completely superimposed and not distinguishable. Like in the cases of polystyrene and acrylic, M1 (green line) is nearer to the latter than M2 (black).

Generalities
It was mentioned in section 3 that the model M4 lends itself to be applied when the plasma is not in non-LTE. In atomic spectroscopy jargon, non-LTE describes the state of an atomic system in quasi-static equilibrium in which the level populations do not follow a Saha-Boltzman distribution at the electron temperature. This can happen when the plasma is thin enough, so that the emitted photons escape the plasma while an external source keeps the temperature constant-on a time scale longer than the average spontaneous radiation lifetime. It has been shown [13] that for an atomic system with many levels, not too far from LTE, the populations can be approximated by a Saha-Bolzman distribution at another temperature, called the ionization temperature T z . For testing M4, we use the algorithm of the RADIOM model [13]. This model defines the ionization temperature T z , such that The necessary ingredients for this algorithm are the charge state Z * (for each element), the electron density N e and the electron temperature T e . These values are available or readily computed for model M4. The RADIOM model may also take into account an ambient radiation field, which is not available in the MIX code. This means that our evaluation of non-LTE effects is valid only for optically thin plasmas. In this work, we give only an estimate of the non-LTE effects. The correct use of the RADIOM model is to implement it in a hydrodynamic code, taking into account the ambient radiation field in every cell at each time step. Then the non-LTE effects change the dynamical evolution of the plasma. Figure 10 shows the ionization temperature T z for polyimide as obtained in model M4 and the RADIOM model in optically thin plasma (i.e. no ambient radiation). The change from T e to T z is most dramatic at low density and high temperature. The curves correspond to different temperatures T e . The black line is for T e = 50 eV, the blue for 100 eV, the green for 150 eV and the red for 200 eV. Figure 11 shows the ratio of the Planck mean opacity computed with the model M4 and the non-LTE option, to the regular M4 (or M3) LTE result. The curves correspond to different temperatures. The black line is for T e = 50 eV, the blue for 100 eV, the green for 150 eV and the red for 200 eV. The effect appears to be quite dramatic at low density and high temperature. However, in a real hydrodynamic simulation, there is always a radiation field that tends to  moderate the non-LTE effects. In addition, the role of Planck opacity in radiation transfer is usually connected to emissivity. The latter involves a correction to the source function. The resulting emissivity is displayed in figure 12.

Non-local thermal equilibrium effects on polyimide
The coloured lines correspond to the same temperature as in figure 11. It can be seen that emissivity can be affected up to a factor of 2 for temperatures in the range 150-200 eV in the density range of 1-10 mg cm −3 .

Discussion
In this work, we compared four models for evaluating the partial densities of the elements composing a mixture. Our STA code gives spectrally resolved opacities for each element on a grid of temperature and densities. The MIX code computes the partial densities using the STA databases. Our procedure assumes that the elements are in the atomic or ionic state.
Model M1 is the simple 'partial density mix' model, in which these partial densities are simply proportional to the mass fraction of each element in the mixture. This model has often been used for astrophysical plasmas [9]. It is also called the 'partial pressure model' because it assumes that the mixture and its components are ideal gases following Dalton's law of mass action. It does not take into account any interaction between atoms. Consequently, for non-ideal gases or plasmas, it is not TC. A modification of the latter is model M2, where we replace the mass fraction of the elements by the requirement of the equality of chemical potential, in an attempt to make the model more realistic. Model M3 is based on the idea of partial volumes, equal pressure and equality of chemical potential. It was shown that this model is TC and minimizes free energy. Because STA is based on an average atom model, the chemical equilibrium implies also electrical equilibrium. Model M4 is also based on partial volumes and equal pressure, but imposes equality of a common free electron density. By the same argument it also implies chemical equilibrium. It is also TC and minimizes free energy. Since M4 is based on electronic density, it can be easily be modified to give estimates of non-LTE effects through the use of the RADIOM model.
We show through some examples that M3 and M4 give essentially identical results in LTE. This shows that our algorithms are consistent. However, the convergence of M4 is slower than that of M3. The difference between the models M1 and M3 can amount to a factor of 10 on the Planck mean, and to even more between M2 and M3. Indeed, although M2 was modified in order to take into account atomic interactions, it turned out worse than M1.
In summary, we justify the use of model M3 in LTE for the group opacities used by the CRASH project. It is TC and minimizes free energy. M4 can be used to give an estimate of non-LTE effects, but its convergence is slower. Other models may yield errors on Rosseland or Planck opacities by one or two orders of magnitude.