Equations of Radiation Hydrodynamics

The basic features of the model used for simulation of the laser-matter interaction are described (Povarnitsyn et al., Reference Povarnitsyn, Andreev, Levashov, Khishchenko and Rosmej2012a). The model takes account of 1D hydrodynamic motion of substance, laser energy absorption, two-temperature nonequilibrium states for electron and ion subsystems, electron thermal conductivity and radiation transport in diffusion approximation. The evolution of material parameters is described using the conservation of mass, momentum and energy of electron and ion subsystems in a single-fluid two-temperature Lagrangian form:

(1)$$\displaystyle{{\partial \lpar 1/{\rm \rho} \rpar } \over {\partial t}} - \displaystyle{{\partial u} \over {\partial m}}=0\comma \;$$

(2)$$\displaystyle{{\partial u} \over {\partial t}}+\displaystyle{{\partial \lpar P_i+P_e \rpar } \over {\partial m}}=0\comma \;$$

(3)$$\eqalign{\displaystyle{{\partial e_e } \over {\partial t}}+P_e \displaystyle{{\partial u} \over {\partial m}} & = - {\rm \gamma} _{ei} \lpar T_e - T_i \rpar /{\rm \rho} +Q_L /{\rm \rho}+\displaystyle{\partial \over {\partial m}}\left({{\rm \rho} {\rm \kappa} _e \displaystyle{{\partial T_e } \over {\partial m}}} \right)- \displaystyle{{\partial S} \over {\partial m}}\comma \; }$$

(4)$$\displaystyle{{\partial e_i } \over {\partial t}}+P_i \displaystyle{{\partial u} \over {\partial m}}={\rm \gamma} _{ei} \lpar T_e - T_i \rpar /{\rm \rho}.$$

Here ρ is the density; t is the time; m is the mass coordinate, dm = ρdz with the space coordinate z; u is the velocity; P, e, T are the pressure, specific energy and temperature, respectively. Indices e and i stand for electron and ion species. The energy exchange between electrons and ions is described by the corresponding term with the coupling factor γ_ei(ρ, T _e, T _i). The electron heat transfer is specified with the aid of the thermal conduction coefficient κ_e(ρ, T _e, T _i). Optical and transport properties of Al and Ag plasma are taken into account by means of wide-range models (Povarnisyn et al., 2012b; Veysman et al., Reference Veysman, Agranat, Andreev, Ashitkov, Fortov, Khishchenko, Kostenko, Levashov, Ovchinnikov and Sitnikov2008). For CH₂ material, we use a similar approach interpolating between dense preionized state (Drude model) and Spitzer ideal plasma (Spitzer & Härm, Reference Spitzer and Härm1953). The key parameter in these models is the frequency of collisions (electron–phonon, electron–ion, and electron–electron). In the metal state, it can be written as follows

(5)$${\rm \nu} _{\rm met}=A_1\, k_B T_i /\hbar+A_2\, k_B T_e^2 /\lpar T_F \hbar \rpar.$$

Here k _B and ħ are the Boltzmann and Planck constants, respectively, and T _F is the Fermi temperature. The frequency is limited from above by the electron free-path between ions so that

(6)$${\rm \nu} _{\rm max}=\displaystyle{{A_3 } \over {r_0 }}\sqrt {{\rm \upsilon} _F^2+k_B T_e /m_e }\comma \;$$

where υ _F is the Fermi speed of electrons, r ₀ is the interatomic distance and m _e is the electron mass.

For hot states T _e ≫ T _F with not very high densities, the plasma model (Spitzer & Härm, Reference Spitzer and Härm1953) is relevant

(7)$${\rm \nu} _{\rm pl}=\displaystyle{4 \over 3}\sqrt {2{\rm \pi} } \displaystyle{{\left\langle Z \right\rangle n_e e^4 \Lambda } \over {\sqrt {m_e } \lpar k_B T_e \rpar ^{3/2} }}\comma \;$$

where 〈Z〉 is the mean charge of ions, n _e is the electron concentration, e is the electron charge, and Λ is the Coulomb logarithm.

The electron thermal conductivity in the metal is calculated according to the Drude formalism as follows

(8)$${\rm \kappa} _{\rm met}=\displaystyle{{{\rm \pi} ^2\, k_B^2 n_e } \over {3m_e {\rm \nu} _{\rm eff} }}T_e\comma \;$$

where ν _eff = min(ν _met, ν _max). The hot plasma limit is

(9)$${\rm \kappa} _{\rm pl}=\displaystyle{{16\sqrt 2\, k_B \lpar k_B T_e \rpar ^{5/2} } \over {{\rm \pi} ^{3/2} \left\langle Z \right\rangle e^4 \sqrt {m_e } \Lambda }}.$$

An interpolation between (8) and (9) in the vicinity of the Fermi temperature gives us a wide-range expression for the thermal diffusivity,

(10)$${\rm \kappa}={\rm \kappa} _{\rm pl}+\lpar {\rm \kappa} _{\rm met} - {\rm \kappa} _{\rm pl} \rpar \exp\left[{ - A_4 T_e /T_F } \right].$$

The dimensionless parameters A ₁, A ₂, A ₃, A ₄ are used to adjust the model to available theoretical and experimental data. For the thermal conduction of Ag, the parameters are A ₁^t = 0.7, A ^t₂ = 0.5, A ^t₃ = 1.0, A ^t₄ = 3.0. In Figure 1, we compare the results of the interpolation (10) with the theoretical data (Apfelbaum, Reference Apfelbaum2011) for normal and 0.1 of normal density. Appropriate coefficients of the model for Al are A ^t₁ = 2.95, A ^t₂ = 0.5, A ^t₃ = 0.16, A ^t₄ = 1.2, see Figure 2 in paper byPovarnitsyn et al. (Reference Povarnitsyn, Andreev, Apfelbaum, Itina, Khishchenko, Kostenko, Levashov and Veysman2012b).

Fig. 1.

(Color online) The thermal conductivity coefficient of Ag. Comparison of interpolation formula (10) (lines) with available theoretical data (Apfelbaum, Reference Apfelbaum2011) (signs) in cases of normal and 0.1 of normal densities, T _i = T _e.

Fig. 2.

(Color online) The coupling factor of Ag according to the model (11) (lines) and theoretical results (Lin et al., Reference Lin, Zhigilei and Celli2008) (signs) for normal density.

The coupling factor is written in the form

(11)$${\rm \gamma} _{ei}=\displaystyle{{3k_B m_e } \over {m_i }}n_e {\rm \nu} _{\rm eff}\comma \;$$

where m _i is the ion mass, and ν _eff = min(ν_met, ν _max, ν _pl) is given by expressions (5), (6), (7). For Ag, the dimensionless parameters are A ^g₁ = 10.0, A ^g₂ = 10.0, A ^g₃ = 0.4. The coupling factor (11) and results of first-principle calculations (Lin et al., Reference Lin, Zhigilei and Celli2008) are presented in Figure 2 for Ag. Appropriate coefficients for Al can be found in paper by Povarnitsyn et al. (Reference Povarnitsyn, Andreev, Apfelbaum, Itina, Khishchenko, Kostenko, Levashov and Veysman2012b).

The laser energy absorption (bremsstrahlung process) is taken into account with the aid of heat source term Q _L (t, z). To find it we solve the Helmholtz wave equation for the laser electric field with proper boundary conditions (Povarnitsyn et al., Reference Povarnitsyn, Andreev, Levashov, Khishchenko and Rosmej2012a)

$$\displaystyle{{\partial ^2 E} \over {\partial z^2 }}+\displaystyle{{{\rm \omega} _L^2 } \over {c^2 }}{\rm \epsilon} \lpar z\rpar E=0.$$

Here the normal incidence is supposed, E is the laser electric field envelope, ω_L is the laser frequency, c is the speed of light, ε(ρ, T _e, T _i) is the complex dielectric function. The resulting laser heating source is

$$Q_L \lpar t\comma \; z\rpar =I\lpar t\rpar \displaystyle{{{\rm \omega} _L } \over c} {\rm Im}\lcub {\rm \epsilon} \lpar t\comma \; z\rpar \rcub \vert E\lpar t\comma \; z\rpar /E_L \lpar t\rpar \vert ^2\comma \;$$

where the incident laser field amplitude is $E_L \lpar t\rpar =\sqrt {8{\rm \pi} I\lpar t\rpar /c} .$

The radiation transport is taken into account by the integral over the spectrum flux $S=\,\vint_{\rm \omega} \, S_{\rm \omega} {\rm d}{\rm \omega} $. The equations of radiation transport are written in diffusion approximation:

(12)$$\displaystyle{{\partial S_{\rm \omega} } \over {\partial z}}=4{\rm \pi} j_{\rm \omega} - {\rm \kappa} _{\rm \omega} cU_{\rm \omega}\comma \;$$

(13)$$S_{\rm \omega}=- \displaystyle{c \over {3{\rm \kappa} _{\rm \omega} }}\displaystyle{{\partial U_{\rm \omega} } \over {\partial z}}\comma \;$$

(14)$$\mathop {\left. {S_{\rm \omega} } \right\vert }\nolimits_\Omega=\mp \displaystyle{{cU_{\rm \omega} } \over 2}\comma \;$$

where S _ω is the radiation flux density, U _ω is the radiant energy density, j _ω and κ _ω are the coefficients of emissivity and absorption, respectively, and Ω is the external boundary of plasma.

Equations of state (EOSs) for Al, CH₂, and Ag determine the relations P _e(ρ, T _e ), P _i(ρ, T _i ), e _e(ρ, T _e ), e _i(ρ, T _i) and are necessary for completeness of system (1)–(4). These EOSs are based on the analytical expression of the Helmholtz free energy that has a form ℱ(ρ, T _i, T _e) = ℱ_i(ρ, T _i) + ℱ_e(ρ, T _e), and is composed of two parts. The first item stands for the ionic part ℱ_i(ρ, T _i) = ℱ_c(ρ) + ℱ_a(ρ, T _i), and consists of the electron-ion interaction term ℱ_c (calculated at T _i = T _e = 0 K) and the contribution of the thermal motion of ions ℱ_a. The analytical form of ℱ_i has different expressions for solid and fluid phases (Khishchenko, Reference Khishchenko2008). The tables of thermodynamic parameters are calculated to take into account phase transitions and metastable regions (Khishchenko et al., Reference Khishchenko, Tkachenko, Levashov, Lomonosov and Vorob'ev2002; Levashov & Khishchenko, Reference Levashov and Khishchenko2007; Oreshkin et al., Reference Oreshkin, Baksht, Ratakhin, Shishlov, Khishchenko, Levashov and Beilis2004). The second term ℱ_e(ρ, T _e) describes the thermal contribution of electrons calculated with the aid of the Thomas–Fermi model (Shemyakin et al., Reference Shemyakin, Levashov, Obruchkova and Khishchenko2010).

Absorption and Emission Coefficients

As was demonstrated in Povarnitsyn et al. (Reference Povarnitsyn, Andreev, Levashov, Khishchenko and Rosmej2012a) the radiation transport is a key effect in the dynamics of thin films illuminated by nanosecond prepulses with intensity ≳10¹¹ W/cm². Transformation of monochromatic laser energy into a wide spectrum of photons results in formation of noticeable re-emission fluxes on both sides of the film. We focus here on the detailed description of the procedure for calculation of absorption and emission coefficients. These coefficients depend on ion concentrations. To obtain concentrations x _ks of k-fold ionized atoms or molecules in the state s the system of level kinetic equations is solved in quasi-stationary approximation:

(15)$$\eqalign{\displaystyle{{{\rm d}x_{ks} } \over {{\rm d}t}} & =- x_{ks} \left({\sum\limits_{s^{\prime}} \, R_{ks \to k - 1\comma s^{\prime}}+\sum\limits_{s^{\prime}} \, I_{ks \to k+1\comma s^{\prime}} } \right. +\left. {\sum\limits_{s^{\prime}} \, {\rm \alpha} _{ks \to ks^{\prime}} } \right)\cr & \quad +\sum\limits_{s^{\prime}} \, x_{k+1\comma s^{\prime}} R_{k+1\comma s^{\prime} \to k\comma s} +\sum\limits_{s^{\prime}} \, x_{k - 1\comma s^{\prime}} I_{k - 1\comma s^{\prime} \to k\comma s}\cr & \quad+\sum\limits_{s^{\prime}} \, x_{ks^{\prime}} {\rm \alpha} _{ks^{\prime} \to ks}=0.}$$

In Eq. (15) the rate of recombination from the ion state ks to the state k − 1, s′ is given by

(16)$$\eqalign{R_{ks \to k - 1\comma s^{\prime}} & ={\rm \alpha} ^{ir} \lpar ks \to k - 1\comma \; s^{\prime}\rpar \cr & \quad +{\rm \alpha} ^{\,phr} \lpar ks \to k - 1\comma \; s^{\prime}\rpar \cr & \quad +{\rm \alpha} ^{dc} \lpar ks \to k - 1\comma \; s^{\prime}\rpar \comma \;}$$

where α^ir(ks → k − 1, s′), α ^phr(ks → k − 1, s′), and α^dc(ks → k − 1, s′) are the rates of three-body recombination, photo-recombination and dielectronic capture, respectively. The rate of ionization from ks to k + 1, s′ state is given by

(17)$$\eqalign{I_{ks \to k+1\comma s^{\prime}} & ={\rm \alpha} ^{ii} \lpar ks \to k+1\comma \; s^{\prime}\rpar \cr & \quad +{\rm \alpha} ^{\,phi} \lpar ks \to k+1\comma \; s^{\prime}\rpar \cr & \quad +{\rm \alpha} ^{ai} \lpar ks \to k+1\comma \; s^{\prime}\rpar \comma \;}$$

where αⁱⁱ(ks → k + 1, s′), α^phi(ks → k + 1, s′), and α^ai(ks → k + 1, s′) are the rates of impact ionization, photo-ionization and auto-ionization, respectively. The rates of transitions from ks to ks′ without a change of ionization stage are written as follows

(18)$${\rm \alpha} _{ks \to ks^{\prime}}=\left\{{\matrix{ {{\rm \alpha} ^{ex} \lpar ks \to ks^{\prime}\rpar +{\rm \alpha} ^{abs} \lpar ks \to ks^{\prime}\rpar \comma \; } \cr \hfill {{\rm if}\, E_{ks}\lt E_{ks^{\prime}}\comma \; } \cr {{\rm \alpha} ^{dex} \lpar ks \to ks^{\prime}\rpar +{\rm \alpha} ^{em} \lpar ks \to ks^{\prime}\rpar \comma \; } \cr \hfill {{\rm if}\, E_{ks}\gt E_{ks^{\prime}}\comma \; }}} \right.$$

where E _ks is the energy of ion state ks, α^ex(ks → ks′) and α^dex(ks → ks′) are the rates of excitation and de-excitation by electron impact, and α^abs(ks → ks′), α^em(ks → ks′) are the rates of radiative excitation (absorption) and radiative emission, respectively. For calculation of these rates the Regemorter, Lotz and Kramers approximations are used (Nikiforov et al., Reference Nikiforov, Novikov and Uvarov2005; Sobel'man et al., Reference Sobel'man, Vainshtein and Youkov1995).

The number of ion states that must be taken into account can be very large. To reduce the system of kinetic equations we apply a radiative unresolved spectra atomic model (RUSAM) (Kim et al., Reference Kim, Novikov, Dolgoleva, Koshelev and Solomyannaya2012; Novikov et al., Reference Novikov, Koshelev, Solomyannaya and Fortov2010). Knowing the concentrations x _ks the tables of spectral absorption coefficients and emissivities are obtained with the code THERMOS (Nikiforov et al., Reference Nikiforov, Novikov and Uvarov2005).

The concentrations x _ks are defined by the radiation field U _ω through the rates of the radiation processes. The photo-ionization of an electron with quantum numbers nℓ (n is the principal quantum number, ℓ is the orbital quantum number) in Kramers approximation for Q _ks-configuration with occupation numbers N ^ks_nℓ has the form:

(19)$${\rm \alpha} _{n\ell }^{\,phi} \lpar ks \to k+1\comma \; s^{\prime}\rpar =N_{n\ell }^{ks} u_{n\ell } \vint\limits_{{\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}} /{\rm \theta} }^\infty \, W\lpar {\rm \xi} \rpar {\rm \xi} ^{ - 1} {\rm d}{\rm \xi}\comma \;$$

where

$$u_{n\ell }=4.45 \times 10^{10} \, Z_k \displaystyle{{\vert {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}} \vert ^{3/2} } \over {2n^2 }}\comma \;$$

ɛ_ks;k+1,s′ = E _ks − E _k+1,s′ is the photo-ionization threshold in atomic units, Z _k is the effective charge, θ is the temperature in atomic units, W(ξ) = C _ξU _ω/ω³ is the radiation energy density in dimensionless units, C _ξ is the dimension factor. For instance, for Planck's radiation field W(ξ) = [exp(ξ) − 1]⁻¹. The photo-recombination rate to the level nℓ in configuration Q _ks with occupation numbers N ^ks_nℓ is given by

(20)$$\eqalign{{\rm \alpha} _{n\ell }^{\,phr} \lpar k+1\comma \; s^{\prime} \to ks\rpar & =C_{n\ell } u_{n\ell } N_{n\ell }^{ks} \ \times \vint\limits_{{\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}} /{\rm \theta} }^\infty \, \displaystyle{{1+W\lpar {\rm \xi} \rpar } \over {{\rm \xi} \, {\rm exp}\lpar {\rm \xi} \rpar }}{\rm d}{\rm \xi}\comma \; }$$

where

$$C_{n\ell }=0.704\, \displaystyle{{N_e } \over {N_A }}\displaystyle{{g_{ks} } \over {g_{k+1\comma s^{\prime}} }}\displaystyle{{\exp\lpar {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}} /{\rm \theta} \rpar } \over {{\rm \theta} ^{3/2} }}\comma \;$$

N _A is the Avogadro constant, N _e is the electron density, g _ks is the statistical weight for the state ks. For given photon energy intervals [ω_i, ω_i+1] (i = 1, N _g), absorption in lines is defined by

(21)$$\eqalign{{\rm \alpha} _{ks\comma ks^{\prime}}^{abs} & =3.2 \times 10^{10} \sum\limits_{i=1}^{N_g } \, \displaystyle{{\mathop {\overline {g_f } }\nolimits_i \lpar ks\comma \; ks^{\prime}\rpar } \over {g_{ks} }}\mathop {\overline {{\rm \omega} _i } ^2\lpar ks\comma \; ks^{\prime}\rpar }\nolimits \times W\left({\displaystyle{{\overline {{\rm \omega} _i } \lpar ks\comma \; ks^{\prime}\rpar } \over {\rm \theta} }} \right)\comma \; }$$

where $\overline {{\rm \omega} _i } \lpar ks\comma \; ks^{\prime}\rpar $ is the center of line group, and $\mathop {\overline {g_f } }\nolimits_i \lpar ks\comma \; ks^{\prime}\rpar $ is the oscillator strength averaged over interval [ω_i, ω_i+1]. The line emission rate is given by

(22)$$\eqalign{{\rm \alpha} _{ks^{\prime}\comma ks}^{em} & =3.2 \times 10^{10} \sum\limits_{i=1}^{N_g } \, \displaystyle{{\mathop {\overline {g_f } }\nolimits_i \lpar ks\comma \; ks^{\prime}\rpar } \over {g_{ks^{\prime}} }}\mathop {\overline {{\rm \omega} _i } ^2\lpar ks\comma \; ks^{\prime}\rpar }\nolimits \cr & \quad\times \left[{1+W\left({\displaystyle{{\overline {{\rm \omega} _i } \lpar ks\comma \; ks^{\prime}\rpar } \over {\rm \theta} }} \right)} \right].}$$

Using RUSAM method the absorption coefficient in cm⁻¹ for a photon with energy ω in the spectral interval [ω_i, ω_i+1] can be written in the form:

(23)$$\eqalign{{\rm \kappa} _i & ={\rm \rho} \displaystyle{{N_A } \over A}\left[{1 - \exp\lpar - {\rm \omega} _i^{{\rm ff}} /{\rm \theta} \rpar } \right]\left\{{\sum\limits_{ks} \, x_{ks} \sum\limits_{ks^{\prime}} \, {\rm \sigma} _i^{{\rm bb}} \lpar ks\comma \; ks^{\prime}\rpar } \right. \cr & \quad\left. {+\sum\limits_{ks} \, x_{ks} \sum\limits_{k+1\comma s^{\prime}} \, {\rm \sigma} _{ks\semicolon k+1\comma s^{\prime}}^{{\rm bf}} \lpar {\rm \omega} _i^{{\rm bf}} \rpar +{\rm \sigma} ^{{\rm ff}} \lpar {\rm \omega} _i^{{\rm ff}} \rpar } \right\}.}$$

The emissivity in TW/(cm³× eV × steradian) is written as follows:

(24)$$\eqalign{j_i & ={\rm \rho} \displaystyle{{N_A } \over A}\displaystyle{1 \over {C_W }}\left[{\displaystyle{{\mathop {\overline {{\rm \omega} _i } }\nolimits^3 \lpar ks\comma \; ks^{\prime}\rpar } \over {4{\rm \pi} ^3 }}{\rm \alpha} ^2 } \right]\left\{{\sum\limits_{ks^{\prime}} \, x_{ks^{\prime}} \sum\limits_{ks} \, {\rm \sigma} _i^{{\rm bb}} \lpar ks^{\prime}\comma \; ks\rpar } \right. \cr & \quad \left. {+\sum\limits_{k+1\comma s^{\prime}} \, x_{k+1\comma s^{\prime}} \sum\limits_{k\comma s} \, {\rm \sigma} _{k+1\comma s^{\prime}\semicolon ks}^{{\rm fb}} \lpar {\rm \omega} _i^{{\rm bf}} \rpar +\exp\lpar \!- \!{\rm \omega} _i^{{\rm ff}} /{\rm \theta} \rpar {\rm \sigma} ^{{\rm ff}} \lpar {\rm \omega} _i^{{\rm ff}} \rpar } \right\}.}$$

Here A is the atomic weight, C _W = 4.23 × 10⁻³, ω^ff_i = 0.5(ω_i+1 + ω_i) is the center of interval,

$${\rm \omega} _i^{{\rm bf}} = \left\{\matrix{0.5\lpar {\rm \omega} _{i+1}+{\rm \omega} _i \rpar \comma \; \, {\rm if}\, {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}}\lt {\rm \omega}_i\comma \; \hfill \cr {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}}\comma \; \, {\rm if}\, {\rm \omega} _i\lt {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}} \leq {\rm \omega}_{i+1}.} \right.$$

The absorption cross-section in spectral lines can be written in the form:

(25)$${\rm \sigma}_i^{\rm bb} \lpar ks\comma \; ks^{\prime}\rpar = 2{\rm \pi}^2 {\rm \alpha} a_0^2 {1 \over {\rm \omega}_{i+1} - {\rm \omega}_i} {\overline{g_f}_i \lpar ks\comma \; ks^{\prime}\rpar \over g_{ks}}\comma \;$$

and the corresponding emission cross-section is given by

(26)$${\rm \sigma}_i^{\rm bb} \lpar ks^{\prime}\comma \; ks\rpar = 2{\rm \pi}^2 {\rm \alpha} a_0^2 {1 \over {\rm \omega}_{i + 1} - {\rm \omega} _i} {\overline{g_f}_i \lpar ks\comma \; ks^{\prime}\rpar \over g_{ks^{\prime}}}\comma \;$$

where α ≈ 1/137, a ₀ ≈ 5.29 × 10⁻⁹ cm is the Bohr radius. For the photo-ionization cross-section, we use the simple Kramers approximation:

(27)$$\eqalign{{\rm \sigma}_{ks\semicolon k + 1\comma s^{\prime}}^{\rm bf} \lpar {\rm \omega} \rpar &= {64 {\rm \pi} {\rm \alpha} a_0^2 \over 3\sqrt{6}} Z_k {{\rm \varepsilon}_{ks\semicolon k + 1\comma s^{\prime}}^{3/2} \over {\rm \omega}^3} \cr &\quad \times {N_{n\ell}^{ks} \over 2n^2} \left[1 - n\lpar {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}} - {\rm \omega}\rpar \right]\comma \; \, {\rm if}\, {\rm \omega}\gt {\rm \varepsilon}_{ks\semicolon k + 1\comma s^{\prime}}\comma \; }$$

σ^bf_{ks;k+1, s′} (ω) = 0, if ω < ɛ_ks;k+1,s′, with photo-ionization threshold ɛ_ks;k+1,s′ = E _ks − E _{k+1, s′} in atomic units. The photo-recombination cross-section to the state ks has the form:

(28)$$\eqalign{{\rm \sigma} _{k+1\comma s^{\prime}\semicolon ks}^{{\rm fb}} \lpar {\rm \omega} \rpar &= {{{64{\rm \pi} {\rm \alpha} a_0^2 } \over {3\sqrt 6 }}}Z_k {{\displaystyle{{\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}}^{3/2} } \over {{\rm \omega} ^3 }}} \cr & \quad \times {{{2\lpar 2\ell+1\rpar - N_{n\ell }^{k+1\comma s^{\prime}} } \over {2n^2 }}}n\lpar {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}} - {\rm \omega} \rpar \comma \; \, {\rm if} {\rm \omega}\gt {\rm \varepsilon} _{ks\semicolon k+1\comma s^{\prime}}\comma \; }$$

σ^fb_k+1,s′;ks(ω) = 0, if ω < ɛ_ks;k+1,s′. Here $n\lpar {\rm \varepsilon} \rpar = \mathop {\left[{1 + \exp\left({ - {\textstyle{{{\rm \varepsilon}+{\rm \mu} } \over {\rm \theta} }}} \right)} \right]}\nolimits^{ - 1}\!$ is the free electrons distribution function μ is the chemical potential. The bremsstrahlung cross section is calculated in the Born-Elwert approximation

(29)$$\eqalign{{\rm \sigma} ^{{\rm ff}} \lpar {\rm \omega} \rpar = 2.384 \times 10^6 {\textstyle{{{\rm \rho} Z_0^2 {\rm \theta} } \over {A{\rm \omega} ^3 }}}\vint\limits_{{\rm \varepsilon} _0 }^\infty \, n\lpar {\rm \varepsilon} \rpar \left[{1 - n\lpar {\rm \varepsilon} ^{\prime}\rpar } \right]g\lpar {\rm \varepsilon} ^{\prime}\comma \; {\rm \varepsilon} \rpar {\rm d}{\rm \varepsilon}\comma \; }$$

where ɛ′ = ɛ+ ω and

$${\rm g}\lpar {\rm \varepsilon}^{\prime}\comma \; {\rm \varepsilon}\rpar = {\sqrt{3} \over {\rm \pi}} \sqrt{{{\rm \varepsilon}^{\prime} \over {\rm \varepsilon}}} \ln\left({\sqrt{{\rm \varepsilon}^{\prime}} + \sqrt{{\rm \varepsilon}} \over \sqrt{{\rm \varepsilon}^{\prime}} - \sqrt{{\rm \varepsilon}}} \right){1 - \exp \left[- 2{\rm \pi} Z_0 /\sqrt{2{\rm \varepsilon}^{\prime}} \right]\over 1 - \exp\left[- 2{\rm \pi} Z_0 /\sqrt{2{\rm \varepsilon}} \right]}.$$

The radiation field and level kinetics are accounted self-consistently by using an approach with interpolation between pre-calculated tables of spectral parameters.

Spectral absorption coefficients K_ω and emissivities j _ω are calculated using a collisional-radiative equilibrium model in quasi-stationary approximation for an optically-transparent plasma (radiation field U _ω = 0) and for plasma in thermodynamic equilibrium with Planck's radiation field (U _ω = U _ω^P) (Novikov & Solomyannaya, Reference Novikov and Solomyannaya1998):

(30)$$U_{\rm \omega} ^P=\displaystyle{{60} \over {c{\rm \pi} ^4 }}{\rm \sigma} \displaystyle{{{\rm \omega} ^3 } \over {\exp\lpar {\rm \omega} /T\rpar - 1}}\comma \;$$

where σ = 1.028 × 10⁻⁷ TW/cm²/eV⁴ is the Stefan–Boltzmann constant (here ω and T are measured in eV, c ≈ 3 × 10¹⁰ cm/s). It should be noted that the transparent plasma approximation comes to thermodynamic equilibrium at high densities (due to collisional processes) and gives the coronal equilibrium at low densities.

The ratio between calculated radiation field and Planck's radiation field is used as an interpolation parameter:

(31)$${\rm \zeta} \lpar z\comma \; t\rpar =\displaystyle{{\vint\limits_0^\infty \, U_{\rm \omega} \lpar z\comma \; t\rpar {\rm d}{\rm \omega} } \over {\vint\limits_0^\infty \, U_{\rm \omega} ^P {\rm d}{\rm \omega} }}.$$

Applying parameter ζ, the absorption coefficients K_ω and emissivities j _ω are calculated as follows:

(32)$$\matrix{ {{\rm \kappa} _{\rm \omega}=\lpar 1 - {\rm \zeta} \rpar {\rm \kappa} _{\rm \omega} ^{\lpar U_{\rm \omega}=0\rpar }+{\rm \zeta} {\rm \kappa} _{\rm \omega} ^{\lpar U_{\rm \omega}=U_{\rm \omega} ^P \rpar }\comma \; } \cr {\,j_{\rm \omega}=\lpar 1 - {\rm \zeta} \rpar j_{\rm \omega} ^{\lpar U_{\rm \omega}=0\rpar }+{\rm \zeta} j_{\rm \omega} ^{\lpar U_{\rm \omega}=U_{\rm \omega} ^P \rpar }\comma \;}}$$

where ${\rm \kappa}_\omega^{(U_\omega = 0)}, j_\omega^{(U_\omega = 0)}$ are the absorption coefficients and emissivities for optically-transparent plasma, and ${\rm \kappa}_\omega^{(U_\omega = U_\omega^P)}, j_\omega^{(U_\omega = U_\omega^P)}$ are the absorption coefficients and emissivities for plasma in thermodynamic equilibrium.

After that two types of mean group absorption coefficients are calculated by using Rosseland R(ξ) = ξ⁴exp(−ξ)/[1 − exp(−ξ)]², and Planck P(ξ) = ξ³exp(−ξ)/[1 − exp(−ξ)] weight functions:

(33)$${\rm \kappa} _{g_i }^R = \displaystyle{{\vint_{{\rm \xi} _i }^{{\rm \xi} _{i + 1} } \, R\lpar {\rm \xi} \rpar {\rm d}{\rm \xi} } \over {\vint_{{\rm \xi} _i }^{{\rm \xi} _{i + 1} } \, R\lpar {\rm \xi} \rpar /{\rm \kappa} _{\rm \omega} {\rm d}{\rm \xi} }}\comma \; {\rm \kappa} _{g_i }^P = \displaystyle{{\vint_{{\rm \xi} _i }^{{\rm \xi} _{i + 1} } \, P\lpar {\rm \xi} \rpar {\rm \kappa} _{\rm \omega} {\rm d}{\rm \xi} } \over {\vint_{{\rm \xi} _i }^{{\rm \xi} _{i + 1} } \, P\lpar {\rm \xi} \rpar {\rm d}{\rm \xi} }}.$$

The group emission coefficient is calculated as follows:

(34)$$\,j_{g_i } = \displaystyle{{\vint_{{\rm \xi} _i }^{{\rm \xi} _{i + 1} } \, j_{\rm \omega} {\rm d}{\rm \xi} } \over {{\rm \xi} _{i + 1} - {\rm \xi} _i }}.$$

Here ξ = ω/T, and the index g _i denotes the heat radiation frequency range ω_i ≤ ω ≤ ω_i+1 for i = 1, … N _g, where N _g is the number of groups (N _g = 96 in the presented below results).

Now, the system (12)–(14) can be written in a group approximation and transformed to the diffusion equation with appropriate boundary conditions so that

(35)$$\displaystyle{\partial \over {\partial z}}\left({ - \displaystyle{c \over {3{\rm \kappa} _{g_i }^R }}\displaystyle{{\partial U_{g_i } } \over {\partial z}}} \right)= 4{\rm \pi} j_{g_i } - {\rm \kappa} _{g_i }^P cU_{g_i }.$$

The tables of coefficients for CH₂ are obtained by applying a mixing procedure of the tables for pure elements C and H. At given temperature the equality of chemical potentials (electron densities) of the mix components is suggested (Nikiforov et al., Reference Nikiforov, Novikov and Uvarov2005).