A space marching method for multidimensional ... - Springer Link

Heat and Mass Transfer 34 (1999) 349±356 Ó Springer-Verlag 1999

A space marching method for multidimensional transient inverse heat conduction problems J. Taler, P. Duda

Abstract A new method for the solution of multidimensional heat conduction problems is formulated. The developed space marching method allows to determine quickly and exactly unsteady temperature distributions in the construction elements of irregular geometry. The method is especially appropriate for determining transient temperature distribution in thick-wall pressure components based on temperature measurements at the outer surface. Two examples are included to demonstrate the capabilities of the new approach. List of symbols Bi1 ˆ h1 b=k; Bi2 ˆ h2 d=k Biot number b length of the plate, m c heat capacity, J/kgK d width of the plate, m d diameter, m ds surface vector f measured value of temperature at time t,  C 2 Fo1 ˆ at=b ; Fo2 ˆ at=d2 Fourier numbers g header thickness, m h heat transfer coef®cient, W/m2 K k thermal conductivity, W/mK n unit outward normal to boundary surface N number of control volumes q heat ¯ux value, W/m2 q heat ¯ux vector r radius, m r in inner radius, m r out outer radius, m S surface of control volume, m2 t time, s T temperature,  C T0 initial temperature,  C T1 ambient temperature,  C V volume of ®nite element, m3 x; y cartesian coordinates, m

Received on 28 October 1998

Prof. dr hab. in_z Jan Taler Dr in_z. Piotr Duda Cracow University of Technology Institute of Process- and Power Engineering PL-31-864 KrakoÂw, ul. Jana Pawla II 37, Poland

…x1 ; y1 †; …x2 ; y2 †; …x3 ; y3 † nodal coordinates of x and y for triangular element, m Greek symbols a thermal diffusivity, m2 =s bi ; c i roots of characteristic equations Dt time step, s Dx; Dy spatial sizes of control volumes, m H dimensionless temperature q density, kg/m3 rf standard deviation of the measurement errors, K Subscripts & Superscripts c calculated ex exact i spatial node index in inner o centre of ®nite element out outer mean

Introduction Space marching methods are widely applied in practice [1±6] which is justi®ed by simple calculation algorithms usually based on ®nite difference method. These methods can only be applied to one dimensional inverse problems i.e. only for bodies of simple shapes. Another signi®cant disadvantage of space marching methods is their high sensibility to random measuring errors in input data. It is caused by strong propagation of error (bound to the temperature measurement) by marching in space from the temperature measurement point to the active surface. A new space marching method for solving the inverse problem of transient heat transfer is presented. The method is appropriate for complex shaped bodies and gives results of high stability. The temperature dependent material properties (k; c and q) are also accounted for. Formulation of the method The equation governing the transient heat conduction problem is given by c…T†q…T†

oT ˆ ÿrq ; ot

…1†

where q is the heat ¯ux vector de®ned by Fourier's law

349

integral formulation corresponding to Eq. (5) can be arrived at by applying the conservation principle for energy The material thermal properties: speci®c heat c, thermal to a control volume V 1aoc , that is ®xed in space. conductivity k and density q are known functions of The resulting integral conservation equation (4), when temperature. The numerical solution of the heat conducapplied to the polygonal control volume surrounding node tion equation will be constructed using a control-volume1 in Fig. 2, can be written as follows: based ®nite-element method. Control-volume-based "Z # Z c 0 methods retrain the full geometric ¯exibility inherent in ÿ q
n ds ‡ q
n ds ®nite-element methods and are receiving continually …6† a 0 increasing usage in the ®nite-element analysis [7]. Consider the control volume V with a bounding surface ˆ Vlaoc
q…T1 †
c…T1 †
dT1 =dt S. The variable-property heat conduction equation may be where n is a unit outward vector normal to the differential integrated over the volume V, to obtain Z Z area element ds. oT The total contribution of element 123 to the conservaq…T†c…T† …3† dV ˆ ÿ rq dV : ot tion equation for the control volume surrounding node 1, V V Using the mean value theorem for integrals on the left and obtain from equation (6) (see Appendix), can be compactly expressed in the following form: the divergence theorem on the right gives

q ˆ ÿk…T† rT :

350

Vq…T†c…T†

dT ˆÿ dt

…2†

Z

S

q
n dS ;

…4†

where the bar indicates a suitable average value over the region. For clarity, the ®nite-element equations are derived in the Cartesian coordinate system in two dimensions. The governing heat conduction equation (1) becomes

    oT o oT o oT ˆ k…T† ‡ k…T† : q…T†c…T† ot ox ox oy oy

2

b c

…5†

o

3

a The domain of interests is ®rst divided into three±node triangular elements. Then the centroids of the elements are joined to the midpoints of the corresponding sides: this 1 creates polygonal control volumes around each node in the calculation domain. A sample domain discretization is shown in Fig. 1. The solid lines denote the domain and element boundaries, the dashed lines represent the control-volume faces, and the shaded areas show the control volumes associated with one internal node and two Fig. 2. Details of the discretization in Fig. 1 and related boundary nodes. In this discretization scheme, curved boundaries are approximated by piecewise-linear curves. nomenclature; an internal node Consider a typical node i in the calculation domain: it could be an internal node, like the one shown in Fig. 2. An 2

b

noc β2

y

β1

c

x 0

α2

noa α1

3

α2 a

Fig. 1. An irregularly shaped calculation domain and its discretization into three-node triangular elements and polygonal volumes

1

Fig. 3. A typical three-node triangular element and related nomenclature; local x, y coordinate system

1  k…Ta †‰…y2 ÿ y3 †T1 ‡ …y3 ÿ y1 †T2 ‡ …y1 ÿ y2 †T3 Šya DET

The temperature distribution in the area V I including the surface temperature S will be evaluated by extrapolation of the solution obtained in the area V D onto the area V I , retaining the condition of equal temperatures and heat ¯ux density on the surface Sm (Fig. 4a). Problem formulated in this way is an inverse, ill-posed heat transfer problem which is very dif®cult to solve correctly. It is caused by extremely high sensitivity of the calculated temperature or the heat ¯ux density on the surface Sin to the random measurement errors connected with time-temperature histories measured on the surface Sm . Even very small temperature measurement errors may cause high oscillations in the parameters being determined [8]. The solution of the inverse problem starts with writing the heat balance equations for the control volumes DN ‡ 1; . . . ; D2N, whose central nodes are located on the surface Sm (Fig. 4b). Because the temperature histories in these nodes are known, the temperature in adjacent nodes 1; . . . ; N in the inverse region V I is evaluated from the heat balance equations (Fig. 4b). Then the heat balance equations are written again for the nodes in the inverse region V I that the temperature was already determined and the temperature in adjacent nodes, closer to the surface Sin is evaluated from these equations. The above procedure is repeated until the temperature in the nodes on the surface Sin is determined. Knowing the temperature distribution in the inverse region V I , the heat ¯ux density on the surface Sin can also be calculated. The method is space marching because of moving in space from the surface Sm where the temperature and heat

ÿ k…Ta †‰…x3 ÿ x2 †T1 ‡ …x1 ÿ x3 †T2 ‡ …x2 ÿ x1 †T3 Šxa ÿ k…Tc †‰…y2 ÿ y3 †T1 ‡ …y3 ÿ y1 †T2 ‡ …y1 ÿ y2 †T3 Šyc ‡ k…Tc †‰…x3 ÿ x2 †T1 ‡ …x1 ÿ x3 †T2 ‡ …x2 ÿ x1 †T3 Šxc ˆ ÿc…T1 †q…T1 †V1aoc

dT1 dt



…7†

where: V1aoc is the ®eld of ®gure laoc. Writing the heat balance equations for all control volumes (Fig. 1) the system of N non-linear ordinary equations is obtained, whose solution are the temperature changes in the central nodes of the control volumes. The number of equations is equal to the number of the control volumes N. The initial condition gives the starting values of nodal temperatures. In this way the initialboundary problem for the heat transfer equation is transformed into the initial problem for N ordinary differential equations. For the direct problem these equations can be easily integrated using one of many available methods for this class of problems, for instance the Runge-Kutta method.

Inverse heat conduction problem Solving the inverse problem is far more dif®cult. This paper presents the new method of solving the inverse transient heat transfer problems of the space-marching kind.

Sout

D2

D1

VD

DN+1

DN+2

DN D2N

Sm S in

D2N-1 N-1

DN+3 3

N+1 N+2 N+3

2N DN-1

VI

2

1

N

Known boundary condition

D3

4

N+4

2N-1

DN+4 D4

S in

Temperature measurement point

N+i i

Unknown boundary condition

a

S in -

Inner surface

Sout -

Outer surface

VD -

Direct region

VI

Inverse region

-

DN+i Di

VD

Sm VI b

into direct V D and inverse V I regions, b triangular ®nite element Fig. 4a, b. Two-dimensional inverse heat conduction problem; Sm ± surface with known boundary conditions at which temper- mesh (solid line) with polygonal ®nite control volumes (dashed ature measurement points are located. a Division of the domain V line)

351

352

¯ux density is known towards the surface Sin where the dT5 c…T5 †
q…T5 †
Dx
Dy boundary condition is unknown (Fig. 4b). dt Detailed way of proceeding is shown using the example k…T5 † ‡ k…T4 † T4 ÿ T5 of the rectangular plate (Fig. 5a). The plate is divided into Dy ˆ Dx 2 the uniform rectangular mesh. In two dimensional probk…T5 † ‡ k…T1 † T1 ÿ T5 lems the temperature distribution must be known on the Dx ‡ surface Sm (Fig. 5a). In this case it is assumed that the Dy 2 temperatures f 9 ; . . . ; f 13 on the back surface are given k…T5 † ‡ k…T6 † T6 ÿ T5 (Fig. 5a). Starting the calculations from the nodes on this Dy ‡ Dx 2 surface and marching in space towards the surface where the boundary conditions are unknown, the inverse probk…T5 † ‡ k…T11 † T11 ÿ T5 Dx ; ‡ lem is solved. Dy 2 The heat balance equations for nodes 9±13 are of the from which the temperature T 1 is calculated: following form:

2
c…T5 †
q…T5 †
…Dy†2 dT5 …n† dt k…T5 † ‡ k…T1 † " …Dy†2 k…T5 † ‡ k…T4 † …T5 ÿ T4 † ‡ …Dx†2 k…T5 † ‡ k…T1…n† † # k…T5 † ‡ k…T6 † …T5 ÿ T6 † ‡ …n† k…T5 † ‡ k…T1 †

…n‡1†

Dy dfi k…fi † ‡ k…fiÿ1 † fiÿ1 ÿ fi Dy c…fi †q…fi † Dx ˆ 2 2 dt Dx 2 k…fi † ‡ k…Tiÿ8 † Tiÿ8 ÿ fi Dx ‡ Dy 2

T1

ˆ T5 ‡

k…fi † ‡ k…fi‡1 † fi‡1 ÿ fi Dy Dx 2 2 …8† ÿ qi Dx ‡

‡

where: i ˆ 9; . . . ; 13. From these equations temperatures T 4 , T 5 and T 6 can be calculated in the iterative way: …n‡1†

Ti

2Dy …n† k…Ti †

‡ k…fi‡8 †

;

…n†

k…T5 † ‡ k…T1 †

…T5 ÿ T11 †; n ˆ 0; 1; 2; . . . …12†

q2 ˆ c…T2 †
q…T2 †

Dy dT2 Dy ÿ 2 dt 4
…Dx†2

 f‰k…T1 † ‡ k…T2 †Š…T1 ÿ T2 † ‡ ‰k…T3 † ‡ k…T2 †Š…T3 ÿ T2 †g k…T6 † ‡ k…T2 † …T6 ÿ T2 † ÿ 2
Dy

i ˆ 4; 5; 6 n ˆ 0; 1; . . . …9†

…0† Ti

k…T5 † ‡ k…T11 †

After having determined the temperature distribution in the inverse region, the heat ¯ux density in nodes on the body surface can be calculated. Writing the heat balance equation for point no. 2 the equation for q2 is obtained:

c…fi‡8 †
q…fi‡8 †
…Dy†2 dfi‡8 …n† dt k…Ti † ‡ k…fi‡8 † " …Dy†2 k…fi‡7 † ‡ k…fi‡8 † …fi‡8 ÿ fi‡7 † ‡ 2…Dx†2 k…Ti…n† † ‡ k…fi‡8 † # k…fi‡9 † ‡ k…fi‡8 † …fi‡8 ÿ fi‡9 † ‡ …n† k…Ti † ‡ k…fi‡8 †

ˆ fi‡8 ‡

‡ qi‡8

…11†

…13†

Additionally the heat transfer coef®cient for the point 2 can be calculated:

h2 …t† ˆ

q…t† ; T1 ÿ T2 …0; t†

…14†

Assuming for examples ˆ f i‡8 …i ˆ 4; . . . ; 6† the iterative procedure is continued until the following condi- where T1 is the temperature of ¯uid ¯owing along the tion is ful®lled: surface 1±3 …y ˆ d†. The precision of calculating the time derivatives in …n‡1† equations (9), (12) and (13) has a strong in¯uence on the …n† ÿ Ti  e; n ˆ 0; 1; 2; 3 . . . ; i ˆ 4; . . . ; 6 Ti quality of obtained results. The Gram polynomials are …10† used to smooth the data. The value of smoothing Gram polynomial and its dewhere e is given tolerance of calculations. rivative are calculated only at the center of the approxiMeasuring time-temperature changes in …2N ‡ 1† mation interval, while for determining of polynomial points on the surface 9±15, the temperature changes in coef®cients future and past data points are used. Using the …2N ÿ 1) nodes placed inside the inverse region and estimated Gram polynomial a moving averaging ®lter is adjacent to the surface 9±15 can be calculated. constructed. To ®t a surface f …y; t† through data on rectThen the heat balance equation for the point no. 5 can angular grid …yi ; tj ; fij † two-way univariate approximation be written: is used. First the ®ltering is applied to all the data, keeping

t j ®xed and varying i, and then by applying it again to the resulting sets of spatially smoothed data points, but this time keeping yi ®xed and varying j. After having determined the temperatures T i in adjacent nodes using the formula (12), the temperatures are smoothed again. The time derivatives are calculated after the temperatures Ti have been smoothed. The process of smoothing is repeated each time after the temperatures on the adjacent surface are calculated. In the case shown in Fig. 5b, ®rst the temperatures measured in nodes …28  35† are smoothed. The smoothing procedure is then repeated to the node sets: 19  26 and 10  16.

solution of the direct problem, the material properties were assumed to be temperature independent:

k ˆ 50 W=mK; c
q ˆ 3768000 J=m3 K : The ``measured data'' for nodes 9±15 are obtained using the superposition method [9]. The calculated temperature distribution on the lower surface with the time interval of Dt ˆ 3 s was taken as an undisturbed ``measured data'' for the inverse method. The ``disturbed measured data'' are obtained by adding the ``measuring errors'' of normal distribution. The average value of the generated error was equal to zero and the standard deviation rf ˆ 1=5 (error  0, 6 K). Figure 6 shows the ``disturbed measured data'' generApplication examples ated for nodes 9±15. The presented space-marching method will ®rst be applied The inverse method calculations were carried out using to solving the inverse heat conduction problem with the control volumes shown in Fig. 5a and 5b. In the case of temperature independent properties. ®ner division (Fig. 5b), both rectangular and polygonal control volumes were used. The results obtained for difTwo dimensional plate ferent rectangular meshes of control volumes are shown in This example shows the calculation of the temperature distribution, heat ¯ux density, and heat transfer coef®cient Fig. 7a. The estimated temperatures of the active upper surface (Fig. 7b) were obtained using polygonal control on the upper surface of the plate (Fig. 5) based on the volumes. Comparing them with the analytical results it can ``measured data'' obtained on the lower surface of the plate. The ``measured data'' were generated from the so- be seen that the presented method is very ef®cient. Adding lution of the direct problem. The plate of initial temper- one layer of elements results in slight improvement of ature T 0 ˆ 20 C was subjected to the step change of the accuracy. It should be noticed, that the error at deter¯uid temperature T 1 ˆ 100 C while heat transfer coef®- mining the temperature (Fig. 7b) is smaller than at cients were: h1 ˆ 1071; 43 W=m2 K and h2 ˆ 500 W=m2 K. determining the heat ¯ux or the heat transfer coef®cient The remaining two sides of the plate are thermally insu- (Fig. 7a). lated. Because the data were generated from the analytical Experimental verification of the presented method High thermal stresses may occur during start-up and shutdown in thick-walled power boiler elements working at high temperatures and pressures due to the non-uniform temperature distribution. In order to determine stresses at the element inner surface adjacent to the highly pressurized working ¯uid it is necessary to know the temperature of this surface, where placing of temperature sensors is very dif®cult. The proposed method allows determining the timespace temperature distribution on this hardly accessible surface based on the temperature measurements on the elements outer surface. Knowing the measured and cal-

Temperature T (°C)

100

60 T9 40 20 0

Fig. 5a, b. Control volume mesh. a Temperature is measured at nodes 9  15, b temperature is measured at nodes 28  36

T15

80

0

36

72 108 Time t (s)

144

180

Fig. 6. Noisy temperature data from the direct calculation used to simulate thermocouple measurements at the nodes 9 to 15 (Fig. 5a) normally distributed random errors between 0:6 K were used to perturbate exact data

353

Heat transfer coefficient h (W/m2K)

∆r

600

21

27 20

14

400

13 7

Exact 3 layers 4 layers

200

3

39

75

111 Time t (s)

r1 1 r6

60

11

r5

r4

r3

r2

18 25

3 10

2

17 9

24

8 16 15

40

23

T1-exact T1-inv.m. T5-exact T5-inv.m.

20 0

5 4

b

80

12

180

147

19 26

6 ∆ϕ

0 r0

Temperature T (°C)

354

28

a

0

36

72

108 Time t (s)

144

22

180

Fig. 8. Cross section of the spray attemperator with control volumes; temperature was measured at nodes 22±28

where temperature T is expressed in  C. The Ni-CrNi sheath thermocouples of outer diameter d ˆ 3 mm were attached on the outer surface in points 22± 28. Measured temperatures and pressure inside the header are shown in Fig. 9a. Using the measured temperature and applying the nonlinear space marching method, the temperatures in all the remaining points 1±21 were calculated. Temperature dependent material properties were accounted for. Fig. 9b shows the difference between the calculated temperatures culated temperatures, thermal stresses can be evaluated. in points 5, 6 and 7 on the inner surface and the measured Furthermore, it is possible to determine the heat ¯ux temperature on the outer surface in points 26, 27 and 28, density and the heat transfer coef®cient on the inner respectively. Knowing the transient temperature distribusurface of the element. To verify the developed method in practice and to test tion, thermal stresses can be calculated. the technique of smoothing the temperature changes, Final Remarks measurements were made in a power plant. Critical elements determining the speed of temperature The presented new space marching method for solving the changes during boiler start-up and shut-down are steam inverse transient heat transfer problems can be applied for boiler drums, superheater headers, and headers of spray complex-shaped bodies with temperature-dependent theattemperators. High thermal stresses occur in thick-walled rmophysical properties. The method gives stable solutions also for random disturbed temperature measurement erattemperator, when injected water does not evaporate rors. The calculation algorithm is simple and easy to code. completely. Apart from this it is working under creep Thanks to its ¯exibility, solving the inverse problems can conditions because of the high steam temperature. The be automated, similarly to the direct problems solved above mentioned factors decide of low durability of the attemperator. This is why this element was used to verify using the Finite Element Method. The applicability of the inverse space marching method was demonstrated the presented method. through comparison with the results obtained from the The measurements were made on the half of the attdirect solution. The examples presented in the paper emperator perimeter because of symmetry. demonstrate high ef®ciency and ¯exibility of the method. The attemperator is perfectly insulated on the outer surface (Fig. 8). It is made of the 10CrMo910 steel and its dimensions are: d0 ˆ 0:355 m; and thickness g ˆ 0:07 m. Appendix The following thermophysical properties were assumed The derivation of algebraic approximations to the integrals in Eq. (6) requires the speci®cation of interpolation [10]: functions for the dependent variable T, conductivity k, density q and speci®c heat c. The scalar dependent variable 1 ÿ2 ÿ5 2 k ˆ 3:53
10 ‡ 2:14
10 T ÿ 4:3
10 T ‰W=mKŠ T is interpolated linearly:

Fig. 7a, b. Comparison between exact input data and inverse results. a Comparison between exact and estimated heat transfer coef®cient at node 2 (Fig. 5a) for two different rectangular control volume meshes (Fig. 5a±3 layers, Fig. 5b±4 layers), b comparison between exact and estimated surface temperature at nodes 1 and 5 for polygonal control volumes shown in Fig. 5b; Ðб analytical direct solution, e, n-inverse analysis at node 1 and 5, respectively

a ˆ 10:2
10ÿ6 ÿ 0:48
10ÿ8 T ÿ 5:06
10ÿ12 T 2

‰m2 =sŠ T ˆ A
x ‡ B
y ‡ C :

…A1†

16

Temperature (°C)

12

23

440

22 24

390

10 8

26

340

0

4

28

a

290

4800

9600

6

p

27

Pressure p (MPa)

14

490

25 14400

19200

24000 28800 Time t (s)

33600

38400

43200

2 47820

355

40

Temperature T (°C)

20 0 -20 -40

5 T5-T26 6 T6-T27 7 T7-T28

-60 -80 b

-100 0

4800

9600

14400

19200

24000 28800 Time t (s)

33600

38400

43200

Fig. 9a, b. Measured temperatures and calculated temperature differences for the spray attemperator. a Temperature history at different thermocouple locations (nodes) at the outside surface of the attemperator during shut-down of the boiler (Fig. 8), b estimated temperature differences (between node pairs 5±26, 6±27 and 7±28) across the wall of the attemperator during shut-down of the boiler

47820

In each element, the constants: A, B and C in this inter- The normal vectors noc and noa needed to evaluate the polation function can be uniquely determined in terms of integrals in the equations (6) are given by the x; y coordinates of the three nodes and the correnoa ˆ cos…a1†
i ‡ cos…a2†
j sponding values of T. Thus with reference to the element ÿya xa 123 and local x; y coordinate system shown in Fig. 3, where ˆ p
i ‡ p
j ; 2 2 2 i ˆ 1, we obtain xa ‡ ya xa ‡ y2a

  A ˆ …y2 ÿ y3 †T1 ‡ …y3 ÿ y1 †T2 ‡ …y1 ÿ y2 †T3 =DET;   B ˆ …x3 ÿ x2 †T1 ‡ …x1 ÿ x3 †T2 ‡ …x2 ÿ x1 †T3 =DET;  C ˆ …x2 y3 ÿ x3 y2 †T1 ‡ …x3 y1 ÿ x1 y3 †T2  ‡ …x1 y2 ÿ x2 y1 †T3 =DET : …A2†

where: DET ˆ …x1 y2 ‡ x2 y3 ‡ x3 y1 ÿ x2 y1 ÿ x3 y2 ÿ x1 y3 †. The derivation of Eq. (6) requires the speci®cation of the heat ¯ux q. In each element, q can be expressed in terms of its components in the x and y directions:

    oT oT q ˆ qx i ‡ qy j ˆ ÿk…T† i ‡ ÿk…T† j ; ox oy

…A3†

noc ˆ cos…b1†
i ‡ cos…b2†
j yc ÿxc ˆ p
i ‡ p
j ; x2c ‡ y2c x2c ‡ y2c Using Eq. (A4) and with reference to element 123 and the local x; y coordinate system in Fig. 3, the integrals in Eq. (6) that represent the diffusion transport of energy can be approximated as follows:

Z

‰ÿya ; xa Š q q
n ds ˆ ‰qx ; qy Š
p x2a ‡ y2a ; x2a ‡ y2a a Z c ‰yc ; ÿxc Š q q
n ds ˆ ‰qx ; qy Š
p x2c ‡ y2c : x2c ‡ y2c 0 0

…A5†

where: i, j are unit vectors in the x and y directions, re- Substituting Eq. (A3) into Eq. (A4) yields: spectively. The interpolation function given in Eq. (A1) is Z 0 used to approximate q: 

q ˆ ÿ ‰k…T†……y2 ÿ y3 †T1 ‡ …y3 ÿ y1 †T2 ‡ …y1 ÿ y2 †T3 †=DETŠi ÿ ‰k…T†……x3 ÿ x2 †T1 ‡ …x1 ÿ x3 †T2 ‡ …x2 ÿ x1 †T3 †=DETŠj :

a

…A4†

q
n ds ˆfk…Ta †‰…y2 ÿ y3 †T1 ‡ …y3 ÿ y1 †T2  a †‰…x3 ÿ x2 †T1 ‡ …y1 ÿ y2 †T3 Šya ÿ k…T 1 ; ‡ …x1 ÿ x3 †T2 ‡ …x2 ÿ x1 †T3 Šxa g DET …A6†

Z

c 0

 c †‰…y2 ÿ y3 †T1 ‡ …y3 ÿ y1 †T2 q
n ds ˆfk…T  c †‰…x3 ÿ x2 †T1 ‡ …y1 ÿ y2 †T3 Šyc ÿ k…T 1 ‡ …x1 ÿ x3 †T2 ‡ …x2 ÿ x1 †T3 Šxc g ; DET …A7†

where

 a †  k…T0 † ‡ k…Ta † ; k…T 2

356 and

 c †  k…T0 † ‡ k…Tc † ; k…T 2 Substituting Eqs. (A6) and (A7) into Eq. (6) the total contribution of element 123 to the conservation equation for the control volume surrounding node 1 is obtained (Eq. 7). References

[1] D'Souza N (1975) Numerical solution of one-dimensional inverse transient heat conduction by ®nite difference method. ASME Paper No. 75-WA/HT-81 [2] Weber Ch F (1981) Analysis and solution of the ill-posed inverse heat conduction problem. Int. J. Heat Mass Transfer 24: 1783±1792 [3] Raynaud M; Bransier J (1986) A new ®nite difference method for nonlinear inverse heat conduction problem. Num. Heat Transfer 9: 27±42 [4] Hensel E H; Hills R G (1986) An initial value approach to the inverse heat conduction problem. ASME Journal of Heat Transfer, 108: 248±256 [5] Vogel J; Sara L; Krejci L (1993) A simple inverse heat conduction method with optimization. Int. J. Heat Mass Transfer 36: 4215±4220 [6] Taler J (1996) A semi-numerical method for solving inverse heat conduction problems. Heat and Mass Transfer 31: 105±111 [7] Baliga B R; Patankar S V Elliptic Systems: Finite Element Method II. Handbook of Numerical Heat Transfer, Wiley 1988 New York [8] Taler J (1995) Theory and Practice of Identi®cation of Heat Transfer Processes. Wroclaw, Ossolineum [9] Taler J; Zima W Solution of inverse heat conduction problems using control volume approach. Int. J. Heat and Mass Transfer (in Press) [10] Richter F (1983) Physikalische Eigenschaften von StaÈhlen und ihre TemperaturabhaÈngigkeit. Mannesmann Forchungsberichte 930/1983, Sonderdruck, Stahleisen-Sonderberichte Heft 10, 1983, Verlag Stahleisen DuÈsseldorf