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

HAL
open science

# Adaptive boundary element methods for regularized combined field integral equations

Théophile Chaumont-Frelet, Gregor Gantner

► To cite this version:

Théophile Chaumont-Frelet, Gregor Gantner. Adaptive boundary element methods for regularized combined field integral equations. 2025. hal-04917182

HAL Id: hal-04917182

https://inria.hal.science/hal-04917182v1

Preprint submitted on 28 Jan 2025

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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# ADAPTIVE BOUNDARY ELEMENT METHODS FOR REGULARIZED COMBINED FIELD INTEGRAL EQUATIONS

## THÉOPHILE CHAUMONT-FRELET* AND GREGOR GANTNER†

**ABSTRACT.** While the exterior Helmholtz problem with Dirichlet boundary conditions is always well-posed, the associated standard boundary integral equations are not if the squared wavenumber agrees with an eigenvalue of the interior Dirichlet problem. Combined field integral equations are not affected by this spurious resonances but are essentially restricted to sufficiently smooth boundaries. For general Lipschitz domains, the latter integral equations are applicable through suitable regularization. Under fairly general assumptions on the regularizing operator, we propose a *posterirori* computable error estimators for corresponding Galerkin boundary element methods of arbitrary polynomial degree. We show that adaptive mesh-refining algorithms steered by these local estimators converge at optimal algebraic rate with respect to the number of underlying boundary mesh elements. In particular, we consider mixed formulations involving the inverse Laplace-Beltrami as regularizing operator. Numerical examples highlight that in the vicinity of spurious resonances the proposed adaptive algorithm is significantly more performant when applied to the regularized combined field equation rather than the standard one.

## 1. INTRODUCTION

1.1. **Model problem.** The propagation of waves is of central importance in a large array of applications in physics and engineering. Many situations of interest may be modeled assuming a fixed frequency $k > 0$, and considering the reflections of an incident field $u$ on a bounded obstacle $\Omega \subset \mathbb{R}^d$, $d \in \{2, 3\}$, surrounded by an infinite homogeneous medium $\Omega^{\text{ext}}$. In this case, reflected waves can be modeled by a complex-valued amplitude $U : \Omega^{\text{ext}} \to \mathbb{C}$ satisfying the Helmholtz equation

$$ \Delta U + k^2 U = 0 \text{ in } \Omega^{\text{ext}}, \qquad (1.1\text{a}) $$

the condition

$$ U = u \text{ on } \Gamma \qquad (1.1\text{b}) $$

on the compact boundary $\Gamma := \partial\Omega$ of the obstacle $\Omega$, and the Sommerfeld radiation condition

$$ \partial_r U - ikU = o(r^{(1-d)/2}) \text{ as } r \to \infty, \qquad (1.1\text{c}) $$

where $r$ denotes the radial coordinate, and $\partial_r$ the radial derivative.

<footnote>*INRIA UNIV. LILLE AND LABORATOIRE PAUL PAINLEVÉ, 59655 VILLENEUVE-D'ASCQ, FRANCE</footnote>

<footnote>†INSTITUTE FOR NUMERICAL SIMULATION, UNIVERSITY OF BONN, FRIEDRICH-HIRZEBRUCH-ALLEE 7, 53115 BONN, GERMANY</footnote>

<footnote>E-mail addresses: theophile.chaumont@inria.fr, gantner@ins-uni.bonn.de.</footnote>

<footnote>Date: January 28, 2025.</footnote>

<footnote n="1991">1991 Mathematics Subject Classification. 35J05, 65N15, 65N30, 65N38, 65N50.</footnote>

<footnote>Key words and phrases. exterior Helmholtz problem; combined field integral equations; boundary element methods; a *posteriori* error estimation; convergence of adaptive algorithms.</footnote>
### 1.2. Standard boundary element methods.
For general obstacles $\Omega$, the solution to (1.1) cannot be found analytically, and needs to be approximated numerically. In the setting considered here, boundary element methods (BEM) are particularly attractive, since they can be directly formulated on the boundary $\Gamma$ and naturally handle the Sommerfeld radiation condition [43, 48, 47]. In this work, we consider adaptive BEM which, given an initial triangulation $\mathcal{T}_0$ of the boundary, allows for an automatically steered localized mesh refinement tailored to the solution. This approach results in accurate and reliable approximations potentially requiring much less degrees of freedom than uniform meshes.

To employ a BEM discretization, a natural direct approach is to consider the first-kind integral equation

$$V_k \phi = (K_k - 1/2)u, \qquad (1.2)$$

where $V_k$ and $K_k$ are single- and double-layer operator for the Helmholtz equation (see, e.g., [43, 48, 47] and Section 2.7 below), and the unknown $\phi$ is the Neumann trace (i.e., the normal derivative on $\Gamma$) of $U$. The unknown $\phi$ lies in the Sobolev space $H^{-1/2}(\Gamma)$, and (1.2) is understood in the corresponding dual space $H^{1/2}(\Gamma)$, leading to a convenient variational framework (see Section 2.4 for a precise definition of $H^{\pm 1/2}(\Gamma)$). Once $\phi$ (or an approximation thereof) is computed, the wave field may be evaluated in $\Omega^{\text{ext}}$ through the representation formula

$$U = \tilde{K}_k u - \tilde{V}_k \phi, \qquad (1.3)$$

where $\tilde{V}_k$ and $\tilde{K}_k$ are the single- and double-layer potential of the Helmholtz equation (see again [43, 48, 47] and Section 2.7).

Similarly, an indirect approach is given by the ansatz $U = \tilde{V}_k \phi$ for some $\phi$ in $H^{-1/2}(\Gamma)$ (which is in general different from $\partial_\nu U$) with corresponding boundary integral equation

$$V_k \phi = u. \qquad (1.4)$$

Both equations (1.2) and (1.4) involve the single-layer operator $V_k$, and a corresponding adaptive BEM is formulated and analyzed in [4].

### 1.3. Spurious resonances.
However, a significant issue arises when using the first-kind boundary integral equations (1.2) and (1.4). Indeed, although the original Helmholtz problem in (1.1) is well-posed for arbitrary frequencies $k > 0$, the operator $V_k$ fails to be invertible if $k$ is a resonant frequency of $\Omega$ (i.e., if $k^2$ is a Dirichlet eigenvalue of the Laplace operator on $\Omega$). This phenomenon is called spurious resonance in the literature. Intuitively, it comes from the fact that (1.2) and (1.4) simultaneously encompass the exterior and the interior Dirichlet problem (with corresponding integral equation $V_k \phi = (K_k + 1/2)u$ or $V_k \phi = u$). While the former are always well-posed, the latter are not for resonant frequencies. Indeed, the (finite-dimensional) kernel of $V_k$ is spanned by the Neumann traces of Dirichlet eigenfunctions with eigenvalue $k^2$ in $\Omega$; see, e.g., [47, Theorem 3.9.1].

This leads to many issues in the numerical discretization at and around spurious resonances, such as ill-conditioning of the resulting matrices, slow convergence of iterative linear solvers, and loss of accuracy of the resulting approximation, cf. the numerical experiments of, e.g., [21, 45]. Closer to the present topic, our numerical experiments in Section 6 demonstrate that the convergence of adaptive BEM is also severely impaired in the vicinity of spurious resonances. These issues arise both for the direct formulation (1.2), but also for the corresponding indirect formulation (1.4).
### 1.4. Combined field integral equations.
These shortcomings have motivated the development of alternative integral equations to be used in place of (1.2) and (1.4). The seminal works [9, 41, 46] simultaneously proposed the indirect ansatz $U = \tilde{V}_k\phi + i\mu\tilde{K}_k\phi$ with real user-defined parameter $\mu \neq 0$ (which typically scales as $k^{-1}$) and corresponding integral equation

$$V_k\phi + i\mu(K_k + 1/2)\phi = u. \qquad (1.5)$$

The key idea is that an impedance boundary condition appears for the interior problem, which is thus always solvable. Similarly, [13] proposed the direct approach

$$V_k\phi + i\mu(K'_k + 1/2)\phi = (K_k - 1/2)u - i\mu W_k u, \qquad (1.6)$$

where $K'_k$ is the adjoint double-layer operator, $W_k$ the hypersingular operator, and the unknown $\phi$ is again the Neumann trace of $U$. It results from exploiting the representation formula (1.3) for the sum $U|_{\Gamma} + i\mu\partial_{\nu}U$. The equations (1.5) and (1.6) are typically called combined field integral equations. These equations are well-posed in the sense that $V_k + i\mu(K_k + 1/2)$ and $V_k + i\mu(K'_k + 1/2)$ are automorphisms on $L^2(\Gamma)$. While this seems not to be widely appreciated in the literature, this is indeed true for general Lipschitz domains; see [18, Section 2]. However, only if the boundary $\Gamma$ is smooth, the operators $K_k, K'_k$ are compact so that the combined integral operators are indeed coercive, which allows for a feasible numerical discretization.

In this work, we therefore consider regularized combined field integral equations, which were already introduced in [46] (for theoretical purposes though). This regularization consists in applying a regularizing operator $M: H^{-1/2}(\Gamma) \to H^{1/2}(\Gamma)$ to restore the variational setting of (1.2) and (1.4). In other words, we consider the integral equations

$$V_k\phi + iM(K_k + 1/2)\phi = u \qquad (1.7)$$

and

$$V_k\phi + iM(K'_k + 1/2)\phi = (K_k - 1/2)u - iMW_k u, \qquad (1.8)$$

where the multiplication by $\mu$ in (1.5) and (1.6) is replaced by an application of $M$

Different numerically feasible choices of $M$ have been suggested in the literature [11, 12, 21, 22, 38]. All of them guarantee that (1.7) and (1.8) are well-posed in the sense that the operators $V_k + iM(K_k + 1/2)$ and $V_k + i(K'_k + 1/2)M$ are coercive isomorphisms between $H^{-1/2}(\Gamma)$ and $H^{1/2}(\Gamma)$ for all $k > 0$.

### 1.5. Contributions of the manuscript.
Our key contribution is to extend algorithms and analysis proposed in [4] for adaptive BEM of the standard integral equations (1.2) and (1.4) to the regularized combined field integral equations (1.7) and (1.8). Specifically, given an initial triangulation $\mathcal{T}_0$, we discretize the latter equations with discontinuous piecewise polynomials of degree $p \in \mathbb{N}_0$ on a sequence of iteratively refined meshes $\mathcal{T}_\ell$.

For every possible refinement $\mathcal{T}_\bullet$ of $\mathcal{T}_0$, we introduce an a posteriori error estimator $\eta_\bullet$ which is reliable and weakly efficient, which is the state of the art for (1.2) and (1.4). This result is stated in Theorem 3.2. We then combine this error estimator with the usual Dörfler marking to an adaptive algorithm that produces a sequence of refined meshes $\mathcal{T}_\ell$ and associated approximations $\phi_\ell$. We show that, if the initial triangulation is sufficiently fine, then the sequence $\phi_\ell$ converges at optimal rate with respect to the number of mesh elements $\#\mathcal{T}_\ell$. This is stated precisely in Theorem 3.6.
The results mentioned above apply under mild assumptions on $M$, which are stated precisely in Section 2.11. In particular, the assumptions are satisfied if $M$ is a positive operator that continuously maps $H^{-1/2}(\Gamma)$ into $H^1(\Gamma)$. Important examples, also considered in [11, 12], that satisfy this assumption are the square of the single-layer operator of the Laplace operator $V_0^2$, the single-layer operator for the bi-Laplacian $\Delta^2$, and the inverse of $\alpha - \Delta_\Gamma$ with $\alpha > 0$ and $\Delta_\Gamma$ the Laplace-Beltrami operator.

Our first set of results is established under the assumption that the operator $M$ is computed exactly. Although this might be a reasonable assumption for the boundary integral operators mentioned above (up to quadrature errors), this is not the case for $M = (\alpha - \Delta_\Gamma)^{-1}$. For this reason, if $M = (\alpha - \Delta_\Gamma)^{-1}$, we introduce as in [11] an additional variable $f := M\phi$ for (1.7) and $f := M[(K'_k + 1/2)\phi + W_k u]$ for (1.8), and consider the resulting mixed formulation. The variable $\phi$ is discretized as above and $f$ is approximated by $H^1(\Gamma)$-conforming piecewise polynomials of degree $p+2$.

Our second set of results therefore concerns these mixed formulations. Specifically, in this case too, we propose a reliable and weakly efficient a posteriori error estimator that leads to an optimally convergent adaptive BEM algorithm. The corresponding results are reported in Theorems 4.4 and 4.7.

We finally present a set of numerical examples where the adaptive BEM based on the regularized mixed formulations is compared to the standard approach based on (1.2) and (1.4). These examples illustrate our theory by highlighting that the proposed estimators for the regularized mixed formulations are reliable and (weakly) efficient and that the adaptive BEM converges at optimal rate. In particular, the quality of the estimator and convergence speed of the adaptive algorithm are not affected by nearby interior resonances. In contrast, our examples showcase that this is not the case for the standard formulations (1.2) and (1.4).

An important comment is that although we do not precisely track the dependency on $k$ in our analysis, all the constants involved remain uniformly bounded whenever $k$ lies in a bounded interval in $(0, \infty)$. This means that although the quality of the estimator and performance of the adaptive algorithm may deteriorate for very small or very large values of $k$ (which is to be expected), it is indeed not affected by spurious resonances.

We finally point out that in this work, we only focus on the so-called sound-soft scattering problems where the boundary condition on $\Gamma$ is of Dirichlet type. Indeed, sound-hard scattering problems involving Neumann boundary conditions are easier to deal with in the context of combined field integral equations, as they do not need any regularization; see, e.g., [38]. As a result, a posteriori error estimation and optimality of the corresponding adaptive BEM follows as for standard boundary integral equations [4].

**1.6. Outline.** The remainder of this work is organized as follows. In Section 2, we make the setting precise, set up notation and collect preliminary results from the literature. Section 3 deals with the indirect formulation (1.7) where the operator $M$ is assumed to be exactly available. We then analyze in Section 4, the mixed formulation of (1.7) where the operator $M = (\alpha - \Delta_\Gamma)^{-1}$ is approximated. Section 5 then briefly discusses how the results established for the indirect formulations can be transferred to their direct counterparts. Finally, we present a set of numerical examples in Section 6.
## 2. PRELIMINARIES

**2.1. General notation.** Throughout and without any ambiguity, $|\cdot|$ denotes the absolute value of scalars, the Euclidean norm of vectors in $\mathbb{R}^n$, or the measure of a set in $\mathbb{R}^n$, e.g., the length of an interval or the area of a surface. The cardinality of a finite set $S$ is denoted by $|S|$. The set of all linear and continuous operators between two normed spaces $X$ and $Y$ is denoted by $\mathcal{L}(X, Y)$.

**2.2. Constants.** In the remainder of this work, we employ the notation $A \lesssim B$ for real numbers $A, B \ge 0$ to mean that there exists a constant $C(d, k, \Gamma, M, \mathcal{T}_0, p)$ depending at most on the dimension $d$, the wavenumber $k$, the boundary $\Gamma$, the regularizing operator $M$, the initial triangulation $\mathcal{T}_0$ of $\Gamma$, and the polynomial degree $p$ of the trial functions such that $A \le C(d, k, \Gamma, M, \mathcal{T}_0, p)B$. Throughout the manuscript, all constants $C$ are additionally uniformly bounded in $k$ on bounded intervals. In particular, this applies to the constants $C_{\text{rel}}$, $C_{\text{eff}}$, $C_{\text{lin}}$, $C_{\text{opt}}$ as well as $\rho_{\text{lin}}$, $\theta_{\text{opt}}$ from Theorems 3.2, 3.6, 4.4 and 4.7. Specifically, for all $K > 1$, the estimate $A \lesssim B$ holds true uniformly for all $k \in [1/K, K]$. If both $A \lesssim B$ and $B \lesssim A$, we also write $A \asymp B$.

**2.3. Domain and boundary.** Throughout this work, we consider a bounded Lipschitz domain $\Omega \subset \mathbb{R}^d$, $d \in \{2, 3\}$; see e.g., [43, Definition 3.28] for a precise definition. We assume that the complement $\Omega^{\text{ext}} := \mathbb{R}^d \setminus \bar{\Omega}$ is connected, and we abbreviate the joint boundary $\Gamma := \partial\Omega = \partial\Omega^{\text{ext}}$.

**2.4. Sobolev spaces on the boundary.** For $\sigma \in [0, 1]$, we define the Hilbert spaces $H^{\pm\sigma}(\Gamma)$ as in [43, page 99] by use of Bessel potentials on $\mathbb{R}^{d-1}$ and liftings via bi-Lipschitz mappings that describe $\Gamma$. For $\sigma = 0$, it holds that $H^0(\Gamma) = L^2(\Gamma)$ with equivalent norms. We thus may define $\|\cdot\|_{H^0(\Gamma)} := \|\cdot\|_{L^2(\Gamma)}$.

For $\sigma \in (0, 1]$, any measurable subset $\omega \subseteq \Gamma$, and all $v \in H^{\sigma}(\Gamma)$, we define the associated Sobolev-Slobodeckij norm

$$ \|v\|_{H^{\sigma}(\omega)}^{2} := \|v\|_{L^{2}(\omega)}^{2} + |v|_{H^{\sigma}(\omega)}^{2} $$

with

$$ |v|_{H^{\sigma}(\omega)}^{2} := \int_{\omega} \int_{\omega} \frac{|v(x) - v(y)|^{2}}{|x - y|^{d-1+2\sigma}} dxdy $$

if $\sigma \in (0, 1)$, and

$$ |v|_{H^1(\omega)} := \|\nabla_{\Gamma} v\|_{L^2(\omega)}. $$

It is well-known that $\|\cdot\|_{H^{\sigma}(\Gamma)}$ provides an equivalent norm on $H^{\sigma}(\Gamma)$. Here, $\nabla_{\Gamma}(\cdot)$ denotes the (weak) surface gradient.

For $\sigma \in (0, 1]$, $H^{-\sigma}(\Gamma)$ is a realization of the dual space of $H^{\sigma}(\Gamma)$. With the extended $L^2$-scalar product $\langle \cdot, \cdot \rangle_{L^2(\Gamma)}$, we define an equivalent norm

$$ \|\phi\|_{H^{-\sigma}(\Gamma)} := \sup_{v \in H^{\sigma}(\Gamma) \setminus \{0\}} \frac{|\langle \phi, v \rangle_{L^2(\Gamma)}|}{\|v\|_{H^{\sigma}(\Gamma)}} \quad \text{for all } \phi \in H^{-\sigma}(\Gamma). $$
2.5. **Model problem.** In the remainder of this work, we fix a wavenumber $k > 0$ as well as a Dirichlet datum

$$ u \in H^1(\Gamma). \qquad (2.1) $$

Our model problem is then to find $U : \Omega^{\text{ext}} \to \mathbb{C}$ such that

$$ \left\{ \begin{array}{ll} \Delta U + k^2 U = 0 & \text{in } \Omega^{\text{ext}}, \\ U = u & \text{on } \Gamma, \\ \partial_r U - ikU = o(r^{(1-d)/2}) & \text{as } r \to +\infty, \end{array} \right. \qquad (2.2) $$

where $r$ denotes the radial coordinate, and $\partial_r$ the radial derivative. While the weak formulation of (2.2) only requires the assumption $u \in H^{1/2}(\Gamma)$, the additional regularity (2.1) is needed for the definition of the considered a *posteriori* error estimators. This mild assumption is satisfied in many applications (e.g., plane wave scattering) and does not preclude $U$ to be highly singular in the vicinity of corners and edges (if $d = 3$) of $\Gamma$.

As alluded to above, instead of directly computing $U$, we will reformulate (2.2) into a boundary integral equation on $\Gamma$.

2.6. **Trace operators.** We denote by $\gamma_0^{\text{int}} : H^1(\Omega) \to H^{1/2}(\Gamma)$ the interior trace operator, which coincides for smooth functions with the usual restriction $(\cdot)|_{\Gamma}$. The exterior trace operator $\tilde{\gamma}_0^{\text{ext}} : H^1(\tilde{\Omega} \setminus \bar{\Omega}) \to H^{1/2}(\Gamma)$ is defined analogously for any bounded Lipschitz domain $\tilde{\Omega} \supset \bar{\Omega}$ as restriction to $\Gamma$.

With $H^1(\Delta; \Omega) := \{u \in H^1(\Omega) : \Delta u \in L^2(\Omega)\}$, we denote by $\gamma_1^{\text{int}} : H^1(\Delta, \Omega) \to H^{-1/2}(\Gamma)$ the normal derivative, which coincides for smooth functions with the usual normal derivative $\partial_\nu$. The exterior normal derivative $\tilde{\gamma}_1^{\text{ext}} : H^1(\Delta; \tilde{\Omega} \setminus \bar{\Omega}) \to H^{-1/2}(\Gamma)$ is defined analogously for any bounded Lipschitz domain $\tilde{\Omega} \supset \bar{\Omega}$ as restriction to $\Gamma$.

If the considered domain is clear from the context, we will simply write $(\cdot)|_{\Gamma}$ and $\partial_\nu$.

2.7. **Layer potentials and boundary integral operators.** For $k > 0$, we recall the Helmholtz kernel

$$ G_k(z) := \frac{i}{4} H_0^{(1)}(k|z|) \quad \text{for } d=2, \qquad G_k(z) := \frac{e^{ik|z|}}{4\pi|z|} \quad \text{for } d=3, $$

where $H_0^{(1)}$ is the first-kind Hankel function of order zero. We further use the fundamental solution of the Laplace operator

$$ G_0(z) := -\frac{1}{2\pi} \log|z| \quad \text{for } d=2, \qquad G_0(z) := \frac{1}{4\pi|z|} \quad \text{for } d=3. $$

For $k \ge 0$ and sufficiently smooth $\phi, u : \Gamma \to \mathbb{C}$, we define the single-layer potential

$$ (\tilde{V}_k \phi)(x) := \int_{\Gamma} G_k(x-y) \phi(y) dy \quad \text{for all } x \in \mathbb{R}^d \setminus \Gamma $$

and the double-layer potential

$$ (\tilde{K}_k u)(x) := \int_{\Gamma} \partial_{\nu(y)} G_k(x-y) u(y) dy \quad \text{for all } x \in \mathbb{R}^d \setminus \Gamma. $$

For any bounded Lipschitz domain $\tilde{\Omega}$, they give rise to bounded linear operators

$$ \tilde{V}_k \in \mathcal{L}(H^{-1/2}(\Gamma), H^1(\tilde{\Omega})) \quad \text{and} \quad \tilde{K}_k \in \mathcal{L}(H^{1/2}(\Gamma), H^1(\tilde{\Omega} \setminus \Gamma)). $$
For $\phi \in H^{-1/2}(\Gamma)$ and $v \in H^{1/2}(\Gamma)$, the potentials $\tilde{V}_k\phi$ and $\tilde{K}_k v$ are both solutions of the exterior Helmholtz problem (2.2) (with suitable Dirichlet data). The exact solution of (2.2) can be represented as

$$ U = \tilde{K}_k(U|_{\Gamma}) - \tilde{V}_k(\partial_{\nu}U). \qquad (2.3) $$

Applying the trace operators defines the single-layer operator

$$ V_k := \gamma_0^{\text{int}} \tilde{V}_k = \gamma_0^{\text{ext}} \tilde{V}_k \in \mathcal{L}(H^{-1/2}(\Gamma), H^{1/2}(\Gamma)), $$

the double-layer operator

$$ K_k := \gamma_0^{\text{int}} \tilde{K}_k + 1/2 = \gamma_0^{\text{ext}} \tilde{K}_k - 1/2 \in \mathcal{L}(H^{1/2}(\Gamma), H^{1/2}(\Gamma)), $$

the adjoint double-layer operator

$$ K'_k := \gamma_1^{\text{int}} \tilde{V}_k - 1/2 = \gamma_1^{\text{int}} \tilde{V}_k + 1/2 \in \mathcal{L}(H^{-1/2}(\Gamma), H^{-1/2}(\Gamma)), $$

and the hypersingular operator

$$ W_k := -\gamma_1^{\text{int}} \tilde{K}_k \in \mathcal{L}(H^{1/2}(\Gamma), H^{-1/2}(\Gamma)). $$

For $k = 0$, the additional mapping properties

$$ V_k \in \mathcal{L}(H^{-1/2+\sigma}(\Gamma) \rightarrow H^{1/2+\sigma}(\Gamma)), \qquad (2.4a) $$

$$ K_k \in \mathcal{L}(H^{1/2+\sigma}(\Gamma), H^{1/2+\sigma}(\Gamma)) \qquad (2.4b) $$

$$ K'_k \in \mathcal{L}(H^{-1/2+\sigma}(\Gamma), H^{-1/2+\sigma}(\Gamma)), \qquad (2.4c) $$

$$ W_k \in \mathcal{L}(H^{1/2+\sigma}(\Gamma), H^{-1/2+\sigma}(\Gamma)) \qquad (2.4d) $$

are well-known for all $\sigma \in [-1/2, 1/2]$. In [4, 32], these mapping properties are also verified for all $k > 0$ (and the critical cases $\sigma \in \{-1/2, 1/2\}$).

The single-layer operator $V_0 : H^{-1/2}(\Gamma) \rightarrow H^{1/2}(\Gamma)$ is symmetric and elliptic, provided that $\text{diam}(\Omega) < 1$ if $d = 2$. If $d = 2$ and $\text{diam}(\Omega) \ge 1$, this is still true up to a compact perturbation, i.e., there is a symmetric and elliptic $V_0^+ : H^{-1/2}(\Gamma) \rightarrow H^{1/2}(\Gamma)$ such that $V_0 - V_0^+$ is compact. If $d = 3$ or $d = 2$ with $\text{diam}(\Omega) < 1$, we abbreviate $V_0^+ := V_0$.

The operator $V_k : H^{-1/2}(\Gamma) \rightarrow H^{1/2}(\Gamma)$ is in general only elliptic up to the compact perturbation $V_0^+ - V_k$ with kernel given by the Neumann traces of Dirichlet eigenfunctions.

For a more detailed introduction to layer potentials and boundary integral operators, we refer to the monographs [43, 48, 47].

## 2.8. Meshes and refinement. Let $T_{\text{ref}}$ denote the reference element defined by

$$ T_{\text{ref}} := (0, 1) \text{ for } d = 2 \text{ and } T_{\text{ref}} := \text{conv}\{(0, 0), (1, 0), (0, 1)\} \text{ for } d = 3. $$

We consider conforming triangulations $\mathcal{T}_\bullet$ of $\Gamma$ characterized similarly as in [1] by the following properties:

* Each $T \in \mathcal{T}_\bullet$ is a relative open subset of $\Gamma$, and there exists a bi-Lipschitz parametrization $\gamma_T : \bar{T}_{\text{ref}} \rightarrow \bar{T}$ such that $\gamma_T \in C^2(T_{\text{ref}}, \mathbb{R}^d)$.

* The union of all elements cover $\Gamma$, i.e., $\Gamma = \bigcup \bar{T}$.

* There are no hanging nodes, i.e., for all $T \ne T' \in \mathcal{T}_\bullet$, the intersection $\bar{T} \cap \bar{T}'$ is either empty, a joint node (for $d \ge 2$), or a joint face (for $d = 3$).
• For $d = 3$, parametrizations of common boundary parts of neighboring elements are compatible, i.e., if $\emptyset \neq \bar{T} \cap \bar{T}' = \gamma_T(E_{\text{ref}}) = \gamma_{T'}(E'_{\text{ref}})$ is a joint edge, then $\gamma_T^{-1} \circ \gamma_{T'} : E'_{\text{ref}} \to E_{\text{ref}}$ is affine.

We define the mesh-size function $h_{\bullet} \in L^{\infty}(\Gamma)$ by $h_{\bullet}|_{T} := h_{T} := |T|^{1/(d-1)}$ for all $T \in \mathcal{T}_{\bullet}$. Moreover, let $\mathcal{N}_{\bullet}$ denote the set of all nodes in $\mathcal{T}_{\bullet}$ with corresponding node patches $\omega_{\bullet}(z) := \bigcup \{T \in \mathcal{T}_{\bullet} : z \in \bar{T}\}$ for $z \in \mathcal{N}_{\bullet}$.

Let $\mathcal{T}_0$ be a conforming initial triangulation meeting the above requirements. (Later, we shall need $\mathcal{T}_0$ to be sufficiently fine; cf. Sections 3.1 and 4.1.) As mesh refinement strategy, we employ the extended 1D bisection from [2] for $d = 2$ and newest vertex bisection from, e.g., [50], for $d = 3$. Given a conforming triangulation $\mathcal{T}_{\bullet}$ and a set of marked elements $\mathcal{M}_{\bullet} \subseteq \mathcal{T}_{\bullet}$, the call $\mathcal{T}_{\circ} = \text{refine}(\mathcal{T}_{\bullet}, \mathcal{M}_{\bullet})$ returns the coarsest refinement $\mathcal{T}_{\circ}$ of $\mathcal{T}_{\bullet}$ such that all $T \in \mathcal{M}_{\bullet}$ have been refined. We define the set $\text{refine}(\mathcal{T}_{\bullet})$ as the set of all conforming $\mathcal{T}_{\circ}$ that can be obtained from $\mathcal{T}_{\bullet}$ by a finite number of refinement steps. In particular, $\mathcal{T}_{\bullet}$ itself is always in the set of its refinements, i.e., $\mathcal{T}_{\bullet} \in \text{refine}(\mathcal{T}_{\bullet})$. For the initial mesh, we call $\mathcal{T} := \text{refine}(\mathcal{T}_0)$ the set of all admissible triangulations. The employed refinement strategies guarantee shape-regularity

$$h_T \approx \operatorname{diam}(T) \quad \text{for all } T \in \mathcal{T}_{\bullet} \in \mathcal{T}$$

and local quasi-uniformity

$$h_T \approx h_{T'} \quad \text{for all } T, T' \in \mathcal{T}_{\bullet} \in \mathcal{T} \text{ with } \bar{T} \cap \bar{T}' \neq \emptyset,$$

with hidden constants depending only on the initial triangulation $\mathcal{T}_0$. Moreover, for sequences of successively refined triangulations $(\mathcal{T}_\ell)_{\ell \in \mathbb{N}_0}$ with $\mathcal{T}_\ell = \text{refine}(\mathcal{T}_{\ell-1}, \mathcal{M}_{\ell-1})$ for some $\mathcal{M}_{\ell-1} \subseteq \mathcal{T}_{\ell-1}$, and $\mathcal{T}_\ell \in \mathcal{T}$, there hold for all $\ell \in \mathbb{N}$ the following three crucial properties:

(R1) Child estimate: $\#\mathcal{T}_\ell \le C_{\text{child}} \#\mathcal{T}_{\ell-1}$ with $C_{\text{child}} = 2$ for $d=2$ and $C_{\text{child}} = 4$ for $d=3$;

(R2) Overlay estimate: $\exists (\mathcal{T}_\ell \oplus \mathcal{T}_\bullet) \in \text{refine}(\mathcal{T}_\ell) \cap \text{refine}(\mathcal{T}_\bullet)$ s.t. $(\mathcal{T}_\ell \oplus \mathcal{T}_\bullet) \le \#\mathcal{T}_\ell + \#\mathcal{T}_\bullet - \#\mathcal{T}_0$;

(R3) Closure estimate: $\#\mathcal{T}_\ell - \#\mathcal{T}_0 \le C_{\text{clo}} \sum_{j=0}^{\ell-1} \#\mathcal{M}_j$ with $C_{\text{clo}} \ge 1$ depending only on $\mathcal{T}_0$.

The child estimate (R1) follows from the fact that in one refinement step, each element is refined at most into 2 and 4 elements, respectively, see, e.g., [2] for $d=2$ and [50] for $d=3$. A triangulation $\mathcal{T}_\ell \oplus \mathcal{T}_\bullet$ with (R2) can indeed be chosen as the literal overlay of two meshes, see, e.g., [2] for $d=2$ and [49] for $d=3$. The closure estimate (R3) is proved in [2] for $d=2$ and in [6, 50] for $d=3$. The latter two references actually require an admissibility condition of $\mathcal{T}_0$, which was shown to be redundant in [40].

We only mention that the results of the present work also apply for other types of meshes, for instance curvilinear quadrangulations; cf. [32, 33, 10] for details in case of the standard integral equation in (1.2).

2.9. Discrete spaces. Let $\mathcal{T}_{\bullet} \in \mathcal{T}$ be a given triangulation. For any polynomial degree $q \in \mathbb{N}_0$, we define the space of all (discontinuous) $\mathcal{T}_{\bullet}$-piecewise polynomials as

$$\mathcal{P}^q(\mathcal{T}_{\bullet}) := \{\psi_h \in L^\infty(\Gamma) : \psi_h|_T \circ \gamma_T \text{ is a polynomial of degree } \le q \text{ for all } T \in \mathcal{T}_{\bullet}\}.$$

For $q \in \mathbb{N}$, we define the space of continuous $\mathcal{T}_{\bullet}$-piecewise polynomials as $\mathcal{S}^q(\mathcal{T}_{\bullet}) := \mathcal{P}^q(\mathcal{T}_{\bullet}) \cap H^1(\Gamma)$. Note that

$$\mathcal{P}^q(\mathcal{T}_{\bullet}) \subset L^2(\Gamma) \subset H^{-1/2}(\Gamma) \text{ for } q \in \mathbb{N}_0 \text{ and } \mathcal{S}^q(\mathcal{T}_{\bullet}) \subset H^1(\Gamma) \subset H^{1/2}(\Gamma) \text{ for } q \in \mathbb{N}.$$
For the remainder of this work, we fix a polynomial degree $p \in \mathbb{N}_0$. We will consider the space $\mathcal{P}^p(\mathcal{T}_\bullet)$ and the space $\mathcal{P}^p(\mathcal{T}_\bullet) \times \mathcal{S}^{p+2}(\mathcal{T}_\bullet)$ as trial space for the discretization of the considered combined field integral equation and its formulation as mixed system, respectively.

We only mention that the results of the present work also apply for splines on quadrangulations; cf. [32, 33, 10] for details in case of standard integral equations (1.2).

2.10. **Inverse estimates.** The following crucial inverse estimates were proved in [1] for $k=0$ and generalized in [4, 32] to arbitrary $k \ge 0$:

$$ \|h_{\bullet}^{1/2} \nabla_{\Gamma} V_k \psi\|_{L^2(\Gamma)} \lesssim \|\psi\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2} \psi\|_{L^2(\Gamma)}, \quad (2.5a) $$

$$ \|h_{\bullet}^{1/2} \nabla_{\Gamma} K_k v\|_{L^2(\Gamma)} \lesssim \|v\|_{H^{1/2}(\Gamma)} + \|h_{\bullet}^{1/2} \nabla_{\Gamma} v\|_{L^2(\Gamma)} \quad (2.5b) $$

$$ \|h_{\bullet}^{1/2} K'_{k} \psi\|_{L^2(\Gamma)} \lesssim \|\psi\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2} \psi\|_{L^2(\Gamma)}, \quad (2.5c) $$

$$ \|h_{\bullet}^{1/2} W_k v\|_{L^2(\Gamma)} \lesssim \|v\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2} v\|_{L^2(\Gamma)} \quad (2.5d) $$

for all triangulations $\mathcal{T}_\bullet \in \mathcal{T}$ and all function $\psi \in L^2(\Gamma)$ and $v \in H^1(\Gamma)$ with hidden constants that do not depend on the wavenumber $k$.

We further recall the discrete inverse estimates

$$ \|h_{\bullet}^{1/2} \psi_{\bullet}\|_{L^2(\Gamma)} \lesssim \|\psi_{\bullet}\|_{H^{-1/2}(\Gamma)} \quad \text{for all } \psi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\bullet}), \quad (2.6a) $$

$$ \|h_{\bullet}^{1/2} \nabla_{\Gamma} v_{\bullet}\|_{L^2(\Gamma)} \lesssim \|v_{\bullet}\|_{H^{1/2}(\Gamma)} \quad \text{for all } v_{\bullet} \in \mathcal{S}^{p+2}(\mathcal{T}_{\bullet}); \quad (2.6b) $$

with hidden constants that do not depend on the polynomial degree $p$; cf. [1, Lemma A.1] and [3, Proposition 5], see also the original works [19, 37].

2.11. **Regularizing operator.** Our regularized combined field integral equations involve a regularizing operator $M: H^{-1/2}(\Gamma) \to H^{1/2}(\Gamma)$. Throughout, we assume that

$$ M \in \mathcal{L}(H^{-1/2}(\Gamma), H^{1/2}(\Gamma)) \text{ is compact,} \quad (M1) $$

that it is positive in the sense that

$$ \mathrm{Re}\langle M\psi, \psi \rangle_{L^2(\Gamma)} > 0 \quad \text{for all } \psi \in H^{-1/2}(\Gamma) \setminus \{0\}, \quad (M2) $$

and that it satisfies the following inverse inequality

$$ \|h_{\bullet}^{1/2} \nabla_{\Gamma} M \psi\|_{L^2(\Gamma)} \lesssim \|\psi\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2} \psi\|_{L^2(\Gamma)} \quad \text{for all } \psi \in L^2(\Gamma), \quad (M3) $$

which is similar to the inverse inequalities in (2.5) for the standard integral operators. Note that (M1) and (M3) particularly imply that

$$ M \in \mathcal{L}(L^2(\Gamma), H^1(\Gamma)). \quad (2.7) $$

Conversely, we point out that (M1) and (M3) are both trivially satisfied provided that $M \in \mathcal{L}(H^{-1/2}(\Gamma), H^1(\Gamma))$, which is the case for all concrete examples of the operator $M$ considered in this work.

**Proposition 2.1.** If $M \in \mathcal{L}(H^{-1/2}(\Gamma), H^1(\Gamma))$, then (M1) and (M3) hold true.

*Proof.* The requirement in (M1) simply follows from the fact that the injection $H^1(\Gamma) \hookrightarrow H^{1/2}(\Gamma)$ is compact. For (M3), the inequality $h_{\bullet} \lesssim \mathrm{diam}(\Gamma)$ and continuity of $M$ show for all $\psi \in H^{-1/2}(\Gamma)$ that

$$ \|h_{\bullet}^{1/2} \nabla_{\Gamma} M \psi\|_{L^2(\Gamma)} \lesssim \|M \psi\|_{H^1(\Gamma)} \lesssim \|\psi\|_{H^{-1/2}(\Gamma)}, $$
and trivially,

$$ \|h_{\bullet}^{1/2} \nabla_{\Gamma} M \psi\|_{L^2(\Gamma)} \lesssim \|\psi\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2} \psi\|_{L^2(\Gamma)} $$

whenever $\psi \in L^2(\Gamma)$.

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

As suggested in [12], the operator $M$ might be realized with an integral operator. For instance, we can use the square of the single-layer operator $V_0^2$ (provided that $\mathrm{diam}(\Omega) < 1$ if $d = 2$) or the single-layer operator of the biharmonic operator $(\alpha - \Delta)^2$ with $\alpha > 0$. Both of them continuously map $H^{-1}(\Gamma)$ to $H^1(\Gamma)$, which implies (M1) and (M3), and also satisfy (M2).

Another interesting choice of $M$ is the inverse of a differential operator. In particular, we will consider, as in [11], $M := (\alpha - \Delta_{\Gamma})^{-1}$ with $\alpha > 0$. Here, $\Delta_{\Gamma} \in \mathcal{L}(H^1(\Gamma), H^{-1}(\Gamma))$ is the Laplace-Beltrami operator defined weakly as

$$ \langle \Delta_{\Gamma} v, w \rangle_{L^2(\Gamma)} := - \langle \nabla_{\Gamma} v, \nabla_{\Gamma} w \rangle_{L^2(\Gamma)} \quad \text{for all } v, w \in H^1(\Gamma). \quad (2.8) $$

Also this choice of $M$ continuously maps $H^{-1}(\Gamma)$ to $H^1(\Gamma)$. Moreover, there holds the ellipticity

$$ \mathrm{Re}\langle M\psi, \psi \rangle_{L^2(\Gamma)} \gtrsim \|\psi\|_{H^{-1}(\Gamma)} \quad \text{for all } \psi \in H^{-1}(\Gamma), $$

which guarantees positivity (M2). However, as already mentioned, such an $M$ is not computable (even up to quadrature errors), and it needs to be approximated. This leads to mixed formulations analyzed in Sections 4 and 5.2.

Similarly, one could consider, as in [11], the inverse of the piecewise Laplace-Beltrami operator for $M$, i.e., $M\psi \in H_0^1(\mathcal{T}_0)$ with

$$ \langle \nabla_{\Gamma} M \psi, \nabla_{\Gamma} v \rangle_{L^2(\Gamma)} = \langle \psi, v \rangle_{L^2(\Gamma)} \quad \text{for all } \psi \in H^{-1}(\Gamma), v \in H_0^1(\mathcal{T}_0), \quad (2.9) $$

where $H_0^1(\mathcal{T}_0) := \{v \in H^1(\Gamma) : v|_{T_0} \in H_0^1(T) \text{ for all } T_0 \in \mathcal{T}_0\}$. This allows to exploit the shift theorem on the regular surfaces $T_0 \in \mathcal{T}_0$, being diffeomorphic images of convex sets.

**Remark 2.2.** The works [38, 21, 22] consider a regularizing operator $M$ that involves the inverse of a stabilized hypersingular operator $\widetilde{W}_0$. For general Lipschitz domains, this operator is not compact. As in Sections 4 and 5.2, where we will consider $M := (\alpha - \Delta_{\Gamma})^{-1}$, this approach naturally leads to a mixed system. Similarly as we do for the Laplace-Beltrami equation, one can use a second weighted-residual estimator for the involved hypersingular integral equation. We expect that the analogous results of Sections 4 and 5.2 can be proved by combining the analysis of this work with the a posteriori and convergence analysis for the standard hypersingular integral equation [27].

### 3. INDIRECT INTEGRAL EQUATION

Given a regularizing operator $M \in \mathcal{L}(H^{-1/2}(\Gamma), H^{1/2}(\Gamma))$ satisfying (M1) and (M2), we follow [11] and make the ansatz

$$ U = \widetilde{V}_k \phi + i \widetilde{K}_k M \phi \quad (3.1) $$

for some suitable density function $\phi \in H^{-1/2}(\Gamma)$ to obtain a solution of the exterior Dirichlet problem (2.2). Applying the (exterior) trace operator to (3.1) and exploiting the Dirichlet condition of (2.2), we see that

$$ u = U|_{\Gamma} = V_k \phi + i(K_k + 1/2) M \phi. $$
With the short-hand notation

$$S := V_k + i(K_k + 1/2)M \in \mathcal{L}(H^{-1/2}(\Gamma), H^{1/2}(\Gamma)),$$

we hence need to find $\phi \in H^{-1/2}(\Gamma)$ such that

$$S\phi = u. \qquad (3.2)$$

Clearly, the operator $S$ is the sum of the symmetric and elliptic $V_0^+$ and the compact perturbation $V_k - V_0^+ + i(K_k + 1/2)M$. According to [11], $S$ is also injective and thus even bijective by the Fredholm alternative. In other words, the problem in (3.2) is well-posed for all $u \in H^{1/2}(\Gamma)$.

**3.1. Discretization.** Given any $\mathcal{T}_\bullet \in \mathbb{T}$, the discrete problem consists in finding a Galerkin approximation $\phi_\bullet \in \mathcal{P}^p(\mathcal{T}_\bullet)$ of $\phi$ satisfying the discrete variational formulation

$$\langle S\phi_\bullet, \psi_\bullet \rangle_{L^2(\Gamma)} = \langle u, \psi_\bullet \rangle_{L^2(\Gamma)} \quad \text{for all } \psi_\bullet \in \mathcal{P}^p(\mathcal{T}_\bullet). \qquad (3.3)$$

To theoretically guarantee well-posedness of (3.3), we assume that the initial triangulation $\mathcal{T}_0$ is sufficiently fine such that

$$\inf_{\chi_\bullet \in \mathcal{P}^p(\mathcal{T}_\bullet)} \sup_{\psi_\bullet \in \mathcal{P}^p(\mathcal{T}_\bullet)} \frac{|\langle S\chi_\bullet, \psi_\bullet \rangle_{L^2(\Gamma)}|}{\|\psi_\bullet\|_{H^{-1/2}(\Gamma)}} \gtrsim 1 \quad \text{for all } \mathcal{T}_\bullet \in \mathbb{T}, \qquad (3.4)$$

where we make the convention that $\frac{0}{0} := 0$; see also Remark 3.5 for a way (which is not fully satisfactory though) to avoid this assumption.

Owing to the fact that $S$ is symmetric and elliptic up to some compact perturbation, it is well known that there exists some $h \gtrsim 1$ such that (3.4) is satisfied whenever $h_T \le h$ for all $T \in \mathcal{T}_0$; see, e.g., [47, Theorem 4.2.9] or [5, Proposition 1]. The dependency on $h$, $k$, $p$, and $\Gamma$ is hard to track, but we refer the reader to [42, 30] for the case without regularization.

**Remark 3.1.** Given any $\mathcal{T}_\bullet \in \mathbb{T}$, the unique Galerkin approximation $\phi_\bullet \in \mathcal{P}^p(\mathcal{T}_\bullet)$ of $\phi$ obtained by solving (3.3) satisfies the Céa lemma

$$\|\phi - \phi_\bullet\|_{H^{-1/2}(\Gamma)} \lesssim \min_{\psi_\bullet \in \mathcal{P}^p(\mathcal{T}_\bullet)} \|\phi - \psi_\bullet\|_{H^{-1/2}(\Gamma)}. \qquad (3.5)$$

**3.2. A posteriori error estimation.** We employ the standard weighted-residual error estimator

$$\eta_\bullet := \|h_\bullet^{1/2} \nabla_\Gamma (u - S\phi_\bullet)\|_{L^2(\Gamma)}. \qquad (3.6)$$

An analogous estimator was introduced for the Laplace single-layer operator $V_0$ instead of $S$ in [16, 17, 14]; see also [4] for the Helmholtz single-layer operator $V_k$. We show that $\eta_\bullet$ is reliable and weakly efficient.

**Theorem 3.2.** Let $M$ be a regularizing operator with (M1), (M2), and (M3). Then, there exist constants $C_{\text{rel}}$, $C_{\text{eff}} > 0$ such that

$$C_{\text{rel}}^{-1} \|\phi - \phi_\bullet\|_{H^{-1/2}(\Gamma)} \le \eta_\bullet \le C_{\text{eff}} \|h_\bullet^{1/2}(\phi - \phi_\bullet)\|_{L^2(\Gamma)} \quad \text{for all } \mathcal{T}_\bullet \in \mathbb{T}. \qquad (3.7)$$

The constants $C_{\text{rel}}$ and $C_{\text{eff}}$ depend only on $d, k, \Gamma, M$, and $\mathcal{T}_0$.
**Remark 3.3.** We note that reliability, i.e., the first inequality of (3.7), does not require the discrete inf-sup stability (3.4) and is in fact true for any (in general non-unique) $\phi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\bullet})$ satisfying (3.3). Moreover, the proof of reliability does not exploit (M3) but rather its implication (2.7). Similarly, without (3.4), weak efficiency, i.e., the second inequality of (3.7), is still satisfied if the right-hand side is replaced by $\|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2}(\phi - \phi_{\bullet})\|_{L^2(\Gamma)}$.

*Proof of reliability.* Since $S$ is an isomorphism, we can readily estimate the discretization error by the norm of the residual, i.e.,

$$ \|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} \approx \|u - S\phi_{\bullet}\|_{H^{1/2}(\Gamma)} $$

for all $\mathcal{T}_{\bullet} \in \mathbb{T}$. While the latter term is in principle computable, it involves a double integral over the boundary $\Gamma$ and does not provide any local information about the error. The localization argument of [23, 24] shows that

$$ \sum_{z \in \mathcal{N}_{\bullet}} |v|_{H^{1/2}(\omega_{\bullet}(z))}^2 \lesssim \|v\|_{H^{1/2}(\Gamma)}^2 \lesssim \sum_{z \in \mathcal{N}_{\bullet}} |v|_{H^{1/2}(\omega_{\bullet}(z))}^2 + \|h_{\bullet}^{-1/2}v\|_{L^2(\Gamma)}^2 \quad (3.8) $$

for all $v \in H^{1/2}(\Gamma)$. If $v$ is orthogonal to the space of piecewise constants $\mathcal{P}^0(\mathcal{T}_{\bullet})$, a Poincaré-type inequality of [23, 24] gives that

$$ \|h_{\bullet}^{-1/2}v\|_{L^2(\Gamma)}^2 \lesssim \sum_{z \in \mathcal{N}_{\bullet}} |v|_{H^{1/2}(\omega_{\bullet}(z))}^2. $$

Choosing $v = u - S\phi_{\bullet}$, we hence see the equivalence

$$ \|u - S\phi_{\bullet}\|_{H^{1/2}(\Gamma)}^2 \approx \sum_{z \in \mathcal{N}_{\bullet}} |u - S\phi_{\bullet}|_{H^{1/2}(\omega_{\bullet}(z))}^2. $$

Owing to (2.1), the mapping properties of $V_k$ and $K_k$ stated in (2.4) with $\sigma = 1/2$, and (2.7), it holds that $u - S\phi_{\bullet} \in H^1(\Gamma)$. As a result, a scaling argument allows to replace the Sobolev–Slobodeckij seminorm on patches by the $H^1$-seminorm on elements. Indeed, it holds that

$$ |v|_{H^{1/2}(\omega_{\bullet}(z))} \lesssim \|h_{\bullet}^{1/2} \nabla_{\Gamma} v\|_{L^2(\omega_{\bullet}(z))} \quad \text{for all } v \in H^1(\Gamma) \text{ and } z \in \mathcal{N}_{\bullet}; \quad (3.9) $$

see, e.g., [31, Proposition 5.2.2] for a detailed proof. Owing to finite overlap of the patches, we get that

$$ \|u - S\phi_{\bullet}\|_{H^{1/2}(\Gamma)}^2 \lesssim \sum_{T \in \mathcal{T}_{\bullet}} \|h_{\bullet}^{1/2} \nabla_{\Gamma}(u - S\phi_{\bullet})\|_{L^2(T)}^2 = \eta_{\bullet}^2, $$

which proves the first inequality in (3.7). $\square$

*Proof of weak efficiency.* We proceed similarly as in [1, Section 3.2.1], where weak efficiency is shown for the Laplace problem. The definition of $S$ and the inverse estimates (2.5a), (2.5b), and (M3) show that

$$ \begin{aligned} \|h_{\bullet}^{1/2} \nabla_{\Gamma}(u - S\phi_{\bullet})\|_{L^2(\Gamma)} &= \|h_{\bullet}^{1/2} \nabla_{\Gamma}[V_k(\phi - \phi_{\bullet}) + i(K_k + 1/2)M(\phi - \phi_{\bullet})]\|_{L^2(\Gamma)} \\ &\lesssim \|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2}(\phi - \phi_{\bullet})\|_{L^2(\Gamma)} + \|M(\phi - \phi_{\bullet})\|_{H^{1/2}(\Gamma)} + \|h_{\bullet}^{1/2} \nabla_{\Gamma} M(\phi - \phi_{\bullet})\|_{L^2(\Gamma)} \\ &\lesssim \|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2}(\phi - \phi_{\bullet})\|_{L^2(\Gamma)}. \end{aligned} $$
Using that the Galerkin projection $\phi \mapsto \phi_{\bullet}$ is $H^{-1/2}(\Gamma)$-stable, one can now argue as in the proof of [1, Corollary 3.3] to see that

$$ \|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} \lesssim \|h_{\bullet}^{1/2}(\phi - \phi_{\bullet})\|_{L^2(\Gamma)}, $$

which concludes the proof. □

**3.3. Adaptive algorithm.** We abbreviate the weighted-residual error indicators by

$$ \eta_{\bullet}(T)^2 := \|h_{\bullet}^{1/2} \nabla_{\Gamma}(u - S\phi_{\bullet})\|_{L^2(T)}^2 \quad \text{for all } T \in \mathcal{T}_{\bullet}. \qquad (3.10) $$

Moreover, we set

$$ \eta_{\bullet}(\mathcal{U}_{\bullet})^2 := \sum_{T \in \mathcal{U}_{\bullet}} \eta_{\bullet}(T)^2 \quad \text{for all } \mathcal{U}_{\bullet} \subseteq \mathcal{T}_{\bullet}. $$

Note that $\eta_{\bullet} = \eta_{\bullet}(\mathcal{T}_{\bullet})$. We can use the weighted-residual error estimator within a standard adaptive algorithm to steer the local mesh refinement of some given sufficiently fine conforming initial triangulation $\mathcal{T}_0$ such that the discrete inf-sup stability (3.4) is satisfied.

**Algorithm 3.4.** Input: Dörfler parameter $0 < \theta \le 1$.

Loop: For each $\ell = 0, 1, 2, \dots$, iterate the following steps:

(i) Compute the Galerkin approximation $\phi_{\ell} \in \mathcal{P}^{p}(\mathcal{T}_{\ell})$ by solving (3.3).

(ii) Compute the refinement indicators $\eta_{\ell}(T)$ from (3.10) for all elements $T \in \mathcal{T}_{\ell}$.

(iii) Determine a minimal<fnref n="1" /> set of marked elements $\mathcal{M}_{\ell} \subseteq \mathcal{T}_{\ell}$ satisfying the Dörfler marking

$$ \theta \eta_{\ell}^2 \le \eta_{\ell}(\mathcal{M}_{\ell})^2. $$

(iv) Generate the refined triangulation $\mathcal{T}_{\ell+1} := \text{refine}(\mathcal{T}_{\ell}, \mathcal{M}_{\ell})$.

**Output:** Refined triangulations $\mathcal{T}_{\ell}$, corresponding Galerkin approximations $\phi_{\ell} \in \mathcal{P}^{p}(\mathcal{T}_{\ell})$, and weighted-residual error estimators $\eta_{\ell}$ for all $\ell \in \mathbb{N}_{0}$.

**Remark 3.5.** In [4], which investigates adaptive BEM for the standard integral equation (1.2), the assumption that $\mathcal{T}_0$ is sufficiently fine is circumvented by an alternative adaptive algorithm: If the discrete system on level $\ell$ is not solvable, then one uniform refinement is performed. Moreover, in each step an extended Dörfler marking is used, which additionally marks the largest elements in the current mesh. Hence, that algorithm eventually guarantees the inf-sup condition (3.4) for all $\mathcal{T}_{\bullet} \in \text{refine}(\mathcal{T}_{L})$ for a sufficiently large $L \in \mathbb{N}_{0}$. The analysis of the present manuscript also covers the analogous version of that algorithm. However, we stress that it is unclear in practice when a discrete system is truly non-solvable and all convergence results only hold starting from level $L$.

**3.4. Optimal convergence.** We set

$$ \mathcal{T}(N) := \{ \mathcal{T}_{\bullet} \in \mathcal{T} : \#\mathcal{T}_{\bullet} - \#\mathcal{T}_{0} \le N \} \quad \text{for all } N \in \mathbb{N}_{0} $$

and, for all $s > 0$,

$$ \|\phi\|_{\mathcal{A}_{s}^{\text{est}}} := \sup_{N \in \mathbb{N}_{0}} \min_{\mathcal{T}_{\bullet} \in \mathcal{T}(N)} (N+1)^s \eta_{\bullet} \in [0, \infty]. \qquad (3.11) $$

<footnote n="1">¹More precisely, we choose $\mathcal{M}_{\ell} \subseteq \mathcal{T}_{\ell}$ with $\theta \eta_{\ell}^2 \le \eta_{\ell}(\mathcal{U}_{\ell})^2$ such that $\#\mathcal{M}_{\ell} \le \#\mathcal{U}_{\ell}$ for all $\mathcal{U}_{\ell} \subseteq \mathcal{T}_{\ell}$ with $\theta \eta_{\ell}^2 \le \eta_{\ell}(\mathcal{U}_{\ell})^2$. Note that this choice might not be unique.</footnote>
We say that the solution $\phi \in H^{-1/2}(\Gamma)$ lies in the approximation class s with respect to the estimator if

$$ \|\phi\|_{A_s^{\text{est}}} < \infty. \qquad (3.12) $$

By definition, $\|\phi\|_{A_s^{\text{est}}} < \infty$ implies that the error estimator $\eta_{\bullet}$ on the optimal triangulations $\mathcal{T}_{\bullet}$ decays at least with rate $O((\#\mathcal{T}_{\bullet})^{-s})$. The following theorem states that each possible rate $s > 0$ is indeed be realized by Algorithm 3.4. The proof is given in Section 3.5, where we will verify the so-called axioms of adaptivity from [15] for the error estimator, which automatically imply the desired result.

**Theorem 3.6.** Suppose that $M$ satisfies (M1), (M2) and (M3), and that the discrete inf-sup stability (3.4) is satisfied. Let $(\mathcal{T}_\ell)_{\ell \in \mathbb{N}_0}$ be the sequence of triangulations generated by Algorithm 3.4. Then, for arbitrary $0 < \theta \le 1$, the estimator converges linearly, i.e., there exist constants $0 < \rho_{\text{lin}} < 1$ and $C_{\text{lin}} \ge 1$ such that

$$ \eta_{\ell+j}^2 \le C_{\text{lin}} \rho_{\text{lin}}^j \eta_{\ell}^2 \quad \text{for all } j, \ell \in \mathbb{N}_0. $$

Moreover, there exists a constant $0 < \theta_{\text{opt}} \le 1$ such that for all $0 < \theta < \theta_{\text{opt}}$, the estimator converges at optimal rate, i.e., for all $s > 0$, there exist constants $c_{\text{opt}}$, $C_{\text{opt}} > 0$ such that

$$ c_{\text{opt}} \|\phi\|_{A_s^{\text{est}}} \le \sup_{\ell \in \mathbb{N}_0} (\#\mathcal{T}_\ell - \#\mathcal{T}_0 + 1)^s \eta_\ell \le C_{\text{opt}} \|\phi\|_{A_s^{\text{est}}}. $$

The constant $\theta_{\text{opt}}$ depends only on $d, k, \Gamma, M, \mathcal{T}_0$, and $p$, while $C_{\text{lin}}$, $\rho_{\text{lin}}$ depend additionally on $\theta$, and $C_{\text{opt}}$ depends furthermore on $s > 0$. The constant $c_{\text{opt}}$ depends only on $d, \#\mathcal{T}_0$, $s$, and if there exists $\ell_0$ with $\eta_{\ell_0} = 0$, then also on $\ell_0$ and $\eta_0$.

**Remark 3.7.** We can slightly lighten the assumptions on the regularizing operator $M$ in Theorem 3.6. It suffices to suppose (M3) for $\psi \in \mathcal{P}^p(\mathcal{T}_{\bullet})$ instead of $\psi \in L^2(\Gamma)$. With the discrete inverse estimate (2.6a), (M3) can thus be replaced by

$$ \|h_{\bullet}^{1/2} \nabla_{\Gamma} M \psi_{\bullet}\|_{L^2(\Gamma)} \lesssim \|\psi_{\bullet}\|_{H^{-1/2}(\Gamma)} \quad \text{for all } \psi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \text{ and all } \mathcal{T}_{\bullet} \in \mathcal{T}. \qquad (3.13) $$

We additionally note that in fact, (3.13) is only required for the axioms (A1) and (A2) below. The other axioms (namely (A3) and (A4)) only require (M1) and (M2).

**3.5. Proof of Theorem 3.6.** As already announced, we only need to verify the axioms of adaptivity from [15]. Then, the abstract result [15, Theorem 4.1] immediately implies linear convergence (3.11) at optimal rate (3.12). We mention that the required axioms for the mesh refinement have already been stated in (R1)–(R3). So we only need to show that there exist uniform constants $C_{\text{stab}}$, $C_{\text{red}}$, $C_{\text{drel}}$, $C_{\text{qo}} > 0$ and $0 < \rho_{\text{red}} < 1$ such that the following axioms of adaptivity for the estimator hold for all $\mathcal{T}_{\bullet} \in \mathcal{T}$, $\mathcal{T}_{\circ} \in \text{refine}(\mathcal{T}_{\bullet})$ with corresponding patch of all refined elements $\mathcal{R}_{\bullet, \circ} := \{T \in \mathcal{T}_{\bullet} : \bar{T} \cap \bigcup(\mathcal{T}_{\bullet} \setminus \mathcal{T}_{\circ}) \neq \emptyset\}$, the sequence $(\mathcal{T}_{\ell})_{\ell \in \mathbb{N}_0}$ of Algorithm 3.4, and all $\ell, n \in \mathbb{N}_0$:

(A1) Stability on non-refined elements: $|\eta_{\circ}(\mathcal{T}_{\circ} \cap \mathcal{T}_{\bullet}) - \eta_{\bullet}(\mathcal{T}_{\bullet} \cap \mathcal{T}_{\circ})| \le C_{\text{stab}} \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)}$;

(A2) Reduction on refined elements: $\eta_{\circ}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\bullet}) \le \rho_{\text{red}} \eta_{\bullet}(\mathcal{T}_{\bullet} \setminus \mathcal{T}_{\circ}) + C_{\text{red}} \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)}$;

(A3) General quasi-orthogonality: $\sum_{j=\ell}^{\ell+n} \|\phi_{j+1} - \phi_j\|_{H^{-1/2}(\Gamma)}^2 \le C_{\text{qo}}^2 \eta_{\ell}^2$;

(A4) Discrete reliability: $\|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} \le C_{\text{drel}} \eta_{\bullet}(\mathcal{R}_{\bullet, \circ})$.
Proof of stability on non-refined elements (A1). Let $\mathcal{T}_{\bullet} \in \mathbb{T}$ and $\mathcal{T}_{\circ} \in \text{refine}(\mathcal{T}_{\bullet})$. The inverse triangle inequality and the fact that $h_{\bullet} = h_{\circ}$ on $\omega := \bigcup(\mathcal{T}_{\bullet} \cap \mathcal{T}_{\circ})$ show that

$$ \|h_{\circ}^{1/2} \nabla_{\Gamma}(u - S\phi_{\circ})\|_{L^2(\omega)} - \|h_{\bullet}^{1/2} \nabla_{\Gamma}(u - S\phi_{\bullet})\|_{L^2(\omega)} \le \|h_{\circ}^{1/2} \nabla_{\Gamma} S(\phi_{\circ} - \phi_{\bullet})\|_{L^2(\Gamma)}. $$

The inverse estimates (2.5a)–(2.5b), (2.6), and (3.13) on the mesh $\mathcal{T}_{\circ}$ imply that

$$ \begin{aligned} \|h_{\circ}^{1/2} \nabla_{\Gamma} S(\phi_{\circ} - \phi_{\bullet})\|_{L^2(\Gamma)} & \lesssim \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|M(\phi_{\circ} - \phi_{\bullet})\|_{H^{1/2}(\Gamma)} \\ & \quad + \|h_{\bullet}^{1/2} \nabla_{\Gamma} M(\phi_{\circ} - \phi_{\bullet})\|_{L^2(\Gamma)} \\ & \lesssim \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)}. \end{aligned} $$

This concludes the proof. $\square$

Proof of reduction on refined elements (A2). Let $\mathcal{T}_{\bullet} \in \mathbb{T}$ and $\mathcal{T}_{\circ} \in \text{refine}(\mathcal{T}_{\bullet})$. The triangle inequality and the fact that $h_{\circ}^{1/2} \le \rho h_{\bullet}^{1/2}$ on $\omega := \bigcup(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\bullet}) = \bigcup(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\circ})$ for some $0 < \rho < 1$ show that

$$ \|h_{\circ}^{1/2} \nabla_{\Gamma}(u - S\phi_{\circ})\|_{L^2(\omega)} \le \rho \|h_{\bullet}^{1/2} \nabla_{\Gamma}(u - S\phi_{\bullet})\|_{L^2(\omega)} + \|h_{\circ}^{1/2} \nabla_{\Gamma} S(\phi_{\circ} - \phi_{\bullet})\|_{L^2(\Gamma)}. $$

The second term has already been estimated in the proof of (A1). $\square$

Proof of general quasi-orthogonality (A3). By discrete inf-sup stability (3.4), [25] is applicable and directly implies this axiom with $C_{qo} = C_{qo}(n) = o(n)$ depending sublinearly on the integer $n$. The post-processing argument [25, Remark 20] or [15, Proposition 4.11] then yield the axiom with a uniform constant $C_{qo} \simeq 1$. $\square$

**Remark 3.8.** Since $S$ is symmetric and elliptic up to some compact perturbation, general quasi-orthogonality (A3) alternatively follows as in [28, Section 3–4] and [15, Proposition 4.11]; see also the abstract framework of [5, Section 3.3]. The drawback of this approach is that the resulting $C_{qo}$ depends on the sequence of Galerkin approximations $(\phi_{\ell})_{\ell \in \mathbb{N}_{0}}$.

Proof of discrete reliability (A4). The proof follows as for the weakly-singular integral equation of the Laplace operator [29, 36] or the standard weakly-singular integral equation for the Helmholtz operator [4], and is only given for the sake of completeness. Discrete inf-sup stability (3.4) and the definition of $\phi_{\circ}$ as well $\phi_{\bullet}$ (see (3.3)) show for all $\psi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\circ})$ that

$$ \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} \lesssim \sup_{\psi_{\circ} \in \mathcal{P}^p(\mathcal{T}_{\circ})} \frac{|\langle S(\phi_{\circ} - \phi_{\bullet}), \psi_{\circ} \rangle_{L^2(\Gamma)}|}{\|\psi_{\circ}\|_{H^{-1/2}(\Gamma)}} = \sup_{\psi_{\circ} \in \mathcal{P}^p(\mathcal{T}_{\circ})} \frac{|\langle S(\phi - \phi_{\bullet}), \psi_{\circ} - \psi_{\bullet} \rangle_{L^2(\Gamma)}|}{\|\psi_{\circ}\|_{H^{-1/2}(\Gamma)}}. $$

We choose

$$ \psi_{\bullet} := \begin{cases} 0 & \text{on } \bigcup(\mathcal{T}_{\bullet} \setminus \mathcal{T}_{\circ}), \\ \psi_{\circ} & \text{on } \bigcup(\mathcal{T}_{\bullet} \cap \mathcal{T}_{\circ}). \end{cases} \qquad (3.14) $$

In particular, the support of $\psi_{\circ} - \psi_{\bullet}$ is contained in $\omega := \bigcup(\mathcal{T}_{\bullet} \setminus \mathcal{T}_{\circ})$. With the hat functions $\varphi_{\bullet,z} \in \mathcal{S}^1(\mathcal{T}_{\bullet})$ defined for all $z \in \mathcal{N}_{\bullet}$ via $\varphi_{\bullet,z}(z') = \delta_{z,z'}$ for all $z' \in \mathcal{N}_{\bullet}$, the continuous cut-off function $\tilde{\chi}_{\omega} := \sum_{z \in \mathcal{N}_{\bullet} \cap \omega} \varphi_{\bullet,z}$ equals 1 on the set $\omega$. This shows that

$$ \langle S(\phi - \phi_{\bullet}), \psi_{\circ} - \psi_{\bullet} \rangle_{L^2(\Gamma)} = \langle \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet}), \psi_{\circ} - \psi_{\circ} \rangle_{L^2(\Gamma)}. $$

We bound $\langle \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet}), \psi_{\circ} \rangle_{L^2(\Gamma)}$ and $\langle \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet}), \psi_{\bullet} \rangle_{L^2(\Gamma)}$ separately. The Cauchy–Schwarz inequality, the bounds $0 \le \tilde{\chi}_{\omega} \le 1$ with $\text{supp}(\tilde{\chi}_{\omega}) \subseteq \bigcup \mathcal{R}_{\bullet,\circ}$, the fact that $h_{\bullet} = h_{\circ}$
and $\psi_{\bullet} = \psi_{\circ}$ on $\bigcup(\mathcal{T}_{\bullet} \cap \mathcal{T}_{\circ})$, the inverse estimate (2.6), and an element-wise Poincaré inequality (see, e.g., [31, Lemma 5.2.1 and 5.3.3]) show for the latter term that

$$ \begin{aligned} |\langle \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet}), \psi_{\bullet} \rangle_{L^2(\Gamma)}| &\le \|h_{\bullet}^{-1/2} \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet})\|_{L^2(\Gamma)} \|h_{\bullet}^{1/2} \psi_{\bullet}\|_{L^2(\bigcup(\mathcal{T}_{\bullet} \cap \mathcal{T}_{\circ}))} \\ &\le \|h_{\bullet}^{-1/2} S(\phi - \phi_{\bullet})\|_{L^2(\bigcup \mathcal{R}_{\bullet, \circ})} \|h_{\circ}^{1/2} \psi_{\circ}\|_{L^2(\Gamma)} \\ &\lesssim \|h_{\bullet}^{1/2} \nabla_{\Gamma} S(\phi - \phi_{\bullet})\|_{L^2(\bigcup \mathcal{R}_{\bullet, \circ})} \|\psi_{\circ}\|_{H^{-1/2}(\Gamma)}. \end{aligned} $$

For the other term, we have that

$$ |\langle \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet}), \psi_{\circ} \rangle_{L^2(\Gamma)}| \le \|\tilde{\chi}_{\omega} S(\phi - \phi_{\bullet})\|_{H^{1/2}(\Gamma)} \|\psi_{\circ}\|_{H^{-1/2}(\Gamma)}. $$

With the Faermann estimate (3.8), we see that

$$ \|\tilde{\chi}_{\omega} S(\phi - \phi_{\bullet})\|_{H^{1/2}(\Gamma)}^2 \lesssim \sum_{z \in \mathcal{N}_{\bullet}} |\tilde{\chi}_{\omega} S(\phi - \phi_{\bullet})|_{H^{1/2}(\omega_{\bullet}(z))}^2 + \|h_{\bullet}^{-1/2} \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet})\|_{L^2(\Gamma)}^2. $$

We have already estimated the $L^2$-term before. The local bound (3.9), the product rule, and the properties of $\tilde{\chi}_{\omega}$ show for the other term that

$$ \begin{aligned} \sum_{z \in \mathcal{N}_{\bullet}} |\tilde{\chi}_{\omega} S(\phi - \phi_{\bullet})|_{H^{1/2}(\omega_{\bullet}(z))}^2 &\lesssim \sum_{z \in \mathcal{N}_{\bullet}} \|h_{\bullet}^{1/2} \nabla_{\Gamma}(\tilde{\chi}_{\omega} S(\phi - \phi_{\bullet}))\|_{L^2(\omega_{\bullet}(z))}^2 \\ &\lesssim \|h_{\bullet}^{1/2} \nabla_{\Gamma} \tilde{\chi}_{\omega} S(\phi - \phi_{\bullet})\|_{L^2(\Gamma)}^2 + \|h_{\bullet}^{1/2} \nabla_{\Gamma} S(\phi - \phi_{\bullet})\|_{L^2(\bigcup \mathcal{R}_{\bullet, \circ})}^2 \\ &\lesssim \|h_{\bullet}^{-1/2} S(\phi - \phi_{\bullet})\|_{L^2(\bigcup \mathcal{R}_{\bullet, \circ})}^2 + \|h_{\bullet}^{1/2} \nabla_{\Gamma} S(\phi - \phi_{\bullet})\|_{L^2(\bigcup \mathcal{R}_{\bullet, \circ})}^2. \end{aligned} $$

Again, the first term has been already estimated before. Overall, we conclude that

$$ |\langle S(\phi - \phi_{\bullet}), \psi_{\circ} - \psi_{\bullet} \rangle_{L^2(\Gamma)}| \le \eta_{\bullet}(\mathcal{R}_{\bullet, \circ}) \|\psi_{\circ}\|_{H^{-1/2}(\Gamma)} $$

which completes the proof. $\square$

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

## 4. MIXED FORM OF INDIRECT INTEGRAL EQUATION

As discussed in Section 2.11, the regularizing operator $M$ can be realized with integral operators, or with the inverse of a differential operator. The former might be cumbersome since the Galerkin matrices associated to the composition of integral operators are typically challenging to assemble. In this section, we consider instead the latter option. Specifically, we consider the choice $M := (\alpha - \Delta_{\Gamma})^{-1} : H^{-1/2}(\Gamma) \to H^1(\Gamma)$, presented in Section 2.11.

While $M$ satisfies all required assumptions of Section 3, the discretization of (3.2) is not feasible, as $M$ cannot be computed exactly (even up to quadrature errors). Instead, we introduce similarly as in [11] the extra variable $f := M\phi$ and rewrite (3.2) as the following mixed system: Find $(\phi, f) \in H^{-1/2}(\Gamma) \times H^1(\Gamma)$ such that

$$ \begin{aligned} V_k \phi + i(K_k + 1/2)f &= u, \\ -\phi + \alpha f - \Delta_{\Gamma} f &= 0, \end{aligned} \quad (4.1) $$

where the equations hold in $H^{1/2}(\Gamma) \times H^{-1}(\Gamma)$. Since $V_k \in \mathcal{L}(H^{-1/2}(\Gamma), H^{1/2}(\Gamma))$ and $\alpha - \Delta_{\Gamma} \in \mathcal{L}(H^1(\Gamma), H^{-1}(\Gamma))$ are symmetric and elliptic at least up to some compact perturbation and $K_k \in \mathcal{L}(H^1(\Gamma), H^1(\Gamma))$, the operator $R = (R_1, R_2) \in \mathcal{L}(H^{-1/2}(\Gamma) \times H^1(\Gamma), H^{1/2}(\Gamma) \times H^{-1}(\Gamma))$ induced by the left-hand side of (4.1) is elliptic up to some compact perturbation. As (3.2) is uniquely solvable, so is (4.1). Hence, $R$ is injective and thus even bijective by the Fredholm alternative. We conclude that problem (4.1) is well-posed.
**Remark 4.1.** With minor modifications, the analysis of the whole Section 4 holds equally true for the (inverse of the) piecewise Laplace-Beltrami operator of (2.9). For brevity though, we do not provide the details.

**4.1. Discretization.** For any $\mathcal{T}_{\bullet} \in \mathbb{T}$, we consider the discrete subspace $\mathcal{P}^p(\mathcal{T}_{\bullet}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\bullet}) \subset H^{-1/2}(\Gamma) \times H^1(\Gamma)$. This leads to the discrete formulation: Find a Galerkin approximation $(\phi_{\bullet}, f_{\bullet}) \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})$ such that

$$ \langle R(\phi_{\bullet}, f_{\bullet}), (\psi_{\bullet}, g_{\bullet}) \rangle_{L^2(\Gamma)} = \langle u, \psi_{\bullet} \rangle_{L^2(\Gamma)} \quad \text{for all } (\psi_{\bullet}, g_{\bullet}) \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\bullet}). \quad (4.2) $$

To theoretically guarantee well-posedness of (4.2), we assume that the initial triangulation $\mathcal{T}_0$ is sufficiently fine such that

$$ \inf_{\substack{\chi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \\ e_{\bullet} \in \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})}} \sup_{\substack{\psi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \\ g_{\bullet} \in \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})}} \frac{|\langle R(\chi_{\bullet}, e_{\bullet}), (\psi_{\bullet}, g_{\bullet}) \rangle_{L^2(\Gamma)}|}{\|(\psi_{\bullet}, g_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)}} \gtrsim 1 \quad \text{for all } \mathcal{T}_{\bullet} \in \mathbb{T}, \quad (4.3) $$

where $\langle \cdot, \cdot \rangle_{L^2(\Gamma)}$ denotes the extended $L^2$-scalar product on $[H^{1/2}(\Gamma) \times H^{-1}(\Gamma)] \times [H^{-1/2}(\Gamma) \times H^1(\Gamma)]$. Owing to the fact that $R$ is symmetric and elliptic up to some compact perturbation, the existence of a maximal mesh size $h \gtrsim 1$ ensuring (4.3) whenever $h_T \le h$ for all $T \in \mathcal{T}_0$ is standard; see again, e.g., [47, Theorem 4.2.9] or [5, Proposition 1].

**Remark 4.2.** Given any $\mathcal{T}_{\bullet} \in \mathbb{T}$, the unique Galerkin approximation $(\phi_{\bullet}, f_{\bullet}) \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})$ of $(\phi, f)$ obtained by solving (4.2) satisfies the Céa lemma

$$ \|(\phi - \phi_{\bullet}, f - f_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)} \lesssim \min_{\substack{\psi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \\ g_{\bullet} \in \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})}} \|(\phi - \psi_{\bullet}, f - g_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)}. \quad (4.4) $$

**Remark 4.3.** For sufficiently smooth $\phi$ and $f$, one has that

$$ \begin{aligned} \min_{\psi_{\bullet} \in \mathcal{P}^p(\mathcal{T}_{\bullet})} \|\phi - \psi_{\bullet}\|_{H^{-1/2}(\Gamma)} &= \mathcal{O}(\|h_{\bullet}\|_{L^{\infty}(\Gamma)}^{p+3/2}) \\ \min_{g_{\bullet} \in \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})} \|f - g_{\bullet}\|_{H^1(\Gamma)} &= \mathcal{O}(\|h_{\bullet}\|_{L^{\infty}(\Gamma)}^{p+2}). \end{aligned} $$

The choice $p+2$ as polynomial degree to discretize $f$ thus does not spoil the convergence rate of $\|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)}$. That being said, the forthcoming analysis holds analogously for $p+2$ replaced by some arbitrary fixed degree $q \in \mathbb{N}$.

**4.2. A posteriori error estimation.** Our error estimator now consists of two parts. For the first equation in (4.1), we introduce

$$ \eta_{1, \bullet} := \|h_{\bullet}^{1/2} \nabla_{\Gamma}(u - V_k \phi_{\bullet} - i(K_k + 1/2) f_{\bullet})\|_{L^2(\Gamma)}, \quad (4.5a) $$

whereas

$$ \eta_{2, \bullet}^2 := \sum_{T \in \mathcal{T}_{\bullet}} h_T^2 \|\phi_{\bullet} - \alpha f_{\bullet} + \Delta_T f_{\bullet}\|_{L^2(T)}^2 + h_T \|[\nabla_{\Gamma} f_{\bullet} \cdot \bar{\nu}]\|_{L^2(\partial T)}^2 \quad (4.5b) $$

controls the residual of the second equation. We then set

$$ \eta_{\bullet}^2 := \eta_{1, \bullet}^2 + \eta_{2, \bullet}^2. \quad (4.5c) $$

Here, $\Delta_T$ denotes the local Laplace-Beltrami operator on $T$ and $[\nabla_{\Gamma}(\cdot) \cdot \bar{\nu}]$ the jump of the normal derivative, which reads on any non-empty face $F = T_{-} \cap T_{+}$, $T_{\pm} \in \mathcal{T}_{\bullet}$, as

$$ [\nabla_{\Gamma}(\cdot) \cdot \bar{\nu}]|_{F} = (\nabla_{\Gamma}(\cdot)|_{T_{-}} \cdot \bar{\nu}_{-} + \nabla_{\Gamma}(\cdot)|_{T_{+}} \cdot \bar{\nu}_{+}) $$
with the tangential normal vectors $\bar{\nu}_{\pm}$ on $\partial T_{\pm}$ pointing outside of $T_{\pm}$ and tangent to $T_{\pm}$; cf. [7, Section 3.1] for a precise definition in terms of the bi-Lipschitz parametrizations $\gamma_{T_{\pm}}$. Note that the jump is actually defined up to a sign (if the roles of $T_{\pm}$ are flipped), which is however not relevant here, as only the magnitude enters the definition of the estimator.

The first part (4.5a) of the estimator is similar to the estimator introduced for the non-mixed formulation in Section 3.2. It is the standard weighted-residual error estimator for integral equations of the first kind; see again [16, 17, 14, 4]. The second part (4.5b) is the standard residual estimator for Poisson problems, here applied on a manifold; see [39, 20, 44, 8, 7].

**Theorem 4.4.** There exists a constant $C_{\text{rel}} > 0$ such that

$$ \|(\phi, f) - (\phi_{\bullet}, f_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)} \le C_{\text{rel}}\eta_{\bullet} \quad \text{for all } \mathcal{T}_{\bullet} \in \mathcal{T}. \qquad (4.6) $$

If the boundary $\Gamma$ is polyhedral and all parametrizations $\gamma_T$, $T \in \mathcal{T}_0$, are chosen to be affine, there further exists a constant $C_{\text{eff}} > 0$ such that

$$ \eta_{\bullet} \le C_{\text{eff}} (\|h_{\bullet}^{1/2}(\phi - \phi_{\bullet})\|_{L^2(\Gamma)} + \|f - f_{\bullet}\|_{H^1(\Gamma)}). \qquad (4.7) $$

The constants $C_{\text{rel}}$ and $C_{\text{eff}}$ depend only on $d, k, \Gamma, \alpha$, and $\mathcal{T}_0$.

**Remark 4.5.** Similarly as in Remark 3.3, the discrete inf-sup stability (4.3) is not required for reliability (4.6), which holds in fact true for any (in general non-unique) $(\phi_{\bullet}, f_{\bullet}) \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})$ satisfying (4.2). This is also the case for weak efficiency (4.7) if the right-hand side is replaced by $\|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|h_{\bullet}^{1/2}(\phi - \phi_{\bullet})\|_{L^2(\Gamma)} + \|f - f_{\bullet}\|_{H^1(\Gamma)}$. Moreover, the assumption of a polyhedral boundary $\Gamma$ is not essential for weak efficiency (4.7), but additional oscillation terms are present in the general case; see [20, 8, 7].

*Proof of reliability.* Since $R$ is an isomorphism, we can readily estimate the discretization error by the norm of the residual, i.e.,

$$ \|(\phi, f) - (\phi_{\bullet}, f_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)} \approx \|(u, 0) - R(\phi_{\bullet}, f_{\bullet})\|_{H^{1/2}(\Gamma) \times H^{-1}(\Gamma)}. $$

As the first component $u - R_1(\phi_{\bullet}, f_{\bullet}) = u - V_k\phi_{\bullet} - i(K_k + 1/2)f_{\bullet}$ is orthogonal to piecewise constants $\mathcal{P}^0(\mathcal{T}_{\bullet})$, its norm can be estimated as in the proof of Theorem 3.2, i.e.,

$$ \|u - R_1(\phi_{\bullet}, f_{\bullet})\|_{H^{1/2}(\Gamma)}^2 \approx \sum_{z \in \mathcal{N}_{\bullet}} |u - R_1(\phi_{\bullet}, f_{\bullet})|^2_{H^{1/2}(\omega_{\bullet}(z))} \lesssim \eta_{1, \bullet}^2. \qquad (4.8) $$

For the second component $R_2(\phi_{\bullet}, f_{\bullet}) = -\phi_{\bullet} + \alpha f_{\bullet} - \Delta_{\Gamma} f_{\bullet}$, standard arguments as for the Poisson problem show that

$$ \|R_2(\phi_{\bullet}, f_{\bullet})\|_{H^{-1}(\Gamma)}^2 \lesssim \eta_{2, \bullet}^2; \qquad (4.9) $$

cf. [39] and [20, 44, 8, 7], which even take some potential geometry approximation into account. $\square$

*Proof of weak efficiency.* In case of a polyhedral boundary $\Gamma$ (and affine parametrizations), standard arguments and the fact that $\phi_{\bullet}$ is discrete also yield the converse estimate of (4.8) without any oscillation term, i.e., $\eta_{2, \bullet} \lesssim \|R_2(\phi_{\bullet}, f_{\bullet})\|_{H^{-1}(\Gamma)}$; cf. [20, 8, 7]. With the Céa lemma (4.4), one then sees as in the proof of Theorem 3.2 the overall weak efficiency (4.7). $\square$
4.3. **Adaptive algorithm.** We abbreviate the weighted-residual error indicators by

$$ \eta_{1,\bullet}(T)^2 := \|h_{\bullet}^{1/2} \nabla_{\Gamma}(u - R_1(\phi_{\bullet}, f_{\bullet}))\|_{L^2(T)}^2, \quad (4.10a) $$

$$ \eta_{2,\bullet}(T)^2 := h_T^2 \|\phi_{\bullet} - \alpha f_{\bullet} + \Delta_T f_{\bullet} \|_{L^2(T)}^2 + h_T \|[\nabla_{\Gamma} f_{\bullet} \cdot \bar{\nu}]\|_{L^2(\partial T)}^2, \quad (4.10b) $$

$$ \eta_{\bullet}(T)^2 := \eta_{1,\bullet}(T)^2 + \eta_{2,\bullet}(T)^2 \quad \text{for all } T \in \mathcal{T}_{\bullet}. \quad (4.10c) $$

Moreover, we set

$$ \eta_{j,\bullet}(\mathcal{U}_{\bullet})^2 := \sum_{T \in \mathcal{U}_{\bullet}} \eta_{j,\bullet}(T)^2, \quad \eta_{\bullet}(\mathcal{U}_{\bullet})^2 := \eta_{1,\bullet}(\mathcal{U}_{\bullet})^2 + \eta_{2,\bullet}(\mathcal{U}_{\bullet})^2 \quad \text{for all } j \in \{1, 2\}, \mathcal{U}_{\bullet} \subseteq \mathcal{T}_{\bullet}. $$

Note that $\eta_{1,\bullet} = \eta_{1,\bullet}(\mathcal{T}_{\bullet})$, $\eta_{2,\bullet} = \eta_{2,\bullet}(\mathcal{T}_{\bullet})$, and $\eta_{\bullet} = \eta_{\bullet}(\mathcal{T}_{\bullet})$. As in Section 3.3, we can use the weighted-residual error estimator within a standard adaptive algorithm to steer the local mesh refinement of some given sufficiently fine conforming initial triangulation $\mathcal{T}_0$ such that the discrete inf-sup stability (4.3) is satisfied.

**Algorithm 4.6.** Input: Dörfler parameter $0 < \theta \le 1$.

Loop: For each $\ell = 0, 1, 2, \dots$, iterate the following steps:

(i) Compute the Galerkin approximation $(\phi_{\ell}, f_{\ell}) \in \mathcal{P}^p(\mathcal{T}_{\ell}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\ell})$ by solving (4.2).

(ii) Compute the refinement indicators $\eta_{\ell}(T)$ from (4.10) for all elements $T \in \mathcal{T}_{\ell}$.

(iii) Determine a minimal set of marked elements $\mathcal{M}_{\ell} \subseteq \mathcal{T}_{\ell}$ satisfying the Dörfler marking

$$ \theta \eta_{\ell}^2 \le \eta_{\ell}(\mathcal{M}_{\ell})^2. $$

(iv) Generate refined triangulation $\mathcal{T}_{\ell+1} := \text{refine}(\mathcal{T}_{\ell}, \mathcal{M}_{\ell})$.

**Output:** Refined triangulations $\mathcal{T}_{\ell}$, corresponding Galerkin approximations $(\phi_{\ell}, f_{\ell}) \in \mathcal{P}^p(\mathcal{T}_{\ell}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\ell})$, and weighted-residual error estimators $\eta_{\ell}$ for all $\ell \in \mathbb{N}_0$.

4.4. **Optimal convergence.** We define the estimator-based approximation norm $\|(\phi, f)\|_{\mathcal{A}_{\text{est}}}$ as in (3.11), and say again that the solution pair $(\phi, f) \in H^{-1/2}(\Gamma) \times H^1(\Gamma)$ lies in the approximation class $s$ with respect to the estimator if $\|(\phi, f)\|_{\mathcal{A}_{\text{est}}} < \infty$. The following theorem is the analogon of Theorem 3.6. It states that each possible rate $s > 0$ is indeed realized by Algorithm 3.4. The proof is given in Section 4.5, where we will verify again the axioms of adaptivity from [15] for the error estimator, which automatically imply the desired result.

**Theorem 4.7.** Suppose that the discrete inf-sup stability (4.3) is satisfied. Let $(\mathcal{T}_{\ell})_{\ell \in \mathbb{N}_0}$ be the sequence of triangulations generated by Algorithm 3.4. Then, for arbitrary $0 < \theta \le 1$, the estimator converges linearly, i.e., there exist constants $0 < \rho_{\text{lin}} < 1$ and $C_{\text{lin}} \ge 1$ such that

$$ \eta_{\ell+j}^2 \le C_{\text{lin}} \rho_{\text{lin}}^j \eta_{\ell}^2 \quad \text{for all } j, \ell \in \mathbb{N}_0. $$

Moreover, there exists a constant $0 < \theta_{\text{opt}} \le 1$ such that for all $0 < \theta < \theta_{\text{opt}}$, the estimator converges at optimal rate, i.e., for all $s > 0$, there exist constants $c_{\text{opt}}, C_{\text{opt}} > 0$ such that

$$ c_{\text{opt}} \|(\phi, f)\|_{\mathcal{A}_{\text{est}}} \le \sup_{\ell \in \mathbb{N}_0} (\#\mathcal{T}_{\ell} - \#\mathcal{T}_0 + 1)^s \eta_{\ell} \le C_{\text{opt}} \|(\phi, f)\|_{\mathcal{A}_{\text{est}}}. $$

The constant $\theta_{\text{opt}}$ depends only on $d, k, \Gamma, \alpha, \mathcal{T}_0$, and $p$, while $C_{\text{lin}}, \rho_{\text{lin}}$ depend additionally on $\theta$, and $C_{\text{opt}}$ depends furthermore on $s > 0$. The constant $c_{\text{opt}}$ depends only on $d, \#\mathcal{T}_0, s$, and if there exists $\ell_0$ with $\eta_{\ell_0} = 0$, then also on $\ell_0$ and $\eta_0$.
4.5. **Proof of Theorem 4.7.** The axioms of adaptivity for the estimator from [15] now take the following form: There exist uniform constants $C_{\text{stab}}$, $C_{\text{red}}$, $C_{\text{drel}}$, $C_{\text{qo}} > 0$ and $0 < \rho_{\text{red}} < 1$ such that for all $\mathcal{T}_{\bullet} \in \mathbb{T}$, $\mathcal{T}_{\circ} \in \text{refine}(\mathcal{T}_{\bullet})$ with corresponding patch of all refined elements $\mathcal{R}_{\bullet, \circ} := \{T \in \mathcal{T}_{\bullet} : T \cap \bigcup(\mathcal{T}_{\bullet} \setminus \mathcal{T}_{\circ}) \neq \emptyset\}$, the sequence $(\mathcal{T}_{\ell})_{\ell \in \mathbb{N}_0}$ of Algorithm 4.6, and all $\ell, n \in \mathbb{N}_0$, there hold:

(A1) Stability on non-refined elements: $|\eta_{\circ}(\mathcal{T}_{\circ} \cap \mathcal{T}_{\bullet}) - \eta_{\bullet}(\mathcal{T}_{\circ} \cap \mathcal{T}_{\circ})| \le C_{\text{stab}} \|(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)}$;

(A2) Reduction on refined elements: $\eta_{\circ}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\bullet}) \le \rho_{\text{red}}\eta_{\bullet}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\circ}) + C_{\text{red}}\|(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)}$;

(A3) General quasi-orthogonality: $\sum_{j=\ell}^{\ell+n} \|(\phi_{j+1} - \phi_j, f_{j+1} - f_j)\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)}^2 \le C_{\text{qo}}^2 \eta_{\ell}^2$;

(A4) Discrete reliability: $\|(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)} \le C_{\text{drel}}\eta_{\bullet}(\mathcal{R}_{\circ, \circ})$.

*Proof of stability on non-refined elements* (A1). Let $\mathcal{T}_{\bullet} \in \mathbb{T}$ and $\mathcal{T}_{\circ} \in \text{refine}(\mathcal{T}_{\bullet})$. By the inverse triangle inequality, it suffices to prove the estimate for $\eta_1$ and $\eta_2$. The inverse triangle inequality and the fact that $h_{\bullet} = h_{\circ}$ on $\omega := \bigcup(\mathcal{T}_{\circ} \cap \mathcal{T}_{\circ})$ show that

$$|\eta_{1,\circ}(\mathcal{T}_{\circ} \cap \mathcal{T}_{\bullet}) - \eta_{1,\bullet}(\mathcal{T}_{\circ} \cap \mathcal{T}_{\circ})| \le \|h_{\circ}^{1/2} \nabla_{\Gamma} R_1(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet})\|_{L^2(\Gamma)}.$$

The inverse estimates (2.5a)–(2.5b), and (2.6) on the mesh $\mathcal{T}_{\circ}$ imply that

$$\begin{aligned} \|h_{\circ}^{1/2} \nabla_{\Gamma} R_1(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet})\|_{L^2(\Gamma)} &\lesssim \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|f_{\circ} - f_{\bullet}\|_{H^{1/2}(\Gamma)} \\ &\quad + \|h_{\bullet}^{1/2} \nabla_{\Gamma}(f_{\circ} - f_{\bullet})\|_{L^2(\Gamma)} \\ &\lesssim \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|f_{\circ} - f_{\bullet}\|_{H^1(\Gamma)}. \end{aligned}$$

Similarly, standard inverse estimates and a trace inequality (in the parameter domain) show that

$$\begin{aligned} |\eta_{2,\circ}(\mathcal{T}_{\circ} \cap \mathcal{T}_{\bullet}) - \eta_{2,\bullet}(\mathcal{T}_{\circ} \cap \mathcal{T}_{\circ})| &\le \|h_{\circ}(\phi_{\circ} - \phi_{\bullet})\|_{L^2(\Gamma)} + \|f_{\circ} - f_{\bullet}\|_{H^1(\Gamma)} \\ &\lesssim \|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|f_{\circ} - f_{\bullet}\|_{H^1(\Gamma)}. \end{aligned}$$

This concludes the proof. $\square$

*Reduction on refined elements* (A2). Let $\mathcal{T}_{\bullet} \in \mathbb{T}$ and $\mathcal{T}_{\circ} \in \text{refine}(\mathcal{T}_{\bullet})$. Again, we can estimate $\eta_1$ and $\eta_2$ separately. The triangle inequality and the fact that $h_{\circ}^{1/2} \le \rho h_{\bullet}^{1/2}$ on $\omega := \bigcup(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\bullet}) = \bigcup(\mathcal{T}_{\bullet} \setminus \mathcal{T}_{\circ})$ for some $0 < \rho < 1$ show that

$$\eta_{1,\circ}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\bullet}) \le \rho \eta_{1,\bullet}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\circ}) + \|h_{\circ}^{1/2} \nabla_{\Gamma} R_1(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet})\|_{L^2(\Gamma)}.$$

The second term has already been estimated in the proof of (A1). Similarly, the fact that $[\nabla_{\Gamma} f_{\circ} \cdot \bar{\nu}] = 0$ on edges introduced within the refinement of $\mathcal{T}_{\bullet}$ to $\mathcal{T}_{\circ}$, standard inverse estimates, and a trace inequality show (as in the proof of (A1)) for some uniform constants $C, C' > 0$ that

$$\begin{aligned} \eta_{2,\circ}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\bullet}) &\le \rho \eta_{2,\bullet}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\circ}) + C (\|h_{\circ}(\phi_{\circ} - \phi_{\bullet})\|_{L^2(\Gamma)} + \|f_{\circ} - f_{\bullet}\|_{H^1(\Gamma)}) \\ &\le \rho \eta_{2,\bullet}(\mathcal{T}_{\circ} \setminus \mathcal{T}_{\circ}) + C' (\|\phi_{\circ} - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)} + \|f_{\circ} - f_{\bullet}\|_{H^1(\Gamma)}). \end{aligned}$$

This concludes the proof. $\square$

*Proof of general quasi-orthogonality* (A3). By discrete inf-sup stability (4.3), [25] is applicable and directly implies this axiom with $C_{\text{qo}} = C_{\text{qo}}(n) = o(n)$ depending sublinearly on the
integer $n$. The post-processing argument [25, Remark 20] or [15, Proposition 4.11] then yield the axiom with a uniform constant $C_{qo} \simeq 1$. $\square$

**Remark 4.8.** As in Remark 3.8, the analysis of [28, 5] allows for an alternative proof of (A3).

*Proof of discrete reliability* (A4). Let $\mathcal{T}_{\bullet} \in \mathcal{T}$ and $\mathcal{T}_{\circ} \in \text{refine}(\mathcal{T}_{\bullet})$. Discrete inf-sup stability (4.3) and the definition of $(\phi_{\circ}, f_{\circ})$ as well $(\phi_{\bullet}, f_{\bullet})$ (see (4.2)) show for all $(\psi_{\bullet}, g_{\bullet}) \in \mathcal{P}^p(\mathcal{T}_{\bullet}) \times \mathcal{S}^{p+2}(\mathcal{T}_{\bullet})$ that

$$ \begin{aligned} \|(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)} &\lesssim \sup_{\substack{\psi_{\circ} \in \mathcal{P}^p(\mathcal{T}_{\circ}) \\ g_{\circ} \in \mathcal{S}^{p+2}(\mathcal{T}_{\circ})}} \frac{|\langle R(\phi_{\circ} - \phi_{\bullet}, f_{\circ} - f_{\bullet}), (\psi_{\circ}, g_{\circ}) \rangle_{L^2(\Gamma)}}{|\|(\psi_{\circ}, g_{\circ})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)}} \\ &= \sup_{\substack{\psi_{\circ} \in \mathcal{P}^p(\mathcal{T}_{\circ}) \\ g_{\circ} \in \mathcal{S}^{p+2}(\mathcal{T}_{\circ})}} \frac{|\langle R(\phi - \phi_{\bullet}, f - f_{\bullet}), (\psi_{\circ} - \psi_{\bullet}, g_{\circ} - g_{\bullet}) \rangle_{L^2(\Gamma)}}{|\|(\psi_{\circ}, g_{\circ})\|_{H^{-1/2}(\Gamma) \times H^1(\Gamma)}}. \end{aligned} $$

Choosing $\psi_{\bullet}$ as in (3.14), we see along the same lines of the non-mixed case that

$$ |\langle R_1(\phi - \phi_{\bullet}, f - f_{\bullet}), \psi_{\circ} - \psi_{\bullet} \rangle_{L^2(\Gamma)}| \lesssim \eta_{1, \bullet}(\mathcal{R}_{\bullet, \circ}) \|\psi_{\circ}\|_{H^{-1/2}(\Gamma)} $$

Choosing $g_{\bullet}$ as the Scott-Zhang projection of $g_{\circ}$, standard arguments as for the Poisson problem show that

$$ |\langle R_2(\phi - \phi_{\bullet}, f - f_{\bullet}), g_{\circ} - g_{\bullet} \rangle_{L^2(\Gamma)}| \lesssim \eta_{2, \bullet}(\mathcal{R}_{\bullet, \circ}) \|g_{\circ}\|_{H^1(\Gamma)}; $$

cf. [8, 7]. This concludes the proof. $\square$

## 5. DIRECT INTEGRAL EQUATIONS

Given again a regularizing operator $M \in \mathcal{L}(H^{-1/2}(\Gamma), H^{1/2}(\Gamma))$ satisfying (M1) and (M2), we consider the representation formula (2.3) with corresponding Calderón system

$$ U|_{\Gamma} = (K_k + 1/2)(U|_{\Gamma}) - V_k(\partial_{\nu}U), \quad (5.1a) $$

$$ \partial_{\nu} U = -W_k(U|_{\Gamma}) - (K'_k - 1/2)(\partial_{\nu} U). \quad (5.1b) $$

We now apply $iM$ to (5.1b), add it to (5.1a), use the Dirichlet condition $U|_{\Gamma} = u$, and obtain for $\phi := \partial_{\nu}U$ that

$$ iM\phi + u = -iMW_k u - iM(K'_k - 1/2)\phi + (K_k + 1/2)u - V_k\phi, $$

or equivalently,

$$ V_k\phi + iM(K'_k + 1/2)\phi = (K_k - 1/2)u - iMW_k u. \quad (5.2) $$

Clearly, the operator

$$ S := V_k + iM(K'_k + 1/2) \in \mathcal{L}(H^{-1/2}(\Gamma), H^{1/2}(\Gamma)) \quad (5.3) $$

on the left-hand side of (5.2) is elliptic and symmetric up to the compact perturbation $V_k - V_0^+ + iM(K'_k + 1/2)$. According to [11], $S$ is also injective and thus even bijective by the Fredholm alternative.

### 5.1. Non-mixed form.
We can follow the lines of Section 3 to show reliability and weak efficiency of the weighted-residual estimator $\eta_{\bullet} := \|h_{\bullet}^{1/2}\nabla_{\Gamma}(u - S\phi_{\bullet})\|_{L^2(\Gamma)}$ and optimal convergence of a corresponding adaptive algorithm. In particular, the mapping property $K'_k \in \mathcal{L}(L^2(\Gamma), L^2(\Gamma))$ of (2.4c) and the inverse estimates (2.5c) and (M3) are exploited.
5.2. **Mixed form.** Similarly as in Section 4, we consider $M := (\alpha - \Delta_{\Gamma})^{-1}$ as mapping from $H^{-1/2}(\Gamma) \subset H^{-1}(\Gamma)$ to $H^1(\Gamma) \subset H^{1/2}(\Gamma)$ and set

$$f := M[(K'_{k} + 1/2)\phi + W_{k}u]. \qquad (5.4)$$

As a matter of fact, $f = 0$ owing to (5.1b). We obtain the mixed system

$$V_k \phi + if = (K_k - 1/2)u, \qquad (5.5)$$

$$-(K'_{k} + 1/2)\phi + \alpha f - \Delta_{\Gamma} f = W_{k}u.$$

in $H^{1/2}(\Gamma) \times H^{-1}(\Gamma)$ with sought $(\phi, f) \in H^{-1/2}(\Gamma) \times H^1(\Gamma)$. Clearly, the operator $R = (R_1, R_2) : H^{-1/2}(\Gamma) \times H^1(\Gamma) \to H^{1/2}(\Gamma) \times H^{-1}(\Gamma)$ induced by the left-hand side of (5.5) is elliptic up to some compact perturbation. Since (5.2) is uniquely solvable, $R$ is also injective and thus even bijective by the Fredholm alternative.

It is well-known that if the Dirichlet datum satisfies $u = U|_{\Gamma} \in H^1(\Gamma)$, then $\phi = \partial_{\nu}U \in L^2(\Gamma)$; see, e.g., [43, Theorem 4.24]. We further recall the mapping properties $K_k \in \mathcal{L}(H^1(\Gamma), H^1(\Gamma))$ and $W_k \in \mathcal{L}(H^1(\Gamma), L^2(\Gamma))$ of (2.4). Thus, we can follow the lines of Section 4 to show reliability and weak efficiency of the weighted-residual estimator $\eta_{\bullet} = (\eta_{1,\bullet}^2 + \eta_{2,\bullet}^2)^{1/2}$ with

$$\eta_{1,\bullet}(T)^2 := \|h_{\bullet}^{1/2} \nabla_{\Gamma}((K_k - 1/2)u - V_k \phi_{\bullet} - if_{\bullet})\|_{L^2(T)}^2,$$

$$\eta_{2,\bullet}(T)^2 := h_T^2 \|W_k u + (K'_k + 1/2)\phi_{\bullet} - \alpha f_{\bullet} + \Delta_T f_{\bullet}\|_{L^2(T)}^2 + h_T \|[\nabla_{\Gamma} f_{\bullet} \cdot \bar{\nu}]\|_{L^2(\partial T)}^2,$$

$$\eta_{\bullet}(T)^2 := \eta_{1,\bullet}(T)^2 + \eta_{2,\bullet}(T)^2 \quad \text{for all } T \in \mathcal{T}_{\bullet}.$$

and optimal convergence of a corresponding adaptive algorithm.

## 6. NUMERICAL EXPERIMENTS

We consider the exterior Dirichlet problem (2.2) with $k^2$ close to an eigenvalue of the interior Dirichlet problem on the following two bounded Lipschitz domains $\Omega$: the circle

$$\Omega = \{x \in \mathbb{R}^2 : |x| < \frac{1}{10}\}$$

and the L-shape

$$\Omega = \text{conv}\{(\frac{1}{10}, 0), (0, \frac{1}{10}), (0, 0), (0, \frac{-1}{10})\} \setminus \text{conv}\{(0, 0), (\frac{-1}{20}, \frac{1}{20}), (\frac{-1}{10}, 0), (\frac{-1}{20}, \frac{-1}{20})\}.$$

Note that both domains satisfy $\text{diam}(\Omega) < 1$ so that the single-layer operator $V_0$ associated to the Laplace operator is elliptic. As exact solution $U$, we prescribe the fundamental solution $U(x_1, x_2) := G_k(x_1, x_2 - \frac{1}{20})$ with singularity in the interior of $\Omega$. We apply the adaptive Algorithm 4.6 with $p = 0$ and $\theta \in \{1, 0.9\}$ for the mixed forms (4.1) and (5.5) of the indirect combined field integral equation (3.2) and the direct combined regularized field integral equation (5.2), respectively. We select the scaling parameter $\alpha = 1$ in the definition of $M$, i.e., $M := (1 - \Delta_{\Gamma})^{-1.2}$. Note that $\theta = 1$ just leads to uniform refinement in each step. As initial mesh $\mathcal{T}_0$ for the boundary of the circle, we choose the four arcs $\{x \in \Gamma : x_1, x_2 \ge 0\}, \{x \in \Gamma : x_1 \le 0, x_2 \ge 0\}, \{x \in \Gamma : x_1, x_2 \le 0\}, \{x \in \Gamma : x_1 \ge 0, x_2 \le 0\}$ and parametrize them using their arclength parametrization (up to scaling to the reference element $T_{\text{ref}} = (0, 1)$). As initial mesh $\mathcal{T}_0$ for the boundary of the L-shape, we simply choose the six straight lines that define the boundary and equip them with affine parametrizations.

<footnote n="2">²This choice is dimensionally consistent as $\text{diam}(\Omega) \sim 1$ in our examples. Setting $\alpha = k^2$ is also appealing, but we expect little difference in the results, since the frequencies considered here are relatively low.</footnote>
We stress that it is in principle unclear whether the initial meshes are indeed sufficiently fine to guarantee the assumed theoretically required discrete inf-sup stability on all refinements of $\mathcal{T}_0$ (see (4.3) in case of the indirect combined field integral formulation).

For comparison, we also consider the standard direct and indirect integral equations (1.2) and (1.4), discretized via

$$ \langle V_k \phi_{\bullet}, \psi_{\bullet} \rangle_{L^2(\Gamma)} = \langle u, \psi_{\bullet} \rangle_{L^2(\Gamma)} \quad \text{for all } \psi_{\bullet} \in \mathcal{P}^0(\mathcal{T}_{\bullet}) \qquad (6.1) $$

and

$$ \langle V_k \phi_{\bullet}, \psi_{\bullet} \rangle_{L^2(\Gamma)} = \langle (K_k - 1/2)u, \psi_{\bullet} \rangle_{L^2(\Gamma)} \quad \text{for all } \psi_{\bullet} \in \mathcal{P}^0(\mathcal{T}_{\bullet}), $$

respectively, with sought $\phi_{\bullet} \in \mathcal{P}^0(\mathcal{T}_{\bullet})$. In this case, we employ the weighted-residual estimator $\|h_{\bullet}^{1/2} \nabla_{\Gamma} V_k(\phi - \phi_{\bullet})\|_{L^2(\Gamma)}$ to steer adaptive refinement. For $k^2$ not being an eigenvalue of the interior Dirichlet problem, [4] shows that this estimator is indeed reliable and the algorithm converges at optimal rate.

In any case, the numerical realization of the adaptive algorithms is based on the MATLAB implementation [35, 34]. The non-discrete Dirichlet datum $u = U|_{\Gamma}$ is replaced by its $L^2$-orthogonal projection $u_{\bullet}$ onto $S^2(\mathcal{T}_{\bullet})$. For the Laplace problem, local data oscillation terms that control the additional error $\|u - u_{\bullet}\|_{H^{1/2}(\Gamma)}$ are investigated in [26], where it is shown that the algorithms still converge at optimal rate. The results of [26] can be readily extended to the considered regularized combined field integral formulations. As a matter of fact, to obtain higher-order oscillations, it would suffice to take $u_{\bullet} \in S^1(\mathcal{T}_{\bullet})$. However, since $f$ is discretized in $S^2(\mathcal{T}_{\bullet})$, the required mass matrix is available anyway.

For the direct approaches, the exact solution $\phi$ is explicitly given as $\phi = \partial_{\nu} U$. The error $\|\phi - \phi_{\bullet}\|_{V_0} := \langle V_0(\phi - \phi_{\bullet}), \phi - \phi_{\bullet} \rangle^{1/2}$ (being equivalent to $\|\phi - \phi_{\bullet}\|_{H^{-1/2}(\Gamma)}$) is approximately computed by replacing $\phi$ with its $L^2$-projection onto $\mathcal{P}^1(\mathcal{T}_{\bullet})$. We expect that adaptivity (i.e., $\theta = 0.9$) yields a convergence rate $O((\#\mathcal{T}_{\bullet})^{-3/2})$ for the approximate error and the weighted-residual estimator.

**6.1. Experiments on circle.** The eigenvalues $k^2$ of the interior Dirichlet problem on the unit circle are given as the squares of the roots of the Bessel functions of the first kind. On the scaled circle

$$ \Omega = \{x \in \mathbb{R}^2 : |x| < \frac{1}{10}\}, $$

an eigenvalue $k^2$ is thus given as 100 times the square of the first root of $J_0$, i.e., $k_O := 24.04825558...$. For $k \in \{k_O + 10^1, k_O + 10^0, \dots, k_O + 10^{-6}, k_O\}$, we plot in Figure 6.1 and 6.2 the (approximate) errors $\|\phi - \phi_{\ell}\|_{V_0}$ and estimators $\eta_{\ell}$ for BEM based on the standard integral equations (1.2) and BEM based on the combined field integral equations (4.1) and (5.5) under uniform ($\theta = 1$) and adaptive ($\theta = 0.9$) refinement.

While for $k$ sufficiently far away from $k_O$, the estimators in Figure 6.1 corresponding to the indirect integral approach converge for each method at optimal rate $O((\#\mathcal{T}_{\ell})^{-3/2})$, the involved multiplicative constant deteriorates quickly for standard adaptive BEM as $k$ gets closer to $k_O$. For $k = k_O$, no convergence at all is observed.

In Figure 6.2 corresponding to the direct integral approach, each method yields the optimal rate $O((\#\mathcal{T}_{\ell})^{-3/2})$ for the error estimator. Unsurprisingly however, although the convergence of the estimator remains unaffected as $k$ gets closer to $k_O$, the multiplicative constant for the error deteriorates for standard BEM under both uniform and adaptive refinement. In particular, the effectivity index, i.e., the ratio between error and estimator, blows up. This is
not the case for BEM based on combined field integral equations, where error and estimator remain close to each other.

Recall that for the combined field integral equations, the estimator actually consists of two parts $\eta_{1,\ell}$ and $\eta_{2,\ell}$. We only mention that the second part $\eta_{2,\ell}$ converges for the indirect and the direct approach at higher rate $O((\#T_\ell)^{-2})$ (not displayed). This higher rate can be motivated by the following arguments: As mentioned in the proof of Theorem 4.4 and Remark 4.5 for the indirect approach, the second part $\eta_{2,\ell}$ is up to potential oscillations equivalent to $\| -\phi_\ell + f_\ell - \Delta_\Gamma f_\ell \|_{H^{-1}(\Gamma)}$ and $\| W_k u + (K'_k + 1/2)\phi_\ell - f_\ell + \Delta_\Gamma f_\ell \|_{H^{-1}(\Gamma)}$, respectively. The triangle inequality and stability of the involved operators show that both terms are bounded by $\|\phi - \phi_\ell\|_{H^{-1}(\Gamma)} + \|f - f_\ell\|_{H^1(\Gamma)}$. From standard approximation theory and $\phi_\ell \in P^0(T_\ell)$ and $f_\ell \in S^2(T_\ell)$, the rate $\eta_{2,\ell} = O((\#T_\ell)^{-2})$ might indeed be expected.

## 6.2. Experiments on L-shape.
An exact eigenvalue $k^2$ of the interior Dirichlet problem on the L-shape

$$ \Omega = \operatorname{conv}\{(\frac{1}{10}, 0), (0, \frac{1}{10}), (0, 0), (0, \frac{-1}{10})\} \setminus \operatorname{conv}\{(0, 0), (\frac{-1}{20}, \frac{1}{20}), (\frac{-1}{10}, 0), (\frac{-1}{20}, \frac{-1}{20})\} $$

is given for $k_L := 20\pi$. For $k \in \{k_L + 10^1, k_L + 10^0, \dots, k_L + 10^{-6}, k_L\}$, we plot in Figure 6.3 and 6.4 the (approximate) errors $\|\phi - \phi_\ell\|_{V_0}$ and estimators $\eta_\ell$ for BEM based on the standard integral equations (1.2) and BEM based on the combined field integral equations (4.1) and (5.5) under uniform ($\theta = 1$) and adaptive ($\theta = 0.9$) refinement.

For $k = k_L + 10$, the estimators in Figure 6.3 corresponding to the indirect integral approach converge at the suboptimal rate $O((\#T_\ell)^{-2/3})$ under uniform refinement. Adaptive refinement regains the optimal rate $O((\#T_\ell)^{-3/2})$. The same holds for BEM based on combined field integral equations as $k$ gets closer to $k_L$. Instead, for standard BEM, we even observe a strong growth of the estimator in the increasing preasymptotic phase. The seemingly optimal convergence under uniform refinement for $k \neq k_L + 10$ also appears to be a preasymptotic behavior as is suggested by the plot for $k = k_L + 10$ (and other values between $k_L + 1$ and $k_L + 10$ which are not displayed).

In Figure 6.4 corresponding to the direct integral approach, we observe a similar behavior as for the circle considered in Section 6.1. We also mention that the second part $\eta_{2,\ell}$ again converges for the indirect and the direct approach at higher rate $O((\#T_\ell)^{-2})$ even for uniform refinement (not displayed).

## ACKNOWLEDGEMENT

GG acknowledges funding by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany's Excellence Strategy - EXC-2047/1 - 390685813.

TCF acknowledges funding by the ANR JCJC project APOWA (research grant ANR-23-CE40-0019-01).

## REFERENCES

[1] M. Aurada, M. Feischl, T. Führer, M. Karkulik, J. M. Melenk, and D. Praetorius. Local inverse estimates for non-local boundary integral operators. *Math. Comp.*, 86(308):2651–2686, 2017.

[2] M. Aurada, M. Feischl, T. Führer, M. Karkulik, and D. Praetorius. Efficiency and optimality of some weighted-residual error estimator for adaptive 2D boundary element methods. *Comput. Methods Appl. Math.*, 13(3):305–332, 2013.

[3] M. Aurada, M. Feischl, T. Führer, M. Karkulik, and D. Praetorius. Energy norm based error estimators for adaptive BEM for hypersingular integral equations. *Appl. Numer. Math.*, 95:15–35, 2015.
<figure><img src="image_5.png" /><figcaption>FIGURE 6.1. Double-logarithmic convergence plots over number of elements #T<sub>l</sub> for circle and indirect integral approach. Estimators η<sub>l</sub> are displayed for BEM based on the standard integral equation (1.2) and BEM based on the combined field integral equation (4.1) under uniform (θ = 1) and adaptive (θ = 0.9) refinement. Reference lines are of order O((#T<sub>l</sub>)<sup>-3/2</sup>).</figcaption></figure>

<footnote>[4] A. Bespalov, T. Betcke, A. Haberl, and D. Praetorius. Adaptive BEM with optimal convergence rates for the Helmholtz equation. *Comput. Methods Appl. Mech. Engrg.*, 346:260–287, 2019.</footnote>

<footnote>[5] A. Bespalov, A. Haberl, and D. Praetorius. Adaptive FEM with coarse initial mesh guarantees optimal convergence rates for compactly perturbed elliptic problems. *Comput. Methods Appl. Mech. Engrg.*, 317:318–340, 2017.</footnote>
<figure><img src="image_6.png" /><figcaption>FIGURE 6.2. Double-logarithmic convergence plots over number of elements #T<sub>l</sub> for circle and direct integral approach. (Approximate) errors ||φ − φ<sub>l</sub>||<sub>V<sub>0</sub></sub> and estimators η<sub>l</sub> are displayed for BEM based on the standard integral equation (1.2) and BEM based on the combined field integral equation (5.5) under uniform (θ = 1) and adaptive (θ = 0.9) refinement. Reference lines are of order O(||T<sub>l</sub>||<sup>−3/2</sup>).</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>FIGURE 6.3. Double-logarithmic convergence plots over number of elements #T<sub>l</sub> for L-shape and indirect integral approach. Estimators η<sub>l</sub> are displayed for BEM based on the standard integral equation (1.2) and BEM based on the combined field integral equation (4.1) under uniform (θ = 1) and adaptive (θ = 0.9) refinement. Reference lines are of order O((#T<sub>l</sub>)<sup>-2/3</sup>) and O((#T<sub>l</sub>)<sup>-3/2</sup>).</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>FIGURE 6.4. Double-logarithmic convergence plots over number of elements #T<sub>l</sub> for L-shape and direct integral approach. (Approximate) errors ||<sub>l</sub> - <sub>0</sub>||<sub>V<sub>0</sub></sub> and estimators η<sub>l</sub> are displayed for BEM based on the standard integral equation (1.2) and BEM based on the combined field integral equation (5.5) under uniform (θ = 1) and adaptive (θ = 0.9) refinement. Reference lines are of order O((#T<sub>l</sub>)<sup>-3/2</sup>).</figcaption></figure>

<footnote>[8] A. Bonito, J. M. Cascón, P. Morin, and R. H. Nochetto. AFEM for geometric PDE: the Laplace-Beltrami operator. In *Analysis and numerics of partial differential equations*, volume 4 of Springer INdAM Ser., pages 257–306. Springer, Milan, 2013.</footnote>
[9] H. Brakhage and P. Werner. Über das Dirichletsche Aussenraumproblem für die Helmholtzsche Schwingungsgleichung. *Arch. Math.*, 16:325–329, 1965.

[10] A. Buffa, G. Gantner, C. Giannelli, D. Praetorius, and R. Vázquez. Mathematical foundations of adaptive isogeometric analysis. *Arch. Comput. Methods Eng.*, 29:4479–4555, 2022.

[11] A. Buffa and R. Hiptmair. Regularized combined field integral equations. *Numer. Math.*, 100:1–19, 2005.

[12] A. Buffa and S. Sauter. On the acoustic single layer potential: stabilization and Fourier analysis. *SIAM J. Sci. Comput.*, 28(5):1974–1999, 2006.

[13] A. J. Burton and G. F. Miller. The application of integral equation methods to the numerical solution of some exterior boundary-value problems. *Proc. Roy. Soc. London Ser. A*, 323:201–210, 1971.

[14] C. Carstensen. An a posteriori error estimate for a first-kind integral equation. *Math. Comp.*, 66(217):139–155, 1997.

[15] C. Carstensen, M. Feischl, M. Page, and D. Praetorius. Axioms of adaptivity. *Comput. Math. Appl.*, 67(6):1195–1253, 2014.

[16] C. Carstensen and E. P. Stephan. A posteriori error estimates for boundary element methods. *Math. Comp.*, 64(210):483–500, 1995.

[17] C. Carstensen and E. P. Stephan. Adaptive boundary element methods for some first kind integral equations. *SIAM J. Numer. Anal.*, 33(6):2166–2183, 1996.

[18] S. N. Chandler-Wilde and S. Langdon. A Galerkin boundary element method for high frequency scattering by convex polygons. *SIAM J. Numer. Anal.*, 45(2):610–640, 2007.

[19] W. Dahmen, B. Faermann, I. G. Graham, W. Hackbusch, and S. A. Sauter. Inverse inequalities on non-quasi-uniform meshes and application to the Mortar element method. *Math. Comp.*, 73(247):1107–1138, 2004.

[20] A. Demlow and G. Dziuk. An adaptive finite element method for the Laplace-Beltrami operator on implicitly defined surfaces. *SIAM J. Numer. Anal.*, 45(1):421–442, 2007.

[21] S. Engleder and O. Steinbach. Modified boundary integral formulations for the Helmholtz equation. *Aust. J. Math. Anal. Appl.*, 331(1):396–407, 2007.

[22] S. Engleder and O. Steinbach. Stabilized boundary element methods for exterior Helmholtz problems. *Numer. Math.*, 110:145–160, 2008.

[23] B. Faermann. Localization of the Aronszajn-Slobodeckij norm and application to adaptive boundary elements methods. Part I. The two-dimensional case. *IMA J. Numer. Anal.*, 20(2):203–234, 2000.

[24] B. Faermann. Localization of the Aronszajn-Slobodeckij norm and application to adaptive boundary element methods. Part II. The three-dimensional case. *Numer. Math.*, 92(3):467–499, 2002.

[25] M. Feischl. Inf-sup stability implies quasi-orthogonality. *Math. Comp.*, 91(337):2059–2094, 2022.

[26] M. Feischl, T. Führer, M. Karkulik, J. M. Melenk, and D. Praetorius. Quasi-optimal convergence rates for adaptive boundary element methods with data approximation, Part I: Weakly-singular integral equation. *Calcolo*, 51(4):531–562, 2014.

[27] M. Feischl, T. Führer, M. Karkulik, J. M. Melenk, and D. Praetorius. Quasi-optimal convergence rates for adaptive boundary element methods with data approximation, Part II: Hyper-singular integral equation. *Electron. Trans. Numer. Anal.*, 44:153–176, 2015.

[28] M. Feischl, T. Führer, and D. Praetorius. Adaptive FEM with optimal convergence rates for a certain class of nonsymmetric and possibly nonlinear problems. *SIAM J. Numer. Anal.*, 52(2):601–625, 2014.

[29] M. Feischl, M. Karkulik, J. M. Melenk, and D. Praetorius. Quasi-optimal convergence rate for an adaptive boundary element method. *SIAM J. Numer. Anal.*, 51(2):1327–1348, 2013.

[30] J. Galkowski and E. Spence. Does the Helmholtz boundary element method suffer from the pollution effect? *SIAM Rev.*, 65(3):806–828, 2023.

[31] G. Gantner. *Optimal adaptivity for splines in finite and boundary element methods*. PhD thesis, TU Wien, 2017.

[32] G. Gantner and D. Praetorius. Adaptive BEM for elliptic PDE systems, part I: abstract framework, for weakly-singular integral equations. *Appl. Anal.*, 101(6):2085–2118, 2022.

[33] G. Gantner and D. Praetorius. Adaptive BEM for elliptic PDE systems, part II: Isogeometric analysis with hierarchical B-splines for weakly-singular integral equations. *Comput. Math. Appl.*, 117:74–96, 2022.

[34] G. Gantner, D. Praetorius, and S. Schimanko. IGABEM2D. Zenodo software 6282998, 2022.
[35] G. Gantner, D. Praetorius, and S. Schimanko. Stable implementation of adaptive IGABEM in 2D in MATLAB. *Comput. Methods Appl. Math.*, 22(3):563–590, 2022.

[36] T. Gantumur. Adaptive boundary element methods with convergence rates. *Numer. Math.*, 124(3):471–516, 2013.

[37] I. G. Graham, W. Hackbusch, and S. A. Sauter. Finite elements on degenerate meshes: inverse-type inequalities and applications. *IMA J. Numer. Anal.*, 25(2):379–407, 2005.

[38] R. Hiptmair. Coercive combined field integral equations. *J. Numer. Math.*, 11(2):115–134, 2003.

[39] M. Holst. Adaptive numerical treatment of elliptic systems on manifolds. *Adv. Comput. Math.*, 15(1):139–191, 2001.

[40] M. Karkulik, D. Pavlicek, and D. Praetorius. On 2D newest vertex bisection: optimality of mesh-closure and $H^1$-stability of $L_2$-projection. *Constr. Approx.*, 38(2):213–234, 2013.

[41] R. Leis. Zur Dirichletschen Randwertaufgabe des Außenraumes der Schwingungsgleichung. *Math. Z.*, 90:205–211, 1965.

[42] M. Löndorf and J. Melenk. Wavenumber-explicit hp-BEM for high frequency scattering. *SIAM J. Numer. Anal.*, 49(6):2340–2363, 2011.

[43] W. McLean. *Strongly elliptic systems and boundary integral equations*. Cambridge University Press, Cambridge, 2000.

[44] K. Mekchay, P. Morin, and R. H. Nochetto. AFEM for the Laplace-Beltrami operator on graphs: design and conditional contraction property. *Math. Comp.*, 80(274):625–648, 2011.

[45] P. E. Meury. *Stable finite element boundary element Galerkin schemes for acoustic and electromagnetic scattering*. PhD thesis, ETH Zurich, 2007.

[46] O. I. Panich. On the question of the solvability of the exterior boundary-value problems for the wave equation and Maxwell's equations. *Russian Math. Surveys*, 20(A):221–226, 1965.

[47] S. A. Sauter and C. Schwab. *Boundary element methods*. Springer, Berlin, 2011.

[48] O. Steinbach. *Numerical approximation methods for elliptic boundary value problems*. Springer, New York, 2008.

[49] R. Stevenson. Optimality of a standard adaptive finite element method. *Found. Comput. Math.*, 7(2):245–269, 2007.

[50] R. Stevenson. The completion of locally refined simplicial partitions created by bisection. *Math. Comp.*, 77(261):227–241, 2008.