A Virtual Boundary Element Method for Three-Dimensional Inverse Heat Conduction Problems in Orthotropic Media

This paper aims to apply a virtual boundary element method (VBEM) to solve the inverse problems of three-dimensional heat conduction in orthotropic media. This method avoids the singular integrations in the conventional boundary element method, and can be treated as a potential approach for solving the inverse problems of the heat conduction owing to the boundary-only discretization and semi-analytical algorithm. When the VBEM is applied to the inverse problems, the numerical instability may occur if a virtual boundary is not properly chosen. The method encounters a highly illconditioned matrix for the larger distance between the physical boundary and the virtual boundary, and otherwise is hard to avoid the singularity of the source point. Thus, it must adopt an appropriate regularization method to deal with the ill-posed systems of inverse problems. In this study, the VBEM and different regularization techniques are combined to model the inverse problem of three-dimensional heat conduction in orthotropic media. The proper regularization techniques not only make the virtual boundary to be allocated freer, but also solve the ill-conditioned equation of the inverse problem. Numerical examples demonstrate that the proposed method is efficient, accurate and numerically stable for solving the inverse problems of three-dimensional heat conduction in orthotropic media.


Introduction
The inverse heat conduction problem has received considerable attention in the past few decades [Weber (1981); Raynaud and Bransier (2007); Babaniyi, Oberai and Barbone (2018); Wang, Hua and Liu (2018)].In the fields of engineering, there are many orthotropic materials, such as woods, sheet metals, cables and new advanced composites [Jin, Zheng and Marin (2010)].It is necessary to consider the thermal conductivity with direction when modelling the heat transfer in these orthotropic media [Mera, Elliott, Ingham et al. (2003); Ohmichi, Noda and Sumi (2017); Somasundharam and Reddy (2018)].Thus, the heat conduction in orthotropic materials has widespread applications in many branches of engineering and its understanding is of great significance.In some engineering problems, it is impractical to directly survey the heat flux or the temperature of the boundary.For instance, it is difficult to directly survey the heat on the surface of a hypersonic vehicle's thermal protection system when the vehicle reenter from the atmosphere [Yang, Yue and Yang (2017)]; the heat flux at the surface of a wall cannot be obtained directly when the surface is subjected to fire [Jarny, Ozisik and Bardon (1991)]; it could be extremely hard to survey the heat flux or the temperature at the tool-work interface of machine cutting [Mohammadiun (2016)].These examples can be treated as the inverse problems of the heat conduction in which the Neumann or Dirichlet data on certain boundaries are unknown [Amirfakhrian, Arghand and Kansa (2016); Cui, Zhao, Xu et al. (2017)].The primary objectives of the inverse heat conduction problems are estimating unknown boundary conditions or the initial conditions with the help of the additional information of measurement [Gu, Chen and Fu (2014)].However, the inverse problems of heat conduction are mathematically illconditioned, the arbitrarily large changes may occur in the numerical solution if there are small deviations in the measured data [Ding and Sun (2015)].The boundary element methods (BEM) have been frequently used to solve these inverse heat conduction problems, it is efficient when the boundary integrals can be obtained accurately [Wang, Jiang and Yang (2016); Mohasseb, Moradi, Sokhansefat et al. (2017)].But the BEM will encounter the issue of the calculation of singular integral if the source point lies close to the integral element or in the integral element, especially for the threedimensional high-order geometry element [Gao (2010)].Like the method of fundamental solution [Sun and He (2017); Wang, Liu and Qu (2018)], the virtual boundary element method (VBEM) is developed to overcome this obstacle [Xu and Yang (2014)].The VBEM introduces the virtual boundary, which can avoid the calculation of singular integral.Additionally, the VBEM is more accurate and efficient than other methods, and the program is also relatively concise [Ling and Yang (2015)].When the VBEM is applied to the inverse problems of three-dimensional heat conduction in orthotropic media, the position of the virtual boundary will influence the precision of the numerical results.Meanwhile, the VBEM discretization of inverse problems will yield the severely ill-conditioned matrix, so it is necessary to solve these problems effectively [Chen, Karageorghis and Li (2016)].The main purpose of this study is to extend the VBEM to solve the inverse problems of three-dimensional heat conduction in orthotropic media.The Gauss elimination method, the conjugate gradient method (CGM) [Chen, Yu, Zhou et al. (2017)], and the Tikhonov regularization technique with the Lcurve and GCV criterion are applied to the solution of the ill-conditioned equations [Berntsson, Kozlov, Mpinganzima et al. (2017)].The results indicate that the CGM and the Tikhonov regularization technique are insensitive to the virtual boundary.Among these methods, the Tikhonov regularization technique combined with the L-curve criterion can effectively improve the accuracy of the results even though the relatively high level of noise is added into the boundary data.The rest of this paper is organized as follows.The second section demonstrates the governing equation and the boundary conditions of the inverse problems of threedimensional heat conduction in orthotropic media.The VBEM for the inverse problems of three-dimensional heat conduction is reviewed in the third section.In the fourth section, the formulations of the VBEM are in conjunction with the regularization techniques to solve the ill-conditioned problems.Three numerical test problems are investigated in the fifth section.Finally, the conclusions are drawn in the sixth section.

Mathematical formulation of the problem
Consider a three-dimensional orthotropic media in a domain 3 R  bounded by a surface  , two parts are assumed to be contained in the surface, which can be defined as 12 =    , where 12 ,     and 12   =  .In this study, the orthotropic steady heat conduction applications are involved in the absence of heat generation and other heat transfer modes.The governing equation of the problem can be described by the partial differential equation comply with the following boundary conditions where ( ) 1,2,3 is the thermal conductivity tensor, () u x denotes the temperature at the  indicates the known boundary on which the temperature and normal heat flux are available, the overline quantities () u x and () q x are the presetting values on the accessible boundary 1  .
The normal heat flux through the boundary  is given by The fundamental solution of governing Eq. ( 1) is where    x represents the collocation point, and 1 2 3 ( , , ) y y y = y is the source point.
In the boundary conditions Eq. (2) and Eq.(3), it can be noted that the boundary 1  is over-specified on which both the temperature and normal heat flux are known, while both the temperature and normal heat flux are unknown on the remaining boundary 2  , and have to be determined.

Virtual boundary element method
The VBEM can be treated as an extension of the BEM, it introduces the infinite domain   to substitute the original area  , and the virtual boundary Γ can be established in the infinite domain   .There is a certain distance between the virtual boundary Γ and the physical boundary  .The virtual boundary Γ makes the distribution density function ()  c as middle variable through continuous way.Then, the method can determine the density function on the virtual boundary Γ through the given physical quantity on the boundary  , so as to eliminate the calculation of singular integral and nearly singular integral in the conventional BEM.A virtual boundary integral equation of the potential problem can be described as follows [Yang, Yue and Yang (2017)]: The integral equation of the interior point of original area  is expressed as If the virtual boundary is discretized by N elements and the constant element is adopted.
The unknown values 12 ( , , ... , ) on the virtual boundary can be obtained by choosing M collocation points on the boundary 1  .Suppose that Eq. ( 5) and Eq. ( 6) are contented with the given boundary conditions, then we can get 2MN  algebra equations, which can be described in the form of matrix as b represents a vector of 2M dimension, it is formed at the collocation points on the physical boundary and based on the given boundary condition.This algebra system is generally ill-conditioned, the traditional Gaussian elimination method and the methods suitable for direct problems cannot efficiently solve the linear system Eq.( 8).However, some proper regularization techniques can be used to their calculation.Two of them, namely the conjugate gradient method (CGM) and the Tikhonov regularization technique, are briefly presented and employed to solve Eq. ( 8) in the following section.

Conjugate gradient method
The CGM is an efficient method of solving the inverse problem of the heat conduction by setting a regularizing stopping criterion; it is a stable scheme to solve the linear algebraic equations [Yang, Chen, Ling et al. (2017)].In this method, Eq. ( 8) can be transformed to a normal linear system: Therefore, the problem has been converted to solve Eq. ( 9), the algorithm is demonstrated as follows: 1. Set an initial guess 0 c , determine the residual 0 0 1 =− r Dc b and then set 00 = tr .
2. For the k th ( 0,1, 2,... k = ) step, the procedure of iterations is formulated as Increase k by one and go to Step 2 until a given stopping criterion where 0   is a constant.In this paper,  is taken as 1e-6.In Eq. ( 11), r denotes the gradient, t is the direction of descent,  is the search step size, and  is the conjugate coefficient.

Tikhonov regularization technique
The Tikhonov regularization technique is an effective and powerful tool to approximate a solution of the ill-posed Eq. ( 8) by minimizing the following function [Gu, Chen and Fu (2014)]: where 2  denotes the Euclidean norm, 0   is the regularization parameter, it usually controls the relative weight of the penalty term in the function.If the parameter  is equal to zero, Eq. ( 12) will degenerate to the least squares problem.More details can be referred to Ref. [Mera, Elliott, Ingham et al. (2003)].The determination of regularization parameter  will be introduced in the next section.
The singular value decomposition (SVD) works well for direct problems but usually difficult to achieve a stable and accurate numerical solution to Eq. ( 8).Therefore, the Tikhonov regularization technique has been developed for solving the ill-conditioning problem with the help of the SVD.The SVD of the matrix A is given by 1 where Based on the SVD, the Tikhonov regularized solution can be obtained by 1 () where the filter factors i f can be expressed as the details of the above regularization technique can be found in Wei et al. [Wei, Hon and Ling (2007)].All the computational codes are written in MATLAB and the code developed by Hansen based on the SVD for solving the ill-posed Eq. ( 8) has been adopted in our computations [Hansen (1994)].

Determination of regularization parameters
An appropriate selection of the regularization parameter  is very important for the efficiency of regularization methods used in the ill-conditioned interpolation matrix equation.But the choice of an optimal regularization parameter is still an open question under the intensive research.To deal with this issue, there are two widely used techniques, which are the L-curve criterion and the generalized cross-validation criterion.

L-curve criterion
The L-curve criterion is often suggested to determine a suitable value of the regularization parameters when solving ill-conditioned problems, which is defined as a continuous curve consisting of all points This curve is usually called the L-curve since it looks like the letter "L" in shape for a large class of problems.According to the basic idea of this criterion, the best regularization parameter should lie on the "corner" of the L-curve where the perturbation and the regularization errors are minimized.This "corner" can be perceived as a point at which the curvature of the L-curve is maximized [Gu, Chen and Fu (2014)].

Generalized cross-validation criterion
As another alternative method, the generalized cross-validation (GCV) criterion can also be adopted to choose the regularization parameters.This method is to minimize the generalized cross-validation function [Wahba (1977)]: for real 0   , here I A denotes a matrix which produces the regularized solution  w , N I represents the identity matrix of order N.

Numerical examples and discussions
To evaluate the performance of the VBEM combined with the CGM and the Tikhonov regularization techniques (VBEM-CGM, VBEM-Tikhonov-GCV and VBEM-Tikhonov-LC), three numerical examples are tested.It is necessary to carefully investigate the effect of the three regularization methods and the scheme stability in relation to the level of noise.In the following test cases, the simulated data with noise are yielded by the following formula: ( ) where i b represents the exact boundary data, T is a standard normally distributed random variable, and  denotes the noise level [Abe, Fujii and Yoshimura (2017)].The random variable T is implemented by using the MATLAB function 'randn ()' in the tests.

Example 1
A central hollow sphere is considered in this example, which is shown in Fig. 1.The inner radius of the spherical shell is 1 1 r = , and the outer radius is 2 2 r = .The accessible boundary 1  is set as the outer surface.
where 1  First, the necessity for introducing the regularization methods should be investigated when an ill-posed Cauchy problem is solved.To do this, the virtual boundary is located at a spherical shell with inner radius of 0.5 and outer radius of 2+d, where d denotes the distance between the physical boundary and the virtual boundary.200 nodes are selected on the outer and inner virtual boundaries.Fig. 2(a) displays the numerical solutions of the temperatures on the circle obtained by using the VBEM combined with the conventional Gauss elimination method under different distance between the virtual and the physical boundary with 1% noisy Cauchy data.A fair comparison is made among the numerical solutions obtained by using the CGM ( 14 e  =− , iter=41), the Tikhonov-LC, and the Tikhonov-GCV regularization method, as illustrated in Fig. 2 1 that the result using the Gauss elimination is highly oscillatory and inaccurate; the result could not be adopted as an approximate solution at all.Additionally, it is reported that the least-squares method and LU factorization could not obtain the accurate and stable solution either.Standard methods could not produce the accurate results for the noisy data and, therefore, the regularization technique should be introduced to solve this problem.However, as can be seen from Fig. 2(b), Fig. 2(c) and Fig. 2(d), the CGM, the Tikhonov-LC, and the Tikhonov-GCV regularization method could yield very accurate results for a wide variation of the distance (d=1, d=10, d=100) between the virtual boundary and the physical boundary.Even though the distance reached 100, these three regularization methods still can reach a high accuracy level.This demonstrates that the VBEM combined with the three regularization methods are the stable and effective tools to solve the inverse problem of the three-dimensional heat conduction.
To compare the CGM, Tikhonov-LC, and Tikhonov-GCV regularization techniques, Fig. 3 presents the root mean square relative errors (RMSREs) of the temperatures at 200 test points uniformly distributed on the inner surface obtained by using the three regularization techniques with various noise levels added into the data.The accuracy of the CGM decreases dramatically as the noise level increases, but the Tikhonov-LC and Tikhonov-GCV are accurate and stable for various noise levels, as shown in Fig. 3. Furthermore, it can also be observed that the Tikhonov-LC is slightly more stable than the Tikhonov-GCV.Therefore, the Tikhonov-LC regularization technique will be employed to deal with the ill-posed problem in the following examples.

Figure 3:
The RMSREs of the temperatures on the inner surface obtained by using the CGM, Tikhonov-LC, Tikhonov-GCV regularization techniques with various noise levels added into the data

Example 2
In this example, the VBEM-Tikhonov-LC is applied to the solution of the inverse heat conduction problem in a complicated domain with the boundary  as shown in Fig. 4, prescribed by the following spherical parametric equation: The accessible boundary 1  is specified as the upper half of the surface, namely, ( ) The distribution of the analytical temperature is given by The orthotropic parameters are 1 To solve the problem numerically, 400 boundary nodes are uniformly placed on a sphere with a radius of 5. Fig. 5 shows the numerical results for the temperatures at 100 test points on the circle obtained by using the VBEM-Tikhonov-LC with 3% noise level.The corresponding results are summarized in Tab. 2. It can be observed from Fig. 5 and Tab. 2 that the numerical solutions agree well with the according analytical solutions, and the maximum absolute error is 2 4.22 10 −  in Fig. 5.To further investigate the accuracy of the proposed scheme for the reconstruction of the temperature distribution on the under-specified boundary, Fig. 6 and Fig. 7 give the distributions of the analytical and numerical results for the temperatures on the lower half of the surface with 5% noisy Cauchy data.It can be seen from these figures that the numerical solutions agree quite well with the exact solutions, indicating that the VBEM, together with the Tikhonov-LC, is accurate and stable for the inverse Cauchy heat conduction problems on the irregular three-dimensional domain even for a relatively high amount of noise (5%) added into the data.In the calculation, the regularization parameter is selected as 0.011723 via the L-curve method, Fig. 8 shows the L-curve for the Tikhonov regularization technique with 3% noisy data.

Example 3
As the last example, a cylinder with non-uniform thickness is considered, as shown in Fig. 9, the inner boundary is a cylindrical surface with the radius 1 1 r = , the outer boundary is hyperboloid of one sheet which is defined parametrically as where , respectively.In this example, the boundary conditions are specified on the inner and outer surfaces of the cylinder, and both the temperature and the heat flux of the lower and upper surfaces are unknown.The exact solution is given by where 1 To tackle this problem numerically, the virtual boundary is taken as a similar geometry with that in Fig. 9.The inner virtual boundary is a cylindrical surface with the radius 1 0.5 r = , the outer virtual boundary is hyperboloid of one sheet which is defined parametrically as ( )

22
( , ) 2.5 1 0.4 cos , 2.5 1 0.4 sin ,  , respectively. 200, 200, 100, and 100 nodes are placed on the outside face, inside face, upper face and lower face, respectively.Fig. 10 shows the numerical results of the temperatures on the circle by using the VBEM combined with the Tikhonov-LC regularization method with varying noise levels.The numerical solutions converge to the exact solution as the noise level decreases.Furthermore, the numerical solutions still agree very well with the analytical solutions even if a relatively high noise level (5%) is added into the data.In addition, the RMSREs of the noise levels 1%, 3%, and 5% are , respectively.This clearly illustrates the accuracy and stability of the proposed numerical scheme.Finally, the temperatures of the upper and lower surfaces are reconstructed via the VBEM combined with the Tikhonov-LC regularization method under 3% noise level.Fig. 11 gives the L-curve plot in the determination of the regularization parameter of the Tikhonov regularization technique.Fig. 12 and Fig. 13 displays the exact solution profile and the numerical solution surfaces of the temperatures on the lower and upper surfaces with 3% noisy Cauchy data, respectively.It can be seen from Fig. 12 and Fig. 13 that the agreement between the numerical and analytical solutions is quite close, indicating the present method is effective, stable and accurate for solving the inverse problems the heat conduction in the complicated three-dimensional orthotropic media.In this paper, it is the first attempt to apply a regularized virtual boundary element method for the inverse problems of three-dimensional heat conduction in orthotropic media.The resulting ill-posed system of algebra equations is regularized by employing the regularization techniques, which are the CGM, the Tikhonov-GCV and the Tikhonov-LC.Numerical solutions for three test examples involving both smooth and segmented smooth geometries have been investigated.It will be more arbitrary to set the location of the virtual boundary when using the three the regularization techniques.Furthermore, compared with the VBEM-CGM and the VBEM-Tikhonov-GCV, the VBEM-Tikhonov-LC is more accurate and stable with respect to the noise level in the input data.
Overall, it can be concluded that the VBEM-Tikhonov-LC is more efficient when involving the same noise level of Cauchy data.The proposed scheme is computationally robust, accurate, stable and efficient, it can be considered as a competitive candidate to the existing methods for solving three-dimensional inverse heat conduction problem in orthotropic media.
cosines of the unit outward normal vector at the boundary point x. 1 matrix with non-negative elements, i  are called the singular values of matrix A .

Figure 1 :Figure 2 :
Figure 1: Steady state heat conduction in a central hollow sphere: (a) Physical boundary; (b) Virtual boundary.

Figure 6 :Figure 7 :
Figure 6: Profile of analytical solution for the temperature on the lower half of the surface

Figure 8 :
Figure 8: The L-curve for the Tikhonov regularization technique with 3% noisy data

Figure 9 :
Figure 9: Steady state heat conduction in a cylinder with non-uniform thickness: (a) physical boundary; (b) virtual boundary.

Table 1 :
The numerical solutions and relative errors of the temperatures at several inner points on the circle

Table 2 :
Comparison between the numerical and analytical solutions of the temperatures on the circle