Project Hephaistos – I. Upper limits on partial Dyson spheres in the Milky Way

ABSTRACT

Dyson spheres are hypothetical megastructures built by advanced extraterrestrial civilizations to harvest radiation energy from stars. Here, we combine optical data from Gaia DR2 with mid-infrared data from AllWISE to set the strongest upper limits to date on the prevalence of partial Dyson spheres within the Milky Way, based on their expected waste-heat signatures. Conservative upper limits are presented on the fraction of stars at G ≤ 21 that may potentially host non-reflective Dyson spheres that absorb 1–90 per cent of the bolometric luminosity of their host stars and emit thermal waste-heat in the 100–1000 K range. Based on a sample of ≈2.7 × 10⁵ stars within 100 pc, we find that a fraction less than ≈2 × 10⁻⁵ could potentially host ∼300 K Dyson spheres at 90 per cent completion. These limits become progressively weaker for less complete Dyson spheres due to increased confusion with naturally occurring sources of strong mid-infrared radiation, and also at larger distances, due to the detection limits of WISE. For the ∼2.9 × 10⁸ stars within 5 kpc in our Milky Way sample, the corresponding upper limit on the fraction of stars that could potentially be ∼300 K Dyson spheres at 90 per cent completion is ≲8 × 10⁻⁴.

Extraterrestrial intelligence, infrared: stars

1 INTRODUCTION

The construction of swarms of artificial, light-absorbing components around stars could in principle allow advanced extraterrestrial civilizations to harvest vast amounts of radiation energy in the form of star light (Dyson 1960). The observational signatures of such ‘Dyson spheres’ include waste-heat from the absorbing sphere, obscured direct star light (resulting in both an apparently underluminous star and potential temporal variations in brightness) and the effects of feedback on the properties of the star from the surrounding sphere (e.g. Dyson 1960; Wright et al. 2016; Wright 2020).

Unlike searches for technologically advanced civilizations based on communication signals, searches based on Dysonian signatures do not strongly hinge on the willingness of other civilizations to contact us. On the nine axes of merits for technosignature searches defined by Sheikh (2020), waste-heat signatures rate favourably on detectability and cost (because searches can largely be done with existing data), ancillary benefits (because candidates that turn out not to be Dyson spheres may still be interesting for mainstream astrophysics), inevitability (because the emission of thermal waste-heat radiation is difficult to evade from thermodynamic considerations), and duration (Dyson spheres can be very long lived).

A number of observational searches for individual Dyson spheres (e.g. Slysh 1985; Jugaku & Nishimura 1991; Timofeev, Kardashev & Promyslov 2000; Jugaku & Nishimura 2004; Carrigan 2009; Zackrisson et al. 2018) and for the large-scale use of Dysonian technology at extragalactic distances have previously been carried out (Annis 1999; Wright et al. 2014a, b; Garrett 2015; Griffith et al. 2015; Zackrisson et al. 2015; Lacki 2016; Olson 2017; Chen & Garrett 2021), although without uncovering any truly compelling Dyson sphere candidates. In terms of sample size, the largest searches for individual Dyson spheres (e.g. Slysh 1985; Jugaku & Nishimura 1991; Timofeev et al. 2000; Carrigan 2009) are based on the IRAS catalogue (Neugebauer et al. 1984), which includes ≈2.5 × 10⁵ point sources. However, modern photometric surveys covering optical, near-infrared (near-IR) and mid-infrared (mid-IR) wavelengths reach object counts of up to ∼10⁹ targets, and should allow for even larger search programmes.

As part of Project Hephaistos¹, we intend to carry out the largest search to date for partial Dyson spheres, using data on ∼10⁸ Milky Way objects from combined Gaia, 2MASS, and WISE catalogues. While Gaia, 2MASS, and WISE all provide photometric data, Gaia crucially also provides parallax-based distances, which allow the spectral energy distributions of the targets to be anchored to an absolute luminosity scale.

In this first paper, we explore the upper limits that can be set on the prevalence of near-complete (capturing |$\gtrsim 10{{\ \rm per\ cent}}$| of the radiation from their host stars) 100–1000 K Dyson spheres at distances up to ≈5 kpc based on their expected waste-heat signatures, using a combination of Gaia and WISE data. In a later paper, we will take a closer look at the targets with the most Dyson sphere-like spectral energy distributions.

While it has been proposed that Dyson spheres and similar radiation-harvesting megastructures could be constructed around a variety of stellar-mass objects, including white dwarfs (Semiz & Oğur 2015), pulsars (Osmanov 2016, 2018) and black holes (Hsiao et al. 2021), we here limit the discussion to Dyson spheres around main-sequence stars. We furthermore assume the Dyson spheres to be non-reflective and that feedback from the Dyson sphere on the star (Huston & Wright 2021) may be neglected.

Our overall method is described in Section 2 and the resulting Dyson sphere limits in Section 3. In Section 4, we discuss the likely nature of some of the Dyson sphere interlopers uncovered in this work and also explore the impact that some of the underlying assumptions in our Dyson sphere models may have on our upper limits. Section 5 summarizes our results.

2 METHODS

In this paper, we combine data from Gaia Data Release 2 (Gaia Collaboration 2016, 2018a) and AllWISE (Cutri & et al. 2014), where Gaia DR2 provides parallaxes and fluxes in three optical bands (G_BP, G, G_RP) and AllWISE mid-IR fluxes at 3.4, 4.6, 12, and 24 µm. The AllWISE program is an extension of the WISE program (Wright et al. 2010) and combines data from different phases of the mission.

To place upper limits on fraction of main-sequence stars that could potentially host Dyson spheres (DS), we define a set of colour–magnitude diagram (CMD) regions where objects compatible with DS are expected to turn up, and count the number of candidates falling into these regions for different distance-based subsamples. In Section 2.1, we describe the models used to predict the positions of DS in these CMDs. In Sections 2.2 and 2.3, we show some examples of how our models modify the spectrum and the photometric broadband fluxes of stars, and how the choice of DS model parameters affect the output. In Section 2.4, we explain the counting method used to impose our upper limits. Finally, in Section 2.5 we demonstrate how the models work on a sample of stars and explain some observed trends.

2.1 Dyson sphere models

When predicting the observational signatures of the composite system formed by a star and its DS, we model the stellar component as an obscured version of its original spectrum and the DS component as a blackbody whose brightness depends on the amount of radiation it collects.

The AGENT formalism described in Wright et al. (2014b) represents a useful way of parametrizing the different energy inputs and outputs of DS. Here, AGENT is mnemonic for αγϵνT_waste, where α represents the energy collected by the DS, ϵ the energy collected by other means (e.g. energy resources from their planet), γ the thermal waste-heat (for a given characteristic temperature T_waste), and ν other losses (e.g. neutrinos, gravitational waves, etc.). Energy conservation is described by equation (1) under this parametrization:

$$\begin{eqnarray} \alpha + \epsilon = \gamma + \nu . \end{eqnarray}$$

(1)

Throughout this paper, we assume that the civilization’s main energy input comes from the stellar radiation captured by the DS (ϵ ≪ α) and that the main form of energy disposal is through thermal photons (ν ≪ γ). Under these assumptions, the energy conservation and steady-state equation (1) becomes equation (2):

$$\begin{eqnarray} \alpha = \gamma . \end{eqnarray}$$

(2)

Additionally, we assume that DS act like grey absorbers and that the thermal loss resembles a blackbody with a temperature in the range of 100–1000 K. This range of temperatures is set by the detectability of the infrared excess radiation from the DS in the mid-IR wavelength range where WISE operates. Temperatures that are lower than 100 K would place this excess in the far-IR, whereas temperatures above 1000 K fall would shift it into the near-IR.

Under these assumptions, only two parameters are required to model the spectra of stars surrounded by DS: the effective temperature T_DS of the DS and the luminosity of the thermal radiation released by the DS. For the sake of simplicity, we re-define γ as the normalized DS energy output:

$$\begin{eqnarray} \gamma = \frac{L_{\rm DS}}{L_{\star }}, \end{eqnarray}$$

(3)

where L_DS is the luminosity of the DS and L_⋆ is the luminosity of the star hosting the DS before being obscured. Under this definition, γ can take values between 0 and 1. In the case of an isotropically radiating star, γ also represent the fractional solid angle of outgoing radiation intercepted by the DS (the covering factor) or, in more casual terms, the level of completion of the DS. Because of this, we use γ and covering factor interchangeably throughout this paper.

2.2 Spectrum

With the assumptions above, we can model the photometric and spectroscopic properties of an obscured star and its DS. In this section, we present a basic spectral model for this compound system.

The overall spectrum of this system corresponds to the sum of the spectral contribution from each component. The intrinsic shape of the star’s spectrum is unchanged, but its specific luminosity dims by a factor (1 − γ) due to the grey absorber assumption. Since the DS spectrum is considered a blackbody at an assumed effective temperature, the specific luminosity of the combined system becomes

$$\begin{eqnarray} L_{\nu } = (1 - \gamma)L_{\nu ,\star } + \gamma BB_{\nu }(T_{\rm DS},L_{\rm DS}), \end{eqnarray}$$

(4)

where L_ν is the specific luminosity of the combined DS + star system, L_{ν, ⋆} is the specific luminosity of the star before being obscured, γ is the covering factor, L_DS is the bolometric luminosity of the DS, T_DS is the temperature of the DS, and BB_ν(T, L) is the specific luminosity for a blackbody-like source with a given temperature (T) and bolometric luminosity (L).

Fig. 1 presents some illustrative examples of how a Sun-like blackbody spectrum (T_eff = 5778 K and L_⋆ = 1 L_⊙) is modified in the presence of a DS, under different assumptions about the DS temperature and covering factor. In both panels, a Sun-like blackbody spectrum is shown. In the top panel, we demonstrate how the compound spectrum changes due to a DS with a fixed temperature of 300 K and a covering factor of either γ = 0.1, 0.5, and 0.9. As seen, the DS causes a boost in the mid-IR and a drop in the luminosity of the stellar component. Both features depend on the DS covering factor. The drop is predominant in the optical and the near-IR. In the bottom panel, we show how the spectrum varies when we apply DS models with a fixed covering factor of 0.5 and DS temperatures of 100, 300, and 600 K. The signatures previously mentioned are recovered: a drop in stellar flux and a boost in the mid-IR. Since we are considering various DS temperatures, the mid-IR blackbody peak shifts in wavelength depending on the temperature of the DS. However, for the temperature range chosen here, the excess remains in the detectable wavelength window of the WISE mission.

Figure 1.

Modified Sun-like blackbody spectra due to the presence of different Dyson spheres. In both panels, the black solid line represents the unmodified blackbody (T_eff = 5778 K, L_⋆ = 1 L_⊙). In the top panel, DS models with a temperature of T = 300 K and covering factors of 0.1, 0.5, and 0.9 are depicted with grey solid, dashed, and dotted lines, respectively. In the bottom panel, DS models with a covering factor of 0.5 and temperatures of 100, 300, and 600 K are depicted with grey dotted, dashed, and solid lines, respectively. The coloured bands represent the wavelength intervals relevant for the Gaia, 2MASS, and WISE missions.

2.3 Photometry

Given our previous assumptions, equation (5) represents the magnitude (M) of the combined star+DS system. This magnitude depends on the DS’s magnitude (M_DS) and the star’s magnitude after being obscured (M_⋆) by the DS. This equation can be generalized to any photometric band, and applies to both apparent and absolute magnitude. For the rest of this section, we use absolute magnitudes. Throughout the rest of the paper, all magnitudes are given in the Vega system.

$$\begin{eqnarray} M = -2.5\log \left(10^{-M_{\star }/2.5} + 10^{-M_{\rm DS}/2.5}\right) . \end{eqnarray}$$

(5)

Here, M_DS represents the absolute magnitude of the DS that is calculated theoretically by integrating the blackbody spectrum. The spectrum of this blackbody depends on the temperature (T_DS) and the bolometric luminosity (L_DS = γL_⋆). Here, T_DS and γ are considered free parameters. However, the bolometric luminosity of the star is required as an input.

On the other hand, the magnitude of the combined system requires deriving the star’s brightness after being obscured by the DS structure. Under the assumption of a grey absorber, the obscured magnitude of the star becomes

$$\begin{eqnarray} M_{\star } = M_{\star ,O} -2.5\log (1 - \gamma) , \end{eqnarray}$$

(6)

where M_⋆ is the magnitude of the star after being obscured, M_⋆,O is the magnitude of the star before being obscured, and γ is the covering factor.

In Figs 2 and 3, we illustrate how the position of stars would change in the (G − W3, G) CMD due to presence of Dyson stars with various γ and T_DS. We use a Sun-like (T = 5778 K, L_bol = 1 L_⊙) and a Proxima Centauri-like star (T = 3042 K, L_bol = 0.0017 L_⊙) as examples.

Figure 2.

Colour–magnitude diagram showing main-sequence stars and fixed-temperature Dyson sphere models. The orange and red stars represent Sun-type and Proxima Centauri-type stars, respectively, while the grey dots represent other main-sequence stars. The coloured dots demonstrate how the Sun and Proxima Centauri shift their location in the diagram when we apply a DS model with T_DS = 300 K and different covering factors.

Figure 3.

Colour–magnitude diagram showing main-sequence stars and Dyson sphere models with a fixed covering factor. The orange and red stars represent a Sun-type star and a Proxima Centauri-type star, while the grey dots are other main-sequence stars. The coloured dots demonstrate how the Sun and Proxima Centauri would be shifted in the diagram when applying a DS model with γ = 0.5 and different temperatures.

2.3.1 Covering factor

In Fig. 2, we study how the covering factor affects our models for a fixed DS temperature of 300 K. First, we notice how the brightness in the G band decreases as γ is increased. According to equation (5), the magnitude of the system in the G band should be |$G = -2.5\log (10^{-G_{\star }/2.5} + 10^{-G_{\rm DS}/2.5})$|⁠. However, a 300 K blackbody does not radiate in the optical at a comparable scale relative to stars. Because of this, the magnitude of the system becomes G ≈ G_⋆ = G_⋆,O − 2.5log (1 − γ). As we can see, the drop in the optical depends only on the covering factor.

Another observation is that stars shift to redder regions (higher G − W3) of the CMD as the covering factor increases. There are two reasons for this phenomenon. First of all, as we increase γ, we increase the drop in the optical flux, which translates as an increase in the G magnitude, hence, a shift to redder regions. Secondly, we have an excess in the mid-IR, which implies a decrease in the value of the W3 magnitude and consequently a further redward displacement. We know that the brightening of the flux in W3 must depend on the covering factor since a higher covering factor denotes a brighter DS by definition. Both phenomena cooperate in shifting the stars to the red regions of the CMD.

2.3.2 Temperature

In Fig. 3, we illustrate the effect in the (G − W3, G) CMD of changing the temperature for a DS with fixed covering factor. In this figure, we model DS with a covering factor of 0.5 and temperature with values of 100, 300, and 600 K.

Since the covering factor is fixed, the drop in the optical is the same for every model, but every G − W3 colour shift is different because the DS temperature modulates the Planck function. In Fig. 3, we can see that the model that has the most significant shift is the one with a temperature of 300 K, for which the blackbody peaks in the W3 band. For the other two temperatures plotted, the peak occurs outside the W3 band which results in smaller G − W3 offsets from the main sequence.

2.3.3 Additional remarks

The analysis made in the previous sections was done by visualizing the behaviour of the stars in a G, G − W3 CMD. However, conclusions are general and apply to all combinations of Gaia–WISE colours for the range of temperatures adopted in this work.

Since we aim to derive upper limits based on WISE data, the range of DS temperatures considered in this paper is limited to 100–1000 K. DS operating outside this range (Wright 2020) would not produce significant excess light in the mid-IR, but instead in the far-IR (for temperatures below 100 K) or in the near-IR (temperatures above 1000 K).

In Section 2.3.1, we neglected the flux contribution of the DS in the optical since we were analysing a 300 K blackbody. We retain this assumption for the rest of the work since blackbodies in the range of temperatures adopted (100–1000 K) do not contribute significantly to the Gaia G_BP, G, or G_RP bands. Hence, we use equation (6) when we estimate how the G magnitude of a star is affected by the DS.

Additionally, an important conclusion from the examples shown in Figs 2 and 3 is how the colour shift produced by the DS depends on the colour of its host star. A DS constructed around a blue star results in a larger G − W3 shift than a DS with the same temperature and covering fraction built around a red star. This happens because a red star already has a low optical to mid-IR flux ratio compared to a blue star, thereby resulting in a smaller relative shift in colour.

2.4 Counting method

To estimate the prevalence of DS in the Milky Way, we first identify the stars whose magnitudes and colours are consistent with those of DS models. For this purpose, we make use of CMDs that include the absolute magnitude in the G band and the colours G − W1, G − W2, G − W3, and G − W4. We study the DS models in these diagrams to check if the stars in our samples fall in these regions for different combinations of temperatures and covering factors. This method is executed in two phases. The first parametrizes the boundaries between what is compatible with a DS and what is not. The second consists of checking if our stars fall in these regions.

The first phase is summarized in Fig. 4. We select main-sequence stars in our sample of stars and then apply DS models of a given temperature and covering factor. Finally, we parametrize the boundary created by these models. For every WISE band (WX) used in the CMD, we get a different parametrization of the border, where f_X symbolizes the function that mathematically describes the boundary for WISE band X. For the purpose of this paper, we say that all the stars falling on the right side of this boundary are sources compatible with a DS of a given temperature and covering factor.

Figure 4.

Summary of the first phase to count stars and DS candidates: parametrization of boundaries. First, we select main-sequence stars, represented with blue dots. Then, we apply DS models of a given T and γ over the main-sequence stars, represented with red dots. We parametrize the boundary of the DS models with a mathematical function f_X, where X is the number 1, 2, 3, or 4 to cover all the bands in the WISE mission.

The second phase consists of counting the stars that fall in the compatibility region of DS for all the combinations of CMDs that we are using. This stage is depicted in Fig. 5. Here, we see representative main-sequence stars and examples of stars exhibiting mid-IR excesses. The boundary is parametrized by the functions f₁, f₂, f₃, and f₄ for the colours G − W1, G − W2, G − W3, G − W4, respectively. A star is only considered fully compatible with a DS of a given T_DS and γ if it falls in the compatibility region for all the diagrams used. In Fig. 5, stars A, B, and C are examples that show how our criterion works. Star A falls in the compatibility region of all the diagrams used; hence, it represents a potential DS candidate for the T_DS and γ used. On the other hand, star B does not fall in the compatibility region of the diagrams G − W1 and G − W4, so we dismiss it as a potential DS. Similarly, star C does not fall in the compatibility region of the diagram G − W4. We discard these stars because the mechanism that produces excesses in specific bands cannot be the smooth blackbody continuum associated with the DS models for this combination of T_DS and γ.

Figure 5.

Summary of the second phase to count stars and DS candidates: location of stars in the compatibility regions for Dyson spheres. Here, we show all the colour–magnitude diagrams in use. In these diagrams blue dots stand for main-sequence stars, while the red lines represent the boundaries of the compatibility regions which are represented with a parametrization f₁, f₂, f₃, and f₄ for the diagrams G − W1, G − W2, G − W3, and G − W4, respectively.

Finally, the upper limit on the fraction of DS in the sample is obtained by dividing the number of Dyson sphere-like sources by the number of sources in the sample as described by

$$\begin{eqnarray} f_\mathrm{DS}(T_\mathrm{DS},\gamma)= \frac{N_\mathrm{DS}}{N_\mathrm{stars}(\gamma) + N_\mathrm{DS}}, \end{eqnarray}$$

(7)

where f_DS is the fraction of stars that resemble DS for a given T_DS and γ, N_DS is the number of stars that fall in the DS region of the CMDs, and N_stars is the number of stars that do not resemble DS.

An important point to consider when estimating our fractions turns out to be the number of stars [N_stars(γ)] in the denominator of equation (7). The instrumental detection limit naturally reduces the number of stars that can produce a detectable model of a certain γ. As an example, if the G-band apparent magnitude detection limit is m_limit, then detectable DS models can only be produced by stars with apparent magnitude m_G < m_limit + 2.5log (1 − γ), whereas fainter stars will have DS models too faint to be detected. We take this phenomenon into account when estimating the fraction of stars that are compatible with DS.

Additionally, we notice that some stars lack detections in the WISE filters. For a star below the WISE detection threshold, we use the S/N < 2 WISE upper limits on its flux, rather than the flux itself, when assessing whether it falls inside or outside a given DS compatibility region in the CMDs.

2.5 First sample

We use combined Gaia DR2–AllWISE data set for this work. The cross-match between Gaia DR2 and AllWISE sources is done by using the AllWISE best neighbour catalogue provided by Gaia. We use as many stars as possible, therefore our samples are just limited by the instruments, especially Gaia, that has a limiting magnitude G = 21. Since distances based on simply inverting the parallax can lead to highly spurious results in cases where the parallax uncertainty is high, we instead adopt the distances provided by Bailer-Jones et al. (2018), which are based on the Gaia DR2 parallaxes and a prior that depends on Galactic longitude and latitude.

As a first approach to estimating upper limits on the fraction of stars that could potentially host DS, we select a sample of stars within 100 pc in Gaia DR2 + AllWISE. In total 265 724 sources are located within 100 pc. Then, we use the method described in Section 2.4 to estimate the number of sources compatible with DS for different T_DS and γ. Since these stars are in the local bubble, we assume zero extinction (the effect of extinction is further discussed in Appendix  C).

Then, we proceed to select the main-sequence stars on which to base the DS models for this sample. Main-sequence stars are selected using their absolute G magnitude and G_BP − G_RP colour. Objects lacking this colour measurements are removed, but this affects only ≈0.8 per cent of the objects in the 100 pc sample.

We illustrate the selection of main-sequence stars in Fig. 6, which features a CMD with colour G_BP − G_RP versus absolute G magnitude, and lines that demonstrate the pruning of white dwarfs, objects with high astrometric excess noise (HAEN) and post-main-sequence stars. HAEN sources are stars for which the astrometric solution has failed and/or stars in multiple systems for which the single-star solution does not work (Gaia Collaboration 2018b). These stars correspond to the cloud between the main-sequence and the white dwarf sequence.

Figure 6.

Stars in the Gaia DR2 + AllWISE sample within 100 pc and the cuts we use. In this figure, all the stars in the sample are depicted with grey dots. We use the red line delimits main-sequence stars from red giants, while we use the blue line separates white dwarfs and high astrometric excess noise stars from main-sequence stars.

The functional form used to filter red giants corresponds to equation (8), whereas the function used to filter white dwarfs and HAEN is given by equation (9). In these equations, M_G corresponds to the absolute magnitude of the star in the G band and G_BP − G_RP its colour.

$$\begin{eqnarray} M_{G}\text{ }\lt \text{ }4 \text{ and } M_G \lt 7\cdot (G_\mathrm{BP} - G_\mathrm{RP}) - 3 , \end{eqnarray}$$

(8)

$$\begin{eqnarray} M_G \gt 3\cdot (G_\mathrm{BP} - G_\mathrm{RP}) + 5. \end{eqnarray}$$

(9)

Even though we attempt to reject post-main-sequence objects and sources with HAEN, these cuts are admittedly not perfect. Nevertheless, we show in Appendix  A that residual members from these classes of objects do not contaminate our models significantly.

2.5.1 Models and their behaviour

Once all main-sequence stars are selected, we proceed to apply DS models to these. However, the bolometric luminosity of each star is required as an input parameter, and Gaia DR2 only provides estimations for a few of them (Andrae et al. 2018). We estimate the missing luminosities by using a polynomial fitting function where we extrapolate from their absolute magnitudes in the G band. This can be done since we are working with main-sequence stars. The details on the luminosity fitting function used are in Appendix  B.

With all the bolometric luminosities, we apply our DS models to all the stars in the main-sequence sample. We do this for a range of temperatures between 100 and 1000 K and covering factors from 10⁻⁴ to 0.9. In Fig. 7, we show some examples of our models in CMDs of G − W2 versus M_G. We identify four trends, which we describe in more detail below.

Figure 7.

Examples of how our models behave in colour–magnitude diagrams for different combinations of DS temperature and covering factor. Main-sequence stars are depicted with grey dots, while blue dots correspond to our DS models. The upper-left panel represents the low covering factor case, where DS are indistinguishable from the MS. The upper-right panel shows the behaviour of cool DS. Only the optical drop is recovered for colours using the W1 and W2 bands. The lower-left panel depicts the intermediate case, when some DS separate from the MS, and some do not. The lower-right panel show the ideal case, when all DS are clearly distinguishable from the main sequence.

2.5.2 Case a – low covering factor

First of all, we have the case of DS with low covering factors, which is illustrated in the upper-left panel of Fig. 7. In this example, the models have a temperature of 300 K and a covering factor of 10⁻⁴. We can see how both main-sequence stars and models overlap almost completely. For a covering factor this small, both the optical dimming and the mid-IR excess are too small to be detected. Consequently, DS of this type are not distinguishable from normal main-sequence stars in the present data set.

2.5.3 Case b – cool Dyson spheres

In the upper-right panel of Fig. 7, we show models with a temperature of 100 K and a covering factor of 0.9. The models behave like the normal main sequence but shift downwards for this particular choice of colour (G − W2). This phenomenon happens because the blackbody spectrum for a 100 K DS gives a negligible contribution to the W2 band, so the only signature recovered is the drop in the optical flux. This phenomenon happens only for cool DS when we try to use the colours G − W1 and G − W2 to measure the infrared excess.

The same effect is not seen for hot DS (T_DS > 600 K) when G − W3 and G − W4 are used as a measure of infrared excess. Even though the blackbody spectrum for hot DS peaks at wavelengths longward of the W3 and W4 passbands, the tail of their spectra still contributes significantly at these wavelengths (as seen in Fig. 1).

2.5.4 Case c – intermediate cases

In the lower-left panel of Fig. 7, we have an example of an intermediate case, when the DS models are only significantly offset from the main sequence in the case of the brightest stars. In this example, the DS temperature is 450 K, and the covering factor is 0.1. The phenomenon described in Section 2.3.3 explains why these DS models are more clearly separated from the main sequence for brighter stars. In general, the range of magnitudes for which the DS are distinguishable from the main sequence depends both on the DS parameters (T, γ) and on the colour used.

2.5.5 Case d – ideal cases

Finally, in the lower-right panel of Fig. 7, we illustrate how the models behave for a DS with a temperature of 600 K and a covering factor of 0.9. In this case, all DS models are clearly separated from the main sequence. This happens when the DS blackbody peak falls in the wavelength range of the WISE bands, and the covering factor is so high that all DSs give rise to a significant colour shift, regardless of whether they are hosted by bright (hot, blue) or faint (cool, red) main-sequence stars.

2.5.6 Boundaries and counting

To estimate the number of objects potentially compatible with DS within a sample, we define four boundaries (one for each of the colours G − W1, G − W2, G − W3, and G − W4) between the DS models and the main-sequence stars. Fig. 8 shows an example of how these boundaries are established. For a given T_DS and γ we derive the DS models (red dots) and pick the leftmost objects among these models for every Δ(G) = 0.5 mag bin (black circles). With these objects, we determine the border by using a linear interpolator (black line) and extrapolate this line for values beyond the models we have. This procedure is carried out for every CMD independently.

Figure 8.

Colour–magnitude diagram that illustrates the border between main-sequence stars and sources compatible with Dyson spheres. In this diagram, we show the distribution of stars (blue dots) and the models (red dots) for a DS with T_DS = 100 K and γ = 0.9. The black dots represent the leftmost stars chosen every 0.5 G magnitude bin that are used to create the boundaries by using a linear interpolator. The red region represents the compatibility region for this model in this specific colour–magnitude diagram.

With every boundary defined, we proceed to count the number of stars on either side. However, objects are not considered DS candidates unless they fall in the DS regions of all four CMDs considered (including G − W1, G − W2, G − W3, and G − W4 on the x-axis).

3 RESULTS

3.1 Exclusion maps for the 100 pc sample

Using the computational machinery outlined in Sections 2.4 and 2.5.6, we now set out to derive upper limits on the fraction of objects that could potentially host DS within 100 pc in the Gaia DR2 + AllWISE combined data set. Since the CMD boundaries for DS candidates depend on the assumed DS temperatures and covering factors, so will the limits. By considering temperatures from 100 to 1000 K and covering factors from 10⁻⁴ to 0.9, we arrive at the exclusion map shown in Fig. 9.

Figure 9.

Exclusion maps for stars within 100 pc: Colour maps that represent the fraction of stars that are compatible with DS of a given temperature and a covering factor. In the top panel, we show the exclusion maps for temperatures between 100 and 1000 K, and covering factors between 10⁻⁴ and 0.9. In the bottom panel, we show the same exclusion map, but for covering factors between 0.1 and 0.9.

From this exclusion map, we can recognize the trends already discussed in Section 2.5.1. For covering factors below log₁₀(γ) = −2, the upper limit on the fraction of potential DS is close to unity, i.e. a very weak limit. This means that we are unable to set meaningful limits on such low-γ DS. In the case of covering factors between log₁₀(γ) = −2 and log₁₀(γ) = −0.25, the upper limits become progressively stronger, as the DS models for the brighter stars in the sample detach from the main sequence and the mid-IR excess becomes recognizable. For covering factors above log₁₀(γ) ≈ −0.125 (⁠|$\gamma \, \approx$| 0.75), the upper limits are below 10⁻⁴ for all temperatures considered.

3.2 Exclusion maps for stars out to 5 kpc

In Fig. 10, we show exclusion maps for the samples containing stars at distances up to 200, 1000, and 5000 pc. As the effects of interstellar dust reddening become more severe at large distances, objects exhibit more scatter in the redward direction across the CMD, which makes the main sequence more difficult to isolate. Instead of using equations (8) and (9) to filter out stars not located on the main sequence, we have in the 5000 pc sample resorted to equation (10) to select main-sequence candidates on which our DS models are based.

$$\begin{eqnarray} M_G \le 1.5\cdot (G_\mathrm{BP} - G_\mathrm{RP}) - 0.1 . \end{eqnarray}$$

(10)

Figure 10.

Exclusion maps at various maximum distances: Colour maps that represent the fraction of stars that are compatible with DS of a given temperature and a covering factor. In the left-hand, middle, and right-hand panels, we show the exclusion maps for the sample of stars within 200, 1000, and 5000 pc, respectively.

In Fig. 10, we limit the constraints plotted to covering factors above 0.1, since the upper limits derived for lower γ are very weak. Overall, the upper limits become progressively weaker when the maximum distance limit of the sample is increased.

The primary reason for this is that the fraction of Gaia stars that fall above the detection threshold of WISE drops with distance. At large distances, intrinsically faint, red main-sequence stars tend to drop out of the Gaia sample, which leaves a larger fraction of intrinsically brighter, bluer main-sequence stars. Because of this, the fraction of stars that lie just above the Gaia detection threshold, yet below the WISE detection threshold, grows with increasing distance.

As all the large star-forming regions of the Milky Way are located at distances beyond 100 pc, a larger number of Dyson sphere interlopers in the form of Young Stellar Objects (YSOs) also start entering the sample. Finally, samples featuring objects at larger distances (smaller parallaxes) will also include more objects with highly uncertain distance estimates, and hence highly unreliable absolute magnitudes. This includes intrinsically very bright objects from outside the distance cut-off (even extragalactic objects, like quasars) contaminating the samples. We stress that none of these samples are complete in a volumetric sense, and instead flux-limited by the Gaia DR2 and AllWISE detection thresholds.

Table 1 summarizes the DS constraints set by each sample for the case of T_DS = 300 K and γ = 0.1, 0.5, 0.9. As seen, the 100 pc sample gives the strongest constraints, and the 5000 pc the weakest. When comparing the upper limits on the Dyson sphere fraction for the case of γ = 0.9 DS, we note that the upper limits are weakened by about one order of magnitude as the sample size is increases by three orders of magnitude (∼2.7 × 10⁵ objects compared to ∼2.9 × 10⁸ objects, when going from a maximum distance of 100 to 5000 pc).

Table 1.

Summary of the samples analysed. In this table, we summarize the samples analysed in this work. We show the distance considered in each sample, the size of the sample, and the percentage of sources compatible with a DS of T = 300 K and covering factors 0.1, 0.5, and 0.9.

Radius (pc)	Size sample	Fraction of stars compatible with a DS with TDS = 300 K
		γ = 0.1	γ = 0.5	γ = 0.9
100	265 724	6.6 × 10−3	1.9 × 10−4	1.8 × 10−5
200	1847 472	2.2 × 10−2	4.7 × 10−4	7.4 × 10−5
300	5243 957	4.8 × 10−2	4.8 × 10−4	8.6 × 10−5
400	10 384 485	6.7 × 10−2	5.6 × 10−4	1.1 × 10−4
500	16 831 302	8.9 × 10−2	7.5 × 10−4	1.9 × 10−4
600	24 162 082	1.0 × 10−1	7.5 × 10−4	1.9 × 10−4
700	32 123 156	1.1 × 10−1	7.8 × 10−4	1.9 × 10−4
800	40 575 873	1.2 × 10−1	8.5 × 10−4	2.0 × 10−4
900	49 472 732	1.2 × 10−1	9.4 × 10−4	2.0 × 10−4
1000	58 769 015	1.3 × 10−1	1.0 × 10−3	2.1 × 10−4
5000	∼2.9 × 108	1.6 × 10−1	7.1 × 10−3	8.3 × 10−4

We stress that, when counting the sources that fall in given DS region for a certain γ₀, we include all the stars on the right-hand side of the boundary in the CMD. This means that we are selecting stars consistent with γ > γ₀ (for a fixed T_DS). However, this γ cannot reach unity, since the DS would then completely obscures these sources and render it undetectable in the optical.

4 DISCUSSION

4.1 Previous upper limits

Prior to the current study, the most recent work on upper limits was that of Carrigan (2009). As a preliminary part of his search, Carrigan studied the photometry provided by the IRAS mission. Based on the fluxes at 12, 25, and 60 µm, he concluded that less than one in 10 000 of the ≈2.5 × 10⁵ IRAS sources could be DS, in the sense of being consistent with pure blackbody spectra in the 150–500 K range across these IRAS bands. This is, for instance, similar in terms of sample size and limit for our 100 pc sample in the case of γ > 0.6 at T_DS = 300 K. Carrigan also searched for DS by analysing the atlas compiled by the Calgary group. This sample contained information from the IRAS low-resolution spectrometer for ≈1.1 × 10⁴ sources. After applying different cuts, he ended up with 16 objects called ‘somewhat interesting sources.’ However, he focused on pure DS, i.e. γ = 1. This assumption implies pure thermal blackbody emission from the DS and complete obscuration of the host starlight. Unlike Carrigan’s work, here we focus on the behaviour of partial DS, where the optical flux of the star remains relevant. Moreover, our method for estimating upper limits relies on the location of stars in different CMDs, whereas Carrigan focused on temperature estimations. The most significant difference between our upper limit estimates turns out to be the characteristics of the data set used, as our largest data set is roughly three orders of magnitude larger and also includes distance estimates for the sources.

4.2 Feedback on the host star

Huston & Wright (2021) recently highlighted the importance of the feedback from DS on the properties of their host stars. They found that the interiors of low-mass stars may be affected under the presence of high feedback values, whereas more massive stars (≥ 1 M_⊙) remain insensitive to internal changes. While our models neglect the impact of feedback on the host stars, the DS temperatures considered in our work correspond to DS radii much larger than the range for which such feedback would have detectable consequences. As an example, let’s take a solar-type star and the highest temperature used in our work (1000 K) with a covering factor of 0.90. The feedback factor, or percentage of the star’s luminosity reflected back on itself, is 0.04 per cent. For a 1 per cent change in a solar-type star’s nuclear luminosity, the feedback factor must be at least 45 per cent, and for a 1 per cent change in a solar-type star’s effective temperature, a feedback factor of roughly 5 per cent. Even for lower mass stars which are more strongly affected by feedback, a feedback level of roughly 1 per cent is required in order to change nuclear luminosity by 1 per cent. The effective temperatures of these stars are not significantly affected for any feedback levels below 50 per cent. Thus, the non-reflective and relatively cool DS models in this work would not produce feedback strong enough to cause noticeable changes in stellar properties.

4.3 Multilayer Dyson spheres

Many variations on the DS theme have been proposed since the seminal paper by Dyson (1960). One of them consists of many DS layers operating at different temperatures. In the case of such objects, our models would still apply for some specific configurations while potentially failing for others. For instance, in a two-layer system, the outer layer could absorb some of the thermal light emitted from the inner sphere. If the radiation of the highest-γ layer completely dominates over the other one, our limits would roughly correspond to the properties (T_DS, γ) of that highest-γ layer. However, in cases where both layers have similar contributions, the situation becomes more complicated. As an example, a double-layer system where the two layers operate at T_DS = 275 and 300 K respectively, and each have γ = 0.7, would result in the same number of compatible sources and upper limits as a single layer with T_DS = 300 and γ = 0.9. A double-layer system with T_DS = 675 and 700 K, with each layer at γ = 0.781 also produces the same result.

4.4 Dyson sphere interlopers

The WISE mission is sensitive to the expected mid-IR signatures of DS in the 100–1000 K range, but dusty sources can produce similar signatures. PAH emission features fall in the 3.4 and 12 µm filters. Additionally, the 4.6 µm filter measures the continuum emission from tiny grains. The 22 µm filter sees both stochastic emission from small grains and the tail of thermal emission from large grains (Wright et al. 2010). Hence, sources surrounded by thick dust shells can exhibit significant mid-IR fluxes and populate the same regions of Gaia–WISE CMDs as our DS models. While bona fide DS are expected to produce a smooth mid-IR continuum, our method of counting sources in separate CMDs does not distinguish between sources with emission peaks in the separate WISE filters and those for which the WISE colours reflect a pure blackbody spectrum. YSOs surrounded by dusty cocoons and evolved stars ejecting dust into the interstellar medium are therefore significant DS interlopers in the current study and likely dominate our upper limits on the fraction of objects exhibiting Dyson sphere-like characteristics. A spectroscopic analysis would be necessary to uncover the true nature of individual DS candidates, but a more detailed analysis of the spectral energy distribution from the optical to the far-IR can also help weed out many types of dusty objects.

To get a sense of the relevance of interlopers in our analysis, we have examined the auxiliary data available for the five objects at distances out to 100 pc that in our diagnostic CMDs are deemed compatible with DS characterized by T_DS = 300 K and γ = 0.9. The available 2MASS photometry reveal that all of these have very red optical–near-IR colours that indicate that they cannot be main-sequence objects, unless they happen to be subject to several magnitudes of optical extinction due to interstellar dust. The latter scenario is implausible for stars within 100 pc (e.g. Lallement et al. 2019), but it should be noted that these objects have highly uncertain distances (as is evident from significantly reduced parallaxes for four out of the five objects in Gaia EDR3 compared to DR2; Gaia Collaboration 2021) and may in reality lie much further away. In either case, they are not good candidates for DS of the type we are considering in this paper.

All of the five objects are located in nebular surroundings, as revealed by available mid-IR and/or far-IR images, indicating that they are likely to be YSOs associated with ongoing star formation. As there are no large star-forming regions known within 100 pc, this highlights the problem of contamination by more distant objects in the sample. For instance, one of these objects lies in the direct vicinity of (and is hence likely associated with) the Taurus molecular cloud, one of the nearest star-forming regions, at ≈140 pc. For the four objects with EDR3 parallaxes, the distance estimates presented by Bailer-Jones et al. (2021) also place them far beyond the 100 pc limit.

In the sample of stars within 100 and 200 pc, 110 sources are compatible with our T_DS = 300 K, γ = 0.9 DS models. Here, we resort to the classifications available in the SIMBAD data base, and bin all the possible labels into five categories: pre-main-sequence stars, stars, post-main-sequence stars, other types, and unclassified. The label ‘stars’ refers to sources that are simply classified as such, with no additional information about their type or evolutionary stage according to SIMBAD. Thirty sources (27.3  per cent) are classified as pre-main-sequence stars, eight sources (7.2  per cent) are classified as stars, five sources (4.5  per cent) correspond to post-main-sequence stars, and eight sources (7.2  per cent) fall in the ’other’ category. We note that most of these sources are uncategorized, and this trend becomes more pronounced for samples at even larger distances.

4.5 Temporal variations

Dyson spheres are usually envisioned as a swarm of small panels or satellites collecting energy from their host star. When hypothetically observing a star surrounded by a partial DS, the number of orbiting panels covering the stellar disc would fluctuate over time, which could give rise to temporal variations in the amount of unobscured star light emitted in our direction. However, whether such effects would be detectable depends on the sizes of the panels used. (Wright et al. 2016) states that for a vast number of small objects, the array of satellites might appear as a translucent screen that cannot be easily detected. Intermediate-size panels would generate time variations that can be mistaken for other phenomena like photospheric noise or asteroseismic variations. Even larger objects/panels might generate light curves characterized by aperiodic events of almost arbitrary depth, duration, and complexity. Hence, while variability could provide an interesting auxiliary diagnostic, we cannot easily dismiss DS candidates based on whether they display variability. Light-curve information should however be able to weed out some types of interlopers, like MIRA variables, that become interlopers when their distances are underestimated.

5 CONCLUSIONS

By analysing the 2.9 × 10⁸ objects in the combined Gaia DR2 + AllWISE data set samples, we have derived a set of conservative upper limits on the fraction of stars that could potentially host near-complete DS with covering fractions in the 0.1–0.9 range. The present analysis is limited to non-reflective DS operating at temperatures from 100–1000 K.

Our conclusions can be summarized as follows:

• The constraints tend to be the strongest for DS temperatures around 300 K (similar to the objects originally envisioned in Dyson 1960), and become progressively weaker at both higher and lower temperatures (Figs 3,  9, and  10).

• Among the ≈2.6 × 10⁵ stars located within 100 pc in the sample, less than one per ∼10³ objects is compatible with a 300 K DS with a covering factor of ≳0.1, less than one per ∼10⁴ sources with a covering factor ≳0.5, and one per ∼10⁵ with a covering factor of ≳0.9 (Fig. 9, Table 1).

• Upper limits derived at distances beyond 100 pc are weaker, as objects with significant mid-IR excess may still lie below the detection threshold of WISE. Among the ≈2.9 × 10⁸ sources in the sample out to 5000 pc, as many as one per ≈3000 objects could in principle be consistent with a 300 K, covering fraction 0.9 DS (Table 1, Fig. 10).

In mid-2022, Gaia DR3 will provide stellar parameters deduced from spectral analysis, which will allow stronger constraints on the nature of the stars with mid-IR excesses associated with DS.

ACKNOWLEDGEMENTS

MS acknowledges funding from the Royal Swedish Academy of Sciences. EZ acknowledges funding from the Magnus Bergvall Foundation and a sabbatical fellowship from AI4Research at Uppsala University. AJK acknowledges funding from the Swedish National Space Agency (SNSA).

This work presents results from the European Space Agency (ESA) space mission Gaia. Gaia data are being processed by the Gaia Data Processing and Analysis Consortium (DPAC). Funding for the DPAC is provided by national institutions, in particular the institutions participating in the Gaia MultiLateral Agreement (MLA). This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, and NEOWISE, which is a project of the Jet Propulsion Laboratory/California Institute of Technology. WISE and NEOWISE are funded by the National Aeronautics and Space Administration.

DATA AVAILABILITY

The data underlying this article will be shared on reasonable request to the corresponding author.

Footnotes

1

https://www.astro.uu.se/~ez/hephaistos/hephaistos.html

APPENDIX A: POST-MAIN-SEQUENCE AND HIGH ASTROMETRIC EXCESS NOISE SOURCES

The process of creating DS models for main-sequence stars used in this paper involves an attempt to first filter out post-main-sequence and HAEN sources from the sample. However, since it is unlikely that this procedure removes every single object of this type, a small fraction of the DS models generated may in fact end up being based on stars that are not on the main sequence. Here, we explore what effects this may have on our overall result, using the 100 pc sample from Section 2.5 as a reference.

First, we study the effect of white dwarfs (WD) and sources with HAEN. We select these stars by using equation (9) that was originally used to filter them.

We apply different DS models to all of these stars and compare their CMD signatures to models based on our sample of main-sequence stars. The result turns out to be quite similar, as illustrated in Fig. A1 for white dwarfs and HAEN stars, under the assumption of a T_DS = 100 K and γ = 0.9 DS model.

Figure A1.

Dyson sphere models (T_DS = 100 K, γ = 0.9) applied to main-sequence stars (MS), white dwarfs (WD), and sources with high astrometric excess noise. Grey dots represent our main-sequence stars (MS), while blue dots are the models applied to these stars (MS DS). On the other hand, orange dots are white dwarfs and stars with high astrometric excess noise (HAEN + WD), while yellow dots are the models applied to these stars (HAEN + WD DS).

For red giant stars, the behaviour is slightly different compared to the other sources. We define a star to be a red giant star if its Gaia magnitudes follow equation (8), originally used to filter them.

In Fig. A2, we demonstrate how DS models based on red giant behave compared to those based on main-sequence stars, using the same DS parameters as in Fig. A1. As seen, the models for main-sequence stars and post-main-sequence stars partially agree. However, the boundary of the DS-compatible region changes in the region of very bright stars. Since we do not find a significant number of stars located in this region of any of the CMDs used, we conclude that the effect on our overall limits must be negligible.

Figure A2.

Dyson sphere models (T_DS = 100 K, γ = 0.9) applied to main-sequence (MS) stars and red giant (RG) stars. Grey dots represent our MSs, while blue dots are the models applied to these stars (MS DS). On the other hand, orange dots are red giants (PMS), while yellow dots are the models applied to these stars (RG DS).

APPENDIX B: LUMINOSITY FITTING FUNCTION

The bolometric luminosity of a star is essential when modelling a system containing a DS. However, Gaia DR2 only provides luminosity estimates for a fraction of the objects in the sample. To compensate for the lack of data, we model luminosities using a third-order polynomial fitting function. We choose main-sequence stars in the sample within 100 pc with luminosity estimations to generate this model. Since Gaia DR2 only provides luminosities for stars with M_G < 10, we use PARSEC models (Bressan et al. 2012) to account for the missing data. PARSEC models correspond to theoretical stellar evolutionary tracks that provide stellar parameters as well as their photometry. We use the version 1.2S that incorporates the improvements presented in Chen et al. (2014), Tang et al. (2014), Chen et al. (2015). From these models, we select isochrones of evolved systems (10¹⁰ yr) to ensure the existence of low mass main-sequence stars. We also take tracks with metallicities between Z = 0.012 to Z = 0.018. This range should be representative of the thin-disc stars within 100 pc. We extend the range of stars with known luminosities by adding these models to G-band absolute magnitudes as faint as ∼ 15.5.

Using both types of luminosity estimates, we proceed to create a fitting function. However, we count with more observational data estimates than theoretical estimations. To avoid a bias favouring the regime of massive stars, we base the fitting function on a pruned sample with equal number of stars in each ΔG = 0.5 bin. In generating this sample, we randomly stars from Gaia DR2 for which the luminosity uncertainties are lower than 1 per cent. We moreover do not consider stars brighter than M_G = 0 since their uncertainties tend to be higher. The resulting fitting function is described by

$$\begin{eqnarray} \log (L/L_{\odot })& = & -4.98\times 10^{-4} M_G^3 + 1.64\times 10^{-2}M_G^2\\ &&-\,\,4.87\times 10^{-1}M_G + 1.98, \end{eqnarray}$$

(B1)

where L is the bolometric luminosity of the star and M_G is its absolute magnitude in the G band. The relation between the data used and the fitting function is shown in Fig. B1.

Figure B1.

Luminosity as a function of the absolute magnitude in the G band. The red dots correspond to the PARSEC model stars, while the blue dots correspond to the stars taken from the Gaia DR2 sample. The black line correspond to the fitting function.

APPENDIX C: EXTINCTION EFFECTS

The DS models used in this paper ignore the effects of interstellar dust. In general, dust causes optical obscuration of sources yet excess emission at infrared wavelengths. Because of this, objects embedded in nebular regions tend to be common interlopers when checking the nature of sources compatible with DS models in the CMDs used throughout this paper.

Here, we explore the impact of dust obscuration on our models by comparing models with and without the extinction corrections provided for a subset of Gaia DR2 stars Andrae et al. (2018).

We use a sample that contains data on stars within 5 kpc. In this sample, 83 per cent of them have A_G values below 1, while 97 per cent of them have values below 2.

While these models differ due to the changes in the input luminosities and magnitudes of the stars, we find that the boundaries used to define DS-compatible regions within the CMDs use remain unchanged. In Fig. C1, we exemplify this for the star with the highest extinction estimate (A_G = 2.3 mag). As seen, the models with and without the extinction corrections both lie very close to the boundary. This trend holds for all the CMDs used in this paper.

Figure C1.

Colour–magnitude diagram showing the behaviour of DS models based on stars with and without extinction corrections. The red star represents a star with A_G = 2.3 prior to extinction correction and the orange star the same object after correction. The red and orange dots represent DS models (T_DS = 300, γ = 0.9) applied to these stars. Grey dots represent uncorrected stars within 100 pc. The blue line is the boundary that separates the DS and non-DS regions, derived from DS models based on stars not corrected for extinction.