Line radiative transfer and statistical equilibrium

Atomic and molecular line emission from protoplanetary disks contains key information of their detailed physical and chemical structures. To unravel those structures, we need to understand line radiative transfer in dusty media and the statistical equilibrium, especially of molecules. I describe here the basic principles of statistical equilibrium and illustrate them through the two-level atom. In a second part, the fundamentals of line radiative transfer are introduced along with the various broadening mechanisms. I explain general solution methods with their drawbacks and also specific difficulties encountered in solving the line radiative transfer equation in disks (e.g. velocity gradients). I am closing with a few special cases of line emission from disks: Radiative pumping, masers and resonance scattering.


Introduction
The previous chapter "Line Observations in Disks" (Dionatos 2015) introduced the different degrees of freedom a molecule has (rotation, vibration, electronic levels) and the respective nomenclature for atomic and molecular lines.Spectra of molecules are particularly powerful for deducing physical conditions in the observed regions, e.g. also in protoplanetary disks.The role of various molecules and specific types of transitions is nicely summarized in Table 5.1 of Stahler & Palla (2004).For disks, interesting diagnostics are CS and CN, which are high density probes, and thus are -especially in their higher rotational lines -less contaminated than CO by emission from the surrounding cloud.
A rule of thumb to remember is that lines preferentially originate close to where the critical density is reached and also close to where the optical depth of the line becomes one.The critical density n cr is defined as the ratio between the Einstein coefficient for spontaneous emission A ul and the collision rate between the upper and lower level of the line transition C ul n cr = A ul C ul . (1) If the density is higher than the critical density, the level populations will follow the Boltzmann distribution due to efficient collisional coupling and the line will be in Local Thermodynamic Equilibrium (LTE).
The interplay between critical density and excitation temperature in defining the optimum emitting conditions for various types of transitions and molecules is shown in figure 1 in the chapter by Dionatos (2015).It nicely illustrates that the various types of CO lines ranging from low rotational lines to high rotational lines and then ro-vibrational transitions have the potential to probe the entire disk from the low density cold outer regions to the very hot, dense inner regions (see Fig. 1).
Figure 1.Disk sketch with the various molecular and atomic line emitting regions.

Statistical Equilibrium
Line fluxes depend on the level populations, or more directly on the column density of an atom/molecule in the upper energy level.We will simplify the radiative transfer for a moment and only consider the probability, β, that a line photon, emitted somewhere in a column of gas that we study, escapes the medium.The escape probability depends on the optical depth τ ul at line center frequency ν ul .We will explain that concept later in more detail (see also chapter by Woitke 2015) and also get back to the full radiative transfer.Then, we can write the relation between the column density of the upper level, N u , and the line flux F ul as where A ul is the Einstein coefficient for spontaneous emission and Ω source the solid angle of the line emitting region.This equation directly illustrates the importance of determining the level populations to calculate the column densities and hence compute line fluxes.In the following, we will discuss how to calculate these level populations in LTE and non-LTE (statistical equilibrium, SE).

Local Thermodynamic Equilibrium
In the case of LTE, the level populations within a specific ionization state of an atom or molecule are given by the Boltzmann equation EPJ Web of Conferences where n u and n l are the level populations of the upper and lower energy level and g u and g l their statistical weights, respectively.ΔE = E u − E l is the energy difference between the two levels, k the Boltzmann constant and T is the gas temperature.Note that even if LTE does not hold, one can still define an excitation temperature T ex that relates the two level populations through

Ionization Equilibrium
The ionization balance of atoms and molecules in LTE follows from the Saha equation where n(X i ) is the level population in the ionization state X i of the atom/molecule, n e is the electron density, χ i,i+1 is the ionization energy from state X i to X i+1 , m is the mass of the atom/molecule, and T is the gas temperature.The partition function Z(X i ) of the ionization state X i is defined as

Non Local Thermodynamic Equilibrium
If the assumption of LTE breaks down, we have to solve the equations for non-LTE statistical equilibrium (SE).This means that we need to take into account all radiative and collisional processes that govern the distribution of level populations.The radiative processes are subject to the quantum mechanical selection rules and comprise spontaneous emission, stimulated emission and absorption, which depend on the respective Einstein coefficients.Collisional coupling can occur between any two levels; however, collision cross sections are generally larger between closely spaced energy levels and dipole allowed transitions.For each energy level, we can balance the rate for processes populating that level and those de-populating that same level.In equilibrium, dn/dt = 0 holds, where n is the vector containing all level populations n i .The complete set of equations for SE becomes then The first term sums over all emissions from higher levels ending up in level i, the second term over all absorptions from lower levels into level i.The third term adds up all collisional rates C ji from lower and higher levels into level i.The fourth to sixth term add up the respective radiative and collisional depopulation of level i. P(ν i j ) denotes the local mean intensity averaged over ray directions n and local line profile functions φ i j (ν, n) as Summer School "Protoplanetary Disks: Theory and Modeling Meet Observations‰ which can cause absorption and stimulated emission.I ν ( n) is the spectral intensity at frequency ν in direction n.The three different Einstein coefficients for spontaneous emission A ji , stimulated emission B ji , and absorption B i j are coupled through the following relations and Note that the exact form of these relations (first term of constants) depends on the units in which the radiation field is given.This set of equations can be written in vector form as where n is an N-dimensional vector composed of the individual level populations n i for all N levels and M is an N × N matrix composed of elements M i j .For example, the column u with u > i contains the elements The largest amount of work lies in the compilation of all the atomic/molecular data, i.e. the level energies, statistical weights, radiative transitions, Einstein A coefficients, collision cross sections.Several databases provide this data such as the LAMDA database (Schöier et al. 2005), NIST, and CHIANTI (Dere et al. 1997).

Solving the Statistical Equilibrium
Once all radiative and collisional data is compiled for a specific atom/molecule, Eq. ( 11) can be inverted and solved for n.We need to respect, however, a very important boundary condition, namely the particle conservation that ensures that the sum of all level populations n i adds up to the total density of that atom/molecule n tot .The way to impose this implicitly in the solution of the SE equation is to replace one of the equations in Eq. ( 11) with the conservation equation The matrix M then becomes This matrix equation can be solved using standard numerical recipes to invert the matrix such as LU decomposition.The problem scales with N 3 , making evident why often limited atoms/molecules are used, i.e. the amount of ionization stages and levels being restricted according to the physical conditions in the region to be investigated.An example is the restriction to pure rotational energy levels of the ground electronic and ground vibrational level of CO in the case of molecular clouds, where densities and temperatures are low (n ∼ 10 4 cm −3 , T ∼ 10 K).

The Two-level Atom
To provide better insight into the SE, we will in the following consider the two-level atom as an example.The atom has two energy levels with populations n 0 and n 1 , which are connected through a single line with frequency ν 10 that has an Einstein coefficient A 10 .The SE in Eq. ( 7) becomes in that case We can re-arrange this equation to isolate the ratio between the two level populations -making also use of the relations between the Einstein A and B coefficient - We see immediately that in the absence of a radiation field, P(ν 10 ) = 0, this simplifies to Using the relation between the upwards and downwards collision rates we find for very large collision rates (C 10 A 10 ) that this equation becomes the Boltzmann equation.Hence, in that case, the two level populations are in LTE according to the local gas temperature T gas .An example of a simple two level atom is the two fine structure levels of the ground state of ionized carbon: Lower level 2 P 3/2 and upper level 2 P 1/2 .They are separated by an energy E 10 of 92 K and the emission line has a wavelength of 157.7 μm.Since the critical density of this line is very low, it is widely observed in the interstellar medium (ISM, e.g.Pineda et al. 2013).
A second case to be considered is a situation where the collision rates are much smaller than the radiative rates (C 10 → 0 and C 01 → 0).In that case, Eq. ( 16) becomes If P(ν 10 ) is a blackbody radiation field of temperature T gas , this becomes Note, that we never imposed LTE in this case, on the contrary, we started with negligible collision rates.Hence, radiation can also impose LTE on the level populations provided that the radiation field is that of a blackbody.Of course, this will be a valid assumption in the optically thick regions of the disk.

Radiative Transfer
Now that we understand how the different rotational, vibrational and electronic levels of an atom/molecule are populated, we can address the second part of understanding line emission: how does radiation cross the medium to reach the observer.From the equations of SE it became already evident that the background radiation field can pump level populations under certain conditions.Hence, the two problems of SE and radiative transfer are closely intertwined, making the combined problem very difficult to solve even numerically.We will come back to that later when we address simplifications made in disk research.The next few paragraphs will first illustrate the basic theory of line emission/absorption and line broadening.

Line Emission Coefficient
The line emission coefficient i j ν is related to the Einstein A coefficient, the level population of the upper level n i and the line profile function φ i j (ν) The line profile function describes how a line at central frequency ν i j is broadened due to the uncertainty principle, collisions, thermal and turbulent motions in the gas, see Sect.3.3.The line can also be Doppler-shifted by systematic motions of the gas in the observers frame, which causes the line profile function to become direction dependent φ i j = φ i j (ν, n) which is, however, not considered in this section.The line profile itself is normalized to one, φ i j (ν) dν = 1, with the effect that the total line emission is smeared out over a certain frequency range symmetric around the rest frequency of the line ν i j .A useful characterization of a line profile -also in observations -is the Full Width Half Maximum (FWHM), which measures the width of the profile at half its peak value.

Line Absorption Coefficient
The line absorption coefficient α i j ν has two components, stimulated emission B i j and absorption B ji .In that case, the sign of the absorption coefficient immediately hints to the presence of laser/maser 1 emission, which will be discussed at the end of this chapter.The line absorption coefficient can be expressed as where n i and n j are the upper and lower level populations respectively and φ i j (ν) is again the line profile function.

Line Broadening
The total line profile is a convolution of profiles from various types of broadening mechanisms: natural broadening, collisional broadening, thermal and turbulent broadening.The first two are related to the lifetime of the upper level of the transition and lead to a Lorentzian profile, while the second two cause Doppler shifts of the frequency and hence produce Gaussian profiles.The convolution of a Lorentzian and a Gaussian profile is a Voigt profile.

Natural Broadening
Natural broadening is caused by the uncertainty principle.The shorter the lifetime of a state, the larger the energy uncertainty of the state and the natural line width thus becomes The width of the profile is related to the spontaneous de-excitation rate from level i into all other levels j, γ = j A i j (Fig. 2).With this, the Lorentzian line profile can be written as (indices i, j omitted)

Collisional Broadening
Atoms and molecules in a gas collide frequently and this reduces the lifetime in the different energy levels (Fig. 2).The broadening Δν is proportional to the collision frequency ν coll .Hence, collisional broadening becomes more important at higher densities.The line profile from the combined effect of natural and collisional broadening is still a Lorentzian with the width Δν tot being the sum of the two widths (γ + 2ν coll ).
In most astrophysical cases, natural broadening is negligible.Collisional broadening plays mostly a role at very high densities, such as the ones reached in stellar atmospheres 10 12 cm −3 .This leaves thermal and turbulent broadening as the dominant mechanisms under ISM, molecular cloud and protoplanetary disk conditions.

Thermal Broadening
The thermal velocities of atoms and molecules show a random Gaussian 3D velocity distribution (Fig. 3), a Maxwell-Boltzmann distribution Summer School "Protoplanetary Disks: Theory and Modeling Meet Observations‰ Here, f v denotes the probability of finding an atom/molecule with a velocity close to v. The other quantities are the mass of the atom/molecule m and the gas temperature T gas .Due to these random velocities, the frequency of the line ν i j is shifted according to the projected radial velocity component v r in the line of sight to an observer The resulting Doppler frequency shift is The most probable velocity of an atom/molecule of mass m in a gas of temperature T gas is A particle with this velocity has an energy of kT gas .A simple expression to remember is the most probable velocity of a hydrogen atom From this, one can easily scale for atoms/molecules of different mass or a gas of different temperature.
The line profile resulting from thermal Doppler broadening is a Gaussian profile with the width Δν D (Eq.29).

Turbulent Broadening
In the presence of turbulent motions v turb , the width of the profile becomes The root mean square of the thermal and turbulent velocities is often denoted by the parameter b in studies of the Interstellar Medium (ISM).Typical turbulent velocities for disks are inferred to be below 1 km/s.The Gaussian profile is normalized to one.Its peak value at ν i j and FWHM can be calculated from the width

Voigt profiles
A Voigt (V) profile is the convolution of a Lorentzian (L, width Δν tot ) and a Gaussian (G, width Δν D ) profile EPJ Web of Conferences line "resonance" region, pumping the levels (Avrett & Hummer 1965).In 1D plane-parallel geometry, the intensity at the surface I ν is related to the intensity at the bottom I ν (0) through Roughly speaking, e −τ gets replaced by the escape probability β(τ), where τ is the line optical depths along the considered ray direction.This approximation has problems if the geometry is multidimensional (no plane-parallel geometry) and the medium becomes optically thin.It will be discussed in more detail in the chapter by Woitke (2015).Finally, the line emission is calculated by performing ray tracing through the computational volume and using the level populations derived from the escape probability approach.

Velocity gradients
A protoplanetary disk has a Keplerian velocity field.This leads to double peaked velocity profiles as described in the seminal paper of Beckwith & Sargent (1993) and detailed in the chapter of Dionatos (2015).If a line photon is emitted in cell r n , it is shifted according to the Keplerian velocity of that cell.If the photon travels a distance s = r n − r (n−1) from one cell in the midplane to the next, this corresponds to a difference in velocities v(r n ) − v(r n−1 ), and hence to a frequency shift of If the frequency shift is larger than the line width, the photon is 'shifted out of the line' and can escape.In Interstellar Medium physics, this is often referred to as the Large Velocity Gradient (LVG) approximation.It resembles the escape probability approximation discussed above since the escape probability can be written as where τ LVG is the total optical depth along the path for any frequency.In disks, however, the LVG approximation does not work since a line photon can interact with the gas at various locations along a ray (e.g.top and bottom of the disk, near-side and far-side, Fig. 5, left panel).This problem is discussed in Pontoppidan et al. (2009) in detail for near-IR lines originating in the inner disk where velocity gradients across a cell can be large.In that case, the opacity needs to be properly sampled despite a finite cell size.To illustrate this, consider a typical turbulent/thermal broadening of < 1 km/s and a grid with a resolution of 0.01 AU at a distance of 1 AU from the star.The velocity gradient across one cell is then ∼ 4 km/s for a solar mass star.One possibility to 'catch' the line within the cell is a subsampling of the source function between the cell boundaries (Fig. 5, right panel).

Special cases
After having discussed the foundations of the statistical equilibrium and radiative transfer, we will turn to a few special cases to illustrate the basic concepts.

Radiative pumping
The molecules in the disk surface are directly exposed to the stellar radiation field.This means that besides their statistical equilibrium being governed by the warm gas temperatures in the surface, radiative pumping can be very strong if the molecule has strong transitions in the UV/optical wavelength range, where the stellar radiation field peaks.An example is the CO molecule, which has electronic bands in the UV (4th positive system: A 1 Π − X 1 Σ + around 1600 Å). Figure 6 illustrates the UV pumping of the upper electronic state and the cascade leading eventually to the emission of the CO ro-vibrational lines from the ground electronic state in the near-IR (∼ 4.6 μm).
The stellar radiation field enters into the equations of SE through the absorption term j<i n j B ji P(ν ji ).The UV fluorescence can pump the populations in the upper vibrational levels of the ground electronic state to values much higher than LTE.Hence, the band ratios between the v =1-0, 2-1, 3-2 bands will be affected and the hotter bands (v =2-1, 3-2) will get stronger with respect to the fundamental band (v = 1-0).The mechanism counterbalancing this fluorescence is collisional Summer School "Protoplanetary Disks: Theory and Modeling Meet Observations‰ 00010-p.13quenching (i.e.collisional de-population of the higher v-levels) due to high densities.The main collision partner in the emitting region is atomic hydrogen (Thi et al. 2013) whose collision rates are two orders of magnitude larger than those of molecular hydrogen.
Another example of radiative pumping is the water molecule.It has key transitions in the mid-IR where the spectral energy distribution of protoplanetary disks peaks.This means that the IR continuum photons (e.g. at 79 and 90 μm) can efficiently pump low rotational levels of water.We saw in Sect.2.1 that we can define an excitation temperature T ex for any two levels.In the most extreme case of radiative pumping, the excitation temperature of the line would approach the local dust temperature (optically thick case) or more precise the radiation temperature of the continuum.This has been illustrated in the study of the water statistical equilibrium and radiative transfer in the disk around TW Hya (Kamp et al. 2013).

Masers
Maser (microwave amplification by the stimulated emission of radiation) emission is frequently observed in star forming regions and is often associated with massive star formation and H ii regions.Maser emission is also seen in the vicinity of AGB stars.Table 1 gives an overview of detected maser transitions from molecules and references to the literature.
In order to observe a maser, we need an inversion of the level populations, so a collisional and/or radiative pumping mechanism.For the case of dusty environments around AGB stars, some possible pumping mechanisms for the 22 GHz water maser are described by Babkovskaia & Poutanen (2006).Characteristics of maser emission are large intensities, very narrow line width and unusual line ratios.The latter indicates very extreme non-LTE conditions.If multiple maser lines of the same molecule are observed from the same physical region, the excitation mechanism and the small scale structure and physical conditions in the gas can be disentangled.

Resonance scattering
An important example for resonance line scattering in disks is Ly α (Bergin et al. 2003;Neufeld 1990).
Figure 7 shows the escape fraction of a Ly α photon from a one dimensional slab model as a function of the atomic hydrogen column density N(H 0 ).The two processes of dust absorption and collisional de-excitation compete in destroying Ly α photons.The solid curves show an increasing ratio of the scattering coefficient for the line photons and the dust absorption coefficient.With less efficient dust absorption, a higher fraction of line photons can escape the slab.The dashed curves illustrate the effect of decreasing collisional de-excitation.Resonance scattering of Ly α photons into the disk can increase the depth into the disk to which molecules can be photo dissociated.This occurs e.g. for molecules that have their dissociation bands overlapping with Ly α such as CO, OH, H 2 O and HCN (Fogel et al. 2011).Other molecules such as CN are unaffected.However, part of this can be compensated by an increased photodesorption rate also due to the enhanced Ly α radiation.
9 th Lecture of the Summer School "Protoplanetary Disks: Theory and Modelling Meet Observations" DOI: 10.1051/ C Owned by the authors, published by EDP Sciences,

1
Laser stands for Light Amplification by Stimulated Emission of Radiation and maser for Microwave Amplification by Stimulated Emission of Radiation.EPJ Web of Conferences 00010-p.6

Figure 3 .
Figure 3. Left panel: Velocity distribution for hydrogen atoms at two different gas temperatures.Shown in color are the most probable, average and root mean square velocity of the distribution.Right panel: Comparison of a Lorentzian, Gaussian and Voigt profile.

Figure 5 .
Figure 5. Left panel: Channel map of a line that is red-shifted by 3 km/s in a disk -visible are the front and back, but also the lower and upper disk part.Right panel: Sketch illustrating the concept of subsampling in the case of large velocity gradients and a finite cell size.Figures from Pontoppidan et al. (2009, c AAS, reproduced with permission).

Figure 6 .
Figure 6.Left panel: Sketch of the CO fluorescence mechanism.Right panel: Sketch of the ortho-and parawater energy levels together with some key lines that have been observed with the Herschel satellite.Figures are from Thi et al. (2013, reproduced with permission c ESO), van Dishoeck et al. (2011).

Figure 7 .
Figure 7. Escape fraction of a Ly α photon from a 1D slab model.Figure adapted from Neufeld (1990, c AAS, reproduced with permission).

Table 1 .
Examples of Maser transitions found in star forming environments and the vicinity of AGB stars.