Heating and Cooling

tempera-31 2.2 Methods ture ofT = 10⁴K (Liffman 2003). It is smaller than the sink size of our computational domain for a 0.5M_⊙ star. Thus, we could miss the mass loss due to photoevaporative flows in the sink region, if excited. However, the base density profile of the ionized flows is anticipated to scale with∝R^−1.5 (Tanaka et al. 2013), where R is the distance in the cylindrical polar coordinates. (We will show below that the profile actually has a similar scaling.) The outer region has a more contribution than the inner region with this scaling. The contribution from the sink region is a tiny fraction of the total photoevaporation rates. This is confirmed by additional test simulations where the inner boundaries are set torinner= 0.1, 0.35, 0.5 au; the simulation results have shown that the mass loss from the inner region (R ≤10 au) is responsible for only about a few percents of the total. This justifies using the inner boundary slightly larger than the effective gravitational radius.

Note that the inner (< 1 au) disk could have an important effect on shielding the direct EUV photons. If the inner disk in the sink region were able to completely attenuate the direct stellar photons, only the diffuse photons would reach the outer disk surface in the computational domain.

We have found, in the simulations with the smaller inner boundaries, that the direct EUV actually reaches the > 1 au region, and heats/ionizes the gas there as in the simulations with the fiducial value of the inner boundary. The resulting photoevaporation rates thus hardly vary with the sink size, at least for 0.1 au≤r_inner≤1 au.

Chapter 2 UV Photoevaporation of PPDs: Metallicity Dependence 32 Eq.(2.19) provides the photoionization rate and the associated heating rate

RIonize=yHI

∫ _∞

ν1

dν σνFν, (2.21)

Γ_EUV= n_HI ρ

∫ _∞

ν₁

dν σ_νh(ν−ν₁)F_ν, (2.22)

respectively. We assume a black body spectrum for the SED of the EUV. The effective tem-perature is set to Teff = 10⁴K, which yields the total stellar EUV emission rate of ΦEUV ≃ 1.5×10⁴¹(R_∗/ R_⊙)²s⁻¹. The frequency-integration range is divided into 81 bins.

In our simulations, the photoelectric heating rate is provided by the analytic formula of Bakes and Tielens (1994) (hereafter BT94). The Mathis-Rumpl-Nordsieck (MRN) distribution (Mathiset al.

1977) is assumed for the dust model to derive the heating rates. The size distributions of small carbon grains including polycyclic aromatic hydrocarbons (PAHs) follow the MRN distribution in the model. We note that observations suggest significantly lower PAH abundances around T Tauri stars by typically several tens of times than the ISM value (Geers et al.2007, Oliveiraet al.2010, Vicenteet al. 2013), although there are large uncertainties in the estimated values. We investigate the effects of the PAH abundance on disk photoevaporation rates in Section 2.4.3. The photoelectric heating function is

Γpe= 10⁻²⁴erg s⁻¹ϵpeGFUV

n_H ρ

Z

Z_⊙, (2.23)

ϵpe=

[ 4.87×10⁻²

1 + 4×10⁻³γ_pe^0.73 +3.65×10⁻²(T /10⁴K)^0.7 1 + 2×10⁻⁴ γpe

]

, (2.24)

whereϵ_pe is the photoelectric effect efficiency andγ_pe is the ratio of the dust/PAH photoionization rate to the dust/PAH recombination rate. The ratio is provided by γ_pe ≡ G_FUV√

T /˜˜ n_e, where T˜ = (T /1 K) and ˜ne = (ne/1 cm⁻³). The photoelectric effect efficiency measures the ratio of the gas heating rate to FUV absorption rate of the grains. In Eq.(2.23), G_FUV is the FUV flux (6 eV< hν <13.6 eV) normalized by the averaged interstellar fluxFISRF= 1.6×10⁻³erg cm⁻²s⁻², and is calculated asG_FUV=L_FUVe^−1.8AV/(4πr²F_ISRF). The last factorZ/Z_⊙in Eq.(2.23) accounts for the reduction of the grain amount with decreasing metallicity.

For cooling processes, we include radiative recombination cooling of H II (Spitzer 1978), Lyα cooling of HI(Anninoset al.1997), fine-structure line cooling of OIand CII(Hollenbach and McKee 1989, Osterbrock 1989, Santoro and Shull 2006), molecular line cooling of H₂ and CO (Galli and Palla 1998, Omukaiet al.2010), and dust-gas collisional cooling (Yorke and Welz 1996). We neglect other collisional excited lines (CELs) which can be coolants in HIIregions. Possible influences this simplification could yield are discussed in Section 2.4.6.

Dust-Gas Collisional Cooling

Heat is transferred via collisions between dust and gas. This effect works as cooling or heating for gas. We use the dust-gas collisional cooling function presented in Yorke and Welz (1996)

Λdust=−4πa²_dustcs

nH

ρ

(ρdust

m_dust )

k(T−Tdust) Z

Z_⊙, (2.25)

wherea_dust, ρ_dust, m_dust, andT_dustare the mean dust size, dust density, mean dust mass, and dust temperature, respectively. The parameters are set toadust= 5×10⁻⁶cm andmdust= 1.3×10⁻¹⁵g (Yorke and Welz 1996).

33 2.2 Methods

Atomic and Molecular Line Cooling

Radiative recombination cooling of HII, Lyαcooling of HI, C IIline cooling, OI line cooling, H₂ line cooling, and CO line cooling are implemented as the coolants of gas.

In HII regions, free electrons recombine to various energy levels of hydrogen. The electrons in excited states quickly go down to lower-energy states with emitting lines (radiative recombination).

This works as cooling in HIIregions. Approximately two-thirds of the electron energy∼kT is lost through this process (e.g., Spitzer 1978). In the simulations, the radiative recombination cooling rate is given as

Λrec= 0.67kT Rk2nenHIIρ⁻¹, (2.26) whereRk2 is the recombination coefficient (the reaction labeled “k2” in Table 2.2).

In a hot environment where atomic hydrogen exists, the energetic collisions excite neutral hydrogen atoms, and then the atoms de-excite through Lyαphoton emission. Hereafter, we simply call it Lyα cooling. We adopt the Lyαcooling function of in Anninoset al.(1997)

ΛLyα=ξLyαnenHIρ⁻¹, (2.27)

ξLyα= 7.5×10⁻¹⁹e^−118348/TK 1 +√

T_K/10⁵ erg cm³s⁻¹, (2.28)

whereTK is gas temperature in the unit of Kelvin.

Table 2.3 Fine-Structure Line Parameters of C II and O I

Species j→i νij(Hz) Aij( s⁻¹) γ_ij^e ( cm³s⁻¹) γ_ij^HI( cm³s⁻¹) Reference C II 2→1 1.9×10¹² 2.4×10⁻⁶ 2.8×10⁻⁷T₂^−0.5 8.0×10⁻¹⁰T₂^0.07 1,2

O I 2→1 4.7×10¹² 8.9×10⁻⁵ 1.4×10⁻⁸ 9.2×10⁻¹¹T2^0.67 1,3 O I 3→1 – 1.0×10⁻¹⁰ 1.4×10⁻⁸ 4.3×10⁻¹¹T2^0.80 1,3

O I 4→1 – 6.3×10⁻³ 1.0×10⁻¹⁰ 1.0×10⁻¹² 1,3

O I 5→1 – 2.9×10⁻⁴ 1.0×10⁻¹⁰ 1.0×10⁻¹² 1,3

O I 3→2 2.1×10¹² 1.7×10⁻⁵ 5.0×10⁻⁹ 1.1×10⁻¹⁰T2^0.44 1,3

O I 4→2 – 2.1×10⁻³ 1.0×10⁻¹⁰ 1.0×10⁻¹² 1,3

O I 5→2 – 7.3×10⁻² 1.0×10⁻¹⁰ 1.0×10⁻¹² 1,3

O I 4→3 4.7×10¹⁴ 7.3×10⁻⁷ 1.0×10⁻¹⁰ 1.0×10⁻¹² 1,3

O I 5→3 – 0 1.0×10⁻¹⁰ 1.0×10⁻¹² 1,3

O I 5→4 5.4×10¹⁴ 1.2 0 0 1,3

Note: Here, T2 ≡ T /100 K. The labels i, j indicate the energy levels of the species, the corresponding frequency of the energy difference between levelsiandjis represented byνij,Aijis an Einstein A coefficient, andγij^λ is the rate of collisions with a speciesλ. We define the labels of energy levels as following:²P1/2 of C II (label 1), and²P3/2 of C II (label 2), respectively. ³P2 of O I (label 1),³P1of O I (label 2),³P0 of O I (label 3),¹D2 of O I (label 4), and¹S0 of O I (label 5), respectively.

Reference: (1) Osterbrock (1989), (2) Santoro and Shull (2006), (3)Hollenbach and McKee (1989).

Spontaneous emission associated with the fine-structure transitions of C II and O I works as coolants. The total cooling rate of each species is given as

ΛX =nX

ρ

∑

j

xj

∑

j>i

Aji∆Eji. (2.29)

Chapter 2 UV Photoevaporation of PPDs: Metallicity Dependence 34 The labelX indicates either of CII or OI,xj is the level population atj, and ∆Eji(∆Eji=hνji) is the energy difference between levelsiandj. The level populations are determined by solving the simultaneous equations of statistical equilibrium

x_i∑

j̸=i

c_ij =∑

i̸=j

x_jc_ji, (2.30)

where

c_ij ≡







Aij+∑

λγ_ij^λnλ for (i > j)

∑

λγ_ij^λnλ for (i < j)

. (2.31)

We treat these line emissions as optically thin for simplicity in this study. Therefore, we omit the calculations of an escaping probability, excitation due to absorption of external radiation, or de-excitation due to induced radiation in Eq.(2.31). Regarding O I cooling, we do not explicitly take into account O I photoionization in our chemistry model. Atomic oxygen might be ionized in H II regions, and this would reduce the O I cooling rate. We approximately incorporate this effect into the O I cooling rate by setting O I abundance to yOI(1−yHII) in H II regions. This approximation could be justified by the similar ionization potentials of atomic hydrogen and atomic oxygen. It would be necessary to treat OIphotoionization in detail to model the chemical structure of oxygen. However, adiabatic cooling is dominant in HIIregions, and thus this simplification would have little effect on our results.

Line emission associated with rovibrational transitions of H₂ and CO molecules also work as coolants for gas. The cooling rates can be calculated in the same manner as C II and O I line cooling with Eq.(2.29), but in order to save computational cost, we give the cooling rates by using the analytic formula of Galli and Palla (1998) for H2 cooling and the tabulated values provided by Omukaiet al.(2010) for CO cooling.