Introduction
The concept of the Galactic Habitable Zone (GHZ), as a part of the Milky Way that is likely to contain habitable planets, is a convenient platform for exploring the habitability of the Galaxy and in the past few years much work has been completed on the subject (Gonzalez et al. Reference Gonzalez, Brownlee and Ward2001; for a fine review see Gonzalez Reference Gonzalez2005). The exact boundaries (and their temporal evolution) of the GHZ are still uncertain mainly due to the lack of observational data (thus far Earth is the only observed planet that hosts life), as well as the dearth of sophisticated theoretical models detailing preconditions for life. Although the basic physical processes that drive planet formation and establish (non-)habitable conditions are clear, there is still a lot of work to be done, preferentially on calibrating computer simulations with the incoming observational data, e.g., those of the space-based missions like Gaia, Terrestrial Planet Finder (TPF) and others. The uncertainties in global Galactic processes and parameters, e.g., star formation rate, initial mass function, and rates of the global regulation mechanisms (GRMs), imply an uncertainty in the parameters of GHZ (namely, the number of Earth-like planetary systems and the probability that some of them host intelligent life). In addition, there is very little knowledge on the origin of intelligence (noogenesis) timescales. The aim of this paper is to present a simple Monte Carlo model of the GRMs in the GHZ and to estimate their influence on the number and temporal evolution of habitats with increased astrobiological complexity. This analysis is a part of the broader work on quantifying and simulating the GHZ. The simple Monte Carlo toy model of the GHZ, presented here, is a refined model from Vukotić & Ćirković (Reference Vukotić and Ćirković2007), henceforth referred to as Paper I. The applicability of the present class of models to such key astrobiological puzzles as Fermi's paradox and Carter's anthropic argument has been presented in other companion papers (Ćirković & Vukotić Reference Ćirković and Vukotić2008; Ćirković et al. Reference Ćirković2009). In contrast, we intend to focus here on the internal mechanics of the model, with particular emphasis on the input data characterizing the GHZ.
Relevant model values
Boundaries of the GHZ
With the lack of complete understanding of the processes that shape its inner and outer boundary, the spatial and temporal extent of the GHZ is still questionable and estimates in the literature vary. Prantzos (Reference Prantzos2008) suggests that the GHZ may encompass the whole thin disk, while the more conservative approach of Lineweaver et al. (Reference Lineweaver, Fenner and Gibson2004) gives the highest probability for Earth-like planets to reside around the stars that formed between 8 and 4 billion years ago and in the annular region between galactocentric radii of 7 and 9 kpc that slowly widens with time due to chemical and dynamical evolution. The outer limit to the GHZ is set by the metallicity gradient of the Milky Way, measured to be ▽Z=0.09 dex kpc−1 (e.g., Tadross Reference Tadross2003). From the purely dynamical (collision) constraints, Ćirković (Reference Ćirković, Knezevic and Milani2005) obtains the value of 1.3 kpc for the inner boundary of the annular ring, but it seems to be too low in comparison to other estimates. The work of Pena-Cabrera & Durand-Manterola (Reference Pena-Cabrera and Durand-Manterola2004) indicates the inner boundary at 4 kpc, on the basis of the assumption that metal-rich planets in the inner Galaxy would grow too fast to be inhabitable; this assumption, however, is highly questionable since another kind of habitable planet could plausibly form under high-Z conditions (Léger et al. Reference Léger2004). Another process that could play some role is the photo-evaporation of protoplanetary disks suppressing terrestrial planet formation, considered in a different context by Adams et al. (Reference Adams, Hollenbach, Laughlin and Gorti2004); modelling of this process in the GHZ has not been done so far.
Given the fact that a certain amount of metallicity is required to build up terrestrial planets, the GHZ should favour thin disk over thick disk stellar populations. Reddy (Reference Reddy, Vazdekis and Peletier2007) argues for the existence of an intermediate population with the thin disk abundances and thick disk kinematics. Bensby et al. (Reference Bensby, Zenn, Oey and Feltzing2007) suggests that thick disk stars have experienced early metallicity enrichment up to the solar metallicity value, ending about 8–9 Gyrs ago. It is generally believed that the thick disk stellar population is mostly older than 10 Gyr, thus making the existence of Sun-like stars unlikely. While the possibility that the thick disk may contain rocky planets around Sun-like stars cannot be rejected at present, we shall not consider the thick disk as a part of the GHZ in the course of the present work. The thin disk exhibited more intense star formation than the thick disk during the past 10 Gyr and is probably a far better target to search for Earth-like life prerequisites.
The number of potential habitats
Table 1 illustrates the mass variations of the GHZ stellar component for the Milky Way synthetic model from Robin et al. (Reference Robin, Reylé, Derrière and Picaud2003). Given the uncertainties in the GHZ boundaries and for simplicity, the corrections for disk warping and flaring were not taken into account. The GHZ is modelled as an annular ring between variable values for the inner (R inn) and outer (R out) radii with the scale height of the thin disk population taken from Robin et al. (Reference Robin, Reylé, Derrière and Picaud2003). The cut-off value for the heights above/below the Galactic plane is taken to be 400 pc. The fraction of the total stellar mass between 3 and 15 kpc is given in brackets. The total stellar mass is more sensitive to the value of R inn, given the considered density profiles, and for R inn>5 kpc it is approximately 109 solar masses. The simple model of the Galactic disk from Naab & Ostriker (Reference Naab and Ostriker2006), in accordance with the results of other authors, gives 2.7×1010 solar masses for the total mass of the disk stellar component out to 26 kpc – this work gives 1.2×1010 solar masses for the thin disk population between 3 and 15 kpc (Table 1).
Table 1. The mass of the GHZ stellar component in 109 solar mass
Rows – values for R inn; columns – values for R out. The value for the fraction of the mass between the minimum given value for R inn and maximum given value for R out is presented in brackets.
Circumstellar habitable zone host stars (CHZHS) can range from early F to late K spectral classes (Kasting et al. Reference Kasting, Whitmire and Reynolds1993), i.e., approximately 0.6–1.3 solar massFootnote 1. The CHZHS population is dominated with low-mass stars and their main sequence life time is larger than the present day age of the thin disk. Kroupa (Reference Kroupa2002) calculated that 29% of the total stellar mass is in stars with masses from 0.5 to 1 solar mass. Using the same mass function, a simple calculation shows that the total stellar mass in the GHZ given in units of solar mass, multiplied by approximately 0.3, will give a rough estimate of the number of CHZHS in the GHZ. The scenarios of planet formation where ‘hot Jupiters’ destroy nascent Earths, combined with the (non-)existence of circumstellar habitable zones around multiple and single stars, should reduce this fractionFootnote 2.
From the above considerations, the present-day number of CHZHS in the Milky Way of around 109 is mostly sensitive to the lower end of the chosen stellar mass interval and the inner boundary of the GHZ. The temporal evolution of these habitats, in the first approximation, reflects the star formation history of the thin disk.
Global regulation mechanisms
GRMs are manifested as events that can partially or completely sterilize habitable planets. Although devastating, they may be a crucial factor for more simple forms of life to advance into becoming more complex life forms. As they damage planetary biospheres, more ecological niches open up and new life forms evolve more quickly from the heritage of the old ones (Dar et al. Reference Dar, Laor and Shaviv1998). After ‘biological tempering’, simple forms of life should advance towards greater astrobiological complexity.
Gamma-ray bursts and supernovae
Although not completely understood, gamma-ray bursts (GRBs) are considered as very luminous events that originate when the cores of massive stellar objects collapse (‘hypernovae’) or during mergers of massive stellar remnants. Since the 1990s there have been many studies concerning the astrobiological impact of GRBs (Thomas Reference Thomas2009 and references therein) and it is argued that such events might have initiated some of the mass extinctions in the history of life on Earth. Scalo & Wheeler (Reference Scalo and Wheeler2002) calculated that GRBs can be lethal for eukaryote life forms to distances greater than 10 kpc, and that GRBs anywhere in the Milky Way can be significant. They also estimated, by modelling the GRBs as a two-dimensional Galactic plane Poisson process, that at least 1000 biologically significant ultraviolet (UV) irradiations have occurred in the 4.4 Gyr long Earth's history. They argued that GRBs are likely to be collimated events (e.g., see the review by Zhang & Mészáros Reference Zhang and Mészáros2004). Using the rate model from the next section, and a nearby GRBs rate of 0.025 Myr−1 galaxy−1, an estimate of approximately 2000 GRBs, in the 10 Gyr long history of the thin disk, is obtained. Although not fully constrained, these collimated high-energy phenomena are likely to be a significant threat to habitability and should be implemented in the astrobiological models of the GHZ.
Scalo & Wheeler (Reference Scalo and Wheeler2002) also imply that the supernovae (SNe) rate is approximately 105 higher than GRB rate. With the SNe fluence distances, smaller than 1 kpc (one of the largest values, suggested by Karam Reference Karam2002), and most probably even less than 10 pc, the largest SNe net astrobiological influence on GHZ should be approximately 102 times the net influence of GRBs. If GRBs are collimated in some preferred direction then this ratio will change, given the geometry of the GHZ. The resulting depletion of GHZ astrobiological complexity caused by SNe can be estimated to no more than 102 times the depletion caused by the GRBs.
The rates of core collapse SNe and GRBs should roughly follow the star formation rate and be consistent with the gas consumption rate in the Galaxy. The history of star formation in the Galaxy is not a smooth function and suffers a great deal of uncertainty. Despite the ‘spikes’ in star formation, on the long timescales relevant for this paper the star formation rate should be decreasing with time as more and more interstellar gas is consumed in the process. This will open longer evolutionary ‘temporal windows’, creating a stable astrobiological environment for simple life forms to evolve and diversify. Li (Reference Li2008) showed that 64 Swift GRBs are in accordance with the cosmic star formation history. The updated star formation rate in Li (Reference Li2008) is decreasing for red shifts smaller than four, similar to that used in Lineweaver (Reference Lineweaver2001). As argued in Paper I, if the GRB rate decreases over time (roughly following the star formation rate), then the GHZ may at some point exhibit an astrobiological phase transition. This transition will manifest itself as a rapid increase in the number of habitats with high astrobiological complexity. In addition, the metallicity buildup over cosmological timescales will result in greater abundance of planet-forming material (for the connection between planet abundance and host star metallicity, see Lineweaver (Reference Lineweaver2001), Fischer & Valenti (Reference Fischer and Valenti2005) and Israelian (Reference Israelian2005)). However, as pointed out by Lineweaver (Reference Lineweaver2001) and Pena-Cabrera & Durand-Manterola (Reference Pena-Cabrera and Durand-Manterola2004), the higher metallicity may result in enhanced Jupiter-like planet formations that can, by means of planetary migration, destroy smaller Earth-like planets plausibly making the planetary system inhabitable.
Solid body collisions
So far, no GRBs and only a few historical SNe (Clark & Stephenson Reference Clark and Stephenson1977) have been observed in the Milky Way, but due to their high luminosity, other galaxies have given us plenty of data. On the contrary, comets and asteroids can be observed only within our own planetary system. This fact gives rise to the selection effect that is a consequence of our very existence as observers (Ćirković Reference Ćirković2007). Stemming from the Bayesian probability analysis, Ćirković argues that there is a plausible chance that the rate of the impact-caused (and other) catastrophes in the Earth's evolutionary history (as inferred from the past records such as craters on Earth, Moon, Mars, etc), has been underestimated.
The characteristics for these type of GRMs (the probability of forming Jovian planet shields and the parameters of Oort clouds and asteroid belts, etc.; Emel'Yanenko et al. (Reference Emel'Yanenko, Asher and Bailey2007)) should be estimated by their dependence on global GHZ parameters: e.g., Matese & Whitmire (Reference Matese and Whitmire1996), Goncharov & Orlov (Reference Goncharov and Orlov2003), among others, implied that there is a plausible connection between the passage of the Sun through the Galactic plane and meteorite impacts. At the moment, simulating and estimating them is a much harder task than doing the same for the SNe. In addition, the GRMs with a small distance of influence are probably, to some extent, spatially biased due to their spiral arm pattern. This makes possible habitats spatially correlated, depending on their spiral arms crossing frequency; according to Fernández et al. (Reference Fernández, Figueras and Torra2001) the Sun is currently between the centre and the outer edge of one of the spiral arms and also very near the co-rotation circle: see also Shaviv (Reference Shaviv2003) and Gies & Helsel (Reference Gies and Helsel2005) for the long-term climatic connection. This, however, might be at least partially due to observation selection effects (Marochnik Reference Marochnik1983). In classical early studies, McCrea (Reference McCrea1975) investigates the possibility and consequences of our Solar System encountering a dense dust cloud during a spiral arm crossing event, while Clark et al. (Reference Clark, McCrea and Stephenson1977) discuss the frequency of nearby SNe as a competing causal mechanism for the increased frequency of climatic and biological catastrophes in relation to the Solar Galactic trajectory.
For future studies in this direction spatial modelling of the GHZ is required. This will demand more computational resources and a detailed knowledge of the kinematics and dynamics of the Galactic disk. For robust quantitative studies, locally acting versatile GRMs should be implemented (such as parent star flares and radiative changes, planet volcanism and others).
Induced correlations and global risk function
The GRMs are global in the sense that they act all over the GHZ, but paraphrasing George Orwell's Animal Farm: Some are more global than others, i.e., the parameters of the Oort cloud and asteroid belts that govern the solid body collisions with habitable planets are different from one planetary system to another. Although of smaller integrated energetic influence than SNe, GRBs are probably more important for the synchronization of astrobiological clocks in GHZ habitats (Paper I). An energetic GRB collimated in the plane of the Galaxy may comprise a very large part of the GHZ. It is very probable that few such events have occurred in several Gyrs of Galactic astrobiological past (e.g., Melott et al. (Reference Melott, Lieberman, Laird, Martin, Medvedev, Thomas, Cannizzo, Gehrels and Jackman2004) suggest that late Ordovician mass extinction was partly caused by a GRB).
A possibility for induced spatio-temporal correlations in the set of habitable planets is likely to be very important for studies that scope the history of the Galactic astrobiological landscape. If some class of highly energetic events (in further discussion referred to as type α GRMs) is capable of simultaneously influencing the large part of the GHZ than it will be possible, in part, to deduce the temporal evolution of Galactic astrobiological complexity from the time rates of such events. The best candidates for such agents are GRBs whose rates tightly correlate with the history of star formation. The smaller scale GRMs (type β), such as off-plane collimated GRBs or local SNe, will induce weaker and more spatially limited correlations between astrobiological complexity and star formation rate, restraining the accuracy of resulting astrobiological history. The total hazard of all GRMs (as well as other biologically adverse events) can be subsumed in the global risk function (GRF), which is unknown at present. The detailed reconstruction of the GRF will yield the astrobiological complexity of the GHZ throughout the history of the thin disk. The toy model presented in the next section implements the type α and type β GRMs in the GRF, but only future research will yield more complex GRFs, creating a more detailed astrobiological landscape.
A stochastic model of the Galactic Habitable Zone
Two major differences between this model and the model presented in Paper I are the somewhat different GRFs (Fig. 1) and a correction for the main sequence stellar life times. As stated in Paper I, GRF is modelled as a series of approximately 200 random events occurring with exponentially decreasing frequency:
with the fixed characteristic timescale tγ=3 Gyr in accordance with the cosmological observations (Bromm & Loeb Reference Bromm and Loeb2002). The probability to reset the astrobiological clocks of possible habitats was presented with the factor Q∊(0, 1). The model from Paper I should represent the type α GRMs. With GRBs being the best candidates for type α GRMs, the number of modelled events should depend on the GRB rate and geometry and the boundaries of the GHZ. Given the uncertainties of the latter it was chosen somewhat arbitrarily, but in accordance with the findings of the previous section. Type β GRMs can be considered as ‘noise’ in the GRF, while type α GRMs will present sharp spikes (Fig. 1). Type β are modelled as approximately 2000 (see previous section) random events with the characteristic timescale from Equation (1). Their reset probability is varied over various runs as a constant fraction of Q.
Fig. 1. Sketched view: GRF with type α (spikes) and type β (‘noise’) random reset events.
The correction for main sequence stellar life times is implemented as follows. The stellar mass is randomly chosen on the linear scale between 0.6 and 1.3 solar mass. The main-sequence life time in Myr is calculated as
where m is in the units of solar mass. Each temporal interval that represents the habitable planet has its own host star. An interval is considered completed if it is finished before the t ms of the host star expires, despite resetting events from the GRF. The length of each interval is chosen on a logarithmic timescale in the range [t ms/10, t ms]. In the absence of observational data and the existence of only arbitrary assumptions on the topic, this kind of modelling for the astrobiological timescales is adopted at present. As modelled in Paper I, the ages of planets are distributed according to the Lineweaver (Reference Lineweaver2001) age distribution for terrestrial planets, but the timescale is different. Instead of 0–13.81 Gyr, the range of 0–10 Gyr is used to represent the lifetime of the thin disk. Accordingly, Lineweaver's (Reference Lineweaver2001) data is sampled from the 3.81–13.81 Gyr range and scaled to the 0–10 Gyr range.
Results and discussion
In every run of the simulation 109 intervals were simulated. The lengths of the intervals are different for each curve presented in Figs 2 and 3, but this will not significantly alter the results and findings of this paper. Other results of the same model are presented in Paper I, Ćirković & Vukotić (Reference Ćirković and Vukotić2008), Vukotić & Ćirković (Reference Vukotić and Ćirković2008) and Ćirković et al. (Reference Ćirković, Vukotić and Dragićević2009).
Fig. 2. Number of habitable planets in the Milky Way with time. Comparative view of results from different models: Orange – without type β GRMs and main sequence correction; sky blue – without type β GRMs. For the model presented in this paper the reset probability of type β GRMs is smaller than for type α GRMs: vermilion – 500 times; bluish green – 300 times; blue – 100 times; black – 5 times smaller. Left: The results for two different random sets of reset events. Right: The same results on logarithmic ordinate.
Fig. 3. Same as in Fig. 2, with the same sets of random reset events.
As expected, in comparison with the model from Paper I (orange curve), the chosen model for the intervals length did not significantly (to an order of magnitude) reduce the number of completed intervals (sky blue curve). The characteristic stairway pattern, a consequence of the exponentially decreasing GRMs rate in Paper I, is present. Even without the type β events (sky blue curve) the same pattern is preserved, implying that type β events have no influence on the general pattern. This can be clearly seen on logarithmic ordinates even for black curves – a tentative support for the temporal correlations and findings of Paper I. A large number of less devastating events, uniformly distributed across the GHZ, will induce the same temporal correlations as a smaller number of more devastating events. In both cases, the resulting pattern will depend on the chosen time rate (Equation (1)). Increasing the lower limit for the length of the intervals from 102 (Paper I), to t ms/10 for the mass range 0.6–1.3, delays the noogenesis for approximately 4–9 Gyr after the formation of the thin disk. The time of Earth's formation then falls somewhere in the middle of this interval. From Figs 2 and 3, the number of habitats with high astrobiological complexity, and possible noogenesis, at the current epoch is up to approximately 102 times smaller than the initial number of habitats. The influence of the type β GRMs is small and it begins to converge for reset probabilities <Q/300 (Figs 2 and 3). This may change for different numbers of type β reset events, since their total influence on the GHZ can be roughly estimated as a product of their number and their reset probability. The last parameter depends on the geometry and boundaries of the GHZ. In addition, it will strongly depend upon the distance of astrobiological influence for the GRM type in question. Estimates of biologically effective distances for GRBs are presented and discussed in Galante & Horvath (Reference Galante and Horvath2007). The authors suggested that effects of nearby GRBs are difficult to model and that there is ample room to study scenarios that address these issues.
Given that the modelled biogenesis timescales are close to the real ones, the GRF model presented here yields, as an upper limit, that about one in a hundred habitats will achieve high astrobiological complexity. According to Figs 2 and 3, a different set of reset events will delay the astrobiological phase transition for varying amounts of time. The time delays seem milder for smaller Q, but for a catastrophic driven evolution (Q~1) the delay difference in Fig. 3 is approximately 1 Gyr. This opens up the possibility that the apparent Search for ExtraTerrestrial Intelligence failure (see Ćirković Reference Ćirković, Vukotić and Dragićević2009) can be understood as an artefact of the stochastic nature of GRMs. The characteristic astrobiological phase transition signature (Ćirković & Vukotić Reference Ćirković and Vukotić2008) is present within the extended (‘noisy’) model of the GRF, which further confirms the findings of Paper I, suggesting that we are currently living in the epoch of astrobiological disequilibrium. Obviously, a lot of further research is necessary in order to obtain more detailed quantitative predictions; as with all young fields, astrobiology is lacking sufficient numerical rigor at present, but the development of GHZ-related computer models presented here and elsewhere (Bjrk Reference Bjrk2007; Forgan Reference Forgan2009; Forgan & Rice Reference Forgan and Rice2010; Iorio Reference Iorio2010), combined with the growing data on existing and newly discovered extrasolar planetary systems, should change this in the near future.
Acknowledgements
The author would like to thank to two anonymous referees for very useful comments that have improved the quality of the paper. Milan M. Ćirković is acknowledged for many helpful discussions and conversations on the topics within the scope of this paper. Many thanks to Andrew Hodges for the careful reading and correction of a previous version of this manuscript. This work has been supported by the Ministry of Science and Technological Development of the Republic of Serbia through the project No. 146012: ‘Gaseous and stellar components of galaxies: interaction and evolution’.
