© The Authors 2024

1 Introduction

During the Big Bang, a highly energetic environment sets the foundation for theories to emerge regarding matter and its composition. Evolved from different cosmological environments, modern theories describe a versatile zoo of dark matter models ranging from, for instance, sterile neutrinos (Dodelson & Widrow 1994; Shi & Fuller 1999), axion-like particles (Georgi et al. 1986; Masso & Toldra 1995; Masso 2003), weakly interacting massive particles (Steigman & Turner 1985), and quark nuggets that typically exhibit nuclear density (Witten 1984; Farhi & Jaffe 1984; de Rujula & Glashow 1984). In addition to the comparatively low-mass dark matter candidates, massive astrophysical compact halo objects (Alcock et al. 1993; Aubourg et al. 1993) with for instance primordial mass black holes could propose a promising alternative as well (Zel’dovich & Novikov 1967; Hawking 1971; Chapline 1975). Constraints on mass, properties, and number densities of these dark matter candidates can be made by various techniques. As for potential dark matter candidates, particle accelerators are instrumental in proposing mass boundaries for elementary particles, while gravitational wave detectors may aid in constraining the distribution of primordial black holes. A dark matter candidate is rarely observationally falsifiable by testing signatures in the electromagnetic spectrum.

Here, we present a detailed analysis of the observational feasibility ofa dark matter model that is proposed to obey properties of cold dark matter, but is still expected to be the origin of observable signatures that can be tested in a wide range of the electromagnetic spectrum, called the axion quark nugget (AQN, Zhitnitsky 2003). These nuggets are not only capable of describing generic cold dark matter features of structure formation but also provide solutions to fundamental cosmological problems. Two of the most prominent mysteries are related to (i) why the abundance of dark matter is not extraordinarily lower or larger than the visible component; and (ii) why we tend to observe more matter over antimatter and where this imbalance comes from.

This paper is structured as follows. In Sect. 2, a short overview is provided by focusing on the formation, structure, and interaction scenario of AQNs in the respective subsections. Furthermore, the foundation of interaction-resulting electromagnetic signatures is set in Sect. 3. A physical description responsible for the electromagnetic counterpart is provided in Sects. 3.1 and 3.2. The implementation method using a smoothed-particle hydrodynamics (SPH) formalism in a constrained cosmological simulation is presented in Sect. 3.3, with more details on the simulation provided in Sect. 3.5. Selection criteria and details on our underlying sample are presented in Sect. 3.6. Our results are provided in Sect. 4, followed by a discussion of our findings in Sect. 5. An outlook in conjunction with a concluding assessment is provided in Sect. 6.

2 General overview of AQNs

2.1 The formation process

In this section, we describe the key steps for the formation of these particles and refer to (Liang & Zhitnitsky 2016; Ge et al. 2017, 2018, 2019) for further details. Quantum chromodynamics (QCD) has a chiral anomaly term characterized by an axial angle θ. This angle θ is physical and observable in the standard model. When θ ≠ 0, the anomaly term induces a charge-conjugation Parity (C𝒫) violation, namely, the laws of physics do not behave in the same way under combined transformations of charge conjugation and parity. This leads to the so-called strong C𝒫 problem, since the C𝒫 violation in QCD is not observed (Abel et al. 2020). The strong C𝒫 problem can be naturally resolved when θ is considered a dynamical field so that the vacuum expectation value ⟨θ⟩ settles to zero after the QCD transition. Consequently, the dynamical field θ produces a hypothetical particle called axion¹. The axion is very light and couples very weakly with ordinary matter.

In the framework of AQN, axions are produced before inflation. In this scenario, a unique physical vacuum occupies the entire Universe, and most types of topological defects [including the N_DW, 1 domain walls (DWs)²] are forbidden. Only the N_DW = 1 DWs, when θ interpolates in a single physical vacuum but different topological branches θ → θ + 2πn (Vilenkin & Everett 1982; Sikivie 1982), can be produced. The N_DW = 1 DWs may form when the axion field starts to oscillate due to the misalignment mechanism at temperature T_osc ∼ 1 GeV. This process continues until the QCD transition T_c ∼ 170 MeV when the quark chiral condensates form. From T_osc to T_c, the axion mass effectively turns on, and the axion oscillation is underdamped. At the beginning of the oscillation T_osc, the axion field is coherent globally. A coherent nonzero θ violates the global C𝒫 symmetry and affects the DW formation. Consequently, a preferred matter or antimatter species of nuggets will be formed – it should be mentioned though, that an equal amount of nuggets would have been formed if θ = 0 at the moment of the DW formation. Only the initial sign of the coherent θ field before oscillation dictates which nugget family is preferably formed.

A small portion of the N_DW = 1 DWs forms closed bubbles and acquires baryon charges from the quark-gluon plasma. Near and after the QCD transition T_c, the DWs on the bubbles can mix the axion with a tilted η′ field (the field of the η^′ meson). Such a tilted substructure boosts the charge accumulation of the DW bubbles. Depending on its inherent topological charge, a bubble acquires either matter or antimatter charges during this phase. Charge accumulation terminates at T_form ∼ 41 MeV and the bubbles form quark nuggets in the form of color superconducting (CS) condensate. Due to the C𝒫 violation effect discussed earlier, more antimatter AQNs are formed compared to matter AQNs. Under the assumption of zero baryon net charge in the Universe, the observed visible-to-dark-matter-density ratio of Ω_dark : Ω_visible ≃ 5:1 implies a baryon charge ratio of $$B_{\overline{\text{AQN}}}:B_{\text{AQN}}:B_{\text{visible }}\simeq 3:2:1,$$(1)

where the subscripts of B correspond to the baryon charges of antimatter AQNs, matter AQNs, and visible matter, respectively. Both, the matter- and the antimatter-AQNs serve as the dark matter component, while the remaining quarks in plasma become the visible component. Unlike conventional dark matter candidates such as weakly interacting massive particles (WIMPs) and axion, the AQN explains the observed cosmological relation Ω_dark ~ Ω_visible without fine-tuning of fundamental parameters, since dark and visible matters now have the same QCD origin. In the AQN framework, freely propagating axions may exist via the misalignment mechanism (Abbott & Sikivie 1983; Dine & Fischler 1983; Preskill et al. 1983) and topological defects (Chang et al. 1998; Kawasaki et al. 2015; Fleury & Moore 2016; Klaer & Moore 2017; Gorghetto et al. 2018). However, they contribute a negligible amount to the dark sector Ω_dark unless very specific fine-tuning applies (Ge et al. 2018). The conclusion (1) is unaffected.

The CS condensate in the AQN is stabilized by the excessive surface tension of the DW bubble. The binding energy in the CS phase is sufficiently large (with a gap Δ ~ 100 MeV) so that AQNs do not modify the basic elements³ of the Big Bang Nucleosynthesis (BBN) at T_BBN ∼ 1 MeV, where elements of metallicity Z ≤ 3 were formed. Quark cores in the CS state are in the lowest energy state possible, and it is assumed that quarks in neutron stars obey the same state (Alford et al. 2008). AQNs remain in this stable state over cosmic time scales.

On the other hand, if the theory of AQNs holds, the asymmetry of matter to antimatter in the observable Universe can be explained by an asymmetric baryon charge separation. A coherent background field θ violates the global $\mathcal{C}\mathcal{P}$ symmetry at the beginning of DW bubble formation. This theory allows an equal amount of antimatter and matter, with more antimatter being hidden in composite nuggets than the ordinary one, causing it to be unobservable. These nuggets are stable over cosmic time, as the quark core is in a dense CS state. It serves as cold dark matter because these macroscopic particles exhibit a low number density due to individual masses of the order of grams (a more specific description will be addressed in the following section). In summary, the AQN model provides intriguing cosmological implications by explaining the matter-antimatter asymmetry and the similarity between dark- and visible-matter densities. This is because dark and visible matters share the same QCD ancestor.

2.2 Structure of AQNs and cosmological properties

It is important to differentiate matter families as (anti-)AQN consist of a core of (anti)quarks surrounded by positrons/electrons that make the electrosphere. The electrosphere is surrounded by the N_DW = 1 axion DW. Interactions of baryonic matter with the axion DW seldom occur. The surface tension of the DW exerts a strong pressure on the quark core causing it to preserve in a CS state. To prevent the DW from collapsing further, Fermi pressure from the AQN counteracts the DW surface tension.

Parameters that describe the properties of an AQN are the baryon number B, the electric charge eQ, where Q typically corresponds to the number of positrons that depleted from an originally neutral AQN and in some scenarios the magnetization ℳ (Santillán & Sempé 2020), which will be neglected in this study for reasons of simplification. It will be shown in later sections that Q highly depends on the temperature of the surrounding gas as positrons can deplete if the conditions suffice. The lower limit on B can be observationally constrained by non-detections from IceCube with 〈B〉 > 3 × 10²⁴ (Lawson et al. 2019). Similar limits are also obtainable from the ANtarctic Impulse Transient Antenna (ANITA) and geothermal data (Gorham 2012). B results in the AQN’s mass M_AQN and size R_AQN, with most recent discussion on their distribution in Raza et al. (2018); Ge et al. (2019). Since mass contributions from the electrosphere are negligible, the baryonic mass component of an AQN can be estimated as follows: $$M_{\text{AQN}}\approx m_{p}B\simeq 16.7\left(\frac{B}{{10}^{25}}\right)\text{g}.$$(2)

In the state of a CS phase, we can assume the nuclear density to be ρ_n = 3.5 × 10¹⁴ gcm^–3 (Zhitnitsky 2018) resulting in an AQN radius of $$R_{\text{AQN}}={\left(\frac{3M_{\text{AQN}}}{4\pi {\rho }_n}\right)}^{1/3}\approx 2.25\times {10}^{-5}\text{cm}.$$(3)

Ionized AQNs are characterized by a charge Q, with Q being a function of T_AQN . A small fraction of positrons in the electrosphere are loosely bound and in the non-relativistic Boltzmann regime (McNeil Forbes & Zhitnitsky 2008b), which rises with the internal AQNs temperature.

Throughout this paper, we fix the parameters to physically describe the structure of AQNs as B = 1025, M_AQN = 16.7 g and R_AQN = 2.25 × 10^–5 cm. A distribution function of these parameters is not well-constrained yet and would drastically increase the complexity of the system.

2.3 Possible interaction scenarios

As pointed out in Sect. 2.2, AQNs consist of two different particle families and will therefore yield different combinations of possible interaction scenarios. Not all of the different interactions are likely to be visible and, therefore, their relevance for our study will be assessed in the following.

1. Antimatter AQN interaction with electrons: electrons can interact with the electrosphere made up of positrons for antimatter AQNs. These e⁺e^– annihilations might possibly explain the 511 keV line emission coming from the galactic center (Oaknin & Zhitnitsky 2005; Zhitnitsky 2007; Forbes et al. 2010; Flambaum & Samsonov 2021) and the COMP- TEL gamma-ray photons in an energy range of 1-20 MeV (Lawson & Zhitnitsky 2008).

2. Antimatter AQN interaction with protons: free protons propagating through the electrosphere can collide with the antiquark core (see Fig. 1). Each collision produces an annihilation energy of 2 GeV. A major fraction (about 90%) of the annihilation energy will thermalize the quark core and transfer to the electrosphere. This heating process leads to the emission of the electrosphere in the form of thermal Bremsstrahlung. Observations that can be explained by thermal AQN emission are, for instance, the diffuse galactic microwave excess (McNeil Forbes & Zhitnitsky 2008b), the 21 cm absorption line (Lawson & Zhitnitsky 2019), and diffuse galactic UV background emission (Zhitnitsky 2022a). A minor fraction (about 10%) of annihilation energy will be emitted as a short pulse of non-thermal photons at the point of impact. Observations, such as Chandra’s diffuse k_BT ≈ 8 keV emission (McNeil Forbes & Zhitnitsky 2008a) might find their nature in non-thermal AQN emission processes.

3. Antimatter AQN annihilation with celestial bodies: as it is discussed in detail by Liang (2022), antimatter AQNs can annihilate with more matter at the same time. Axions, residing off-shell in the DW can be emitted, when AQNs lose potential energy after annihilation processes in the core took place. To retain a minimum total energy of the system and therefore maintain the stability of the AQN, axions will be emitted to lower the DW’s mass. In addition to axions, antimatter AQNs produce electromagnetic radiation as described in (1) and (2). As a consequence, phenomena such as coronal heating, extreme ultraviolet (EUV), and X- ray emissions might be attributed to the outcomes of AQN interactions in the sun (Zhitnitsky 2017, 2018; Raza et al. 2018; Ge et al. 2020). The nature of fast radio bursts may be explained by AQNs by triggering magnetic reconnection events in magnetars as well (Van Waerbeke & Zhitnitsky 2019). Furthermore, according to Budker et al. (2022), infrasound, acoustic, and seismic waves might be the direct outcome of AQN-earth encounters. Non-inverted polarity detections by ANITA were discussed by Liang & Zhitnitsky (2022a), and short time bursts detected by the Telescope Array (TA) were analyzed by Zhitnitsky (2021b); Liang & Zhitnitsky (2022b) in the framework of AQN annihilations. Cosmic ray (CR)-like events detected by the Pierre Auger Observatory might be explained by AQNs inducing lightning strikes with subsequent direct emission (Zhitnitsky 2022b). Multi-modal events detected by Horizon-T that were discussed using AQNs by Zhitnitsky (2021a) should be distinct from conventional cosmic-ray showers.

4. Matter AQN interaction with antimatter: technically, one could expect the same emission processes from points (1) and (2) with matter AQNs with positrons and antiprotons. However, because of the low abundance of antimatter CRs (e.g. Hillas 2006; Blum et al. 2017) and the slightly lower abundance of matter AQNs over antimatter AQNs, it is expected that matter AQNs remain dark for most of the time and electromagnetic signals originating from matter AQNs can be neglected.

5. Matter AQN-antimatter AQN interaction: even though the electromagnetic outcome of matter-antimatter AQN interactions would result in an impressive electromagnetic signature, it is because of the low number density of n_AQN ~ 10^–29 cm^–3 that this interaction scenario can be entirely neglected (see comments on this for example in Zhitnitsky 2022b). And although both of the AQNs can have an attractive influence on each other if ionized, the higher effective cross-section would not drastically increase interaction rates.

6. AQN-AQN interaction from same particle family: Since AQNs from the same particle family would only collide by chance with no interaction enhancement by their ionization state, collisions are expected to be even rarer than described in point (5). An electromagnetic signature would at most be visible by thermalization of a heated AQN after collision.

As we have assumed the intracluster medium (ICM) to consist of ionized hydrogen, this work will only focus on emission phenomena coming from proton-antimatter AQN interactions. Line emission from e⁺e^– annihilation, axion emission by strong AQN emission, and AQN-AQN self-interactions will be neglected in this study as thermal and non-thermal AQN emission from proton interactions is presumed to dominate the spectrum.

With the fixed parameters describing a single AQN, we can infer the antimatter AQN’s number density by assuming the mass ratio relation from Eq. (1). According to the mass ratio relation, a fraction of 3/5 of the entire dark matter accounts for anti- AQNs. Since ionized gas in the ICM mainly consists of normal matter, the highest rate of annihilation processes – and hence visible interactions – will be expected from antimatter AQNs. One third of the total AQN mass comes from the N_DW = 1 axion DW (Ge et al. 2018). Therefore – to account for the antimatter AQN number density n_AQN from the total “observable” dark matter density p_DM - we would further have to multiply a factor of 2/3 by excluding the DW contribution to the AQN mass. We therefore obtain an anti-AQN number density estimate of $$n_{\text{AQN}}=\frac{2}{3}\times \frac{3}{5}\times \frac{{\rho }_{\text{DM}}}{M_{\text{AQN}}},$$(4)

which typically scales in a galaxy cluster environment as n_AQN ~ 10^–29 cm^–3. It is because of the low number density that AQNs mostly preserve cold dark matter features, provided they are in a suitable environment. However, the interactions that occur may even yield different emission processes with respect to AQNs. For the sake of simplicity, we further refer to antimatter axion quark nuggets by using the abbreviation “AQN”.

Fig. 1 Schematic sketch for the antimatter-AQN’s internal structure (not true to scale) with an indication of typical sizes (left panel) and the processes, which lead to the different radiation signatures (middle and right panel). The antimatter-AQN consists of a core of anti-quarks surrounded by positrons. When colliding with the ambient gas, the antimatter-AQN heats up and radiation occurs, which depends on the speed at which the collision occurs, the gas temperature and density as well as the ionization state of the gas. The right panel shows the two key emissions, namely thermal and non-thermal emissions. From an incoming proton (green, small circle), 2 GeV of energy will be available after collision with the AQN. Through annihilation, energetic quarks and gluons are produced and stream either deep within the nugget (around 90% of the 2 GeV) or toward the nugget surface (around 10% of the 2 GeV). Quarks and gluons, which remain inside, will contribute to the thermalization of the nugget. Heat will be transferred from the core to the electrosphere which cools down via thermal Bremsstrahlung emission. Quarks and gluons, which proceed to the surface without thermalization, can transfer energy to positrons and produce non-thermal emission.

3 Methodology

3.1 Internal AQN temperatures through ambient gas interactions

The internal temperature of AQNs is influenced by annihilation processes of protons impacting a nugget with anti-quarks from the inner core. The internal temperature of a nugget can be estimated from the conservation of energy. In thermal equilibrium, if a proton is captured by an AQN, the injected energy must be equal to the radiative output. Protons can be captured more efficiently if the anti-quark nugget is ionized, but loses efficiency under certain conditions of the environment. Therefore, we have to introduce an effective cross-section σ_eff to model collision rates. The energy conversion efficiency during a collision depends on a proxy for the probability of quantum reflection for neutral gas called (1 – f), and the fraction of thermal photon production (1 – g) ≈ 0.9. Because of the increased Coulomb potential in an ionized gas, the quantum reflection probability $\mathcal{P}=\left(1-f\right)$ is highly suppressed, and we therefore obtain f ~ 1. The Coulomb potential causes the hydrogen ion to be trapped around an AQN until the final collision and therefore annihilation occurs. The injected energy per unit time can be estimated as $$\frac{\text{d}E}{\text{d}t}=(1-g)f{\sigma }_{\text{eff}}{n}_{\text{gas}}\text{Δ}E\text{Δν}.$$(5)

∆E is the energy that will be released by an annihilation event of a single proton amounting up to ∆E ~ 2m_pc² ≈ 2 GeV and ∆v = |v_AQN − v_gas| is the relative AQN-proton speed. The radiative output is the luminosity of a nugget with radius R_AQN, namely, $$L=4\pi R_{\text{AQN}}^2{F}^{\text{th}}.$$(6)

If we set dE/dt = L and after re-expressing the effective cross-section ${\sigma }_{\text{eff}}=\kappa \pi R_{\text{AQN}}^2$ using the scaling parameter κ, we obtain the expression $$F^{\text{th}}=2\text{GeV}\frac{\kappa }{4}(1-g)f{n}_{\text{gas}}\text{Δν}.$$(7)

As a consequence of the internal temperature increase due to the annihilation process, the positrons in the electrosphere will be heated and emit Bremsstrahlung. In natural units, namely, c = k_B = ħ = 1 and h = 2π, the total surface emissivity of an AQN is calculated by McNeil Forbes & Zhitnitsky (2008b) and given as $$F^{\text{th}}\approx \frac{16{\alpha }^{5/2}}{3\pi }{T}_{\text{AQN}}^4{\left(\frac{T_{\text{AQN}}}{m_e}\right)}^{1/4}.$$(8)

Here, α is the fine-structure constant, and T_AQN is the internal temperature of the AQN. After equating Eqs. (7) with (8), we can solve for T_AQN and obtain $$T_{\text{AQN}}={\left(\frac{3\pi }{4}\frac{2\text{\hspace{0.17em}GeV\hspace{0.17em}}}{16{\alpha }^{5/2}}(1-g)f\kappa m_e^{1/4}{n}_{\text{gas}}\text{Δv}\right)}^{4/17}.$$(9)

Equation (9) can be represented more intuitively if we use galaxy cluster-typical gas densities of n_gas,cluster ~ 10⁻³cm⁻³ and a relative velocity of ∆ν ~ 10⁻³ c: $$\begin{array}{c}{T}_{\text{AQN}}\approx 6.95\times {10}^{-2}\text{eV}{\left(\frac{n_{\text{gas }}}{{10}^{-3}{\text{cm}}^{-3}}\right)}^{4/17}{\left(\frac{\text{Δν}}{{10}^{-3}c}\right)}^{4/17}\\ \times {\left(\frac{(1-g)f}{0.9}\right)}^{4/17}{\kappa }^{4/17}.\end{array}$$(10)

The numerical factor κ boosts the internal AQN temperature if the effective cross-section is larger than its geometrical cross-section. We thus use the following condition for each AQN: $$\kappa =\{\begin{array}{c}{\left(\frac{R_{\text{eff}}}{R_{\text{AQN}}}\right)}^2,R_{\text{eff}}>R_{\text{AQN}}\\ 1,R_{\text{eff}}\le R_{\text{AQN}}.\end{array}$$(11)

To find an expression for R_eff, two important contributions have to be taken into account. First, the charge of an AQN Q, which is a proxy for the ionization stage, and second, the temperature of the surrounding plasma. In galaxy clusters, the ICM consists of highly ionized gas with temperatures of T_ICM ~ keV. In an entirely ionized gas, Zhitnitsky (2023) shows that the cross-section of AQNs can scale with R_eff instead of R_AQN, leading to a potential increase in the collision efficiency with the surrounding gas. It is important to note that geometrical radii of A_QNs are of the order of R_AQN ~ 0.1 µm and – given that the gas number density in galaxy clusters is of the order of n_gas ~ 10⁻³ cm⁻³ – collision rates with surrounding particles can be quite low if R_eff < R_AQN. It is not sensible to incorporate effective radii smaller than their geometrical size, and therefore for R_eff < R_AQN we set R_eff = R_AQN. It is possible to derive R_eff from the physical properties of the nugget and its environment.

In the following Zhitnitsky (2023) proposes that the potential energy of attraction for AQNs in an ionized gas scales with the thermal energy of the plasma in the environment, namely, $$\frac{\alpha Q}{R_{\text{cap }}}\sim \frac{m_p{\nu }_p^2}{2}\sim T_{\text{gas}},$$(12)

with the fine-structure constant, α, the charge, Q, the proton mass, m_p and proton velocity, v_p. Equation (12) shows the radius at which an ionized AQN can capture protons from the gas with thermal velocities of $T_{\text{gas}}~m_p{v}_p^2$. What follows from this scaling relation is a capture radius which decreases for increasing gas temperatures. Zhitnitsky (2023) used the capture radius instead of the effective radius with $R_{\text{cap}}=\sqrt{\pi }{\mathrm{R}}_{\text{eff}}$. In this paper, however, we adopt the convention of using R_eff instead of R_cap. After re-expressing Eq. (12), we obtain the following relation: $$R_{\text{eff}}=\frac{\alpha Q}{T_{\text{gas}}\sqrt{\pi }}.$$(13)

It is important to note that Q is a function of T_AQN. The higher T_AQN, the higher the thermal motion of positrons in the electrosphere with positrons following a non-relativistic Boltzmann regime (McNeil Forbes & Zhitnitsky 2008b). If kinetic positron energies exceed the potential energy of the AQN, positrons can deplete, leading to a negative charge for AQNs. The positrons are rather weakly bound and the number of positrons Q which are likely to evaporate from the AQN can be estimated by using the mean field approximation (McNeil Forbes & Zhitnitsky 2008b) of the local density of positrons n(ɀ, T_AQN) at a distance ɀ from the surface of an AQN. According to Zhitnitsky (2023), one obtains for typical T_AQN in cluster environments $$Q\simeq \frac{4\pi R_{\text{AQN}}^2}{\sqrt{2\pi \alpha }}\left(m_e{T}_{\text{AQN}}\right){\left(\frac{T_{\text{AQN}}}{m_e}\right)}^{1/4}$$(14) $$\approx 4.62\times {10}^3{\left(\frac{T_{\text{AQN}}}{{10}^{-2}\text{eV}}\right)}^{5/4}{\left(\frac{R_{\text{AQN}}}{2.25\times {10}^{-5}\text{cm}}\right)}^2,$$(15)

with m_e and units of cm from Eq. (14) being re-expressed in units of eV using natural units (i.e. c = ħ = 1 and h = 2π). Now, after using Eqs. (9), (13), (15), and κ = (R_eff /R_AQN)², we find $$\begin{array}{c}{R}_{\text{eff}}\simeq 3.89\times {10}^{-7}\text{cm}{\left(\frac{n_{\text{gas }}}{{10}^{-3}{\text{cm}}^{-3}}\right)}^{5/7}{\left(\frac{\text{Δv}}{{10}^{-3}\text{c}}\right)}^{5/7}\\ \times {\left(\frac{R_{\text{AQN}}}{2.25\times {10}^{-5}\text{cm}}\right)}^{24/7}{\left(\frac{1\text{\hspace{0.17em}keV}}{T_{\text{gas}}}\right)}^{17/7}.\end{array}$$(16)

Equation (16) immediately implies that in a galaxy cluster environment, the majority of AQNs will be treated with κ = 1, since R_AQN > R_eff.

3.2 Radiation processes of AQNs

AQNs can emit both thermal and non-thermal radiation. The major fraction (1 − 𝑔) ≈ 0.9 of the annihilation energy of 2 GeV is emitted by thermal Bremsstrahlung. The remaining fraction 𝑔 ≈ 0.1 is radiated non-thermally. The total surface emissivity was already introduced in Eq. (8), which was derived as McNeil Forbes & Zhitnitsky (2008b) $$F_{\nu }^{\text{th}}\sim \frac{8}{45}{\alpha }^{5/2}{T}_{\text{AQN}}^3{\left(\frac{T_{\text{AQN}}}{m_e}\right)}^{1/4}\times f(x),$$(17)

with f (x) being: $$f(x)=(1+x)e^{-x}\{\begin{array}{c}(17-12\ln (x/2)),x<1\hfill \\ (17+12\ln (2)),x\ge 1,\hfill \end{array}$$(18)

and x = 2πν/T_AQN. Majidi et al. (2024) rewrote Eq. (17) in physical units. This also changes x to $\tilde{x}=2\pi \hslash \nu /\left(k_B{T}_{\text{AQN}}\right)$, such that we get: $$F_{\nu }^{\text{th}}=\frac{8}{45}\frac{{\alpha }^{5/2}{k}_B^3{T}_{\text{AQN}}^3}{{ħ}^2{c}^2}{\left(\frac{k_B{T}_{\text{AQN}}}{m_e{c}^2}\right)}^{1/4}\times f(\tilde{x}).$$(19)

Figure 2 shows a plot of dF(v)/dν for a set of internal AQN temperatures for ν ∊ 2.42 × [10⁸,1 0²⁰] Hz. Temperatures that fall into the cluster regime are mostly abundant in the radio and microwave bands. There is another approach (Flambaum & Samsonov 2022a, b) that assumes an emission spectrum completely different from ours. As the total flux intensity must be the same in both approaches, we do not implement this model since it is beyond the scope of this work.

In the case of non-thermal emission, McNeil Forbes & Zhitnitsky (2008a) proposed the non-thermal emission per frequency to scale with: $$F_{\nu }^{\text{non-th }}=\frac{\nu }{{\nu }_c}{{{\displaystyle \int }}^{\text{}}}_{\nu /{\nu }_c}^{\infty }{K}_{3/5}(x)\text{d}x\text{, }$$(20)

with K_3/5(x) being the second modified Bessel function. The critical frequency, ν_c = ω_c/(2π) ≈ 30 keV/(2π), was approximated by McNeil Forbes & Zhitnitsky (2008a). It is important to note that ν_c = 30 keV/(2π), whose numerical value is set to 30 for conventional reasons, is not strictly constrained. For the expression of the synchrotron function, various approximations can be estimated by different approaches and integration methods (see for example Crusius & Schlickeiser 1986; Rybicki & Lightman 1986; Weniger & Cížek 1990; Fouka & Ouichaoui 2013, 2014; Yang & Chu 2017; Palade & Pomârjanschi 2023). For our study, we have adapted the approximation of the integrated modified Bessel function from Majidi et al. (2024): $$F_{\nu }^{\text{non-th }}\approx C(\beta )x^{\beta }{e}^{-x},$$(21)

with C(1/3) ≈ 1.81. It is important to note that β is usually not fixed and can vary depending on the energy distribution of positrons. Additionally, the exponent depends on the emission mechanism at the nugget’s surface, which can be free-free-like emission rather than synchrotron-like, resulting in β = 0. Given the various analytical approximations and the choice of utilizing the synchrotron function to model the non-thermal spectrum requiring more physical evidence, one has to keep in mind that the non-thermal AQN spectrum is not constrained well enough to propose definite predictions. With this in mind, the non-thermal AQN emission has to be understood with caution. Nevertheless, the bolometric surface emissivity of the non-thermal component can then be calculated by substituting x with ν/ν_c in Eq. (21) and integrating over ν ∊ [0, ∞]. The final expression of $F_{\nu }^{\text{non-th}}$ is derived in great detail in Majidi et al. (2024) and is directly adopted from the cited literature: $$F_{\nu }^{\text{non}-\text{th}}=\frac{2\text{\hspace{0.17em}GeV\hspace{0.17em}}fg}{4\text{Γ}(4/3)}\frac{n_{\text{gas}}\text{Δv}}{{\nu }_c}{\left(\frac{R_{\text{eff}}}{R_{\text{AQN}}}\right)}^2{\left(\frac{\nu }{{\nu }_c}\right)}^{1/3}{e}^{-\nu /{\nu }_c},$$(22)

and expressed in terms of T_AQN: $$F_{\nu }^{\text{non}-\text{th}}=\frac{16{\alpha }^{5/2}}{3\pi \text{Γ}(4/3)}\frac{k_B^{17/4}}{{ħ}^3{c}^{5/2}}\frac{g}{1-g}\frac{1}{{\nu }_c}\frac{T_{\text{AQN}}^{17/4}}{m_e^{1/4}}{\left(\frac{\nu }{{\nu }_c}\right)}^{1/3}{e}^{-\nu /{\nu }_c}.$$(23)

A plot of the non-thermal spectrum for typical R_eff, n_gas = 10⁻³cm⁻³ and ∆v = 10⁷cm s⁻¹ in a galaxy clusters is shown in Fig. 3. In most cases, R_eff = R_AQN, and when comparing this to Fig. 2 for typical AQN temperatures in galaxy clusters (T_AQN ∊ [10⁻², 10⁻¹] eV), it can be seen that the non-thermal emission is expected to be a few orders of magnitude lower than the thermal emission.

In a system of multiple AQNs, a more appropriate representation of radiation is the emissivity j_ν, which depends on the number density of AQNs n_AQN, and it is therefore crucial to evaluate individual local number densities instead of assuming a global n_AQN (the denser the environment, the more emission we can expect). It is therefore important to evaluate the emissivity of each AQN tracer particle individually. With the individual emissivity, we can obtain a cluster spectrum or emission maps after summing over all individual j_ν. A more detailed description of how local properties in a particle simulation are inferred will be presented in Sect. 3.3. We adopted Eq. (4) for describing the number density, while we have to keep in mind that ρ_DM,i should be evaluated at the ith AQN tracer particle in the simulation. The emissivity is represented as $$j_{\nu ,i}=4\pi R_{\text{AQN}}^2{n}_{\text{AQN},i}{F}_{\nu ,i}.$$(24)

Fig. 2 Spectral range for thermal emission from AQNs. The magnitude and width of a thermal spectrum increase for increasing TAQN and can reach from radio to X-ray frequencies for extraordinary hot AQNs.

3.3 Mapping AQN properties onto halo tracer particles

Cosmological simulations in the SPH formalism use tracer particles to represent environmental properties in a specific region. Depending on the resolution in our simulation, a tracer has a mass of M_part ~ 10⁹–10¹⁰ h⁻¹ M_⊙. In the underlying AQN model, a particle represents properties of temperature, density, and velocity relative to the environment within the region it exists in. In our analysis of the simulation, particle families are split into gas and halo tracer particles. The motion of halo tracer particles is defined only by their gravitational interactions with their surroundings. It is because of this property that halo tracer particles can serve as good tracers for AQNs in this analysis. As proposed in Eq. (10), physical properties of AQNs rely on the cluster environment – more specifically, on the underlying gas properties of the ICM. It is not trivial to find that gas properties of SPH tracer particles can be mapped onto a single halo tracer particle. Therefore, we aim to describe in greater detail how the mapping process can be numerically implemented.

Initially, the number of neighboring gas particles that can influence the AQN properties is defined to be N_neigh = 200. To reduce the computational costs, we use the package NearestNeighbors.jl (Carlsson et al. 2022) and sort particles by the k-d tree algorithm to find the nearest gas neighbors surrounding an AQN particle.

Once all the gas particles have been localized for each AQN, the relative velocities in each of the x, y, and z-coordinates to all neighboring gas particles can be calculated to infer the speed, which is then saved in an array. Depending on their distance relative to the gas particles, AQNs are characterized by different weights, 𝒲, for each neighboring gas particle. Kernel weights are usually a function of the corresponding kernel function, k_K, the AQN distance to the jth gas particle, d_i,j, and the smoothing length, h_AQN,i = max(d) if d ≠ 0 of an AQN with respect to the jth gas particle; namely, 𝒲_AQN(k_K, d_i,j, h_AQN,i), and can be calculated by using the package SPHKernels.jl (Böss 2023). A 3D Wendland C2 kernel is chosen for AQNs and a 3-dimensional Wendland C4 for gas particles (see Dehnen & Aly (2012) for choosing kernels). When considering the ·th AQN particle, any arbitrary SPH property A_j of the jth gas particle can be calculated according to Dolag et al. (2008) via $$\langle A_i\rangle ={\displaystyle \sum _{j=1}^{N_{\text{neigh }}}\frac{m_{\text{gas},j}}{{\rho }_{\text{gas},j}}}{A}_j{\mathcal{W}}_{\text{AQN}}\left(k_{\text{C}2},d_{i,j},h_{\text{AQN},i}\right).$$(25)

In this case, this is how 〈∆v_i〉 is calculated. However, to properly infer gas properties at positions of AQN particles relative to their environment, it is crucial to define a kernel that computes weights only from the given gas HSML property and their corresponding distances to the AQN. For the gas properties, we therefore obtain the mapped gas particle SPH property onto the AQN via $$\langle A_i\rangle ={\displaystyle \sum _{j=1}^{N_{\text{neigh }}}\frac{m_{\text{gas},j}}{{\rho }_{\text{gas},j}}}{A}_j{\mathcal{W}}_{\text{AQN}}\left(k_{\text{C}4},d_{i,j},h_{\text{gas},j}\right).$$(26)

This finally gives us the following set of equations that have to be solved numerically: $$\langle \Delta {\text{v}}_i\rangle ={\displaystyle \sum _{j=1}^{N_{\text{neigh }}}\frac{m_{\text{gas},j}}{{\rho }_{\text{gas},j}}}\Delta {\text{v}}_j{\mathcal{W}}_{\text{AQN}}\left(k_{\text{C}2},d_{i,j},h_{\text{AQN},i}\right),$$(27) $$\langle n_{\text{gas},i}\rangle ={\displaystyle \sum _{j=1}^{N_{\text{neigh}}}\frac{m_{\text{gas},j}}{{\rho }_{\text{gas},j}}}{n}_{\text{gas},j}{\mathcal{W}}_{\text{AQN}}\left(k_{\text{C}4},d_{i,j},h_{\text{gas},j}\right),$$(28) $$\langle T_{\text{gas},i}\rangle ={\displaystyle \sum _{j=1}^{N_{\text{ne}igh}}\frac{m_{\text{gas},j}}{{\rho }_{\text{gas},j}}}{T}_{\text{gas},j}{\mathcal{W}}_{\text{AQN}}\left(k_{\text{C}4},d_{i,j},h_{\text{gas},j}\right).$$(29)

Here, h_gas,j can be directly obtained from the HSML block of each gas particle. To ensure that numerical values do not fall below sensible SPH values, the conditions $$\langle T_{\text{gas},i}\rangle =\min \left(T_{\text{gas}}\right)\text{ for }\langle T_{\text{gas},i}\rangle <\min \left(T_{\text{gas}}\right)$$(30) $$\langle n_{\text{gas},i}\rangle =\min \left(n_{\text{gas}}\right)\text{ for }\langle n_{\text{gas},i}\rangle <\min \left(n_{\text{gas}}\right)$$(31)

are applied. For reasons of simplification, only the absolute relative velocity for each particle was taken to calculate an SPH representative for the AQN, which can be calculated via $$\Delta {\text{v}}_{ij}=\sqrt{{\left({\nu }_{j,\text{gas}}^x-{\nu }_{i,\text{AQN}}^x\right)}^2+{\left({\nu }_{j,\text{gas}}^y-{\nu }_{i,\text{AQN}}^y\right)}^2+{\left({\nu }_{j,\text{gas}}^z-{\nu }_{i,\text{AQN}}^z\right)}^2}.$$(32)

It is a direct consequence of this approach that 〈∆vi〉 is obtained independently of the direction of the velocity vector of a gas particle relative to the AQN. Only a single scalar is saved to describe the smoothed speed of an AQN with respect to its environment. It is therefore possible that particles are strongly weighted because of their relative distance to the AQN particle but have already passed the AQN and are moving away. This circumstance, on the other hand, is not too big of an issue, since tracer particles are representations of particle distributions, and direct collisions are therefore not desired in the setup of an SPH simulation.

Fig. 3 Non-thermal AQN emission spectrum for different TAQN serving as a validation to reproduce Fig. 3 from Majidi et al. (2024).

3.4 Additional radiation sources in simulated galaxy clusters

3.4.1 ICM Bremsstrahlung

It is important to compare the AQN emission to gas emission to infer AQN excesses and to identify potential emission windows. Amongst other sources of radiation in the X-ray continuum, Bremsstrahlung coming from the hot ICM is the most dominant emission process (Felten et al. (1969), and Sarazin (1986) for an overview). Following Sect. 3.3, two different SPH gas properties are accessible: the gas properties at each AQN position and the gas properties at each gas particle. For gas emission, the latter will be chosen, since the spatial distribution of gas particles differs from the ones of halo tracer particles. It is specifically important to analyze different regions where gas emission is more prominent and where AQN positions might differ. The X-ray emissivity j_X-ray,ν in the energy band [ν₀, ν₁] depends on the hydrogen mass fraction, X_H, the gas temperature, T_gas, the gas density, ρ_gas, on the gaunt factor, 𝑔(ν, T_gas), and the frequency, ν. It takes the following form: $$\begin{array}{l}{j}_{\text{X-ray }}=4C_{j}g\left(\nu ,T_{\text{gas}}\right){\left(1+X_H\right)}^2{\left(\frac{\tilde{n}}{\tilde{f}{m}_p}\right)}^2{\rho }_{\text{gas}}^2\left[\text{g}\hspace{0.17em}{\text{cm}}^{-3}\right]\sqrt{T_{\text{gas}}[\text{keV}]}\\ \times \left(\exp \left(-\frac{h{\nu }_0[\text{keV}]}{T_{\text{gas }}[\text{keV}]}\right)-\exp \left(-\frac{h{\nu }_1[\text{keV}]}{T_{\text{gas }}[\text{keV}]}\right)\right)\text{erg} {\text{cm}}^{-3} {\text{s}}^{-1},\end{array}$$(33)

with $$\tilde{n}=\frac{X_H+0.5\left(1-X_H\right)}{2X_H+0.75\left(1-X_H\right)},$$(34)

and $$\tilde{f}=\frac{4}{5X_H+3},$$(35)

where X_H = 0.76 is a parameter of the simulation, and C_j = 2.42 × 10⁻²⁴ is a numerical factor provided by Bartelmann & Steinmetz (1996). Since a detailed comparison of the AQN and gas spectra requires the emissivity per unit frequency and Eq. (33) only provides the emissivity in an integrated form, we must take the derivative per frequency bin to obtain the units of erg cm⁻³Hz⁻¹ s⁻¹. In the numerical implementation, a frequency array of 1000 elements was defined in the frequency range of ν ∊ [10⁻⁵, 10⁸]/(2π)eV given in log-space. Each consecutive element has an increased factor of a = ν_i+1/ν_i, which is constant for any i. For the factor a, we can therefore choose i = 0 for simplification, namely, the first element in the frequency array. Then, ∆ν_i = av_i − ν_i is then the ith frequency band, providing a sufficient approximation for the derivation with: $$j_{\text{X}-\text{ray},\nu }(\nu )=\frac{j_{\text{X}-\text{ray}}(\nu )}{\Delta {\nu }_i}.$$(36)

3.4.2 Estimating the synchrotron emission of the ICM

There are various sources of synchrotron emission within galaxy clusters. For a detailed, recent review on radio emission from galaxy clusters, see (van Weeren et al. 2019). Most importantly for our work, some clusters show large-scale volume filling radio emission that roughly follows the distribution of the ICM, so-called giant radio halos. Their exact nature is still under debate and theoretical models involve turbulent re-acceleration processes and/or secondary electrons from hadronic interactions. Other sources like peripheral radio emissions, so-called relics that are related to shocks or radio emission related to individual galaxies, or active galactic nuclei and their jets will be neglected for this first study. Given the large uncertainty involved when directly modeling the radio emission from first principles from our simulation (see discussion in Böss et al. 2023), we decided to follow to more conservative approach to producing the possible synchrotron signal from the ICM when comparing to the predicted AQN signal. Here, we made use of the fact that observations of a large sample of galaxy clusters suggest that, compared to the Coma galaxy cluster, other clusters typically show less or even no diffuse radio emission (Hanisch 1982; Cuciti et al. 2021). In addition, due to the quasi-near linear correlation observed between the radio signal and the Sunyaev-Zeldovich effect (Planck Collaboration Int. X 2013) we can directly relate the monochromatic radio emissivity to the pressure of the ICM. The radio halo of Coma is also observed at various frequencies so that an overall spectrum can be used (Thierbach et al. 2003). Making use of having a very good counterpart of Coma in our constrained simulation (Dolag et al. 2023; Hernández-Martínez et al. 2024) we, therefore, assume a scaling of the pressure within the simulation with the radio emissivity at 352 MHz, which reproduces the scaling presented in Planck Collaboration Int. X (2013) and modified at different frequencies to reproduce the spectrum measured in Thierbach et al. (2003). This reproduces the observed scaling relation between mass and radio power at 1.4 GHz very well (Cuciti et al. 2021); although with this approach, every cluster shows a giant radio halo. Thus, it gives an upper limit of the expected synchrotron emission from the ICM in the simulations.

3.5 The constrained simulation of the local universe “SLOW”

The underlying simulation is a constrained cosmological simulation of the local universe, called Simulating the LOcal Web (SLOW, Dolag et al. 2023). Derived from the CosmicFlows-2 catalog (Tully et al. 2013), the peculiar velocity field can help to constrain initial conditions for a simulation of our local universe, when tracing back the trajectories up to a defined redshift. Tracers, such as the velocity and density fields, are direct proxies for the underlying gravitational fields which are important for the initial conditions of the simulation. This study uses the constrained cosmological simulated simulation, covering a volume of (500 h⁻¹Mpc)³ and including 3072³ gas and dark matter particles. It assumes a cosmology based on Planck measurements (Planck Collaboration XVI 2014), with a Hubble constant H₀ = 67.77 km s⁻¹ Mpc⁻¹, a total matter density Ω_M = 0.307115, a cosmological constant Ω_Λ = 0.692885, a baryon fraction Ω_B = 0.0480217, a normalization of the power spectrum σ₈ = 0.829 and a slope of the primordial fluctuation spectra n = 0.961 as described in Dolag et al. (2023).

3.6 The sample

In this work, we tested the spectral emission of AQN quantitatively for a large sample and qualitatively for a small sample. A large sample of 150 galaxy clusters (which we call “sample 𝒜”) ordered in five corresponding mass bins was used to test scaling relations and mass dependencies for the search of emission excesses in the spectra. Each mass bin was defined to contain 30 galaxy clusters with bin sizes of M_{vir, 1} ∊ [0.8, 0.9] × 10¹⁴M_⊙, M_{vir, 2} ∊ [1.1, 2.0] × 10¹⁴M_⊙, M_vir,1 ∊ [4.0, 5.0] × 10¹⁴M_⊙, M_vir,1 ∊ [7.0, 7.9] × 10¹⁴M_⊙ and M_vir,1 ∊ [10.7, 31.7] × 10¹⁴M_⊙. Mass bins were arbitrarily chosen to maintain a relatively equal mass distribution and conserve a well-distributed population of virial masses. The lowest mass range was chosen to be 0.8 × 10¹⁴M_⊙ since it becomes tricky at lower masses to separate galaxy clusters from galaxy groups. As an upper mass range, the 31.7 × 10¹⁴M_⊙ were obtained by taking the most massive galaxy clusters above the boundary of 10¹⁴M_⊙. We used the subfind algorithm (Springel et al. 2001; Dolag et al. 2009) implemented into the GadgetI0.jl package (Böss & Valenzuela 2023) to receive the necessary halo IDs to read out the halo and gas tracer particles within a certain box.

“Sample ℬ” consists of 11 out of 45 cross-identified galaxy clusters (Hernández-Martínez et al. 2024) that resemble digital twins to their real counterparts, which was then utilized to identify the most promising galaxy cluster candidates that exhibit AQN signatures. According to Hernández-Martínez et al. (2024), galaxy clusters from observations were cross-identified with simulated galaxy clusters by comparing their mass, X-ray luminosity, temperature, and Compton-y signal. A probability of cross-identification can be established by measuring the probability of obtaining each cross-match in a random, unconstrained simulation. By computing the significance, the authors could assess how well the cross-identified clusters are associated with each other. The most massive galaxy cluster from the digital twin-sample is Coma (M_vir = 18.82 × 10¹⁴ M_⊙) and the least massive is Fornax (M_vir = 0.61 × 10¹⁴ M_⊙).

Galaxy clusters from both sample 𝒜 and sample ℬ were read out from a sphere of radius r = b × r_vir, where b = 1.5 is a chosen multiplier used to extend the sphere.

4 Results

4.1 Environmental impact onto physical properties of AQNs

The histograms presented in Fig. 4 show the distribution of T_AQN, Δv, n_AQN, R_eff, n_gas, and T_gas, respectively. Out of the distributions of all 150 galaxy clusters, the median distribution in each mass bin (namely for 30 galaxy clusters each) was taken and plotted. It can be seen that the larger the range gets, all properties start to show a relatively wider distribution with fewer outliers. Since the cluster mass directly scales with the particle numbers, a general trend of increased total counts can be seen in all histograms; however it is interesting to see how the shape of the histograms changes for increasing mass bins.

• Relative velocity: in the Δv-histogram, most of the particles are distributed towards Δv ~ 10⁸cms⁻¹ ~ 0.3%c. Intriguingly, an additional low-velocity population component of Δv ≲ 10⁶cms⁻¹ seems to be present, especially in massive galaxy clusters. A Δv peak shift towards higher velocities for larger cluster masses can also be observed.

• Internal AQN temperature: the T_AQN-histogram is heavily weighted for values of T_AQN ∈ [10⁻²,10⁻¹]eV and occasionally high T_AQN ≳ 1 eV can be possible for massive clusters. Peaks in the distribution shift from approximately 4 × 10⁻² eV to 6 × 10⁻² eV for increasing cluster masses. Here, Δv is likely to play the dominating factor (as discussed before). On the low T_AQN -side, fewer counts are caused by non-optimal combinations of Δv, n_gas and T_gas. The higher T_AQN -side is limited by an increasing amount of high counts of high T_gas .

• Effective AQN radius: the R_eff-histogram directly shows that most of the R_eff are set to R_AQN = 2.25 × 10⁻⁵ cm, since κ = 1 in most environments. Depending on the cluster mass, the general shape of the R_eff distribution does not drastically change.

• AQN number density: following Equation (4), values in the n_AQN -histogram are generally low, since AQNs are large and massive compared to elementary particles. Therefore, a smaller abundance of AQN is required, which results in low n_AQN.

• Gas number density: n_gas shows a reasonable distribution with a heavier weight in number densities of n_gas < 10⁻³cm⁻³ because of the extended readout-radius of 1.5 r_vir.

• Gas temperature: For increasing M_vir in galaxy clusters the peak in the T_gas-histograms shifts to larger values. Higher gas temperatures for heavier galaxy clusters can be expected due to a steeper potential that in return attracts more substructures. A constant infall of new galaxies continuously heats the ICM and provides a general increase in T_gas (see Sarazin 1986, for a review).

A better understanding of the spatial distribution of AQN properties in galaxy clusters can be obtained when considering their radial profiles. Figure 5 shows how Δv, T_AQN, and R_eff, mapped onto the position of AQNs scale with increasing distances from the cluster’s center of mass. Radial profiles were taken from each galaxy cluster in sample 𝒜 and the median radial profile was taken for each mass bin. SPH properties were inferred according to the proposed methodology in Sect. 3.3.

• Relative velocity: An interesting outcome of the radial Δv profiles is that for the most massive galaxy clusters, the profile increases towards normalized distances of r/r_vir ~ 0.2 and decreases at larger distances. All mass bins show relatively flat radial profiles until r/r_vir ~ 0.1.

• Effective AQN radius: R_eff represents an intriguing property of AQNs since (due to the high ICM temperature in the central regions) the gas properties do not permit effective radii to be larger than R_AQN . At larger distances, individual regional properties enable small increases in R_eff that highly depend on the degenerate AQN properties present in the surrounding gas. Single substructures might be capable of obtaining large effective radii, and it is subject to Sect. 4.3 to identify these regions in spatially resolved maps of individual galaxy clusters to properly infer their physical properties. It is important to note that a slight trend of increasing R_eff can already be observed at r ≤ 1.5 r_vir. Even though the increase is almost negligible with a maximum increase of R_eff hosted by the lowest cluster mass bin of ΔR_eff ≃ 3 × 10⁻⁷ cm (namely, an increase of ~1.3%), the phenomenon of higher R_eff in the peripheries becomes more significant the further the distance is. This is a feature of the underlying AQN model that was anticipated for this paper, since AQNs are expected to be surrounded by a fully ionized gas. In regions of r ≫ 1.5 r_vir, smaller T_gas values consequently yield large R_eff, as ionized particles can be captured over larger distances if the thermal motion of the ions is sufficiently low. This will lead to a runaway effect for increasing R_eff that is not physically meaningful. It is, therefore, crucial to constrain the sample to galaxy cluster regions with an extent of approximately r_max ~ 1.5 r_vir at maximum for the underlying AQN model.

• Internal AQN temperature: T_AQN are most likely increasing for different masses with higher Δv. For larger radii and in all mass bins, T_AQN seems to converge towards a single radially normalized value. Given the underlying AQN model, it is assumed that T_AQN will artificially rise for radii r ≫ 1.5 r_vir, because R_eff → ∞ for T_gas → 0 K. Moreover, this means that at a certain relative distance, T_AQN , will mainly be influenced by the properties of filaments and voids.

4.2 Spectral energy distributions

A spectral energy distribution can be obtained by using Eqs. (19), (23), and (33) for thermal and non-thermal AQN emission, and thermal Bremsstrahlung emission from the ICM, respectively. For synchrotron emission, we followed the steps described in Sect. 3.4.2. We note that we have to apply Eqs. (24)–(19), and Eq. (23), to obtain an emissivity in units of erg cm⁻³Hz⁻¹s⁻¹.

The emissivity is then calculated for each particle in each frequency in a range of ν ∈ [10⁻⁵,10⁸]/(2π)eV with a spectral resolution of 10³ bins for AQN and gas. For the CR synchrotron emission, frequency ranges of ν_synch ∈ [10⁻⁵,10⁻²]/(2π) eV were chosen in a spectral resolution of 100 bins. Figure 6 displays the spectral energy distribution of the five mass bins from sample 𝒜: here, the median for each frequency in each mass bin was taken out of the 30 respective galaxy clusters. AQN thermal emission, non-thermal emission, ICM Bremsstrahlung, and CR synchrotron emission are represented by solid, dash-dot-dotted, dashed, and dash-dotted lines, respectively.

It is noticeable that higher cluster masses result in an increment of the general spectral energy. In this figure, our attention is drawn towards the regions where both of the AQN emission become comparable to the ICM Bremsstrahlung emission in the low- and high-energy regimes, as it is possible to obtain small frequency windows where AQN emission would overtake the ICM Bremsstrahlung. In the low-energy regime of ν_low ≲ 10⁻³eV, AQN and ICM emission roughly increase equally for higher M_vir. This opens the possibility that thermal AQN radiation has the potential to outshine ICM Bremsstrahlung independently of cluster masses. The fact that AQNs can deliver a Bremsstrahlung intensity similar to hot X- ray gas in the radio regime is particularly interesting since in contrast to T_gas, T_AQN is not as strictly related to n_gas, which highly depends on the equation of state and the equilibrium condition.

On the other hand, we would like to point out that the spectrum was extended towards minimum frequencies of ν_min = 10⁻⁵/(2π) eV to provide a low-frequency coverage by the synchrotron emission. For AQNs, however, it can become challenging to provide a physically meaningful interpretation for frequencies lower than ν_min, because of the finite size effect (see McNeil Forbes & Zhitnitsky 2008b). Here, the plane wave approximation breaks down, as the positron wavelength becomes comparable to the size of the nugget.

This becomes a bit trickier when comparing non-thermal AQN emission to the high energy Bremsstrahlung regime in frequencies of ν_high ≳ 10 keV. Non-thermal AQN emission will increment not as strongly as X-ray emission of the hot cluster gas. Since $j_{\text{X}-\text{ray}}~{\rho }_{\text{gas}}^2{T}_{\text{gas}}^{1/2}$ (see Eq. (33)) and T_gas increases with M_vir (see Fig. 4), j_X-ray,ν will especially increase drastically for increasing M_vir in the high energy regime, leading to the circumstance that non-thermal AQN is expected to only be visible in sufficiently low cluster masses.

The contribution in synchrotron emission in Fig. 6 is increasing for higher M_vir because of the Compton-y parameter. Synchrotron emission is likely to interfere with the thermal AQN emission, especially in the regions where thermal AQN emission would dominate the ICM Bremsstrahlung emission. It is therefore worth noting that one has to account for synchrotron emission and ICM Bremsstrahlung in the low-energy regime when searching for an optimal frequency window for thermal AQN emission.

In the following, we focus on the promising frequency bands where AQN emission is expected to dominate. First, to measure how strongly AQN emission outshines the background emission a metric has to be defined that includes the frequency and the magnitude of the increment in emission. Depending on the mass and other physical properties of the galaxy cluster, thermal and non-thermal AQN emission will distinctively be higher or lower in individual frequencies for different galaxy clusters. Therefore, a transition frequency ν_T has to be defined, where the absolute summed emissivity of all AQN tracer particles per cluster per frequency (j_AQN,ν) transitions from higher to lower values compared to the thermal gas emission. From a minimum frequency ν_min to the transition frequency all background emission sources (including synchrotron emission) will be integrated. In the simulation, the numerical integration is conducted by applying the trapezoidal rule which is implemented in the NumericalIntegration.jl package⁴.

Fig. 4 Histograms of physical properties relevant to the AQN emission. Median values were taken for each mass bin, i.e. for 30 galaxy clusters, respectively.

Fig. 5 Radial profiles of Δv, Reff, and TAQN mapped onto AQN tracer particles. The median radial profile was taken for each mass bin (i.e. 30 galaxy clusters each) and compared to other mass bins in a region of r ∈ [0,1.5] rvir.

Fig. 6 Median cluster spectra with lines for each mass bin out of sample 𝒜 over a frequency range of ν ∈ [10−5,108]/(2π)eV. The solid, dash-dot-dotted, dashed, and dash-dotted lines represent AQN thermal emission, non-thermal emission, ICM Bremsstrahlung, and synchrotron emission, respectively.

4.2.1 Low-frequency regime

For the low-energy regime, we define ν_min as the lowest possible frequency that was used to calculate the spectra, namely, ν_min = 10⁻⁵/(2π)eV because of the finite size effect of AQNs. After integration from ν_min to ν_T, all background sources will be summed to represent a total galaxy cluster emission. The integrated spectral ratio of AQNs $\tilde{J}$ significantly contributing to the background emission can be calculated as $$j_{\text{tot }}=\underset{=:j_{\text{AQN}}}{\underbrace{{\displaystyle {\int }_{{\nu }_{\min }}^{{\nu }_T}{j}_{\nu ,\text{AQN}}} \text{d}\nu }}+{\displaystyle {\int }_{{\nu }_{\min }}^{{\nu }_T}{j}_{\nu ,\text{gas}}} \text{d}\nu +{\displaystyle {\int }_{{\nu }_{\min }}^{{\nu }_T}{j}_{\nu ,\text{ synch }}}\text{d}\nu ,$$(37)

and $$\tilde{J}=\frac{j_{\text{AQN}}}{j_{\text{tot}}}.$$(38)

In the following analysis, $\tilde{J}$ was determined in both the low- and high-energy regimes for sample 𝒜 and sample ℬ. After applying the condition that thermal AQN emission must exceed thermal gas emission at frequencies below ν_T, approximately 40.99% of all galaxy clusters combined from sample 𝒜 and sample ℬ exhibit stronger AQN emission compared to the thermal background gas emission from ν_min to ν_T.

Using this subset of galaxy clusters and applying Eq. (38), produces Fig. 7, which displays the integrated spectral ratio as a function of the transition frequency for both frequency windows where AQN emission exhibits an excess. On the left side of the figure, the thermal regime of the AQNs dominates; whereas on the right side, at high-energy frequencies, the non-thermal AQN regime dominates. An important outcome is that galaxy clusters with relatively low M_vir are still capable of reaching high $\tilde{J}$ and follow relatively well-constrained tangents indicated by the colormap for different cluster masses.

It is important to find out how cross-identified galaxy clusters populate in this relation to identify promising real-world cluster candidates. Therefore, it is shown in Fig. 7 how sample ℬ is distributed in the ${\nu }_T-\tilde{J}$ diagram. It can be seen in the low-energy regime that Fornax and Virgo appear to be the best candidates for the strongest AQN emission offset for the background gas emission while extending to the highest ν_T. Even though it is the most massive galaxy cluster in sample ℬ, the Coma cluster does not meet the condition to exhibit stronger AQN emission compared to the background emission in the energy window from ν_min to ν_T, and is therefore not displayed in Fig. 7. Furthermore, even if synchrotron emission is excluded in Eq. (37), AQN emission does not dominate the background emission at any frequency for Coma. For a comparison of how AQN emission from other crossidentified galaxy clusters performs with and without synchrotron emission included in the background emission, see Table 1.

Fig. 7 Transition frequency-integrated ratio-relation, color-coded by Mvir with frequency ranges of positive AQN contributions. The low-energy frequency range shows in which the thermal AQN emission dominates. Low cluster masses tend to reside on a tangent, slightly offset towards higher values. The high-energy frequency range shows in which the nonthermal AQN emission dominates. Only AQN emission, originating from the lowest Mvir shows an excess over the thermal gas emission in the high-energy regime. We removed one cluster from the plot where for numerical reasons we could not compute νT and the associated $\tilde{J}$ reliably.

4.2.2 High-frequency regime

Since non-thermal AQN emission is expected to show two intersections in the high energy regime, one has to define ν_min to be the first transition frequency where j_{AQN,nonthermal,ν} > j_{background,ν} and ν_T is again where j_{AQN,nonthermal,ν} < j_{background,ν} with the condition that ν_min < ν_T. Non-thermal AQN emission causing an excess over the background emission in the high energy regime is mostly dictated by the individual X-ray brightness of a galaxy cluster (see Fig. 6).

It can be seen on the right-hand side of the broken axis in Fig. 7 that only a small fraction of galaxy clusters could potentially dominate background gas emission with AQN emission. When applying the conditions to determine an excess in the high energy regime by AQNs, 9.32%, namely, a total number of 15 out of both samples was identified. Figure A.1 depicts the frequency window in which thermal gas emission is expected to intersect with the non-thermal AQN emission. In this image, dashed lines represent the Bremsstrahlung emission from the ICM and dash-dot-dotted lines are the non-thermal AQN emission. If there is an excess occurring, the signature is only very subtle and many low-mass galaxy clusters are still too bright in X-rays to permit a window for non-thermal AQN emission to dominate. When comparing both regimes in Fig. 7, it can be seen that in the high-energy regime, an AQN excess can only be visible for low cluster masses, which then results in a dimmer gas emission. Galaxy clusters in a higher dynamical state of dark matter particles to the surrounding gas particles could additionally positively contribute to a non-thermal AQN-excess, since j_{AQN,nonthermal,ν}~n_AQNn_gas∆v.

It should be noted that the non-thermal emission model is not as robust as the thermal emission model of the AQNs. In future studies, the exponent of 1/3 in Eqs. (22) and (23) for the frequency, will most likely be replaced by a free parameter β, which would be allowed to remain adaptive. The exponent of 1/3 is an approximate solution to solve the synchrotron function. It is crucial to note that particular solutions for the synchrotron function to model the non-thermal spectrum have to be motivated in greater detail and new precise models should aim for a more detailed physical implementation. In addition, the critical frequency as proposed in Sect. 3.2 is set to 30 keV by convention and should be understood in an “order of magnitude estimate” (McNeil Forbes & Zhitnitsky 2008a). To provide improved estimates, a refined theory of the non-thermal physics in the structure of AQNs is desired. The current model utilizes an approximate treatment in the mean-field of the Maxwell equations. Therefore, it is not yet possible to discard the possibility of detecting AQN signatures in the high energy regime immediately and a final conclusion shall be made once a better non-thermal AQN model is developed.

Fig. 8 Ratio map expressed in $\tilde{J}$, shown in the left column. The thermal gas, AQN, and synchrotron emission of Coma integrated over 1.40 GHz to 4.85 GHz, projected 3 rvir along the line of sight are given in the right column (from top to bottom). In this frequency window, the synchrotron emission of Coma is well constrained through radio emission (Thierbach et al. 2003). Due to the strong radio emission of Coma, here the contribution to the diffuse emission by AQNs is in the sub-percent level.

4.3 Cluster maps

It is important to note that one loses information when only considering spectral features. Since AQNs spatially follows the halo distribution, positions of dark matter with respect to gas might slightly differ from each other. Extremely low emissions do not contribute to the spectrum, because j_ν was summed over all particles in the corresponding frequency bin, so spatial regions of low emission are not resolved in a spectrum. It is especially important for the gas emission to identify regions of low emissivity where AQN emission would be proficient to outshine thermal gas emission. It can even be possible to find an AQN excess in small regions in spatially resolved maps where one did not expect an excess only by considering the spectrum (see, for example, the Coma cluster presented in Fig. 8). On the other hand, it is not possible to look for excess features in all 161 galaxy clusters by finding their individual optimal region where AQN features dominate over the background emission when scanning through a spectral resolution of 1000 frequency bins. We, therefore, stick to the methodology presented in Eq. (37) to analyze cluster regions, where AQN radiation certainly dominates in the spectrum by inferring the emission offset and only consider promising candidates of cross-identified galaxy clusters. Note that for simplicity the maps for all clusters are done in super- galactic x/y coordinates and not rotated to match a real observer in the local universe.

Given this reasoning, in this section, the SPH maps of gas and AQNs are represented using their most contributing physical properties and fluxes integrated over an adaptive best possible frequency range. Figures 9 and 10 display the Fornax and Virgo clusters in a selection of relevant properties.

Relative pressure maps were obtained by generating SPH maps by utilizing the ideal gas relation, p_gas ~ n_gas T_gas. Hence, the product of n_gas and T_gas will represent the scaling of the pressure in the gas particles, and SPH maps were constructed this way. With the obtained pressure map, a Gaussian filter was applied to the projected bitmap of the galaxy cluster with σ = 5 in the Gaussian kernel. The relative pressure change is therefore obtained via $$\delta p=\frac{p_{\text{original }}-p_{\text{gaussian }}}{p_{\text{gaussian }}}.$$(39)

In Fornax (see Fig. 9), fewer filamentary structures with strong emission are visible in thermal AQN emission compared to thermal gas emission because dark matter is not affected by environmental friction effects that would perturb the particle distribution. In the map depicting $\tilde{J}$, the strongest contributions originate from a larger infalling substructure, which is visible in the effective radius map as well. The $\tilde{J}\text{-map}$ specifically shows the presence of the glowing axion quark nugget component, which seems to show strong emissions in regions of the ICM.

Prominent regions of R_eff > R_AQN are only possible if the gas environment provides the necessary parameter combination. Compared to other cross-identified galaxy clusters, it is rather rare to see large regions with high R_eff values. The effective radius increments with increasing n_gas and ∆v and decrements for increasing T_gas. T_gas scales with the largest exponent and therefore shows the strongest influence on R_eff with especially high gas temperatures in the central regions. Fornax is a relatively small cluster with a comparable low thermal gas emission coming from a low T_gas. R_eff can therefore reach values of R_eff > R_AQN in central regions, too. However, when comparing R_eff maps in Fornax to Virgo and Coma, we can see that R_eff > R_AQN values are more likely to be found in the peripheries of galaxy clusters.

An intriguing physical property is that regions of R_eff > 10⁻⁴ cm are typically not abundant in relative velocity maps, and it is especially in these regions where ∆v seem to be lower compared to their nearest surroundings. Following $R_{eff}\sim n_{\text{gas}}^{5/7}\Delta {\text{v}}^{5/7}{T}_{\text{gas}}^{-17/7}$ (cf. Eq. (16)), gas overdensities in conjunction with low T_gas yield higher R_eff. In the case of Fornax, the prominent R_eff > R_AQN-region is too close to the cluster center and will therefore be hotter than in the outskirts. In this case, n_gas would be the only free parameter that could directly influence R_eff.

Regions of high T_AQN do not necessarily follow the n_AQN distribution, since T_AQN is also influenced by the surrounding gas. Shocked regions are not significantly embedded in the AQN emission map and vice versa. A further interesting feature is that the AQN emission does not seem to show a significantly stronger abundance at extended regions. This is not trivial since dark matter particles do not suffer from friction and do not interact other than with gravity in the simulation. Therefore, it would have been sensible to expect additional emission contributions in the periphery of Fornax that are not strongly present in the gas emission maps.

The Virgo cluster (see Fig. 10) differs strongly from Fornax as it is more massive, hotter, and in a more relaxed state. No signs of a current major merger are visible in the relative pressure maps (other than visible in Fornax with its prominent northern bow shock). And yet, both of these clusters independently show the strongest AQN excess. It is not obvious that shock features are no tracers for AQN emission as shocked regions are always accompanied by significant offsets in thermodynamic properties (e.g., density or temperature). Even though it was expected that at least shocked regions would play a diminishing or reinforcing role on the AQN emission, this could not be directly confirmed for any of the galaxy clusters in this study. Consequently, shock maps of other cross-identified galaxy clusters do not show any signs of shock signatures in the AQN maps either.

In the $\tilde{J}$-map of the Virgo cluster, we can see that substructures show a mix of stronger and weaker contributions from the AQN emission to the ambient background emission. High ∆v in the outskirts of Virgo have an impact on larger values of $\tilde{J}$ in the peripheries. This direct influence, however, is not obvious in the case of Fornax. In the Virgo cluster, the strongest contribution of $\tilde{J}$ seems to originate from the regions of the ICM instead of the infalling substructures. This suggests we ought to look for AQN signatures where no galaxies are located when searching for signatures outside the center of the galaxy cluster.

In both clusters, we can see that high ∆v can also be reached in the clusters’ peripheries. The distribution of strong ∆v appears to be anisotropic and AQN emission does not necessarily follow the distribution of high relative velocities on large scales.

To learn how AQN signatures evolve in bands of higher energies, individual energy ranges of four different instruments were taken. AQN and gas emission were integrated within the respective energy ranges and AQN emission was superimposed by its thermal and non-thermal components. For the sake of versatile energy coverage, the instruments Wilkinson Microwave Anisotropy Probe (WMAP; Bennett et al. 2013), Planck (Tauber et al. 2010), Euclid (Laureijs et al. 2011), and the X-Ray Imaging and Spectroscopy Mission (XRISM; Mori et al. 2022) were considered. Table 2 shows the corresponding maximum energy range with upper and lower limits denoted by ν_min and ν_max for each instrument in units of GHz.

In the following, the two most promising galaxy clusters, determined by the top figure in Fig. 7 were selected for a multiband analysis. Figures 11 and 12 show the gas and AQN emission with the corresponding ratio image in the selected bands. To compare how the brightness of the emission evolves, all frames in their corresponding emission feature share the same colorbar limits.

While the gas emission shows a continuing brightening for increasing band energies, the AQN signatures show a signiflcant decrease in brightness in the Euclid band. This decrement originates from the spectral feature that thermal AQN emission transitions into the non-thermal regime. Euclid operates in the energy range, where the thermal emission experiences its cutoff (cf. Fig. 6). Even though AQN signatures show a re-brightening once the non-thermal regime takes over, AQN emission cannot compete with the thermal gas emission in the X-ray regime. It is evident in the ratio-images that AQN emission would only be able to dominate in the low-energy regime. A re-brightening in the higher energy bands is only visible in infalling substructures, however, the ambient glow is no longer as abundant.

We may conclude, thus, that Euclid and XRISM are the least promising instruments for a potential AQN signature detection as the thermal gas emission dominates over the AQN emission in this energy range. Furthermore, WMAP and Planck, exhibit AQN signatures in the ratio-images, even when combining thermal gas emission with synchrotron emission.

Fig. 9 Cross-identified Fornax cluster represented in different relevant properties.

Fig. 10 Cross-identified Virgo cluster represented in different relevant properties.

Fig. 11 Fornax cluster AQN, gas, and synchrotron emission generated from WMAP, Planck, Euclid, and XRISM energy ranges. Synchrotron emission is predominantly abundant in WMAP and Planck bands. Thermal gas emission shows the most significant contribution to the background emission in the Euclid and XRISM energy bands.

Fig. 12 Virgo cluster AQN, gas, and synchrotron emission generated from WMAP, Planck, Euclid, and XRISM energy ranges. Synchrotron emission is predominantly abundant in WMAP and Planck bands. Thermal gas emission shows the most significant contribution to the background emission in Euclid and XRISM energy bands.

4.4 Observational feasibility of promising cluster candidates

It is clear from the discussion in Sect. 4.3 that the spatial features of strong AQN signatures are hard to identify and differ from cluster to cluster and their dynamical states. Even though AQN emission excess is abundant in galaxy clusters, it becomes challenging to pinpoint regions of clusters that are expected to be the most promising.

In observations, galaxy cluster emissions are typically analyzed by radial profiles in a given frequency band. This approach enables the identification of important spatial features that are dependent on distance and the corresponding frequency band. In this section, the radial profiles are shown for the cross-identified galaxy clusters that show an excess signature of thermal AQN emission integrated over their individual frequency bands from ν_min to ν_T. Radial profiles in 50 bins were taken for each emission property and a final radial profile was calculated after summing the AQN, gas, and synchrotron maps.

Figure A.2 shows how the $\tilde{J}\left(r\right)$ evolves for the cross-identified galaxy clusters over a radius of r ∊ [10^-2,1.5] r_vir. An interesting feature in the radial profiles is that galaxy clusters exhibiting high $\tilde{J}$ in Fig. 7 also show high ratios in their central regions. In comparison to the remaining cross-identified galaxy clusters, Fornax and Virgo stand out in Fig. A.2.

In all radial profiles, the most significant contribution of thermal AQN emission comes from the peripheral regions, where thermal gas and synchrotron emissions have sufficiently decreased.

5 Discussion

Numerical analysis of AQNs in galaxy clusters in a cosmological simulation has shown that direct signatures of this new dark matter model are possible, but can be challenging to pinpoint. The spatial distributions of significant AQN emission contributions depend on factors, like the dynamical cluster state and thermal AQN emission is spectrally superimposed with non-thermal synchrotron and thermal gas emission.

Given that the spectral AQN signature is not exactly Bremsstrahlung and differs from synchrotron, component separation (where templates are used and spectral signatures help disentangle sources of emission) could be a reasonable strategy to look for weak AQN contributions. Multifrequency data would especially improve this approach (see Majidi et al. 2024 for a discussion).

In Sect. 2.3, tracers of AQN signatures, and further unexplained observational phenomena were pointed out and analyzed by various studies. The study by Lawson & Zhitnitsky (2019) proposed that thermal AQN emission might be responsible for a stronger absorption feature in the 21 cm line at a redshift of z = 17 during the early universe. Clusters of galaxies have not formed at this time yet and early protoclusters – observationally confirmed to be progenitors of present-day galaxy clusters – are detected at z ~ 6–7 (Harikane etal. 2019). Cosmological, hydro-dynamical simulations suggest that a hot ICM is already present at z = 4.3, at least within the core of such a protocluster (Remus et al. 2023). We therefore cannot assume the same conditions for the AQN increment which was obtained in the framework of this paper. However, it is worth mentioning that protoclusters at z ~ 7 could be capable of influencing the radio background if the AQN signal would be similarly significant as in galaxy clusters at present time. An increment in AQN emission was concluded in Majidi et al. (2024) after constructing a light cone to a redshift of z = 5.4 using the Magneticum⁵ simulation (Dolag et al. 2016). The sum of the AQN emission increment in the radial profiles from the cross-identified galaxy clusters of our work at z = 0 can be attributed to $\sum \tilde{J}\left(r\right)\approx 4.80\%$ thermal AQN emission.

It is, however, important to mention that early clusters possess different physical properties (Chiang et al. 2013; Overzier 2016), which would likely yield different AQN signatures. Results from this study only provided information of AQN signatures from cluster properties at z = 0. It would be sensible to implement AQN features already in snapshots at higher redshift and develop a dynamical AQN-treatment, to grasp a full evolution of AQN signatures from early to present-time galaxy clusters that could contribute to the background radio emission.

The study by Zhitnitsky (2022a) analyzed the diffuse galactic UV radiation and suggested that hot AQNs with T_AQN ~ 5 eV could reproduce such observations – however since our study focuses on scales of galaxy clusters and not galactic scales, we cannot verify these predictions. The same argumentation is applied to the proposal to the contribution of Chandra’s diffuse 8 keV emission by AQNs (McNeil Forbes & Zhitnitsky 2008a). However, it is worth commenting that in a galaxy cluster, it is rather unlikely to find AQNs with temperatures of 5 eV because of the suppression of capturing protons in the hot ICM.

Galactic environments differ strongly from galaxy cluster environments with different T_gas, n_gas, n_DM, and Δv as the interstellar medium is cooler and denser than the ICM. In addition, galaxies exhibit a significantly lower gas emission, and especially non-thermal AQN emission in the high energy regime might indeed be capable of dominating highly energetic galactic background emission.

It is not possible to comment predictions on the 511 keV, as predicted by Oaknin & Zhitnitsky (2005); Zhitnitsky (2007); Forbes et al. (2010); Flambaum & Samsonov (2021) as we did not consider e⁺e^- annihilation lines in our approach.

In Sect. 4.3, AQN emission was studied for different instruments in their corresponding bands. While this paper focused on AQN properties in galaxy cluster environments rather than properties within the Milky Way, a direct comparison to McNeil Forbes & Zhitnitsky (2008b) cannot be conducted. Nevertheless, it is striking to observe the abundance of AQN signatures in the WMAP maps depicted in Figs. 11 and 12. AQN signatures in the WMAP band-pass can be identified even after including synchrotron background emission. While it is important to interpret the results with caution given the significant differences between cluster and galactic properties, the possibility of an excess in microwaves in the galactic core due to AQN emission may not be entirely unfounded.

6 Conclusions

We studied the glow of axion quark nugget dark matter as signatures in a large sample of 161 simulated galaxy clusters using the constrained cosmological simulation of the local universe, SLOW. We provide lower-limit predictions on electromagnetic counterparts of AQNs in the environment of galaxy clusters by inferring their thermal and non-thermal emission spectrum originating from axion quark nugget-cluster gas interactions. Dividing our cluster sample into a sub-sample ordered in five mass bins ranging from 0.8 to 31.7 × 10¹⁴ M_⊙ allows us to specify the contribution of the AQN signal to the emission of the ICM as a function of cluster mass. Having 11 cross-identified galaxy clusters within this sample allows us in addition to propose promising galaxy cluster candidates, where the AQN signal might show the strongest contribution to the overall emission.

6.1 Key results

Throughout this analysis, we fixed the parameters of the underlying AQNs model to have M_AQN = 16.7 g. We then computed the effective radius (R_eff) and temperatures (T_AQN) of the AQNs based on the local gas and dark matter properties within the simulations, which allowed us to compute the expected electromagnetic signal of the AQNs as well as the thermal and non-thermal emission of the ICM, where for the later we assumed the same ratio of between CR electrons and thermal energy as observed in the Coma galaxy cluster. Our key findings can be summarized as follows:

1. Even though more massive galaxy clusters typically exhibit higher ICM temperatures, the largest population of high R_eff can be observed for the most massive mass bin in sample 𝒜 as Δv scales with cluster mass as well. This is an interesting feature, since T_AQN directly scales with R_eff and the AQN emission highly depends on T_AQN. Consequently, strong AQN emission can be expected for massive galaxy clusters. The increased AQN emission, however, is always accompanied by high cluster gas emission, as thermal gas emission of the ICM scales with cluster gas temperatures as well. To some degree, relative emission domination is therefore equilibrated for higher cluster masses.

2. A strong AQN emission contribution scales with the mass of a galaxy cluster in the low energy regime. Galaxy clusters with different M_vir can be populated at a fixed ν_T and differences in the absolute value of $\tilde{J}$ are attributed to a combination of T_gas, n_gas, and Δv.

3. In the low-energy frequency range, approximately 40.99% of galaxy clusters show a stronger emission from ν_min to ν_T in the AQN emission over thermal gas emission. After including synchrotron emission, values of the integrated spectral ratio $\tilde{J}$ are distributed in the range of 0.33–53.76%. However, as not all galaxy clusters feature giant radio halos, this opens a window of opportunity when the synchrotron emission of individual galaxy clusters is sufficiently low.

4. Only low-mass galaxy clusters showed an excess of nonthermal AQN emission attributing approximately 9.32% to the 161 galaxy clusters out of both samples, with values of $\tilde{J}$ distributed between 50.71 and 73.65%.

5. A small amount of galaxy clusters exhibit large values of ν_T. This allows for follow-up observations in a wider range of possible surveys to be conducted. In particular with respect to the Fornax and Virgo host environments, which are responsible for extremely high transition frequencies of ν_T = 3842.19 GHz and 382.45 GHz, respectively. To make use of this feature, which is naturally embedded in the thermal emission property of AQNs, improvements in pinpointing regions in galaxy clusters that can only be attributed to an AQN emission excess are required.

6. Shock features are no tracers for regions of strong AQN emission.

7. The spatial distribution of Δv in the SPH maps is not strongly correlated with AQN emission signatures.

8. Fornax and Virgo show in the WMAP and Planck $\tilde{J}$-maps significant AQN signatures. Due to their different dynamical states, regions of strong AQN emission vary in their spatial distribution. Fornax, as the dynamical cluster has strong AQN emission contributions in the infalling substructures; whereas Virgo predominantly shows high $\tilde{J}$-values originating from locations of the ICM and distributed in both central and peripheral regions of the cluster.

9. Even though a general increase of emission in the radial profiles can be attributed to AQN signatures, it is hard to identify specific regions that only trace AQN signatures, since radial AQN emission profiles share a similar morphology as gas emission profiles. We assumed the model to be true and showed that AQN emission affects the spectrum and is visible in the emission maps. However, the AQN-background emission offset is so small that it is not possible to verify or discard this dark matter model observationally, even though it might be one of the most promising ones.

6.2 Future work

First of all, our study was constrained to regions enclosing a radius of 1.5 r_vir. This provides the first valuable insights on expectable AQN signatures from the innermost regions. However, it is expected that especially in the peripheries of galaxy clusters - more specifically in the warm-hot intergalactic medium – the environmental gas temperature of T_gas ∊ [10⁵,10⁷] K (Davé et al. 2001) may permit a better detection of AQN signatures. Lower gas temperatures and number densities would yield a lower thermal gas emission. On the other hand, a lower number density of gas and AQNs would play a decrementing role in peripheral AQN excess detections. Simultaneously, the background emission polluting the AQN signature would decrease as well and both of the emissions will probably fall below the detection limit. However, it is still important to study the outskirts of galaxy clusters for tracers of AQNs, for example in simulations, to confirm or discard possible non-detection predictions.

Furthermore, we did not use a mass and size distribution for AQNs throughout this study and set every AQN to a size of R_AQN = 2.25 × 10^-5 cm and a mass of M_AQN = 16.7 g. Of course, this is a simplification that is not expected in the real world (see for example Ge et al. 2019). It is reasonable to assume the utilized R_AQN and M_AQN as a rough reference value in an analysis conducted for a single snapshot at z = 0 from the cosmological simulation. However, it would also be sensible to implement an on-the-fly model for AQNs directly evolving in a cosmological simulation. This would imply different halo distributions due to dynamical friction of differently massive AQN particles, leading to segregation of more massive AQNs in the central regions of galaxy clusters and probably influencing the radial dark matter density profiles, too. Consequently, not only density profiles would be altered by a constant influence of AQNs on the cluster environment throughout time. The annihilation and heating processes of AQNs would result in feedback mechanisms that act directly on the ICM. Therefore, temperature profiles in galaxy clusters would be influenced too, if the AQN model in the simulation would be treated dynamically. An additional consequence would be that the heated ICM will become visible as tail-like traces behind substructures that fall towards the center of the galaxy cluster. Tails of galaxies are already detected in radio observations (Vallee 1988; Sun et al. 2005; Terni de Gregory et al. 2017; Chen et al. 2020; de Vos et al. 2021; Hu et al. 2021; Müller et al. 2021; Pal & Kumari 2023) and are mostly attributed to AGN feedback mechanisms and ram-pressure stripping. Constant gas interactions with an infalling subhalo would cause similar effects that are recommended to be verified in a dynamical AQN treatment in the simulation. A suitable AQNs mass distribution would therefore be reasonable if on-the-fly AQN models are implemented in the simulations.

Throughout this paper, the similarity between AQNs and neutron stars was pointed out multiple times, as both objects are expected to host cores in a CS state - the same connection can be drawn to magnetars. It is assumed that strong magnetic fields of magnetized neutron stars can be established if the core is ferromagnetic with nuclear density. A solid foundation of the model for magnetized quark nuggets (MQNs) and predictions were developed over the years in for example VanDevender et al. (2020a,b, 2021b); Sloan et al. (2021); VanDevender et al. (2021a, 2024). If MQNs exist, interesting implications on features for galaxy clusters and observations would follow from their properties. First, MQNs would propose an additional source for magnetic fields in galaxy clusters as multiple MQNs could possibly form (when aligned) a large magnetic field and MQNs as additional magnetic field sources are not far-fetched. Second, if MQNs propagated through the ICM, electrons will be accelerated by their strong magnetic fields due to an additional Lorentz force component. Since MQNs could be aligned to form large-scale magnetic fields, they could further act as a Faraday rotating medium of polarised synchrotron emission from the ICM. Since no studies were found addressing these large-scale implications of MQNs, these features are purely speculative, and it is important to verify if these implications hold true or not by analyzing simulations and observations.

The physical effect of magnetization in quark composites serving as dark matter is that their effective radius would be influenced by the surrounding magnetic field. Furthermore, studies show that the strong magnetic field in MQNs influences the interaction with the surrounding plasma as each MQN would develop a magnetopause that causes them to lose kinetic energy while moving through the plasma (VanDevender et al. 2017, 2020a). This would have the effect that dark matter halos of substructures in galaxy clusters would experience a frictional component. By estimating an on-the-fly stopping power effect for a dark matter halo, it would be possible to observe a change in the radial dark matter profiles. On the other hand, the aforementioned studies focused on magnetized quark nuggets with a different composition in the quark core and slightly different properties than axion quark nuggets. Santillán & Sempé (2020) studied magnetic effects specifically for quark nuggets that contain an axion domain wall in their structure and suggested that once an axion domain wall is present, a ferromagnetic state cannot be established. They concluded, however, that this state might be possible for ordinary matter nuggets.

In conclusion, the general theory of AQNs and similar derivations of this model yield physically interesting outcomes. Remarkably, the AQN model was initially introduced to resolve the fine-tuning problem of the cosmological relation Ω_dark ~ Ω_visible and the matter-antimatter asymmetry. Long-lived cosmo-logical problems are naturally resolved by this model. Certain observable features suggest predictions that can be validated by numerical methods, as we show throughout this study. This work serves as a first foundation to show that observable AQN features can be studied in cluster environments using cosmological simulations. We show that it might be possible to detect AQN signatures, especially in galaxy clusters with low synchrotron emission from the ICM. Our results also show that to disentangle pure AQN contributions from the radiation background in galaxy clusters, we need improved strategies and more sophisticated models to predict the non-thermal emission of galaxy clusters.

Acknowledgements

We thank the Center for Advanced Studies (CAS) of LMU Munich for hosting the collaborators of the LOCALIZATION project for a week-long workshop. This work was supported by the grant agreements ANR-21-CE31-0019 / 490702358 from the French Agence Nationale de la Recherche / DFG for the LOCALIZATION project. K.D. acknowledges support by the Excellence Cluster ORIGINS which is funded by the Deutsche Forschungsge-meinschaft (DFG, German Research Foundation) under Germany’s Excellence Strategy – EXC-2094 – 390783311 and funding for the COMPLEX project from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program grant agreement ERC-2019-AdG 882679. X.L., L.V.W., A.Z., and F.M. acknowledge the support from NSERC. The calculations for the hydro-dynamical simulations were carried out at the Leibniz Supercomputer Center (LRZ) under the project pn68na.

Appendix A Supplementary figures

In Fig. A.1, we show the high energy part of the spectra for 11 cross-identified clusters from our sample. The different amounts of the contribution of the AQN emission to the X-ray signal can be clearly seen for the individual clusters. For Fornax, the AQN signal in our simulation even outshines the thermal emission in the hard X-ray regime.

In Fig. A.2, we show the radial profiles of $\tilde{J}\left(r\right)$ (see Eq. 38) for seven of the cross-identified clusters. They are the same ones that showed a thermal AQN emission excess in the spectral features. Note that for some of the systems (especially Virgo, Fornax, and A2256) the contribution within the core of the clusters by the AQNs emission is reaching a 10% level.

Fig. A.1 Spectral plot for non-thermal AQN emission and thermal Bremsstrahlung emission of the hot ICM incorporating the spectra from sample 𝓑. Fornax, being the only galaxy cluster, shows a small region where AQN emission dominates gas emission.

Fig. A.2 Radial profiles of $\tilde{J}\left(r\right)$ (see Eq. 38) for cross-correlated galaxy clusters that showed a thermal AQN emission excess in the spectral features. Galaxy clusters with high $\tilde{J}$ in Fig. 7 exhibit higher $\tilde{J}\left(r\right)$ in central regions of the galaxy cluster.

References

All Tables

All Figures

Fig. 1 Schematic sketch for the antimatter-AQN’s internal structure (not true to scale) with an indication of typical sizes (left panel) and the processes, which lead to the different radiation signatures (middle and right panel). The antimatter-AQN consists of a core of anti-quarks surrounded by positrons. When colliding with the ambient gas, the antimatter-AQN heats up and radiation occurs, which depends on the speed at which the collision occurs, the gas temperature and density as well as the ionization state of the gas. The right panel shows the two key emissions, namely thermal and non-thermal emissions. From an incoming proton (green, small circle), 2 GeV of energy will be available after collision with the AQN. Through annihilation, energetic quarks and gluons are produced and stream either deep within the nugget (around 90% of the 2 GeV) or toward the nugget surface (around 10% of the 2 GeV). Quarks and gluons, which remain inside, will contribute to the thermalization of the nugget. Heat will be transferred from the core to the electrosphere which cools down via thermal Bremsstrahlung emission. Quarks and gluons, which proceed to the surface without thermalization, can transfer energy to positrons and produce non-thermal emission.

In the text

	Fig. 2 Spectral range for thermal emission from AQNs. The magnitude and width of a thermal spectrum increase for increasing TAQN and can reach from radio to X-ray frequencies for extraordinary hot AQNs.
In the text

	Fig. 3 Non-thermal AQN emission spectrum for different TAQN serving as a validation to reproduce Fig. 3 from Majidi et al. (2024).
In the text

	Fig. 4 Histograms of physical properties relevant to the AQN emission. Median values were taken for each mass bin, i.e. for 30 galaxy clusters, respectively.
In the text

	Fig. 5 Radial profiles of Δv, Reff, and TAQN mapped onto AQN tracer particles. The median radial profile was taken for each mass bin (i.e. 30 galaxy clusters each) and compared to other mass bins in a region of r ∈ [0,1.5] rvir.
In the text

	Fig. 6 Median cluster spectra with lines for each mass bin out of sample 𝒜 over a frequency range of ν ∈ [10−5,108]/(2π)eV. The solid, dash-dot-dotted, dashed, and dash-dotted lines represent AQN thermal emission, non-thermal emission, ICM Bremsstrahlung, and synchrotron emission, respectively.
In the text

Fig. 7 Transition frequency-integrated ratio-relation, color-coded by Mvir with frequency ranges of positive AQN contributions. The low-energy frequency range shows in which the thermal AQN emission dominates. Low cluster masses tend to reside on a tangent, slightly offset towards higher values. The high-energy frequency range shows in which the nonthermal AQN emission dominates. Only AQN emission, originating from the lowest Mvir shows an excess over the thermal gas emission in the high-energy regime. We removed one cluster from the plot where for numerical reasons we could not compute νT and the associated $\tilde{J}$ reliably.

In the text

Fig. 8 Ratio map expressed in $\tilde{J}$, shown in the left column. The thermal gas, AQN, and synchrotron emission of Coma integrated over 1.40 GHz to 4.85 GHz, projected 3 rvir along the line of sight are given in the right column (from top to bottom). In this frequency window, the synchrotron emission of Coma is well constrained through radio emission (Thierbach et al. 2003). Due to the strong radio emission of Coma, here the contribution to the diffuse emission by AQNs is in the sub-percent level.

In the text

	Fig. 9 Cross-identified Fornax cluster represented in different relevant properties.
In the text

	Fig. 10 Cross-identified Virgo cluster represented in different relevant properties.
In the text

	Fig. 11 Fornax cluster AQN, gas, and synchrotron emission generated from WMAP, Planck, Euclid, and XRISM energy ranges. Synchrotron emission is predominantly abundant in WMAP and Planck bands. Thermal gas emission shows the most significant contribution to the background emission in the Euclid and XRISM energy bands.
In the text

	Fig. 12 Virgo cluster AQN, gas, and synchrotron emission generated from WMAP, Planck, Euclid, and XRISM energy ranges. Synchrotron emission is predominantly abundant in WMAP and Planck bands. Thermal gas emission shows the most significant contribution to the background emission in Euclid and XRISM energy bands.
In the text

	Fig. A.1 Spectral plot for non-thermal AQN emission and thermal Bremsstrahlung emission of the hot ICM incorporating the spectra from sample 𝓑. Fornax, being the only galaxy cluster, shows a small region where AQN emission dominates gas emission.
In the text

	Fig. A.2 Radial profiles of $\tilde{J}\left(r\right)$ (see Eq. 38) for cross-correlated galaxy clusters that showed a thermal AQN emission excess in the spectral features. Galaxy clusters with high $\tilde{J}$ in Fig. 7 exhibit higher $\tilde{J}\left(r\right)$ in central regions of the galaxy cluster.
In the text