Analytical Polarizable Continuum Model for Wavelets on NURBS Patches

Maharavo Randrianarivony^*

Pappelweg 7, Zimmer 21, Sankt Augustin, Germany.

DOI: 10.4236/am.2017.88081                

Abstract

This article concerns the application of wavelet techniques on molecular surfaces constituted of four-sided patches. The Polarizable Continuum Model, which is governed by the Poisson-Boltzmann equation, is treated by means of boundary integral equations. The media inside and outside the molecular surface consist respectively of the solute and the solvent. For a given electrically charged molecule, the principal unknown is the electrostatic solvation energy when the permittivity is specified. The wavelet basis functions are constructed on the unit square which are subsequently mapped onto the patches that are assumed to be isotropically shaped and to admit similar surface areas. The initial transmission problem is recast as an integral equation in term of both the single and the double layers. Domain decomposition preconditioner serves as acceleration of the linear solver of the single layer which is badly conditioned.

Keywords

Polarizable Continuum Model, Wavelet, Poisson-Boltzmann, Patch, Electrostatic Solvation, Energy

1. Introduction

Potential applications of solute-solvent interactions include synthetic medical design, charge transfer, simulation of membranes and nanotubes as well as studies of biological process. The molecular surface $\Gamma$ is surrounded by the solvent media consisting of mobile ions admitting a specified permittivity. The domain enclosed within the molecular surface consists of the solute composed of fixed atoms admitting their respective charges [1] [2] . We focus on the Boundary Element Method (BEM) [3] [4] [5] [6] using wavelets on NURBS (Non-Uniform Rational B-Splines) patches [7] [8] [9] [10] to solve the Polarizable Continuum

Model (PCM) which treats a particular case of the Poisson-Boltzmann equation when the coefficient related to the modified Debye-Hückel parameter vanishes. Approaches based on BEM [1] [2] [11] become a preferred tool for solving linearized solvent problems because it is inconvenient to apply 3D-solvers such as FEM (Finite Element Method) [12] [13] [14] for several reasons. For one, the related transmission problem is solved by FEM in the entire 3D-space [15] [16] [17] while only the solution on the molecular surface is required. In addition, most FEM programs and analysis assume that the right hand side (RHS) is a sufficiently smooth function whereas one considers in the solvent problem a RHS which is a sum of nonsmooth Diracs defined only in the sense of distribution [18] leading to very fine adaptive mesh refinements. In the opposite, FEM appears to be the only recourse to treat nonlinear Poisson-Boltzmann [19] [20] . The principal BEM unknown is the apparent surface charge that is not a physical entity but a fictive value which behaves very much like a density of charge distributed on the molecular surface. In the application, the entire solution to the original transmission problem of the PCM is not needed as only the apparent surface charge suffices to deduce the electrostatic reaction potential. The initial transmission partial differential equations are recast as a PCM integral equation involving the single layer and the double layer operators. That will be decoupled into a pair of equations which are solved separately. First, one solves a second kind integral equation governed by the double layer operator only. The solution to that first problem serves subsequently as a RHS to another equation involving the single layer operator. The linear systems from BEM [21] [22] [23] [24] are fully populated in the case that the classical polynomial basis is utilized. Wavelet basis [25] [26] serves well as an efficient technique to obtain quasi-sparse matrices [27] . The construction of the wavelet basis [21] [28] [29] starts from univariate basis on the unit interval. Taking the tensor product in the unit square and mapping onto the patches result in the wavelet basis on the whole manifold. We assume the absence of anisotropic patches where some curved edges are relatively much longer in comparison to other ones. In addition, the surface areas of the NURBS patches ought not be considerably dissimilar. Accomplishing such a geometric property is not straightforward if the only available information consists of the nuclei coordinates and the Van-der-Waals radii of the atoms. Formerly, we have made efforts to achieve a geometric structure resulting in patches admitting comparable surface areas as well as isotropically shaped structure. The inputs are usually acquired from the format PDB or PQR for which some conversion techniques exist [30] . The system for the double layer operator produces a bounded condition number and no preconditioner is required. In the opposite, the single layer potential is severely badly conditioned and a domain decomposition technique [31] is used as a preconditioner for the reduction of the condition number which substantially decreases the required number of iterations to solve the linear system iteratively. Before proceeding further, we briefly survey our previous papers as well as computer implementations. A splitting method for CAD surfaces has been

proposed in [32] [33] for BEM simulation. Additionally, methods for checking the regularity of the mappings have been proved in [34] . The surface structure, which is required by the wavelet-BEM, is unfortunately very complicated to construct in contrast to the standard mesh generation [35] . While approximations are required to obtain global continuity in [34] for CAD objects, it can be achieved exactly for molecular surfaces in [36] . Furthermore, some quantum simulation by using wavelet BEM on single processor is described in [37] . A wavelet BEM simulation using domain decomposition techniques was described in [38] which treats the case of ASM (Additive Schwarz Method). By using Sobolev norms with negative orders [39] , it was utilized as an efficient preconditioner for the wavelet single layer potential which is badly conditioned. Recently, we gained experience [13] in elasticity nanosimulation where one has used nanotube fibers immersed in polymer matrices as quantum models. Separate studies dedicated exclusively for the wavelet double layer potential are documented in [31] which describes as well the modeling techniques to achieves the molecular NURBS patching. In this current paper, we mainly compare computational results with analytical solutions. We display numerical outcomes pertaining to the single and double layer operators separately for complete molecules. The only well-known molecular model that has explicit analytical formula for the solvation energy appears to be the single atom, which is also known as Born ion. More precisely, it consists of an atom admitting specified charge and radius which is surrounded by the solvent media with given permittivity. More details pertinent to chemical solvation energy for more complex molecules will be deferred to a subsequent paper.

2. Transmission Problem on NURBS Patches

We present in this section the required geometric structure needed by the wavelet formulation of the PCM. The problem setting of the transmission problem which is related to the PCM will also be developed. The pertaining integral equation formulation is introduced because we need only the values on the interface $\Gamma$ . We consider $M$ atoms ${\left\{x_k\right\}}_{k=1}^M$ which are located in the solute region ${\Omega }^{\text{solute}}$ and we suppose the geometry $\Gamma$ satisfies the following conditions:

・ We have a covering of the surface by four-sided patches $\Gamma ={\displaystyle {\cup }_{p=1}^P}{\Gamma }_p$ ,

・ The intersection of two different patches ${\Gamma }_p$ and ${\Gamma }_q$ is supposed to be either empty, a common curvilinear edge or a common vertex,

・ Each patch ${\Gamma }_p$ where $p=1,2,\cdots ,P$ is the image by ${\gamma }_p\mathrm{<mo>\; :\; </mo>}\square \mathrm{<mo>\; :\; </mo>}={\left[\mathrm{0,1}\right]}^2\to {\Gamma }_p$ which is described by a bivariate function that is bijective, sufficiently smooth and admitting bounded Jacobians,

・ The patch decomposition has a global continuity: for each pair of patches ${\Gamma }_p$ , ${\Gamma }_q$ sharing a curvilinear edge, the parametric representation is subject to a matching condition. That is, a bijective affine mapping $\Xi \mathrm{<mo>\; :\; </mo>}\square \to \square$ exists such that for all $x={\gamma }_p\left(s\right)$ on the common curvilinear edge, one has ${\gamma }_p\left(s\right)=\left({\gamma }_q\circ \Xi \right)\left(s\right)$ . In other words, the images of the functions ${\gamma }_p$ and ${\gamma }_q$ agree pointwise at common edges after some reorientation,

・ The manifold $\Gamma$ is orientable and the normal vector $n\left(x\right)$ is consistently pointing outward for any $x\in \Gamma$ ,

・ The patches admit similar surface areas.

An illustration of the above surface structure is depicted in Figure 1. The CAD representation of the former mappings ${\gamma }_p$ uses the concept of B-spline and NURBS [7] [9] [32] . Consider two integers $n\mathrm{,}k$ such that $n\ge k\ge 1$ . The interval [0, 1] is subdivided by a knot sequence $\tau ={\left({\tau }_i\right)}_{i=0}^{n+k}$ such that ${\tau }_i<{\tau }_{i+1}$ for $i=k-1,\cdots ,n-1$ and such that the initial and the final entries of the knot sequence are clamped ${\tau }_0=\cdots ={\tau }_{k-1}=0$ and ${\tau }_n=\cdots ={\tau }_{n+k}=1$ . One defines the B-splines [7] [8] [10] basis functions as

$N_i^{\tau ,k}\left(t\right):=\left({\tau }_{i+k}-{\tau }_i\right)\left[{\tau }_i,\cdots ,{\tau }_{i+k}\right]{\left(\cdot -t\right)}_{+}^{k-1}\text{for}i=0,\cdots ,n\text{and}t\in \left[0,1\right]$

where we employ the divided difference $\left[{\tau }_i\mathrm{,}{\tau }_{i+1}\mathrm{,}\cdots \mathrm{,}{\tau }_p\right]f$ in which we use the truncated power functions ${\left(\cdot -t\right)}_{+}^k$ given by ${\left(x-t\right)}_{+}^k:={\left(x-t\right)}^k$ if $x\ge t$ , while it is zero otherwise. The integer k controls the polynomial degree $k-1$ of the B-spline which admits an overall smoothness of ${\mathcal{C}}^{k-2}$ while the integer n controls the number of B-spline functions for which each B-spline basis $N_i^{\tau \mathrm{,}k}$ is supported by $\left[{\tau }_i\mathrm{,}{\tau }_{i+k}\right]$ . The NURBS patch ${\gamma }_p$ admitting the control points $d_{i\mathrm{,}j}\in {ℝ}^3$ and weights $w_{i\mathrm{,}j}\in {ℝ}^{+}$ is expressed as

${\gamma }_p\left(u,v\right)=\frac{{{\displaystyle \sum }}_{i=0}^n{{\displaystyle \sum }}_{j=0}^n{w}_{i,j}{d}_{i,j}{N}_i^{\tau ,k}\left(u\right)N_j^{\tau ,k}\left(v\right)}{{{\displaystyle \sum }}_{i=0}^n{{\displaystyle \sum }}_{j=0}^m{w}_{i,j}{N}_i^{\tau ,k}\left(u\right)N_j^{\tau ,k}\left(v\right)}\in {ℝ}^3\forall \left(u,v\right)\in \square$ (1)

which describes thus a parametric function from the unit square $\square =\left[0,1\right]\times \left[0,1\right]$ onto the four-sided patch ${\Gamma }_p\subset {ℝ}^3$ . The similarity of the surface areas are needed in practice for the wavelet threshold to function in a unified manner. We

Figure 1. Patch representation of a Water Cluster with 1089 NURBS.

will the wavelet threshold to function in a unified manner. We will consider only geometries which are globally smooth and which admit moderate curvature. For each patch ${\Gamma }_p$ , the Gram determinant is denoted by

$G_p\left(t\right)=G_p\left(t_1\mathrm{,}{t}_2\right)\mathrm{<mo>\; :\; </mo>}=\Vert \frac{\partial {\gamma }_p}{\partial t_1}\left(t\right)\times \frac{\partial {\gamma }_p}{\partial t_2}\left(t\right)\Vert \forall t=\left(t_1\mathrm{,}{t}_2\right)\in \square$ (2)

which quantifies the norm of the normal vector at the image ${\gamma }_p\left(t\right)\in {\Gamma }_p$ of any point $t$ belonging to the unit square $\square =\left[0,1\right]\times \left[0,1\right]$ by the parametric function ${\gamma }_p$ . The Gram determinant $G_p$ and its partial derivatives are assumed to be bounded

$0<c\le \underset{p=1,\cdots ,M}{min}\underset{t\in \square }{inf}{G}_p\left(t\right)<\underset{p=1,\cdots ,M}{max}\underset{t\in \square }{sup}{G}_p\left(t\right)\le C<\infty \mathrm{,}$ (3)

$\underset{t\in \square }{sup}\left|{\partial }_{\alpha }{G}_p\left(t\right)\right|=\underset{t\in \square }{sup}\left|\frac{{\partial }^{\left|\alpha \right|}}{{\partial }^{{\alpha }_1}{t}_1{\partial }^{{\alpha }_2}{t}_2}{G}_p\left(t\right)\right|\le C<\infty \mathrm{,}t=\left(t_1\mathrm{,}{t}_2\right)$ (4)

for $\alpha =\left({\alpha }_1,{\alpha }_2\right)$ where $\left|\alpha \right|={\alpha }_1+{\alpha }_2\le \eta$ for $\eta$ sufficiently large. We consider the transmission problem for a parameter $\epsilon \ge 1$ related to the dielectric coefficient:

$-{\nabla }^{2}u\left(x\right)={\displaystyle \underset{k=1}{\overset{M}{\sum }}}{q}_k\delta \left(x-x_k\right)=\mathrm{<mo>\; :\; </mo>}\rho \left(x\right)\text{for}x\in {\Omega }^{\text{solute}}\mathrm{,}$ (5)

$-\epsilon {\nabla }^{2}u\left(x\right)=0\text{for}x\in {\Omega }^{\text{solvent}}\mathrm{,}$ (6)

$\underset{\begin{array}{l}x\to x_0\\ x\in {\Omega }^{\text{solute}}\end{array}}{lim}u\left(x\right)=\underset{\begin{array}{l}x\to x_0\\ x\in {\Omega }^{\text{solute}}\end{array}}{lim}u\left(x\right)\text{for}{x}_0\in \Gamma \mathrm{,}$ (7)

$\underset{\begin{array}{l}x\to x_0\\ x\in {\Omega }^{\text{solute}}\end{array}}{lim}\langle \nabla u\left(x\right)\mathrm{,}n\left(x_0\right)\rangle =\epsilon \underset{\begin{array}{l}x\to x_0\\ x\in {\Omega }^{\text{solute}}\end{array}}{lim}\langle \nabla u\left(x\right)\mathrm{,}n\left(x_0\right)\rangle \text{for}{x}_0\in \Gamma \mathrm{,}$ (8)

$u\left(x\right)\to 0\text{as}\Vert x\Vert \to \infty$ (9)

where $q_k$ denotes the electric charge of the k-th atom while $\delta \left(\cdot -x_k\right)$ is the Dirac distribution centered at the coordinates of the nucleus $x_k$ . It consists of a special case of the Poisson-Boltzmann equation in the situation that the effect of the Debye-Hückel parameter is neglected. We are not solving those equations directly, rather we solve only the pertaining integral equations which are located on the interface surface $\Gamma$ . Introduce the single layer $V$ , the double layer $K$ and the adjoint double layer $K^{\mathrm{<mo>\; *\; </mo>}}$ operators for $x\in \Gamma$ as follows

$\left(Vv\right)\left(x\right)\mathrm{<mo>\; :\; </mo>}=\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{1}{\Vert x-y\Vert }v\left(y\right){\text{d}}_{\Gamma }y\mathrm{,}$

$\left(Kv\right)\left(x\right)\mathrm{<mo>\; :\; </mo>}=\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{\partial }{\partial n\left(y\right)}\left[\frac{1}{\Vert x-y\Vert }\right]v\left(y\right){\text{d}}_{\Gamma }y=\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{\langle n\left(y\right)\mathrm{,}x-y\rangle }{{\Vert x-y\Vert }^3}v\left(y\right){\text{d}}_{\Gamma }y\mathrm{,}$

$\left(K^{\mathrm{<mo>\; *\; </mo>}}v\right)\left(x\right)\mathrm{<mo>\; :\; </mo>}=\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{\partial }{\partial n\left(x\right)}\left[\frac{1}{\Vert x-y\Vert }\right]v\left(y\right){\text{d}}_{\Gamma }y=\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{\langle n\left(x\right)\mathrm{,}x-y\rangle }{{\Vert x-y\Vert }^3}v\left(y\right){\text{d}}_{\Gamma }y\mathrm{.}$

Consider the components $\left(u_i\mathrm{,}{u}_e\right)$ of the solution u inside the solute ${\Omega }^{\text{solute}}$ and the solvent ${\Omega }^{\text{solvent}}$ respectively. We recall very briefly the transformation of the transmission problem (5)-(9) into integral equations which follow the procedure of [1] [2] where one replaces $\left(u_i\mathrm{,}{u}_e\right)$ by $\left(w_i\mathrm{,}{w}_e\right)$ such that

#Math_82# (10)

The apparent surface charge is defined as $\sigma \mathrm{<mo>\; :\; </mo>}=V^{-1}{w}_i$ so that one obtains

$u_i\left(x\right)=\frac{1}{4\text{π}}{\displaystyle {\int }_{\Gamma }}\frac{\sigma \left(x\right)}{\Vert x-y\Vert }{\text{d}}_{\Gamma }y+\frac{1}{4\text{π}}{\displaystyle \underset{k=1}{\overset{M}{\sum }}}\frac{q_k}{\Vert x-x_k\Vert }$ (11)

which shows that $\sigma$ behaves as a charge distribution over $\Gamma$ while the second term is a sum over the actual charges $\left\{q_k\right\}$ . The representation formula [4] yields

$V\left(\partial w_i\right)-K{w}_i=\frac{1}{2}{w}_i\mathrm{,}\frac{1}{\epsilon }V\left(\partial w_e\right)-K{w}_e=-\frac{1}{2}{w}_i$ (12)

where $\partial$ represents the conormal derivative. By combining them with the transmission jump conditions $\mathrm{<mo\; stretchy="false">\; [\; </mo>}\text{}\mathrm{<mo\; stretchy="false">\; [\; </mo>}w\mathrm{<mo\; stretchy="false">\; ]\; </mo>}\text{}\mathrm{<mo\; stretchy="false">\; ]\; </mo>\; <mo>\; :\; </mo>}=w_e-w_i=f_e-f_i$ and $\mathrm{<mo\; stretchy="false">\; [\; </mo>}\text{}\mathrm{<mo\; stretchy="false">\; [\; </mo>}\partial w\mathrm{<mo\; stretchy="false">\; ]\; </mo>}\text{}\mathrm{<mo\; stretchy="false">\; ]\; </mo>\; <mo>\; :\; </mo>}=\partial w_e-\partial w_i=\partial f_e-\partial f_i$ , one obtains

$V\left(\partial w_i-\partial f_e+\partial f_i\right)-\epsilon K\left(w_i-f_e+f_i\right)=-\frac{\epsilon }{2}\left(w_i-f_e+f_i\right)$ (13)

of which the full detail is described in [1] [2] . The fact that $KV=V{K}^{\mathrm{<mo>\; *\; </mo>}}$ , $w_i=V\sigma$ and the Newton potential ${\mathcal{N}}_{\rho }\mathrm{<mo>\; :\; </mo>}=f_i$ lead to

$\left[\left(\frac{1}{2}I-K\right)V+\frac{1}{\epsilon }V\left(\frac{1}{2}I+K^{\mathrm{<mo>\; *\; </mo>}}\right)\right]\sigma =-\left(\frac{1}{2}I-K\right){\mathcal{N}}_{\rho }-\frac{1}{\epsilon }V\left(\partial {\mathcal{N}}_{\rho }\right)\mathrm{.}$ (14)

Apply the representation formula to the Newton potential which is harmonic to obtain $\left(\frac{1}{2}I-K\right){\mathcal{N}}_{\rho }+V\left(\partial {\mathcal{N}}_{\rho }\right)=0$ which results in

$\left[\left(\frac{1}{2}I-K\right)+\frac{1}{\epsilon }\left(\frac{1}{2}I+K\right)\right]V\sigma =-\left(\frac{1}{2}I-K\right){\mathcal{N}}_{\rho }+\frac{1}{\epsilon }\left(\frac{1}{2}I-K\right){\mathcal{N}}_{\rho }$ (15)

$\left[\frac{1}{2}\left(\epsilon +1\right)I+\left(1-\epsilon \right)K\right]V\sigma ={\mathcal{N}}_{\rho }-\left[\frac{1}{2}\left(\epsilon +1\right)I+\left(1-\epsilon \right)K\right]{\mathcal{N}}_{\rho }$ (16)

$V\sigma ={\left[\frac{1}{2}\left(\epsilon +1\right)I+\left(1-\epsilon \right)K\right]}^{-1}{\mathcal{N}}_{\rho }-{\mathcal{N}}_{\rho }\mathrm{.}$ (17)

In practice, we are not interested in the entire solution $u$ of the transmission problem but only in the electrostatic reaction potential

$P^{REAC}\mathrm{<mo>\; :\; </mo>}=\frac{1}{2}{\displaystyle \underset{k=1}{\overset{M}{\sum }}}\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{\sigma \left(y\right)}{\Vert x_k-y\Vert }{\text{d}}_{\Gamma }y\mathrm{.}$ (18)

Thus, our objective is to seek for the apparent surface charge $\sigma$ satisfying

$V\sigma =\lambda {\left(I-\mu K\right)}^{-1}m+g$ (19)

where $\lambda$ and $\mu$ are constant parameters while m and g are given functions which are sufficiently smooth. For the practical applications related to (5)-(9), the coefficients $\lambda$ and $\mu$ and the functions m and g are as follows

$\lambda =2\epsilon +\mathrm{1,}\mu =2\left(\epsilon -1\right)\epsilon +\mathrm{1,}m\left(x\right)=-g\left(x\right)={\displaystyle \underset{k=1}{\overset{M}{\sum }}}\frac{q_k}{\Vert x-x_k\Vert }\text{for}x\in \Gamma \mathrm{.}$ (20)

In order to simplify the presentation, we assume henceforth the parameters $\lambda$ and $\mu$ are unity. Everything remains valid with little modifications for other constant parameters. The above problem (19) is decoupled as two integral equations:

$V\sigma =f+g$ (21)

where

$f\mathrm{<mo>\; :\; </mo>}={\left(I-K\right)}^{-1}m\text{or}\text{equivalently}\text{solve}\left(I-K\right)f=m\mathrm{.}$ (22)

That is, the solution to (22) is used as a RHS in (21). To numerically treat those last problems by means of the wavelet technique, several approximations are involved:

・ The Galerkin error related to the single layer $V$ by using a finite dimensional space ${\mathcal{V}}_L\left(\Gamma \right)$ ,

・ Similar Galerkin error but for the double layer potential $K$ ,

・ Replacing the RHS for $V\sigma =f$ by the approximate solution from the double layer equation,

・ Replacing the operator $V$ by its truncated version $\tilde{V}$ ,

・ Similar truncation but for substituting $K$ by $\tilde{K}$ .

We will consider the space of piecewise constants ${\mathcal{V}}_L\left(\Gamma \right)$ whose construction will be specified later on. Concerning the discrete Galerkin variational formulation in ${\mathcal{V}}_L\left(\Gamma \right)$ for discretizing the integral Equations ((21) and (22)), one searches for ${\sigma }_L\mathrm{,}{f}_L\in {\mathcal{V}}_L\left(\Gamma \right)$ such that for all $v_L\in {\mathcal{V}}_L\left(\Gamma \right)$

${\displaystyle {\int }_{\Gamma }}{\displaystyle {\int }_{\Gamma }}{\mathcal{K}}^{SL}\left(x\mathrm{,}y\right){\sigma }_L\left(x\right)v_L\left(y\right){\text{d}}_{\Gamma }x{\text{d}}_{\Gamma }y={\displaystyle {\int }_{\Gamma }}{f}_L\left(x\right)v_L\left(x\right){\text{d}}_{\Gamma }x+{\displaystyle {\int }_{\Gamma }}g\left(x\right)v_L\left(x\right){\text{d}}_{\Gamma }x\mathrm{,}$ (23)

${\displaystyle {\int }_{\Gamma }}{f}_L\left(x\right)v_L\left(x\right)d_{\Gamma }x-{\displaystyle {\int }_{\Gamma }}{\displaystyle {\int }_{\Gamma }}{\mathcal{K}}^{DL}\left(x\mathrm{,}y\right)f_L\left(x\right)v_L\left(y\right){\text{d}}_{\Gamma }x{\text{d}}_{\Gamma }y={\displaystyle {\int }_{\Gamma }}m\left(x\right)v_L\left(x\right){\text{d}}_{\Gamma }x$ (24)

where ${\mathcal{K}}^{SL}\left(x\mathrm{,}y\right)$ and ${\mathcal{K}}^{DL}\left(x\mathrm{,}y\right)$ are respectively the kernel functions for the single layer $V$ and double layer $K$ . Most of the following theoretical works are inspired from different sources including [1] [2] [21] [22] [27] [40] [41] [42] [43] .

3. A-Priori Estimate for Multi-Wavelet PCM

We consider in this section the multi-wavelet Galerkin formulation of the PCM problem. We recall several important properties which are used in the deduction of the a-priori error estimate when no wavelet. For the unidimensional basis functions, we introduce the next knot sequence on level $\mathcal{l}=0,1,\cdots ,L^{1D}$ ,

${\zeta }^{\mathcal{l}}=\left\{{\zeta }_0^{\mathcal{l}}\mathrm{,}{\zeta }_1^{\mathcal{l}}\mathrm{,}\cdots \mathrm{,}{\zeta }_{2^{\mathcal{l}}}^{\mathcal{l}}\right\}\subset \left[\mathrm{0,1}\right]\text{where}{\zeta }_i^{\mathcal{l}}=i{2}^{-\mathcal{l}}\mathrm{.}$ (25)

The internal knots on the next level $\left(\mathcal{l}+1\right)$ are obtained by inserting a new knot in the middle of two consecutive knots on the lower level $\mathcal{l}$ . Introduce the piecewise constant linear space in the unit interval $\left[\mathrm{0,1}\right]$ on level $\mathcal{l}$ :

${\mathbb{V}}_{\mathcal{l}}\left[\mathrm{0,1}\right]\mathrm{<mo>\; :\; </mo>}=\text{span}\left\{{\varphi }_i^{\mathcal{l}}={\chi }_{\left[{\zeta }_{i-1}^{\mathcal{l}}\mathrm{,}{\zeta }_i^{\mathcal{l}}\right]}\mathrm{,}i=\mathrm{1,}\cdots {\mathrm{,2}}^{\mathcal{l}}\right\}$ (26)

where ${\chi }_D$ designates the characteristic function with respect to $D$ . By using the two scale relation

${\varphi }_1^0\left(t\right)={\varphi }_1^0\left(2t\right)+{\varphi }_1^0\left(2t-1\right)\text{for}\text{all}t\in \left[0,1\right]$ (27)

and the inclusion ${\zeta }^{\mathcal{l}}\subset {\zeta }^{\mathcal{l}+1}$ , the spaces ${\mathbb{V}}_{\mathcal{l}}\left[\mathrm{0,1}\right]$ form a nested sequence of subspaces:

${\mathbb{V}}_0\left[\mathrm{0,1}\right]\subset {\mathbb{V}}_1\left[\mathrm{0,1}\right]\subset \cdots \subset {\mathbb{V}}_{L^{1D}}\left[\mathrm{0,1}\right]\subset {\mathbb{V}}_2\left[\mathrm{0,1}\right]\mathrm{.}$ (28)

On account of the nestedness (28), the space ${\mathbb{V}}_{\mathcal{l}}\left[\mathrm{0,1}\right]$ is expressed as an orthogonal sum

${\mathbb{V}}_{\mathcal{l}}\left[\mathrm{0,1}\right]={\mathbb{V}}_{\mathcal{l}-1}\left[\mathrm{0,1}\right]\oplus {\mathbb{W}}_{\mathcal{l}}\left[\mathrm{0,1}\right]$ (29)

with respect to the ${\mathbb{L}}^2$ -scalar product where ${\mathbb{W}}_{\mathcal{l}}\left[\mathrm{0,1}\right]$ is the complementary wavelet space

${\mathbb{W}}_{\mathcal{l}}\left[\mathrm{0,1}\right]=\text{span}\left\{{\psi }_i^{\mathcal{l}}\in {\mathbb{V}}_{\mathcal{l}}\left[\mathrm{0,1}\right]\mathrm{,}{\langle {\psi }_i^{\mathcal{l}}\mathrm{,}\varphi \rangle }_{{\mathbb{L}}^2\left[\mathrm{0,1}\right]}=\mathrm{0,}\forall \varphi \in {\mathbb{V}}_{\mathcal{l}-1}\left[\mathrm{0,1}\right]\right\}\mathrm{.}$ (30)

For the explicit expression of the wavelet functions ${\psi }_i^{\mathcal{l}}$ , we use the Haar wavelet defined on $\left[\mathrm{0,1}\right]$ by

$\{\begin{array}{l}{\psi }^{Haar}\left(t\right)\mathrm{<mo>\; :\; </mo>}=+1\text{for}t\in \left[\mathrm{0,}1/2\right)\\ {\psi }^{Haar}\left(t\right)\mathrm{<mo>\; :\; </mo>}=-1\text{for}t\in \left[1/2\mathrm{,1}\right]\end{array}\text{such}\text{that}{\displaystyle {\int }_0^1}{\psi }^{Haar}\left(t\right)\text{d}t=0$ (31)

whose relation with the single scale basis is ${\psi }^{Haar}={\chi }_{\left[\mathrm{0,}1/2\right)}-{\chi }_{\left[1/2\mathrm{,1}\right]}$ . By using dilation and shift, one defines for $\mathcal{l}=1,\cdots ,L^{1D}$ and $i=1,\cdots ,2^{\mathcal{l}-1}$

#Math_156# (32)

The wavelet functions ${\psi }_i^{\mathcal{l}}$ constitute an orthonormal basis

${\displaystyle {\int }_0^1}{\psi }_{i_1}^{{\mathcal{l}}_1}\left(t\right){\psi }_{i_2}^{{\mathcal{l}}_2}\left(t\right)\text{d}t={\delta }_{{\mathcal{l}}_1,{\mathcal{l}}_2}{\delta }_{i_1,i_2}$ (33)

where the first Dirac ${\delta }_{{\mathcal{l}}_1\mathrm{,}{\mathcal{l}}_2}$ comes from the inter-level orthogonality while the second Dirac ${\delta }_{i_1\mathrm{,}{i}_2}$ is justified by the non-overlapping of $\text{Support}\left({\psi }_{i_1}^{\mathcal{l}}\right)$ and $\text{Support}\left({\psi }_{i_2}^{\mathcal{l}}\right)$ on the same level $\mathcal{l}$ . By applying the decomposition (29) recursively, one obtains on the maximal level $L^{1D}$

${\mathbb{V}}_{L^{1D}}\left[\mathrm{0,1}\right]={\mathbb{V}}_0\left[\mathrm{0,1}\right]\oplus \left(\underset{\mathcal{l}=1}{\overset{L^{1D}}{{\displaystyle \oplus }}}{\mathbb{W}}_{\mathcal{l}}\left[\mathrm{0,1}\right]\right)=\underset{\mathcal{l}=0}{\overset{L^{1D}}{{\displaystyle \oplus }}}{\mathbb{W}}_{\mathcal{l}}\left[\mathrm{0,1}\right]$ (34)

where

${\mathbb{W}}_0\left[\mathrm{0,1}\right]\mathrm{<mo>\; :\; </mo>}={\mathbb{V}}_0\left[\mathrm{0,1}\right]\text{and}{\psi }_1^0\mathrm{<mo>\; :\; </mo>}={\varphi }_1^0$ (35)

so that we have the dimensionalities

${\omega }_{\mathcal{l}}\mathrm{<mo>\; :\; </mo>}=\text{dim}\left({\mathbb{W}}_{\mathcal{l}}\left[\mathrm{0,1}\right]\right)\text{where}{\omega }_0=1\text{and}{\omega }_{\mathcal{l}}=2^{\mathcal{l}-1}\forall \mathcal{l}=\mathrm{1,}\cdots \mathrm{,}{L}^{1D}\mathrm{.}$ (36)

Thus, a function $u\in {\mathbb{V}}_{L^{1D}}\left[\mathrm{0,1}\right]$ has two representations: in the single-scale basis and in the multiscale basis, we have respectively

$u\left(t\right)={\displaystyle \underset{i=1}{\overset{2^{L^{1D}}}{\sum }}}{u}_i^{S.Sc.}{\varphi }_i^L\left(t\right)\text{where}{u}_i^{S.Sc.}\in ℝ\mathrm{,}$ (37)

$u\left(t\right)={\displaystyle \underset{\mathcal{l}=0}{\overset{L^{1D}}{\sum }}}{\displaystyle \underset{k=1}{\overset{{\omega }_{\mathcal{l}}}{\sum }}}{u}_{\mathcal{l}\mathrm{,}k}^{M.Sc.}{\psi }_k^{\mathcal{l}}\left(t\right)\text{where}{u}_{\mathcal{l}\mathrm{,}k}^{M.Sc.}\in ℝ\mathrm{.}$ (38)

From here onward, we use the notation $X\precsim Y$ if there is a constant c such that $X\le cY$ in which $c$ is independent on the level $\mathcal{l}$ and the maximal level. Besides, the notation $X\simeq Y$ amounts to $Y\precsim X\precsim Y$ . Moreover, we denote

${{\displaystyle \exp }}_2\left\{X\right\}:=2^X$ (39)

when the expression of $X$ is complicated. It has the properties:

${{\displaystyle \exp }}_2\left\{X_1+X_2\right\}={{\displaystyle \exp }}_2\left\{X_1\right\}{{\displaystyle \exp }}_2\left\{X_2\right\},{{\displaystyle \exp }}_2\left\{-X\right\}=1/{{\displaystyle \exp }}_2\left\{X\right\}$ (40)

${\left[{{\displaystyle \exp }}_2\left\{X\right\}\right]}^r={{\displaystyle \exp }}_2\left\{rX\right\},{{\displaystyle \exp }}_2\left\{0\right\}=1,{{\displaystyle \exp }}_2\left\{1\right\}=2.$ (41)

The next norm equivalences related to the coefficients are valid [21] [29]

${\Vert u\Vert }_{{\mathbb{L}}^2\left[\mathrm{0,1}\right]}\simeq {\Vert {\left\{u_i^{S.Sc.}\right\}}_i\Vert }_{{\mathcal{l}}^2}\simeq {\Vert {\left\{u_{\mathcal{l}\mathrm{,}k}^{M.Sc.}\right\}}_{\mathcal{l}\mathrm{,}k}\Vert }_{{\mathcal{l}}^2}\mathrm{.}$ (42)

Due to the property (33) and ${\displaystyle {\int }_0^1}{\varphi }_1^0\left(t\right)\text{d}t=1$ , the orthogonal projection of any $u\in {\mathbb{V}}_{L^{1D}}\left[\mathrm{0,1}\right]$ onto ${\mathbb{V}}_q\left[\mathrm{0,1}\right]$ verifies

$Q_{q}u={\displaystyle \underset{\mathcal{l}=0}{\overset{q}{\sum }}}{\displaystyle \underset{k=1}{\overset{{\omega }_{\mathcal{l}}}{\sum }}}{\langle u\mathrm{,}{\psi }_k^{\mathcal{l}}\rangle }_{{\mathbb{L}}^2\left[\mathrm{0,1}\right]}{\psi }_k^{\mathcal{l}}\mathrm{.}$ (43)

The 2D-wavelet space on the unit square , is defined for any level $\mathcal{l}=0,1,\cdots ,L:=2L^{1D}$ in term of tensor product as follows

${\mathbb{W}}_{\mathcal{l}}(\square )\mathrm{<mo>\; :\; </mo>}=\underset{{\mathcal{l}}_1+{\mathcal{l}}_2=\mathcal{l}}{{\displaystyle \oplus }}\left\{{\mathbb{W}}_{{\mathcal{l}}_1}\left[\mathrm{0,1}\right]\otimes {\mathbb{W}}_{{\mathcal{l}}_2}\left[\mathrm{0,1}\right]\right\}$ (44)

$=\text{span}\left\{{\tilde{\psi }}_{\left(\mathcal{l}\mathrm{,}i\right)}={\psi }_{i_1}^{{\mathcal{l}}_1}\otimes {\psi }_{i_2}^{{\mathcal{l}}_2}\mathrm{,}{\mathcal{l}}_1+{\mathcal{l}}_2=\mathcal{l}\right\}\mathrm{,}$ (45)

${\mathbb{V}}_{\mathcal{l}}(\square )\mathrm{<mo>\; :\; </mo>}=\underset{k=0}{\overset{\mathcal{l}}{{\displaystyle \oplus }}}{\mathbb{W}}_k(\square )={\mathbb{V}}_{\mathcal{l}-1}(\square )\oplus {\mathbb{W}}_{\mathcal{l}}(\square )\mathrm{.}$ (46)

On each patch ${\Gamma }_p$ ( $p=1,\cdots ,P$ ) and on the whole surface $\Gamma ={\displaystyle {\cup }_{p=1}^P{\Gamma }_p}$ , we define for the level $\mathcal{l}=0,1,\cdots ,L$

${\mathcal{V}}_{\mathcal{l}}\left({\Gamma }_p\right)\mathrm{<mo>\; :\; </mo>}=\text{span}\left\{{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\mathrm{<mo>\; :\; </mo>}={\tilde{\psi }}_{\left(\mathcal{l}\mathrm{,}i\right)}\circ {\gamma }_p^{-1}\mathrm{<mo>\; :\; </mo>}{\tilde{\psi }}_{\left(\mathcal{l}\mathrm{,}i\right)}\in {\mathbb{V}}_{\mathcal{l}}(\square )\right\}\mathrm{,}$ (47)

${\mathcal{W}}_{\mathcal{l}}\left({\Gamma }_p\right)\mathrm{<mo>\; :\; </mo>}=\text{span}\left\{{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\mathrm{<mo>\; :\; </mo>}={\tilde{\psi }}_{\left(\mathcal{l}\mathrm{,}i\right)}\circ {\gamma }_p^{-1}\mathrm{<mo>\; :\; </mo>}{\tilde{\psi }}_{\left(\mathcal{l}\mathrm{,}i\right)}\in {\mathbb{V}}_{\mathcal{l}}(\square )\right\}\mathrm{,}$ (48)

${\mathcal{V}}_{\mathcal{l}}\left(\Gamma \right)\mathrm{<mo>\; :\; </mo>}=\underset{p=1}{\overset{P}{{\displaystyle \oplus }}}{\mathcal{V}}_{\mathcal{l}}\left({\Gamma }_p\right)\mathrm{,}{\mathcal{W}}_{\mathcal{l}}\left(\Gamma \right)\mathrm{<mo>\; :\; </mo>}=\underset{p=1}{\overset{P}{{\displaystyle \oplus }}}{\mathcal{W}}_{\mathcal{l}}\left({\Gamma }_p\right)\mathrm{.}$ (49)

The space ${\mathcal{V}}_L\left(\Gamma \right)$ corresponds to the piecewise constant functions on a tensor product mesh admitting a step size of $h_L=\mathcal{O}\left(2^{-L^{1D}}\right)=\mathcal{O}\left(2^{-L/2}\right)$ . It is deduced from the above construction that the next nested inclusions are valid

${\mathcal{V}}_0\left(\Gamma \right)\subset \cdots \subset {\mathcal{V}}_{\mathcal{l}}\left(\Gamma \right)\subset {\mathcal{V}}_{\mathcal{l}+1}\left(\Gamma \right)\subset \cdots \subset {\mathcal{V}}_L\left(\Gamma \right)\mathrm{.}$ (50)

Some immediate properties for all wavelet functions ${\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\in {\mathcal{W}}_{\mathcal{l}}\left(\Gamma \right)$ are as follows:

$\text{meas}\left[S\left(\mathcal{l}\mathrm{,}i\right)\right]=\mathcal{O}\left(2^{-\mathcal{l}}\right)\text{where}S\left(\mathcal{l}\mathrm{,}i\right)\mathrm{<mo>\; :\; </mo>}=\text{support}\left[{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\right]\mathrm{,}$ (51)

${\Delta }_{\mathcal{l}}\mathrm{<mo>\; :\; </mo>}=\text{diam}\left[S\left(\mathcal{l}\mathrm{,}i\right)\right]=\mathcal{O}\left(2^{-\mathcal{l}/2}\right)\mathrm{,}$ (52)

${\max }_{u,v}{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\left(u\mathrm{,}v\right)=\mathcal{O}\left(2^{\mathcal{l}/2}\right)\mathrm{,}$ (53)

${\Delta }_{\mathcal{l}}=\mathcal{O}\left(2^{-\mathcal{l}/2}\right)\succsim {\displaystyle \underset{\left|\alpha \right|=1}{\sum }}\frac{1}{\alpha \mathrm{<mo>\; !\; </mo>}}\underset{u_1\mathrm{,}{u}_2\in S\left(\mathcal{l}\mathrm{,}i\right)}{max}{\left|u_1-u_2\right|}^{\alpha }\mathrm{.}$ (54)

They are easily derived from the 1D definitions where $\mathcal{l}={\mathcal{l}}_1+{\mathcal{l}}_2$ and the boundedness of the mappings ${\gamma }_p$ .

The space of square integrable functions is

${\mathbb{L}}^2\left(\Gamma \right)\mathrm{<mo>\; :\; </mo>}=\left\{f\mathrm{<mo>\; :\; </mo>}\Gamma \to ℝ\mathrm{,}{\displaystyle {\int }_{\Gamma }}{\left|f\left(x\right)\right|}^2{\text{d}}_{\Gamma }x<\infty \right\}$ (55)

which is equipped with the following scalar product and norm after transformation onto $\square$ ,

${\langle u\mathrm{,}v\rangle }_{{\mathbb{L}}^2\left(\Gamma \right)}\mathrm{<mo>\; :\; </mo>}={\displaystyle {\int }_{\Gamma }}u\left(x\right)v\left(x\right){\text{d}}_{\Gamma }x={\displaystyle \underset{p=1}{\overset{P}{\sum }}}{\displaystyle {\int }_{\square }}u\left({\gamma }_p\left(t\right)\right)v\left({\gamma }_p\left(t\right)\right)G_p\left(t\right)\text{d}t\mathrm{,}$ (56)

${\Vert v\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\mathrm{<mo>\; :\; </mo>}={\left[{\langle v\mathrm{,}v\rangle }_{{\mathbb{L}}^2\left(\Gamma \right)}\right]}^{1/2}\mathrm{.}$ (57)

The Sobolev space on $\Gamma$ for a non-negative integer $k$ is

${ℍ}^k\left(\Gamma \right)\mathrm{<mo>\; :\; </mo>}=\left\{f\in {\mathbb{L}}^2\left(\Gamma \right)\mathrm{<mo>\; :\; </mo>}{\Vert {\partial }^{\alpha }f\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}<\infty \text{for}\text{all}\left|\alpha \right|\le k\right\}$ (58)

where the differentiation ${\partial }^{\alpha }f$ is interpreted in the sense of distribution [18] such that ${\langle {\partial }^{\alpha }f\mathrm{,}g\rangle }_{{\mathbb{L}}^2\left(\Gamma \right)}={\left(-1\right)}^{\left|\alpha \right|}{\langle f\mathrm{,}{\partial }^{\alpha }g\rangle }_{{\mathbb{L}}^2\left(\Gamma \right)}$ for all compactly supported smooth functions $g$ . The Sobolev space ${ℍ}^k\left(\Gamma \right)$ is endowed with the norm

${\Vert f\Vert }_{{ℍ}^k\left(\Gamma \right)}^2\mathrm{<mo>\; :\; </mo>}={\displaystyle \underset{\left|\alpha \right|\le k}{\sum }}{\Vert {\partial }^{\alpha }f\Vert }_{{ℍ}^2\left(\Gamma \right)}^2\mathrm{.}$ (59)

Concerning a real positive order $p=k+\theta$ such that $k$ is an integer and $\theta \in \left]\mathrm{0,1}\right[$ , the Sobolev space ${ℍ}^p\left(\Gamma \right)$ consists of the functions such that their norms with respect to the next Slobodeckij norm are finite

${\Vert f\Vert }_{{ℍ}^p\left(\Gamma \right)}^2={\Vert f\Vert }_{{ℍ}^k\left(\Gamma \right)}^2+{\displaystyle \underset{\left|\alpha \right|=k}{\sum }}{\displaystyle {\int }_{\Gamma \times \Gamma }}\frac{{\left|{\partial }^{\alpha }f\left(x\right)-{\partial }^{\alpha }f\left(y\right)\right|}^2}{{\Vert x-y\Vert }^{2+2\theta }}{\text{d}}_{\Gamma }x{\text{d}}_{\Gamma }y\mathrm{.}$ (60)

For negative orders, one employs the dual spaces ${ℍ}^{-p}\left(\Gamma \right)={\left[{ℍ}^p\left(\Gamma \right)\right]}^{\mathrm{<mo>\; *\; </mo>}}$ equipped with the dual norm

${\Vert u\Vert }_{-p}={\Vert u\Vert }_{{ℍ}^{-p}\left(\Gamma \right)}=\underset{0\ne v\in {ℍ}^p\left(\Gamma \right)}{sup}\frac{{\langle u\mathrm{,}v\rangle }_{{ℍ}^2\left(\Gamma \right)}}{{\Vert v\Vert }_{{ℍ}^p\left(\Gamma \right)}}\mathrm{.}$ (61)

We will denote the orthogonal projection with respect to the ${\mathbb{L}}_2$ -scalar product onto ${\mathcal{V}}_{\mathcal{l}}\left(\Gamma \right)$ by ${\mathcal{Q}}_{\mathcal{l}}$ such that

${\langle {\mathcal{Q}}_{\mathcal{l}}v\mathrm{,}w\rangle }_{{\mathbb{L}}_2\left(\Gamma \right)}={\langle v\mathrm{,}w\rangle }_{{\mathbb{L}}_2\left(\Gamma \right)}\forall w\in {\mathcal{V}}_{\mathcal{l}}\left(\Gamma \right)\mathrm{.}$ (62)

Theorem 1. On the 2D level $L$ , consider the discrete Galerkin PCM problems for seeking ${\sigma }_L\mathrm{,}{f}_L\in {\mathcal{V}}_L\left(\Gamma \right)$ such that

$\langle V{\sigma }_L\mathrm{,}{v}_L\rangle =\langle f_L\mathrm{,}{v}_L\rangle +\langle g\mathrm{,}{v}_L\rangle \forall v_L\in {\mathcal{V}}_L\left(\Gamma \right)$ (63)

$\langle \left(I-K\right)f_L\mathrm{,}{v}_L\rangle =\langle m\mathrm{,}{v}_L\rangle \forall v_L\in {\mathcal{V}}_L\left(\Gamma \right)\mathrm{.}$ (64)

Then, the accuracy between the apparent surface charge $\sigma$ from (21) and ${\sigma }_L$ from (63) for $r\ge 1/2$ and $t\ge -1/2$ satisfies

${\Vert \sigma -{\sigma }_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim \frac{1}{2^{L\left(r+1/2\right)/2}}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+\frac{1}{2^{L\left(t-1/2\right)/2}}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (65)

In addition, the accuracy of the reaction potential $P^{REAC}$ from (18) and

$P_L^{REAC}\mathrm{<mo>\; :\; </mo>}=\frac{1}{2}{\displaystyle \underset{k=1}{\overset{M}{\sum }}}\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{{\sigma }_L\left(y\right)}{\Vert x_k-y\Vert }{\text{d}}_{\Gamma }y\text{fulfills}$ (66)

$\left|P^{REAC}-P_L^{REAC}\right|\precsim \frac{1}{2^{Lr/2}}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+\frac{1}{2^{L\left(t-1\right)/2}}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (67)

Proof. The single layer operator $V$ can be decomposed [4] as $V=V_0+L$ where the principal part $V_0$ is an elliptic operator $\langle V_{0}u\mathrm{,}u\rangle \succsim {\Vert u\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}^2$ and $L$ is a compact operator. Therefore, it is a Fredholm operator of zero index which is in addition injective [44] . As a consequence, Equations ((21) and (63)) are uniquely solvable. On account of the Strang lemma [12] with perturbed right hand side, the orthogonal projection ${\mathcal{Q}}_L\mathrm{<mo>\; :\; </mo>}{ℍ}^{-1/2}\left(\Gamma \right)\to {\mathcal{V}}_L\left(\Gamma \right)$ leads to

$\begin{array}{c}{\Vert \sigma -{\sigma }_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim \underset{w_L\in {\mathcal{V}}_L\left(\Gamma \right)}{min}{\Vert \sigma -w_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}+\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f+g\mathrm{,}{v}_L\rangle -\langle f_L+g\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\\ \le {\Vert \sigma -{\mathcal{Q}}_L\sigma \Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}+\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f-f_L\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\mathrm{.}\end{array}$

On the 2D-level $L$ for piecewise constant setup, one has

$\begin{array}{c}{\Vert \sigma -{\mathcal{Q}}_L\sigma \Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}=\underset{{\Vert \varphi \Vert }_{{ℍ}^{1/2}\left(\Gamma \right)}=1}{sup}\langle \left(I-{\mathcal{Q}}_L\right)\sigma \mathrm{,}\varphi \rangle =\underset{{\Vert \varphi \Vert }_{{ℍ}^{1/2}\left(\Gamma \right)}=1}{sup}\langle {\left(I-{\mathcal{Q}}_L\right)}^2\sigma \mathrm{,}\varphi \rangle \\ \precsim \underset{{\Vert \varphi \Vert }_{{ℍ}^{1/2}\left(\Gamma \right)}=1}{sup}{\Vert \left(I-{\mathcal{Q}}_L\right)\sigma \Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}{\Vert \left(I-{\mathcal{Q}}_L\right)\varphi \Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\\ \precsim \underset{{\Vert \varphi \Vert }_{{ℍ}^{1/2}\left(\Gamma \right)}=1}{sup}{2}^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}{2}^{-L/4}{\Vert \varphi \Vert }_{{ℍ}^{1/2}\left(\Gamma \right)}=2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}\mathrm{.}\end{array}$

In addition, apply the inverse estimate ${\Vert v_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\precsim 2^{L/4}{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}$ to obtain

$\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f-f_L\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\le \underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{{\Vert f-f_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}{\Vert v_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}$ (68)

$\le \underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{{\Vert f-f_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}{\Vert v_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}}{2^{-L/4}{\Vert v_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}}$ (69)

$\precsim 2^{L/4}{\Vert f-f_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\mathrm{.}$ (70)

A combination of the above relations yields therefore

${\Vert \sigma -{\sigma }_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim 2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L/4}{\Vert f-f_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\mathrm{.}$ (71)

On the other hand, for the estimation of ${\Vert f-f_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}$ , one uses the same reasoning as above where $K$ is now the compact operator. Hence, one has the estimation

${\Vert f-f_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\precsim \underset{q_L\in {\mathcal{V}}_L\left(\Gamma \right)}{min}{\Vert f-q_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\le {\Vert f-{\mathcal{Q}}_{L}f\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}$ (72)

$\precsim 2^{-Lt/2}{\Vert f\Vert }_{{ℍ}^t\left(\Gamma \right)}\precsim 2^{-Lt/2}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (73)

The a-priori error estimate for the apparent surface charge follows

${\Vert \sigma -{\sigma }_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim 2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{-L\left(t-1/2\right)/2}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (74)

Concerning the reaction potential, one obtains on account of the Cauchy- Schwarz inequality

$\left|P^{REAC}-P_L^{REAC}\right|=\left|\langle \sigma -{\sigma }_L\mathrm{,}\frac{1}{2}{\displaystyle \underset{k=1}{\overset{M}{\sum }}}\frac{1}{\text{4π}}\frac{1}{\Vert x_k-\cdot \Vert }\rangle \right|$ (75)

$\le {\Vert \sigma -{\sigma }_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}{\Vert \frac{1}{2}{\displaystyle \underset{k=1}{\overset{M}{\sum }}}\frac{1}{\text{4π}}\frac{1}{\Vert x_k-\cdot \Vert }\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\mathrm{.}$ (76)

We introduce the minimal distance between the nuclei coordinates ${\left\{x_k\right\}}_{k=1}^M$ and the molecular surface $\Gamma$ ,

$D=\underset{k=\mathrm{1,}\cdots \mathrm{,}M}{min}\underset{y\in \Gamma }{min}\Vert x_k-y\Vert >0.$ (77)

By using the orthogonal projection ${\mathcal{Q}}_L$ and an inverse inequality, one obtains

$\begin{array}{c}\left|P^{REAC}-P_L^{REAC}\right|\precsim M{D}^{-1}\left({\Vert \sigma -{\mathcal{Q}}_L\sigma \Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}+{\Vert {\sigma }_L-{\mathcal{Q}}_L\sigma \Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\right)\\ \precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+{\Vert {\sigma }_L-{\mathcal{Q}}_L\sigma \Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\\ \precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L/4}{\Vert {\sigma }_L-{\mathcal{Q}}_L\sigma \Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\\ \precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L/4}\left\{{\Vert {\mathcal{Q}}_L\sigma -\sigma \Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}+{\Vert \sigma -{\sigma }_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\right\}\\ \precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L/4}\left\{2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+{\Vert \sigma -{\sigma }_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\right\}\\ \precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L/4}{\Vert \sigma -{\sigma }_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\mathrm{.}\end{array}$

As a consequence, one concludes the following accuracy by applying (74),

$\begin{array}{c}\left|P^{REAC}-P_L^{REAC}\right|\precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L\mathrm{<mo>\; /\; </mo>4}}\left(2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{-L\left(t-1/2\right)/2}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\right)\\ \precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{-L\left(t-1\right)/2}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}\end{array}$

,

4. Multiscale for Solvation Energy

This section will be occupied by the treatment of the PCM by means of the wavelet technique. The matrix entries where the distance between the supports is beyond a level dependent threshold are annihilated. For a matrix $M={\left[M_{ij}\right]}_{i=1,\cdots ,P}^{j=1,\cdots ,Q}$ , define

$\mathcal{N}\left(M\right)\mathrm{<mo>\; :\; </mo>}=\text{card}\left\{M_{ij}\ne \mathrm{0\; <mo>\; :\; </mo>}i=\mathrm{1,}\cdots \mathrm{,}P\mathrm{<mo>\; ;\; </mo>}j=\mathrm{1,}\cdots \mathrm{,}Q\right\}\mathrm{,}$ (78)

${\Vert M\Vert }_{\infty }\mathrm{<mo>\; :\; </mo>}=\underset{j=1,\cdots ,Q}{max}\left\{{\displaystyle \underset{i=1}{\overset{P}{\sum }}}\left|M_{ij}\right|\right\},$ (79)

${\Vert M\Vert }_1\mathrm{<mo>\; :\; </mo>}=\underset{j=1,\cdots ,P}{max}\left\{{\displaystyle \underset{j=1}{\overset{Q}{\sum }}}\left|M_{ij}\right|\right\}\mathrm{,}$ (80)

${\Vert M\Vert }_2\mathrm{<mo>\; :\; </mo>}=\underset{x\ne 0}{sup}\Vert Mx\Vert /\Vert x\Vert \mathrm{.}$ (81)

For two levels $\mathcal{l},{\mathcal{l}}^{\prime }=0,1,\cdots ,L$ , a threshold ${\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ whose expression will be specified subsequently is used. We consider the matrix coefficients

${\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}={\displaystyle {\int }_{\square \times \square }}{\tilde{\mathcal{K}}}_{p\mathrm{,}q}\left(u\mathrm{,}v\right){\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\left(u\right){\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\left(v\right)\text{d}u\text{d}v$ (82)

where

${\tilde{\mathcal{K}}}_{p\mathrm{,}q}\left(u\mathrm{,}v\right)\mathrm{<mo>\; :\; </mo>}=\mathcal{K}\left[{\gamma }_p\left(u\right)\mathrm{,}{\gamma }_p\left(v\right)\right]G_p\left(u\right)G_q\left(v\right)\mathrm{.}$ (83)

They correspond to the matrix entries in (23) and (24) for the single layer kernel $\mathcal{K}={\mathcal{K}}^{SL}$ and the double layer kernel $\mathcal{K}={\mathcal{K}}^{DL}$ after transformation onto $\square$ . The following level depended truncation procedure for each pair of levels $\left(\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }\right)$ is applied:

${\tilde{\mathcal{A}}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\mathrm{<mo>\; :\; </mo>}=\{\begin{array}{ll}{\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\hfill & \text{if}\text{dist}\left[S\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)\right]\le {\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\hfill \\ 0\hfill & \text{otherwise}\mathrm{.}\hfill \end{array}$ (84)

The block matrices ${\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ and ${\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ have respectively the entries ${\left[{\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\right]}_{i\mathrm{,}{i}^{\prime }}$ and ${\left[{\tilde{\mathcal{A}}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\right]}_{i\mathrm{,}{i}^{\prime }}$ such that we have blockwise:

$\mathcal{A}=\left[\begin{array}{ccc}{\mathcal{A}}_{\mathrm{0,0}}& \cdots & {\mathcal{A}}_{\mathrm{0,}L}\\ ⋮& & ⋮\\ {\mathcal{A}}_{L\mathrm{,0}}& \cdots & {\mathcal{A}}_{L\mathrm{,}L}\end{array}\right]\mathrm{,}\tilde{\mathcal{A}}=\left[\begin{array}{ccc}{\tilde{\mathcal{A}}}_{\mathrm{0,0}}& \cdots & {\tilde{\mathcal{A}}}_{\mathrm{0,}L}\\ ⋮& & ⋮\\ {\tilde{\mathcal{A}}}_{L\mathrm{,0}}& \cdots & {\tilde{\mathcal{A}}}_{L\mathrm{,}L}\end{array}\right]\mathrm{.}$ (85)

For the single layer and the double layer, define respectively

${\theta }^{SL}\mathrm{<mo>\; :\; </mo>}=1{\theta }^{DL}\mathrm{<mo>\; :\; </mo>}=2{\theta }_0^{SL}\mathrm{<mo>\; :\; </mo>}=-1{\theta }_0^{DL}\mathrm{<mo>\; :\; </mo>}=0.$ (86)

For the involved kernels, one has on account of the Calderon-Zygmund inequality:

$\left|\frac{{\partial }^{\left|\alpha \right|}}{\partial x^{\alpha }}\frac{{\partial }^{\left|\beta \right|}}{\partial y^{\beta }}\mathcal{K}\left(x\mathrm{,}y\right)\right|\le \frac{C_{\alpha ,\beta }}{{\Vert x-y\Vert }^{\theta +\left|\alpha \right|+\left|\beta \right|}}\mathrm{.}$ (87)

Lemma 1. Fix any real parameters $\left(s\mathrm{,}p\mathrm{,}q\right)$ and consider two levels $\mathcal{l},{\mathcal{l}}^{\prime }=0,\cdots ,L$ . Suppose (see (39), (40),(41))

${\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\mathrm{<mo>\; :\; </mo>}=max\{\begin{array}{l}{\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\mathrm{<mo>\; :\; </mo>}={{\displaystyle exp}}_2\left\{\frac{\left(L-\mathcal{l}-{\mathcal{l}}^{\prime }\right)s}{2\left({\theta }_0+2\right)}-\frac{\mathcal{l}p}{2\left({\theta }_0+2\right)}-\frac{{\mathcal{l}}^{\prime }q}{2\left({\theta }_0+2\right)}\right\}\hfill \\ \lambda \left({\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}\right)\hfill \end{array}$ (88)

where $\lambda >1$ and ${\Delta }_{\mathcal{l}}$ denotes the diameter of the supports of any wavelet ${\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}$ on level $\mathcal{l}$ . Then, the following error estimates in matrix norms are fulfilled:

${\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_{\infty }=\mathcal{O}\left(L{2}^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-{\mathcal{l}}^{\prime }}\right)$ (89)

${\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_1=\mathcal{O}\left(L{2}^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-\mathcal{l}}\right)$ (90)

${\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_2=\mathcal{O}\left(L{2}^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-\mathcal{l}/2}{2}^{-{\mathcal{l}}^{\prime }/2}\right)\mathrm{.}$ (91)

By using the thresholds $\left({\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\right)$ for all $\mathcal{l},{\mathcal{l}}^{\prime }=0,\cdots ,L$ as in (88) such that the parameters $\left(s\mathrm{,}p\mathrm{,}q\right)$ satisfy

$s,p,q>0,s+p+q={\theta }_0+2,$ (92)

one has for $u\in {ℍ}^{\gamma }\left(\Gamma \right)$ and $v_L\in {\mathcal{V}}_L\left(\Gamma \right)$ in the case of the single layer operator

$\left|\langle \left(V-\tilde{V}\right){\mathcal{Q}}_{L}u\mathrm{,}{v}_L\rangle \right|\precsim L{\left(2^{-L/2}\right)}^{\gamma +1/2}{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}$ (93)

while the double layer operator yields

$\left|\langle \left(K-\tilde{K}\right){\mathcal{Q}}_{L}u\mathrm{,}{v}_L\rangle \right|\precsim L{\left(2^{-L/2}\right)}^{\gamma }{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}{\Vert v_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\mathrm{.}$ (94)

Proof. The value of ${\tilde{\mathcal{K}}}_{p\mathrm{,}q}\left(u\mathrm{,}v\right)$ from (83) will be examined on a pair of patches ${\Gamma }_p\times {\Gamma }_q$ . The indices $\left(p\mathrm{,}q\right)$ are dropped to simplify the notation. Consider a tensor product wavelet basis ${\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\otimes {\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}$ and fix $\left(\tilde{u}\mathrm{,}\tilde{v}\right)=\left[\left({\tilde{u}}_1\mathrm{,}{\tilde{u}}_2\right)\mathrm{,}\left({\tilde{v}}_1\mathrm{,}{\tilde{v}}_2\right)\right]$ as the midpoints of $S\left(\mathcal{l}\mathrm{,}i\right)$ and $S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)$ . By considering any $\left(u\mathrm{,}v\right)=\left[\left(u_1\mathrm{,}{u}_2\right)\mathrm{,}\left(v_1\mathrm{,}{v}_2\right)\right]\in S\left(\mathcal{l}\mathrm{,}i\right)\times S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)$ , the Taylor expansion leads to

#Math_316# (95)

for some $\left(u^{\mathrm{<mo>\; *\; </mo>}}\mathrm{,}{v}^{\mathrm{<mo>\; *\; </mo>}}\right)\in S\left(\mathcal{l}\mathrm{,}i\right)\times S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)$ . For the first term, by taking the integration over $\square \times \square$ , one obtains

${\displaystyle {\int }_{\square \times \square }}\tilde{\mathcal{K}}\left(\tilde{u}\mathrm{,}\tilde{v}\right){\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\left(u\right){\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\left(v\right)\text{d}u\text{d}v=\tilde{\mathcal{K}}\left(\tilde{u}\mathrm{,}\tilde{v}\right){\displaystyle {\int }_{\square }}{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\left(u\right)\text{d}u{\displaystyle {\int }_{\square }}{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\left(v\right)\text{d}v$ (96)

which is zero due to (33). We will consider now the value of ${\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}$ in the case ${\Delta }_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\mathrm{<mo>\; :\; </mo>}=\text{dist}\left[S\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)\right]\ge {\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ . For the summation over $\left|\alpha \right|=\left|\beta \right|=1$ from (95), supposing that the Jacobians of the mappings ${\gamma }_p$ and ${\gamma }_q$ are bounded, one has for $x^{\mathrm{<mo>\; *\; </mo>}}={\gamma }_p\left(u^{\mathrm{<mo>\; *\; </mo>}}\right)$ and $y^{\mathrm{<mo>\; *\; </mo>}}={\gamma }_q\left(v^{\mathrm{<mo>\; *\; </mo>}}\right)$

$\left|\frac{{\partial }^{\left|\alpha \right|}}{\partial u^{\alpha }}\frac{{\partial }^{\left|\beta \right|}}{\partial v^{\beta }}\tilde{\mathcal{K}}\left(u^{\mathrm{<mo>\; *\; </mo>}}\mathrm{,}{v}^{\mathrm{<mo>\; *\; </mo>}}\right)\right|\le \frac{C_{\alpha ,\beta }}{{\Vert x^{\mathrm{<mo>\; *\; </mo>}}-y^{\mathrm{<mo>\; *\; </mo>}}\Vert }^{\theta +2}}\precsim {\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-\left({\theta }_0+2+1/L\right)}\frac{1}{{\Vert x^{\mathrm{<mo>\; *\; </mo>}}-y^{\mathrm{<mo>\; *\; </mo>}}\Vert }^{\theta -{\theta }_0-1/L}}\mathrm{.}$ (97)

Since ${\Delta }_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\ge {\delta }_{\left(\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }\right)}\ge \lambda \left({\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}\right)$ , one has the next inequalities for $x={\gamma }_p\left(u\right)$ and $y={\gamma }_q\left(v\right)$ :

$\Vert x-y\Vert \le \Vert x-x^{\mathrm{<mo>\; *\; </mo>}}\Vert +\Vert x^{\mathrm{<mo>\; *\; </mo>}}-y^{\mathrm{<mo>\; *\; </mo>}}\Vert +\Vert y-y^{\mathrm{<mo>\; *\; </mo>}}\Vert \le \Vert x^{\mathrm{<mo>\; *\; </mo>}}-y^{\mathrm{<mo>\; *\; </mo>}}\Vert +{\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}\mathrm{,}$ (98)

${\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}\le \frac{1}{\lambda }{\Delta }_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\le \frac{1}{\lambda }\Vert x-y\Vert \mathrm{,}$ (99)

$\Vert x^{\mathrm{<mo>\; *\; </mo>}}-y^{\mathrm{<mo>\; *\; </mo>}}\Vert \ge \left(1-1/\lambda \right)\Vert x-y\Vert \succsim \Vert x-y\Vert$ (100)

which, when inserted into (97), leads to

$\left|\frac{{\partial }^{\left|\alpha \right|}}{\partial u^{\alpha }}\frac{{\partial }^{\left|\beta \right|}}{\partial v^{\beta }}\tilde{\mathcal{K}}\left(u^{\mathrm{<mo>\; *\; </mo>}}\mathrm{,}{v}^{\mathrm{<mo>\; *\; </mo>}}\right)\right|\precsim {\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-\left({\theta }_0+2+1/L\right)}\frac{1}{{\Vert x-y\Vert }^{\theta -{\theta }_0-1/L}}\mathrm{.}$ (101)

A combination of (95), (96) and (101) leads to

$\begin{array}{c}{\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\le {\displaystyle \underset{\left|\alpha \right|=1,\left|\beta \right|=1}{\sum }}\frac{1}{\alpha !}\frac{1}{\beta !}{\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-\left({\theta }_0+2+1/L\right)}max\left|{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\right|\underset{u_1\mathrm{,}{u}_2\in S\left(\mathcal{l}\mathrm{,}i\right)}{max}{\left|u_1-u_2\right|}^{\alpha }\\ \times max\left|{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\right|\underset{v_1\mathrm{,}{v}_2\in S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}{max}{\left|v_1-v_2\right|}^{\beta }{\displaystyle {\int }_{S\left(\mathcal{l}\mathrm{,}i\right)\times S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}}\frac{{\text{d}}_{\Gamma }x{\text{d}}_{\Gamma }y}{{\Vert x-y\Vert }^{\theta -{\theta }_0-1/L}}\mathrm{.}\end{array}$

On account of the properties (51)-(54), one deduces

${\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\precsim {\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-\left({\theta }_0+2+1/L\right)}{\displaystyle {\int }_{S\left(\mathcal{l}\mathrm{,}i\right)\times S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}}\frac{{\text{d}}_{\Gamma }x{\text{d}}_{\Gamma }y}{{\Vert x-y\Vert }^{\theta -{\theta }_0-1/L}}\mathrm{.}$ (102)

The function $\left(x\mathrm{,}y\right)\mapsto ex{p}_2\left[-\frac{\left(1-x-y\right)s}{2\left({\theta }_0+2\right)}+\frac{xp}{2\left({\theta }_0+2\right)}+\frac{yq}{2\left({\theta }_0+2\right)}\right]$ is continuous on the compact $\left[\mathrm{0,1}\right]\times \left[\mathrm{0,1}\right]$ . Hence, the value

${\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-1/L}=ex{p}_2\left[-\frac{\left(1-\mathcal{l}/L-{\mathcal{l}}^{\prime }/L\right)s}{2\left({\theta }_0+2\right)}+\frac{\left(\mathcal{l}/L\right)p}{2\left({\theta }_0+2\right)}+\frac{\left({\mathcal{l}}^{\prime }/L\right)q}{2\left({\theta }_0+2\right)}\right]$ (103)

is bounded independently of $\mathcal{l},{\mathcal{l}}^{\prime }=0,\cdots ,L$ . As a consequence, use the expression of ${\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ to obtain

${\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-\left({\theta }_0+2+1/L\right)}={\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-\left({\theta }_0+2\right)}{\tilde{\delta }}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{-1/L}\precsim 2^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}$ (104)

and hence,

${\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\precsim 2^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{\displaystyle {\int }_{S\left(\mathcal{l}\mathrm{,}i\right)\times S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}}\frac{{\text{d}}_{\Gamma }x{\text{d}}_{\Gamma }y}{{\Vert x-y\Vert }^{\theta -{\theta }_0-1/L}}\mathrm{.}$ (105)

As a consequence, one deduces the norm estimate

$\begin{array}{c}{\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_{\infty }=\underset{i^{\prime }\in Z_{{\mathcal{l}}^{\prime }}}{max}\left\{{\displaystyle \underset{i\in Z_{\mathcal{l}}}{\sum }}\left|{\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}-{\tilde{\mathcal{A}}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\right|\right\}=\underset{i^{\prime }\in Z_{{\mathcal{l}}^{\prime }}}{max}\left\{{\displaystyle \underset{\begin{array}{l}i\in Z_{\mathcal{l}}\\ {\Delta }_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\ge {\delta }_{\mathcal{l},{\mathcal{l}}^{\prime }}\end{array}}{\sum }\left|{\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\right|}\right\}\\ \precsim 2^{\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)\frac{s}{2}}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}\underset{i^{\prime }\in Z_{{\mathcal{l}}^{\prime }}}{max}\left\{{\displaystyle {\int }_{S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}}\left[{\displaystyle \underset{\begin{array}{l}i\in Z_{\mathcal{l}}\\ {\Delta }_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\ge {\delta }_{\mathcal{l},{\mathcal{l}}^{\prime }}\end{array}}{\sum }{\displaystyle {\int }_{S\left(\mathcal{l}\mathrm{,}i\right)}}\frac{{\text{d}}_{\Gamma }y}{{\Vert x-y\Vert }^{\theta -{\theta }_0-1/L}}}\right]{\text{d}}_{\Gamma }x\right\}\end{array}$

For any $x\in \Gamma$ , we have

${\displaystyle \underset{\begin{array}{c}i\in Z_{\mathcal{l}}\\ {\Delta }_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\ge {\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\end{array}}{\sum }}{\displaystyle {\int }_{S\left(\mathcal{l}\mathrm{,}i\right)}}\frac{{\text{d}}_{\Gamma }y}{{\Vert x-y\Vert }^{\theta -{\theta }_0-1/L}}\le {\displaystyle {\int }_{\Gamma }}\frac{{\text{d}}_{\Gamma }y}{{\Vert x-y\Vert }^{\theta -{\theta }_0-1/L}}={\displaystyle {\int }_{\Gamma }}\frac{{\text{d}}_{\Gamma }y}{{\Vert x-y\Vert }^{2-1/L}}=\mathcal{O}\left(L\right)$

where the last equality uses a local mapping and polar coordinates. By using the measure property (51), obtain ${\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_{\infty }\precsim L{2}^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-{\mathcal{l}}^{\prime }}$ . The 1-norm is processed in an analogous manner and therefore

${\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_1=\underset{i\in Z_{\mathcal{l}}}{max}\left\{{\displaystyle \underset{i^{\prime }\in Z_{{\mathcal{l}}^{\prime }}}{\sum }}\left|{\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}-{\tilde{\mathcal{A}}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\right|\right\}\precsim L{2}^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-\mathcal{l}}$

${\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_2^2\le {\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_{\infty }{\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_1\precsim {\left[L{2}^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-\mathcal{l}/2}{2}^{-{\mathcal{l}}^{\prime }/2}\right]}^2\mathrm{.}$

Consider the orthogonal projection ${\mathcal{Q}}_{\mathcal{l}}$ onto ${\mathcal{V}}_{\mathcal{l}}\left(\Gamma \right)$ and introduce

$S\mathrm{<mo>\; :\; </mo>}=\langle \left(\mathcal{A}-\tilde{\mathcal{A}}\right){\mathcal{Q}}_{L}u\mathrm{,}{\mathcal{Q}}_{L}v\rangle$ (106)

in which

${\mathcal{Q}}_{L}u={\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}\left({\mathcal{Q}}_{\mathcal{l}}-{\mathcal{Q}}_{\mathcal{l}-1}\right)u\text{where}{\mathcal{Q}}_{-1}\mathrm{<mo>\; :\; </mo>}=0.$ (107)

Blockwise according to the levels

$\mathcal{A}-\tilde{\mathcal{A}}=\left[\begin{array}{ccc}\left({\mathcal{A}}_{\mathrm{0,0}}-{\tilde{\mathcal{A}}}_{\mathrm{0,0}}\right)& \cdots & \left({\mathcal{A}}_{\mathrm{0,}L}-{\tilde{\mathcal{A}}}_{\mathrm{0,}L}\right)\\ ⋮& & ⋮\\ \left({\mathcal{A}}_{L\mathrm{,0}}-{\tilde{\mathcal{A}}}_{L\mathrm{,0}}\right)& \cdots & \left({\mathcal{A}}_{L\mathrm{,}L}-{\tilde{\mathcal{A}}}_{L\mathrm{,}L}\right)\end{array}\right]$ (108)

$S=\langle \left(\mathcal{A}-\tilde{\mathcal{A}}\right){\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}\left({\mathcal{Q}}_{\mathcal{l}}-{\mathcal{Q}}_{\mathcal{l}-1}\right)u,{\displaystyle \underset{{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}\left({\mathcal{Q}}_{{\mathcal{l}}^{\prime }}-{\mathcal{Q}}_{{\mathcal{l}}^{\prime }-1}\right)v\rangle$ (109)

$\simeq {\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{\displaystyle \underset{{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}\langle \left(\mathcal{A}-\tilde{\mathcal{A}}\right){\displaystyle \underset{i\in Z_{\mathcal{l}}}{\sum }}\langle u,{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\rangle {\psi }_{\left(\mathcal{l}\mathrm{,}i\right)},{\displaystyle \underset{i^{\prime }\in Z_{{\mathcal{l}}^{\prime }}}{\sum }}\langle v,{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\rangle {\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\rangle$ (110)

$={\displaystyle \underset{\mathcal{l},{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{\displaystyle \underset{i,i^{\prime }}{\sum }}\langle \left(\mathcal{A}-\tilde{\mathcal{A}}\right){\psi }_{\left(\mathcal{l}\mathrm{,}i\right)},{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\rangle \langle u,{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\rangle \langle v,{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\rangle$ (111)

$={\displaystyle \underset{\mathcal{l},{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{\displaystyle \underset{i,i^{\prime }}{\sum }}\left[{\mathcal{A}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}-{\tilde{\mathcal{A}}}_{\left(\mathcal{l}\mathrm{,}i\right)\mathrm{,}\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\right]\langle u,{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\rangle \langle v,{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\rangle .$ (112)

By utilizing the next estimates

${\Vert {\left\{\langle u,{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\rangle \right\}}_i\Vert }_{{\mathcal{l}}^2}\precsim 2^{-\mathcal{l}\gamma /2}{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}\mathrm{,}{\Vert {\left\{\langle v,{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\rangle \right\}}_{i^{\prime }}\Vert }_{{\mathcal{l}}^2}\precsim 2^{-{\mathcal{l}}^{\prime }{\gamma }^{\prime }/2}{\Vert v\Vert }_{{ℍ}^{{\gamma }^{\prime }}\left(\Gamma \right)}\mathrm{,}$ (113)

one obtains, based on the Cauchy-Schwarz inequality and (91),

$\begin{array}{c}\left|S\right|\precsim {\displaystyle \underset{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{\Vert {\mathcal{A}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}-{\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\Vert }_2{\Vert {\left\{\langle u,{\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}\rangle \right\}}_i\Vert }_{{\mathcal{l}}^2}{\Vert \left\{\langle v,{\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}\rangle \right\}\Vert }_{{\mathcal{l}}^2}\\ \precsim L{\displaystyle \underset{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{2}^{\frac{1}{2}\left(\mathcal{l}+{\mathcal{l}}^{\prime }-L\right)s}{2}^{\mathcal{l}p/2}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-\mathcal{l}/2}{2}^{-{\mathcal{l}}^{\prime }/2}{2}^{-\mathcal{l}\gamma /2}{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}{2}^{-{\mathcal{l}}^{\prime }{\gamma }^{\prime }/2}{\Vert v\Vert }_{{ℍ}^{{\gamma }^{\prime }}\left(\Gamma \right)}\\ \precsim L{2}^{-Ls/2}{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}{\Vert v\Vert }_{{ℍ}^{{\gamma }^{\prime }}\left(\Gamma \right)}\left[{\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{2}^{\frac{1}{2}\mathcal{l}s}{2}^{\mathcal{l}p/2}{2}^{-\mathcal{l}/2}{2}^{-\mathcal{l}\gamma /2}\right]\left[{\displaystyle \underset{{\mathcal{l}}^{\prime }\mathrm{<mo>\; =\; </mo>0}}{\overset{L}{\sum }}}{2}^{\frac{1}{2}{\mathcal{l}}^{\prime }s}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-{\mathcal{l}}^{\prime }/2}{2}^{-{\mathcal{l}}^{\prime }{\gamma }^{\prime }/2}\right]\mathrm{.}\end{array}$

The use of the sum of geometric sequence leads to

${\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{2}^{\frac{1}{2}\mathcal{l}s}{2}^{\mathcal{l}p/2}{2}^{-\mathcal{l}/2}{2}^{-\mathcal{l}\gamma /2}=\mathcal{O}\left(2^{\frac{1}{2}Ls}{2}^{Lp/2}{2}^{-L/2}{2}^{-L\gamma /2}\right)\mathrm{,}$

${\displaystyle \underset{{\mathcal{l}}^{\prime }\mathrm{<mo>\; =\; </mo>0}}{\overset{L}{\sum }}}{2}^{\frac{1}{2}{\mathcal{l}}^{\prime }s}{2}^{{\mathcal{l}}^{\prime }q/2}{2}^{-{\mathcal{l}}^{\prime }/2}{2}^{-{\mathcal{l}}^{\prime }{\gamma }^{\prime }/2}=\mathcal{O}\left(2^{\frac{1}{2}Ls}{2}^{Lq/2}{2}^{-L/2}{2}^{-L{\gamma }^{\prime }/2}\right)$

and hence

$\left|S\right|\precsim L{2}^{-Ls/2}{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}{\Vert v\Vert }_{{ℍ}^{{\gamma }^{\prime }}\left(\Gamma \right)}{2}^{Ls}{2}^{Lp/2}{2}^{Lq/2}{2}^{-L}{2}^{-L\gamma /2}{2}^{-L{\gamma }^{\prime }/2}\mathrm{.}$

For piecewise constant setup, obtain for the constant parameters $s\mathrm{,}p\mathrm{,}q$

$\left|\langle \left(\mathcal{A}-\tilde{\mathcal{A}}\right){\mathcal{Q}}_{L}u\mathrm{,}{\mathcal{Q}}_{L}v\rangle \right|\precsim L{\left(2^{-L/2}\right)}^{2-s-p-q+\gamma +{\gamma }^{\prime }}{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}{\Vert v\Vert }_{{ℍ}^{{\gamma }^{\prime }}\left(\Gamma \right)}\mathrm{.}$ (114)

For the single layer $V$ , set ${\gamma }^{\prime }=-1/2$ and use $s+p+q={\theta }_0^{SL}+2$ to obtain from (86)

$2-s-p-q+\gamma +{\gamma }^{\prime }=2-\left({\theta }_0^{SL}+2\right)+\gamma -1/2=\gamma +1/2\mathrm{,}$

$\frac{\left|\langle \left(V-\tilde{V}\right){\mathcal{Q}}_{L}u\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}=\frac{\left|\langle \left(V-\tilde{V}\right){\mathcal{Q}}_{L}u\mathrm{,}{v}_L\rangle \right|}{{\Vert {\mathcal{Q}}_L{v}_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\precsim L{\left(2^{-L/2}\right)}^{\gamma +1/2}{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}\mathrm{.}$

Similarly for the double layer, setting ${\gamma }^{\prime }=0$ results in

$2-s-p-q+\gamma +{\gamma }^{\prime }=2-\left({\theta }_0^{DL}+2\right)+\gamma =\gamma$ (115)

$\left|\langle \left(K-\tilde{K}\right){\mathcal{Q}}_{L}u\mathrm{,}{v}_L\rangle \right|{\Vert v_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\precsim L{\left(2^{-L/2}\right)}^{\gamma }{\Vert u\Vert }_{{ℍ}^{\gamma }\left(\Gamma \right)}\mathrm{.}$ (116)

,

Theorem 2. Consider the truncated PCM equations on level $L$ :

$\langle \tilde{V}{\tilde{\sigma }}_L\mathrm{,}{v}_L\rangle =\langle {\tilde{f}}_L\mathrm{,}{v}_L\rangle +\langle g\mathrm{,}{v}_L\rangle ,\forall v_L\in {\mathcal{V}}_L\left(\Gamma \right)$ (117)

$\langle \left(I-\tilde{K}\right){\tilde{f}}_L\mathrm{,}{v}_L\rangle =\langle m\mathrm{,}{v}_L\rangle ,\forall v_L\in {\mathcal{V}}_L\left(\Gamma \right)$ (118)

where the operators $\tilde{V}$ and $\tilde{K}$ are obtained from $V$ and $K$ by using the threshold ${\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{SL}$ and ${\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^{DL}$ respectively.

Suppose the constant parameters $\left(s\mathrm{,}p\mathrm{,}q\right)$ are chosen such that

$s,p,q>0,s+p+q={\theta }_0+2.$ (119)

For $P$ patches, the numbers of nonzero matrix coefficients are then reduced from $\mathcal{O}\left[{\left(P{2}^L\right)}^2\right]$ to $\mathcal{N}\left(\tilde{V}\right)=\mathcal{O}\left[L\left(P{2}^L\right)\right]$ and $\mathcal{N}\left(\tilde{K}\right)=\mathcal{O}\left[L\left(P{2}^L\right)\right]$ .

In addition, the accuracy between $\sigma$ in (21) and ${\tilde{\sigma }}_L$ in (117) for $r\ge 1/2$ and $t\ge -1/2$ verifies

${\Vert \sigma -{\tilde{\sigma }}_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim \frac{L}{2^{L\left(r+1/2\right)/2}}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+\frac{L}{2^{L\left(t-1/2\right)/2}}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (120)

Further, the error between the reaction potential $P^{REAC}$ from (18) and

${\tilde{P}}_L^{REAC}\mathrm{<mo>\; :\; </mo>}=\frac{1}{2}{\displaystyle \underset{k=1}{\overset{M}{\sum }}}\frac{1}{\text{4π}}{\displaystyle {\int }_{\Gamma }}\frac{{\tilde{\sigma }}_L\left(y\right)}{\Vert x_k-y\Vert }{\text{d}}_{\Gamma }y\text{fulfills}$ (121)

$\left|P^{REAC}-{\tilde{P}}_L^{REAC}\right|\precsim \frac{L}{2^{Lr/2}}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+\frac{L}{2^{L\left(t-1\right)/2}}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (122)

Proof. We consider first the number of nonzero entries of the compressed matrices $\tilde{V}$ and $\tilde{K}$ . For computing $\mathcal{N}\left({\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\right)$ for each fixed level pair $\left(\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }\right)$ , consider a basis function ${\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}$ on level $\mathcal{l}$ . According to [40] , the number of basis functions ${\psi }_{\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)}$ on level ${\mathcal{l}}^{\prime }$ whose support $S\left({\mathcal{l}}^{\prime }\mathrm{,}{i}^{\prime }\right)$ is within a distance of ${\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}+{\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ from the support $S\left(\mathcal{l}\mathrm{,}i\right)$ is in the worst case $2^{{\mathcal{l}}^{\prime }}{\left({\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}+{\delta }_{\mathcal{l},{\mathcal{l}}^{\prime }}\right)}^2$ . For that, Lemma 3.2. of [40] uses an argument using a sphere of radius ${\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}+{\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ centered at a point located upon the support $S\left(\mathcal{l}\mathrm{,}i\right)$ . There are $2^{\mathcal{l}}$ basis functions ${\psi }_{\left(\mathcal{l}\mathrm{,}i\right)}$ on level $\mathcal{l}$ for each patch ${\Gamma }_p$ where $p=1,\cdots ,P$ . The number of the interlevel nonzero entries on all patches is thus

$\mathcal{N}\left({\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\right)\precsim P{2}^{\mathcal{l}}{2}^{{\mathcal{l}}^{\prime }}{\left({\Delta }_{\mathcal{l}}+{\Delta }_{{\mathcal{l}}^{\prime }}+{\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\right)}^2\mathrm{.}$ (123)

By summing over all levels

$\begin{array}{c}\mathcal{N}\left(\tilde{\mathcal{A}}\right)={\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{\displaystyle \underset{{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}\mathcal{N}\left({\tilde{\mathcal{A}}}_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}\right)\precsim {\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{\displaystyle \underset{{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}P{2}^{\left(\mathcal{l}+{\mathcal{l}}^{\prime }\right)}\left({\Delta }_{\mathcal{l}}^2+{\Delta }_{{\mathcal{l}}^{\prime }}^2+{\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^2\right)\\ \precsim \left\{P{\displaystyle \underset{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{2}^{\left(\mathcal{l}+{\mathcal{l}}^{\prime }\right)}\left(2^{-\mathcal{l}}+2^{-{\mathcal{l}}^{\prime }}\right)\right\}+\left\{P{\displaystyle \underset{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{2}^{\left(\mathcal{l}+{\mathcal{l}}^{\prime }\right)}{\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}^2\right\}\\ ={\mathcal{N}}_1+{\mathcal{N}}_2\mathrm{.}\end{array}$

The first sum on all patches verifies

${\mathcal{N}}_1=P{\displaystyle \underset{{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}L{2}^{{\mathcal{l}}^{\prime }}+P{\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}L{2}^{\mathcal{l}}=\mathcal{O}\left[L\left(P{2}^L\right)\right].$ (124)

For the second sum ${\mathcal{N}}_2$ , the expression of the threshold ${\delta }_{\mathcal{l}\mathrm{,}{\mathcal{l}}^{\prime }}$ leads to

$\begin{array}{c}{\mathcal{N}}_2=P{\displaystyle \underset{\mathcal{l},{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{2}^{\left(\mathcal{l}+{\mathcal{l}}^{\prime }\right)}{\exp }_2\left\{\frac{\left(L-\mathcal{l}-{\mathcal{l}}^{\prime }\right)s}{\left({\theta }_0+2\right)}-\frac{\mathcal{l}p}{\left({\theta }_0+2\right)}-\frac{{\mathcal{l}}^{\prime }q}{\left({\theta }_0+2\right)}\right\}\\ =P{2}^{Ls/\left({\theta }_0+2\right)}{\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{\exp }_2\left\{\mathcal{l}\left[1-\frac{s}{{\theta }_0+2}-\frac{p}{{\theta }_0+2}\right]\right\}{\displaystyle \underset{{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{\exp }_2\left\{{\mathcal{l}}^{\prime }\left[1-\frac{s}{{\theta }_0+2}-\frac{q}{{\theta }_0+2}\right]\right\}.\end{array}$

By using sum of geometric sequence, deduce from $s+p+q={\theta }_0+2$

$\begin{array}{l}{\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{\exp }_2\left\{\mathcal{l}\left[1-\frac{s}{{\theta }_0+2}-\frac{p}{{\theta }_0+2}\right]\right\}\\ \precsim {\displaystyle \underset{\mathcal{l}=0}{\overset{L}{\sum }}}{\exp }_2\left\{\mathcal{l}\left[\frac{q}{{\theta }_0+2}\right]\right\}=\mathcal{O}\left({\exp }_2\left\{L\left[\frac{q}{{\theta }_0+2}\right]\right\}\right).\end{array}$

In a similar manner, one obtains

${\displaystyle \underset{{\mathcal{l}}^{\prime }=0}{\overset{L}{\sum }}}{\exp }_2\left\{{\mathcal{l}}^{\prime }\left[1-\frac{s}{{\theta }_0+2}-\frac{q}{{\theta }_0+2}\right]\right\}=\mathcal{O}\left({\exp }_2\left\{L\left[\frac{p}{{\theta }_0+2}\right]\right\}\right).$

Hence, obtain for the second sum by using $s+p+q={\theta }_0+2$

${\mathcal{N}}_2\precsim Pex{p}_2\left\{L\left[\frac{s+p+q}{{\theta }_0+2+1/L}\right]\right\}=\mathcal{O}\left(P{2}^L\right)$ (125)

Hence, $\mathcal{N}\left(\tilde{\mathcal{A}}\right)\le {\mathcal{N}}_1+{\mathcal{N}}_2=\mathcal{O}\left[L\left(P{2}^L\right)\right]+\mathcal{O}\left(P{2}^L\right)=\mathcal{O}\left[L\left(P{2}^L\right)\right]$ for all $P$ patches. Concerning the first Equation (117), the next Strang relation [12] holds

$\begin{array}{c}{\Vert \sigma -{\tilde{\sigma }}_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim \underset{w_L\in {\mathcal{V}}_L\left(\Gamma \right)}{min}\left\{{\Vert \sigma -w_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}+\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle \left(V-\tilde{V}\right)w_L\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\right\}\\ +\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f+g\mathrm{,}{v}_L\rangle -\langle {\tilde{f}}_L+g\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\\ \precsim {\Vert \sigma -{\mathcal{Q}}_L\sigma \Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}+\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle \left(V-\tilde{V}\right){\mathcal{Q}}_L\sigma \mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\\ +\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f-{\tilde{f}}_L\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\end{array}$

$\begin{array}{c}\precsim 2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle \left(V-\tilde{V}\right){\mathcal{Q}}_L\sigma \mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\\ +\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f-{\tilde{f}}_L\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\mathrm{.}\end{array}$

By inserting (115) in the last estimate, the next inequality is fulfilled

${\Vert \sigma -{\tilde{\sigma }}_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim 2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+L{2}^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f-{\tilde{f}}_L\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\mathrm{.}$

In the same manner as (70), one has

$\underset{v_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle f-{\tilde{f}}_L\mathrm{,}{v}_L\rangle \right|}{{\Vert v_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}}\precsim 2^{L/4}{\Vert f-{\tilde{f}}_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}$ (126)

and hence

${\Vert \sigma -{\tilde{\sigma }}_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim L{2}^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L/4}{\Vert f-{\tilde{f}}_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\mathrm{.}$ (127)

On the other hand, for the second Equation (118) by applying the Strang lemma [12] for the second time

$\begin{array}{c}{\Vert f-{\tilde{f}}_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\precsim \underset{q_L\in {\mathcal{V}}_L\left(\Gamma \right)}{min}\left\{{\Vert f-q_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}+\underset{p_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle \left[\left(I-K\right)-\left(I-\tilde{K}\right)\right]q_L\mathrm{,}{p}_L\rangle \right|}{{\Vert p_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}}\right\}\\ \precsim {\Vert f-{\mathcal{Q}}_{L}f\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}+\underset{p_L\in {\mathcal{V}}_L\left(\Gamma \right)}{sup}\frac{\left|\langle \left(K-\tilde{K}\right){\mathcal{Q}}_{L}f\mathrm{,}{p}_L\rangle \right|}{{\Vert p_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}}\mathrm{.}\end{array}$

By using (116), deduce

${\Vert f-{\tilde{f}}_L\Vert }_{{\mathbb{L}}^2\left(\Gamma \right)}\precsim 2^{-Lt/2}{\Vert f\Vert }_{{ℍ}^t\left(\Gamma \right)}+L{2}^{-Lt/2}{\Vert f\Vert }_{{ℍ}^t\left(\Gamma \right)}$ (128)

$\precsim L{2}^{-Lt/2}{\Vert f\Vert }_{{ℍ}^t\left(\Gamma \right)}\precsim L{2}^{-Lt/2}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (129)

Based on the last estimate and (127), deduce

${\Vert \sigma -{\tilde{\sigma }}_L\Vert }_{{ℍ}^{-1/2}\left(\Gamma \right)}\precsim L{2}^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+L{2}^{-L\left(t-1/2\right)/2}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}$ (130)

Concerning the accuracy with respect to the electrostatic reaction potential ${\tilde{P}}^{REAC}$ , one proceeds as in (75) and (76) in order to obtain

$\begin{array}{c}\left|P^{REAC}-{\tilde{P}}_L^{REAC}\right|\precsim 2^{-Lr/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{L/4}L\left(2^{-L\left(r+1/2\right)/2}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+2^{-L\left(t-1/2\right)/2}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\right)\\ \precsim \frac{L}{2^{Lr\mathrm{<mo>\; /\; </mo>2}}}{\Vert \sigma \Vert }_{{ℍ}^r\left(\Gamma \right)}+\frac{L}{2^{L\left(t-1\right)/2}}{\Vert m\Vert }_{{ℍ}^t\left(\Gamma \right)}\mathrm{.}\end{array}$

,

5. Practical Implementation and Numerical Results

We want to present in this section some results of the former method. We will describe first the process of obtaining the computational results. We distinguish

Figure 2. Relative errors for several RHS in function of levels: (a) double layer, (b) single layer.

separate tests for the single layer operator, the double layer operator and the electrostatic reaction potential. The exact solutions which have zero Laplacians are specified by the user in the case of the double/single layer potentials. The corresponding right hand sides are then obtained by applying the Lapacian operator to the exact solutions as illustrated in Figure 2(a) and Figure 2(b). Transforming Laplace problems into integral equations using double/single layers is found in standard textbooks of integral equations [4] [45] . The assembly of the basis functions consists in constructing piecewise constant wavelets on the unit interval. After taking the tensor product, one obtains basis functions on the unit square which are mapped by the NURBS functions onto the patches that are embedded in the 3D-space. The assembly of the matrix entries uses hierarchical integrations which are described in our former paper [31] . That procedure can be achieved on arbitrary geometries. The linear system obtained from the double layer operator does not require any preconditioner as a direct GMRes linear solver suffices to solve the system. For the single layer operator however, we use a domain decomposition preconditioner to reduce the number of iterations. The error computation consists in comparing the chosen exact solution with the outputs which are acquired by applying the BEM-simulation. The reason for displaying separate results pertaining to the double layer potential and the single layer potential is that if they are solved individually, analytical results for arbitrary boundary $\Gamma$ exist for comparisons. As for the computation of the electrostatic reaction potential, we make use of the Born ion. That is, the only well-known analytical result for the reaction potential $P^{REAC}$ is a simple geometry composed of a single atom. For that case, the process consists in fixing some configurations by setting the values of the radius, the electric charges as well as the dielectric parameter. The actual process consists therefore in applying the PCM simulation on a resolution which is specified by the wavelet level. The electrostatic potential is obtained by traversing the patches and by computing some integrals related to the apparent surface charge as specified in (18). For the numerical computation of those integrals, a standard Gaussian quadrature for smooth integrand suffices because the nuclei ${\left\{x_k\right\}}_k$ are not located on the NURBS patches ${\left\{{\Gamma }_p\right\}}_p$ . Different parameters are varied in order to validate the accuracy of the PCM-simulation.

5.1. Molecular Models

We employ different sorts of quantum models including molecules which are acquired from PDB/PQR files [46] as well as water clusters which are obtained from a former MD (Molecular Dynamics) simulation. When the MD-iteration attains its equilibrium state where the total energy becomes stable, a water cluster is obtained by extracting the water molecules which are contained in some given large sphere whose radius controls the final size of the water cluster. The formerly proposed method assumes that the sizes of the patches are proportional. Efforts were done to obtain patches whose shapes related to the normal derivatives are as good as possible. In order to measure the quality of the patches, we examined the proportionality of the area, length of the curved edges and the norm of the normal vectors. That is important because we use the same wavelet threshold for all patches. We avoid a patch where some curved boundary edges are very long while others are very short. For each patch ${\Gamma }_q$ , the quality gauge for the area is

$max\left\{\frac{\text{area}\left({\Gamma }_q\right)}{\text{averagearea}}\mathrm{,}\frac{\text{averagearea}}{\text{area}\left({\Gamma }_q\right)}\right\}\mathrm{.}$ (131)

The average of the above value is collected on Table 1 which shows that although the areas are different, their ratios with the average area are not excessive. Similar quality measurements as (131) have been utilized for the length of curved edges separating the patches [47] as well as the normal vectors which are computed at some tensor product samples. It is observed from Table 1 that all involved entities in the patches are practically proportional.

5.2. Double Layer Potential

For the double layer potential, one can solve the linear system iteratively without recourse to any preconditioner because the system is well conditioned. A GM-Res method serves as a solver of the system which is non-symmetric. We

Table 1. Quality of the NURBS patches for several molecules.

Table 2. BEM accuracy for the double layer and the single layer in function of the levels.

shall examine first the error reductions in term of the multiscale level. The results are collected in Table 2 which contains both the absolute errors and the relative errors. The absolute errors are obtained by comparing with the exact solution at some fixed samples in the interior $\Omega$ . A division by the largest value of the exact solution provides the relative error.

In the first part of Table 2, we collect the convergence of the BEM simulation in function of the level ranging from 1 to 4. Those results have been obtained from a water cluster molecule containing 386 patches. The ratio between the nonzero counts of the compressed and uncompressed matrices are also exhibited. We examine now the general linear characteristic of the relative errors. We display in Figure 2(a) the error evolution where we consider five levels applied to a propane molecule using several right hand sides. We consider five exact solutions which are respectively ${\mathcal{U}}_1\left(x\right)=0.2x^2-0.15y^2-0.05z^2$ and ${\mathcal{U}}_2\left(x\right)=$ $exp\left(0.5x\right)cos\left(0.5y\right)$ , ${\mathcal{U}}_3\left(x\right)={\left[{\left(x-5\right)}^2+{\left(y-5\right)}^2+{\left(z-5\right)}^2\right]}^{-1/2}$ , ${\mathcal{U}}_4\left(x\right)=0.1x$ and ${\mathcal{U}}_5\left(x\right)=-\exp \left(0.01y\right)\sin \left(0.01x\right)+\exp \left(0.2z\right)\cos \left(0.2x\right)-\exp \left(0.1z\right)\sin \left(0.1y\right)$ $+\exp \left(0.05y\right)\cos \left(0.05x\right)$ that all have vanishing Laplacian. The right hand side $g\left(x\right)$ is the restriction of the function $\mathcal{U}$ on the boundary $\Gamma$ . The propane molecule admits 75 patches. The error reduction is lightly affected by the used right hand side but in general the errors reduce satisfactorily in function of the wavelet levels. In fact, they decrease linearly in logarithmic scale in function of the BEM levels.

5.3. Single Layer Potential

The single layer yields a linear system which is badly conditioned [48] . We use an additive domain decomposition preconditioner to remedy that bad conditioning. We briefly summarize the procedure whose comprehensive description is detailed in our former work [38] . In term of geometric structure, the overlapping domain decomposition is as follows

$\Gamma ={\displaystyle \underset{p=1}{\overset{M}{\cup }}}{\Omega }_p\text{such}\text{that}{\Omega }_p={\displaystyle \underset{i\in {\mathcal{N}}_p}{\cup }}{\Gamma }_i$ (132)

where

${\Omega }_p\cap {\Omega }_q\text{is}\text{not}\text{necessarily}\text{empty}\text{for}p\ne q\mathrm{.}$ (133)

In term of linear spaces, this leads to the decomposition

$\mathcal{V}\left(\Gamma \right)=\mathcal{V}\left({\Omega }_1\right)+\cdots +\mathcal{V}\left({\Omega }_M\right)\text{where}\mathcal{V}\left({\Omega }_p\right)\mathrm{<mo>\; :\; </mo>}=\underset{i\in {\mathcal{N}}_p}{{\displaystyle \oplus }}{\mathcal{V}}_L\left({\Gamma }_i\right)\mathrm{.}$ (134)

Denote the orthogonal projection onto $\mathcal{V}\left({\Omega }_p\right)$ with respect to the bilinear form $\langle V\cdot \mathrm{,}\cdot \rangle$ related to the single layer operator $V$ by

$P_{{\Omega }_p}\mathrm{<mo>\; :\; </mo>}\mathcal{V}\left(\Gamma \right)\to \mathcal{V}\left({\Omega }_p\right)\mathrm{,}\langle V{P}_{{\Omega }_p}v\mathrm{,}\varphi \rangle =\langle Vv\mathrm{,}\varphi \rangle ,\forall \varphi \in \mathcal{V}\left({\Omega }_p\right)\mathrm{.}$ (135)

The ASM operator is defined by

$P_{\Gamma }^{ASM}\mathrm{<mo>\; :\; </mo>}\mathcal{V}\left(\Gamma \right)\to \mathcal{V}\left(\Gamma \right)\mathrm{,}{P}_{\Gamma }^{ASM}\mathrm{<mo>\; :\; </mo>}=P_{{\Omega }_1}+\cdots +P_{{\Omega }_M}\mathrm{.}$ (136)

The ASM procedure consists in solving a local problem which is projected by $P_{{\Omega }_p}$ onto the subdomain ${\Omega }_p$ . The results of the number of iterations are summarized in Figure 3(a) where it is observed that the ASM method needs significantly less iterations than the direct method. A larger view for the iterations less than 50 is exhibited in Figure 3(b). Like the double layer, the same test for the compression is detailed in the second part of Table 2 in the case of the single layer potential. Further, the decrease of the BEM error in function of the levels for different RHS is presented in Figure 2(b) where a propane molecule is used as in the double layer case.

5.4. Apparent Surface Charge

We want now to examine the values of electrostatic reaction potential $P^{REAC}$ which are computed with the PCM-wavelet technique. It is beyond the scope of

Figure 3. (a) Preconditioning the single layer potential for the entire iterations. (b) Zoom for preconditioning the single layer potential for iterations 0 - 60.

this article to provide a detailed comparison with other softwares as all results here are only computed for analytical comparison where we do not use any physical units. The explicit values for the reaction potential are unknown except for some simple geometries. We consider the Born ion that consists of a single sphere of radius $R>0$ which is centered at the origin and to which an electric charge q is assigned while the dielectric coefficient $\epsilon$ of the solute is adjusted. The analytical expression [11] for the reaction potential reads

$P^{REAC}=-\left(q^2/8\text{π}R\right)\left(1-1/\epsilon \right)$ . All those parameters are unitless in the sense that no physical units are used to measure them. Our next test consists in examining the effect of the dielectric coefficient $\epsilon \ge 1$ . The corresponding result is depicted in Figure 4(a) in which we consider three configurations corresponding to $\left(R,q\right)=\left(1.2,2.0\right)$ , $\left(R,q\right)=\left(2.2,2.0\right)$ and $\left(R,q\right)=\left(1.5,1.0\right)$ respectively. We execute the wavelet-PCM simulation on multiscale level 3 while the dielectric coefficient varies from $\epsilon =1$ to $\epsilon =200$ . We observe for all three configurations that the exact reaction potential and the computed PCM results align well in Figure 4(a) where a zoom is also provided for $\epsilon \in \left[\mathrm{1,20}\right]$ in which the reaction potential drops down very quickly. As a next test, we investigate the variation of

Figure 4. Unitless comparison for the Born ion: (a) dielectric coefficients, (b) electric charge, (c) radius.

the electrostatic potential in function of the electric charge q. The corresponding outcome is depicted in Figure 4(b) where we consider four configurations corresponding to $R=4$ , $R=3$ , $R=2$ , $R=1$ while the dielectric coefficient is constantly $\epsilon =200$ . The electric charge is allowed to vary in the range $q\in \left[\mathrm{1.5,4.25}\right]$ . One observes in Figure 4(b) that the computed reaction potential varies quadratically in function of the electric charge. In fact, a parabolic dependence on the charge q is exhibited. The purpose of the next test displayed in Figure 4(c) is to highlight the agreement between the computed PCM result and the exact solution when the ion radius varies from $R=0.5$ to $R=3.5$ . One considers three configurations where the electric charge is respectively $q=1.5$ , $q=2.0$ and $q=2.5$ while the dielectric coefficient is consistently $\epsilon =12.5$ . For all three configurations, the electrostatic reaction potential is inversely proportional with regard to the ion radius. For all cases, the results computed by means of the PCM-wavelet align well with the theoretical expectations. As a final test, we examine the accuracy of the electrostatic potential in function of the BEM-levels which vary from 1 to 6. We consider two configurations where the radii and electric charges are $\left(R,q\right)=\left(1.2,2.0\right)$ and $\left(R,q\right)=\left(2.3,2.5\right)$ respectively. The results of the BEM accuracy are summarized in Table 3 where we consider four dielectric values $\epsilon =15$ , $\epsilon =5$ , $\epsilon =100$ and $\epsilon =45$ . The exact and the computed reaction potentials which are exhibited there imply that the BEM-approximation becomes satisfactorily precise as the multiscale levels increase.

Table 3. Computed PCM-reaction potential in function of the multi-scale levels for the Born ion.

6. Conclusion

We consider the transmission problem occurring in the interaction between the solute and the solvent media. Our method is based on an integral equation formulation which is solved by means of the wavelet Boundary Element Method. The entire molecular surface is supposed to be structured in four-sided NURBS patches onto which wavelets constructed on the unit square are embedded. For the Born ion, the exact reaction potential is known explicitly without using physical units. Our results align well with the exact solutions when we vary various parameters including the electric charge, the radius and the dielectric coefficient. The current approach is limited to the linear case where the coefficients in the transmission problem are constant parameters. Prospective future works could comprise the nonlinear Poisson-Boltzmann admitting nonsmooth Dirac load function. Further, non-constant matrix-valued coefficients which are difficult to treat by BEM might be examined. These challenging works need to be approached and presented step by step.

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Cancès, E., Mennucci, B. and Tomasi, J. (1998) A New Integral Equation Formalism for the Polarizable Continuum Model: Theoretical Background and Applications to Isotropic and Anisotropic Dielectrics. Journal of Chemical Physics, 107, 3032-3041. https://doi.org/10.1063/1.474659
[2]	Cancès, E. (1998) Simulation Moléculaire et Effets d’Environnement. Une Perspective Mathématique et Numérique. Ph.D. Thesis, école Nationale des Ponts et Chau- ssées, Paris.
[3]	Maischak, M., Stephan, E. and Tran, T. (2004) A Multiplicative Schwarz Algorithm for the Galerkin Boundary Element Approximation of the Weakly Singular Integral Operator in Three Dimensions. International Journal of Pure and Applied Mathematics, 12, 1-21.
[4]	McLean, W. (2000) Strongly Elliptic Systems and Boundary Integral Equations. Cambridge University Press, Cambridge.
[5]	Pechstein, C. (2013) Special Lecture on Boundary Element Methods, Lecture Notes. Institute of Computational Mathematics, Johannes Kepler University Linz, Linz.
[6]	Qiu, T. and Sayas, F. (2016) The Costabel-Stephan System of Boundary Integral Equations in the Time Domain. Mathematics of Computation, 85, 2341-2364. https://doi.org/10.1090/mcom3053
[7]	Piegl, L. and Tiller, W. (1995) The NURBS Book. Springer, Berlin. https://doi.org/10.1007/978-3-642-97385-7
[8]	De Boor, C. and Fix, G. (1973) Spline Approximation by Quasi-Interpolants. Journal of Approximation Theory, 8, 19-45. https://doi.org/10.1016/0021-9045(73)90029-4
[9]	Hoschek, J. and Lasser, D. (1996) Fundamentals of Computer Aided Geometric Design. Taylor & Francis, Abingdon.
[10]	Brunnett, G. (1995) Geometric Design with Trimmed Surfaces. Computer Supply, 10, 101-115. https://doi.org/10.1007/978-3-7091-7584-2_7
[11]	Bajaj, C., Xu, G. and Zhang, Q. (2009) A Fast Variational Method for the Construction of Adaptive Resolution C2 Smooth Molecular Surfaces. Computer Methods in Applied Mechanics and Engineering, 198, 1684-1690. https://doi.org/10.1016/j.cma.2008.12.042
[12]	Ciarlet, P. (1987) The Finite Element Method for Elliptic Problems. Society for Industrial and Applied Mathematics, Philadelphia.
[13]	Diedrich, C., Dijkstra, D., Hamaekers, J., Henniger, B. and Randrianarivony, M. (2015) A Finite Element Study on the Effect of Curvature on the Reinforcement of Matrices by Randomly Distributed and Curved Nanotubes. Journal of Computational and Theoretical Nanoscience, 12, 2108-2116. https://doi.org/10.1166/jctn.2015.3995
[14]	Randrianarivony, M. (2004) Anisotropic Finite Elements for the Stokes Problem: A- Posteriori Error Estimator and Adaptive Mesh. Journal of Computational and Applied Mathematics, 169, 255-275. https://doi.org/10.1016/j.cam.2003.12.025
[15]	Baker, N., Sept, D., Holst, M. and McCammon, J. (2001) The Adaptive Multilevel Finite Element Solution of the Poisson-Boltzmann Equation on Massively Parallel Computers. IBM Journal of Research and Development, 45, 427-438. https://doi.org/10.1147/rd.453.0427
[16]	Baker, N., Sept, D., Joseph, S., Holst, M. and McCammon, J. (2001) Electrostatics of Nanosystems: Application to Microtubes and the Ribosome. Proceedings of the National Academy of Sciences of the United States of America, 98, 10037-10041. https://doi.org/10.1073/pnas.181342398
[17]	Bank, R. and Holst, M. (2003) A New Paradigm for Parallel Adaptive Meshing Algorithms. SIAM Review, 45, 291-323. https://doi.org/10.1137/S003614450342061
[18]	Schwarz, L. (1966) Théorie des Distributions. Hermann, Paris.
[19]	Holst, M. (2001) Adaptive Numerical Treatment of Elliptic System on Manifolds. Advances in Computational Mathematics, 15, 139-191. https://doi.org/10.1023/A:1014246117321
[20]	Holst, M. and Saied, F. (1993) Multigrid Solution of the Poisson-Boltzmann Equation. Journal of Computational Chemistry, 14, 105-113. https://doi.org/10.1002/jcc.540140114
[21]	Dahmen, W. and Schneider, R. (1999) Wavelets on Manifolds I: Construction and Domain Decomposition. SIAM Journal on Mathematical Analysis, 31, 184-230. https://doi.org/10.1137/S0036141098333451
[22]	Harbrecht, H. and Randrianarivony, M. (2010) From Computer Aided Design to Wavelet BEM. Computing & Visualization in Science, 13, 69-82. https://doi.org/10.1007/s00791-009-0129-1
[23]	Hsiao, G., Steinbach, O. and Wendland, W. (2000) Domain Decomposition Methods via Boundary Integral Equations. Journal of Computational and Applied Mathematics, 125, 521-537. https://doi.org/10.1016/S0377-0427(00)00488-X
[24]	Heuer, N., Stephan, E. and Tran, T. (1998) Multilevel Additive Schwarz Method for the H-P Version of the Galerkin Boundary Element Method. Mathematics of Computation, 67, 501-518. https://doi.org/10.1090/S0025-5718-98-00926-0
[25]	Cohen, A., Daubechies, I. and Feauveau, J. (1992) Biorthogonal Bases of Compactly Supported Wavelets. Communications on Pure and Applied Mathematics, 45, 485- 560. https://doi.org/10.1002/cpa.3160450502
[26]	Lyche, T., Morken, K. and Quak, E. (2001) Theory and Algorithms for Non-Uni- form Spline Wavelets, Multivariate Approximation and Applications. Cambridge University Press, Cambridge, 152-187.
[27]	Beylkin, G. (1992) On the Representation of Operators in Bases of Compactly Supported Wavelets. SIAM Journal on Numerical Analysis, 29, 1716-1740. https://doi.org/10.1137/0729097
[28]	Lage, C. and Schwab, C. (1999) Wavelet Galerkin Algorithms for Boundary Integral Equations. SIAM Journal on Scientific Computing, 20, 2195-2222. https://doi.org/10.1137/S1064827597329989
[29]	Petersdorff, T. and Schwab, C. (1996) Wavelet Approximations for First Kind Boundary Integral Equations on Polygons. Numerische Mathematik, 74, 479-519. https://doi.org/10.1007/s002110050226
[30]	Dolinsky, T., Nielsen, J., McCammon, J. and Baker, N. (2004) PDB2PQR: An Automated Pipeline for the Setup, Execution, and Analysis of Poisson-Boltzmann Electrostatics Calculations. Nucleic Acids Research, 32, 665-667. https://doi.org/10.1093/nar/gkh381
[31]	Randrianarivony, M. (2016) Multiwavelet Boundary Element Method for Cavities Admitting Many NURBS Patches. Modeling and Numerical Simulation, 6, 69-93. https://doi.org/10.4236/mnsms.2016.64007
[32]	Randrianarivony, M. (2006) Geometric Processing of CAD Data and Meshes as Input of Integral Equation Solvers. Ph.D. Thesis, Technical University of Chemnitz, Chemnitz.
[33]	Randrianarivony, M. and Brunnett, G. (2008) Preparation of CAD and Molecular Surfaces for Meshfree Solvers. Lecture Notes in Computational Science and Engineering Description, 65, 231-245. https://doi.org/10.1007/978-3-540-79994-8_14
[34]	Randrianarivony, M. (2009) On Global Continuity of Coons Mappings in Patching CAD Surfaces. Computer-Aided Design, 41, 782-791. https://doi.org/10.1016/j.cad.2009.04.012
[35]	Randrianarivony, M. (2008) Harmonic Variation of Edge Size in Meshing CAD Geometries from IGES Format. Lecture Notes in Computer Science, 5102, 56-65. https://doi.org/10.1007/978-3-540-69387-1_7
[36]	Randrianarivony, M. and Brunnett, G. (2008) Molecular Surface Decomposition Using Geometric Techniques. Proc. Conf. Bildverarbeitung für die Medizine, Berlin, 6-8 April 2008, 197-201.
[37]	Weijo, V., Randrianarivony, M., Harbrecht, H. and Frediani, L. (2010) Wavelet For- mulation of the Polarizable Continuum Model. Journal of Computational Chemistry, 31, 1469-1477.
[38]	Randrianarivony, M. (2016) Domain Decomposition for Wavelet Single Layer on Geometries with Patches. Applied Mathematics, 7, 1798-1823. https://doi.org/10.4236/am.2016.715151
[39]	Oswald, P. (1998) Multilevel Norms for . Computing, 61, 235-255. https://doi.org/10.1007/BF02684352
[40]	Micchelli, C., Xu, Y. and Zhao, Y. (1997) Wavelet Galerkin Methods for Second- Kind Integral Equations. Journal of Computational and Applied Mathematics, 86, 251-270. https://doi.org/10.1016/S0377-0427(97)00160-X
[41]	Beylkin, G., Coifman, R. and Rokhlin, V. (1991) Fast Wavelet Transforms and Numerical Algorithms I. Communications on Pure and Applied Mathematics, 44, 141- 183. https://doi.org/10.1002/cpa.3160440202
[42]	Graham, I., Hackbusch, W. and Sauter, S. (2000) Discrete Boundary Element Methods for General Meshes in 3D. Numerische Mathematik, 86, 103-137. https://doi.org/10.1007/PL00005399
[43]	Rathsfeld, A. (1995) A Wavelet Algorithm for the Solution of the Double Layer Potential Equation over Polygonal Boundaries. Journal of Integral Equations and Applications, 7, 47-98. https://doi.org/10.1216/jiea/1181075850
[44]	Nedelec, J. and Planchard, J. (1973) Une Méthode Variationelle d’éléments Finis pour la Résolution Numérique d’un Problème Extérieur dans . Revue Franaise d’Automatique, Informatique, Recherche Opérationelle, Mathématique, 3, 105-129.
[45]	Hsiao, G. and Wendland, W. (2008) Boundary Integral Equations. Springer, Berlin, 164. https://doi.org/10.1007/978-3-540-68545-6
[46]	Dolinsky, T., Czodrowski, P., Li, H., Nielsen, J., Jensen, J., Klebe, G. and Baker, N. (2007) PDB2PQR: Expanding and Upgrading Automated Preparation of Biomolecular Structures for Molecular Simulations. Nucleic Acids Research, 35, 522-525. https://doi.org/10.1093/nar/gkm276
[47]	Dahlke, S. and Weimar, M. (2015) Besov Regularity for Operator Equations on Patchwise Smooth Manifolds. Foundations of Computational Mathematics, 15, 1533-1569. https://doi.org/10.1007/s10208-015-9273-9
[48]	Gemmrich, S., Gopalakrishnan, J. and Nigam, N. (2012) Convergence Analysis of a Multigrid for the Acoustic Single Layer Equation. Applied Numerical Mathematics, 62, 767-786. https://doi.org/10.1016/j.apnum.2012.02.003