# Influence of Cosmic Voids on the propagation of TeV Gamma Rays and the
Puzzle of GRB 221009A

Hassan Abdalla1,2,3, Soebur Razzaque4,5,6, Markus Böttcher3, Justin Finke7 and
Alberto Domínguez1
1 IPARCOS and Department of EMFTEL, Universidad Complutense de Madrid, E-28040
Madrid, Spain
2 Department of Astronomy and Meteorology, Omdurman Islamic University,
Omdurman 382, Sudan
3 Centre for Space Research, North-West University, Potchefstroom 2520, South
Africa
4 Centre for Astro-Particle Physics (CAPP) and Department of Physics,
University of Johannesburg, PO Box 524, Auckland Park 2006, South Africa
5 Department of Physics, The George Washington University, Washington, DC
20052, USA
6 National Institute for Theoretical and Computational Sciences (NITheCS),
South Africa
7 U.S. Naval Research Laboratory, Code 7653, 4555 Overlook Ave. SW,
Washington, DC, 20375-5352, USA E-mail<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

The recent detection of gamma-ray burst GRB 221009A has attracted attention
due to its record brightness and first-ever detection of $\gtrsim 10$ TeV
$\gamma$ rays from a GRB. Despite being the second-nearest GRB ever detected,
at a redshift of $z=0.151$, the distance is large enough for severe
attenuation of $\gamma$-ray flux at these energies due to $\gamma\gamma\to
e^{\pm}$ pair production with the extragalactic background light (EBL). Here,
we investigate whether the presence of cosmic voids along the line of sight
can significantly impact the detectability of very-high energy (VHE, $>$ 100
GeV) gamma rays from distant sources. Notably, we find that the gamma-gamma
opacity for VHE gamma rays can be reduced by approximately 10% and up to 30%
at around 13 TeV, the highest-energy photon detected from GRB 221009A, for
intervening cosmic voids along the line-of-sight with a combined radius of 110
Mpc, typically found from voids catalogs, and 250 Mpc, respectively. This
reduction is substantially higher for TeV photons compared to GeV photons,
attributable to the broader target photon spectrum that TeV photons interact
with. This finding implies that VHE photons are more susceptible to variations
in the EBL spectrum, especially in regions dominated by cosmic voids. Our
study sheds light on the detection of $\gtrsim 10$ TeV gamma rays from GRB
221009A in particular, and on the detection of extragalactic VHE sources in
general.

###### keywords:

Radiation mechanisms: non-thermal — Galaxies: active, GRB 221009A — Galaxies:
jets — Cosmology: Extragalactic Background Light, cosmic void

††pubyear: 2024††pagerange: Influence of Cosmic Voids on the propagation of
TeV Gamma Rays and the Puzzle of GRB 221009A–Influence of Cosmic Voids on the
propagation of TeV Gamma Rays and the Puzzle of GRB 221009A

##  1 Introduction

Gamma-rays from sources located at a cosmological distance and with energies
greater than the threshold energy for electron-positron ($e^{\pm}$) pair
production can be annihilated due to $\gamma\gamma$ absorption by low-energy
extragalactic background photons. These intergalactic $\gamma\gamma$
absorption signatures have attracted great interest in astrophysics and
cosmology due to their potential to indirectly measure the Extragalactic
Background Light (EBL, e.g., Finke & Razzaque, 2009; Sánchez-Conde et al.,
2009; Domínguez et al., 2011a, b; Sinha et al., 2014; Abdalla et al., 2017,
2021), and thereby probe the cosmic star-formation history (e.g., Dwek &
Krennrich, 2013; Abdollahi et al., 2018). Extensive research has been
conducted to study the $\gamma\gamma$ absorption imprints on the gamma-ray
fluxes from distant sources (e.g., Aharonian et al., 2006; Razzaque et al.,
2009; Kneiske & Dole, 2010; Finke et al., 2010; Gilmore et al., 2012;
Domínguez & Ajello, 2015; Biteau & Williams, 2015; Tavecchio, F. & Bonnoli,
G., 2016; Dzhatdoev et al., 2017; Martínez-Huerta et al., 2018; Korochkin et
al., 2019; Dzhatdoev et al., 2019; Lang et al., 2019; Long et al., 2020;
Troitsky, 2021; Li & Ma, 2021; Terzić et al., 2021; Finke et al., 2022; Azzam
& Hasan, 2023; Pérez de los Heros & Terzić, 2023; Qin et al., 2023; Alves
Batista et al., 2023). Observations have revealed that the spectra of some
distant blazars at very high energies (VHE; E $>$ 100 GeV) appear to be
unexpectedly harder after correcting for $\gamma\gamma$ absorption by the EBL
(e.g., Furniss et al., 2013, 2015; Abdalla & Böttcher, 2017). One potential
explanation for the observed spectral hardening could be a lower density of
the EBL than what is currently estimated by the existing models. The EBL might
exhibit an inhomogeneous distribution, with the presence of cosmic voids along
the line of sight to distant blazars. These voids could feasibly lessen the
attenuation of gamma rays emitted by these blazars, thereby rendering them
detectable at greater redshifts (e.g., Furniss et al., 2013; Abdalla &
Böttcher, 2017; Kudoda & Faltenbacher, 2017, 2018).

The recent detection of VHE $\gamma$ rays from the gamma-ray burst GRB 221009A
at a redshift of $z=0.151$ (Malesani et al., 2023; Williams et al., 2023;
Frederiks et al., 2023; Lesage et al., 2023; Burns et al., 2023; Klinger et
al., 2023; Wang et al., 2023; Mbarubucyeye et al., 2023; Zhao et al., 2023;
Aharonian et al., 2023; Xia et al., 2024; Zheng et al., 2024; Zhang et al.,
2024a; Khangulyan et al., 2024; Banerjee et al., 2024; Ror et al., 2024;
Castro-Tirado et al., 2024) by the Large High Altitude Air Shower Observatory
(LHAASO, Huang et al., 2022; Cao et al., 2023b) has increased focus on EBL
models and fundamental physics, since the $\gamma\gamma$ absorption effect is
expected to reduce the source flux by a factor as large as $\sim 2.2\times
10^{4}$ at 18 TeV (Finke & Razzaque, 2023), the energy of the highest-energy
photon detected by LHAASO. In its recent publication Cao et al. (2023a), the
LHAASO Collaboration conservatively estimated the highest-energy photon to be
a 13 TeV or 18 TeV photon, depending on whether the intrinsic spectrum is a
power-law function with an exponential cutoff or a log-parabola, respectively,
after considering the systematic uncertainties in the photon-energy
determination. In order to rigorously evaluate various EBL models, the study
detailed in Cao et al. (2023a) adopts a methodological approach. This involves
initially characterizing the intrinsic spectrum of Gamma-Ray Bursts (GRBs)
utilizing a log-parabolic (LP) function, $dN/dE=J_{0}(E/1\rm
TeV)^{-a-b\log(E/1\rm TeV)}$, or a power-law function with exponential cutoff
(PLEC), $dN/dE=J_{0}(E/1\rm TeV)^{-\alpha}e^{-E/E_{\rm cut}}$, and
subsequently, obtaining the parameters of the functions from fitting LHAASO
data. This systematic procedure facilitates a comprehensive and quantitative
analysis of the different EBL models under consideration. Precise parameters
were obtained by Cao et al. (2023a) from the fitting process for the time
intervals 230-300 s and 300-900 s after $T_{90}$.

Although there is a non-zero probability that the highest-energy photon
originated from a Galactic HAWC source (e.g., Fraija et al., 2022).
Nevertheless, explaining even an $\gtrsim 10$ TeV photon from this source is
difficult. One possibility is via Ultra-high energy cosmic rays (UHECRs)
emitted from GRB 221009A which produce $\gamma$ rays by interacting with CMB
and EBL photons along our line-of-sight (e.g., Das & Razzaque, 2023; Alves
Batista, 2022; Huang et al., 2023b; He et al., 2024; Kalashev et al., 2024).
Exotic physics involving Lorentz invariance violation (LIV, e.g., Dzhappuev et
al., 2022; Baktash et al., 2022; Li & Ma, 2023; Finke & Razzaque, 2023; Zheng
et al., 2023; Abdalla et al., 2024; Yang et al., 2024; The LHAASO
collaboration, 2024) or mixing of $\gamma$ rays with Axion-like particles
(ALPs, e.g., Baktash et al., 2022; Troitsky, 2022; Galanti et al., 2023; Wang
& Ma, 2023; Gonzalez et al., 2023; Nakagawa et al., 2023; Gao et al., 2024;
Troitsky, 2024; Rojas et al., 2024) or radiative decay during the propagation
of a hypothetical heavy neutrino produced by GRB 221009A (e.g., Cheung, 2022;
Smirnov & Trautner, 2023; Huang et al., 2023a; Guo et al., 2023; Brdar & Li,
2023; Aiello et al., 2024; Dhuria, 2024; ShivaSankar K. et al., 2023) have
also been proposed to evade $e^{\pm}$ pair-production with EBL photons.

In this work, we seek to explain the detection of tens of TeV $\gamma$ rays
from GRB 221009A as a result of a reduced EBL intensity, and hence weaker
$\gamma\gamma$ absorption, due to presence of a cosmic void, regions devoid of
stars and dust, along the line-of-sight. In previous works Abdalla & Böttcher
(2017, 2018a, 2018b), we investigated the effect of such a cosmic void along
the line-of-sight to a distant $\gamma$-ray source on the resulting angle-
averaged EBL energy density by considering the inhomogeneity of the direct
starlight component of the EBL model by Razzaque et al. (2009) in the optical-
ultraviolet. In the present work, we also consider the inhomogeneity of the
dust component of the EBL in the infrared-optical and accordingly modify the
model by Finke et al. (2010). We also considered the potential existence of
cosmic voids, incorporating Domínguez et al. (2011a), Saldana-Lopez et al.
(2021) and Finke et al. (2022) EBL models. Specifically, we examined the
effect of a reduced $\gamma$-ray opacity due to the presence of such voids on
the intrinsic spectra of GRB 221009A. We organize the paper in the following
way. In Section 2, we summarize observations of voids from large-scale surveys
and perform a Monte Carlo study to estimate the potential amount of void that
could exist between us and a source located at the same redshift as GRB
221009A. In Section 3, we describe the model setup and the method used to
evaluate the EBL characteristics and the resulting $\gamma\gamma$ opacity. The
results and conclusions are presented in Sections 4 and 5, respectively.
Throughout this paper the following cosmological parameters are assumed:
$H_{0}=70$ km s-1 Mpc-1, $\Omega_{m}=0.3$, $\Omega_{\Lambda}=0.7$.

##  2 The Presence of Cosmic Voids in Large-Scale Cosmic Surveys

###  2.1 Observational considerations

In recent years, there has been an intensive study of the large-scale
structure and dynamics of galaxy distribution. The role of cosmic voids as
dynamically distinct elements in the universe’s large-scale structure has been
firmly established through both theoretical and observational frameworks
(e.g., Hoffman & Shaham, 1982; Abazajian et al., 2009; Pan et al., 2012;
Tavasoli et al., 2013; Sutter et al., 2015; Wojtak et al., 2016; Douglass et
al., 2023). The most extensive void catalog currently available is designed to
enable precise cosmological experiments using voids and to facilitate more
accurate testing of galaxy formation theories is SDSS DR7 (e.g., Sutter et
al., 2015; Douglass et al., 2023). Characterizing cosmic voids provides
distinct insights into the large-scale distribution of galaxies, their
evolution, and the prevailing cosmological paradigm (Micheletti et al., 2014).
Traditionally, the quest to identify voids within galactic distributions has
employed one of two methodologies: firstly, by Identifying empty spherical
regions amidst the galactic expanse, or secondly, through the implementation
of a watershed algorithm designed to interconnect regions of low density. The
preeminent catalogs of cosmic voids, utilized extensively in contemporary
astrophysical research, have been derived primarily from the seventh data
release of the Sloan Digital Sky Survey SDSS DR7 (see, Abazajian et al., 2009;
Sutter et al., 2015). Recently, Douglass et al. (2023) developed a series of
cosmic void catalogs, utilizing data from the SDSS DR7. Their approach
incorporated three sophisticated void-finding algorithms: VoidFinder and two
algorithms based on the ZOBOV framework (Neyrinck, 2008), namely VIDE and
REVOLVER. This comprehensive analysis led to the identification of over 500
cosmic voids, each exhibiting radii in excess of 10 Mpc/h. The research delved
into various characteristics of these voids, including their radial density
profiles and the spatial distribution of galaxies within them. Furthermore,
this study’s results were cross-referenced with simulated data from the
Horizon Run 4 N-body simulation’s mock galaxy catalogs, demonstrating a
remarkable congruence that substantiated the efficacy of their methodologies
and the validity of their findings.

Figure 1:  This illustration, based on the updated VoidFinder SDSS DR7 from
Douglass et al. (2023), shows the distribution of 2,303 cosmic voids with
effective radii between 14.53 to 43.91 Mpc. The colors on the map reflect the
voids’ relative comoving effective radii, weighted by radius and plotted using
HEALPix grid cells to ensure consistent solid angle coverage across the sky.
This method accurately represents spatial distribution by assigning weights to
corresponding pixels on the HEALPix map. This is done by using a HEALPix map
with mollview for a Mollweide projection in Celestial Equatorial coordinates.
Figure 2: The density of voids within a specified solid angle on the sky is
quantified by calculating the total comoving effective radii within each grid
cell shown in Figure 1. Each cell covers a 10-degree square in both Right
Ascension and Declination. These cells are weighted based on the cumulative
comoving effective radius encompassed within. The color bars represent the
relative amount of void content in each cell compared to others, visually
illustrating variations in void density across the sky. This is done by using
a HEALPix map with mollview for a Mollweide projection in Celestial Equatorial
coordinates.

To provide a clear visual representation, we have plotted the locations of
voids of different sizes in the sky in Figure 1 that demonstrates the quantity
of voids within each grid cell. These cells specifically fall within the
ranges of 140∘ to 240∘ in Right Ascension (RA) and 0∘ to 60∘ in Declination
(Dec). Each cell depicted in Figure 1, encompasses the cumulative voids (total
effective radius) up to a redshift of $z=0.114$. This value represents the
maximum redshift encompassed by the VoidFinder updated Void Catalogs from
Douglass et al. (2023). It is important to note that the sky region
encompassing the extraordinarily luminous GRB 221009A, positioned at RA =
288.264° and Dec = 19.768°, is not included in these catalogs. Consequently,
accurately estimating the extent of voids along our line of sight to GRB
221009A is currently unfeasible. To visualize the variations of cumulative
void sizes in different locations of the sky, we present Figure 1, which
illustrates the distribution of relative void sizes within a region covered by
the updated Void Catalogs SDSS DR7 (Douglass et al., 2023) using the
VoidFinder algorithm. In Figure 2, we computed the accumulated comoving
effective radii contained in the volume of the solid angle subtended by each
cell within the celestial coordinates of 140 degrees to 240 degrees in Right
Ascension and 0 degrees to 60 degrees in Declination. A notable finding is
that there are regions significantly under-dense, with some cells containing
around three times more voids than others. Our choice of a 10-degree square
grid is arbitrary. At smaller scales, such as within 5-degree square cells,
variations become more pronounced. We found that the minimum accumulated
effective radii within a grid cell can be as low as zero, while the maximum,
for the grid cell containing the most voids, can reach up to approximately 250
Mpc. The goal of using a 10-degree grid in Figure 2 is to demonstrate that
significant variations in void quantities can be found even at such large
scales. Therefore, it is conceivable that GRB 221009A might have been situated
in a direction where an under-dense region exists between us and the source.
Such a condition could result in reduced EBL energy density, potentially
facilitating the arrival of TeV photons. In Section 3 we present detailed
calculation by assuming an arbitrary possible void sizes, using the Finke et
al. (2010) EBL model. We characterize the voids in the rest of this Section.

###  2.2 Monte Carlo study on the average radii of cosmic voids

The objective of Figures 1 and 2 is to provide a general understanding based
on observations that the distribution of cosmic voids varies; certain areas of
the sky possess a higher concentration of voids in comparison to others. Here,
we proceed to obtain a more precise quantitative estimate of the cumulative
void radius within a spherical region in the direction towards a source, from
an observational perspective. This involves assessing the likelihood of
encountering such amount of voids in a specific area of the sky where a
particular source might be situated behind these voids. To estimate the
potential cumulative void radius of voids between us and a source at redshift
$0.151$, such as GRB 221009A, we initially extrapolated observed data beyond
$z=0.114$, the limit covered by the most recent Void Catalogs from the
VoidFinder SDSS DR7, extending up to $z=0.151$. This was done by assuming that
the distribution of voids follows the same pattern up to this redshift. The
properties of voids are not expected to change significantly over such low
redshifts (see, Sheth & van de Weygaert, 2004; van de Weygaert & Platen, 2011;
Paz et al., 2013). The primary goal of this work is to determine if voids
exist along the line-of-sight to a very high energy (VHE) gamma-ray source,
potentially leading to a reduction in the EBL density and, consequently, the
optical depth for TeV photons. In Section 3, for simplicity, we assumed the
possibility of multiple voids with typical radii distributed in a large
spherical region between the Earth and the source. This section explores the
feasibility of this assumption from an observational consideration. However,
the voids catalogs of SDSS do not contain the portion of sky where GRB221009A
is located; therefore, we cannot obtain an exact realistic estimation of the
total effective radius of voids towards the source. Therefore, we estimate the
likelihood of such an amount of voids being present in a large spherical
region between us and the source at $z=0.151$. To achieve this, we conducted a
comprehensive Monte Carlo simulation using the following approach:

Figure 3: Results of Monte Carlo simulations of the cumulative effective radii
of voids within spherical regions, situated between us and various random
coordinates. These coordinates fall within the Right Ascension (RA) range of
140 to 230 degrees and the Declination (Dec) range of 10 to 60 degrees. The
Gaussian distribution, depicted by the red line, has a mean of 113 Mpc and a
standard deviation of 50 Mpc.

We randomly selected points within specified celestial coordinates, ranging
from 140 to 230 degrees in RA and 10 to 60 degrees in Dec at redshift z =
0.15. For each iteration, we draw a spherical region with its diameter
extending from the observer to one of these points. Within this spherical
comoving volume, we calculated the cumulative comoving effective radius
($\Sigma R_{i}$), which is found in the VoidFinder SDSS DR7 void catalog and
our extrapolated data beyond redshift z = 0.114. This procedure was repeated
across a significant sample size of $10^{5}$ iterations, creating an extensive
dataset. Subsequently, we calculated the mean and standard deviation of these
cumulative effective radii to assess the distribution and variability across
various directions. Our analysis indicated an average $\Sigma R_{i}=113$ Mpc,
with a standard deviation of 50 Mpc. Furthermore, we evaluated the
probabilities of $\Sigma R_{i}\geq$ 150 Mpc, 200 Mpc, and 250 Mpc, which were
found to be $23.6~{}\%$, $4.6~{}\%$, and $0.6~{}\%$, respectively. Figure 3
visually represents the possible cumulative effective radii of voids within
spherical regions across specific areas of the sky. It features a histogram
that illustrates the probability of encountering varying amounts of voids in
these regions, The red line depicts the probability density function (PDF)
smoothly summarizing the distribution trend. This detailed visualization aids
in understanding the distribution, variation, and overall density of cosmic
voids. Additionally, we conducted a Monte Carlo simulation using the same
region of the sky to calculate the fraction of the line of sight that
traverses a void. Our results showed that the average cumulative effective
radius is 12.97 Mpc, with a standard deviation of 19.33 Mpc. Therefore, on
average, approximately 2% of the distance between us and the source is
traversed through a void.

##  3 EBL in the Presence of Cosmic Voids

Our calculations of the inhomogeneous EBL are based on a modified version of
the formalism presented in Finke et al. (2010), considering the direct
starlight and the re-processing of starlight by dust. For the purpose of a
generic study of the effects of cosmic voids along the line of sight to the
$\gamma$-ray source (i.e., GRB 221009A), we assume that a spherical cosmic
void is located with its center at redshift $z_{v}$ and radius $R$ between the
observer and the $\gamma$-ray source at redshift $z_{s}$. The geometry is
illustrated in Figure 4. We calculate the angle- and photon-energy-dependent
EBL energy density at each point between the observer and the source by using
co-moving coordinates, converting redshifts $z$ to distances $l(z)$. However,
the voids are vast areas that are almost empty of luminous matter (see, Foster
& Nelson, 2009; Ricciardelli et al., 2014). The matter density at the center
of the void can be less by up to 80% compared to its surrounding areas (see,
e.g., Ricciardelli et al., 2014). Therefore, for simplicity, we assumed that
the contribution to the EBL energy density produced within cosmic voids is
negligible. In this work, we set the total emissivity (including both stellar
and dust components) within these voids to zero.

$\ell_{v}$$R$$R$$\theta~{}~{}$$\ell_{1}$$\ell_{2}$$\tilde{\ell}$$~{}~{}~{}~{}~{}\theta_{max}$$R$$z$$z_{v}$$z_{s}$
Figure 4: Illustration of an underdense region (void) situated between the
observer at redshift $0$ and the source at redshift $z_{s}$, where a
$\gamma$-ray photon is emitted from the source and travelling towards the
observer. The void is assumed to have a radius of $R$, and the redshift at its
center is denoted as $z_{v}$. The distances $\ell_{1,2}$ represent the points
where the direction of travel of the EBL photon intersects the boundaries of
the void. $\ell_{v}$ is the distance from the photon’s position at redshift
$z$ to the center of the void, and $\theta_{\rm max}$ represents the maximum
angle at which the direction of travel of the EBL photon still intersects the
boundary of the void at one point.

From Figure 4 we can find the distances $\ell_{1,2}$ where the EBL photon
travel direction $\Omega=(\theta,\phi)$ crosses the boundaries of the void, in
the case where the $\gamma$-ray photon is not inside the void, as

$\ell_{1,2}(\theta)=\ell_{v}~{}\mu~{}{\pm}~{}\sqrt{R^{2}-\ell_{v}^{2}\sin^{2}\theta},$
(1)

where $\mu=\cos\theta$, $\theta$ being the angle between the line of sight and
the EBL photon travel direction, and $\ell_{v}$ is the distance from the
position of the photon at redshift $z$ to the center of the void. The maximum
angle $\theta_{\rm max}$ at which the EBL photon travel direction still
crosses the boundary at one point, is given by:

$\sin\theta_{\rm max}=\frac{R}{\ell_{v}}.$ (2)

and the corresponding distance to the tangential point, $\tilde{\ell}_{\rm
max}$, is given by

$\tilde{\ell}_{\rm max}=\ell_{v}\,\cos\theta_{\rm max}.$ (3)

To calculate the EBL energy density inside the void, there is exactly one
point where an EBL photon crosses the boundary of the void, at a distance
$\ell_{1}$ from $z$. Note that $z$ represents an arbitrary redshift along the
line of sight between the observer and the source, where it equals zero at the
observer’s location and $z_{s}$ at the source. To account for the EBL
inhomogeneity, we introduce a step function that sets the emissivities inside
the void to zero. For a point at a distance $\tilde{\ell}$ from the location
of $z$ (as a function of $\theta$), this step function can be written as

$f_{\rm void}^{\rm outside}~{}=~{}\bigg{\\{}\begin{tabular}[]{cc}0&~{} ~{} if
~{} $\ell_{1}(\theta)~{}<~{}\tilde{\ell}~{}<~{}\ell_{2}(\theta)$\\\ 1&~{} ~{}
otherwise\end{tabular}$ (4)

for points $z$ along the line of sight that are located outside the void, and

$f_{\rm void}^{\rm inside}~{}=~{}\bigg{\\{}\begin{tabular}[]{cc}0&~{} ~{} if
~{} $\tilde{\ell}~{}<~{}\ell_{1}(\theta)$\\\ 1&~{} ~{} otherwise\end{tabular}\
.$ (5)

We adapt the expression presented in Finke et al. (2010) to modify the
comoving luminosity density (luminosity per unit comoving volume), also termed
as emissivity. This expression allows us to investigate the dependence of the
luminosity density on comoving energy $\epsilon$ at a specific redshift
$z_{e}$:

$\begin{split}\epsilon
j^{stars}(\epsilon,z_{e})=m_{e}c^{2}\epsilon^{2}f_{esc}(\epsilon)\int_{m_{\rm
min}}^{m_{\rm max}}dm\,\xi(m)\\\ \times\int_{{z_{e}}}^{z_{\rm
max}}dz_{1}\left|\frac{dt_{*}}{dz_{1}}\right|\psi(z_{1})\dot{N}(\epsilon,m,t_{*}[z_{e},z_{1}])\
;\end{split}$ (6)

where the initial mass function (IMF) is represented by $\xi(m)$, the comoving
star formation rate (SFR) density (i.e., the SFR per unit comoving volume) is
denoted as $\psi(z_{1})$, $f_{esc}(\epsilon)$ is the fraction of photons that
manage to escape a galaxy without being absorbed by interstellar dust, and
$\dot{N}(\epsilon,m,t_{*}(z_{e},z_{1}))$ is the total number of photons
emitted per unit energy per unit time. The connection between cosmic time and
redshift can be expressed as

$\left|\frac{dt_{*}}{dz_{1}}\right|=\frac{1}{H_{0}(1+z_{1})\sqrt{\Omega_{m}(1+z_{1})^{3}+\Omega_{\Lambda}}}.$
(7)

To incorporate the effects of dust and account for the inhomogeneity in the
comoving luminosity density, we have made modifications to the expression
introduced in Finke et al. (2010),

$\begin{split}\epsilon j^{dust}(\epsilon,z_{e})=\frac{15}{\pi^{4}}\int
d\epsilon\left[\frac{1}{f_{esc}(\epsilon)}-1\right]j^{stars}(\epsilon,z_{e})\\\
\times\sum_{n=1}^{3}\frac{f_{n}}{\Theta^{4}_{n}}\frac{\epsilon^{4}}{\exp{\left(\frac{\epsilon}{\Theta_{n}}\right)-1}},\end{split}$
(8)

where $f_{n}$ represents the fraction of the absorbed emissivity reradiated
within a specific dust component, $\Theta_{n}=k_{B}T_{n}/m_{e}c^{2}$ denotes
the dimensionless temperature of the respective dust component with subscripts
n = 1, 2, and 3 corresponding to the warm dust, hot dust, and PAHs,
respectively (for more details, see, Finke et al., 2010). The differential
photon number density of the EBL at a given redshift $z$ along the travel
direction of the $\gamma$-ray photon consists of the direct contribution from
stars throughout cosmic history and the indirect contribution from the
absorption and re-emission by dust. The energy density $u_{\rm EBL}$ of the
EBL in the comoving frame per unit comoving volume as a function of $\epsilon$
can be written as:

$\epsilon u_{\rm EBL}(\epsilon,z)=\int_{z}^{z_{\rm
max}}dz_{e}\left|\frac{dt}{dz_{e}}\right|\frac{\epsilon^{\prime}j(\epsilon^{\prime,},z_{e})}{(1+z_{e})}$
(9)

The quantity $j(\epsilon^{\prime},\tilde{z})$ is the stellar and dust
emissivity function defined as $j(\epsilon^{\prime},\tilde{z})=j^{\rm
stars}(\epsilon^{\prime},\tilde{z})+j^{\rm dust}(\epsilon^{\prime},\tilde{z})$
given by Eqs. (6) and (8). $dt/d\tilde{z}$ is evaluated using Eq. (7). An
observer at a redshift $z>0$ would perceive a universe that is $(1+z)^{3}$
times smaller than its current size. Photons that we observe today with energy
$\epsilon$ would have been seen by them with a proper energy
$\epsilon_{p}=(1+z)\epsilon$. Taking into account the cosmic voids described
earlier, the energy density expected to be observed can be expressed as:

$\begin{split}\epsilon_{p}u_{{\rm
EBL,}p}({\epsilon_{p},z,\Omega})=(1+z)^{4}\int_{z}^{z_{\rm
max}}d\tilde{z}\left|\frac{dt}{d\tilde{z}}\right|\\\
\frac{\epsilon^{\prime}j(\epsilon^{\prime},\tilde{z},\Omega)}{(1+\tilde{z})}f_{\rm
void}(z,\tilde{z},\Omega),\end{split}$ (10)

where $f_{\rm void}(z,\tilde{z},\Omega)$ is the step function set to be zero
inside the void (Eqs. 4 and 5), and

$\epsilon^{\prime}=\epsilon(1+\tilde{z})=\frac{(1+\tilde{z})}{(1+z)}\epsilon_{p}.$

When including the effects of voids, the EBL energy density exhibits an
angular dependence. Consequently, we calculate the optical depth for a
$\gamma$-ray photon from a source at redshift $z_{s}$ with observed energy
$E$, using the angle-dependent EBL photon energy density from Equation (10) as
(see, e.g., Razzaque et al., 2009):

$\begin{split}\tau_{\gamma\gamma}(E,z_{s})=c~{}~{}{\int_{0}^{z_{s}}}d{z_{1}}\left|\frac{dt}{d{z_{1}}}\right|\int_{0}^{\infty}d\epsilon_{1}\\\
\oint\
d\Omega\frac{\mu_{\epsilon_{1}}(\epsilon_{1},z_{1},\Omega)}{\epsilon_{1}}~{}(1-\mu)\sigma_{\gamma\gamma}(s).\end{split}$
(11)

The $\gamma~{}\gamma$ pair-production cross section $\sigma_{\gamma\gamma}(s)$
can be written as:

$\sigma_{\gamma\gamma}(s)=\frac{1}{2}\pi
r_{e}^{2}(1-\beta_{cm}^{2})(3-\beta_{cm}^{4})\ln\left(\frac{1+\beta_{cm}}{1-\beta_{cm}}\right)-2\beta_{cm}(2-\beta_{cm}^{2})$

where $r_{e}$ is the classical electron radius and
$\beta_{cm}=(1-\frac{1}{s})^{1/2}$ is the electron-positron velocity in the
center-of-momentum (c.m.) frame of the $\gamma\gamma$ interaction,
$s=\frac{s_{0}}{2}(1-\cos\theta)$ is the c.m. frame electron/positron energy
squared and $s_{0}=\frac{\epsilon E}{m_{e}^{2}c^{4}}$. In the case of a
homogeneous and isotropic EBL (with which we will compare our results for the
inhomogeneous EBL case), equation (11) can be simplified using the
dimensionless function $\bar{\varphi}$ defined by Gould & Schreder (1967):

$\bar{\varphi}[s_{0}(\epsilon)]=\int_{1}^{s_{0}(\epsilon)}s\bar{\sigma}(s)ds,$

where $\bar{\sigma}(s)=\frac{2\sigma(s)}{\pi r_{e}^{2}}$ and
$s_{0}(\epsilon)=E(1+z)\epsilon/m_{e}^{2}c^{4}$, so that equation (11) can be
written as

$\begin{split}\tau_{\gamma\gamma}^{\rm hom}(E,z)=~{}c~{}\pi
r_{e}^{2}\left(\frac{m_{e}^{2}c^{4}}{E}\right)^{2}{\int_{0}^{z_{s}}}\frac{d{z_{1}}}{(1+z_{1})^{2}}\left|\frac{dt}{d{z_{1}}}\right|\\\
\times\int_{{\frac{m_{e}^{2}c^{4}}{E(1+z1)}}}^{\infty}d\epsilon_{1}~{}\frac{\mu_{\epsilon_{1}}}{\epsilon_{1}^{3}}~{}\bar{\varphi}[s_{0}(\epsilon)].\end{split}$
(12)

Knowing the optical depth $\tau_{\gamma\gamma}$, we can calculate the
attenuation of the intrinsic photon flux $F_{\nu}^{int}$ as

$F_{\nu}^{\rm obs}=F_{\nu}^{\rm int}e^{-\tau_{\gamma\gamma}(E,z)},$ (13)

where $F_{\nu}^{\rm obs}$ is the observed spectrum. In Fig. 5 we show the
$\gamma\gamma$ opacity for $z=0.151$, the distance to GRB 221009A, and
contributions to the opacity by EBL photons of different wavelengths. It is
clear that there would be significant attenuation of $\gamma$-ray flux at
energies beyond 10 TeV. In Section 4, we will explore to what extent the
presence of voids may decrease the impact on the EBL ${\gamma\gamma}$
absorption in this energy range. The primary objective is to examine whether
the presence of a cosmic void along the line-of-sight to GRB 221009A may
explain the detection of $\gtrsim 10$ TeV photons.

##  4 Results and discussion

In this work we hypothesize that the region along our line of sight to GRB
221009A is dominated by voids of various sizes. We assume that the center of
the cumulative void is located at a redshift of $z_{v}=0.075$, given that the
redshift of GRB 221009A is $z=0.151$. This assumption allows us to simplify
the calculation using a single void where the center of the hypothetical large
sphere, is at the midpoint between us and GRB 221009A. In Abdalla & Böttcher
(2017), we found that the EBL deficit is proportional to the void’s size. The
impact of several small voids, each with radius $R_{i}$, is roughly equivalent
to the effect of a single large void with a radius R equal to the sum of all
the small voids radii $\Sigma R_{i}$. Therefore, the single, arbitrarily void
radius $R$ used in our calculations in Section 3 is effectively equivalent to
the cumulative radius $\Sigma R_{i}$ used in Section2. We consider cases with
void diameters between 15% (corresponding to $R=50$ Mpc) and 80%
(corresponding to $R=250$ Mpc) of the distance to GRB 221009A.

In the top panel of Figure 6, we show the effect of voids of size $50$ Mpc,
$100$ Mpc, $150$ Mpc, $200$ Mpc and $250$ Mpc on the $\gamma\gamma$ optical
depth compared to the homogeneous case (black solid line, Finke et al., 2010).
The bottom panel of Figure 6 shows that the optical depth
$\tau_{\gamma\gamma}$ can be reduced by up to around 10% and 30% at 13 TeV for
the cases of $R=110$ Mpc and $R=250$ Mpc, respectively. Note that the average
combined radius that we obtain from our calculations in Section 3 is R $\sim$
113 Mpc. The deficit due to cosmic voids is much higher for TeV photons than
for GeV photons. This can be explained by the fact that the TeV photons
interact with a broader target photon spectrum, which means that the deficit
for very high energy photons is expected to be due to the deficit in the
entire EBL spectrum, while the deficit for GeV photons is only due to the
deficit in the optical-ultraviolet EBL component. This explanation is clearly
depicted in Figure 5, where we show $\tau_{\gamma\gamma}$ calculated for EBL
photons of different wavelengths. For the same reason, higher energy photons
are more sensitive to the values of the Hubble constant when EBL
$\gamma\gamma$ absorption is used to constrain cosmology (e.g., Domínguez et
al., 2019, 2024).

Figure 5:  Optical depth $\tau_{\gamma\gamma}$, for VHE photons traversing
from a source located at a redshift $z=0.151$ attributable to interactions
with the EBL, evaluated across various EBL energy bands. For additional
details see Figure 10 in Acciari et al. (2019).

Figure 6: Top panel: $\gamma\gamma\to e^{+}e^{-}$ optical depth due to the EBL
as a function of the $\gamma$-ray energy for a source located at redshift
$z=0.151$. The black solid line indicates the optical depth for a homogeneous
EBL distribution. The dashed lines show the opacities in the presence of voids
of various sizes, as indicated by the legend. The vertical red dot-dashed line
in both panels indicates a photon energy of 13 TeV. Bottom panel: Relative
deficit of the optical depth with and without voids as a function of the
$\gamma$-ray energy for the same cases as in the top panel.

###  4.1 Cosmic voids using different EBL models

In Section 3, we formulated a methodology for incorporating the concept of EBL
inhomogeneity, or the potential existence of cosmic voids, within the
framework of the EBL model proposed by Finke et al. (2010). This involved
conducting meticulous calculations, assuming that the cumulative emissivity
within the void regions could be considered negligible. Incorporating the
notion of cosmic voids or under-dense regions within various EBL models is
crucial for conducting a thorough analysis and facilitating comparisons with
the EBL models featured in Cao et al. (2023a). In the top panel of Figure 6 we
show the optical depth calculated for the Finke et al. (2010) EBL model in the
homogeneous case and in the case of cosmic voids of various sizes for a source
at $z=0.151$, the redshift of GRB 221009A. The bottom panel shows the relative
deficit of optical depth without and with voids. It is noteworthy that the
anticipated effects of cosmic voids are expected to remain consistent across
various EBL models.

To incorporate the effect of cosmic voids into other EBL models we define
$\tau_{Finke_{10}}^{hom}$ and $\tau_{Finke_{10}}^{void}$ as the optical depths
calculated under two scenarios: one based on the standard homogeneous model
and the other considering the presence of cosmic voids, both utilizing Finke
et al. (2010) EBL model. Similarly, $\tau_{model_{x}}^{hom}$ and
$\tau_{model_{x}}^{void}$ represent analogous optical depths when employing an
alternative EBL model denoted as “$\text{model}_{x}$” under the same two
scenarios (homogeneous and void-inclusive). Therefore, the optical depth
considering voids for any EBL model can be calculated as,

$\tau_{model_{x}}^{void}=\frac{\tau_{Finke_{10}}^{void}}{\tau_{Finke_{10}}^{hom}}\tau_{model_{x}}^{hom}\,.$
(14)

Based on the hypothesis that the expected optical depth deficit remains
consistent regardless of the chosen EBL model, we can introduce EBL
inhomogeneity into various EBL models using Eq. (14). Figure 7 displays the
optical depths for photons originating from GRB 221009A as calculated using
the Domínguez et al. (2011a), Saldana-Lopez et al. (2021) and Finke et al.
(2022) EBL models. The amount by which the cosmic voids lower the optical
depth depends on the amount of potential voids that exist along the line-of-
sight to the source. Next we study how this change in the optical depth could
potentially impact the observed flux for each EBL model, which amount of voids
could provide an explanation for the LHAASO data, and how much this could
increase the detectability of VHE $\gamma$-rays from these sources.

Figure 7: The optical depth for $\gamma\gamma\to e^{+}e^{-}$ interactions in
the EBL, as a function of gamma-ray energy, is presented for a source located
at redshift $z=0.151$. In the top-left panel, different lines correspond to
various EBL models. We consider both the standard homogeneous case and
introduce various void sizes, as illustrated in the corresponding figures.
Cosmic voids have been incorporated into the EBL models proposed by Domínguez
et al. (2011a), Saldana-Lopez et al. (2021), and Finke et al. (2022). The
vertical red dot-dashed line indicates a photon energy of 13 TeV, which
represents the highest energy photon detected from GRB 221009A by LHAASO, as
reported in Cao et al. (2023a).

###  4.2 GRB 221009A: Observed and intrinsic fluxes

The LHAASO-KM2A detector has observed an unexpected abundance of 140 gamma-ray
events with energies exceeding 3 TeV within a time span of 230 to 900 seconds
following the initial trigger. A log-parabola function (referred to as LP) and
a power-law with an exponential cutoff function (referred to as PLEC) has been
employed to describe the intrinsic energy spectrum of these gamma rays,
incorporating essential corrections to account for the effects of EBL
absorption Cao et al. (2023a). A detailed analysis of GRB 221009A, accounting
for the systematic uncertainties in photon-energy determination, indicates
that the maximum energy observed from this source is, depending on the assumed
spectrum, approximately 18 TeV (LP) or 13 TeV (PLEC), as reported by Cao et
al. (2023a).

In this study, we adopted an arbitrary power-law intrinsic flux and
incorporated EBL attenuation, represented by the equation
$dN/dE=N(E/TeV)^{-\alpha}\times e^{-\tau_{\gamma\gamma}}$ utilizing the EBL
models from Finke et al. (2010), Domínguez et al. (2011a), Saldana-Lopez et
al. (2021) and Finke et al. (2022). Specifically, we employed a fitting
procedure to match the attenuated flux from the adopted EBL model to the
observed data from Cao et al. (2023a), thereby effectively constraining the
values of the normalization and index of the power-law function for the
intrinsic flux. For each EBL model, we calculated the fitting parameters in
the standard homogeneous case first. We then incorporated the cosmic voids of
different sizes for each EBL model and identified the best-fit parameters for
each case. The optimal values of the power-law index $\alpha$ for each
scenario are presented in Table 1. It can be noted that $\alpha$ becomes
smaller (the spectrum becomes harder) with increasing void size, as expected
from a decreasing opacity. However, the value is not significantly different
from the spectral index in the homogeneous case ($\boldsymbol{\alpha}_{hom}$).

Table 1: The spectral fitting analyses of GRB 221009A, conducted over two time intervals $T_{0}+230$ s to 300 s and $T_{0}+300$ s to 900 s, employed a power-law functional form, $dN/dE=N(E/TeV)^{-\alpha}$. The best fitted value for $\alpha$ in the context of homogeneous scenarios, denoted by $\alpha_{hom}$, and in the models accounting for varying cosmic void sizes, represented by $\alpha_{R}$, where $R$ is the void radius. Our analysis revealed that the normalization parameter $N$ for the time interval $T_{0}+230$ s, the values are as follows: $N=(2.9,~{}2.0,~{}2.1,~{}2.1)\times 10^{-6}$ $TeV^{-1}cm^{-2}s^{-1}$ for Finke et al. (2010), Domínguez et al. (2011a), Finke et al. (2022) and Saldana-Lopez et al. (2021) respectively. For the time interval $T_{0}+300$ s, the values are $N=(4.6,~{}3.7,~{}3.3,~{}3.3)\times 10^{-7}$ $TeV^{-1}cm^{-2}s^{-1}$ for Finke et al. (2010), Domínguez et al. (2011a), Finke et al. (2022) and Saldana-Lopez et al. (2021) respectively. Time Interval | EBL Models | $\alpha_{hom}$ | $\alpha_{50}$ | $\alpha_{100}$ | $\alpha_{150}$ | $\alpha_{200}$ | $\alpha_{250}$
---|---|---|---|---|---|---|---
230 - 300 s | Finke et al. (2010) | 2.20$\pm$0.07 | 2.19$\pm$0.07 | 2.18$\pm$0.08 | 2.17$\pm$0.09 | 2.16$\pm$0.11 | 2.15$\pm$0.11
Domínguez et al. (2011a) | 2.32$\pm$0.08 | 2.31$\pm$0.07 | 2.31$\pm$0.06 | 2.30$\pm$0.06 | 2.29$\pm$0.05 | 2.28$\pm$0.05
Saldana-Lopez et al. (2021) | 2.40$\pm$0.11 | 2.36$\pm$0.11 | 2.36$\pm$0.09 | 2.35$\pm$0.09 | 2.34$\pm$0.08 | 2.33$\pm$0.07
Finke et al. (2022) | 2.41$\pm$0.08 | 2.40$\pm$0.07 | 2.39$\pm$0.06 | 2.38$\pm$0.06 | 2.37$\pm$0.06 | 2.36$\pm$0.06
300 - 900 s | Finke et al. (2010) | 2.15$\pm$0.06 | 2.14$\pm$0.07 | 2.13$\pm$0.08 | 2.12$\pm$0.09 | 2.11$\pm$0.09 | 2.09$\pm$0.10
Domínguez et al. (2011a) | 2.30$\pm$0.05 | 2.29$\pm$0.05 | 2.28$\pm$0.05 | 2.27$\pm$0.05 | 2.26$\pm$0.05 | 2.24$\pm$0.05
Saldana-Lopez et al. (2021) | 2.36$\pm$0.06 | 2.35$\pm$0.06 | 2.34$\pm$0.05 | 2.33$\pm$0.05 | 2.32$\pm$0.04 | 2.31$\pm$0.05
Finke et al. (2022) | 2.38$\pm$0.04 | 2.37$\pm$0.04 | 2.36$\pm$0.04 | 2.35$\pm$0.04 | 2.34$\pm$0.05 | 2.33$\pm$0.05

Figure 8: The red data points are LHAASO observations for the interval from
$T_{0}+230$ s to $T_{0}+300$ s, and the intrinsic spectrum obtained through a
power-law (PL) model is represented by the solid blue line, employing the EBL
models described in Finke et al. (2010), Domínguez et al. (2011a), Saldana-
Lopez et al. (2021), and Finke et al. (2022). The anticipated EBL-induced
absorption, as estimated by the indicated EBL models, is represented by the
solid black line. Additionally, we explore the influence of cosmic voids by
incorporating various void sizes, depicted as dashed lines.

Figures 8 and 9 showcase the intrinsic spectra of GRB 221009A, which have been
derived using the power-law model based on the fitted parameters detailed in
Table 1. These figures also present the anticipated observable flux as
estimated by various EBL models. These spectra are shown for two distinct time
intervals: $T_{0}+230$ s to $T_{0}+300$ s (Figure 8) and $T_{0}+300$ s to
$T_{0}+900$ s (Figure 9) as observed by LHAASO (see, Cao et al., 2023a). The
expected EBL-induced absorption, as estimated by the EBL models in the
homogeneous case, is presented as a solid black line. The potential influence
of cosmic voids of various sizes (dashed lines) is investigated. The observed
flux by LHAASO, especially in the $T_{0}+300$ s to $T_{0}+900$ s (Figure 9)
exhibits a hardening of the spectrum above 10 TeV and the fit seems better
with some EBL models in the presence of voids. Although, we notice that no
single void model can explain the whole spectrum.

Figure 9: The blue data points are LHAASO observations for the interval from
$T_{0}+300$ s to $T_{0}+900$ s, and the intrinsic spectrum obtained through a
power-law (PL) model is represented by the solid blue line, employing the EBL
models described in Finke et al. (2010), Domínguez et al. (2011a), Saldana-
Lopez et al. (2021), and Finke et al. (2022). The anticipated EBL-induced
absorption, as estimated by the indicated EBL models, is represented by the
solid black line. Additionally, we explore the influence of cosmic voids by
incorporating various void sizes, depicted as dashed lines.

For voids with accumulated size $R=250$ Mpc, the EBL-attenuated flux at 13 TeV
is can be greater than a factor of ten than the expected one in the case of a
homogeneous EBL. Although illustrative, we note that the LHAASO flux is an
optimistic estimate of the source-intrinsic $\gtrsim 10$ TeV flux from GRB
221009A, but the reduction of the optical depth is independent of the
$\gamma$-ray emission model and intrinsic $\gamma$-ray spectrum. When combined
with our scenario of cosmic voids along the line-of-sight, physical models
such as those presented by Das & Razzaque (2023); Cao et al. (2023b); Zhang et
al. (2024b, 2023) or proton synchrotron emission Zhang et al. (2023) would
require a lower total energy output from GRB 221009A to explain the detection
of >10 TeV photons by LHAASO.

###  4.3 Likelihood Analysis:

Our goal is to determine how well various EBL models, featuring different void
sizes, align with LHAASO’s observational data as presented in Figures 8 and 9.
We calculated the reduced chi-square statistic for each model with different
void sizes, and convert these values to likelihoods for a probabilistic
evaluation using the formula
$\exp\left(-\frac{\chi_{\text{reduced}}^{2}}{2}\right)$.

Figure 10:  The likelihood for the dataset from both intervals, $T_{0}+230$ s to 300 s and $T_{0}+300$ s to 900 s, obtained by fitting a Power-Law (PL) model under various EBL models while assuming different quantities of voids in each scenario, is depicted by the solid line. The maximum likelihood value for each model is marked by vertical dashed lines, with colors corresponding to each specific EBL model represented by the solid lines.  Table 2: The cumulative voids radii $\Sigma R_{i}$ in Mpc necessary for each EBL model get the best fit with LHASSO data. Time Interval | EBL Models |  $\Sigma R_{i}$ in Mpc
---|---|---
230 - 300 s | Finke et al. (2010) | 55±52
Domínguez et al. (2011a) | 215±36
Saldana-Lopez et al. (2021) | 152±51
Finke et al. (2022) | 113±58
300 - 900 s | Finke et al. (2010) | 94±60
Domínguez et al. (2011a) | 250±36
Saldana-Lopez et al. (2021) | 209±46
Finke et al. (2022) | 147±58

Based on our estimates from the Monte Carlo study (see, 3), all these various
quantities of voids, as required by each model, are possible. However, in
certain scenarios, the quantity of voids needed to get the best fit with the
LHASSO data appears less likely. Specifically, for Domínguez et al. (2011a)
and Saldana-Lopez et al. (2021) EBL models, the accumulated radii of voids
$\Sigma R$ necessary to achieve the best fit can exceed 200 Mpc. The
likelihood of such a mount of voids existing for an object at a redshift of
$z=0.151$ is less than 5%.

##  5 Summary and Conclusions

Cosmic voids remain an active area of research, and their precise sizes remain
uncertain from an observational perspective. Determining the minimum size
required for a region to be classified as a void, as well as quantifying the
matter density compared to the average, poses significant challenges.
Estimating the exact reduction in matter density within a spherical region,
assumed to be a void, is particularly complex (e.g., Douglass et al., 2023).
Further in-depth analysis will be crucial in the future to rigorously test our
hypothesis and refine our understanding of cosmic voids and their implications
In this work, we have explored implications of cosmic voids on the $\gamma$
ray propagation along the line-of-sight and through a region of of reduced EBL
opacity, using TeV flux detected by the LHAASO experiment from GRB 221009A. We
have computed the optical depth for $\geq 0.1$ TeV photons propagating to the
earth from $z=0.151$, the redshift of GRB 221009A. We found that the presence
of cosmic voids can reduce the optical depth significantly in the $\gtrsim 1$
TeV range. Correspondingly, assuming the same intrinsic power-law spectrum,
the observable $\gamma$-ray flux from GRB 221009A can be greater than a factor
of ten for the cases with the larger combined radius compared to the
homogeneous EBL models. This can help explain detection of $\gtrsim 10$ TeV
photons from GRB 221009A.
We show that the optical depth may be reduced by about 10% and 30% at 13 TeV,
in the cases of $R=110$ Mpc and $R=250$ Mpc, respectively. Furthermore, our
findings suggest that the fit of all EBL models is more accurate within the
context of possible voids present along our line-of-sight to GRB 221009A. The
quantity of voids needed to achieve optimal fitting falls within the realm of
feasibility, as suggested by our Monte Carlo study employing an updated void
catalog from VoidFinder SDSS DR7. Our results, using EBL de-absorption in the
presence of void shows a harder intrinsic spectrum than for the homogeneous
EBL. Therefore, prospects for future detection of TeV from extragalactic
sources such as GRBs or blazars are enhanced, if those take place in the
directions where cosmic voids are present.

## Acknowledgements

We thank the anonymous referee for a careful reading of the manuscript and
helpful suggestions, which have led to significant improvements in the
manuscript. The work of H.A. is supported through the María Zambrano
researcher fellow funded by the European Union, the AfAS Seed Research Grant,
SA-GAMMA ST/S002952/1 and by the University of Johannesburg (UJ) URC grant.
S.R. was supported by grants from the National Research Foundation (NRF),
South Africa, NITheCS, UJ URC. J.D.F. was support by NASA through contract
S-15633Y and the Office of Naval Research. A.D. is thankful for the support of
the Ramón y Cajal program from the Spanish MINECO, Proyecto
PID2021-126536OA-I00 funded by MCIN / AEI / 10.13039/501100011033.

## Data Availability

The codes used to produce the results featured in this article are available
upon reasonable request to the corresponding author.

## References

  * Abazajian et al. (2009) Abazajian K. N., et al., 2009, Astrophys. J. Suppl., 182, 543
  * Abdalla & Böttcher (2017) Abdalla H., Böttcher M., 2017, The Astrophysical Journal, 835, 237
  * Abdalla & Böttcher (2018a) Abdalla H., Böttcher M., 2018a, The Astrophysical Journal, 865, 159
  * Abdalla & Böttcher (2018b) Abdalla H., Böttcher M., 2018b, PoS, HEASA2017, 028
  * Abdalla et al. (2017) Abdalla H., et al., 2017, Astron. Astrophys., 606, A59
  * Abdalla et al. (2021) Abdalla H., et al., 2021, JCAP, 02, 048
  * Abdalla et al. (2024) Abdalla H., Cotter G., Backes M., Kasai E., Böttcher M., 2024, Classical and Quantum Gravity, 41, 015022
  * Abdollahi et al. (2018) Abdollahi et al., 2018, Science, 362, 1031
  * Acciari et al. (2019) Acciari V. A., et al., 2019, Mon. Not. Roy. Astron. Soc., 486, 4233
  * Aharonian et al. (2006) Aharonian F., et al., 2006, Nature, 440, 1018
  * Aharonian et al. (2023) Aharonian F., et al., 2023, ApJ, 946, L27
  * Aiello et al. (2024) Aiello S., et al., 2024, arXiv e-prints, p. arXiv:2404.05354
  * Alves Batista (2022) Alves Batista R., 2022, arXiv e-prints, p. arXiv:2210.12855
  * Alves Batista et al. (2023) Alves Batista R., et al., 2023, arXiv e-prints, p. arXiv:2312.00409
  * Azzam & Hasan (2023) Azzam W., Hasan A., 2023, Journal of Applied Mathematics and Physics, 11, 2179
  * Baktash et al. (2022) Baktash A., Horns D., Meyer M., 2022, arXiv e-prints, p. arXiv:2210.07172
  * Banerjee et al. (2024) Banerjee B., et al., 2024, arXiv e-prints, p. arXiv:2405.15855
  * Biteau & Williams (2015) Biteau J., Williams D. A., 2015, The Astrophysical Journal, 812, 60
  * Brdar & Li (2023) Brdar V., Li Y.-Y., 2023, Phys. Lett. B, 839, 137763
  * Burns et al. (2023) Burns E., et al., 2023, Astrophys. J. Lett., 946, L31
  * Cao et al. (2023a) Cao Z., et al., 2023a, Sci. Adv., 9, adj2778
  * Cao et al. (2023b) Cao Z., et al., 2023b, Science, 380, adg9328
  * Castro-Tirado et al. (2024) Castro-Tirado A. J., et al., 2024, A&A, 683, A55
  * Cheung (2022) Cheung K., 2022, arXiv e-prints, p. arXiv:2210.14178
  * Das & Razzaque (2023) Das S., Razzaque S., 2023, A&A, 670, L12
  * Dhuria (2024) Dhuria M., 2024, Phys. Rev. D, 109, 063007
  * Domínguez & Ajello (2015) Domínguez A., Ajello M., 2015, ApJ, 813, L34
  * Domínguez et al. (2011a) Domínguez A., et al., 2011a, MNRAS, 410, 2556
  * Domínguez et al. (2011b) Domínguez A., Sánchez-Conde M. A., Prada F., 2011b, J. Cosmology Astropart. Phys., 2011, 020
  * Domínguez et al. (2019) Domínguez A., et al., 2019, ApJ, 885, 137
  * Domínguez et al. (2024) Domínguez A., et al., 2024, MNRAS, 527, 4632
  * Douglass et al. (2023) Douglass K. A., Veyrat D., BenZvi S., 2023, Astrophys. J. Suppl., 265, 7
  * Dwek & Krennrich (2013) Dwek E., Krennrich F., 2013, Astroparticle Physics, 43, 112
  * Dzhappuev et al. (2022) Dzhappuev D. D., et al., 2022, The Astronomer’s Telegram, 15669, 1
  * Dzhatdoev et al. (2017) Dzhatdoev T. A., Khalikov E. V., Kircheva A. P., Lyukshin A. A., 2017, A&A, 603, A59
  * Dzhatdoev et al. (2019) Dzhatdoev T., Khalikov E., Podlesnyi E., Telegina A., 2019, Journal of Physics: Conference Series, 1181, 012049
  * Finke & Razzaque (2009) Finke J. D., Razzaque S., 2009, ApJ, 698, 1761
  * Finke & Razzaque (2023) Finke J. D., Razzaque S., 2023, Astrophys. J. Lett., 942, L21
  * Finke et al. (2010) Finke J. D., Razzaque S., Dermer C. D., 2010, ApJ, 712, 238
  * Finke et al. (2022) Finke J. D., Ajello M., Domínguez A., Desai A., Hartmann D. H., Paliya V. S., Saldana-Lopez A., 2022, ApJ, 941, 33
  * Foster & Nelson (2009) Foster C., Nelson L. A., 2009, ApJ, 699, 1252
  * Fraija et al. (2022) Fraija N., Gonzalez M., HAWC Collaboration 2022, The Astronomer’s Telegram, 15675, 1
  * Frederiks et al. (2023) Frederiks D., et al., 2023, Astrophys. J. Lett., 949, L7
  * Furniss et al. (2013) Furniss A., et al., 2013, ApJ, 768, L31
  * Furniss et al. (2015) Furniss A., Sutter P. M., Primack J. R., Domínguez A., 2015, Monthly Notices of the Royal Astronomical Society, 446, 2267
  * Galanti et al. (2023) Galanti G., Nava L., Roncadelli M., Tavecchio F., Bonnoli G., 2023, Phys. Rev. Lett., 131, 251001
  * Gao et al. (2024) Gao L.-Q., Bi X.-J., Li J., Yao R.-M., Yin P.-F., 2024, JCAP, 01, 026
  * Gilmore et al. (2012) Gilmore R. C., Somerville R. S., Primack J. R., Domínguez A., 2012, MNRAS, 422, 3189
  * Gonzalez et al. (2023) Gonzalez M. M., Rojas D. A., Pratts A., Hernandez-Cadena S., Fraija N., Alfaro R., Araujo Y. P., Montes J. A., 2023, Astrophys. J., 944, 178
  * Gould & Schreder (1967) Gould R. J., Schreder G. P., 1967, Phys. Rev., 155, 1404
  * Guo et al. (2023) Guo S.-Y., Khlopov M., Wu L., Zhu B., 2023, Phys. Rev. D, 108, L021302
  * He et al. (2024) He H.-N., Zhang B. T., Fan Y.-Z., 2024, Astrophys. J., 963, 109
  * Hoffman & Shaham (1982) Hoffman Y., Shaham J., 1982, ApJ, 262, L23
  * Huang et al. (2022) Huang Y., Hu S., Chen S., Zha M., Liu C., Yao Z., Cao Z., Experiment T. L., 2022, GRB Coordinates Network, 32677, 1
  * Huang et al. (2023a) Huang J., Wang Y., Yu B., Zhou S., 2023a, JCAP, 04, 056
  * Huang et al. (2023b) Huang J.-K., Huang X.-L., Cheng J.-G., Ren J., Zhang L.-L., Liang E.-W., 2023b, ApJ, 947, 84
  * Kalashev et al. (2024) Kalashev O., Aharonian F., Essey W., Inoue Y., Kusenko A., 2024, arXiv e-prints, p. arXiv:2405.05402
  * Khangulyan et al. (2024) Khangulyan D., Aharonian F., Taylor A. M., 2024, ApJ, 966, 31
  * Klinger et al. (2023) Klinger M., Taylor A. M., Parsotan T., Beardmore A., Heinz S., Zhu S. J., 2023, Mon. Not. Roy. Astron. Soc., 529, L47
  * Kneiske & Dole (2010) Kneiske T. M., Dole H., 2010, A&A, 515, A19
  * Korochkin et al. (2019) Korochkin A., Rubtsov G., Troitsky S., 2019, JCAP, 12, 002
  * Kudoda & Faltenbacher (2017) Kudoda A. M., Faltenbacher A., 2017, Monthly Notices of the Royal Astronomical Society, 467, 2896
  * Kudoda & Faltenbacher (2018) Kudoda A. M., Faltenbacher A., 2018, Mon. Not. Roy. Astron. Soc., 481, 405
  * Lang et al. (2019) Lang R. G., Martínez-Huerta H., de Souza V., 2019, Physical Review D, 99
  * Lesage et al. (2023) Lesage S., et al., 2023, ApJ, 952, L42
  * Li & Ma (2021) Li H., Ma B.-Q., 2021, Journal of High Energy Astrophysics, 32, 1
  * Li & Ma (2023) Li H., Ma B.-Q., 2023, Astroparticle Physics, 148, 102831
  * Long et al. (2020) Long G. B., Lin W. P., Tam P. H. T., Zhu W. S., 2020, Phys. Rev. D, 101, 063004
  * Malesani et al. (2023) Malesani D. B., et al., 2023, arXiv e-prints, p. arXiv:2302.07891
  * Martínez-Huerta et al. (2018) Martínez-Huerta H., Lang R. G., de Souza V., 2018, PoS, BHCB2018, 010
  * Mbarubucyeye et al. (2023) Mbarubucyeye J. D., et al., 2023, PoS, ICRC2023, 705
  * Micheletti et al. (2014) Micheletti D., et al., 2014, Astron. Astrophys., 570, A106
  * Nakagawa et al. (2023) Nakagawa S., Takahashi F., Yamada M., Yin W., 2023, Phys. Lett. B, 839, 137824
  * Neyrinck (2008) Neyrinck M. C., 2008, Mon. Not. Roy. Astron. Soc., 386, 2101
  * Pan et al. (2012) Pan D. C., Vogeley M. S., Hoyle F., Choi Y.-Y., Park C., 2012, Monthly Notices of the Royal Astronomical Society, 421, 926
  * Paz et al. (2013) Paz D., Lares M., Ceccarelli L., Padilla N., Lambas D. G., 2013, Mon. Not. Roy. Astron. Soc., 436, 3480
  * Pérez de los Heros & Terzić (2023) Pérez de los Heros C., Terzić T., 2023, Lect. Notes Phys., 1017, 241
  * Qin et al. (2023) Qin L., Wang J., Gao Q., Na W., Li H., Wang A., Yang C., Yang J., 2023, Mon. Not. Roy. Astron. Soc., 521, 6219
  * Razzaque et al. (2009) Razzaque S., Dermer C. D., Finke J. D., 2009, ApJ, 697, 483
  * Ricciardelli et al. (2014) Ricciardelli E., Quilis V., Varela J., 2014, Mon. Not. Roy. Astron. Soc., 440, 601
  * Rojas et al. (2024) Rojas D. A., Hernández-Cadena S., González M. M., Pratts A., Alfaro R., Serna-Franco J., 2024, The Astrophysical Journal, 966, 114
  * Ror et al. (2024) Ror A. K., Gupta R., Aryan A., Bhushan Pandey S., Oates S. R., Castro-Tirado A. J., Kumar S., 2024, arXiv e-prints, p. arXiv:2406.01220
  * Saldana-Lopez et al. (2021) Saldana-Lopez A., Domínguez A., Pérez-González P. G., Finke J., Ajello M., Primack J. R., Paliya V. S., Desai A., 2021, Monthly Notices of the Royal Astronomical Society, 507, 5144
  * Sánchez-Conde et al. (2009) Sánchez-Conde M. A., Paneque D., Bloom E., Prada F., Domínguez A., 2009, Phys. Rev. D, 79, 123511
  * Sheth & van de Weygaert (2004) Sheth R. K., van de Weygaert R., 2004, Mon. Not. Roy. Astron. Soc., 350, 517
  * ShivaSankar K. et al. (2023) ShivaSankar K. A., Das A., Lambiase G., Nomura T., Orikasa Y., 2023, arXiv e-prints, p. arXiv:2308.14483
  * Sinha et al. (2014) Sinha A., Sahayanathan S., Misra R., Godambe S., Acharya B. S., 2014, ApJ, 795, 91
  * Smirnov & Trautner (2023) Smirnov A. Y., Trautner A., 2023, Phys. Rev. Lett., 131, 021002
  * Sutter et al. (2015) Sutter P., et al., 2015, Astronomy and Computing, 9, 1
  * Tavasoli et al. (2013) Tavasoli S., Vasei K., Mohayaee R., 2013, Astron. Astrophys., 553, A15
  * Tavecchio, F. & Bonnoli, G. (2016) Tavecchio, F. Bonnoli, G. 2016, A&A, 585, A25
  * Terzić et al. (2021) Terzić T., Kerszberg D., Strišković J., 2021, Universe, 7, 345
  * The LHAASO collaboration (2024) The LHAASO collaboration 2024, arXiv e-prints, p. arXiv:2402.06009
  * Troitsky (2021) Troitsky S., 2021, Eur. Phys. J. C, 81, 264
  * Troitsky (2022) Troitsky S. V., 2022, Pisma Zh. Eksp. Teor. Fiz., 116, 745
  * Troitsky (2024) Troitsky S., 2024, JCAP, 01, 016
  * Wang & Ma (2023) Wang L., Ma B.-Q., 2023, Phys. Rev. D, 108, 023002
  * Wang et al. (2023) Wang K., Ma Z.-P., Liu R.-Y., Zou Y.-C., Li Z., Dai Z.-G., 2023, Science China Physics, Mechanics, and Astronomy, 66, 289511
  * Williams et al. (2023) Williams M. A., et al., 2023, ApJ, 946, L24
  * Wojtak et al. (2016) Wojtak R., Powell D., Abel T., 2016, Mon. Not. Roy. Astron. Soc., 458, 4431
  * Xia et al. (2024) Xia Z.-Q., Wang Y., Yuan Q., Fan Y.-Z., 2024, arXiv e-prints, p. arXiv:2210.13052
  * Yang et al. (2024) Yang Y.-M., Bi X.-J., Yin P.-F., 2024, J. Cosmology Astropart. Phys., 2024, 060
  * Zhang et al. (2023) Zhang B. T., Murase K., Ioka K., Song D., Yuan C., Mészáros P., 2023, Astrophys. J. Lett., 947, L14
  * Zhang et al. (2024a) Zhang Y.-Q., et al., 2024a, arXiv e-prints, p. arXiv:2403.12851
  * Zhang et al. (2024b) Zhang Z., Lin H.-X., Li Z., Xiong S.-L., Zhang Y.-Q., Zhang Q.-y., Yi S.-X., 2024b, arXiv e-prints, p. arXiv:2405.12977
  * Zhao et al. (2023) Zhao Z.-C., Zhou Y., Wang S., 2023, Eur. Phys. J. C, 83, 92
  * Zheng et al. (2023) Zheng Y. G., Kang S. J., Zhu K. R., Yang C. Y., Bai J. M., 2023, Phys. Rev. D, 107, 083001
  * Zheng et al. (2024) Zheng C., et al., 2024, Astrophys. J. Lett., 962, L2
  * van de Weygaert & Platen (2011) van de Weygaert R., Platen E., 2011, in International Journal of Modern Physics Conference Series. pp 41–66 (arXiv:0912.2997), doi:10.1142/S2010194511000092