Cédric Herzet
email: cedric.herzet@inria.fr

Mamadou Diallo

• M Diallo

Performance guarantees for a variational "multi-space" decoder

Keywords: Model-order reduction, Multi-space prior information, Performance guarantees

.

Introduction

Many approximation methods encountered in the domain of model-order reduction rely on the following ingredients: 1) some "prior" knowledge about the location of the target solution h in the ambient space H. 1 This prior information typically takes the form of a subset M prior ⊂ H such that h ∈ M prior . 2) a set of partial observations of h , say {y j } m j=1 . In this paper, we assume that y j = w j , h where {w j } m j=1 is a basis of some m-dimensional subspace W m .

Examples of problems exploiting these two ingredients include "projectionbased" reduction of linear parametrized partial differential equations (PPDE) [START_REF] Quarteroni | Reduced Basis Methods for Partial Differential Equations[END_REF], some approximation methods for non-linear operators (e.g., GEIM [START_REF] Maday | The generalized empirical interpolation method: Stability theory on hilbert spaces with an application to the stokes equation[END_REF], DEIM [START_REF] Chaturantabut | Nonlinear model reduction via discrete empirical interpolation[END_REF], Gappy POD [START_REF] Everson | Karhunen-Loève procedure for gappy data[END_REF]) or data-assimilation schemes with reduced-order models [START_REF] Binev | Data assimilation in reduced modeling[END_REF].

Both M prior and {y j } m j=1 provide some information about the location of h in H. In particular, in the case of linear measurements, observing {y j } m j=1 indicates that h must belong to some m-dimensional affine subspace of H: h ∈ P y {h ∈ H : w j , h = y j for j = 1 . . . m}

= h + h : h ∈ W ⊥ m , (1)
where W ⊥ m is the orthogonal complement of W m in H. The goal of an approximation method then consists in combining effectively (in some sense) these two sources of information to provide a good approximation ĥ of h .

The choice of M prior and W m usually plays a crucial role in the performance of the approximation methods since they determine a trade-off between accuracy and complexity. A standard option consists in choosing M prior as:

M prior = {h ∈ H : dist(h, V n ) ≤ ˆ n } (2)
where V n is an n-dimensional subspace, dist(h, V n ) min h ∈Vn hh and ˆ n ≥ 0. Prior model [START_REF] Babuska | Error-bounds for finite element method[END_REF] indicates that the sought solution h is close (up to some error bounded by ˆ n ) to an n-dimensional subspace V n . This simple model has been exploited in most contributions of the literature, see e.g., [START_REF] Chaturantabut | Nonlinear model reduction via discrete empirical interpolation[END_REF][START_REF] Everson | Karhunen-Loève procedure for gappy data[END_REF][START_REF] Maday | The generalized empirical interpolation method: Stability theory on hilbert spaces with an application to the stokes equation[END_REF][START_REF] Quarteroni | Reduced Basis Methods for Partial Differential Equations[END_REF].

Recently, triggered by the work by Binev et al. [START_REF] Binev | Data assimilation in reduced modeling[END_REF], more refined definitions of M prior have been considered in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF]. In these works, M prior is built from a sequence of embedded approximation subspaces, that is

M prior = ∩ n k=0 {h ∈ H : dist(h, V k ) ≤ ˆ k } (3)
where

ˆ k ≥ 0, dim(V k ) = k and V 0 ⊂ V 1 ⊂ . . . ⊂ V n . (4)
Prior ( 3) is usually referred to as "multi-space" prior. It indicates that h is located at some bounded distances of a sequence of embedded subspaces {V k } n k=0 . This kind of information may for example stem from prior physical knowledge about the process under study or come as a by-product of the procedures computing the approximation subspaces, see e.g., [START_REF] Quarteroni | Reduced Basis Methods for Partial Differential Equations[END_REF].

Two different decoders exploiting multi-space information were proposed in [START_REF] Binev | Data assimilation in reduced modeling[END_REF] and [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF]. In [START_REF] Binev | Data assimilation in reduced modeling[END_REF], the authors suggested to approximate h by some feasible point of the set M prior ∩ P y . They showed theoretically that this particular decoder may outperform standard single-space procedures in specific settings. In [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF] a different approach was considered: an approximation ĥ of h was computed as the solution of an optimization problem with constraints ĥ ∈ M prior and ĥ ∈ V n . 2 The authors showed empirically that this procedure may have superior approximation performance than single-space approaches in different applicative domains. The goal of the present work is to give some theoretical ground to these observations. More specifically, we derive theoretical guarantees relating the approximation performance of the decoder considered in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF] to the distances between the target solution h and the approximation subspaces {V k } n k=0 . 3 We show morerover in several examples that, for particular choices of subspaces {V k } n k=0 and W m , this decoder can lead to superior approximation guarantees than its standard counterpart [START_REF] Babuska | Error-bounds for finite element method[END_REF].

The rest of this paper is organized as follows. In Section 2, we present several applicative domains where the general framework described in this introduction may apply. In Section 3, we recall the standard reconstruction guarantees holding in the single-space setup, i.e., when the prior information exploited by the approximation procedure only involves one n-dimensional subspace V n . Section 4 is dedicated to the multi-space setup. In Section 4.1, we give a precise definition of the multi-space decoder considered in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF] and provide some theoretical guarantees on its performance (see Theorem 2). In Section 4.2, we particularize this result to different setups and show that the multi-space decoder proposed in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF] has provably more favorable performance than the single-space approach in these cases. We also provide a graphical interpretation of the geometry underlying the single and multi-space problems, which leads to some further insights into the scenarios where multispace priors may be of interest. Finally, in Section 4.3 we discuss the similarities and differences between the multi-space decoders proposed in [START_REF] Binev | Data assimilation in reduced modeling[END_REF] and [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF]. In particular we show that when {w j } m j=1 is an orthonormal basis, the theoretical guarantees derived in this paper are slightly more favorable than those in [START_REF] Binev | Data assimilation in reduced modeling[END_REF] but that no general ordering of the decoder performances is possible. The proof of our main result is finally detailed in Section 5.

Standard approximation problems

In this section, we present several particularizations of the general approximation problem introduced in Section 1. This allows us to replace the different contributions of the literature dealing with multi-space priors in their context and motivate our theoretical results in Section 4 from an applicative perspective.

Reduced-order modeling

In the context of reduced-order modeling of PPDEs, see e.g., [START_REF] Quarteroni | Reduced Basis Methods for Partial Differential Equations[END_REF], the goal consists in approximating the elements of a set of solutions:

M = {h ∈ H : PDE(h , θ) = 0 for some θ ∈ Θ} ( 5
)
where Θ is some set of parameters and PDE(h , θ) = 0 is an abstract notation for the PPDE. In the case of elliptic PPDEs, the weak formulation of the differential problem reads:

Find h ∈ H such that a θ (h , h) = b θ (h) ∀h ∈ H, (6)
where a θ : H × H → R and b θ : H → R are some problem-specific bilinear and linear operators depending on θ.

If H is infinite (or very high) dimensional, solving (6) may be computationally expensive. In order to alleviate this problem, model reduction combines a subset of the constraints defined in [START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF] with some prior information on h to compute an approximated solution ĥ. For example, in the context of singlespace prior information (2), a typical approach consists in imposing " ĥ ∈ V n " where V n is some proper n-dimensional subspace. Moreover, the subset of constraints in [START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF] used in the reduced model is usually defined as " ĥ must be such that a θ ( ĥ,

z j ) = b θ (z j ) ∀1 ≤ j ≤ m", (7)
where {z j } m j=1 is a family of linearly-independent elements of H. We note that [START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF] is equivalent to imposing a constraint of the form (1). Indeed, letting w j ∈ H be the Riesz representer of the linear form a θ (•, z j ), we have that

w j , h = a θ (h , z j ) = b θ (z j ), (8)
where the last equality follows from [START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF]. Hence, defining y j b θ (z j ) and W m span {w j } m j=1 , it is easy to see that ( 7) is tantamount to enforcing " ĥ ∈ P y ". We note that in this particular example, the definition of w j depends on θ and thus varies for each problem instance.

Another example of approximation problem closely-connected to model reduction of PPDEs is the time-integration procedure considered in [START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF]. In this contribution, the use of a semi-implicit integration scheme requires the resolution of a problem similar to (6) at each time step. In this particular case, the bilinear operator appearing in the approximation problem (and thus the definition of the subspace W m ) depends on the solution of the integration scheme at the previous time step. The authors exploit a particular multispace procedure (see Section 4.1) to find a good approximated solution of this problem.

Generalized empirical interpolation and data assimilation

The "Generalized Empirical Interpolation" and "Data Assimilation" problems, considered respectively in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Maday | The generalized empirical interpolation method: Stability theory on hilbert spaces with an application to the stokes equation[END_REF] and [START_REF] Binev | Data assimilation in reduced modeling[END_REF], are closely connected to the modelreduction setup described above. In these frameworks, the set of target solutions M takes the form of a compact subset of H and the goal is to recover h ∈ M (up to some accuracy) from a set of m linear measurements of h , say { w j , h } m j=1 . As a typical example of such a setup, one may mention (see e.g., [START_REF] Maday | The generalized empirical interpolation method: Stability theory on hilbert spaces with an application to the stokes equation[END_REF]) the case where M corresponds to a subset of the square-integrable functions on a bounded domain Ω ⊂ R d (e.g., defining the response of a physical system to different operating conditions), and {w j } m j=1 characterizes the behavior of some sensing devices probing the state of the physical system at different locations, e.g.,4

w j : Ω → R x → exp - x-mj 2 2 2σ 2 (9)
where m j ∈ Ω and σ 2 > 0 are some parameters. Going back to the generic approximating problem introduced in Section 1, we see that the "observation" subspace W m is here trivially given by W m = span({w j } m j=1 ). Different priors have been considered in the literature to complement the collected observations. The standard GEIM is based on single-space information, where the estimate ĥ is contrained to belong some n-dimensional subspace V n . Recently, GEIM has been reformulated in the context of multi-space information (3) in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF]. The decoder considered in this contribution is presented in Section 4.1. In [START_REF] Binev | Data assimilation in reduced modeling[END_REF] both single-space and multi-space prior information were considered in the context of data assimilation. The multi-space proposed in this contribution is different from the one considered in the present paper and its definition is recalled in Section 4.3.

The single-space approximation

In this section, we review some of the guarantees of performance holding in the single-space setup, that is when prior information of the form (2) is available.

The material presented in this section is not novel but is intended to support our discussion in the next section.

When dealing with single-space approximations, a common option to compute an approximation of h consists in solving the following problem:

Find ĥSS ∈ V n such that w j , ĥSS = y j for all j = 1 . . . m.

(

) 10
In a nutshell, ( 10) is equivalent to finding an element of V n satisfying the observation constraints. Problem [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] is for example the standard approach used in reduced modeling or GEIM. The setup of most interest in this case is n = m, since a unique solution to [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] then exists under mild conditions (see Theorem 1 below). We note that, in the general case h / ∈ V n and therefore ĥSS = h . Nevertheless, depending on the problem at stake, ĥSS may "close" to h in some situations. In fact, the quality of the approximation obtained from the singlespace problem [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] can be characterized as a function of the subspaces W m and V n involved in this problem. More specifically, typical recovery results in the single-space setup take the form

ĥ -h ≤ C(W m , V n ) dist(h , V n ) ( 11
)
where C(W m , V n ) is some constant depending on W m and V n . For example, in the context of Hilbert spaces the following result holds:

Theorem 1 Let m = n and σ min inf h∈Vn sup h ∈Wm h , h h h . ( 12
)
If σ min > 0, (10) admits a unique solution ĥSS and

h -ĥSS ≤ σ -1 min dist(h , V n ). (13)
Moreover, (13) is tight in the following sense:

∀ n > 0, ∃h ∈ H such that dist(h , V n ) = n and h -ĥSS = σ -1 min n . (14)
A proof of this result can be found in [START_REF] Maday | The generalized empirical interpolation method: Stability theory on hilbert spaces with an application to the stokes equation[END_REF]Theorem 2.4] for the first part and [3, Theorem 2.9] for the second. We note that (13) is only one particular example of guarantee valid for single-space approximations in Hilbert spaces. Different guarantees exist in other setups. For example, in Banach spaces [START_REF] Quarteroni | Reduced Basis Methods for Partial Differential Equations[END_REF] holds with a slightly larger constant

C(W n , V n ) = 1 + σ -1 min , see [9, Theorem 1.3].
In the context of reduced-order model (cf. Section 2.1), similar results exist where the constant C(W n , V n ) is related to the "inf-sup stability" of the operator characterizing the PPDE, see e.g., [START_REF] Babuska | Error-bounds for finite element method[END_REF]. Other works also characterize the link between the worst-case approximation error and the Kolmogorov width of a set of target solutions, see e.g., [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF]. In what follows, we will restrict our attention to (13), since it is valid in the context of Hilbert spaces considered in this paper, and is more amenable to comparisons with our main result stated in Theorem 2.

Before concluding this section, let us mention that the authors in [3, Theorem 2.9] showed that the approximation computed in [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] is in fact the "worst-case" optimal solution when single-space prior information (2) is available. More specifically, if m = n and σ min > 0, (10) can be rewritten as

ĥSS = arg min h∈H sup h ∈Mprior∩Py h -h , (15)
where M prior is defined in [START_REF] Babuska | Error-bounds for finite element method[END_REF]. Hence, as long as prior ( 2) is considered, ĥSS corresponds to the element of H leading to the best worst-case approximation error.

The multi-space approximation

In this section, we present a theoretical result supporting the performance of the "multi-space" decoder considered in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF]. We state our main result in Theorem 2 and provide two examples of scenarios in which the multi-space decoder has provably better performance guarantees than the standard "single space" decoder [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF]. We then provide a graphical illustration of the geometry underlying the single-space and the multi-space problems. Finally, we elaborate on the differences and similarities between the results presented in this paper and those in [START_REF] Binev | Data assimilation in reduced modeling[END_REF].

4.1 Recovery guarantee for the multi-space decoder considered in [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF] Before stating our result, we recall the definition of the multi-space decoder considered in [1, 6-8]:5

Find ĥMS ∈ arg min

h∈Vn m j=1 (y j -w j , h ) 2 (16) subject to dist(h, V k ) ≤ ˆ k , k = 0 . . . n.
In the sequel, we assume without loss of generality that the vectors {w j } m j=1

are linearly independent (but not necessarily orthonormal). We also suppose that the subspaces {V k } n k=0 are embedded, that is, obey (4).

We note that (16) can be seen as an extension of the single-space decoder (10) discussed in Section 3. More specifically, if one removes the constraints in (16), the problem becomes Find ĥ ∈ arg min

h∈Vn m j=1 (y j -w j , h ) 2 . ( 17
)
It is then quite easy to see that the unique solution of this problem is given by [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] as soon as m = n and σ min > 0. On the other hand, we note that the constraints in (16) define a set of feasible points of the form (3). Hence, if the subspaces {V k } n k=0 and widths {ˆ k } n k=0 are properly chosen, these constraints add some valuable information about the position of the sought solution h in H. We may thus expect the multi-space decoder (16) to lead to enhanced performance in some specific situations. In this section, we provide a mathematical support to this intuition.

We remind the reader that we assume that h ∈ M prior . When multi-space prior (3) is considered, this entails that

dist(h , V k ) ≤ ˆ k for all k = 0 . . . n, (18)
that is the target solution h is a feasible point of the optimization problem (16). Under this assumption, we provide a mathematical characterization of the performance achievable by the multi-space decoder (16). More specifically, we derive an upper bound on the approximation error ĥMSh depending on the distance between h and the different approximation subspaces

{V k } n k=0 .
Our result is presented in Theorem 2 below. In order to state our result we need to introduce the following quantities. We first define the short-hand notations

k = dist(h , V k ) for all k = 0 . . . n, (19)
and

γ = sup h∈V ⊥ n , h =1   m j=1 w j , h 2   1 2 , ( 20
)
where V ⊥ n is the orthogonal complement of V n in H. We let {v j } n j=1 be an orthonormal basis of V n such that

V k = span {v j } k j=1 . (21)
We note that such a basis always exists since we assume that the sequence of subspaces {V k } n k=0 obeys (4). We define the Gram matrix G as

G = [ w i , v j ] i,j ∈ R m×n , (22)
and let

δ j = n k=1 |x kj |(ˆ k-1 + k-1 ) j = 1 . . . n, (23)
where x kj are the elements of matrix X appearing in the singular value decomposition of G, that is G = UΛX T , where U ∈ R m×m , X ∈ R n×n are orthogonal matrices and Λ ∈ R m×n is the diagonal matrix of singular values {σ j } min(m,n) j=1

. In the sequel we will consider the extended set {σ j } n j=1 by using the following convention: if n > m, we define σ j = 0 for all j > m. Without loss of generality, we assume that the singular values {σ j } n j=1 are sorted by decreasing order of magnitude.

Finally, we introduce the function

B {δ j } n j=1 , τ =      n j= +1 δ 2 j + ρ δ 2 + 2 n 1 2 if n j=1 σ 2 j δ 2 j ≥ τ 2 , n j=1 δ 2 j + 2 n 1 2 otherwise, (24)
where is the largest integer such that

n j= σ 2 j δ 2 j ≥ τ 2 , (25)
and ρ ∈ [0, 1] is defined as

ρσ 2 δ 2 + n j= +1 σ 2 j δ 2 j = τ 2 . ( 26
)
Using these notations, our result reads:

Theorem 2 Assume h verifies (18) and let y j = w j , h for j = 1 . . . m.

Then any solution ĥMS of (16) verifies

h -ĥMS ≤ B {δ j } n j=1 , 2γ n . ( 27
)
Moreover, if σ n > 0, (16) admits a unique solution.

The proof of Theorem 2 is presented in Section 5. We can make the following remarks about the result presented in Theorem 2. First, we note that all the quantities appearing in (24)-( 26) are easily computable so that the upper bound B({δ j } n j=1 , 2γ n ) can be straightforwardly evaluated numerically. 6 We also notice that the performance guarantee stated in Theorem 2 is valid for any values of m and n. In particular, if ˆ k < ∞ for all k = 0 . . . n, we always have B({δ j } n j=1 , 2γ n ) < ∞ since all the terms in the right-hand side of (24) are finite. This means that the multi-space decoder (16) always leads to a bounded error even if m < n. This is in good accordance with the fact that both h and ĥMS belong to the bounded set M prior ∩ P y . We note however that when m < n, we have σ n = 0 and the solution of ( 16) is usually no longer unique.

Finally, let us point out that the upper bound ( 27) is not sharp. This can for example be seen when m = n, ˆ k → ∞ for k = 0 . . . n, and {w j } m j=1 is an orthonormal basis. In this case, problem (16) reduces7 to (17) which, in turn, is equivalent to the single-space problem [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] as soon as σ n > 0. 8 Considering the latter setup, we thus have ĥMS = ĥSS and a sharp upper bound on the performance of the multi-space decoder is given by (13) in Theorem 1:

h -ĥSS ≤ σ -1 n n . ( 28
)
On the other hand, particularizing (27) to this setup, we easily obtain that

h -ĥSS ≤ 4σ -2 n + 1 1 2 n , (29)
which is slightly more pessimistic than (28). In the next section, we nevertheless emphasize that the result stated in Theorem 2 is sufficiently accurate to show the superiority of the multi-space decoder over its single-space counterpart in some specific setups.

Comparison with the single-space approximation

We now turn our attention to the comparison between the performance achievable by the standard single-space decoder [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] and its multi-space version (16). More specifically, we show that the recovery guarantees obtained in the multispace setup may be much more favorable than in the single-space setup in some specific settings. In order to allow a comparison between the single and multispace decoders, we consider the case where m = n and assume that {w j } m j=1 is an orthonormal basis. We note that, in such a case, we have σ min = σ n and γ ≤ 1. The specific settings considered hereafter are inspired from [START_REF] Binev | Data assimilation in reduced modeling[END_REF] and are described in Examples 1 and 2. They correspond to different choices of matrix X appearing in the singular value decomposition of the Gram matrix G defined in (22). Since the latter matrix directly depends on the bases {w j } m j=1 and {v j } n j=1 , these examples thus correspond to some particular choices of the observation and approximation subspaces W m and {V k } n k=0 .

Example 1 We first assume that X = I n in the singular-value decomposition of G. We set ˆ j = j and assume that

j =    1 j = 0 . . . n -3, 1 2 j = n -2, n -1, j = n, (30)
for some 1. Moreover, we let

σ j =    1 j = 1 . . . n -3, 1 2 j = n -2, n -1, j = n. (31)
In this setup, the upper bound (13) of Theorem 1 becomes:

ĥSS -h ≤ σ -1 n dist(h , V n ) = -1 = 1. (32)
We note that this is sharp. In particular, there exists some h ∈ H such that dist(h , V n ) = n and the equality holds in (32) (see second part of Theorem 1). On the other hand, because X = I n and ˆ j = j , we have

δ j = ˆ j-1 + j-1 = 2 j-1 . (33)
The index appearing in the definition of the function B in ( 24) is then smaller or equal to n -1 since

σ 2 n δ 2 n = σ 2 n (2 n-1 ) 2 = 4 3 4 2 , σ 2 n-1 δ 2 n-1 = σ 2 n-1 (2 n-2 ) 2 = 4 2 ,
and thus

σ 2 n-1 δ 2 n-1 + σ 2 n δ 2 n ≥ 4 2 ≥ 4γ 2 2 (34) since γ ≤ 1. The upper bound in Theorem 2 becomes h -ĥMS ≤ δ 2 n-1 + δ 2 n + 2 n 1 2 , = 4 + 4 + 2 1 2 , ≤ 3
1 2 . ( 35
)
Hence the bound in the multi-space setup (35) can be arbitrarily small as compared to (32) when → 0.

Example 2 We now consider X = n -1 2 1 n×n where 1 n×n is an n × n matrix of 1's. We set ˆ j = j and assume that

j =    1 2 j = 0, 1 2(n-1) j = 1 . . . n -1, j = n, (36)
for some n -1 (Note that we must have: ≤ 1 2(n-1) by definition). Moreover, we let

σ j = σ j = 1 . . . n -1, 2 j = n, (37)
for some 1 ≥ σ > whose value will be specified below. With these choices, the upper bound (13) of Theorem 1 becomes:

ĥSS -h ≤ σ -1 n dist(h , V n ) = -2 = -1 . ( 38
)
On the other hand, we have

δ j = n k=1 |x kj |(ˆ k-1 + k-1 ), = 2n -1 2 n k=1 k-1 , = 2n -1 2 . ( 39
)
By choosing σ such that (we remind the reader that σ n-1 = σ by definition (37))

σ 2 n-1 δ 2 n-1 + σ 2 n δ 2 n = 4 2 , ( 40
)
we obtain that index appearing in the definition of the function B in ( 24) is smaller or equal to n -1 since γ ≤ 1. The upper bound in Theorem 2 then reads

h -ĥMS ≤ δ 2 n-1 + δ 2 n + 2 n 1 2 , = 4n -1 + 4n -1 + 2 1 2 , ≤ 3n -1 2 , ( 41
)
where the last inequality follows from our initial assumption n -1 .

We conclude this section by providing a graphical representation of the geometry of problem (16). Fig. 1 gives an illustration of the feasible set and the iso-contours of the cost function appearing in (16); these quantities are plotted in the plane V n for m = n = 2.

If σ n > 0, it can be seen (see Appendix A) that the cost function f (h) n j=1 (y jw j , h ) 2 can also be rewritten for any h ∈ V n as

f (h) = n j=1 σ 2 j v * j , ĥSS -v * j , h 2 , ( 42
)
where

v * j = n i=1 x ij v i . ( 43
)
Here, the elements x ij 's correspond to the components of the orthonormal matrix X appearing in the singular value decomposition of G. From (42), we thus see that the iso-contours of f (h) in V n correspond to n-dimensional ellipsoids with center equal to ĥSS and principal axes equal to v * j n j=1 . Moreover, the elongation of the ellipsoids along each axis v * j is inversely proportional to the singular value σ j .

From a geometric point of view, bad recovery guarantees in the singlespace setup corresponds to situations where the "iso-contour" ellipsoids is very elongated along (at least) one direction (i.e., σ n is small): the center of the ellipsoid ĥSS may then be quite distant from the optimal orthogonal projection P Vn (h ).

The prior information used in the multi-space decoder (16) may provide a solution to this problem by constraining ĥMS to belong to some prespecified feasible set. Fig. 1 gives an illustration of such a situation. The feasible set defined by the constraints in ( 16) is represented by the gray shaded area. In the simplified setup considered here, it corresponds to the intersection of two sets: the constraint associated to V 0 imposes ĥMS to belong to a ball of radius ˆ 0 ; the constraint corresponding to V 1 requires that ĥMS does not deviate from the line passing through v 1 by more than ˆ 1 . The multi-space estimate ĥMS then corresponds to the element of the feasible set leading to the smallest value of the cost function. We see in Fig. 1 that constraining the estimate ĥMS to belong to the feasible set prevents it from deviating too far from P Vn (h ). In particular, in the simple example described in Fig. 1, the multi-space estimate leads to better approximation performance than its single-space counterpart.

The same type of conclusions can in fact be drawn in more general setups: the multi-space decoder ( 16) is able to enhance the performance of the singlespace approach [START_REF] Maday | Convergence analysis of the generalized empirical interpolation method[END_REF] as soon as the prior information used in (16) can compensate for large deviations of the "iso-contour" ellipsoids. The gain achievable in the multi-space setup thus depends on the relative configuration of the "isocontour ellipsoids" and the feasible set. We note that the shape of the ellipsoids depends on the basis {v * j } n j=1 and the singular value {σ j } n j=1 . Similarly, the feasible set is fully defined by the basis {v j } n j=1 and widths {ˆ j } n j=0 . Moreover, {v j } n j=1 and {v * j } n j=1 only differ up to an orthogonal transformation X, see (43). This explains why the parameters {σ j } n j=1 , {ˆ j } n j=0 and {δ j } n j=1 (which depend on X) play a crucial role in the characterization of the performance of the multi-space decoder in Theorem 2.

ˆ 1 ˆ 0 • ĥSS • ĥMS v * 2 v * 1 •P V n (h ) v 2 v 1 V n "iso-contours" of n j=1 (y j -w j , h ) 2 {h : dist(h, V 0 ) ≤ ˆ 0 } {h : dist(h, V 1 ) ≤ ˆ 1 }

Comparison with the result by Binev et al. in [3]

In [START_REF] Binev | Data assimilation in reduced modeling[END_REF], the authors considered a different multi-space decoder:

Find ĥMSU ∈ M prior ∩ P y , (44)
where M prior is defined in [START_REF] Binev | Data assimilation in reduced modeling[END_REF]. In this section, we elaborate on the connections and the differences between decoders (16) and (44).

Let us first mention that, although exploiting some "multi-space" prior, the two decoders differ in their definition. Whereas (16) addresses a variational problem where the estimate ĥMS is the solution of a convex optimization problem, (44) focusses on a "feasility" problem, where one element in the intersection of M prior and P y must be found. The two decoders also differ in the constraints they impose on the solution: (44) enforces ĥMSU to belong to the affine subspace P y , whereas (16) imposes ĥMS ∈ V n but minimizes some distance (depending on {w j } m j=1 ) between ĥMS and P y . Although these choices may lead to some differences in the implementation of these decoders, we will restrict ourselves hereafter to a discussion of their performance.

In [START_REF] Binev | Data assimilation in reduced modeling[END_REF], the authors derived an upper bound on the approximation error achievable by the multi-space decoder (44). In order to state their result, we introduce the following notations. We first define

G = [ wi , v j ] i,j ∈ R m×n . ( 45
)
where { wi } m i=1 is an orthonormal basis of W m . Moreover, we let

δj = n k=1 |x kj |ˆ k-1 j = 1 . . . n, (46)
where xkj are the elements of matrix X appearing in the singular value decomposition of G, that is

G = Ũ Λ XT ( Ũ ∈ R m×m , X ∈ R n×n are orthogonal matrices and Λ ∈ R m×n is the diagonal matrix of singular values {σ j } min(m,n) j=1
). We define the extended set {σ j } n j=1 by using the following convention: if n > m, we set σj = 0 for all j > m. We assume that the singular values {σ j } n j=1 are sorted by decreasing order of magnitude. The upper bound derived in [START_REF] Binev | Data assimilation in reduced modeling[END_REF] then reads:

Theorem 3 ([3, Section 3]) Assume h verifies (18). Then any solution ĥMSU of (44) verifies

h -ĥMSU ≤ 4B { δj } n j=1 , ˆ n . ( 47
)
We note that, due to the close connection between problems ( 16) and (44), the results stated in Theorems 2 and 3 are unsurprinsingly similar although involving different quantities and constants. On the one hand, the quantities appearing in Theorem 2 (e.g., {σ j } n j=1 , {δ j } n j=1 or γ) depends on the basis {w j } m j=1 used to make the measurements {y j } m j=1 . On the other hand, {σ j } n j=1 and { δj } n j=1 used in Theorem 3 are only functions of the subspace W m but not specifically on the generating vectors {w j } m j=1 .9 This is in good accordance with the definitions of the decoders in ( 16) and ( 44) as discussed at the beginning of this section. Despite of their differences, the upper bounds stated in Theorems 2 and 3 can be compared when {w j } m j=1 is an orthonormal basis and ˆ k = k ∀k. In this case, simple calculations show that

σ j = σj , δ j = 2 δj , γ ≤ 1, (48)
and ( 27) particularizes to

h -ĥMS ≤ 2B { δj } n j=1 , ˆ n . (49)
Comparing ( 49) to (47), we then see that our result in Theorem 2 is more favorable than Theorem 3 by a factor 2. However, no conclusion on the relative performance of the multi-space decoders can be drawn from Theorems 2 and 3 since the bounds presented in these results are generally not sharp. In fact, we emphasize herefafter that no general ordering of the performance of decoders ( 16) and ( 44) can be made from a "worst-case" perspective, i.e., by considering the following figure of merit

E( ĥ) sup h∈Mprior∩Py h -ĥ , (50)
where M prior is defined in [START_REF] Binev | Data assimilation in reduced modeling[END_REF]. In particular, we provide two examples of scenarios where either decoder ( 16) outperforms (44) or vice-versa. First, we note that since M prior ∩ P y is a bounded, closed, convex set, E( ĥ) admits a unique minimizer which necessarily belongs to M prior ∩ P y , see [START_REF] Binev | Data assimilation in reduced modeling[END_REF]Lemma 2.1 & Remark 2.4]. We denote this minimizer ĥWC in what follows. From this remark, the scenario where (44) achieves better performance than ( 16) is as follows: since ĥMSU ∈ M prior ∩ P y , it may happen that ĥMSU is (luckily) equal to ĥWC . In such a case, (44) thus attains the best-achievable worst-case performance. On the other hand, if P y ∩V n = ∅, we have necessarily that ĥMS / ∈ M prior ∩ P y and, from the unicity of ĥWC , E( ĥMSU ) < E( ĥMS ). A scenario where the reverse inequality holds is as follows. Assume that P y ∩ V n = ∅. In this case, ĥWC ∈ V n as soon as M prior ∩ P y is symmetric around some point in V n (see [START_REF] Binev | Data assimilation in reduced modeling[END_REF]Remark 2.2]). Moreover, ĥWC is a solution 10of the constrained problem (16) since it satisfies the problem's constraints and leads to a value of the cost function equal to zero. In this case, (16) thus attains the optimal worst-case performance. On the other hand, since (44) may output any point of M prior ∩ P y , we do not have necessarily that ĥMSU = ĥWC . In the latter case, we thus have E( ĥMS ) < E( ĥMSU ).

Proof of Theorem 2

In this section, we provide a proof of the result stated in Theorem 2. We first note that problem (16) is equivalent to finding the minimum of a quadratic where U ∈ R m×m and X ∈ R n×n are the orthonormal matrices appearing in (52). We note that {v * j } n j=1 is an orthonormal basis whereas {w * j } m j=1 is not necessarily orthonormal. By definition, {v * j } n j=1 and {w * j } m j=1 enjoy the following desirable property: Let us define D 1 as

w * i , v * j = σ j if i = j 0 otherwise. ( 55
D 1 =    d = n j=1 β j v * j : n j=1 σ 2 j β 2 j ≤ 4γ 2 2 n    , ( 56
)
where γ is defined in (20). We show hereafter that P Vn (h ) -ĥMS ∈ D 1 .

Let us first consider the intermediate set

S = h : f (h) ≤ γ 2 2 n , (57)
where f (h) m j=1 (y jw j , h ) 2 is the cost function appearing in the variational formulation of multi-space decoder (16).

Clearly P Vn (h ) ∈ S because

f (P Vn (h )) = m j=1 y j -w j , P Vn (h ) 2 = m j=1 w j , h -w j , P Vn (h ) 2 = m j=1 w j , P ⊥ Vn (h ) 2 ≤ γ 2 P ⊥ Vn (h ) 2 ≤ γ 2 2 n . (58)
Moreover, we also have ĥMS ∈ S. This can be seen from the following arguments. First, P Vn (h ) is a feasible point for problem (16), that is

dist(P Vn (h ), V k ) ≤ ˆ k for k = 0 . . . n. (59)
Indeed, rewriting h as

h = n j=1 v j , h v j + z, (60)
where z ∈ V ⊥ n , we have

ˆ k ≥ dist(h , V k ) = P ⊥ V k (h ) = n j=k+1 v j , h v j + z = n j=k+1 v j , h v j 2 + z 2 ≥ n j=k+1 v j , h v j = P ⊥ V k P Vn (h ) = dist(P Vn (h ), V k ). (61)
The first inequality follows from our initial assumption (18). The third equality is true because z ∈ V ⊥ n . Now, since ĥMS is a minimizer of f (h) over the set of feasible points, we have f ( ĥMS ) ≤ f (P Vn (h )) ≤ γ 2 2 n and therefore ĥMS ∈ S. Finally, we show that ĥMS ∈ S and P Vn (h ) ∈ S implies P Vn (h ) -ĥMS ∈ D 1 . Let us first note that, if h ∈ V n , the cost function f (h) can be rewritten as:

f (h) = m j=1 ( w j , h -w j , h ) 2 , = m j=1 w * j , h -w * j , h 2 , = min(m,n) j=1 w * j , h -σ j v * j , h 2 + m j=n+1 w * j , h 2 , (62)
where the second equality follows from the fact that {w j } m j=1 and {w * j } m j=1 differ up to an orthonormal transformation; the last equality is a consequence of (55) and the fact that h ∈ V n by hypothesis.

We note that, since {v * j } n j=1 is an orthonormal basis of V n , P Vn (h ) -ĥMS can be written as

n j=1 β j v * j by setting β j = v * j , P Vn (h ) -v * j , ĥMS .
In particular, this implies that

| v k+1 , h | ≤ ˆ k .
Using the fact that

v * j = n k=1
x kj v k , we obtain (66). In a similar way, we can find that

v * j , P Vn (h ) ≤ η j , (67)
by using the fact that dist(P Vn (h ), V k ) ≤ k from (61).

Let us now show that P Vn (h )-ĥMS ∈ D 2 . Since {v * j } n j=1 is an orthonormal basis of V n , P Vn (h ) -ĥMS can be written as n j=1 β j v * j by setting β j = v * j , P Vn (h )v * j , ĥMS . This leads to

|β j | = v * j , P Vn (h ) -v * j , ĥMS , ≤ v * j , P Vn (h ) + v * j , ĥMS , ≤ ηj + η j = δ j ,
where the last inequality follows from (66) and (67).

If n j=1 σ 2 j δ 2 j ≥ 4γ 2 2 n , the solution of (68) is given by

sup d∈D d 2 = n j= +1 δ 2 j + ρ δ 2 , ( 70
)
where is the largest integer such that

n j= σ 2 j δ 2 j ≥ 4γ 2 2 n , (71)
and ρ ∈ [0, 1] is defined as

ρσ 2 δ 2 + n j= +1 σ 2 j δ 2 j = 4γ 2 2 n . (72)
This can be seen by verifying the optimality condition of problem (68). We note that problem (68) is the same (up to some constants) to the one considered in [START_REF] Binev | Data assimilation in reduced modeling[END_REF]Section 3.1]. The solution (70) is therefore similar, up to some different constants, to the one obtained in that paper.

Conclusions

In this paper, we provide a mathematical characterization of the performance of some particular decoder exploiting a set of approximation subspaces. This decoder was previously shown to lead to good empirical results in several contributions [START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF], although no proof of its theoretical performance was provided in these works. Our result shows how the performance of this "multispace" decoder is related to the parameters defining the approximation problem: the observation W m and prior subspaces {V k } n k=0 , the quality of the prior information (i.e., the widths {ˆ k } n k=0 ) and the distance between the target solution h and the approximation subspaces (i.e., { k } n k=0 ). Based on this result, we show that the "multi-space" decoder can have provably better reconstruction guarantees than its standard ("single-space") counterpart in some situations.

A Proof of (42)

In this appendix, we show that the cost function f (h) n j=1 (y j -w j , h ) 2 can be rewritten as in (42) when h ∈ Vn and σn > 0. 11 First, using the definition of y j , we have

f (h) = n j=1 ( w j , h -w j , h ) 2 .
(73)

Moreover, using the particular bases introduced in Section 5.1.1, we obtain

f (h) = n j=1 w * j , h -w * j , h 2 , = n j=1 w * j , h -σ j v * j , h 2 , (74)
where the first equality follows from the fact that {w j } n j=1 and {w * j } n j=1 differ up to an orthogonal transformation; the second is a consequence of (55) and our hypothesis h ∈ Vn.

Since ĥSS corresponds to the minimum of f (h) over Vn (see ( 17)), we simply have

v * j , ĥSS = w * j , h σ j , (75)
if σn > 0. Hence, under this assumption, (74) can also be rewritten as in (42).

Fig. 1

 1 Fig. 1 Graphical representation of the geometry of the multi-space problem (16).

) 5 . 1 . 2

 512 Definition of D 1

5. 2 2 Weσ 2 j δ 2 j < 4γ 2 2 n

 2222 Expression of sup d∈D d , the first constraint in (68) is always inactive and the solution simply reads sup d∈D

INRIA, Centre Rennes Bretagne-Atlantique, Campus de Beaulieu, 35000 Rennes, France

In this paper, we assume that H is a Hilbert space with inner product •, • and induced norm • .

We give a precise formulation of the decoders considered in[START_REF] Binev | Data assimilation in reduced modeling[END_REF] and[START_REF] Argaud | Stabilization of (G)EIM in presence of measurement noise: application to nuclear reactor physics[END_REF][START_REF] Fick | A Reduced Basis Technique for Long-Time Unsteady Turbulent Flows[END_REF][START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF][START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF] in Section 4 of this paper.

We note that the main result of this paper (see Theorem 2 in Section

[START_REF] Chaturantabut | Nonlinear model reduction via discrete empirical interpolation[END_REF] was presented in the scientific conferences Enumath'17[START_REF] Herzet | Beyond petrov-galerkin projection by using multi-space prior[END_REF] and Morepas'18[START_REF] Herzet | Beyond Petrov-Galerkin projection by using multispace prior[END_REF] but not published elsewhere with a complete description of the proofs.

In this example, w j models the fact that the state of the system can be sensed (with some precision σ 2 ) around m j ∈ Ω.

In this paper we assume that the constraints are defined ∀k ∈ {0 . . . n}. All the derivations presented hereafter may nevertheless be easily extended to the case where the constraints in (16) are only available for some k ∈ {0 . . . n}.

If { k = dist(h , V k )} n k=1is known, the result provides an upper bound on the reconstruction error for a particular h . If not, setting k = ˆ k provides an "a priori" upper bound on the reconstruction error of any h verifying (18).

In particular, all the constraints in (16) become inactive.

We note that σ min = σn because {w j } m j=1 is an orthonormal basis.

Indeed, we note that the singular values {σ j } n j=1 and matrix X only depend on the subspace Wm but not on the orthonormal basis { wi } m i=1 used to define it.

It is in fact the unique solution as soon as n ≤ m and the vectors {w j } m j=1 are linearly independent.

The authors thank the "Agence nationale de la recherche" for its financial support through the Geronimo project (ANR-13-JS03-0002).

function over a closed bounded subset of V n . A minimizer thus always exists. Moreover, the unicity of the minimizer stated at the end of Theorem 2 follows from the strict convexity of the cost function over V n when σ n > 0.

In the rest of this section, we thus mainly focus on the derivation of the upper bound (27). Our proof is based on the following steps. First, since ĥMS ∈ V n , we have that

where P Vn (•) (resp. P ⊥ Vn (•)) denotes the orthogonal projector onto V n (resp. V ⊥ n ). We then derive an upper bound on P Vn (h ) -ĥMS 2 as follows:

• We identify a set D such that P Vn (h ) -ĥMS ∈ D in Section 5.1. This implies in particular that P Vn (h ) -ĥMS 2 ≤ sup d∈D d 2 .

• We derive the analytical expression of sup d∈D d 2 as a function of the parameters

Combining these results, we obtain (27).

Definition of D

We express D as the intersection of two sets D 1 and D 2 that we define in Sections 5.1.2 and 5.1.3 respectively. In order to properly define these quantities, we introduce some particular bases for V n and W m in Section 5.1.1.

Some particular bases for V n and W m

Let

be the singular value decomposition of the Gram matrix defined in (22), where U ∈ R m×m and X ∈ R n×n are orthonormal matrices and Λ ∈ R m×n is the diagonal matrix of singular values. We denote by {σ j } n j=1 the set of singular values of G sorted in their decreasing order of magnitude. We remind the reader that if m < n, we adopt the convention σ j = 0 ∀j > m.

We define the following bases for V n and W m :

Therefore, we have

where the first inequality follows from the standard inequality (a + b) 2 ≤ 2(a 2 + b 2 ), the second from (62), and the last one from the fact that ĥMS ∈ S and P Vn (h ) ∈ S.

Definition of D 2

Let

where

and the x ij 's are the elements of the matrix X appearing in the SVD decomposition (52). We define D 2 as

We show hereafter that P Vn (h ) -ĥMS ∈ D 2 .

We first note that if h is feasible for problem (16), we must have

Indeed, if h is feasible, the constraint dist(h, V k ) ≤ ˆ k simply writes as


Aldo Auditore

Matteo Conforti

Costantino De Angelis

Alejandro B Aceves

Dark-antidark solitons in waveguide arrays

Keywords:

I. INTRODUCTION

In the last years, discrete optical systems and waveguide arrays have been a very active research area in optics [START_REF] Lederer | Discrete solitons in optics[END_REF]. Binary waveguide arrays in particular have been studied because their intrinsic two bands structure can be very helpful in order to control wave propagation in the linear and nonlinear regimes [START_REF] Longhi | Multiband diffraction and refraction control in binary arrays of periodically curved waveguides[END_REF][START_REF] Guasoni | Peculiar properties of photonic crystal binary waveguide array[END_REF][START_REF] Efremidis | Wave propagation in waveguide arrays with alternating positive and negative couplings[END_REF]. More recently, the interplay between plasmonic waveguiding and periodicity has been also considered, inasmuch as plasmonic confinement offers an extra degree of freedom to be usefully exploited in all-optical devices [START_REF] Conforti | Subwavelength diffraction management[END_REF][START_REF] Guasoni | Light propagation in nonuniform plasmonic subwavelength waveguide arrays[END_REF][START_REF] Nam | Deep subwavelength surface modes in metal-dielectric metamaterials[END_REF].

In this framework solitons represent an important class of solutions, as their particle like behaviour can be very useful for switching applications and their peculiar features often represent an unvaluable tool to understand the overall dynamics of the system in the nonlinear regime. This certainly explains the huge effort that the scientific community has put in finding soliton solutions in these systems [START_REF] Efremidis | Wave propagation in waveguide arrays with alternating positive and negative couplings[END_REF][START_REF] Kivshar | Optical solitons: From Fibers to Photonic Crystals[END_REF][START_REF] Sukhorukov | Discrete gap solitons in modulated waveguide arrays[END_REF][START_REF] Morandotti | Observation of discrete gap solitons in binary waveguide arrays[END_REF][START_REF] Sukhorukov | Soliton control and Bloch-wave filtering in periodic photonic lattices[END_REF][START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF][START_REF] Conforti | Modulational stability and gap solitons of gapless systems: Continuous versus discrete limits[END_REF].

It is well known that in nonlinear discrete systems as those describing light propagation in waveguide arrays, bright localized modes may exist in the form of gap solitons when a gap opens in the linear dispersion relation [START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF][START_REF] Marini | Coupledmode approach to surface plasmon polaritons in nonlinear periodic structures[END_REF]. On the other hand, without the band gap in the linear dispersion relation, a forbidden frequency region can still exist for a nonuniform nonlinear response and solitons sitting on a pedestal can be found [START_REF] Conforti | Modulational stability and gap solitons of gapless systems: Continuous versus discrete limits[END_REF][START_REF] Kivshar | Self-induced gap solitons[END_REF][START_REF] Flytzanis | Gap solitons in gapless coupled chains[END_REF]. Among the situations where the nouniform nonlinearity can be exploited it is certainly worth quoting the case of the linear-nonlinear interlaced waveguide arrays [18].

In previous papers we have used a continuous approximation to exploit bright solitary wave solutions of this system [START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF] and we have then used a fully discrete model to explore the existence and stability of solitons sitting on a non zero background [START_REF] Conforti | Modulational stability and gap solitons of gapless systems: Continuous versus discrete limits[END_REF]. In this paper we extend the continuous approximation to the case of a non zero background and we demonstrate that this continuous approximation can capture many features of the discrete system; * Electronic address: aldo.auditore@ing.unibs.it remarkably we also show that the continuum approximation gives a reasonable description of discrete states even when they are confined to a very small number of sites.

II. PHYSICAL SETTINGS AND THEORETICAL ANALYSIS

According to coupled mode theory and taking into account third-order nonlinearities in the form of a pure Kerr effect, the field amplitude propagation in a binary waveguide array can be described by the following two sets of coupled equations with constant coefficients: [START_REF] Sukhorukov | Discrete gap solitons in modulated waveguide arrays[END_REF][START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF]

iA nz + ∆β 2 A n + C 1 B n + C 2 B n+1 + γ 1 |A n | 2 A n = 0 iB nz -∆β 2 B n + C 2 A n-1 + C 1 A n + γ 2 |B n | 2 B n = 0(1)
where A n , B n are respectively the mode amplitudes in the even and odd waveguides, ∆β represents the difference between the propagation constants in even and odd waveguides, γ 1 (γ 2 ) is the site-dependent nonlinear coefficient for n even (odd), C 1 and C 2 are the coupling coefficients when n is odd and even respectively. Finally, without loss of generality we can set C 2 = 1.

In close proximity of the band edge (i.e. around k x = 0 for C 1 < 0 and for k x around π for C 1 > 0) a very useful equivalent continuous model can be derived by performing a Taylor expansion to obtain (as a first order approximation):

iu z + ∆β 2 u + w x + ǫw + γ 1 |u| 2 u = 0 iw z - ∆β 2 w -u x + ǫu + γ 2 |w| 2 w = 0 (2)
where we have also defined C 1 = ±1 + ǫ, with the + (-) sign that has to be used for k x = 0 (k x = π). To look for both stationary and walking self confined solutions of the system defined by Eqs. (2), we use the following trial functions [START_REF] Aceves | Self-induced transparency solitons in nonlinear refractive periodic media[END_REF]:

u(x, z) = 1 2 (K 1 g 1 (ξ) + iK 2 g 2 (ξ)) exp(i cos(Q)ψ) w(x, z) = 1 2i (K 1 g 1 (ξ) -iK 2 g 2 (ξ)) exp(i cos(Q)ψ) ξ = x + vz √ 1 -v 2 , ψ = vx + z √ 1 -v 2 K 1 = ( 1 + v 1 -v ) 1/4 , K 2 = ( 1 -v 1 + v ) 1/4 (3)
with g 1,2 two arbitrary complex functions, -1 < v < 1.

Although not necessary, for the sake of clarity, from now on we set ∆β = 0. Substituting the ansatz (3) into Eqs. ( 2) and following the procedure as in [START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF], it is straightforward to obtain an Hamiltonian form for the equations, observing that

P = |g 1 | 2 -|g 2 | 2 is a constant of motion for the dynami- cal system. Indeed, setting g 1,2 (ξ) = f 1,2 (ξ)exp[iθ 1,2 (ξ)], η = f 2
2 and µ = θ 1 -θ 2 , η and µ obey the following Hamiltonian system:

η = - ∂H ∂µ μ = ∂H ∂η H = 2η cos Q + 2 η(η + P )ǫ cos µ - s 8 η η K 4 1 2 + K 4 2 2 + 2 -cos(2µ) + P K 4 1 + 2 -cos(2µ) - d 4 η(η + P ) η(K 2 1 + K 2 2 ) + P K 2 1 sin µ. ( 4
)
where we also set

s = γ 1 + γ 2 and d = γ 1 -γ 2 .
It is straightforward to show that this Hamiltonian has the following symmetry: mapping P and v into -P and -v induces only a non nonessential shift by a constant into the Hamiltonian's value. For this reason, from now on we set P ≥ 0. Note also that the Hamiltonian system described by Eqs.(4) reduces obviously to the one considered in [START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF] for P = 0. Moreover in the situation considered in [START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF] the condition ǫ = 0 was necessary in the quest for bright soliton solutions; here, on the contrary, we are interested in discussing solitons with a non zero background and their existence is not related to the presence of a band gap in the linear spectrum, i.e. they exist even in the case ǫ = 0 as we have discussed in [START_REF] Conforti | Modulational stability and gap solitons of gapless systems: Continuous versus discrete limits[END_REF]. As a matter of fact, the key properties of these solutions sitting on a non zero background do not depend on the presence of a bandgap in the linear spectrum; for these reasons from now on we consider the case ǫ = 0, corresponding to the presence of a Dirac point at zero transverse momentum in the linear spectrum. Soliton solutions of Eqs. 2 correspond to separatrix trajectories emanating from and sinking into unstable fixed points of the dynamical system described by Eqs. 4. Obviously the validity of the continuum approximation becomes more questionable as the degree of localization of the solitons increases; in order to get a feeling of how far one can pushes the use of the continuum model while still having a reasonable description of the discrete system, we performed a thorough comparison between the approximated soliton solutions obtained in the continuum limit and exact soliton states obtained numerically using the Newton conjugate-gradient method [START_REF] Yang | Newton-conjugate-gradient methods for solitary wave computations[END_REF].

To summarize our findings we report in figure 1 the comparison between the modulus of the approximated soliton states obtained from the continuum limit (blue crosses and dashed line) and exact numerical results obtained using the Newton conjugate-gradient method (red open circles). Note that despite the very strong degree of localization (the soliton structure is basically confined to 5 sites only), the continuum approximation still shows an excellent agreement with the exact solution and this holds true both for the modulus (reported in figure 1) and the phase (not shown here). If we further increase the degree of confinement we gradually loose the validity of the continuum approximation and for soliton states confined to 3 sites only the continuum approximation is not able to capture closely the features of the discrete solutions. For such an high level of confinement one has to go back to the discrete model and find there solutions as was recently done by using asymptotic expansions in [START_REF] Conforti | Modulational stability and gap solitons of gapless systems: Continuous versus discrete limits[END_REF]. It is remarkable to note that by using these two different approximations (asymptotic expansion and continuum model) one can describe in a simple and accurate fashion the entire spectrum of the dark-antidark soliton states of this system. The derived hamiltonian system thus represents a valuable tool in describing and in understanding the features of the solutions of the problem; the goal of the rest of the paper is to prove that thanks to this hamiltonian system we are able to introduce new soliton solutions for the problem at hand. Moreover, we will always test the validity of our results, obtained in the framework of a continuum approximation, in the original (truly discrete) system. The above goals will be pursued using some representative examples; we will thus focus our attention on three different binary arrays: the first case we will discuss is the case of uniform nonlinearity in the array; the second case we will face is that of the linear-nonlinear interlaced binary array and the third case will be the case of nonlinearities with different signs along the array, i.e. the focusing-defocusing interlaced binary waveguide array.

III. EXAMPLES

In this section we present some results derived from the analysis of the Hamiltonian system (Eqs. 4). We thus first look for unstable fixed points of the dynamical system and then obtain the separatrices corresponding to solitary wave solutions. To assess the validity of our approach we then propagate the obtained waveforms in the truly discrete system described by Eqs.1. If we were looking for bright solitons, as we did in [START_REF] Conforti | Energy localization and transport in binary waveguide arrays[END_REF], we would pick P = 0 and this in turn would simplify considerably the algebra required in the analysis of the dynamical system described by Eqs.4; here, on the contrary, we can not limit ourselves to the simplest case P = 0 and thus we focus on the more general case P = 0.

As a first example, we consider an array with uniform nonlinearity (i.e. d = 0). In figure 2a we report the bifurcation diagram of the amplitude η of the fixed points as a function of Q; unstable fixed points and thus dark solitons do exist only for Q < arccos sP (3 + k 2 1 )/16 ; note that the unstable fixed points here correspond to two different branches with different generalized phase (µ = ±π/2). In figure 2b we report the phase plane analysis of the system and we see there two saddle points at µ = ±π/2 and η = (16 cos(Q) -sP (K 4 1 + 3))/(s(K 4 1 + K 4 2 + 6)); the separatrices emanating from and sinking into the saddles turn around the center located at µ = 0 and η = (16 cos Q -sP (K 4 1 + 1))/(s(K 4 1 + K 4 2 + 2)). The existence of these saddle points is possible if the following constraints on Q and P are satisfied: 0 < Q < arccos(sP (K 4 1 + 3)/16) and 0 < P < 16/(s(K 4 1 + 3)). To test the validity of these results we use the obtained solutions as the initial condition (i.e. at z=0) in the numerical integration of the equations describing propagation in the binary waveguide array (i.e. Eqs. 1); the results are reported in figure 3a where the initial condition corresponds to the trajectory enlightened with a thicker line in the inset of the figure. Note that propagation of this solution into the array reveals only a small amount of radiation (almost non visible) thus proving the validity of the continuum approximation; moreover these solutions are strong enough to survive for non zero transverse velocities; in figure 3b we report one example of non zero transverse velocity and we observe once again that the continuum system is able to capture the nature of the discrete one. Obviously the solution cannot move without radiating some energy; however the amount of radiation is almost negligible and the soliton can propagate with high velocity along the array. As a second example we consider the case of an array of alternating linear-nonlinear waveguides, i.e. s = d. The situation is quite different with respect to the d = 0 case; two fixed points (one center and one saddle) exist for Q < π/2; only one stable center exists for Q > π/2 (see figure 4a). In figure 5 we report the beam propagation along the array 6a, for Q > π/2 we find one center and one saddle. In the corresponding phase plane in figure 6b we report the heteroclinic trajectory emanating from and sinking into the saddle points at µ = -π/2 and µ = 2π -π/2. In figure 7 we show the field evolution along the array using as initial condition the field profiles obtained from Eqs. 3 and corresponding to the separatrix enlightened in figure 6b. Once again we can observe that propagation of this solution into the array reveals only an almost invisible small amount of radiation thus proving the validity of our analytical approach.

CONCLUSIONS

In this work we have obtained dark and antidark soliton solutions in a binary waveguide array with alternating positive and negative linear couplings between adjacent waveguides and in the presence of focusing and/or defocusing Kerr nonlinearity. These solutions do exist also in a linear-nonlinear interlaced array and they even survive in focusing-defocusing interlaced arrays. We have also numerically verified the soundness of our approach by a detailed comparison with exact results obtained by numerically solving the discrete system; remarkably our results, obtained in the framework of a continuum approximation, retain their validity also for very strong degrees of localization.

FIG. 1 :

 1 FIG. 1: Comparison between approximated soliton states obtained from the continuum limit (blue crosses and dashed line) and exact numerical results obtained using the Newton conjugate-gradient method (red open circles). With reference to Eqs. 1 here we set γ1 = 1, γ2 = 0, C1 = -1, C2 = 1, ∆β = 0 and limn→±∞ |An| = √ 2. a): even sites; b): odd sites.

FIG. 2 :

 2 FIG. 2: s = 2, d = 0, P = 1, v = 0. a) Bifurcation diagram of the Hamiltonian system: continuous lines for stable centers, dashed lines for unstable saddle points. b) Phase plane analysis (Q = 1): crosses correspond to the unstable points; the dot shows the stable center. The thicker lines are the heteroclinic separatrices connecting the two unstable saddles and correspond to solitons on a non zero background.

FIG. 3 :

 3 FIG. 3: Field evolution along the array for s = 2, d = 0,P = 1: the initial condition is the dark soliton solution corresponding to the thicker line in the inset where the phase plane of system 4 is reported. a): Q = 1, v = 0; b): Q = 0.7, v = 0.5

FIG. 4 :FIG. 5 :

 45 FIG. 4: s = 2, d = 2, P = 0.6, v = 0. a) Bifurcation diagram of the Hamiltonian system: continuous lines for stable centers, dashed lines for unstable saddle points. b) Phase plane analysis (Q = 0.6): the cross corresponds to the unstable point; the dot shows the stable center. The thicker line is the homoclinic separatrix corresponding to solitons on a non zero background.

FIG. 6 :

 6 FIG. 6: s = 2, d = 2.1, P = 1.0, v = 0. a) Bifurcation diagram of the Hamiltonian system: continuous lines for stable centers, dashed lines for unstable saddle points. b) Phase plane analysis (Q = 1.7): the crosses correspond to the unstable points; the dots show the stable centers. The thicker line is the heteroclinic separatrix corresponding to solitons on a non zero background.

ACKNOWLEDGMENTS

The authors wish to thank Professor T. R. Akylas for very helpful discussions. AA, CDA and MC acknowl-edge financial support from CARIPLO Foundation under grant no. 2010-0595 and US ARMY under grant no. W911NF-12-1-0202.