# Universal role of curvature in vacuum entanglement

Hari K<EMAIL_ADDRESS>Subhajit Barman
<EMAIL_ADDRESS>Dawood Kothawala<EMAIL_ADDRESS>Centre
for Strings, Gravitation and Cosmology, Department of Physics, Indian
Institute of Technology Madras, Chennai 600 036, India

###### Abstract

We highlight some universal features concerning the role of spacetime
curvature in the entanglement induced between quantum probes coupled to a
quantum field in a suitable vacuum state. The probes are initially causally
disconnected and non-entangled. We explore the parameter space
$\\{{\omega},d_{0},\bm{v}_{0}\\}$ spanned by the energy gap $\omega$ of the
detectors, and the initial values of separation distance $d_{0}$ and relative
velocity $\bm{v}_{0}$, both covariantly defined in arbitrary curved spacetime.
We also obtain numerical results in de Sitter spacetimes and use these to
explore strong curvature regime, while also corroborating our perturbative
results in arbitrary curved spacetime. Our analysis shows that curvature can
induce entanglement features in certain regions of the above parameter space,
in a manner which facilitates using entanglement as a probe of spacetime
curvature.

###### Contents

  1. I Introduction
  2. II Geometric aspects of quantum probes coupled to a quantum field in curved spacetime
    1. II.1 Geodesic interval between two detectors in curved spacetime
    2. II.2 Green’s function in Hadamard representation
    3. II.3 Causal structure and the nature of entanglement induced between the probes
  3. III Mathematical details for computing the entanglement
    1. III.1 de Sitter background
    2. III.2 General curved spacetime
  4. IV Characteristics of $\mathcal{N}^{(2)}$ in de Sitter spacetime
  5. V Entanglement as a probe of curvature
  6. VI Concluding remarks
  7. A Geodesic interval between events on two timelike trajectories
  8. B Evaluation of $\mathcal{I}_{\varepsilon}$ in general curved spacetime: complete perturbative analysis
    1. B.0.1 Evaluation of $\mathcal{I}_{\varepsilon}$ for two initially static detectors
    2. B.0.2 Evaluation of $\mathcal{I}_{\varepsilon}$ for two geodesic detectors
  9. C 3D plots - dependence of negativity on velocity and angle

##  I Introduction

Entanglement – one of the key features of quantum mechanics – is also of
fundamental importance in studying phenomena at the interface of gravity and
quantum physics. Such an interface has conventionally been associated with the
domain of quantum gravity, but of course it does not have to be. Interesting
aspects associated with the study of low energy quantum systems in non-
singular curved backgrounds can in fact provide significant insights into what
one might expect in the extreme, quantum gravitational conditions. Work along
these lines has gained considerable attention over the past decade,
particularly in the light of universal role that gravity might play in quantum
processes such as decoherence and entanglement Pikovski et al. (2015), and the
question of whether this has any implication for the quantum nature of gravity
Bose et al. (2017); Marletto and Vedral (2017). In this work, we will focus on
identifying how entanglement between quantum probes, mediated through their
coupling to a quantum field, is affected by spacetime curvature. A deeper
question we would like to ask is whether a measure of entanglement can be used
to re-construct the curvature of spacetime. This would be in the spirit of the
broader formalism, developed by one of us and collaborators, which aims to
describe spacetime in terms of non-local, bi-tensorial objects such as Synge
world function and correlators; for recent reviews and relevant references,
see Kothawala (2023a, b).

Unruh-de Witt (UdW) detectors provide the simplest and the most operational
way to study such effects. It is well known that individual UdW detectors
probe vacuum fluctuations of the quantum field which they are coupled to. When
a pair of such detectors are employed, entanglement in the vacuum of the
quantum field can be transferred to that between the two detectors. What makes
this latter a more operationally and conceptually appropriate probe of the
structure of quantum field vacuum is that the entanglement induced in
detectors can be non-zero even for two inertial detectors coupled to a field
in Minkowski vacuum, provided they have some initial relative velocity Koga et
al. (2018). (A single such detector would not respond at all.) Using two
detectors as a probe of vacuum entanglement is an idea first introduced by
Reznik Reznik (2003), and has since been extensively studied in the
literature. (It is also referred to as “entanglement harvesting” in the
literature Koga et al. (2018, 2019); Martin-Martinez (2015); Tjoa and Martín-
Martínez (2021).) Many of these works have studied the nature of entanglement
induced in the detectors for various scenarios, such as the effects of motion
of the detectors Cong et al. (2020); Zhang and Yu (2020); Barman et al.
(2022a, b); Hu and Yu (2015); Menezes (2018), the effects of presence of a
thermal bath Brown (2013); Simidzija and Martín-Martínez (2018); Barman et al.
(2021), the dependence on spacetime dimensions Pozas-Kerstjens and Martin-
Martinez (2015), and the effect of a passing gravitational wave Xu et al.
(2020); Gray et al. (2021); Barman et al. (2023). In a related but different
context, radiative processes involving entangled detectors have also been
extensively studied Ostapchuk et al. (2012); Hu et al. (2012); Menezes and
Svaiter (2015, 2016); Arias et al. (2016); Menezes (2016); Liu et al. (2018);
Picanço et al. (2020); Barman and Majhi (2021); Barman et al. (2022c).

However, the effect of background curvature in all its generality, explicitly
in terms of the Riemann tensor, has not been studied so far. This is important
since any entanglement dynamics involving the detector-field system will be
affected by the background curvature, both explicitly, via additional terms in
the vacuum two-point correlators, and implicitly, through the geometrical
quantities appearing in the correlators pulled back onto the detector
trajectories. One might expect to gain some insight by restricting to a
special class of background spacetimes, taken as some exact solution of
Einstein equations which is simple enough for the analysis to be tractable and
yield analytical results Cliche and Kempf (2011); Kukita and Nambu (2017);
Tjoa and Mann (2020). Evidently, such attempts are limited by the fact that
the choice of a specific solution hides the manner in which Riemann curvature
explicitly appears and affects the final results, thereby yielding only
limited insights. We will comment further on these works vis-a-vis our results
in the final section. In this work, we will set-up the problem directly in
terms of Riemann tensor components, thereby highlighting explicitly their role
in various measures of entanglement, and uncovering interesting universal
features even from a perturbative analysis.

We will briefly describe here the key features of the setup conventionally
used to examine the possibility of two uncorrelated detectors getting
entangled over time and then quantify the measure of this entanglement. The
simplest such probes of vacuum structure of a quantum field are the Unruh-
DeWitt detectors, originally (and still widely) used to demonstrate and
understand various aspects of the Unruh effect Unruh (1976). A study of
entanglement involves two such detectors, 1 and 2, interacting with a massless
background scalar field $\phi(x^{k})$. We will assume these detectors to be
two-level systems with ground state $\lvert g\rangle$, an excited state
$\lvert e\rangle$, with energy difference $\omega$ between them. The
interaction action of this detector-field system is taken as

$\displaystyle\mathcal{A}_{\rm int}=\int
c_{0}\,\chi(\tau)\,m_{1}(\tau)\,\phi(x_{1}(\tau))\,d\tau+\int
c_{0}\,\chi(s)\,m_{2}(s)\,\phi(x_{2}(s))\,ds~{}.$ (1)

Here, $\chi$ is the switching function, and $m_{1,2}$ represent the monopole
operators of the detectors. Let $\lvert\rm in\rangle$ be the product state in
the asymptotic past, $\lvert{\rm in}\rangle=\lvert 0\rangle\lvert
g\rangle_{1}\lvert g\rangle_{2}$. Here, $\lvert 0\rangle$ is the vacuum state
of the field. In the asymptotic future, the state will evolve into,
$\lvert{\rm out}\rangle=T\big{\\{}e^{i\mathcal{A}_{\rm{int}}}\lvert{\rm
in}\rangle\big{\\}}$, with $T$ being the time ordering operator. The reduced
density matrix corresponding to the detector system is obtained by tracing
over all the final field states, i.e., $\rho_{12}={\rm Tr}_{\phi}\,\lvert{\rm
out}\rangle\langle{\rm out}\rvert$. Introducing $c_{0}$ as a book-keeping
variable for the interaction strength, one can express the reduced density
matrix in the bases $\left\\{\;\lvert g\rangle_{1}\lvert g\rangle_{2},\;\lvert
e\rangle_{1}\lvert g\rangle_{2},\;\lvert g\rangle_{1}\lvert
e\rangle_{2},\;\lvert e\rangle_{1}\lvert e\rangle_{2}\;\right\\}$ (ref. Koga
et al. (2018)) as

$\displaystyle\rho_{12}=\begin{bmatrix}1-c_{0}^{2}\left(\mathcal{P}_{1}+\mathcal{P}_{2}\right)&0&0&c_{0}^{2}\mathscr{E}^{*}\\\
0&c_{0}^{2}\mathcal{P}_{1}&c_{0}^{2}\mathcal{P}_{12}&0\\\
0&c_{0}^{2}\mathcal{P}_{12}&c_{0}^{2}\mathcal{P}_{2}&0\\\
c_{0}^{2}\mathscr{E}&0&0&0\\\ \end{bmatrix}+\mathcal{O}(c_{0}^{4})~{}.$ (2)

For understanding entanglement, the quantities we will need are
$\mathcal{P}_{1},\mathcal{P}_{2}$ and $\mathscr{E}$ Koga et al. (2018), and
these are given by

$\displaystyle\mathcal{P}_{1}$ $\displaystyle=$
$\displaystyle\left\lvert\langle e\rvert m_{1}(0)\lvert
g\rangle_{1}\right\rvert^{2}\,\mathcal{I}_{1}\hskip 28.45274pt({\rm
similarly\;for\;2})$ (3a) $\displaystyle\mathscr{E}$ $\displaystyle=$
$\displaystyle\langle e\rvert m_{2}(0)\lvert g\rangle_{2}\;\langle e\rvert
m_{1}(0)\lvert g\rangle_{1}\,\mathcal{I}_{\varepsilon}$ (3b)

where

$\displaystyle\mathcal{I}_{1}$ $\displaystyle=$
$\displaystyle\int_{-\infty}^{\infty}d\tau^{\prime}\int_{-\infty}^{\infty}d\tau\;\chi(\tau)\chi(\tau^{\prime})e^{i\omega(\tau-\tau^{\prime})}G_{W}(x_{1}(\tau^{\prime}),x_{1}(\tau))$
(4a) $\displaystyle i\,\mathcal{I}_{\varepsilon}$ $\displaystyle=$
$\displaystyle\int_{-\infty}^{\infty}ds\int_{-\infty}^{\infty}d\tau\;\chi(\tau)\chi(s)e^{i\omega(\tau+s)}G_{F}(x_{2}(s),x_{1}(\tau))~{}.$
(4b)

Here, $G_{W}(x^{\prime},x)$ is the Wightmann function, and
$G_{F}(x_{2},x_{1})$ is the Feynman propagator. The quantities
$\mathcal{P}_{1,2}$ are local to each detector, and correspond to individual
detector responses, whereas $\mathscr{E}$ represents the non-local correlation
between the detectors.

Now according to Peres (1996), a state is entangled if the partial
transposition of the relevant density matrix has at least one negative
eigenvalue. With the reduced detector density matrix from Eq. (2), this
requirement leads to the condition:

$\displaystyle\mathcal{P}_{1}\mathcal{P}_{2}<|\mathscr{E}|^{2}~{}$ (5)

which can be rewritten as

$\mathcal{I}_{1}\mathcal{I}_{2}<|\mathcal{I}_{\varepsilon}|^{2}$ (6)

where, $\mathcal{I}_{1}$ and $\mathcal{I}_{2}$ are the individual detector
responses. The above condition can be used to analyze the possibility of two
Unruh-DeWitt detectors coupled to a free scalar field getting entangled. It is
worth highlighting that the non-local term (RHS) above depends on the Feynman
propagator, while the individual responses (LHS) depend on the Wightman
function.

Several measures quantify the entanglement induced between the detectors by
providing a bound on this entanglement in the system. For example,
entanglement negativity $\mathcal{N}(\rho_{12})$ is defined as the sum of all
negative eigenvalues of the partial transpose of the density matrix
$\rho_{12}$. One can also define the logarithmic negativity
$\log_{2}{(\mathcal{N}(\rho_{12})+1)}$, see Koga et al. (2018), which provides
an upper bound on the induced entanglement. Another very relevant measure of
entanglement is concurrence $\mathcal{C}(\rho_{12})$ Koga et al. (2018). With
the reduced density matrix from Eq. (2), we have the expressions for the
negativity and the concurrence (to lowest order in perturbation, up to
$c_{0}^{2}$ terms), respectively as

$\displaystyle\mathcal{N}(\rho_{12})$ $\displaystyle=$ $\displaystyle
c^{2}_{0}\,\Big{[}\sqrt{(\mathcal{P}_{1}-\mathcal{P}_{2})^{2}+4\,|\mathscr{E}|^{2}}-\mathcal{P}_{1}-\mathcal{P}_{2}\Big{]}/2$
(7a) $\displaystyle~{}\mathcal{C}(\rho_{12})$ $\displaystyle=$ $\displaystyle
2\,c^{2}_{0}\,\Big{[}|\mathscr{E}|-\sqrt{\mathcal{P}_{1}\,\mathcal{P}_{2}}\Big{]}~{}.$
(7b)

The expression of the monopole moment operator (for each detector) is given by

$\displaystyle m(0)=\lvert e\rangle\langle g\rvert+\lvert g\rangle\langle
e\rvert~{}.$ (8)

Then one can observe that the expectation values of these operators as used in
Eq. (3) is $\langle e\rvert m(0)\lvert g\rangle=1$. Utilizing this
observation, in the special scenario of
$\mathcal{I}_{1}=\mathcal{I}_{2}=\mathcal{I}$, one can express the negativity
and the concurrence as
$\mathcal{N}(\rho_{12})=\mathcal{C}(\rho_{12})/2=\big{[}|\mathcal{I}_{\varepsilon}|-\mathcal{I}\big{]}c^{2}_{0}$,
taking terms upto $\mathcal{O}(c_{0}^{2})$. This is the case that we shall
encounter in our system. Therefore, it is convenient to define negativity
specifically as $\mathcal{N}^{(2)}=|\mathcal{I}_{\varepsilon}|-\mathcal{I}$ to
denote a measure of the entanglement in our system.

The rest of the manuscript is organised as follows: In Sec. II, we provide a
description for some of the quantities in curved spacetime, like Green’s
functions in terms of the geodesic distance, that are essential in realizing
the phenomenon of entanglement. In the subsequent Sec. III.1, we consider the
maximally symmetric de Sitter background and investigate the characteristics
of negativity. Following in Sec. III.2, we study the nature of negativity in a
generalized curved spacetime with a relatively small curvature compared to the
distance between the two detectors, i.e., the curvature length scale is larger
than the initial separation between the detectors. We discuss our observed
outcomes regarding negativity in a de Sitter background and elucidate their
physical interpretations in Sec. IV. In Sec. V, we describe how our results
quantifying the effect of curvature on entanglement can be utilized as a tool
to probe curvature of arbitrary spacetime, in the spirit of the comment made
at the end of the first paragraph of this section. Finally, we summarize with
a concluding discussion in Sec. VI.

##  II Geometric aspects of quantum probes coupled to a quantum field in
curved spacetime

As elucidated in the previous section, the criteria for entanglement from Eq.
(6) and its quantification depends on two specific kinds of terms in the
expression for negativity: the local terms involving $\mathcal{I}_{1,2}$ that
denote individual detector transition probabilities, and the non-local term
$\mathcal{I}_{\varepsilon}$ that signifies the correlation between the two
detectors. To obtain explicit expressions of these quantities, an
understanding of different types of Green’s functions is necessary. These
Green’s functions will again contain the information about the detector
trajectories and the background in them. In this section, we intend to provide
general expressions of Green’s functions in curved background. In particular,
we shall start by providing a general expression for the geodesic interval
between two geodesic detectors and delineate the necessary Green’s functions
in terms of this geodesic distance. In this section, we shall also discuss the
curved background’s causal structure and compare it to the Minkowski
spacetime. As we will see, it has a close correspondence with the nature of
Green’s functions.

###  II.1 Geodesic interval between two detectors in curved spacetime

The geometric set-up and kinematical parameters associated with two detectors
(timelike curves) is described in Fig. 1. To obtain the geodesic interval
between the two detectors, we construct Fermi Normal Coordinates(FNC) along
the curve $\mathscr{C}_{1}$ and identify the point on curve $\mathscr{C}_{2}$
by sending out space-like geodesics orthogonal to curve $\mathscr{C}_{1}$. Let
the first detector follows the curve $\mathscr{C}_{1}$ with proper time $\tau$
and be at $A$ and the second detector follows the curve $\mathscr{C}_{2}$ with
proper time $s$ be at $B$. Let $\tau_{2}^{\prime}$ be the time component of
the FNC on $\mathscr{C}_{1}$ of the detector $B$ at $B^{\prime}$. Now, the
geodesic interval between $A$ and $B$ can be expressed as a series expansion
in terms of the spacelike geodesic interval between $B$ and $B^{\prime}$, and
the geodesic interval between $A$ and $B^{\prime}$. The point $O$ with
$\tau=0$ provides the Cauchy data. Employing a generalized variant of an
expansion scheme proposed by Synge, we can obtain the geodesic interval
between arbitrary events $A$ and $B$ on the two probe trajectories using
coincidence limits of covariant derivatives of Synge’s world function. This
rather long calculation is described in Appendix A. We introduce the unit
vector $\xi^{\alpha}$ (representing direction cosines) defined by
${X}^{\alpha}(0):=d_{0}\,\xi^{\alpha}$ where, $\xi^{\alpha}$ (the directional
cosines) satisfy $\eta_{\alpha\beta}\xi^{\alpha}\xi^{\beta}=1$, and the
parallel and perpendicular components of initial $v_{0}^{\alpha}$ w.r.t.
$\xi^{\alpha}$ by $v_{0}^{\alpha}=v_{||}\,\xi^{\alpha}+v_{\perp}\,n^{\alpha}$,
with $n^{\alpha}$ an arbitrary unit vector in the 2-plane orthogonal to
$\bm{\xi}$ (and, of course, $\bm{u}_{1}$). (Note that $v_{||}=v_{0}\cos\theta$
and $v_{\perp}=v_{0}\sin\theta$.)

Figure 1: Schematic diagrams for two point-like detectors in geodesic
trajectories in general curved spacetime. In left we depict the two detectors
as denoted by $u_{1}^{i}$ and $u_{2}^{i}$. On the $\tau=0$ hypersurface, the
shortest distance between these detectors is $d_{0}$. Note the angle between
this distance and the velocity vector of the second detector is $\theta$. On
the right, we depict the same scenario, which now corresponds to the
construction of the Fermi coordinates for the two detectors with respect to
the trajectory of detector $A$. Here, $\Omega_{AB}\equiv\Omega_{\mathcal{G}}$,
$\Omega_{AB^{\prime}}\equiv\Omega_{\mathcal{G}^{\prime}}$, and
$\Omega_{B^{\prime}B}\equiv\Omega_{\mathcal{G}s}$. Also note that
$\Omega_{\mathcal{G}s}(\tau=0)=d_{0}$ denotes the separation between the two
detectors at the $\tau_{1}=0$ hypersurface.

The final expression in terms of these quantities can be written in a compact
form

$\displaystyle-\sigma^{2}_{\mathcal{G}}\;=\;\Delta\tau_{\mathcal{G}}^{2}\;=\;\left(\tau_{2}^{\prime}-\tau\right)^{2}\mathpzc{f}(\tau_{2}^{\prime})-\mathpzc{g}(\tau_{2}^{\prime})+\mathcal{O}(\nabla
R,R^{2})~{}$ (9)

where the auxiliary functions $\mathpzc{f}(\tau_{2}^{\prime})$ and
$\mathpzc{g}(\tau_{2}^{\prime})$ are given by

$\displaystyle\mathpzc{f}(\tau_{2}^{\prime})$ $\displaystyle=$ $\displaystyle
1+\frac{1}{3}\left[R_{\xi 0\xi 0}d_{0}^{2}+2d_{0}\left(R_{\xi 0\xi
0}v_{||}+R_{\xi 0n0}v_{\perp}\right)\tau_{2}^{\prime}+\left(R_{\xi 0\xi
0}v_{||}^{2}+2R_{\xi
0n0}v_{||}v_{\perp}+R_{n0n0}v_{\perp}^{2}\right){\tau_{2}^{\prime}}^{2}\right]\hskip
28.45274pt$ (10a) $\displaystyle\mathpzc{g}(\tau_{2}^{\prime})$
$\displaystyle=$
$\displaystyle\left(d_{0}^{2}+v_{0}^{2}{\tau_{2}^{\prime}}^{2}+2d_{0}v_{||}\tau_{2}^{\prime}\right)+d_{0}^{2}\,\left(\frac{2}{3}R_{n\xi
0\xi}v_{\perp}-R_{\xi 0\xi 0}+\frac{1}{3}R_{n\xi
n\xi}v_{\perp}^{2}\right){\tau_{2}^{\prime}}^{2}$ (10b)
$\displaystyle-\frac{4}{3}d_{0}\left(R_{\xi 0\xi 0}v_{||}+R_{\xi
0n0}v_{\perp}\right){\tau_{2}^{\prime}}^{3}-\frac{1}{3}\left(R_{\xi 0\xi
0}v_{||}^{2}+2R_{n0\xi
0}v_{||}v_{\perp}+R_{n0n0}v_{\perp}^{2}\right){\tau_{2}^{\prime}}^{4}$

To improve readability, we have introduced the following compact notation:
$R_{\xi 0\xi 0}:=R_{\alpha i\beta j}\xi^{\alpha}\xi^{\beta}u^{i}u^{j},R_{\xi
0n0}:=R_{\alpha i\beta j}\xi^{\alpha}n^{\beta}u^{i}u^{j},R_{n\xi
0\xi}:=R_{\alpha\mu i\nu}n^{\alpha}\xi^{\mu}\xi^{\nu}u^{i}$ and so on for the
other Riemann tensor components.

It is easy to recover the following results from the above expression:

Flat spacetime:

$\Delta\tau_{\mathcal{G}}^{2}=\left(\tau_{2}^{\prime}-\tau\right)^{2}-\Biggl{[}\left(d_{0}+v_{\parallel}\tau_{2}^{\prime}\right)^{2}+v_{\perp}^{2}{\tau_{2}^{\prime}}^{2}\Biggl{]}$
(11)

Initially static detectors; $v_{0}=0$:

When $v_{0}=0$, the above equation simplifies to

$\Delta\tau_{\mathcal{G}}^{2}=\left(\tau_{2}^{\prime}-\tau\right)^{2}\Biggl{[}1+\frac{1}{3}R_{\xi
0\xi 0}d_{0}^{2}\Biggl{]}-d_{0}^{2}\Biggl{[}1+R_{\xi 0\xi
0}{\tau_{2}^{\prime}}^{2}\Biggl{]}$ (12)

To compute the non-local term in entanglement measures, we also need to map
the proper time $s$ of trajectory $\mathscr{C}_{2}$ to the Fermi time
$\tau_{2}^{\prime}$ based on $\mathscr{C}_{1}$’s trajectory. This is best done
using the Fermi metric, and the derivation is given in Appendix A. The final
result is:

$\displaystyle s$ $\displaystyle=$
$\displaystyle\frac{\tau_{2}^{\prime}}{\gamma_{0}}+\frac{\gamma_{0}}{2}\tau_{2}^{\prime}\left[\mathpzc{h}_{1}(v_{0},d_{0})+\tau_{2}^{\prime}\left\\{d_{0}\,\mathpzc{h}_{2}(v_{0},d_{0})+\frac{\tau_{2}^{\prime}}{3}\mathpzc{h}_{3}(v_{0},d_{0})\right\\}\right]+\mathcal{O}(R^{2},\nabla
R)~{}.$ (13)

where, $\gamma_{0}=1/\sqrt{1-v_{0}^{2}}$ and
$\mathpzc{h}_{1},\mathpzc{h}_{2},\mathpzc{h}_{3}$ are given in Appendix A.
When the initial velocity of the second observer measured with respect to the
first one’s FNC vanishes, we have $v_{0}=0$ and also $v_{||}=0=v_{\perp}$. In
that scenario, the relation between $s$ and $\tau_{2}^{\prime}$ boils down to
$s=\tau_{2}^{\prime}(1+R_{\xi 0\xi 0}\,d_{0}^{2}/2)$. We also mention that to
obtain the expression of Eq. (13), we have resorted to expanding $\gamma$ from
(49) perturbatively, assuming small curvature terms. Let us also consider
$s=\tau_{2}^{\prime}/\bar{\gamma}$, where the expression of $\bar{\gamma}$ is
evident from the previous equation.

###  II.2 Green’s function in Hadamard representation

The Hadamard representation of the Feynman propagatorHadamard (2003); Decanini
and Folacci (2006) for KG equation of massless scalar field in 4 dimensions is
given by,

$G_{F}\left(x,x^{\prime}\right)=\frac{i}{4\pi^{2}}\left[\frac{\Delta^{1/2}}{\sigma^{2}(x,x^{\prime})+i\epsilon}+V_{n}(x,x^{\prime})\ln\left(\sigma^{2}(x,x^{\prime})+i\epsilon\right)+W_{n}(x,x^{\prime})\right]~{}.$
(14)

Here, the functions, $V_{n}(x,x^{\prime})$ and $W_{n}(x,x^{\prime})$ are
expressed as series,
$V_{n}(x,x^{\prime})=\sum_{n=0}^{\infty}v_{n}(x,x^{\prime})\sigma^{2n}(x,x^{\prime})$,
$W_{n}(x,x^{\prime})=\sum_{n=0}^{\infty}w_{n}(x,x^{\prime})\sigma^{2n}(x,x^{\prime})$,
and $-(1/2)\sigma^{2}(x,x^{\prime})\equiv\Omega(x,x^{\prime})$ is half the
square of geodesic distance also known as Synge’s world functionSynge (1960).

Similarly, for Wightman function we use the formKay and Wald (1991),

$\displaystyle
G^{+}\left(x,x^{\prime}\right)=\frac{1}{4\pi^{2}}\left[\frac{\Delta^{1/2}}{\sigma^{2}_{\epsilon}(x,x^{\prime})}+V_{n}(x,x^{\prime})\ln\left(\sigma^{2}_{\epsilon}(x,x^{\prime})\right)+W_{n}(x,x^{\prime})\right]~{}.$
(15)

Here,
$\sigma^{2}_{\epsilon}:=\sigma^{2}(x,x^{\prime})+i\epsilon\left[T(x)-T(x^{\prime})\right]+\epsilon^{2}$,
with $T(x)$ being any global time function that increases towards the future
and other parameters are similar to that defined for Feynman propagator.

The Hadamard coefficient, $V_{n}(x,x^{\prime})$ obeys certain recursion
relation with specific boundary condition and can be solved, whereas the
boundary condition $W_{0}(x,x^{\prime})$ for $W_{n}(x,x^{\prime})$ is not
restrained by the recursion relations. But the Hadamard coefficients,
$W_{n}(x,x^{\prime})$ can be determined once $W_{0}(x,x^{\prime})$ is
specified. For our evaluation of the integral, we are interested in the
leading $1/\sigma^{2}$ term since this term determines the pole structure of
the Feynman propagator, and $\Delta^{1/2}$ will not contribute to the poles.
We also approximate, $\Delta^{1/2}\approx 1$. Hence, the geodesic interval
between the two detectors is required to evaluate the Feynman propagator.

_The van Vleck determinant_ :

For both single detector probability and non-local term in the density matrix,
the Hadamard form of the Wightmann function and the Feynman propagator have
van Vleck determinant in the expression. The van Vleck determinant can be
expanded in terms of geodesic distance as in Ref. Christensen (1976),

$\Delta^{1/2}=1+\frac{1}{12}R_{ab}\sigma^{a}\sigma^{b}-\frac{1}{24}\nabla_{c}R_{ab}\sigma^{a}\sigma^{b}\sigma^{c}+\mathcal{O}(\sigma^{4})$
(16)

Note that $\sigma^{a}(x,x^{\prime})=-(\Delta\tau_{\rm geod})u^{a}$, where
$\Delta\tau$ is the geodesic interval and $u^{a}$ is the tangent vector of the
geodesic. When this expansion is used in the Hadamard form,

$\displaystyle\frac{\Delta^{1/2}}{4\pi^{2}\sigma^{2}_{\epsilon}}\approx$
$\displaystyle\frac{1}{4\pi^{2}\sigma^{2}_{\epsilon}}+\frac{(1/12)R_{ab}\sigma^{a}\sigma^{b}}{4\pi^{2}\sigma^{2}_{\epsilon}}+\mathcal{O}(\sigma^{3})$
$\displaystyle\approx$
$\displaystyle\frac{1}{4\pi^{2}\sigma^{2}_{\epsilon}}+\frac{(1/12)R_{uu}\sigma^{2}}{4\pi^{2}\sigma^{2}_{\epsilon}}+\mathcal{O}(\sigma^{3})$

Hence the contribution from the curvature through the van Vleck determinant
will be only an additive constant. The additional corrections from derivatives
of curvature and higher-order terms in the curvature are ignored.

###  II.3 Causal structure and the nature of entanglement induced between the
probes

The non-local entangling term from Eq. (4b) will depend on different
parameters such as $d_{0},v_{0},\theta,$ and different curvature components.
To compare the contribution of genuine entanglement and entanglement through
field-mediated communication channel, we express the non-local term in the
form

$\displaystyle\mathcal{I}_{\varepsilon}=-\mathcal{I}_{\varepsilon}^{+}-\mathcal{I}_{\varepsilon}^{-}~{},$
(17)

where $\mathcal{I}_{\varepsilon}^{+}$ and $\mathcal{I}_{\varepsilon}^{-}$
depend exclusively on the expectation values of the field anti-commutator and
field commutators. The part $\mathcal{I}_{\varepsilon}^{-}$ depending on the
field commutator, is considered in the literature to signify the entanglement
due to the causal connection. Whereas $\mathcal{I}_{\varepsilon}^{+}$
corresponds to the contribution that signifies entanglement in the true sense.
In particular, these quantities can be explicitly expressed in terms of the
Wightman functions $G_{W}(\bar{u},\bar{v})$ as

$\displaystyle\mathcal{I}_{\varepsilon}^{+}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\int_{-\infty}^{\infty}d\bar{u}\int_{-\infty}^{\infty}d\bar{v}\,e^{i\,\omega\,v}\chi_{A}(\bar{u},\bar{v})\,\chi_{B}(\bar{u},\bar{v})\,\text{Re}[G_{W}(\bar{u},\bar{v})]~{},$
(18a) $\displaystyle\mathcal{I}_{\varepsilon}^{-}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\int_{-\infty}^{\infty}d\bar{u}\int_{-\infty}^{\infty}d\bar{v}\,e^{i\,\omega\,v}\chi_{A}(\bar{u},\bar{v})\,\chi_{B}(\bar{u},\bar{v})\,\left[1+2\,\Theta(\bar{u})\right]\,\text{Im}[G_{W}(\bar{u},\bar{v})]~{}.$
(18b)

Here, $\Theta(\bar{u})$ denotes the Heaviside step function, which arrives
from the decomposition of the Feynman propagator in terms of time-ordered
Wightman functions, e.g.,
$i\,G_{F}(x,x^{\prime})=\Theta(t-t^{\prime})G_{W}(x,x^{\prime})+\Theta(t^{\prime}-t)G_{W}(x^{\prime},x)$.
In the special case, when under the transformation of $\bar{v}\to-\bar{v}$ the
quantities $\chi_{A}(\bar{u},\bar{v})$, $\chi_{B}(\bar{u},\bar{v})$, and
$G_{W}(\bar{u},\bar{v})$ remain unchanged, the above two quantities
$\mathcal{I}_{\varepsilon}^{\pm}$ are simply given by

$\displaystyle\mathcal{I}_{\varepsilon}^{+}$ $\displaystyle=$
$\displaystyle\text{Re}[\mathcal{I}_{\varepsilon}]~{},~{}~{}\textup{and}~{}~{}\mathcal{I}_{\varepsilon}^{-}=\text{Im}[\mathcal{I}_{\varepsilon}]~{}.$
(19)

The Wightman function in Minkowski spacetime with the Gaussian switching
satisfies this condition. Furthermore, we shall see that the Wightman function
is de Sitter background also satisfies this condition.

On the other hand, in situations when at least one of the quantities
$\chi_{A}(\bar{u},\bar{v})$, $\chi_{B}(\bar{u},\bar{v})$, or
$G_{W}(\bar{u},\bar{v})$ is not invariant under the transformation of
$\bar{v}\to-\bar{v}$, one cannot simply use the expressions from Eq. (19) for
the evaluation of $\mathcal{I}_{\varepsilon}^{\pm}$. In this scenario, one
should use the general expressions from Eq. (18). We should mention that the
Green’s functions in a generalized curved background with the geodesic
interval from Eq. (II.1) do not remain invariant under the aforementioned
transformation. Therefore, when investigating the entanglement in a
generalized curved spacetime, one has to use the expressions from Eq. (18) to
understand the natures of $\mathcal{I}_{\varepsilon}^{\pm}$.

Furthermore, as discussed in Martin-Martinez (2015); Tjoa and Martín-Martínez
(2021) one can define suitable quantities denoting the spacelike and
communication-based entanglement as

$\displaystyle\mathcal{N}^{\pm}=|\mathcal{I}_{\varepsilon}^{\pm}|-\mathcal{I}~{}.$
(20)

Here also, the superscript ${}^{\prime}+^{\prime}$ or ${}^{\prime}-^{\prime}$
respectively correspond to the spacelike and communication channel-mediated
entanglement. In the later part of our present analysis, we shall observe how
the light-cone diagrams from Fig. 3 are closely related and give a physical
interpretation to the quantitative values of these specific components
$\mathcal{N}^{\pm}$.

Figure 2: The overlap of the causal diamond formed between the switch on and
switch off time for each detector determines the contribution from genuine
entanglement and field-mediated communication. We will choose our parameter
values for the detector setup such that this overlap is minimal but the value
of negativity is large enough to see the effects due to curvature.

The commutator and anti-commutator can be used to characterize the
entanglement through field-mediated communication channels and the
entanglement through non-local correlations. The contribution from these
correlators will depend on the causal structure of the spacetime and the given
detector configuration. The overlap of the causal diamonds $\mathcal{D}_{1,2}$
determined by the switching functions (see Figs. 2, 3) can be used to estimate
these contributions. Any non-zero intersection of the causal diamonds would
indicate possible entanglement due to the communication channel will be
present. Hence, if the detectors are placed far apart to ensure there is
little to no overlap between causal diamonds, the genuine entanglement can be
obtained.

Figure 3: The above figure represents schematic diagrams for the causal
structure in a general curved background with two different detector
trajectories. We have used Fermi normal coordinates to represent the curves
and light cones to compare it with Minkowski. As the distance and time
progresses, the curved spacetime light cones deviate from the Minkowski light
cones. We also observe that the two detector trajectories move away from each
other due to curvature. With finite switching, the first detector $A$ can
access a lesser portion of the lightcone on the trajectory of detector $B$ as
it moves away.

##  III Mathematical details for computing the entanglement

###  III.1 de Sitter background

We start our investigation on entanglement in a background that is a maximally
symmetric solution of Einstein’s equations with a positive cosmological
constant, i.e., in the de Sitter spacetime. The motivation behind considering
the de Sitter spacetime is twofold. First, the de Sitter era signifies an
accelerated expansion of the universe among the general homogeneous and
isotropic descriptions given by the Friedmann–Lemaître–Robertson–Walker (FLRW)
metric. Our universe, as observed, is in accelerated expansion in its current
state, and the study of the de Sitter background becomes relevant from this
point of view. Second, obtaining an exact expression of the geodesic interval
between two geodesic observers is straightforward in de Sitter spacetime due
to its maximal symmetry.

In de Sitter spacetime, the line element Deser and Levin (1997) in standard
coordinates is expressed as

$\displaystyle ds^{2}$ $\displaystyle=$ $\displaystyle-
dt^{2}+R^{2}\cosh^{2}{(t/R)}\Big{[}d\chi^{2}+\sin^{2}{(\chi)}\big{\\{}d\theta^{2}+\sin^{2}{(\theta)}\,d\phi^{2}\big{\\}}\Big{]}~{}.$
(21)

For static observers in this coordinates of de Sitter spacetime, the
coordinate time $t$ acts as a proper time. We utilize the embedding
coordinates for the de Sitter spacetime to perceive a geodesic distance
between two static detectors. We use these embedding coordinates along with
the considerations of constant $\theta$ and constant $\phi$ plane to obtain
the geodesic distance. In the de Sitter background, with the identification of
$\Delta\chi=d_{0}/R$, this geodesic distance is obtained to be

$\displaystyle\sigma^{2}_{DS}=\Lambda^{-1}~{}\cos^{-2}\Big{\\{}-\sinh(\tau\sqrt{\Lambda})\sinh(s\sqrt{\Lambda})+\cosh(\tau\sqrt{\Lambda})\cosh(s\sqrt{\Lambda})\cos(d_{0}\sqrt{\Lambda})\Big{\\}}~{},$
(22)

where $\Lambda=1/R^{2}$ and the Ricci scalar in this background is
$\mathcal{R}=12\,\Lambda$. In a de Sitter background, considering a massless
scalar field with conformal coupling, one can obtain the Green’s function, see
Garbrecht and Prokopec (2004), as,

$\displaystyle G(x,x^{\prime})$ $\displaystyle=$
$\displaystyle-\frac{i}{4\pi^{2}}\,\frac{\Lambda}{4\sin^{2}{(\sqrt{\Lambda}\,\sigma_{DS}/2)}}~{}.$
(23)

For the coincidence limit, Green’s function reduces to Hadamard form and with
$d_{0}\to 0$ in Eq. (22), the Wightman function that characterizes the single
detector transition can be obtained with appropriate $i\epsilon$ prescription.

The local and non-local terms in negativity

We use the expression of $\sigma^{2}_{DS}$ from Eq. (22) to get the expression
of the Wightman function (23). One can notice that for a single detector
$t_{A}=t_{B}$ and $d_{0}=0$. In that scenario, with the assumption
$\tau-\tau^{\prime}=u$ and $\tau+\tau^{\prime}=v$, we can express the local
individual detector transition energy denoting terms $\mathcal{I}$ in the de
Sitter background as

$\displaystyle\mathcal{I}=-\frac{\Lambda}{32\pi^{2}}\int_{-\infty}^{\infty}dv\int_{-\infty}^{\infty}du\,\frac{e^{-i\,\omega\,u}e^{-\frac{u^{2}+v^{2}}{2\bar{T}^{2}}}}{\sinh^{2}{\big{\\{}(\sqrt{\Lambda}\,u/2)-i\,\epsilon\big{\\}}}}~{}.$
(24)

Here one can perform the $v$ integration analytically, which results in
$\int_{-\infty}^{\infty}dv\,e^{-v^{2}/(2\bar{T}^{2})}=\sqrt{2\pi}\,\bar{T}$.
To perform the $u$ integration, we have used numerical methods. In Fig. 4, we
have plotted these local terms denoting the individual detector transition
probabilities. We elaborately discuss the characteristics of these local terms
in a later section of our manuscript.

On the other hand, we use the relation
$\sin^{2}\big{\\{}{(\cos^{-1}{z})/2}\big{\\}}=(1-z)/2$ to get a simplified
expression of the Feynman propagator from Eq. (23). Using the expression of
the Feynman propagator, and with a change of variables to $\bar{u}=s-\tau$ and
$\bar{v}=s+\tau$, we get the non-local entangling term as,

$\displaystyle\mathcal{I}_{\varepsilon}$ $\displaystyle=$
$\displaystyle\frac{\Lambda}{8\pi^{2}}\int_{-\infty}^{\infty}d\bar{v}\int_{-\infty}^{\infty}d\bar{u}\,\frac{e^{i\,\omega\,\bar{v}}e^{-\frac{\bar{u}^{2}+\bar{v}^{2}}{2\bar{T}^{2}}}}{2-\alpha_{+}\,\cosh{(\sqrt{\Lambda}\,\bar{u})}-\alpha_{-}\,\cosh{(\sqrt{\Lambda}\,\bar{v})}+i\,\epsilon\,\Lambda}~{},$
(25)

where, we have the expressions of
$\alpha_{\pm}=\cos{(d_{0}\,\sqrt{\Lambda})}\pm 1$ . Here, we evaluate both the
$u$ and $v$ integrals numerically.

In Fig. 5 we plot the negativity in a de Sitter background, which is the
subtraction of the modulus of this non-local entangling term and the local
terms. In a later section, we discuss elaborately the characteristics of these
plots. Subsequently, in Fig. 6 we have plotted the difference in negativity
between the de Sitter and the Minkowski case for a very small curvature length
scale. The discussion on its characteristics and its similarity with the
general curved case is discussed in the same later section. We also plot the
different components of negativity and understand their characteristics
through the plots from Figs. 7 and 8.

###  III.2 General curved spacetime

As discussed previously, the formulation to understand the effects of tidal
curvature on entanglement constrains the interaction time between the
detectors and the background field. Then it will be suitable to consider a
specific switching function $\chi(\tau)$, through which we can implement this
constraint. In particular, introducing a nontrivial switching function will
complicate the computations, but it will provide physically more relevant and
consistent results with the mathematical formulations. In this regard, we
consider the Gaussian switching $\chi(\tau)=e^{-\tau^{2}/\bar{T}^{2}}$, and
shall evaluate the local $\mathcal{I}$ and non-local
$\mathcal{I}_{\varepsilon}$ terms in entanglement.

The local and non-local terms in negativity:

Evaluation of $\mathcal{I}$:

Let us first evaluate the local terms $\mathcal{I}$ in the concurrence, which
corresponds to the individual detector transition probabilities. We take the
expression of this integral from Eq. (4a), and shall be using the Gaussian
switching functions in its evaluation. One can notice that the integral looks
simplified with a change of variables to $u=\tau^{\prime}-\tau$ and $\tau$. In
particular, with this change of variables, this integral looks like

$\displaystyle\mathcal{I}=-\frac{1}{4\pi^{2}}\int_{-\infty}^{\infty}d\tau\,\int_{-\infty}^{\infty}du\,\frac{e^{-i\omega\,u}}{(u-i\epsilon)^{2}}\,e^{-[\tau^{2}+(\tau+u)^{2}]/\bar{T}^{2}}~{}.$
(26)

Here the $\tau$ integral is easily doable and results in

$\displaystyle\mathcal{I}=-\frac{1}{4\pi^{2}}\sqrt{\frac{\pi}{2}}\bar{T}\int_{-\infty}^{\infty}du\,\frac{e^{-\left[i\omega\,u+u^{2}/(2\bar{T}^{2})\right]}}{(u-i\epsilon)^{2}}~{}.$
(27)

To analytically perform this last integral, we specifically Fourier transform
the factor
$e^{-u^{2}/(2\bar{T}^{2})}=\bar{T}/(\sqrt{2\pi})\int_{-\infty}^{\infty}d\kappa\,e^{i\kappa
u-\xi^{2}\bar{T}^{2}/2}$. Then we shall get a factor of
$e^{i(\kappa-{\omega})u}$ in the numerator and $(u-i\epsilon)^{2}$ in the
denominator. The integral will have a pole of order two at $u=i\epsilon$. This
integral can be performed through contour integration, with non-vanishing
contributions when the contour is taken in the upper half complex plane. When
we take the contour in the upper half complex plane, we must satisfy
$\kappa-{\omega}>0$ for the damping of other integration. This leads to the
expression:

$\displaystyle\mathcal{I}$ $\displaystyle=$
$\displaystyle-\frac{\bar{T}^{2}}{8\pi^{2}}\int_{{\omega}}^{\infty}d\kappa\,e^{-\kappa^{2}\bar{T}^{2}/2}\,(i(\kappa-{\omega})e^{{\omega}\epsilon-\kappa\epsilon})2\pi
i~{}$ (28) $\displaystyle\simeq$
$\displaystyle\frac{\bar{T}^{2}}{4\pi}\int_{{\omega}}^{\infty}d\kappa\,e^{-\kappa^{2}\bar{T}^{2}/2}\,(\kappa-{\omega})~{}.$

After carrying out this integration over $\xi$, we get

$\displaystyle\mathcal{I}$ $\displaystyle=$
$\displaystyle\frac{1}{4\pi}\bigg{[}e^{-{\omega}^{2}\bar{T}^{2}/2}-\sqrt{\frac{\pi}{2}}{\omega}\,\bar{T}\,\mathcal{E}\text{rfc}\left(\frac{{\omega}\,\bar{T}}{\sqrt{2}}\right)\bigg{]}~{}.$
(29)

Here one can notice that with the use of the identity
$\mathcal{E}\text{rfc}(z)=\Gamma(1/2,z^{2})/\sqrt{\pi}$, this expression
matches exactly with the one provided in Sriramkumar and Padmanabhan (1996).
Thus we have obtained that the integrals $\mathcal{I}$ are non-vanishing with
finite switching, and they depend on the form of the window function. This
expression will be used in the expression Eq. (7a) along with the non-local
terms $\mathcal{I}_{\varepsilon}$ to investigate the entanglement measure
negativity. It should be noted that in the infinite interaction time limit,
i.e., when $\bar{T}\to\infty$, the entire $\mathcal{I}$ from Eq. (29) vanishes
(as $\bar{T}\mathcal{E}\text{rfc}\left({\omega}\,\bar{T}/\sqrt{2}\right)\sim
e^{-{\omega}^{2}\bar{T}^{2}/2}$ when $\bar{T}\to\infty$), providing the usual
eternal interaction result for the inertial detectors.

Evaluation of the non-local term $\mathcal{I}_{\varepsilon}$:

Now let us focus on evaluating the non-local entangling term
$\mathcal{I}_{\varepsilon}$ using the Gaussian switching function
$\chi(\tau)=e^{-\tau^{2}/\bar{T}^{2}}$. In this regard, the Hadamard form of
the Feynman propagator from Eq. (14) with the expression of the geodesic
distance (II.1), is used in the integral of Eq. (4b). We should mention that
here providing a final analytical expression does not seem possible. However,
we have tried to obtain as much analytical results as possible and then
numerically evaluate the final form of the integral
$\mathcal{I}_{\varepsilon}$. For the simplicity of calculation, we consider
two identical detectors, i.e., ${\omega}^{1}={\omega}^{2}={\omega}$, and then
the concerned integral can be represented as

$\displaystyle\mathcal{I}_{\varepsilon}=\frac{1}{4\pi^{2}}\int_{-\infty}^{\infty}\frac{d\tau_{2}^{\prime}}{\gamma}\int_{-\infty}^{\infty}d\tau\frac{e^{i{\omega}\left(\tau+\tau_{2}^{\prime}/\bar{\gamma}\right)}}{\left[\sigma_{\mathcal{G}}^{2}(\tau_{2}^{\prime}-\tau,\tau_{2}^{\prime})+i\epsilon\right]}e^{-\tau^{2}/\bar{T}^{2}}e^{-{\tau_{2}^{\prime}}^{2}/(\bar{\gamma}^{2}\bar{T}^{2})}~{}.$
(30)

We take the geodesic distance between the detectors
$\sigma_{\mathcal{G}}^{2}(\tau_{2}^{\prime}-\tau,\tau_{2}^{\prime})$ from Eq.
(II.1), which was expressed in a form,
$\sigma_{\mathcal{G}}^{2}=-(\tau_{2}^{\prime}-\tau)^{2}\mathpzc{f}(\tau_{2}^{\prime})+\mathpzc{g}(\tau_{2}^{\prime})$.
We consider a change of variables $\tau_{2}^{\prime}-\tau=\bar{u}$ and
$\tau_{2}^{\prime}=\bar{\tau}$, similar to the one done in the previous case
for evaluating $\mathcal{I}_{j}$. To perform the integration over $\bar{u}$ we
express the factor $e^{-(\bar{u}-\bar{\tau})^{2}/\bar{T}^{2}}$ in the integral
in terms of Fourier transform, i.e.,
$e^{-(\bar{u}-\bar{\tau})^{2}/\bar{T}^{2}}=\bar{T}/(2\sqrt{\pi})\int_{-\infty}^{\infty}d\kappa\,e^{i\kappa(\bar{u}-\bar{\tau})-\kappa^{2}\bar{T}^{2}/4}$.
The resulting integral looks like

$\displaystyle\mathcal{I}_{\varepsilon}$ $\displaystyle=$
$\displaystyle\frac{\bar{T}}{8\pi^{2}\sqrt{\pi}}\int_{-\infty}^{\infty}\frac{d\bar{\tau}\,e^{-\bar{\tau}^{2}/(\bar{\gamma}^{2}\bar{T}^{2}})}{\gamma}\int_{-\infty}^{\infty}d\kappa\,e^{-\kappa^{2}\bar{T}^{2}/4}e^{i\left\\{{\omega}\left(1+1/\bar{\gamma}\right)-\kappa\right\\}\bar{\tau}}\,\int_{-\infty}^{\infty}d\bar{u}\,\frac{e^{i(\kappa-{\omega})\bar{u}}}{\mathpzc{f}(\bar{\tau})\,\bar{u}^{2}-\mathpzc{g}(\bar{\tau})-i\epsilon}~{}.$
(31)

This integral has poles at
$\bar{u}=\pm\sqrt{\mathpzc{g}(\bar{\tau})+i\epsilon}/\sqrt{\mathpzc{f}(\bar{\tau})}$
of order one, which is respectively in the upper and lower half complex plane
depending on the ‘+’ or ‘-’ sign. When the $\kappa>{\omega}$, one should
consider a contour in the upper half complex plane and the pole
$\bar{u}=\sqrt{\mathpzc{g}(\bar{\tau})+i\epsilon}/\sqrt{\mathpzc{f}(\bar{\tau})}$
would contribute to the non-vanishing integral. Let us denote this pole as
$\sqrt{\mathpzc{g}(\bar{\tau})+i\epsilon}/\sqrt{\mathpzc{f}(\bar{\tau})}=\bar{u}_{0}(\bar{\tau})$.
Whereas, when $\xi<{\omega}$, one should consider a contour in the lower half
complex plane and the pole $\bar{u}=-\bar{u}_{0}(\bar{\tau})$ contributes to
the non-zero integral. Therefore we will have the expression of the previous
integral as

$\displaystyle\mathcal{I}_{\varepsilon}$ $\displaystyle=$
$\displaystyle\frac{2\pi\,i\,\bar{T}}{8\pi^{2}\sqrt{\pi}}\int_{-\infty}^{\infty}\frac{d\bar{\tau}\,e^{-\bar{\tau}^{2}/(\bar{\gamma}^{2}\bar{T}^{2}})}{\gamma}e^{i{\omega}\left(1+1/\bar{\gamma}\right)\bar{\tau}}\Bigg{[}\int_{{\omega}}^{\infty}d\kappa\,e^{-\kappa^{2}\bar{T}^{2}/4}e^{-i\,\kappa\bar{\tau}}\,\frac{e^{i(\kappa-{\omega})\bar{u}_{0}(\bar{\tau})}}{2\,\mathpzc{f}(\bar{\tau})\,\bar{u}_{0}(\bar{\tau})}~{}$
$\displaystyle~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}+\,\int_{-\infty}^{{\omega}}d\kappa\,e^{-\kappa^{2}\bar{T}^{2}/4}e^{-i\,\kappa\bar{\tau}}\,\frac{e^{-i(\kappa-{\omega})\bar{u}_{0}(\bar{\tau})}}{2\,\mathpzc{f}(\bar{\tau})\,\bar{u}_{0}(\bar{\tau})}\Bigg{]}\,.$
(32)

After carrying out these integrations over $\kappa$, the entire non-local term
$\mathcal{I}_{\varepsilon}$ looks like

$\displaystyle\mathcal{I}_{\varepsilon}$ $\displaystyle=$
$\displaystyle\frac{i\,\bar{T}}{4\pi\sqrt{\pi}}\int_{-\infty}^{\infty}\frac{d\bar{\tau}\,e^{-\bar{\tau}^{2}/(\bar{\gamma}^{2}\bar{T}^{2}})}{\gamma}\frac{e^{i{\omega}\left(1+1/\bar{\gamma}\right)\bar{\tau}}}{\bar{T}\,\mathpzc{f}(\bar{\tau})\,\bar{u}_{0}(\bar{\tau})}\Bigg{[}i\pi^{3/2}\exp\left(-\frac{\bar{u}_{0}(\bar{\tau})\left(2\bar{\tau}-i\bar{T}^{2}{\omega}\right)+\bar{\tau}^{2}+\bar{u}_{0}^{2}(\bar{\tau})}{\bar{T}^{2}}\right)$
$\displaystyle\times$
$\displaystyle~{}\left\\{\mathcal{E}\text{rf}\left(\frac{2i\bar{u}_{0}(\bar{\tau})+2i\bar{\tau}+\bar{T}^{2}{\omega}}{2\bar{T}}\right)+\mathcal{E}\text{rfc}\left(\frac{-2i\bar{u}_{0}(\bar{\tau})+2i\bar{\tau}+\bar{T}^{2}{\omega}}{2\bar{T}}\right)\exp\left(\frac{2\,\bar{u}_{0}(\bar{\tau})\left(2\bar{\tau}-i\bar{T}^{2}{\omega}\right)}{\bar{T}^{2}}\right)+1\right\\}\Bigg{]}\,.$
(33)

We could not find an analytical result for this integral. However, one can
always do a numerical analysis, and here also, we have numerically evaluated
this integral.

_Initially static detectors $(v_{0}=0)$_ :

Let us first consider the case when the initial coordinate velocity of the
detector $B$ in the frame of detector $A$ is zero. In flat spacetime, this
corresponds to two static detectors separated by a distance, $d_{0}$. Since
the detectors are static, the separation distance will not change for a flat
spacetime configuration. But when these detectors are placed in curved
spacetime, the curvature of the spacetime will produce deviation between the
initially static detectors. For $v_{0}=0$, the only curvature tensor component
in the expression for geodesic distance (II.1) will be $R_{\xi 0\xi 0}$.

We define $R_{\xi 0\xi 0}:=k_{1}$ and to make the approximations valid, we are
choosing the parameter range in the order, $\bar{T}\sim
d_{0},\,\bar{T}\ll\ell_{R}$ where $\ell_{R}$ is the curvature length scale
corresponding to the tidal curvature. Scaling all the length scales associated
with the problem by $\bar{T}$ leads to the domain of validity of our
calculations for $d_{0}/\bar{T}\sim 1,\,\ell_{R}/\bar{T}\gg 1$. Furthermore,
for $v_{0}=0$ we have the expressions of the auxiliary functions from Eq. (10)
as $\mathpzc{f}_{0}(\tau_{2}^{\prime})=1+k_{1}d_{0}^{2}/3$, and
$\mathpzc{g}_{0}(\tau_{2}^{\prime})=d_{0}^{2}(1-k_{1}\,{\tau_{2}^{\prime}}^{2})$.
Whereas we have from Eqs. (49) and (13) the expressions
$\gamma=(1+k_{1}\,d_{0}^{2})^{-1/2}$ and
$\bar{\gamma}=(1+k_{1}\,d_{0}^{2}\,\tau_{2}^{\prime}/2)^{-1}$.

For initially static detectors, $\mathcal{I}_{\varepsilon}$ depends on the
energy gap of the detectors, initial separation, and curvature length scale.
In Fig. 9 we have plotted the local terms $\mathcal{I}$ and the non-local
entangling term $|\mathcal{I}_{\varepsilon}|$ as a function of the
dimensionless detector transition energy $({\omega}\,\bar{T})$ in the case of
$v_{0}=0$. One should note that the local terms $\mathcal{I}$ from Eq. (29)
remain the same in both flat and curved spacetimes and even in the scenario of
$v_{0}\neq 0$. The mentioned figure asserts that both the local and the non-
local terms decrease with increasing detector transition energy. We also
observe that the non-local term is $\sim 10^{2}$ times larger than the local
terms. Therefore, the characteristics of negativity, which is expressed as
$\mathcal{N}^{(2)}=|\mathcal{I}_{\varepsilon}|-\mathcal{I}$, will also be
dictated by the non-local entangling term.

On the other hand, in Fig. 10, we have plotted the difference in negativity
$\Delta\mathcal{N}^{(2)}$ between the curved and flat spacetimes as functions
of the dimensionless detector transition energy $({\omega}\,\bar{T})$, the
initial separation $d_{0}/\bar{T}$, and the curvature length scale
$\ell_{R}/\bar{T}$. This difference decreases as the curvature length scale
increases, i.e., as the curvature itself decreases. The difference
$\Delta\mathcal{N}^{(2)}$ also decreases with increasing initial separation
$d_{0}/\bar{T}$ between the detectors.

_General geodesic detectors $(v_{0}\neq 0)$_ :

Let us now talk about the scenario when the initial coordinate velocity of the
detector $B$ in the frame of detector $A$ is not zero, i.e., $v_{0}\neq 0$. In
that scenario, there are six distinct curvature terms that affect the non-
local entangling term $\mathcal{I}_{\varepsilon}$. For instance, let us
consider $R_{\xi 0\xi 0}=k_{1}$, $R_{\xi 0n0}=k_{2}$, $R_{n0n0}=k_{3}$,
$R_{\xi 0\xi n}=k_{4}$, $R_{\xi n\xi n}=k_{5}$, and $R_{\xi nn0}=k_{6}$. Then
the quantities $\mathpzc{f}(\tau_{2}^{\prime})$ and
$\mathpzc{g}(\tau_{2}^{\prime})$ will be given by

$\displaystyle\mathpzc{f}(\tau_{2}^{\prime})$ $\displaystyle=$ $\displaystyle
1+\frac{1}{3}\left[k_{1}d_{0}^{2}+2d_{0}\left(k_{1}v_{||}+k_{2}v_{\perp}\right)\tau_{2}^{\prime}+\left(k_{1}v_{||}^{2}+2k_{2}v_{||}v_{\perp}+k_{3}v_{\perp}^{2}\right){\tau_{2}^{\prime}}^{2}\right]~{},$
(34a) $\displaystyle\mathpzc{g}(\tau_{2}^{\prime})$ $\displaystyle=$
$\displaystyle\left(d_{0}^{2}+v_{0}^{2}{\tau_{2}^{\prime}}^{2}+2d_{0}v_{||}\tau_{2}^{\prime}\right)+d_{0}^{2}\,\left(\frac{2}{3}k_{4}v_{\perp}-k_{1}+\frac{1}{3}k_{5}v_{\perp}^{2}\right){\tau_{2}^{\prime}}^{2}$
(34b)
$\displaystyle-\frac{4}{3}d_{0}\left(k_{1}v_{||}+k_{2}v_{\perp}\right){\tau_{2}^{\prime}}^{3}-\frac{1}{3}\left(k_{1}v_{||}^{2}+2k_{2}v_{||}v_{\perp}+k_{3}v_{\perp}^{2}\right){\tau_{2}^{\prime}}^{4}~{}.$

While the auxiliary functions
$\mathpzc{h}_{1}(v_{0},d_{0}),\,\mathpzc{h}_{2}(v_{0},d_{0})$ and
$\mathpzc{h}_{3}(v_{0},d_{0})$ necessary to define $\gamma$ and $\bar{\gamma}$
will be

$\displaystyle\mathpzc{h}_{1}(v_{0},d_{0})$ $\displaystyle:=$ $\displaystyle
k_{1}d_{0}^{2}+\frac{2}{3}k_{4}d_{0}^{2}v_{\perp}$ (35a)
$\displaystyle\mathpzc{h}_{2}(v_{0},d_{0})$ $\displaystyle:=$ $\displaystyle
k_{1}v_{||}+k_{2}v_{\perp}+\frac{2}{3}v_{\perp}\left(k_{4}v_{||}+k_{6}v_{\perp}\right)$
(35b) $\displaystyle\mathpzc{h}_{3}(v_{0},d_{0})$ $\displaystyle:=$
$\displaystyle k_{1}v_{||}^{2}+2k_{2}v_{||}v_{\perp}+k_{3}v_{\perp}^{2}~{}.$
(35c)

We use these expressions to evaluate the entangling term
$\mathcal{I}_{\varepsilon}$. The behavior of this term for inertial probes in
flat spacetime with the Gaussian switching function is similar to that
reported in Ref. Koga et al. (2018). This is expected as we have considered
the curvature terms to be perturbatively small. In Fig. 12, we have plotted
the difference in negativity $\Delta\mathcal{N}^{(2)}$ between the curved and
flat spacetime in this case with respect to different system parameters.

##  IV Characteristics of $\mathcal{N}^{(2)}$ in de Sitter spacetime

In this section, we will describe several observations about negativity in the
de Sitter background, utilizing the expression of Sec. III.1. In Sec. III.1,
we have estimated the necessary quantities to investigate the entanglement in
a de Sitter background. We will now use numerical plots to extract key
features from the resultant expressions.

Individual detector response:

In Fig. 4, we have plotted the local terms $\mathcal{I}$ that denote
individual detector transition probabilities. One can observe that these
probabilities decrease with increasing energy gap and with increasing
curvature length scale. However, for any curvature scale $\Lambda\,\bar{T}\ll
1$, it remains almost the same and quantitatively similar to the ones obtained
in Minkowski background. In general curved spacetimes, the characteristics of
the local terms are expected to be qualitatively similar to the de Sitter ones
for large curvature length scales (for geodesic detectors, these will come
entirely from the van Vleck determinant).

Figure 4: The local term $\mathcal{I}$ in entanglement, which denotes
individual detector transition probability, is plotted as functions of the
detector transition energy and curvature for static detector in de Sitter
background. The transition probability depends on the curvature constant. As
the curvature increases, the transition probability increases. However, this
increment is prominent only in the $\Lambda\,\bar{T}^{2}>1$ regime, as
observed from the right plot.

Entanglement between the detectors:

In Fig. 5 we have plotted the negativity in a de Sitter background as
functions of the dimensionless detector transition energy
$({\omega}\,\bar{T})$. The left plot here depicts curves with different values
of the initial separation, and the right plot gives curves for a fixed initial
separation, but in different curved backgrounds. One can observe that, for a
fixed energy gap, as the separation increases, at first the negativity
decreases (from black $\to$ green $\to$ blue curve), and then increases (blue
$\to$ red curve). Therefore, the behaviour of negativity is not monotonous
with respect to the initial separation. A similar feature also holds when one
looks at negativity as a function of curvature, with fixed initial separation
(the plot on the right). From the right plot, we observe that in the low
detector transition energy regime, the negativity decreases monotonically as
the background curvature increases. However, there exists a band of energy gap
(roughly $1\lesssim\omega{\bar{T}}\lesssim 3$ in the given plot), over which
negativity is no longer a monotonous function of curvature. We should
highlight that, our results are consistent with those in Kukita and Nambu
(2017), when restricted to the parameter range used in that paper.

Figure 5: In both the above plots we have plotted the negativity
$\mathcal{N}^{(2)}$ as a function of the dimensionless detector energy gap
${\omega}\,\bar{T}$ in a de Sitter background. In particular, on the left, we
have depicted negativity $\mathcal{N}^{(2)}$ for different initial separations
between the detectors. On the right, we have plotted the negativity in the de
Sitter background for different curvature length scales. We should mention
that the local terms $\mathcal{I}$, along with the non-local terms
$\mathcal{I}_{\varepsilon}$, are also affected by high curvature, which
results in significant distortion in the negativity profiles. The right plot
suggests that one can observe entanglement due to curved backgrounds using
detectors with a sufficiently large energy gap.

Imprints of curvature on negativity:

We can combine the plots and discussion from previous subsections to explore
the possibility of using entanglement as a probe of background curvature,
which we now discuss in the context of de Sitter before turning to general
curved spacetime (below). One of the main objectives behind this analysis is
to provide qualitative, and to some extent, quantitative, similarities between
de Sitter and general curved backgrounds with small curvature. Since our de
Sitter results do not assume small $\Lambda$, we can also use them to validate
our perturbative technique.

In Fig. 6 we have plotted the difference in negativity between the de Sitter
and Minkowski backgrounds. Unlike in Fig. 5, here we have chosen small values
for curvature which will facilitate comparison with the perturbative results
for general curved backgrounds in the later section.

Figure 6: The difference in the negativity between the de Sitter and flat
spacetimes $\Delta\mathcal{N}^{(2)}$ is plotted as a function of the
dimensionless detector transition energy ${\omega}\,\bar{T}$. In both plots,
we have considered the curvature length scale to be very large, i.e., the
curvature itself is very small so that one can compare these results with the
ones from our perturbative approach. Comparing these plots with the
perturbative approach of generic curved spacetime, we observe that the plots
from the two cases are qualitatively and quantitatively similar.

Fig. 7 shows the contour plot for negativity components that particularly
result from spacelike separation $(\mathcal{N}^{+})$, causal communication
$(\mathcal{N}^{-})$, and the total negativity $(\mathcal{N}^{(2)})$ in a de
Sitter background. These different negativity components with their
significance were mentioned in Sec. II.3. These contour diagrams are plotted
as functions of the dimensionless initial separation and detector energy gap.
The two different rows of the figure correspond to two different curvature
length scales. One can observe that the negativity due to spacelike
separation, often associated with genuine entanglement effects Martin-Martinez
(2015); Tjoa and Martín-Martínez (2021); Barman et al. (2022a), is positive
only in the regime of large energy gaps and initial separations. At the same
time, the negativity due to the communication channel is positive in the
entire region except for the very large initial. These phenomena are true for
both of the considered curvature length scales.

The physical reason behind the above observations is straightforward: for a
finite switching time scale, as the initial separation between the detectors
increases, the second detector moves away from the accessible causal lightcone
patch of the first detector, see Fig. 3 for illustration. Nevertheless, one
observes that there exists a range of parameters that will yield the so-called
genuine entanglement effects.

Figure 7: The contribution in negativity due to spacelike separation $N^{+}$
(left plot), purely due to the causal communication $N^{-}$ (middle plot), and
the entire negativity $N^{(2)}$ (right plot) are plotted for $\Lambda=10^{-4}$
(top row) and $\Lambda=0.1$ (bottom row), see Eq. (20) for the definitions of
$N^{\pm}$. These plots correspond to detectors in the de Sitter background
with initial separation $d_{0}$. As is evident from the plots, entanglement
without causal communication is possible only in the high detector transition
energy $\omega\,\bar{T}$ regime. Whereas, entanglement from the causal
communication is prominent in the low separation $d_{0}/\bar{T}$ and low
transition energy regime ${\omega}\,\bar{T}$.

Finally, in Fig. 8, we have plotted the spacelike negativity component
$(\mathcal{N}^{+})$, the causal communication-based negativity component
$(\mathcal{N}^{-})$, and their difference as functions of the initial detector
separation. For these plots, we have considered a suitable parameter space,
taking inspiration from Fig. 7, a parameter space where, for large detector
separation, the spacelike negativity can become larger compared to the causal
negativity. From these plots, one can observe that for small initial
separation $d_{0}$, the communication-based negativity $(\mathcal{N}^{-})$ is
always greater compared to the spacelike negativity. However, for large
initial detector separation, $\mathcal{N}^{+}$ is becoming larger compared to
$\mathcal{N}^{-}$. From the difference between these two quantities we observe
that, for small background curvature, the transition from
$\mathcal{N}^{+}<\mathcal{N}^{-}$ to $\mathcal{N}^{+}>\mathcal{N}^{-}$ happens
near $d_{0}/\bar{T}\simeq 1$, and this value keeps on decreasing as the
curvature increases. The physical reasoning could be as the curvature
increases, the past and the future-directed light cones depart more from each
other with increasing separation, see Fig. 3. Therefore, as the curvature
increases, the scope of causal communication for a certain detector separation
is also diminished for finite detector switching.

Figure 8:  The contribution in the non-local entangling term
$\mathcal{I}_{\varepsilon}$ with and without causal communication, which are
denoted as $\mathcal{I}_{\varepsilon}^{-}$ and
$\mathcal{I}_{\varepsilon}^{+}$, are plotted in the top row of this figure as
functions of the dimensionless distance $d_{0}/\bar{T}$ for different
curvature scales of the background. We observe that, as the distance
decreases, the $|\mathcal{I}_{\varepsilon}^{-}|$ part keeps increasing, while
$|\mathcal{I}_{\varepsilon}^{+}|$ seems to reach saturated values. We have
also plotted the modulus of the entire entangling term
$\mathcal{I}_{\varepsilon}$ that mostly resembles
$|\mathcal{I}_{\varepsilon}^{-}|$. In the below two rows of the figure, we
have plotted the difference $(\mathcal{N}^{+}-\mathcal{N}^{-})$ as functions
of the dimensionless initial separation $d_{0}/\bar{T}$ between the detectors.
Each of the plots from these two rows corresponds to a different background
curvature. From these plots, we observe that as the separation increases,
there is a certain value of $d_{0}/\bar{T}$ when $\mathcal{N}^{+}$ overtakes
$\mathcal{N}^{-}$. We further observe that this crossing point in
$d_{0}/\bar{T}$ decreases with increasing curvature. In the plots from the
lower two rows, the dashed lines correspond to the respective values of
$1/(\sqrt{\Lambda}\,\bar{T})$, which indicate the curvature length scales.

##  V Entanglement as a probe of curvature

In this part, we shall discuss our observations on entanglement in a general
curved spacetime, for which we have used the perturbative technique as
described in Sec. III.2.

Figure 9: The local term $\mathcal{I}$ and the non-local term
$|\mathcal{I}_{\varepsilon}|$ (in $v_{0}=0$ scenario) in the negativity are
plotted as functions of the dimensionless detector transition energy
$\bar{T}\,{\omega}$ for the Gaussian switching function
$\chi(\tau)=e^{-\tau^{2}/\bar{T}^{2}}$. The local terms $\mathcal{I}$ signify
individual detector transition probabilities and are identical in the flat and
the curved backgrounds, at least for the geodesic detectors. From the right
figure, we observe that the non-local term $\mathcal{I}_{\varepsilon}$ is
indistinguishable in the curved background compared to the flat spacetime.
However, there is a perceivable difference between them that will be clear in
our subsequent analysis. The epilogue in the figure signifies the appropriate
value for the energy we have to choose, such that
$|\mathcal{I}_{\varepsilon}|>I$ and hence leads to entanglement between
detectors.

Initially static detectors, $v_{0}=0$:

Figure 10: The difference in the negativity
$\Delta\mathcal{N}^{(2)}=\big{[}\mathcal{N}^{(2)}\big{]}_{\text{Curv}}-\big{[}\mathcal{N}^{(2)}\big{]}_{\text{Flat}}$
is plotted as functions of the dimensionless detector transition energy
${\omega}\,\bar{T}$ in the two left figures. In the upper left figure,
different curves correspond to different initial separations $d_{0}/\bar{T}$,
and fixed curvature length scale $\ell_{R}/\bar{T}=50$. While in the lower
left plot, different curves correspond to different curvature length scales
$\ell_{R}/\bar{T}$, and fixed initial separation $d_{0}/\bar{T}=2$. In the
upper right figure, we have plotted $\Delta\mathcal{N}^{(2)}$ as a function of
the initial separation $d_{0}/\bar{T}$. In this plot, the curvature length
scale is fixed at $\ell_{R}/\bar{T}=50$. While in the lower right subfigure,
we have plotted $\Delta\mathcal{N}^{(2)}$ as a function of the dimensionless
curvature length scale, and the separation between the detectors is fixed at
$d_{0}/\bar{T}=2$. This plot signifies that the correction due to the
curvature becomes smaller as the curvature length scale increases, i.e., as
the curvature decreases. In all of these plots the initial velocity is
$v_{0}=0$.

Figure 11:  The difference in the negativity
$\Delta\mathcal{N}^{(2)}=\big{[}\mathcal{N}^{(2)}\big{]}_{\text{Curv}}-\big{[}\mathcal{N}^{(2)}\big{]}_{\text{Flat}}$
is plotted as functions of the dimensionless detector transition energy
${\omega}\,\bar{T}$ for negative curvature. One can observe that these plots
are flipped upside down compared to the ones from Fig. 10, i.e., the curved
background now reports more entanglement compared to the flat background,
unlike the previous case. However, certain characteristics remain the same,
like as the initial separation between the detectors increases and the
absolute value of curvature decreases, the difference in negativity
$\Delta\mathcal{N}^{(2)}$ from the flat background decreases.

In Fig. 9, we have plotted the local individual detector transition
probability denoting term $\mathcal{I}$ and the modulus of the non-local
entangling term $|\mathcal{I}_{\varepsilon}|$. As we are dealing with a nearly
flat spacetime background, the individual detector transition probabilities
$\mathcal{I}$ from Eq. (29) are the same as Minkowski. Previously, in de
Sitter background, we have seen that this quantity depends on the curvature,
but for a very small curvature, like the one considered here, it does not
change and remains the same as the one obtained in Fig. 9. From the right plot
of Fig. 9, it seems that $|\mathcal{I}_{\varepsilon}|$ from the curved and
Minkowski backgrounds overlap. However, there is a subtle difference in these
curves, which will be prominent in the low transition energy regime, see Fig.
10.

In Fig. 10, we have plotted the difference in negativity
$\Delta\mathcal{N}^{(2)}$ between the curved and flat backgrounds as a
function of the dimensionless detector energy gap for initial static
detectors. In the top row, we have depicted the variation of
$\Delta\mathcal{N}^{(2)}$ with respect to different initial detector
separations. Whereas, in the lower plots, the variation of
$\Delta\mathcal{N}^{(2)}$ is shown for different curvature. These lower plots
signify that as the curvature length scale increases, i.e., the curvature
decreases, the negativity profile from the curved background tends to resemble
the flat spacetime negativity, which is expected and consistent in our
analysis.

From the top row of Fig. 10, we observe that the negativity becomes more
Minkowski-like as the initial separation $d_{0}$ between the detectors
increases. To elaborate on this, for a smaller initial separation between the
detectors, the curvature effects in the entanglement become prominent. At
first glance, this outcome may seem a bit counter-intuitive, as some of the
curvature corrections in our perturbative analysis, see Eq. (10) and Eq. (50),
arrive with multiplicative $d_{0}$ terms. Thus curvature corrections should
increase with increasing $d_{0}$, which is not observed from our plots.
However, with a deeper investigation, one can perceive that this is not
unphysical. The reason can be found in the fact that the correlation between
the detectors degrades with increasing distance. If this degradation is more
than the increasing curvature effects, then the negativity profile from the
curved background tends to become Minkowski-like as the separation increases.
Further, for very small initial separations, the entanglement due to field-
mediated communication channels will dominate and if the detectors are placed
farther apart, the total negativity itself will decrease. Even though the
genuine entanglement starts to dominate after, $d_{0}\geq 2c\,\bar{T}$, the
net negativity is small.

With an increasing distance between inertial detectors, the decay of
correlation is reported in Koga et al. (2018) for eternal switching, i.e.,
$\chi(\tau)=1$. Moreover, when one considers a finite switching like the
Gaussian switching considered in our present work, this decay of correlation
with increasing detector separation becomes rapid. To support this argument,
one can look into the expressions of the non-local terms
$\mathcal{I}^{M}_{\varepsilon}({\omega})$ and
$\mathcal{I}^{C}_{\varepsilon}({\omega})$ from Appendix B (Eqs. 54 and B.0.1),
which respectively have $e^{-d_{0}^{2}/\bar{T}^{2}}$ and
$e^{-d_{0}^{2}/(2\bar{T}^{2})}$ terms multiplied in them. To further elucidate
this, we refer the reader to Fig. 3, where one can observe that as two
geodesic trajectories move away from each other, the detectors moving on these
curves interact less with each other with finite switching. Therefore, with
increasing initial separation between the detectors, one should get diminished
negativity and negativity like near flat spacetime.

Furthermore, in Fig. 11, we have plotted the difference in negativity
$\Delta\mathcal{N}^{(2)}$ for two initial static detectors for opposite sign
of curvature than Fig. 10. This plot shows that the entanglement from the
curved spacetime is now greater than the flat background. Therefore, if in the
previous case (Fig. 10) the curvature moved the geodesic trajectories apart,
in the present scenario (Fig. 11) it can be considered to get closer.
Moreover, we observed that the overall behavior of negativity to the variation
of initial separation and curvature remains the same. As the initial
separation between the detectors increases and the absolute value of curvature
decreases, the difference in negativity $\Delta\mathcal{N}^{(2)}$ between the
curved and flat backgrounds decreases.

General geodesic detectors, $v_{0}\neq 0$:

Figure 12: Top row: Behavior of $\mathcal{N}^{(2)}$ as functions of
${\omega}\,\bar{T}$, initial detector velocity $v_{0}/c$ and the angle
$\theta$ between the initial separation and the velocity vector. In the
extreme left plot, we have considered different initial velocities of the
second detector. The angle between the initial velocity and the initial
separation is taken as $\theta=\pi/4$, and the curvature length scale for all
the components of the Riemann tensor is $\ell_{R}/\bar{T}=50$, and the initial
separation is $d_{0}/\bar{T}=2$. In the right of both the top and bottom rows,
the initial velocity is taken as $v_{0}=0.01c$ for the second detector, and
the energy gap as $\omega\,\bar{T}=2.5$.
Bottom row: Behavior of
$\Delta\mathcal{N}^{(2)}=\big{[}\mathcal{N}^{(2)}\big{]}_{\text{Curv}}-\big{[}\mathcal{N}^{(2)}\big{]}_{\text{Flat}}$
as functions of ${\omega}\,\bar{T}$, initial detector velocity $v_{0}/c$ and
the angle $\theta$ between the initial separation and the velocity vector. We
have chosen the same parameter values as the top row. Figure 13: In both the
above plots, we have considered the initial configuration such that the
relative velocity of the detectors is perpendicular to their relative
separation; that is, $\theta=\pi/2$. This facilitates isolating a combination
of Riemann tensor components. When $v_{0}=0$, the only Riemann component that
contributes is $R_{\xi 0\xi 0}$. For $\theta=\pi/2$, if we choose $\bm{\xi}$
and $\bm{n}$ along the eigen-directions of the electric part of Riemann, the
difference from the initially static $(v_{0}=0)$ case arises solely due to
magnetic parts of the Riemann tensor. This is depicted in the left plot. In
the right plot, all the electric parts of the Riemann tensor are considered in
the same order of magnitude.

Figure 14: The difference in negativity between the curved and flat spacetimes
$\Delta\mathcal{N}^{(2)}$ is plotted as functions of different electric
components of the curvature tensor, i.e., with respect to $R_{\xi 0\xi 0}$,
$R_{n0n0}$, and $R_{n\xi n\xi}$. The initial separation is chosen as
$d_{0}/\bar{T}=2$ and the energy gap as $\omega\,\bar{T}=2.5$. In the upper
left plot, $R_{\xi 0\xi 0}$ is kept fixed, while in the upper right plot,
$R_{n0n0}$ is fixed. On the other hand, in the lower plot, $R_{n\xi n\xi}$ is
kept fixed. The upper two plots indicate no visible variation of
$\Delta\mathcal{N}^{(2)}$ with respect to the particular curvature term
$R_{n\xi n\xi}$. While there is a finite visible change in
$\Delta\mathcal{N}^{(2)}$ when $R_{\xi 0\xi 0}$ or $R_{n0n0}$ is varied. From
the lower plot, it is clear that the change in $\Delta\mathcal{N}^{(2)}$ due
to the change in $R_{\xi 0\xi 0}$ is greater compared to $R_{n0n0}$.

In Fig. 12, in the left column, we have plotted $\mathcal{N}^{(2)}$ and
$\Delta\mathcal{N}^{(2)}$ as functions of the dimensionless detector
transition energy for two detectors with different initial velocities
respectively. From Eq. (10) and Eq. (50), one can observe that now, with
velocity, many more different curvature terms will contribute to the
negativity profile. In the first two figures on the top row of Fig. 12,
negativity decreases with an increase in velocity, and the same is true if one
varies the angle between initial separation and the initial velocity of the
two detectors. The change in negativity, as given in the bottom row of Fig. 12
shows that $\Delta\mathcal{N}^{(2)}$ degrades as initial velocity increases.
The change in negativity decreases and then increases when the angle changes
from small angles to $\theta=\pi/2$. As the velocity increases, the distance
between detectors will increase rapidly and hence the correlation degrades. If
the angle is chosen very small, the detectors can come close and entanglement
due to communication can be enhanced. Hence choosing larger angles makes sure
there is no causal contact between the detectors.

In Figs. 10, 11, 12, 5, and 6 we have plotted either the difference in
negativity or the difference in the modulus of the non-local entangling term.
In these particular figures, we observed that as the dimensionless detector
energy exceeds certain values, in most cases as ${\omega}\,\bar{T}>1$, there
is a sudden change in behavior in the entanglement profiles. For instance, in
Fig. 10, we observe that for ${\omega}\,\bar{T}<1$, the curved space
entanglement is lower than the flat space, but for ${\omega}\,\bar{T}>1$ the
curved space entanglement becomes higher than the flat space entanglement. It
is not unusual to believe that there could be some interesting physics hidden
behind these phenomena.

Figure 13 portrays the effects due to the magnetic part of the Riemann tensor.
The choice of orientation can be used to isolate certain combinations of
Riemann tensor components in the expression for negativity. This can be
employed as a probe of curvature, especially those components which are
otherwise difficult to determine through other experiments.

In Fig. 14, the horizontal axes correspond to different curvature terms.
Interestingly in 14, we observed that the dependence of the negativity profile
on different curvature terms is not the same. For instance, with respect to
$R_{n\xi n\xi}$ the negativity profiles do not alter at all. While with
respect to $R_{n0n0}$ and $R_{\xi 0\xi 0}$ the negativity profiles change, but
again the change is larger when $R_{\xi 0\xi 0}$ varies compared to
$R_{n0n0}$. This may provide a means to distinguish between different curved
backgrounds in terms of entanglement.

##  VI Concluding remarks

Given the already exhaustive discussion and commentary in preceding sections,
we will keep the discussions in this concluding section brief, focusing
largely on some general observations extracted from the results we have
established.

Let us first summarize the key steps of our analysis: We have extracted the
role of curvature in the entanglement induced between two Unruh-DeWitt
detectors, interacting solely with a background scalar field in an arbitrary
background spacetime. Assuming the field state to be Hadamard, the key
curvature contribution comes from the relative kinematics of the trajectories.
While the curvature will also affect the dynamics of the field, this effect is
suppressed compared to that arising from the leading terms in the Hadamard
expansion, and this is nicely corroborated by the de Sitter case, where the
exact known expressions yield results for entanglement measure that closely
matches the one obtained from our perturbative scheme. The key geometric
reason for why this works, in terms of the causal structure of spacetime, was
already identified in Sec. II.

We now list some general observations based on this work. The list is only
representative, to aid in navigating through the bulk of the paper; we hope
the reader will be able to extract more insights from the various analytical
and graphical results presented here.

  1. 1.

Contribution to entanglement measure (negativity here) from spacelike
correlations, versus that due to timelike correlations, can be identified by
decomposing the Feynman propagator into its anti-commutator and commutator
parts. For smaller distances, contribution of commutator dominates, while at
larger separations, the anti-commutator contribution dominates, which is what
one would want to attribute genuine quantum entanglement with. Hence, choosing
a suitable $d_{0}/\bar{T}$ can lead to genuine entanglement (in the suitable
range of energy gap), and this ratio is altered by presence of Riemann tensor.

  2. 2.

For lower energy gaps, the local term corresponding to the individual detector
transition due to vacuum polarisation will dominate over the non-local term
for initial separation, $d_{0}\geq\bar{T}$. As the energy gap becomes large
enough, i.e., $\omega\geq 2\hbar/\bar{T}$, the non-local term starts to
dominate over the individual detector transition probabilities.

  3. 3.

For parameters chosen such that $d_{0}\approx 2\bar{T}$ and $\omega\approx
2\hbar/\bar{T}$, the genuine entanglement between detectors is present.
Presence of curvature enhances the entanglement for the above parameter range.

  4. 4.

When the detectors have an initial relative velocity ${\bm{v}_{0}}$, the
entanglement depends on this relative velocity and its angle $\theta$ with
respective to the initial separation $\bm{\xi}$. The genuine entanglement
decreases with an increase in speed. Overall, the entanglement also increases
with increase in $\theta$.

The key upshot of the above dependences, and the interplay between them, is
that these may allow one to measure spacetime curvature using quantum probes
with suitable initial conditions.

  5. 5.

Our investigation on the role of curvature in entanglement between the probes
indicates that the entanglement is most sensitive to the eigenvalue of tidal
part of the Riemann tensor – $R_{a0c0}=R_{abcd}u^{b}u^{d}$ – along the
direction of separation of the detectors $\bm{\xi}$; that is,
$R_{a0b0}\xi^{a}\xi^{b}$. This contribution is present even for initially
static detectors. Entanglement is also sensitive to $R_{a0b0}n^{a}n^{b}$, the
eigenvalue of the tidal tensor in the direction $n^{a}$ of the relative
initial velocity of the detectors. This dependence, however, is weaker than
the previous one. Interestingly, the entanglement does not depend on the
components such as $R_{abcd}n^{a}\xi^{b}n^{c}\xi^{d}$, which denotes the
purely spatial part of the curvature.

The initial relative configuration of the detectors can be used as a probe to
determine the curvature of spacetime. An explicit example of one such set-up
is demonstrated in Fig. 13, which can be used to extract the contribution of
the magnetic part of Riemann tensor. It would be interesting to compare the
viability of probing such components of Riemann using entanglement, with
classical experiments.

  6. 6.

There exists in literature a couple of works Cliche and Kempf (2011); Kukita
and Nambu (2017) which attempt to study entanglement between detectors in
specific backgrounds with non-zero curvature. While helpful, a major drawback
of working with specific solutions of Einstein equations is that it hides the
explicit role of curvature and hence does not yield physical insights which
can be claimed to be generic. Our motivation, set-up, as well overall analysis
in this work are much more general and hence have a much wider scope. We can,
however, compare our results with the ones in these works in cases when the
choice of initial conditions and parameters coincide. In Ref. Kukita and Nambu
(2017), the authors study entanglement in a de Sitter background and conclude
that entanglement at the super horizon scale is possible only with multiple
quantum probes. Our present work finds that the parameter values for
entanglement extraction are in the same range as those obtained in Ref. Kukita
and Nambu (2017) for the de Sitter background. On the other hand, in Ref.
Cliche and Kempf (2011), entanglement generation is studied outside of weakly
gravitating, spherically symmetric, bodies such as the Earth and the Sun.
Their result indicating enhancement of entanglement (compared to Minkowski)
for large energy gaps is also observed in our case.

  7. 7.

Our work should be relevant in observationally viable scenarios such as for
tests of the violation of Bell’s inequality. In this context, we would like to
specifically mention approaches to understand entanglement between two non-
intersecting spatial regions in de Sitter spacetime, which do not suggest the
presence of entanglement Espinosa-Portalés and Vennin (2022). It will be
interesting to relate our analysis with these and identify the points of
conflicts and similarities between the same. It should not be difficult to
construct a set-up similar to Martin and Vennin (2021a); Espinosa-Portalés and
Vennin (2022); Martin and Vennin (2021b), using the tools introduced here, to
discern entanglement induced in a general curved background. This
investigation is currently under progress.

  8. 8.

Finally, we would like to mention a broader context in which this work will
serve as the first step. This concerns a known subtlety in application of the
equivalence principle to deduce the role of gravity in quantum phenomenon by
analyzing these phenomenon in non-inertial frames in Minkowski spacetime.
While doing so is usually easy since one is working in Minkowski spacetime, it
has its pitfalls, and can serve as a hurdle to obtaining insights which a
perturbative analysis can never give. For instance, in Ref. Hari and Kothawala
(2021, 2023), it is shown that even a single accelerated probe can pick up
non-trivial effects due to tidal curvature which cannot be seen in a
perturbative expansion. Our analysis in this paper, though still perturbative
in curvature, develops the necessary mathematical tools that, along with some
of the exact results we obtain in de Sitter, have the potential to yield
stronger results beyond a perturbative expansion. This is under investigation.

###### Acknowledgements.

We would like to thank Vincent Vennin, Alessandro Pesci, and L. Sriramkumar
for useful discussions and comments on the manuscript. H.K. thanks the Indian
Institute of Technology Madras and the Ministry of Human Resources and
Development (MHRD), India for financial support. S.B. would like to thank the
Science and Engineering Research Board (SERB), Government of India (GoI), for
supporting this work through the National Post Doctoral Fellowship (N-PDF,
File number: PDF/2022/000428).

##  Appendix A Geodesic interval between events on two timelike trajectories

We start with the following expansion (ref. Fig. (1)):

$\Omega_{\mathcal{G}}=[\Omega_{\mathcal{G}}]+\Omega_{\mathcal{G}_{s}}[t^{i}\nabla_{i}\Omega_{\mathcal{G}}]+\frac{1}{2}\Omega^{2}_{B^{\prime}B}\left[t^{i}\nabla_{i}\left(t^{k}\nabla_{k}\Omega_{\mathcal{G}}\right)\right]+...$
(36)

Here, $\Omega_{xx^{\prime}}:=\Omega(x,x^{\prime})$ represents Synge’s world
function which is half the square of the geodesic distance between $x$ and
$x^{\prime}$ and $[\,]$ represents the limit,
$\Omega_{\mathcal{G}s}\rightarrow 0$. The FNC at $B^{\prime}$ is defined by,
$x^{i^{\prime}}=\left(\tau_{2}^{\prime},X^{\alpha^{\prime}}(\tau=\tau_{2}^{\prime})\right)$
and
$X^{\alpha^{\prime}}(\tau=\tau_{2}^{\prime}):=\sigma_{\mathcal{G}s}(\tau_{2}^{\prime})t^{i^{\prime}}(\tau_{2}^{\prime}){\rm
e}^{\alpha^{\prime}}_{\phantom{\alpha^{\prime}}i^{\prime}}(\tau_{2}^{\prime})$
where,
$\sigma_{\mathcal{G}s}(\tau_{2}^{\prime})t^{i^{\prime}}=-\sigma^{i^{\prime}}_{\mathcal{G}s}(\tau_{2}^{\prime}):=-g^{i^{\prime}k^{\prime}}\nabla_{k^{\prime}}\Omega_{\mathcal{G}s}(\tau_{2}^{\prime})$.
Using this definition, the expansion for $\Omega_{\mathcal{G}}$ can be
expressed in terms of the Fermi coordinates.

$\Omega_{\mathcal{G}}=\Omega_{\mathcal{G}^{\prime}}+[\nabla_{\alpha^{\prime}}\Omega_{\mathcal{G}}]X^{\alpha^{\prime}}(\tau_{2}^{\prime})+\frac{1}{2}\left[\nabla_{\alpha^{\prime}}\nabla_{\beta^{\prime}}\Omega_{\mathcal{G}}\right]X^{\alpha^{\prime}}(\tau_{2}^{\prime})X^{\beta^{\prime}}(\tau_{2}^{\prime})+...$
(37)

Here, for a reference geodesic curve,
$\Omega_{\mathcal{G}^{\prime}}=-\sigma^{2}_{\mathcal{G}^{\prime}}/2\equiv-1/2\left(\tau_{2}^{\prime}-\tau\right)^{2}$.
Now, the expansions of $[\,]$ quantities and
$X^{\alpha^{\prime}}(\tau_{2}^{\prime})$ should be expressed in terms of data
given at $\tau_{1}=0$.

The $[\,]$ quantities need to be evaluated carefully. The derivatives of
$\Omega_{\mathcal{G}}$ are expanded in double Taylor series of
$\Omega_{\mathcal{G}^{\prime}}$ with $\tau$ and $\tau_{2}^{\prime}$ as
variables.

$\displaystyle\nabla_{\alpha^{\prime}}\Omega_{\mathcal{G}}$
$\displaystyle=\nabla_{\alpha^{\prime}}\Omega_{\mathcal{G}^{\prime}}+\tau
u^{i}\nabla_{i}\left(\nabla_{\alpha^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right)+\tau_{2}^{\prime}u^{i^{\prime}}\nabla_{i^{\prime}}\left(\nabla_{\alpha^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right)+...$
(38a)
$\displaystyle\nabla_{\alpha^{\prime}}\nabla_{\beta^{\prime}}\Omega_{\mathcal{G}}$
$\displaystyle=\nabla_{\alpha^{\prime}}\nabla_{\beta^{\prime}}\Omega_{\mathcal{G}^{\prime}}+\tau
u^{i}\nabla_{i}\left(\nabla_{\alpha^{\prime}}\nabla_{\beta^{\prime}}\Omega_{\mathcal{G}}\right)+\tau_{2}^{\prime}u^{i^{\prime}}\nabla_{i^{\prime}}\left(\nabla_{\alpha^{\prime}}\nabla_{\beta^{\prime}}\Omega_{\mathcal{G}}\right)+...$
(38b)

Similarly, all higher derivatives are expressed and then the coincidence limit
is correctly taken. The calculations are cumbersome and hence we have used
CADABRAPeeters (2007a, b) for the calculations.

Note that, $\nabla_{\alpha^{\prime}}\Omega_{\mathcal{G}^{\prime}}\equiv{\rm
e}^{i^{\prime}}_{\phantom{i^{\prime}}{\alpha^{\prime}}}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}$
and note that $u^{i}\nabla_{i}{\rm
e}^{i^{\prime}}_{\phantom{i^{\prime}}\alpha^{\prime}}=0$ since the tetrad is
at $x^{\prime}$ and derivative is taken at $x$. Similarly,
$u^{i^{\prime}}\nabla_{i^{\prime}}{\rm
e}^{i^{\prime}}_{\phantom{i^{\prime}}\alpha^{\prime}}=0$ since the tetrad is
Fermi-Walker transported along a geodesic and will not generally be zero for
other curves. Substituting these observations into Eq. (38b), the RHS is
expressed in terms of covariant derivatives of $\Omega_{\mathcal{G}^{\prime}}$
at $x^{\prime}$ and $x$. Then the coincidence limit,
$\tau_{2}^{\prime},\tau\rightarrow 0$ is taken. Some covariant derivatives are
as follows:

$\displaystyle\left[\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]$
$\displaystyle=0$
$\displaystyle\left[\nabla_{j}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]$
$\displaystyle=-g_{ij}$
$\displaystyle\left[\nabla_{j^{\prime}}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]$
$\displaystyle=g_{ij}$
$\displaystyle\left[\nabla_{k^{\prime}}\nabla_{j^{\prime}}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]=\left[\nabla_{k}\nabla_{j^{\prime}}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]$
$\displaystyle=\left[\nabla_{k^{\prime}}\nabla_{j}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]=\left[\nabla_{k}\nabla_{j}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]=0$
$\displaystyle\left[\nabla_{l^{\prime}}\nabla_{k^{\prime}}\nabla_{j^{\prime}}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]=\left[\nabla_{l}\nabla_{k}\nabla_{j}\nabla_{i}\Omega_{\mathcal{G}^{\prime}}\right]$
$\displaystyle=-\left[\nabla_{l}\nabla_{k^{\prime}}\nabla_{j^{\prime}}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]=\left[\nabla_{l}\nabla_{k}\nabla_{j^{\prime}}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]=-\frac{1}{3}\left(R_{ikjl}+R_{iljk}\right)$
$\displaystyle\left[\nabla_{l}\nabla_{k}\nabla_{j}\nabla_{i^{\prime}}\Omega_{\mathcal{G}^{\prime}}\right]$
$\displaystyle=-\frac{1}{3}\left(R_{ijkl}+R_{ikjl}\right)$
$\displaystyle.....$

Substituting the coincidence limits to the covariant expansion of
$\Omega_{\mathcal{G}}$ in Eq. (37) results in the following equation,

$\displaystyle\Omega_{\mathcal{G}}$
$\displaystyle=-\frac{1}{2}\left(\tau_{2}^{\prime}-\tau\right)^{2}+\frac{1}{2}\left(g_{\alpha\beta}-\frac{1}{3}\left(\tau_{2}^{\prime}-\tau\right)^{2}R_{ijkl}{\rm
e}^{i}_{\phantom{i}\alpha}{\rm
e}^{k}_{\phantom{k}\beta}u^{j}u^{l}\right)X^{\alpha}\left(\tau_{2}^{\prime}\right)X^{\beta}\left(\tau_{2}^{\prime}\right)+\mathcal{O}(\tau^{4},X^{5},R^{2})$
(39)

Notice that, we have assumed that the derivatives of the Riemann tensor are
negligible compared to the Riemann tensor itself for the above expression. The
FNC needs to be expanded in terms of $\tau_{2}^{\prime}$. For FNC, the
expansions are simple Taylor expansion in terms of $\tau_{2}^{\prime}$ given
by,

$X^{\alpha^{\prime}}(\tau_{2}^{\prime})=X^{\alpha^{\prime}}(0)+\left.\frac{dX^{\alpha}}{d\tau}\right|_{O}\tau_{2}^{\prime}+\frac{1}{2}\left.\frac{d^{2}X^{\alpha}}{d\tau^{2}}\right|_{O}t^{2}_{B^{\prime}}+...$
(40)

The coefficients of the Taylor expansion can be evaluated using the definition
of deviation vector in terms of Synge’s world function as given by Ref.Vines
(2015). Let the deviation vector be, $\xi^{i}=X^{\alpha}{\rm
e}^{i}_{\phantom{i}\alpha}$

$\displaystyle X^{\alpha}$ $\displaystyle=-\sigma^{i}_{\mathcal{G}s}(\tau){\rm
e}^{\alpha}_{\phantom{\alpha}i}(\tau)$ (41)
$\displaystyle\frac{\mathrm{d}X^{\alpha}}{\mathrm{d}\tau}$
$\displaystyle=-{\rm
e}^{\alpha}_{\phantom{\alpha}i}(\tau)u^{k}\nabla_{k}\sigma^{i}_{\mathcal{G}s}(\tau)-{\rm
e}^{\alpha}_{\phantom{\alpha}i}(\tau)v^{\tilde{k}}\nabla_{\tilde{k}}\sigma^{i}_{\mathcal{G}s}(s)$
(42) $\displaystyle...$

Here,
$v^{\tilde{k}}:=\left({\mathrm{d}s}/{\mathrm{d}\tau}\right)u^{\tilde{k}}_{B}$,
where $u^{\tilde{k}}_{B}$ is the 4-velocity of the detector $B$ along its
trajectory. Using the expansion for derivatives of Synge world function
givesVines (2015),

$\displaystyle\nabla_{k}\sigma^{i}_{\mathcal{G}s}(\tau)$
$\displaystyle=\delta^{i}_{\phantom{i}k}-\frac{1}{3}R^{i}_{\phantom{i}\alpha
k\beta}X^{\alpha}X^{\beta}+\mathcal{O}(R^{2},\nabla R,X^{3})$ (43)
$\displaystyle\nabla_{\tilde{k}}\sigma^{i}_{\mathcal{G}s}(s)$
$\displaystyle=-g^{j}_{\phantom{j}\tilde{k}}\left(\delta^{i}_{\phantom{i}j}+\frac{1}{6}R^{i}_{\phantom{i}\alpha
j\beta}X^{\alpha}X^{\beta}+\mathcal{O}(R^{2},\nabla R,X^{3})\right)$ (44)

Here, $g^{j}_{\phantom{j}\tilde{k}}$ is the parallel propagator. The
coefficients for the Taylor expansion of FNC can be evaluated using Eq. (A)
and Eq. (44) and can be obtained as,

$\displaystyle\frac{\mathrm{d}X^{\alpha}}{\mathrm{d}\tau}$
$\displaystyle=v^{\alpha}_{0}+\frac{1}{3}R^{\alpha}_{\phantom{\alpha}\mu
0\nu}X^{\mu}X^{\nu}+\frac{1}{6}R^{\alpha}_{\phantom{\alpha}\mu
k\nu}X^{\mu}X^{\nu}v^{k}_{0}+\mathcal{O}(R^{2},\nabla R,X^{3})$
$\displaystyle\frac{\mathrm{d}^{2}X^{\alpha}}{\mathrm{d}\tau^{2}}$
$\displaystyle=-R^{\alpha}_{\phantom{\alpha}0\beta
0}X^{\beta}+\mathcal{O}(R^{2},\nabla R,X^{2})$
$\displaystyle\frac{\mathrm{d}^{3}X^{\alpha}}{\mathrm{d}\tau^{3}}$
$\displaystyle\approx-R^{\alpha}_{\phantom{\alpha}0\beta
0}v_{0}^{\beta}+\mathcal{O}(R^{2},\nabla R,X^{2})$

$X^{\alpha}(\tau)=\left(X^{\alpha}+v^{\alpha}_{0}\tau_{2}^{\prime}+\frac{1}{3}R^{\alpha}_{\phantom{\alpha}\mu
0\nu}X^{\mu}X^{\nu}\tau_{2}^{\prime}+\frac{1}{6}R^{\alpha}_{\phantom{\alpha}\mu
k\nu}X^{\mu}X^{\nu}v^{k}_{0}\tau_{2}^{\prime}-\frac{1}{2}R^{\alpha}_{\phantom{\alpha}0\mu
0}X^{\mu}{\tau_{2}^{\prime}}^{2}-\frac{1}{6}R^{\alpha}_{\phantom{\alpha}0\beta
0}v_{0}^{\beta}{\tau_{2}^{\prime}}^{3}\right)+\mathcal{O}(R^{2},\nabla R)$
(45)

Using the above simplifications in the Eq. (39), the geodesic distance to
$R^{2},\nabla R$ order will be,

$\displaystyle-\sigma^{2}_{\mathcal{G}}=\Delta\tau_{\mathcal{G}}^{2}=$
$\displaystyle\left(\tau_{2}^{\prime}-\tau\right)^{2}-\left[\left(X^{\alpha}(0)+v^{\alpha}_{0}\tau_{2}^{\prime}\right)^{2}+\frac{2}{3}R_{\alpha\mu
0\nu}v_{0}^{\alpha}X^{\mu}X^{\nu}{\tau_{2}^{\prime}}^{2}-R_{\alpha 0\beta
0}X^{\alpha}(0)X^{\beta}(0){\tau_{2}^{\prime}}^{2}-\frac{4}{3}R_{\alpha 0\beta
0}X^{\alpha}(0)v_{0}^{\beta}{\tau_{2}^{\prime}}^{3}\right.$
$\displaystyle\quad\quad\quad\quad\quad\quad\quad\left.+\frac{1}{3}R_{\alpha\mu\beta\nu}X^{\alpha}X^{\beta}v_{0}^{\mu}v_{0}^{\nu}{\tau_{2}^{\prime}}^{2}-\frac{1}{3}R_{\mu
0\nu 0}v_{0}^{\mu}v_{0}^{\nu}{\tau_{2}^{\prime}}^{4}\right]$
$\displaystyle+\frac{1}{3}\left(\tau_{2}^{\prime}-\tau\right)^{2}R_{\alpha
0\beta
0}\left[X^{\alpha}(0)X^{\beta}(0)+2X^{\alpha}(0)v^{\beta}_{0}\tau_{2}^{\prime}+v_{0}^{\alpha}v_{0}^{\beta}{\tau_{2}^{\prime}}^{2}\right]+\mathcal{O}(\nabla
R,R^{2})$ (46)

The initial separation is ${X}^{\alpha}(0):=d_{0}\,\xi^{\alpha}$ where,
$\xi^{\alpha}$ is the directional cosines and
$\eta_{\alpha\beta}\xi^{\alpha}\xi^{\beta}=1$. The initial velocity can be
decomposed into components along the initial separation and the component
along the normal to the initial separation denoted by unit vector $n^{\alpha}$
that is in the plane orthogonal to $X^{\alpha}(0)$,
$v_{0}^{\alpha}=v_{||}\,\xi^{\alpha}+v_{\perp}\,n^{\alpha}$. Note that
$v_{||}=v_{0}\cos\theta$ and $v_{\perp}=v_{0}\sin\theta$. Using this
decomposition in the above expression, we finally obtain Eq. (II.1).

The mapping between proper times of the detectors: Apart from geodesic
interval, another quantity we will need to characterize entanglement is the
mapping between proper times of the two detectors. The relevant quantitiy here
is $\gamma=\left({\mathrm{d}s}/{\mathrm{d}\tau}\right)$, which can be easily
found using the Fermi metric based on $\mathcal{C}_{1}$, assumed to bbe
geodesic (see, for e.g., Poisson (2004)).

$\displaystyle\mathrm{d}S^{2}$ $\displaystyle=$
$\displaystyle-~{}\mathrm{d}{\tau_{2}^{\prime}}^{2}\left(1+R_{0\alpha
0\beta}X^{\alpha}X^{\beta}\right)-2\mathrm{d}X^{\mu}\mathrm{d}\tau_{2}^{\prime}\left(\frac{2}{3}R_{0\alpha\mu\beta}X^{\alpha}X^{\beta}\right)$
$\displaystyle~{}~{}~{}+~{}\mathrm{d}X^{\mu}\mathrm{d}X^{\nu}\left(\delta_{\mu\nu}-\frac{1}{3}R_{\mu\alpha\nu\beta}X^{\alpha}X^{\beta}\right)+\mathcal{O}(X^{3},\nabla
R,R^{2})~{}.$ (47)

Since, $\mathrm{d}S^{2}=-\mathrm{d}s^{2}$ the ratio is,

$\displaystyle\gamma(\tau^{\prime}_{2}):=\frac{\mathrm{d}\tau_{2}^{\prime}}{\mathrm{d}s}=\left[1-\left(\frac{dX^{\alpha}}{d\tau_{2}^{\prime}}\right)^{2}+R_{0\alpha
0\beta}X^{\alpha}X^{\beta}+\frac{4}{3}R_{0\alpha\mu\beta}X^{\alpha}X^{\beta}\frac{dX^{\mu}}{d\tau_{2}^{\prime}}+\frac{1}{3}R_{\mu\alpha\nu\beta}X^{\alpha}X^{\beta}\frac{dX^{\mu}}{d\tau_{2}^{\prime}}\frac{dX^{\nu}}{d\tau_{2}^{\prime}}\right]^{-1/2}~{}.$
(48)

Using Eq. (45) and resolving the components as done for the geodesic distance,
in the equation for $\gamma$, to linear order in curvature the equation
becomes,

$\displaystyle\gamma(\tau^{\prime}_{2})=\frac{\mathrm{d}\tau_{2}^{\prime}}{\mathrm{d}s}$
$\displaystyle=$
$\displaystyle\left[1-v^{2}_{0}+\mathpzc{h}_{1}(v_{0},d_{0})+2d_{0}\mathpzc{h}_{2}(v_{0},d_{0})\tau_{2}^{\prime}+\mathpzc{h}_{3}(v_{0},d_{0}){\tau_{2}^{\prime}}^{2}\right]^{-1/2}+~{}\mathcal{O}(R^{2},\nabla
R)~{}.$ (49)

Here, $\mathpzc{h}_{1}(v_{0},d_{0}),\,\mathpzc{h}_{2}(v_{0},d_{0})$ and
$\mathpzc{h}_{3}(v_{0},d_{0})$ are auxiliary functions defined by,

$\displaystyle\mathpzc{h}_{1}(v_{0},d_{0})$ $\displaystyle:=$ $\displaystyle
R_{\xi 0\xi 0}d_{0}^{2}+\frac{2}{3}R_{\xi 0\xi n}d_{0}^{2}v_{\perp}$ (50a)
$\displaystyle\mathpzc{h}_{2}(v_{0},d_{0})$ $\displaystyle:=$ $\displaystyle
R_{\xi 0\xi 0}v_{||}+R_{\xi 0n0}v_{\perp}+\frac{2}{3}v_{\perp}\left(R_{\xi
0\xi n}v_{||}+R_{n\xi 0n}v_{\perp}\right)$ (50b)
$\displaystyle\mathpzc{h}_{3}(v_{0},d_{0})$ $\displaystyle:=$ $\displaystyle
R_{\xi 0\xi 0}v_{||}^{2}+2R_{\xi 0n0}v_{||}v_{\perp}+R_{n0n0}v_{\perp}^{2}$
(50c)

The proper time of $B$ can be expressed as, $s=\int
d\tau_{2}^{\prime}/\gamma(\tau_{2}^{\prime})$. (When the coordinate velocity
of second detector is zero, $\gamma=\sqrt{1+R_{\xi 0\xi 0}d_{0}^{2}}$.)

##  Appendix B Evaluation of $\mathcal{I}_{\varepsilon}$ in general curved
spacetime: complete perturbative analysis

We observed that the individual detector transition probabilities
$\mathcal{I}$ and non-local entangling terms $\mathcal{I}_{\varepsilon}$
depend on the geodesic distance that in a general curved spacetime was
estimated in a perturbative manner in (II.1). We have already analytically
estimated the individual local terms $\mathcal{I}$ for small curvatures in
III.2. Now we shall evaluate the non-local entangling term
$\mathcal{I}_{\varepsilon}$ by expressing it perturbatively around the
Minkowski scenario, i.e., by assuming the corrections due to the curvature to
be very small compared to other system parameters specifically the initial
separation between the detectors. Let us proceed this way and express
$\mathcal{I}_{\varepsilon}$ for two initially static and general geodesic
detectors.

####  B.0.1 Evaluation of $\mathcal{I}_{\varepsilon}$ for two initially static
detectors

Let us first consider the case of two static detectors. In this scenario, the
only quantity that depends on curvature terms and which is present in the
geodesic distance between the two detectors is $R_{\alpha 0\beta
0}X^{\alpha}(0)X^{\beta}(0):=k_{1}d_{0}^{2}$. Let us term this quantity
$c_{1}$, i.e., $c_{1}=k_{1}d_{0}^{2}$. We take the geodesic distance between
the detectors
$\sigma_{\mathcal{G}}^{2}(\tau^{\prime}_{2}-\tau,\tau^{\prime}_{2})$ from Eq.
(II.1). We consider a change of variables $\tau_{2}^{\prime}-\tau=\Delta\tau$
and $\tau_{2}^{\prime}=v$, similar to the one done in the case for evaluating
$\mathcal{I}_{j}$. First, the Feynman propagator in the non-local entangling
term can be expanded considering $c_{1}\ll 1$. Second, the non-local term
contains other factors involving the switching functions and the exponential
with the detector transition energy. These factors will also contain the
curvature term $c_{1}$ due to their representation in terms of Fermi-normal
coordinates with respect to a certain geodesic, see Eq. (13). Finally, the
non-local entangling term can be expressed as

$\mathcal{I}_{\varepsilon}=\mathcal{I}^{M}_{\varepsilon}+\mathcal{I}^{C}_{\varepsilon}~{},$
(51)

where with the Gaussian switching function we have

$\displaystyle\mathcal{I}^{M}_{\varepsilon}({\omega})$ $\displaystyle=$
$\displaystyle-\frac{1}{4\pi^{2}}\int_{-\infty}^{\infty}dv\,\int_{-\infty}^{\infty}d\Delta\tau\,\,\frac{e^{-\frac{(v-\Delta\tau)^{2}+v^{2}}{\bar{T}^{2}}+i{\omega}(2v-\Delta\tau)}}{\left(\Delta\tau^{2}-d_{0}^{2}-i\epsilon\right)}~{},$
(52a) $\displaystyle\mathcal{I}^{C}_{\varepsilon}({\omega})$ $\displaystyle=$
$\displaystyle-\frac{1}{4\pi^{2}}\int_{-\infty}^{\infty}dv\,\int_{-\infty}^{\infty}d\Delta\tau\,\,\frac{e^{-\frac{(v-\Delta\tau)^{2}+v^{2}}{\bar{T}^{2}}+i{\omega}(2v-\Delta\tau)}}{\left(\Delta\tau^{2}-d_{0}^{2}-i\epsilon\right)}\Bigg{[}\frac{(\Delta\tau^{2}/3)-v^{2}}{\left(\Delta\tau^{2}-d_{0}^{2}-i\epsilon\right)}-\left\\{\frac{i\,v\,{\omega}+1}{2}-\frac{v^{2}}{\bar{T}^{2}}\right\\}\Bigg{]}\,c_{1}\,$
$\displaystyle~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}~{}\hskip
284.52756pt+~{}\mathcal{O}(c_{1}^{2})~{}.$ (52b)

Let us first evaluate the Minkowski part
$\mathcal{I}^{M}_{\varepsilon}({\omega})$. Here the integration over $\tau$
can be readily carried out as the denominators of $G^{M,C}_{F}(\Delta\tau,v)$
are independent of $v$. Using the Gaussian integration formula
$\int_{-\infty}^{\infty}dv\,e^{-\alpha(v-\beta)^{2}}=\sqrt{\pi/\alpha}$ one
can express

$\displaystyle\mathcal{I}^{M}_{\varepsilon}({\omega})=-\frac{\bar{T}e^{-\frac{{\omega}^{2}\bar{T}^{2}}{2}}}{4\sqrt{2}\pi^{3/2}}\int_{-\infty}^{\infty}d\Delta\tau\,\frac{e^{-\frac{\Delta\tau^{2}}{2\bar{T}^{2}}}}{\left(\Delta\tau^{2}-d_{0}^{2}-i\epsilon\right)}~{}.$
(53)

To evaluate this integral, we express the involved Gaussian factor as
$e^{-\frac{\Delta\tau^{2}}{2\bar{T}^{2}}}=(\bar{T}/\sqrt{2\pi})\int_{-\infty}^{\infty}d\kappa\,e^{-(\zeta^{2}\bar{T}^{2}/2)+i\Delta\tau\,\kappa}$,
i.e., in terms of its Fourier transform. The integral of the previous equation
has poles at $\Delta\tau=\pm\sqrt{d_{0}^{2}+i\epsilon}$, which are in the
upper and the lower half complex planes respectively, depending on the ‘+ve’
and the ‘-ve’ signs. Therefore, we must consider the contour in the upper or
lower semicircle to dampen the necessary parts of the contour integral. After
finding the residues, the integration limit of $\kappa$ will also be
restricted. For example, for the upper pole, we have $\kappa\in[0,\infty)$,
and for the lower pole, we have $\kappa\in(-\infty,0]$. Then, carrying out the
integration over $\kappa$, we have the final expression

$\displaystyle\mathcal{I}^{M}_{\varepsilon}({\omega})=-\frac{i\,\bar{T}\left(\mathcal{E}\text{rf}\left(\frac{i\,d_{0}}{\sqrt{2}\,\bar{T}}\right)+1\right)}{4\,\sqrt{2\pi}\,d_{0}}\,e^{-\frac{d_{0}^{2}}{2\bar{T}^{2}}-\frac{{\omega}^{2}\bar{T}^{2}}{2}}\,.$
(54)

Let us now evaluate the integral $\mathcal{I}^{C}_{\varepsilon}({\omega})$
that corresponds to the contribution entirely from the curved spacetime due to
the terms up to $\mathcal{O}(c_{1})$. It is to be noted that several
numerators in the term $\mathcal{I}^{C}_{\varepsilon}({\omega})$ also depend
on $\tau$. If we integrate over this variable $\tau$, we shall get

$\displaystyle\mathcal{I}^{C}_{\varepsilon}({\omega})$ $\displaystyle=$
$\displaystyle\frac{c_{1}e^{-\frac{{\omega}^{2}\bar{T}^{2}}{2}}}{48\sqrt{2}\pi^{3/2}\bar{T}}\int_{-\infty}^{\infty}\frac{d\Delta\tau\,e^{-\frac{\Delta\tau^{2}}{2\bar{T}^{2}}}}{\left(\Delta\tau^{2}-d_{0}^{2}-i\epsilon\right)^{2}}\,\Bigg{[}\bigg{\\{}3\Delta\tau^{2}\left(\Delta\tau^{2}-i\epsilon\right)-3d_{0}^{2}\left(\Delta\tau^{2}+\bar{T}^{2}(-1+i{\omega}\Delta\tau)\right)-3{\omega}^{2}\bar{T}^{6}~{}$
$\displaystyle~{}~{}\hskip
113.81102pt+~{}\bar{T}^{4}(3+6i{\omega}\Delta\tau)+\bar{T}^{2}\left(\Delta\tau^{2}(4+3i{\omega}\Delta\tau)+3\epsilon({\omega}\Delta\tau+i)\right)\bigg{\\}}\Bigg{]}~{}.$
(55)

Here also the poles are at the same places
$\Delta\tau=\pm\sqrt{d_{0}^{2}+i\epsilon}$ but with order two, and we proceed
in a similar fashion to integrate over the variable $\Delta\tau$. The final
expression for $\mathcal{I}^{C}_{\varepsilon}({\omega})$ will be given by

$\displaystyle\mathcal{I}^{C}_{\varepsilon}({\omega})$ $\displaystyle=$
$\displaystyle\frac{c_{1}}{192\,\pi\,d_{0}^{3}\,\bar{T}}\,e^{-\frac{{\omega}^{2}\bar{T}^{2}}{2}-\frac{d_{0}^{2}}{2\bar{T}^{2}}}\bigg{[}\sqrt{2\pi}\,\mathcal{E}\text{rf}\text{i}\left(\frac{d_{0}}{\sqrt{2}\bar{T}}\right)\left\\{d_{0}^{4}+d_{0}^{2}\left(2\bar{T}^{2}-3{\omega}\,^{2}\bar{T}^{4}\right)-3{\omega}\,^{2}\bar{T}^{6}+3\bar{T}^{4}\right\\}-i\,\sqrt{2\pi}\,d_{0}^{4}$
$\displaystyle+~{}i\sqrt{2\pi}\,d_{0}^{2}\,\bar{T}^{2}\left(3{\omega}\,^{2}\bar{T}^{2}-2\right)+6d_{0}\bar{T}^{3}e^{\frac{d_{0}^{2}}{2\bar{T}^{2}}}\left({\omega}\,^{2}\bar{T}^{2}-1\right)-14d_{0}^{3}\bar{T}e^{\frac{d_{0}^{2}}{2\bar{T}^{2}}}+3i\sqrt{2\pi}\bar{T}^{4}\left({\omega}\,^{2}\bar{T}^{2}-1\right)\bigg{]}.$
(56)

Therefore, we have the analytical expression for both the
$\mathcal{I}^{M}_{\varepsilon}({\omega})$ and
$\mathcal{I}^{C}_{\varepsilon}({\omega})$, which will lead one to the entire
expression of $\mathcal{I}_{\varepsilon}({\omega})$ from Eq. (51). We have
verified that the plots of this non-local entangling term are similar in
nature compared to the ones from Figs. 9, 10, and 11. It further validates our
previous procedure.

####  B.0.2 Evaluation of $\mathcal{I}_{\varepsilon}$ for two geodesic
detectors

For two detectors in curved spacetime with an initial separation and relative
velocity between them, one can express the non-local entangling term like the
previous case as
$\mathcal{I}_{\varepsilon}=\mathcal{I}^{M}_{\varepsilon}+\mathcal{I}^{C}_{\varepsilon}$.
Here $\mathcal{I}^{M}_{\varepsilon}$ corresponds to the contribution from the
Minkowski spacetime and is given by

$\displaystyle\mathcal{I}^{M}_{\varepsilon}=-\frac{1}{4\pi^{2}\gamma_{0}}\int_{-\infty}^{\infty}dv\,\int_{-\infty}^{\infty}d\Delta\tau\,\frac{e^{-\frac{(v-\Delta\tau)^{2}}{\bar{T}^{2}}-\frac{v^{2}}{\bar{T}^{2}\gamma_{0}^{2}}+i{\omega}(v(1+1/\gamma_{0})-\Delta\tau)}}{\left\\{\Delta\tau^{2}-(d_{0}^{2}+2d_{0}\,v_{0}\,\cos{(\theta)}\,v+v_{0}^{2}v^{2})-i\epsilon\right\\}}~{}.$
(57)

Whereas $\mathcal{I}^{C}_{\varepsilon}$ corresponds to the term solely
dependent on the presence of curvature. This quantity is given by

$\displaystyle\mathcal{I}^{C}_{\varepsilon}$ $\displaystyle=$
$\displaystyle\frac{c_{1}}{4\pi^{2}\gamma_{0}}\int_{-\infty}^{\infty}dv\,\int_{-\infty}^{\infty}d\Delta\tau\,\frac{e^{-\frac{(v-\Delta\tau)^{2}}{\bar{T}^{2}}-\frac{v^{2}}{\bar{T}^{2}\gamma_{0}^{2}}+i{\omega}(v(1+1/\gamma_{0})-\Delta\tau)}}{\left\\{\Delta\tau^{2}-(d_{0}^{2}+2d_{0}\,v_{0}\,\cos{(\theta)}\,v+v_{0}^{2}v^{2})-i\epsilon\right\\}}\,$
$\displaystyle\times$
$\displaystyle\Bigg{[}-\frac{\bar{\mathpzc{f}}(v)\,\Delta\tau^{2}+\bar{\mathpzc{g}}(v)}{\left\\{\Delta\tau^{2}-(d_{0}^{2}+2d_{0}\,v_{0}\,\cos{(\theta)}\,v+v_{0}^{2}v^{2})-i\epsilon\right\\}}+\bigg{\\{}\frac{\mathcal{H}_{1}(v)\,\gamma_{0}^{2}+iv\,\gamma_{0}{\omega}\,\mathcal{H}_{2}(v)}{2}-\frac{v^{2}\,\mathcal{H}_{2}(v)}{\bar{T}^{2}}\bigg{\\}}\Bigg{]}\,.$
(58)

Here $\bar{\mathpzc{f}}(v)=\mathpzc{f}(v)-1$, and
$\bar{\mathpzc{g}}(v)=\mathpzc{g}(v)-\left(d_{0}^{2}+v_{0}^{2}v^{2}+2d_{0}v_{||}v\right)$
can be obtained from Eq. (10). Whereas,
$\mathcal{H}_{1}(v)=\left\\{\mathpzc{h}_{1}(v_{0},d_{0})+2d_{0}\mathpzc{h}_{2}(v_{0},d_{0})+\mathpzc{h}_{3}(v_{0},d_{0})v^{2}\right\\}$,
and
$\mathcal{H}_{2}(v)=\left\\{\mathpzc{h}_{1}(v_{0},d_{0})+d_{0}\mathpzc{h}_{2}(v_{0},d_{0})+\mathpzc{h}_{3}(v_{0},d_{0})v^{2}/3\right\\}$
are obtained utilizing the expressions from Eq. (50). Therefore, one can
easily understand that the first part inside the brackets in Eq. (B.0.2)
appears from the curvature corrections to the geodesic interval. On the other
hand, the second part inside the brackets in Eq. (B.0.2) appears from the
transformation from $s$ to $\tau_{2}^{\prime}$. The above integrals can be
easily performed using the Fourier transform
$e^{-\frac{(\Delta\tau-v)^{2}}{\bar{T}^{2}}}=(\bar{T}/\sqrt{2\pi})\int_{-\infty}^{\infty}d\kappa\,e^{i\,\kappa(\Delta\tau-v)-\kappa^{2}\bar{T}^{2}/4}$.
One can first perform the integral over $\Delta\tau$ through the help of
contour integration and then over $\kappa$. Both of these integrations can be
performed analytically. Finally, one is left with the integration over $v$,
for which we could not find an analytical result. We have sought the help of
numerical integration to perform this final integration. We have compared the
results obtained from here with the ones from our previous calculation from
Fig. 12, and confirmed similar outcomes.

##  Appendix C 3D plots - dependence of negativity on velocity and angle

Figure 15: The negativity as a function of both initial velocity and angle
between initial separation and initial velocity is shown in the left. The
difference in negativity between the curved and flat spacetimes
$\Delta\mathcal{N}^{(2)}$ is plotted as functions of the initial velocity $v$
of the detector $B$ and the angle $\theta$ between the initial velocity and
the separation on the right. In both of the plots, the curvature length scales
are of the order of $\ell_{R}=50\bar{T}$, the dimensionless detector energy
gap is chosen as, $\omega\,\bar{T}=2.5$ and the initial separation distance is
$d_{0}/T=2$.

## References

  * Pikovski et al. (2015) I. Pikovski, M. Zych, F. Costa, and Č. Brukner, Nature Physics 11, 668 (2015).
  * 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, Phys. Rev. Lett. 119, 240401 (2017).
  * Marletto and Vedral (2017) C. Marletto and V. Vedral, Phys. Rev. Lett. 119, 240402 (2017).
  * Kothawala (2023a) D. Kothawala (2023a), eprint arXiv:2305.09265.
  * Kothawala (2023b) D. Kothawala, in _Spacetime, Matter, Quantum Mechanics_ (2023b), eprint arXiv:2304.01995.
  * Koga et al. (2018) J.-i. Koga, G. Kimura, and K. Maeda, Phys. Rev. A 97, 062338 (2018), eprint arXiv:1804.01183.
  * Reznik (2003) B. Reznik, Found. Phys. 33, 167 (2003), eprint arXiv:quant-ph/0212044.
  * Koga et al. (2019) J.-i. Koga, K. Maeda, and G. Kimura, Phys. Rev. D 100, 065013 (2019), eprint arXiv:1906.02843.
  * Martin-Martinez (2015) E. Martin-Martinez, Phys. Rev. D 92, 104019 (2015), eprint arXiv:1509.07864.
  * Tjoa and Martín-Martínez (2021) E. Tjoa and E. Martín-Martínez, Phys. Rev. D 104, 125005 (2021), eprint arXiv:2109.11561.
  * Cong et al. (2020) W. Cong, C. Qian, M. R. R. Good, and R. B. Mann, JHEP 10, 067 (2020), eprint arXiv:2006.01720.
  * Zhang and Yu (2020) J. Zhang and H. Yu, Phys. Rev. D 102, 065013 (2020), eprint arXiv:2008.07980.
  * Barman et al. (2022a) S. Barman, D. Barman, and B. R. Majhi, JHEP 09, 106 (2022a), eprint arXiv:2112.01308.
  * Barman et al. (2022b) D. Barman, S. Barman, and B. R. Majhi, Phys. Rev. D 106, 045005 (2022b), eprint arXiv:2205.08505.
  * Hu and Yu (2015) J. Hu and H. Yu, Phys. Rev. A 91, 012327 (2015), eprint arXiv:1501.03321.
  * Menezes (2018) G. Menezes, Phys. Rev. D97, 085021 (2018), eprint arXiv:1712.07151.
  * Brown (2013) E. G. Brown, Phys. Rev. A 88, 062336 (2013), eprint arXiv:1309.1425.
  * Simidzija and Martín-Martínez (2018) P. Simidzija and E. Martín-Martínez, Phys. Rev. D 98, 085007 (2018), eprint arXiv:1809.05547.
  * Barman et al. (2021) D. Barman, S. Barman, and B. R. Majhi, JHEP 07, 124 (2021), eprint arXiv:2104.11269.
  * Pozas-Kerstjens and Martin-Martinez (2015) A. Pozas-Kerstjens and E. Martin-Martinez, Phys. Rev. D 92, 064042 (2015), eprint arXiv:1506.03081.
  * Xu et al. (2020) Q. Xu, S. A. Ahmad, and A. R. H. Smith, Phys. Rev. D 102, 065019 (2020), eprint arXiv:2006.11301.
  * Gray et al. (2021) F. Gray, D. Kubiznak, T. May, S. Timmerman, and E. Tjoa, JHEP 11, 054 (2021), eprint arXiv:2105.09337.
  * Barman et al. (2023) S. Barman, I. Chakraborty, and S. Mukherjee, JHEP 09, 180 (2023), eprint arXiv:2305.17735.
  * Ostapchuk et al. (2012) D. C. M. Ostapchuk, S.-Y. Lin, R. B. Mann, and B. L. Hu, JHEP 07, 072 (2012), eprint arXiv:1108.3377.
  * Hu et al. (2012) B. L. Hu, S.-Y. Lin, and J. Louko, Class. Quant. Grav. 29, 224005 (2012), eprint arXiv:1205.1328.
  * Menezes and Svaiter (2015) G. Menezes and N. Svaiter, Phys. Rev. A 92, 062131 (2015), eprint arXiv:1508.04513.
  * Menezes and Svaiter (2016) G. Menezes and N. Svaiter, Phys. Rev. A 93, 052117 (2016), eprint arXiv:1512.02886.
  * Arias et al. (2016) E. Arias, J. Dueñas, G. Menezes, and N. Svaiter, JHEP 07, 147 (2016), eprint arXiv:1510.00047.
  * Menezes (2016) G. Menezes, Phys. Rev. D94, 105008 (2016), eprint arXiv:1512.03636.
  * Liu et al. (2018) X. Liu, Z. Tian, J. Wang, and J. Jing, Phys. Rev. D 97, 105030 (2018), eprint arXiv:1805.04470.
  * Picanço et al. (2020) G. Picanço, N. F. Svaiter, and C. A. Zarro, JHEP 08, 025 (2020), eprint arXiv:2002.06085.
  * Barman and Majhi (2021) S. Barman and B. R. Majhi, JHEP 03, 245 (2021), eprint arXiv:2101.08186.
  * Barman et al. (2022c) S. Barman, B. R. Majhi, and L. Sriramkumar (2022c), eprint arXiv:2205.01305.
  * Cliche and Kempf (2011) M. Cliche and A. Kempf, Phys. Rev. D 83, 045019 (2011), eprint arXiv:1008.4926.
  * Kukita and Nambu (2017) S. Kukita and Y. Nambu, Entropy 19, 449 (2017), eprint arXiv:1708.01359.
  * Tjoa and Mann (2020) E. Tjoa and R. B. Mann, JHEP 08, 155 (2020), eprint arXiv:2007.02955.
  * Unruh (1976) W. Unruh, Phys.Rev. D14, 870 (1976).
  * Peres (1996) A. Peres, Phys. Rev. Lett. 77, 1413 (1996), eprint arXiv:quant-ph/9604005.
  * Hadamard (2003) J. Hadamard, _Lectures on Cauchy’s Problem in Linear Partial Differential Equations_ , Dover phoenix editions (Dover Publications, 2003), ISBN 9780486495491.
  * Decanini and Folacci (2006) Y. Decanini and A. Folacci, Phys. Rev. D 73, 044027 (2006), eprint arXiv:gr-qc/0511115.
  * Synge (1960) J. Synge, _Relativity: The General Theory_ , North-Holland series in physics (North-Holland Publishing Company, 1960), ISBN 9780444102799.
  * Kay and Wald (1991) B. S. Kay and R. M. Wald, Phys. Rept. 207, 49 (1991).
  * Christensen (1976) S. M. Christensen, Phys. Rev. D 14, 2490 (1976).
  * Deser and Levin (1997) S. Deser and O. Levin, Class. Quant. Grav. 14, L163 (1997), eprint arXiv:gr-qc/9706018.
  * Garbrecht and Prokopec (2004) B. Garbrecht and T. Prokopec, Class. Quant. Grav. 21, 4993 (2004), eprint arXiv:gr-qc/0404058.
  * Sriramkumar and Padmanabhan (1996) L. Sriramkumar and T. Padmanabhan, Classical and Quantum Gravity 13, 2061 (1996).
  * Espinosa-Portalés and Vennin (2022) L. Espinosa-Portalés and V. Vennin, JCAP 07, 037 (2022), eprint arXiv:2203.03505.
  * Martin and Vennin (2021a) J. Martin and V. Vennin, Phys. Rev. D 104, 085012 (2021a), eprint arXiv:2106.14575.
  * Martin and Vennin (2021b) J. Martin and V. Vennin, JCAP 10, 036 (2021b), eprint arXiv:2106.15100.
  * Hari and Kothawala (2021) Hari K. and D. Kothawala, Phys. Rev. D 104, 064032 (2021), eprint arXiv:2106.14496.
  * Hari and Kothawala (2023) Hari K. and D. Kothawala (2023), eprint arXiv:2307.16413.
  * Peeters (2007a) K. Peeters, Computer Physics Communications 176, 550 (2007a), eprint arXiv:cs/0608005.
  * Peeters (2007b) K. Peeters, arXiv:hep-th/0701238 (2007b).
  * Vines (2015) J. Vines, Gen. Rel. Grav. 47, 59 (2015), eprint arXiv:1407.6992.
  * Poisson (2004) E. Poisson, Living Rev. Rel. 7, 6 (2004), eprint arXiv:gr-qc/0306052.