Habitable Zones in Binary Star Systems: A Zoology

Several concepts have been brought forward to determine where terrestrial planets are likely to remain habitable in multi-stellar environments. Isophote-based habitable zones, for instance, rely on insolation geometry to predict habitability, whereas radiative habitable zones take the orbital motion of a potentially habitable planet into account. Dynamically informed habitable zones include gravitational perturbations on planetary orbits, and full scale, self consistent simulations promise detailed insights into the evolution of select terrestrial worlds. All of the above approaches agree that stellar multiplicity does not preclude habitability. Predictions on where to look for habitable worlds in such environments can differ between concepts. The aim of this article is to provide an overview of current approaches and present simple analytic estimates for the various types of habitable zones in binary star systems.


Introduction
Defining habitable zones in binary star systems has its challenges. Sharing a system with two host stars means that the amount and spectral composition of the light arriving at a potentially habitable planet can vary on relatively short timescales [Eggl et al., 2012, Kane and Hinkel, 2013, Cuntz, 2014, Forgan, 2016. The multitude of interactions between binary stars and their potentially habitable worlds has led to the development of a variety of approaches that tackle the problem where one may find habitable worlds in such environments. Here, we discuss a subset of the various concepts used to determine habitable zones in binary star systems found in current literature. For the sake of transparency, we shall group the various approaches into the following categories: The remainder of this article is dedicated to discussing the various approaches and provide analytic approximations to calculate them where possible.

Single Star Habitable Zones
One of the first to introduce the concept of the Habitable Zone (HZ) was Huang [1959], who, under the assumption that the amount of radiative flux (insolation) is the main driver for climate evolution, proposed that a planet should be habitable if it orbits its host star at a distance where L is the luminosity of the star cast into a steradian (L sr), S is the insolation in solar constants (L au −2 sr = s ≈ 1361 W m −2 ) and r the distance between the planet and its host star measured in astronomical units (au). Equation (1) states that for a sun-like star (L = 1 L ) a planet receives one solar constant (S = 1 L au −2 sr) if it is on a circular orbit with a semi-major axis of one astronomical unit (r = 1 au). The Earth itself, however, is not on a perfectly circular orbit. Together with slight changes in the sun's luminosity this fact causes variations in the amount of light our planet receives. Since the Earth is still habitable, it is not unreasonable to think that the Earth's climate remains robust for a range of higher (S I ) and lower (S O ) insolation values. Hence, the single star habitable zone (SSHZ) can be defined as a circumstellar shell between The subscripts I and O stand for the inner and outer limit of the habitable zone. Huang [1959Huang [ , 1960 proposed that a terrestrial planet could remain habitable for insolation values between S I = 5 s and S O = 0.1 s , from well inside Venus' orbit all the way to the asteroid belt. Instead of using best guesses for habitable zone insolation limits later works used the principle that the presence of liquid water on the surface of a planet can regulate its climate [Rasool and de Bergh, 1970, Hart, 1978, 1979, Kasting, 1988, Kasting et al., 1993. The borders of the habitable zone could then be defined as follows: Too much insolation and oceans evaporate leading to a steep increase in surface temperatures. Too little insolation and the planet falls into a cold-trap.
To calculate habitable zones around main sequence stars Kasting et al. [1993] determined effective insolation thresholds (S I,O ) that would lead to climatic runaway states. Once found, those insolation limits could be translated to habitable zone borders using Equation (2). Kasting et al. [1993] also discovered that not only the magnitude of radiative flux but also its spectral distribution is relevant to planetary climates. Kasting et al. [1993] found that given the same amount of insolation energy, light originating from an M-class star is more potent in heating an Earth-like world than that of an F-class star. Given the effective temperature of a star (T e f f ), its spectral energy distribution can often be approximated with an energy equivalent black body emission profile. The insolation thresholds S I,O have, therefore, been parametrized as a function of T e f f of the star [Kopparapu et al., 2014] S I = a I + b I T + c I T 2 + d I T 3 + e I T 4 ,

The Trouble with Two Stars
One implicit assumption that goes into calculations of classical habitable zone limits for single star systems as given by Kasting et al. [1993] and Kopparapu et al. [2014] is that potentially habitable worlds move  Kopparapu et al. [2014]. about their host stars on circular orbits. In single star systems circular orbital motion traces stellar isophotes, i.e., regions of equal insolation around the host star. This is convenient, because it means that insolation can be considered constant over time. Near constant insolation is a good approximation in the case of the Earth, as the value of our 'solar constant' changes very little at present [Laskar et al., 2004]. Exoplanetary system architectures are diverse and can differ from the one of our Solar System, however. Many authors have investigated the effect of planetary eccentricity, and hence variable insolation, on potentially habitable worlds [e.g., , Williams and Pollard, 2002, Kane and Gelino, 2012, Bolmont et al., 2016, Kane and Torres, 2017, Way and Georgakarakos, 2017, Méndez and Rivera-Valentín, 2017, Deitrick et al., 2018, Adams et al., 2019, Guendelman and Kaspi, 2020. In such cases the assumption of constant insolation can break down entirely.
Similar situations arise when dealing with binary star systems where the distances between the planet and the stars change continuously. This means that the amount of starlight received by the planet can vary substantially. The situation is often more complex than having a single planet on an eccentric orbit around a single star because the variations of the planetary orbital elements can occur on different timescales. Huang [1960] realized this in his study on binary star systems published not long after his work on single star habitable zones. Huang also understood that, on top of a second source of radiation, the system has to allow for stable orbital motion 1 of the planet to be considered habitable. Being expelled from the binary star system is likely not conducive to planetary habitability. Since then, a substantial body of work has been dedicated to investigating the stability of planets in binary star systems [e.g., Dvorak, 1986, Rabl and Dvorak, 1988, Whitmire et al., 1998, Mardling and Aarseth, 2001, Holman and Wiegert, 1999, Pilat-Lohinger and Dvorak, 2002, Pilat-Lohinger et al., 2003, Pichardo et al., 2005, Doolin and Blundell, 2011, Jaime et al., 2012, Georgakarakos, 2013, Quarles and Lissauer, 2016). Among many other important results, it has been found that stable orbits are possible in the vicinity of circumstellar and circumbinary habitable zones, but that much depends on the exact setup of the systems involved. To have an approximate idea on where one expects stable and unstable systems one can resort to numerically generated fit functions [Dvorak, 1986, Rabl and Dvorak, 1988, Holman and Wiegert, 1999, Mardling and Aarseth, 2001, Quarles and Lissauer, 2016. One example of such an orbital stability fit function is reproduced here [Holman and Wiegert, 1999]. For planets orbiting a star in an S-type configuration, the stability limit in terms of the semi-major axis of the planet is given by where a p is the critical initial semi-major axis of the planetary orbit, a b and e b the semi-major axis and eccentricity of the binary star orbit and µ = m B /(m A + m B ) the stellar mass ratio. Here, A is the star that hosts the planet and B is the perturbing star. For a circumbinary planet we find The coefficients j-v and their respective uncertainties for circumstellar and circumbinary orbits are given in Table 2. Due to the vast parameter space of multi-body systems generic stability studies often use mass-less test-particles as proxies for planets. Given the intricacies of resonant interactions in N-body systems detailed numerical simulations are always recommended in order to test where a particular system permits stable orbits [Wittenmyer et al., 2012, Makarov andBerghea, 2014]. Later studies have shown that even planets on stable, initially circular orbits may become uninhabitable due to the gravitational perturbations the two stars cause in the orbit of the planet [e.g., Eggl et al., 2012, Andrade-Ines and Eggl, 2017. The complex interaction of gravity and stellar radiation in binary star systems is what makes the definition of habitable zones challenging.  Table 2: Orbital stability fit parameters and their respective uncertainties as given in Holman and Wiegert [1999] for prograde circumstellar and circumbinary motion of a massless test-planet.

Isophote Based Habitable Zones
Orbital dynamics aside we can define regions in binary star systems that would yield the correct amount of insolation to make a planet habitable. As we have two stars in the system, contour curves of insolation are not necessarily circular with one star at the center. Instead, the insolation geometry in the system has to be investigated in detail in order to derive isophote-based habitable zones [Kane and Hinkel, 2013, Müller and Haghighipour, 2014. Figures 1 and 2 show the instantaneous insolation values in S-type binary star systems akin to α Centauri (α Cen A, B) on close, circular orbits. Stellar parameters for α Cen A and B are given in Table  3. The continuous black lines in Figure 1 trace lines of constant insolation (isophotes). The four curves are drawn using spectrally weighted insolation values that correspond to the inner and outer habitable zone limits for each star. We shall refer to the area around each star enclosed by those isophotes as "isophote-based habitable zone" (IHZ). Figure 1: A binary star system with circumstellar habitable zones similar to α Centauri on a compact orbit. The plot shows the system at a distance of a b = 5 au. Single star habitable zones (green) are shown on top of the larger isophote-based habitable zones (red). If the stars have different spectral types we also have to account for the spectral energy distribution of the combined insolation. Kane and Hinkel [2013] proposed to solve this issue by superimposing the spectra of both stars weighted with the respective star-to-planet distances. Several other studies such as Eggl et al. [2012], , , Cuntz [2014Cuntz [ , 2015 have adopted a slightly different approach. Instead of combining the spectra and investigating the impact on the planetary atmosphere each star is assumed to heat the potentially habitable world independently. If we assume that star A and B have the same spectral energy distribution each star contributes a flux inversely proportional to the square of its distance to the planet. Let a and b be the distances between the planet and star A and the planet and star B, respectively. We can sum up both stellar flux contributions on the planet and require that L A a 2 to guarantee a planet is habitable. Luminosities per unit area L A,B are approximately constant for stars as long as they are on the main sequence. The distances a and b, however, are not. Similar to the implicit equation for an ellipse (x/u) 2 + (y/w) 2 = 1, Equation (7) trace the inner (I) and outer (O) habitable zone borders when left and right hand sides are equal. The underlying assumption that both stars have to have the same spectral energy distribution is a major shortcoming of Equation (7). A significantly better approach is not to use the same S I,O values for both stars, but to weight each star with its proper effective insolation constant instead, so that where ) are the habitable flux limits for the inner and outer edges of the single star habitable zone using the effective temperatures of stars A and B, respectively [Eggl, 2018]. We are using the effective insolation limits for each individual star as weights to model how much star A affects the climate of a terrestrial planet compared to star B. We shall, therefore, refer to SA I,O and SB I,O as "spectral weights". In a binary star system with two stars identical to our Sun, a planet would receive exactly one solar constant's worth of insolation with a solar spectral energy distribution (SA = SB = 1s ) where 1L sr Hence, for a = b, the effective insolation from both stars at a star-planet distance of a = √ 2 au is equal to one solar constant. We can rewrite Equation (8) in a more concise form by introducing spectrally weighted luminosities of the isophote-based habitable zone in binary star systems. Inserting in Equation (11), where x, y, and z are coordinates with respect to the origin of a convenient coordinate center that has both stars along the x axis at a distance of d from each other, we see that the isophote-based habitable zone borders depend on the mutual distance d between the two stars. Since d changes over time for binary star systems on eccentric orbits, the isophote-based habitable zone also varies with time [Müller and Haghighipour, 2014]. For co-planar systems, i.e., z = 0, one can find analytic solutions to Equation (11) by expressing y = F(x, d, A, B). This leads to a quartic equation that, although unwieldy, can be solved Cuntz [2014Cuntz [ , 2015. Alternatively, one can use either numerical methods or analytic approximations based on fixed-point iterations to solve Equation (11). The latter approach will be discussed shortly.
Figure 1 compares single star habitable zone insolation limits, where the contribution of the second star is ignored, to isophote-based habitable zones. In close binary star systems with α Centauri-like stellar components there is a clear difference between single star habitable zones and isophote-based habitable zones. The combined flux of the two stars causes the isophote-based habitable zones around both stars to extend toward each other. In S-type systems the largest displacement of the isophotes is registered along the line connecting the centers of the two stars. For binary stars on elliptic orbits, the isophote displacement is a function of their mutual distance d and, thus, time. Generally speaking, the stars influence each other's habitable zone most during their closest approach, i.e., near the pericenter where a b is the semi-major axis in [au] and e b the orbital eccentricity of the binary. Circumstellar habitable zones will be least affected by the respective companion stars when they are at apocenter Q b = a b (1 + e b ). For circular orbits, we have e b = 0 and, thus, To quantify the largest shift in isophote-based habitable zone borders we can substitute b = d − a in Equation (11). We shall only consider distances along line connecting the star centers (line of centers) for the moment.
Multiplying both sides of the equation by a we can interpret the resulting expression as a fixed-point iteration, where j is an integer and a j+1 denotes the distance of the (spectrally weighted) isophotes centered around star A after the j + 1 iteration step. As starting points we use the classical habitable zone limits, Choosing the positive roots of Equation (14) and stopping after the first iteration, we find the following approximation for the new isophote positions between the two stars  (14) to calculate the habitable zone borders on the opposite side of star A.
(16) Figure 3 illustrates the above points showing a zoom on star A in Figure 1. We see that the distances a +,− I,O are intersections of the isophotes corresponding to the habitable zone limits with the line of centers. This means that the isophote-based habitable zone borders along the line connecting the two stars are given by in the direction of the second star B (+), and in the opposite direction (−), respectively. Equation (17) (17) and (18) characterize the extent of the isophote-based habitable zone. The dimension of isophote-based habitable zone borders is that of a length, in our case [au], as can be seen from Equations (11), (15) and (16).
The difference between the single star habitable zone and the corresponding isophotes is largest in the direction of star B, e.g., The asymmetry between the two directions (+,−) leads to a deformation of the single star habitable zone into the tear-shaped isophote-based habitable zones. It is evident from Figure 1 that in some cases the single star habitable zone provides a reasonable enough approximation to the isophote-based habitable zone. As a rule of thumb that is the case for two sun-like stars that have a pericenter distance of no less than 10 au. To determine where the single star habitable zone is no longer a good proxy for the actual isophote-based habitable zone in S-type configurations we can use Equations (13) and (14) to calculate the relative displacement ∆a of isophote-based habitable zone borders with respect to the single star habitable zone borders. This yields Figure 4 illustrates the displacement of the isophote-based habitable zone borders as a function of the binary star's pericenter distance (q). Three cases are shown, one for α Cen A, one for α Cen B and one for a system of solar twins. One can see that the displacement is always larger for the outer isophote-based habitable zone borders than for the inner ones. Absent dynamical constraints this means that the presence of the second star leads to a net growth of the isophote-based habitable zone compared to the single star habitable zone. For very close binaries, on the other hand, the individual isophote-based habitable zones of both stars merge into a single circumbinary isophote-based habitable zone as shown in Figure 2. In order to calculate circumbinary isophote-based habitable zone borders we insert Equation (12) into Equation (11) and define δ := d/2. Still assuming a coplanar configuration (z = 0) this yields For very small separations of the binary δ → 0, Equation (20) gives where c = x 2 + y 2 is the distance of the planet to the origin of the coordinate system. For compact, equal mass binary stars c is the distance to the barycenter of the two stars. A crude approximation to the isophote-based habitable zone in such systems can be constructed if we assume that both stars have exactly the same location at the origin of our coordinate system. The circumbinary isophote-based habitable zone then resembles a classical habitable zone around a 'hybrid-star' featuring the combined spectrally weighted luminosities of both stars. The inner and outer borders of the circumbinary isophote-based habitable zone are then approximated by are denoted by empty and full squares, respectively. The radiative habitable zone as defined by Cuntz [2014] is the area shown in gray. The radiative habitable zone constitutes a hollow sphere with inner radius a + I and outer radius a − O .
This approximation is only reasonable for systems where d c I as Figure 2 reveals. If the distance between the stars is comparable to the inner limit of the 'hybrid-habitable zone', then isophote-based habitable zone and radiative habitable zone limits can no longer be calculated using Equation (21). Since the assumption that the single star habitable zone is a good starting point for the fixed-point iteration no longer holds we have to construct another fixed-point iteration with c I,O = A I,O + B I,O as initial guess. Restricting Equation (20) to the line of centers between the binary stars (y = 0) we find as approximants for the circumbinary isophote-based habitable zone borders. Once more, Equation (25) represents four equations as the corresponding spectral weights need to be accounted for in A I,O and B I,O . Figure 5 shows the evolution of the isophote-based habitable zone size as a function of the binary star pericenter distance and orbital eccentricity. For two sun-like stars on circular orbits, the isophote-based habitable zone borders can expand up to ≈25% before orbital instability starts to chew away on the outer isophote-based habitable zone limit. No stable circumstellar orbits are left in the isophote-based habitable zone when the sun-like binary stars orbit each other at a distance closer than 4 au. Systems similar to α Centauri allow for the isophote-based habitable zone to grow around the primary for as little as 10% before orbital instability sets in. For binaries on highly eccentric orbits the isophote-based habitable zone is substantially truncated even at large pericenter distances. Such an example can be found in Georgakarakos and Eggl [2019]. (19) showing the maximum displacement of the inner (I, dashed) and outer (O, continous) border of the isophote-based habitable zone as a function of the binary orbit pericenter distance q. The displacement is given in relative to the original single star habitable zone borders. A ∆a = 100% means that the new border is twice as far from its host star than the single star habitable zone pendant. We consider α Centauri-like systems and a binary consisting of two sun-like (G2V) stars.  Figure 4 lead to an expansion of the isophote-based habitable zone compared to the single star habitable zone. At the same time orbital stability restrictions shrink the size of the isophote-based habitable zone. A shrinkage of 50% means that the outer half of the isophote-based habitable zone is dynamically unstable whereas a shrinkage of 100% means that no planet can remain on a stable orbit in the isophote-based habitable zone. Expansion and shrinkage of the isophote-based habitable zone are computed with respect to the primary of an α Centauri like system and a binary consisting of two sun-like G2V stars.

Figure 4: A visualization of Equation
We conclude that if orbital stability of the planet is required, the maximum radiative contribution of the second star to the extent and location of the isophote-based habitable zone is relatively small. Which systems are expected to have a circumstellar habitable zone that is not truncated due to orbital instability? Combining Equation (15) with a simplified form of Equation (5) we find that the entire circumstellar habitable zone is stable if the binary pericenter distance is larger than where q b = a b (1 − e b ) represents the pericenter distance of the stellar binary and A O is the spectrally weighted luminosity for the outer edge of the single star habitable zone. Choosing values for A O and B O of the actual α Centauri system Equation (26) predicts a the stability of the entire circumstellar habitable zone around α Centauri A for q b > 11.8 au. This is in good agreement with the results presented in Quarles and Lissauer [2016], that find the orbital stability limit around α Centauri A to be close to 2 au which is also the outer limit of the single star habitable zone [Kopparapu et al., 2014]. α Centauri has a pericenter distance of q b ≈ 11.3 au. Conversely, Equation (26) suggests that the circumstellar habitable zone around α Centauri B would be stable even if q b ≈ 7.1 au. The more generous stability limit is due to the fact that α Centauri B is only half as luminous as our Sun. The outer border of the corresponding single star habitable zone is, therefore, at 1.24 au well within the dynamically stable region [Quarles and Lissauer, 2016].
For circumbinary planets, Equation (6) dictates that the distance between the planet and the center of mass of the binary must be several times the semi-major axis of the stellar orbit in order to allow for stable configurations. Using Equations (6) and (21) we can define an approximate limit condition for the existence of dynamically stable circumbinary isophote-based habitable zones, namely where Q b = a b (1 + e b ) is the apocenter of the binary star orbit, and A I and B I are the spectrally weighted luminosities for the inner edge of the respective single star habitable zones, see Equation (10). If the binary star orbit is tight enough, circumbinary isophote-based habitable zones are dynamically stable. As we shall explore in the next sections, dynamical stability and isophote-based habitable zone constraints may not be enough, however, to determine where Earth-like planets can be habitable in binary star systems.
Moreover, Equations (15) and (16) show that the deformation of the single star habitable zone is a function of the binary star distance d. The latter is, however, time dependent for all but systems where the stars have a circular orbit with respect to their common center of mass. To cope with this issue, Müller and Haghighipour [2014] introduced rotating, pulsating isophote-based habitable zones. Those can be thought of as analogous to pulsating coordinate systems or zero velocity curves that are sometimes used to study the elliptic restricted three body problem [Szebehely and Grebenikov, 1969]. Having time-varying habitable zone borders means that isophote-based habitable zones sweep over planets on relatively short timescales. This leads to the problem of determining to which degree planets that are only partly inside habitable zones are actually habitable -a topic of an ongoing investigation [Williams and Pollard, 2002, Bolmont et al., 2016.

Radiative Habitable Zones
One of the issues of using isophote-based habitable zones is that planetary orbits in binary star systems do not follow isophotes. As a consequence Cuntz [2014] introduced the so-called "radiative habitable zone" (RHZ). In analogy to single star habitable zones the radiative habitable zone is based on the assumption that planets are moving on circular orbits either around one or both of the stars forming the binary. The radiative habitable zone is then defined as the largest spherical shell that can be inscribed in the isophote-based habitable zone. As such its limits can be approximated through: Due to its symmetry and in contrast to Equation (25) the radiative habitable zone has only one equation for its inner and one equation for its outer edge. For the exact expressions we would like to refer the reader to Cuntz [2014Cuntz [ , 2015. Figure 3 shows the radiative habitable zone around star A of an α Centauri-like system with d = 5 au. The radiative habitable zone is smaller than the isophote-based habitable zone, but both start to coincide for large distances between the primary and secondary star. In fact, in the limit d → ∞, the radiative habitable zone, isophote-based habitable zone and single star habitable zone are identical.
Similar to the circumstellar (S-type) case, we can define a circumbinary (P-type) radiative habitable zone as Contrary to circumstellar radiative habitable zones, the existence of circumbinary radiative habitable zones is not guaranteed. For certain stellar distances d the circumbinary isophote-based habitable zone starts to deform and finally separate into individual stellar isophote-based habitable zones. Thus, at a certain limit distance d the shape of the isophote-based habitable zone does not permit to inscribe spherical shells anymore (see Figure 6). To determine when that is the case, let us investigate the isophote intersections with the ordinate of our coordinate system. Letting x = 0 in Equation (20) we find that Two important constraints follow from the above expression, one for the existence of the radiative habitable zone and the other for the merging of the individual isophote-based habitable zones of star A and B into a single circumbinary I HZ AB Equation (33) tells us that we can only inscribe spherical shells, if the inner radiative habitable zone border does not intersect the outer isophote-based habitable zone border. On the other hand, Equation (34) predicts that isophote-based habitable zones merge whenever the distance between the stars is smaller than twice the square root of the sum of the spectrally weighted luminosities.  (15) and (16), whereas the ones on the right result from Equation (25), respectively. The radiative habitable zone vanishes for the system in the left panel. The closer system on the right has a radiative habitable zone shown in gray. Figure 7: Dynamically informed habitable zones. The permanently habitable zone in blue represents the most conservative limits. A planet started in the permanently habitable zone will never exceed habitable insolation limits throughout its orbital evolution. In contrast, a planet in the averaged habitable zone (yellow) can experience ample excursions beyond habitable insolation limits as long as the insolation average permits liquid water near the surface of the planet. The green curve represents an orbit of the extended habitable zone as defined in Eggl et al. [2012].

Dynamically Informed Habitable Zones
Radiative habitable zones are largely sufficient to provide an idea as to where planets can be located in a binary star system and still retain liquid water near their surface. If we take into account that planets move on perturbed Keplerian orbits we have to acknowledge the fact that the amount of light a planet receives from the binary star can change drastically with time [Eggl et al., 2012, Eggl, 2018. How planetary atmospheres react to such changes is a matter of ongoing investigation [Williams and Pollard, 2002, Forgan, 2016, Popp and Eggl, 2017, Moorman et al., 2019. To reduce the complexity of the problem Eggl et al. [2012] introduced the concept of climate inertia on planetary habitability. An analogy would be the apparent thermal inertia [Price, 1985] for celestial objects, such as the Moon [Williams et al., 2017] and the Earth [Ermida et al., 2019]. If the climate of a planet reacts to changes in insolation with very little latency (low climate inertia), then the planet must remain inside the so-called permanently habitable zone (PHZ) to allow for liquid water to exist near its surface. For a planet to reside in the permanently habitable zone, maximum and minimum values of the insolation function must not exceed habitable limits at any time. The permanently habitable zone corresponds to a strict, "classical" definition of a habitable zone when orbital dynamics are taken into account. A more relaxed definition allows some parts of the planetary orbit to lie outside the permanently habitable zone. Following the argument of Williams and Pollard [2002] that the atmosphere and oceans of a planet can buffer insolation variability we can define an averaged habitable zone (AHZ). Insolation extrema are ignored as long as the time averaged insolation stays within habitable bounds. Formally the above habitable zones are defined as PHZ : max(S I ) ≤ 1 and min(S O ) ≥ 1 AHZ : is the combined spectrally weighted insolation on the planet, a function of time, and S denotes the timeaveraged combined stellar insolation. Again, a and b are the distances between the planet and star A and the planet and star B, respectively, and subscripts I, O represent the inner and outer edge of the habitable zone. Note that all dynamically informed habitable zones neither depend on angular variables nor on time.
Consequently, permanent habitable zone and averaged habitable zone form concentric rings around the center of reference, just as the classical habitable zone, see Figure 7.

Circumstellar Habitable Zones
Permanently habitable zones in S-type systems are regions where the following conditions hold: . The maximum eccentricity the planetary orbit attains during its evolution is called e max p = max t (e p (t)). The maximum eccentricity e max p of the planet depends on the orbital elements of the planet and the binary in a non-trivial way (Eggl et al. [2012], see the Appendix therein). Equation (36) can be solved numerically with respect to a p , however. The two solutions then correspond to the inner and outer permanently habitable zone limits. If the initial orbit of a potentially habitable world supports a semi-major axis in the range PHZ I ≤ a p ≤ PHZ O , orbital evolution will never carry the planet beyond habitable insolation limits. In systems where the binary does not influence the orbit of the planet strongly, i.e., e max p 1, we can estimate the permanent habitable zone borders analytically via where we can use a p = (A I,O ) 1/2 as initial guesses to calculate q p and Q p , respectively. When formulating the maximum insolation condition for the inner edge of the permanent habitable zone we have assumed that the radiative contribution of star B does not overpower that of star A. This condition reads If the planet receives more light at its apocenter than at its pericenter and Equation (37) must be adapted accordingly and If the climate of the planet has a high capacity to buffer changes in insolation, averaged habitable zone limits can be derived from insolation averages. To simplify the calculation of planetary insolation averages [Eggl et al., 2012] have introduced the so-called equivalent radii. Equivalent radii are constant distances with respect to the host star that yield the same average amount of insolation a planet would receive, were it on an elliptic orbit. In other words, the equivalent radius of a star-planet system corresponds to the semi-major axes of a circular orbit where a planet receives the same average insolation as it would on the original elliptic orbit. In contrast to Eggl et al. [2012] the equivalent radiir p andr b are here chosen so as to be consistent with two body insolation averages. In other words, where we usedr p ≈ (A I,O ) 1/2 . Secular orbit evolution theory yields relatively compact expressions for e max p and e 2 p , see, for instance, Andrade-Ines and . For planets on initially circular orbits those expressions read e max p = 2 , e 2 p = 2 2 (44) For dynamically less excited states, i.e., e p (0) = , we have A word of caution: Eccentricity estimates such as the one above are based on secular orbit evolution theory. As they lack short period and resonant terms, they may not always provide accurate estimates. More elaborate estimates can be found in Eggl et al. [2012], Georgakarakos and Eggl [2015], Andrade-Ines et al.
[2016], Georgakarakos et al. [2016]. Dynamically informed habitable zones for the α Centauri system are presented in Figure 8. The difference in the permanent and averaged habitable zones suggests a strong dependence of the extent of the habitable region on the climate inertia of a potentially habitable planet. Worlds that cannot effectively buffer variations in incoming radiation would not be able to retain liquid water near their surface over roughly 50% of the classical habitable zone around α Cen A. The situation is similar for potentially habitable planets orbiting α Cen B. If α Centauri had a higher orbital eccentricity, its circumstellar habitable zones would become dynamically unstable. Figure 8: Dynamically informed habitable zones around α Centauri A and B. The top panels show how the extent of dynamically informed habitable zones changes with the orbital eccentricity of the system. Permanently habitable zones (blue), averaged habitable zones (yellow), dynamically unstable zones (purple) and non-habitable regions (red) are presented. The vertical, dashed lines describe single star habitable zone limits. The horizontal line represents the actual α Centauri system. The bottom panels represent a top down view on the actual system. Dynamically stable circumstellar zones around star A and B are colored green.

Circumbinary Habitable Zones
The maximum insolation configuration in P-type systems occurs when the planet comes closest to the brightest star. The minimum insolation configuration is reached when the brightest star is farthest from the planet. Hence, we can derive permanent habitable zone limits via PHZ I : where µ = m B /(m A + m B ) and A I > B I . Here, as in S-type systems, the pericenter (q p ) and apocenter (Q p ) distances of the planet evolve with time. The maximum insolation is, thus, always related to the maximum in the orbital eccentricity attained by the planet with respect to time. To explicitly calculate the borders, Equation (47) can be solved numerically for a p . The corresponding analytic estimates for small planetary orbital eccentricities read with q p = c I = √ A I + B I and a p = q p /(1 − e max p ) for the inner border of the circumbinary habitable zone, and for the outer border. Maximum eccentricities are evaluated at q p and Q p for the inner and outer border, respectively. While Equations (44) and (46) still hold for e max p and e 2 p , the expression for the forced eccentricity in circumbinary systems is different from that in circumstellar systems. Following Moriwaki and Nakagawa [2004] we find that for circumbinary planets orbiting in the same plane as the stars. To calculate circumbinary averaged habitable zone borders, we make use of equivalent radii. Defininḡ the insolation averaged over a planetary orbit becomes The circumbinary averaged habitable zone then reads (53) Figure 9 shows dynamically informed habitable zones for a Kepler-35-like system. We assume that no other planets except a potentially habitable world are present. Circumbinary habitable zones in systems with similar stars do not depend as strongly on planetary climate inertia as S-type systems. This is evident form the small difference in the extent of permanently and averaged habitable zones for a broad range of stellar orbital eccentricities and a consequence of the fact that the forced eccentricity vanishes for mass ratios µ ≈ 0.5. Different stellar types in a binary, however, will cause the habitability of such systems to become more dependent to the climate inertia of the planet.

Self-Consistent Models
Combining analytic insolation estimates with precomputed spectral weights allows for a quick assessment of where habitable worlds can be expected in binary star systems. This approach does have its limits, however. Orbit evolution models based on the three body problem, for instance, do not account for additional perturbers. Other planets in the system can alter the orbit of a potentially habitable planet. More complete dynamical models are required to account for such effects [Bazsó et al., 2016, Forgan, 2016. Precomputed climate collapse criteria, such as runaway greenhouse or freeze-out limits of atmospheric greenhouse gases, may not always accurately reflect the response of a planetary climate to insolation forcing, either. Self-consistent simulations of climate and orbital evolution of a planet are more suitable to study such phenomena in detail including resonant responses. By coupling orbit propagators to climate models fully self consistent approaches can provide a more detailed picture of the climate evolution of a planet in binary star systems. Such examples are the works of Moorman et al. [2019], Forgan [2016], Popp and Eggl [2017], which used 1D climate models, longitudinally averaged energy balance models (LEBMs) and general circulation models (GCMs), respectively. Self-consistent habitable zone calculations are time-consuming, however, and tuned to a specific climate model. Hence, results have to be interpreted with care [Forgan, 2013, Popp and.

Comparing Habitable Zones
Comparing single star habitable zones, radiative habitable zones and dynamically informed habitable zones for the actual α Centauri system, as well as several circumbinary systems we find that analytically derived habitable zone estimates tend to coincide for systems with planets that can buffer insolation variations to a high degree. This can be seen in Table 4 when comparing radiative habitable zone to averaged habitable zone limits. Since perturbations on the orbit of the planet do not impact its habitability significantly in such cases, even results obtained from single star approximations differ very little from radiative habitable zone and averaged habitable zone values. On the other hand, if potentially habitable worlds have a lower climate inertia, which is likely the case near the outer rim of the single star habitable zone [Popp and Eggl, 2017], habitable zones in binary star systems could be substantially smaller. This can be seen in Figures 8 and 9 when the permanent habitable zone is compared to the averaged habitable zone. More details on this subject can be found in Eggl [2018]. Table 4: Habitable zone borders for the α Centauri and Kepler-35 system. All habitable zone values are given in [au]. The habitable zone limits for S-type A systems are given with respect to star A, and for S-type B systems with respect to star B, respectively. Circumbinary habitable zone borders are given with respect to the barycenter of the binary. Permanently habitable zones (PHZ) are derived assuming that the planet started on an initially circular orbit, whereas PHZ * values are derived for planets on orbits with forced eccentricity (e p = ). (i) the corresponding habitable zone borders may be affected by orbital instability. ( * ) In the case of Kepler-35 the single star habitable zone is derived using the combined flux of Kepler-35 A and B originating from the barycenter of the system. Stellar parameters for Kepler-35 A and B are given in Table  3. The gravitational effect of the exoplanet Kepler-35ABb has been neglected. ( * * ) same configuration as in Figure 2 right panel a b = 0.5 au, e b = 0, ( * * * ) same configuration as above only with e b = 0.5.

Summary and Conclusions
The aim of this article was to grant the inclined reader some insight into the various species of habitable zones in binary star systems. Understanding the rationale behind different concepts is crucial in choosing which concept is best suited for predictions of where to look for habitable worlds in a multi star environment.
Isophote-based habitable zones may be one of the most readily accessible concepts, but care has to be taken in interpreting the results. A large isophote-based habitable zone, for instance, does not necessarily imply that a system has a greater chance of hosting habitable worlds. An increase in isophote-based habitable zones compared to classical circumstellar habitable zones can be counteracted by dynamical instability and gravitational perturbations distorting the orbit of the planet. The combination of radiative habitable zone and orbital stability provides a better framework that only starts to break down when gravitational perturbations on planets with a low climate inertia become non-negligible. Dynamically informed habitable zones can deal with the later cases, but they require more knowledge of the properties of the system. Dynamically informed habitable zones based on analytic estimates presented in this work start to become inaccurate when resonances come into play. Full scale simulations can be used to study resonant phenomena in climate and orbital dynamics of exoplanetary systems. However, detailed simulations with climate models are computationally expensive and it can be difficult to generalize results. We, thus recommend that the model complexity is adapted to the specific use case.

Open Source Software
All materials, data, and computer code associated with this article are publicly accessible upon request to the corresponding author. Python codes for calculating and visualizing dynamically informed habitable zones are publicly available at https://github.com/eggls6/dihz.