Constraints on metal oxide and metal hydroxide abundances in the winds of AGB stars - Potential detection of FeO in R Dor

Using ALMA, we observed the stellar wind of two oxygen-rich Asymptotic Giant Branch (AGB) stars, IK Tau and R Dor, between 335 and 362 GHz. One aim was to detect metal oxides and metal hydroxides (AlO, AlOH, FeO, MgO, MgOH), some of which are thought to be direct precursors of dust nucleation and growth. We report on the potential first detection of FeO (v=0, Omega=4, J=11-10) in RDor (mass-loss rate, Mdot, ~1e-7 Msun/yr). The presence of FeO in IK Tau (Mdot~5e-6 Msun/yr) cannot be confirmed due to a blend with 29SiS, a molecule that is absent in R Dor. The detection of AlO in R Dor and of AlOH in IK Tau was reported earlier by Decin et al. (2017). All other metal oxides and hydroxides, as well as MgS, remain undetected. We derive a column density N(FeO) of 1.1+/-0.9e15 cm^{-2} in R Dor, or a fractional abundance [FeO/H]~1.5e-8 accounting for non-LTE effects. The derived fractional abundance [FeO/H] is a factor ~20 larger than conventional gas-phase chemical kinetic predictions. This discrepancy may be partly accounted for by the role of vibrationally excited OH in oxidizing Fe, or may be evidence for other currently unrecognised chemical pathways producing FeO. Assuming a constant fractional abundance w.r.t. H_2, the upper limits for the other metals are [MgO/H_2]<5.5e-10 (R Dor) and<7e-11 (IK Tau), [MgOH/H_2]<9e-9 (R Dor) and<1e-9 (IK Tau), [CaO/H_2]<2.5e-9 (R Dor) and<1e-10 (IK Tau), [CaOH/H_2]<6.5e-9 (R Dor) and<9e-10 (IK Tau), and [MgS/H_$]<4.5e-10 (R Dor) and<6e-11 (IK Tau). The retrieved upper limit abundances for these latter molecules are in accord with the chemical model predictions.


INTRODUCTION
The gas-phase elements Ca, Fe, Mg, Si and Ti are depleted w.r.t. the solar abundances in diffuse clouds. The formation of metal oxides and metal hydroxides and of dust species is suggested as major cause for this depletion. Indeed, a variety of metal oxides and hydroxides are prominent in a wide range of temperature and density environments. The metal oxides TiO, VO, CrO, YO, ZrO are present in the atmospheres of cool M stars (see, e.g., Scalo & Ross 1976;Castelaz et al. 2000). SiO, TiO, TiO 2 , AlO, and AlOH are detected in the winds of oxygen-rich Asymptotic Giant Branch (AGB) stars (e.g. Schöier et al. 2004;Decin et al. 2010;De Beck et al. 2015;Kamiński et al. 2016;De Beck et al. 2017;Kamiński et al. 2017;Decin et al. 2018). Other metal oxides and hydroxides such as CaO, CaOH, FeOH, MgO, and MgOH have been searched in molecular clouds and stars (e.g. Hocking et al. 1979;Kamiński et al. 2013;Sánchez Contreras et al. 2015;Quintana-Lacaci et al. 2016;Velilla Prieto et al. 2017) without success. Laboratory measurements show that FeO can be formed at high temperatures (≥1000 K, Cheung et al. 1982) and hence could be abundant in the atmospheres and inner winds of AGB stars. However, until now, FeO has only been detected in interstellar space in absorption along the line of sight toward the galactic center HII region Sagittarius B2 Main (Sgr B2 M) (Walmsley et al. 2002;Furuya et al. 2003). This detection is interpreted as due to shocks associated with star formation which might liberate some fraction of gas-phase elements from the refractory grains. FeO has remained, however, undetected in stellar atmospheres and stellar winds. Here we report the first potential detection of FeO in the stellar wind of the low mass-loss rate AGB star R Dor (see Sect. 2). No spectral features of other metal oxides and hydroxides (CaO, CaOH, MgO, MgOH) and MgS have been seen in the winds of the two oxygen-rich AGB stars R Dor and IK Tau surveyed with ALMA. In Sect. 3 we show that FeO only accounts for a tiny fraction of the solar iron abundance and derive upper limit abundances for the undetected metal species.In Sect. 4 we discuss the derived abundances in the framework of local thermodynamic equilibrium (LTE) and pulsation shock induced non-equilibrium chemistry models for the stellar winds.

OBSERVATIONS
We used ALMA to observe the high mass-loss rate AGB star IK Tau (Ṁ ∼5×10 −6 M /yr) and the low mass-loss rate AGB star R Dor (Ṁ ∼1×10 −7 M /yr). Data were obtained in August-September 2015 in Band 7 (335-362 GHz) with a spatial resolution of ∼150 mas (proposal 2013.1.00166.S, PI L. Decin). Data reduction was done using CASA (McMullin et al. 2007) and is described in detail in Decin et al. (2018). The spectral restoring beam parameters are in the range of 120−180 mas for IK Tau and 130−180 mas for R Dor. The channel σ rms noise varies between spectral windows and is in the range of 3−9 mJy for IK Tau and 2.7−5.7 mJy for R Dor. The velocity resolution is 1.6−1.7 km/s for IK Tau and 0.8−0.9 km/s for R Dor.
Some two hundred spectral features from 15 molecules were identified. Detected species include the gaseous precursors of dust grains such as SiO, AlO, AlOH, TiO, and TiO 2 (Decin et al. , 2018. 66 lines remain unidentified, some of which may belong to OH and H 2 O (Decin et al. 2018) or higher excitation rotational transitions not included in the current spectral line catalogues of the Jet Propulsion Laboratory (JPL, Pickett et al. 1998) and the Cologne Database for Molecular Spectroscopy (CDMS, Müller et al. 2001Müller et al. , 2005Endres et al. 2016). The rest frequencies of the unidentified features were carefully compared to the predicted line frequencies of various metal oxides and hydroxides. Rotational transitions of CaO, CaOH, MgO, MgOH (and MgS) do not correspond to any of the unidentified lines 1 . However, one of the unidentified spectral features in R Dor has a central frequency around 336.815 GHz, with a minor blend at the blue side due to a TiO 2 line at 366.8241 GHz (see Fig. 1). We attribute this feature to the FeO (v=0, Ω = 4, J=11-10) transition in the ground electronic 5 ∆ i state with rest frequency 336 816.030±0.05 MHz (Allen et al. 1996). This is the only (potential) FeO line detected in our ALMA data (see also Sect. 3). We can not confirm if this FeO transition is present in IK Tau due to a blend with the strong 29 SiS (v = 1, J = 19−18) line at 336.815 GHz. However, SiS (and CS) are 1 The rotational spectrum of FeOH in the X 6 A i ground state has not been measured in the laboratory. Detailed quantum mechanical calculations indicate the dipole moment of 1.368 Debye is favourable, but FeOH is quasilinear with a small barrier to linearity of less than 300 cm −1 and its spectrum may be complex (Hirano et al. 2010). absent 2 in R Dor (see Decin et al. 2018), enabling the possibility of detecting and identifying this rotational transition of FeO.
3. ANALYSIS AND RESULTS 3.1. The FeO (v = 0, Ω = 4, J = 11−10) transition The ground state of FeO is 5 ∆ i in Hund's case a (Cheung et al. 1981(Cheung et al. , 1982. As such there are five spin-orbit components separated by intervals of 190 cm −1 , which are labelled by the quantum number Ω = Λ + Σ. With Λ and Σ both equal to 2, vector addition gives possible values for Ω between 4 and 0, with the 5 ∆ 4 component (Ω = 4) lying lowest in energy. The (sub)millimeter wave spectrum of FeO in the X 5 ∆ i state (for v = 0) has been measured up to 400 GHz in the laboratory by Allen et al. (1996).
The rotational transitions for J = 11−10 of all five Ω components lie in the range of the observed ALMA frequencies (see Table 1 in Allen et al. 1996). Only the J = 11−10 in the lowest energy Ω = 4 component at 336.816030 GHz corresponds to a spectral feature in the ALMA spectrum of R Dor (see Fig. 1 and the channel map in App. A). The higher excitation Ω components remain undetected. This is unfortunate, since the detection of more than one transition of FeO would strengthen its identification. However, this outcome is not completely unexpected since transitions between the spin-orbit components are highly forbidden due to the very strong case a coupling ) and hence do not support the argument that radiative transitions could cause significant transfer of population between the Ω-ladders of the X 5 ∆ ground state. In this high density inner wind region, collisions might however pump the population to the higher energy Ω-levels.
Other rotational transitions in the Ω = 4 spin-orbit component lie outside the observed ALMA frequency range 3 . The only other transition of FeO hitherto detected in interstellar space belongs to the (v = 0, Ω = 4) spin-orbit component as well, and is the lowest (J = 5−4) rotational transition at 153.135 GHz towards the galactic center HII region Sgr B2 M (Walmsley et al. 2002;Furuya et al. 2003). FeO has been searched for in stellar winds since the early 1980s  without success. If the spectral feature near 336.815 GHz is indeed caused by FeO (and a minor blend with a TiO 2 line), this would be the first detection of FeO in the wind of an evolved star.
The ALMA data here offer the possibility of estimating the column density of FeO, N (FeO). We therefore need to subtract the contribution of the TiO 2 (23 8,16 − 23 7,17 ) transition which is slightly blended with the FeO line in the blue wing. We therefore have selected three other TiO 2 transitions with almost equal quantum numbers, excitation energies, line strengths, and channel maps: TiO 2 (24 2,22 −23 3,21 ) Figure 1. Continuum subtracted ALMA spectrum for R Dor around 336.8 GHz in black extracted for a circular beam with aperture of 300 mas. Two lines of previously identified molecules are observed: SO2 (ν2=1,201,18), and TiO2 (238,16 − 237,17) indicated by the vertical dotted lines in black. The feature at 336.81603 GHz, tentatively identified as the FeO (v = 0, Ω = 4, J = 11−10) transition, is blended with the line of TiO2 at 336.82407 GHz. The flux of the TiO2 line 8 MHz higher in frequency than FeO was estimated by referring to three lines of TiO2 with similar excitation energies and channel maps that were observed in the same program (see inset in the right hand side of the panel). From these a red synthetic profile for the TiO2 line was derived. Plotted in blue is the observed profile minus the synthetic TiO2 profile attributed to FeO. at 347.788 GHz, TiO 2 (26 0,26 − 25 1,25 ) at 350.399 GHz, and TiO 2 (25 2,24 − 24 1,23 ) at 350.708 GHz (see Fig. 1). An average flux density of these three lines scaled to the same peak flux was calculated and fitted using a 'soft-parabola' function (see Eq. 1 in Decin et al. (2018); see red line in Fig. 1 here), that was then subtracted from the ALMA data. The resulting spectrum (for a circular aperture with beam of 300 mas) is shown in blue in Fig. 1. The peak flux of the FeO line is 0.093 Jy and the integrated line flux is 0.79 Jy km/s. Correcting for a local standard of rest velocity, v LSR , of 7 km/s (Decin et al. 2018), the half-width of the line at zero intensity, ∆v, is 7.5 km/s. The FeO emission is essentially unresolved with the ALMA beam of ∼150 mas (see the channel map in Fig. 6 and discussion in App. B) giving us an upper limit for the detected FeO emission of radius 2.5 R .

FeO column density
We use a population diagram analysis to estimate the column density of FeO, N (FeO). Assuming the emission is optically thin, we derive that the column density in the upper (J = 11) state, N thin u , is 1.7×10 13 cm −2 (see App. B). Assuming local thermodynamic equilibrium (LTE), N (FeO) can be calculated as with g u the degeneracy of the upper rotational state, Q the partition function, and E u the energy level of the upper state. The energy of the J = 11 rotational level in the Ω = 4 ladder (E u = 57.203 cm −1 ) was calculated with the spectroscopic constants in Allen et al. (1996). For a constant excitation temperature, T x , the partition function for a linear diatomic molecule is given by (Tielens 2010) with B the rotational constant of FeO in units of Hz; B = 15493.63255 MHz (Allen et al. 1996). A first estimate on the excitation temperature, T x , can be obtained from calculating the (upper limit) of the crossladder temperature using the fact that the FeO (v = 0, Ω = 3, J = 11−10) at 338.844 GHz is undetected in our survey. For an energy difference between the Ω-ladders of 190 cm −1 ) and a σ rms of 4 mJy (hence 3σ rms of 12 mJy), we derive that the cross-ladder temperature, T CL x , is <130 K. If the cross-ladder populations were controlled by collisions, T CL x would be a direct measure of the kinetic temperature, T kin , provided cross-ladder radiative transitions are negligible, and a lower limit if not (Thaddeus et al. 1984). The estimated kinetic temperature in the region between 1 and 2.5 R ranges between 1300-2400 K . Since cross-ladder transitions are only weakly permitted (see Sect. 3.1), the cross-ladder temperature is expected to be higher than the temperature within the ladders because populations (within a ladder) rapidly decay by emission of millimeter-wave photons (Thaddeus et al. 1984). Since other rotational transitions in the Ω = 4 ladder are not in the frequency window of the ALMA data, we can not calculate the excitation temperature within a ladder. Thaddeus et al. (1984) were able to derive the 'within' and 'cross'ladder temperature for the X 1 A 1 SiCC molecule -whose ∆K a = 2 electric dipole transition moments are smallin the carbon-rich AGB star CW Leo, being 10 K and 140 K respectively. Using a lower limit of 10 K for the excitation temperature of FeO seems unreasonably low, since FeO is detected in the inner wind of R Dor (r ≤ 2.5 R ) in contrast to detection of SiCC in the outer wind of CW Leo. We henceforth assume a lower limit for the excitation temperature of 80 K. Assuming that the excitation temperature can be as high as 2000 K, the derived column density of FeO, N (FeO), varies between 2 × 10 14 and 2 × 10 15 cm −2 , or Using the equation of mass conservationṀ = 4πr 2 ρ(r)v(r), with ρ(r) the gas density and v(r) the gas velocity, one can calculate the H 2 density, n(H 2 ), assuming all hydrogen to be locked in H 2 . The H 2 column density, N (H 2 ), is dependent on the gas velocity v(r). As shown by Decin et al. (2018), our knowledge of the gas velocity in the inner wind of R Dor is limited. Using the velocity β-laws described by Decin et al. (2018) (their Eqs.(2)-(4)), we derive that N (H 2 ) is ∼2×10 22 cm −2 for a column with length between 1 and 2.5 R . Hence, the ratio N (FeO)/N (H 2 ) is estimated to be ∼5.5±4.5 × 10 −8 .
The FeO level populations might, however, violate the assumption of a Boltzmann distribution. There are strong transitions between the ground X 5 ∆ electronic state of FeO to the 5 Π and 5 Φ excited electronic states near 10 000 cm −1 (Cheung et al. 1982), and owing to the large spin-orbit coupling between these three states the X 5 ∆ state has some excited state character that may enhance ∆Ω ± 1 transitions between spin-components. Also rotational levels within one Ω ladder might be subject to radiative excitation effects. To check for the impact of the latter, we have calculated the frequencies, upper state energies, and Einstein A coefficients for the rotational transitions in the (v = 0, Ω = 4) spin component (see App. C To derive the abundance of FeO and the upper limit abundance of the undetected species, we use the same procedure as for the determination of the AlO, AlOH and AlCl abundances in R Dor and IK Tau outlined in Decin et al. (2017). In short, we modelled the ALMA data using a non-LTE radiative transfer model based on the Accelerated Lambda Iteration (ALI) method (Maercker et al. 2008) that allows us to derive the global mean molecular density assuming a 1D ge-ometry. The gas kinetic temperature and velocity have been approximated by a power law distribution. The gas density, ρ(r), is calculated from the equation of mass conservation. The fractional abundance of the Al-species w.r.t. H 2 was assumed either to be constant up to a certain maximum radius, R max , or to decline according to a Gaussian profile centred on the star for both targets f (r) = f 0 exp(−(r/R e ) 2 ), with f 0 the initial abundance and R e the e-folding radius.
Collisional excitation rates have not been published for these six molecules (as was the case for AlO, AlOH, and AlCl). Hence we have used the values from other molecules as substitutes, scaling for the difference in molecular weight. The HCN-H 2 system (Green & Thaddeus 1974) was used to extract the collisional rates for MgOH, CaOH (and AlOH) and the SiO-H 2 system (Dayou & Balança 2006) for FeO, CaO, MgO, MgS (and AlO, AlCl). The collisional rates are used in the form that appears in the LAMDA database 4 (Schöier et al. 2005). Einstein-A coefficients were calculated from the quantum-mechanical line strength, S, as given in CDMS (MgOH, CaOH, CaO), JPL (MgS, MgO), or App. C (FeO). The ALMA R Dor FeO channel map does not allow us to properly deconvolve the ALMA beam (see App. B). Assuming the ALMA emission to be essentially unresolved leads to a maximum extent of 150 mas in diameter (or 2.5 R in radius). Using the moment-0 maps would lead to a larger FeO source size with radius ∼6 R . Setting R e to 2.5 R (by analogy with AlCl; Decin et al. 2017) and R max to 6 R , the fractional abundance of FeO can be derived. However, for a turbulent velocity, v turb , of 1 km/s , the predicted FeO line profile is smaller than observed. As described by Decin et al. (2018), our understanding of the gas velocity in the inner wind region is limited. Pulsation-induced shocks might result in a radial velocity amplitude of a few km/s (Nowotny et al. 2010). Allowing the turbulent velocity to be 3 km/s permits the predicted line profile to reach the observed ∆v of 7.5 km/s (see Fig. 2). Using these parameters for the undetected species, the (upper limit) abundance of the metal species is derived (see Table 1). Since for IK Tau all metal species remain undetected, the upper limit abundances were only calculated for a constant abundance profile with v turb = 1 km/s and R max = 40 R (as determined from AlOH and AlCl, Decin et al. 2017).
A principal uncertainty in the abundance calculations concerns the unknown collisional rates. Changing the collisional rates by one order of magnitude only alters the retrieved abundances by 10% or less. The only exception is a lowering by 60% of the calculated FeO abundance if the collisional rates were a factor 10 lower.  (Cristallo et al. 2015); dust formation is not accounted for since we focus on the region where the bulk of the dust has not yet formed.
The gas-phase reaction rate coefficients are taken from the literature where available, and extrapolated to the high temperatures of an outflow using Transition State Theory (TST, Atkins 1998) with molecular constants (vibrational frequencies, rotational constants) calculated using quantum theory. The rate coefficients for reverse reactions were then calculated assuming detailed balance. The list of all chemical reactions involved is given in Table 3 in App. D.
The oxides of Ca, Fe, and Mg are produced by reactions with O 2 , CO 2 , H 2 O and OH releasing O, CO, H 2 , and H, respectively, and the hydroxides of Ca, Fe, and Mg are formed by reactions with H 2 O and OH. Moreover, the metal oxides (CaO, FeO, MgO) are linked to the hydroxides by reaction with molecular hydrogen H 2 . Generally, small (reduced) networks might introduce oversimplifications compared to extensive, complete reaction networks. However, the metallic Ca-Mg-Fe chemistry is largely decoupled from the remaining gas phase chemical families (e.g. sulphur, nitrogen, silicon). We also compared the modelled OH (and H 2 O) abundance with the study of Gobrecht et al. (2016) who used an extensive chemical network with 100 species and 424 reactions (including the N, S and Si chemistry). We find similar trends and absolute values of the OH abundance in both models. The reactions R7, R8 and R11-R16 have the largest impact on the OH chemistry and determine the H 2 O-OH balance. In addition, the abundances of the prevalent species CO, CO 2 , H 2 O, and OH agree with observations. The physical conditions experienced by the upper atmosphere of R Dor are described by a parcel of gas which is initially at rest at the photosphere and is in thermodynamic (thermal, chemical, radiative and mechanical) equilibrium. We assume that the stellar pulsation, originating from the interior of the star, has steepened in a shock and hits the gas parcel. As a consequence, gas in the cube is compressed, heated and accelerated outwards. The temperature and density profiles are calculated following Bertschinger & Chevalier (1985) for a 10 km/s shock and a diatomic gas with preshock conditions of T 0 = 2400 K and n 0 = 1×10 14 cm −3 (see, e.g., Fig. 10 in Freytag et al. 2017). Hydrodynamic calculations (Nowotny et al. 2010) have shown that the amplitude of the velocity variation, and hence the shock velocity, is slightly larger than the terminal wind velocity (being ∼5.5 km/s; see discussion in Decin et al. 2018). At the shock front, gas temperature and density take peak values of T = 3500 K and n = 6×10 14 cm −3 , respectively, and subsequently decrease with an exponential decay in the post-shock gas (see Fig. 3). We ran the gas-phase chemistry model over a full pulsation period (being 332 days, Bedding et al. 1998) and followed the change in the atomic and molecular abundance profiles over one pulsation phase (with the phase defined as the decimal part of (t − T 0 )/P , with T 0 the epoch of the start of the pulsation cycle and P the period); see Fig. 4. We note that the abundance variations within a pulsation period are much larger than the cycle-to-cycle variations in the periodic pulsation model.
To probe the accuracy of the chemical kinetics code and reaction network we performed runs at constant temperature (T = 2400 K) and density (1×10 14 cm −3 ), and compared the results to equilibrium abundances with the same conditions. For the large majority of molecules, the differences are small and within a factor of 2. Exceptions are O 2 , FeOH, MgO and MgOH that differ by factors of up to 6-9. We conclude that our small network describes the chemical behaviour in the inner wind of R Dor with sufficient accuracy. The relatively small differences between equilibrium and chemicalkinetic abundances may arise from the incomplete network or from an insufficient characterisation of some molecules (FeOH, MgOH) and their related reaction rate coefficients.   The predicted abundances for the metal oxides and metal hydroxides are shown in Fig. 4. In the immediate postshock region (pulsation phase, φ, 0.0-0.2) the dissociated molecules start to reform and reach peak abundances between 1.1×10 −12 ([MgO/H]) and 7.4×10 −10 ([FeO/H]) in the cooling post-shock gas (φ = 0.2-1.0). While the upper limit abundances for the undetected species in R Dor (Table 1) are in accord with the model results, the predicted FeO abundance is a factor 20 lower than derived from the ALMA data. The main processes leading to the formation of FeO are Fe + OH (R42 in Table 3) and Fe + H 2 O (R58). The first reaction dominates at early phases φ < 0.5, whereas the latter is only important at later phases φ >0.5. The main FeO destruction channel is the reaction FeO + H → Fe + OH (R43).
Provided the FeO identification is correct, a number of suggestions can be put forward to explain the discrepancy between observed and predicted FeO abundance. It might be that the chemical network is not complete or that the use of detailed balance to estimate some of the rate coefficients is not correct in this environment where molecular vibrational models may not be thermally equilibrated. Another possibility is the sputtering of dust grains, although this seems unlikely since the grains close to the star should be Fe-free silicates or alumina (Khouri et al. 2016) and sputtering products such as O will actually decrease FeO (R39). Although fresh molecular O 2 might react with Fe, which is abundantly present (3.1 × 10 −5 relative to H), R38 has a very large activation energy (see Table 3). We here propose an alternative scenario, following the idea of Elitzur et al. (1976), that the Fe + OH reaction (R42) might occur from vibrationally excited OH, where R42 would no longer be endothermic. We account for this possibility by reducing the activation barrier of R42 to zero. As a result, the FeO fractional abundance increases by a factor ∼4 (see dashed brown line in Fig. 4, denoted as FeO ).
The fraction of vibrationally excited OH in the inner wind of R Dor is, however, unknown. A first-order estimate could come from the assumption of a Boltzmann distribution of states, but this would not represent vibrational disequilibrium. The impact of the amount of vibrationally excited OH can be gauged by reducing the activation barrier, E a , in reaction (R42) stepwise from 3348 K to 0 K, the latter situation assuming all OH is vibrationally excited hence representing an upper limit (see Fig. 5). As expected, the [FeO/H] maximum increases, and the maximum value reached is 2.96 × 10 −9 for E a /R = 0 K (see dashed brown line in Fig. 4 and full brown line in Fig. 5), which is a factor ∼5 lower than observed. Accounting for the uncertainties of the thermodynamics properties of the inner wind region, we conclude that this may be a viable route for the formation of gaseous FeO. However, the results also showcase that using the best available chemical kinetics, FeO is hard to make at the level tentatively observed and presents an important challenge for future chemical models. Phase Φ E a /R= 0 K E a /R= 500 K E a /R=1000 K E a /R=1500 K E a /R=2000 K E a /R=2500 K E a /R=3000 K E a /R=3348 K Figure 5. Predicted FeO abundances with respect to the total gas number density as a function of pulsation phase φ at 1 R for different values of the activation barrier, Ea/R, in reaction R42.
Council (ST/P00041X/1), and CAG and KLKL from NSF grant AST-1615847. This paper uses the ALMA data ADS/JAO.ALMA2013.1.00166.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada) and NSC and ASIAA (Taiwan), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. This paper makes use of the CASA data reduction package: http://casa.nra.edu -Credit: International consortium of scientists based at the National Radio Astronomical Observatory (NRAO), the European Southern Observatory (ESO), the National Astronomical Observatory of Japan (NAOJ), the CSIRO Australia Telescope National Facility (CSIRO/ATNF), and the Netherlands Institute for Radio Astronomy (ASTRON) under the guidance of NRAO.
which in the Rayleigh-Jeans regime simplifies to with k the Boltzmann constant and λ the wavelength. The flux density S F ν of the source is defined as Ω here being the solid angle and S denoting that one integrates over the source solid angle S dΩ = Ω S . For a source of uniform brightness, T b (Ω) can be taken out of the integral and Eq. (B4) becomes The source size Ω S is dependent on the observed frequency ν as can be seen in Fig. 6. To derive Ω S for each frequency, we fitted a single 2D Gaussian component to the emission in each channel from the observed frequencies between 336.7006 and 336.7455 GHz. The emission is, however, too compact, patchy and irregular to deconvolve the beam reliably from the apparent angular size. This was close to the restoring beam (0. 175×0. 127) size and showed large uncertainties. We hence assume that the emission is unresolved, and use the beam size (Ω beam = Ω S = 0.175 × 0.127 × π/(4 ln(2)) = 0.025 arcsec 2 ) to estimate the brightness temperature 5 .
Defining W as with v the velocity (corrected for the v LSR ) and using Eq (B5) for an integrated flux density of the FeO (v = 0 Ω = 4 J = 11−10) line of 0.79 Jy km/s, yields W = 330 K km/s (for a source angular size of 150 mas in diameter or 2.5 R in radius). Assuming the emission is optically thin, one can write the column density in the upper state N u as (Goldsmith & Langer 1999) N thin with k the Boltzmann constant, c the velocity of light, and A i,j the Einstein-A coefficient for the transition. For a linear molecule, the Einstein-A coefficient for a rotational transition J → J − 1 is given by (Tielens 2010) with µ i,j the transition moment and h the Planck constant. For a transition J → J − 1, the transition moment and the quantum mechanical line strength S is calculated with the standard expression for a symmetric top by replacing K, the angular momentum along the symmetry axis of the symmetric top, by Ω the projection of the angular momentum along the molecular axis ) where µ is the permanent electric dipole moment. This yields A J,J−1 = 1.16 × 10 −11 µ 2 ν 3 J 2 − Ω 2 J(2J + 1) (B11) for ν in GHz and µ in Debye. Hence for W in units of Jy km/s. Using a permanent electric dipole moment of 4.7 Debye (Steimle et al. 1989), we obtain that N thin u = 1.7 × 10 13 cm −2 . We have calculated the rotational transition frequencies and upper state energies, E u , in the FeO (v = 0, Ω = 4) spin component by using the leading spectroscopic constants in Allen et al. (1996) and following Merer et al. (1982). The Einstein-A coefficients are calculated from Eq. (B11); see Table 2. By comparison to the results of Allen et al. (1996), who measured the FeO spectrum for frequencies lower than 400 GHz, the accuracy of the calculated frequencies is about 1 MHz. This is sufficient for our radiative transfer calculations, but we note that for spectroscopic identifications the required accuracy of the calculated frequencies in Table 2 CaO, CaOH, CaH, MgO, MgOH, MgH, FeO, FeOH, and FeH) species that take part in 59 reactions. The gas-phase reaction rate coefficients are taken from the literature where available, and extrapolated to the high temperatures of an outflow using Transition State Theory (TST, Atkins 1998) with molecular constants (vibrational frequencies, rotational constants) calculated using quantum theory. The rate coefficients for reverse reactions were then calculated assuming detailed balance.
Most of the uncertainty in the TST calculations arises from the accuracy of the energy barrier. For most of this work we have calculated the equilibrium constant and then used this for detailed balance where either the forward or backward reaction has been measured (i.e., barrier not needed). In some cases, a TST expression is fitted to a measured rate constant, with the barrier calculated from quantum theory adjusted to optimise the fit. This greatly improves the accuracy of the estimated rate constants. An example of our applied methodology and resulting accuracy can, e.g., be found in Self & Plane (2003) for reactions R38 and R39.
A special note concerns reactions R42 and R43. we have calculated the potential energy surface for these reactions. Fig. 7 illustrates the surface for a fixed Fe-O-H angle of 150 • . Note there are no barriers in the entry channels of either R42 (Fe + OH) or R43 (FeO + H). We have therefore set the barrier of R43 to zero in the expression estimated by Rumminger et al. (1999). R42 is endothermic by 25.5±5.8 kJ mol −1 , using the measured FeO bond energy (0 K) of 398.5±5.8 kJ mol −1 (Metz et al. 2005). The rate coefficient for R42 can then be calculated by detailed balance, yielding k 42 (T ).
The list of all chemical reactions involved is given in Table 3.  a The rates are given in the Arrhenius form k(T ) = A × T 300 n × exp(−Ea/T ), where T is the gas temperature, A the Arrhenius coefficient in cm 3 s −1 or cm 6 s −1 for a bimolecular or termolecular process, respectively, n the temperature dependence of the rate coefficient, and Ea is the activation energy barrier in K. b NIST: National Institute for Standards and Technology (http:// kinetics.nist.gov), TST estimate: estimation of the rate by Transition State Theory. Figure 7. Diagram of the potential energy surface for the Fe + OH → FeO + H reaction (R42), calculated for a fixed Fe-O-H angle of 150 • at the b3lyp/6-311+g(2d,p) level of theory using the GAUSSIAN 16 suite of programs (Frisch et al. 2016). Note the absence of barriers in the entrance channels for R42 (Fe + OH) or R43 (FeO + H). The deep well in the surface is due to the formation of FeOH.