# The Power Spectrum of Climate Change

Albert Sneppen  e-mail<EMAIL_ADDRESS>(corresponding author) Niels Bohr
Institute, University of Copenhagen, Blegdamsvej 17, København Ø 2200,
Denmark.  Cosmic Dawn Center (DAWN)

###### Abstract

Both global, intermediate and local scales of Climate Change have been studied
extensively, but a unified diagnostic framework for examining all spatial
scales concurrently has remained elusive. Here we present a new tool-set using
spherical harmonics to examine climate change through surface temperature
anomalies from 1850 to 2021 on spatial scales ranging from planetary to 50 km
scales. We show that the observed temperature anomalies are accurately
decomposed in spherical harmonics typically within 0.05 K. This decomposition
displays a remarkably simple dependence on spatial scale with a universal
shape across seasons and decades. The decomposition separates two distinct
regimes by a characteristic turnover-length of approximately 3000 km. The
largest scales confirm established trends, while local fluctuations are
consistent with 2-dimensional turbulence. We observe a downward cascade, from
which it follows that climate change feeds increasing volatility on all
spatial scales from $2.000$ to $50$ km. This increase is primarily driven by
growing volatility along longitudes.

##  1 Introduction

Climate Change affects the Earth on a global scale, but its impact exhibits
inhomogeneities over many spatial scales. Both large scale features of climate
change [1, 2, 3] and the evolution of local fluctuations in weather and
climate [4, 5] have been extensively studied. However, the multi-scale nature
of climate-change emphasises the need for a unified diagnostic examining
correlations on all spatial scales [6, 7, 8, 9]. Historically, Fourier
decomposition of temperature-fluctuations has been utilised to determine the
relative importance of different spatial scales [10, 11]. Such spectral domain
analyses has revealed atmospheric fields (such as absolute temperature, wind
velocity, pressure, etc.) are commonly scale-free over a wide dynamical range
[12, 13, 14, 15]. However, the evolution of spectral coefficients remain
unexplored across the time-scales of climate change.

Therefore, this article introduces temperature anomalies to the spatial
spectral domain and uniquely constrains the temporal evolution of this
spectral regime. The analysis extends across a wide range of spatial scales
from planetary to 50 km resolution with a monthly temporal resolution from
1850 to 2022. This is achieved by presenting the climate and its evolution on
a range spatial scales through a collection of spherical harmonics. This basis
is chosen across many fields of physics, because the visual nature of the
different modes of spherical harmonics allows an interpretable decomposition
of any scalar field, such as temperature fluctuations, on a sphere.

The significance of different spatial scales is described by the power
spectrum of spherical harmonics. Therefore, we determine the spectrum’s shape
(see § 3.2), where large spatial scales confirm previous results within
climatology and small scale features follow the analytical predictions of
turbulence confined to two dimensions. Furthermore, we show that temperature
anomalies on these different scales are causally connected because variations
cascade from large to small spatial scales.

We extend the analysis with a temporal evolution, where the contribution of
each mode is detailed across time. Thus, in § 3.1 we explore the evolution of
the lowest modes with key results conforming to established observational
trends within climatology. Furthermore, the time-evolution of the spectrum
(see § 3.3) quantifies how climate change results in an increased magnitude
and frequency of local temperature fluctuations. This increase in volatility
is not isotropic, but predominantly due to increased fluctuations along
longitudes not latitudes (see § 3.5). Ultimately, spherical harmonics prove an
intuitive framework in decomposing the local effects of global warming.

Figure 1: Projected visualisation of first three degrees ($l={0,1,2}$) of
spherical harmonics with order $m=0$. Left: The 0th spherical harmonics is an
average value around the sphere. Center: illustrates the different behaviour
of the Northern and Southern hemisphere. Right: emphasises the difference
between the poles and the equator.

##  2 Methods

###  2.1 An Introduction to Spherical Harmonics

Spherical harmonics are ortho-normal base functions on the surface of a
sphere. They are in this analysis used to decompose the temperature anomalies
across the planet, $T(\theta,\varphi,t)$ at any given time, $t$:

$T(\theta,\varphi,t)=\sum_{l=1}^{\infty}\sum_{m=-l}^{l}a_{l,m}(t)Y_{l}^{m}(\theta,\varphi)$
(1)

where $a_{l,m}$ is the coefficient of the spherical harmonic
$Y_{l}^{m}(\theta,\varphi)$, which is given by the Legendre-polynomial
$P_{l}^{m}(\cos\theta)$:

$Y_{l}^{m}(\theta,\varphi)=(-1)^{m}\bigg{[}\frac{2l+1}{4\pi}\frac{(l-m)!}{(l+m)!}\bigg{]}^{1/2}P_{l}^{m}(\cos\theta)e^{im\varphi}$
(2)

The shape of $Y_{l}^{m}$ is described by the degree, $l$, which determines the
number of oscillations around the sphere. In this context, the l’th harmonic
is a description of fluctuations on angular scales
$\theta\approx\frac{180}{l}$, which correspond to spatial scales of
$\frac{20.000km}{l}$. For instance, the monopole, $l=0$, gives the average
value around the sphere; The dipole, $l=1$, shows the behaviour of the
opposing hemispheres with different signs; The quadruple, $l=2$, describes the
behaviour of the opposing hemispheres with the same sign (see Fig. 1).

Furthermore, the order, $m$ (ranging from $-l$ to $l$), determines whether the
waves are oriented along longitudes ($m=0$) or latitudes ($m=\pm l$).
$a_{l,m}$ is the coefficient of the $Y_{l}^{m}$’th harmonic and are determined
from the temperature field $T(\theta,\varphi,t)$:

$a_{l,m}(t)=\int_{-1}^{1}\int_{0}^{2\pi}[Y_{l}^{m}(\theta,\phi)]^{*}\,T(\theta,\phi,t)\,d\phi\;d(\cos\theta)$
(3)

Here, $|a_{l,m}|^{2}$ represents the contribution of any spherical harmonic to
the decomposed function, so the significance of a degree is given by the
power, $C_{l}$:

$C_{l}(t)=\frac{1}{2l+1}\sum_{m=-l}^{l}|a_{l,m}(t)|^{2}$ (4)

Thus, investigating the power spectrum is a simple diagnostic tool for
characterizing the magnitude of temperature fluctuations as a function of
angular scale.

Figure 2: Left: Observed temperature anomalies, Center: Spherical Harmonical
Reconstruction using all degrees $l<90$, Right: Residual between observed and
reconstructed for May 2021. Note, we emphasise discrepancies by increasing the
contrast on residuals. The magnitude of residuals is typically small with 50%
of the surface being within 0.05 degrees.

###  2.2 Data-sets and Potential Bias

The data was collected from The Berkeley Earth Land/Ocean Temperature Record
[3]. The data set is a combination of the Berkeley Earth land-surface
temperature field [16] and a reinterpolated version of the HadSST ocean
temperature field [17]. It provides the monthly temperature anomalies along
with their uncertainties across the globe in a 360$\times$180 grid (longitude
vs. latitude), spanning a period of 171 years between Jan 1850 and Dec 2021.
The anomalies are reported relative to the locally inferred average
temperature for that month over the period Jan 1951-Dec 1980 [for further
details see [3]]. Note, that both absolute temperature and temperature
anomalies can be decomposed in spherical harmonics, but the former is affected
by systematic shifts [for instance the dependence of temperatures with
elevation]. For completeness the power spectrum of absolute temperature is
computed in Appendix B. However, subtracting the underlying mean, we decouple
these systematic shifts and thereby examine only the structure of dynamical
temperature fluctuations caused by kinematics of the atmosphere. Thus, we
follow the convention of utilizing anomalies, which has been shown to more
accurately describe climate variability over larger areas [3].

Importantly, the coverage of the grid is a dynamic quantity ranging from 57 %
in 1850 to 99 % coverage in 2015. Limited coverage can induce a bias in the
estimated spherical harmonical decomposition [18, 19]. However, the
decomposition in this analysis is robust within this range in coverage, as the
overall trends remain unchanged when limiting the analysis to the coverage
used in 1850 (see Appendix C). Additionally, the significant increase of
active recording stations could bias any temporal trend of large modes as more
stations may resolves features on finer scales. From 1850 to 1960 the number
of active land-temperature stations and sea-temperature observations increased
60-fold and 20-fold respectively (see Appendix D); Over the same period only a
slight evolution of the power-spectrum is observed. In contrast, from 1960 to
2020 the number of active land-temperature stations and sea-temperature
observations stayed relatively consistent, with an increase of 15% and 20%
respectively (see Appendix D). However, in this period a drastic growth is
observed within the power-spectrum. We emphasise, that increasing the number
of observations by orders of magnitude still yields consistent power-spectra,
so the slightly increased statistics of the last 60 years does explain the
recent and drastic evolution of the power-spectrum. For completeness the paper
presents the decomposition over all 171 years, but the emphasis remains on the
recent evolution from 1960 to 2021, where the coverage remains above 95 % at
all times. Furthermore, this limited period mitigates the potential historical
biases in SST from the early 20th century [20].

The uniform gridding of the Berkeley Earth Land/Ocean Temperature Record is
achieved through an interpolated field with a kriging-based approach [for
detail on implementation see [3]]. Such interpolation will potentially
introduce a bias by smoothing structures on the spatial scales of
interpolation. However, the coverage grows drastically and in extension
interpolation-lengths diminish from 1850-1960. If this bias was significant
within this data-set the power of intermediate and high modes should increase
from 1850 to 1960 as the suppression of large interpolation-lengths
diminishes. Nevertheless, intermediate and large modes [eg. $C_{20}(t)$,
$C_{50}(t)$ or $C_{80}(t)$] show no consistent growth throughout this period.
From 1960-2020 the sampling over land is dense containing approximately 20.000
active measurement sites with typically around 100 km between stations. A bias
from smoothing on scales below 100 km will not significantly affect the
planetary or mesoalpha scales analysed below (see Appendix E). In contrast,
the sampling over sea is relatively sparse with up to 500 km between grid-
cells in HadSST; These interpolation-scales could theoretically suppress modes
of degree $l>40$. However, note this suppression is not observationally
suggested from the power-spectrum itself with modes $l=90$ following the
temporal trends and spatial structure of degrees well-outside the potentially
suppressed regime (eg. $l=20$). This article provides all spherical harmonics
to spatial scales of 200 km, but we caution against any interpretation of
structures (especially over sea) below 500 km.

In addition to the data-resolution differences between land and sea, there
exist well established kinematic and dynamical differences between land and
sea [16, 17]. Therefore, we also present the relative contribution of sea and
land temperatures (see Appendix F). Here we show that land-temperatures
(despite the smaller surface area) are the dominant contributor to the shape
of the power-spectrum.

Lastly, it is emphasised that the recorded findings presented remain robust
across any coarser grid and any temporal resolution from the scale of decades
to days. The monthly time-resolution was chosen due to its global coverage
(including both land and ocean temperatures). Furthermore, the presented
results are insensitive to the specific choice of analysing The Berkeley Earth
Land/Ocean Temperature Record with other coincident catalogues yielding a
similar functional form of the power spectrum with a similar evolution across
seasons and decades. In Appendix A, we present a similar power-spectral
decomposition of the ERA5 data-set.

###  2.3 Decomposition Accuracy

The spherical harmonical decomposition is computed with the well-documented
and readily accessible pyshtools-software [19]. With increasingly detailed
coverage the estimated Spherical Harmonics decomposition is ensured to
converge towards the observed temperature-landscape. However, for any limited
spatial resolution the spherical harmonics decomposition remains an
approximation of the underlying landscape. In Fig. 2 we show an example of a
decomposition of temperature anomalies to illustrate the high accuracy of the
approach. Utilizing all spherical harmonics up to degree $l=90$, does provide
an excellent estimation of the underlying temperature-field. Typical
deviations between observations and the spherical harmonic reconstructions
(using up to degree $l=90$) are 0.05 K. Furthermore, the observed and
reconstructed temperature anomalies are not dominated by systematic trends
with altitude or latitude. Such effects are decoupled from the investigated
signal because we utilise the temperature anomalies and not the absolute
temperature.

##  3 Results

###  3.1 The Lowest Modes

With spatially resolved temperatures, we can determine the weights,
$a_{l,m}(t)$, of the spherical harmonic functions for the scalar temperature
field. To exemplify the information available in the time-series of spectral
weights, we show the lowest degree coefficients from 1850 to 2021 in Fig. 3.
Each of these time-series contain unique perspectives on the spatial
distribution of temperature anomalies, which will now be shortly summarized.

Figure 3:  $a_{0,0},a_{1,0},a_{2,0}$ coefficients of the spherical harmonics
derived from global temperature anomalies resolved from Jan 1850 - Jan 2021.
Figure 4:  Right: Power of select modes from 1950-2020 with linear fit
overlayed. No systematic temporal evolution of power is seen for the mode l=5,
but the power on smaller spatial scales is increasing towards the present.
Left: Power spectrum of temperature distribution from 1850s (in black) to
2010s (in red). Dashed-dotted purple line indicates a generic powerlaw with
powerlaw index $\alpha=-3$. The power of modes $l\geq 10$ (ie. on spatial
scales smaller than 2000 km) is increasing with time especially from 1960 to
2020. The identical power-spectrum of the ERA5-catalogue (which has finer
spatial resolution but less temporal range) is shown in Appendix A.

Global Temperature $(a_{0,0})$: First, the lowest degree (i.e. the spherically
symmetric mode) increases towards the present. This mode’s evolution
demonstrates the increase of the global average temperature [3].

The North-South Inequality $(a_{1,0})$: The $Y_{1}^{0}(\theta,\varphi)$
spherical harmonic describes variations between northern and southern regions
of the sphere. Here $a_{1,0}<0$ implies that the increase in temperature is
larger for the Northern Hemisphere than for the Southern Hemisphere.
Therefore, the evolution of $a_{1,0}(t)$ reaffirms the well-known
interhemispheric temperature asymmetry [21]. In the new millennium the
Northern Hemisphere has experienced a majority of the temperature increase
with the annual interhemispheric temperature asymmetry being approximately 0.4
°C in 2010 [3].

Increased Heating of Polar Regions $(a_{2,0})$: The
$Y_{2}^{0}(\theta,\varphi)$ spherical harmonic describes variations between
polar and equatorial regions of the sphere. Thus, the growth of $a_{2,0}$ in
recent decades validates the well-known polar amplification [22, 23]. Note,
the north-south symmetry of $Y_{2}^{0}(\theta,\varphi)$ means that $a_{2,0}$
does not constrain the asymmetry of polar amplification [24].

###  3.2 The Power Spectrum

While the lowest degrees of spherical harmonics remain easily interpretable,
the entire power spectrum is required for a comprehensive understanding of the
spatial structure. In Fig. 4 the power-spectrum is shown for each decade from
1850 to 2020.

The overall shape of the power-spectrum is unchanged across decades; the power
is constant on the largest of scales, but at a turnover length-scale
transitions to a powerlaw-like decline. These two different regimes are
dominated by different processes.

First, the observed power is approximately flat for $2\leq l\leq 7$
corresponding to a spatial wavelength from $10.000$ to $3.000$ km. These
largest scales contain dynamics on the size of global Rossby waves, large
scale tropical waves and the Madden-Julian Oscillation [25]. Thus, white noise
is observed for temperature anomalies on these planetary scales.

Second, after the turnover higher degree spherical harmonics become
progressively less important following a powerlaw like decline, where the
power $C_{l}\propto l^{\alpha}$. For $l>10$ the power of each harmonic
$Y_{l}^{m}$ decreases with a powerlaw index which is consistent with
$\alpha=-3$ for all parts of the time-series. The powerlaw is observed from
spatial scales of 2000 km to the spatial resolution of the data-set (ie. order
50 km - see Fig. A1).

Notable, several atmospheric fields including humidity, temperature, zonal-,
meridional- and vertical winds at varying altitudes of the ECMWF interim
analysis have empirically shown a similar scale-free spectra with varying
powerlaw slopes [26]. However, the specific slope $\alpha=-3$ is remarkable in
accordance with the analytical prediction of turbulent flows constrained to 2
dimensions. Under such conditions, the expectation is that the energy spectrum
scales as $E(k)\propto k^{-3}$, where the wave-number, $k$, is the 2D-spatial
analog of the angular scales of $l$ [27]. Note, this characteristic scaling of
$k=-3$ has been observed previously for potential temperature near the
tropopause across planetary and synoptic scales [10].

An important test of the classical models of 2D/3D turbulence is to examine
the spectral energy transfer. For 2D the energy cascades from small to larger
scales, while enstrophy is characterised by a downscale transfer. Previous
empirical determinations have yielded little conclusive constraints, although
3rd order velocity structure corroborate the downscale nature of atmospheric
cascades [28, 29]. Another analytical perspective on the cause of powerlaw-
tails in power spectra is discussed in Ahlers2014 [30]. Here Boltzmann’s
equation is used to derive a coupled set of differential equations relating
the power of each mode to the power of neighbour modes, assuming a simple non-
linear coupling associated to turbulence. These coupled equations are solved
assuming hierarchical generation, where large spatial fluctuations generate
fluctuations on smaller scales. This formalism, suggests that a downward
cascade of fluctuations on a sphere would produce a power-spectrum of the form
$C_{l}\propto l^{-3}$.

Regardless of the causal direction of the cascade a generic prediction of
cascades is that an increasing amplitude on large scale variations must have a
corresponding increase on all smaller spatial scales. Given the observed
growth in the power of the lowest modes (eg. $a_{1,0},a_{2,0}$) we must
therefore expect an increasing power across the spectrum with time.
Furthermore, in § 3.4, § 3.5 and in Appendix G we present evidence of a
downward cascade with large scale features being driven externally (such as by
differential solar heating and variations over latitudes) which turbulence
feeds into a range of fluctuations on all smaller spatial scales. Importantly,
this downward cascade predicts that the increasing amplitude of fluctuations
on large spatial scales due to climate change, also feeds more volatility for
weather on more local scales.

Figure 5: Power spectrum marginalized over month on northern hemisphere. The
amplitude of temperature fluctuations on all spatial scales are larger in
winter than in summer. Naturally, for the southern hemisphere the order is
inverted with most power in June and July.

###  3.3 Yearly Evolution of The Power Spectrum

From 1850 to 1960 the power spectrum shows variability but no systematic
coherent growth, but a few noteworthy features emerge after 1960 (see Fig. 4).
First, the power of the lowest degrees is increasing towards the present as
discussed in § 3.1. Second, the flat part of the power-spectrum (degrees
$2\leq l\leq 8$) remains unchanged up to statistical fluctuations. Third, the
turnover length-scale increases, ie. the powerlaw cascade begins at lower
modes (for elaboration see appendix Fig. A10). An analogous shift of the
energy spectrum to larger spatial scales is often correlated with increased
eddy activity and weather volatility [11]. Lastly, the power of higher degrees
($10\leq l$) are increasing on all scales from 2000 km to the resolution of
the data-set (ie. order 50 km for ERA5) from 1960 to the present. Local
temperature fluctuations are growing in magnitude and the growth is
accelerating over the past century.

###  3.4 Seasonal Evolution of The Power Spectrum

It is observationally well-established that the kinetic energy in the
atmosphere is larger in winter than in summer [11]. In summer the meridional
temperature gradient is both weaker and shifted poleward, resulting in a lower
eddy activity and a shift to smaller spatial scales. These seasonal features
are also seen in the power-spectrum of temperature fluctuations (see Fig. 5).
Summer-months shows decreased fluctuations on all scales (ie. smaller eddy
activity), while the transition to the powerlaw decline is shifted to higher
wave-numbers (for further elaboration see Appendix H). Importantly, this
indicates that variations on small spatial scales follow the growth and decay
of fluctuations on the largest scale (ie. that we are observing a downward
cascade).

Figure 6: Color indicates power of the $Y_{l}^{m}$’th spherical harmonic for
the period 1960-1970 [left] and 2010-2020 [middle]. For both periods the power
is calculated each month and then averaged over a decade. Right shows
residual, suggesting that the increasing power of higher order modes is not
symmetric across $m$, but dominated by the increasing power of modes from
north-to-south [ie. $m=-l/2$ to $m=l/2$, indicated with dashed lines].

Beyond the yearly fluctuation no other periodicity is detected across the
power-spectrum. Therefore, the power of modes on smaller spatial scales only
oscillate seasonally (combined with the previously mentioned growth over
decades). While the mean global temperature is a complex convolution of
anthropogenic climate change and non-anthropogenic oscillations [31], the
fluctuations from the mean (described by smaller modes) suggests a singular
temporal trend. Thus, the power-spectrum yields a diagnostic with a simple and
predictive evolution across time.

###  3.5 The Horizontal Anisotropy of Climate Change

The expectation from atmospheric reanalysis fields is that horizontal
anisotropies are common, with fluxes elongated more zonally than meriodionally
[26]. The easily interpretable decomposition of spherical harmonics does
provide a straightforward framework for examining horizontal anisotropies, but
by summing over $m$ we have yet to discuss the spatial orientation of these
fluctuations. In spherical harmonics, spherical symmetry implies every $m$
should contribute equally to the power of the $l$’th mode. This is elaborated
through simulated temperature fields in Appendix I. Here it is shown that when
generating temperature fluctuations with isotropic variances along latitudes
and longitudes the power is uniform across $m$. However, if the length-scale
of temperature fluctuations is smaller along longitudes than latitudes the
modes from $m\approx-\frac{l}{2}$ to $m\approx\frac{l}{2}$ are shown to be the
dominant contribution to the power.

Indeed, anisotropy is observed throughout the timeseries as seen in Fig. 6
[left and middle panel]. Evidently, the climate generically has a larger
magnitude of temperature fluctuations meridionally than zonally, as expected
for a rotating sphere [26]. However, in contrast to previous analyses the
detailed timeseries of global temperatures allows computing the change in
power across decades [right panel]. We see that on all spatial scales
resolvable it is mainly fluctuations oriented north-to-south which are
increasing in amplitude. This coherence is a prediction of downward cascades,
where the observed growing fluctuation on large scales (eg. $a_{1,0},a_{2,0}$,
etc.) feed increasing amplitudes on all subsequent scales. Because the large
scale fluctuations are oriented north-to-south the cascade displays a similar
anisotropy. In contrast, an inverse cascade should not display any horizontal
anisotropy on small scales, as there are no local differences between latitude
and longitude. Climate change therefore results in rapid weather transitions
which are growing in frequency and magnitude, but crucially the increase is
primarily driven by increased volatility along longitudes.

##  4 Conclusion

Spherical harmonics provide a succinct description of the Earth’s climate.
Within this paper we prove that the lowest 90 modes can reconstruct global
temperature anomalies with a typically accuracy of 0.05 K. This highly
accurate decomposition is the result of a simple combination of oscillations
on different spatial scales, which is succinctly summarized by the shape of
the power-spectrum. The largest power is coming from modes $l\leq 8$, with
progressively smaller scales contributing less and less. Importantly, the
decline in power on spatial scales from $2.000$ to $50$ km follows a powerlaw
with exponent $\alpha=-3$, which is characteristic of 2 dimensional turbulence
and other theoretical models of cascading fluctuations. Furthermore, mapping
the seasonal fluctuations and anisotropies between longitudes and latitudes we
show that large scale temperature variations cascades into fluctuations on all
smaller spatial scales.

As previously observed large scale fluctuations (which are oriented from
north-to-south) are growing in magnitude. Therefore, the observed cascade
dictates that on any scale the amplitude of fluctuations along longitudes has
to increase. From the planetary scale to that of weather fronts, temperature
oscillations are becoming more volatile due to climate change. Ultimately, the
methodology provides an illuminating perspective on our evolving climate by
bridging the gap between known large scale structures of the climate and small
scale volatility.

Furthermore, the methodology introduced within the paper can be extended with
significant theoretical and observational benefits. First, larger spatial
resolution will allow improvements on the accuracy of reconstructions by
decomposing the structure at even smaller spatial scales. The power structure
of higher modes could reveal the extent and the small-scale limit of the
observed cascade; a more detailed grid could even examine the regime of
3-dimensional turbulence. Second, larger temporal resolution allows a closer
examination of the flow of fluctuations, i.e. how the spatial scales are
causally connected. Thereby, probing the underlying physics connecting the
cascade (discussed in Appendix G). Lastly, we emphasise that the observed
power-spectrum is remarkable simple. Three parameters, a normalisation, a
turnover-length and a powerlaw-slope, characterize the entire power spectrum
of temperature anomalies for over 170 years of observations from the planetary
scale through the mesoalpha. This remarkable universality should be examined,
reproduced and predicted by future climate models in order to further
constrain the spatial scales of climate change.

## Acknowledgements

The author would like to thank Mogens Høgh Jensen, Peter Ditlevsen, Markus
Ahlers and Rune Thinggard Hansen for useful discussions on drafts and
insightful feedback. Furthermore, Kevin Kumar, Gergely Friss, Gustav Ernst
Madsen, Rasmus Damgaard Nielsen and Jason Koskinen all deserve praise for kind
encouragement with the conceptual formulation. Lastly, Julie Kiel Holm
provided the audacity required for publication. The Cosmic Dawn Center (DAWN)
is funded by the Danish National Research Foundation under grant No. 140.

## Data Availability Statement

The data was collected from The Berkeley Earth Land/Ocean Temperature Record,
which can be accessed at http://berkeleyearth.org/data/.

## References

  * Screen and Simmonds [2014] J. Screen and I. Simmonds, Amplified mid-latitude planetary waves favour particular regional weather extremes, Nature Climate Change 4, 704 (2014).
  * Kornhuber  _et al._ [2019] K. Kornhuber, D. Coumou, E. Vogel, C. Lesk, J. F. Donges, J. Lehmann, and R. M. Horton, Amplified Rossby waves enhance risk of concurrent heatwaves in major breadbasket regions, Nature Climate Change 10, 48 (2019).
  * Rohde and Hausfather [2020] R. A. Rohde and Z. Hausfather, The berkeley earth land/ocean temperature record, Earth System Science Data 12, 3469 (2020).
  * Sutton  _et al._ [2015] R. Sutton, E. Suckling, and E. Hawkins, What does global mean temperature tell us about local climate?, Philosophical Transactions of the Royal Society of London Series A 373, 20140426 (2015).
  * Vincze  _et al._ [2017] M. Vincze, I. D. Borcia, and U. Harlander, Temperature fluctuations in a changing climate: an ensemble-based experimental approach, Scientific Reports 7, 254 (2017), arXiv:1702.07048 [physics.flu-dyn] .
  * Donges  _et al._ [2009] J. F. Donges, Y. Zou, N. Marwan, and J. Kurths, Complex networks in climate dynamics. Comparing linear and nonlinear network construction methods, European Physical Journal Special Topics 174, 157 (2009).
  * Donges  _et al._ [2015] J. F. Donges, I. Petrova, A. Loew, N. Marwan, and J. Kurths, How complex climate networks complement eigen techniques for the statistical analysis of climatological data, Climate Dynamics 45, 2407 (2015).
  * Di Cecco and Gouhier [2018] G. Di Cecco and T. Gouhier, Increased spatial and temporal autocorrelation of temperature under climate change, Scientific Reports 10.1038/s41598-018-33217-0 (2018).
  * Zhang and Moore [2015] Z. Zhang and J. C. Moore, Chapter 6 - empirical orthogonal functions, in  _Mathematical and Physical Fundamentals of Climate Change_, edited by Z. Zhang and J. C. Moore (Elsevier, Boston, 2015) pp. 161–197.
  * Nastrom and Gage [1985] G. D. Nastrom and K. S. Gage, A Climatology of Atmospheric Wavenumber Spectra of Wind and Temperature Observed by Commercial Aircraft., Journal of Atmospheric Sciences 42, 950 (1985).
  * Straus and Ditlevsen [1999] D. M. Straus and P. Ditlevsen, Two-dimensional turbulence properties of the ECMWF reanalyses, Tellus Series A 51, 749 (1999).
  * Lovejoy  _et al._ [2008] S. Lovejoy, D. Schertzer, M. Lilley, K. Strawbridge, and A. Radkevich, Scaling turbulent atmospheric stratification. i: Turbulence and waves, Quarterly Journal of the Royal Meteorological Society 134, 277 (2008).
  * Lovejoy  _et al._ [2010] S. Lovejoy, A. F. Tuck, and D. Schertzer, Horizontal cascade structure of atmospheric fields determined from aircraft data, Journal of Geophysical Research (Atmospheres) 115, D13105 (2010).
  * Chen  _et al._ [2016] W. Chen, S. Lovejoy, and J.-P. Muller, Mars’ atmosphere: The sister planet, our statistical twin, Journal of Geophysical Research (Atmospheres) 121, 11,968 (2016).
  * Cavanaugh  _et al._ [2017] N. R. Cavanaugh, T. A. O’Brien, W. D. Collins, and W. C. Skamarock, Spherical Harmonic Spectral Estimation on Arbitrary Grids, Monthly Weather Review 145, 3355 (2017).
  * Muller  _et al._ [2013] R. Muller, R. Rohde, R. Jacobsen, E. Muller, and C. Wickham, A new estimate of the average earth surface land temperature spanning 1753 to 2011, Geoinformatics & Geostatistics: An Overview 01 (2013).
  * Kennedy  _et al._ [2019] J. J. Kennedy, N. A. Rayner, C. P. Atkinson, and R. E. Killick, An ensemble data set of sea surface temperature change from 1850: The met office hadley centre hadsst.4.0.0.0 data set, Journal of Geophysical Research: Atmospheres 124, 7719 (2019), https://agupubs.onlinelibrary.wiley.com/doi/pdf/
10.1029/2018JD029867 .

  * Ahlers [2016] M. Ahlers, Deciphering the Dipole Anisotropy of Galactic Cosmic Rays, Physical Review Letters 117, 151103 (2016), arXiv:1605.06446 [astro-ph.HE] .
  * Wieczorek and Meschede [2018] M. Wieczorek and M. Meschede, Shtools: Tools for working with spherical harmonics, Geochemistry, Geophysics, Geosystems 19 (2018).
  * Chan  _et al._ [2019] D. Chan, E. C. Kent, D. I. Berry, and P. Huybers, Correcting datasets leads to more homogeneous early-twentieth-century sea surface warming, Nature (London) 571, 393 (2019).
  * Friedman  _et al._ [2013] A. R. Friedman, Y.-T. Hwang, J. C. H. Chiang, and D. M. W. Frierson, Interhemispheric Temperature Asymmetry over the Twentieth Century and in Future Projections, Journal of Climate 26, 5419 (2013).
  * Budyko [1969] M. I. Budyko, The effect of solar radiation variations on the climate of the earth., Tellus 21, 611 (1969).
  * Goosse  _et al._ [2018] H. Goosse, J. E. Kay, K. C. Armour, A. Bodas-Salcedo, H. Chepfer, D. Docquier, A. Jonko, P. J. Kushner, O. Lecomte, F. Massonnet, H.-S. Park, F. Pithan, G. Svensson, and M. Vancoppenolle, Quantifying climate feedbacks in polar regions, Nature Communications 9, 1919 (2018).
  * Salzmann [2017] M. Salzmann, The polar amplification asymmetry: role of Antarctic surface height, Earth System Dynamics 8, 323 (2017).
  * Zhang [2013] C. Zhang, Madden–julian oscillation: Bridging weather and climate, Bulletin of the American Meteorological Society 94, 1849 (2013).
  * Lovejoy and Schertzer [2011] S. Lovejoy and D. Schertzer, Space-time cascades and the scaling of ECMWF reanalyses: Fluxes and fields, Journal of Geophysical Research (Atmospheres) 116, D14117 (2011).
  * Kraichnan [1967] R. H. Kraichnan, Inertial Ranges in Two-Dimensional Turbulence, Physics of Fluids 10, 1417 (1967).
  * Lindborg [1999] E. Lindborg, Can the atmospheric kinetic energy spectrum be explained by two-dimensional turbulence?, Journal of Fluid Mechanics 388, 259 (1999).
  * Cho and Lindborg [2001] J. Y. N. Cho and E. Lindborg, Horizontal velocity structure functions in the upper troposphere and lower stratosphere. 1. Observations, Journal of Geophysical Research 106, 10,223 (2001).
  * Ahlers [2014] M. Ahlers, Anomalous Anisotropies of Cosmic Rays from Turbulent Magnetic Fields, Physical Review Letter 112, 021101 (2014), arXiv:1310.5712 [astro-ph.HE] .
  * Wang [2019] J.-L. Wang, Key Modes for Time-Space Evolutions of ENSO and PDO by ESMD Method, arXiv e-prints , arXiv:1912.02013 (2019), arXiv:1912.02013 [physics.ao-ph] .
  * Lilly [1989] D. K. Lilly, Two-Dimensional Turbulence Generated by Energy Sources at Two Scales., Journal of Atmospheric Sciences 46, 2026 (1989).
  * Wieczorek and Simons [2005] M. A. Wieczorek and F. J. Simons, Localized spectral analysis on the sphere, Geophysical Journal International 162, 655 (2005), https://academic.oup.com/gji/article-pdf/162/3/655/6040509/162-3-655.pdf .

Figure A1: Left: Power spectrum of temperature distribution from 1950s (in
black) to 2010s (in red) from ERA5 catalogue. Dashed-dotted purple line
indicates a generic powerlaw with powerlaw index $\alpha=-3$. The power of
modes $l\geq 10$ (ie. on spatial scales smaller than 2000-50 km) is increasing
with time.

## Appendix

###  A Power Spectral Analysis of ERA5 Dataset

While the above analysis centers on The Berkeley Earth Land/Ocean Temperature
Record - other coincident catalogues tell a similar story. In Fig. A1, the
power spectrum for each decade from temperature anomalies in the ERA5
reanalysis (from 1979 to present) and the preliminary ERA5 reanalysis (from
1950 to 1978) is illustrated. The anomalies are defined relative to the
locally inferred average temperature for that month over the period Jan
1951-Dec 1980 (ie. similar to the Berkeley Earth Land/Ocean Temperature
Record). Crucially, this shows an identical functional form of the power
spectrum with a similar normalization, turnover length-scale and powerlaw-
decline (see Fig. 4). Additionally, the evolution with time displays a
coherent growth of power for modes within and above the turbulent cascade.

Note, the analysis of the 30km spatial resolution of ERA5 allows an
examination of small modes ($90<l<360$). The previous observed powerlaw-regime
with $\alpha=-3$ is observed within the entire range from 2000 to 50 km. The
scale-height of the atmosphere (order 10 km) still remains significantly
smaller than the horizontal resolution, so at first glance one might believe
we are still well above the scales of 3-dimensional turbulence, but several
previous analysis (especially of wind spectra within the Global Assimilation
and Prognosis System) indicate $\alpha=-5/3$ regime up to several hundreds of
kilometers [10, 32]. Thus, the singular powerlaw extending over several orders
of magnitude for temperature anomalies on the earths surface is seemingly in
contrast with some atmospheric fields at higher altitudes. Nevertheless, we
emphasise that ERA5 reanalysis is not an empirical field - it is obtained
through variational data assimilation techniques, which require smoothing over
small spatial scales. Thus, the spectral power on the smallest scales remains
ill-constrained is subject to the potential systematics effects of smoothing
(see Appendix E).

###  B Power Spectrum of Average Temperature

Temperature anomalies are reported relative to the baseline average
temperature for each month from Jan. 1951 to Dec. 1980. The power spectrum for
this baseline temperature is seen in Fig. A2. Firstly, the power of all modes
is larger than observed for temperature anomalies; this is expected as
absolute temperature varies far more drastically - ranging from $-70\degree$C
to $40\degree$C with elevation, climate and latitude. For instance the
systematic variation with elevation will increase the amplitude of temperature
fluctuations on spatial scales ranging from continents to mountains. Secondly,
the shape of the power-spectrum is different than temperature anomalies, but
still consistent with a powerlaw for large degrees. The powerlaw slope
$\alpha$ is for all months within the range $-3.1$ to $-3.4$, which remains
significantly steeper than for temperature anomalies. The cause of the
increased steepness remains (as of yet) unresolved. Thirdly, the spread from
the powerlaw has increased with large systematic bias observed across months
for specific degrees (ie. $C_{3}$ is remarkably small). Such rugged power-
spectra are commonly observed when examining altitude [33], which illustrates
the apparent systematic bias in decomposing absolute temperature at the
surface.

Spherical harmonics have not previously been used as a standard descriptor of
climate change, because structures like elevation, latitude and differences
between land and ocean dominate the variation of absolute temperature.
However, by limiting the analysis to temperature anomalies we focus only on
the underlying fluctuations due to atmospheric processes. This effectively
decouples the effect of climate change from the larger irrelevant structures.

Figure A2: Power spectrum of average temperature for each month from Jan. 1951
to Dec. 1980. The power for all modes is larger than for temperature
anomalies.

###  C Power Spectrum Characteristics are robust to limiting Coverage

While the specific numerical value of the power for any mode is dependent on
the coverage, the characteristic features of the power-spectrum (discussed in
§ 3.2 and § 3.3) are unchanged when limiting the coverage to that 1850 [See
Fig. A3]. The turnover-length remains at $\approx 3000$ km, the powerlaw slope
remains consistent with $\alpha=-3$ for the entire time-series and the
relative growth of small to large modes is unchanged.

Figure A3: Power spectrum of temperature distribution from 1850s (in black) to
2010s (in red) in log-log plot limiting coverage to that of 1850. The power
for modes $l\geq 10$ is increasing with time especially from 1960 to 2020. For
comparison see the full coverage in Fig. 4.

###  D Increasing Resolution

From 1850-2020 the number of surface temperature measurements over land and
sea increased drastically. In Fig. A4 we show the evolution in number of
active temperature stations over land (reported in [16]) and number of
”superobservations” of temperatures over sea. Super-observations are the
Winsorised mean of temperature anomalies used in HadSST [17].

We emphasise that most of the evolution of the power-spectrum is from
1960-2020, where no significant change in number of either land stations nor
super-observations over sea takes place. On the other hand, the methodology
yields consistent power-spectra from 1850-1960 even though the number-of-
observations increased by orders of magnitude. Thus, increasing the number of
measurements will improve the statistical uncertainty on the monthly
temperature-estimates, but improving accuracy does not correlate with a bias
in the spherical harmonics decomposition.

Figure A4: Top: power of degree l relative to power in 1850. Bottom: Number of
active temperature stations over land (brown line) and number of
”superobservations” of temperatures over sea (blue dots) from 1850-2021. Grey
shaded region indicates the period 1960-2021. The significant growth in power-
spectrum does not correlate with an increasing number of observing stations.

###  E Interpolation Bias

Fig. A5 shows the underestimation of power when generated temperature
anomalies are interpolated over a characteristic length-scale, $L$. As
expected larger interpolation length results in a larger bias. For
interpolation over $0.5\degree$, $1\degree$ and $2\degree$ latitude and
longitude the estimated power remains accurate within respectively 1%, 3% and
10% up to degree $l=50$. Therefore, the typical bias remains significantly
smaller, than the factor 2 growth in power of these modes from 1960-2020.
Furthermore, the bias has a clear dependence on angular scale. In contrast,
from 1960-2020 all power within the cascade grew coherently maintaining the
powerlaw-slope, $\alpha=-3$ (see Fig. A6).

Over land the typical scale of interpolation from 1960-2020 is less than 100
km so any potential bias would remains negligible for modes $l\leq 90$ [16].
Note, the growth observed in the power-spectrum from 1960-2020 is also
detected in a when only examines temperature anomalies over land using a
multi-tapered power-spectrum.

The HadSST grid-cells are separated by up to 500 km, which is re-interpolated
into a $1\degree\times 1\degree$ lat-long grid [3]. Thus, interpolation
lengths are typically of order 250 km (which suggests underestimation of
around 20% for degree $l=50$). Thus, sparse sea coverage may bias the power-
spectrum on small scales at all times, but neither the magnitude nor the
functional shape of the bias matches the growing power from 1960-2020.

Figure A5: Inferred power for smoothed temperature anomalies (linearly
interpolated over 0.5°, 1°, 2°and 5 °latitude and longitude) relative to power
of unsmoothed temperature anomalies. The smoothed power is underestimated with
most suppression at small angular scales. Larger interpolation lengths result
in smoothing structure at larger scales. Figure A6: Inferred power for
smoothed temperature anomalies (linearly interpolated over $2\degree\times
2\degree$ lat-lon) relative to power of unsmoothed temperature anomalies (red
line). Observed power in 1960 to 1970 relative to power in 2010 to 2020 (blue
error-bars). The angular dependence of any interpolation bias contrast the
coherent growth of modes $20\leq l\leq 90$.

###  F Difference between Land And Sea

The analysis presented thus far has averaged over the entire globe, thereby
neglecting any kinematic or dynamical differences between land and sea.
However, if we limit the spatial coverage, we can determine the contribution
of different parts of the globe. In general, the procedure of obtaining a
multitaper power spectrum estimate is well established. [33, 19] As seen in
Fig. A7, the contribution of earth and sea temperatures are not
interchangeable with slightly different functional shapes and powerlaw slopes.
While the typical temperature fluctuations over land and sea have a similar
turnover length, the fluctuation on most length scales have larger amplitudes
on land than at sea. For instance the dipole power, $C_{1}$, is substantially
larger over land than sea - as extensively observed the seasonal difference is
larger over land.

Figure A7: Power spectrum from 2000-2020 when only calculated from Land or Sea
temperatures. Importantly, the temperature fluctuations are not consistent
with identical within 5$\sigma$. The power spectrum over land is similar to
the entire surface (indicated by grey shading) with similar powerlaw slopes,
turnover length-scale and overall power on large scales. Sea temperatures
display fewer large scale temperature fluctuations, but with a shallower
powerlaw decline than land.

We again emphasise that limiting coverage will introduce a bias for degrees
close to the spectral bandwidth of the aliased structure [18, 19]. This is a
direct consequence of the signal from low degrees leaking into high degrees
over the spectral bandwidth of the excluded regions. However, previous studies
specifically dividing land and sea have highlighted 1) that any leakage of
signal is small between these domains and 2) any bias remains small for large
degrees [19]. Thus, we conclude that the spectral differences of land and sea
is not merely caused by aliasing the spatial structure. Furthermore, Fig. A7
clearly indicates that the land-surface temperature field remains the dominant
contribution to the observed power spectrum. The overall power, the turnover
length and the powerlaw decline for the entire globe traces the power spectrum
from land.

###  G The Flow of Temperature Fluctuations

To determine whether the observed cascade is downwards (from large spatial
scales to smaller) or inverse (from small to large) we compute the Pearson
Correlation Coefficient between the power of each mode;

$r_{C_{l}(t),C_{l}(t+\Delta t)}=\frac{Cov(C_{l}(t),C_{l}(t+\Delta
t))}{\sigma_{C_{l}(t)}\sigma_{C_{l}(t+\Delta t)}}$ (1)

For $\Delta t=0$, we can determine the correlation between the power of
different modes. In Fig. A8 we see that all modes within the cascade ($l\geq
10$) show large internal correlation with $r>0.9$. The large scale spatial
structures are far more uncorrelated with $r<0.6$.

Figure A8: Correlation Coefficient between power of modes. Typically modes of
similar spatial scales are more correlated while all modes within the cascade
show large correlation. Figure A9: Correlation Coefficient between power of
modes in neighbouring months (left) and residual correlation coefficient [ie.
difference in correlation between power in prior and subsequent months]
(right). We still see that modes of similar spatial scales are more
correlated, while all modes within the cascade show large correlation.
However, the correlation between powers across subsequent months is not
symmetric. There is more correlation between lower modes in previous months
with higher modes in subsequent months then between higher modes in previous
months with lower modes in ensuing months.

For $\Delta t=1$ Month, we can determine the correlation between the power of
different modes in neighbouring months, see Fig. A9 (left). Over time
correlations in general decay, as the former temperature fluctuations diffuse
and dissipate. Nevertheless, the relative distribution of correlations remain
unchanged with for instance the highest correlation remaining within the
cascade. However, the correlation is no longer symmetrical; for instance
$C_{20}(t+\Delta t)$ is more correlated with the power of higher modes in
subsequent months, than the correlation between $C_{20}(t)$ and higher modes
of prior months. This implies, that increasing power on large spatial scales
is often followed by increased power on smaller spatial scales, but the
inverse does not hold. Thus, the direction of the cascade of temperature
oscillations is downwards.

In Fig. A9 (right) this contrast is emphasised by taking the difference in
correlation between previous and following month. All modes from $10\leq l\leq
40$ the flow is from large to smaller spatial scales (ie. the correlations
from large to small are forward in time). Ultimately, this resolution across
time yields observational evidence, that temperature oscillations follow a
downward cascade from large to small spatial scales. For modes $l<8$ the
inverse relationship seems to hold, smaller spatial scales precede
oscillations on larger scales. This tentatively suggests, that temperature
oscillations are generally driven around the characteristic transition within
the power-spectrum (ie. $l\approx 8$), which then feeds temperature
fluctuations on both smaller and larger spatial scales.

Figure A10: Left: the turnover length scale as a function of months for the
northern and southern hemisphere [averaged over the period 1950-2020]. Right:
the turnover length scale from 1950-2020 for entire globe. We see the smallest
turnover-lengths in summer months, with the larger turnover-length in winter.
Furthermore, in recent decades the turnover length has increased.

###  H The Temporal Evolution of Turnover Length

The emphasis on the coherent temporal evolution of the power within the
cascade underlines the growing weather volatility. However, another
interesting temporal trend is the changing spatial scale for the transition
from flat to cascade. We can determine the intersection between the powerlaw
(fitted from modes $10\leq l\leq 50$) with the average power of modes in the
flat regime ($2\leq l\leq 8$), which yields an estimate of the turnover
length. As seen in Fig. A10 this definition proves that the characteristic
length scale is dynamically evolving. A noteworthy feature is that the
turnover length has increased from $2960\pm 30$ km (for 1950-1960) to $3230\pm
60$ km (for 2010-2020). The discrepancy in the length scales for these two
periods is more than $3\sigma$ significant.

Figure A11: Left panels show simulated temperature field composed from 250
isotropically distributed Gaussian temperature fluctuations with right panels
indicating corresponding power spectrum for any spherical harmonic
$Y_{l}^{m}$. Top illustrates equal gaussian variance along longitudes and
latitudes, with bottom showing 5 times larger variance along longitudes than
latitudes. Evidently anisotropies of the scale of temperature fluctuations
causes an asymmetry in power across $m$.

###  I Simulated Anisotropic Temperature Fields

To simulate a global temperature field with a prescribed scale of
fluctuations, we randomly generate 250 points sampled isotropically across the
sphere. Each point represents the center of a Gaussian with the variance
setting the scale of fluctuations. Importantly, the variance along latitudes
and longitudes can be set independently. The center and spread of the $i$’th
Gaussian is denoted respectively ($\phi_{i}$, $\theta_{i}$) and
($\sigma_{\phi},\sigma_{\theta}$). Thus, the total temperature field is
determined as a superposition of all individual fluctuations:

$T(\theta,\phi)=\sum_{i}{\frac{1}{2\pi\sqrt{\sigma_{\phi}\sigma_{\theta}}}\exp{\left(-\frac{(\phi-\phi_{i})^{2}\sin(\theta)^{2}}{2\sigma_{\phi}^{2}}-\frac{(\theta-\theta_{i})^{2}}{2\sigma_{\theta}^{2}}\right)}}$
(2)

As seen in Fig. A11 varying the scale of fluctuations between longitudes and
latitudes varies the power across $m$ at fixed $l$.


# Spin entanglement in two-proton emission from 6Be

Tomohiro Oishi<EMAIL_ADDRESS>RIKEN Nishina Center for
Accelerator-Based Science, Wako 351-0198, Japan

###### Abstract

This paper presents an evaluation of coupled-spin entanglement in the two-
proton ($2p$) radioactive emission. The three-body model of 6Be with the
proton-proton interaction, which is adjusted to reproduce the experimental
energy release, is utilized. Time-dependent calculation is performed to
compute the coupled-spin state of the emitted two protons. The spin-
correlation function $S$ as Clauser-Horne-Shimony-Holt (CHSH) indicator is
obtained as $\left|S\right|\cong 2.72$. Namely, the $2p$-spin entanglement
beyond the limit of local-hidden-variable (LHV) theory is concluded. This
entanglement is sensitive to the proton-proton interaction. The short-lived
(broad-width) $2p$ state has the weaker spin entanglement. In parallel, the
core-proton interactions do not harm this entanglement during the time-
dependent decaying process. The CHSH measurement can be a probe into the
effective nuclear interaction inside finite systems.

Introduction \- Clauser-Horne-Shimony-Holt (CHSH) inequality is one essential
property of the quantum-entangled state Clauser  _et al._ (1969). This is one
variant of Bell inequality introduced by John Clauser et. al. for the proof of
Bell theorem, which claims that certain consequence of quantum mechanics
cannot be reproduced by the local-hidden-variable (LHV) theory Bell (1964,
2004). By using the CHSH indicator $S$, the limit of LHV theory is
symbolically given as $\left|S\right|\leq 2$. Throughout the history of Bell-
CHSH examinations Clauser and Shimony (1978); Aspect  _et al._ (1982); Weihs
_et al._ (1998); Pan  _et al._ (2000); Rowe  _et al._ (2001); Scheidl  _et
al._ (2010); Giustina  _et al._ (2015); Shalm  _et al._ (2015); Vasilyev  _et
al._ (2017); Storz  _et al._ (2023), the violation of LHV-theory limit has
been confirmed. In these examinations, the entangled states of photons,
electrons, and atoms are populated and measured to satisfy $\left|S\right|>2$.
I emphasize that, for this purpose, many efforts have been devoted to close
the loopholes, including those of detection Rowe  _et al._ (2001); Giustina
_et al._ (2015); Shalm  _et al._ (2015), locality Aspect  _et al._ (1982);
Weihs  _et al._ (1998); Scheidl  _et al._ (2010), and memory Barrett  _et al._
(2002).

In the nuclear physics, the quantum entanglement plays a role in various
scenes Lamehi-Rachti and Mittig (1976); Sakai  _et al._ (2006); Kanada-En’yo
(2015); Bulgac (2023); Pazy (2023); Gu  _et al._ (2023); Hengstenberg  _et
al._ (2023); Pérez-Obiol  _et al._ (2023); Kou  _et al._ (2024); Kirchner  _et
al._ (2024); Robin  _et al._ (2021); Kruppa  _et al._ (2021); Miller (2023a,
b); Tichai  _et al._ (2023); Sun  _et al._ (2023). In Ref. Lamehi-Rachti and
Mittig (1976) by Lamehi-Rachti and Mittig, the low-energy proton-proton
scattering was measured for testing the Bell inequality. In this pioneer work,
some extra assumptions were required to certify the violation of LHV-theory
limit. In Ref. Sakai  _et al._ (2006), Sakai et. al. performed a novel
measurement to demonstrate that a strong entanglement is realized between two
protons. This experiment measures the spin-singlet two protons made in the
reaction of 2H+$p\longrightarrow^{2}$He+$n$. The spin-correlation function as
the CHSH indicator is deduced as $S_{\rm expt}=2.83\pm 0.24_{\rm stat}\pm
0.07_{\rm sys}$, which is in agreement with the non-local quantum mechanics
and beyond the LHV-theory limit. Recently, several theoretical works have been
devoted to compute the entanglement entropy in atomic nuclei Kanada-En’yo
(2015); Bulgac (2023); Pazy (2023); Gu  _et al._ (2023); Hengstenberg  _et
al._ (2023); Pérez-Obiol  _et al._ (2023); Kou  _et al._ (2024); Kirchner  _et
al._ (2024). In Ref. Bulgac (2023), due to the nuclear short-range
correlations, the occupation probabilities of nuclear orbits change to
increase the entanglement entropy. In Ref. Bai (2024), the evaluation of
nuclear spin entanglement with the quantum-state tomography is suggested to be
feasible. In Refs. Johnson and Gorton (2023); Miller (2023b), a finite proton-
neutron entanglement is suggested.

Another possible example to observe the nuclear quantum entanglement is the
two-proton ($2p$) radioactive emission Bertulani  _et al._ (2008); Grigorenko
(2009); Pfützner  _et al._ (2012); Blank and Ploszajczak (2008); Blank and
Borge (2008); Pfützner  _et al._ (2023); Qi  _et al._ (2019). In this
radioactive decay, the parent nuclei spontaneously decay by emitting two
protons. Especially in so-called “prompt” $2p$ emission Grigorenko (2009);
Goldansky (1960, 1961), the two protons are expected to have the diproton-like
clustering and/or the dominant spin-singlet configuration. This is
attributable to the effective $2p$ interaction inside finite systems, being in
a contrast to the vacuum $2p$ interaction supporting no bound state.

In this paper, I evaluate the quantum entaglement in the $2p$ emission of 6Be.
Main motivation is to utilize this entanglement as a probe into the effective
interaction. The 6Be is the lightest $2p$ emitter and well approximated with
the simple three-body model with time dependence Oishi  _et al._ (2014, 2017);
Wang and Nazarewicz (2021). For the two protons spontaneously emitted, the
measurement of their coupled-spin correlation is assumed Bertulani  _et al._
(2008). For quantitizing the entanglement, the CHSH indicator $S$, which was
originally introduced for testing the CHSH inequality Clauser  _et al._
(1969); Sakai  _et al._ (2006), is computed. Because of the three-body
problem, the $2p$ entanglement is under the effect from the third particle,
namely, the daughter alpha nucleus. Whether this effect destroys or not the
entanglement is investigated. The sensitivity of CHSH indicator to the proton-
proton interaction is also discussed.

Formalism and Model \- I consider the coupling of two protons, i.e., identical
spin-$1/2$ fermions. In such a case, the CHSH indicator is represented as
follows. First the four options of measurement by the two observers, so-called
“Alice” and “Bob” by convenance, are introduced:

$\displaystyle\hat{A}_{1,\theta}(1)\otimes\hat{B}_{1,\theta}(2),~{}\hat{A}_{2,\theta}(1)\otimes\hat{B}_{1,\theta}(2),$
$\displaystyle\hat{A}_{1,\theta}(1)\otimes\hat{B}_{2,\theta}(2),~{}\hat{A}_{2,\theta}(1)\otimes\hat{B}_{2,\theta}(2),$
(1)

Namely, Alice observes the first fermion with one choosen from the two
options, $\hat{A}_{1,\theta}$ and $\hat{A}_{2,\theta}$. Bob does the second
fermion with one choosen from $\hat{B}_{1,\theta}$ and $\hat{B}_{2,\theta}$.
Those operators including the parameter angle $\theta$ are given as

$\hat{A}_{1}=\hat{\sigma}_{z},~{}~{}~{}~{}\hat{A}_{2,\theta}=\hat{\sigma}_{z}\cos
2\theta+\hat{\sigma}_{x}\sin 2\theta,$ (2)

for Alice, whereas

$\begin{array}[]{l}\hat{B}_{1,\theta}=\hat{\sigma}_{z}\cos\theta+\hat{\sigma}_{x}\sin\theta,\\\
\hat{B}_{2,\theta}=\hat{\sigma}_{z}\cos\theta-\hat{\sigma}_{x}\sin\theta,\end{array}$
(3)

for Bob. For an arbitrary two-fermion state $\ket{\Psi(1,2)}$, their
expectation values are obtained:

$\Braket{A_{i}B_{j}}=\Braket{\Psi(1,2)\mid\hat{A}_{i,\theta}(1)\otimes\hat{B}_{j,\theta}(2)\mid\Psi(1,2)}.$
(4)

Then the CHSH indicator $S$ is determined as

$S=\max\left\\{\left|S_{-+++}\right|,\left|S_{+-++}\right|,\left|S_{++-+}\right|,\left|S_{+++-}\right|\right\\},$
(5)

where

$\displaystyle S_{-+++}$ $\displaystyle=$
$\displaystyle-\Braket{A_{1}B_{1}}+\Braket{A_{2}B_{1}}+\Braket{A_{1}B_{2}}+\Braket{A_{2}B_{2}},$
$\displaystyle S_{+-++}$ $\displaystyle=$
$\displaystyle\Braket{A_{1}B_{1}}-\Braket{A_{2}B_{1}}+\Braket{A_{1}B_{2}}+\Braket{A_{2}B_{2}},$
$\displaystyle S_{++-+}$ $\displaystyle=$
$\displaystyle\Braket{A_{1}B_{1}}+\Braket{A_{2}B_{1}}-\Braket{A_{1}B_{2}}+\Braket{A_{2}B_{2}},$
$\displaystyle S_{+++-}$ $\displaystyle=$
$\displaystyle\Braket{A_{1}B_{1}}+\Braket{A_{2}B_{1}}+\Braket{A_{1}B_{2}}-\Braket{A_{2}B_{2}}.$
(6)

In FIG. 1, the CHSH indicator for the spin-singlet state,
$\ket{\Psi(1,2)}=\ket{\uparrow\downarrow-\downarrow\uparrow}/\sqrt{2}$, is
presented. At $\theta=\pi/4,~{}3\pi/4,~{}5\pi/4$, and $7\pi/4$, this indicator
has the maximum value, $S=2\sqrt{2}$ (Tsirelson’s bound) Cirel’son (1980).
This situation resembles Sakai’s experiment Sakai  _et al._ (2006). In the
following sections, $\theta=\pi/4$, except when modified.

FIG. 1: CHSH indicator of the pure spin-singlet state.

As well known, when the initial state is one of the Bell states,
$\ket{\Psi(1,2)}=\ket{e_{n}}$, the CHSH indicator satisfies $S=2\sqrt{2}$ when
$\theta=\pi/4$. Here the Bell states read Nielsen and Chuang (2000)

$\displaystyle\ket{e_{1}}=\frac{1}{\sqrt{2}}\ket{\uparrow\uparrow+\downarrow\downarrow},~{}~{}~{}~{}~{}\ket{e_{2}}=\frac{i}{\sqrt{2}}\ket{\uparrow\uparrow-\downarrow\downarrow},$
$\displaystyle\ket{e_{3}}=\frac{i}{\sqrt{2}}\ket{\uparrow\downarrow+\downarrow\uparrow},~{}~{}~{}~{}~{}\ket{e_{4}}=\frac{1}{\sqrt{2}}\ket{\uparrow\downarrow-\downarrow\uparrow}.$
(7)

These states can be used as basis to represent an arbitrary coupling of two
spin-$1/2$ fermions. In addition, there can be also coordinate degrees of
freedom. Thus, an arbitrary two-fermion state is generally expanded as

$\ket{\Psi(1,2)}=\sum_{i=1}^{4}G_{i}(\bm{r}_{1},\bm{r}_{2})\ket{e_{i}},$ (8)

where $G_{i}(\bm{r}_{1},\bm{r}_{2})$ is the coordinate part.

For the coupled spin $\hat{S}_{12}$, its eigenvalues are given as
$\hat{S}^{2}_{12}\ket{d_{k}}=S(S+1)\ket{d_{k}}$ and
$\hat{S}_{12,z}\ket{d_{k}}=V\ket{d_{k}}$, where

$\displaystyle\ket{d_{1}}$ $\displaystyle=$
$\displaystyle\ket{S=1,V=+1}=\ket{\uparrow\uparrow},$
$\displaystyle\ket{d_{2}}$ $\displaystyle=$
$\displaystyle\ket{S=1,V=-1}=\ket{\downarrow\downarrow},$
$\displaystyle\ket{d_{3}}$ $\displaystyle=$
$\displaystyle\ket{S=1,V=0}=\frac{1}{\sqrt{2}}\ket{\uparrow\downarrow+\downarrow\uparrow},$
$\displaystyle\ket{d_{4}}$ $\displaystyle=$
$\displaystyle\ket{S=0,V=0}=\frac{1}{\sqrt{2}}\ket{\uparrow\downarrow-\downarrow\uparrow}.$
(9)

Then, the unitary transformation $U$ is formulated:

$\begin{array}[]{l}\ket{e_{1}}\\\ \ket{e_{2}}\\\ \ket{e_{3}}\\\
\ket{e_{4}}\end{array}=U\begin{array}[]{l}\ket{d_{1}}\\\ \ket{d_{2}}\\\
\ket{d_{3}}\\\
\ket{d_{4}}\end{array}=\left(\begin{array}[]{llll}\frac{1}{\sqrt{2}}&\frac{1}{\sqrt{2}}&0&0\\\
\frac{i}{\sqrt{2}}&\frac{-i}{\sqrt{2}}&0&0\\\ 0&0&i&0\\\
0&0&0&1\end{array}\right)\begin{array}[]{l}\ket{d_{1}}\\\ \ket{d_{2}}\\\
\ket{d_{3}}\\\ \ket{d_{4}}\end{array}.$ (10)

In the following sections, the total-spin basis is also utilized:
$\ket{\Psi(1,2)}=\sum_{i}F_{i}(\bm{r}_{1},\bm{r}_{2})\ket{d_{i}}$. Thus,
$G_{i}(\bm{r}_{1},\bm{r}_{2})=\sum_{k}\left(U^{\dagger}\right)_{ik}F_{k}(\bm{r}_{1},\bm{r}_{2})$
to convert to Bell basis. Notice that, for computing the expectation values in
Eq. (6), these coordinate parts can be simply integrated.

I employ the three-body model, which has been developed and utilized in Refs.
Suzuki and Ikeda (1988); Bertsch and Esbensen (1991); Esbensen  _et al._
(1997); Hagino and Sagawa (2005); Oishi  _et al._ (2014, 2017). The system
contains an alpha particle as the rigid core with mass $m_{C}$ and two valence
protons. The two valence protons feel the spherical mean field $V(\bm{r})$
generated by the alpha core. Thus, the three-body Hamiltonian reads

$\displaystyle\hat{H}_{3B}=\hat{h}(\bm{r}_{1})+\hat{h}(\bm{r}_{2})+v_{pp}(\bm{r}_{1},\bm{r}_{2})+\frac{\bm{p}_{1}\cdot\bm{p}_{2}}{m_{C}},$
(11)
$\displaystyle\hat{h}(r_{i})=-\frac{\hbar^{2}}{2\mu}\frac{d^{2}}{dr_{i}^{2}}+\frac{\hbar^{2}}{2\mu}\frac{l(l+1)}{r_{i}^{2}}+V(r_{i}),$
(12)

where $\mu=m_{p}m_{C}/(m_{p}+m_{C})$ and $m_{p}=938.272$ MeV$/c^{2}$ for
protons. Here $\hat{h}(\bm{r}_{k})$ is of the $k$th proton-alpha subsystem.
For its interaction $V(r_{i})$, the same Woods-Saxon and Coulomb potentials in
the previous work Oishi  _et al._ (2017) are employed. The proton-proton
interaction reads

$v_{pp}(\bm{r}_{1},\bm{r}_{2})=v_{nucl}(d)+\frac{e^{2}}{4\pi\epsilon_{0}}\frac{1}{d}+v_{add}(\bm{r}_{1},\bm{r}_{2}),$
(13)

where $d=\left|\bm{r}_{2}-\bm{r}_{1}\right|$. The nuclear-force term is
described by the spin-dependent Gaussian potential Thompson  _et al._ (1977):

$\displaystyle v_{nucl}(d)$ $\displaystyle=$
$\displaystyle\left[V_{R}e^{-a_{R}d^{2}}+V_{S}e^{-a_{S}d^{2}}\right]\hat{P}_{S=0}$
(14)
$\displaystyle+\left[V_{R}e^{-a_{R}d^{2}}+V_{T}e^{-a_{T}d^{2}}\right]\hat{P}_{S=1}.$

The operator $\hat{P}_{S=0}$ ($\hat{P}_{S=1}$) indicates the projection into
the spin-singlet (spin-triplet) channel of the proton-proton subsystem.
Parameters are given as $V_{R}=200$ MeV, $V_{S}=-91.85$ MeV, $V_{T}=-178$ MeV,
$a_{R}=1,487$ fm-2, $a_{S}=0.465$ fm-2, and $a_{T}=0.639$ fm-2. These
parameters correctly reproduce the experimental vacuum-scattering length of
two protons Thompson  _et al._ (1977). In addition, the surface-dependent term
$v_{add}$ is employed:

$v_{add}(\bm{r}_{1},\bm{r}_{2})=w_{0}e^{-(R-R_{0})^{2}/B_{0}^{2}}\delta(\bm{r}_{1}-\bm{r}_{2}),$
(15)

where $R=\left|(\bm{r}_{1}+\bm{r}_{2})/2\right|$, $R_{0}=1.68$ fm,
$B_{0}=0.6R_{0}$, and $w_{0}=-470$ MeV$\cdot$fm3. This additional term is
necessary to reproduce the experimental energy: the $2p$ energy and width are
obtained as $E_{2p}=1.356$ MeV and $\Gamma_{2p}=0.055$ MeV, respectively,
whereas the experimental data read $E_{2p}=1.372(5)$ MeV and
$\Gamma_{2p}=0.092(6)$ MeV Ajzenberg-Selove (1988, 1991). Notice that this
additional potential vanishes when one of the three particles is infinitely
separated. Thus, the vacuum properties of two-body subsystems can be conserved
in the time-development calculations. I especially notify that the emitted two
protons should be unbound in vacuum.

FIG. 2:  Time-dependent decaying probability of $2p$ emission from 6Be. Its
spin-singlet ($S_{12}=0$) component is also plotted. The factor $f_{pp}=1$
indicates that the proton-proton interaction is adjusted to reproduce the
experimental $2p$-energy.

Results and Discussions \- The $2p$-emitting process is simulated with the
time-dependent method Oishi  _et al._ (2014, 2017):

$\ket{\Psi(t)}=\exp\left[-it\frac{\hat{H}_{3B}}{\hbar}\right]\ket{\Psi(0)},$
(16)

where the initial state is solved as the confined $2p$ state inside the
Coulomb barrier Oishi  _et al._ (2014). The decaying state
$\ket{\Psi_{d}(t)}$, which describes the emitted component outside the
barrier, is determined as

$\ket{\Psi_{d}(t)}=\ket{\Psi(t)}-\beta(t)\ket{\Psi(0)},$ (17)

where $\beta(t)$ is the survival coefficient,
$\beta(t)=\Braket{\Psi(0)}{\Psi(t)}$.

In FIG 2, the time-dependent decaying probability is displayed:

$P_{decay}(t)=\Braket{\Psi_{d}(t)\mid\Psi_{d}(t)}=1-\left|\beta(t)\right|^{2}.$
(18)

One clearly finds that the decaying probability increases in time development.
The spin-singlet state,
$\ket{d_{4}}=\ket{e_{4}}=\ket{\uparrow\downarrow-\downarrow\uparrow}/\sqrt{2}$,
is always dominant. I confirmed that the survival probability,

$P_{surv}(t)\equiv\Braket{\Psi(0)\mid\Psi(t)}^{2}=\left|\beta(t)\right|^{2}=1-P_{decay}(t),$
(19)

is well approximated by the exponential damping. Namely, $P_{surv}(t)\cong
e^{-t/\tau}$, where the width can be evaluated as $\Gamma_{2p}=\hbar/\tau\cong
0.055$ MeV by numerical fitting. Also, by observing the decaying density
distribution,
$\rho_{d}(t,\bm{r}_{1},\bm{r}_{2})=\left|\Psi_{d}(t,\bm{r}_{1},\bm{r}_{2})\right|^{2}$,
I confirmed that the present $2p$ emission can be interpreted as the diproton-
correlating emission Oishi  _et al._ (2014, 2017).

FIG. 3:  CHSH indicator of the time-dependent decaying state
$\ket{\Psi_{d}(t)}$ of the $2p$ emission of 6Be $\longrightarrow\alpha+p+p$.
Here $\theta=\pi/4$. The factor $f_{pp}$ indicates the strength of the proton-
proton interaction.

In Fig. 3, the CHSH indicator is evaluated for the decaying state by using the
expansion,

$\ket{\Psi_{d}(t,\bm{r}_{1},\bm{r}_{2})}=\sum_{i=1}^{4}G_{i}(t,\bm{r}_{1},\bm{r}_{2}))\ket{e_{i}}.$
(20)

One clearly finds that the CHSH indicator becomes larger than $2$, i.e. the
LHV-theory limit, during the time evolution. For $ct\geq 1500$ fm, I obtain
$S\cong 2.72$ with the default setting of $v_{pp}$ ($f_{pp}=1.0$) for
reproducing $E_{2p}=1.356$ MeV of 6Be Ajzenberg-Selove (1988, 1991). This
coupled-spin entanglement is of course attributable to the dominant spin-
singlet component in FIG 2.

Notice that the present time evolution in Eq. (16) also contains the effect of
core-proton interactions $V(r_{i})$. These interactions do not harm the $2p$
entanglement.

For deeper investigation, the sensitivity of CHSH indicator to the proton-
proton interaction is studied. For this purpose, an effective tuning factor is
employed:

$v_{pp}(\bm{r}_{1},\bm{r}_{2})\longrightarrow
f_{pp}v_{pp}(\bm{r}_{1},\bm{r}_{2}).$ (21)

With $f_{pp}=1$, the $v_{pp}$ approximates the experimental energy,
$E_{2p}=1.356$ MeV Ajzenberg-Selove (1988, 1991). From FIG. 3, one can read
that this CHSH indicator becomes smaller when the interaction is weakened. In
the present case of 6Be, it goes below the LHV-theory limit, $S<2$, with
$f_{pp}\leq 0.2$. Consequently, the present entanglement is a product of the
proton-proton interaction. With the $v_{pp}$ weakened, the $2p$ energy and
width also change: $E_{2p}=2.864$ MeV and $\Gamma_{2p}=0.8579$ MeV with
$f_{pp}=0.5$; $E_{2p}=2.925$ MeV and $\Gamma_{2p}\cong 2.9$ MeV with
$f_{pp}=0.2$. Thus, the short-lived (broad-width) $2p$ state has the weaker
spin entanglement.

Summary \- In this paper, the CHSH indicator for the $2p$ emission from 6Be is
evaluated in order to measure its coupled-spin entanglement. As a product of
proton-proton interaction, $S\cong 2.72$ is obtained in the time-dependent
decaying state: an entanglement beyond the limit of LHV theory is suggested.
This entanglement is not harmed by the core-proton interactions. The CHSH
measurement can be a probe into the effective nuclear interaction inside
finite systems.

This work is limited to the lightest 6Be nucleus. Whether the spin
entanglement exists commonly in other $2p$ emitters or not is an open
question. The dominance of spin-singlet Bell state is not trivial in other
systems. For evaluating the quantum entanglement, the von-Neumann entropy has
been considered as an essential quantity Kanada-En’yo (2015); Bulgac (2023);
Pazy (2023); Gu  _et al._ (2023); Hengstenberg  _et al._ (2023); Pérez-Obiol
_et al._ (2023); Kou  _et al._ (2024); Kirchner  _et al._ (2024), whereas this
work focuses only on the CHSH indicator. Evaluation and discussion of the von-
Neumann entropy for the time-dependent state is in progress. Whether the
entanglement is apparent or not may depend on the observable quantity of
interest, e.g. the coupled spin, momenta, and energy distribution. These
topics will be addressed in forthcoming studies.

In the experimental side, mass production of $2p$ emitters, including the 6Be,
for a sufficient statistics is still challenging. For measuring the two
protons by independent detectors after decay, an advanced design of experiment
will be necessary. Closing all the loopholes of detection Rowe  _et al._
(2001); Giustina  _et al._ (2015); Shalm  _et al._ (2015), locality Aspect
_et al._ (1982); Weihs  _et al._ (1998); Scheidl  _et al._ (2010), and memory
Barrett  _et al._ (2002) should require another lot of efforts. On the other
side, the experimental survey of $2p$ emitters is rapidly in progress, e.g.
that in RIKEN RIBF. Confirmation of spin entanglement in $2p$ radioactivity
can be one landmark in the nuclear physics as well as quantum many-body
science.

Acknowledgments \- I especially thank Yutaka Shikano, Masaki Sasano, Tomoya
Naito, Tokuro Fukui, and Masaaki Kimura for fruitful discussions. Numerical
calculations are supported by the cooperative project of supercomputer
Yukawa-21 in Yukawa Institute for Theoretical Physics, Kyoto University. I
appreciate the Multi-disciplinary Cooperative Research Program (MCRP)
2023-2024 by Center for Computational Sciences, University of Tsukuba (project
ID wo23i034), allocating computational resources of supercomputer
Wisteria/BDEC-01 (Odyssey) in Information Technology Center, University of
Tokyo.

## References

  * Clauser  _et al._ (1969) J. F. Clauser, M. A. Horne, A. Shimony, and R. A. Holt, Phys. Rev. Lett. 23, 880 (1969).
  * Bell (1964) J. Bell, Physics 1, 195 (1964).
  * Bell (2004) J. S. Bell,  _Speakable and Unspeakable in Quantum Mechanics_ , 2nd ed., Collected Papers on Quantum Philosophy (Cambridge University Press, Cambridge, UK, 2004).
  * Clauser and Shimony (1978) J. F. Clauser and A. Shimony, Reports on Progress in Physics 41, 1881 (1978).
  * Aspect  _et al._ (1982) A. Aspect, P. Grangier, and G. Roger, Phys. Rev. Lett. 49, 91 (1982).
  * Weihs  _et al._ (1998) G. Weihs, T. Jennewein, C. Simon, H. Weinfurter, and A. Zeilinger, Phys. Rev. Lett. 81, 5039 (1998).
  * Pan  _et al._ (2000) J.-W. Pan, D. Bouwmeester, M. Daniell, H. Weinfurter, and A. Zeilinger, Nature 403, 515 (2000).
  * Rowe  _et al._ (2001) M. A. Rowe, D. Kielpinski, V. Meyer, C. A. Sackett, W. M. Itano, C. Monroe, and D. J. Wineland, Nature 409, 791 (2001).
  * Scheidl  _et al._ (2010) T. Scheidl, R. Ursin, J. Kofler, S. Ramelow, X.-S. Ma, T. Herbst, L. Ratschbacher, A. Fedrizzi, N. K. Langford, T. Jennewein, and A. Zeilinger, Proceedings of the National Academy of Sciences 107, 19708 (2010), https://www.pnas.org/doi/pdf/10.1073/pnas.1002780107 .
  * Giustina  _et al._ (2015) M. Giustina, M. A. M. Versteegh, S. Wengerowsky, J. Handsteiner, A. Hochrainer, K. Phelan, F. Steinlechner, J. Kofler, J.-A. Larsson, C. Abellán, W. Amaya, V. Pruneri, M. W. Mitchell, J. Beyer, T. Gerrits, A. E. Lita, L. K. Shalm, S. W. Nam, T. Scheidl, R. Ursin, B. Wittmann, and A. Zeilinger, Phys. Rev. Lett. 115, 250401 (2015).
  * Shalm  _et al._ (2015) L. K. Shalm, E. Meyer-Scott, B. G. Christensen, P. Bierhorst, M. A. Wayne, M. J. Stevens, T. Gerrits, S. Glancy, D. R. Hamel, M. S. Allman, K. J. Coakley, S. D. Dyer, C. Hodge, A. E. Lita, V. B. Verma, C. Lambrocco, E. Tortorici, A. L. Migdall, Y. Zhang, D. R. Kumor, W. H. Farr, F. Marsili, M. D. Shaw, J. A. Stern, C. Abellán, W. Amaya, V. Pruneri, T. Jennewein, M. W. Mitchell, P. G. Kwiat, J. C. Bienfang, R. P. Mirin, E. Knill, and S. W. Nam, Phys. Rev. Lett. 115, 250402 (2015).
  * Vasilyev  _et al._ (2017) D. Vasilyev, F. O. Schumann, F. Giebels, H. Gollisch, J. Kirschner, and R. Feder, Phys. Rev. B 95, 115134 (2017).
  * Storz  _et al._ (2023) S. Storz, J. Schär, A. Kulikov, P. Magnard, P. Kurpiers, J. Lütolf, T. Walter, A. Copetudo, K. Reuer, A. Akin, J.-C. Besse, M. Gabureac, G. J. Norris, A. Rosario, F. Martin, J. Martinez, W. Amaya, M. W. Mitchell, C. Abellan, J.-D. Bancal, N. Sangouard, B. Royer, A. Blais, and A. Wallraff, Nature 617, 265 (2023).
  * Barrett  _et al._ (2002) J. Barrett, D. Collins, L. Hardy, A. Kent, and S. Popescu, Phys. Rev. A 66, 042111 (2002).
  * Lamehi-Rachti and Mittig (1976) M. Lamehi-Rachti and W. Mittig, Phys. Rev. D 14, 2543 (1976).
  * Sakai  _et al._ (2006) H. Sakai, T. Saito, T. Ikeda, K. Itoh, T. Kawabata, H. Kuboki, Y. Maeda, N. Matsui, C. Rangacharyulu, M. Sasano, Y. Satou, K. Sekiguchi, K. Suda, A. Tamii, T. Uesaka, and K. Yako, Phys. Rev. Lett. 97, 150405 (2006).
  * Kanada-En’yo (2015) Y. Kanada-En’yo, Progress of Theoretical and Experimental Physics 2015, 043D04 (2015).
  * Bulgac (2023) A. Bulgac, Phys. Rev. C 107, L061602 (2023).
  * Pazy (2023) E. Pazy, Phys. Rev. C 107, 054308 (2023).
  * Gu  _et al._ (2023) C. Gu, Z. H. Sun, G. Hagen, and T. Papenbrock, Phys. Rev. C 108, 054309 (2023).
  * Hengstenberg  _et al._ (2023) S. M. Hengstenberg, C. E. P. Robin, and M. J. Savage, The European Physical Journal A 59, 231 (2023).
  * Pérez-Obiol  _et al._ (2023) A. Pérez-Obiol, S. Masot-Llima, A. M. Romero, J. Menéndez, A. Rios, A. Garcia-Sáez, and B. Juliá-Diaz, The European Physical Journal A 59, 240 (2023).
  * Kou  _et al._ (2024) W. Kou, J. Chen, and X. Chen, Physics Letters B 849, 138453 (2024).
  * Kirchner  _et al._ (2024) T. Kirchner, W. Elkamhawy, and H.-W. Hammer, Few-Body Systems 65, 29 (2024).
  * Robin  _et al._ (2021) C. Robin, M. J. Savage, and N. Pillet, Phys. Rev. C 103, 034325 (2021).
  * Kruppa  _et al._ (2021) A. T. Kruppa, J. Kovács, P. Salamon, and O. Legeza, Journal of Physics G: Nuclear and Particle Physics 48, 025107 (2021).
  * Miller (2023a) G. A. Miller, Phys. Rev. C 108, L041601 (2023a).
  * Miller (2023b) G. A. Miller, Phys. Rev. C 108, L031002 (2023b).
  * Tichai  _et al._ (2023) A. Tichai, S. Knecht, A. Kruppa, O. Legeza, C. Moca, A. Schwenk, M. Werner, and G. Zarand, Physics Letters B 845, 138139 (2023).
  * Sun  _et al._ (2023) Z. H. Sun, G. Hagen, and T. Papenbrock, Phys. Rev. C 108, 014307 (2023).
  * Bai (2024) D. Bai, Phys. Rev. C 109, 034001 (2024).
  * Johnson and Gorton (2023) C. W. Johnson and O. C. Gorton, Journal of Physics G: Nuclear and Particle Physics 50, 045110 (2023).
  * Bertulani  _et al._ (2008) C. Bertulani, M. Hussein, and G. Verde, Physics Letters B 666, 86 (2008).
  * Grigorenko (2009) L. V. Grigorenko, Physics of Particles and Nuclei 40, 674 (2009).
  * Pfützner  _et al._ (2012) M. Pfützner, M. Karny, L. V. Grigorenko, and K. Riisager, Rev. Mod. Phys. 84, 567 (2012).
  * Blank and Ploszajczak (2008) B. Blank and M. Ploszajczak, Reports on Progress in Physics 71, 046301 (2008).
  * Blank and Borge (2008) B. Blank and M. Borge, Progress in Particle and Nuclear Physics 60, 403 (2008).
  * Pfützner  _et al._ (2023) M. Pfützner, I. Mukha, and S. Wang, Progress in Particle and Nuclear Physics 132, 104050 (2023).
  * Qi  _et al._ (2019) C. Qi, R. Liotta, and R. Wyss, Progress in Particle and Nuclear Physics 105, 214 (2019).
  * Goldansky (1960) V. I. Goldansky, Nucl. Phys. 19, 482 (1960).
  * Goldansky (1961) V. I. Goldansky, Nucl. Phys. 27, 648 (1961).
  * Oishi  _et al._ (2014) T. Oishi, K. Hagino, and H. Sagawa, Phys. Rev. C 90, 034303 (2014).
  * Oishi  _et al._ (2017) T. Oishi, M. Kortelainen, and A. Pastore, Phys. Rev. C 96, 044327 (2017).
  * Wang and Nazarewicz (2021) S. M. Wang and W. Nazarewicz, Phys. Rev. Lett. 126, 142501 (2021).
  * Cirel’son (1980) B. S. Cirel’son, Letters in Mathematical Physics 4, 93 (1980).
  * Nielsen and Chuang (2000) M. Nielsen and I. Chuang,  _Quantum Computation and Quantum Information_ (Cambridge University Press, Cambridge, UK, 2000).
  * Suzuki and Ikeda (1988) Y. Suzuki and K. Ikeda, Phys. Rev. C 38, 410 (1988).
  * Bertsch and Esbensen (1991) G. Bertsch and H. Esbensen, Annals of Physics 209, 327 (1991).
  * Esbensen  _et al._ (1997) H. Esbensen, G. F. Bertsch, and K. Hencken, Phys. Rev. C 56, 3054 (1997).
  * Hagino and Sagawa (2005) K. Hagino and H. Sagawa, Phys. Rev. C 72, 044321 (2005).
  * Thompson  _et al._ (1977) D. Thompson, M. Lemere, and Y. Tang, Nuclear Physics A 286, 53 (1977).
  * Ajzenberg-Selove (1988) F. Ajzenberg-Selove, Nuclear Physics A 490, 1 (1988), note: several versions with the same title has been published.
  * Ajzenberg-Selove (1991) F. Ajzenberg-Selove, Nuclear Physics A 523, 1 (1991), with revised version at http://www.tunl.duke.edu/nucldata/.