# Chiral spin symmetry and hot QCD ††thanks: Presented at Excited QCD 2024;
Benasque, 14 - 20.01.2024

L. Ya. Glozman
Institute of Physics, University of Graz, A-8010 Graz, Austria

###### Abstract

In this talk we overview main results indicating existence in QCD of three
qualitatively different regimes connected by smooth crossovers upon heating: a
hadron gas, a stringy fluid and a quark-gluon plasma. In the combined large
$N_{c}$ and chiral limit these regimes likely become distinct phases separated
by phase transitions: a chiral restoration phase transition around $T_{ch}\sim
130$ MeV and a deconfinement phase transition around $T_{d}\sim 300$ MeV. It
should be an important task to verify this issue on the lattice. We will
introduce a chiral spin symmetry, which is a symmetry of the electric part of
electrodynamics and of QCD with light quarks. It is realized approximately in
QCD above the chiral restoration crossover and disappears in the QGP regime.
The center symmetry of the pure glue action and the chiral spin symmetry of
the electric part of the QCD Lagrangian with light quarks are complementary to
distinguish the confining regime and its disappearance. We also address other
lattice evidences for stringy fluid: hadron resonances extracted from the
lattice correlators; breakdown of the thermal perturbation theory at
$T<~{}600$ MeV and fluctuations of conserved charges that point out the
$N_{c}$ scaling above $T\sim 155$ MeV.

##  1 Chiral spin symmetry

Consider Maxwell equations that describe evolution of the electric and
magnetic fields in a given reference frame:

$\displaystyle div\vec{E}$ $\displaystyle=$ $\displaystyle 4\pi\rho$
$\displaystyle rot\vec{B}-\frac{1}{c}\frac{\partial\vec{E}}{\partial t}$
$\displaystyle=$ $\displaystyle\frac{4\pi}{c}\vec{j}$ $\displaystyle
rot\vec{E}+\frac{1}{c}\frac{\partial\vec{B}}{\partial t}$ $\displaystyle=$
$\displaystyle 0$ $\displaystyle div\vec{B}$ $\displaystyle=$ $\displaystyle
0.$ (1)

We define $\vec{E}$ and $\vec{B}$ in a given Lorentz frame in a gauge-
invariant way through their action on charge and current,

$\vec{F}=q\vec{E}+q\frac{\vec{v}}{c}\times\vec{B}.$ (2)

It is possible to measure directly $\vec{F}$ in electrodynamics but not
possible in quantum chromodynamics. Is there another method to distinguish
$\vec{E}$ and $\vec{B}$?

Consider charges to be massless particles with $s=1/2$. They can be either
right- or left-handed:

$\left(\begin{array}[]{c}R\\\ L\end{array}\right).$

Consider a $SU(2)_{CS}$ chiral spin transformation that mixes $R$ and $L$:

$\left(\begin{array}[]{c}R\\\
L\end{array}\right)\;\rightarrow\left(\begin{array}[]{c}R^{\prime}\\\
L^{\prime}\end{array}\right)=\exp\left(i\frac{\varepsilon^{n}\sigma^{n}}{2}\right)\left(\begin{array}[]{c}R\\\
L\end{array}\right).$ (3)

What happens with the charge density $\rho$?

$R^{\prime\dagger}R^{\prime}+L^{\prime\dagger}L^{\prime}=R^{\dagger}R+L^{\dagger}L,$
(4)

i.e.

$\rho^{\prime}=\rho.$ (5)

The charge density is invariant under the chiral spin transformation. However,
upon the chiral spin transformation the current density $\vec{j}$ and
$\vec{v}$ change. We conclude that the interaction of a charge with the
electric field is invariant under the chiral spin transformation, while
interaction of a current with the magnetic field is not. This means that the
chiral spin symmetry allows one to distinguish the electric and magnetic
fields in a given reference frame. The electric part of the EM theory is more
symmetric than the magnetic part.

The chromoelectric field of QCD is defined via interaction with the color
charge

$\vec{F}=Q^{a}\vec{E}^{a};~{}~{}~{}~{}Q^{a}=\int
d^{3}x~{}q^{\dagger}(x)T^{a}q(x),~{}~{}~{}~{}a=1,...,8.$ (6)

The color charge, which is Lorentz-invariant, is invariant under $SU(2)_{CS}$:

$~{}~{}~{}q\rightarrow
q^{\prime}=\exp\left(i\frac{\varepsilon^{n}\Sigma^{n}}{2}\right)q,\;~{}~{}~{}~{}\Sigma=\\{\gamma_{k},-i\gamma_{5}\gamma_{k},\gamma_{5}\\}.$
(7)

In a given Lorentz frame interaction of quarks with the electric part of the
gluonic field is chiral spin invariant like in electrodynamics. The
$SU(2)_{CS}$ symmetry includes $U(1)_{A}$ as a subgroup.

We can extend the $SU(2)_{CS}$ symmetry to $SU(2N_{F})$: $SU(2)_{CS}\times
SU(N_{F})\subset SU(2N_{F})$ and $SU(N_{F})_{R}\times SU(N_{F})_{L}\times
U(1)_{A}\subset SU(2N_{F})$. The $SU(2N_{F})$ is also a symmetry of the color
charge and of the electric part of the QCD Lagrangian. I.e., the color charge
and electric part of the theory have a $SU(2N_{F})$ symmetry that is larger
than the chiral symmetry of QCD as a whole. The fundamental vector of
$SU(2N_{F})$ at $N_{F}=2$ is

$\Psi^{T}=(u_{R}~{}u_{L}~{}d_{R}~{}d_{L}).$ (8)

Notice that the Dirac Lagrangian is not invariant under $SU(2)_{CS}$ and
$SU(2N_{F})$ above: The chiral spin symmetry and its flavor extension are
explicitly broken by magnetic interactions and by the quark kinetic terms.

Since the confining interaction in QCD is associated with the emergence of the
color-electric flux tubes that connect static quarks we can consider the
chiral spin symmetry as a symmetry of the confining interaction with the
ultrarelativistic light quarks. The necessary and sufficient conditions for
emergence of the approximate chiral spin symmetry are: (i) both chiral
symmetries in QCD must be at least approximately restored and (ii) the quark-
electric interaction must strongly dominate over the quark-magnetic
interaction and over quark kinetic terms. The latter condition means that the
physics must be dominated by a confining electric field.

The discussed symmetries were introduced in Refs. [1, 2]. For review on
symmetries and their implications for hadrons and for hot QCD see Ref. [3].

##  2 Emergence of the chiral spin symmetry above chiral restoration crossover
and its implications

The most detailed information about QCD is encoded in correlation functions

$C_{\Gamma}(t,x,y,z)=\langle
O_{\Gamma}(t,x,y,z)\,O_{\Gamma}(0,\mathbf{0})^{\dagger}\rangle\;.$ (9)

They carry the full spectral information.

In Fig. 1 we show spatial correlators of all possible $J=0,1$ isovector
bilinears $O_{\Gamma}$ at different temperatures obtained within $N_{F}=2$ QCD
with chirally symmetric Dirac operator [4].

Figure 1: Spatial correlators of all possible $J=0,1$ bilinears. Source: From
Ref. [4].

We see a distinct multiplet structure of the correlators. The multiplet $E1$
evidences restored $U(1)_{A}$ symmetry (at least approximately); the multiplet
$E2$ demonstrates approximate chiral spin and $SU(4)$ symmetries; the
multiplet $E3$ is consistent with chiral symmetry alone and with chiral spin
(and $SU(4)$) symmetry and hence could be ignored. The chiral spin symmetry
and its flavor extension persist up to $T\sim 500$ MeV. At larger temperatures
two distinct multiplets $E1$ and $E2$ disappear because the QCD correlators
approach correlators obtained with free quarks. Notice that the latter
correlators do not have the chiral spin symmetry because the Dirac Lagrangian
is not CS-symmetric; for analytic and lattice results for all possible ”free
correlators” see Ref. [4]. An apparent mergence of the $E1$ and $E2$
multiplets at very high temperatures is due to the fact that at
$T\rightarrow\infty$ all $J=0,1$ screening masses from $E1$ and $E2$ approach
their limiting value $2\pi T$, which is not related to any dynamical symmetry.

In Fig. 2 we show temporal correlators for all possible $J=0,1$ isovector
bilinears [5].

Figure 2: Temporal correlators for $12\times 48^{3}$ lattices. Left:
correlators obtained with free quarks. Right: full QCD correlators at $T=220$
MeV. Source: From Ref. [5]

The full QCD correlators demonstrate $U(1)_{A}$ and $SU(2)_{L}\times
SU(2)_{R}$ multiplets as well as approximate $SU(2)_{CS}$ and $SU(4)$
multiplets. In contrast, qualitatively different correlators obtained with
noninteracting quarks exhibit only $U(1)_{A}$ and $SU(2)_{L}\times SU(2)_{R}$
multiplets.

We conclude that above the chiral restoration crossover QCD partition function
is not only chiral symmetric, but is also approximately $SU(2)_{CS}$ and
$SU(4)$ symmetric. This suggests that QCD is still in the confining regime
until roughly $500$ MeV where confining chromoelectric field gets screened and
$SU(2)_{CS}$ and $SU(4)$ symmetries smoothly disappear. Similar results have
recently been obtained in $2+1+1$ QCD [6].

##  3 Breakdown of thermal perturbation theory below $500-600$ MeV

Lattice results for pseudoscalar and vector screening masses obtained at very
large temperatures $T\sim$ 1 - 160 GeV [7] are shown in Fig. 3, left panel.

Figure 3: Left:Temperature dependence of the pseudoscalar and vector screening
masses at very large T. Source: From Ref. [7]. Right: The same but for lower
temperatures. Source: From Ref. [8]

They can be parameterized [7] over two orders of magnitude by

$\displaystyle\frac{m_{PS}}{2\pi T}$ $\displaystyle=$ $\displaystyle
1+p_{2}\,\hat{g}^{2}(T)+p_{3}\,\hat{g}^{3}(T)+p_{4}\,\hat{g}^{4}(T)\;,$
$\displaystyle\frac{m_{V}}{2\pi T}$ $\displaystyle=$
$\displaystyle\frac{m_{PS}}{2\pi T}+s_{4}\,\hat{g}^{4}(T)\;,$

where $p_{2}$ is fixed by EQCD calculations and $p_{3},p_{4},s_{4}$ are
numbers that are fitted to lattice results. A perturbative description of
screening masses suggests partonic degrees of freedom, which is a signal of
QGP.

Screening masses at lower temperatures, above chiral restoration crossover up
to $T\sim 2.5$ GeV are shown in the right panel [8]. While the temperature
dependence of screening masses above 500 - 600 MeV is flat and consistent with
those shown in the left panel, at lower temperatures the screening masses
demonstrate a steep increase. The temperature dependence in the partonic (QGP)
regime is only in coupling constants. Given that the coefficients
$p_{2},p_{3},p_{4},s_{4}$ are fixed by the high T behavior of screening
masses, the perturbative expansion cannot explain screening masses below
500-600 MeV. We observe an apparent change of dynamics at temperatures below
500-600 MeV. The behavior of meson screening masses from 12 different channels
provides an independent demonstration of the existence of the regime below
500-600 MeV where chiral symmetry is restored but the dynamics is inconsistent
with the partonic description. The discussed behavior of screening masses must
be reflected in the equation of state. Indeed, a very steep increase of
$p/T^{4}$ with temperature in the same temperature interval is observed [9].
The analysis presented in this section was done in Ref. [10].

##  4 Pion spectral function

Direct evidence of the hadron-like degrees of freedom in stringy fluid should
be observation of states in spectral functions. Using a finite T
generalization of the Källen-Lehmann spectral representation [11] Refs. [12,
13] reconstructed pion and kaon spectral functions from the spatial lattice
correlators above chiral restoration. The pion spectral function is shown in
Fig. 4 [12]. The spectral function demonstrates two distinct peaks
corresponding to pion and its first radial excitations. These peaks get
broader with temperature and eventually melt above 500-600 MeV out. The so
obtained spectral function can be controlled because it predicts the temporal
correlators that can be compared with lattice results. The comparison shows a
good agreement at large Euclidean times.

Figure 4: Pion spectral function at different temperatures reconstructed from
spatial lattice correlators presented in Fig. 1. Source: From Ref. [12].

##  5 Conserved charges and their fluctuations in hadron gas and at higher
temperatures. Temperature evolution of the Polyakov loop.

The hadron resonance gas model assumes a dilute system of point-like
structureless hadrons that do not interact. For $T<155$ MeV the hadron
resonance gas model reproduces fluctuations of conserved charges calculated on
the lattice, see e.g. [14]. At larger temperatures the lattice results
radically deviate from the HRG model predictions.

Here [15] we focus on charges associated with the net number of up, down and
strange quarks:

$N_{q}\equiv\int d^{3}x~{}n_{q}(x)\;\;\;{\rm
with}\;\;\;n_{q}(x)=\bar{q}(x)\gamma^{0}q(x),\;\;\;q=u,d,s$ (10)

Each quark can be in one of the $N_{c}$ color states. This means that the
conserved flavor charges $N_{q}$, scale as $N_{c}^{1}$. The key point is that
the fluctuations of quark number densities scale as $N_{c}^{1}$ above chiral
crossover and that in the stringy fluid phase at large $N_{c}$, they would be
expected to differ from their vacuum and hadron gas values (of order
$N_{c}^{0}$) by $N_{c}$.

In Fig. 5 we show typical results for fluctuations of quark numbers of $u,d$
quarks taken from Ref. [14] and their comparison with the hadron resonance gas
model. We see that the fluctuations of the quark numbers deviate from the HRG
just at the chiral restoration temperature 155 MeV.

Figure 5: Left:Fluctuations of conserved net $u,d$ quark numbers in 2+1 QCD at
physical quark masses. Source: From Ref. [14]. Right: Temperature evolution of
the properly renormalized Polyakov loop in 2+1 QCD at physical quark masses.
Source: From Ref. [19]

Since the quark numbers scale as $N_{c}^{1}$ the above behavior of
fluctuations of conserved charges is consistent with the crossover from the
hadron gas to the stringy fluid regime, where all main thermodynamical
quantities scale as $N_{c}^{1}$ (it will be discussed in the next section). To
demonstrate a transition to the QGP regime one needs an observable that is
sensitive to presence of $\sim N_{c}^{2}$ deconfined gluons. A natural
observable of this kind is Polyakov loop.

The expectation value of the Polyakov loop in a pure glue theory is the order
parameter for ${\cal Z}_{N_{c}}$ center symmetry and for confinement. In the
center symmetric confining phase the expectation value of the Polyakov loop
identically vanishes while above the first order deconfinement phase
transition at $T_{d}\sim 270-300$ MeV the center symmetry gets spontaneously
broken and the expectation value of the Polyakov loop jumps to 0.5-0.6. An
important issue is that the Polyakov loop in the deconfined phase, where it is
not zero, is explicitly sensitive to $N_{c}^{2}-1$ gluons. Confining
properties of QCD with light quarks and of pure Yang-Mills theory are
identical in the large $N_{c}$ limit. Then one expects a similar deconfinement
temperature in QCD with light quarks at large $N_{c}$.

In the real world $N_{c}=3$ in QCD with light quarks the center symmetry of
the action is explicitly broken by quark loops and the deconfinement first
order phase transition is replaced by a very smooth crossover. The
renormalized Polyakov loop informs us about the deconfinement crossover
region. The temperature evolution of the properly renormalized Polyakov loop,
taken from Ref. [19], is shown in Fig. 5. We observe that indeed the
deconfinement crossover is nothing else but a smeared broad transition around
$T_{d}\sim 300$ MeV. Above the chiral restoration temperature around
$T_{ch}\sim 155$ MeV the Polyakov loop is very small which suggests that here
QCD is in the confining regime. At the same time at these temperatures
fluctuations of conserved charges demonstrate that the hadron gas picture does
not work. Both these facts are consistent with the crossover from the hadron
gas to the stringy fluid. The Polyakov loop reaches the value around 0.5 at a
temperature roughly $400-500$ MeV, in agreement with the temperature of smooth
disappearance of the chiral spin symmetry.

##  6 Large $N_{c}$ QCD phase diagram at $\mu_{B}=0$

In previous sections we have demonstrated numerous lattice evidences that upon
increasing temperature there are three qualitatively different regimes in QCD
that are connected by smooth crossovers: a hadron gas below $T\sim 155$ MeV,
an intermediate stringy fluid regime at $155<T<\sim 500$ MeV and a quark-gluon
plasma at higher temperatures. These regimes reflect different approximate
symmetries and effective degrees of freedom. In the hadron gas the degrees of
freedom are well separated hadrons and chiral symmetry is spontaneously
broken. In the stringy fluid chiral symmetry is restored and approximate
chiral spin symmetry emerges, which is a symmetry of confining interaction.
One could view this medium as a densely packed system of the color-singlet
clusters. The quark interchanges between clusters required by Pauli principle
should be significant. Still it is a system with confinement. The QGP matter
at high temperatures consists of decofined quark and gluon quasiparticles.

The large $N_{c}$ limit of QCD with massless quarks should clarify issues
since in this case QCD with light quarks is manifestly center symmetric. This
allows one to define unambiguously possible phases with confinement or
deconfinement and with spontaneously broken or restored chiral symmetry.
Standard large $N_{c}$ analysis [16] suggests that three regimes at $N_{c}=3$
and small but nonzero quark masses might become distinct phases separated by
phase transitions with different scaling of main thermodynamical quantities:
energy density $\epsilon$, pressure $P$ and entropy density $s$:

$\epsilon_{\rm HG}\sim N_{c}^{0}\;\;,\;\;P_{\rm HG}\sim N_{c}^{0}\;,\;s_{\rm
HG}\sim N_{c}^{0}\;,$ (11)

$\epsilon_{\rm str}\sim N_{c}^{1}\;\;,\;\;P_{\rm str}\sim N_{c}^{1}\;,\;s_{\rm
str}\sim N_{c}^{1}\;,$ (12)

$\epsilon_{\rm QGP}\sim N_{c}^{2}\;\;,\;\;P_{\rm QGP}\sim N_{c}^{2}\;,\;s_{\rm
QGP}\sim N_{c}^{2}\;.$ (13)

A peculiar feature of the intermediate phase is that it should contain a gas
of noninteracting glueballs, as in the hadron gas, while the degrees of
freedom containing quarks are radically changed. Lattice large $N_{c}$ studies
could clarify this picture.

It is known that the deconfinement temperature in pure Yang-Mills theory is
practically $N_{c}$-independent [17]. At the same time confining properties of
QCD with light quarks are identical with those in Yang-Mills at infinite
$N_{c}$. Then one expects a deconfinement temperature in large $N_{c}$ QCD
with light quarks to be around $T_{d}\sim 300$ MeV. It is also found on the
lattice that at $T=0$ the quark condensate scales practically exactly as
$N_{c}^{1}$ and the $1/N_{c}$ corrections are negligibly small [18]. This
strongly suggests that the physics of chiral symmetry breaking in QCD at
$N_{c}=3$ and at large $N_{c}$ is the same. Then one expects a chiral
restoration phase transition in large $N_{c}$ QCD around the same temperature
as at $N_{c}=3$. The latter temperature was established to be around
$T_{ch}\sim 130$ MeV [20]. If correct, one anticipates two different phase
transitions at large $N_{c}$, as was discussed above.

##  7 Outlook

The next important step to clarify the QCD phase diagram would be to perform
lattice measurements of the chiral restoration temperature at large $N_{c}$.
It is important to stress that the chiral and large $N_{c}$ limits do not
commute and a proper sequence of limits should be taken. First the large
$N_{c}$ limit and then the chiral limit. This sequence of limits significantly
simplifies the numerical problem. As the first step the pure glue (quenched)
configurations should be generated at large $N_{c}$. There are two strategies
of doing this. The first strategy would be to employ gluonic ensembles at
large volumes at reasonably large $N_{c}\sim 10-20$ as was done in Ref. [17].
Then the eigenvalue problem for the Dirac operator should be solved to
determine the quark condensate via the Banks-Casher relation. Upon increasing
temperature the temperature of the chiral restoration phase transition could
be obtained. The alternative strategy would be to follow the Eguchi-Kawai
volume independence to address the problem at very large $N_{c}$ but in a not
large volume, as was done in Ref. [18]. The results for chiral restoration
temperatures should give in principle the same answer in both approaches.

## References

  * [1] L. Y. Glozman, Eur. Phys. J. A 51 (2015) no.3, 27 doi:10.1140/epja/i2015-15027-x [arXiv:1407.2798 [hep-ph]].
  * [2] L. Y. Glozman and M. Pak, Phys. Rev. D 92 (2015) no.1, 016001 doi:10.1103/PhysRevD.92.016001 [arXiv:1504.02323 [hep-lat]].
  * [3] L. Y. Glozman, Prog. Part. Nucl. Phys. 131 (2023), 104049 doi:10.1016/j.ppnp.2023.104049 [arXiv:2209.10235 [hep-lat]].
  * [4] C. Rohrhofer, Y. Aoki, G. Cossu, H. Fukaya, C. Gattringer, L. Y. Glozman, S. Hashimoto, C. B. Lang and S. Prelovsek, Phys. Rev. D 100 (2019) no.1, 014502 doi:10.1103/PhysRevD.100.014502 [arXiv:1902.03191 [hep-lat]].
  * [5] C. Rohrhofer, Y. Aoki, L. Y. Glozman and S. Hashimoto, Phys. Lett. B 802 (2020), 135245 doi:10.1016/j.physletb.2020.135245 [arXiv:1909.00927 [hep-lat]].
  * [6] T. W. Chiu, Phys. Rev. D 107 (2023) no.11, 114501 doi:10.1103/PhysRevD.107.114501 [arXiv:2302.06073 [hep-lat]].
  * [7] M. Dalla Brida, L. Giusti, T. Harris, D. Laudicina and M. Pepe, JHEP 04 (2022), 034 doi:10.1007/JHEP04(2022)034 [arXiv:2112.05427 [hep-lat]].
  * [8] A. Bazavov, S. Dentinger, H. T. Ding, P. Hegde, O. Kaczmarek, F. Karsch, E. Laermann, A. Lahiri, S. Mukherjee and H. Ohno, et al. Phys. Rev. D 100 (2019) no.9, 094510 doi:10.1103/PhysRevD.100.094510 [arXiv:1908.09552 [hep-lat]].
  * [9] A. Bazavov, P. Petreczky and J. H. Weber, Phys. Rev. D 97 (2018) no.1, 014510 doi:10.1103/PhysRevD.97.014510 [arXiv:1710.05024 [hep-lat]].
  * [10] L. Y. Glozman, O. Philipsen and R. D. Pisarski, Eur. Phys. J. A 58 (2022) no.12, 247 doi:10.1140/epja/s10050-022-00895-4 [arXiv:2204.05083 [hep-ph]].
  * [11] J. Bros and D. Buchholz, Z. Phys. C 55 (1992), 509-514 doi:10.1007/BF01565114
  * [12] P. Lowdon and O. Philipsen, “Pion spectral properties above the chiral crossover of QCD,” JHEP 10 (2022), 161 doi:10.1007/JHEP10(2022)161 [arXiv:2207.14718 [hep-lat]].
  * [13] D. Bala, O. Kaczmarek, P. Lowdon, O. Philipsen and T. Ueding, [arXiv:2310.13476 [hep-lat]].
  * [14] R. Bellwied, S. Borsanyi, Z. Fodor, S. D. Katz, A. Pasztor, C. Ratti and K. K. Szabo, Phys. Rev. D 92 (2015) no.11, 114505 doi:10.1103/PhysRevD.92.114505 [arXiv:1507.04627 [hep-lat]].
  * [15] T. D. Cohen and L. Y. Glozman, [arXiv:2401.04194 [hep-ph]].
  * [16] T. D. Cohen and L. Y. Glozman, [arXiv:2311.07333 [hep-ph]].
  * [17] B. Lucini and M. Panero, Phys. Rept. 526 (2013), 93-163 doi:10.1016/j.physrep.2013.01.001 [arXiv:1210.4997 [hep-th]].
  * [18] C. Bonanno, P. Butti, M. García Peréz, A. González-Arroyo, K. I. Ishikawa and M. Okawa, [arXiv:2309.15540 [hep-lat]].
  * [19] P. Petreczky and H. P. Schadler, Phys. Rev. D 92 (2015) no.9, 094517 doi:10.1103/PhysRevD.92.094517 [arXiv:1509.07874 [hep-lat]].
  * [20] H. T. Ding et al. [HotQCD], Phys. Rev. Lett. 123 (2019) no.6, 062002 doi:10.1103/PhysRevLett.123.062002 [arXiv:1903.04801 [hep-lat]].


#  The extended binary quadratic residue code of length $42$ holds a
$3-$design

A. Bonnecaze, P. Solé<EMAIL_ADDRESS>corresponding author,
<EMAIL_ADDRESS>AB and PS are with Aix Marseille Univ, CNRS, Centrale Marseille,
I2M, Marseille, France.

###### Abstract

The codewords of weight $10$ of the $[42,21,10]$ extended binary quadratic
residue code are shown to hold a design of parameters $3-(42,10,18).$ Its
automorphism group is isomorphic to $PSL(2,41).$ Its existence can be
explained neither by a transitivity argument, nor by the Assmus-Mattson
theorem.

Keywords: designs, quadratic residue codes, Assmus-Mattson theorem

AMS Math Sc. Cl. (2010): 94 B15, 62K10

##  1 Introduction

The quadratic residue codes have been known to hold designs of large strength
since the seminal paper [1], which constructs in particular the Witt designs
by coding theoretic techniques. The argument in that famous paper (see also
[2] for background) is the so-called Assmus-Mattson theorem which derives the
existence of designs from certain facts on the weight distribution of the dual
code [13, Chap. 6, Th. 29].

In the present note we produce, by electronic calculation, a $3-(42,10,18)$
design which is not in the table [7, pp. 82–83], cannot be derived from the
Assmus-Mattson Theorem, and does not follow by the standard transitivity
argument. Its blocks are the supports of the codewords of weight $10$ in the
$[42,21,10]$ extended binary quadratic residue code. Its automorphism group
coincides with that of that code. Thus it is isomorphic to $PSL(2,41)$ which
is not $3$-homogeneous. Note that $PSL(2,q)$ is $3$-homogeneous in general for
$q\equiv 3\pmod{4}$ [8]. Thus our design cannot be obtained by the methods of
[6]. Neither can it be obtained from the orbits glueing methods of [3, Th. 5],
or [10, Th. 5.1], since $41\equiv 9\pmod{16}.$ Even though the code we are
using is a formally self-dual even code, the results of [9] do not apply,
because the minimum distance is not high enough.

##  2 Background material

###  2.1 Permutation groups

A permutation group $G$ acting on a set $X$ is transitive on $X$ if it has a
single orbit on $X.$ Given an integer $s$ it is $s$-fold transitive if it is
transitive on ordered $s$-uples of distinct elements of $X.$ It is $s$-fold
homogeneous if it is transitive on unordered $s$-uples of distinct elements of
$X.$ Let $p$ denote an odd prime. The group $PSL(2,p)$ is defined in its
action on the projective line $\mathbb{F}_{p}\cup\infty$ as the set of all
transforms

$y\mapsto\frac{a+by}{c+dy}$

with $ad-bc=1,$ and the conventions $\frac{0}{0}=0$ and $\frac{\neq
0}{0}=\infty.$ See [13, Chap. 16, §5], or [14, Chap 17] for background.

###  2.2 Designs

A block design of parameters $t-(b,v,k,r,\lambda)$ ( shortly a
$t-(v,k,\lambda)$ design) is an incidence structure
$(\mathcal{P},\mathcal{B},I)$ satisfying the following axioms.

  1. 1.

The set $\mathcal{P}$ has $v$ points.

  2. 2.

The set $\mathcal{B}$ has $b$ blocks.

  3. 3.

Each block is incident to $k$ points.

  4. 4.

Each point is incident to $r$ blocks.

  5. 5.

Each $t$-tuple of points is incident in common with $\lambda$ blocks.

###  2.3 Codes

Let $\mathbb{F}_{2}=\\{0,1\\}$ denote the finite field of order $2.$ A binary
code of length $n$ is a $\mathbb{F}_{2}$ subspace of $\mathbb{F}_{2}^{n}.$ The
weight of a vector of $\mathbb{F}_{2}^{n}$ is the number of its nonzero
coordinates. The weight distribution of a code $C$ is the sequence $A_{w}$ of
number of codewords of $C$ of weight $w.$ It is written in Magma [12] notation
as the list with generic element $\langle w,A_{w}\rangle$ where $w$ ranges
over the weights of $C.$ A binary code is cyclic if it is invariant under the
cyclic shift. Cyclic codes are in one to one correspondence with ideals of the
residue class ring $\mathbb{F}_{2}[x]/(x^{n}-1).$ The generator polynomial of
a cyclic code is then the generator of the corresponding ideal. The Quadratic
residue codes are the cyclic codes of length $p,$ with $p$ an odd prime,
defined for $p\equiv\pm 1\pmod{8}$ by the generator polynomial of degree
$\frac{p-1}{2}$

$\prod_{r=\Box}(x-\alpha^{r}),$

with $\alpha$ a primitive root of order $p.$ Since $2$ is a quadratic residue
modulo an odd prime $p$ this polynomial is indeed in $\mathbb{F}_{2}[x].$ As a
simple example consider the case of $p=7$ when the polynomial

$(x-\alpha)(x-\alpha^{2})(x-\alpha^{4})=x^{3}+x+1$

generates the Hamming code $[7,4,3].$ See [13, Chap.16] for background.

The automorphism group a binary code of length $n$ is the subgroup of the
symmetric group on the $n$ coordinate places that leaves the code wholly
invariant.

The Assmus-Mattson Theorem is as follows for binary codes [13, Chap. 6, Th.
29] .

###### Theorem 1

(Assmus-Mattson) If $C$ is a binary $[n,k,d]$ code such that the weight
distribution of its dual code contains at most $d-t$ nonzero weights $\leq
n-t$ then the codewords of weight $d$ of $C$ form a $t$-design.

A folklore theorem is that if the automorphism group of a binary code is
$t$-homogeneous then the codewords of any given weight $\geq t$ hold
$t$-designs [14, p. 308]. .

##  3 Construction

Let $Q$ be the binary quadratic residue code of prime length $41,$ and denote
by $C$ its extension by an overall parity-check. The code $C$ is not self-
dual, but it is equivalent to its dual [13, Chap. 16, Fig. 16.3]. Designs in
even formally self-dual codes have been explored in [9], but the codes there
are extremal, which would require a $[42,21,12].$

The weight distribution of $C$ is also its dual weight distribution

$[\langle 0,1\rangle,\langle 10,1722\rangle,\langle 12,10619\rangle,\langle
14,49815\rangle,\langle 16,157563\rangle,\langle 18,341530\rangle,$ $\langle
20,487326\rangle,\langle 22,487326\rangle,\langle 24,341530\rangle,\langle
26,157563\rangle,\langle 28,49815\rangle,\langle 30,10619\rangle,$ $\langle
32,1722\rangle,\langle 42,1\rangle].$

Thus the Assmus-Mattson theorem, which requires at most $10-3=7$ non zero dual
weights $\leq 42-3=39,$ cannot apply here.

###### Theorem 2

The $1722$ codewords of weight $10$ of $C$ hold a $3-(1722,42,10,410,18)$
block design $\mathcal{D}.$

Proof. Magma computation [12].

Its automorphism group is a permutation group of order $34440$ with generators

$(3,30,29,31)(4,9,18,7)(5,24,25,17)(6,22,38,42)(8,34,11,28)$
$(10,36,16,33)(12,32,23,21)(13,15,14,41)(19,20,39,27)(26,35,40,37),$
$(3,8,6,33,15,29,11,38,36,41)(4,35,27,21,5,18,37,20,32,25)$
$(7,26,39,23,17,9,40,19,12,24)(10,14,31,28,42,16,13,30,34,22),$
$(1,32)(2,21)(3,36)(5,10)(6,26)(7,38)(8,20)(9,22)(11,25)(12,23)(13,33)(14,28)$
$(15,39)(16,27)(17,30)(18,37)(19,31)(24,41)(29,34)(40,42),$
$(2,32,23,12,36,29,37,14,24,10,8,15,40,6,4,31,41,18,26,19)$
$(3,42,11,34,9,28,20,17,30,33,39,7,16,22,25,38,35,27,5,21),$

This group can be characterized as an abstract group as follows.

###### Proposition 1

The automorphism group $G$ of $\mathcal{D}$ is isomorphic to $PSL(2,41).$

Proof. The automorphism group of $C$ is known to be isomorphic to $PSL(2,41)$
by [14, Chap. 17]. Since the blocks of $\mathcal{D}$ consist of all the
codewords of weight $10$ of $C$, they are wholly invariant under the
automorphism group of $C.$ A Magma computation [12] shows that the order of
the automorphism group $\mathcal{D}$ is $34440=\frac{1}{2}41\times(41^{2}-1),$
the known order of $PSL(2,41)$ [13, Chap. 16, Th. 9 (b)]. The result follows.

The group $G$ is not $3$-homogeneous.

###### Proposition 2

The group $G$ partitions the set of triples into two orbits of equal size.

Proof. A Magma computation shows that the orbits of $\\{1,2,3\\},$ and
$\\{1,3,8\\},$ are disjoint and both of size $5740={42\choose 3}/2.$

Remarks:

  1. 1.

The codewords of weights $9$ and $10$ in the quadratic residue code of length
$41$ hold $2$-designs of parameters $2-(41,9,18)$ and $2-(41,10,72)$ They
cannot be explained by the Assmus-Mattson theorem or group action. The
$2-(41,9,18)$ is a derived design of $\mathcal{D}.$

  2. 2.

It can be shown that the linear span of $\mathcal{D}$ is $C.$ This gives an
alternative proof of Proposition 1.

  3. 3.

The codewords of minimum weight in the extended quadratic residue code of
length $74$ do not hold 3-designs.

  4. 4.

The codewords of weight 10 and 12 in the duadic codes of length $74$ of [11]
do not hold 3-designs.

##  4 Conclusion

In this work, we have constructed a new $3$-design in the minimum weight
codewords of a quadratic residue codes. Its existence cannot be explained by
the Assmus-Mattson Theorem, nor by group action.

The arguments are based on machine calculation. It would be more satisfying to
have a conceptual derivation. A method by multivariate weight enumerators and
invariant theory in the spirit of [5], would be more conceptual, but certainly
not computer-free. The fact that the low weight codewords in the extended
quadratic residue code of length $74$ or in the duadic codes of that length do
not hold 3-designs suggests that the design in this note is exceptional enough
that a general conjecture is impossible to formulate.

Acknowledgement: The authors thank Professor Jenny Key for helpful
discussions.

## References

  * [1] E.F. Assmus jr, H.F. Mattson jr, New 5-designs, J. of Combinatorial Theory A,6, (1969), 122–151.
  * [2] E.F. Assmus jr, H.F. Mattson jr, R. Turyn, Research to develop the algebraic theory of codes, Sylvania report AFCRL-67-0365, (1967), available from https://apps.dtic.mil/dtic/tr/fulltext/u2/656783.pdf.
  * [3] N. Balachandran, D. Ray-Chaudhuri, Simple $3$-designs and $PSL(2,q)$ with $q\equiv 1\pmod{4}$. Des. Codes Cryptogr. 44 (2007), 263–274.
  * [4] T. Beth, D. Jungnickel, H. Lenz, Design Theory, 2nd edition, vol. I, (1999), Cambridge University Press.
  * [5] A. Bonnecaze, B. Mourrain, P. Solé, Jacobi Polynomials, Type II Codes, and Designs. Des. Codes Cryptogr. 16, (1999), 215–234.
  * [6] P.J. Cameron, H.R. Maimanid, G.R. Omidic, B. Tayfeh-Rezaie, 3-Designs from $PSL(2,q),$ Discrete Mathematics 306 (2006) 3063–3073.
  * [7] C.J. Colbourn, J.H. Dinitz (2007), Handbook of Combinatorial Designs (2nd ed.), Boca Raton: Chapman & Hall/ CRC, ISBN 1-58488-506-8.
  * [8] C. Ding, Designs from linear codes, World Scientific (2018) Singapore.
  * [9] G.T. Kennedy and V. Pless, On designs and formally self-dual codes, Designs Codes Cryptography 4, (1994) 43–55.
  * [10] W. Liu, J-X. Tang, Y-X. Wu, Some new $3$-designs from $PSL(2,q)$ with $q\equiv 1\pmod{4},$ Science in China, 55,(2012), 1901–1911.
  * [11] J.S. Leon, J. M. Masley, V. Pless, On duadic codes, IEEE Trans. on Information Theory, (1984), 709–714.
  * [12] http://magma.maths.usyd.edu.au/calc/
  * [13] F.J. MacWilliams, N.J. A. Sloane, The theory of Error Correcting Codes, North-Holland, Amsterdam (1981).
  * [14] V.S. Pless, W.C. Huffman, Handbook of Coding theory, North Holland (1998), Amsterdam.


# Sensitive detection of local magnetic field changes with atomic
interferometry by using superconducting Meissner effects

Y. Q. Chai  School of Physical Science and Technology, Southwest Jiaotong
University, Chengdu 610031, China  Information Quantum Technology Laboratory,
International Cooperation Research Center of China Communication and Sensor
Networks for Modern Transportation, School of Information Science and
Technology, Southwest Jiaotong University, Chengdu 610031, China  M.
<EMAIL_ADDRESS>School of Physical Science and Technology,
Southwest Jiaotong University, Chengdu 610031, China  L. F. Wei  Information
Quantum Technology Laboratory, International Cooperation Research Center of
China Communication and Sensor Networks for Modern Transportation, School of
Information Science and Technology, Southwest Jiaotong University, Chengdu
610031, China

###### Abstract

Sensitive detection of magnetic field is one of the open problem in metrology.
Here, we propose an Mach-Zehnder atomic interferometry to sensitively detect
the very weak local magnetic field, which is expelled by the superconductor
(as the “testing magnet”) due to the Meissner effect. The induced magnetic
field gradient near the superconductor provides a centripetal acceleration of
the atomic motion in the interferometry and thus can be detected by using the
atomic interferences. Given gravity acceleration of the atoms have been
measured at the accuracy of $10^{-12}$ g, the measured sensitivity of the
expelled local field could reach $10^{-14}$ T.

##  I introduction

High precision measurements of physical quantities take the important roles in
modern science and technology. As a powerful and precise measuring tool,
atomic interferometry has been widely used for various metrologies 1 ; Zhu ; 3
, e.g., gyroscope 4 , gravimeter 5 ; 6 ; 7 , the measurement of fine structure
constant 8 , the test of equivalence principle EEP , the demonstration of
electrodynamics Aharonov-Casher (AC) effect AC , the gravitational Aharonov-
Bohm (AB) effect AB , the gravitational wave detection gravitational wave ,
and the test of chameleon theories of dark energy Dark Energy . Particularly,
the accuracy of measuring gravitational acceleration can reach $\delta
g=10^{-12}g$ EEP  by using the Mach-Zehnder (MZ) atom interferometer. In this
atom interferometer, three laser pulses are applied to split and recombine the
atomic paths. The gravitational potential between the two separated paths of
atom results in a signal $kgt^{2}$ of atomic matter-wave interference, with
$k$ being the effective wave number of laser and $t$ the freely falling time
of atom in vacuum. The high precision measurement on gravity is attributed to
the large value of optical wave number $k$, the ensemble of many cold atoms,
and the atom-detection at single atom level.

It is well known that, in addition to the gravity field, there is also a
magnetic field called usually as the geomagnetic field around the earth. The
precise measurement of geomagnetic field has important significance and
applications, such as for the geophysical exploration, geomagnetism,
archaeology and so on 10 ; 11 ; 12 ; 13 . Usually, the tool for measuring
magnetic field is called as the magnetometer 14 ; 15 ; 16 ; 17 ; 18 ; 19 .
Given the superconducting quantum interferometer (SQUID) has been developed
well for the various sensitive magnetometer applications 15 , and also the
atomic interferometer possesses the high-precision measurement ability for the
gravitational acceleration, in this paper we discuss how to combine these two
metrology techniques for the implementation of the sensitivity detection of
the earth’s magnetic field. Basically, the earth’s magnetic field is too weak
to be detected directly by the atomic interferometer, which is originally
designed to implement the sensitive measurement of the gravitational
acceleration. However, if a superconductor Meissner0 ; Meissner1  served as
the “testing mass” is introduced, then the relevant earth’s magnetic field is
expelled outside (due to the Meissner effect), yielding a gradient of the
geomagnetic field be generated. Furthermore, an atomic interferometer is
designed to measure such a geomagnetic field gradient. As a consequence, the
information of geomagnetic field in the local location of the conductor (which
becomes the superconductor under the low temperature) can be extracted.
Physically, the measurement accuracy of the local magnetic field is related to
the one of the local gravitational acceleration.

This article is organized as follows. In Sec. II, by specifically solving the
Maxwell (London) equation outside (inside) the superconductor, we show how a
geomagnetic field gradient is generated. In Sec. III, an atomic interferometer
is designed to implement the sensitive detection of such a magnetic gradient.
Its achievable accuracy is also discussed. Finally, we summarize our work in
Sec. IV.

##  II The geomagnetic field gradient generated by using superconducting
Meissner effect

First, we discuss how to the local geomagnetic field gradient from the very
weak geomagnetic field background. Let us introduce a cylindrical conductor to
serve as the “test mass”, which can be cooled to become a superconductor.
Before the cooling, the local geomagnetic field around it is assumed to be
uniformly distributed and thus served as the very weak geomagnetic field
background. By cooling the conductor to become superconductivity, the local
background will be expelled due to the Meissner effect. The configuration is
shown in Fig. 1.

Figure 1: The coordinate system is established with the axis of the
superconducting cylinder as the $z$ axis. Considering the length $L$ of
cylindrical superconductor is much larger than its radius $R_{0}$, then the
magnetic field gradient along $z$-direction is negligible within the uniform
geomagnetic background.

###  II.1 The magnetic field outside the superconductor

Basically, the magnetic field B outside the superconductor and its vector
potential A satisfies the relation:

$\displaystyle\nabla\times\textbf{B}=\nabla\times(\nabla\times\textbf{A})=-\nabla^{2}\textbf{A}=0\,,$
(1)

where $\nabla$ is the laplace operator and the Coulomb gauge
$\nabla\cdot\textbf{A}=0$ is used. In the cylindrical coordinate system
$(r,\theta,z)$ with $r=\sqrt{x^{2}+y^{2}}$ and $\theta=\arctan(y/x)$, the
magnetic field reads

B
$\displaystyle=B_{r}(r,\theta)\textbf{e}_{r}+B_{\theta}(r,\theta)\textbf{e}_{\theta}+B_{z}\textbf{e}_{z}\,,$
(2)

which can be further expressed as

$\displaystyle B_{r}(r,\theta)=\frac{1}{r}\frac{\partial
A_{z}(r,\theta)}{\partial\theta}\,,$ (3) $\displaystyle
B_{\theta}(r,\theta)=-\frac{\partial A_{z}(r,\theta)}{\partial r}\,.$ (4)

Note that, the $B_{z}$ is the $z$-directional component of the magnetic field,
which is invariant due to the symmetry of the long cylindrical superconductor.
In Eqs. (3) and (4), the $A_{z}$ is the $z$-directional component of the
magnetic vector potential.

In the cylindrical coordinate system, the laplace equation (1) for $A_{z}$
outside the superconductor reads

$\displaystyle\frac{1}{r^{2}}\frac{\partial^{2}A_{z}}{\partial\theta^{2}}+\frac{\partial^{2}A_{z}}{\partial
r^{2}}+\frac{1}{r}\frac{\partial A_{z}}{\partial r}=0\,.$ (5)

Writing $A_{z}=f(r)\Phi(\theta)$, the laplace equation becomes

$\displaystyle r^{2}\frac{{\rm d}^{2}f}{{\rm d}r^{2}}+r\frac{{\rm d}f}{{\rm
d}r}-\beta^{2}f=0\,,$ (6)

and

$\displaystyle\frac{{\rm d}^{2}\Phi}{{\rm d}\theta^{2}}+\beta^{2}\Phi=0\,,$
(7)

with $\Phi(\theta)=\Phi(\theta+2\pi)$ and $\beta=0,1,2,\cdots$. The solutions
of Eqs. (6) and (7) are well-known

$\displaystyle
f(r)=c_{0}+d_{0}\ln(r)+\sum_{n=1}^{\infty}\left(c_{n}r^{n}+d_{n}r^{-n}\right)\,,$
(8)
$\displaystyle\Phi(\theta)=\sum_{n=0}^{\infty}a_{n}\cos(n\theta)+b_{n}\sin(n\theta)\,,$
(9)

with $a_{n}$, $b_{n}$, $c_{n}$, $d_{n}$ being the arbitrary real numbers.

For $r\rightarrow\infty$, the magnetic field reduces to that of Earth, i.e.,

$\displaystyle\textbf{B}_{0}$
$\displaystyle=B_{0x}\textbf{e}_{x}+B_{0y}\textbf{e}_{y}$ (10)
$\displaystyle=\left[B_{0x}\cos(\theta)+B_{0y}\sin(\theta)\right]\textbf{e}_{r}+\left[B_{0y}\cos(\theta)-B_{0x}\sin(\theta)\right]\textbf{e}_{\theta}\,.$

This means that,

$\left\\{\begin{aligned} &\frac{1}{r}\frac{\partial
A_{z}(r,\theta)}{\partial\theta}=B_{0x}\cos(\theta)+B_{0y}\sin(\theta)\,,\,\,\,\,{\rm
with}\,\,r\rightarrow\infty\,,\\\ \\\ &-\frac{\partial
A_{z}(r,\theta)}{\partial
r}=B_{0y}\cos(\theta)-B_{0x}\sin(\theta)\,,\,\,\,\,{\rm
with}\,\,r\rightarrow\infty\,.\end{aligned}\right.$ (11)

According to Eqs. (8), (9), and (11), we have

$\displaystyle
A_{z}=C_{0}+r^{-1}\left[C_{1}\cos(\theta)+C_{1}\sin(\theta)\right]+r\left[B_{0x}\sin(\theta)-B_{0y}\cos(\theta)\right]\,,$
(12)

and here $C_{0}$ and $C_{1}$ are the regularized coefficient.

###  II.2 The magnetic field inside the superconductor

The magnetic field $\tilde{\textbf{B}}$ inside the superconductor can be
described by the London equation London0 ; London1 ; London2 :

$\displaystyle\nabla\times\textbf{j}_{s}=-\frac{1}{\mu_{0}\lambda^{2}}\tilde{\textbf{B}}\,,$
(13)

with $\textbf{j}_{s}$, $\lambda$, and $\mu_{0}$ being the supercurrent
density, London penetration depth, and the vacuum permeability, respectively.
The above London equation can be further written as

$\displaystyle\nabla^{2}\tilde{\textbf{A}}-\frac{1}{\lambda^{2}}\tilde{\textbf{A}}=0\,,$
(14)

with $\textbf{j}_{s}=-\tilde{\textbf{A}}/\lambda^{2}$. Inside the
superconductor, the magnetic potential
$\tilde{A}_{z}=\tilde{f}(r)\tilde{\Phi}(\theta)$ is described by the London
equation (14), i.e.,

$\displaystyle\frac{1}{r^{2}}\frac{\partial^{2}\tilde{A}_{z}}{\partial\theta^{2}}+\frac{\partial^{2}\tilde{A}_{z}}{\partial
r^{2}}+\frac{1}{r}\frac{\partial\tilde{A}_{z}}{\partial
r}=\frac{1}{\lambda^{2}}\tilde{A}_{z}\,.$ (15)

Obviously, the solution of $\tilde{\Phi}$ is similar to Eq. (9), i.e.,

$\displaystyle\tilde{\Phi}(\theta)=\sum_{n=0}^{\infty}\tilde{a}_{n}\cos(n\theta)+\tilde{b}_{n}\sin(n\theta)\,,$
(16)

with the arbitrary constant coefficients $\tilde{a}_{n}$ and $\tilde{b}_{n}$.
While, the function $\tilde{f}(r)$ obeys the following equation,

$\displaystyle r^{2}\frac{{\rm d}^{2}\tilde{f}}{{\rm d}r^{2}}+r\frac{{\rm
d}\tilde{f}}{{\rm
d}r}-\left(\frac{r^{2}}{\lambda^{2}}+n^{2}\right)\tilde{f}=0\,,$ (17)

which is solved as

$\displaystyle\tilde{f}(r)=\sum_{n=0}^{\infty}I_{n}\left(\frac{r}{\lambda}\right)\,.$
(18)

Here,

$\displaystyle
I_{n}(x)=\left(\frac{x}{2}\right)^{n}\sum_{k=0}^{\infty}\frac{1}{\Gamma(n+k+1)}\left(\frac{x}{2}\right)^{2k}\,,$
(19)

is the deformed Bessel function of order $n$. As consequence, the solution of
$\tilde{A}_{z}$ reads

$\displaystyle\tilde{A}_{z}=\sum_{n,m}^{\infty}I_{m}\left(\frac{r}{\lambda}\right)\left[\tilde{a}_{n}\cos(n\theta)+\tilde{b}_{n}\sin(n\theta)\right]\,.$
(20)

According the boundary condition $A_{z}(R_{0})=\tilde{A}_{z}(R_{0})$ at the
surface of the cylindrical superconductor, we have

$\displaystyle\tilde{A}_{z}=I_{0}(\frac{r}{\lambda})\tilde{a}_{0}+I_{1}(\frac{r}{\lambda})\left[\frac{B_{0x}(R_{0}+C_{1}R_{0}^{-1})}{I_{1}(\frac{R_{0}}{\lambda})}\sin(\theta)-\frac{B_{0y}(R_{0}+C_{1}R_{0}^{-1})}{I_{1}(\frac{R_{0}}{\lambda})}\cos(\theta)\right]\,.$
(21)

Furthermore, using the continuity condition
$\textbf{e}_{\theta}\cdot\textbf{B}(R_{0})=\textbf{e}_{\theta}\cdot\tilde{\textbf{B}}(R_{0})$
of magnetic fields, i.e.,

$\displaystyle\frac{\partial A_{z}}{\partial
r}\bigg{|}_{r=R_{0}}=\frac{\partial\tilde{A}_{z}}{\partial
r}\bigg{|}_{r=R_{0}}\,,$ (22)

we have

$\left\\{\begin{aligned} &\tilde{a}_{0}=0\,,\\\ \\\
&C_{1}=2R_{0}\lambda\frac{I_{1}(\frac{R_{0}}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}-R_{0}^{2}\,.\end{aligned}\right.$
(23)

According to Eqs. (12), (21), and (23), the magnetic field outside and inside
the superconductor can be obtained as,

$\displaystyle\textbf{B}=$
$\displaystyle\left[1-\frac{R_{0}^{2}}{r^{2}}+\frac{2R_{0}\lambda}{r^{2}}\frac{I_{1}(\frac{R_{0}}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\right][B_{0y}\sin(\theta)+B_{0x}\cos(\theta)]\textbf{e}_{r}$
(24)
$\displaystyle+\left[1+\frac{R_{0}^{2}}{r^{2}}-\frac{2R_{0}\lambda}{r^{2}}\frac{I_{1}(\frac{R_{0}}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\right][B_{0y}\cos(\theta)-B_{0x}\sin(\theta)]\textbf{e}_{\theta}$
$\displaystyle=$ $\displaystyle
B_{x}(x,y)\textbf{e}_{x}+B_{y}(x,y)\textbf{e}_{y}\,,$

and

$\displaystyle\tilde{\textbf{B}}=$ $\displaystyle\frac{2\lambda
I_{1}(\frac{r}{\lambda})}{rI_{0}(\frac{R_{0}}{\lambda})}\left[B_{0y}\sin(\theta)+B_{0x}\cos(\theta)\right]\textbf{e}_{r}$
(25) $\displaystyle+$
$\displaystyle\frac{2}{I_{0}(\frac{R_{0}}{\lambda})}\left[I_{0}(\frac{r}{\lambda})-\frac{\lambda}{r}I_{1}(\frac{r}{\lambda})\right]\left[B_{0y}\cos(\theta)-B_{0x}\sin(\theta)\right]\textbf{e}_{\theta}$
$\displaystyle=$
$\displaystyle\tilde{B}_{x}(x,y)\textbf{e}_{x}+\tilde{B}_{y}(x,y)\textbf{e}_{y}\,.$

In the rectangular coordinate system, the magnetic field inside and outside
the superconductor are expressed as

$\left\\{\begin{aligned}
&B_{x}=B_{0x}+\frac{R_{0}(y^{2}-x^{2})}{r^{4}}\left[R_{0}-2\lambda\frac{I_{1}(\frac{R_{0}}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\right]B_{0x}-\frac{R_{0}xy}{r^{4}}\left[R_{0}-2\lambda\frac{I_{1}(\frac{R_{0}}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\right]B_{0y}\,,\\\
\\\
&B_{y}=B_{0y}+\frac{R_{0}(x^{2}-y^{2})}{r^{4}}\left[R_{0}-2\lambda\frac{I_{1}(\frac{R_{0}}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\right]B_{0y}-\frac{R_{0}xy}{r^{4}}\left[R-2\lambda\frac{I_{1}(\frac{R_{0}}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\right]B_{0x}\,,\\\
\\\ &\tilde{B}_{x}=\frac{2\lambda
x^{2}}{r^{3}}\frac{I_{1}(\frac{r}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\left(\frac{y}{x}B_{0y}+B_{0x}\right)+\frac{2x^{2}}{r^{2}I_{0}(\frac{R_{0}}{\lambda})}\left[I_{0}(\frac{r}{\lambda})-\frac{\lambda}{r}I_{1}(\frac{r}{\lambda})\right]\left[\frac{y^{2}}{x^{2}}B_{0x}-\frac{y}{x}B_{0y}\right]\,,\\\
\\\ &\tilde{B}_{y}=\frac{2\lambda
x^{2}}{r^{3}}\frac{I_{1}(\frac{r}{\lambda})}{I_{0}(\frac{R_{0}}{\lambda})}\left(\frac{y^{2}}{x^{2}}B_{0y}+\frac{y}{x}B_{0x}\right)+\frac{2x^{2}}{r^{2}I_{0}(\frac{R_{0}}{\lambda})}\left[I_{0}(\frac{r}{\lambda})-\frac{\lambda}{r}I_{1}(\frac{r}{\lambda})\right]\left[B_{0y}-\frac{y}{x}B_{0x}\right]\,.\end{aligned}\right.$
(26)

Following Eqs. (26), the Meissner effect of superconductor within the Earth’s
magnetic field is resolvable, as it showed in Fig. 2.

Figure 2: Numerical solution of magnetic field (26), with the radius $R_{0}=1$
cm of cylindrical superconductor, the London penetration depth
$\lambda=0.01R_{0}$, and Earth’s magnetic field $B_{0x}=B_{0y}=0.00005$ T. The
figure shows the Meissner effect that the superconductor placed in an uniform
magnetic field will affect the magnetic field distribution, resulting in a
gradient magnetic field around the superconductor.

##  III Sensitive measurement of the generated geomagnetic field gradient by
using the atomic interferometry

###  III.1 The motion of the cooled atom in a gradient magnetic field

As shown in Sec. II, the inhomogeneous magnetic field can be generated by the
Meissner effect of superconductor in the uniform geomagnetic background. Such
a gradient magnetic field influences the motion of the cooled free atoms with
atomic internal magnetic dipoles. The quantized Hamiltonian of atom reads

$\displaystyle\hat{H}$
$\displaystyle=\frac{\hat{\textbf{p}}^{2}}{2m}+\bm{\mu}\cdot\textbf{B}(\hat{x},\hat{y},\hat{z})+mg_{z}\hat{z}\,,$
(27)

with atomic kinetic energy $\hat{\textbf{p}}^{2}/(2m)$, magnetic moment
$\bm{\mu}$, the gravitational acceleration $g=g_{z}\approx 9.8~{}{\rm
m/s^{2}}$ along $z$-direction. To get the analytic solution of Hamiltonian
(27), we consider the lowest order gradient of magnetic field (26).

For that, we denote the atomic horizontal position
$(x=x_{0}+x^{\prime},y=y_{0}+y^{\prime})$, with $(x_{0},y_{0})$ being the
initial position of atom, and $(x^{\prime},y^{\prime})$ the dynamically
variable for horizontal motion. Furthermore, considering the London
penetration depth $\lambda\rightarrow 0$, and the radius $R_{0}$ of
cylindrical superconductor is much larger than $|x^{\prime}|$ and
$|y^{\prime}|$, we can approximately write the magnetic field (26) as follows
by neglecting the high orders of $x^{\prime}/R_{0}$ and $y^{\prime}/R_{0}$,

$\left\\{\begin{aligned}
B_{x}&\approx\left[1+\frac{(y_{0}^{2}-x_{0}^{2})}{(x_{0}^{2}+y_{0}^{2})^{2}}R_{0}^{2}+\frac{x^{\prime}}{R_{0}}\gamma_{1}+\frac{y^{\prime}}{R_{0}}\gamma_{2}\right]B_{0x}-\left[\frac{x_{0}y_{0}}{(x_{0}^{2}+y_{0}^{2})^{2}}R_{0}^{2}+\frac{x^{\prime}}{R_{0}}\gamma_{3}+\frac{y^{\prime}}{R_{0}}\gamma_{4}\right]B_{0y}\,,\\\
\\\
B_{y}&\approx\left[1+\frac{x_{0}^{2}-y_{0}^{2}}{(x_{0}^{2}+y_{0}^{2})^{2}}R_{0}^{2}-\frac{x^{\prime}}{R_{0}}\gamma_{1}-\frac{y^{\prime}}{R_{0}}\gamma_{2}\right]B_{0y}-\left[\frac{x_{0}y_{0}}{(x_{0}^{2}+y_{0}^{2})^{2}}R_{0}^{2}+\frac{x^{\prime}}{R_{0}}\gamma_{3}+\frac{y^{\prime}}{R_{0}}\gamma_{4}\right]B_{0x}\,.\end{aligned}\right.$
(28)

Here,

$\left\\{\begin{aligned}
&\gamma_{1}=-\frac{2R_{0}^{3}x_{0}}{(x_{0}^{2}+y_{0}^{2})^{2}}\left[1+2\frac{y_{0}^{2}-x_{0}^{2}}{x_{0}^{2}+y_{0}^{2}}\right]\,,\\\
\\\
&\gamma_{2}=\frac{2R_{0}^{3}y_{0}}{(x_{0}^{2}+y_{0}^{2})^{2}}\left[1-2\frac{y_{0}^{2}-x_{0}^{2}}{x_{0}^{2}+y_{0}^{2}}\right]\,,\\\
\\\
&\gamma_{3}=\frac{R_{0}^{3}y_{0}}{(x_{0}^{2}+y_{0}^{2})^{2}}\left[1-\frac{4x_{0}^{2}}{x_{0}^{2}+y_{0}^{2}}\right]\,,\\\
\\\
&\gamma_{4}=\frac{R_{0}^{3}x_{0}}{(x_{0}^{2}+y_{0}^{2})^{2}}\left[1-\frac{4y_{0}^{2}}{x_{0}^{2}+y_{0}^{2}}\right]\,.\end{aligned}\right.$
(29)

Following the magnetic field (28), the Hamiltonian (27) can be approximately
written as

$\displaystyle\hat{H}$
$\displaystyle\approx\frac{\hat{\textbf{p}}^{2}}{2m}+mg_{x}\hat{x}+mg_{y}\hat{y}+mg_{z}\hat{z}+{\rm
constant}\,,$ (30)

with

$\left\\{\begin{aligned}
g_{x}&=\frac{\mu_{b}}{mR_{0}}[(\gamma_{1}-\gamma_{3})B_{0x}-(\gamma_{1}+\gamma_{3})B_{0y}]\,,\\\
\\\
g_{y}&=\frac{\mu_{b}}{mR_{0}}[(\gamma_{2}-\gamma_{4})B_{0x}-(\gamma_{2}+\gamma_{4})B_{0y}]\,,\end{aligned}\right.$
(31)

being the horizontal accelerations induced by the lowest order magnetic
gradient, and where $\mu_{b}$ is the well-known Bohr magneton. Immediately, we
have Heisenberg operator of atomic position,

$\displaystyle\hat{\textbf{r}}(t)$
$\displaystyle=U^{\dagger}(t)\hat{\textbf{r}}U(t)=\hat{\textbf{r}}+\frac{\hat{\textbf{p}}}{m}t-\frac{1}{2}\textbf{g}t^{2}\,,$
(32)

with the well-known evolution operator $\hat{U}(t)=\exp(-it\hat{H}/\hbar)$ and
the acceleration vector $\textbf{g}=(g_{x},g_{y},g_{z})$.

Worth of note that, the Schrödinger equation with Hamiltonian (30) is also
solvable. For example, considering $x$-directional wave function of atom is
initially cooled in the Gaussian one, i.e.,

$\displaystyle\langle
x|\Psi(0)\rangle=\left(\frac{1}{2\pi\sigma_{x}^{2}}\right)^{-\frac{1}{4}}\exp\left(-\frac{x^{2}}{4\sigma_{x}^{2}}\right)\,,$
(33)

then the time-dependent wave function reads

$\displaystyle\Psi(x,t)$ $\displaystyle=\langle x|\hat{U}(t)|\Psi(0)\rangle$
(34) $\displaystyle=\psi(t)\exp\left[\frac{\hbar m}{\hbar
m+2it\sigma_{p}^{2}}\left(\frac{4xg_{x}t^{2}-x^{2}-g_{x}^{2}t^{4}}{16\sigma_{x}^{2}}-im\frac{g_{x}^{2}t^{3}-2g_{x}tx}{2\hbar}\right)\right]\,,$

with

$\displaystyle\psi(t)=\left(\frac{1}{2\pi}\right)^{-\frac{1}{4}}\sqrt{\frac{2m\sigma_{p}}{\hbar
m+2it\sigma_{p}^{2}}}\exp\left[-\frac{m^{2}g_{x}^{2}t^{2}}{4\sigma_{p}^{2}}-\frac{5img_{x}^{2}t^{3}}{6\hbar}+\frac{\hbar
m^{3}g_{x}^{2}t^{2}}{(\hbar m+2it\sigma_{p}^{2})\sigma_{p}^{2}}\right]\,.$
(35)

Above, $\sigma_{p}$ and $\sigma_{x}$ are atomic initial momentum and position
uncertainties, respectively. The solutions of wave function $\Psi(y,t)$ and
$\Psi(z,t)$ describing the $y$\- and $z$-directional motions of atom have the
similar forms as (34).

###  III.2 The MZ atom interferometer

Following Hamiltonian (30), we propose using the horizontal atom interference
to measure the lowest order gradient magnetic field generated by Meissner
effect within the background of Earth’s magnetic field
$\textbf{B}_{0}=\textbf{B}_{0x}\textbf{e}_{x}+\textbf{B}_{0y}\textbf{e}_{y}$.

Figure 3: Sketch for the MZ atom interferometer near the cross-section of
cylindrical superconductor. Three laser pulses act as the atomic beam
splittings, namely, BS1, BS2, and BS3, for realizing the matter-wave
interference and measuring the magnetic fields gradient in $(x,y)$ plane
generated by superconducting Meissner effect and the geomagnetic background.
The gravitational acceleration of atom is along $z$ direction.

The MZ atom interferometer is dramatically showed by Fig. 3, where the atomic
path is split and recombined by three short laser pulses, namely BS1, BS2,
BS3. In experiments, the BS (atom-beam splitting) can be realized by the well-
known Kapitza-Dirac (KD) scattering of atom in the pulsed standing wave of
light RMP-AKD ; KD-E . The light is large detuning with atomic internal
states, generates the AC stark effect between two internal levels of atom, but
not excite the transition between them. The AC stark effect is associated with
the external state of atomic center-of-mass motion, and changes atomic
momentum in terms of $\pm n2\hbar k$, with $k$ being the wave number of the
laser beam, and $n=0,1,2,\cdots$ the order of KD scattering. The probability
of high-order KD scattering is small, so we consider just only the lowest-
order KD scattering of atoms, i.e., $n=1$. As a consequence, the separated
distance between atom beams is on the order of $(2\hbar k/m)\times t\approx 6$
mm with $k=2\pi/(780~{}{\rm nm})$ and $t=0.5$ s.

The state of an atom moving along the two selected paths in Fig. 3 can be
written as 27

$\displaystyle|\psi_{1}\rangle=\hat{U}(t)e^{-i2\textbf{k}\cdot\textbf{r}}\hat{U}(t)e^{i2\textbf{k}\cdot\textbf{r}}|\psi_{0}\rangle\,,$
(36)

and

$\displaystyle|\psi_{2}\rangle=e^{-i2\textbf{k}\cdot\textbf{r}}\hat{U}(t)e^{i2\textbf{k}\cdot\textbf{r}}\hat{U}(t)|\psi_{0}\rangle\,.$
(37)

Here, $|\psi_{0}\rangle$ is atomic initial state, $\exp(\pm
i2\textbf{k}\cdot\textbf{r})$ is the two-photon momentum recoil operator
caused by the applied laser pulse, and $\hat{U}(t)$ the time evolution
operator of atom between the laser pulses. The interference of atomic matter-
wave is described by

$\displaystyle\langle\psi_{1}|\psi_{2}\rangle$
$\displaystyle=\langle\psi_{0}|\hat{U}^{\dagger}(t)e^{-i2\textbf{k}\cdot\textbf{r}}\hat{U}^{\dagger}(t)e^{i2\textbf{k}\cdot\textbf{r}}\hat{U}(t)e^{-i2\textbf{k}\cdot\textbf{r}}\hat{U}(t)e^{i2\textbf{k}\cdot\textbf{r}}|\psi_{0}\rangle$
(38)
$\displaystyle=\langle\psi_{0}|e^{-i2\textbf{k}\cdot\textbf{r}(t)}e^{i2\textbf{k}\cdot\textbf{r}(2t)}e^{-i2\textbf{k}\cdot\textbf{r}(t)}e^{i2\textbf{k}\cdot\textbf{r}}|\psi_{0}\rangle$
$\displaystyle=e^{i2\textbf{k}\cdot\textbf{g}t^{2}}\,.$

In the second line of this equation, we have used the property
$\hat{U}^{{\dagger}}(t)\hat{U}(t)=\hat{U}(t)\hat{U}^{{\dagger}}(t)=1$ of
unitary evolution, and used the Heisenberg operator (32) of atom moving within
the gradient magnetic field. For the third line of (38), we have used the so-
called Zassenhaus formula 28

$e^{\hat{A}+\hat{B}}=e^{\hat{A}}e^{\hat{B}}e^{-\frac{1}{2}[\hat{A},\hat{B}]}e^{\frac{1}{6}[\hat{A},[\hat{A},\hat{B}]]+\frac{1}{3}[\hat{B},[\hat{A},\hat{B}]]}\cdots\,.$
(39)

The magnetic field induced phase-shift $2\textbf{k}\cdot\textbf{g}t^{2}$in Eq.
(38) has the form same as that of the standard MZ atom interferometer for
gravitational acceleration measurements Zhu . Considering the applied laser
beam is horizontal, i.e.,
$2\textbf{k}\cdot\textbf{g}t^{2}=2(k_{x}g_{x}+k_{y}g_{y})t^{2}$, we can
directly use the experimental accuracy $\delta g_{z}$ of the gravitational
acceleration detections to make an estimation on the accuracy of magnetic
field measurements, for example, $\delta g_{x}=\delta g_{y}=\delta
g_{z}\approx 10^{-8}~{}{\rm m/s^{2}}$ with the freely falling time $t=160$ ms
Zhu . According to magnetic field induced acceleration (31), we have

$\left\\{\begin{aligned} \delta B_{0x}&=\frac{mR_{0}\delta
g_{x}}{\mu_{b}(\gamma_{1}-\gamma_{3})}+\frac{mR_{0}(\gamma_{1}+\gamma_{3})[\delta
g_{y}(\gamma_{1}-\gamma_{3})-\delta
g_{x}(\gamma_{2}-\gamma_{4})]}{2\mu_{b}(\gamma_{2}\gamma_{3}-\gamma_{1}\gamma_{4})(\gamma_{1}-\gamma_{3})}\,,\\\
\\\ \delta B_{0y}&=\frac{mR_{0}[\delta g_{y}(\gamma_{1}-\gamma_{3})-\delta
g_{x}(\gamma_{2}-\gamma_{4})]}{2\mu_{b}(\gamma_{2}\gamma_{3}-\gamma_{1}\gamma_{4})}\,.\end{aligned}\right.$
(40)

Numerically, considering the initial position $(x_{0},y_{0})=(2,2)$ cm of
atom, then $\gamma_{1}=-\gamma_{2}=\gamma_{3}=\gamma_{4}\approx 0.06$, with
atomic mass $m=10^{-25}$ kg and the radius $R_{0}=1$ cm of cylindrical
superconductor. As a consequence, the measuring accuracy for Earth’s magnetic
field is on the order of $B_{0x}\approx B_{0y}\approx 0.00005\pm 10^{-11}$ T.
Note that, the measurement-accuracy of atomic acceleration can be
significantly improved by increasing the freely falling time $t$. For example,
the experiment EEP  showed the relative acceleration between the isotope atoms
at the level of $10^{-12}g_{z}$. This means that the accuracy of geomagnetic
field measurements may reach to a level of $B_{0x}\approx B_{0y}\approx
0.00005\pm 10^{-14}$ T. In the ground-based laboratory, the freely falling
time of atoms is ultimately limited by the sizes of practical vacuum
installation. There are already plans for testing the matter-wave interference
by launching the space-mission 29 . The advantages of space experiments are
the long freely falling time of test particles and the small nongravity
noises. Thus, the future measurements for magnetic fields could reach at the
unprecedented accuracy by using atomic matter-wave interferometry.

##  IV conclusion

We have showed an approach to measure the local magnetic field changes by
using atomic matter-wave interferometry and the superconducting Meissner
effects. Firstly, we computed the magnetic field distributions inside and
outside a cylindrical superconductor within the geomagnetic field background.
Secondly, we designed a horizontal MZ atom interferometer to measure the
nonuniform magnetic field generated by the superconducting Meissner effects
and the uniform geomagnetic field background. This atom interferometer is
realized by three laser pulses which act as the atomic beam splittings but do
not excite the atomic internal states. We considered the lowest order magnetic
gradient and showed that the induced phase shift has the form same as that of
the standard MZ atom interferometer for measuring the gravitational
acceleration. Compared with the standard MZ atom interferometer, the proposed
measurement does not significantly change the experimental installation, just
place a superconductor near the path of freely falling atoms and use the
horizontal laser pulses. The scheme relies on the Meissner effect of
superconductors, so it is not only suitable for measuring the geomagnetic
field, but also suitable for measuring any magnetic field and even the
superconductivity itself such as the London penetration depth.

## References

  * (1) M. Kasevich and S. Chu, Measurement of the gravitational acceleration of an atom with a light-pulse atom interferometer, Appl. Phys. B 54, 321 (1992).
  * (2) A. Peters, K. Y. Chung, and S. Chu, Measurement of gravitational acceleration by dropping atoms, Nature (London) 400, 849 (1999).
  * (3) A. M. Nobili, A. Anselmi, and R. Pegna, Systematic errors in high-precision gravity measurements by light-pulse atom interferometry on the ground and in space, Phys. Rev. Research 2, 012036(R) (2020).
  * (4) T. L. Gustavson, P. Bouyer, and M. A. Kasevich, Precision Rotation Measurements with an Atom Interferometer Gyroscope, Phys. Rev. Lett. 78, 2046 (1997).
  * (5) Z. K. Hu, B. L. Sun, X. C. Duan, M. K. Zhou, L. L. Chen, S. Zhan, Q. Z. Zhang, and J. Luo, Demonstration of an ultrahigh-sensitivity atom-interferometry absolute gravimeter, Phys. Rev. A 88, 043610 (2013).
  * (6) J. M. McGuirk, G. T. Foster, J. B. Fixler, M. J. Snadden, and M. A. Kasevich, Sensitive absolute-gravity gradiometry using atom interferometry, Phys. Rev. A 65, 033608 (2002).
  * (7) J. B. Fixler, G. T. Foster, J. M. Mcguirk, and M. A. Kasevich, Atom Interferometer Measurement of the Newtonian Constant of Gravity, Science 315, 74 (2007).
  * (8) R. H. Parker, C. Yu, W. Zhong, B. Estey, and H. Müller, Measurement of the fine-structure constant as a test of the Standard Model, Science 360, 191 (2018).
  * (9) P. Asenbaum, C. Overstreet, M. Kim, J. Curti, and M. A. Kasevich, Atom-Interferometric Test of the Equivalence Principle at the $10^{-12}$ Level, Phys. Rev. Lett. 125, 191101 (2020).
  * (10) K. Sangster, E. A. Hinds, S. M. Barnett, and E. Riis, Measurement of the Aharonov-Casher phase in an atomic system, Phys. Rev. Lett. 71, 3641 (1993).
  * (11) M. A. Hohensee, B. Estey, P. Hamilton, A. Zeilinger, and H. Müller, Force-Free Gravitational Redshift: Proposed Gravitational Aharonov-Bohm Experiment, Phys. Rev. Lett. 108, 230404 (2012).
  * (12) S. Dimopoulos, P. W. Graham, J. M. Hogan, M. A. Kasevich, and S. Rajendran, Atomic gravitational wave interferometric sensor, Phys. Rev. D 78, 122002 (2008).
  * (13) P. Hamilton, M. Jaffe, P. Haslinger, Q. Simmons, H. Müller, J. Khoury, Atom-interferometry constraints on dark energy, Science 349, 849 (2015).
  * (14) N. Linford, P. Linford, L. Martin, A. Payne, Recent results from the English Heritage caesium magnetometer system in comparison with recent fluxgate gradiometers, Archaeol. Prospect. 14, 151 (2007).
  * (15) C. Walter, A. Braun, G. Fotopoulos, High-resolution unmanned aerial vehicle aeromagnetic surveys for mineral exploration targets, Geophys. Prospect. 68, 334 (2019).
  * (16) V. Paoletti, A. Buggi and R. Pašteka, UXO Detection by Multiscale Potential Field Methods, Pure and Applied Geophysics 176, 4363 (2019).
  * (17) S. D. Billings, C. Pasion, S. Walker, and L. Beran, Magnetic models of unexploded ordnance, IEEE Transactions on Geoscience and Remote Sensing 44, 2115 (2006).
  * (18) P. Carelli, and M. G. Castellano, High-sensitivity DC-SQUID measurements, Physica B: Condensed Matter 280, 537 (2000).
  * (19) R. L. Fagaly, Superconducting quantum interference device instruments and applications, Review of Scientific Instruments 77, 101101 (2006).
  * (20) A. Johansson, G. Sambandamurthy, D. Shahar, N. Jacobson, and R. Tenne, Nanowire Acting as a Superconducting Quantum Interference Device, Phys. Rev. Lett. 95, 116805 (2005).
  * (21) H. M. Florez, and T. Pyragius, Floquet description of optically pumped magnetometers, Phys. Rev. A 103, 033113 (2021).
  * (22) G. Oelsner, R. IJsselsteijn, T. Scholtes, A. Krüger, V. Schultze, G. Seyffert, G. Werner, M. Jäger, A. Chwala, and R. Stolz, Integrated Optically Pumped Magnetometer for Measurements within Earth’s Magnetic Field, Phys. Rev. Applied 17, 024034 (2022).
  * (23) C. Abel, et al., Optically pumped Cs magnetometers enabling a high-sensitivity search for the neutron electric dipole moment, Phys. Rev. A 101, 053419 (2020).
  * (24) W. Meissner, and R. Ochsenfeld, Ein neuer Effekt bei eintritt der Supraleitfähigkeit, Naturwissenschaften 21, 787 (1933).
  * (25) J. Bardeen, Theory of the Meissner Effect in Superconductors, Phys. Rev. 97, 1724 (1955).
  * (26) H. London, and F. London, The Electromagnetic Equations of the Supraconductor, Proceedings of the Royal Society A, 149, 71 (1935).
  * (27) R. Prozorov, R. W. Giannetta, A. Carrington, and F. M. Araujo-Moreira, Meissner-London state in superconductors of rectangular cross section in a perpendicular magnetic field, Phys. Rev. B 62, 115 (2000).
  * (28) M. C. N. Fiolhais, and H. Essén, Magnetic field expulsion from an infinite cylindrical superconductor, Physica C 497, 54 (2014).
  * (29) A. D. Cronin, J. Schmiedmayer, and D. E. Pritchard, Optics and interferometry with atoms and molecules, Rev. Mod. Phys. 81, 1051 (2009).
  * (30) R. E. Sapiro, R. Zhang, and G. Raithel, Atom interferometry using Kapitza-Dirac scattering in a magnetic trap, Phys. Rev. A 79, 043630 (2009).
  * (31) M. Zhang and Y. Wang, Nonclassical effects in geodesic motion, Phys. Rev. D 101, 026005 (2020).
  * (32) F. Casas, A. Murua, and M. Nadinic, Efficient computation of the Zassenhaus formula, Comput. Phys. Commun. 183, 2386 (2012).
  * (33) J. Williams, S. Chiow, N. Yu, and H. Müller, Quantum test of the equivalence principle and space-time aboard the International Space Station, New J. Phys. 18, 025018 (2016).


# An Empirical Determination of the Intergalactic Background Light Using NIR
Deep Galaxy Survey Data out to 5 $\mu$m and the Gamma-ray Opacity of the
Universe

Sean T. Scully  Department of Physics and Astronomy, James Madison University
Harrisonburg, VA 22807<EMAIL_ADDRESS>Los Angeles, CA90095-1547  Matthew
A. Malkan  Department of Physics and Astronomy, University of California, Los
Angeles  Los Angeles, CA90095-1547<EMAIL_ADDRESS>Harrisonburg, VA
22807<EMAIL_ADDRESS>Floyd W. Stecker  Astrophysics Science Division,
NASA/Goddard Space Flight Center  Greenbelt, MD 20771;
<EMAIL_ADDRESS>Department of Physics and Astronomy, University of
California, Los Angeles

###### Abstract

We extend our previous model-independent determination of the intergalactic
background light, based purely on galaxy survey data, out to a wavelength of 5
$\mu$m. Our approach enables us to constrain the range of photon densities,
based on the uncertainties from observationally determined luminosity
densities and colors. We further determine a 68% confidence upper and lower
limit on the opacity of the universe to $\gamma$-rays up to energies of
$1.6/(1+z)$ TeV. A comparison of our lower limit redshift-dependent opacity
curves to the opacity limits derived from the results of both ground-based air
Čerenkov telescope and Fermi LAT observations of PKS 1424+240 allows us to
place a new upper limit on the redshift of this source, independent of IBL
modeling.

diffuse radiation – galaxies: observations – gamma-rays: theory – galaxies:
PKS 1424+240

##  1 Introduction

Past work on estimating the spectral and redshift characteristics of the
intergalactic photon density (IBL) (with the $z=0$ IBL usually referred to as
the EBL), has depended on various assumptions as to the evolution of stellar
populations and dust absorption in galaxies. There have also been attempts to
probe the EBL using studies of blazar $\gamma$-ray spectra (Ackermann et al.
2012; Abramowski et al. 2013), an approach originally suggested by Stecker, De
Jager & Salamon (1992).

We previously have pursued an independent, fully empirical approach to
calculating the IBL deep using galaxy survey data (Stecker, Malkan & Scully
2012 – hereafter SMS; see also Helgason & Kashlinsky 2012). In this paper we
extend the previous SMS results on the $\gamma$-ray opacity of the Universe to
TeV energies accessible to ground-based air Čerenkov telescopes, here
including galaxy survey data extending beyond the photometric $I$ band to the
$M$ band (i.e., 0.8 $\mu$m to $\sim$ 5 $\mu$m ). We stress that the SMS
approach, being completely observationally based, is model independent,
relying only on published luminosity densities. Such an approach has a
particular advantage over model-based methods, as it enables a determination
of both the IBL and its observational uncertainties without making any
assumptions about how the galaxy luminosity functions (LFs) evolve, as some
approaches require (e.g. Malkan & Stecker 1998, 1; Kneiske, Mannheim &
Hartmann, 2002: Stecker, Malkan & Scully 2006; Franceschini et al. 2008;
Finke, Razzaque, & Dermer, 2010; Kneiske & Dole 2010) or by using semi-
analytic models that make assumptions concerning galaxy evolution, stellar
population synthesis models, star formation rates, and dust attenuation (e.g.,
Gilmore et al. 2009; Somerville et al. 2011). Domínguez et al. (2011) used the
redshift evolution of the $K$-band LFs of Cirasuolo et al. (2010), together
with model templates based on AEGIS data to model the IBL. A recent detailed
review of the problem has been given by Dwek & Krennrich (2013).

Direct integration of the galaxy survey data in the wavelength range from the
FUV to $\sim$ 5 $\mu$m has now become possible. Observational data in this
wavelength range, even out to redshifts $2<z<4$, have become sufficiently
complete to eliminate the necessity of using model-based approaches. The high-
redshift galaxies are sufficiently sampled to enable us to interpolate between
observationally determined luminosity densities obtained for specific
wavebands. We use observationally determined galaxy colors and their
associated errors to fill in the gaps in wavelength bands where there are less
data.

Here we apply the technique of SMS to determine the IBL to a wavelength of
$\sim$ 5 $\mu$m, thus extending our previous range of $\gamma$-ray opacities
of the Universe to energies up to 1.6 TeV. This range is not only relevant for
the sources with the highest redshifts observed by the Fermi-LAT telescope,
but also for results in the sub-TeV and TeV energy range made using ground-
based air Čerenkov telescopes. Our technique to determine the expected
$\gamma$-ray opacity of the Universe is thus complementary to the technique of
using $\gamma$-ray observations directly, because the intrinsic (unabsorbed)
spectra of the $\gamma$-ray sources are uncertain.

To further show the utility of our approach, we apply our lower limit
opacities to the results of both ground-based air Čerenkov telescope
observations and Fermi-LAT observations of the BL Lac object PKS 1424+240,
which allows us to place a true upper limit on the redshift of this source,
independent of IBL models. Previous observations of the ultraviolet absorption
from the IGM of this object place a lower limit on its redshift of
$z\geq.6035$ (Furniss et al. 2013).

Our results give the $\gamma$-ray opacity as a function of energy and redshift
to within an observationally determined 68% confidence band. Therefore, a
direct comparison with $\gamma$-ray spectral data will allow an independent
determination of the previously suggested effects on such spectra. These
effects include secondary $\gamma$-ray production generated by: (1) cosmic-ray
interactions along the line of sight to the source (Essey et al. 2010; Essey &
Kusenko 2012) and (2) line-of-sight photon-axion oscillations during
propagation (e.g., De Angelis et al. 2007). Both of these processes can
produce an apparent reduction of the pair-production opacity effect as derived
from $\gamma$-ray spectra alone. Thus our results can be critical in analyzing
the implications of present and future $\gamma$-ray spectral data.

##  2 Determination of the IBL from Galaxy Photon Emissivities

Galaxy luminosity functions (LFs), $\Phi_{\nu}(L)$, can be determined by
properly sampling galaxies in a survey for a given wavelength and accounting
for biases. As in SMS, in this work we chose to use only those references that
give values for the integrated LF, i.e., the luminosity density (LD),
$\rho_{L_{\nu}}(z)$. This is because there is generally a lack of knowledge
about the covariance of the errors in the Schechter function parameters used
to determine the dominant statistical errors in their analyses. Sufficient
numbers of galaxies need to be observed up and down the luminosity function
for a reliable constraint on their total luminosity density, i.e., the
integral of the LF. It is not even guaranteed that a Schechter function will
always provide the best analytic fit to the LF data. The co-moving radiation
energy density, $u_{\nu}(z)$, is derived from the equivalent co-moving
specific emissivity ${\cal E}_{\nu}(z)$ = $\rho_{L_{\nu}}(z)$.

###  2.1 Luminosity Densities

In this work, we use the same set of galaxy survey data and resulting
luminosity densities that covered rest-frame wavelengths from the far-
ultraviolet to the $I$ band as in SMS; however we now extend our calculations
to the $J,{\rm and}~{}K$ bands. We have excluded the $H$ band due to a paucity
of observational data but we have checked that our approach is consistent with
what data exists in this band. We further use the colors derived in the next
subsection to extend the calculation into the $L$ and $M$ bands (out to 4.8
$\mu$m).

The co-moving radiation energy density $u_{\nu}(z)$ is the time integral of
the co-moving specific emissivity ${\cal E}_{\nu}(z)$,

$u_{\nu}(z)=\int_{z}^{z_{\rm max}}dz^{\prime}\,{\cal
E}_{\nu^{\prime}}(z^{\prime})\frac{dt}{dz}(z^{\prime})$ (1)

where $\nu^{\prime}=\nu(1+z^{\prime})/(1+z)$ and $z_{\rm max}$ is the redshift
corresponding to initial galaxy formation (Salamon & Stecker 1998), and

$\frac{dt}{dz}{(z)}={[H_{0}(1+z)\sqrt{\Omega_{\Lambda}+\Omega_{m}(1+z)^{3}}}]^{-1},$
(2)

with $\Omega_{\Lambda}=0.72$ and $\Omega_{m}=0.28$. We note that we have
scaled all of the observational data for ${\cal E}_{\nu}(z)$ to a value of h =
0.7 for consistency.

###  2.2 Average Colors

The continuum emission from galaxies between 0.8 $\mu$m and 5 $\mu$m arises
predominantly from stellar photospheres. At these near infrared (NIR)
wavelengths, the light is mostly emitted by red giant stars. Within a few
million years of the first generation of star formation, massive stars will
have left the main sequence to begin populating the red giant branch. At later
times, red giant branch stars do not vary greatly from one stellar population
to another. The resulting spectral energy distributions (SEDs) of galaxies are
therefore much more similar to each other in the NIR range than they are at
shorter wavelengths. The scatter in red and NIR SEDs of galaxies is especially
small, since all of them contain old, red giant branch stellar populations,
which dominate the stellar mass and the continuum emission around 1 $\mu$m.

Dai et al. (2009) determine galaxy luminosity functions at 3.6, 4.5, 5.8, and
8.0 $\mu$m for a sample of 3800 – 5800 galaxies utilizing combined photometry
from the Spitzer/IRAC Shallow Survey with redshifts from the AGN and Galaxy
Evolution Survey of the NOAO Deep Wide-Field Survey Boötes field. They
obtained well defined complete luminosity functions in the local redshift bin
of $z\leq$ 0.2. They then derived galaxy luminosity densities and generated a
best-fit SED (their figure 13). We take these data to be representative of low
redshift galaxies ($z\leq$ 0.5).

To estimate the red and NIR SEDs of galaxies at higher redshifts, we have used
the template SEDs derived by Kriek et al. (2010,2011)111see also
http://astro.berkeley.edu/ mariska/comp/. Those authors utilized extensive
multi-band photometry obtained for a sample of 3500 $K$-band selected
galaxies, at redshifts between $z=0.5$ and $z=2.0$, close to a mass-limited
sample of galaxies. For redshifts between $z=0.5$ and $z=1.0$, the galaxy
stellar masses are uniformly distributed between the range of 109.7 and 1011.3
M☉. At redshifts from 1 to 2, the mass range covered by these SEDs shrinks to
1010.3 to 1011.3. Kriek et al. (2010) grouped their galaxy SEDs into 30
average templates spanning the full range from the bluest to the reddest
galaxies. These typical $L_{*}$ galaxies are the ones that produce most of the
cosmic red and NIR emission. Therefore we take the unweighted average of the
Kriek et al. (2010,2011) SEDs to represent the full sample of galaxies. This
should be reasonably accurate as there is only a 0.5 magnitude $Y-L$ color
variation (from 1.05 to 3.5 $\mu$m) between the bluest and reddest galaxies.
Therefore, over this full wavelength range, the average colors we have adopted
must be quite close to the colors of any galaxy, at least at these redshifts.
We express the average colors we have adopted as flux ratios normalized to $J$
band in Table 1.

Table 1:  Average Colors

Waveband | Relative $\nu L_{\nu}$ | Std. Deviation | Min Value | Max Value
---|---|---|---|---
$J$ (1.2 $\mu$m) | 1∗ | – | – | –
$H$ (1.6 $\mu$m) | 1.19 | 0.08 | 1.09 | 1.37
$K$ (2.2 $\mu$m) | 0.97 | 0.12 | 0.83 | 1.23
$L$ (3.5 $\mu$m) | 0.70 | 0.15 | 0.43 | 0.99
$M$ (4.8 $\mu$m) | 0.50 | 0.16 | 0.3 | 0.84

Note. — The bias in the colors towards bluer averages occurs because there are
not many galaxies in the very red categories, so the mean color stays on the
blue side of (Max + Min)/2. ∗All differences are relative to the $J$ band.

The peak in the galaxy SEDs that we used occurs around a rest wavelength of
1.6 $\mu$m. We note that there is a subtle trend for the red and NIR SEDs to
become bluer at higher redshift. This is because the younger stellar
populations were more important in the early universe. This effect produces a
small shift of the SED peak to slightly shorter rest wavelengths at higher
redshifts. This trend may continue to redshifts higher than 2. However,
because this effect is small, we assume that the SEDs of Kriek et al. (2010)
apply to all higher redshifts. We assume that these average colors also apply
to the less massive galaxies that Kriek et al. could not study. Fortunately,
the total NIR emission of galaxies tends to be dominated by the more massive
galaxies (e.g.., Ly et al 2011). Furthermore, galaxies of masses smaller by an
order of magnitude are only slightly bluer.

###  2.3 Photon Density Calculations

As in SMS, we derive a luminosity confidence band in each of our additional
wavebands by using a robust rational fitting function characterized by

$\rho_{L_{\nu}}={\cal E}_{\nu}(z)={{ax+b}\over{cx^{2}+dx+e}}$ (3)

where $x=\log(1+z)$ and $a$, $b$, $c$, $d$, and $e$ are free parameters. We
compute the 68% confidence band from Monte Carlo simulations by finding
100,000 realizations of the data and then fitting the rational function given
by equation (3). As in SMS, we model symmetric error bars with a Gaussian
distribution, while choosing a skew normal distribution to model asymmetric
errors.

The $K$ band data do not extend beyond a redshift of $z$ $\sim$ 2, while the
$J$ band data extend to $z\sim$ 3\. Thus, the fits determined using the above
method can not be trusted beyond these limits. To handle this issue, we
truncate our fits at a redshift of $z=2$, and then color transform our
previous $I$ band from SMS to fit at the higher redshifts. This is a
reasonable assumption owing to the expected similarity in the SEDs of galaxies
at these wavelengths (see previous sections). Indeed, the overlapping $J$ band
data beyond a redshift of 2 are in good agreement with this color-shifted
band.

Figure 1 shows the redshift evolution of the luminosity density for the
various wavebands based on the survey data published in the literature for the
wavelength bands that we use to extend the previous work of SMS.222References
for the values of ${\cal E}_{\nu}(z)$ used to construct Figure 1 are listed in
the legend. The upper and lower limits of the bands correspond to the highest
and lowest IBL consistent with 68% confidence errors. The error bars on the
$J$ and $K$ band data taken from Pozzetti et al. (2003) are small in
comparison to other data used. Their uncertainty was estimated strictly by
considering the range of acceptable Schechter parameters when fitting their
LFs. As other authors considered systematic errors such as cosmic variance and
extrapolation to the faint end, we have chosen to assign a linear error of 20%
to these data to bring their errors more in line with what is typically quoted
so that these data do not inordinately influence the fitting routine. The
upper and lower limits of the bands correspond to the highest and lowest IBL
consistent with 68% confidence errors.

Our method of dealing with the confidence band for redshifts beyond where
there are any data for $J$ and $K$ is not critical for the opacity
calculation. There is very little contribution to the $\gamma$-ray opacities
from photons in these wavebands beyond a redshift of 2 because of the short
time interval of the emission from galaxies at higher redshifts (see equation
(2)). To verify this, we set an upper limit by assuming that the value of the
photon density at $z=2$ remains constant out to higher redshifts. Since we
expect that the LD will drop off at the higher redshifts, this assumption
gives an upper limit. We found no appreciable differences in the calculated
opacities under this assumption (see next section) as expected.

Having determined the colors between the $J,K,L,{\rm and}~{}M$ bands, we use
them to transform the $J$ and $K$ LDs into the $L$ and $M$ bands and again
apply our fitting technique to produce confidence limits for those
wavelengths. The color-transformed data along with the confidence bands for
$L$ and $M$ can also be seen in Figure 1. We note that the wider resulting
confidence bands for $L$ and $M$ occur because we have propagated the error
determined from the colors to the $J$ and $K$ data (reflected in the size of
the error bars) before running the Monte Carlo simulations and fits. This is
as expected, since the band of acceptable photon densities derived from the
data should reflect both errors. As in SMS, we then interpolate between the
upper limit and lower limit of the confidence bands between the various
wavebands separately to find the upper and lower limit rest frame luminosity
densities.

##  3 Determination of the Optical Depth to $\gamma$-rays and resulting
$\gamma$-ray Horizon

Using the results of the previous section, the specific emissivity can then be
derived for the highest and lowest IBL LDs. The co-moving radiation energy
density is then determined from equation (1). The photon densities derived
thereby are given by

$n(\epsilon,z)=u(\epsilon,z)/\epsilon\ \ ,$ (4)

with $\epsilon=h\nu$, and with $u(\epsilon,z)$ given by equation (1).

With the co-moving photon density $n(\epsilon,z)$ evaluated, the optical depth
for $\gamma$-rays owing to electron-positron pair production interactions with
photons of the stellar radiation background can be determined from the
expression (Stecker, De Jager & Salamon 1992)

$\tau(E_{0},z_{e})=c\int_{0}^{z_{e}}dz\,\frac{dt}{dz}\int_{0}^{2}dx\,\frac{x}{2}\int_{\epsilon_{th}}^{\infty}d\epsilon\,(1+z)^{3}n(\epsilon,z)\sigma_{\gamma\gamma}[s(z)],$
(5)

In equation (5), $dt/dz$ is given by equation (2), $E_{0}$ is the observed
$\gamma$-ray energy at redshift zero, $\epsilon$ is evaluated at redshift $z$,
$z_{e}$ is the redshift of the $\gamma$-ray source, $x=(1-\cos\theta)$,
$\theta$ being the angle between the $\gamma$-ray and the soft background
photon, and the pair production cross section $\sigma_{\gamma\gamma}$ is zero
for the rest system center-of-mass energy $\sqrt{s}<2m_{e}c^{2}$, $m_{e}$
being the electron mass. Above this threshold, the pair production cross
section is given by

$\sigma_{\gamma\gamma}(s)=\frac{3}{16}\sigma_{\rm
T}(1-\beta^{2})\left[2\beta(\beta^{2}-2)+(3-\beta^{4})\ln\left(\frac{1+\beta}{1-\beta}\right)\right],$
(6)

where $\sigma_{T}$ is the Thompson scattering cross section and
$\beta=(1-4m_{e}^{2}c^{4}/s)^{1/2}$ (Breit & Wheeler 1934; Jauch & Rohrlich
1955).

As derived in SMS, the pair-production cross section energy has a threshold at
$\lambda=4.8\ \mu{\rm m}\cdot E_{\gamma}({\rm TeV})$, determined from the
energy required to produce twice the electron rest mass in the center of mass
frame. Since the maximum $\lambda$ is in the rest frame M band at 4.8 $\mu$m
at redshift $z$, the resulting energy at threshold is $\sim$ 1.6 TeV at $z=0$.
The energy at interaction in the rest frame is given by $(1+z)E_{\gamma}$
meaning the maximum $\gamma$-ray energy affected by the photon range that we
consider is $\sim 1.6(1+z)^{-1}$ TeV.

The 68% opacity ranges for $z=0.1,0.5,1,3~{}$and $5$, calculated using the SMS
methods as described above, are plotted in Figure 2. The increasing
uncertainties in the $\gamma$-ray opacity towards higher redshifts are a
reflection of the increasing widths of the uncertainty bands in the luminosity
densities shown in Figure 1, as follows from equations (4) and (5).

Figure 3 shows an energy-redshift plot of the highest energy photons from
extragalactic sources at various redshifts from Fermi as given by Abdo et al.
(2010) along with our 68% confidence band for $\tau=1$, extending our result
from SMS down to a redshift of $z=0.2$. We stress that our 68% confidence band
as shown in Figure 3 is the range of optical depths allowed as derived from
our observationally determined IBL results.

##  4 Application of our results to PKS 1424+240

To show the utility of our results, we consider redshift limits obtained from
$\gamma$-ray observations of the distant blazar PKS 1424+240. BL Lac objects
such as this one typically display bright and featureless continunous spectra
making definite spectroscopic determinations of their redshifts challenging. A
recent determination of a lower limit of the redshift of PKS 1424+240 has been
obtained by Furniss et al. (2013). They find $z\geq$ 0.6035, inferred from
Hubble observations of the measured positions of Lyman lines from intervening
Lyman $\alpha$ absorbers at lower redshifts. PKS 1424+240 has been observed at
$\gamma$-ray energies lower than 100 GeV by Fermi-LAT and at energies above
500 GeV by $VERITAS$ (Acciari et al. 2010).

Using our opacity results we can place an upper limit on the redshift of this
source. A power-law spectrum with index $\sim 1.8$ has been determined from
Fermi for energies below 100 GeV (Furniss et al. 2013). In this energy range
there should be little $\gamma$-ray attenuation. Therefore, extending this
power-law spectrum out to energies measured by $VERITAS$ should represent an
upper limit on its intrinsic spectrum. An upper limit for the opacity for each
of the measured $\gamma$-ray energies can then be determined from the flux
difference between the extended power-law spectrum and the measured fluxes.
This upper limit can then be compared with our redshift dependent opacity
curves generated from our lower limit IBL, thereby determining an upper limit
on the redshift of PKS 1424+240. Our upper limit on $z$ is found to be 1.0;
our lower limit opacities for redshifts greater than $z=1.0$ exceed the upper
limits on the opacity determined from the flux differences. Figure 4 shows the
upper limits on the opacities determined from the VERITAS data along with our
lower limit $z=1.0$ opacity curve. Aleksić et al. (2014) have used MAGIC data
to place an upper limit on the redshift of PKS 1424+240 of 0.81. However, the
limit obtained by the MAGIC collaboration depends on assuming absolute
accuracy of the IBL model of Franceschini et al. (2008). The loosening of this
limit to a slightly larger upper limit of $z=1.0$ is a result of using our
observationally derived uncertainty band rather than a theoretical model.

##  5 Comparison with Other Results

Helgason & Kashlinsky (2012) (HK) have used a method similar to ours, using
galaxy LFs to determine the $\gamma$-ray opacity of the universe with
uncertainties over the same energy range as considered here. Their results are
also model independent. However, their treatment differs from our ours in two
respects: (1) The HK uncertainty band is determined by fitting to Schechter
function parameterizations with upper and lower limits determined by fitting
to faint end slopes of the Schecter functions; we use observationally
determined LDs with observational errors combined with our Monte Carlo
treatment as described in Section 2.3. Our use of observationally determined
errors takes account of systematic effects such as cosmic variance. (2) We use
color data as described in Section 2.3 to extrapolate our LDs to higher
redshifts, allowing us to make use of more observational data, particularly in
the $L$ and $M$ bands; HK postulate an exponential cutoff in their LDs at
higher redshifts. The result of these differences is that our opacities are
somewhat larger than those of HK.

Our uncertainty band for the IBL at $z=0$ (the EBL) is in full agreement with
that obtained by Abramowski et al. (2013) using $\gamma$-ray spectra from
H.E.S.S. (see their Figure 5). Our results are consistent with the work of the
Fermi collaboration (2010) as shown in our Figure 3. Our results for $z=1$ are
in fair agreement with those of the analyses of the Fermi collaboration
(Ackermann et al. 2012). Our $\tau=1$ (horizon) band is consistent with the
result of Domínguez et al. (2013).

##  6 Conclusions

We have presented an extension of our previous determination of a 68%
confidence bands giving upper and lower limits on the IBL out to 4.8 $\mu$m.
This model-independent determination is entirely based on observationally
derived luminosity functions from local and deep galaxy survey data and color
data. This has enabled us to directly derive both the $\gamma$-ray opacity and
its observational uncertainties as a function of both energy and redshift out
to an energy of $1.6/(1+z)$ TeV for $z\leq 5$. We have applied our lower limit
opacities to the results of both ground-based air Čerenkov telescope and Fermi
LAT observations of PKS 1424+240, allowing us to place a new upper limit on
the redshift of this source, viz., $z\leq 1.0$, independent of IBL modeling.

We find no direct evidence in our spectral analysis for a required
modification of our predicted $\gamma$-ray opacities as shown in Fig. 4, and
for the opacity in the direction of PKS1424+240, either owing to axion-photon
mixing as suggested Meyer & Horns (2013) or by the existence of secondary
production effects as suggested by Essey & Kusenko (2013). However, such
effects may need to be invoked should a subsequently determined redshift of
PKS 1424+240 turn out to be greater than 1.

## Acknowledgments

We wish to thank Amy Furniss for supplying us with the VERITAS $\gamma$-ray
data on PKS1424+240 that we used in our analysis. We also thank Mariska Kriek
for supplying us with galaxy template SEDs.

## References

  * (1)
  * (2) Abdo, A., Ackermann, M., Ajello, M. et al. 2010, ApJ 723, 1082
  * (3)
  * (4) Abramowski, A., Acero, F., Aharonian, F. et al. 2013, A&A 550, A4
  * (5)
  * (6) Acciari, V. A., Aliu, E., Arlen, et al. 2010, ApJ 708, L100
  * (7)
  * (8) Ackermann, M., et al. 2012, Science 338, 1190
  * (9)
  * (10) Aleksić, S., Ansoldi, L., Antonelli, A., Antoranz, P. et al. 2014, e-print arXiv:1401.0464
  * (11)
  * (12) Arnouts, S. and Walcher, C. J., Le Fèvre, O. et al. 2007, A&A 476, 137
  * (13)
  * (14) Bell, E. F., McIntosh, D. H., Katz, N. & Weinberg, M.D. 2003, ApJ 585, 117
  * (15)
  * (16) Breit, G and Wheeler, J. A. 1934, Phys. Rev. 46, 1087
  * (17)
  * (18) Dai, X., Assef, R. J., Kochanek, C. S. et al. 2009, ApJ 697, 506
  * (19)
  * (20) De Angelis, A., Roncadelli, M. & Mansutti, O. 2007, Phys. Rev. D 76, 121301
  * (21)
  * (22) Domínguez, A., Primack, J. R., Rosario, D. J. et al. 2011, MNRAS 410, 2556
  * (23)
  * (24) A. Dom´nguez, A. J. D. Finke, J. D., Prada, F. et al. 2013, ApJ 770:77
  * (25)
  * (26) Dwek, E. & Krennrich, F. 2013, Astropart. Phys. 43, 112
  * (27)
  * (28) Essey, W., Kalashev, O. E., Kusenko, A. and Beacom, J. F. 2010, Phys. Rev. Letters 104, 141102
  * (29)
  * (30) Essey, W. & Kusenko, A. 2012, ApJ 751, L11
  * (31)
  * (32) Essey, W. & Kusenko, A. 2013, e-print arXiv:1310.3440
  * (33)
  * (34) Fazio, G. G. & Stecker, F. W. 1970, Nature 226, 135
  * (35)
  * (36) Feulner, G., Bender, R., Drory, N., Hopp, U. et al. 2003, MNRAS 342, 605
  * (37)
  * (38) Finke, J. D., Razzaque, S. & Dermer, C. D. 2010, ApJ 712, 238
  * (39)
  * (40) Franceschini, A., Rodighiero, G. & Vaccari, M. 2008, A&A, 487, 837
  * (41)
  * (42) Feulner, G., Bender, R., Drory, N., Hopp, U. et al. 2003, MNRAS 342, 605
  * (43)
  * (44) Furniss, A., Williams, D. A., Danforth, C. et al. 2013, ApJ 768, L31
  * (45)
  * (46) Gilmore, R. C. and Madau, P., Primack, J. R. et al. 2009, MNRAS 399, 1694
  * (47)
  * (48) Heath-Jones, D. Peterson, B.A., Colless, M. & Saunders, W. 2006, MNRAS 369, 25
  * (49)
  * (50) Hill, D. T., Driver, S. P., Cameron, E., et al. 2010, MNRAS 404, 1215
  * (51)
  * (52) Helgason, K. & Kashlinsky, A. 2012, ApJ 758, L13
  * (53)
  * (54) Jauch, J. M. & Rohrlich, F. 1955, The Theory of Photons and Electrons (Cambridge, MA: Addison-Wesley)
  * (55)
  * (56) Kochanek, C. S., Pahre, M. A., & Falco, E. E. 2001, ApJ 560, 566
  * (57)
  * (58) Kneiske, T. M., Mannheim, K. & Hartmann, D. H. 2002, A&A 386, 1
  * (59)
  * (60) Kneiske, T. M. & Dole, H. 2010, A&A 515, A19
  * (61)
  * (62) Kriek, M, Labbé, I., Convoy, C. et al. 2010, ApJ 722, 64
  * (63)
  * (64) Kriek, M., van Dokkum, P. G., Whitaker, K. E. et al. 2011, ApJ 743, 168
  * (65)
  * (66) Ly, C., Malkan, M. A., Hayashi, M. et al. , ApJ 735, L91
  * (67)
  * (68) Malkan, M. A. & Stecker, F. W. 1998, ApJ 496, 13
  * (69)
  * (70) Malkan, M. A. & Stecker, F. W. 2001, ApJ 555, 641
  * (71)
  * (72) Meyer, M. & Horns, D. 2013, e-print arXiv:1310.2058
  * (73)
  * (74) Pozzetti, L., Cimatti, A., Zamorani, G., Daddi, E. et al. 2003, A&A 402, 837
  * (75)
  * (76) Salamon, M. H. & Stecker, F. W. 1998, ApJ 493, 547
  * (77)
  * (78) Somerville, R. S., Gilmore, R. C., Primack, J. R. & Domínguez, A. 2012, MNRAS 423, 1992
  * (79)
  * (80) Stecker, F. W., De Jager, O.C. & Salamon, M. H. 1992, ApJ 390, L49
  * (81)
  * (82) Stecker, F.W., Malkan, M. A. & Scully, S. T. 2006, ApJ 648, 774
  * (83)
  * (84) Stecker, F.W., Malkan, M. A. & Scully, S. T. 2012, ApJ 761, 128
  * (85)
  * (86) Stefanon, M. and Marchesini, D. 2013, MNRAS 429, 881
  * (87)

Figure Captions

Figure 1: The specific emissivities for $J,K,L,{\rm and}~{}M$ wavebands. In
the $L$ and $M$ panels, the $J$ and $K$ band data have been shifted using the
color relations given in the text in order to fully determine the specific
emissivities in these wavebands. The grey shading represents the 68%
confidence bands (see text).

Figure 2: The empirically determined opacities for redshifts of 0.1, 0.5, 1,
3, 5 extended from SMS. The dashed lines indicate the opacities $\tau=1$ and
$\tau=3$.

Figure 3: A $\tau=1$ energy-redshift plot (Fazio & Stecker 1970) showing our
uncertainty band results compared with the Fermi plot of their highest energy
photons from FSRQs (red), BL Lacs (black) and and GRBs (blue) vs. redshift
(from Abdo et al. 2010).

Figure 4: Determination of the upper limit redshift for PKS 1424+240. The
points represent the upper limits to the $\gamma$-ray opacity derived from the
$VERITAS$ data. The dashed curve corresponds to our lower limit opacity for
$z=.6035$ while the solid curve is our best fit lower limit opacity to the
opacity limits derived from the results shown in Fig. 3, corresponding to an
upper limit redshift of $z=1.0$.

Figure 1.

Figure 2.

Figure 3.

Figure 4.