# Solution to Two-Dimensional Steady Inverse Heat Transfer Problems with Interior Heat Source Based on the Conjugate Gradient Method

Authors: Shoubin Wang; Li Zhang; Xiaogang Sun; Huangchao Jia
Journal: Mathematical Problems in Engineering
 (2017)

---

## Abstract

The compound variable inverse problem which comprises boundary temperature distribution and surface convective heat conduction coefficient of two-dimensional steady heat transfer system with inner heat source is studied in this paper applying the conjugate gradient method. The introduction of complex variable to solve the gradient matrix of the objective function obtains more precise inversion results. This paper applies boundary element method to solve the temperature calculation of discrete points in forward problems. The factors of measuring error and the number of measuring points zero error which impact the measurement result are discussed and compared with L-MM method in inverse problems. Instance calculation and analysis prove that the method applied in this paper still has good effectiveness and accuracy even if measurement error exists and the boundary measurement points’ number is reduced. The comparison indicates that the influence of error on the inversion solution can be minimized effectively using this method.


---

## Body

## 1. Introduction

Inverse heat transfer problem (Inverse Heat Transfer Problems, IHTP) means inversion unknown characteristic parameters of heat conduction objects using the internal or surface local heat measurement [1, 2], such as thermal physical parameters, thermal boundary conditions, geometric boundary shape, and heat conduction coefficient. Inverse heat conduction problem has a broad application prospect in nondestructive testing, geometry optimization, aerospace engineering, power engineering, mechanical engineering, construction engineering, biological engineering, metallurgical engineering, materials processing, biomedical and food engineering, and other fields [1–12]. Focusing on this problem, domestic and foreign scholars have done a lot of research. Zhu et al. studied the disadvantages of optimization algorithms and inherent space distribution characteristics of the measurement in inverse heat conduction problems. They proposed a decentralized fuzzy reasoning mechanism and established a decentralized fuzzy inference system of two-dimensional steady inverse heat conduction problems [3, 4]. Cui et al. gave an improved conjugate gradient method. They introduced the complex variable derivation into the traditional conjugate gradient method, calculated the sensitivity coefficients accurately, and identified the boundary conditions [5]. Yu combined the boundary element method with complex variables derivation to inverse the inhomogeneous material coefficients of thermal conductivity [6]. Wang et al. applied particle swarm optimization algorithm and the least square method for solving the inverse problems. That improved the precision and shorted the solving time at the same time [13]. Yaparova studied the heat conduction boundary value inverse problem by solving the stability boundary based on Laplace and Fourier transformation [14]. Tian combined SPSO algorithm with the conjugate gradient method and the fast convergence of the traditional regularized gradient algorithm and global convergence of the intelligent optimization algorithm were highlighted [15]. Using the boundary element method, Zhou et al. analyzed two-dimensional transient conduction problems, introduced the conjugate gradient method to find the heat conduction coefficient, and verified the validity and stability of the method [16]. A stable differential method was proposed to solve the inverse heat conduction problems by Baranov et al. on the basis of the differential transform [17]. In this paper, the two-dimensional steady forward problems with internal heat source are solved using the boundary element method. The existing surface convective thermal dissipation is analyzed, and the conjugate gradient method is applied to solve the inverse problems. And the complex variable derivation method is introduced to increase the sensitivity coefficient and the accuracy in the inversion.


## 2. Thermal Steady State Forward Problems

### 2.1. The Boundary Integral Equation

The mathematical model of two-dimensional steady heat conduction problems with internal heat source is(1)∂2T∂x2+∂2T∂y2=f∈Ω,T=T-∈Γ1,∂T∂n=q-∈Γ2,∂T∂n+hλT=q′′-∈Γ3,where Γ1, Γ2, and Γ3 are the boundary of the domain Ω, Γ=Γ1+Γ2+Γ3, f=-s/λ, s is the internal heat source, h is the heat conduction coefficient between solid and fluid, λ=HTC, T= temp, and n is the normal vector outside the boundary. q′′-=h/λTf, Tf is environment temperature, and q- is heat flux.Weight functionT∗ is introduced into weighted residual of control equation. (2)∫Ω∇2T-fT∗dΩ=∫Γ3∂T∂n+hλT-q′′-T∗dΓ+∫Γ2∂T∂n-q-T∗dΓ-∫Γ1T-T-∂T∗∂ndΓ.Green’s function of Laplace equation is∬Dv∇2u-u∇2vdΩ=∫cv∂u/∂n-u∂v/∂nds, where C is the boundary curve of plane closed area D and ds is arc differential [18].To decompose the left side of the equation, it becomes(3)∫Ω∇2T-fT∗dΩ=∫Ω∇2TT∗dΩ-∫ΩfT∗dΩ.According to Laplace Green’s function, we can have(4)∫Ω∇2TT∗dΩ=∫ΩT∗∇2T-T∇2T∗dΩ+∫ΩT∇2T∗dΩ.Further reduction gives(5)∫Ω∇2TT∗dΩ=∫ΓT∗∂T∂n-T∂T∗∂ndΓ+∫ΩT∇2T∗dΩ.Substituting (5) for the same term in (3), (2) becomes(6)∫ΓT∗∂T∂n-T∂T∗∂ndΓ+∫ΩT∇2T∗dΩ-∫ΩfT∗dΩ=∫Γ3∂T∂n+hλT-q′′-T∗dΓ+∫Γ2∂T∂n-q-T∗dΓ-∫Γ1T-T-∂T∗∂ndΓ.Further reduction gives(7)∫ΩT∇2T∗dΩ=∫Γ3∂T∂n+hλT-q′′-T∗dΓ+∫Γ2∂T∂n-q-T∗dΓ-∫Γ1T-T-∂T∗∂ndΓ-∫ΓT∗∂T∂n-T∂T∗∂ndΓ+∫ΩfT∗dΩ.T
∗ is the basic solution of Laplace equation and it can make ∇T∗+δr-ri=0, where the point “i” refers to the unit point concentration. Because ∂T∗/∂n=q∗ and ∂T/∂n=q, (7) can be further reduced.The left side of (7) becomes(8)∫ΩT∇2T∗dΩ=∫ΩT∇2T∗+δr-ri-Tδr-ridΩ.Calculating furthermore, it is shown as follows:(9)∫ΩT∇2T∗+δr-ri-Tδr-ridΩ=0-∫ΩTδr-ridΩ.For the selectivity ofδ function, ∫-∞+∞fxδx-x0dr=f(x0). Then the reduced result on the left side of (7) is(10)-∫ΩTδr-ridΩ=-Ti.Substituting (10) into (7),(11)-Ti=∫Γ3q+hλT-q′′-T∗dΓ+∫Γ2q-q-T∗dΓ-∫Γ1T-T-q∗dΓ-∫ΓT∗q-Tq∗dΓ+∫ΩfT∗dΩ.BecauseΓ=Γ1+Γ2+Γ3, ∫ΓT∗q-Tq∗dΓ is decomposed: (12)∫ΓT∗q-Tq∗dΓ=∫Γ1T∗q-Tq∗dΓ+∫Γ2T∗q-Tq∗dΓ+∫Γ3T∗q-Tq∗dΓ.Substituting (12) for the same term in (11), (13)-Ti=∫Γ3hλTT∗-q′′-T∗+Tq∗dΓ+∫Γ2Tq∗-q-T∗dΓ+∫Γ1T-q∗-T∗qdΓ+∫ΩfT∗dΩ.Further reduction gives(14)-Ti=∫Γ3hλTT∗dΓ+∫Γ3Tq∗dΓ+∫Γ2Tq∗dΓ+∫Γ1T-q∗dΓ-∫Γ3q′′-T∗dΓ+∫Γ2q-T∗dΓ+∫Γ1qT∗dΓ+∫ΩfT∗dΩ.The second and the third terms are combined, which makes(15)Ti+∫Γq∗TdΓ+∫Γ3hλTT∗dΓ=∫ΓT∗qdΓ-∫ΩfT∗dΩ.The two-dimensional steadyT∗=1/2πln⁡1/r, r=x-xi2-y-yi2.Moving the source pointi to the boundary, the integral equation at any point on the boundary can be obtained. Thus, the basic solution T∗ and the singularity of the integral ∂T∗/∂n need to be considered. Assuming that a circular arc centering on the point i near the border and ε is the small radius of the circle, (15) can be written as(16)Ti+∫Γ-εq∗TdΓ+∫Γεq∗TdΓ+∫Γ3hλTT∗dΓ=∫Γ-εT∗qdΓ+∫ΓεT∗qdΓ-∫ΩfT∗dΩ,where Γε is the new circular border and Γ-ε is the boundary outside the new circular boundary in Γ. Considering the smooth boundary, the integral mean value theorem is introduced. (17)∫ΓεTq∗dΓ=∫ΓεT∂∂n14πε2dΓ=∫ΓεT-14πε2dΓ=Tξ-14πε22πε2,where ξ is the point on Γε. Let ε→0T(ξ)→Ti; then(18)limε→0⁡∫ΓεTq∗dΓ=-12Ti,limε-0⁡∫Γ-εq∗TdΓ=∫Γq∗TdΓ,limε-0⁡∫Γ-εqT∗dΓ=∫Γq∗TdΓ.Synthesizing the above limit results, we have(19)12Ti+∫Γq∗TdΓ+∫Γ3hλTT∗dΓ=∫ΓT∗qdΓ-∫ΩfT∗dΩ.Considering (15) and (19) simultaneously, boundary integral equation of any point in the domain Ω or on the boundary can be found out:(20)CiTi+∫Γq∗TdΓ+∫Γ3hλTT∗dΓ=∫ΓT∗qdΓ-∫ΩfT∗dΩ,Ci=1∈Ω,12∈Γ.


### 2.2. Calculation of the Unknown Value Based on Boundary Element Method

Discreteness of boundary integral equations: because the quantitiesT and q on the boundary are known, the only thing to do is divide the boundary discretely. Boundary Γ is divided into N boundary units Γj(j=1,…,N), where there are N1 units on the boundary Γ1, N2 units on the boundary Γ2, and N3 units on the boundary Γ3. N=N1+N2+N3.The unknown points in the boundary units are called “nodes.” For the requirements of different difference, the expressions ofT and q to be found out can be written as(21)T=∑i=1NφiξTi,q=∑i=1Nφiξqi,where N is the interpolation node number in a boundary unit and φ is the interpolation function. When j is in Γ1orΓ2, (19) can be written as(22)12Ti+∑j=1N1+N2∫Γjq∗TdΓ=∑j=1N1+N2∫ΓjT∗qdΓ-∫ΩfT∗dΩ.This paper uses linear interpolation. Because φ1ξ=1-ξ/2, φ2ξ=1+ξ/2 is the interpolation function of linear unit, the straight line is used to approximate the boundary curve, and the two linearized terminal values are used to approximate the values of T and q. So the node values can be converted into the following:(23a)Tξ=φ1ξTj+φ2ξTj+1=φ1ξ,φ2ξTjTj+1,(23b)qξ=φ1ξqj+φ2ξqj+1=φ1ξ,φ2ξqjqj+1.Substituting (23a) and (23b) into (22), (24)12Ti+∑j=1N1+N2hij1Tj+hij2Tj+1=∑j=1N1+N2gij1qj+gij2qj+1-∫ΩfT∗dΩ,where(25)hije=∫Γjφeq∗dΓ,gije=∫ΓjφeT∗dΓ,e=1,2.Let (26)Hij=hij-12+hij1,Gij=gij-12+gij1,and we have(27)12Ti+∑j=1N1+N2HijTj=∑j=1N1+N2Gijqj-∫ΩfT∗dΩ.Similarly, when the point locates in Γ3, (24) becomes(28)12Ti+∑j=1N3hij1+hλgij1Tj+hij2+hλgij2Tj+1=∑j=1N3gij1qj+gij2qj+1-∫ΩfT∗dΩ.Let(29)Hij=hij-12+hλgij-12+hij1+hλgij1,Gij=gij-12+gij1,and the boundary integral equation is(30)12Ti+∑j=1N3HijTj=∑j=1N3Gijqj-∫ΩfT∗dΩ.Rearranging (30), we have(31)∑j=1NH^ijTj=∑j=1NGijqj-∫ΩfT∗dΩ,where (32)H^ij=Hij,i≠j,Hij+12,i=jconsidering Bi=∫ΩfT∗dΩ, (31) can be expressed as a matrix form HT=GQ-Bi. The boundary conditions T and Q are written as T=T-1T2, Q=Q1Q-2, where T-1, Q-2 are known temperature and heat flux. Equation (31) becomes(33)H1H2T-1T2=G1G2Q1Q-2-Bi.Put the known quantities on the left.(34)H2-G1T2Q1=G2-H1Q-2T-1-Bi.Appling the radial integration method, domain integration caused by heat source is converted to boundary integration [19–21].(35)Bi=∫ΩfT∗dΩ=∫ΓFR∂r∂ndΓ,where R is the distance between the source and the field point and F is radial integration which is expressed as(36)F=∫0RfT∗rdr.Becausef is a known function and the function form is simple, T∗ is the basic function, the radial integration can be found out through (36), and the boundary integration can be obtained from (35).Thus, (34) can be written as AX=K to solve the above equations and to find T and Q.


## 2.1. The Boundary Integral Equation

The mathematical model of two-dimensional steady heat conduction problems with internal heat source is(1)∂2T∂x2+∂2T∂y2=f∈Ω,T=T-∈Γ1,∂T∂n=q-∈Γ2,∂T∂n+hλT=q′′-∈Γ3,where Γ1, Γ2, and Γ3 are the boundary of the domain Ω, Γ=Γ1+Γ2+Γ3, f=-s/λ, s is the internal heat source, h is the heat conduction coefficient between solid and fluid, λ=HTC, T= temp, and n is the normal vector outside the boundary. q′′-=h/λTf, Tf is environment temperature, and q- is heat flux.Weight functionT∗ is introduced into weighted residual of control equation. (2)∫Ω∇2T-fT∗dΩ=∫Γ3∂T∂n+hλT-q′′-T∗dΓ+∫Γ2∂T∂n-q-T∗dΓ-∫Γ1T-T-∂T∗∂ndΓ.Green’s function of Laplace equation is∬Dv∇2u-u∇2vdΩ=∫cv∂u/∂n-u∂v/∂nds, where C is the boundary curve of plane closed area D and ds is arc differential [18].To decompose the left side of the equation, it becomes(3)∫Ω∇2T-fT∗dΩ=∫Ω∇2TT∗dΩ-∫ΩfT∗dΩ.According to Laplace Green’s function, we can have(4)∫Ω∇2TT∗dΩ=∫ΩT∗∇2T-T∇2T∗dΩ+∫ΩT∇2T∗dΩ.Further reduction gives(5)∫Ω∇2TT∗dΩ=∫ΓT∗∂T∂n-T∂T∗∂ndΓ+∫ΩT∇2T∗dΩ.Substituting (5) for the same term in (3), (2) becomes(6)∫ΓT∗∂T∂n-T∂T∗∂ndΓ+∫ΩT∇2T∗dΩ-∫ΩfT∗dΩ=∫Γ3∂T∂n+hλT-q′′-T∗dΓ+∫Γ2∂T∂n-q-T∗dΓ-∫Γ1T-T-∂T∗∂ndΓ.Further reduction gives(7)∫ΩT∇2T∗dΩ=∫Γ3∂T∂n+hλT-q′′-T∗dΓ+∫Γ2∂T∂n-q-T∗dΓ-∫Γ1T-T-∂T∗∂ndΓ-∫ΓT∗∂T∂n-T∂T∗∂ndΓ+∫ΩfT∗dΩ.T
∗ is the basic solution of Laplace equation and it can make ∇T∗+δr-ri=0, where the point “i” refers to the unit point concentration. Because ∂T∗/∂n=q∗ and ∂T/∂n=q, (7) can be further reduced.The left side of (7) becomes(8)∫ΩT∇2T∗dΩ=∫ΩT∇2T∗+δr-ri-Tδr-ridΩ.Calculating furthermore, it is shown as follows:(9)∫ΩT∇2T∗+δr-ri-Tδr-ridΩ=0-∫ΩTδr-ridΩ.For the selectivity ofδ function, ∫-∞+∞fxδx-x0dr=f(x0). Then the reduced result on the left side of (7) is(10)-∫ΩTδr-ridΩ=-Ti.Substituting (10) into (7),(11)-Ti=∫Γ3q+hλT-q′′-T∗dΓ+∫Γ2q-q-T∗dΓ-∫Γ1T-T-q∗dΓ-∫ΓT∗q-Tq∗dΓ+∫ΩfT∗dΩ.BecauseΓ=Γ1+Γ2+Γ3, ∫ΓT∗q-Tq∗dΓ is decomposed: (12)∫ΓT∗q-Tq∗dΓ=∫Γ1T∗q-Tq∗dΓ+∫Γ2T∗q-Tq∗dΓ+∫Γ3T∗q-Tq∗dΓ.Substituting (12) for the same term in (11), (13)-Ti=∫Γ3hλTT∗-q′′-T∗+Tq∗dΓ+∫Γ2Tq∗-q-T∗dΓ+∫Γ1T-q∗-T∗qdΓ+∫ΩfT∗dΩ.Further reduction gives(14)-Ti=∫Γ3hλTT∗dΓ+∫Γ3Tq∗dΓ+∫Γ2Tq∗dΓ+∫Γ1T-q∗dΓ-∫Γ3q′′-T∗dΓ+∫Γ2q-T∗dΓ+∫Γ1qT∗dΓ+∫ΩfT∗dΩ.The second and the third terms are combined, which makes(15)Ti+∫Γq∗TdΓ+∫Γ3hλTT∗dΓ=∫ΓT∗qdΓ-∫ΩfT∗dΩ.The two-dimensional steadyT∗=1/2πln⁡1/r, r=x-xi2-y-yi2.Moving the source pointi to the boundary, the integral equation at any point on the boundary can be obtained. Thus, the basic solution T∗ and the singularity of the integral ∂T∗/∂n need to be considered. Assuming that a circular arc centering on the point i near the border and ε is the small radius of the circle, (15) can be written as(16)Ti+∫Γ-εq∗TdΓ+∫Γεq∗TdΓ+∫Γ3hλTT∗dΓ=∫Γ-εT∗qdΓ+∫ΓεT∗qdΓ-∫ΩfT∗dΩ,where Γε is the new circular border and Γ-ε is the boundary outside the new circular boundary in Γ. Considering the smooth boundary, the integral mean value theorem is introduced. (17)∫ΓεTq∗dΓ=∫ΓεT∂∂n14πε2dΓ=∫ΓεT-14πε2dΓ=Tξ-14πε22πε2,where ξ is the point on Γε. Let ε→0T(ξ)→Ti; then(18)limε→0⁡∫ΓεTq∗dΓ=-12Ti,limε-0⁡∫Γ-εq∗TdΓ=∫Γq∗TdΓ,limε-0⁡∫Γ-εqT∗dΓ=∫Γq∗TdΓ.Synthesizing the above limit results, we have(19)12Ti+∫Γq∗TdΓ+∫Γ3hλTT∗dΓ=∫ΓT∗qdΓ-∫ΩfT∗dΩ.Considering (15) and (19) simultaneously, boundary integral equation of any point in the domain Ω or on the boundary can be found out:(20)CiTi+∫Γq∗TdΓ+∫Γ3hλTT∗dΓ=∫ΓT∗qdΓ-∫ΩfT∗dΩ,Ci=1∈Ω,12∈Γ.


## 2.2. Calculation of the Unknown Value Based on Boundary Element Method

Discreteness of boundary integral equations: because the quantitiesT and q on the boundary are known, the only thing to do is divide the boundary discretely. Boundary Γ is divided into N boundary units Γj(j=1,…,N), where there are N1 units on the boundary Γ1, N2 units on the boundary Γ2, and N3 units on the boundary Γ3. N=N1+N2+N3.The unknown points in the boundary units are called “nodes.” For the requirements of different difference, the expressions ofT and q to be found out can be written as(21)T=∑i=1NφiξTi,q=∑i=1Nφiξqi,where N is the interpolation node number in a boundary unit and φ is the interpolation function. When j is in Γ1orΓ2, (19) can be written as(22)12Ti+∑j=1N1+N2∫Γjq∗TdΓ=∑j=1N1+N2∫ΓjT∗qdΓ-∫ΩfT∗dΩ.This paper uses linear interpolation. Because φ1ξ=1-ξ/2, φ2ξ=1+ξ/2 is the interpolation function of linear unit, the straight line is used to approximate the boundary curve, and the two linearized terminal values are used to approximate the values of T and q. So the node values can be converted into the following:(23a)Tξ=φ1ξTj+φ2ξTj+1=φ1ξ,φ2ξTjTj+1,(23b)qξ=φ1ξqj+φ2ξqj+1=φ1ξ,φ2ξqjqj+1.Substituting (23a) and (23b) into (22), (24)12Ti+∑j=1N1+N2hij1Tj+hij2Tj+1=∑j=1N1+N2gij1qj+gij2qj+1-∫ΩfT∗dΩ,where(25)hije=∫Γjφeq∗dΓ,gije=∫ΓjφeT∗dΓ,e=1,2.Let (26)Hij=hij-12+hij1,Gij=gij-12+gij1,and we have(27)12Ti+∑j=1N1+N2HijTj=∑j=1N1+N2Gijqj-∫ΩfT∗dΩ.Similarly, when the point locates in Γ3, (24) becomes(28)12Ti+∑j=1N3hij1+hλgij1Tj+hij2+hλgij2Tj+1=∑j=1N3gij1qj+gij2qj+1-∫ΩfT∗dΩ.Let(29)Hij=hij-12+hλgij-12+hij1+hλgij1,Gij=gij-12+gij1,and the boundary integral equation is(30)12Ti+∑j=1N3HijTj=∑j=1N3Gijqj-∫ΩfT∗dΩ.Rearranging (30), we have(31)∑j=1NH^ijTj=∑j=1NGijqj-∫ΩfT∗dΩ,where (32)H^ij=Hij,i≠j,Hij+12,i=jconsidering Bi=∫ΩfT∗dΩ, (31) can be expressed as a matrix form HT=GQ-Bi. The boundary conditions T and Q are written as T=T-1T2, Q=Q1Q-2, where T-1, Q-2 are known temperature and heat flux. Equation (31) becomes(33)H1H2T-1T2=G1G2Q1Q-2-Bi.Put the known quantities on the left.(34)H2-G1T2Q1=G2-H1Q-2T-1-Bi.Appling the radial integration method, domain integration caused by heat source is converted to boundary integration [19–21].(35)Bi=∫ΩfT∗dΩ=∫ΓFR∂r∂ndΓ,where R is the distance between the source and the field point and F is radial integration which is expressed as(36)F=∫0RfT∗rdr.Becausef is a known function and the function form is simple, T∗ is the basic function, the radial integration can be found out through (36), and the boundary integration can be obtained from (35).Thus, (34) can be written as AX=K to solve the above equations and to find T and Q.


## 3. The Steady Inverse Problem

### 3.1. The Objective Function of Inverse Problem

For aforementioned heat conduction system, the heat transfer coefficient between solid and fluidh and the thermal conductivity λ are known. The boundary temperature to be found out is unknown which can be determined based on the known measured boundary temperature and the known condition of forward problem.The following objective function is defined as(37)Jq=∑i=1NTiq-Timea2.In the objective function,q is a parameter vector whose temperature needs to be inversed and N is the number of temperature measuring points on the boundary.(38)q=r1,r2,…,rm=T11,T21,…,Tm1.T
i
(
q
) is the calculation of measuring point i in the forward problem and Tmea is the measurement of temperature. The minimum of the objective function J(q) is the parameter vector q of inverse problem.


### 3.2. Conjugate Gradient Method of Inverse Problem

Conjugate gradient method is a method which combines the conjugacy and the steepest descent method. It derives from perturbation principle and the inverse problem is converted into three questions, such as forward problem, sensitivity problem, and adjoint problem. In order to solve the effects of these three questions, this paper introduces the complex variable derivation method into the traditional conjugate gradient method, which makes the calculation of the sensitivity coefficient accurate.During the calculation, unconstrained optimization algorithm is achieved by iteration. Considering that thekth iteration point qk has been available, the (k+1)th iteration calculates according to the following formula.The iterative equation solving the inverse problem by conjugate gradient method includes(39)qk+1=qk+αkdk,qk+1=qk-αkdk,dk=-gk,k=0,-gk+βkdk-1,k≥1,where αk is step length, obtained by searching some one-dimensional line; dk is the searching direction in which gk=∇J(qk) and βk is a scalar. Different βk corresponds to different nonlinear conjugate gradient methods.The searching step lengthαk, the conjugate coefficient βk, and gradient ∇Jqk are needed to be found out.Conjugate coefficientβk equals the ratio of the square of normal form between the current and the previous-step-gradient paradigm:(40)βk=∇Jqk-∇Jqk-1∇JTqk∇Jqk-1∇JTqk-1,(41)∇Jqk=∂J∂r1,∂J∂r2,…,∂J∂rm+1,(42)∂Jqk∂rj=2∑i=1NTiqk-Ti∗∂Tiqk∂rjj=1,2,…,m+1.If the iteration step is k+1 and searching step length is αk, (37) Jqk+1=∑i=1NTiqk+αkdk-T∗2 gives(43)αk=∑i=1NTiqk-Ti∗∇Tidk∑i=1N∇Tidk2,where ∇Ti=∂Ti/∂r1,∂Ti/∂r2,…,∂Ti/∂rm+1.


### 3.3. Complex Variable Derivation Method [22–24]

For any real functionf(x), a very small imaginary part h is added to the real variable x. It is expressed in complex function fx+ih and its Taylor series is(44)fx+ih=fx+ihf′x-h22f′′x+oh3.As h is very small,(45)dfdx=Im⁡fx+ihh,d2fdx2=2fx-Re⁡fx+ihh2,where ∂Ti(qk)/∂rj is sensitivity coefficient and its matrix form is(46)∂T∂r=∂T1r∂r1∂T1r∂r2⋯∂T1r∂rm∂T2r∂r1∂T2r∂r2⋯∂T2r∂rm⋮⋮⋮⋮∂TNr∂r1∂TNr∂r2⋯∂TNr∂rm.The complex variable derivation method is used to calculate the matrix of sensitivity coefficient.(47)∂T∂r=Im⁡T1r1+hhIm⁡T1r2+hh⋯Im⁡T1rm+hhIm⁡T2r1+hhIm⁡T2r2+hh⋯Im⁡T2rm+hh⋮⋮⋮⋮Im⁡TNr1+hhIm⁡TNr2+hh⋯Im⁡TNrm+hh.The complex variable derivation method is introduced into the traditional conjugate gradient method, which makes the calculation of the sensitivity coefficient accurate and avoids the sensitivity problem and the adjoint problem.


### 3.4. Solving Process of the Inverse Problem

Initialization:k=0, q0,r0=0, and let ε be a small positive number.Solve the forward problem given by (1), calculate temperature Ti, and judge if the following condition is true:(48)Jq=∑i=1NTiq-Timea2<ε.The equation means the condition required by the instruction DO WHILE has not been satisfied. Calculating the gradient∇J using (41)Calculating the searching-down directiondk using (27)Calculating the searching step lengthαk using (43)Calculating the new estimation using (34)Solving the forward problem given by (1) and obtaining Ti:(49)k=k+1End do


## 3.1. The Objective Function of Inverse Problem

For aforementioned heat conduction system, the heat transfer coefficient between solid and fluidh and the thermal conductivity λ are known. The boundary temperature to be found out is unknown which can be determined based on the known measured boundary temperature and the known condition of forward problem.The following objective function is defined as(37)Jq=∑i=1NTiq-Timea2.In the objective function,q is a parameter vector whose temperature needs to be inversed and N is the number of temperature measuring points on the boundary.(38)q=r1,r2,…,rm=T11,T21,…,Tm1.T
i
(
q
) is the calculation of measuring point i in the forward problem and Tmea is the measurement of temperature. The minimum of the objective function J(q) is the parameter vector q of inverse problem.


## 3.2. Conjugate Gradient Method of Inverse Problem

Conjugate gradient method is a method which combines the conjugacy and the steepest descent method. It derives from perturbation principle and the inverse problem is converted into three questions, such as forward problem, sensitivity problem, and adjoint problem. In order to solve the effects of these three questions, this paper introduces the complex variable derivation method into the traditional conjugate gradient method, which makes the calculation of the sensitivity coefficient accurate.During the calculation, unconstrained optimization algorithm is achieved by iteration. Considering that thekth iteration point qk has been available, the (k+1)th iteration calculates according to the following formula.The iterative equation solving the inverse problem by conjugate gradient method includes(39)qk+1=qk+αkdk,qk+1=qk-αkdk,dk=-gk,k=0,-gk+βkdk-1,k≥1,where αk is step length, obtained by searching some one-dimensional line; dk is the searching direction in which gk=∇J(qk) and βk is a scalar. Different βk corresponds to different nonlinear conjugate gradient methods.The searching step lengthαk, the conjugate coefficient βk, and gradient ∇Jqk are needed to be found out.Conjugate coefficientβk equals the ratio of the square of normal form between the current and the previous-step-gradient paradigm:(40)βk=∇Jqk-∇Jqk-1∇JTqk∇Jqk-1∇JTqk-1,(41)∇Jqk=∂J∂r1,∂J∂r2,…,∂J∂rm+1,(42)∂Jqk∂rj=2∑i=1NTiqk-Ti∗∂Tiqk∂rjj=1,2,…,m+1.If the iteration step is k+1 and searching step length is αk, (37) Jqk+1=∑i=1NTiqk+αkdk-T∗2 gives(43)αk=∑i=1NTiqk-Ti∗∇Tidk∑i=1N∇Tidk2,where ∇Ti=∂Ti/∂r1,∂Ti/∂r2,…,∂Ti/∂rm+1.


## 3.3. Complex Variable Derivation Method [22–24]

For any real functionf(x), a very small imaginary part h is added to the real variable x. It is expressed in complex function fx+ih and its Taylor series is(44)fx+ih=fx+ihf′x-h22f′′x+oh3.As h is very small,(45)dfdx=Im⁡fx+ihh,d2fdx2=2fx-Re⁡fx+ihh2,where ∂Ti(qk)/∂rj is sensitivity coefficient and its matrix form is(46)∂T∂r=∂T1r∂r1∂T1r∂r2⋯∂T1r∂rm∂T2r∂r1∂T2r∂r2⋯∂T2r∂rm⋮⋮⋮⋮∂TNr∂r1∂TNr∂r2⋯∂TNr∂rm.The complex variable derivation method is used to calculate the matrix of sensitivity coefficient.(47)∂T∂r=Im⁡T1r1+hhIm⁡T1r2+hh⋯Im⁡T1rm+hhIm⁡T2r1+hhIm⁡T2r2+hh⋯Im⁡T2rm+hh⋮⋮⋮⋮Im⁡TNr1+hhIm⁡TNr2+hh⋯Im⁡TNrm+hh.The complex variable derivation method is introduced into the traditional conjugate gradient method, which makes the calculation of the sensitivity coefficient accurate and avoids the sensitivity problem and the adjoint problem.


## 3.4. Solving Process of the Inverse Problem

Initialization:k=0, q0,r0=0, and let ε be a small positive number.Solve the forward problem given by (1), calculate temperature Ti, and judge if the following condition is true:(48)Jq=∑i=1NTiq-Timea2<ε.The equation means the condition required by the instruction DO WHILE has not been satisfied. Calculating the gradient∇J using (41)Calculating the searching-down directiondk using (27)Calculating the searching step lengthαk using (43)Calculating the new estimation using (34)Solving the forward problem given by (1) and obtaining Ti:(49)k=k+1End do


## 4. The Instance Calculation and Analysis

To testify the availability of the aforementioned method, the simulation experiment having two sets of data was designed. The experiment discussed the effect of the testing points’ number and the measuring error on the inversion results and compared CGM and L-MM. The schematic diagram is shown as Figure1.Figure 1
Model of the instance.In the experiment, the concrete slab’s thickness isδ=0.3 m, and the thickness of casting concrete slab between two wooden templates is δ1=δ2=0.02 m, and coefficient of thermal conductivity is λ1=λ2=0.16 w/(mK), the air temperature is Tf=293 K, andD1 andD3 are the convective heat conduction boundary conditions. The heat transfer coefficients are hi=8 W/(m2K) and he=20 W/(m2K). The thermal conductivity coefficient is λ=1.5 W/(mK) and the power density of inner heat source is s=1000 W/m3 (similar to even-distributed and constant in short time). The real temperature distribution of the boundaryD4 to be solved is(50)y=311+48x-1000x23.According to the real temperatures, the calculated temperaturesTi of testing points on the known boundaryD2 are counted through the forward question.(51)Timea=Ti+ωσ,ε=∑i=1Nσ2,where ω is standard-normal-distributed random number and σ is the standard deviation of measuring temperatures.


### 4.1. The Influence of Measuring Points’ Number on Inversion

When the measurement error is zero,Tf=293 K, hi=8 W/(m2K), and he=20 W/(m2K). The thermal conductivity coefficient is λ=1.5 W/(mK) and the power density of inner heat source is s=1000 W/m3. The influence of measuring points’ number on inversion is shown as Figure 2. When the measuring point N = 5, 9, or 12, the inversion results reveal that the more points are measured, the more accuracy of inversion can be achieved.Figure 2
The influence of measuring points’ number on inversion.


### 4.2. The Influence of the Measuring Error on Inversion

N
=
12, Tf=293 K, hi=8 W/(m2K), and he=20 W/(m2K). The thermal conductivity coefficient is λ=1.5 W/(mK), and the power density of inner heat source is s=1000 W/m3. The influence of different measuring error on inversion is shown as Figure 3. As the measuring standard deviation σ=0,0.2,0.4, the inversion results reveal that the method can obtain better inversion results for a relatively small measurement standard deviation.Figure 3
The influence of different measuring error on inversion.


### 4.3. The Comparison between CGM and L-MM

This paper applies CGM and L-MM to calculate the inversion and compare with each other under three conditions:σ=0, N=9, and σ=0, N=12, and N=12, σ=0.2 (other conditions are the same as the aforementioned).The inversions under the conditionsσ=0, N=9, and σ=0, N=12 are shown as Figures 4 and 5. Figures 4 and 5 show that the inversions using CGM and L-MM are similar and both have high accuracy. Both of them are satisfied.Figure 4
The inversion under the condition:σ=0 and N=9.Figure 5
The inversion under the condition:σ=0 and N=12.The inversions under the conditionsσ=0.2, N=9 are shown as Figure 6. Comparing with Figure 5, the inversions using L-MM cause obvious fluctuation when the measuring error increases, while the inversions using CGM are satisfied. That means CGM has the better stability.Figure 6
The inversion under the condition:σ=0.2 and N=12.


## 4.1. The Influence of Measuring Points’ Number on Inversion

When the measurement error is zero,Tf=293 K, hi=8 W/(m2K), and he=20 W/(m2K). The thermal conductivity coefficient is λ=1.5 W/(mK) and the power density of inner heat source is s=1000 W/m3. The influence of measuring points’ number on inversion is shown as Figure 2. When the measuring point N = 5, 9, or 12, the inversion results reveal that the more points are measured, the more accuracy of inversion can be achieved.Figure 2
The influence of measuring points’ number on inversion.


## 4.2. The Influence of the Measuring Error on Inversion

N
=
12, Tf=293 K, hi=8 W/(m2K), and he=20 W/(m2K). The thermal conductivity coefficient is λ=1.5 W/(mK), and the power density of inner heat source is s=1000 W/m3. The influence of different measuring error on inversion is shown as Figure 3. As the measuring standard deviation σ=0,0.2,0.4, the inversion results reveal that the method can obtain better inversion results for a relatively small measurement standard deviation.Figure 3
The influence of different measuring error on inversion.


## 4.3. The Comparison between CGM and L-MM

This paper applies CGM and L-MM to calculate the inversion and compare with each other under three conditions:σ=0, N=9, and σ=0, N=12, and N=12, σ=0.2 (other conditions are the same as the aforementioned).The inversions under the conditionsσ=0, N=9, and σ=0, N=12 are shown as Figures 4 and 5. Figures 4 and 5 show that the inversions using CGM and L-MM are similar and both have high accuracy. Both of them are satisfied.Figure 4
The inversion under the condition:σ=0 and N=9.Figure 5
The inversion under the condition:σ=0 and N=12.The inversions under the conditionsσ=0.2, N=9 are shown as Figure 6. Comparing with Figure 5, the inversions using L-MM cause obvious fluctuation when the measuring error increases, while the inversions using CGM are satisfied. That means CGM has the better stability.Figure 6
The inversion under the condition:σ=0.2 and N=12.


## 5. Conclusion

This paper applies CGM based on complex variable derivation to study the multivariable inverse problem which combines boundary distribution with convection coefficient of two-dimensional steady system with inner heat source. The inversions applying CGM and L-MM are compared. The influence of measuring points’ number and measuring error on inversion is tested. The simulation reveals that applying boundary element method and conjugate gradient method to solve the inverse heat conduction problem is successful and has good stability.



---

*Source: 2861342-2017-07-11.xml*