The AtomictoMolecular Transition in Galaxies. III. A New Method for Determining the Molecular Content of Primordial and Dusty Clouds
Abstract
Understanding the molecular content of galaxies is a critical problem in star formation and galactic evolution. Here we present a new method, based on a Strömgrentype analysis, to calculate the amount of H I that surrounds a molecular cloud irradiated by an isotropic radiation field. We consider both planar and spherical clouds, and H formation either in the gas phase or catalyzed by dust grains. Under the assumption that the transition from atomic to molecular gas is sharp, our method gives the solution without any reference to the photodissociation cross section. We test our results for the planar case against those of a PDR code, and find typical accuracies of about 10%. Our results are also consistent with the scaling relations found in Paper I of this series, but they apply to a wider range of physical conditions. We present simple, accurate analytic fits to our results that are suitable for comparison to observations and to implementation in numerical and semianalytic models.
Subject headings:
ISM: clouds — ISM: molecules — molecular processes — radiative transfer — stars: formation1. Introduction
1.1. Description of the Problem
Stars form in molecular gas. It is therefore of central importance to understand how molecular gas forms in the interstellar medium of galaxies, which are pervaded by far ultraviolet (FUV) radiation that destroys the molecules by photodissociation. In the diffuse interstellar medium (ISM) of galaxies like the Milky Way, this radiation is sufficient to maintain the gas occupying most of the volume in an atomic or ionized state. Sufficient concentrations of gas are able to exclude the dissociating radiation, both by selfshielding and by dust absorption, and thereby become molecular; these are the molecular clouds. Regions in and around molecular clouds in which the thermal and chemical state of the gas is determined by FUV radiation are termed photodissociation regions (PDRs; e.g. Hollenbach & Tielens, 1999). The surface density of atomic gas in a PDR, , sets the threshold for gas to become molecular and thus to be able to form stars.
The study of PDRs is a vast field, but the majority of the work that has been done is not suited to the problem of making general statements about the thicknesses of HI layers. Most of it is highly numerical (Federman et al., 1979; van Dishoeck & Black, 1986; Black & van Dishoeck, 1987; Draine & Bertoldi, 1996; Neufeld & Spaans, 1996; Stoerzer et al., 1996; Spaans & Neufeld, 1997; Liszt & Lucas, 2000; Liszt, 2002; Browning et al., 2003), which is useful for trying to match detailed properties of particular regions, but makes it difficult to draw general conclusions that allow extrapolation beyond the particular numerical parameters selected for a calculation. The same problem of how to draw general conclusions from particular results arises with galacticscale simulations that model H formation and dissociation as one of many processes occurring in a galaxy (Robertson & Kravtsov, 2008; Gnedin et al., 2009). In contrast, the analytic work preceeding this series of papers has been limited to the comparatively simple onedimensional cases of a unidirectional beam of radiation impinging on a slab (Draine & Bertoldi, 1996; Sternberg, 1988) or a radiallyconverging radiation field striking the surface of a sphere (Elmegreen, 1993); the beam and slab case is a reasonable approximation for a PDR irradiated by a single star, but not for a cloud immersed in the interstellar radiation field produced by the ensemble of many stars. Since these authors assumed that the radiation is beamed, they were able to introduce a shielding function that describes the attenuation of the dissociating radiation with depth; this does not appear to be possible for isotropic incident radiation, particularly for spherical clouds. Moreover, these treatments did not identify the important dimensionless numbers that characterize the problem of a spherical cloud immersed in an isotropic radiation field. Perhaps the paper that is closest in spirit to the present series of papers is that of Sternberg (1988), who determined scaling laws for infrared fluorescent emission lines of H emitted by an irradiated slab. He assumed that the radiation was beamed, so that he could use a shielding function.
We determined the dimensionless quantities that characterize irradiated spherical clouds for the first time in Krumholz et al. (2008, hereafter Paper I), where we also gave an approximate solution for the atomictomolecular ratio of an irradiated gas cloud in terms of these numbers. In this paper we extend this analysis by introducing a new and more accurate method for solving the equations of radiative transfer and photodissociation. We also extend our treatment to the case of primordial H formation, where there is no dust and H instead forms via the H channel. Finally, we apply our results to the detemination of the molecular fraction in atomicmolecular complexes.
1.2. The Strömgren Analysis
PDRs are one example of a general class of problem in which one must find the chemical state of a gas cloud that simultaneously obeys coupled conditions of chemical and radiative equilibrium. The most famous example of this type of problem is that of determining the structure of an HII region, and the structure of PDRs and similar problems can be solved by adapting the classic solution to the HII region problem by Strömgren (1939). In the absence of dust, the thickness of the layer of ionized gas created by an incident flux of ionizing photons is
(1) 
where the superscript indicates photon quantities, is the density of electrons, and is the radiative recombination coefficient. Similarly, for a point source of ionizing radiation, the radius of the ionized gas (the Strömgren radius) is given by
(2) 
Note that these results are independent of the value of the cross section, , for absorption of the ionizing radiation, which is the process that actually produces the ionized gas. However, in order for this simple Strömgren analysis to be useful in determining , it is necessary for the thickness of the layer in which the gas switches from primarily ionized to primarily atomic to be thin: , which is generally satisfied for HII regions around OB stars (Strömgren, 1939). If this condition is not satisfied, the Strömgren analysis still provides the total volume emission measure, . The effects of dust absorption on Strömgren layers and spheres can be readily included (Petrosian et al., 1972).
The same type of analysis can be applied to the atomic gas in a PDR (Paper I). This is particularly useful because, whereas photoionization in an HII region is a relatively simple continuum absorption process, photodissociation is a line absorption process; sophisticated calculations (e.g., Draine & Bertoldi, 1996) are needed to follow the absorption as the lines become opaque and even overlap. Only a fraction of the line absorptions result in dissociation of the H molecule; the rest result in the emission of fluorescence radiation. ( is a weak function of frequency and of the level populations of the H molecules; the average value is in the range 0.12—Draine & Bertoldi, 1996—to 0.11—Browning et al., 2003, although it may be higher in regions within a PDR where a large fraction of the H molecules are in excited states.) The process that is analogous to radiative recombination is molecule formation, which occurs at a rate , where is the number density of hydrogen nuclei, is the H I fraction, and is the rate coefficient for H formation. In the Galaxy, H forms on dust grains and the rate coefficient is cm s (Draine & Bertoldi, 1996). In primordial, metalfree gas, H forms via the H process or the threebody process, which can be described similarly (§3). In the absence of dust absorption, a point source with a photon luminosity in the LymanWerner bands between 91.2 nm and 110.8 nm that is embedded in a medium of constant density, , can create a sphere of atomic gas with a total number of neutral hydrogen atoms
(3) 
In the absence of dust, the transition between the atomic and molecular gas is not sharp (Paper I); nonetheless, we can define a characteristic radius of the atomic gas by setting , obtaining
(4) 
which has the same form as the Strömgren radius. Similarly, a dissociating photon flux incident on the surface of a slab of gas can maintain an HI column density
(5) 
Setting gives the characteristic size of this layer,
(6) 
An important limitation of the Strömgren analysis is that the bandwidth of the dissociating radiation must be known in advance. As shown by Sternberg (1988), the effective bandwidth begins to shrink for metallicities , where dust absorption is sufficiently strong that the individual absorption lines no longer overlap. For the cases of greatest interest, however, and the results presented here are applicable.
1.3. Summary of Previous Papers in this Series
In Paper I, we gave an approximate determination of the thickness of the HI layer in a sphere irradiated from the outside by isotropic radiation. For such a sphere, the incident flux is unknown, since the flux due to radiation striking the surface from outside the sphere can be compensated in part by the flux due to radiation that passes through the sphere. It is therefore advantageous to characterize the radiation field by the ambient value of the mean density of dissociating photons. If we let be the energy density of photodissociating radiation, then is the corresponding photon density; let be the corresponding ambient value. We define the characteristic value of the thickness of the HI layer as with replaced by :
(7) 
For an opaque slab embedded in an isotropic radiation field, the energy density just outside the slab is half the ambient value and the flux is half that, or ; as a result, in this case. The importance of dust is measured by the characteristic dust optical depth,
(8) 
where is one of the two dimensionless parameters introduced in Paper I. Photodissociation regions are often characterized by the ratio , where is the ratio of the dissociating radiation field, , to the typical value in the Milky Way, which is photons cm s from Draine (1978). In terms of this ratio, is
(9) 
Physically, is the ratio of the number of LW photons absorbed by dust to the number absorbed by H molecules. Thus, dust has a significant impact on the structure of the HI layer for . The parameter is of order unity in the Milky Way (Krumholz et al., 2009a, hereafter Paper II; see eq. 13 below). Clumps in the PDRs discussed by Hollenbach & Tielens (1999) have , while the interclump gas in those PDRs has values about 100 times larger. Observe that for the case in which H forms on dust grains, both and are effectively measures of the total surface area of dust grains mixed into the gas; measures the area available for absorbing photons, and the area available for adsorbing hydrogen atoms. Thus their ratio should be independent of metallicity, dusttogas ratio, density, and most other quantities. For an opaque slab, we have seen that the incident flux , so that and
(10) 
keep in mind that these simple estimates are based on neglecting absorption of the incident FUV radiation by dust. For future reference, we note that this implies that the total HI column density in a dustfree slab (eq. 5) is . Finally, we note that the ratio of to the cloud radius, , is
(11) 
where , the dust optical depth associated with the cloud radius, is the other dimensionless parameter in Paper I.
As in the case of HII regions, a Strömgrentype analysis is particularly useful if the transition from atomic to molecular gas is sharp, . In Paper I we showed that this transition is indeed relatively sharp for , so that in that case it is possible to infer the size of the atomic gas by setting . For the case with dust, we make this approximation here; this is the principal approximation in our work. In the absence of dust, we can determine the total mass of H I (in the spherical case) or column of H I (in the slab case) without making this approximation.
The advance made in Paper II was to recognize that the density in the atomic gas in the transition region is set by the requirement that it be in twophase equilibrium (Wolfire et al., 2003). The latter authors showed that the minimum density of cold HI (the CNM) in galaxies is
(12) 
where and are the radiation field and metallicity normalized to the local Milky Way value. If the density is several times the minimum value ( with ), then
(13) 
depends primarily on metallicity, and weakly at that; here cm), etc. The surface density of H I in a galaxy, assuming that it is in a slab and neglecting dust absorption, is then
(14)  
for ; here g is the mass per hydrogen for a gas of cosmic abundances. Paper II showed that dust absorption reduces this by about a factor of 2. On the other hand, if the slab of gas is illuminated from both sides, then the total column of H I is increased by a factor 2. Paper II generalized these results to the case of finite clouds and demonstrated that the results were in good agreement with existing observations. Indeed, since the theory applies to individual cloud complexes, predictions were made on how the results would change as the resolution of the observations improves.
Krumholz et al. (2009b) applied the results of Papers I and II together with the earlier work of Krumholz & McKee (2005) to determine the star formation rate in galaxies as a function of metallicity. For low surface densities ( pc), the rate is dominated by the transition from H I to H; for intermediate columns ( pc pc, it is determined by the properties of giant molecular clouds; and for larger column densities it is determined by the pressure of the galactic ISM. Similarly, Krumholz et al. (2009c) showed how the results of Papers I and II can be used to explain the observed absence of damped Lyman systems with high column densities and metallicities. The treatment of the H I—H transition in both these papers relied on an analysis that is more broadly applicable than the one developed in Paper I. We provide that analysis here.
In §2, we present the general formalism for determining the thickness of the H I shielding layer for spherical molecular clouds. This formalism is based on a Strömgren–type analysis. In §3, we present the results for both slabs and spheres in the absence of dust; this is relevant to the study of primordial clouds. §4 gives the results including dust absorption, §5 shows how these results can be generalized to the case of clouds in which the atomic and molecular regions have differing densities, and the conclusions are given in §6.
2. Formalism
2.1. Basic Equations
Consider a spherical cloud of radius and uniform density exposed to an isotropic UV radiation field. We describe the radiation field in terms of the specific intensity in photon units, . Let be the opacity due to absorption by H molecules, where is the density of H molecules and is the H absorption cross section at frequency . Let be the opacity due to absorption by dust, where is the density of H nuclei and is the dust absorption cross section per H in the photodissociation part of the spectrum, from 91.2 nm to 110.8 nm; we ignore the weak frequency dependence of this cross section. The equation of transfer is then
(15) 
for the intensity in the direction. The opacity depends on position, but we have suppressed the argument for clarity. For a cloud bathed in a uniform radiation field of intensity , the intensity at a point inside the cloud is
(16) 
where the optical depths are proportional to the column density: the dust optical depth is
(17) 
and the H optical depth is
(18) 
where the range of integration extends along the ray from where the ray enters the cloud to the point in question and is the corresponding H column density along the ray.
We next integrate over the range of frequencies that can photodissociate H, , corresponding to wavelengths nm nm. The frequencyintegrated equation of transfer is
(19) 
where
(20) 
We include the argument in to emphasize that this opacity depends on angle through its dependence on ; also depends on position. The intensityweighted H opacity is then a function of ,
(21) 
The dependence of on is governed by two countervailing factors. On one hand, as increases more and more photons at frequencies near line centers are absorbed. Since decreases away from line center, the mean the opacity per H molecule is reduced. On the other hand, is also proportional to the density of H molecules , which increases with . The intensityweighted optical depth that determines the attenuation of is , so that
(22) 
In terms of the first three angular moments of the radiation field,
(23) 
the first two angular moments of the equation of radiative transfer are
(24)  
(25) 
Integrating over frequency gives
(26)  
(27) 
where
(28)  
We have assumed that the direction of is independent of frequency; as a result, is given by an expression like that for with replaced by and with a factor in the angular integrals. Writing out equations (26) and (27) in the case of spherical symmetry, we obtain
(29)  
(30) 
where and is the component of the radiation pressure tensor.
Now define the characteristic flux
(31) 
which is the flux required to photodissociate a layer of thickness , as discussed in §1. The solution is determined by the two quantities defined in Paper I,
(32)  
(33) 
where is the number density of dissociating photons far from the cloud. In terms of these parameters, we have
(34) 
We define the normalized flux and photon density by
(35) 
The normalized photon density far from the cloud is then
(36) 
We now use the assumption that the formation and destruction of H are in balance in order to determine the opacity ; this is the essence of the Strömgren argument. Observe that it is difficult to calculate explicitly, since it involves first a frequency integration over the cross section and the specific intensity, both of which are complicated functions of frequency due to line absorption, and then an average over angle. Balancing the rates of formation and destruction of H gives
(37)  
(38) 
where is the atomic fraction and is the rate coefficient for H formation. The principal approximation in our analysis is to assume that the transition from atomic to molecular gas is sharp, so that in the atomic gas. The balance equation then implies that
(39) 
Our seond approximation is to relate the fluxweighted opacity to by
(40) 
where is a constant. The results for the 1D case below verify that this is a good approximation. We expect to exceed unity since the rays that carry the flux are closer to the radial direction and have therefore undergone less attenuation and interact with a larger cross section, as mentioned below equation (21).
To close the moment equations of radiative transfer, we introduce a variable Eddington factor,
(41) 
Expressing equations (29) and (30) in dimensionless form and using equation (39) to eliminate , we obtain
(42)  
(43) 
where . These equations must be solved subject to the boundary conditions that both the photon density and the flux vanish at , the transition point from atomic to molecular gas; i.e., . The values of the photon density and flux at the surface are not known, except in the case in which the cloud is completely opaque.
2.2. Approximate Determination of the Eddington Factor
Many astrophysical problems are characterized by a central source of radiation in an opaque region, where the radiation is nearly isotropic, surrounded by a less opaque region in which the radiation is beamed outwards. Shu (1991) writes a general form for a relation to close the moment equations of radiative transfer as
(44) 
He points out that the choice makes a smooth transition from the isotropic case () to the beamed case (). Our problem is quite different, however: the radiation is isotropic at large distances from the cloud, but close to the cloud there is a “hole” in the radiation field since radiation cannot pass through the molecular core. We therefore need to develop a different approach to closure.
As remarked above, we close the moment equations by using variable Eddington factor, . Our approximation for at any point inside the cloud is based on representing the main qualititive feature of the radiation field—that the intensity decreases for rays at larger angles to the outward normal and vanishes for rays that have intersected the molecular core—by a step function decrease at some angle . (We adopt the convention in Paper I that is the angle between the outward normal and the direction from which the photon originates – see Figure 1.) We thus approximate the radiation field as being of constant intensity for (corresponding to , where ) and vanishing for (corresponding to ). We then obtain
(45) 
It remains to adopt an ansatz for . First consider the dustfree case. For points close to the molecular core , we assume that is the angle subtended by the molecular core, so that . In the opposite limit, in which the molecular core is assumed to be small () and the point is far from the core, we assume that a ray at intersects the midplane at a distance from the center of the cloud, so that . We combine these two limits by assuming that is the sum of the squares of the tangents in the two limiting cases,
(46) 
The corresponding cosine is , where the minus sign reflects the fact that the angle is in the second quadrant:
(47) 
The parameter is an eigenvalue of the problem, which is adjusted to give the correct net flux at the surface.
How does this compare with Shu’s (1991) closure approximation? Using the model that the radiation field is uniform for and zero for implies that the flux is . Solving for the closure factor in equation (44), we find . Since is negative, it follows that the closure factor is negative, which is qualitatively different than the standard case. The fact that implies that in equation (45) is always less than ; it is in the range .
Next, consider the case in which dust is included. For a dusty slab, the dissociating photon density at an optical depth from the surface is
(48) 
Evaluation of the integral gives
(49) 
where is an exponential integral and where the approximation is accurate to about 25%. The advantage of the approximation is that it shows that the photon density is as if the intensity along the normal filled a solid angle corresponding to with , so that
(50) 
For the spherical case, we take , where is the dust optical depth from the cloud center to edge. To join the case where dust dominates the opacity onto the case where molecular absorption dominates it, we use the simple ansatz
(51) 
with . Note that whereas is always negative in the dustfree case, it can be positive when dust is important.
2.3. Boundary Conditions for the Spherical Case
As in Paper I, we use the solution of the radiative transfer equation at the cloud surface to constrain the solution. At the surface (), photons at are unattenuated, those that traverse the cloud outside the molecular core, , with
(52) 
are attenuated, and those that strike the molecular core () are absorbed. Integration of equation (22) then gives
(53) 
We then introduce our next approximation: we replace , with (see eq. 28). Note that the error introduced by this approximation is limited, since the term in which it is made contributes at most half the total. Proceeding in the same fashion with the evaluation of the flux at the surface, we obtain
(54) 
where the factor 2 in front arises since the unattenuated intensity is whereas the unattenuated flux is . The optical depth is given by
(55) 
where is the value of of closest approach to the center of the cloud for a ray at an angle .
The equations we have now written down constitute a complete set that fully determines the solution. For a given choice of and , a solution consists of two unknown numbers, (or equivalently , since the two are related by equation 51) and , and two unknown functions, and . The two functions are constrained by the ordinary differential equations (42) and (43), while the two numbers are constrained by the algebraic equations and the consistency condition equation (54). We must choose a fixed value of , since we lack an equation to determine an additional parameter. We discuss how to choose this value in § 2.4, and we defer discussion of how to obtain the solution numerically once we have chosen until § 3.2 and § 4.2.
2.4. The SemiInfinite Slab Limit
The normalization we have used for the spherical case, which is based on , breaks down for a onedimensional, semiinfinite slab, which corresponds to the limit . We therefore normalize with respect to the ambient radiation field, by defining
(56) 
thus, and at the surface of the slab. The first two moments of the radiative transfer equation become:
(57)  
(58) 
where , with being zero at the surface and increasing inward. In the 1D limit (), the angle that enters the Eddington factor approaches zero (eq. 47). The Eddington factor thus becomes
(59) 
For small (i.e., close to the surface), as it should, since the radiation is approximately isotropic in the halfspace . For large , we have , and ; this too is as it should be, since highly extincted radiation is beamed so that .
As with the spherical case, we have now written down enough equations to close the system. The independent parameter that determines the solution is the dimensionless radiation intensity . For a given , a solution consists of two numbers, the ratio of the fluxmean opacity to the energymean opacity and the dust optical depth at which the molecular transition occurs, plus two functions, and . These quantities are constrained by two algebraic equations, , and the two ODEs (57) and (58). These ODEs are subject to the boundary conditions and . This set of constraints therefore fully specifies the system. We defer discussion of our solution algorithm to § 3.1 and 4.1. We pause here to point out an important distinction between the spherical and slab cases: for the spherical case, the angle at which we transition from zero flux to finite flux represents an extra parameter to be determined that is not present in the slab case. As a result, for the slab case we can use our constraint equations to determine , and we can then use the resulting value of in the spherical case. As we shall see, a single value of works quite well over a very broad range of radiation fields.
3. The DustFree Case: Primordial Clouds
Although the motivation for this series of papers was to understand the atomicmolecular transition in galaxies today, where H is formed on grains, our methodology applies to primordial gas clouds also, where H is formed via gas phase reactions. The effective values of for gasphase production of H are given in the Appendix. The results obtained here are also useful for comparison with the results for dusty clouds in the following section.
3.1. SemiInfinite Slab
In the limit that , equations (57) and (58) for the 1D case become
(60)  
(61) 
where we used the fact that (see below eq. 10) and defined
(62) 
which varies from 0 at the surface of the slab to 1 at the point that the slab becomes fully molecular. Note that we have not assumed anything about the spatial variation of the HI fraction, . As pointed out in §2.4, the boundary conditions are and ; we also pointed out that for small , the Eddington factor approaches . In the dustfree case, the use of at large optical depths in the H lines is justified by the fact that the photons that dominate the photodissociation are in the line wings; they are not highly extincted and therefore are not beamed. Integration of the equations then gives
(63)  
(64) 
Since must vanish at the same point that does (i.e., at ), it follows that
(65) 
We shall see in §4.1 that this remains a good approximation in the dusty case for all physically relevant values of .
3.2. Spherical Clouds
The problem of the dustfree spherical cloud does not have an exact solution with the Strömgren method because the net incident flux, , depends on the structure of the molecular hydrogen in the transition zone; some of the rays that strike the cloud will penetrate all the way through, reducing the incident flux by an amount that depends on the unknown spatial distribution of the H. We therefore consider two complementary approximations: In the first method, we assume that the gas is entirely atomic in the transition zone () and calculate the incident flux selfconsistently. In the second method, we allow for the fact that molecular gas exists in the transition zone (), but assume that the transition zone is thin enough that the incident flux has the value appropriate for an opaque cloud, . We shall see that the results overlap, giving a selfconsistent determination of the region of validity of our results.
In general, both and go to zero for dustfree spherical clouds, since both are proportional to the dust cross section . However, the ratio
(66) 
remains finite, and this becomes the independent parameter describing the normalized radiation intensity in the dustfree case. We anticipate that the cloud will be fully atomic at high values of , whereas the atomic gas will be confined to a thin shell at low .
3.2.1 Method 1:
With and , equation (42) reduces to
(67) 
which may be integrated analytically, subject to the boundary condition , to obtain
(68) 
Substituting this into equation (43) gives
(69) 
We can now obtain the relationship between and via the following algorithm:

Select a value for . Our algorithm will determine the corresponding .

Given , determine the function by integrating the ordinary differential equation (69) from to , using the boundary condition that .
Applying this algorithm yields the curve for versus shown in Figure 2. Note that in this approximation, goes to zero exactly for . In an exact calculation with , the gas would not be fully atomic at this point, but this result shows that clouds with are substantially atomic. The formula
(70) 
matches the numerical solution to within 11% in the range , and may be used in place of a numerical evaluation for most practical calculations.
3.2.2 Method 2:
If we retain the dependence on but assume that the incident flux is known, an analytic solution is possible. In this case, equation (42) reduces to
(71) 
The boundary condition becomes
(72) 
with the aid of equation (36). Integration of equation (71) gives
(73) 
where is the number of HI atoms outside and is the total number of H nuclei in the cloud. To relate this result to that from the first method, we define an effective value of such that
(74) 
where , the total number of HI atoms, is evaluated at . Equation (73) then implies
(75) 
This has the same form as equation (70), which is based on the first method, for . Since the first method properly accounts for the flux that penetrates through the cloud, we conclude that the total number of HI atoms in a spherical, dustfree cloud is
(76) 
for . Since Method 1 showed that the cloud becomes fully atomic for and both methods agree for , we conclude that the results of this section are valid everywhere except where the cloud is close to being fully atomic.
4. Dusty Clouds ()
4.1. SemiInfinite Slab
We now consider the case of clouds with dust, , beginning with the semiinfinite slab. As discussed in §1, the presence of dust sharpens the transition between atomic and molecular gas, enabling us to make the approximation that