Multilevel Monte Carlo on a high-dimensional parameter space for transmission problems with geometric uncertainties
MMultilevel Monte Carlo on a high-dimensional parameterspace for transmission problems with geometricuncertainties
Laura Scarabosio ∗ June 27, 2017
Abstract
In the framework of uncertainty quantification, we consider a quantity of interest whichdepends non-smoothly on the high-dimensional parameter representing the uncertainty.We show that, in this situation, the multilevel Monte Carlo algorithm is a valid optionto compute moments of the quantity of interest (here we focus on the expectation), asit allows to bypass the precise location of discontinuities in the parameter space. Weillustrate how such lack of smoothness occurs for the point evaluation of the solution to a(Helmholtz) transmission problem with uncertain interface, if the point can be crossed bythe interface for some realizations. For this case, we provide a space regularity analysisfor the solution, in order to state converge results in the L ∞ -norm for the finite elementdiscretization. The latter are then used to determine the optimal distribution of samplesamong the Monte Carlo levels. Particular emphasis is given on the robustness of ourestimates with respect to the dimension of the parameter space. Keywords: multilevel Monte Carlo, shape uncertainty, interface problem, L ∞ -estimates, uncer-tainty quantification. In many engineering applications, the behavior of a physical system depends on a parametervector y belonging to a high-dimensional parameter space P J ⊆ R J with J ∈ N large. The vec-tor y ∈ P J may represent, for instance, random variations in material or geometrical propertiesof the physical system. Equipping P J with a σ -algebra A J and a probability measure µ J , weobtain the probability space ( P J , A J , µ J ). In such cases, it is of interest to compute statistics,with respect to the parameter, of a quantity q ( y ; u ( y )) (quantity of interest, Q.o.I. for short)depending on the solution u to a partial differential equation (PDE):Find u s.t. : u ( y ) ∈ X for every y ∈ P J , D ( y ; u ( y )) = 0 for every y ∈ P J . (1.1)In (1.1), X denotes a separable Banach space. For every y ∈ P J and every J ∈ N , q ( y ; · ) : X → Y , that is, every realization of q belongs to a separable Hilbert space Y . For instance,if q is the solution u itself, then Y = X , if q is some linear output functional, then Y = R or Y = C . ∗ Technical University of Munich, Germany, Department of Mathematics, Institute for Nu-merical Mathematics,
Boltzmannstrasse 3, 85748 Garching b. M¨unchen.
Email : [email protected] a r X i v : . [ m a t h . NA ] J un ntroducing the quantity Q : P J → Y such that Q ( y ) := q ( y ; u ( y )) for every y ∈ P J , thepresent work focuses on the case when Q is non-smooth , with respect to the high-dimensionalparameter y , across a submanifold P Γ J ⊂ P J which is not easy to track. Here, by non-smoothwe mean ‘not analytic in y ’, and in our treatment we allow Q to have jumps across P Γ J . Moreprecisely, in this paper (1.1) is a (acoustic) transmission problem, where the shape of thescatterer is subject to random variations modeled by the high-dimensional parameter y , andthe Q.o.I. is the point evaluation of the solution in locations that, depending on the realization,may be either inside or outside the scatterer. We focus on the computation of the mean E µ J [ Q ] := (cid:90) P J Q ( y ) d µ J ( y ) , (1.2)and aim at numerical methods which are robust with respect to the dimension J of the param-eter space, that is, whose convergence rates do not deteriorate for large J , possibly tending toinfinity. Related work.
We first review the literature on the computation of moments of a Q.o.I.,and then, in view of our application to a transmission problem with random interface, theliterature in shape uncertainty quantification.If the randomness in the system consists of deviations from a deterministic quantity that aresmall enough, it is possible to apply a perturbation approach [17], and approximate momentsof the Q.o.I. exploiting its Taylor expansion centered at the deterministic quantity. Otherwise,we have to compute (1.2) directly (or analogous expression for higher order moments), whichmeans employing quadrature formulas on the parameter space. In this work we focus on thissecond option. Quadrature rules on a (high-dimensional) parameter space can be classifiedin two main cathegories: Monte Carlo-like rules and deterministic rules. As with quadraturerules for functions of one real variable, there is a compromise between speed of convergencewith respect to the number of function evaluations and smoothness required on the integrand.The Monte Carlo approach to compute (1.2), consisting of random sampling [9], convergesalmost surely to the exact mean provided the Q.o.I. is Lebesgue integrable with respect to theparameter. This is ensured by the strong law of large numbers [9, Sect. 2]. If the Q.o.I. has alsofinite variance, then the Monte Carlo quadrature converges with rate M − / , where M is thenumber of samples [9, Thm. 2.1]. The high computational effort due to the slow convergencerate can be reduced using the multilevel Monte Carlo (MLMC) method [26, 27, 35, 36] or othervariance reduction techniques [28]. To converge, MLMC requires square integrability of theQ.o.I., and details are provided in Section 4 of this paper. Deterministic quadrature rulescomprise quasi-Monte Carlo (QMC) methods and spectral methods. We refer to [20] for acomprehensive treatment of QMC. It is possible to construct QMC sequences of quadraturepoints such that the speed of convergence is M − ε , for any (cid:15) > M the number ofquadrature points) [21, Prop. 2.18, Thm. 3.20 and Sect. 3.4], under the assumption that theintegrand has continuous first order mixed derivatives. If the integrand has higher regularity,then higher order QMC quadrature rules can be constructed, with convergence rates that arerobust with respect to the dimension of the parameter space [19]. Spectral methods can bedivided in stochastic Galerkin [2, 52, 57] and stochastic collocation [1, 45] approaches. Theyprovide high order convergence rates if the Q.o.I. admits an analytic extension to the complexplane: for finite-dimensional parameter spaces, the rate is exponential with respect to number ofevaluation points, but it depends on the dimension and deteriorates as the latter increases [1,5];the dimension-independent convergence rate, which still holds in infinite-dimensional parameterspaces, is algebraic, and it depends only on the ‘sparsity class of the unknown’ [13,51,52]. If theQ.o.I. is not globally smooth with respect to the parameter, but it is piecewise smooth, thenone possibility is to employ discontinuity detection methods (as, for instance, the one suggestedin [59]) to detect the surfaces of non-smoothness, and then apply a high order quadrature rule2eparately on each subdomain on which the Q.o.I. is smooth. However, this approach is notapplicable for complicated surfaces of discontinuity. This issue is discussed in more details insubsection 3.4 of this work, which then motivates why MLMC is a valid option when non-smoothness occurs across manifolds that are not easy to track.In the model problem that we consider, the randomness stems from uncertain variations ofthe scatterer boundary. Several approaches are possible to tackle shape uncertainty quantifi-cation: perturbation techniques [12, 31, 33] (analogous to [17] using shape calculus to constructthe Taylor expansions), level set methods [46, 47], the fictitious domain approach [10] andthe mapping technique [54, 58]. Recently, a new approach has been suggested in [34] in theframework of a Helmholtz scattering problem, where a boundary integral formulation is used toreconstruct the expansion of the solution in spherical or cylindrical harmonics; however, explicitformulas for the coefficients seem to be available, for the moment, only when the parameterspace is low-dimensional. In our paper, we adopt the mapping technique, because it allows todeal with not small perturbations, it provides a natural way of resolving the interface for thespatial discretization [32, Sect. 5], and, transforming a PDE on a random domain to a PDEon a deterministic domain with stochastic coefficients, it simplifies both theoretical analysisand practical implementation. The regularity of the solution to a PDE with respect to thehigh-dimensional parameter describing the shape variations has been studied in [11, 15, 32, 37]and [39]. The authors of these papers prove holomorphic dependence, with respect to the high-dimensional parameter, of the solution on the nominal, deterministic domain introduced by thedomain mapping. The work [11] deals with an elliptic boundary value problem, [32] tackles alsoan elliptic interface problem, and [37] treats the same Helmholtz transmission problem as theone addressed in the present paper (and considers also some linear output functionals). Thepaper [15] provides, in the framework of the stationary Navier-Stokes equations, a unified math-ematical treatment of the mapping method, independent of the domain parametrization, andintroduces the concept of ‘shape holomorphy’. The techniques presented in [15] have been ap-plied, in [39], to the Maxwell equations in frequency domain. However, the smooth dependenceon the parameter breaks down for point evaluations of the solution to an interface problem onthe physical domain, where the interface changes for every realization [49, Ch. 8]. This is thecase treated in this paper. For an application of the mapping technique to the inverse problemsetting, we refer to [23] and [38], where the inverse problem in electrical impedance tomographyis considered. In particular, in [23] the authors prove the regularity of the posterior measurewith respect to the high-dimensional parameter associated to the shape variations. Scope and outline of the paper.
One goal of this paper is to highlight the presence andimpact of the non-smooth dependence on the stochastic parameter in the case of an importantclass of transmission problems with stochastic interface, namely Helmholtz transmission prob-lems. We prove, and confirm by numerical experiments, that MLMC offers a robust treatmentfor such class of problems, allowing to bypass the precise location of discontinuities in the pa-rameter space. The methodology used clearly conveys that MLMC is a viable approach alsofor other problems lacking smoothness with respect to the stochastic parameter. The secondgoal is to provide a full numerical analysis for point evaluation in (Helmholtz) transmissionproblems with geometric uncertainties, including the regularity of the solution with respectto the parameter and to the spatial coordinate, and their implications in the convergence ofMLMC.The paper is organized as follows. In Section 2 we introduce our model transmission problem.Sections 3 and 5 are the core of this paper. In Section 3, we show that the point value of thesolution in locations that might be crossed by the random interface is a Q.o.I. which does notdepend smoothly on the parameter describing the shape variations. The main contributionthere is Proposition 3.2, where we state the regularity of the Q.o.I. with respect to the high-dimensional parameter. In the same section, we discuss possible ways to handle the non-smooth3 r ( y ; ϕ ) R out D in ( y ) D out,R out ( y )Γ( y ) ∂D R out u i Figure 2.1: Geometry of our model problem. parameter dependence, and provide our motivation for choosing the MLMC method. Thelatter is reviewed in Section 4, with focus on Q.o.I.s depending on the solution of a partialdifferential equation. In Section 5, we first analize the space regularity of the solution tothe model transmission problem, and then state convergence results for MLMC when using afinite element discretization. Finally, in Section 6, we show numerical experiments matchingthe theoretical predictions. For ease of presentation, technical details on the space regularityof the solution, used in the proofs of Proposition 3.2 and Theorem 5.1, have been moved toAppendices A and B, since they consist in the adaptation of already existing results to ourHelmholtz transmission problem.
As a model problem, we address the Helmholtz transmission problem in R , describing thescattering of an incoming wave u i from a penetrable object whose shape is subject to randomvariations. We formally define Γ( y ), y ∈ P J , to be the boundary of the scatterer, and denote by D in ( y ) the domain enclosed inside Γ( y ). We consider a circle of fixed radius R out containing allrealizations of the scatterer in its interior, and indicate and by D out,R out ( y ) the part of the outer,unbounded domain contained inside this circle. Finally, D R out := D in ( y ) ∪ Γ( y ) ∪ D out,R out ( y ).Geometry and notation are clarified in Fig. 2.1.The transmission problem for the Helmholtz equation reads: − ∇ · ( α (Γ( y ); x ) ∇ u ) − κ (Γ( y ); x ) u = 0 in D in ( y ) ∪ D out,R out ( y ) , (cid:74) u (cid:75) Γ( y ) = 0 , (cid:74) α (Γ( y ); x ) ∇ u · n (cid:75) Γ( y ) = 0 ,∂∂ n out ( u − u i ) = DtN( u ) − DtN( u i ) on ∂D R out , for every y ∈ P J , (2.1a)(2.1b)(2.1c)where we consider real-valued, piecewise-constant coefficients α (Γ( y ); x ) = (cid:40) x ∈ D out,R out ( y ) ,α if x ∈ D in ( y ) , κ (Γ( y ); x ) = (cid:40) κ if x ∈ D out,R out ( y ) ,α κ if x ∈ D in ( y ) . (2.2)We assume u i to be a plane wave, that is u i ( x ) = e jκ d · x , where d is a direction vector with (cid:107) d (cid:107) = 1 and j = √−
1. The unknown u = u ( y ; x ) represents the total field, whereas κ , κ > α is a positive4oefficient. In equation (2.1b), (cid:74) · (cid:75) Γ( y ) denotes the jump across the random interface Γ( y ).Equation (2.1c) is the exact boundary condition on the disc of radius R out , and correspondsto the radiation condition in free space (Sommerfeld radiation condition). Such boundarycondition is stated in terms of the Dirichlet-to-Neumann map (DtN) on the scattered wave,see [44, Sect. 6.2.3] for its definition.We work in the large wavelength regime , assuming the wavelength to be large enough com-pared to the size of the scatterer (see Assumption 3.1). Mathematically, this means that weaddress the case when the bilinear form associated to (2.1) is coercive.We consider here an explicit description for the interface. We assume the scatterer to bestar-shaped with respect to the origin, and set, in polar coordinates, Γ( y ) := { ( r, ϕ ) : r = r ( y ; ϕ ) , ϕ ∈ [0 , π ) } , where r is a stochastic, angle-dependent radius (see Fig.2.1). We expressthe latter as: r ( y ; ϕ ) = r ( ϕ ) + J (cid:88) j =1 β j y j ψ j ( ϕ ) , ψ j ( ϕ ) = (cid:40) sin( j +12 ϕ ) for j odd , cos( j ϕ ) for j even , (2.3)for every ϕ ∈ [0 , π ), J ∈ N and y = ( y , . . . , y J ) ∈ P J . The quantity r ∈ C k,βper ([0 , π )), forsome k ≥ < β <
1, is an approximation to the mean radius. The real parameters { y j } Jj =1 are the images of independent, identically distributed (i.i.d.) uniform random variables Y j ∼ U ([ − , ≤ j ≤ J . Thus, P J = [ − , J and µ J is the product measure µ J = (cid:0) (cid:1) J . Forevery J ∈ N , the radius (2.3) is a well-defined random variable on the closure of the subspacespan (cid:110) , ( ψ j ) Jj =1 (cid:111) in the C k,βper ([0 , π ))-norm. For the expansion (2.3), we require: Assumption 2.1.
The sequence ( β j ) j ≥ in (2.3) has a monotonic majorant in (cid:96) p ( N ) with < p < . Furthermore, (cid:80) j ≥ | β j | ≤ r − , with r − = inf ϕ ∈ [0 , π ) r ( ϕ ) > . The bounds on (cid:80) j ≥ | β j | and r − ensure positivity and boundedness of the radius for everyrealization. The condition on the decay of ( β j ) j ≥ , instead, is a regularity assumption (withrespect to ϕ ) on the radius: the smaller the p , the smoother the radius [49, Lemma 2.1.6]. Inparticular, p < ensures that every realization of r ( y ; · ) ∈ C ,β (0 , π ) for some β ∈ (0 , J and y -independent bound (cf. proof of Lemma 2.1.6 in [49]).Our Q.o.I. is Q ( y ) = u ( y ) := { u ( y ; x i ) } N − i =0 ⊂ C N , y ∈ P J , the value of the solution u to (2.1) at fixed points { x , . . . , x N − } ⊂ R . In particular, we are interested in the case thatthese evaluation points are close to the interface, so that they may lie on different sides of Γ( y )for different realizations of y . The aim of this section is to highlight the non-smooth dependence of u ( y ) = { u ( y ; x i ) } N − i =0 onthe high-dimensional parameter y ∈ P J . To better explain the non-smooth behavior, we firstconsider, in subsection 3.1, a one-dimensional transmission problem. Then, in subsection 3.2,we move to the model problem introduced in the previous section. Due to the failure of highorder quadrature methods to compute E µ J [ u ], illustrated in subsection 3.3, in subsection 3.4we discuss how this issue can be overcome, and motivate why we opt for MLMC. We consider the one-dimensional problem 5 .25 0.3 0.35 0.4 0.45 0.5 0.55 0.6 0.65 0.7 0.750.60.620.640.660.680.70.720.740.760.78 y u ( y , x = . ) Point value at x=0.5 y u ( y , x = . ) Point value at x=0.3 y u ( y , x = . ) Point value at x=0.2
Figure 3.1: Point evaluations of the solution to (3.1) for α l = 3 and α r = 1, in dependence of y ∈ (cid:2) , (cid:3) .For the evaluation at x = 0 . x = 0 . x ∈ (cid:2) , (cid:3) , while the dependence is smooth for x = 0 . − ( α ( y, x ) u (cid:48) ( y, x )) (cid:48) = 0 , x ∈ (0 , ,u (0) = 1 , u (1) = 0 , for every y ∈ (cid:2) , (cid:3) , with α ( y, x ) = (cid:40) α l if x ∈ (0 , y ) ,α r if x ∈ ( y, , (3.1)where (cid:48) denotes the derivative with respect to x , and α l , α r ∈ R + \ { } , α l (cid:54) = α r . The exactsolution to (3.1) is u ( y, x ) = (cid:40) − α r α l (1 − y )+ α r y x + 1 if x ∈ (0 , y ) , α l α l (1 − y )+ α r y (1 − x ) if x ∈ ( y, , and presents a kink at the interface y . Consequently, the evaluation of the solution at a point x ∈ (0 ,
1) is only continuous (and in particular not analytic) as a function of y if x is a pointthat can be crossed by the interface, that is if x ∈ (cid:2) , (cid:3) . The y -dependence of point values ofthe solution for the points x = 0 . x = 0 . x = 0 . By analogy, we can expect a behavior similar to the one-dimensional case when considering themodel problem (2.1) with coefficients (2.2) which are discontinuous across the interface.To better understand, in the general case, the regularity of the point evaluation with respectto the parameter, we introduce a reference interface ˆΓ := Γ( y = ) = { r ( ϕ ) , ϕ ∈ [0 , π ) } , and aso-called nominal configuration , corresponding to the domain configuration when the interface isˆΓ (see e.g. [11,32,37]). We denote by ˆ D in and ˆ D out,R out , respectively, the inner and outer domainin the nominal configuration. In the following, we use the terminology actual configuration todenote the domain configuration when the interface is Γ( y ), y ∈ P J , and distinguish it from thenominal configuration. The nominal configuration can be mapped to the actual configurationby a parameter-dependent diffeomorphism Φ( y ) : D R out → D R out , y ∈ P J , J ∈ N . This isthe so-called mapping approach, first introduced in [54] and [58]. It is natural, from (2.3),to consider the domain mapping as a perturbation of the identity (see also [53, Sect. 2.8]and [15, Sect. 5.2]). Here we consider x ( y ) = Φ( y ; ˆ x ) = ˆ x + χ ( ˆ x ) ( r ( y ; ˆ ϕ ˆ x ) − r ( ˆ ϕ ˆ x )) ˆ x (cid:107) ˆ x (cid:107) , (3.2)where ˆ x denotes the coordinates in the nominal configuration, ˆ ϕ ˆ x := arg( ˆ x ), and χ : D R out → [0 ,
1] is a mollifier. From now on, we assume the following on χ : it acts on the radial componentof ˆ x , it is supported in [ r − , ˜ R ] for some ˜ R ≤ R out , it is strictly decreasing in [ r − , r ] and6trictly increasing in [ r , ˜ R ], it has at least the same smoothness as the nominal radius r , andmax (cid:110) (cid:107) χ (cid:107) C ( ˆ D in ) , (cid:107) χ (cid:107) C ( ˆ D out,Rout ) (cid:111) < √ r − . These assumptions guarantee that Φ as in (3.2) is anorientation preserving diffeomorphism with the same spatial smoothness as the the radius r in(2.3), and the singular values of its Jacobian matrix have J -, y - and ˆ x -independent lower andupper bounds σ min and σ max , see Sect. 3.2 and Appendix E in [37] for details. Since in thefollowing we will use H¨older norms of Φ( y ; · ) and its inverse, here we note that this domainmapping is well-defined as random variable taking values on H¨older spaces: although thesespaces are not separable, it is clear from (2.3) and (3.2) that Φ and its inverse take values inseparable subspaces of these spaces.With the mapping Φ at hand, we can rewrite the variational formulation of (2.1) on thenominal configuration asFind ˆ u ( y ) ∈ H ( D R out ) : (cid:90) D Rout ˆ α ( y ; ˆ x ) ˆ ∇ ˆ u ( y ) · ˆ ∇ ˆ v − ˆ κ ( y ; ˆ x )ˆ u ( y )ˆ v d ˆ x − (cid:90) ∂D Rout
DtN(ˆ u ( y ))ˆ v d S = (cid:90) ∂D Rout (cid:18) ∂u i ∂ n out − DtN( u i ) (cid:19) ˆ v d S, for every ˆ v ∈ H ( D R out ) and y ∈ P J , (3.3)where ˆ u := Φ ∗ ( u ) (Φ ∗ denoting the pullback with respect to Φ), ˆ n is the normal to ˆΓ pointingto ˆ D out,R out , ˆ ∇ denotes differentiation with respect to ˆ x andˆ α ( y ; ˆ x ) = D Φ( y ) − D Φ( y ) −(cid:62) det D Φ( y ) α ( y ; Φ − ( y ; x )) , ˆ κ ( y ; ˆ x ) = det D Φ( y ) κ ( y ; Φ − ( y ; x )) , (3.4)with D Φ the Jacobian matrix of Φ. Since Φ( y ; · ) and its inverse take values in separablesubspaces of H¨older spaces, the coefficients (3.4) also do, and they are well-defined as H¨older-space valued random variables.We have mentioned in Section 2 that we work in the large-wavelength regime. In quantitativeterms, this means that we assume the following: Assumption 3.1 (Large wavelength assumption) . The wavenumbers in (2.2) satisfy the con-dition: κ , κ ≤ τ inf w ∈ H ( D Rout ) | w | H ( D Rout ) + (cid:107) w (cid:107) L ( ∂D Rout ) (cid:107) w (cid:107) L ( D Rout ) , (3.5) with τ < σ min σ max min (cid:110) α , α (cid:111) (where σ min and σ max are the J -, y - and ˆ x -independent upper andlower bounds on the singular values of D Φ ). In [37, Sect. 5.3], we have shown that, provided Assumption 3.1 holds, ˆ u depends smoothlyon y (it admits a holomorphic extension to polyellipses in the complex plane). This is possiblebecause, on the nominal configuration, the interface ˆΓ is fixed, for every parameter realization.However, this is not the case when considering the solution on the actual configuration.Since we are interested in the evaluation of the solution u to (2.1) at points that may belocated on either side of the interface, we introduce the set P Γ J ( x ) := { y ∈ P J : x ∈ Γ( y ) } = (cid:110) y ∈ P J : Φ − ( y ; x ) ∈ ˆΓ( y ) (cid:111) . (3.6)Due to the affine parametrization (2.3) for the interface, for every x ∈ R the set P Γ J ( x ), ifnot empty, is a hyperplane, affine to a ( J − P J .7n the following proposition we show that, for the point evaluation in the actual config-uration, the smoothness with respect to y is C or C , and in general not C k for k ≥ C the continuity constant of the bilinear form a p ( v, w ) := (cid:104) ˆ ∇ ˆ v, ˆ ∇ ˆ w (cid:105) − (cid:104) ˆ ∇ ˆ v · n out , ˆ w (cid:105) (cid:104) H − ( ∂D Rout ) ,H ( ∂D Rout ) (cid:105) on H ( D R out ), and by γ p the coercivity constant of a p restrictedto functions that satisfy the radiation condition. Proposition 3.2.
Let u be the solution to (2.1) with coefficients (2.2) . Let Assumptions2.1 and 3.1 hold, and let us assume that we can build the mapping Φ in (3.2) such that σ min σ max min (cid:110) α , α (cid:111) ≥ − γ p C . Consider x ∈ D R out such that P Γ J ( x ) is not empty.For every J ∈ N , the map y (cid:55)→ u ( y ; x ) from P J to C is continuous, and (cid:107) u ( · ; x ) (cid:107) C ( P J ) has a J -independent upper bound.If α = 1 in (2.2) , then the map y (cid:55)→ u ( y ; x ) is of class C , and (cid:107) u ( · ; x ) (cid:107) C ( P J ) has a J -independent upper bound.Proof. We first consider the general case α (cid:54) = 1. Using the mapping from the nominal config-uration, we can write:max y ∈P J | u ( y ; x ) | = max y ∈P J | ˆ u ( y ; Φ − ( y ; x )) | ≤ max y ∈P J (cid:107) ˆ u ( y ; · ) (cid:107) C ( D Rout ) . (3.7)This means that it is sufficient to show that the mapping y (cid:55)→ ˆ u is continuous from P J to C ( D R out ), with a J -independent bound on the last term in (3.7).Assumption 2.1 ensures that, for every y ∈ P J , (cid:107) r ( y ; · ) (cid:107) C ,βper ([0 , π )) has a J - and y -independentbound for some β ∈ (0 , α ( y ; · ) and ˆ κ ( y ; · )belong to C β ( ˆ D in ) ∪ C β ( ˆ D out,R out ) for every y ∈ P J , with J -independent bounds on the norms.Then, using Assumption 3.1 to ensure coercivity, we can apply Lemma 2 in [40] and a slightgeneralization of Theorem 3.1 in [7] to conclude that every realization of the scattered waveˆ u s ( y ; · ) := ˆ u ( y ; · ) − u i (Φ( y ; · )) has a J - and y -independent bound on its H β (cid:48) ( D R out )-norm, for0 < β (cid:48) < β sufficiently small. This implies that the right-hand side in (3.7) has a J -independentbound, thanks to the Sobolev embedding theorem [25, Thm. 7.26]. The adaptation of Theorem3.1 in [7] to our case is reported in Appendix A, and it requires that σ min σ max min (cid:110) α , α (cid:111) ≥ − γ p C .The smoothness with respect to y is limited by the spatial smoothness of ˆ u , due to theapplication of the chain rule on ˆ u ( y ; Φ − ( y ; x )). Thus, if α = 1 in (2.2), we obtain higherregularity with respect to the parameter.In particular, for a generic α , the J - and y -independent upper bound on (cid:107) ˆ u ( y ; · ) (cid:107) C ( D Rout ) and the H¨older regularity of the PDE coefficients (with J - and y -independent norm bounds)also imply a J - and y -upper bound on (cid:107) ˆ u ( y ; · ) (cid:107) C β ( ˆ D in ) and (cid:107) ˆ u ( y ; · ) (cid:107) C β ( ˆ D out,Rout ) [25, Thm.8.33]. We specify that the result in [25] is for a boundary value problem, but it can be adaptedto a interface problem proceeding as elaborated in Appendix B. If α = 1, then the transmissionconditions at ˆΓ ensure that ˆ u ∈ C ,β ( D R out ). It is shown in [49, Lemma 4.3.8] that r = r ( y ),as C -valued map, is analytic with respect to y even in the case that P J has infinite dimension( J → ∞ ). Then, thanks to the regularity of Φ in (3.2) with respect to y [49, Lemma 4.3.9],the claim for α = 1 follows from chain rule. Remark 3.3 (Alternative proof) . The assumption that σ min σ max min (cid:110) α , α (cid:111) ≥ − γ p C in the aboveproposition is due to a technicality when adapting the proof of Theorem 3.1 in [7] to our case,because of the presence of the DtN map (see Appendix A for details). Such requirement can bedropped if, instead of using the result in [7], we show H -regularity of ˆ u in each subdomain. Inthe latter case, though, we would have a J -independent bound on the right-hand side in (3.7) only when, in Assumption 2.1, p < . Moreover, the requirement σ min σ max min (cid:110) α , α (cid:111) ≥ − γ p C s not needed at all if instead of boundary conditions with the DtN map we have Dirichlet orNeumann boundary conditions [7].
The proof of Proposition 3.2 shows that the smoothness of the point evaluation with respect tothe high-dimensional parameter depends on the spatial smoothness of the solution u to (2.1)across the interface. In particular, as the solution is in general only continuous (or C ) across theinterface, we cannot expect, in general, a holomorphic dependence of the point evaluation withrespect to y . Consequently, it is likely that high order quadrature methods such as stochasticGalerkin [52], stochastic collocation [1, 51, 56] or high order quasi-Monte Carlo rules [18], willnot show full convergence rates when trying to compute statistics of u ( y ; x ) for a point x ∈ R that can be crossed by the interface. We show this effect on sparse grids.We have run the Smolyak adaptive algorithm described in [51] using R -Leja quadraturepoints. We have set J = 16 and β j − = β j = r j − for j = 1 , . . . ,
8, where r = 0 .
01 isthe nominal radius (we work with non-dimensional quantities), chosen to be constant. Thecoefficients in (2.1) have been chosen as α = 4, κ = κ and κ = 4 κ , with k = 209 . Thequantity of interest is Re u h := { Re u h ( y ; x i ) } N − i =0 , for some N ∈ N , where u h ( y ; x i ) denotesthe discrete approximation to u ( y ; x i ), i = 0 , . . . , N −
1. In the case of holomorphic parameterdependence of the Q.o.I., we would expect a convergence rate of s = 2 − ε with respect to thecardinality of the index set Λ, for ε > estimated error (cid:80) ν ∈N (Λ) (cid:107) ∆ Qν (Re u h ) (cid:107) ∞ computed by the algorithm at each iteration; here N (Λ) is the set of neighbors of the index set Λ, and ∆ Qν are the difference operators (seee.g. [51] for their definition). The left plot in Figure 3.2 shows the convergence of the algorithmwhen applied to one point evaluation ( N = 1), for different points on the horizontal axis.We can see that the curve saturates if the point is crossed by the interface Γ( y ) for manyparameter realizations ( x = ( r , r = 0 . x = (0 . ,
0) and x = (0 . , x = (0 . ,
0) and x = (0 . , N ≥ N = 8 point evaluations, we observe a convergence rate of 0 . Computational details: Smolyak algorithm run on a discrete solution obtained from finite element dis-cretization of (3.3) with piecewise linear, globally continuous ansatz functions; the DtN map (2.1c) has beenapproximated using a circular PML, starting at R out = 0 .
055 and ending at R (cid:48) out = 0 . .
5; the mesh is quasi-uniform and consists of 558705 nodes; domain mapping (3.2) with χ ( ˆ x ) = (cid:107) ˆ x − r (cid:107) r − r for (cid:107) ˆ x (cid:107) ∈ (cid:2) r , r (cid:3) , χ ( ˆ x ) = R out −(cid:107) ˆ x (cid:107) R out − r for (cid:107) ˆ x (cid:107) ∈ [ r , R out ] and χ ( ˆ x ) = 0 elsewhere (the interface (cid:107) ˆ x (cid:107) = r is resolved, so that the discontinuity in the mapping does not affect the finite element convergence); software:NGSolve finite element library (http://sourceforge.net/apps/mediawiki/ngsolve) coupled with the MKL versionof PARDISO (https://software.intel.com/en-us/intel-mkl) for the direct solver. − − − − − − (cid:93) Λ E s t i m a t e d e rr o r Mean of Re u , J = 16 x = 0 . x = 0 . x = 0 . x = 0 . x = 0 .
015 10 − − − − − − (cid:93) Λ E s t i m a t e d e rr o r Mean of Re u , J = 16 − − − − − − (cid:93) PDE solves E s t i m a t e d e rr o r Mean of Re u , J = 16 Figure 3.2: Convergence of adaptive Smolyak algorithm for one point evaluation at different points(left) and more point evaluations simultaneously, the latter with respect to the cardinality of the indexset Λ (center) and to the number of PDE solves (right). In the first plot, x denotes the first coordinateof x . The error reported is the one estimated by the algorithm, (cid:80) ν ∈N (Λ) (cid:107) ∆ Qν (Re u h ) (cid:107) ∞ , see [51] forits definition. Two strategies can be identified in order to handle the loss of smoothness of the Q.o.I. with re-spect to the parameter: detect the surface of non-smoothness and apply a high order quadraturemethod in each subdomain of P J where the Q.o.I. is smooth, or adopt a low order quadraturemethod requiring less smoothness of the Q.o.I.The problem with the first strategy is that, in our case, the discontinuities are not easy totrack. For a single point evaluation, we could apply already existing discontinuity detectiontechniques, see [59] and references therein. For multiple evaluations, the complexity of thesurface of non-smoothness increases with the number of points, as we have a hyperplane ofdiscontinuity for each of them. In such a case, the method proposed in [59] cannot be appliedanymore; we do not exclude that the algorithm in [59] could be adapted to tackle multiplediscontinuity detection, but its complexity would probably grow with the number of hyperplanesof discontinuity. Another possibility to pursue the first strategy is to adopt an approach basedon X-FEM in the parameter space, as proposed in [46] in the framework of a level set approach todescribe the uncertain geometry (and named X-SFEM by the authors). However, the algorithmproposed there to track the uncertain boundary seems to be applicable only to a low-dimensionalparameter space (cf. in particular Sect. 6.2.3 in [46]).Our choice is therefore the second strategy, namely to use a quadrature rule that doesnot suffer from the so-called ‘curse of dimensionality’ and requires weak assumptions on theregularity of the Q.o.I.. Namely, a quadrature rule that does not need information aboutthe location of the surfaces of non-smoothness to provide the full convergence rate. We opttherefore for a Monte Carlo approach and in particular, in order to reduce the computationaleffort, to its multilevel version (MLMC). The latter requires only square integrability of theQ.o.I., which is a much weaker smoothness assumption than those for high order quadraturemethods. Moreover, since the convergence rate of the MLMC integration does not depend on10he dimension of the parameter space, such quadrature rule is well suited for high-dimensionalproblems. In this section we provide a brief overview of MLMC, in particular for a Q.o.I. depending onthe solution of an elliptic PDE. Our survey is based on [26], [14] and [4], and we use the samenotation as in the introduction.Due to the need to solve a PDE to compute the Q.o.I., usually we do not have at ourdisposal the quantity Q itself, but an approximation to it. We consider a sequence ( X l ) l ≥ offinite-dimensional subspaces of XX ⊂ X ⊂ . . . ⊂ X l ⊂ . . . ⊂ X , (4.1)with X l associated to the discretization parameter h l ∈ R , l ∈ N = N ∪ { } . Thinking of h l asthe meshsize at level l , we can assume, without loss of generality, that h ≥ h ≥ h ≥ . . . . Wedenote by u l the discrete solution to (1.1) on the level l :Find u l s.t. : u l ( y ) ∈ X l for every y ∈ P J , D l ( y ; u l ( y )) = 0 for every y ∈ P J , (4.2)where the subscript in D l denotes the discretization of D in (1.1) at level l . With u l at ourdisposal, for some l ∈ N , we can compute the approximation of Q at the level l for every y ∈ P J , which we denote by Q l ( y ) := q l ( y ; u l ( y )). Note that the discretization error in Q l might be due not only to the replacement of u by u l , but also to an error coming from thecomputation of q on u l (for example, if q is an output functional defined as an integral quantitythat needs numerical integration).The MLMC method is a modification to the single-level Monte Carlo (MC) algorithm inorder to improve the computational efficiency. In single-level Monte Carlo, the quantity E µ J [ Q ]is estimated by E M [ Q L ] := 1 M M (cid:88) i =1 Q iL ∈ X L , (4.3)where we have assumed that u is approximated by the solution u L to (4.2) at a fixed level L ∈ N ,and Q iL := q L ( y i ; u L ( y i )), i = 1 , . . . , M , M ∈ N , are independent, identically distributedrealizations of Q L ( y ), y ∈ P J . Note that the definition (4.3) is independent of J ∈ N . Theapproximation error of the MC estimator can be decomposed as (cf. [4, Sect. 4.2]) (cid:107) E µ J [ Q ] − E M [ Q L ] (cid:107) L ( P J , Y ) ≤ E µ J ( (cid:107) Q − Q L (cid:107) Y ) + 1 √ M Var µ J [ Q L ] , (4.4)provided Q and Q L have finite variance. The norm on the left-hand side is defined as (cid:107) v (cid:107) L p ( P J , Y ) = (cid:16)(cid:82) P J (cid:107) v ( y ) (cid:107) p Y d µ J ( y ) (cid:17) p if p < ∞ , esssup y ∈P J (cid:107) v ( y ) (cid:107) Y if p = ∞ (4.5)(analogous definition holds when replacing Y by X or any other separable Banach space). Thefirst summand in (4.4) measures the bias of Q L with respect to Q , and depends on the spatialdiscretization error, that is on the accuracy with which Q L approximates Q for every parameterrealization. The second summand is the so-called sampling error, depending on the variance11ar µ J of Q L and the number of samples M . To balance the two error contributions for acertain threshold on the total error, √ M has to be chosen to be inversely proportional to thediscretization error. For fine meshes, this can be very expensive.The idea of the multilevel version of Monte Carlo is to reduce the error contribution fromthe second summand in (4.4) with a lower computational effort than the MC method. Settingby convention Q − := 0 and exploiting that E µ J [ Q L ] = L (cid:88) l =0 E µ J [ Q l − Q l − ] , the classical MLMC method consists in estimating E µ J [ Q ] by E L [ Q ] := L (cid:88) l =0 E M l [ Q l − Q l − ] , (4.6)for a given L ∈ N (again, note that the definition (4.6) is independent of J ). In this way, thevariance is reduced at each level l ∈ N estimating the mean of the tail Q l − Q l − , and this is thereason why the MLMC is also said to be a variance reduction technique. If Var µ J [ Q l − Q l − ]decreases as l increases, as we will see to be usually the case, then it is possible to savecomputational effort with respect to MC, taking more samples on the coarser grids and onlyfew samples on the finer ones. In other words, the advantage of MLMC consists in balancingthe two opposite effects, as l increases, of the decay of Var µ J [ Q l − Q l − ] and the increase ofWork l , the computational cost to compute a sample of Q l − Q l − . This balancing is achieved bydetermining the optimal number of samples M l for each of the levels l = 0 , . . . , L (and possiblyalso the optimal maximal level L ) in order to achieve a certain accuracy for the total error atminimal computational cost. Theorem 4.1 (Theorem 1 in [14]) . Suppose that, for every J ∈ N , there exist positive constants α, β, γ, C , C and C , independent of l ∈ N , such that α ≥ min( β, γ ) , and(i) (cid:107) E µ J [ Q l − Q ] (cid:107) Y ≤ C h αl ,(ii) Var µ J [ Q l − Q l − ] ≤ C h βl ,(iii) Work l ≤ C h − γl .Then, for every ε < e − , there exist a value L and a sequence ( M l ) Ll =0 such that (cid:107) E L [ Q ] − E µ [ Q ] (cid:107) L ( P J , Y ) < ε, (4.7) for every J ∈ N , and there exists a positive constant C such that the total computational cost Work tot ( E L ) is bounded by Work tot ( E L ) ≤ C C ε − if β > γ,C C ε − (log ε ) if β = γ,C C ε − − γ − βα if β < γ. (4.8) If the constants C and C are independent of J ∈ N , then C is independent of J ∈ N . The proof can be found in Appendix A of [14], where it can be checked that the constant C is dependent on C and C but not on C . In general, however, we cannot expect the cost C of a single solve Work l to be independent of J .12he previous theorem indicates that, in order to compute the optimal distribution of sampleson each level, it is necessary to determine the values of the exponents α, β and γ . The valueof γ depends on the method used to discretize (1.1) and on the quantity of interest Q . Forexample, if Q is the solution ifself and (4.2) correspond to linear finite element discretizations,then a multigrid solver for the linear system has linear complexity with respect to the numberof degrees of freedom, and we can set γ = d , with d = 1 , , α and β can be determined, instead, from the convergence estimateof Q l to Q as l → ∞ . Proposition 4.2.
Assume that Q ∈ L ( P J , Y ) for every J ∈ N , and that there exists a constant C A > , independent of h l , l ∈ N , of y ∈ P J and of J ∈ N , and a positive real number t ,independent of l , such that the approximations of Q fulfill (cid:107) q ( y ; u ( y )) − q l ( y ; u l ( y )) (cid:107) Y ≤ C A h tl (cid:107) u ( y ) (cid:107) W , for every y ∈ P J \ N P J and every J ∈ N , (4.9) for a subspace W ⊂ X , and N P J any nondense subset of P J with measure zero. Moreover, let u ∈ L ( P J , W ) with a J -independent norm bound.Then, given a geometric sequence ( h l ) l ≥ of discretization parameters: • the inequality ( i ) in Theorem 4.1 holds with a constant C independent of J and α = t ; • the inequality ( ii ) in Theorem 4.1 holds with a constant C independent of J and β = 2 t .Proof. We first address the bound ( i ) in Theorem 4.1. From the properties of the Bochnerintegral and (4.9) we have, for every l ∈ N : (cid:107) E µ J [ Q l − Q ] (cid:107) Y ≤ E µ J [ (cid:107) Q l − Q (cid:107) Y ] ≤ C A h tl (cid:107) u (cid:107) L ( P J , W ) . Thus we obtain the bound ( i ) in Theorem 4.1 with C = C A sup J ∈ N (cid:107) u (cid:107) L ( P J , W ) and α = t .For the bound ( ii ) in Theorem 4.1, we have:Var µ J [ Q l − Q l − ] ≤ E µ J [ (cid:107) Q l − Q l − (cid:107) Y ] = (cid:107) Q l − Q l − (cid:107) L ( P J , Y ) ≤ (cid:107) Q l − Q (cid:107) L ( P J , Y ) + 2 (cid:107) Q l − − Q (cid:107) L ( P J , Y ) . (4.10)(4.11)Owing to (4.9), the first summand is bounded by (cid:107) Q l − Q (cid:107) L ( P J , Y ) ≤ C A h tl sup J ∈ N (cid:107) u (cid:107) L ( P J , W ) . The analogous holds for the second summand in (4.11) replacing h l by h l − . We remind that weassume a geometric sequence of discretization parameters, that is h l − h l ≤ C H for every l ∈ N and some C H >
0. Then, substituting the above bounds in (4.11), we obtain (cid:107) Q l − Q l − (cid:107) L ( P J , Y ) ≤ C A ( h tl + h tl − ) (cid:107) u (cid:107) L ( P J , W ) ≤ C A (1 + C tH ) h tl (cid:107) u (cid:107) L ( P J , W ) , that is the bound ( ii ) in Theorem 4.1 holds with C = 2 C A (1 + C tH ) sup J ∈ N (cid:107) u (cid:107) L ( P , W ) and β = 2 t .It is clear from the proof that we could slightly relax the assumption on the constant C A ,requiring it to belong to L ( P J , R ) with a J -independent bound, instead of being y -independent.We also note that in (4.10) we have shown that (cid:107) Q l − Q l − (cid:107) L ( P ,Y ) ≤ C h βl , which is astronger requirement that Var µ [ Q l − Q l − ] ≤ C h βl . Such choice gives automatically α ≥ β inTheorem 4.1 (cf. [27, Sect. 2.1]). 13 emark 4.3 (log factors in convergence rates) . The statement of Proposition 4.2 can be adaptedeasily to the case when the convergence rate in (4.9) is of the kind h tl | log h l | ¯ t (cid:107) u ( y ) (cid:107) W , for somepower ¯ t > of log h l . Indeed, as l → ∞ , h tl | log h l | ¯ t ≤ h t (cid:48) for any t (cid:48) < t and any ¯ t > , and onecan use the result of Proposition 4.2 with t (cid:48) in place of t . For the estimate (4.8) we have thenthe following situations: if t > γ , we can choose t (cid:48) such that t > t (cid:48) > γ , and still obtain Work tot ( E L ) ≤ ε − ; if t = γ , then, using t (cid:48) , we switch from the second to the third case, with Work tot ( E L ) ≤ ε − δ , for any δ > ; if t < γ , then Work tot ( E L ) ≤ ε − − γ − tt + δ for any δ > . In the next section we establish under which conditions the assumptions of Proposition4.2 are fulfilled for our model transmission problem when discretized using finite elements.For the point evaluation, in (4.9) we have q = ( y ; u ( y )) = u ( y ; x ) = ˆ u ( y ; ˆ x ( y )), withˆ x = Φ − ( y ; ˆ x ), and Y = R , and we need to determine the exponent t , the space W , and showthat u ∈ L ( P J , W ) with a J -independent bound. Note that the assumption Q ∈ L ( P J , Y )holds as we have shown in Proposition 3.2 the continuity of the map y (cid:55)→ u ( y ; x ) from P J to C . To determine the space W in (4.9), we have to understand which is the proper convergenceestimate for the point evaluation. Once this has been settled, we can address under whichconditions the solution u to (2.1) belongs to L ( P J , W ) with a J -independent bound.A first observation is that, for every J , P Γ J as defined in (3.6) is a zero measure set in P J (as it is a hyperplane). Thus, in (4.9) we can set N P J = P Γ J , and it is sufficient to determinethe convergence estimate in the case that x ∈ D in or x ∈ D out,R out .A second observation is that (4.9) has to be established for every y ∈ P J \ N P J fixed. There-fore, instead of studying the convergence for u ( y ; x ), we can work in the reference configurationand study the convergence estimate for ˆ u ( y ; ˆ x ) with ˆ x = Φ − ( y ; x ). Note, however, thattruncating (2.3) at the term with index J means setting to zero all the entries of y in positiongrater than J , and, for a fixed realization, ˆ x depends on J , as the mapping Φ − ( y ; · ) does.For this reason, we need a convergence estimate for ˆ u ( y ; ˆ x ) which is uniform in the secondargument (that is, independent of ˆ x ).An option to determine the convergence rate would be to consider the point evaluationas the Dirac delta functional δ ˆ x (ˆ v ) = ˆ v ( ˆ x ), which is a bounded on the space H ε ( ˆ D in ) ∪ H ε ( ˆ D out,R out ), for any ε >
0. If s is the finite element convergence rate in H ( D R out ) (withrespect to the meshwidth) for the solution to (3.3), then we would infer the convergence rate s − ε for the point evaluation, for ε > L ∞ -norm, it is possible to achieve a convergence rate of h s +1 (log | h | ) ¯ s for h → s = 1 if s = 1, and ¯ s = 0 if s ≥ W s +1 , ∞ -regularity. Wewill state these convergence estimates rigorously in subsection 5.2, after having established, inthe next subsection, the regularity of the solution to (3.3) in the space W = C k ( ˆ D in ) ∪ C k ( ˆ D out,R out ) , (5.1)for some k ≥
2, equipped with the norm (cid:107)·(cid:107) W := max (cid:110) (cid:107)·(cid:107) C k ( ˆ D in ) , (cid:107)·(cid:107) C k ( ˆ D out,Rout ) (cid:111) .14 .1 Space regularity of the solution To obtain upper bounds on C k ( ˆ D in ) ∪ C k ( ˆ D out,R out ) for some k ≥
2, we consider Schauderestimates (see [25, Ch. 6] and [55, Ch. 6]), as they require milder space regularity of thecoefficients (3.4) than Sobolev estimates [25, Ch. 8] followed by an application of the Sobolevembedding theorem [25, Thm. 7.26].Starting from k = 2, we notice that it is not possible to obtain bounds on the norm of ˆ u in C ( ˆ D in ) ∪ C ( ˆ D out,R out ), because estimates in this last norm are in general not well defined(cf. p.52 and Problem 4.9 in [25]). For this reason, we state estimates in the H¨older spaces C k,βpw ( D R out ) := C k,β ( ˆ D in ) ∪ C k,β ( ˆ D out,R out ), β ∈ (0 , Theorem 5.1.
Let β ∈ (0 , and k ≥ , and let ˆΓ and ∂D R out be simple closed curves of class C k,β . Let the coefficients in (3.3) be such that, for every y ∈ P J : ˆ α ( y ; · ) has J -, y - and ˆ x -uniform upper and lower bounds Λ min and Λ max on its singular values, (cid:107) ˆ α ( y ; · ) (cid:107) C k − ,βpw ( D Rout ) ≤ C α and (cid:107) ˆ κ ( y ; · ) (cid:107) C k − ,βpw ( D Rout ) ≤ C κ , with C α and C κ independent of J ∈ N and y ∈ P J . Thenthe solution ˆ u to (3.3) is such that (cid:107) ˆ u ( y ) (cid:107) C k,βpw ( D Rout ) ≤ C (cid:16) (cid:107) ˆ u ( y ) (cid:107) C ( D Rout ) + (cid:107) u i (cid:107) C k,β ( D Rout ) (cid:17) , (5.2) with a constant C = C ( k, β, C α , C κ , Λ min , Λ max ) independent of J ∈ N and y ∈ P J .Proof. This result is a slight modification for interface problems of the Schauder estimatesin [25, Ch. 6] and [55, Ch. 6], taking care of J - and y -independence in the norm bounds. Werefer to Appendix B for details.We have already seen in the proof of Proposition 3.2 that the smoothness of the PDEcoefficients in (3.3) derives from the smoothness of the radius. More precisely, if the sequence( β j ) j ≥ in (2.3) fulfills Assumption 2.1 and the nominal radius r is sufficiently smooth, then,for every J ∈ N and every y ∈ P J , r ( y ) ∈ C k,βper ([0 , π )), with J - and y -independent normbound and (cid:40) k = (cid:98) p − (cid:99) , and β < p − − k if p − ,k = p − , and any β ∈ (0 ,
1) otherwise . (5.3)Using the expression (3.2) for the domain mapping and (3.4) for the PDE coefficients, we have k ≥ p < in Assumption 2.1. The bounds on the singular values of ˆ α ( y ; · )hold if Assumption 3.1 does.To bound (cid:107) ˆ u ( y ) (cid:107) C ( D Rout ) , we note that, if Assumption 3.1 holds (and r is sufficientlysmooth), and if p < in Assumption 2.1, then (cid:107) ˆ u ( y ) (cid:107) H ( D in ) ∪ H ( D out,Rout ) for every y ∈ P J ,and the norm has a J -independent bound [49, Thm. 6.1.7]. Then the Sobolev embeddingtheorem [25, Thm. 7.26] and the continuity of ˆ u across ˆΓ imply (cid:107) ˆ u ( y ) (cid:107) C ( D Rout ) ≤ C (cid:18) (cid:107) u i (cid:107) H ( ∂D Rout ) + (cid:13)(cid:13)(cid:13) ∂u i ∂ n out (cid:13)(cid:13)(cid:13) H ( ∂D Rout ) (cid:19) , (5.4)for every J ∈ N and y ∈ P J , with a constant C = C ( γ − ) independent of J ∈ N and y ∈ P J ,but dependent on the coercivity constant γ − of the bilinear form in (3.3) (which is J - and y -independent, see [49, Lemma 3.2.5]).We arrive then to the following important corollary to Theorem 5.1.15 orollary 5.2. Let the sequence ( β j ) j ≥ in (2.3) fulfill Assumption 2.1 with p < , let r ∈ C k,βper ([0 , π )) with k ≥ and β ∈ (0 , , and let Assumption 3.1 hold. Then the solution ˆ u to (3.3) belongs to C k,βpw ( D R out ) with k ≥ and β ∈ (0 , as in (5.3) , and (cid:107) ˆ u ( y ) (cid:107) C k,βpw ( D Rout ) ≤ C (cid:107) u i (cid:107) C k,β ( D Rout ) . (5.5) The constant C = C ( k, β, γ − , C α , C κ , σ min , σ max ) is independent of J ∈ N and y ∈ P J (here σ min and σ max are the J - and y -independent bounds on the singular values of D Φ , and theother constants are as defined in this subsection).In particular, ˆ u ∈ L ( P J , W ) with a J -independent bound and W as in (5.1) (with k from (5.3) ). We consider the finite element space of globally continuous ansatz functions which are polyno-mials of degree s on each element of a quasi-uniform mesh with meshsize h l > nominal configuration. We denote this space by S sh l ( D R out ). Setting X l := S sh l ( D R out ) and consideringa nested sequence of meshes and thus a geometric sequence of meshsize parameters ( h l ) l ≥ , weare in the framework for MLMC as in (4.1).Our starting point is the L ∞ -estimate for finite element solutions to elliptic boundary valueproblems. Theorem 5.3 (Theorem 2.1 in [50]) . For a domain D ⊂ R n , n ≥ , we consider the bilinearform a bvp ( y ; ˆ w, ˆ v ) := (cid:90) D ˆ α A ( y ; ˆ x ) ˆ ∇ ˆ w · ˆ ∇ ˆ v + ˆ β A ( y ; ˆ x ) · ˆ ∇ ˆ w ˆ v + ˆ κ A ( y ; ˆ x ) ˆ w ˆ v d ˆ x , y ∈ P J , J ∈ N , for every ˆ w, ˆ v ∈ H ( D ) , with ˆ α A ( y ; ˆ x ) ∈ R n × n , ˆ β A ( y ; ˆ x ) ∈ R n and ˆ κ A ( y ; ˆ x ) ∈ R for every ˆ x ∈ D , J ∈ N , y ∈ P J . For s ≥ , let the following assumptions be satisfied:(i) ∂ D is of class C s +4 ;(ii) for every J ∈ N and every y ∈ P J , ˆ α A ∈ C s +3 ( D ) and ˆ β A , ˆ κ A ∈ C s +2 ( D ) , with J - and y -independent bounds on the norms ;(iii) a bvp ( · , · ) has a J - and y -uniform lower, positive bound on the coercivity constant;(iv) the matrix ˆ α A has a J - and y - and ˆ x -uniform lower, positive bound on the ellipticityconstant.Let ˆ w ( y ; · ) ∈ C ( D ) and ˆ w h l ( y ; · ) ∈ S sh l ( D ) satisfy a bvp ( y ; ˆ w ( y ) − ˆ w h l ( y ) , ˆ v h l ) = 0 for all ˆ v h l ∈ S sh l ( D ) . Then there exists a constant C , independent of ˆ w , ˆ w h l , l ∈ N , of J ∈ N and of y ∈ P J such that (cid:107) ˆ w ( y ) − ˆ w h l ( y ) (cid:107) L ∞ ( D ) ≤ Ch l (cid:18) log 1 h l (cid:19) ¯ s inf χ ∈S shl ( D ) (cid:107) ˆ w ( y ) − χ (cid:107) C ( D ) , (5.6) for every l ∈ N , J ∈ N and y ∈ P J , with ¯ s = 1 if s = 1 and ¯ s = 0 if s ≥ . According to Remark 1.1 in [50], we would need ˆ α A ∈ C s +2 ( D ). However, the reference provided therefor this claim is [41], according to which (see p.107) we need the higher order coefficient in C s +2 ( D ) if theoperator is not in divergence form, and thus we need the higher order coefficient in C s +3 ( D ) when consideringthe operator in divergence form. roof. Repeating the proof of Theorem 2.1 in [50], it is easy to check that, under the assumptionof J - and y -uniform bounds on the norms of the coefficients and on the coercivity and ellipticityconstants of the bilinear form, the constant C in (5.6) is J - and y -independent, too.To be more precise, Theorem 2.1 in [50] provides a sharper estimate using a weighted W , ∞ -norm instead of the C -norm on the right-hand side. However, what we are interested in is theconvergence rate rather than a quantitative estimate, and for this the C -norm is sufficient.Moreover, an extension of L ∞ -estimates to the case that D is a convex polygon can be foundin [29] (although in the case of constant coefficients).Going back to our model problem, in the variational formulation (3.3), differently from theassumptions of Theorem 5.3, the coefficients are smooth in ˆ D in and in ˆ D out,R out , but in generalthey are not smooth across ˆΓ. We can expect that, if the interface ˆΓ is resolved ‘well enough’(in a sense to be made precise), then we still achieve the same convergence rates as in Theorem5.3 when discretizing our interface problem. Finite element estimates taking into account theresolution of the interface have been proven in [42] for the convergence in the H - and L -norms.It is plausible that similar results hold for the convergence in the L ∞ -norm, but, to the author’sknowledge, they seem not to be available in the literature. Also in more recent applications of L ∞ -estimates to interface problems [30], the issue of the approximation of ˆΓ is not addressed.Since proving it goes far beyond the scope of this paper, we formulate the following assumption,and test numerically its plausibility for our model problem in the next subsection. Assumption 5.4.
If as domain D we consider D R out = ˆ D in ∪ ˆΓ ∪ ˆ D out,R out , if ˆ α A ∈ C s +3 pw ( D R out ) and ˆ b A , ˆ κ A ∈ C s +2 pw ( D R out ) with J - and y -independent norm bounds, and if every finite elementmesh provides a piecewise s th -order polynomial approximation for ˆΓ , then the result of Theorem5.3 still holds, in the sense that, for ˆ w ∈ C pw ( D ) and ˆ w h l ∈ S sh l ( D ) satisfying a bvp ( y ; ˆ w ( y ) − ˆ w h l ( y ) , ˆ v h l ) = 0 for all ˆ v h l ∈ S sh l ( D ) : (cid:107) ˆ w ( y ) − ˆ w h l ( y ) (cid:107) L ∞ ( D Rout ) ≤ Ch l (cid:18) log 1 h l (cid:19) ¯ s inf χ ∈S shl ( D Rout ) (cid:107) ˆ w ( y ) − χ (cid:107) C ( ˆ D in ) ∪ C ( ˆ D out,Rout ) , (5.7) with ¯ s as in Theorem 5.3 and C a J - and y -independent constant. If we set ˆ α A = ˆ α , ˆ b A ≡ κ A = ˆ κ , then (5.7) gives us the convergence rate for thesolution to (3.3) (the boundary condition with the DtN map is smooth).As in Corollary 5.2, we can deduce the regularity of the coefficients ˆ α and κ in (3.4) fromthe decay of the coefficients in the radius expansion (2.3). Combining this with Corollary 5.2itself, we obtain Theorem 5.5.
Let the sequence ( β j ) j ≥ in (2.3) fulfill Assumption 2.1 with p < s +5 , s ∈ N ,and let the wavenumbers fulfill Assumption 3.1. Let Assumption 5.4 hold and let the finiteelement meshes provide a piecewise s th -order polynomial approximation to ˆΓ . Then the finiteelement solutions ˆ u h l ∈ S sh l ( D R out ) to (3.3) , l ∈ N , satisfy: (cid:107) ˆ u ( y ) − ˆ u h l ( y ) (cid:107) L ∞ ( D ) ≤ Ch s +1 l (cid:18) log 1 h l (cid:19) ¯ s (cid:107) ˆ u ( y ) (cid:107) C s +1 ( ˆ D in ) ∪ C s +1 ( ˆ D out,Rout ) , (5.8) with ¯ s as in Theorem 5.3 and a constant C independent of l ∈ N , of J ∈ N and of y ∈ P J (butdependent on the mesh regularity parameters, on some J - and y -independent bounds on thenorms of the coefficients in (3.3) and on a J - and y -independent lower bound on the coercivityconstant).Moreover, the norm on the right-hand side in (5.8) is bounded independently of J ∈ N and y ∈ P J . roof. The decay of the sequence ( β j ) j ≥ ensures that, for every y ∈ P J and every J ∈ N ,the radius (2.3) belongs to C s +4 ,βper ([0 , π ) for some β ∈ (0 , J - and y -independentnorm bound, see subsection 5.1. Proceeding as in the proof of Corollary 5.2, the smoothnessof the mapping Φ ensures that, for every y ∈ P J and every J ∈ N , ˆ α ( y ; · ) and ˆ κ ( y ; · ) belongto C s +3 ,βpw ( D R out ), with J - and y -independent norm bounds. The J - and y -independent lowerand upper bounds on the singular values of D Φ ensure a J - and y -independent lower boundon the ellipticity constant of ˆ α , which, together with Assumption 3.1, implies a J - and y -uniform lower bound on the coercivity constant of the bilinear form in (3.3) [49, Lemma 3.2.5].Then Theorem 5.3, together with Assumption 5.4, implies the estimate (5.7). Finally, theinterpolation properties of the spaces S sh l ( D R out ) ensure thatinf χ ∈S shl ( D Rout ) (cid:107) ˆ u ( y ) − χ (cid:107) C ( ˆ D in ) ∪ C ( ˆ D out,Rout ) ≤ C (cid:48) h sl (cid:107) ˆ u ( y ) (cid:107) C s +1 ( ˆ D in ) ∪ C s +1 ( ˆ D out,Rout ) , (5.9)for a constant C (cid:48) dependent on the mesh regularity parameters but clearly not on J ∈ N and y ∈ P J . The norm on the right-hand side has a J - and y -independent bound thanks toCorollary 5.2.In Theorem 5.5 we have not formulated any regularity assumption on ∂D R out as we assumeit to be a circle, and thus of class C ∞ . Corollary 5.6.
Under the assumptions of Theorem 5.5, assumption ( ii ) of Proposition 4.2holds with t = 2 − ε and any ε > for linear finite elements, and with t = s + 1 for Lagrangianfinite elements of degree s with s ≥ . Remark 5.7 (Regularity of coefficients) . In order to have a J - and y -independent bound on (cid:107) ˆ u ( y ) (cid:107) C s +1 ( ˆ D in ) ∪ C s +1 ( ˆ D out ) , it is sufficient that the decay parameter p for the sequence ( β j ) j ≥ satisfies p < s +2 , see Theorem 5.1. The stronger requirement that p < s +5 is due to a techni-cality in the proof of the L ∞ -estimate (5.6) presented in [50], requiring stronger smoothness onthe PDE coefficients. In particular, it is needed for the decay estimate of the Green’s functionassociated to (3.3) . One might ask whether such stronger requirement is necessary.The decay estimate on the Green’s function and the space regularity required on the coeffi-cients is reported Lemma 1.1 and Remark 1.1 of [50], which refer to [41] (whose assumptions onthe coefficients can be found on p. 107). It might be, however, that the estimate reported in [41]still holds on milder assumptions on the regularity of the boundary and of the coefficients (cf.estimate (8.3) and Theorem 19.VII in [43], and Theorem 8.1.11, Corollary 8.1.12 and Remark8.1.13 in [8]). In this subsection we show numerical results to validate the convergence estimates of the previ-ous subsections. We address the case s = 1 in Theorem 5.5, because in the MLMC simulationswe will use linear finite elements.As in subsection 3.3, we work with non-dimensional quantities. In (2.2), we set α = 4, α = 1, κ = κ and κ = 2 κ , where κ = 209 .
44 denotes the wavenumber in free space.The incident wave u i ( x ) = e jκ d · x is coming from the left, that is d = (1 , R out = 0 . .
02 and absorption coefficient (or dampingparameter) 0 . r = 0 .
01. In (2.3), weconsider β j − = β j = . r j p , j = 1 , . . . , J , with three decays p = , , , and four dimensiontruncations J = 8 , , ,
64. The case J = 8 will not be used in the numerical experiments for18LMC, but we consider it here in order to better investigate the dependence of the convergenceestimates of Theorem 5.5 on the dimension of the parameter space. The domain mapping is(3.2) with mollifier χ ( ˆ x ) = (cid:107) ˆ x (cid:107) ≤ r , (cid:107) ˆ x (cid:107)− r r − r if r < (cid:107) ˆ x (cid:107) ≤ r , R out −(cid:107) ˆ x (cid:107) R out − r if r ≤ (cid:107) ˆ x (cid:107) ≤ R out . (5.10)The non-smoothness of this mollifier at (cid:107) ˆ x (cid:107) = r can be easily handled treating the circle ofradius r as an additional interface resolved by the finite element meshes, cf. Assumption 5.4,Theorem 5.5 and [42].We consider six nested, unstructed quasi-uniform meshes on the reference configuration,with 581, 2250, 8855, 35133, 139961 and 558705 degrees of freedom, respectively, and use anadditional refinement, with 2232545 degrees of freedom, to obtain reference solutions. Thecircles of radius r and r have been approximated by piecewise linear curves.Each finite element solution has been obtained using the NGSolve finite element library (version 5.1), coupled to the MKL PARDISO direct solver to solve the algebraic system re-sulting from the discretization.We study the convergence of the point evaluation of Re u ( y ; · ), the real part of the solutionto (2.1), for two points in the actual configuration: x = ( r ,
0) and x = (0 , r ). For x ,we consider the realization y with all entries set to 1, so that, for every J , ˆ x = Φ − ( y ; x ) islocated in ˆ D in (although the coordinates of ˆ x depend on J ). For x , we consider the realization y with all entries set to −
1, so that, for every J , ˆ x = Φ − ( y ; x ) ∈ ˆ D out,R out .The results are reported in Figures 5.1 and 5.2. From Theorem 5.5, we expect a convergencerate close to 2 with respect to meshwidth, and thus a rate close to 1 with respect to the numberof degrees of freedom N dof , for every decay p and every dimension J . However, we expectthe constant multiplying the rate in (5.7) (incorporating the norm of the solution) to have a J -independent upper bound only for p < , and thus in none of our test cases. Taking intoaccount Remark 5.7, we could expect J -independence of the constant for p < . Figures 5.1and 5.2 show that the convergence rate predicted by the theory is correct, but the constantseems to have a J -independent upper bound for all values of p considered.The last observation can indicate two things. A possibility is that our theory of subsections5.1 and 5.2 is not sharp and can be improved. Another possible interpretation is that, due tothe decay of the coefficient sequence ( β j ) j ≥ , there is a ‘natural’ dimension truncation from themesh, that does not allow to track the high frequency perturbations. Furthermore, because ofthe nonlinear dependence of the Q.o.I. on the high-dimensional parameter, it could be that, alsowhen the mesh is able to capture some high-frequency shape variations, they contribute to avariation in the Q.o.I. which is smaller than the discretization error. To give an idea about thesize of the shape perturbations, the maximum shape variation for p = is around 0 . r for J = 16 and 0 . r for J = 32, which means that the harmonics added from J = 16 to J = 32contribute for 1 . · − r to the maximum shape variation. The meshsize around r is insteadof the order of 1 . · − on the finest mesh. In Figures 5.1 and 5.2, for p = , we see indeeda slight difference in the convergence curves at the finest level, but it is negligible. Passingfrom J = 32 to J = 64, the contribution of the higher order shape variations is even smallerthan from J = 16 to J = 32, and the convergence curves are indistinguishable. To furtherinvestigate the influence of shape variations, we may ask ourselves how far are the solutionscorresponding to J = 16, J = 32 and J = 64, for a fixed decay p of the coefficient sequence.The fact that the convergence lines are very close to each other gives us no information about http://sourceforge.net/apps/mediawiki/ngsolve − − − − − N dof | R e u ( x ) − R e u l ( x ) | J=8J=16J=32J=64 10 − − − − − N dof | R e u ( x ) − R e u l ( x ) | J=8J=16J=32J=64 10 − − − − − N dof | R e u ( x ) − R e u l ( x ) | J=8J=16J=32J=64
Figure 5.1: Finite element convergence for the point evaluation at x = (0 . , y j = 1for j = 1 , . . . , J , using linear finite elements. Coefficient sequence β j − = β j = j − p with p = 2 (left), p = 3 (center) and p = 4 (right) and j = 1 . . . J . this. We have performed a crossed comparison for each of the cases p = 2 and p = 3: we haveconsidered as reference solution the one obtained on the finest grid for J = 16, and studied theconvergence to this value for the solutions corresponding to J = 32 and J = 64. The outcomefor the evaluation at x = (0 , r ) and with all entries of y set to − p = 2 and p = 3, the solution for J = 32 converges to a value that differs from the exact solution for J = 16 by a quantity thatis some orders of magnitude smaller than the finite element error on the last mesh considered.The right plot in Figure 5.3 shows instead that, for p = 2, the exact solution for J = 16and the exact solution for J = 64 differ by a quantity of the order of 10 − , and this affectsonly the convergence on the last two meshes. Returning to the left plot in Figure 5.2, we seethat the line for J = 64 slightly departs from the line for J = 16. This does not happen for J = 64 and the faster decay p = 3, and in the correponding line in the right plot of Figure5.3 we observe convergence until the last mesh considered. From these last experiments we canconclude that the high frequency perturbations of the shape can be observed only when goingto very fine meshes, supporting the hypothesis of ‘natural’ dimension truncation coming fromthe discretization.Finally, we mention that the achievement of the full convergence rate prescribed by Theorem5.3 when using a piecewise linear approximation for ˆΓ supports the validity of Assumption 5.4. In this section we report the numerical results for the estimation of E µ [Re u ], where u ( y ) = { u ( y ; x i ) } N − i =0 is a set of N point evaluations of the solution u to (2.1). We consider the cases of N = 1 , , , N fixed, x i = r (cos ϕ i , sin ϕ i ) with ϕ i = 2 π iN , i = 0 , . . . , N − β j − =20 − − − − − − N dof | R e u ( x ) − R e u l ( x ) | J=8J=16J=32J=64 10 − − − − − − N dof | R e u ( x ) − R e u l ( x ) | J=8J=16J=32J=64 10 − − − − − − N dof | R e u ( x ) − R e u l ( x ) | J=8J=16J=32J=64
Figure 5.2: Finite element convergence for the point evaluation at x = (0 , . y j = − j = 1 , . . . , J , using linear finite elements. Coefficient sequence β j − = β j = j − p with p = 2 (left), p = 3 (center) and p = 4 (right) and j = 1 . . . J . − − − − − − N dof | R e u ( x ) − R e u l ( x ) | p = 2 p = 3 10 − − − − − − N dof | R e u ( x ) − R e u l ( x ) | p = 2 p = 3 Figure 5.3: Point evaluation at x = (0 , . y j = − j = 1 , . . . , J , using linearfinite elements: convergence of solution for J = 32 to solution for J = 16 (left) and convergence ofsolution for J = 64 to solution for J = 16 (right). In each of the two cases, the coefficient sequences β j − = β j = j − p , j = 1 , . . . , J , with p = 2 and p = 3 are considered. aximal level M M M M M M L = 0 1 L = 1 31 6 L = 2 570 107 20 L = 3 9075 1697 305 54 L = 4 134460 25144 4513 790 136 L = 5 1923719 359729 64557 11293 1943 331 Table 1: Number of samples ( M l ) Ll =0 for the numerical experiments of this section. β j = 0 . r j − p , j = 1 , . . . , J , with p = 2 , ,
4, and dimensions J = 16 , ,
64 of the parameterspace.The physical and geometrical parameters and the domain mapping are as in subsection 5.3.For the MLMC levels, we consider the first five meshes used in the finite element convergencestudies of the previous section, that is unstructed, quasi-uniform meshes with 581, 2250, 8855,35133, and 139961 degrees of freedom, corresponding to L = 0 , . . . ,
4, respectively. The finiteelement setting is as in the previous section (same PML parameters, first order elements, samefinite element solver).The MLMC estimators have been computed using the gMLQMC library [24], with dis-tribution of the samples among the levels determined by solving the optimization problem ofminimizing the total error for a given amount of total computational cost. The work per sam-ple has been estimated as Work l = N dof,l · J , l = 0 , . . . L , where N dof,l is the number of finiteelement degrees of freedom at level l , and J the dimension of the parameter space. The totalwork is calculated as Work tot = (cid:80) Ll =0 Work l . To compute the total error, we have taken intoaccount the logarithmic factor in the convergence rate as from Theorem 5.5. The distributionof the samples among the levels used in all our experiments is reported in Table 1.The error (cid:107) E µ [Re u ] − E L [Re u ] (cid:107) L ( P J , R N ) has been approximated by the average over 10realizations of it, considering, on R N , the Euclidean norm. As reference solution for E µ [Re u ],we use the MLMC estimator E L [Re u ] for L = 5, where the mesh at the fifth level consists of558705 degrees of freedom.Figure 6.1 shows the error versus work for one point evaluation, that is when u = u ( x ) with x = ( r , L = 5, but also with respect to the solution obtained by the Smolyak algorithmwith R -Leja quadrature points before the estimated error saturates (cf. Fig. 3.2). The dashedline reports the theoretical rate of error versus work estimated when running the optimizationalgorithm to choose the number of samples at each level. In Figure 6.2, we compare, for the caseof a 16-dimensional parameter space, the performance of the MLMC estimator with the singlelevel estimator when samples chosen as M = N dof (log N dof ) − (for the single level estimatorthe error is computed over 15 repetitions).Figures 6.3, 6.4 and 6.5 show the performance of MLMC when considering, respectively, 2,4 and 8 point evaluations.From Figures 6.1, 6.2, 6.3, 6.4 and 6.5 we can draw the following conclusions: • the convergence rate of error versus work predicted when running the optimization al-gorithm (dashed line with slope − .
45) is achieved, in all experiments; for low errorthresholds, significant cost savings can be observed when comparing MLMC with singlelevel MC; • the right shift of the error curves as the dimension J of the parameter space increases isonly due to the fact that we compute the work of a single solve as Work l = N dof,l · J , l = 0 , . . . L ; this increase of the computational cost with respect to J is inevitable unless https://gitlab.math.ethz.ch/gantnerr/gMLQMC − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=16, Smkd=32, selfd=32, Smkd=64, selfd=64, Smk 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=16, Smkd=32, selfd=32, Smkd=64, selfd=64, Smk 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=16, Smkd=32, selfd=32, Smkd=64, selfd=64, Smk Figure 6.1: MLMC convergence for 5 mesh levels ( L = 4) and one point evaluation ( u = u ( x ) with x = ( r , β j = ( j (cid:48) ) − p with p = 2 (left), p = 3 (center) and p = 4 (right).Reference solution computed with MLMC on 6 levels ( L = 5, label ‘self’) and with the adaptiveSmolyak algorithm (label ‘Smk’). The dashed line corresponds to the theoretical error versus workrate estimated when running the optimization algorithm to choose the number of samples at eachlevel. − − − − . − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, SLd=16, MLMC 10 − − − − . − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, SLd=16, MLMC 10 − − − − . − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, SLd=16, MLMC Figure 6.2: Comparison of MLMC and single level MC for d = 16, 5 mesh levels and one pointevaluation ( u = u ( x ) with x = ( r , β j = ( j (cid:48) ) − p with p = 2 (left), p = 3(center) and p = 4 (right). Reference solution computed with MLMC on 6 levels. The dashed linecorresponds to the theoretical error versus work rates: 0 .
45 for MLMC and 0 .
33 for single level MC. − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self Figure 6.3: MLMC convergence for 5 mesh levels ( L = 4) and two point evaluations ( u = { u ( x i ) } i =1 ,with x = ( r ,
0) and x = ( − r , β j = ( j (cid:48) ) − p with p = 2 (left), p = 3(center) and p = 4 (right). Reference solution computed with MLMC on 6 levels ( L = 5). The dashedline corresponds to the theoretical error versus work rate estimated when running the optimizationalgorithm to choose the number of samples at each level. − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self Figure 6.4: MLMC convergence for 5 mesh levels ( L = 4) and four point evaluations ( u = { u ( x i ) } i =1 ,with x = ( r , x = (0 , r ), x = ( − r , x = (0 , − r )). Coefficient sequence β j = ( j (cid:48) ) − p with p = 2 (left), p = 3 (center) and p = 4 (right). Reference solution computed with MLMC on 6 levels( L = 5). The dashed line corresponds to the theoretical error versus work rate estimated when runningthe optimization algorithm to choose the number of samples at each level. − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self 10 − − − − . Work (cid:107) E µ [ R e u ] − E L [ R e u ] (cid:107) L ( P J , R ) d=16, selfd=32, selfd=64, self Figure 6.5: MLMC convergence for 5 mesh levels ( L = 4) and eight point evaluations ( u = { u ( x i ) } i =1 ,with x i = (cos( ϕ i ) , sin( ϕ i )), and ϕ i = π ( i − i = 1 , . . . , β j = ( j (cid:48) ) − p with p = 2 (left), p = 3 (center) and p = 4 (right). Reference solution computed with MLMC on 6 levels( L = 5). The dashed line corresponds to the theoretical error versus work rate estimated when runningthe optimization algorithm to choose the number of samples at each level. an algorithm to adapt J to the discretization level is considered, as also suggested in theconclusions in [14]; • the rate of convergence of MLMC is dimension robust; as already observed for the finiteelement convergence, and thus not surprisingly here, the results are even better thanpredicted by theory, in the sense that dimension robustness occurs also for p = 2 ( p = 3is a the limit case, see Theorem 5.5); • requiring only square integrability of the Q.o.I., MLMC is robust with respect to thenumber of singularities in the parameter space, and provides full convergence rate for N = 2 , , u does. We have shown that the MLMC method is effective in computing statistics (in particular themean) of a Q.o.I. whose dependence on the parameter is non-smooth, with discontinuities whichare not easy to track. As model we have considered the computation of point values of thesolution to a Helmholtz transmission problem with stochastic interface. For this case, we haveanalyzed the convergence rate of the finite element discretization and shown how it can be usedto compute the optimal distribution of samples in the MLMC algorithm. Particular attentionhas been dedicated to the robustness of the convergence rates with respect to the dimension ofthe parameter space. The numerical experiments confirm the theoretical results, and show thatmaybe the result on the J -independence of the finite element convergence rate for the pointevaluation can be improved. 25oncerning the application to the point evaluation, we highlight that the results are notconfined to our model problem. The affine parametrization of the stochastic interface does notneed to be in polar coordinates and with respect to the Fourier basis: a more general expansionfor a stochastic interface is possible, as long as C ,β -smoothness is guaranteed. Moreover, theanalysis on the space regularity of the solution to the PDE carries over to any other ellipticPDE associated to a coercive bilinear form, with a parameter-independent lower bound onthe coercivity constant. Finally, the methodology presented in this paper still holds for three-dimensional problems.The results of this work open the way to further investigations. As observed in Proposition3.2, if the highest order PDE coefficient is continuous across the interface (i.e. α = 1 inour model problem), then the solution has C -dependence on the parameter, and it wouldbe interesting to analyze the performance of quasi-Monte Carlo quadrature rules in this case.Another interesting aspect is the possibility to adapt the truncation dimension J to the meshlevels in the MLMC algorithm. As observed in subsection 5.3, indeed, it is likely that on coarserlevels the high-frequency perturbations of the domain cannot be captured by the discretization,and this could be exploited to save computational effort and have a J -dependence of the cost ofone solve which is milder than Work l = N dof,l · J (for l = 0 , . . . , L ). This observation can alsobe found in [14]; a first step in this direction has been done [22], where the truncation levelshave been chosen empirically. Acknowledgements
The author would like to thank Prof. Ralf Hiptmair and Prof. Christoph Schwabfor their suggestions during the development of this work, and for their feedback ondrafts of this paper. She would also like to thank Dr. Vanja Nikoli´c for pointing out thereference [40], and Robert Gantner for his support with the MLMC library. This workhas been funded by ETH under CHIRP Grant CH1-02 11-1 and partially by the TechnicalUniversity of M¨unich.
Appendix A. H β (cid:48) -regularity of the solution In the proof to Proposition 3.2 we have used the fact that, if the coefficients ˆ α ( y ; · ) and ˆ κ ( y ; · )are piecewise H¨older continuous, then the scattered wave ˆ u s ( y ; · ) is in H β (cid:48) ( D R out ) for some β (cid:48) >
0, and from this continuity of ˆ u ( y ; · ) follows. Here we present the H β (cid:48) -regularity resulton the scattered wave, and we do it slightly modifying the proofs contained in [40] and [7, Sect.3]. The scattered wave fulfills the variational formulationFind ˆ u s ( y ) ∈ H ( D R out ) : a ( y ; ˆ u s , ˆ v ) = (cid:104) f s ( y ) , ˆ v (cid:105) (cid:104) ( H ( D Rout )) (cid:48) ,H ( D Rout ) (cid:105) , for every ˆ v ∈ H ( D R out ) , J ∈ N , y ∈ P J , (A.1)where ( H ( D R out )) (cid:48) denotes the dual space of ( H ( D R out )), a ( y ; ˆ u s , v ) := (cid:90) D Rout ˆ α ( y ; ˆ x ) ˆ ∇ ˆ u s ( y ) · ˆ ∇ ˆ v − ˆ κ ( y ; ˆ x )ˆ u s ( y )ˆ v d ˆ x − (cid:90) ∂D Rout
DtN(ˆ u s ( y ))ˆ v d S,f s ( y ) := ˆ ∇ · (cid:16) ˆ α ( y ; ˆ x ) ˆ ∇ ˆ u i ( y ) (cid:17) + ˆ κ ( y ; ˆ x )ˆ u i ( y ) , (A.2)and ˆ u i ( ˆ x ) := u i (Φ( y ; ˆ x )), for all ˆ x ∈ D R out . Our goal is to prove the following:26 heorem A.1. If ˆ α ( y ; · ) ∈ C βpw ( D R out ) , β ∈ (0 , ) , and ˆ κ ( y ; · ) ∈ C pw ( D R out ) , with J - and y -independent bounds on the norms, if Assumption 3.1 holds and σ min σ max min (cid:110) α , α (cid:111) ≥ − γ p C (with γ p , C as in Proposition 3.2), then there exists < β (cid:48) < β such that (cid:107) ˆ u s (cid:107) H β (cid:48) ( D Rout ) ≤ C (cid:107) u i (cid:107) C ( D Rout ) , (A.3) with a constant C independent of J ∈ N and y ∈ P J . We first note that showing the above result for (A.1) is equivalent to showing the result forthe variational formulationFind ˆ u s ( y ) ∈ H ( D R out ) : a p ( y ; ˆ u s , v ) = (cid:104) ˜ f s ( y ) , ˆ v (cid:105) (cid:104) ( H ( D Rout )) (cid:48) ,H ( D Rout ) (cid:105) , for every ˆ v ∈ H ( D R out ) , y ∈ P J , J ∈ N , (A.4)with the low order term of the bilinear form moved to the right-hand side: a p ( y ; ˆ u s , v ) := (cid:90) D Rout ˆ α ( y ; ˆ x ) ˆ ∇ ˆ u s ( y ) · ˆ ∇ ˆ v d ˆ x − (cid:90) ∂D Rout
DtN(ˆ u s ( y ))ˆ v d S, ˜ f s ( y ) := ˆ ∇ · (cid:16) ˆ α ( y ; ˆ x ) ˆ ∇ ˆ u i ( y ) (cid:17) + ˆ κ ( y ; ˆ x )ˆ u ( y ) . (A.5)From now on, we use bold symbols for Sobolev spaces of vector-valued functions; for in-stance, L ( D R out ) := ( L ( D R out )) . Using the notation of [40] and [7], we define the operators J : ( H ( D R out )) (cid:48) → H ( D R out ) and S : L ( D R out ) → ( H ( D R out )) (cid:48) by: (cid:104) ˆ ∇ ( J f ) , ˆ ∇ ˆ v (cid:105) (cid:104) L ( D Rout ) , L ( D Rout ) (cid:105) − (cid:104) DtN( J f ) , ˆ v (cid:105) (cid:104) H − ( ∂D Rout ) ,H ( ∂D Rout ) (cid:105) = (cid:104) f, ˆ v (cid:105) (cid:104) ( H ( D Rout )) (cid:48) ,H ( D Rout ) (cid:105) , (cid:104)S F , ˆ v (cid:105) (cid:104) ( H ( D Rout )) (cid:48) ,H ( D Rout ) (cid:105) = (cid:104) F , ˆ ∇ ˆ v (cid:105) (cid:104) L ( D Rout ) , L ( D Rout ) (cid:105) − (cid:104) F · n out , ˆ v (cid:105) (cid:104) H − ( ∂D Rout ) ,H ( ∂D Rout ) (cid:105) , (A.6)(A.7)for all ˆ v ∈ H ( D R out ), f ∈ ( H ( D R out )) (cid:48) and F ∈ L ( D R out ).On the lines of Lemmas 3.1 and 3.2 in [7], we prove the mapping properties of the operators J and S . Lemma A.2 (Analogous to Lemma 3.1 in [7]) . For all s ∈ [0 , and all F ∈ H s ( D R out ) , S F ∈ H s − ( D R out ) and (cid:107)S F (cid:107) H s − ( D Rout ) ≤ C − s (cid:107) F (cid:107) H s ( D Rout ) , (A.8) with C as in Proposition 3.2.Proof. For s = 0: (cid:104)S F , ˆ v (cid:105) (cid:104) ( H ( D Rout )) (cid:48) ,H ( D Rout ) (cid:105) ≤ C(cid:107) F (cid:107) L ( D Rout ) (cid:107) ˆ v (cid:107) H ( D Rout ) , for every ˆ v ∈ H ( D R out ). For s = 1: (cid:104)S F , ˆ v (cid:105) (cid:104) L ( D Rout ) ,L ( D Rout ) (cid:105) = −(cid:104) ˆ ∇ · F , ˆ v (cid:105) (cid:104) L ( D Rout ) ,L ( D Rout ) (cid:105) ≤ (cid:107) F (cid:107) H ( D Rout ) (cid:107) ˆ v (cid:107) L ( D Rout ) , for every ˆ v ∈ L ( D R out ). The claim follows then from the Riesz-Thorin Theorem.27 emma A.3 (Analogous to Lemma 3.2 in [7]) . For all q ∈ [0 , ) , there exists K = K ( D R out , q, γ p ) such that, for all f ∈ H q − ( D R out ) , J f ∈ H q ( D R out ) and (cid:107)J f (cid:107) H q ( D Rout ) ≤ K (cid:107) f (cid:107) H q − ( D Rout ) , (A.9) and, for all s ∈ [0 , q ] and all f ∈ H s − ( D R out ) , J f ∈ H s ( D R out ) and (cid:107)J f (cid:107) H s ( D Rout ) ≤ (cid:18) γ p (cid:19) − sq K sq (cid:107) f (cid:107) H s − ( D Rout ) , (A.10) with γ p as in Proposition 3.2.Proof. For s = 0 we have: γ p (cid:107)J f (cid:107) H ( D Rout ) ≤ (cid:104) ˆ ∇ ( J f ) , ˆ ∇ ( J f ) (cid:105) (cid:104) L ( D Rout ) , L ( D Rout ) (cid:105) − (cid:104) DtN( J f ) , J f (cid:105) (cid:104) H − ( ∂D Rout ) ,H ( ∂D Rout ) (cid:105) = (cid:104) f, J f (cid:105) (cid:104) ( H ( D Rout )) (cid:48) ,H ( D Rout ) (cid:105) ≤ (cid:107) f (cid:107) ( H ( D Rout )) (cid:48) (cid:107)J f (cid:107) H ( D Rout ) , and thus (cid:107)J f (cid:107) H ( D Rout ) ≤ γ p (cid:107) f (cid:107) ( H ( D Rout )) (cid:48) . For s = q , the inequality follows from Theorem 4 and Remark 4.5 in [48]. The latter ensuresthe existence of a constant K = K ( D R out , q, γ p ) such that (cid:107)J f (cid:107) H q ( D Rout ) ≤ K (cid:107) f (cid:107) H q − ( D Rout ) , for all f ∈ H q − ( D R out ). The estimate (A.10) is then obtained by interpolation.In the next lemma, we use the symbol E ν to denote the multiplier associated to a tensor ν [7, 40]. Lemma A.4 (Similar to Prop. 2.1 in [7]) . If a tensor ν belongs to C βpw ( D R out ) , then, for q ∈ [0 , β ) , there exists a constant C > such that (cid:107)E ν (cid:107) H q → H q ≤ ν max N ν,q , with N ν,q = max (cid:40) , C (cid:107) ν (cid:107) C βpw ( D Rout ) ν max (cid:41) (A.11) and ν max the maximum singular value of ν . Moreover, for s ∈ [0 , q ) , (cid:107)E ν (cid:107) H s → H s ≤ ν max N sq ν,q . (A.12) Proof.
For s = 0, (cid:107)E ν v (cid:107) L ≤ ν max (cid:107) v (cid:107) L , for every v ∈ L ( D R out ). Equation (A.11) is a directconsequence of Lemma 2 in [40]. Then (A.12) is obtained by interpolation.We are now ready to address the proof of Theorem A.1. For this, we proceed on the linesof the proof of Theorem 3.1 in [7], with our modified definition of the operators J and S as in(A.6)-(A.7).For a positive number k >
0, we can write:˜ f s ( y ) = S ( ˆ α ( y ; · ) ˆ ∇ ˆ u s ( y )) = S ˆ ∇ ( k ˆ u s ( y )) − S (cid:18)(cid:18) I − k ˆ α ( y ; · ) (cid:19) ˆ ∇ ( k ˆ u s ( y )) (cid:19) in (cid:0) H ( D R out ) (cid:1) (cid:48) , J ∈ N and every y ∈ P J , where I ∈ R × denotes the identity matrix (for the equationabove, we remind that ˆ α ( y ; ˆ x ) = I for ˆ x ∈ ∂D R out , for every J ∈ N and every y ∈ P J ).We set ¯ α ( y ; · ) := I − k ˆ α ( y ; · ) and ˆ w ( y ) := k ˆ u s ( y ). Since ˆ w ( y ) ∈ H ( D R out ) fulfills theradiation condition, we have that J S ˆ ∇ ˆ w = ˆ w . Thus, we can write:ˆ w ( y ) − Q ( y ; ˆ w ( y )) = J ˜ f s ( y ) , for every J ∈ N and y ∈ P J , with Q ( y ) := J S ( E ¯ α ( y ) ˆ ∇ ). If we can show that Q ( y ) ∈ L (cid:0) H β (cid:48) +1 ( D R out ) , H β (cid:48) +1 ( D R out ) (cid:1) and (cid:107)Q(cid:107) H β (cid:48) +1 → H β (cid:48) +1 ≤ C Q < J ∈ N , every y ∈ P J and some β (cid:48) ∈ (0 , k (cid:107) ˆ u s ( y ) (cid:107) H β (cid:48) +1 ( D Rout ) ≤ (cid:107)J (cid:107) − C Q (cid:107) ˜ f s ( y ) (cid:107) H β (cid:48)− ( D Rout ) , (A.13)for every y . If k is chosen to be independent of J and y , the claim of Theorem A.1 followsonce we prove that (cid:107) ˜ f s ( y ) (cid:107) H β (cid:48)− ( D Rout ) has a J - and y -independent bound.We now show that, for every y , Q ( y ) is a contraction from H β (cid:48) +1 ( D R out ) to H β (cid:48) +1 ( D R out ).For ˆ v ∈ H β (cid:48) +1 ( D R out ), we have that ˆ ∇ ˆ v ∈ H β (cid:48) ( D R out ). Since, by assumption, ˆ α ( y ; · ) and thus¯ α ( y ; · ) are piecewise H¨older continuous, Lemma A.4 ensures that ¯ α ( y ) ˆ ∇ ˆ v ∈ H β (cid:48) ( D R out ) for any0 < β (cid:48) < β , and every J ∈ N , y ∈ P J . Finally, using Lemmas A.2 and A.3, we obtain that Q ( y ; ˆ v ) ∈ H β (cid:48) +1 ( D R out ). Moreover, for every J ∈ N , every y ∈ P J and for 0 < β (cid:48) < q < β : (cid:107)Q ( y ) (cid:107) H β (cid:48) +1 → H β (cid:48) +1 ≤ (cid:107)J (cid:107) H β (cid:48)− → H β (cid:48) +1 (cid:107)S(cid:107) H β (cid:48) → H β (cid:48)− (cid:107)E ¯ α ( y ) (cid:107) H β (cid:48) → H β (cid:48) ≤ (cid:18) γ p (cid:19) − β (cid:48) q K β (cid:48) q C − β (cid:48) ¯ α max ( y ) N β (cid:48) q ¯ α ( y ) ,q ≤ C γ p ( K C − q N ¯ α ( y ) ,q ) β (cid:48) q ¯ α max ( y ) , where for the last inequality we have used that C γ p ≥
1, and ¯ α max ( y ) denotes the maximum sin-gular value of ¯ α ( y ; · ). Let us denote by Λ min and Λ max , respectively, the J - and y -independentlower and upper bounds on the eigenvalues of ˆ α (these bounds exist thanks to the boundson the singular values of D Φ). Then, if we choose k such that 1 − k Λ max >
0, we have0 < − k Λ max ≤ ¯ α max ( y ) ≤ − k Λ min for every J and every y , and thus ¯ α max ( y ) has a J - and y -independent upper bound. For the same reason and because of the J - and y -uniform upperbound on the H¨older norm of ˆ α , if we choose k independent of J and y then N ¯ α ( y ) ,q has a J -and y -independent upper bound N q .The norm (cid:107)Q ( y ) (cid:107) H β (cid:48) +1 → H β (cid:48) +1 has a J - and y -independent bound which is smaller than oneif β (cid:48) < q min , log (cid:18) γ p C ( − k Λ min ) (cid:19) log( K C − q N q ) . (A.14)A β (cid:48) > (cid:0) − k Λ min (cid:1) < γ p C , which, combined with the requirement that 0 < − k Λ max ,implies that we must choose Λ max < k < Λ min − γp C . Such a k exists and can be chosen independentlyof J and y if Λ max < Λ min − γp C , that is Λ min Λ max > − γ p C and this is ensured by the requirement σ min σ max min (cid:110) α , α (cid:111) ≥ − γ p C . Note that N q ≥ K C − q ≥ J S ˆ ∇ ˆ w =ˆ w for all ˆ w in H ( D R out ) that satisfy the radiation condition (see Remark 3.2 in [7]).29o complete the proof of Theorem A.1, we have to show a J - and y -uniform bound on (cid:107) ˜ f s ( y ) (cid:107) H β (cid:48)− ( D Rout ) . We have: (cid:107) ˜ f s ( y ) (cid:107) H β (cid:48)− ( D Rout ) ≤ (cid:107)E ˆ α ( y ) ˆ ∇ ˆ u i ( y ) (cid:107) H β (cid:48) ( D Rout ) + ¯ C (cid:107)E ˆ κ ( y ) ˆ u ( y ) (cid:107) L ( D Rout ) ≤ C (cid:107) ˆ α ( y ; · ) (cid:107) C βpw ( D Rout ) (cid:107) ˆ u i ( y ) (cid:107) H β (cid:48) +1 ( D Rout ) + ¯ C (cid:107) ˆ κ ( y ; · ) (cid:107) C pw ( D Rout ) (cid:107) ˆ u ( y ) (cid:107) L ( D Rout ) with C (the same as in (A.11)) and ¯ C independent of J and y . The norms of the coefficientshave a J - and y -independent bound thanks to the assumptions of Theorem A.1. The norm (cid:107) ˆ u ( y ) (cid:107) L ( D Rout ) can be bounded as (cid:107) ˆ u ( y ) (cid:107) L ( D Rout ) ≤ (cid:107) ˆ u ( y ) (cid:107) H ( D Rout ) ≤ C (cid:18) (cid:107) u i (cid:107) H ( ∂D Rout ) + (cid:13)(cid:13)(cid:13) ∂u i ∂ n out (cid:13)(cid:13)(cid:13) H − ( ∂D Rout ) (cid:19) , where C is a J - and y -independent constant, thanks to Assumptions 3.1 and the properties ofΦ [49, Cor. 3.2.6]. Finally, there exists C independent of J and y such that (cid:107) ˆ u i ( y ) (cid:107) H β (cid:48) +1 ( D Rout ) ≤ C (cid:107) ˆ u i ( y ) (cid:107) C ,β (cid:48) pw ( D Rout ) ≤ C max (cid:110) , (cid:107) D Φ( y ) (cid:107) C β (cid:48) pw ( D Rout ) (cid:111) (cid:107) u i (cid:107) C ,β (cid:48) ( D Rout ) , and (cid:107) D Φ( y ) (cid:107) C pw ( D Rout ) isuniformly bounded with respect to J and y thanks to the uniform bounds on the radius. Appendix B. Schauder estimates for the transmission problem
We present here the proof to Theorem 5.1. We adapt the results of [25, Ch. 6] and [55, Ch.6], stated for boundary value problems, to the transmission problem (3.3), with particularemphasis on having constants which are independent of J ∈ N and y ∈ P J . We first addressthe local regularity at the interface, then the interior regularity, and finally the global regularityestimate. Also, we prove these estimates for k = 2 in (5.2), and extend them for any k ≥ n the spatial dimension, that in Theorem 5.1 is n = 2.We use the following abbreviations for norms and seminorms: Notation
Let v ( x ) be a function on Ω ⊂ R n , n ≥ . For β ∈ (0 , and k ∈ N , we use thefollowing notation for the norms and seminorms in C k (Ω) : (cid:107) v (cid:107) = sup x ∈ Ω | v ( x ) | , | v | k ;Ω = (cid:88) | ν | = k (cid:13)(cid:13)(cid:13)(cid:13) ∂ | ν | v∂ x ν (cid:13)(cid:13)(cid:13)(cid:13) , (cid:107) v (cid:107) k ;Ω = (cid:88) | ν |≤ k (cid:13)(cid:13)(cid:13)(cid:13) ∂ | ν | v∂ x ν (cid:13)(cid:13)(cid:13)(cid:13) , and for the norms and seminorms in C k,β (Ω) : | v | β ;Ω = sup x , x ∈ Ω , x (cid:54) = x | v ( x ) − v ( x ) || x − x | β , | v | k,β ;Ω = (cid:88) | ν | = k (cid:12)(cid:12)(cid:12)(cid:12) ∂ | ν | v∂ x ν (cid:12)(cid:12)(cid:12)(cid:12) β ;Ω , (cid:107) v (cid:107) β ;Ω = (cid:107) v (cid:107) + | v | β ;Ω , (cid:107) v (cid:107) k,β ;Ω = (cid:88) | ν |≤ k (cid:12)(cid:12)(cid:12)(cid:12) ∂ | ν | v∂ x ν (cid:12)(cid:12)(cid:12)(cid:12) β ;Ω . If Ω = Ω in ∪ Γ ∪ Ω out , where Γ is an interface separating the two subdomains, then we denote (cid:12)(cid:12) v ± (cid:12)(cid:12) β ;Ω ± = | v | β ;Ω in ∪ Γ + | v | β ;Ω out ∪ Γ , (cid:13)(cid:13) v ± (cid:13)(cid:13) β ;Ω ± = (cid:107) v (cid:107) β ;Ω in ∪ Γ + (cid:107) v (cid:107) β ;Ω out ∪ Γ , (cid:12)(cid:12) v ± (cid:12)(cid:12) β ;Ω ± = | v | β ;Ω in + | v | β ;Ω out , (cid:13)(cid:13) v ± (cid:13)(cid:13) β ;Ω ± = (cid:107) v (cid:107) β ;Ω in + (cid:107) v (cid:107) β ;Ω out , and analogously for the piecewise- C k and piecewise- C k,β norms and seminorms. .1 Local estimates at the interface ˆΓ Without loss of generality, we assume ˆΓ to be the boundary of the upper half-plane, as every C ,β boundary is C ,β -diffeomorphic to the upper half-plane (with J - and y -independent continuityconstants).The standard technique to prove Schauder estimates for the solution to (3.3) is the methodof solidifying coefficients (see [55, Sect. 6.3.2] and [25, Proof of Thm. 6.2]).Fixed ˆ x ˆΓ ∈ ˆΓ and a ball B R ( ˆ x ˆΓ ) of radius R centered in ˆ x ˆΓ , we can write (3.3) restrictedto B R ( ˆ x ˆΓ ) as − ˆ ∇ · (cid:16) ˆ α ( y ; ˆ x ˆΓ ) ˆ ∇ ˆ u (cid:17) = ˆ F ( y ; ˆ x ) in B − R ( ˆ x ˆΓ ) ∪ B + R ( ˆ x ˆΓ ) , (cid:74) ˆ u (cid:75) B R (ˆ x ˆΓ ) ∩ ˆΓ = 0 , (cid:114) ˆ α ( y ; ˆ x ˆΓ ) ˆ ∂ ˆ u ˆ ∂ ˆ n (cid:122) B R (ˆ x ˆΓ ) ∩ ˆΓ = ˆ g ( y ; ˆ x ) , (B.1a)(B.1b)with ˆ F ( y ; ˆ x ) := ˆ κ ( y ; ˆ x )ˆ u + ˆ ∇ · (cid:16) ( ˆ α ( y ; ˆ x ) − ˆ α ( y ; ˆ x ˆΓ )) ˆ ∇ ˆ u (cid:17) , ˆ g ( y ; ˆ x ) := (cid:114) ( ˆ α ( y ; ˆ x ˆΓ ) − ˆ α ( y ; ˆ x )) ˆ ∂ ˆ u ˆ ∂ ˆ n (cid:122) B R (ˆ x ˆΓ ) ∩ ˆΓ , (B.2)(B.3)and with B + R ( ˆ x ˆΓ ) := B R ( ˆ x ˆΓ ) ∩ ˆ D out,R out and B − R ( ˆ x ˆΓ ) := B R ( ˆ x ˆΓ ) ∩ ˆ D in . We develop our analysistaking ˆ F ( y ; ˆ x ) = ˆ f ( y ; ˆ x ) + ( ˆ α ( y ; ˆ x ) − ˆ α ( y ; ˆ x ˆΓ )) ◦ ˆD ˆ u, (B.4)where ( A ◦ B ) ij = A ij B ij , i, j = 1 , . . . , n , is the Hadamard product for matrices, and ˆD ˆ u denotesthe Hessian matrix of ˆ u . The term ˆ f ( y ; ˆ x ) is a generic right-hand side, possibly including lowerorder terms; in our case, ˆ f ( y ; ˆ x ) = ˆ κ ( y ; ˆ x )ˆ u + (cid:16) ˆ ∇ · ˆ α ( y ; ˆ x ) (cid:17) · ˆ ∇ ˆ u .Since, by assumption, the constant matrix ˆ α ( y ; ˆ x ˆΓ ) is symmetric positive definite, for every y ∈ P J and every J ∈ N there exists an orthonormal matrix M y , dependent on y and J , suchthat M (cid:62) y ˆ α ( y ; ˆ x ˆΓ ) M y = Λ y , (B.5)with Λ y a diagonal matrix dependent of y ∈ P J , J ∈ N . The entries of Λ y have J - and y -uniform lower and upper bounds, because, by assumption, ˆ α ( y ; · ) has J -, y - and ˆ x -uniformlower and upper bounds on its singular values. We denote these bounds by Λ min and Λ max ,respectively.Introducing the change of coordinates ˜ x y = M (cid:62) y ˆ x for ˆ x ∈ B R ( ˆ x ˆΓ ), and using the symbols˜ ∇ y and ˜ ∂ y to denote differentiation with respect to ˜ x y , (B.1) becomes: − ˜ ∇ y · (cid:16) Λ y ˜ ∇ y ˜ u y (cid:17) = ˜ F y in B − R ( ˜ x ˜Γ ) ∪ B + R ( ˜ x ˜Γ ) , (cid:74) ˜ u y (cid:75) ˜Γ y = 0 , (cid:114) Λ y ˜ ∂ y ˜ u y ˜ ∂ y ˜ n y (cid:122) ˜Γ y = ˜ g y , (B.6a)(B.6b)for every y ∈ P J and every J ∈ N . We have denoted ˜ F y := ˆ F ( y ; M y ˜ x y ) and ˜ g y := ˆ g ( y ; M y ˜ x y ).Since M y is orthonormal, B R ( ˆ x ˆΓ ) is mapped to another ball with the same radius R and just adifferent center ˜ x ˜Γ . In (B.6a), B − R ( ˜ x ˜Γ ) is the preimage of B − R ( ˆ x ˆΓ ) under M y , corresponding toa half-ball with radius R and center in ˜ x ˜Γ = M (cid:62) y ˆ x ˆΓ ; the same convention applies for B + R ( ˜ x ˜Γ ).In (B.6b), ˜Γ is a short notation for the preimage of B R ( ˆ x ˆΓ ) ∩ ˆΓ under M y .In the following, C k,β ˜ pw ( B R ( ˜ x ˜Γ )) := C k,β ( B + R ( ˜ x ˜Γ ) ∪ ˜Γ) ∪ C k,β ( B − R ( ˜ x ˜Γ ) ∪ ˜Γ), k ∈ N .In subsection B.4, we provide the proof to the following lemma:31 emma B.1. Let
R > and let ˜ u y be a solution to (B.6) in B R ( ˜ x ˜Γ ) , ˜ x ˜Γ ∈ ˜Γ y . If ˜ F y ∈ C ,β ˜ pw ( B − R ( ˜ x ˜Γ )) and ˜ g y ∈ C ,β (˜Γ y ) , then (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R/ (˜ x ˜Γ ) ≤ C (cid:18) R β (cid:13)(cid:13) ˜ u ± y (cid:13)(cid:13) B ± R (˜ x ˜Γ ) + 1 R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (˜ x ˜Γ ) + 1 R β (cid:107) ˜ g y (cid:107) y + | ˜ g y | ,β ;˜Γ y (cid:19) . (B.7) The constant C = C ( n, Λ max , Λ min , β ) is independent of the center ˜ x ˜Γ of the ball B R , of J ∈ N and of y ∈ P J . Inserting in (B.7) the expressions obtained from (B.4) and (B.3), and denoting ˜ α y :=ˆ α ( y ; M y ˜ x y ) and ˜ f y := ˆ f ( y ; M y ˜ x y ), we obtain: Lemma B.2.
Let < R ≤ and let ˜ u y be a solution to (B.6) in B R ( ˜ x ˜Γ ) , ˜ x ˜Γ ∈ ˜Γ y , with ˆ F given by (B.4) . Let the assumptions of Theorem 5.1 hold. Then (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R/ (˜ x ˜Γ ) ≤ C (cid:18) R β (cid:107) ˜ u y (cid:107) B R (˜ x ˜Γ ) + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (˜ x ˜Γ ) + 1 R β (cid:13)(cid:13)(cid:13) ˜ f ± y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) + (cid:12)(cid:12)(cid:12) ˜ f ± y (cid:12)(cid:12)(cid:12) β ; B ± R (˜ x ˜Γ ) (cid:19) . (B.8) The constant C = C (cid:16) n, Λ min , Λ max , (cid:13)(cid:13) ˜ α ± y (cid:13)(cid:13) ,β ; B ± R (˜ x ˜Γ ) (cid:17) is independent of the center ˜ x ˜Γ of the ball B R , of J ∈ N and of y ∈ P J .Proof. In this proof we denote B R := B R ( ˜ x ˜Γ ) and use the symbol ˜D for the Hessian withrespect to ˜ x y .We have ˜ F y = ˜ f y ( y ; ˜ x y ) + (cid:0) M (cid:62) y ( ˜ α y ( y ; ˜ x y ) − ˜ α y ( y ; ˜ x ˜Γ )) M y (cid:1) : ˜D y ˜ u y and ˜ g y = (cid:114) ( ˜ α y ( y ; ˜ x ˜Γ ) − ˜ α y ( y ; ˜ x y )) M y ˜ ∂ y ˜ u y ˜ ∂ y ˜ n y (cid:122) ˜Γ y . Using interpolation inequalities (cf. [55, Cor. 1.2.1]), together with thefact that M y is orthonormal for every y , and that 0 < R ≤
1, we obtain: (cid:12)(cid:12)(cid:12) M (cid:62) y ( ˜ α y ( y ; ˜ x y ) − ˜ α y ( y ; ˜ x ˜Γ )) M y ◦ ˜D y ˜ u ± y (cid:12)(cid:12)(cid:12) β ; B ± R ≤ (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) B ± R (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R + (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R R β ≤ (cid:13)(cid:13) ˜ α ± y (cid:13)(cid:13) β ; B ± R (cid:16)(cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) B ± R + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:17) ≤ (cid:13)(cid:13) ˜ α ± y (cid:13)(cid:13) β ; B ± R (cid:18) R (cid:13)(cid:13) ˜ u ± y (cid:13)(cid:13) B ± R + 2 R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) ≤ (cid:13)(cid:13) ˜ α ± y (cid:13)(cid:13) β ; B ± R (cid:18) R β (cid:13)(cid:13) ˜ u ± y (cid:13)(cid:13) B ± R + 2 R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) , (cid:13)(cid:13)(cid:13) M (cid:62) y ( ˜ α y ( y ; ˜ x y ) − ˜ α y ( y ; ˜ x ˜Γ )) M y ◦ ˜D y ˜ u ± y (cid:13)(cid:13)(cid:13) B ± R ≤ R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) B ± R (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R ≤ R β (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R (cid:18) R β (cid:13)(cid:13) ˜ u ± y (cid:13)(cid:13) B ± R + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) | ˜ g y | ,β ;˜Γ y ≤ (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) B ± R (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) ,β ; B ± R + (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) B ± R + (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) B ± R (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R ≤ (cid:13)(cid:13) ˜ α ± y (cid:13)(cid:13) ,β ; B ± R (cid:18)(cid:18) R + 1 R β + 1 R (cid:19) (cid:13)(cid:13) ˜ u ± y (cid:13)(cid:13) B ± R + (cid:0) R β + R + R β (cid:1) (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) ≤ (cid:13)(cid:13) ˜ α ± y (cid:13)(cid:13) ,β ; B ± R (cid:18) R β (cid:13)(cid:13) ˜ u ± y (cid:13)(cid:13) B ± R + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) , (cid:107) ˜ g y (cid:107) y ≤ R β (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) B ± R ≤ R β (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R (cid:18) R (cid:107) ˜ u y (cid:107) B R + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) = R β (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R (cid:18) R (cid:107) ˜ u y (cid:107) B R + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) ≤ R β (cid:12)(cid:12) ˜ α ± y (cid:12)(cid:12) β ; B ± R (cid:18) R β (cid:107) ˜ u y (cid:107) B R + R β (cid:12)(cid:12) ˜ u ± y (cid:12)(cid:12) ,β ; B ± R (cid:19) . Inserting these estimates in (B.7), we obtain (B.8).We now return to the variable ˆ x . We notice that M y ( B R ( ˜ x ˜Γ )) = B R ( ˆ x ˆΓ ), and, thanks tothe orthonormality of M y , the H¨older norms in the ˆ x -space and in the ˜ x y -space do coincide,for every ˆ x ˆΓ ∈ ˆΓ, every y ∈ P J and every J ∈ N . Let us denote by C α the J - and y -uniformupper bound on the C ,βpw ( D R out )-norm of ˆ α ( y ; · ). From (B.8), we have the following estimate: (cid:12)(cid:12) ˆ u ± (cid:12)(cid:12) ,β ; B ± R/ (ˆ x ˆΓ ) ± ≤ C (cid:18) R β (cid:107) ˆ u (cid:107) D Rout + R β (cid:12)(cid:12) ˆ u ± (cid:12)(cid:12) ,β ; D Rout ± + 1 R β (cid:13)(cid:13)(cid:13) ˆ f ± (cid:13)(cid:13)(cid:13) D Rout ± + (cid:12)(cid:12)(cid:12) ˆ f ± (cid:12)(cid:12)(cid:12) β ; D Rout ± (cid:19) , (B.9)with C = C ( n, Λ min , Λ max , C a ) independent of ˆ x ˆΓ , of J ∈ N and of y ∈ P J .To take into account the lower order terms in (B.2), we proceed as following: • we write, in (B.9), ˆ f ( y ; ˆ x ) = ˆ κ ( y ; ˆ x )ˆ u + (cid:16) ˆ ∇ · ˆ α ( y ; ˆ x ) (cid:17) · ˆ ∇ ˆ u ; • use the interpolation inequalities (cf. [55, Cor. 1.2.1]) to obtain the bounds | ˆ u ± | D Rout ± ≤ R β | ˆ u ± | ,β ; D Rout ± + R (cid:107) ˆ u (cid:107) D Rout and | ˆ u ± | ,β ; D Rout ± ≤ R | ˆ u ± | ,β ; D Rout ± + R β (cid:107) ˆ u (cid:107) D Rout ; • exploit that 0 < R ≤
1, as we have done in the proof of Lemma B.2.Summarizing the last steps, the local estimate at ˆΓ reads:
Theorem B.3.
Let the assumptions of Theorem 5.1 be fulfilled, and let us denote by C α and C κ the J - and y -independent upper bounds on the C ,βpw ( D R out ) -norm of ˆ α ( y ; · ) and on the C ,βpw ( D R out ) -norm of ˆ κ ( y , · ) , respectively. If ˆ u ∈ C ,β ˆ pw ( D R out ) is a solution to (3.3) , then, forevery ˆ x ˆΓ ∈ ˆΓ : (cid:12)(cid:12) ˆ u ± (cid:12)(cid:12) ,β ; B ± R/ (ˆ x ˆΓ ) ± ≤ C (cid:18) R β (cid:107) ˆ u (cid:107) D Rout + R β (cid:12)(cid:12) ˆ u ± (cid:12)(cid:12) ,β ; D Rout ± (cid:19) , (B.10) for a radius < R < min { , dist ( ˆ x ˆΓ , ∂D R out ) } such that B + R ( ˆ x ˆΓ ) ⊂ D out,R out ∪ ˆΓ and B − R ( ˆ x ˆΓ ) ⊂ D in ∪ ˆΓ .The constant C = C ( n, Λ min , Λ max , C α , C κ ) in (B.10) is independent of ˆ x ˆΓ , of J ∈ N andof y ∈ P J . .2 Local interior estimates Proceeding as for the local estimate at ˆΓ, it is easy to verify that analogous estimates hold inthe interior of D in and D out,R out : Theorem B.4.
Let the assumptions of Theorem 5.1 be fulfilled. If ˆ u ∈ C ,β ˆ pw ( D R out ) is a boundedsolution to (3.3) , then, for every ˆ x ∈ ˆ D in ∪ ˆ D out,R out : | ˆ u | ,β ; B R/ (ˆ x ) ≤ C (cid:18) R β (cid:107) ˆ u (cid:107) D Rout + R β (cid:12)(cid:12) ˆ u ± (cid:12)(cid:12) ,β ; D Rout ± (cid:19) , (B.11) for a radius < R < min { , dist ( ˆ x , ∂D R out ) } if ˆ x ∈ ˆ D out,R out , < R < min (cid:110) , dist (cid:16) ˆ x , ˆΓ (cid:17)(cid:111) if ˆ x ∈ ˆ D in . The constant C = C ( n, Λ min , Λ max , C α , C κ ) in (B.11) is independent of ˆ x , of J ∈ N and of y ∈ P J (with C α and C κ as in Theorem B.3).Proof. We refer to [55], Sections 6.2.2, 6.2.3, 6.2.6 and 6.3.2. We remark that the J - and y -uniform ellipticity condition on ˆ α , together with J - and y -independence of the norms of ˆ α andˆ κ , ensure the J - and y -independence of the constant C in (B.11). B.3 Global estimates
The local estimates at ∂D R out are very similar to the local estimates at ˆΓ. Therefore, wedo not present them explicitly, and refer to [25, Sect. 6.7] for details. What we obtain isthat, under the assumptions of Theorems B.3 and B.4, for every ˆ x ∂D ∈ ∂D out and 0 < R < min (cid:110) , dist (cid:16) ˆ x ∂D , ˆΓ (cid:17)(cid:111) : | ˆ u | ,β ; B − R/ (ˆ x ∂D ) ≤ C (cid:18) R β (cid:107) ˆ u (cid:107) D Rout + R β (cid:12)(cid:12) ˆ u ± (cid:12)(cid:12) ,β ; D Rout ± (cid:19) + C (cid:107) u i (cid:107) ,β ; D Rout , (B.12)with B − R/ ( ˆ x ∂D ) := B R/ ( ˆ x ∂D ) ∩ ˆ D out,R out and the constants C and C independent of ˆ x ∂D , of J ∈ N and y ∈ P J ( C is possibly depending on R ).For the global estimate, we recall that,owing to the interpolation inequalities [55, Cor. 1.2.1],in order to bound (cid:107) ˆ u ± (cid:107) ,β ; D Rout ± it is sufficient to bound | ˆ u ± | ,β ; D Rout ± and (cid:107) ˆ u (cid:107) D Rout . Usinga finite covering argument on D R out together with Theorems B.3, B.4 and equation (B.12), weobtain: Theorem B.5.
Let the assumptions of Theorem 5.1 be fulfilled. If ˆ u ∈ C ,βpw ( D R out ) is a boundedsolution to (3.3) , then (cid:13)(cid:13) ˆ u ± (cid:13)(cid:13) ,β ; D Rout ± ≤ C (cid:16) (cid:107) ˆ u (cid:107) D Rout + (cid:107) u i (cid:107) ,β ; D Rout (cid:17) , (B.13) with a constant C = C ( n, Λ min , Λ max , C α , C κ ) independent J ∈ N and of y ∈ P J (with C α and C κ the J - and y -uniform bounds on the norms of the coefficients as in Theorem B.3). To obtain the estimate on (cid:107) ˆ u ± (cid:107) k,β ; D Rout ± for k >
2, one proceeds considering the differencequotient for k = 3 and then, for k >
3, proceeds by induction. The J - and y -independenceon the constants is preserved, provided the assumptions of Theorem 5.1 are fulfilled. We referto [25, Thm. 6.17] and [25, Thm. 6.19] for details.34 emark B.6. It is clear that the regularity results reported in this section are not restrictedto the Helmholtz transmission problem. In particular, they still hold true if the elliptic operatorcontains a transport term b ( y ; ˆ x ) · ˆ ∇ ˆ u , where (cid:107) ˆ b ( y ; · ) (cid:107) C k − ,βpw ( D Rout ) is bounded independently of J and y . Indeed, the results in [25] and [55] (our guidelines throughtout this section) are statedfor an elliptic operator containing a tranport term. A nonzero right-hand side in (3.3) can betreated adding it to ˆ F in (B.2) and including it in ˆ f in (B.4) . An extension to nonhomogeneoustransmission conditions at ˆΓ is also possible: a jump in the Dirichlet trace can be treatedsimilarly to nonhomogeneous Dirichlet boundary conditions, and a jump in the Neumann tracecan be added to ˆ g in (B.3) . B.4 Proof of Lemma B.1
We present here the proof to the Schauder estimate of Lemma B.1. Schauder estimates can beproved either using Green’s representation formula for the solution to (B.6), as done in [25, Ch.6], or using Campanato norms, as in [55, Ch. 6]. Here we follow the latter approach.We consider the solution to the Poisson equation (B.6). In this section, we denote B R := B R ( ˜ x ˜Γ ), B + R := B + R ( ˜ x ˜Γ ) and B − R := B − R ( ˜ x ˜Γ ). Without loss of generality, we assume that˜Γ y = { ˜ x y ∈ R n : ˜ x y,n = 0 } , where ˜ x y,i denotes the i th component of ˜ x y , i = 1 , . . . , n . If ˜Γ y is ageneric hyperplane in R n , the estimates we will obtain still work if we substitute derivatives withrespect to the cartesian coordinates by derivatives with respect to the normal and tangentialdirections with respect to ˜Γ y , see [55, Rmk. 6.2.8]. Also, we can assume the solution ˜ u y to besufficiently smooth, see Proposition 6.2.1 in [55] (the latter still holds true if we consider thetransmission problem (B.6)).Analogously to [25, Sect. 6.7], we first assume that ˜ g y ≡
0, and only at the end return tothe general case of nonzero ˜ g y . B.4.1 Preliminaries
This subsection contains some technical lemmas that will be used in the next subsections.
Lemma B.7 (p. 174 in [55]) . For every w ∈ L (Ω) , Ω ⊂ R n , the function p ( λ ) := (cid:90) Ω ( w ( z ) − λ ) d z , λ ∈ R is strictly convex, and attains its minimum at λ = w Ω := | Ω | (cid:82) Ω w ( z ) d z .Proof. By trivial calculations one sees that d pdλ ( λ ) ≡ | Ω | > dpdλ ( λ ) = 0 ⇔ λ = w Ω . Lemma B.8 (Thm. 6.1.1 and Rmk. 6.1.2 in [55]) . Let B R ∈ R n be a ball with radius R and,for any < λ < , consider the quantity | w | ( λ ) p,µ ; B R := sup x ∈ B R , <ρ<λ diam B R (cid:32) ρ − µ (cid:90) B ρ ( x ) ∩ B R | w ( x ) − w x ,ρ | p d z (cid:33) p , (B.14) where w x ,ρ := | B ρ ( x ) ∩ B R | (cid:82) B ρ ( x ) ∩ B R w ( z ) d z .Then, for β = µ − np ∈ (0 , , (B.14) is a seminorm equivalent to the H¨older seminorm | w | β ; B R ,that is, there exist positive constants C , C , depending only on n, R, p, µ and λ , such that C | w | β ; B R ≤ | w | ( λ ) p,µ ; B R ≤ C | w | β ; B R . emma B.9 (Iteration lemma, Sect. 6.2.5 in [55]) . Assume that ψ ( R ) is a nonnegative andnondecreasing function on [0 , R ] , satisfying ψ ( ρ ) ≤ A (cid:16) ρR (cid:17) α ψ ( R ) + BR α , < ρ < R ≤ R , where α , α are constants with < α < α . Then there exists a constant C , depending onlyon A, α and α , such that ψ ( ρ ) ≤ C (cid:16) ρR (cid:17) α ( ψ ( R ) + BR α ) , < ρ < R ≤ R . In the next lemma, we consider the matrix Λ y as defined in (B.5). Lemma B.10.
For every w ∈ R n , Λ min (cid:107) w (cid:107) ≤ (cid:107) Λ y w (cid:107) ≤ Λ max (cid:107) w (cid:107) , where Λ min , Λ max > are, respectively, the J - and y -independent lower and upper bounds forthe eigenvalues of Λ y .Proof. The proof is trivial. We just remark that Λ y is well defined thanks to the assumptionthat Λ min , Λ max > Notation
Given an open domain
Ω = Ω + ∪ ˜Γ y ∪ Ω − ⊂ R n divided into two parts, Ω + and Ω − ,by ˜Γ y , and given a function h ∈ L (Ω + ) ∪ L (Ω − ) , we introduce the short notation (cid:90) Ω ± h ± ( z ) d z := (cid:90) Ω + h ( z ) | Ω + d z + (cid:90) Ω − h ( z ) | Ω − d z . Also, ( h ± ) := (cid:40) ( h | Ω + ) , in Ω + , ( h | Ω − ) , in Ω − , , and, if h ∈ H (Ω + ) ∪ H (Ω − ) , ( ˜ ∇ y h ) ± := (cid:40) ˜ ∇ y h | Ω + in Ω + , ˜ ∇ y h | Ω − in Ω − . Analogous notations with the symbol ± as exponent will follow the same rule.Furthermore, we use the symbol ˜ ∂ y,i := ˜ ∂ y ˜ ∂ y ˜ x y,i , i = 1 , . . . , n , to denote partial differentiation, ˜D for the Hessian in the ˜ x y -coordinates, and in general, ˜D j , j ∈ N , to denote the tensorcontaining all the partial derivatives of order j . B.4.2 Cacciopoli’s inequalitiesTheorem B.11.
Let ˜ u y be a solution to (B.6) in B R = B R ( ˜ x ˜Γ ) with ˜ g y ≡ . Then, for every < ρ < R and every λ ∈ R , it holds: (cid:90) B ± ρ (cid:12)(cid:12)(cid:12) ( ˜ ∇ y ˜ u y ( z )) ± (cid:12)(cid:12)(cid:12) d z ≤ C (cid:34) R − ρ ) (cid:90) B R (˜ u y ( z ) − λ ) d z + ( R − ρ ) (cid:90) B ± R ( ˜ F ± y ) ( z ) d z (cid:35) , (cid:90) B ± ρ (cid:12)(cid:12)(cid:12) ( ˜ ∇ y ˜ w y ( z )) ± (cid:12)(cid:12)(cid:12) d z ≤ C (cid:34) R − ρ ) (cid:90) B R ˜ w y ( z ) d z + (cid:90) B ± R (cid:16) ˜ F ± y ( z ) − ˜ F ± y,R (cid:17) d z (cid:35) , (B.15)(B.16) where ˜ w y := ˜ ∂ y,i ˜ u y , i = 1 , . . . , n − , B ρ = B ρ ( ˜ x ˜Γ ) and ˜ F + y,R := 1 | B + R | (cid:90) B + R ˜ F y ( z ) d z , ˜ F − y,R := 1 | B − R | (cid:90) B − R ˜ F y ( z ) d z . The constants C = C ( n, Λ min , Λ max ) , C = C ( n, Λ min , Λ max ) are independent of the center ˜ x ˜Γ of B R and B ρ , and, overall, they are independent of J ∈ N and y ∈ P J . roof. (On the lines of the proof of Thm. 6.2.2 in [55].) We first proof (B.15). Let η ∈ C ∞ ( B R )be a cut-off function such that:0 ≤ η ( ˜ x y ) ≤ , η ( ˜ x y ) = 1 in B ρ , max ˜ x y ∈ B R (cid:107) ˜ ∇ y η ( ˜ x y ) (cid:107) ≤ C η ( R − ρ )Λ max , (B.17)for some J - and y -independent constant C η > n ). Multiplying (B.6a)by η (˜ u y − λ ), integrating by parts on B + R and B − R and using (B.6b), we obtain: (cid:90) B ± R η (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z = − (cid:90) B ± R η (cid:16) Λ y ˜ ∇ y ˜ u y (cid:17) ± · (cid:16) ˜ ∇ y η (cid:17) (˜ u y − λ ) d z + (cid:90) B ± R η (˜ u y − λ ) ˜ F ± y d z . From this, applying Cauchy’s inequality on the right-hand side with ε = for the first termand ε = R − ρ ) for the second term: (cid:90) B ± R η (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ± R η (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + 2 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y η (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) (˜ u y − λ ) d z + ( R − ρ ) (cid:90) B ± R η (cid:16) ˜ F ± y (cid:17) d z + 12( R − ρ ) (cid:90) B ± R η (˜ u y − λ ) d z . Exploiting the properties (B.17) of η in the equation above and Lemma B.10, we have: (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:18) C η + 14 (cid:19) R − ρ ) (cid:90) B R (˜ u y ( z ) − λ ) d z + ( R − ρ ) (cid:90) B ± R ( ˜ F ± y ) ( z ) d z . The estimate (B.15) follows then just using Lemma B.10 to have the lower bound Λ min (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) ≤ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) on the integrand at the left-hand side. It is clear then that C in (B.15) dependsonly on n , Λ min and Λ max and it does not depend on the center ˜ x ˜Γ of B R .To prove (B.16), we differentiate (B.6) by ˜ x y,i , i = 1 , . . . , n −
1. Then ˜ w y = ˜ ∂ y,i ˜ u y , i =1 , . . . , n −
1, satisfies − ˜ ∇ y · (cid:16) Λ y ˜ ∇ y ˜ w y (cid:17) = ˜ ∂ y,i (cid:16) ˜ F y − ˜ F ± y,R (cid:17) , in B + R ∪ B − R , (cid:74) ˜ w y (cid:75) ˜Γ y = 0 , (cid:114) Λ y ˜ ∂ y ˜ w y ˜ ∂ y ˜ n y (cid:122) ˜Γ y = 0 . (B.18a)(B.18b)Then, multiplying (B.18a) by η ˜ w y and integrating by parts: (cid:90) B ± R η (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ w y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z = − (cid:90) B ± R η ˜ w y (cid:16) Λ y ˜ ∇ y ˜ w y (cid:17) ± · ˜ ∇ y η d z + (cid:90) B ± R ˜ ∂ y,i (cid:16) ˜ F ± y − ˜ F ± y,R (cid:17) η ˜ w y d z . If we integrate by parts (on B + R and B − R separately) the last term on the right-hand side,then (cid:74) η ˜ w y (cid:16) ˜ F ± y − ˜ F ± y,R (cid:17) ˜ n y,i (cid:75) ˜Γ y = 0 for i = 1 , . . . , n −
1, due to the fact that the tangentialcomponents (with respect to ˜Γ y ) of ˜ n y are zero. Thus (cid:90) B ± R η (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ w y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z = − (cid:90) B ± R η ˜ w y (cid:16) Λ y ˜ ∇ y ˜ w y (cid:17) ± · ˜ ∇ y η d z − (cid:90) B ± R ˜ ∂ y,i ( η ˜ w y ) (cid:16) ˜ F ± y − ˜ F ± y,R (cid:17) d z . Splitting the derivative in the last integral, and using Cauchy’s inequality on each term with ε = ε (Λ min ) sufficiently small, the properties (B.17) of η lead to (B.16). As for the constant in(B.15), it is clear that C = C ( n, Λ min , Λ max ) in (B.16) is independent of the center ˜ x ˜Γ of B R and in general of J ∈ N and y ∈ P J . 37 emark B.12 (Analogous to Rmk. 6.2.5 in [55]) . Since (B.6a) can be rewritten as (Λ y ) nn ˜ ∂ y,nn ˜ u y = − (cid:80) n − j =1 (Λ y ) jj ˜ ∂ y,jj ˜ u y − ˜ F y , we obtain, using (B.16) : (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:32) n − (cid:88) j =1 (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y (cid:16) ˜ ∂ y,j ˜ u y (cid:17)(cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + (cid:90) B ± ρ (cid:16) ˜ F ± y (cid:17) d z (cid:33) ≤ C (cid:48) (cid:34) R − ρ ) (cid:90) B R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) ( z ) d z + (cid:90) B ± R (cid:16) ˜ F ± y ( z ) − ˜ F ± y,R (cid:17) d z + (cid:90) B ± R (cid:16) ˜ F ± y (cid:17) ( z ) d z (cid:35) , (B.19) where, thanks to Lemma B.10 and Theorem B.11, C and C (cid:48) depend only on n , Λ min and Λ max . If we apply (B.19) with ρ = R and R , and (B.15) with ρ = R and R (and λ = 0), weobtain [55, Cor. 6.2.4]: Corollary B.13.
Let ˜ u y be a solution to (B.6) in B R = B R ( ˜ x ˜Γ ) with ˜ g y ≡ . Then it holds: (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:20) R (cid:90) B R ˜ u y ( z ) d z + R n (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R + R n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (cid:21) , (B.20) where C = C ( n, Λ min , Λ max ) is independent of the center ˜ x ˜Γ of B R and, overall, is independentof J ∈ N and y ∈ P J . Corollary B.14.
Let ˜ u y be a solution to (B.6) in B = B ( ˜ x ˜Γ ) . If ˜ F y ≡ and ˜ g y ≡ , then,for any positive integer k ∈ N : (cid:107) ˜ u y (cid:107) H k ( B +1 / ) + (cid:107) ˜ u y (cid:107) H k ( B − / ) ≤ C (cid:16) (cid:107) ˜ u y (cid:107) L ( B +1 / ) + (cid:107) ˜ u y (cid:107) L ( B − / ) (cid:17) , (B.21) with C = C ( n, Λ min , Λ max , k ) independent of the center ˜ x ˜Γ of B and overall of J ∈ N and y ∈ P J .Proof. (On the lines of the proof of Corollary 6.2.5 in [55].) For k = 1, the claim follows directlyfrom Cacciopoli’s inequality (B.15) with λ = 0. For k = 2, it follows from (B.19). For k > Corollary B.15.
Let ˜ u y be a solution to (B.6) in B = B ( ˜ x ˜Γ ) . If ˜ F y ≡ in and ˜ g y ≡ , then sup B + R/ ∪ B − R/ | ˜ u y | ≤ C (cid:32) R n (cid:90) B + R ∪ B − R ˜ u y ( z ) d z (cid:33) , (B.22) with C = C ( n, Λ min , Λ max ) independent of the center ˜ x ˜Γ of B R and overall of J ∈ N and y ∈ P J .Proof. (Analogous to proof of Corollary 6.2.6 in [55]) We first establish the estimate assuming R = 1. The Sobolev embedding theorem [25, Thm. 7.26] applied in B +1 / and in B − / impliesthat, for k > n , sup B +1 / ∪ B − / | ˜ u y | ≤ C (cid:16) (cid:107) ˜ u y (cid:107) H k ( B +1 / ) + (cid:107) ˜ u y (cid:107) H k ( B − / ) (cid:17) , where the constant C depends on n only. The claim for R = 1 follows then from CorollaryB.14.For a generic radius R , the result follows from a scaling argument, defining v ( z ) := ˜ u y ( R z ), z ∈ B , and applying the estimate for R = 1 to the function v .38 .4.3 Interface estimate for the Laplace equationTheorem B.16. Let ˜ u y be a solution to (B.6) , with ˜ F y ≡ and ˜ g y ≡ . Then, for every < ρ ≤ R and every i ∈ N : (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D i ˜ u y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:16) ρR (cid:17) n (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D i ˜ u y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z , (B.23) where C = C ( n, Λ min , Λ max ) is independent of the center ˜ x ˜Γ of B R and overall of J ∈ N and y ∈ P J .Proof. (On the lines of the proof of Thm. 6.2.4 in [55].) Case i = 0 . If 0 < ρ < R , then (cid:90) B ρ (˜ u ± y ) d z ≤ | B ρ | sup B + ρ ∪ B − ρ ˜ u y = C n ρ n sup B + ρ ∪ B − ρ ˜ u y ≤ C (cid:16) ρR (cid:17) n (cid:90) B ρ (˜ u ± y ) d z , where the last constant is given by the product of C n with the constant from Corollary B.15,and C n is a constant depending on n only.If R ≤ ρ ≤ R , then, trivially, (cid:90) B ρ (˜ u ± y ) d z ≤ (cid:90) B R (˜ u ± y ) d z ≤ n (cid:16) ρR (cid:17) n (cid:90) B R (˜ u ± y ) d z . Case i = 1 . We first consider 0 < ρ < R . For k − > n , the Sobolev embedding theorem [25, Thm.7.26] ensures that (cid:107) ˜ ∇ y ˜ u y (cid:107) L ∞ ( B + R/ ) ≤ C k (cid:107) ˜ ∇ y ˜ u y (cid:107) H k − ( B + R/ ) , for a constant C k dependent on n only (and analogously in B − R/ ). Exploiting this fact, we have: (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C n ρ n sup B + R/ ∪ B − R/ (cid:12)(cid:12)(cid:12) ˜ ∇ y ˜ u y (cid:12)(cid:12)(cid:12) ≤ C n C k ρ n k − (cid:88) l =0 R l − n (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D l ( ˜ ∇ y ˜ u y ) (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z = C n C k ρ n k (cid:88) j =1 R j − − n (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D j ˜ u y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z , where C n is a constant depending only on n . The factors R l − n , l = 0 , . . . , k −
1, in the secondinequality are due to a scaling argument as in Corollary B.15. We note that k depends on n only. Denoting ˜ u y,R := | B R | (cid:82) B R ˜ u y ( z ) d z , and observing that (˜ u y − ˜ u y,R ) fulfills (B.6) with˜ F y ≡ g y ≡
0, we can apply Corollary B.14 and derive: (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ Cρ n k (cid:88) j =1 R j − − n R − j (cid:90) B ± R (cid:12)(cid:12) ˜ u ± y − ˜ u y,R (cid:12)(cid:12) d z = kC (cid:16) ρR (cid:17) n R − (cid:90) B ± R (cid:12)(cid:12) ˜ u ± y − ˜ u y,R (cid:12)(cid:12) d z ≤ kCC p (cid:16) ρR (cid:17) n R − R (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z . C in the previous inequalities is the product of C n C k by the constant of CorollaryB.14. The constant C p is instead a scalar factor, independent of R , coming from application ofthe Poincar´e inequality for balls (that we could apply being ˜ u y in H ( B R )).For R ≤ ρ ≤ R , the inequality (B.23) follows trivially taking C ≥ n . Case i = 2 . For j = 1 , . . . , n − (cid:74) ˜ ∂ y,j ˜ u y (cid:75) ˜Γ y = 0 and (cid:114) Λ y ˜ ∂ y ( ˜ ∂ y,j ˜ u y )˜ ∂ y ˜ n y (cid:122) ˜Γ y = 0. Thus, the case for i = 1implies that (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ( ˜ ∂ y,j ˜ u y ) (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:16) ρR (cid:17) n (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ( ˜ ∂ y,j ˜ u y ) (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z , (B.24)with C = C ( n, Λ min , Λ max ). For j = n , we can use that (Λ y ) nn ˜ ∂ y,nn ˜ u y = − (cid:80) n − j =1 (Λ y ) jj ˜ ∂ y,jj ˜ u y to obtain: (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,nn ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C n Λ max Λ min n − (cid:88) j =1 (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,jj ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ CC n Λ max Λ min (cid:16) ρR (cid:17) n n − (cid:88) j =1 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ( ˜ ∂ y,j ˜ u y ) (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ CC n Λ max Λ min (cid:16) ρR (cid:17) n (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z , where the constant C n depends on n only, and the constant C = C ( n, Λ min , Λ max ) is theconstant in (B.24). The latter inequality together with (B.24) imply finally (B.23). Case i = 3 . The result follows similarly as for the case i = 2: the estimates associated to ˜D ( ˜ ∂ y,j ˜ u y ), for j = 1 , . . . , n −
1, follow from the case i = 2; for ˜ ∂ y,nnn ˜ u y , we observe that, by differentiation,(B.6a) implies that (Λ y ) nn ˜ ∂ y,nnn ˜ u y = − (cid:80) n − j =1 (Λ y ) jj ˜ ∂ y,jjn ˜ u y , and we can proceed as we did inthe case i = 2 for ˜ ∂ y,nn ˜ u y .The case i > i = 3. Theorem B.17.
Let ˜ u y be a solution to (B.6) , with ˜ F y ≡ and ˜ g y ≡ . Then, for every < ρ ≤ R : (cid:90) B ± ρ (cid:0) ˜ u ± y ( z ) − ˜ u ± y,ρ (cid:1) d z ≤ C (cid:16) ρR (cid:17) n +2 (cid:90) B R ˜ u y ( z ) d z , (B.25) where ˜ u + y,ρ = | B + ρ | (cid:82) B + ρ ˜ u y ( z ) d z and ˜ u − y,ρ = | B − ρ | (cid:82) B − ρ ˜ u y ( z ) d z . The constant C = C ( n, Λ min , Λ max ) is independent of the center ˜ x ˜Γ of B R and, overall, of J ∈ N and y ∈ P J .Proof. (On the lines of proof of Thm. 6.2.5 in [55].) If 0 < ρ < R : (cid:90) B ± ρ (cid:0) ˜ u ± y − ˜ u ± y,ρ (cid:1) d z ≤ C p ρ (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ CC p ρ (cid:16) ρR (cid:17) n (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48) ρ (cid:16) ρR (cid:17) n R (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z . In the first step we have applied the Poincar´e inequality in B + ρ and B − ρ separately, and denotedby C p the ρ -independent scalar factor in the Poincar´e constant. For the second inequality,40e have used Theorem B.16, and the constant C corresponds to the constant in (B.23). Fi-nally, the last line follows from Cacciopoli’s inequality (B.15), and we have denoted by C (cid:48) themultiplication of CC p with the constant in (B.15).If R ≤ ρ ≤ R , (B.25) follows simply taking C ≥ n +2 . B.4.4 Interface estimate for the Poisson equationTheorem B.18.
Let ˜ u y be a solution to (B.6) in B R , with ˜ g y ≡ , and let ˜ w y = ˜ ∂ y,i ˜ u y , i = 1 , . . . , n − . Then, for any < ρ ≤ R ≤ R : ρ n +2 β (cid:90) B ± ρ n (cid:88) j =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C R n +2 β (cid:90) B ± R n − (cid:88) j =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z + C R n +2 β (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± R (cid:12)(cid:12)(cid:12)(cid:12) d z + C (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R , where (cid:16) ˜ ∂ y,j ˜ w y (cid:17) + ρ := | B + ρ | (cid:82) B + ρ ˜ ∂ y,j ˜ w y d z and (cid:16) ˜ ∂ y,j ˜ w y (cid:17) − ρ := | B − ρ | (cid:82) B − ρ ˜ ∂ y,j ˜ w y d z , and similarly for (cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± R . The constants C = C ( n, Λ min , Λ max , β ) and C = C ( n, Λ min , Λ max ) are indepen-dent of the center ˜ x ˜Γ of B R and, overall, they are independent of J ∈ N and y ∈ P J .Proof. (On the lines of the proof of Thm. 6.2.9 in [55].) We decompose ˜ w y as ˜ w y = ˜ w (cid:48) y + ˜ w (cid:48)(cid:48) y ,where − ˜ ∇ y · (cid:16) Λ y ˜ ∇ y ˜ w (cid:48) y (cid:17) = 0 , in B + R ∪ B − R , (cid:74) ˜ w (cid:48) y (cid:75) ˜Γ y = 0 , (cid:114) Λ y ˜ ∂ y ˜ w (cid:48) y ˜ ∂ y ˜ n y (cid:122) ˜Γ y = 0 , ˜ w (cid:48) y | ∂B R = ˜ w y , − ˜ ∇ y · (cid:16) Λ y ˜ ∇ y ˜ w (cid:48)(cid:48) y (cid:17) = ˜ ∂ y,i (cid:16) ˜ F y − ˜ F ± y,R (cid:17) , in B + R ∪ B − R , (cid:74) ˜ w (cid:48)(cid:48) y (cid:75) ˜Γ y = 0 , (cid:114) Λ y ˜ ∂ y ˜ w (cid:48)(cid:48) y ˜ ∂ y ˜ n y (cid:122) ˜Γ y = 0 , ˜ w (cid:48)(cid:48) y | ∂B R = 0 , and ˜ F ± y,R = ˜ F + y = | B + R | (cid:82) B + R ˜ F y ( z ) d z in B + R , ˜ F − y = | B − R | (cid:82) B − R ˜ F y ( z ) d z in B − R . We first consider j = 1 , . . . , n −
1. In this case, ˜ ∂ y,j ˜ w (cid:48) y solves (B.6) with ˜ F y ≡ g y ≡ (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48) y (cid:17) ± − (cid:16) ˜ ∂ y,j ˜ w (cid:48) y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:16) ρR (cid:17) n +2 (cid:90) B R (cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y j ˜ w (cid:48) y (cid:17)(cid:12)(cid:12)(cid:12) d z , (B.26)where C = C ( n, Λ min , Λ max ) is the constant in (B.25). Then, for ˜ ∂ y,j ˜ w y , using Lemma B.7 wecan write: (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± − (cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + 2 (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48) (cid:16) ρR (cid:17) n +2 (cid:90) B R (cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17)(cid:12)(cid:12)(cid:12) d z + C (cid:48)(cid:48) (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z , (B.27)41or C (cid:48) , C (cid:48)(cid:48) independent of ˜ x ˜Γ , J and y . The last summand in the above inequality can bebounded as (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ min (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) Λ y ˜ ∇ y ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z = − min (cid:90) B ± R (cid:16) ˜ ∇ y · (cid:16) Λ y ˜ ∇ y ˜ w (cid:48)(cid:48) y (cid:17) ˜ w (cid:48)(cid:48) y (cid:17) ± d z = 1Λ min (cid:90) B ± R ˜ ∂ y,i (cid:16) ˜ F ± y − ˜ F ± y,R (cid:17) (cid:0) ˜ w (cid:48)(cid:48) y (cid:1) ± d z = − min (cid:90) B ± R (cid:16) ˜ F ± y − ˜ F ± y,R (cid:17) (cid:16) ˜ ∂ y,i ˜ w (cid:48)(cid:48) y (cid:17) ± d z ≤ min (cid:32) ε (cid:90) B ± R (cid:12)(cid:12)(cid:12) ˜ F ± y − ˜ F ± y,R (cid:12)(cid:12)(cid:12) d z + ε (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z (cid:33) , where for the third and fifth line we have used integration by parts, and the jump terms on ˜Γvanished because of the transmission conditions in the first case, and because of null tangentialcomponents of ˜ n y in the second case; the boundary terms vanished because of the Dirichletboundary conditions. In the last step, we have applied Cauchy’s inequality for a generic ε > ε = ε (Λ min ) sufficiently small, we finally obtain (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48)(cid:48)(cid:48) R n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R , (B.28)for a positive constant C (cid:48)(cid:48)(cid:48) = C (cid:48)(cid:48)(cid:48) (Λ min ). Combining the last estimate with (B.27), we infer (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48) (cid:16) ρR (cid:17) n +2 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z + ( C (cid:48)(cid:48) C (cid:48)(cid:48)(cid:48) ) R n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R , (B.29)(B.30)for j = 1 , . . . , n − j = n and 0 < ρ < R . Applying the Poincar´e inequality, Lemma B.7,Theorem B.16 and equation (B.28) (which holds for j = n , too): (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± − (cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w (cid:48) y (cid:17) ± − (cid:16) ˜ ∂ y,n ˜ w (cid:48) y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z + 2 (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w (cid:48)(cid:48) y (cid:17) ± − (cid:16) ˜ ∂ y,n ˜ w (cid:48)(cid:48) y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C p ρ (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y (cid:16) ˜ ∂ y,n ˜ w (cid:48) y (cid:17)(cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + 2 (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C p Cρ (cid:16) ρR (cid:17) n (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ w (cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + C (cid:48)(cid:48)(cid:48) R n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R , with C p a scalar factor, independent of R , coming from application of the Poincar´e inequality,42nd C is the constant in Theorem B.16. We can bound the integral on the right-hand side by: (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ w (cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,nn ˜ w (cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + 2 n − (cid:88) j =1 (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ( ˜ ∂ y,j ˜ w (cid:48) y ) (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C n − (cid:88) j =1 (cid:90) B ± R/ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ( ˜ ∂ y,j ˜ w (cid:48) y ) (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ CC c R n − (cid:88) j =1 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ CC c R (cid:32) n − (cid:88) j =1 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + n − (cid:88) j =1 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w (cid:48)(cid:48) y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z (cid:33) ≤ CC c R (cid:32) n − (cid:88) j =1 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + ( n − C (cid:48)(cid:48)(cid:48) R n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (cid:33) . In the second line we have used the equality (Λ y ) nn ˜ ∂ y,nn ˜ w (cid:48) y = − (cid:80) n − j =1 (Λ y ) jj ˜ ∂ y,jj ˜ w (cid:48) y and LemmaB.10, and thus C = C ( n, Λ min , Λ max ) is J - and y -independent. In the third line we have ex-ploited the Cacciopoli inequality (B.15), and, in the last step, the bound (B.28). Summarizing,for j = n and 0 < ρ < R we have (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± − (cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48) (cid:16) ρR (cid:17) n +2 n − (cid:88) j =1 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) d z + C (cid:48) R n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R , for two constants C (cid:48) = C (cid:48) ( n, Λ min , Λ max ) and C (cid:48) = C (cid:48) ( n, Λ min , Λ max ) independent of the centerof B R , of J ∈ N and of y ∈ P J .Combining the latter estimate with the estimate (B.29), we finally obtain: (cid:90) B ± ρ n (cid:88) j =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:16) ρR (cid:17) n +2 (cid:90) B ± R n − (cid:88) j =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z + C (cid:16) ρR (cid:17) n +2 (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± R (cid:12)(cid:12)(cid:12)(cid:12) d z + C R n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R , for two positive constants C = C ( n, Λ min , Λ max ) and C = C ( n, Λ min , Λ max ). The claim for0 < ρ < R follows then by application of Lemma B.7 and the Iteration Lemma B.9.If instead R ≤ ρ ≤ R , the claim holds simply by choosing C ≥ n +2 and using LemmaB.7. Theorem B.19.
Let ˜ u y be a solution to (B.6) in B R , with ˜ g y ≡ , and let ˜ w y = ˜ ∂ y,i ˜ u y , i = 1 , . . . , n . Then, for any < ρ ≤ R : (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∇ y ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ Cρ nβ M R , (B.31) with M R = 1 R β (cid:107) ˜ u y (cid:107) B R + 1 R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R , (B.32) and C = C ( n, Λ min , Λ max , β ) independent of the center ˜ x ˜Γ of B R and B ρ , and, overall, of J ∈ N and y ∈ P J . The term (cid:16) ˜ ∇ y ˜ w y (cid:17) ± ρ has the same meaning as in Theorem B.18. roof. (On the lines of the proof of Thm. 6.2.10 in [55].) For i = 1 , . . . , n −
1, we can applyTheorem B.18, which, together with Lemma B.7 and Corollary B.13, brings: (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∇ y ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∇ y ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z = (cid:90) B ± ρ n (cid:88) j =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48) ρ n +2 β (cid:32) R n +2 β (cid:90) B ± R/ n − (cid:88) j =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z + 1 R n +2 β (cid:90) B ± R (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,n ˜ w y (cid:17) ± R (cid:12)(cid:12)(cid:12)(cid:12) d z + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (cid:19) ≤ C (cid:48) ρ n +2 β (cid:32) R n +2 β (cid:90) B ± R/ n (cid:88) j =1 (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,j ˜ w y (cid:17) ± ( z ) (cid:12)(cid:12)(cid:12)(cid:12) d z + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (cid:33) ≤ C (cid:48)(cid:48) ρ n +2 β (cid:18) R n +2 β +4 (cid:90) B R ˜ u y d z + 1 R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (cid:19) ≤ C (cid:48)(cid:48) ρ n +2 β (cid:18) R β +4 (cid:107) ˜ u y (cid:107) R + 1 R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (cid:19) , with C (cid:48) and C (cid:48)(cid:48) depending only on n , Λ min , Λ max and β .For i = n , we observe that, from (B.6a):(Λ y ) nn ˜ ∂ y,nn ˜ u y + n − (cid:88) j =1 (Λ y ) jj (cid:16) ˜ ∂ y,jj ˜ u y (cid:17) ± ρ + ˜ F ± y,ρ = − n − (cid:88) j =1 (Λ y ) jj (cid:18) ˜ ∂ y,jj ˜ u y − (cid:16) ˜ ∂ y,jj ˜ u y (cid:17) ± ρ (cid:19) − (cid:16) ˜ F y − ˜ F ± y,ρ (cid:17) in B + R ∪ B − R , where (cid:16) ˜ ∂ y,jj ˜ u y (cid:17) ± ρ denotes (cid:16) ˜ ∂ y,jj ˜ u y (cid:17) + ρ = | B + ρ | (cid:82) B + ρ ˜ ∂ y,jj ˜ u y d z in B + ρ and (cid:16) ˜ ∂ y,jj ˜ u y (cid:17) − ρ = | B − ρ | (cid:82) B − ρ ˜ ∂ y,jj ˜ u y d z in B − ρ , for j = 1 , . . . , n . Analogous notation holds for ˜ F ± y,ρ . Setting λ ± := − (cid:80) n − j =1 (cid:16) ˜ ∂ y,jj ˜ u y (cid:17) ± ρ − ˜ F ± y,ρ , we can apply Lemma B.7 and the result for j = 1 , . . . , n − (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,nn ˜ w y (cid:17) ± ( z ) − (cid:16) ˜ ∂ y,nn ˜ w y (cid:17) ± ρ (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ± ρ (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ ∂ y,nn ˜ w y (cid:17) ± ( z ) − λ ± (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48)(cid:48)(cid:48) ρ n +2 β M R + (cid:90) B ± ρ (cid:12)(cid:12)(cid:12) ˜ F ± y − ˜ F ± y,ρ (cid:12)(cid:12)(cid:12) d z ≤ C (cid:48)(cid:48)(cid:48) ρ n +2 β M R + ρ n +2 β (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± ρ ≤ ( C (cid:48)(cid:48)(cid:48) + 1) ρ n +2 β M R , with C (cid:48)(cid:48)(cid:48) = C (cid:48)(cid:48)(cid:48) ( n, Λ min , Λ max , β ). Theorem B.20.
Let ˜ u y be a solution to (B.6) in B R = B R ( ˜ x ˜Γ ) , with ˜ g ≡ . Then (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) β ; B ± R/ (˜ x ˜Γ ) ≤ C (cid:18) R β (cid:107) ˜ u y (cid:107) B R (˜ x ˜Γ ) + 1 R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (˜ x ˜Γ ) (cid:19) , (B.33) where C = C ( n, Λ min , Λ max , β ) is independent of the center ˜ x ˜Γ of B R , and overall it is inde-pendent of J ∈ N and y ∈ P J . roof. (On the lines of the proof of Thm. 6.2.11 in [55].) For ˜ x y ∈ ˜Γ y ∩ B R/ ( ˜ x ˜Γ ) and 0 < ρ ≤ R ,Theorem B.19 implies: (cid:90) B ± ρ (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± − (cid:16) ˜D ˜ u y (cid:17) ± B ρ (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48) ρ n +2 β M R/ ≤ C (cid:48) ρ n +2 β M R , (B.34)with C (cid:48) = C (cid:48) ( n, Λ min , Λ max , β ) (where again ( · ) ± B ± ρ (˜ x y ) denotes the mean on B + ρ and B − ρ ).If instead ˜ x y ∈ B R/ ( ˜ x ˜Γ ) but ˜ x y / ∈ ˜Γ y , we denote ˜ x (cid:48) y := (˜ x y, , . . . , ˜ x y,n − , ∈ ˜Γ y , anddistinguish two cases: 0 < ˜ x y,n < R and R ≤ ˜ x y,n < R . In the first case, we consider twosubcases: ˜ x y,n ≤ ρ ≤ R and 0 < ρ < ˜ x y,n .If 0 < ˜ x y,n < R and ˜ x y,n ≤ ρ ≤ R , then B ρ ( ˜ x y ) ∩ B R/ ( ˜ x ˜Γ ) ⊂ B ρ ( ˜ x (cid:48) y ), and, using (B.34),we can write: (cid:90) ( B ρ (˜ x y ) ∩ B R/ (˜ x ˜Γ ) ) ± (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± − (cid:16) ˜D ˜ u y (cid:17) ± B ρ (˜ x y ) ∩ B R/ (˜ x ˜Γ ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ρ (˜ x (cid:48) y ) ± (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± − (cid:16) ˜D ˜ u y (cid:17) ± B ρ (˜ x (cid:48) y ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48) ρ n +2 β M R . (B.35)If 0 < ˜ x y,n < R and 0 < ρ < ˜ x y,n , then B ρ ( ˜ x y ) is in the interior, meaning it does not crossthe interface ˜Γ y . Therefore, using the analogue of Theorem B.18 for the interior [55, Thm.6.2.6] (where it can be checked, as for the interface case, that the constants in the bound areindependent of the center of the ball, of J ∈ N and y ∈ P J ), we have: (cid:90) B ρ (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12) ˜D ˜ u y − (cid:16) ˜D ˜ u y (cid:17) B ρ (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C (cid:48)(cid:48) ρ n +2 β (cid:32) x n +2 βy,n (cid:90) B ˜ xy,n (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) − (cid:16) ˜D ˜ u y (cid:17) B ˜ xy,n (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12) d z + (cid:12)(cid:12)(cid:12) ˜ F y (cid:12)(cid:12)(cid:12) β ; B ˜ xy,n (˜ x y ) (cid:33) ≤ C (cid:48) C (cid:48)(cid:48) ρ n +2 β M R , with C (cid:48)(cid:48) = C (cid:48)(cid:48) ( n, Λ min , Λ max , β ). In the last step we have used (B.35) with ρ = ˜ x y,n (as B ˜ x y,n ( ˜ x y ) ⊂ B x y,n ( ˜ x (cid:48) y )).Altogether, if 0 < ˜ x y,n < R , then, for 0 < ρ < R : (cid:90) ( B ρ (˜ x y ) ∩ B R/ (˜ x ˜Γ ) ) ± (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± − (cid:16) ˜D ˜ u y (cid:17) ± B ρ (˜ x y ) ∩ B R/ (˜ x ˜Γ ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C ρ n +2 β M R , (B.36)with C = C ( n, Λ min , Λ max , β ).If instead R ≤ ˜ x y,n < R and 0 < ρ ≤ R , then either B ρ ( ˜ x y ) ⊂ B R/ ( ˜ x y ) ⊂ B +3 R/ ( ˜ x ˜Γ ) ⊂ B + R ( ˜ x ˜Γ ), or B ρ ( ˜ x y ) ⊂ B R/ ( ˜ x y ) ⊂ B − R/ ( ˜ x ˜Γ ) ⊂ B − R ( ˜ x ˜Γ ). In the first case, the analogous ofTheorem B.19 for the interior (where again it can be checked that the constants in the boundare independent of the center of the ball, of J ∈ N and y ∈ P J ) implies: (cid:90) B ρ (˜ x y ) ∩ B + R/ (˜ x ˜Γ ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ˜D ˜ u y − (cid:16) ˜D ˜ u y (cid:17) B ρ (˜ x y ) ∩ B + R/ (˜ x ˜Γ ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) d z ≤ (cid:90) B ρ (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12) ˜D ˜ u y − (cid:16) ˜D ˜ u y (cid:17) B ρ (˜ x y ) (cid:12)(cid:12)(cid:12)(cid:12) d z ≤ C ρ n +2 β M R , with C = C ( n, Λ min , Λ max , β ). An analogous estimate holds in the second case. Consideringthis last bound together with (B.36), and using Lemma B.8 with µ = n + 2 β , p = 2 and λ = R ,45e finally obtain, for C = C ( n, β ): (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ D ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) β ; B ± R/ (˜ x ˜Γ ) ≤ C (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ D ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) (cid:16) (cid:17) ,n +2 β ; B ± R/ (˜ x ˜Γ ) ≤ CM R , from which (B.33) follows, with C = C ( n, Λ min , Λ max , β ).The previous result states the local estimate for the Poisson equation in case of homogeneoustransmission conditions. We are now in the position to consider the case ˜ g (cid:54) = 0: Theorem B.21.
Let ˜ u y be a solution to (B.6) in B R = B R ( ˜ x ˜Γ ) . Then (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) β ; B ± R/ (˜ x ˜Γ ) ≤ C (cid:18) R β (cid:107) ˜ u y (cid:107) B R (˜ x ˜Γ ) + 1 R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) + (cid:12)(cid:12)(cid:12) ˜ F y (cid:12)(cid:12)(cid:12) β ; B ± R (˜ x ˜Γ ) + 1 R β (cid:107) ˜ g y (cid:107) y + | ˜ g y | β ;˜Γ y (cid:19) , where C = C ( n, Λ min , Λ max , β ) is independent of the center ˜ x ˜Γ of B R , of J ∈ N and y ∈ P J .Proof. For this proof we use a similar argument as in [25, pp. 124–125].Consider a nonnegative function η ∈ C ( R n − ), such that (cid:82) R n − η ( z (cid:48) ) d z (cid:48) = 1. Lemma 6.38in [25] ensures that ˜ g y can be extended outside ˜Γ y in such a way that its extension belongs to C ,β ( R n − ). With some abuse of notation, we still denote by ˜ g y this extension.We define ˜ ψ and ˜ ψ as the functions fulfilling the following equalities:(Λ y ) + nn ˜ ψ ( ˜ x y ) = 12 ˜ x y,n (cid:90) R n − ˜ g y ( ˜ x (cid:48) y − ˜ x y,n z (cid:48) ) η ( z (cid:48) ) d z (cid:48) , (Λ y ) − nn ˜ ψ ( ˜ x y ) = −
12 ˜ x y,n (cid:90) R n − ˜ g y ( ˜ x (cid:48) y − ˜ x y,n z (cid:48) ) η ( z (cid:48) ) d z (cid:48) , (B.37)(B.38)where ˜ x (cid:48) y = (˜ x (cid:48) y, , . . . , ˜ x (cid:48) y,n − ). It can be checked (see (B.41)) that ˜ ψ , ˜ ψ ∈ C ,β ( R n ), and that:˜ ψ ( ˜ x (cid:48) y ,
0) = ˜ ψ ( ˜ x (cid:48) y ,
0) = 0 , (Λ y ) + nn ˜ ∂ y ˜ ∂ y ˜ x y,n ˜ ψ ( ˜ x (cid:48) y , − (Λ y ) − nn ˜ ∂ y ˜ ∂ y ˜ x y,n ˜ ψ ( ˜ x (cid:48) y ,
0) = ˜ g y ( ˜ x (cid:48) y ) . The solution ˜ u y to (B.6) can be decomposed as ˜ u y = ˜ v y + ˜ ψ y , where˜ ψ y | B + R = ˜ ψ , ˜ ψ y | B − R = ˜ ψ . Then ˜ v y fulfills − ˜ ∇ y · (cid:16) Λ y ˜ ∇ y ˜ v y (cid:17) = ˜ F y + ˜ ∇ y · (cid:16) Λ y ˜ ∇ y ˜ ψ y (cid:17) , in B + R ∪ B − R , (cid:74) ˜ v y (cid:75) ˜Γ y = 0 , (cid:114) Λ y ˜ ∂ y ˜ v y ˜ ∂ y ˜ n y (cid:122) ˜Γ y = 0 . Applying Theorem B.20 to ˜ v y (and with the help of Lemma B.10), we infer: (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜D ˜ v y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) β ; B ± R/ (˜ x ˜Γ ) ≤ C (cid:48) R β (cid:107) ˜ v y (cid:107) B R (˜ x ˜Γ ) + C (cid:48) R β (cid:18)(cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) + (cid:13)(cid:13)(cid:13) ˜ ψ y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) (cid:19) + C (cid:48) (cid:18)(cid:12)(cid:12)(cid:12) ˜ F y (cid:12)(cid:12)(cid:12) β ; B ± R (˜ x ˜Γ ) + (cid:12)(cid:12)(cid:12) ˜ ψ y (cid:12)(cid:12)(cid:12) ,β ; B ± R (˜ x ˜Γ ) (cid:19) , (B.40)46ith C (cid:48) = C (cid:48) ( n, Λ min , Λ max , β ).Denoting by ˜ ∇ y, ˜Γ the gradient with respect to the first n − ∂ y,i the derivative with respect to the i th component, in B + R ( ˜ x ˜Γ ) we have:(Λ y ) + nn ˜ ∂ y,ij ˜ ψ y = 12 (cid:90) R n − ˜ ∂ y,i ˜ g y ( ˜ x (cid:48) y − ˜ x y,n z (cid:48) ) ˜ ∂ y,j η ( z (cid:48) ) d z (cid:48) , for i, j (cid:54) = n, (Λ y ) + nn ˜ ∂ y,in ˜ ψ y = − (cid:90) R n − z (cid:48) · ˜ ∇ y, ˜Γ ˜ g y ( ˜ x (cid:48) y − ˜ x y,n z (cid:48) ) ˜ ∂ y,i η ( z (cid:48) ) d z (cid:48) , for i (cid:54) = n, (Λ y ) + nn ˜ ∂ y,nn ˜ ψ y = 12 (cid:90) R n − z (cid:48) · ˜ ∇ y, ˜Γ ˜ g y ( ˜ x (cid:48) y − ˜ x y,n z (cid:48) ) (cid:104) ( n − η ( z (cid:48) ) + z (cid:48) · ˜ ∇ y, ˜Γ η ( z (cid:48) ) (cid:105) d z (cid:48) , and thus (cid:13)(cid:13)(cid:13) D ˜ ψ y (cid:13)(cid:13)(cid:13) B + R (˜ x ˜Γ ) ≤ C η | ˜ g y | y , (cid:12)(cid:12)(cid:12) ˜ ψ y (cid:12)(cid:12)(cid:12) ,β ; B + R (˜ x ˜Γ ) ≤ C η | ˜ g y | ,β ;˜Γ y . (B.41a)(B.41b)Analogous results hold for the norms on B − R ( ˜ x ˜Γ ). The constant C η depends on the norms of η on ˜Γ y , and thus, in principle, it could depend on y ∈ P J and J ∈ N . However, if, for every˜ x ˜Γ ∈ ˜Γ y considered, we use, in B R ( ˜ x ˜Γ ), the same function η translated so that it is centeredin ˜ x ˜Γ , then C η is independent of J ∈ N and of y ∈ P J . Combining (B.41) with (B.40), andusing the interpolation inequalities (cf. [55, Cor. 1.2.1]) to bound (cid:107) ˜ g (cid:107) y , we gather the desiredestimate: (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) ˜ D ˜ u y (cid:17) ± (cid:12)(cid:12)(cid:12)(cid:12) β ; B ± R/ (˜ x ˜Γ ) ≤ C (cid:18) R β (cid:107) ˜ u y (cid:107) B R (˜ x ˜Γ ) + 1 R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) + 1 R β | ˜ g y | y (cid:19) + C (cid:18) | ˜ g y | ,β ;˜Γ y + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (˜ x ˜Γ ) (cid:19) ≤ C (cid:18) R β (cid:107) ˜ u y (cid:107) B R (˜ x ˜Γ ) + 1 R β (cid:107) ˜ g y (cid:107) y + | ˜ g y | ,β ;˜Γ y (cid:19) + C (cid:18) R β (cid:13)(cid:13)(cid:13) ˜ F ± y (cid:13)(cid:13)(cid:13) B ± R (˜ x ˜Γ ) + (cid:12)(cid:12)(cid:12) ˜ F ± y (cid:12)(cid:12)(cid:12) β ; B ± R (˜ x ˜Γ ) (cid:19) , with C = C ( n, Λ min , Λ max , β, η ) independent of ˜ x ˜Γ , of J ∈ N and of y ∈ P J . References [1]
I. Babuˇska, F. Nobile, and R. Tempone , A stochastic collocation method for ellipticpartial differential equations with random input data , SIAM Journal on Numerical Analysis,45 (2007), pp. 1005–1034.[2]
I. Babuska, R. Tempone, and G. E. Zouraris , Galerkin finite element approxi-mations of stochastic elliptic partial differential equations , SIAM Journal on NumericalAnalysis, 42 (2004), pp. 800–825.[3]
C. Bacuta, J. H. Bramble, and J. E. Pasciak , New interpolation results and applica-tions to finite element methods for elliptic boundary value problems , Journal of NumericalMathematics, 9 (2001), pp. 179–198.[4]
A. Barth, C. Schwab, and N. Zollinger , Multi-level Monte Carlo finite elementmethod for elliptic PDEs with stochastic coefficients , Numerische Mathematik, 119 (2011),pp. 123–161. 475]
J. Beck, F. Nobile, L. Tamellini, and R. Tempone , Convergence of quasi-optimalstochastic Galerkin methods for a class of PDEs with random coefficients , Computers &Mathematics with Applications, 67 (2014), pp. 732–751.[6]
J.-P. Berenger , A perfectly matched layer for the absorption of electromagnetic waves ,Journal of Computational Physics, 114 (1994), pp. 185–200.[7]
A. Bonito, J.-L. Guermond, and F. Luddens , Regularity of the Maxwell equationsin heterogeneous media and Lipschitz domains , Journal of Mathematical Analysis andapplications, 408 (2013), pp. 498–512.[8]
S. C. Brenner and R. Scott , The mathematical theory of finite element methods ,vol. 15, Springer Science & Business Media, 2008.[9]
R. E. Caflisch , Monte carlo and quasi-monte carlo methods , Acta numerica, 7 (1998),pp. 1–49.[10]
C. Canuto and T. Kozubek , A fictitious domain approach to the numerical solutionof PDEs in stochastic domains , Numerische mathematik, 107 (2007), pp. 257–293.[11]
J. E. Castrillon-Candas, F. Nobile, and R. F. Tempone , Analytic regularity andcollocation approximation for PDEs with random domain deformations , Comput. Math.Appl., 71 (2016), pp. 1173–1197.[12]
A. Chernov and C. Schwab , First order k-th moment finite element analysis of non-linear operator equations with stochastic data , Mathematics of Computation, 82 (2013),pp. 1859–1888.[13]
A. Chkifa, A. Cohen, and C. Schwab , Breaking the curse of dimensionality in sparsepolynomial approximation of parametric PDEs , Journal de Math´ematiques Pures et Ap-pliqu´ees, 103 (2015), pp. 400–428.[14]
K. A. Cliffe, M. B. Giles, R. Scheichl, and A. L. Teckentrup , Multilevel MonteCarlo methods and applications to elliptic PDEs with random coefficients , Computing andVisualization in Science, 14 (2011), pp. 3–15.[15]
A. Cohen, C. Schwab, and J. Zech , Shape Holomorphy of the stationary Navier-StokesEquations , Report 2016-45, Seminar for Applied Mathematics, ETH Z¨urich, Switzerland.[16]
F. Collino and P. Monk , The perfectly matched layer in curvilinear coordinates , SIAMJournal on Scientific Computing, 19 (1998), pp. 2061–2090.[17]
M. D. Dettinger and J. L. Wilson , First order analysis of uncertainty in numericalmodels of groundwater flow part: 1. Mathematical development , Water Resources Research,17 (1981), pp. 149–161.[18]
J. Dick, Q. T. L. Gia, and C. Schwab , Higher Order Quasi–Monte Carlo Integrationfor Holomorphic, Parametric Operator Equations , SIAM/ASA Journal on UncertaintyQuantification, 4 (2016), pp. 48–79.[19]
J. Dick, F. Y. Kuo, Q. T. Le Gia, D. Nuyens, and C. Schwab , Higher order QMCPetrov–Galerkin discretization for affine parametric operator equations with random fieldinputs , SIAM Journal on Numerical Analysis, 52 (2014), pp. 2676–2702.[20]
J. Dick, F. Y. Kuo, and I. H. Sloan , High-dimensional integration: the quasi-MonteCarlo way , Acta Numerica, 22 (2013), pp. 133–288.4821]
J. Dick and F. Pillichshammer , Digital nets and sequences: Discrepancy Theory andQuasi–Monte Carlo Integration , Cambridge University Press, 2010.[22]
T. J. Dodwell, C. Ketelsen, R. Scheichl, and A. L. Teckentrup , A hierarchicalmultilevel Markov chain Monte Carlo algorithm with applications to uncertainty quantifi-cation in subsurface flow , SIAM/ASA Journal on Uncertainty Quantification, 3 (2015),pp. 1075–1108.[23]
R. Gantner and M. Peters , Higher Order Quasi-Monte Carlo for Bayesian ShapeInversion , Report 2015-31, Seminar for Applied Mathematics, ETH Z¨urich, Switzerland.[24]
R. N. Gantner , A Generic C++ Library for Multilevel Quasi-Monte Carlo , in Proceed-ings of the Platform for Advanced Scientific Computing Conference, PASC ’16, New York,NY, USA, 2016, ACM, pp. 11:1–11:12.[25]
D. Gilbarg and N. S. Trudinger , Elliptic partial differential equations of secondorder , Springer, 2015.[26]
M. B. Giles , Multilevel monte carlo path simulation , Operations Research, 56 (2008),pp. 607–617.[27] ,
Multilevel Monte Carlo methods , Acta Numerica, 24 (2015), p. 259.[28]
P. Glasserman , Monte Carlo methods in financial engineering , vol. 53, Springer Science& Business Media, 2013.[29]
J. Guzm´an, D. Leykekhman, J. Rossmann, and A. H. Schatz , H¨older estimatesfor Green’s functions on convex polyhedral domains and their applications to finite elementmethods , Numerische Mathematik, 112 (2009), pp. 221–243.[30]
J. Guzm´an, M. S´anchez, and M. Sarkis , On the accuracy of finite element ap-proximations to a class of interface problems , Mathematics of Computation, 85 (2016),pp. 2071–2098.[31]
H. Harbrecht and J. Li , First order second moment analysis for stochastic interfaceproblems based on low-rank approximation , ESAIM: Mathematical Modelling and Numer-ical Analysis, 47 (2013), pp. 1533–1552.[32]
H. Harbrecht, M. Peters, and M. Siebenmorgen , Analysis of the domain map-ping method for elliptic diffusion problems on random domains , Numerische Mathematik,(2016), pp. 1–34.[33]
H. Harbrecht, R. Schneider, and C. Schwab , Sparse second moment analysis forelliptic problems in stochastic domains , Numerische Mathematik, 109 (2008), pp. 385–414.[34]
Y. Harness , Low-Dimensional Spatial Embedding Method for Shape Uncertainty Quan-tification in Acoustic Scattering , arXiv preprint arXiv:1704.07727, (2017).[35]
S. Heinrich , Monte Carlo complexity of global solution of integral equations , Journal ofComplexity, 14 (1998), pp. 151–175.[36]
S. Heinrich and E. Sindambiwe , Monte Carlo complexity of parametric integration ,Journal of Complexity, 15 (1999), pp. 317–341.4937]
R. Hiptmair, L. Scarabosio, C. Schillings, and C. Schwab , Large deforma-tion shape uncertainty quantification in acoustic scattering , Report 2015-31, Seminar forApplied Mathematics, ETH Z¨urich, Switzerland. .[38]
N. Hyvo¨onen, V. Kaarnioja, L. Mustonen, and S. Staboulis , Polynomial colloca-tion for handling an inaccurately known measurement configuration in electrical impedancetomography , SIAM Journal on Applied Mathematics, 77 (2017), pp. 202–223.[39]
C. Jerez-Sanchez, C. Schwab, and J. Zech , Electromagnetic Wave Scattering byRandom Surfaces: Shape Holomorphy , Report 2016-49, Seminar for Applied Mathematics,ETH Z¨urich, Switzerland.[40]
F. Jochmann , An Hs-regularity result for the gradient of solutions to elliptic equationswith mixed boundary conditions , Journal of mathematical analysis and applications, 238(1999), pp. 429–450.[41]
J. P. Krasovski , Isolation of singularities of the Green’s function , Mathematics of theUSSR-Izvestiya, 1 (1967), p. 935.[42]
J. Li, J. M. Melenk, B. Wohlmuth, and J. Zou , Optimal a priori estimates forhigher order finite elements for elliptic interface problems , Applied numerical mathematics,60 (2010), pp. 19–37.[43]
C. Miranda , Partial differential equations of elliptic type , vol. 2, Springer Science &Business Media, 2012.[44]
J.-C. N´ed´elec , Acoustic and electromagnetic equations: integral representations for har-monic problems , vol. 144, Springer Science & Business Media, 2001.[45]
F. Nobile, R. Tempone, and C. G. Webster , A sparse grid stochastic collocationmethod for partial differential equations with random input data , SIAM Journal on Numer-ical Analysis, 46 (2008), pp. 2309–2345.[46]
A. Nouy, A. Clement, F. Schoefs, and N. Mo¨es , An extended stochastic finiteelement method for solving stochastic partial differential equations on random domains ,Computer Methods in Applied Mechanics and Engineering, 197 (2008), pp. 4663–4682.[47]
A. Nouy, F. Schoefs, and N. Mo¨es , X-SFEM, a computational technique based on X-FEM to deal with random shapes , European Journal of Computational Mechanics/RevueEurop´eenne de M´ecanique Num´erique, 16 (2007), pp. 277–293.[48]
G. Savar´e , Regularity results for elliptic equations in Lipschitz domains , Journal of Func-tional Analysis, 152 (1998), pp. 176–201.[49]
L. Scarabosio , Shape uncertainty quantification for scattering transmission problems ,PhD thesis, ETH Z¨urich, 2016. Diss. No. 23574. http://e-collection.library.ethz.ch/eserv/eth:49652/eth-49652-02.pdf .[50]
A. Schatz , Pointwise error estimates and asymptotic error expansion inequalities forthe finite element method on irregular grids: Part I. Global estimates , Mathematics ofComputation of the American Mathematical Society, 67 (1998), pp. 877–899.[51]
C. Schillings and C. Schwab , Sparse, adaptive Smolyak quadratures for Bayesianinverse problems , Inverse Problems, 29 (2013), p. 065011.5052]
C. Schwab and C. J. Gittelson , Sparse tensor discretizations of high-dimensionalparametric and stochastic PDEs , Acta Numerica, 20 (2011), pp. 291–467.[53]
J. Sokolowski and J.-P. Zolesio , Introduction to Shape Optimization , Springer, 1992.[54]
D. M. Tartakovsky and D. Xiu , Stochastic analysis of transport in tubes with roughwalls , Journal of Computational Physics, 217 (2006), pp. 248–259.[55]
Z. Wu, J. Yin, and C. Wang , Elliptic & parabolic equations , World Scientific, 2006.[56]
D. Xiu and J. S. Hesthaven , High-order collocation methods for differential equationswith random inputs , SIAM Journal on Scientific Computing, 27 (2005), pp. 1118–1139.[57]
D. Xiu and G. E. Karniadakis , The Wiener–Askey polynomial chaos for stochasticdifferential equations , SIAM journal on scientific computing, 24 (2002), pp. 619–644.[58]
D. Xiu and D. M. Tartakovsky , Numerical methods for differential equations in ran-dom domains , SIAM Journal on Scientific Computing, 28 (2006), pp. 1167–1185.[59]