{'fraction_non_alphanumeric': 0.0351637928114096, 'fraction_numerical': 0.012026152647447406, 'mean_word_length': 4.5198204936424835, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Multiple copy distillation and purification of phase diffused squeezed states

We provide a detailed theoretical analysis of multiple copy purification and distillation protocols for phase diffused squeezed states of light. The standard iterative distillation protocol is generalized to a collective purification of an arbitrary number of N copies. We also derive a semi-analytical expression for the asymptotic limit of the iterative distillation and purification protocol and discuss its properties.


I. INTRODUCTION
Continuous variable (CV) quantum information processing [1] is complementary to discrete variable quantum information based on qubits. It employs the encoding of information into field modes of light [2,3,4] or into the collective spin state of many particles such as in a cloud of atoms [5,6,7]. Encoding quantum information into continuous variables of optical modes exhibits several advantages over approaches based on single-photon qubits. For example, it allows for the deterministic realization of several important protocols, such as the generation of entanglement [8,9,10,11] and entanglement swapping [12], quantum teleportation [2,3], quantum cloning [13], and quantum dense coding [14], using only linear optics, optical parametric oscillators as squeezed light sources and balanced homodyne detection. It is also possible to store continuous variables of light in atomic memories, either by means of direct interaction [6] or by teleportation [7].
In the description of CV systems, an important role is played by probability quasidistributions, which represent each field mode by a real two-dimensional function [15]. When this function is of Gaussian shape, the corresponding state is referred to as Gaussian. Due to the ease of both experimental and theoretical treatment of such states, they are, together with Gaussian operations preserving the Gaussian nature of the quasidistributions, a significant part of CV quantum information. However, attempts to estabilish entanglement between two distant parties, which is one of the fundamental tasks quantum information is facing, are hampered by decoherence and thus, purification and distillation of entanglement are needed to remove the added noise and restore, at least approximately, the initial entangled state [16,17]. Indeed, the entanglement distillation and purification techniques allow to extract from many copies of weakly entangled mixed state a single copy of highly entangled almost pure state by means of local quantum operations on each party's side and classical communication between the two parties sharing the states. In general, distillation refers to process of increasing entanglement (or squeezing), while purification is associated with reduction of the state's mixedness. Within the framework of this paper both effects involve each other and the corresponding terms will be used interchangeably. Unfortunately, purification and distillation of entanglement are impossible to accomplish for Gaussian states by using feasible Gaussian operations alone [18,19,20]. Note, that a similar no-go theorem also holds for single mode squeezed states [21]. Namely, it is impossible to use passive Gaussian operations, balanced homodyne detection and feedforward to transform an arbitrary number of Gaussian states into a single state with squeezing better than the initial one, where the squeezing is characterized by the lowest eigenvalue of the covariance matrix. On the other hand, Gaussian operations are sufficient to improve purity of a squeezed state if loss of squeezing is accepted [22].
To allow for increasing the entanglement of Gaussian states by local means, one has to at one point step out of the Gaussian domain, as was utilized in [23,24,25]. The complete distillation protocol, proposed by [25,26,27], consists of two steps: de-Gaussification by a single photon subtraction [24,28,29,30,31] followed by Gaussification by means of interference on beamsplitters, measurement and conditioning. However, if the initial states are already non-Gaussian, for example due to the effect of non-Gaussian noise, it is possible to distill and purify them by employing only linear optics, homodyne detection and post-selection [32,33,34,35]. Recently, we have proposed [32] and experimentally demonstrated [33,35] distillation of squeezing from states disturbed by phase diffusion noise. This noise, caused by random fluctuations of optical phase, commonly occurs in optical communication links. The distillation protocol is based on the interference of two copies of the state at a balanced beamsplitter, followed by homodyne detection of one mode triggering the acceptance or rejection of the other mode that contains the purified state [25,26,32]. The procedure can be straightforwardly extended to accommodate two-mode squeezed states and can also serve as the Gaussification part of the universal CV distillation protocol based on de-Gaussification and subsequent Gaussification.
In this paper we provide detailed theoretical numerical treatment of the distillation and purification of phase diffused single-mode squeezed states. We consider the standard iterative purification protocol [32] as well as the simultaneous collective purification of an arbitrary number of copies of the decohered state. We also discuss the possibility of conditioning on the measurement of an arbitrary quadrature operator, up to abandoning the concept of definite quadrature completely and using phase-randomized homodyne detection. Finally we derive a semi-analytical expression for the state obtained in the limit of infinite number of iterations of the distillation procedure and discuss its properties.

II. ITERATIVE PURIFICATION
The single step of the purification protocol for singlemode squeezed states is depicted in Fig. 1. It consists of interference on a balanced beamsplitter (BS) followed by homodyne measurement of a single mode, detecting a value of the quadrature operator that was, prior to effects of the noise, squeezed [32,33]. If the absolute value of the measured quadrature x is lower than a certain preselected threshold X then the purification succeeded and otherwise it failed. The purified state is present in the other output port of the beamsplitter and is available for further applications or for the next iteration of the purification procedure, whereby two purified copies of the state are purified again. The success of the procedure can be verified by suitable homodyne measurements on the purified state, as schematically shown in Fig. 1. We will treat the purification in the Fock basis which allows us to calculate the full density matrix of the purified state after k iterations of the purification protocol, where |n denotes the n-photon Fock state of the field mode. The initial density matrix of a phase-diffused squeezed vacuum state can be expressed as Here ρ G denotes density matrix of the initial Gaussian state, Φ(φ) denotes the probability distribution of the random phase shift φ and U (φ) is a unitary phaseshift operator with matrix elements m|U |n = e inφ δ mn , where the δ mn stands for the Kronecker delta. We can rewrite Eq. (2) as follows, where ρ G,m,n = m|ρ G |n and is a Fourier transformation of Φ(φ). Throughout this paper we assume that the phase fluctuations have Gaussian distribution, On inserting this probability distribution into Eq. (4) we obtain f n = e −n 2 σ 2 /2 . The Wigner function corresponding to the squeezed vacuum state ρ G reads where V x and V p denote the variances of the x and p quadratures, respectively. For a vacuum state one has V x = V p = 1 2 and the state is squeezed when V x < 1 2 or V p < 1 2 . Note that due to the proper normalisation, all variables used throughout the paper are dimensionless. Presently, these values could go as far as V x ≈ 0.06 -0.05 (−9dB to −10dB), V p ≈ 16 -25 (15dB -17dB) [36,37], but values around V x ≈ 0.16 -0.22 (−5dB to −3dB) are common in contemporary experiments [33,34,35].
The density matrix elements ρ G,m,n can be evaluated by noting that the Husimi Q-function defined as Q(α) = 1 π α|ρ|α , where |α denotes a coherent state with amplitude α, is a generating function of the density matrix elements in the Fock state basis, The Q-function corresponding to the Wigner function (6) is also Gaussian, On inserting function (8) into expression (7) we obtain after some algebra if m + n is even and ρ G,m,n = 0 if m + n is odd. The distillation and purification scheme shown in Fig. 1 produces with certain probability from two copies of the state ρ (k−1) a single copy of a purified state ρ (k) and we can formally write the purification map as As the formula suggests, the protocol can be iterated and the outputs of each step can be used as inputs for the next iteration. As we shall show, each iteration increases the squeezing and purity of the state and Gaussifies it. The interference of two copies of the state ρ (k−1) on a balanced beamsplitter is governed by a unitary transformation. Since this transformation preserves the total photon number, we can write it in the Fock basis as follows, where After mixing on a BS the first output mode is measured in a balanced homodyne detector where it is projected on the eigenstate of the x quadrature |x . The positive operator valued measure (POVM) element corresponding to conditioning on |x| ≤ X reads C = X −X |x x| and its matrix elements in the Fock basis can be expressed as This can be easily evaluated numerically by recalling that where H n (x) denotes the Hermite polynomial. For any realistic detector with efficiency η < 1 we can use the model where the inefficient detector is replaced by a beamsplitter with transmittance η whose auxiliary input port is in the vacuum state and which is followed by an ideal perfect detector [15]. This yields where We are now in a position to combine all the above expressions and write down the purification map (11) in the Fock state basis, ×C m1+a,n1+b (η) |m 2 − a n 2 − b|. The density matrix (18) is not normalized and its trace is equal to the probability of success (success rate) of k iterations of the purification protocol, Typical numerical results are shown in Fig. 2. We can see that each step of the iterative purification decreases the variance of x, i.e. increases the squeezing. Also, the purity of the state P = Tr[ρ 2 ]/(Tr[ρ]) 2 increases after each iteration of the protocol and approaches some asymptotic value, which is generally less than unity, c.f. Fig 2(b). So our protocol purifies the state but it does not generally distill perfectly pure states from the initially mixed states. One important parameter determining the practical feasibility and usefulness of the scheme is the probability of success of the protocol, which is plotted in Fig. 2(c). We can see that for the chosen acceptance window X = 0.45 the success rate P (N ) is very high, the probability of success of each iteration is about 50% and the total probability of success for a protocol involving four iterations is still about 10%. The trade-off between the performance of the protocol and the success probability is illustrated in Fig. 3 which shows the final variance and purity of the state after four iterations as a function of the total success rate P . Higher purity and stronger squeezing are achievable at the expense of reduced P . Our numerical calculations suggest that it is suitable to choose the threshold X ≈ √ V x , which achieves good purification and squeezing enhancement at a reasonably high probability of success.
The purification procedure also Gaussifies the state. To quantify this, we evaluate the fidelity of the state ρ (k) after k iterations of the protocol with the Gaussian stateρ G that has the same mean quadrature values and covariance matrix as the state ρ (k) . The fidelity of two mixed states is defined as follows, and it holds that F = 1 if and only ifρ G = ρ (k) . A value F < 1 is a clear signature of a non-Gaussian character of the state ρ (k) . The results are shown in Fig. 4 which confirms that the present protocol indeed Gaussifies the state. After four iterations, we have F > 0.999 so the state is almost perfectly Gaussified. So far we have considered only conditioning on the measurement of the initially squeezed x quadrature. However, the method works even for a measurement of an arbitrary quadrature q(θ) = x cos θ + p sin θ [35], where the conditioning is again of the form |q(θ)| < Q. This can be straightforwardly implemented into our Fock-state basis formalism by re-defining the matrix elements of the POVM C(η) as follows, Note, however, that there are two prominent quadratures, x = q(0) and p = q(π/2). Conditioning on one of these quadratures yields, for separate ranges of σ, the best improvement of the purified state. Furthermore, since improvement in the purified state can be seen no matter which quadrature q(θ) is measured [35], it is possible to consider also a phase-randomized homodyning. In this case, Eq. (21) has to be averaged over the random phase θ and only diagonal elements survive, C rand m,n (η) = C m,n (η)δ m,n . The results of numerical simulations of the iterative purification are given in Fig. 5. Note that the advantage of conditioning on p grows with the number of iterations. For example, if σ = 0.5 then after the first iteration the variance is more reduced by conditioning on |x| < Q than on |p| < Q. However, for a four-step iterative procedure this is reversed and it is slightly better to condition on |p| < Q rather than on |x| < Q. The phase-randomised purification can be seen as a compromise between conditioning on |x| < Q and |p| < Q -for large (small) values of σ it yields results better (worse) than conditioning on p (x) but worse (better) than conditioning on x (p). In fact, values obtained by phase-randomised purification are almost identical to those produced by random switching between x and p detections.

III. COLLECTIVE PURIFICATION
The iterative purification scheme analyzed in the previous section requires at least 2 k copies for k iterations of the protocol. Here we analyze a more general collective purification scheme that works for any number N of available copies of the state. In particular, this protocol can be employed when three copies of the state are available and it reduces to the previous scheme when N = 2 k . The proposed setup is shown in Fig. 6. The N copies of the phase-diffused squeezed state are combined on an array of N − 1 beamsplitters with transmittances t j and reflectances r j , t 2 j + r 2 j = 1. All output modes except for the last one are monitored by balanced homodyne detectors which measure the x quadrature of each mode. The purification is successful if the absolute value of each measurement outcome is below the corresponding threshold, |x j,out | ≤ X j . The choice of the thresholds would depend on the transmittances of the beam splitters. Here we shall consider ideal homodyne detectors with unit efficiency, η = 1, and conditioning on |x j,out | = 0 for which the protocol achieves the best performance. In practice, this would correspond to choosing sufficiently narrow acceptance windows, i.e. very small thresholds X j .
Let us for a moment fix the random phase shifts φ j of each mode. Then the variance of the quadrature x j of the j-th mode reads and in what follows we shall not explicitly display the dependence of V j on φ j for the sake of brevity. The joint probability distribution of the x quadratures of all N input modes factors into the product of marginal Gaussian probability distributions, In the Heisenberg picture, the output x quadratures after interference on the array of beamsplitters can be expressed as linear combinations of the input x quadratures, where U j,l denote elements of a real orthogonal matrix which describes the action of the interferometer composed of the beamsplitters. For the scheme shown in Fig. 6 it holds that Here we formally define r 0 = 1. The input-output relation (25) can be inverted and we have The joint distribution of the output quadratures can be obtained by inserting the x j expressed in terms of x l,out in Eq. (24). Taking into account that the Jacobian of the transformation (25) is det U = 1, we obtain (28) The un-normalized distribution of x N,out conditional on projections x l,out = 0, l = 1, . . . , N − 1 then reads The variance of x N,out can be calculated by averaging properly weightedṼ over the N random phase shifts, Here dφ = dφ 1 · · · dφ N and N N is a normalization constant, Two choices of the beamsplitting ratios are of particular interest. Using balanced beamsplitters is experimentally most straightforward, since these beamsplitters are most commonly used in the lab. In this case we have t j = r j = 2 −1/2 which implies U N,1 = 2 −(N −1)/2 and U N,l = 2 −(N −l+1)/2 , l ≥ 2. Notice that the output quadrature x N,out is an unbalanced combination of the input quadratures. The second choice is to use unbalanced beamsplitters which yield a balanced superposition of the input quadratures, U N,l = 1/ √ N . This happens if t j = j/(j + 1). This latter scheme becomes fully equivalent to the k-step iterative purification procedure discussed in the previous section when N = 2 k . Note that in the iterative scheme we can also first let all the beams interfere on the balanced beamsplitters and then measure on all output ports except one. If we condition on x l,out = 0 then the two schemes become immediately fully equivalent. For realistic finite thresholds X j the schemes become equivalent if conditioning on appropriate linear combinations of x l,out is adopted. The results of numerical simulations of the collective purification procedure are shown in Fig. 7. We can see that each additional copy helps to further suppress the noise and reduce the quadrature fluctuations. The difference between the performances of the two-copy and three-copy protocols is sufficiently large so that it should be possible to observe it in present-day experiments. The numerics also reveal that the difference between the schemes using balanced and unbalanced beamsplitters is very small and the scheme is quite robust and insensitive to changes in reflectances or transmittances of the beamsplitters of the order of a few percent, which is an important practical advantage.
A conditional measurement in the form of a detection of arbitrary quadrature q(θ) can be also considered for the collective purification. To treat in a unified way conditioning on measurements of arbitrary quadratures of the first N − 1 output modes, we shall present a more detailed treatment of the scenario, based on the formalism of covariance matrices (CM). Let ξ = (x 1 , p 1 , x 2 , p 2 , . . . , x N , p N ) denote the vector of quadrature operators. The covariance matrix Σ jk = 1 2 {∆ξ j , ∆ξ k } comprises the second moments. The CM of the single-mode squeezed state is diagonal, For fixed random phase shifts φ j , the covariance matrix of the N-mode Gaussian state impinging on the array of N − 1 beamsplitters is given by The covariance matrix of the output modes can be obtained from Σ in by a symplectic transformation, where the symplectic matrix S BS governs the linear mixing of the quadrature components in the passive linear interferometer consisting of the array of N − 1 beamsplitters, see Fig. 6. To make our treatment simple, we shall again consider the limiting case of conditioning on q j = 0, j = 1, . . . , N − 1. The quadratures which are being observed are specified by N − 1 relative phases θ j and we assume that all detectors have the same efficiency η. We introduce a new covariance matrix where I stands for the identity matrix, S(θ) = N j=1 R(θ j ) and we set θ N = 0. The single-mode phase shifts in the above formula map all measured quadratures onto the x quadratures. We want to calculate the variance of the quadrature x N conditional on all measurement outcomes being equal to 0. We first construct a submatrix of Σ aux which comprises variances and covariances of the x quadratures only, Here S x is a N ×2N matrix defined as S j,2j−1 = 1 and all other elements are equal to zero. The joint probability distribution of the N quadratures x = (x 1 , x 2 , . . . , x N ) can be expressed as where the |.| denotes determinant of the matrix. The probability distribution of x N conditional on is the conditional variance of x N . Since we assume that the mean values of all quadratures are initially zero and the purification scheme preserves this property, the resulting variance of x N after purification can be evaluated as properly weighted average of V N , Here dφ = dφ 1 · · · dφ N and N N is a normalization constant, Note that V out determined in this way corresponds to the variance of the purified beam that would be measured by balanced homodyne detector with efficiency η. The expressions (40) and (41) generalize the formulas (31) and (32). The numerical results for three-copy purification are plotted in Fig. 8. This figure confirms the general trend that for weak phase fluctuations it is advantageous to condition on measurements of p while for strong noise it is better to measure the quadrature x [35]. The intermediate strategy where the first balanced homodyne detector (BHD) measures x while the second BHD measures p does not bring any advantage so we can conclude that the optimum strategy that provides maximum reduction of quadrature variance consists of measuring the same quadrature (either x or p) by both homodyne detectors.

IV. ASYMPTOTIC LIMIT
It is possible to find a state to which the purification procedure (11) converges in the limit of infinite number of iterations, k → ∞. This asymptotic state can serve as a reference for the performance of the purification, giving values that can never be surpassed by finite number of iterations. Due to the nature of the result we seek, we shall consider the post-selection threshold X in Eq. (14) to be zero. That is, given two copies of a noisy state ρ in,j , j = 1, 2, the state after a single step of the purification is given by the map where U 12 denotes a unitary operator describing balanced beamsplitter coupling modes 1 and 2 and |x = 0 2 x = 0| is a projector on the eigenstate of the x quadrature. In the following we shall exploit the possibility of expressing the x eigenstate as where the operator S r = exp[r(a 2 − a †2 )/2] stands for finite squeezing with parameter r and |0 is the vacuum state [27]. Since the operator S r,1 ⊗ S r,2 commutes with the beamsplitter transformation U 12 , we can rewrite the relation (42) as with ρ ′ in,jr = S † r ρ in,j S r and where the purification map E |0 is the same as in Ref. [26]. This is schematically demonstrated in Fig. 9. Note that the inverse squeezing operation performed on the output state cancels out with the initial squeezing of the next purification step. Thus, the iterative purification scheme, using projection onto the quadrature eigenvectors in all steps, is equivalent to iterative purification with projection onto vacuum states, accompanied by squeezing of all initial states and inverse squeezing of the purification outcome [27], see Fig. 10. This allows us to find the asymptotic limit, analogously as in Ref. [26]. Due to this similarity, we shall mention here only the most important parts of the derivation and the reader can find more details in Ref. [26].
The key realization in finding the limit is that stationary points of the purification map E |0 (ρ ⊗ ρ), i.e. states that remain unchanged by the map, are Gaussian states with zero mean value. The initial non-Gaussian state will, over the course of purification, converge to one of these states. As another important observation, let us note that the coefficients σ 11 , σ 20 and σ 02 defined as σ ij = ρ ij /ρ 00 are also unchanged by the purification map E |0 and there is a one-to-one correspondence between these three coefficients and the set of covariance matrices, which in turn completely describe any singlemode Gaussian state with zero mean.
The procedure of finding the limit is the following: first, the initial states are described by the Wigner func- tion where the vector ξ = (x, p) comprises the two conjugate quadratures and and R(φ) is defined in Eq. (33). The Q-function corresponding to this Wigner function reads with Γ r,φ = (Σ r,φ + I/2) −1 , where I stands for the identity matrix. Now, with the substitution x = (α+α * )/ √ 2, p = (α − α * )/i √ 2 we can use relation (7) and find the density matrix elements where φ denotes averaging over the distribution Φ(φ). From the density matrix elements we can obtain the coefficients σ 11 , σ 20 and σ 02 , σ ij = ρ ij /ρ 00 , which are invariant over the course of the purification and, therefore, describe the Gaussian limit to which the state converges.
It is now a simple matter of finding the inverse relations After some algebra we find that the formulas (48) and (49) can be simplified to The final covariance matrix is given by Σ = S −1 r (Γ −1 − I/2)S −1 r . Due to the symmetric nature of the phase noise, the matrix Σ is diagonal and fully specified by values of variances of the squeezed and antisqueezed quadratures, V x,lim and V p,lim . We can directly find the output variances by taking the limit r → ∞, obtaining where Note that V x,lim V p,lim = V x V p so the purification procedure asymptotically suppresses all added noise and the asymptotic Gaussian state has the same purity P = 1/(2 V x V p ) as the initial state before de-phasing. This, however, holds only for idealized conditioning on x = 0. For finite acceptance window, X > 0, the purity of the asymptotic distilled state would generally be lower than the purity of the initial state. It is possible to modify the method to incorporate other scenarios discussed in previous sections. Finding the limit state for iterations incorporating conditioning measurements along an arbitrary quadrature q(θ) is straightforward, requiring only the use of an alternative squeezing operation described by S r,θ = R(θ)S r R T (θ). Adapting the calculations to incorporate for imperfect detectors is a bit more tricky. Inefficient detection with efficiency η can be simulated by a beamsplitter with transmittance η followed by a perfect detector. The attenuation can be represented by a map, M(ρ), acting on the density matrix and transforming the corresponding covariance matrix of the Gaussian state as, The inverse operation, M −1 (ρ) acting as is, of course, nonphysical. We can, however, formally use it and in analogy with Eq. (44) represent a single step of the purification procedure by the relation where the ε is given by Now again, starting from the modified input states S † r M(ρ)S r we can find the asymptotic limit for the map E |0 and with the use of relation (55) it can be transformed into the desired result. However, if we are interested in results that can be observed by imperfect detection, the final operation M −1 is not necessary.
Some numerical results are shown in Fig. 11. It can be seen that the iterations converge to the limit given by Eq. (51) in both the resulting variance and purity. Moreover, the purity of the limit is independent of the actual amount of initial phase fluctuations and equals to the purity of the original state (modified by imperfect detectors), although it is less when the states are compared directly, without the imperfect verification measurement. Figures 12 and 13 illustrate the behavior of the asymptotic limit relative to the angle θ of the conditioning measurement. Although all choices of θ lead to a state purified to some extent, the two most interesting choices are θ = 0, corresponding to the measurement of the quadrature operator x and θ = π/2, corresponding to the measurement of p. Only these two choices can lead to a purity of the final state to be maximal, that is, equal to the purity of the initial state (modified by detection efficiency). This is not very surprising, since these two choices coincide with the basis in which the covariance matrix of the initial state is diagonal. We can see that for weak phase noise (small σ) it is advantageous to employ conditioning on measurements of the p quadrature while for strong phase fluctuations the x based post-selection is optimal.

V. SUMMARY
We have provided a detailed analysis of multiple copy purification protocols for phase diffused single-mode squeezed states. We performed most calculations in the Fock basis, which allowed for efficient numerical treatment. We have also discussed the performance of collective purification, employing more than two copies of the initial state in one step. The collective purification/distillation uses an arbitrary number of copies to improve the purity and non-classicality of the state, and therefore was found to be a generalization of iterative purification and distillation which requires 2 k copies. In particular we have shown that collective distillation based on three copies of phase diffused squeezed states is possible, and already provides a significantly improved degree of squeezing when compared with two-copy distillation. We have also investigated the convergence of the iterative purification protocol and derived semianalytical expressions for variances of the asymptotic Gaussian state. Finally, let us note that in the spirit of ref. [32], the results presented here can be straightforwardly extended to describe the distillation of phase diffused two-mode squeezed states.