Pooneh Maghoul
email: pooneh.maghoul@gmail.com

Behrouz Gatmiri
email: behrouz.gatmiri@andra.fr

Murray Fredlund
email: murray@soilvision.com

Denis Duhamel
email: denis.duhamel@enpc.fr

A Boundary Element Formulation for the Wave Propagation in the Unsaturated Soils

This paper aims at obtaining an advanced formulation of the time-domain Boundary Element Method (BEM) for two-dimensional dynamic analysis of unsaturated soil. Unlike the usual time-domain BEM the present formulation applies a Convolution Quadrature which requires only the Laplace-domain instead of the time-domain fundamental solutions. The coupled equations governing the dynamic behavior of unsaturated soils ignoring contributions of the inertia effects of the fluids (water and air) are derived based on the poromechanics theory within the framework of the suction-based mathematical model. In this formulation, the solid skeleton displacements, water pressure and air pressure are presumed to be independent variables. As there is no analytical solution for the 2D wave propagation in unsaturated soils in the literature, to verify the accuracy of this implementation, the displacement response obtained by the boundary element formulation is partially verified by comparison with the elastodynamics problem.

INTRODUCTION

It has long been recognized that surface topographies can have crucial influences on damage severity and its spatial distribution during strong earthquakes. Site response analysis of topographical structures could only be solved accurately, economically and under realistic conditions, with the aid of numerical methods such as the Boundary Element Method (BEM).

Unsaturated soils are encountered near the earth's surface where most engineering structures are ultimately supported. Saturated and dry soils can become unsaturated due to seasonal variations. In the current state of the art, it could be claimed that behavior of the saturated porous media has been well understood. Conversely, the study of the dynamic behavior of unsaturated porous media is a relatively new area in the field of geotechnical earthquake engineering.

The BEM is a very effective numerical tool for dynamic analysis of linear elastic bounded and unbounded media. The method is very attractive for wave propagation problems, because the discretization is done only on the boundary, yielding smaller meshes and systems of equations. Another advantage is that this method represents efficiently the outgoing waves through infinite domains, which is very useful when dealing with scattered waves by topographical structures.

This paper aims at obtaining an advanced formulation of the time-domain Boundary Element Method (BEM) for two-dimensional dynamic analysis of unsaturated soil. Unlike the usual time-domain BEM the present formulation applies a Convolution Quadrature developed by [START_REF] Lubich | Convolution quadrature and discretized operational calculus, I[END_REF] which requires only the Laplace-domain instead of the time-domain fundamental solutions.

In this paper first of all, a set of fully coupled governing differential equations of a porous medium saturated by two compressible fluids (water and air) subjected to dynamic loadings is obtained. These phenomenal formulations are presented based on the experimental observations and with respect to the poromechanics theory within the framework of the suction-based mathematical model presented by [START_REF] Gatmiri | UDAM: a powerful finite element software for the analysis of unsaturated porous media[END_REF]. In this model, the effect of deformations on the suction distribution in the soil skeleton and the inverse effect are included in the formulation via a suctiondependent formulation of state surfaces of void ratio and degree of saturation. In this formulation, the solid skeleton displacements , water pressure and air pressure are presumed to be independent variables. Secondly, the Boundary Integral Equation (BIE) is developed directly from those equations via the use of the weighted residuals method in a way that permits an easy discretization and implementation in a numerical code. The associated fundamental solution obtained by [START_REF] Maghoul | Boundary integral formulation and two-dimensional fundamental solutions for dynamic behaviour analysis of unsaturated soils[END_REF] is used in the BIE.

Finally, as it seems that there is no analytical solution for the 2D wave propagation in unsaturated soils in the literature, to verify the accuracy of this implementation, the displacement response obtained by the boundary element formulation is partially verified by comparison with the elastodynamic problem.

GOVERNING EQUATIONS

Governing differential equations consist of mass conservation equations of liquid and gaseous phases, the equilibrium equation of the skeleton associated with water and air flow equations and constitutive relation. The assumption of infinitesimal transformation and incompressibility of solid matrix is considered. Solid Skeleton. The equilibrium equation and the constitutive law for the soil's solid skeleton including the effect of suction are written:

( ) (1)
( ) ( ) (2)
where are Lame coefficients, = is the water or air pressure, is the Kronecker delta and is the suction modulus matrix:

(3) in which is a vector obtained from the state surface of void ratio which is a function of the independent variables of and :

( ) [ ] (4)
The elasticity matrix ( ) can be presented by using the bulk modulus and the tangent modulus:

(5)

where is tangent elastic modulus which can be evaluated as: (6) in which is the elastic modulus in absence of suction and ( 7) being a constant, represents the effect of suction on the elastic modulus. is the bulk modulus of an open system and evaluated from the surface state of void ratio:

( ) (8)
Mass Conservation of Water. The mass of the water in a representative elementary volume can be written as:

where and . In this equation, = is the displacement of water or air relative to solid, = is the degree of saturation relative to water or air, = is the compressibility of water or air = and . Mass Conservation of Air. With the same approach presented before, the mass conservation of the air can be written as:

(10) where and . Flow Equation for the Water. Based on generalized Darcy's law for describing the balance of the forces acting on the liquid phase of the representative elementary volume, the water velocity in the unsaturated soil takes the following form: (11) in which denotes the water permeability in an unsaturated soil. Flow Equation for the Air. With the same approach presented for the water based on generalized Darcy's law, the air velocity in the unsaturated soil takes the following form: (12) in which denotes the air permeability in an unsaturated soil. Summary of the Governing Differential Equations in Laplace Domain. By introducing (2) into (1), ( 11) into ( 9) and ( 12) into (10) and by applying the Laplace transform to eliminate the time variable of partial differential equations, assuming = = = and = = , we obtain the final set of governing equations in Laplace transform domain:

[ ̃ ̃ ̃] [ ̃ ] (13)
with the not self-adjoint operator :

[ ] (14)
where and . In equations ( 13) and ( 14), ̅̅̅̅ in two dimensional problems. Also in ( 14), the partial derivative , is denoted by and is the Laplacian operator.

BOUNDARY INTEGRAL EQUATION

The boundary integral equations for this problem will be derived by taking the fundamental solution as the weighted function and using the method of weighted residuals, which is essentially an integration by parts technique. In this method, the integral equation is derived directly by equating the inner product of eq. ( 13) and the matrix of the adjoint fundamental solutions ̃ implying that ̃ (15)

to a null vector, i.e. ∫ Ω [ ̃ ̃ ̃] ̃ Ω with ̃ [ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ] [ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ] ( 16
)
where the integration is performed over a domain Ω with boundary and vanishing body forces and sources are assumed. In eq. (15), is the adjoint operator to used for the deduction of fundamental solutions [START_REF] Maghoul | Boundary integral formulation and two-dimensional fundamental solutions for dynamic behaviour analysis of unsaturated soils[END_REF].

After integrating by parts twice over the domain according to the theory of Green's formula and using partial integration, the operator is transformed from acting on the vector of unknowns [ ̃ ̃ ̃] to the matrix of fundamental solutions ̃. After some algebraic manipulations and using the property of Dirac's delta function, we reach the transformed dynamic unsaturated poroelastic boundary integral representation for the transformed internal displacements and pressures given in matrix form, i.e.,

[ ] [ ̃ ̃ ̃ ] ∫ [ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ] [ ̃ ̃ ̃ ] ∫ [ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ̃ ] [ ̃ ̃ ̃ ] (17)
where the traction vector, the normal water flux and the normal air flux are respectively:

̃ [( ̃ ̃ ̃ ̃) ( ̃ ̃ )] (18) ̃ ( ̃ ̃ ) (19) ̃ ( ̃ ̃ ) (20)
The coefficient has a value for points inside Ω and zero outside Ω. The value of for points on the boundary is determined from the Cauchy principal value of the integrals.

The corresponding 2D fundamental solutions ̃ as well as the ̃ , ̃ and ̃ in Eq. ( 17) which can be interpreted as the adjoint terms to the traction vector ̃ , the water flux ̃ and the air flux ̃ were already derived in Maghoul et el. (2011).

The time dependent boundary integral equation for the unsaturated soil is obtained by a transformation to time domain:

∫ ∫[ ] (21)

BOUNDARY ELEMENT FORMULATION

Eq. ( 21) is an exact represent of the dynamic response of a multiphase porous medium, involving integrations over the surface as well as the time history. For the practical problem, suitable approximations are needed for both spatial and temporal variations of field variables. As will be shown, temporal integrations of the time functions involved will be performed numerically using an operational convolution quadrature method (CQM), as like as the spatial integration which be evaluated using numerical techniques. The salient features of the temporal and spatial integrations are outlined below.

Temporal integration. Since it seems too difficult to obtain the time-dependent fundamental solution in an explicit analytical form by an inverse transformation of the frequency domain results, the convolution quadrature method (CQM) proposed by Lubich (1988 a, b) is used. In this formulation, the convolution integral is numerically approximated by a quadrature formula whose weights are determined by the Laplace transform of the fundamental solution and a linear multistep method [START_REF] Schanz | A new visco-and elastodynamic time domain boundary element formulation[END_REF]. By applying this method, the convolution integrals between the fundamental solutions and the nodal values in eq. ( 21) are approximated by:

∫ ∑ ( ) = (22) ∫ ∑ ( ) = (23)
in which ( ) and ( ) are the influence function which are defined by:

( ) ℛ ∑ ̃ = (24) ( ) ℛ ∑ ̃ = (25)
where is given by (ℛ ) . By substituting Eqs. ( 22) and (23) into Eq. ( 21), the time-convoluted boundary element equation is:

∑ ∫( ) = ∫( ) (26)
Spatial integration. In this section in order to numerically evaluate the integrals on the surface Eq. , a spatial discretization will be performed. Using isoparametric quadratic elements and assuming a quadratic variation over both geometry and field variables, the functions (displacements and tractions) at any point over an element can be expressed in terms of the nodal values as , , , , and where for 2D and for a quadratic element and are the shape functions in the local intrinsic coordinates of the element. Once the spatial discretization process described above has been accomplished, the nodal quantities can be brought outside the surface integrals of Eq. ( 26), since now the integrands contain only known functions. Therefore, the discretized BE equation corresponding to point can be written as:

∑ ∑ ∑ [( ) = = = = ( ) ] (27)
in which the integrals which have to be evaluated over the isoparametric element , can be written in intrinsic coordinates. Then:

∫ ( ) | | (28) ∫ ( ) | | (29)
where | | is the Jacobian of transformation. The usual point collocations scheme, i.e. by allowing point to coincide sequentially with all the nodal points of the boundary, is used to establish a set of integral equations in order to obtain unknown boundary values. Since there are four degrees of freedom for each node of the problem [ ] , there are integral equations per location . Eq. ( 27) can be written in matrix form by replacing the double sum Σ = = Σ = = by a matrix multiplication, where coefficients are assembled in a similar way as in the FEM. Then,

̅ ∑ [ ] = (30)

VERIFICATION

As explained above, this research is the first attempt to develop BEM formulations for the problem of 2D wave propagation in unsaturated soils. The complexity of the coupled hydro-mechanical model in unsaturated porous media, subjected to dynamic loadings, leads to the fact that no analytical solution is known in the literature.

Therefore, the results obtained by BEM for the dynamic behavior of unsaturated soils will be partially verified by comparison with elastodynamic result which is related to a purely mechanical simulation. In this case, the nodal degrees of freedom related to water pressure and air pressure are blocked.

A 3m high soil column, meshed by 32 quadratic elements on the boundary and a total of 64 nodes is considered (Fig. 1). At the upper surface of the column a sudden and constant longitudinal pressure, which is distributed uniformly is applied. On the side walls the movements in the horizontal direction are restricted, while at the bottom surface they are blocked in the vertical direction. In addition, all the boundaries are assumed to be permeable, i.e. the water pressure is zero. Also, the air pressure is assumed to be zero. The hydraulic and mechanical properties of the soil are presented in the table 1.

The result of the numerical calculation by BEM at the central point of the upper surface is compared to elastodynamic results (Fig. 2). There is a good agreement between the elastodynamic solution and the numerical results obtained by the proposed formulation.

Figure 1 .

 1 Figure 1. Unsaturated column of soil.

  , an advanced formulation of the time-domain Boundary Element Method (BEM) for two-dimensional dynamic analysis of unsaturated soil is obtained. Unlike the usual time-domain BEM the present formulation applies a Convolution Quadrature which requires only the Laplace-domain instead of the time-domain fundamental solutions.

Figure 2 .

 2 Figure 2. Vertical displacements obtained at the central point at the upper surface

Table 1 . Hydraulic and mechanical properties of unsaturated column

 1

	. ×	. ×	.	.		.	. ×
	. ×	.	.	. ×	.	×


Gabriel Raúl Barrenechea

AND Lyonell Boulton
email: l.boulton@hw.ac.uk

Nabile Boussaid
email: nboussai@univ-fcomte.fr

enclosures. 2013 Nabile Boussaid Eigenvalue

Hal- 00837475v3

Gabriel R Barrenechea
email: gabriel.barrenechea@strath.ac.uk

Nabile Boussa Ïd

Eigenvalue enclosures

Keywords: eigenvalue enclosures, spectral pollution, finite element method, Maxwell equation

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

 [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF], and a method developed more recently by Davies and Plum [START_REF] Davies | Spectral pollution[END_REF]. These turn out to be highly robust and they can be applied to a wide variety of settings with minimal implementation difficulties.

The approach of Zimmermann and Mertins is based on an extension of the Lehmann-Maehly-Goerisch method [START_REF] Goerisch | The convergence of a new method for calculating lower bounds to eigenvalues[END_REF][START_REF] Weinberger | Variational Methos for Eigenvalue Approximation[END_REF] and it has proved to be highly successful in concrete numerical implementations. These include the computation of bounds for eigenvalues of the radially reduced magnetohydrodynamics operator [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF][START_REF] Boulton | Eigenvalue enclosures for the MHD operator[END_REF], the study of complementary eigenvalue bounds for the Helmholtz equation [START_REF] Behnke | Bounds for eigenvalues with the use of finite elements[END_REF] and the calculation of sloshing frequencies in the left definite case [START_REF] Behnke | Lower and upper bounds for sloshing frequencies[END_REF].

The method of Davies and Plum on the other hand, is based on a notion of approximated spectral distance which is highly geometrical in character. Its original formulation dates back to [START_REF] Davies | Spectral enclosures and complex resonances for general self-adjoint operators[END_REF][START_REF]A hierarchical method for obtaining eigenvalue enclosures[END_REF][START_REF] Davies | Spectral pollution[END_REF], and it is yet to be tested properly on models of dimension other than one. Our main motivation for the analysis conducted below, initiated with the results presented in [START_REF] Davies | Spectral pollution[END_REF]Section 6] where it is shown that both techniques are equivalent. Below we determine in a more precise manner the nature of this equivalence and examine their convergence properties.

In Section 2 we extend various canonical results from [START_REF] Davies | Spectral pollution[END_REF]. Notably, we include multiplicity counting (propositions 1 and 4) and a description of how eigenfunctions are approximated (Proposition 6). The method of Zimmermann and Mertins, on the other hand, is introduced in Section 3. We derive the latter in a self-contained manner independently from the work [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF]. See Theorem 9 and Corollary 10.

Section 4 addresses the questions of convergence and upper bounds for residuals in both methods. The main statements in this respect are Theorem 13, Corollary 14 and Theorem 15, where we formulate general convergence estimates with explicit bounds for a finite group of contiguous eigenvalues.

Section 5 is devoted to a concrete computational application in the spectral pollution regime. For this purpose, we consider the model of the resonant cavity, for which it has been well-documented that nodal elements lead to spurious eigenvalues. Remarkably the present approach on nodal elements allows estimation of sharp eigenvalue bounds. A companion Comsol Multiphysics v4.3b Livelink code which was employed to produce some of the results presented in Section 5 as well as further numerical experiments on this model, is available in the appendix.

Approximated local counting functions

Let A : D(A) -→ H be a self-adjoint operator acting on a Hilbert space H. Decompose the spectrum of A in the usual fashion, as the union of discrete and essential spectrum, σ(A) = σ disc (A) ∪ σ ess (A). Let J be any Borel subset of R. The spectral projector associated to A is denoted by 1 J (A) = J dE λ . Hence Tr 1 J (A) = dim 1 J (A)H. We write E J (A) = ⊕ λ∈J ker(A -λ) with the convention E λ (A) = E {λ} (A). Generally E J (A) ⊆ 1 J (A)H, however there is no reason for these two subspaces to be equal.

Let t ∈ R. Let q t : D(A) × D(A) -→ C be the closed bi-linear form

(1) q t (u, w) = (A -t)u, (A -t)w ∀u, w ∈ D(A).

For any u ∈ D(A) we will constantly make use of the following t-dependant seminorm, which is a norm if t is not an eigenvalue,

(2) |u| t = q t (u, u) 1/2 = (A -t)u .

By virtue of the min-max principle, q t characterizes the spectrum which lies near the origin of the positive operator (A -t) 2 . In turn, this gives rise to a notion of local counting function at t for the spectrum of A as we will see next. Let

d j (t) = inf dim V =j V ⊂D(A) sup u∈V |u| t u so that 0 ≤ d j (t) ≤ d k (t) for j < k. Then d 1 (t) is the Hausdorff distance from t to σ(A), (3)
d 1 (t) = min{λ ∈ σ(A) : |λ -t|} = inf u∈D(A) |u| t u .
Similarly d j (t) are the distances from t to the jth nearest point in σ(A) counting multiplicity in a generalized sense. That is, stopping when the essential spectrum is reached. Moreover

d j (t) = d j-1 (t) ⇐⇒      either dim E [t-dj-1(t),t+dj-1(t)] (A) > j -1 or t + d j-1 (t) ∈ σ ess (A) or t -d j-1 (t) ∈ σ ess (A).
Without further mention, below we will always count spectral points of A relative to t, regarding multiplicities in this generalized sense. We now show how to extract certified information about σ(A) in the vicinity of t from the action of A onto finite-dimensional trial subspaces L ⊂ D(A), see [START_REF] Davies | Spectral enclosures and complex resonances for general self-adjoint operators[END_REF]Section 3]. For j ≤ n = dim L, let (4)

F j L (t) = min dim V =j V ⊂L max u∈V |u| t u . Then 0 ≤ F 1 L (t) ≤ . . . ≤ F n L (t) and F j L (t) ≥ d j (t) for all j = 1, 2, . . . , n. Since [t -d j (t), t + d j (t)] ⊆ [t -F j L (t), t + F j L (t)]
, there are at least j spectral points of A in the segment t -F j L (t), t + F j L (t) including, possibly, the essential spectrum. That is [START_REF] Bernhard | On a theorem of Danskin with an application to a theorem of von Neumann-Sion[END_REF] Tr

1 [t-F j L (t),t+F j L (t)] (A) ≥ j ∀j = 1, . . . , n.
Hence F j L (t) is an approximated local counting function for σ(A). As a consequence of the triangle inequality,

F j L is a Lipschitz continuous function such that (6) |F j L (t) -F j L (s)| ≤ |t -s| ∀s, t ∈ R and j = 1, . . . , n.
Moreover, F j L (t) is the jth smallest eigenvalue µ of the non-negative weak problem:

(7) find (µ, u) ∈ [0, ∞) × L\{0} such that q t (u, v) = µ 2 u, v ∀v ∈ L. Hence (8) F j L (t) = max dim V =j-1 V ⊂L min u∈L V |u| t u = max dim V =j-1 V ⊂H min u∈L V |u| t u .
2.1. Optimal setting for detection of the spectrum. As we show next, it is possible to detect the spectrum of A to the left/right of t by means of F j L in an optimal setting. This turns out to be a crucial ingredient in the formulation of the strategy proposed in [START_REF] Davies | Spectral enclosures and complex resonances for general self-adjoint operators[END_REF][START_REF]A hierarchical method for obtaining eigenvalue enclosures[END_REF][START_REF] Davies | Spectral pollution[END_REF].

The following notation simplifies various statements below. Let

n - j (t) = sup{s < t : Tr 1 (s,t] (A) ≥ j} and n + j (t) = inf{s > t : Tr 1 [t,s) (A) ≥ j}. Then n ∓ j (t)
is the jth point in σ(A) to the left(-)/right(+) of t counting multiplicities. Here t ∈ σ(A) is allowed and neither t nor n ∓ 1 (t) have to be isolated from the rest of σ(A). Note that n - j (t) = -∞ for Tr 1 (-∞,t] (A) < j and n + j (t) = +∞ for Tr 1 [t,+∞) (A) < j. Without further mention, all statements below regarding bounds on n ∓ j (t) will be void (hence redundant) in either of these two cases. Proposition 1. Let t -< t < t + . Then

(9) F j L (t -) ≤ t -t - ⇒ t --F j L (t -) ≤ n - j (t) F j L (t + ) ≤ t + -t ⇒ t + + F j L (t + ) ≥ n + j (t). Moreover, let t - 1 < t - 2 < t < t + 2 < t + 1 . Then (10) F j L (t - i ) ≤ t -t - i for i = 1, 2 ⇒ t - 1 -F j L (t - 1 ) ≤ t - 2 -F j L (t - 2 ) ≤ n - j (t) F j L (t + i ) ≤ t + i -t for i = 1, 2 ⇒ t + 1 + F j L (t + 1 ) ≥ t + 2 + F j L (t + 2 ) ≥ n + j (t).
Proof. We firstly show [START_REF] Bonito | Approximation of the eigenvalue problem for the time harmonic Maxwell system by continuous Lagrange finite elements[END_REF]. Suppose that t ≥ F j L (t -) + t -. Then Tr

1 [t --F j L (t -),t] (A) ≥ j. Since n - j (t) ≤ . . . ≤ n - 1 (t) are the only spectral points in the segment [n - j (t), t], then necessarily n - j (t) ∈ [t --F j L (t -), t].
The bottom of ( 9) is shown in a similar fashion.

The second statement follows by observing that the maps t → t ± F j L (t) are monotonically increasing as a consequence of [START_REF] Birman | The self-adjoint Maxwell operator in arbitrary domains[END_REF].

The structure of the trial subspace L determines the existence of t ± satisfying the hypothesis in [START_REF] Bonito | Approximation of the eigenvalue problem for the time harmonic Maxwell system by continuous Lagrange finite elements[END_REF]. If we expect to detect σ(A) at both sides of t, a necessary requirement on L should certainly be the condition [START_REF] Bramble | The approximation of the Maxwell eigenvalue problem using a least-squares method[END_REF] min

u∈L Au, u u, u < t < max u∈L Au, u u, u .
By virtue of lemmas 7 and 8 below, for j = 1, the left hand side inequality of [START_REF] Bramble | The approximation of the Maxwell eigenvalue problem using a least-squares method[END_REF] implies the existence of t -and the right hand side inequality implies the existence of t + , respectively.

Remark 1. From Proposition 1 it follows that optimal lower bounds for n - j (t) are achieved by finding tj ≤ t, the closer point to t, such that F j L ( tj ) = t -tj . Indeed, by virtue of [START_REF] Boulton | Eigenvalue enclosures for the MHD operator[END_REF]

, t --F j L (t -) ≤ t- j -F j L ( t- j ) ≤ n - j (t)
for any other t -as in [START_REF] Bonito | Approximation of the eigenvalue problem for the time harmonic Maxwell system by continuous Lagrange finite elements[END_REF]. Similarly, optimal upper bounds for n + j (t) are found by analogous means. This observation will play a crucial role in Section 3.

The main result of this section is Proposition 1, which is central to the hierarchical method for finding eigenvalue inclusions examined a few years ago in [START_REF] Davies | Spectral enclosures and complex resonances for general self-adjoint operators[END_REF][START_REF]A hierarchical method for obtaining eigenvalue enclosures[END_REF]. For fixed L this method leads to bounds for eigenvalues which are far sharper than those obtained from the obvious idea of estimating local minima of F 1 L (t). From an abstract perspective, Proposition 1 provides an intuitive insight on the mechanism for determining complementary bounds for eigenvalues (in the left definite case, for example). The method proposed in [START_REF] Davies | Spectral enclosures and complex resonances for general self-adjoint operators[END_REF][START_REF]A hierarchical method for obtaining eigenvalue enclosures[END_REF][START_REF] Davies | Spectral pollution[END_REF]] is yet to be explored more systematically in the practical setting, however in most circumstances the technique described in [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF] is easier to implement.

Geometrical properties of the first approximated counting function.

We now determine further geometrical properties of F 1 L and its connection to the spectral distance. Let the Hausdorff distances from t ∈ R to σ(A) \ (-∞, t] and σ(A) \ [t, ∞), respectively, be given by [START_REF] Buffa | Solving electromagnetic eigenvalue problems in polyhedral domains[END_REF] δ

+ (t) = inf{µ -t : µ ∈ σ(A), µ > t} and δ -(t) = inf{t -µ : µ ∈ σ(A), µ < t}.
In general, t -n - 1 (t) ≤ δ -(t) and n + 1 (t) -t ≤ δ + (t). In fact, |n ± 1 (t) -t| = δ ± (t) for t ∈ σ(A). However, these relations can be strict whenever t ∈ σ(A). Indeed, n + 1 (t) -t = δ + (t) iff there exists a decreasing sequence t + n ∈ σ(A) such that t + n ↓ t, whereas t -n - 1 (t) = δ -(t) iff there exists an increasing sequence t - n ∈ σ(A) such that t - n ↑ t. An emphasis in distinguishing |n ± 1 (t) -t| from δ ± (t) seems unnecessary at this stage. However, this distinction in the notation will be justified later on. Without further mention below we write δ ± (t) = ±∞ to indicate that either of the sets on the right side of (12) is empty.

Let λ ∈ σ(A) be an isolated point. If there exists a non-vanishing u ∈ L ∩ E λ (A), then

|u| s u = |λ -s| = d 1 (s) ∀s ∈ λ - δ -(λ) 2 , λ + δ + (λ) 2 .
According to the convergence analysis carried out in Section 4, the smaller the angle between L and the spectral subspace E λ (A), the closer the

F 1 L (t) is to d 1 (t) for t ∈ λ -δ -(λ) 2 , λ + δ + (λ) 2
. The special case of this angle being zero is described by the following lemma.

Lemma 2. For λ ∈ σ(A) isolated from the rest of the spectrum, the following statements are equivalent.

a) There exists a minimizer u ∈ L of the right side of (4) for j = 1, such that

|u| t = d 1 (t) for a single t ∈ λ -δ -(λ) 2 , λ + δ + (λ) 2 , b) F 1 L (t) = d 1 (t) for a single t ∈ λ -δ -(λ) 2 , λ + δ + (λ) 2 , c) F 1 L (s) = d 1 (s) for all s ∈ [λ -δ -(λ) 2 , λ + δ + (λ) 2 ], d) L ∩ E λ (A) = {0}.
Proof. Since L is finite-dimensional, a) and b) are equivalent by the definitions of d 1 (t), F 1 L (t) and q t . From the paragraph above the statement of the lemma it is clear that d) ⇒ c) ⇒ b). Since |u| t / u is the square root of the Rayleigh quotient associated to the operator (A -t) 2 , the fact that λ is isolated combined with the Rayleigh-Ritz principle, gives the implication a)⇒d).

As there can be a mixing of eigenspaces, it is not possible to replace b) in this lemma by an analogous statement including t = λ ± δ ± (λ)

2 . If λ = λ + δ + (λ) is an eigenvalue, for example, then

F 1 L λ+λ 2 = d 1 λ+λ 2
ensures that L contains elements of E λ (A) ⊕ E λ (A). However it is not guaranteed to be orthogonal to either of these two subspaces.

∂ ± t f (t) = lim τ →0 + ± f (t ± τ ) -f (t) τ ,
the one-side derivatives of a function f : J -→ R. Let V be a compact topological space. For given J : J × V -→ R we write

J (t) = max v∈V J (t, v) and Ṽ(t) = ṽ ∈ V : J (t) = J (t, ṽ) . Lemma 3. If the map J is upper semi-continuous and ∂ ± t J (t, v) exist for all (t, v) ∈ J × V, then also ∂ ± t J (t) exist for all t ∈ J and (13) ∂ ± t J (t) = max ṽ∈ Ṽ(t) ∂ ± t J (t, ṽ).
In the statement of this lemma, note that the left and right derivatives of both J and J might possibly be different. Proposition 4. Let j = 1, . . . , n and t ∈ R be fixed. The following assertions are equivalent.

a)

|F j L (t) -F j L (s)| = |t -s| for some s = t. b)
There exists an open segment J ⊂ R containing t in its closure, such that

|F j L (t) -F j L (s)| = |t -s| ∀s ∈ J. c) There exists an open segment J ⊂ R containing t in its closure, such that ∀s ∈ J, either L ∩ E s+F j L (s) = {0} or L ∩ E s-F j L (s) (A) = {0}. Proof. a) ⇒ b). Assume a). Since r → r ± F j L (r)
are continuous and monotonically increasing, then they have to be constant in the closure of

J = {τ t + (1 -τ )s : 0 < τ < 1}. This is precisely b). b) ⇒ c). Assume b). Then s → F j L (s)
is differentiable in J and its one-side derivatives are equal to 1 or -1 in the whole of this interval. For this part of the proof, we aim at applying [START_REF] Chatelin | Spectral Approximation of Linear Operators[END_REF], in order to get another expression for these derivatives.

Let F j be the family of (j -1)-dimensional linear subspaces of L. Identify an orthonormal basis of L with the canonical basis of C n . Then any other orthonormal basis of L is represented by a matrix in O(n), the orthonormal group. By picking the first (j-1) columns of these matrices, we cover all possible subspaces V ∈ F j . Indeed we just have to identify

(v 1 | . . . |v j-1 ) for [v kl ] n kl=1 ∈ O(n) with V = Span{v k } j-1 k=1 . Let K j = (v 1 , . . . , v j-1 ) : [v kl ] n kl=1 ∈ O(n) ⊂ C n × . . . × C n j-1
.

Then K j is a compact subset in the product topology of the right hand side. According to (8),

F j L (s) = max (v 1 ,...,v j-1 )∈Kj g(s; v 1 , . . . , v j-1 ) where g(s; v 1 , . . . , v j-1 ) = min (a1,...,aj-1)∈C j-1 |a k | 2 =1 a k ṽk s .
Here we have used the correspondence between v k ∈ C n and ṽk ∈ L in the orthonormal basis set above. We write

g(r, V ) = g(r; v 1 , . . . , v j-1 ) for V = Span{ṽ k } j-1 k=1 ∈ F j . The map g : J × K j -→ R + is
the minimum of a differentiable function, so the hypotheses of Lemma 3 are satisfied by J = -g. Hence, by virtue of ( 13),

∂ ± s g(s, V ) = min u∈L V, u =1 |u|s=g(s,V ) Re l s (u, u) |u| s .
As minima of continuous functions, g(s, V ) and ∂ ± s g(s, V ) are upper semi-continuous. Therefore, a further application of Lemma 3 yields

∂ ± s F j L (s) = max (v 1 ,...,v j-1 )∈Kj g(s;v 1 ,...,v j-1 )=F j L (s) ∂ ± s g(s, v 1 , . . . , v j-1 ) = max V ∈Fj g(s,V )=F j L (s) min u∈L V, u =1 |u|s=g(s,V ) Re l s (u, u) |u| s . Now, this shows that max V ∈Fj g(s,V )=F j L (s) min u∈L V, u =1 |u|s=g(s,V ) Re l s (u, u) |u| s = 1.
As L is finite dimensional, there exists a vector u ∈ L satisfying

|u| s = F j L (s) such that | Re l s (u, u)| |u| s = 1. Thus | Re (A -s)u, u | = (A -s)u, (A -s)u = F j L (s)
. Hence, according to the "equality" case in the Cauchy-Schwarz inequality, u must be an eigenvector of A associated with either s + F j L (s) or s -F j L (s). This is precisely c). c) ⇒ a). Under the condition c), there exists an open segment J ⊆ J, possibly smaller, such that t ∈ J and F j L (s) = d j (s) for all s ∈ J. As |d j (s) -d j (r)| = |s -r|, then either a) is immediate, or it follows by taking r → t.

As a consequence of this statement, we find the following extension of Proposition 1 for t an eigenvalue.

Corollary 5. Let t ∈ σ(A) be an eigenvalue of multiplicity m. Let t -< t < t + . If E t (A) ∩ L = {0}, then (14) F j L (t -) ≤ t -t - ⇒ t --F j L (t -) ≤ n - j+m (t) F j L (t + ) ≤ t + -t ⇒ t + + F j L (t + ) ≥ n + j+m (t). Proof. According to (5), Tr 1 [t --F j L (t -),t -+F j L (t -)] (A) ≥ j. Thus, if t > F j L (t -) + t -, there is nothing to prove. Consider now the case t = F j L (t -) + t -. If there exists τ < t -such that t = F j L (τ ) + τ , then from Proposition 4 there exists an open segment J ⊂ R containing (τ, t -) such that ∀s ∈ J, either L ∩ E s+F j L (s) = {0} or L ∩ E s-F j L (s) (A) = {0}.
From the assumption, only the second alternative takes place, and necessarily

∀s ∈ (τ, t -), s -F j L (s) ∈ σ p (A). Hence, as s -F j L (s)
is continuous and H is separable, this function should be constant in the segment (τ, t -). We also notice that due to monotonicity for any

s ∈ (τ, t -), s + F j L (s) = t -. Hence if s ∈ (τ, t -) → s -F j L (s)
is constant, and equal to some value (say v), then s is the midpoint between t and v for any s ∈ (τ, t -), which is a contradiction with the fact that τ = t -. Hence t > F j L (τ ) + τ, ∀τ < t - and so τ -F j L (τ ) ≤ n - j+m (t), for all τ < t -. Thus, by continuity, also

t --F j L (t -) ≤ n - j+m (t)
. The bottom of ( 14) is shown in a similar fashion.

2.4. Approximated eigenspaces. We conclude this section by examining extensions of the implications b) ⇒ d) of Lemma 2 into a more general context. In combination with the results of Section 3, the next proposition shows how to obtain certified information about spectral subspaces.

Here and below {u t j } n j=1 ⊂ L will denote an orthonormal family of eigenfunctions associated to the eigenvalues µ = F j L (t) of the weak problem [START_REF] Boffi | Finite element approximation of eigenvalue problems[END_REF]. In a suitable asymptotic regime for L, the angle between these eigenfunctions and the spectral subspaces of |A -t| in the vicinity of the origin is controlled by a residual which is

as small as O F j L (t) -d j (t) for F j L (t) -d j (t) → 0.
Assumption 1. Unless otherwise specified, from now on we will always fix the parameter m ≤ n = dim L and suppose that

(15) [t -d m (t), t + d m (t)] ∩ σ(A) ⊆ σ disc (A). Set δ j (t) = dist t, σ(A) \ {t ± d k (t)} j k=1
. By virtue of (15), δ j (t) > d j (t) for all j ≤ m.

Remark 2. If t = n - j (t)+n + j (t) 2
for a given j, the vectors φ t j introduced in Proposition 6 and invoked subsequently, might not be eigenvectors of A despite of the fact that |A -t|φ t j = d j (t)φ t j . However, in any other circumstance φ t j are eigenvectors of A.

Proposition 6. Let t ∈ R and j ∈ {1, . . . , m}. Assume that F j L (t) -d j (t) is small enough so that 0 < ε j < 1 holds true for the residuals constructed inductively as follows,

ε 1 = F 1 L (t) 2 -d 1 (t) 2 δ 1 (t) 2 -d 1 (t) 2 ε j = F j L (t) 2 -d j (t) 2 δ j (t) 2 -d j (t) 2 + j-1 k=1 ε 2 k 1 -ε 2 k 1 + d j (t) 2 -d k (t) 2 δ j (t) 2 -d j (t) 2 .
Then, there exists an orthonormal basis

{φ t j } m j=1 of E [t-dm(t),t+dm(t)] (A) such that φ t j ∈ E {t-dj (t),t+dj (t)} (A), u t j -u t j , φ t j φ t j ≤ ε j and (16)
|u t j -u t j , φ t j φ t j | t ≤ F j L (t) 2 -d j (t) 2 + d j (t) 2 ε 2 j . ( 17
)
Proof. As it is clear from the context, in this proof we suppress the index t on top of any vector. We write Π S to denote the orthogonal projection onto the subspace S with respect to the inner product •, • .

Let us first consider the case j = 1. Let S 1 = E {t-d1(t),t+d1(t)} (A), and decompose 16) holds immediately and ( 17) is achieved by clearing (A -t)u ⊥ 1 2 from (18). We define the needed basis, and show ( 16) and ( 17), for j up to m inductively as follows. Set

u 1 = Π S1 u 1 + u ⊥ 1 where u ⊥ 1 ⊥ S 1 . Since A is self-adjoint, (18) F 1 L (t) 2 = (A -t)u 1 2 = d 1 (t) 2 Π S1 u 1 2 + (A -t)u ⊥ 1 2 . Hence F 1 L (t) 2 ≥ d 1 (t) 2 (1 -u ⊥ 1 2 ) + δ 1 (t) 2 u ⊥ 1 2 . Since δ 1 (t) > d 1 (t), clearing from this identity u ⊥ 1 2 yields u ⊥ 1 ≤ε 1 . Hence Π S1 u 1 2 ≥ 1 -ε 2 1 > 0. Let φ 1 = 1 Π S1 u 1 Π S1 u 1 so that Π S1 u 1 = | u 1 , φ 1 |. Then (
φ j = 1 Π Sj u j Π Sj u j where S j = E {t-dj (t),t+dj (t)} (A) Span{φ l } j-1 1
and Π Sj u j = 0, all this for 1 ≤ j ≤ k -1. Assume that ( 16) and ( 17) hold true for j up to k -

1. Define S k = E {t-d k (t),t+d k (t)} (A) Span{φ l } k-1 1
. We first show that Π S k u k = 0, and so we can define [START_REF] Girault | Finite element methods for Navier-Stokes equations[END_REF]

φ k = 1 Π S k u k Π S k u k ensuring φ k ⊥ Span{φ l } k-1 l=1
. After that we verify the validity of ( 16) and ( 17) for j = k.

Decompose

u k = Π S k u k + 1 l=k-1 u k , φ l φ l + u ⊥ k where u ⊥ k ⊥ Span{φ l } k-1 l=1 ⊕ S k . Then F k L (t) 2 = d k (t) 2 Π S k u k 2 + 1 l=k-1 d l (t) 2 | u k , φ l | 2 + (A -t)u ⊥ k 2 ≥ d k (t) 2 Π S k u k 2 + 1 l=k-1 d l (t) 2 | u k , φ l | 2 + δ k (t) 2 u ⊥ k 2 = d k (t) 2 (1 -u ⊥ k 2 ) + 1 l=k-1 (d l (t) 2 -d k (t) 2 )| u k , φ l | 2 + δ k (t) 2 u ⊥ k 2 .
The conclusion [START_REF]A hierarchical method for obtaining eigenvalue enclosures[END_REF]

up to k -1, implies | u l , φ l | 2 ≥ 1 -ε 2 l for l = 1, . . . , k -1. Since u k , u l = 0 for l = k, | u l , φ l || u k , φ l | = | u k , u l -u l , φ l φ l |.
Then, the Cauchy-Schwarz inequality alongside with ( 16) yield

(20) | u k , φ l | 2 ≤ ε 2 l 1 -ε 2 l . Hence, since d l (t) ≤ d k (t), F k L (t) 2 ≥ d k (t) 2 + 1 l=k-1 (d l (t) 2 -d k (t) 2 ) ε 2 l 1 -ε 2 l + (δ k (t) 2 -d k (t) 2 ) u ⊥ k 2 .
Clearing u ⊥ k 2 from this inequality and combining with the validity of ( 20) and ( 16) up to k -1, yields Π S k u k = 0.

Let φ k be as in [START_REF] Girault | Finite element methods for Navier-Stokes equations[END_REF]. Then ( 16) is guaranteed for j = k. On the other hand, (16) up to j = k, (20) and the identity

F k L (t) 2 = d k (t) 2 | u k , φ k | 2 + (A -t)(u k -u k , φ k φ k ) 2 , yield (17) up to j = k.

Local bounds for eigenvalues

Let t ∈ R and L ⊂ D(A) be a specified trial subspace as above. Recall that q t is given by [START_REF] Amrouche | Vector potentials in three-dimensional non-smooth domains[END_REF]. Let l t : D(A) × D(A) -→ C be the (generally not closed) bi-linear form associated to (A -t), l t (u, w) = (A -t)u, w ∀u, w ∈ D(A).

Our next purpose is to characterize the optimal parameters t ± in Proposition 1 as described in Remark 1 by means of the following weak eigenvalue problem,

(Z L t ) find u ∈ L \ {0} and τ ∈ R such that τ q t (u, v) = l t (u, v) ∀v ∈ L.
This problem is central to the method of eigenvalue bounds calculation examined in [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF]. Let

τ - 1 (t) ≤ . . . ≤ τ - n -(t) < 0 and 0 < τ + n + (t) ≤ . . . ≤ τ + 1 (t),
be the negative and positive eigenvalues of (Z L t ) respectively. Here and below n ∓ (t) is the number of these negative and positive eigenvalues, which are both locally constant in t. Below we will denote eigenfunctions associated with τ ∓ j (t) by u ∓ j (t). Assumption 2. For the purpose of clarity of exposition and without further mention, below we write most statements only for the case of "lower bounds for the eigenvalues of A which are to the left of t". As the position of t relative to the essential spectrum is irrelevant here, evidently this assumption does not restrict generality. The corresponding results regarding "upper bounds for the eigenvalues of A which are to the right of t" can be recovered by replacing A by -A.

The left side of the hypotheses [START_REF] Bramble | The approximation of the Maxwell eigenvalue problem using a least-squares method[END_REF] ensures the existence of τ - 1 (t). A more concrete connection with the framework of Section 2 is made precise in the following lemma. Its proof is straightforward, hence omitted. Lemma 7. The following conditions are equivalent, a -) F 1 L (s) > t -s for all s < t b -) Au,u u,u > t for all u ∈ L c -) all the eigenvalues of (Z L t ) are positive.

Remark 3. Let L = Span{b j } n j=1 . The matrix [q t (b j , b k )] n jk=1 is singular if and only if E t (A) ∩ L = {0}. On the other hand, the kernel of (Z L t ) might be nonempty. If n 0 (t) is the dimension of this kernel and n

∞ (t) = dim(E t (A) ∩ L), then n = n ∞ (t) + n 0 (t) + n -(t) + n + (t). Assumption 3. Note that n ∞ (t) ≥ 1 if and only if F j L (t) = 0 for j = 1, . . . , n ∞ (t).
In this case the conclusions of Lemma 8 and Theorem 9 below become void. In order to write our statements in a more transparent fashion, without further mention from now on we will suppose that

(21) L ∩ E t (A) = {0}.
By virtue of the next three results, finding the negative eigenvalues of (Z L t ) is equivalent to finding s = tj ∈ R such that [START_REF] Weinberger | Variational Methos for Eigenvalue Approximation[END_REF] t -s = F j L (s), and in this case t-

j = t + 1 2τ - j (t)
. It then follows from Remark 1 that (Z L t ) encodes information about the optimal bounds for the spectrum around t, achievable by [START_REF] Boulton | Eigenvalue enclosures for the MHD operator[END_REF] in Proposition 1.

3.1.

The eigenvalue immediately to the left. We begin with the case j = 1, see [START_REF] Davies | Spectral pollution[END_REF]Theorem 11].

Lemma 8. Let t ∈ R. The smallest eigenvalue τ = τ - 1 (t) of (Z L t )
is negative if and only if there exists s < t such that (22) holds true. In this case

s = t + 1 2τ - 1 (t)
and

F 1 L (s) = - 1 2τ - 1 (t) = |u - 1 (t)| s u - 1 (t)
for u = u - 1 (t) ∈ L the corresponding eigenvector.

Proof. For all u ∈ L and s ∈ R,

q s (u, u) -F 1 L (s) 2 u, u = q t (u, u) + 2(t -s)l t (u, u) + (t -s) 2 -F 1 L (s) 2 u, u . Suppose that F 1 L (s) = t -s. Then q s (u, u) -F 1 L (s) 2 u, u = q t (u, u) + 2F 1 L (s)l t (u, u).
As the left side of this expression is non-negative,

l t (u, u) q t (u, u) ≥ - 1 2F 1 L (s)
for all u ∈ L \ {0} and the equality holds for some u ∈ L. Hence -1 2F 1 L (s) is the smallest eigenvalue of (Z L t ), and thus necessarily equal to τ - 1 (t). In this case s -

F 1 L (s) = t -2F 1 L (s) = t + 1 τ - 1 (t)
. Here the vector u for which equality is achieved is exactly u = u - 1 (t). Conversely, let τ - 1 (t) and u - 1 (t) be as stated. Then

τ - 1 (t) ≤ l t (u, u) q t (u, u)
for all u ∈ L with equality for u = u - 1 (t). Re-arranging this expression yields

q t (u, u) - 1 τ - 1 (t) l t (u, u) ≥ 0 for all u ∈ L with equality for u = u - 1 (t). The substitution t = s -1 2τ - 1 (t) then yields q t (u, u) - 1 (2τ - 1 (t)) 2 u, u ≥ 0
for all u ∈ L. The equality holds for u = u - 1 (t). This expression further re-arranges as

|u| 2 s u 2 ≥ 1 (2τ - 1 (t)) 2 . Hence F 1 L (s) 2 = 1 (2τ - 1 (t)) 2 , as needed.
3.2. Further eigenvalues. An extension to j = 1 is now found by induction. Theorem 9. Let t ∈ R and 1 ≤ j ≤ n be fixed. The number of negative eigenvalues n -(t) in (Z L t ) is greater than or equal to j if and only if Au, u u, u < t for some u ∈ L Span{u - 1 (t), . . . , u - j-1 (t)}.

Assuming this holds true, then τ = τ - j (t) and u = u - j (t) are solutions of (Z L t ) if and only if

F j L t + 1 2τ - j (t) = - 1 2τ - j (t) = u - j (t) t+ 1 2τ - j (t) u - j (t)
.

Proof. For j = 1 the statements are Lemma 8 taking into consideration [START_REF] Bramble | The approximation of the Maxwell eigenvalue problem using a least-squares method[END_REF]. For j > 1, due to the self-adjointness of the eigenproblem (Z L t ), it is enough to apply again Lemma 8 by fixing L = L Span{u - 1 (t), . . . , u - j-1 (t)} as trial spaces. Note that the negative eigenvalues of (Z L t ) are those of (Z L t ) except for τ - 1 (t), . . . , τ - j-1 (t).

A neat procedure for finding certified spectral bounds for A, as described in [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF], can now be deduced from Theorem 9. By virtue of Proposition 1 and Remark 1, this procedure turns out to be optimal in the context of the approximated counting functions discussed in Section 2, see [START_REF] Davies | Spectral pollution[END_REF]Section 6]. We summarize the core statement as follows.

Corollary 10. For all t ∈ R and j ∈ {1, . . . , n ± (t)}, [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF] t

+ 1 τ - j (t) ≤ n - j (t) and n + j (t) ≤ t + 1 τ + j (t)
.

In recent years, numerical techniques based on this statement have been designed to successfully compute eigenvalues for the radially reduced magnetohydrodynamics operator [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF][START_REF] Boulton | Eigenvalue enclosures for the MHD operator[END_REF], the Helmholtz equation [START_REF] Behnke | Bounds for eigenvalues with the use of finite elements[END_REF] and the calculation of sloshing frequencies in the left definite case [START_REF] Behnke | Lower and upper bounds for sloshing frequencies[END_REF]. We will explore the case of the Maxwell operator in sections 5.

Convergence and error estimates

Our first goal in this section will be to show that, if L captures an eigenspace of A within a certain order of precision O(ε) as specified below, then the bounds which follow from Proposition 1 are a) at least within O(ε) from the true spectral data for any t ∈ R, b) within O(ε 2 ) for t ∈ σ(A). This will be the content of theorems 12 and 13, and Corollary 14. We will then show that, in turns, the estimates [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF] have always residual of size O(ε 2 ) for any t ∈ R. See Theorem 15. In the spectral approximation literature this property is known as optimal order of convergence/exactness, see [START_REF] Chatelin | Spectral Approximation of Linear Operators[END_REF]Chapter 6] or [START_REF] Weinberger | Variational Methos for Eigenvalue Approximation[END_REF].

Recall Remark 2, and the assumptions 1 and 3. Below {φ t j } m j=1 denotes an orthonormal set of eigenvectors of E [t-dm(t),t+dm(t)] (A) which is ordered so that |A -t|φ t j = d j (t)φ t j for j = 1, . . . , m.

Whenever 0 < ε j < 1 is small, as specified below, the trial subspace L ⊂ D(A) will be assumed to be close to Span{φ t j } m j=1 in the sense that there exist w t j ∈ L such that w t j -φ t j ≤ ε j and (A 0 )

|w t j -φ t j | t ≤ ε j . (A 1 )
We have split this condition into two, in order to highlight the fact that some times only (A 1 ) is required. Unless otherwise specified, the index j runs from 1 to m.

From [START_REF] Davies | Spectral enclosures and complex resonances for general self-adjoint operators[END_REF] it follows that the family {φ s j } m j=1 ⊂ E [t-dm(t),t+dm(t)] (A) and the family {w s j } m j=1 ⊂ L above can always be chosen piecewise constant for s in a neighbourhood of t. Moreover, they can be chosen so that jumps only occur at s ∈ σ(A). Assumption 4. Without further mention all t-dependant vectors below will be assumed to be locally constant in t with jumps only at the spectrum of A.

A set {w t j } m j=1 subject to (A 0 )-(A 1 ) is not generally orthonormal. However, according to the next lemma, it can always be substituted by an orthonormal set, provided ε j is small enough. Lemma 11. There exists a constant C > 0 independent of L ensuring the following.

If {w t j } m j=1 ⊂ L is such that (A 0 )-(A 1 ) hold for all ε j such that ε = m j=1 ε 2 j < 1 √ m ,
then there is a set {v t j } m j=1 ⊂ L orthonormal in the inner product •, • such that

|v t j -φ t j | t + v t j -φ t j < Cε.
Proof. As it is clear from the context, in this proof we suppress the index t on top of any vector. The desired conclusion is achieved by applying the Gram-Schmidt procedure. Let G = [ w k , w l ] m kl=1 ∈ C m×m be the Gram matrix associated to {w j }. Set

v j = m k=1 (G -1/2 ) kj w k . Then G -I ≤ m k,l=1 | w k , w l -φ k , φ l | 2 ≤ 2 m k,l=1 w k -φ k 2 ( w l + φ l ) 2 ≤ √ 2(2 + ε)ε. Since v j -w j 2 = m k=1 (G -1/2 -I) kj w k 2 = m k,l=1 (G -1/2 -I) kj (G -1/2 -I) lj w k , w l = m k=1 (G -1/2 -I) kj m l=1 G kl (G -1/2 -I) lj = m k=1 (G -1/2 -I) kj (G 1/2 -G) jk = (I -G 1/2 ) 2 jj then v j -w j ≤ I -G 1/2 .
As G 1/2 is a positive-definite matrix, for every v ∈ C m we have

(G 1/2 + I)v 2 = G 1/2 v 2 + 2 G 1/2 v, v + v 2 ≥ v 2 .
Then det(I + G 1/2 ) = 0 and

(I + G 1/2 ) -1 ≤ 1. Hence (24) v j -w j ≤ (I -G)(I + G 1/2 ) -1 ≤ I -G (I + G 1/2 ) -1 ≤ (2 + ε)ε . Now, identify v = (v 1 , . . . , v m ) ∈ C m with v = m k=1 v k φ k . As G 1/2 v = m j=1 v, φ j w j ≥ v - m j=1 v, φ j (w j -φ j ) ≥ (1 -ε) v then G -1/2 ≤ 1 1 -ε .
Hence

|v j -w j | t ≤ m k=1 |(G -1/2 -I) jk ||w k | t ≤ m k=1 |(G -1/2 -I) jk |(ε k + d k (t)) ≤ m k,l=1 |(G -1/2 ) kl ||(G 1/2 -I) lj |(ε k + d k (t)) ≤ √ m(ε + d m (t))(2 + ε) 1 -ε ε. (25)
The desired conclusion follows from ( 24) and (25).

Convergence of the approximated local counting function.

The next theorem addresses the claim a) made at the beginning of this section. According to Lemma 11, in order to examine the asymptotic behaviour of F j L (t) as ε j → 0 under the constraints (A 0 )-(A 1 ), we can assume without loss of generality that the trial vectors w t j form an orthonormal set in the inner product •, • .

Theorem 12. Let {w t j } m j=1 ⊂ L be a family of vectors which is orthonormal in the inner product •, • and satisfies (A 1 ). Then

F j L (t) -d j (t) ≤ j k=1 ε 2 k 1/2 ∀j = 1, . . . , m.
Proof. From the min-max principle we obtain

F j L (t) ≤ max |c k | 2 =1 j k=1 c k w k t ≤ max |c k | 2 =1 j k=1 c k (w k -φ k ) t + max |c k | 2 =1 j k=1 c k φ k t = max |c k | 2 =1 j k=1 c k (w k -φ k ) t + d j (t).
This gives

F j L (t) -d j (t) ≤ max |c k | 2 =1 j k=1 |c k ||w k -φ k | t ≤ max |c k | 2 =1 j k=1 |c k | 2 1/2 j k=1 |w k -φ k | 2 t 1/2 ≤ j k=1 ε 2 k 1/2
as needed.

In terms of order of approximation, Theorem 12 will be superseded by Theorem 13 for t ∈ σ(A). However, if t ∈ σ(A), the trial space L can be chosen so that

F 1 L (t)-d 1 (t) is linear in ε 1 . Indeed, fixing any non-zero u ∈ D(A) and L = Span{u}, yields F 1 L (t) -d 1 (t) = F 1 L (t) = ε 1 .
This shows that Theorem 12 is optimal, upon the presumption that t is arbitrary.

The next theorem addresses the claim b) made at the beginning of this section. Its proof is reminiscent of that of [START_REF] Strang | An Analysis of the Finite Element Method[END_REF]Theorem 6.1].

Theorem 13. Let t ∈ σ(A). Suppose that the ε j in (A 1 ) are such that

(26) m j=1 ε 2 j < d 1 (t) 2 6 .
Then,

(27) F j L (t) -d j (t) ≤ 3 d j (t) d 1 (t) 2 j k=1 ε 2 k ∀j = 1, . . . , m.
Proof. Since t ∈ σ(A), then (D(A), q t (•, •)) is a Hilbert space. Let P L : D(A) -→ L be the orthogonal projection onto L with respect to the inner product q t (•, •), so that q t (u

-P L u, v) = 0 ∀v ∈ L. Then |u| 2 t = |P L u| 2 t + |u -P L u| 2 t for all u ∈ D(A) and |u -P L u| t ≤ |u -v| t for all v ∈ L. Hence (28) |φ k -P L φ k | t ≤ ε k ∀k = 1, . . . , m. Let E j = Span φ k } j k=1 . Define F j = {φ ∈ E j : φ = 1 and µ j L (t) = max φ∈Fj 2 Re φ, φ -P L φ -φ -P L φ 2 .
Here µ j L depends on t, as P L does. We first show that, under hypothesis (26),

µ j L (t) < 1 2 . Indeed, given φ ∈ F j we decompose it as φ = j k=1 c k φ k . Then | φ, φ -P L φ | = j k=1 c k φ k , φ -P L φ = j k=1 c k d k (t) 2 q t (φ k , φ -P L φ) = q t j k=1 c k d k (t) 2 φ k , φ -P L φ = q t j k=1 c k d k (t) 2 (φ k -P L φ k ), φ -P L φ ≤ j k=1 c k d k (t) 2 (φ k -P L φ k ) t j k=1 c k (φ k -P L φ k ) t . ( 29
)
For each multiplying term in the latter expression, the triangle and Cauchy-Schwarz's inequalities yield (take

α k = c k or α k = c k d k (t) 2 ) j k=1 α k (φ k -P L φ k ) t ≤ j k=1 |α k | |φ k -P L φ k | t ≤ j k=1 |α k | 2 1/2 j k=1 |φ k -P L φ k | 2 t 1/2 . (30) Then (31) |2 Re φ, φ -P L φ | ≤ 2 j k=1 |c k | 2 d k (t) 4 1/2 j k=1 |c k | 2 1/2 j k=1 ε 2 k ≤ 2 d 1 (t) 2 j k=1 ε 2 k for all φ ∈ F j .
The other term in the expression for µ j L (t) has an upper bound found as follows. According to the min-max principle

(32) φ -P L φ 2 ≤ 1 d 1 (t) 2 q t (φ -P L φ, φ -P L φ).
Therefore, by repeating analogous steps as in ( 29) and (30), we get

φ -P L φ 2 ≤ 1 d 1 (t) 2 j k=1 c k q t (φ k -P L φ k , φ -P L φ) = q t j k=1 c k d 1 (t) 2 (φ k -P L φ k ), φ -P L φ = q t j k=1 c k d 1 (t) 2 (φ k -P L φ k ), j l=1 c l (φ l -P L φ l ) ≤ 1 d 1 (t) 2 j k=1 ε 2 k . (33)
Hence, from (31) and ( 33), (34)

µ j L (t) ≤ 3 d 1 (t) 2 j k=1 ε 2 k < 1 2
as a consequence of (26). Next, observe that dim(P L E j ) = j. Indeed P L ψ = 0 for ψ = 1 would imply

µ j L (t) ≥ 2 Re ψ, ψ -P L ψ -ψ -P L ψ 2 = ψ 2 = 1,
which would contradict the fact that µ j L (t) < 1. Then,

F j L (t) 2 ≤ max u∈P L Ej |u| 2 t u 2 = max φ∈Ej |P L φ| 2 t P L φ 2 = max φ∈Fj |P L φ| 2 t P L φ 2 .
As

P L φ 2 = φ 2 -2 Re φ, φ -P L φ + φ -P L φ 2 ≥ 1 -µ j L (t), we get (35) F j L (t) 2 ≤ max φ∈Fj |φ| 2 t 1 -µ j L (t) = max |c k | 2 =1 j k=1 |c k | 2 d k (t) 2 1 -µ j L (t) = d j (t) 2 1 -µ j L (t)
.

Finally, (35) and (34) yield

F j L (t) 2 -d j (t) 2 ≤ µ j L (t) 1 -µ j L (t) d j (t) 2 ≤ 2µ j L (t)d j (t) 2 ≤ 2 3 d 1 (t) 2 d j (t) 2 j k=1 ε 2 k . ( 36
)
The proof is completed by observing that F j L (t) + d j (t) ≥ 2d j (t).

As the next corollary shows, a quadratic order of decrease for F j L (t) -d j (t) is prevented for t ∈ σ(A) in the context of theorems 12 and 13, only for j up to dim E t (A).

Corollary 14. Let t ∈ σ disc (A), = 1 + dim E t (A) and k ∈ { , . . . , m}. Let α k (t) = 1 4 min {|d l (t) -d l-1 (t)| : d l (t) = d l-1 (t), l = , ..., k} > 0.
There exists ε > 0 independent of k ensuring the following. If (A 1 ) holds true for m j=1 ε 2 j < ε, then

F k L (t) -d k (t) ≤ 3 d k (t) α k (t) 2 k j=1 ε 2 j .
Proof. Without loss of generality we assume that t

+ d k (t) ∈ σ(A). Otherwise t -d k (t) ∈ σ(A)
and the proof is analogous to the one presented below. Let t = t + α k (t). Then t ∈ σ(A) and t + d k (t) = t + d k ( t). Since the map s → s + F j L (s) is non-decreasing as a consequence of Proposition 1, Theorem 13 applied at t yields

F k L (t) -d k (t) = t + F k L (t) -(t + d k (t)) ≤ t + F k L ( t) -( t + d k ( t)) = F k L ( t) -d k ( t) ≤ 3 d k ( t) d 1 ( t) 2 k j=1 ε 2 k ≤ 3 d k (t) α k (t) 2 k j=1 ε 2 j as needed.
4.2. Convergence of local bounds for eigenvalues. For the final part of this section, we formulate precise statements on the convergence of the method described in Section 3. Theorem 15 below improves upon two crucial aspects of a similar result established in [10, Lemma 2]. It allows j > 1 and it allows t ∈ σ(A). These two improvements are essential in order to obtain sharp bounds for those eigenvalues which are either degenerate or form a tight cluster.

Remark 4. The constants εt and C ± t below do have a dependence on t that may be determined explicitly from Theorem 13, Corollary 14 and the proof of Theorem 15. Despite of the fact that they can deteriorate as t approaches the isolated eigenvalues of A and they can have jumps precisely at these points, they may be chosen locally independent of t in compacts outside the spectrum.

Set

ν - j (t) = sup{s < t : Tr 1 (s,t) (A) ≥ j} ν + j (t) = inf{s > t : Tr 1 (t,s) (A) ≥ j}.
Note that these are the spectral points of A which are strictly to the left and strictly to the right of t respectively. The inequality ν ± j (t) = n ± j (t) only occurs when t is an eigenvalue. In view of (12), δ ± (t) = |t -ν ± 1 (t)|. Remark 5. By virtue of Corollary 10 and Corollary 5, 1 τ - j (t) ≤ ν - j (t) -t and

1 τ + j (t) ≥ ν + j (t) -t. Then t- j = t + 1 2τ - j (t) ≤ t + ν - j (t) 2 ≤ ν + j (t) + ν - j (t) 2 ≤ ν + j (t) + t 2 ≤ t + 1 2τ + j (t) = t+ j .
The following is one of the main results of this paper.

Theorem 15. Let J ⊂ R be a bounded open segment such that J ∩σ(A) ⊆ σ disc (A).

Let {φ k } m k=1 be a family of eigenvectors of A such that Span{φ k } m k=1 = E J (A). For fixed t ∈ J, there exist constants εt > 0 and C - t > 0 independent of the trial space L, ensuring the following. If there are {w j } m j=1 ⊂ L such that Let j be such that ν - j (t) ∈ J. Since ν - j (t) -(α + t) ≤ (t + α) -ν - 1 (t) for all α such that

(37)   m j=1 w j -φ j 2 + |w j -φ j | 2 t   1/2 ≤ ε < εt , then 0 < ν - j (t) -t + 1 τ - j (t) ≤ C - t ε 2 for all j ≤ n -(t) such that ν - j (t) ∈ J. Proof.
ν - j (t)+ν - 1 (t) 2 -t ≤ α ≤ 0, then d j (s) = s -ν - j (t) ∀s ∈ ν - 1 (t) + ν - j (t) 2 , t + ν - j (t) 2 .
Let g(α) = F j L (t + α) + α. Then g is an increasing function of α and g(0) = F j L (t) > 0. For the strict inequality in the latter, recall Assumption 3. Moreover, according to (38)

g ν - j (t) + ν - 1 (t) 2 -t = F j L ν - j (t) + ν - 1 (t) 2 -t + ν - 1 (t) - ν - 1 (t) -ν - j (t) 2 = F j L ν - j (t) + ν - 1 (t) 2 -t + ν - 1 (t) -d j ν - j (t) + ν - 1 (t) 2 ≤ t -ν - 1 (t) 2 -(t -ν - 1 (t)) < 0 .
Hence, the Mean Value Theorem ensures the existence of α ∈ ν - 1 (t)+ν - j (t) 2 -t, 0 such that α = F j L (t + α). According to Theorem 9, α is unique and α = 1 2τ - j (t) .

The proof is now completed as follows. By virtue of Remark 5, t-

j (t) = t + 1 2τ - j (t) ∈ ν - 1 (t) + ν - j (t) 2 , t + ν - j (t) 2 and F j L ( t- j (t)) = 1 2τ - j (t)
.

Then, Theorem 13 or Corollary 14, as appropriate, ensure the existence of C - t > 0 yielding

ν - j (t) -t + 1 τ - j (t) = F j L ( t- j ) -d j ( t- j ) ≤ C - t j k=1 ε 2 k < C - t ε 2 ,
as needed.

4.3. Convergence for eigenfunctions. We conclude this section with a statement on convergence of eigenfunctions.

Corollary 16. Let J ⊂ R be a bounded open segment such that J ∩σ(A) ⊆ σ disc (A).

Let {φ k } m k=1 be a family of eigenvectors of A such that Span{φ k } m k=1 = E J (A). For fixed t ∈ J, there exist constants εt > 0 and C ± t > 0 independent of the trial space L, ensuring the following. If there are {w j } m j=1 ⊂ L guaranteeing the validity of (37), for all j ≤ n ± (t) such that ν ± j (t) ∈ J we can find

ψ ε± j ∈ E {ν - j (t),ν + j (t)} (A) satisfying |u ± j (t) -ψ ε± j | t + u ± j (t) -ψ ε± j ≤ C ± t ε.
Proof. Fix t ∈ J. By virtue of Theorem 9, u ± j (t) = u t± j j in the notation for eigenvectors employed in Proposition 6. The claimed conclusion is a consequence of the latter combined with Theorem 13 or Corollary 14, as appropriate.

The finite element method for the Maxwell eigenvalue problem

Let Ω ⊂ R 3 be a polyhedron which is open, bounded, simply connected and Lipschitz in the sense of [START_REF] Amrouche | Vector potentials in three-dimensional non-smooth domains[END_REF]Notation 2.1]. Let ∂Ω be the boundary of Ω and denote by n its outer normal vector. The physical phenomenon of electromagnetic oscillations in a resonator filled with a homogeneous medium is described by the isotropic Maxwell eigenvalue problem, (39)

     curl E = iωH in Ω curl H = -iωE in Ω E × n = 0 on ∂Ω.
Here the angular frequency ω ∈ R and the field phasor (E, H) = 0 is restricted to the solenoidal subspace, characterized by the Gauss law (40) div(E) = 0 = div(H).

The orthogonal complement of this subspace is the gradient space, which has infinite dimension and it lies in the kernel of the eigenvalue equation (39). In turns, this means that (39)-(40) and the unrestricted problem (39), have the same non-zero spectrum and the same corresponding eigenspace. Let

H(curl; Ω) = u ∈ [L 2 (Ω)] 3 : curl u ∈ [L 2 (Ω)] 3 equipped with the norm (41) u 2 curl,Ω = u 2 0,Ω + curl u 2 0,Ω .
Let R max denote the operator defined by the expression "curl" acting on the domain D(R max ) = H(curl; Ω), the maximal domain. Let

R min = R * max = R max [D(Ω)] 3 . The domain of R min is D(R min ) = H 0 (curl; Ω) = {u ∈ H(curl; Ω) : curl u, v Ω = u, curl v Ω ∀v ∈ H(curl; Ω)}.
By virtue of Green's identity for the rotational [START_REF] Girault | Finite element methods for Navier-Stokes equations[END_REF]Theorem I.2.11],

H 0 (curl; Ω) = {u ∈ H(curl; Ω) : u × n = 0 on ∂Ω}.
The linear operator associated to (39) is then,

M = 0 iR max -iR min 0 on the domain (42) D(M) = D(R min ) × D(R max ) ⊂ [L 2 (Ω)] 6 . Note that M : D(M 1 ) -→ [L 2 (Ω)] 6 is self-adjoint, as R max and R min are mutually adjoints [6, Lemma 1.2].
The numerical estimation of the eigenfrequencies of (39)-( 40) is known to be extremely challenging in general. The operator M does not have a compact resolvent and it is strongly indefinite. The self-adjoint operator associated to (39)-(40) has a compact resolvent but it is still strongly indefinite. By considering the square of M on the solenoidal subspace, one obtains a positive definite eigenvalue problem (involving the bi-curl) which can in principle be discretized via the Galerkin method. A serious drawback of this idea for practical computations is the fact that the standard finite element spaces are not solenoidal. Usually, spurious modes associated to the infinite-dimensional kernel appear and give rise to spectral pollution. This has been well documented and it is known to be a manifested problem whenever the underlying mesh is unstructured, [START_REF] Arnold | Finite element exterior calculus: from hodge theory to numerical stability[END_REF][START_REF] Boffi | Finite element approximation of eigenvalue problems[END_REF] and references therein.

Various ingenious methods, e.g. [START_REF] Bonito | Approximation of the eigenvalue problem for the time harmonic Maxwell system by continuous Lagrange finite elements[END_REF][START_REF] Bramble | The approximation of the Maxwell eigenvalue problem using a least-squares method[END_REF][START_REF] Buffa | Solving electromagnetic eigenvalue problems in polyhedral domains[END_REF][START_REF] Boffi | Computational models of electromagnetic resonators: analysis of edge element approximation[END_REF][START_REF] Boffi | Finite element approximation of eigenvalue problems[END_REF], capable of approximating the eigenvalues of (39) by means of the finite element method have been documented in the past. Let us apply the framework of Section 3 for finding eigenvalue bounds for M employing Lagrange finite elements on unstructured meshes. Convergence and absence of spectral pollution are guaranteed, as a consequence of Corollary 10 and Theorem 15.

Let {T h } h>0 be a family of shape-regular triangulations of Ω [START_REF] Ern | Theory and practice of finite elements[END_REF], where each

element K ∈ T h is a simplex with diameter h K such that h = max K∈T h h K . For r ≥ 1, let V r h = {v h ∈ [C 0 (Ω)] 3 : v h | K ∈ [P r (K)] 3 ∀K ∈ T h }, V r h,0 = {v h ∈ V r h : v h × n = 0 on ∂Ω} and set L h = V r h,0 × V r h ⊂ D(M).
Let ω 1 ≤ ω 2 ≤ . . . be the positive eigenvalues of M. The upper bounds ω + j and lower bounds ω - j reported below are found by fixing t ∈ R, solving (Z L h t ) numerically, and then applying [START_REF] Zimmermann | Variational bounds to eigenvalues of self-adjoint eigenvalue problems with arbitrary spectrum[END_REF].

The only hypothesis in the analysis carried out above ensuring that the ω ± j are close to ω j , is for the trial space to capture well the eigenfunctions in the graph norm of D(M). Therefore, as we have substantial freedom to choose these spaces and they constitute the simplest alternative, we have picked the Lagrange nodal elements. A direct application of Theorem 12 and classical interpolation estimates e.g. [START_REF] Ciarlet | The finite element method for elliptic problems[END_REF]Theorem 3.1.6], leads to convergence of the approximated eigenvalues and eigenspaces. Moreover, if the eigenspaces are regular, then the optimal convergence rates of order h 2r for eigenvalues and h r for eigenspaces can be proved.

This regularity assumption on the corresponding vector spaces can be formulated in different ways in order to suit the chosen algorithm. For the one we have employed here, if we wish to obtain a lower/upper bound for the j-eigenvalue to the left/right of a fixed t (and consequently obtain approximate eigenvectors) all the vectors of the sum of all eigenspaces up to j have to be regular. If by some misfortune, an intermediate eigenspace does not fullfill this requirement, then the algorithm will converge slowly. To circumvent this difficulty, the computational procedure can be modified in many ways. For instance, it can be allowed to split iteratively the initial interval, once it is clear that some accuracy can not be achieved after a fixed number of steps. 5.1. Orders of convergence on a cube. The eigenfunctions of (39) are regular in the interior of a convex domain. In this case, the Zimmermann-Mertins method for the resonant cavity problem achieves an optimal order of convergence in the context of finite elements.

Let Ω = Ω c = (0, π) Here {l, m, n} ⊂ N ∪ {0} and not two indices are allowed to vanish simultaneously. The vector α determines the multiplicity of the eigenvalue for a given triplet (l, m, n). That is, for example, ω 1 = √ 2 (the first positive eigenvalue) has multiplicity 3 corresponding to indices {(1, 1, 0), (0, 1, 1), (1, 0, 1)} each one of them contributing to one of the dimensions of the eigenspace. However, ω 2 = √ 3 (the second positive eigenvalue) corresponding to index {(1, 1, 1)} has multiplicity 2 determined by α on a plane.

In Figure 1 we have depicted the decrease in enclosure width and exact residual, for the lower bounds. 5.2. Benchmark eigenvalue bounds for the Fichera domain. In this next experiment we consider the region Ω = Ω F = (0, π) 3 \ [0, π/2] 3 . Some of the eigenvalues can be obtained by domain decomposition and the corresponding eigenfunctions are regular. For example, eigenfunctions on the cube of side π/2 can be assembled in the obvious fashion, in order to build eigenfunctions on Ω F . Therefore the set {±2 √ l 2 + m 2 + n 2 } where not two indices vanish simultaneously certainly lies inside σ(M). The first eigenvalue in this set is 2 √ 2. We conjecture that there are exactly 15 eigenvalues in the interval (0, 2 √ 2). Furthermore, we conjecture that the multiplicity counting of the spectrum in this interval is 1, 2, 3, 2, 1, 2, 1, 3.

ω + 2 -ω -
The table on the right of Figure 2 shows a numerical estimation of these eigenvalues. We have considered a mesh refined along the re-entrant edges as shown on the left side of this figure.

The slight numerical discrepancy shown in the table for the seemingly multiple eigenvalues appears to be a consequence of the fact that the meshes employed are not entirely symmetric with respect to permutation of the spacial coordinates.

3 ⊂ R 3 .

 33 The non-zero eigenvalues are ω = ± l 2 + m 2 + n 2 and the corresponding eigenfunctions are E(x, y, z) =   α 1 cos(lx) sin(my) sin(nz) α 2 sin(lx) cos(my) sin(nz) α 3 sin(lx) sin(my) cos(nz)

2 and ω + 2 -ω 2 ,

 222 for the computed bounds of the eigenvalue ω 2 = √ 3 by means of Lagrange elements of order r = 1, 2, 3. In this experiment we have chosen a sequence of unstructured tetrahedral mesh. The values for the slopes of the straight lines indicates that the enclosures obey the estimate of the form (43)|ω ± -ω| ≤ ch 2r , which is indeed the optimal convergence rate.

Figure 1 .

 1 Figure 1. Log-log graph associated to Ω c and ω 2 = √ 3. Vertical axis: enclosure width. Horizontal axis: maximum element size h. Here we have chosen Lagrange elements of order r = 1, 2, 3 on a sequence of unstructured meshes. Here we have chosen t = √ 2+ √ 3 2

  Convergence of the approximated local counting function 4.2. Convergence of local bounds for eigenvalues 4.3. Convergence for eigenfunctions 5. The finite element method for the Maxwell eigenvalue problem 5.1. Orders of convergence on a cube 5.2. Benchmark eigenvalue bounds for the Fichera domain Appendix A. A Comsol v4.3 LiveLink code operators: a method considered a few years ago by Zimmermann and Mertins

	1. Introduction
	2. Approximated local counting functions
	2.1. Optimal setting for detection of the spectrum
	2.2. Geometrical properties of the first approximated counting function
	2.3. Geometrical properties of the subsequent approximated counting
	functions
	2.4. Approximated eigenspaces
	3. Local bounds for eigenvalues
	3.1. The eigenvalue immediately to the left
	3.2. Further eigenvalues
	4. Convergence and error estimates
	4.1.

  The hypotheses ensure that the number of indices j ≤ n -(t) such that ν - j (t) ∈ J never exceeds m. Therefore this condition in the conclusion of the theorem is consistent. Letm(t) = max{m ∈ N : [t -d m (t), t + d m (t)] ⊂ J}.The hypothesis on L guarantees that (A 0 )-(A 1 ) hold true for m = m(t) and with

	m j=1 ε 2 j	1/2	< ε. By combining Lemma 11 and Theorem 12 and the fact that
	we can pick {w t j } m(t) j=1 ⊆ {w k } m k=1 , there exists εt > 0 small enough, such that (37)
	yields
	(38)	F j L (s) -d j (s) ≤	t -ν -1 (t) 2	∀j = 1, . . . , m and s ∈ J.

  'std1').feature.create('eigv', 'Eigenvalue'); model.study('std1').feature('eigv').set('shift', num2str(searchtau)); model.study('std1').feature('eigv').set('neigs', neval); std1.run; % STORING SOLUTION FOR POST PROCESSING [SZ,NDOFS,DATA,NAME,TYPE]= mphgetp(model,'solname','sol1');

	cons2.set('R', 3, 'E3');
	cons2.selection.set([2 5 7]);
	cons3=model.physics('w').feature.create('cons3', 'Constraint');
	cons3.set('R', 1, 'E1');
	cons3.set('R', 2, 'E2');
	cons3.selection.set([3 4 6]);
	% HOW MANY EIGENVALUES TO LOOK FOR AROUND t
	neval=3;
	% SOLVING THE MODEL
	std1=model.study.create('std1');
	model.study(% DISPLAYING SOLUTION
	for inde=1:neval,
	tauinv=(real(DATA(inde)));
	bd=parat+tauinv;
	if tauinv<0, disp(['lower= ',num2str(bd,10)]);
	else disp(['upper= ',num2str(bd,10)]);
	end
	disp(['DOF= ',num2str(NDOFS)])
	end

Acknowledgements References 1. Introduction

Below we examine in close detail the equivalence between two pollution-free techniques for numerical computation of eigenvalue enclosures for general self-adjoint

Acknowledgements

We kindly thank Michael Levitin and Stefan Neuwirth for their suggestions during the preparation of this manuscript. We kindly thank Université de Franche-Comté, University College London and the Isaac Newton Institute for Mathematical Sciences, for their hospitality. Funding was provided by MOPNET, the British-French project PHC Alliance (22817YA), the British Engineering and Physical Sciences Research Council (EP/I00761X/1) and the French Ministry of Research (ANR-10-BLAN-0101).

hex1.set('p',{'0' '0' '0' '0' 'pi' 'pi' 'pi' 'pi'; '0' '0' 'pi' 'pi' '0' '0' 'pi' 'pi'; '0' 'pi' 'pi' '0' '0' 'pi' 'pi' '0'}); hex2=geom1.feature.create('hex2', 'Hexahedron'); hex2.set('p',{'0' '0' '0' '0' 'pi/2' 'pi/2' 'pi/2' 'pi/2'; '0' '0' 'pi/2' 'pi/2' '0' '0' 'pi/2' 'pi/2'; '0' 'pi/2' 'pi/2' '0' '0' 'pi/2' 'pi/2' '0'}); dif1 = geom1.feature.create('dif1', 'Difference'); dif1.selection('input').set({'hex1'}); dif1.selection('input2').set({'hex2'}); geom1.run; %CREATING THE GEOMETRY model.mesh('mesh1').automatic(false); model.mesh('mesh1').feature('size').set('custom', 'on'); model.mesh('mesh1').feature('size').set('hmax', '. % FINITE ELEMENTS TO USE AND ORDER w.prop('ShapeProperty').set('shapeFunctionType', 'shlag'); w.prop('ShapeProperty').set('order', 3); % PHYSICS w.feature('wfeq1').set('weak',1 ,'(H3y-H2z)*(H3y_test-H2z_test)i*2*tt*(H3y-H2z)*E1_test+tt^2*E1*E1_test+(i*(H3y-H2z)-tt*E1)*E1t_test'); w.feature('wfeq1').set('weak',2 ,'(H1z-H3x)*(H1z_test-H3x_test)i*2*tt*(H1z-H3x)*E2_test+tt^2*E2*E2_test+(i*(H1z-H3x)-tt*E2)*E2t_test'); w.feature('wfeq1').set('weak',3 ,'(H2x-H1y)*(H2x_test-H1y_test)i*2*tt*(H2x-H1y)*E3_test+tt^2*E3*E3_test+(i*(H2x-H1y)-tt*E3)*E3t_test'); w.feature('wfeq1').set('weak',4 ,'(E3y-E2z)*(E3y_test-E2z_test)+ i*2*tt*(E3y-E2z)*H1_test+tt^2*H1*H1_test+((-i)*(E3y-E2z)-tt*H1)*H1t_test'); w.feature('wfeq1').set('weak',5 ,'(E1z-E3x)*(E1z_test-E3x_test)+ i*2*tt*(E1z-E3x)*H2_test+tt^2*H2*H2_test+((-i)*(E1z-E3x)-tt*H2)*H2t_test'); w.feature('wfeq1').set('weak',6 ,'(E2x-E1y)*(E2x_test-E1y_test)+ i*2*tt*(E2x-E1y)*H3_test+tt^2*H3*H3_test+((-i)*(E2x-E1y)-tt*H3)*H3t_test'); % BOUNDARY CONDITIONS cons1=model.physics('w').feature.create('cons1', 'Constraint'); cons1.set('R', 2, 'E2'); cons1.set('R', 3, 'E3'); cons1.selection.set([1 8 9]); cons2=model.physics('w').feature.create('cons2', 'Constraint'); cons2.set('R', 1, 'E1');