* Dimastrogiovanni  _et al._ (2023) E. Dimastrogiovanni, M. Fasiello, A. Malhotra, and G. Tasinato, “Enhancing gravitational wave anisotropies with peaked scalar sources,” JCAP 01, 018 (2023), arXiv:2205.05644 [astro-ph.CO] .
  * Zheng  _et al._ (2023) Y. Zheng, N. Kouvatsos, J. Golomb, M. Cavaglià, A. I. Renzini, and M. Sakellariadou, “The angular power spectrum of gravitational-wave transient sources as a probe of the large-scale structure,” (2023), arXiv:2305.02652 [astro-ph.CO] .
  * Ballmer (2006) S. W. Ballmer, “A Radiometer for stochastic gravitational waves,” Class. Quant. Grav. 23, S179–S186 (2006), arXiv:gr-qc/0510096 .
  * Thrane  _et al._ (2009) E. Thrane, S. Ballmer, J. D. Romano, S. Mitra, D. Talukder, S. Bose, and V. Mandic, “Probing the anisotropies of a stochastic gravitational-wave background using a network of ground-based laser interferometers,” Phys. Rev. D 80, 122002 (2009), arXiv:0910.0858 [astro-ph.IM] .
  * Abbott  _et al._ (2021d) R. Abbott _et al._ (LIGO, Virgo, KAGRA), “Search for anisotropic gravitational-wave backgrounds using data from Advanced LIGO and Advanced Virgo’s first three observing runs,” Phys. Rev. D 104, 022005 (2021d), arXiv:2103.08520 [gr-qc] .
  * Tsukada  _et al._ (2023) L. Tsukada, S. Jaraba, D. Agarwal, and E. Floden, “Bayesian parameter estimation for targeted anisotropic gravitational-wave background,” Phys. Rev. D 107, 023024 (2023), arXiv:2208.14421 [astro-ph.IM] .
  * Banagiri  _et al._ (2021) S. Banagiri, A. Criswell, T. Kuan, V. Mandic, J. D. Romano, and S. R. Taylor, “Mapping the gravitational-wave sky with LISA: a Bayesian spherical harmonic approach,” mnras 507, 5451–5462 (2021), arXiv:2103.00826 [astro-ph.IM] .
  * Abadie  _et al._ (2011) J. Abadie _et al._ (LIGO), “Directional limits on persistent gravitational waves using LIGO S5 science data,” Phys. Rev. Lett. 107, 271102 (2011), arXiv:1109.1809 [astro-ph.CO] .
  * Abbott  _et al._ (2017e) B. P. Abbott _et al._ (LIGO, Virgo), “Directional Limits on Persistent Gravitational Waves from Advanced LIGO’s First Observing Run,” Phys. Rev. Lett. 118, 121102 (2017e), arXiv:1612.02030 [gr-qc] .
  * Thrane and Talbot (2019) E. Thrane and C. Talbot, “An introduction to Bayesian inference in gravitational-wave astronomy: Parameter estimation, model selection, and hierarchical models,” pasa 36, e010 (2019), arXiv:1809.02293 [astro-ph.IM] .
  * Abbott  _et al._ (2021e) R. Abbott _et al._ (LIGO, Virgo, KAGRA), “Upper limits on the isotropic gravitational-wave background from Advanced LIGO and Advanced Virgo’s third observing run,” Phys. Rev. D 104, 022004 (2021e), arXiv:2101.12130 [gr-qc] .
  * Ade  _et al._ (2016) P. A. R. Ade _et al._ (Planck), “Planck 2015 results. XIII. Cosmological parameters,” Astron. Astrophys. 594, A13 (2016), arXiv:1502.01589 [astro-ph.CO] .
  * Romano and Cornish (2017) J. D. Romano and N. J. Cornish, “Detection methods for stochastic gravitational-wave backgrounds: a unified treatment,” Living Rev. Rel. 20, 2 (2017), arXiv:1608.06889 [gr-qc] .
  * Allen and Ottewill (1997) B. Allen and A. C. Ottewill, “Detection of anisotropies in the gravitational wave stochastic background,” Phys. Rev. D 56, 545–563 (1997), arXiv:gr-qc/9607068 .
  * Allen and Romano (1999) B. Allen and J. D. Romano, “Detecting a stochastic background of gravitational radiation: Signal processing strategies and sensitivities,” Phys. Rev. D 59, 102001 (1999), arXiv:gr-qc/9710117 .
  * Maggiore (2008) M. Maggiore,  _Gravitational Waves: Volume 1: Theory and Experiments_, Gravitational Waves (OUP Oxford, 2008).
  * Abbott  _et al._ (2019c) B. P. Abbott _et al._ (LIGO, Virgo), “Directional limits on persistent gravitational waves using data from Advanced LIGO’s first two observing runs,” Phys. Rev. D 100, 062001 (2019c), arXiv:1903.08844 [gr-qc] .
  * Biwer  _et al._ (2016) C. Biwer _et al._ (LIGO), “Validating gravitational-wave detections: The advanced ligo hardware injection system,” https://dcc.ligo.org/LIGO-P1600285/public (2016).
  * Lazzarini  _et al._ (2004) A. Lazzarini, S. Bose, P. Fritschel, M. McHugh, T. Regimbau, K. Reilly, J. D. Romano, J. T. Whelan, S. E. Whitcomb, and B. F. Whiting, “Optimal combination of signals from co-located gravitational wave interferometers for use in searches for a stochastic background,” Phys. Rev. D 70, 062001 (2004), arXiv:gr-qc/0403093 .
  * Lazzarini and Romano (2004) A. Lazzarini and J. Romano (LIGO), “Use of overlapping windows in the stochastic background search,” https://dcc.ligo.org/T040089/public (2004).
  * Lazzarini (2004) A. Lazzarini (LIGO), “Bias from power spectrum measurement in parameter estimation for the stochastic gravitational wave background,” https://dcc.ligo.org/LIGO-T040128/public (2004).
  * Chung  _et al._ (2022) A. K.-W. Chung, A. C. Jenkins, J. D. Romano, and M. Sakellariadou, “Targeted search for the kinematic dipole of the gravitational-wave background,” Phys. Rev. D 106, 082005 (2022), arXiv:2208.01330 [gr-qc] .
  * Matas and Romano (2021) A. Matas and J. D. Romano, “Frequentist versus Bayesian analyses: Cross-correlation as an approximate sufficient statistic for LIGO-Virgo stochastic background searches,” Phys. Rev. D 103, 062003 (2021), arXiv:2012.00907 [gr-qc] .
  * Bernstein (2005) D. Bernstein,  _Matrix Mathematics: Theory, Facts, and Formulas with Application to Linear Systems Theory_ (Princeton University Press, 2005).
  * Mitra  _et al._ (2008) S. Mitra, S. Dhurandhar, T. Souradeep, A. Lazzarini, V. Mandic, S. Bose, and S. Ballmer, “Gravitational wave radiometry: Mapping a stochastic gravitational wave background,” Phys. Rev. D 77, 042002 (2008), arXiv:0708.2728 [gr-qc] .
  * Biscoveanu  _et al._ (2020a) S. Biscoveanu, C. Talbot, E. Thrane, and R. Smith, “Measuring the primordial gravitational-wave background in the presence of astrophysical foregrounds,” Phys. Rev. Lett. 125, 241101 (2020a), arXiv:2009.04418 [astro-ph.HE] .
  * Floden  _et al._ (2022a) E. Floden, V. Mandic, A. Matas, and L. Tsukada, “Angular resolution of the search for anisotropic stochastic gravitational-wave background with terrestrial gravitational-wave detectors,” Phys. Rev. D 106, 023010 (2022a), arXiv:2203.17141 [astro-ph.IM] .
  * Panda  _et al._ (2019) S. Panda, S. Bhagwat, J. Suresh, and S. Mitra, “Stochastic gravitational wave background mapmaking using regularized deconvolution,” Phys. Rev. D 100, 043541 (2019), arXiv:1905.08276 [gr-qc] .
  * Suresh  _et al._ (2021) J. Suresh, A. Ain, and S. Mitra, “Unified mapmaking for an anisotropic stochastic gravitational wave background,” Phys. Rev. D 103, 083024 (2021), arXiv:2011.05969 [gr-qc] .
  * Agarwal  _et al._ (2021) D. Agarwal, J. Suresh, S. Mitra, and A. Ain, “Upper limits on persistent gravitational waves using folded data and the full covariance matrix from Advanced LIGO’s first two observing runs,” Phys. Rev. D 104, 123018 (2021), arXiv:2105.08930 [gr-qc] .
  * Floden  _et al._ (2022b) E. Floden, V. Mandic, A. Matas, and L. Tsukada, “Angular resolution of the search for anisotropic stochastic gravitational-wave background with terrestrial gravitational-wave detectors,” Phys. Rev. D 106, 023010 (2022b), arXiv:2203.17141 [astro-ph.IM] .
  * Agarwal  _et al._ (2023) D. Agarwal, J. Suresh, S. Mitra, and A. Ain, “Angular power spectra of anisotropic stochastic gravitational wave background: developing statistical methods and analyzing data from ground-based detectors,” (2023), arXiv:2302.12516 [gr-qc] .
  * Bartolo  _et al._ (2022) N. Bartolo _et al._ (LISA Cosmology Working Group), “Probing anisotropies of the Stochastic Gravitational Wave Background with LISA,” JCAP 11, 009 (2022), arXiv:2201.08782 [astro-ph.CO] .
  * Chowdhury  _et al._ (2022) D. Chowdhury, G. Tasinato, and I. Zavala, “On the response of the Einstein Telescope to Doppler anisotropies,” (2022), arXiv:2209.05770 [gr-qc] .
  * Zonca  _et al._ (2019) A. Zonca, L. Singer, D. Lenz, M. Reinecke, C. Rosset, E. Hivon, and K. Gorski, “healpy: equal area pixelization and spherical harmonics transforms for data on the sphere in python,” Journal of Open Source Software 4, 1298 (2019).
  * Górski  _et al._ (2005) K. M. Górski, E. Hivon, A. J. Banday, B. D. Wandelt, F. K. Hansen, M. Reinecke, and M. Bartelmann, “HEALPix: A Framework for High-Resolution Discretization and Fast Analysis of Data Distributed on the Sphere,” Astrophys. J. 622, 759–771 (2005), arXiv:astro-ph/0409513 .
  * Hall and Evans (2019) E. D. Hall and M. Evans, “Metrics for next-generation gravitational-wave detectors,” arXiv e-prints , arXiv:1902.09485 (2019), arXiv:1902.09485 [astro-ph.IM] .
  * Ain  _et al._ (2018) A. Ain, J. Suresh, and S. Mitra, “Very fast stochastic gravitational wave background map making using folded data,” Phys. Rev. D 98, 024001 (2018), arXiv:1803.08285 [gr-qc] .
  * Romero  _et al._ (2021b) A. Romero, K. Martinovic, T. A. Callister, H.-K. Guo, M. Martínez, M. Sakellariadou, F.-W. Yang, and Y. Zhao, “Implications for first-order cosmological phase transitions from the third ligo-virgo observing run,” Phys. Rev. Lett. 126, 151301 (2021b).
  * Reitze  _et al._ (2019) D. Reitze _et al._ , “Cosmic Explorer: The U.S. Contribution to Gravitational-Wave Astronomy beyond LIGO,” Bull. Am. Astron. Soc. 51, 035 (2019), arXiv:1907.04833 [astro-ph.IM] .
  * Punturo  _et al._ (2010) M. Punturo _et al._ , “The Einstein Telescope: A third-generation gravitational wave observatory,” Class. Quant. Grav. 27, 194002 (2010).
  * Abbott  _et al._ (2004) B. Abbott _et al._ (LIGO), “Analysis of first LIGO science data for stochastic gravitational waves,” Phys. Rev. D 69, 122004 (2004), arXiv:gr-qc/0312088 .
  * Abbott  _et al._ (2019d) R. Abbott _et al._ (LIGO), “Data for ‘Search for the isotropic stochastic background using data from Advanced LIGO’s second observing run’,” https://dcc.ligo.org/T1900058/public (2019d).
  * Abbott  _et al._ (2021f) R. Abbott _et al._ (LIGO, Virgo, KAGRA), “Data products and supplemental information for O3 stochastic directional paper,” https://dcc.ligo.org/LIGO-G2002165/public (2021f).
  * Abbott  _et al._ (2016d) B. P. Abbott _et al._ (LIGO, Virgo), “Characterization of transient noise in Advanced LIGO relevant to gravitational wave signal GW150914,” Class. Quant. Grav. 33, 134001 (2016d), arXiv:1602.03844 [gr-qc] .
  * Thrane  _et al._ (2013) E. Thrane, N. Christensen, and R. Schofield, “Correlated magnetic noise in global networks of gravitational-wave interferometers: observations and implications,” Phys. Rev. D 87, 123009 (2013), arXiv:1303.2613 [astro-ph.IM] .
  * Thrane  _et al._ (2014) E. Thrane, N. Christensen, R. M. S. Schofield, and A. Effler, “Correlated noise in networks of gravitational-wave detectors: subtraction and mitigation,” Phys. Rev. D 90, 023013 (2014), arXiv:1406.2367 [astro-ph.IM] .
  * Smith and Thrane (2018) R. Smith and E. Thrane, “Optimal Search for an Astrophysical Gravitational-Wave Background,” Phys. Rev. X 8, 021019 (2018), arXiv:1712.00688 [gr-qc] .
  * Biscoveanu  _et al._ (2020b) S. Biscoveanu, C. Talbot, E. Thrane, and R. Smith, “Measuring the primordial gravitational-wave background in the presence of astrophysical foregrounds,” Phys. Rev. Lett. 125, 241101 (2020b), arXiv:2009.04418 [astro-ph.HE] .


# Probing Modified Gravity with Entanglement of Microspheres

Ankit Kumar<EMAIL_ADDRESS>Department of Physics, Indian
Institute of Technology Roorkee, Roorkee 247667, India  International Centre
for Theory of Quantum Technologies, University of Gdańsk, 80-308 Gdańsk,
Poland  Yen-Kheng Lim  School of Mathematics and Physics, Xiamen University
Malaysia, 43900 Sepang, Malaysia  P. Arumugam  Department of Physics, Indian
Institute of Technology Roorkee, Roorkee 247667, India  Centre for Photonics
and Quantum Communication Technology, Indian Institute of Technology Roorkee,
Roorkee 247667, India  Tom Złośnik  Institute of Theoretical Physics and
Astrophysics, Faculty of Mathematics, Physics and Informatics, University of
Gdańsk, 80-308 Gdańsk, Poland  Tomasz Paterek  School of Mathematics and
Physics, Xiamen University Malaysia, 43900 Sepang, Malaysia  Institute of
Theoretical Physics and Astrophysics, Faculty of Mathematics, Physics and
Informatics, University of Gdańsk, 80-308 Gdańsk, Poland

###### Abstract

While a wide variety of astrophysical and cosmological phenomena suggest the
presence of Dark Matter, all evidence remains via its gravitational effect on
the known matter. As such, it is conceivable that this evidence could be
explained by a modification to gravitation and/or concepts of inertia. Various
formulations of modified gravity exist, each giving rise to several non-
canonical outcomes. This motivates us to propose an experiment searching for
departures from (quantum) Newtonian predictions in a bipartite setting with
gravitational accelerations $\lesssim 10^{-10}$ m/s2, i.e., where the
effective force needs to be stronger than Newtonian to account for the Dark
Matter effects. Since quantum particles naturally source weak gravitation,
their non-relativistic dynamics offers opportunities to test this small
acceleration regime. We show that two nearby mesoscopic quantum masses
accumulate significantly larger entanglement in modified gravity models, such
as the Modified Newtonian Dynamics. Our calculations include Casimir-Polder
forces as well as tidal effects next to the surface of the earth, and confirm
that entanglement is observable within the limits imposed by environmental
decoherence. We demonstrate how the temperature can be fine-tuned such that
modified gravity is certified simply by witnessing the entanglement generated
from uncorrelated thermal states, eliminating the need for precise noise
characterization. Overall, the required parameters could be realized in a
tabletop experiment.

_Introduction._ The Newtonian limit of General Relativity is very successful
on the scale of the solar system. For example, balancing the centrifugal and
the gravitational forces of objects in approximately circular orbits around
the Sun implies that orbital velocity falls as the square root of the
distance. This is famously known as Keplerian decline and has been observed to
hold for all planets [1]. Spiral galaxies have a lot in common with the solar
system. Most of their mass is also concentrated towards the center, but the
stars do not show any asymptotic Keplerian decline [2]. Their orbital speeds
generally do not fall, and the rotation curves saturate [3]. Consequently, the
stars in the outer regions appear to be orbiting so fast that they should not
be gravitationally bound. This is not happening, and hence there seems to be
more gravity than expected based on the known visible mass at the center of
spiral galaxies. This is a prime example of the Dark Matter (DM) effect. The
name originates in the proposal of the existence of an invisible matter
distributed throughout the galaxies [4], generating an extra gravitational
pull that balances the centrifugal force. Despite being the most widely
accepted explanation and with evidence appearing even on the largest
cosmological scales, existence of DM has not been directly detected or
confirmed by any experiment so far [5], and hence the continued interest in
alternative solutions. A plausible route involves modifications to our present
understanding of gravity. Modified Newtonian Dynamics (MOND) is one such
proposal where, without invoking DM, Newton’s second law and/or the law of
universal gravitation is modified to account for DM effects in galaxies [6].
While the experiment we propose is independent of any concrete formulation of
alternative gravity models, for quantitative statements we follow the
parameters present in the MOND theory. We therefore begin with more details
about it.

_Modified Newtonian Dynamics._ A general form of MOND, which encompasses a
wide variety of proposed variants of the model [7, 8, 9, 10, 11, 12] ,
generalizes Newton’s second law to

$\nu\bigg{(}\frac{|\vec{a}|}{a_{0}}\bigg{)}\vec{a}=-\nabla\Phi+\vec{F}/m,$ (1)

where $\vec{F}$ represents the sum of all non-gravitational forces on an
object of mass $m$, $\vec{a}$ is its acceleration, and $a_{0}$ is the
acceleration scale where the generalization is to take place. The potential
$\Phi$ belongs to the set of potentials $\\{\phi_{i}\\}$, or it could be their
linear combination. Collectively, the potentials obey a system of potentially
nonlinear and coupled Poisson equations:

$\sum_{j}\nabla\cdot\bigg{[}\mu_{ij}\bigg{(}\frac{{\nabla}\phi_{1}}{a_{0}},\frac{{\nabla}\phi_{2}}{a_{0}},\dots\bigg{)}\nabla\phi_{j}\bigg{]}=4\pi
G_{i}\rho,$ (2)

where $\rho$ is the density of matter, $G_{i}$ is the generalized Newton’s
constant quantifying the coupling of matter with the potential $\phi_{i}$, and
$(\nu,\mu_{ij})$ are the interpolating functions that govern the transition
from Newtonian dynamics to the modified regime. Many existing models are
special cases of the above expressions (see Supplemental Material for various
existing formulations). Proposals where $d\nu(x)/dx\neq 0$ are referred to as
models of _modified inertia_. It is not known whether a viable action
principle for these models exists, and so it is not clear whether the degrees
of freedom describing an object can be canonically quantized. Since we explore
the possibility of entanglement accumulation as the signal revealing non-
canonical gravity, we effectively put to test only the _modified gravity_
formulations of MOND ($\nu(x)=1$), as the existence of a Hamiltonian allows us
to employ standard quantum mechanical methods. The exact functional forms for
$\mu_{ij}$ are still debated, but consistency with astronomical observations
fixes their asymptotics. For example, in the Lagrangian-based AQUAL model
$\nu=1$ and there exists a single potential $\Phi$ with interpolating function
$\tilde{\mu}$ [8]. In the limiting cases $\tilde{\mu}(x\ll 1)\to x$ to explain
galactic rotation curves, and $\tilde{\mu}(x\gg 1)\to 1$ to recover Newtonian
dynamics in regimes of stronger gravity. Note that even a single massive
object may source multiple potentials $\phi_{i}$. Notably, although each model
is non-relativistic, some exist as limits of relativistic completions. A
recent example of such a completion, building on antecedent work [9], was
shown to be consistent with data from the large-scale cosmic structure and
inhomogeneities of the cosmic microwave background [13]. Although there does
not exist a unique model for the MOND paradigm, we will see that many of them
and potentially other modified gravity theories predict deviations from
Newtonian gravity that may be observable in bipartite quantum mechanical
experiments.

A wide array of DM phenomenology in galaxies is associated with acceleration
due to gravity below $a_{0}\approx 1.2\times 10^{-10}$ m/s2 [6, 14, 15]. Such
small accelerations are readily accessible between quantum masses. We show
here how the entanglement dynamics between two microspheres senses the force
gradient between them, and design tests where a simple act of entanglement
witnessing reveals non-canonical interaction. The calculations include
Casimir-Polder forces, tidal effects next to the surface of the earth, and
estimates for the effects of environmental decoherence. The modified gravity
models have a stronger force gradient than Newtonian in the discussed limit,
thereby leading to a stronger entanglement that is robust against noisy
measurements. Such an experiment could be done using the masses recently
cooled in Vienna [16] which, when separated by a distance of a few times their
radius, generate an acceleration deep into the discussed relevant regime. Our
calculations at first assume the particles are in free space devoid of all
external gravitational fields. Although this might seem unrealistic, there may
be places in the solar system where the net gravitational acceleration is in
the MOND regime [17, 18]. Furthermore, it is also instructive to perform this
experiment either with freely falling microspheres in a large orbit around the
earth, where pressure falls to $\sim 10^{-15}$ Pa and the decoherence effects
are negligible, or in a laboratory on earth. In the latter case, the
acceleration along the line between the masses, sourced by the tidal force, is
two orders of magnitude larger than the acceleration due to mutual gravity.
Yet, we show that it has a negligible effect on generated entanglement and
hence opens the way towards a terrestrial experiment. An observation of non-
canonical effects in agreement with our calculations would indicate that the
modified gravity adheres to the Strong Equivalence Principle (SEP), i.e., the
internal dynamics stays the same in all uniformly accelerated frames. See also
Ref. [19], which discusses probing MOND with the Satellite Test of the
Equivalence Principle. Note that SEP is certainly violated in the MOND model.
Many other proposals have been put forward to test the predictions of MOND.
Torsion pendula experiments show an agreement with Newton’s second law down to
accelerations $\approx a_{0}/2400$ [20], and with Newton’s law of universal
gravity down to $\approx a_{0}/60$ [21]. The accelerations mentioned in these
experiments are with respect to the local laboratory, and the net acceleration
due to heavenly sources is well above $a_{0}$. Therefore, it remains unclear
to what extent departures from Newtonian dynamics are to be expected.

The driving idea behind modified gravity models is the following: just like
Newtonian gravity is an approximation of General Relativity when the
gravitational field is not too strong, it might also be an approximation of an
underlying theory when the Newtonian accelerations are not too small. In order
to match the stellar rotation curves, the modification has to involve
accelerations smaller than $a_{0}\approx 1.2\times 10^{-10}$ m/s2, and the
simplest solution is to demand that gravitational force in this regime scales
inversely to the distance. Accordingly, in the deep MOND regime the modified
gravitational potential at a position $\vec{r}_{2}$ due to a mass $m$ located
at position $\vec{r}_{1}$ is given by $\Phi=\sqrt{Gma_{0}}\
\ln(\absolutevalue{\vec{r}_{2}-\vec{r}_{1}})$. Note that this violates the law
of equal and opposite action and reaction: the force on a mass $m_{2}$ due to
$m_{1}$ is $\sim m_{2}\sqrt{m_{1}}$, whereas the force on mass $m_{1}$ due to
$m_{2}$ is $\sim m_{1}\sqrt{m_{2}}$. This issue is rectified if one solves the
full nonlinear Poisson-like equation governing the bipartite dynamics in
MONDian gravity. The resultant gravitational potential energy of two identical
particles of mass $m$ is given by [22, 23, 24]:

$V_{M}(\vec{r}_{1},\vec{r}_{2})=\frac{4}{3}\quantity(\sqrt{2}-1)\
m\sqrt{Gma_{0}}\ \ln(L+r),$ (3)

where $r$ is the _relative displacement_ of the two masses from their initial
separation $L$. This is very different from the usual Newtonian potential and
underlies the differences in observable quantities. In particular, we focus on
quantum entanglement between small objects.

_Entanglement Dynamics._

Figure 1: Gravitational and Casimir entanglement between two identical
Platinum spheres of radius 10 $\mu$m. The two masses are initially separated
by a distance of 100 $\mu$m, and their initial quantum states are prepared by
cooling in harmonic traps of frequency $10$ MHz. Entanglement is measured by
logarithmic negativity (note also logarithmic vertical scale), $\bar{n}$ is
the average phonon number (temperature), and the grey shaded region
(logarithmic negativity $<10^{-4}$) signifies entanglement sensitivity beyond
the limits achievable with near-future technology [25, 26]. The dashed lines
present calculations with only Newtonian potential, only gravity due to
Modified Newtonian Dynamics (MOND) and taking into account only Casimir
attraction. The solid lines describe the simultaneous presence of (modified)
gravity and Casimir. The regions marked as “Experiment” correspond to
entanglement values indicating modified gravity. Detecting any entanglement
before the time when the red curve comes out of the grey region provides the
simplest check of non-canonical interaction.

Entanglement dynamics of two nearby quantum masses in empty space, starting
with uncorrelated thermal states, has been studied in detail in Refs. [27, 28,
29]. The system gets entangled because of the position dependence of the
gravitational force: the parts of the wave packet closer to each other are
attracted more than the parts further away. Accordingly, different momenta are
generated at different positions, producing entanglement as time passes.
Notably, the methods of Ref. [29] apply to any central interactions without
the need for an explicit form of the Hamiltonian. The displacement-to-
separation ratio $r/L$ is always small in gravitational experiments, and hence
we expand Eq. (3) in a Maclaurin series. Note that the cubic term is relevant
in entanglement dynamics only when there is a significant relative motion
between the two particles [29], and hence we truncate at the quadratic term:

$V_{M}(r)\approx\frac{4}{3}\quantity(\sqrt{2}-1)\
m\sqrt{Gma_{0}}\quantity(\ln(L)+\frac{r}{L}-\frac{r^{2}}{2L^{2}}).$ (4)

For quantum particles this interaction is quantized by replacing variable $r$
with operator $\hat{r}$. The system admits Gaussianity at all times, and we
quantify entanglement through the logarithmic negativity of the bipartite
covariance matrix [30, 31, 32, 33]. The covariance matrix starting from the
uncorrelated ground states of individual harmonic traps, each with frequency
$\omega_{0}$, is derived in exact closed form for any quadratic and central
potential in Ref. [29], and it reads:

$\displaystyle\bm{\sigma}_{00}$ $\displaystyle=$
$\displaystyle\frac{\hbar}{4m\omega_{0}}\quantity[2+\omega_{0}^{2}t^{2}+\quantity(1+\frac{\omega_{0}^{2}}{\omega^{2}})\sinh^{2}(\omega
t)],$ $\displaystyle\bm{\sigma}_{02}$ $\displaystyle=$
$\displaystyle\frac{\hbar}{4m\omega_{0}}\quantity[\omega_{0}^{2}t^{2}-\quantity(1+\frac{\omega_{0}^{2}}{\omega^{2}})\sinh^{2}(\omega
t)],$ $\displaystyle\bm{\sigma}_{11}$ $\displaystyle=$
$\displaystyle\frac{m\hbar\omega_{0}}{4}\quantity[2+\quantity(1+\frac{\omega^{2}}{\omega_{0}^{2}})\sinh^{2}(\omega
t)],$ $\displaystyle\bm{\sigma}_{13}$ $\displaystyle=$
$\displaystyle-\frac{m\hbar\omega_{0}}{4}\quantity(1+\frac{\omega^{2}}{\omega_{0}^{2}})\sinh^{2}(\omega
t),$ $\displaystyle\bm{\sigma}_{01}$ $\displaystyle=$
$\displaystyle\frac{\hbar}{8}\quantity[2\omega_{0}t+\quantity(\frac{\omega_{0}}{\omega}+\frac{\omega}{\omega_{0}})\sinh(2\omega
t)],$ $\displaystyle\bm{\sigma}_{03}$ $\displaystyle=$
$\displaystyle\frac{\hbar}{8}\quantity[2\omega_{0}t-\quantity(\frac{\omega_{0}}{\omega}+\frac{\omega}{\omega_{0}})\sinh(2\omega
t)],$ (5)

where $\omega$ encodes the mutual interaction between the two particles as
follows. For an arbitrary central potential expanded in a series in terms of
the displacement-to-separation ratio, the parameter $\omega$ is found by
equating the coefficient of the quadratic term with $-m\omega^{2}/4$ [29]. For
more than one central interaction present simultaneously, the “total” $\omega$
is given by a Pythagoras-like sum for the individual interactions:
$\omega^{2}=\sum_{i}\omega_{i}^{2}$ [29]. Note that the covariance matrix is
entirely determined by frequencies $\omega_{0}$ and $\omega$. For the standard
Newtonian gravity one finds $\omega_{N}^{2}=4Gm/L^{3}$, and for the MONDian
gravity described by Eq. 3 we obtain

$\omega^{2}_{M}=\frac{8}{3}(\sqrt{2}-1)\frac{\sqrt{Gma_{0}}}{L^{2}}.$ (6)

In an actual experiment, the gravitational force is accompanied by non-
gravitational interactions, which also admit an entangling effect. The most
important of these is the Casimir attraction between the two masses. Following
the algorithm described above, we find that contribution to entanglement is
parameterized by [34, 35] (see Supplemental Material for complete derivation)

$\omega_{C}^{2}=\frac{2\hbar cR_{0}^{6}}{m\pi
L^{9}}\sum_{n=0}^{9}C_{n}(n+7)(n+8)\frac{R_{0}^{n}}{L^{n}},$ (7)

where the constants $\quantity{C_{n}}$ were determined in Ref. [34] and
$R_{0}$ stands for the radius of each metallic sphere. Placing the two objects
in close proximity increases gravitational coupling, but such a setup would
also encounter an overwhelmingly large Casimir attraction, practically
suppressing all gravitational effects. On the other hand, increasing the
separation to an extent that the Casimir forces are negligible also diminishes
the gravitational coupling drastically, leading to a large amount of time
required to accumulate a detectable entanglement. We therefore propose a
configuration where the two masses are separated by a moderate distance such
that the force gradients of gravitational and Casimir interactions are
comparable.

In Fig. 1 we show the entanglement gain between two Platinum
spheres111Although Osmium is the densest naturally occurring material, it is
hard and brittle, making it difficult to shape as a sphere. Platinum is almost
as dense and malleable, and it can be forged into a spherical shape of a
predetermined size. of radius $R_{0}=10$ $\mu$m separated by $L=10R_{0}$. The
configuration admits an internal _relative acceleration_ of $\approx
a_{0}/100$, readily in the regime of our interest and beyond the limits tested
with torsion pendulums experiments [21]. The entanglement accumulation due to
the MOND, Newtonian gravity, and the Casimir interaction are independently
characterized by $(\omega_{M},\omega_{N},\omega_{C})=(3.06,1.55,3.53)\times
10^{-4}$. The density of Platinum is $21.45$ g/cm3, which implies the
considered objects have a mass $m\approx 90$ ng. We propose to prepare their
initial (Gaussian) states by cooling the masses in harmonic traps of frequency
$10$ MHz, which results in an initial ground state localization of
$\sqrt{\hbar/2m\omega_{0}}\approx 0.1$ fm, achievable with currently available
optical tweezers. The left panel of Fig. 1 shows the ideal case of
$\bar{n}=0$, where the masses are cooled to the ground state. The middle panel
shows a realistic scenario where the initial state is a mixed Gaussian with a
no- zero phonon number $\bar{n}=0.0015$. The covariance matrix for initially
uncorrelated thermal states, characterized by the average phonon number
$\bar{n}$, is given by $\quantity(2\bar{n}+1)$ times the one obtained for the
initial ground states. The solid lines represent the entanglement generated by
the simultaneous action of gravitational and Casimir forces, while dashed
lines serve as a reference, illustrating the partial contributions of
individual interactions. The MONDian entanglement is much higher than what is
accumulated by Newtonian gravity. Note that the shaded region within the two
curves corresponds to other possible approaches to modified gravity predicting
a force gradient between Newtonian and MONDian. Hence, the experimenter is not
required to quantify entanglement, and a mere entanglement witnessing within
the time window $110-125$ milliseconds will certify departures from canonical
gravity. The right panel of Fig. 1 shows this experimental time window for
initial states with different phonon numbers.

_Discussion._ In order to experimentally estimate the amount of entanglement
one needs to measure elements of the covariance matrix, which in principle
requires finding position and momentum correlations [36]. The question emerges
as to whether a simpler method could exist to detect the stronger-than-
Newtonian force, e.g., by measurements of position only. This is of course
possible, but it turns out that the average displacements of the masses are
about three orders of magnitude smaller than their standard deviations, at
least for the parameters considered here. While the alternative models give
rise to larger than Newtonian variances, from a practical perspective one
needs to ensure that the predicted increment is not mimicked by the noise in
the system. The added value of entanglement witnessing is its independence
from noise characterization, whereas the added complexities are not high since
the measurements of different mechanical quadratures can be accomplished,
e.g., by choice of local phase in homodyne detection [36].

_Environmental Decoherence:_ We now estimate the effects of decoherence. The
interaction of the considered masses with ambient thermal photons and air
molecules suppresses coherences over long distances. For the gravitational
entanglement proposals, such decoherence mechanisms have been studied in Refs.
[37, 38, 27, 39, 40]. In the regime where the superposition is much smaller
than the wavelength of the scattering particles, the coherence time for a
quantum superposition across a distance $\Delta x$ is given by
$\tau=\frac{1}{\Lambda(\Delta x^{2})}$, where $\Lambda$ depends on the
characteristics of environment. The interactions with ambient thermal photons
is quantified by $\Lambda_{\text{ph}}\approx
10^{36}R_{0}^{6}T_{\text{env}}^{9}$ [41]. In the discussed setup, we take the
final width of the wave packet at time 125 ms (upper limit of time window at
$\bar{n}=0.0015$), i.e., $\Delta x\approx 0.76$ nm. An ambient temperature of
2 mK (achievable with 3He/4He dilution refrigerator [42]) implies that the
coherence time is $\tau_{\text{ph}}\sim 10^{36}$ s. This is not surprising as
the intensity of the thermal background peaks at $\approx 1.45$ m, which is
nine orders of magnitude larger than the position spread. Interactions with
ambient air molecules is still an issue. At 2 mK ambient temperature the
thermal de Broglie wavelength of air molecules is an order of magnitude larger
than $\Delta x$, and hence the decoherence is quantified by [41]

$\Lambda_{\text{air}}=\frac{8n_{\text{air}}}{3\hbar^{2}}\sqrt{2\pi
m_{\text{air}}}(k_{B}T_{\text{env}})^{3/2}R_{0}^{2},$ (8)

where $m_{\text{air}}\approx 0.5\times 10^{-25}$ kg is the average mass of an
air molecule, and $n_{\text{air}}$ is the air density in the vacuum chamber.
In an ultra-high vacuum with a particle density $\approx 3\times 10^{4}$ cm-3
[43], we get a coherence time $\approx$ 0.94 s, well above the requirements of
this proposal. Note the short experimental time requires a free fall inside a
vacuum chamber of height $\lesssim 10$ cm, adding to the ease of creating and
maintaining low pressures and temperatures. All this suggests a plausible,
though demanding, terrestrial experiment. Accordingly, in the following
section we further address the effects of Earth’s proximity.

_Tabletop Experiment:_

Figure 2: Origin of tidal forces in a bipartite system of two identical masses
$m$ separated by a distance $L$, next to the surface of the earth (extremely
not to scale). Here $g=|\vec{g}|=9.81$ m/s2 is the terrestrial acceleration,
and $R_{E}=6371$ km denotes the radius of the earth. The distances are over-
exaggerated as in fact $L\sim 10^{-11}R_{E}$. Yet, the tidal components along
the line joining the masses are two orders of magnitude larger than mutual
gravity.

Next to the surface of the Earth each particle experiences a strong pull
towards the center of the planet, which contributes non-zero components along
the line connecting the centers of the masses. These tidal forces arising due
to the non-uniform gravity of the earth are illustrated in Fig. 2. Since that
the acceleration due to gravity is $\sim 10^{11}a_{0}$, we treat the tidal
forces as horizontal components of the Newtonian force $mg$, which add
linearly to the mutual gravity. As seen in the figure, the magnitude of the
tidal forces is $mg\sin\theta$, and they drift the masses towards each other.
Replacing $\sin\theta$ using geometry in Fig. 2 we find:

$\displaystyle\vec{F}^{\text{(Tidal)}}_{A}$ $\displaystyle=$
$\displaystyle-m\frac{g}{R_{E}}\quantity(-\frac{L}{2}+x_{A})\ \hat{e}_{x},$
$\displaystyle\vec{F}^{\text{(Tidal)}}_{B}$ $\displaystyle=$
$\displaystyle-m\frac{g}{R_{E}}\quantity(\frac{L}{2}+x_{B})\ \hat{e}_{x},$ (9)

where $\hat{e}_{x}$ is the unit vector along $x$ axis (origin at the center
between the masses), and $x_{A}(x_{B})$ is the displacement of mass $A(B)$
from its initial position. The tidal potential is therefore

$\displaystyle V^{\text{(Tidal)}}$ $\displaystyle=$
$\displaystyle\frac{1}{2}m\frac{g}{R_{E}}\left[\quantity(-\frac{L}{2}+x_{A})^{2}+\quantity(\frac{L}{2}+x_{B})^{2}\right],$
(10)

which appears to confine each mass in a harmonic trap of frequency
$\omega_{E}=\sqrt{g/R_{E}}\approx 0.00124$. Even though $\omega_{E}$ is an
order of magnitude larger than the characteristic frequencies of mutual
gravity (and the tidal acceleration about hundred times that of mutual
attraction), it is a _local_ interaction and as such it cannot generate any
entanglement. Yet, it does influences the wave function locally, which in turn
can speed up / slow down the entanglement accumulation due to mutual
interactions. To check if this effect is relevant we solve the time evolution
of covariance matrix with mutual interactions truncated at the quadratic term
and with the tidal potential. The solution reads:

$\displaystyle\bm{\sigma}_{00}(\bm{\sigma}_{02})$ $\displaystyle=$
$\displaystyle\frac{\hbar}{4m\omega_{0}}\quantity{\quantity(\cos^{2}(\omega_{E}t)+\frac{\omega_{0}^{2}}{\omega_{E}^{2}}\sin^{2}(\omega_{E}t))+\\!(\\!-\\!)\
\quantity[\cos^{2}\quantity(\sqrt{\omega_{E}^{2}-\omega^{2}}\,t)+\frac{\omega_{0}^{2}}{\omega_{E}^{2}-\omega^{2}}\sin^{2}\quantity(\sqrt{\omega_{E}^{2}-\omega^{2}}\,t)]},$
$\displaystyle\bm{\sigma}_{11}(\bm{\sigma}_{13})$ $\displaystyle=$
$\displaystyle\frac{1}{4}m\hbar\omega_{0}\quantity{\quantity(\cos^{2}(\omega_{E}t)+\frac{\omega_{E}^{2}}{\omega_{0}^{2}}\sin^{2}(\omega_{E}t))+\\!(\\!-\\!)\
\quantity[\cos^{2}\quantity(\sqrt{\omega_{E}^{2}-\omega^{2}}\,t)+\frac{\omega_{E}^{2}-\omega^{2}}{\omega_{0}^{2}}\sin^{2}\quantity(\sqrt{\omega_{E}^{2}-\omega^{2}}\,t)]},$
$\displaystyle\bm{\sigma}_{01}(\bm{\sigma}_{03})$ $\displaystyle=$
$\displaystyle\frac{\hbar}{8}\Bigg{[}\Big{(}\frac{\omega_{0}}{\omega_{E}}-\frac{\omega_{E}}{\omega_{0}}\Big{)}\sin(2\omega_{E}t)+\\!(\\!-\\!)\
\quantity(\frac{\omega_{0}}{\sqrt{\omega_{E}^{2}-\omega^{2}}}-\frac{\sqrt{\omega_{E}^{2}-\omega^{2}}}{\omega_{0}})\sin\quantity(2\sqrt{\omega_{E}^{2}-\omega^{2}}\,t)\Bigg{]},$
(11)

where again $\omega_{0}$ is the frequency of harmonic trap used for preparing
initial Gaussian states, and $\omega$ encodes all mutual interactions between
the two masses. While this solution looks distinct from that in Eq. (5), the
entanglement negativity in Fig. 1 decreases relatively by a negligible amount
$\sim 10^{-9}$, many orders of magnitude smaller than the achievable
laboratory precision. In short, for the parameters considered in this work,
tidal forces are insignificant in entanglement dynamics.

_Summary._ We proposed an entanglement experiment with the potential outcome
revealing non-Newtonian gravity at small accelerations. It was motivated by
modified gravity models, and the concrete calculations were performed using
MOND parameters. Quantum particles naturally source weak gravity and small
accelerations and provide an alternative route for testing modified gravity.
The entanglement measurement brings independence from detailed noise
characterization while adding only moderate complexity. Calculations including
Casimir forces, tidal effects, and environmental decoherence reveal the
possibility of a tabletop experiment.

_Acknowledgments._ This work is jointly supported by (i) Xiamen University
Malaysia, via projects XMUMRF/2022-C10/IPHY/0002 and XMUMRF/2021-C8/IPHY/0001,
(ii) SERB-DST, Govt. of India, via project no. CRG/2022/009359, and (iii)
project 2021/43/P/ST2/02141 co-funded by the NCN and the EU’s Horizon 2020
under the Marie Skłodowska-Curie grant agreement 945339. A.K. is supported by
the (i) Foundation for Polish Science through the IRAP project “International
Centre for Theory of Quantum Technologies” (contract no. 2018/MAB/5) co-
financed by EU’s Smart Growth Operational Program, and (ii) QuantERA II
project “Experiment and Theory of Resources in Quantum Technologies” (contract
no. 2021/03/Y/ST2/00178) that has received funding from EU’s Horizon 2020. We
gratefully thank the anonymous reviewer for their thorough reading and
attention to to detail; their insightful comments resulted in a significant
improvement in the clarity, correctness, and completeness of this work. A.K.
thanks Rishabh Duhan (TIFR, India) for explaining the current technological
limits for creating and maintaining low pressures and temperatures. We
acknowledge National Supercomputing Mission (NSM) for providing computing
resources of ‘PARAM Ganga’ at Indian Institute of Technology Roorkee, which is
implemented by C-DAC and supported by the Ministry of Electronics and
Information Technology (MeitY) and Department of Science and Technology (DST),
Government of India.

## References

  * [1] D. R. Williams, NASA Planetary Fact Sheet, https://nssdc.gsfc.nasa.gov/planetary/factsheet/.
  * Sanders [2010] R. H. Sanders,  _The Dark Matter Problem: A Historical Perspective_ (Cambridge University Press, 2010).
  * Babcock [1939] H. W. Babcock, Lick Observatory Bulletin 498, 41 (1939).
  * Zwicky [1937] F. Zwicky, The Astrophysical Journal 86, 217 (1937).
  * Bertone and Hooper [2018] G. Bertone and D. Hooper, Reviews of Modern Physics 90, 045002 (2018).
  * Milgrom [1983] M. Milgrom, The Astrophysical Journal 270, 365 (1983).
  * Milgrom [1994] M. Milgrom, Annals of Physics 229, 384 (1994).
  * Bekenstein and Milgrom [1984] J. D. Bekenstein and M. Milgrom, The Astrophysical Journal 286, 7 (1984).
  * Bekenstein [2004] J. D. Bekenstein, Physical Review D 70, 083509 (2004).
  * Milgrom [2010] M. Milgrom, Monthly Notices of the Royal Astronomical Society 403, 886 (2010).
  * Milgrom [2023a] M. Milgrom, Physical Review D 108, 063009 (2023a).
  * Milgrom [2023b] M. Milgrom, Physical Review D 108, 084005 (2023b).
  * Skordis and Złośnik [2021] C. Skordis and T. Złośnik, Physical Review Letters 127, 161302 (2021).
  * Famaey and Binney [2005] B. Famaey and J. Binney, Monthly Notices of the Royal Astronomical Society 363, 603 (2005).
  * Gentile  _et al._ [2011] G. Gentile, B. Famaey, and W. J. G. de Blok, Astronomy & Astrophysics 527, A76 (2011).
  * Delić  _et al._ [2020] U. Delić, M. Reisenbauer, K. Dare, D. Grass, V. Vuletić, N. Kiesel, and M. Aspelmeyer, Science 367, 892 (2020).
  * Bekenstein and Magueijo [2006] J. Bekenstein and J. Magueijo, Physical Review D 73, 103513 (2006).
  * Magueijo and Mozaffari [2012] J. Magueijo and A. Mozaffari, Physical Review D 85, 043527 (2012).
  * Pereira  _et al._ [2016] J. P. Pereira, J. M. Overduin, and A. J. Poyneer, Physical Review Letters 117, 071103 (2016).
  * Gundlach  _et al._ [2007] J. H. Gundlach, S. Schlamminger, C. D. Spitzer, K.-Y. Choi, B. A. Woodahl, J. J. Coy, and E. Fischbach, Physical Review Letters 98, 150801 (2007).
  * Little and Little [2014] S. Little and M. Little, Classical and Quantum Gravity 31, 195008 (2014).
  * Milgrom [2014] M. Milgrom, Physical Review D 89, 024016 (2014).
  * Zhao  _et al._ [2010] H. Zhao, B. Li, and O. Bienaymé, Physical Review D 82, 103001 (2010).
  * Milgrom [1986] M. Milgrom, The Astrophysical Journal 302, 617 (1986).
  * Palomaki  _et al._ [2013] T. A. Palomaki, J. D. Teufel, R. W. Simmonds, and K. W. Lehnert, Science 342, 710 (2013).
  * Krisnanda  _et al._ [2023] T. Krisnanda, T. Paterek, M. Paternostro, and T. C. H. Liew, Physical Review D 107, 086014 (2023).
  * Krisnanda  _et al._ [2020] T. Krisnanda, G. Y. Tham, M. Paternostro, and T. Paterek, npj Quantum Information 6, 12 (2020).
  * Qvarfort  _et al._ [2020] S. Qvarfort, S. Bose, and A. Serafini, Journal of Physics B: Atomic, Molecular and Optical Physics 53, 235501 (2020).
  * Kumar  _et al._ [2023] A. Kumar, T. Krisnanda, P. Arumugam, and T. Paterek, Quantum 7, 1008 (2023).
  * Simon [2000] R. Simon, Physical Review Letters 84, 2726 (2000).
  * Adesso and Illuminati [2005] G. Adesso and F. Illuminati, Physical Review A 72, 032334 (2005).
  * Vidal and Werner [2002] G. Vidal and R. F. Werner, Physical Review A 65, 032314 (2002).
  * Adesso  _et al._ [2004] G. Adesso, A. Serafini, and F. Illuminati, Physical Review A 70, 022318 (2004).
  * Emig  _et al._ [2007] T. Emig, N. Graham, R. L. Jaffe, and M. Kardar, Physical Review Letters 99, 170403 (2007).
  * Casimir and Polder [1948] H. B. G. Casimir and D. Polder, Physical Review 73, 360 (1948).
  * Aspelmeyer  _et al._ [2014] M. Aspelmeyer, T. J. Kippenberg, and F. Marquardt, Reviews of Modern Physics 86, 1391 (2014).
  * Bose  _et al._ [2017] S. Bose, A. Mazumdar, G. W. Morley, H. Ulbricht, M. Toroš, M. Paternostro, A. A. Geraci, P. F. Barker, M. S. Kim, and G. Milburn, Physical Review Letters 119, 240401 (2017).
  * Marletto and Vedral [2017] C. Marletto and V. Vedral, Physical Review Letters 119, 240402 (2017).
  * Rijavec  _et al._ [2021] S. Rijavec, M. Carlesso, A. Bassi, V. Vedral, and C. Marletto, New Journal of Physics 23, 043040 (2021).
  * Datta and Miao [2021] A. Datta and H. Miao, Quantum Science and Technology 6, 045014 (2021).
  * Schlosshauer [2007] M. A. Schlosshauer,  _Decoherence and the Quantum-To-Classical Transition_, 1st ed. (Springer-Verlag, Germany, 2007).
  * Zu  _et al._ [2022] H. Zu, W. Dai, and A. de Waele, Cryogenics 121, 103390 (2022).
  * [43] VAC AERO International Inc., The Fundamentals of Vacuum Theory, https://vacaero.com/information-resources/vac-aero-training/170466-the-fundamentals-of-vacuum-theory.html.

SUPPLENMENTAL MATERIAL

## Different MOND formulations

The generalizations of Newton’s second law of motion and Poisson’s equation
describe a modification to the acceleration $\vec{a}$ of an object of mass $m$
under the influence of a non-gravitational force $\vec{F}$ and gravitational
force $-m{\nabla}\Phi$, where $\Phi$ is determined via the solution of a set
of potentially modified Poisson equations for the set of potentials
$\quantity{\phi_{i}}$ [see main text for details]. Examples of several
specific proposals are given in Table 1, where non-zero $\mu_{ij}$ and $G_{i}$
are detailed. We note that there can additionally exist formulations of MOND
where the quantity $a_{0}$ itself may depend on $\Phi$ or $|\nabla\Phi|$ [23],
where $\Phi$ appears as a ‘mass term’ in the field equations [13], or where
terms involving higher derivatives of potentials appear in the field equations
[12].

Table 1: Different formulations of the MOND theory. Model | Potentials | Properties
---|---|---
Modified Inertia [7] | $\phi_{1}=\Phi$ | $\nu=\nu(|\vec{a}|/a_{0})$
|  | $\mu_{11}=1$
AQUAL [8] | $\phi_{1}=\Phi$ | $\nu=1$
|  | $\mu_{11}=\mu_{11}(|\nabla\Phi|/a_{0})$
|  | $G_{1}=G$
TeVeS [9] | $(\phi_{1},\phi_{2})=$ | $\nu=1$
| $(\Phi-\varphi,\varphi)$ | $\mu_{11}=1$
|  | $\mu_{22}=\mu_{22}(|\nabla\varphi|/a_{0})$
|  | $G_{1}=G_{2}=G$
QUMOND [10] | $(\phi_{1},\phi_{2})=$ | $\nu=1$
| $(\Phi,\varphi)$ | $\mu_{11}=1$
|  | $\mu_{12}=\mu_{12}(|\nabla\varphi|/a_{0})$
|  | $\mu_{22}=1$
|  | $G_{2}=G$
TRIMOND [11] | $(\phi_{1},\phi_{2},\phi_{3})=$ | $\nu=1$
| $(\Phi,\varphi,\psi)$ | $\mu_{13}=1$
|  | $\mu_{22}=\mu_{22}(\nabla\varphi,\nabla\psi)$
|  | $\mu_{23}=\mu_{23}(\nabla\varphi,\nabla\psi)$
|  | $\mu_{31}=1$
|  | $\mu_{32}=\mu_{32}(\nabla\varphi,\nabla\psi)$
|  | $\mu_{33}=\mu_{33}(\nabla\varphi,\nabla\psi)$
|  | $G_{1}=G$

## Bipartite covariance matrix

The covariance matrix formalism is based on the first two statistical moments
of a quantum state. Given a bipartite system $AB$ with the displacement and
momentum quadratures
$\hat{u}=(\hat{x}_{A},\hat{p}_{A},\hat{x}_{B},\hat{p}_{B})^{T}$, the
(symmetric) covariance matrix is defined by [32, 33, 31]:

$\bm{\sigma}_{jk}=\bm{\sigma}_{kj}=\frac{1}{2}\expectationvalue{\anticommutator{\hat{u}_{j}}{\hat{u}_{k}}}-\expectationvalue{\hat{u}_{j}}\expectationvalue{\hat{u}_{k}}.$
(12)

In the block form it reads

$\bm{\sigma}\equiv\matrixquantity(\bm{\alpha}&\bm{\gamma}\\\
\bm{\gamma}^{T}&\bm{\beta}),$ (13)

where $\bm{\alpha}(\bm{\beta})$ contains the local mode correlation for
$A(B)$, and $\bm{\gamma}$ describes the intermodal correlation. In this work
we deal with two identical masses prepared in identical initial states, which
implies that the local modes are same at all times: $\bm{\alpha}=\bm{\beta}$.

## The Casimir interaction

The Casimir interaction is central, and in the relevant limit of perfect
metals the energy of attraction between two identical spheres is given by [34,
35]:

$V_{C}(r)=-\frac{\hbar
cR_{0}^{6}}{\pi(L+r)^{7}}\sum_{n=0}^{\infty}C_{n}\quantity(\frac{R_{0}}{L+r})^{n},$
(14)

where $R_{0}$ is the radius of the sphere and the constants for the first ten
terms were determined in Ref. [34]:

$\displaystyle C_{0}$ $\displaystyle=$ $\displaystyle+143/16,$ $\displaystyle
C_{1}$ $\displaystyle=$ $\displaystyle 0,$ $\displaystyle C_{2}$
$\displaystyle=$ $\displaystyle+7947/160,$ $\displaystyle C_{3}$
$\displaystyle=$ $\displaystyle+2065/32,$ $\displaystyle C_{4}$
$\displaystyle=$ $\displaystyle+27\ 705\ 347/100\ 800,$ $\displaystyle C_{5}$
$\displaystyle=$ $\displaystyle-55\ 251/64,$ $\displaystyle C_{6}$
$\displaystyle=$ $\displaystyle+1373\ 212\ 550\ 401/144\ 506880,$
$\displaystyle C_{7}$ $\displaystyle=$ $\displaystyle-7\ 583\ 389/320,$
$\displaystyle C_{8}$ $\displaystyle=$ $\displaystyle-2\ 516\ 749\ 144\ 274\
023/44\ 508\ 119\ 040,$ $\displaystyle C_{9}$ $\displaystyle=$
$\displaystyle+274\ 953\ 589\ 659\ 739/275\ 251\ 200.$ (15)

While the series, in principle, does not converge for any value of $R_{0}/L$,
the results presented in Fig. 1 of Ref. [34] show that the leading term alone
is sufficient to describe the configuration considered in this work.
Nevertheless, for possible implementation in other configurations, we expand
all available terms in binomial series of $r/L$:

$\displaystyle V_{C}(r)$ $\displaystyle\approx$ $\displaystyle-\frac{\hbar
cR_{0}^{6}}{\pi}\sum_{n=0}^{9}C_{n}\frac{R_{0}^{n}}{L^{n+7}}\quantity(1+\frac{r}{L})^{-(n+7)}$
(16) $\displaystyle=$ $\displaystyle-\frac{\hbar
cR_{0}^{6}}{\pi}\sum_{n=0}^{9}C_{n}\frac{R_{0}^{n}}{L^{n+7}}\sum_{k=0}^{\infty}\binom{-(n+7)}{k}\frac{r^{k}}{L^{k}}$
$\displaystyle=$ $\displaystyle-\frac{\hbar cR_{0}^{6}}{\pi
L^{7}}\sum_{n=0}^{9}C_{n}\frac{R_{0}^{n}}{L^{n}}\sum_{k=0}^{\infty}(-1)^{k}\binom{n+k+6}{k}\frac{r^{k}}{L^{k}}.\hskip
28.45274pt$

Using the algorithm developed in Ref. [29] and discussed in the main text, we
find that the Casimir energy truncated at quadratic term ($k=2$) gives the
contribution to entanglement parameterized by

$\omega_{C}^{2}=\frac{2\hbar cR_{0}^{6}}{m\pi
L^{9}}\sum_{n=0}^{9}C_{n}(n+7)(n+8)\frac{R_{0}^{n}}{L^{n}}.$ (17)