# A finite element method on composite grids based on Nitsche's method

Anita Hansbo, Peter Hansbo, Mats G. Larson

► To cite this version:

Anita Hansbo, Peter Hansbo, Mats G. Larson. A finite element method on composite grids based on Nitsche's method. ESAIM: Mathematical Modelling and Numerical Analysis, 2003, 37 (3), pp.495-514. 10.1051/m2an:2003039. hal-01353508

HAL Id: hal-01353508

https://hal.science/hal-01353508

Submitted on 11 Aug 2016

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à 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.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A FINITE ELEMENT METHOD ON COMPOSITE GRIDS BASED ON NITSCH'S METHOD

ANITA HANSBO<fnref n="1" />, PETER HANSBO<fnref n="2" /> AND MATS G. LARSON<fnref n="3" />

**Abstract.** In this paper we propose a finite element method for the approximation of second order elliptic problems on composite grids. The method is based on continuous piecewise polynomial approximation on each grid and weak enforcement of the proper continuity at an artificial interface defined by edges (or faces) of one the grids. We prove optimal order *a priori* and energy type *a posteriori* error estimates in 2 and 3 space dimensions, and present some numerical examples.

*Keywords and phrases.* Nitsche's method, overlapping grids.

## 1. INTRODUCTION

Composite overlapping grids are commonly used with finite difference and difference-related finite volume methods, e.g., as a tool for local mesh refinement [7]. In finite element and related finite volume methods, which are inherently unstructured and thus allow for local mesh refinement, the need for overlapping mesh methods is less obvious. However, a general finite element/volume methodology for handling overlapping meshes would be a useful tool to deal with the often complicated mesh generation problem. Examples of specific applications include: (a) construction of a global mesh for a complex geometry by using overlapping meshes of elementary parts; (b) coupling of unstructured and structured meshes; and (c) coupling of boundary fitted meshes to structured or unstructured meshes, see Figure 1.

The purpose of this paper is to introduce and analyze such a method for a model second order elliptic problem in 2 and 3 space dimensions. Unlike composite grid methods where interpolation is performed on the boundary of the overlap, cf. [3, 7], our approach is based on weak enforcement of the proper continuity across an artificial interface defined by edges (faces) of one of the meshes. The weak enforcement proposed here is constructed in such a way that the resulting scheme is stable and 'arbitrary order consistent' in the sense the exact solution satisfies the discrete equation. Hence we are able to prove optimal *a priori* error estimates for arbitrary order of polynomial approximation under weak mesh conditions; in particular the meshes may overlap in quite an arbitrary fashion.

<footnote n="1">¹ Department of Informatics and Mathematics, University of Trollhättan Uddevalla, Box 957, 461 39 Trollhättan, Sweden.
e-mail: hansbo@solid.chalmers.se</footnote>

<footnote n="2">² Department of Applied Mechanics, Chalmers University of Technology, 412 96 Göteborg, Sweden.</footnote>

<footnote n="3">³ Department of Mathematics, Chalmers University of Technology, 412 96 Göteborg, Sweden.</footnote>
<figure><img src="image_2.png" /><figcaption>FIGURE 1. Overlapping mesh FEM may be used to couple a boundary fitted mesh to an unstructured mesh (top), or to construct a global mesh by using overlapping meshes of elementary parts (bottom).</figcaption></figure>

In composite grid methods, approximations on the overlap are usually computed on all meshes. This may be done by projections [1,6] or interpolation [7] on the interfaces. Alternatively, one may use integration of products of test functions living on both meshes, as in the finite element method proposed by Brezzi, Lions, and Pironneau [5]. In contrast, in our method one only computes an approximation on *one* of the meshes on the overlap; in fact, we do not require the meshes to overlap at all even though this is the situation we have in mind.

Our method stems from the work of Nitsche [12], where a method for consistent weak enforcement of Dirichlet boundary conditions was introduced, and is related, in particular, to Becker, Hansbo and Stenberg [2], where the meshes were assumed to be non-matching on the interface but shape regular on both sides of this interface. Here this last condition is relaxed, which, *e.g.*, makes possible to use highly structured meshes in irregularly shaped subdomains. Related work includes also the mixed penalty approach analyzed in Lazarov *et al.* [10,11].

An outline of the paper is as follows: in Section 2 we formulate the second order elliptic model problem and in Section 3 we state the mesh assumptions and define the numerical method used for the approximation. In Sections 4.1–4.2 we demonstrate the stability of the method and derive the approximation properties of its (non-standard) finite element spaces. Optimal order *a priori* error estimates in a discrete energy norm, as well
as in $L_2$-norm, are shown in Section 4.3. Our a priori analysis is in parts akin to Hansbo and Hansbo [9], where optimal order convergence was shown for a method where a material discontinuity interface was allowed to cut through the elements in an arbitrary fashion. Using similar lines of arguments as in [9], a posteriori error estimates for the control of linear functionals of the error may be derived for the present method. In this work we instead focus our attention in Section 5 on two variants of residual based a posteriori energy norm error estimates, where the element indicators of the second one are designed for ease of implementation, reducing the complications due to the geometry of the mesh. These estimates do not presuppose the saturation assumption used in [2]. Finally, in Section 6, we discuss some implementation details and present numerical examples, including a convergence study using quadratic elements in 2D as well as examples using the a posteriori error estimates as a basis for the implementation of an adaptive algorithm.

## 2. PROBLEM FORMULATION AND PRELIMINARIES

As a model problem, we consider Poisson's equation in a bounded domain $\Omega$ in $\mathbb{R}^n$, $n = 2$ or $n = 3$, with, for simplicity, a convex polygonal boundary $\partial\Omega$. Find $u: \Omega \to \mathbb{R}$ such that

$$ -\Delta u = f \quad \text{in } \Omega, \qquad (1) $$

$$ u = 0 \quad \text{on } \partial\Omega, \qquad (2) $$

with $f \in L^2(\Omega)$.

Consider two given triangulations $T_i^h$, $i = 1, 2$, where $h$ is a mesh size parameter. Assume that they together cover $\Omega$, so that $\bar{\Omega} = \bar{\Omega}_1^* \cup \bar{\Omega}_2^*$ where $\bar{\Omega}_i^* = \cup_{K \in T_i^h} K$. The meshes may overlap in an arbitrary fashion; further assumptions are given below. We then choose an (artificial) internal interface $\Gamma$ composed of edges from the triangles in $T_i^h$ and dividing $\Omega$ into two open disjoint sets $\Omega_i$, $i = 1, 2$, such that $\Omega_i \subset \Omega_i^*$ and $\Omega = \Omega_1 \cup \Omega_2 \cup \Gamma$. We assume that the interface $\Gamma$ does not depend on $h$, by, i.e., assuming that the mesh family $T_i^h$ is obtained by refinement from a single coarse mesh, or by remeshing of a region $\Omega_1$ defined from a selected mesh.

For any sufficiently regular function $u$ in $\Omega_1 \cup \Omega_2$ we define the jump of $u$ on $\Gamma$ by $[u] := u_1|_{\Gamma} - u_2|_{\Gamma}$, where $u_i = u|_{\Omega_i}$ is the restriction of $u$ to $\Omega_i$. Conversely, for $u_i$ defined in $\Omega_i$ we identify the pair $\{u_1, u_2\}$ with the function $u$ which equals $u_i$ on $\Omega_i$.

Our model problem may now, due to the (artificial) interface, be written as:

$$ -\Delta u = f \quad \text{in } \Omega_1 \cup \Omega_2, \qquad (3) $$

$$ u = 0 \quad \text{on } \partial\Omega, \qquad (4) $$

$$ [u] = 0 \quad \text{on } \partial\Omega, \qquad (5) $$

$$ [\nabla_n u] = 0 \quad \text{on } \partial\Omega. \qquad (6) $$

Here $n$ denotes the exterior unit normal to $\Omega_1$ and $\nabla_n v = n \cdot \nabla v$ is the normal flux.

For a bounded open connected domain $D$ we shall use standard Sobolev spaces $H^r(D)$ with norm $||\cdot||_{r,D}$ and spaces $H_0^r(D)$ with zero trace on $\partial D$. The inner products in $H^0(D) = L_2(D)$ is denoted $(\cdot, \cdot)_D$. For a bounded open set $G = \cup_{i=1}^2 D_i$, where $D_i$ are open disjoint components of $G$, we let $H^k(D_1 \cup D_2)$ denote the Sobolev space of functions in $G$ such that $u|_{D_i} \in H^k(D_i)$ with norm

$$ \| \cdot \|_{k, D_1 \cup D_2} = \left( \sum_{i=1}^{2} \| \cdot \|_{k, D_i}^2 \right)^{1/2} . $$
## 3. FINITE ELEMENT SPACES AND METHOD

### 3.1. The finite element spaces

We will use the following notation for mesh related quantities. Let $h_K$ be the diameter of an element $K \in T_i^h$ and $h = \max_{K \in T_i^h, i=1,2} h_K$. To distinguish elements from the two meshes, we will sometimes use indexed element notation $K_i \in T_i^h$ for clarity.

The nodes on $\Gamma$ of the elements in $T_1^h$, together with the points of intersection between elements in $T_2^h$ and $\Gamma$, define a partition of $\Gamma$, $\Gamma = \cup_{j \in J_h} \Gamma^j$. Note that each part $\Gamma^j$ belongs to two elements, one from each mesh. We denote these elements by $K_1^j$ and $K_2^j$, respectively. A local meshsize on $\Gamma$ is defined by

$$h(x) = h_{K_1^j}, \quad x \in \Gamma^j. \tag{7}$$

For any element $K \in T_i^h$, let $P_K = K \cap \Omega_i$ denote the part of $K$ in $\Omega_i$.

We make the following assumptions regarding the meshes:

A1) The triangulations are non-degenerate, i.e.,

$$h_K/\rho_K \le C \quad \forall K \in T_i^h, \quad i = 1,2,$$

where $h_K$ is the diameter of $K$ and $\rho_K$ is the diameter of the largest ball contained in $K$.

A2) The meshes have locally compatible meshsize over $\Gamma$. More precisely, let $K_1^j \in T_1^h$ and $K_2^j \in T_2^h$ be the elements which contain a specific part $\Gamma^j$ of $\Gamma$. We assume that

$$ch_{K_1^j} \le h_{K_2^j} \le Ch_{K_1^j} \quad \forall j \in J_h.$$

Here and below, $C$ and $c$ denote generic constants.

We shall seek a discrete solution $U = (U_1, U_2)$ in the space $V^h = V_1^h \times V_2^h$, where

$$V_i^h = \{ \phi \in H^1(\Omega_i) : \phi|_{K \cap \Omega_i} \text{ is a polynomial of degree } p \quad \forall K \in T_i^h, \phi|_{\partial \Omega_i} = 0 \}.$$

Note that functions in $V^h$ are, in general, discontinuous across $\Gamma$. As for the nodal representation of polynomials on the parts in $\Omega_2$, see Figure 2.

### 3.2. The finite element method

The method is defined by the variational problem: find $U \in V^h$ such that

$$a^h(U, \phi) = l(\phi), \quad \forall \phi \in V^h, \tag{8}$$

where

$$a^h(U, \phi) = (\nabla U, \nabla \phi)|_{\Omega_1 \cup \Omega_2} - (\langle \nabla_n U \rangle, [\phi])_\Gamma - ([U], \langle \nabla_n \phi \rangle)_\Gamma + (\lambda h^{-1}[U], [\phi])_\Gamma, \\ l(\phi) = (f, \phi)|_{\Omega_1}$$

with

$$\langle \nabla_n v \rangle = \nabla_n v_1 \quad \text{on } \Gamma,$$

and where $h$ is the local meshsize (7). The continuity conditions of $u$ and $\nabla_n u$ at $\Gamma$ are satisfied weakly by means of a variant of Nitsche's method [12] for consistent weak enforcement of Dirichlet boundary conditions. To ensure stability, the parameter $\lambda$ has to be taken sufficiently large and we return to this issue in Lemma 4.4 below.
<figure><img src="image_3.png" /><figcaption>FIGURE 2. The interface $\Gamma$ consists of element edges from elements in $T_1^h$. Each part $\Gamma^j$ belongs to two triangles, $K_1^j \in T_1^h$ and $K_2^j \in T_2^h$. Nodes for representing a quadratic polynomial on the element $K_2^j$ are indicated. The same nodes are used in the implementation to represent a polynomial on the part $P_{K_2^j} = K_2^j \cap \Omega_2$.</figcaption></figure>

Here we use a one sided approximation of the normal flux on the interface instead of the usual symmetric average

$$ \langle \nabla_n v \rangle = (\nabla_n v_1 + \nabla_n v_2)/2 \quad \text{on } \Gamma, $$

commonly used in the discontinuous Galerkin method and in the context of interfaces between meshes that are shape regular on both sides of the interface by Becker, Hansbo and Stenberg [2]. The latter situation rules out the possibility of arbitrary small elements in $\Omega_2$ close to the interface. As has been noted by Stenberg [14], any convex combination of the fluxes yields consistent methods. In Hansbo and Hansbo [9], this fact was exploited to allow for internal discontinuities along an interface in shape regular elements, choosing convex combinations that take into account the size of the parts of the cut element to ensure stability. In fact, a one-sided flux approximation could have been used there too, however with negligible gain in implementation complexity. Likewise, in the case of the present work a stable two-sided variant may be defined. However, in this case such a method becomes more complicated to construct and implement since the element parts on each side of the interface now stems from two different meshes. This is an important practical point, and the main reason why we use one-sided fluxes here.

With these definitions, we have the following consistency relation.

**Proposition 3.1.** The discrete problem (8) is consistent in the sense that, for $u$ solving (1, 2) there holds

$$ a^h(u, \phi) = l(\phi), \quad \forall \phi \in V^h, $$

or, equivalently,

$$ a^h(u - U, \phi) = 0, \quad \forall \phi \in V^h. \qquad (9) $$

*Proof.* Let $u$ be the solution of the Poisson problem (1, 2). Then $\langle \nabla_n u \rangle = n_1 \nabla \cdot u_1 = -n_2 \nabla \cdot u_2 = \nabla_n u$ and $[u] = 0$. By Green's formula it follows that

$$ \begin{aligned} a^h(u, \phi) &= (\nabla u, \nabla \phi)_{\Omega_1 \cup \Omega_2} - (\langle \nabla_n u \rangle, [\phi])_{\Gamma} \\ &= (-\Delta u, \phi)_{\Omega_1 \cup \Omega_2} + (\nabla_n u - \langle \nabla_n u \rangle, [\phi])_{\Gamma} = (f, \phi)_{\Omega_1 \cup \Omega_2}, \end{aligned} $$

which proves the result.

<figure><img src="image_4.png" /></figure>
# 4. A PRIORI ANALYSIS

## 4.1. Interpolation error estimates

In the error analysis, we shall use the following mesh dependent norms:

$$ \|v\|_{1/2,h,\Gamma}^2 := \|h(x)^{-1/2}v\|_{0,\Gamma}^2 = \sum_{j \in J_h} h_{K_1^j}^{-1} \|v\|_{0,\Gamma^j}^2, $$

$$ \|v\|_{-1/2,h,\Gamma}^2 := \|h(x)^{1/2}v\|_{0,\Gamma}^2 = \sum_{j \in J_h} h_{K_1^j}^{-1} \|v\|_{0,\Gamma^j}^2, $$

and

$$ \|v\|^2 := \|\nabla v\|_{0,\Omega_1 \cup \Omega_2}^2 + \|\langle \nabla n v \rangle\|_{-1/2,h,\Gamma}^2 + \|[v]\|_{1/2,h,\Gamma}^2. \quad (10) $$

Note for future reference that

$$ (u,v)_{\Gamma} \le \|v\|_{1/2,h,\Gamma} \|v\|_{-1/2,h,\Gamma}. \quad (11) $$

To show that functions in $V^h$ approximates functions $v \in H_0^1(\Omega) \cap H^{p+1}(\Omega)$ to the order $h^p$ in the norm $\|\cdot\|$, we define an interpolant $I^h v \in V^h$ of $v$ by $I^h v = I_i^h v$ on $\Omega_i$, $i = 1, 2$. Here $I_i^h$ is the standard Lagrange nodal interpolant on the mesh $T_i^h$ of $\Omega_i^*$. The following local interpolation error estimates holds, see, e.g., [4].

$$ \|v - I_i^h v\|_{m,K} \le Ch_K^{p+1-m} |v|_{p+1,K}, \quad m = 0, 1, p \ge 1, K \in T_h^i, i = 1, 2. \quad (12) $$

One may note that a node of interpolation used to define $I_2^h v$ lies in $\bar{\Omega}_2^* = \cup_{K \in T_2^h} \bar{K}$ but not necessarily in $\bar{\Omega}_2$. The following interpolation error estimate holds.

**Lemma 4.1.** Let $I^h v = I_i^h v$ on $\Omega_i$, $i = 1, 2$, where $I_i^h$ is the Lagrange nodal interpolant on the mesh $T_i^h$ of $\Omega_i^*$. Then, for $p \ge 1$,

$$ \|v - I^h v\| \le Ch^p |v|_{p+1,\Omega}, \quad \forall v \in H_0^1(\Omega) \cap H^{p+1}(\Omega). \quad (13) $$

In the proof of this result, we need to estimate the interpolation error at the interface. We recall the well known trace inequality

$$ \|w\|_{0,\partial \tilde{K}}^2 \le C \|w\|_{0,\tilde{K}} \|w\|_{1,\tilde{K}}, \quad \forall w \in H^1(\tilde{K}). \quad (14) $$

on a reference element $\tilde{K}$. The following Lemma provides a scaled version of (14) which can be used to estimate traces not only on the boundary but on arbitrary lines (planes) intersecting the element.

**Lemma 4.2.** Let $L$ be the intersection between a line (plane) and an element $K$. Then

$$ \|w\|_{0,L}^2 \le Ch_K^{-1} \|w\|_{0,K}^2 + h_K \|w\|_{1,K}^2, \quad \forall w \in H^1(K), \quad (15) $$

where the constant $C$ is independent of $L$.

*Proof.* Map the element by an affine map to a reference element $\tilde{K}$ and denote by $\tilde{L}$ the image of $L$. For the plane case, let $(\xi, \eta)$ denote local coordinates on $\tilde{K}$ such that $\eta = 0$ on $\tilde{L}$. If $\tilde{L}$ divides $\tilde{K}$ into two subsets let $\tilde{K}_1$ denote one of these sets, and else set $\tilde{K}_1 = \tilde{K}$. Let $n$ denote the outward pointing unit normal of $\tilde{K}_1$ and note that we may assume that $n_\eta = 1$ on $\tilde{\Gamma}$. By the divergence theorem,

$$ \begin{aligned} 2 \int_{\tilde{K}_1} w \frac{\partial w}{\partial \zeta} dV &= \int_{\tilde{K}_1} \operatorname{div} (0, w^2) dV = \int_{\partial \tilde{K}_1} n \cdot (0, w^2) dA \\ &= \int_{\tilde{L}} w^2 dA + \int_{\partial \tilde{K}_1 \setminus \tilde{L}} n_\eta w^2 dA. \end{aligned} $$
We thus find, using Cauchy–Schwarz’ inequality and (14), that

$$
\begin{aligned}
\|w\|_{0,\bar{L}}^2 &\le 2 \|w\|_{0,\bar{K}_1} \|w\|_{1,\bar{K}_1} + \|w\|_{0,\partial\bar{K}_1\setminus\bar{L}}^2 \\
&\le 2 \|w\|_{0,\bar{K}_1} \|w\|_{1,\bar{K}_1} + \|w\|_{0,\partial\bar{K}}^2 \le C \|w\|_{0,\bar{K}} \|w\|_{1,\bar{K}}.
\end{aligned}
$$

The result of the lemma now follows by scaling. The proof in three dimensions is similar.

<figure><img src="image_5.png" /></figure>

Proof of Lemma 4.1. Starting from the definition of the norm (10) we have three terms to estimate. Beginning with the interior contributions we find that

$$
\|\nabla(v - I^h v)\|_{0,\Omega_i} \le \|\nabla(v - I_i^h v)\|_{0,\Omega_i^*} \le C h^p |v|_{p+1,\Omega_i^*} \le C h^p |v|_{p+1,\Omega}, \quad (16)
$$

for $i = 1, 2$. Here we have used the fact that $\Omega_i \subset \Omega_i^* \subset \Omega$ and the interpolation error estimate (12).

Next, for the contribution from the jump at the interface, we note that

$$
\|[v - I^h v]\|_{1/2, h, \Gamma} \le \|v - I_1^h v\|_{1/2, h, \Gamma} + \|v - I_2^h v\|_{1/2, h, \Gamma}. \quad (17)
$$

and consider first the second term on the right. Let $L(\Gamma_j)$ be the line segment (plane domain) obtained by extending $\Gamma_j$ to the boundary of $K_2^j \in T_2^h$.

By the definition of the discrete 1/2-norm,

$$
\|v - I_2^h v\|_{1/2, h, \Gamma} \le \sum_{K \in T_h^2} \sum_{\Gamma_j \subset K} h_{K_1^j}^{-1/2} \|v - I_2^h v\|_{0, L(\Gamma_j)}. \quad (18)
$$

By Lemma 4.2, we have that

$$
h_{K_1^j}^{-1} \|v - I_2^h v\|_{0, L(\Gamma_j)}^2 \le C h_{K_2^j} h_{K_1^j}^{-1} \left( h_{K_2^j}^{-2} \|v - I_2^h v\|_{0, K_2^j}^2 + \|v - I_2^h v\|_{1, K_2^j}^2 \right). \quad (19)
$$

It follows from assumption A2 that $h_{K_2^j} h_{K_1^j}^{-1} \le C$. Hence by (19) and the interpolation estimate (13) we obtain

$$
h_{K_1^j}^{-1} \|v - I_2^h v\|_{0, L(\Gamma^j)}^2 \le C h^{2p} |v|_{p+1, K_2^j}^2. \quad (20)
$$

Combining (18) and (20) we arrive at the desired estimate

$$
\|v - I_2^h v\|_{1/2, h, \Gamma} \le C h^p |v|_{p+1, \Omega_2^*} \le C h^p |v|_{p+1, \Omega}. \quad (21)
$$

Here we have used that by assumption A2, the number of terms in the inner sum in (18) is bounded, uniformly with respect to the mesh size. The estimate for the first term on the right-hand side in (17) is readily found by similar arguments.

Finally, by Lemma 4.2 with $w = \langle \nabla_n (v - I^h v) \rangle = \nabla_n (v - I_1^h v)$, we have that

$$
h_{K_1^j} \| \langle \nabla_n (v - I^h v) \rangle \|_{0, L(\Gamma^j)}^2 \le C \left( \| \nabla (v - I_1^h v) \|_{1, K_1^j}^2 + h_{K_1^j}^2 \| \nabla (v - I_1^h v) \|_{2, K_1^j}^2 \right).
$$

We find in the same way as above, using the interpolation error estimate (13) and summing the contributions over the interface, that

$$
\|\langle \nabla_n (v - I^h v) \rangle\|_{-1/2, h, \Gamma} \le C h^p |v|_{p+1, \Omega}. \quad (22)
$$

The lemma now follows from (16, 21) and (22).

<figure><img src="image_6.png" /></figure>
## 4.2. Coercivity

To prove coercivity of the bilinear form we need the following known inverse inequality. We include its proof for completeness.

**Lemma 4.3.** For $\phi \in V^h$, the following inverse inequality holds:

$$ \|\langle \nabla_n \phi \rangle\|_{-1/2, h, \Gamma}^2 \le C_I \|\nabla \phi\|_{0, \Omega_1}^2. $$

*Proof.* Recall the definition of $\langle \nabla_n \phi \rangle = \nabla_n \phi_1$. Note that on a reference triangle (tetrahedron) $\tilde{K}$ which is the image of $K \in T_1^h$ under an affine map, we have

$$ \|\langle \nabla_n \phi \rangle\|_{\Gamma}^2 \le C \|\nabla \phi\|_{0, \tilde{K}}^2. $$

since if the right-hand side is zero so is the left hand side, and since the space of polynomials of degree $p-1$ is finite dimensional. The result then follows by scaling, using the inverse of the affine map, and summation over all elements with an edge on the interface $\Gamma$. □

**Lemma 4.4.** The discrete form $a^h(\cdot, \cdot)$ is coercive on $V^h$, i.e.,

$$ a^h(v, v) \ge C \|v\|^2 \quad \forall v \in V^h, $$

provided $\lambda$ is chosen sufficiently large. It is also continuous, i.e.,

$$ a^h(u, v) \le C \|u\| \|v\| \quad \forall u, v \in V. $$

*Proof.* Continuity of the discrete form follows directly from the definitions. To prove coercivity, we use (11) to estimate the form from above:

$$ \begin{aligned} a^h(v, v) &= \|\nabla v\|_{0, \Omega_1}^2 + \|\nabla v\|_{0, \Omega_2}^2 - 2([v], \langle \nabla_n v \rangle)_\Gamma + \|\lambda^{1/2} [v]\|_{1/2, h, \Gamma}^2 \\ &\ge \|\nabla v\|_{0, \Omega_1}^2 + \|\nabla v\|_{0, \Omega_2}^2 - 2 \|\langle \nabla_n v \rangle\|_{-1/2, h, \Gamma} \|[v]\|_{1/2, h, \Gamma} + \lambda \|[v]\|_{1/2, h, \Gamma}^2. \end{aligned} $$

It follows from the inverse inequality in Lemma 4.3 that

$$ \begin{aligned} -2 \|\langle \nabla_n v \rangle\|_{-1/2, h, \Gamma} \|[v]\|_{1/2, h, \Gamma} &\ge -\epsilon \|\{\nabla_n v\}\|_{-1/2, h, \Gamma}^2 - \epsilon^{-1} \|[v]\|_{1/2, h, \Gamma}^2 \\ &\ge -\epsilon C_I \|\nabla v\|_{0, \Omega_1}^2 - \epsilon^{-1} \|[v]\|_{1/2, h, \Gamma}^2. \end{aligned} $$

Combining these estimates we obtain

$$ a^h(v, v) \ge (1 - \epsilon C_I) \|\nabla v\|_{0, \Omega_1}^2 + \|\nabla v\|_{0, \Omega_2}^2 + (\lambda - \epsilon^{-1}) \|[v]\|_{1/2, h, \Gamma}^2. $$

Taking $\epsilon = 1/(2C_I)$, coercivity follows if $\lambda \ge 1/2 + 1/\epsilon$. □

## 4.3. A priori error estimates

**Theorem 4.5.** For $U$ solving (8) and $u$ solving (1, 2), the following a priori error estimates hold

$$ \|u - U\| \le C h^p |u|_{p+1, \Omega}, \qquad (23) $$

and

$$ \|u - U\|_{0, \Omega} \le C h^{p+1} |u|_{p+1, \Omega}. \qquad (24) $$
Proof. For any $v \in V^h$, $\|u - U\| \le \|u - v\| + \|v - U\|$. Further, by Lemma 4.4 and orthogonality, we have that

$$ \begin{aligned} \|U - v\|^2 &\le Ca^h(U - v, U - v) \\ &= Ca^h(u - v, U - v) \\ &\le C \|u - v\| \|U - v\|, \end{aligned} $$

and it follows that

$$ \|u - U\| \le C \|u - v\| \quad \forall v \in V^h. $$

Taking $v = I^h u$ and invoking the interpolation result of Theorem 4.1, (23) follows.

For (24) we use a duality argument. Let $z$ be defined by

$$ \begin{aligned} -\Delta z &= e \quad \text{in } \Omega, \\ z &= 0 \quad \text{on } \partial\Omega, \end{aligned} \quad (25) $$

with $e = u - U$. Multiplying (25) with $e$ and using Green's formula gives

$$ \begin{aligned} \|e\|_{0,\Omega}^2 &= -(\Delta z, e)_\Omega \\ &= (\nabla z, \nabla e)_{\Omega_1} + (\nabla z, \nabla e)_{\Omega_2} - (\nabla_n z, [e])_\Gamma \\ &= a^h(z, e), \end{aligned} $$

since $[z] = 0$ and $\nabla_n z = \langle \nabla_n z \rangle$. Now, using the symmetry of $a^h(\cdot, \cdot)$, the orthogonality relation (9), and Theorem 4.1, we find that

$$ \|e\|_{0,\Omega}^2 = a^h(z - I^h z, e) \le C \|z - I^h z\| \|e\| \le Ch \|z\|_{2,\Omega} \|e\|. \quad (26) $$

Finally, by elliptic regularity, we have

$$ \|z\|_{2,\Omega} \le C \|e\|_{0,\Omega}, $$

and thus the estimate (24) follows from (23) and (26). $\square$

## 5. A POSTERIORI ERROR ESTIMATES

We first introduce an interpolation operator suitable for the a posteriori error analysis. From Scott and Zhang [13] we deduce the existence of Clement type interpolation operators $r_{h,i}$, $i = 1, 2$, defined on $H^1(\Omega_i^*)$, which preserve Dirichlet boundary conditions on $\partial\Omega \cap \bar{\Omega}_i^*$ and satisfy the following local interpolation error estimates.

$$ h_K^{m-1} \|r_{h,i} v - v\|_{m,K} \le C \|\nabla v\|_{0,\Delta K}, \quad m = 0, 1, K \in T_{h,i}. \quad (27) $$

Here $\Delta K$ denotes a patch of elements which are neighbors of $K$. We then define

$$ \pi v := (r_{h,1}v_1, (r_{h,2}\mathbb{E}_2 v_2)|_{\Omega_2}) \quad \text{for } v \in H^1(\Omega_1 \cup \Omega_2). \quad (28) $$

Here we have used an extension operator $\mathbb{E}_2: H^1(\Omega_2) \to H^1(\Omega)$ such that $(\mathbb{E}_2 w)|_{\Omega_2} = w$ and

$$ |\mathbb{E}_2 w|_{m,\Omega} \le C |w|_{m,\Omega_i} \quad \forall w \in H^m(\Omega_2), \quad m = 0, 1. \quad (29) $$

In the following lemma we collect some useful estimates for the interpolation operator $\pi$.
**Lemma 5.1.** Let $\pi$ be the interpolation operator defined in (28). Under assumptions A1 and A2 the following estimates hold for $v \in H^1(\Omega_1 \cup \Omega_2)$:

$$ \sum_{K \in T_i^h} h_K^{-2} \|v_i - (\pi v)_i\|_{K \cap \Omega_i}^2 \le C \|\nabla v_i\|_{0, \Omega_i}^2, \quad i = 1, 2, \qquad (30) $$

$$ \|v_i - (\pi v)_i\|_{1/2, h, \Gamma} \le C \|\nabla v_i\|_{0, \Omega_i}, \quad i = 1, 2, \qquad (31) $$

$$ \sum_{K \in T_i^h} h_K^{-1} \|v - \pi v\|_{0, \partial(K \cap \Omega_i) \setminus \Gamma}^2 \le C \|\nabla v_i\|_{0, \Omega_i}^2, \quad i = 1, 2, \qquad (32) $$

$$ \|\nabla(\pi v)\|_{0, \Omega_i} \le C \|\nabla v\|_{0, \Omega_i}, \quad i = 1, 2. \qquad (33) $$

*Proof.* We consider the case $i = 2$ as the proof for $i = 1$ is similar. Let $v_2^* = E_2 v_2$ denote the extension of $v_2$ to $\Omega$ and recall that $v_2 - (\pi v)_2 = v_2^* - r_{h,2}v_2^*$ on $\Omega_2$. Using (27) with $m = 0$ we obtain, for $K \in T_2^h$, that

$$ h_K^{-1} \|v_2 - (\pi v)_2\|_{0, K \cap \Omega_2} \le \|v_2^* - r_{h,2}v_2^*\|_{0, K} \le C \|\nabla v_2^*\|_{0, \Delta K}. $$

As the number of elements in the patches $\Delta K$ is uniformly bounded with respect to the mesh size by assumption A1, it follows that

$$ \sum_{K \in T_i^h} h_K^{-2} \|v_i - (\pi v)_i\|_{K \cap \Omega_2}^2 \le C \|\nabla v_2^*\|_{0, \Omega_2^*}^2, $$

whence (30) follows by the bound for the extentions operator (29). Turning to the second inequality (31) of the lemma, it follows from trace inequality (15) and interpolation estimate (30) that

$$ \begin{aligned} h_{K_1^j}^{-1} \|v_2 - (\pi v)_2\|_{0, \Gamma^j}^2 &\le C h_{K_2^j}^{-2} \|v_2 - (\pi v)_2\|_{0, K_2}^2 + \|v_2 - (\pi v)_2\|_{1, K_2}^2 \\ &\le C \|\nabla v_2\|_{0, \Delta K_2}^2, \end{aligned} $$

where we also have used that $h_{K_2^j} h_{K_1^j}^{-1} \le C$ by assumption A2. Summing over the elements and using (29) gives

$$ \|v_2 - (\pi v)_2\|_{1/2, h, \Gamma} \le C \|\nabla v_2^*\|_{0, \Omega_2^*} \le C \|\nabla v\|_{0, \Omega_2}. $$

The third estimate (32) is readily shown by similar arguments. Finally, for the fourth inequality of the lemma, it follows from (27) and (29) that

$$ \|\nabla(\pi v)\|_{0, \Omega_2} \le \|\nabla(v_2^* - r_{h,2}v_2^*)\|_{0, \Omega_2^*} + \|\nabla v\|_{0, \Omega_2} \le C \|\nabla v\|_{0, \Omega_2}, $$

and the proof is complete.

<figure><img src="image_7.png" /></figure>

We shall also need the following trace inequality.

**Lemma 5.2.** Under our mesh assumptions there holds, for all $K \in T_1^h$ and $\Gamma_j \subset \partial K$,

$$ \|\nabla n v_1\|_{-1/2, \Gamma^j} \le C \left( \|\nabla v_1\|_{0, K} + h_K \|\Delta v_1\|_{0, K} \right), \quad v \in H^2(K). \qquad (34) $$

*Proof.* On a reference element $\tilde{K}$ there holds (cf. [8], Th. 2.2)

$$ \|n \cdot w\|_{-1/2, \partial \tilde{K}} \le C \left( \|w\|_{0, \tilde{K}} + \|\nabla \cdot w\|_{0, \tilde{K}} \right), \quad w \in L_2(\tilde{K})^n : \nabla \cdot w \in L_2(\tilde{K}). $$

The result follows from this estimate, scaled by the map from the reference element, with $w = \nabla v$. □

<figure><img src="image_8.png" /></figure>
We are now ready to show an a posteriori error estimate in a discrete energy norm, using the following notation. At an edge of an element K that is common with a neighbouring element L we let $[w] = w|_K - w|_L$ denote the jump of w over the edge. Further, $n_P$ denotes the outward pointing unit normal of the boundary of a domain P.

**Theorem 5.3.** Assume A1, A2, and $\lambda \ge 1$. For U solving (8) and u solving (3), the following a posteriori error estimate holds:

$$ \|\nabla e\|_{0,\Omega_1 \cup \Omega_2}^2 + \|[e]\|_{1/2,h,\Gamma}^2 \le C \sum_{i=1}^2 \sum_{K \in T_i^h} \rho_{K,i}^2. \quad (35) $$

Here the element error indicators $\rho_{K,i}$ are defined by

$$ \rho_{K,i}^2 = h_K^2 \|f + \Delta U\|_{0,P_K}^2 + h_K \|[n_{P_K} \cdot \nabla U]\|_{0,\partial P_K}^2 + h_K^{-1} \|[U]\|_{0,\partial P_K \cap \Gamma}^2 + \sum_{\Gamma_j \subset \bar{K}} \|[U]\|_{1/2,\Gamma_j}^2, $$

where $P_K = K \cap \Omega_i$ for $K \in T_i^h$.

*Proof.* Using the definition of the method (8) we have the identity

$$ I = (\nabla e, \nabla e)_{\Omega_1 \cup \Omega_2} = (\nabla e, \nabla(e - \pi e))_{\Omega_1 \cup \Omega_2} + (\langle \nabla n e \rangle, [\pi e])_{\Gamma} + ([e], \langle \nabla n \pi e \rangle)_{\Gamma} - (\lambda h^{-1}[e], [\pi e])_{\Gamma}. \quad (36) $$

Integration by parts gives

$$ \begin{aligned} I = & \sum_{i=1}^{2} \sum_{K \in T_i^h} \left( (-\Delta e, (e - \pi e))_{P_K} + \frac{1}{2} ([n_K \cdot \nabla e], e - \pi e)_{\partial P_K \setminus \Gamma} \right) \\ & + (\langle \nabla n e \rangle, [e - \pi e])_{\Gamma} + ([\nabla n e], e_2 - \pi e_2)_{\Gamma} \\ & + (\langle \nabla n e \rangle, [\pi e])_{\Gamma} + ([e], \langle \nabla n \pi e \rangle)_{\Gamma} - (\lambda h^{-1}[e], [\pi e])_{\Gamma}. \end{aligned} \quad (37) $$

Note that we may write

$$ -(\lambda h^{-1}[e], [\pi e])_{\Gamma} = (\lambda h^{-1}[e], [e - \pi e])_{\Gamma} - \|[e]\|_{1/2,h,\Gamma}^2 - ((\lambda - 1)h^{-1}[e], [e])_{\Gamma}. $$

Using that u and its derivatives have zero jumps over the element edges and the interface, and also that $-\Delta e = f + \Delta U$ in the interior of the elements, we obtain the error representation formula

$$ \begin{aligned} \|\nabla e\|_{0,\Omega_1 \cup \Omega_2}^2 + \|[e]\|_{1/2,h,\Gamma}^2 = & \sum_{i=1}^2 \sum_{K \in T_i^h} \left( (f + \Delta U, (e - \pi e))_{P_K} + \frac{1}{2} ([n_K \cdot \nabla U], e - \pi e)_{\partial P_K \setminus \Gamma} \right) \\ & + ([\nabla_n U], e_2 - \pi e_2)_{\Gamma} + ([U], \langle \nabla_n e \rangle)_{\Gamma} + ([U], \langle \nabla_n \pi e \rangle)_{\Gamma} \\ & + (\lambda h^{-1}[U], [e - \pi e])_{\Gamma} - ((\lambda - 1)h^{-1}[e], [e])_{\Gamma}. \end{aligned} \quad (38) $$

Here the last term is non-positive under our assumptions. We now proceed to estimate the other terms on the right-hand side. For $\epsilon > 0$ to be chosen below we have that

$$ \begin{aligned} (f + \Delta U, (e - \pi e))_{P_K} & \le \|(f + \Delta U)_{0,P_K} \|_{0,P_K} \|e - \pi e\|_{0,P_K} \\ & \le \frac{1}{2\epsilon} h_K^2 \|f + \Delta U\|_{0,P_K}^2 + \frac{\epsilon}{2} h_K^{-2} \|e - \pi e\|_{0,P_K}^2. \end{aligned} $$
Treating the other term in the sum over the elements in the same way, we find by Lemma 5.1 that the first term in (38) may be estimated as follows.

$$ 
\begin{aligned}
\sum_{K \in T_i^h} & \left( (f + \Delta U, (e - \pi e))_{P_K} + \frac{1}{2} ([n_K \cdot \nabla U], e - \pi e)_{\partial P_K \setminus \Gamma} \right) \\
& \le \frac{C}{\epsilon} \sum_{K \in T_i^h} \left( h_K^2 \|f + \Delta U\|_{0,P_K}^2 + h_K \|[n_K \cdot \nabla U]\|_{0,\partial P_K \setminus \Gamma}^2 \right) \\
& \qquad + C\epsilon \sum_{K \in T_i^h} \left( h_K^{-2} \|e - \pi e\|_{0,P_K}^2 + h_K^{-1} \|e - \pi e\|_{0,\partial P_K \setminus \Gamma}^2 \right) \\
& \le \frac{C}{\epsilon} \sum_{K \in T_i^h} \left( h_K^2 \|f + \Delta U\|_{0,P_K}^2 + h_K \|[n_K \cdot \nabla U]\|_{0,\partial P_K \setminus \Gamma}^2 \right) \\
& \qquad + C\epsilon \|\nabla e\|_{0,\Omega_i}^2.
\end{aligned} 
\tag{39} $$

Using again Lemma 5.1, we obtain for the second term $([\nabla_n U], e_2 - \pi e_2)_{\Gamma}$ that

$$ 
\begin{aligned}
([\nabla_n U], e_2 - \pi e_2)_{\Gamma} &\le C \|[\nabla_n U]\|_{-1/2,h,\Gamma} \|e_2 - \pi e_2\|_{1/2,h,\Gamma} \\
&\le C \|[\nabla_n U]\|_{-1/2,h,\Gamma} \|\nabla e\|_{0,\Omega_2} \\
&\le \frac{C}{\epsilon} \|[\nabla_n U]\|_{-1/2,h,\Gamma}^2 + C\epsilon \|\nabla e\|_{0,\Omega_2}^2 \\
&\le \frac{C}{2\epsilon} \sum_{i=1}^2 \sum_{K \in T_i^h} h_K \|[\nabla_n U]\|_{\partial P_K \cap \Gamma}^2 + C\epsilon \|\nabla e\|_{0,\Omega_2}^2.
\end{aligned} 
\tag{40} $$

In the last step above, we have merely distributed the error indicator at $\Gamma$ over the elements intersected by $\Gamma$, using that the element sizes are compatible by assumption A2. For the third term $([U], \langle \nabla_n e \rangle)_{\Gamma}$, we begin by noting that

$$ 
\begin{aligned}
([U], \langle \nabla_n e \rangle)_{\Gamma} &\le \sum_j \|[U]\|_{1/2, \Gamma_j} \|\nabla_n e_1\|_{-1/2, \Gamma_j} \\
&\le \sum_j \frac{1}{2\epsilon} \|[U]\|_{1/2, \Gamma_j}^2 + \frac{\epsilon}{2} \|\nabla_n e_1\|_{-1/2, \Gamma_j}^2.
\end{aligned} 
$$

By the trace inequality of Lemma 5.2 we have that

$$ 
\begin{aligned}
\|\nabla_n e_1\|_{-1/2, \Gamma_j}^2 &\le C(\|\nabla e\|_{0, K_1^j}^2 + h_{K_1^j}^2 \|\nabla \cdot \nabla e\|_{0, K_1^j}^2) \\
&= C(\|\nabla e\|_{0, K_1^j}^2 + h_{K_1^j}^2 \|f + \Delta U\|_{0, K_1^j}^2),
\end{aligned} 
$$

and hence, since the number of parts $\Gamma_j$ in each element is uniformly bounded,

$$ 
\begin{aligned}
([U], \langle \nabla_n e \rangle)_{\Gamma} &\le \frac{C}{\epsilon} \sum_j \left( h_{K_1^j}^2 \|f + \Delta U\|_{0, K_1^j}^2 + \|[U]\|_{1/2, \Gamma_j}^2 \right) + C\epsilon \|\nabla e\|_{0, \Omega_1}^2 \\
&\le \frac{C}{\epsilon} \sum_{i=1}^2 \sum_{K \in T_i^h} \left( h_K^2 \|f + \Delta U\|_{0, K_1^j}^2 + \sum_{\Gamma_j \subset \bar{K}} \|[U]\|_{1/2, \Gamma_j}^2 \right) + C\epsilon \|\nabla e\|_{0, \Omega_1}^2.
\end{aligned} 
\tag{41} $$
The fourth term $([U], \langle \nabla_n \pi e \rangle)_\Gamma$ is bounded as follows

$$([U], \langle \nabla_n \pi e \rangle)_\Gamma \le \|[U]\|_{1/2,h,\Gamma} \| \langle \nabla_n \pi e \rangle \|_{-1/2,h,\Gamma}.$$

Further using the inverse inequality of Lemma 4.3 and (33) we obtain

$$\| \langle \nabla_n \pi e \rangle \|_{-1/2, h, \Gamma} \le C \| \nabla \pi e \|_{0, \Omega_1}^2 \le C \| \nabla e \|_{0, \Omega_1}^2.$$

Likewise, for the fifth and last term we find by Lemma 5.1 that

$$(h^{-1}[U], [e - \pi e])_\Gamma \le \|[U]\|_{1/2,h,\Gamma} \|[e - \pi e]\|_{1/2,h,\Gamma} \le C \|[U]\|_{1/2,h,\Gamma} \|\nabla e\|_{0,\Omega_1 \cup \Omega_2}.$$

Thus we get for the fourth and fifth terms that

$$([U], \langle \nabla_n \pi e \rangle)_\Gamma + (h^{-1}[U], [e - \pi e])_\Gamma \le C \|[U]\|_{1/2,h,\Gamma} \|\nabla e\|_{0,\Omega_1 \cup \Omega_2} \\ \le \frac{C}{\epsilon} \sum_{i=1}^2 \sum_{K \in T_i^h} h_K^{-1} \|[U]\|_{0,\partial P_K \cap \Gamma}^2 + C\epsilon \|\nabla e\|_{0,\Omega_1 \cup \Omega_2}^2. \quad (42)$$

Choosing $\epsilon$ small enough, the theorem now follows from (38), (39), (40), (41), and (42). $\square$

The error indicator of Theorem 5.3 contains a term $D_K := h_K^{-1} \|[U]\|_{0,\partial P_K \cap \Gamma}^2$ corresponding to a discrete 1/2-norm over $\Gamma$, as well as a term

$$S_K := \sum_{\Gamma_j \subset \bar{K}} \|[U]\|_{1/2, \Gamma_j}^2$$

with continuous 1/2-norms over the parts $\Gamma_j$. We shall consider the computation of $S_K$ in the special case of linear elements in two dimensions in Section 5 below, using that for a one-dimensional interface (see, e.g., [4])

$$\|[U]\|_{1/2, \Gamma_j}^2 := \|[U]\|_{0, \Gamma_j}^2 + \int_{\Gamma_j} \int_{\Gamma_j} \frac{|[U](\xi) - [U](\eta)|^2}{|\xi - \eta|^2} d\xi d\eta. \quad (43)$$

In a general case, however, $S_K$ is somewhat complicated to compute and one would like to simplify the error indicators. Indeed, $S_K$ is bounded by (but not equal to) $\|[U]\|_{1/2,\partial P_K \cap \Gamma}^2$, and it is therefore natural to ask if an inverse inequality can be found which would make possible to remove $S_K$ from the error indicator. We note that even though we have the inverse inequality $\|[U]\|_{1/2,\Gamma_j}^2 \le C|\Gamma_j|^{-1}\|[U]\|_{0,\Gamma_j}^2$, the corresponding elementwise inverse inequality $\|[U]\|_{1/2,\partial P_K \cap \Gamma}^2 \le C h_K^{-1} \|[U]\|_{0,\partial P_K \cap \Gamma}^2$ does not follow since the parts $\Gamma_j$ may become arbitrarily small even when the meshsize is bounded away from zero. Nevertheless, $S_K$ and $D_K$ may, as we shall show, both be bounded by a third quantity of the same order of magnitude. In Theorem 5.4 below we use this and some further simplifications to obtain a less sharp but more implementation-friendly a posterori error estimate. Integration over the parts $P_K$ is not required to compute these error indicators, nor is integration over the parts $\Gamma_j$; all terms are integrals of single polynomials over the original elements or its edges.

**Theorem 5.4.** For any piecewise polynomial $w$ on the partition $\{\Gamma_j\}$ of $\Gamma$, let $w_j$ denote the polynomial which defines $w$ on $\Gamma_j$. Theorem 5.3 holds with $\rho_{K,i}^2$ replaced by

$$\theta_{K,i}^2 = h_K^2 \|f + \Delta U\|_{0,K}^2 + h_K \|[n_K \cdot \nabla U]\|_{0,\partial K \cap \Omega_i^*}^2 \\ + \sum_{j: \Gamma_j \subset \bar{K}} \left( h_K \|[\nabla_n U]_j\|_{0,\partial K_1^j \cap \Gamma}^2 + h_K^{-1} \|[U]_j\|_{0,\partial K_1^j \cap \Gamma}^2 \right).$$
Proof. We shall show that $\rho_{K,i} \le C\theta_{K,i}$. Obviously,

$$ \|f + \Delta U\|_{0,P_K}^2 \le \|f + \Delta U\|_{0,K}^2. \qquad (44) $$

Further, for $K \in T_i^h$ we have that

$$ \begin{aligned} h_K \|[n_{P_K} \cdot \nabla U]\|_{0, \partial P_K}^2 &= h_K \|[n_K \cdot \nabla U]\|_{0, \partial P_K \setminus \Gamma}^2 + \sum_{\Gamma_j \subset \bar{K}} h_K \|[n_{P_K} \cdot \nabla U]\|_{0, \Gamma_j}^2 \\ &\le h_K \|[n_K \cdot \nabla U]\|_{0, \partial K \cap \Omega_i^*}^2 + h_K \sum_{j: \Gamma_j \subset \bar{K}} \|[n \cdot \nabla U]_j\|_{0, \partial K_1^j \cap \Gamma}^2. \end{aligned} \qquad (45) $$

We now turn to the estimate of $S_K = \sum_{\Gamma_j \subset \bar{K}} \|[U]\|_{1/2, \Gamma_j}^2$. For $K \in T_i^h$ and $\Gamma_j \subset \bar{K}$ we also have $\Gamma_j \subset \bar{K}_1^j$. Consider a reference element $\tilde{K}$ for $K_1^j$ with $\Gamma_j$ mapped onto $\tilde{\Gamma}_j$ in the edge (side) $\tilde{E}$. By equivalence of norms there holds for all polynomials $q$ of degree $p$

$$ \|q\|_{1/2, \tilde{\Gamma}_j}^2 \le \|q\|_{1/2, \tilde{E}}^2 \le C \|q\|_{0, \tilde{E}}^2. $$

Mapping back to $K_1^j$ and taking $q = [U]_j$ we find that

$$ \|[U]_j\|_{1/2, \Gamma_j}^2 \le C h_{K_1^j}^{-1} \|[U]_j\|_{0, \partial K_1^j \cap \Gamma}^2, $$

where the constant is independent of $j$. Hence, using assumption A2 to replace $h_{K_1^j}$ by $h_K$ if $i=2$,

$$ S_K = \sum_{\Gamma_j \subset \bar{K}} \|[U]\|_{1/2, \Gamma_j}^2 \le C h_K^{-1} \sum_{j: \Gamma_j \subset \bar{K}} \|[U]_j\|_{0, \partial K_1^j \cap \Gamma}^2. \qquad (46) $$

Since $D_K = h_K^{-1}\|[U]\|_{0, \partial P_K \cap \Gamma}^2$ is indeed bounded by the right-hand side above, the theorem now follows from (44, 45) and (46). $\square$

# 6. NUMERICAL EXAMPLES

## 6.1. Implementation

For the numerical examples to be presented, we chose to make the following simplifying assumption: the interface is assumed to be made up of straight lines that are so long that each element on the cut grid is intersected by just one corner of the interface. We further assumed that the area enclosed by the interface lay completely within the cut mesh. We then needed to consider only seven cases, depicted in Figures 3–5.

For the numerical integration, we divided any non-triangular cut part of an element into triangles, using sufficiently high integration for the stiffness matrix to be exactly evaluated. On the interface we performed exact (numerical) integration of all terms using the union of the intersection points on the cut side and the nodes on the uncut side to define the intervals of integration.
<figure><img src="image_9.png" /><figcaption>FIGURE 3. Elements containing one corner of the interface: zero nodes, one node, and two nodes on the opposite side of the interface.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>FIGURE 4. Elements cut by a straight segment: one node and two nodes on the opposite side of the interface.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>FIGURE 5. Overlapped elements: two sides and three sides cut.</figcaption></figure>

## 6.2. Convergence study

The example was solved on the domain $(-4, 4) \times (-4, 4)$ with an overlapping domain according to Figure 6. The locations of the corners of the overlapping domain were determined by

$$ \begin{bmatrix} x \\ y \end{bmatrix} = \begin{bmatrix} \cos \phi & -\sin \phi \\ \sin \phi & \cos \phi \end{bmatrix} \begin{bmatrix} -2.5 & 0.5 & 0.5 & -2.5 \\ -2.5 & -2.5 & 0.5 & 0.5 \end{bmatrix} + \begin{bmatrix} 0.75 \\ 1 \end{bmatrix}, $$

with $\phi = 1$.

We imposed zero Dirichlet boundary conditions and applied a forcing term

$$ f = 64 - 2x^2 - 2y^2, $$

corresponding to the exact solution

$$ u = (x - 4)(x + 4)(y - 4)(y + 4). $$

We solved this problem numerically using finite elements with both a linear and a quadratic polynomial *ansatz*. Elevations of the different solutions on the same mesh can be seen in Figure 7. In Figure 8, we show the
<figure><img src="image_12.png" /><figcaption>FIGURE 6. Coarse triangulation of the domains.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>FIGURE 7. Elevation of the linear and the quadratic solutions on the coarse mesh.</figcaption></figure>

convergence in $L_2$-norm of on a sequence of refined meshes. As expected, we obtain second and third order convergence for the linear and quadratic approximations, respectively.

## 6.3. Adaptive refinement

In order to evaluate the *a posteriori* estimate (35) we need to compute the continuous half-norm (43) of the jump. We have implemented the error estimator for linear elements, whence we compute the half-norm as follows: since $[U]$ is linear on each $\Gamma_j$, $[U](\xi) = a\xi + b$, say, we can write

$$ \frac{|[U](\xi) - [U](\eta)|^2}{|\xi - \eta|^2} = \frac{|a(\xi - \eta)|^2}{|\xi - \eta|^2} = a^2 $$

and thus

$$ \|[U]\|_{1/2, \Gamma_j}^2 = \|[U]\|_{0, \Gamma_j}^2 + a^2 |\Gamma_j|^2. $$
<figure><img src="image_14.png" /><figcaption>FIGURE 8. Convergence of the different approximations.</figcaption></figure>

We considered an example with exact solution

$$u = \frac{1}{256}e^{-(x^2+y^2)}(4-x)(4+x)(4-y)(4+y),$$

on the same domain as in the previous example, and used the two a posteriori estimates of Theorem 5.3, called estimate 1, and Theorem 5.4, called estimate 2. These estimates were used to control the local meshsize by a fixed fraction refinement strategy (refining the elements containing the largest contributions). No attempt was made to calibrate the unknown constants appearing in the estimates; instead we computed the effectivity index, defined as

$$\text{effectivity index} = \frac{\text{estimated error}}{\text{exact error}},$$

in order to numerically verify that it remained constant as the mesh was refined. In Figure 9, we show the mesh obtained when using estimate 1, and in Figure 10 the corresponding result when using estimate 2. The solution on this mesh is shown as an elevation in Figure 11. Note that the meshes in both cases has a tendency to refine more at the interface. This is because the local error is the largest there, as has been noted previously, cf. [2,9]. This phenomenon is more noticeable using estimate 2, as expected since this estimate is more conservative at the interface. Finally, in Figure 12, we show the effectivity index on the sequence of meshes for both estimators. As can be seen, the effectivity index is almost constant for both estimates, which indicates that the refinement on the interface does not much affect the global norm. We also show the ratio of the errors on the sequence of meshes obtained with estimate 2 to that obtained with estimate 1. We note that this ratio is slightly above 1, indicating that the degrees of freedom are put to slightly better use using estimate 1, as expected.
<figure><img src="image_15.png" /><figcaption>FIGURE 9. Final adapted mesh in a sequence using Estimate 1.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>FIGURE 10. Final adapted mesh in a sequence using Estimate 2.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>FIGURE 11. Elevation of the approximate solution.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>FIGURE 12. Effectivity index obtained on a sequence of adaptively refined meshes.</figcaption></figure>

REFERENCES

[1] Y. Achdou and Y. Maday, The mortar element method with overlapping subdomains. *SIAM J. Numer. Anal.* **40** (2002) 601–628.

[2] R. Becker, P. Hansbo and R. Stenberg, A finite element method for domain decomposition with non-matching grids. *ESAIM: M2AN* **37** (2003) 209–225.

[3] M.J. Berger, On conservation at grid interfaces. *SIAM J. Numer. Anal.* **24** (1987) 967–984.

[4] S.C. Brenner and L.R. Scott, *The Mathematical Theory of Finite Element Methods*. Springer-Verlag, Berlin (1994).

[5] F. Brezzi, J.-L. Lions and O. Pironneau, Analysis of a Chimera method. *C. R. Acad. Sci. Paris Sér. I Math.* **332** (2001) 655–660.

[6] X.-C. Cai, M. Dryja and M. Sarkis, Overlapping nonmatching grid mortar element methods for elliptic problems. *SIAM J. Numer. Anal.* **36** (1999) 581–606.
[7] G. Chesshire and W.D. Henshaw, Composite overlapping meshes for the solution of partial-differential equations. *J. Comput. Phys.* **90** (1990) 1–64.

[8] V. Girault and P.-A. Raviart, *Finite Element Approximation of the Navier-Stokes Equations*. Springer-Verlag, Berlin (1979).

[9] A. Hansbo and P. Hansbo, An unfitted finite element method, based on Nitsche's method, for elliptic interface problems. *Comput. Methods Appl. Mech. Engrg.* **191** (2002) 5537–5552.

[10] R.D. Lazarov, J.E. Pasciak, J. Schöberl and P.S. Vassilevski, Almost optimal interior penalty discontinuous approximations of symmetric elliptic problems on non-matching grids. Technical Report, ISC-01-05-MATH (2001).

[11] R.D. Lazarov, S.Z. Tomov and P.S. Vassilevski, Interior penalty discontinuous approximations of elliptic problems. *Comput. Meth. Appl. Math.* **1** (2001) 367–382.

[12] J. Nitsche, Über ein Variationsprinzip zur Lösung von Dirichlet-Problemen bei Verwendung von Teilräumen, die keinen Randbedingungen unterworfen sind. *Abh. Math. Sem. Univ. Hamburg* **36** (1971) 9–15.

[13] L.R. Scott and S. Zhang, Finite element interpolation of nonsmooth functions satisfying boundary conditions. *Math. Comp.* **190** (1990) 483–493.

[14] R. Stenberg, Mortaring by a method of J.A. Nitsche, in *Computational Mechanics: New Trends and Applications*, S. Idelsohn, E. Onate and E. Dvorkin Eds., CIMNE, Barcelona (1998).