Issue |
A&A
Volume 684, April 2024
|
|
---|---|---|
Article Number | A48 | |
Number of page(s) | 13 | |
Section | Astrophysical processes | |
DOI | https://doi.org/10.1051/0004-6361/202346533 | |
Published online | 03 April 2024 |
Effective equation of state of a radiatively cooling gas
Self-similar solution of spherical collapse
1
Department of Earth Sciences, National Taiwan Normal University, Taipei 116, Taiwan
e-mail: ynlee@ntnu.edu.tw
2
Center of Astronomy and Gravitation, National Taiwan Normal University, Taipei 116, Taiwan
3
Physics Division, National Center for Theoretical Sciences, Taipei 106, Taiwan
Received:
29
March
2023
Accepted:
10
January
2024
Context. The temperature of the interstellar medium (ISM) is governed by several physical processes, including radiative cooling, external UV/cosmic-ray heating, and mechanical work due to compression and expansion. In regimes where the dynamical effect is important, the temperature deviates from that derived by simply balancing the heating and cooling functions. This renders the expression of the gas energy evolution with a simple equation of state (EOS) less straightforward.
Aims. Given a cooling function, the behavior of the gas is subject to the combined effect of dynamical compression and radiative cooling. The goal of the present work is to derive the effective EOS of a collapsing gas within a full fluid solution.
Methods. We solved the Navier-Stokes equations with a parametric cooling term in spherical coordinates, and looked for a self-similar collapse solution.
Results. We present a solution that describes a cloud that is contracting while losing energy through radiation. This yields an effective EOS that can be generally applied to various ISM contexts, where the cooling function is available from first principles and is expressed as a power-law product of the density and temperature.
Conclusions. Our findings suggest that a radiatively cooling gas under self-gravitating collapse can easily manifest an effective polytropic EOS, even isothermal in many scenarios. The present model provides theoretical justification for the simplifying isothermal assumptions of simulations at various scales, and can also provide a more realistic thermal recipe without additional computation costs.
Key words: ISM: clouds / evolution
© The Authors 2024
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1. Introduction
In a thermally isolated system, that is, when the variation of gas internal energy results only from mechanical work without heat exchange, the gas adiabatic index γ describes the relation between pressure P and density ρ, such that P ∝ ργ. However, in astrophysical systems, cooling and heating mechanisms often come into play and act on a timescale that is sufficiently small compared to the dynamical timescale. The thermal behavior of the gas therefore deviates from adiabaticity. For convenience, an effective polytropic equation of state (EOS) is very often employed to avoid explicitly solving the energy equation in many contexts. An EOS following
where Γ is the polytropic index, is used to replace the energy conservation equation in a hydrodynamic system. This practice is widely accepted, given numerous supporting evidence, either observational (Myers 1978) or theoretical (Tarafdar et al. 1985; Boland & de Jong 1984; Larson 1973b).
The molecular interstellar medium (ISM) is known to be efficiently cooled due to the numerous molecular and atomic lines. The ISM cools from ∼1000 K at 1 cm−3 to ∼10 K at 104 cm−3, which can be described with a polytropic index Γ ∼ 0.73 (Larson 1985). In the density range between 104 and 1010 cm−3, the temperature is maintained at roughly 10 K. At higher densities, the dust component, which takes up only one percent of the mass fraction, starts to contribute significantly to the opacity, and the radiative heat loss becomes very inefficient. Therefore, many studies, without loss of generality, assume an isothermal behavior for the star-forming ISM, which only becomes adiabatic at high densities, corresponding to the first hydrostatic core formation.
The polytropic index is important for understanding the gas dynamics because it determines whether the energy balance allows gravitational collapse to occur. For collapse solution to exist, the critical value Γcrit = 4/3 for spherical geometry, and is equal to unity for cylindrical geometry. Given a polytropic index larger than the critical value, gravitational contraction leads to increased gas internal energy and results in a thermally supported structure.
In the diffuse ISM, the heating mostly comes from the cosmic rays or background UV radiation, which can be regarded as almost constant. The volumetric heating thus scales linearly with the density. On the other hand, cooling usually results from the deexcitation of the collisionally excited molecules, and therefore scales quadratically with the density and has a steep nonlinear dependance on the temperature. When the dynamical timescale of the system is significantly longer than the heating and cooling timescales, the gas temperature is simply determined through the balance between the two latter mechanisms. The present work considers the regime where the cooling timescale is greater than, or at least comparable to, the dynamic or thermal sound crossing time. Therefore, the cooling modifies the thermal behavior during the collapse, instead of reaching instantaneous equilibrium with other heating mechanisms.
The purpose of the present study is to derive the effective polytropic index, Γ, from a parametric cooling function that is a power-law product of density and temperature, Λ ∝ ρmTn. We solve the full hydrodynamic equations, adding a cooling term to the energy equation, and look for self-similar solutions. The solutions are then used to understand the apparent polytropic behavior of the collapsing gas. We successfully demonstrate that a gas element evolves along a polytrope with the material polytropic index, ΓM = 1 + (3/2 − m)/(n − 1), in the flat center of a cloud under homologous collapse, and with ΓM = 1 + (2 − m)/(n − 1) during the accretion phase onto the central mass singularity. However, when profiling the cloud as an ensemble, the thus derived apparent polytropic index, ΓA, has different values in different cloud regions. We clarify in this work that such a polytropic index measured with the temperature-density relation of instantaneous profiles does not reflect directly the thermodynamic properties of the cloud, and that their usage should not be confused.
2. Formalism
2.1. The physical system
The compressible Navier-Stokes equations with cooling writes
where t is time, u is the velocity vector, ϕ is the gravitational potential, G is the gravitational constant, ϵ is the specific internal energy, and Λ the volumetric cooling rate. The heat exchange includes heating and cooling, while here we only discuss the regime where cooling is more important than heating and there is a net energy loss. The gas follows the ideal gas law,
where kB is the Boltzmann constant, μ is the mean particle (molecular) weight, mp is the proton mass, and T is the temperature. By parameterizing the cooling function as
the equation set can be rewritten for a one-dimensional system in spherical symmetry:
2.2. Cooling models in the literature
The radiative cooling functions were calculated from first principles (e.g., Hollenbach & McKee 1979; Neufeld & Kaufman 1993; Neufeld et al. 1995; Omukai 2001), and were applied to and tested in many studies1. In general, the cooling due to collision-induced radiations scales quadratically with the density (m = 2) and linearly with the abundance of the cooling agent in question. On the other hand, with increased optical depth reaching local thermal equilibrium (LTE), cooling scales linearly with the density of the cooling agent (total density multiplied by the abundance; m = 1). However, we can still assume that the photons can escape efficiently once they exit the local region where they are emitted due to the large velocity gradients inside the cloud.
There have been many theoretical models of the cooling gas temperature, while a full collapse solution has not been proposed due to the complexity of the problem. For example, Goldsmith (2001) considered a static one-zone model and studied the effect of depletion on molecular cooling. They concluded that m ≈ 2 and 1.4 < n < 3.8 for the major coolants at low densities (between 102 and 107). Hunter (1969), Disney et al. (1969), and Larson (1973a) considered simplified, stationary models of balancing heating by cosmic rays and compression against cooling by atomic lines and dust emission, and assuming uniform density for the estimated optical depth. This simplified calculation already allowed the characteristic behavior of the cloud to be inferred. They all reached very similar conclusions, that is, isothermality during collapse, while they did not specify whether the center is in free fall or not. In the density range between 105 and 1011 cm−3, the cooling process is strongly coupled to dust grains and is very sensitive to temperature (∝T7; Larson 1973b). As we demonstrate in Sect. 3.1, the apparent EOS is indeed close to isothermality when the temperature dependence is steep.
2.3. Self-similar equations
Following the transformation of Murakami et al. (2004),
we look for a self-similar system. For self-similar solution to exist, the system must satisfy
We chose the similarity ansatz:
Here A, B ∈ ℛ+ are free parameters, the physical meanings of which are explained later. The time zero t = 0 corresponds to the epoch of formation of the central singular mass M0, or Ω0 in dimensionless form. Solutions should exist for t < 0 and for t > 0, corresponding to the collapse phase and accretion phase, respectively. The two phases are connected at ξ = ∞.
The dimensionless ODE set has the form
where the lower and upper signs correspond to t < 0 and t > 0, respectively. From the mass conservation law ∂M/∂t = −4πr2ρu, one derives
The constants
Alternatively,
When the absolute value of the cooling function is lower, the physical extent of the dimensional system is correspondingly larger. On the other hand, the central density of the system is fixed by the values of m and n. We note that a physical solution exists only for K1, Ω ∈ ℛ+, which signify density and mass positivity. Net cooling (Λ0 > 0) requires K2 to be of the same sign as a.
3. Results
We adopt the nomenclature by Whitworth & Summers (1985), and identify three regimes: early interior path (t < 0 and small ξ), exterior path (large ξ), and late path (t > 0 and small ξ). There is no analytical solution for this nonlinear system (except for one trivial solution at t < 0; see Appendix A), while we can derive asymptotic solutions in limiting conditions and perform numerical integration to obtain the full solutions.
3.1. Asymptotic behaviors of the self-similar system
The asymptotic solutions allow us to analytically evaluate the polytropic indices in the characteristic regimes of the collapsing cloud, and also serve as boundary conditions for numerical integration. The asymptotes are presented directly here and we refer to Appendix B for detailed derivations. We show the thus-derived apparent (instantaneous) polytropic indices
as functions of m and n for the three regimes in Fig. 1.
![]() |
Fig. 1. Apparent polytropic indices Γi, Γe, and Γl as functions of m and n. The critical values of 1, 4/3, and 5/3 are marked with white, yellow, and green contours, respectively. In regimes shown in red, the gas temperature increases with increasing density. In regimes shown in blue, the opposite happens. The values of Γ are shown irrespectively of the existence of a collapse solution, while the physical domain is shown in Fig. 3. |
3.1.1. Early interior path, homologous collapse at t < 0, ξ ≪ 1
Time t < 0 corresponds to the epoch before mass singularity formation at the center of the collapsing region. The central region is almost uniformly contracting, and has flat density and temperature profiles. The asymptotes have the form
where the coefficients are derived by inserting the asymptotic forms into the governing equation set. The effective polytropic index writes
Some values of Γi are given in Table 1 (see also left panel of Fig. 1). For reasonably large values of m and n, Γi is close to unity and can be either larger or smaller. We note that this polytropic index value is independent of the constants A and B.
Effective Γi of the early interior path with γ = 5/3.
3.1.2. Exterior path, stationary solution at ξ ≫ 1
The asymptotic solution at large distances or close to the mass singularity formation epoch has the form
where the coefficients are to be determined by numerical integration from ξ ≪ 1. This leads to an effective polytropic index
The value of Γe is displayed as a function of m and n in the middle panel of Fig. 1.
3.1.3. Late path, runaway collapse at t > 0, ξ ≪ 1
After the formation of the central point mass at t = 0, the density and velocity diverge when ξ → 0. The asymptotes are
where
The convergence toward the asymptotic solution can be slow when ω is significantly negative. The constants Ω0 and τ+ must be found by numerical integration from ξ ≫ 1. The right panel of Fig. 1 displays the resulting effective polytropic index as function of m and n, which writes
3.2. Full self-similar collapse solutions
The system was solved numerically by starting the integration at ξ ≪ 1 and t < 0. The value of K1 that allows the existence of a continuous solution was found with the shooting method. The solution with t < 0 was then mapped to the solution with t > 0 at ξ ≫ 1 and then integrated toward small ξ. Example profiles are shown in Fig. 2 for three selected combination of m and n. We show the profiles of dimensionless density, velocity, temperature, and enclosed mass. Solutions for t < 0 and t > 0 are plotted with dashed and solid lines, respectively. Some comparisons with numerical results in the literature are presented in Appendix C.
![]() |
Fig. 2. Complete self-similar solutions for three cases (from left to right): (m, n) = (1.3, 2.3),(1.5, 3), and (2, 5). The dimensionless density, velocity, temperature, and mass are shown in blue, orange, green, and red, respectively. The solutions for t < 0 and t > 0 are shown with dashed and solid lines respectively. |
3.3. Effective polytropic indices and parameter domain
The effective polytropic indices of the three asymptotes are shown in Fig. 1. We only show n > 1, where the cooling rate scales at least linearly with the internal energy. The red and blue areas correspond to Γ greater and less than unity, respectively. The characteristic values of Γ = 1, Γ = Γcrit = 4/3, and Γ = γ = 5/3 are marked with white, yellow, and green lines, respectively.
It is possible to infer that for m ≳ 1.5 the central flat region can have slightly decreasing temperature toward the central region with higher density (Γi < 1). The freefalling late path always has Γl greater than that of the exterior region Γe. This can be physically interpreted with the decreased dynamical time, leading to less efficient cooling. There are some special cases: the exterior path is isothermal when m = 1.5, and so is the late path with m = 2.
When cooling is inefficient (i.e., small m and small n values), some combinations can lead to Γe > γ. More specifically, this happens when
which corresponds to the region below the green line in the middle panel of Fig. 1. The cooling in the inner region with higher density and higher temperature is not sufficiently efficient compared to the outer region. This situation leads to an apparent polytropic index larger than the adiabatic index that seems to be unphysical at first glance. However, this is derived from the apparent EOS that describes an ensemble of gas instead of the thermal evolution of an individual fluid element. In other words, the actual thermal properties of the gas should be studied along a material line instead of using instantaneous snapshots of density and temperature profiles. This will be detailed further in Sect. 3.5. A similar situation can lead to apparent superadiabaticity Γl > γ, which happens when
This corresponds to the region below the green line in the right panel of Fig. 1.
Not all combinations of m and n allow the existence of a self-similar solution throughout the whole domain. Constraints such as negative velocity, net cooling, and positive mass, need to be satisfied for physical collapse solutions. These constraints are illustrated in the left panel of Fig. 3, where forbidden zones are color-shaded. The derivations for these conditions are detailed in Appendix D.1. The above-mentioned criteria for solutions at t < 0 and t > 0 are independent, while they must be satisfied simultaneously for solutions to be connected across time zero. If we compare Figs. 1 and 3, we see that superadiabaticity (Γl > γ) never occurs in the physically allowed domain of the late path. Γi > γ can occur in a narrow band, although this solution is not part of a full solution.
![]() |
Fig. 3. Parameter space analyses. Left: forbidden domain of self-similar solution existence. A complete solution only exists for combinations of m and n in the white area. Right: large-scale velocity behavior. The velocity converges to zero at infinite distance in the shaded zone. The parameters used for the example cases are marked with blue stars. |
On the other hand, apparent superadiabaticity can easily occur for the exterior path since there are no physical constraints, although such a solution cannot continue smoothly to ξ ≪ 1. The collapsing cloud is connected to larger scales, where the flow could be converging, static, or expanding. Therefore, conditions in the exterior path do not impose criteria for the existence of solutions. The cloud that has vanishing velocity at infinity can be interpreted as isolated, while that with nonzero velocity can be regarded as embedded in a converging flow. These regimes with vanishing velocity are illustrated in the right panel of Fig. 3 (see derivation in Appendix D.2). It is not possible to tell whether the large-scale velocity is converging or expanding without performing direct numerical simulation, while it is unlikely to have an expanding envelope if the center of the self-similar solution is contracting.
3.4. Evolution in physical coordinates
We transform the self-similar solution back to physical coordinates by assigning a value for Λ0. The evolution of the profiles are shown in Fig. 4 for m = 1.3 and n = 2.3, while the other two cases can be found in Appendix E. The behavior is basically very similar to what was suggested in previous isothermal and polytropic calculations. The inner flat region shrinks in size and increases in density before time zero, and the freefall onto the central singularity proceeds inside-out. The exterior profile is stationary.
![]() |
Fig. 4. Profiles in physical coordinate with Λ(106 cm−3, 10 K) = 10−26 erg cm−3 s−1, m = 1.3, and n = 2.3. Solutions for t < 0 and t > 0 are presented with dashed and solid lines, respectively. The increasing line thickness corresponds to increasing absolute time. |
To obtain some physical insights of the characteristic values, we express
where ρref and Tref are arbitrary reference values. The temporal evolution of central density is solely determined by the combination of m and n, while the temperature, mass, and the physical extent of the inner solution (early interior path and late path) decrease with more effective cooling (larger Λ0). Due to the power-law dependences on m and n, it is not straightforward to derive characteristic reference values. Nonetheless, under the limits that m is small and n is large, we recover the similarity normalization of the isothermal case (Shu 1977; Whitworth & Summers 1985), where
with . It is also evident from Fig. 1 that when n ≫ 1, all three ΓA values approach unity irrespective of m.
3.5. Material equation of state
The above-mentioned polytropic indices are relevant for single time frames of a collapsing sphere, which should be compared to resolved observations of a region or a snapshot from a simulation. The results represent an ensemble EOS, or apparent EOS, which is defined with Eq. (29). The effective EOS, or more precisely the temperature-density relation, is shown in Fig. 5 for the three example cases. The early interior path corresponds to a very short segment on these plots since both density and temperature are almost constant. The solution at low density and that at positive time tend very quickly to the power-law asymptotes (described in Sect. 3.1). In particular, the temperature approaches a constant value on the late path for m = 1.5 and on the exterior path for m = 2.
![]() |
Fig. 5. Temperature–density relation for three cases (from top to bottom): (m, n) = (1.3, 2.3), (1.5, 3), and (2, 5). Solutions for t < 0 and t > 0 are shown with dashed and solid lines, respectively. The asymptotic polytropic indices are (Γi, Γe, Γl) = (0.96, 1.15, 1.534), (0.93, 1.0, 1.25), and (0.87, 0.88, 1.0), respectively. |
On the other hand, in order to understand the thermodynamical behavior of one single fluid parcel, it is necessary to discuss the material EOS, which is in fact the true EOS,
where the material derivative
We can then derive
Inserting the asymptotic solutions previously derived, we can obtain the material polytropic index for the asymptotes. On the early interior path, a fluid parcel follows
on the T − ρ plane, which is parallel to the exterior path. On the exterior path, the material polytropic index
which falls somewhere between Γe and Γl. We can verify that ΓM, e < γ in the physical domain of m and n. On the late path, the fluid parcel evolves along the polytrope with index
which is identical to that of the ensemble EOS. It is sometimes not explicitly stated in the literature whether the apparent or material EOS is used. We discuss and clarify the usage in Appendix F.
The material equation of state can be inferred by tracing the density and temperature evolution at a fixed mass coordinate. Figure 6 shows some examples of the material evolution following enclosed mass M = 1, 10, 100, 1000, 10 000 M⊙ (gray lines with increasing thickness). Each gray line represents the time evolution of a fluid parcel at a given mass coordinate, in contrast to the spatial profiles at given instants presented with colored lines. When t < 0 the material EOSs at any mass fall on the same line, parallel to the exterior path, and move along ΓM, i. When a fluid parcel is swept by the inside-out collapse wave, its trajectory turns upward in the T − ρ plane. Finally when entering the freefall at t > 0, the material EOS has the same slope, ΓM, l, as that of the late path, but with varied absolute values for each mass coordinate.
![]() |
Fig. 6. Temperature–density relation snapshots of the three example cases shown in Figs. 4, E.1, and E.2. Solutions for t < 0 and t > 0 are presented with dashed and solid lines, respectively. The colored lines with increasing thickness correspond to increasing absolute time. The gray curves with increasing thickness trace the time evolution of fluid parcels at enclosed mass coordinate M = 1, 10, 100, 1000, 10 000 M⊙. |
4. Discussion
4.1. Characteristic timescales
There are several mechanisms that act simultaneously. According to the relative value of the timescales, the physical behavior of the system is thus different. To discuss the various regimes, we evaluate the relevant timescales:
-
Freefall timescale:
;
-
Cooling timescale:
;
-
Dynamical timescale:
;
-
Pressure timescale:
;
Nishida (1968) has discussed the collapse of a prestellar core by comparing the timescales. When the cooling or heating timescale is significantly smaller than the freefall timescale, the gas contracts or expands very quickly toward the cooling-heating equilibrium until the cooling timescale becomes comparable to the freefall timescale. Subsequently, the gas contracts along the equilibrium tcool ≈ tff.
Consistently found in our work, the collapsing gas always tends to reach a state where tdyn ∼ tcool. In a medium with relatively short cooling and heating timescales, this stationary state should occur very close to the equilibrium line between cooling and heating. This temperature is still slightly offset with respect to the equilibrium temperature, with a small net cooling that has a timescale comparable to the dynamical time.
Having obtained the collapse solutions, we can then calculate the relation between the timescales:
On the early interior path, tP/tff scales linearly with the radius, being small near the center, while all other timescales are of the same order of magnitude as the freefall timescale. The small pressure timescale maintains the flat pressure and density profiles of the central region. On the exterior path, all the above-mentioned timescales are maintained with constant ratios. The system therefore reaches a stationary state. On the late path,
which increases with decreasing radius (cf. the physical parameter domain in Fig. 3). The thermal pressure becomes unimportant as the freefall proceeds. Both tcool and tdyn remain in a constant ratio with the freefall time. The difference among the three regimes comes mainly from the different scaling behaviors of the pressure timescale.
4.2. Convective stability
The Schwarzschild instability criterion can be translated into2
The interpretation is straightforward: if the environmental temperature gradient is too steep compared to the material variation (governed by radiation in this context), a perturbation cannot be restored by buoyancy and will lead to convection. This is valid when the local dynamical timescale of the perturbation is larger than the radiation timescale. Otherwise, radiation will not have enough time to act, and the gas polytropic index ΓM should be replaced with γ. In the physical regime of m and n combinations (cf. Fig. 3), the material polytropic index in always larger than the apparent polytropic index, as can be seen from the expressions of ΓM in Sect. 3.5 and Fig. 1. If we look at the points where a gray line intercepts a colored line in Fig. 6, their slopes indicate the material and apparent (environment) polytropic indices, respectively. The slope of the gray line is always greater than the colored line, implying that the stability against convection and a perturbed gas parcel always tends to move back to its original position.
4.3. New perspective on the collapse criteria
From energy arguments it can be demonstrated that the effective polytropic index, Γ, of a gas must be smaller than a critical value for collapse to happen. However, instead of defining a critical value of Γ for collapse to occur, we demonstrate in this work that it is physically more relevant to define a cooling function that allows the existence of a collapse solution. The mass positivity condition (Eq. (D.7)) coincides with Γe < 4/3, while it is evident that this collapsing envelope could be connected to an inner freefalling region with Γl > 4/3. The existence of runaway collapse only requires that Γl < γ = 5/3, not 4/3.
5. Conclusions
In this study, we discussed a spherically collapsing gas cloud that follows a parametric cooling law Λ ∝ ρmTn. The effective EOS of a cooling and self-gravitating gas can be divided into different regimes according to the relative values of the dynamical and thermal timescales. A self-similar solution is found for such a system, which consists of a central region that is either flat or freefalling before and after the mass singularity formation, and an outer part that is undergoing stationary contraction. This is very similar to other existing collapse solutions.
The distinct feature of the current work is that the effective gas EOS, for the three regimes of the solution, can be directly derived given a parametric cooling function. We derived both the ensemble (apparent) EOS of a collapsing profile snapshot and the material (real) EOS of the fluid parcels within. Instead of examining the effective Γ as the collapse criteria, the cooling function should be used to determine whether the energy can be efficiently evacuated during the collapse process.
The spherically symmetric configuration was discussed in the present manuscript. Nonetheless, on the scales of cloud structure formation, cylindrical collapse, corresponding to filament formation, and plane-parallel geometry, corresponding to colliding flows, should also be of physical interest and worth further investigation. We will continue the calculations in the upcoming work.
The sign of Eq. (66) in Murakami et al. (2004) should be reversed.
Acknowledgments
The author would like to thank the anonymous referee for constructive comments that helped to make the manuscript more concise and comprehensive. This work received financial support of the NSTC (project number 111-2636-M-003-002), and the MOE Yushan Young Fellow program. This work was initiated during the JSPS fellowship in 2020 at Nagoya University.
References
- Abel, T., Anninos, P., Zhang, Y., & Norman, M. L. 1997, New Astron., 2, 181 [NASA ADS] [CrossRef] [Google Scholar]
- Boland, W., & de Jong, T. 1984, A&A, 134, 87 [NASA ADS] [Google Scholar]
- Clark, P. C., Glover, S. C. O., Klessen, R. S., & Bromm, V. 2011, ApJ, 727, 110 [NASA ADS] [CrossRef] [Google Scholar]
- Disney, M. J., McNally, D., & Wright, A. E. 1969, MNRAS, 146, 123 [NASA ADS] [CrossRef] [Google Scholar]
- Flower, D. R., Le Bourlot, J., & Pineau des Forêts, G., & Roueff, E., 2000, MNRAS, 314, 753 [NASA ADS] [CrossRef] [Google Scholar]
- Galli, D., & Palla, F. 1998, A&A, 335, 403 [NASA ADS] [Google Scholar]
- Goldsmith, P. F. 2001, ApJ, 557, 736 [Google Scholar]
- Greif, T. H., Springel, V., White, S. D. M., et al. 2011, ApJ, 737, 75 [NASA ADS] [CrossRef] [Google Scholar]
- Hollenbach, D., & McKee, C. F. 1979, ApJS, 41, 555 [Google Scholar]
- Hunter, C. 1977, ApJ, 218, 834 [NASA ADS] [CrossRef] [Google Scholar]
- Hunter, J. H. Jr. 1969, MNRAS, 142, 473 [NASA ADS] [CrossRef] [Google Scholar]
- Larson, R. B. 1973a, ARA&A, 11, 219 [NASA ADS] [CrossRef] [Google Scholar]
- Larson, R. B. 1973b, Fund. Cosmic. Phys., 1, 1 [NASA ADS] [Google Scholar]
- Larson, R. B. 1985, MNRAS, 214, 379 [NASA ADS] [Google Scholar]
- Lesaffre, P., Belloche, A., Chièze, J. P., & André, P. 2005, A&A, 443, 961 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- McGreer, I. D., & Bryan, G. L. 2008, ApJ, 685, 8 [CrossRef] [Google Scholar]
- Murakami, M., Nishihara, K., & Hanawa, T. 2004, ApJ, 607, 879 [NASA ADS] [CrossRef] [Google Scholar]
- Myers, P. C. 1978, ApJ, 225, 380 [NASA ADS] [CrossRef] [Google Scholar]
- Neufeld, D. A., & Kaufman, M. J. 1993, ApJ, 418, 263 [NASA ADS] [CrossRef] [Google Scholar]
- Neufeld, D. A., Lepp, S., & Melnick, G. J. 1995, ApJS, 100, 132 [Google Scholar]
- Nishida, M. 1968, PASJ, 20, 162 [NASA ADS] [Google Scholar]
- Omukai, K. 2001, ApJ, 546, 635 [NASA ADS] [CrossRef] [Google Scholar]
- O’Shea, B. W., & Norman, M. L. 2007, ApJ, 654, 66 [CrossRef] [Google Scholar]
- Pagani, L., Lesaffre, P., Jorfi, M., et al. 2013, A&A, 551, A38 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
- Shu, F. H. 1977, ApJ, 214, 488 [Google Scholar]
- Suto, Y., & Silk, J. 1988, ApJ, 326, 527 [CrossRef] [Google Scholar]
- Tarafdar, S. P., Prasad, S. S., Huntress, W. T. Jr., Villere, K. R., & Black, D. C. 1985, ApJ, 289, 220 [NASA ADS] [CrossRef] [Google Scholar]
- Whitworth, A., & Summers, D. 1985, MNRAS, 214, 1 [Google Scholar]
- Yahil, A. 1983, ApJ, 265, 1047 [NASA ADS] [CrossRef] [Google Scholar]
- Yoshida, N., Omukai, K., Hernquist, L., & Abel, T. 2006, ApJ, 652, 6 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Trivial solution for t < 0
A trivial solution can be identified when K1 = K0 = d2/6:
This solution describes the homologous collapse of a flat structure that extends to infinity, and is valid only for t < 0. For this solution Γ is undetermined since both g and τ are constant. It is of less interest in the current context, and thus we do not discuss it in detail.
Appendix B: Derivation of the asymptotic solutions and the sonic point
B.1. Early interior path, t < 0, ξ ≪ 1
We presume the polynomial series for the centrally flat profiles before the mass singularity formation:
Since A and B are free parameters, we can thus choose a scaling such that the leading terms of g and τ are both unity. Physical boundary conditions require
which implies the absence of singularity at the origin. Substituting into the equations, it is easily found that
In short, all odd-order terms of g and τ and even-order terms of v should vanish for symmetry reasons. The asymptotic solutions simplify to
Inserting into the governing equation set, we obtain
It should be noted that K1 determines the sign and amplitude of g2. Fixing the temperature leading term τ0 = 1 is equivalent to fixing the value of K2. The variation of τ0 and K2 leads to rescaling of the system, which corresponds to exactly the same solution when transformed back into physical units.
The effective polytropic index can be derived by inserting the asymptotes into Eq. (29), which yields
There are two conditions that give Γi = γ. One is when (n − m)/(n − 1) = γ, and the other is when (c − (γ − 1)d) = 0. The conditions translate respectively to
For γ = 5/3, the first condition gives 3m + 2n = 5 and the second 6m + 4n = 13. The denominator vanishes on the curve
where the sign of Γi changes from positive to negative infinity on the two sides. We note that this polytropic index value is independent of the choice of A and B. When K1 = K0 = d2/6, the trivial solution is recovered.
B.2. Exterior path, ξ ≫ 1
The asymptotic solution at large distances, or close to the mass singularity formation epoch, can be expanded in power-law series
where the highest order constants are to be determined by numerical integration if one wishes to connect the solutions smoothly to ξ ≪ 1. The lower order terms follow
B.3. Alternative asymptotes for ξ ≫ 1
The asymptotic solution at large distances or close to the mass singularity formation time has an alternative form with nonvanishing velocity
where
This asymptotic solution exists for opposite domains of m-n combination for t → 0− and t → 0+. Therefore, the solutions cannot be connected smoothly across time zero.
B.4. The sonic point (surface)
When writing the equation set as linear combinations of the derivatives, some algebraic manipulation yields
where
is the determinant of the system, and
The criterion Δ(g, v, τ) = 0 defines a sonic surface. The solutions with t > 0 does not pass through through the sonic surface. On the other hand, for the solutions with t < 0 to pass smoothly through the sonic surface, the conditions
must be satisfied simultaneously on the sonic point, where the solution intercepts the sonic surface, and
also holds through linear dependance.
There are two types of solution that can pass through the sonic point smoothly. The first type is the trivial solution where the density and the temperature both remain constant (Appendix A). The other valid value of K1 must be inferred by performing the shooting method toward the sonic point with numerical integration. Depending on m and n, sometimes K1 < K0 and leads to increasing density with increasing distance from the center. In this case the exterior solution crosses the sonic surface again, while not satisfying Eq. (B.33). If we somehow relax the continuity condition at the sonic point, it is then still possible to find a pair of early interior path and exterior path solutions that are connected through a shock. Due to its complexity, we do not study the detailed properties of the sonic point in the present work.
B.5. Time reversibility
The gas might be subject to heating at low density and low temperature, while the current study does not consider the heating mechanism. The effect of heating can be easily obtained by changing the sign of the cooling term. With a transformation of the time coordinate t → −t, this results in an expansion solution that has exactly the same effective polytropic indices.
For the monoatomic gas with γ = 5/3 considered in this study, a self-similar collapse solution can only exist when cooling is considered; we recall that the critical value in spherical geometry is 4/3. If heating is present instead, the thermal energy must exceed the gravitational energy at some moment, and thus the collapse cannot continue infinitely. On the other hand, the present self-similar solution remains valid if the sign of the entire system is reversed.
Appendix C: Thermal behavior of gas in simulations with cooling
We compared our results to the numerical results in the literature. Since the current model is a generic model that only requires knowing the cooling function, it can therefore be easily applied to various cooling mechanisms.
C.1. Regular ISM
Tarafdar et al. (1985) simulated the thermal evolution of a collapsing cloud with a cooling chemical network. Lesaffre et al. (2005) also performed 1D spherical collapse simulation coupled to a chemical network and radiation, while solving the dust temperature explicitly. The evolution of density, temperature, and velocity profiles obtained in these two studies is very similar to our self-similar solution with m ≈ 2 before time zero (see Fig. 2 in Tarafdar et al. 1985). The central region has almost flat density and temperature profiles and undergoes homologous collapse. The size of this central region shrinks and the temperature decreases with time, while the envelope is stationary and the velocity vanishes at infinity.
Lesaffre et al. (2005) found a density profile with slope shallower than −2 and a decreasing temperature profile toward the center when cooling is considered. This is exactly reproduced with the current solution with m > 3/2. However, our solution is not applicable to their results when the central density reaches ∼109cm−3, where the cooling becomes optically thick. Since some of the major chemical reactions in a star-forming core have comparable timescales to that of the collapse, correctly calculating the temperature is important for out-of-equilibrium chemistry (e.g., Pagani et al. 2013). However, a full 3D hydrodynamical simulation with a complete chemical network and radiative transfer is still numerically challenging at present.
C.2. Primordial ISM
Yoshida et al. (2006) used the cooling model for primordial gas (Abel et al. 1997; Galli & Palla 1998), with cooling dominated by hydrogen molecules. Molecular hydrogen being the major coolant, m ≈ 2 for the optically thin low density gas and m ≈ 1 at higher densities, and a steep temperature dependence gives n > 1. Our model explains why the effective Γ switches from below unity to above at density ≈103 − 104 cm−3 in their Fig. 3. At densities above ≈1014 cm−3, the cooling is dominated by the collision-induced emission and follows m ≈ 1, giving Γ > 1.
At densities of around 109 − 1011 cm−3, the rapidly increasing molecular fraction leads to values of m that are significantly larger than unity since the cooling rate depends on the H2 number density rather than the total number density, which is dominated by hydrogen atoms. Also shown in the simulations by Clark et al. (2011, their Fig. 11) is that the molecular fraction increases roughly linearly with density. Consequently, m ≈ 2, and the effective Γ is slightly lower than unity.
Greif et al. (2011) used the molecular hydrogen cooling with slightly different prescriptions. They performed cosmological simulations and analyzed some minihaloes formed within. At very low density < 1 cm−3 the gas heats up adiabatically. After the cooling becomes effective, the behavior at density < 108 cm−3 is very similar to that in Yoshida et al. (2006). Interestingly, all minihaloes do not behave identically, and the prolonged HD cooling maintains the temperature at ≈200K in the density range 104 − 108 cm−3. This is probably the result of the combination between HD cooling with intrinsically m ≈ 1.5 at density and temperature slightly below values where LTE can be reached (Flower & Le Bourlot 2000) and the HD enhancement at low temperature (McGreer & Bryan 2008).
O’Shea & Norman (2007) performed cosmological simulations for an ensemble of haloes. They found that at higher redshifts, a higher environmental temperature leads to a higher H2 fraction. This results in haloes with lower central temperature and lower mass accretion rates at comparable central density. These findings are exactly accounted for by our solution, as shown in Eqs. (15) and (16), due to larger Λ0 value in Eq. (27). Our model also yields a smaller central density plateau under the condition shown in Eq. (14), which is seen in their Fig. 11.
Appendix D: Parameter space
D.1. Physical regime
On the early interior path, the contracting condition for the center, v1/a < 0, is automatically satisfied. The cooling condition requires K2/a > 0 (cf. Eq. (26)). This is equivalent to
or alternatively,
The forbidden parameter space is shaded in blue with dots in the left panel of Fig. 3.
On the late path, the existence of the runaway solution requires the positivity of Eq. (37), which means that
which yields
or alternatively,
This forbidden parameter space is shaded in green with circles in Fig. 3 (left panel). Another mass positivity criterion from Eqs. (24), (34), and (35) for a collapsing late path solution requires
which translates into
The forbidden parameter space is shaded in pink with slashes in Fig. 3 (left panel). The late path solution should satisfy Eqs. (D.4) and (D.7) at the same time.
D.2. Large-scale behavior
There are no physical constraints on the cooling parameters at large scales. However, the parameter space is separated into regimes, where the velocity vanishes at infinity or not. Vanishing velocity at infinite distance requires c < 0, which yields
which can be interpreted as isolated or accreting clouds. Depending on the large-scale behavior, the cloud can be interpreted to be isolated or connected to a converging flow.
Appendix E: Evolution in physical coordinates
Profiles of density, velocity, temperature, and enclosed mass are plotted at several time frames. They are presented in Figs. E.1 and E.2 for two cases: m = 1.5 and n = 3, and m = 2 and n = 5.
![]() |
Fig. E.1. Profiles in physical coordinates with Λ(106cm−3, 10K) = 10−26 erg cm−3s−1, m = 1.5, and n = 3. Solutions for t < 0 and t > 0 are shown with dashed and solid lines, respectively. The increasing line thickness corresponds to increasing absolute time. |
![]() |
Fig. E.2. Profiles in physical coordinates with Λ(106cm−3, 10K) = 10−26 erg cm−3s−1, m = 2, and n = 5. Solutions for t < 0 and t > 0 are shown with dashed and solid lines, respectively. The increasing line thickness corresponds to increasing absolute time. |
Appendix F: Ambiguity in the usage of the polytropic Γ
As an extension to the self-similar solution of an isothermal gas (Shu 1977; Hunter 1977; Whitworth & Summers 1985), Yahil (1983) studied the self-similar solution of a polytropic gas with3
where K is a constant. This EOS fixes the trajectory of all fluid elements onto a single straight line on the log T − log ρ plane (i.e., Γ ≡ ΓA = ΓM). Suto & Silk (1988) further generalized this calculation to a gas that follows
where K is a function of time. This solution describes an ensemble of gas that falls on a straight line on the log T − log ρ plane with fixed slope, while the exact position of this line evolves with time. A constant K value being a special case, the actual trajectory of gas on the log T − log ρ plane is generally different from the gas ensemble distribution: Γ ≡ ΓA ≠ ΓM. In some of their solutions, the gas has to be either cooled or heated in different regions to maintain this overall polytropic appearance.
Murakami et al. (2004) discussed a similar problem to the present study, where the gas has an intrinsic adiabatic index γ and is cooled explicitly with a prescribed cooling function. They used the radiative diffusion as the relevant cooling mechanism in the optically thick regime. They defined
is the Péclet number. This is the definition of the material polytropic index, which describes the thermal behavior of a fluid element that heats up during compression, while some of the energy is lost through radiative diffusion. The thus derived polytropic index, Γ ≡ ΓM, can be compared to, for example, the evolution of the central region of the first Larson core. While not explicitly stated in their work, they also derived the apparent polytropic index, which is smaller than the material polytropic index.
In the present context, we studied the radiative cooling process in the optically thin regime, with a parameterization as a power-law function of density and temperature. This regime is more relevant for the collapse and structure formation at molecular cloud scales. We demonstrated that the instantaneous temperature and density profiles of a collapsing cloud can exhibit an apparent polytropic index different from that derived from the temporal evolution of the overall (or local) cloud properties.
All Tables
All Figures
![]() |
Fig. 1. Apparent polytropic indices Γi, Γe, and Γl as functions of m and n. The critical values of 1, 4/3, and 5/3 are marked with white, yellow, and green contours, respectively. In regimes shown in red, the gas temperature increases with increasing density. In regimes shown in blue, the opposite happens. The values of Γ are shown irrespectively of the existence of a collapse solution, while the physical domain is shown in Fig. 3. |
In the text |
![]() |
Fig. 2. Complete self-similar solutions for three cases (from left to right): (m, n) = (1.3, 2.3),(1.5, 3), and (2, 5). The dimensionless density, velocity, temperature, and mass are shown in blue, orange, green, and red, respectively. The solutions for t < 0 and t > 0 are shown with dashed and solid lines respectively. |
In the text |
![]() |
Fig. 3. Parameter space analyses. Left: forbidden domain of self-similar solution existence. A complete solution only exists for combinations of m and n in the white area. Right: large-scale velocity behavior. The velocity converges to zero at infinite distance in the shaded zone. The parameters used for the example cases are marked with blue stars. |
In the text |
![]() |
Fig. 4. Profiles in physical coordinate with Λ(106 cm−3, 10 K) = 10−26 erg cm−3 s−1, m = 1.3, and n = 2.3. Solutions for t < 0 and t > 0 are presented with dashed and solid lines, respectively. The increasing line thickness corresponds to increasing absolute time. |
In the text |
![]() |
Fig. 5. Temperature–density relation for three cases (from top to bottom): (m, n) = (1.3, 2.3), (1.5, 3), and (2, 5). Solutions for t < 0 and t > 0 are shown with dashed and solid lines, respectively. The asymptotic polytropic indices are (Γi, Γe, Γl) = (0.96, 1.15, 1.534), (0.93, 1.0, 1.25), and (0.87, 0.88, 1.0), respectively. |
In the text |
![]() |
Fig. 6. Temperature–density relation snapshots of the three example cases shown in Figs. 4, E.1, and E.2. Solutions for t < 0 and t > 0 are presented with dashed and solid lines, respectively. The colored lines with increasing thickness correspond to increasing absolute time. The gray curves with increasing thickness trace the time evolution of fluid parcels at enclosed mass coordinate M = 1, 10, 100, 1000, 10 000 M⊙. |
In the text |
![]() |
Fig. E.1. Profiles in physical coordinates with Λ(106cm−3, 10K) = 10−26 erg cm−3s−1, m = 1.5, and n = 3. Solutions for t < 0 and t > 0 are shown with dashed and solid lines, respectively. The increasing line thickness corresponds to increasing absolute time. |
In the text |
![]() |
Fig. E.2. Profiles in physical coordinates with Λ(106cm−3, 10K) = 10−26 erg cm−3s−1, m = 2, and n = 5. Solutions for t < 0 and t > 0 are shown with dashed and solid lines, respectively. The increasing line thickness corresponds to increasing absolute time. |
In the text |
Current usage metrics show cumulative count of Article Views (full-text 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 48-96 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.