Moving boundary analysis in heat conduction with multilayer composites by finite block method

Engineering Analysis with Boundary Elements 89 (2018) 36–44

Contents lists available at ScienceDirect

Engineering Analysis with Boundary Elements journal homepage: www.elsevier.com/locate/enganabound

Moving boundary analysis in heat conduction with multilayer composites by finite block method M. Lei a,b, M. Li b,∗, P.H. Wen c,∗, C.G. Bailey c a

Department of Mathematics, City University of Hong Kong, Hong Kong, China College of Mathematics, Taiyuan University of Technology, Taiyuan, China c School of Engineering and Materials Science, Queen Mary University of London, London E14NS, UK b

a r t i c l e

i n f o

Keywords: Moving boundary Inverse problem of heat transfer Inhomogeneous heat transfer Finite block method Lagrange interpolation Roots of Chebyshev polynomial

a b s t r a c t An inverse reconstruction investigation is presented to determine the inner boundary location (corrosion points) for the heat transfer in composite walls from measurement data on exterior boundary. Finite Block Method (FBM) is utilized in this paper to deal with transient heat problems across the multilayered composite walls. Starting from one-dimensional problems, Lagrange interpolation with equally spaced nodes is applied to create first order differential matrices and thereafter the higher order differential matrices are obtained. Then combining with mapping technique, physical domain is mapped into a normalized domain for two-dimensional or three-dimensional problems with 8 seeds or 20 seeds respectively. Both time-spatial approach and Laplace transform technique with Durbin’s inversion method are employed in the simulating procedure. In addition, roots of Chebyshev polynomial of first kind are considered in FBM for the first time, which can improve the degree of convergence significantly. Three numerical examples are presented to validate the accuracy of FBM. Comparisons between Finite Element Method (FEM), FBM and Point Collocation Method (PCM) are demonstrated respectively. Numerical observation indicates that FBM has much higher degree of accuracy even with few collocation points. © 2018 Elsevier Ltd. All rights reserved.

1. Introduction Detecting the corrosion boundary of multilayer materials occurs extensively in engineering applications such as metallurgy. In the procedure of iron-steel-making, corrosion will take place on the inner surface of steel-smelting furnace which is always made of multilayer materials with different physical properties. We need to monitor the corrosion timely to prevent accidents. The corrosion degree inside steel-smelting furnace can be simulated through outermost layer information such as temperature or heat flux which may be measured easily. This kind of problem is always ill-posed and belongs to inverse problems. The identification of corrosion boundary in heat transfer problems has been studied numerically by Aparicio and Atkinson [1], Bryan and Caudill [2,3] and Fredman [4]. Lots of numerical algorithms such as Finite Element Method (FEM), Finite Difference Method (FDM) as well as Boundary Element Method (BEM) are available for inverse problems in engineering, see Liu and Zhang [5], Raynaud and Bransier [6], Aliabadi and Wen [7]. Although FEM is considered as the most powerful and well developed tool in numerical engineering, advanced numerical methods have become more attractive recently such as meshless approaches. In the last decade, meshless methods based on interpolation techniques

∗

like Radial Basis Functions (RBF) and Moving Least Square (MLS) provide an efficient tool in engineering analysis, see Powell [8], Jackson [9], and Reinhard Farwig [10]. Method of Fundamental Solutions (MFS) developed as a useful technique in numerical computation by Mathon and Johnston [11] and then Golberg and Chen [12] extended its application further. Hon and Wei [13] utilized MFS to solve inverse heat transfer problems. Wei and Yamamoto [14] employed this method to identify moving boundary for heat conduction problems with Cauchy conditions. Li and Wei [15] proved the uniqueness of corresponding solutions. However, it is difficult to make further convergence study. Same problems with Stefan–Boltzmann conditions had been studied by Hu and Chen [16] and the existence of the corresponding solutions had been proved. Niu et al. [17,18] employed time-spatial method as well as Kansa method to solve inverse heat transfer problems in inhomogeneous media. In their study, time is treated as an extra coordinate and the method of discrete Tikhonov regularization is used due to ill-posed characteristics of inverse problems. Differential Quadrature Method (DQM) developed by Bellman [19] is a powerful numerical method for lots of boundary problems. It is applied extensively in engineering and sciences by Bert et al. [20,21]. Liew et al. [22–26] proposed a lot of investigations on DQM in plate vibration as well as bending. When they use same number of discrete

Corresponding authors. E-mail addresses: [email protected] (M. Li), [email protected] (P.H. Wen).

https://doi.org/10.1016/j.enganabound.2018.01.009 Received 4 November 2016; Received in revised form 12 January 2018; Accepted 13 January 2018 0955-7997/© 2018 Elsevier Ltd. All rights reserved.

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

nodes, it can be found that DQM has better convergent results than FEM. Bert and Malik [27] gave a comprehensive review on DQM as well as its application. Recently Wen et al. [28], Li et al. [29,30] proposed a new meshless method called Finite Block Method (FBM) which is based on point collocation concept. It has been applied to solve heat transfer and elastodynamic problems in both 2D and 3D in Functionally Graded Materials (FGMs) with “star” performance. The novelty of FBM is that both first order and high order partial differential matrices in physical domain are based on the first order differential matrices constructed by Lagrange polynomial in normalized domain. Meanwhile, partial differential matrices of higher order are not essential in the governing equations with non-constant coefficients such as thermo-elasticity problems in inhomogeneous media. At the interface between two adjacent blocks in similar material, all stress components are found to be continuous [30]. The aim of this paper is to propose an inverse reconstruction procedure to determine the inner boundary (corrosion points) location in the problem of heat conduction for composite walls from measurement data of both temperature and heat flux on exterior boundary. FBM is extended to reconstruct the shape of corrosion boundary line (or surface) for 2D (or 3D) problems. First, with Lagrange interpolation where collocation points are distributed uniformly along a straight line, the first order differential matrix can be constructed for 1D problem, after which differential matrices of higher order can be obtained directly. Partial differential matrices for 2D and 3D problems can be evaluated similarly. Object with irregular boundary configuration in physical domain should be divided into quadratic blocks. With several quadratic shape functions, each block is mapped into a normalized square. In this procedure, elements with 8 (or 20) seeds are employed for 2D (or 3D) problems. Thereafter, the continuity conditions of heat flux and temperature along interface between two blocks should be satisfied and the partial differential operators in each block should satisfy a strong form too. Two approaches are proposed to deal with time dependent behavior in this paper, i.e. time-space approach (TS) where we do not distinguish time and space variables, and Laplace transform approach (LT), where the dependence on time variable t is reduced to the dependence on the transform variable s . For the second approach, in physical domain, the Durbin’s inversion technique [31] for Laplace transformation is applied to get the nodal values. The paper is organized as follows. Fundamentals of FBM which can be used directly to solve time independent problems are given in Sections 2 and 3. Time-space approach and Laplace transform approach are shown in Section 4 to deal with time dependent behavior. Three numerical examples are analyzed to demonstrate the procedure of FBM in identifying the moving boundary. Comparisons among FBM (TS), FBM (LT), FEM and RBF have been shown in Section 5. Finally a conclusion is presented in Section 6.

Fig. 1. Uniformly distributed nodes in normalized domain.

where uj = u(𝜉 j ) is the nodal value. Through (2), the first order derivative of u(𝜉) can be obtained 𝑀 𝑀 𝑀 𝑀 ∏ )−1 ∑ ( ) 𝑑𝑢 ∑ ∏ ( 𝑢𝑗 𝜉𝑗 − 𝜉𝑘 𝜉 − 𝜉ℎ . = 𝑑𝜉 𝑗=1 𝑘=1, 𝑖≠𝑗 ℎ=1,ℎ≠𝑗,ℎ≠𝑖

(3)

𝑘≠𝑗

For all nodes, the first order derivative can be expressed, in matrix form, as 𝐔𝜉 = 𝐃0 𝐮,

(4)

where u = [u1 ,u2 ,....uM ]T is vector of nodal values, 𝐔𝜉 = [𝑢′1 , 𝑢′2 , ..., 𝑢′𝑀 ]𝑇 is vector of the first order derivative of u(𝜉) at each node. D0 is the first order differential matrix. For m-th order derivative, we have 𝐔(𝜉𝑚) = 𝐃𝑚 𝐮, 𝑚 > 0, 0

(5)

𝑚) 𝑇 where 𝐔(𝜉𝑚) = [𝑢(1𝑚) , 𝑢(2𝑚) , ..., 𝑢(𝑀 ] is vector of m-th order derivative of nodal value.

2.2. Multi-dimensions 2. Lagrange interpolations for FBM

Similar to 1D problems in normalized domain, we have 𝑗−1) 𝜂𝑗 = −1 + 2(𝑁−1 , 𝑗 = 1, 2....𝑁,where N indicates the number of discrete points along 𝜂-axis. With Lagrange interpolation, smooth function u(𝜉, 𝜂) can be approximated as

2.1. One dimension problems For 1D problem in normalized domain, collocation points are equally spaced along 𝜉 axis, as follows 𝜉𝑖 = −1 +

2(𝑖 − 1) , 𝑖 = 1, 2....𝑀, 𝑀 −1

𝑢(𝜉, 𝜂) =

𝑖=1 𝑗=1

(1)

𝑢(𝜉) =

𝑀 ∏

𝑗=1 𝑘=1,𝑘≠𝑗

𝜉 − 𝜉𝑘 𝑢 , 𝜉𝑗 − 𝜉𝑘 𝑗

𝐹 (𝜉, 𝜉𝑖 )𝐺(𝜂, 𝜂𝑗 )𝑢𝑘 =

𝑄 ∑ 𝑘=1

𝜓𝑘 (𝜉, 𝜂)𝑢𝑘 ,

(6)

where

where M denotes the number of discrete points on 𝜉-axis. By Lagrange polynomial interpolation, smooth function u(𝜉)( − 1 ≤ 𝜉 ≤ 1) can be approximated by 𝑀 ∑

𝑀 ∑ 𝑁 ∑

𝐹 (𝜉, 𝜉𝑖 ) =

𝑀 𝑁 ∏ ∏ (𝜉 − 𝜉𝑚 ) (𝜂 − 𝜂𝑛 ) , 𝐺(𝜂, 𝜂𝑗 ) = , ( 𝜉 − 𝜉 ) ( 𝜂𝑗 − 𝜂𝑛 ) 𝑖 𝑚 𝑚=1 𝑛=1 𝑚≠𝑖

(7)

𝑛≠𝑗

For 2D problems, collocation nodes in normalized domain are shown in Fig. 1. uk = u(𝜉 k ,𝜂 k ) indicates the nodal value at points Pk which stands for P(𝜉 k ,𝜂 k ), where k = (j − 1) × M + i, i = 1, 2, ..., M and

(2) 37

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

Fig. 2. Seeds distribution for 2D and 3D mapping procedure: (a) seeds for 2D mapping; (b) seeds for 3D mapping.

into square blocks in normalized domain with 𝜉o𝜂 coordinate system. The coordinate transform with 8 seeds shown in Fig. 2(a) is defined as

j = 1, 2, ..., N. The number of nodes in total is Q = M × N. Shape function 𝜓 k (𝜉,𝜂) can be used in the following analysis as 𝜓𝑘 (𝜉, 𝜂) = 𝜓𝑖𝑗 (𝜉, 𝜂) = 𝐹 (𝜉, 𝜉𝑖 )𝐺(𝜂, 𝜂𝑗 ) =

𝑀 𝑁 ∏ (𝜉 − 𝜉𝑚 ) ∏ (𝜂 − 𝜂𝑛 ) . (𝜉𝑖 − 𝜉𝑚 ) 𝑛=1 (𝜂𝑗 − 𝜂𝑛 ) 𝑚=1 𝑚≠𝑖

𝑥=

(8)

𝑖=1

𝑛≠𝑗

𝜕𝜉

𝜕 𝜓𝑖𝑗 (𝜉, 𝜂) 𝜕 𝐺(𝜂, 𝜂𝑗 ) 𝜕𝐹 (𝜉, 𝜉𝑖 ) = 𝐺(𝜂, 𝜂𝑗 ), = 𝐹 (𝜉, 𝜉𝑖 ) , 𝜕𝜉 𝜕𝜂 𝜕𝜂

𝜕𝐹 = 𝜕𝜉

𝑀 ∏

𝑙=1 𝑖=1,𝑖≠𝑖,𝑖≠𝑙

(𝜉 − 𝜉𝑖 )∕

𝑀 ∏ 𝑚=1,𝑚≠𝑖

(9)

𝜕𝑢 𝜕𝑥 𝜕𝑢 𝜕𝑦 𝜕𝑢 + , = 𝜕𝜉 𝜕𝜉 𝜕𝑥 𝜕𝜉 𝜕𝑦

(10)

𝜕𝑢 𝜕𝑥 𝜕𝑢 𝜕𝑦 𝜕𝑢 + . = 𝜕𝜂 𝜕𝜂 𝜕𝑥 𝜕𝜂 𝜕𝑦

By solving equations above we can obtain ( ( ) ) 𝜕𝑢 1 𝜕𝑦 𝜕𝑢 𝜕𝑦 𝜕𝑢 𝜕𝑢 1 𝜕𝑥 𝜕𝑢 𝜕𝑥 𝜕𝑢 = , = − , − + 𝜕𝑥 𝐽 𝜕𝜂 𝜕𝜉 𝜕𝜉 𝜕𝜂 𝜕𝑦 𝐽 𝜕𝜂 𝜕𝜉 𝜕𝜉 𝜕𝜂 where | 𝜕𝑥 | 𝜕𝜉 𝐽 = || 𝜕𝑦 | 𝜕𝜉 |

where } { 𝜕𝑢(𝑃𝑄 ) T { }T 𝜕𝑢(𝑃1 ) 𝜕𝑢(𝑃2 ) 𝐔𝜈 = , 𝐮 = 𝑢(𝑃1 ), 𝑢(𝑃2 ), … , 𝑢(𝑃𝑄 ) , (12) , ,…, 𝜕𝜈 𝜕𝜈 𝜕𝜈

𝜕𝑥 | | 𝜕𝜂 | 𝜕𝑦 |, 𝛾11 𝜕𝜂 ||

=

𝜕𝑦 𝜕𝑦 𝜕𝑥 𝜕𝑥 , 𝛾 = − , 𝛾21 = − , 𝛾22 = . 𝜕𝜂 12 𝜕𝜉 𝜕𝜂 𝜕𝜉

(17)

(18)

(19)

𝐔𝑥 = 𝐁11 𝐔𝜉 + 𝐁12 𝐔𝜂 = 𝐁11 𝐃𝜉 𝐮 + 𝐁12 𝐃𝜂 𝐮 = 𝐃𝑥 𝐮, 𝐔𝑦 = 𝐁21 𝐔𝜉 + 𝐁22 𝐔𝜂 = 𝐁21 𝐃𝜉 𝐮 + 𝐁22 𝐃𝜂 𝐮 = 𝐃𝑦 𝐮,

(20)

where Dx and Dy are differential matrices in physical domain, B𝛽μ are diagonal matrices written as

(13)

(1)

where D𝜉 and D𝜂 are first order partial differential matrices in (11). Similarly, for higher order partial derivative in 3D problems with respect to 𝜉,𝜂,𝜍, we have 𝑚𝑛𝑙) 𝑛 𝑙 𝐔(𝜉𝜂𝜍 = 𝐃𝑚 𝜉 𝐃𝜂 𝐃𝜍 𝐮, 𝑚, 𝑛, 𝑙 > 0

(16)

So the first order partial derivatives are obtained as

Higher order partial derivative in 2D problems with respect to both 𝜉 and 𝜂 can be calculated by 𝜕 𝑚+𝑛 𝑢 𝑛 = 𝐃𝑚 𝜉 𝐃𝜂 𝐮, 𝑚, 𝑛 > 0 𝜕𝜉𝑚 𝜕𝜂𝑛

(15)

For smooth function u(x, y), we have

The nodal value of the first order partial derivative U𝜈 can be obtained in vector form from (6), in terms of nodal value u(Pk ), as } { 𝜕 𝜓𝑖𝑗 𝐔𝜈 = 𝐃𝜈 𝐮, 𝐃𝜈 = , (𝑖, 𝑗 = 1, 2, ..., 𝑄; 𝜈 = 𝜉, 𝜂), (11) 𝜕𝜈 𝑄×𝑄

𝐔(𝜉𝜂𝑚𝑛) =

𝑁𝑖 (𝜉, 𝜂)𝑦𝑖 ,

𝑖=1

1 (1 + 𝜉𝑖 𝜉)(1 + 𝜂𝑖 𝜂)(𝜉𝑖 𝜉 + 𝜂𝑖 𝜂 − 1), 𝑖 = 1, 2, 3, 4 4 1 𝑁𝑖 (𝜉, 𝜂) = (1 − 𝜉 2 )(1 + 𝜂𝑖 𝜂), 𝑖 = 5, 7 2 1 𝑁𝑖 (𝜉, 𝜂) = (1 − 𝜂 2 )(1 + 𝜉𝑖 𝜉), 𝑖 = 6, 8 2

(𝜉𝑖 − 𝜉𝑚 ),

𝑁 𝑁 𝑁 ∑ ∏ ∏ 𝜕𝐺 (𝜂 − 𝜂𝑗 )∕ (𝜂𝑗 − 𝜂𝑛 ), = 𝜕𝜂 𝑙=1 𝑗 =1,𝑗 ≠𝑖,𝑗 ≠𝑙 𝑛=1,𝑛≠𝑗

8 ∑

𝑁𝑖 (𝜉, 𝜂) =

where 𝑀 ∑

𝑁𝑖 (𝜉, 𝜂)𝑥𝑖 , 𝑦 =

where shape functions Ni (𝜉,𝜂) are defined as follows

So the first order partial derivative of 𝜓 ij (𝜉,𝜂) can be given as 𝜕 𝜓𝑖𝑗 (𝜉, 𝜂)

8 ∑

𝐁𝛽𝜇

(14)

⎛ 𝛾𝛽𝜇 ⎜ 𝐽 (1) [ ] ⎜ (𝑘) = 𝑑 𝑖𝑎𝑔 𝛾𝛽𝜇 ∕𝐽 = ⎜0 ⎜⋯ ⎜ ⎜ ⎝0

𝛽, 𝜇 = 1, 2; 𝑘 = 1, 2 … , 𝑄.

3. Mapping procedure

0 (2) 𝛾𝛽𝜇 𝐽 (2)

⋯

⋯

⋯ ⋯

0

⋯

⎞ ⎟ ⎟ 0 ⎟, ⋯ ⎟⎟ (𝑄) 𝛾𝛽𝜇 ⎟ ⎠ 𝐽 (𝑄) 0

(21)

where 𝛾 𝛽μ (𝛽,μ = 1, 2) are given by (19). Higher order partial derivative with respect to both x and y for 2D problems is

For 2D problems, physical domain with irregular boundary is subdivided into several quadratic blocks. Mapping technology is used to transform quadratic blocks in physical domain with xoy coordinate system

𝑚𝑛) 𝑛 𝐔(𝑥𝑦 = 𝐃𝑚 𝑥 𝐃𝑦 𝐮, 𝑚, 𝑛 > 0.

38

(22)

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

The 1D time dependent problems in space coordinates are called 2D problems in time-space approach in the following numerical examples because there are two variables time coordinate t and space coordinate x in time-spatial domain. First, for simplicity, attention is restricted to 2D non-stationary problems in multilayer materials in time-spatial domain. In a three-layered wall shown in Fig. 3, temperatures u𝛼 (x,t), 𝛼 = 1, 2, 3 for each layer obey the following governing equations 𝜕 𝑢𝛼 (𝑥, 𝑡) 𝜕 2 𝑢𝛼 (𝑥, 𝑡) − 𝜆𝛼 = 𝑓𝛼 (𝑥, 𝑡), 𝜕𝑡 𝜕 𝑥2

𝑥, 𝑡 ∈ Ω𝛼 ,

𝛼 = 1, 2, 3,

(28)

where subscript 𝛼 indicates the number of layer, u𝛼 (x,t) is temperature for each layer at time t, f𝛼 (x,t) indicates the heat source. Coefficients 𝜆𝛼 are defined as 𝜆𝛼 = 𝜅 𝛼 /𝜌𝛼 c𝛼 , where 𝜅 𝛼 ,𝜌𝛼 ,c𝛼 are thermal conductivity coefficients, density and specific heat capacity. Sub-domains are Ω1 = [0, l1 ] × [0, T], Ω2 = [l1 ,l2 ] × [0, T], Ω3 = [l2 ,l3 ] × [0, T]. In addition, boundary conditions at the left end of layer 1 are specified as 𝑢1 (0, 𝑡) = 𝑢0 (𝑡), 0 ≤ 𝑡 ≤ 𝑇 , 𝜕 𝑢1 (0, 𝑡) 𝜅1 = 𝑞0 (𝑡), 0 ≤ 𝑡 ≤ 𝑇 , 𝜕𝑥

Fig. 3. 2D time-spatial domain with blocks and the location of smelting points xs (t).

Consider 3D problems with 20 seeds and the coordinate transform 𝑥=

20 ∑ 𝑖=1

𝑁𝑖 (𝜉, 𝜂, 𝜍)𝑥𝑖 , 𝑦 =

20 ∑ 𝑖=1

𝑁𝑖 (𝜉, 𝜂, 𝜍)𝑦𝑖 , 𝑧 =

20 ∑ 𝑖=1

with initial conditions in each sub-domain 𝑁𝑖 (𝜉, 𝜂, 𝜍)𝑧𝑖 ,

(23)

𝑢1 (𝑥, 0) = 𝑈1 (𝑥), 0 ≤ 𝑥 ≤ 𝑙1 , 𝑢2 (𝑥, 0) = 𝑈2 (𝑥), 𝑙1 ≤ 𝑥 ≤ 𝑙2 ,

in which the shape functions Ni (𝜉,𝜂, 𝜍) are given by 𝑁𝑖 (𝜉, 𝜂, 𝜍) = 𝑖= 𝑁𝑖 (𝜉, 𝜂, 𝜍) = 𝑁𝑖 (𝜉, 𝜂, 𝜍) = 𝑁𝑖 (𝜉, 𝜂, 𝜍) =

1 (1 + 𝜉𝑖 𝜉)(1 + 𝜂𝑖 𝜂)(1 + 𝜍𝑖 𝜍)(𝜉𝑖 𝜉 + 𝜂𝑖 𝜂 + 𝜍𝑖 𝜍 − 2), 8 1, 2, 3, 4, 5, 6, 7, 8 1 (1 − 𝜉 2 )(1 + 𝜂𝑖 𝜂)(1 + 𝜍𝑖 𝜍), 𝑖 = 9, 11, 17, 19 4 1 (1 − 𝜂 2 )(1 + 𝜉𝑖 𝜉)(1 + 𝜍𝑖 𝜍), 𝑖 = 10, 12, 18, 20 4 1 (1 − 𝜍 2 )(1 + 𝜉𝑖 𝜉)(1 + 𝜂𝑖 𝜂), 𝑖 = 13, 14, 15, 16 4

𝑢3 (𝑥, 0) = 𝑈3 (𝑥), 𝑙2 ≤ 𝑥 ≤ 𝑙3 ,

𝑢1 (𝑙1 , 𝑡) = 𝑢2 (𝑙1 , 𝑡), 𝑢2 (𝑙2 , 𝑡) = 𝑢3 (𝑙2 , 𝑡), 𝜕 𝑢 (𝑙 , 𝑡) 𝜕 𝑢1 (𝑙1 , 𝑡) 𝜕 𝑢 (𝑙 , 𝑡) 𝜕 𝑢 (𝑙 , 𝑡) 𝜅1 = 𝜅2 2 1 , 𝜅2 2 2 = 𝜅3 3 2 , 0 ≤ 𝑡 ≤ 𝑇 . 𝜕𝑥 𝜕𝑥 𝜕𝑥 𝜕𝑥

(24)

=

𝑛 𝑙 𝐃𝑚 𝑥 𝐃𝑦 𝐃𝑧 𝐮, 𝑚, 𝑛, 𝑙

>0

[𝐃𝑡𝛼 − 𝜆𝛼 𝐃2𝑥𝛼 ]𝐮𝛼 = 𝐟 (𝛼) ,𝑥 ∈ Ω𝛼 , 𝛼 = 1, 2, 3

(26)

(32)

where 𝐟 (𝛼) = {𝑓1(𝛼) , 𝑓2(𝛼) , ..., 𝑓𝑄(𝛼) }𝑇 (Q𝛼 = M𝛼 × N𝛼 ). Q𝛼 is the number of 𝛼 collocation points for sub-domain 𝛼. As the time is treated as the second dimension (y) in time-space approach, we have Dt𝛼 = Dy𝛼 . Then the boundary conditions and initial conditions become

4. Detecting corrosion boundary 4.1. Time-space approach (TS)

𝑢(1) 1+𝑀

In metallurgy engineering, the inner part of steel-smelting furnace is always unreachable and corroded. Therefore, we need to rebuild and simulate the corrosion boundary (or surface). In the dynamic model, three layers of composite walls are assumed to simulate the structure of furnace. We assume that xs (t) is the location of corrosion points in zone Ω3 as shown in Fig. 3. This location can be predicted via exterior boundary inputs at the left end of the layer 1 in zone Ω1 . Both homogeneous and inhomogeneous inverse heat transfer problems with Cauchy conditions are considered. Further perturbation analysis of xs (t) will be studied in following numerical examples. In the first time-space approach named TS-I, a series of equidistant nodes are used in Lagrange interpolation, while the second time-space approach (TS-II) utilizes roots of Chebyshev polynomial of first kind as (2𝑘−1) 𝜉̄𝑘 = cos 𝜋,𝑘 = 1, 2, 3......𝑛, 2𝑛

(31)

Temperature u3 (xs ,t) at the corrosion points can be defined as u3 (xs ,t) = us (t), where xs is corrosion location and us (t) is material property (smelting) and set to be zero in this paper for the sake of simplicity. Applying the differential matrices of FBM to the governing equations for each layer yields

(25)

where Dx ,Dy ,Dz are differential matrices in the physical domain and higher order partial derivative with respect to x, y and z can be obtained by 𝑚𝑛𝑙) 𝐔(𝑥𝑦𝑧

(30)

where u0 (t) and q0 (t) are temperature and heat flux, U𝛼 (x)(𝛼 = 1, 2, 3) are initial temperatures at t = 0 for each sub-domain. For the nodes collocated on the interfaces between two layers, the continuity conditions read

Similar to 2D cases, we have 𝐔𝑥 = 𝐃𝑥 𝐮, 𝐔𝑦 = 𝐃𝑦 𝐮, 𝐔𝑧 = 𝐃𝑧 𝐮,

(29)

𝜅1

1 ×(𝑗−1)

𝑄1 [( ∑ 𝑘=1

= 𝑢0 (𝑡𝑗 ),

(𝑥1) 𝑑1+ 𝑀

𝑢(𝑖𝛼) = 𝑈𝛼 (𝑥𝑖 ),

)] 1 ×(𝑗−1),𝑘

𝑥 ∈ Ω𝛼 ,

= 𝑞0 (𝑡𝑗 ), 𝑢(1) 𝑘 𝛼 = 1, 2, 3;

0 ≤ 𝑡 ≤ 𝑇, 𝑖 = 1, 2, ..., 𝑀𝛼 ,

(33)

where the subscripts i and j are shown in Fig. 1. The elements in the differential matrix are 𝐃𝑥𝛼 = [𝑑𝑖𝑗(𝑥𝛼) ]𝑄 ×𝑄 and vector 𝐮𝛼 = [𝑢(𝑗𝛼) ] . In 𝛼 𝛼 addition, the continuity conditions between two layers in (31) should be considered. 4.2. Laplace transform approach (LT) Applying Laplace transform over both sides of (28) gives

(27)

where the two ends are relocated at 𝜉̄1 = −1 and 𝜉̄𝑛 = 1 respectively.

𝑠𝑢̃ 𝛼 (𝑥, 𝑠)−𝑈𝛼 (𝑥) − 𝜆𝛼 39

𝜕 2 𝑢̃ 𝛼 (𝑥, 𝑠) 𝜕 𝑥2

= 𝑓̃𝛼 (𝑥, 𝑠),

𝑥 ∈ Ω′𝛼 ,

𝛼 = 1, 2, 3,

(34)

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

where s stands for the Laplace transform parameter, Ω’1 = [0, 𝑙1 ], Ω’2 = [𝑙1 , 𝑙2 ], Ω’3 = [𝑙2 , 𝑙3 ]. The transformed variable in the Laplace space is defined by 𝑔̃ (𝑠) =

∞

∫0

𝑔(𝑡)𝑒−𝑠𝑡 𝑑𝑡.

The perturbation relative root mean-square errors PRMS(xs (ti )) is defined as √ √ 𝑁3 𝑁3 ∑[ ( ) ( ( )) √ ]2 ∑ [ ]2 𝑃 𝑅𝑀𝑆 𝑥̂ 𝑠 𝑡𝑖 = √ (43) 𝑥̂ 𝑠 𝑡𝑖 − 𝑥∗𝑠 (𝑡𝑖 ) ∕ 𝑥∗𝑠 (𝑡𝑖 ) ,

(35)

𝑖=1

where 𝑥̂ 𝑠 (𝑡𝑖 ) indicates the numerical results of the corrosion coordinate at time ti with noise level 𝛿. In following analysis, we chose M𝛼 = N𝛼 .

Introducing differential matrices over (34) yields 𝑠𝐮̃ 𝛼 − 𝐔𝛼 − 𝜆𝛼 𝐃2𝛼 𝐮̃ 𝛼 = 𝐟̃𝛼 ,

𝑥 ∈ Ω𝛼 ,

𝛼 = 1, 2, 3,

(36)

Example 5.1. An inhomogeneous heat problem across multilayer furnace

with boundary conditions at x = 0

First, we consider an inhomogeneous inverse heat conduction problem across multilayer furnace with three-component walls shown in Fig. 3. Corrosion occurs in Ω3 and corrosion location xs (t) is to be determined by solving equation u3 (x,t) = 0. Parameters in (31) and (32) are selected as 𝜅 1 = 2, 𝜅 2 = 3, 𝜅 3 = 1, 𝜌𝛼 c𝛼 = 1, 𝛼 = 1, 2, 3 and heat sources for each layer are specified as f1 (x,t) = −12x + 1, f2 (x,t) = −12x + 4/3, f3 (x,t) = −12x. In addition, we set T = 1, l3 = 1 in this example. It is easy to verify that the analytical solution of (28) is

𝑢̃ (1) (𝑠) = 𝑢̃ 0 (𝑠), 1 𝜅1

𝑀 ∑ [(𝑑1(𝑥𝑗 1) )]𝑢̃ (1) 𝑗 = 𝑞̃0 (𝑠).

(37)

𝑗=1

The continuity conditions at interfaces in the Laplace domain become 𝑢̃ (1) (𝑙1 , 𝑠) = 𝑢̃ (2) (𝑙1 , 𝑠), 𝑢̃ (2) (𝑙2 , 𝑠) = 𝑢̃ (3) (𝑙2 , 𝑠), 𝜅1

𝜕 𝑢̃ (1) (𝑙1 , 𝑠) 𝜕 𝑢̃ (2) (𝑙1 , 𝑠) 𝜕 𝑢̃ (2) (𝑙2 , 𝑠) 𝜕 𝑢̃ (3) (𝑙2 , 𝑠) = 𝜅2 , 𝜅2 = 𝜅3 . 𝜕𝑥 𝜕𝑥 𝜕𝑥 𝜕𝑥

( ) ( ) ( ) 1 5 2 103 5 81 52853 𝑥+ 𝑥+ +2 𝑡− + + , 4 4 120 4 32 72000 ( )2 ( ) ( ) 2 1 5 11 5 81 121 𝑥+ 𝑥+ +2 𝑡− + 𝑢2 (𝑥, 𝑡) = 𝑥3 + + , 3 9 4 15 4 32 100 ( )2 ( ) 5 81 32 𝑢3 (𝑥, 𝑡) = 2𝑥3 + 𝑥 + +2 𝑡− . − 4 32 375 𝑢1 (𝑥, 𝑡) = 𝑥3 +

(38)

A powerful method proposed by Durbin [31] is adopted for the inverse procedure of Laplace transformation. By picking K + 1 samples in the Laplace space sk ,k = 0, 1, ..., K, the transformed 𝑔̃ (𝑠𝑘 ) variables are obtained following the procedure in the chapter above. Then time dependent function g(t) can be obtained from ) ( 𝐾 ] 2𝑒𝜎𝑡∕𝑇 1 ( ) ∑ [ 𝑔(𝑡) = − 𝑔̃ 𝑠0 + Re 𝑔̃ (𝑠𝑘 )𝑒2𝑘𝜋𝑡𝑖∕𝑇 , (39) 𝑇0 2 𝑘=0

𝑢3 (𝑥, 𝑡) =

𝑁3

∑ ∏ 𝑗=1 𝑘=1,𝑘≠𝑗

𝑥 − 𝑥𝑘 𝑢 (𝑡). 𝑥𝑗 − 𝑥𝑘 𝑗

𝑢0 (𝑡) = 2𝑡 − 𝑞0 (𝑡) =

1 2 89 3223 𝑥 + 𝑥− , 𝑥 ∈ Ω1 , 4 60 1125 2 1 91 1243 𝑈2 (𝑥) = 𝑥3 + 𝑥2 + 𝑥− , 𝑥 ∈ Ω2 , 3 9 90 450 5 2689 𝑈3 (𝑥) = 2𝑥3 + 𝑥2 + 𝑥 − , 𝑥 ∈ Ω3 . 2 750

(45)

(46)

Two time-space approaches (TS-I, TS-II) and Laplace transform approach (LT) with uniformly distributed nodes are implemented respectively. In FBM (LT), the free parameters in Durbin’s inverse method are selected as 𝜎 = 5, T0 = 6. Initial conditions in Ω𝛼 ,𝛼 = 1, 2, 3 are given in (46) and boundary conditions at x = 0 can be written, from (45) in transformed domain, as

5. Numerical examples

2 3223 − , 𝑠2 1125𝑠 𝜕 𝑢 (0, 𝑠) 89 𝑞̃0 (𝑠) = 𝜅1 1 = . 𝜕𝑥 30𝑠

𝑢̃ 0 (𝑠) = 𝑢1 (0, 𝑠) =

In this section, the computation of the moving boundary with corrosion points xs (t) in the sub-domain Ω3 is fulfilled. Three layers with l1 = 0.2 and l2 = 0.4 are shown in Fig. 3. Sensitivity has to be investigated due to the ill-posed characteristics of inverse problems. A random noise rand(i) is introduced at the left end of layer 1 to test the stability of the FBM by

(47)

In the Point Collocation Method (PCM) with radial basis functions, smooth function u(x, t) can be approximated as follows

𝑢̂ 0 (𝑡) = 𝑢0 (𝑡)[1 + 𝛿 ⋅ 𝑟𝑎𝑛𝑑(𝑖)],

𝑢(𝑥, 𝑡) =

𝑄𝛼 ∑ 𝑖=1

(41)

𝑎𝑖 𝜙𝑖 (𝑟), 𝛼 = 1, 2, 3,

(48)

where ai are coefficients and the radial basis function is chosen as √ inverse multi-quadrics, i.e. 𝜙𝑖 (𝑟) = 1∕ 𝑟2𝑖 + 𝜔2 . Optimal selection of parameter 𝜔 can be found in [17,32]. The relative errors RMS (xs ) for FBM (TS-I), FBM (TS-II) and FBM (LT) are shown in Table 1. RMS (xs ) is calculated after all corrosion points xs (t) are found in the corrosion procedure, which occurs from 0 to T. It is obvious that the degree of accuracy by FBM (TS-I)

where u0 andq0 are defined in (29), rand(i) ∈ [0, 1] is random variable and 𝛿 is the level of noise. The relative root mean-square errors presenting the average error between numerical solution (xs ) and analytic solution (𝑥∗𝑠 ) are defined as

𝑖=1

89 , 0 ≤ 𝑡 ≤ 𝑇, 30

𝑈1 (𝑥) = 𝑥3 +

By using the Newton–Downhill method, we can obtain the roots (xs (t)) of equation u3 (x,t) = 0 numerically for each time step.

√ √ 𝑁3 𝑁3 √∑ [ ]2 ∑ [ ]2 𝑅𝑀𝑆(𝑥𝑠 ) = √ 𝑥𝑠 (𝑡𝑖 ) − 𝑥∗𝑠 (𝑡𝑖 ) ∕ 𝑥∗𝑠 (𝑡𝑖 ) .

3223 , 0 ≤ 𝑡 ≤ 𝑇, 1125

with initial conditions at t = 0.

(40)

𝑞̂0 (𝑡) = 𝑞0 (𝑡)[1 + 𝛿 ⋅ 𝑟𝑎𝑛𝑑(𝑖)],

(44)

The exact location xs (t) can be calculated by solving u3 (x,t) = 0 in (44) and the initial corrosion location is 𝑥𝑠 (0) = 1599 . The Dirichlet as 2024 well as Neumann boundary conditions at x = 0 in Ω1 are

where the Laplace transform parameters are set as: 𝑠𝑘 = √ (𝜎 + 2𝑘𝜋𝑖)∕𝑇0 , (𝑖 = −1) in which 𝜎 and T0 are two free parameters. The choosing of parameters T0 depends on the observing period in time domain. In numerical examples, all variables are normalized for analysis convenience. In order to obtain the location of smelting points xs (t) as shown in Fig. 3 at time t, where the temperature u3 [xs ,t] = 0 , the Lagrange polynomial can be used to interpolate the temperature u3 [x,t] = 0 in the domain Ω3 as follows 𝑁3

𝑖=1

(42)

𝑖=1

40

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

Table 1 RMS (xs ) for different M and different methods. M

FBM(TS-I)

FBM(LT)

FBM(TS-II)

6 7 8 9 10 16 22 28

2.3914 × 10-9 5.8233 × 10-8 6.5488 × 10-9 7.7152 × 10-9 6.7089 × 10-9 4.8959 × 10-7 0.8231 0.2567

4.2468 × 10-3 3.8877 × 10-3 7.1411 × 10-3 3.0017 × 10-3 2.7385 × 10-3 5.6828 × 10-2 3.2341 × 10-2 2.4594 × 10-2

2.2333 × 10-8 3.5227 × 10-9 6.7788 × 10-9 2.0761 × 10-9 5.7843 × 10-9 1.8289 × 10-7 8.1548 × 10-4 9.5428 × 10-4

Table 2 PRMS (xs ) with different levels of noise 𝛿. 𝛿\

0.01 0.05 0.1

M = 10

PCM

FBM(TS-I)

FBM(TS-II)

1.5207 × 10-4 1.5074 × 10-4 1.0622 × 10-4

1.6193 × 10-4 7.9764 × 10-4 1.5675 × 10-3

0.3914 × 10-2 0.3475 × 10-2 0.6170 × 10-2

Fig. 5. Relative errors RMS (xs ) for FBM(TS-I), FBM(TS-II) and FBM(LT) with different M. Table 3 RMS (xs ) with different M by FBM (TS-I and TS-II). M 8 13 18 23 28 33

FBM(TS-I)

FBM(TS-II) -4

9.7268 × 10 2.0050 × 10-5 5.3579 × 10-4 4.5832 × 10-2 0.5668 0.2587

8.4624 × 10-4 6.4900 × 10-5 1.2560 × 10-5 3.9025 × 10-5 3.0073 × 10 − 5 8.0353 × 10-4

Consider a 2D homogeneous inverse heat transfer problem though a multilayer furnace with moving boundary in time-spatial domain. In this case, both FEM and FBM (TS) are employed with different level of noise. In FBM (TS), the maximum length and observing time are assumed as l3 = 1.0, T = 0.5 Material properties are assumed as 𝜅 1 = 6, 𝜅 2 = 3, 𝜅 3 = 2, 𝜌1 c1 = 2/3, 𝜌2 c2 = 3/4, 𝜌3 c3 = 2. Zero heat sources are assumed with following initial conditions Fig. 4. Comparisons between analytic and numerical solutions for corrosion location when M = 14, 𝛿 = 0.1.

1 4 𝑈 1 ( 𝑥 ) = − 𝑥 + , 𝑥 ∈ Ω1 , 3 5 2 13 𝑈2 (𝑥) = − 𝑥 + , 𝑥 ∈ Ω2 , 3 15 𝑈3 (𝑥) = 1 − 𝑥, 𝑥 ∈ Ω3 .

(about 5.8233 × 10-8 ) is much higher than that by FBM (LT) (nearly 3.0017 × 10-3 ). Besides, FBM (TS-I) and FBM (TS-II) nearly stand at the same level when M < 16. However, the numerical results by FBM (TS-I) become divergent with the increasing of collocation points, whilst the solutions by both FBM (TS-II) and FBM (LT) are still stable. Table 2 shows the relative errors of PRMS (xs ) versus different levels of noise 𝛿 by FBM (TS-I), FBM (TS-II) and PCM respectively. In the case of M = 10, PRMS (xs ) by FBM (TS-I) are maintained at the level of 1.5074 × 10-4 , while FBM (TS-II) are at the level of 1.5675 × 10-3 , which means that TS-I and TS-II are of higher degree of accuracy than PCM (about 0.3475 × 10-2 ) for different levels of noise. When 0.01 < 𝛿 < 0.05, the relative errors by FBM (TS-I) and FBM (TS-II) are less than 0.002 which is less than input noise. It indicates the noises 𝛿 from the measurements have limited effect on numerical results. Comparisons with analytic solutions are presented in Figs. 4 and 5 which shows the relative errors versus different M for FBM (TS-I), FBM (TS-II) and FBM (LT) respectively. It can be concluded that FBM (TS-II) is the best among them and FBM (TSI) performs better than FBM (LT). Due to the numerical solutions match the exact solutions very well even when 𝛿 = 0.1, we can conclude that the numerical methods proposed in this paper are stable and accurate.

(49)

Neumann boundary is given by 𝜕 u1 (0,t)/𝜕 x = −1/3, 0 < t < T while the Dirichlet boundary condition is unknown at x = 0. However, we can determine it from the corresponding direct problem where the boundary conditions at x = 0 andx = xs (t) can be given as ( ) 𝜕 𝑢1 (0, 𝑡)∕𝜕𝑥 = −1∕3, 𝑢3 𝑥𝑠 (𝑡), 𝑡 = 0, 0 ≤ 𝑡 ≤ 𝑇 . (50) Moving boundary is supposed as xs (t) = 2.4t2 − 1.2t + 1.0 in the direct problem. So we can get the temperature distribution at x = 0, which can be employed in former inverse procedure to reconstruct and testify the moving boundary. Three blocks I, II and III shown in Fig. 3 are used in this direct problem in time-spatial domain. For more complicated moving boundary, domain Ω3 needs to be subdivided into two or more blocks. In addition, the continuity conditions in terms of heat flux as well as temperature should be considered. RMS (xs ) versus different M by FBM (TS-I) and FBM (TS-II) are shown in Table 3. It is clear that RMS (xs ) by FBM (TS-I) is of small oscillation around the level of 5.3579 × 10-4 and FBM (TS-II) almost stands at 1.2560 × 10-5 when 8 ≤ M ≤ 18. When 23 ≤ M ≤ 33, FBM (TSII) is more accurate than FBM (TS-I). Numerical solutions of PRMS (xs )

Example 5.2. A homogeneous heat problem across multilayer furnace 41

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

Table 4 PRMS (xs ) with different 𝛿 when M = 10.

Table 5 RMS (xs ) for different M and Q.

𝛿

FBM(TS-I)

FEM

FBM(TS-II)

M

Q

FBM(TS-I)

FBM(TS-II)

0 0.01 0.03 0.05 0.07 0.1

8.6458 × 10-5 9.4384 × 10-4 2.6960 × 10-3 4.3778 × 10-3 6.0006 × 10-3 8.3423 × 10-3

3.0619 × 10-4 2.6343 × 10-3 8.1065 × 10-3 1.4463 × 10-2 2.2212 × 10-2 3.8096 × 10-2

1.9311 × 10-4 1.3111 × 10-3 3.5661 × 10-3 5.6736 × 10-3 7.6609 × 10-3 1.0466 × 10-2

7 8 9 10 11 15 18

343 512 729 1000 1331 3375 5832

4.4944 × 10-9 2.0064 × 10-8 5.1311 × 10-8 5.1509 × 10-5 3.1872 × 10-5 2.1328 × 10-3 1.9623

8.2311 × 10-9 6.8738 × 10-10 3.8424 × 10-9 5.7316 × 10-7 1.3778 × 10-7 1.2018 × 10-6 6.5527 × 10-5

Example 5.3. Three-dimensional homogeneous heat problem across monolayer plate For simplicity, we only consider a 3D homogeneous inverse heat transfer problem in monolayer plate in time-spatial domain as shown in Fig. 8(a). The aim of this example is to determine the corrosion surface P[xs (t),ys (t)]. In FBM (TS), we only consider the outermost layer Ω = [0, 1] × [0, 1] × [0, T]. Where the maximum time is T = 3/16. The governing equation is given by [ 2 ] 𝜕𝑢(𝑥, 𝑦, 𝑡) 𝜕 𝑢(𝑥, 𝑦, 𝑡) 𝜕 2 𝑢(𝑥, 𝑦, 𝑡) −𝜆 + = 𝑓 (𝑥, 𝑦, 𝑡), 𝑥, 𝑦, 𝑡 ∈ Ω (51) 𝜕𝑡 𝜕 𝑥2 𝜕 𝑦2 where 𝜆 = 2, f(x, y, t) = −9. Analytical solution is assumed as 𝑢(𝑥, 𝑦, 𝑡) = (𝑥 − 1)2 + 𝑦2 − 𝑡 −

1 . 16

(52)

Similar with Example 5.2, Dirichlet boundary condition 𝑢0 (0, 𝑦, 𝑡) at x = 0 is unknown. So a direct problem is also needed here to determine the temperature distribution at x = 0 first. Two blocks I and II shown in Fig. 8(b) are used in this direct problem, boundary conditions are given as

Fig. 6. Comparisons between analytic solutions and numerical solutions of corrosion location xs (t) when M = 10, 𝛿 = 0.0.

𝑞0 (0, 𝑦, 𝑡) = 4, 𝑞0 (1, 𝑦, 𝑡) = 0, 𝑞0 (𝑥, 0, 𝑡) = 0, 𝑞0 (𝑥, 1, 𝑡) = 4, 𝑢(𝑥𝑠 , 𝑦𝑠 , 𝑡) = 0,

𝑡 ∈ [0, 𝑇 ]

(53) 1 16

in which (𝑥𝑠 − 1)2 + 𝑦2𝑠 − 𝑡 − = 0 and initial conditions are given from the analytical solution (52) as 𝑈 (𝑥, 𝑦, 0) = (𝑥 − 1)2 + 𝑦2 −

1 . 16

(54)

𝑢0 (0, 𝑦, 𝑡) can be obtained numerically by solving the direct problem. Thereafter, in the inverse problem, boundary conditions are amended by 𝑢(0, 𝑦, 𝑡) = 𝑢0 , 𝑞0 (0, 𝑦, 𝑡) = 4, 𝑞0 (1, 𝑦, 𝑡) = 0, 𝑞0 (𝑥, 0, 𝑡) = 0, 𝑞0 (𝑥, 1, 𝑡) = 4,

𝑡 ∈ [0, 𝑇 ]

(55)

initial conditions are given in (54). The unknown moving surface P[xs (t),ys (t)] can be identified with application of FBM. Random noise rand(i) ∈ [0, 1] is added on Dirichlet data and Neumann data at x = 0 as follows 𝑢̂ (0, 𝑦, 𝑡) = 𝑢0 (0, 𝑦, 𝑡)(1 + 𝛿 ⋅ 𝑟𝑎𝑛𝑑(𝑖)),

Fig. 7. Comparisons between analytic solutions and numerical solutions of corrosion location xs (t) when M = 10, 𝛿 = 0.05.

𝑞̂(0, 𝑦, 𝑡) = 𝑞0 (0, 𝑦, 𝑡)(1 + 𝛿 ⋅ 𝑟𝑎𝑛𝑑(𝑖)).

by FBM (TS-I), FBM (TS-II) and FEM are shown in Table 4. In order to study the efficiency, 10 quadratic elements for FEM and 10 collocation points for FBM (TS-I and TS-II) are used. It is observed that PRMS (xs ) by FBM (TS-I) is better than FEM and FBM (TS-II) when 𝛿 = 0.1. It is also clear that both FBM (TS-I) and FBM (TS-II) are stable with different 𝛿. Comparisons between analytic solutions and numerical solutions of corrosion location without noise and with noise 𝛿 = 0.05 when M = 10 are demonstrated in Figs. 6 and 7 respectively.

Similarly 𝛿 ∈ [0, 1] is the level of noise. Numerical results are presented in Tables 5, 6 and Fig. 9. In order to observe the effect of node density, the relative error RMS (xs ) with different M by FBM (TS-I) and FBM (TS-II) without noise are shown in Fig. 9. It indicates that FBM (TS-II) is much better than FBM (TS-I). Table 5 shows the relative error RMS (xs ) of corrosion location for different numbers M and Q by FBM (TS-I) (uniform) and FBM (TS-II) (Chebyshev). It can be found that the error of FBM (TS-II) is less than FBM (TS-I), especially while M is bigger than 15. When M = 18, for FBM (TS-II) with roots 42

(56)

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

Fig. 8. (a) 3D time-spatial domain (b) direct problem domain.

Table 6 PRMS (xs ) for different levels of noise whenM = 8.

6. Conclusion

𝛿

FBM(TS-I)

FBM(TS-II)

0.01 0.03 0.05

4.6592 × 10-2 4.3119 × 10-2 6.5802 × 10-2

3.7503 × 10-2 5.6833 × 10-2 8.1237 × 10-2

An inverse computational procedure to determine the inner corrosion location of heat conduction composite walls from Cauchy measurement data of temperature and heat flux on the exterior boundary is proposed in the paper. FBM is employed to simulate the moving boundary of the inner surface of multilayer composite materials. First, transient heat transform for both 1D inverse problems in Laplace space in space coordinates and 2D/3D inverse problems in time-spatial domain are investigated. Then Lagrange interpolation and Newton–Downhill method are applied to obtain the location of the corrosion points. Comparisons between FBM, PCM (RBF) and FEM are given to study the stability, convergence and efficiency of FBM. Compared with Radial Basis Functions (RBF), there are no free parameters in FBM (TS-I, TS-II) and FBM (LT) due to the use of Lagrange interpolation. Roots of Chebyshev polynomial of first kind are utilized in the FBM for the first time, leading to a marked improvement in the convergence degree. Three different working cases with various corrosion boundaries (surfaces) are investigated. All the results indicate that FBM is well preformed to detect the corrosion location with its unique efficiency and stability. FBM (TS-II) is of much higher degree of convergence over FBM (TS-I) for large number of nodes. Acknowledgments The work of this paper was partially supported by a grant from the National Youth Science Foundation of China (grant no. 11401423) and National Natural Science Foundation of China (grant no. 11472184). Supplementary materials Supplementary material associated with this article can be found, in the online version, at doi:10.1016/j.enganabound.2018.01.009.

Fig. 9. RMS (xs ) for different M without noise.

References [1] Aparicio ND, Atkinson C. An inverse boundary problem for the steady-diffusion equation with moving boundaries with applications to the pearlite-austenite transformation in steel. Inverse Prob 2002;18(4):1013–23. [2] Bryan K, Caudill L. Reconstruction of an unknown boundary portion from Cauchy data in n dimensions. Inverse Prob 2005;21(1):239–55. [3] Bryan K, Caudill LF Jr. Stability and reconstruction for an inverse problem for the heat equation. Inverse Prob 1998;14(6):1429–53. [4] Fredman TP. A boundary identification method for an inverse heat conduction problem with an application in ironmaking. Heat Mass Transfer 2004;41:95–103.

of Chebyshev polynomial of first kind, RMS (xs ) is obtained around 6.5527 × 10-5 while FBM (TS-I) is diverged. Table 6 shows the PRMS (xs ) with different levels of noise 𝛿 when M = 8. It can be noticed that the PRMS (xs ) for FBM (TS-I) and FBM (TS-II) almost have very little oscillation on the level of 5.6833 × 10-2 . So it can be concluded that the method proposed is stable against different levels of noise in this case. 43

M. Lei et al.

Engineering Analysis with Boundary Elements 89 (2018) 36–44

[5] Liu GR, Zhang GY. A normed G space and weakened weak(W-2)formulation of a cell-based smoothed point interpolation method. Int J Comput Methods 2009;6(1):147–79. [6] Raynaud M, Bransier J. A new finite difference method for the nonlinear inverse heat conduction problem. Numer Heat Transfer 1986;9:27–42. [7] Aliabadi MH, Wen PH. Boundary element methods in engineering and sciences. London: Imperial College Press; 2010. [8] Powell MJD. Radial basis functions for multivariable interpolation: a review. Algorithms Approximation 1987:143–67. [9] Jackson IRH. Convergence properties of radial basis functions. Construct Approximation 1988;4(1):243–64. [10] Farwig R. Multivariate interpolation of arbitrarily spaced data by moving least squares methods. J Comput Appl Math 1986;16(1):79–93. [11] Mathon R, Johnston RL. The approximate solution of elliptic boundary-value problems by fundamental solutions. Siam J Numer Anal 1977;14(14):638–50. [12] Golberg MA, Chen CS. The theory of radial basis functions applied to the BEM for inhomogeneous partial differential equations. Boundary Elem Commun 1994;5:57–61. [13] Hon YC, Wei T. The method of fundamental solution for solving multidimensional inverse heat conduction problems. Comput Model Eng Sci 2005;7(2):119–32. [14] Wei T, Yamamoto M. Reconstruction of a moving boundary from Cauchy data in one dimensional heat equation. Inverse Prob Sci Eng 2009;17(4):551–67. [15] Li YS, Wei T. Identification of a moving boundary for a heat conduction problem in a multilayer medium. Heat Mass Transfer 2010;46:779–89. [16] Hu XY, Chen WB. Determination of unknown boundary in the composite materials with Stefan-Boltzmann conditions. Chin Ann Math 2010;31B(2):145–62. [17] Niu RP, Liu GR, Li M. Reconstruction of moving boundary of a multilayer composite wall as a heat conduction media. Eng Anal Boundary Elem 2014;42:92–8. [18] Niu RP, Liu GR, Li M. Inverse analysis of heat transfer across a multilayer composite wall with Cauchy boundary conditions. Int J Heat Mass Transfer 2014;79:727–35. [19] Bellman RE, Casti J. Differential quadrature and long-term integration. J Math Anal Appl 1971;34:235–8.

[20] Bert CW, Wang X, Striz AG. Differential quadrature for static and free vibration analysis of anisotropic plate. Int J Solids Struct 1993;30:1737–44. [21] Bert CW, Wang X, Striz AG. Convergence of the DQ method in the analysis of anisotropic plates. J Sound Vib 1994;170:140–4. [22] Liew KM, Han JB, Xiao ZM, Du H. Differential quadrature method for Mindlin plates on Winkler foundations. Int J Mech Sci 1996;38:405–21. [23] Liew KM, Han J-B, Xiao ZM. Differential quadrature method for thick symmetric cross-ply laminates with first-order shear flexibility. Int J Solids Struct 1996;33:2647–58. [24] Liew KM, Han J-B, Xiao ZM. Vibration analysis of circular Mindlin plates using the differential quadrature method. J Sound Vib 1997;205:603–17. [25] Liew KM, Han J-B. Bending analysis of simply supported shear deformable skew plates. J Eng Mech ASCE 1997;123:214–21. [26] Liew KM, Teo TM, Han JB. Comparative accuracy of DQ and HDQ methods for three-dimensional vibration analysis of rectangular plates. Int J Numer Methods Eng 2014;45:1831–48. [27] Bert CW, Malik M. Differential quadrature method in computational mechanics: a review. Appl Mech Rev 1996;49(1):1–28. [28] Wen PH, Cao P, Korakianitis T. Finite block method in elasticity. Eng Anal Boundary Elem 2014;46:116–25. [29] Li M, Wen PH. Finite block method for transient heat conduction analysis in functionally graded media. Int J Numer Methods Eng 2014;99(5):372–90. [30] Li M, Lei M, Munjiza A, Wen PH. Finite block method for contact analysis in functionally graded materials. Int J Numer Methods Eng 2015;5(3):668–88. [31] Durbin F. Numerical inversion of Laplace transforms: an efficient improvement to Dubner and Abate’s method. Comput J 1975;17:371–6. [32] Niu RP, Liu GR, Li M. A technique with a correction term for unsteady state heat transfer problem with moving boundary. Int J Comput Methods 2016;13(02).

44