Issue 
A&A
Volume 611, March 2018



Article Number  A20  
Number of page(s)  8  
Section  Interstellar and circumstellar matter  
DOI  https://doi.org/10.1051/00046361/201731871  
Published online  21 March 2018 
A simple approach to CO cooling in molecular clouds
School of Physics and Astronomy, Cardiff University,
Cardiff
CF24 3AA, UK
email: ant@astro.cf.ac.uk
Received:
31
August
2017
Accepted:
14
December
2017
Carbon monoxide plays an important role in interstellar molecular clouds, both as a coolant, and as a diagnostic molecule. However, a proper evaluation of the cooling rate due to CO requires a determination of the populations of many levels, the spontaneous and stimulated radiative deexcitation rates between these levels, and the transfer of the emitted multiline radiation; additionally, this must be done for three isotopologues. It would be useful to have a simple analytic formulation that avoided these complications and the associated computational overhead; this could then be used in situations where CO plays an important role as a coolant, but the details of this role are not the main concern. We derive such a formulation here, by first considering the two asymptotic forms that obtain in the limits of (a) low volumedensity and optical depth, and (b) high volumedensity and optical depth. These forms are then combined in such a way as to fit the detailed numerical results from Goldsmith & Langer (1978, ApJ, 222, 881; hereafter GL78). The GL78 results cover low temperatures, and a range of physical conditions where the interplay of thermal and subthermal excitation, opticaldepth effects, and the contributions from rare isotopologues, are all important. The fit is obtained using the MetropolisHastings method, and reproduces the results of GL78 well. It is a purely local and analytic function of state — specifically a function of the density, ρ, isothermal sound speed, a, CO abundance, X_{CO}, and velocity divergence, ∇⋅υ. As an illustration of its use, we consider the cooling layer following a slow steady nonmagnetic planar Jshock. We show that, in this idealised configuration, if the postshock cooling is dominated by CO and its isotopologues, the thickness of the postshock cooling layer is very small and approximately independent of the preshock velocity, υ_{o}, or preshock isothermal sound speed, a_{o}.
Key words: hydrodynamics / molecular processes / radiation mechanisms: thermal / shock waves / ISM: clouds
© ESO 2018
1 Introduction
Carbon monoxide, CO, is a critical molecule in the physics of molecular clouds and star formation. First, it is believed to be the next most abundant molecule in the Universe after molecular hydrogen, H_{2}. Second, as compared with H_{2}, it has a relatively high moment of inertia, and a permanent dipole moment, so it emits readily at the low temperatures found in molecular clouds, whereas H_{2} does not; therefore CO and its isotopologues are often used to trace the structure and dynamics of molecular clouds (e.g. RomanDuval et al. 2016). Third, because it emits readily, it plays an important role in the thermal balance of molecular clouds, helping to maintain their pervasive low temperatures. However, the physics underlying the formation/destruction of CO, and the net line emission from CO, is complicated.
The pathways to CO formation are very diverse and uncertain; they depend on the formation of several other species, on reaction rates that are only known approximately, and on the agency of cosmic rays to produce molecular ions. Likewise, the destruction of CO is difficult to model, largely because it involves line radiation; as a consequence, evaluating selfshielding is difficult, and is made more difficult because selfshielding has to compete with some lines also being blocked by H_{2}, with dust attenuation, and, at high densities, with CO freezeout. A variety of schemes has been proposed to estimate the abundance of CO without considering all these details (e.g. Nelson & Langer 1997, 1999), and indeed these approximate schemes are the ones used in many interstellar chemistry codes (e.g. Glover & Clark 2012a). In the Solar vicinity, it seems that gasphase CO is usually the dominant form of carbon at volume densities in the range and columndensities . At lower volume and columndensities, CO gives way to atomic and ionic carbon (C^{o} and C^{+}), due to the slow rate of the twobody reactions leading to CO formation, and the strong ultraviolet irradiation which rapidly destroys CO. At higher volumedensities, CO appears to freeze out onto dust; parenthetically, CO’s role as a coolant also becomes less important at these higher densities, because the gas starts to couple thermally to the dust, and hence to cool by continuum emission from dust.
Line emission from CO is complicated by the fact that, in the range of volume and columndensity where CO tends to be abundant, transfer of population between the different levels on the rotational ladder of CO involves a balance between collisional and radiative excitation and deexcitation, and the details of this balance shift with changing volumedensity, columndensity, temperature and velocity dispersion. CO line emission is seldom in either of the asymptotic limits of very low volume and columndensity (hence low linecentre optical depth), or very high volume and columndensity (hence high linecentre optical depth), where its collective line emission can be described by simple algebraic equations (e.g. Goldreich & Kwan 1974).
Here we derive an approximate analytic formulation for the net cooling rate from CO and its isotopologues. The formulation is simple in the sense that it is a purely local function of state, depending only on four parameters: (i) the massdensity, ρ, or equivalently the volumedensity of molecular hydrogen, ; (ii) the isothermal sound speed, a, or equivalently the gaskinetic temperature, T; (iii) the relative abundance of CO, ; and (iv) the velocity divergence, ∇ ⋅υ. The formulation is obtained by first deriving the simple algebraic equations that describe the dependence of the CO cooling rate on density, temperature, CO abundance and velocity divergence, in the asymptotic limit of low volumedensity and optical depth, and in the asymptotic limit of high volumedensity and optical depth; hereafter we refer to these equations as the asymptotic forms. Their derivation is based on fundamental physical arguments, and forms the subject of Sect. 2. Then, in Sect. 3 the coefficients in front of the asymptotic forms, and the variation between them, are fit by reproducing the detailed results of GL78, using an ad hoc mathematical function. In Sect. 4 we illustrate the application of the approximate analytic formulation by using it to evaluate the thickness of the postshock cooling layer behind a lowvelocity steady Jshock in a nonmagnetic molecular cloud. In Sect. 5 we summarise our conclusions.
For mathematical convenience, we use the massdensity, ρ and the isothermal sound speed, (in place of the numberdensity of molecular hydrogen, , and the gaskinetic temperature T) in much of the analysis. Here is Boltzmann’s constant and is the mean gasparticle mass. For the purpose of illustration, we use a reference temperature To = 10 K and a mean gasparticle mass (appropriate for molecular gas with elemental composition X = 0.70, Y = 0.28, Z = 0.02);^{1} hence the reference isothermal sound speed is ao = 0.187 km s^{−1}. On the assumption that virtually all the hydrogen is in the form of H_{2} , we also define the mean mass per hydrogen molecule, , and hence .
2 The CO cooling rate: asymptotic theory
In this section we derive, using basic physical arguments, the dependence of the CO line cooling rate on the density, on the isothermal sound speed (or equivalently the temperature), on the gasphase CO abundance, and on the local velocity divergence. We treat the two asymptotic limits of (i) very low volumedensity and optical depth (hereafter the “lo” limit), and (ii) very high volumedensity and optical depth (hereafter the “hi” limit). We do not derive the absolute physical value of the cooling rate in these limits, but only the dependence on density, isothermal sound speed, abundance and velocity divergence. The coefficients converting these dependences into actual physical values are obtained in the following section (Sect. 3) by comparing the theoretical predictions with the detailed results of GL78.
2.1 CO cooling in the lo limit
For a linear molecule with moment of inertia , the rotational levels are characterised by quantum number J, and have energy therefore, if the gaskinetic temperature is T, they are significantly excited up to (1)
where is a factor greater than, but of order, unity. Radiative deexcitations from level J to level J − 1 release photons with energy . Therefore the mean energy of emitted photons is (2)
Here, we have replaced a sum over levels with an integral. This is reasonable for higher temperatures, T≳20 K, but may introduce a small systematic error at lower temperatures where only a few levels are involved.
At sufficiently low volumedensity, i.e. well below the critical density, the level populations are not thermalised. Most molecules sit in their ground state, and most collisional excitations are followed by radiative deexcitations. For example, at T ~ 10 K, the first rotationally excited level of CO, J = 1, has a critical density of order 10^{3} H_2 cm^{−3}; for higherJ levels, higher temperatures and densities are required before they are thermalised. At sufficiently low optical depth (i.e. sufficiently low columndensity and/or high velocity divergence), the lines are optically thin, so most of the emitted photons escape directly; consequently the thermal kinetic energy that caused the initial excitation is lost, and the gas is cooled. For example, at T ~ 10 K and with purely thermal velocity dispersion, i.e. no turbulence and no velocity gradient, the first rotational transition of CO (J = 1 → 0) becomes optically thick for columndensities ; this limit is increased if there is turbulence and/or a velocity gradient delivering extra velocity dispersion. The rate of collisional excitation per unit volume depends on the product of the numberdensities of the molecule and the exciting collider (here presumed to be H_{2} , see below) and their relative speed, so it is approximately proportional to . The number of lines that are excited is approximately proportional to T^{1∕2} (see Eq. (1)). And the mean energy of the emitted photons is approximately proportional to T^{1∕2} (see Eq. (2)). Therefore, in the limit of very low volumedensity and optical depth, the cooling rate per unit volume can be approximated by (2)
In treating cooling by CO, we presume that H_{2} is the dominant agent of collisional excitation, since, in regions where there is CO, H_{2} is likely to be by far the most abundant species; GL78 make the same assumption.
On the basis of Eq. (3), we posit that, in the lo limit, the cooling rate due to CO is given by (4)
The coefficient will be estimated in Sect. 3, but the above dependence of on , ρ and a is now fixed.
2.2 CO cooling in the hi limit
At sufficiently high volumedensity, i.e. well above the critical density, the level populations are approximately thermalised, and radiative deexcitations account for only a small fraction of deexcitations. At sufficiently high optical depth (i.e. sufficiently high columndensity and/or low velocity divergence), the intensity in the emission line cores is approximately equal to the Planck function, and the net frequency width of all the excited lines is (5)
Here Δυ is the velocity spread along the line of sight and is the factor by which the effective width of an optically thick line with central frequency νo differs from^{2} νo Δυ∕c. The effective width of the Planck spectrum is where = 14.4, and so the total cooling rate per unit volume in this limit can be approximated by
In deriving Eq. (6), we have adopted the most generic geometry we could imagine, viz. a spherical cloud of diameter D expanding or collapsing homologously, so that Δυ∕D ≡∇⋅υ∕3; GL78 adopt the same configuration. However, in obtaining Eq. (7), this choice becomes immaterial, since we only retain the dependence on a and ∇ ⋅ υ. The effect of geometry is encapsulated entirely in the ∇⋅υ term, which measures how far, in different directions on the sky, linephotons have to travel before the Dopplershift has separated their frequency from the frequencies that can be absorbed by the molecules they are passing.
In effect we are adopting the large velocity gradient (LVG) approximation to treat optically thick cooling radiation. The LVG approximation was introduced by Sobolev (1960), and developed further by Castor (1970) and Lucy (1971), in the context of stellar outflows, but it can also be applied to the cooling of molecular clouds (e.g. Goldreich & Kwan 1974). Basically it assumes that line radiation only experiences selfabsorption in a local region whose linear size is of order σ ∇ ⋅υ^{−1}, where σ is the local velocity dispersion (thermal plus turbulent) and υ is the local bulk velocity. Any material outside this region is moving at a sufficiently different velocity (bulk plus or minus the dispersion) from the emitting region that it cannot absorb its line radiation. Strictly, the LVG approximation requires that the largescale spatial variation of the bulk velocity is monotonic along any line of sight.
On the basis of Eq. (7), we posit that, in the hi limit, the CO cooling rate is given by (8)
the coefficient will be estimated in Sect. 3, but the above dependence of on a and ∇ ⋅ υ is now fixed — and is independentof ρ and .
We note parenthetically that the only parameters in Eq. (6) that depend on the specific molecule are and , and that this dependence is in general rather weak, so Eq. (6) is approximately valid for many other linear molecules in this high volumedensity, high optical depth limit. Thus, for example, the results presented by GL78 for cooling by O_{2} are very similar to those for CO (at the same density and temperature).
2.3 In between the lo and hi limits
In the next section (Sect. 3) we show how Eqs. (4) and (8) can be combined to obtain an approximate analytic formulation for the total CO cooling rate in the intermediate regime (between the lo and hi limits) using (9)
with β = β(ρ, a).
3 The CO cooling rate: calibration
In the preceding section (Sect. 2) we have obtained expressions for the CO cooling rate (in the lo and hi limits) in terms of the massdensity, ρ, the isothermal sound speed, a, the abundance of CO, , and the velocity divergence, ∇⋅υ. This is because these are the most convenient variables to use when treating the equations of hydrodynamics. Here, we replace ρ with the density of molecular hydrogen , and a with the temperature, T, since these are the variables used by GL78. It is straightforward to switch between the two, using and .
3.1 The detailed results of Goldsmith & Langer (1978)
GL78 consider a uniformdensity spherical cloud with a linear radial velocity field, and determine the net cooling rate, , from CO and its isotopologues. To do this, they solve in detail (i) the equations of statistical equilibrium, to determine CO level populations; then (ii) the equations of radiative emission (spontaneous and stimulated), to determine the net CO emissivity; finally (iii) they use an escape probability in lieu of solving the equation of radiative transfer. In their Fig. 2, which is the input we use to calibrate our approximate analytic formulation, they plot the quantity against , for different values of T and different values of^{3} , specifically T = 10, 20, 40, and 60K and . We note that is the cooling rate per CO molecule. It depends on , because – all other things being equal – the optical depth increases with increasing and decreases with increasing ∇ ⋅ υ. In the lo limit, where the emission is optically thin, is independent of .
Given the complexity of the physics underlying the net CO cooling rate, and the consequent computational cost of calculating it properly, it would be useful to have an approximate analytic formulation that captured the dependence on local variables but avoided the associated computational cost. This formulation could then be used in situations where CO cooling played an important role, but was not the principal interest, as for example, contracting prestellar cores, or the accretion shock at the boundary of an assembling filament.
3.2 The fitting procedure
To fit the GL78 results analytically, we have read from their plots the values of for each treated combination of T, i.e. T = 10, 20, 40, and 60K, and , i.e. and − 4, at the densities , thus a total of 140 discrete values. We estimate that the uncertainty which derives from our reading these values of GL78’s Fig. 2 by eye is <0.03 in i.e. at worst ± 7%, and less than the thickness of the lines in Fig. 1.
We then use a MetropolisHastings Markov Chain Monte Carlo algorithm to find the five fitting parameters in the formulation
that give the best fit to these discrete values.
We note that the exponents from Eqs. (4) and (8) are not being allowed to vary, they retain the values derived in Sect. 2. Because we are here fitting (rather than ), the ρ^{2} term in Eq. (4) has become n_{H2}^{1}, and the ρ^{0} term implicit in Eq. (8) has become . The a^{3} term in Eq. (4) has become T^{3∕2} and the a^{8} term in Eq. (4) has become T^{4}.
Fig. 1 Continuous curves give the total CO cooling rate, as predicted by the approximate analytic formulation derived here (Eqs. (15) through (18)). The different panels correspond to (a) T = 10K ; (b) T = 20K; (c) T = 40K; and (d) T = 60K (these temperatures are given in the top lefthand corner of each panel). The different curves represent different values of (these values are given down the righthand margin of each panel). The rates obtained by the detailed computations of GL78 are given at discrete values of the density, ; plus signs, stars, open circles, crosses and open squares correspond respectively to . 
3.3 The best fit and its accuracy
The best fit to the results from GL78 is obtained with the following values for the fitting parameters: (14)
In otherwords, the best fit is
The continuous lines in Fig. 1 show the predictions of Eqs. (15) through (18) for T = 10, 20, 40, and 60K (respectively panels a, b, c, d), and (separate lines labelled on the righthand margin of each panel). Like GL78 we limit the plots to densities in the range . The continuous curves should be compared with the discrete symbols, which give the values that we have read from Fig. 2 of GL78 and that we used in the MetropolisHastings fitting procedure.
In general, the agreement is good, with the magnitude of the fractional offset being everywhere <0.18 in , and on average ≲0.07 in . The fractional offset is worst in the intermediate regime. This is as expected, since the analytic forms that obtain in the asymptotic limits (lo and hi) are physically motivated, and therefore in some sense absolute. In contrast, the intermediate regime between these limits is being fit with an algebraic expression that has no physical motivation beyond the fact that it relaxes to the asymptotic forms in the corresponding limits; it is therefore unable to produce an exact fit. The error in reading the discrete values of GL78’s Fig. 1 by eye ( ± 3%) is much smaller than the fractional offsets cited above, and can therefore be ignored.
There are some additional uncertainties that should be noted. First, GL78 used the coefficients for collisional deexcitation of CO by H_{2} computed by Green & Thaddeus (1976). More recent computations by Yang et al. (2010) have increased these coefficients somewhat, at the lowest temperatures, T≲20 K. Therefore the GL78 cooling rates, and our approximate analytic formulation based on them, may be a little low at these low temperatures. Second, the device of replacing with integrals, summations over the discrete contributions from individual levels, will be most inaccurate at low temperatures where only a few levels are involved.
We note that the rarer isotopologues of CO make little contribution to the cooling in the lo limit, but do contribute once the ^{12}C^{16}O lines start to become optically thick, because the rarer isotopologues remain optically thin to much higher columndensities. They also make an important contribution in the hi limit, by providing extra bandwidth at high columndensities; in fact, in the hi limit, each isotopologue makes essentially the same contribution to the net cooling rate.
3.4 Use and range of applicability of the approximate analytic formulation
Given this level of accuracy, the approximate analytic formulation (i.e. Eqs. (15) through (18)) provides a convenient and computationally very inexpensive way to estimate the net CO cooling rate: convenient because it is a local function of state, and computationally inexpensive because it entails only a small number of arithmetic operations. It can therefore be used to greatly speed up numerical and/or analytic integrations of the energy equation in hydrodynamic or hydrostatic simulations, provided only that the conditions are within the range explored by GL78.
The numerical results of GL78 cover the density range , the temperature range 10 K ≤ T ≤ 60 K, and the abundance/velocitydivergence range , and therefore this is the range over which our approximate analytic formulation is most secure. Given the regularity of the curves in GL78, and their apparent selfsimilarity, it may be safe to extrapolate further. It is unlikely that gasphase CO is the dominant form of CO outside the density range defined above, so we only consider extrapolation of the temperature and abundance/velocitydivergence ranges. Specifically, we presume that the approximate analytic formulation can be extended to somewhat higher temperatures, T ~100 K, and somewhat lower abundances or higher velocity divergences, .
4 Slow steady planar Jshocks in nonmagnetic molecular clouds
As a simple example of the application of our approximate analytic formulation, we consider the role of CO in postshock cooling layers in molecular clouds. There are many other possible applications, for example the role of CO cooling in collapsing and fragmenting prestellar cores, and the thermodynamics of the accretion shock bounding a growing filament. We plan to investigate these in future papers.
Molecular clouds are observed to be highly turbulent, and supersonic velocity dispersions are the norm on scales ≳0.1 pc (Larson 1981; Solomon et al. 1987; Goodman et al. 1998; Caselli et al. 2002; Heyer & Brunt 2004). Shocks are therefore endemic in molecular clouds, and a critical issue is then how quickly the postshock gas cools back to something like its preshock value. In particular, if postshock cooling is very quick, shocks can for many purposes be treated as isothermal discontinuities, and this may greatly simplify analysis of the dynamics on larger scales.
In the interests of simplicity, we focus on slow, steady, planar Jshocks (e.g. Brand 1989; Flower et al. 2003). We note (i) that, if the shock velocity is high, υo≳15 km s^1, the associated kinetic energy, and resulting postshock thermal energy, may be sufficient to drive significant chemical changes; (ii) that these chemical changes may influence the dynamical coupling between the gas, the dust and the magnetic field; and (iii) that the cooling radiation from the hot gas immediately behind the shock may have important effects on the chemistry and thermal state of the gas flowing into the shock, as will the ambient radiation field. However, these three considerations are probably not important for the lowvelocity shocks, υo ≲1.5 km s^1, considered here.
4.1 The hydrodynamics of a Jshock plus cooling layer
We consider a slow steady planar Jshock at x = 0, with gas flowing in from x < 0, having density ρ(x < 0) = ρo, velocity parallel to the x axis υ(x < 0) = υo, and isothermal sound speed a(x < 0) = ao. There is no magnetic field, so we are ignoring the possibility that the shock is cushioned by a lateral magnetic field; in that case, the shock would be a Cshock, and the immediate postshock temperature would be somewhat lower, but, as shown by Pon et al. (2012), the cooling would still very likely be dominated by CO. We are also ignoring the possibility that the preshock gas is significantly heated by the radiation from the postshock cooling layer; this is a reasonable assumption, since gas in the postshock cooling layer is necessarily at a very different velocity to the gas flowing into the shock, and therefore the CO cooling radiation cannot be absorbed by the inflowing gas.
We assume that the gas flowing through the shock is and remains molecular, with constant adiabatic exponent γ = 5∕3. This is strictly only true in lowtemperature regions, T≲100 K, where the rotational degrees of freedom of molecular hydrogen are not significantly excited (Whitworth & Clarke 1997; Boley et al. 2007). This requirement is naturally met, since the GL78 results (and hence our approximate analytic formulation) are limited to temperatures T≲60 K. If the preshock temperature is ~10 K, this means we must restrict our model to shocks with υo≲1.0 km s^1. If we extrapolate the approximate analytic formulation to T≲100 K, we can treat shocks with υo≲1.5 km s^1.
The gas emerging from the shock has density, ρs, velocity, υs, and isothermal sound speed, as, given by the Rankine–Hugoniot conditions:
where the second expressions apply provided υo^{2} ≫ ao^{2}, and we recall that here ao is the isothermal sound speed.
In the postshock region (x > 0), conservation of mass requires constant ρ(x)υ(x), hence (22)
conservation of momentum requires constant , hence (23)
and conservation of energy requires (24)
where Λ(ρ, a, dυ∕dx) is the cooling rate per unit volume, and in planar geometry ∇⋅υ →dυ∕dx.
Since, in the postshock cooling layer, υ(x)≲υo∕4, Eqs. (23) and (24) approximate to
where we have substituted dυ∕dx = −dυ∕ddx, since the gas here is decelerating.
4.2 The thermodynamics of a Jshock plus cooling layer
We assume that the abundance of CO is fixed at . In other words, we assume that essentially all the carbon is in CO, that the abundance of C is 1.5 × 10^{−4} (cf. Sembach et al. 2000), and that CO is not destroyed in the shock, nor has it had time to freeze out significantly (e.g. Goldsmith 2001).
Substituting from Eq. (9) into Eq. (26), we obtain (27)
the second expression is obtained by substituting for ρ and a from Eqs. (22) and (25).
In Eq. (27), the leading term (outside the braces) represents deceleration due to cooling in the lowdensity optically thin regime. The term in braces represents the correction due to the cooling trending towards the highdensity optically thick regime, and this can in principle cause the postshock cooling to stall. Since the gas must decelerate to υ = ao (see discussion following Eq. (29) below), it follows that postshock cooling will not stall due to the CO cooling becoming optically thick, provided that (28)
where is the numberdensity of molecular hydrogen in the preshock gas. If this inequality is well satisfied, we can neglect the term in braces, and this is generally the case. For example, if the preshock density is high, say , and the temperature is low, say To = 10 K, we require υo ≫ 0.012 km s^{−1}, which is a rather weak constraint on υo. In most cases of interest, the preshock density is lower than 10^{5} cm^{−3}, and the temperature cannot be much below ~10 K, in which case the minimum velocity is even lower. Furthermore, if the preshock density is any higher than 10^{5} cm^{−3}, the gas is probably so well thermally coupled to the dust that molecular line cooling is redundant (e.g. Glover & Clark 2012b).
4.3 The thickness of the postshock cooling layer
The thickness of the postshock cooling layer is (29)
The lower limit on this integral is the immediate postshock velocity, given by Eq. (20). If we were to follow the gas until it cooled right back down to , the upper limit on the integral would be υ = ao^{2}∕υo. However, by this stage the velocity would be highly subsonic and the LVG approximation would no longer be valid (in the sense that it would be seriously underestimating the cooling rate). Therefore we set the upper limit to υ = ao; in other words, we follow the deceleration until the velocity becomes subsonic. We believe this is justified, since the integral in Eq. (29) is dominated by the lower limit (υ = υo∕4).
If we substitute from Eq. (27) in Eq. (29), neglecting the correction term for optical thickness, we obtain (30)
We have replaced ≃ with ≲ in Eq. (30) – and also in Eq. (32) below – because (a) we have neglected the contribution from the upper limit in the integral of Eq. (29), and this contribution can be quite large for relatively slow shocks, thereby reducing further; (b) for fast shocks there will be a significant additional cooling contribution from dust (Glover & Clark 2012b), and possibly also other molecules (Neufeld et al. 1995), due to the high postshock density, reducing still further.
The gas passes through in a time
and the integrated flux in CO cooling lines from the shock is
4.4 Caveats on shock model
As already stated, gasphase CO appears to be the dominant form of carbon at volume densities in the range , and therefore our shock model can only be applied to shocks involving preshock gas with density, , in this range. Many shocks arising in turbulent molecular clouds will involve gas with density in this range.
Jshocks in turbulent molecular clouds will only approximate to being steady and planar if the counterflows that form them have coherence lengths significantly larger than the thickness of the postshock cooling layer, . If we adopt coherence lengths from Larson’s relations, i.e. L ~ 1600 (Larson 1981), this condition is very well satisfied. For example, if we set , we have . This suggests that the assumption of a steady planar shock is reasonable.
If we assume that the approximate analytic formulation can be extended to T ~ 100 K, and that the preshock gas has T ~ 10 K, then our shock model can only be applied to shocks with preshock velocity υo ≲1.5 km s^1, otherwise the immediate postshock gaskinetic temperature (i.e. pre cooling) is too high. This is a significant, but not critical, restriction, since there are likely to be many shocks satisfying this constraint.
For example, columndensity maps of lowmass cores (Könyves et al. 2015) and filaments (Palmeirim et al. 2013) have centrally condensed profiles. However, the columndensity contrast between the background and the line of sight through the centre of the core or the spine of the filament is small, typically ≲10. This indicates that the ram pressure of the gas accreting onto the core or filament is low, otherwise the profile would be flat. It also suggests that the background volumedensity is not hugely lower than that inside the core or filament, hence the Mach number of the approximately isothermal accretion shock at the boundary of the core or filament must be low, , and the preshock velocity must also be low, υo≲1 km s^1. In addition we note that gravitational acceleration of the material accreting onto a marginally Jeansunstable core or Ostrikerunstable filament can only generate marginally transsonic velocities, so gravitational acceleration does not alter this conclusion.
To first order, we can approximate (35)
using Eq. (30). If again we set , then the limits on , i.e. yield (36)
Since we already have the constraint υo≲1.5 km s^{−1}, this reduces to , which is only slightly more restrictive than the upper limit on derived at the start of this section.
Our shock model assumes that the carbon in the preshock gas is mainly in CO, and that the CO is not destroyed in the shock, nor does it freeze out. If this is not the case, cooling will be provided by other species, for example C^{+} (see Glover & Clark 2012a), but the gas will not cool so fast, and will probably not be able to get down to ~ 10 K, unless the postshock density becomes high enough for the gas to couple thermally to the dust.
4.5 Discussion of shock results
The thickness of the postshock cooling layer, , and the postshock cooling time, , are both likely to be very small, as compared with other length and timescales in the molecular cloud. With typical preshock densities in a molecular cloud, , and setting , we have and . These values are similar to the values derived by Pon et al. (2014, 2016), although they treated magnetically cushioned Cshocks, in which the compression and heating are less abrupt and less extreme.
One reason why CO is an effective postshock coolant is that it has a low dipole moment, and therefore its lines remain optically thin up to quite high columndensities. A second reason is that the steady deceleration of the postshock gas ensures that there is a broad range of velocities, and hence wavelengths, over which each line can be emitted, and therefore the lines do not readily become optically thick.
When treating the dynamics of turbulent interstellar clouds – for example, their assembly, cloud/cloud collisions, the formation of sheets and filaments and cores – Eqs. (30) and (32) can be used to assess rather quickly whether it is acceptable to treat the associated shocks as isothermal discontinuities, or alternatively, what resolution is needed to resolve the postshock cooling layers. That is, provided that the inherent assumptions in our analysis are valid, viz. steady slow nonmagnetic planar Jshocks, in which CO dominates the cooling and is not destroyed or removed from the gasphase (as discussed in Sect. 4.4).
Pon et al. (2012) and Pon et al. (2016) estimate that the volume fraction of a turbulent molecular cloud that is involved in shock dissipation (i.e. the fraction that is in postshock cooling layers) is . Given the short cooling times estimated above (Eq. (32)), it follows that a representative fluid element in such a cloud should be shocked at least once every Myr. Indeed, Pan & Padoan (2009) identify such shocks as an important overall heating mechanism for molecular clouds, although our estimates indicate that this heating should be very localised and transient.
5 Conclusions
We have derived a simple approximate analytic formulation for the cooling rate due to rotational transitions of the CO molecule, (37)
An equivalent formulation, but giving in terms of and T (instead of ρ and a), is given in Eqs. (15) through (18). Our approximate analytic formulation extends from the lowdensity optically thin limit, to the highdensity optically thick limit, and includes the contributions from the different isotopologues of CO. It is based on physical considerations, but is calibrated against the detailed numerical results of GL78, and reproduces those results well (see Fig. 1). It should be usable for
Using this formulation, we have derived estimates of the thickness of postshock cooling layers, , and postshock coolingtimes, , under the assumptions (i) that the shocks are steady slow nonmagnetic planar Jshocks, and (ii) that CO dominates the cooling, and is not destroyed or removed from the gasphase in the shock or the postshock cooling layer. With these caveats, and given the constraints on shock velocity, preshockdensity and preshock temperature detailed in Sect. 4.4, our estimates can be used to justify treating shocks as isothermal discontinuities when the primary concern is dynamics on larger scales; to evaluate the resolution required to model a shock; and to estimate the integrated CO flux from a postshock cooling region. The values of suggest that a typical fluid element in a molecular cloud should be shocked at least once every Myr.
Acknowledgements
We thank the anonymous referee for a very thoughtful and constructive report, which greatly improved the original version of this paper. APW gratefully acknowledges the support of a consolidated grant (ST/K00926/1), and SEJ a PhD studentship, both from the UK Science and Technology Funding Council. This work was performed using the computational facilities of the Advanced Research Computing at Cardiff (ARCCA) Division, Cardiff University.
References
 Boley, A. C., Hartquist, T. W., Durisen, R. H., & Michael, S. 2007, ApJ, 656, L89 [NASA ADS] [CrossRef] [Google Scholar]
 Brand, P. W. J. L. 1989, in IAU Colloq. 120: Structure and Dynamics of the Interstellar Medium, eds. G. TenorioTagle, M. Moles, & J. Melnick (Berlin: SpringerVerlag), Lecture Notes in Physics, 350, 38 [NASA ADS] [CrossRef] [Google Scholar]
 Caselli, P., Benson, P. J., Myers, P. C., & Tafalla, M. 2002, ApJ, 572, 238 [NASA ADS] [CrossRef] [Google Scholar]
 Castor, J. I. 1970, MNRAS, 149, 111 [NASA ADS] [Google Scholar]
 Flower, D. R., Le Bourlot, J., Pineau des Forêts, G., & Cabrit, S. 2003, MNRAS, 341, 70 [NASA ADS] [CrossRef] [Google Scholar]
 Glover, S. C. O., & Clark, P. C. 2012a, MNRAS, 421, 116 [NASA ADS] [Google Scholar]
 Glover, S. C. O., & Clark, P. C., 2012b, MNRAS, 421, 9 [NASA ADS] [Google Scholar]
 Goldreich, P., & Kwan, J. 1974, ApJ, 189, 441 [NASA ADS] [CrossRef] [Google Scholar]
 Goldsmith, P. F. 2001, ApJ, 557, 736 [NASA ADS] [CrossRef] [MathSciNet] [Google Scholar]
 Goldsmith, P. F., & Langer, W. D. 1978, ApJ, 222, 881 [NASA ADS] [CrossRef] [Google Scholar]
 Goodman, A. A., Barranco, J. A., Wilner, D. J., & Heyer, M. H. 1998, ApJ, 504, 223 [NASA ADS] [CrossRef] [Google Scholar]
 Green, S., & Thaddeus, P. 1976, ApJ, 205, 766 [NASA ADS] [CrossRef] [Google Scholar]
 Heyer, M. H., & Brunt, C. M. 2004, ApJ, 615, L45 [NASA ADS] [CrossRef] [Google Scholar]
 Könyves, V., André, P., Men’shchikov, A., et al. 2015, A&A, 584, A91 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Larson, R. B. 1981, MNRAS, 194, 809 [NASA ADS] [CrossRef] [Google Scholar]
 Lucy, L. B. 1971, ApJ, 163, 95 [NASA ADS] [CrossRef] [Google Scholar]
 Nelson, R. P., & Langer, W. D. 1997, ApJ, 482, 796 [NASA ADS] [CrossRef] [Google Scholar]
 Nelson, R. P., & Langer, W. D. 1999, ApJ, 524, 923 [NASA ADS] [CrossRef] [Google Scholar]
 Neufeld, D. A., Lepp, S., & Melnick, G. J. 1995, ApJS, 100, 132 [NASA ADS] [CrossRef] [Google Scholar]
 Palmeirim, P., André, P., Kirk, J., et al. 2013, A&A, 550, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pan, L., & Padoan, P. 2009, ApJ, 692, 594 [NASA ADS] [CrossRef] [Google Scholar]
 Pon, A., Johnstone, D., & Kaufman, M. J. 2012, ApJ, 748, 25 [NASA ADS] [CrossRef] [Google Scholar]
 Pon, A., Johnstone, D., Kaufman, M. J., Caselli, P., & Plume, R. 2014, MNRAS, 445, 1508 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Pon, A., Johnstone, D., Caselli, P., et al. 2016, A&A, 587, A96 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 RomanDuval, J., Heyer, M., Brunt, C. M., et al. 2016, ApJ, 818, 144 [NASA ADS] [CrossRef] [Google Scholar]
 Sembach, K. R., Howk, J. C., Ryans, R. S. I., & Keenan, F. P. 2000, ApJ, 528, 310 [NASA ADS] [CrossRef] [Google Scholar]
 Sobolev, V. V. 1960, Moving envelopes of stars (Cambridge: Harvard Univ. Press) [Google Scholar]
 Solomon, P. M., Rivolo, A. R., Barrett, J., & Yahil, A. 1987, ApJ, 319, 730 [NASA ADS] [CrossRef] [Google Scholar]
 Whitworth, A. P., & Clarke, C. J. 1997, MNRAS, 291, 578 [NASA ADS] [Google Scholar]
 Yang, B., Stancil, P. C., Balakrishnan, N., & Forrey, R. C. 2010, ApJ, 718, 1062 [NASA ADS] [CrossRef] [Google Scholar]
If the equivalent width of an optically thick emission line falls on the flat portion of the curve of growth, this factor is approximately constant, and that is what we assume here. If the equivalent width falls on the squareroot portion of the curve of growth, depends on the columndensity (as ). However, this requires extremely large columndensities. By the time such high columndensities are reached, molecules like CO are likely to have frozen out — and, even if they have not, the density is so high that molecularline cooling has given way to dust cooling.
All Figures
Fig. 1 Continuous curves give the total CO cooling rate, as predicted by the approximate analytic formulation derived here (Eqs. (15) through (18)). The different panels correspond to (a) T = 10K ; (b) T = 20K; (c) T = 40K; and (d) T = 60K (these temperatures are given in the top lefthand corner of each panel). The different curves represent different values of (these values are given down the righthand margin of each panel). The rates obtained by the detailed computations of GL78 are given at discrete values of the density, ; plus signs, stars, open circles, crosses and open squares correspond respectively to . 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.