A Bespalov
email: a.bespalov@bham.ac.uk

S Nicaise
email: serge.nicaise@univ-valenciennes.fr

The BEM with graded meshes for the electric field integral equation on polyhedral surfaces

Keywords: electromagnetic scattering, electric field integral equation, Galerkin discretisation, boundary element method, Raviart-Thomas interpolation, anisotropic elements, graded mesh AMS Subject Classification: 65N38, 65N12, 78M15

We consider the variational formulation of the electric field integral equation on a Lipschitz polyhedral surface Γ. We study the Galerkin boundary element discretisations based on the lowest-order Raviart-Thomas surface elements on a sequence of anisotropic meshes algebraically graded towards the edges of Γ. We establish quasioptimal convergence of Galerkin solutions under a mild restriction on the strength of grading. The key ingredient of our convergence analysis are new componentwise stability properties of the Raviart-Thomas interpolant on anisotropic elements.

Introduction

In this paper, we study the Galerkin boundary element method (BEM) on graded meshes for numerical solution of the electric field integral equation (EFIE) on a Lipschitz polyhedral surface Γ in R 3 (i.e., Γ = ∂Ω, where Ω ⊂ R 3 is a Lipschitz polyhedron). The EFIE models the scattering of time-harmonic electromagnetic waves at a perfect conductor, and the Galerkin BEM is widely used in engineering practice for simulation of this physical phenomenon.

The Galerkin BEM considered in this paper employs div Γ -conforming lowest-order Raviart-Thomas surface elements to discretise the variational formulation of the EFIE (known as Rumsey's principle). This approach is referred to as the natural BEM for the EFIE (there exist other approaches, e.g., based on a stable mixed reformulation of Rumsey's principle, see [START_REF] Buffa | Boundary element methods for Maxwell's equations on non-smooth domains[END_REF]). Non-coercivity of the bilinear form in Rumsey's principle (due to the infinite-dimensional kernel of div Γ , cf. (2.1)) significantly complicates the convergence analysis of Galerkin schemes. This problem can be overcome by using appropriate decompositions of vector fields in order to isolate the kernel of div Γ (we refer to discussion in [18, Section 3], to an abstract theory in [START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF][START_REF] Buffa | Remarks on the discretization of some noncoercive operator with applications to heterogeneous Maxwell equations[END_REF], and we outline available techniques for constructing such decompositions in Section 4). These ideas have led to major advances in the convergence analysis and a priori error analysis of the BEM for the EFIE on (open and closed) Lipschitz surfaces, see [START_REF] Hiptmair | Natural boundary element methods for the electric field integral equation on polyhedra[END_REF][START_REF] Buffa | Boundary element methods for Maxwell's equations on non-smooth domains[END_REF][START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF][START_REF] Buffa | Boundary element methods for Maxwell transmission problems in Lipschitz domains[END_REF][START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF] for the h-version of the BEM and [START_REF]Natural p-BEM for the electric field integral equation on screens[END_REF][START_REF] Bespalov | Convergence of the natural hp-BEM for the electric field integral equation on polyhedral surfaces[END_REF][START_REF]The hp-BEM with quasi-uniform meshes for the electric field integral equation on polyhedral surfaces: a priori error analysis[END_REF][START_REF]Natural hp-BEM for the electric field integral equation with singular solutions[END_REF] for high-order methods (p-and hp-BEM). All these results, however, assume shape-regularity of the underlying meshes on Γ.

It is well-known that convergence rates of the h-BEM with quasi-uniform and shaperegular meshes are bounded by the poor regularity of solutions to the EFIE on non-smooth surfaces. For example, on a closed polyhedral surface Γ = ∂Ω, the solution may be only H ε (Γ)-regular (with a small ε > 0 in the case of non-convex polyhedron Ω, cf. [START_REF] Costabel | Singularities of electromagnetic fields in polyhedral domains[END_REF]Section 4.4.2]), and convergence rate of the h-BEM is only 1 2 + ε in this case, whereas in the case of smooth solutions the lowest-order h-BEM converges with the optimal rate of 3 2 (see [START_REF] Hiptmair | Natural boundary element methods for the electric field integral equation on polyhedra[END_REF]Theorem 8.2] and [START_REF]The hp-BEM with quasi-uniform meshes for the electric field integral equation on polyhedral surfaces: a priori error analysis[END_REF]Theorem 2.2]). Taking the cue from the h-BEM results for the Laplacian (see [START_REF] Petersdorff | Regularity of mixed boundary value problems in R 3 and boundary element methods on graded meshes[END_REF][START_REF] Petersdorff | Randwertprobleme der Elastizitätstheorie für Polyeder -Singularitäten und Approximation mit Randelementmethoden[END_REF]), we expect that an optimal convergence rate of the h-BEM for the EFIE can be recovered on the non-smooth surface Γ, if one employs the meshes that are appropriately graded towards the edges of Γ. These meshes contain highly anisotropic elements along the edges of Γ, and none of the results mentioned above is applicable in this case. Moreover, to the best of our knowledge, the quasi-optimality of the Galerkin h-BEM with graded meshes for the EFIE has not been studied in the literature, and with this paper we fill this theoretical gap.

In the next section, we introduce necessary notation and formulate the EFIE in its variational form. In Section 3, we construct graded meshes on Γ, introduce the boundary element space, and formulate the main result of the paper-Theorem 3.1-that establishes quasi-optimal convergence of Galerkin solutions on graded meshes. The proof of Theorem 3.1 follows the approach suggested in [START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF][START_REF] Buffa | Boundary element methods for Maxwell transmission problems in Lipschitz domains[END_REF], summarised in [START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF]Section 9.1], and extended to a general class of operators in [12, Section 3]. At the heart of this approach is the decomposition technique described in Section 4. Section 5 is instrumental in the construction of the corresponding discrete decomposition: here we establish new stability properties of the Raviart-Thomas interpolant of low-regular vector fields on anisotropic elements. In Section 6, we introduce the discrete decomposition and complete the proof of Theorem 3.1. An essential ingredient here is the projection operator Q h with enhanced approximation properties (see Proposition 6.1). The proof of Proposition 6.1 is given in Section 7.

The electric field integral equation

The variational formulation of the EFIE is posed on the Hilbert space

X = H -1/2 (div Γ , Γ) := {u ∈ H -1/2 (Γ); div Γ u ∈ H -1/2 (Γ)}.
Here, div Γ denotes the surface divergence operator, H -1/2 (Γ) is the dual space of H 1/2 (Γ)

(the tangential trace space of H 1 (Ω) on Γ, see [START_REF] Buffa | On traces for functional spaces related to Maxwell's equations, Part I: An integration by parts formula in Lipschitz polyhedra[END_REF][START_REF] Buffa | On traces for H(curl, Ω) in Lipschitz domains[END_REF]), and H -1/2 (Γ) is the dual space of H 1/2 (Γ). The space X is equipped with its graph norm • X . We refer to [START_REF] Buffa | On traces for functional spaces related to Maxwell's equations, Part I: An integration by parts formula in Lipschitz polyhedra[END_REF][START_REF]On traces for functional spaces related to Maxwell's equations, Part II: Hodge decompositions on the boundary of Lipschitz polyhedra and applications[END_REF][START_REF] Buffa | On traces for H(curl, Ω) in Lipschitz domains[END_REF][START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF] for definitions and properties of H -1/2 (div Γ , Γ) and other involved trace spaces. We also recall from [START_REF] Buffa | On traces for functional spaces related to Maxwell's equations, Part I: An integration by parts formula in Lipschitz polyhedra[END_REF][START_REF] Buffa | On traces for H(curl, Ω) in Lipschitz domains[END_REF] that X is the natural tangential trace space of H(curl, Ω).

In the present article, we use the same notation as in [START_REF]The hp-BEM with quasi-uniform meshes for the electric field integral equation on polyhedral surfaces: a priori error analysis[END_REF], where we recalled definitions of the full range of Sobolev spaces and differential operators needed for convergence analysis of the BEM for the EFIE (see Section 3.1 therein). In particular, we use a traditional notation for the Sobolev spaces (of scalar functions) H s , Hs (s ∈ [-1, 1]), H s 0 (s ∈ (0, 1]) and their norms on Lipschitz domains and surfaces (see [START_REF] Lions | Non-Homogeneous Boundary Value Problems and Applications I[END_REF][START_REF] Mclean | Strongly Elliptic Systems and Boundary Integral Equations[END_REF]). The norm and inner product in L 2 (D) = H 0 (D) on a domain or surface D will be denoted by • 0,D and (•, •) 0,D , respectively. The notation (•, •) 0,D will be used also for appropriate duality pairings extending the L 2 (D)-pairing for functions on D.

For vector fields we will use boldface symbols (e.g., u = (u 1 , u 2 )), and the spaces (or sets) of vector fields are also denoted in boldface (e.g.,

H s (D) = (H s (D)) 2 with D ⊂ R 2 ).
The norms and inner products in these spaces are defined componentwise. The notation for the Sobolev spaces of tangential vector fields on Γ follows [START_REF] Buffa | On traces for functional spaces related to Maxwell's equations, Part I: An integration by parts formula in Lipschitz polyhedra[END_REF][START_REF]On traces for functional spaces related to Maxwell's equations, Part II: Hodge decompositions on the boundary of Lipschitz polyhedra and applications[END_REF][START_REF] Buffa | On traces for H(curl, Ω) in Lipschitz domains[END_REF]. In particular, L 2 t (Γ) denotes the space of two-dimensional, tangential, square integrable vector fields on Γ. The norm and inner product in this space will be denoted by • 0,Γ and (•, •) 0,Γ , respectively, and we will also use (•, •) 0,Γ for appropriate duality pairings extending the L 2 t (Γ)-pairing for tangential vector fields on Γ. The similarity of this notation with the one for scalar functions should not lead to any confusion, as the meaning will always be clear from the context .

For a fixed wave number k > 0 and for a given source functional f ∈ X ′ , the variational formulation for the EFIE reads as: find a complex tangential field u ∈ X such that

a(u, v) := Ψ k div Γ u, div Γ v -k 2 Ψ k u, v = f, v ∀v ∈ X. (2.1)
Here, Ψ k (resp., Ψ k ) denote the scalar (resp., vectorial) single layer boundary integral operator on Γ for the Helmholtz operator -∆-k 2 , see [16, Section 4.1] (resp., [START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF]Section 5]).

To ensure the uniqueness of the solution to (2.1) we always assume that k 2 is not an electrical eigenvalue of the interior problem in Ω.

First, let us describe the construction of graded meshes on individual faces of Γ. Here, we follow [START_REF] Petersdorff | Regularity of mixed boundary value problems in R 3 and boundary element methods on graded meshes[END_REF]Section 3]. For simplicity, we can assume that all faces of Γ are triangles. On general polygonal faces the construction is similar, or one can first subdivide the polygon into triangles. On a triangular face F ⊂ Γ, we first draw three lines through the centroid and parallel to the sides of F . This makes F divided into three parallelograms and three triangles (see Figure 1). Each of the three parallelograms can be mapped onto the unit square Q = (0, 1) 2 by a linear transformation such that the vertex (0, 0) of Q is the image of a vertex of F . Analogously, each of the three sub-triangles can be mapped onto the unit triangle

T = {x = (x 1 , x 2 ); 0 < x 1 < 1, 0 < x 2 < x 1 } ⊂ Q such that the vertex (1, 1)
of T is the image of the centroid of F . Next, the graded mesh on Q (and hence on T ) is generated by the lines

x 1 = i N β , x 2 = j N β , i, j = 0, 1, . . . , N.
Here, β ≥ 1 is the grading parameter, and N ≥ 1 corresponds to the level of refinement.

Mapping each cell of these meshes back onto the face F , we obtain a graded mesh of triangles and parallelograms on F (see Figure 1). Note that the diameter of the largest element of this mesh is proportional to βN -1 . Hence, h = 1/N defines the mesh parameter, and we will denote by T = {∆ β h } a family of graded meshes ∆ β h = {K; ∪ K = Γ} generated on Γ by following the procedure described above.

Let us now introduce the boundary element space X h . It is known that Raviart-Thomas surface elements provide an affine equivalent family of div Γ -conforming finite elements under the Piola transformation, see [START_REF] Brezzi | Mixed and Hybrid Finite Element Methods[END_REF]Section III.3]. We will write RT 0 (K) for the local lowest-order Raviart-Thomas space on a generic (triangular or quadrilateral) element K, and we denote by X h = RT 0 (∆ β h ) the corresponding space of div Γ -conforming boundary elements over the graded mesh ∆ β h . The following theorem states the unique solvability and quasi-optimal convergence of the Galerkin BEM on graded meshes for the EFIE. Theorem 3.1 There exists h 0 < 1 such that for any f ∈ X ′ and for any graded mesh ∆ β h with h ≤ h 0 and β ∈ [1, 3), the Galerkin boundary element discretisation of (2.1) admits a unique solution u h ∈ X h and the h-version of the Galerkin BEM on graded meshes ∆ β h converges quasi-optimally, i.e.,

u -u h X ≤ C inf v∈X h u -v X , (3.1)
where the constant C may depend only on the geometry of Γ and the grading parameter β.

The proof of Theorem 3.1 relies on an abstract theory for analysing convergence of Galerkin discretisations for non-coercive variational problems like (2.1). This theory was developed in [START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF], [START_REF] Buffa | Boundary element methods for Maxwell transmission problems in Lipschitz domains[END_REF], [START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF]Section 9.1], and in [START_REF] Buffa | Remarks on the discretization of some noncoercive operator with applications to heterogeneous Maxwell equations[END_REF]Section 3]. In particular, it follows from the latter article that in order to prove Theorem 3.1 we need to establish the following properties:

Figure 1: Graded mesh on the triangular face F ⊂ Γ. The triangular (resp., parallelogram) block of elements T F (resp., Q F ) is the image of the graded mesh on the unit triangle T (resp., the unit square Q).

(A) the existence of a stable direct decomposition X = V ⊕ W such that a| V×V and -a| W×W are both X-coercive, and a| V×W and a| W×V are both compact;

(B) the existence of the corresponding discrete decomposition

X h = V h +W h , W h ⊂ W,
that is uniformly stable with respect to the mesh parameter h;

(C) the gap property

sup v h ∈V h inf v∈V v -v h X v h X ≤ ε(h) with ε(h) → 0 as h → 0. (3.2)
We will prove Theorem 3.1 by verifying these properties in Sections 4 and 6 below.

Remark 3.1 Theorem 3.1 remains valid if Γ is a piecewise plane orientable open surface (see [START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF] for the problem formulation and the underlying tangential trace spaces in this case).

The proof repeats the arguments in Sections 4 and 6 below by using a specific construction of the decomposition X = V ⊕ W as described in [9, Section 3].

Remark 3.2 If some information about the regularity of the solution u to (2.1) is available, then convergence result of Theorem 3.1 translates into an a priori error estimate in the natural X-norm. For scattering problems with sufficiently smooth source functional f (e.g., with f representing the excitation by an incident plane wave), the regularity of the solution depends only on the geometry of Γ. In particular, nonsmoothness of Γ leads to singularities in the solution of the EFIE, severely affecting convergence rates of the h-BEM on shaperegular meshes. However, similar to the case of the Laplacian in [START_REF] Petersdorff | Regularity of mixed boundary value problems in R 3 and boundary element methods on graded meshes[END_REF][START_REF] Petersdorff | Randwertprobleme der Elastizitätstheorie für Polyeder -Singularitäten und Approximation mit Randelementmethoden[END_REF], by employing the graded meshes with sufficiently large grading parameter β (depending on the strength of singularities in u) one may hope to recover the optimal convergence rate (i.e., the rate of the h-BEM on quasi-uniform meshes in the case of a smooth solution). The main question here is whether the restriction on the grading parameter β that guarantees quasi-optimality of the Galerkin BEM in Theorem 3.1, is sufficient for recovering this optimal convergence. We will address this issue in the forthcoming article [START_REF] Bespalov | A priori error analysis of the BEM with graded meshes for the electric field integral equation on polyhedra[END_REF].

Throughout the paper, C, C 1 , etc. denote generic positive constants that are independent of the mesh parameter h and involved functions but may depend on the geometry of Γ and the grading parameter β. We will also write a b and A ≃ B, which means the existence of generic positive constants C, C 1 , C 2 such that a ≤ Cb and C 1 B ≤ A ≤ C 2 B, respectively.

Decomposition technique

Let us address property (A) in Section 3. One way to obtain a suitable decomposition is to employ the L 2 t (Γ)-orthogonal Hodge decomposition of X, cf. [START_REF]On traces for functional spaces related to Maxwell's equations, Part II: Hodge decompositions on the boundary of Lipschitz polyhedra and applications[END_REF]. In the context of the h-BEM on shape-regular meshes, this idea was successfully exploited in [START_REF] Buffa | Boundary element methods for Maxwell's equations on non-smooth domains[END_REF][START_REF] Hiptmair | Natural boundary element methods for the electric field integral equation on polyhedra[END_REF][START_REF] Buffa | Boundary element methods for Maxwell transmission problems in Lipschitz domains[END_REF][START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF]. However, in the case of non-smooth surfaces, the regularity of surface gradients in the V-component of the decomposition may be poor, and this causes substantial technical difficulties in the analysis of the p-and hp-BEM. For the p-BEM on plane open screens, a modification of the above strategy was suggested in [START_REF]Natural p-BEM for the electric field integral equation on screens[END_REF], where we consistently used the H-1/2 -inner product and proved H-1/2 -orthogonality of the Hodge decomposition. Unfortunately, these ideas do not generalise immediately to closed polyhedral surfaces or piecewise plane open screens, neither to the hp-BEM with quasi-uniform meshes, cf. [START_REF] Bespalov | On the convergence of the hp-BEM with quasiuniform meshes for the electric field integral equation on polyhedral surfaces[END_REF]. When attempting to use the Hodge decomposition for convergence analysis of the h-BEM on graded meshes, poor regularity of vector fields in the V-component leads to severe restrictions on the grading parameter β.

An alternative technique employs a regularising projection R : X → X to construct a decomposition of X with enhanced regularity of the V-component (see [START_REF] Hiptmair | Coupling of finite elements and boundary elements in electromagnetic scattering[END_REF], [START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF]Section 3], and [START_REF] Buffa | Remarks on the discretization of some noncoercive operator with applications to heterogeneous Maxwell equations[END_REF]Section 4.3.1]). The projection R is defined by employing the H 1 (Ω)-regular vector potentials from the following lemma.

: curl H(curl, Ω) → H 1 (Ω) such that curl L Φ = Φ for all Φ ∈ curl H(curl, Ω).
Let u ∈ X. The regularising projection R : X → X is defined as follows:

R u := (L Φ) × ν | Γ ,
where L is from Lemma 4.1, ν denotes the unit outward normal to Ω, Φ := ∇w, and w ∈ H 1 (Ω) is the solution to the problem

-∆w = 0 in Ω, ∇w • ν = div Γ u on Γ.
The fact that Σ div Γ u dS = 0 for each connected component Σ of Γ guarantees Φ ∈ curl H(curl, Ω), and Lemma 4.1 can be applied. By using elliptic lifting theorems, trace theorems, and the continuity of L, we conclude:

∃C = C(Γ) > 0 such that R u H 1/2 ⊥ (Γ) ≤ C div Γ u H -1/2 (Γ) ∀ u ∈ X, (4.1)
where [START_REF] Buffa | On traces for functional spaces related to Maxwell's equations, Part I: An integration by parts formula in Lipschitz polyhedra[END_REF][START_REF] Buffa | On traces for H(curl, Ω) in Lipschitz domains[END_REF]. By construction of R, we have on

H 1/2 ⊥ (Γ) ⊂ X is the rotated tangential trace space of H 1 (Ω) on Γ = ∂Ω, see
Γ div Γ R u = div Γ u for any u ∈ X, (4.2)
that is R 2 = R. Now we can define the decomposition

X = V ⊕ W with V := R(X) ⊂ H 1/2 ⊥ (Γ) and W := (Id -R)(X). (4.3)
Decomposition (4.3) was used in [START_REF] Bespalov | Convergence of the natural hp-BEM for the electric field integral equation on polyhedral surfaces[END_REF] to prove unique solvability and quasi-optimal convergence of the hp-BEM with locally variable polynomial degrees on shape-regular meshes. As we will see in this paper, the same decomposition technique can be used effectively in the analysis of the h-BEM on graded meshes. By (4.2) we conclude that W comprises div Γ -free vector fields. Stability of decomposition (4.3) follows from inequality (4.1) and the continuous embedding

H 1/2 ⊥ (Γ) ֒→ X. Furthermore, the embedding V ֒→ L 2 t (Γ) is compact by (4.1
) and Rellich's theorem. Thus, thanks to the [START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF]Lemmas 8,[START_REF]-conforming p-interpolation operator in two dimensions[END_REF]), the X-coercivity of a| V×V and a| W×W is proved by the same arguments as in [START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF]proof of Theorem 3.4]. The compactness of a| V×W and a| W×V is due to the continuity of [START_REF] Buffa | Galerkin boundary element methods for electromagnetic scattering[END_REF]Lemma 9]). This proves (A).

H -1/2 (Γ)-coercivity (resp., H -1/2 ⊥ (Γ)-coercivity) of Ψ k (resp. Ψ k ) (see
Ψ k : H -1/2 ⊥ (Γ) → H 1/2 ⊥ (Γ) and the compactness of the embedding V ֒→ H -1/2 ⊥ (Γ) (see
Before we can define the discrete counterpart of decomposition (4.3), we need to find a suitable projector onto the space of Raviart-Thomas surface elements. Stability of the discrete decomposition will follow from stability properties of the Raviart-Thomas interpolation on anisotropic elements, which is the subject of the next section.

Raviart-Thomas interpolation on anisotropic elements

In this section, we establish new stability properties of the Raviart-Thomas interpolant on anisotropic elements. We will also prove the corresponding interpolation error estimates.

In the context of the finite element method, the analysis of interpolation operators on anisotropic elements can be found in [START_REF] Apel | Anisotropic Finite Elements: Local Estimates and Applications[END_REF][START_REF] Acosta | Error estimates for Raviart-Thomas interpolation of any order on anisotropic tetrahedra[END_REF]. For the Raviart-Thomas interpolation (see, e.g., [1, Section 3]), the main idea is to study componentwise stability of the interpolant on a reference element K. For sufficiently regular vector fields, this study relies on the fact that the standard (scalar) trace operator is well defined for functions in W 1,p ( K) for any p > 1. In our BEM application, however, the stability result is needed for low-regular vector fields living in fractional Sobolev spaces

H s ( K) ∩ H(div, K) with 0 < s ≤ 1/2.
For such vector fields, the trace of the normal component only exists in a weak sense. Therefore, instead of the standard trace argument used in [START_REF] Acosta | Error estimates for Raviart-Thomas interpolation of any order on anisotropic tetrahedra[END_REF], we use Green's formula (5.3) and consistently employ the anisotropic seminorms defined below. More precisely, for s ∈ (0, 1/2] we introduce the H s -seminorms of anisotropic type. On the reference square Q = (0, 1) 2 these are defined as follows:

|u| 2 AH s 1 ( Q) = 1 0 |u(•, x 2 )| 2 H s (0,1) dx 2 , |u| 2 AH s 2 ( Q) = 1 0 |u(x 1 , •)| 2 H s (0,1) dx 1 .
These definitions are meaningful for all u ∈ H s ( Q) due to [25, Theorem 10.2] which yields that

u H s ( Q) ≃ u 0, Q + |u| AH s 1 ( Q) + |u| AH s 2 ( Q) . (5.1)
On the reference triangle T = {(x 1 , x 2 ); 0 < x 1 < 1, 0 < x 2 < x 1 }, the following seminorms

|u| 2 AH s 1 ( T ) = 1 0 |u(•, x 2 )| 2 H s (x 2 ,1) dx 2 , |u| 2 AH s 2 ( T ) = 1 0 |u(x 1 , •)| 2 H s (0,x 1 ) dx 1
are also well defined for all u ∈ H s ( T ). Indeed, by Theorem 1.4.3.1 of [START_REF] Grisvard | Elliptic Problems in Nonsmooth Domains[END_REF] there exists a continuous linear operator (called the extension operator) E :

H s ( T ) → H s ( Q) such that E u| T = u, ∀u ∈ H s ( T ).
Hence, applying (5.1) to E u, we have

|u| AH s 1 ( T ) + |u| AH s 2 ( T ) u H s ( T ) . (5.2)
We recall that the Raviart-Thomas interpolant Π RT u is well-defined for any u ∈ H s (K) (s > 0) such that div u ∈ L 2 (K), where K is any triangle or rectangle. Indeed, for such vector fields the following Green's formula has a meaning (see, e.g., [7, Lemma 2.1])

(u, ∇ϕ) 0,K + K div u ϕ = (u • n, ϕ) 0,∂K , ∀ϕ ∈ H 1-ε (K), (5.3)
with ε ∈ (0, s) and n denoting the outward unit normal to ∂K. Hence, taking ϕ ∈ H 1-ε (K) such that ϕ = 1 on the edge e ⊂ ∂K and ϕ = 0 on ∂K\e, we can define

(u • n, 1) 0,e := (u • n, ϕ) 0,∂K .
In what follows, we will denote by Π RT the Raviart-Thomas interpolation operator on the reference element K ( K = Q or T ).

The reference square

On the reference square Q, we denote by e 1 and e 3 the edges parallel to the x 1 -axis, and by e 2 and e 4 the edges parallel to the x 2 -axis (see Figure 2). We recall from [START_REF] Brezzi | Mixed and Hybrid Finite Element Methods[END_REF] that the lowest order Raviart-Thomas elements on Q are defined as

RT 0 ( Q) = (a + cx 1 , b + dx 2 ) ⊤ ; a, b, c, d ∈ R ,
(i) For all u ∈ H s ( Q) with s > 1/2, we have for l = 1, 2 ( Π RT u) l 0, Q u l H s ( Q) . (5.4) (ii) For all u ∈ H s ( Q) ∩ H(div, Q) with 0 < s ≤ 1/2, we have for l = 1, 2 ( Π RT u) l 0, Q u l H s ( Q) + |u l+1 | AH s l+1 ( Q) + div u 0, Q (mod. 2).
(5.5)

Proof. By symmetry, it will suffice to prove both statements for l = 1. One has

Π RT u = (a + cx 1 , b + dx 2 ) ⊤ with a, b, c, d ∈ R. Hence ( Π RT u) 1 0, Q = a + cx 1 0, Q |a| + |c|,
and it remains to estimate |a| and |c|. Observe that

a = -(u • n, 1) 0, e 4 and a + c = (u • n, 1) 0, e 2 .
Again by symmetry, it will suffice to estimate

(u • n, 1) 0, e 4 . If u ∈ H s ( Q) with s > 1/2
, then the trace of u 1 on e 4 is well-defined, and

|a| = |(u • n, 1) 0, e 4 | = |(u 1 , 1) 0, e 4 | u 1 H s-1/2 ( e 4 ) u 1 H s ( Q) .
This proves statement (i) (for l = 1). Now, let us consider u ∈ H s ( Q) ∩ H(div, Q) for 0 < s ≤ 1/2. In order to estimate (u • n, 1) 0, e 4 in this case, we fix a function ϕ ∈ H 1-ε ( Q), ε ∈ (0, s), such that ϕ = 1 on e 4 and ϕ = 0 on ∂ Q\ e 4 . Then by Green's formula (5.3) we have

(u • n, 1) 0, e 4 = (u, ∇ϕ) 0, Q + Q divu ϕ = (u 1 , ∂ 1 ϕ) 0, Q + (u 2 , ∂ 2 ϕ) 0, Q + Q divu ϕ. (5.6) For the term (u 1 , ∂ 1 ϕ) 0, Q , we first use a standard duality argument, |(u 1 , ∂ 1 ϕ) 0, Q | ≤ u 1 H ε ( Q) ∂ 1 ϕ H -ε ( Q) ,
and by the continuity property of

∂ 1 : H 1-ε ( Q) → H -ε ( Q) (see Theorem 1.4.4.6 of [21]) we find |(u 1 , ∂ 1 ϕ) 0, Q | u 1 H s ( Q) ϕ H 1-ε ( Q) .
(5.7)

The term (u 2 , ∂ 2 ϕ) 0, Q requires more subtle analysis. First by (5.1) and Fubini's theorem, we can write (hereafter, I = (0, 1))

(u 2 , ∂ 2 ϕ) 0, Q = 1 0 (u 2 (x 1 , •), ∂ 2 ϕ(x 1 , •)) 0, I dx 1 . (5.8)
Now we use a density argument to show that

1 0 h(x 1 )(1, ∂ 2 ϕ(x 1 , •)) 0, I dx 1 = 0, (5.9)
for all h ∈ L 2 ( I). Indeed, fix a sequence of smooth function ϕ n such that

ϕ n → ϕ in H 1-ε ( Q) as n → ∞.
Then for all x 1 ∈ I = (0, 1), we have

(1, ∂ 2 ϕ n (x 1 , •)) 0, I = 1 0 ∂ 2 ϕ n (x 1 , x 2 ) dx 2 = ϕ n (x 1 , 1) -ϕ n (x 1 , 0).
For any h ∈ L 2 ( I), multiplying this identity by h and integrating the result in x 1 ∈ (0, 1), we obtain

1 0 h(x 1 )(1, ∂ 2 ϕ n (x 1 , •)) 0, I dx 1 = e 3 hϕ n - e 1 hϕ n .
Hence, as n → ∞ we find that

1 0 h(x 1 )(1, ∂ 2 ϕ(x 1 , •)) 0, I dx 1 = e 3 hϕ - e 1 hϕ,
which proves (5.9) by recalling that ϕ = 0 on e 1 and e 3 .

Coming back to (5.8) and using (5.9), we have

(u 2 , ∂ 2 ϕ) 0, Q = 1 0 u 2 (x 1 , •) -M I (u 2 (x 1 , •)), ∂ 2 ϕ(x 1 , •) 0, I dx 1 ,
where

M I (u 2 (x 1 , •)) = 1 0 u 2 (x 1 , x 2 ) dx 2 is the mean of u 2 (x 1 , •) on I = (0, 1) (clearly M I (u 2 (x 1 , •)) ∈ L 2 ( I))
. At this stage we first use a duality argument and then again Theorem 1.4.4.6 of [START_REF] Grisvard | Elliptic Problems in Nonsmooth Domains[END_REF] to obtain

|(u 2 , ∂ 2 ϕ) 0, Q | ≤ 1 0 u 2 (x 1 , •) -M I (u 2 (x 1 , •) H ε (0,1) ∂ 2 ϕ(x 1 , •) H -ε (0,1) dx 1 1 0 u 2 (x 1 , •) -M I (u 2 (x 1 , •) H ε (0,1) ϕ(x 1 , •) H 1-ε (0,1) dx 1 . (5.10)
We use Friedrichs' inequality to estimate

u 2 (x 1 , •) -M I (u 2 (x 1 , •) H ε (0,1) ≤ u 2 (x 1 , •) -M I (u 2 (x 1 , •) H s (0,1) |u 2 (x 1 , •)| H s (0,1) .
Using this estimate in (5.10) and applying the Cauchy-Schwarz inequality we arrive at

|(u 2 , ∂ 2 ϕ) 0, Q | |u 2 | AH s 2 ( Q) 1 0 ϕ(x 1 , •) 2 H 1-ε (0,1) dx 1 1 2 .
(5.11)

The last term on the right-hand side of (5.6) is estimated by applying the Cauchy-Schwarz inequality:

Q div uϕ ≤ div u 0, Q ϕ 0, Q .
Using this estimate and inequalities (5.7), (5.11) in (5.6), we obtain (5.5) (for l = 1). ✷

Corollary 5.1 (i) For all u ∈ H s ( Q) with s > 1/2, we have for l = 1, 2 u l -( Π RT u) l 0, Q |u l | H s ( Q) . (ii) For all u ∈ H s ( Q) ∩ H(div, Q) with 0 < s ≤ 1/2, we have for l = 1, 2 u l -( Π RT u) l 0, Q |u l | H s ( Q) + |u l+1 | AH s l+1 ( Q) + div u 0, Q (mod. 2).
Proof. It is sufficient to prove only statement (ii). For l = 1, we take

u = u-(M Q u 1 , 0) ⊤ , where M Q u 1 = Q u 1 . One has u -Π RT u = u -Π RT u, div u = div u,
and

|u 1 | H 1/2 ( Q) = |ũ 1 | H 1/2 ( Q) ≃ ũ1 H 1/2 ( Q) .
The assertion then follows by applying estimate (5.5) to u. The proof is analogous for

l = 2. ✷ Corollary 5.2 Let s ∈ (0, 1/2]. For all u ∈ H s ( Q) such that div u ∈ R, we have for l = 1, 2 ( Π RT u) l 0, Q u l H s ( Q) + |u l+1 | AH s l+1 ( Q) (mod. 2) (5.12) and u l -( Π RT u) l 0, Q |u l | H s ( Q) + |u l+1 | AH s l+1 ( Q) (mod. 2).
(5.13)

Proof. Note that the function ϕ ∈ H 1-ε ( Q) in the proof of Theorem 5.1 (ii) can be chosen to have zero average on Q, i.e., Q ϕ = 0. Indeed, if this is not the case then we

fix ψ ∈ C ∞ 0 ( Q) such that Q ψ = 1, and consider ϕ -Q ϕ ψ ∈ H 1-ε ( Q)
that has the same values on ∂ Q as ϕ and also has zero average on Q. Then, for div u ∈ R, the last term on the right-hand side of (5.6) vanishes, and (5.12) (resp., (5.13)) follows by the same arguments as in the proof of Theorem 5.1 (ii) (resp., Corollary 5.1 (ii)). ✷

The reference triangle

On the references triangle T , we denote by e 1 (resp. e 2 ) the edge on the x 1 -axis (resp., parallel to the x 2 -axis), and by e 3 the oblique edge (see Figure 2). We recall from [START_REF] Brezzi | Mixed and Hybrid Finite Element Methods[END_REF] that the lowest order Raviart-Thomas elements on T are defined as

RT 0 ( T ) = (a, b) ⊤ + c(x 1 , x 2 ) ⊤ ; a, b, c ∈ R ,
and that the associated degrees of freedom are given by

e i u • n ds, i = 1, 2, 3. Theorem 5.2 (i) For all u ∈ H s ( T ) ∩ H(div, T ) with s > 1/2, we have for l = 1, 2 ( Π RT u) l 0, T u l H s ( T ) + div u 0, T . (5.14) (ii) For all u ∈ H s ( T ) ∩ H(div, T ) with 0 < s ≤ 1/2, we have for l = 1, 2 ( Π RT u) l 0, T u l H s ( T ) + |u l+1 | AH s l+1 ( T ) + div u 0, T (mod. 2). (5.15)
Proof. We will prove both statements for l = 1. The proof is analogous in the case l = 2.

One has

Π RT u = (a, b) ⊤ + c(x 1 , x 2 ) ⊤ with a, b, c ∈ R. Hence ( Π RT u) 1 0, T = a + cx 1 0, T |a| + |c| ≤ |a + c| + 2|c|. Observe that 2c = div Π RT u = 2 T div u.
By the Cauchy-Schwarz inequality this yields

|c| div u 0, T .
Therefore in order to estimate ( Π RT u) 1 0, T one needs to bound |a

+ c| = |(u • n, 1) 0, e 2 |.
If u ∈ H s ( T ) with s > 1/2, then we can use trace theorem to estimate

|a + c| = |(u • n, 1) 0, e 2 | = |(u 1 , 1) 0, e 2 | u 1 H s ( T ) .
This proves statement (i) (for l = 1). Now, let us consider u ∈ H s ( T ) ∩ H(div, T ) for 0 < s ≤ 1/2. In order to estimate (u • n, 1) 0, e 2 in this case, we fix a function ϕ ∈ H 1-ε ( T ), ε ∈ (0, s), such that ϕ = 1 on e 2 and ϕ = 0 on ∂T \ e 2 . Then by Green's formula (5.3) we have

(u • n, 1) 0, e 2 = (u, ∇ϕ) 0, T + T divu ϕ = (u 1 , ∂ 1 ϕ) 0, T + (u 2 , ∂ 2 ϕ) 0, T + T divu ϕ. (5.16)
As before (cf. (5.7)), it is easy to show that

|(u 1 , ∂ 1 ϕ) 0, T | u 1 H s ( T ) ϕ H 1-ε ( T ) .
(5.17)

The term (u 2 , ∂ 2 ϕ) 0, T is treated as in the case of the reference square except that the interval I = (0, 1) in the x 2 -variable will be replaced by the interval I(x 1 ) := (0, x 1 ), but the arguments remain mainly unchanged because x 1 < 1. More precisely, first by (5.2) and Fubini's theorem, we may write

(u 2 , ∂ 2 ϕ) 0, T = 1 0 (u 2 (x 1 , •), ∂ 2 ϕ(x 1 , •)) 0,I(x 1 ) dx 1 . (5.18)
With a property similar to (5.9) we deduce that

(u 2 , ∂ 2 ϕ) 0, T = 1 0 u 2 (x 1 , •) -M I(x 1 ) (u 2 (x 1 , •)), ∂ 2 ϕ(x 1 , •) 0,I(x 1 ) dx 1 ,
where

M I(x 1 ) (u 2 (x 1 , •)) = 1 x 1 x 1 0 u 2 (x 1 , x 2 ) dx 2 is the mean of u 2 (x 1 , •) on I(x 1 ) = (0, x 1 ). Then a duality argument yields |(u 2 , ∂ 2 ϕ) 0, T | ≤ 1 0 u 2 (x 1 , •) -M I(x 1 ) (u 2 (x 1 , •)) H ε (I(x 1 )) ∂ 2 ϕ(x 1 , •) H -ε (I(x 1 )) dx 1 . (5.19)
Using a scaling argument and Friedrichs' inequality we estimate

u 2 (x 1 , •) -M I(x 1 ) (u 2 (x 1 , •)) H ε (I(x 1 )) ≤ u 2 (x 1 , •) -M I(x 1 ) (u 2 (x 1 , •)) H s (I(x 1 )) ≤ C|u 2 (x 1 , •)| H s (I(x 1 )) , (5.20)
with C > 0 independent of x 1 ∈ (0, 1).

For the second factor in the integrand in (5.19), in order to apply Theorem 1.4.4.6 of [START_REF] Grisvard | Elliptic Problems in Nonsmooth Domains[END_REF] on a fixed domain, we first notice that ϕ(x 1 , •) ∈ H 1-ε 0 (I(x 1 )) a. e. in (0, 1) ∋ x 1 . Therefore, for almost all x 1 ∈ (0, 1) there exists a sequence of functions ϕ n ∈ C ∞ 0 (I(x 1 )), n = 1, 2, . . ., such that ϕ n → ϕ(x 1 , •) in H 1-ε 0 (I(x 1 )) as n → ∞ and

x 1 0 ∂ϕ n (x 2 ) dx 2 = 0 ∀n = 1, 2, . . . .

Using a scaling argument, we deduce that

∂ϕ n H -ε (I(x 1 )) = sup v∈H ε (I(x 1 )) v =0 x 1 0 ∂ϕ n v dx 2 v H ε (I(x 1 )) ≤ x -1 2 +ε 1 sup v∈H ε (0,1) v =0 1 0 ∂ ϕ n v d x 2 | v| H ε (0,1) .
As ∂ ϕ n has zero average, we can estimate

∂ϕ n H -ε (I(x 1 )) ≤ x -1 2 +ε 1 sup v∈H ε (0,1) v =0, 1 0 v=0 1 0 ∂ ϕ n v d x 2 | v| H ε (0,1)
x

-1 2 +ε 1 sup v∈H ε (0,1) v =0, 1 0 v=0 1 0 ∂ ϕ n v d x 2 v H ε (0,1) ≤ x -1 2 +ε 1 ∂ ϕ n H -ε (0,1) .
Hence, by Theorem 1.4.4.6 of [START_REF] Grisvard | Elliptic Problems in Nonsmooth Domains[END_REF] we prove that

∂ϕ n H -ε (I(x 1 )) x -1 2 +ε 1 ϕ n H 1-ε (0,1) x -1 2 +ε 1 | ϕ n | H 1-ε (0,1) .
Mapping back to the interval I(x 1 ) = (0, x 1 ) we have

∂ϕ n H -ε (I(x 1 )) |ϕ n | H 1-ε (I(x 1 )) .
As n → ∞ we find

∂ 2 ϕ(x 1 , •) H -ε (I(x 1 )) ≤ C 1 |ϕ(x 1 , •)| H 1-ε (I(x 1 )) a. e. on (0, 1) ∋ x 1 , (5.21)
with C 1 > 0 independent of x 1 .

Using estimates (5.20) and (5.21) in (5. [START_REF] Buffa | Boundary element methods for Maxwell transmission problems in Lipschitz domains[END_REF]) we arrive at

(u 2 , ∂ 2 ϕ) 0, T 1 0 |u 2 (x 1 , •)| H s (0,x 1 ) ϕ(x 1 , •) H 1-ε (0,x 1 ) dx 1 .
Then the Cauchy-Schwarz inequality yields

(u 2 , ∂ 2 ϕ) 0, T |u 2 | AH s 2 ( T ) 1 0 ϕ(x 1 , •) 2 H 1-ε (0,x 1 ) dx 1 1 2 . ( 5.22)
The last term on the right-hand side of (5.16) is estimated by using the Cauchy-Schwarz inequality:

T divu ϕ ≤ div u 0, T ϕ 0, T .

Using this estimate and inequalities (5.17), (5.22) in (5.16), we obtain (5.15) (for l=1). ✷ Corollary 5.3 (i) For all u ∈ H s ( T ) ∩ H(div, T ) with s > 1/2, we have for l = 1, 2

u l -( Π RT u) l 0, T |u l | H s ( T ) + div u 0, T . (ii) For all u ∈ H s ( T ) ∩ H(div, T ) with 0 < s ≤ 1/2, we have for l = 1, 2 u l -( Π RT u) l 0, T |u l | H s ( T ) + |u l+1 | AH s l+1 ( T ) + div u 0, T (mod. 2
). The proof of this statement is similar to the proof of Corollary 5.1.

Counterexample 5.1 Here we provide a counterexample which demonstrates that for lowregular vector fields the terms |u l+1 | AH s l+1 ( Q) in (5.5) and |u l+1 | AH s l+1 ( T ) in (5.15) cannot be omitted. This is in contrast to the case of sufficiently-regular fields in (5.4) and (5.14) (see also Lemma 3.3 in [1]). In particular, if we assume that

( Π RT u) 2 0, T u 2 H 1/2 ( T ) + div u 0, T (5.23)
for all u ∈ H 1 2 ( T ) ∩ H(div, T ), then we will arrive at a contradiction. Indeed, inspired by Example 2.6 of [3], we define on T × (0, 1)

v ε (x 1 , x 2 , x 3 ) = (x 1 -1)w ε (x 2 , x 3 ) with w ε (x 2 , x 3 ) = min 1, ε log log e r
for any ε > 0.

Here, r = (x 2 2 +x 2 3 )

1 2 , and e is the Euler number. Taking ∇v ε ×ν on T (here ν = (0, 0, -1)), we find a divergence-free vector field

u ε (x 1 , x 2 ) = (1 -x 1 )∂ 2 w ε (x 2 , 0), w ε (x 2 , 0) ⊤ .
Simple calculations show that Π RT u ε = (0, 1) ⊤ , and by trace theorem we have

(u ε ) 2 H 1/2 ( T ) = w ε (x 2 , 0) H 1/2 ( T ) w ε H 1 ( T ×(0,1)) . Since w ε H 1 ( T ×(0,1)) → 0 as ε → 0 (see [3, Example 2.6]), we conclude that (u ε ) 2 → 0 in H 1 2 (
T ) as ε → 0. This seems to contradict (5.23) but not directly because the first component of u ε is not in H 1 2 ( T ). Hence, in order to arrive at a contradiction, we need to show that if (5.23) holds for all u ∈ H 1 2 ( T ) ∩ H(div, T ), then u ε satisfies (5.23) (with a constant indepedent of ε). Indeed, for a fixed ε > 0, as w ε (•, 0) ∈ H 1 (0, 1), we can consider a sequence of smooth functions

w n ∈ C ∞ ([0, 1]) such that w n → w ε (•, 0) in H 1 (0, 1) as n → ∞.
Then we define

u n (x 1 , x 2 ) = ((1 -x 1 )∂w n (x 2 ), w n (x 2 )) ⊤ .
One has u n ∈ H 1 2 ( T ) and div u n = 0. Moreover, Π RT u n → Π RT u ε as n → ∞. Therefore, applying estimate (5.23) to u n and letting n → ∞, we conclude that u ε satisfies (5.23).

Remark 5.1 By Counterexample 5.1 we can easily show that a result similar to Lemma 3.3 in [START_REF] Acosta | Error estimates for Raviart-Thomas interpolation of any order on anisotropic tetrahedra[END_REF] for the Nédélec interpolant on the tetrahedron T 3 = {(x 1 , x 2 , x 3 ) ∈ R 3 ; x i > 0, i = 1, 2, 3 and 0 < x 1 + x 2 + x 3 < 1} is not valid. In other words, the anisotropic estimate

( Π Ned v) l 0, T 3 v l H 1 ( T 3 ) + curl v 0, T 3 , l = 1, 2, 3 does not hold for all v ∈ H 1 ( T 3 ) = H 1 ( T 3 ) 3 .

Anisotropic elements

In this subsection we will denote the functions on the elements K and K by u and u, respectively. Analogous notation will be used for coordinates (e.g., x ∈ K and x ∈ K) and for differential operators (e.g., div and div).

First, let us prove the following auxiliary result.

Lemma 5.1 Let K be the image of the reference element K ( K = T or K = Q) under diagonal scaling with matrix B = h 1 0 0 h 2 , where h l > 0. Then for any u ∈ H -1/2 (K) there holds

u H -1/2 ( K) max {h 1/2 1 , h 1/2 2 } h 1 h 2 u H -1/2 (K) , (5.24)
where u( x) = u(B x), x = ( x 1 , x 2 ) ∈ K.

Proof. By the definition of the dual norm

u H -1/2 ( K) = sup v∈ H1/2 ( K) ( u, v) 0, K v H1/2 ( K) .
(5.25)

We now estimate the norm v H1/2 ( K) . If v ∈ H 1 ( K), then diagonal scaling yields

v 2 0, K = (h 1 h 2 ) -1 v 2 0,K , ∂ 1 v 2 0, K ≃ h 1 h -1 2 ∂ 1 v 2 0,K , ∂ 2 v 2 0, K ≃ h -1 1 h 2 ∂ 2 v 2 0,K . Therefore, v 2
H 1 0 ( K) = ∂ 1 v 2 0, K + ∂ 2 v 2 0, K min {h 2 1 , h 2 2 } h 1 h 2 v 2 H 1 0 (K)
, and by interpolation between L 2 and H 1 0 we find that

v 2 H1/2 ( K) min {h 1 , h 2 } h 1 h 2 v 2 H1/2 (K) = 1 max {h 1 , h 2 } v 2 H1/2 (K) ∀ v ∈ H1/2 ( K). (5.26)
Since ( u, v) 0, K = (h 1 h 2 ) -1 (u, v) 0,K , we use (5.26) in (5.25) to obtain inequality (5.24). ✷ Now, we are in a position to prove the stability result and the corresponding error estimate for the Raviart-Thomas interpolation on anisotropic elements. Theorem 5.3 Let K be either the triangle T with vertices (0, 0), (h 1 , 0), (h 1 , h 2 ), or the rectangle Q with vertices (0, 0), (h 1 , 0), (0, h 2 ), (h 1 , h 2 ), where h l > 0. Denote h max := max {h 1 , h 2 }. Then for any u ∈ H 1/2 (K) with div u ∈ R there holds for l = 1, 2

(Π RT u) l 2 0,K u l 2 0,K + h 3 max h 1 h 2 |u l | 2 H 1/2 (K) + |u l+1 | 2 AH 1/2 l+1 (K) + div u 2 H -1/2 (K) (mod. 2)
(5.27) and

u l -(Π RT u) l 2 0,K h 3 max h 1 h 2 |u l | 2 H 1/2 (K) + |u l+1 | 2 AH 1/2 l+1 (K) + div u 2 H -1/2 (K) (mod. 2).
(5.28)

Proof. We consider only the case of the triangle, K = T . The proof is similar in the case K = Q.

We use the Piola transformation to define u ∈ H 1/2 ( T ) ∩ H(div, T ) on the reference triangle T as follows:

u( x) = h 1 h 2 B -1 u(B x) with B = h 1 0 0 h 2 .
Then we have

u 1 2 0, T = h 2 2 (h 1 h 2 ) -1 u 1 2 0,T = h 2 h -1 1 u 1 2 0,T , | u 1 | 2 H 1/2 ( T ) = T T | u 1 ( x) -u 1 ( y)| 2 (| x 1 -y 1 | 2 + | x 2 -y 2 | 2 ) 3/2 d x d y = h 2 2 (h 1 h 2 ) -2 T T |u 1 (x) -u 1 (y)| 2 (h -2 1 |x 1 -y 1 | 2 + h -2 2 |x 2 -y 2 | 2 ) 3/2 dx dy ≤ h -2 1 h 3 max |u 1 | 2 H 1/2 (T ) , and
| u 2 | 2 AH 1/2 2 ( T ) = 1 0 | u 2 ( x 1 , •)| 2 H 1/2 (0, x 1 ) d x 1 = h -1 1 h 1 0 h 2 x 1 /h 1 0 h 2 x 1 /h 1 0 h 2 1 |u 2 (x 1 , x 2 ) -u 2 (x 1 , y 2 )| 2 h -2 2 |x 2 -y 2 | 2 dx 2 dy 2 h 2 2 dx 1 = h 1 |u 2 | 2 AH 1/2
2 (T ) . Furthermore, the standard properties of the Piola transformation yield

Π RT u(x) = 1 h 1 h 2 B Π RT u( x) = 1/h 2 0 0 1/h 1 Π RT u( x) and div u = h 1 h 2 div u ∈ R.
Therefore, applying Theorem 5.2 (ii) and Lemma 5.1 and using the fact that div u ∈ R (hence, div u 0, T ≃ div u H -1/2 ( T ) ), we obtain

(Π RT u) 1 2 0,T = h -2 2 h 1 h 2 ( Π RT u) 1 2 0, T h 1 h -1 2 u 1 2 0, T + | u 1 | 2 H 1/2 ( T ) + | u 2 | 2 AH 1/2 2 ( T ) + div u 2 0, T u 1 2 0,T + h 3 max h 1 h 2 |u 1 | 2 H 1/2 (T ) + h 2 1 h 2 |u 2 | 2 AH 1/2 2 (T ) + h 1 h 2 h max div u 2 H -1/2 (T ) .
Recalling that h max = max {h 1 , h 2 }, it is easy to see that

h 2 1 h 2 ≤ h 1 h 2 h max ≤ h 3 max h 1 h 2 ,
and then inequality (5.27) follows (for l = 1 and K = T ). Arguing as above and using Corollary 5.3 (ii) instead of Theorem 5.2 (ii) we establish the error estimate in (5.28) for l = 1. The proof is analogous in the case l = 2. ✷

We can now estimate the L 2 -error of the Raviart-Thomas interpolation on the graded mesh ∆ β h on Γ. The specific estimate that we need is for H 1/2 -regular vector fields with discrete divergence.

Lemma 5.2 For any u ∈ H 1/2 ⊥ (Γ) such that div Γ u ∈ div Γ X h there holds u -Π RT u 0,Γ h 1-β/2 u H 1/2 ⊥ (Γ) . (5.29)
Proof. Let F be a face of Γ, and let T F ⊂ F be a triangular block of elements, see Figure 1 (the arguments are analogous for the parallelogram block of elements Q F ). The triangle T F is mapped onto the unit triangle T by the affine transformation which is independent of the mesh parameter h. Let us first establish the error estimate for the Raviart-Thomas interpolation on the unit triangle T partitioned into elements as shown in Figure 1.

The graded mesh on T comprises the quadrilaterals K ij = I i ×I j (i, j = 1, . . . , N, i > j) isomorphic to (0, h i ) × (0, h j ) with h i ≥ h j and the triangles K ii isomorphic to the triangle with vertices (0, 0), (0, h i ), (h i , h i ). Applying error estimates from Theorem 5.3 on each element K ij (i ≥ j), we have for l = 1, 2:

u l -(Π RT u) l 2 0,K ij h 2 i h -1 j |u l | 2 H 1/2 (K ij ) + |u l+1 | 2 AH 1/2 l+1 (K ij ) + div u 2 H -1/2 (K ij ) (mod. 2).
Summing these estimates over all elements in T and recalling that

h 2 i h -1 j h 2-β for 1 ≤ i, j ≤ N, we obtain u -Π RT u 2 0, T h 2-β N i,j=1 i≥j |u| 2 H 1/2 (K ij ) + 2 l=1 |u l | 2 AH 1/2 l (K ij ) + div u 2 H -1/2 (K ij ) . (5.30)
Note that by standard superposition argument

N i=1 i j=1 |u 2 | 2 AH 1/2 2 (K ij ) N i=1 I i |u 2 (x 1 , •)| 2 H 1/2 (0,x 1 ) dx 1 = 1 0 |u 2 (x 1 , •)| 2 H 1/2 (0,x 1 ) dx 1 = |u 2 | 2 AH 1/2 2 ( T ) (5.2) u 2 2 H 1/2 ( T ) , (5.31)
and similarly for u 1 . Hence, using standard superadditivity properties of the H 1/2 -seminorm and the H -1/2 -norm, we deduce from (5.30) the following error estimate on T :

u -Π RT u 0, T h 1-β/2 u H 1/2 ( T ) + div u H -1/2 ( T ) .
Applying now the Piola transformation associated with the mapping T F → T , patching together all individual blocks of elements on all faces of Γ, and using the superadditivity of H 1/2 -and H -1/2 -norms (as the functions of subdomains), we obtain

u -Π RT u 0,Γ h 1-β/2 u H 1/2 -(Γ) + div Γ u H -1/2 (Γ) (5.32)
(here and below we use the space H s -(Γ), s > 0, which is defined in a piecewise fashion by localisation to each face of Γ, with the norm u 2

H s -(Γ) := F ⊂Γ u| F 2 H s (F ) ).
Inequality (5.29) follows from (5.32) due to the continuity property of div Γ : [START_REF]On traces for functional spaces related to Maxwell's equations, Part II: Hodge decompositions on the boundary of Lipschitz polyhedra and applications[END_REF]Section 4.2]). ✷

H 1/2 ⊥ (Γ) → H -1/2 (Γ) (see

Discrete decomposition and the gap property

Following the ideas from [START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF] and [18, Section 9.1], we can use the Raviart-Thomas interpolation operator Π RT to define the discrete counterparts of V and W in (4.3) (e.g., we can set

V h := Π RT (R(X h ))
, where R is the regularised projector introduced in Section 4). However, as it follows from the results in Section 5, the Raviart-Thomas interpolation of low-regular vector fields on graded meshes ∆ β h is only stable (with respect to the L 2 -norm) when β < 2. Since the definition of the energy space X for the EFIE involves the dual space H -1/2 (Γ) with a weaker norm than • 0,Γ , we can relax the restriction on the grading parameter β by employing a different projection onto the boundary element space and using a duality argument on individual faces of Γ. This approach was successfully used by Hiptmair and Schwab in [START_REF] Hiptmair | Natural boundary element methods for the electric field integral equation on polyhedra[END_REF]Section 8] and by Buffa and Christiansen in [START_REF] Buffa | The electric field integral equation on Lipschitz screens: definitions and numerical approximation[END_REF]Section 4.2.2] in the context of the h-BEM with shape-regular meshes for the EFIE (see [START_REF]Natural p-BEM for the electric field integral equation on screens[END_REF][START_REF]The hp-BEM with quasi-uniform meshes for the electric field integral equation on polyhedral surfaces: a priori error analysis[END_REF] for applications of these ideas to the analysis of the p-BEM and the hp-BEM with quasi-uniform meshes). We will demonstrate below that using these techniques together with stability properties and error estimates for the Raviart-Thomas interpolation on anisotropic elements, one can design a stable discrete decomposition of the boundary element space on ∆ β h and prove the corresponding gap property (3.2) for any β < 3.

The construction of the desired projection operator is technically involved. Therefore, we formulate here the final result relevant to our discussion and give a detailed proof in the next section. In the Proposition below, Π 0 denotes the L 2 (Γ)-projection onto the space of piecewise constant functions over the mesh ∆ β h , and

H -1/2 - (Γ) denotes the dual space of H 1/2 -(Γ) (with L 2
t (Γ) as pivot space).

Proposition 6.1 There exists an operator

Q h : H s -(Γ) ∩ H(div Γ , Γ) → X h (s > 0) such that div Γ • Q h = Π 0 • div Γ , (6.1)
and for any ε > 0

u -Q h u H -1/2 - (Γ) h 1/2-ε u -Π RT u H(div Γ ,Γ) ∀ u ∈ H s -(Γ) ∩ H(div Γ , Γ). (6.2)
Thus, the operator Q h inherits the crucial commuting diagram property (6.1) of the classical RT-interpolation operator and, at the same time, allows to gain an extra power of h when estimating the error (u -Q h u) in the dual norm. Corollary 6.1 For any u ∈ H

1/2 ⊥ (Γ) such that div Γ u ∈ div Γ X h , one has Q h u ∈ X h , div Γ Q h u = div Γ u,
and for any ε > 0 there holds

u -Q h u X = u -Q h u H -1/2 (Γ) h 3/2-β/2-ε u H 1/2 ⊥ (Γ) . (6.3) Since • H -1/2 (Γ) • H -1/2 - (Γ)
, this result immediately follows from Proposition 6.1 and Lemma 5.2 due to the commuting diagram property for Π RT .

Since R X h ⊂ H 1/2 ⊥ (Γ) (see (4.1)) and div Γ R X h = div Γ X h (see (4.2)), the following definitions are valid thanks to Proposition 6.1:

V h := (Q h • R)X h , W h := (Id -Q h • R)X h .
Using the commuting diagram property (6.1), we have

div Γ Q h R u h = Π 0 (div Γ R u h ) (4.2) = Π 0 (div Γ u h ) = div Γ u h ∀ u h ∈ X h . (6.4) Therefore, R Q h R = R on X h , (6.5)
and hence Q h • R : X h → X h is a projection. This fact confirms that X h = V h ⊕ W h . Property (6.4) also implies W h ⊂ W, and Corollary 6.1 yields stability of the discrete decomposition in the following sense: there exists C = C(Γ, β) such that

Q h R u h X ≤ C R u h H 1/2 ⊥ (Γ) (4.1) ≤ C div Γ u h H -1/2 (Γ) ≤ C u h X ∀ u h ∈ X h ,
provided that β < 3. This verifies property (B) from Section 3. It remains to establish the gap property (C). Inequality (3.2) in (C) is another consequence of Corollary 6.1:

sup v h ∈V h inf v∈V v -v h X v h X ≤ sup v h ∈V h R v h -v h X v h X (6.5) = sup v h ∈V h R v h -Q h R v h X v h X (6.3) h 3/2-β/2-ε sup v h ∈V h R v h H 1/2 ⊥ (Γ) v h X (4.1) h 3/2-β/2-ε .
This completes the proof of Theorem 3.1.

7 Proof of Proposition 6.1

In this section, we give a constructive proof of Proposition 6.1. For any u ∈ H s -(Γ) ∩ H(div Γ , Γ) we construct Q h u in the Raviart-Thomas spaces on individual faces of Γ. Let F be a single face of Γ, and let ∆ β h (F ) denote the restriction of the graded mesh ∆ β h onto F . For the sake of simplicity of notation we will omit the subscript F for differential operators over this face, e.g., we will write div for div F . We will also write (•, •) for the L 2 (F )-and L 2 (F )-inner products, and similarly • for the corresponding norms of scalar functions and vector fields. First, let us prove the following auxiliary result. Lemma 7.1 For any s > 1/2, the Raviart-Thomas interpolation operator Π RT :

H s (F ) ∩ H(div, F ) → RT 0 (∆ β h (F )), is L 2 (F )-stable, i.e.
, there exists a constant C > 0 independent of h such that

Π RT u 0,F ≤ C u H s (F ) + div u 0,F ∀u ∈ H s (F ) ∩ H(div, F ). (7.1)
Proof. Similarly to the proof of Lemma 5.2, it is sufficient to establish (7.1) for the unit triangle T partitioned into elements as shown in Figure 1 (this is because the affine transformations that map triangular blocks of elements T F ⊂ F onto T are independent of h).

The graded mesh on T (see Figure 1) comprises anisotropic quadrilaterals K ij = I i × I j (i, j = 1, . . . , N, i > j) isomorphic to (0, h i ) × (0, h j ) with h i ≥ h j and shape-regular triangles K ii isomorphic to the triangle with vertices (0, 0), (0, h i ), (h i , h i ). Using the Piola transform associated with the mapping K → K ij (i ≥ j; K = T or Q), we define u ∈ H s ( K) ∩ H(div, K). Then, by the standard properties of the Piola transform we have

u 1 2 0, K ≃ h -1 i h j u 1 2 0,K ij , div u 2 0, K ≃ h i h j div u 2 0,K ij , (Π RT u) 1 2 0,K ij ≃ h i h -1 j ( Π RT u) 1 2
0, K . The application of the scaling argument yields:

| u 1 | 2 H s ( T ) = T T | u 1 ( x) -u 1 ( y)| 2 (| x 1 -y 1 | 2 + | x 2 -y 2 | 2 ) 1+s d x d y ≃ h 2 i (h i h i ) -2 K ii K ii |u 1 (x) -u 1 (y)| 2 (h -2 i |x 1 -y 1 | 2 + h -2 i |x 2 -y 2 | 2 ) 1+s dx dy = h 2s i |u 1 | 2 H s (K ii ) , | u 1 | 2 AH s 1 ( Q) = 1 0 | u 1 (•, x 2 )| 2 H s (0,1) d x 2 = h -1 j h j 0 h i 0 h i 0 h 2 j |u 1 (x 1 , x 2 ) -u 1 (y 1 , x 2 )| 2 h -1-2s i |x 1 -y 1 | 1+2s dx 1 dy 1 h 2 i dx 2 ≃ h 2s-1 i h j |u 1 | 2 AH s 1 (K ij ) (i > j),
and analogously,

| u 1 | 2 AH s 2 ( Q) ≃ h -1 i h 1+2s j |u 1 | 2 AH s 2 (K ij ) (i > j).
Therefore, applying Theorem 5.2 (i), we obtain

(Π RT u) 1 2 0,K ii ≃ ( Π RT u) 1 2 0, T u 1 2 H s ( T ) + div u 0, T ≃ u 1 2 0,K ii + h 2s i |u 1 | 2 H s (K ii ) + h 2 i |div u| 2 0,K ii . (7.2)
Similarly, applying Theorem 5.1 (i) and recalling (5.1), we have for i > j

(Π RT u) 1 2 0,K ij ≃ h i h -1 j ( Π RT u) 1 2 0, Q h i h -1 j u 1 2 H s ( Q) ≃ h i h -1 j u 1 2 0, Q + | u 1 | 2 AH s 1 ( Q) + | u 1 | 2 AH s 2 ( Q) ≃ u 1 2 0,K ij + h 2s i |u 1 | 2 AH s 1 (K ij ) + h 2s j |u 1 | 2 AH s 2 (K ij ) . (7.3)
The estimates analogous to (7.2) and ( 7.3) are also valid for (Π RT u) 2 0,K ij with i ≥ j.

Combining the estimates for both components of Π RT u over all elements in T and then using the superposition argument as in (5.31) for anisotropic seminorms and the superadditivity property of the H 1/2 -norm, we arrive at the desired result. ✷

Our construction of the operator Q h follows the technique used by Hiptmair and Schwab in the proof of Lemma 8.1 in [START_REF] Hiptmair | Natural boundary element methods for the electric field integral equation on polyhedra[END_REF] but relies on stability properties of the Raviart-Thomas interpolation on graded meshes over individual faces of Γ. Given u ∈ H s (F ) ∩ H(div, F ), s > 0, we consider the following mixed problem:

Find (z, f ) ∈ H(div, F ) × L 2 * (F ) such that (z, v) + (div v, f ) = (u, v) ∀v ∈ H 0 (div, F ), (div z, g) = (div u, g) ∀g ∈ L 2 * (F ), z • ñ = u • ñ on ∂F . (7.4)
Here, L 2 * (F ) := v ∈ L 2 (F ); (v, 1) = 0 , ñ is the unit outward normal vector to ∂F , and

H 0 (div, F ) := {v ∈ H(div, F ); v • ñ| ∂F = 0}.
The unique solvability of (7.4) is proved by standard techniques (see [11, Chapter II]). In fact, it is clear that the pair (u, 0) solves (7.4).

A conforming Galerkin approximation of problem (7.4) with Raviart-Thomas elements on the graded mesh ∆ β h (F ) reads as:

Find (z h , f h ) ∈ X h (F ) × R h (F ) such that (z h , v) + (div v, f h ) = (u, v) ∀v ∈ X h (F ) ∩ H 0 (div, F ), (div z h , g) = (div u, g) ∀g ∈ R h (F ), z h • ñ = Π RT u • ñ on ∂F . (7.5)
Here, X h (F ) denotes the restriction of X h onto the face F , and R h (F ) := {g ∈ L 2 (F ); g| K = const, ∀ K ∈ ∆ β h (F ) and (g, 1) = 0}.

Note that the third equation in (7.5) implies (div(uz h ), 1) = 0. Hence, the second identity in (7.5) holds for any piecewise constant function g ∈ div X h (F ). Thus, div z h is the L 2 (F )-projection of div u onto div X h (F ). In particular, if div u ∈ div X h (F ) then div z h = div u.

We now prove the unique solvability of (7.5). First, for any g h ∈ R h (F ) we find a function φ ∈ H 1 * (F ) := {φ ∈ H 1 (F ); (φ, 1) = 0} solving the variational problem

(∇φ, ∇ φ) = (g h , φ) ∀ φ ∈ H 1 * (F ). (7.6)
Applying the standard regularity result for problem (7.6) (see, e.g., [22, p. 82]), we conclude that φ ∈ H 1+r (F ) with some r ∈ ( 1 2 , π ω ) (here, ω < 2π denotes the maximal internal angle at the vertices of F ), and

∇ φ H r (F ) φ H 1+r (F ) g h . (7.7)
Therefore, ∇ φ ∈ H r (F ) ∩ H 0 (div, F ), r > 1 2 , and the interpolant Π RT ∇ φ ∈ X h (F ) ∩ H 0 (div, F ) is well defined and stable, due to Lemma 7.1. Moreover, div (Π RT ∇ φ) = Π 0 (div ∇ φ) (7.6) = g h . Hence, using (7.1) and (7.7) we prove the discrete inf-sup condition:

sup v h ∈X h (F )∩H 0 (div,F ) v h =0 (div v h , g h ) v h H(div,F ) ≥ (div (Π RT ∇ φ), g h ) Π RT ∇ φ H(div,F ) ≥ g h 2 C ∇ φ H r (F ) + div ∇ φ + div (Π RT ∇ φ) ≥ C g h ∀ g h ∈ R h (F ).
This condition along with the property div X h (F ) ∩ H 0 (div, F ) = R h (F ) ensures existence, uniqueness, and quasi-optimality of the solution (z h , f h ) to (7.5) (see [START_REF] Brezzi | Mixed and Hybrid Finite Element Methods[END_REF]). In particular, using the quasi-optimality and recalling that z = u, f = 0, we estimate

u -z h H(div,F ) inf v h ∈X h (F ) (v h -Π RT u)•ñ| F =0 u -v h H(div,F ) + inf g h ∈R h (F ) f -g h u -Π RT u H(div,F ) . (7.8)
We now estimate uz h H-1/2 (F ) . One has for any ε ∈ (0, 1 2 )

uz h H-1/2 (F ) ≤ uz h H-1/2+ε (F ) = sup Similarly to (7.7), the regularity result for ϕ reads as

ϕ ∈ H 3/2-ε (F ), ϕ H 3/2-ε (F ) f (H 1/2+ε (F )) ′ w H 1/2-ε (F ) , (7.11)
where f ∈ (H 1/2+ε (F )) ′ is defined by f (φ) = -(w, ∇ φ), ∀ φ ∈ H 1/2+ε (F ).

Then we set q := w + ∇ ϕ ∈ H 1/2-ε (F ) ∩ H 0 (div, F ). (7.12)

It also follows from (7.10) that div q = div w + div ∇ ϕ = 0. Furthermore, we have by (7.11)-(7.12) that

q H 1/2-ε (F ) w H 1/2-ε (F ) + ϕ H 3/2-ε (F ) w H 1/2-ε (F ) . (7.13)
We now use (7.12) and integration by parts to represent the numerator in (7.9) as (uz h , w) = (uz h , q) -(uz h , ∇ ϕ) = (uz h , q) + (div (uz h ), ϕ) -((uz h ) • ñ, ϕ) 0,∂F .

Hence, using (7.4), (7.5) and recalling that z = u, f = 0, we find for any q h ∈ X h (F ) ∩ H 0 (div, F ) and arbitrary ϕ h ∈ R h (F )

|(u -z h , w)| = |(u -z h , q -q h ) + (u -z h , q h ) + (div (u -z h ), ϕ -ϕ h ) -((u -Π RT u) • ñ, ϕ) 0,∂F |
= |(uz h , qq h ) + (div q h , f h ) + (div (uz h ), ϕ -ϕ h ) -((u -Π RT u) • ñ, ϕ) 0,∂F | ≤ uz h qq h + |(div q h , f h )| + div (uz h ) ϕ -ϕ h

+ (u -Π RT u) • ñ H -1+ε (∂F ) ϕ H 1-ε (∂F ) . (7.14)
Let Π q/u RT denote the Raviart-Thomas interpolation operator on the 'coarse' quasi-uniform and shape-regular mesh ∆ q/u h (F ) obtained from the graded mesh ∆ β h (F ) by patching together long and thin elements (see Figure 3). We also denote by Π q/u 0 the L 2 (F )-projector onto the space of piecewise constant functions on ∆ q/u h (F ). Then we set q h := Π q/u RT q ∈ X h (F ) ∩ H 0 (div, F ) and ϕ h := Π q/u 0 ϕ ∈ R h (F ).

By the standard properties of the Raviart-Thomas interpolation and the L 2 -projection on quasi-uniform and shape-regular meshes, we have div q h = Π q/u 0 div q = 0, (7.15)

qq h h 1/2-ε q H 1/2-ε (F ) (7.13)

h 1/2-ε w H 1/2-ε (F ) , (7.16)
ϕ -ϕ h h ϕ H 3/2-ε (F ) (7.11) h w H 1/2-ε (F ) . (7.17) Figure 3: An example of shape-regular quasi-uniform mesh (thicker lines) obtained by patching together elements of the graded mesh (thinner lines).

To estimate (u -Π RT u) • ñ H -1+ε (∂F ) we recall that e h (u -Π RT u) • ñ = 0 for any element edge e h ⊂ ∂F . Therefore, we can use a standard duality argument to prove (cf. [13, p. 259])

(u -Π RT u) • ñ H -1+ε (∂F ) max e h ⊂∂F |e h | 1-ε (u -Π RT u) • ñ L 2 (∂F ) .
Then by interpolation we obtain (u -Π RT u) • ñ H -1+ε (∂F ) h 1/2-ε (u -Π RT u) • ñ H -1/2 (∂F ) h 1/2-ε u -Π RT u H(div,F ) (7.18) (here, we also used the continuity of the normal trace operator v → v • ñ| ∂F as a mapping H(div, F ) → H -1/2 (∂F )).

Furthermore, one has ϕ H Using this estimate in (7.9) we obtain uz h H-1/2 (F ) h 1/2-ε u -Π RT u H(div,F ) . (7.20)

Lemma 4 . 1 [ 2 ,

 412 Section 3] For any bounded Lipschitz domain Ω ⊂ R 3 there exists a continuous mapping L

  and that the associated degrees of freedom are given by e i u • n ds, i = 1, 2, 3, 4.

Figure 2 :

 2 Figure 2: The reference square Q and the reference triangle T .

w∈H 1 /

 1 2-ε (F )\{0} |(uz h , w)| w H 1/2-ε (F ) .(7.9)For a given w ∈ H 1/2-ε (F ), we solve the following problem:Find ϕ ∈ H 1 * (F ) such that (∇ ϕ, ∇ φ) = -(w, ∇ φ) ∀ φ ∈ H 1 * (F ). (7.10)

  z h , w)| h 1/2-ε u -Π RT u H(div,F ) w H 1/2-ε (F ) .

		(7.11)
	1-ε (∂F )	ϕ H 3/2-ε (F )	w H 1/2-ε (F ) .	(7.19)
	Now, using (7.15)-(7.19) in (7.14) and recalling (7.8) we find
	|(u

Galerkin BEM on graded meshes. The main result.For approximate solution of (2.1) we apply the natural BEM based on Galerkin discretisations with lowest-order Raviart-Thomas spaces on graded meshes.

Acknowledgement. A significant part of this work has been done while A.B. was visiting Laboratoire de Mathématiques et ses Applications de Valenciennes, Université de Valenciennes et du Hainaut-Cambrésis (Valenciennes, France). This author is grateful to the colleagues in that department for their hospitality and stimulating research atmosphere.

Now we can prove the desired result.

Proof of Proposition 6.1. For any u ∈ H s -(Γ) ∩ H(div Γ , Γ), we define Q h u ∈ X h face by face as Q h u| F := z h for any face F ⊂ Γ, where z h is a unique (vectorial) solution to (7.5). Then the commuting diagram property (6.1) follows from the second identity in (7.5), and inequality (7.20) yields estimate (6.2). ✷