# Asymmetries in the simulated ozone distribution on TRAPPIST-1e due to
orography

Anand Bhongade  School of Earth and Environment, University of Leeds, Leeds,
UK  Daniel R Marsh  School of Physics and Astronomy, University of Leeds,
Leeds, UK  School of Chemistry, University of Leeds, Leeds, UK  Felix
Sainsbury-Martinez  School of Physics and Astronomy, University of Leeds,
Leeds, UK  Gregory Cooke  School of Physics and Astronomy, University of
Leeds, Leeds, UK  Institute of Astronomy, University of Cambridge, UK

###### Abstract

TRAPPIST-1e is a tidally locked rocky exoplanet orbiting the habitable zone of
an M dwarf star. Upcoming observations are expected to reveal new rocky
exoplanets and their atmospheres around M dwarf stars. To interpret these
future observations we need to model the atmospheres of such exoplanets. We
configured CESM2-WACCM6, a chemistry climate model, for the orbit and stellar
irradiance of TRAPPIST-1e assuming an initial Earth-like atmospheric
composition. Our aim is to characterize the possible ozone (O3) distribution
and explore how this is influenced by the atmospheric circulation shaped by
orography, using the Helmholtz wind decomposition and meridional mass
streamfunction. The model included Earth-like orography and the substellar
point was located over the Pacific Ocean. For such a scenario, our analysis
reveals a North-South asymmetry in the simulated O3 distribution. The O3
concentration is highest below 10 hPa level (below $\sim$30 km) near the South
pole. This asymmetry results from the presence of land masses on the night
side that cause drag in near-surface flows and lead to an asymmetric
meridional overturning circulation. Catalytic species were roughly
symmetrically distributed and were not found to be primary driver for the O3
asymmetry. The total ozone column (TOC) density was higher for TRAPPIST-1e
compared to Earth, with 8000 Dobson Units (DU) near the South pole and 2000 DU
near the North pole. The results emphasize the sensitivity of O3 to model
parameters, illustrating how incorporating Earth-like orography can affect
atmospheric dynamics and O3 distribution. This link between surface features
and atmospheric dynamics underlines the importance of how changing model
parameters used to study exoplanet atmospheres can influence the
interpretation of observations.

Exoplanet atmospheres (487); Exoplanet atmospheric composition (2021);
Exoplanet atmospheric dynamics (2307); Transmission spectroscopy (2133); James
Webb Space Telescope (2291)

##  1 Introduction

Scientists have long wondered about the existence of life on other planets and
this curiosity has motivated them to explore celestial bodies beyond our solar
system. To date, over 5600 confirmed exoplanets have been identified by
instruments like the Kepler Space Telescope, the Transiting Exoplanet Survey
Satellite, and various space-based and ground-based
telescopes.111https://exoplanetarchive.ipac.caltech.edu/docs/counts_detail.html
\- Date: 15/04/2024 Studies utilizing Kepler mission data have found that the
presence of small rocky exoplanets around M dwarf stars exceeds that around
Sun-like stars (Howard et al., 2012; Dressing & Charbonneau, 2015; Mulders et
al., 2015; Gaidos et al., 2016), although this is likely because of the
observational bias of telescopes for detecting exoplanets around smaller,
cooler and dimmer stars. To date, 200 confirmed rocky exoplanets have been
detected.222https://exoplanets.nasa.gov/what-is-an-exoplanet/planet-
types/terrestrial/ \- Date: 15/04/2024

M dwarf stars, constituting approximately 70% of all known stars in our
galaxy, have garnered scientific attention due to their abundance and compact
planetary systems (Bochanski et al., 2010). The location of the habitable zone
around a star depends upon its stellar properties (Huang, 1959). M dwarf
stars, characterized by their low temperature and flux, have habitable zones
situated closer to them. Detection methods such as transit photometry and
radial velocity rely on the planet-to-star mass and size ratio, making the
search for rocky exoplanets in habitable zones around M dwarf stars more
feasible due to the star’s relatively small size (Shields et al., 2016;
Reiners et al., 2018; Gould et al., 2003; Nutzman & Charbonneau, 2008).

Exoplanets orbiting M dwarfs can be in such proximity to their stars that they
become tidally locked, resulting in a permanent day side and night side. Tidal
locking occurs due to the gravitational force exerted by the star, which
distorts the planet into an elongated shape. This results in synchronous
rotation, where the planet’s rotational period equals its orbital period
(Barnes, 2017).

In addition to stellar properties, the habitable zone is significantly
influenced by atmospheric composition, particularly the presence or absence of
greenhouse gases. Ozone (O3), acting as a greenhouse gas on Earth, influences
the vertical temperature structure and humidity of the atmosphere. In absence
of O3, the average surface temperature of an Earth-like planet would be 7K
cooler (Gómez-Leal et al., 2019). Furthermore, O3 serves as a filter for
incoming UV radiation, safeguarding life on Earth’s surface. Therefore, O3
also plays a role in determining the potential habitability of an exoplanet.

To determine whether an exoplanet harbors life, scientists are keen on
detecting biosignatures in its atmosphere. This detection can be achieved
through direct imaging or transit spectroscopy. Direct imaging is particularly
suitable for larger exoplanets in wider orbits around bright and massive
stars. However, detecting small rocky planets through direct imaging poses a
challenge due to the overshadowing effect of their host stars, especially with
current technology, and also because habitable zones for M-dwarf stars are
very close to them.

Scientists have considered using O3 as a proxy to detect molecular oxygen (O2)
in the atmosphere of an exoplanet. However, a study conducted by Kozakis et
al. (2022) has highlighted the challenge in accurately determining the amount
of O2 through O3 measurements, especially when the UV flux of the host star is
uncertain (Cooke et al., 2023). Nonetheless, if we have knowledge about the
host star’s UV spectrum and certain other information, O3 measurements might
offer insights into the potential habitability of an exoplanet.

On Earth O3 is primarily produced in the tropical stratosphere via the Chapman
mechanism (Chapman, 1930). During daytime, UV radiation breaks the O2 molecule
to form two O atoms which then reacts with O2 molecules and a third body M
through a three-body process to form O3. The Brewer-Dobson circulation
distributes O3 in the stratosphere (Dobson, 1956; Brewer, 1949; Newell, 1963).
In this circulation pattern, air near the tropics rises and then moves toward
the poles. The Brewer-Dobson circulation slightly differs between the two
hemispheres due to variances in land and ocean structure in the Southern and
the Northern hemispheres (Remsberg, 2015).

This O3 then either gets photodissociated by UV radiation to form O and O2, or
it reacts with O to form two O2 molecules. O3 also gets catalytically
destroyed by NOx, HOx, Br, Cl, etc (Portmann et al., 2012). If X is a
catalytic species, then the O3 is catalytically destroyed by the following
reactions:

$\displaystyle\qquad\text{X + O}_{3}\rightarrow\text{XO + O}_{2}$
$\displaystyle\qquad\text{XO + O}\rightarrow\text{X + O}_{2}$
$\displaystyle\text{Net:}\quad\text{O}_{3}+\text{O}\rightarrow 2\text{O}_{2}$

In 2016, three Earth-sized planets (TRAPPIST-1b, c, and d) were detected
orbiting an ultra-cool M dwarf star known as TRAPPIST-1 or 2MASS J23062928 -
0502285, through observations made by the TRAPPIST (TRAnsiting Planets and
PlanetIsimals Small Telescope) (Gillon et al., 2016). Subsequently, in 2017,
Spitzer revealed the existence of four more exoplanets around TRAPPIST-1
(TRAPPIST-1e, f, g, and h), establishing it as the first system with seven
Earth-sized planets (Gillon et al., 2017). All the planets in the TRAPPIST-1
system are closer to their host star than Mercury is to the Sun.

Future observations with JWST (Gardner et al., 2006) and new telescopes, such
as ELT (Extremely Large Telescope) (Hook, 2009) are expected to reveal more
about small rocky exoplanets and their atmospheres. Looking ahead, future
concepts like the Habitable Worlds Observatory (National Academies of
Sciences, Engineering, and Medicine, 2021), which was proposed by combining
two earlier concepts - the HabEx (Habitable Exoplanet Observatory Mission)
(Gaudi et al., 2020) and the LUVOIR (Large UV/Optical/IR Surveyor) (LUVOIR
Team, 2019) \- are anticipated to directly image and characterize the
atmospheres of Earth-like exoplanets. Hence, it becomes crucial to model rocky
exoplanet atmospheres to aid in interpreting observations by extracting the
atmospheric properties and explaining the underlying physical processes
occurring in these exoplanetary atmospheres.

In previous studies, general circulation models (GCMs) have been employed to
investigate the atmospheres of Earth-like planets, including studies designed
to understand the atmospheric circulation and O3 chemistry of tidally locked
Earth-like exoplanets. For example, the TRAPPIST-1 Habitable Atmosphere
Intercomparison (THAI) project (Turbet et al., 2022; Sergeev et al., 2022)
compared the results from four GCMs which included slab oceans for both dry
and moist N2-dominated and CO2-dominated atmospheres. Carone et al. (2018)
used the MITgcm to study the stratospheric circulation of a tidally locked
ExoEarth scenario for TRAPPIST-1b, TRAPPIST-1d, Proxima Centauri b and GJ 667
C f. Yates et al. (2020) employed the Met Office Unified Model to explore the
O3 chemistry of the tidally locked exoplanet Proxima Centauri b around an M
dwarf, while Proedrou & Hocke (2016) used the CESM1(WACCM) model to simulate
the 3D O3 distribution of a tidally-locked Earth-like planet around a Sun-like
star. A recent study conducted by Braam et al. (2023) utilized a slab ocean
model of a tidally locked exoplanet around an M dwarf star with Proxima
Centauri b parameters to study O3 spatial distribution. They found that O3
accumulates on the night side, demonstrating a day side-night side hemispheric
asymmetry in O3 distribution.

In our study we incorporated an Earth-like land-ocean structure to examine for
the first time how orography might alter the atmospheric dynamics and
chemistry of a TRAPPIST-1e exoplanetary atmosphere scenario. Our aim is to
characterize the possible O3 distribution and the influence of atmospheric
circulation on it. We focus on O3 because it is affected by photochemistry,
catalytic cycles, atmospheric transport, and has strong spectral features from
the UV to the mid-infrared. Additionally, the presence of O2, a potential
biosignature, can be inferred from a detection of O3.

##  2 Methods

###  2.1 Model description and setup

For this study, we modeled TRAPPIST-1e using the Whole Atmosphere Community
Climate Model version 6 (WACCM6)(Gettelman et al., 2019). WACCM6 is an
atmospheric model that operates as a configuration of the Community Earth
System Model version 2 (CESM2). CESM2 is an Earth system model consisting of
submodels simulating the atmosphere, ocean, land, sea ice, land ice, river
runoff, and surface waves (Danabasoglu et al., 2020). The WACCM6 configuration
comprises 70 vertical levels, starting from the surface at 1000 hPa and
extending up to 140 km at 4.5 × 10-6 hPa (lower thermosphere), with a
horizontal resolution of 1.875° latitude × 2.5° longitude.

The chemistry applied in WACCM6 is based on the Model of Ozone and Related
Chemical Tracers (MOZART). MOZART serves as a global chemical transport model
encompassing physical and chemical processes that span the troposphere,
stratosphere, mesosphere, and lower thermosphere (Emmons et al., 2020; Conley
et al., 2012). WACCM6 employs the Rapid Radiative Transfer Model (RRTMG)
radiation code for solving longwave (lower energy radiation emitted from the
Earth’s surface and atmosphere) and shortwave (higher energy wavelengths such
as UV, visible light, and a part of the near-infrared spectrum, associated
with the host star) radiative transfer equations.

WACCM6 has been previously utilized to study the climate and O3 chemistry of
prehistoric Earth and Earth-like exoplanets. For example, Cooke et al. (2022)
demonstrated how the total O3 column at O2 concentrations between 0.1 - 50%
the present atmospheric level, may have been lower than previously thought.
The possible reasons behind this result are investigated in Ji et al. (2023).
Cooke et al. (2023) showed how for an Earth-like exoplanet, observations and
spectral signatures of chemical species could be affected by the line of
sight, albedo, clouds, and chemistry. Liu et al. (2023) studied how changing
the eccentricity of an Earth-like exoplanet affects the abundance and loss of
water present in the atmosphere.

The model considered here used the BWma1850
compset333https://docs.cesm.ucar.edu/models/cesm2/config/2.1.3/compsets.html
of WACCM6, which included a pre-industrial Earth-like atmosphere and
orography, and was modified to allow for synchronous
rotation.444https://github.com/exo-
cesm/CESM2.1.3/tree/main/Tidally_locked_exoplanets/cases The model composition
was with O2 set to present atmospheric level (PAL), i.e., 21% by volume, and
N2 with 78% by volume. The volume mixing ratios of CH4 (0.8 ppmv), CO2 (280
ppmv), N2O (270 ppbv), and H2 (500 ppbv), are fixed at the surface. We ran the
simulation for 300 years, of which we utilised the last 40 years of data in
order to eliminate any effects associated with the model adjusting to tidally-
locked conditions. The substellar point was fixed at 180° longitude and 0°
latitude over the Pacific Ocean (see Figure 1). Table 1 shows the parameters
used to model TRAPPIST-1e. For this study, we have used a stellar spectrum
based on the work of Peacock et al. (2019) who modelled the stellar energy
distribution of TRAPPIST-1 and produced model 1A, 2A and 2B of which we used
model 1A, which best matched the TRAPPIST-1 Ly$\alpha$ reconstruction from
Bourrier et al. (2017). Figure 2 compares the top-of-the-atmosphere
irradiation of TRAPPIST-1e and Earth.

Figure 1: Stellar irradiation received by TRAPPIST-1e in our model. The black dot represent the substellar point which is fixed at 180° longitude and 0° latitude over the Pacific Ocean. The black outlines represent the land masses. Parameter | Units | Value
---|---|---
Semimajor axis | au | 0.029
Orbital period | Earth days | 6.1
Rotation period | Earth days | 6.1
Obliquity |  | 0
Eccentricity |  | 0
Instellation | W m-2 | 900
Planet radius | km | 5797
Gravity | m s-2 | 9.14
Table 1: Planetary parameters used for model of TRAPPIST-1e (Grimm et al.,
2018). Figure 2: Stellar irradiance at top of the atmosphere of TRAPPIST-1e
(orange) and Earth (blue). The orange shaded region represents the range of UV
wavelengths.

###  2.2 Data analysis

For our analysis, we calculated the time average of the 40 years of data and
examined various parameters, including surface temperature, horizontal wind
velocities ($u$ and $v$), vertical wind velocity ($w$), reaction rates, and
volume mixing ratios of O3, O2, OH, HO2, NO, NO2, Br, and Cl. To calculate the
concentration of chemical species in units of molecules m-3, we multiplied the
air number density by the volume mixing ratio of the respective chemical
species.

We calculated the Ox production rate using the formula $J_{O_{2}}$*[$O_{2}$],
where $J_{O_{2}}$ is the oxygen photolysis rate constant and [$O_{2}$] is the
O2 concentration. The catalytic destruction rates of O3 by various chemical
species (X) were calculated using the formula $k$*[X]*[O3], where $k$ is the
rate constant for X + O3 reactions, [X] is the concentration of the catalytic
species, and [O3] is the concentration of O3. In this context, X can be OH,
HO2, NO, NO2, Br, or Cl.

To study the horizontal wind circulation we used the Helmholtz wind
decomposition, which has been previously used to study the atmospheric
circulations of Earth, tidally locked terrestrial exoplanets and hot Jupiters
(Hammond & Lewis, 2021). The Helmholtz wind decomposition breaks down the
total horizontal wind into rotational ($\boldsymbol{u_{r}}$) and divergent
components ($\boldsymbol{u_{d}}$), allowing us to study the winds responsible
for the overall circulation in more detail. For a tidally locked planet, the
rotational flow consists of two parts: an equatorial jet (zonal-mean
rotational component) which is a narrow belt of winds moving around the
equator of a planet and stationary waves (eddy rotational component). The
stationary waves drive the equitorial jet and can accelerate it resulting in a
superrotating jet which moves faster than the planet’s rotation rate (Showman
et al., 2013). The divergent flow consists of a global overturning circulation
which is a result of thermally driven circulation where the air rises on the
day side of the planet and sinks on the night side. Mathematically, the
Helmholtz wind decomposition divides the total circulation given by
$\boldsymbol{u}$ = ($u$,$v$) into the two aforementioned components (Hammond &
Lewis, 2021):

$\boldsymbol{u}=\boldsymbol{u_{d}}+\boldsymbol{u_{r}}=-\boldsymbol{\nabla}\chi+\boldsymbol{\mathrm{k}}\times\boldsymbol{\nabla}\psi$
(1)

$\nabla{{}^{2}}\chi=\delta$ (2) $\nabla{{}^{2}}\psi=\zeta$ (3)

Here, $\boldsymbol{u}$ represents the horizontal wind vector, where $u$ and
$v$ denote the zonal (along latitude) and meridional (along longitude)
velocities, respectively. The velocity potential function, $\chi$, is derived
from the divergence ($\delta$ \- Equation 2), while the streamfunction,
$\psi$, is derived from the vorticity ($\zeta$ \- Equation 3) of the wind.

To examine the meridional overturning circulation, we used the meridional mass
streamfunction. This is commonly used to study atmospheric or oceanic
circulation patterns and helps in visualizing mass transport in the meridional
plane (North-South direction). Mathematically, the relationship between mass
streamfunction $\psi$ and both vertical velocity ($w$) and meridional velocity
($v$) is expressed by the equations (Sidorenko et al., 2020)

$\frac{1}{R}\frac{\partial\psi}{\partial\theta}=w,\quad\frac{\partial\psi}{\partial
z}=-v$ (4)

Here, $\theta$ is the latitude in radians, $z$ is the altitude and $R$ is the
planet’s radius.

##  3 Results and Discussion

###  3.1 Ozone distribution

Figure 3 a) and b) show the O3 concentration on the day and the night side of
our TRAPPIST-1e model. O3 is predominantly present in the lower atmosphere,
below 10 hPa level (below $\sim$ 30 km) and O3 concentrations are higher near
the poles than near the equator. The highest O3 concentrations are found near
the South pole between 600 and 70 hPa levels ($\sim$5 and $\sim$15 km).
Further, there is little difference in O3 concentrations between the day-side
and the night-side, much less than the North-South asymmetry, suggesting that
the O3 produced on the day side is transported to the night side by horizontal
winds.

Figure 4 a) and b) display the O3 concentrations at 103 hPa ($\sim$15 km) and
609 hPa ($\sim$5 km) levels respectively. At both of these pressures, the
North-South asymmetry in O3 distribution is clearly visible. On Earth, O3 is
distributed in a roughly symmetrical manner about the equator, with the
majority located between 300 and 10 hPa levels. Between 30 and 10 hPa levels,
O3 is densely packed around low latitudes, while between 300 and 30 hPa
levels, O3 is more abundant around the poles (see Ejzak et al. (2007)).
Therefore, the vertical structure and the latitudinal distribution of O3 in
our TRAPPIST-1e model differs significantly from that of Earth.

Figure 3: Cross section of the O3 number density on the meridional plane
passing through the a) substellar point and b) antistellar point for our
TRAPPIST-1e model. Figure 4: Horizontal slices of the O3 number density for
our TRAPPIST-1e model at a pressure (altitude) of a) 103 hPa ($\sim$15 km) and
b) 609 hPa ($\sim$5 km). The substellar point is at 180° longitude and 0°
latitude. N represents the North pole and S represents the South pole. Lower
panels show the polar projections.

Figure 5 shows the Ox production rate (O + O3 production rate) on the
meridional plane passing through the substellar point. The Ox production rate
is symmetric, with peak production occurring in the upper atmosphere (above 10
hPa level or above $\sim$30 km) due to the large amount of UV radiation
received here. UV radiation does not penetrate the lower atmosphere because it
becomes denser, absorbing most of the incoming UV radiation, and it is weakest
at high latitudes due to geometric effects. Hence, in the lower atmosphere,
there is less Ox production near the equator and no production at higher
latitudes (as seen in figure 5, where white areas indicate regions of zero
production). When comparing the pattern of O3 concentration and Ox production
rate, it is evident that the region of peak Ox production rate does not
coincide with the region of peak O3 concentration. There is relatively low O3
concentrations where the Ox production is high, whereas the region with the
highest O3 concentration (near the South pole between 600 and 70 hPa levels)
has either low or no Ox production. This indicates that O3 is primarily
produced high up in the atmosphere and either all the O3 is transported from
the upper atmosphere to the lower atmosphere, or some of it is transported
while the remainder is catalytically destroyed. However, this does not explain
the North-South asymmetry in O3 concentration. To explain this asymmetry, we
need to examine if there are asymmetries in catalytic species distribution and
atmospheric circulation.

Figure 5: Cross section of Ox production rate on the meridional plane passing
through the substellar point for our TRAPPIST-1e model. The white region on
the plot indicates zero production rate.

###  3.2 Role of catalytic species

Figure 6 shows the distribution of catalytic species (OH, HO2, NO, NO2, Cl,
and Br) on the meridional plane passing through the substellar point for our
TRAPPIST-1e model. For regions below 10 hPa level ($\sim$30 km), our initial
hypothesis, after observing the symmetric Ox production and the asymmetric
distribution of O3, was that there might be an uneven distribution of
catalytic species. Such an imbalance could potentially be responsible for O3
depletion being more pronounced near the North pole compared to the South
pole. However, our analysis suggests that this is not the case, with the
distribution of catalytic species being roughly symmetric and the
concentrations low when compared to the O3 concentration. In the case of Cl
(Figure 6 e) and Br (Figure 6 f) we find that below 100 hPa level the
concentrations of these species is slightly higher near the South pole as
compared to the North pole, since this is aligned with the peak in ozone
concentration, it cannot explain the high concentrations of O3 we find below
10 hPa level near the South pole.

Figure 7 shows the time taken by the catalytic species to destroy O3 on the
day side. We observe that in the regions with the highest concentrations of
O3, catalytic species require more time to destroy O3 as compared to the
regions where O3 is produced. This observation implies that below the 10 hPa
level, catalytic destruction proceeds at a slower rate, allowing sufficient
time for atmospheric circulation to influence the distribution of O3. As such,
this suggests that catalytic species may not be the primary factor driving the
asymmetry in the distribution of O3. Therefore, we next look at the transport
of O3 throughout the atmosphere to see if this might have a more significant
effect.

Figure 6: Catalytic species concentration on the meridional plane passing
through the substellar point for a) OH, b) HO2, c) NO, d) NO2, e) Cl and f) Br
for our TRAPPIST-1e model. The catalytic species are roughly symmetrically
distributed and their concentrations are lower than the ozone concentrations.
Figure 7: Time taken by the catalytic species to destroy O3 on the meridional
plane passing through the substellar point for our TRAPPIST-1e model.

###  3.3 Role of atmospheric circulation

Figure 8 a) and b) shows the meridional mass streamfunction on a plane passing
through the substellar point and the antistellar point, respectively. Here,
red indicates clockwise transfer of mass and blue indicates anticlockwise
transfer of mass. On the day side, without orography we would expect a
symmetric circulation with upwelling at the equator, meridional flow towards
the poles, downwelling near the poles and meridional flow from the poles to
the equator (see Sergeev et al. (2022)). In our simulation, we find that there
is an asymmetry in this circulation: near the North pole the downwelling is
interrupted at 200 hPa level and does not reach the surface, rather there is a
meridional transport from the North pole towards mid latitudes, downwelling at
mid latitude and then meridional flow from mid latitude towards the equator
(Figure 8 a). This asymmetry seems to be caused by the presence of a small
anticlockwise (blue) circulation cell near the North pole. Similarly, on the
night side, without orography we would expect to find upwelling near the
poles, meridional transport from the poles to the equator, downwelling at the
equator and meridional transport from the equator towards the poles near
surface, thus closing the global overturning circulation (see Sergeev et al.
(2022)). However, figure 8 b) reveals there is an asymmetry in this
circulation. Near the North pole there are two circulation cells: below 200
hPa level, there is upwelling from the North pole, meridional transport from
the North pole to the equator, downwelling near the equator until 200 hPa
level, and meridional transport back to the North pole. Above 200 hPa level,
there is downwelling near the North pole, meridional transport from the North
pole, across the equator, to the South pole, and then upwelling near the South
pole.

So, both on the day side and night side we find an asymmetry in the meridional
overturning circulation, near the North pole due to the presence of a small
anticlockwise circulation cell on the day side and a large clockwise
circulation cell on the night side. This asymmetry is caused by the presence
of Earth-like orography. Both on the day side and the night side the Northern
hemisphere has higher fraction of land mass as compared to the Southern
hemisphere. This landmass helps to reshape the wind, leading to the asymmetric
meridional overturning circulation near the North pole.

The horizontal flow is complex, but the Helmholtz wind decomposition gives us
an easier way to understand it. Figure 9 and Figure 10 show the divergent
component and the eddy rotational component of the horizontal winds in the
stratosphere (20 hPa level), near the tropopause (103 hPa level), and near the
surface (800 hPa level). When comparing the mean wind speeds of both
components at each pressure level, we find that the eddy rotational component
has a higher mean wind speed than the divergent component. This implies that
the eddy rotational component plays a major role in meridional transport.

We can explain the O3 concentrations North-South asymmetry using both the
meridional overturning circulation and the Helmholtz wind decomposition. O3
formed on the day side between 200 and 10 hPa levels ($\sim$10 and $\sim$30
km) is evenly transported toward the poles by meridional flow and then from
the day side to the night side by the rotational component of horizontal
winds, which includes both the eddy and zonal-mean components. On the night
side, in the lower atmosphere (below 100 hPa level), the meridional flow (see
Figure 10 c) forces large amounts of O3 from the North pole to the South pole.
The North-South asymmetry on the night side in the meridional flow, near the
surface, occurs due to the presence of Earth-like orography. Land-ocean
boundaries cause wave breaking (not shown), resulting in a flow from the North
pole towards the South pole disrupting the symmetric day/night transport found
by Braam et al. (2023). Also, near the North pole, below 100 hPa level, the
presence of a small anticlockwise circulation cell on the day side and large
clockwise circulation cell on the night side in the lower atmosphere, prevents
O3 from accumulating at higher latitudes.

Considering, Figure 5 and Figure 10 a) and b), we find that the O3 formed at
the substellar point between 200 and 10 hPa levels ($\sim$10 and $\sim$30 km)
is first pushed towards the 120° longitude, then towards the poles, and from
the day side to the night side by the rotationally driven winds. Below 200 hPa
level ($\sim$30 km) on the night side (300° – 360° longitude), the flow from
the North pole to the South pole (see Figure 10 c) forces O3 towards the South
pole which then gets well mixed by the rotational component of the horizontal
winds resulting in high O3 concentrations in the lower atmosphere near the
South pole on both the day side and the night side.

Figure 8: Meridional mass streamfunction on the meridional plane passing
through the a) substellar point and the b) antistellar point for our
TRAPPIST-1e model. The positive values (red) indicates clockwise circulation
and the negative values (blue) indicates anticlockwise circulation. The values
are in log-scale. Figure 9: Helmholtz divergent winds at a) 20 hPa, b) 103 hPa
and c) 800 hPa levels. Arrows represents winds velocity vectors and point in
the direction of flow, with the length of the arrows representing the
magnitude of the wind speed (with the mean wind speed shown in the top right
corner of each plot). Figure 10: Helmholtz rotational eddy winds at a) 20 hPa,
b) 103 hPa and c) 800 hPa. Arrows represent winds velocity vectors and point
in the direction of flow, with the length of the arrows representing the
magnitude of the wind speed (with the mean wind speed shown in the top right
corner of each plot).

###  3.4 Total ozone column density of TRAPPIST-1e and Earth

Total ozone column (TOC) density is reported in Dobson Units (DU), where 1 DU
is equivalent to a 10 $\mu$m thick layer of pure O3 at 273 K and 1 atm
pressure. For Earth the values of TOC are 200-250 DU near the South pole and
300-350 DU near the North pole.555https://ozonewatch.gsfc.nasa.gov/SH.html \-
Date: 24/08/2023 Figure 11 shows the latitudinal variation in zonal mean TOC
density for TRAPPIST-1e and Earth. From this figure, we clearly observe the
North-South asymmetry in O3 distribution for our TRAPPIST-1e model, whereas O3
on Earth is roughly symmetrically distributed about the equator. It is also
interesting to note that the total O3 concentration on TRAPPIST-1e is much
higher compared to that of Earth. The TOC density on TRAPPIST-1e near the
South pole is 8000 DU, which is around 28 times that of Earth ($\sim$250 DU),
and near the North pole, it is approximately 2000 DU, which is 7 times that of
Earth ($\sim$300 DU). There is slight asymmetry in TOC on Earth with slightly
higher TOC near the North pole as compared to the South pole, but for
TRAPPIST-1e this asymmetry is reversed, and significantly larger.

From our results we can say that TRAPPIST-1e has the potential to exhibit
large O3 concentrations, which could be potentially observable. Furthermore,
if O3 is detected, this would be indicative of atmospheric oxygen. Chen et al.
(2019) demonstrated through simulated transmission spectra of M-dwarf planets
that JWST could potentially detect O3 features during primary transit.
However, detecting the prominent O3 features would require over 100 transits
in conditions of zero cloud coverage, and twice that number in instances of
100% cloud coverage (Lustig-Yaeger et al., 2019). There are numerous factors
which influence the O3 spectral signature, such as the line of sight of the
telescope and the presence of other chemical species (Proedrou & Hocke, 2016).
The brightness of an exoplanet changes as it orbits its host star and is
depicted by the full-phase light curve. This curve shows the entire range of
phases from minimum to maximum illumination (Seager, 2010). The presence of
land and ocean masses can lead to non-uniform O3 distributions as found in our
model and hence as the planet orbits its host star, the amounts of O3 that is
visible might also change. Hence, the detection of O3 might depend on the
phase of the exoplanet we observe. The presence of other chemical constituents
in the atmosphere could also hinder O3 detection by affecting its spectral
signature. Chemical constituents which have similar spectral features to O3
could overlap with an O3 spectral features making it difficult to detect O3.
Finally, TRAPPIST-1 is a highly active star with dark spots and bright faculae
on its surface which could make O3 detection difficult by interfering with the
transmission spectrum of the planet (Lim et al., 2023).

Figure 11: Zonal mean latitudinal variation of total ozone column density for
our TRAPPIST-1e model (blue) and Earth (orange).

###  3.5 Limitations and Future Work

The limitation of our study is that we considered an Earth-like orography and
atmospheric composition, which is likely not to be the case in reality. Slab
ocean are also unlikely to be physical. In a follow up study we would like to
look into detail how the land structures disrupts the near surface winds. We
would also like to study how changing the location of the substellar point
(with respect to the land-mass distribution) affects the distribution and
overall concentration of O3, atmospheric dynamics and chemistry and how would
these might influence observations.

##  4 Summary

In this study we used the CESM2-WACCM6 GCM which includes Earth-like orography
to simulate the tidally locked exoplanet TRAPPIST-1e and investigate the O3
distribution and how the atmospheric circulations shape it. The atmospheric
composition was set to pre-industrial Earth-like composition (O2-N2
dominated). The substellar point was positioned over the ocean at 180°
longitude and 0° latitude, and the stellar irradiance was based on the
TRAPPIST-1 spectrum of Peacock et al. (2019). We analyzed the time average of
40 years of steady state simulation data for our study.

Our results reveal a North-South asymmetry in O3 distribution. O3 is highly
concentrated near the South pole below 10 hPa level (below $\sim$30 km). This
does not relate to Ox production, which is symmetric at low pressures and is
least prevalent in the deep atmosphere, especially near the poles. On the day
side, Ox generation is symmetrical and reaches its peak above 10 hPa level
(above $\sim$30 km). Below 10 hPa level (below $\sim$30 km), Ox production is
minimal in low latitudes and lacking or minimal in high latitudes. This
suggests that where high O3 concentrations exist, Ox production is either
significantly reduced or entirely absent. On examining the distribution of
catalytic species (OH, HO2, NO, NO2, Br, and Cl) we found a relatively
symmetrical distribution. This suggests that catalytic destruction is not
primarily responsible for asymmetric O3 distribution.

Consequently we next investigated if O3 transport is responsible for the O3
asymmetry. Indeed after examining the atmospheric circulation patterns, we
identified an asymmetry in the meridional overturning circulation near the
North pole. This discrepancy is caused by wave breaking due to land-ocean
boundaries, impacting near surface winds and resulting in a flow from the
North pole to the South pole on the night side. Analyzing the Helmholtz wind
decomposition revealed that the rotational component, specifically the eddy
wind component dominated over the divergent component of the wind i.e. the
wave driven transport dominates the global overturning circulation.

Upon comparing meridional circulation on the day side with Ox production, we
observed that O3 generated above 10 hPa level (above $\sim$30 km) near low
latitudes is carried towards the poles by this circulation. Subsequently,
horizontal winds (linked with the rotational component) carry O3 from the day
side to the night side below 100 hPa level (below $\sim$10 km) and on the
night side, a North-to-South flow near the surface forces O3 near the South
pole which gets distributed zonally by a wind gyre at the South pole. The eddy
winds which drives the meridional transport is what leads to asymmetry in the
observed meridional overturning circulation and O3 distribution.

The asymmetry in O3 distribution for our TRAPPIST-1e model is intriguing. On
Earth, O3 is distributed by the Brewer-Dobson circulation, leading to an O3
distribution that is approximately symmetric about the equator. However, in
our TRAPPIST-1e model, the combination of tidal locking and the inclusion of
the Earth-like orography generates an asymmetric circulation pattern, leading
to O3 asymmetry. Our TRAPPIST-1e model was found to exhibit a higher TOC
density than Earth with 8000 DU, nearly 28 times higher than Earth near the
South pole, and 2000 DU, which is roughly 7 times higher than Earth near the
North pole. It is interesting how tidal locking and orography affects the
global overturning circulation and O3 chemistry.

F. Sainsbury-Martinez would like to thank UK Research and Innovation for
support under grant number MR/T040726/1. G.J.C. acknowledges the studentship
funded by the Science and Technology Facilities Council of the United Kingdom
(STFC; grant number ST/T506230/1). This work was undertaken on ARC4, part of
the High Performance Computing facilities at the University of Leeds, UK.

## References

  * Barnes (2017) Barnes, R. 2017, Celestial Mechanics and Dynamical Astronomy, 129, 509
  * Bochanski et al. (2010) Bochanski, J. J., Hawley, S. L., Covey, K. R., et al. 2010, The Astronomical Journal, 139, 2679, doi: 10.1088/0004-6256/139/6/2679
  * Bourrier et al. (2017) Bourrier, V., Ehrenreich, D., Wheatley, P., et al. 2017, Astronomy & Astrophysics, 599, L3
  * Braam et al. (2023) Braam, M., Palmer, P., & Decin, L. 2023, in EGU General Assembly Conference Abstracts, EGU General Assembly Conference Abstracts, EGU–15501, doi: 10.5194/egusphere-egu23-15501
  * Brewer (1949) Brewer, A. 1949, Quarterly Journal of the Royal Meteorological Society, 75, 351
  * Carone et al. (2018) Carone, L., Keppens, R., Decin, L., & Henning, T. 2018, MNRAS, 473, 4672, doi: 10.1093/mnras/stx2732
  * Chapman (1930) Chapman, S. 1930, Mem. Roy. Meteor. Soc., 3, 103
  * Chen et al. (2019) Chen, H., Wolf, E. T., Zhan, Z., & Horton, D. E. 2019, The Astrophysical Journal, 886, 16, doi: 10.3847/1538-4357/ab4f7e
  * Conley et al. (2012) Conley, A. J., Garcia, R., Kinnison, D., et al. 2012, NCAR technical note, 3
  * Cooke et al. (2022) Cooke, G., Marsh, D., Walsh, C., Black, B., & Lamarque, J.-F. 2022, A revised lower estimate of ozone columns during Earth’s oxygenated history.  https://arxiv.org/abs/2102.11675
  * Cooke et al. (2023) Cooke, G., Marsh, D., Walsh, C., & Youngblood, A. 2023, Degenerate interpretations of O3 spectral features in exoplanet atmosphere observations due to stellar UV uncertainties: a 3D case study with TRAPPIST-1e
  * Cooke et al. (2023) Cooke, G. J., Marsh, D. R., Walsh, C., Rugheimer, S., & Villanueva, G. L. 2023, MNRAS, 518, 206, doi: 10.1093/mnras/stac2604
  * Danabasoglu et al. (2020) Danabasoglu, G., Lamarque, J.-F., Bacmeister, J., et al. 2020, Journal of Advances in Modeling Earth Systems, 12, e2019MS001916
  * Dobson (1956) Dobson, G. M. B. 1956, Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences, 236, 187
  * Dressing & Charbonneau (2015) Dressing, C. D., & Charbonneau, D. 2015, The Astrophysical Journal, 807, 45
  * Ejzak et al. (2007) Ejzak, L. M., Melott, A. L., Medvedev, M. V., & Thomas, B. C. 2007, The Astrophysical Journal, 654, 373, doi: 10.1086/509106
  * Emmons et al. (2020) Emmons, L. K., Schwantes, R. H., Orlando, J. J., et al. 2020, Journal of Advances in Modeling Earth Systems, 12.  https://api.semanticscholar.org/CorpusID:216525929
  * Gaidos et al. (2016) Gaidos, E., Mann, A., Kraus, A., & Ireland, M. 2016, Monthly Notices of the Royal Astronomical Society, 457, 2877
  * Gardner et al. (2006) Gardner, J. P., Mather, J. C., Clampin, M., et al. 2006, Space Science Reviews, 123, 485, doi: 10.1007/s11214-006-8315-7
  * Gaudi et al. (2020) Gaudi, B. S., Seager, S., Mennesson, B., et al. 2020, The Habitable Exoplanet Observatory (HabEx) Mission Concept Study Final Report.  https://arxiv.org/abs/2001.06683
  * Gettelman et al. (2019) Gettelman, A., Mills, M., Kinnison, D., et al. 2019, Journal of Geophysical Research: Atmospheres, 124, 12380
  * Gillon et al. (2016) Gillon, M., Jehin, E., Lederer, S. M., et al. 2016, Nature, 533, 221, doi: 10.1038/nature17448
  * Gillon et al. (2017) Gillon, M., Triaud, A. H. M. J., Demory, B.-O., et al. 2017, Nature, 542, 456, doi: 10.1038/nature21360
  * Gómez-Leal et al. (2019) Gómez-Leal, I., Kaltenegger, L., Lucarini, V., & Lunkeit, F. 2019, Icarus, 321, 608, doi: 10.1016/j.icarus.2018.11.019
  * Gould et al. (2003) Gould, A., Pepper, J., & DePoy, D. 2003, The Astrophysical Journal, 594, 533
  * Grimm et al. (2018) Grimm, Demory, Brice-Olivier, Gillon, Michaël, et al. 2018, A&A, 613, A68, doi: 10.1051/0004-6361/201732233
  * Hammond & Lewis (2021) Hammond, M., & Lewis, N. T. 2021, Proceedings of the National Academy of Sciences, 118, doi: 10.1073/pnas.2022705118
  * Hook (2009) Hook, I. 2009, in Science with the VLT in the ELT Era, Springer, 225–232
  * Howard et al. (2012) Howard, A. W., Marcy, G. W., Bryson, S. T., et al. 2012, The Astrophysical Journal Supplement Series, 201, 15
  * Huang (1959) Huang, S.-S. 1959, Publications of the Astronomical Society of the Pacific, 71, 421
  * Ji et al. (2023) Ji, A., Kasting, J. F., Cooke, G. J., Marsh, D. R., & Tsigaridis, K. 2023, Roy. Soc. Open Sci., 10, 230056, doi: 10.1098/rsos.230056
  * Kozakis et al. (2022) Kozakis, T., Mendonça, J. M., & Buchhave, L. A. 2022, Astronomy and Astrophysics, 665, A156, doi: 10.1051/0004-6361/202244164
  * Lim et al. (2023) Lim, O., Benneke, B., Doyon, R., et al. 2023, ApJ, 955, L22, doi: 10.3847/2041-8213/acf7c4
  * Liu et al. (2023) Liu, B., Marsh, D. R., Walsh, C., & Cooke, G. 2023, MNRAS, 524, 1491, doi: 10.1093/mnras/stad1828
  * Lustig-Yaeger et al. (2019) Lustig-Yaeger, J., Meadows, V. S., & Lincowski, A. P. 2019, The Astronomical Journal, 158, 27, doi: 10.3847/1538-3881/ab21e0
  * LUVOIR Team (2019) LUVOIR Team. 2019, The LUVOIR Mission Concept Study Final Report.  https://arxiv.org/abs/1912.06219
  * Mulders et al. (2015) Mulders, G. D., Pascucci, I., & Apai, D. 2015, The Astrophysical Journal, 814, 130
  * National Academies of Sciences, Engineering, and Medicine (2021) National Academies of Sciences, Engineering, and Medicine. 2021, Decadal survey on astronomy and astrophysics 2020, The National Academies Press Washington, DC
  * Newell (1963) Newell, R. E. 1963, Quarterly Journal of the Royal Meteorological Society, 89, 167
  * Nutzman & Charbonneau (2008) Nutzman, P., & Charbonneau, D. 2008, Publications of the Astronomical Society of the Pacific, 120, 317
  * Peacock et al. (2019) Peacock, S., Barman, T., Shkolnik, E. L., Hauschildt, P. H., & Baron, E. 2019, The Astrophysical Journal, 871, 235
  * Portmann et al. (2012) Portmann, R., Daniel, J., & Ravishankara, A. 2012, Philosophical Transactions of the Royal Society B: Biological Sciences, 367, 1256
  * Proedrou & Hocke (2016) Proedrou, E., & Hocke, K. 2016, Earth, Planets and Space, 68, 1
  * Reiners et al. (2018) Reiners, A., Zechmeister, M., Caballero, J., et al. 2018, Astronomy & Astrophysics, 612, A49
  * Remsberg (2015) Remsberg, E. 2015, Atmospheric Chemistry and Physics, 15, 3739, doi: 10.5194/acp-15-3739-2015
  * Seager (2010) Seager, S. 2010, Exoplanets
  * Sergeev et al. (2022) Sergeev, D. E., Fauchez, T. J., Turbet, M., et al. 2022, The Planetary Science Journal, 3, 212
  * Shields et al. (2016) Shields, A. L., Ballard, S., & Johnson, J. A. 2016, Physics Reports, 663, 1
  * Showman et al. (2013) Showman, A. P., Wordsworth, R. D., Merlis, T. M., & Kaspi, Y. 2013, in Comparative Climatology of Terrestrial Planets (University of Arizona Press), doi: 10.2458/azu_uapress_9780816530595-ch12
  * Sidorenko et al. (2020) Sidorenko, D., Danilov, S., Koldunov, N., Scholz, P., & Wang, Q. 2020, Geoscientific Model Development, 13, 3337, doi: 10.5194/gmd-13-3337-2020
  * Turbet et al. (2022) Turbet, M., Fauchez, T. J., Sergeev, D. E., et al. 2022, The Planetary Science Journal, 3, 211
  * Yates et al. (2020) Yates, J. S., Palmer, P. I., Manners, J., et al. 2020, Monthly Notices of the Royal Astronomical Society, 492, 1691


# High-fidelity simulation of pebble beds: Toward an improved understanding of
the wall channeling effect ††thanks: Citation: Authors. Title. Pages….
DOI:000000/11111.

David Reger, Elia Merzari, and Saya Lee
Pennsylvania State University
University Park, PA 16803
<EMAIL_ADDRESS><EMAIL_ADDRESS>
&Paolo Balestra
Idaho National Laboratory
Idaho Falls, ID 83415
<EMAIL_ADDRESS>
&Yassin Hassan
Texas A&M University
College Station, TX 77843
<EMAIL_ADDRESS>
###### Abstract

Wall channeling is a phenomena of interest for Pebble Bed Reactors (PBRs)
where flow is diverted into high-porosity regions near the wall. This
diversion of flow can have a significant impact on maximum fuel temperatures
and core bypass flow. Porous media models that are currently used to model
PBRs for design scoping and transient simulation are lacking in their
capabilities to model the wall channel effect. Recent efforts at Penn State
have produced an improved porous media pressure drop equation that is more
capable of modeling the velocity variations caused by the wall channel effect
in a porous media model. Several pebble beds were divided into concentric
rings of $0.05D_{peb}$, and average flow quantities and porosities were
extracted for the ring. A correlation between the form loss coefficient and
the local ring porosity was found, allowing for the addition of a correction
factor to the form loss term of the KTA equation. The developed correlation
was purely empirical, and thus a more thorough understanding of the underlying
flow phenomena is desired. This study investigates geometric and flow features
that can explain the observed correlation between the form coefficient and the
local porosity that was used to generate the improved pressure drop equation.
The solid surface area to volume ratio $S_{v}$ along with the production of
Turbulent Kinetic Energy (TKE) is analyzed. A relationship between $S_{v}$ and
the local porosity and an inverse relationship between the negative TKE
production and the local porosity were found, pointing to the idea that
inertial effects caused by different pore geometry in each ring contribute to
the variation of the form constant with the local porosity.

_K_ eywords Wall Channel Effect $\cdot$ Porous Media $\cdot$ Pebble Beds
$\cdot$ CFD

##  1 Introduction

The Pebble Bed Reactor (PBR) design has seen a resurgence in interest in
recent years. PBRs are currently being developed in the United States and
China is currently constructing several PBRs. As these systems are approaching
more widespread deployment, fast and accurate simulation tools are necessary
for design scoping and simulation of accident scenarios.

Porous media modeling is the current state-of-the-art for intermediate-
fidelity simulation of packed beds. Especially in the case of randomly
organized porous media like pebble beds, resolution of the complex void
regions and all fluid-solid interfaces is incredibly expensive, requiring
billions of gridpoints [1]. It is worth noting, however, that these pebble-
resolved computations are not impossible, as recent increases in computing
power have made such calculations feasible [2]. Regardless, the immense cost
involved with pebble-resolved simulations makes them unrealistically expensive
for design scoping and plant-level simulation. Porous media models mend this
issue by homogenizing the porous media with spatial averaging, reducing the
computational cost by several orders of magnitude. Closure models are also
required to approximate the effects of small flow features on the macroscale
behavior of the flow. These closure models provide estimations of drag
coefficients, effective conductivities, and interphase heat transfer
coefficients among other parameters [3]. The accuracy of the closure models
will directly influence the accuracy of the porous media model as a whole, and
thus it is important to ensure that the models effectively represent real-
world physics.

One area where current closure models are insufficient is in the near-wall
region of PBRs. In this region, the presence of the wall influences the
packing of the pebbles, causing them to form more orderly structures and
increasing the porosity. This effect can be seen in the projection of the
pebble centers in a PBR found in Figure 1. Many existing closure models have
been developed to predict the average behavior of the pebble bed as a whole,
leading to inaccuracies when applied to the near-wall region that differs
greatly from the bed interior. Accurate prediction of the flow in the near-
wall region is critical, as it can have a significant impact on bypass flow
where coolant is diverted through gaps between reflector blocks. Additionally,
depending on the specific PBR design, there may be a peak in the fuel power
near the wall where neutrons are reflected and thermalized by the graphite
reflector. Accurate modeling of the near-wall region is therefore critical to
ensure accurate prediction of fuel temperature maximums and provide confidence
in predicted safety margins. For these reasons, improving the understanding
and modeling capabilities of the near-wall region has been identified by the
United States Nuclear Regulatory Commission (NRC) as an issue of high
importance [4].

Figure 1: (left) Instantaneous velocity field for a bed of 7,000 pebbles
($D/d_{peb}=30$). (right) Projection of the 7,000 pebble centers in a PBR onto
an axial plane. The wall-channeling effect is visible in the organized ring of
pebbles near the wall.

The near-wall region of pebble beds has been a topic of research interest for
many years, with many researchers looking to improve understanding of this
region to enhance modeling of the near-wall porosity, pressure drop, and heat
transfer. The phenomena has been studied experimentally by Amini [5]. Their
study used hot wire anemometry probes to measure the near-wall flow velocities
in several differently-shaped near-wall gaps. Their results experimentally
confirmed the existence of the high-velocity flow channels that form in the
gaps near the wall. They then examined the different behaviors observed in the
two different near-wall gap geometries. Nguyen has also experimentally studied
the flow behavior in the wall region through the use of Particle Image
Velocimetry [6]. They examined the cross flows between pebbles and the bypass
flow in the near-wall region through the use of Proper Orthogonal
Decomposition (POD). Computationally, the near-wall region has been studied by
Fick and Merzari [7]. Their study employed a regular arrangement of pebbles
with a confining wall. They performed a Direct Numerical Simulation (DNS) and
examined the 2nd and 3rd order flow statistics near the wall. Their study
began to reveal some of the defining flow characteristics present in the near-
wall region. Modeling of the near-wall region has been improved through the
use of porosity, pressure drop, and heat transfer equations. De Klerk
developed a porosity correlation that is capable of accurately modeling the
oscillatory porosity variation near the wall [8]. Reichelt was one of the
first to investigate the effects of wall-channeling on the bed pressure drop.
He found significant errors in the Ergun equation, one of the widely used
pressure drop equations at the time, when applied to slender beds where the
near-wall effects make up a large portion of the bed [9]. He then suggested a
new correlation that accounts for the ratio between pebble diameter and bed
diameter. Eisefeld and Schnitzlein [10] identified the effects that wall-
channeling may have on the pressure drop by comparing experimental results to
correlations available at the time. They identified the approach by Reichelt
as being promising in effectively modeling the influence of the near-wall
region and developed their own improved correlation based on the findings from
their study. With regards to heat transfer, work by Achenbach has studied the
near-wall region and developed a correlation to model the near-wall heat
transfer coefficient in packed beds [11].

This work reviews the current findings by our group at Penn State with regards
to near-wall flow behavior in a pebble bed reactor [12, 13]. Additional
investigation into the local geometry and Turbulent Kinetic Energy (TKE)
production is then presented to better explain the observations that have been
made thus far. Section 2 details the methods used to study the near-wall flow,
and results and analysis are then presented in Section 3.

##  2 Codes and Methods

###  2.1 Introduction to NekRS

Argonne National Laboratory’s spectral-element CFD code NekRS [14] was chosen
as the high-fidelity code for this study. NekRS is a GPU-oriented variant of
the well-established open-source code Nek5000 [15]. It demonstrates excellent
scalability [16], and is capable of linking to Nek5000 to utilize its existing
pre- and post-processing features.

The simulations performed in this work use the incompressible, constant-
properties Navier-Stokes equations in dimensionless form:

$\frac{\partial{\vec{v_{i}}}}{\partial{t}}+{\vec{v_{i}}}\cdot\nabla{\vec{v_{i}}}=-\nabla{P}+\frac{1}{Re}\nabla^{2}{\vec{v_{i}}}$
(1) $\nabla\cdot{\vec{v_{i}}}=0$ (2)

where $v$ is the fluid velocity, $P$ is the fluid pressure, and $Re$ is the
Reynolds number based on the pebble diameter and inlet velocity
($\frac{\rho{v_{inlet}}{D_{peb}}}{\mu}$). The nondimensionalization of
variables follows the following scheme:

$x^{*}=\frac{x}{D_{peb}}$ (3) $v^{*}=\frac{v}{v_{inlet}}$ (4)
$t^{*}=\frac{{t}{v_{inlet}}}{D_{peb}}$ (5)
$P^{*}=\frac{P}{\rho{v_{inlet}}^{2}}$ (6)

where $D_{peb}$ is the pebble diameter, $v_{inlet}$ is the inlet velocity, $t$
is the time, and $\rho$ is the fluid density. The simulations performed in
NekRS are wall-resolved LES, where an explicit filter was used to approximate
the effects of dissipation on the subgrid scales [17]. Simulations in NekRS
were run with a polynomial order of 7.

###  2.2 Mesh Creation

The pebble beds created for high-fidelity simulation in this work were
generated using the Discrete Element Method (DEM) in the open-source code
Project Chrono [18].

The DEM simulations used to generate the beds for this work used the material
properties of graphite, found in Table 1 [19, 20, 21]. Additional information
about the contact model used in Project Chrono can be found in [22].

Table 1: Graphite material properties used in the DEM simulation. Properties were used for both sphere-sphere and sphere-wall contacts. [19, 20, 21] Property | Value
---|---
Density | 2260 $\mathrm{kg/m^{3}}$
Elastic Modulus | 8 GPa
Poisson’s Ratio | 0.12
Coefficient of Restitution | 0.6
Sliding Friction Coefficient | 0.3
Rolling Friction Coefficient | 0.1
Simulation Timestep |  $5\times 10^{-5}$ s

The beds were generated by randomly sampling sheets of pebbles separated by
$2D_{peb}$ of axial distance. These sheets were then dropped into a
cylindrical vessel to randomly pack. Several thousand pebbles were used for
the packing, and then a section of pebbles was extracted from the center of
the resulting bed to avoid any influence of the cylinder bottom or top.

A meshing script is then used to generate an all-hexahedral mesh for
simulation with NekRS. Developed as part of the Cardinal multiphysics project
[2], the script uses a novel Voronoi-cell approach. It receives the pebble
center coordinates along with the pebble and cylinder diameters as inputs and
begins by generating a voronoi cell for the void region around each pebble.
The faces of each voronoi cell are equidistant between pebbles. At the top and
bottom of the bed, additional pebble centers are provided to determine the top
and bottom voronoi faces, although these pebbles are not included in the
resulting mesh. The voronoi faces are then modified to improve the resulting
mesh quality by collapsing small edges or dividing very long edges. Quad
elements are then generated on the faces and are projected down onto the
pebble surfaces, producing an all-hexahedral mesh. Pebble-pebble and pebble-
wall contacts are handled by adding a small chamfer at the point of contact.
This simply inserts a small cylinder through the contact point to widen it
slightly. This method has been shown to have a minimal effect on the resulting
porosity and pressure drop compared to other methods, such as shrinking each
pebble to avoid contacts [23]. An example of the resulting mesh is pictured in
Figure 2.

We note that we designed the mesh to have resolution sufficient to resolve the
Taylor micro-scales based on estimates from a previous work [24] and
additional calculations computed for this work with RANS. We also design the
mesh to resolve the boundary layer and have the first grid point at $y^{+}<1$
and a sufficient number of points in the viscous sub-layer.

A flat-profile inlet velocity condition and a stabilized outflow condition
[25] are then applied along with no-slip conditions on the pebble and cylinder
walls.

Figure 2: Example of the high-fidelity model meshes used for this work showing
the chamfer between pebbles.

###  2.3 Turbulent Kinetic Energy Budgets

The Turbulent Kinetic Energy (TKE) represents the mean kinetic energy that is
carried by the eddies found in turbulent flow. TKE is typically produced
through shear, forcing, or friction. It is then transferred donw the energy
cascade where it is eventually dissipated by viscous forces in the smallest
eddies. The TKE can be represented as the sum of the individual contributing
mechanisms to produce the TKE equation:

$\underbrace{\frac{\partial{k}}{\partial{t}}}_{\text{TKE
Derivative}}+\underbrace{\overline{u_{j}}\frac{\partial{k}}{\partial{x_{j}}}}_{\text{Advection}}=\underbrace{-\frac{1}{\rho}\frac{\partial\overline{u^{\prime}_{i}{p^{\prime}}}}{1}}_{\text{Pressure
Diffusion}}-\underbrace{\frac{1}{2}\frac{\partial\overline{{u^{\prime}_{j}}{u^{\prime}_{j}}{u^{\prime}_{i}}}}{\partial{x_{i}}}}_{\text{Turbulent
Transport}}+\underbrace{\nu\frac{\partial^{2}{k}}{\partial{x}^{2}_{j}}}_{\text{Molecular
Transport}}-\underbrace{\overline{{u^{\prime}_{i}}{u^{\prime}_{j}}}\frac{\partial\overline{u_{i}}}{\partial{x_{j}}}}_{\text{Production}}-\underbrace{\nu\overline{\frac{\partial{u^{\prime}_{i}}}{\partial{x_{j}}}\frac{\partial{u^{\prime}_{i}}}{\partial{x_{j}}}}}_{\text{Dissipation}}$
(7)

These terms are commonly referred to as the budgets of the TKE. Analysis of
the production term for a packed bed is carried out in section 3.

###  2.4 Validation of NekRS

Validation of NekRS’s ability to reproduce the velocity and pressure drop in a
packed bed has been performed with experimental data from Texas A&M
University. Velocity validation was performed on a bed of 146 pebbles in work
by Yildiz [24]. Additional validation has been performed for a bed of 67
pebbles, a comparison of velocity profiles between experiment and simulation
can be seen in Figure 3. Validation of the pressure drop was performed by
comparing results between experimental and simulation pressure gradients over
five Reynolds numbers for a second bed of 789 pebbles. This comparison is
shown in Figure 4, where it can be seen that NekRS falls within the range of
experimental uncertainty at all points of comparison.

Figure 3: In-Plane velocity magnitude over one sampling line for the 67-pebble
validation study. The shaded gray region indicates the location of a pebble.
Retrieved from [26] Figure 4: Comparison of pressure gradients between
experiment and corresponding simulation.

##  3 Results

###  3.1 Near-Wall Form Coefficient Trends

Current work by our group has investigated the effect of the cylinder wall on
the form loss coefficients [12, 13]. An understanding of this effect is
critical to properly model the radial variation in the streamwise velocity of
a PBR. The goal of this study was to improve the capabilities of the KTA drag
correlation[27] to more accurately model the near-wall velocity profile in a
porous media model. NekRS simulation was used to produce an LES flow dataset
for two pebble beds of roughly 1,568 and 1,700 pebbles at aspect ratios
($D_{bed}/D_{peb}$) of 13 and 14. This analysis was then performed by
separating the beds into concentric rings of $0.05D_{peb}$ width. The ring-
volume-average porosity, fluid velocity, and pressure drop could then be
calculated for each ring. This information can then be used to calculated the
form loss coefficient in each ring. An overview of the data extraction method
is shown in Figure 5 with more information available in Reference [13].

Figure 5: Data extraction method used to generate the improved pressure drop
equation from [13].

The form loss terms were then plotted against the Reynolds number in each
ring, and the constant in the numerator of the KTA form term was calculated
for each ring. These constants were then plotted against the ring porosity,
shown in Figure 6. It was determined that there is a correlation between the
form constant and the ring porosity. A fourth-order polynomial was fit to the
data to describe this correlation:

$f(\epsilon)=253.9\epsilon^{4}-499.3\epsilon^{3}+364.7\epsilon^{2}-115.6\epsilon+14.21$
(8)

where $\epsilon$ is the local ring porosity. This correction term may be added
to the KTA equation to produce an improved drag correlation:

$\frac{\Delta{P}}{L}=\left(\frac{320}{{Re}_{m}}+\frac{6{f(\epsilon)}}{{{Re}_{m}}^{0.1}}\right)\left(\frac{1-\epsilon}{\epsilon^{3}}\right)\left(\frac{\rho^{2}{v_{s}}^{2}}{{D}_{p}}\right)\left(\frac{1}{2\rho}\right)$
(9)

It should be noted that this correlation can currently be stated to be valid
for $10^{2}<Re_{m}<10^{4}$ and $0.2<\epsilon<0.9$. Additional data may help to
further quantify the uncertainties of this correlation at the high and low
limits of the porosity where the data is currently sparse.

Figure 6: Form loss term constants plotted against the porosity of each
respective ring. A fourth order fit is applied to the data.

After the form analysis was performed and a new correlation was determined,
several knowledge gaps still persisted. Most importantly, it was not
immediately clear as to why this correlation between the region porosity and
the form loss existed. Similarly, it was also unclear as to whether this
variation is truly a correlation with porosity, or rather with wall distance.
The next section presents some preliminary analysis into possible explanations
for the dramatic increase in the form loss constant as the ring porosity
increases.

###  3.2 Toward an improved understanding of the wall-channeling effect

Among previously derived equations for the pressure loss in a packed bed of
spheres, the pebble surface to volume ratio $S_{v}$ is commonly used. Carman
and Kozeny [28] used this definition in their derivation of the viscous and
inertial losses in a packed bed. It was also used by Ergun [29] who built off
of the work of Carman and Kozeny. In these previous derivations, $S_{v}$ of a
sphere is used, which can easily be calculated as $6/D_{p}$. The actual
surface to volume ratio was calculated for concentric rings of $0.05D_{peb}$
width for several computational pebble beds of different sizes. These values
were then plotted against the porosity in each ring, shown if Figure 7. From
this plot, it can be seen that $S_{v}$ is roughly equal to $6/D_{peb}$ in low
and medium porosity regions. The average porosity of most beds falls around
0.35-0.5 depending on the bed aspect ratio, meaning that the $6/{D_{p}}$ value
that many correlations have assumed for $S_{v}$ is fairly accurate when
applied to an averaged bed. As the porosity increases, however, this value
increases significantly. Also shown in Figure 7 is the trend of the form
constant $6f(\epsilon)$ for comparison. The trend in $S_{v}$ is significantly
more linear than that of $f(\epsilon)$, although this increase is the surface
to volume ratio may be a contributor to the additional dependency of the form
constant on the porosity.

Pressure drop equations such as the Ergun, Carman-Kozeny, KTA, and our
improved KTA equation separate the porosity effects from the viscous and
inertial loss coefficients. This has previously allowed for the representation
of the losses with a quantity that represents the losses per pore and a
porosity-dependent quantity that represents the number of pores. The
multiplication of these two terms then yields a drag coefficient. The result
from Figure 7 suggests that this previous approach may not be accurate, as the
loss per pore changes as a function of the pore shape (which is itself a
function of the porosity). This result suggests higher losses per pore as the
porosity is increased. This at first may seem counter intuitive, as one would
expect a more open and uniform flow geometry to have lower pressure loss. This
intuition remains correct even with the findings from Figure 7, as although
the losses per pore are suggested to be higher at high porosities, the total
loss coefficients will still remain lower than a lower-porosity region, as the
$\frac{1-\epsilon}{\epsilon^{3}}$ term that represents the number of pores
will be small for high porosity regions.

Figure 7: Ratio between pebble surface area and pebble volume $S_{v}$
calculated for rings of $0.05D_{peb}$ width versus ring porosity for several
pebble beds. The form constant term with correction factor, $6f(\epsilon)$ is
also shown. All beds are nondimensionalized to $D_{peb}=1$ with the bed aspect
ratio shown in parentheses in the legend.

####  3.2.1 Investigation of the Turbulent Kinetic Energy production

A better understanding of the near-wall flow physics can be obtained by
examining the budgets of the TKE in the near-wall region. The budget terms
were calculated from a Direct Numerical Simulation (DNS) of a small bed of 67
pebbles at Re = 1,460 that was previously used for experimental velocity
validation in section 2.4. Previous work examined line samples of the fluid
velocity at three locations [26], one of which is shown in Figure 3. For the
DNS in this work, the bed was simulated at a polynomial order of 9 for 50
convective units to converge the 3rd-order statistics. The bed is shown
alongside a centerplane slice of the TKE in Figure 8. Some of the notable
features are the peaks in the TKE, typically in the wake regions behind the
pebbles. There is also a noticeable left-side bias in the TKE at the bed
outlet. This is can be attributed to the distribution of the top layer of
pebbles, causing a turbulent plume in this region. Analysis of the TKE
production term, as was done in previous work [7], can reveal additional
information on the inertial effects throughout the bed.

Figure 8: 67 pebble bed used for the DNS simulation (left) and centerplane
slice of the TKE (right).

Figure 9 shows the isosurfaces of highly positive and negative TKE production
in the bed. The areas of negative production are particularly interesting. In
many canonical cases, such as channel flow, the production term is positive
and acts as a source of TKE. A complex case such as a pebble bed, however, can
see this term change sign and become negative, suppressing the TKE rather than
strengthening it. There are large areas of negative production on the bottom
surfaces of the pebbles. Flow accelerates around the pebble in these regions,
causing negative production and suppression of the TKE. Additionally, there
are some elongated negative-production structures in the large void regions
near the wall. The areas of high TKE production exist in the wake regions
behind pebbles where the velocity gradients and covariances are high.

Figure 9: 3D Visualization of isosurfaces enclosing areas where the TKE
production is less than or equal to -5 and greater than or equal to 30.

Further investigation into the TKE production was also performed to determine
if there is any relationship to the observed form correction factor. It was
theorized that differences in the amount of negative production could
potentially explain the correlation between the form constant and the porosity
that was observed previously. As has been done with the other investigations,
the 67-pebble bed was separated into rings of $0.05D_{peb}$ width, and the
average negative TKE production ($\langle|P|\rangle$ for $P<0$) was calculated
for each ring. An inverse relationship between the negative TKE production and
the the porosity was found that closely matches the relationship between the
form constant and the porosity. Figure 10 presents this relationship, along
with the correction factor $f(\epsilon)$ for reference. It can be seen that
the trends of $1/\langle|P|\rangle>$ with the porosity exhibit many
similarities to the trends of $f(\epsilon)$, as there are increases in both
the high and lower porosity regions with a slight linear increase in the
medium-porosity regions. This points to the idea that inertial flow effects
vary greatly based on the porosity as a result of the different void
geometries. These inertial effects will have an effect on the form loss
coefficient which can perhaps explain the trends described with the
$f(\epsilon)$ correction factor. Regions with high negative TKE production
experience more laminarization and a lower form loss per pore. Meanwhile,
areas with lower negative TKE production exhibit the opposite behavior, with
less laminarization and higher form loss per pore, aligning with the observed
trend in the form constant found in Figure 6. The data point nearest to the
wall in Figure 10, however, still requires additional investigation, as it
does not follow the trend seen with the rest of the data. The bed used to
extract this data is fairly small, leading to small averaging volumes and poor
averaging statistics. Although the results presented from the investigation of
the negative production suggest a possible explanation for the trend in the
form constant, it remains difficult to decisively confirm this explanation.
Additional data for larger beds to gather additional data points is necessary
to further reinforce this claim.

Figure 10: Average of negative production ($1/\langle|P|\rangle$ for $P<0$)
versus local ring porosity. The correction factor $f(\epsilon)$ is also
included on a separate axis for reference.

##  4 Conclusion

This study presented a review of current efforts at Penn State University to
gain a better understanding of the near-wall region of pebble bed reactors.
Work has thus far studied the flow in a PBR by dividing the bed into multiple
concentric regions and examining temporally and spatially-averaged flow
characteristics in each region. Examination of the form-loss coefficients has
revealed a correlation between the form loss term and the porosity of each
ring region. This information has been used to generate an improved pressure
drop correlation that is capable of more accurately reproducing the radial
velocity profile in a PBR with a porous media code [13].

Although this improved correlation has proven to be a promising result,
additional work is currently being conducted to better understand the near-
wall flow phenomena that may be behind the observed correlation between the
form term and the porosity. This study investigated the geometry of the
concentric rings, particularly the ratio between the solid surface area and
the solid volume. This ratio (denoted $S_{v}$) has been used as part of the
derivation of many past pressure drop equations, where it has been assumed to
be the theoretical value of $6/D_{p}$. Calculation of $S_{v}$ for several
computational beds revealed that $S_{v}$ actually varies with the porosity,
ranging from $3/D_{p}$ to greater than $30/D_{p}$ based on the porosity. This
helps to explain the increase in the form constant term with the porosity,
although the relationship between $S_{v}$ and $\epsilon$ is fairly linear.
This indicates that although $S_{v}$ may play a role in influencing the
increasing form constant, it is likely not the only contributor.

The TKE production in each ring was then investigated to further examine
potential flow phenomena that may contribute to the varying form constant. A
relationship was found between the inverse of the negative TKE production and
the porosity that closely matches the trend observed in the form constant.
This suggests that the observed trend in the form constant with the porosity
may be related to inertial effects that are caused by the different pore
geometries at different porosities. The bed used to generate the TKE budgets
was rather small, at only 67 pebbles, and thus further data for larger beds is
needed to reinforce the hypothesis on the effect of the negative production.

## ACKNOWLEDGMENTS

This material is based upon work supported under an Integrated University
Program Graduate Fellowship.

This research used resources of the Oak Ridge Leadership Computing Facility at
the Oak Ridge National Laboratory, which is supported by the Office of Science
of the U.S. Department of Energy under Contract No. DE-AC05-00OR22725.

This research made use of Idaho National Laboratory computing resources which
are supported by the Office of Nuclear Energy of the U.S. Department of Energy
and the Nuclear Science User Facilities under Contract No. DE-AC07-05ID14517.

## References

  * [1] C. Wu, Y. Feng, C. Chieng, and C. Liu.  Investigating the advantages and disadvantages of realistic approach and porous approach for closely packed pebbles in cfd simulation.  Nuclear Engineering and Design, 240:1151–1159, 2010.
  * [2] E. Merzari, H. Yuan, M. Min, D. Shaver, R. Rahaman, P. Shriwise, P. Romano, A. Talamo, Y. Lan, D. Gaston, R. Martineau, P. Fischer, and Y. A. Hassan.  Cardinal: A lower length-scale multiphysics simulator for pebble-bed reactors.  Nuclear Technology, 2021.
  * [3] A. J. Novak, R. W. Carlsen, S. Schunert, P. Balestra, D. Reger, R. N. Slaybaugh, and R. C. Martineau.  Pronghorn: A multidimensional coarse-mesh application for advanced reactor thermal hydraulics.  Nuclear Technology, 207(7):1015–1046, 2021.
  * [4] S.J. Ball.  Next Generation Nuclear Plant Phenomena Identification and Ranking Tables (PIRTs).  Technical report, NRC, 2007.
  * [5] N. Amini and Y. Hassan.  Experimental study of bypass flow in near wall gaps of a pebble bed reactor using hot wire anemometry technique.  Annals of Nuclear Energy, 65:60–71, 2014.
  * [6] T. Nguyen, R. Muyshondt, and Y. Hassan.  Experimental investigation of cross flow mixing in a randomly packed bed and streamwise vortex characteristics using particle image velocimetry and proper orthogonal decomposition analysis.  Physics of Fluids, 31(2), 2019.
  * [7] L. Fick and E. Merzari.  Inertial effects and anisotropy for the flow in a domain of close packed spheres with a bounding wall.  Nuclear Technology, 208:539–561, 2020.
  * [8] A. de Klerk.  Voidage variation in packed beds at small column to particle diameter ratio.  AIChE, 49(8):2022–2029, 2003.
  * [9] W. Reichelt.  Zur berechnung des druckverlustes einphasig durch- stromter kugel- und zylinderschuttungen.  Rheological Acta, 17:676–692, 1978.
  * [10] B. Eisfeld and K. Schnitzlein.  The influence of confining walls on the pressure drop in packed beds.  Chemical Engineering Science, 56:4321–4329, 2001.
  * [11] W. Hahn and E. Achenbach.  Determination of hte wall heat transfer coefficient in pebble beds.  Kernforschungsanlage Juelich, 1986.
  * [12] D. Reger, E. Merzari, P. Balestra, S. Schunert, Y. Hassan, H. Yuan, P. Fischer, M. Min, and Y. Lan.  Pressure drop correlation improvement for the near-wall region of pebble-bed reactors.  Nuclear Technology, 2022.
  * [13] D. Reger, E. Merzari, P. Balestra, S. Schunert, Y.A. Hassan, and H. Yuan.  An improved pressure drop correlation for modeling localized effects in a pebble bed reactor.  Nuclear Engineering and Designed, 403, 2023.
  * [14] P. Fischer, S. Kerkemeier, M. Min, Y. H. Lan, M. Phillips, T. Rathnayake, E. Merzari, A Tomboulides, A. Karakus, N Chalmers, and T. Warburton.  Nekrs, a gpu-accelerated spectral element navier-stokes solver.  2021\.
  * [15] Paul Fischer, J Lottes, S Kerkemeier, O Marin, K Heisey, E Obabko, E Merzari, and Y Peet.  Nek5000 user documentation.  Argonne National Laboratory, Lemont, IL, Report No. ANL/MCS-TM-351, 2016.
  * [16] Elia Merzari, Paul Fischer, Misun Min, Stefan Kerkemeier, Aleksandr Obabko, Dillon Shaver, Haomin Yuan, Yiqi Yu, Javier Martinez, Landon Brockmeyer, et al.  Toward exascale: overview of large eddy simulations and direct numerical simulations of nuclear reactor flows with the spectral element method in nek5000.  Nuclear Technology, 206(9):1308–1324, 2020.
  * [17] Saleh Rezaeiravesh, Ricardo Vinuesa, and Philipp Schlatter.  On numerical uncertainties in scale-resolving simulations of canonical wall turbulence.  Computers & Fluids, page 105024, 2021.
  * [18] A. Tasora, R. Serban, H. Mazhar, A. Pazouki, D. Melanz, J. Fleischmann, M. Taylor, H. Sugiyama, and D. Negrut.  Chrono: An open source multi-physics dynamics engine.  High Performance Computing in Science and Engineering, 2015.
  * [19] A. Toptan, W. Jiang, S. Novascone, and J. Hales.  Matrix graphite material models in pebbles and compacts for bison.  Technical report, Idaho National Laboratory, 2021.
  * [20] L. Xiaowei, L. Xiaotian, and Y. Suyuan.  Nuclear graphite friction properties and the influence of friction properties on the pebble bed.  Nuclear Engineering and Design, 240, 2010.
  * [21] H. Qi, N. Gui, X. Yang, J. Tu, and S. Jiang.  Effects of restitution coefficient on pebble motion in a thin pebble bed - lagrangian analysis.  Annals of Nuclear Energy, 145, 2020.
  * [22] Jonathan Fleischmann.  Dem-pm contact model with multi-step tangential contact displacement history, 2015.
  * [23] Y. Lan, P. Fischer, E. Merzari, and M. Min.  All-hex meshing strategy for densely packed spheres.  2016\.
  * [24] M. A. Yildiz, G. Botha, H. Yuan, E. Merzari, R. Kurwitz, and Y. A. Hassan.  Direct numerical simulation of the flow through a randomly packed pebble bed.  Journal of Fluids Engineering, 142(4), 2020.
  * [25] Suchuan Dong, George E Karniadakis, and Chryssostomos Chryssostomidis.  A robust and accurate outflow boundary condition for incompressible flow simulations on severely-truncated unbounded domains.  Journal of Computational Physics, 261:83–105, 2014.
  * [26] D. Reger, E. Merzari, P. Balestra, S. Schunert, Y. Hassan, and S. King.  Direct numerical simulation and large eddy simulation of a 67–pebble bed experiment.  Nuclear Technology, 2022.
  * [27] K. Ausschusses and KTA.  Reactor core design of high-temperature gas-cooled reactors - part 3: Loss pressure through friction pebble bed cores.  1981\.
  * [28] P. Carman.  Fluid flow through granular beds.  Transactions of the Institution of Chemical Engineering, 75, 1937\.
  * [29] S. Ergun.  Fluid flow through packed columns.  Journal of Chemical Engineering Progress, 48(2):89–94, 1952.


# Complex homothetic sections and projections
through a Helly type Theorem for cosets of $\mathbb{S}^{1}$

Jorge Luis Arocha  ,  Javier Bracho  and  Luis Montejano  UNAM, Mexico
<EMAIL_ADDRESS>
###### Abstract.

We prove that two closed subsets of complex space $\mathbb{C}^{n}$ with
corresponding complex homothetic sections (projections) are complex
homothetic. The proof uses a new Helly-type theorem for cosets of closed
subgroups of $\mathbb{S}^{1}$.

##  1\. Introduction

Let $V$ be an $n$-dimensional vector space and let $K_{1},K_{2}$ be two closed
sets, $n\geq 3.$ Suppose that for every 2-dimensional subspace (hyperplane)
$\Lambda$ the corresponding sections

$\Lambda\cap K_{1}\mbox{ and }\Lambda\cap K_{2}$

are “equivalent”. Is it true that

$K_{1}\mbox{ is ``equivalent'' to }K_{2}?$

Where the meaning of “equivalent” will be clarified in the sequel.

Suppose for example that $K_{1},K_{2}$ are convex bodies containing the origin
in euclidean $n$-space, $n\geq 3$, and the meaning of “equivalent” is
congruence. That is, suppose that for every hyperplane $\Lambda$ through the
origin there is a linear isometry $g\in O(\Lambda)$ and a vector $v\in\Lambda$
(depending on $\Lambda$) such that

$g(\Lambda\cap K_{1})+v=\Lambda\cap K_{2}.$

Is it true that $K_{1}$ is either a translated copy of $K_{2}$ or $-K_{2}$?

V.P. Golubyatnikov [3] and D. Ryagobin [4], [5], [6], [1], have studied this
problem extensively. Perhaps the most representative result in this direction
is Ryagobin’s Theorem in [4]. Suppose $K_{1},K_{2}$ are convex bodies in
euclidean $n$-space such that for every plane $\Lambda$ through the origin
there is a rotation $r$ of $\Lambda$ about the origin (depending on $\Lambda$)
for which $r(\Lambda\cap K_{1})=(\Lambda\cap K_{2})$, then either
$K_{1}=K_{2}$ or $K_{1}=-K_{2}$.

Ryagobin’s proof rests heavily on the fact that there are only a countable
number of closed subgroups of $SO(2)$. The analogous is not true for $SO(k)$,
$k>2$, and this is a reason why the result is unknown for $k$-subspaces in
$\mathbb{R}^{n}$.

In the charming paper [7], Rogers proved the following. Suppose $K_{1},K_{2}$
are convex bodies in euclidean $n$-space, $n\geq 3$, such that for every
hyperplane $\Lambda$ through the origin

$\Lambda\cap K_{1}\mbox{ is positively homothetic to }\Lambda\cap K_{2},$

then $K_{1}$ is positively homothetic to $K_{2}.$

Note that the center of each homothecy between $\Lambda\cap K_{1}$ and
$\Lambda\cap K_{2}$ is not fixed and of course depends on the hyperplane
$\Lambda$. This result was crucial in many characterizations of ellipsoids in
the second half of the 20th century.

Studying complex ellipsoids in [2], it was necessary to prove a complex
version of Rogers’ Theorem for real homotheties. More precisely, suppose
$K_{1},K_{2}$ are convex bodies in complex $n$-space, $\mathbb{C}^{n}$, $n\geq
3$, such that for every complex hyperplane $\Lambda$ through the origin

$\Lambda\cap K_{1}\mbox{ is positively homothetic to }\Lambda\cap K_{2},$

then $K_{1}$ is positively homothetic to $K_{2}.$

It would be interesting to know if the following complex version of Ryagobin’s
Theorem is true. Suppose two convex bodies $K_{1},K_{2}$ in complex $n$-space
$\mathbb{C}^{n}$ have the property that for every complex plane $\Lambda$,
there is $g\in U(\Lambda)$ such that

$g(\Lambda\cap K_{1})=\Lambda\cap K_{2}\,.$

Is there $z\in\mathbb{S}^{1}$ with $|z|=1$ such that $zK_{1}=K_{2}$?

It is precisely in this direction that we obtain the following theorem as a
special case of the above. We note that our result is true not only for convex
subsets but also for closed subsets of $\mathbb{C}^{n}$ (including of course
finite sets).

######  Theorem 1.

Let $K_{1},K_{2}$ be two closed subsets of $\mathbb{C}^{n}$, $n\geq 3$.
Suppose that for every complex plane $\Lambda$ through the origin there is a
complex number $z_{\Lambda}\in\mathbb{C}$ such that

$z_{\Lambda}(\Lambda\cap K_{1})=\Lambda\cap K_{2}.$

Then, there is $z\in\mathbb{C}$, such that

$zK_{1}=K_{2}.$

The analogue of Theorem 1 in the real setting is essentially Lemma 1 of [4].

The proof of Theorem 1 makes use of the following Helly-type theorem for
cosets of closed subgroups of $\mathbb{S}^{1}$.

######  Theorem 2.

Let $\mathcal{F}$ be a pairwise intersecting family of cosets of closed
subgroups of $\mathbb{S}^{1}$. Then $\mathcal{F}$ has non-empty intersection.

It would be interesting to know which groups have this Helly property. That
is, for which groups $G$ there exists an integer $m>2$ with the property that
if $\mathcal{F}$ is a family of m+1 cosets of closed subgroups of $G$ and each
$m$ members de $\mathcal{F}$ is intersecting, then $\mathcal{F}$ is
intersecting. Indeed, the integers $\mathbb{Z}$ have this property. We will
prove that every three pairwise intersecting arithmetic sequences is
intersecting (see Theorem 4).

Similarly, using Theorem 2, we prove the corresponding result for projections.
Let $\Lambda$ be a complex subspace of $\mathbb{C}^{n}$, we denote by
$\pi_{\Lambda}:\mathbb{C}^{n}\to\Lambda$, the orthogonal projection.

######  Theorem 3.

Let $K_{1},K_{2}\subset\mathbb{C}^{n}$ be either two closed convex sets or two
finite sets, $n\geq 3$. Suppose that for every complex plane $\Lambda$ through
the origin there is a complex number $z_{\Lambda}\in\mathbb{C}$ such that

$z_{\Lambda}(\pi_{\Lambda}(K_{1}))=\pi_{\Lambda}(K_{2}).$

Then, there is $z\in\mathbb{C}$, such that

$zK_{1}=K_{2}.$

Regarding projections of convex bodies in $\mathbb{R}^{5}$ whose projections
are $SU(2)$-congruent, Ryagobin obtained interesting results in [6].

##  2\. The group of complex homotheties

From now on, unless otherwise stated, all our lines, planes, hyperplanes pass
through the origin.

Let $V$ be a complex vector space and let $K\subset V$ be a closed subset.
Define the group of homotheties of  $K$ as:

$H(K)=\\{\,z\in\mathbb{S}^{1}\,|\,zK=K\,\\},$

which is a closed subgroup of $\mathbb{S}^{1}$. For every subspace
$\Delta\subset\mathbb{C}^{n}$,

(1) $H(K)\mbox{ is a subgroup of }H(\Delta\cap K),$ (2) $H(K)\mbox{ is a
subgroup of }H(\pi_{\Delta}(K)).$

Note that (1) is trivially true. On the other hand, if
$\pi_{\Delta}:\mathbb{C}^{n}\to\Delta$ is the orthogonal projection, then
$z\pi_{\Delta}(x)=\pi_{\Delta}(zx)$, for every $z\in\mathbb{C}$ and every
$x\in\mathbb{C}^{n}$. Hence, $zK=K$ implies $z\pi_{\Delta}K=\pi_{\Delta}K$ and
consequently (2) is true.

Now, if $K_{1},K_{2}$ are closed subsets of $V$, let us define

$H(K_{1},K_{2})=\\{z\in\mathbb{S}^{1}|zK_{1}=K_{2}\\}.$

For every $z_{0}\in H(K_{1},K_{2})$,

(3) $z_{0}H(K_{1})=H(K_{1},K_{2})\,.$

Therefore, if $H(K_{1},K_{2})$ is nonempty, it is a coset of the closed
subgroup $H(K_{1})$ of $\mathbb{S}^{1}$.

Let $\mathcal{L}$ denote the set of complex lines through the origin in
$\mathbb{C}^{n}$. Of course, it is the complex projective space of dimensión
$n-1$, but we need to refer to its elements as complex lines.

######  Lemma 1.

Let $K\subset\mathbb{C}^{n}$, $n\geq 3$, be a closed subset. Then

$H(K)=\bigcap_{L\in\mathcal{L}}H(L\cap K).$

###### Proof.

By (1), $H(K)\subset\bigcap_{L\in\mathcal{L}}H(L\cap K)$. Suppose now,
$z\in\bigcap_{L\in\mathcal{L}}H(L\cap K)$, that is, $z(L\cap K)=L\cap K$, for
every line $L\in\mathcal{L}$. Consequently, $zK=K$, which implies that $z\in
H(K)$ as we wished. ∎

##  3\. Helly type results

For the proof of our main Helly type theorem (Theorem 2), we need two lemmas.

######  Lemma 2.

Let $a,b\in\mathbb{Z}$, then $a\mathbb{Z}$, $b\mathbb{Z}$ and
$a+(b-a)\mathbb{Z}$ have non empty intersection.

###### Proof.

Let $d$ be the greatest common divisor of $a$ and $b$. There exist
$a^{\prime},b^{\prime},x,y\in\mathbb{Z}$ such that

$a=da^{\prime}\quad\text{,}\quad b=db^{\prime}$

and

$xa^{\prime}+yb^{\prime}=1\,.$

Adding and subtracting $xb^{\prime}$

$x(a^{\prime}-b^{\prime})+(x+y)b^{\prime}=1\,.$

From which

$(x+y)b^{\prime}=1+x(b^{\prime}-a^{\prime})\,,$

and multiplying by $da^{\prime}=a$, we get

$(x+y)da^{\prime}b^{\prime}=a+xa^{\prime}(b-a)\,.$

The left hand side is in $a\mathbb{Z}\cap b\mathbb{Z}$; the right hand side is
in $a+(b-a)\mathbb{Z}$. ∎

######  Lemma 3.

Let $A,B,C$ be cosets of closed subgroups of $\mathbb{S}^{1}$ such that each
pair is intersecting, then they intersect.

###### Proof.

If one of the given cosets is $\mathbb{S}^{1}$, we are done, so we may assume
that the three cosets are finite. Without loss of generality, we can suppose
that $A$ and $B$ are subgroups of $\mathbb{S}^{1}$ (so that ${1\in A\cap B}$),
for otherwise, take $z\in A\cap B$, and consider $z^{-1}A,z^{-1}B,z^{-1}C$.

By hypothesis, there exist $\alpha\in A\cap C$ and $\beta\in B\cap C$. Observe
that $\alpha\langle\beta\alpha^{-1}\rangle\subset C$, where we denote by
$\langle z\rangle$ the subgroup of $\mathbb{S}^{1}$ generated by
$z\in\mathbb{S}^{1}$.

Let $p=order(\alpha)$ and $q=order(\beta)$. So that
$\langle\alpha\rangle\subset A$ and $\langle\beta\rangle\subset B$ are the
cyclic groups of order $p$ and $q$ respectively.

Let $m$ be the minimum common multiple of $p$ and $q$. And let

$\displaystyle\varphi:\mathbb{Z}$ $\displaystyle\to\mathbb{S}^{1}$
$\displaystyle\varphi(n)$ $\displaystyle=exp(n\frac{2\pi}{m})$

where $exp:\mathbb{R}\to\mathbb{S}^{1}$ is the standard exponential map. By
construction, there exist $a,b\in\mathbb{Z}$ such that $\varphi(a)=\alpha$ and
$\varphi(b)=\beta$.

By Lemma 2, there exists $c\in a\mathbb{Z}\cap
b\mathbb{Z}\cap(a+(b-a)\mathbb{Z})$, so that

$\varphi(c)\in\langle\alpha\rangle\cap\langle\beta\rangle\cap\alpha\langle\beta\alpha^{-1}\rangle\subset
A\cap B\cap C\,.$

∎

Proof of Theorem 2. Let $\mathcal{F}$ be a pairwise intersecting family of
cosets of closed subgroups of $\mathbb{S}^{1}$.

First, we prove the case where $\mathcal{F}$ is finite by induction on the
cardinality, $r$, of $\mathcal{F}$. For $r=3$ it is Lemma 3. Suppose the
Theorem true for families of size $r-1$. Consider $A\in\mathcal{F}$ and let
$\mathcal{F}^{\prime}=\\{\,A\cap X\,|\,X\in\mathcal{F}\setminus\\{A\\}\,\\}$.
By Lemma 3, every pair of elements of $\mathcal{F}^{\prime}$ intersect and by
induction they all have non-empty intersection which is also the intersection
of all the elements of $\mathcal{F}$.

Finally, consider the general case. Suppose that $\mathcal{F}$ is not
intersecting. Then the family of its complements cover $\mathbb{S}^{1}$. They
are open sets, so that the compactness of $\mathbb{S}^{1}$ gives us a finite
open subcover. By the preceding paragraph, in the corresponding finite family
of closed cosets of groups of $\mathbb{S}^{1}$ there must be a non
intersecting pair, which is a contradiction. ∎

Again using Lemma 2, it is easy to prove that not only the group
$\mathbb{S}^{1}$ has this Helly property but also the group $\mathbb{Z}$ for
finite families.

######  Theorem 4.

Let ${\mathcal{F}}$ be a finite pairwise intersecting family of cosets of
subgroups of $\mathbb{Z}$. Then ${\mathcal{F}}$ has non empty intersection.

###### Proof.

Following the proof of Theorem 2, it is enough to consider the case when
${\mathcal{F}}=\\{A,B,C\\}$ consists of three cosets of subgroups of
$\mathbb{Z}$. Therefore $A,B$ and $C$ are arithmetic sequences and without
loss of generality $A=a\mathbb{Z}$, and $B=b\mathbb{Z}$. Assume $a^{\prime}\in
C\cap A$ and $b^{\prime}\in C\cap B$, Consequently,
$a^{\prime}\mathbb{Z}\subset A$, $b^{\prime}\mathbb{Z}\subset B$ and
$a^{\prime}+(b^{\prime}-a^{\prime})\mathbb{Z}\subset C$. By Lemma 2, we have
that ${\mathcal{F}}=\\{A,B,C\\}$ is intersecting. This concludes the proof of
the theorem. ∎

##  4\.  Complex homothetic sections and projections

Let $K_{1},K_{2}$ be two non-homothetic complex ellipsoids centered at the
origin in $\mathbb{C}^{n}$. Hence, for every line $L$, $L\cap K_{1}$ and
$L\cap K_{2}$ are disks centered at the origin and consequently there is
$z_{L}$ such that

$z_{L}(L\cap K_{1})=L\cap K_{2}.$

Nevertheless, $zK_{1}\not=K_{2}$, for $z\in\mathbb{C}$. This shows that
Theorem 1 is false if we replace complex planes with complex lines.

For the proof of Theorem 1, we will start by noting that in its hypothesis we
can assume, without loss of generality, that the module of all $z_{\Lambda}$
is 1.

######  Lemma 4.

Let $K_{1},K_{2}$ be closed subsets of $\mathbb{C}^{n},n\geq 3$ and suppose
that for every complex plane $\Lambda$ there is $z_{\Lambda}\in\mathbb{C}$,
such that

$z_{\Lambda}(\Lambda\cap K_{1})=\Lambda\cap K_{2}.$

Then, we may assume that $|z_{\Lambda}|$ is constant.

###### Proof.

Let us call a subset of $\mathbb{C}^{n}$ _nontrivial_ if it is neither empty
nor the origin, and observe that for a nontrivial set
$A\subset\mathbb{C}^{n}$, if $z_{1}A=z_{2}A$ then $|z_{1}|=|z_{2}|$. We may
clearly assume $K_{1}$ is nontrivial, so that there exists a line $L$ such
that $L\cap K_{1}$ is nontrivial. If $\Lambda_{1}$ and $\Lambda_{2}$ are two
complex planes which contain the complex line $L$, then by hypothesis,

$z_{\Lambda_{1}}(L\cap K_{1})=L\cap K_{2}=z_{\Lambda_{2}}(L\cap K_{1})$

and then, $|z_{\Lambda_{1}}|=|z_{\Lambda_{2}}|=r_{L}$. Now consider a plane
$\Lambda$ that does not contain $L$. If $\Lambda\cap K_{1}$ is trivial, we may
choose $z_{\Lambda}=r_{L}$. Otherwise, there is a line
$L^{\prime}\subset\Lambda$ such that $L^{\prime}\cap K_{1}$ is nontrivial. Let
$\Lambda^{\prime}$ be the plane generated by $L$ and $L^{\prime}$. Then,

$z_{\Lambda}(L^{\prime}\cap K_{1})=L^{\prime}\cap
K_{2}=z_{\Lambda^{\prime}}(L^{\prime}\cap K_{1})\,,$

and therefore, $|z_{\Lambda}|=|z_{\Lambda^{\prime}}|=r_{L}$. ∎

Proof of Theorem 1. First, by Lemma 4 we may assume that $|z_{\Lambda}|=1$ for
every complex plane $\Lambda$. Let $L$ be a complex line and
$L\subset\Lambda$, where $\Lambda$ is a complex plane. Then $z_{\Lambda}(L\cap
K_{1})=(L\cap K_{2}).$ Therefore $z_{\Lambda}$ lies in the coset $H(L\cap
K_{1},L\cap K_{2})\subset\mathbb{S}^{1}$ which is not empty. Let

$\mathcal{F}=\\{\,H(L\cap K_{1},L\cap K_{2})\,\mid\,L\in\mathcal{L}\,\\}\,.$

We shall prove that $\mathcal{F}$ is a pairwise intersecting family of cosets
of closed subgroups of $\mathbb{S}^{1}$. Let $L_{1},L_{2}$ be two different
complex lines and let $\Lambda$ be the complex plane generated by $L_{1}$ and
$L_{2}$. Hence

$z_{\Lambda}\in H(L_{1}\cap K_{1},L_{1}\cap K_{2})\cap H(L_{2}\cap
K_{1},L_{2}\cap K_{2})\not=\emptyset\,.$

By Theorem 2 and Lemma 1, there exists

$z\in\bigcap_{L\in\mathcal{L}}H(L\cap K_{1},L\cap K_{2})=H(K_{1},K_{2})\,.$

Therefore, $zK_{1}=K_{2}$ as we wished. ∎

For the proof of Theorem 3 we need two lemmas.

######  Lemma 5.

If $K_{1},K_{2}\subset\mathbb{C}^{n}$ are either two closed convex sets or two
finite sets such that $\pi_{L}(K_{1})=\pi_{L}(K_{2})$ for every complex line
$L$, then $K_{1}=K_{2}$.

###### Proof.

Suppose first $K_{1},K_{2}\subset\mathbb{C}^{n}$ are two closed convex sets.
If $K_{1}\not=K_{2}$, then without loss of generality there is a point $x\in
K_{1}$ and a real affine hyperplane $\Gamma$ such that $x$ is contained in the
open semispace generated by $\Gamma$ and $K_{2}$ is contained in its
complement. Let $\Delta$ be the unique complex affine hyperplane contained in
$\Gamma$ and let $L$ be the complex line through the origin orthogonal to
$\Delta$. Consider $\pi_{L}:\mathbb{C}^{n}\to L$ be the orthogonal projection.
Then $\pi_{L}(x)\in\pi_{L}(K_{1})$ is contained in the open half plane
generated by the real line $\pi_{L}(\Gamma)$ of $L$ and $\pi_{L}(K_{2})$ is
contained in the complement, contradicting the assumption that
$\pi_{L}(K_{1})=\pi_{L}(K_{2})$.

Suppose now $K_{1},K_{2}\subset\mathbb{C}^{n}$ are two finite sets. If they
are different, then without loss of generality we may chose $x\in
K_{1}\setminus K_{2}$. Since $K_{1}\cup K_{2}$ is finite there exists a
complex line through the origin $L$ such that $\pi_{L}\mid_{K_{1}\cup K_{2}}$
is injective. Then $\pi_{L}(x)\notin K_{2}$, contradicting the assumption that
$\pi_{L}(K_{1})=\pi_{L}(K_{2})$. ∎

######  Lemma 6.

Let $K,K_{1},K_{2}\subset\mathbb{C}^{n}$, $n\geq 3$, be either closed convex
sets or finite sets. Then

  1. (1)

$H(K)=\bigcap_{L\in\mathcal{L}}H(\pi_{L}(K))\,,$

  2. (2)

$H(K_{1},K_{2})=\bigcap_{L\in\mathcal{L}}H(\pi_{L}(K_{1}),\pi_{L}(K_{2})).$

###### Proof.

We know that $H(K)\subset\bigcap_{L\in\mathcal{L}}H(\pi_{L}(K))$ from (2).
Suppose now, ${z\in\bigcap_{L\in\mathcal{L}}H(\pi_{L}(K))}$, that is,
$z\pi_{L}(K)=\pi_{L}(K)$, for every line $L\in\mathcal{L}$. Consequently,
$\pi_{L}(zK)=\pi_{L}(K)$, for every line $L$ and therefore, by Lemma 5,
$zK=K$, which implies that $z\in H(K)$ as we wished.

Moreover, for $z_{0}\in H(K_{1},K_{2})$, by (3) and item (1) of this lemma, we
have:

$\displaystyle H(K_{1},K_{2})$
$\displaystyle=z_{0}H(K_{1})=z_{0}\bigcap_{L\in\mathcal{L}}H(\pi_{L}(K_{1}))=$
$\displaystyle=\bigcap_{L\in\mathcal{L}}z_{0}H(\pi_{L}(K_{1}))=\bigcap_{L\in\mathcal{L}}H(\pi_{L}(K_{1}),\pi_{L}(K_{2})).$

∎

Proof of Theorem 3. As in Lemma 4, let us first prove that we may assume that
$|z_{\Lambda}|$ is constant and does not depend on $\Lambda$. Remember from
Section 2 that if $A$ is a closed subset of $\mathbb{C}^{n}$, $\Delta$ is a
subspace of $\mathbb{C}^{n}$, $z\in\mathbb{C}$ and $zA=A$, then
$z\pi_{\Delta}(A)=\pi_{\Delta}(A)$.

Now, call a line $L$ _nontrivial_ if $\pi_{L}(K_{1})$ is nontrivial. If $L$ is
a nontrivial complex line contained in the complex planes $\Lambda_{1}$ and
$\Lambda_{2}$, then
$z_{\Lambda_{1}}(\pi_{L}(K_{1}))=\pi_{L}(K_{2})=z_{\Lambda_{2}}(\pi_{L}(K_{1}))$
and so, $|z_{\Lambda_{1}}|=|z_{\Lambda_{2}}|.$ If $\Lambda_{1}$ and
$\Lambda_{2}$ are two complex planes with nontrivial
$\pi_{\Lambda_{i}}(K_{1})$, there is a third complex plane $\Lambda_{3}$ such
that $\Lambda_{i}\cap\Lambda_{3}$ is a nontrivial complex line, $i=1,2$, so we
obtain that $|z_{\Lambda_{1}}|=|z_{\Lambda_{3}}|=|z_{\Lambda_{2}}|$.
Consequently, we may assume without loss of generality that $|z_{\Lambda}|=1$,
for every complex plane.

This time, let

$\mathcal{F}=\\{H(\pi_{L}(K_{1}),\pi_{L}(K_{2}))\mid L\in\mathcal{L}\\}\,.$

We shall prove that $\mathcal{F}$ is a pairwise intersecting family of cosets
of closed subgroups of $\mathbb{S}^{1}$. Let $L_{1},L_{2}$ be two different
complex lines and let $\Lambda$ be the complex plane generated by $L_{1}$ and
$L_{2}$. Hence,

$z_{\Lambda}\in H(\pi_{L_{1}}(K_{1}),\pi_{L_{1}}(K_{2}))\cap
H(\pi_{L_{2}}(K_{1}),\pi_{L_{2}}(K_{2}))\,.$

By Theorem 2, there is

$z\in\bigcap_{L\in\mathcal{L}}H(\pi_{L}(K_{1}),\pi_{L}(K_{2}))\,.$

Therefore, by Lemma 6(2), $z\in H(K_{1},K_{2})$. ∎

Acknowledgments. Luis Montejano acknowledges support from CONACyT under
project 166306 and from PAPIIT-UNAM under project IN112614. Javier Bracho
acknowledges support from PAPIIT-UNAM under project IN109023.

## References

  * [1] M. Alfonseca, M. Cordier and D. Ryabogin, _On bodies with directly congruent projections and sections_ , Israel J. Math., 215 (2016), no. 2, 765-799
  * [2] J. L. Arocha, J. Bracho, and L. Montejano, _Complex ellipsoids and complex symmetry_ , Aequationes mathematicae, 97, 295–306 (2023)
  * [3] V. P. Golubyatnikov, _Uniqueness Questions in Reconstruction of Multidimensional Objects from Tomography Type Projection_ Data, Inverse and Posed Problems Series. Utrecht, Boston (2000).
  * [4] D. Ryabogin, _On the continual Rubik’s cube_ Advances in Mathematics 231 (2012) 3429–3444
  * [5] D. Ryabogin, _On symmetries of projections and sections of convex bodies_ , Discrete Geometry and Symmetry. In: Marston, D.E., Conder, A.D. and Weiss, A.I. (eds.) Honor of Károly Bezdek’s and Egon Schulte’s 60th Birthdays. Springer Proceedings in Mathematics and Statistics, 2017.
  * [6] D. Ryabogin, _On a functional equation related to convex bodies with SU(2)-congruent projections_ To appear in Geometriae Dedicata.
  * [7] C. A. Rogers, _Sections and projections of convex bodies._ Port. Math. 24, 99–103 (1965)