<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
11institutetext: Departamento de Física, Universidade do Estado do Rio Grande
do Norte, 59610–210, Mossoró – RN, Brazil 22institutetext: Departamento de
Ciências Exatas e Naturais, Universidade Federal Rural do Semi–Árido,
59900–000, Pau dos Ferros – RN, Brazil

# Constraint on the chemical potentials of hydrogen and proton in
recombination

L. L. Salese1,addr1  F. C. Carvalhoe2,addr1  H. T. C. M. Souzae4,addr3

(Received: date / Accepted: date)

###### Abstract

In this paper, we revisit the hydrogen recombination history from a novel
perspective: the evolution of chemical potentials. We derive expressions for
the chemical potentials, which depend on the thermal bath temperature and the
ionization degree of the universe. Our main finding reveals a constraint
between the chemical potentials of hydrogen and proton at $z\approx 1200$ when
the free electron fraction is $X_{e}\approx 1/3$. Furthermore, we present
important data on the chemical potentials during recombination, highlighting
the differences between the predictions of the Peebles’ and CosmoRec code
solutions. Finally, we discuss a particular case related to the chemical
potential of hydrogen.

###### Keywords:

Early universe Atomic hydrogen recombination Chemical potential

††journal: Eur. Phys. J. C

##  1 Introduction

The chemical potential is an abstract measure that relates to the propensity
for physical or chemical transformation of a species gibbs1878equilibrium . In
the field of thermodynamics, the chemical potential of a species stands for
the energy that can either be absorbed or released due to a change in the
particle number resulting from a chemical reaction or phase transition. This
physical quantity is defined as the rate of change of the free energy of a
thermodynamic system concerning the alteration in the number of atoms or
molecules of the species introduced to the system atkins2014atkins .

Recombination is a crucial phase transition that occurred over approximately
339,000 years, during which electrons and protons combined to form the neutral
hydrogen atom. As the universe expanded and cooled, the primordial plasma
underwent a process of neutralization, giving rise to the first hydrogen
atoms. At the end of this era, known as the last scattering surface, nearly
99% of the matter was combined with a significant amount of neutral hydrogen
peebles1993principles . With the matter in this state almost entirely neutral,
the remaining photons traveled freely. They dominated the energy content
beyond our Milky Way galaxy, manifesting as the Cosmic Microwave Background
Radiation (CMBR) that was first detected by radio astronomers Penzias and
Wilson in 1965 penzias1965measurement .

One systematic way of understanding the evolution of physical processes
throughout the thermal history of the universe is to study it by the Boltzmann
equation for annihilation. The reduced form of this equation is given by
dodelson2003modern

$\displaystyle\frac{df}{dt}=C[f]\;,$ (1)

with $f$ representing the distribution function and $C[f]$ referring to a
functional containing all possible collision terms between species (fermions
or bosons). In the case of cosmological recombination, with some
simplifications in $C[f]$, the equilibrium condition, known as the Saha
equation, is established.

Saha’s ionization equation is built under chemical equilibrium and assumes
that recombination occurs directly from the ground state. However, since the
recombination history is characterized by physical processes that are out of
equilibrium, the Saha equation is an inadequate tool for fully describing
these processes. In contrast, Peebles’ approach to recombination considers the
exact solution of the Boltzmann equation, which takes into account non-
equilibrium states. Peebles assumes that recombination occurs through the
capture of electrons in excited states, which then cascade down to the ground
state via the Lyman-$\alpha$ or 2s-1s two-photon transition. Peebles’s final
expression includes the probability of an excited hydrogen atom decaying via
two-photon emission peebles1968recombination , thereby releasing photons with
insufficient energy to ionize the newly formed hydrogen atoms.

This paper focuses on describing the evolution of primordial chemical
potentials, particularly during hydrogen recombination in the early universe.
Our main objective is to analyze recombination from a new perspective, in
order to identify results that are not predicted by traditional approaches. To
do so, we shall derive the chemical potentials as a function of the thermal
bath temperature and the ionization degree.

The paper is organized as follows. In section 2, we describe the universe
model adopted to study the recombination period. Details of Saha and Peebles’
approaches are presented in section 3. In section 4, we show how primordial
chemical potentials are ascertained. A special case is discussed in section 5.
Lastly, the conclusion of the paper is presented in section 6.

##  2 Universe model

In a Friedmann-Robertson-Walker (FRW) universe, where the metric is given by

$ds^{2}=c^{2}dt^{2}-a^{2}(t)\left[\frac{dr^{2}}{1-kr^{2}}-r^{2}(d\theta^{2}+sen^{2}\theta
d\phi^{2})\right]\;,$ (2)

Einstein’s equations read as

$\frac{8\pi G}{3}\rho=\frac{\dot{a}^{2}}{a^{2}}+\frac{c^{2}k}{a^{2}}\;,$ (3)
$\frac{8\pi
G}{c^{2}}p=-\frac{2{\ddot{a}}}{a}-\frac{\dot{a}^{2}}{a^{2}}-\frac{c^{2}k}{a^{2}}\,,$
(4)

where $a$ and $\dot{a}$ are the scale factor and the time derivative of the
scale factor, respectively. Here, $\rho$ is the total energy density of the
universe. Pressure $p$ is related to energy density by the state equation
$p=\omega\rho$. That being the case, the best-fit results straightforwardly to
$\omega=-1$, which corresponds to the cosmological constant, i.e., the vacuum
energy aghanim2020planck . The nonrelativistic component, baryonic and dark
matter, has pressures zero (i.e., $\omega=0$).

Supposing the various components of energy are conserved separately, from the
energy conservation equation,

$\dot{\rho}+\frac{3\dot{a}}{a}(\rho+p)=0\;,$ (5)

we found the relation between the energy density of gravitating matter and $a$
written as

$\rho_{m}=\rho_{b}+\rho_{d}\propto a^{-3}\;,$ (6)

where $\rho_{b}$ is the baryonic matter density and $\rho_{d}$ the dark matter
density. The relation between the dark energy density and the scalar factor is

$\rho_{x}\propto a^{-3(1+\omega)}\;.$ (7)

The matter density parameter,

$\Omega_{m}=\frac{8\pi G}{3H_{0}^{2}}\rho_{m}=\frac{\rho_{m}}{\rho_{0c}}\;,$
(8)

dark energy parameter,

$\Omega_{x}=\frac{8\pi G}{3H_{0}^{2}}\rho_{x}=\frac{\rho_{x}}{\rho_{0c}}\;,$
(9)

and the curvature,

$\Omega_{k}=-\frac{c^{2}k}{R_{0}^{2}H_{0}^{2}}\;,$ (10)

(where $\rho_{0c}\equiv 3H_{0}^{2}/8\pi G$ is the current critical density of
the universe) are related from Eq. (3) as:

$\Omega_{k}=1-\Omega_{m}-\Omega_{x}\;.$ (11)

From equations (3), (6) and (7) we obtain the age of the universe as follows:

$t={1\over
H_{o}}\int^{(1+z)^{-1}}_{0}\frac{dx}{(\Omega_{m}x^{-1}+\Omega_{x}x^{-(1+3\omega)}+\Omega_{k})^{1/2}}\;,$
(12)

where $x=a/a_{0}=(1+z)^{-1}$, with $a_{0}$ being the current scalar factor and
$z$ the cosmological redshift.

For estimation purposes, we adopt the latest data from the Planck Space
Telescope: $\Omega_{b}=0.0492$, $\Omega_{d}=0.264$, $h=0.674$ and $H_{0}=67.4$
km s-1 Mpc-1 aghanim2020planck . We are assuming the flat case, i.e,
$\Omega_{k}=0$.

##  3 Description of recombination

###  3.1 Saha’s approach

Developed by the Indian astrophysicist Meghnad Saha in 1920 saha1920liii , the
Saha equation plays an important role in understanding several physical
processes in the early universe. The expression relates the ionization state
of a gas in thermal equilibrium to the temperature and pressure. In order to
start Saha’s approach to recombination, let us consider the following physical
reaction:

$\displaystyle e^{-}+p\leftrightarrow H+\gamma\;,$ (13)

where $e^{-}$, $p$, $H$, and $\gamma$ stand for the electron, proton, neutral
hydrogen, and photon, respectively. Assuming that the reaction (13) remains in
equilibrium, the thermal ionization equation can be conveniently written as

$\displaystyle\frac{n_{e}n_{p}}{n_{H}}=\frac{n_{e}^{(0)}n_{p}^{(0)}}{n_{H}^{(0)}}\;,$
(14)

and each species $i$ comes from the integral over momentum space

$\displaystyle
n_{i}=g_{i}e^{\beta\mu_{i}}\int\frac{d^{3}p}{(2\pi\hbar)^{3}}e^{-\beta
E_{i}}\;,$ (15)

and

$\displaystyle n_{i}^{(0)}\equiv
g_{i}\int\frac{d^{3}p}{(2\pi\hbar)^{3}}e^{-\beta E_{i}}\;,$ (16)

with $\beta=1/k_{B}T$, where $k_{B}$ and $T$ are the Boltzmann constant and
the thermal bath temperature, respectively. Also, $m_{i}$ and $\mu_{i}$ are
the mass and chemical potential of species $i$, respectively, $\hbar$ the
reduced Planck’s constant and $c$ the speed of light in vacuum. The
superscript “($0$)” means that the particle number density is defined at
equilibrium, i.e., by configuring the chemical potential equal to zero. The
factor $g_{i}$ is the spin states number and $E_{i}=m_{i}c^{2}+p^{2}/2m_{i}$
corresponds to the energy for nonrelativistic particles ($m_{i}c^{2}\gg
k_{B}T$).

The Saha equation for recombination is given by dodelson2003modern ;
peebles1993principles ; weinberg2008cosmology

$\displaystyle\frac{X_{e}^{2}}{1-X_{e}}=\frac{1}{n_{b}}\left(\frac{m_{e}k_{B}T}{2\pi\hbar^{2}}\right)^{3/2}e^{-\epsilon_{0}/k_{B}T}~{},$
(17)

where $\epsilon_{0}=2.18\times 10^{-11}\textrm{erg}$ represents the ground
state binding energy of the hydrogen atom, and the free electron fraction is
defined as follows:

$\displaystyle X_{e}=\frac{n_{e}}{n_{b}}\;,$ (18)

with $n_{b}=n_{e}+n_{H}$ being the baryons density, ascertained from baryon to
photon ratio by

$\displaystyle n_{b}=1.12\times
10^{-5}\Omega_{b0}h^{2}(1+z)^{3}~{}\mathrm{cm}^{-3}\;.$ (19)

In the equation above, $\Omega_{b0}$ is the baryonic matter density parameter
at the present time and $h$ the reduced Hubble constant.

Saha’s equation accounts for recombination directly from the ground state.
However, during this process, highly energetic photons are produced, which
implies that newly formed neutral hydrogen atoms become ionized. As a result,
as $X_{e}$ decreases, the recombination rate decreases, making it more
difficult to maintain equilibrium dodelson2003modern ; peebles1993principles .
This leads to a rapid decrease in the ionization degree, reaching zero
quickly. Nevertheless, an asymptotic behavior around $10^{-4}$ is expected and
necessary for the emergence of the next generation of bound systems, namely
molecules peebles2009finding . Hence, the Saha equation is only useful for
estimating outcomes while maintaining chemical equilibrium. Still, as we will
conclude later, this condition will prove to be an exciting indicator for
constraining the beginning of the atomic recombination of hydrogen to
$z<1400$.

###  3.2 Peebles’ approach

The ionization rate informs us about how the recombination epoch progresses.
We use the recombination rate found by Peebles peebles1993principles ;
peebles1968recombination , which takes into account the expansion of the
universe, the photo-ionization due to the CMBR, and the photons that arise
from the recombination process:

$\displaystyle\frac{dX_{e}}{dt}=\left[\beta_{e}(1-X_{e})e^{-[(B_{1}-B_{2})/k_{B}T]}-a_{\mathrm{rec}}X_{e}^{2}n_{b}\right]C\;.$
(20)

In this expression, $B_{1}$ and $B_{2}$ are the binding energies of the ground
state and the first excited state of the hydrogen atom, respectively. The
parameter $a_{\mathrm{rec}}$ is the recombination coefficient given by

$\displaystyle a_{\mathrm{rec}}\equiv\langle\sigma v\rangle=2.84\times
10^{-11}T_{m}^{-1/2}~{}\mathrm{cm}^{3}~{}\mathrm{s}^{-1}\;,$ (21)

where $T_{m}$ is the matter temperature, ruled by the equation

$\displaystyle\frac{dT_{m}}{dt}=T_{m}\left[-2\frac{\dot{a}}{a}-\frac{\dot{X}_{e}}{3(1+X_{e})}\right]-\frac{8}{3}\frac{\sigma_{T}bT^{4}X_{e}}{m_{e}c(1+X_{e})}(T_{m}-T)~{}.$

The photo-ionization rate reads as

$\displaystyle\beta_{e}=a_{\mathrm{rec}}\left(\frac{m_{e}k_{B}T}{2\pi\hbar^{2}}\right)^{3/2}e^{-(B_{2}/k_{B}T)}\;.$
(23)

Besides, the factor $C$ is the probability of an excited hydrogen atom decay
via two photons emission ${2s}$-${1s}$, given by

$\displaystyle C=\frac{\Lambda_{2s,1s}}{\Lambda_{2s,1s}+\beta_{e}}\;\;,$ (24)

where $\Lambda_{2s,1s}=8.23~{}\mathrm{s}^{-1}$ is the decay rate of two
photons from the metastable level $2s$ to the ground state $1s$
peebles1993principles ; jones1985ionisation .

The most accurate results for the ionization history are given by numerical
codes, such as CosmoRec111Available at
http://www.jb.man.ac.uk/~jchluba/Science/CosmoRec/CosmoRec.html, which take
into account several important corrections to the cosmological hydrogen
recombination history chluba2011towards . In this study, we assume that
recombination begins when $90\%$ ($X_{e}=0.9$) of the universe’s matter is
still in an ionized state, and ends when $99\%$ has already combined, i.e.,
when $X_{e}=0.01$. Fig. 1 shows the ionization history during recombination
obtained through numerical integration of Eq. (20) and CosmoRec.

Figure 1: Degree of ionization in Peebles’ approach and with the support from
CosmoRec code

As shown in Fig. 1, Peebles’ solution agrees reasonably well with the CosmoRec
code for $z>800$. For comparison purposes, Table 1 displays important
numerical results for recombination obtained through both Peebles’ approach
and the CosmoRec code.

Table 1: Some data on the recombination epoch. The _beg_ and _end_ indices denote the beginning and end of recombination, respectively. Additionally, $\Delta z_{\mathrm{rec}}$ represents the duration of recombination in redshift, while the redshift of decoupling ($X_{e}=0.5$) is denoted by $z_{\mathrm{dec}}$ Approach | $z_{\mathrm{beg}}$ | $z_{\mathrm{end}}$ | $z_{\mathrm{dec}}$ | $\Delta z_{\mathrm{rec}}$
---|---|---|---|---
Peebles | 1444.9 | 887.1 | 1264.9 | 557.8
CosmoRec | 1453.4 | 885.3 | 1277.8 | 568.1

##  4 Evolution of primordial chemical potentials

In this section, we shall ascertain the chemical potentials of species $i$ as
a function of the ionization degree and the thermal bath temperature. This
will be performed through the utilization of the definition of fugacity. We
shall assess the evolution of primordial chemical potentials by employing both
Peebles’ and CosmoRec solutions. To commence this evaluation, first consider
the following definition:

$\frac{n_{i}}{n_{i}^{(0)}}=e^{\beta\mu_{i}}\;.$ (25)

Using the ionization fraction given by Eq. (18) together with Eq. (25), we
find the electron chemical potential

$\mu_{e}=k_{B}T\left[\ln(n_{b}X_{e})-\ln(n_{e}^{(0)})\right]\;,$ (26)

where

$n_{e}^{(0)}=2\left(\frac{m_{e}k_{B}T}{2\pi\hbar^{2}}\right)^{3/2}e^{-\beta
m_{e}c^{2}}\;.$ (27)

Employing the universe neutrality property ($n_{e}=n_{p}$) in Eq. (18), the
proton chemical potential is obtained from

$\mu_{p}=k_{B}T\left[\ln(n_{b}X_{e})-\ln(n_{p}^{(0)})\right]\;,$ (28)

where

$n_{p}^{(0)}=2\left(\frac{m_{p}k_{B}T}{2\pi\hbar^{2}}\right)^{3/2}e^{-\beta
m_{p}c^{2}}\;.$ (29)

In order to determine the hydrogen chemical potential, we resort the
definition of baryons density ($n_{b}=n_{e}+n_{H}$) in Eq. (25), so that

$\frac{n_{e}}{n_{H}^{(0)}}=\frac{n_{b}}{n_{H}^{(0)}}-e^{\beta\mu_{H}}\;\;,$
(30)

where

$n_{H}^{(0)}=4\left(\frac{m_{H}k_{B}T}{2\pi\hbar^{2}}\right)^{3/2}e^{-\beta
m_{H}c^{2}}\;.$ (31)

Using Eq. (30) in Eq. (18), we get the hydrogen chemical potential

$\mu_{H}=k_{B}T\left[\ln(1-X_{e})+\ln\alpha\right]\;\;,$ (32)

with $\alpha=n_{b}/n_{H}^{(0)}$.

For the sake of simplicity, equations (26), (28) and (32) can be written
respectively as

$\displaystyle\beta(\mu_{e}-m_{e}c^{2})=\ln(X_{e}n_{b})+\ln\left(\frac{h^{3}}{2}\right)-\frac{3}{2}\ln(2\pi
m_{e}k_{B}T)\;,$
$\displaystyle\beta(\mu_{p}-m_{p}c^{2})=\ln(X_{e}n_{b})+\ln\left(\frac{h^{3}}{2}\right)-\frac{3}{2}\ln(2\pi
m_{p}k_{B}T)\;,$

and

$\displaystyle\beta(\mu_{H}-m_{H}c^{2})$ $\displaystyle=$
$\displaystyle\ln((1-X_{e})n_{b})+\ln\left(\frac{h^{3}}{4}\right)$ (35)
$\displaystyle-$ $\displaystyle\frac{3}{2}\ln(2\pi m_{H}k_{B}T)\;.$

Based on the equations presented above, we can investigate how the chemical
potentials of the species evolve throughout the recombination history. It is
important to note that these potentials depend on the thermal bath temperature
and the free electron fraction.

To incorporate the results of the Peebles’ and CosmoRec approaches into our
analysis, we use the values of $X_{e}$ provided by the solution of Eq. (20)
and the CosmoRec code in equations (4), (4) and (35). Fig. 2 depicts the
evolution of chemical potentials as a function of redshift during the
recombination history. Note that the chemical potential lines generated via
Peebles’ solution practically overlap with those yielded through the CosmoRec
code. As shown in Fig. 1, for $z<800$, there is a slight discrepancy between
the forecasts of the two approaches. However, from the perspective of chemical
potentials, there is no significant difference observed, even for redshifts as
low as $z<800$. It is worth stressing that no proton or electron is created
and their particle chemical potentials are conserved at these late stages. One
aspect to be highlighted is that the quantity $\mu_{i}-m_{i}c^{2}$ is of the
order of hydrogen binding energy.

Figure 2: Evolution of chemical potentials during recombination. The dashed,
dotted, and solid lines stand for the electron, proton, and hydrogen,
respectively. The color black represents Peebles’ solution, and the color blue
denotes the CosmoRec solution

The main finding of our approach is underscored in Fig. 2, which shows that
during the recombination history, the curves for $\mu_{H}-m_{H}c^{2}$ and
$\mu_{p}-m_{p}c^{2}$ intersect for both Peebles’ and CosmoRec solutions. This
result imposes a noteworthy constraint on the chemical potentials of hydrogen
and proton, which can be expressed by

$\mu_{H}-\mu_{p}=8.988\times 10^{-7}\;\mathrm{erg}\;,$ (36)

where we use the constants $m_{H}=1.674\times 10^{-24}~{}\mathrm{g}$,
$m_{p}=1.673\times 10^{-24}~{}{\mathrm{g}}$ and $c=2.998\times
10^{10}~{}\mathrm{cm/s}$. By matching equations (4) and (35), we might be able
to show that the intersection occurs when the ionization degree of the
universe is

$X_{e}=\frac{1}{1+2\left(\frac{m_{H}}{m_{p}}\right)^{3/2}}~{}.$ (37)

Considering $m_{p}\approx m_{H}$, we obtain $X_{e}\approx 1/3$, corresponding
to a redshift of about $1200$. In terms of cosmic time, this corresponds to
when the universe was approximately 413,890 years old. It should be mentioned
that the traditional approach to the recombination history does not predict
this result.

Table 2 presents important data for primordial chemical potentials during the
recombination history using Peebles’ and CosmoRec code solutions. The first
column displays the two methods utilized. Columns two to seven show the
electron, proton, and hydrogen values of $y_{i}^{j}$, where we define
$y_{i}^{j}=\mu_{i}-m_{i}c^{2}$ with
$j=\mathrm{beg},\>\mathrm{end},\;\mathrm{dec}$. The eighth to last column
shows the results for the decoupling. The results presented in Table 2 only
confirm what we had already observed in Fig. 2: the difference between the
predictions according to Peebles and CosmoRec for chemical potentials is very
small.

Table 2: Data on chemical potentials in recombination. The approximate results
for $y_{i}^{j}$ have an order of magnitude of $10^{-11}$

Approach | $y_{e}^{\mathrm{beg}}$ | $y_{e}^{\mathrm{end}}$ | $y_{p}^{\mathrm{beg}}$ | $y_{p}^{\mathrm{end}}$ | $y_{H}^{\mathrm{beg}}$ | $y_{H}^{\mathrm{end}}$ | $y_{e}^{\mathrm{dec}}$ | $y_{p}^{\mathrm{dec}}$ | $y_{H}^{\mathrm{dec}}$
---|---|---|---|---|---|---|---|---|---
Peebles | -2.28 | -1.58 | -2.90 | -1.95 | -3.06 | -1.82 | -2.04 | -2.57 | -2.61
CosmoRec | -2.30 | -1.57 | -2.91 | -1.95 | -3.07 | -1.82 | -2.06 | -2.60 | -2.63

##  5 A special case: $\mu_{H}-m_{H}c^{2}$

Examining the recombination history through the lens of chemical potentials
offers distinct advantages, as it enables the exploration of aspects that may
not be perceptible using traditional approaches. In this section, we discuss
the behavior of the chemical potential of hydrogen regarding equilibrium and
non-equilibrium states. For equilibrium states, we consider Saha’s solution
for $X_{e}$, and for non-equilibrium states, we employ $X_{e}$ from the
CosmoRec code solution.

Fig. 3 portrays the evolution of chemical potentials under both scenarios.
Notably, the curves for the hydrogen chemical potential exhibit reasonable
agreement, particularly for $z<1200$. As is widely recognized, Saha’s solution
diverges from the CosmoRec solution at the onset of recombination, so that
$X_{e}$ quickly goes to zero. Thus, one might expect the graph of
$\mu_{H}-m_{H}c^{2}$ yielded by Saha’s solution to diverge from the one
provided by the CosmoRec solution soon after recombination begins, mainly
towards its end. However, as depicted in Fig. 3, this is not the case.
Therefore, we can establish the following result:

$(\mu_{H}-m_{H}c^{2})_{\mathrm{CosmoRec}}\sim(\mu_{H}-m_{H}c^{2})_{\mathrm{Saha}}~{}.$
(38)

This finding implies that the $\mu_{H}-m_{H}c^{2}$ quantity can does not
discern between equilibrium and non-equilibrium states.

Figure 3: Behavior of chemical potentials for equilibrium and non-equilibrium
states. The CosmoRec solution is represented by solid lines, while Saha’s
solution is depicted by dashed lines. Electron, proton, and hydrogen are
represented by black, blue, and red colors, respectively

On the other hand, this does not turn out for the electron and proton chemical
potentials. Specifically, $\mu_{e}-m_{e}c^{2}$ and $\mu_{p}-m_{p}c^{2}$
effectively distinguish between equilibrium and non-equilibrium states,
particularly for $z<1400$. The point at which the chemical potential curves of
the electron and proton begin to diverge identifies the threshold of the phase
transition from equilibrium to non-equilibrium regimes. This threshold is
useful for estimating the start of primordial hydrogen recombination.

##  6 Conclusion

In this paper, we have examined the hydrogen recombination history from the
perspective of chemical potentials. Our analysis has uncovered a constraint
between the hydrogen and proton chemical potentials at $z\approx 1200$ when
the free electron fraction is approximately $1/3$. Additionally, we have
presented outcomes on the chemical potentials during recombination,
highlighting the differences between the Peebles’ and CosmoRec code
predictions. It was also shown that the $\mu_{H}-m_{H}c^{2}$ measure does not
discern between equilibrium and non-equilibrium states. This finding offers a
potential avenue to simplify the description of recombination without
resorting directly to kinetic equations. This is relevant and might be further
addressed in future studies. In conclusion, our investigation has yielded
meaningful results not previously observed via traditional approaches.

###### Acknowledgements.

The authors are grateful to the Brazilian agency CAPES for financial support.
FCC was supported by CNPq/FAPERN/PRONEM.

## References

  * (1) J.W. Gibbs, Am. J. Sci. 3, 441 (1878).  DOI https://doi.org/10.2475/ajs.s3-16.96.441
  * (2) P. Atkins, P.W. Atkins, J. de Paula, _Atkins’ physical chemistry_ (Oxford university press, 2014)
  * (3) P.J.E. Peebles, _Principles of physical cosmology_ (Princeton University Press, 1993)
  * (4) A.A. Penzias, R.W. Wilson, ApJ 142, 419 (1965).  DOI 10.1086/148307
  * (5) S. Dodelson, _Modern cosmology_ (Elsevier, 2003)
  * (6) P.J.E. Peebles, ApJ 153, 1 (1968).  DOI 10.1086/149628
  * (7) N. Aghanim, Y. Akrami, M. Ashdown, J. Aumont, C. Baccigalupi, M. Ballardini, A. Banday, R. Barreiro, N. Bartolo, S. Basak, et al., A&A 641, A6 (2020).  DOI https://doi.org/10.1051/0004-6361/201833910
  * (8) M.N. Saha, Philos. Mag. 40, 472 (1920).  DOI https://doi.org/10.1080/14786441008636148
  * (9) S. Weinberg, _Cosmology_ (OUP Oxford, 2008)
  * (10) P.J.E. Peebles, L.A. Page Jr, R.B. Partridge, _Finding the big bang_ (Cambridge University Press, 2009)
  * (11) B.J. Jones, R.F. Wyse, Astron. Astrophys. 149, 144 (1985)
  * (12) J. Chluba, R.M. Thomas, MNRAS 412, 748 (2011).  DOI https://doi.org/10.1111/j.1365-2966.2010.17940.x


# SpatialCodec: Neural Spatial Speech Coding

###### Abstract

In this work, we address the challenge of encoding speech captured by a
microphone array using deep learning techniques with the aim of preserving and
accurately reconstructing crucial spatial cues embedded in multi-channel
recordings. We propose a neural spatial audio coding framework that achieves a
high compression ratio, leveraging single-channel neural sub-band codec and
SpatialCodec. Our approach encompasses two phases: (i) a neural sub-band codec
is designed to encode the reference channel with low bit rates, and (ii), a
SpatialCodec captures relative spatial information for accurate multi-channel
reconstruction at the decoder end. In addition, we also propose novel
evaluation metrics to assess the spatial cue preservation: (i) spatial
similarity, which calculates cosine similarity on a spatially intuitive
beamspace, and (ii), beamformed audio quality. Our system shows superior
spatial performance compared with high bitrate baselines and black-box neural
architecture. Demos are available at https://xzwy.github.io/SpatialCodecDemo.
Codes and models are available at https://github.com/XZWY/SpatialCodec.

Index Terms—  microphone array, audio codec, spatial audio

##  1 Introduction

Audio and speech codec aims at compressing the signals into low bitrate codes
for efficient storage or network streaming applications. Traditionally, these
coding schemes usually take advantage of some signal models and psycho-
acoustics. Speech codecs like CELP [1], SPEEX [2], OPUS [3], and EVS [4] all
use linear predictive modeling for signal analysis. Likewise, audio or music
codecs like MP3 [5] and OPUS [3], apply classic perceptual coding technique
[6] inspired by the perceptual masking effect of human hearing. In addition,
conventional quantization and entropy coding methods [7] are applied for
discretization and efficient coding respectively in these classic codecs.
However, the performance of these conventional methods suffers with very low
bit rates, e.g., at 6 kps.

While traditional codec(s) struggle to achieve high-quality, perceptually
accurate reconstruction at extremely low bitrates, neural codecs are capable
of overcoming this limitation. SoundStream [8] uses time-domain CNNs as
encoding and decoding blocks with a residual vector quantization (RVQ) to
compress the intermediate latent representations. TF-Codec [9] uses a temporal
linear predictive coding technique to further remove temporal redundancies.
More recently, Encodec [10] has been proposed which employs transformer-based
network to model the code distribution, an approach originally intended for
arithmetic coding, in order to achieve improved compression. AudioDec [11]
further proposes a two-stage training scheme such that the encoder and decoder
can be more flexible and easily switched for different applications. HiFi-
Codec [12] substitutes the RVQ block with grouped residual vector quantization
(GRVQ) and achieves better performance. These aforementioned neural approaches
are primarily derived or adapted from VQ-VAE [13] and GAN-based Vocoders [14,
15]. Thus, the encoded information from these codec(s) is also treated as
learned representations for generative tasks [16, 17, 18].

Besides single-channel codecs, spatial audio codec aims to compress multi-
channel audio while preserving the spatial information [19]. Such multi-
channel codecs are commonly designed for playback systems, or multi-speaker
systems, e.g., stereo coding [20], MPEG-3D Audio [21], MPEG-Surround [22].
Typically, a spatial audio codec adheres to a pipeline consisting of the
following steps: (i) downmix the multi-channel audio into mono or stereo, and
code with a traditional audio codec, (ii) some sub-band spatial parameters are
extracted from the multi-channel audio and coded channel-wise and band-wise,
(iii) The decoder then resynthesizes the multi-channel audio from the previous
two components. Since these codecs are designed for specific playback systems,
they do not consider microphone array spatial recording systems. Bao et al.
[23] applies this scheme to microphone array recorded speech but without any
reverberation (only the direct path exists). These methods do not use any
neural networks and also do not fully exploit inter-channel or inter-band
correlations. This means for a decent reconstruction, the system needs to have
reasonably large bands coded separately for each channel, resulting in high
coding bit rates.

Fig. 1: An Overview of the proposed SpatialCodec framework.

Our SpatialCodec aims to address the high coding rate challenge with neural
networks. Similar to the conventional methods, our design also has two
branches: the first branch codes the reference channel audio, while the second
branch codes the spatial information. Then on the decoder side, the first
branch’s decoder outputs the reconstructed reference channel. Then the second
branch’s decoder’s output and the reconstructed reference channel are used
jointly to synthesize all non-reference channels. We train and test our
proposed codec on a synthesized multi-channel spatially rich reverberant
dataset with speech from a single speaker.

We also propose several novel metrics to evaluate spatial cue preservation.
One is spatial similarity, which calculates the cosine similarity between the
estimated and ground-truth spatial features. Spatial features are designed by
beamforming towards a few fixed directions. We believe this metric is a more
intuitive metric because the spatial features are directly related to real-
world directions. We also propose to use beamforming performance as a metric
to validate our SpatialCodec’s ability to preserve both the spectral quality
and the main directivity. Our proposed SpatialCodec with 12 kbps of bitrate
performs significantly better than 96 kbps (8 channels x 12 kbps/channel) OPUS
and other channel-independent neural codecs. We also designed one black-box
model for comparison.

##  2 Problem Formulation

We consider the challenge of compressing spatial audio recordings from a
$M$-channel microphone array in a reverberant environment. Let $s(t)$ and
$h_{i}(t)$ represent the clean speech from the speaker and room impulse
response (RIR) from the speaker to the $i$-th microphone. The signals captured
by an $M$-channel microphone array, $\mathbf{x}(t)$ (termed as “mixture”) at
time ‘$t$’ is defined as:

$\mathbf{x}(t)=\big{[}h_{1}(t),....,h_{M}(t)\big{]}\ast s(t)$ (1)

where ‘$\ast$’ denotes the convolution operation, and
$\mathbf{x}(t)=\big{[}x_{1}(t),....,x_{M}(t)\big{]}$ includes speech captured
from the speaker’s direct path as well as early and late reflections using the
$M$-channel array. The goal of the proposed model is to compress $\mathbf{x}$
to a low-bitrate representation ‘$\mathbf{C}$’ such that the reconstructed
multi-channel audio $\hat{\mathbf{x}}$ at the decoder preserves all spatial
cues. Our proposed model consists of an encoder ($\Psi_{\mathrm{Enc}}$), a
quantizer ($\Psi_{\mathrm{Quant}}$), and a decoder ($\Psi_{\mathrm{Dec}}$)
which are jointly optimized such that $\hat{\mathbf{x}}$ is approximating
$\mathbf{x}$ from both spectral (perceptual quality) and spatial (direct path,
early, and late reflections) perspectives.

$\hat{\mathbf{x}}=\Psi_{\mathrm{Dec}}\big{(}\mathbf{C}\big{)};\quad\text{where}\quad\mathbf{C}=\Psi_{\mathrm{Quant}}\Big{(}\Psi_{\mathrm{Enc}}\big{(}\mathbf{x}\big{)}\Big{)}$
(2)

##  3 Method

The overall architecture of the proposed SpatialCodec is depicted in Fig.1
which comprises two main branches: (i) a single-channel sub-band codec pre-
trained to code the reference channel of the microphone array, and (ii) a
SpatialCodec that codes spatial information to reconstruct multi-channel audio
signals.

###  3.1 Single-Channel Sub-band Codec (First Branch)

The reference channel sub-band codec is a neural frequency domain sub-band
codec. The reason we are designing this instead of using existing time-domain
codecs is that the structure aligns with the SpatialCodec in the frequency
domain, which will be discussed in 3.2. The input
$x_{\text{ref}}\in\mathbb{R}^{2\times T\times F}$ is the STFT of the reference
channel audio, where $2$ corresponds to real and imaginary components. Then
the whole encoder-decoder architectures are 2D-CNNs with residual blocks
treating real-imaginary as the channel dimension. The architecture is similar
to HiFi-Codec [12] except 1-D CNN becomes 2-D CNN and downsampling time
becomes downsampling frequency.

Encoder and Decoder all have six convolutional layers, each followed by a
residual unit. For the encoder, The six layers’ kernel and stride for the time
dimension are always 3 and 1, respectively. The kernels and strides for the
frequency dimension are $[5,3,3,3,3,4]$ and $[2,2,2,2,1]$, respectively. The
output channel dimensions for all layers are $[16,32,64,128,128,256]$. We use
640-point FFT, which means the encoder would compress the frequency dimension
from 321 to 6 convolutional sub-bands. Then we code these 6 sub-bands
independently using residual vector quantization. The decoder is just the
opposite transpose convolution version of the encoder.

Each residual unit contains two residual blocks. Each block contains three 2-D
time-dilated CNN layers with skip connections. The first block’s three layers’
kernel and dilation sizes are $[(3,3),(3,5),(3,5)]$ and $[(1,1),(3,1),(5,1)]$,
respectively, in order of (time, freq). The second block’s corresponding
configurations are $[(7,3),(7,5),(7,5)]$ and $[(1,1),(3,1),(5,1)]$. Details
can be checked in the source code.

###  3.2 SpatialCodec (Second Branch)

The SpatialCodec has the same structure as the reference channel codec, except
input, output, and channel dimensions of the six convolutional layers. The
input of SpatialCodec is the reference channel STFT and spatial covariance
matrix concatenated in the channel dimension. Given $M$-channel STFT
$\mathbf{X}(t,f)\in\mathbb{C}^{M\times 1}$, the spatial covariance matrix
$\mathbf{\Phi}(t,f)\in\mathbb{C}^{M\times M}$ is defined to be:

$\mathbf{\Phi}(t,f)=\mathbf{X}(t,f)\mathbf{X}(t,f)^{\text{H}}\vspace{-4pt}$
(3)

Then the real and imaginary parts of $\mathbf{\Phi}(t,f)$ are concatenated
with the real and imaginary part of the reference channel STFT, which gives a
$2(M^{2}{+}1)$ dimensional real feature for each time-frequency bin. The
feature is then treated as the channel dimension when fed into the
SpatialCodec. The output channel dimensions for all the layers in the encoder
are $[128,128,128,128,256,256]$. Otherwise, the encoder and quantizer are the
same as the reference channel codec. For the spatial decoder, the output is
$M{-}1$ complex ratio filters (CRFs) [24]
$W_{m}(t,f)\in\mathbb{C}^{2L{+}1,2K{+}1},m\in[1,...,M{-}1]$ for all the non-
reference channels. CRFs encode the spatial relative transfer functions.
Assume the output of the reference channel STFT is
$\hat{X}_{ref}\in\mathbb{C}^{T\times F}$, then all non-reference channels for
$m\in[1,...,M{-}1]$ are:

$\hat{X}^{\text{non\\_ref}}_{m}(t,f)=\sum_{l=-L}^{L}\sum_{k=-K}^{K}W_{m}(t,f,l,k)\hat{X}_{\text{ref}}(t{+}l,f{+}k)\vspace{-3pt}$
(4)

Thus the last layer of the spatial decoder’s output channel dimension is
$2\times(2L+1)\times(2K+1)\times(M-1)$, where the first 2 means real and
imaginary.

###  3.3 Training and Loss

The training loss of the first branch mostly follows the HiFi-Codec’s [12]
strategy, using reconstruction loss, adversarial loss, and codebook learning
loss. We also add an additional time-domain SNR loss with weighting
$\lambda$=5. The SNR loss is defined as:

$L_{\text{SNR}}(x,\hat{x})\triangleq-10\text{log}_{10}(\frac{||x||^{2}}{||x-\hat{x}||^{2}})$
(5)

The second branch is trained separately with respect to the first branch.
During training, the complex ratio filter (CRF) [24] is applied to the
original reference channel audio instead of the reconstructed reference
channel from the first branch. The reason is that even after pretraining, the
first branch can only output reconstructed speech that’s perceptually
equivalent to the original speech. The underlying spectrogram or waveform does
not have an exact match. This means if we apply the CRF to the reconstructed
speech, the original non-reference channel audio cannot be used as learning
targets, because of the mismatching problem of the first branch. Thus during
training, in contrast to Eq. 4, we have

$\hat{X}^{\text{non\\_ref}}_{m}(t,f)=\sum_{l=-L}^{L}\sum_{k=-K}^{K}W_{m}(t,f,l,k)X_{\text{ref}}(t{+}l,f{+}k)\vspace{-6pt}$
(6)

where $\hat{X}_{\text{ref}}$ is substituted as the original reference channel
signal $X_{\text{ref}}$. Then we simply use the time domain SNR loss averaged
over all non-reference channels:

$\displaystyle
L_{\text{all}}=\frac{1}{M-1}\sum_{m=1}^{M{-}1}L_{\text{SNR}}(\text{ISTFT}(X^{\text{non\\_ref}}_{m}),\text{ISTFT}(\hat{X}^{\text{non\\_ref}}_{m}))$
(7)

###  3.4 Inference

In inference, we still get $\hat{X}^{\text{non\\_ref}}_{m}(t,f)$ using Eq.4.
The first branch’s sub-band codec reconstructs the reference channel audio.
Then the second branch’s SpatialCodec reconstructs $M{-}1$ complex ratio
filters, which are applied to the reconstructed reference channel audio to get
$M{-}1$ non-reference channel audio reconstructed.

##  4 Experiments and Dataset

We use the AISHELL-2 speech dataset [25] as clean speech sources to synthesize
the multi-channel reverberant dataset. Our array configuration is from a real-
world microphone array designed for the meeting scenario, which is an
8-channel linear non-uniform array where the distances between each pair of
neighbor microphones are $[2,2,2,14,2,2,2]$ in centimeters. A total of 10k
multi-channel RIRs are generated with random room characteristics using the
image-source method. The reverberation time lies in the range [0, 0.7s]. A
total of 90k, 7.5k, and 2k utterances are generated for the ‘Train’, ‘Dev’,
and ‘Test’ datasets. We use 16kHz sampling rate.

Table 1: Evaluation results for 8-channel reverberant audio reconstruction on the test set. CI denotes Channel Independent. SS denotes spatial similarity. SB-CODEC denotes our sub-band codec. DoA Error is in degrees while RTF Error is in radians. SIG, BAK, and OVRL are from DNSMOS [26] raw scores. MIMO E2E corresponds to our black-box end-to-end MIMO model as in Section 6. SYSTEM/METRICS | KBPS | CI | BEAMFORMING | SPATIAL
---|---|---|---|---
Intrusive | Non-Intrusive
SNR$\mathbf{\uparrow}$ | PESQ$\mathbf{\uparrow}$ | STOI$\mathbf{\uparrow}$ | SIG$\mathbf{\uparrow}$ | BAK$\mathbf{\uparrow}$ | OVRL$\mathbf{\uparrow}$ | DoA Err$\mathbf{\downarrow}$ | RTF Err$\mathbf{\downarrow}$ | SS$\mathbf{\uparrow}$
reverberant clean (groundtruth) | – | – | $\infty$ | 4.5 | 1 | 2.83 | 3.50 | 2.38 | 11.79 | 0 | 1
OPUS12 |  8*12=96 | Yes | 7.13 | 2.53 | 0.79 | 2.44 | 3.33 | 2.03 | 22.05 | 0.90 | 0.90
OPUS6 | 8*6=48 | 0.91 | 2.21 | 0.73 | 2.12 | 2.16 | 1.77 | 56.02 | 1.29 | 0.85
ENCODEC [10] | 3.91 | 2.56 | 0.76 | 2.50 | 3.37 | 2.08 | 56.88 | 1.24 | 0.85
HIFI-CODEC [12] | 2.36 | 2.62 | 0.77 | 2.51 | 3.39 | 2.09 | 56.88 | 1.25 | 0.86
SB-CODEC(Prop.) | 9.66 | 2.71 | 0.79 | 2.57 | 3.30 | 2.15 | 27.56 | 0.86 | 0.91
MIMO E2E(Prop.) | 12 | No | 6.14 | 3.02 | 0.85 | 2.99 | 3.80 | 2.56 | 21.67 | 0.95 | 0.86
ENCODEC+SPATIALCODEC(Prop.) | 6+6=12 | 5.01 | 2.80 | 0.80 | 2.73 | 3.61 | 2.30 | 13.78 | 0.75 | 0.95
HIFI-CODEC+SPATIALCODEC(Prop.) | 4.32 | 2.87 | 0.82 | 2.81 | 3.76 | 2.39 | 14.03 | 0.75 | 0.95
SB-CODEC+SPATIALCODEC(Prop.) | 8.28 | 2.92 | 0.82 | 2.86 | 3.72 | 2.43 | 14.83 | 0.73 | 0.95

For our models, STFT is applied with 640-point hanning window size and
320-point hop size, resulting in 50 frames per second. For the complex ratio
filter mentioned in Section 3.2, we use K=1 and L=4. For the vector
quantization in the sub-band codec and SpatialCodec in Section 3, all the
codebooks have 1024 entries. We use 2 residual vector quantization layers so
that both sub-band codec and SpatialCodec have 50*6*2*10 = 6kbps bitrate,
12kbps in total.

For training the reference channel codec, we follow the training scripts in
HiFi-Codec [12]111https://github.com/yangdongchao/AcademiCodec and set the
batch size to be 8 and the segment length to be 4 seconds. The training takes
two million steps. The SpatialCodec is trained in a similar way. We use ADAM
optimizer with $1\mathrm{e}{-4}$ learning rate.

##  5 Evaluation Metrics

This work focuses on spatial information preservation and reconstruction for
multi-channel audio. For the reference channel codec, we use HiFi-Codec [12]
and Encodec [10] as benchmark systems that can guarantee the performance of a
single-channel audio codec. We mainly evaluate our methods using several
spatial metrics and beamforming performance. We first use RTF (Relative
Transfer Function) error [27], and directional of arrival (DoA) estimation
error with the classic MUSIC [28] algorithm. Then we propose a more intuitive
spatial metric called spatial similarity. Besides the spatial metrics, we also
evaluate the target DoA’s beamforming performance to double-check both the
signal quality and the main directivity of the target source. For all metrics
in this section, we use 2048-point FFT, 512-point hop size, and 2048-point
hanning window. Below are the detailed descriptions of each metric.

###  5.1 RTF Error

The RTF (Relative Transfer Function) error is defined in [27] as the angle
between the groundtruth RTF $a(f)$ and estimated RTF $\hat{a}(f)$ averaged
over all frequency bins:

$\text{RTF\ Error}=\frac{1}{F}\sum\text{arccos\
Re}(\frac{\hat{a}(f)^{\text{H}}a(f)}{|\hat{a}(f)||a(f)|})$ (8)

where ground-truth RTF $a(f)$ and estimated RTF $\hat{a}(f)$ are extracted
from the first principal component of the $M$-channel STFT matrix
$X(f)\in\mathbb{C}^{M\times T}$ and $\hat{X}(f)\in\mathbb{C}^{M\times T}$,
respectively.

###  5.2 MUSIC DoA Error

MUSIC DoA error is the DoA estimation error corresponding to the ground-truth
DoA. This is a very intuitive metric to roughly measure the ability of the
algorithm to capture the direct path component. We use the classic MUSIC [28]
algorithm’s implementation in pyroomacoustics [29].

###  5.3 Spatial Similarity

We provide a spatially more intuitive metric which we call spatial similarity
(SS). We first define a $B$ dimensional spatial feature vector
$\textbf{S}(f)$, which are just $B$ fixed super-directive beamforming
responses’ magnitude averaged over time. Given a multi-channel input $X$, the
spatial feature is defined as:

$\displaystyle\vspace{-14pt}Y_{k}(f)$
$\displaystyle=\frac{1}{T}\sum_{t}|\text{Beamformer}_{b}(X(t,f))|,b\in[1,...,B]$
(9) $\displaystyle\boldsymbol{S}(f)$
$\displaystyle=[Y_{1}(f),Y_{2}(f),Y_{3}(f),...,Y_{B}(f)]$ (10)

The spatial feature is similar to the beamspace-domain concept in [30, 31],
which uses a few spatially sampled beamformers to transform the multi-channel
STFT from RTF domain to the beamspace domain, which matches our spatial
intuition. We uniformly sample beamformers’ directions in the inter-channel
time difference domain, which means all pairs of neighboring beamform
directions have the same inter-channel time difference from the perspective of
the linear array. The $b$-th beamformer’s direction $\theta_{b}$ follows:

$\vspace{-3pt}\theta_{b}=\text{arccos}(1-\frac{2b}{B})$ (11)

For the super-directive beamforming, we set diagonal loading to be 1e-2.
Details can be checked here [30] or our code.

Then we define the spatial similarity based on the defined spatial features.
After we get the spatial features of the original and reconstructed multi-
channel audio, $\boldsymbol{S}(f)$ and $\boldsymbol{\hat{S}}(f)$ accordingly,
the spatial similarity is defined as:

$\text{Spatial\
Similarity}=\frac{1}{F}\sum_{f}\frac{\boldsymbol{S}(f)^{T}\boldsymbol{\hat{S}}(f)}{||\boldsymbol{S}(f)||\cdot||\boldsymbol{\hat{S}}(f)||}\vspace{-5pt}$
(12)

We believe this spatial similarity is a more intuitive metric to measure
spatial cue preservation for all directions. In our case, it can capture the
direct path and early reflections’ direction of arrival. We also set B=50 when
calculating the spatial similarity.

###  5.4 Beamforming Performance

Besides beamforming towards a few fixed directions, we can also beamform
towards the ground truth direct path’s direction and check the beamforming
result. This also provides a perceptual listening opportunity to check if the
direct path’s signal is preserved well or distorted in the source direction.
We report PESQ, SNR, STOI, and the non-intrusive DNSMOS [26] raw scores (SIG,
BAK, OVRL).

##  6 Results and Discussions

Our baselines include a few channel-independent coding methods (Use single
channel codec to code each channel independently). They are OPUS [3] with two
versions (6 and 12 kbps per channel), HiFi-Codec [12], Encodec [10], and our
proposed sub-band codec (all 6 kbps per channel). HiFi-Codec and Encodec are
retrained on our dataset. We also design a black-box multi-channel in multi-
channel out end-to-end (MIMO E2E) model for comparison and analysis, which has
the same architecture as the SpatialCodec except the last decoder layer
directly outputs the multi-channel STFT reconstruction. The number of residual
vector quantization layers is set to 4 so the bitrate is 12kbps. For the
SpatialCodec, We also substitute the first branch’s sub-band codec with HiFi-
Codec and Encodec as two other baselines. Note that HiFi-Codec and Encodec are
time-domain models, so STFT is needed before applying the complex ratio
filter.
Overall Comparisons: Table 1 shows the evaluation results of all proposed and
baseline models. Generally, we can observe that our proposed SpatialCodec and
MIMO E2E model achieves much better performance in all metrics than any other
higher bitrate baselines, while maintaining only 12 kbps of bitrate in total.
Also, we can see that channel-independent coding baselines are not able to
preserve spatial information well even though each channel is coded with
decent coding rates. Last, the black-box MIMO E2E model is worse than our
proposed two-branch approach in spatial performance, e.g., spatial similarity
(SS) 0.86 vs. 0.95, which shows it’s hard to directly learn spatial
preservation through such a black-box network. However, it performs better
than SpatialCodec in beamforming performance because the main lobe of the
beamformer is generally wide and not sharply focused on the target direction.
Spatial performance: For DoA error, SpatialCodec achieves the best result when
combined with Enodec, which differs from the DoA error of the ground-truth
reverberant clean signal by only 2 degrees. When combined with sub-band codec
or Hifi-Codec, SpatialCodec still achieves similar performance. All channel-
independent baselines have over 22 degrees of DoA error, which is much worse
than our proposed SpatialCodec methods. For RTF error, our SpatialCodec has
the best performance when combined with our sub-band codec. Overall, it
achieves an RTF error of 0.73, which exceeds the OPUS12 (8*12=96 kbps)
baseline by over 18% and exceeds the OPUS6 (8*6=48 kbps) baseline by over 43%.
Again, channel-independent methods all have RTF errors higher than 0.8,
showing a comparably large error. The black-box MIMO E2E model can only
achieve 0.95 error which is even worse than OPUS12. For spatial similarity
(SS), our SpatialCodec achieves a superior score of 0.95 when combined with
any reference channel neural codec. For baseline methods, only OPUS12 and sub-
band codec have a spatial similarity of more than 0.9. The black-box MIMO E2E
model again shows inferior results here. Figure 2 shows a visualization of
B=50 normalized spatial features of four methods: original multi-channel
audio, OPUS with two bitrate versions, and our proposed SpatialCodec (with
sub-band codec). We can see that our SpatialCodec very closely aligns with the
ground-truth feature in both 1 and 3kHz frequencies, while OPUS12 and MIMO E2E
show quite some deviations from the ground-truth pattern. Overall, for spatial
performance, the proposed SpatialCodec method exceeds other baselines when
it’s combined with any reference channel codecs. Note that both the black-box
MIMO E2E model and high bitrate channel-independent coding methods are unable
to achieve similar performance with our two-branch methods.
Beamforming Performance: Our SpatialCodec and MIMO E2E all have promising
performance in both intrusive and non-intrusive beamforming metrics, showing
decent spectral reconstruction. The reconstructed multi-channel audio demos
are also given on our website. Overall, MIMO E2E performs the best
performance, probably because it incorporates all channels to code spectral
information. Then the SpatialCodec performs slightly worse but still much
better than all other baselines. For PESQ, MIMO E2E has a score of 3.02, while
the best SpatialCodec can achieve 2.92 (coupled with sub-band codec). OPUS12
can only achieve a score of 2.53 despite its high bitrate. Other neural
channel-independent baselines achieve scores around 2.65, while sub-band codec
performs a bit better. The situation is similar for STOI, with MIMO E2E
scoring 0.85 and SpatialCodec scoring 0.8+, all other baselines score below
0.8. For SNR, the sub-band codec achieves the highest 9.66 dB, while our best
SpatialCodec can achieve 8.28 dB. Then comes OPUS12’s 7.13 dB and MIMO E2E’s
6.14 dB. SNR is probably not the best metric here to show perceptual spectral
performance but it’s a good sanity check. For non-intrusive DNSMOS raw scores,
the overall score is a bit low for all systems because the target signal is
reverberant. MIMO E2E again performs the best. Together with MIMO E2E,
SpatialCodec + sub-band Codec exceeds the performance of the ground-truth
reverberant clean, because the original AISHELL2 [25] speech data is not pure
clean. However, this shows our model has promising spectral performance.
Furthermore, SpatialCodec performs significantly better than channel-
independent methods.

Fig. 2: Spatial Features Visualization (1kHz and 3kHz).

##  7 Conclusion and Future Work

This paper proposes a novel two-branch codec framework for neural spatial
speech coding. The first branch is a reference channel codec to code the
source spectral information. The second branch aims to extract and code the
spatial information for multi-channel resynthesis. We also propose a few novel
metrics to measure spatial and spectral quality, including spatial similarity
and beamforming performance. Our 12 kbps approach performs much better than
all baselines including 96 kbps OPUS12. We also designed a black-box MIMO E2E
model for comparison. Although we only test our system on single-speaker
reverberant speech signals, this approach should be able to generalize to more
complicated scenarios like multi-speaker, music sources, and moving sources.
We leave these more challenging directions for future research.

## References

  * [1] M. Schroeder and B. Atal,  “Code-excited linear prediction(celp): High-quality speech at very low bit rates,”  in ICASSP ’85. IEEE International Conference on Acoustics, Speech, and Signal Processing, 1985, vol. 10, pp. 937–940.
  * [2] Jean-Marc Valin,  “Speex: A free codec for free speech,” 2016.
  * [3] Jean-Marc Valin, Gregory Maxwell, Timothy B. Terriberry, and Koen Vos,  “High-quality, low-delay music coding in the opus codec,” 2016.
  * [4] Martin Dietz, Markus Multrus, Vaclav Eksler, Vladimir Malenovsky, Erik Norvell, Harald Pobloth, Lei Miao, Zhe Wang, Lasse Laaksonen, Adriana Vasilache, Yutaka Kamamoto, Kei Kikuiri, Stephane Ragot, Julien Faure, Hiroyuki Ehara, Vivek Rajendran, Venkatraman Atti, Hosang Sung, Eunmi Oh, Hao Yuan, and Changbao Zhu,  “Overview of the EVS codec architecture,”  in 2015 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), 2015, pp. 5698–5702.
  * [5] H.G. Musmann,  “Genesis of the MP3 audio coding standard,”  IEEE Transactions on Consumer Electronics, vol. 52, no. 3, pp. 1043–1049, 2006.
  * [6] Ning Xiang and G.M. Sessler,  Acoustics, Information, and Communication: Memorial Volume in Honor of Manfred R. Schroeder,  01 2015.
  * [7] Allen Gersho and Robert M. Gray,  Vector Quantization and Signal Compression,  Springer New York, NY, 11 1991.
  * [8] Neil Zeghidour, Alejandro Luebs, Ahmed Omran, Jan Skoglund, and Marco Tagliasacchi,  “Soundstream: An end-to-end neural audio codec,” 2021.
  * [9] Xue Jiang, Xiulian Peng, Huaying Xue, Yuan Zhang, and Yan Lu,  “Latent-domain predictive neural speech coding,”  IEEE/ACM Transactions on Audio, Speech, and Language Processing, vol. 31, pp. 2111–2123, 2023.
  * [10] Alexandre Défossez, Jade Copet, Gabriel Synnaeve, and Yossi Adi,  “High fidelity neural audio compression,” 2022.
  * [11] Yi-Chiao Wu, Israel D. Gebru, Dejan Marković, and Alexander Richard,  “Audiodec: An open-source streaming high-fidelity neural audio codec,”  in ICASSP 2023 - 2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). jun 2023, IEEE.
  * [12] Dongchao Yang, Songxiang Liu, Rongjie Huang, Jinchuan Tian, Chao Weng, and Yuexian Zou,  “HiFi-Codec: Group-residual vector quantization for high fidelity audio codec,” 2023.
  * [13] Aaron van den Oord, Oriol Vinyals, and Koray Kavukcuoglu,  “Neural discrete representation learning,” 2017.
  * [14] Kundan Kumar, Rithesh Kumar, Thibault de Boissiere, Lucas Gestin, Wei Zhen Teoh, Jose Sotelo, Alexandre de Brebisson, Yoshua Bengio, and Aaron Courville,  “MelGAN: Generative adversarial networks for conditional waveform synthesis,” 2019.
  * [15] Jungil Kong, Jaehyeon Kim, and Jaekyoung Bae,  “HiFi-GAN: Generative adversarial networks for efficient and high fidelity speech synthesis,” 2020.
  * [16] Adam Polyak, Yossi Adi, Jade Copet, Eugene Kharitonov, Kushal Lakhotia, Wei-Ning Hsu, Abdelrahman Mohamed, and Emmanuel Dupoux,  “Speech resynthesis from discrete disentangled self-supervised representations,” 2021.
  * [17] Jade Copet, Felix Kreuk, Itai Gat, Tal Remez, David Kant, Gabriel Synnaeve, Yossi Adi, and Alexandre Défossez,  “Simple and controllable music generation,” 2023.
  * [18] Chengyi Wang, Sanyuan Chen, Yu Wu, Ziqiang Zhang, Long Zhou, Shujie Liu, Zhuo Chen, Yanqing Liu, Huaming Wang, Jinyu Li, Lei He, Sheng Zhao, and Furu Wei,  “Neural codec language models are zero-shot text to speech synthesizers,” 2023.
  * [19] D.T. Yang, C. Kyriakakis, and C.C.J. Kuo,  High-fidelity Multichannel Audio Coding,  EURASIP book series on Signal Processing and Communications: European Association for Speech, Signal and Image Processing. Hindawi Publishing Corporation, 2005.
  * [20] Jeroen Breebaart, Steven Par, Armin Kohlrausch, and Erik Schuijers,  “Parametric coding of stereo audio,”  EURASIP Journal on Advances in Signal Processing, vol. 2005, 06 2005\.
  * [21] Jürgen Herre, Johannes Hilpert, Achim Kuntz, and Jan Plogsties,  “MPEG-H 3D audio—the new standard for coding of immersive spatial audio,”  IEEE Journal of Selected Topics in Signal Processing, vol. 9, no. 5, pp. 770–779, 2015.
  * [22] Johannes Hilpert and Sascha Disch,  “The MPEG surround audio coding standard [standards in a nutshell],”  IEEE Signal Processing Magazine, vol. 26, no. 1, pp. 148–52, 2009\.
  * [23] Qiang Ji and Changchun Bao,  “Multi-channel speech coding combining spatial information and speex codec,”  in 2020 15th IEEE International Conference on Signal Processing (ICSP), 2020, vol. 1, pp. 136–140.
  * [24] Wolfgang Mack and Emanuel A. P. Habets,  “Deep filtering: Signal extraction and reconstruction using complex time-frequency filters,”  IEEE Signal Processing Letters, vol. 27, pp. 61–65, 2020.
  * [25] Jiayu Du, Xingyu Na, Xuechen Liu, and Hui Bu,  “AISHELL-2: Transforming mandarin asr research into industrial scale,” 2018.
  * [26] Chandan K A Reddy, Vishak Gopal, and Ross Cutler,  “DNSMOS p.835: A non-intrusive perceptual objective speech quality metric to evaluate noise suppressors,” 2022.
  * [27] Ryan Corey and Andrew Singer,  Relative Transfer Function Estimation from Speech Keywords, pp. 238–247,  06 2018.
  * [28] Ralph O. Schmidt,  “Multiple emitter location and signal parameter estimation,”  1986\.
  * [29] Robin Scheibler, Eric Bezzam, and Ivan Dokmanic,  “Pyroomacoustics: A python package for audio room simulation and array processing algorithms,”  in 2018 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP). apr 2018, IEEE.
  * [30] Wenzhe Liu, Andong Li, Xiao Wang, Minmin Yuan, Yi Chen, Chengshi Zheng, and Xiaodong Li,  “A neural beamspace-domain filter for real-time multi-channel speech enhancement,”  Symmetry, vol. 14, no. 6, 2022.
  * [31] Andong Li, Guochen Yu, Zhongweiyang Xu, Cunhang Fan, Xiaodong Li, and Chengshi Zheng,  “TaBE: Decoupling spatial and spectral processing with taylor’s unfolding method in the beamspace domain for multi-channel speech enhancement,”  Information Fusion, vol. 101, pp. 101976, 2024.


# SOME RECENT DEVELOPMENTS IN 5TH FORCE SEARCHES

Ephraim Fischbacha,b Dennis E. Krausec,a Megan H. McDuffiea and Michael J.
Mueterthiesa aDepartment of Physics and Astronomy, Purdue University West
Lafayette, IN 47907 USA
bSnare, Inc., West Lafayette, IN 47906, USA
cDepartment of Physics, Wabash College, Crawfordsville, IN 47933, USA

###### Abstract

Recently we drew attention to the fact that most recent 5th force searches and
tests of the weak equivalence principle (WEP) utilize only one or two pairs of
test samples. We argue that, despite the great precision of these experiments,
the lack of diversity of samples may mean they are unable to detect new
composition-dependent forces, which requires the observation of a pattern in
the data. Such a pattern was observed in the experiment by Eötvös, Pekár, and
Fekete (EPF), the last high precision test of the WEP which used a significant
number of different samples. We advocate for new experiments utilizing a
sufficient number and range of samples to either confirm or refute the pattern
found in the EPF experimental data.

##  1 Introduction

Beginning with Newton and continuing through the early 20th century,
experiments searching for a composition-dependence of the force of gravity,
i.e., violations of the weak equivalence principle (WEP), utilized a rich
variety of materials.[1, 2, 3, 4] However, this is not true for modern, more
precise experiments[5, 6, 7, 8] as seen in Table 1. The signature of a
composition-dependent force would be a pattern in the acceleration data of a
WEP experiment which can only be observed when a sufficiently large number of
different samples is used, as was done in the earlier tests of the WEP. Such a
pattern, a coupling to baryon number was, in fact, observed in the data of the
most sensitive experiment which used more than three samples,[9] the torsion
balance experiment by Eötvös, Pekár, and Fekete (EPF).[3] Yet, despite a
significant subsequent experimental effort, no compelling additional evidence
for a new “fifth force” coupling to baryon number has been observed. However,
as we pointed out recently,[10] the more contemporary experiments have only
utilized one or two pairs of different materials.

Here we briefly review the work by Fischbach, et al., [9, 11] which discovered
the pattern in the EPF data, and what the data can tell us about the force
which might have caused it. We also discuss how subsequent experiments might
have failed to observed this force. As an example, we show that the EPF
experiment, and experiments by the Eöt-Wash group,[7] have different
sensitivities to new forces, this despite the fact they both use torsion
balances. Together, these observations build a compelling case for conducting
new WEP experiments with sensitivities comparable to the most recent
experiments, but which compare accelerations for a significantly broader range
of samples.[10]

Table 1: Tests of the WEP with the method, the number of different samples and pairs of samples used, and their sensitivity $\Delta\kappa=\Delta a/g$. (The number of pairs for pendulum experiments are not given since it is the same as the number of samples.) Experiment |  | Method |  | # of Samples |  | Sensitivity ($\Delta\kappa$)
---|---|---|---|---|---|---
Newton [1] (1687) |  | Pendulum |  | 9 |  | $\sim 10^{-3}$
Bessel [2] (1832) |  | Pendulum |  | 10 |  | $\sim 2\times 10^{-5}$
Eötvös [3] (1922) |  | Torsion Balance |  | 11 (10 pairs) |  | $\sim 10^{-9}$
Potter [4] (1923) |  | Pendulum |  | 7 |  | $\sim 3\times 10^{-6}$
Dicke [5] (1964) |  | Torsion Balance |  | 2 (1 pair) |  | $\sim 10^{-11}$
Braginskii [6] (1972) |  | Torsion Balance |  | 2 (1 pair) |  | $\sim 10^{-12}$
Eöt-Wash [7] (2012) |  | Torsion Balance |  | 3 (2 pairs) |  | $\sim 10^{-13}$
MICROSCOPE [8] (2017) |  | Satellite |  | 2 (1 pair) |  | $\sim 10^{-15}$

##  2 The Eötvös Pattern and Paradox

The last high-precision WEP experiment which used a significant number of
different samples was conducted by Eötvös, Pekár, and Fekete (EPF), the
results of which were published in 1922.[3] Motivated by a number of
theoretical and experimental considerations, Fischbach and colleagues
reanalyzed the EPF data, discovering evidence of a new Yukawa “fifth force”
coupling to baryon number $B$ with range $\lambda$, [9]

$V_{5}(r)\propto\frac{B_{1}B_{2}}{r}e^{-r/\lambda}.$ (1)

Specifically, when the acceleration differences relative to gravity
($\Delta\kappa=\Delta a/g$) of the samples were plotted versus their baryon-
to-mass ratio differences $\Delta(B/\mu$) ($\mu=m/m_{H}$, the sample mass
relative to atomic hydrogen), one observes that the data fall along a line
with a vanishing intercept (Fig. 1). That is, if we write

$\Delta\kappa=a\,[\Delta(B/\mu)]+b,$ (2)

a weighted least-squares fit to the EPF data finds that the slope $a$ and
intercept $b$ are given by[9, 10]

$a=(5.65\pm 0.71)\times 10^{-6},\phantom{space}b=(4.83\pm 6.44)\times
10^{-10}.$ (3)

Artificially inflating the uncertainties to force the data to agree with the
WEP only worsens the $\chi^{2}$ of the fit without affecting the slope.[10] We
are then led to the Eötvös Paradox: There is no evidence of any fault in the
experiment by EPF, or in the subsequent re-analysis of their data by
Fischbach, et al., which revealed the pattern shown in Fig. 1. Yet, no
convincing evidence of a Yukawa force coupling to baryon number has been found
in subsequent experiments.[7, 12, 13, 14, 15, 16] How is this possible?

Figure 1: The Eötvös Pattern: Data from the EPF experiment exhibit a linear
dependence of $\Delta\kappa=\Delta a/g$ on the baryon-to-mass ratio difference
$\Delta(B/\mu)$ which passes through the origin.

##  3 Towards a Resolution of the Eötvös Paradox

We take the first steps towards resolving the Eötvös paradox by noting that
the pattern shown in Fig. 1 does not give any indication of the spatial
dependence of the interaction which may be causing the pattern. While recent
experiments have set very stringent limits on new Yukawa interactions of the
form given by Eq. (1), this does not necessarily exclude other non-Yukawa
forces coupling to baryon number. The Eötvös pattern can be produced by any
force proportional to baryon number. If we write the formula for the force on
a test body as $\vec{F}_{B}=B\vec{\mathbb{F}}$, where $\vec{\mathbb{F}}$ is a
baryon force field, the difference in accelerations of two test bodies
relative to their gravitational accelerations is

$\frac{\Delta\vec{a}}{g}=\Delta\\!\left(\frac{B}{\mu}\right)\frac{\vec{\mathbb{F}}}{m_{H}g}\equiv\Delta\\!\left(\frac{B}{\mu}\right)\vec{\cal
F},$ (4)

where $\vec{\cal F}$ is the dimensionless force field. In the original fifth
force hypothesis by Fischbach et al, the source of the Eötvös pattern was the
Earth with the Yukawa interaction given by Eq. (1). Then the dimensionless
fifth force field would be [9]

${\cal F}={\cal
F}_{5}=\frac{f^{2}\epsilon(R_{\oplus}/\lambda)}{Gm^{2}_{H}}\left(\frac{B_{\oplus}}{\mu_{\oplus}}\right),$
(5)

where $f$ is the fifth force charge, $B_{\oplus}/\mu_{\oplus}\simeq 1$ is the
average baryon-to-mass ratio of the Earth, and $\epsilon(R_{\oplus}/\lambda)$
is a constant arising from integrating over the volume of the Earth.

While a large number of experiments have found no evidence for a Yukawa fifth
force, this does not exclude the possibility that another model for a fifth
force may exist that can produce the Eötvös pattern, and yet would not have
been seen in these more recent experiments. There are several reasons for
this. First, as shown above, any force with a linear coupling in $B$ will
produce the Eötvös pattern, not just a Yukawa potential. Second, none of the
newer experiments utilize more than two pairs of different samples and so
cannot convincingly reproduce the Eötvös pattern. By investing great effort in
maximizing precision at the price of using only two or three samples, it is
possible that an unexpected new force might be unintentionally zeroed-out by
misinterpreting it as an irrelevant background effect. Finally, one usually
optimizes an experiment to maximize its sensitivity to an expected signal, but
possibly at the expense of reducing its sensitivity to an unexpected signal.

To investigate the latter effect, Mueterthies compared the sensitivity of the
EPF experiment with the Eöt-Wash (EW) group’s experiment,[7] both of which use
torsion balances, to a very general dimensionless fifth force.[17]
Specifically, he assumed both EPF and EW torsion balances where acted upon by
a force whose $i$th component was given by
$F_{i}=B\left(\mathbb{F}_{i}+\sum_{j}\mathbb{D}_{ij}r_{j}\right)$, where
$\mathbb{F}_{i}$ and $\mathbb{D}_{ij}$ are components of vector and tensor
force fields, and $r_{j}$ is the $j$th component of the sample’s position
relative to the torsion balance’s pivot point. Mueterthies then showed that
the differential accelerations produced by this force on the EW and EPF
apparatuses were

$\displaystyle\Delta\kappa_{EW}$ $\displaystyle\equiv$
$\displaystyle\left(\frac{\Delta a}{g}\right)_{\rm
EW}=\sqrt{2}\Delta\\!\left(\frac{B}{\mu}\right)\sqrt{\left({\cal F}_{x}+{\cal
F}_{z}\sin\beta\right)^{2}+{\cal F}_{y}^{2}},$ (6)
$\displaystyle\Delta\kappa_{EPF}$ $\displaystyle\equiv$
$\displaystyle\left(\frac{\Delta a}{g}\right)_{\rm
EPF}=\Delta\\!\left(\frac{B}{\mu}\right)\left[\left({\cal
F}_{z}\sin\beta+{\cal F}_{x}+\frac{\nu}{\sigma}{\cal F}_{y}\right)\right.$ (7)
$\displaystyle\mbox{}-L\frac{\nu}{2\delta}({\cal D}_{yy}-{\cal
D}_{xx})+h\left({\cal D}_{xz}+\frac{\nu}{\sigma}{\cal D}_{yz}\right)$
$\displaystyle\left.-h\sin\beta\left({\cal D}_{xx}+\frac{\nu}{\sigma}{\cal
D}_{yx}-{\cal D}_{zz}\right)\right],$

respectively. Here ${\cal F}_{i}=\mathbb{F}_{i}/m_{H}g$ and ${\cal
D}_{ij}=\mathbb{D}_{ij}/m_{H}g$ are the components of the normalized force
fields, and $\beta$ is the angle of the torsion fiber relative to the
vertical. The other parameters $\nu$, $\sigma$, $L$, $h$ characterize the
dimensions of the EPF experimental setup.[17]

The differences between $\Delta\kappa_{\rm EW}$ and $\Delta\kappa_{\rm EPF}$
given by Eqs. (6) and (7) are clearly evident. While both experiments are
sensitive to the fifth force components ${\cal F}_{i}$, their functional
dependence differs significantly. Furthermore, the EPF experiment is sensitive
to a fifth force that depends on ${\cal D}_{ij}$ (force gradients), while
these terms are absent from $\Delta\kappa_{\rm EW}$. This vividly illustrates
an important point: while the EPF and EW experiments both use torsion
balances, their different configurations of test masses and methodologies make
them potentially sensitive to different signals. It is thus clearly possible
for a new force to be detected by one experiment and be missed by the other,
as demonstrated by a simple model given by Mueterthies.[17]

##  4 Conclusions

We have argued that, despite their great precision, modern WEP tests utilize
too few pairs of samples to detect a pattern such as seen in the EPF
experimental data. New high-precision WEP experiments using pairs of samples
with the widest possible range of $B/\mu$ are needed to search for, or
exclude, the Eötvös pattern.

## Acknowledgments

We thank J. T. Gruenwald and C. Y. Scarlett for useful discussions on this
work.

## References

## References

  * [1] I. Newton, The Principia (University of California Press, 1999).
  * [2] F. Bessel, Ann. Phys. (Leipzig) 101, 401 (1832).
  * [3] R. V. Eötvös, D. Pekár, and E. Fekete, Ann. Phys. (Leipzig) 68, 11 (1922); Z. Szabó, ed., Three Fundamental Papers of Loránd Eötvös (Eötvös Loránd Geophysical Institute of Hungary, 1998); É. Kilényi, The Eötvös Experiment in its Historical Context (Unicus Mühely, 2019).
  * [4] H. Potter, Proc. R. Soc. A 104, 588 (1923).
  * [5] P. G. Roll, R. V. Krotkov, R. H. Dicke, Ann. Phys. (NY) 26, 442 (1964).
  * [6] V. B. Braginskii, V. I. Panov, Sov. Phys. JETP 34, 463 (1972).
  * [7] T. A. Wagner, S. Schlamminger, J. H. Gundlach, E. G. Adelberger, Class. Quantum Grav. 29, 184002 (2012).
  * [8] P. Touboul et al., Phys. Rev. Lett. 119, 231101 (2017); P. Touboul et al., Class. Quantum Grav. 36, 225006 (2019).
  * [9] E. Fischbach, D. Sudarsky, A. Szafer, C. Talmadge, S. H. Aronson, Phys. Rev. Lett. 53, 3 (1986).
  * [10] E. Fischbach, J. T. Gruenwald, D. E. Krause, M. H. McDuffie, M. J. Mueterthies, C. Y. Scarlett, Phys. Lett. A 399, 127300 (2021).
  * [11] E. Fischbach, D. Sudarsky, A. Szafer, C. Talmadge, S. H. Aronson, Ann. Phys. (NY) 182, 1 (1988).
  * [12] G. M. Tino et al., Prog. Part. Nucl. Phys. 112, 103772 (2020).
  * [13] A. Franklin, E. Fischbach, The Rise and Fall of the Fifth Force, 2nd ed. (Springer, 2016).
  * [14] E. Fischbach, C. Talmadge, The Search for Non-Newtonian Gravity (Springer, 1999).
  * [15] E. G. Adelberger, J. H. Gundlach, B. R. Heckel, S. Hoedl, and S. Schlamminger, Prog. Part. Nucl. Phys. 62, 102 (2009).
  * [16] C. M. Will, Theory and Experiment in Gravitational Physics, 2nd ed. (Cambridge University Press, 2018).
  * [17] M. J. Mueterthies, et al., arXiv:2001.03752; M. J. Mueterthies, Ph.D. thesis, Purdue University, 2019.