The estimation of the strength of the heat source in the heat conduction problems

(1)

The estimation of the strength of the heat source

in the heat conduction problems

David T.W. Lin

a,*

, Ching-yu Yang

b

a

Department of Mechanical Engineering, Far East University, No. 49, Chung Hua Road, Hsin-Shih, Tainan Country 744, Taiwan, ROC

b

Department of Mold and Die Engineering, National Kaohsiung University of Applied Sciences, Kaohsiung 807, Taiwan, ROC Received 1 September 2005; accepted 11 October 2006

Available online 21 December 2006

Abstract

A general method is proposed to determine the strength of the heat source in the Fourier and non-Fourier heat con-duction problems. A finite difference method, the concept of the future time and a modified Newton–Raphson method are adopted in the problem. The undetermined heat source at each time step is formulated as an unknown variable in a set of equations from the measured temperature and the calculated temperature. Then, an iterative process is used to solve the set of equations. No selected function is needed to represent the undetermined function in advance. Three exam-ples are used to demonstrate the characteristics of the proposed method. The validity of the proposed method is confirmed by the numerical results. The results show that the proposed method is an accurate and stable method to determine the strength of the heat source in the inverse hyperbolic heat conduction problems. Furthermore, the result shows that more future times are needed in the hyperbolic equation than that of parabolic equation. Moreover, the robustness and the accu-racy of the estimated results in the non-Fourier problem are not as well as those of the Fourier problem.

1. Introduction

The inverse source problems deal with the determination of the strength of the heat source in analysis such as the internal energy source, or the quantity of the energy generation in a computer chip, or in a microwave heat-ing process, or in a chemical reaction process, et al. They have been widely applied in many design and man-ufacturing problems especially when the direct measurements for the problem are difficult. Wide attention has been called to this problem, and most studies investigate the Fourier heat model [1–11] to determine the unknown sources in the inverse problems. The Fourier heat model expresses that the heat flux is directed pro-portional to temperature gradient in that the model has an infinite speed of heat propagation. Huang and Ozisik[1]used the conjugate gradient method combined with the adjoint equation to estimate the strength of heat source in an internal plate. Yang[2–5]solved the source strength problems by a reverse matrix method, a symbolic method and a numerical sequential method based on a linear least-squares error method. Silva and

* _{Corresponding author.}

E-mail address:david@cc.fec.edu.tw(D.T.W. Lin).

(2)

Ozisik[6]investigated the two sources problem by the conjugate gradient method. Niliot and Lefevre[7–10]

used a boundary element method to formulate the problem and lead to a non-linear optimization to identify the source. Khachfe and Jarny[11]used a finite element method combined the conjugate gradient method to determine the heat transfer coefficient and the heat source. However, in many practical engineering problems, the heat model is described as the non-Fourier law such as rapid heating process, or slow conduction process. Many studies of the non-Fourier model have shown in references [12–29]. The non-Fourier heat model describes that the heat transfer propagates at a finite speed. Lor and Chu[12]analyzed the problem with the interface thermal resistance. Antaki[13]discussed the heat transfer in solid-phase reactions. Sanderson et al.

[14]and Liu et al.[15]discussed the laser generated ultrasound models. Mullis[16,17]examined the rapid solid-iﬁcation problems. Roetzel and Ranong[18]and Roetzel et al.[19,20]calculated the heat exchanger problems. The non-Fourier ﬁn problems under the periodic thermal conditions are computed by Lin[21], Abdel-Hamid

Nomenclature

C speciﬁc heat capacity

g(t) the strength of the heat source J error function

k thermal conductivity

Nt number of the temporal measurements

p number of grids at spatial coordinate r number of the future time

t temporal coordinate

T temperature

Tm value of temperature at the mth time step

Xm sensitivity function of T with respect to the undetermined condition at m-time step

x0 vector of the initial guess

x spatial coordinate Y measured temperature

q density

b relaxation time l eigenvalue of matrix U vector constructed from U

U calculated temperature minus measured temperature Uc calculated temperature

Umeas measured temperature

b

Uu component of vector U /q;iq

m unknown heat ﬂux condition at iqth grid and mth time step

W sensitivity matrix Wu,v component of vector W

D increment of the search step e, d value of the stopping criterion

r standard deviation of measurement error ki,j random number

Subscripts

i, j, m, u, v indices Superscripts

exact exact temperature meas measured temperature

(3)

[22], Tang and Araki[23]. As well, the thin film problems are also explored by Tan and Yang[24,25]. Moreover, some works have been done in the inverse hyperbolic heat problem[26–29]. Weber[26]used a 2nd-order explicit difference equation to discrete the inverse problem domain. Al-Khalidy [27,28] adopted a control volume method to discrete the spatial domain and the space marching method to solve the inverse problem. Chen et al.[29]used Laplace transform and control volume method combined with non-linear least-square scheme to estimate the boundary condition. Yang [30,31] developed a forward difference method to estimated the boundary condition of the one-dimensional problems. However, the inverse source problem with the non-Fourier effect is not investigated by the past research. Therefore, it is necessary to develop an accurate and sta-ble method to deal with the inverse source prosta-blem with the non-Fourier effect. Furthermore, the estimated results of the Fourier and non-Fourier problems are also made to compare each other.

In this paper, a sequential method combined with the concept of the future time[32]is proposed to solve the problems progressively. As well, a modified Newton–Raphson method [33,34]is used to search the inverse solution at each time step. In the proposed approach, the determination of the source’s strength at each time step includes two phases, one is the process of direct analysis and the other is the process of inverse analysis. In the forward analysis process, the strength of the heat source is assumed as the known values and then directed to solve the temperature profile of the hyperbolic heat conduction equation through a finite difference method. Solutions from the above process are substituted into the sensitivity analysis and integrated with the available temperature measured at the sensor’s location. Thus, a set of non-linear equation is formulated for the process of the inverse estimation. An iterative method is used to guide the exploring points systematically to approach to the undetermined source in the inverse analysis process. Then, the intermediate value is substituted for the unknown source in the following analysis. As such, several iterations are needed for obtaining the undeter-mined strength of the heat source. In the present research, the proposed method formulates the problem from the difference between the calculated temperature and the one measured directly. Therefore, the inverse for-mulation derived from the proposed method is simpler than that from the non-linear least-squares method.

This paper includes five sections. In the first section, the current researches in the fields of inverse estimation of the heat source are introduced and the feature of using the proposed method in the problem is also stated. Next, a finite difference formulation for the hyperbolic equations is stated and the stable condition for the algorithm is introduced in the second section. In the third section, the characteristics of solving the inverse problem are delineated and the content of the concept of the future time, the direct problem, the sensitivity problem, and the algorithm of the proposed method are presented. Meanwhile, the criterion to stop the iter-ative process is illustrated. In the fourth section, the computational algorithm of the proposed method is shown. Three examples are employed to demonstrate and discuss the results of the proposed method in the fourth section. Finally, the overall contribution of this research to the field of inverse heat conduction problem is discussed in the final section.

2. Problem statement

The inverse heat source problem in one spatial dimensional consists of ﬁnding the strength at the medium while the temperature measurements at the front are given. Consider a homogeneous slab with l thickness and constant thermal properties. This slab originally has a uniformly distributed temperature. The adiabatic con-dition is applied to x = 0 and x = l. At a speciﬁc time, a heat source g(t) is applied to x = xs. A mathematical

formation of the heat conduction is presented as follows: ko 2_{Tðx; tÞ} ox2 þ gðtÞdðx xsÞ ¼ bqC o2Tðx; tÞ ot2 þ qC oTðx; tÞ ot t >0; 0 < x < l; ð1Þ Tðx; 0Þ ¼ T0 0 6 x 6 l; ð2Þ oTðx; tÞ ot ¼ 0 at t¼ 0 and 0 6 x 6 l; ð3Þ oTðx; tÞ ox ¼ 0 at x¼ 0 and t >0; ð4Þ koTðx; tÞ ox ¼ 0 at x ¼ 0 and x¼ l; t > 0; ð5Þ

(4)

where T represent the temperature ﬁeld T(x, t). k is the thermal conductivity and qC is the heat capacity per unit volume. b is the relaxation time that is non-negative. The value of b is not vanished when the problem is a non-Fourier type.

The inverse problem is given the temperature measured at x = 0 to estimate the strength of the heat source g(t).

3. The direct solution of the heat conduction equations

The proposed method uses a finite difference method with the equidistant grids in the spatial and temporal coordinates. The finite difference method has been implemented in the researches of Weber[26]and Carey and Tsai[35]. However, the stable condition is not clearly defined by the above researches and the relaxation time is limited to a small value. Therefore, the following derivation investigates the stable condition for solving the hyperbolic condition equation through eigenvalue analysis[30,31]. In this study, the spatial step size is Dx and the temporal-step size is Dt. The differential terms oTðx;tÞ_ot , o2T_otðx;tÞ2 and

o2_Tðx;tÞ

ox2 can be approached by the Taylor series in x = xiand t = tjas follows:

oT otðxi; tjÞ ¼ Tðxi; tjþ DtÞ T ðxi; tjÞ Dt Dt 2 o2T ot2ðxi;gjÞ; ð6Þ o2T ot2 ðxi; tjÞ Tðxi; tj DtÞ 2T ðxi; tjÞ þ T ðxi; tjþ DtÞ Dt2 Dt2 12 o4T ot4ðxi;jjÞ; ð7Þ o2T ox2ðxi; tjÞ ¼ Tðxi Dx; tjÞ 2T ðxi; tjÞ þ T ðxiþ Dx; tjÞ Dt2 Dx2 12 o4T ox4ðmi; tjÞ; ð8Þ

where gj2 (tj, tj+ Dt), jj2 (tj Dt, tj+ Dt) and mi2 (xi Dx, xi+ Dx) Accordingly, Eq.(1)can be discretized

as the follows: bqCTðxi; tj DtÞ 2T ðxi; tjÞ þ T ðxi; tjþ DtÞ Dt2 þ qC Tðxi; tjþ DtÞ T ðxi; tjÞ Dt kTðxi Dx; tjÞ 2T ðxi; tjÞ þ T ðxiþ Dx; tjÞ Dx2 gðtjþ1Þdðxi xsÞ ¼ si;j; ð9Þ where si;j¼ qCDt2 o2T ot2ðxi;gjÞ þ bqCDt 2 12 o4T ot4ðxi;jjÞ kDx 2 12 o4T

ox4ðmi; tjÞ and si,jis the error term of Taylor

approxima-tion. After neglect the error term si,jand the diﬀerence equation is shown as

bqCTi;j1 2Ti;jþ Ti;jþ1

Dt2 þ qC

Ti;jþ1 Ti;j

Dt k

Ti1;j 2Ti;jþ Tiþ1;j

Dx2 gjþ1dðxi xsÞ ¼ 0 ð10Þ and Ti;jþ1¼ kTi1;jþ 2 1 Dt 2ak k Ti;jþ kTiþ1;j k b aTi;j1þ Dt2 bqCgjþ1dðxi xsÞ; ð11Þ where a¼ kDt2 qcDx2 and k¼ a bþDt¼ 1 bþDt kDt2 qcDx2. T_1;jþ1 T_2;jþ1 .. . T_p1;jþ1 2 6 6 6 4 3 7 7 7 5¼ 2 1Dt 2ak k k 0 0 k 2 1Dt 2ak k k 0 0 k 0 0 k 2 1Dt 2ak k 2 6 6 6 6 4 3 7 7 7 7 5 ðp1Þðp1Þ T1;j T2;j .. . T_p1;j 2 6 6 6 4 3 7 7 7 5 kb a T1;j1 T2;j1 .. . Tp1;j1 2 6 6 6 4 3 7 7 7 5þ kT0;j 0 .. . 0 kTp;j 2 6 6 6 6 4 3 7 7 7 7 5þ Dt2 bqC g_jþ1dðx1 xsÞ g_jþ1dðx2 xsÞ .. . g_jþ1dðxp1 xsÞ 2 6 6 6 4 3 7 7 7 5 ð12Þ where p is the grid number of spatial coordinate.

(5)

The ith eigenvalue of the matrix is l_i¼ 2 Dt bþ Dt 4k sin ip 2p 2 : ð13Þ

Therefore, the condition for stability is

max 2 Dt bþ Dt 4k sin ip 2p 2 61 where i¼ 1; 2; . . . ; p 1; ð14Þ 1 4 1 Dt bþ Dt 6_k _sin ip 2p 2 61 4 3 Dt bþ Dt : ð15Þ

The stability requires that this inequality condition hold as Dx! 0, i.e. p ! 1 lim p!1 sin 2 p 1 2p p ¼ 1: ð16Þ

Consequently, that stability region is conﬁned in 1 4 1 Dt bþ Dt 6 1 bþ Dt kDt2 qcDx26 1 4 3 Dt bþ Dt : ð17Þ

4. The method to determine the strength of the heat source

In each time step, an iterative algorithm is used to estimate the strength of the heat source while the tem-perature measured at the front. Some treatments are needed in the process of solving the inverse problem. There are the forward problem, the sensitivity problem, the operational algorithm, and the stopping criterion. The forward problem is used to determine the temperature profile and the sensitivity problem is used to find the search step in the inverse problem. The operational algorithm is used to fulfill the process of the inverse analysis when the solutions of the forward problem and the sensitivity problem are available. Finally, the stop-ping criterion is an important index to decide the termination of the iterative process.

4.1. The forward problem

The proposed method is based on a sequential algorithm and the inverse solution is solved at each time step. Hence, Eqs. (1)–(5) are limited to only one time step and the transient problem at t = tmis governed

by the following equations: ko 2_{Tðx; t} mÞ ox2 þ gmþ1dðx xsÞ ¼ bqC o2Tðx; tmÞ ot2 þ qC oTðx; tmÞ ot at t¼ tm and 0 < x < l; ð18Þ Tðx; tm1Þ ¼ Tm1 0 6 x 6 l; ð19Þ oTðx; tÞ ot ¼ 0 at t¼ tm1 and 0 6 x 6 l; ð20Þ oTðx; tÞ ox ¼ 0 at x¼ 0 and t¼ tmþ1; ð21Þ koTðx; tÞ ox ¼ 0 at x ¼ 0 and x¼ l and t¼ tmþ1; ð22Þ where m = 1, 2, 3, . . ..

Eqs.(18)–(20) are solved by the forward finite difference method and Eq.(18) spans the temporal index from t = tm1 to t = tm+1 (see Eq. (10)). However, the temperature field at t = tm+1 is evaluated based on

the known temperature ﬁeld at t = tm1and t = tm(i.e., Eqs. (19) and (20)). Therefore, the temporal index

(6)

The inverse solution of the above problem is ill-posed and it is often unstable when the measured data has a slight variation in the experimental measurements. The concept of future time is used to improve the stability of the estimation in this research. The future time is included in the measurement to estimate the present state and aﬀects the behavior of the experimental data at future time steps. In the present research, the proposed method formulates the problem from the diﬀerence between the calculated temperature and the one measured directly. As well, equation solver instead of optimization algorithm solves the inverse problem.

When the computed time is t = tm, the estimated condition between t = tm1and t = tmhas been evaluated

and the problem is to estimate the strength of the heat source at t = tm+1. To stabilize the estimated results in

the inverse algorithms, the sequential procedure is assumed temporally that several future values of the esti-mation are constant[30]. Then, the unknown conditions at the future time are assumed to be constants equal to their present values, i.e.

g_mþ2¼ ¼ gmþr1¼ gmþr¼ gmþ1 ð23Þ

Here r is the number of the future time.

The forward problem Eqs. (18)–(22) are solved in r steps (from t = tm+1to tm+r) and the undetermined

boundaries are set by Eq.(23). 4.2. The sensitivity problem

In this paper, the sensitivity analysis of the inverse problem is examined numerically by the modiﬁed New-ton–Raphson method to decide the search step in the each iteration. The derivative o

ogmþ1is taken at both sides of Eqs. (18)–(22). Then, we have

ko 2_X_{ðx; t} mÞ ox2 dðx xsÞ ¼ bqC o2Xðx; tmÞ ot2 þ qC oXðx; tmÞ ot at t¼ tm and 0 < x < l; ð24Þ Xðx; tm1Þ ¼ 0 0 6 x 6 l; ð25Þ oXðx; tÞ ot ¼ 0 at t¼ tm1 and 0 6 x 6 l; ð26Þ oXðx; tÞ ox ¼ 0 at x¼ 0 and t¼ tmþ1; ð27Þ koXðx; tÞ ox ¼ 0 at x ¼ 0 and x¼ l and t¼ tmþ1; ð28Þ where Xm¼oTðx_ogi;tmÞ mþ1 :

Eqs. (24)–(28) describe the mathematical equations for sensitivity coeﬃcient Xm that can be explicitly

found. The equations are the linear equations and the dependent variable Xmis with respect to the

indepen-dent variables and Therefore, the sensitive data can be determined directly through a finite difference method. 4.3. A modified Newton–Raphson method

The Newton–Raphson method[36] has been widely adopted to solve a set of non-linear equations. This method is applicable to solve the non-linear problem when the number of the equations and the number of the unknown variables are the same. In the inverse problem, the number of equations is usually larger than the number of variables; hence a modiﬁed version of the Newton–Raphson method is necessary to deal with the inverse problem.

In the study, the problem is formulated by the proposed method from the comparison between the calcu-lated temperature and the one measured directly. For that reason, the calcucalcu-lated temperature Ucði; jÞ and the

measured temperature Umeasði; jÞ at the i-grid of the spatial coordinate and j at grid of the temporal coordinate

is necessary to be evaluated ﬁrst. Then, the estimation of the unknown source at each time step can be recast as the solution of a set of non-linear equations:

Uði; jÞ ¼ Ucði; jÞ Umeasði; jÞ ¼ 0; ð29Þ

(7)

The number of equations is the number of the future time r. This detail procedure can be shown as follows. Substitute the index j from m + 1 to m + r and the index i¼ 0, we have

U¼ ½Uð0; m þ 1Þ; Uð0; m þ 2Þ; Uð0; m þ 3Þ; ; Uð0; m þ rÞT ¼ f bUug ð30Þ

where bUuis the component of vector U

The undetermined coeﬃcient is set as v. The derivative of bUuwith respect to v is solved through Eqs.(24)–

(28)and it can be expressed as follows:

Wu¼

o bUu

ov : ð31Þ

The sensitivity matrix W can be deﬁned as follows:

W¼ fWug; ð32Þ

where u = 1, 2, 3, . . . , r and Wuis the element W at uth row. With the starting vector v0and the above

deriva-tions from Eqs.(30)–(32), we have the following equation:

v_kþ1¼ vkþ Dk: ð33Þ

Dkis a linear least-squares solution for a set of over-determined linear equations and it can be derived as

follows:

WðvkÞDk¼ UðvkÞ; ð34Þ

Dk¼ ½WTðvkÞWðvkÞ 1_WT

ðvkÞUðvkÞ: ð35Þ

The above derivation is applied at each time step. This method is implemented in the multi-sensors’ mea-surement. Under this condition, the number of the elements in Eq.(30)is based on the number of measured locations and the number of future times.

4.4. The stopping criteria

The modiﬁed Newton–Raphson method (Eqs.(33)–(35)) is used to determine the unknown value v. The step size Dkgoes from vkto vk+1and it is determined from Eq.(33). Once Dkis calculated, the iterative process

to determine vk+1is executed until the stopping criterion is satisﬁed.

The discrepancy principle[37]is widely used to evaluate the value of the stopping criterion in the inverse technique. Nevertheless, the convergence of the inverse solution is not guaranteed by the stopping criterion. Accordingly, two criteria used by Frank and Wolfe[38]are chosen to assure the convergence and to stop the iteration: kvkþ1 vkk 6 dkvkþ1k; ð36Þ kJðvkþ1Þ JðvkÞk 6 ekJðvkþ1Þ; ð37Þ wherekJðvkþ1Þk ¼¼ Xp i¼1 Xr j¼1

½Ucði; jÞ Umði; jÞ 2

; ð38Þ

where d and e are small positive values.

The values of d and e are the converge tolerances. 5. Computational algorithm

The iterative procedure for the proposed method can be summarized as follows: First, we choose the num-ber of future times r, the mesh conﬁguration of the problem domain, the temporal size Dt, the measured grid and the estimated grid. Given overall convergence tolerance d and e and the initial guess v0. The value of vkis

(8)

Step 1. Let j = m and the temperature distribution at {Tj1} and {Tj} are known.

Step 2. Collect the measurement Umeasði; jÞ which are Yijl; Y il

jþl; L; Y il

jþrl.

Step 3. Assume the initial guess v0.

Step 4. Solve the forward problem (Eqs.(18)–(22)), and compute the calculated temperature Ucði; jÞ.

Step 5. Integrate the calculated temperature Ucði; jÞ with the measured temperature Umeasði; jÞ to

construct U.

Step 6. Calculate the sensitivity matrix W through Eqs.(24)–(28). Step 7. Knowing W and U, compute the step size Dkfrom Eq.(35).

Step 8. Knowing Dkand vk, compute vk+1from Eq.(33).

Step 9. Terminate the process if the stopping criterion (Eqs. (36,37)) is satisﬁed. Otherwise return to Step 4.

Step 10. Terminate the process if the ﬁnal time step is attached. Otherwise, let j = m + 1 return to Step 2.

6. Results and discussion

Our simulations deﬁne from Eqs.(1)–(5)that estimate the strength of the heat source. Three diﬀerent source functions over temporal domain; namely, a triangular function, a sinusoidal function, and a quarter sinusoidal function are adopted to illustrate the numerical modeling. The stability and the accuracy of the estimation are discussed and compared to the solutions of Yang’s[3]approach that is the Fourier heat problem. The exact temperature and the source strength used in the following examples are selected so that these functions can satisfy Eqs. (1)–(5). The accuracy is assessed by the comparison between the estimated and preset source strength. Meanwhile, the measurement of the simulated temperature is generated from the exact temperature in each problem and it is supposed to have the measurement errors. In other words, the random errors of mea-surement are added to the exact temperature. It can be shown in the following equation:

Tmeas i;j ¼ T

exact

i;j þ ki;jr ð39Þ

where the subscripts i and j are the grid number of spatial-coordinate and temporal-coordinate respectively. The Texact

i;j and T meas

i;j in Eq.(39)are the exact temperature and the measured temperature separately.

Further-more, r and ki;j are the standard deviation of measurement errors and the random number, respectively. The

value of ki;jis calculated by the IMSL subroutine DRNNOR[39]and located within2.576–2.576, this range

represents the 99% conﬁdence bound for the measured temperature.

The time evolution in all cases is from 0 to 1.8 with 0.05 increments. As well, the increment of spatial coor-dinate is 0.1. The heat sources are applied at x = 0.5 and the temperature measurements are taken from x = 1 in all examples. The values of l, k, qC and b are set to unity in Eq. (1). Two kinds of random noise level r= 0.001 and r = 0.005 are adopted. The application of the proposed approach is demonstrated by the fol-lowing examples in this present. The strength of the heat source is presented as the time-varying function. Detailed descriptions for the examples are shown as follows:

Example 1 gðtÞ ¼ 0:3 þ7 9t when 0 < t < 0:9 and gðtÞ ¼ 1:5 þ 5 9t when 0:9 6 t 6 1:8: ð40Þ Example 2 gðtÞ ¼ sin 5pt 6 when 0 < t < 1:8: ð41Þ Example 3

(9)

The inverse solution of the problems is to identify the magnitudes of the time-dependent heat source. The feature of the hyperbolic equation is that the heat wave propagates with a ﬁnite speed. Accordingly, it is unex-pected to estimate the input source from the measured temperature immediately. Moreover, the concept of the future time is adopted to resolve the problem in order to recover the strength of the heat source from a ‘‘delay-temperature’’ measurement. The estimated results without the measurement errors are shown inFigs. 1, 4 and 7and four future time steps are used. All examples have excellent approximations when measurement errors are free. When measurement errors are included, the increasing future times are used to stabilize the estimated results in the example problems. As the random noise level r = 0.001, the results are shown inFig. 2for Exam-ple 1,Fig. 5forExample 2, andFig. 8forExample 3. A close examination of these ﬁgures discloses that they

0 0.2 0.4 0.6 0.8 1 1.2 0 0.3 0.6 0.9 1.2 1.5 1.8 Exact Estimated r = 4

Estimated heat source, g(t)

Temporal-coordinate σ = 0

Fig. 1. Estimation of the strength of the heat source inExample 1(measurement error r = 0).

0 0.2 0.4 0.6 0.8 1 1.2 0 0.3 0.6 0.9 1.2 1.5 1.8 Exact Estimated r = 6 Estimated r = 8

Temporal-coordinate σ = 0.001

(10)

have more stable outcome when r = 8. In r = 0.005, the results are shown inFig. 3forExample 1,Fig. 6for

Example 2, andFig. 9for Example 3. The results show that they have more stable estimation when r = 10. Overall inspection on these examples discloses that more future times are needed to stabilize the estimation when error level is raised.

To investigate the deviation of the estimated results from the exact solution, the relative average errors for the estimated solutions are deﬁned as follows:

l¼ 1 Nt XNt j¼1 f ^f ^ f ; ð43Þ 0 0.2 0.4 0.6 0.8 1 1.2 0 0.3 0.6 0.9 1.2 1.5 1.8 Exact r = 8 r = 10

Fig. 3. Estimation of the strength of the heat source inExample 1(measurement error r = 0.005).

-1.5 -1 -0.5 0 0.5 1 1.5 0 0.3 0.6 0.9 1.2 1.5 1.8 Exact Estimated r = 4

(11)

where f is the estimated result with measurement errors, ^f is the exact result, and Ntis the number of the

tem-poral steps. It is clear that a smaller value of l indicates a better estimation and vice versa.

When measurement errors are considered, the relative average errors of the present estimation are shown in

Table 1. It is seen that the result with the larger measurement error is less accurate than that with smaller error. For example, when r = 6 inExample 1, the value of relative errors is 0.250946 and 1.339404 when r = 0.0001 and r = 0.005 respectively. It is also interesting to investigate the relationship among the average relative error and the number of the future time. No matter what the value of r is 0.001 or 0.005, the estimation appears better as the number of future time increases. For example in all cases, when r = 0.005 inExample 3, the best estimation is r = 10 and the worst estimation is r = 6. It reveals that the proposed method is robust and stable when the measurement error is included in the estimation. The similar examples with the Fourier eﬀect are

-1.5 -1 -0.5 0 0.5 1 1.5 0 0.3 0.6 0.9 1.2 1.5 1.8 Exact Estimated r = 6 Estimated r = 8

(12)

done by Yang[3]. Yang adopts r = 0 for the inverse estimation when 0 and the exact solution can be approx-imated. However, in this research, the future time step is set to four r = 4 in that the exact solution can be approximated. In other words, the non-Fourier problem needs more future time steps to estimate the strength of the heat source than those of the Fourier model. Additionally, the relative average errors of the Yang’s results [3] are also shown in Table 2. The result shows that the large values of the relative average errors are computed in the non-Fourier heat conduction model. It represents the estimation of the non-Fourier heat conduction needs more computed strategies than that of the Fourier model.

To be concluded, the numerical results show that the exact solution of the examples can be found when future time is r = 4. Yet, this is under the condition that the measurement errors are neglected. When measure-ment errors are included, it is suggested that the r = 8 for and r = 10 for r = 0.005 are needed for the better estimations in the examples. In addition, the numerical results show that the estimation of the non-Fourier

0 0.2 0.4 0.6 0.8 1 1.2 0 0.3 0.6 0.9 1.2 1.5 1.8 Exact Estimated r = 4

Fig. 7. Estimation of the strength of the heat source inExample 3(measurement error r = 0).

(13)

heat conduction is more complicated than that of the Fourier model. That is to say, more future time steps are needed to improve the stability and accuracy of the estimated results in the inverse hyperbolic heat conduction problem.

7. Conclusions

The forward problem computed by a finite difference method within a stable interval is used to predict the inverse heat source problems with non-Fourier effect based on the sequential approach. As well, the inverse solution at each time step is solved by a modified Newton–Raphson method. The non-linear least-squares error is not adopted to formulate the inverse problem, but it is employed a direct comparison of the measured temperature and calculated temperature. Special features about this method are that no preselect functional form for the unknown function is necessary and non-linear least squares is not necessary in the algorithm. Three examples have been used to show the usage of the proposed method. The results show that the exact solution can be found through the proposed method without measurement errors. When measurement errors are increased, it is suggested that the increasing of the future times to stabilize the fluctuation of the estimation

Table 1

The relative average errors of example one and three in the non-Fourier heat conduction

Measurement error r= 0 r= 0.001 r= 0.001 r= 0.001 r= 0.005 r= 0.005 r= 0.005

r = 4 r = 4 r = 6 r = 8 r = 6 r = 8 r = 10

Relative error ofExample 1 0 63.49735 0.250946 0.015417 1.339404 0.059275 0.02817 Relative error ofExample 3 0 76.31025 0.190392 0.045507 1.536017 0.076226 0.060945

Table 2

The relative average errors of example one and three in the Fourier heat conduction[3]

Measurement error r= 0 r= 0.001 r= 0.001 r= 0.005 r= 0.005 r= 0.005

r = 1 r = 2 r = 4 r = 4 r = 6 r = 8

Relative error ofExample 1 3.44877E07 0.041932 0.010743 0.048317 0.024803

(14)

from the exact solution (i.e., r = 8 for r = 0.0001 and r = 10 for r = 0.005). It can be concluded that more future time steps are needed to improve the stability and accuracy of the estimated results, likewise this result can be refereed to adopt the number of future times in the inverse source problem in the future researches. Moreover, the proposed method is also applicable to the other kinds of inverse problems such as boundary estimation in the one- or multi-dimensional heat transfer problems.

References

[1] C.H. Huang, M.N. Ozisik, Inverse problem of determining the unknown strength of an internal plate heat source, J. Franklin Inst. 329 (1992) 751–764.

[2] C.Y. Yang, Noniterative solution of inverse heat conduction problems in one dimension, Commun. Numer. Method Eng. 13 (1997) 419–427.

[3] C.Y. Yang, A sequential method to estimate the strength of the heat source based on symbolic computation, Int. J. Heat Mass Transfer 41 (1998) 2245–2252.

[4] C.Y. Yang, The determination of two heat sources in an inverse heat conduction problem, Int. J. Heat Mass Transfer 42 (1999) 345– 356.

[5] C.Y. Yang, Solving the two dimensional inverse heat source problem through the linear least-squares error method, Int. J. Heat Mass Transfer 41 (1998) 393–398.

[6] A.J. Silva Neto, M.N. Ozisik, Inverse problem of simultaneously estimating the time-varying strengths of two plane heat sources, J. Appl. Phys. 73 (1993) 2132–2137.

[7] C. Le Niliot, F. Lefevre, A method for multiple steady line heat sources identiﬁcation in a diﬀusive system: application to an experimental 2D problem, Int. J. Heat Mass Transfer 44 (2001) 1425–1438.

[8] C. Le Niliot, F. Lefevre, Multiple transient point heat sources identiﬁcation in heat diﬀusion: application to numerical 2D and 3D problems, Numer. Heat Transfer, Part B 39 (2001) 277–301.

[9] F. Lefevre, C. Le Niliot, The BEM for point heat source estimation: application to multiple static sources and moving sources, Int. J. Therm. Sci. 41 (2002) 536–545.

[10] F. Lefevre, C. Le Niliot, Multiple transient point heat sources identiﬁcation in heat diﬀusion: application to experimental 2D problems, Int. J. Heat Mass Transfer 45 (2002) 1951–1964.

[11] R.A. Khachfe, Y. Jarny, Determination of heat sources and heat transfer coeﬃcient for two-dimensional heat ﬂow-numerical and experimental study, Int. J. Heat Mass Transfer 44 (2001) 1309–1322.

[12] W.B. Lor, H.S. Chu, Eﬀect of interface thermal resistance on heat transfer in a composite medium using the thermal wave model, Int. J. Heat Mass Transfer 43 (2000) 653–663.

[13] P.J. Antaki, Importance of nonFourier heat conduction in solid-phase reactions, Combust. Flame 112 (1998) 329–341.

[14] T. Sanderson, C. Ume, J. Jarzynski, Hyperbolic heat equations in laser generated ultrasound models, Ultrasonics 33 (1995) 415–421. [15] L.H. Liu, H.P. Tan, T.W. Tong, Non-Fourier eﬀects on transient temperature response in semitransparent medium caused by laser

pulse, Int. J. Heat Mass Transfer 44 (2001) 3335–3344.

[16] A.M. Mullis, Rapid solidification and a finite velocity for the propagation of heat, Mater. Sci. Eng. A226–A228 (1997) 28–32. [17] A.M. Mullis, Rapid solidification within the framework of a hyperbolic conduction model, Int. J. Heat Mass Transfer 40 (1997) 4085–

4094.

[18] K.S. Ranjit, R. Wilfried, Hyperbolic axial dispersion model for heat exchangers, Int. J. Heat Mass Transfer 45 (2002) 1261–1270. [19] R. Wilfried, K.D. Sarit, Hyperbolic axial dispersion model: concept and its application to a plate heat exchanger, In. J. Heat Mass

Transfer 38 (1995) 3065–3076.

[20] R. Wilfried, N.R. Chakkrit, Consideration of maldistributation in heat exchangers using the hyperbolic dispersion model, Chem. Eng. Process. 38 (1999).

[21] J.Y. Lin, The non-Fourier effect on the fin performance under periodic thermal conditions, Appl. Math. Modell. 22 (1998) 629–640. [22] A.H. Bishri, Modelling non-Fourier heat conduction with periodic thermal oscillation using the finite integral transform, Appl. Math.

Modell. 23 (1999) 899–914.

[23] D.W. Tang, N. Araki, Non-Fourier heat conduction in a ﬁnite medium under periodic surface thermal disturbance, Int. J. Heat Mass Transform 39 (1996) 1585–1590.

[24] Z.M. Tan, W.J. Yang, Propagation of thermal waves in transient heat conduction in a thin ﬁlm, J. Franklin Inst. 336B (1999) 185– 197.

[25] Z.M. Tan, W.J. Yang, Heat transfer during asymmetrical collision of thermal waves in a this ﬁlm, Int. J. Heat Mass Transfer 40 (1997) 3999–4006.

[26] C. Weber, Analysis and solution of the ill-posed inverse heat conduction problem, Int. Heat Mass Transfer 24 (1981) 783–1791. [27] A.K. Nehad, On the solution of parabolic and hyperbolic inverse heat conduction problems, Int. J. Heat Mass Transfer 41 (1998)

3731–3740.

[28] A.K. Nehad, Analysis of boundary inverse heat conduction problems using space marching with Savitzky–Gollay digital ﬁlter, Int. Commun. Heat Mass Transfer 26 (1999) 199–208.

[29] H.T. Chen, S.Y. Peng, L.C. Fang, Numerical method for hyperbolic inverse heat conduction problems, Int. Commun. Heat Mass Transfer 28 (2001) 847–856.

(15)

[30] C.Y. Yang, Direct and inverse solutions of the hyperbolic heat conduction problems, AIAA J. Thermophys. Heat Transfer 19 (2005) 217–225.

[31] C.Y. Yang, Estimation of the periodic thermal conditions on the non-Fourier ﬁn problem, Int. J. Heat Mass Transfer 48 (2005) 3506– 3515.

[32] J.V. Beck, B. Blackwell, C.R.St. Clair, Inverse Heat Conduction-Ill Posed Problem, ﬁrst ed., Wiley, New York, 1985.

[33] C.Y. Yang, Determination of temperature dependent thermophysical properties from temperature responses measured ad medium’s boundaries, Int. J. Heat Mass Transfer 43 (2000) 1261–1270.

[34] C.Y. Yang, Estimation of boundary conditions in nonlinear inverse heat conduction problems, AIAA J. Thermophys. Heat Transfer 17 (2003) 389–395.

[35] G.F. Carey, M. Tsai, Hyperbolic heat transfer with reﬂection, Numer. Heat Transfer 5 (1982) 309–327.

[36] B. Carnahan, H.A. Luther, J.O. Wilkes, Applied Numerical Methods, ﬁrst ed., John Wiley & Sons Inc, New York, 1977. [37] A.N. Tikhonov, V.Y. Arsenin, Solutions of Ill-posed Problems, Winston, Washington, DC, 1977.

[38] N. Frank, P. Wolfe, An algorithm for quadratic programming, Naval Res. Log. Quart. 3 (1956) 95–110. [39] User’s Manual: Math Library Version 1.0, IMSL Library Edition 10.0, IMSL, Houston, TX, 1987.