Convergence analysis of a high-order Nyström integral-equation method for surface scattering problems

Abstract

In this paper we present a convergence analysis for the Nyström method proposed in [J Comput Phys 169 (1):80–110, 2001] for the solution of the combined boundary integral equation formulations of sound-soft acoustic scattering problems in three-dimensional space. This fast and efficient scheme combines FFT techniques and a polar change of variables that cancels out the kernel singularity. We establish the stability of the algorithms in the $L^2$ norm and we derive convergence estimates in both the $L^2$ and $L^\infty$ norms. In particular, our analysis establishes theoretically the previously observed super-algebraic convergence of the method in cases in which the right-hand side is smooth.

1 Introduction

In this paper we present a full convergence analysis of the high-order Nyström method introduced by Bruno and Kunyansky in [7, 8] for boundary integral equations (BIE) related to the scattering of time-harmonic acoustic waves by arbitrary (smooth) surfaces in three-dimensional space. In particular, our analysis establishes theoretically the previously observed super-algebraic convergence of the method for smooth right-hand sides. To the best of the authors’ knowledge, this proof constitutes one of the few instances of analysis of a Nyström type method for a boundary integral equation in three dimensions.

As is well known, Galerkin methods for BIE of the first kind have enjoyed thorough theoretical analyses since their inception—on the basis of ellipticity properties and discrete Fredholm theory. Compactness arguments can also be used to establish convergence of Galerkin methods for equations of the second kind. Few results exist, on the other hand, concerning convergence for three-dimensional BIE collocation methods—in which finite-element bases are used for approximation, but testing relies on point sampling. We refer to the excellent text-book [3, Chapter 9] for a brief introduction on this topic.

This state of affairs has led to the widespread perception that, being even more “discrete” than collocation schemes, Nyström methods for BIEs of the second kind with weakly singular kernels could not be easily analyzed. This paper will hopefully help dispell this view and lend additional credibility to Nyström methods—whose qualities can be very attractive in practice [7].

One of the main difficulties in the design of three-dimensional integral solvers concerns development of high-quality quadrature rules for approximation of singular integral terms. Wienert [25] constructed a singular integration rule on the basis of spherical-harmonic transforms for surfaces for which a diffeomorphism to the sphere can be constructed, (see also [13], Sect. 3.6). A Galerkin version of this approach was introduced and analyzed in [17] where, in particular, the superalgebraic convergence of the method was established. Reference [16] presents a collocation method, with corresponding analysis, for the Laplace equation, which shares the good convergence properties of the method by Wienert but which again is limited to surfaces for which a smooth mapping from the sphere is available. We emphasize that such limitations, which are highly restrictive in practice, are not imposed by the Nyström method studied in this paper.

We thus consider the numerical method [7] for the second-kind combined field integral equation associated to the problem of sound-soft scattering by a three dimensional obstacle with smooth boundary $S$. The method relies on a series of geometric constructions: (1) Representation of the surface by means of a set of $J$ overlapping smooth charts (parametrizations); (2) A smooth partition of unity subordinated to these charts which decomposes both the overall integral operators as well as the solution of the BIE as the sum of contributions defined on the parametrized patches; (3) A family of floating smooth cut-off functions $\eta _\delta : S \times S \rightarrow \mathbb R$ (see (2.10) below) that is used to isolate the singularity of the kernel function, and thus produces a splitting of the integral operator of the form $K=K_{\mathrm{reg},\delta }+K_{\mathrm{sing},\delta }$—whereby the regular operator is an integral operator with a smooth kernel, and the singular part which enjoys a reduced domain of integration but contains the weak singularity of the integral kernel.

The Nyström method under consideration is based on approximation of the integral operator by a quadrature rule which treats the regular and singular parts $K_{\mathrm{reg},\delta }$ and $K_{\mathrm{sing},\delta }$ separately. The $J$ local charts mentioned above are used to push forward $J$ uniform grids from the unit square to the surface; these grids are used for both, approximation and integration. The quadrature rule used for the regular operator $K_{\mathrm{reg},\delta }$ is based on application of the two-dimensional trapezoidal rule in the parameter space $[0,1]\times [0,1]$ for each one of the $J$ parametric domains: since the corresponding integrands are smooth with compact support strictly contained in the unit square, these trapezoidal quadratures give rise to super-algebraically accurate approximations. For the singular part, a change of variables to polar coordinates around each integration point is applied. This procedure results in a smooth integrand to which, upon necessary Fourier-based interpolations, the trapezoidal rule is applied for radial and angular integration—once again yielding super-algebraically accurate approximations.

As a result of these constructions we obtain an algorithm that evaluates the action of an approximating discrete operator on a continuous function, using only the values of the function on the selected quadrature points.

Our theoretical treatment relies on use of both existing and new analytical tools. In a first key step of the analysis the problem is re-expressed as a system of integral equations in a space of periodic functions. This is accomplished by means of yet another set of local cut-off functions, whose presence does not affect the actual numerical scheme. The use of periodic Sobolev spaces allows us to take full advantage of numerous results for approximation of Fourier series and interpolation operators [23], as well as the theory of collectively compact operators by Anselone [2]. Recasting the numerical scheme of quadrature type as a discretization method in $L^2$-type (Sobolev) spaces gives rise to a number of difficulties. In particular, for the sake of the analysis we introduce Fourier projection operators, bi-periodic trigonometric interpolation operators, and a discrete operator [11] that produces a linear combination of Dirac delta distributions on a uniform grid from a smooth function input. The convergence analysis for the operators arising from the regular part of the original boundary integral operator follows the lines of the theory [11, 15] on periodic integral equations in one variable. The final (rather technical) element of our convergence proof is a detailed analysis of the integration error arising from the numerical polar coordinate integration algorithm for products of smooth functions, “shrinking” floating cut-off functions and bi-variate trigonometric polynomials, in terms of Sobolev norms of the latter polynomials.

We point out that, for efficiency, a variety of acceleration techniques were used in [7, 8] in conjunction with the Nyström algorithm we consider. Some of these algorithmic components have been taken into account in our analysis. In the formulation considered in the present paper, for example, the computation of the singular part requires the evaluation of bivariate trigonometric polynomials that approximate the unknown at points on a polar grid. A deliberate choice of the radial quadrature nodes for this integral makes it possible to reduce this process to 1–dimensional trigonometric interpolation problems (see Sect. 2.3) on the horizontal and vertical lines of the grid. In [7], such trigonometric polynomials are then approximated by means of piecewise Hermite interpolation on dyadic grids, which can be evaluated much more rapidly than the either of the underlying trigonometric polynomials. We have analyzed the effect of these additional approximations on the full convergence of the Nyström method. The corresponding results can be found in [5, Appendix B]; briefly, upon correct parameter selections, the resulting (more efficient) method retains the super-algebraic convergence of the original approach. Additional, more sophisticated acceleration techniques which, based on use of equivalent sources and Fast Fourier Transforms, provide a means to reduce the solution cost for high-frequency problems (but on which the Nyström method itself does not depend, and whose use is not advantageous for problems of lower frequencies) were introduced in [7]. The impact of such equivalent-source acceleration methods on convergence are not considered either in the present paper or in reference [5].

This paper is structured as follows. In Sect. 2 we describe, in a compact form, the Nyström method under consideration, and we state the main convergence results of this paper. In Sect. 3 we then recast both the continuous and discrete problems as systems of equations in spaces of biperiodic functions, we derive bounds, on various norms, of the main continuous integral operators in our periodic formulation, we establish unique solvability of the continuous system of periodic integral equations as well as the equivalence of this system to the original BIE, and, finally, we state the main approximation results in the biperiodic frame: norm convergence of the discrete operators to the continuous ones together with corresponding error estimates. Finally, in Sect. 4 we present key estimates on the singular operators that appear in the biperiodic formulation, and in Sect. 5, in turn, we provide the proofs of the main results stated in Sects. 3 and 2—in that order.

We conclude our introduction with two remarks concerning notation.

Remark 1

To make a clearly visible distinction between points on the surface $S$ and coordinates for their parametrization, we use boldface letters (e.g. $\mathbf{r}$) for points on the scattering surface $S\subset \mathbb R ^3$, and underlined letters (such as $\underline{u}$) for points in $\mathbb R ^2$. The coordinates of such points will be denoted according to $\underline{u}=(u_1,u_2)$.

Remark 2

Throughout this paper the letter $C$ denotes a positive constant independent of the parameters $h$ and $\delta$ and any other variable quantities appearing in the equation. When necessary a subscripted letter $C$ is used—either to avoid confusion or to explicitly signify dependence on parameters other than $h$ and $\delta$.

2 The Nyström method

2.1 The boundary integral equation

We consider the problem of time-harmonic acoustic scattering by a sound-soft obstacle with smooth boundary $S$ in three-dimensional space:

$$\begin{aligned} \left[ \begin{array}{rcll} \varDelta U+\kappa ^2 U&{}=&{}0 &{}\quad \text{ in}\,\mathrm{ext}(S),\\ U&{}=&{}-U^\mathrm{inc}&{}\quad \text{ on}\,S,\\ \partial _r U-\mathrm{i} \kappa U&{}=&{} o(r^{-1})&{}\quad \text{ as} r\rightarrow \infty . \end{array} \right. \end{aligned}$$

(2.1)

Here $\kappa >0$ is the wave number, $U^\mathrm{inc}$ is the incident wave, $r := |\mathbf{x} |$ and $\partial _r$ denotes the radial derivative. Letting $\varPhi _\kappa$ denote the fundamental solution of the Helmholtz equation,

$$\begin{aligned} \varPhi _{\kappa } (\mathbf{r},\mathbf{r}^{\prime }):=\frac{\exp (\mathrm{i}\kappa |\mathbf{r}-\mathbf{r}^{\prime }|)}{4\pi \,|\mathbf{r}-\mathbf{r}^{\prime }|}, \end{aligned}$$

then the solution $U$ of (2.1) can be expressed as the combined (or Brakhage-Werner [4]) potential

$$\begin{aligned} U(\mathbf{r})=\int _S \frac{\partial \varPhi _{\kappa }(\mathbf{r},\mathbf{r}^{\prime })}{\partial \varvec{\nu }(\mathbf{r}^{\prime })} \psi (\mathbf{r}^{\prime })\mathrm{d}\mathbf{r}^{\prime }-\mathrm{i}\eta \int _S \varPhi _{\kappa }(\mathbf{r},\mathbf{r}^{\prime })\psi (\mathbf{r}^{\prime })\mathrm{d}\mathbf{r}^{\prime }, \end{aligned}$$

where $\frac{\partial }{\partial \varvec{\nu }}$ denotes the outward normal derivative on $S$ and $\eta > 0$ is a coupling parameter. The density $\psi$ is the unique solution of the integral equation

$$\begin{aligned} {\frac{1}{2}}\psi (\mathbf{r})+A\psi (\mathbf{r})=-U^\mathrm{inc}(\mathbf{r}) \qquad \forall \mathbf{r}\in S, \end{aligned}$$

(2.2)

where

$$\begin{aligned} A\psi (\mathbf{r})=\int _S K(\mathbf{r},\mathbf{r}^{\prime })\psi (\mathbf{r}^{\prime })\,\mathrm{d}\mathbf{r}^{\prime }:=\int _S \Big (\frac{\partial \varPhi _{\kappa }(\mathbf{r},\mathbf{r}^{\prime })}{\partial \varvec{\nu }(\mathbf{r}^{\prime })}-\mathrm{i}{\eta } \varPhi _{\kappa }(\mathbf{r},\mathbf{r}^{\prime })\Big )\psi (\mathbf{r}^{\prime })\mathrm{d}\mathbf{r}^{\prime }.\quad \end{aligned}$$

(2.3)

Various choices of the coupling parameter $\eta$ have been proposed for accuracy and numerical efficiency; see e.g. [7, 9, 12, 14]. Note that the kernel $K$ can be expressed in the form $K = K_0+K_1$, where

$$\begin{aligned} K_0(\mathbf{r},\mathbf{r^{\prime }})&:= \eta \frac{\sin (\kappa |\mathbf{r}-\mathbf{r^{\prime }}|)}{4\pi |\mathbf{r}-\mathbf{r^{\prime }}|}\nonumber \\&+\mathrm{i}\left( \kappa \cos (\kappa |\mathbf{r}-\mathbf{r^{\prime }}|) -\frac{\sin (\kappa |\mathbf{r}-\mathbf{r^{\prime }}|)}{|\mathbf{r}-\mathbf{r^{\prime }}|}\right) \frac{(\mathbf{r}-\mathbf{r^{\prime }})\cdot \varvec{\nu }({\mathbf{r^{\prime }}})}{4\pi |\mathbf{r}-\mathbf{r^{\prime }}|^2} \end{aligned}$$

(2.4)

and

$$\begin{aligned} K_1({\mathbf{r}},\mathbf{r^{\prime }})&:= -\mathrm{i}\eta \frac{\cos (\kappa |{\mathbf{r}} -\mathbf{r^{\prime }}|)}{4\pi |\mathbf{r}-\mathbf{r^{\prime }}}|\nonumber \\&- \frac{\kappa |\mathbf{r}-\mathbf{r^{\prime }}|\sin (\kappa |\mathbf{r}-\mathbf{r^{\prime }}|)+\cos (\kappa |\mathbf{r}-\mathbf{r^{\prime }}|)}{4\pi |\mathbf{r}-\mathbf{r^{\prime }}|}\frac{(\mathbf{r}-\mathbf{r^{\prime }})\cdot \varvec{\nu }(\mathbf{r^{\prime }})}{|\mathbf{r}-\mathbf{r^{\prime }}|^2}.\quad \quad \end{aligned}$$

(2.5)

Remark 3

It is easy to check that $K_0\in \mathcal C ^\infty (S\times S)$. The kernel $K_1$, on the other hand, while not $\mathcal C ^\infty (S\times S)$, can be integrated with super-algebraic accuracy by means of a combination of a polar change of variables and the trapezoidal rule; see [7] and Eq. (2.23).

In this paper we focus on the Brakhage-Werner formulation presented above; a similar analysis can be used to treat the closely related Burton-Miller formulation [10].

2.2 Geometry

The numerical method studied in this work relies heavily on the use of a system of local charts for description of the surface $S$. We will thus use a set of a number $J$ of open overlapping coordinate patches $\{{ S}^j\}_{j=1}^J$ that cover $S$,

$$\begin{aligned} S=\bigcup _{j=1}^J { S}^j, \end{aligned}$$

(2.6)

each one of which is the image of a $\mathcal C ^\infty$ parametrization

$$\begin{aligned} \begin{array}{rccl} \mathbf{r}^j:&{}D^j&{}\longrightarrow &{}S^j\\ &{}\underline{u}:=(u_1,u_2)&{}\longmapsto &{}\mathbf{r}^j(u_1,u_2), \end{array} \end{aligned}$$

where $\overline{D^j}\subset I_2:=(0,1)\times (0,1)$. We assume that $\mathbf{r}^j$ can be extended to a $\mathcal C ^\infty$ bijective diffeomorphism between $\overline{D^j}$ and $\overline{S^j}$ so that, in particular, the Jacobians

$$\begin{aligned} a^j(\underline{u}):= \bigg |\frac{\partial \mathbf{r}^j(u_1,u_2)}{\partial u_1}\times \frac{\partial \mathbf{r}^j(u_1,u_2)}{\partial u_2}\bigg |>{ 0},\qquad j=1,\ldots ,J \end{aligned}$$

(2.7)

($a^j:D^j \rightarrow \mathbb R$) are $\mathcal C ^\infty$ functions of $\underline{u}$.

The method requires explicit use of a smooth partition of unity $\{\omega ^j\}_{j=1}^J\subset \mathcal{C }^\infty (S)$ on $S$, subordinated to the covering (2.6), that is

$$\begin{aligned} \omega ^j \ge 0, \qquad \mathrm{supp}\,\omega ^j \subset {S}^j\quad \text{ and}\quad \sum _{j=1}^J \omega ^j \equiv 1. \end{aligned}$$

(see Fig. 1).

Fig. 1

Sketch of the geometric layout of a single chart and its associated cut-off function

The hypotheses on availability of local charts and a partition of unity is not restrictive in practice: such parameterizations can be constructed for smooth arbitrary geometries (see e.g. [6]). We will also assume that the boundary of $\mathrm{supp}\,\omega ^j$ is the finite union of Lipschitz arcs, a restriction which, again, is easy to accommodate [6].

For any $\delta >0$ and $\mathbf{r}\in S$ we let

$$\begin{aligned} B({\mathbf{r}},\delta ):=\{ {\mathbf{r}^{\prime }}\in S\, :\, |\mathbf{r}^{\prime }-\mathbf{r}|<\delta \} \end{aligned}$$

and, selecting parameters $0< \epsilon _0<\epsilon _1\le 1$ that will otherwise be fixed throughout this paper, we define

$$\begin{aligned} S^j_{\delta }:=\overline{\bigcup _{\mathbf{r}\in \mathrm{supp}\,\omega ^j} B(\mathbf{r},2\epsilon _1\delta )}=\bigcup _{\mathbf{r}\in \mathrm{supp}\omega ^j}\overline{B(\mathbf{r},2\epsilon _1\delta )}. \end{aligned}$$

(2.8)

Clearly there exists $\delta _0>0$ such that for all $j$

$$\begin{aligned} \overline{B(\mathbf{r},\epsilon _1\delta _0)}\cap S^j_{2\delta _0}=\emptyset \qquad \forall \mathbf{r}\in S\setminus S^j \end{aligned}$$

(2.9)

—as it can be checked by considering a pair of points that realize the distance between the boundaries of $\mathrm{supp}\,\omega ^j$ and $S^j$. In particular, this implies that $S^j_{2\delta _0}\subset S^j$. The final element in our geometric constructions is a $\mathcal C ^\infty$ function $\upsilon :S\times [0,\infty )\rightarrow [0,1]$ such that

$$\begin{aligned} \begin{array}{l} \upsilon (\mathbf{r},\,\cdot \,) \equiv 1 \quad \text{ in}\,\,\, [0,\epsilon _0\delta _0]\\ \upsilon (\mathbf{r},\,\cdot \,) \equiv 0 \quad \text{ in}\,\,\, [\epsilon _1\delta _0,\infty ) \end{array} \qquad \forall \mathbf{r}\in S. \end{aligned}$$

Given $\delta \in (0,\delta _0]$ we now define the functions $\eta _\delta :S \times S \rightarrow [0,1]$

$$\begin{aligned} \eta _\delta ({\mathbf{r}},{\mathbf{r^{\prime }}}):=\upsilon \Big ({\mathbf{r}},\frac{\delta _0|{\mathbf{r^{\prime }}}-{\mathbf{r}}|}{\delta }\Big ). \end{aligned}$$

(2.10)

Clearly

$$\begin{aligned} \mathrm{supp}\,\eta _\delta (\mathbf{r},\,\cdot \,) \subset \overline{B(\mathbf{r},\epsilon _1\delta )}, \qquad \forall \mathbf{r}\in S. \end{aligned}$$

(2.11)

In view of the definition of the sets $S_{\delta }^j$ and (2.9) we also have

$$\begin{aligned} \mathrm{supp}\,\eta _{\delta }(\mathbf{r},\cdot )\cap S^j_{3\delta /2}=\emptyset \qquad \forall \mathbf{r}\in S\setminus S^j_{2\delta }. \end{aligned}$$

(2.12)

Given $\delta$ and considering the decomposition $K=K_0+K_1$ in term of the kernel functions given in Eqs. (2.4) and (2.5), we define the regular part of the kernel of the integral operator (2.3) as

$$\begin{aligned} K_{\mathrm{reg},\delta }(\mathbf{r},\mathbf{r}^{\prime }):= K_0(\mathbf{r},\mathbf{r}^{\prime })+ \big (1-\eta _\delta (\mathbf{r}, \mathbf{r}^{\prime })\big ) \, K_1(\mathbf{r},\mathbf{r}^{\prime }). \end{aligned}$$

(2.13)

Clearly $K_{\mathrm{reg},\delta }\in \mathcal C ^\infty (S\times S)$. The remainder is the singular part of the kernel,

$$\begin{aligned} K_{\mathrm{sing},\delta }:= K - K_{\mathrm{reg},\delta }= \eta _\delta \, K_1, \end{aligned}$$

(2.14)

which, like the kernel $K_1$, can be integrated accurately by means of a polar change of variables; see Remark 3. The parameter $\delta$, which controls the support of the kernel $K_\mathrm{sing,\delta }(\mathbf{r},\,\cdot \,)$, plays an essential role in both, the performance of the algorithm and its theoretical analysis; see Remark 6 for details.

We next introduce

$$\begin{aligned} K^{ij}(\underline{u},\underline{v}):= K(\mathbf{r}^i(\underline{u}),\mathbf{r}^j(\underline{v})) \,a^j(\underline{v}), \end{aligned}$$

[see Eq, (2.7)] with corresponding definitions of $K^{ij}_{\mathrm{sing},\delta }$ and $K^{ij}_{\mathrm{reg},\delta }$ (cf. Eqs. (2.13) and (2.14)). Noting that $K^{ij}, K^{ij}_{\mathrm{sing},\delta }$ and $K^{ij}_{\mathrm{reg},\delta }$, are defined in $D^i \times D^j$, we extend these functions by zero (possibly thereby introducing discontinuities on the boundary of $D^j$) to the full product of squares $I_2\times I_2$. Clearly,

$$\begin{aligned} \int _S K(\mathbf{r}^i(\underline{u}),\mathbf{r}^{\prime })\xi (\mathbf{r}^{\prime })\,\mathrm{d}\mathbf{r}^{\prime }=\sum _{j=1}^J \int _{D^j} K^{ij}(\underline{u},\underline{v})\omega ^j(\mathbf{r}^j(\underline{v}))\xi (\mathbf{r}^j(\underline{v}))\,\mathrm{d}\underline{v} \qquad \forall \underline{u} \in D^i. \end{aligned}$$

Therefore, if $\psi$ is the exact solution of (2.2), then $\psi ^i:=\psi \circ \mathbf{r}^i :D^i\rightarrow \mathbb C$ ($i=1,\ldots ,J$), is a solution of the system

$$\begin{aligned}&\frac{1}{2}\psi ^i(\underline{u})+ \sum _{j=1}^J \int _{I_2} K^{ij}_{\mathrm{reg},\delta }(\underline{u},\underline{v}) \omega ^j(\mathbf{r}^j(\underline{v}))\psi ^j(\underline{v})\,\mathrm{d}\underline{v}\nonumber \\&\quad +\sum _{j=1}^J \int _{I_2} K^{ij}_{\mathrm{sing},\delta }(\underline{u},\underline{v}) \omega ^j(\mathbf{r}^j(\underline{v})) \psi ^j(\underline{v})\,\mathrm{d}\underline{v}= - U^\mathrm{inc}(\mathbf{r}^i(\underline{u}))\nonumber \\&\qquad \forall \underline{u} \in D^i,\quad i=1,\ldots , J.\qquad \qquad \end{aligned}$$

(2.15)

Remark 4

Note that the functions $\psi ^j$ that appear in the integrals over $I_2$ in Eq. (2.15) are only defined in $D^j$. However, they are multiplied by the cutoff function $\omega ^j\circ \mathbf{r}^j$ which vanishes outside $D^j$ and thus provides a natural extension for the product throughout $I_2$.

The solution of this system is used in Sect. 3 to reconstruct the solution of the original problem (2.2).

We can clearly distinguish two different types of integral operators in (2.15), namely, integral operators with smooth and singular kernels. The discretizations of these operators are produced, accordingly, by means of two different strategies—as discussed in the following section.

2.3 Discretization of the integral operators in Eq. (2.15)

For each patch $S^j$, we select a positive integer $N_j$, we let $h_j:=1/N_j$, and we introduce the grid points

$$\begin{aligned} \underline{x}^j_{\underline{m}} := h_j \underline{m}= h_j\, (m_1,m_2), \qquad 0\le m_1,m_2 \le N_j-1. \end{aligned}$$

We assume that these grids are quasiuniform: letting

$$\begin{aligned} N:=\min _{j=1,\ldots ,J} N_j, \qquad h:=1/N=\max _{j=1,\ldots ,J} h_j, \end{aligned}$$

we assume there exists $c>0$ such that

$$\begin{aligned} h<c h_j,\qquad j=1,\ldots , J. \end{aligned}$$

(2.16)

This is not a restrictive assumption in view of the assumed smoothness of $S$ and $U^\mathrm{inc}$, and, therefore, of the solution $(\psi ^j)$.

The algorithm under consideration is a Nyström method which produces point-wise values of the unknown $(\psi ^j)$ at the discretization points $\underline{x}_{\underline{m}}^j$. The quadrature rules that ultimately define the method are described in what follows; for notational simplicity our approximate quadrature formulae use the set

$$\begin{aligned} \mathbb Z _{N_j}:=\big \{\underline{m}:=(m_1,m_2)\in \mathbb Z ^2 : 0\le m_1,m_2\le N_j-1\big \} \end{aligned}$$

(2.17)

of two-dimensional summation indices $\underline{m}$.

The approximate integration method used by the algorithm to treat the regular portion of the kernel is, simply, based on the trapezoidal rule,

$$\begin{aligned} \int _{I_2} K^{ij}_{\mathrm{reg},\delta }(\,\cdot \,, \underline{v}) \xi (\underline{v})\mathrm{d}\underline{v} \approx h_j^2\sum _{\underline{m}\in \mathbb Z _{N_j}} \! K^{ij}_{\mathrm{reg},\delta }(\,\cdot \,, \underline{x}^j_{\underline{m}}) \xi (\underline{x}^j_{\underline{m}}), \end{aligned}$$

(2.18)

with the convention that $K^{ij}_{\mathrm{reg},\delta }(\,\cdot \,, \underline{v})\equiv 0$ for $\underline{v}\not \in D^j$. Notice that for a $\mathcal C ^\infty$ function $\xi$ compactly supported in

$$\begin{aligned} \varOmega ^j:=(\mathbf{r}^{j})^{-1}(\mathrm{supp}\,\omega ^j)\subset D^j, \end{aligned}$$

(2.19)

(such as $\xi =(\omega ^j\psi )\circ \mathbf{r}^j$ with $\psi \in \mathcal{C }^\infty (S)$) the rule embodied in Eq. (2.18) gives rise to super-algebraic convergence.

Fig. 2

The overlapping of two charts and the associated domains

To approximate integrals of the form

$$\begin{aligned} \int _{I_2} K^{ij}_{\mathrm{sing},\delta }(\underline{u}, \underline{v}) \xi (\underline{v}) \mathrm{d} \underline{v}, \qquad \xi \in \mathcal C ^\infty (I_2), \quad \mathrm{supp}\,\xi \subseteq \varOmega ^j, \end{aligned}$$

(2.20)

that include the singular part of the kernel, the algorithm uses polar coordinates around the singularity. To properly account for contributions arising from various regions delineated by local charts and overlaps we let

$$\begin{aligned} \begin{array}{c} \underline{r}^{ji}:=(\mathbf{r}^{j})^{-1} \circ \mathbf{r}^i:(\mathbf{r}^i)^{-1}(S^i\cap S^j)\subset D^i\rightarrow D^j \\ S^{ij}_{\delta }:=S^i \cap S^j_{\delta } \subset S^i\cap S^j,\qquad \varOmega ^{ij}_{\delta }:=(\mathbf{r}^{i})^{-1}(S^i \cap S^j_{\delta })=(\mathbf{r}^i)^{-1} (S^{ij}_\delta ); \end{array}\nonumber \\ \end{aligned}$$

(2.21)

see Fig. 2 and Eq. (2.8). A close inspection of the definition of the sets $S^j_\delta$ shows that $\overline{B(\mathbf{r},\epsilon _1\delta )}\cap \mathrm{supp}\,\omega ^j=\emptyset$ if $\mathbf{r}\not \in S^j_{\delta /2}$. Since, additionally,

$$\begin{aligned} \mathrm{supp}\,K^{ij}_{{\mathrm{sing}},\delta }(\underline{u}, \,\cdot \,) \subset \mathrm{supp}\,\eta _\delta (\mathbf{r}^i(\underline{u}), \mathbf{r}^j(\,\cdot \,))\subset (\mathbf{r}^{j})^{-1}\big ( \overline{B(\mathbf{r}^i(\underline{u}),\epsilon _1\delta )}\cap S^j\big ), \end{aligned}$$

we conclude that the integral (2.20) with $\mathrm{supp}\,\xi \subseteq \varOmega ^j$, vanishes for all $\underline{u}$ such that $\mathbf{r}^i(\underline{u})\not \in S^j_{\delta /2}$, i.e., for $\underline{u}\not \in \varOmega _{\delta /2}^{ij}$. The algorithm therefore only produces approximations of (2.20) for $\underline{u}\in \varOmega _{\delta /2}^{ij}$. To do this consider the relations

$$\begin{aligned}&\int _\mathbb{R ^2} K^{ij}_{\mathrm{sing},\delta }( \underline{u}, \underline{v}) \xi (\underline{v}) \mathrm{d} \underline{v}\!\nonumber \\&\quad = \int _\mathbb{R ^2} K_{\mathrm{sing},\delta }^{ij}\left( \underline{u}, \underline{r}^{ji}(\underline{u}) + \underline{w}\right) \xi (\underline{r}^{ji}(\underline{u})+\underline{w})\, \mathrm{d} \underline{w}\nonumber \\&\quad =\frac{1}{2}\int _{0}^{2\pi }\bigg (\int _{-\infty }^\infty |\rho |\,K_{\mathrm{sing},\delta }^{ij}\left( \underline{u}, \underline{r}^{ji}(\underline{u}) + \rho \underline{e}(\theta )\right) \,\xi (\underline{r}^{ji}(\underline{u})+ \rho \underline{e}(\theta ))\,{\mathrm{d}} \rho \bigg ) \mathrm{d} \theta \qquad \quad \end{aligned}$$

(2.22)

where $\underline{e}(\theta ):=(\cos \theta ,\sin \theta )$ and note that the function

$$\begin{aligned} (\rho ,\theta ) \mapsto |\rho | K_{\mathrm{sing},\delta }^{ij}\left( \underline{u}^i,\, \underline{r}^{ji}( \underline{u} ) + \rho \underline{e}(\theta )\right) \end{aligned}$$

(2.23)

is smooth (as shown in Sect. 4, cf. [7, 8]), $2\pi$-periodic in $\theta$ and compactly supported as a function of the variable $\rho$. In what follows we temporarily let $\underline{z}=(z_1,z_2):=\underline{r}^{ji}(\underline{u})$, so that the dependences on the patch indices $(i,j)$ and parametric coordinate $\underline{u}$ are not displayed.

The integral in the angular variable is approximated using a trapezoidal rule on a uniform partition of $[0,2\pi ]$ in $\varTheta _j$ subintervals of length $k_j=2\pi /\varTheta _j$. Therefore, we consider the angles

$$\begin{aligned} \theta _p^j:= p\, k_j, \qquad p=0,\ldots , \varTheta _j-1 \end{aligned}$$

and approximate (2.22) by

$$\begin{aligned} \frac{ k_j}{2}\sum _{p=0}^{\varTheta _j-1} \Bigg ( \int _{-\infty }^\infty K^{ij}_{\mathrm{sing},\delta } (\underline{u},\underline{z}+ \rho \,\underline{e}(\theta _p^j) )\, |\rho |\, \xi (\underline{z}+\rho \, \underline{e}(\theta _p^j) ) \mathrm{d} \rho \Bigg ). \end{aligned}$$

(2.24)

Recall that a uniform Cartesian grid in $I_2$ with mesh length $h_j$ has been introduced in the approximation of the regular part (2.18). We will now use points of this grid to approximate the radial integrals in (2.24).

If $|\cos \theta _p^j|\ge \sqrt{2}/2$, that is, modulo $2\pi$, $\theta _p^j \in [-\pi /4,\pi /4] \cup [3\pi /4, 5\pi /4]$, we look for the intersections

$$\begin{aligned} \{\underline{x}= \underline{z}+\rho \, e(\theta _p^j)\, :\, \rho \in \mathbb R \}\cap \{ \underline{x}= ( q\, h_j,\xi )\, :\, \xi \in \mathbb R \}, \qquad q=0,\ldots , N_j-1, \end{aligned}$$

i.e., the intersections of the double rays stemming from $\underline{z}$ with angle $\pm \theta _p^j$ with the vertical lines of the uniform grid: the points of intersections are given by

$$\begin{aligned} \big ( q\, h_j, z_2+( q\,h_j-z_1) \tan \theta _p^j\big ),\qquad q=0,1,\ldots , N_j-1; \end{aligned}$$

(2.25)

clearly, the distance between two consecutive points is $h_j/|\cos \theta _p^j|\le \sqrt{2}h_j.$ Alternatively we can express these points in the form

$$\begin{aligned} \underline{z} + \rho _q^{h_j} (\underline{z},\theta _p^j)\, \underline{e}(\theta _p^j), \qquad \text{ where}\qquad \rho _q^h (\underline{z},\theta ):= \frac{q\,h-z_1}{\cos \theta }. \end{aligned}$$

For these values of $\theta _p^j$, the corresponding integral in (2.24) is approximated by

$$\begin{aligned} h_j\,c(\theta _p^j) \sum _{q=0}^{N_j-1} K^{ij}_{\mathrm{sing},\delta } \big (\underline{u},\underline{z}+ \rho _q^{ h_j} (\underline{z},\theta _p^j)\,\underline{e}(\theta _p^j)\big ) \, |\rho _q^{h_j} (\underline{z},\theta _p^j)| \, \xi \left( \underline{z}+\rho _q^{ h_j} (\underline{z},\theta _p^j)\, \underline{e}(\theta _p^j) \right) ,\nonumber \\ \end{aligned}$$

(2.26)

where $c(\theta )=1/|\cos \theta |$. For the remaining angles, intersections are computed with the corresponding horizontal lines:

$$\begin{aligned} \{\underline{x}= \underline{z}+\rho \, \underline{e}(\theta _p^j)\, :\, \rho \in \mathbb R \}\cap \{ \underline{x}= (\xi , q\, h_j)\, :\, \xi \in \mathbb R \}, \qquad q=0,\ldots , N_j-1.\nonumber \\ \end{aligned}$$

(2.27)

(We show in Fig. 3 a picture illustrating this procedure).

Fig. 3

Points used in the integration of the singular part

The quadrature points are deliberately selected to lie on lines with one of the coordinates equal to $q \, h_j$. As discussed in the introduction, such a selection was introduced in [7] to enable fast interpolation of the function $\xi$ on the radial quadrature points [that is, to produce a fast algorithm for evaluation of the operator $R_{N_j}$ in (2.33)].

Using the angle-dependent weight

$$\begin{aligned} c(\theta ):= \min \Big \{\frac{1}{|\cos \theta |},\frac{1}{|\sin \theta |}\Big \}, \end{aligned}$$

(2.28)

as well as the nodal points radii

$$\begin{aligned} \rho ^h_q(\underline{z},\theta ) := \left\{ \begin{array}{ll} \frac{qh-z_1}{|\cos \theta |}, &{}\quad \text{ if}\,|\cos \theta | \sqrt{2}/2,\\ \frac{qh-z_2}{|\sin \theta |}, &{}\quad \text{ otherwise}, \end{array}\right. \end{aligned}$$

expression (2.26) provides approximation of all the integrals in (2.24).

In sum, we define our discrete singular operator (which depends on $h_j$ and $k_j$, or, equivalently, on $N_j$ and $\varTheta _j$) by

$$\begin{aligned} (\mathrm{L}^{ij}_{\delta ,h}\xi )(\underline{u}):= \left\{ \begin{array}{ll} \displaystyle \frac{h_j k_j}{2} \!\!\sum _{p=0}^{\varTheta _j-1} \sum _{q=0}^{N_j-1}\!\! c(\theta _p^j)\bigg (K^{ij}_\mathrm{sing,\delta } \Big (\underline{u},\underline{r}^{ji}(\underline{u})+ \rho _q^{h_j} (\underline{r}^{ji}(\underline{u}),\theta _p^j)\,\underline{e}(\theta _p^j) \Big )&{}\\ \displaystyle \times |\rho _q^{h_j} (\underline{r}^{ji}(\underline{u}),\theta _p^j)|\, \xi \Big (\underline{r}^{ji}(\underline{u})+\rho _q^{h_j} ( \underline{r}^{ji}(\underline{u}),\theta _p^j)\, \underline{e}(\theta _p^j) \Big )\bigg ),&{}\text{ if}\ \underline{u} \in \varOmega ^{ij}_{\delta /2},\nonumber \\ 0,&{} \text{ otherwise}, \end{array}\right. \nonumber \\ \end{aligned}$$

(2.29)

where $\varOmega ^{ij}_{\delta /2}=(\mathbf{r}^i)^{-1}(S^{ij}_{\delta /2})= (\mathbf{r}^i)^{-1}(S^{i}\cap S^{j}_{\delta /2})$.In what follows we use the scaling

$$\begin{aligned} \varTheta _j\approx N_j^{1+\alpha },\qquad \alpha >0, \end{aligned}$$

(2.30)

i.e., we assume that there exist $C,c>0$ such that for all $j$, $cN_j^{1+\alpha }\le \varTheta _j\le C N_j^{1+\alpha }$. We point out that $\rho _q^h(\underline{z},\,\cdot \,)$ is discontinuous at $\{\pi /4,3\pi /4,5\pi /4, 7\pi /4\}$ but both side limits exist. It is immaterial which value (the right or left limit or even the average of both) is set as value to $\rho _q^{h}(\underline{z},\,\cdot \,)$ at these points.

2.4 The numerical method

We are now in a position to lay down a fully discrete version of Eq. (2.15). The unknowns are taken to be approximate values

$$\begin{aligned} \varphi ^j_{\underline{\ell }}, \qquad \underline{\ell }\in \varOmega ^j_h:=\big \{\underline{m}\in \mathbb Z _{N_j}\,:\,\mathbf{r}^j(\underline{x} ^j_{\underline{m}})\in \mathrm{int}(\mathrm{supp}\,\omega ^j)\big \}, \qquad j=1,\ldots ,J\nonumber \\ \end{aligned}$$

(2.31)

of the true function values $\psi ^j(\underline{x}^j_{\underline{\ell }})$. In what follows we denote by $\varvec{\varphi }^j\in \mathbb C (\varOmega ^j_h)$ the array whose $\underline{\ell }= (\ell _1,\ell _2)$ entry is $\varphi ^j_{\underline{\ell }}$ for all indices $\underline{\ell }\in \varOmega ^j_h$. We consider the set of indices

$$\begin{aligned} \mathbb Z _{N_j}^* := \big \{\underline{m}:=(m_1,m_2)\in \mathbb Z ^2\, :\, -N_j/2\le m_1,m_2 < N_j/2\big \}, \end{aligned}$$

(2.32)

(note that the cardinality of $\displaystyle \mathbb Z _{N_j}^*$ is $\# \displaystyle \mathbb Z ^*_{N_j}=\# \mathbb Z _{N_j}=N_j^2$), the space of trigonometric polynomials

$$\begin{aligned} \mathbb T _{N_j}:= \mathrm{span}\, \,\big \{ \exp ( 2\pi \mathrm{i}\, \underline{m}\cdot \underline{u})\, \,:\,\, \underline{m} \in \mathbb Z _{N_j}^*\big \}, \end{aligned}$$

the interpolation operator $\mathrm{R}_{N_j}: \mathbb C (\varOmega ^j_h) \rightarrow \mathbb T _{N_j}$ defined by the equations

$$\begin{aligned} {\mathrm{R}}_{N_j}\varvec{\varphi }(\underline{x}^j_{\underline{\ell }})= \left\{ \begin{array}{ll} \varphi ^j_{\underline{\ell }}, &{}\quad \underline{\ell }\in \varOmega ^j_h,\\ 0, &{}\quad \underline{\ell }\in \mathbb Z _{N_j}\setminus \varOmega ^j_h,\end{array}\right. \end{aligned}$$

and the vectors ${\varvec{\omega }}^j\in \mathbb C (\varOmega ^j_h)$ with components $\omega ^j_{\underline{\ell }}:=\omega ^j(\mathbf{r}^j( \underline{x}^j_{\underline{\ell }}))$. Using these notations, the discrete Nyström equations for the system (2.15) are given by

$$\begin{aligned}&{\frac{1}{2}} \varphi ^{i}_{\underline{\ell }} + \sum _{j=1}^J h_j^2 \sum _{\underline{m}\in \varOmega ^j_h} K^{ij}_{\mathrm{reg},\delta }(\underline{x}^i_{\underline{\ell }}, \underline{x}^j_{\underline{m}}) \omega ^j_{\underline{m}}\, \varphi ^j_{\underline{m}} + \sum _{j=1}^J \Big (\mathrm{L}^{ij}_{\delta ,h} {\mathrm{R}}_{N_j}(\varvec{\omega }^j \odot \varvec{\varphi }^j) \Big ) (\underline{x}_{\underline{\ell }}^i)\nonumber \\&\quad =-U^\mathrm{inc}(\mathbf{r}^i(\underline{x}^i_{\underline{\ell }})) \end{aligned}$$

(2.33)

for $\underline{\ell }\in \mathbb C (\varOmega ^j_h)$ and $j=1,\ldots ,N$, where the binary operator $\odot$ denotes componentwise product of arrays: for example, $\varvec{\omega }^j \odot \varvec{\varphi }^j$ is the array whose $\underline{\ell }$-th entry is given by $\omega ^j_{\underline{\ell }}\varphi ^j_{\underline{\ell }}$. Once these point values $\varphi ^j_{\underline{\ell }}$ have been computed, the Nyström method provides the reconstruction formula

$$\begin{aligned} \psi ^{i}_h =-2\sum _{j=1}^J \bigg ( h_j^2 \sum _{\underline{m}\in \varOmega ^j_h} K^{ij}_{\mathrm{reg},\delta }(\,\cdot \,, \underline{x}^j_{\underline{m}}) \omega ^j_{\underline{m}}\, \varphi ^j_{\underline{m}} + \mathrm{L}^{ij}_{\delta ,h} \mathrm{R}_{N_j}(\varvec{\omega }^j\odot \varvec{\varphi }^j) \bigg )-2U^\mathrm{inc}\circ \mathbf{r}^i, \nonumber \\ \end{aligned}$$

(2.34)

($i=1,\ldots ,J$). With these definitions, we clearly have

$$\begin{aligned} \psi ^j_h(\underline{x}^j_{\underline{\ell }})=\varphi ^j_{\underline{\ell }} \qquad \underline{\ell }\in \varOmega ^j_h. \end{aligned}$$

Finally, the $J$ parameter-space continuous functions $\psi ^j_h:I_2\rightarrow \mathbb C$ can be assembled into a single continuous approximate solution $\psi _h(\mathbf{r})$ defined on $S$:

$$\begin{aligned} \psi _h(\mathbf{r})\!:=\!\sum _{j \in \mathcal I (\mathbf{r})} \omega ^j(\mathbf{r}) \psi ^j_h((\mathbf{r}^j)^{-1}(\mathbf{r})), \qquad \text{ where}\,\mathcal I (\mathbf{r})\!:=\!\{ j\,:\, \mathbf{r}\in \mathrm{supp}\,\omega ^j \}.\qquad \quad \quad \end{aligned}$$

(2.35)

The main convergence result of this paper can now be stated; its proof is given in Sect. 3 through 5. The regularity estimates given in the statement of the theorem are expressed in terms of standard Sobolev norms on the surface $S$ (see for instance [1, 22, 24] or any standard text on Sobolev spaces).

Remark 5

The notation $\delta \approx h^\beta$ means that

$$\begin{aligned} c h^\beta \le \delta \le Ch^\beta , \end{aligned}$$

for some constants $c>0$ and $C>0$, independent of $h$ and $\beta$. Hereafter, the parameter $\beta \in (0,1)$ will be taken to be fixed. Dependence of constants on this parameter will not be shown explicitly.

Theorem 1

Let $\psi$ be the solution of (2.2), and assume $\delta \approx h^\beta$ with $\beta \in (0,1)$. Then for sufficiently large values of $N$ the Nyström equations (2.33) admit a unique solution. Letting $\psi _h$ denote the reconstructed function on $S$ as defined by (2.34) and (2.35), we have the error estimate

$$\begin{aligned} \Vert \psi _h-\psi \Vert _{L^2(S)}\le C_{t} \big (|\log h|^{\frac{3}{2}}h^{t(1-\beta )}+h^{t-1}\big )\Vert \psi \Vert _{H^t(S)}, \end{aligned}$$

(2.36)

for all $t\ge 2$. Finally, for all $\varepsilon >0$ and $t\ge 2$

$$\begin{aligned} \Vert \psi _h-\psi \Vert _{L^\infty (S)}\le C_{t,\varepsilon } \big (h^{t(1-\beta )-(1+\varepsilon ) \beta }+h^{t-1}\big )\Vert \psi \Vert _{H^{t}(S)}. \end{aligned}$$

(2.37)

Theorem 1 tells us that for a smooth surface $S$ and a $\mathcal C ^\infty$ right-hand side $U^\mathrm{inc}$ (from which it follows that $\psi \in \mathcal C ^\infty (S)$) the Nyström algorithm under consideration converges super-algebraically fast.

Remark 6

As mentioned in Sect. 2.2, the parameter $\delta \approx h^\beta$ plays central roles in both the theory and the actual performance of the the algorithm under consideration. With regards to the latter we briefly mention here that use of the floating cut-off function (2.10), whose support is controlled by the parameter $\delta$, helps restrict the use of the costly polar integration scheme to a small region around the singular point (thus reducing the overall computing time required by the algorithm), and, further, it enables acceleration via a sparse, parallel-face FFT-based equivalent-source technique; see [7] for details. (In particular we note that the value $\beta = 1/3$ is used in [7, 8] for optimal speed of the equivalent-source accelerated Nyström method.) The parameter $\delta \approx h^\beta$ also has a significant impact on our theoretical treatment. Indeed, one of the most delicate points in our stability analysis concerns the convergence in norm of a discrete singular operator to a corresponding continuous singular operator. One of the terms in our estimate of the norm of the difference of these operators (Eq. (5.18) in Proposition 9) is bounded by $|\log h|^{3/2} h^{\beta /2}$. By taking $\beta >0$ we ensure that this terms also tends to zero, from which the desired convergence in norm results.

3 Biperiodic framework

3.1 Continuous equations

The analysis of the method will be carried out by recasting the system (2.15) as a system of periodic integral equations with all unknown functions defined on the unit square $I_2$. To introduce our periodic formulation we utilize a second family of cut-off functions, $\widetilde{\omega }^j_\delta \in \mathcal{C }^\infty (I_2)$, that depends on $\delta$ and satisfies the following assumptions:

$$\begin{aligned}&\mathrm{supp}\,\widetilde{\omega }^j_\delta \subset \varOmega ^{jj}_{3\delta /2}=(\mathbf{r}^j)^{-1}(S^{j}_{3\delta /2}),\end{aligned}$$

(3.1a)

$$\begin{aligned}&\widetilde{\omega }_\delta ^j\ \equiv \ 1 \quad \text{ in} \quad \varOmega ^{jj}_{\delta }=(\mathbf{r}^j)^{-1}(S^j_\delta )\supset \varOmega ^j,\end{aligned}$$

(3.1b)

$$\begin{aligned}&\Vert \partial _{\underline{\alpha }}\widetilde{\omega }_\delta ^j \Vert _{\infty ,I_2}\le C_{\underline{\alpha }}\delta ^{-|\underline{\alpha }|}\qquad \forall \underline{\alpha }\ge \underline{0}, \end{aligned}$$

(3.1c)

where for a given non-negative bi-index $\underline{\alpha }=(\alpha _1,\alpha _2)$,

$$\begin{aligned} \partial _{\underline{\alpha }}\xi (\underline{u})=\partial ^{\alpha _1}_{ {u}_1 } \partial ^ { \alpha _2 } _ { {u}_2}\xi (u_1,u_2) \end{aligned}$$

and $\Vert \cdot \Vert _{\infty ,I_2}$ is the $L^\infty (I_2)$ norm. We will use the characteristic function

$$\begin{aligned} \widetilde{\omega }_0^j(\underline{u}):=\left\{ \begin{array}{ll} 1&{}\quad \text{ if} \quad \underline{u}\in \varOmega ^j,\\ 0&{}\quad \text{ otherwise},\\ \end{array} \right. \end{aligned}$$

(3.2)

which can be viewed as the limit of $\widetilde{\omega }^j_{\delta }$ as $\delta \rightarrow 0$.

Using these functions, we define the following periodic integral operators

$$\begin{aligned} \mathrm{A}^{ij}_{\mathrm{reg},\delta }\xi&:= (\omega ^i\circ \mathbf{r}^i) \int _{I_2} K^{ij}_{\mathrm{reg},\delta }(\,\cdot \,,\underline{u})\widetilde{\omega }^j_\delta (\underline{u})\xi (\underline{u}) \mathrm{d}\underline{u},\end{aligned}$$

(3.3)

$$\begin{aligned} \mathrm{A}^{ij}_{\mathrm{sing},\delta }\xi&:= (\omega ^i\circ \mathbf{r}^i) \int _{I_2} K^{ij}_\mathrm{sing,\delta }(\,\cdot \,,\underline{u})\widetilde{\omega }^j_\delta (\underline{u})\xi (\underline{u})\mathrm{d}\underline{u},\end{aligned}$$

(3.4)

$$\begin{aligned} \mathrm{A}^{ij}_\delta \xi&:= (\mathrm{A}^{ij}_{\mathrm{reg},\delta }+ \mathrm{A}^{ij}_\mathrm{sing,\delta })\xi , \end{aligned}$$

(3.5)

as well as the right-hand sides $U^i:=-(\omega ^iU^\mathrm{inc})\circ \mathbf{r}^i$, properly extended by zero to the full unit square $I_2$. If $\underline{u} \in I_2 \setminus \varOmega ^{ij}_{2\delta }$, then $\mathbf{r}^i(\underline{u}) \not \in S^{j}_{2\delta }$ and

$$\begin{aligned} \mathrm{supp}\,(K^{ij}_{\mathrm{sing},\delta }(\underline{u},\,\cdot \,) \,\widetilde{\omega }^j_\delta )&\subset \mathrm{supp}\, \eta _\delta (\mathbf{r}^i(\underline{u}), \mathbf{r}^j(\,\cdot \,) )\cap \mathrm{supp}\widetilde{\omega }^j_\delta \\&\subset \mathrm{supp}\,\eta _\delta (\mathbf{r}^i(\underline{u}), \mathbf{r}^j(\,\cdot \,) ) \cap \varOmega ^{jj}_{3\delta /2}=\emptyset \end{aligned}$$

by (2.12) and (3.1a), which implies that

$$\begin{aligned} (\mathrm{A}^{ij}_{\mathrm{sing},\delta }\xi )(\underline{u})=0\qquad \forall \underline{u}\in I_2\setminus \varOmega ^{ij}_{2\delta }. \end{aligned}$$

(3.6)

Letting $\psi ^j:=\psi \circ \mathbf{r}^j$ be as in Eq. (2.15) and since $\widetilde{\omega }_\delta ^j\,(\omega ^j\circ \mathbf{r}^j)\equiv \omega ^j\circ \mathbf{r}^j$ (see (3.1b) and note that $\mathrm{supp}\,(\omega ^j\circ \mathbf{r}^j)\subset \varOmega ^{jj}_\delta$), we see that the functions

$$\begin{aligned} \phi ^j:=(\omega ^j\psi )\circ \mathbf{r}^j=(\omega ^j\circ \mathbf{r}^j)\psi ^j:I_2\rightarrow \mathbb C \qquad j=1,\ldots ,J \end{aligned}$$

(3.7)

extended by zero outside $\varOmega ^j=\mathrm{supp}\,(\omega ^j\circ \mathbf{r}^j)$ constitute a solution of the system

$$\begin{aligned} {\frac{1}{2}} \phi ^i+\sum _{j=1}^J \mathrm{A}^{ij}_{\mathrm{reg},\delta }\phi ^j+\sum _{j=1}^J \mathrm{A}^{ij}_{\mathrm{sing},\delta }\phi ^j=U^i, \qquad i=1,\ldots ,J. \end{aligned}$$

(3.8)

Equation (3.8) amounts to a system of equations for the vector $(\phi ^j)\in \left( L^2(I_2)\right) ^J$ for a given right-hand side $(V^i)\in \left( L^2(I_2)\right) ^J$. With this understanding, Theorem 3 shows that the system (3.8) has a unique solution for any right-hand side $(U^i)\in \left( L^2(I_2)\right) ^J$. [It follows that for the particular right-hand side (3.8) the solution is (3.7)—which, clearly, is independent of $\delta$ in spite of the $\delta$-dependence of the system of Eq. (3.8)]. Theorem 2 shows that for the case $U^i:=-(\omega ^iU^\mathrm{inc})\circ \mathbf{r}^i$, the solution of (2.2) can be reconstructed from the solution of (3.8).

3.2 Analysis of the continuous system

Theorem 2

Let $(\phi ^1,\ldots ,\phi ^J)$ be a solution of (3.8) with $U^i:=-(\omega ^i U^{\mathrm{inc}})\circ \mathbf{r}^i$. Then the solution of (2.2) can be expressed in the form

$$\begin{aligned} \psi (\mathbf{r}):= \sum _{j \in \mathcal I (\mathbf{r})}\phi ^j((\mathbf{r}^j)^{-1}(\mathbf{r})), \end{aligned}$$

see Eq. (2.35).

Proof

Because of the particular form of the right-hand side as well as the presence of the factor $\omega ^i \circ \mathbf{r}^i$ in the operator $\mathrm{A}^{ij}_\delta$ [see Eqs. (3.3), (3.4) and  (3.8)], it follows that $\mathrm{supp}\quad \phi ^i\subseteq \mathrm{supp}(\omega ^i \circ \mathbf{r}^i) =\varOmega ^i$. Therefore, by (3.1b), $\widetilde{\omega }^i_\delta \phi ^i \equiv \phi ^i$ for all $i$.

Consider now the functions $\psi ^i: D^i \rightarrow \mathbb C$ given by

$$\begin{aligned} \psi ^i:= -2 \sum _{j=1}^J \int _{I_2} K^{ij} (\,\cdot \,, \underline{v}) \phi ^j(\underline{v}) \mathrm{d }\underline{v}-2\,U^{\mathrm{inc}}\circ \mathbf{r}^i. \end{aligned}$$

These functions are constructed so that $\phi ^i= (\omega ^i \circ \mathbf{r}^i)\, \psi ^i$ for all $i$. Although the functions $\psi ^i$ are infinitely differentiable only up to the boundary of $D^i$ (the domain of the chart $\mathbf{r}^i$), the function $\phi ^i$ can be smoothly extended by zero to $\mathbb R ^2$. Similarly, we can consider the functions $\varPhi ^i:S \rightarrow \mathbb C$ such that $\varPhi ^i(\mathbf{r})= \omega ^i(\mathbf{r}) \psi ^i ((\mathbf{r}^i)^{-1}(\mathbf{r}))=\phi ^i((\mathbf{r}^i)^{-1}(\mathbf{r}))$ if $\mathbf{r}\in S^i$ and are zero otherwise. These are $\mathcal C ^\infty$ functions on the surface $S$ and $\psi = \sum _{j=1}^J \varPhi ^j$.

Then, if $\mathbf{r}\in S$ and $i \in \mathcal I (\mathbf{r})$, we can write $\mathbf{r}=\mathbf{r}^i(\underline{u}^i)$ with $\underline{u}^i:=(\mathbf{r}^i)^{-1}(\mathbf{r})$, and

$$\begin{aligned} \left( {\frac{1}{2}}\psi +A\psi \right) (\mathbf{r})&= {\frac{1}{2}}\sum _{i \in \mathcal I (\mathbf{r})} \phi ^i(\underline{u}^i)+ \sum _{i=1}^J\omega ^i(\mathbf{r})(A\psi )(\mathbf{r})\\&= \sum _{i\in \mathcal I (\mathbf{r})}\Bigg ({\frac{1}{2}}\phi ^i(\underline{u}^i)+ \omega ^i(\mathbf{r}) \int _S K(\mathbf{r},\mathbf{r}^{\prime }) \Big (\sum _{j=1}^J \varPhi ^j(\mathbf{r}^{\prime })\Big ) \mathrm{d} \mathbf{r^{\prime }}\Bigg )\\&= \sum _{i\in \mathcal I (\mathbf{r})}\Bigg ({\frac{1}{2}}\phi ^i(\underline{u}^i) + \sum _{j=1}^J \omega ^i(\mathbf{r}) \int _{\mathrm{supp}\omega ^j} K(\mathbf{r},\mathbf{r}^{\prime }) \varPhi ^j(\mathbf{r}^{\prime }) \mathrm{d} \mathbf{r}^{\prime }\Bigg )\\&= \sum _{i\in \mathcal I (\mathbf{r})}\Bigg ({\frac{1}{2}}\phi ^i(\underline{u}^i) + \sum _{j=1}^J \omega ^i(\mathbf{r}) \int _{I_2} K(\mathbf{r},\mathbf{r}^j(\underline{v}))\,a^j(\underline{v})\,\phi ^j(\underline{v})\,\mathrm{d} \underline{v}\Bigg )\\&= \sum _{i\in \mathcal I (\mathbf{r})}\Bigg ({\frac{1}{2}}\phi ^i(\underline{u}^i) + \sum _{j=1}^J \omega ^i(\mathbf{r}^i(\underline{u}^i)) \int _{I_2} K^{ij}(\underline{u}^i,\underline{v}) \widetilde{\omega }^j_\delta (\underline{v}) \phi ^j(\underline{v}) \mathrm{d} \underline{v}\Bigg )\\&= \sum _{i\in \mathcal I (\mathbf{r})}\Bigg ({\frac{1}{2}}\phi ^i(\underline{u}^i)+ \sum _{j=1}^J \mathrm{A}^{ij}_\delta \phi ^j (\underline{u}^i)\Bigg )=\sum _{i\in \mathcal I (\mathbf{r})}U^i(\underline{u}^i) \!=\!-U^{\mathrm{inc}}(\mathbf{r}). \end{aligned}$$

Note that we have used the fact that $\widetilde{\omega }^j_\delta \phi ^j \equiv \phi ^j$. This finishes the proof. $\square$

The product space

$$\begin{aligned} \mathcal{H }^s:=\underbrace{{H}^s\times {H}^s\times \cdots \times {H}^s}_{{J\,\mathrm{times}}}, \end{aligned}$$

will be endowed with the product norm, also denoted by $\Vert \cdot \Vert _s$. We can then consider the matrices of operators

$$\begin{aligned} \mathcal{A }^\mathrm{reg}_\delta :=\left( \mathrm{A}^{ij}_{\mathrm{reg},\delta }\right) _{i,j=1}^J,\quad \mathcal{A }^\mathrm{sing}_\delta :=\left( \mathrm{A}^{ij}_{\mathrm{sing},\delta }\right) _{i,j=1}^J,\qquad \mathcal{A }_\delta := \mathcal{A }^\mathrm{reg}_\delta + \mathcal{A }_\delta ^\mathrm{sing}, \end{aligned}$$

as well as the identity operator $\mathcal{I }:\mathcal{H }^s\rightarrow \mathcal{H }^s$. With this notation, (3.8) can be written in operator form as

$$\begin{aligned} \mathcal{B }_\delta \varvec{\phi }:=\left( {\frac{1}{2}}\mathcal{I }+ \mathcal{A }_\delta \right) \varvec{\phi }=\mathbf{U} \end{aligned}$$

(3.9)

where $\varvec{\phi } :=(\phi ^1,\phi ^2,\ldots ,\phi ^J)$ and

$$\begin{aligned} \mathbf{U}:=\left( U^1,\ldots , U^J\right) \qquad U^j=-(\omega ^jU^{\mathrm{inc}})\circ \mathbf{r}^j. \end{aligned}$$

(3.10)

Remark 7

Note that, while the operator $\mathcal{A }_\delta$ depends on $\delta$, Eqs. (3.3)–(3.5) show that, for elements $\varvec{\phi }=(\phi ^1,\ldots ,\phi ^J)$ such that $\mathrm{supp}\,\phi ^i\subseteq \mathrm{supp}(\omega ^i\circ \mathbf{r}^i)$, $\mathcal A _\delta \varvec{\phi }$ is independent of $\delta$. The following theorem shows that the operator in Eq. (3.9) is invertible, and, thus, in view of Theorem 2, for right-hand sides of the form (3.10), the solution of Eq. (3.9) is independent of $\delta$ as well.

Theorem 3

For all  $0\le \delta \le \delta _0$, the operators $\mathcal B _\delta =\frac{1}{2}\mathcal I +\mathcal A _\delta :\mathcal H ^0\rightarrow \mathcal H ^0$ are invertible. Moreover

$$\begin{aligned} \Vert \mathcal{B }_\delta \Vert _{\mathcal{H }^0\rightarrow \mathcal{H }^0}+\Vert \mathcal{B }_\delta ^{-1}\Vert _{\mathcal{H }^0\rightarrow \mathcal{H }^0}\le C. \end{aligned}$$

Proof

We will first show that $\mathcal{B }_\delta :\mathcal{H }^0\rightarrow \mathcal{H }^0$ is invertible. Propositions 1 and 3 and the compact injection $H^{1}\subset H^{0}$ prove that $\mathcal{A }_\delta :\mathcal{H }^0\rightarrow \mathcal{H }^0$ is compact (and uniformly bounded in $\delta$). Therefore, $\mathcal{B }_\delta :\mathcal{H }^0\rightarrow \mathcal{H }^0$ is bounded and Fredholm of index zero. Thus, it suffices to prove the injectivity of the operator.

Let then $\varvec{\phi }\in \mathrm{Ker}\,\mathcal{B }_\delta$, that is,

$$\begin{aligned} {\frac{1}{2}}\phi ^i +\sum _{j=1}^J \mathrm{A}^{ij}_\delta \phi ^j=0,\qquad i=1,\ldots ,J. \end{aligned}$$

Arguing as in the proof of Theorem 2, it is clear that $\widetilde{\omega }_\delta ^j\phi ^j\equiv \phi ^j$ for all $j$, which means that $\mathrm{Ker}\,\mathcal B _\delta$ is independent of $\delta$. Since $\mathrm{A}^{ij}_{\delta }:H^s\rightarrow H^{s+1}$ are continuous for all $s$, then $\phi ^i\in \mathcal{C }^\infty (\overline{I}_2)$. Following the argumentation in the proof of Theorem 2, we define

$$\begin{aligned} \psi ^i:=-2\sum _{j=1}^J \int _{I_2} K^{ij}(\,\cdot \,,\underline{v}) \widetilde{\omega }_\delta ^j (\underline{v})\phi ^j(\underline{v})\,\mathrm{d}\underline{v}=-2\sum _{j=1}^J \int _{I_2} K^{ij}(\,\cdot \,,\underline{v}) \phi ^j(\underline{v})\mathrm{d} \underline{v},\quad \qquad \end{aligned}$$

(3.11)

so that $\phi ^i=(\omega ^i\circ \mathbf{r}^i)\psi ^i$. We now set

$$\begin{aligned} \psi (\mathbf{r}):=\sum _{j \in \mathcal I (\mathbf{r})}\phi ^j((\mathbf{r}^j)^{-1}(\mathbf{r}))= \sum _{j=1}^J\varPhi ^j(\mathbf{r}), \end{aligned}$$

where the functions $\varPhi ^j:S\!\rightarrow \! \mathbb C$ are defined by extending $\omega ^j \, \psi ^j \circ (\mathbf{r}^j)^{-1}\!=\!\phi ^j \circ (\mathbf{r}^j)^{-1}$ by zero to $S \setminus S^j$. Proceeding as before, we prove that $\frac{1}{2}\psi \!+\! A\psi \!=\!0$ and therefore $\psi \!=\!0$.

From (3.11) we deduce that for $\underline{u}\in D^i$,

$$\begin{aligned} \psi ^i(\underline{u})&= -2\sum _{j=1}^J \int _{I_2} K^{ij}(\underline{u},\underline{v})\,\phi ^j(\underline{v})\,{ \mathrm d}\underline{v}=-2\sum _{j=1}^J \int _{S^j} K(\mathbf{r}^i(\underline{u}),\mathbf{r}^{\prime })\varPhi ^j(\mathbf{r}^{\prime })\,\mathrm{d}\mathbf{r}^{\prime }\\&= -2 \int _{S} K(\mathbf{r}^i(\underline{u}),\mathbf{r}^{\prime })\psi (\mathbf{r}^{\prime })\,{ \mathrm d}\mathbf{r}^{\prime }=0 \end{aligned}$$

since $\psi$ itself is zero and the functions $\varPhi ^j$ are each supported in $S^j$. Therefore $\phi ^i=(\omega ^i\circ \mathbf{r}^i)\psi ^i=0$ and the injectivity of $\mathcal B _\delta$ is proven.

To prove the uniform boundedness of $\mathcal{B }_\delta ^{-1}:\mathcal{H }^0\rightarrow \mathcal{H }^0$ for $\delta \ge 0$ we proceed as in the proof of Theorem 10.9 in [20]. Recall that $\mathcal{A }_\delta : \mathcal{H }^0\rightarrow \mathcal{H }^1$ is uniformly bounded in $\delta$ by Proposition 1 with $s=0$. Since the injection $\mathcal{H }^1\subset \mathcal{H }^0$ is compact, the set $\{\mathcal{A }_\delta \}_\delta$ turns out to be collectively compact. Moreover, it is easy to verify that $\lim _{\delta \rightarrow 0} \mathcal{A }_\delta \varvec{\xi }=\mathcal A _0\varvec{\xi }$ in $\mathcal H ^0$ for all $\varvec{\xi }\in \mathcal H ^0$. Applying that pointwise convergence is uniform on compact sets and the collective compactness of the set of operators $\{ \mathcal A _\delta \}$ it follows (see [20, Corollaries 10.5 and 10.8]) that

$$\begin{aligned} \Vert \mathcal{B }_0^{-1}(\mathcal{A }_0-\mathcal{A }_{\delta })\mathcal{A }_{\delta }\Vert _{\mathcal{H }^0\rightarrow \mathcal{H }^0}\rightarrow 0. \end{aligned}$$

(3.12)

Consider now $\mathcal{C }_\delta :=2(\mathcal I -\mathcal{B }_0^{-1}\mathcal{A }_\delta )$, which is uniformly bounded. Straightforward computations show that

$$\begin{aligned} \mathcal{C }_{\delta }\mathcal{B }_\delta = \mathcal{I }+2\mathcal{B }_0^{-1}(\mathcal{B }_0\mathcal{A }_\delta -\mathcal{B }_\delta \mathcal{A }_{\delta }) =\mathcal{I }+2\mathcal{B }_0^{-1}(\mathcal{A }_0-\mathcal{A }_{\delta })\mathcal{A }_{\delta }. \end{aligned}$$

Therefore (3.12) proves that for $\delta$ small enough $\mathcal C _\delta \mathcal B _\delta$ is invertible with uniformly bounded inverse and therefore, since $\mathcal C _\delta$ is uniformly bounded, so is $\mathcal B _\delta ^{-1}$. $\square$

3.3 Discrete system

In order to obtain a continuous system of equations from the fully discrete system (2.33) we introduce an interpolation operator $\mathrm{Q}_{N_j}:\mathcal{C }^0(\overline{I}_2)\rightarrow \mathbb T _{N_j}$ given by

$$\begin{aligned} (\mathrm{Q}_{N_j}\xi )(\underline{x}^j_{\underline{\ell }})=\xi (\underline{x}^j_{\underline{\ell }})\qquad \forall \underline{\ell }\in \mathbb Z _{N_j}. \end{aligned}$$

(3.13)

The discrete counterparts of the functions $\phi ^j$ are

$$\begin{aligned} \phi ^j_{h}:=\mathrm{Q}_{N_j}((\omega ^j\circ {\mathbf{r}}^j)\psi ^j_h), \end{aligned}$$

(3.14)

where $\psi ^j_h$ is defined by (2.34) with $\varphi ^j_{\underline{m}}$ obtained as the solution of the system (2.33). Note that

$$\begin{aligned} \phi ^j_h(\underline{x}^j_{\underline{\ell }})\!=\!\mathrm{Q}_{N_j}((\omega ^j\circ {\mathbf{r}}^j)\, \psi ^j_h)(\underline{x}_{\underline{\ell }}^j)\!=\! \mathrm{R}_{N_j}({\varvec{\omega }}^j\odot \varvec{\varphi }^j) (\underline{x}_{\underline{\ell }}^j)\!=\! \left\{ \begin{array}{ll} \omega ^j_{\underline{\ell }}\, \varphi ^j_{\underline{\ell }},\qquad &{}\quad \text{ if} \quad \underline{\ell }\in \varOmega ^j_h,\\ 0,&{}\quad \mathrm{otherwise}. \end{array} \right. \nonumber \\ \end{aligned}$$

(3.15)

On the other hand, $\big (\widetilde{\omega }^j_{\delta } \phi ^j_{h}\big )(\underline{x}_{\underline{\ell }}^j)= \phi ^j_{h}(\underline{x}_{ \underline{\ell }}^j)$ for all $\underline{\ell }$, since $\phi ^j_{h}(\underline{x}_{\underline{\ell }}^j)=0,$ for all $\underline{\ell }\not \in \varOmega ^j_h$ and inside $\varOmega ^j$ the functions $\widetilde{\omega }^j$ do not have any influence by (3.1b). Therefore (2.33) can be equivalently re-expressed as an equation for continuous biperiodic functions $(\phi ^1_h,\ldots , \phi ^J_h )$ such that for all $i\in \{1,\ldots , J\}$

$$\begin{aligned}&{\frac{1}{2}} \phi ^i_h + \sum _{j=1}^J \mathrm{Q}_{N_i} \bigg ( h_j^2 \sum _{\underline{m} \in \mathbb Z _{N_j}} \left( \omega ^i\circ \mathbf{r}^i\right) K^{ij}_{\mathrm{reg},\delta } (\,\cdot \,, \underline{x}_{\underline{m}}^j) \widetilde{\omega }^j_\delta (\underline{x}^j_{\underline{m}}) \phi _{h}^j (\underline{x}_{\underline{m}}^j) \bigg )\nonumber \\&\quad + \sum _{j=1}^J \mathrm{Q}_{N_i}\Big ( (\omega ^i\circ \mathbf{r}^i) \mathrm{L}^{ij}_{\delta ,h} \phi ^j_{h} \Big ) = \mathrm{Q}_{N_i} U^i. \end{aligned}$$

(3.16)

Note that if the system (3.16) is uniquely solvable, the solution belongs necessarily to $\mathbb T _{N_1}\times \cdots \times \mathbb T _{N_J}$. The solutions of (2.33) and (3.16) are related by the formula (3.15).

For the sake of our analysis we recast the system (3.16) in an equivalent operator form. To do this, we introduce the discrete operators

$$\begin{aligned} \mathrm{A}^{ij}_{\mathrm{sing},\delta ,h}:= (\omega ^i\circ \mathbf{r}^i) \mathrm L^{ij}_{\delta ,h}, \end{aligned}$$

(3.17)

and, for $\xi \in \mathcal C ^0(\overline{I}_2)$, we define (cf. [11] and [15])

$$\begin{aligned} \mathrm{D}_{N_j} \xi := h_j^2 \sum _{\underline{m}\in \mathbb Z _{N_j}} \xi (\underline{x}^j_{\underline{m}}) \delta _{\underline{x}_{\underline{m}}^j}, \end{aligned}$$

(3.18)

where $\delta _{\underline{u}}$ denotes here the Dirac delta distribution at the point $\underline{u} \in \overline{I}_2$: $\delta _{\underline{u}}\psi :=\psi (\underline{u})$. Further, we note that since the operators $\mathrm{A}^{ij}_{\mathrm{reg},\delta }$ have continuous kernels, they may be applied to delta distributions:

$$\begin{aligned} \mathrm{A}^{ij}_{\mathrm{reg},\delta } \delta _{\underline{u}}= \left( \omega ^i\circ \mathbf{r}^i\right) K^{ij}_{\mathrm{reg},\delta } (\,\,\cdot \,\,, {\underline{u}}) \widetilde{\omega }^j_\delta ({\underline{u}}). \end{aligned}$$

In view of these definitions, we clearly have

$$\begin{aligned} \mathrm{A}^{ij}_{\mathrm{reg},\delta } \mathrm{D}_{N_j}\xi = h_j^2 \sum _{\underline{m} \in \mathbb Z _{N_j}} \left( \omega ^i\circ \mathbf{r}^i\right) K^{ij}_{\mathrm{reg},\delta } (\,\,\cdot \,\,, \underline{x}_{\underline{m}}^j) \widetilde{\omega }^j_\delta (\underline{x}^j_{\underline{m}})\xi (\underline{x}_{\underline{m}}^j), \end{aligned}$$

and (3.16) is equivalent to the equation

$$\begin{aligned} {\frac{1}{2}} \phi ^i_{h} \!+\! \sum _{j=1}^J \mathrm{Q}_{N_i} \Big (\mathrm{A}^{ij}_{\mathrm{reg},\delta } \mathrm{D}_{N_j} \phi ^j_{h} \!+\! \mathrm{A}^{ij}_{\mathrm{sing},\delta ,h} \phi ^j_{h}\Big ) \!=\! \mathrm{Q}_{N_i} U^i, \qquad i\!=\!1,\ldots , J,\quad \quad \end{aligned}$$

(3.19)

for the unknowns $(\phi ^1_{h},\ldots ,\phi ^J_{h})\in \mathbb T _{N_1}\times \cdots \times \mathbb T _{N_J}$.

In order to recast this system of operator equations as a system defined in $L^2$ spaces, we insert the orthogonal projections $\mathrm{P}_{N_j}:L^2(I_2)\rightarrow \mathbb T _{N_j}$ and write (3.16) in the form

$$\begin{aligned} {\frac{1}{2}} \phi ^i_h + \sum _{j=1}^J \mathrm{Q}_{N_i} \Big (\mathrm{A}^{ij}_{\mathrm{reg},\delta } \mathrm{D}_{N_j} {\mathrm{P}}_{N_j} \phi ^j_h+ \mathrm{A}^{ij}_{\mathrm{sing},\delta ,h} \mathrm{P}_{N_j}\phi ^j_h\Big ) = \mathrm{Q}_{N_i} U^i, \qquad i=1,\ldots , J.\nonumber \\ \end{aligned}$$

(3.20)

for the unknowns $(\phi ^1_h,\ldots ,\phi ^J_h )\in L^2(I_2)\times \cdots \times L^2(I_2)$.

Remark 8

Note that the projection operator $\mathrm{P}_{N_j}$, which maps $L^2$ to the space of trigonometric polynomials, makes it possible to recast the fully discrete equation (3.16) [whose unknowns are approximate point values of the continuous solution of Eq. (2.15)] via its version (3.19), posed in the space of trigonometric polynomials, as an Eq. (3.20) in the complete space $L^2$.

3.4 Mapping properties of the operators introduced in Sect. 3.1

As a result of the work in Sects. 3.1 and 3.3, our overall integral equation and its Nyström discretization have been re-expressed as Eqs. (3.8) and (3.20), respectively, which are posed in terms of unknowns defined (and compactly supported) on the open unit square $I_2$. Such functions can naturally be extended to period-1 biperiodic functions defined in $\mathbb R ^2$, on which our analysis is based. In this section, in particular, we study the mapping properties of the operators introduced in Sect. 3.1, when viewed as operators defined on spaces of periodic functions.

To do this we first consider the biperiodic Sobolev spaces, see e.g. [23], and we study their connections (that result through our use of local charts and partitions of unity) with the classical Sobolev spaces on the surface $S$. The Fourier coefficients of any locally integrable biperiodic function $\xi :\mathbb R ^2\rightarrow \mathbb C$ are given by

$$\begin{aligned} \widehat{\xi }(\underline{m}):=\int _{I_2} \xi (\underline{u}) e_{\underline{m}}(\underline{u})\mathrm{d}\underline{u}, \qquad e_{\underline{m}}(\underline{u}):=\exp (2\pi \mathrm{i}\underline{m} \cdot \underline{u}), \quad \underline{m}\in \mathbb Z ^2. \end{aligned}$$

For arbitrary $s\in \mathbb R$, the Sobolev norm

$$\begin{aligned} \Vert \xi \Vert _{s}:=\Bigg (|\widehat{\xi }(\underline{0})|^2+\sum _{\underline{m} \ne \underline{0} }(|m_1|^2+|m_2|^2)^s|\widehat{\xi }(\underline{m})|^2\Bigg )^{1/2} \end{aligned}$$

(3.21)

is well defined for all $\xi$ in the space $\mathbb T :=\mathrm{span}\{ e_{\underline{m}}\,:\, \underline{m} \in \mathbb Z ^2\}=\cup _N \mathbb T _N$ of trigonometric polynomials. The Sobolev space $H^s$ is defined as the completion of $\mathbb T$ under the norm $\Vert \,\cdot \,\Vert _s$. Note that $H^0$ is the space of biperiodic extensions of functions in $L^2(I_2)$. By a simple density argument, we can define the Fourier coefficients for any element of $H^s$ and any $s<0$. It is also possible to define partial derivatives of any order: for a given non-negative bi-index $\underline{\alpha }=(\alpha _1,\alpha _2)$, using the notation $|\underline{\alpha }|=\alpha _1+\alpha _2$, we see that $\partial _{\underline{\alpha }}$ is a bounded operator from $H^s$ to $H^{s-|\underline{\alpha }|}$ for all $s$.

The atlas introduced in Sect. 2 for representation of the surface $S$ together with the $H^s$ norm just defined gives rise to a definition of the Sobolev norm

$$\begin{aligned} \Vert U\Vert _{H^s(S)}:=\bigg (\sum _{i=1}^J\Vert (\omega ^i U)\circ \mathbf{r}^i\Vert ^2_{s}\bigg )^{1/2} \end{aligned}$$

(3.22)

for any $U\in \mathcal C ^\infty (S)$ and any $s \in \mathbb R$. (In this formula it has been implicitly assumed that, even though $\mathbf{r}^i$ is only defined on $D^i\subset I_2$, the function $(\omega ^i U)\circ \mathbf{r}^i$ can be extended by zero to a function in $\mathcal C ^\infty (I_2)$—since the support of $\omega ^i\circ \mathbf{r}^i$ is contained in $D^i$—and the result can be extended as a $\mathcal C ^\infty$ biperiodic function to all of $\mathbb R ^2$.) The space $H^s(S)$ can then be then defined, for instance, as the completion of $\mathcal{C }^\infty (S)$ in the above norm: this definition is equivalent to the classical definitions given in e.g. [1], [18, §1.3.3], [24, §2.4] and [22, Chapter 3].

Lemma 1

For all $s\in \mathbb R$ and $j\in \{1,\ldots ,J\}$

$$\begin{aligned} \Vert (\widetilde{\omega }_\delta ^j\xi )\circ (\mathbf{r}^j)^{-1}\Vert _{H^s(S)} \le C_s\delta ^{-|s|}\Vert \xi \Vert _{s} \qquad \forall \xi \in H^s, \qquad \forall \delta , \end{aligned}$$

(3.23)

where $C_s$ is a constant independent of $\delta$.

Proof

Because of (3.22) we only need to prove that the maps

$$\begin{aligned} \xi \longmapsto \mathrm{T}^{ji}_\delta \xi := \underbrace{(\omega ^i\circ \mathbf{r}^i)(\widetilde{\omega }^j_\delta \circ \underline{r}^{ji})}_{=:\varrho _\delta ^{ji}}\, (\xi \circ \underline{r}^{ji}) \end{aligned}$$

(3.24)

map $H^s$ to $H^s$ for every $s$ and that their norms are bounded by a multiple of $\delta ^{-|s|}$. To do this we first note that

$$\begin{aligned} \mathrm{supp}\,\varrho ^{ji}_\delta \subset (\mathbf{r}^i)^{-1} (S^j\cap S^i), \end{aligned}$$

(3.25)

(which is the domain where $\underline{r}^{ji}$ is defined), and that, in view of Eq. (3.1c), we have

$$\begin{aligned} \Vert \partial _{\underline{\alpha }}\varrho ^{ji}_\delta \Vert _{\infty ,I_2}\le C_{\underline{\alpha }}\delta ^{-|\alpha |}. \end{aligned}$$

(3.26)

Since the operator $\mathrm{T}^{ji}_\delta$ is given by multiplication by the $\mathcal{C }^\infty$ function $\varrho ^{ji}_\delta$ preceded by application of the smooth ($\delta$-independent) diffeomorphism $\underline{r}^{ji}:(\mathbf{r}^i)^{-1} (S^j\cap S^i) \rightarrow (\mathbf{r}^j)^{-1}(S^j\cap S^i)$, using the differential form of the Sobolev norms $H^s$ for positive integers $s$ together with Eqs. (3.25) and (3.26), we see that

$$\begin{aligned} \Vert \mathrm{T}^{ji}_\delta \xi \Vert _s \le C_s \delta ^{-s}\Vert \xi \Vert _s \qquad \forall \xi \in H^s \qquad s=0,1,2,\ldots \end{aligned}$$

Together with the interpolation properties [23] of the Sobolev spaces $H^s$, this bound establishes the result for all $s\ge 0$. The result for $s<0$ follows from a duality argument. $\square$

Here and in the sequel, $\Vert A\Vert _{X\rightarrow Y}$ denotes the operator norm of the bounded operator $A:X\rightarrow Y$ between the Hilbert spaces $X$ and $Y$.

Proposition 1

For all $s\in \mathbb R$ and all indices $i,j$, the integral operators $\mathrm{A}^{ij}_\delta :H^s\rightarrow H^{s+1}$ given by Eq. (3.5) are continuous and we have

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\delta }\Vert _{H^s\rightarrow H^{s+1}}\le C_s\delta ^{-|s|}. \end{aligned}$$

Proof

We can write $\mathrm{A}^{ij}_\delta$ as the composition (from right to left) of the maps $\xi \mapsto (\widetilde{\omega }^j_\delta \xi )\circ (\mathbf{r}^j)^{-1}$ with $A$ [the combined integral operator given in (2.3)] and then $\varphi \mapsto (\omega ^i \varphi )\circ \mathbf{r}^i$. Lemma 1 provides a bound for the norm of first of these maps, whereas the norm of the last function as a map from $H^s(S)$ to the periodic Sobolev space $H^s$ is clearly bounded in view of the definition (3.22) of the surface Sobolev norms. The result therefore follows from the fact that $A$ is a bounded operator from $H^s(S)$ to $H^{s+1}(S)$ for all $s \in \mathbb R$—as it follows from standard results concerning pseudodifferential operators on smooth surfaces (see e.g. [19, Chapters 6–9]). $\square$

The next result studies the mapping properties of the regular and singular part of $\mathrm{A}^{ij}_\delta$ in the frame of the periodic Sobolev spaces and gives some estimates for the continuity constants in terms of $\delta$. In some of the arguments, it is convenient to use the space $\mathcal D :=\cap _{s\in \mathbb R }H^s$, consisting of $\mathcal C ^\infty$ biperiodic functions.

Proposition 2

For all $i,j=1,\ldots ,J$, and $s,t\in \mathbb R$ the operators

$$\begin{aligned} \mathrm{A}_{\mathrm{reg},\delta }^{ij}:{H}^s\rightarrow {H}^t,\qquad \mathrm{A}_{\mathrm{sing},\delta }^{ij}: {H}^s\rightarrow {H}^{s+1} \end{aligned}$$

are continuous. Moreover, for all $t\ge 0\ge s$ and all $r\in \mathbb R$ there exist positive constants $C_{s,t}$ and $C_r$ such that

$$\begin{aligned} \Vert \mathrm{A}_{\mathrm{reg},\delta }^{ij}\Vert _{{H}^{s}\rightarrow {H}^t}&\le C_{s,t}\delta ^{s-t},\end{aligned}$$

(3.27)

$$\begin{aligned} \Vert \mathrm{A}_\mathrm{sing,\delta }^{ij}\Vert _{{H}^{r}\rightarrow {H}^{r+1}}&\le C_r\delta ^{-\max \{-r,r+1,1\}}. \end{aligned}$$

(3.28)

Finally

$$\begin{aligned} \Vert \mathrm{A}_\mathrm{sing,\delta }^{ij}\Vert _{{H}^{0}\rightarrow {H}^{0}} \le C_0. \end{aligned}$$

(3.29)

Proof

Let

$$\begin{aligned} F_{\mathrm{reg},\delta }^{ij}(\underline{u},\underline{v}):=(\omega ^i\circ \mathbf{r}^i)(\underline{u}) K^{ij}_{\mathrm{reg},\delta }(\underline{u},\underline{v}) \widetilde{\omega }^j_\delta (\underline{v}) \end{aligned}$$

be the integral kernel of the operator $\mathrm{A}_{\mathrm{reg},\delta }^{ij}$. This function is well defined on $D^i\times D^j$ and can be extended by zero to the rest of $I_2\times I_2$ thanks to the cut-off functions that appear in its definition. In particular, $F_{\mathrm{reg},\delta }^{ij}$ admits a $\mathcal C ^\infty$ biperiodic extension to $\mathbb R ^2\times \mathbb R ^2$, and, thus $\mathrm{A}_{\mathrm{reg},\delta }^{ij}$ is a continuous operator from $H^s$ to ${H}^t$ for all $s,t\in \mathbb R$ (see [23, Theorem 6.1.1]). Since $\mathrm{A}_{\mathrm{sing},\delta }^{ij}=\mathrm{A}^{ij}_\delta -\mathrm{A}_{\mathrm{reg},\delta }^{ij}$, using Proposition 1, the claimed mapping properties of $\mathrm{A}_{\mathrm{sing},\delta }^{ij}$ follow directly.

To establish (3.27) we consider the operators

$$\begin{aligned} D_{\underline{u}}^{m}\xi := \partial _{u_1}^{2m}\xi + \partial _{u_2}^{2m}\xi ,\qquad \varLambda _m \xi :=D_{\underline{u}}^{m}\xi +\widehat{\xi }(\underline{0}). \end{aligned}$$

(3.30)

(Note that $\varLambda _m:H^s\rightarrow H^{s-2m}$ is a bounded isomorphism for all $s$.) Letting

$$\begin{aligned} G(\underline{u},\underline{v}):=D_{\underline{u}}^{m} F_{\mathrm{reg},\delta }^{ij}(\underline{u},\underline{v})+\int _{I_2} F_{\mathrm{reg},\delta }^{ij}(\underline{w},\underline{v})\,\mathrm{d}{\underline{w}}, \end{aligned}$$

we note that $G\in \mathcal{C }^\infty (\mathbb R ^2\times \mathbb R ^2)$, that $G$ is periodic in each of its variables, and that, when restricted to $I_2\times I_2$, $G$ has compact support. Now, for all $\eta \in \mathcal D$ we have

$$\begin{aligned} \varLambda _m \mathrm{A}_{\mathrm{reg},\delta }^{ij}\varLambda _n \eta&= \int _{I_2} G(\cdot ,\underline{v})D_{\underline{v}}^{n}\eta (\underline{v}) \,\mathrm{d}\underline{v}+\widehat{\eta }(\underline{0}) \int _{I_2} G(\cdot ,\underline{v})\,\mathrm{d}\underline{v} \\&= \int _{I_2} D_{\underline{v}}^{n} G(\cdot ,\underline{v})\eta (\underline{v}) \,\mathrm{d}\underline{v}+\widehat{\eta }(\underline{0}) \int _{I_2} G(\cdot ,\underline{v})\,\mathrm{d}\underline{v}. \end{aligned}$$

Further, bounds of the form

$$\begin{aligned} \sup _{\underline{u},\underline{v}\in I_2} |D_{\underline{u}}^{m}D_{\underline{v}}^{n} F_{\mathrm{reg},\delta }^{ij}(\underline{u},\underline{v})|\le C_{m,n} \delta ^{-2m-2n}, \end{aligned}$$

(3.31)

follow from the definition the kernel functions $K^{ij}_{\mathrm{reg},\delta }$ [see (2.13) and subsequent lines] and from the assumptions (3.1c) on $\widetilde{\omega }^{j}_\delta$. It follows that, for $\xi \in \mathcal{D }$ we have

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\xi \Vert _{2m}&\le C_{2m}\Vert \varLambda _m \mathrm{A}^{ij}_{\mathrm{reg},\delta }\xi \Vert _{0}=C_{2m}\Vert \varLambda _m \mathrm{A}^{ij}_{\mathrm{reg},\delta }\varLambda _n \varLambda _n^{-1}\xi \Vert _{0}\\&\le C_{2m,2n}\delta ^{-2m-2n}\Vert \varLambda _n^{-1}\xi \Vert _{0}\le C_{2m,2n}\delta ^{-2m-2n}\Vert \xi \Vert _{-2n}, \end{aligned}$$

where we have used the integro-differential form of $\varLambda _m \mathrm{A}^{ij}_{\mathrm{reg},\delta }\varLambda _n$ and (3.31). This inequality establishes (3.27) in the particular case $t=2m, s=-2n$; the result for general $t\ge 0\ge s$ then follows by interpolation [23].

To establish (3.28), in turn, we first note that (3.27) implies

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\Vert _{H^r\rightarrow H^{r+1}}\le \left\{ \begin{array}{ll} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\Vert _{H^r\rightarrow H^{0}}\le C_{r,0} \delta ^{r},&{}\quad \text{ if}\,r<-1,\\ \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\Vert _{H^r\rightarrow H^{r+1}}\le C_{r,r+1} \delta ^{-1},&{}\quad \text{ if}\,-1\le r<0,\\ \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\Vert _{H^{0}\rightarrow H^{r+1}}\le C_{0,r+1} \delta ^{-(r+1)},&{}\quad \text{ if}\,r 0, \end{array} \right. \end{aligned}$$

or, equivalently

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\Vert _{H^r\rightarrow H^{r+1}} \le C^{\prime }_{r} \delta ^{-\max \{-r,r+1,1\}}. \end{aligned}$$

Since $\mathrm{A}^{ij}_{\mathrm{sing},\delta }=\mathrm{A}^{ij}_\delta -\mathrm{A}^{ij}_{\mathrm{reg},\delta }$, and in view of Proposition 1, we have

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }\Vert _{H^r\rightarrow H^{r+1}}&\le \Vert \mathrm{A}^{ij}_\delta \Vert _{H^r\rightarrow H^{r+1}}+\Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\Vert _{H^r\rightarrow H^{r+1}}\\&\le C_r \delta ^{-|r|} +C_r^{\prime }\delta ^{-\max \{-r,r+1,1\}} \le C_r^{\prime \prime }\delta ^{-\max \{-r,r+1,1\}}. \end{aligned}$$

which establishes (3.28).

Inequality (3.29), finally, follows from (3.27) with $s=t=0$ and Proposition 1:

$$\begin{aligned} \Vert \mathrm{A}^{ij}_\delta \Vert _{H^{0}\rightarrow H^{0}}\le \Vert \mathrm{A}^{ij}_\delta \Vert _{H^{0}\rightarrow H^{1}} \le C, \end{aligned}$$

with $C$ independent of $\delta .$ $\square$

The limiting case $\delta =0$ is studied next.

Proposition 3

For all $i,j$, the operators $\mathrm{A}^{ij}_0:H^0\rightarrow H^1$ are bounded.

Proof

Clearly, $\mathrm{A}^{ij}_{0}\xi =\mathrm{A}^{ij}_{\delta _0}\,(\widetilde{\omega }_0^j\xi )$, where $\widetilde{\omega }_0^j$ coincides with the characteristic function of set $\varOmega ^j$ (see (3.2)). Since the operator of multiplication by $\widetilde{\omega }_0^j$ defines a continuous operator in $H^0$, the result follows from Proposition 1. $\square$

3.5 Convergence estimates in the biperiodic framework

We now state a convergence theorem (Theorem 5) in the biperiodic framework; as shown in Sect. 5, this theorem is equivalent to our main convergence result, Theorem 1. The proofs of the two results presented in this section (Theorems 4 and 5), which require certain analytical tools that are developed in the next section, are given in sect. 5.3.

Let us introduce the operators

$$\begin{aligned} \mathcal Q _h \varvec{\xi }&:= \left( \mathrm{Q}_{N_1}\xi _1,\ldots ,\mathrm{Q}_{N_J}\xi _J\right) ,\\ \mathcal P _h \varvec{\xi }&:= \left( \mathrm{P}_{N_1}\xi _1,\ldots ,\mathrm{P}_{N_J}\xi _J\right) ,\\ \mathcal D _h \varvec{\xi }&:= \left( \mathrm{D}_{N_1}\xi _1,\ldots ,\mathrm{D}_{N_J}\xi _J\right) , \end{aligned}$$

see Sect. 3.3, and let

$$\begin{aligned} \mathcal{A }^\mathrm{sing}_{\delta , h}:=\big (\mathrm{A}_{\mathrm{sing},\delta ,h}^{ij}\big )_{i,j=1}^J. \end{aligned}$$

Clearly, Eq. (3.20) [and, thus, in view of Remark 8, the main Nyström system presented in Eq. (2.33)], can be re-expressed in the form

$$\begin{aligned} \mathcal{B }_{\delta ,h}\varvec{\phi }_h:=\left( {\frac{1}{2}}\mathcal{I }+\mathcal Q _h \mathcal{A }_\delta ^\mathrm{reg}\mathcal D _h\mathcal P _h+ \mathcal Q _h \mathcal{A }^\mathrm{sing}_{\delta ,h} \mathcal P _h\right) \varvec{\phi }_h=\mathcal Q _h \mathbf{U} \end{aligned}$$

(3.32)

where $\varvec{\phi }_h:=(\phi _{h}^1,\phi _{h}^2,\ldots ,\phi _{h}^J)$. (Note that, although not explicit in the notation, the unique solution $\varvec{\phi }_h$ of Eq. (3.32) does depend on the parameter $\delta$.)

Theorem 4

Let $\delta \approx h^{\beta }$ with $\beta \in (0,1)$. Then

$$\begin{aligned} \lim _{h \rightarrow 0}\Vert \mathcal{B }_{\delta ,h}-\mathcal{B }_\delta \Vert _{\mathcal{H }^0\rightarrow \mathcal{H }^0}= 0. \end{aligned}$$

Thus, for all $h$ small enough, $\mathcal{B }_{\delta ,h}:\mathcal{H }^0\rightarrow \mathcal{H }^0$ is invertible, with $h$-uniformly bounded inverse.

Hence, for $h$ small enough the numerical scheme admits a unique solution which depends continuously on the right-hand side.

Theorem 5

Let $\delta \approx h^{\beta }$ with $\beta \in (0,1)$ and let $\mathbf U$ be given by (3.10). Let $\varvec{\phi }$ and $\varvec{\phi }_h$ be the respective solutions of $\mathcal{B }_\delta \varvec{\phi }=\mathbf{U}$ and $\mathcal{B }_{\delta ,h}\varvec{\phi }_h=\mathcal Q _h\mathbf{U}$. Then, there exist constants $C_{s,t}>0$ for all $t\ge s\ge 0$ and $t>1$ such that

$$\begin{aligned} \Vert \varvec{\phi }-\varvec{\phi }_h\Vert _{s}\le C_{s,t} {h^{t(1-\beta )-s}}\Vert \varvec{\phi }\Vert _{t}. \end{aligned}$$

Note that $\varvec{\phi }$ does not depend on $\delta$ (see Remark 7), although $\varvec{\phi }_{h}$ is a $\delta$ dependent quantity.

4 Estimates for the singular kernel and associated singular operator

In this section we discuss the regularity properties of the singular kernel, and we present a non-standard estimate involving $L^\infty$ and $L^2$ norms over spaces of trigonometric polynomials (Proposition 5) for the continuity constants of the associated singular operator.

Lemma 2

The functions

$$\begin{aligned} (\underline{u},\rho ,\theta ) \longmapsto |\rho |K_1(\mathbf{r}^i(\underline{u}),\mathbf{r}^j(\underline{r}^{ji}(\underline{u})+\rho \,\underline{e}(\theta ))\,) \end{aligned}$$

are $\mathcal C ^\infty$ in their domain of definition.

Proof

The kernel function $K_1$ can be decomposed (see (2.5)) as

$$\begin{aligned} K_1(\mathbf{r},\mathbf{r}^{\prime })=\frac{F_1(\mathbf{r},\mathbf{r}^{\prime })}{|\mathbf{r}-\mathbf{r}^{\prime }|}+ \frac{F_2(\mathbf{r},\mathbf{r}^{\prime })}{|\mathbf{r}-\mathbf{r}^{\prime }|^3}\,(\mathbf{r}-\mathbf{r}^{\prime })\cdot \varvec{\nu }(\mathbf{r}^{\prime }), \end{aligned}$$

where $F_1, F_2\in \mathcal C ^\infty (S\times S)$. Let

$$\begin{aligned} f_\ell (\underline{u},\rho ,\theta )&:= F_\ell (\mathbf{r}^i(\underline{u}),\mathbf{r}^j(\underline{r}^{ji}(\underline{u})+\rho \,\underline{e}(\theta )))\qquad \ell =1,2,\\ g_1(\underline{u},\rho ,\theta )&:= \rho ^{-2} |\mathbf{r}^i(\underline{u})- \mathbf{r}^j (\underline{r}^{ji}(\underline{u})+\rho \underline{e}(\theta ))|^2\\ g_2(\underline{u},\rho ,\theta )&:= \rho ^{-2}\Big (\mathbf{r}^i(\underline{u})- \mathbf{r}^j (\underline{r}^{ji}(\underline{u})+\rho \underline{e}(\theta ))\Big ) \cdot \varvec{\nu }(\underline{r}^{ji}(\underline{u})+\rho \,\underline{e}(\theta )). \end{aligned}$$

It is clear that $f_1$ and $f_2$ are $\mathcal C ^\infty$. Noticing that the functions

$$\begin{aligned} \underline{z}&\longmapsto&|\mathbf{r}^i(\underline{u})-\mathbf{r}^j(\underline{r}^{ji}(\underline{u})+\underline{z})|^2\end{aligned}$$

(4.1)

$$\begin{aligned} \underline{z}&\longmapsto&\Big (\mathbf{r}^i(\underline{u})-\mathbf{r}^j(\underline{r}^{ji}(\underline{u})+\underline{z})\Big )\cdot \varvec{\nu }(\underline{r}^{ji}(\underline{u})+\underline{z}) \end{aligned}$$

(4.2)

satisfy the hypotheses of Lemma 3, it follows that $g_1$ and $g_2$ are $\mathcal C ^\infty$. It is also clear that $g_1$ is positive for $\rho \ne 0$. The Hessian matrix at $\underline{z}=\underline{0}$ of the function defined in (4.1) is the matrix with elements

$$\begin{aligned} 2 \partial _{z_{n_1}} \mathbf{r}^j(\underline{u})\cdot \partial _{z_{n_2}} \mathbf{r}^j(\underline{u}) \qquad 1\le n_1, n_2\le 2. \end{aligned}$$

Clearly this matrix is positive semidefinite, and since its determinant is

$$\begin{aligned} 4 |\partial _{z_1}\mathbf{r}^j(\underline{u})\times \partial _{z_2}\mathbf{r}^j(\underline{u})|^2 = 4a^j(\underline{u})^2, \end{aligned}$$

(see (2.7)), it is positive definite. Using the equality (4.3), it follows that $g_1(\underline{u},0,\theta )\ne 0$. The mapping mentioned in the statement of the present lemma can be formulated in terms of the expression

$$\begin{aligned} g_1^{-1/2} f_1+ g_1^{-3/2} f_2\,g_2; \end{aligned}$$

the previous arguments show this mapping is infinitely differentiable, and the lemma thus follows. $\square$

Lemma 3

Let $F$ be a $\mathcal C ^\infty$ function in a neighborhood of the origin in $\mathbb R ^2$. If $F(\underline{0})=0$ and $\nabla F(\underline{0})=\underline{0}$, then the function

$$\begin{aligned} f(\rho ,\theta ):=\rho ^{-2}F(\rho \,\underline{e}(\theta )) \end{aligned}$$

is $\mathcal C ^\infty$. Moreover,

$$\begin{aligned} f(0,\theta )={\frac{1}{2}}\, \underline{e}(\theta )\cdot \mathrm HF(\underline{0})\underline{e}(\theta ), \end{aligned}$$

(4.3)

where $\mathrm{H}F(\underline{0})$ is the Hessian matrix of $F$ at the origin.

Proof

The result follows from an application of the Taylor formula

$$\begin{aligned} h(\rho ,\theta )= h(0,\theta ) + \rho \partial _\rho h(0,\theta ) + \rho ^2 \int _0^1 \partial _\rho ^2 h(\rho \,u,\theta )\,(1-u)\mathrm{d} u. \end{aligned}$$

to the function $h(\rho ,\theta ):= F(\rho \,\underline{e}(\theta ))$. $\square$

The forthcoming analysis relies heavily on use of families $\{ \chi _\delta \,:\,0<\delta \le \delta _0\}$ of functions for which there exist constants $c_0>0$ and $C_{\underline{\alpha }}$ for $\underline{\alpha }\ge \underline{0}$ such that for all $\delta$ ($0<\delta \le \delta _0$)

$$\begin{aligned}&\chi _\delta \in \mathcal C ^\infty (\mathbb R ^2), \end{aligned}$$

(4.4a)

$$\begin{aligned}&\chi _\delta (\rho ,\theta )=\chi _\delta (\rho ,\theta +2\pi ) \qquad \forall (\rho ,\theta ) \in \mathbb R ^2,\end{aligned}$$

(4.4b)

$$\begin{aligned}&\mathrm{supp}\,\chi _\delta \subset (-c_0\,\delta ,c_0\,\delta )\times \mathbb R ,\end{aligned}$$

(4.4c)

$$\begin{aligned}&\Vert \partial _{\underline{\alpha }} \chi _\delta \Vert _{L^\infty (\mathbb R ^2)}\le C_{\underline{\alpha }} \delta ^{-|\underline{\alpha }|}\ \quad \forall \underline{\alpha }\ge \underline{0}. \end{aligned}$$

(4.4d)

Proposition 4

For all $\underline{u}\in {\varOmega ^{ij}_{2\delta _0}}={(\mathbf{\mathbf{r}}^i)^{-1}(S^i\cap S_{2\delta _0}^j)}$ and $\delta \le \delta _0$, let us define:

$$\begin{aligned} \chi _{\underline{u},\delta }(\rho ,\theta )\!:=\!{\frac{1}{2}} \omega ^i(\mathbf{r}^i(\underline{u}))|\rho | \, K^{ij}_{\mathrm{sing},\delta }(\underline{u},\underline{r}^{ji}(\underline{u})\!+\!\rho \underline{e}(\theta ))\, \widetilde{\omega }_\delta ^j (\underline{r}^{ji}(\underline{u})\!+\!\rho \underline{e}(\theta ))\qquad \end{aligned}$$

(4.5)

Then the sequence $\{\chi _{\underline{u},\delta }\}$ satisfies conditions (4.4) with $c_0$ and $C_{\underline{\alpha }}$ independent of $\underline{u}$.

Proof

We start by considering the simpler functions

$$\begin{aligned} \psi _{\underline{u},\delta }(\rho ,\theta ):=|\rho |\, K^{ij}_\mathrm{sing,\delta }(\underline{u},\underline{r}^{ij}(\underline{u})+\rho \underline{e}(\theta )). \end{aligned}$$

(4.6)

Since $\varOmega ^{ij}_{2\delta _0}\subset (\mathbf{r}^{i})^{-1}(S^i\cap S^j)$, the mapping $\underline{r}^{ji}:=(\mathbf{r}^j)^{-1}\circ \mathbf{r}^i:\varOmega ^{ij}_{2\delta _0} \subset I_2\rightarrow I_2$ is well defined. Therefore, the functions

$$\begin{aligned} \eta _{\underline{u},\delta }(\underline{z})&:= \eta _\delta ({\mathbf{r}}^i(\underline{u}), {\mathbf{r}}^j(\underline{r}^{ji}(\underline{u})+ \underline{z})),\\ f_{\underline{u}}(\underline{z})&:= |\underline{z}|\, K_{1}(\mathbf{r}^i(\underline{u}),\mathbf{r}^j(\underline{r}^{ji} (\underline{u})+\underline{z})) \end{aligned}$$

are well defined in $\varOmega _{\underline{u}}:=\{\underline{z}\,:\, \underline{r}^{ji}(\underline{u})+\underline{z}\in D^j\}$ (that is $\varOmega _{\underline{u}} = D^j-\underline{r}^{ji}(\underline{u})$).

As a simple consequence of (2.9), if $\mathbf{r}\in S^j_{2\delta _0}$, then $\overline{B(\mathbf{r},\epsilon _1\delta )} \subset S^j$. Applying (2.11) and the fact that $\mathbf{r}^i(\underline{u})\in S^{j}_{2\delta _0}$, it follows that

$$\begin{aligned} \mathrm{supp}\,\eta _\delta (\mathbf{r}^i (\underline{u}),\,\cdot \,) \subset \overline{B(\mathbf{r}^i(\underline{u}),\epsilon _1\delta )} \subset S^j. \end{aligned}$$

(4.7)

Consequently

$$\begin{aligned} \mathrm{supp}\,\eta _{\underline{u},\delta } \subset (\mathbf{r}^j)^{-1}( \overline{B(\mathbf{r}^i(\underline{u}),\epsilon _1\delta ) )} -\underline{r}^{ji}(\underline{u}) \subset (\mathbf{r}^j)^{-1}(S^j)-\underline{r}^{ji}(\underline{u})=\varOmega _{\underline{u}}. \end{aligned}$$

Therefore $\mathrm{supp}\quad \eta _{\underline{u},\delta }$ is strictly contained in $\varOmega _{\underline{u}}$ and $\eta _{\underline{u},\delta }$ can be extended by zero to a function in $\mathcal C ^\infty (\mathbb R ^2)$.

On the other hand, for $\underline{z} \in \mathrm{supp}\,\eta _{\underline{u},\delta }$ we have

$$\begin{aligned} |\underline{z}|&= |\underline{r}^{ji}(\underline{u})+\underline{z}-\underline{r}^{ji}(\underline{u})| \le C_j |\mathbf{r}^j(\underline{r}^{ji}(\underline{u})+\underline{z})- \mathbf{r}^j (\underline{r}^{ji}(\underline{u}))|\\&= C_j |\mathbf{r}^j(\underline{r}^{ji}(\underline{u})+\underline{z})-\mathbf{r}^i(\underline{u})|\le C_j \epsilon _1\delta \end{aligned}$$

[we have used (4.7) in the final inequality]. Therefore we can fix $c>0$ independent of the particular chart number $j$ so that $\mathrm{supp}\,\eta _{\underline{u},\delta }\subset B(\underline{0},c\delta )$ for all $\underline{u} \in \varOmega ^{ij}_{2\delta _0}$. Moreover, from (2.10) we see that

$$\begin{aligned} \Vert \partial _{\underline{\alpha }} \eta _{\underline{u},\delta }\Vert _{L^\infty (\mathbb R ^2)}\le C_{\underline{\alpha }}\delta ^{-|\underline{\alpha }|} \qquad \forall \underline{u}\in \varOmega ^{ij}_{2\delta _0}. \end{aligned}$$

Thus, $\eta _{\underline{u},\delta }(\rho \underline{e}(\theta ))$ satisfies the conditions (4.4).

Finally, $\psi _{\underline{u},\delta }(\rho ,\theta )= (f_{\underline{u}}\eta _{{\underline{u}},\delta })(\rho \cos \theta ,\rho \sin \theta )$ and since, by Lemma 2, $f_{\underline{u}}(\rho \underline{e}(\theta ))$ is $\mathcal C ^\infty$ on the support of $\eta _{\underline{u},\delta }(\rho \underline{e}(\theta ))$, it follows that $\{\psi _{\underline{u},\delta }\}$ satisfies (4.4) with $c_0$ and $C_{\underline{\alpha }}$ independent of $\underline{u}$. In view of the inequalities (3.1) the same is true of the family $\{\chi _{\underline{u},\delta }\}$, and the result follows. $\square$

Proposition 5

There exists $C>0$ such that for all $N\in \mathbb N$ and all $\delta$, $0 < \delta \le \delta _0$,

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }\xi _N\Vert _{\infty ,I_2}\le C|\log N|^{3/2} \Vert \xi _N\Vert _0\qquad \forall \xi _N\in \mathbb T _N. \end{aligned}$$

Proof

By (3.6), it suffices to bound the values of $(\mathrm{A}^{ij}_{\mathrm{sing},\delta }\xi _N)(\underline{u})$ for $\underline{u}\in \varOmega ^{ij}_{2\delta }$. The polar coordinate form (2.22) gives, for $\underline{u} \in \varOmega ^{ij}_{2\delta }$,

$$\begin{aligned} (\mathrm{A}^{ij}_{\mathrm{sing},\delta }\xi _N)(\underline{u})&= \int _{-c_0\delta }^{c_0\delta } \int _0^{2\pi } \chi _{\underline{u},\delta }(\rho ,\theta ) \xi _N(\underline{r}^{ji}(\underline{u})+\rho \,\underline{e}(\theta ))\,\mathrm{d}\rho \mathrm{d} \theta \\&= \sum _{\underline{m} \in \mathbb Z _N^* } \widehat{\xi }_N(\underline{m}) e_{\underline{m}} (\underline{r}^{ji}(\underline{u})) \, \int _{-c_0\delta }^{c_0 \delta } \int _0^{2\pi } \chi _{\underline{u},\delta }(\rho ,\theta ) e_{\underline{m}}(\rho \,\underline{e}(\theta ))\mathrm{d} \rho \mathrm{d} \theta , \end{aligned}$$

where $\chi _{\underline{u},\delta }$ is given by (4.5). Therefore, by the Cauchy-Schwarz inequality we have

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }\xi _N\Vert _{\infty ,I_2} \le \Vert \xi _N\Vert _0 \Big ( \sum _{\underline{m}\in \mathbb Z _N^*} D_{\underline{m},\delta }^2\Big )^{1/2}, \end{aligned}$$

(4.8)

where

$$\begin{aligned} D_{\underline{m},\delta }:=\max _{\underline{u} \in \varOmega ^{ij}_{2\delta }} \int _0^{2\pi }\left| \,\,\int _{-c_0\delta }^{c_0\delta } \chi _{\underline{u},\delta }(\rho ,\theta ) \exp ( {2\pi }\mathrm{i}\rho \underline{m} \cdot \underline{e}(\theta ))\, \mathrm{d}\rho \right| \mathrm{d}\theta . \end{aligned}$$

As a direct consequence of Proposition 4 it follows that

$$\begin{aligned} \left| \,\, \int _{-c_0\delta }^{c_0\delta } \chi _{\underline{u}, \delta }(\rho ,\theta ) \exp ( {2\pi }\mathrm{i}\rho \,\underline{m}\cdot \underline{e}(\theta ) )\,\mathrm{d}\rho \right| \le C \delta , \end{aligned}$$

and, integrating by parts, that

$$\begin{aligned} \left| \,\,\int _{-c_0\delta }^{c_0\delta } \!\!\!\!\chi _{\underline{u}, \delta }(\rho ,\theta ) \exp ( {2\pi }\mathrm{i}\rho \,\underline{m}\cdot \underline{e}(\theta ) )\,\mathrm{d}\rho \right|&= \left| \mathrm{i} \int _{-c_0\delta }^{c_0\delta } \partial _\rho \chi _{\underline{u}, \delta }(\rho ,\theta )\, \frac{\exp ({2\pi }\mathrm{i}\rho \,\underline{m} \cdot \underline{e}(\theta ) )}{\underline{m} \cdot \underline{e}(\theta )}\mathrm{d}\rho \right| \\&\le C \frac{1}{|\underline{m}\cdot \underline{e}(\theta )|}, \end{aligned}$$

where the constant in both inequalities is independent of $\underline{u}$. Therefore

$$\begin{aligned} D_{\underline{m},\delta } \le C \int _0^{2\pi }\!\!\min \Big \{\delta ,\frac{1}{|m_1\cos \theta +m_2\sin \theta |}\Big \} \mathrm{d} \theta \le C^{\prime } \frac{\log N}{\sqrt{ m_1^2+m_2^2}},\quad \forall \underline{m} \in \mathbb Z _N \setminus \{ \underline{0}\}, \end{aligned}$$

where we have applied Lemma 4 for the last bound. Inserting this bound in (4.8) and using the fact that

$$\begin{aligned} \sum _{\underline{m}\in \mathbb Z _N^*\setminus \{ \underline{0}\}} \frac{1}{m_1^2+m_2^2} \le C \log N, \end{aligned}$$

(this can be proved by comparison with the integral of $1/(x^2+y^2)$) the result follows readily. $\square$

Lemma 4

For all $\underline{m}\in \mathbb Z ^2\setminus \{ \underline{0}\}$

$$\begin{aligned} \int _0^{2\pi }\min \Big \{ 1, \frac{1}{|m_1\cos \theta +m_2\sin \theta |}\Big \} \mathrm{d}\theta \le \frac{2\pi \big (1+\log \sqrt{m_1^2+m_2^2}\big )}{\sqrt{m_1^2+m_2^2}}. \end{aligned}$$

Proof

With some simple trigonometric arguments we prove that

$$\begin{aligned} \int _0^{2\pi }\!\!\min \{ 1, |m_1\cos \theta +m_2\sin \theta |^{-1}\}\mathrm{d} \theta&= \!\! \int _0^{2\pi }\!\!\min \{ 1,|\mathrm{Re}( (m_1-\mathrm{i}m_2) \exp (\mathrm{i}\theta )) |^{-1}\}\mathrm{d} \theta \\&= 4 \int _0^{\pi /2}\!\! \min \left\{ 1, \left( \sqrt{m_1^2+m_2^2}\,\sin \theta \right) ^{-1}\right\} \mathrm{d}\theta . \end{aligned}$$

If we now shorten $c:=\sqrt{m_1^2+m_2^2}$, we can estimate

$$\begin{aligned} \int _0^{\pi /2}\min \{ 1,(c\,\sin \theta )^{-1}\} \mathrm{d} \theta&= c^{-1} \int _0^{\pi /2}\min \{ c,(\sin \theta )^{-1}\}\mathrm{d} \theta \\&\le c^{-1}\int _0^{\pi /2}\min \left\{ c,\frac{\pi }{2\theta }\right\} \mathrm{d}\theta =\frac{\pi }{2c}(1+\log c), \end{aligned}$$

which finishes the proof. $\square$

5 Proofs of the main results

5.1 Inverse inequalities. Auxiliary approximation properties

In this subsection we collect some properties concerning the bivariate trigonometric polynomials $\xi _N\in \mathbb T _N$, which are needed for the analysis of the Nyström method under consideration.

From the definition of Sobolev norms (3.21) it is easy to establish the inverse inequalities

$$\begin{aligned} \Vert \xi _N\Vert _t \le (\sqrt{2}\, h)^{s-t}\Vert \xi _N\Vert _s \qquad \forall \xi _N \in \mathbb T _N, \qquad s \le t \end{aligned}$$

(5.1)

and

$$\begin{aligned} \Vert \xi _N\Vert _{\infty ,I_2}\!\le \! \sum _{\underline{m}\in \mathbb Z _N^*} |\widehat{\xi }_N(\underline{m})|\!\le \! N \bigg (\sum _{\underline{m}\!\in \!\mathbb Z _N^*} |\widehat{\xi }_N(\underline{m})|^2\bigg )^{ 1/2}\!=\! h^{-1 }\Vert \xi _N\Vert _{0} \quad \forall \xi _N \in \mathbb T _N.\qquad \end{aligned}$$

(5.2)

The operator $\mathrm{P}_N:\cup _s H^s\rightarrow \mathbb T _N$ that cuts off the tail of the Fourier series and at the same time gives the best $H^s$ approximation in $\mathbb T _N$ for all $s$ is given by

$$\begin{aligned} {\mathrm{P}}_N\xi :=\sum _{\underline{m}\in \mathbb Z _{N}^*} \widehat{\xi }(\underline{m}) e_{\underline{m}}. \end{aligned}$$

Recalling that $h=1/N$, it is easy to check that

$$\begin{aligned} \Vert \mathrm{P}_N\xi -\xi \Vert _{s}\le (2\,h)^{t-s}\Vert \xi \Vert _t \qquad \forall \xi \in H^t, \qquad s\le t. \end{aligned}$$

(5.3)

The interpolation operators introduced in (3.13) satisfy (cf. [23, Theorem 8.5.3])

$$\begin{aligned} \Vert \mathrm{Q}_N\xi -\xi \Vert _{s}\le C_{s,t} h^{t-s}\Vert \xi \Vert _{t}\qquad \forall \xi \in H^t, \qquad 0 \le s \le t,\qquad t> 1. \end{aligned}$$

(5.4)

The following lemma studies the uniform boundedness of $\mathrm{Q}_N$ as an operator from the space of continuous bivariate periodic functions to $H^0$. (See [21, Lemma 11.5] for a different proof of this result in the univariate case.)

Lemma 5

For all $N$

$$\begin{aligned} \Vert \mathrm{Q}_N\xi \Vert _{0} = \frac{1}{N}\bigg (\sum _{\underline{n}\in \mathbb Z _N} \left| \xi (\underline{x}_{\underline{n}}) \right| ^2\bigg )^{1/2}\le \Vert \xi \Vert _{\infty , I_2 } \qquad \forall \xi \in \mathcal C ^0(\overline{I}_2). \end{aligned}$$

(5.5)

Proof

Given $\xi _N \in \mathbb T _N$ we can write

$$\begin{aligned} \xi _N\left( {\textstyle \frac{2\pi \mathrm{i}}{N}}\,\underline{n}\right) = \sum _{\underline{m} \in \mathbb Z _N^*} \widehat{\xi }_N(\underline{m}) \exp ({\textstyle \frac{2\pi \mathrm{i}}{N}}\,(\underline{n}\cdot \underline{m})), \qquad \underline{n} \in \mathbb Z _N, \end{aligned}$$

and therefore, using the Parseval Theorem for the 2-dimensional discrete finite Fourier transform, we obtain

$$\begin{aligned} \Vert \xi _N\Vert _0^2=\sum _{\underline{m}\in \mathbb Z _N^*}|\widehat{\xi }_N(\underline{m})|^2 = \frac{1}{N^2} \sum _{\underline{n}\in \mathbb Z _N}|\xi _N\left( {\frac{2\pi \mathrm{i}}{N}\,\underline{n}}\right) |^2. \end{aligned}$$

When $\xi _N=\mathrm{Q}_N\xi$, this gives the equality in (5.5), whereas the inequality is a simple consequence of the fact that $\# \mathbb Z _N=N^2$. $\square$

We finally give an approximation result for the operator $\mathrm{D}_N$ defined in (3.18):

Lemma 6

There exists $C_{s,t}>0$ independent of $h$ such that

$$\begin{aligned} \Vert \mathrm{D}_N \mathrm{P}_N \xi -\xi \Vert _{s}\le C_{s,t} h^{t-s}\Vert \xi \Vert _t \qquad \forall \xi \in H^t, \qquad s\le t \le 0, \qquad s<-1. \end{aligned}$$

Proof

The result is proven in [11, Lemma 6] (see also [15]) for univariate 1-periodic functions (with the restriction $s<-1/2$ instead). The proofs can be easily adapted to the current case. $\square$

5.2 Error analysis of the polar coordinate quadrature rules

The results presented in this subsection, which lie at the heart of the main convergence proof presented in this paper, provides estimates on the quadrature error

$$\begin{aligned} \mathcal E _{h,k,\upgamma }(\psi ):=\left| \,\, \int _{-\infty }^\infty \int _0^{2\pi } \psi (\rho ,\theta ) \,\mathrm{d}\rho \,\mathrm{d}\theta - Q_{h,k,\upgamma }\psi \right| , \end{aligned}$$

(5.6)

(which are used in Sect. 5.3) for the family of polar-integration quadrature rules

$$\begin{aligned} Q_{h,k,\upgamma } \psi := k\sum _{p=0}^{\varTheta -1} \bigg [c(\theta _p)h \sum _{q=-\infty }^{\infty } \psi (\upgamma (\theta _p)+q c(\theta _p)h,\theta _p)\bigg ]\approx \int _{-\infty }^{\infty }\int _{0}^{2\pi } \psi (\rho ,\theta )\,\mathrm{d}\rho \,\mathrm{d}\theta \qquad \quad \end{aligned}$$

(5.7)

with integrands arising from certain trigonometric polynomials. The parameters in this quadrature formula are the positive integers $N$ and $\varTheta$ and the mesh-sizes $h=1/N$ and $k=2\pi /\varTheta$; the angular and radial quadrature nodes, in turn, are given by $\theta _p=pk$ for $p= 0,\dots ,\varTheta -1$ and $\upgamma (\theta _p)+q c(\theta _p)h$ for $q=-\infty ,\dots , \infty$ with $c(\theta ):=\min \{1/|\cos \theta |$, $1/| \sin \theta |\}$ respectively (see (2.28)). Finally, the function $\upgamma$ is a 2$\pi -$periodic piecewise continuous function. Equation (5.7) embodies trapezoidal quadrature rules in the variables $\theta$ and $\rho$, where the grid points used for integrating in $\rho$ depend on $\theta$. The error estimate (5.6) is applied in Propositions 7 and 8 on the $j$-th coordinate patch for each $j=1,\dots , J$—with $N=N_j$, $\varTheta =\varTheta _j$ and $h=h_j$.

Our main results concerning this family of rules are collected in the next theorem. For brevity our proof of this result is omitted here; all details in these regards can be found in [5, Appendix A].

Theorem 6

Let $\varTheta \approx N^{1+\alpha }$ with $\alpha >0$ and assume the sequence $\{ \chi _\delta \}_{\delta >0}$ satisfies the conditions listed in (4.4). Then, there exists $C>0$, independent of $\delta$, $h$ and $\upgamma$ such that, for all $\delta$ and $h$ satisfying $c_0\delta < 1/2$ and $h<\delta$, we have

$$\begin{aligned} \mathcal E _{h,k,\upgamma }(\chi _\delta \,\widetilde{\xi }_N) \le C \Big ( h\delta ^{-1}|\log h|+\delta \Big )\Vert \xi _N\Vert _{0}\qquad \forall \xi _N \in \mathbb T _N. \end{aligned}$$

(5.8)

In addition, for each $r \ge 2$ there exists $C_r>0$ such that

$$\begin{aligned} \mathcal E _{h,k,\upgamma }(\chi _\delta \,\widetilde{\xi }_N) \le C_r h^r\delta ^{1-r} \Vert \xi _N\Vert _r \qquad \forall \xi _N\in \mathbb T _N. \end{aligned}$$

(5.9)

5.3 Proofs of the results of Sect. 3.5

Recall that $N$ has been taken to be the largest of the discrete parameters $N_j$ and therefore $h_j\le h$ for all $j=1,\ldots ,J$. Also, by (2.16), we can bound $h_j^{-1}\le C h^{-1}$ and $|\log h_j|\le C |\log h|$ whenever needed.

Proposition 6

For all $t> 2$, there exists $C_t>0$ such that

$$\begin{aligned} \left\| \mathcal{A }^\mathrm{reg}_\delta -\mathcal Q _h\mathcal{A }^\mathrm{reg}_{\delta }\mathcal D _h\mathcal P _h\right\| _{\mathcal{H }^0\rightarrow \mathcal{H }^0}\le C_t \delta ^{-t} h ^{t}. \end{aligned}$$

Proof

Let us first consider the decomposition

$$\begin{aligned} \mathrm{A}^{ij}_{\mathrm{reg},\delta }-\mathrm{Q}_{N_i} \mathrm{A}^{ij}_{\mathrm{reg},\delta }\mathrm{D}_{N_j} \mathrm{P}_{N_j}&= \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\mathrm{I}-\mathrm{D}_{N_j} \mathrm{P}_{N_j})+(\mathrm{I}-\mathrm{Q}_{N_i})\mathrm{A}^{ij}_{\mathrm{reg},\delta }(\mathrm{D}_{N_j} \mathrm{P}_{N_j}-\mathrm{I})\\&+ (\mathrm{I}-\mathrm{Q}_{N_i})\mathrm{A}^{ij}_{\mathrm{reg},\delta }. \end{aligned}$$

By Proposition 2 and Lemma 6, for all $t>2$,

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\mathrm{I}-\mathrm{D}_{N_j}\mathrm{P}_{N_j})\xi \Vert _{0}\le C_t \delta ^{-t}\Vert (\mathrm{I}-\mathrm{D}_{N_j}\mathrm{P}_{N_j})\xi \Vert _{-t}\le C^{\prime }_t h^t\delta ^{-t}\Vert \xi \Vert _0. \end{aligned}$$

The second term is bounded using (5.4), Proposition 2 and Lemma 6:

$$\begin{aligned} \Vert (\mathrm{I}-\mathrm{Q}_{N_i}) \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\mathrm{D}_{N_j} \mathrm{P}_{N_j}-\mathrm{I}) \xi \Vert _{0}&\le C_t h^{t/2} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\mathrm{D}_{N_j}\mathrm{P}_{N_j}-\mathrm{I}) \xi \Vert _{t/2}\\&\le C_t^{\prime } h^{t/2} \delta ^{-t}\Vert (\mathrm{D}_{N_j}\mathrm{P}_{N_j}-\mathrm{I}) \xi \Vert _{-t/2}\le C_t^{\prime \prime }\delta ^{-t} h^t\Vert \xi \Vert _{0}. \end{aligned}$$

Finally, by (5.4), Proposition 2 and Lemma 6,

$$\begin{aligned} \Vert (\mathrm{I}-\mathrm{Q}_{N_i})\mathrm{A}^{ij}_{\mathrm{reg},\delta }\xi \Vert _{0}&\le C_t h^t\Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\xi \Vert _{t}\le C^{\prime }_t h^t\delta ^{-t}\Vert \xi \Vert _{0}, \end{aligned}$$

and the proof is finished. $\square$

Proposition 7

There exists $C>0$ such that for all $\xi _{N_j}\in \mathbb T _{N_j}$ and $\delta >h$

$$\begin{aligned} \sup _{\underline{u}\in \varOmega ^{ij}_{\delta /2}}\Big | (\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{\mathrm{sing},\delta ,h})\xi _{N_j}(\underline{u}) \Big |&\le C\Big (\delta ^{-1} h_j|\log h_j|+\delta \Big )\Vert \xi _{N_j}\Vert _{0},\qquad \end{aligned}$$

(5.10)

$$\begin{aligned} \sup _{\underline{u}\in \varOmega ^{ij}_{2\delta }\setminus \varOmega ^{ij}_{\delta /2}}\Big | (\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{\mathrm{sing},\delta ,h})\xi _{N_j}(\underline{u}) \Big |&\le C|\log h_j|^{3/2}\Vert \xi _{N_j}\Vert _{0},\end{aligned}$$

(5.11)

$$\begin{aligned} \sup _{\underline{u}\in I_2\setminus \varOmega ^{ij}_{2\delta }}\Big |( \mathrm{A}^{ij}_{\mathrm{sing},\delta }-\mathrm{A}^{ij}_{\mathrm{sing},\delta ,h})\xi _{N_j} (\underline{u})\Big |&= 0. \end{aligned}$$

(5.12)

Proof

Recalling the definition of the discrete operator $\mathrm{A}^{ij}_{\mathrm{sing},\delta ,h}=(\omega ^i\circ \mathbf{r}^i) \mathrm{L}^{ij}_{\delta ,h}$ [the operator $\mathrm{L}^{ij}_{\delta ,h}$ is defined in (2.29)], it is clear that

$$\begin{aligned} (\mathrm{A}^{ij}_{\mathrm{sing},\delta ,h}\xi _{N_j})(\underline{u})=0 \qquad \underline{u}\in I_2\setminus \varOmega ^{ij}_{\delta /2}. \end{aligned}$$

(5.13)

Using also (3.6), it is obvious that both operators in (5.12) vanish for $\underline{u} \in I_2\setminus \varOmega ^{ij}_{2\delta }$. Also (5.11) is a simple consequence of (5.13) and Proposition 5.

For $\underline{u}\in \varOmega ^{ij}_{\delta /2}$, recalling that $\widetilde{\xi }_{N_j}(\rho ,\theta )=\xi _{N_j}(\underline{r}^{ji}(\underline{u})+ \rho \underline{e}(\theta ))$,

$$\begin{aligned} (\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{\mathrm{sing},\delta ,h}) \xi _{N_j} (\underline{u})= \int _{-\infty }^{\infty }\int _0^{2\pi }\! (\chi _{\underline{u},\delta } \widetilde{\xi }_{N_j})(\rho ,\theta )\mathrm{d}\rho \,\mathrm{d}\theta -Q_{h,k,\gamma }(\chi _{\underline{u},\delta }\widetilde{\xi }_{N_j}), \end{aligned}$$

where $\chi _{\underline{u},\delta }$ is defined by (4.5) and the quadrature rule $Q_{h,k,\gamma }$ is given in (5.7). By Proposition 4 we can apply Theorem 6, which estimates the above quadrature error in terms of the constants that appear in (4.4). Since these constants can be taken to be independent of $\underline{u} \in \varOmega ^{ij}_{\delta /2}$ [this is part of the assertion of Proposition 4], then (5.10) follows readily. $\square$

In the following sequence of results we will use the geometric cut-off operator

$$\begin{aligned} \mathcal G \varvec{\xi }= (\mathrm G_1\xi _1,\ldots ,\mathrm G_J \xi _J):= \big ( (1-\widetilde{\omega }_0^1) \xi _1,\ldots ,(1-\widetilde{\omega }_0^J)\xi _J\big ), \end{aligned}$$

(5.14)

where, as a reminder, $\widetilde{\omega }_0^j$ is the characteristic function of the domain $\varOmega ^j$.

Proposition 8

For all $t\ge 2$, there exists $C_t>0$ such that for all $\xi _{N_j}\in \mathbb T _{N_j}$ and $\delta > h$

$$\begin{aligned} \Vert (\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\xi _{N_j}\Vert _{\infty ,I_2}\le C_t\delta ^{1-t}h_j^{t}\Vert \xi _{N_j}\Vert _{t}+C \delta \Vert \mathrm{G}_j \xi _{N_j}\Vert _{\infty ,I_2}. \end{aligned}$$

(5.15)

Proof

Following the argument of the preceding proof, but using (5.9) of Theorem 6, it follows that for any integer $t\ge 2$,

$$\begin{aligned} \sup _{\underline{u}\in \varOmega ^{ij}_{\delta /2}}\Big | (\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\xi _{N_j}(\underline{u}) \Big |\le C_t \delta ^{1-t}h^t\Vert \xi _{N_j}\Vert _t. \end{aligned}$$

(5.16)

From the definition (2.8), it follows that if $\mathbf{r}\not \in S^j_{\delta /2}$, then $\overline{B(\mathbf{r},\epsilon _1\delta )}\cap \mathrm{supp}\omega ^j =\emptyset$. Therefore, if $\underline{u}\not \in \varOmega ^{ij}_{\delta /2}$,

$$\begin{aligned} \mathrm{supp}\,K^{ij}_{{\mathrm{sing}},\delta }(\underline{u},\,\cdot \,)\cap \varOmega ^j&\subset \mathrm{supp}\,\eta _{\delta }(\mathbf{r}^i(\underline{u}),\mathbf{r}^j(\,\cdot \,)) \cap \varOmega ^j\\&\subset (\mathbf{r}^j)^{-1}\big (\overline{B(\mathbf{r}^i(\underline{u}),\epsilon _1\delta )}\cap \mathrm{supp}\omega ^j)=\emptyset , \end{aligned}$$

where we have used (2.11). This means that if $\underline{u}\in \varOmega ^{ij}_{2\delta }\setminus \varOmega ^{ij}_{\delta /2}$, only the value of $\xi _{N_j}$ on $I_2\setminus \varOmega ^j$ (where $\xi _{N_j}=\mathrm{G}_j \xi _{N_j}$) is relevant. We then change to polar coordinates as in the proof of Proposition 5 to obtain:

$$\begin{aligned} |(\mathrm{A}^{ij}_{{\mathrm{sing}},\delta }-\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\xi _{N_j}(\underline{u})|&= |(\mathrm{A}^{ij}_{{\mathrm{sing}},\delta }\xi _{N_j})(\underline{u})|\nonumber \\ \nonumber&\le \Vert \xi _{N_j}\Vert _{\infty ,I_2\setminus \varOmega ^j} \int _{-c_0\delta }^{c_0\delta } \int _0^{2\pi }|\chi _{\underline{u},\delta }(\rho ,\theta )|\,\mathrm{d}\rho \,\mathrm{d}\theta \nonumber \\&\le C \delta \Vert \mathrm{G}_j\xi _{N_j}\Vert _{\infty ,I_2}. \end{aligned}$$

(5.17)

(In the last inequality we have applied Proposition 4, according to which the constants $c_0$ and $C$ do not depend on $\underline{u}$). Equation (5.15) now follows from Eqs. (5.12), (5.16) and (5.17). $\square$

Proposition 9

For all $\varepsilon >0$ there exists $C_\varepsilon >0$ independent of $h$ and $\delta >h$ such that

$$\begin{aligned} \left\| \mathcal{A }^{\mathrm{sing}}_\delta -\mathcal Q _h\mathcal{A }^{{\mathrm{sing}}}_{\delta ,h}\mathcal P _h\right\| _{\mathcal{H }^0\rightarrow \mathcal{H }^0}\le C_\varepsilon \Big (\delta ^{-1-\varepsilon }h+ \delta ^{-1}h\,|\log h| +\delta ^{1/2}|\log h|^{3/2} \Big ).\quad \quad \quad \end{aligned}$$

(5.18)

Furthermore, for any integer $r\ge 2$,

$$\begin{aligned}&\big \Vert \big (\mathcal{A }^{\mathrm{sing}}_\delta -\mathcal Q _h\mathcal{A }^{{\mathrm{sing}}}_{\delta ,h}\mathcal P _h\big )\varvec{\xi }\big \Vert _{0}\nonumber \\&\quad \le C_r \big (\delta ^{-r} h^r\Vert \varvec{\xi }\Vert _{r} +\delta \Vert \mathcal P _h\varvec{\xi }-\varvec{\xi }\Vert _{\infty ,I_2}+\delta \Vert \mathcal G \varvec{\xi }\Vert _{\infty ,I_2}\Big ) \quad \forall \varvec{\xi }\in \mathcal H ^0. \end{aligned}$$

(5.19)

Proof

We start with the decomposition

$$\begin{aligned} \mathrm{A}_{\mathrm{sing},\delta }^{ij}-\mathrm{Q}_{N_i} \mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h}\mathrm{P}_{N_j}&= \mathrm{A}^{ij}_{\mathrm{sing},\delta }(\mathrm{I}-\mathrm{P}_{N_j})+({\mathrm{I}}-\mathrm{Q}_{N_i})\mathrm{A}^{ij}_{\mathrm{sing},\delta }\mathrm{P}_{N_j}\nonumber \\&+\mathrm{Q}_{N_i}(\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\mathrm{P}_{N_j}. \end{aligned}$$

(5.20)

For all $t\ge 1$, Proposition 2 with $r=-1$ and (5.3) imply

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }(\mathrm{I}-\mathrm{P}_{N_j})\xi \Vert _{0}\le C \delta ^{-1} \Vert \xi -\mathrm{P}_{N_j}\xi \Vert _{-1}\le C_t \delta ^{-1}h^{t}\Vert \xi \Vert _{t-1}. \end{aligned}$$

(5.21)

For $t>1$, by (5.4), Proposition 2 with $r=t-1$ and the inverse inequality (5.1), we can bound

$$\begin{aligned} \Vert (\mathrm{I}-\mathrm{Q}_{N_i})\mathrm{A}^{ij}_{\mathrm{sing},\delta }\mathrm{P}_{N_j}\xi \Vert _{0}&\le C_t h^t\Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }\mathrm{P}_{N_j}\xi \Vert _{t}\le C_t^{\prime } \delta ^{-t}\,h^t \Vert \mathrm{P}_{N_j}\xi \Vert _{t-1}\end{aligned}$$

(5.22)

$$\begin{aligned}&\le C^{\prime }_t \delta ^{-t}h\Vert \mathrm{P}_{N_j}\xi \Vert _{0}\le C^{\prime }_t \delta ^{-t}h\Vert \xi \Vert _0. \end{aligned}$$

(5.23)

For fixed $\xi$ we introduce the auxiliary function $g:=(\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\mathrm{P}_{N_j}\xi$, so that by Lemma 5,

$$\begin{aligned} \Vert \mathrm{Q}_{N_i}(\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\mathrm{P}_{N_j}\xi \Vert ^2_{0} = \Vert \mathrm{Q}_{N_i}g\Vert _0^2= h_i^{2}\sum _{\underline{m}\in \mathbb Z _{N_i}}|g(\underline{x} ^ i_{\underline{m}})|^2. \end{aligned}$$

(5.24)

We can then split $\mathbb Z _{N_i}$ as follows:

$$\begin{aligned} A_{N_i}&:= \big \{\underline{m}\in \mathbb Z _{N_i}\, :\, \underline{x}_{\underline{m}}^i\not \in \varOmega ^{ij}_{2\delta } \big \},\\ B_{N_i}&:= \big \{\underline{m}\in \mathbb Z _{N_i}\, :\, \underline{x}_{\underline{m}}^i\in \varOmega ^{ij}_{2\delta } \setminus \varOmega ^{ij}_{\delta /2} \big \},\\ C_{N_i}&:= \big \{\underline{m}\in \mathbb Z _{N_i}\, :\, \underline{x}_{\underline{m}}^i\in \varOmega ^{ij}_{\delta /2} \big \}. \end{aligned}$$

A key point is the fact that $\# B_{N_i}\le C \delta N_i^2$, which is proved in Lemma 7. From (5.12) it follows that $g(\underline{x}_{\underline{m}}^i)=0$ for every $\underline{m}\in A_{N_i}$. Also, by (5.10)

$$\begin{aligned} \Bigg ( h_i^2\sum _{\underline{m}\in C_{N_i}}|g(\underline{x}^i_{\underline{m}})|^2\Bigg )^{1/2}\le C\Big (\delta ^{-1} h_j|\log h_j|+\delta \Big )\Vert \mathrm{P}_{N_j}\xi \Vert _{0}. \end{aligned}$$

Finally, from (5.11), we obtain

$$\begin{aligned} \Bigg ( h_i^2\sum _{\underline{m}\in B _{N_i}}|g(\underline{x}^i_{\underline{m}})|^2\Bigg )^{1/2}\le C \delta ^{1/2} |\log h_j|^{3/2} \Vert \mathrm{P}_{N_j}\xi \Vert _{0}. \end{aligned}$$

Going back to (5.24), we have proved that

$$\begin{aligned} \Vert \mathrm{Q}_{N_i}(\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\mathrm{P}_{N_j}\xi \Vert _{0}\le C\Big (\delta ^{-1} h|\log h|+\delta ^{1/2}|\log h|^{3/2} \Big )\Vert \mathrm{P}_{N_j}\xi \Vert _{0}.\quad \quad \quad \end{aligned}$$

(5.25)

A direct estimate from (5.24), using (5.15) now, gives

$$\begin{aligned} \Vert \mathrm{Q}_{N_i}(\mathrm{A}^{ij}_{\mathrm{sing},\delta } -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h})\mathrm{P}_{N_j}\xi \Vert _{0}\le C_r\delta ^{1-r}h^{r}\Vert \xi \Vert _{r}+C \delta \Vert \mathrm{G}_j \mathrm{P}_{N_j}\xi \Vert _{\infty ,I_2}. \end{aligned}$$

(5.26)

Using the decomposition (5.20) and the bounds (5.21) with $t=1$, (5.23) with $t=1+\varepsilon$ and (5.25) we can easily prove (5.18). At the same time, (5.19) follows from the same decomposition, using now (5.21) with $t=r+1$, (5.22) with $t=r$ and (5.26). In a last step we use the fact that the cut-off operators $\mathrm G^j$ are bounded in $L^\infty (I_2)$. This finishes the proof of the Proposition. $\square$

Lemma 7

There exists $C$ independent of $\delta >h$ such that

$$\begin{aligned} \# \big \{\underline{m} \in \mathbb Z _{N_i}\, :\, \underline{x}_{\underline{m}}^i \in \varOmega ^{ij}_{2\delta } \setminus \varOmega ^{ij}_{\delta /2}\big \} \le C\delta N_i^2. \end{aligned}$$

Proof

We start by setting $\varGamma ^j:=\partial (\mathrm{supp}\omega ^j)$, which has been assumed to be the finite union of Lipschitz arcs. Therefore, for every $\delta >0$, we can pick up a finite set $\varGamma _\delta ^j \subset \varGamma ^j$, with the following properties

$$\begin{aligned} \varGamma \subset \bigcup _{\mathbf{r}\in \varGamma _\delta ^j} \overline{B(\mathbf{r},\epsilon _1\delta )} \qquad \text{ with}\,\, \#\varGamma _\delta ^j \le C \delta ^{-1}. \end{aligned}$$

The constant $C$ that controls the size of the discrete set $\varGamma _\delta ^j$ depends on the Lipschitz constants related to $\varGamma ^j$. On the other hand, for $\delta$ small enough

$$\begin{aligned} S^{j}_{2\delta }\setminus (\mathrm{supp}\omega ^j)\subset \bigcup _{{\mathbf{r}}\in \varGamma } \overline{B({\mathbf{r}},4\epsilon _1\delta )} \end{aligned}$$

and therefore

$$\begin{aligned} S^j_{2\delta }\setminus S^j_{\delta /2}\subset S^j_{2\delta }\setminus (\mathrm{supp}\omega ^j) \subset \bigcup _{\mathbf{r}\in \varGamma ^j_\delta } \overline{B(\mathbf{r}, 5\epsilon _1\delta )}. \end{aligned}$$

(5.27)

We now go back to the parametric domain using $(\mathbf{r}^i)^{-1}$ to define the sets

$$\begin{aligned} \mathcal O _{\mathbf{r}}:= (\mathbf{r}^i)^{-1} (S^i\cap \overline{B(\mathbf{r}, 5\epsilon _1\delta )}) \qquad \mathbf{r}\in \varGamma _\delta ^j. \end{aligned}$$

Since $|\underline{u}-\underline{v}|\le C^i|\mathbf{r}^i(\underline{u})-\mathbf{r}^i(\underline{v})|$ for all $\underline{u},\underline{v} \in D^i$, each set $\mathcal O _{\mathbf{r}}$ is either empty or can be surrounded by a closed ball of radius $5 C^i\epsilon _1\delta$. This gives a collection of points $\{ \underline{u}_1^\delta ,\ldots ,\underline{u}_{k_\delta }^\delta \}\subset \mathbb R ^2$ with $k_\delta \le C \delta ^{-1}$ such that

$$\begin{aligned} \varOmega ^{ij}_{2\delta }\setminus \varOmega ^{ij}_{\delta /2}&= (\mathbf{r}^i)^{-1} (S^{ij}_{2\delta }\setminus S^{ij}_{\delta /2}) = (\mathbf{r}^i)^{-1} (S^i \cap (S^j_{2\delta }\setminus S^j_{\delta /2}))\\&\subset \bigcup _{\mathbf{r}\in \varGamma ^j_\delta } \mathcal O _\mathbf{r}\subset \bigcup _{\ell =1}^{k_\delta } \big \{ \underline{v} \in \mathbb R ^2\, :\, |\underline{v} -\underline{u}_\ell ^\delta |\le 5 C^i\epsilon _1\delta \big \}. \end{aligned}$$

The number of points of the uniform grid $\{ \underline{x}_{\underline{m}}^i\,:\, \underline{m} \in \mathbb Z _{N_i}\}$ that fit in a closed disk of radius $r$ is bounded by $(2 r N_i+1)^2$. Therefore, the cardinal of the intersection of the uniform grid with $\varOmega ^{ij}_{2\delta }\setminus \varOmega ^{ij}_{\delta /2}$ can be bounded by

$$\begin{aligned} k_\delta (10 C^i \epsilon _1\delta \,N_i+1)^2 \le C \delta ^{-1} (10 C^i \delta \,N_i+1)^2 \le C^{\prime } \delta N_i^2, \end{aligned}$$

where we have applied that $1\le \delta N_i$ (which is implied by $\delta > h$), and the result follows. $\square$

Proof of Theorem 4

This result follows now from an adequate choice of parameters in Propositions 6 and 9. Since $\delta \approx h^\beta$, then $\delta ^{-1} h \approx h^{1-\beta }$. We choose $\varepsilon =(1-\beta )/(2\beta )$, so that $\delta ^{-1-\varepsilon } h \approx h^{(1-\beta )/2}$. We then apply Proposition 6 with $t=3$ and the first bound of Proposition 9 with the above $\varepsilon$. This yields the bound

$$\begin{aligned} \Vert \mathcal{B }_{\delta ,h}-\mathcal{B }_\delta \Vert _{\mathcal{H }^0\rightarrow \mathcal{H }^0}\le C_\varepsilon (h^{3(1-\beta )}+ h^{(1-\beta )/2}+ h^{1-\beta }|\log h| +h^{\beta /2}|\log h|^{3/2}), \end{aligned}$$

which ensures convergence as $h \rightarrow 0$. The uniform bound for the inverses of $\mathcal B _\delta$ (Theorem 3) ensures then a uniform bound for the inverses of the operators $\mathcal B _{\delta ,h}$ when $h$ is small enough and $\delta \approx h^\beta$. $\square$

Proof of Theorem 5

For $\psi \in H^r(S)$, let us define $\varvec{\phi }:= ((\omega ^1\psi )\circ \mathbf{r}^1,\ldots ,(\omega ^J\psi )\circ \mathbf{r}^J) \in \mathcal H ^r$. By construction (cf. (5.14)), $\mathcal G \varvec{\phi }=0$. Applying Proposition 6 and the second bound of Proposition 9 it follows that for any integer $r\ge 3$

$$\begin{aligned} \Vert (\mathcal{B }_{\delta ,h}-\mathcal{B }_{\delta } )\varvec{\phi }\Vert _0 \le C_r\big ( h^{r(1-\beta )} \Vert \varvec{\phi }\Vert _r+ h^\beta \Vert \mathcal P _h\varvec{\phi }-\varvec{\phi }\Vert _ {\infty ,I_2} \big ). \end{aligned}$$

(5.28)

Taking now

$$\begin{aligned} r\ge \max \{1/\beta ,\beta /(1-\beta )\}, \end{aligned}$$

(5.29)

so that $1+\beta \le 1+r(1-\beta )\le r$, from Sobolev’s embedding theorem and the approximation estimate (5.3) we obtain

$$\begin{aligned} h^\beta \Vert \mathcal P _h\varvec{\phi }-\varvec{\phi }\Vert _{\infty ,I_2}&\le C_\beta h^\beta \Vert \mathcal P _h\varvec{\phi }-\varvec{\phi }\Vert _{1+\beta }\le C_{r,\beta } h^{r(1-\beta )}\Vert \varvec{\phi }\Vert _{1+r(1-\beta )}\\&\le C_{r,\beta } h^{r(1-\beta )}\Vert \varvec{\phi }\Vert _r. \end{aligned}$$

Applying this estimate in (5.28) we deduce that for every integer $r\ge 3$ satisfying (5.29), we have

$$\begin{aligned} \Vert (\mathcal{B }_{\delta ,h}-\mathcal{B }_{\delta } )\varvec{\phi }\Vert _0 \le C_r h^{r(1-\beta )}\Vert \varvec{\phi }\Vert _r. \end{aligned}$$

(5.30)

Since $\beta$ is fixed, the dependence of the various bounding constants on the parameter $\beta$ is dropped from the notation in what follows (recall Remark 5).

For $\psi \in H^r(S)$ we now define $\mathcal M _{\delta ,h}\psi :=(\mathcal B _{\delta ,h}-\mathcal B _\delta )\varvec{\phi }$ with $\varvec{\phi }$ as above; clearly $\mathcal M _{\delta ,h}: H^r(S)\rightarrow \mathcal H ^0$ is a continuous map. In view of (3.22) Eq. (5.30) can be re-expressed in the form

$$\begin{aligned} \Vert \mathcal M _{\delta ,h}\psi \Vert _0 \le C_r h^{r(1-\beta )}\Vert \psi \Vert _{H^r(S)} \qquad \forall \psi \in H^r(S). \end{aligned}$$

(5.31)

Since $\Vert \mathcal B _{\delta ,h}-\mathcal B _\delta \Vert _\mathcal{H ^0\rightarrow \mathcal H ^0}\rightarrow 0$, it follows that the sequence $\mathcal M _{\delta ,h}$ is uniformly bounded,

$$\begin{aligned} \Vert \mathcal M _{\delta ,h}\psi \Vert _0 \le C \Vert \psi \Vert _{H^0(S)}, \end{aligned}$$

and therefore, using Sobolev interpolation theory  [22, Appendix B] for the operator $\mathcal M _{\delta ,h}$, it follows that (5.31), and therefore (5.30), hold for all $r\ge 0$.

Let now $\varvec{\phi }$ and $\varvec{\phi }_h$ be the respective solutions of $\mathcal B _\delta \varvec{\phi }=\mathbf{U}$ and $\mathcal B _{\delta ,h}\varvec{\phi }_h=\mathcal Q _h\mathbf{U}$. From Theorem 4 it follows that for $h$ small enough, the inverse of $\mathcal{B }_{\delta ,h}:\mathcal{H }^0\rightarrow \mathcal{H }^0$ is uniformly bounded with respect to $h$. Hence

$$\begin{aligned} \Vert \varvec{\phi }-\varvec{\phi }_h\Vert _0&\le C \Vert \mathcal{B }_{\delta ,h}(\varvec{\phi }-\varvec{\phi }_h)\Vert _0\le C \Vert \mathcal{B }_{\delta } \varvec{\phi }-\mathcal{B }_{\delta ,h}\varvec{\phi }_h\Vert _0+ C\Vert (\mathcal{B }_{\delta ,h}-\mathcal{B }_{\delta } )\varvec{\phi }\Vert _0\nonumber \\&= C \Vert \mathbf{U} -\mathcal Q _h\mathbf{U}\Vert _0 +C \Vert (\mathcal{B }_{\delta ,h}-\mathcal{B }_{\delta } )\varvec{\phi }\Vert _0. \end{aligned}$$

(5.32)

But, from (5.4), it follows that for $t>1$

$$\begin{aligned} \Vert \mathbf{U} -\mathcal Q _h\mathbf{U}\Vert _0\le C_t h^t\Vert \mathbf{U}\Vert _{t}=C_t h^t\Vert \mathcal{B }_{\delta _0} \varvec{\phi }\Vert _{t}\le C^{\prime }_t h^t\Vert \varvec{\phi }\Vert _t. \end{aligned}$$

(5.33)

Thus, substituting (5.30) and (5.33) in (5.32), we conclude that for all $t>1$ there exists $C_t>0$ such that

$$\begin{aligned} \Vert \varvec{\phi }-\varvec{\phi }_h\Vert _0\le C_t h^{t(1-\beta )}\Vert \varvec{\phi }\Vert _{t}. \end{aligned}$$

(5.34)

Error estimates in stronger Sobolev norms can now be obtained by means of the inverse inequalities (5.1) together with the fact that $\mathcal P _h$ provides the best approximation, for all periodic Sobolev norms, in the discrete subspaces of trigonometric polynomials. Indeed, for $s>0$ we have

$$\begin{aligned} \Vert \varvec{\phi }\!-\!\varvec{\phi }_h\Vert _s \le \Vert \varvec{\phi }\!-\! \mathcal P _h\varvec{\phi }\Vert _s + C h^{-s} \Vert \mathcal P _h\varvec{\phi }\!-\!\varvec{\phi }_h\Vert _0 \le \Vert \varvec{\phi }\!-\! \mathcal P _h\varvec{\phi }\Vert _s + C h^{-s} \Vert \varvec{\phi }-\varvec{\phi }_h\Vert _0.\nonumber \\ \end{aligned}$$

(5.35)

The proof now follows by substituting (5.3) and (5.34) in (5.35). $\square$

5.4 Proof of Theorem 1

In view of the reconstruction formula of Theorem 2 and since $\mathcal B _\delta \varvec{\phi }=\mathbf{U}$, the exact solution of (2.2) can be expressed in the form

$$\begin{aligned} \psi (\mathbf{r})&= \sum _{i\in \mathcal I (\mathbf{r})} (\omega ^i \psi ^i)\left( (\mathbf{r}^i)^{-1}(\mathbf{r})\right) =\sum _{i\in \mathcal I (\mathbf{r})} \phi ^i\left( (\mathbf{r}^i)^{-1}(\mathbf{r})\right) \\&= \sum _{i\in \mathcal I (\mathbf{r})} 2\bigg (U^i-\sum _{j=1}^J\mathrm{A}^{ij}_{\mathrm{reg},\delta }\phi ^j-\sum _{j=1}^J \mathrm{A}^{ij}_\mathrm{{\mathrm{sing}},\delta }\phi ^j\bigg )\left( (\mathbf{r}^i)^{-1}(\mathbf{r})\right) . \end{aligned}$$

Similarly, for the discrete Nyström solution we have the relation

$$\begin{aligned} \psi _h(\mathbf{r})&= \sum _{i\in \mathcal I (\mathbf{r})} (\omega ^i \psi ^i_h)\left( (\mathbf{r}^i)^{-1}(\mathbf{r})\right) \\&= \sum _{i\in \mathcal I (\mathbf{r})} 2\Bigg (U^i-\sum _{j=1}^J\mathrm{A}^{ij}_{\mathrm{reg},\delta }\mathrm{D}_{N_j}\phi _h^j-\sum _{j=1}^J \mathrm{A}^{ij}_{ {\mathrm{sing}},\delta ,h}\phi _h^j\Bigg )\left( (\mathbf{r}^i)^{-1}(\mathbf{r})\right) , \end{aligned}$$

which follows by re-expressing Eqs. (2.34) and (2.35) in terms of $\phi _h^j=\mathrm{Q}_{N_j}((\omega ^j\circ \mathbf{r}^j)\psi ^j_h)$.

Since $\Vert \psi \Vert _{H^s(S)}=\Vert \varvec{\phi }\Vert _s$, by Propositions 10 and 11 below we thus have

$$\begin{aligned} \Vert \psi -\psi _h\Vert _{L^2(S)}&\le C\max _{i,j} \Big ( \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\mathrm{D}_{N_j}\phi ^j_h-\mathrm{A}^{ij}_{\mathrm{reg},\delta }\phi ^j\Vert _0+ \Vert \mathrm{A}^{ij}_{{\mathrm{sing}},\delta }\phi ^j_h-\mathrm{A}^{ij}_{{\mathrm{sing}},\delta }\phi ^j\Vert _0\Big )\\&\le C_t \Big ( h^{t-1}+|\log h|^{3/2} h^{t(1-\beta )}\Big ) \Vert \varvec{\phi }\Vert _t\\&= C_t \big ( h^{t-1}+|\log h|^{3/2} h^{t(1-\beta )}\big ) \Vert \psi \Vert _{H^t(S)}, \end{aligned}$$

and (2.36) follows. Note that the constants $C_t$ depend on $\beta$ (see Remark 5). The $L^\infty (S)$ bound (2.37) follows from similar arguments, and the proof of Theorem 1 is thus complete. $\square$

Proposition 10

For all $t>1$ there exists $C_t>0$ such that

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\mathrm{D}_{N_j}\phi _h^j-\mathrm{A}^{ij}_{\mathrm{reg},\delta }\phi ^j\Vert _0 \le C_t h^{t(1-\beta )}\Vert \varvec{\phi }\Vert _{t}. \end{aligned}$$

(5.36)

In addition, for all $\varepsilon >0$ and $t\ge 2$ there exists $C_{\varepsilon ,t}>0$ such that

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\mathrm{D}_{N_j}\phi _h^j-\mathrm{A}^{ij}_{\mathrm{reg},\delta }\phi ^j\Vert _{\infty ,I_2}\le C_{\varepsilon ,t}h^{t(1-\beta )-(1+\varepsilon )\beta } \Vert \varvec{\phi }\Vert _{t}. \end{aligned}$$

(5.37)

Proof

Using successively Proposition 2 (recall that $\delta \approx h^{\beta }$), Lemma 6 and Theorem 5 we derive the first estimate:

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\mathrm{D}_{N_j}\phi _h^j-\mathrm{A}^{ij}_{\mathrm{reg},\delta }\phi ^j\Vert _0&\le \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\mathrm{D}_{N_j}\phi _h^j-\phi _h^j)\Vert _{0}+\Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\phi _h^j-\phi ^j)\Vert _{0}\nonumber \\&\le C_t h^{-t\beta }\Vert \mathrm{D}_{N_j}\phi _h^j-\phi _h^j\Vert _{-t}+C_0 \Vert \phi _h^j-\phi ^j\Vert _{0}\nonumber \\&\le C_t^{\prime } h^{t(1-\beta )}(\Vert \phi _h^j\Vert _0 +\Vert \phi ^j\Vert _{t})\le C_t^{\prime \prime }h^{t(1-\beta )}\Vert \varvec{\phi }\Vert _t. \quad \quad \quad \end{aligned}$$

(5.38)

To derive an estimate in the $L^\infty (I_2)$ norm, we proceed similarly:

$$\begin{aligned}&\Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }\mathrm{D}_{N_j}\phi _h^j-\mathrm{A}^{ij}_{\mathrm{reg},\delta }\phi ^j\Vert _{\infty ,I_2}\nonumber \\&\quad \le C_{\varepsilon }\left( \Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\mathrm{D}_{N_j}\phi _h^j-\phi _h^j)\Vert _{1+\varepsilon }+\Vert \mathrm{A}^{ij}_{\mathrm{reg},\delta }(\phi _h^j-\phi ^j)\Vert _{1+\varepsilon }\right) \nonumber \\&\quad \le C_{\varepsilon ,t}\Big ( h^{-(1+\varepsilon +t)\beta }\Vert \mathrm{D}_{N_j}\phi _h^j-\phi _h^j\Vert _{-t}+ h^{-(1+\varepsilon )\beta } \Vert \phi _h^j-\phi ^j\Vert _{0}\Big )\nonumber \\&\quad \le C_{\varepsilon ,t}^{\prime }\Big ( h^{t(1-\beta )-(1+\varepsilon )\beta } \Vert \phi _h^j\Vert _0 + h^{t(1-\beta )-(1+\varepsilon )\beta } \Vert \phi ^j\Vert _{t}\Big )\nonumber \\&\quad \le C_{\varepsilon ,t}^{\prime \prime }h^{t(1-\beta )-(1+\varepsilon )\beta } \Vert \varvec{\phi }\Vert _t, \end{aligned}$$

(5.39)

where we applied sequentially Proposition 2, Lemma 6 and Theorem 5. $\square$

Proposition 11

For all $t\ge 2$, there exists $C_t>0$ such that

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h}\phi _h^j-\mathrm{A}^{ij}_{ {\mathrm{sing}},\delta }\phi ^j\Vert _{0}\le C_t\big (|\log h|^{3/2}\ h^{t(1-\beta )}+ h^{t-1}\big )\Vert \varvec{\phi }\Vert _{t}. \end{aligned}$$

(5.40)

Further, for all $\varepsilon >0$ and $t\ge 2$, there exists $C_{\varepsilon ,t}>0$ such that

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h}\phi _h^j-\mathrm{A}^{ij}_{ {\mathrm{sing}},\delta }\phi ^j\Vert _{\infty ,I_2}\le C_{\varepsilon ,t}\big ( h^{t(1-\beta )-(1+\varepsilon )\beta }+h^{t-1}\big )\Vert \varvec{\phi }\Vert _{t}. \end{aligned}$$

(5.41)

Proof

We start by considering the decomposition

$$\begin{aligned} \mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h}\phi _h^j-\mathrm{A}^{ij}_{{\mathrm{sing}},\delta }\phi ^j&= (\mathrm{A}^{ij}_{{\mathrm{sing} },\delta ,h} -\mathrm{A}^{ij}_{{\mathrm{sing}},\delta })(\phi _h^j-{\mathrm{P}}_{N_j}\phi ^j)\nonumber \\&+\, (\mathrm{A}^{ij}_{{\mathrm{sing} },\delta ,h}-\mathrm{A}^{ij}_{{\mathrm{sing}},\delta }) \mathrm{P}_{N_j}\phi ^j+\mathrm{A}^{ij}_{{\mathrm{sing}},\delta }(\phi _h^j-\phi ^j).\qquad \qquad \end{aligned}$$

(5.42)

In order to establish Eqs. (5.40) and (5.41) we first estimate the $L^\infty (I_2)$ norm (and, thus, the $H^0$ norm) for the first two terms of the right hand side of (5.42). For the first term, from Proposition 7 and (5.3) we have

$$\begin{aligned} \Vert ( \mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h}\!-\!\mathrm{A}^{ij}_{{\mathrm{sing}},\delta })(\phi _h^j\!-\!\mathrm{P}_{N_j}\phi ^j)\Vert _{\infty ,I_2}\!\!&\le \!\!C|\log h|^{3/2}\Vert \phi _h^j-\mathrm{P}_{N_j}\phi ^j\Vert _0\nonumber \\&\le \!\! C|\log h|^{3/2}\big (\Vert \phi _h^j\!-\!\phi ^j\Vert _0\!+\! \Vert \mathrm{P}_{N_j}\phi ^j\!-\!\phi ^j\Vert _0\big )\nonumber \\&\le \!\! C_t|\log h|^{3/2}\big ( h^{t(1-\beta )} \!+\!h^{t}\big )\Vert \varvec{\phi }\Vert _t. \end{aligned}$$

(5.43)

For the second term, on the other hand, Proposition 8 and the fact that $\mathrm{G}_j\phi ^j=0$ show that for $t\ge 2$,

$$\begin{aligned} \Vert ( \mathrm{A}^{ij}_{{\mathrm{sing}},\delta ,h}-\mathrm{A}^{ij}_{{\mathrm{sing}},\delta }) \mathrm{P}_{N_j}\phi ^j \Vert _{\infty ,I_2}&\le C_t h^{t(1-\beta )+\beta }\Vert \mathrm{P}_{N_j}\phi ^j\Vert _{t}+C h^\beta \Vert \mathrm{G}_j\mathrm{P}_{N_j}\phi ^j\Vert _{\infty ,I_2} \nonumber \\&\le C_t h^{t(1-\beta )+\beta }\Vert \phi ^j\Vert _t+ C h^\beta \Vert \mathrm{P}_{N_j}\phi ^j-\phi ^j\Vert _{\infty ,I_2} \nonumber \\&\le C_t h^{t(1-\beta )+\beta }\Vert \phi ^j\Vert _t + C^{\prime }_t h^{t-1} \Vert \phi ^j\Vert _t, \end{aligned}$$

(5.44)

where in the last step we have applied the Sobolev embedding theorem $H^{\beta +1}\subset L^\infty (I_2)$ and (5.3). It thus remains to estimate the $L^\infty (I_2)$ and $H^0$ norms of the last term on the right hand side of Eq. (5.42).

We establish the $H^0$ bound first: in view of Proposition 2 and Theorem 5 we have

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }(\phi _h^j-\phi ^j)\Vert _0\le C \Vert \phi _h^j-\phi ^j\Vert _{0}\le C_t h^{t(1-\beta )}\Vert \varvec{\phi }\Vert _{t}; \end{aligned}$$

(5.45)

substituting this bound together with (5.43) and (5.44) in (5.42) yields (5.40). To obtain the $L^\infty (I_2)$ bound we proceed similarly, but using the bound

$$\begin{aligned} \Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }(\phi _h^j-\phi ^j)\Vert _{\infty , I_2}&\le C_r \Vert \mathrm{A}^{ij}_{\mathrm{sing},\delta }(\phi _h^j-\phi ^j)\Vert _{1+r}\le C_{r}^{\prime } h^{-(1+r)\beta } \Vert \phi _h^j-\phi ^j\Vert _{r}\nonumber \\&\le C_{r,t} h^{t(1-\beta )-(1+\varepsilon )\beta }\Vert \varvec{\phi }\Vert _{t}, \end{aligned}$$

(5.46)

instead of (5.45). The first inequality in Eq. (5.46) follows from the Sobolev embedding theorem. In view of the relation $\delta \approx h^\beta$, in turn, the second inequality results from (3.28) with $r:=\beta \,\varepsilon /(1+\beta )>0$. The third inequality, finally, follows from Theorem 5 and the fact that $r+(1+r)\beta =(1+\varepsilon )\beta$. The proof is now complete. $\square$