Julien Laurat

Gaëlle Keller

José Augusto Oliveira-Huguenin

Claude Fabre

Thomas Coudreau
email: coudreau@spectro.jussieu.fr

Alessio Serafini♮♯

Gerardo Adesso♮

Fabrizio Illuminati♮

Entanglement of two-mode Gaussian states: characterization and experimental production and manipulation

Keywords: numbers: 03.67.Mn, 42.65.Yj, 42.50.Dv, 42.50.Lc

A powerful theoretical structure has emerged in recent years on the characterization and quantification of entanglement in continuous-variable systems. After reviewing this framework, we will illustrate it with an original setup based on a type-II OPO with adjustable mode coupling. Experimental results allow a direct verification of many theoretical predictions and provide a sharp insight into the general properties of two-mode Gaussian states and entanglement resource manipulation.

Introduction

Quantum information aims at the treatment and transport of information using the laws of quantum physics. For these goals, continuous variables (CV) of the electromagnetic field have emerged as a powerful tool [START_REF]Quantum information with Continuous Variables[END_REF][START_REF] Braunstein | Quantum information with continuous variables[END_REF][START_REF]Quantum Information with Continuous Variables of Atoms and Light[END_REF]. In this context, entanglement is an essential ressource. The purpose of this paper is to link the very powerful mathematical description of gaussian states based on covariance matrices and the experimental production and manipulation of entanglement.

Experimentally, entanglement can be obtained directly by type-II parametric interaction below the oscillation deamplifying either the vacuum fluctuations as was demonstrated in the seminal experiment by Ou et al. [START_REF] Ou | Realization of the Einstein-Podolsky-Rosen paradox for continuous variables[END_REF] (or in recent experiments [START_REF] Laurat | Compact source of Einstein-Podolsky-Rosen entanglement and squeezing at very low noise frequencies[END_REF]) or the fluctuations of a weak injected beam [START_REF] Zhang | Experimental generation of bright two-mode quadrature squeezed light from a narrow-band nondegenerate optical parametric amplifier[END_REF]. It can also be obtained above the oscillation threshold in a non frequency degenerate configuration [START_REF] Villar | Generation of Bright Two-Color Continuous Variable Entanglement[END_REF]. It can also be obtained indirectly by mixing on a beam splitter two independent squeezed beams. The required squeezing can be produced by Kerr effects -using optical fibers [START_REF] Ch | Generation of continuous variable Einstein-Podolsky-Rosen entanglement via the Kerr nonlinearity in an optical fiber[END_REF] or cold atoms in an optical cavity [START_REF] Josse | Continuous Variable Entanglement using Cold Atoms[END_REF] -or by type-I parametric interaction in a cavity [START_REF] Wu | Generation of squeezed states by parametric down conversion[END_REF][START_REF] Bowen | Experimental investigation of continuous-variable quantum teleportation[END_REF]. Single-pass type-I interaction in a non-colinear configuration can also generate directly entangled beams as demonstrated recently by Wenger et al. in the pulsed regime [START_REF] Wenger | Time-resolved homodyne characterization of individual quadrature-entangled pulses[END_REF]. All these methods produce a symmetric entangled state enabling dense coding, the teleportation of coherent [START_REF] Furusawa | Unconditional quantum teleportation[END_REF][START_REF] Bowen | Experimental investigation of continuous-variable quantum teleportation[END_REF][START_REF] Takei | High-Fidelity Teleportation beyond the No-Cloning Limit and Entanglement Swapping for Continuous Variables[END_REF] or squeezed states [START_REF] Takei | Experimental demonstration of quantum teleportation of a squeezed state[END_REF] or entanglement swapping [START_REF] Jia | Experimental demonstration of unconditional entanglement swapping for continuous variables[END_REF][START_REF] Takei | High-Fidelity Teleportation beyond the No-Cloning Limit and Entanglement Swapping for Continuous Variables[END_REF]. These experiments generate an entangled state with a covariance matrix in the so-called 'standard form' [START_REF] Duan | Inseparability criterion for continuous variable systems[END_REF][START_REF] Simon | Peres-Horodecki separability criterion for continuous variable systems[END_REF], without having to apply any local linear unitary transformations such as beam-splitting or phase-shifts to exploit it optimally in quantum information protocols.

However, it has been recently shown in [START_REF] Laurat | Effects of mode coupling on the generation on Einstein-Podolsky-Rosen entanglement in a type-II OPO below threshold[END_REF] that, when a birefringent plate is inserted inside the cavity of a type-II optical parametric oscillator, i.e. when mode coupling is added, the generated two-mode state remains symmetric but entanglement is not observed on orthogonal quadratures: the state produced is not in the standard form. The entanglement of the two emitted modes in this configuration is not optimal: it is indeed possible by passive "global" operations to select modes that are more entangled. Our original system provides thus a good insight into the quantification and manipulation of the entanglement resources of two-mode Gaussian states. In particular, as just anticipated, it allows to confirm experimentally the theoretical predictions on the entangling capacity of passive optical elements and on the selection of the optimally entangled bosonic modes [START_REF] Wolf | Entangling power of passive optical elements[END_REF].

The paper is organized as follows: we start by giving a general overview of gaussian states and defining the covariance matrix formalism, focusing in particular on two-mode states. Over the last years, a great deal of attention has been devoted to defining not only criterion of entanglement but also quantification of this entanglement. Section 3 focuses on different such measures which are interpreted in the covariance matrix formalism. Effects of mode coupling on the generated entanglement is then discussed in detail (section 4). The experimental setup used to generate and quantify the entanglement is presented in section 5. In section 6, experimental measures of entanglement are realized and we discuss on the effect of noise of the covariance matrix on the entanglement measures. Finally, we show how optimization of the resource is obtained by a passive operation -namely a polarization adjustment -operated on the two-mode state (section 7).

Gaussian states: general overview

A continuous variable (CV) system is described by a Hilbert space H = n i=1 H i resulting from the tensor product of infinite dimensional Fock spaces H i 's. Let âi and â † i be the annihilation and creation operators acting on H i (ladder operators), and xi = (â i + â † i ) and pi = (â iâ † i )/i be the related quadrature phase operators. The corresponding phase space variables will be denoted by x i and p i . Let X = (x 1 , p1 , . . . , xn , pn ) denote the vector of the operators xi and pi . The canonical commutation relations for the Xi can be expressed in terms of the symplectic form Ω

[ Xi , Xj ] = 2iΩ ij , with Ω ≡ n i=1 ω , ω ≡ 0 1 -1 0 .
The state of a CV system can be equivalently described by a positive trace-class operator (the density matrix ̺) or by quasi-probability distributions. In the following, we shall adopt the Wigner quasi-probability function W (R) defined, for any density matrix, as the Fourier transform of the symmetrically ordered characteristic function [START_REF] Barnett | Methods in theoretical quantum optics[END_REF]:

W (R) ≡ 1 π 2 Ê 2n Tr [̺D Y ] e iY T ΩR d 2n Y , R ∈ Ê 2n ,
where the displacement operators D Y (describing the effect on the field of a classical driving current) are defined as

D Y = e iY T Ω X , Y ∈ Ê 2n . (1)
The 2n-dimensional space of definition of the Wigner function, associated to the quadratic form Ω, is referred to as phase space, in analogy with classical Hamiltonian dynamics. In Wigner phase space picture, the tensor product H = H i of the Hilbert spaces H i 's of the n modes results in the direct sum Γ = Γ i of the phase spaces Γ i 's.

States with Gaussian characteristic functions and quasi-probability distributions are referred to as Gaussian states. Such states are at the heart of information processing in CV systems [START_REF] Braunstein | Quantum information with continuous variables[END_REF] and are the subject of our analysis. By definition, a Gaussian state ̺ is completely characterised by the first and second statistical moments of the field operators, which will be denoted, respectively, by the vector of first moments X ≡ X1 , X2 , . . . , X2n-1 , X2n and the covariance matrix

(CM) Γ of elements Γ ij ≡ 1 2 Xi Xj + Xj Xi -Xi Xj , (2)
where, for any observable ô, the expectation value ô ≡ Tr(̺ô). Notice that, according to our definition of the quadrature operators in terms of the ladder operators, the entries of the CM are real numbers. Being the variances and covariances of quantum operators, such entries are obtained by noise variance and noise correlation measurements. They can be expressed as energies by multiplying them by the quantity ω, where ω is the frequency of the considered mode. In fact, for any n-mode state the quantity ωTr (Γ/4) is just the contribution of the second moments to the average of the "free" Hamiltonian

n i=1 (â † i âi + 1/2).
Coherent states, resulting from the application of displacement operators D X to the vacuum state |0 , are Gaussian states with CM Γ = ½ and first statistical moments X = X. First moments can be arbitrarily adjusted by local unitary operations (displacements), which cannot affect any property related to entropy or entanglement. In the present experimental case, the fields are produced by an OPO below the oscillation threshold, and the means are actually zero. Let us note that this can be done without any loss of generality.

The canonical commutation relations and the positivity of the density matrix ̺ imply

Γ + iΩ ≥ 0 , (3)
meaning that all the eigenvalues of the (hermitian) matrix Γ + iΩ have to be greater or equal than zero. Inequality (3) is the necessary and sufficient constraint the matrix Γ has to fulfill to be a CM corresponding to a physical Gaussian state [START_REF] Simon | Gaussian-Wigner distributions in quantum mechanics and optics[END_REF][START_REF] Simon | Quantum noise matrix for multimode systems: U(n) invariance, squeezing and normal forms[END_REF]. Note that the previous condition is necessary for the CM of any (generally non Gaussian) state, as it generalises to many modes the Robertson-Schrödinger uncertainty relation, reducing to the familiar Heisenberg principle for pure, uncorrelated states. We mention that such a constraint implies Γ ≥ 0. A major role in the theoretical and experimental manipulation of Gaussian states is played by unitary operations which preserve the Gaussian character of the states on which they act. Such operations are all those generated by terms of the first and second order in the field operators. So, beside the already discussed displacements, also the unitary operations generated by terms of the second order are Gaussian. As a consequence of the Stone-Von Neumann theorem, any such operation at the Hilbert space level corresponds, in phase space, to a symplectic transformation, i.e. to a linear transformation S which preserves the symplectic form Ω, so that Ω = S T ΩS, i.e. it preserves the commutators between the different operators. Symplectic transformations on a 2n-dimensional phase space form the (real) symplectic group, denoted by Sp (2n,Ê) . Such transformations act linearly on first moments and "by congruence" on covariance matrices (i.e. so that Γ → S T ΓS). One has Det S = 1, ∀ S ∈ Sp (2n,Ê) .

Ideal beam splitters, phase shifters and squeezers are described by symplectic transformations. In fact single and two-mode squeezings occurring, respectively, in degenerate and non degenerate parametric down conversions, are described by the operators

U ij,r,ϕ = e 1 2 (εâ † i â † j -ε * âi âj ) with ε = r e i2ϕ , r ∈ Ê, ϕ ∈ [0, 2π], (4)
resulting in single-mode squeezing of mode i for i = j. The representation in phase space of the operation U ij,r,ϕ for i = j is given by the linear transformation S ij,r,ϕ

S ij,r,ϕ =     c -hs 0 ks 0 0 c + hs 0 -ks ks 0 c + hs 0 0 -ks 0 c -hs     for i = j , (5)
where c = cosh(2r), s = sinh(2r), h = cos(2ϕ), k = sin(2ϕ) (r and ϕ are same as in Eq. ( 4)) and the matrix is understood to act on the couple of modes i and j. Beam splitters are described by the operators

O ij,θ = e θâ † i âj -θâiâ † j , θ ∈ [0, 2π], (6)
corresponding to symplectic rotations R ij,θ in phase space

R ij,θ =     cos(θ) 0 -sin(θ) 0 0 cos(θ) 0 -sin(θ) sin(θ) 0 cos(θ) 0 0 sin(θ) 0 cos(θ)     for i = j . (7)
The angle θ is defined by Eq. ( 6) and, again, the matrix is understood to act on the i and j modes. Single-mode symplectic operations are easily retrieved as well, being just combinations of single mode (two dimensional) rotations and of single mode squeezings of the form diag ( e r , e -r ) for r > 0. Now, symplectic transformations in phase space are generated by exponentiation of matrices written as JΩ, where J is antisymmetric [START_REF] Arvind | The Real Symplectic Groups in Quantum Mechanics and Optics[END_REF]. Such generators can be symmetric or antisymmetric. The operations R ij,θ , generated by antisymmetric operators are orthogonal and, acting by congruence on the CM Γ, preserve the value of Tr Γ. Since Tr Γ gives the contribution of the second moments to the average of the Hamiltonian i a † i a i , these transformations are said to be 'passive', or 'energy preserving' (they belong to the compact subgroup of Sp (2n,Ê) ). Instead, operations S ij,r,ϕ , generated by symmetric operators, are not orthogonal and do not preserve Tr Γ (they belong to the non compact subgroup of Sp (2n,Ê) ). This mathematical difference between squeezers and phase space rotations accounts, in a quite elegant way, for the difference between 'active' (i.e. energy consuming) and 'passive' (i.e. energy preserving) optical transformations.

Finally, let us recall that, due to Williamson theorem [26], the CM of a n-mode Gaussian state can always be written as

Γ = S T νS , (8)
where S ∈ Sp (2n,R) and ν is the CM

ν = diag(ν 1 , ν 1 , . . . , ν n , ν n ) , (9)
corresponding to a tensor product of thermal states with diagonal density matrix ̺ ⊗ given by

̺ ⊗ = i 2 ν i + 1 ∞ k=0 ν i -1 ν i + 1 |k ii k| , (10)
|k i being the k-th number state of the Fock space H i . The dual (Hilbert space) formulation of Eq. ( 8) then reads: ̺ = U † ̺ ⊗ U , for some unitary U . The quantities ν i 's form the symplectic spectrum of the CM Γ and can be computed as the eigenvalues of the matrix |iΩΓ| [27]. Such eigenvalues are in fact invariant under the action of symplectic transformations on the matrix Γ.

The symplectic eigenvalues ν i encode essential informations on the Gaussian state Γ and provide powerful, simple ways to express its fundamental properties. For instance, in terms of the symplectic eigenvalues ν i , the uncertainty relation (3) simply reads

ν i ≥ 1 . ( 11
)
Also the entropic quantities of Gaussian states can be expressed in terms of their symplectic eigenvalues and invariants. Notably, the purity Tr ̺ 2 of a Gaussian state ̺ is simply given by the symplectic invariant Det Γ = n i=1 ν i , being Tr

̺ 2 = (Det Γ) -1/2 [28].

Two-mode states

Since this work is focused on two-mode Gaussian states, we briefly review here some of their basic properties. The expression of the two-mode CM Γ in terms of the three 2 × 2 matrices α, β, γ will be useful

Γ ≡ α γ γ T β . ( 12
)
For any two-mode CM Γ there exist local symplectic operations S 1 and S 2 (each S j acting on one of the two modes), such that their direct sum S l = S 1 ⊕S 2 (corresponding to the tensor product of local unitary operations) brings Γ to the so called standard form Γ sf [START_REF] Simon | Peres-Horodecki separability criterion for continuous variable systems[END_REF][START_REF] Duan | Inseparability criterion for continuous variable systems[END_REF]]

S T l ΓS l = Γ sf ≡     a 0 c + 0 0 a 0 c - c + 0 b 0 0 c -0 b     . (13)
States whose standard form fulfills a = b are said to be symmetric. Let us recall that any pure state is symmetric and fulfills

c + = -c -= √ a 2 -1.
The correlations a, b, c + , and c -are determined by the four local symplectic invariants DetΓ = (ab

-c 2 + )(ab -c 2 -), Detα = a 2 , Detβ = b 2 , Detγ = c + c -.
Therefore, the standard form corresponding to any CM is unique (up to a common sign flip in c - and c + ).

For two-mode states, the uncertainty principle Ineq. ( 3) can be recast as a constraint on the Sp (4,R) invariants DetΓ and ∆(Γ) = Detα + Detβ + 2 Detγ:

∆(Γ) ≤ 1 + DetΓ . (14)
The symplectic eigenvalues of a two-mode Gaussian state will be named ν -and ν + , with ν -≤ ν + , with the Heisenberg uncertainty relation reducing to

ν -≥ 1 . (15)
A simple expression for the ν ± can be found in terms of the two Sp (4,R) invariants (invariants under global, two-mode symplectic operations) [START_REF] Vidal | Computable measure of entanglement[END_REF][START_REF] Serafini | Symplectic invariants, entropic measures and correlations of Gaussian states[END_REF] 2ν

2 ± = ∆(Γ) ± ∆(Γ) 2 -4 Det Γ . (16)
A subclass of Gaussian states with a major interest in experimental quantum optics and in the practical realization of CV quantum information protocols is constituted by the nonsymmetric two-mode squeezed thermal states. Let S r = S 12,r,π/4 be the two mode squeezing operator with real squeezing parameter r, and let ̺ ⊗ νi be a tensor product of thermal states with CM Γ ν± = ½ 1 ν -⊕ ½ 2 ν + , where ν ± is the symplectic spectrum of the state. A nonsymmetric two-mode squeezed thermal state ξ νi,r is defined as

ξ νi,r = S r ̺ ⊗ νi S † r , corresponding to a standard form with a = ν -cosh 2 r + ν + sinh 2 r , b = ν -sinh 2 r + ν + cosh 2 r , (17)
c ± = ± ν -+ ν + 2 sinh 2r .
In the symmetric instance (with ν -= ν + = ν) these states reduce to two-mode squeezed thermal states. The covariance matrices of these states are symmetric standard forms with a = ν cosh 2r , c ± = ±ν sinh 2r .

These are the states which occur in most realistic parametric down conversion processes, like the one which will be discussed later on in the paper. In the pure case, for which ν = 1, one recovers the two-mode squeezed vacuum. Such states encompass all the standard forms associated to pure states: any two-mode Gaussian state can be reduced to a squeezed vacuum by means of unitary local operations.

Entanglement of two-mode Gaussian states

This section aims at reviewing the main results on the qualification and quantification of entanglement for Gaussian states of CV systems, which will be exploited in the following.

The positivity of the partially transposed state (PPT criterion) is necessary and sufficient for the separability of two-mode Gaussian states [START_REF] Simon | Peres-Horodecki separability criterion for continuous variable systems[END_REF] (and, more generally, of all (1+n)-mode Gaussian states under 1×n-mode bipartitions [START_REF] Werner | Bound entangled Gaussian states[END_REF] and of bisymmetric (m + n)-mode Gaussian states under m × n-mode bipartitions [START_REF] Serafini | Unitarily localizable entanglement of Gaussian states[END_REF]). In general, the partial transposition ̺ of a bipartite quantum state ̺ is defined as the result of the transposition performed on only one of the two subsystems in some given basis. It can be promptly inferred from the definition of the Wigner function W (X) that the action of partial transposition amounts, in phase space, to a mirror reflection of one of the four canonical variables. In terms of

Sp (2,R) ⊕ Sp (2,R) invariants, this reduces to a sign flip in Det γ. Therefore the invariant ∆(Γ) is changed into ∆(Γ) = ∆( Γ) = Det α + Det β -2 Det γ. Now, the symplectic eigenvalues ν∓ of Γ read ν∓ = ∆(Γ) ∓ ∆(Γ) 2 -4 Det Γ 2 . ( 19
)
The PPT criterion thus reduces to a simple inequality that must be satisfied by the smallest symplectic eigenvalue ν-of the partially transposed state

ν-≥ 1 , (20)
which is equivalent to

∆(Γ) ≤ Det Γ + 1 . (21)
The above inequalities imply Det γ = c + c -< 0 as a necessary condition for a twomode Gaussian state to be entangled. Therefore, the quantity ν-encodes all the qualitative characterization of the entanglement for arbitrary (pure or mixed) twomodes Gaussian states.

Let us now briefly focus on the entanglement qualification of symmetric states, which will be the subject of the experimental investigations presented in the paper. It is immediately apparent that, because a = b, the partially transposed CM in standard form Γ (obtained by flipping the sign of c -) is diagonalized by the orthogonal, symplectic transformation R 12,π/4 of Eq. ( 7), resulting in a diagonal matrix with entries a ∓ |c ∓ |. The symplectic eigenvalues of such a matrix are then easily retrieved by applying local squeezings. In particular, the smallest eigenvalue ν-is simply given by

ν-= (a -|c + |)(a -|c -|) . ( 22
)
Note that also the original standard CM Γ with a = b could be diagonalized (not symplectically, since the four diagonal entries are generally all different) by the same beam splitter transformation R 12,π/4 , with the same orthogonal eigenvalues a∓|c ∓ |. It is immediate to verify that ν-is just given by the geometric average between the two smallest of such orthogonal eigenvalues of Γ. The two quadratures resulting from the previous beam splitter transformation select orthogonal directions in phase space with respect to the original ones, so they will be referred to as 'orthogonal' quadratures. Notice that, in the experimental practice, this allows for the determination of the entanglement through the measurement of diagonal entries (noise variances) of the CM after the application of a balanced beam splitter (which embodies the transformation R 12,π/4 ).

To explore further consequences of this fact, let us briefly recall some theoretical results on the generation of entanglement under passive (energy-preserving) transformations, which will be precious in the following. As shown in Ref. [START_REF] Wolf | Entangling power of passive optical elements[END_REF], the minimum value for ν-(i.e. the maximal entanglement) attainable by passive transformations is given by ν2 -= λ 1 λ 2 , where λ 1 and λ 2 are the two smallest eigenvalues of Γ. Therefore, the entanglement of symmetric states in standard form cannot be increased through energy preserving operations, like beam splitter and phase shifters. On the other hand, as it will be carefully discussed in great detail in the next section, the insertion of a birefringent plate in a type-II optical parametric oscillator results in states symmetric but not in standard form. In such a case the entanglement can be optimized by the action of a (passive) phase shifter.

A measure of entanglement which can be computed for general Gaussian states is provided by the negativity N , first introduced in Refs. [START_REF] Życzkowski | Volume of the set of separable states[END_REF]34], later thoroughly discussed and extended in Ref. [START_REF] Vidal | Computable measure of entanglement[END_REF] to CV systems. The negativity of a quantum state ̺ is defined as

N (̺) = ̺ 1 -1 2 , ( 23
)
where ̺ is the partially transposed density matrix and ô 1 = Tr|ô| stands for the trace norm of the hermitian operator ô. The quantity N (̺) is equal to | i λ i |, the modulus of the sum of the negative eigenvalues of ̺, quantifying the extent to which ̺ fails to be positive. Strictly related to N is the logarithmic negativity E N , defined as E N ≡ log 2 ̺ 1 , which constitutes an upper bound to the distillable entanglement of the quantum state ̺. Both the negativity and the logarithmic negativity have been proven to possess the crucial property of being monotone under LOCC (local operations and classical communications) [START_REF] Vidal | Computable measure of entanglement[END_REF][START_REF] Plenio | The logarithmic negativity is a full entanglement monotone under measuring LOCC and PPT-operations[END_REF].

For any two-mode Gaussian state ̺ the negativity is a simple decreasing function of ν-, which is thus itself an inverse quantifier of entanglement:

̺ 1 = 1 ν- ⇒ N (̺) = max 0, 1 -ν- 2ν - , (24)
E N (̺) = max [0, -log 2 ν-] . (25)
This expression quantifies the amount by which PPT inequality ( 20) is violated.

The symplectic eigenvalue ν-thus completely qualifies and quantifies the quantum entanglement of a Gaussian state Γ.

For symmetric Gaussian states one can also compute [START_REF] Giedke | Entanglement of formation for symmetric Gaussian states[END_REF] the entanglement of formation E F [START_REF] Bennett | Mixed-state entanglement and quantum error correction[END_REF]. We recall that the entanglement of formation E F of a quantum state ̺ is defined as

E F (̺) = min {pi,|ψi } i p i E(|ψ i ) , (26)
where the minimum is taken over all the pure states realizations of ̺:

̺ = i p i |ψ i ψ i | .
The quantity E F satisfies all the requirements of a proper entanglement measure.

The asymptotic regularization of the entanglement of formation is equal to the entanglement cost E C (̺), defined as the minimum number of singlets (maximally entangled antisymmetric two-qubit states) which is needed to prepare the state ̺ through local operations and classical communication [START_REF] Hayden | The asymptotic entanglement cost of preparing a quantum state[END_REF]. In formulae:

E C (̺) = lim n→∞ E F (̺ ⊗n ) n . ( 27
)
The optimal convex decomposition of Eq. ( 26) can be found for symmetric states and turns out to be Gaussian, allowing for the determination of the entanglement of formation E F :

E F = max [0, h(ν -)] , (28)
with

h(x) = (1 + x) 2 4x log 2 (1 + x) 2 4x - (1 -x) 2 4x log 2 (1 -x) 2 4x .
Such a quantity is, again, a decreasing function of ν-, thus providing a quantification of the entanglement of symmetric states equivalent to the one provided by the logarithmic negativity E N .

In the nonsymmetric case, an important result is that for any entangled two-mode Gaussian state ̺, the symplectic eigenvalue ν-(and, consequently, the logarithmic negativity) can be estimated with remarkable accuracy by only determining the global purity Tr ̺ 2 of the state, and the two local purities Tr ̺ 2 1,2 of each of the two reduced single-mode states ̺ i = Tr j ̺. For the aforementioned class of nonsymmetric thermal squeezed states, the estimate becomes actually an exact quantification, since the logarithmic negativity of these states is a function of the three purities only [START_REF] Adesso | Determination of continuous variable entanglement by purity measurements[END_REF]. These states are indeed the maximally entangled two-mode Gaussian states at fixed global and local purities, and thus they are the states one would like to produce and exploit in any continuous-variable quantum information processing. On the other hand, the symmetric instance, which carries the highest possible entanglement [START_REF] Adesso | Extremal entanglement and mixedness in continuous variable systems[END_REF] (and which is the experimental product of the present paper), is the one that enables continuousvariable teleportation of an unknown coherent state [START_REF] Furusawa | Unconditional quantum teleportation[END_REF][START_REF] Bowen | Experimental investigation of continuous-variable quantum teleportation[END_REF][START_REF] Takei | High-Fidelity Teleportation beyond the No-Cloning Limit and Entanglement Swapping for Continuous Variables[END_REF] with a fidelity arbitrarily close to 1 even in the presence of noise (mixedness), provided that the state is squeezed enough (ideally, a unit fidelity requires infinite squeezing). Actually, the fidelity of such an experiment, if the squeezed thermal states employed as shared resource are optimally produced, turns out to be itself a measure of entanglement and provides a direct, operative quantification of the entanglement of formation present in the resource [START_REF] Adesso | Equivalence between entanglement and the optimal fidelity of continuous variable teleportation[END_REF].

Effects of mode coupling on the entanglement generation

As mentioned in the introduction, entanglement is very often produced by mixing on a beam splitter two squeezed modes. In the general case, the squeezed quadratures have an arbitrary phase difference. We denote θ + π/2 the phase difference between the two squeezed quadratures. The CM of the squeezed modes is then

Γ A+A-=     a 0 0 0 0 1/a 0 0 0 0 b c 0 0 c b ′     , (29)
while the CM of the two modes after the beam-splitter is

Γ A1A2 = R T π/4 .Γ A+A-.R π/4 =     n 1 k ′ k k ′ k ′ n 2 k ′ -k k k ′ n 1 k ′ k ′ -k k ′ n 2     (30) where b = cos 2 θ a + a sin 2 θ , b ′ = a cos 2 θ + sin 2 θ a , c = a - 1 a sin θ cos θ , n 1 = cos 2 θ + a 2 (sin 2 θ + 1) 2a , n 2 = a 2 cos 2 θ + sin 2 θ + 1 2a , k = 1 -a 2 2a cos 2 θ , k ′ = a 2 -1 2a sin θ cos θ .
Let us first note that expression 30 can be brought back to the expression given on fig. 5 of [START_REF] Laurat | Effects of mode coupling on the generation on Einstein-Podolsky-Rosen entanglement in a type-II OPO below threshold[END_REF] via local unitary operations which do not modify the entanglement.

The CM of the squeezed (A ± ) modes gives a good insight into the properties of the two-mode state. One can see that the intermodal blocks are zero, meaning that the two modes are uncorrelated. Consequently, they are the two most squeezed modes of the system (no further passive operation can select more squeezed quadratures). But one can also note that the two diagonal blocks are not diagonal simultaneously. This corresponds to the tilt angle of the squeezed quadrature. In order to maximize the entanglement, the two squeezed quadratures have to be made orthogonal, which can be done by a phase-shift of one mode relative to the other.

It is easy in fact to compute the logarithmic negativity quantifying entanglement between the entangled modes A 1 and A 2 , when the two squeezed quadratures are rotated of π/2 + θ. One has

E N (Γ A1A2 ) = -(1/2) log ν2 , with ν2 = 1 4a 2 2 a 4 + 1 cos 2 (θ) + 4a 2 sin 2 (θ) - √ 2 a 2 -1 cos 2 (θ) a 4 + 6a 2 + (a 2 -1) 2 cos(2θ) + 1 . ( 31
)
The symplectic eigenvalue ν is obviously a periodic function of θ, and it is globally minimized for θ = k π, with k ∈ . The entanglement, in other words, is maximized for orthogonal modes in phase space, as already predicted in Ref. [START_REF] Wolf | Entangling power of passive optical elements[END_REF]. Notice that this results holds for general nonsymmetric states, i.e. also in the case when the two modes A 1 and A 2 possess different individual squeezings. For symmetric states, the logarithmic negativity is depicted as a function of the single-mode squeezing a and the tilt angle θ in figure 1.

In the experiment we will discuss below, the entanglement is produced by a single device, a type-II OPO operated below threshold. When no coupling is present in the optical cavity, the entangled modes are along the neutral axis of the crystal while the squeezed modes corresponds to the ±45 • linear polarization basis. However, we have shown theoretically and experimentally in [START_REF] Laurat | Effects of mode coupling on the generation on Einstein-Podolsky-Rosen entanglement in a type-II OPO below threshold[END_REF] that a coupling can be added via a birefringent plate which modifies the quantum properties of this device: the most squeezed quadratures are non-orthogonal with an angle depending on the plate angle.

When the plate angle increases, the squeezed (A -) quadrature rotates of a tilt angle θ and the correlations are degraded. The evolution is depicted on figure 2 through the noise ellipse of the superposition modes.

A + A - q

A + A - q A + A - Figure 2.
Fresnel representation of the noise ellipse of the ±45 • rotated modes when the coupling is increased. The noise ellipse of the -45 • mode rotates and the noise reduction is degraded when the coupling increases while the +45 • rotated mode is not affected.

Eqn. 31 shows that when coupling is present, it is necessary to perform an operation on the two modes in order to optimize the available entanglement. Before

IN Trig H 3 l/2 KTP PZT PZT - Q 0 - Ramp IN Trig SA 1 H 2
Frequency-Doubled Nd:YAG Laser developing experimental measures of entanglement and optimization of the available resource in our system, let us detail our experimental setup.

Ä 12 MHz Servo Servo Isolator { Q 1 H 1 Q 2 Q 4 H 4 H 5 - IN Trig

Experimental setup and homodyne measurement

The experimental scheme is depicted on figure 3 and relies on a frequency-degenerate type-II OPO below threshold. The system is equivalent to the one of the seminal experiment by Ou et al. but a λ/4 birefringent plate has been inserted inside the optical cavity (Q 0 ). When this plate is rotated, it results in a linear coupling between the signal and idler modes which induces above threshold a phase locking effect at exact frequency degeneracy [START_REF] Mason | Observation of two distinct phase states in a self-phase-locked type II phase-matched optical parametric oscillator[END_REF][START_REF] Longchambon | Non-linear and quantum optics of a type II OPO containing a birefringent element, 1. classical operation[END_REF]. This triply-resonant OPO is pumped below threshold with a continuous frequency-doubled Nd:YAG laser. The input flat mirror is directly coated on one face of the 10mm-long KTP crystal. The reflectivities for the input coupler are 95% for the pump (532nm) and almost 100% for the signal and idler beams (1064nm). The output coupler (R=38mm) is highly reflecting for the pump and its transmission is 5% for the infrared. At exact triple resonance, the oscillation threshold is less than 20 mW. The OPO is actively locked on the pump resonance by the Pound-Drever-Hall technique. The triple resonance is reached by adjustment of both the crystal temperature and the frequency of the pump laser. Under these conditions, the OPO can operate stably during more than one hour without modehopping. The birefringent plate Q 0 inserted inside the cavity is exactly λ/4 at 1064 nm and almost λ at the 532 nm pump wavelength. Very small rotations of this plate around the cavity axis can be operated thanks to a piezoelectric actuator.

Measurements of the quantum properties of arbitrary quadratures of light mode are generally made using homodyne detection [START_REF] Yuen | Noise in homodyne and heterodyne detection[END_REF]. When an intense local oscillator is used, one obtains a photocurrent which is proportional to the quantum noise of the light in a quadrature defined by the phase-shift between the local oscillator and the beam measured. This photocurrent can be either sent to a spectrum analyzer which calculates the noise power spectrum or numerized for further treatments like tomographic measurements of the Wigner function [START_REF] Vogel | Determination of quasiprobability distributions in terms of probability distributions for the rotated quadrature phase[END_REF] or selection [START_REF] Laurat | Conditional preparation of a quantum state in the continuous variable regime: generation of a sub-Poissonian state from twin beams[END_REF]. As mentioned above, one can also characterize the entanglement by looking at linear combinations of the optical modes as opposed to linear combinations of the photocurrents [START_REF] Duan | Inseparability criterion for continuous variable systems[END_REF][START_REF] Simon | Peres-Horodecki separability criterion for continuous variable systems[END_REF]. The two modes which form the entangled state must be transformed via the beam splitter operation (R π/4 ), that is they are mixed on a 50/50 beam splitter or a polarizing beam-splitter preceded by λ/2 waveplate (H 2 ) into two modes which will be both squeezed if the original state is entangled.

Homodyne detection allows for a simple and direct measurement of the 2 × 2 diagonal blocks of the 4×4 CM. In order to measure the 2×2 off-diagonal blocks, linear combinations of the fields can be used as we will show below. In order to characterize simultaneously two modes a single phase reference is needed. To implement this, we have built a simultaneous double homodyne detection (Fig. 3, in box). The difference photocurrents are sent into two spectrum analyzers triggered by the same signal (SA 1 and SA 3 ). Two birefringent plates (Q 4 , H 3 ) inserted in the local oscillator path are rotated in order to compensate residual birefringence. A λ/4 (Q 3 ) plate can be added on the beam exiting the OPO in order to transform the in-phase detections into inquadrature ones (making the transformation (x + , p+ , x-, p-) → (x + , p+ , p+ , x+ )). In such a configuration, two states of light with squeezing on orthogonal quadratures give in phase squeezing curves on the spectrum analyzers. This has two goals: firstly, it simplifies the measurements of the phase-shift between the two homodyne detections and secondly it is necessary for the measurement of the off-diagonal blocks of the CM as we will show now.

Let us describe precisely the procedure used to extract the values of the CM from the homodyne detection signals. These signals consist in an arbitrary pair of spectrum analyzer traces which are represented in figure 4. The horizontal axis is the local oscillator (LO) phase which is scanned as a function of time via a piezoelectric transducer while the vertical axis gives the noise power relative to the shot noise expressed in decibels (dB).

We make no assumption on the form of the CM which is written in the general case

Γ = Γ + C ± C T ± Γ - =     a b c d b e f g c f h i d g i j    
When the LO phase is chosen so that zero corresponds to the long axis of the noise ellipse of the first mode, the CM is written in the form

Γ =     a ′ 0 c ′ d ′ 0 e ′ f ′ g ′ c ′ f ′ h ′ i ′ d ′ g ′ i ′ j ′     = Γ ′ + C ′ ± C ′T ± Γ ′ -
a ′ and c ′ correspond respectively to the maximum and minimum noise levels measured in a linear scale on the spectrum analyzer for the first mode, which we will choose arbitrarily to be A + . One can also easily determine h ′ , i ′ and j ′ from the spectrum analyzer traces for A -: when the LO phase is chosen so that zero corresponds to the long axis of the noise ellipse of the second mode, the CM is written in the form

Γ ′′ -=
h ′′ 0 0 j ′′ and Γ ′ -can be easily deduced from Γ ′′ -by applying a rotation. The angle of this rotation is given by the phase shift ϕ between the two traces (see Fig. 4). This operation is performed numerically. We have now measured both diagonal blocks.

In order to measure the non-diagonal blocks, one records on an additional spectrum analyzer a third signal, the difference between the two homodyne detection signals (these signals being themselves the difference between their respective photodiodes photocurrents, see Fig. 3). Let us consider the case where the waveplate Q 3 is not present, for a given LO phase ψ 1 , the homodyne detections will give photocurrents which are proportional to the amplitude noise for the A + beam, x+ and to the phase noise for the A -beam, p -. The signal recorded on spectrum analyzer SA 2 is, in this case, proportional ŝ = x+p-whose variance is î2 = x2 + + p2 --2 x+ p-= a ′ + j ′ -2d ′ . a ′ and j ′ being already known, it is easy to extract d ′ from this measurement. For a LO phase ψ 1 + π/2, one will get using a similar procedure e ′ , h ′ and f ′ . Let us know add the wave plate Q 3 . For a LO phase ψ 1 , one will get a ′ , h ′ and c ′ and for ψ 1 + π/2 e ′ , j ′ and g ′ thus completing the measurement of the CM.

Experimental measures of entanglement by the negativity

As described previously, the entanglement is deduced experimentally from the CM. As all experimental measurements, the measurement of the CM is subject to noise. It is thus critical to evaluate the influence of this noise on the entanglement. A quantitative analysis, relating the errors on the measured CM entries (in the A ± basis) to the resulting error in the determination of the logarithmic negativity (the latter quantifying entanglement between the corresponding A 1 and A 2 modes) has 2 × 2 blocks of the measured CM Γ in the A ± basis, given by Eq. [START_REF] Serafini | Unitarily localizable entanglement of Gaussian states[END_REF]. In plot (a): the solid red curve refers to equal errors (of value δΓ) on the eight entries of the diagonal blocks, the dotted blue curve refers to equal errors on the four diagonal entries of the diagonal blocks while the dashed green curve refers to equal errors on the off-diagonal entries of the diagonal blocks (non standard form entries). At δΓ 0.16 some of the considered states get unphysical. In plot (b): the solid red curve refers to equal errors on the four entries of the off-diagonal block, the dotted blue curve refers to equal errors on the two off-diagonal entries of the off-diagonal block (non standard form entries), while the dashed green curve refers to equal errors on the diagonal entries of the off-diagonal block (standard form entries). At δΓ 0.19 some of the considered states get unphysical.

been carried out and is summarized in Fig. 5 in absence of mode coupling. In general, the determination of the logarithmic negativity is much more sensitive to the errors on the diagonal 2 × 2 blocks α and β (referring to the reduced states of each mode, see Eq. ( 12)) of the CM Γ than on the off-diagonal ones (γ, and its transpose γ T , whose expectations are assumed to be null). Let us remark that the relative stability of the logarithmic negativity with respect to the uncertainties on the off-diagonal block adds to the reliability of our experimental estimates of the entanglement. Notice also that, concerning the diagonal blocks, the errors on diagonal (standard form) entries turn out to affect the precision of the logarithmic negativity more than the errors on off-diagonal (non standard form) entries. This behavior is reversed in the off-diagonal block, for which errors on the off-diagonal (non standard form) entries affect the uncertainty on the entanglement more than errors on the diagonal (standard form) entries. Experimentally, we have measured the noise on the CM elements to be at best on the order of a few percents of the measured values for the diagonal blocks, corresponding to a fraction of a dB [START_REF]The noise expressed in decibels (dB) is obtained from the covariance matrix elements via[END_REF]. This is the case for the diagonal blocks which are well-known since they are directly related to the noise measurements of A + and A -. The situation is less favorable for the off-diagonal blocks which require the two homodyne detections to be very precisely in phase since signals from these two detections are used. As a consequence, these blocks show a large experimental noise. As shown on Fig. 5 (b), this may lead in some cases to unphysical CM, yielding for instance a negative determinant and complex values for the logarithmic negativity. Typical values are 0.07 ± 0.05 dB or 1 ± 1 dB for the diagonal elements of the offdiagonal blocks. A phase difference between the two homodyne detections of a few degrees can easily explain such values, especially when the squeezing is large as in our case. In the following, we will set these terms to zero in agreement with the form of the CM of Eq. 29.

Let us first give an example of entanglement determination from measurements of CM elements, in the absence of mode coupling. Without the plate, the squeezing of the two superposition modes is expected on orthogonal quadratures (the ideal CM is then in the form 29 with θ = 0). Spectrum analyzer traces while scanning the local oscillator phase are shown on figure 6: the rotated modes are squeezed on orthogonal quadratures.

The state is produced directly in the standard form and the CM in the A ± basis can be extracted from this measurement:

Γ(ρ = 0) =     0.33 0 (0) (0) 0 7.94 (0) (0) (0) (0) 7.94 0 (0) (0) 0 0.33     (32)
The resulting smallest symplectic eigenvalue is the geometric average of the two minimal diagonal elements : ν-= 0.33, yielding a logarithmic negativity E N = log 2 (ν -) = 1.60 between the modes A 1 and A 2 .

Experimental non standard form and optimization by linear optics

As discussed previously, when the plate angle is increased, the state produced is not anymore in the standard form but rather similar to eqn. 29.

In this instance, the logarithmic negativity between A 1 and A 2 is much lower than the previous value: E N = 1.13.

As stated above, entanglement can be increased via passive operations performed simultaneously on the two modes. Such operations include phase-shifts and beam splitters, which can be readily performed on co-propagating, orthogonally polarized beams [START_REF] Simon | Quantum noise matrix for multimode systems: U(n) invariance, squeezing and normal forms[END_REF]. The minimal combination of waveplates can be shown to consist in three waveplates : two λ/4 waveplates (Q 1 and Q 2 on Fig. 3) and one λ/2 waveplate (H 1 ). When using a combination of these three plates, the state can be put back into standard form which will maximize the entanglement. This operation consists in a phase-shift of the rotated modes. Figure 7 gives the normalized noise variances before and after this operation. The CM is changed to: Γ(ρ = 0.

giving a logarithmic negativity E N = 1.32 between A 1 and A 2 , larger than the value before the operation. It is also the maximal value than can be obtained considering the available entanglement. This value is smaller than the one obtained when the cavity waveplate (Q 0 ) is not tilted since tilting this waveplate reduces the correlations Let us remark again that this transformation is non-local in the sense of the EPR argument: it has to be performed before spatially separating the entangled modes for a quantum communication protocol for instance.

Conclusion

We have described the powerful tools underlying the description of continuous variables systems in quantum optics. These tools allow for a nice pictorial view of two-mode Gaussian entangled states. Specifically, we have illustrated their properties through the description and manipulation of the entanglement produced experimentally with an original device, a type-II optical parametric oscillator containing a birefringent plate. We have demonstrated the capabilities of this system through entanglement measurements and manipulation of entangled states, showing in particular how entanglement can be maximized using purely passive operations.

We have also studied quantitatively the influence on the entanglement of the noise affecting the measurement of the elements of the CM. We thus showed that the most significant covariances (exhibiting the highest stability against noise) for an accurate entanglement quantification are the diagonal terms of the diagonal singlemode blocks, and the off-diagonal terms of the intermodal off-diagonal block, the latter being the most difficult to measure with high precision. Alternative methods have been devised to tackle this problem [START_REF] Fiurasek | How to measure squeezing and entanglement of Gaussian states without homodyning[END_REF][START_REF] Adesso | Determination of continuous variable entanglement by purity measurements[END_REF] based on direct measurements of global and local invariants of the CM. Such techniques have been implemented in the case of pulsed beams [START_REF] Wenger | Pulsed squeezed vacuum characterization without homodyning[END_REF] but no experiment to date has been performed for continuous-wave beams.

2 ΘFigure 1 .

 21 Figure 1. Logarithmic negativity as a function of the single-mode squeezing a and the tilt angle θ between the two non-orthogonal quadratures in presence of mode coupling.

Figure 3 .

 3 Figure 3. Experimental setup. A continuous-wave frequency-doubled Nd:YAG laser pumps below threshold a type II OPO with a λ/4 plate inserted inside the cavity (Q 0 ). The generated two-mode vacuum state is characterized by two simultaneous homodyne detections. The infrared output of the laser is used as local oscillator after filtering by a high-finesse cavity. SA 1,2,3 : spectrum analyzers. Q 1,2,3,4 and H 1,...,5 : respectively quarter and half waveplates. PD Lock: FND-100 photodiode for locking of the OPO. PD Split: split two-element InGaAs photodiode for tilt-locking of the filtering cavity.

Figure 4 .

 4 Figure 4. Spectrum analyzer traces as a function of the local oscillator phase.

Figure 5 .

 5 Figure 5. Error δE N on the logarithmic negativity between modes A 1 and A 2 , as a function of the error δΓ on the entries of the diagonal (a) and off-diagonal (b) 2 × 2 blocks of the measured CM Γ in the A ± basis, given by Eq.[START_REF] Serafini | Unitarily localizable entanglement of Gaussian states[END_REF]. In plot (a): the solid red curve refers to equal errors (of value δΓ) on the eight entries of the diagonal blocks, the dotted blue curve refers to equal errors on the four diagonal entries of the diagonal blocks while the dashed green curve refers to equal errors on the off-diagonal entries of the diagonal blocks (non standard form entries). At δΓ 0.16 some of the considered states get unphysical. In plot (b): the solid red curve refers to equal errors on the four entries of the off-diagonal block, the dotted blue curve refers to equal errors on the two off-diagonal entries of the off-diagonal block (non standard form entries), while the dashed green curve refers to equal errors on the diagonal entries of the off-diagonal block (standard form entries). At δΓ 0.19 some of the considered states get unphysical.

Figure 6 .

 6 Figure 6. Normalized noise variances at 3.5 MHz of the ±45 • modes while scanning the local oscillator phase. The first plot corresponds to in-phase homodyne detections and the second one in-quadrature. Squeezing is well observed on orthogonal quadratures. (RBW 100 kHz, VBW 1 kHz)

Figure 7

 7 gives the normalized noise variances at 3.5 MHz of the rotated modes while scanning the local oscillator phase for an angle of the plate ρ = 0.3 • . The first plot shows that the squeezing is not obtained on orthogonal quadratures. The CM takes the following form:

Figure 7 .

 7 Figure 7. Normalized noise variances at 3.5 MHz of the rotated modes while scanning the local oscillator phase for an angle of the plate of 0.3 • , before and after the non-local operation. The homodyne detections are in-quadrature. After this operation, squeezing is observed on orthogonal quadratures.

Acknowledgement Laboratoire Kastler-Brossel, of the Ecole Normale Supérieure and the Université Pierre et Marie Curie, is associated with the Centre National de la Recherche Scientifique (UMR 8552). Laboratoire Matériaux et Phénomènes Quantiques, of Université Denis Diderot, is associated with the Centre National de la Recherche Scientifique (UMR 7162). This work has been supported by the European Commission project QUICOV (IST-1999-13071) and ACI Photonique (Ministère de la Recherche et de la Technologie). JAOH acknowledges financial support from CAPES/COFECUB. AS acknowledges financial support from EPSRC through the grant QIP-IRC. GA and FI acknowledge INFM, INFN, and MIUR under national project PRIN-COFIN 2002 for financial support.


Mathilde Bouvel
email: mbouvel@liafa.jussieu.fr

Dominique Rossin
email: rossin@liafa.jussieu.fr

Longest Common Pattern between two Permutations

In this paper, we give a polynomial (O(n 8 )) algorithm for finding a longest common pattern between two permutations of size n given that one is separable. We also give an algorithm for general permutations whose complexity depends on the length of the longest simple permutation involved in one of our permutations.

Introduction and basic concepts

The study of patterns in permutations has blossomed these last years: from a combinatorial point of view with the recent proof of the Stanley-Wilf conjecture by Marcus and Tardös, and from an algorithmic one with the development of algorithms for pattern and general pattern involvement. Although the general pattern involvement problem is N P -hard, some polynomial solutions exist for special kinds of patterns like the separable ones [START_REF] Bose | Pattern matching for permutations[END_REF][START_REF] Ibarra | Finding pattern matchings for permutations[END_REF]. In this article we study the problem of finding a longest common pattern between two permutations σ 1 and σ 2 i.e. a permutation σ which is involved in both permutations σ 1 and σ 2 . This is a generalization of the pattern involvement problem since finding if the longest pattern between σ 1 and σ 2 is equal to σ 1 is equivalent to the pattern involvement problem.

First, we give a polynomial algorithm based on the work of [START_REF] Bose | Pattern matching for permutations[END_REF] for finding the longest common pattern if one permutation is separable. Then we generalize this algorithm for general permutations. The complexity of our algorithm is highly based on the length of the longest simple permutation [Bri] involved in our permutations.

Permutations

A permutation σ of an interval I of N is a bijective map from I to itself. We denote by σ i the image of i by σ. The permutation σ could either be seen as a function or a word σ i σ i+1 . . . σ j , where I = {h : i ≤ h ≤ j}. For example the permutation σ = 1 4 2 5 6 3 is the bijective function such that σ(1) = 1, σ(2) = 4, σ(3) = 2, σ(4) = 5 . . .. In the following, when we consider permutations without giving explicitly the interval I, we mean that I = {1, . . . , n} for some n.

Definition 1.1. A permutation π = π 1 . . . π k is called a pattern of the permutation σ = σ i+1 . . . σ i+n of I = {h : i + 1 ≤ h ≤ i + n}, with k ≤ n, if and only if there exist integers i + 1 ≤ i 1 < i 2 < . . . < i k ≤ i + n such that σ i ℓ < σ im whenever π ℓ < π m . We will also say that π is involved in σ or that σ contains π. The subsequence σ i1 σ i2 . . . σ i k is called an occurrence of π in σ.

A permutation σ that does not contain π as a pattern is said to avoid π.

Example 1.2. For example σ = 1 4 2 5 6 3 contains the pattern 1 3 4 2 and 1 5 6 3, 1 4 6 3, 2 5 6 3 and 1 4 5 3 are the occurrences of this pattern in σ. But σ avoids the pattern 3 2 1 as no subsequence of length 3 of σ is isomorphic to 3 2 1 i.e. is decreasing.

A number of enumerative results has been proved on classes of pattern avoiding permutations for patterns of length 3, 4 and multiple patterns. More recently results about the algebricity of the generating function of general classes of permutations have been given [START_REF] Brignall | Simple permutations and algebraic generating functions[END_REF].

Another field of study of these permutations is from the point of view of pattern involvement. The problem of deciding if a permutation π is a pattern of a permutation σ is N P -hard but this problem is proved to be polynomial if the pattern is separable [START_REF] Bose | Pattern matching for permutations[END_REF][START_REF] Ibarra | Finding pattern matchings for permutations[END_REF].

Definition 1.3. A permutation σ of size n is called separable if it avoids the patterns 3 1 4 2 and 2 4 1 3 or equivalently if it has a binary separating tree.

Definition 1.4. A binary separating tree is a binary ordered tree with n leaves such that each internal vertex is labeled by + or -.

For each such tree, there is a unique way [START_REF] Bose | Pattern matching for permutations[END_REF] to decorate its leaves (considering them from left to right) by σ 1 , σ 2 , . . . , σ n such that:

1. σ 1 . . . σ n is a permutation of {1 . . . n}.

Each node (internal or leaf

) is decorated by a permutation of an interval.

3. Each internal node labeled with + (resp. -) is decorated by a permutation of {i . . . j} such that its left child is decorated by a permutation of {i . . . h-1} (resp. {h . . . j}) and its right child is decorated by a permutation of {h . . . j} (resp. {i . . . h -1}), for some h ∈ {i + 1, . . . , j}.

It is easy to prove [START_REF] Bose | Pattern matching for permutations[END_REF] that to each separable permutation one can associate a binary separating tree (see Figure 1). Note that this tree is not uniquely defined as shown in figure 1. However, one can associate a unique tree to each separable permutation by taking arbitrary ordered trees instead of binary ordered trees. These trees are the contraction of the binary ones by contracting every edge between two nodes with the same label + or -. In these trees, the signs of the internal nodes are thus ranked: if the root of the tree has label + (resp. -), then every node at odd depth has label + (resp. -) and every node at even depth label -(resp. +), so that all labels are determined from the label of the root.

Definition 1.5. The unique contracted separating tree associated to a separable permutation σ is obtained from any binary separating tree of σ by contracting every edge between nodes of the same sign.

Modular decomposition of graphs ; Common interval decomposition of permutations

The contracted separating trees we introduced in Definition 1.5 also appear in graph theory. Namely, those trees are a special case of common interval decomposition trees (on which we however need to add a labeling). The common interval decomposition trees are an equivalent for permutations of the modular decomposition trees for graphs [START_REF] De Montgolfier | Décomposition modulaire des graphes : théories extensions et algorithmes[END_REF][START_REF] Bui Xuan | Revisiting t. uno and m. yagiura's algorithm[END_REF]. Before we come to the pattern matching problem in permutations, we need to introduce the common interval decomposition trees, the labeled decomposition trees, and finally the expanded decomposition trees that are the key structure we use in our algorithms.

The common interval decomposition of a permutation σ of size n (common meaning common with Id n ) is defined as follows. First consider all the common intervals of σ that is to say all the subsequences σ j σ j+1 . . . σ k of consecutive entries of σ such that {σ j , σ j+1 , . . . , σ k } is an interval of N. Among the common intervals, the strong common intervals (or strong intervals for short) are those that do not overlap any other common interval1 . Figure 2 illustrates the notion of common intervals and strong intervals. The inclusion ordering yields a tree-like ordering on the set of strong intervals. This ordering is represented by a tree whose leaves are σ 1 , σ 2 , . . ., σ n from left to right in this order, whose root is σ, and such that each internal node is the union of its children.

Note that there are two different types of internal nodes in the tree of strong intervals. For some nodes, say V , with k children V 1 , . . . , V k from left to right, there do not exist (i, j) = (1, k) such that 1 ≤ i < j ≤ k and the union of V i , V i+1 , . . . , V j is an interval. These nodes are called prime nodes and are of type P .

The other nodes V are such that every union of consecutive children form an interval. Those node are called linear nodes and are of type L.

The tree along with the types P and L of the internal nodes (see Figure 3) is called the common interval decomposition tree of σ. In this tree, the order of the children of a node depends on σ so that we have an ordered tree, unlike the modular decomposition trees for graphs.

Note that nodes of arity 2 verify both linear and prime definitions. We choose to consider them of type L. This choice will be explained later. For further explanation on these trees, like the proof that a node is either linear or prime, refer to [BXHP]. In [BCdMR, BXHP], the authors show how to compute these trees in linear time.

For our algorithmic use of common interval decomposition tree, the types of the internal nodes are not sufficient and we need to label the internal nodes of the common interval decomposition tree of a permutation. 1 It is easy to see that for each node V of type L, the intervals of values corresponding to the children of V are ordered either by increasing order or by decreasing order, when considering these children from left to right. The linear nodes are subsequently labeled + or -respectively.

The labeling of a prime node with d children consists of a permutation σ of size d that does not have any common intervals except {1}, . . . , {d} and σ. Such permutations are known as simple permutations [AA05, AAK03, Bri, BHV06]. σ represents the ordering of the children V 1 , . . . , V d of V between them with respect to the values in the common intervals corresponding to the V i 's. Namely, σ i < σ j if and only if the common interval corresponding to V i contains values that are smaller than those contained in the common interval corresponding to V j . For example, the simple permutation labeling the root of the tree on figure 3 is 3 1 4 2.

This common interval decomposition tree along with labels +,-and simple permutations can be computed easily in O(n 2 ln n) which is sufficient for our purpose. Remind that the common interval decomposition tree can be computed in linear time. Thus it remains to add a label on each internal node V . This can be done by sorting the intervals corresponding to the children of V .

When a common interval decomposition tree is labeled, the common intervals corresponding to the nodes can be deduced from this labeling, like in the case of binary separating trees (see definition 1.4). So that the common intervals can be seen as a decoration on the nodes of the common interval decomposition tree with labels on internals nodes.

Definition 1.6. The labeled decomposition tree of a permutation σ is the common interval decomposition tree of σ, where we add the labeling on internal nodes described above, and where we forget the decoration.

Figure 4 gives the labeled decomposition tree of σ = 5 1 10 9 6 7 8 11 2 4 3.

We have the following nice characterization of separable permutations in terms of labeled decomposition trees: Proposition 1.1. The separable permutations are exactly those having a common interval decomposition tree with no prime nodes. For any separable permutation, its contracted separating tree and its labeled decomposition tree are equal.

The proof of this proposition is straightforward from the definition of separable permutations. Note that it is important to give the type L to binary internal nodes for stating that proposition.

Proposition 1.1 states that labeled decomposition trees are a generalization to all permutations of the contracted separating tree defined only for separable permutations. For binary separating trees, we have an analogous generalization, called expanded decomposition trees. In our algorithms, it is easier to work on binary nodes so that the natural representation used later is the expanded decomposition tree.

To transform a labeled decomposition tree into an expanded decomposition tree, take each linear node V with children V 1 , . . . , V k and note that we can represent it by (. . . ((V 1 , V 2 ), V 3 ), . . .), V k ) as shown in figure 4. Each positive (resp. negative) internal node of arity k gives k -1 positive (resp. negative) binary internal nodes.

A consequence of proposition 1.1 is:

Proposition 1.2. The expanded decomposition tree of a separable permutation is one of its binary separating trees.

Longest common pattern between two permutations, including one separable

In this section, we describe a polynomial time algorithm for finding a longest common pattern between two permutations σ and τ provided that σ is separable. This algorithm uses the same technique as the one of Bose, Buss and Lubiw [START_REF] Bose | Pattern matching for permutations[END_REF] for finding an occurrence of a separable pattern in a general permutation. Namely, it computes a binary separating tree T σ of σ and uses it as a guide in the search of the longest common pattern with τ . To avoid a complexity blow-up, the key point is to use dynamic programming, the initial problem being decomposed into sub-problems according to the structure of T σ . First recall that it can be decided in linear time whether a permutation is separable or not. If it is, its binary separating tree can also be computed in linear time, as described in [START_REF] Bose | Pattern matching for permutations[END_REF]. This result is actually a special case of a more general one, proved in [BCdMR, BXHP], and stating that the common interval decomposition tree of any permutation can be computed in linear time.

Instead of two permutations σ and τ , our algorithm takes as an input a binary separating tree T σ for a separable permutation σ of size k, and a permutation τ in the usual representation τ 1 τ 2 . . . τ n . Notice that a binary separating tree for σ has O(k) nodes.

More precisely, the algorithm fills in the array

M = {M (V, i, j, a, b) : V a node of T σ , 1 ≤ i ≤ j ≤ n, 1 ≤ a ≤ b ≤ n}.
For any node V in T σ , let us denote by σ(V ) the subpermutation of σ corresponding to the subtree of T σ rooted at V . With the notations of section 1, σ(V ) is the permutation decorating the node V . The cell M (V, i, j, a, b) of the array M contains a longest common pattern π between σ(V ) on one hand, and the subpermutation τ i . . . τ j of τ on the other hand, with the additional restriction that the occurrence of π in τ i . . . τ j uses only entries of τ whose values are between a and b. The empty pattern, of size 0, will be denoted ǫ.

Example 2.1. If the node V represents the pattern 2 1 (i.e. σ(V ) = (i + 1) i), and τ = 6 4 2 5 3 1 then we have for instance M (V, 2, 4, 3, 5) = 1, M (V, 2, 5, 3, 4) = 2 1 and M (V, 4, 5, 1, 2) = ǫ.

The algorithm works as follows. To start the computation, we fill in the subarrays M (V, , , , ) for all the leaves V of T σ . Then, we compute M (V, , , , ) for any internal node V using the subarrays M (V L , , , , ) and M (V R , , , , ) corresponding to the left child (V L ) and the right child (V R ) of V . In order to combine the patterns found in M (V L , , , , ) and M (V R , , , , ) with the intention of filling M (V, , , , ), we need a definition of pattern concatenation. Definition 2.2. Consider two patterns π and π ′ of respective lengths k and k ′ . The positive and negative concatenations of π and π ′ are defined respectively by: It is clear from definition 2.2 that the positive (resp. negative) concatenation of two patterns produces again a pattern.

π ⊕ π ′ = π 1 • • • π k (π ′ 1 + k) • • • (π ′ k ′ + k) and π ⊖ π ′ = (π 1 + k ′ ) • • • (π k + k ′ )π ′ 1 • • • π ′ k ′ Example 2.3.
Algorithm 1 Longest common pattern between two permutations, first one separable 1: Input: A binary separating tree T σ of a separable permutation σ of size k and a permutation τ of size n 2: Create an array M : 3: for any node V of T σ and any integers i, j, a and b between 1 and n do 4:

M (V, i, j, a, b) ← ǫ 5: end for 6: Fill in M for the leaves of T σ : 7: for any leaf V of T σ do 8:

for any integers i, j, a and b between 1 and n, i ≤ j, a ≤ b do 9:

if there exists some h ∈ {i, i + 1, . . . , j} such that a ≤ τ h ≤ b then 10:

M (V, i, j, a, b) ← 1 11: end if 12:
end for 13: end for 14: Fill in the rest of M : 15: for any internal node V of T σ , considering the nodes in the postfix ordering do 16:

if V is a positive node then 17:

for any integers i, j, a and b between 1 and n, i ≤ j, a ≤ b do 18:

M (V, i, j, a, b) ← Longest({M (V L , i, h-1, a, c-1)⊕M (V R , h, j, c, b) : i ≤ h ≤ j + 1, a ≤ c ≤ + 1}) 19:
end for 20:

else 21:

/* V is a negative node /*

22:

for any integers i, j, a and b between 1 and n, i ≤ j, a ≤ b do 23:

M (V, i, j, a, b) ← Longest({M (V L , i, h-1, c, b)⊖M (V R , h, j, a, c-1) : i ≤ h ≤ j + 1, a ≤ c ≤ b + 1}) 24:
end for 25:

end if 26: end for

27: Output: M (root of T σ , 1, n, 1, n)
The detailed dynamic programming algorithm is given in Algorithm 1. In the previous, for any set S of patterns, Longest(S) returns a longest pattern in the set S. If we want Longest to be a uniquely defined application, we can choose Longest(S) as the smallest pattern with respect to the lexicographical ordering among the patterns in S having maximal length.

Proposition 2.1. Algorithm 1 is correct: it outputs a longest common pattern between the two permutations σ and τ given in input.

Proof. The proof is by induction. We show that the algorithm stores in M (V, i, j, a, b) a longest common pattern between σ(V ) and τ i . . . τ j whose occurrence in τ i . . . τ j uses only values between a and b.

If V is a leaf, the above statement is clearly true. If V is an internal node, with two children V L (its left child) and V R (its right child), then let i, j, a and b be integers such that 1 ≤ i ≤ j ≤ n and 1 ≤ a ≤ b ≤ n. We assume in the rest of the proof that V is a positive node, the case of a negative node being very similar. To begin with, it is easy to see that M (V, i, j, a, b) contains a common pattern between σ(V ) and τ i . . . τ j using only values between a and b in τ i . . . τ j . Indeed, by induction hypothesis, we infer that every pattern in the set

S = {M (V L , i, h -1, a, c -1) ⊕ M (V R , h, j, c, b) : i ≤ h ≤ j + 1, a ≤ c ≤ b + 1}
, and a fortiori Longest(S), is a common pattern between σ(V ) and τ i . . . τ j using only values between a and b in τ i . . . τ j .

Example 2.4. Consider the permutations σ = 1 4 2 3 6 5 7 8 and τ = 4 1 3 2 5 6 8 9 7. A separating tree T σ of σ is represented on the right. For example, we can choose V to be the right child of the root of T σ . Then we have σ(V ) = 4 2 3 6 5 7 8, σ(V L ) = 4 2 3 and σ(V R ) = 6 5 7 8. Now choose i = 2, j = 7, a = 2 and b = 8. We want to show that for any h ∈ {i, . . . j + 1} and c ∈ {a, . . By induction hypothesis, M (V L , i, h -1, a, c -1) contains a longest common pattern between σ(V L ) and τ i . . . τ h-1 using only values between a and c -1 in τ i . . . τ h-1 . Here, M (V L , i, h -1, a, c -1) = 2 1, an occurrence of 2 1 in τ i . . . τ h-1 = 1 3 2 using values between 2 and 3 being 1 3 2, and an occurrence of 2 1 in σ(V L ) = 4 2 3 being 4 2 3. Similarly, we have M (V R , h, j, c, b) = M (V R , 5, 7, 4, 8) = 1 2 3, as shown by the occurrences 5 6 8 in τ h . . . τ j = 5 6 8 and 6 5 7 8 in σ(V R ) = 6 5 7 8. The occurrence of M (V L , i, h -1, a, c -1) ⊕ M (V R , h, j, c, b) = 2 1 3 4 5 in τ i . . . τ j using values between a and b is thus obtained by considering simultaneously the two occurrences in τ enlightened before. Namely, an occurrence of 2 1 3 4 5 in τ i . . . τ j = 1 3 2 5 6 8 using values between 2 and 8 is 1 3 2 5 6 8. Notice also that the occurrence 4 2 3 6 5 7 8 of 2 1 3 4 5 in σ(V ) = 4 2 3 6 5 7 8 is again obtained considering simultaneously the occurrences of 2 1 and 1 2 3 in σ(V L ) and σ(V R ) respectively.

. b + 1}, M (V L , i, h -1, a, c -1) ⊕ M (V R , h, j, c
It remains to show that:

Lemma 2.5. Longest(S) is of maximal length among all the common patterns between σ(V ) and τ i . . . τ j using only values between a and b in τ i . . . τ j .

Proof. We denote by π a longest common pattern between σ(V ) and τ i . . . τ j , using only values between a and b in τ i . . . τ j . As shown on Figure 5, there exist integers h ∈ {i, . . . , j+1} and c ∈ {a, . . . , b+1} such that π is decomposed into π = π 1 ⊕ π 2 , with π 1 a common pattern between σ(V L ) and τ i . . . τ h-1 , using only values between a and c -1 in τ i . . . τ h-1 , and π 2 a common pattern between σ(V R ) and τ h . . . τ j , using only values between c and b in τ h . . . τ j .

+ node V σ(V L ) σ(V R ) σ(V ) = π 1 π 2 ⊕ τ h . . . τ j τ i . . . τ j = τ i . . . τ h-1 π = Figure 5: Proof of Lemma 2.5
Notice that in this decomposition π 1 or π 2 might be the empty pattern. It can be easily seen that if π 1 (respectively π 2 ) were not a pattern of maximal length for the given intervals of indices and values, then π would not be of maximal length either, contradicting the definition of π. Consequently, π 1 and π 2 are longest common patterns between σ(V L ) and σ(V R ) respectively, and τ in the given intervals of indices and values. By induction hypothesis,

|M (V L , i, h -1, a, c -1)| = |π 1 | and |M (V R , h, j, c, b)| = |π 2 |. The pattern stored in M (V, i, j, a, b) by the algorithm is of length at least |M (V L , i, h -1, a, c -1) ⊕ M (V R , h, j, c, b)| = |π 1 ⊕ π 2 | = |π|.
π being of maximal length by assumption, we conclude that M (V, i, j, a, b) is also of maximal length.

Finally, M (V, i, j, a, b) contains a longest common pattern between σ(V ) and τ i . . . τ j whose occurrence in τ i . . . τ j uses only values between a and b.

When V is a negative node, we decompose π into π 1 ⊖ π 2 , with π 1 a common pattern between σ(V L ) and τ i . . . τ h-1 , using only values between c and b in τ i . . . τ h-1 , and π 2 a common pattern between σ(V R ) and τ h . . . τ j , using only values between a and c -1 in τ h . . . τ j , and the proof follows the same steps.

Proposition 2.2. Algorithm 1 has a time complexity in O(min(k, n)kn 6 ).

Proof. Algorithm 1 handles an array M of size O(kn 4 ), where each cell contains a pattern of length at most min(n, k), so that the total space complexity is O(min(n, k)kn 4 ). For filling in the subarrays M (V, _, _, _, _) for all the leaves V of T σ (line 6 to 13 of algorithm 1), the total time complexity is O(kn 5 ). And for any internal node V , filling in one entry of the subarray M (V, _, _, _, _) costs O(min(k, n)n 2 ), since at line 18 (or 23) of algorithm 1, we search for an element of maximal length among O(n 2 ) elements, each of size O(min(k, n)). Consequently, filling in completely this subarray M (V, _, _, _, _) requires a time complexity O(min(k, n)n 6 ). Since there are O(k) internal nodes in T σ , we have the announced result. This complexity can be improved to O(kn 6 ), storing an integer, a label (⊕ or ⊖), and two pointers in M (V, i, j, a, b) (when V is an internal node) instead of a pattern. Namely, if algorithm 1 fills in M (V, i, j, a, b) with the pattern ρ = M (V L , i, h -1, a, c -1) ⊕ M (V R , h, j, c, b), it is sufficient to store in M (V, i, j, a, b) the length of ρ, the label ⊕, and two pointers pointing to the entries M (V L , i, h -1, a, c -1) and M (V R , h, j, c, b) of the array M . At the end of the algorithm, this system of pointers gives a binary separating tree of a longest common pattern π between σ and τ . From this tree, π can be computed in linear time [START_REF] Bose | Pattern matching for permutations[END_REF].

A consequence of properties 2.1 and 2.2 is:

Theorem 2.6. The problem of finding a longest common pattern between two permutations, the first one being separable, is in P .

Longest common pattern between two permutations

The result of theorem 2.6 can be easily extended to classes of permutations that are less restricted than separable permutations. Using the common interval decomposition tree introduced in section 1.2, we will see that a longest common pattern between two permutations σ and τ can be computed in polynomial time as soon as the arity of any prime node in the common interval decomposition tree of σ is bounded by a constant d chosen independently.

Proposition 1.2 states that expanded decomposition trees are a generalization to all permutations of binary separating trees, defined only for the separable ones. From this remark, it becomes natural to try and use expanded decomposition trees in an algorithm for finding a longest common pattern between two general permutations. In the following, we describe such an algorithm and analyze its complexity: it is not a polynomial time algorithm, but the complexity analysis reveals classes of permutations for which the algorithm runs in polynomial time.

First, we notice that the expanded decomposition tree of a permutation can be computed in polynomial time, and even in linear time if we have an independent bound on the arity of the prime nodes. [BCdMR, BXHP] provide a O(n)-time algorithm for computing the common interval decomposition tree T σ of any permutation σ of size n. This tree can easily be labeled in time O(d log d • n) where d is a bound on the arity of prime nodes: for every internal node V of T σ -that are in number O(n) -it is enough to sort its children to find the label of V . Finally, the vertical expansion necessary to obtain the expanded decomposition tree of σ requires again a linear time.

Algorithm 2 takes as an input an expanded decomposition tree T σ of a permutation σ and a permutation τ . It outputs a longest common pattern between σ and τ . It works just like algorithm 1, except for the case of prime nodes in T σ . The procedure in this additional case is described in algorithm 2. /* V is a prime node /* 10:

Let ρ be the simple permutation labeling V 11:

Let d be the arity of V , and V 1 , . . . , V d the children of V , from left to right 12:

for any integers i, j, a and b between 1 and n, i ≤ j, a ≤ b do 13:

M (V, i, j, a, b) ← Longest(S) where

S = ⊙ ρ ( M (V k , h k-1 , h k -1, c ρ k -1 , c ρ k -1) 1≤k≤d ) : i = h 0 ≤ h 1 ≤ . . . ≤ h d = j + 1, a = c 0 ≤ c 1 ≤ . . . ≤ h d = b + 1 14:
end for 15:

end if 16: end for 17: Output : M (root of T σ , 1, n, 1, n) Algorithm 2 uses a more general kind of pattern concatenation than just ⊕ of M (V, _, _, _, _), algorithm 1 computes a longest pattern in a set containing O(n 2 ) elements, whereas in algorithm 2, the set from which we have to extract a longest pattern contains O(n 2d-2 ) elements, if d is the arity of V (see line 13 of algorithm 2).

With no hypothesis on a permutation σ of size k, the only bound we can give on the maximal arity d of a prime node in the expanded decomposition tree of σ is d ≤ k. This bound is optimal since the equality d = k is achieved when σ is a simple permutation. The total time complexity of algorithm 2 is consequently O(min(n, k)kn 2k+2 ), and it is not polynomial. However, if we consider classes of permutations such that the arity of any prime node in their expanded decomposition tree is bounded by a constant d, algorithm 2 has a time complexity O(min(n, k)kn 2d+2 ). In this particular case, algorithm 2 runs in polynomial time.

This can be summarized in the following theorem:

Theorem 3.3. Let d be a integer. Consider the class R of permutations having an expanded decomposition tree with all prime nodes of arity smaller than d.

Then the problem of finding a longest common pattern between a permutation in R and another unrestricted permutation is in P .

Conclusion and open problems

We generalize the algorithm given in [START_REF] Bose | Pattern matching for permutations[END_REF] for the problem longest common pattern. Yet our algorithm seems far from optimal. For example, for separable permutations, our work is based on [START_REF] Bose | Pattern matching for permutations[END_REF] but Ibarra [START_REF] Ibarra | Finding pattern matchings for permutations[END_REF] give a faster (O(n 5 )) algorithm for the pattern involvement problem. Could this algorithm be adapted to the longest common pattern problem? Yet a lower bound is given by the edit distance problem [START_REF] Zhang | Simple fast algorithms for the editing distance between trees and related problems[END_REF] as the edit distance problem between two trees is a special case of the longest common pattern problem as shown in [START_REF] Micheli | Average edit distance between ordered trees[END_REF].

Figure 1 :

 1 Figure 1: Two decorated binary separating trees and the decorated contracted separating tree of σ = 4 2 3 1 6 5 8 9 7

Figure 2 :

 2 Figure 2: Common interval decomposition of σ = 5 1 10 9 6 7 8 11 2 4 3. Common intervals are represented by horizontal lines, with strong intervals corresponding to bold lines.

LFigure 3 :

 3 Figure 3: Common interval decomposition tree for σ = 5 1 10 9 6 7 8 11 2 4 3

Figure 4 :

 4 Figure 4: Labeled decomposition tree and expanded decomposition tree for σ = 5 1 10 9 6 7 8 11 2 4 3

  , b) is a common pattern between σ(V ) and τ i . . . τ j using only values between a and b in τ i . . . τ j . Take for example h = 5 and c = 4.

Algorithm 2

 2 Longest common pattern between two permutations 1: Input : An expanded decomposition tree T σ of a permutation σ of size k and a permutation τ of size n 2: Create an array M : proceed as in algorithm 1 3: Fill in M for the leaves of T σ : proceed as in algorithm 1 4: Fill in the rest of M : 5: for any internal node V of T σ , considering the nodes in the postfix ordering do 6:if V is a positive or a negative node then

The definition of overlapping common intervals follows the intuition: we say that two common intervals I and J are overlapping when I \ J = ∅, J \ I = ∅ and I ∩ J = ∅

By A ≺ B, we mean that ∀a ∈ A, ∀b ∈ B, a < b.

and ⊖. The ρ-concatenation, or pattern concatenation according to ρ, is defined as follows: Definition 3.1. Given ρ = ρ 1 . . . ρ n a permutation of size n, and n patterns π 1 , . . . , π n of respective size k 1 , . . . , k n , the ρ-concatenation of the (π i ) 1≤i≤n is

) for all x between 1 and k i Example 3.2. We can also notice that ⊙ 12 = ⊕ and ⊙ 21 = ⊖.

The idea behind algorithm 2 is quite simple. When filling in M (V, i, j, a, b) for a prime node V labeled by ρ and having d children V 1 . . . V d , we "slice" the intervals {i, . . . , j} and {a, . . . , b} into I 1 . . . I d and A 1 . . . A d respectively, such that I p ≺ I k and A p ≺ A k 2 whenever p < k. Then we ρ-concatenate longest common patterns between the σ(V k ) and τ in the intervals I k of indices and A ρ k of values. With the notation of algorithm 2,

Proposition 3.1. Algorithm 2 is correct: it outputs a longest common pattern between the two permutations σ and τ given in input.

Proof. Similar to the proof of proposition 2.1.

With the notations of the proof of proposition 2.1, in the case of a prime node V labeled by ρ, with children V 1 . . . V d , there exist integers i = h 0 ≤ h 1 ≤ . . . ≤ h d = j + 1 and a = c 0 ≤ c 1 ≤ . . . ≤ c d = b + 1, such that we can decompose π into π = ⊙ ρ (π 1 , . . . , π d ), with π k a common pattern between σ(V k ) and

Using this decomposition of π, we can use the induction hypothesis on the nodes (V k ) 1≤k≤d and finish the proof as before.

In this proof, the trick relies on the fact that a common pattern between σ(V ) and τ i . . . τ j is always a concatenation of common patterns between the children of V and "slices" of τ i . . . τ j . This stability when going from parents to children in the expanded decomposition tree also appears in a paper of Albert and Atkinson [START_REF] Albert | Simple permutations and pattern restricted permutations[END_REF], for example in their lemma 15.

The main difference between algorithms 1 and 2 lies in the complexity analysis. Those two algorithms deal with dynamic programming arrays of the same size, but the cost for computing one entry can be very superior in algorithm 2 than in algorithm 1. Indeed, for any internal node V , in order to fill in one entry


Claire Lacour

Rates of convergence for nonparametric deconvolution

This Note presents original rates of convergence for the deconvolution problem. We assume that both the estimated density and noise density are supersmooth and we compute the risk for two kinds of estimators.

Résumé

Vitesses de convergence en déconvolution nonparamétrique. Cette Note présente des vitesses de convergence originales pour le problème de déconvolution. On suppose que la densité estimée ainsi que la densité du bruit sont (( supersmooth )) et on calcule le risque pour deux types d'estimateurs.

Version française abrégée

On considère le problème de déconvolution suivant :

Y i = X i + ε i i = 1, . . . , n
où les X i sont des variables aléatoires indépendantes et identiquement distribuées de densité inconnue g et les variables aléatoires ε i sont i.i.d de densité connue f ε . Les suites (X i ) et (ε i ) sont de plus supposées indépendantes. L'objectif est d'estimer g à partir des données Y 1 , . . . , Y n . Le cadre d'hypothèses est le suivant. Notons, pour toute fonction u, u * la transformée de Fourier de u : u * (x) = e ixt u(t)dt. On suppose que le bruit est tel que pour tout x de R, f * ε (x) = 0 et qu'il satisfait l'hypothèse suivante :

N1. Il existe s ≥ 0, b ≥ 0, γ ∈ R (γ > 0 si s = 0) et k 0 , k 1 > 0 tels que k 0 (x 2 + 1) -γ/2 exp(-b|x| s ) ≤ |f * ε (x)| ≤ k 1 (x 2 + 1) -γ/2 exp(-b|x| s )
Email address: lacour@math-info.univ-paris5.fr (Claire LACOUR).

Preprint submitted to Elsevier Science

On suppose de plus que g appartient à l'espace

A δ,r,a (L) = {g densité sur R et |g * (x)| 2 (x 2 + 1) δ exp(2a|x| r ) ≤ L} avec r ≥ 0, a ≥ 0, δ ∈ R (δ > 1/2 si r = 0), L > 0.
Lorsque r est strictement positif, la fonction est dite superrégulière et ordinairement régulière sinon, la terminologie étant la même pour le bruit. Ce problème a été largement étudié pour une fonction g appartenant à un espace de Sobolev ou de Hölder (i.e. r = 0) : on peut citer entre autres [START_REF] Carrol | Optimal rates of convergence for deconvolving a density[END_REF], [START_REF] Devroye | Consistent deconvolution in density estimation[END_REF], [START_REF] Fan | On the optimal rates of convergence for nonparametric deconvolution problem[END_REF][START_REF] Fan | Adaptatively local one-dimensional subproblems with application to a deconvolution problem[END_REF], [START_REF] Liu | A consistent non-parametric density estimator for the deconvolution problem[END_REF], [START_REF] Stefanski | Rates of convergence of some estimators in a class of deconvolution problems[END_REF]. Les vitesse de convergences médiocres (en puissances de ln n) obtenues lorsque le bruit est superrégulier (et en particulier pour les distributions gaussiennes) ont conduit à considérer des fonctions g également superrégulières. En premier [START_REF] Pensky | Adaptative wavelet estimator for nonparametric density deconvolution[END_REF] et plus récemment [START_REF] Butucea | Deconvolution of supersmooth densites with smooth noise[END_REF], [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF] et [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF] ont étudié des estimateurs dans ce contexte.

Dans cette Note, nous fournissons des vitesses de convergence exactes et explicites, même dans le cas r > 0, s > 0 pour lequel jusqu'à maintenant les vitesses n'étaient données que de façon implicite, excepté dans des cas très particuliers. Ces vitesses sont calculées pour deux types d'estimateurs.

L'estimateur à noyau classique est le suivant :

ĝn (x) = 1 nh n i=1 K x -Y i h
où K est défini comme la transformée de Fourier inverse de K * (x) = ½ {|x|≤1} /f * ε (x/h). Le risque ponctuel (noté MSE) et le risque intégré (notée MISE) sont établis dans [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF] (voir la proposition 2.1 ci-dessous). Cet estimateur a l'avantage d'être optimal au sens minimax asymptotique fin (voir [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF]) mais ne fournit d'estimateur adaptatif que dans des cas particuliers. C'est pourquoi nous nous intéressons également à l'estimateur de projection introduit dans [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF].

Soit

ϕ(x) = sin(πx)/(πx) et ϕ m,j (x) = √ L m ϕ(L m x -j). En posant u t (x) = (1/2π)(t * (-x)/f * ε (x)), l'estimateur de projection est défini par ĝm (x) = |j|≤Kn âm,j ϕ m,j où âm,j = 1 n n i=1 u * ϕm,j (Y i )
Pour cet estimateur, le calcul du rique intégré est présenté dans [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF] (voir la proposition 2.2 ci-dessous).

Un estimateur adaptatif peut ensuite être défini en utilisant une technique de sélection de modèles (voir [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF] pour les détails).

On peut observer que les deux estimateurs ont les mêmes vitesses de convergence, qui s'obtiennent en minimisant l'ordre du risque. On est alors ramené à résoudre l'équation suivante exp( 2b

h s + 2a h r )h α = O(n) où α = r -2δ -2γ -1 lorsque l'on considère l'erreur intégrée et α = -2δ -2γ + (s -1)
+ lorsque l'on considère l'erreur ponctuelle. Dans la plupart des cas, la solution de cette équation est bien connue (voir Tables 1 et 2). Seul le cas r > 0, s > 0 n'a pas été complètement résolu. C'est ce cas qui est étudié ici : les vitesses sont données explicitement dans le théorème 3.1.

Trois cas sont à différencier : r = s, r < s et r > s. Si r est égal à s, les vitesses (ponctuelles et intégrées) sont en n -a/(a+b) modifiées par un facteur logarithmique. Si r est différent de s, on observe un phénomène original : les vitesses dépendent de l'intervalle ]k/(k + 1), (k + 1)/(k + 2)], k ∈ N auquel appartient r/s ou s/r. Si r est strictement inférieur à s (biais dominant), le terme principal est d'ordre exp[b 0 (ln n) r/s ] avec b 0 = -2a/(2b) r/s , et si r est strictement supérieur à s (variance dominante), le terme principal est d'ordre exp[-d 0 (ln n) s/r ]/n avec d 0 = -2b/(2a) s/r .

Ainsi ces vitesses originales décroissent plus vite que n'importe quelle fonction logarithmique. De plus elles sont optimales lorque les bornes inférieures correspondantes sont connues, c'est-à-dire r = s = 1 (voir [START_REF] Tsybakov | On the best rate of adaptive estimation in some inverse problems[END_REF]) et r < s (voir [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF]). Il est à noter qu'étant donné la complexité de ces vitesses, il est réellement intéressant de construire des estimateurs adaptatifs, c'est-à-dire des estimateurs dont le risque atteint automatiquement les vitesses minimax.

Introduction

We consider the following deconvolution problem:

Y i = X i + ε i i = 1, . . . , n
where the X i 's are independent and identically distributed random variables with an unknown density g and the random variables ε i are i.i.d with known density f ε . Moreover (X i ) and (ε i ) are independent.

The aim is to estimate g from data Y 1 , . . . , Y n . The hypothesis framework is the following. Denote, for all function u, u * the Fourier transform of u: u * (x) = e ixt u(t)dt. We suppose that noise is such that for all x in R, f * ε (x) = 0 and that it satisfies the following assumption:

N1. There exist s ≥ 0, b ≥ 0, γ ∈ R (γ > 0 if s = 0) and k 0 , k 1 > 0 such that k 0 (x 2 + 1) -γ/2 exp(-b|x| s ) ≤ |f * ε (x)| ≤ k 1 (x 2 + 1) -γ/2 exp(-b|x| s )
We assume that g belongs to the space

A δ,r,a (L) = {g is a probability density on R and |g * (x)| 2 (x 2 + 1) δ exp(2a|x| r ) ≤ L} with r ≥ 0, a ≥ 0, δ ∈ R (δ > 1/2 if r = 0), L > 0.
When r > 0 the function is known as supersmooth, and as ordinary smooth else. The terminology is the same for noise. This problem has been extensively studied for a function g belonging to a Sobolev or Hölder class (i.e. r = 0): see among others [START_REF] Carrol | Optimal rates of convergence for deconvolving a density[END_REF], [START_REF] Devroye | Consistent deconvolution in density estimation[END_REF], [START_REF] Fan | On the optimal rates of convergence for nonparametric deconvolution problem[END_REF][START_REF] Fan | Adaptatively local one-dimensional subproblems with application to a deconvolution problem[END_REF], [START_REF] Liu | A consistent non-parametric density estimator for the deconvolution problem[END_REF], [START_REF] Stefanski | Rates of convergence of some estimators in a class of deconvolution problems[END_REF]. The bad rates of convergence (power of ln n) for supersmooth noise (and then in particular for Gaussian distributions) lead to consider supersmooth functions. First [START_REF] Pensky | Adaptative wavelet estimator for nonparametric density deconvolution[END_REF] and more recently [START_REF] Butucea | Deconvolution of supersmooth densites with smooth noise[END_REF], [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF] and [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF] studied estimators in this context.

The contribution of this paper is to provide exact and explicit rates of convergence, even in the case r > 0 and s > 0 where up to now the rates were not explicitly available except in very particular cases.

Estimators and preliminar results

Estimators

The classical kernel estimator is the following:

ĝn (x) = 1 nh n i=1 K x -Y i h
where K is the function defined as the inverse Fourier transform of K * (x) = ½ {|x|≤1} /f * ε (x/h). The pointwise mean squared error (denoted by MSE) and mean integrated squared error (denoted by MISE) are established in [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF]: Proposition 2.1 If g belongs to A δ,r,a (L), then under Assumption N1,

M ISE = E g -ĝn 2 2 = O h 2δ exp(-2a/h r ) + h s-1-2γ
n exp(2b/h s ) and

M SE = E|g(x) -ĝn (x)| 2 = O h 2δ+r-1 exp(-2a/h r ) + min(1, h s-1 ) h s-1-2γ n exp(2b/h s )
This estimator has the advantage to be optimal in the sharp asymptotic minimax sense (see [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF]) but provides an adaptive estimator only in particular cases. That is why we present the projection estimator introduced in [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF]. Let ϕ(x) = sin(πx)/(πx) and ϕ m,j

(x) = √ L m ϕ(L m x -j). Consider u t (x) = (1/2π)(t * (-x)/f * ε (x)).
Then the projection estimator is defined by

ĝm (x) = |j|≤Kn âm,j ϕ m,j
where âm,j = 1

n n i=1 u * ϕm,j (Y i )
For this estimator, the following result is proved in [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF]: Proposition 2.2 Assume that f ε ∈ L 2 (i.e. γ > 1/2 when s = 0) and that g is a L 2 function which verifies x 2 g 2 (x)dx ≤ M . If g belongs to A δ,r,a (L), then under Assumption N1,

M ISE = E g -ĝm 2 2 = O L -2δ m exp(-2aπ r L r m ) + L 2γ+1-s m n exp(2bπ s L s m )
Then, an adaptive estimator can be defined using a model selection method (see [START_REF] Comte | Penalized contrast estimator for adaptive density deconvolution[END_REF] for details).

Rates of convergence

We can observe that both estimators have the same L 2 rate of convergence (take h -1 = πL m ). To compute this rate, we have to minimize the risk orders in h (or L m ). By setting to zero the derivative of this quantity we obtain the equation exp( 2b

h s + 2a h r )h α = O(n) (1)
where α = r -2δ -2γ -1 if we consider the integrated error and α = -2δ -2γ + (s -1) + if we consider the pointwise error. In most cases, the solution of this equation is well known and leads to the following tables where different regularities for g and f ε are examined:

s = 0 s > 0 r = 0 n - 2δ 2δ+2γ+1 (ln n) -2δ s r > 0 (ln n) 2γ+1 r n Theorem 3.1
Table 1 Rates of convergence for the MISE.

s = 0 s > 0 r = 0 n 1-2δ 2δ+2γ (ln n) 1-2δ s r > 0 (ln n) 2γ+1 r n Theorem 3.1 Table 2
Rates of convergence for the MSE.

Except for the bottom right cells (to be completed in the next section), these rates are known to be optimal minimax rates: see [START_REF] Fan | On the optimal rates of convergence for nonparametric deconvolution problem[END_REF] and [START_REF] Butucea | Deconvolution of supersmooth densites with smooth noise[END_REF] for the lower bounds.

Results

The rates of convergence in the case (r > 0, s > 0) depend on the integer k such that r/s or s/r belongs to the interval (k/(k + 1), (k + 1)/(k + 2)]: Theorem 3.1 We assume r > 0 and s > 0. Let k ∈ N and λ = µ -1 = r/s. Then

• if r = s, if ξ = [2δb + (s -2γ -1)a]/[(a + b)s] M ISE = O n -a/(a+b) (ln n) -ξ ; M SE = O n -a/(a+b) (ln n) -ξ+ (1-s) + b (a+b)s • if r < s and k k + 1 < λ ≤ k + 1 k + 2
, there exist reals b i such that

M ISE = O (ln n) -2δ/s exp[ k i=0 b i (ln n) (i+1)λ-i ] ; M SE = O (ln n) (-2δ-r+1)/s exp[ k i=0 b i (ln n) (i+1)λ-i ] • if r > s and k k + 1 < µ ≤ k + 1 k + 2 , there exist reals d i such that, M ISE = O (ln n) (1+2γ-s)/r n exp[- k i=0 d i (ln n) (i+1)µ-i ] ; M SE = O (ln n) (1+2γ-s-(s-1)+)/r n exp[- k i=0 d i (ln n) (i+1)µ-i ]
The coefficients b i and d i are computable, see Section 4 for the exact form of reals b i . Notice that these original rates have the property to decrease faster than any logarithmic function. Moreover, they are optimal in the cases where the corresponding lower bounds are known, i.e r = s = 1 (see [START_REF] Tsybakov | On the best rate of adaptive estimation in some inverse problems[END_REF]) and r < s (see [START_REF] Butucea | Sharp optimality for density deconvolution with dominating bias[END_REF]). We can also remark that, given the complexity of these rates, it is woth finding adaptive estimators, i.e. estimators whose risk automatically achieves the minimax rates.

Proof of Theorem 3.1

• If r = s, we check that h * = (2a+2b) 1/s ln n + α s ln ln n -1/s satisfies Equation (1). The corresponding risks are easily obtained. 1) is verified. The risks are then computed by using the same tools.

• If r < s and k k + 1 < λ ≤ k + 1 k + 2 , let h * = (2b)
• The reasoning is the same in the case r > s

k

  i=0 b i (ln n) (i+1)λ-i ] -1/s exp[ 2B h * s + 2A h * r + α ln h * ] = Kn exp[ 2a (2b) λ (ln n) λ (1 + λu n + .. + λ(λ -1)..(λ -k) ln n) (i+1)λ-iα s ln(1 + o(1))] with u n := α s ln ln n ln n + k i=0 b i (ln n) (i+1)λ-(i+1) By noting that for 2 ≤ j ≤ k + 1, u j n = k+1 i=j p1+..+pj=i b p1-1 ..b pj -1 (ln n) (λ-1)i + o((ln n) (λ-1)(k+1) ), we compute exp[ 2b h * s + 2a h * r + α ln h * ] = Kn exp[o(1) + k+1 i=0 M i (ln n) (i+1)λ-i + o((ln n) (k+2)λ-(k+1) )] with M 0 = b 0 + 2a (2b) λ and for 1 ≤ i ≤ k + 1, M i = b i + 2a (2b) λ i j=1 λ . . . (λ -j + 1) j! p1+..+pj =i b p1-1 ..b pj -1by denoting b k+1 = 0. The b i 's are recursively defined by M 0 = ... = M k = 0 and thus Equation (

  1/s [ln n + α s ln ln n +