\section{Introduction}

In general spacetimes the notion of particles loses its invariant meaning. One consequence of this observation is the detection of thermal particles in the Minkowski vacuum by a uniformly accelerated observer. This phenomenon is known as the Unruh effect \cite{Unruh:1976db}. It has immediate consequences for quantum states observed by accelerated observers. Their description becomes observer-dependent. In particular, the entanglement present in quantum states cannot be assumed to be independent of the observer \cite{Peres:2002wx, FuentesSchuller:2004xp}.

In the past, several investigations of entanglement in uniformly accelerating systems used nonlocalized Fock states, see e.g. \cite{FuentesSchuller:2004xp, Bruschi:2010mc, Richter:2015wha} and references therein. Therefore, it was not clear how to obtain an operational meaning of entanglement in these states. Furthermore, as pointed out in \cite{Dragan:2012hy, Ahmadi:2016fbd}, it is questionable that these simplified models capture the consequences of acceleration correctly. Another approach that is commonly used to study entanglement in the Minkowski vacuum is employing accelerated Unruh-de Witt detectors \cite{birrell1984quantum}. In these models, a two-level system is linearly coupled to a scalar field and the effect of finite time interactions is studied \cite{benatti2004entanglement}. In particular, it was shown that entanglement can be extracted from vacuum correlations \cite{Salton:2014jaa, Lorek:2014dwa}. The process of extracting entanglement from the vacuum sometimes is referred to as ``entanglement harvesting'' \cite{MartinMartinez:2012sg}. Since these effects are typically small, a direct observation is a challenging task. However, analog systems to study these and related phenomena were proposed \cite{delRey:2011jt, garcia2016entanglement, richter2017}. Another approach to localize quantum states was taken in \cite{Dragan:2012hy, Ahmadi:2016fbd, Kubicki:2016ngz, Dragan:2013dna}. Instead of making use of Unruh-de Witt detectors, the modes constituting the state itself were constructed to be localized. That is also the approach we are following in this work. 

The goal of this work is to give a fully analytic treatment of the effect of acceleration on general localized fermionic Gaussian states. These states are of great relevance, as the class of fermionic Gaussian states contains a broad variety of states such as vacuum states of quadratic Hamiltonians, thermal states and Bell states. Our approach is based on a method that was developed to study bosonic two-mode Gaussian states of localized wave packets \cite{Ahmadi:2016fbd}. Similarly to the bosonic case, we are able to formulate the transformation connecting inertial and accelerated observers as the action of a Gaussian quantum channel.

The framework developed in this work is then used to study the entanglement of the vacuum. In particular, we quantify the amount of entanglement that could possibly be harvested from the vacuum by local detectors and show that it increases with increasing acceleration. Furthermore, we are able to obtain the entanglement in Bell states, a problem that attracted a lot of attention \cite{FuentesSchuller:2004xp, Bruschi:2010mc, Richter:2015wha}, and show that this effect is due to an inevitable mismatch of the wave functions in Minkowski and Rindler spacetime.


The structure of this work is the following. In Sec.\ \ref{framework} we introduce the setting that we study in this work. In Sec.\ \ref{fermionicnoninertial} we calculate the transformation of general fermionic Gaussian states for the case of two uniformly accelerated observers with arbitrary accelerations. With these results in hand, we study the vacuum entanglement in Sec.\ \ref{vacuumcase}. In particular, we compare our findings to the case of bosonic Gaussian states. In Sec.\ \ref{bellstates}, we study the entanglement in Bell states. Finally, we give the conclusions in Sec.\ \ref{conclusions}.



\section{Framework}\label{framework}

\subsection{Outline}

As anticipated above, we study the effect of acceleration on spatially localized fermionic modes. These modes are solutions of the equation of motion for the Dirac field $\hat{\Psi}$. The Dirac equation is obtained by varying the action
\begin{equation}\label{actiondirac}
S_D=\int \text{d}^2 x \sqrt{-g}\left(i\bar{\hat{\Psi}}g_{\mu\nu}\gamma^\mu \partial^\nu \hat{\Psi}-m \bar{\hat{\Psi}}\hat{\Psi}\right),
\end{equation}
where $g$ denotes the determinant of the metric $g_{\mu\nu}$, $\gamma^\mu$ are the gamma matrices, $m$ is the mass of the field and $\bar{\hat{\Psi}}=\hat{\Psi}^\dagger \gamma^0$. The field can be expanded in a complete set of solutions of the Dirac equation.  Since we are interested in modes that are sufficiently localized, we choose to expand $\hat{\Psi}$ in terms of localized wavepackets $\phi_k^\pm$. This gives the expansion
\begin{equation}\label{expansionminkowski}
\hat{\Psi}=\sum_k \left(\phi_k^+ \hat{f}_k +\phi_k^- \hat{g}_k^\dagger \right)
\end{equation}
in Minkowski space, where the $\hat{f}_k^\dagger/\hat{g}_k^\dagger$ are the creation operators for the wavepackets $\phi_k^\pm$ of particles/antiparticles. The Dirac equation dictates the anticommutativity of the creation and annihilation operators; $\{\hat{f}_k , \hat{f}_l\}= 0$, $\{\hat{f}_k^\dagger , \hat{f}_l^\dagger\}= 0$  and $\{\hat{f}_k,\hat{f}_l^\dagger\}= \delta_{kl}$, where the same holds for $\hat{g}_k$ and $\hat{g}_k^\dagger$ with vanishing mixed anticommutators. Alternatively, we can carry out an equivalent expansion in Rindler space and obtain in terms of the respective creation and annihilation operators $\hat{d}_k^\dagger/\hat{e}_k^\dagger$ and $\hat{d}_k/\hat{e}_k$
\begin{equation}\label{expansionrindler}
\hat{\Psi}= \sum_k \left(\psi_k^+ \hat{d}_k +\psi_k^- \hat{e}_k^\dagger\right).
\end{equation}
Here the $\psi_k^\pm$ denote the respective wavepackets and the $\hat{d}_k$ and $\hat{e}_k$ satisfy the same algebra as the $\hat{f}_k$ and $\hat{g}_k$.

Let Alice and Bob be two observers that are each in possession of a (localized) mode of the fermion field $\hat{\Psi}$. These modes have negligible overlap and we denote these by $\phi_\text{I}^+$ and $\phi_\text{II}^+$, respectively. It follows that the corresponding operators commute; $[\hat{f}_\text{I}, \hat{f}_\text{II}^\dagger]=0$. Initially, these modes are prepared in a fermionic Gaussian state $\rho$ that is completely characterized by its first and second moments (covariance matrix). Then to describe the state of the shared pair of modes from an accelerated perspective,  Rindler space offers a convenient reference frame. The modes in both spacetimes are related to each other by a Bogolyubov transformation. Hence, the Gaussianity of state $\rho$ is preserved under such a transformation. The resulting state of Alice's and Bob's modes is thus again a Gaussian of the transformed modes $\psi_\text{I}^\pm$ and $\psi_\text{II}^\pm$. As for bosons \cite{Ahmadi:2016fbd, PhysRevA.63.032312}, the map transforming the covariance matrix of the inertial modes, into the state of the accelerated modes, is a trace preserving Gaussian map $\sigma^{(f)}\rightarrow\sigma^{(d)}$, which takes the form
\begin{equation}\label{trafocovmatrix}
\sigma^{(d)}=M\sigma^{(f)} M^T +N,
\end{equation}
where $M$ and $N$ are $8\times 8$ matrices \cite{Bravyi05} and  we use the superscripts $(f)$ and $(d)$ to indicate the modes that are used to calculate the covariance matrix. One of the main results of the present work is the exact analytic calculation of the matrices $M$ and $N$ for the generic case of two uniformly accelerated observers.






\subsection{Fermionic Gaussian states}\label{fergaustates}


While bosonic Gaussian states are intensively studied in the context of quantum information and extensive literature exists \cite{RevModPhys.77.513, RevModPhys.84.621}, fermionic Gaussian states are much less investigated. Therefore, we insert a brief discussion of these states. The most prominent members of the family of fermionic Gaussian states are vacuum states of quadratic Hamiltonians, thermal states and the Bell states. Even so fermionic Gaussian states share some similarities with their bosonic counterparts, there are crucial differences due to the anticommutativity of fermions. Starting from the creation and annihilation operators $\hat{f}_k$ and $\hat{f}_k^\dagger$ that satisfy the CAR algebra, i.e, $\{\hat{f}_k , \hat{f}_l\}= 0$ and $\{\hat{f}_k,\hat{f}_l^\dagger\}= \delta_{kl}$, we define the Majorana fermion operators  $\hat{c}_k$, as
\begin{equation}\label{defofmajorana}
\hat{c}_{2j-1}=\frac{1}{\sqrt{2}}(\hat{f}_j^\dagger+\hat{f}_j), \hspace{5mm} \hat{c}_{2j}=\frac{1}{i\sqrt{2}} (\hat{f}_j^\dagger-\hat{f}_j).
\end{equation}
These form a Clifford algebra 
\begin{equation}
\{\hat{c}_k,\hat{c}_l\}= \delta_{kl}.
\end{equation}
By definition, the density matrix of an even fermionic Gaussian state can be written as
\begin{equation}
\rho=C \,e^{-\frac{i}{2}\hat{c}^T A \hat{c}},
\end{equation}
where $A$ is a real, antisymmetric matrix and $C$ is the normalization \cite{Bravyi05}. That is, these states are thermal states of a quadratic Hamiltonian $H=\frac{i}{2}\hat{c}^T A \hat{c}$. Upon an $SO(2n)$ transformation $O$ with $n$ being the number of Majorana fermions $\hat{c}_i$, we can write the density matrix in the form
\begin{equation}
\rho=\frac{1}{2^n} \prod_{k=1}^n \left(1+i\lambda_k \hat{c}_{2k-1}' \hat{c}_{2k}'\right),
\end{equation}
where the $\lambda_i$ live in $[-1,1]$ and $\hat{c}'=O^T\hat{c}$. Further, using the Bloch-Messiah reduction \cite{bloch1962canonical}, we can find a basis of modes ${\hat{h}_k}$  such that every pure fermionic Gaussian state can be written in the form 
\begin{equation}
\rho=\prod_{k=1}^n \left(u_k+v_k \hat{h}_{k}^\dagger \hat{h}_{-k}^\dagger\right),
\end{equation}
where $|u_k|^2+|v_k|^2=1$. Even Gaussian states are completely characterized by the corresponding covariance matrix $\sigma_{kl}$ \cite{Bravyi05}.  The real antisymmetric covariance matrix $\sigma_{kl}$ for a Gaussian state $\rho$ is given by
\begin{equation}\label{covmatrix}
\sigma_{kl}= i \,Tr(\rho[\hat{c}_k,\hat{c}_l])= i \langle \hat{c}_k \hat{c}_l - \hat{c}_l \hat{c}_k \rangle.
\end{equation}
The diagonal elements in (\ref{covmatrix}) are vanishing due to the commutator and the covariance matrix is completely characterized by the elements
\begin{equation}
\sigma_{kl}= 2i \langle \hat{c}_k \hat{c}_l \rangle,\hspace{10mm}\text{with} \,\, k>l.
\end{equation}
All higher moments of $\rho$ can be obtained by Wick's theorem \cite{Bravyi05}. So far we concentrated on even states. While these are completely characterized by $\sigma_{kl}$, the description of odd Gaussian states requires, in general, also the knowledge of the first moments
\begin{equation}\label{firstmoments}
 Tr(\rho \hat{c}_k)=  \langle \hat{c}_k \rangle.
\end{equation}
These are naturally vanishing for even states due to the vanishing commutator with the parity operator $P$, $[P, \rho]=0$ \cite{Bravyi05}. More information regarding fermionic Gaussian states and Gaussian linear maps can, for example, be found in \cite{kraus2009pairing, de2013power, bravyi2005classical,  Bradler:2010st, greplova2016degradability}.



\subsection{Dirac field in Rindler spacetime}

We briefly introduce solutions of the Dirac equation that is obtained by varying action (\ref{actiondirac}) with respect to $\bar{\hat{\Psi}}$. For simplicity, we restrict our attention to the $1+1$ dimensional case. However, the generalization to $3+1$ dimensions is straightforward. Furthermore, we use units such that $c=\hbar=1$ throughout the entire work. We start by deriving the solutions of the Dirac equation in Minkowski coordinates $(t,x)$. The Dirac equation takes the form 
\begin{equation}
i\partial_t\hat{\Psi}=\left(-i\alpha_3\partial_x+m\beta\right)\hat{\Psi},
\end{equation}
where $m$ denotes the mass of the field and the matrices $\beta$ and $\alpha_3$ are given by
\begin{equation}\label{matrices}
\beta=\left( \begin{array}{cc}
1 & 0  \\
0& -1   \end{array} \right), \hspace{10mm} \alpha_3=\left( \begin{array}{cc}
0 & 1  \\
1& 0   \end{array}\right).
\end{equation} 
We find the solutions in Minkowski coordinates to be given by
\begin{align}
u_{k,\pm}=&\frac{1}{\sqrt{4\pi\omega_k}}  \left(
\begin{array}{c}
\sqrt{\omega_k\pm m}\\
\pm\sqrt{\omega_k\mp m}
\end{array}
\right)e^{\mp i\omega_k t+ikx}  ,
\end{align} 
where $\omega_k=\sqrt{k^2+m^2}$. These are normalized as given in (\ref{orthonormalization}). For further details, see \cite{birrell1984quantum}. Next, to describe the physics of an accelerated observer \cite{Crispino:2007eb, takagi1986vacuum}, we work in Rindler coordinates 
\begin{subequations}\label{coordtrafo}
\begin{align}
t=& \chi \sinh(a\eta),\\
x=&\chi \cosh(a\eta),
\end{align} 
\end{subequations}
 where $a$ is a parameter that does not have an immediate physical meaning. However, if we consider the worldline of a particle with a proper acceleration $\cal A$, namely $\chi=\frac{1}{\cal A}$, its proper time $\tau$ can be related
to the coordinate time $\eta$ as ${\cal A} \tau =a \eta$. In Rindler coordinates, the equation of motion obtained from (\ref{actiondirac}) is given by
\begin{equation}\label{diracinrindler}
i \frac{1}{a}\partial_\eta\hat{\Psi}=\left(m\chi \beta-\frac{i}{2} \alpha_3-i\chi\alpha_3\partial_{\chi}\right)\hat{\Psi},
\end{equation}
where $m$ denotes the mass of the field and the matrices $\beta$ and $\alpha_3$ are given by (\ref{matrices}). Solving the Dirac equation using the ansatz $w^\pm_{\text{I}\Omega}(\eta,\chi)= e^{\mp i\Omega\eta}w^\pm(\chi)$ for particles and antiparticles respectively, one finds the normalized solutions to be given by
\begin{widetext}
\begin{align}
w^\pm_{\text{I}\Omega}=&\sqrt{\frac{m \cosh(\frac{\pi \Omega}{a})}{2\pi^2 a}} \left(
\begin{array}{c}
K_{\pm i\frac{\Omega}{a}+\frac{1}{2}}(m\chi)+i K_{\pm i\frac{\Omega}{a}-\frac{1}{2}}(m\chi)\\
- K_{\pm i\frac{\Omega}{a}+\frac{1}{2}}(m\chi)+i K_{\pm i\frac{\Omega}{a}-\frac{1}{2}}(m\chi)
\end{array}
\right)  e^{\mp i\Omega\eta}, \hspace{13mm} \text{in Rindler wedge $I$},\\
w^\pm_{\text{II}\Omega}=&\sqrt{\frac{m \cosh(\frac{\pi \Omega}{a})}{2\pi^2 a}}  \left(
\begin{array}{c}
K_{\pm i\frac{\Omega}{a}+\frac{1}{2}}(-m\chi)+i K_{\pm i\frac{\Omega}{a}-\frac{1}{2}}(-m\chi)\\
- K_{\pm i\frac{\Omega}{a}+\frac{1}{2}}(-m\chi)+i K_{\pm i\frac{\Omega}{a}-\frac{1}{2}}(-m\chi)
\end{array}
\right) e^{\pm i\Omega\eta}, \hspace{6mm} \text{  in Rindler wedge $II$},
\end{align}
\end{widetext}
where $K_{ y}(x)$ is the modified Bessel function of the second kind of order $y$.  The Dirac inner product of the mode functions $\omega_1$ and $\omega_2$ is defined as $\left(\omega_1(x),\omega_2(x)\right)_\Sigma=\int \text{d}\Sigma^\mu \omega_1^\dagger(x)\gamma^0\gamma_\mu \omega_2(x)$, where $\text{d}\Sigma^\mu=n^\mu \text{d}s$ is space-like with normal vector $n^\mu$ and volume element $\text{d}s$. It has the properties $\left(\omega_1(x),\omega_2(x)\right)^*_\Sigma=\left(\omega_1^*(x),\omega_2^*(x)\right)_\Sigma =\left(\omega_2(x),\omega_1(x)\right)_\Sigma$. Here, in Rindler coordinates, we make the particular choice
\begin{equation}
\left(\omega_1,\omega_2\right)=\int \text{d}\chi\, \omega_1^\dagger\omega_2,
\end{equation}
where we suppressed the coordinate dependence of the mode functions. The positive and negative energy solutions satisfy the relations
\begin{equation}\label{orthonormalization}
(w^\pm_{\text{I}\Omega},w^\pm_{\text{I}\Omega'})= \delta( \Omega-\Omega'),\hspace{5mm}(w^\pm_{\text{I}\Omega},w^\mp_{\text{I}\Omega'})=0.
\end{equation}
The decomposition of $\hat{\Psi}$ in terms of localized modes is given by
\begin{equation}\label{fieldexpansion}
\hat{\Psi}=\sum_k \left(\phi_k^+ \hat{f}_k +\phi_k^- \hat{g}_k^\dagger \right)= \sum_k \left(\psi_k^+ \hat{d}_k +\psi_k^- \hat{e}_k^\dagger\right),
\end{equation}
where $\hat{f}_k/\hat{g}_k$ are the annihilation operators of particles/antiparticles in the Minkowski case and $\hat{d}_k/\hat{e}_k$ are the annihilation operators of particles/antiparticles in the Rindler case. These are related to the Rindler particle and antiparticle operators, $\hat{b}_\Omega$ and $\hat{a}_\Omega$,  by 
\begin{subequations}\label{relationdb}
\begin{align}
\hat{d}_\text{I}=&\int \text{d}\Omega (\psi_\text{I}^+, w^+_{\text{I}\Omega})\hat{b}_{\text{I}\Omega},\\
\hat{d}_\text{II}=&\int \text{d}\Omega (\psi_\text{II}^+, w^+_{\text{II}\Omega})\hat{b}_{\text{II}\Omega},\\
\hat{e}_\text{I}^\dagger=&\int \text{d}\Omega (\psi_\text{I}^-, w^-_{\text{I}\Omega})\hat{a}^{\dagger}_{\text{I}\Omega},\\
\hat{e}_\text{II}^\dagger=&\int \text{d}\Omega (\psi_\text{II}^-, w^-_{\text{II}\Omega})\hat{a}^{\dagger}_{\text{II}\Omega}.
\end{align}
\end{subequations}
The Minkowski vacuum $|0\rangle_\text{M}$ is related to the Rindler vacuum $|0\rangle_\text{R}$ by a squeezing operator $S$ as
\begin{equation}
|0\rangle_\text{M}=S |0\rangle_\text{R},
\end{equation}
where $S$ acts on $\hat{b}_{\text{I}\Omega}$ and $\hat{a}_{\text{II}\Omega}$  as
\begin{subequations}\label{squeezing}
\begin{align}
S^\dagger \hat{b}_{\text{I}\Omega} S=& \cos(r_\Omega) \hat{b}_{\text{I}\Omega}-\sin(r_\Omega) \hat{a}_{\text{II}\Omega}^{\dagger},\\
S^\dagger \hat{a}_{\text{II}\Omega} S=& \cos(r_\Omega) \hat{a}_{\text{II}\Omega}+\sin(r_\Omega) \hat{b}_{\text{I}\Omega}^{\dagger}
\end{align}
\end{subequations}
with $\tan(r_\Omega)=e^{-\frac{\pi |\Omega|}{a}}$. To obtain the transformations of the operators $\hat{b}_{\text{II}\Omega}$ and $\hat{a}_{\text{I}\Omega}$, one simply has to interchange $\hat{a}\leftrightarrow \hat{b}$ in (\ref{squeezing}). 


\subsection{Modes}\label{modes}
  
  \begin{figure}[t]

\includegraphics[scale=1]{figure1.pdf}
\caption{ Comparison between the real parts of the spatial modes $\phi^+$ and $\psi^+$ for the following choice of parameters: $x_0^{-1} = {\cal A}=0.1$, $L=2$, $\Omega_0=4.71$, $m=0.1$. These modes are localized and, therefore, a single proper acceleration ${\cal A}$ can be assigned.}
\label{picmodes}
\end{figure}

Most of our calculations are general and do not depend on the particular choice of modes. However, when considering concrete examples, namely in Sec.\ \ref{vacuumcase} and  Sec.\ \ref{bellstates}, we have to make a specific choice. Physically, cavity modes are suitable candidates for localized  solutions. In particular, motivated by the solutions obtained in \cite{Friis:2013eva}, for $\Lambda\in \{\text{I},\text{II}\}$, we consider the Minkowski modes
\begin{equation}\label{inercavity}
\phi_\Lambda^+(x,0)= C_\phi \, \xi_\phi(k) \, e^{-2\left(\frac{x_0}{L}\log(\frac{x}{x_0})\right)^2+ik x},
\end{equation}
where   $L$ is the width of the wave packet that is centred at $x_0$ and $\xi_\phi(k)$ is 
\begin{equation}
\xi_\phi(k)=\left(
\begin{array}{c}
\cos\left(\frac{\kappa}{2}\right)+\sin\left(\frac{\kappa}{2}\right)\\
\cos\left(\frac{\kappa}{2}\right)-\sin\left(\frac{\kappa}{2}\right)
\end{array}
\right)  e^{-\frac{i}{2}\kappa}e^{-i k x_0}
\end{equation}
with $\kappa=\arctan(\frac{m}{k})$ and $C_\phi$ being a constant of normalization. As in \cite{Ahmadi:2016fbd}, we use a Gaussian profile to sufficiently localize the mode.


For the accelerated case, we choose wave packets of particle solutions of the Dirac equation in Rindler space, (\ref{diracinrindler}). Specifically, the localized modes we work with are 
\begin{equation}\label{acccavity}
 \psi_\Lambda^+(\chi,0)= C_\psi \, \xi_\psi(\chi) \, e^{-2\left(\frac{\chi_0}{L}\log(\frac{\chi}{\chi_0})\right)^2},
\end{equation}
where  $\xi_\psi(\chi)$ is given by the spinor
\begin{equation}
\xi_\psi(\chi)=\left(
\begin{array}{c}
I_{i\frac{ \Omega }{a}-\frac{1}{2}}(m \chi )+i I_{i\frac{\Omega }{a}+\frac{1}{2}}(m \chi )\\
I_{i\frac{ \Omega }{a}-\frac{1}{2}}(m \chi )-i I_{i\frac{ \Omega }{a}+\frac{1}{2}}(m \chi )
\end{array}
\right)
\end{equation}
and $C_\psi$ is 
\begin{equation}
C_\psi=\tilde{C}_\psi\left( I_{-i\frac{ \Omega }{a}-\frac{1}{2}}\left(m\chi_0\right)-I_{-i\frac{ \Omega }{a}+\frac{1}{2}}\left(m\chi_0\right)\right)
\end{equation}
with the constant of normalization, $\tilde{C}_\psi$ and $I_{i\mu}(x)$ being the modified Bessel functions of the first kind. The modes are plotted in Fig.\ \ref{picmodes}. 
 
The physical motivation for choosing localized cavity modes is that these modes can possibly be detected by a local detector (in our case a cavity). Therefore, the choice of localized modes enables us to assign an operational meaning to quantities like entanglement. Furthermore, the modes (\ref{inercavity}) and (\ref{acccavity}) are purely positive frequency solutions. This guarantees, that the average particle and antiparticle numbers are zero in the Rindler vacuum, as expected. In consequence, the modes do necessarily have noncompact support. However, due to the Gaussian profile they are still sufficiently localized.




\section{Accelerated fermionic Gaussian states}\label{fermionicnoninertial}

In this section, we study the Bogolyubov transformation of fermionic Gaussian states  relating modes in Minkowski space to those of Rindler space. We give the analytical solution for the general case of two uniformly accelerated observers.


\subsection{Transformation of the Minkowski vacuum}

To study the effect of acceleration on the Minkowski vacuum, we restrict our attention to two localized fermion modes with respective creation operators $\hat{f}_1^\dagger$ and $\hat{f}_2^\dagger$. Then, using (\ref{defofmajorana}), we can  write the corresponding Majorana operators $\hat{c}_i$  as
\begin{equation}
\left(
\begin{array}{c}
\hat{c}_1\\
\hat{c}_2\\
\hat{c}_3\\
\hat{c}_4\\
\end{array}
\right)
=\frac{1}{\sqrt{2}}
\left(
\begin{array}{c}
\hat{f}_\text{I}^\dagger+\hat{f}_\text{I}\\
\frac{1}{i} (\hat{f}_\text{I}^\dagger-\hat{f}_\text{I})\\
\hat{f}_\text{II}^\dagger+\hat{f}_\text{II}\\
\frac{1}{i} (\hat{f}_\text{II}^\dagger-\hat{f}_\text{II})\\
\end{array}
\right).
\end{equation}
From (\ref{covmatrix}) it is clear that in the case of two modes, there are six independent entries in the fermionic covariance matrix. The same is true for the corresponding antiparticles. Therefore, we find the covariance matrix  of the Minkowski vacuum 
\begin{align}\label{covmatrixinertial}
\sigma^{(f)}_\text{M}=&\left(\begin{array}{cccccccc}
0 & 1& 0 & 0  & 0 & 0& 0 & 0     \\
-1 & 0& 0 & 0  & 0 & 0& 0 & 0    \\
0 & 0& 0 & 1 & 0 & 0& 0 & 0     \\
0 & 0 & -1 & 0 & 0 & 0& 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 1& 0 & 0 \\
0 & 0 & 0 & 0 & -1 & 0& 0 & 0 \\
0 & 0 & 0 & 0 & 0 & 0& 0 &1 \\
0 & 0 & 0 & 0 & 0 & 0& -1 & 0 \\
\end{array}\right),
\end{align}
where this expression describes two particle modes and two antiparticle modes. However, it generalizes trivially to an arbitrary number of particles and antiparticle modes; see Appendix \ref{appminkvac} for details. As it is known for pure states \cite{Bravyi05}, $\sigma^{(f)^T}_\text{M} \sigma^{(f)}_\text{M}=1$ and the first moments vanish, as the vacuum is an even state.  To obtain  the covariance matrix $\sigma^{(d)}_\text{M}$ of the Minkowski vacuum in the Rindler frame, we need to transform the operators $\hat{f}$ and $\hat{g}$ (for details see Appendix \ref{appminkvactrafo}). We find the covariance matrix of the Minkowski vacuum in the Rindler basis to be of the form
\begin{equation}\label{covmatrixnoninertial}
\sigma^{(d)}_{}=\left(\begin{array}{cc}
\sigma_{+}^{(d)} & \sigma_{c}^{(d)}    \\
\tilde{\sigma}_{c}^{(d)}  & \sigma_{-}^{(d)}  
\end{array}\right),
\end{equation}
where the block $\sigma_{+}^{(d)}$ captures the effect of acceleration on the respective mode and displays the thermal character of the Unruh effect ($\sigma_{-}^{(d)}$ is the analog for antiparticles). $\sigma_{c}^{(d)} $ and $\tilde{\sigma}_{c}^{(d)}$ describe the correlations between particles and antiparticles that arise due to the acceleration. Explicitly, (\ref{covmatrixnoninertial}) reads
\begin{widetext}
\begin{equation}\label{covmatrixnoninertialdntriv}
\sigma^{(d)}_{}=\left(\begin{array}{cccccccc}
0 & N_\text{I}^+& 0 & 0   & 0 & 0& \Im[N_\text{I,II}^+] & \Re[N_\text{I,II}^+]   \\
-N_\text{I}^+ & 0& 0 & 0   & 0 & 0& \Re[N_\text{I,II}^+] & -\Im[N_\text{I,II}^+]  \\
0 & 0 & 0 & N_\text{II}^+  &   -\Im[N_\text{I,II}^-] & -\Re[N_\text{I,II}^-]   & 0 & 0  \\
0& 0& -N_\text{II}^+ & 0  &  -\Re[N_\text{I,II}^-] & \Im[N_\text{I,II}^-]  & 0 & 0 \\
0 & 0& \Im[N_\text{I,II}^-] & \Re[N_\text{I,II}^-]  &  0 & N_\text{I}^-& 0 & 0  \\
0 & 0& \Re[N_\text{I,II}^-] & -\Im[N_\text{I,II}^-]  &  -N_\text{I}^- & 0& 0 & 0 \\
-\Im[N_\text{I,II}^+] & -\Re[N_\text{I,II}^+]   & 0 & 0  &  0 & 0 & 0 & N_\text{II}^- \\
-\Re[N_\text{I,II}^+] & \Im[N_\text{I,II}^+]  & 0 & 0& 0& 0& -N_\text{II}^- & 0 
\end{array}\right),
\end{equation}
\end{widetext}
where the entries of the matrix are given by
\begin{subequations}\label{covmatelefermi}
\begin{align}
N_\text{I}^\pm=& 1-2\int \text{d}\Omega\, \frac{|(\psi_\text{I}^\pm , w^\pm_{\text{I}\Omega})|^2}{1+e^{\frac{2\pi\Omega}{a}}} \, ,\label{diagonalelement1}\\
N_\text{II}^\pm=& 1-2\int \text{d}\Omega\, \frac{|(\psi_\text{II}^\pm , w^\pm_{\text{II}\Omega})|^2}{1+e^{\frac{2\pi\Omega}{a}}}\,,\label{diagonalelement2}\\
N_\text{I,II}^\pm=&-2 \int \text{d}\Omega\, \frac{(\psi_\text{I}^\pm, w^\pm_{\text{I}\Omega}) (\psi_\text{II}^\mp, w^\mp_{ \text{II}\Omega})}{1+e^{\frac{2\pi\Omega}{a}}}\,e^{\frac{\pi\Omega}{a}}\label{diagonalelement3}
\end{align}
\end{subequations}
and we used $(\psi_\text{I}^\pm, w^\pm_{\text{I} \Omega})^*=(\psi_\text{I}^\mp, w^\mp_{\text{I} \Omega})$. We note that (\ref{covmatrixnoninertial}) reduces to  (\ref{covmatrixinertial}) in the limit of vanishing accelerations, $a\to 0$. That is, the correlated noise is absent in this limit.





\subsection{Transformation of general fermionic Gaussian states}\label{trafoofgeneralstates}


Using the results of the previous section, we derive the transformation for general fermionic Gaussian states. We recall that the transformed covariance matrix $\sigma^{(d)}$ can be obtained using (\ref{trafocovmatrix}).  Therefore, we can completely describe the transformation of general Gaussian states by the two matrices $M$ and $N$ that we calculate in the following. To obtain $M$, we consider the transformation of the first moments (\ref{firstmoments}). The Bogolyubov transformations can be obtained from (\ref{fieldexpansion}) and read
\begin{subequations}\label{gernaltrafo1}
\begin{align}
\hat{d}_l=&\sum_k (\psi_l^+, \phi_k^+) \hat{f}_k + (\psi_l^+, \phi_k^-) \hat{g}_k^\dagger ,\\
\hat{e}_l^\dagger=&\sum_k (\psi_l^-, \phi_k^+) \hat{f}_k +(\psi_l^-, \phi_k^-)  \hat{g}_k^\dagger  ,
\end{align}
\end{subequations}
where, in the following, we denote the overlaps by $\alpha_\text{I}^+=(\psi_\text{I}^+, \phi_\text{I}^+)$, $\beta_\text{I}^+=-(\psi_\text{I}^+, \phi_\text{I}^-)$, $\alpha_\text{I}^-=(\psi_\text{I}^-, \phi_\text{I}^-)^*$, $\beta_\text{I}^-=-(\psi_\text{I}^-, \phi_\text{I}^+)$, and accordingly for modes $\text{II}$. We note that both the Minkowski wave packets $\phi_{\text{I}/\text{II}}^\pm$ as well as the Rindler wave packets $\psi_{\text{I}/\text{II}}^\pm$ are composed of either only  particles or only  antiparticles. Therefore, the overlaps $\beta_{\text{I}/\text{II}}^\pm$ are vanishing and we can write (\ref{gernaltrafo1}) as
\begin{equation}\label{gernaltrafo2}
\hat{d}_l=\alpha_\text{I}^+ \hat{f}_l,\hspace{10mm} \hat{e}_l=\alpha_\text{I}^-  \hat{g}_l\, .
\end{equation}
The transformation of the first moments can be written as
\begin{equation}
\langle \hat{c}^{(d)}\rangle= M\langle \hat{c}^{(f)}\rangle
\end{equation}
with the block-diagonal $8\times 8$-matrix
\begin{equation}\label{Mmatrix}
M=\left(\begin{array}{cccc}
M(\alpha_\text{I}^+) & 0 &0&0     \\
0 & M(\alpha_\text{II}^+)  &0&0\\
0&0&M(\alpha_\text{I}^-)&0    \\
0&0&0&M(\alpha_\text{II}^-)    
\end{array}\right),
\end{equation}
where we defined $M(\cdot)$ as the matrix given by
\begin{equation}
M(\cdot)=\left(\begin{array}{cc}
\Re[\cdot] & \Im[\cdot]      \\
-\Im[\cdot] & \Re[\cdot]      \\
\end{array}\right).
\end{equation}
We now gathered all quantities that we need to specify the transformation of a general Gaussian state characterized by $\langle \hat{c}^{(f)}\rangle$ and $\sigma^{(f)}$. First, we  explicitly calculate the noise matrix $N$. The matrix describing the noise arising due to acceleration is given by
\begin{equation}\label{noisematrix}
N=\sigma^{(d)}-M \sigma^{(f)}_\text{M} M^T,
\end{equation}
where $\sigma^{(d)}$ is given in (\ref{covmatrixnoninertial}). Therefore, in total we obtain
\begin{equation}\label{explnoise}
N =\left(\begin{array}{cc}
N_{+} & N_{\text{vac}}    \\
\tilde{N}_{\text{vac}}   & N_- 
\end{array}\right),
\end{equation}
where $N_{\text{vac}}$ and $\tilde{N}_{\text{vac}}$ are given by the off-diagonal blocks of (\ref{covmatrixnoninertial}) and
\begin{widetext}
\begin{equation}
N_{\pm}=\left(\begin{array}{cccc}
0 & N_\text{I}^\pm-|\alpha_{\text{I}}^\pm|^2& 0 & 0     \\
-N_\text{I}^\pm+|\alpha_{\text{I}}^\pm|^2 & 0& 0 & 0    \\
0 & 0 & 0 & N_{\text{II}}^\pm-|\alpha_{\text{II}}^\pm|^2     \\
0& 0& -N_{\text{II}}^\pm+|\alpha_{\text{II}}^\pm|^2 & 0
\end{array}\right).
\end{equation}
\end{widetext}
With the matrices $M$ and $N$ we now have all data that is necessary to completely characterize the effect of acceleration to fermionic Gaussian states. In the following, we use the developed formalism to quantify vacuum entanglement and, as a example of a nonvacuum state, we study a maximally entangled state, the Bell state $|B\rangle=\frac{1}{\sqrt{2}}(|00\rangle+|11\rangle)$.





\section{Vacuum entanglement}\label{vacuumcase}


\subsection{Vacuum entanglement for fermions}

In this section, we study the entanglement in the vacuum as seen by an accelerated observer. As a measure of entanglement we employ the logarithmic negativity. In \cite{eisler2015partial}, it was shown that the partial transpose of a fermionic  Gaussian state is, in general, not Gaussian and therefore it is difficult to calculate the negativity exactly. However, it was shown that a lower bound can be obtained \cite{eisert2016entanglement, eisler2015partial}. Building on this construction, we derive the respective bound  $\tilde{\mathcal{E}}_\mathcal{N}$ for our case in Appendix \ref{lognegativity}. It is given by
\begin{align}\label{lowerboundfermiond0}
\tilde{\mathcal{E}}_\mathcal{N}=&\ln(\frac{1}{2}(1+N_\text{I}^+N_\text{II}^-+|N_\text{I,II}^+|^2  +\nonumber\\
+& \Re[\sqrt{(N_\text{I}^+-N_\text{II}^-)^2-4|N_\text{I,II}^+|^2}]+\nonumber\\
+& \Im[\sqrt{(N_\text{I}^+-N_\text{II}^-)^2-4|N_\text{I,II}^+|^2}] ))
\end{align}
where $N_\text{I}^+$, $N_\text{II}^-$ and $N_\text{I,II}^+$ are the elements of the covariance matrix (\ref{covmatrixnoninertial}) that are given by (\ref{diagonalelement1}) - (\ref{diagonalelement3}). 


\begin{figure}[t]
\includegraphics[scale=1.05]{figure2.pdf}
\caption{Negativity $\tilde{\mathcal{E}}_\mathcal{N}$ of the Minkowski vacuum  as a function of the proper accelerations. The parameters are $m=0.1$, $L=2$ and $\Omega_0=4.71$. The vacuum entanglement increases with increasing acceleration.}
\label{picnegvac}
\end{figure}


To obtain the entanglement of the vacuum, we choose the modes as described in Sec.\ \ref{modes} and solve (\ref{lowerboundfermiond0}) numerically. We show the results for the entanglement between two modes of different accelerations in Fig.\ \ref{picnegvac}.

Firstly, it is important to emphasize that $\tilde{\mathcal{E}}_\mathcal{N}$ is a lower bound for the negativity and therefore it is difficult to make strong quantitative statements. However, it is reasonable to assume that the actual value of the negativity is close to the lower bound $\tilde{\mathcal{E}}_\mathcal{N}$. Therefore, in the following, we refer to $\tilde{\mathcal{E}}_\mathcal{N}$ as negativity, keeping this subtlety in mind.

The entanglement between the state of a mode in Rindler wedge I and its respective counterpart in wedge II increases with increasing acceleration, see Fig.\ \ref{picnegvac}. For the parameter regime we are interested in, the entanglement is of the order $10^{-12}$, as measured by the logarithmic negativity. That is, the particles that are produced due to the Unruh effect are correlated across the acceleration horizon. Interestingly, the vacuum entanglement has an operational meaning, as it can, in principle, be extracted by suitable detectors \cite{Salton:2014jaa, Lorek:2014dwa}.



\subsection{Comparison to bosons}\label{bosons}

The aim of this section is to compare our results for fermionic states to the results in the  case of localized modes of a massive scalar field studied in \cite{Ahmadi:2016fbd}. For bosonic Gaussian states the covariance matrix $\sigma^{(d)}_{\text{bos}}$ of the vacuum takes the form 
\begin{equation}\label{covmatrixnoninertialbosons}
\sigma^{(d)}_{\text{bos}}=
\left(\begin{array}{cccc}
N_\text{I}^b & 0& \Re[N_\text{I,II}^b] & \Im[N_\text{I,II}^b]      \\
0 & N_\text{I}^b& \Im[N_\text{I,II}^b] & -\Re[N_\text{I,II}^b]     \\
\Re[N_\text{I,II}^b] & \Im[N_\text{I,II}^b] & N_\text{II}^b &       \\
\Im[N_\text{I,II}^b] & \-Re[N_\text{I,II}^b] & 0 & N_\text{II}^b
\end{array}\right),
\end{equation}
where the matrix elements are given by
\begin{subequations}\label{covmatbosonele}
\begin{align}
N_{\text{I}/\text{II}}^b=& 1+2\int \text{d}\Omega \frac{|(\psi_\text{I/II}, w_{\text{I}/\text{II} \Omega})|^2 }{e^{\frac{2\pi\Omega}{a}}-1},\\
N_\text{I,II}^b=& 2\int \text{d}\Omega \frac{(\psi_\text{I}, w_{\text{I} \Omega}) (\psi_\text{II}, w_{\text{II} \Omega})}{e^{\frac{2\pi\Omega}{a}}-1}\,e^{\frac{\pi\Omega}{a}},
\end{align}
\end{subequations}
where $\psi_\text{I}$, $\psi_\text{II}$ are localized modes and the $w_{\text{I}/\text{II} \Omega}$ are solutions of the Klein-Gordon equation spanning the Hilbert space for regions I and II, respectively. The effect of acceleration on bosons is manifest in the covariance matrix (\ref{covmatrixnoninertialbosons}). While the Minkowski vacuum is described by the covariance matrix $\sigma^{(\text{bos})}_{M}=\mathbb{I}$, for nonvanishing acceleration, correlations build up and the off-diagonal elements increase.  As expected, the occupations are governed by the Bose-Einstein distribution that are explicit in the matrix elements (\ref{covmatbosonele}). This contrasts with the case of fermions where the occupation is characterized by the Fermi-Dirac distribution, (\ref{covmatelefermi}). Furthermore, due to the presence of antiparticles for fermions, the correlations that are build up in fermionic Gaussian states are between particles and antiparticles, while for bosons particles and antiparticles are identical. 

Besides the differences outlined above, the results for fermions qualitatively agree with the findings for bosons. Quantitatively,  the entanglement we observed for fermions is less than in the bosonic case. In the parameter regime we studied, bosons develop ten times more entanglement \cite{Ahmadi:2016fbd}. However, one should keep in mind that we were calculating a lower bound for the negativity and that the amount of entanglement also depends on the explicit choice of modes. 


In the next section, we go beyond vacuum entanglement and study the degradation of entanglement for Bell states.


\section{Entanglement in Bell states}\label{bellstates}

Interestingly, Bell states of fermions are fermionic Gaussian states, while maximally entangled states of two bosons are nonGaussian. This fact, enables us to apply our formalism to Bell states. Therefore, as an example for a nonvacuum state, we consider a maximally entangled state $|B\rangle$ of two particles in the following. We choose
\begin{equation}\label{bellstate}
|B\rangle=\frac{1}{\sqrt{2}}(|00\rangle+|11\rangle)=\frac{1}{\sqrt{2}}(1+\hat{f}_1^\dagger \hat{f}_2^\dagger)|0\rangle.
\end{equation}
This state is an even Gaussian state and it is described by the following covariance matrix (cf.\ Appendix \ref{appminkbell})
\begin{equation}
\sigma^{(f)}_{\text{Bell}}=\left(\begin{array}{cccccccc}
0 & 0&0 &  1 & 0& 0& 0& 0    \\
0 & 0& 1& 0 & 0& 0& 0& 0   \\
0 & - 1 & 0 & 0  & 0& 0& 0& 0  \\
- 1 & 0 & 0 & 0& 0& 0& 0& 0   \\
0 & 0 & 0 & 0& 0& 1& 0& 0   \\
0 & 0 & 0 & 0& -1& 0& 0& 0   \\
0 & 0 & 0 & 0& 0& 0& 0& 1   \\
0 & 0 & 0 & 0& 0& 0& -1& 0   
\end{array}\right).
\end{equation}
The transformed state is obtained according to equation (\ref{trafocovmatrix}) with the matrix $M$ given in (\ref{Mmatrix}) and $N$ given in (\ref{explnoise}). Then the transformed covariance matrix, neglecting some terms much smaller than 1, is given by
\begin{widetext}
\begin{equation}\label{transformedcovmatbell}
\sigma^{(d)}_\text{Bell}=\left(\begin{array}{cccc}
0 & N_\text{I}^+-|\alpha_\text{I}^+|^2 & (\sigma^{(d)}_\text{Bell})_{13} & (\sigma^{(d)}_\text{Bell})_{14}   \\
_-N_\text{I}^++|\alpha_\text{I}^+|^2& 0&   (\sigma^{(d)}_\text{Bell})_{23}   &     (\sigma^{(d)}_\text{Bell})_{24}     \\
- (\sigma^{(d)}_\text{Bell})_{13} & -  (\sigma^{(d)}_\text{Bell})_{23} & 0 & N_\text{II}^+ -|\alpha_\text{II}^+|^2    \\
-(\sigma^{(d)}_\text{Bell})_{14}  & -  (\sigma^{(d)}_\text{Bell})_{24} & -N_\text{II}^++|\alpha_\text{II}^+|^2 & 0
\end{array}\right),
\end{equation}
\end{widetext}
where 
\begin{subequations}
\begin{align}
(\sigma^{(d)}_\text{Bell})_{13}=&\Im[\alpha_\text{II}^+] \Re[\alpha_\text{I}^+] + \Im[\alpha_\text{I}^+] \Re[\alpha_\text{II}^+],\\
(\sigma^{(d)}_\text{Bell})_{14}=& \frac{1}{2} \left(\alpha_\text{I}^+ \alpha_\text{II}^+ + (\alpha_\text{I}^+\alpha_\text{II}^+)^*\right), \\
(\sigma^{(d)}_\text{Bell})_{23}=&\frac{1}{2} \left(\alpha_\text{I}^+ \alpha_\text{II}^+ + (\alpha_\text{I}^+\alpha_\text{II}^+)^*\right),\\
(\sigma^{(d)}_\text{Bell})_{24}=&-\Im[\alpha_\text{II}^+] \Re[\alpha_\text{I}^+] - \Im[\alpha_\text{I}^+] \Re[\alpha_\text{II}^+]
\end{align}
\end{subequations}
with $\alpha_{\text{I}/\text{II}}^+$ as defined in (\ref{gernaltrafo2}). In (\ref{transformedcovmatbell}), we neglected contributions from the vacuum noise matrix $N$, as these are insignificant corrections to $\sigma^{(d)}_\text{Bell}=M\sigma^{(f)}_\text{Bell} M^T$. 

To quantify the entanglement we employ a lower bound for the negativity given by 
\begin{align}\label{lowerboundfermionbell}
\tilde{\mathcal{E}}_\mathcal{N}=\ln(\frac{1}{2}+&\frac{\Re[\nu^+\nu^- +\nu^+-\nu^-]}{2}+\nonumber\\
+&\frac{\Im[\nu^+\nu^-+\nu^+-\nu^-]}{2}),
\end{align}
where the $\nu^\pm$ can be calculated from (\ref{transformedcovmatbell}); cf.\ (\ref{lowerboundneg}) in Appendix \ref{lognegativity}.

\begin{figure}
\centering
\includegraphics[scale=1]{figure3.pdf}
\caption{Negativity $\tilde{\mathcal{E}}_\mathcal{N}$ of a Bell state as a function of the proper acceleration $\mathcal{A}=\mathcal{A}_\text{I}=\mathcal{A}_\text{II}$. The parameters are $m=0.1$, $\Omega_0=4.71$ and $L=2$. We normalized $\tilde{\mathcal{E}}_\mathcal{N}$ such that for vanishing accelerations $\tilde{\mathcal{E}}_\mathcal{N}=1$. The blue continuous line gives the numerical results, while the gray dotted line is an extrapolation to $\mathcal{A}=0$. From an inertial perspective Bell states are maximally entangled. However, for an accelerated observer, entanglement gets degraded with increasing acceleration.}
\label{picbell}
\end{figure}

We show the results for the entanglement of state $|B\rangle$ in Fig.\ \ref{picbell}. For vanishing acceleration the modes are maximally entangled, while the entanglement decreases with increasing acceleration. The constraint $\mathcal{A}L\lesssim 1$ that ensures that a single acceleration can be associated to a localized mode prevents us from studying arbitrarily large accelerations. 

The physical reason for the degradation of entanglement is that the overlap between the inertial wavepackets $\phi^+_k$ and the accelerated ones $ \psi_\Omega^+$ decreases with increasing acceleration. This is an inevitable effect of the Bogolyubov transformation connecting Minkowski and Rindler solutions. Therefore, with increasing acceleration the mismatch between the modes increases and the entanglement in the initially maximally entangled Bell state decreases. Interestingly, the effect of Unruh radiation is negligibly small. That is, the production of thermal particles due to the Unruh effect has only a very small effect to entanglement degradation of states of localized modes.






\section{Conclusions}\label{conclusions}


In this work, we developed a general framework to describe the effect of acceleration on arbitrary fermionic Gaussian states of two localized modes. We formulated the transformation between the quantum state seen by initial observers and the state seen by uniformly accelerated ones as the action of a fermionic Gaussian channel, and completely characterized this channel. This enabled us to study how the entanglement of the vacuum and the entanglement in Bell states is affected by acceleration.

We found that vacuum entanglement is enhanced by acceleration. That is, the correlations in the vacuum that build up due to the Unruh effect lead to entanglement between particles and antiparticles. In particular, as our framework is employing localized modes, entanglement has an operational meaning and can be extracted by suitable local detectors; a process that is sometimes referred to as ``entanglement harvesting''.  We also compare these findings  with the case of bosons \cite{Ahmadi:2016fbd} and find a qualitative agreement.

Furthermore, we quantified the entanglement degradation in a localized maximally entangled state of fermions. The degradation is due to an increasing mismatch between the initial modes and the accelerated modes that causes their overlap to decrease and, therefore, entanglement decreases as well. The effect of particle creation on entanglement degradation is negligible.

We emphasize that due to the localization of the modes, the framework presented in this work can be applied to quantum information protocols, in scenarios where local gravitational effects or effects due to acceleration are not negligible. In particular, the entanglement can be exploited as a resource.

For the future, we are interested in extending our studies to localized modes of fermions with arbitrary mutual separation, as well as studying localized Gaussian states in higher dimensional spacetimes.





\begin{acknowledgements}

B.R. thanks the Institute of Theoretical Physics at the University of Warsaw for hospitality during his visit. B.R. and Y.O. acknowledge support from Funda\c{c}\~{a}o para a Ci\^{e}ncia e a Tecnologia (Portugal), namely through programmes PTDC/POPH/POCH and projects UID/EEA/50008/2013, IT/QuSim, IT/QuNet, ProQuNet, partially funded by EU FEDER, from the EU FP7 project PAPETS (GA 323901), and from the JTF project NQN (ID 60478). Furthermore, B.R. acknowledges the support from the DP-PMI and FCT through scholarship SFRH/BD/52651/2014. K.L. and A.D. would like to acknowledge financial support from the National Science Center, Sonata BIS Grant No. DEC-2012/07/E/ST2/01402.
 
\end{acknowledgements}


