Radiative transfer of hydrogen lines from supernova remnant shock waves: contributions of 2s-state hydrogen atoms

TL;DR

This paper investigates how the 2s-state hydrogen atoms influence radiative transfer and emission line profiles in supernova remnant shock waves, revealing effects on line asymmetry and temperature diagnostics.

Contribution

It introduces a detailed model of hydrogen line transfer considering the 2s-state, highlighting its impact on line profiles and interpretation of shock conditions.

Findings

Hα line profiles become asymmetric due to 2s-state effects.

The broad-to-narrow Hα ratio varies with upstream ionization degree.

Line scattering occurs in dense H I clouds, affecting emission diagnostics.

Abstract

Radiative transfer in hydrogen lines in supernova remnant (SNR) shock waves is studied taking into account the population of the hydrogen atom 2s-state. Measurements of Balmer line emission, especially of H~$\alpha$, are often relied upon to derive physical conditions in the SNR shock. On the other hand, Lyman series photons, especially Ly~$\beta$, are mostly absorbed by upstream hydrogen atoms. As a result, atoms are excited to the 3p state, and then emit H~$\alpha$ by the spontaneous transition from 3p to 2s. Thus, the nature of H~$\alpha$ depends on how many Ly~$\beta$ photons are converted to H~$\alpha$ photons. Moreover, the Balmer lines can be scattered by the 2s-state hydrogen atoms, which are excited not only by collisional excitation but also by the Lyman-Balmer conversion. It is shown for example that the H~$\alpha$ photons are scattered if the shock propagates into an H~$_{\rm I}$ cloud with a density of $\sim30~{\rm cm^{-3}}$ and a size of $\sim 1$~pc. We find that the line profile of H~$\alpha$ becomes asymmetric resulting from the difference between line centre frequencies among the transitions from 3s to 2p, from 3p to 2s and from 3d to 2p. We also find that the broad-to-narrow ratio of H~$\alpha$, which is often used to estimate the ion-electron temperature equilibrium, varies at most $\simeq 10$ per cent depending on the ionization degree of the upstream medium because of incomplete conversion of Lyman lines to Balmer lines.

Equations

$$n_{\mathrm{H},1s}{C}_{1s,2s}+n_{\mathrm{H},3p}{A}_{3p,2s}-n_{\mathrm{H},2s}{A}_{2s,1s}=0,$$

$$n_{\mathrm{H},1s}\left(C_{1s,3p}+{\int }_0^{\infty }\frac{4\pi {\sigma }_{\nu }^{1s,3p}}{h\nu }{I}_{\nu }\mathrm{d}\nu \right)-n_{\mathrm{H},3p}\left(A_{3p,1s}+A_{3p,2s}\right)=0,$$

$$I_{\nu }=S_{\nu }(1-{\mathrm{e}}^{-{\tau }_{\nu }})=\frac{\frac{h\nu }{4\pi }{A}_{3p,1s}{n}_{\mathrm{H},3p}}{{\sigma }^{\prime }{n}_{\mathrm{H},1s}}(1-{\mathrm{e}}^{-{\tau }_{\nu }}),$$

$$n_{\mathrm{H},2s}=\frac{C_{1s,2s}}{A_{2s,1s}}\left[1+\frac{A_{3p,2s}}{{\mathrm{e}}^{-{\tau }_0}{A}_{3p,1s}+A_{3p,2s}}\frac{C_{1s,3p}}{C_{1s,2s}}\right]n_{\mathrm{H},1s},$$

$$k_0(\mathrm{H}\alpha )$$

$$n_{\mathrm{p}}(0)$$

$$u_2$$

$$k_{\mathrm{B}}{T}_{\mathrm{p}}$$

$$T_{\mathrm{e}}$$

$${\mu }^{\prime }=1-\left(1-{\mu }_{\odot }^{\prime }\right)\frac{\beta -\frac{m_{\mathrm{e}}}{m_{\mathrm{p}}}}{{\mu }_{\odot }^{\prime }+(1-{\mu }_{\odot }^{\prime })\beta -\frac{m_{\mathrm{e}}}{m_{\mathrm{p}}}},$$

$$\frac{\partial n_{\mathrm{H},1s}^{\mathrm{N}}}{\partial z}=-n_{\mathrm{H},1s}^{\mathrm{N}}\frac{C_{\mathrm{I},\mathrm{N}}+C_{\mathrm{CX},\mathrm{N}}}{V_{\mathrm{sh}}},$$

$$C_{\mathrm{I},\mathrm{N}}$$

$$C_{\mathrm{CX},\mathrm{N}}$$

$$f_{\mathrm{H}}^{\mathrm{N}}$$

$$f_q$$

$$f_{\mathrm{H}}^{\mathrm{B}}$$

$$\frac{\partial n_{\mathrm{H},1s}^{\mathrm{B}}}{\partial z}=\frac{n_{\mathrm{H},1s}^{\mathrm{N}}{C}_{\mathrm{CX},\mathrm{N}}-n_{\mathrm{H},1s}^{\mathrm{B}}{C}_{\mathrm{I},\mathrm{B}}}{u_2},$$

$$C_{\mathrm{I},\mathrm{B}}$$

$$\frac{\partial n_{\mathrm{p}}}{\partial z}=\frac{n_{\mathrm{H},1s}^{\mathrm{N}}{C}_{\mathrm{I},\mathrm{N}}+n_{\mathrm{H},1s}^{\mathrm{B}}{C}_{\mathrm{I},\mathrm{B}}}{u_2}.$$

$$\frac{\mathrm{d}{n}_{\mathrm{H},j}}{\mathrm{d}t}=\sum _k\left\{n_{\mathrm{H},k}\left(C_{k,j}+P_{k,j}\right)-n_{\mathrm{H},j}\left(C_{j,k}+P_{j,k}\right)\right\},$$

$$C_{j,k}=\sum _{q=\{\mathrm{e},\mathrm{p}\}}{n}_q\int f_{\mathrm{H}}({\mathbf{v}}_{\mathrm{H}})f_q({\mathbf{v}}_{\mathbf{q}})\Delta v_q{\sigma }_q^{j,k}{\mathrm{d}}^3{\mathbf{v}}_{\mathrm{H}}{\mathrm{d}}^3{\mathbf{v}}_{\mathbf{q}},$$

$$P_{k,j}=A_{k,j}-\int \frac{4\pi }{h\nu }\mathrm{d}\nu {\int }_{-1}^1\frac{{\sigma }_{\nu ,\mu }^{k,j}{I}_{\nu ,\mu }}{2}\mathrm{d}\mu ,$$

$$P_{j,k}=\int \frac{4\pi }{h\nu }\mathrm{d}\nu {\int }_{-1}^1\frac{{\sigma }_{\nu ,\mu }^{j,k}{I}_{\nu ,\mu }}{2}\mathrm{d}\mu ,$$

$$\frac{\mathrm{d}{I}_{\nu ,\mu }}{\mathrm{d}s}=-k_{\nu ,\mu }{I}_{\nu ,\mu }+j_{\nu ,\mu },$$

$$I_{\nu ,\mu }({\tau }_{\nu ,\mu })={\int }_0^{{\tau }_{\nu ,\mu }}{S}_{\nu ,\mu }({\tilde{\tau }}_{\nu ,\mu }){\mathrm{e}}^{-\frac{{\tilde{\tau }}_{\nu ,\mu }-{\tau }_{\nu ,\mu }}{\mu }}\frac{\mathrm{d}{\tilde{\tau }}_{\nu ,\mu }}{\mu },$$

$${\tau }_{\nu ,\mu }={\int }_{z_{\mathrm{out}}}^z{k}_{\nu ,\mu }\frac{\mathrm{d}z}{\mu }(\mathrm{for}\mu >0),$$

$${\tau }_{\nu ,\mu }={\int }_{z_{\mathrm{in}}}^z{k}_{\nu ,\mu }\frac{\mathrm{d}z}{\mu }(\mathrm{for}\mu <0).$$

$$n_{\mathrm{H},1s}{C}_{1s,j}+\sum _k\left\{n_{\mathrm{H},k}{P}_{k,j}-n_{\mathrm{H},j}{P}_{j,k}\right\}=0.$$

$${\sigma }_{\nu ,\mu }^{j,k}=\pi r_{\mathrm{e}}c{f}_{j,k}{\varphi }_{\nu ,\mu }^{j,k},$$

$${\varphi }_{\nu ,\mu }^{j,k}$$

$$\Delta {\nu }_{\mathrm{D}}^{\mathrm{N}}={\nu }_{j,k}^{\prime }\sqrt{\frac{2k_{\mathrm{B}}{T}_0}{m_{\mathrm{H}}{c}^2}},$$

$$\Delta {\nu }_{\mathrm{D}}^{\mathrm{B}}={\nu }_{j,k}^{\prime }\sqrt{\frac{2k_{\mathrm{B}}{T}_{\mathrm{p}}}{m_{\mathrm{H}}{c}^2}}.$$

Full text

††pubyear: 2018††pagerange: Radiative transfer of hydrogen lines from supernova remnant shock waves: contributions of 2s-state hydrogen atoms–Radiative transfer of hydrogen lines from supernova remnant shock waves: contributions of 2s-state hydrogen atoms

1 Introduction

Balmer line emissions from supernova remnant (SNR) shock waves are relied upon as a probe of the physics of collisionless shocks. The shock transition occurs on a length scale much shorter than that associated with a particle mean free path to Coulomb scattering, so thermal equilibrium is much less strongly enforced. Balmer lines can be used to diagnose the resulting departures from equilibrium, such as the effects on ion-electron temperature ratio, the nature of the shock precursor and the acceleration of non-thermal particles (see Raymond, 1991; Heng, 2010, and Section 2 for reviews). Such a shock is often called a Balmer dominated shock (BDS).

The physics of particle acceleration may be the most important issue because it is crucial to the origin of cosmic-rays. Moreover, accelerated particles in collisionless shocks are often considered to be responsible for the radiation from high-energy astrophysical sources from the radio to the TeV band. An important concern is to specify the density of accelerated particles in SNR shocks, a necessary step towards confirming SNR shocks as the main sites producing Galactic cosmic-rays. In addition, it allows us to quantify the back reaction of accelerated particles on the background shock structure. Note that in the standard model of particle acceleration (i.e. diffusive shock acceleration), the distribution function of accelerated particles depends on the shock structure (e.g. Berezhko & Ellison, 1999, and references therein). If the amount of accelerated particles is significant, in other words, if the kinetic energy of the shock consumed due to the particle acceleration is a large fraction of the shock energy, the downstream temperature becomes considerably lower than the case of an adiabatic shock, that is, there is some missing thermal energy. This energy loss from the shock has been widely investigated (e.g. Hughes et al., 2000; Tatischeff & Hernanz, 2007; Helder et al., 2009; Morlino et al., 2013a, b, 2014; Shimoda et al., 2015; Hovey et al., 2018). In the latest development, Shimoda et al. (2018) showed that the linear polarization degree of H $\alpha$ observed perpendicularly to the shock velocity vector depends on the energy-loss from the shock. Such polarized H $\alpha$ was originally predicted by Laming (1990) to estimate the ion-electron temperature equilibrium for the adiabatic shock and was recently discovered in SN 1006 by Sparks et al. (2015).

As well as being directly excited from the ground state, Balmer lines (intensity, line profile, polarization and so on) are affected by the conversion of Lyman lines to Balmer lines. For example, the absorption of Ly $\beta$ by a hydrogen atom results in radiative excitation from 1s to 3p, and the excited atom can emit H $\alpha$ by the spontaneous transition from 3p to 2s. Simultaneously, the conversion yields the 2s-state hydrogen atom, which creates the two-photon continuum by the spontaneous transition from 2s to 1s. Thus, the Ly $\beta$ to H $\alpha$ conversion impacts the total intensity, line profile and net polarization of H $\alpha$. Moreover, an adequate density of 2s-sate atoms can further scatter H $\alpha$ photons. Although such fundamental physics is well known, it has not been well studied in SNR shocks. In fact, it is usually assumed that the Ly $\beta$ photons are at the limits of either completely optically thick or optically thin at SNR shocks, that is, they are completely converted to H $\alpha$ photons or not at all (e.g. Heng & McCray, 2007; van Adelsberg et al., 2008; Morlino et al., 2012, 2013b; Shimoda et al., 2018). Contrary to this, Ghavamian et al. (2001) studied the conversion of Ly $\beta$ and Ly $\gamma$ to H $\alpha$ and H $\beta$ by Monte Carlo simulations and claimed that intermediate conversion occurs. However, they and previous studies did not consider the population of 2s-state hydrogen atoms. In this paper, we provide a formulation of the radiative line transfer with the rate equation of atomic population and study the nature of Balmer line emissions from SNR shocks. In this paper we do not consider the polarization, deferring that instead to a later work. Note that as a first step, our model makes several simplifications in the treatment of the SNR shock, handling the hydrogen atoms as fluids and supposing no particles leaking back upstream (e.g. cosmic-rays). Our calculation of radiative transfer is based on so-called the ray-tracing method and uses updated atomic data from the literature (e.g. Heng & Sunyaev, 2008; Tseliakhovich et al., 2012). Moreover, we consider only hydrogen line emissions and ignore bremsstrahlung radiation, thermal emissions from the SNR ejecta, and external radiation sources. Thus, our model possibly predicts somewhat smaller population of 2s-state hydrogen atoms than would be the case in a realistic SNR shock.

This paper is organized as follows. In Section 2, we briefly review the BDSs as a probe of the collisionless shock physics and give a simple estimate of the occupation number of the 2s-state hydrogen atoms. In Section 3, we formulate the radiative line transfer problem for the SNR shock. In Section 4, we present the results for atomic populations. In Section 5, we consider how the hydrogen lines are observed based on the results of the atomic population computations. Finally, we summarize our results.

2 Diagnostics of Balmer Dominated Shocks

In this section, we will briefly review the diagnostics of BDSs and estimate the previously neglected population of 2s-state hydrogen atoms.

The basic theory of the Balmer line emission from SNR shocks was described by Chevalier et al. (1980). They pointed out that hydrogen atoms in the interstellar medium (ISM) do not suffer shock-heating across the shock front because the SNR shock is ‘collisionless’. The shock transition is formed by the interaction between charged particles and plasma waves rather than by Coulomb collisions. Then, the hydrogen atoms entering the downstream region collide with the shock-heated, charged particles. The collisions result in several atomic reactions such as ionization, direct excitation and charge-transfer between the hydrogen atoms and the shock-heated protons. Atoms which have not experienced any charge-transfer reactions emit a ‘narrow’ line with width of $\sim 10~{}{\rm km~{}s^{-1}}$, while the atoms that have undergone a charge-transfer reaction emit a ‘broad’ line with a width characteristic of the post-shock ion temperature; $\sim 1000~{}{\rm km~{}s^{-1}}$. Thus, the profile of hydrogen lines observed in the SNR shocks consists of at least these two components. The width of each component corresponds to the upstream temperature (narrow) and the downstream proton temperature (broad), respectively. Thus, the downstream temperature can be derived from the width of broad component. The intensity ratio of the broad component to the narrow component is often relied on to estimate the ion-electron temperature ratio (Ghavamian et al., 2013). The Balmer decrement, which is the intensity ratio of the narrow H $\alpha$ to the narrow H $\beta$, depends on the ionization degree of the ambient gas around the SNR shock. The H $\alpha$ (H $\beta$) intensity is enhanced via the conversion from Ly $\beta$ (Ly $\gamma$) to H $\alpha$ (H $\beta$). Since the absorption cross-section of Ly $\beta$ is larger than that of Ly $\gamma$, the Ly $\beta$ to H $\alpha$ conversion may occur more than the Ly $\gamma$ to H $\beta$ conversion. Thus, the Balmer decrement reflects the optical properties of BDSs, which depends on the ionization structure of hydrogen. This basic model did not consider the existence of particles leaking back upstream (e.g. cosmic-rays and fast-neutral particles produced by the charge-transfer reaction) and was further developed in several papers (e.g. Ghavamian et al., 2001; Heng & McCray, 2007; van Adelsberg et al., 2008).

The existence of particles leaking back to the upstream region is implied by observations of the full width at half maximum (FWHM) of the narrow component of H $\alpha$ (30–50 km s*-1*) and H $\alpha$ emission from the upstream region; i.e. the shock-precursor emission (e.g. Smith et al., 1994; Ghavamian et al., 2000; Sollerman et al., 2003; Lee et al., 2007, 2010; Medina et al., 2014; Katsuda et al., 2016; Knežević et al., 2017). Note that the FWHM of 30–50 km s*-1* implies an upstream temperature of 2.5–5.6 eV, which would be too high for neutral hydrogen atoms to exist if it was the equilibrium temperature of the ISM. Therefore, non-thermal pre-heating by the leaking particles in the upstream region is expected. Moreover, Raymond et al. (2010) showed that the H $\alpha$ line profile observed in Tycho’s SNR can be fitted by three Gaussian functions. This implies that there is an intermediate component with temperature between the narrow and broad components (i.e. a non-thermal velocity distribution of hydrogen atoms undergoing the charge-transfer reaction) or the velocity distribution of downstream protons undergoing charge-transfer reaction deviates from the Maxwellian distribution due to a non-thermal wing (Raymond et al., 2010, 2017). Semi-analytical models (Morlino et al., 2012, 2013b) and hybrid simulations (Ohira, 2016a) suggest that the cosmic-rays and/or the fast-neutral particles emerging from the charge-transfer reaction lead to such shock-precursor emission accompanied by the intermediate component and the anomalous width of narrow component. Note that the shock-precursor emission due to the absorption of Lyman photons may be inevitable. In this paper, we refer to such precursor emission as ‘shock-precursor-like emission’ or ‘photo-precursor emission’ to distinguish it from the cosmic-ray/fast-neutral precursor.

Here we give a simple estimate for the population of 2s-state hydrogen atoms by considering a 3-level system (1s, 2s and 3p). We set the rate equation for the 2s population as

[TABLE]

where $n_{{\rm H},j}$ is the number density of hydrogen atoms in the state $j$, $C_{j,k}$ and $A_{j,k}$ are the collisional excitation rate and spontaneous decay rate for the transition from $j$ to $k$, respectively. For the bound states, we use the notation $j=n_{j}l_{j}$, where $n_{j}$ is the principal quantum number of the state $j$. Similarly, $l_{j}=0,1,2,3,...,n_{j}-1$ (equivalently: s, p, d, f,…) is the orbital angular-momentum quantum number of the state $j$. Here we suppose that the depopulation term of 2s-state atoms is dominated by the spontaneous transition at the rate of $A_{\rm 2s,1s}\simeq 8.2~{}{\rm s^{-1}}$. In reality, the collisional transition from 2s to 2p can be a subdominant process for depopulation. For a collision at a velocity $\sim 10^{8}~{}{\rm cm~{}s^{-1}}$, which is a typical velocity scale for young SNR shocks, the cross-section is $\sim 10^{-13}~{}{\rm cm^{2}}$, giving a reaction rate $\sim 10^{-5}~{}{\rm cm^{3}~{}s^{-1}}$ (e.g. Janev et al., 1987; Sahal-Brechot et al., 1996). Thus, if the density is $\sim 10^{6}~{}{\rm cm^{-3}}$, the collisional depopulation becomes important. Note that we assume no strong radiation field inducing the radiative transition from 2s to any other state.111 The cross-section of 2s-2p collisional excitation has the maximum value $\sim 10^{-10}~{}{\rm cm^{2}}$ around the relative velocity $\sim 10^{6}\mathchar 45\relax 10^{7}~{}{\rm cm~{}s^{-1}}$ for both proton and electron impacts (Janev et al., 1987), leading to a rate of $10^{-4}-10^{-3}$ cm3s*-1* and critical densities of $10^{4}-10^{5}$ cm*-3* at slower shocks. The occupation number of 3p, $n_{\rm H,3p}$, depends on the absorption of Ly $\beta$. Here we assume an isotropic radiation field for Ly $\beta$. Then, we obtain the rate equation for 3p as

[TABLE]

where $h$, $\nu$, $\sigma_{\nu}^{\rm 1s,3p}$ and $I_{\nu}$ are the Planck constant, frequency, absorption cross-section for the transition from 1s to 3p and the specific intensity, respectively. The intensity is set to be

[TABLE]

where $S_{\nu}$ and $\tau_{\nu}$ are the source function and optical depth, respectively. $\sigma^{\prime}$ is a combination of physical constants relevant to the radiative absorption cross-section. Thus, we derive the occupation number of 2s as

[TABLE]

where $\tau_{0}$ is the optical depth at the line centre. Here we assume a narrow line profile function $\phi_{\nu}$ for which we can approximate as $\int_{0}^{\infty}(1-{\rm e}^{-\tau_{\nu}}){\rm d}\nu\approx 1-{\rm e}^{-\tau_{0}}$. The terms in the brackets […] indicate the contribution of the combination of the absorption and cascades. Note that roughly say, the ratios are $A_{\rm 3p,2s}/(A_{\rm 3p,1s}+A_{\rm 3p,2s})\simeq 0.118$, $C_{\rm 1s,3p}/C_{\rm 1s,2s}\sim 2\mathchar 45\relax 10$ and $C_{\rm 1s,2s}/A_{\rm 2s,1s}\sim 10^{-9}n_{\rm p}$, where $n_{\rm p}$ is the proton number density. Thus, if Ly $\beta$ is in the optically thick limit, $n_{\rm H,2s}$ is enhanced roughly at most ten times compared with the optically thin case. The absorption coefficient of H $\alpha$ at the line centre becomes

[TABLE]

where $\sigma_{0}({\rm H\alpha})$ is the radiative cross-section of H $\alpha$ at the line centre for given temperature $T_{0}$. Thus, if the SNR shock interacts with somewhat dense clump with a density of $\sim 30~{}{\rm cm^{-3}}$ and a size of $\sim 1~{}{\rm pc}$, the H $\alpha$ emission can be scattered. Note that the H $\beta$ emission can also be scattered but its absorption coefficient is about quarter of the H $\alpha$ coefficient. The interaction between the shock and a dense clump is implied by the ripple of an SNR shock with a length-scale of $\sim 10$ per cent of SNR radius (e.g. Ishihara et al., 2010; Williams et al., 2013, 2016; Miceli et al., 2014; Sano et al., 2017; Tsubone et al., 2017, and see the discussion of Shimoda et al. 2015). Note that according to magnetohydrodynamic simulations performed by Inoue et al. (2009); Inoue et al. (2012), even if the shock propagates into a simulated ISM having density contrast ranging in $\sim 1\mathchar 45\relax 30~{}{\rm cm^{-3}}$ as a consequence of thermal instability, the scale length of rippling is $\sim 10$ per cent of the length of sides of simulation box.

The temperature $6000$ K we assumed is often taken for the warm neutral medium of ISM (e.g. Ferrière, 2001). If we suppose the temperatures implied by the measured H $\alpha$ widths of Sollerman et al. (2003), we obtain the absorption coefficients around half those assumed at $6000$ K. Note that for SNR Cygnus Loop, Medina et al. (2014) pointed out that the pre-shock gas is photoionized and heated up to $\sim 17000$ K by the emissions from post-shock region.

In an actual SNR shock, the optical depth of Ly $\beta$ may be intermediate (Ghavamian et al., 2001). Moreover, the validity of this simple estimate is still unclear because of the many complexities of BDSs, especially the ionization structure of hydrogen. Therefore, we consider more sophisticated formulas for the line transfer problem and solve them numerically.

3 Formulation of line transfer

The line transfer problem is reviewed in several papers (e.g. Castor, 2004). We apply their formulation to the problem for SNR shocks propagating into pure atomic hydrogen plasma, which consists of hydrogen atoms (denoted ’H’), protons (’p’) and electrons (’e’). The shock is set to be stationary, axially symmetric about the $z$-axis, plane-parallel to $x\mathchar 45\relax y$ plane and located at $z=0$ (see, Fig. 1a). We set two free-escape boundaries for photons upstream ($z=z_{\rm out}$) and downstream ($z=z_{\rm in}$) of the shock. For simplicity, we assume that there are no particles leaking to the upstream region and that the radiation field consists of only the H line emissions (i.e. bremsstrahlung radiation, emission from the SNR ejecta and any other external radiation sources are neglected). Moreover, we assume temperature equilibrium for the upstream plasma and fix the upstream temperature at $T_{\rm 0}=6000$ K for simplicity.

Firstly, we describe the ionization structure of hydrogen. Let $n_{{\rm H},j}^{\rm N}$ be the number density of ’narrow’ (i.e. cold) hydrogen atoms, which have not experienced charge-exchange reactions, while $n_{{\rm H},j}^{\rm B}$ is the number density of ’broad’ (i.e. hot) hydrogen atoms emerging from charge-exchange reactions. Obviously, we have the relation $n_{{\rm H},j}=n_{{\rm H},j}^{\rm N}+n_{{\rm H},j}^{\rm B}$. Fig. 1b is a schematic illustration of the spatial distribution of particles. We consider that the partially ionized plasma flows from the far upstream region ($z<z_{\rm out}$) and presume that it is in ionization equilibrium. Hence, we set the boundary conditions as $n_{\rm H,1s}^{\rm N}(z<0)=n_{\rm H,1s}^{\rm N}(z_{\rm out})$, $n_{\rm H,1s}^{\rm B}(z<0)=0$ and $n_{\rm p}(z<0)=n_{\rm p}(z_{\rm out})$, where $n_{\rm p}$ is the number density of protons. At the shock ($z=0$), we assume the strong shock jump conditions,

[TABLE]

where $V_{\rm sh}$ is the shock velocity, $k_{\rm B}$ is the Boltzmann constant and $m_{\rm p}$ is the proton mass. $T_{\rm p}$ and $T_{\rm e}$ are the downstream temperatures of protons and electrons, respectively. The effective mean molecular weight is $\mu^{\prime}$, which is defined as

[TABLE]

where $\mu_{\odot}^{\prime}=0.62$ and $m_{\rm e}$ is the electron mass (see Shimoda et al., 2018, for details). Note that the number density of downstream protons is function of $z$, while $u_{2}$, $T_{\rm p}$ and $T_{\rm e}$ are kept constant in the model.222SNR shocks propagating into a dense medium with solar metallicity can be radiative. Hollenbach & McKee (1979) give the cooling length of shock heated gas (as post shock column density) by using the cooling function of Raymond et al. (1976), $N_{\rm c}\simeq 2\times 10^{17}\beta^{0.6}V_{\rm sh,7}{}^{4.2}~{}{\rm cm^{-2}}$, where $V_{{\rm sh},m}\equiv 10^{m}~{}{\rm cm~{}s^{-1}}$. Here we regard the cooling function as determined by the electron temperature. In this paper, we consider $V_{\rm sh}\ga 10^{8}~{}{\rm cm~{}s^{-1}}$ for which $N_{\rm c}\ga 10^{21.2}~{}{\rm cm^{-2}}$. Thus, the adiabatic shock approximation can be valid for the atomic transition layer from $z=0$ to $z=z_{\rm in}$ at which the column density $\sim 4\times 10^{16}~{}{\rm cm^{-2}}(n/1~{}{\rm cm^{-3}})$ is lower than $N_{\rm c}$. In the following, we neglect the radiative recombination rate $\sim 10^{-13}~{}{\rm cm^{3}~{}s^{-1}}$, which is much smaller than any other rates. Moreover, we assume $n_{{\rm H},j\neq{\rm 1s}}\ll n_{\rm H,1s}$ (see Eq. (4)). We will address these assumptions later. Then, the spatial distribution of narrow hydrogen atoms is given by

[TABLE]

where we define the collisional ionization rate,

[TABLE]

and the charge-exchange rate,

[TABLE]

Here $n_{q}$ is the number density of particle $q$, the symbols $\bm{v_{q}}$ and $\bm{v}_{\rm H}$ denote the velocity vectors of particle $q$ and the hydrogen atom, respectively, \Delta v_{q}\equiv\big{|}\bm{v_{\rm H}}-\bm{v_{q}}\big{|} is the relative velocity between the hydrogen atom and particle $q$, $\sigma_{q}^{\rm I}$ is the ionization cross-section by collision with particle $q$ and $\sigma_{\rm CX}$ is the total cross-section of charge-exchange reactions. We omit the rates of collisions between hydrogen atoms, which are small compared with the rates by proton/electron collisions because of the lack of shock compression or heating. Note that the collisional rates of He2+ impacts, which are also omitted in this paper, would make a moderate contribution (Laming et al., 1996). We assume the velocity distribution function of narrow hydrogen atoms to be

[TABLE]

where $m_{\rm H}$ is the hydrogen atom mass and $\bm{V_{\rm sh}}=(0,0,V_{\rm sh})$. Similarly, the distribution functions of the downstream protons/electrons are given by

[TABLE]

where $\bm{u_{2}}=(0,0,u_{2})$ and $q=\{\rm p,e\}$. For simplicity, we assume that the broad atoms have the same mean velocity and temperature as the downstream protons. Note that the velocity distribution of broad atoms can substantially deviate from the proton’s due to the velocity dependence of the cross-section of charge-exchange reaction (e.g. Heng & McCray, 2007), but the nature of radiative line transfer depends mainly on the ionization structure of ‘narrow’ atoms, which can be well approximated by Eq. (11). Then, for the broad hydrogen atoms, we obtain the distribution function,

[TABLE]

and the differential equation of their spatial distribution,

[TABLE]

where

[TABLE]

Accordingly, the differential equation for the spatial distribution of downstream protons is

[TABLE]

Note that the electron number density is given by the charge neutrality condition $n_{\rm e}=n_{\rm p}$. We consider the radiative line transfer and the population of bound-state ($j\neq{\rm 1s}$) atoms under the ionization structure given by the above formulae.

Here we consider the population of excited hydrogen atoms. The rate equation for the excited-state $j$ is

[TABLE]

where $P_{k,j}$ and $C_{k,j}$ ($P_{j,k}$ and $C_{j,k}$) are radiative and collisional rates per unit time for the transition from $k$ to $j$ ($j$ to $k$), respectively. The collisional rate is

[TABLE]

where $\sigma_{q}^{k,j}$ is a sum of all kinds of collisional cross-sections between particle $q$ and hydrogen atom resulting in the transition from $k$ to $j$. Here we do not distinguish between the broad and narrow hydrogen atoms so that $f_{{\rm H},j}=f_{\rm H}^{\rm N}+f_{\rm H}^{\rm B}$. For $n_{k}>n_{j}$ (henceforth, we refer $k>j$), the radiative rates are

[TABLE]

and

[TABLE]

where $A_{k,j}$ is the rate of the spontaneous transition from $k$ to $j$ and $\mu\equiv\cos\theta$ indicates the direction of ray making an angle $\theta$ to the $z$-axis (see Fig. 1a). Here $\sigma_{\nu,\mu}^{k,j}$ is the radiative cross-section for the ray in the direction $\mu$ at the frequency $\nu$ resulting in the transition from $k$ to $j$, and $I_{\nu,\mu}$ is the specific intensity of the ray directed in $\mu$ (in unit erg cm*-2* s*-1* Hz*-1str-1*). In order to evaluate the radiative rate $P_{j,k}$, we need to simultaneously solve the radiation transfer equation

[TABLE]

where $s$ is the unit length measured along the path of ray, and $k_{\nu,\mu}$ and $j_{\nu,\mu}$ are the absorption coefficient and emission coefficient (or emissivity), respectively. Note that these coefficients depend on the occupation number of the atomic states. We give their formulae later. According to the two free-escape boundaries at $z=z_{\rm out}$ and $z=z_{\rm in}$, we derive the formal solution of the intensity as

[TABLE]

where $S_{\nu,\mu}\equiv j_{\nu,\mu}/k_{\nu,\mu}$ is the source function. The optical depth $\tau_{\nu,\mu}$ is defined as

[TABLE]

and

[TABLE]

We can solve the set of equations iteratively.333 We do not solve the momentum gain of atoms due to the absorption of photons for simplicity. For example, the velocity change of atom for each scattering of Ly ${\rm\alpha}$ is $\sim 3\times 10^{-3}~{}{\rm km~{}s^{-1}}$, which is much smaller than the thermal velocity of the hydrogen atoms $\ga 10~{}{\rm km~{}s^{-1}}$.

Several of the terms of Eq. (20) can be omitted because of the enormous difference in time scales. The orders of magnitude of the time scales are $t_{\rm r}\sim 10^{-8}\mathchar 45\relax 10^{-1}~{}{\rm s}$ for the decay time of excited atoms due to the spontaneous transition, $t_{\rm c}\sim 10^{4}\mathchar 45\relax 10^{8}~{}{\rm s}~{}(n/1~{}{\rm cm^{-3}})^{-1}$ for the collisions between the hydrogen atoms and shock heated protons/electrons, and the mean collision time for the upstream medium with temperature of $\sim 10^{4}$ K is $\sim 10^{11}$ s. Note that the radiative rates due to photon absorption P_{j,k}\big{|}_{j<k} may be at most comparable with the maximum value of collisional excitation rates $\sim 1/t_{\rm c}$ (in the sense of the orders of magnitude estimate) because the radiation field consists only of the line emissions. The ratio of the excitation rate $\sim 1/t_{\rm c}$ to the decay rate $\sim 1/t_{\rm r}$ indicates a very small occupation number of excited state atoms $n_{{\rm H},j\neq 1s}\sim 10^{-9}\mathchar 45\relax 10^{-16}n_{\rm H,1s}$ (see Eq. (4)). In addition, the recombination rate of hydrogen atoms, $\sim 10^{-14}~{}{\rm cm^{3}~{}s^{-1}}$, is negligibly small. Assessing the above factors, we can ignore the collisional deexcitation, recombination and any other collisional processes in the upstream region. Thus, for excited states ($j\neq$ 1s), we obtain the statistical equilibrium condition

[TABLE]

Here we approximate ${\rm d}n_{{\rm H},j\neq{\rm 1s}}/{\rm d}t\approx 0$ and consider direct excitation and charge-exchange as the collisional excitation processes.

We can derive the ratio of $n_{{\rm H},j}$ to $n_{\rm H,1s}$ from Eq. (28) if the radiative rates $P_{j,k}$ are given. To do this, the absorption coefficient $k_{\nu}$ and the emissivity $j_{\nu}$ are required. In this paper, we consider only the resonant scattering, whose cross-section is typically $\sim 10^{12}$ times the Thomson cross-section. Moreover, we neglect the coherence of all atomic (quantum) processes for the line profile function, because its contribution is usually very small and appears at a frequency far from the line centre (e.g. the Lorentzian wing). Furthermore, we ignore any overlaps in the frequency among each line. This is true for the lines we are interested in (e.g. Ly $\alpha$, Ly $\beta$, Ly $\gamma$, H $\alpha$, H $\beta$, Pa $\alpha$ and so on). Thus, we individually treat the specific intensities, absorption coefficients and emissivities of each line induced by the transition from $n_{k}$ to $n_{j}$. The radiative cross-section for the ray in the direction $\mu$ resulting in the transition from $j$ to $k$ ($j<k$) is given by

[TABLE]

where $r_{\rm e}$ is the classical electron radius, $c$ is the speed of light, $f_{j,k}$ is the oscillator strength, and $\phi_{\nu,\mu}^{j,k}$ is the line profile function for the transition from $j$ to $k$:

[TABLE]

where the normalization condition $\int\phi_{\nu,\mu}^{j,k}{\rm d}\nu{\rm d}\mu=1$ is satisfied. $\Delta\nu_{\rm D}^{\rm N}$ is the Doppler frequency for the narrow hydrogen atoms,

[TABLE]

where $\nu_{j,k}^{\prime}$ is the frequency at the line centre measured in the atom rest frame, and $\Delta\nu_{\rm D}^{\rm B}$ is the Doppler frequency for the broad hydrogen atoms,

[TABLE]

The centroid frequencies are, respectively,

[TABLE]

and

[TABLE]

Then, the absorption coefficient for the transition from $j$ to $k$ is

[TABLE]

where $g_{l_{j}}$ is the statistical weight of the state $j$ and we use the relation $f_{k,j}=-\frac{g_{l_{j}}}{g_{l_{k}}}f_{j,k}$. Note that $\sigma_{\nu}^{k,j}n_{{\rm H},k}$ indicates the stimulated emission. In the following, we consider only dipole transitions ($|l_{j}-l_{k}|=1$) for the absorption coefficient. The net absorption coefficient for each line resulting in the transition from $n_{j}$ to $n_{k}$ is given by

[TABLE]

where for fixed $n_{j}$ and $n_{k}$, we take the summation of $k_{\nu,\mu}^{j,k}$ for $l_{j}$ and $l_{k}$ under the constraint $|l_{j}-l_{k}|=1$. For the emissivity, we consider not only the dipole transition but also the transition from 2s to 1s (i.e. 2$\gamma$–decay). This 2s–1s transition yields two photons to satisfy the conservation of net angular momentum before and after the transition. The frequencies of two photons range $0<\nu<\nu_{\rm Ly\alpha}$ and their sum is equal to the frequency of Ly $\alpha$, $\nu_{\rm Ly\alpha}$, because of energy conservation. The spontaneous 2s-1s transition rate depends on the frequency. The net transition rate is given by

[TABLE]

where $\varphi_{\nu}$ is equivalently the line profile function. This has a peak at a frequency of $\nu=0.5\nu_{\rm Ly\alpha}$ and has the half-maximum values at $\nu\simeq 0.1\nu_{\rm Ly\alpha}$ and $\nu\simeq 0.9\nu_{\rm Ly\alpha}$ (e.g. Chluba & Sunyaev, 2008). Thus, the 2$\gamma$-decay has a much wider profile compared with $\phi_{\nu,\mu}^{j,k}$ and therefore behaves as continuum emission with respect to each line with frequency lower than $\nu_{\rm Ly\alpha}$ (i.e. $\varphi_{\nu}=0$ for $\nu>\nu_{\rm Ly\alpha}$). Hence, the emissivity for the ray directed along $\mu$ at the frequency $\nu$ for each line is

[TABLE]

where we take the Doppler shift into account for the line profile function $\phi_{\nu,\mu}^{k,j}$ but neglect it for the 2$\gamma$-decay profile $\varphi_{\nu}$. Note that the line profile function for the emissivity, $\phi_{\nu,\mu}^{k,j}$, has the same form as Eq. (30) because the direction $\mu$ is fixed and any coherences are ignored. Using the absorption coefficient $k_{\nu,\mu}$ (Eq. (36)) and the emissivity $j_{\nu,\mu}$ (Eq. (38)), we derive the specific intensity of each line from its formal solution Eq. (25). Then, the radiative rates, Eqs. (22) and (23), are calculated. Hence, we obtain the ratio $n_{{\rm H},j}$ to $n_{\rm H,1s}$ from Eq. (28). Using the ratio of $n_{{\rm H},j}$ to $n_{\rm H,1s}$, we obtain newly $k_{\nu,\mu}$ and $j_{\nu,\mu}$, solving the radiation field, and eventually obtain the ratio again. We iterate the above procedure until the ratio converges. Note that $n_{\rm H,1s}=n_{\rm H,1s}^{\rm N}+n_{\rm H,1s}^{\rm B}$ and $n_{\rm p}$ are solved separately from any excited states $n_{{\rm H},j\neq{\rm 1s}}$ by Eqs. (11),  (17) and (19). 444 Here we neglect the Lorentz transformations for the intensity $I_{\nu,\mu}=(\nu/\nu^{\prime})^{3}I^{\prime}_{\nu^{\prime},\mu^{\prime}}$, the emissivity $j_{\nu,\mu}=(\nu/\nu^{\prime})^{2}j^{\prime}_{\nu^{\prime},\mu^{\prime}}$, and the absorption coefficient $k_{\nu,\mu}=(\nu^{\prime}/\nu)k^{\prime}_{\nu^{\prime},\mu^{\prime}}$, where prime indicates the rest frame of hydrogen atom, i.e. $I_{\nu,\mu}\approx I^{\prime}_{\nu^{\prime},\mu^{\prime}}$, $j_{\nu,\mu}\approx j^{\prime}_{\nu^{\prime},\mu^{\prime}}$ and $k_{\nu,\mu}\approx k^{\prime}_{\nu^{\prime},\mu^{\prime}}$.

For the numerical calculation of $n_{{\rm H},j}$, we use data for collisional cross-sections and their fitting functions provided in the literature (Barnett et al., 1990; Janev & Smith, 1993; Bray & Stelbovics, 1995; Heng & Sunyaev, 2008; Tseliakhovich et al., 2012). Note that since the available data are limited, we calculate states up to 4f. Therefore, we solve Ly $\alpha$, Ly $\beta$, Ly $\gamma$, H $\alpha$, H $\beta$ and Pa $\alpha$. Moreover, the data of direct-collisional excitation to $n_{j}=4$ level by proton impact are unavailable for the range $v_{\rm H}\leq 1000~{}{\rm km~{}s^{-1}}$ (Tseliakhovich et al., 2012). We treat them to be zero in this range (see also Shimoda et al., 2018). We also refer to the data table for the spontaneous transition rates, the oscillator strengths and the centroid frequencies $\nu^{\prime}_{j,k}$ provided by Wiese & Fuhr (2009). Table 1 presents the $A_{k,j}$ and $f_{j,k}$ we are interested in. The centroid frequencies measured in the atom rest frame, $\nu^{\prime}_{j,k}$, are calculated from the centroid wavelengths ($\lambda_{\rm vac}$ in Wiese & Fuhr, 2009) divided by the speed of light. Here we assume $\nu^{\prime}_{\rm 2s,1s}=\nu^{\prime}_{\rm 2p,1s}$. Note that the offset of centroid frequencies, for example, $1-\nu^{\prime}_{\rm 2s,3p}/\nu^{\prime}_{\rm 2p,3s}\simeq 1.8\times 10^{-5}$ is comparable with the Doppler shift of hydrogen atoms with velocity of $\simeq 5~{}{\rm km~{}s^{-1}}$. Thus, the widths of Balmer and Paschen lines will be broader than the width given by only the thermal Doppler shift. Moreover, we use the fitting function of $\varphi_{\nu}$ provided by Chluba & Sunyaev (2008). Note that the statistical weights for each state are $g_{\rm s}=2$, $g_{\rm p}=6$, $g_{\rm d}=10$ and $g_{\rm f}=14$, respectively.

4 Atomic population of hydrogen atoms in supernova remnant shocks

In this section, the results of the atomic population calculations are exhibited. We parameterize the shock by the total number density in the upstream region $n_{\rm tot,0}\equiv n_{\rm H,out}+n_{\rm p,out}$, the upstream ionization degree (or the upstream proton fraction) $\chi_{0}\equiv n_{\rm p,out}/n_{\rm tot,0}$, the downstream proton temperature $T_{\rm p}$ and the electron temperature $T_{\rm e}$ (i.e. $\beta=T_{\rm e}/T_{\rm p}$). The two escape-boundaries are set at $z_{\rm out}=-5\times 10^{16}~{}{\rm cm}\left(\frac{n_{\rm tot,0}}{1~{}{\rm cm^{-3}}}\right)^{-1}$ and $z_{\rm in}=4u_{2}/(n_{\rm tot,0}C_{\rm I,B})$.

Firstly, we give the ionization structure of hydrogen. By rewriting Eq. (11) as

[TABLE]

we can easily recognize that each density fraction $n_{\rm H}^{\rm N}/n_{\rm tot,0}$, $n_{\rm H}^{\rm B}/n_{\rm tot,0}$ and $n_{\rm p}/n_{\rm tot,0}$ is a self-similar function with respect to $n_{\rm tot,0}z$ (having the same value at the same $n_{\rm tot,0}z$) for fixed $T_{\rm p}$, $\beta$ and $\chi_{0}$. Fig. 2 shows $n_{\rm H,1s}^{\rm N}/n_{\rm tot,0}$, $n_{\rm H,1s}^{\rm B}/n_{\rm tot,0}$, and $n_{\rm p}/n_{\rm tot,0}$ for given $k_{\rm B}T_{\rm p}=12~{}{\rm keV}$ and $\beta=0.1$ with fixed values $\chi_{0}=$(0.1, 0.3 and 0.5).

The spatial distribution of $n_{{\rm H},j}$ is characterized by the ratio of the mean free path of Lyman lines $l\propto 1/n_{\rm H,1s}$ to the ionization length of hydrogen atoms $L_{\rm I}\propto 1/n_{\rm p}$. Other series hardly contribute to the occupation number because $n_{{\rm H},j\neq{\rm 1s}}/n_{\rm H,1s}\ll 1$. The ratio,

[TABLE]

indicates that the spatial distribution of $n_{{\rm H},j\neq{\rm 1s}}$ has also self-similarity with respect to $n_{\rm tot,0}z$. In addition, from Eq. (4) or Eq. (28), we recognize the scaling relation $n_{{\rm H},j\neq{\rm 1s}}\propto n_{\rm tot,0}{}^{2}$. Thus, for fixed $T_{\rm p}$, $\beta$ and $\chi_{0}$, $n_{{\rm H},j\neq{\rm 1s}}/n_{\rm tot,0}{}^{2}$ has the same value at the same $n_{\rm tot,0}z$. Note that this scalability of the excited state may be kept up to $n_{\rm tot,0}\sim 10^{4}~{}{\rm cm^{-3}}$ in which the occupation number of 2s state is comparable with the ground state (i.e. our formulation would not be valid).

Fig. 3 represents the spatial distribution of $n_{{\rm H},j}/n_{\rm tot,0}{}^{2}$ for given $k_{\rm B}T_{\rm p}=12~{}{\rm keV}$ and $\beta=0.1$ with fixed values $\chi_{0}=$(0.1, 0.5 and 0.9), while Fig. 4 represents the case of given $k_{\rm B}T_{\rm p}=12~{}{\rm keV}$ and $\chi_{0}=0.5$ with fixed values $\beta=$(0.01, 0.1, 0.5 and 1). In this case, the maximum value of each occupation number varies $\sim 10$ per cent due to the variation of $\beta$. Note that although there are no energetic particles leaking and no extra radiation sources (e.g. thermal emission from the SNR ejecta), the shock-precursor-like structure is formed by Lyman lines (see the bottom panel of Fig. 3).

Fig. 5 shows the collisional excitation rates, $n_{\rm H,1s}C_{j,k}/n_{\rm tot,0}{}^{2}$, and radiative excitation rates, $n_{\rm H,1s}P_{j,k}/n_{\rm tot,0}{}^{2}$, for the transitions from 1s to 2p (Ly $\alpha$), from 1s to 3p (Ly $\beta$) and from 1s to 4p (Ly $\gamma$), respectively. The radiative rates are comparable with the collisional rates even in the upstream region.

Fig. 6 shows a comparison of $n_{\rm H,2s}$ from numerical calculation and $n_{\rm H,2s}$ from the simple estimate for the optically thick limit (${\rm e}^{-\tau_{0}}=0$ in Eq. (4)). The numerical value of $n_{\rm H,2s}$ less than the estimated $n_{\rm H,2s}$ in the optically thick limit implies the incomplete conversion of Lyman lines to Balmer lines. Fig. 7 represents intensities of narrow Ly $\alpha$, Ly $\beta$ and Ly $\gamma$ in the direction $\mu=\pm 1$,

[TABLE]

where

[TABLE]

Lyman photons escape toward the downstream region ($\mu=1$) from the shock. We define the fraction of escaping photons as

[TABLE]

where

[TABLE]

is the energy flux in the narrow line. Note that the sum of $|F^{\rm N}(z_{\rm out})|$ and $|F^{\rm N}(z_{\rm in})|$ indicates the net energy taken away from the shock by the lines.

The top panel of Fig. 8 shows the escape fraction of the Lyman lines for a shock with $k_{\rm B}T_{\rm p}$ and $\beta=0.1$. The fraction of Ly $\beta$ depends somewhat strongly on the ionization degree $\chi_{0}$ compared to the others. This is because the difference in optical thickness. The bottom panel of Fig. 8 shows the optical depth of Ly $\beta$ at the line centre in the direction of $\mu=1$. The depths of Ly $\alpha$ and Ly $\gamma$ are $(f_{\rm 1s,2p}/f_{\rm 1s,3p})\simeq 5.26$ and $(f_{\rm 1s,4p}/f_{\rm 1s,3p})\simeq 0.367$ times the Ly $\beta$ depth, respectively. The value of Ly $\beta$ optical-thickness decreases from $\simeq 5.2$ to $\simeq 0.14$ with increasing of $\chi_{0}$, while the thickness of Ly $\alpha$ and Ly $\gamma$ vary from $\simeq 27$ to $\simeq 0.74$ and from $\simeq 1.9$ to $\simeq 0.05$, respectively. The variations of attenuation, $e^{-\tau}$, are from $5.5\times 10^{-3}$ to $0.87$ for Ly $\beta$, from $1.9\times 10^{-12}$ to $0.47$ for Ly $\alpha$, and from $0.15$ to $0.95$ for Ly $\gamma$. Thus, the Ly $\beta$ attenuation varies the most drastically, giving the strongest dependence of the escape fraction on $\chi_{0}$. Note that the fraction would substantially increase roughly as the velocity width of narrow hydrogen atoms $\propto T_{0}{}^{\frac{1}{2}}$ because the optical depth is function of $T_{0}{}^{-\frac{1}{2}}$. Thus, the conversion of narrow Lyman lines to narrow Balmer lines is incomplete.

Fig. 9 represents the absorption coefficient of H $\alpha$ at the line centre of the ray directed along $\mu=0.01$ for given $k_{\rm B}T_{\rm p}=12~{}{\rm keV}$ and $\beta=0.1$ with fixed values $\chi=$( 0.1, 0.2, 0.3, 0.5, 0.7 and 0.9). The value of the coefficient ranges from $10^{-23}$ to $10^{-22}~{}{\rm cm^{-1}}$, which is consistent with our simple estimate using Eq. (2). Note that the absorption coefficient of H $\beta$ is about quarter of $k_{\nu,\mu}({\rm H\alpha})$; $k_{\nu,\mu}({\rm H\beta})\simeq(f_{\rm 2s,4p}/f_{\rm 2s,3p})k_{\nu,\mu}({\rm H\alpha})$.

5 synthetic observation

Here we consider how lines are observed at an SNR shock based on the calculated atomic populations. We suppose that the SNR shock propagates into a realistic ISM, which consists of diffuse gas and clumpy gas (e.g. Heiles & Troland, 2003), and the interaction between the shock and the inhomogeneous medium results in H $\alpha$ filaments on the sky (e.g. Hester, 1987; Shimoda et al., 2015, see also Fig. 10 of Inoue & Inutsuka 2012). Then, we presume that the shock becomes partly plane-parallel due to the density contrast, and that the length of the sides of the plane corresponds to the typical length scale of the density contrast. Note that the length scale may be at least $10$ per cent of the SNR radius as indicated by the rippling of the observed H $\alpha$ filaments. Moreover, we regard each plane as isolated for simplicity. In the following, our line of sight is fixed orthogonally to the shock normal (i.e. to along the $y$-axis and $\mu=0$).

Let $L_{\rm cl}$ be the extent of shock along our line of sight. Thus, the observed intensity (surface intensity) of each line is written as

[TABLE]

Note that roughly speaking, the source function hardly depends on the frequency but the absorption coefficient has sharp peak around the line centre. Therefore, for the optically thick case, the line shape is somewhat flattened from the centroid frequency to the critical frequency at which $k_{\nu,0}L_{\rm cl}\sim 1$, that is, $I_{\nu,0}\approx S_{\nu,0}=j_{\nu,0}/k_{\nu,0}$. For the frequency corresponding to $k_{\nu,0}L_{\rm cl}\la 1$, i.e. in the optically thin case, the line shape follows the line profile function; $I_{\nu,0}\approx S_{\nu,0}k_{\nu,0}L_{\rm cl}$. Note that we have the scaling relations $S_{\nu,0}\propto~{}n_{\rm tot,0}$ for the Lyman lines, $S_{\nu,0}\propto~{}n_{\rm tot,0}{}^{0}$ for the Balmer/Paschen lines, and $S_{\nu,0}k_{\nu,0}=j_{\nu,0}\propto n_{\rm tot,0}{}^{2}$. In the following, we fix $L_{\rm cl}=2~{}{\rm pc}$. Note that for SN 1006, Raymond et al. (2007) estimated the ambient density $n_{\rm tot,0}=0.25\mathchar 45\relax 0.4~{}{\rm cm^{-3}}$, the ionization degree $\chi_{0}\sim 0.1$ and the length of line of sight $\sim 2\times 10^{18}~{}{\rm cm}$, thus the optically thin condition is implied. It is possibly consistent with the observed width of narrow component at SN 1006, $\simeq 21~{}{\rm km~{}s^{-1}}$, that is exceptionally narrow compared to the other SNRs (Sollerman et al., 2003).

The observed H $\alpha$ spectra,

[TABLE]

are displayed in Fig. 10 in hydrogen velocity space for given $k_{\rm B}T_{\rm p}=12$ keV and $\beta=0.1$ with fixed values $n_{\rm tot,0}=$(1 cm*-3*, 20 cm*-3*, 40 cm*-3*, 60 cm*-3* and 80 cm*-3*). Note that we take the interval of spatial integration from $z_{\rm out}$ to $z_{\rm in}$ for simplicity. For shocks propagating into a tenuous medium ($\sim 0.1~{}{\rm cm^{-3}}$), the length of the precursor-like emission of H $\alpha$ becomes several times $10^{16}~{}{\rm cm}$ (see Figs. 5 or 7). Therefore, in typical spectroscopic observations of Galactic SNRs with a slit width of $\simeq 1~{}{\rm arcsec}\simeq 10^{16}~{}{\rm cm}(d/1~{}{\rm kpc})$, the photo-precursor emissions would be missed. Here we define the line centre ($v_{\rm H}=0$) by the mean value of relevant centroid frequencies $\nu_{j,k}^{\prime}$, $\nu_{\rm centre}=(\nu^{\prime}_{\rm 2s,3p}+\nu^{\prime}_{\rm 2p,3s}+\nu^{\prime}_{\rm 2p,3d})/3\simeq 4.56680\times 10^{14}~{}{\rm Hz}$. The two–photon emission dominates at velocity $|v_{\rm H}|\ga 3000~{}{\rm km~{}s^{-1}}$ in this case. Note that in reality, the broad component would be becoming fainter at $|v_{\rm H}|\ga 2000~{}{\rm km~{}s^{-1}}$ because of the velocity dependence of the charge-exchange cross-section, that is not accounted for in our model. The bottom panel of Fig. 10 shows the narrow component. The line shape evolves asymmetrically with increasing $n_{\rm tot,0}$. This evolution is enhanced (suppressed) for the case of lower (higher) $\chi_{0}$. Fig. 11 shows the cases of $\chi_{0}=0.1$ and $\chi_{0}=0.9$.

These modifications of the line profile from the distribution function of hydrogen atom come from the existence of the 2s-state hydrogen atoms and the radiation transfer effects. The asymmetry of the line profile arises from the offset of centroid frequencies, for instance, $1-\nu_{\rm 2s,3p}^{\prime}/\nu_{\rm centre}\simeq-0.8\times 10^{-5}$ (equivalently, $v_{\rm H}\simeq+2.3~{}{\rm km~{}s^{-1}}$).

The top panel of Fig. 12 shows the absorption and emission coefficients of H $\alpha$ just behind the shock front for a given $k_{\rm B}T_{\rm p}=12$ keV, $\beta=0.1$, $\chi_{0}=0.5$ and $n_{\rm tot}=1~{}{\rm cm^{-3}}$. The absorption coefficient has a peak at $v_{\rm H}\simeq 2.3~{}{\rm km~{}s^{-1}}$ (i.e. $\nu\simeq\nu^{\prime}_{\rm 2s,3p}$). It is obvious because almost all of the H $\alpha$ is scattered by the 2s-state hydrogen atoms. On the other hand, the emissivity $j_{\nu,0}$ has the peak at $v_{\rm H}\simeq 0$ ($\nu\simeq\nu_{\rm centre}$). Note that the value of peak frequency depends on the relative occupation number among 3s, 3p and 3d states. Thus, with increasing the optical depth, the line shape becomes more asymmetric due to the efficient scattering around $v_{\rm H}\simeq 2.3~{}{\rm km~{}s^{-1}}$, which is not aligned with the peak frequency of emissivity. In the optically thick limit, the line shape follows the source function that is shown in the bottom panel of Fig 12. The asymmetry of source function comes from the offset between the peak frequencies of $j_{\nu,0}$ and $k_{\nu,0}$. Note that as a consequence of the modification of line profile, the line width is apparently broadened.

The modifications of other narrow Balmer lines are modest compared with the H $\alpha$ because of the small absorption coefficient. Fig. 13 shows the case of H $\beta$.

Fig. 14 shows the frequency-integrated intensity of each line,

[TABLE]

for given $k_{\rm B}T_{\rm p}=12~{}{\rm keV}$, $\beta=0.1$, $\chi_{0}=0.5$ and $n_{\rm tot,0}=1~{}{\rm cm^{-3}}$. Here we also display the intensities of narrow component,

[TABLE]

and broad component,

[TABLE]

The narrow Lyman lines are mostly absorbed, while the broad Lyman lines are transparent. The Balmer lines are also optically thin in this tenuous-gas case ($n_{\rm tot,0}\ll 30~{}{\rm cm^{-3}}$). The Pa $\alpha$ is always in the optically thin limit. There is the shock-precursor-like emission due to the leaking of Lyman photons (see the bottom panel of Fig. 14). Obviously, the length of photo-precursor emission corresponds to the mean free path of the Lyman photons, which depends on $n_{\rm tot,0}{}^{-1}$. We display the intensity ratio of the photo-precursor emission to the downstream emission,

[TABLE]

for H $\alpha$, H $\beta$ and Pa $\alpha$ in Fig. 15.

The relative intensity of photo-precursor emission is increasing with increasing ionization degree $\chi_{0}$ because the Lyman series lines emerging in the downstream region tend to leak more to the upstream region. Thus, the shock-precursor-like emission is ubiquitously observed even if there are no leaking particles.

We define the total observed intensity of the narrow component as

[TABLE]

Similarly, that of broad component is defined as

[TABLE]

The intensity ratio of broad H $\alpha$ to narrow H $\alpha$ is often relied on to estimate $\beta$. Note that since our model simplifies the broad hydrogen atoms, we predict only how the broad-to-narrow ratio depends on $\chi_{0}$ and $\beta$. The ratio is shown in Fig. 16 for given $\chi_{0}=0.5$ and $n_{\rm tot,0}=1~{}{\rm cm^{-3}}$ with fixed values $\beta_{0}=$(0.01, 0.1, 0.5 and 1), while Fig. 17 shows the ratio for given $\beta=0.1$ and $n_{\rm tot,0}=1~{}{\rm cm^{-3}}$ with fixed values of $\chi_{0}=$(0.1,0.5 and 0.9).

The broad-to-narrow ratio of H $\alpha$ also depends on $\chi_{0}$ due to the conversion of narrow Ly $\beta$ to narrow H $\alpha$ but the dependence on $\chi_{0}$ is modest compared with the case of $\beta$ at lower $T_{\rm p}$.

Because of the mildly high opacity of H $\alpha$ for $n_{\rm tot,0}\ga 10~{}{\rm cm^{-3}}$, the peak of narrow H $\alpha$ is reduced in comparison of the case of $n_{\rm tot,0}=1~{}{\rm cm^{-3}}$. Therefore, the intensity ratio of broad H $\alpha$ to narrow H $\alpha$ depends on the ambient density.

Fig. 18 shows the broad-to-narrow ratio of H $\alpha$ as function of $n_{\rm tot,0}$ for given $k_{\rm B}T_{\rm p}=12~{}{\rm keV}$ and $\beta=0.1$ with fixed values $\chi_{0}=$(0.1, 0.5, 0.7 and 0.9).

Fig. 19 shows the ratio of the total narrow intensity of H $\beta$ to that of H $\alpha$, $\zeta_{0}^{\rm N}({\rm H\beta})/\zeta_{0}^{\rm N}({\rm H\alpha})$, i.e. the Balmer decrement. Note that the ratio of photon counts shown in Shimoda et al. (2018) is given by $(\nu_{\rm 3p,2s}^{\prime}/\nu_{\rm 4p,2s}^{\prime})\zeta_{0}^{\rm N}({\rm H\beta})/\zeta_{0}^{\rm N}({\rm H\alpha})\simeq 0.74\zeta_{0}^{\rm N}({\rm H\beta})/\zeta_{0}^{\rm N}({\rm H\alpha})$555 The Balmer decrement shown in Shimoda et al. (2018) was implicitly defined by $I_{\nu}/h\nu$ to compare with the observation by Sparks et al. (2015). . The ratio depends on how many Ly $\beta$ and Ly $\gamma$ photons are converted to H $\alpha$ or H $\beta$ photons. Since the escape fraction of Ly $\beta$, $\varepsilon_{\rm esc}({\rm Ly\beta})$, depends strongly on the ionization degree $\chi_{0}$ compared to the fraction of $\rm Ly\gamma$ (see Fig. 8), in other words $\left|\frac{\partial\zeta_{0}^{\rm N}({\rm H\alpha})}{\partial\chi_{0}}\right|>\left|\frac{\partial\zeta_{0}^{\rm N}({\rm H\beta})}{\partial\chi_{0}}\right|$, a larger $\chi_{0}$ results in a less conversion of Ly $\beta$ to H $\alpha$ than the case of Ly $\gamma$ to H $\beta$ conversion. Thus, the observed ratio $\zeta_{0}^{\rm N}({\rm H~{}\beta})/\zeta_{0}^{\rm N}({\rm H~{}\alpha})$ would decrease with increasing $\chi_{0}$. On the other hand, the ratio also depends on the interval of spatial integration of Eq. (46). For example, when we take the interval from $z=0$ to $z=z_{\rm in}$ (i.e. observe the downstream region), the ratio increases with increasing $\chi_{0}$, with values around $0.4\mathchar 45\relax 0.5$.

Fig. 20 shows the Balmer decrement for given $\chi_{0}=0.5$ and $n_{\rm tot,0}=1~{}{\rm cm^{-3}}$ with fixed values of $\beta=$(0.01, 0.1, 0.5 and 1). The ratio varies with various $\chi_{0}$ and $\beta$ with $\sim 10$ per cent. Note that the values of ratio around $0.4\mathchar 45\relax 0.5$ seem to be consistent with observed value in SN 1006, $\sim 0.5\pm 0.1$ (Raymond et al., 2017). On the other hand, the ratios $\la 0.3$ are observed by Ghavamian et al. (2001, 2002) for Tycho’s SNR and SN 1006 (at different position from Raymond et al. (2017)). The small value of ratio could be explained if the SNR shock suffers extreme energy losses (Shimoda et al., 2018, see the results of $\eta>0.5$).

The difference of absorption coefficients between H $\alpha$ and H $\beta$ leads to a variation of the Balmer decrement $\zeta_{0}^{\rm N}({\rm H\beta})/\zeta_{0}^{\rm N}({\rm H\alpha})$ with the ambient density $n_{\rm tot,0}$.

Fig. 21 shows the Balmer decrement as function of $n_{\rm tot,0}$ for the given $k_{\rm B}T_{\rm p}=12~{}{\rm keV}$ and $\beta=0.1$ with fixed values $\chi_{0}=$(0.1, 0.5, 0.7 and 0.9). Note that both variations of $\zeta_{0}^{\rm N}({\rm H\beta})/\zeta_{0}^{\rm N}({\rm H\alpha})$ (see Fig. 18) and $\zeta_{0}^{\rm B}({\rm H\alpha})/\zeta_{0}^{\rm N}({\rm H\alpha})$ on $n_{\rm tot,0}$ result from the scattering of narrow H $\alpha$. Therefore, they resemble each other.

6 Summary and Discussion

We have studied the radiative transfer of hydrogen lines taking into account the atomic population. The absorption of Lyman photons results in the emission of Balmer photons and yields 2s-state hydrogen atoms. We have shown that a fraction of the Lyman photons escape from the shock toward the far downstream region and therefore the Lyman-Balmer conversion is incomplete, which is consistent with the results of Ghavamian et al. (2001). As a result, the observed intensity ratio of narrow H $\beta$ to narrow H $\alpha$ (the Balmer decrement) and the observed broad-to-narrow ratio of H $\alpha$ vary $\sim 10$ per cent depending on the ionization degree of upstream medium. Note that the broad-to-narrow ratio of H $\alpha$ is often relied upon to derive the ion-electron temperature equilibrium in the downstream region rather than the ionization degree. The intensity ratios depend also on the ambient density because of the scattering of H $\alpha$ by the 2s-state hydrogen atoms if the shock propagates into the medium with density of $\ga 10~{}{\rm cm^{-3}}$. In the case of optically thick H $\alpha$, the line shape of narrow H $\alpha$ becomes asymmetric as a consequence of the differences in atomic energy levels. The degree of asymmetry increases with increasing the optical depth. Since the narrow H $\alpha$ can be scattered, the observed Balmer decrement and broad-to-narrow depend also on the ambient density.

We have supposed that the H I cloud in the path of the shock is clearly separated into regions of dense clumps and a diffuse, tenuous component and that the shock-fronts propagating into each part are isolated. Moreover, we have fixed the temperature for H I gas. However, the temperature depends on the density in reality (e.g. Field, 1965). The clumpy gas (‘cold neutral medium’) is usually supposed to have the temperature of $\sim 50\mathchar 45\relax 100$ K with the density of $\sim 10\mathchar 45\relax 100~{}{\rm cm^{-3}}$, while the diffuse gas (‘warm neutral medium’) is supposed to have the temperature of $\sim 6000\mathchar 45\relax 10000$ K with the density of $\sim 0.2\mathchar 45\relax 0.5~{}{\rm cm^{-3}}$ (e.g. Ferrière, 2001, and references therein). On the other hand, motivated by investigations of the molecular cloud formation via accretion flows of H I clouds (i.e. star formation), recent numerical simulations show that during the formation, the density and temperature of H I cloud range over two to three orders of magnitude (e.g. Inoue & Inutsuka, 2012, 2016, see also Figures 5 and 6 of Fukui et al. 2018). The probability distribution function in the density-temperature plane presented in Fukui et al. (2018) shows the existence of H I clouds with density of $\sim 10\mathchar 45\relax 100~{}{\rm cm^{-3}}$ and temperature of $\sim 10^{3}\mathchar 45\relax 10^{4}~{}{\rm K}$. The size of such a cloud is $\sim 1~{}{\rm pc}$, while the size of denser clouds ($\ga 100~{}{\rm cm^{3}}$ and $\sim 100~{}{\rm K}$) is $\sim 0.1$ pc. The separation of these clouds is typically $\sim 1~{}{\rm pc}$, but sometimes $\sim 0.1~{}{\rm pc}$. Thus, for the case of actual SNR shocks, the densest clump, tenuous gas and intermediately dense cloud can be co-existing on our line of sight. Note that such H I clouds associated with the molecular cloud formation may be a minor component of the ISM, but it may be possible that the SNR shock propagates in such H I regions because of the very long dynamical time of ISM, $\sim(10~{}{\rm pc}/10~{}{\rm km~{}s^{-1}})\sim 10^{6}~{}{\rm yr}$.

We have not considered any other radiation sources, especially the radiation from supernova ejecta. Moreover, we neglect the leaking of broad hydrogen atoms and the existence of cosmic-rays. The cosmic-rays or neutral leakage precursors are accelerated relative to the pre-shock gas, giving a velocity offset as an extra parameter in the radiative transfer calculation (Boulares & Cox, 1988). Note that the broad hydrogen atoms leaking to the upstream region may lead to the modification of shock structure as well as that arising in a cosmic-ray modified shock (Blasi et al., 2012; Ohira, 2012, 2013, 2016a, 2016b). Furthermore, to quantify the production of non-thermal particles, polarization measurements are additionally required (Shimoda et al., 2018). We will address these issues (realistic ISM, radiation from supernova ejecta, leaking particles and polarization) in forthcoming papers. In principle, the 2s1/2 population may also be quenched by mixing with the 2p1/2 level in a motional electric field (e.g. Drake, 1988). In an electric field of $E=10^{-6}$ StatVolt cm*-1*, arising from motion with a velocity 3000 km s*-1* through a magnetic field of $10^{-4}$ G, the mixing amplitude evaluates to $2\times 10^{-6}\left(E/10^{-6}\right)$, leading to a negligible increase in the decay rate of $\left(2\times 10^{-6}\right)^{2}\times 6.25\times 10^{8}=0.0025$ s*-1*. A factor of 60 increase in electric field is required to make this new contribution comparable with the two-photon decay rate of $\sim 8$ s*-1*, which for now seems to be out of reach for plausible cosmic ray generated magnetic fields and SNR shock velocities (e.g. Vink & Laming, 2003).

Although our model is still too simple to compare line profiles quantitatively with actual observations, it can predict different line profiles between H $\alpha$ and H $\beta$. This qualitative prediction is based on the atomic physics (i.e. oscillator strength $f_{j,k}$) rather than the property of BDSs. In particular, if the anomalous width of H $\alpha$ reflects the intrinsic velocity dispersion of hydrogen atoms, both of H $\alpha$ and H $\beta$ should have the same width. On the other hand, in the case of scattering, the width should depend on the direction of our line of sight (i.e. optical depth). For $\mu\simeq 0$ ($\mu\simeq 1$), the optical depth tends to be maximum (minimum). Interestingly, the width of narrow H $\alpha$ observed in SNR 0509 67.5 depends on the line of sight (Smith et al., 1994). The width observed at the west rim of SNR ($\mu\simeq 0$) is $\simeq 30\pm 2~{}{\rm km~{}s^{-1}}$, while at the centre of SNR ($\mu\simeq 1$) is $\simeq 25\pm 2~{}{\rm km~{}s^{-1}}$. It is qualitatively consistent with the difference of optical depth. Note that Smith et al. (1994) suggested the existence of intermediate component for the west rim. This scenario can be tested by future observations of H $\beta$ with high-resolution spectroscopy.

Long et al. (1992) reported the two-photon continuum emission from SNR Cygnus Loop at which $k_{\rm B}T_{\rm p}\sim 0.1$ keV (Medina et al., 2014). Unfortunately, since there are no cross-section data on proton impact excitation to $n_{k}=4$ at a relative velocity less than $1000~{}{\rm km~{}s^{-1}}$, our present calculation is limited at $k_{\rm B}T_{\rm p}\ga 5$ keV. According to Figs. 10 and 13, the two-photon continuum can be a diagnostic of BDSs. To do this, we should treat the broad hydrogen atoms self-consistently (e.g., see Raymond et al., 2008; Blasi et al., 2012; Ohira, 2012) and consider any other radiation sources that yield continuum components. We will extend our model to provide a diagnostic from the two-photon continuum in future work.

Acknowledgements

We thank Dr. Makito Abe for valuable comments that helped us complete this work. We also thank the referee, John Raymond, for his comments further improve the paper. This work is partially supported by JSPS KAKENHI grant no. JP18H01245. JML was supported by the Guest Investigator Grant HST-GO-13435.001 from the Space Telescope Science Institute and by the NASA Astrophysics Theory Program (80HQTR18T0065), as well by Basic Research Funds of the CNR.