A Numerical Method for Filtering the Noise in the Heat Conduction Problem

Abstract

In this paper, we give an effective numerical method for the heat conduction problem connected with the Laplace equation. Through the use of a single-layer potential approach to the solution, we get the boundary integral equation about the density function. In order to deal with the weakly singular kernel of the integral equation, we give the projection method to deal with this part, i.e., using the Lagrange trigonometric polynomials basis to give an approximation of the density function. Although the problems under investigation are well-posed, herein the Tikhonov regularization method is not used to regularize the aforementioned direct problem with noisy data, but to filter out the noise in the corresponding perturbed data. Finally, the effectiveness of the proposed method is demonstrated using a few examples, including a boundary condition with a jump discontinuity and a boundary condition with a corner. Whilst a comparative study with the method of fundamental solutions (MFS) is also given.

1. Introduction

It is well-known that numerous problems are described by boundary value problems governed by the partial differential equations, which arises in many areas of science, such as wave propagation, vibration, electromagnetic scattering, nondestructive testing, geophysics, and cardiology. The steady state heat conduction equation or the Laplace equation is considered to describe the temperature in heat conducting material. This phenomenon is usually considered in the engineering problem. In practice, the boundary values are the measured noisy data, thus the solution may not be unique. Thus some regularization method is needed to get a stable solution.

Many problems can be solved by the boundary element method or boundary integral equation method, such as for heat conduction problems [1,2,3,4], acoustic problem [5,6,7,8], elasticity [9], potential problem [10], etc. The direct boundary element method can be obtained from Green’s third identity or Betti’s law as well as Somigliana’s identity, while the indirect boundary element method was motivated by the superposition of singular solutions. It is interesting that the roles of source and field points are reversed in the two systems. A clear link is established in Hong and Chen’s paper [11]. Sometimes, the boundary element method will result in a degenerate scale [10,12,13]. Thus, many authors have enriched this method. In [14], Alves and Valtchev investigates the boundary conditions with jump discontinuities of the Laplace equation by enriching the MFS approximation basis with a set of particular solutions. We can refer [15] for the invariant method of fundamental solutions and [16] for the normal derivatives of the fundamental solutions. In [17], Lin, et al. give a multiple/scale/direction Trefftz method for the numerical solution of the direct problem as well as the Cauchy problem of the multi-dimensional Laplace equation in an arbitrary domain. The method can largely reduce the number of unknown coefficients appeared in the series expansion of the numerical solution. In [18], Liu, et al. introduce a novel concept to construct Trefftz energy bases based on the MFS for the numerical solution of the Cauchy problem in an arbitrary star plane domain. The Trefftz energy bases used for the solution not only satisfy the Laplace equation but also preserve the energy. In [19], Fu, et al. present a strong-form boundary collocation method, the boundary knot method, in conjunction with Laplace transform to solve the heat conduction equations. In [20], Feng, et al. present a new transformation technique from multi-domain integrals into equivalent boundary integrals, by utilizing this technique, for solving heat conduction in multi-non-homogeneous media with multiple heat sources. In [21], Fu, et al. investigate the thermal behavior inside skin tissue with the presence of a tumor by solving the 3D Pennes bioheat equation with the method of approximate particular solutions. To ensure the interface conditions were exactly satisfied, an affine space decomposition technique is adopted. Recently, Gu, et al. [22] investigate the application of the generalized finite difference method to three-dimensional transient heat conduction in anisotropic composite materials. In the algorithm, the Krylov deferred correction method is introduced to perform temporal discretization in time-domain. In [23], a modified dual-level fast multipole boundary element method is investigated for large-scale three-dimensional potential problems, and the main idea is based on a dual-level structure to handle the excessive storage requirement and ill-conditioned problems resulting from the fully-populated interpolation matrix of the boundary element method. Another effective method, finite element methods, is also widely used for the numerical solution, see [24,25]. Another popular method for the direct and inverse heat conduction problems is the singular boundary method, we can see [26,27,28,29,30] for further reading.

The method of fundamental solutions (MFS) can give a good approximation of the solution [16,31,32], but it will result in an ill-conditioned system. In this paper, a novel numerical scheme is proposed to give a stable numerical solution of the Dirichlet problem for the Laplace equation, and it can reduce the condition number of the resulting system. This will display in the last two examples. In practice, the boundary data contains the noise. From the classical theorem of the partial differential equations [33], we know that the numerical solution in the solution domain is continuously dependent on the boundary data. Thus we are only to filter out the noise in the corresponding perturbed data on the boundary, and we will get the accurate solution in the solution domain. Here we use the Tikhonov regularization method [34] to get a stable solution of the boundary data. In fact, it is enough to study the calculation results on the boundary.

This paper is organized as follows. In Section 2, we give the mathematical formulation and the solution method. In Section 3, we give a brief review of the regularization method. At last, some examples are given to show the effectiveness of the method even for a boundary condition with a jump discontinuity.

2. Mathematical Formulation and the Solution Method

Let $D\subset {\mathbb{R}}^2$ be a simply connected bounded domain with piecewise smooth boundary $\partial D$. Consider the following boundary value problem: given the temperature $f\left(\mathbf{x}\right),\mathbf{x}\in \partial D$, on the boundary $\partial D$, find the solution $u\left(\mathbf{x}\right)$ satisfying

$$\Delta u\left(\mathbf{x}\right)=0,\mathbf{x}\in D,$$

(1)

$$u\left(\mathbf{x}\right)=f\left(\mathbf{x}\right),\mathbf{x}\in \partial D.$$

(2)

Moreover, we also assume in this paper that the boundary condition associated with the heat conduction problem is the noisy data denoted by $f^{\delta }$, and this corresponds to real life problems encountered in practice.

From [35], we know that, for some density function $\phi \left(\mathbf{y}\right),\mathbf{y}\in \partial D$, the single layer potential

$$u\left(\mathbf{x}\right)={\int }_{\partial D}\mathsf{\Phi }(\mathbf{x},\mathbf{y})\phi \left(\mathbf{y}\right)\mathrm{d}S\left(\mathbf{y}\right),\mathbf{x}\in D,$$

(3)

is a solution of the Equation (1), since the fundamental solution $\mathsf{\Phi }(\mathbf{x},\mathbf{y})=-\frac{1}{2\pi }ln|\mathbf{x}-\mathbf{y}|$ satisfy the Equation (1).

If the single layer potential (3) satisfies the boundary condition (2), i.e.,

$${\int }_{\partial D}\mathsf{\Phi }(\mathbf{x},\mathbf{y})\phi \left(\mathbf{y}\right)\mathrm{d}S\left(\mathbf{y}\right)=f\left(\mathbf{x}\right),\mathbf{x}\in \partial D,$$

(4)

then the single layer potential (3) is the solution of the problem consisting Equations (1) and (2). It means that we have

$$-\frac{1}{2\pi }{\int }_{\partial D}\phi \left(\mathbf{y}\right)ln|\mathbf{x}-\mathbf{y}|\mathrm{d}S\left(\mathbf{y}\right)=f\left(\mathbf{x}\right),\mathbf{x}\in \partial D.$$

(5)

In this paper, we consider the solution domain with the boundary $\partial D$, which has a parameter expression $\gamma \left(t\right),0\le t\le 2\pi$. We denote the boundary point $\gamma \left(t\right)=(x_1\left(t\right),x_2\left(t\right))$ as the parametrization of the boundary with semiaxis 1 and 2. With the simple calculation, the above equation, combining the parameter expression of $\partial D$, changes to

$$-\frac{1}{2\pi }{\int }_0^{2\pi }\varphi \left(s\right)ln|\gamma \left(t\right)-\gamma \left(s\right)|\mathrm{d}s=f\left(\gamma \left(t\right)\right),$$

(6)

where $\varphi \left(s\right)=\phi \left(\gamma \left(s\right)\right)\sqrt{{\left[x_1^{\prime }\left(s\right)\right]}^2+{\left[x_2^{\prime }\left(s\right)\right]}^2}$, and $t,s$ denote the parametrization value corresponding to the points $\mathbf{x},\mathbf{y}$.

It is well-known that the Equation (4) has a unique solution when D is not a unit circle, see ([36] theorem 7.38). In fact, the restrict on $\partial D$ can be removed by adding a constant c in (3) (see also [10]), i.e.,

$$u\left(x\right)={\int }_{\partial D}\mathsf{\Phi }(x,y)\phi \left(y\right)\mathrm{d}S\left(y\right)+c.$$

(7)

Now, let us consider the numerical solution of the Equation (6). Since the function $ln\left|\gamma \right(t)-\gamma (s\left)\right|$ has a singular point at $t=s$. First, we split the kernels as follows

$$-\frac{1}{2\pi }ln|\gamma \left(t\right)-\gamma \left(s\right)|=-\frac{1}{4\pi }ln\left(4{sin}^2\frac{t-s}{2}\right)+K(t,s).$$

(8)

When $t=s$, there is

$$\begin{array}{lll}K(t,t)& =& {\displaystyle \underset{s\to t}{lim}}\left\{-\frac{1}{2\pi }ln|\gamma \left(t\right)-\gamma \left(s\right)|+\frac{1}{4\pi }ln\left(4{sin}^2\frac{t-s}{2}\right)\right\}\\ & =& -\frac{1}{2\pi }{\displaystyle \underset{s\to t}{lim}}\left\{ln|\gamma \left(t\right)-\gamma \left(s\right)|-\frac{1}{2}ln\left(4{sin}^2\frac{t-s}{2}\right)\right\}\\ & =& -\frac{1}{2\pi }{\displaystyle \underset{s\to t}{lim}}ln\frac{\left|\gamma \right(t)-\gamma (s\left)\right|}{2|sin\frac{t-s}{2}|}\\ & =& -\frac{1}{2\pi }{\displaystyle \underset{s\to t}{lim}}ln\frac{\left|\gamma \right(t)-\gamma (s\left)\right|}{|t-s|}\\ & =& -\frac{1}{2\pi }ln\sqrt{{\left[x_1^{\prime }\left(t\right)\right]}^2+{\left[x_2^{\prime }\left(t\right)\right]}^2}.\end{array}$$

(9)

Thus the integral Equation (6) can be reduced to

$$-\frac{1}{4\pi }{\int }_0^{2\pi }\varphi \left(s\right)ln\left(4{sin}^2\frac{t-s}{2}\right)\mathrm{d}s+{\int }_0^{2\pi }\varphi \left(s\right)K(t,s)\mathrm{d}s=f\left(\gamma \left(t\right)\right).$$

(10)

In the collocation method, we choose $2n$ interpolation conditions distributed in $[0,2\pi )$, denoted by $t_j=\frac{j\pi }{n},j=0,1,2\dots 2n-1,$ and the second part of (4) turns to

$${\int }_0^{2\pi }\varphi \left(s\right)K(t,s)\mathrm{d}s\approx \frac{\pi }{n}{\displaystyle \sum _{j=0}^{2n-1}}K(t,t_j)\varphi \left(t_j\right).$$

(11)

For the weakly singular kernel $ln\left(4{sin}^2\frac{t-s}{2}\right)$, we use the Lagrange trigonometric polynomials basis ([36] Chapter 11.3)

$$X_n=\left\{{\displaystyle \sum _{j=0}^{n-1}}{a}_{j}cosjt+{\displaystyle \sum _{j=1}^n}{b}_{j}sinjt\right\}$$

and write $\varphi \left(t\right)\approx {\sum }_{j=0}^{2n-1}\varphi \left(t_j\right)L_j\left(t\right)$. The Lagrange interpolation basis ([36] Chapter 11.3) functions $L_j,j=0,1,2\dots 2n-1,$ are given by

$$L_j\left(t\right)=\frac{1}{2n}sinn(t-t_j)cot\frac{t-t_j}{2}.$$

(12)

We have

$$\begin{array}{ll}& -\frac{1}{4\pi }{\int }_0^{2\pi }\varphi \left(s\right)ln\left(4{sin}^2\frac{t-s}{2}\right)\mathrm{d}s\\ \approx & -\frac{1}{4\pi }\sum _{j=0}^{2n-1}{\int }_0^{2\pi }{L}_j\left(s\right)ln\left(4{sin}^2\frac{t-s}{2}\right)\mathrm{d}s\\ =& \sum _{j=0}^{2n-1}\varphi \left(t_j\right)R_j\left(t\right),\end{array}$$

(13)

where

$$\begin{array}{lll}{R}_j\left(t\right)& =& -\frac{1}{4\pi }{\int }_0^{2\pi }{L}_j\left(s\right)ln\left(4{sin}^2\frac{t-s}{2}\right)\mathrm{d}s\end{array}$$

(14)

$$\begin{array}{ll}=& \frac{1}{2n}\left\{\frac{1}{2n}cosn(t-t_j)+{\displaystyle \sum _{m=1}^{n-1}}cosm(t-t_j)\right\}.\end{array}$$

(15)

Now the Equation (2) is discrete to

$${\displaystyle \sum _{j=0}^{2n-1}}\varphi \left(t_j\right)\left(R_j\left(t\right)+\frac{\pi }{n}K(t,t_j)\right)=f\left(\gamma \left(t\right)\right).$$

In the collocation method, the values of $\varphi \left(t_j\right),j=0,1,2,\dots ,2n-1,$ are determined by solving a system of $2n$ equations consisting of the boundary conditions at the $2n$ collocation points

$$A\varphi =b.$$

(16)

The matrix A and the right-hand side b are given by the following:

$$A_{i,j}=R_{|i-j|}+\frac{\pi }{n}K(t_{i-1},t_{j-1}),b_i=f\left(\gamma \left(t_{i-1}\right)\right),i,j=1,\cdots ,2n,$$

(17)

where

$$R_{|i-j|}=\frac{1}{2n}\left\{\frac{{(-1)}^{i-j}}{2n}+{\displaystyle \sum _{m=1}^{n-1}}cos\frac{m(i-j)\pi }{n}\right\}.$$

3. Regularization Method

As following, the Tikhonov regularization method combining Morozov discrepancy principle is used to solve the system (16). In general, the right hand-side vector b of the system (16) is a perturbed vector $b^{\delta }$, thus we should consider the following perturbed equations

$$A{\varphi }^{\delta }=b^{\delta }.$$

(18)

More precisely, $b^{\delta }$ is the measured noisy data satisfying $b_i^{\delta }=b_i+\delta \mathrm{rand}\left(i\right)b_i,$ where $\delta$ is the percentage noise and the number rand(i) is a pseudo-random number drawn from the standard uniform distribution on the interval $[-1,1]$ generated by the Matlab commend $-1+2\mathrm{rand}(1,i)$.

The standard Tikhonov regularized solution of the system (18) is given by

$$(\alpha I+A^{\top }A){\varphi }_{\alpha }=A^{\top }{b}^{\delta },$$

(19)

where $\alpha$ is the regularization parameter chosen by the Morozov’s discrepancy principle. Although this problem is a direct problem. Here the Tikhonov regularization method is not used to regularize the aforementioned direct problem with noisy data, but to filter out the noise in the corresponding perturbed data.

We can see that the solution accuracy depends on the regularization parameter, $\alpha$, from the Equation (19). Thus how to choose an optimal regularization parameter $\alpha$ is crucial. We will choose the regularization parameter $\alpha$ by the Morozov’s discrepancy principle, which was developed in [4]. The computation of $\alpha$ can be carried out with Newton’s method as follows:

• Set $n=0$, and give an initial regularization parameter ${\alpha }_0>0$;
• Get ${\varphi }_{{\alpha }_n}^{\delta }$ from ${\displaystyle (A^{*}A+{\alpha }_{n}I){\varphi }_{{\alpha }_n}^{\delta }=A^{*}{b}^{\delta }}$;
• Get $\frac{\mathrm{d}}{\mathrm{dff}}{\varphi }_{{\alpha }_n}^{\delta }$, from $({\alpha }_{n}I+A^{*}A)\frac{\mathrm{d}}{\mathrm{dff}}{\varphi }_{{\alpha }_n}^{\delta }=-{\varphi }_{{\alpha }_n}^{\delta }$;
• Get ${\displaystyle F\left({\alpha }_n\right)}$ and ${\displaystyle F^{{}^{\prime }}\left({\alpha }_n\right)}$ by

  $${\displaystyle F\left({\alpha }_n\right)={∥A{\varphi }_{{\alpha }_n}^{\delta }-b^{\delta }∥}^2-{\delta }^2}$$

  and

  $${\displaystyle F^{{}^{\prime }}\left({\alpha }_n\right)=2{\alpha }_n{∥A\frac{\mathrm{d}}{\mathrm{dff}}{\varphi }_{{\alpha }_n}^{\delta }∥}^2+2{\alpha }_n^2{∥\frac{\mathrm{d}}{\mathrm{dff}}{\varphi }_{{\alpha }_n}^{\delta }∥}^2,}$$

  respectively.
• Set ${\alpha }_{n+1}={\alpha }_n-\frac{F\left({\alpha }_n\right)}{F^{{}^{\prime }}\left({\alpha }_n\right)}$. If $∥{\alpha }_{n+1}-{\alpha }_n∥<\epsilon (\epsilon \ll 1)$, end. Else, set $n=n+1$ and return to 2.

When the regularization parameter $\alpha$ is fixed with the value ${\alpha }^{*}$, we can obtain the corresponding regularized solution.

4. Numerical Examples and Discussion

In this section, we report some examples to demonstrate the effectiveness of our algorithm. The implementation of the algorithm is based on the Matlab software. From the classical theorem of the partial differential equations [33], we know that

$${∥u-u_{num}∥}_{H^1\left(D\right)}\le C{∥f-f^{\delta }∥}_{L^2(\partial D)}.$$

Thus we are only to filter out the noise in the corresponding perturbed data on the boundary, and it is enough to study the calculation results on the boundary.

In all examples, we choose $n=100,$ and denote the analytical and numerical solution by $f^{\left(an\right)}\left(\mathit{x}\right)$ and $f^{\left(num\right)}\left(\mathit{x}\right)$, respectively. Define the normalised relative root mean square (RMS) error on $\partial D$ as following

$$E_{\partial D}\left(f\right)=\frac{{\left\{\frac{1}{2n}{\sum }_{l=0}^{2n-1}{\left|f^{\left(an\right)}\left({\mathit{x}}_l\right)-f^{\left(num\right)}\left({\mathit{x}}_l\right)\right|}^2\right\}}^{\frac{1}{2}}}{{max}_{l\in \{0,1,\dots ,2n-1\}}\left|f^{\left(an\right)}\left({\mathit{x}}_l\right)\right|}.$$

Example 1.

D is a triangle-shaped domain (see Figure 1a), and $\partial D$ has a parameter expression [37]

$$\gamma \left(t\right)=(2+0.3cos3t)(cost,sint).$$

The temperature on the boundary is given by $f(x_1,x_2)=sin{x}_1{e}^{x_2}.$

Example 2.

D is an apple-shaped domain (see Figure 1b), and $\partial D$ has a parameter expression [16]

$$\gamma \left(t\right)=\left(1+\frac{sin2t}{5+4cost}\right)(cost,sint).$$

The temperature on the boundary is given by

$$f(x_1,x_2)=x_1^3-3x_1{x}_2^2+sin2x_1{e}^{2x_2}-cos{x}_2{e}^{x_1}.$$

Figure 2 shows the numerical results corresponding to Examples 1 and 2, and Figure 3 shows the relative errors of f. It can be seen that the numerical solutions are stable approximations of the analytic solution, and the numerical solution converges to the exact solution as the level of noise decreases.

Table 1. The values of the regularization parameter, $\alpha$, and the corresponding errors, $E_{\partial D}\left(f\right)$, obtained using the presented method, and various amounts of noise added in the boundary data, for Examples 1–4.

	Noise Level	$\mathit{\alpha }$	${\mathit{E}}_{\mathit{\partial }\mathit{D}}\left(\mathit{f}\right)$
Example 1	$1\%$	$2.13\times {10}^{-5}$	$5.78\times {10}^{-3}$
	$5\%$	$8.91\times {10}^{-5}$	$2.35\times {10}^{-2}$
	$10\%$	$2.42\times {10}^{-4}$	$6.05\times {10}^{-2}$
Example 2	$1\%$	$2.98\times {10}^{-5}$	$5.99\times {10}^{-3}$
	$5\%$	$1.38\times {10}^{-4}$	$2.70\times {10}^{-2}$
	$10\%$	$3.02\times {10}^{-4}$	$5.66\times {10}^{-2}$
Example 3	$1\%$	$5.11\times {10}^{-3}$	$5.89\times {10}^{-3}$
	$5\%$	$4.61\times {10}^{-2}$	$2.29\times {10}^{-2}$
	$10\%$	$1.16\times {10}^{-1}$	$5.77\times {10}^{-2}$

gives the values of the regularization parameter, $\alpha$, and the corresponding errors, $E_{\partial D}\left(f\right)$, obtained using the presented method, and various amounts of noise added in the boundary data, for Examples 1–3.

Example 3.

In the previous examples, the exact solutions are given, and the solutions are analytic. In this example, we consider D as an ellipse domain, and $\partial D$ has a parameter expression $\gamma \left(t\right)=(3cost,2sint)$. The temperature on the boundary is given by

$$f\left(t\right)=\left\{\begin{array}{cc}1-\frac{t}{\pi }& ,0\le t<\pi ,\\ \frac{t}{\pi }-1& ,\pi \le t<2\pi .\end{array}\right.$$

(20)

We can see that there is not an exact solution for this boundary value problem. Figure 4 shows the numerical results and the relative errors of f obtained using the presented method, and various amounts of noise added in the boundary data. It can be seen that the numerical solutions are stable approximations of the temperature on the boundary. We can see that the numerical solutions on the corner are less accurate than the case in which an analytic solution is available, especially when the noise level is becoming larger. We think it is reasonable and meets expectations.

In order to examine the efficiency of the presented method, we give a comparative study between the present method and the MFS.

Example 4.

D is a non-convex kite-shaped [16], and $\partial D$ described by the parametric representation (see Figure 5)

$$\gamma \left(t\right)=(cost+0.65cos2t-0.65,1.5sint).$$

The temperature on the boundary is given by $f(x_1,x_2)=x_1^2+x_2^2.$ We can see that there is not an exact solution for this boundary value problem.

Figure 6a,b show the numerical results obtained by the present method and the MFS. Figure 7a,b show the relative errors of f obtained with various amounts of noise added in the boundary data on $\partial D$. From these figures, we can see that the numerical solutions obtained by the present method are stable approximations of the temperature on the boundary. It also should be noticed that the presented method yields more accurate, convergent and stable numerical results than the MFS.

Example 5.

D is an apple-shaped domain, and $\partial D$ has a parameter expression

$$\gamma \left(t\right)=\left(1+\frac{sin2t}{5+4cost}\right)(cost,sint).$$

The temperature on the boundary is given by

$$f\left(t\right)=\left\{\begin{array}{cc}1& ,0\le t<\pi ,\\ -1& ,\pi \le t<2\pi .\end{array}\right.$$

(21)

In this example, the temperature on the boundary has a jump. Figure 8a,b show the numerical results obtained by the present method and the MFS. In order to check the efficiency of the present method (see Figure 9), we give a comparative study between the present method (PM) and the MFS.

Table 2. The values of the regularization parameter, $\alpha$, and the corresponding errors, $E_{\partial D}\left(f\right)$, obtained using the presented method, and various amounts of noise added in the boundary data, for Examples 4 and 5.

	Noise Level	${\mathit{\alpha }}_{\mathit{M}\mathit{F}\mathit{S}}$	Err(MFS)	${\mathit{\alpha }}_{\mathit{P}\mathit{M}}$	Err(PM)
Example 4	$1\%$	$1.59\times {10}^{-3}$	$5.72\times {10}^{-3}$	$8.19\times {10}^{-5}$	$5.55\times {10}^{-3}$
	$5\%$	$1.10\times {10}^{-1}$	$2.98\times {10}^{-2}$	$4.75\times {10}^{-4}$	$2.90\times {10}^{-2}$
	$10\%$	$2.81\times {10}^{-1}$	$5.73\times {10}^{-2}$	$1.00\times {10}^{-3}$	$5.43\times {10}^{-2}$
	$20\%$	$8.66\times {10}^{-1}$	$1.20\times {10}^{-1}$	$2.40\times {10}^{-3}$	$1.01\times {10}^{-1}$
Example 5	$1\%$	$3.55\times {10}^{-16}$	$8.47\times {10}^{-2}$	$3.48\times {10}^{-5}$	$5.78\times {10}^{-3}$
	$5\%$	$3.55\times {10}^{-16}$	$8.47\times {10}^{-2}$	$1.85\times {10}^{-4}$	$2.93\times {10}^{-2}$
	$10\%$	$3.55\times {10}^{-16}$	$8.47\times {10}^{-2}$	$4.06\times {10}^{-4}$	$6.03\times {10}^{-2}$
	$20\%$	$1.08\times {10}^{-9}$	$1.17\times {10}^{-1}$	$8.67\times {10}^{-4}$	$1.10\times {10}^{-1}$

shows the values of the regularization parameter, $\alpha$, and the corresponding errors, $E_{\partial D}\left(f\right)$, obtained using the presented method and MFS, and various amounts of noise added in the boundary data, for Examples 4 and 5. The condition numbers of the MFS are $5.74\times {10}^{18}$ and $1.22\times {10}^{19}$, but the present method can reduce the condition numbers to $1.66$ and $1.49$. From this figure, we can see that the numerical solutions obtained by the present method are stable approximations of the temperature on the boundary. It also should be noticed that the presented method yields more accurate, convergent and stable numerical results than the MFS. Especially in the jump points, the presented method gives more stable numerical results than the MFS.

Especially from these figures, we can see that the present method also gives a accurate approximation when the boundary value is not a harmonic function, even with jump discontinuities.

5. Conclusions

In this paper, a numerical algorithm is presented to solve the isotropic heat conduction equation with noisy boundary conditions. In the numerical examples, both convex and non-convex geometries with smooth boundaries are considered. For these cases, the present method can well filter the noise in the boundary data. In this process, the Tikhonov regularization method combining Morozov’s discrepancy principle is used to select an appropriate regularization parameter. The results are very accurate for perturbed/noisy data. In the future work, we will report this method for the acoustic radiation [38,39,40], seismic wave scattering [41] and more.