Detecting Inverse Boundaries by Weighted High-Order Gradient Collocation Method

Abstract

The weighted reproducing kernel collocation method exhibits high accuracy and efficiency in solving inverse problems as compared with traditional mesh-based methods. Nevertheless, it is known that computing higher order reproducing kernel (RK) shape functions is generally an expensive process. Computational cost may dramatically increase, especially when dealing with strong-from equations where high-order derivative operators are required as compared to weak-form approaches for obtaining results with promising levels of accuracy. Under the framework of gradient approximation, the derivatives of reproducing kernel shape functions can be constructed synchronically, thereby alleviating the complexity in computation. In view of this, the present work first introduces the weighted high-order gradient reproducing kernel collocation method in the inverse analysis. The convergence of the method is examined through the weights imposed on the boundary conditions. Then, several configurations of multiply connected domains are provided to numerically investigate the stability and efficiency of the method. To reach the desired accuracy in detecting the outer boundary for two special cases, special treatments including allocation of points and use of ghost points are adopted as the solution strategy. From four benchmark examples, the efficacy of the method in detecting the unknown boundary is demonstrated.

1. Introduction

The origin of the reproducing kernel (RK) approximation can be dated back to the early development of meshfree methods in the 1970s. Among these methods, smooth particle hydrodynamics (SPH) [1] was the beginning, followed by diffuse element methods (DEM) [2] and lastly by element-free Galerkin methods (EFG) [3] and reproducing kernel particle methods (RKPM) [4,5]. The aforementioned weak-form-based meshfree methods have been successfully applied to solve partial differential equations (PDEs). Nevertheless, integration in the domain and indirect imposition of boundary conditions make these methods cumbersome. Along the direction of solving PDEs, the strong-form collocation methods were proposed almost at the end of the twentieth century [6]. In particular, the reproducing kernel collocation method (RKCM) has become one promising method with the following merits: localized feature in the approximation and well-conditioned systems in solving PDEs in comparison with global functions adopted in the approximation. The applications of RKCM include solving direct problems such as Poisson problems and elasticity problems [7,8], where a least-squares formulation is introduced to derive the weighted version of RKCM by minimizing errors on the boundary.

As it is known that the inverse problems are generally ill-posed and unstable [9], how to effectively ensure the accuracy of solutions while overcoming the disturbance arising in the measurement is of great importance in computational mechanics [10,11]. In view of the convergence and stability of the weighted RKCM in solving direct problems, Yang and her co-workers first introduced the weighted RKCM to the inverse analyses including inverse Cauchy problems [12], inverse problems with singularity [13], and inverse elasticity problems [14]. For inverse problems of Cauchy type and elasticity, the weights to be imposed on Dirichlet and Neumann boundary conditions were estimated from the error analysis by using the least-squares functional and equivalent quadrature rule via minimization; furthermore, it was shown that the weights on Dirichlet boundary conditions are different for these two types of problems, and such a difference of weights is consistent with those estimated in the direct problems of Cauchy type and elasticity. Nevertheless, computing the second derivatives of RK shape function in solving second-order PDEs is computationally expensive due to sub-functions composed in the RK shape function. To this end, the gradient reproducing kernel approximation was proposed to avoid calculating direct derivatives twice in solving second-order PDEs [15]. Such an efficient gradient RKCM was recently introduced to investigate multiply connected inverse Cauchy problems [16]; it was shown that the accuracy of approximation can be reached if proper weights are imposed on the boundary conditions, which makes similar comments on the boundary weights as given in the analysis of direct problems. Still, the gradient approximation of RK shape functions requires computing the derivatives once in solving the second-order PDEs. As pointed out in the literature [17], by using the gradient reproducing conditions, the higher-order gradient RK shape functions can be constructed implicitly without taking direct derivatives of RK shape functions, which benefits the second- and fourth-order PDEs. To overcome discretization restriction, the harmonic-enriched reproducing kernel approximation is proposed to solve highly oscillatory PDEs, in which the idea of implicit derivatives of RK shape functions is also adopted [18]. The present work further introduces the higher-order gradient reproducing kernel collocation method (HG-RKCM) in the inverse analysis. Under the least-squares framework, the weights on the boundary conditions are investigated. Then, several complex configurations of multiply connected inverse problems are studied to demonstrate the ability of the method in detecting the unknown boundary.

The arrangement of the paper is described as follows: Section 2 presents the mathematical model of the multiply connected inverse problem. In Section 3, the RK approximation and high-order gradient RK approximation are introduced in the weighted collocation method. The numerical examples are provided in Section 4; numerical issues including convergence, stability, efficiency, and strategy for enhancing accuracy are discussed. Section 5 concludes the present work.

2. Strong Form of Inverse Cauchy Problems

A multiply connected domain is depicted in Figure 1. In the direct problems, the domain $\Omega$ is enclosed completely by the boundary $\partial \Omega$ with known boundary conditions, i.e. $\partial \Omega =\partial {\Omega }_1\cup \partial {\Omega }_2\cup \partial {\Omega }_3\cup \partial {\Omega }_4$. In contrast, the inverse problems may have an incomplete boundary. For illustration purposes, consider the configuration of the inverse problem shown in Figure 1 with the following general expression:

$$\begin{array}{r}a\left(x\right)\frac{{\partial }^{2}u}{\partial x^2}+b\left(x\right)\frac{{\partial }^{2}u}{\partial x\partial y}+c\left(x\right)\frac{{\partial }^{2}u}{\partial y^2}+d\left(x\right)\frac{\partial u}{\partial x}+e\left(x\right)\frac{\partial u}{\partial y}+f\left(x\right)u=m\left(x\right)\\ \forall \hspace{0.17em}x\in \Omega \cup \partial {\Omega }_2\cup \partial {\Omega }_4\end{array}$$

(1)

$$u=g\left(x\right)\forall \hspace{0.17em}x\in \partial {\Omega }_1$$

(2)

$$u_n=h\left(x\right)\forall \hspace{0.17em}x\in \partial {\Omega }_1$$

(3)

$$u=n\left(x\right)\forall \hspace{0.17em}x\in \partial {\Omega }_3$$

(4)

where the coefficients $a\left(x\right)$, $b\left(x\right)$, $c\left(x\right)$, $d\left(x\right)$, $e\left(x\right)$, and $f\left(x\right)$ are known; $m\left(x\right)$, $n\left(x\right)$, $g\left(x\right)$, and $h\left(x\right)$ are given functions in the domain and on the boundary. As a variety of multiply connected inverse Cauchy problems will be investigated in this work, it is assumed that both Neumann and Dirichlet boundary conditions are imposed on the boundary $\partial {\Omega }_1$ while $\partial {\Omega }_2$ and $\partial {\Omega }_4$, respectively, denote the unknown outer and inner boundary conditions, i.e., the missing boundary. In addition, $\partial {\Omega }_3$ is the known boundary prescribed by the Dirichlet boundary condition.

3. Weighted High-Order Gradient Reproducing Kernel Collocation Method

3.1. Reproducing Kernel Approximation

To approximate an unknown function $u(x)$ by using the collocation method, the reproducing kernel (RK) shape function evaluated at $N_s$ source points can be used in conjunction with the generalized coefficients $a_I$ given by

$$u\left(x\right)\approx v\left(x\right)={\displaystyle \sum _{I=1}^{N_s}{\Psi }_I\left(x\right)} a_I$$

(5)

where the RK shape function ${\Psi }_I\left(x\right)$ has the following form:

$${\Psi }_I\left(x\right)=H^{\mathrm{T}}\left(x-x_I\right)b\left(x\right){\phi }_a\left(x-x_I\right)$$

(6)

which consists of a monomial $H\left(x-x_I\right)$, coefficient vector $b\left(x\right)$, and kernel function ${\phi }_a\left(x-x_I\right)$. The influence of kernel function is related to the RK support size as denoted by the subscript a. A general form of a two-dimensional monomial basis is expressed as

$$H^{\mathrm{T}}\left(x-x_I\right)={\left\{{\left(x-x_I\right)}^{\alpha }\right\}}_{\left|\alpha \right|\le p}$$

(7)

where the multi-index $\alpha$ is defined by $\alpha =\left({\alpha }_1,{\alpha }_2\right)$ with length defined as $\left|\alpha \right|={\alpha }_1+{\alpha }_2$; the superscript $p$ denotes the order of approximation to reproduce a polynomial. As the inverse Cauchy problems considered in this work are governed by the second-order differential equations, a quintic B-spline kernel function is adopted and is given in the previous work [12]. With regard to the coefficient vector, it is determined by the following reproducing conditions:

$${\displaystyle \sum _{I=1}^{N_s}{\Psi }_I\left(x\right)x_I^{\alpha }}=x^{\alpha },\left|\alpha \right|\le p$$

(8)

from which the coefficient vector is found to be

$$b\left(x\right)=M^{-1}\left(x\right)H\left(0\right)$$

(9)

where $M\left(x\right)$ is the moment matrix given by

$$M\left(x\right)={\displaystyle \sum _{I=1}^{N_s}H\left(x-x_I\right)H^{\mathrm{T}}\left(x-x_I\right){\phi }_a\left(x-x_I\right)}$$

(10)

By substituting Equation (9) into Equation (6), the RK shape function reaches the following form:

$${\Psi }_I\left(x\right)=H^{\mathrm{T}}\left(0\right)M^{-1}\left(x\right)H\left(x-x_I\right){\phi }_a\left(x-x_I\right)$$

(11)

3.2. High-Order Gradient Reproducing Kernel Approximation

As an extension of the gradient reproducing kernel approximation [15], the second derivatives of RK shape functions are obtained by using the high-order gradient reproducing kernel (HG-RK) approximation [17]. Unlike the direct derivatives of RK shape functions, the implicit derivatives of RK shape functions in two dimensions are constructed as follows:

$$\begin{array}{l}{\Psi }_I^x\left(x\right)=-H_{,x}^{\mathrm{T}}\left(0\right)M^{-1}\left(x\right)H\left(x-x_I\right){\phi }_a\left(x-x_I\right)\\ {\Psi }_I^y\left(x\right)=-H_{,y}^{\mathrm{T}}\left(0\right)M^{-1}\left(x\right)H\left(x-x_I\right){\phi }_a\left(x-x_I\right)\end{array}$$

(12)

and

$$\begin{array}{l}{\Psi }_I^{xx}\left(x\right)=-H_{,xx}^{\mathrm{T}}\left(0\right)M^{-1}\left(x\right)H\left(x-x_I\right){\phi }_a\left(x-x_I\right)\\ {\Psi }_I^{yy}\left(x\right)=-H_{,yy}^{\mathrm{T}}\left(0\right)M^{-1}\left(x\right)H\left(x-x_I\right){\phi }_a\left(x-x_I\right)\end{array}$$

(13)

where $H_{,x}^{\mathrm{T}}$ and $H_{,y}^{\mathrm{T}}$ are the first derivatives of $H^{\mathrm{T}}\left(x-x_I\right)$; $H_{,xx}^{\mathrm{T}}$ and $H_{,yy}^{\mathrm{T}}$ are the second derivatives of $H^{\mathrm{T}}\left(x-x_I\right)$. Similarly, $H^{\mathrm{T}}\left(x-x_I\right)$ is the $q$th order monomial basis given by

$$H^T\left(x-x_I\right)={\left\{{\left(x-x_I\right)}^{\alpha }\right\}}_{\left|\alpha \right|\le q}$$

(14)

It is noted that $q\ne p$ in general. Nevertheless, $q\ge 2$ is the convergent requirement for solving second-order differential equations in the weighted collocation method.

3.3. Weighted Discrete Collocation Equations by HG-RK Approximation

According to the schematic diagram in Figure 1, the following strong form by using high-order gradient reproducing kernel approximation is considered:

$$L^1{u}_{,xx}+L^2{u}_{,yy}=m \mathrm{in} \Omega \cup \partial {\Omega }_2\cup \partial {\Omega }_4$$

(15)

$$B_{g}u=g \mathrm{on} \partial {\Omega }_1\cup \partial {\Omega }_3$$

(16)

$$B_h^1{u}_{,x}+B_h^2{u}_{,y}=h \mathrm{on} \partial {\Omega }_1$$

(17)

where $L^1$ and $L^2$ denote the differential operators in the domain domain $\Omega \cup \partial {\Omega }_2\cup \partial {\Omega }_4$; $B_h^1$ and $B_h^2$ denote the Neumann boundary operators on $\partial {\Omega }_1$; $B_g$ denotes the Dirichlet boundary operator on $\partial {\Omega }_1\cup \partial {\Omega }_3$. As the detailed derivations of weighted collocation equations for solving inverse problems are presented in several studies [12,14], this work mainly focuses on the application of weighted RKCM in the inverse analysis. A brief description is given herein.

Based on Equations (15)–(17), a weighted least-squares functional is established with the following expression:

$$\begin{array}{l}E\left(u\right)=\frac{1}{2}{\displaystyle {\int }_{\Omega \cup \partial {\Omega }_2\cup \partial {\Omega }_4}{\left(L^1{u}_{,xx}+L^2{u}_{,yy}-m\right)}^{\mathrm{T}}\left(L^1{u}_{,xx}+L^2{u}_{,yy}-m\right)d\Omega }\\ \hspace{1em}\hspace{1em}+\frac{{\alpha }_h}{2}{\displaystyle {\int }_{\partial {\Omega }_1}{\left(B_h^1{u}_{,x}+B_h^2{u}_{,y}-h\right)}^{\mathrm{T}}\left(B_h^1{u}_{,x}+B_h^2{u}_{,y}-h\right)d{\Omega }_1}\\ \hspace{1em}\hspace{1em}+\frac{{\alpha }_g}{2}{\displaystyle {\int }_{\partial {\Omega }_1\cup \partial {\Omega }_3}{\left(B_{g}u-g\right)}^{\mathrm{T}}\left(B_{g}u-g\right)d{\Omega }_1}\end{array}$$

(18)

where these three terms correspond to least-squares residuals in the domain and on the boundary. To balance the errors of approximation in the domain and on the boundary, the weights ${\alpha }_h$ and ${\alpha }_g$ are further introduced.

In the collocation method, the collocation points can be classified into different groups according to their locations in the problem domain. In this way, the integral form of the functional can be evaluated equivalently by treating the collocation points as the integration points in the discrete form. Then, introducing the RK approximation together with HG-RK approximation to this discrete form and using the stationary condition of the functional lead to the following weighted discrete collocation system:

$$\left(\begin{array}{c}{A}^1\\ \sqrt{{\alpha }_h}{A}^2\\ \sqrt{{\alpha }_g}{A}^3\end{array}\right)u=\left(\begin{array}{c}{b}^1\\ \sqrt{{\alpha }_h}{b}^2\\ \sqrt{{\alpha }_g}{b}^3\end{array}\right)$$

(19)

where the sub-matrices are given by

$$A^1=\left(\begin{array}{c}{L}^1\left({\Psi }^{x{x}^{\mathrm{T}}}\left(x_{p1}\right)\right)+L^2\left({\Psi }^{y{y}^{\mathrm{T}}}\left(x_{p1}\right)\right)\\ L^1\left({\Psi }^{x{x}^{\mathrm{T}}}\left(x_{p2}\right)\right)+L^2\left({\Psi }^{y{y}^{\mathrm{T}}}\left(x_{p2}\right)\right)\\ ⋮\\ L^1\left({\Psi }^{x{x}^{\mathrm{T}}}\left(x_{p{N}_p}\right)\right)+L^2\left({\Psi }^{y{y}^{\mathrm{T}}}\left(x_{p{N}_p}\right)\right)\end{array}\right)$$

(20)

$$A^2=\left(\begin{array}{c}{B}_h^1\left({\Psi }^{x^{\mathrm{T}}}\left(x_{q1}\right)\right)+B_h^2\left({\Psi }^{y^{\mathrm{T}}}\left(x_{q1}\right)\right)\\ B_h^1\left({\Psi }^{x^{\mathrm{T}}}\left(x_{q2}\right)\right)+B_h^2\left({\Psi }^{y^{\mathrm{T}}}\left(x_{q2}\right)\right)\\ ⋮\\ B_h^1\left({\Psi }^{x^{\mathrm{T}}}\left(x_{q{N}_q}\right)\right)+B_h^2\left({\Psi }^{y^{\mathrm{T}}}\left(x_{q{N}_q}\right)\right)\end{array}\right)$$

(21)

$$A^3=\left(\begin{array}{c}{B}_g\left({\Psi }^{\mathrm{T}}\left(x_{r1}\right)\right)\\ B_g\left({\Psi }^{\mathrm{T}}\left(x_{r2}\right)\right)\\ ⋮\\ B_g\left({\Psi }^{\mathrm{T}}\left(x_{r{N}_r}\right)\right)\end{array}\right)$$

(22)

$$b^1=\left(\begin{array}{c}m\left(x_{p1}\right)\\ m\left(x_{p2}\right)\\ ⋮\\ m\left(x_{p{N}_p}\right)\end{array}\right), b^2=\left(\begin{array}{c}h\left(x_{q1}\right)\\ h\left(x_{q2}\right)\\ ⋮\\ h\left(x_{q{N}_q}\right)\end{array}\right), b^3=\left(\begin{array}{c}g\left(x_{r1}\right)\\ g\left(x_{r2}\right)\\ ⋮\\ g\left(x_{r{N}_r}\right)\end{array}\right)$$

(23)

In Equation (19), $u$ is the vector containing generalized coefficients in both RK and HG-RK approximations, i.e., $u={\left[a_1,a_2,\dots \dots ,a_{N_s}\right]}^{\mathrm{T}}$. Besides, $x_p$ contains $N_p$ points in the domain; $x_q$ contains $N_q$ points on the Neumann boundary; $x_r$ contains $N_r$ points on the Dirichlet boundary.

4. Numerical Implementation

The convergence of the weighted HG-RKCM in solving inverse Cauchy problems is first demonstrated. Then, several complicated multiply connected inverse problems are provided to demonstrate the efficacy of the method.

4.1. Convergence Study

The convergence of the weighted collocation method is closely related to the weights imposed on the boundary conditions. According to the weights estimated in the inverse analysis of Cauchy problems, the following weights are suggested [12]:

$$\sqrt{{\alpha }_h}\approx O\left(1\right),\sqrt{{\alpha }_g}\approx O\left(\sqrt{N_s}\right)$$

(24)

Nevertheless, in the gradient approximation, the following weights are suggested [15]:

$$\sqrt{{\alpha }_h}\approx O\left(1\right),\sqrt{{\alpha }_g}\approx O\left(\frac{1}{a}\right)$$

(25)

where $p=q=2$ is assumed; $a$ is the RK support size. To determine the most suitable weights in the higher-order gradient approximation in the inverse analysis, an inhomogeneous Helmholtz problem is considered, which is given by

$${\nabla }^{2}u+2u=\cos x \mathrm{in} \Omega =\left(0,1\right)\times \left(0,1\right)$$

(26)

The boundary conditions are

$$\begin{array}{ccc}u=\left[\sin \left(y\right)+1\right]\left(1+\frac{s}{100}\times rand\right)& \mathrm{on}& x=0\\ u=\cos \left(1\right)\left[\sin \left(y\right)+1\right]\left(1+\frac{s}{100}\times rand\right)& \mathrm{on}& x=1\\ u=\cos \left(x\right)\left[\sin \left(1\right)+1\right]\left(1+\frac{s}{100}\times rand\right)& \mathrm{on}& y=1\end{array}$$

(27)

and

$$\begin{array}{ccc}{u}_{,x}=0+\frac{s}{100}\times rand& \mathrm{on}& x=0\\ u_{,x}=-\sin \left(1\right)\left[\sin \left(y\right)+1\right]\left(1+\frac{s}{100}\times rand\right)& \mathrm{on}& x=1\\ u_{,y}=\cos \left(x\right)\cos \left(1\right)\left(1+\frac{s}{100}\times rand\right)& \mathrm{on}& y=1\end{array}$$

(28)

where $s$ is the noise applied to the boundary conditions, and $rand$ is a random number in the interval $0\sim 1$. The analytical solution along $y=0$ is given by the infinite-order polynomial by setting $s=0$:

$$u\left(x,y\right)=\cos x\sin y+\cos x$$

(29)

In this square domain, four sets of uniformly discretized source points and collocation points are adopted in the convergence study: $N_s=N_c=\left\{{10}^2,{20}^2,{30}^2,{40}^2\right\}$. The support size of the RK shape function is $4.6d$ with $d$ representing the average nodal distance between two adjacent source points. The error norms are defined as

$${\Vert v-u\Vert }_0=\sqrt{{\displaystyle {\int }_{\Omega }{\left(v-u\right)}^{2}d\Omega }}$$

(30)

$${\Vert v_{\alpha }-u_{,\alpha }\Vert }_0=\sqrt{{\displaystyle {\int }_{\Omega }\left[{\left(v_x-u_{,x}\right)}^2+{\left(v_y-u_{,y}\right)}^2\right]d\Omega }}$$

(31)

in which $v$ denotes the approximation of $u$. As depicted Figure 2, a similar trend of convergence behavior is observed for $\sqrt{{\alpha }_g}=10$ and $\sqrt{{\alpha }_g}=\frac{1}{a}$ while $\sqrt{{\alpha }_g}=\frac{1}{a}$ exhibits higher convergence rates among the two weights. Since the expression of $\sqrt{{\alpha }_g}\approx O\left(\frac{1}{a}\right)$ is implicitly related to the discretization via RK support size, it is adopted in the following study. For uniform discretization with $N_s=N_c=40\times 40$, is adopted in both uniform and non-uniform discretization while 5% disturbance is added to the domain points to obtain non-uniform discretization as shown in Figure 3. For uniform discretization, the contour plot of absolute error (%) shows that large absolute error occurs along the missing boundary y = 0 with a maximal value of 0.11506% ($s=0$) as shown in Figure 4a; similarly, the non-uniform discretization has a maximal absolute error 0.14779% as shown in Figure 4b. In the uniform discretization, by adding some noise s = 1,2,3, on the boundary conditions, the corresponding maximal absolute errors for $s=0,1,2,3$ are 1.01703%, 1.41618%, and 2.61759%, respectively; referring to Figure 5a, good agreement between the analytical solution and approximations is observed except $s=3$. For non-uniform discretization, the comparison of the analytical solution and approximations are shown in Figure 5b; although the deviation from the analytical solution is more pronounced than the results shown in Figure 5a, the maximal absolute error for $s=3$ is 5.24307%, which is still small and acceptable. From the above investigation, the stability of the method is demonstrated. Finally, the CPU time recorded for RKCM, G-RKCM, and HG-RKCM is given in Figure 6, and HG-RKCM is the most efficient one among the three methods.

4.2. Annular Domains

Consider three annular domains as discretized in Figure 7a–c, where these annular domains have different boundary conditions: half unknown outer boundary, half unknown inner boundary, and unknown inner boundary. In particular, the star denotes the symbol for points in the domain; the circle denotes the symbol for points imposed by Neumann and Dirichlet boundary conditions; the triangle denotes the symbol for the inner boundary points, as prescribed by the Dirichlet boundary condition. The governing equation is described by the following Laplace equation:

$${\nabla }^{2}u=0,R_i\le r\le R_o, 0\le \theta <2\pi$$

(32)

where $R_i=0.5$ and $R_o=1$. The general expressions for the boundary conditions in these three configurations are given by

$$\begin{array}{l}u\left(R_o,\theta \right)=R_o^2\cos \left(2\theta \right)\left(1+\frac{s}{100}\times rand\right), {\beta }_1\le \theta \le {\beta }_2\\ u_n\left(R_o,\theta \right)=2R_o\cos \left(2\theta \right)\left(1+\frac{s}{100}\times rand\right), {\beta }_1\le \theta \le {\beta }_2\end{array}$$

(33)

$$u\left(R_i,\theta \right)=R_i^2\cos \left(2\theta \right), {\beta }_3\le \theta \le {\beta }_4$$

(34)

where ${\beta }_i$ denotes the angle corresponding to the prescribed boundary conditions as shown in Figure 7a–c, which satisfies $0\le {\beta }_i\le \pi$. Along the missing boundary, $u$ should satisfy the analytical solution given below:

$$u=R_o^2\cos \left(2\theta \right)$$

(35)

for the outer boundary or

$$u=R_i^2\cos \left(2\theta \right)$$

(36)

for the inner boundary.

By using $N_s=N_c=386$ and $a=3d$, the contour plots of absolute errors (%) are depicted in Figure 8a–c. Obviously, large errors occur along the missing boundary with an extremely small magnitude of absolute errors. Upon using the quadratic basis to construct the RK and HG-RK shape functions, high accuracy of approximation to the quadratic polynomial is retained. Then, the outer boundary conditions are disturbed by some noise $s=0,1,2,3$. The comparison between numerical solutions along the missing boundary and analytical solution is made in Figure 9a–c, from which great agreement is achieved even when $s=3$. The maximal absolute errors (%) for $u$, $u_{,x}$, and $u_{,y}$ obtained by various $s$ are summarized

Table 1. Comparison of maximal absolute errors (%) obtained by various s: (a) half unknown outer boundary; (b) half unknown inner boundary; (c) unknown inner boundary.

s		0	1	2	3
(a)	u	6.2683 × 10−10	1.0422	2.3082	3.7493
	u,x	6.2619 × 10−10	2.9022	4.7829	5.5098
	u,y	1.2626 × 10−9	2.7713	4.9511	6.7354
(b)	u	3.4001 × 10−12	0.1038	0.2400	0.4209
	u,x	2.8721 × 10−11	0.9011	1.8392	2.6904
	u,y	2.2637 × 10−11	0.8725	1.7428	1.7441
(c)	u	3.6526 × 10−12	0.2722	0.3525	0.5416
	u,x	1.8374 × 10−11	1.2760	1.7217	1.9987
	u,y	1.6998 × 10−11	1.1475	1.5912	2.0693

, in which high accuracy is observed. The stability of the method in detecting the unknown boundary is demonstrated. As for the CPU calculation, HG-RKCM exhibits the most efficiency, as shown in Figure 10, while G-RKCM is more efficient than RKCM.

4.3. Multiply Connected Domains: Type One

Referring to Figure 11a–d, four multiply connected domains equipped with various boundary conditions are in consideration: half unknown inner boundary, two half unknown inner boundaries, unknown inner boundary, and half unknown outer boundary. The domain is enclosed by a circular perimeter with radius $R_o=1$; each inner hole has a radius equal to 0.2. Again, the domain is governed by the Laplace equation while a third-order approximation is adopted. The general expressions for the boundary conditions on the outer boundary in the four configurations are described by

$$\begin{array}{l}u\left(R_o,\theta \right)=\left[-R_o^3\sin \left(3\theta \right)+R_o^2\cos \left(2\theta \right)\right]\left(1+\frac{s}{100}\times rand\right), {\beta }_1\le \theta \le {\beta }_2\\ u_n\left(R_o,\theta \right)=\left[-3R_o^2\sin \left(3\theta \right)+2R_o\cos \left(2\theta \right)\right]\left(1+\frac{s}{100}\times rand\right), {\beta }_1\le \theta \le {\beta }_2\end{array}$$

(37)

where ${\beta }_i$ has the same definition as stated previously. For $s=0$, $u$ on the missing boundary should satisfy the analytical solution in the inverse analysis, which is

$$u=-R^3\sin \left(3\theta \right)+R^2\cos \left(2\theta \right)$$

(38)

As depicted in Figure 11a–d, the same set of source points and collocation points are adopted in the discretization of four configurations, i.e., $N_s=N_c=470$. For Figure 11a–c, RK support size is $a=3.3d$. For Figure 11d, RK support size is $a=3.5d$; further, it is noted that the points on the circumference have been spun with a small angle to avoid overlaps of points, as can be observed in the configuration. The contour plots of absolute errors (%) are shown in Figure 12a–d; it is observed that large errors occur along the missing boundary, and the distribution of errors exhibits symmetry as can be expected due to symmetric configurations of the problems. Except the case shown in Figure 12d, the absolute errors are quite small. The use of quadratic basis in the formulation of the RK and HG-RK shape functions is able to yield the desired accuracy in the approximation of the cubic polynomial. Concerning the percentage of disturbance s added on the outer boundary conditions, Figure 13a–d compares the numerical solutions with analytical ones, and good agreement is observed even when $s=3$, except the difficult case of Figure 11d.

Table 2. Comparison of maximal absolute errors (%) obtained by various s: (a) half unknown inner boundary; (b) two half unknown inner boundaries; (c) unknown inner boundaries; (d) half unknown outer boundary.

s		0	1	2	3
(a)	u	0.3942	0.5251	0.5947	0.8998
	u,x	4.1349	4.5853	4.9640	5.2201
	u,y	3.0719	3.3959	3.9454	4.1871
(b)	u	0.3899	0.7809	0.9303	1.1415
	u,x	3.4197	5.2031	5.5692	8.8528
	u,y	3.1640	3.7078	4.0539	8.7244
(c)	u	0.3932	0.6791	1.0372	1.6399
	u,x	4.0919	4.6996	5.0791	6.5832
	u,y	3.0315	3.3222	3.9394	4.4157
(d)	u	4.6230	9.7452	11.6293	-
	u,x	29.5079	32.3523	39.8299	-
	u,y	27.6338	35.6646	41.2359	-

summarizes the maximal absolute errors (%) for $u$, $u_{,x}$, and $u_{,y}$ obtained by various $s$ added on the outer boundary conditions. Again, high accuracy in the approximation is reached, excluding the difficult case of Figure 11d. As such, the stability of the method is demonstrated. Referring to Figure 14, it is observed explicitly that HG-RKCM is the most efficient method.

4.4. Multiply Connected Domains: Type Two

As shown in Figure 15a–d, four multiply connected domains having various boundary conditions are considered: unknown upper inner boundary, three half unknown inner boundaries, unknown inner boundaries, and half unknown outer boundary. The outer radius $R_o$ of the domain is 1.0 while the inner radius in each hole is 0.2. The governing equation and boundary conditions are the same as the ones described in Section 4.3. For the cases shown in Figure 15a–c, $N_s=N_c=467$ and RK support size is $a=3d$. For the case shown in Figure 15d, a half layer of ghost points is arranged outside the upper boundary with a distance $\left(1+1.25d\right)$ from the point $\left(0,0\right)$ in the domain in order to reach desired accuracy; thus, $N_s=N_c=498$ and RK support size is $a=3.9d$. The contour plots of absolute errors (%) are shown in Figure 16a–d, in which the locations of large errors are consistent with the locations of unknown boundaries, and symmetric distribution of errors is explicitly observed. The overall trend of absolute errors is small. Then, the disturbance $s=0,1,2,3$ is added to the outer boundary conditions as given in Equation (37). Figure 17a–d compares the numerical solutions with analytical solutions, and good agreement is observed even when s = 3, except the case shown in Figure 17d. The maximal absolute errors (%) for $u$, $u_{,x}$, and $u_{,y}$ obtained by various $s$ on the unknown boundaries are summarized in

Table 3. Comparison of maximal absolute errors (%) obtained by various $s$: (a) unknown upper inner boundary; (b) three half unknown inner boundaries; (c) unknown inner boundaries; (d) half unknown outer boundary.

s		0	1	2	3
(a)	u	0.9392	0.9703	1.1766	1.2071
	u,x	4.3020	3.8733	3.8337	4.9856
	u,y	6.5409	6.4914	7.1475	7.2616
(b)	u	1.2561	1.3938	1.4540	1.6904
	u,x	5.5221	6.3356	7.0172	8.0663
	u,y	8.7736	7.5323	7.7360	7.5732
(c)	u	0.9973	1.0513	1.2115	1.4096
	u,x	4.8008	5.1810	5.6641	7.2072
	u,y	7.4692	8.0126	5.3486	8.0440
(d)	u	5.6311	6.5935	10.9980	-
	u,x	25.7791	34.7108	37.6828	-
	u,y	30.7985	34.9769	44.3546	-

. It is observed that Figure 15d is the most difficult case, where $u_{,x}$ and $u_{,y}$ are sensitive to the disturbance. Apart from this, the stability of the method is demonstrated. As shown in Figure 18, HG-RKCM is the most efficient method among the three methods.

5. Conclusions

In this work, the weighted HG-RKCM is introduced to solve inverse Cauchy problems having multiply connected domains. In particular, for two multiply connected domains of different geometry, four different configurations have been provided to investigate the ability of the method in detecting the unknown boundary. Firstly, the weights for Neumann and Dirichlet boundary conditions have been numerically studied. Based on this, the stability of the method is tested by adding disturbance on both Neumann and Dirichlet boundary conditions. From the numerical results, high accuracy in the approximation is observed. For the cases of missing boundary on the circumference, two possible treatments have been adopted: spin of points on the circumference to avoid overlapping points in the domain and arrangement of one layer of ghost points on the known boundary. Lastly, the comparison of CPU time for solving various problems has come to the conclusion that the HG-RKCM is the most efficient method in the inverse analysis.