Asymptotic analysis for elliptic equations with small perturbations on domains in high-contrast medium

Abstract

We provide a comprehensive study on the asymptotic solutions of an interface problem corresponding to an elliptic partial differential equation with Dirichlet boundary condition and transmission condition, subject to the small geometric perturbation and/or the high contrast ratio of the conductivity. All asymptotic terms can be solved in the unperturbed reference domains, which significantly reduces computations in practice, especially for random perturbations. Our setting is quite general and allows two types of elliptic problems: the perturbation of the domain boundary where the Dirchlet condition is imposed and the perturbation of the interface where the transmission condition is imposed. As the perturbation size and the ratio of the conductivities tends to zero, the two-parameter asymptotic expansions on the reference domain are derived to any order after the single parameter expansions are solved beforehand. The results suggest the emergence of the Neumann or Robin boundary condition, depending on the relation of the two asymptotic parameters. Our method is the classic asymptotic analysis techniques but in a new unified approach to both problems.

Keywords

Asymptotic analysis interface problem high-contrast ratio two-parameter expansion

1. Introduction

Let $D \subset R^{d}$ , $d = 1, 2, 3$ , be a simply-connected Lipschitz continuous domain. Consider the perturbation of the domain D given by the perturbed boundary $\partial D_{ε}$ defined as $\begin{matrix} (1.1) & \partial D_{ε} = {x^{'} : x^{'} = x + ε h (x) n (x) : x \in \partial D}, \end{matrix}$ where $ε \in (0, ε_{0}]$ for a fixed small number $ε_{0} ≪ 1$ represents the small characteristic size of the perturbation, $h (x)$ is a continuous function defined on $\partial D$ , and $n (x)$ is the (outward) normal direction of D. For sufficiently small ε, the boundary $\partial D_{ε}$ uniquely defines a perturbed domain $D_{ε}$ . If h is non-negative, then $D_{ε}$ contains D. We assume that h is a sufficiently smooth function.

The main problem of our concern is related to the following Dirichlet boundary value elliptic problem imposed in the perturbed domain $D_{ε}$ : $\begin{matrix} (1.2) & \{\begin{matrix} L u_{ε} = f & in D_{ε}, \\ u_{ε} = g & on \partial D_{ε}, \end{matrix} \end{matrix}$ where $L$ is the second order elliptic operator, having the divergence form $\begin{matrix} (1.3) & L u = - \sum_{i, j = 1}^{d} \partial_{x_{j}} (a^{i j} (x) \partial_{x_{i}} u) + \sum_{i = 1}^{d} b^{i} (x) \partial_{x_{i}} u + c (x) u . \end{matrix}$ The second order coefficient functions $a^{i j}$ , $i, j = 1, \dots, d$ , form a non-degenerate positive definite matrix $a = (a^{i j})$ , i.e., $a^{j i} = a^{i j}$ , and $\begin{matrix} (1.4) & \sum_{i, j = 1}^{d} a^{i j} (x) ξ_{i} ξ_{j} > 0, \end{matrix}$ for every $x \in \overline{D \cup D_{ε_{0}}}$ and non-zero vector $(ξ_{1}, \dots, ξ_{d}) \in R^{d}$ . The coefficients $b^{i}$ , $i = 1, \dots, d$ , and c are assumed smooth in $R^{d}$ . The boundary value function g is also assumed smooth in an open neighbourhood of $\partial D$ .

If the coefficient $a$ is assumed to be continuous everywhere, then the solution $u_{ε}$ is the perturbation of a classic elliptic equation with uncertainty in characterizing the domain. How to quantify the uncertainty in the solution due to the geometric perturbation, particularly when h is a random function, is an interesting and important topic in uncertainty quantification. The more challenging case is that $a$ is not continuous across some interface. Then the transmission condition should be specified on the jump interface. In such cases, the interface may also be subject to small perturbations.

There are two scenarios of the geometric perturbations in the transmission problems. The first one is to consider the previous domain perturbation setup with a non-negative h, then $D \subset D_{ε}$ and the interface is $Γ = \partial D$ , which is fixed and separates the domain D and the thin layer $\begin{matrix} L_{ε} = {x^{'} : x^{'} = x + t h (x) n (x), 0 < t < ε, x \in \partial D, h (x) \neq 0} . \end{matrix}$ We call this model the thin layer problem. The second scenario is to partition a fixed domain D into two subdomains: $D = D_{ε}^{+} \cup D_{ε}^{-} \cup Γ_{ε}$ , where $Γ_{ε}$ is the dividing interface, which is assumed as a perturbation on a fixed interface Γ. The difference between $Γ_{ε}$ and Γ can be also described by a function h, as in the first scenario. The detailed definitions of $D_{ε}^{+}, D_{ε}^{-}, Γ_{ε}$ will be specified later. We call this model the perturbed interface problem. In the first problem, we attach a thin layer $L_{ε}$ to encircle the fixed domain D and the layer thickness vanishes as ε tends to zero. The interface there is fixed. In the second problem, we partition a fixed domain D into two subdomains $D_{ε}^{\pm}$ by a perturbed interface $Γ_{ε}$ and the two subdomains have comparable size. Refer to Fig. 1 and Fig. 2 later for illustrations of geometry.

All these perturbations can be either deterministic or random, depending on whether h is a deterministic function or a random field. For the latter case, after h is represented in random space by Karhunen–Loève theorem $h (x, ω) = \sum h_{i} (x) ϕ_{i} (ω)$ , or by the Monte Carlo samples $h (x, ω) \sim h_{i} (x)$ , the problem is then transformed to a set of deterministic perturbations if the correlation length of h is not vanishing. So, we only focus on the deterministic h here; the application to the random case may follow the standard approaches used in many literatures such as [5,8,15,24].

There is a distinctive class of perturbations of the domain for the PDE (1.2): the so called “rough boundary/rough domain”, in which the spatial scale of the profile h also depends on ε, for instance, D is perturbed by the form $ε h_{ε} (x) = ε \bar{h} (x / ε)$ for a periodic function $\bar{h}$ (see [19] and references therein). When the boundary condition itself also involves the similar multiscale feature, the multiscale finite element method was applied and analyzed by [21].

To explicitly show the transmission condition and to introduce our second asymptotic parameter other than the perturbation size ε, we take the simplest case of the thin layer problem corresponding to the first scenario mentioned above. In this case, $D \subset D_{ε}$ , $Γ = \partial D$ is the interface, separating the domain D and the thin layer $L_{ε} = D_{ε} ∖ \bar{D}$ . Assume that the coefficients b and c vanish and that $a$ is scalar-valued and is piecewisely homogeneous in D and $L_{ε}$ . Then the corresponding transmission problem takes the form $\begin{matrix} (1.5) & \{\begin{matrix} - Δ u_{ε} = f & in D, \\ - σ Δ u_{ε} = f & in L_{ε}, \\ u_{int} = u_{ext}, \partial_{n} u_{int} = σ \partial_{n} u_{ext} & on Γ, \\ u_{ε} = 0 & on \partial D_{ε}, \end{matrix} \end{matrix}$ where σ is a constant parameter representing the ratio of conductivity in two different domains. $u_{int}$ and $u_{ext}$ are the restrictions of the solution $u_{ε}$ on two subdomains D and $L_{ε}$ , respectively. The similar form of the transmission condition will be specified later for the general problems. If the material property across the interface has a significant difference, then the value of σ can take a very small value or a very large value. The resulted transmission problem in this high-contrast media is an important subject in multiscale analysis and computation.

The elliptic model (1.2) and the transmission problem such as (1.5) originate from many applications such as diffusion processes, electrostatics, porous media and heat conduction. One of our motivating examples is the diffusion model of exciton in organic semiconductors [6,14,18]. For the discontinuous coefficient model (1.5), a well-known problem is the electromagnetic model for bodies coated with a dielectric layer $L_{ε}$ with distinctive material coefficients. In porous media applications, the permeability of subsurface regions is described as a quantity with high-contrast and multiscale features.

We here mainly concern the asymptotic analysis in terms of the two different parameters, ε and σ, where ε represents the amplitude of the geometric perturbation on the domain or the interface, and σ represents the ratio of different material coefficients. In this paper, we shall first consider the asymptotic effect of each parameter separately and then work on the more complicated two-parameter expansions.

Many theories and methods have been developed and used to study the above elliptic problems and the interface problems. We review some general methodologies on the asymptotic study for the solution $u_{ε}$ subject to the geometric perturbations. The first idea to handle the irregular domain $D_{ε}$ is the domain mapping, which is to find a smooth mapping to change the irregular domain to a fixed reference domain. See the reference [5,16,24] for the applications and the analysis of this method. This method works for any irregular domain as long as a diffeomorphism can be found regardless it is a small perturbation or not. By applying the diffeomorphism transformation, all geometric information is transformed into a new differential operator and a new boundary condition, which are both more complicated than the original form on irregular domain. The second method, particularly for the perturbed interface problem, is a generalization of calculus of variation to the geometric setting – the shape derivative [15,17]. The method of shape derivatives is widely used for the sensitivity analysis of the geometry of the boundary and shape optimization. Although it is quite easy to obtain the first order derivative, the calculation is very complicated for the higher order derivatives. The last method, which is also our main tool here, is the asymptotic expansion, which actually refers to a collection of problem-specific methods and relies on the correct use of the ansätz [1–4,8,10,23]. In this method, by using a good regularity of the solution in the correct (sub)domains, one can apply certain ansätz in the form of the series expansion to approximate the boundary conditions on the fixed domain. More details on the application of this method to our problems of concerns will be reviewed and commented in subsequent sections.

The main motivation of this article is to give a comprehensive study on the (formal) asymptotic expansions of the solutions to the above various elliptic problems, including the thin layer problem and the interface problem, up to an arbitrary order in theory. Specifically, we shall address the following four problems.

The first task is that for the elliptic model (1.2) with smooth $a$ , we want to have in D $\begin{matrix} (1.6) & u_{ε} = u_{0} + ε u_{1} + ε^{2} u_{2} + ε^{3} u_{3} + \dots, \end{matrix}$ in certain sense, where all terms $u_{i}$ are independent of ε explicitly. Then we want to construct a sequence ${u^{[n]}}$ of functions satisfying the following properties: (i) Each $u^{[n]}$ is the solution to a boundary value problem defined only on the fixed domain D; (ii) The error between the restriction of $u_{ε}$ to D and $u^{[n]}$ is limited to the order $O (ε^{n + 1})$ ; (iii) The numerical computation of $u^{[n]}$ (the detailed implementation is not our objective in this paper) should be easier than directly solving the original equation (1.2). Note that $u^{[n]}$ is not simply the partial sum $\sum_{i = 0}^{n} ε^{n} u_{n}$ , because the latter may not satisfy a closed boundary value problem.

The second task is to generalize the results in (I) to the thin layer problem (1.5) for the case of the discontinuous coefficient a.

The third one is the generalization of (II) to the high-contrast material, i.e., σ, the ratio of material coefficients across the interface Γ, is very large or very small. We want to derive the two-parameter expansions when the limits of both ε and σ are considered. We are concerned with the three scaling regimes for ε and σ: $ε / σ \to 0$ , $ε / σ \to \infty$ , and $ε / σ \to c \in (0, \infty)$ . The final result is the boundary value problem for each term in the two-parameter asymptotic expansions $u_{int} (x) = \sum_{m, n} u_{m, n} (x) ε^{m} μ^{n}$ , where $m, n$ are integers, and μ is linked to the ratio of ε and μ, whose specific form depends on the asymptotic regimes. We shall show that the three scalings will give arise to the Dirichlet, Neumann or Robbin boundary condition for $u_{m, n}$ , respectively.

The last one is on the perturbed interface problem where the interface $Γ_{ε}$ is not fixed as in (II) and (III), but is associated with a perturbed domain partition $D = D_{ε}^{+} \cup D_{ε}^{-} \cup Γ_{ε}$ . Meanwhile, the high-contrast ratio limit is also considered, and we derive the two-parameter asymptotic expansions, where we find there is no special dependence on the scaling of ε and σ.

From Section 2 to Section 5, we solve each of these four problems in each section. The main contributions of our works are summarized below. As far as the authors know, the specific techniques we used for (I) and (II) are different from the classic dilation methods in the existing literature (e.g. [8]) and the two-parameter asymptotic expansions for (III) and (IV) in this paper are completely new results. The main techniques we apply here for all four problems are the Taylor expansion applied in various contents, which all requires a good regularity of the underlying function. Despite of the simplicity of this Taylor expansion idea, one of our important contributions is to rigorously justify its applicability in particular for discontinuous $a$ . For the thin layer problem or the interface problem, where the solution $u_{ε}$ apparently does not posses such smoothness on the interface, our idea is first to extend each smooth component of the solution $u_{ε}$ restricted on each subdomain onto some ε-independent domains. This is achieved by imposing certain Cauchy problems on the interface when interpreting the elliptic equation as a time-evolution equation in which the normal direction of the interface is the time marching direction. The second important idea we used is to apply the inverse Lax–Wendroff procedure [22] to convert the high order derivatives in the normal direction on the interface to those along the tangent directions and the first order normal derivative, for which the original transmission condition on the interface is utilized.

Convergence results for the asymptotic expansions are given in the Appendix, based on interior estimates from traces. Besides, some explicit examples are detailed in order to show that a truncated series may provide a better approximation than an impedance condition. From the practical perspective, the main advantage of a truncated series concerns the approximation of $u_{ε}$ for many values of ε. In this context, the terms constituting the series can be computed once and for all over a regular domain, and then any new value of ε only requires a linear combination, while the original impedance problem needs to be solved over irregular domains independently for each ε. Interested readers may find some application of the current work to design an efficient numerical method for the estimation of exciton diffusion length in organic semiconductors [7].

To end this introduction, we review several existing works which are closely related to the problems we considered here. The work in [8] considered the thin layer problem (1.5) with a fixed σ as $ε \to 0$ . The main idea in [8] is to write the differential operator $L$ in terms of local coordinate in the thin layer $L_{ε}$ , and apply the ansätz $L = \sum_{n ⩾ - 2} ε^{n} L_{n}$ to derive a system of (infinitely number of) recursive equations for the expansion of the solution in this dilated layer. Then with the aid of the transmission condition on the interface Γ, the boundary conditions of these equations in the layer $L_{ε}$ are linked to the solutions in the interior (fixed) domain D. In [3], to assist the construction of local solutions in the multiscale finite element methods for the elliptic equations in high-contrast media, the authors derived asymptotic expansions for the solutions of the elliptic problems with high contrast ratio, i.e., σ tends to 0 or ∞. But their analysis is for the fixed domain and interface.

2. The elliptic problem with smooth coefficients

In this section, we study the equation (1.2) on $D_{ε}$ by assuming that $a (x)$ is sufficiently smooth everywhere and $h (x)$ in (1.1) is also sufficiently smooth on $\partial D$ . This means that the Taylor expansion for these two functions are available up to any order. The signs of $h (x)$ can be arbitrary at different $x \in \partial D$ and the operator $L$ in (1.3) is not limited to the Laplace operator.

Recall that the perturbed thin layer $L_{ε}$ is defined by $\begin{matrix} L_{ε} = {x^{'} : x^{'} = x + t h (x) n (x), 0 < t < ε, x \in \partial D, h (x) \neq 0} . \end{matrix}$ The condition $h (x) \neq 0$ ensures that $L_{ε}$ is also a domain (open set). Depending on the sign of the function h, we can decompose the thin layer $L_{ε}$ into the interior layer $L_{ε, int}$ and the external layer $L_{ε, ext}$ : $\begin{matrix} L_{ε} = L_{ε, int} \cup L_{ε, ext}, \end{matrix}$ where $\begin{array}{l} L_{ε, int} : = L_{ε} \cap D = {x^{'} : x^{'} = x + t h (x) n (x), 0 < t < ε, x \in \partial D, h (x) < 0}, \\ L_{ε, ext} : = L_{ε} ∖ D = {x^{'} : x^{'} = x + t h (x) n (x), 0 < t < ε, x \in \partial D, h (x) > 0} . \end{array}$ $L_{ε, int} \subset D$ and $L_{ε, ext} \cap D = \emptyset$ . Then $D_{ε}$ is the interior of $\overline{(D ∖ L_{ε, int}) \cup L_{ε, ext}}$ . Refer to the schematic illustration in Fig. 1.

Fig. 1.

Schematic illustration of the domain perturbation. The regular domain D is in the “ball” shape and the perturbed domain $D_{ε}$ is in the “ellipse” shape.

2.1. Approximate expansions

The problem (1.2) is defined on the ε-dependent domain $D_{ε}$ . We extend it to a fixed domain $D \cup D_{ε_{0}}$ and justify this extension in Section 2.1.1. Then in Section 2.1.2, we use the Taylor expansion near $\partial D$ to derive the asymptotic expansion $u_{ε} = \sum_{n ⩾ 0} ε^{n} u_{n}$ , for which the inverse Lax–Wendroff procedure is applied to convert the high order normal derivatives into the first order normal derivative and the tangential derivatives along the boundary $\partial D$ .

[8] already derived the first three terms, $u_{0}$ , $u_{1}$ and $u_{2}$ . But the method we give below seems simpler and does not require the dilation technique and any asymptotic form for the differential operator $L$ used in [8]. Actually, that kind of singular perturbation suits for the case that the solution itself develops a sharp peak in the thin layer, such as the traditional boundary layer analysis in fluid mechanics. However, the problem here does not have this feature and the solutions on D and $D_{ε}$ both behave very normally at the order $O (1)$ . We find that the direct expansion for the boundary condition of $u_{ε}$ in an appropriate way is sufficient to derive the boundary condition of $u_{n}$ . To present our main technique, we start with the smooth $a$ case in this section and then show how to generalize to the discontinuous $a$ in Section 3.

2.1.1. The extension of the solution to the fixed domain

Note that $D \cup D_{ε}$ is increasing in ε since $L_{ε, ext}$ always expands as ε increases. So it is convenient to make the extension to the whole domain $D \cup D_{ε_{0}}$ since we only consider $ε \in (0, ε_{0}]$ . On this fixed domain $D \cup D_{ε_{0}}$ , the solution is known on the part $\overline{D_{ε}}$ ; we thus consider the difference $△_{ε}$ which consists of the disjoint thin layers: $\begin{matrix} △_{ε} : = (D \cup D_{ε_{0}}) ∖ \overline{D_{ε}} = L_{ε, int} \cup N_{ε}, where N_{ε} : = L_{ε_{0}, ext} ∖ \overline{L_{ε, ext}} . \end{matrix}$ Denote the solution extended on $△_{ε}$ by ${\tilde{u}}_{ε}$ , and assume that ${\tilde{u}}_{ε}$ and $u_{ε}$ have the same values and the same normal derivatives on the common boundary $\partial D_{ε}$ . Specifically, ${\tilde{u}}_{ε}$ is the unique solution to the following Cauchy problem posed in the thin layers $L_{ε, int}$ and $N_{ε}$ : $\begin{matrix} (2.1) & \{\begin{matrix} L {\tilde{u}}_{ε} = f & in L_{ε, int} \cup N_{ε}, \\ {\tilde{u}}_{ε} = u_{ε} = g, \partial_{n} {\tilde{u}}_{ε} = \partial_{n} u_{ε} & on \partial D_{ε}, \end{matrix} \end{matrix}$ where $u_{ε}$ , the solution to equation (1.2), is presumably given, $n$ is the outward normal of $D_{ε}$ on $\partial D_{ε}$ . Note that $\partial D_{ε}$ is a proper subset of the boundaries of $L_{ε, int}$ and $N_{ε}$ . The problem (2.1) is actually a Cauchy problem of ${\tilde{u}}_{ε}$ , not a boundary-valued elliptic problem, because the value and the “velocity” of ${\tilde{u}}_{ε}$ are specified on $\partial D_{ε}$ – a part of its complete boundary. The boundary $\partial D_{ε}$ satisfies the noncharacteristic condition $\sum_{i, j = 1}^{d} a^{i j} n_{i} n_{j} \neq 0$ trivially since $L$ is, by assumption, an elliptic operator satisfying (1.4). Thus by the Cauchy–Kovalevskaya theorem [11], the solution on $\partial D$ can propagate to the boundary $\partial △_{ε}$ and the above Cauchy problem (2.1) is well-posed for sufficiently small $ε_{0}$ .

Remark 2.1.
The above method of extending the solution to a larger (and ε-indepedent) domain can also preserve the regularity of the solution and helps clarify the rigorous meaning of the Taylor expansion we shall apply. This extension idea by the use of the Cauchy problem of a time-evolution equation will be applied repeatedly in this paper, especially for the interface problem so that each smooth component of the solution on each subdomain may be approximated by the Taylor expansion along some interface.

Now it is clear that we can define a function $w_{ε}$ piecewisely on the whole (fixed) domain $\overline{D} \cup \overline{D_{ε_{0}}} = \overline{D_{ε}} \cup \overline{△_{ε}}$ as follows: $\begin{matrix} (2.2) & w_{ε} (x) : = \{\begin{matrix} u_{ε} (x) & in \overline{D_{ε}}, \\ {\tilde{u}}_{ε} (x) & in \overline{△_{ε}} . \end{matrix} \end{matrix}$ This definition is justified by the boundary condition in (2.1) which dictates that $u_{ε}$ and ${\tilde{u}}_{ε}$ coincide on the common boundary $\partial D_{ε}$ . Then $w_{ε}$ satisfies the equation on the fixed domain $\begin{matrix} (2.3) & L w_{ε} = f in D \cup D_{ε_{0}}, \end{matrix}$ and on the ε-dependent boundary. $\begin{matrix} (2.4) & w_{ε} = g, on \partial D_{ε} . \end{matrix}$ Note that (2.4) does not serve as a boundary condition to the equation (2.3). $w_{ε}$ is simply a combination of $u_{ε}$ from the boundary value problem (1.2) and ${\tilde{u}}_{ε}$ from the Cauchy problem (2.1). The above argument of extension ensures that $w_{ε}$ has the same regularity of $u_{ε}$ , but on $D \cup D_{ε_{0}}$ .
2.1.2. Asymptotic expansion on the whole domain

By the above extension, we can assume the following ansätz for $w_{ε}$ , $\begin{matrix} (2.5) & w_{ε} (x) = \sum_{n = 0}^{\infty} ε^{n} w_{n} (x) for x \in \overline{D \cup D_{ε_{0}}} . \end{matrix}$ Plug this ansätz into the equation (2.3), and match the terms at the same order of ε, then we obtain the following equations for $w_{n}$ in $D \cup D_{ε_{0}}$ : $\begin{matrix} (2.6) & L w_{n} = δ_{0, n} f . \end{matrix}$ Here $δ_{i, j} = 1$ if $i = j$ and $δ_{i, j} = 0$ if $i \neq j$ .

For the condition (2.4), $w_{ε} = g$ on $\partial D_{ε}$ , by noticing the fact that $x + ε h (x) n (x) \in \partial D_{ε}$ for all $x \in \partial D$ , we have $\begin{matrix} (2.7) & w_{ε} (x + ε h (x) n (x)) = \sum_{n = 0}^{\infty} ε^{n} w_{n} (x + ε h (x) n (x)) = g (x + ε h (x) n (x)) . \end{matrix}$ The Taylor expansions in ε on the right-hand side read $\begin{array}{l} (2.8) & w_{n} (x + ε h (x) n (x)) = \sum_{k = 0}^{\infty} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} w_{n} (x), \\ (2.9) & g (x + ε h (x) n (x)) = \sum_{k = 0}^{\infty} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} g (x), \end{array}$ where for any vector field $n (x) = (n_{1} (x), \dots, n_{d} (x))$ , the kth directional derivative along $n$ at $x = x_{0}$ is defined by $\begin{matrix} \partial_{n}^{k} |_{x = x_{0}} : = {(\sum_{i = 1}^{d} n_{i} (x_{0}) \partial_{x_{i}} |_{x = x_{0}})}^{k} . \end{matrix}$ Then (2.7), (2.8) and (2.9) together lead to $\begin{matrix} \sum_{n = 0}^{\infty} ε^{n} \sum_{k = 0}^{\infty} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} w_{n} (x) = \sum_{k = 0}^{\infty} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} g (x), \end{matrix}$ which, by a change of the indices $m = k + n$ , is equivalent to $\begin{matrix} \sum_{m = 0}^{\infty} ε^{m} \sum_{k = 0}^{m} \frac{{(h (x))}^{k}}{k!} \partial_{n}^{k} w_{m - k} (x) = \sum_{m = 0}^{\infty} \frac{ε^{m} {(h (x))}^{m}}{m!} \partial_{n}^{m} g (x) . \end{matrix}$ Then by matching the terms with the same order of ε, we obtain that $\begin{matrix} \sum_{k = 0}^{m} \frac{{(h (x))}^{k}}{k!} \partial_{n}^{k} w_{m - k} (x) = \frac{{(h (x))}^{m}}{m!} \partial_{n}^{m} g (x), \end{matrix}$ i.e., $\begin{matrix} (2.10) & \{\begin{matrix} w_{0} (x) = g (x), \\ w_{m} (x) = \frac{{(h (x))}^{m}}{m!} \partial_{n}^{m} g (x) - \sum_{k = 1}^{m} \frac{{(h (x))}^{k}}{k!} \partial_{n}^{k} w_{m - k} (x), \forall m ⩾ 1 . \end{matrix} \end{matrix}$ This provides a recursive expression of the boundary condition on $\partial D$ for the mth order term $w_{m}$ .

Define $u_{n}$ as the restriction of $w_{n}$ to D. Then $u_{ε} = \sum_{n = 0}^{\infty} ε^{n} u_{n}$ . By (2.6) and (2.10), $u_{n}$ satisfies the following sequence of boundary value problems on D where the boundary conditions on $\partial D$ are defined recursively: $\begin{matrix} (2.11) & \{\begin{matrix} L u_{0} = f & in D, \\ u_{0} = g & on \partial D, \end{matrix} \end{matrix}$ and for $n ⩾ 1$ , $\begin{matrix} (2.12) & \{\begin{matrix} L u_{n} = 0 & in D, \\ u_{n} (x) = \frac{{(h (x))}^{n}}{n!} \partial_{n}^{n} g (x) - \sum_{k = 1}^{n} \frac{{(h (x))}^{k}}{k!} \partial_{n}^{k} u_{n - k} (x) & on \partial D . \end{matrix} \end{matrix}$ In particular, for $n = 1, 2, 3$ , the above boundary conditions on $\partial D$ are $\begin{array}{l} (2.13) & u_{1} (x) = h (x) \partial_{n} g (x) - h (x) \partial_{n} u_{0} (x), \\ (2.14) & u_{2} (x) = \frac{{(h (x))}^{2}}{2} \partial_{n}^{2} g (x) - h (x) \partial_{n} u_{1} (x) - \frac{{(h (x))}^{2}}{2} \partial_{n}^{2} u_{0} (x), \\ (2.15) & u_{3} (x) = \frac{{(h (x))}^{3}}{6} \partial_{n}^{3} g (x) - h (x) \partial_{n} u_{2} (x) - \frac{{(h (x))}^{2}}{2} \partial_{n}^{2} u_{1} (x) - \frac{{(h (x))}^{3}}{6} \partial_{n}^{3} u_{0} (x) . \end{array}$

Remark 2.2.
Using the shape calculus method, one may also derive a “shape-Taylor expansion” of $u_{ε}$ on any compact set $K \subset D \cap D_{ε}$ (see [17] and the references therein), $\begin{matrix} u_{ε} (x) = u_{0} (x) + ε d [U] (x) + \frac{ε^{2}}{2} d^{2} [U, U] (x) + O (ε^{3}), \end{matrix}$ where $u_{0}$ is the solution to (2.11), $d [U]$ is the first order shape derivative on the boundary variation $U$ , which is given by the Dirichlet problem $\begin{matrix} \{\begin{matrix} L d [U] = 0 in D, \\ d [U] = U \cdot n \partial_{n} (g - u_{0}) on \partial D . \end{matrix} \end{matrix}$ $d^{2} [U, U^{'}]$ is the second order shape derivative, i.e., the “shape Hessian”, on the pair $(U, U^{'})$ of boundary variations, which is given by the Dirichlet problem $\begin{matrix} \{\begin{matrix} L d^{2} [U, U^{'}] = 0 & in D, \\ d^{2} [U, U^{'}] = \partial_{U} \partial_{U^{'}} (g - u_{0}) - \partial_{U} d [U^{'}] - \partial_{U^{'}} d [U] & on \partial D . \end{matrix} \end{matrix}$ It is easy to see that when the boundary variation $U (x)$ is given by $U (x) = h (x) n (x)$ for $x \in \partial D$ , then $d [U] = u_{1}$ and $d^{2} [U, U] = 2 u_{2}$ . Therefore the shape calculus method produces the same result as our method.

The right-hand side of the boundary condition (2.12) for each $u_{n}$ involves the normal derivatives of all lower order terms. The inverse Lax–Wendroff procedure, which is used to construct high order numerical methods such as in [22], enables us to convert the high order normal derivatives into the first order normal derivative and the tangential derivatives on the boundary $\partial D$ . See Lemma 2.3 below. This conversion procedure here seems only optional in theory, but as we shall show in Section 3, for piecewisely smooth coefficients, this step is essential for the use of transmission conditions on the interface to link the interior solution and the exterior solution.
Lemma 2.3.
Let u satisfy $L u = f$ where $L$ is the elliptic operator in (1.3). Then all the normal derivatives $\partial_{n}^{k} u$ on a smooth surface Γ with order $k ⩾ 2$ can be expressed in terms of the boundary Γ, the restrictions of the function u and its normal derivative $\partial_{n} u$ on Γ, and the coefficient functions $a^{i j}$ , $b^{i}$ , c, $i, j = 1, \dots, d$ . Therefore for every $k ⩾ 2$ , every smooth surface Γ, every elliptic operator $L$ and every smooth function f, there exists an operator $\begin{matrix} F_{k, Γ, L, f} [\cdot, \cdot] \end{matrix}$ acting on a pair of functions defined on Γ such that for any smooth function u satisfying $L u = f$ , its kth normal derivative $\partial_{n}^{k} u$ on Γ is given by $F_{k, Γ, L, f} [u, \partial_{n} u]$ .

In addition, it is easy to see the following properties of the operator $F_{k, Γ, L, f} [\cdot, \cdot]$ from the linearity of $L$ : $\begin{array}{l} F_{k, Γ, L, f} [u, \partial_{n} u] + F_{k, Γ, L, φ} [v, \partial_{n} v] = F_{k, Γ, L, f + φ} [u + v, \partial_{n} u + \partial_{n} v], \\ c F_{k, Γ, L, f} [u, \partial_{n} u] = F_{k, Γ, L, c f} [c u, c \partial_{n} u], \forall c \in R, \end{array}$ where u and v solve $L u = f$ and $L v = φ$ respectively. In particular, taking $c = 0$ in the last equality yields $F_{k, Γ, L, 0} [0, 0] = 0$ .

For the proof of this lemma, refer to Theorem 1 in Section 4.6 of [11]. The crucial assumption for the proof is the noncharacteristic condition of Γ, which is automatically guaranteed by the ellipticity of $L$ . This lemma will be used later multiple times and the dependency on Γ and $L$ in the notation of the mapping F may be dropped out if they are self-explanatory.

With this notation F, the boundary condition for $u_{n} (n ⩾ 1)$ in (2.12) can be formally written as $\begin{matrix} u_{n} (x) = & \frac{{(h (x))}^{n}}{n!} \partial_{n}^{n} g (x) - h (x) \partial_{n} u_{n - 1} (x) \\ - \sum_{k = 2}^{n} \frac{{(h (x))}^{k}}{k!} F_{k, \partial D, L, δ_{k, n} f} [u_{n - k}, \partial_{n} u_{n - k}] (x) . \end{matrix}$

To demonstrate the above theory and show how the conversion of the higher order normal derivatives works, in Appendix A, we present two examples in 2D. The first is our motivating example of exciton diffusion and the second is the Poisson equation. Furthermore, in Appendix A, we demonstrate how to generalize our method to the Neumann boundary condition and the reaction-diffusion equation with nonlinear terms.
2.2. The partial sums

We have formally derived the hierarchic systems of the boundary value problems for the expansion terms ${u_{n}}$ in Section 2.1. We next derive the closed boundary value problems which the partial sums approximately satisfy. The procedure is the same as in [8]. Define the partial sums $\begin{matrix} v^{[n]} (x) : = \sum_{k = 0}^{n} ε^{k} u_{k} (x), n ⩾ 0 . \end{matrix}$ On the boundary $\partial D$ , by using (2.12), we have $\begin{array}{l} v^{[n]} (x) & = \sum_{k = 0}^{n} ε^{k} u_{k} (x) \\ = \sum_{k = 0}^{n} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} g (x) - \sum_{k = 0}^{n} ε^{k} \sum_{j = 1}^{k} \frac{{(h (x))}^{j}}{j!} \partial_{n}^{j} u_{k - j} (x) \\ (2.16) & = \sum_{k = 0}^{n} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} g (x) - \sum_{j = 1}^{n} \frac{ε^{j} {(h (x))}^{j}}{j!} \partial_{n}^{j} v^{[n - j]} (x) . \end{array}$ It is worth pointing out that the system of the boundary value problems for $v^{[n]}$ is defined recursively. To obtain $v^{[n]}$ , one needs to solve the boundary value problems from $v^{[0]}$ (i.e., $u_{0}$ ) up to $v^{[n - 1]}$ . Thus, in total, $(n + 1)$ Dirichlet boundary value problems have to be solved. However, it is possible to directly solve one boundary value problem to obtain the approximation with the same order as $v^{[n]}$ by replacing the $v^{[n - j]}$ terms on the right-hand side of (2.16) by $v^{[n]}$ . Then one obtains the following closed boundary value problem, whose solution is denoted by $u^{[n]}$ : $\begin{matrix} (2.17) & \{\begin{matrix} L u^{[n]} = f & in D, \\ \sum_{k = 0}^{n} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} u^{[n]} (x) = \sum_{k = 0}^{n} \frac{ε^{k} {(h (x))}^{k}}{k!} \partial_{n}^{k} g (x) & on \partial D . \end{matrix} \end{matrix}$ In particular, the boundary value problems for $u^{[1]}$ and $u^{[2]}$ are $\begin{array}{l} \{\begin{matrix} L u^{[1]} = f & in D, \\ u^{[1]} + ε h \partial_{n} u^{[1]} = g + ε h \partial_{n} g & on \partial D, \end{matrix} \\ \{\begin{matrix} L u^{[2]} = f & in D, \\ u^{[2]} + ε h \partial_{n} u^{[2]} + \frac{ε^{2} h^{2}}{2} \partial_{n}^{2} u^{[2]} = g + ε h \partial_{n} g + \frac{ε^{2} h^{2}}{2} \partial_{n}^{2} g & on \partial D . \end{matrix} \end{array}$

The following theorem gives the approximation error of $v^{[n]}$ , whose proof is given in Appendix B.1.

Assumption 2.4.
Assume $D \subset D_{ε} \subset D_{ε_{0}}$ and $\partial D \in C^{\infty}$ . Let the operator $L$ given by (1.3) be strictly elliptic in $D_{ε_{0}}$ and have the coefficients $a^{i j}$ , $b^{i}$ , c belong to $C^{\infty} (\overline{D_{ε_{0}}})$ and $c ⩾ 0$ . Also assume $f, g \in C^{\infty} (\overline{D_{ε_{0}}})$ and $h \in C^{\infty} (\partial D)$ .
Theorem 2.5.
Under the Assumption 2.4, $\forall n, m ⩾ 0$ , $\begin{matrix} (2.18) & {‖ v^{[n]} - u_{ε} ‖}_{H^{m} (D)} = O (ε^{n + 1}) . \end{matrix}$

The following approximation error of $u^{[n]}$ has been proved in [8] for $n = 0, 1, 2$ , $\begin{matrix} {‖ u^{[n]} - u_{ε} ‖}_{H^{1} (D)} = O (ε^{n + 1}) . \end{matrix}$ Note that although $u^{[n]}$ and $v^{[n]}$ have the same approximation order, there might still be a considerable difference in the accuracy of their approximation errors due to the effects of the prefactors. The numerical results in [8] show that the approximation $u^{[n]}$ produces much less accurate results than $v^{[n]}$ for $n = 1, 2$ . This can be easily confirmed by the following simple one-dimensional example: $\begin{matrix} \{\begin{matrix} u_{ε}^{″} = 2 in D_{ε} = (0, 1 + ε), \\ u_{ε} (0) = u_{ε} (1 + ε) = 0 . \end{matrix} \end{matrix}$ The true solution is $u_{ε} (x) = x^{2} - (1 + ε) x$ . The equation for $u_{0}$ reads $\begin{matrix} \{\begin{matrix} u_{0}^{″} = 2 in D = (0, 1), \\ u_{0} (0) = u_{0} (1) = 0, \end{matrix} \end{matrix}$ with the solution $u_{0} (x) = x^{2} - x$ . Then the equation for $u_{1}$ is $\begin{matrix} \{\begin{matrix} u_{1}^{″} = 0 in D = (0, 1), \\ u_{1} (0) = 0, u_{1} (1) = - u_{0}^{'} (1) = - 1 . \end{matrix} \end{matrix}$ So $u_{1} (x) = - x$ , and then the partial sum $v^{[1]} (x) = u_{0} (x) + ε u_{1} (x) = x^{2} - x - ε x$ . Hence $\begin{matrix} v^{[1]} (x) - u_{ε} (x) = 0 . \end{matrix}$ The equation for $u^{[1]}$ is $\begin{matrix} \{\begin{matrix} {(u^{[1]})}^{″} = 2 in D = (0, 1), \\ u^{[1]} (0) = 0, u^{[1]} (1) + ε {(u^{[1]})}^{'} (1) = 0 . \end{matrix} \end{matrix}$ We find $u^{[1]} (x) = x^{2} - \frac{1 + 2 ε}{1 + ε} x$ , which is a worse approximation than $v^{[1]}$ since $\begin{matrix} u^{[1]} (x) - u_{ε} (x) = \frac{ε^{2}}{1 + ε} x = O (ε^{2}) . \end{matrix}$ To attain the zero error as $v^{[1]}$ , one needs to proceed to the next order $u^{[2]}$ by solving $\begin{matrix} \{\begin{matrix} {(u^{[2]})}^{″} = 2 in D = (0, 1), \\ u^{[2]} (0) = 0, u^{[2]} (1) + ε {(u^{[2]})}^{'} (1) + ε^{2} = 0 . \end{matrix} \end{matrix}$ It turns out $u^{[2]} (x) = v^{[1]} (x) = x^{2} - x - ε x$ .
3. The thin layer problem

Next, we generalize the above method from the continuous material coefficients $a (x)$ to the transmission problem associated with the piecewisely smooth coefficients. The Taylor expansion used in Section 2.1 is still applicable since we essentially apply the expansion on each subdomain where $a$ is smooth. The next step is to use Lemma 2.3 (the inverse Lax–Wendroff procedure) to convert the high order normal derivatives on the interface to the first order normal derivative and the tangential derivatives. This critical step facilitates the transmission condition given on the interface to build the connection between the solutions on each subdomain.

For ease of exposition, we only deal with the outward perturbation where $h (x) ⩾ 0$ for all $x \in \partial D$ . So D is a (proper) subset of $D_{ε}$ and the difference $D_{ε} ∖ \overline{D}$ is the thin layer $L_{ε}$ . The transmission condition is thus imposed on $\begin{matrix} Γ : = {x \in \partial D : h (x) > 0} \subset \partial D . \end{matrix}$ Note that $\overline{Γ} = \partial D \cap \partial L_{ε}$ and $\partial D ∖ Γ = {x \in \partial D : h (x) = 0} = \partial D \cap \partial D_{ε}$ . Assume that the second order coefficients $a^{i j}$ , $i, j = 1, \dots, d$ , are piecewisely smooth and have jumps only across the transmission interface Γ. In addition, the term f on the right-hand side of the equation is also allowed (but not necessarily) to have jumps on Γ. Specifically, we assume for $i, j = 1, \dots, d$ , $\begin{matrix} a^{i j} (x) = \{\begin{matrix} a_{int}^{i j} (x) & for x \in D \\ a_{ext}^{i j} (x) & for x \in L_{ε_{0}}, \end{matrix} f (x) = \{\begin{matrix} f_{int} (x) & for x \in D \\ f_{ext} (x) & for x \in L_{ε_{0}}, \end{matrix} \end{matrix}$ where $a_{int}^{i j}$ and $f_{int}$ are smooth functions on $\overline{D}$ while $a_{ext}^{i j}$ and $f_{ext}$ smooth on $\overline{L_{ε_{0}}}$ , and in general, $a_{int}^{i j} (x) \neq a_{ext}^{i j} (x)$ for $x \in Γ$ .

Write $\begin{matrix} u (x) = \{\begin{matrix} u_{int} (x) & for x \in D, \\ u_{ext} (x) & for x \in L_{ε}, \end{matrix} \end{matrix}$ then the transmission problem of our concern takes the form: $\begin{matrix} (3.1) & \{\begin{matrix} L u_{int} = f_{int} in D, \\ L u_{ext} = f_{ext} in L_{ε}, \\ u_{int} = u_{ext}, \sum_{i, j = 1}^{d} a_{int}^{i j} n_{i} \partial_{x_{j}} u_{int} = \sum_{i, j = 1}^{d} a_{ext}^{i j} n_{i} \partial_{x_{j}} u_{ext} on Γ, \\ u_{int} = g on \partial D_{ε} \cap \partial D, \\ u_{ext} = g on \partial D_{ε} \cap \partial L_{ε} . \end{matrix} \end{matrix}$

3.1. Asymptotic expansions in D and $L_{ε}$

Conceptually, we may first extend the domain of $u_{ext}$ to a fixed larger domain $L_{ε_{0}}$ , as in Section 2.1.1, and for simplicity we still use $u_{ext}$ for its extension. Assume the following two ansätze for $u_{int}$ and $u_{ext}$ respectively: $\begin{array}{l} (3.2) & u_{int} (x) = \sum_{n = 0}^{\infty} ε^{n} u_{int, n} (x) for x \in D, \\ (3.3) & u_{ext} (x) = \sum_{n = 0}^{\infty} ε^{n} u_{ext, n} (x) for x \in L_{ε_{0}} . \end{array}$ Plug these ansätze into (3.1), and match the terms at the same order of ε, then we obtain the following equations for $u_{int, n}$ and $u_{ext, n}$ , $\begin{matrix} L u_{int, n} = δ_{0, n} f_{int} in D, and q L u_{ext, n} = δ_{0, n} f_{ext} in L_{ε_{0}}, \end{matrix}$ and the transmission conditions on Γ for $u_{int, n}$ and $u_{ext, n}$ , $\begin{array}{l} (3.4) & u_{int, n} = u_{ext, n}, \\ (3.5) & \sum_{i, j = 1}^{d} a_{int}^{i j} n_{i} \partial_{x_{j}} u_{int, n} = \sum_{i, j = 1}^{d} a_{ext}^{i j} n_{i} \partial_{x_{j}} u_{ext, n} . \end{array}$ The boundary conditions on $\partial D_{ε} \cap \partial D$ for $u_{int, n}$ is $u_{int, n} = δ_{0, n} g$ and $u_{ext, n}$ share the same condition on $\partial D_{ε} \cap \partial L_{ε}$ .

Our goal is to derive the correct boundary conditions on $\partial D$ for $u_{int, n}$ . Note that we already have these conditions on $\partial D \cap \partial D_{ε}$ , thus it remains to find the boundary conditions on Γ for $u_{int, n}$ . To this end, we actually first derive the boundary conditions on Γ for $u_{ext, n}$ , and then convert $u_{ext, n}$ to $u_{int, n}$ by the transmission conditions (3.4) and (3.5).

To work on the exterior solution $u_{ext}$ , which behaves nicely in $L_{ε_{0}}$ , we apply the Taylor expansion method used in Section 2.1.2 to the ansätz (3.3) with the boundary condition $u_{ext} = g$ on $\partial D_{ε} \cap \partial L_{ε}$ . The obtained result is the following recursive expression of the boundary conditions on Γ for $u_{ext, n}$ : $\begin{matrix} (3.6) & \{\begin{matrix} u_{ext, 0} = g, \\ u_{ext, n} = \frac{h^{n}}{n!} \partial_{n}^{n} g - \sum_{k = 1}^{n} \frac{h^{k}}{k!} \partial_{n}^{k} u_{ext, n - k} \\ = \frac{h^{n}}{n!} \partial_{n}^{n} g - h \partial_{n} u_{ext, n - 1} \\ - \sum_{k = 2}^{n} \frac{h^{k}}{k!} F_{k, δ_{k, n} f_{ext}} [u_{ext, n - k}, \partial_{n} u_{ext, n - k}], \forall n ⩾ 1, \end{matrix} \end{matrix}$ where the operator $F_{k, f}$ is the operator $F_{k, Γ, L, f}$ introduced in Lemma 2.3 and the subindices Γ and $L$ are dropped for simplicity.

To handle the terms $\partial_{n} u_{ext, n - k}$ on the right-hand side of (3.6), we need the following lemma, proven in Appendix B.2.

Lemma 3.1.
For any integer $n ⩾ 0$ and any $x \in Γ$ , one can uniquely determine the value of the normal derivative $\partial_{n} u_{ext, n} (x)$ on Γ from the information of $u_{int, n}$ by using (3.4) and (3.5). More precisely, $\partial_{n} u_{ext, n} (x)$ for $x \in Γ$ only depends on
the normal vector $n = n (x)$ and

the value of $u_{int, n} (x^{'})$ for all $x^{'} \in Γ$ and

$\partial_{n} u_{int, n} (x)$ and

the second order coefficients $a_{int}^{i j} (x)$ and $a_{ext}^{i j} (x)$ , $i, j = 1, \dots, d$ .

Now the transmission conditions (3.4) and (3.5), serve the bridge from $u_{ext, n}$ to $u_{int, n}$ , with the aid of Lemma 3.1. Then the calculation following the procedure in the proof of Lemma 3.1 shows that (3.6) leads to the following final results for the boundary condition of ${u_{int, n}}$ on Γ: $\begin{matrix} (3.7) & \{\begin{matrix} u_{int, 0} = g, \\ u_{int, n} = \frac{h^{n}}{n!} \partial_{n}^{n} g - h q_{1} - \sum_{k = 2}^{n} \frac{h^{k}}{k!} F_{k, δ_{k, n} f_{ext}} [u_{int, n - k}, q_{k}], \forall n ⩾ 1, \end{matrix} \end{matrix}$ where for any $x \in \partial D$ , $\begin{array}{l} q_{k} (x) : = \frac{Q_{int} (n) \partial_{n} u_{int, n - k} + \sum_{i, j = 1}^{d} (a_{int}^{i j} - a_{ext}^{i j}) n_{i} \partial_{τ_{j}} u_{int, n - k}}{Q_{ext} (n)}, \\ Q_{ext} (n) : = \sum_{i, j = 1}^{d} a_{ext}^{i j} n_{i} n_{j}, Q_{int} (n) : = \sum_{i, j = 1}^{d} a_{int}^{i j} n_{i} n_{j} . \end{array}$
Remark 3.2.
Since we have $h = 0$ and $u_{int, n} = δ_{0, n} g$ on the boundary $\partial D ∖ Γ$ , the boundary conditions (3.7) also holds on $\partial D ∖ Γ$ and thus on the whole boundary $\partial D$ .

As an illuminating example, let us consider the elliptic operator $L = - \nabla \cdot (σ (x) \nabla)$ with a discontinuous $σ (x)$ , which has been studied in Example A.2 when $σ (x)$ is a smooth function. Example 3.3.
Set $g = 0$ and $L = - \nabla \cdot (σ (x) \nabla)$ . Assume $\begin{matrix} (3.8) & σ (x) = \{\begin{matrix} σ_{int} (x) & for x \in D, \\ σ_{ext} (x) & for x \in L_{ε_{0}}, \end{matrix} \end{matrix}$ where $σ_{int}$ and $σ_{ext}$ are smooth functions on $\overline{D}$ and $\overline{L_{ε_{0}}}$ respectively, and in general, they are distinct on the common boundary. To ensure the ellipticity of $L$ , we assume that $σ_{int}$ and $σ_{ext}$ are both positive everywhere in their domain. Then the transmission condition (3.5) reads $\begin{matrix} σ_{int} \partial_{n} u_{int, n} = σ_{ext} \partial_{n} u_{ext, n} . \end{matrix}$ Thus we deduce $\begin{matrix} (3.9) & \partial_{n} u_{ext, n} = \frac{σ_{int} \partial_{n} u_{int, n}}{σ_{ext}} . \end{matrix}$ Next, we compute explicitly the boundary conditions on Γ for the first three orders $u_{int, n}$ .

Order $n = 0$ . The boundary condition (3.7) on Γ for $u_{int, 0}$ is simply $\begin{matrix} (3.10) & u_{int, 0} = 0 . \end{matrix}$

Order $n = 1$ . The boundary condition (3.6) on Γ for $u_{ext, 1}$ reads $\begin{matrix} u_{ext, 1} (x) = - h (x) \partial_{n} u_{ext, 0} (x) . \end{matrix}$ Then by (3.4) and (3.9), we obtain the boundary condition on Γ for $u_{int, 1}$ : $\begin{matrix} (3.11) & u_{int, 1} = - \frac{h σ_{int}}{σ_{ext}} \partial_{n} u_{int, 0} . \end{matrix}$

Order $n = 2$ . Applying the boundary condition (A.8) for $u_{2}$ in Example A.2 to $u_{ext, 2}$ yields $\begin{matrix} u_{ext, 2} = - h \partial_{n} u_{ext, 1} + \frac{h^{2}}{2} (\frac{\partial_{n} σ_{ext}}{σ_{ext}} + κ) \partial_{n} u_{ext, 0} + \frac{h^{2} f_{ext}}{2 σ_{ext}} . \end{matrix}$ κ is the curvature of $\partial D$ , defined in Example A.2. Then substituting (3.4) and (3.9) into the last equation gives the boundary condition on Γ for $u_{int, 2}$ : $\begin{matrix} (3.12) & u_{int, 2} = - \frac{h σ_{int}}{σ_{ext}} \partial_{n} u_{int, 1} + \frac{h^{2} σ_{int}}{2 σ_{ext}} (\frac{\partial_{n} σ_{ext}}{σ_{ext}} + κ) \partial_{n} u_{int, 0} + \frac{h^{2} f_{ext}}{2 σ_{ext}} . \end{matrix}$

3.2. The approximate boundary conditions for the partial sums

Define the partial sums $\begin{matrix} v^{[n]} (x) = \{\begin{matrix} v_{int}^{[n]} (x) : = \sum_{k = 0}^{n} ε^{k} u_{int, k} (x) & for x \in D, \\ v_{ext}^{[n]} (x) : = \sum_{k = 0}^{n} ε^{k} u_{ext, k} (x) & for x \in L_{ε_{0}} . \end{matrix} \end{matrix}$ As in Section 2.2, the goal here is to derive the recursive boundary condition for the partial sums and to find the closed boundary value problems for the approximations $u^{[n]}$ .

To derive the boundary conditions that the partial sums $v_{int}^{[n]}$ satisfy, we have two equivalent approaches. The first one is to directly derive the boundary conditions for $v_{int}^{[n]}$ from the boundary conditions for $u_{int, n}$ which are already obtained above; the second approach is to apply (2.16) to $v_{ext}^{[n]}$ and then transfer to $v_{int}^{[n]}$ via the following transmission conditions $\begin{array}{l} v_{int}^{[n]} = v_{ext}^{[n]}, \\ \sum_{i, j = 1}^{d} a_{int}^{i j} n_{i} \partial_{x_{j}} v_{int}^{[n]} = \sum_{i, j = 1}^{d} a_{ext}^{i j} n_{i} \partial_{x_{j}} v_{ext}^{[n]}, \end{array}$ which can be easily deduced from (3.4) and (3.5). Let us continue to work on Example 3.3 to illustrate the first approach.

Order $n = 1$ . On the boundary $\partial D$ , (3.10) and (3.11) give $\begin{matrix} (3.13) & v_{int}^{[1]} = - ε \frac{h σ_{int}}{σ_{ext}} \partial_{n} u_{int, 0} = - \frac{ε h σ_{int}}{σ_{ext}} \partial_{n} v_{int}^{[0]} . \end{matrix}$ Then we are motivated to introduce the following Robin boundary value problem for $u^{[1]}$ : $\begin{matrix} (3.14) & \{\begin{matrix} L u^{[1]} = f & in D, \\ u^{[1]} + \frac{ε h σ_{int}}{σ_{ext}} \partial_{n} u^{[1]} = 0 & on \partial D . \end{matrix} \end{matrix}$

Order $n = 2$ . On the boundary $\partial D$ , (3.12) and (3.13) show $\begin{matrix} v_{int}^{[2]} & = - \frac{ε h σ_{int}}{σ_{ext}} \partial_{n} v_{int}^{[0]} - \frac{ε^{2} h σ_{int}}{σ_{ext}} \partial_{n} u_{int, 1} + \frac{ε^{2} h^{2} σ_{int}}{2 σ_{ext}} (\frac{\partial_{n} σ_{ext}}{σ_{ext}} + κ) \partial_{n} u_{int, 0} + \frac{ε^{2} h^{2} f_{ext}}{2 σ_{ext}} \\ = - \frac{ε h σ_{int}}{σ_{ext}} \partial_{n} v_{int}^{[1]} + \frac{ε^{2} h^{2} σ_{int}}{2 σ_{ext}} (\frac{\partial_{n} σ_{ext}}{σ_{ext}} + κ) \partial_{n} v_{int}^{[0]} + \frac{ε^{2} h^{2} f_{ext}}{2 σ_{ext}}, \end{matrix}$ thus the closed Robin boundary value problem for $u^{[2]}$ can be imposed as: $\begin{matrix} (3.15) & \{\begin{matrix} L u^{[2]} = f & in D, \\ u^{[2]} + \frac{ε h σ_{int}}{σ_{ext}} \partial_{n} u^{[2]} - \frac{ε^{2} h^{2} σ_{int}}{2 σ_{ext}} (\frac{\partial_{n} σ_{ext}}{σ_{ext}} + κ) \partial_{n} u^{[2]} = \frac{ε^{2} h^{2} f_{ext}}{2 σ_{ext}} & on \partial D . \end{matrix} \end{matrix}$

To compare with the results derived in [8] where the coefficient σ is piecewise constant, we set $σ_{int} = σ_{0}$ and $σ_{ext} = 1$ and $h (x) > 0$ for all $x \in \partial D$ . Then (3.14) becomes $\begin{matrix} \{\begin{matrix} L u^{[1]} = f & in D, \\ u^{[1]} + ε h σ_{0} \partial_{n} u^{[1]} = 0 & on \partial D, \end{matrix} \end{matrix}$ which is the same as that in [8]. The equation (3.15) becomes $\begin{matrix} \{\begin{matrix} L u^{[2]} = f & in D, \\ u^{[2]} + ε h σ_{0} (1 - \frac{ε κ h}{2}) \partial_{n} u^{[2]} = \frac{ε^{2} h^{2} f}{2} & on \partial D . \end{matrix} \end{matrix}$ Multiplying the boundary condition for $u^{[2]}$ by $(1 + \frac{ε κ h}{2})$ yields $\begin{matrix} (1 + \frac{ε κ h}{2}) u^{[2]} + ε h σ_{0} (1 - \frac{ε^{2} κ^{2} h^{2}}{2}) \partial_{n} u^{[2]} = (1 + \frac{ε κ h}{2}) \frac{ε^{2} h^{2} f}{2}, \end{matrix}$ that is, $\begin{matrix} (1 + \frac{ε κ h}{2}) u^{[2]} + ε h σ_{0} \partial_{n} u^{[2]} = \frac{ε^{2} h^{2} f}{2} + O (ε^{3}) . \end{matrix}$ By neglecting the third order term $O (ε^{3})$ , we have the same equation in [8] for $u^{[2]}$ .

4. The thin layer problem with high-contrast ratio

From this section, we take into account of the contrast ratio parameter σ together with the geometric perturbation parameter ε. This section considers the following transmission problem on $D_{ε}$ : $\begin{matrix} (4.1) & \{\begin{matrix} - Δ u_{int} = f_{int} & in D, \\ - σ Δ u_{ext} = f_{ext} & in L_{ε}, \\ u_{int} = u_{ext} \partial_{n} u_{int} = σ \partial_{n} u_{ext}, & on Γ, \\ u_{int} = g & on \partial D_{ε} \cap \partial D, \\ u_{ext} = g & on \partial D_{ε} \cap \partial L_{ε}, \end{matrix} \end{matrix}$ where σ is a positive constant. The geometry of the domains are exactly the same as in Section 3, i.e., $D \subset D_{ε}$ and $L_{ε} = D_{ε} ∖ \overline{D}$ . Γ is the interface separating two materials with different conductivity. A large σ means a large conductivity in the thin layer $L_{ε}$ and a small σ means a (relatively) large conductivity in the interior D.

We want to investigate the limiting behavior, as well as the asymptotic expansions, of the interior solution $u_{int}$ as $ε \to 0$ and $σ \to 0$ or $σ \to \infty$ . Before we present the abstract analysis, let us first heuristically show how three scaling regimens can appear by considering a simple 1D example.

Example 4.1.
Let $D = (0, 1)$ , $D_{ε} = (- ε h_{0}, 1 + ε h_{1})$ with two numbers $h_{0}, h_{1} ⩾ 0$ , and take $f_{int} = f_{ext} = - 2$ and $g = 1$ . Then it is easy to find the interior solution is $\begin{matrix} (4.2) & \begin{matrix} u_{int} (x) = & x^{2} - A x + B - \frac{h_{0} h_{1} ε}{σ} C, \end{matrix} \end{matrix}$ and the exterior solution is $\begin{matrix} (4.3) & u_{ext} (x) = \{\begin{matrix} \frac{x^{2} - A x}{σ} + B - \frac{h_{0} h_{1} ε}{σ} C, & - h_{0} ε ⩽ x ⩽ 0, \\ \frac{x^{2} - A x + A - 1}{σ} + (1 - A + B) - \frac{h_{0} h_{1} ε}{σ} C, & 1 ⩽ x ⩽ 1 + h_{1} ε, \end{matrix} \end{matrix}$ where $\begin{matrix} A = \frac{σ + 2 h_{1} ε + (h_{1}^{2} - h_{0}^{2}) ε^{2}}{σ + (h_{0} + h_{1}) ε}, B = \frac{σ + h_{1} ε - h_{0}^{2} ε^{2}}{σ + (h_{0} + h_{1}) ε}, C = \frac{2 ε + (h_{0} + h_{1}) ε^{2}}{σ + (h_{0} + h_{1}) ε} . \end{matrix}$

The limiting behavior of the interior solution (4.2) and exterior solution (4.3) for this example is different in the following three cases

$ε / σ \to 0$ ,

$σ / ε \to 0$ ,

$ε / σ \to c$ , where $c \in (0, \infty)$ .

In Case (i), as ε and $μ : = ε / σ$ tend to 0, we have the interior solution (4.2) $u_{int} (x) \to x^{2} - x + 1 \sim O (1)$ , and the exterior solution $u_{ext} (x) \sim O (σ^{- 1}) + O (1) + O (μ) = O (μ / ε) + O (1)$ .

In Case (ii), introduce $λ : = σ / ε$ , then both ε and λ go to 0, If $h_{0} h_{1} > 0$ , i.e., the domain perturbation is applied to the whole boundary $\partial D$ , then $u_{int}$ is at the order $O (λ^{- 1})$ ; otherwise, one has $h_{0} = 0$ or $h_{1} = 0$ , and so $u_{int} (x) \sim O (1)$ . In both circumstances, $u_{ext}$ is at the order $O (σ^{- 1}) = O (ε^{- 1} λ^{- 1})$ .

In Case (iii), as $ε \to 0$ and $ε / σ \to c$ , $u_{int} (x) \sim O (1)$ and $u_{ext}$ is at the order $O (σ^{- 1}) = O (ε^{- 1})$ .

For general problems, the scalings of the magnitudes of $u_{int}$ and $u_{ext}$ behave exactly the same as in the above example. In the next, we develop the two-parameter asymptotic analysis for the general transmission problem (4.1) by discussing the above three cases. The results we obtained below are written recursively up to any order in an abstract way. The readers can find explicit boundary conditions and solvability conditions for some lower order terms for each case in Appendix C.
4.1. Case (i): $ε / σ \to 0$ , $ε \to 0$

We now treat ε and $\begin{matrix} μ = ε / σ \end{matrix}$ as independent small parameters. Introduce the rescaled exterior solution ${\tilde{u}}_{ext} = ε u_{ext}$ , then rewrite the original equation (4.1) in terms of $u_{int}$ and ${\tilde{u}}_{ext}$ : $\begin{matrix} (4.4) & \{\begin{matrix} - Δ u_{int} = f_{int} & in D, \\ - Δ {\tilde{u}}_{ext} = μ f_{ext} & in L_{ε}, \\ {\tilde{u}}_{ext} = ε u_{int}, \partial_{n} {\tilde{u}}_{ext} = μ \partial_{n} u_{int} & on Γ, \\ u_{int} = g & on \partial D_{ε} \cap \partial D, \\ {\tilde{u}}_{ext} = ε g & on \partial D_{ε} \cap \partial L_{ε} . \end{matrix} \end{matrix}$ Assume $u_{int}$ and ${\tilde{u}}_{ext}$ have double asymptotic expansions $\begin{matrix} u_{int} (x) = \sum_{m, n = 0}^{\infty} u_{int, m, n} (x) ε^{m} μ^{n}, x \in D, \\ {\tilde{u}}_{ext} (x) = \sum_{m, n = 0}^{\infty} {\tilde{u}}_{ext, m, n} (x) ε^{m} μ^{n}, x \in L_{ε_{0}} . \end{matrix}$ After substituting these into (4.4) and equating terms of each pair of powers of ε and μ, we get the following results: $\begin{array}{l} - Δ u_{int, m, n} = δ_{0, m} δ_{0, n} f_{int} in D, \\ - Δ {\tilde{u}}_{ext, m, n} = δ_{0, m} δ_{1, n} f_{ext} in L_{ε_{0}}, \\ {\tilde{u}}_{ext, m, n} = \{\begin{matrix} u_{int, m - 1, n}, & m ⩾ 1 \\ 0, & m = 0 \end{matrix} on Γ, \\ \partial_{n} {\tilde{u}}_{ext, m, n} = \{\begin{matrix} \partial_{n} u_{int, m, n - 1}, & n ⩾ 1 \\ 0, & n = 0 \end{matrix} on Γ, \\ u_{int, m, n} = δ_{0, m} δ_{0, n} g on \partial D_{ε} \cap \partial D . \end{array}$ For the boundary condition ${\tilde{u}}_{ext} = ε g$ on $\partial D_{ε} \cap \partial L_{ε}$ , applying the Taylor expansion method as in Section 2 and Section 3 yields the following recursive boundary conditions on Γ for ${\tilde{u}}_{ext, m, n}$ : $\begin{matrix} (4.5) & \{\begin{matrix} {\tilde{u}}_{ext, 0, n} = 0, \\ {\tilde{u}}_{ext, m, n} = δ_{0, n} \frac{h^{m - 1}}{(m - 1)!} \partial_{n}^{m - 1} g - \sum_{k = 1}^{m} \frac{h^{k}}{k!} \partial_{n}^{k} {\tilde{u}}_{ext, m - k, n} \\ = δ_{0, n} \frac{h^{m - 1}}{(m - 1)!} \partial_{n}^{m - 1} g - h \partial_{n} {\tilde{u}}_{ext, m - 1, n} \\ - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, δ_{k, m} δ_{0, n} f_{ext}} [{\tilde{u}}_{ext, m - k, n}, \partial_{n} {\tilde{u}}_{ext, m - k, n}], \forall m ⩾ 1 . \end{matrix} \end{matrix}$

Next, we transform these boundary conditions on Γ for ${\tilde{u}}_{ext, m, n}$ into those for $u_{int, m, n}$ . One has on the interface Γ $\begin{matrix} u_{int, m, n} = {\tilde{u}}_{ext, m + 1, n} = & δ_{0, n} \frac{h^{m}}{m!} \partial_{n}^{m} g - h \partial_{n} {\tilde{u}}_{ext, m, n} \\ - \sum_{k = 2}^{m + 1} \frac{h^{k}}{k!} F_{k, δ_{k, m + 1} δ_{1, n} f_{ext}} [{\tilde{u}}_{ext, m + 1 - k, n}, \partial_{n} {\tilde{u}}_{ext, m + 1 - k, n}] . \end{matrix}$ Thus for $m = 0$ , we have on Γ $\begin{matrix} u_{int, 0, n} = δ_{0, n} g - h \partial_{n} {\tilde{u}}_{ext, 0, n} = \{\begin{matrix} g, & n = 0, \\ - h \partial_{n} u_{int, 0, n - 1}, & n ⩾ 1; \end{matrix} \end{matrix}$ and for $m ⩾ 1$ and $n = 0$ , on Γ $\begin{matrix} u_{int, m, 0} = & \frac{h^{m}}{m!} \partial_{n}^{m} g - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, 0} [u_{int, m - k, 0}, 0] - \frac{h^{m + 1}}{(m + 1)!} F_{m + 1, 0} [0, 0] \\ = & \frac{h^{m}}{m!} \partial_{n}^{m} g - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, 0} [u_{int, m - k, 0}, 0] . \end{matrix}$ Note that here we used the trivial fact $F_{m + 1, 0} [0, 0] = 0$ by definition.

For $m, n ⩾ 1$ , on Γ $\begin{matrix} u_{int, m, n} = & - h \partial_{n} u_{int, m, n - 1} - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, 0} [u_{int, m - k, n}, \partial_{n} u_{int, m + 1 - k, n - 1}] \\ - \frac{h^{m + 1}}{(m + 1)!} F_{m + 1, δ_{1, n} f_{ext}} [0, \partial_{n} u_{int, 0, n - 1}] . \end{matrix}$

4.2. Case (ii): $σ / ε \to 0$ , $ε \to 0$

Now both ε and $\begin{matrix} λ : = σ / ε \end{matrix}$ are small parameters. Introduce ${\hat{u}}_{ext} = σ u_{ext} = λ ε u_{ext}$ . Then (4.1) becomes $\begin{matrix} \{\begin{matrix} - Δ u_{int} = f_{int} & in D, \\ - Δ {\hat{u}}_{ext} = f_{ext} & in L_{ε}, \\ {\hat{u}}_{ext} = ε λ u_{int} \partial_{n} {\hat{u}}_{ext} = \partial_{n} u_{int} & on Γ, \\ u_{int} = g & on \partial D_{ε} \cap \partial D, \\ {\hat{u}}_{ext} = ε λ g & on \partial D_{ε} \cap \partial L_{ε} . \end{matrix} \end{matrix}$ We have to further study two subcases and treat them separately.

4.2.1. Case (ii)₁: $Γ \neq \partial D$ , or $\partial D \cap \partial D_{ε} \neq \emptyset$

This means the domain perturbation is only applied to a proper subset Γ of the boundary $\partial D$ .

Assume the double asymptotic expansions $\begin{matrix} u_{int} (x) = \sum_{m, n = 0}^{\infty} u_{int, m, n} (x) ε^{m} λ^{n}, x \in D, \\ {\hat{u}}_{ext} (x) = \sum_{m, n = 0}^{\infty} {\hat{u}}_{ext, m, n} (x) ε^{m} λ^{n}, x \in L_{ε_{0}} . \end{matrix}$ Substituting these into (4.1) and equating terms of each pair of powers of ε and λ, we find that $\begin{array}{l} - Δ u_{int, m, n} = δ_{0, m} δ_{0, n} f_{int} in D, \\ - Δ {\hat{u}}_{ext, m, n} = δ_{0, m} δ_{0, n} f_{ext} in L_{ε_{0}}, \\ {\hat{u}}_{ext, m, n} = \{\begin{matrix} u_{int, m - 1, n - 1}, & m, n ⩾ 1 \\ 0, & otherwise \end{matrix} on Γ, \\ \partial_{n} {\hat{u}}_{ext, m, n} = \partial_{n} u_{int, m, n} on Γ, \\ u_{int, m, n} = δ_{0, m} δ_{0, n} g on \partial D_{ε} \cap \partial D . \end{array}$ Applying the Taylor expansion method to the boundary condition ${\hat{u}}_{ext} = ε λ g$ on $\partial D_{ε} \cap \partial L_{ε}$ , we obtain the following recursive boundary conditions on Γ for ${\hat{u}}_{ext, m, n}$ : $\begin{matrix} (4.6) & \{\begin{matrix} {\hat{u}}_{ext, 0, n} = 0, \\ {\hat{u}}_{ext, m, n} = δ_{1, n} \frac{h^{m - 1}}{(m - 1)!} \partial_{n}^{m - 1} g - \sum_{k = 1}^{m} \frac{h^{k}}{k!} \partial_{n}^{k} {\hat{u}}_{ext, m - k, n} \\ = δ_{1, n} \frac{h^{m - 1}}{(m - 1)!} \partial_{n}^{m - 1} g - h \partial_{n} {\hat{u}}_{ext, m - 1, n} \\ - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, δ_{k, m} δ_{0, n} f_{ext}} [{\hat{u}}_{ext, m - k, n}, \partial_{n} {\hat{u}}_{ext, m - k, n}], \forall m ⩾ 1 . \end{matrix} \end{matrix}$

Next, we convert these boundary conditions for ${\hat{u}}_{ext, m, n}$ into those for $u_{int, m, n}$ . It turns out that the Neumann boundary condition on Γ appears in this case. For $m = 0$ , on Γ $\begin{matrix} u_{int, 0, n} = {\hat{u}}_{ext, 1, n + 1} = δ_{0, n} g - h \partial_{n} {\hat{u}}_{ext, 0, n + 1} = δ_{0, n} g - h \partial_{n} u_{int, 0, n + 1}, \end{matrix}$ thus we get for $n ⩾ 1$ , $\begin{matrix} \partial_{n} u_{int, 0, n} = \frac{δ_{1, n} g}{h} - \frac{1}{h} u_{int, 0, n - 1}; \end{matrix}$ moreover, we have for $m = n = 0$ , on Γ $\begin{matrix} \partial_{n} u_{int, 0, 0} = \partial_{n} {\hat{u}}_{ext, 0, 0} = - \frac{1}{h} {\hat{u}}_{ext, 1, 0} = 0 . \end{matrix}$ For $m ⩾ 1$ and $n = 0$ , one has on Γ, $\begin{matrix} \partial_{n} u_{int, m, 0} = & \partial_{n} {\hat{u}}_{ext, m, 0} \\ = & - \frac{1}{h} {\hat{u}}_{ext, m + 1, 0} - \sum_{k = 2}^{m + 1} \frac{h^{k - 1}}{k!} F_{k, δ_{k, m + 1} f_{ext}} [{\hat{u}}_{ext, m + 1 - k, 0}, \partial_{n} {\hat{u}}_{ext, m + 1 - k, 0}] \\ = & - \sum_{k = 2}^{m + 1} \frac{h^{k - 1}}{k!} F_{k, δ_{k, m + 1} f_{ext}} [0, \partial_{n} u_{int, m + 1 - k, 0}]; \end{matrix}$ and for $m, n ⩾ 1$ , on Γ $\begin{matrix} \partial_{n} u_{int, m, n} = & \partial_{n} {\hat{u}}_{ext, m, n} \\ = & δ_{1, n} \frac{h^{m - 1}}{m!} \partial_{n}^{m} g - \frac{1}{h} {\hat{u}}_{ext, m + 1, n} - \sum_{k = 2}^{m + 1} \frac{h^{k - 1}}{k!} F_{k, 0} [{\hat{u}}_{ext, m + 1 - k, n}, \partial_{n} {\hat{u}}_{ext, m + 1 - k, n}] \\ = & δ_{1, n} \frac{h^{m - 1}}{m!} \partial_{n}^{m} g - \frac{1}{h} u_{int, m, n - 1} - \sum_{k = 2}^{m} \frac{h^{k - 1}}{k!} F_{k, 0} [u_{int, m - k, n - 1}, \partial_{n} u_{int, m + 1 - k, n}] \\ - \frac{h^{m}}{(m + 1)!} F_{m + 1, 0} [0, \partial_{n} u_{int, 0, n}] . \end{matrix}$

Note that the boundary conditions on $\partial D$ for $u_{int, m, n}$ are the mixture of the Neumann conditions on Γ and the Dirichlet conditions $u_{int, m, n} = δ_{0, m} δ_{0, n} g$ on $\partial D ∖ Γ$ .

4.2.2. Case (ii)₂: $Γ = \partial D$ , or $\partial D \cap \partial D_{ε} = \emptyset$

In this case, the domain perturbation is applied to the whole boundary $\partial D$ . It turns out that $u_{int}$ is at the order $O (λ^{- 1})$ . So we assume $\begin{matrix} u_{int} (x) = \sum_{m = 0}^{\infty} \sum_{n = - 1}^{\infty} u_{int, m, n} (x) ε^{m} λ^{n}, x \in D . \end{matrix}$ Consequently, the transmission conditions on $Γ = \partial D$ become $\begin{array}{l} {\hat{u}}_{ext, m, n} = \{\begin{matrix} u_{int, m - 1, n - 1}, & m ⩾ 1 \\ 0, & m = 0 \end{matrix} on \partial D, \\ \partial_{n} u_{int, m, n} = \{\begin{matrix} \partial_{n} {\hat{u}}_{ext, m, n}, & n ⩾ 0 \\ 0, & n = - 1 \end{matrix} on \partial D . \end{array}$ In addition, (4.6) still holds. We already have $\partial_{n} u_{int, m, - 1} = 0$ ; on $\partial D$ , and for $n ⩾ 0$ , on $\partial D$ $\begin{matrix} \partial_{n} u_{int, m, n} = & \partial_{n} {\hat{u}}_{ext, m, n} \\ = & δ_{1, n} \frac{h^{m - 1}}{m!} \partial_{n}^{m} g - \frac{1}{h} {\hat{u}}_{ext, m + 1, n} \\ - \sum_{k = 2}^{m + 1} \frac{h^{k - 1}}{k!} F_{k, δ_{k, m + 1} δ_{0, n} f_{ext}} [{\hat{u}}_{ext, m + 1 - k, n}, \partial_{n} {\hat{u}}_{ext, m + 1 - k, n}] . \end{matrix}$ Thus on $\partial D$ , one has for $m = 0$ , $n ⩾ 0$ , $\begin{matrix} \partial_{n} u_{int, 0, n} = \frac{δ_{1, n} g}{h} - \frac{1}{h} {\hat{u}}_{ext, 1, n} = \frac{δ_{1, n} g}{h} - \frac{1}{h} u_{int, 0, n - 1}; \end{matrix}$ and for $m ⩾ 1$ , $n ⩾ 0$ , $\begin{matrix} \partial_{n} u_{int, m, n} = & δ_{1, n} \frac{h^{m - 1}}{m!} \partial_{n}^{m} g - \frac{1}{h} u_{int, m, n - 1} - \frac{h^{m}}{(m + 1)!} F_{m + 1, δ_{0, n} f_{ext}} [0, \partial_{n} u_{int, 0, n}] \\ - \sum_{k = 2}^{m} \frac{h^{k - 1}}{k!} F_{k, 0} [u_{int, m - k, n - 1}, \partial_{n} u_{int, m + 1 - k, n}] . \end{matrix}$

The above Neumann boundary value problems for $u_{int, m, n}$ are not well-posed, since the solution to the Poisson equation with pure Neumann boundary condition $\begin{matrix} \{\begin{matrix} - Δ u = f, & in D, \\ \partial_{n} u = g, & on \partial D, \end{matrix} \end{matrix}$ can only be determined up to constant. However, note that a necessary condition for the existence of a solution to the Neumann problem is $\begin{matrix} \int_{\partial D} g = - \int_{D} f . \end{matrix}$ Applying this solvability condition to the Neumann problem for $u_{int, m, n + 1}$ leads to an additional boundary integral condition for $u_{int, m, n}$ . Specifically, the following solvability conditions can uniquely determine $u_{int, m, n}$ : $\begin{matrix} \int_{\partial D} \frac{u_{int, 0, n}}{h} = \int_{D} δ_{- 1, n} f_{int} + \int_{\partial D} \frac{δ_{0, n} g}{h}, n ⩾ 0, \end{matrix}$ and for $m ⩾ 1$ , $n ⩾ 0$ , $\begin{matrix} \int_{\partial D} \frac{u_{int, m, n}}{h} = & \int_{\partial D} δ_{0, n} \frac{h^{m - 1}}{m!} \partial_{n}^{m} g \\ - \sum_{k = 2}^{m} \int_{\partial D} \frac{h^{k - 1}}{k!} F_{k, 0} [u_{int, m - k, n}, \partial_{n} u_{int, m + 1 - k, n + 1}] \\ - \int_{\partial D} \frac{h^{m}}{(m + 1)!} F_{m + 1, δ_{- 1, n} f_{ext}} [0, \partial_{n} u_{int, 0, n + 1}] . \end{matrix}$

4.3. Case (iii): $ε / σ \to c \in (0, \infty)$ , $ε \to 0$

For this case, we introduce the small parameter $\begin{matrix} θ : = \frac{ε}{σ} - c, \end{matrix}$ and also rescale the exterior solution ${\tilde{u}}_{ext} = ε u_{ext}$ as in Case (i). Plugging the ansätz $\begin{matrix} u_{int} (x) = \sum_{m, n = 0}^{\infty} u_{int, m, n} (x) ε^{m} θ^{n} x \in D, \\ {\tilde{u}}_{ext} (x) = \sum_{m, n = 0}^{\infty} {\tilde{u}}_{ext, m, n} (x) ε^{m} θ^{n} x \in L_{ε_{0}} \end{matrix}$ into (4.4) yields the following $\begin{array}{l} (4.7) & - Δ u_{int, m, n} = δ_{0, m} δ_{0, n} f_{int} in D, \\ (4.8) & - Δ {\tilde{u}}_{ext, m, n} = δ_{0, m} δ_{0, n} c f_{ext} + δ_{0, m} δ_{1, n} f_{ext} in L_{ε_{0}}, \\ (4.9) & {\tilde{u}}_{ext, m, n} = \{\begin{matrix} u_{int, m - 1, n}, & m ⩾ 1 \\ 0, & m = 0 \end{matrix} on Γ, \\ (4.10) & \partial_{n} {\tilde{u}}_{ext, m, n} = \{\begin{matrix} \partial_{n} u_{int, m, n - 1} + c \partial_{n} u_{int, m, n}, & n ⩾ 1 \\ c \partial_{n} u_{int, m, n}, & n = 0 \end{matrix} on Γ, \\ (4.11) & u_{int, m, n} = δ_{0, m} δ_{0, n} g on \partial D_{ε} \cap \partial D . \end{array}$ From the boundary condition ${\tilde{u}}_{ext} = ε g$ on $\partial D_{ε} \cap \partial L_{ε}$ , the recursive boundary conditions on Γ for ${\tilde{u}}_{ext, m, n}$ are derived in (4.5). The derivation of the boundary conditions of $u_{int, m, n}$ from those of ${\tilde{u}}_{ext, m, n}$ is below.

For $m = 0$ , one has on Γ $\begin{matrix} u_{int, 0, n} = {\tilde{u}}_{ext, 1, n} = & δ_{0, n} g - h \partial_{n} {\tilde{u}}_{ext, 0, n} \\ = & \{\begin{matrix} g - c h \partial_{n} u_{int, 0, 0}, & n = 0, \\ - c h \partial_{n} u_{int, 0, n} - h \partial_{n} u_{int, 0, n - 1}, & n ⩾ 1 . \end{matrix} \end{matrix}$ Thus we obtain the following Robin boundary conditions $\begin{matrix} u_{int, 0, n} + c h \partial_{n} u_{int, 0, n} = \{\begin{matrix} g, & n = 0, \\ - h \partial_{n} u_{int, 0, n - 1}, & n ⩾ 1 . \end{matrix} \end{matrix}$

For $m ⩾ 1$ and $n = 0$ , on Γ, $\begin{matrix} u_{int, m, 0} = & {\tilde{u}}_{ext, m + 1, 0} \\ = & \frac{h^{m}}{m!} \partial_{n}^{m} g - h \partial_{n} {\tilde{u}}_{ext, m, 0} - \sum_{k = 2}^{m + 1} \frac{h^{k}}{k!} F_{k, δ_{k, m + 1} c f_{ext}} [{\tilde{u}}_{ext, m + 1 - k, 0}, \partial_{n} {\tilde{u}}_{ext, m + 1 - k, 0}] \\ = & \frac{h^{m}}{m!} \partial_{n}^{m} g - c h \partial_{n} u_{int, m, 0} - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, 0} [u_{int, m - k, 0}, c \partial_{n} u_{int, m + 1 - k, 0}] \\ - \frac{h^{m + 1}}{(m + 1)!} F_{m + 1, c f_{ext}} [0, c \partial_{n} u_{int, 0, 0}], \end{matrix}$ hence the Robin boundary condition on Γ is $\begin{matrix} u_{int, m, 0} + c h \partial_{n} u_{int, m, 0} = & \frac{h^{m}}{m!} \partial_{n}^{m} g - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, 0} [u_{int, m - k, 0}, c \partial_{n} u_{int, m + 1 - k, 0}] \\ - \frac{h^{m + 1}}{(m + 1)!} F_{m + 1, c f_{ext}} [0, c \partial_{n} u_{int, 0, 0}] . \end{matrix}$

For $m, n ⩾ 1$ , on Γ, we have $\begin{matrix} u_{int, m, n} = & {\tilde{u}}_{ext, m + 1, n} \\ = & - h \partial_{n} {\tilde{u}}_{ext, m, n} - \sum_{k = 2}^{m + 1} \frac{h^{k}}{k!} F_{k, δ_{k, m + 1} δ_{1, n} f_{ext}} [{\tilde{u}}_{ext, m + 1 - k, n}, \partial_{n} {\tilde{u}}_{ext, m + 1 - k, n}] \\ = & - c h \partial_{n} u_{int, m, n} - h \partial_{n} u_{int, m, n - 1} \\ - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, 0} [u_{int, m - k, n}, c \partial_{n} u_{int, m + 1 - k, n} + \partial_{n} u_{int, m + 1 - k, n - 1}] \\ - \frac{h^{m + 1}}{(m + 1)!} F_{m + 1, δ_{1, n} f_{ext}} [0, c \partial_{n} u_{int, 0, n} + \partial_{n} u_{int, 0, n - 1}], \end{matrix}$ and thus the Robin boundary condition on Γ is $\begin{matrix} u_{int, m, n} + c h \partial_{n} u_{int, m, n} \\ = - h \partial_{n} u_{int, m, n - 1} - \frac{h^{m + 1}}{(m + 1)!} F_{m + 1, δ_{1, n} f_{ext}} [0, \partial_{n} u_{int, 0, n - 1} + c \partial_{n} u_{int, 0, n}] \\ - \sum_{k = 2}^{m} \frac{h^{k}}{k!} F_{k, 0} [u_{int, m - k, n}, \partial_{n} u_{int, m + 1 - k, n - 1} + c \partial_{n} u_{int, m + 1 - k, n}] . \end{matrix}$

To summarize the above three cases, we find that the limit $ε / σ \to c$ is quite important: the value of c determines the type of the boundary conditions in the asymptotic series. $c = 0$ means ε decays faster than σ or σ is not a small value, and our result shows that the boundary conditions for the asymptotic expansions remain the Dirichlet type. $c = \infty$ corresponds to a very small conductivity in the exterior layer, and in this case, it is interesting to see the Neumann conditions on $\partial D$ for all terms in the asymptotic expansions. The case of $c \in (0, \infty)$ that leads to the Robin boundary conditions can be regarded as between the above two extreme cases.

5. Asymptotic expansion for the perturbed interface problem

The previous sections on the interface problem assume that the interface is the boundary of the fixed domain D. The geometric perturbation is only applied to the outside layer. In this section, we focus on the situation where the interface is perturbed. The setting is the following. Assume D is a smooth bounded domain and is partitioned into two subdomains separated by an interface $Γ_{ε}$ : $\begin{matrix} (5.1) & D = D_{ε}^{+} \cup D_{ε}^{-} \cup Γ_{ε} . \end{matrix}$ $Γ_{ε} = \partial D_{ε}^{-} \cap \partial D_{ε}^{+}$ is assumed smooth. The interface $Γ_{ε}$ is modelled in a perturbative way. Assume there is a fixed interface Γ and let $n (x)$ be the unit normal vector on Γ pointing outward of $D^{-}$ . That is, the whole domain D has a fixed decomposition $D = D^{-} \cup D^{+} \cup Γ$ . Then we define $Γ_{ε}$ for $ε < ε_{0}$ $\begin{matrix} (5.2) & Γ_{ε} = {x^{'} : x^{'} = x + ε h (x) n (x), x \in Γ} . \end{matrix}$

We consider the following interface problem on D with transmission condition on the interface $Γ_{ε}$ : $\begin{matrix} (5.3) & \{\begin{matrix} - \nabla \cdot (σ^{\pm} (x) \nabla u_{ε}^{\pm} (x)) = f (x) in D_{ε}^{\pm}, \\ u_{ε}^{+} (x) = u_{ε}^{-} (x), σ^{+} (x) \partial_{n_{ε}} u_{ε}^{+} (x) = σ^{-} (x) \partial_{n_{ε}} u_{ε}^{-} (x) on Γ_{ε}, \\ u_{ε}^{\pm} = g on \partial D \cap \partial D_{ε}^{\pm}, \end{matrix} \end{matrix}$ where $σ^{\pm} (x) > 0$ for every $x \in D_{ε}^{\pm}$ , and $n_{ε} (x)$ is the unit normal vector on $Γ_{ε}$ pointing outward of $D_{ε}^{-}$ . Denote $u_{ε}$ restricted on $D_{ε}^{+}$ and $D_{ε}^{-}$ by $u_{ε}^{+}$ and $u_{ε}^{-}$ , respectively. For this interface problem (5.3), the variational formulation reads as follows: Seek $u_{ε} \in H_{0}^{1} (D)$ such that $\begin{matrix} (5.4) & \int_{D_{ε}^{+}} σ^{+} \nabla u_{ε}^{+} \cdot \nabla v d x + \int_{D_{ε}^{-}} σ^{-} \nabla u_{ε}^{-} \cdot \nabla v d x = \int_{D} f v d x, \forall v \in H_{0}^{1} (D) . \end{matrix}$ We assume that $σ^{\pm}$ are defined on sufficiently large domains such that for every sufficiently small ε, $σ^{\pm} \in C^{\infty} (D_{ε}^{\pm})$ . We also assume $f \in C^{\infty} (D)$ . Then, $u_{ε}^{\pm} \in C^{\infty} (\overline{D_{ε}^{\pm}})$ .

Assume the coefficient $σ (x)$ is the piecewise homogeneous case: $\begin{matrix} (5.5) & σ (x) = \{\begin{matrix} 1, & x \in D_{ε}^{+}, \\ σ, & x \in D_{ε}^{-}, \end{matrix} \end{matrix}$ where σ is a positive constant. We are interested in the high-contrast ratio limit, which corresponds to a very small or very large value of σ.

[15] has studied the first order and second order perturbations to the problem (5.3) by the method of shape calculus for small ε. The second order approximation was obtained by considering the Hessian with respect to the perturbation function h on the reference interface Γ. We shall show how to derive the expansions for small ε up to any order by the method of Taylor expansion. The main tool used here is similar to our previous work in [13] to calculate the first order derivative. After deriving the ε-expansion, we proceed to the two-parameter expansion.

5.1. Asymptotic expansions in ε

5.1.1. The extension of $u_{ε}^{\pm}$

Fig. 2.

Schematic illustration of the interface perturbations and extensions for two different cases of the interface problem. The unperturbed interfaces Γ are the vertical diameter (left) and the inner circle (right) respectively, while the perturbed interfaces $Γ_{ε}$ are the dashed lines for both cases. The unperturbed, ε-perturbed and $ε_{0}$ -perturbed subdomains are respectively $D^{+} = a \cup b \cup c$ , $D^{-} = d \cup e \cup f$ ; $D_{ε}^{+} = a \cup b \cup f$ , $D_{ε}^{-} = d \cup c \cup e$ ; $D_{ε_{0}}^{+} = a \cup e \cup f$ , $D_{ε_{0}}^{-} = d \cup b \cup c$ . $u_{ε}^{\pm}$ are extended to the sufficiently large fixed domains $D^{+} \cup D_{ε_{0}}^{+} = a \cup b \cup c \cup e \cup f$ and $D^{-} \cup D_{ε_{0}}^{-} = d \cup b \cup c \cup e \cup f$ respectively. The Cauchy problems for ${\tilde{u}}_{ε}^{\pm}$ are imposed in the thin layers $△_{ε}^{+} = c \cup e$ and $△_{ε}^{-} = b \cup f$ respectively.

The first technical issue when applying the Taylor expansion is how to extend the solutions $u_{ε}^{\pm}$ of (5.4) from their own subdomains $D_{ε}^{\pm}$ onto the larger and fixed domains which both include the interface $Γ_{ε}$ for all $ε \in [0, ε_{0}]$ . Such domains are chosen as $D^{\pm} \cup D_{ε_{0}}^{\pm}$ . On these fixed domains $D^{\pm} \cup D_{ε_{0}}^{\pm}$ , $u_{ε}^{\pm}$ are known on the parts $\overline{D_{ε}^{\pm}}$ ; we thus consider the differences $△_{ε}^{\pm}$ which consist of the disjoint thin layers: $\begin{matrix} △_{ε}^{\pm} : = (D^{\pm} \cup D_{ε_{0}}^{\pm}) ∖ \overline{D_{ε}^{\pm}} = (D^{\pm} ∖ \overline{D_{ε}^{\pm}}) \cup (D_{ε_{0}}^{\pm} ∖ \overline{D_{ε}^{\pm}}) . \end{matrix}$ Refer to Fig. 2. Denote the solution extended on $△_{ε}^{\pm}$ by ${\tilde{u}}_{ε}^{\pm}$ , and assume that ${\tilde{u}}_{ε}^{\pm}$ and $u_{ε}^{\pm}$ have the same values and the same normal derivatives on the common boundary $Γ_{ε}$ . Specifically, ${\tilde{u}}_{ε}^{\pm}$ are constructed as the unique solutions to the following Cauchy problems posed in the thin layers $D^{\pm} ∖ \overline{D_{ε}^{\pm}}$ and $D_{ε_{0}}^{\pm} ∖ \overline{D_{ε}^{\pm}}$ respectively: $\begin{matrix} \{\begin{matrix} - \nabla \cdot (σ (x) \nabla {\tilde{u}}_{ε}^{\pm} (x)) = f (x) & in △_{ε}^{\pm} = (D^{\pm} ∖ \overline{D_{ε}^{\pm}}) \cup (D_{ε_{0}}^{\pm} ∖ \overline{D_{ε}^{\pm}}), \\ {\tilde{u}}_{ε}^{\pm} = u_{ε}^{\pm}, \partial_{n} {\tilde{u}}_{ε}^{\pm} = \partial_{n} u_{ε}^{\pm} & on Γ_{ε}, \end{matrix} \end{matrix}$ where $u_{ε}^{\pm}$ , the solution to equation (5.4), are presumably given. As in Section 2.1.1, the Cauchy–Kovalevskaya theorem [11] guarantees that such extensions can be realized analytically for sufficiently small $ε_{0}$ so that the Taylor expansion can be applied in a neighbourhood of Γ.

5.1.2. Asymptotic expansions on the fixed subdomains

D^{\pm}

For ease of notation, we will still use $u_{ε}^{\pm}$ to denote their extensions defined above. Let us consider $\begin{matrix} (5.6) & u_{ε}^{\pm} = u_{0}^{\pm} + ε u_{1}^{\pm} + ε^{2} u_{2}^{\pm} + \dots . \end{matrix}$

First, for the transmission condition $u_{ε}^{+} = u_{ε}^{-}$ on $Γ_{ε}$ in (5.3), by (5.2) we have $\begin{matrix} u_{ε}^{+} (x + ε h (x) n (x)) = u_{ε}^{-} (x + ε h (x) n (x)), x \in Γ, \end{matrix}$ and by (5.6), we have $\begin{matrix} \sum_{n = 0}^{\infty} ε^{n} u_{n}^{+} (x + ε h (x) n (x)) = \sum_{n = 0}^{\infty} ε^{n} u_{n}^{-} (x + ε h (x) n (x)), x \in Γ . \end{matrix}$ Then the Taylor expansions in ε on both sides as before can yield $\begin{matrix} \sum_{k = 0}^{n} \frac{{(h (x))}^{k}}{k!} \partial_{n}^{k} u_{n - k}^{+} (x) = \sum_{k = 0}^{n} \frac{{(h (x))}^{k}}{k!} \partial_{n}^{k} u_{n - k}^{-} (x), x \in Γ, \end{matrix}$ thus $\begin{matrix} (5.7) & [u_{n} (x)] = - \sum_{k = 1}^{n} \frac{{(h (x))}^{k}}{k!} [\partial_{n}^{k} u_{n - k} (x)] on Γ, \end{matrix}$ where $[v] : = v^{+} - v^{-}$ denotes the jump across the subdomains from $D^{-}$ to $D^{+}$ .

On the other hand, from the variational form (5.4), we obtain $\begin{matrix} (5.8) & \begin{matrix} \int_{D^{+}} σ^{+} \nabla u_{ε}^{+} \cdot \nabla v d x + \int_{D^{-}} σ^{-} \nabla u_{ε}^{-} \cdot \nabla v d x \\ - \int_{δ D_{ε}} σ^{+} \nabla u_{ε}^{+} \cdot \nabla v d x + \int_{δ D_{ε}} σ^{-} \nabla u_{ε}^{-} \cdot \nabla v d x = \int_{D} f v d x, \end{matrix} \end{matrix}$ where $δ D_{ε} = (D_{ε}^{+} ∖ D^{+}) \cup (D^{+} ∖ D_{ε}^{+}) = (D^{-} ∖ D_{ε}^{-}) \cup (D_{ε}^{-} ∖ D^{-})$ , and the integrand on $δ D_{ε}$ is taken with a minus sign over $D_{ε}^{+} ∖ D^{+} = D^{-} ∖ D_{ε}^{-}$ and a plus sign over $D^{+} ∖ D_{ε}^{+} = D_{ε}^{-} ∖ D^{-}$ . Note that $\begin{array}{l} D^{+} ∖ D_{ε}^{+} = D_{ε}^{-} ∖ D^{-} = {x + t n (x) : x \in Γ, h (x) > 0, 0 ⩽ t < ε h (x)}, \\ D_{ε}^{+} ∖ D^{+} = D^{-} ∖ D_{ε}^{-} = {x + t n (x) : x \in Γ, h (x) < 0, ε h (x) < t ⩽ 0} . \end{array}$

To handle the integration $\int_{δ D_{ε}}$ in (5.8), we introduce the curvilinear coordinates $(ξ, t)$ in a sufficiently small tubular neighborhood of Γ, which are defined by $\begin{matrix} x = θ (ξ) + t n (ξ), \end{matrix}$ where $ξ = (ξ_{1}, \dots, ξ_{d - 1}) \in Ω ⟼ θ (ξ) \in Γ$ is a parametrization of the interface Γ and $n (ξ) : = n (θ (ξ))$ . Then for any smooth function g, by making use of a change of variables, we have $\begin{matrix} (5.9) & \begin{matrix} \int_{δ D_{ε}} g (x) d x = & \int_{Ω} (\int_{0}^{ε h (θ (ξ))} g (θ (ξ) + t n (ξ)) | J (ξ, t) | d t) d ξ, \end{matrix} \end{matrix}$ where $J (ξ, t)$ is the Jacobian determinant of the mapping $(ξ, t) ⟼ x$ . By Appendix B.3, (5.9) is equivalent to $\begin{matrix} (5.10) & \int_{δ D_{ε}} g (x) d x = \int_{Γ} (\int_{0}^{ε h (θ)} \tilde{g} (θ, t) det (I + t W) d t) d S_{Γ} (θ), \end{matrix}$ where $\tilde{g} (θ, t) : = g (θ + t n (θ))$ , $θ \in Γ$ , I denotes the identity matrix, $W = (W_{i}^{j})$ is the matrix representation of the Weingarten map, and $d S_{Γ} (θ)$ denotes the surface area element on the hypersurface Γ.

The equality (5.10) is the major foundation to apply the asymptotic expansion. We show how to proceed this task by considering the first two orders $u_{0}^{\pm}$ and $u_{1}^{\pm}$ . Since $\begin{matrix} (5.11) & \tilde{g} (θ, t) det (I + t W) = \tilde{g} (θ, 0) + O (t) = g (θ) + O (t), \end{matrix}$ then $\begin{matrix} (5.12) & \begin{matrix} \int_{δ D_{ε}} g (x) d x = ε \int_{Γ} h (θ) g (θ) d S_{Γ} (θ) + O (ε^{2}) . \end{matrix} \end{matrix}$ Note that on Γ, we have the following orthogonal decomposition of the gradient operator: $\begin{matrix} \nabla = \nabla_{Γ} + n \partial_{n}, \end{matrix}$ where $\nabla_{Γ}$ denotes the surface gradient operator. Then $\begin{matrix} (5.13) & [σ \nabla u_{k}] \cdot \nabla v = [σ \nabla_{Γ} u_{k}] \cdot \nabla_{Γ} v + [σ \partial_{n} u_{k}] \partial_{n} v on Γ . \end{matrix}$ Substituting (5.6) into (5.8), and applying (5.12) and (5.13), we are led to $\begin{matrix} \sum_{n = 0}^{1} ε^{n} \int_{D^{+}} σ^{+} \nabla u_{n}^{+} \cdot \nabla v d x + \sum_{n = 0}^{1} ε^{n} \int_{D^{-}} σ^{-} \nabla u_{n}^{-} \cdot \nabla v d x \\ - ε \int_{Γ} h ([σ \nabla_{Γ} u_{0}] \cdot \nabla_{Γ} v + [σ \partial_{n} u_{0}] \partial_{n} v) d S_{Γ} + O (ε^{2}) = \int_{D} f v d x . \end{matrix}$ Now we collect terms with equal powers of ε and obtain: $\begin{array}{l} (5.14) & \int_{D^{+}} σ^{+} \nabla u_{0}^{+} \cdot \nabla v d x + \int_{D^{-}} σ^{-} \nabla u_{0}^{-} \cdot \nabla v d x = \int_{D} f v d x; \\ \int_{D^{+}} σ^{+} \nabla u_{1}^{+} \cdot \nabla v d x + \int_{D^{-}} σ^{-} \nabla u_{1}^{-} \cdot \nabla v d x \\ (5.15) & - \int_{Γ} h ([σ \nabla_{Γ} u_{0}] \cdot \nabla_{Γ} v + [σ \partial_{n} u_{0}] \partial_{n} v) d S_{Γ} = 0 \forall v \in H_{0}^{1} (D) . \end{array}$ These weak formulations together with (5.7) with $n = 0, 1$ , lead to the following two PDEs for $u_{0}$ and $u_{1}$ , respectively: $\begin{matrix} (5.16) & \{\begin{matrix} - \nabla \cdot (σ (x) \nabla u_{0} (x)) = f (x) & in D^{+} \cup D^{-}, \\ [u_{0} (x)] = [σ (x) \partial_{n} u_{0} (x)] = 0 & on Γ, \\ u = g & on \partial D . \end{matrix} \end{matrix}$ and $\begin{matrix} (5.17) & \{\begin{matrix} - \nabla \cdot (σ (x) \nabla u_{1} (x)) = 0 & in D^{+} \cup D^{-}, \\ [u_{1}] = - h [\partial_{n} u_{0}], [σ \partial_{n} u_{1}] = \nabla_{Γ} \cdot (h [σ] \nabla_{Γ} u_{0}) & on Γ, \\ u_{1} = 0 & on \partial D . \end{matrix} \end{matrix}$

The equations for higher order terms, $u_{n}, n ⩾ 2$ , can be derived in the same way by considering the higher order Taylor approximations for (5.11).

5.2. Two-parameter expansions

The expansion of high-contrast ratio without interface perturbation is derived in [3]. In the sequel, we show the two-parameter expansion results by combining our ε expansion and the σ-expansion in [3]. We need to consider the following two different cases:

$ε \to 0, σ \to \infty$ ,

$ε \to 0, σ \to 0$ .

Note that σ is defined on the subdomain

D_{ε}^{-}

by (5.5). For Case (i), the solution is still bounded; for Case (ii), the solution on

D_{ε}^{-}

behaves at the order

O (1 / σ)

. The difference between these two cases is mainly a scaling factor

1 / σ

. We focus on Case (i) here. The derivation for Case (ii) can be found in Appendix C.

In Case (i), we have $ε \to 0$ and $σ \to + \infty$ . We introduce $μ = 1 / σ$ and treat ε and μ as independent small parameters. Assume $u_{ε}^{+}$ and $u_{ε}^{-}$ have double asymptotic expansions $\begin{matrix} u_{ε}^{\pm} (x) = \sum_{m, n = 0}^{\infty} u_{m, n}^{\pm} (x) ε^{m} μ^{n}, x \in D^{\pm}, \end{matrix}$ and introduce the notation $\begin{matrix} (5.18) & u_{m, \cdot}^{\pm} (x) = \sum_{n = 0}^{\infty} u_{m, n}^{\pm} (x) μ^{n}, x \in D^{\pm} . \end{matrix}$

From the condition (5.7) we obtain $\begin{matrix} (5.19) & [u_{0, \cdot}] = 0, [u_{1, \cdot}] = - h [\partial_{n} u_{0, \cdot}] on Γ, \end{matrix}$

Substituting (5.18) into (5.19) and matching the terms with the same order of μ yield that $\begin{array}{l} (5.20) & u_{0, n}^{-} = u_{0, n}^{+}, n ⩾ 0, \\ (5.21) & u_{1, n}^{-} = u_{1, n}^{+} + h [\partial_{n} u_{0, n}] n ⩾ 0 . \end{array}$ For the piecewisely homogeneous case of σ considered here, (5.14) and (5.15) become that for all $v \in H_{0}^{1} (D)$ , $\begin{array}{l} (5.22) & \int_{D^{+}} \nabla u_{0, \cdot}^{+} \cdot \nabla v d x + \int_{D^{-}} σ \nabla u_{0, \cdot}^{-} \cdot \nabla v d x = \int_{D} f v d x; \\ \int_{D^{+}} \nabla u_{1, \cdot}^{+} \cdot \nabla v d x + \int_{D^{-}} σ \nabla u_{1, \cdot}^{-} \cdot \nabla v d x \\ (5.23) & - \int_{Γ} h [(\nabla_{Γ} u_{0, \cdot}^{+} - σ \nabla_{Γ} u_{0, \cdot}^{-}) \cdot \nabla_{Γ} v + (\partial_{n} u_{0, \cdot}^{+} - σ \partial_{n} u_{0, \cdot}^{-}) \partial_{n} v] d S_{Γ} = 0 . \end{array}$ Substituting (5.18) into (5.22) and (5.23), we have that for all $v \in H_{0}^{1} (D)$ , $\begin{array}{l} (5.24) & \int_{D^{-}} \nabla u_{0, 0}^{-} \cdot \nabla v d x = 0, \\ (5.25) & \int_{D^{-}} \nabla u_{1, 0}^{-} \cdot \nabla v d x = \int_{Γ} h [\nabla_{Γ} u_{0, 0}^{-} \cdot \nabla_{Γ} v + \partial_{n} u_{0, 0}^{-} \partial_{n} v] d S_{Γ}, \end{array}$ and for $n ⩾ 0$ , $\begin{array}{l} (5.26) & \int_{D^{+}} \nabla u_{0, n}^{+} \cdot \nabla v d x + \int_{D^{-}} \nabla u_{0, n + 1}^{-} \cdot \nabla v d x = δ_{0, n} \int_{D} f v d x; \\ \int_{D^{+}} \nabla u_{1, n}^{+} \cdot \nabla v d x + \int_{D^{-}} \nabla u_{1, n + 1}^{-} \cdot \nabla v d x \\ (5.27) & - \int_{Γ} h [(\nabla_{Γ} u_{0, n}^{+} - \nabla_{Γ} u_{0, n + 1}^{-}) \cdot \nabla_{Γ} v + (\partial_{n} u_{0, n}^{+} - \partial_{n} u_{0, n + 1}^{-}) \partial_{n} v] d S_{Γ} = 0 . \end{array}$

From (5.20) and the weak formulation (5.24), (5.26), we have the following PDEs for each term: $\begin{matrix} (5.28) & \{\begin{matrix} - Δ u_{0, 0}^{-} = 0 & in D^{-}, \\ \partial_{n} u_{0, 0}^{-} = 0 & on Γ, \\ u_{0, 0}^{-} = g & on \partial D^{-} \cap \partial D; \end{matrix} \end{matrix}$ and for $n ⩾ 0$ , $\begin{matrix} (5.29) & \{\begin{matrix} - Δ u_{0, n}^{+} = δ_{0, n} f & in D^{+}, \\ u_{0, n}^{+} = u_{0, n}^{-} & on Γ, \\ u_{0, n}^{+} = δ_{0, n} g & on \partial D^{+} \cap \partial D; \end{matrix} \end{matrix}$ and for $n ⩾ 1$ , $\begin{matrix} (5.30) & \{\begin{matrix} - Δ u_{0, n}^{-} = δ_{1, n} f & in D^{-}, \\ \partial_{n} u_{0, n}^{-} = \partial_{n} u_{0, n - 1}^{+} & on Γ, \\ u_{0, n}^{-} = 0 & on \partial D^{-} \cap \partial D . \end{matrix} \end{matrix}$

We also list the PDEs for the terms with $m = 1$ : $\begin{matrix} \{\begin{matrix} - Δ u_{1, 0}^{-} = 0 & in D^{-}, \\ \partial_{n} u_{1, 0}^{-} = - \nabla_{Γ} \cdot (h \nabla_{Γ} u_{0, 0}^{-}) & on Γ, \\ u_{1, 0}^{-} = 0 & on \partial D^{-} \cap \partial D; \end{matrix} \end{matrix}$ and for $n ⩾ 0$ , $\begin{matrix} \{\begin{matrix} - Δ u_{1, n}^{+} = 0 & in D^{+}, \\ u_{1, n}^{+} = u_{1, n}^{-} - h (\partial_{n} u_{0, n}^{+} - \partial_{n} u_{0, n}^{-}) & on Γ, \\ u_{1, n}^{+} = 0 & on \partial D^{+} \cap \partial D; \end{matrix} \end{matrix}$ and $\begin{matrix} \{\begin{matrix} - Δ u_{1, n}^{-} = 0 & in D^{-}, \\ \partial_{n} u_{1, n}^{-} = \partial_{n} u_{1, n - 1}^{+} + \nabla_{Γ} \cdot [h (\nabla_{Γ} u_{0, n - 1}^{+} - \nabla_{Γ} u_{0, n}^{-})] & on Γ, \\ u_{1, n}^{-} = 0 & on \partial D^{-} \cap \partial D . \end{matrix} \end{matrix}$

Here we need to pay attention to a special situation that $\partial D^{-} = Γ$ , or equivalently, $\partial D^{-} \cap \partial D = \emptyset$ . Refer to the right panel in Fig. 2. The boundary value problems above then may become Neumann problems, which are uniquely solvable only up to an arbitrary constant. To determine those constants, as we have done in Section 4.2.2, we need the solvability condition from the next order. For $m = 0$ , and $n ⩾ 1$ , the solvability condition for (5.30) reads $\begin{matrix} (5.31) & \int_{Γ} \partial_{n} u_{0, n - 1}^{+} d S_{Γ} = - δ_{1, n} \int_{D^{-}} f d x . \end{matrix}$ (5.28) shows that $u_{0, 0}^{-} \equiv C_{0}$ . To determine $C_{0}$ , we need to look at $u_{0, 0}^{+}$ , which satisfies $\begin{matrix} (5.32) & \{\begin{matrix} - Δ u_{0, 0}^{+} = f & in D^{+}, \\ u_{0, 0}^{+} = C_{0} & on Γ, \\ u_{0, 0}^{+} = g & on \partial D, \end{matrix} \end{matrix}$ by (5.29). By the solvability condition (5.31), we have $\int_{Γ} \partial_{n} u_{0, 0}^{+} d S_{Γ} = - \int_{D^{-}} f d x$ , which uniquely determines the constant $\begin{matrix} C_{0} = - \frac{\int_{D^{-}} f + \int_{Γ} (\partial_{n} ϕ_{2} + \partial_{n} ϕ_{3})}{\int_{Γ} \partial_{n} ϕ_{1}}, \end{matrix}$ where $ϕ_{i}$ solve the following equations, respectively, $\begin{matrix} (5.33) & \{\begin{matrix} - Δ ϕ_{1} = 0 & in D^{+}, \\ ϕ_{1} = 1 & on Γ, \\ ϕ_{1} = 0 & on \partial D, \end{matrix} \{\begin{matrix} - Δ ϕ_{2} = f & in D^{+}, \\ ϕ_{2} = 0 & on Γ, \\ ϕ_{2} = 0 & on \partial D, \end{matrix} \{\begin{matrix} - Δ ϕ_{3} = 0 & in D^{+}, \\ ϕ_{3} = 0 & on Γ, \\ ϕ_{3} = g & on \partial D . \end{matrix} \end{matrix}$

Footnotes

Examples and generalizations for Section 2

Collection of proofs

Explicit formula of boundary conditions for low order terms of two-parameter expansion in Section 4

${u_{int, m, n}}$ satisfy the equations $\begin{matrix} - Δ u_{int, m, n} = δ_{0, m} δ_{0, n} f_{int} in D . \end{matrix}$ Their boundary conditions for a few lower order are listed below.

Case (ii) in Section 5

For Case (ii), we treat ε and σ as independent small parameters. Therefore we assume $u_{ε}^{+}$ and $u_{ε}^{-}$ have double asymptotic expansions $\begin{array}{l} u_{ε}^{+} (x) = \sum_{m, n = 0}^{\infty} u_{m, n}^{+} (x) ε^{m} σ^{n}, x \in D^{+}, \\ u_{ε}^{-} (x) = \sum_{m = 0}^{\infty} \sum_{n = - 1}^{\infty} u_{m, n}^{-} (x) ε^{m} σ^{n}, x \in D^{-} . \end{array}$ Note that for $u_{ε}^{-} (x)$ , the terms $σ^{n}$ start from $n = - 1$ . Define $\begin{array}{l} (D.1) & u_{m, \cdot}^{+} (x) = \sum_{n = 0}^{\infty} u_{m, n}^{+} (x) σ^{n}, x \in D^{+}, \\ (D.2) & u_{m, \cdot}^{-} (x) = \sum_{n = - 1}^{\infty} u_{m, n}^{-} (x) σ^{n}, x \in D^{-} . \end{array}$

Inserting (D.1), (D.2) into the last three equations and collecting terms with equal powers of σ, we obtain on Γ $\begin{array}{l} (D.3) & u_{0, - 1}^{-} = 0, \\ (D.4) & u_{0, n}^{-} = u_{0, n}^{+}, n ⩾ 0, \\ (D.5) & u_{1, - 1}^{-} = - h \partial_{n} u_{0, - 1}^{-}, \\ (D.6) & u_{1, n}^{-} = u_{1, n}^{+} + h [\partial_{n} u_{0, n}] n ⩾ 0; \end{array}$ and for $n ⩾ 0$ , and all $v \in H_{0}^{1} (D)$ $\begin{array}{l} (D.7) & \int_{D^{+}} \nabla u_{0, n}^{+} \cdot \nabla v d x + \int_{D^{-}} \nabla u_{0, n - 1}^{-} \cdot \nabla v d x = δ_{0, n} \int_{D} f v d x; \\ \int_{D^{+}} \nabla u_{1, n}^{+} \cdot \nabla v d x + \int_{D^{-}} \nabla u_{1, n - 1}^{-} \cdot \nabla v d x \\ (D.8) & - \int_{Γ} h [(\nabla_{Γ} u_{0, n}^{+} - \nabla_{Γ} u_{0, n - 1}^{-}) \cdot \nabla_{Γ} v + (\partial_{n} u_{0, n}^{+} - \partial_{n} u_{0, n - 1}^{-}) \partial_{n} v] d S_{Γ} = 0 . \end{array}$ Therefore we deduce from (D.3), (5.20) and (D.7) that: $u_{0, - 1}^{-}$ satisfies $\begin{matrix} \{\begin{matrix} - Δ u_{0, - 1}^{-} = f & in D^{-}, \\ u_{0, - 1}^{-} = 0 & on Γ, \\ u_{0, - 1}^{-} = 0 & on \partial D^{-} \cap \partial D . \end{matrix} \end{matrix}$ For $n ⩾ 0$ , $u_{0, n}^{+}$ satisfies $\begin{matrix} \{\begin{matrix} - Δ u_{0, n}^{+} = δ_{0, n} f & in D^{+}, \\ \partial_{n} u_{0, n}^{+} = \partial_{n} u_{0, n - 1}^{-} & on Γ, \\ u_{0, n}^{+} = δ_{0, n} g & on \partial D^{+} \cap \partial D; \end{matrix} \end{matrix}$ and $u_{0, n}^{-}$ satisfies $\begin{matrix} \{\begin{matrix} - Δ u_{0, n}^{-} = 0 & in D^{-}, \\ u_{0, n}^{-} = u_{0, n}^{+} & on Γ, \\ u_{0, n}^{-} = δ_{0, n} g & on \partial D^{-} \cap \partial D . \end{matrix} \end{matrix}$ Similarly, from (D.5), (D.6) and (D.8), we obtain that $u_{1, - 1}^{-}$ satisfies $\begin{matrix} \{\begin{matrix} - Δ u_{1, - 1}^{-} = 0 & in D^{-}, \\ u_{1, - 1}^{-} = - h \partial_{n} u_{0, - 1}^{-} & on Γ, \\ u_{1, - 1}^{-} = 0 & on \partial D^{-} \cap \partial D . \end{matrix} \end{matrix}$ For $n ⩾ 0$ , $u_{1, n}^{+}$ satisfies $\begin{matrix} \{\begin{matrix} - Δ u_{1, n}^{+} = 0 & in D^{+}, \\ \partial_{n} u_{1, n}^{+} = \partial_{n} u_{1, n - 1}^{-} - \nabla_{Γ} \cdot [h (\nabla_{Γ} u_{0, n}^{+} - \nabla_{Γ} u_{0, n - 1}^{-})] & on Γ, \\ u_{1, n}^{+} = 0 & on \partial D^{+} \cap \partial D; \end{matrix} \end{matrix}$ and $u_{1, n}^{-}$ satisfies $\begin{matrix} \{\begin{matrix} - Δ u_{1, n}^{-} = 0 & in D^{-}, \\ u_{1, n}^{-} = u_{1, n}^{+} + h (\partial_{n} u_{0, n}^{+} - \partial_{n} u_{0, n}^{-}) & on Γ, \\ u_{1, n}^{-} = 0 & on \partial D^{-} \cap \partial D . \end{matrix} \end{matrix}$

In this Case (ii), there is no emergence of the pure Neumann boundary value problem, even for the situation in the right panel in Fig. 2.

Acknowledgements

J. Chen acknowledges support from National Natural Science Foundation of China grant 21602149, the Young Thousand Talents Program of China, and the Innovation and entrepreneurial talent program in Jiangsu. L. Lin and X. Zhou acknowledge the financial support of Hong Kong GRF (109113, 11304314, 11304715). Z. Zhang acknowledges the financial support of Hong Kong RGC grants (27300616, 17300817) and National Natural Science Foundation of China via grant 11601457. J. Chen would like to thank the hospitality of Department of Mathematics, City University of Hong Kong where part of the work was done.

References

Bendali and

Lemrabet, The effect of a thin coating on the scattering of a time-harmonic wave for the Helmholtz equation, SIAM J. Appl. Math. 56 (1996), 1664–1693. doi:10.1137/S0036139995281822.

Bendali and

J.-R.

Poirier, Scattering by a highly oscillating surface, Math Methods Appl. Sci. 38 (2015), 2785–2802. doi:10.1002/mma.3261.

V.M.

Calo,

Efendiev and

Galvis, Asymptotic expansions for high-contrast elliptic equations, Math. Models Methods Appl. Sci. 24 (2014), 465–494. doi:10.1142/S0218202513500565.

Caloz,

Costabel,

Dauge and

Vial, Asymptotic expansion of the solution of an interface problem in a polygonal domain with thin layer, Asymptot. Anal. 50 (2006), 121–173.

J.E.

Castrillón-Candás,

Nobile and

R.F.

Tempone, Analytic regularity and collocation approximation for elliptic PDEs with random domain deformations, Comput. Math. Appl. 71 (2016), 1173–1197. doi:10.1016/j.camwa.2016.01.005.

Chen,

J.D.A.

Lin and

T.-Q.

Nguyen, Towards a unified macroscopic description of exciton diffusion in organic semiconductors, Commun. Comput. Phys. 20 (2016), 754–772. doi:10.4208/cicp.050615.010216a.

Chen,

Lin,

Zhang and

Zhou, Estimation of exciton diffusion lengths of organic semiconductors in random domains, J. Comp. Phys. 376 (2019), 894–912. doi:10.1016/j.jcp.2018.10.008.

Dambrine,

Greff,

Harbrecht and

Puig, Numerical solution of the Poisson equation on domains with a thin layer of random thickness, SIAM J. Numer. Anal. 54 (2016), 921–941. doi:10.1137/140998652.

M.P.

do Carmo, Differential Geometry of Curves and Surfaces, Prentice-Hall, 1976.

10.

Engquist and

J.C.

Nédélec, Effective boundary conditions for acoustic and electromagnetic scattering in thin layers, rapport interne du C.pM.A.P, Ecole Polytechnique 278 (1993).

11.

L.C.

Evans, Partial Differential Equations, American Mathematical Society, 1998.

12.

Gilbarg and

N.S.

Trudinger, Elliptic Partial Differential Equations of Second Order, Springer, 2001.

13.

Gu,

Lin and

Zhou, Sensitivity analysis and optimization of reaction rate, Commun. Math. Sci. 15 (2017), 1507–1525. doi:10.4310/CMS.2017.v15.n6.a2.

14.

Guide,

J.D.A.

Lin,

C.M.

Proctor,

Chen,

Garcia-Cervera and

T.-Q.

Nguyen, Effect of copper metalation of tetrabenzoporphyrin donor material on organic solar cell performance, J. Mater. Chem. A 2 (2014), 7890–7896. doi:10.1039/c4ta01068j.

15.

Harbrecht and

Li, First order second moment analysis for stochastic interface problems based on low-rank approximation, ESAIM: M2AN 47 (2013), 1533–1552. doi:10.1051/m2an/2013079.

16.

Harbrecht,

Peters and

Siebenmorgen, Analysis of the domain mapping method for elliptic diffusion problems on random domains, Numer. Math. 134 (2016), 823–856. doi:10.1007/s00211-016-0791-4.

17.

Harbrecht,

Schneider and

Schwab, Sparse second moment analysis for elliptic problems in stochastic domains, Numer. Math. 109 (2008), 385–414. doi:10.1007/s00211-008-0147-9.

18.

J.D.A.

Lin,

O.V.

Mikhnenko,

Chen,

Masri,

Ruseckas,

Mikhailovsky,

R.P.

Raab,

Liu,

P.W.M.

Blom,

M.A.

Loi,

C.J.

Garcia-Cervera,

I.D.W.

Samuel and

T.-Q.

Nguyen, Systematic study of exciton diffusion length in organic semiconductors by six experimental methods, Mater. Horiz. 1 (2014), 280–285. doi:10.1039/C3MH00089C.

19.

A.L.

Madureira, Modeling PDEs in domains with rough boundaries, in: Numerical Methods and Analysis of Multiscale Problems,

A.L.

Madureira, ed., Springer International Publishing, Cham, 2017, pp. 67–84. doi:10.1007/978-3-319-50866-5_5.

20.

Marschall, The trace of Sobolev–Slobodeckij spaces on Lipschitz domains, Manuscripta Math. 58 (1987), 47–65. doi:10.1007/BF01169082.

21.

P.-B.

Ming and

Xu, A multiscale finite element method for oscillating Neumann problem on rough domain, Multiscale Model. Simul. 14 (2016), 1276–1300. doi:10.1137/15M1044709.

22.

Tan and

C.-W.

Shu, Inverse Lax–Wendroff procedure for numerical boundary conditions of conservation laws, J. Comput. Phys. 229 (2010), 8144–8166. doi:10.1016/j.jcp.2010.07.014.

23.

Vial, Efficiency of approximate boundary conditions for corner domains coated with thin layers, C. R. Acad. Sci. Paris, Ser. I 340 (2005), 215–220. doi:10.1016/j.crma.2004.12.011.

24.

Xiu and

D.M.

Tartakovsky, Numerical methods for differential equations in random domains, SIAM J. Sci. Comput. 28 (2006), 1167–1185. doi:10.1137/040613160.

Asymptotic analysis for elliptic equations with small perturbations on domains in high-contrast medium

Abstract

Keywords

1. Introduction

2. The elliptic problem with smooth coefficients

2.1.1. The extension of the solution to the fixed domain

3.1. Asymptotic expansions in D and L ε

4. The thin layer problem with high-contrast ratio

4.2. Case (ii): σ / ε → 0 , ε → 0

4.2.1. Case (ii)1: Γ ≠ ∂ D , or ∂ D ∩ ∂ D ε ≠ ∅

4.2.2. Case (ii)2: Γ = ∂ D , or ∂ D ∩ ∂ D ε = ∅

4.3. Case (iii): ε / σ → c ∈ ( 0 , ∞ ) , ε → 0

5. Asymptotic expansion for the perturbed interface problem

5.1. Asymptotic expansions in ε

5.1.1. The extension of u ε ±

5.2. Two-parameter expansions

Footnotes

Examples and generalizations for Section 2

Collection of proofs

Explicit formula of boundary conditions for low order terms of two-parameter expansion in Section 4

Case (ii) in Section 5

Acknowledgements

References

3.1. Asymptotic expansions in D and $L_{ε}$

4.2. Case (ii): $σ / ε \to 0$ , $ε \to 0$

4.2.1. Case (ii)₁: $Γ \neq \partial D$ , or $\partial D \cap \partial D_{ε} \neq \emptyset$

4.2.2. Case (ii)₂: $Γ = \partial D$ , or $\partial D \cap \partial D_{ε} = \emptyset$

4.3. Case (iii): $ε / σ \to c \in (0, \infty)$ , $ε \to 0$

5.1.1. The extension of $u_{ε}^{\pm}$