Detection of mode-intrinsic quantum entanglement
Carlos E. Lopetegui Mathieu Isoard Nicolas Treps Mattia Walschaers

Abstract

Quantum correlations are at the core of the power of quantum information and are necessary to reach a quantum computational advantage. In the context of continuous-variable quantum systems, another necessary ressource for quantum advantages is non-Gaussianity. In this work, we propose a witness, based on previously known relations between metrological power and quantum correlations, to detect a strong form of entanglement that only non-Gaussian states possess and that cannot be undone by passive optical operations, i.e., entanglement in all mode bases. The strength of our witness is two-fold: it only requires measurements in one basis to check entanglement in any arbitrary mode basis; it can be made applicable experimentally using homodyne measurements and without requiring a full tomography of the state.

I Introduction

Major efforts are being conducted worldwide to harness the power of quantum mechanics for technological applications. The use of these systems for performing computations beyond the capabilities of classical processors stands as one of the most encouraging perspectives. A wide variety of physical platforms are being considered for this purpose Steffen et al. (2011); Acharya et al. (2023); Henriet et al. (2020); Evered et al. (2023); Kielpinski et al. (2002); Gouzien et al. (2023); Author Bartolucci et al. (2023); Bourassa et al. (2021), among which optical systems prove to be a promising framework given its potential for scalability Larsen et al. (2019); Asavanant et al. (2019).

The full power of optical quantum technologies is attained when we can coherently control the collective state of many modes. These collective states manifest entanglement, a unique feature of quantum systems that has challenged the common conceptions of physicists over the last century Einstein et al. (1935); Schrödinger (1935); Horodecki et al. (2009). Entanglement and quantum correlations in general, reflect in measurement statistics that are not consistent with classical correlations. Mathematically this implies that they cannot be explained using a local hidden variable model Bell (1964). In the context of quantum information processing and quantum computing entanglement has revealed to be an essential resource Acín et al. (2018); Jozsa and Linden (2003). However, in continuous variables (CV) systems, entanglement is not enough if the state is fully Gaussian, or if its associated Wigner function – which describes all the statistical properties of the quantum state – is positive. Indeed, it has been shown that negativity in the Wigner function is a necessary resource for realizing sampling protocols that cannot be efficiently simulated with classical resources Mari and Eisert (2012). This result significantly boosted the already ongoing efforts aimed at producing non-Gaussian states Ra et al. (2020); Endo et al. (2023); Neergaard-Nielsen et al. (2006); Wakui et al. (2007), with remarkable results achieved in microwave cavities and superconducting circuits in the production of exotic non-Gaussian states Kudra et al. (2022); Touzard et al. (2018).

Yet, recent results have pointed out that not any kind of quantum entanglement (even if the state is Wigner negative) is sufficient to render computations with CV quantum systems hard to simulate with a classical device. In particular, in Ref. Chabaud and Walschaers (2023) it was shown that in a certain family of sampling protocols a strong form of quantum entanglement is required: not passive separability, i.e., the fact that entanglement cannot be undone with optical passive transformations (beamsplitters and phase shifters). This kind of entanglement was previously introduced in Sperling et al. (2019), from a different perspective. They refer to it as mode-independent quantum entanglement and point out its potential usefulness for communication protocols. We refer to this kind of quantum correlations as mode-intrinsic entanglement . Finding and characterizing non-Gaussian states that exhibit this refined quantum feature and that can be produced in the lab are thus important steps in the domain of CV quantum information.

Many witnesses and criteria have been derived for detecting entanglement in Gaussian states Simon (2000); Duan et al. (2000); Giovannetti et al. (2003); van Loock and Furusawa (2003); Abiuso et al. (2021). All of them rely in some way in the knowledge of the covariance matrix. Once measured in one basis, such a matrix can be obtained in any mode basis by means of an orthogonal transformation. For this reason these methods are quite versatile for studying entanglement in every mode basis. However, Williamson’s decomposition Williamson (1936); Braunstein (2005) implies that there always exists a basis in which the covariance matrix is consistent with a separable Gaussian state, independently of the nature of the quantum state itself. Therefore, these methods are not suitable for demonstrating mode-intrinsic entanglement. Furthermore, it shows that mode-intrinsic entanglement is inherently a non-Gaussian feature.

Over the years many protocols have been proposed to detect entanglement in non-Gaussian states. If one has access to the density matrix of a given state it is possible to use the positive partial transpose (PPT) criterion Peres (1996) to detect quantum correlations. However, the reconstruction of the density matrix requires a full tomography of the state Lvovsky and Raymer (2009), a very demanding and time-consuming experimental task. Efforts to circumvent this issue have lead to the introduction of different methods to witness quantum correlations Shchukin and Vogel (2005); Rodó et al. (2008); Walborn et al. (2011); Saboia et al. (2011); Gessner et al. (2016); Yadin et al. (2021); Gärttner et al. (2023); Gao et al. (2024); Valido et al. (2014).

Among those there is a subset Gessner et al. (2016); Yadin et al. (2021) that rely on connections to quantum metrology Giovannetti et al. (2011); Pezzè et al. (2018); Tóth and Apellaniz (2014). These witnesses are based on the Fisher Information Liu et al. (2019), a quantity that reflects the highest sensitivity with which we can estimate a parameter of interest using measurements of an observable. These methods rely on the idea that quantum correlations can lead to metrological sensitivities that could not be achieved with only classical correlations. CV-specialized protocols based on these witnesses Lopetegui et al. (2022); Barral et al. (2023) have been shown to succeed in detecting quantum correlations in relevant non-Gaussian states where usual Gaussian criteria fail. This suggests that non-Gaussian entanglement is not only crucial in the framework of quantum computation but can also potentially offer advantages in sensing applications. This idea relates naturally to recent results Lupu-Gladstein et al. (2022) that link the negativity of some quasiprobability distributions to enhancements in quantum phase estimation in quantum optics.

All the methods we have mentioned so far to detect entanglement in non-Gaussian states are developed to detect entanglement once we have fixed a mode basis of interest. On the other hand, to study entanglement in every mode basis, a method that retains the versatility of the Gaussian protocols, but potentially succeeds where the former fails, is necessary. In this paper, we propose such a protocol based on metrological witnesses. We show that our protocol works for relevant non-Gaussian states that can be produced with state-of-the-art experiments, namely photon subtracted states Ra et al. (2020); Endo et al. (2023). Likewise, our witness also demonstrates effectiveness in multimode systems such as (non-Gaussian) cluster states – that can be obtained by de-gaussifying readily available large CV cluster states Larsen et al. (2019); Asavanant et al. (2019); Ra et al. (2020), which have been widely studied for their scalability in the pursuit of large-scale quantum computing. While sensitive to losses, the mode-intrinsic-entanglement witness that we suggest here can be applied experimentally through homodyne detection and the estimation of the corresponding classical Fisher information – using the Hellinger distance as proposed in Ref. Barral et al. (2023); it thus avoids full tomography of the state, making metrological witnesses a very efficient tool in the detection of mode-intrinsic entanglement and non-Gaussian features. Moreover, it proves useful for the less specific task of studying the entanglement properties of a given state beyond the basis in which measurements are performed, with higher accuracy than Gaussian methods.

The paper is organized as follows: in Sec. II we introduce the metrological entanglement witness derived in Ref. Gessner et al. (2016) and extend it to the detection of mode-intrinsic entanglement in Sec. III. Then, in Sec. IV, we apply our witness to experimentally relevant multimode states, namely photon subtracted states. In the last section (Sec. V), we consider a lower bound of the mode-intrinsic entanglement witness derived in the previous sections; we show that this witness is experimentally accessible through homodyne measurements; in particular, we test this lower bound witness on a two-photon subtracted state using simulated experimental data.

II Metrology based entanglement witnesses

In this section we introduce metrological witnesses of entanglement in which the protocol that we present is based Gessner et al. (2016). These witnesses rely on the analysis of a specific metrological task. We consider two parties, Alice and Bob, who share a state ρ^^𝜌\hat{\rho}over^ start_ARG italic_ρ end_ARG.

They are interested in jointly estimating a parameter κ𝜅\kappaitalic_κ that they encode through a unitary transformation induced by the sum of local generators H^=H^A+H^B^𝐻subscript^𝐻𝐴subscript^𝐻𝐵\hat{H}=\hat{H}_{A}+\hat{H}_{B}over^ start_ARG italic_H end_ARG = over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT + over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT (see Fig.1). In this case, the density matrix of the state transforms as ρ^κ=exp⁡(−i⁢κ⁢H^)⁢ρ^⁢exp⁡(i⁢κ⁢H^)subscript^𝜌𝜅𝑖𝜅^𝐻^𝜌𝑖𝜅^𝐻\hat{\rho}_{\kappa}=\exp(-i\kappa\hat{H})\,\hat{\rho}\,\exp(i\kappa\hat{H})over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_κ end_POSTSUBSCRIPT = roman_exp ( - italic_i italic_κ over^ start_ARG italic_H end_ARG ) over^ start_ARG italic_ρ end_ARG roman_exp ( italic_i italic_κ over^ start_ARG italic_H end_ARG ). The sensitivity with which they can estimate the parameter, given that they can measure some (possibly non local) observables M^A⁢(B)subscript^𝑀𝐴𝐵\hat{M}_{A(B)}over^ start_ARG italic_M end_ARG start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT, is upper bounded by the Fisher information associated to the probability distribution of outcomes of these observables.

The Fisher information is given by

F⁢(p⁢(MA,MB|κ))=∫MA,MB𝑑MA⁢𝑑MB(∂ℒ⁢(p⁢(MA,MB|κ))∂κ)2×p⁢(MA,MB|κ),𝐹𝑝subscript𝑀𝐴conditionalsubscript𝑀𝐵𝜅subscriptsubscript𝑀𝐴subscript𝑀𝐵differential-dsubscript𝑀𝐴differential-dsubscript𝑀𝐵superscriptℒ𝑝subscript𝑀𝐴conditionalsubscript𝑀𝐵𝜅𝜅2𝑝subscript𝑀𝐴conditionalsubscript𝑀𝐵𝜅\begin{split}F(p(M_{A},M_{B}|\kappa))=\int_{M_{A},M_{B}}\!\!\!\!\!\!dM_{A}\,dM% _{B}&\left(\frac{\partial\mathcal{L}(p(M_{A},M_{B}|\kappa))}{\partial\kappa}% \right)^{2}\\ &\times p(M_{A},M_{B}|\kappa),\end{split}start_ROW start_CELL italic_F ( italic_p ( italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT | italic_κ ) ) = ∫ start_POSTSUBSCRIPT italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_d italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT italic_d italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT end_CELL start_CELL ( divide start_ARG ∂ caligraphic_L ( italic_p ( italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT | italic_κ ) ) end_ARG start_ARG ∂ italic_κ end_ARG ) start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL × italic_p ( italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT | italic_κ ) , end_CELL end_ROW (1)

where ℒ⁢(p⁢(MA,MB|κ))=log⁡(p⁢(MA,MB|κ))ℒ𝑝subscript𝑀𝐴conditionalsubscript𝑀𝐵𝜅𝑝subscript𝑀𝐴conditionalsubscript𝑀𝐵𝜅\mathcal{L}(p(M_{A},M_{B}|\kappa))=\log(p(M_{A},M_{B}|\kappa))caligraphic_L ( italic_p ( italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT | italic_κ ) ) = roman_log ( italic_p ( italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT | italic_κ ) ), stands for the logarithmic likelihood, and p⁢(MA,MB|κ)=Tr⁢(ρ^κ⁢Π^A⁢B)𝑝subscript𝑀𝐴conditionalsubscript𝑀𝐵𝜅Trsubscript^𝜌𝜅subscript^Π𝐴𝐵p(M_{A},M_{B}|\kappa)={\rm Tr}\left(\hat{\rho}_{\kappa}\,\hat{\Pi}_{AB}\right)italic_p ( italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT | italic_κ ) = roman_Tr ( over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_κ end_POSTSUBSCRIPT over^ start_ARG roman_Π end_ARG start_POSTSUBSCRIPT italic_A italic_B end_POSTSUBSCRIPT ) is the joint probability distribution conditioned by the parameter κ𝜅\kappaitalic_κ, where Π^A⁢B=|MA,MB⟩⁢⟨MA,MB|subscript^Π𝐴𝐵ketsubscript𝑀𝐴subscript𝑀𝐵brasubscript𝑀𝐴subscript𝑀𝐵\hat{\Pi}_{AB}=|M_{A},M_{B}\rangle\langle M_{A},M_{B}|over^ start_ARG roman_Π end_ARG start_POSTSUBSCRIPT italic_A italic_B end_POSTSUBSCRIPT = | italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT ⟩ ⟨ italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT | is the projector into the outcomes (MA,MB)subscript𝑀𝐴subscript𝑀𝐵(M_{A},M_{B})( italic_M start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , italic_M start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT ) of M^A⊗M^Btensor-productsubscript^𝑀𝐴subscript^𝑀𝐵\hat{M}_{A}\otimes\hat{M}_{B}over^ start_ARG italic_M end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT ⊗ over^ start_ARG italic_M end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT. The ultimate sensitivity achievable in the estimation of κ𝜅\kappaitalic_κ is given by the quantum Fisher information (QFI) of the state, defined as the maximum Fisher information, optimized over all possible measurements, i.e.,

FQ⁢(ρ^,H^)=maxM^A,M^B⁡F⁢[Tr⁢(ρ^κ⁢Π^A⁢B)].subscript𝐹𝑄^𝜌^𝐻subscriptsubscript^𝑀𝐴subscript^𝑀𝐵𝐹delimited-[]Trsubscript^𝜌𝜅subscript^Π𝐴𝐵F_{Q}(\hat{\rho},\hat{H})=\max_{\hat{M}_{A},\hat{M}_{B}}F\left[{\rm Tr}\left(% \hat{\rho}_{\kappa}\,\hat{\Pi}_{AB}\right)\right].italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG , over^ start_ARG italic_H end_ARG ) = roman_max start_POSTSUBSCRIPT over^ start_ARG italic_M end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , over^ start_ARG italic_M end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT end_POSTSUBSCRIPT italic_F [ roman_Tr ( over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_κ end_POSTSUBSCRIPT over^ start_ARG roman_Π end_ARG start_POSTSUBSCRIPT italic_A italic_B end_POSTSUBSCRIPT ) ] . (2)

It can be shown that for pure states

FQ⁢(|ψ⟩⁢⟨ψ|,H^)=4⁢Var⁢(|ψ⟩⁢⟨ψ|,H^),subscript𝐹𝑄ket𝜓bra𝜓^𝐻4Varket𝜓bra𝜓^𝐻F_{Q}\left(|\psi\rangle\langle\psi|,\,\hat{H}\right)=4\text{Var}\left(|\psi% \rangle\langle\psi|,\hat{H}\right),italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( | italic_ψ ⟩ ⟨ italic_ψ | , over^ start_ARG italic_H end_ARG ) = 4 Var ( | italic_ψ ⟩ ⟨ italic_ψ | , over^ start_ARG italic_H end_ARG ) , (3)

where Var⁢(|ψ⟩⁢⟨ψ|,H^)Varket𝜓bra𝜓^𝐻\text{Var}\left(|\psi\rangle\langle\psi|,\hat{H}\right)Var ( | italic_ψ ⟩ ⟨ italic_ψ | , over^ start_ARG italic_H end_ARG ) is the variance of the generator H^^𝐻\hat{H}over^ start_ARG italic_H end_ARG in the state |ψ⟩ket𝜓|\psi\rangle| italic_ψ ⟩. In the case of mixed states there is also an analytical expression for the quantum Fisher information, but for the moment it suffices to say that FQ⁢(ρ^,H^)≤4⁢Var⁢(ρ^,H^)subscript𝐹𝑄^𝜌^𝐻4Var^𝜌^𝐻F_{Q}(\hat{\rho},\hat{H})\leq 4\text{Var}(\hat{\rho},\hat{H})italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG , over^ start_ARG italic_H end_ARG ) ≤ 4 Var ( over^ start_ARG italic_ρ end_ARG , over^ start_ARG italic_H end_ARG ). Another useful property for our purpose is the aditivity of the QFI, namely

FQ⁢(ρ^A⊗ρ^B,H^A+H^B)=FQ⁢(ρ^A,H^A)+FQ⁢(ρ^B,H^B),subscript𝐹𝑄tensor-productsubscript^𝜌𝐴subscript^𝜌𝐵subscript^𝐻𝐴subscript^𝐻𝐵subscript𝐹𝑄subscript^𝜌𝐴subscript^𝐻𝐴subscript𝐹𝑄subscript^𝜌𝐵subscript^𝐻𝐵F_{Q}(\hat{\rho}_{A}\otimes\hat{\rho}_{B},\hat{H}_{A}+\hat{H}_{B})=F_{Q}(\hat{% \rho}_{A},\hat{H}_{A})+F_{Q}(\hat{\rho}_{B},\hat{H}_{B}),italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT ⊗ over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT , over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT + over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT ) = italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT ) + italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT , over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT ) , (4)

where ρ^A⁢(B)=TrB⁢(A)⁢(ρ^)subscript^𝜌𝐴𝐵subscriptTr𝐵𝐴^𝜌\hat{\rho}_{A(B)}={\rm Tr}_{B(A)}(\hat{\rho})over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT = roman_Tr start_POSTSUBSCRIPT italic_B ( italic_A ) end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG ) is the reduced density matrix. Combining these two properties and the convexity of the quantum Fisher information, one can show that that if Alice and Bob share a separable state ρ^^𝜌\hat{\rho}over^ start_ARG italic_ρ end_ARG it leads to the inequality Gessner et al. (2016)

FQ⁢(ρ^,H^=H^A+H^B)≤4⁢(Var⁢(ρ^A,H^A)+Var⁢(ρ^B,H^B)).subscript𝐹𝑄^𝜌^𝐻subscript^𝐻𝐴subscript^𝐻𝐵4Varsubscript^𝜌𝐴subscript^𝐻𝐴Varsubscript^𝜌𝐵subscript^𝐻𝐵F_{Q}(\hat{\rho},\hat{H}=\hat{H}_{A}+\hat{H}_{B})\leq 4\left(\text{Var}(\hat{% \rho}_{A},\hat{H}_{A})+\text{Var}(\hat{\rho}_{B},\hat{H}_{B})\right).italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG , over^ start_ARG italic_H end_ARG = over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT + over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT ) ≤ 4 ( Var ( over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT ) + Var ( over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT , over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT ) ) . (5)

This inequality can only be violated by entangled states and can be thus used as a witness. In practical terms this implies that if they share an entangled state they better join forces than trying to estimate independently the parameter. One would be tempted to say that this implies that entanglement offers a metrological advantage in this setup; yet the situation is a bit more subtle than that, because we are comparing the estimation they can do jointly with respect to what they can estimate with their noisy reduced states. However, this does not imply that there is no better pure single mode state they could independently use and get each a better Fisher information.

The full power of this kind of inequalities can be exploited when we optimize over a set of allowed choices of generators. For any H^A⁢(B)subscript^𝐻𝐴𝐵\hat{H}_{A(B)}over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT, the violation of the metrological inequality provides a sufficient condition for entanglement. Yet, if a violation is not attained it might easily be that there is a different choice of generators that produces one. If Alice and Bob each have access to a set of local generators 𝒜A⁢(B)={H^A⁢(B)(1),…,H^A⁢(B)(NA⁢(B))}subscript𝒜𝐴𝐵superscriptsubscript^𝐻𝐴𝐵1…superscriptsubscript^𝐻𝐴𝐵subscript𝑁𝐴𝐵\mathcal{A}_{A(B)}=\{\hat{H}_{A(B)}^{(1)},...,\hat{H}_{A(B)}^{(N_{A(B)})}\}caligraphic_A start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT = { over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( 1 ) end_POSTSUPERSCRIPT , … , over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_N start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT }, where NA⁢(B)subscript𝑁𝐴𝐵N_{A(B)}italic_N start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT is the number of generators in each set, they can choose to implement generators given by H^A⁢(B)=𝐜A⁢(B)T⁢𝒜A⁢(B)subscript^𝐻𝐴𝐵superscriptsubscript𝐜𝐴𝐵𝑇subscript𝒜𝐴𝐵\hat{H}_{A(B)}=\mathbf{c}_{A(B)}^{T}\mathcal{A}_{A(B)}over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT = bold_c start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_T end_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT. We group 𝐜A⁢(B)subscript𝐜𝐴𝐵\mathbf{c}_{A(B)}bold_c start_POSTSUBSCRIPT italic_A ( italic_B ) end_POSTSUBSCRIPT in a single vector 𝐜={𝐜A,𝐜B}𝐜subscript𝐜𝐴subscript𝐜𝐵\mathbf{c}=\{\mathbf{c}_{A},\mathbf{c}_{B}\}bold_c = { bold_c start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , bold_c start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT }, and similarly 𝒜={𝒜A,𝒜B}𝒜subscript𝒜𝐴subscript𝒜𝐵\mathcal{A}=\{\mathcal{A}_{A},\mathcal{A}_{B}\}caligraphic_A = { caligraphic_A start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT , caligraphic_A start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT }. This allows to write the Quantum Fisher Information in a matrix form as Gessner et al. (2016)

FQ⁢(ρ^,𝐜T⁢𝒜)=𝐜T⁢𝒬ρ^𝒜⁢𝐜,subscript𝐹𝑄^𝜌superscript𝐜𝑇𝒜superscript𝐜𝑇superscriptsubscript𝒬^𝜌𝒜𝐜F_{Q}(\hat{\rho},\mathbf{c}^{T}\mathcal{A})=\mathbf{c}^{T}\mathcal{Q}_{\hat{% \rho}}^{\mathcal{A}}\mathbf{c},italic_F start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG , bold_c start_POSTSUPERSCRIPT italic_T end_POSTSUPERSCRIPT caligraphic_A ) = bold_c start_POSTSUPERSCRIPT italic_T end_POSTSUPERSCRIPT caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT bold_c , (6)

where 𝒬ρ^𝒜superscriptsubscript𝒬^𝜌𝒜\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}}caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT is the Quantum Fisher Information matrix and its elements are given, provided the spectral decomposition ρ^=∑kpk⁢|ψk⟩⁢⟨ψk|^𝜌subscript𝑘subscript𝑝𝑘ketsubscript𝜓𝑘brasubscript𝜓𝑘\hat{\rho}=\sum_{k}p_{k}|\psi_{k}\rangle\langle\psi_{k}|over^ start_ARG italic_ρ end_ARG = ∑ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT italic_p start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT | italic_ψ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT ⟩ ⟨ italic_ψ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT |, by

(𝒬ρ^𝒜)i,j(m,n)=∑k,l(pk−pl)2pk+pl⁢⟨ψk|H^i(m)|ψl⟩⁢⟨ψl|H^j(n)|ψk⟩,superscriptsubscriptsuperscriptsubscript𝒬^𝜌𝒜𝑖𝑗𝑚𝑛subscript𝑘𝑙superscriptsubscript𝑝𝑘subscript𝑝𝑙2subscript𝑝𝑘subscript𝑝𝑙quantum-operator-productsubscript𝜓𝑘superscriptsubscript^𝐻𝑖𝑚subscript𝜓𝑙quantum-operator-productsubscript𝜓𝑙superscriptsubscript^𝐻𝑗𝑛subscript𝜓𝑘\left(\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}}\right)_{i,j}^{(m,n)}=\sum_{k,l}% \frac{(p_{k}-p_{l})^{2}}{p_{k}+p_{l}}\langle\psi_{k}|\hat{H}_{i}^{(m)}|\psi_{l% }\rangle\langle\psi_{l}|\hat{H}_{j}^{(n)}|\psi_{k}\rangle,( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT ) start_POSTSUBSCRIPT italic_i , italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_m , italic_n ) end_POSTSUPERSCRIPT = ∑ start_POSTSUBSCRIPT italic_k , italic_l end_POSTSUBSCRIPT divide start_ARG ( italic_p start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT - italic_p start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT end_ARG start_ARG italic_p start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT + italic_p start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT end_ARG ⟨ italic_ψ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT | over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_m ) end_POSTSUPERSCRIPT | italic_ψ start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT ⟩ ⟨ italic_ψ start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT | over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_n ) end_POSTSUPERSCRIPT | italic_ψ start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT ⟩ , (7)

where the sum extends over all pairs k,l𝑘𝑙k,litalic_k , italic_l such that pk+pl≠0subscript𝑝𝑘subscript𝑝𝑙0p_{k}+p_{l}\neq 0italic_p start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT + italic_p start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT ≠ 0. The state ρ^^𝜌\hat{\rho}over^ start_ARG italic_ρ end_ARG is typically considered to be the one after the implementation of the evolution. Yet here we consider only small values of the parameters, so that the spectral decomposition considered is that of the probe state, and all higher orders on the value of the parameters are omited.

The right hand side of inequality (5) on the other hand corresponds to the variance of the generator in the product state of the local reduced states, i.e., Π⁢(ρ^)=ρ^A⊗ρ^BΠ^𝜌tensor-productsubscript^𝜌𝐴subscript^𝜌𝐵\Pi(\hat{\rho})=\hat{\rho}_{A}\otimes\hat{\rho}_{B}roman_Π ( over^ start_ARG italic_ρ end_ARG ) = over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT ⊗ over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT. This term can be expressed also in a matrix form by means of the covariance matrix of the generators defined by

Γi,j(m,n)=C⁢o⁢v⁢(H^i(m),H^j(n))Π⁢(ρ^).superscriptsubscriptΓ𝑖𝑗𝑚𝑛𝐶𝑜𝑣subscriptsuperscriptsubscript^𝐻𝑖𝑚superscriptsubscript^𝐻𝑗𝑛Π^𝜌\Gamma_{i,j}^{(m,n)}=Cov(\hat{H}_{i}^{(m)},\hat{H}_{j}^{(n)})_{\Pi(\hat{\rho})}.roman_Γ start_POSTSUBSCRIPT italic_i , italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_m , italic_n ) end_POSTSUPERSCRIPT = italic_C italic_o italic_v ( over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_m ) end_POSTSUPERSCRIPT , over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_n ) end_POSTSUPERSCRIPT ) start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT . (8)

The variance of the generator for a particular choice of 𝐜𝐜\mathbf{c}bold_c is given by

Var⁢(𝐜AT⁢𝒜A)+Var⁢(𝐜BT⁢𝒜B)=𝐜T⁢ΓΠ⁢(ρ^)𝒜⁢𝐜.Varsuperscriptsubscript𝐜𝐴Tsubscript𝒜𝐴Varsuperscriptsubscript𝐜𝐵Tsubscript𝒜𝐵superscript𝐜𝑇superscriptsubscriptΓΠ^𝜌𝒜𝐜\text{Var}(\mathbf{c}_{A}^{\rm\scriptscriptstyle T}\mathcal{A}_{A})+\text{Var}% (\mathbf{c}_{B}^{\rm\scriptscriptstyle T}\mathcal{A}_{B})=\mathbf{c}^{T}\Gamma% _{\Pi(\hat{\rho})}^{\mathcal{A}}\mathbf{c}.Var ( bold_c start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_A end_POSTSUBSCRIPT ) + Var ( bold_c start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_B end_POSTSUBSCRIPT ) = bold_c start_POSTSUPERSCRIPT italic_T end_POSTSUPERSCRIPT roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT bold_c . (9)

Combining (6) and (9) into (5) we get that, for any separable state

𝐜T⁢(𝒬ρ^𝒜−4⁢ΓΠ⁢(ρ^)𝒜)⁢𝐜≤0,∀𝐜.superscript𝐜𝑇superscriptsubscript𝒬^𝜌𝒜4superscriptsubscriptΓΠ^𝜌𝒜𝐜0for-all𝐜\mathbf{c}^{T}(\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}}-4\Gamma_{\Pi(\hat{\rho})% }^{\mathcal{A}})\mathbf{c}\leq 0,\qquad\forall\mathbf{c}.bold_c start_POSTSUPERSCRIPT italic_T end_POSTSUPERSCRIPT ( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT - 4 roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT ) bold_c ≤ 0 , ∀ bold_c . (10)

Given that this condition is satisfied for all choices of 𝐜𝐜\mathbf{c}bold_c, it can be rewritten as

λmax⁢(𝒬ρ^𝒜−4⁢ΓΠ⁢(ρ^)𝒜)≤0,subscript𝜆maxsuperscriptsubscript𝒬^𝜌𝒜4superscriptsubscriptΓΠ^𝜌𝒜0\lambda_{\rm\scriptscriptstyle max}\left(\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}% }-4\Gamma_{\Pi(\hat{\rho})}^{\mathcal{A}}\right)\leq 0,italic_λ start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT ( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT - 4 roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT ) ≤ 0 , (11)

where λmax⁢(M)subscript𝜆max𝑀\lambda_{\rm\scriptscriptstyle max}(M)italic_λ start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT ( italic_M ) stands for the maximum eigenvalue of M𝑀Mitalic_M. This inequality provides us a witness of entanglement

Eρ𝒜≡λmax⁢(𝒬ρ^𝒜−4⁢ΓΠ⁢(ρ^)𝒜).superscriptsubscript𝐸𝜌𝒜subscript𝜆maxsuperscriptsubscript𝒬^𝜌𝒜4superscriptsubscriptΓΠ^𝜌𝒜E_{\rho}^{\mathcal{A}}\equiv\lambda_{\rm\scriptscriptstyle max}\left(\mathcal{% Q}_{\hat{\rho}}^{\mathcal{A}}-4\Gamma_{\Pi(\hat{\rho})}^{\mathcal{A}}\right).italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT ≡ italic_λ start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT ( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT - 4 roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT ) . (12)

If this witness takes possitive values we can certify the presence of entanglement. Notably, if the set of generators is a complete set and the state is a pure state the latter becomes not just a sufficient but also a necessary condition Gessner et al. (2016).

All the above results can be straightforwardly generalized to multimode systems and to multipartite entanglement detection. In this case, the m𝑚mitalic_m parties have access to a set of local generators 𝒜m={H^m(1),…,H^m(Nm)}subscript𝒜𝑚subscriptsuperscript^𝐻1𝑚…subscriptsuperscript^𝐻subscript𝑁𝑚𝑚\mathcal{A}_{m}=\{\hat{H}^{(1)}_{m},\ldots,\hat{H}^{(N_{m})}_{m}\}caligraphic_A start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT = { over^ start_ARG italic_H end_ARG start_POSTSUPERSCRIPT ( 1 ) end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT , … , over^ start_ARG italic_H end_ARG start_POSTSUPERSCRIPT ( italic_N start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT ) end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT }. The elements of the Quantum Fisher Information matrix and the covariance matrix can be then computed from equations (7) and (8), with Π⁢(ρ^)=ρ^𝒫1⊗…⊗ρ^𝒫pΠ^𝜌tensor-productsubscript^𝜌subscript𝒫1…subscript^𝜌subscript𝒫𝑝\Pi(\hat{\rho})=\hat{\rho}_{\scriptscriptstyle\mathcal{P}_{1}}\otimes\ldots% \otimes\hat{\rho}_{\scriptscriptstyle\mathcal{P}_{p}}roman_Π ( over^ start_ARG italic_ρ end_ARG ) = over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT caligraphic_P start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT end_POSTSUBSCRIPT ⊗ … ⊗ over^ start_ARG italic_ρ end_ARG start_POSTSUBSCRIPT caligraphic_P start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT end_POSTSUBSCRIPT, where we split the m𝑚mitalic_m-multimode system into p𝑝pitalic_p partitions 𝒫isubscript𝒫𝑖\mathcal{P}_{i}caligraphic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT to detect p𝑝pitalic_p-partite entanglement.

III Witness mode-intrinsic entanglement

A key aspect of non-Gaussian entanglement lies in its robustness against mode basis change. Contrary to Gaussian states whose entanglement can always be destroyed through a passive linear optics transformation, no such transformation can entirely remove the entanglement present in some non-Gaussian states; these states are said to be not passively separable and we call this type of entanglement mode-intrinsic entanglement. To detect such entanglement we propose to use the metrological witness (11) introduced in Sec.II as follows: a state will be mode-intrinsic entangled if in every mode basis Eρ𝒜>0superscriptsubscript𝐸𝜌𝒜0E_{\rho}^{\mathcal{A}}>0italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT > 0.

In the following we focus on mutlimode quantum optical systems Fabre and Treps (2020). In CV quantum optics the natural and fundamental observables are the amplitude and phase quadratures q^ksubscript^𝑞𝑘\hat{q}_{k}over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT and p^ksubscript^𝑝𝑘\hat{p}_{k}over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT, which are related to creation and annihilation operators of the corresponding bosonic field by

a^k=q^k+i⁢p^k2,a^k†=q^k−i⁢p^k2,∀k∈{1,…,m},formulae-sequencesubscript^𝑎𝑘subscript^𝑞𝑘𝑖subscript^𝑝𝑘2formulae-sequencesuperscriptsubscript^𝑎𝑘†subscript^𝑞𝑘𝑖subscript^𝑝𝑘2for-all𝑘1…𝑚\hat{a}_{k}=\frac{\hat{q}_{k}+i\hat{p}_{k}}{2},\quad\hat{a}_{k}^{\dagger}=% \frac{\hat{q}_{k}-i\hat{p}_{k}}{2},\quad\forall k\in\{1,\ldots,m\},over^ start_ARG italic_a end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT = divide start_ARG over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT + italic_i over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT end_ARG start_ARG 2 end_ARG , over^ start_ARG italic_a end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT start_POSTSUPERSCRIPT † end_POSTSUPERSCRIPT = divide start_ARG over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT - italic_i over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT end_ARG start_ARG 2 end_ARG , ∀ italic_k ∈ { 1 , … , italic_m } , (13)

where a^ksubscript^𝑎𝑘\hat{a}_{k}over^ start_ARG italic_a end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT and a^k†subscriptsuperscript^𝑎†𝑘\hat{a}^{\dagger}_{k}over^ start_ARG italic_a end_ARG start_POSTSUPERSCRIPT † end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT are, respectively, the creation and annihilation operators corresponding to mode kt⁢hsuperscript𝑘𝑡ℎk^{th}italic_k start_POSTSUPERSCRIPT italic_t italic_h end_POSTSUPERSCRIPT out of m𝑚mitalic_m available modes. Quadrature operators satisfy the canonical commutation relation [q^k,p^l]=2⁢i⁢δk,lsubscript^𝑞𝑘subscript^𝑝𝑙2𝑖subscript𝛿𝑘𝑙[\hat{q}_{k},\hat{p}_{l}]=2i\,\delta_{k,l}[ over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT ] = 2 italic_i italic_δ start_POSTSUBSCRIPT italic_k , italic_l end_POSTSUBSCRIPT. The statistics of these operators can be recovered experimentally through homodyne measurements Leonhardt and Paul (1995), making them a judicious choice to construct a set of generators that the m𝑚mitalic_m parties involved in the multimode system can experimentally implement. To that aim, we introduce the following generators of order N𝑁Nitalic_N

H^i𝒮=𝒮⁢(q^1k1(i)⁢…⁢q^mkm(i)⁢p^1km+1(i)⁢…⁢p^mk2⁢m(i)),with⁢∑j=12⁢mkj(i)=N,formulae-sequencesuperscriptsubscript^𝐻𝑖𝒮𝒮superscriptsubscript^𝑞1superscriptsubscript𝑘1𝑖…superscriptsubscript^𝑞𝑚superscriptsubscript𝑘𝑚𝑖superscriptsubscript^𝑝1superscriptsubscript𝑘𝑚1𝑖…superscriptsubscript^𝑝𝑚superscriptsubscript𝑘2𝑚𝑖withsuperscriptsubscript𝑗12𝑚superscriptsubscript𝑘𝑗𝑖𝑁\hat{H}_{i}^{\mathcal{S}}=\mathcal{S}\left(\hat{q}_{1}^{k_{1}^{(i)}}\ldots\hat% {q}_{m}^{k_{m}^{(i)}}\,\hat{p}_{1}^{k_{m+1}^{(i)}}\ldots\hat{p}_{m}^{k_{2m}^{(% i)}}\right),\,\text{with}\,\sum_{j=1}^{2m}k_{j}^{(i)}=N,over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT = caligraphic_S ( over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT … over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_m + 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT … over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT 2 italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) , with ∑ start_POSTSUBSCRIPT italic_j = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 italic_m end_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT = italic_N , (14)

where 𝒮⁢(H^)𝒮^𝐻\mathcal{S}(\hat{H})caligraphic_S ( over^ start_ARG italic_H end_ARG ) is the symmetrised version of H^^𝐻\hat{H}over^ start_ARG italic_H end_ARG such that H^i𝒮superscriptsubscript^𝐻𝑖𝒮\hat{H}_{i}^{\mathcal{S}}over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT is Hermitian. As one can see from the above expression, H^i𝒮superscriptsubscript^𝐻𝑖𝒮\hat{H}_{i}^{\mathcal{S}}over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT contains products of N𝑁Nitalic_N quadrature operators and is entirely described by the combination (k1(i),…,k2⁢m(i))superscriptsubscript𝑘1𝑖…superscriptsubscript𝑘2𝑚𝑖(k_{1}^{(i)},\ldots,k_{2m}^{(i)})( italic_k start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT , … , italic_k start_POSTSUBSCRIPT 2 italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT ), where each kj(i)superscriptsubscript𝑘𝑗𝑖k_{j}^{(i)}italic_k start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT corresponds to the power of the quadrature operators in the decomposition (14) and whose sum equals N𝑁Nitalic_N. Likewise, this type of generators can be non-local. Here, locality is defined as follows: local generators only belong to one subsystem, while non-local generators can be shared between multiple subsystems; for instance the quadrature operator p^12superscriptsubscript^𝑝12\hat{p}_{1}^{2}over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT is a local generator, while p^1⁢p^2subscript^𝑝1subscript^𝑝2\hat{p}_{1}\,\hat{p}_{2}over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT is a non-local generator shared between subsystems 1 and 2. Since the metrological protocol described in the previous section is only based on local generators, it may seem unnecessary to include non-local generators. However, they will play a crucial role to correctly implement the mode basis change as we will see below.

Then, we denote by 𝒜N,m(NL)={H^1𝒮,…⁢H^ℓ⁢(N,m)𝒮}subscriptsuperscript𝒜NL𝑁𝑚superscriptsubscript^𝐻1𝒮…superscriptsubscript^𝐻ℓ𝑁𝑚𝒮\mathcal{A}^{\rm\scriptscriptstyle(NL)}_{\scriptscriptstyle N,% \scriptscriptstyle m}=\{\hat{H}_{1}^{\mathcal{S}},\ldots\hat{H}_{\ell(N,m)}^{% \mathcal{S}}\}caligraphic_A start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT = { over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT , … over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT roman_ℓ ( italic_N , italic_m ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT } the set of all (local and non-local) symmetrised generators of the form (14) up to order N𝑁Nitalic_N, in a m𝑚mitalic_m-mode system. The superscript (NL) stands for non-local and

ℓ⁢(N,m)=∑k=1N(2⁢m+k−1k)ℓ𝑁𝑚superscriptsubscript𝑘1𝑁matrix2𝑚𝑘1𝑘\ell(N,m)=\sum_{k=1}^{N}\begin{pmatrix}2m+k-1\\ k\end{pmatrix}roman_ℓ ( italic_N , italic_m ) = ∑ start_POSTSUBSCRIPT italic_k = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT ( start_ARG start_ROW start_CELL 2 italic_m + italic_k - 1 end_CELL end_ROW start_ROW start_CELL italic_k end_CELL end_ROW end_ARG ) (15)

is the length of the set. As an illustration, for m=2𝑚2m=2italic_m = 2 and N=2𝑁2N=2italic_N = 2 this set of operators reads

𝒜2,2(NL)={q^1,p^1,q^2,p^2,q^12,12(q^1p^1+p^1q^1),p^12,q^22,12(q^2p^2+p^2q^2),p^22,q^1q^2,p^1p^2,q^1p^2,p^1q^2}.subscriptsuperscript𝒜NL22subscript^𝑞1subscript^𝑝1subscript^𝑞2subscript^𝑝2superscriptsubscript^𝑞1212subscript^𝑞1subscript^𝑝1subscript^𝑝1subscript^𝑞1superscriptsubscript^𝑝12superscriptsubscript^𝑞2212subscript^𝑞2subscript^𝑝2subscript^𝑝2subscript^𝑞2superscriptsubscript^𝑝22subscript^𝑞1subscript^𝑞2subscript^𝑝1subscript^𝑝2subscript^𝑞1subscript^𝑝2subscript^𝑝1subscript^𝑞2\begin{split}\mathcal{A}^{\rm\scriptscriptstyle(NL)}_{2,2}=&\left\{\hat{q}_{1}% ,\hat{p}_{1},\hat{q}_{2},\hat{p}_{2},\hat{q}_{1}^{2},\frac{1}{2}(\hat{q}_{1}% \hat{p}_{1}+\hat{p}_{1}\hat{q}_{1}),\hat{p}_{1}^{2},\hat{q}_{2}^{2},\right.\\ &\left.\frac{1}{2}(\hat{q}_{2}\hat{p}_{2}+\hat{p}_{2}\hat{q}_{2}),\hat{p}_{2}^% {2},\hat{q}_{1}\hat{q}_{2},\hat{p}_{1}\hat{p}_{2},\hat{q}_{1}\hat{p}_{2},\hat{% p}_{1}\hat{q}_{2}\right\}.\end{split}start_ROW start_CELL caligraphic_A start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT start_POSTSUBSCRIPT 2 , 2 end_POSTSUBSCRIPT = end_CELL start_CELL { over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT , divide start_ARG 1 end_ARG start_ARG 2 end_ARG ( over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ) , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT , over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT , end_CELL end_ROW start_ROW start_CELL end_CELL start_CELL divide start_ARG 1 end_ARG start_ARG 2 end_ARG ( over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT + over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT , over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT } . end_CELL end_ROW (16)

Following the above notation, we also denote by 𝒜N,m(L)subscriptsuperscript𝒜L𝑁𝑚\mathcal{A}^{\rm\scriptscriptstyle(L)}_{\scriptscriptstyle N,m}caligraphic_A start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT the set of symmetrised generators of order N𝑁Nitalic_N containing only local generators, in an m𝑚mitalic_m-mode system. In the case of the set (16), this amounts to remove the four generators q^1⁢q^2,p^1⁢p^2,q^1⁢p^2subscript^𝑞1subscript^𝑞2subscript^𝑝1subscript^𝑝2subscript^𝑞1subscript^𝑝2\hat{q}_{1}\hat{q}_{2},\hat{p}_{1}\hat{p}_{2},\hat{q}_{1}\hat{p}_{2}over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT , over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT and p^1⁢q^2subscript^𝑝1subscript^𝑞2\hat{p}_{1}\hat{q}_{2}over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT.

Now that we have introduced the set of local and non-local generators ordered by N𝑁Nitalic_N, we can use them to test mode-intrinsic entanglement of a given quantum state of light. To proceed, we probe entanglement using the metrological witness in every mode bases. It is worth stressing here that it is not necessary to measure the quantum Fisher information matrix 𝒬ρ^𝒜superscriptsubscript𝒬^𝜌𝒜\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}}caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT as well as the covariance matrix of the local reduced state ΓΠ⁢(ρ^)𝒜superscriptsubscriptΓΠ^𝜌𝒜\Gamma_{\Pi(\hat{\rho})}^{\mathcal{A}}roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT in all bases. Indeed, once these matrices are computed or measured in one basis, one can calculate them in another basis in post-processing. To illustrate this important point, let us consider the unitary matrix U⁢(ϑ)𝑈bold-italic-ϑU(\bm{\vartheta})italic_U ( bold_italic_ϑ ) parameterized by some parameters contained in a vector ϑbold-italic-ϑ\bm{\vartheta}bold_italic_ϑ and which enables to go from one set of generators 𝒜N,m(NL)superscriptsubscript𝒜𝑁𝑚NL\mathcal{A}_{\scriptscriptstyle N,m}^{\rm\scriptscriptstyle(NL)}caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT in a given basis to another set of generators (𝒜N,m(NL))′superscriptsuperscriptsubscript𝒜𝑁𝑚NL′(\mathcal{A}_{\scriptscriptstyle N,m}^{\rm\scriptscriptstyle(NL)})^{\prime}( caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT in another basis, i.e., (𝒜N,m(NL))′=U⁢(ϑ)⁢𝒜N,m(NL)superscriptsuperscriptsubscript𝒜𝑁𝑚NL′𝑈bold-italic-ϑsuperscriptsubscript𝒜𝑁𝑚NL(\mathcal{A}_{\scriptscriptstyle N,m}^{\rm\scriptscriptstyle(NL)})^{\prime}=U(% \bm{\vartheta})\,\mathcal{A}_{\scriptscriptstyle N,m}^{\rm\scriptscriptstyle(% NL)}( caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT = italic_U ( bold_italic_ϑ ) caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT. The exact form of U⁢(ϑ)𝑈bold-italic-ϑU(\bm{\vartheta})italic_U ( bold_italic_ϑ ) will be detailed below after Eq. (21). At the level of the quantum Fisher information matrix the change of basis gives

(𝒬ρ^𝒜N,m(NL))′=U⁢(ϑ)⁢𝒬ρ^𝒜N,m(NL)⁢UT⁢(ϑ).superscriptsuperscriptsubscript𝒬^𝜌superscriptsubscript𝒜𝑁𝑚NL′𝑈bold-italic-ϑsuperscriptsubscript𝒬^𝜌superscriptsubscript𝒜𝑁𝑚NLsuperscript𝑈𝑇bold-italic-ϑ\left(\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}\right)^{\prime}=U(\bm{\vartheta})\,\mathcal{Q}_{\hat{% \rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm\scriptscriptstyle(NL)}}U^{T}% (\bm{\vartheta}).( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT = italic_U ( bold_italic_ϑ ) caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT italic_U start_POSTSUPERSCRIPT italic_T end_POSTSUPERSCRIPT ( bold_italic_ϑ ) . (17)

Therefore, knowing how the non-local set of generators transforms under a mode basis change allows to find the transformation at the level of the matrix 𝒬ρ^𝒜superscriptsubscript𝒬^𝜌𝒜\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}}caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT (and similarly for ΓΠ⁢(ρ^)𝒜superscriptsubscriptΓΠ^𝜌𝒜\Gamma_{\Pi(\hat{\rho})}^{\mathcal{A}}roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT) used in Eq. (10), and thus allows to find the metrological witness in the new basis.

There are some subtleties though: the basis change must include non-local operators to be correct. Indeed, let us take a simple example to illustrate this important point: the generator p^12superscriptsubscript^𝑝12\hat{p}_{1}^{2}over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT for a two-mode subsystem will be mapped onto (cos⁡θ⁢p^1′+sin⁡θ⁢p^2′)2superscript𝜃superscriptsubscript^𝑝1′𝜃superscriptsubscript^𝑝2′2(\cos\theta\hat{p}_{1}^{\prime}+\sin\theta\hat{p}_{2}^{\prime})^{2}( roman_cos italic_θ over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT + roman_sin italic_θ over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT when applying a basis change, and it clearly contains the non-local term p^1′⁢p^2′superscriptsubscript^𝑝1′superscriptsubscript^𝑝2′\hat{p}_{1}^{\prime}\,\hat{p}_{2}^{\prime}over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT. Therefore, we need to keep such a non-local term to ensure that matrix 𝒬ρ^𝒜superscriptsubscript𝒬^𝜌𝒜\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}}caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A end_POSTSUPERSCRIPT in Eq. (17) correctly transforms under the mode basis change U⁢(ϑ)𝑈bold-italic-ϑU(\bm{\vartheta})italic_U ( bold_italic_ϑ ).

However, once each matrix 𝒬ρ^𝒜N,m(NL)superscriptsubscript𝒬^𝜌superscriptsubscript𝒜𝑁𝑚NL\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT and ΓΠ⁢(ρ^)𝒜N,m(NL)superscriptsubscriptΓΠ^𝜌superscriptsubscript𝒜𝑁𝑚NL\Gamma_{\Pi(\hat{\rho})}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT in Eq. (10) have been found in the new basis, one should then drop all the rows and columns of these matrices associated to non-local generators. Indeed, the metrological protocol as described in the previous section consists of taking one local generator for each party of the multimode state, thus excluding non-local generators shared between subsystems. Keeping this important point in mind, the overall procedure then reads as follows:

1.

Compute the quantum Fisher information matrix 𝒬ρ^𝒜N,m(NL)superscriptsubscript𝒬^𝜌superscriptsubscript𝒜𝑁𝑚NL\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT for the whole set of operators 𝒜N,m(NL)superscriptsubscript𝒜𝑁𝑚NL\mathcal{A}_{\scriptscriptstyle N,m}^{\rm\scriptscriptstyle(NL)}caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT. When the state is pure this matrix equals four times the extended covariance matrix Γρ^𝒜N,m(NL)superscriptsubscriptΓ^𝜌superscriptsubscript𝒜𝑁𝑚NL\Gamma_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}roman_Γ start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT (whose definition has been given in Eq. (8) with Π⁢(ρ^)Π^𝜌\Pi(\hat{\rho})roman_Π ( over^ start_ARG italic_ρ end_ARG ) replaced by ρ^^𝜌\hat{\rho}over^ start_ARG italic_ρ end_ARG).

Then, compute Γρ^𝒜N,m(NL)superscriptsubscriptΓ^𝜌superscriptsubscript𝒜𝑁𝑚NL\Gamma_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}roman_Γ start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT. Both matrices have dimension ℓ⁢(N,m)×ℓ⁢(N,m)ℓ𝑁𝑚ℓ𝑁𝑚\ell(N,m)\times\ell(N,m)roman_ℓ ( italic_N , italic_m ) × roman_ℓ ( italic_N , italic_m ).

2.

Make the change of basis according to Eq. (17) to obtain (𝒬ρ^𝒜N,m(NL))′superscriptsuperscriptsubscript𝒬^𝜌superscriptsubscript𝒜𝑁𝑚NL′\left(\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}\right)^{\prime}( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT and (Γρ^𝒜N,m(NL))′superscriptsuperscriptsubscriptΓ^𝜌superscriptsubscript𝒜𝑁𝑚NL′\left(\Gamma_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(NL)}}\right)^{\prime}( roman_Γ start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_NL ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT.

3.

Remove the rows and columns related to non-local operators. The final matrices denoted by (𝒬ρ^𝒜N,m(L))′superscriptsuperscriptsubscript𝒬^𝜌superscriptsubscript𝒜𝑁𝑚L′\left(\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(L)}}\right)^{\prime}( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT and (Γρ^𝒜N,m(L))′superscriptsuperscriptsubscriptΓ^𝜌superscriptsubscript𝒜𝑁𝑚L′\left(\Gamma_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(L)}}\right)^{\prime}( roman_Γ start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT will have dimension ℓloc⁢(N,m)×ℓloc⁢(N,m)subscriptℓloc𝑁𝑚subscriptℓloc𝑁𝑚\ell_{\rm\scriptscriptstyle loc}(N,m)\times\ell_{\rm\scriptscriptstyle loc}(N,m)roman_ℓ start_POSTSUBSCRIPT roman_loc end_POSTSUBSCRIPT ( italic_N , italic_m ) × roman_ℓ start_POSTSUBSCRIPT roman_loc end_POSTSUBSCRIPT ( italic_N , italic_m ), with ℓloc⁢(N,m)=m⁢N⁢(N+3)/2subscriptℓloc𝑁𝑚𝑚𝑁𝑁32\ell_{\rm\scriptscriptstyle loc}(N,m)=m\,N\,(N+3)/2roman_ℓ start_POSTSUBSCRIPT roman_loc end_POSTSUBSCRIPT ( italic_N , italic_m ) = italic_m italic_N ( italic_N + 3 ) / 2.

4.

Compute (ΓΠ⁢(ρ^)𝒜N(L))′superscriptsuperscriptsubscriptΓΠ^𝜌superscriptsubscript𝒜𝑁L′\left(\Gamma_{\Pi(\hat{\rho})}^{\mathcal{A}_{\scriptscriptstyle N}^{\rm% \scriptscriptstyle(L)}}\right)^{\prime}( roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT. This computation amounts to keep the elements of (Γρ^𝒜N(L))′superscriptsuperscriptsubscriptΓ^𝜌superscriptsubscript𝒜𝑁L′\left(\Gamma_{\hat{\rho}}^{\mathcal{A}_{\scriptscriptstyle N}^{\rm% \scriptscriptstyle(L)}}\right)^{\prime}( roman_Γ start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT which only involve two generators H^i𝒮superscriptsubscript^𝐻𝑖𝒮\hat{H}_{i}^{\mathcal{S}}over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT and H^l𝒮superscriptsubscript^𝐻𝑙𝒮\hat{H}_{l}^{\mathcal{S}}over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT belonging to the same subsystem and put zero elsewhere.

5.

Compute Eρ𝒜N,m⁢(ϑ)=λmax⁢[(𝒬ρ^𝒜N,m(L))′−4⁢(ΓΠ⁢(ρ^)𝒜N,m(L))′]superscriptsubscript𝐸𝜌subscript𝒜𝑁𝑚bold-italic-ϑsubscript𝜆maxdelimited-[]superscriptsuperscriptsubscript𝒬^𝜌superscriptsubscript𝒜𝑁𝑚L′4superscriptsuperscriptsubscriptΓΠ^𝜌superscriptsubscript𝒜𝑁𝑚L′E_{\rho}^{\mathcal{A}_{\scriptscriptstyle N,m}}(\bm{\vartheta})=\lambda_{\rm% \scriptscriptstyle max}\left[\left(\mathcal{Q}_{\hat{\rho}}^{\mathcal{A}_{% \scriptscriptstyle N,m}^{\rm\scriptscriptstyle(L)}}\right)^{\prime}-4\,\left(% \Gamma_{\Pi(\hat{\rho})}^{\mathcal{A}_{\scriptscriptstyle N,m}^{\rm% \scriptscriptstyle(L)}}\right)^{\prime}\right]italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( bold_italic_ϑ ) = italic_λ start_POSTSUBSCRIPT roman_max end_POSTSUBSCRIPT [ ( caligraphic_Q start_POSTSUBSCRIPT over^ start_ARG italic_ρ end_ARG end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT - 4 ( roman_Γ start_POSTSUBSCRIPT roman_Π ( over^ start_ARG italic_ρ end_ARG ) end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT ], i.e., the entanglement criterion in the new basis. Then, if

𝒲Q⁢(ρ^,𝒜N,m)≡minϑ⁢Eρ𝒜N,m⁢(ϑ)>0,subscript𝒲𝑄^𝜌subscript𝒜𝑁𝑚bold-italic-ϑminsuperscriptsubscript𝐸𝜌subscript𝒜𝑁𝑚bold-italic-ϑ0\mathcal{W}_{Q}(\hat{\rho},\mathcal{A}_{N,m})\equiv\underset{\bm{\vartheta}}{% \text{min}}\,E_{\rho}^{\mathcal{A}_{\scriptscriptstyle N,m}}(\bm{\vartheta})>0,caligraphic_W start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG , caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT ) ≡ underbold_italic_ϑ start_ARG min end_ARG italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , italic_m end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( bold_italic_ϑ ) > 0 , (18)

the state will be not passively separable.

To apply the above procedure we now need to derive the matrix transformation U⁢(ϑ)𝑈bold-italic-ϑU(\bm{\vartheta})italic_U ( bold_italic_ϑ ). We recall that a generator H^i𝒮superscriptsubscript^𝐻𝑖𝒮\hat{H}_{i}^{\mathcal{S}}over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT of type (14) is described by the combination (k1(i),…,k2⁢m(i))superscriptsubscript𝑘1𝑖…superscriptsubscript𝑘2𝑚𝑖(k_{1}^{(i)},\ldots,k_{2m}^{(i)})( italic_k start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT , … , italic_k start_POSTSUBSCRIPT 2 italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT ), where each kj(i)superscriptsubscript𝑘𝑗𝑖k_{j}^{(i)}italic_k start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT corresponds to the power of each quadrature operator involved in the decomposition (14). The operators in the new basis are obtained from the operators in the old basis as

(H^i𝒮)′=∑lUi⁢l⁢(ϑ)⁢H^l𝒮,∀i∈{1,…,ℓ⁢(N,m)}formulae-sequencesuperscriptsuperscriptsubscript^𝐻𝑖𝒮′subscript𝑙subscript𝑈𝑖𝑙bold-italic-ϑsuperscriptsubscript^𝐻𝑙𝒮for-all𝑖1…ℓ𝑁𝑚\left(\hat{H}_{i}^{\mathcal{S}}\right)^{\prime}=\sum_{l}U_{il}(\bm{\vartheta})% \,\hat{H}_{l}^{\mathcal{S}},\quad\forall i\in\{1,\ldots,\ell(N,m)\}( over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT ) start_POSTSUPERSCRIPT ′ end_POSTSUPERSCRIPT = ∑ start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT italic_U start_POSTSUBSCRIPT italic_i italic_l end_POSTSUBSCRIPT ( bold_italic_ϑ ) over^ start_ARG italic_H end_ARG start_POSTSUBSCRIPT italic_l end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_S end_POSTSUPERSCRIPT , ∀ italic_i ∈ { 1 , … , roman_ℓ ( italic_N , italic_m ) } (19)

where

Ui⁢l⁢(ϑ)=1∏j=12⁢mkj(l)!⁢∑σ∈SN∏s=1NOαs⁢βσ⁢(s)T⁢(ϑ),subscript𝑈𝑖𝑙bold-italic-ϑ1superscriptsubscriptproduct𝑗12𝑚superscriptsubscript𝑘𝑗𝑙subscript𝜎subscript𝑆𝑁superscriptsubscriptproduct𝑠1𝑁superscriptsubscript𝑂subscript𝛼𝑠subscript𝛽𝜎𝑠Tbold-italic-ϑU_{il}(\bm{\vartheta})=\frac{1}{\prod_{j=1}^{2m}k_{j}^{(l)}!}\sum_{\sigma\in S% _{N}}\prod_{s=1}^{N}O_{\alpha_{s}\,\beta_{\sigma(s)}}^{\rm\scriptscriptstyle T% }(\bm{\vartheta}),italic_U start_POSTSUBSCRIPT italic_i italic_l end_POSTSUBSCRIPT ( bold_italic_ϑ ) = divide start_ARG 1 end_ARG start_ARG ∏ start_POSTSUBSCRIPT italic_j = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 italic_m end_POSTSUPERSCRIPT italic_k start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_l ) end_POSTSUPERSCRIPT ! end_ARG ∑ start_POSTSUBSCRIPT italic_σ ∈ italic_S start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT end_POSTSUBSCRIPT ∏ start_POSTSUBSCRIPT italic_s = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N end_POSTSUPERSCRIPT italic_O start_POSTSUBSCRIPT italic_α start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT italic_β start_POSTSUBSCRIPT italic_σ ( italic_s ) end_POSTSUBSCRIPT end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT ( bold_italic_ϑ ) , (20)

with

α={1,…⁢1⏟k1(l),…,2⁢m,…,2⁢m⏟k2⁢m(l)},β={1,…⁢1⏟k1(i),…,2⁢m,…,2⁢m⏟k2⁢m(i)}.formulae-sequence𝛼subscript⏟1…1superscriptsubscript𝑘1𝑙…subscript⏟2𝑚…2𝑚superscriptsubscript𝑘2𝑚𝑙𝛽subscript⏟1…1superscriptsubscript𝑘1𝑖…subscript⏟2𝑚…2𝑚superscriptsubscript𝑘2𝑚𝑖\begin{split}\alpha&=\{\underbrace{1,\ldots 1}_{k_{1}^{(l)}},\ldots,% \underbrace{2m,\ldots,2m}_{k_{2m}^{(l)}}\},\\ \beta&=\{\underbrace{1,\ldots 1}_{k_{1}^{(i)}},\ldots,\underbrace{2m,\ldots,2m% }_{k_{2m}^{(i)}}\}.\end{split}start_ROW start_CELL italic_α end_CELL start_CELL = { under⏟ start_ARG 1 , … 1 end_ARG start_POSTSUBSCRIPT italic_k start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_l ) end_POSTSUPERSCRIPT end_POSTSUBSCRIPT , … , under⏟ start_ARG 2 italic_m , … , 2 italic_m end_ARG start_POSTSUBSCRIPT italic_k start_POSTSUBSCRIPT 2 italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_l ) end_POSTSUPERSCRIPT end_POSTSUBSCRIPT } , end_CELL end_ROW start_ROW start_CELL italic_β end_CELL start_CELL = { under⏟ start_ARG 1 , … 1 end_ARG start_POSTSUBSCRIPT italic_k start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT end_POSTSUBSCRIPT , … , under⏟ start_ARG 2 italic_m , … , 2 italic_m end_ARG start_POSTSUBSCRIPT italic_k start_POSTSUBSCRIPT 2 italic_m end_POSTSUBSCRIPT start_POSTSUPERSCRIPT ( italic_i ) end_POSTSUPERSCRIPT end_POSTSUBSCRIPT } . end_CELL end_ROW (21)

In expression (20), SNsubscript𝑆𝑁S_{N}italic_S start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT represents the permutation group of degree N𝑁Nitalic_N and O⁢(ϑ)𝑂bold-italic-ϑO(\bm{\vartheta})italic_O ( bold_italic_ϑ ) is the (linear) mode basis change transformation at the level of 𝒜N=1,m={q^1,…,q^m,p^1,…,p^m}subscript𝒜𝑁1𝑚subscript^𝑞1…subscript^𝑞𝑚subscript^𝑝1…subscript^𝑝𝑚\mathcal{A}_{N=1,m}=\{\hat{q}_{1},\ldots,\hat{q}_{m},\hat{p}_{1},\ldots,\hat{p% }_{m}\}caligraphic_A start_POSTSUBSCRIPT italic_N = 1 , italic_m end_POSTSUBSCRIPT = { over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , … , over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_m end_POSTSUBSCRIPT } (note that here there are only local operators). For a given number of modes m𝑚mitalic_m the most general linear transformation O𝑂Oitalic_O that we are concerned with can be represented in a compact way, using, for example, Clements’s decomposition Clements et al. (2016), by the combination of m⁢(m−1)/2𝑚𝑚12m(m-1)/2italic_m ( italic_m - 1 ) / 2 real beamsplitters and m⁢(m−1)/2𝑚𝑚12m(m-1)/2italic_m ( italic_m - 1 ) / 2 local rotations (phase shifts), as represented in Fig.2. The transformation depends on m⁢(m−1)𝑚𝑚1m(m-1)italic_m ( italic_m - 1 ) parameters that we gather in the vector ϑ=(𝜽,𝝋)bold-italic-ϑ𝜽𝝋\bm{\vartheta}=(\bm{\theta},\bm{\varphi})bold_italic_ϑ = ( bold_italic_θ , bold_italic_φ ), where 𝜽𝜽\bm{\theta}bold_italic_θ contains the m⁢(m−1)/2𝑚𝑚12m(m-1)/2italic_m ( italic_m - 1 ) / 2 angles of the beamspliters, and 𝝋𝝋\bm{\varphi}bold_italic_φ contains the m⁢(m−1)/2𝑚𝑚12m(m-1)/2italic_m ( italic_m - 1 ) / 2 phases of the phase shifters. To reproduce the full set of passive transformations and additional set of m𝑚mitalic_m local phase shifters should be included, nevertheless it is irrelevant for our task, given that these do not alter the entanglement property of the states we consider. More details on how to derive Eq. (20) can be found in Appendix [REF].

To illustrate the quite complicated expression (20), let us take an example with m=2𝑚2m=2italic_m = 2 and N=3𝑁3N=3italic_N = 3: consider the generator (2,0,1,0)2010(2,0,1,0)( 2 , 0 , 1 , 0 ) [we decide to order the set of generators such that it corresponds to i=3𝑖3i=3italic_i = 3] which corresponds to the symmetrised version of q^12⁢p^1superscriptsubscript^𝑞12subscript^𝑝1\hat{q}_{1}^{2}\,\hat{p}_{1}over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT. We then apply for simplicity the following transformation

O≡R⁢(θ)=(cos⁡θsin⁡θ00−sin⁡θcos⁡θ0000cos⁡θsin⁡θ00−sin⁡θcos⁡θ).𝑂𝑅𝜃matrix𝜃𝜃00𝜃𝜃0000𝜃𝜃00𝜃𝜃O\equiv R(\theta)=\begin{pmatrix}\cos\theta&\sin\theta&0&0\\ -\sin\theta&\cos\theta&0&0\\ 0&0&\cos\theta&\sin\theta\\ 0&0&-\sin\theta&\cos\theta\\ \end{pmatrix}.italic_O ≡ italic_R ( italic_θ ) = ( start_ARG start_ROW start_CELL roman_cos italic_θ end_CELL start_CELL roman_sin italic_θ end_CELL start_CELL 0 end_CELL start_CELL 0 end_CELL end_ROW start_ROW start_CELL - roman_sin italic_θ end_CELL start_CELL roman_cos italic_θ end_CELL start_CELL 0 end_CELL start_CELL 0 end_CELL end_ROW start_ROW start_CELL 0 end_CELL start_CELL 0 end_CELL start_CELL roman_cos italic_θ end_CELL start_CELL roman_sin italic_θ end_CELL end_ROW start_ROW start_CELL 0 end_CELL start_CELL 0 end_CELL start_CELL - roman_sin italic_θ end_CELL start_CELL roman_cos italic_θ end_CELL end_ROW end_ARG ) . (22)

Let us derive the term U31⁢(θ)subscript𝑈31𝜃U_{31}(\theta)italic_U start_POSTSUBSCRIPT 31 end_POSTSUBSCRIPT ( italic_θ ), where l=1𝑙1l=1italic_l = 1 corresponds to the generator (1,1,1,0)1110(1,1,1,0)( 1 , 1 , 1 , 0 ), i.e., the symmetrised version of q^1⁢p^1⁢q^2subscript^𝑞1subscript^𝑝1subscript^𝑞2\hat{q}_{1}\,\hat{p}_{1}\,\hat{q}_{2}over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT. In this case the vectors defined in Eq. (21) are α=(1,2,3)𝛼123\alpha=(1,2,3)italic_α = ( 1 , 2 , 3 ) and β=(1,1,3)𝛽113\beta=(1,1,3)italic_β = ( 1 , 1 , 3 ), such that it leads to

U31⁢(θ)=21!⁢O11T⁢O21T⁢O33T=2⁢cos2⁡θ⁢sin⁡θ,subscript𝑈31𝜃21superscriptsubscript𝑂11Tsuperscriptsubscript𝑂21Tsuperscriptsubscript𝑂33T2superscript2𝜃𝜃U_{31}(\theta)=\frac{2}{1!}O_{11}^{\rm\scriptscriptstyle T}\,O_{21}^{\rm% \scriptscriptstyle T}\,O_{33}^{\rm\scriptscriptstyle T}=2\cos^{2}\theta\,\sin\theta,italic_U start_POSTSUBSCRIPT 31 end_POSTSUBSCRIPT ( italic_θ ) = divide start_ARG 2 end_ARG start_ARG 1 ! end_ARG italic_O start_POSTSUBSCRIPT 11 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT italic_O start_POSTSUBSCRIPT 21 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT italic_O start_POSTSUBSCRIPT 33 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT roman_T end_POSTSUPERSCRIPT = 2 roman_cos start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT italic_θ roman_sin italic_θ ,

since all other terms are zero (it comes with a factor 2 since we need to sum over all permutations). By computing (cos⁡θ⁢q^1+sin⁡θ⁢q^2)2⁢(cos⁡θ⁢p^1+sin⁡θ⁢p^2)superscript𝜃subscript^𝑞1𝜃subscript^𝑞22𝜃subscript^𝑝1𝜃subscript^𝑝2(\cos\theta\hat{q}_{1}+\sin\theta\hat{q}_{2})^{2}\,(\cos\theta\hat{p}_{1}+\sin% \theta\hat{p}_{2})( roman_cos italic_θ over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + roman_sin italic_θ over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT ( roman_cos italic_θ over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + roman_sin italic_θ over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ), one can indeed check that the coefficient in front of q^1⁢p^1⁢q^2subscript^𝑞1subscript^𝑝1subscript^𝑞2\hat{q}_{1}\,\hat{p}_{1}\,\hat{q}_{2}over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT is indeed given by 2⁢cos2⁡θ⁢sin⁡θ2superscript2𝜃𝜃2\cos^{2}\theta\,\sin\theta2 roman_cos start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT italic_θ roman_sin italic_θ.

Expression (20) is fully general for any transformation O𝑂Oitalic_O and any number of modes m𝑚mitalic_m and any set of generators 𝒜Nsubscript𝒜𝑁\mathcal{A}_{\scriptscriptstyle N}caligraphic_A start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT. However, the computation time increase rapidly when considering high order generators given the factorial scaling with N𝑁Nitalic_N.

IV Application of the witness to experimentally relevant states

To explore the capability of the witness derived in the previous section to detect mode intrinsic entanglement for bipartite and multimode systems, we use it on a specific class of non-Gaussian states: photon subtracted states. These states are produced by subtracting one or more photons from a Gaussian state. The removal of photons is a non-Gaussian operation, thus rendering the state non-Gaussian. Therefore, it makes it possible to explore how quantum properties are affected by the non-Gaussianity using a class of states which falls within a well-established theoretical framework Walschaers (2021); Walschaers et al. (2019) and has the advantage to be already produced in the lab with state-of-the-art experiments using different spatial and temporal modes of indistinguishable sources Neergaard-Nielsen et al. (2006); Wakui et al. (2007); Endo et al. (2023), or in temporal-frequency domain Ra et al. (2020).

The general expression of the resulting states is given by

|ψ⟩∝(∏i=1NSa^i)⁢O^⁢(∏j=1MS^⁢(rj,ϕ¯j))⁢|𝟎⟩,proportional-toket𝜓superscriptsubscriptproduct𝑖1subscript𝑁𝑆subscript^𝑎𝑖^𝑂superscriptsubscriptproduct𝑗1𝑀^𝑆subscript𝑟𝑗subscript¯italic-ϕ𝑗ket0|\psi\rangle\propto\left(\prod_{i=1}^{N_{S}}\hat{a}_{i}\right)\hat{O}\left(% \prod_{j=1}^{M}\hat{S}(r_{j},\bar{\phi}_{j})\right)|\mathbf{0}\rangle,| italic_ψ ⟩ ∝ ( ∏ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT end_POSTSUPERSCRIPT over^ start_ARG italic_a end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ) over^ start_ARG italic_O end_ARG ( ∏ start_POSTSUBSCRIPT italic_j = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_M end_POSTSUPERSCRIPT over^ start_ARG italic_S end_ARG ( italic_r start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT , over¯ start_ARG italic_ϕ end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ) ) | bold_0 ⟩ , (23)

up to a normalization factor, where NSsubscript𝑁𝑆N_{S}italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT is the number of photons subtracted and M𝑀Mitalic_M is the number of modes in the cluster. S^⁢(rj,ϕ¯j)^𝑆subscript𝑟𝑗subscript¯italic-ϕ𝑗\hat{S}(r_{j},\bar{\phi}_{j})over^ start_ARG italic_S end_ARG ( italic_r start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT , over¯ start_ARG italic_ϕ end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ) is a single mode squeezing operation that squeezes the quadrature ξ^ϕ¯j=q^j⁢cos⁡ϕ¯j+p^j⁢sin⁡ϕ¯jsubscript^𝜉subscript¯italic-ϕ𝑗subscript^𝑞𝑗subscript¯italic-ϕ𝑗subscript^𝑝𝑗subscript¯italic-ϕ𝑗\hat{\xi}_{\bar{\phi}_{j}}=\hat{q}_{j}\cos\bar{\phi}_{j}+\hat{p}_{j}\sin\bar{% \phi}_{j}over^ start_ARG italic_ξ end_ARG start_POSTSUBSCRIPT over¯ start_ARG italic_ϕ end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT end_POSTSUBSCRIPT = over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT roman_cos over¯ start_ARG italic_ϕ end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT + over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT roman_sin over¯ start_ARG italic_ϕ end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT, of mode j𝑗jitalic_j, with a squeezing level rjsubscript𝑟𝑗r_{j}italic_r start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT. O^^𝑂\hat{O}over^ start_ARG italic_O end_ARG is a general passive transformation, and a^isubscript^𝑎𝑖\hat{a}_{i}over^ start_ARG italic_a end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT are the annihilation operator corresponding to mode i𝑖iitalic_i. In Fig. 3 we illustrate the case Ns=2subscript𝑁𝑠2N_{s}=2italic_N start_POSTSUBSCRIPT italic_s end_POSTSUBSCRIPT = 2. We also point out that photon subtractions can happen in a superposition of different modes of the cluster state; this case will only be considered for two-mode states, as in any case it can be casted into the form above by including the mixing operation in the unitary O^^𝑂\hat{O}over^ start_ARG italic_O end_ARG (see Sec. IV.2).

We split this section in two subsections: we start with the case of multimode cluster states to show that the witness can be applied and successfully detect mode-intrinsic entanglement. In particular, we apply the metrological witness to three-mode, four-mode and five-mode cluster states. Then we move to two-mode photon subtracted states. In some configurations discussed below this leads to mode-intrinsic entangled states that the witness can detect. For such a two-mode system the unitary transformation U⁢(ϑ)𝑈bold-italic-ϑU(\bm{\vartheta})italic_U ( bold_italic_ϑ ) to go from one basis to another (as represented in Fig. 2) only depends on two parameters ϑ=(θ,φ)bold-italic-ϑ𝜃𝜑\bm{\vartheta}=(\theta,\varphi)bold_italic_ϑ = ( italic_θ , italic_φ ). Therefore, one can easily visualize the metrological witness Eρ𝒜N,2⁢(θ,φ)superscriptsubscript𝐸𝜌subscript𝒜𝑁2𝜃𝜑E_{\rho}^{\mathcal{A}_{\scriptscriptstyle N,2}}(\theta,\varphi)italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N , 2 end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( italic_θ , italic_φ ) defined in Eq. (12) in all possible bases in a two-dimensional parameter space. Likewise, as we discuss in Sec. V, the mode-intrinsic entanglement witness for two-mode photon subtracted states could be measured experimentally using the Fisher information matrix, making it a very interesting and promising study case towards the experimental detection of non-Gaussian entanglement.

IV.1 Multimode CV-cluster states

A CV-cluster state is made of a set of modes which are coupled together by applying entangling CZ=exp⁡(i⁢q^i⊗q^j)subscript𝐶𝑍tensor-product𝑖subscript^𝑞𝑖subscript^𝑞𝑗C_{Z}=\exp(i\,\hat{q}_{i}\otimes\hat{q}_{j})italic_C start_POSTSUBSCRIPT italic_Z end_POSTSUBSCRIPT = roman_exp ( italic_i over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ⊗ over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_j end_POSTSUBSCRIPT ) gates. Overall, they form a graph structure (𝒱,ℰ)𝒱ℰ(\mathcal{V},\mathcal{E})( caligraphic_V , caligraphic_E ), where 𝒱𝒱\mathcal{V}caligraphic_V is the set of nodes, i.e., the set of modes in the multimode state, while ℰℰ\mathcal{E}caligraphic_E represents the set of edges between the nodes, i.e., the quantum correlations that connect the modes together after applying the CZsubscript𝐶𝑍C_{Z}italic_C start_POSTSUBSCRIPT italic_Z end_POSTSUBSCRIPT gates. The graph formed by the cluster state can be described by the so-called adjacency matrix V𝑉Vitalic_V whose coefficients Vi⁢jsubscript𝑉𝑖𝑗V_{ij}italic_V start_POSTSUBSCRIPT italic_i italic_j end_POSTSUBSCRIPT are 1 if an edge exists between nodes i𝑖iitalic_i and j𝑗jitalic_j (in other words if modes i𝑖iitalic_i and j𝑗jitalic_j are entangled to each other), 0 otherwise. We show in Fig. 4 the three cluster states that we will consider in this section. As an example we give the adjacency matrix of the three-mode cluster state:

V3−modes=(010101010).subscript𝑉3modesmatrix010101010V_{\rm\scriptscriptstyle 3-modes}=\begin{pmatrix}0&1&0\\ 1&0&1\\ 0&1&0\end{pmatrix}.italic_V start_POSTSUBSCRIPT 3 - roman_modes end_POSTSUBSCRIPT = ( start_ARG start_ROW start_CELL 0 end_CELL start_CELL 1 end_CELL start_CELL 0 end_CELL end_ROW start_ROW start_CELL 1 end_CELL start_CELL 0 end_CELL start_CELL 1 end_CELL end_ROW start_ROW start_CELL 0 end_CELL start_CELL 1 end_CELL start_CELL 0 end_CELL end_ROW end_ARG ) .

The red nodes in Fig. 4 indicate in which mode we subtract a photon to create our photon subtracted cluster state.

A cluster state associated with a given adjacency matrix V𝑉Vitalic_V is obtained from a set of single-mode squeezed states by applying a unitary matrix

UV=(𝟙+i⁢V)⁢(V2+𝟙)−1/2⁢𝒪,subscript𝑈𝑉1𝑖𝑉superscriptsuperscript𝑉2112𝒪U_{\scriptscriptstyle V}=(\mathds{1}+i\,V)\,(V^{2}+\mathds{1})^{-1/2}\,% \mathcal{O},italic_U start_POSTSUBSCRIPT italic_V end_POSTSUBSCRIPT = ( blackboard_1 + italic_i italic_V ) ( italic_V start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT + blackboard_1 ) start_POSTSUPERSCRIPT - 1 / 2 end_POSTSUPERSCRIPT caligraphic_O , (24)

where 𝟙1\mathds{1}blackboard_1 is the unitary matrix and 𝒪𝒪\mathcal{O}caligraphic_O can be any real orthogonal matrix van Loock et al. (2007); Ferrini et al. (2015); Sansavini and Parigi (2020). The choice of 𝒪𝒪\mathcal{O}caligraphic_O can be optimized to target a specific cluster configuration as discussed in Ref. Ferrini et al. (2015). Indeed, a perfect cluster state where correlations between modes are exactly given by the adjacency matrix can only be obtained in the limit of infinite squeezing. In this case, the variances of certain combinations of quadratures, called nullifiers, are exactly zero; these nullifiers are defined as

δ^i≡p^i−∑k=1mVi⁢k⁢q^k,∀i∈{1,…,m},formulae-sequencesubscript^𝛿𝑖subscript^𝑝𝑖superscriptsubscript𝑘1𝑚subscript𝑉𝑖𝑘subscript^𝑞𝑘for-all𝑖1…𝑚\hat{\delta}_{i}\equiv\hat{p}_{i}-\sum_{k=1}^{m}V_{ik}\hat{q}_{k},\quad\forall i% \in\{1,\ldots,m\},over^ start_ARG italic_δ end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT ≡ over^ start_ARG italic_p end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT - ∑ start_POSTSUBSCRIPT italic_k = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_m end_POSTSUPERSCRIPT italic_V start_POSTSUBSCRIPT italic_i italic_k end_POSTSUBSCRIPT over^ start_ARG italic_q end_ARG start_POSTSUBSCRIPT italic_k end_POSTSUBSCRIPT , ∀ italic_i ∈ { 1 , … , italic_m } , (25)

where m𝑚mitalic_m is the number of modes in the cluster. In the case of finite squeezing, the variances of the nullifiers are not zero but one can minimize them to tend to the limit of the perfect cluster state by optimizing the orthogonal matrix 𝒪𝒪\mathcal{O}caligraphic_O in Eq. (24). We thus follow the procedure described in Ref. Ferrini et al. (2015) using a genetic algorithm to minimize the variances of the nullifiers and obtain the best cluster states given the three graph structures of Fig. 4. In order to consider realistic squeezing parameters we take the ones measured in the 4-mode cluster states of Ref. Medeiros de Araújo et al. (2014): −77-7- 7 dB, −66-6- 6 dB, −44-4- 4 dB, 0 dB. For the other cluster states considered here we use similar squeezing parameters: −77-7- 7 dB, −44-4- 4 dB, 0 dB (3-mode cluster state); −77-7- 7 dB, −66-6- 6 dB, −44-4- 4 dB, −22-2- 2 dB, 0 dB (5-mode cluster state).

Before we subtract a photon, the cluster states are fully Gaussian. Therefore, we do not expect to witness any mode-intrinsic entanglement since there will always exist a basis where the states are separable. The subtraction of a photon in the first mode (the red nodes in Fig. 4) renders the state non-Gaussian and will affect the modes in its vicinity; we thus anticipate to detect mode-intrinsic entanglement. In the case of a m𝑚mitalic_m-multimode system we have the freedom to choose how we split it into p𝑝pitalic_p different partitions 𝒫1⁢|…|⁢𝒫psubscript𝒫1…subscript𝒫𝑝\mathcal{P}_{1}|\ldots|\mathcal{P}_{p}caligraphic_P start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT | … | caligraphic_P start_POSTSUBSCRIPT italic_p end_POSTSUBSCRIPT and eventually detect p𝑝pitalic_p-partite mode-intrinsic entanglement. It is important to note that, in our case, the labelling of the modes is irrelevant; only the number of modes in each partition matters. Indeed, for any given mode splitting, we can always find a linear passive operation that changes the partitioning while keeping the same number of modes in each partition 𝒫isubscript𝒫𝑖\mathcal{P}_{i}caligraphic_P start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT (but with different mode labels). Therefore, as the witness covers all possible mode bases, any two mode splittings differing only by a reordering of the mode labels are completely equivalent regarding mode-intrinsic entanglement.

Tables 1, 2, 3 gather the values of the mode-intrinsic entanglement witness 𝒲Q⁢(ρ^,𝒜N=1,m)subscript𝒲𝑄^𝜌subscript𝒜𝑁1𝑚\mathcal{W}_{Q}(\hat{\rho},\mathcal{A}_{N=1,m})caligraphic_W start_POSTSUBSCRIPT italic_Q end_POSTSUBSCRIPT ( over^ start_ARG italic_ρ end_ARG , caligraphic_A start_POSTSUBSCRIPT italic_N = 1 , italic_m end_POSTSUBSCRIPT ) – using first-order generators only – obtained for the three clusters and for all possible mode splittings (up to a reordering of the mode labels as discussed above). Given the number of parameters that we need to optimize to find the minimum of the entanglement witness Eρ𝒜N⁢(ϑ)superscriptsubscript𝐸𝜌subscript𝒜𝑁bold-italic-ϑE_{\rho}^{\mathcal{A}_{\scriptscriptstyle N}}(\bm{\vartheta})italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( bold_italic_ϑ ) (see Fig. 2 where each blue and red block of the unitary transformation adds one degree of freedom), we also use here a genetic algorithm to solve this optimization problem. As we can clearly see from the results, for many mode splittings, the subtraction of the photon renders the cluster states mode-intrinsic entangled. On the other hand, it is very instructive to notice that only in the case where we split the cluster states into a bipartition where one of the subsystems contain m−1𝑚1m-1italic_m - 1 modes, the value of the witness becomes zero. In this case, while we cannot guarantee that the cluster is not entangled in at least one basis since we only use a witness of entanglement, the non-Gaussian correlations can be probably localized in the subsytem containing m−1𝑚1m-1italic_m - 1 modes leaving only Gaussian correlations between both subsystems.

While it is difficult to draw any systematic conclusion from the results above given that we only tested the witness in three specific cluster configurations, it shows that the mode-intrinsic entanglement witness with first-order generators is able to detect mode-intrinsic entanglement for experimentally relevant cluster states. It also opens a new toolbox to explore and better understand the entanglement structure of large multimode non-Gaussian states such as those considered in Ref. Walschaers et al. (2023). A deeper exploration of this tool for the study of large clusters is left for future works. In the rest of this paper we rather focus on the study of the smallest version of these systems, i.e., two mode, two photon subtracted states. In this more controlled setting we can fully explore the capabilities and limitations of the method we develop, while also restricting to a scenario with more immediate potential to be tested experimentally.

IV.2 Two-mode photon subtracted states

As explained above, we now turn our attention to two-mode photon subtracted states defined by

|ψ⟩=∏i=1NS(cos⁡Θi⁢a^1+sin⁡Θi⁢a^2)⁢S^⁢(r1)⁢S^⁢(r2)⁢|00⟩,ket𝜓superscriptsubscriptproduct𝑖1subscript𝑁𝑆subscriptΘ𝑖subscript^𝑎1subscriptΘ𝑖subscript^𝑎2^𝑆subscript𝑟1^𝑆subscript𝑟2ket00|\psi\rangle=\prod_{i=1}^{N_{\scriptscriptstyle S}}(\cos\Theta_{i}\,\hat{a}_{1% }+\sin\Theta_{i}\,\hat{a}_{2})\,\hat{S}(r_{1})\,\hat{S}(r_{2})\,|00\rangle,| italic_ψ ⟩ = ∏ start_POSTSUBSCRIPT italic_i = 1 end_POSTSUBSCRIPT start_POSTSUPERSCRIPT italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( roman_cos roman_Θ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT over^ start_ARG italic_a end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT + roman_sin roman_Θ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT over^ start_ARG italic_a end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) over^ start_ARG italic_S end_ARG ( italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ) over^ start_ARG italic_S end_ARG ( italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) | 00 ⟩ , (26)

where S^⁢(r1)^𝑆subscript𝑟1\hat{S}(r_{1})over^ start_ARG italic_S end_ARG ( italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ) and S^⁢(r2)^𝑆subscript𝑟2\hat{S}(r_{2})over^ start_ARG italic_S end_ARG ( italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ) are single-mode squeezing operators, equivalent to those in equation (23), where we removed the phase dependence ϕ¯isubscript¯italic-ϕ𝑖\bar{\phi}_{i}over¯ start_ARG italic_ϕ end_ARG start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT and allow risubscript𝑟𝑖r_{i}italic_r start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT to be negative. Angles ΘisubscriptΘ𝑖\Theta_{i}roman_Θ start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT control mode superposition in which photons are subtracted, which effectively modifies the probability that the subtracted photon came from either mode for each subtraction j=1,…,NS𝑗1…subscript𝑁𝑆j=1,\ldots,N_{\scriptscriptstyle S}italic_j = 1 , … , italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT. In general we could have considered the subtraction to happen in an already entangled Gaussian state, nevertheless, those can always be brought to the form in equation (26) by the action of some passive unitary. In this subsection we consider two cases: NS=1subscript𝑁𝑆1N_{\scriptscriptstyle S}=1italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 1 and NS=2subscript𝑁𝑆2N_{\scriptscriptstyle S}=2italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 2. In both cases we consider the single mode squeezed states to be either squeezed in position (ri>0subscript𝑟𝑖0r_{i}>0italic_r start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT > 0) or in momentum (ri<0subscript𝑟𝑖0r_{i}<0italic_r start_POSTSUBSCRIPT italic_i end_POSTSUBSCRIPT < 0). One-photon subtracted states (NS=1subscript𝑁𝑆1N_{\scriptscriptstyle S}=1italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 1) are not mode-intrinsic entangled in the two following cases: (i) r1=r2subscript𝑟1subscript𝑟2r_{1}=r_{2}italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT and ∀Θ1for-allsubscriptΘ1\forall\,\Theta_{1}∀ roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT, (ii) Θ1=n⁢π/2subscriptΘ1𝑛𝜋2\Theta_{1}=n\pi/2roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = italic_n italic_π / 2, n∈ℤ𝑛ℤn\in\mathbb{Z}italic_n ∈ blackboard_Z and ∀r1,r2for-allsubscript𝑟1subscript𝑟2\forall r_{1},r_{2}∀ italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT; and, similarly, two-photon subtracted (NS=2subscript𝑁𝑆2N_{\scriptscriptstyle S}=2italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 2) are not mode-intrinsic entangled for: (i) r1=r2subscript𝑟1subscript𝑟2r_{1}=r_{2}italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT and ∀Θ1for-allsubscriptΘ1\forall\,\Theta_{1}∀ roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT, (ii) Θ1=n⁢π/2subscriptΘ1𝑛𝜋2\Theta_{1}=n\pi/2roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = italic_n italic_π / 2, n∈ℤ𝑛ℤn\in\mathbb{Z}italic_n ∈ blackboard_Z and ∀r1,r2for-allsubscript𝑟1subscript𝑟2\forall r_{1},r_{2}∀ italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT (NS=1subscript𝑁𝑆1N_{\scriptscriptstyle S}=1italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 1); (i) r1=r2subscript𝑟1subscript𝑟2r_{1}=r_{2}italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT, Θ1subscriptΘ1\Theta_{1}roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT and Θ2subscriptΘ2\Theta_{2}roman_Θ start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT orthogonal to each other, (ii) Θ1=n⁢π/2subscriptΘ1𝑛𝜋2\Theta_{1}=n\pi/2roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = italic_n italic_π / 2 and Θ2=m⁢π/2subscriptΘ2𝑚𝜋2\Theta_{2}=m\pi/2roman_Θ start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT = italic_m italic_π / 2 n,m∈ℤ𝑛𝑚ℤn,m\in\mathbb{Z}italic_n , italic_m ∈ blackboard_Z and ∀r1,r2for-allsubscript𝑟1subscript𝑟2\forall r_{1},r_{2}∀ italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT , italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT. In all other cases one can prove that these states are entangled in all mode bases. To illustrate how the metrological witness (18) can detect mode intrinsic entanglement we now consider two specific examples: NS=1subscript𝑁𝑆1N_{\scriptscriptstyle S}=1italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 1, r1=−r2=0.2subscript𝑟1subscript𝑟20.2r_{1}=-r_{2}=0.2italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = - italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT = 0.2, Θ1=π/4subscriptΘ1𝜋4\Theta_{1}=\pi/4roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = italic_π / 4; NS=2subscript𝑁𝑆2N_{\scriptscriptstyle S}=2italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 2, r1=−r2=0.2subscript𝑟1subscript𝑟20.2r_{1}=-r_{2}=0.2italic_r start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = - italic_r start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT = 0.2, Θ1=−Θ2=π/4subscriptΘ1subscriptΘ2𝜋4\Theta_{1}=-\Theta_{2}=\pi/4roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT = - roman_Θ start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT = italic_π / 4.

In the case of a two-mode state, the mode basis change transformation at the level of quadrature operators using Clement’s decomposition can be written as

O⁢(θ,φ)=R⁢(θ)⁢P⁢h⁢(φ),𝑂𝜃𝜑𝑅𝜃𝑃ℎ𝜑O(\theta,\varphi)=R(\theta)\,Ph(\varphi),italic_O ( italic_θ , italic_φ ) = italic_R ( italic_θ ) italic_P italic_h ( italic_φ ) , (27)

where R⁢(θ)𝑅𝜃R(\theta)italic_R ( italic_θ ) is given in Eq. (22) and

P⁢h⁢(φ)=(cos⁡φ0sin⁡φ00100−sin⁡φ0cos⁡φ00001).𝑃ℎ𝜑matrix𝜑0𝜑00100𝜑0𝜑00001Ph(\varphi)=\begin{pmatrix}\cos\varphi&0&\sin\varphi&0\\ 0&1&0&0\\ -\sin\varphi&0&\cos\varphi&0\\ 0&0&0&1\\ \end{pmatrix}.italic_P italic_h ( italic_φ ) = ( start_ARG start_ROW start_CELL roman_cos italic_φ end_CELL start_CELL 0 end_CELL start_CELL roman_sin italic_φ end_CELL start_CELL 0 end_CELL end_ROW start_ROW start_CELL 0 end_CELL start_CELL 1 end_CELL start_CELL 0 end_CELL start_CELL 0 end_CELL end_ROW start_ROW start_CELL - roman_sin italic_φ end_CELL start_CELL 0 end_CELL start_CELL roman_cos italic_φ end_CELL start_CELL 0 end_CELL end_ROW start_ROW start_CELL 0 end_CELL start_CELL 0 end_CELL start_CELL 0 end_CELL start_CELL 1 end_CELL end_ROW end_ARG ) .

We can then compute the mode basis transformation U⁢(θ,φ)𝑈𝜃𝜑U(\theta,\varphi)italic_U ( italic_θ , italic_φ ) using Eq. (20) and apply the metrological witness Eρ𝒜N⁢(θ,φ)superscriptsubscript𝐸𝜌subscript𝒜𝑁𝜃𝜑E_{\rho}^{\mathcal{A}_{\scriptscriptstyle N}}(\theta,\varphi)italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( italic_θ , italic_φ ) [see Eq. (18)] using first order and second order generators 𝒜N=1(L)subscriptsuperscript𝒜L𝑁1\mathcal{A}^{\scriptscriptstyle\rm(L)}_{N=1}caligraphic_A start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_N = 1 end_POSTSUBSCRIPT and 𝒜N=2(L)subscriptsuperscript𝒜L𝑁2\mathcal{A}^{\scriptscriptstyle\rm(L)}_{N=2}caligraphic_A start_POSTSUPERSCRIPT ( roman_L ) end_POSTSUPERSCRIPT start_POSTSUBSCRIPT italic_N = 2 end_POSTSUBSCRIPT, respectively. Notice that we have omitted the subscript m=2𝑚2m=2italic_m = 2 from the definition of the set of generators, as from this point on we will only be dealing with two-mode states.

The results are shown in Fig. 5 for NS=1subscript𝑁𝑆1N_{\scriptscriptstyle S}=1italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 1 and Fig. 6 for NS=2subscript𝑁𝑆2N_{\scriptscriptstyle S}=2italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 2. In both figures, the top graphs labeled as (a) correspond to the results using generators in the set 𝒜N=1subscript𝒜𝑁1\mathcal{A}_{N=1}caligraphic_A start_POSTSUBSCRIPT italic_N = 1 end_POSTSUBSCRIPT, while we used the set 𝒜N=2subscript𝒜𝑁2\mathcal{A}_{N=2}caligraphic_A start_POSTSUBSCRIPT italic_N = 2 end_POSTSUBSCRIPT to produce the bottom graphs labeled as (b). The black dots indicate the minimal value of Eρ𝒜N⁢(θ,φ)superscriptsubscript𝐸𝜌subscript𝒜𝑁𝜃𝜑E_{\rho}^{\mathcal{A}_{\scriptscriptstyle N}}(\theta,\varphi)italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( italic_θ , italic_φ ) when spanning the parameter space (θ,φ)𝜃𝜑(\theta,\varphi)( italic_θ , italic_φ ). One sees that for the one-photon subtracted state the metrological witness is able to detect mode intrinsic entanglement both at first order and second order on the generators: the minimal value over the parameter space is always strictly greater than zero, thus demonstrating entanglement in all mode bases. The robustness of the witness against optical losses is adressed below (see in particular Figs. 7 and 8). In the case of the two-photon subtracted state considered here the situation is different; as shown in Fig. 6 (a), first order fails to detect entanglement in all bases: the minimum of Eρ𝒜N=1⁢(θ,φ)superscriptsubscript𝐸𝜌subscript𝒜𝑁1𝜃𝜑E_{\rho}^{\mathcal{A}_{\scriptscriptstyle N=1}}(\theta,\varphi)italic_E start_POSTSUBSCRIPT italic_ρ end_POSTSUBSCRIPT start_POSTSUPERSCRIPT caligraphic_A start_POSTSUBSCRIPT italic_N = 1 end_POSTSUBSCRIPT end_POSTSUPERSCRIPT ( italic_θ , italic_φ ) reaches zero for a range of values θ=0𝜃0\theta=0italic_θ = 0 (∀ϕ)for-allitalic-ϕ(\forall\phi)( ∀ italic_ϕ ), θ=π/2𝜃𝜋2\theta=\pi/2italic_θ = italic_π / 2 (∀ϕ)for-allitalic-ϕ(\forall\phi)( ∀ italic_ϕ ) and ϕ=π/2italic-ϕ𝜋2\phi=\pi/2italic_ϕ = italic_π / 2 (∀θ)for-all𝜃(\forall\theta)( ∀ italic_θ ) (see the black lines at the borders of the figure). Therefore, the first-order protocol cannot be used to assess the mode intrinsic entanglement of the state. On the contrary, as shown in Fig. 6 (b), once we use second order generators, the minimal value when spanning all mode-bases become strictly positive (it equals 0.98 in this precise example), thus meaning that the state is mode-intrinsic entangled. This illustrates the interest to use higher generators and to not restrict the metrological witness to first order generators. More generally one can show that for two photon subtracted states (NS=2subscript𝑁𝑆2N_{\scriptscriptstyle S}=2italic_N start_POSTSUBSCRIPT italic_S end_POSTSUBSCRIPT = 2) first order always fails when Θ1subscriptΘ1\Theta_{1}roman_Θ start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT and Θ2subscriptΘ2\Theta_{2}roman_Θ start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT are orthogonal to each other.

So far we did not include optical losses inherent to any experimental setup. We can thus reproduce the previous results and apply the metrological witness to lossy one and two-photon subtracted states. To simulate such losses we apply two single-mode loss channels L^1⁢(η)subscript^𝐿1𝜂\hat{L}_{1}(\eta)over^ start_ARG italic_L end_ARG start_POSTSUBSCRIPT 1 end_POSTSUBSCRIPT ( italic_η ) and L^2⁢(η)subscript^𝐿2𝜂\hat{L}_{2}(\eta)over^ start_ARG italic_L end_ARG start_POSTSUBSCRIPT 2 end_POSTSUBSCRIPT ( italic_η ) to each mode as described