# Exploring Bell inequalities and quantum entanglement in vector boson
scattering

R. A. Morales

###### Abstract

Quantum properties of vector boson scattering $V^{\prime}_{1}V^{\prime}_{2}\to
V_{1}V_{2}$, related to entanglement and violation of Bell inequalities, are
explored in this paper. The analysis is based on the construction of the
polarization density matrix associated to the final state $V_{1}V_{2}$ by
means of the computation of the corresponding tree level amplitudes within the
Standard Model. The aim of this work is to determine the regions of the phase
space where the final vector bosons after the scattering result entangled and
if is it possible to test the Bell inequalities in those regions. We found
that in all cases the entanglement is present. The amount of it depends on the
process and the Maximally Entangled state is reached in some particular
channels. Concerning the Bell inequality, it could be also tested in certain
kinematical regions for some of these processes. This work is a first step in
the analysis of these quantum properties for this kind of processes and it is
postponed for future studies the reconstruction of the polarization density
matrix and the related quantum parameters from experimental data through
Monte-Carlo simulations using quantum tomography techniques.

##  1 Introduction

Quantum entanglement plays a fundamental role in communication between
particles separated by macroscopic distances [1]. It is due to the presence of
correlations in quantum systems which cannot be replicated by classical
systems. One of the most relevant consequences is the violation of Bell
inequality [2], a nonviable fact in any theory consistent with the classical
concepts of realism and locality. This special feature of quantum systems was
the key for discovering quantum technologies as cryptography [3],
teleportation [4] and quantum computation [5]. Furthermore, theoretical
physicists make remarkable progresses on this topic in the context of quantum
field theories, see for instance the review [6].

The study of quantum entanglement and Bell inequalities in the context of High
Energy Physics, from both phenomenological and experimental point of views,
has received a very recent attention. In particular, it provides the
possibility of investigate quantum entanglement in a relativistic limit with
the highest possible energies at colliders. Concretely, testing Bell
inequalities have been proposed in $e^{-}e^{+}$ collisions [7], charmonium [8,
9] and positronium [10] decays, neutrino oscillations [11, 12] and more
recently for Higgs decays into gauge bosons [13, 14, 15, 16, 17]. Also
proposals to test Bell inequalities have been made in systems composed by
mesons [18, 19, 20, 21], top-quarks [22, 23, 24, 25, 26, 27, 28, 29], tau-
lepton and photon pairs [30, 31] and two massive gauge bosons [32, 17].

Violation of Bell inequalities has entanglement as a necessary but not
sufficient condition, then both concepts deserve dedicated attention. In
particular, the violation of such inequalities probes that there is no hidden
variable theory for accomplished the generated entanglement. The aim of this
work is to perform a systematic study of them through vector boson scattering
(VBS) processes in the SM, following closely [15, 16, 17, 33]. The latter of
these previous works has presented a complete theoretical analysis of
entanglement in QED processes. The others were focused on testing both
entanglement and Bell inequalities for Higgs boson decays into massive gauge
bosons and diboson production at the LHC and future lepton colliders by means
of Monte-Carlo simulations. In particular, a general method for the
experimental reconstruction of the density matrix and related entanglement
measurements from angular decay data was developed in [16]. This technique is
known as quantum tomography and, in the context of High Energy Physics, was
previously applied to lower-dimension systems, see for instance the $t\bar{t}$
case [22]. Now the VBS processes allow to study a variety of bipartite systems
with more intricate density matrices compared to those studied in
$t\bar{t}$-pair and Higgs boson decays.

Undoubtedly, the relevance of the VBS processes for examining the deepest
structure of the electroweak (EW) interactions in the SM is well established
in the community. In particular, the precise cancellation of potentially large
contributions among diagrams with only trilinear gauge self-couplings, quartic
gauge self-couplings and with the Higgs boson, is responsible for unitarity
restoration in this kind of processes. Then probing these interactions reveal
the dynamics behind the Higgs mechanism. In this line, there is an active
program of experimental searches by ATLAS and CMS Collaborations, see for
instance the recent reviews [34, 35]. Moreover, VBS is a suitable observable
for new physics in the EW sector looking for anomalous triple (aTGC) and
quartic (aQGC) gauge couplings and also with Higgs boson interactions. Taking
these considerations into account, the goal of this paper is to perform
theoretical predictions of entanglement and violation of Bell inequalities for
different VBS processes, which are not explored in the literature as far as we
know. In particular, this analysis allows us to locate kinematical regions of
the phase space where interesting quantum mechanical measurements might be
performed. The determination of the related quantities from Monte-Carlo
simulations in the complete collider events, i.e. a quantum tomography
analysis, is beyond the scope of this work and it is postponed for a future
study.

The paper is organized as follows: Section 2 presents the theoretical
framework connecting the scattering amplitudes and the density matrix
formalism’s. The quantifiers related to entanglement and Bell inequalities are
also introduced in this section. The main results of this work for the VBS
processes are collected in Section 3. This section is closed with a brief
discussion of a possible quantum tomography analysis for this kind of
processes at colliders. We summarize the main findings and future perspectives
in Section 4. The appendices contain the details for the amplitude
computation, additional plots and analytical expressions for $W^{\pm}\gamma\to
W^{\pm}\gamma$ and $\gamma\gamma\to W^{+}W^{-}$ processes.

##  2 Formalism

The traditional approach to QFT in Particle Physics focuses on scattering
amplitudes by means of Feynman diagrams. In this work, we will consider
2$\to$2 particle scattering processes among vector bosons within the SM,

$V^{\prime}_{1}(p^{\prime}_{1},s^{\prime}_{1})+V^{\prime}_{2}(p^{\prime}_{2},s^{\prime}_{2})\to
V_{1}(p_{1},s_{1})+V_{2}(p_{2},s_{2})$ (2.1)

where $V_{i}$ ($V^{\prime}_{i}$) denotes final (initial) photons, $W^{\pm}$ or
$Z$ bosons with momentum $p_{i}$ ($p^{\prime}_{i}$) and polarization $s_{i}$
($s^{\prime}_{i}$), for $i=1,2$. In particular, ‘$\pm$’ denote the transverse
polarizations and ‘$0$’ is used for the longitudinal one (only for the massive
gauge bosons).

In the language of Quantum Information Theory (QIT), knowledge about the
quantum system is represented by the density matrix $\rho$. Seminal works [36,
37, 38] relate the density matrix formalism with the S-matrix operator for the
scattering processes. In this work, we will focus on the quantum entanglement
of the vector boson polarizations $s_{i}$. The bipartite final state system
$|f\rangle$ is defined in the Hilbert space
$\mathcal{H}_{1}\otimes\mathcal{H}_{2}$ where each $\mathcal{H}_{i}$ has
dimension $d_{i}=dim(\mathcal{H}_{i})$ equals to 2 for photons or to 3 for
massive $W^{\pm}$ and $Z$ gauge bosons. In QIT language, photons correspond to
qubits and massive gauge bosons are qutrits. In addition, the chosen basis
corresponds to the spin-eigenvalue assignment for the outcomes $\\{+,-\\}$ and
$\\{+,0,-\\}$, respectively. In this sense, we write
$|f\rangle=|s_{1}\rangle\otimes|s_{2}\rangle$ where momenta is omitted for
simplicity. Similar definitions for the initial state
$|i\rangle=|s^{\prime}_{1}\rangle\otimes|s^{\prime}_{2}\rangle$ in the Hilbert
space $\mathcal{H}^{\prime}_{1}\otimes\mathcal{H}^{\prime}_{2}$ are
implemented. The unitary time evolution of two incoming particles state
$|i\rangle$ into two outgoing particles state $|f\rangle$ is given by the
scattering amplitude $\mathcal{M}(V^{\prime}_{1}V^{\prime}_{2}\to V_{1}V_{2})$
which is related to the S-matrix elements by momentum conservation relation
$S_{fi}=i(2\pi)^{4}\delta^{(4)}(p^{\prime}_{1}+p^{\prime}_{2}-p_{1}-p_{2})\mathcal{M}(V^{\prime}_{1}V^{\prime}_{2}\to
V_{1}V_{2})$. For the present calculation, these amplitudes were computed in
the SM at tree level and relevant details on the kinematics are summarized in
the Appendix A.

The density matrix at $t_{0}=-\infty$ for unpolarized initial state is

$\rho(t_{0}=-\infty)=|i\rangle\langle
i|=\sum_{s^{\prime}_{1},\,s^{\prime}_{2}}\frac{1}{d^{\prime}_{1}d^{\prime}_{2}}|s^{\prime}_{1}\,s^{\prime}_{2}\rangle\langle
s^{\prime}_{1}\,s^{\prime}_{2}|$ (2.2)

where $d^{\prime}_{i}=dim(\mathcal{H}^{\prime}_{i})$ and
$|s^{\prime}_{1}\,s^{\prime}_{2}\rangle=|s^{\prime}_{1}\rangle\otimes|s^{\prime}_{2}\rangle$
for a compact notation. The unitary time-evolution of this density matrix is
performed through the S-matrix operator [38, 33] and the resulting filtered
state is

$\rho(t=+\infty)=\frac{1}{N}\sum_{f,\tilde{f}}|f\rangle\langle
f|S\rho(t_{0}=-\infty)S^{\dagger}|\tilde{f}\rangle\langle\tilde{f}|=\frac{1}{N}\sum_{f,\tilde{f}}|f\rangle
S_{fi}S_{\tilde{f}i}^{*}\langle\tilde{f}|$ (2.3)

where the normalization factor $N$ is fixed by the condition
Tr$[\rho(t=+\infty)]=1$. From now on, we will write $\rho$ instead of
$\rho(t=+\infty)$ since entanglement of the final state polarizations is
computed here.

Finally, the corresponding bipartite density matrix elements are

$\langle
s_{1}\,s_{2}|\rho|\tilde{s}_{1}\,\tilde{s}_{2}\rangle=\frac{\mathcal{M}_{s_{1},s_{2}}\mathcal{M}^{\dagger}_{\tilde{s}_{1},\tilde{s}_{2}}}{|\overline{\mathcal{M}}|^{2}}$
(2.4)

where $\mathcal{M}_{s_{1},s_{2}}$ is the averaged amplitude
$\mathcal{M}(V^{\prime}_{1}V^{\prime}_{2}\to V_{1}V_{2})$ over all initial
state polarizations $s^{\prime}_{i}$. Also, $|\overline{\mathcal{M}}|^{2}$ is
the total unpolarized square amplitude summing over final state which accounts
for the normalization factor $N$ in Eq. (2.3). It is important to stress that
these amplitudes depends on two kinematical variables, the energy and the
scattering angle, then the $\rho$ matrix does too.

In order to reconstruct the density matrix of each scattering process from
experimental data, i.e. the goal of quantum tomography, it is useful to
introduce the following parametrization [16]

$\rho=\frac{1}{d_{1}d_{2}}I_{d_{1}d_{2}}+\frac{1}{2d_{2}}\sum_{i=1}^{d_{1}^{2}-1}A_{i}\lambda_{i}^{(d_{1})}\otimes
I_{d_{2}}+\frac{1}{2d_{1}}\sum_{j=1}^{d_{2}^{2}-1}B_{j}I_{d_{1}}\otimes\lambda_{j}^{(d_{2})}+\frac{1}{4}\sum_{i=1}^{d_{1}^{2}-1}\sum_{j=1}^{d_{2}^{2}-1}C_{ij}\lambda_{i}^{(d_{1})}\otimes\lambda_{j}^{(d_{2})}$
(2.5)

where $I_{n}$ is the dimension-$n$ identity matrix and the $\lambda^{(d)}$’s
are the dimension-$d$ generalized Gell-Mann matrices. The relevant ones for
this work are the three Pauli matrices
$\vec{\sigma}=\left(\lambda^{(2)}_{1},\lambda^{(2)}_{2},\lambda^{(2)}_{3}\right)$
for the qubits and the eight Gell-Mann matrices
$\vec{\lambda}=\left(\lambda^{(3)}_{1},...\,,\lambda^{(3)}_{8}\right)$ for the
qutrits (see Appendix B). The coefficients $A_{i}={\rm
Tr}\big{[}\rho\cdot\lambda_{i}^{(d_{1})}\otimes I_{d_{2}}\big{]}$ and
$B_{j}={\rm Tr}\big{[}\rho\cdot I_{d_{1}}\otimes\lambda_{j}^{(d_{2})}\big{]}$
provide the polarization of the gauge bosons. Moreover, the matrix
$C_{ij}={\rm
Tr}\big{[}\rho\cdot\lambda_{i}^{(d_{1})}\otimes\lambda_{j}^{(d_{2})}\big{]}$
represents the correlation among them.

###  2.1 Entanglement vs separability

The knowledge of the density matrix allows the computation of different
entanglement quantifiers in order to determine the level of the correlations
in the corresponding system. Furthermore, since entanglement is a necessary
but not sufficient condition for the violation of Bell inequalities, it is
pertinent to test both phenomena in an energy regime never explored before. In
particular, both concepts impose restrictions on $C_{ij}$, i.e. assess
particular relations between coefficients of the correlation matrix. A
‘separable’ bipartite state can be expressed as

$\rho_{\rm sep}=\sum_{n}p_{n}\rho_{n}^{(V_{1})}\otimes\rho_{n}^{(V_{2})}$
(2.6)

where $\rho_{n}^{(V_{1})}$ and $\rho_{n}^{(V_{2})}$ are the density matrices
of the subsystems $V_{1}$ and $V_{2}$, respectively, and $p_{n}$ are classical
probabilities adding 1. These separable states can be created by classical
communications and local operations. On the contrary, a bipartite state is
defined as ‘entangled’ if it is not separable, then its density matrix cannot
be written as in Eq. (2.6). In practice it is hard to decide if a given state
is separable or entangled based on the previous deﬁnition and this is called
the ‘separability problem’. There are no known necessary and sufficient
conditions for evaluating quantum entanglement of a bipartite system with
arbitrary dimension [1]. The VBS processes allow us to explore different cases
corresponding to two-qubits ($2\otimes 2$ dimension), qutrit$\otimes$qubit
($3\otimes 2$) and two-qutrits ($3\otimes 3$).

For two-qubits and qutrit$\otimes$qubit cases, the Positivity of the Partial
Transpose (PPT), also called Peres-Horodecki criterion, gives the necessary
and sufficient conditions for entanglement [39, 40]. Concretely, the partially
transpose matrix given by

$\langle s_{1}\,s_{2}|\rho^{\rm
T_{2}}|\tilde{s}_{1}\,\tilde{s}_{2}\rangle=\langle
s_{1}\,\tilde{s}_{2}|\rho|\tilde{s}_{1}\,s_{2}\rangle$ (2.7)

has at least one negative eigenvalue if and only if the density matrix in Eq.
(2.4) represents an entangled state. Denoting by $\lambda^{\rm T_{2}}_{k}$ the
eigenvalues of $\rho^{\rm T_{2}}$, the amount of the system entanglement is
computed by the Negativity [41]

$\mathcal{N}(\rho)=\sum_{k}\frac{|\lambda^{\rm T_{2}}_{k}|-\lambda^{\rm
T_{2}}_{k}}{2}$ (2.8)

Notice that the Negativity vanishes if and only if all the eigenvalues of the
partially transpose density matrix $\rho^{\rm T_{2}}$ are non-negative.

The PPT criterion is just a sufficient condition for entanglement in the two-
qutrits case. Only in some special cases it is also a necessary condition. For
example, it is shown in [15, 14] that PPT criterion is also a necessary
condition for the Higgs boson decays into $ZZ$ and $W^{+}W^{-}$. Furthermore,
the general two-qutrits processes considered in this work lead to density
matrices in Eq. (2.4) which correspond to pure states, i.e. $\rho$ is
idempotent and acts as a projector. For these bipartite pure states, the
Entropy of Entanglement $\mathcal{S}_{\rm EE}$ [42], that is the von Neumann
Entropy of either subsystem $V_{1}$ or $V_{2}$, can be computed as

$\displaystyle\mathcal{S}_{\rm EE}(\rho)$ $\displaystyle=$ $\displaystyle-{\rm
Tr}[\rho_{red\,1}\log\rho_{red\,1}]=-{\rm Tr}[\rho_{red\,2}\log\rho_{red\,2}]$
(2.9) $\displaystyle=$
$\displaystyle-\sum_{l}\lambda^{red}_{l}\log(\lambda^{red}_{l})$

where the partial trace over each subsystem yields to the reduced density
matrices $\rho_{red\,1}={\rm Tr}_{2}[\rho]=\sum_{s_{2}}\langle
s_{2}|\rho|s_{2}\rangle$ and $\rho_{red\,2}={\rm
Tr}_{1}[\rho]=\sum_{s_{1}}\langle s_{1}|\rho|s_{1}\rangle$. Both reduced
matrices have the same non-vanishing eigenvalues $\lambda^{red}_{l}$ and the
entropy of entanglement can be computed by the sum over the non-zero
eigenvalues of any reduced matrix, as in the second line of the previous
equation. The Entropy of Entanglement vanishes if and only if the pure state
is separable.

Another relevant entanglement quantifier is the ‘Concurrence’ of the system,
which for a two-qutrits pure state is computed as [43]

$\mathcal{C}[\rho]=\sqrt{2(1-{\rm Tr}[(\rho_{red\,1})^{2}])}=\sqrt{2(1-{\rm
Tr}[(\rho_{red\,2})^{2}])}$ (2.10)

which is a generalization of the two-qubit Concurrence case [44]. Therefore
the corresponding state is separable if and only if the Concurrence is equal
to zero.

The Maximally Entangled two-qubits ($d=2$) and two-qutrits ($d=3$) states are
defined as

$|\Psi_{\rm
MaxEnt}\rangle=\frac{1}{\sqrt{d}}\sum_{k=1}^{d}|k\rangle\otimes|k\rangle$
(2.11)

where $|k\rangle$ correspond to the orthonormal basis of the Hilbert spaces
$\mathcal{H}_{1,2}$. For these particular states, the Negativity, Entropy of
Entanglement and Concurrence reach their theoretical maximum values
$\frac{d-1}{2}$, $\log(d)$ and $\sqrt{\frac{2(d-1)}{d}}$, respectively [45].

###  2.2 Testing Bell inequalities

Besides the entanglement due to correlations in a quantum system, a stronger
requirement is the violation of Bell inequalities. The quantifier
$\mathcal{I}$ is introduced in order to discriminate among predictions coming
from deterministic local theories from those coming from quantum mechanics. In
particular, local realist models must obey the so-called Bell inequality

$\mathcal{I}\leq 2$ (2.12)

which can be violated in QFT. For the two-qubits cases, the optimal Clauser-
Horne-Shimony-Holt (CHSH) operator [46] is defined as

$\mathcal{B}_{CHSH}=\hat{a}_{1}\otimes\hat{b}_{1}-\hat{a}_{1}\otimes\hat{b}_{2}+\hat{a}_{2}\otimes\hat{b}_{1}+\hat{a}_{2}\otimes\hat{b}_{2}$
(2.13)

with $\hat{a}_{i}=\vec{a}_{i}\cdot\vec{\sigma}$ and
$\hat{b}_{i}=\vec{b}_{i}\cdot\vec{\sigma}$ are Hermitian operators acting on
the Hilbert spaces $\mathcal{H}_{1}$ and $\mathcal{H}_{2}$, respectively. In
practical computations, the quantifier $\mathcal{I}_{2\otimes
2}\equiv\mathcal{I}_{2}$ is defined as the maximal expectation value of this
Bell operator over the unit vectors $\vec{a}_{i}$ and $\vec{b}_{i}$ in
$\mathbb{R}^{3}$

$\displaystyle\mathcal{I}_{2}$ $\displaystyle=$
$\displaystyle\stackbin[\vec{a}_{i},\vec{b}_{i}]{}{\rm Max}\big{\\{}{\rm
Tr}[\rho\cdot\mathcal{B}_{CHSH}]\big{\\}}$ (2.14) $\displaystyle=$
$\displaystyle\stackbin[\vec{a}_{i},\vec{b}_{i}]{}{\rm
Max}\big{\\{}\vec{a}_{1}^{\rm T}C(\vec{b}_{1}-\vec{b}_{2})+\vec{a}_{2}^{\rm
T}C(\vec{b}_{1}+\vec{b}_{2})\big{\\}}$ $\displaystyle=$ $\displaystyle
2\sqrt{r_{1}+r_{2}}$

where the decomposition of the density matrix in Eq. (2.5) was used in the
second line (notice that only the correlation matrix $C$ is relevant for this
expectation value). In the third line, $r_{1}$ and $r_{2}$ are the two largest
eigenvalues of the matrix $C^{\rm T}C$ [47]. Hence the Bell inequality in Eq.
(2.12) can be violated if and only if $r_{1}+r_{2}$ is larger than 1.
Moreover, the maximum theoretical value is the Cirelson bound $2\sqrt{2}$
[48], corresponding to the Maximally Entangled state.

For the qutrit$\otimes$qubit cases, and in a similar way than in [49], we
consider a generalization to the previous CHSH. Now the 3$\otimes$2 Bell
operator is

$\mathcal{B}_{CHSH}^{gen}=\vec{n}_{1}\cdot\vec{S}\otimes(\vec{n}_{2}-\vec{n}_{4})\cdot\vec{\sigma}+\vec{n}_{3}\cdot\vec{S}\otimes(\vec{n}_{2}+\vec{n}_{4})\cdot\vec{\sigma}$
(2.15)

where $\vec{n}_{i}$ are unit vectors in $\mathbb{R}^{3}$ and the dimension-3
spin-1 matrices $\vec{S}=(S_{x},S_{y},S_{z})$ are shown in Eq. (B.4). The
resulting quantifier is

$\displaystyle\mathcal{I}_{3\otimes 2}$ $\displaystyle=$
$\displaystyle\stackbin[\vec{n}_{i}]{}{\rm Max}\big{\\{}{\rm
Tr}[\rho\cdot\mathcal{B}_{CHSH}^{gen}]\big{\\}}$ (2.16) $\displaystyle=$
$\displaystyle 2\sqrt{\tilde{r}_{1}+\tilde{r}_{2}}$

where, following the two-qubits case, $\tilde{r}_{1,2}$ are the two largest
eigenvalues of $\tilde{C}^{\rm T}\tilde{C}$. In this case, the spin
correlation matrix $\tilde{C}$ is defined from the decomposition of the
density matrix in Eq. (2.5) as

$\tilde{C}_{1j}=\frac{1}{\sqrt{2}}(C_{1j}+C_{6j})\,,\qquad\tilde{C}_{2j}=\frac{1}{\sqrt{2}}(C_{2j}+C_{7j})\,,\qquad\tilde{C}_{3j}=\frac{1}{2}(C_{3j}+\sqrt{3}C_{8j})$
(2.17)

For the two-qutrits cases, the optimal quantifier $\mathcal{I}_{3\otimes
3}\equiv\mathcal{I}_{3}$ correspond to Collins-Gisin-Linden-Massar-Popescu
(CGLMP) [50]. This $\mathcal{I}_{3}$ is constructed as the expectation value
of an appropriate Bell operator $\hat{\mathcal{B}}_{CGLMP}$ as usual, ${\rm
Tr}[\rho\cdot\hat{\mathcal{B}}_{CGLMP}]$. A clever election of the
$\hat{\mathcal{B}}_{CGLMP}$ operator improves the violation of Bell inequality
in Eq. (2.12). It was shown in [13] that a suitable operator for Higgs boson
decaying into gauge boson pair corresponds to

$\hat{\mathcal{B}}_{CGLMP}^{xy}=-\frac{2}{\sqrt{3}}\left(S_{x}\otimes
S_{x}+S_{y}\otimes
S_{y}\right)+\lambda_{4}\otimes\lambda_{4}+\lambda_{5}\otimes\lambda_{5}$
(2.18)

where the spin and the Gell-Mann matrices in dimension-3 are collected in the
Appendix B. In turn, local unitary changes of basis of the states
$|s_{1}\rangle\otimes|s_{2}\rangle$, which define the density matrix, allow an
optimization of the $\mathcal{I}_{3}$ values for the VBS processes considered
here. In particular, a maximization procedure for each $\rho$ density matrix
is performed

$\mathcal{I}_{3}=\stackbin[\vec{\alpha}_{i},\vec{\beta}_{i}]{}{\rm
Max}\big{\\{}{\rm Tr}[\rho\cdot(U_{1}\otimes
U_{2})^{\dagger}\cdot\hat{\mathcal{B}}_{CGLMP}^{xy}\cdot(U_{1}\otimes
U_{2})]\big{\\}}$ (2.19)

where the maximization is over the domain $0\leq\alpha_{1},\,\alpha_{2}\leq
2\pi$ and $0\leq\beta_{1},\,\beta_{2}\leq\pi$ for the angles defining the
dimension-3 rotation matrices

$U_{k}=\exp(-iS_{z}\alpha_{k})\exp(-iS_{y}\beta_{k})$ (2.20)

The maximum for $\mathcal{I}_{3}$ is $1+\sqrt{11/3}\approx\,2.915$ [51] and
notoriously it is not achieved for the Maximally Entangled state which has
$\mathcal{I}_{3}^{\rm MaxEnt}=4(6\sqrt{3}+9)/27\approx\,2.873$.

In summary, the quantifiers related to entanglement detection and test Bell
inequalities require the full knowledge of the quantum state which is
collected in the density matrix. This $\rho$ matrix is determined using Eq.
(2.4) by computing the scattering amplitudes for the corresponding VBS
process. With this analysis, we can explore kinematical regions relevant for
quantum mechanical measurements at colliders. Depending on the final state,
the following quantifiers will be derived,

  * •

for the two-qubits case: the only VBS process of this kind is
$W^{+}W^{-}\to\gamma\gamma$ and its Negativity $\mathcal{N}$ of Eq. (2.8), the
Entropy of Entanglement $\mathcal{S}_{\rm EE}$ of Eq. (2.9) and
$\mathcal{I}_{2}$ of Eq. (2.14) will be determined.

  * •

for the qutrit$\otimes$qubit case: also $\mathcal{N}$, $\mathcal{S}_{\rm EE}$,
and $\mathcal{I}_{3\otimes 2}$ of Eq. (2.16) will be calculated. The
considered VBS processes of this kind are $W^{\pm}\gamma\to W^{\pm}\gamma$,
$W^{+}W^{-}\to Z\gamma$ and $W^{\pm}Z\to W^{\pm}\gamma$.

  * •

for the two-qutrits case: the Entropy of Entanglement $\mathcal{S}_{\rm EE}$,
Concurrence $\mathcal{C}$ and $\mathcal{I}_{3}$ of Eq. (2.9), Eq. (2.10) and
Eq. (2.19) will be computed, respectively. In this case, the considered VBS
processes are $\gamma\gamma\to W^{+}W^{-}$, $W^{\pm}\gamma\to W^{\pm}Z$,
$W^{+}W^{-}\to W^{+}W^{-}$, $W^{\pm}W^{\pm}\to W^{\pm}W^{\pm}$, $W^{+}W^{-}\to
ZZ$, $W^{\pm}Z\to W^{\pm}Z$, $Z\gamma\to W^{+}W^{-}$, $ZZ\to W^{+}W^{-}$ and
$ZZ\to ZZ$.

##  3 Numerical Results

The SM $2\to 2$ scattering amplitudes are computed at tree level in the
center-of-mass frame (see Appendix A for details) using FeynArts [52] and
FormCalc [53]. Some comments are in order: firstly, the radiative corrections
can be relevant near the extreme theoretical values of the entanglement
quantifiers in order to decide if the final state is separable or maximally
entangled. Also for testing Bell inequality when $\mathcal{I}\sim 2$, however
the generation of entanglement at the lowest order in perturbation theory is
interesting by itself and we work on it. Secondly, the amplitudes depend on
two kinematical quantities corresponding to the scattering angle $\theta$ and
the Mandelstam variable $S$ (related to the invariant mass of the final gauge
bosons). The numerical predictions of the mentioned entanglement quantifiers
were computed in the plane $[\cos(\theta),\sqrt{S}]$ in the range -1 to 1 for
the cosine function and from the scattering threshold up to 3 TeV for the
energy.

###  3.1 two-qubits case

The only tree level VBS with photon-pair in the final state is
$W^{+}W^{-}\to\gamma\gamma$. This lowest dimension case can be treated
analytically and clarify the density matrix formalism from scattering
amplitudes introduced in the previous section. For a given photon-pair
polarization $|s_{1}\,s_{2}\rangle$, the averaged amplitude over the initial
state polarizations is

$\displaystyle\mathcal{M}_{s_{1},s_{2}}$ $\displaystyle=$
$\displaystyle\frac{1}{9}\left(\mathcal{M}(W^{+}_{+}W^{-}_{+}\to\gamma_{s_{1}}\gamma_{s_{2}})+\mathcal{M}(W^{+}_{+}W^{-}_{L}\to\gamma_{s_{1}}\gamma_{s_{2}})+\mathcal{M}(W^{+}_{+}W^{-}_{-}\to\gamma_{s_{1}}\gamma_{s_{2}})\right.$
(3.1)
$\displaystyle\left.+\mathcal{M}(W^{+}_{L}W^{-}_{+}\to\gamma_{s_{1}}\gamma_{s_{2}})+\mathcal{M}(W^{+}_{L}W^{-}_{L}\to\gamma_{s_{1}}\gamma_{s_{2}})+\mathcal{M}(W^{+}_{L}W^{-}_{-}\to\gamma_{s_{1}}\gamma_{s_{2}})\right.$
$\displaystyle\left.+\mathcal{M}(W^{+}_{-}W^{-}_{+}\to\gamma_{s_{1}}\gamma_{s_{2}})+\mathcal{M}(W^{+}_{-}W^{-}_{L}\to\gamma_{s_{1}}\gamma_{s_{2}})+\mathcal{M}(W^{+}_{-}W^{-}_{-}\to\gamma_{s_{1}}\gamma_{s_{2}})\right)$

where the sub-indices stand for the polarizations. The explicit computation of
these amplitudes in terms of $S$ and $c=\cos(\theta)$ result in

$\displaystyle\mathcal{M}_{\pm,\pm}$ $\displaystyle=$
$\displaystyle\frac{1}{9}\frac{2e^{2}}{S(1-c^{2})+4c^{2}m_{\rm
W}^{2}}\left(-12\left(1-c^{2}\right)m_{\rm
W}^{2}+8\sqrt{2}ic\sqrt{1-c^{2}}m_{\rm
W}\sqrt{S}+\left(1+3c^{2}\right)S\right)$ $\displaystyle\mathcal{M}_{\pm,\mp}$
$\displaystyle=$
$\displaystyle\frac{1}{9}\frac{8e^{2}}{S(1-c^{2})+4c^{2}m_{\rm
W}^{2}}\left(S-3m_{\rm W}^{2}\right)$ (3.2)

and the corresponding total unpolarized square amplitude is

$\displaystyle|\overline{\mathcal{M}}|^{2}$ $\displaystyle=$
$\displaystyle|\mathcal{M}_{++}|^{2}+|\mathcal{M}_{+-}|^{2}+|\mathcal{M}_{-+}|^{2}+|\mathcal{M}_{--}|^{2}$
(3.3) $\displaystyle=$
$\displaystyle\left(\frac{1}{9}\frac{2e^{2}}{S(1-c^{2})+4c^{2}m_{\rm
W}^{2}}\right)^{2}D_{\gamma\gamma}$

with

$D_{\gamma\gamma}=2\left(144\left(c^{4}-2c^{2}+2\right)m_{\rm
W}^{4}-8\left(7c^{4}-10c^{2}+15\right)m_{\rm
W}^{2}S+\left(9c^{4}+6c^{2}+17\right)S^{2}\right)$ (3.4)

Therefore in the basis $\\{|++\rangle,|+-\rangle,|-+\rangle,|--\rangle\\}$,
the density matrix $\rho$ and its partially transpose $\rho^{\rm T_{2}}$ have
a compact form

$\rho=\begin{pmatrix}\rho_{11}&\rho_{12}&\rho_{12}&\rho_{11}\\\
\rho_{12}^{*}&\rho_{22}&\rho_{22}&\rho_{12}^{*}\\\
\rho_{12}^{*}&\rho_{22}&\rho_{22}&\rho_{12}^{*}\\\
\rho_{11}&\rho_{12}&\rho_{12}&\rho_{11}\\\ \end{pmatrix}\qquad{\rm
and}\qquad\rho^{\rm
T_{2}}=\begin{pmatrix}\rho_{11}&\rho_{12}^{*}&\rho_{12}&\rho_{22}\\\
\rho_{12}&\rho_{22}&\rho_{11}&\rho_{12}^{*}\\\
\rho_{12}^{*}&\rho_{11}&\rho_{22}&\rho_{12}\\\
\rho_{22}&\rho_{12}&\rho_{12}^{*}&\rho_{11}\\\ \end{pmatrix}$ (3.5)

where the three independent entries are:

$\displaystyle\rho_{11}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{\gamma\gamma}}\left(144\left(1-c^{2}\right)^{2}m_{\rm
W}^{4}-8\left(7c^{4}-10c^{2}+3\right)m_{\rm
W}^{2}S+\left(1+3c^{2}\right)^{2}S^{2}\right)$ $\displaystyle\rho_{12}$
$\displaystyle=$ $\displaystyle\frac{4(S-3m_{\rm
W}^{2})}{D_{\gamma\gamma}}\left(-12\left(1-c^{2}\right)m_{\rm
W}^{2}+8i\sqrt{2}c\sqrt{1-c^{2}}m_{\rm W}\sqrt{S}+(1+3c^{2})S\right)$
$\displaystyle\rho_{22}$ $\displaystyle=$
$\displaystyle\frac{16}{D_{\gamma\gamma}}\left(S-3m_{\rm W}^{2}\right)^{2}$
(3.6)

The eigenvalues of the compact $\rho^{\rm T_{2}}$ in Eq. (3.5) can be written
as

$\displaystyle\lambda^{\rm T_{2}}_{1}$ $\displaystyle=$
$\displaystyle\rho_{11}+\rho_{22}-\rho_{12}-\rho_{12}^{*}\,,\qquad\lambda^{\rm
T_{2}}_{2}=\rho_{11}+\rho_{22}+\rho_{12}+\rho_{12}^{*}\,,$
$\displaystyle\lambda^{\rm T_{2}}_{3}$ $\displaystyle=$
$\displaystyle-\lambda^{\rm
T_{2}}_{4}=-\sqrt{(\rho_{11}-\rho_{22})^{2}-(\rho_{12}-\rho_{12}^{*})^{2}}$
(3.7)

Using Eq. (3.6), the previous expressions for the $W^{+}W^{-}\to\gamma\gamma$
scattering are

$\lambda^{\rm
T_{2}}_{1}=\frac{U_{1}^{\gamma\gamma}}{D_{\gamma\gamma}}\,,\quad\lambda^{\rm
T_{2}}_{2}=\frac{U_{2}^{\gamma\gamma}}{D_{\gamma\gamma}}\,,\quad\lambda^{\rm
T_{2}}_{3}=-\frac{\sqrt{U_{1}^{\gamma\gamma}U_{2}^{\gamma\gamma}}}{D_{\gamma\gamma}}\,,\quad\lambda^{\rm
T_{2}}_{4}=\frac{\sqrt{U_{1}^{\gamma\gamma}U_{2}^{\gamma\gamma}}}{D_{\gamma\gamma}}\\\
$

in terms of the functions

$\displaystyle U_{1}^{\gamma\gamma}$ $\displaystyle=$ $\displaystyle
144c^{4}m_{\rm W}^{4}+56c^{2}\left(1-c^{2}\right)m_{\rm
W}^{2}S+9\left(1-c^{2}\right)^{2}S^{2}$ $\displaystyle U_{2}^{\gamma\gamma}$
$\displaystyle=$ $\displaystyle 144(2-c^{2})^{2}m_{\rm
W}^{4}-8(30-13c^{2}+7c^{4})m_{\rm W}^{2}S+(5+3c^{2})^{2}S^{2}$ (3.8)

Notice that $U_{1}^{\gamma\gamma}+U_{2}^{\gamma\gamma}=D_{\gamma\gamma}$ since
Tr[$\rho^{\rm T_{2}}$]=Tr[$\rho$]=1.

In the relevant phase space of $[\cos(\theta),\sqrt{S}]$, the functions
$D_{\gamma\gamma}$, $U_{1}^{\gamma\gamma}$ and $U_{2}^{\gamma\gamma}$ are
positives. Therefore only the third eigenvalue $\lambda^{\rm T_{2}}_{3}$ is
negative and the rest are positives. In consequence the Negativity is
analytically given by

$\mathcal{N}_{W^{+}W^{-}\to\gamma\gamma}=\frac{\sqrt{U_{1}^{\gamma\gamma}U_{2}^{\gamma\gamma}}}{U_{1}^{\gamma\gamma}+U_{2}^{\gamma\gamma}}$
(3.9)

The left panel of Fig. 1 shows the behaviour of the Negativity for this
process as a function of the energy and the scattering angle. The minima of
this quantifier is located at the upper corners, i.e., for $\cos(\theta)=\pm
1$ and large energies. From the analytical expression of Eq. (3.9), in this
regime the Negativity behaves as $\sim 1.5m_{\rm W}^{2}/S$ and reaches the
minimum value $\sim 10^{-3}$ for $\sqrt{S}=3$ TeV. In particular, it never
vanishes. On the other hand, the arithmetic and geometric means inequality
provides the maxima $1/2$ for the Negativity in Eq. (3.9) when
$U_{1}^{\gamma\gamma}=U_{2}^{\gamma\gamma}$. This maxima is the theoretical
maximum expected for the Negativity corresponding to the Maximally Entangled
pure states. The curve of Negativity equals to $1/2$ is given by

$S|_{\rm MaxEnt}=\frac{12(1-\cos^{2}(\theta))m_{\rm
W}^{2}}{1+3\cos^{2}(\theta)}$ (3.10)

and it is represented by the solid black line in the red region of Fig. 1.

|
---|---
Figure 1: Negativity (left) and $\mathcal{I}_{2}$ quantifier (right) for
$W^{+}W^{-}\to\gamma\gamma$ in the plane $[\cos(\theta),\sqrt{S}]$. Dashed
contour lines are shown for an easy comparison of the numerical values. The
solid contour line corresponds to the maximal Negativity equals to $1/2$
associated to Maximally Entangled states of Eq. (3.10).

Furthermore, the Entropy of Entanglement of Eq. (2.9) can also be computed
analytically. The reduced matrices respect to each photon coincide and are
given by

$\rho_{red}=\begin{pmatrix}\frac{1}{2}&\frac{8}{D_{\gamma\gamma}}(S-3m_{\rm
W}^{2})((1+3c^{2})S-12(1-c^{2})m_{\rm W}^{2})\\\
\frac{8}{D_{\gamma\gamma}}(S-3m_{\rm W}^{2})((1+3c^{2})S-12(1-c^{2})m_{\rm
W}^{2})&\frac{1}{2}\end{pmatrix}$ (3.11)

with eigenvalues equal to

$\lambda^{red}_{1}=\frac{U^{\gamma\gamma}_{1}}{D_{\gamma\gamma}}\quad{\rm
and}\quad\lambda^{red}_{2}=\frac{U^{\gamma\gamma}_{2}}{D_{\gamma\gamma}}$
(3.12)

Therefore, the resulting Entropy of Entanglement is

$\mathcal{S}_{\rm
EE}=\frac{U^{\gamma\gamma}_{1}}{U_{1}^{\gamma\gamma}+U_{2}^{\gamma\gamma}}\log\left(1+\frac{U^{\gamma\gamma}_{2}}{U^{\gamma\gamma}_{1}}\right)+\frac{U^{\gamma\gamma}_{2}}{U_{1}^{\gamma\gamma}+U_{2}^{\gamma\gamma}}\log\left(1+\frac{U^{\gamma\gamma}_{1}}{U^{\gamma\gamma}_{2}}\right)$
(3.13)

The behaviour of this quantifier in the plane $[\cos(\theta),\sqrt{S}]$ is
very similar to the Negativity and the resulting plot is relegated to the
upper-left corner of Fig. 9 in Appendix C for saving space here. The minima
are located for $\cos(\theta)=\pm 1$ and large energies decreasing as
$\sim(2.25m_{\rm W}^{4}/S^{2})(1-\log(2.25m_{\rm W}^{4}/S^{2}))$. On the other
hand, Eq. (3.13) for $U^{\gamma\gamma}_{1}=U^{\gamma\gamma}_{2}$ reaches the
maximum theoretical value $\log(2)$ corresponding to the Maximally Entangled
state described by the curve in Eq. (3.10) for which the reduced matrix is the
half of the identity matrix.

Finally, the decomposition of Eq. (2.5) corresponding to this density matrix
have the non-vanishing coefficients

$\displaystyle A_{1}=B_{1}=\frac{16}{D_{\gamma\gamma}}\left(S-3m_{\rm
W}^{2}\right)\left((1+3c^{2})S-12(1-c^{2})m_{\rm W}^{2}\right)\,,$
$\displaystyle C_{11}=1\,,\qquad
C_{23}=C_{32}=-\frac{128\sqrt{2}}{D_{\gamma\gamma}}\left(c\sqrt{1-c^{2}}(S-3m_{\rm
W}^{2})m_{\rm W}\sqrt{S}\right)\,,$ $\displaystyle
C_{22}=-C_{33}=\frac{2}{D_{\gamma\gamma}}\left((15-6c^{2}-9c^{4})S^{2}-8(9+10c^{2}-7c^{4})m_{\rm
W}^{2}S+144c^{2}(2-c^{2})m_{\rm W}^{4}\right)$ (3.14)

Hence the two largest eigenvalues of $C^{\rm T}C$ are

$r_{1}=1\quad{\rm and}\quad
r_{2}=1-\frac{256}{D_{\gamma\gamma}^{2}}\left(S-3m_{\rm
W}^{2}\right)^{2}\left((1+3c^{2})S-12(1-c^{2})m_{\rm W}^{2}\right)^{2}$ (3.15)

resulting in the $\mathcal{I}_{2}$ quantifier of Eq. (2.14) for this process
as

$\mathcal{I}_{2}=2\sqrt{2-\frac{256}{D_{\gamma\gamma}^{2}}\left(S-3m_{\rm
W}^{2}\right)^{2}\left((1+3c^{2})S-12(1-c^{2})m_{\rm W}^{2}\right)^{2}}$
(3.16)

The behaviour of this function in the plane $[\cos(\theta),\sqrt{S}]$ is shown
in the right panel of Fig. 1. The minima correspond to large energies in the
direction $\cos(\theta)=\pm 1$ and can be written as $\sim 2+9m_{\rm
W}^{4}/S^{2}$. In particular, this quantifier is always greater than 2
signaling a violation of Bell inequality in the whole kinematical plane for
this observable. Conversely, for the Maximally Entangled state described by
the curve in Eq. (3.10), this quantifier reaches the theoretical maximum
$2\sqrt{2}$ corresponding to the Cirelson bound [48]. These remarks show that
$W^{+}W^{-}\to\gamma\gamma$ could be an ideal laboratory for a Bell inequality
test among the VBS processes but it requires polarization measurements of the
final photons111Although the polarization of high-energy photons is not
currently measure in ATLAS and CMS, in contrast to the case of massive gauge
bosons, the LHCb Collaboration performed analysis for photon polarization in
$b$-baryon decays [54]. There are also proposals to study CP properties of the
Higgs boson through the di-photon decay [55, 56]. (see for instance a related
discussion in [30]).

|
---|---
|
|
Figure 2: Negativity (left) and $\mathcal{I}_{3\otimes 2}$ quantifier (right)
for $W^{\pm}\gamma\to W^{\pm}\gamma$ (first row), $W^{+}W^{-}\to Z\gamma$
(second row) and $W^{\pm}Z\to W^{\pm}\gamma$ (third row) in the plane
$[\cos(\theta),\sqrt{S}]$. Contour lines are shown for an easy comparison of
the numerical values.

###  3.2 qutrit$\otimes$qubit case

The considered VBS processes of this kind are $W^{\pm}\gamma\to
W^{\pm}\gamma$, $W^{+}W^{-}\to Z\gamma$ and $W^{\pm}Z\to W^{\pm}\gamma$. The
separability of these 3$\otimes$2 final states is also determined by the PPT
criterion as in the previous section. Just the $W^{\pm}\gamma\to
W^{\pm}\gamma$ process was treated analytically and the coefficients of the
decomposition in Eq. (2.5) are collected in Appendix D. Remember the relevance
of these coefficients due to their relation with the experimental data through
quantum tomography. In addition, the analytical expression corresponding to
the Entropy of Entanglement for this process is shown in this appendix.

The corresponding Negativity for $W^{\pm}\gamma\to W^{\pm}\gamma$ is presented
in the left of the first row of Fig. 2. In this process, the Negativity
vanishes at the points $[\cos(\theta),m_{\rm W}]$ and $[1,\sqrt{S}]$ of the
plane $[\cos(\theta),\sqrt{S}]$. For both kinematical regions, the six
eigenvalues of the partial transpose density matrix are $\\{1,0,0,0,0,0\\}$.
In consequence, the final state at the threshold and in the forward direction
is separable. This fact can be also understood with the analytical expression
of $\mathcal{S}_{\rm EE}$ in Eq. (D.5). In turn, the maximum value of the
Negativity is $\sim$0.345 which indicates that the Maximally Entangled state
never occurs in these scattering. On the other hand, the Negativity for
$W^{+}W^{-}\to Z\gamma$ and $W^{\pm}Z\to W^{\pm}\gamma$ is shown in the left
of the second and third rows of Fig. 2. Similar to the two-qubits case, the
minima are located at the upper corners with values $\sim$8$\cdot 10^{-4}$.
Also, the theoretical maximum $1/2$ is achieved for some points in the red
region but the corresponding analytical curve cannot be computed.

The Entropy of Entanglement for these processes follows the same pattern of
the Negativity, then corresponding plots are omitted here since no additional
conclusions can be extracted. They are relegated to the Fig. 9 in Appendix C
for saving space here.

Regarding the Bell inequality, the quantifier $\mathcal{I}_{3\otimes 2}$ of
Eq. (2.16) is shown in the right column of Fig. 2 for each process. This
quantifier never exceeds 2 over the whole kinematical plane. Concretely,
$\mathcal{I}_{3\otimes 2}$ varies in the range $\sim(3\cdot 10^{-3},1.88)$ for
$W^{\pm}\gamma\to W^{\pm}\gamma$, in the range $\sim(1.38,2)$ for
$W^{+}W^{-}\to Z\gamma$ and $\sim(0.94,2)$ for $W^{\pm}Z\to W^{\pm}\gamma$.
The maximum value 2 for the last two processes is reached for the points with
maximal Nagativity. Some comments are in order: firstly, the plots show that
the entanglement is not a sufficient condition for Bell inequality violation.
In particular, we have non-vanishing values of the Negativity but the
corresponding points never exceeds 2 for the $\mathcal{I}_{3\otimes 2}$
quantifier. Secondly, these results are expected [13, 49] since the
generalized CHSH operator of Eq. (2.15) is diminished by the vanishing outcome
of the spin operator for massive gauge bosons. As far as we know, there is no
optimization for Bell operator in the $3\otimes 2$ case and it deserves a
further study which also applied to single-top processes.

###  3.3 two-qutrits case

|
---|---
|

Figure 3: Entropy of Entanglement for different VBS processes in the plane
$[\cos(\theta),\sqrt{S}]$. Contour lines are shown for an easy comparison of
the numerical values.

The VBS processes corresponding to bipartite 3$\otimes$3 states have massive
gauge bosons in the final state. Concretely, they are $\gamma\gamma\to
W^{+}W^{-}$, $W^{\pm}\gamma\to W^{\pm}Z$, $W^{+}W^{-}\to W^{+}W^{-}$,
$W^{\pm}W^{\pm}\to W^{\pm}W^{\pm}$, $W^{+}W^{-}\to ZZ$, $W^{\pm}Z\to
W^{\pm}Z$, $Z\gamma\to W^{+}W^{-}$, $ZZ\to W^{+}W^{-}$ and $ZZ\to ZZ$. As
presented in Section 2, the level of entanglement of these pure states is
given by non-vanishing values of the Entropy of Entanglement in Eq. (2.9) and
Concurrence in Eq. (2.10). The analytical expressions for the coefficients of
the decomposition in Eq. (2.5) and the Entropy of Entanglement are given in
Appendix D just for $\gamma\gamma\to W^{+}W^{-}$ process.

The Fig. 3 collects the $\mathcal{S}_{\rm EE}$ for $\gamma\gamma\to
W^{+}W^{-}$, $W^{+}W^{-}\to ZZ$, $Z\gamma\to W^{+}W^{-}$, $ZZ\to W^{+}W^{-}$
and $W^{\pm}\gamma\to W^{\pm}Z$. The green regions match with values lower
than $\sim$0.1 and are located in both directions $\cos(\theta)\sim\pm 1$. As
the energy increases, this entropy diminishes reaching values $\sim 10^{-5}$
but never vanishes. The red regions in the direction $\cos(\theta)\sim 0$ and
near the threshold correspond to the maximal entropy values. In general, they
are between 0.7 and 0.8 except for the $ZZ\to W^{+}W^{-}$ (right panel of
second row) having 1.04 which is closer to the theoretical maximum $\log(3)$
corresponding to the Maximally Entangled state.

|
---|---
|
Figure 4: Entropy of Entanglement for the rest of VBS processes in the plane
$[\cos(\theta),\sqrt{S}]$. Contour lines are shown for an easy comparison of
the numerical values.

The Fig. 4 gathers the rest of the VBS processes $W^{+}W^{-}\to W^{+}W^{-}$,
$W^{\pm}W^{\pm}\to W^{\pm}W^{\pm}$, $W^{\pm}Z\to W^{\pm}Z$ and $ZZ\to ZZ$. As
before, none of these processes have vanishing entropy. The first two (showed
in the first row) have the minima at the lower corners, i.e. near the
threshold and not at large energies as for processes of Fig. 3. The maxima are
also near the threshold but in the direction $\cos(\theta)\sim 0$. The
$W^{\pm}Z\to W^{\pm}Z$ (left panel of the second row) exhibits a strong
asymmetry in the scattering angle, with lower values of the entropy for
$\cos(\theta)\sim-1$ and minima located at large energies. On the contrary,
the maxima is $\sim 1$ situated in the lower-right corner. The later process
has values between 0.1 and 0.15 in the whole plane in contrast to the large
variations showed in the other VBS processes.

The behaviour of the Entanglement Entropy is also manifest in a very similar
way for the Concurrence in each VBS process. The corresponding plots are
relegated to Figs. 10-11 in Appendix C. As before, this quantifier never
reaches the zero value, i.e. the final states are entangled in the whole
kinematical plane. The maximal value of the Concurrence is 1.12 also for
$ZZ\to W^{+}W^{-}$ which is close to the theoretical maximum is $2/\sqrt{3}$.

|
---|---

Figure 5: $\mathcal{I}_{3}$ quantifier for $\gamma\gamma\to W^{+}W^{-}$ (left of the first row), $W^{+}W^{-}\to W^{+}W^{-}$ (right of the first row) and $ZZ\to ZZ$ (second row) in the plane $[\cos(\theta),\sqrt{S}]$. The solid black contour line corresponds to $\mathcal{I}_{3}$=2, limiting the region for Bell inequality violation. Dashed lines are shown for an easy comparison of the numerical values. |
---|---
|
|
Figure 6: $\mathcal{I}_{3}$ quantifier for the rest of VBS processes in the
plane $[\cos(\theta),\sqrt{S}]$. The solid black contour line corresponds to
$\mathcal{I}_{3}$=2, limiting the region for Bell inequality violation. Dashed
lines are shown for an easy comparison of the numerical values.

Finally, the theoretical predictions for the violation of Bell inequalities by
means of the $\mathcal{I}_{3}$ parameter in Eq. (2.19) are also presented in
the plane $[\cos(\theta),\sqrt{S}]$ for each VBS process in this two-qutrits
case. For each kinematical point, the rotation matrices $U_{1}$ and $U_{2}$
are determined in order to get maximal values for this quantifier. The first
row of Fig. 5 shows the predictions for $\gamma\gamma\to W^{+}W^{-}$ (left)
and $W^{+}W^{-}\to W^{+}W^{-}$ (right). The solid black line corresponds to
$\mathcal{I}_{3}=2$ and delimits the orange-red region having values greater
than 2. In both cases, the achieved maximum is $\sim$2.38 For energies above
300 GeV there is violation of Bell inequality when
$-0.5\lesssim\cos(\theta)\lesssim 0.5$ ($\gamma\gamma\to W^{+}W^{-}$)

and $-0.6\lesssim\cos(\theta)\lesssim 0.05$ ($W^{+}W^{-}\to W^{+}W^{-}$). On
the contrary, the $ZZ\to ZZ$ (second row of Fig. 5) never exceeds 2 in the
whole plane and reaches the maximum $\sim$1.76 in the upper corners.

For the rest of the VBS processes, the region corresponding to
$\mathcal{I}_{3}>2$ has energy below 600 GeV. Then they are collected in Fig.
6 for this low energy region in which the solid black line also delimits the
region for Bell inequality violation.

The first row contains the $ZZ\to W^{+}W^{-}$ and $W^{\pm}Z\to W^{\pm}Z$
having maximum $\sim$2.82 and $\sim$2.72, respectively. The maximum for
$W^{\pm}W^{\pm}\to W^{\pm}W^{\pm}$ and $W^{+}W^{-}\to ZZ$ (second row) are
$\sim$2.71 and $\sim$2.62. In addition, $W^{\pm}\gamma\to W^{\pm}Z$ and
$Z\gamma\to W^{+}W^{-}$ are in the last row and both have maximum $\sim$2.21
All these maxima are located in $\cos(\theta)\sim 0$ and energies between 230
GeV and 290 GeV depending on the process, except for $W^{\pm}Z\to W^{\pm}Z$
for which is located at $[\cos(\theta),\sqrt{S}]\sim$[1,210 GeV].

As discussed in the qutrit$\otimes$qubit processes, entanglement is just a
necessary but not sufficient condition for Bell inequality violation. We
showed that all the final states result entangled after the scattering process
in the whole $[\cos(\theta),\sqrt{S}]$ plane but $\mathcal{I}_{3}$ greater
than 2 is just achieved in small regions of the phase space (even worse for
$ZZ\to ZZ$ which never violates the Bell inequality).

Figure 7: Experimental signature for VBS at sub-process level in the LHC
(future electron-positron colliders). In the latter case, the companion
leptons $l_{1,2}$ can be electrons, positrons or neutrinos depending on the
gauge bosons $V^{\prime}_{1,2}$. The leptonic decays of the $V_{1,2}$ are
represented by the final solid lines.

###  3.4 Prospects at colliders

The analyzed VBS channels correspond to a sub-process at collider, yielding to
the following experimental signatures:

$\displaystyle pp$ $\displaystyle\to$ $\displaystyle(V_{1}V_{2})_{\rm
leptonic\,decay}+j_{1}+j_{2}\quad{\rm for\,\,the\,\,LHC,}$ $\displaystyle
e^{-}e^{+}$ $\displaystyle\to$ $\displaystyle(V_{1}V_{2})_{\rm
leptonic\,decay}+l_{1}+l_{2}\quad{\rm for\,\,future\,\,electron-
positron\,\,colliders,}$ (3.17)

which are generically represented in Fig. 7.

In this figure, the gray circle corresponds to the tree level VBS processes
described in Appendix A. This kind of processes are purely EW and very rare at
the LHC, however they allow to probe the core of the EW Symmetry Breaking
mechanism due to the interplay among the triple and quartic gauge self-
coupling with the Higgs boson interactions. Nowadays, measurements of fiducial
and total cross-sections as well as polarization of the outgoing gauge bosons
are performed by ATLAS and CMS Collaborations [57]. The reported results with
the 13 TeV data for all222$W^{+}W^{-}\to\gamma\gamma$ process has not been
measured yet at the LHC, due to the immense QCD multijet background [34].
However recent sensibility studies can be found in [58]. the considered final
vector boson pairs in this work, in association with two jets, are:
$W^{\pm}\gamma$ [59], $Z\gamma$ [60, 61], $W^{\pm}W^{\pm}$ [62, 63],
$W^{\pm}Z$ [64, 65, 66, 67], $W^{+}W^{-}$ [68, 69] and $ZZ$ [70, 71].

The LHC signal has a very characteristic kinematics at detector level given by
two energetic jets in the forward region with large invariant mass
$m_{j_{1}j_{2}}$ and a wide pseudo-rapidity separation
$\Delta\eta_{j_{1}j_{2}}$. A similar topology for the companion leptons in the
future electron-positron colliders is expected [35]. The relevant decay modes
for the final gauge bosons are the leptonic ones since the momenta of these
leptons provide a measurement of the gauge boson polarizations, see for
instance [72] and also the references therein. In particular, the
corresponding angular distributions of the decays lead to a reconstruction of
the density matrix from the experimental data as it was developed in [16].

The aim of this work was to perform a theoretical determination of the density
matrix and the corresponding entanglement quantifiers for different VBS
processes, which are not explored in the literature as far as we know. Notice
that at the complete process level in Eq. (3.17), the resulting density matrix
represents a mixed-state at collider [25] via a convex combination of the
density matrices for each unpolarized gauge boson $V^{\prime}_{1,2}$ cases
considering them as ‘partons’ inside the initial fermions. In particular, this
mixed density matrix could be computed by means of the Weizsacker-Williams
Approximation [73, 74] for photons and the Effective $W$ Approximation [75,
76] for massive vector bosons. In that case just a lower bound of the
Concurrence can be determined [17] in order to decide if the system
$V_{1}\otimes V_{2}$ is entangled or separable. In addition, the experimental
challenges in each production mechanism and for each collider are different
then, dedicated studies are mandatory. This computation and the corresponding
reconstruction of the entanglement quantifiers from Monte-Carlo simulation is
beyond the scope of this work and it is postponed for a future study.

##  4 Summary and perspectives

In this work, the quantum properties of vector boson scattering were explored
by the computation of the quantifiers associated to entanglement/separability
and violation of Bell inequalities. In particular, the Negativity
$\mathcal{N}$, Entropy of Entanglement $\mathcal{S}_{\rm EE}$ and Concurrence
$\mathcal{C}$ can be treated as entanglement/separability parameters, in a
sense that they characterize a degree of entanglement, whereas the violation
of Bell inequalities is determined by the $\mathcal{I}$ parameter
corresponding to CHSH or CGLMP operators. This kind of processes allows to
examine pure bipartite systems, associated to the spin of the final vector
bosons, with density matrix of dimensions 2$\otimes$2 (for
$W^{+}W^{-}\to\gamma\gamma$), 3$\otimes$2 (for $W^{\pm}\gamma\to
W^{\pm}\gamma$, $W^{+}W^{-}\to Z\gamma$ and $W^{\pm}Z\to W^{\pm}\gamma$) and
3$\otimes$3 (for $\gamma\gamma\to W^{+}W^{-}$, $W^{\pm}\gamma\to W^{\pm}Z$,
$W^{+}W^{-}\to W^{+}W^{-}$, $W^{\pm}W^{\pm}\to W^{\pm}W^{\pm}$, $W^{+}W^{-}\to
ZZ$, $W^{\pm}Z\to W^{\pm}Z$, $Z\gamma\to W^{+}W^{-}$, $ZZ\to W^{+}W^{-}$ and
$ZZ\to ZZ$). The corresponding amplitudes were computed at tree level in the
context of the SM. Analytical expressions of density matrix and entanglement
quantifiers were presented for $W^{+}W^{-}\to\gamma\gamma$, $W^{\pm}\gamma\to
W^{\pm}\gamma$ and $\gamma\gamma\to W^{+}W^{-}$.

The goal of this paper was to determine the kinematical region in the plane
$[\cos(\theta),\sqrt{S}]$ where the resulting final vector bosons after the
scattering were entangled and then, if is it possible to test the Bell
inequality in that region. We found that all the final states are entangled
after the scattering, except for $W^{\pm}\gamma\to W^{\pm}\gamma$ in the
forward direction $\cos(\theta)=1$ or at energy equal to $\sqrt{S}=m_{\rm W}$.
Also, for $W^{+}W^{-}\to\gamma\gamma$, $W^{+}W^{-}\to Z\gamma$ and
$W^{\pm}Z\to W^{\pm}\gamma$ processes, the Maximally Entangled state is
reached in particular kinematical configuration since the entanglement
quantifiers achieve the maximum theoretical values there. On the other hand,
compared to the rest of the VBS processes, the $ZZ\to ZZ$ has the less
entangled final state in the whole kinematical plane.

Regarding the Bell inequality, we conclude that in the whole range of the
scattering angle and up to energy of 3 TeV, the CHSH $\mathcal{I}_{2}$
parameter for $W^{+}W^{-}\to\gamma\gamma$ is greater than 2 and the Maximally
Entangled state of Eq. (3.10) provides the maximal violation, corresponding to
the Cirelson bound 2$\sqrt{2}$. In addition, Bell inequality violation is
expected for the rest of the VBS processes but it occurs in small regions of
the phase space, except for the qutrit$\otimes$qubit processes nor $ZZ\to ZZ$
which the corresponding Bell parameters are always lower than 2. Therefore, it
is manifest that entanglement is not a sufficient condition for Bell
inequality violation. It is important to stress that the implemented
$\mathcal{I}_{3\otimes 2}$ quantifier was not optimized for these VBS
processes in the sense that the corresponding Bell operator is the
generalization of the CHSH to higher dimension. Analogously, $\mathcal{I}_{3}$
was optimized for Higgs boson decays into massive gauge bosons. Therefore, the
determination of more appropriate Bell operators for qutrit$\otimes$qubit and
two-qutrits VBS processes is an interesting improvement for future works.

The previous theoretical predictions intend to guide the experimental search
for the quantum properties of VBS in the kinematical plane since different
reconstruction techniques must be implemented depending on how boosted are the
final particles or if the available energy is near the threshold or much
higher. This study is a first step in that direction for this kind of
processes. At collider level, they correspond to sub-process in a more complex
scattering yielding to a mixed system for which the computation of the
entanglement quantifiers results in a non-trivial maximization problem over a
convex sum of the analyzed bipartite systems. The next step is to perform this
computation using the Effective W Approximation in order to get predictions
for the LHC and future lepton colliders. In addition, a quantum tomography
analysis with Monte-Carlo simulations will be developed to estimate the
significances to these observables, as it was done for Higgs boson decay and
diboson production in [15, 16, 17].

Other lines for further studies are related to BSM physics. Similar to EFT
analysis in [26, 30, 28, 77] or demanding maximal entanglement as fundamental
principle [78], quantum information measurements on VBS processes could
constrain new physics operators in the EWSB sector not explored yet.

## Acknowledgments

I am grateful to Bianca Polari who encouraged the abidance of this study. I
also thank to Alejandro Szynkman, Ernesto Arganda and Francisco Alonso for
their helpful discussions and suggestions. The present work has received
financial support from CONICET and ANPCyT under projects PICT 2017-2751, PICT
2018-03682 and PICT-2021-I-INVI-00374.

## Appendices

##  Appendix A VBS kinematics and amplitudes

This appendix is devoted to summarize the relevant details for the computation
of the amplitudes corresponding to VBS processes in Eq. (2.1). Without loss of
generality, the center-of-mass (CM) frame is chosen with the incoming
particles traveling along the $z$ axis and then are scattered into the $x-z$
plane with angle $\theta$. In that case, the 2$\to$2 amplitudes can be written
in terms of two kinematical variables: the CM energy $\sqrt{S}$ and the
scattering angle $\cos(\theta)$. Explicitly, the momenta $p^{\prime}_{i}$,
$p_{i}$ and the polarization vectors
$\varepsilon^{\prime}_{i}(s^{\prime}_{i})$, $\varepsilon_{i}(s_{i})$ with the
usual normalizations are

$\displaystyle p^{\prime}_{1}$ $\displaystyle=$
$\displaystyle\left(E^{\prime}_{1},0,0,-p_{in}\right)\hskip
56.9055ptp_{1}=\left(E_{1},-p_{out}\sin(\theta),0,-p_{out}\cos(\theta)\right)$
$\displaystyle\varepsilon^{\prime}_{1}(\pm)$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}\left(0,-i,\mp 1,0\right)\hskip
34.1433pt\varepsilon^{*}_{1}(\pm)=\frac{1}{\sqrt{2}}\left(0,i\cos(\theta),\mp
1,-i\sin(\theta)\right)$ $\displaystyle\varepsilon^{\prime}_{1}(0)$
$\displaystyle=$
$\displaystyle\frac{1}{m_{V^{\prime}_{1}}}\left(p_{in},0,0,-E^{\prime}_{1}\right)\hskip
22.76219pt\varepsilon^{*}_{1}(0)=\frac{1}{m_{V_{1}}}\left(p_{out},-E_{1}\sin(\theta),0,-E_{1}\cos(\theta)\right)$
$\displaystyle p^{\prime}_{2}$ $\displaystyle=$
$\displaystyle\left(E^{\prime}_{2},0,0,p_{in}\right)\hskip
65.44133ptp_{2}=\left(E_{2},p_{out}\sin(\theta),0,p_{out}\cos(\theta)\right)$
$\displaystyle\varepsilon^{\prime}_{2}(\pm)$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2}}\left(0,-i,\mp 1,0\right)\hskip
34.1433pt\varepsilon^{*}_{2}(\pm)=\frac{1}{\sqrt{2}}\left(0,i\cos(\theta),\mp
1,-i\sin(\theta)\right)$ $\displaystyle\varepsilon^{\prime}_{2}(0)$
$\displaystyle=$
$\displaystyle\frac{1}{m_{V^{\prime}_{2}}}\left(-p_{in},0,0,-E^{\prime}_{2}\right)\hskip
17.07164pt\varepsilon^{*}_{2}(0)=\frac{1}{m_{V_{2}}}\left(-p_{out},-E_{2}\sin(\theta),0,-E_{2}\cos(\theta)\right)$
(A.1)

with the trimomentum and energies given by

$\displaystyle p_{in}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\sqrt{\lambda(S,m_{V^{\prime}_{1}}^{2},m_{V^{\prime}_{2}}^{2})/S}\,,\hskip
22.76219ptE^{\prime}_{i}=\sqrt{m_{V^{\prime}_{i}}^{2}+p_{in}^{2}}\,,$
$\displaystyle p_{out}$ $\displaystyle=$
$\displaystyle\frac{1}{2}\sqrt{\lambda(S,m_{V_{1}}^{2},m_{V_{2}}^{2})/S}\,,\hskip
22.76219ptE_{i}=\sqrt{m_{V_{i}}^{2}+p_{out}^{2}}$ (A.2)

where the Kallen function is $\lambda(a,b,c)=a^{2}+b^{2}+c^{2}-2ab-2bc-2ca$.
The polarizations of the gauge bosons define the conventional basis,
corresponding to the third component of the spin, as
$\\{|+\rangle,|-\rangle\\}$ for photons and
$\\{|+\rangle,|0\rangle,|-\rangle\\}$ for $W^{\pm}$ and $Z$ bosons.

The SM amplitudes were computed in the Feynman-’t Hooft gauge at tree level.
The $W^{\pm}\gamma\to W^{\pm}\gamma$ and $W^{\pm}\gamma\to W^{\pm}Z$ have the
same Feynman diagrams, which are presented in the first row of Fig. 8 as
$W^{\pm}\gamma\to W^{\pm}V$. Diagrams with $W^{\pm}$ and charged Goldstone
$\varphi^{\pm}$ bosons in the $U$-channel are omitted for simplicity. Notice
that the Higgs boson is not present in these processes. Of course,
$W^{+}W^{-}\to\gamma V$ and $W^{\pm}Z\to W^{\pm}\gamma$ are related by
crossing symmetry and have similar diagrams. For $W^{\pm}Z\to W^{\pm}Z$, and
the related processes by crossing symmetry $W^{+}W^{-}\to ZZ$ and $ZZ\to
W^{+}W^{-}$, there is an additional diagram with the Higgs boson as mediator
respect to the previous ones as can be seen in the second row. On the other
hand, the $W^{+}W^{-}\to W^{+}W^{-}$ and $W^{\pm}W^{\pm}\to W^{\pm}W^{\pm}$
processes in the third row of Fig. 8 have two diagrams with the Higgs boson as
mediator. Finally, the $ZZ\to ZZ$ process only has 3 diagrams with the Higgs
boson as mediator in $S$-, $T$\- and $U$-channels, which are shown in the last
row.

Figure 8: Representative Feynman diagrams for some of the considered VBS
channels ($V$ stands for photon and $Z$ boson). Processes related by crossing-
symmetry are described in the text.

##  Appendix B Generalized Gell-Mann and spin matrices for dimension 2 and 3

For completeness, this appendix gathers the explicit form of the generalized
Gell-Mann matrices that enter in the decomposition of the density matrix in
Eq. (2.5). An important relation is the trace orthogonality of these matrices:

${\rm Tr}\big{[}\lambda^{(d)}_{i}\big{]}=0\quad\text{and}\quad{\rm
Tr}\big{[}\lambda^{(d)}_{i}\lambda^{(d)}_{j}\big{]}=2\delta_{ij}$ (B.1)

For qubits, the dimension-2 Pauli matrices are

$\sigma_{1}=\lambda_{1}^{(2)}=\begin{pmatrix}0&1\\\ 1&0\end{pmatrix}\,,\hskip
22.76219pt\sigma_{2}=\lambda_{2}^{(2)}=\begin{pmatrix}0&-i\\\
i&0\end{pmatrix}\,,\hskip
22.76219pt\sigma_{3}=\lambda_{3}^{(2)}=\begin{pmatrix}1&0\\\
0&-1\end{pmatrix}$ (B.2)

and also the dimension-3 representation of the eight Gell-Mann matrices (for
simplicity, the superscript ‘(3)’ is omitted) are

$\displaystyle\lambda_{1}$ $\displaystyle=$
$\displaystyle\begin{pmatrix}0&1&0\\\ 1&0&0\\\ 0&0&0\end{pmatrix}\,,\hskip
14.22636pt\lambda_{2}=\begin{pmatrix}0&-i&0\\\ i&0&0\\\
0&0&0\end{pmatrix}\,,\hskip 14.22636pt\lambda_{3}=\begin{pmatrix}1&0&0\\\
0&-1&0\\\ 0&0&0\end{pmatrix}\,,\hskip
14.22636pt\lambda_{4}=\begin{pmatrix}0&0&1\\\ 0&0&0\\\ 1&0&0\end{pmatrix}\,,$
$\displaystyle\lambda_{5}$ $\displaystyle=$
$\displaystyle\begin{pmatrix}0&0&-i\\\ 0&0&0\\\ i&0&0\end{pmatrix}\,,\hskip
14.22636pt\lambda_{6}=\begin{pmatrix}0&0&0\\\ 0&0&1\\\
0&1&0\end{pmatrix}\,,\hskip 14.22636pt\lambda_{7}=\begin{pmatrix}0&0&0\\\
0&0&-i\\\ 0&i&0\end{pmatrix}\,,\hskip
14.22636pt\lambda_{8}=\frac{1}{\sqrt{3}}\begin{pmatrix}1&0&0\\\ 0&1&0\\\
0&0&-2\end{pmatrix}$ (B.3)

In addition, the spin-1 matrices in Eq. (2.18) are

$\displaystyle
S_{x}=\frac{1}{\sqrt{2}}(\lambda_{1}+\lambda_{6})=\frac{1}{\sqrt{2}}\begin{pmatrix}0&1&0\\\
1&0&1\\\ 0&1&0\end{pmatrix}\,,\hskip
19.91692ptS_{y}=\frac{1}{\sqrt{2}}(\lambda_{2}+\lambda_{7})=\frac{1}{\sqrt{2}}\begin{pmatrix}0&-i&0\\\
i&0&-i\\\ 0&i&0\end{pmatrix}\,,$ $\displaystyle\hskip
85.35826ptS_{z}=\frac{1}{2}(\lambda_{3}+\sqrt{3}\lambda_{8})=\begin{pmatrix}1&0&0\\\
0&0&0\\\ 0&0&-1\end{pmatrix}$ (B.4)

##  Appendix C Additional plots of entanglement quantifiers

This appendix collects the plots of the quantifiers that were omitted in the
main text for saving space. Fig. 9 corresponds to the Entropy of Entanglement
for the two-qubits and qutrit$\otimes$qubit processes. As for the Negativity
of these VBS, $\mathcal{S}_{\rm EE}$ never vanishes in the kinematical plane
then the final states are entangled333Except for $W^{\pm}\gamma\to
W^{\pm}\gamma$ at the threshold and in the forward direction as can be seen
from Eq. (D.5).. The maximum theoretical values for this quantifier,
corresponding to the Maximally Entangled states, are achieved in
$W^{+}W^{-}\to\gamma\gamma$ (showed in solid black line in the upper-left
plot), and in both $W^{+}W^{-}\to Z\gamma$ and $W^{\pm}Z\to W^{\pm}\gamma$
(second row of this figure).

|
---|---
|
Figure 9: Entropy of Entanglement for $W^{+}W^{-}\to\gamma\gamma$ (upper-
left), $W^{\pm}\gamma\to W^{\pm}\gamma$ (upper-right), $W^{+}W^{-}\to Z\gamma$
(lower-left) and $W^{\pm}Z\to W^{\pm}\gamma$ (lower-right) in the plane
$[\cos(\theta),\sqrt{S}]$. Contour lines are shown for an easy comparison of
the numerical values.

Figs. 10-11 show the Concurrence for all the processes with massive gauge
bosons in final state, i.e. bipartite 3$\otimes$3 system. The same scale for
all the plots is used, making more transparent the comparison among them. Now,
the lower values are denoted in violet colour and never reach the zero value,
i.e. the final states are entangled in the whole kinematical plane. The
maximal value of the Concurrence is 1.12 for $ZZ\to W^{+}W^{-}$ which is close
to the theoretical maximum is $2/\sqrt{3}$.

|
---|---
|

Figure 10: Concurrence for different 3$\otimes$3 VBS processes in the plane
$[\cos(\theta),\sqrt{S}]$. Contour lines are shown for an easy comparison of
the numerical values.

need space

|
---|---
|
Figure 11: Concurrence for the rest of 3$\otimes$3 VBS processes in the plane
$[\cos(\theta),\sqrt{S}]$. Contour lines are shown for an easy comparison of
the numerical values.

need space

##  Appendix D Analytical expressions for $W^{\pm}\gamma\to W^{\pm}\gamma$ and
$\gamma\gamma\to W^{+}W^{-}$

In this appendix, the analytical expressions of the non-vanishing coefficients
in the decomposition of Eq. (2.5) for the density matrix corresponding to
$W^{\pm}\gamma\to W^{\pm}\gamma$ and $\gamma\gamma\to W^{+}W^{-}$ are
presented. In addition, the eigenvalues of the reduced matrices that enter in
the computation of the Entanglement Entropy in Eq. (2.9) are also shown. These
observables are written in terms of the kinematical variables $\sqrt{S}$ and
$c=\cos(\theta)$ and allow us to understand the numerical results.

###  D.1 $W^{\pm}\gamma\to W^{\pm}\gamma$ process

For this qutrit$\otimes$qubit case, the decomposition of the 6$\times$6
density matrix is

$\rho=\frac{1}{6}I_{6}+\frac{1}{4}\sum_{i=1}^{8}A_{i}\lambda_{i}\otimes
I_{2}+\frac{1}{6}\sum_{j=1}^{3}B_{j}I_{3}\otimes\sigma_{j}+\frac{1}{4}\sum_{i=1}^{8}\sum_{j=1}^{3}C_{ij}\lambda_{i}\otimes\sigma_{j}$
(D.1)

For a compact notation, we define the quantity

$\displaystyle D_{W\gamma}$ $\displaystyle=$ $\displaystyle(1-c)^{4}m_{\rm
W}^{8}-2(1-c)^{4}m_{\rm W}^{6}S+(1-c)^{2}\left(2c^{2}+7\right)m_{\rm
W}^{4}S^{2}$ (D.2) $\displaystyle-2\left(c^{4}-3c^{2}-4c+6\right)m_{\rm
W}^{2}S^{3}+\left(c^{4}+3c^{2}+2c+6\right)S^{4}$

The resulting non-vanishing coefficients $A_{i}$ and $B_{j}$ are

$\displaystyle A_{1}=A_{6}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{W\gamma}}S\left((1-c)^{3}m_{\rm
W}^{6}+(1-c)^{2}(3c+1)m_{\rm W}^{4}S\right.$ $\displaystyle\left.\hskip
22.76219pt-\left(3c^{3}-9c^{2}+c+5\right)m_{\rm
W}^{2}S^{2}+\left(c^{3}+c^{2}+3c+3\right)S^{3}\right)$ $\displaystyle
A_{2}=-A_{7}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{W\gamma}}\sqrt{2(1-c^{2})S}m_{\rm
W}\left(-(1-c)^{3}m_{\rm W}^{6}+(1-c)^{2}(c+4)m_{\rm W}^{4}S\right.$
$\displaystyle\left.\hskip 22.76219pt+\left(c^{3}+c^{2}+3c-5\right)m_{\rm
W}^{2}S^{2}+\left(c^{3}+c+2\right)S^{3}\right)$ $\displaystyle
A_{3}=-\sqrt{3}A_{8}$ $\displaystyle=$
$\displaystyle\frac{1}{2D_{W\gamma}}\left((1-c)^{4}m_{\rm
W}^{8}-4(1-c)^{3}(c+2)m_{\rm W}^{6}S\right.$ $\displaystyle\left.\hskip
28.45274pt+2(1-c)^{2}\left(c^{2}+6c+8\right)m_{\rm
W}^{4}S^{2}+4\left(c^{4}+2c-3\right)m_{\rm
W}^{2}S^{3}+\left(c^{4}-4c+3\right)S^{4}\right)$ $\displaystyle A_{4}$
$\displaystyle=$
$\displaystyle\frac{2}{D_{W\gamma}}S\left(-(1-c)^{3}(c+1)m_{\rm
W}^{6}+2\left(c^{3}-3c+2\right)m_{\rm W}^{4}S\right.$
$\displaystyle\left.\hskip 22.76219pt+\left(c^{4}+4c-5\right)m_{\rm
W}^{2}S^{2}+2\left(c^{2}+1\right)S^{3}\right)$ $\displaystyle B_{1}$
$\displaystyle=$
$\displaystyle\frac{1}{D_{W\gamma}}S^{2}\left(5(1-c)^{2}m_{\rm
W}^{4}+2\left(c^{2}+4c-5\right)m_{\rm
W}^{2}S+\left(5c^{2}+2c+5\right)S^{2}\right)$ (D.3)

and the non-vanishing coefficients $C_{ij}$ of the correlation matrix are

$\displaystyle C_{11}=C_{61}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{W\gamma}}S\left((1-c)^{3}m_{\rm
W}^{6}+(1-c)^{2}(3c+1)m_{\rm W}^{4}S\right.$ $\displaystyle\left.\hskip
22.76219pt-\left(3c^{3}-9c^{2}+c+5\right)m_{\rm
W}^{2}S^{2}+\left(c^{3}+c^{2}+3c+3\right)S^{3}\right)$ $\displaystyle
C_{12}=C_{23}=-C_{62}=C_{73}$ $\displaystyle=$
$\displaystyle-\frac{1}{D_{W\gamma}}(1-c)\sqrt{2(1-c^{2})S}m_{\rm W}(S-m_{\rm
W}^{2})\left((1-c)^{2}m_{\rm W}^{4}\right.$ $\displaystyle\left.\hskip
34.1433pt-(1-c)m_{\rm W}^{2}S+c(c+1)S^{2}\right)$ $\displaystyle
C_{13}=C_{22}=-C_{63}=-C_{72}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{W\gamma}}(1-c)(S-m_{\rm
W}^{2})S\left((1-c)^{2}m_{\rm W}^{4}-(1+c)^{2}S^{2}\right)$ $\displaystyle
C_{21}=-C_{71}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{W\gamma}}\sqrt{2(1-c^{2})S}m_{\rm
W}\left(-(1-c)^{3}m_{\rm W}^{6}+(1-c)^{2}(c+4)m_{\rm W}^{4}S\right.$
$\displaystyle\left.\hskip 22.76219pt+\left(c^{3}+c^{2}+3c-5\right)m_{\rm
W}^{2}S^{2}+\left(c^{3}+c+2\right)S^{3}\right)$ $\displaystyle
C_{31}=-\sqrt{3}C_{81}$ $\displaystyle=$
$\displaystyle-\frac{1}{D_{W\gamma}}(1-c)S\left(3(1-c)^{2}(c+1)m_{\rm
W}^{6}+\left(6c^{2}+c-7\right)m_{\rm W}^{4}S\right.$
$\displaystyle\left.\hskip 34.1433pt+\left(3c^{3}+3c^{2}+c+5\right)m_{\rm
W}^{2}S^{2}-(1-c)S^{3}\right)$ $\displaystyle
C_{32}=\frac{1}{2}C_{53}=\frac{1}{\sqrt{6}}C_{82}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{W\gamma}}(1-c)^{2}\sqrt{2(1-c^{2})S}m_{\rm
W}\left(m_{\rm W}^{4}S-2m_{\rm W}^{2}S^{2}+S^{3}\right)$ $\displaystyle
C_{33}=-\frac{1}{2}C_{52}=\frac{1}{\sqrt{3}}C_{83}$ $\displaystyle=$
$\displaystyle-\frac{1}{2D_{W\gamma}}(1-c)(S-m_{\rm
W}^{2})\left(-(1-c)^{3}m_{\rm W}^{6}+(1-c)^{2}(c+3)m_{\rm W}^{4}S\right.$
$\displaystyle\left.\hskip 34.1433pt+\left(c^{3}+c^{2}+3c-5\right)m_{\rm
W}^{2}S^{2}+\left(c^{3}+c^{2}+3c+3\right)S^{3}\right)$ $\displaystyle C_{41}$
$\displaystyle=$ $\displaystyle\frac{1}{D_{W\gamma}}\left((1-c)^{4}m_{\rm
W}^{8}-4(1-c)^{3}m_{\rm W}^{6}S+2\left(c^{4}+2c^{2}-8c+5\right)m_{\rm
W}^{4}S^{2}\right.$ (D.4) $\displaystyle\left.\hskip
22.76219pt+4\left(c^{2}+2c-3\right)m_{\rm
W}^{2}S^{3}+\left(c^{4}+2c^{2}+5\right)S^{4}\right)$

Finally, the two eigenvalues of the reduced density matrix respect to the
$W^{\pm}$ boson are

$\displaystyle\lambda_{1}^{red}$ $\displaystyle=$
$\displaystyle\frac{1}{2D_{W\gamma}}(1-c)^{2}\left(S-m_{\rm
W}^{2}\right)^{2}\left((1-c)^{2}m_{\rm W}^{4}+(c+1)^{2}S^{2}\right)$
$\displaystyle\lambda_{2}^{red}$ $\displaystyle=$
$\displaystyle\frac{1}{2D_{W\gamma}}\left((1-c)^{4}m_{\rm
W}^{8}-2(1-c)^{4}m_{\rm W}^{6}S+2(1-c)^{2}\left(c^{2}+6\right)m_{\rm
W}^{4}S^{2}\right.$ (D.5) $\displaystyle\left.\hskip
28.45274pt-2\left(c^{4}-4c^{2}-8c+11\right)m_{\rm
W}^{2}S^{3}+\left(c^{4}+8c^{2}+4c+11\right)S^{4}\right)$

The normalization of the reduced density matrix yields to
$\lambda_{1}^{red}+\lambda_{2}^{red}=1$. Also, in the analyzed phase space,
$D_{W\gamma}>0$ and it is easy to see that $\lambda_{1}^{red}\geq 0$ with
equality if and only if $\cos(\theta)=1$ or $S=m_{\rm W}^{2}$. In that
kinematical configurations, $\lambda_{2}^{red}=1$ and the resulting Entropy of
Entanglement of Eq. (2.9) vanishes, then we conclude that the $W^{\pm}\gamma$
final state is separable. This result was discussed in the text from the plot
of the Negativity in the upper-left panel of Fig. 2. The behaviour of the
resulting $\mathcal{S}_{\rm EE}$ in the whole kinematical plane is shown in
the upper-right panel of Fig. 9

###  D.2 $\gamma\gamma\to W^{+}W^{-}$ process

For this two-qutrits case, the decomposition of the 9$\times$9 density matrix
is

$\rho=\frac{1}{9}I_{9}+\frac{1}{6}\sum_{i=1}^{8}A_{i}\lambda_{i}\otimes
I_{3}+\frac{1}{6}\sum_{j=1}^{8}B_{j}I_{3}\otimes\lambda_{j}+\frac{1}{4}\sum_{i=1}^{8}\sum_{j=1}^{8}C_{ij}\lambda_{i}\otimes\lambda_{j}$
(D.6)

For a compact notation, we define the quantity

$D_{WW}=48\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}-8\left(3c^{4}-5c^{2}+6\right)m_{\rm
W}^{2}S+\left(3c^{4}+2c^{2}+11\right)S^{2}$ (D.7)

The resulting non-vanishing coefficients $A_{i}$ and $B_{j}$ are

$\displaystyle A_{2}=-A_{7}=B_{2}=-B_{7}$ $\displaystyle=$
$\displaystyle\frac{16}{D_{WW}}c\sqrt{2(1-c^{2})}m_{\rm W}S^{3/2}$
$\displaystyle A_{3}=-\sqrt{3}A_{8}=B_{3}=-\sqrt{3}B_{8}$ $\displaystyle=$
$\displaystyle\frac{4}{D_{WW}}S\left(4\left(2c^{2}-3\right)m_{\rm
W}^{2}+\left(c^{2}+1\right)S\right)$ $\displaystyle A_{4}=B_{4}$
$\displaystyle=$ $\displaystyle\frac{8}{D_{WW}}S\left(-4m_{\rm
W}^{2}+(1+c^{2})S\right)$ (D.8)

and the non-vanishing coefficients $C_{ij}$ of the correlation matrix are

$\displaystyle C_{11}=C_{66}$ $\displaystyle=$
$\displaystyle-\frac{2}{D_{WW}}S\left(4\left(c^{4}-c^{2}+2\right)m_{\rm
W}^{2}+(1-c^{4})S\right)$ $\displaystyle C_{15}=C_{51}$ $\displaystyle=$
$\displaystyle=-C_{56}=-C_{65}$ $\displaystyle=$
$\displaystyle-\frac{4}{D_{WW}}c\sqrt{2(1-c^{2})S}m_{\rm
W}\left(4\left(c^{2}-2\right)m_{\rm W}^{2}+(1-c^{2})S\right)$ $\displaystyle
C_{16}=C_{61}$ $\displaystyle=$
$\displaystyle\frac{4}{D_{WW}}\left(8\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}-2\left(c^{4}-c^{2}+2\right)m_{\rm
W}^{2}S-\left(1-c^{2}\right)S^{2}\right)$ $\displaystyle C_{22}=C_{77}$
$\displaystyle=$
$\displaystyle\frac{2}{D_{WW}}S\left(\left(-12c^{4}+12c^{2}+8\right)m_{\rm
W}^{2}+(1-c^{4})S\right)$ $\displaystyle C_{23}=C_{32}$ $\displaystyle=$
$\displaystyle\frac{8}{D_{WW}}c\sqrt{2(1-c^{2})S}m_{\rm
W}\left(2\left(c^{2}-2\right)m_{\rm W}^{2}+c^{2}S\right)$ $\displaystyle
C_{24}=C_{42}$ $\displaystyle=$ $\displaystyle=-C_{47}=-C_{74}$
$\displaystyle=$ $\displaystyle\frac{4}{D_{WW}}c\sqrt{2(1-c^{2})S}m_{\rm
W}\left(4\left(c^{2}-2\right)m_{\rm W}^{2}+\left(c^{2}+3\right)S\right)$
$\displaystyle C_{27}=C_{72}$ $\displaystyle=$
$\displaystyle\frac{4}{D_{WW}}\left(8\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}+2\left(3c^{4}-3c^{2}-2\right)m_{\rm
W}^{2}S-\left(1-c^{2}\right)S^{2}\right)$ $\displaystyle C_{28}=C_{82}$
$\displaystyle=$
$\displaystyle-\frac{8}{D_{WW}}c\sqrt{\frac{2}{3}(1-c^{2})S}m_{\rm
W}\left(6\left(c^{2}-2\right)m_{\rm W}^{2}+S\right)$ $\displaystyle C_{33}$
$\displaystyle=$
$\displaystyle\frac{2}{D_{WW}}\left(8\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}+4\left(3c^{4}-5c^{2}+2\right)m_{\rm
W}^{2}S+\left(c^{4}+1\right)S^{2}\right)$ $\displaystyle C_{34}=C_{43}$
$\displaystyle=$ $\displaystyle=-\sqrt{3}C_{48}=-\sqrt{3}C_{84}$
$\displaystyle=$
$\displaystyle\frac{4}{D_{WW}}S\left(\left(6c^{4}-6c^{2}-4\right)m_{\rm
W}^{2}+\left(c^{2}+1\right)S\right)$ $\displaystyle C_{37}=C_{73}$
$\displaystyle=$ $\displaystyle-\frac{4}{D_{WW}}c\sqrt{2(1-c^{2})S}m_{\rm
W}\left(8\left(c^{2}-2\right)m_{\rm W}^{2}+\left(c^{2}+1\right)S\right)$
$\displaystyle C_{38}=C_{83}$ $\displaystyle=$
$\displaystyle\frac{4}{\sqrt{3}D_{WW}}\left(-12\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}+2\left(-3c^{4}+c^{2}+6\right)m_{\rm
W}^{2}S+\left(c^{2}-2\right)S^{2}\right)$ $\displaystyle C_{44}$
$\displaystyle=$
$\displaystyle\frac{1}{D_{WW}}\left(8\left(2\left(c^{2}-2\right)m_{\rm
W}^{2}+S\right)^{2}+2\left(c^{2}+1\right)^{2}S^{2}\right)$ $\displaystyle
C_{55}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{WW}}\left(8\left(2\left(c^{2}-2\right)m_{\rm
W}^{2}+S\right)^{2}-2\left(c^{2}+1\right)^{2}S^{2}\right)$ $\displaystyle
C_{78}=C_{87}$ $\displaystyle=$
$\displaystyle\frac{4}{D_{WW}}c(-1+3c^{2})\sqrt{\frac{2}{3}(1-c^{2})}m_{\rm
W}S^{3/2}$ $\displaystyle C_{88}$ $\displaystyle=$
$\displaystyle\frac{1}{3D_{WW}}\left(-48\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}+8\left(3c^{4}-13c^{2}+18\right)m_{\rm W}^{2}S\right.$ (D.9)
$\displaystyle\left.\hskip
34.1433pt+2\left(3c^{4}+4c^{2}-5\right)S^{2}\right)$

Finally, the three eigenvalues of the reduced density matrix that enter in Eq.
(2.9) are

$\displaystyle\lambda_{1}^{red}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{WW}}\left(4\left(c^{2}-2\right)m_{\rm
W}^{2}+(1-c^{2})S\right)^{2}$ $\displaystyle\lambda_{2}^{red}$
$\displaystyle=$
$\displaystyle\frac{1}{D_{WW}}\left(16\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}-8\left(c^{4}-c^{2}+2\right)m_{\rm W}^{2}S+(5+2c^{2}+c^{4})S^{2}\right.$
$\displaystyle\left.\hskip 22.76219pt-4S\sqrt{16\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}-8\left(c^{4}-c^{2}+2\right)m_{\rm
W}^{2}S+\left(c^{2}+1\right)^{2}S^{2}}\right)$
$\displaystyle\lambda_{3}^{red}$ $\displaystyle=$
$\displaystyle\frac{1}{D_{WW}}\left(16\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}-8\left(c^{4}-c^{2}+2\right)m_{\rm W}^{2}S+(5+2c^{2}+c^{4})S^{2}\right.$
(D.10) $\displaystyle\left.\hskip
22.76219pt+4S\sqrt{16\left(c^{2}-2\right)^{2}m_{\rm
W}^{4}-8\left(c^{4}-c^{2}+2\right)m_{\rm
W}^{2}S+\left(c^{2}+1\right)^{2}S^{2}}\right)$

The behaviour of the resulting $\mathcal{S}_{\rm EE}$ in the whole kinematical
plane is shown in the upper-left panel of Fig. 3

## References

  * Horodecki  _et al._ [2009] R. Horodecki, P. Horodecki, M. Horodecki, and K. Horodecki, Rev. Mod. Phys. 81, 865 (2009), arXiv:quant-ph/0702225
  * Bell [1964] J. S. Bell, Physics Physique Fizika 1, 195 (1964)
  * Ekert [1991] A. K. Ekert, Phys. Rev. Lett. 67, 661 (1991)
  * Bennett  _et al._ [1993] C. H. Bennett, G. Brassard, C. Crépeau, R. Jozsa, A. Peres, and W. K. Wootters, Phys. Rev. Lett. 70, 1895 (1993)
  * Raussendorf and Briegel [2001] R. Raussendorf and H. J. Briegel, Phys. Rev. Lett. 86, 5188 (2001)
  * Casini and Huerta [2023] H. Casini and M. Huerta, PoS TASI2021, 002 (2023), arXiv:2201.13310 [hep-th]
  * Tornqvist [1981] N. A. Tornqvist, Found. Phys. 11, 171 (1981)
  * Baranov [2008] S. P. Baranov, Journal of Physics G: Nuclear and Particle Physics 35, 075002 (2008)
  * Chen  _et al._ [2013] S. Chen, Y. Nakaguchi, and S. Komamiya, Progress of Theoretical and Experimental Physics 2013 (2013), 10.1093/ptep/ptt032, 063A01, https://academic.oup.com/ptep/article-pdf/2013/6/063A01/19300324/ptt032.pdf
  * Acín  _et al._ [2001] A. Acín, J. I. Latorre, and P. Pascual, Phys. Rev. A 63, 042107 (2001)
  * Blasone  _et al._ [2009] M. Blasone, F. Dell’Anno, S. D. Siena, and F. Illuminati, Europhysics Letters 85, 50002 (2009)
  * Banerjee  _et al._ [2015] S. Banerjee, A. K. Alok, R. Srikanth, and B. C. Hiesmayr, The European Physical Journal C 75 (2015), 10.1140/epjc/s10052-015-3717-x
  * Barr [2022] A. J. Barr, Phys. Lett. B 825, 136866 (2022), arXiv:2106.01377 [hep-ph]
  * Aguilar-Saavedra [2023] J. A. Aguilar-Saavedra, Phys. Rev. D 107, 076016 (2023), arXiv:2209.14033 [hep-ph]
  * Aguilar-Saavedra  _et al._ [2023] J. A. Aguilar-Saavedra, A. Bernal, J. A. Casas, and J. M. Moreno, Phys. Rev. D 107, 016012 (2023), arXiv:2209.13441 [hep-ph]
  * Ashby-Pickering  _et al._ [2023] R. Ashby-Pickering, A. J. Barr, and A. Wierzchucka, JHEP 05, 020 (2023), arXiv:2209.13990 [quant-ph]
  * Fabbrichesi  _et al._ [2023a] M. Fabbrichesi, R. Floreanini, E. Gabrielli, and L. Marzola, (2023a), arXiv:2302.00683 [hep-ph]
  * Bertlmann and Hiesmayr [2001] R. A. Bertlmann and B. C. Hiesmayr, Physical Review A 63 (2001), 10.1103/physreva.63.062112
  * Banerjee  _et al._ [2016] S. Banerjee, A. K. Alok, and R. MacKenzie, The European Physical Journal Plus 131 (2016), 10.1140/epjp/i2016-16129-0
  * Takubo  _et al._ [2021] Y. Takubo, T. Ichikawa, S. Higashino, Y. Mori, K. Nagano, and I. Tsutsui, Physical Review D 104 (2021), 10.1103/physrevd.104.056004
  * Fabbrichesi  _et al._ [2023b] M. Fabbrichesi, R. Floreanini, E. Gabrielli, and L. Marzola, “Bell inequality is violated in $b^{0}\to j/\psi\,k^{\star}(892)^{0}$ decays,” (2023b), arXiv:2305.04982 [hep-ph]
  * Afik and de Nova [2021] Y. Afik and J. R. M. n. de Nova, Eur. Phys. J. Plus 136, 907 (2021), arXiv:2003.02280 [quant-ph]
  * Fabbrichesi  _et al._ [2021] M. Fabbrichesi, R. Floreanini, and G. Panizzo, Phys. Rev. Lett. 127, 161801 (2021), arXiv:2102.11883 [hep-ph]
  * Severi  _et al._ [2022] C. Severi, C. D. E. Boschi, F. Maltoni, and M. Sioli, Eur. Phys. J. C 82, 285 (2022), arXiv:2110.10112 [hep-ph]
  * Afik and de Nova [2022] Y. Afik and J. R. M. n. de Nova, Quantum 6, 820 (2022), arXiv:2203.05582 [quant-ph]
  * Aoude  _et al._ [2022] R. Aoude, E. Madge, F. Maltoni, and L. Mantani, Phys. Rev. D 106, 055007 (2022), arXiv:2203.05619 [hep-ph]
  * Aguilar-Saavedra and Casas [2022] J. A. Aguilar-Saavedra and J. A. Casas, Eur. Phys. J. C 82, 666 (2022), arXiv:2205.00542 [hep-ph]
  * Severi and Vryonidou [2023] C. Severi and E. Vryonidou, JHEP 01, 148 (2023), arXiv:2210.09330 [hep-ph]
  * Dong  _et al._ [2023] Z. Dong, D. Gonçalves, K. Kong, and A. Navarro, (2023), arXiv:2305.07075 [hep-ph]
  * Fabbrichesi  _et al._ [2023c] M. Fabbrichesi, R. Floreanini, and E. Gabrielli, Eur. Phys. J. C 83, 162 (2023c), arXiv:2208.11723 [hep-ph]
  * Altakach  _et al._ [2023] M. M. Altakach, P. Lamba, F. Maltoni, K. Mawatari, and K. Sakurai, Phys. Rev. D 107, 093002 (2023), arXiv:2211.10513 [hep-ph]
  * Barr  _et al._ [2023] A. J. Barr, P. Caban, and J. Rembieliński, Quantum 7, 1070 (2023), arXiv:2204.11063 [quant-ph]
  * Fedida and Serafini [2023] S. Fedida and A. Serafini, Phys. Rev. D 107, 116007 (2023), arXiv:2209.01405 [quant-ph]
  * Covarelli  _et al._ [2021] R. Covarelli, M. Pellen, and M. Zaro, Int. J. Mod. Phys. A 36, 2130009 (2021), arXiv:2102.10991 [hep-ph]
  * Buarque Franzosi  _et al._ [2022] D. Buarque Franzosi _et al._ , Rev. Phys. 8, 100071 (2022), arXiv:2106.01393 [hep-ph]
  * Fano [1957] U. Fano, Rev. Mod. Phys. 29, 74 (1957)
  * Stanton [1971] L. Stanton, Molecular Physics 20, 655 (1971)
  * Hagston and Roberts [1980] W. E. Hagston and M. Roberts, Journal of Physics A: Mathematical and General 13, 2695 (1980)
  * Peres [1996] A. Peres, Phys. Rev. Lett. 77, 1413 (1996), arXiv:quant-ph/9604005
  * Horodecki  _et al._ [1996] M. Horodecki, P. Horodecki, and R. Horodecki, Phys. Lett. A 223, 1 (1996), arXiv:quant-ph/9605038
  * Vidal and Werner [2002] G. Vidal and R. F. Werner, Phys. Rev. A 65, 032314 (2002)
  * Bennett  _et al._ [1996] C. H. Bennett, H. J. Bernstein, S. Popescu, and B. Schumacher, Phys. Rev. A 53, 2046 (1996)
  * Rungta  _et al._ [2001] P. Rungta, V. Buž ek, C. M. Caves, M. Hillery, and G. J. Milburn, Physical Review A 64 (2001), 10.1103/physreva.64.042315
  * Hill and Wootters [1997] S. Hill and W. K. Wootters, Phys. Rev. Lett. 78, 5022 (1997), arXiv:quant-ph/9703041
  * Eltschka  _et al._ [2015] C. Eltschka, G. Tóth, and J. Siewert, Physical Review A 91 (2015), 10.1103/physreva.91.032327
  * Clauser  _et al._ [1969] J. F. Clauser, M. A. Horne, A. Shimony, and R. A. Holt, Phys. Rev. Lett. 23, 880 (1969)
  * Horodecki  _et al._ [1995] R. Horodecki, P. Horodecki, and M. Horodecki, Physics Letters A 200, 340 (1995)
  * Cirelson [1980] B. S. Cirelson, Lett. Math. Phys. 4, 93 (1980)
  * Caban  _et al._ [2008] P. Caban, J. Rembielinski, and M. Wlodarczyk, Phys. Rev. A 77, 012103 (2008), arXiv:0801.3200 [quant-ph]
  * Collins  _et al._ [2002] D. Collins, N. Gisin, N. Linden, S. Massar, and S. Popescu, Phys. Rev. Lett. 88, 040404 (2002)
  * Acín  _et al._ [2002] A. Acín, T. Durt, N. Gisin, and J. I. Latorre, Phys. Rev. A 65, 052325 (2002)
  * Hahn [2001] T. Hahn, Comput. Phys. Commun. 140, 418 (2001), arXiv:hep-ph/0012260
  * Hahn and Perez-Victoria [1999] T. Hahn and M. Perez-Victoria, Comput. Phys. Commun. 118, 153 (1999), arXiv:hep-ph/9807565
  * Aaij  _et al._ [2022] R. Aaij _et al._ (LHCb), Phys. Rev. D 105, L051104 (2022), arXiv:2111.10194 [hep-ex]
  * Bishara  _et al._ [2014] F. Bishara, Y. Grossman, R. Harnik, D. J. Robinson, J. Shu, and J. Zupan, JHEP 04, 084 (2014), arXiv:1312.2955 [hep-ph]
  * Gritsan  _et al._ [2022] A. V. Gritsan _et al._ , (2022), arXiv:2205.07715 [hep-ex]
  * Covarelli [2021] R. Covarelli (ATLAS, CMS), PoS LHCP2021, 126 (2021)
  * Forster [2022] F. Z. Forster, “Sensitivity studies for the observation of electroweak diphoton+jets production at the LHC,”  (2022), presented 2022
  * Tumasyan  _et al._ [2023a] A. Tumasyan _et al._ (CMS), Phys. Rev. D 108, 032017 (2023a), arXiv:2212.12592 [hep-ex]
  * Tumasyan  _et al._ [2021] A. Tumasyan _et al._ (CMS), Phys. Rev. D 104, 072001 (2021), arXiv:2106.11082 [hep-ex]
  * Aad  _et al._ [2023a] G. Aad _et al._ (ATLAS), (2023a), arXiv:2305.19142 [hep-ex]
  * Aaboud  _et al._ [2019a] M. Aaboud _et al._ (ATLAS), Phys. Rev. Lett. 123, 161801 (2019a), arXiv:1906.03203 [hep-ex]
  * Sirunyan  _et al._ [2021a] A. M. Sirunyan _et al._ (CMS), Phys. Lett. B 812, 136018 (2021a), arXiv:2009.09429 [hep-ex]
  * Aaboud  _et al._ [2019b] M. Aaboud _et al._ (ATLAS), Phys. Lett. B 793, 469 (2019b), arXiv:1812.09740 [hep-ex]
  * Aad  _et al._ [2019] G. Aad _et al._ (ATLAS), Phys. Rev. D 100, 032007 (2019), arXiv:1905.07714 [hep-ex]
  * Sirunyan  _et al._ [2020] A. M. Sirunyan _et al._ (CMS), Phys. Lett. B 809, 135710 (2020), arXiv:2005.01173 [hep-ex]
  * Tumasyan  _et al._ [2022] A. Tumasyan _et al._ (CMS), Phys. Lett. B 834, 137438 (2022), arXiv:2112.05259 [hep-ex]
  * Aad  _et al._ [2021] G. Aad _et al._ (ATLAS), Phys. Lett. B 816, 136190 (2021), arXiv:2010.04019 [hep-ex]
  * Tumasyan  _et al._ [2023b] A. Tumasyan _et al._ (CMS), Phys. Lett. B 841, 137495 (2023b), arXiv:2205.05711 [hep-ex]
  * Aad  _et al._ [2023b] G. Aad _et al._ (ATLAS), Nature Phys. 19, 237 (2023b), arXiv:2004.10612 [hep-ex]
  * Sirunyan  _et al._ [2021b] A. M. Sirunyan _et al._ (CMS), Phys. Lett. B 812, 135992 (2021b), arXiv:2008.07013 [hep-ex]
  * Rahaman and Singh [2022] R. Rahaman and R. K. Singh, Nucl. Phys. B 984, 115984 (2022), arXiv:2109.09345 [hep-ph]
  * von Weizsacker [1934] C. F. von Weizsacker, Z. Phys. 88, 612 (1934)
  * Williams [1934] E. J. Williams, Phys. Rev. 45, 729 (1934)
  * Dawson [1985] S. Dawson, Nucl. Phys. B 249, 42 (1985)
  * Johnson  _et al._ [1987] P. W. Johnson, F. I. Olness, and W.-K. Tung, Phys. Rev. D 36, 291 (1987)
  * Fabbrichesi  _et al._ [2023d] M. Fabbrichesi, R. Floreanini, E. Gabrielli, and L. Marzola, (2023d), arXiv:2304.02403 [hep-ph]
  * Cervera-Lierta  _et al._ [2017] A. Cervera-Lierta, J. I. Latorre, J. Rojo, and L. Rottoli, SciPost Phys. 3, 036 (2017), arXiv:1703.02989 [hep-th]