Abundance ratios of OH/CO and HCO+/CO as probes of the cosmic ray ionization rate in diffuse clouds

TL;DR

This paper introduces a new method to estimate the cosmic-ray ionization rate in diffuse clouds using the abundance ratios of OH/CO and HCO+/CO, providing an alternative to traditional tracers that are difficult to detect.

Contribution

It proposes an analytical approach linking molecular abundance ratios to the cosmic-ray ionization rate, validated with ALMA observations, enhancing modeling of interstellar chemistry.

Findings

Abundance ratios correlate with cosmic-ray ionization rate $\

Derived CRIR values are higher than those at higher extinction, indicating CR attenuation.

Analytical expressions accurately predict OH abundance for $\

Abstract

The cosmic-ray ionization rate (CRIR, $\zeta_2$) is one of the key parameters controlling the formation and destruction of various molecules in molecular clouds. However, the current most commonly used CRIR tracers, such as H$_3^+$, OH$^+$, and H$_2$O$^+$, are hard to detect and require the presence of background massive stars for absorption measurements. In this work, we propose an alternative method to infer the CRIR in diffuse clouds using the abundance ratios of OH/CO and HCO$^+$/CO. We have analyzed the response of chemical abundances of CO, OH, and HCO$^+$ on various environmental parameters of the interstellar medium in diffuse clouds and found that their abundances are proportional to $\zeta_2$. Our analytic expressions give an excellent calculation of the abundance of OH for $\zeta_2$ $\leq$10$^{-15}$ s$^{-1}$, which are potentially useful for modelling chemistry in hydrodynamical simulations. The abundances of OH and HCO$^+$ were found to monotonically decrease with increasing density, while the CO abundance shows the opposite trend. With high-sensitivity absorption transitions of both CO (1--0) and (2--1) lines from ALMA, we have derived the H$_2$ number densities ($n_{\rm H_2}$) toward 4 line-of-sights (LOSs); assuming a kinetic temperature of $T_{\rm k}=50\,{\rm K}$, we find a range of (0.14$\pm$0.03--1.2$\pm$0.1)$\times$10$^2$ cm$^{-3}$}. By comparing the observed and modelled HCO$^+$/CO ratios, we find that $\zeta_2$ in our diffuse gas sample is in the { range of $1.0_{-1.0}^{+14.8}$ $\times$10$^{-16}- 2.5_{-2.4}^{+1.4}$ $\times$10$^{-15}$ s$^{-1}$. This is $\sim$2 times higher than the average value measured at higher extinction, supporting an attenuation of CRs as suggested by theoretical models.

Tables

Table 1. Table 1: Optical depths of CO (1–0) and (2–1); H 2 density derived from MCMC runs.

Sources	Velocity	${\tau }_{\mathrm{CO}}$			$n_{H_2}$/102 cm-3
	km s-1	1–0	2–1		$T_k$=10	50	100 K
	-10.46$\pm$0.04	0.055$\pm$0.002	0.027$\pm$0.002		2.5$\pm$0.7	0.7$\pm$0.2	0.5$\pm$0.1
3C454.3	-9.47$\pm$0.01	0.419$\pm$0.004	0.248$\pm$0.007		3.6$\pm$0.3	1.0$\pm$0.1	0.68$\pm$0.04
	-8.94$\pm$0.01	0.149$\pm$0.003	0.092$\pm$0.004		4.4$\pm$0.4	1.2$\pm$0.1	0.81$\pm$0.07
J0609-1542	7.37$\pm$0.02	0.12$\pm$0.02	0.032$\pm$0.007		0.3$\pm$0.3	0.2$\pm$0.1	0.18$\pm$0.09
J1733-1304	5.03$\pm$0.01	1.55$\pm$0.04	0.553$\pm$0.004		0.46$\pm$0.10	0.14$\pm$0.03	0.11$\pm$0.01
J1743-0350	3.89$\pm$0.05	0.17$\pm$0.01	0.06$\pm$0.02		0.49$\pm$0.53	0.26$\pm$0.18	0.21$\pm$0.12
	5.6$\pm$0.5	0.065$\pm$0.007	0.039$\pm$0.008		3.2$\pm$2.6	0.9$\pm$0.6	0.6$\pm$0.4

Table 2. Table 2: The E 𝐸 E (B-V) values, column densities of CO and HCO + , fit parameters from 3d-pdr models (see Section 6.3 ).

Sources	$E$(B-V)	$N_{\mathrm{CO}}$	$N_{{\mathrm{HCO}}^{+}}$	Fit parameters
	mag	1014 cm-2	1011 cm-2	$n_H$/102 cm-3	${\zeta }_2$/10-15 s-1
3C454.3 V1	0.106$\pm$0.002	0.84$\pm$0.03	1.41$\pm$0.03	${1.0}_{-0.7}^{1.0}$	${1.85}_{-1.66}^{0.30}$(*)
3C454.3 V2	0.106$\pm$0.002	3.35$\pm$0.03	3.54$\pm$0.09	${1.0}_{0.0}^{0.0}$	${1.58}_{-0.26}^{0.26}$(*)
3C454.3 V3	0.106$\pm$0.002	1.03$\pm$0.02	1.09$\pm$0.14	${1.3}_{-0.6}^{1.9}$	${0.12}_{-0.02}^{2.81}$(*)
J0609-1542	0.203$\pm$0.004	0.52$\pm$0.07	6.82$\pm$0.21	${0.1}_{0.0}^{0.0}$	${0.14}_{-0.02}^{0.05}$
J1733-1304	0.513$\pm$0.009	11.40$\pm$0.33	17.41$\pm$0.49	${0.3}_{0.0}^{0.0}$	${0.34}_{-0.06}^{0.06}$
J1743-0350 V1	0.530$\pm$0.008	2.79$\pm$0.19	7.34$\pm$1.23	${0.3}_{-0.2}^{0.4}$	${0.25}_{-0.18}^{0.38}$(*)
J1743-0350 V2	0.530$\pm$0.008	2.28$\pm$0.24	4.97$\pm$0.56	${0.5}_{-0.4}^{0.3}$	${0.25}_{-0.20}^{0.75}$(*)
3C120	0.265$\pm$0.006	1.14$\pm$0.24	2.71$\pm$0.58	${0.8}_{-0.7}^{0.8}$	${1.17}_{-1.14}^{0.42}$(*)
J1745-0753 V1	0.672$\pm$0.008	2.18$\pm$0.89	1.27$\pm$0.87	${3.2}_{-3.1}^{28.5}$	${0.86}_{-0.86}^{9.14}$(*)
J1745-0753 V2	0.672$\pm$0.008	8.06$\pm$1.26	5.60$\pm$2.46	${1.0}_{-0.9}^{2.2}$	${0.86}_{-0.81}^{3.12}$(*)
J1745-0753 V3	0.672$\pm$0.008	0.11$\pm$0.32	2.22$\pm$1.20	${0.2}_{-0.1}^{1.8}$	${0.10}_{-0.10}^{1.48}$(*)
J0211+1051	0.097$\pm$0.006	5.71$\pm$0.40	8.43$\pm$0.32	${0.5}_{0.0}^{0.0}$	${0.54}_{-0.09}^{0.09}$
J0325+2224	0.239$\pm$0.008	14.92$\pm$0.81	11.21$\pm$0.19	${0.3}_{0.0}^{0.0}$	${0.34}_{-0.06}^{0.06}$
J0356+2903	0.159$\pm$0.003	25.71$\pm$0.67	16.64$\pm$1.00	${31.6}_{-15.8}^{0.0}$	${0.34}_{-0.34}^{0.06}$
J0401+0413	0.230$\pm$0.009	3.02$\pm$0.76	5.44$\pm$0.23	${0.2}_{-0.1}^{0.6}$	${0.34}_{-0.26}^{0.52}$(*)
J0403+2600	0.141$\pm$0.002	9.84$\pm$1.06	6.66$\pm$0.32	${0.6}_{0.0}^{0.2}$	${0.74}_{-0.10}^{0.12}$(*)
J0406+0637	0.248$\pm$0.006	10.00$\pm$0.63	5.99$\pm$0.57	${0.8}_{-0.2}^{0.2}$	${1.00}_{-0.37}^{0.17}$(*)
J0407+0742	0.159$\pm$0.006	3.02$\pm$0.67	5.88$\pm$0.34	${0.4}_{-0.2}^{0.2}$	${0.63}_{-0.47}^{0.10}$(*)
J0426+2327	0.292$\pm$0.008	76.82$\pm$0.71	28.51$\pm$0.63	${31.6}_{-15.8}^{0.0}$	${0.18}_{-0.10}^{0.11}$
J0427+0457	0.292$\pm$0.017	6.19$\pm$0.71	6.88$\pm$0.27	${0.3}_{-0.2}^{0.3}$	${0.54}_{-0.38}^{0.19}$(*)
J0437+2037	0.451$\pm$0.006	10.63$\pm$0.41	17.09$\pm$0.81	${0.3}_{0.0}^{0.0}$	${0.34}_{-0.06}^{0.06}$
J0431+1731	0.354$\pm$0.004	11.11$\pm$0.78	11.21$\pm$1.22	${0.5}_{-0.1}^{0.0}$	${0.54}_{-0.14}^{0.09}$(*)
J0440+1437	0.442$\pm$0.005	13.17$\pm$0.46	13.42$\pm$0.34	${0.4}_{-0.1}^{0.0}$	${0.46}_{-0.12}^{0.07}$
J0449+1121	0.398$\pm$0.009	3.65$\pm$0.52	7.21$\pm$0.24	${0.2}_{-0.1}^{0.2}$	${0.29}_{-0.19}^{0.17}$(*)
J0502+1338	0.371$\pm$0.003	39.04$\pm$0.49	20.08$\pm$0.65	${31.6}_{0.0}^{0.0}$	${10.00}_{-1.61}^{1.61}$
J0510+1800	0.212$\pm$0.008	2.70$\pm$0.57	1.44$\pm$0.06	${0.8}_{-0.5}^{2.4}$	${2.51}_{-2.35}^{1.47}$(*)

Table 3. Table 3: The gas phase elemental abundances relative to hydrogen in the PDR models.

Elements	Abundance	Elements	Abundance
C+	1.4$\times$10-4	H2	3$\times$10-1
He	1$\times$10-1	H	4$\times$10-1
O	3$\times$10-4

Equations

$$\ln \mathrm{p}=-\frac{1}{2}\sum _{\mathrm{i}}\left[\frac{{\left({\tau }_{\mathrm{obs}}^{\mathrm{i}}-{\tau }_{\mathrm{model}}^{\mathrm{i}}\right)}^2}{{{\sigma }_{\mathrm{obs}}^{\mathrm{i}}}^2}+\ln \left(2\pi {{\sigma }_{\mathrm{obs}}^{\mathrm{i}}}^2\right)\right],$$

$$N_{tot}=\frac{3h}{8{\pi }^3{\left|{\mu }_{lu}\right|}^2}\frac{Q_{rot}}{g_u}\frac{e^{\frac{E_u}{k{T}_{ex}}}}{e^{\frac{h\nu }{k{T}_{ex}}}-1}\int {\tau }_{\nu }d\upsilon ,$$

$$Q_{tot}=\frac{k{T}_{\mathrm{ex}}}{h{B}_0}{e}^{\frac{h{B}_0}{3k{T}_{\mathrm{ex}}}},$$

$${\mathrm{H}}^{+}+\mathrm{O}\to {\mathrm{O}}^{+}+\mathrm{H},$$

$${\mathrm{O}}^{+}+{\mathrm{H}}_2\to \mathrm{O}{\mathrm{H}}^{+}+H.$$

$$\mathrm{O}+{\mathrm{H}}_3^{+}\to \mathrm{O}{\mathrm{H}}^{+}+{\mathrm{H}}_2.$$

$$\mathrm{O}{\mathrm{H}}^{+}+{\mathrm{H}}_2\to$$

$${\mathrm{H}}_2{\mathrm{O}}^{+}+{\mathrm{H}}_2\to$$

$${\mathrm{H}}_2{\mathrm{O}}^{+}+\mathrm{e}\to$$

$$\to$$

$${\mathrm{H}}_3{\mathrm{O}}^{+}+\mathrm{e}\to$$

$$\to$$

$${\mathrm{H}}_2\mathrm{O}+\mathrm{h}\nu \to \mathrm{OH}+H.$$

$${\mathrm{H}}_2\mathrm{O}+{\mathrm{C}}^{+}\to \mathrm{HC}{\mathrm{O}}^{+}+H.$$

$$\mathrm{OH}+\mathrm{h}\nu \to \mathrm{O}+\mathrm{H},$$

$$\mathrm{OH}+{\mathrm{C}}^{+}\to \mathrm{C}{\mathrm{O}}^{+}+\mathrm{H},$$

$$\mathrm{OH}+{\mathrm{C}}^{+}\to \mathrm{CO}+{\mathrm{H}}^{+},$$

$$\mathrm{OH}+{\mathrm{H}}^{+}\to \mathrm{O}{\mathrm{H}}^{+}+H.$$

$$x(\mathrm{OH})=$$

$$\times [1-\frac{n({\mathrm{H}}_2){\mathrm{k}}_{R5}\delta +n(\mathrm{e})({\mathrm{k}}_{R7}+{\mathrm{k}}_{R8})}{n({\mathrm{H}}_2){\mathrm{k}}_{R5}+n(\mathrm{e})({\mathrm{k}}_{R6}+{\mathrm{k}}_{R7}+{\mathrm{k}}_{R8})}]$$

$$=$$

$$x(\mathrm{OH})=$$

$$\mathrm{C}{\mathrm{O}}^{+}+{\mathrm{H}}_2\to$$

$${\mathrm{C}}^{+}+{\mathrm{H}}_2\mathrm{O}\to$$

$$\mathrm{CH}+\mathrm{O}\to$$

$$\mathrm{HC}{\mathrm{O}}^{+}+\mathrm{e}\to \mathrm{CO}+H.$$

$$x(\mathrm{HC}{\mathrm{O}}^{+})=$$

$$+\frac{x(\mathrm{O})}{x({\mathrm{C}}^{+})}x(\mathrm{CH}){\mathrm{k}}_{R20}].$$

$$x(\mathrm{HC}{\mathrm{O}}^{+})=$$

$$+\frac{x(\mathrm{O})}{x({\mathrm{C}}^{+})}\frac{x(\mathrm{CH})}{x(\mathrm{OH})}{\mathrm{k}}_{R20}],$$

$$ϵ=$$

$$\{{\mathrm{k}}_{\mathrm{pd}}(R12)+[{\mathrm{k}}_{\mathrm{pd}}(R12)+n(\mathrm{H})x({\mathrm{C}}^{+}){\mathrm{k}}_{R13}]\times$$

$$[\frac{{\mathrm{k}}_{R9}+{\mathrm{k}}_{R10}}{{\mathrm{k}}_{R11}}+\frac{{\mathrm{k}}_{R9}+{\mathrm{k}}_{R10}+{\mathrm{k}}_{R11}}{{\mathrm{k}}_{R11}}\frac{x(\mathrm{e}){\mathrm{k}}_{R6}}{x({\mathrm{H}}_2){\mathrm{k}}_{R5}}]{\}}^{-1}.$$

$$x(\mathrm{HC}{\mathrm{O}}^{+})\approx$$

$$\mathrm{CO}+{\mathrm{H}}_3^{+}\to \mathrm{HC}{\mathrm{O}}^{+}+{\mathrm{H}}_2.$$

Full text

Interstellar medium(847) — Interstellar molecules(849) — Chemical abundances(224) — Molecular clouds(1072)

††facilities: ALMA††software: CASA (McMullin et al., 2007), Astropy (Astropy Collaboration et al., 2013, 2018), radex (van der Tak et al., 2007), emcee (Foreman-Mackey et al., 2013), 3d-pdr (Bisbas et al., 2012)

1 Introduction

In the environments of the cold neutral component of the interstellar medium (ISM), where stellar photons cannot penetrate, the low-energy (0.1 $<$ $E$ $<$ 1 GeV) cosmic-rays (CRs) play a critical role in determining the ionization degree, in controlling the thermal balance and for initiating the chemistry at high optical depths (Dalgarno, 2006; Padovani et al., 2013; Vaupré et al., 2014; Grenier et al., 2015). The cosmic-ray ionization rate (CRIR; $\zeta_{2}$111Throughout the text, $\zeta_{2}$ denotes the CRIR of molecular hydrogen.) has been measured during the last few decades (Spitzer & Tomasko, 1968; Webber, 1998; Dalgarno, 2006) with the measured value to strongly depend on the adopted methodology (Dalgarno, 2006; Indriolo & McCall, 2012; Indriolo et al., 2015; Bacalla et al., 2019). Various observations found that the CRIR is higher toward the Galactic center (e.g., CMZ, $\zeta_{2}$ = 10*-15* $\sim$ 10*-14* s*-1*, Oka et al., 2005; Le Petit et al., 2016a) and supernova remnants (e.g., IC443, W49B, $\zeta_{2}$ $\sim$ 2$\times$10*-15* s*-1*, Indriolo et al., 2010; Zhou et al., 2022) than that of nearby molecular clouds (a few 10*-18* $\sim$ 10*-16* s*-1*, Caselli et al., 1998; Indriolo & McCall, 2012; Neufeld & Wolfire, 2017). Due to the interaction between CR particles and the ISM, the low-energy CRs attenuate while propagating at higher column densities (Strong et al., 2007; Padovani et al., 2018, 2020). In low-density diffuse clouds (e.g., $n_{\rm H}$ $\sim$ a few hundred cm*-2*), CRIR is -on average- almost an order of magnitude higher than that in dense clouds (Indriolo & McCall, 2012).

While $\zeta_{2}$ cannot be easily observed directly, the use of tracers is favored instead. H${}^{+}_{3}$ is one of the most commonly used tracers of the CRIR. It is produced by the CR ionization of H2 and it is destroyed through reactions with abundant neutral species (e.g., CO, O) and electrons (McCall et al., 1999; Geballe et al., 1999; Dalgarno, 2006; Indriolo & McCall, 2012). The CRIR can be derived once the gas temperature, volume density ($n_{\rm H}$), and column densities (e.g. H2, CO) are known. Other molecules that are directly relevant to the H${}^{+}_{3}$ chemistry are considered as potential probes of the CRIR, such as HCO+ and DCO+ (Guelin et al., 1982; van der Tak & van Dishoeck, 2000; Caselli et al., 1998), OH+ (Indriolo et al., 2015; Bacalla et al., 2019), and H2O+ (Gerin et al., 2010; Neufeld & Wolfire, 2017; Bialy et al., 2019). Measuring $\zeta_{2}$ with ions such as H${}^{+}_{3}$, OH+, and H2O+, requires the presence of bright massive stars in the background, which is not very common. Furthermore, the deuterium species can only be detected in high extinction regions due to its relatively low abundance. Although the above methodology can provide reasonable measurements, it is difficult and somehow impossible for general use.

However, oxygen-bearing molecules have -in principle- the potential to constrain the CRIR; most of the formation of oxygen-bearing species (e.g., OH, HCO+) starts from the hydrogenation of ionized oxygen (O+). Since the ionization potential of atomic oxygen (13.62 eV) is very close to that of atomic hydrogen (13.6 eV), the majority of oxygen in the cold neutral medium (CNM) is in its atomic form. Thus, in such FUV-shielded regions, O+ and oxygen-bearing species are the indirect products of CRs.

OH has long been proposed to be an alternative tracer of molecular gas due to its fairly constant abundance in the ISM (Liszt & Lucas, 1998, 2002; Xu & Li, 2016; Li et al., 2018). The thermal emission and absorption lines of OH 18 cm have been detected extensively toward the CNM, which is extended to the outskirts of the molecular clouds and where the CO emission is faint or undetectable (Turner, 1979; Magnani et al., 1988; Wannier et al., 1993; Cotten et al., 2012; Li et al., 2018; Busch et al., 2021). HCO+, in the traditional sense, is a dense gas tracer in the ISM. However, recent observations found that HCO+ is ubiquitous in diffuse and translucent clouds (Pety et al., 2017; Luo et al., 2020), especially with absorption measurements against strong continuum sources (e.g., quasars, H ii regions).

In this paper, we combine absorption observations of HCO+, and absorption and emission observations of CO to calculate their column densities along each LOS. We attempt to find the chemical connection between these oxygen-bearing molecules (CO, OH, and HCO+) that are ubiquitous in diffuse clouds. We investigate the variance of chemical abundances of oxygen-bearing molecules under different environmental parameters (e.g., gas volume density, FUV intensity, CRIR), especially the potential connection of their molecular abundance ratios as a probe of CRIR.

This paper is organized as follows. The observations and archival data used in this work are presented in Section 2. We present the results of column densities and constraints on the gas density in Section 3. In Section 4, we perform photodissociation region (PDR) modelling of CO, OH, and HCO+ under various environmental parameters. We discuss the chemistry of CO, OH, and HCO+ in the diffuse cloud and derive the chemical abundances in chemical equilibrium in Section 5. We derive the abundance ratios of OH/CO and HCO+/CO and constrain the CRIR by combining the observed HCO+-to-CO abundance ratio and chemical models in the low-density ($n_{\mathrm{H}}$ $\sim$ 102 cm*-3*) diffuse cloud in Section 6. The main results and conclusions are summarized in Section 7.

2 Observations and Data

2.1 HCO+ (1–0) and CO (1–0)

The HCO+ (1–0) and CO (1–0) absorption observations toward 13 strong continuum sources were carried out during Apr. 2016 to May. 2016 with ALMA (project ID: 2015.1.00503.S, PI: L. Bronfman). The calibration of the raw data was performed using the Common Astronomy Software Applications (McMullin et al., 2007). Self-calibration was performed toward 2 sources (3C454.3 and 3C120) to increase the signal-to-noise ratio and eliminate the spectral contamination from bandpass calibrators. For each line-of-sight (LOS), we decompose the absorption spectra ($\tau_{\nu}$) into different Gaussian components to derive the column density at each velocity component. A detailed description of observations and data reduction can be found in Luo et al. (2020).

The HCO+ (1–0) integrated optical depth (in units of km s*-1*) and CO (1–0) emission toward another 15 sources are taken from Table 1 in Liszt & Gerin (2023). The original data of HCO+ absorption was presented in Lucas & Liszt (1996), Liszt & Lucas (2000), and Liszt & Gerin (2018). We exclude sources in which the brightness temperature of CO is bright ($T_{\rm mb}$ $\geq$ 5 K) since we only focus on diffuse LOSs.

2.2 CO (2–1) data

The CO (2–1) absorption spectra toward 4 sources (3C454.3, 0607-157, 1730-130, and 1741-038) are taken from a blind survey of the absorption lines of bandpass calibrators in ALMA regular observations (Luo et al. 2023, in prep). The data reduction follows the same procedures as that of CO (1–0).

2.3 Reddening $E$(B-V) data

The $E$(B-V) data is taken from Green et al. (2019), in which the values were derived by combining stellar photometry from Pan-STARRS 1 and 2MASS, and parallaxes from Gaia. We use the $E$(B-V) values to derive the total gas column density along the LOS.

3 Results

3.1 Constraints on the gas densities

For those components where both CO (1–0) and (2–1) transitions are available, we use the non-LTE radiative transfer code radex (van der Tak et al., 2007) to constrain the gas densities along the LOSs. We account for H2 as the main colliding partner of CO. We perform several models covering H2 volume densities of 101 $\leq$ $n_{\rm H_{2}}$ $\leq$ 104 cm*-3*, CO column densities of 1011 $\leq$ $N_{\rm CO}$ $\leq$ 1016 cm*-2* and kinetic temperatures of $10^{1}\leq T_{\rm k}\leq 10^{2}\,{\rm K}$. We then find the optimum solutions by maximizing the likelihood function using the Markov chain Monte Carlo (MCMC) method, which is encoded in $emcee$ (Foreman-Mackey et al., 2013). The likelihood function is defined as:

[TABLE]

where the ${\rm\tau_{obs}^{i}}$ and ${\rm\sigma^{i}_{obs}}$ are the optical depth and uncertainty of the observed i-th transition, respectively. ${\rm\tau_{model}^{i}}$ is the modelled optical depth by radex.

The representative optimum results of $n_{\rm H_{2}}$ under different $T_{\rm k}$ are shown in the left panel of Fig. 1. As $T_{\rm k}$ increases by a factor of 10, the resultant $n_{\rm H_{2}}$ decreases by a factor of 2$\sim$5. However, the column densities are less influenced by $T_{\rm k}$ and vary by $<$4% (for CO) and $<$0.1% (for HCO+) (with respect to the values at $T_{\rm k}$=50 K, as shown in the middle and right panels of Fig. 1.)

Table 1 summarises the derived $n_{\rm H_{2}}$ and $n_{\rm H}$ of each velocity component at $T_{\rm k}$ = 10, 50, and 100 K. The derived $n_{\rm H_{2}}$ is in the range of (0.11$\pm$0.01–4.4$\pm$0.4)$\times$102 cm*-3*.

3.2 Calculation of column densities

For the rest of the components where only CO (1–0) and HCO+ (1–0) are available, we consider an excitation temperature of $T_{\rm ex}=4\,{\rm K}$ for CO and $T_{\rm ex}=2.73\,{\rm K}$ for HCO+ (e.g. Liszt & Lucas, 1996; Godard et al., 2010; Luo et al., 2020).

The column densities are calculated with (Mangum & Shirley, 2015):

[TABLE]

where $\mathrm{\left|\mu_{lu}\right|^{2}}$ is the dipole matrix element, $\mathrm{Q_{rot}}$ is the rotational partition function, $\mathrm{g_{u}}$ is the degeneracy of the upper energy level, and $\mathrm{E_{u}}$ is the energy of the upper energy level. For each transition, the $\mathrm{\left|\mu_{lu}\right|^{2}}$, $\mathrm{E_{u}}$, and the rest frequency $\nu$ are taken from the CDMS database (Müller et al., 2001, 2005). For linear molecules, the partition function is given by (McDowell, 1987):

[TABLE]

where $B_{0}$ is the rigid rotor rotation constant. The calculated column densities of CO and HCO+ are summarised in Table 2.

4 Photodissociation region modelling

Variations in chemical abundances can be used as a diagnostic tool for estimating the environmental parameters of the ISM. To better understand the response of the abundances of the molecular species we examine (CO, OH, HCO+) and which are ubiquitously detected in diffuse and translucent clouds, we perform chemical simulations under a range of different environmental parameters (e.g. varying densities, FUV intensities, and cosmic-ray ionization rates). For these calculations, we use the publicly available code 3d-pdr222https://uclchem.github.io/3dpdr/ (Bisbas et al., 2012).

In our modelling, we use a suite of one-dimensional uniform slabs with densities of $n_{\rm H}=10^{m}\,{\rm cm}^{-3}$ where $m=2.0,2.5,3.0,3.5$, interacting with four FUV intensities $\chi$/$\chi_{0}$=1,10,30,100 (where $\chi$ is normalized to the spectrum of Draine 1978), and four CRIR $\zeta_{2}=(5,10,50,100)\times$10${}^{-17}\,{\rm s}^{-1}$. The FUV radiation field is considered as plane-parallel that impinges from one direction. The diffuse component of radiation is not accounted for. The maximum value of $A_{\mathrm{V}}$ in our simulations is 20 mag. The gas-phase element abundances relative to hydrogen are shown in Table 3. We set all of the element carbon in the ionized phase with an abundance of 1.4$\times$10*-4* and 60% of the total hydrogen in the molecular phase as the initial conditions in our PDR models (see Röllig et al., 2007, for further details). Throughout the text, the diffuse cloud is referred to as molecular gas with $n_{\rm H}$ $\leq$500 cm*-3* and $A_{\rm V}$ $\leq$1 mag, the translucent cloud is referred to 500$\leq$ $n_{\rm H}$ $\leq$5000 cm*-3* and $A_{\rm V}$ $\leq$5 mag.

Figure 2 shows our simulation results for a fixed $\zeta_{2}=10^{-16}\,{\rm s}^{-1}$, in which the relative abundance of each species to H2 is defined as the ratio of the corresponding column densities along the LOS. The variance of CO abundance is shown in the first row of Fig. 2. CO is mostly photodissociated, especially at lower extinction regions (e.g., $A_{\rm V}$ $<$1 mag). As a consequence, the abundance of CO at low $A_{\rm V}$ decreases by two orders of magnitude (from $\sim$ 10*-7* to $\sim$ 10*-9*) as the FUV intensity increases from $\chi$/$\chi_{0}$=1 to 100. The abundance of CO increases with increasing $A_{\rm V}$ (e.g., 0.6$\leq$ $A_{\rm V}$ $\leq$5 mag) due to FUV shielding by dust and also self-shielding. However, the abundance of CO slightly decreases by a factor of a few from low to intermediate extinctions ($A_{\rm V}$ $\leq$0.6 mag), which is due to the decrease of the abundance of its precursors (see Section 5). For a fixed FUV intensity, the abundance of CO increases with increasing density. The divergence of CO abundance at different densities is smaller (within a factor of a few, from 5$\times$10*-8* to 5$\times$10*-7* for $\chi$/$\chi_{0}$ = 1) at low extinction ($A_{\rm V}$ $<$1 mag) but larger, over an order of magnitude, at high $A_{\rm V}$. At high-density or low FUV intensity models, the abundance of CO reaches the canonical value (e.g., $\sim$10*-4* as in the solar neighborhood) at lower extinction ($A_{\rm V}$ = 2 mag for $\chi$/$\chi_{0}$ = 1 and $n_{\rm H}$ = 103.5 cm*-3*) than that of low-density or high FUV intensity models ($A_{\rm V}$ = 5 mag for $\chi$/$\chi_{0}$ = 10 and $n_{\rm H}$ = 103.5 cm*-3*).

On the contrary, the abundance of OH has a different pattern to that of CO, as shown in the second row of Fig. 2. At low $A_{\rm V}$ (e.g., $\lesssim$2 mag), the abundance of OH increases with the increasing FUV intensity, especially for high-density simulations (e.g., $n_{\rm H}$ = 103.5 cm*-3*). It increases from $\sim$10*-9* to 10*-8* as the FUV intensity increases from $\chi$/$\chi_{0}$=1 to 100. The abundance of OH decreases with the increasing $A_{\rm V}$ when $A_{\rm V}$ $<$2 mag. The slope of the decreasing trend of OH abundance with $A_{\rm V}$ becomes larger for lower densities and higher FUV intensities. As can be seen from Fig. 2, the abundance of OH decreases by a factor of 2 for the high-density ($n_{\rm H}$ = 103.5 cm*-3*) and low FUV intensity ($\chi$/$\chi_{0}$=1) model, while it decreases by two orders of magnitude for the low-density ($n_{\rm H}$ = 102 cm*-3*) and high FUV intensity ($\chi$/$\chi_{0}$=100) model.

The knee point where the OH abundance trend turns, is approximately at an $A_{\rm V}$ of 2 mag, after which it monotonically decreases with the increase of density below the $A_{\rm V}$ value, while it converges above the $A_{\rm V}$ value (except for $n_{\rm H}$ = 103.5 cm*-3*). Particularly, the divergence of OH abundance at low FUV intensities is larger than that of at high FUV intensity. With the density increase from 102 cm*-3* to 103 cm*-3*, the abundance of OH decreases by nearly two orders of magnitude at $\chi$/$\chi_{0}$=1, while it decreases only less than one order of magnitude at $\chi$/$\chi_{0}$=100.

Furthermore, the abundance of HCO+ follows a similar trend to that of OH at low FUV intensity; HCO+ decreases with the increase of density and $A_{\rm V}$ ($A_{\rm V}$ $<$2 mag). While this trend is maintained for the high-density (e.g., $n_{\rm H}$ $>$103 cm*-3*) and low FUV intensity models (e.g., $\chi$/$\chi_{0}$ $<$100), the abundance of HCO+ declines at low extinction for low-density and high FUV intensity models. At high FUV intensities, HCO+ does not form as efficiently as OH does. The difference between HCO+ and OH is most likely due to the interruption of the formation of HCO+ through the OH channel at high FUV intensities (see Section 5).

Figure 3 shows our simulations with a constant FUV intensity ($\chi$/$\chi_{0}$=1). The abundances of CO, OH, and HCO+ have a similar trend as $\zeta_{2}$ varies: they increase with increasing CRIR at low extinction (e.g., $A_{\rm V}$ $<$ 1 mag for CO at $n_{\rm H}$=103 cm*-3*), while they have an inverse relation at high extinctions. The abundance of CO increases by an order of magnitude as the CRIR increase from 10*-16* to 10*-15* s*-1* at low extinctions for low-density models (e.g., $n_{\rm H}$ $\leq$103 cm*-3*). The $A_{\rm V}$ value of the “inverse point” where the dependence of CO abundance on the CRIR turns opposite, is higher ($A_{\rm V}$ $\sim$5 mag) for low-density models and lower ($A_{\rm V}$ $\sim$0.3 mag) for high-density models. For the $n_{\rm H}$ = 103.5 cm*-3* models, an order of magnitude increase in CRIR (from $\zeta_{2}=10^{-16}\,{\rm s}^{-1}$ to $10^{-15}\,{\rm s}^{-1}$) affects the CO abundance from a factor of $\sim 2$ at low $A_{\rm V}$ to $\gtrsim 10^{2}$ at higher $A_{\rm V}$.

As an important ISM environmental parameter, the CRIR has a high impact on the abundance of OH. In particular, it is almost proportional to $\zeta_{2}$ at low $A_{\rm V}$. As seen in the second row of Fig. 3, the abundance of OH increases by an order of magnitude as $\zeta_{2}$ increases from 10*-16* to 10*-15* s*-1* at low $A_{\rm V}$ (e.g., $A_{\rm V}$ $<$1 mag) for all four densities explored. The $A_{\rm V}$ value of the “inverse point” is larger than that of CO, which decreases with increasing density. At high extinction ($A_{\rm V}$ $>$2 mag) and high-densities ($n_{\rm H}$ $>$103 cm*-3*), the abundance of OH decreases with the increasing $\zeta_{2}$.

The dependence of HCO+ abundance on the CRIR is similar to that of OH. With the CRIR increase by an order of magnitude, the abundance of HCO+ increases by an order of magnitude at low $A_{\rm V}$ for $n_{\rm H}$ $\leq$103 cm*-3*. However, at low-density models ($n_{\rm H}$ = 102 cm*-3*), an increase of CRIR from 5$\times$10*-16* s*-1* to 10*-15* s*-1* only results in a slight increase ($\leq$30%) on the abundance of HCO+, meaning that the abundance of HCO+ would hit a maximum value in this case (see Section 5.2 for more discussion). At high-densities ($n_{\rm H}\sim$103.5 cm*-3*), the increase of CRIR by an order of magnitude increases HCO+ abundance by a factor of $\sim$5. The $A_{\rm V}$ value of the “inverse point” is the same as that of OH.

5 Analysis and Discussion

5.1 The abundance of OH in diffuse cloud

In diffuse clouds, the formation of OH can be traced back from two channels that form OH+. The first one starts from the charge transfer reaction:

[TABLE]

with O+ forming OH+ through reaction:

[TABLE]

In addition to the reactions R1 and R2, atomic oxygen can directly react with H${}_{3}^{+}$ to form OH+ (the second channel):

[TABLE]

Once OH+ has formed, it can hydrogenate to form the precursors of OH, namely H2O+ and H3O+:

[TABLE]

in which reaction R4 is the main destruction process for OH+. OH is formed through electron recombination reactions:

[TABLE]

and

[TABLE]

From reaction R11, H2O also contributes to the formation of OH by photodissociation:

[TABLE]

Reaction R12 is also the main destruction path of H2O. The minor destruction path of H2O is through C+ in diffuse clouds:

[TABLE]

The destruction of OH in diffuse clouds includes photodissociation:

[TABLE]

reaction with C+:

[TABLE]

and reaction with H+:

[TABLE]

The OH abundance ($x$(OH)=$n$(OH)/$n_{\rm H}$) in chemical equilibrium can be solved when the formation and destruction processes reach a balance (see detailed derivation in Appendix A). In diffuse clouds, photodissociation dominates the destruction of OH, and reactions R15 and R16 play a minor role. Since the abundance of C+ is an order of magnitude higher than H+ when the $\zeta_{2}$ is not high (e.g., $\zeta_{2}$ $\lesssim$10*-15* s*-1*, Le Petit et al., 2016b), we ignore here the contribution from R17 (the last term in equation A17). Thus, the abundance of OH can be written333In the following, $\rm k_{Rx}$ represents the reaction rate of Reaction Rx, and $\rm k_{pd}(Rx)$ represents the photodissociation rate of Reaction Rx. as:

[TABLE]

where $\phi$ is the term in the square brackets of equation 4 and $\theta$ is in the form of equation A15. Since the term $\delta$ is small ($\sim$0.1) in average diffuse ISM conditions (e.g., T =30 K, $A_{\rm V}$ = 1 mag), $\phi$ $\approx$1. To a good approximation, equation 4 can be simplified as:

[TABLE]

Since there is an anti-correlation between x(OH) and $n_{\mathrm{H}}$ as density increases, the abundance of OH decreases monotonically ( similar to that shown in Fig. 2). On the other hand, if the UV radiation field increase or the extinction decreases, the gas temperature will increase. Therefore, all reaction rates in equation 5 will decrease and, therefore, the abundance of OH will also increase.

As can be seen from equation 5, the abundance of OH is proportional to $\zeta_{2}$ in the diffuse cloud. Figure 4 shows the resultant abundance of OH ($A_{\rm V}$ = 1 mag, $n_{\rm H}$ = 102 cm*-2*) from both equations 5 as well as the 3d-pdr models, in which the reaction rates are taken from McElroy et al. (2013). It should be emphasized that equation 5 can only be used under average ISM conditions. At high FUV radiation fields or very low $A_{\rm V}$, the main formation pathway of OH is contributed by neutral-neutral reactions through atomic O and H, which are not included in our analysis444The reaction rate of O and H2 is $<<$10*-30* cm3 s*-1* at a gas temperature below 100 K, which can be ignored. However, once the temperature increases above 200 K, the formation rate increases to 10*-20* cm3 s*-1*, which is comparable to that in equation A16.. Equation 5 represents a good approximation ( with $<$60% uncertainty) for $\zeta_{2}$ $<$10*-15* s*-1* when compared with 3d-pdr models.

However, in the PDR simulations, the abundance of OH reaches the maximum at $\zeta_{2}$ $\sim$ 2$\times$10*-15* s*-1*. For higher $\zeta_{2}$, the OH abundance decreases. This is because the destruction of OH at high CRIR is dominated by H+ (Reaction R17) rather than by photodissociation and C+ (Reactions R14$\sim$R16), which has been ignored in deriving equation 5. Similar work using isothermal simulations has been reported by Bisbas et al. (2017, their Fig. 10), in which the OH abundance peaks at $\zeta_{2}$/$n_{\rm H}$ $\approx$ 2$\times$10*-17* cm3 s*-1*.

Though there are limitations, we still highlight the usage of equation 5 at $\zeta_{2}$ $<$10*-15* s*-1* instead of running chemical networks, especially when coupled with hydrodynamical simulations, as this dramatically reduces the overall computational expense.

5.2 The abundance of HCO+ in diffuse cloud

The formation of HCO+ starts from ion-neutral reactions:

[TABLE]

where CO+ is the result of reaction R15. Reaction R18 is very efficient as almost all CO+ forms HCO+555This assumption may not be true if the UV radiation field is much higher than the usual ISM and the molecular fraction is small. In such a case, the destruction of CO+ would involve atomic hydrogen and free electrons instead of H2..

The destruction of HCO+ is always dominated by free electrons:

[TABLE]

Thus, in chemical equilibrium, the HCO+ abundance ($\rm{\it x}(HCO^{+})$) can be written as:

[TABLE]

By substituting equations A1-A17 to the above and relating $\rm{\it x}(H_{2}O)$ with $\rm{\it x}(OH)$, we obtain:

[TABLE]

where $\epsilon$ is:

[TABLE]

For low values of $\zeta_{2}$ (e.g., $\zeta_{2}$ $\lesssim$ 10*-15* s*-1*), the production of electrons results primarily from the ionization of atomic carbon. The electron abundance can be, thus, approximated with the C+ abundance, e.g. $\rm{\it x}(e)$ = $\rm{\it x}(C^{+})$ (Goldsmith, 2001). If the ionization degree of the cloud is large or the gas density is low, x(O)/x(C+) $\approx$ 2, and the abundance of CH is an order of magnitude lower than that of OH. Thus, the contribution from the last term of equation 7 is small. Equation 7 can, then, be simplified as:

[TABLE]

For higher $\zeta_{2}$ (e.g., $\zeta_{2}$ $\gg$ 10*-15* s*-1*), a significant fraction of electrons are produced by H+ (Le Petit et al., 2016b). In this case, equation 7 should be applied.

However, we should still keep in mind that at low extinctions and high FUV intensities, the destruction of CO+ by H and free electrons becomes more efficient than that of H2 (Reaction R18). The formation of HCO+ is, therefore, interrupted. The response of HCO+ abundance does not follow a similar trend to that of OH at low $A_{\rm V}$ and high FUV intensities (as shown in Fig. 2). In such a case, equation 9 is not valid.

Overall, equation 9 shows that $x$(HCO+) is proportional to x(OH) in diffuse clouds. High-sensitivity absorption observations found that the integrated optical depth of HCO+ has a tight linear relation with that of OH ($N$(HCO+)/$N$(OH) = 0.03, Liszt & Lucas, 1996). Considering a diffuse cloud with gas temperature ranging from 50 to 100 K, at an extinction range of $0.5-1.5\,{\rm mag}$, the resultant HCO+/OH ratio is nearly constant ($\sim$3.5$\times$10*-3*, see Fig. 5). Note that while this value is an order of magnitude lower than the above observations, the abundance of OH measured in diffuse gas has large uncertainty between different LOSs (see Section 5.4 for more discussion). Additionally, most of the HCO+ observations of Liszt & Lucas (1996) have high optical depths ($\tau>$1). This means that the formation and destruction of HCO+ may be different from what we considered here. Furthermore, the CH abundance in translucent clouds is likely to be comparable to that of OH ($\sim$4$\times$10*-8*, Liszt & Lucas, 2002; Sheffer et al., 2008), meaning that the HCO+/OH ratio could be underestimated by equation 9.

Finally, in the case of a high CO abundance, the latter becomes the precursor of HCO+ through reaction:

[TABLE]

In this case, the abundance of HCO+ may be underestimated.

5.3 The abundance of CO in diffuse cloud

There are three major formation paths of CO: Reactions R16, R21, and

[TABLE]

The destruction of CO in the diffuse and translucent clouds is dominated by photodissociation:

[TABLE]

and by the interaction with He+

[TABLE]

The CO abundance ($\rm{\it x}(CO)$) in chemical equilibrium can be written as:

[TABLE]

$\rm\eta_{CO}$ denotes the destruction rate of CO by He+, which is proportional to the CRIR. If we substitute x(HCO+) as x(OH) using equation 7, we obtain:

[TABLE]

The CO abundance is proportional to the gas density. In low-density diffuse gas where photodissociation dominates the destruction of CO, the OH channels dominate the formation of CO (Luo et al., 2023). Equation 11 can be simplified as :

[TABLE]

Since $\rm{\it x}(OH)$ is proportional to the CRIR, $\rm{\it x}(CO)$ increases with the increasing CRIR.

5.4 Comparison between the observed molecular abundances and model predictions

Radio emission line observations toward high latitude diffuse clouds (0.4$\lesssim$ $A_{\rm V}$ $\lesssim$1.1 mag) have reported an abundance of OH between 1.6$\times$10*-7* and 4$\times$10*-6* (Magnani et al., 1988). Absorption measurements toward quasars at radio wavelengths suggest a relative abundance ratio to total H column density ($N_{\rm OH}$/$N_{\rm H}$) of $2.5-5\times 10^{-8}$ at $A_{\rm V}$ of 1 mag (Crutcher, 1979; Liszt & Lucas, 1996). Considering that the molecular fraction is 10%-20% at such LOSs (Lucas & Liszt, 1996), the abundance of OH with respect to H2 is in the range of $\sim$10*-7* to 10*-6*. Recent radio observations by Tang et al. (2021) covering a broad range of $A_{\rm V}$ (0.2$-$60 mag) suggest that the abundance of OH is higher ($\sim$10*-6*) at low $A_{\rm V}$ and lower ($\sim$10*-7*) at higher $A_{\rm V}$. The high spatial resolution ($\sim$0.12 pc) observations by Xu et al. (2016) in the Taurus boundary also found a decreasing trend of OH abundance with the $A_{\rm V}$.

Our models show that the abundance of OH strongly depends on the density distribution and the ISM environmental parameters (see Section 4, Fig. 2 and 3) and can vary by more than two orders of magnitude ($\sim$ 10*-9* $-$ 10*-7* at the same $A_{\rm V}$). As can be seen from equation 5, x(OH) is anti-correlated with ${\rm k_{pd}(R22)}$ and $n_{\rm H}$. This is consistent with the OH survey in nearby molecular clouds of Tang et al. (2021). Considering the aforementioned large scatter of OH abundance in diffuse clouds, we adopt an OH abundance in the range of 10*-7* to 10*-6* as the “typical value”. Thus, for any given FUV intensity (1$\leq$ $\chi$/$\chi_{0}$ $\leq$100) and density (102 $<$ $n_{\rm H}$ $\leq$ 103.5 cm*-3*), the model underestimates the abundance of OH at an $A_{\rm V}\simeq 0.2-2\,{\rm mag}$ if a $\zeta_{2}=10^{-16}\,{\rm s}^{-1}$ is adopted (Fig. 2). As seen from Fig. 3, the CRIR should be no less than 10*-16* s*-1* if we are to reproduce the observed abundance of OH in diffuse clouds.

Due to the sub-thermal excitation of HCO+ transitions in low-density gas (Godard et al., 2010; Luo et al., 2020), its abundance can vary over an order of magnitude without precise measurement of both the excitation temperature and optical depth (Luo et al., 2023). The abundance of HCO+ in diffuse clouds has been measured frequently through absorption observations against strong continuum sources (e.g., quasars, H ii regions). Despite the different methods in obtaining the column density of H2, the abundance of HCO+ is fairly constant in diffuse gas ((1.7–3.1)$\times$10*-9*, Lucas & Liszt, 1996; Liszt & Gerin, 2016; Gerin et al., 2019; Luo et al., 2020). At a CRIR of 10*-16* s*-1*, our models underestimate the abundance of HCO+ in diffuse gas in all density and FUV intensity ranges explored (Fig. 2). Therefore, to reproduce the observed abundance of HCO+, the gas density should be approximately in the range of 102 $\leq$ $n_{\rm H}$ $\leq$ 102.5 cm*-3* and the CRIR should be $\zeta_{2}>10^{-16}\,{\rm s}^{-1}$ (Fig. 3). The inferred density is consistent with the inferred density in Section 3.1 and those by radio and UV absorptions (e.g., 80$-$160 cm*-3*, Goldsmith, 2013; Liszt & Gerin, 2016; Luo et al., 2020).

In diffuse clouds, most of the carbon is in the form of C+ or C0 due to insufficient shielding from UV photons. Similar to HCO+, the CO low-$J$ transitions are usually sub-thermally excited in the low-density diffuse gas (Goldsmith et al., 2008; Luo et al., 2020). The abundance of CO in diffuse clouds can vary by two orders of magnitude in different environments (e.g., 2.6$\times$10*-8*–2$\times$10*-5* from UV absorption measurements). In particular, it increases with increasing $N_{\rm H_{2}}$ (or $A_{\rm V}$) (Burgh et al., 2007; Liszt, 2007; Sheffer et al., 2008). The abundance of CO through CO (J=1–0) emission line measurements in low extinction regions ($A_{\rm V}\sim 1\,{\rm mag}$) in Taurus666The value is obtained by averaging the pixels without CO detection. is approximately (1.2–7)$\times$10*-6* (Goldsmith et al., 2008; Pineda et al., 2010), which is over an order of magnitude lower than the canonical value in well-shielded regions (e.g., $\sim$10*-4*, Frerking et al., 1982). This value is similar to that of absorption observations in diffuse LOSs ($A_{\rm V}$ = 0.19–2.08 mag, $f_{\rm CO}$ = (0.2$\pm$0.1–5$\pm$4)$\times$10*-6*, Luo et al., 2020). Our models shown in Fig. 3 can reasonably explain the observed CO abundances.

Considering all the above, we find that the abundances of OH, HCO+, and CO in diffuse clouds suggest an ISM environment with low gas densities ($n_{\rm H}$ $\lesssim$102.5 cm*-3*) and high CRIRs ($\zeta_{2}$ $\gtrsim$10*-16* s*-1*). A more quantitative analysis is discussed below (§6.3).

6 The abundance ratios and CRIR in diffuse clouds

As can be seen from equations 5 and 9, the abundances of OH and HCO+ are proportional to the CRIR in diffuse gas, which implies they can be used to constrain CRIR with a given density. However, molecular hydrogen does not emit radiation that can be observed from radio telescopes due to the lack of permanent dipole moment, the measurement of an accurate H2 column density in the CNM, as well as its exact abundance, is very hard. Instead, it is possible to put constraints on the CRIR by combining the observed molecular column densities and their ratios with chemical models.

6.1 The abundance ratio of OH/CO

The abundance ratio of OH/CO can be obtained from equation 11:

[TABLE]

The second term in the square brackets can be safely ignored in the diffuse gas777The term ${\rm k_{pd}(R24)}$+${\rm\eta_{CO}}$ is comparable to ${\rm k_{R20}+k_{R23}}$ at $A_{\rm V}$ = 1 mag, while 1/${\rm{\it n}(H)}$ is apparently a few orders of magnitude higher than ${\rm\frac{{\it x}(CH)}{{\it x}(CO)}{\it x}(O)}$, thus, the second term can be ignored.:

[TABLE]

Figure 6 shows the predicted abundance ratio of OH/CO with the column density of CO ($N_{\rm CO}$) from 1D slab model simulations ($n_{\rm H}$ = 102 cm*-3*). The abundance ratio of OH/CO decreases as the increasing $N_{\rm CO}$. The abundance ratio of OH/CO will increase with the increasing FUV intensity by a factor of a few at low $N_{\rm CO}$, while it remains approximately constant at high $N_{\rm CO}$. This is because, at a low FUV intensity (e.g., $\chi$/$\chi_{0}\sim 1$), the gas temperature is significantly lower than that of a moderate FUV intensity (e.g., $\chi$/$\chi_{0}\sim 10$). The decreasing gas temperature would increase the reaction rates of R15, R16, and R19, leading to a lower OH/CO ratio. However, at high $N_{\rm CO}$, CO molecules exist mainly in well-shielded regions, in which the gas temperature does not significantly increase even when high external FUV intensities exist.

As shown in Fig. 6, the abundance ratio of OH/CO monotonically increases with the increasing CRIR at a given $N_{\rm CO}$. Since high CRIR heats the gas, the reaction rates of R15, R16, and R19 will decrease. At the meantime, $\eta_{\rm CO}$ will increase as the increasing CRIR, leading to a higher OH/CO ratio.

The results shown in Fig. 6 indicate that once the gas density can be constrained from line ratios (e.g., C i/CO, CO rotational line ladder), the CRIR can be inferred from OH and CO observations without knowledge of their exact abundance relative to H2. This is potentially useful in high spatial resolution radio spectral line observations, especially when the column densities of H2 cannot be determined.

However, the main challenge of generalizing the use of OH and CO as a probe of CRIR is that the excitation temperature of OH is usually within a few K above the Galactic synchrotron background (Li et al., 2018). This will lead to an order of magnitude uncertainty in $N_{\rm OH}$ through emission lines if $T_{\rm ex}$ varies between 0.1–1.0 K above $T_{\rm bg}$ (note that $N_{\rm OH}$ $\propto$ $T_{\rm ex}$/($T_{\rm ex}$-$T_{\rm bg}$)). Absorption measurements are less suffered from the above issues ($N_{\rm OH}$ $\propto$ $T_{\rm ex}$), while the optical depth of OH is usually between a few to tens of 10*-3* (more than an order of magnitude lower than that of HCO+). This affects the feasibility of detecting the absorption of OH in a relatively short integration time even with the JVLA. With the proposed capability of high-sensitivity radio telescopes such as FAST, and SKA in the future (McClure-Griffiths et al., 2015), it is possible to constrain the CRIR through OH and CO in both the Milky Way and external galaxies.

6.2 The abundance ratio of HCO+/CO

Contrary to OH, HCO+ can be easily detected in diffuse clouds through absorption measurements against strong continuum sources (e.g., a few to tens of minutes with ALMA), and has less uncertainty in calculating its column density.

Following the same way, if we replace equation 9 with equation 14, the abundance ratio of HCO+/CO can be written as:

[TABLE]

Figure 7 shows the predicted abundance ratio of HCO+/CO with $N_{\rm CO}$ from 1D slab model simulations ($n_{\rm H}$ = 102 cm*-3*), overlaid with the observed values from absorption measurements. Red dots are the observed values, and different line curves represent PDR models under different CRIRs. As seen from Fig. 7, the majority (22 out of the total 26) of the observational values are within $10^{-16}\lesssim\zeta_{2}\lesssim 10^{-15}\,{\rm s}^{-1}$.

At a given $N_{\rm CO}$, the abundance ratio of HCO+/CO will first increase with increasing CRIR. This is similar to the variance of OH/CO as we have explained in Section 6.1. With the increasing photodissociation rate (R24) and the gas temperatures, reaction rates R15, R16, R19, and R23 decrease while the reaction rates of electron recombination (R21) increase. Thus, the increasing trend of HCO+/CO is slowing down because of the increasing destruction rate of HCO+ by free electrons. Consequently, there is a plateau of the HCO+/CO ratio (as well as the abundance of HCO+), where the increasing CRIR no longer increases the HCO+/CO ratio.

However, the HCO+/CO ratio decreases at an even higher CRIR. This is because the most important precursors of HCO+ (R18), such as CO+, are destroyed by atomic hydrogen and free electrons at high temperature, reducing the formation efficiency of HCO+. The abundance ratio of HCO+/CO reaches a maximum when the CRIR is approximately 10*-15* s*-1* for $n_{\rm H}=10^{2}\,{\rm cm}^{-3}$.

6.3 The CRIR inferred from HCO+ and CO absorptions

In order to quantitatively constrain $\zeta_{2}$ toward each source, we have run PDR models with ${\rm log_{10}\,(\zeta_{2}/s^{-1})}=-18$ to $-14$ with a step size of 0.07 dex and ${\rm log_{10}\,({\it n}_{\rm H}/cm^{-3})}$ = 1 – 3.5 with a step size of 0.1 dex, interacting with three FUV intensities ($\chi/\chi_{0}=1,5,10$). We vary $\zeta_{2}$ and $n_{\rm H}$ to find the optimum model by minimizing the reduced $\chi^{2}$ function:

[TABLE]

where $\rm{f_{obs}^{i}}$ and ${\rm\sigma^{i}_{obs}}$ are the observed molecular column densities and their corresponding uncertainties of the $i$th species, respectively. $\rm{f_{model}^{i}}$ is the column density from PDR models and N = 1 is the degree of freedom.

We show two representative ${\rm\chi^{2}_{red}}$ distributions from our fitting in Figure 8. The higher values of ${\rm\chi^{2}_{red}}$ are colored with bright blue and the lower colored with dark blue. White color indicates values beyond the maximum in this colour-bar. Red contours represent the boundary of the best-fit parameter where the deviation between modeled values and observations is 1$\sigma$. We treat a fitting result with a minimal value of ${\rm\chi^{2}_{red}}\lesssim 1$ (upper panel of Figure 8) as a “high confidence” fit, and a fitting result with a minimal value of ${\rm\chi^{2}_{red}}\gg 1$ (lower panel of Figure 8) as a “low confidence” fit.

As seen from Figure 8, the best-fit values drift to both higher $\zeta_{2}$ and $n_{\rm H}$ as the increasing FUV intensity. However, at high FUV intensity, the best-fit value of $n_{\rm H}$ from PDR models is much higher than the allowable density range by MCMC runs (Table 1. On the other hand, increasing $\chi$/$\chi_{0}$ would increase $T_{\rm k}$. In that case, the inferred gas density should be lower (as seen in Figure 1). Thus, it is not likely that FUV intensity is high for our sources.

Here, we made a rough estimation of the FUV intensity in diffuse clouds. The interstellar radiation field is a function of dust temperature ($T_{\rm d}$, Beuther et al., 2014):

[TABLE]

where $\gamma(A_{\rm V})$ is the attenuation factor. We consider $\gamma(A_{\rm V})$ = 0.35 at $A_{\rm V}$ = 1 mag (Glover & Clark, 2012). For $T_{\rm d}$ = 15 K, the derived FUV intensity is $\chi/\chi_{0}$ = 1.4. A variation on $T_{\rm d}$ of $\pm$2 K would only result in a difference on $\chi/\chi_{0}$ by a factor of 2. In the following, we only consider the fitting results at $\chi$/$\chi_{0}$ = 1.

The fitting results of all sources are summarised in Table 2, in which high confidence fit parameters are labelled with a “” at the end of each row. The $\zeta_{2}$ in our sample (high confidence results) is in the range of $\zeta_{2}\sim 1.0_{-1.0}^{+14.8}$ $\times$10${}^{-16}-2.5_{-2.4}^{+1.4}$ $\times$10-15* s*-1*. The $n_{\rm H}$ is in the range of (0.1–3.2)$\times$102 cm*-3*, which is consistent with that obtained from MCMC runs in Section 3.1.

In Fig. 9, we plot the Luo et al. (2023) $\zeta_{2}$ toward IC 348 measured from HCO+ and CO measurements, and the Indriolo & McCall (2012) $\zeta_{2}$ from H${}_{3}^{+}$ measurements. For comparison, we also include the Padovani et al. (2022) CR attenuation models $\mathscr{L}$, $\mathscr{H}$, and with low-energy spectral slope = $-$1.2. To plot $\zeta_{2}$ as a function of $E$(B-V), we convert N${}_{\rm H_{2}}$ to $E$(B-V) by assuming that the atomic fraction $f_{\rm atomic}$ follows the power law: $f_{\rm atomic}$ = 0.64/$A_{\rm V}$ at $A_{\rm V}>$1 mag (Luo et al., 2023), and the $f_{\rm atomic}$ = 0.64888The atomic fraction is between 0.6-0.8 at low extinction LOSs (Luo et al., 2020). at $A_{\rm V}\leq$1 mag999Note that the discontinuity seen in models at $E$(B-V) = 0.32 mag is not a peculiarity of the models themselves, but it is due to the assumption on $f_{\rm atomic}$.. The average value of CRIR in our samples (red dots) is (8.0$\pm$6.4)$\times$10*-16* s*-1*. This is $\sim 2$ times higher than the value measured with H${}^{+}_{3}$ toward nearby massive stars ($\zeta_{2}$ = 3.5${}^{+5.3}_{-3.0}$ $\times$10*-16* s*-1*) and toward IC 348 (4.7$\pm$1.5 $\times$10*-16* s*-1*, see Fig. 9). Considering that our sightlines are located at much lower $E$(B-V)101010We convert $A_{\rm V}$ to $E$(B-V) by adopting $A_{\rm V}$ = 3.1$E$(B-V) (Schlafly & Finkbeiner, 2011). (0.1 – 0.7 mag), the values are reasonably consistent with the low $E$(B-V) portion by Indriolo & McCall (2012). Our results are consistent with the model $\mathscr{H}$ by Padovani et al. (2018), in which CRs (and consequently the CRIR) are attenuated as a function of $A_{\rm V}$ (Padovani et al., 2018, 2022; Gaches et al., 2022). However, model $\mathscr{L}$, which corresponds to a “low” cosmic-ray spectrum based on Voyager-1 data, is found to underestimate the CRIR by almost an order of magnitude.

7 Conclusions

We analyze the abundances of CO, OH, and HCO+, and their abundance ratios in chemical equilibrium. We present a new approach to constrain the CRIR using these species. We calculated the column densities of HCO+ and CO toward diffuse LOSs against quasars and compare them with 3d-pdr models. Our inferred values of $\zeta_{2}$ show good consistency with previous measurements using H${}^{+}_{3}$. The main conclusions are as follows:

The gas volume densities ($n_{\rm H_{2}}$) obtained from CO (1–0) and (2–1) transitions toward 4 LOSs are in the range of (0.14$\pm$0.03 – 1.2$\pm$0.1)$\times$102 cm*-3* (at $T_{\rm k}$ = 50 K), which suggests a sub-thermal excitation environment for CO and HCO+. 2. 2.

Analyzing the chemical response of different molecules, we found that the abundance of CO increases with gas density and decreases with increasing FUV intensity, while the abundances of OH and HCO+ mostly have an opposite trend to that of CO. 3. 3.

Our analytic expressions give an excellent abundance of OH when $\zeta_{2}$ $<$10*-15* s*-1*. This is potentially useful for hydrodynamical simulations as it reduces the computational expense of chemical networks. In the diffuse gas, the abundance of OH is proportional to the CRIR. 4. 4.

At a given $N_{\rm CO}$, the abundance ratio of OH/CO is monotonically increasing with the increase $\zeta_{2}$ in the diffuse cloud, while the abundance ratio of HCO+/CO increases reaching a local maximum value at $\zeta_{2}$ $\approx$ 10*-15* s*-1*, before it decreases again. The downward trend of HCO+/CO at high $\zeta_{2}$ is caused due to the increased destruction of the precursor of HCO+–CO+ by electrons and atomic H at higher gas temperatures. 5. 5.

By comparing the observational values of HCO+/CO and chemical models, we find that the average $\zeta_{2}$ in our sample is (8.0$\pm$6.4)$\times$10*-15* s*-1*. This value is $\sim$2 times higher than that of higher extinction regions, which is consistent with the hypothesis of decreasing $\zeta_{2}$ as the increasing $A_{\rm V}$ in theoretical studies.

With high-sensitivity measurements from HCO+ and CO, we have demonstrated the possibility of using HCO+/CO ratios to constrain the CRIR in diffuse gas without knowing the exact molecular abundances relative to H2. We propose that the abundance ratios of OH/CO and HCO+/CO can be used to constrain the CRIR, especially with interferometry observations where high-resolution H2 information is inaccessible. Due to the difficulties in obtaining accurate excitation and optical depth of OH with current facilities, we foresee that future instruments, such as SKA, will produce large samples of data sets for which our approach will be potentially useful.

We thank Marco Padovani for the useful discussions on the CR attenuation models and for providing the latest data of these models. This work has been supported by the National Natural Science Foundation of China (grant No. 12041305), China Postdoctoral Science Foundation (grant No. 2021M691533), the Program for Innovative Talents, Entrepreneur in Jiangsu, the science research grants from the China Manned Space Project with NO.CMS-CSST-2021-A08. P. Z. acknowledges support from the National Natural Science Foundation of China (grant No. 12273010). This paper makes use of the following ALMA data: ADS/JAO.ALMA$\#$2015.1.00503.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) and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ.

Appendix A Derivation of abundance of OH in chemical equilibrium

The formation of OH is given by:

[TABLE]

where $\rm{\it n}(H_{2}O^{+})$ can be written as:

[TABLE]

and $\rm{\it n}(H_{2}O^{+})$ can be written as:

[TABLE]

The destruction of H2O is dominated by photodissociation and C+ in the diffuse cloud (note that the destruction of H2O will be dominated by H+ in dense cloud), thus, we have

[TABLE]

Then, equation A1 can be written as:

[TABLE]

where $\delta$ is:

[TABLE]

Considering chemical equilibrium of $\rm{\it n}(OH^{+})$, the term $\rm{\it n}(OH^{+}){\it n}(H_{2})k_{R4}$ can be written as:

[TABLE]

where $\zeta_{\rm O}$ $\approx$ 2.83$\zeta_{2}$ is the cosmic-ray ionization rate of atomic oxygen.

H+ is produced by cosmic-rays (CRs):

[TABLE]

and removed by atomic oxygen (reaction R1) and electron recombination reaction:

[TABLE]

In chemical equilibrium, $\rm{{\it n}(H^{+})}$ is in the form of:

[TABLE]

in which H${}_{2}^{+}$ is formed through:

[TABLE]

and removed by H through reaction RA3 or by H2 through:

[TABLE]

Thus, equation A8 can be written as:

[TABLE]

where $\rm{k_{RA1}}$ = $\zeta_{1}$ = 1/2$\zeta_{2}$, $\rm{k_{RA2}}$ = 0.02 $\zeta_{2}$, $\rm{k_{RA5}}$ = 0.88 $\zeta_{2}$, and $\xi$ is in the form of:

[TABLE]

Thus, if $\rm{{\it f}_{mol}}$ is close to 1, $\xi$=0. Otherwise, if $\rm{{\it f}_{mol}}$ is close to 0, $\xi$=1.

However, equation A9 does not consider the formation of H+ by He+. The formation of He+ is mainly due to the He ionization due to CRs:

[TABLE]

The destruction of He+ is either dominated by H or H2 in diffuse cloud:

[TABLE]

In equilibrium, n(He+) can be written as:

[TABLE]

where $\rm{k_{RA7}}$ = 0.5 $\zeta_{2}$.

If we take reactions RA8 and RA9 into account, equation A9 can be rewritten as:

[TABLE]

H${}_{3}^{+}$ is formed through reaction RA6 and removed by electrons:

[TABLE]

Thus, $\rm{{\it n}(H_{3}^{+})}$ is in the form of:

[TABLE]

Note that equation A13 does not consider the destruction of H${}^{+}_{3}$ by CO and O (see also equation 18 in Indriolo & McCall (2012)), however, equation A13 is a good approximation in diffuse and translucent clouds since the destruction rate by electrons is more than two orders of magnitude higher than that of CO or O.

If we substitute equation A7 with equation A9 and A13, we have:

[TABLE]

where $\theta$ is in the form of:

[TABLE]

which denotes the last term in equation A14. Combining equations A5 and A14, we thus obtain:

[TABLE]

The destruction of OH is given by:

[TABLE]

In equilibrium, the formation and destruction of OH reach a balance ($\rm F[{\it n}(OH)]$ = $\rm D[{\it n}(OH)]$).