Evaluating 2D domain integrals by sinh transformation for transient heat conduction problem

The time-dependent boundary element method (BEM) is widely used to solve transient heat conduction problem. The integrands of the domain integrals in the time-dependent BEM are close to singular when small time step is used. A straightforward computation of these integrals using Gaussian quadrature will produce large errors. To evaluate domain integrals accurately and efficiently, the sinh transformation combined with element subdivision method is proposed in this paper. With the sinh transformation, the integrands of the domain integrals become smoother. Meanwhile, the Gaussian points are shifted towards the source point due to the element subdivision. Thus, the 2D domain integrals in the time-dependent BEM can be computed accurately. Numerical examples have demonstrated the accuracy and efficiency of the proposed method.


Introduction
The numerical method is a successful approach for solving complicated transient heat conduction problems. The boundary element method (BEM) [1][2][3][4][5][6][7][8][9][10] is one of the most used numerical method for transient heat conduction analysis. However, when using small time step, the numerical results of the BEM will be unstable. Dong et al [11][12][13] have analyzed the unstable problem and found that inaccurate evaluation of the domain integrals was the main reason for this problem. They have proposed (α, β) transformation and (α, β, γ) transformation for 2D and 3D domain integrals, respectively.
The integrand of the domain integral in the time-dependent BEM is regular and there is no singularity when the source point is on the integration element. But when the time step is very small, the integrand of the domain integral will vary dramatically around the source point. Thus, a straightforward use of ordinary numerical integration method will produce large errors.
In order to compute domain integrals accurately, the sinh transformation combined with element subdivision method is proposed in this paper. The sinh transformation was used for computing the nearly singular integrals [14][15][16]. Johnston et al [17] firstly introduced the sinh transformation for the nearly singular integrals in detail. Then Gu et al [18] improve this method. The sinh transformation can make the integrand of the domain integral smoother and shift the Gaussian points towards the source point. Thus, more accurate numerical results for the domain integral can be obtained. In this paper, the minimum distance from the filed point to the integration element is replaced by the time step. Combined with element subdivision method, the integral accuracy can be further improved. The paper is outlined as follows. The boundary integral equation and the domain integral are presented in Section 2. The sinh transformation is introduced in Section 3. Two numerical examples are shown in Section 4. The conclusion is given in Section 5.

Time-dependent boundary integral equation
The time-dependent boundary integral equation for transient heat conduction problem is where ξ and x are the source point and the field point, respectively. c(ξ) is a constant depending on geometrical shape at point ξ. k is the diffusion coefficient. u * and q * are the time-dependent fundamental solution and its derivative.

2D domain integral
The 2D domain integral in Eq. (1) can be written as where u 0 is the initial temperature and it is a constant in the above domain integral. The expression of the time-dependent fundamental solution u * is where r is distance between the source point and the field point. τ is the time step. Figure 1 shows the variation of u * with respect to r for different time steps. It can be seen that the time-dependent fundamental solution u * vary dramatically around the source point as the time step decreases. Therefore, a straightforward use of Gaussian quadrature cannot lead to accurate numerical results.

Sinh transformation
The sinh transformation is mainly used to evaluate nearly singular integrals. In this paper, the minimum distance in traditional sinh transformation is replaced by the time step τ. The expression of the sinh transformation can be written as 1 The Jacobian of the sinh transformation is as follow: Figure 2 shows the variation of u * with respect to r after sinh transformation. One can see that compared with figure 1, the time-dependent fundamental solution u * vary gently around the source point even for very small-time step. It can be illustrated that the 2D domain integrals can be computed accurately using the sinh transformation.
The diffusion coefficient k is assumed to be 1. The sinh transformation with 5×10 Gaussian points is used on the sub-triangles and 5×5 point Gaussian quadrature is used on the sub-quadrangles.
It can be seen from table 1 that a straightforward Gaussian quadrature would produce large relative errors when small time steps are used. However, with the proposed method, high integral accuracy can be obtained even for τ = 0.000001.

Example 2
In this example, the integration element is the same with that in previous example. The relative errors for different position of source points are shown in Table 2. '(0.0, 0.0)' is the coordinate of the source point.
One can see that for different position of source points, high integral accuracy can be obtained by the proposed method and the numerical results are stable.