University of Hail

Faculty of Engineering
DEPARTMENT OF MECHANICAL ENGINEERING

ME 315 Heat Transfer

Lecture notes

Chapter 4
Part 2
Numerical analysis of conduction

Prepared by : Dr. N. Ait Messaoudene

Based on:
Introduction to Heat Transfer
Incropera, DeWitt, Bergman, and Lavine,
5th Edition, John Willey and Sons, 2007.
2nd semester 2011-2012

The Finite-Difference Method

An approximate method for determining temperatures at discrete
(nodal) points of the physical system.
Procedure:
Represent the physical system by a nodal network (discrete number of points
representing the center of regular regions). Each node represents a certain region,
and its temperature is a measure of the average temperature of the region.
The selection of nodal points is rarely arbitrary, depending often on matters
such as geometric convenience and the desired accuracy.
The numerical accuracy of the calculations depends strongly on the number of
designated nodal points. If this number is large (a fine mesh), accurate solutions
can be obtained.

Use the energy balance method to obtain a finite-difference equation for each
node of unknown temperature.
Solve the resulting set of algebraic equations for the unknown nodal temperatures.
Special treatment for points on the boundaries where conditions are prescribed

The Nodal Network

The nodal network identifies discrete
points at which the temperature is to be
determined and uses an m,n notation to
designate their location.
The aggregate of points is termed a nodal
network, grid, or mesh.

A finite-difference approximation
is used to represent temperature
gradients in the domain.

Finite-Difference Form of the Heat Equation

the value of the second derivative in the x direction at the (m, n) nodal point may be
approximated as

Use the
expressions
in the figure

Similarly:

If a regular grid is used (x= y), substitution of the above expressions in the 2D
steady state heat equation (with no heat generation) yields:

Or simply that the temperature of an interior node be equal to the average of the
temperatures of the four neighboring nodes.

The Energy Balance Method

As the direction of heat flow is often unknown, assume all heat flows are into the
nodal region of interest, and express all heat rates accordingly. Hence, the energy
balance (with heat generation) becomes:
(4.30)

Consider application to an interior nodal point (one that exchanges heat by

conduction with four, equidistant nodal points):

Per unit length in the z direction

where, for example,

(4.31)

The assumption of heat flow into (m, n) requires assuming Tm.n less then all 4 adjoining
temperatures so that it is substracted from these temperatures in the 4 expressions of heat rates
to (m, n).
In fact, the final values of temperatures will insure the correct direction of heat flow since it is
impossible for all 4 to be into (m, n) in steady state (unless we have a heat sink).

The remaining conduction

rates may be expressed as

Substituting the heat rates into the energy balance and if x = y, it follows that the finitedifference equation for an interior node (for a regular square grid) with generation is:

4.35

Which is the same as Equ. 4.29 when there is no generation.

So what is the interest or the advantage of this method?
In fact, this method is most powerful in more difficult cases as in:
Boundaries of irregular shapes and mixed conditions,
Problems involving multiple materials,
embedded heat sources.

To illustrate this, consider the node corresponding to an internal corner that exchanges
energy by convection with an adjoining fluid at T.
The 4 conduction heat rates may be expressed as

The total convection rate (also supposed into (m, n)may be expressed as:
Implicit in this expression is the assumption that the
exposed surfaces of the corner are at a uniform
temperature corresponding to the nodal temperature Tm,n.

In the case of no generation, steady-state and regular square grid, the heat balance yields
(sum of all rates=0):

 A summary of finite-difference equations for common nodal regions is provided

in Table 4.2.
Consider an external
corner with convection
heat transfer (Case 4).

q m 1, n m, n q m, n 1 m, n q m, n 0
y Tm 1, n Tm, n
x Tm, n 1 Tm, n
k

k
x
y
2
2
x
y
h
T Tm, n h
T Tm, n 0
2
2

or, with x ,y
Tm 1, n Tm, n 1 2

hx
hx
T 2
1 Tm, n 0
k
k

(4.43)

where

If the energy balance is

satisfied, the left-hand
side of this equation will
be identically
equal to zero. Substituting
values, we obtain

The inability to precisely satisfy

the energy balance can be
attributable to temperature
measurement errors, the
approximations employed in
developing the finite-difference
equations, and the use of a
relatively coarse mesh.

 Note potential utility of using thermal resistance concepts to express rate

equations. E.g., conduction between adjoining dissimilar materials with
an interfacial contact resistance.

thermal contact resistance

(in m2.K/W)

(4.46)

Solutions Methods
Matrix Inversion: Expression of system of N finite-difference equations for
N unknown nodal temperatures as:

AT C
Coefficient
Matrix (NxN)

Solution

Solution Vector
(T1,T2, TN)

(4.48)
Right-hand Side Vector of Constants
(C1,C2CN)

T A1 C
(4.49)

Inverse of Coefficient Matrix

Gauss-Seidel Iteration: Each finite-difference equation is written in explicit

form, such that its unknown nodal temperature appears alone on the lefthand side:
N aij
C i 1 aij k
k
Ti i Tj
Tj( k 1)
(4.51)
aii

j 1 aii

j i 1 aii

where i =1, 2,, N and k is the level of iteration.

Iteration proceeds starting from initial assumed values for each temperature
Ti until a prescribed convergence criterion is satisfied for all nodes:
k
k 1
Ti Ti

Acceptable error in
temperature values

Verifying the Accuracy of the Solution

It is good practice to verify that a numerical solution has been correctly formulated by
performing an energy balance on a control surface surrounding all nodal regions
whose temperatures have been evaluated.
The temperatures should be substituted into the energy balance equation, and if the
balance is not satisfied to a high degree of precision, the finite difference equations
should be checked for errors. If no errors are found, the grid should be refined.

Node 5 : T3 T6 T7 500 4T5 0

EXAMPLE 4.4
A major objective in advancing gas turbine engine technologies is to increase the temperature limit
associated with operation of the gas turbine blades. This limit determines the permissible turbine gas inlet
temperature, which, in turn, strongly influences overall system performance. In addition to fabricating
turbine blades from special, high-temperature, high-strength superalloys, it is common to use internal
cooling by machining flow channels within the blades and routing air through the channels. We wish to
assess the effect of such a scheme by approximating the blade as a rectangular solid in which rectangular
channels are machined. The blade, which has a thermal conductivity of k = 25 W/m.K, is 6 mm thick, and
each channel has a 2mm x 6 mm rectangular cross section, with a 4-mm spacing between adjoining
channels.

Problem 4.39: Finite-difference equations for (a) nodal point on a diagonal

surface and (b) tip of a cutting tool.

(a) Diagonal surface

(b) Cutting tool.

Schematic:

ASSUMPTIONS: (1) Steady-state, 2-D conduction, (2) Constant properties

ANALYSIS: (a) The control volume about node m,n is triangular with sides x and y and diagonal
(surface) of length 2 x.
The heat rates associated with the control volume are due to conduction, q1 and q2, and to convection,
qc. An energy balance for a unit depth normal to the page yields
Ein 0
q1 q2 qc 0
Tm,n-1 Tm,n
Tm+1,n Tm,n
k x 1
k y 1
h
y
x

2 x 1 T Tm,n 0.

With x = y, it follows that

hx
hx

T 2 2
Tm,n 0.
k
k

(b) The control volume about node m,n is triangular with sides x/2 and y/2 and a lower diagonal
surface of length 2 x/2 .
Tm,n-1 Tm+1,n 2

The heat rates associated with the control volume are due to the uniform heat flux, qa, conduction, qb,
and convection qc. An energy balance for a unit depth yields

Ein =0
qa q b qc 0
x
x
y Tm+1,n Tm,n

qo 1 k 1
h 2
T Tm,n 0.
x
2
2
2

or, with x = y,

Tm+1,n 2

hx
x
hx
T qo
1 2
Tm,n 0.
k
k
k

Problem 4.76: Analysis of cold plate used to thermally control IBM multi-chip,
thermal conduction module.
Features:
Heat dissipated in the chips is transferred
by conduction through spring-loaded
aluminum pistons to an aluminum cold
plate.

Nominal operating conditions may be

assumed to provide a uniformly
5
2
distributed heat flux of qo 10 W/m
at the base of the cold plate.
Heat is transferred from the cold
plate by water flowing through
channels in the cold plate.

Find: (a) Cold plate temperature distribution

for the prescribed conditions. (b) Options
for operating at larger power levels while
remaining within a maximum cold plate
temperature of 40C.

Schematic:

ASSUMPTIONS: (1) Steady-state conditions, (2) Two-dimensional conduction, (3) Constant properties.

ANALYSIS: Finite-difference equations must be obtained for each of the 28 nodes. Applying the energy
balance method to regions 1 and 5, which are similar, it follows that

x T2 x y T6 y x x y T1 0
Node 5: y x T4 x y T10 y x x y T5 0
Node 1:

Nodal regions 2, 3 and 4 are similar, and the energy balance method yields a finite-difference equation of
the form
Nodes 2,3,4:

x Tm 1, n Tm 1, n 2 x y Tm, n 1 2 y x x y Tm, n 0

Energy balances applied to the remaining combinations of similar nodes yield the following
finite-difference equations.

Nodes 6, 14:

Nodes 7, 15:

y T1 y x T7 x y y x hx k T6 hx k T

x
x

y T19 y x T15 x y y x hx k T14 hx k T

x T6 T8 2 x y T2 2 y x x y hx k T7 2hx k T

x T14 T16 2 x y T20 2 y x x y hx k T15 2hx k T

Nodes 8, 16:

Node 11:

x T7 2 y x T9 x y T11 2 x y T3 3 y x 3 x y

x T15 2 y x T17 x y T11 2 x y T21 3 y x 3 x y

h k x y T8 h k x y T

h k x y T16 h k x y T

y T8 x y T16 2 y x T12 2 x y y x hy k T11 2hy k T

Nodes 9, 12, 17, 20, 21, 22:

x Tm 1, n y x Tm 1, n x y Tm, n 1 x y Tm, n 1 2 x y y x Tm, n 0

Nodes 10, 13, 18, 23:

y Tn 1, m x y Tn 1, m 2 y x Tm 1, n 2 x y y x Tm, n 0

Node 19:
Nodes 24, 28:

x
x

Nodes 25, 26, 27:

y T14 x y T24 2 y x T20 2 x y y x T19 0

y T19 y x T25 x y y x T24 qo x k

y T23 y x T27 x y y x T28 qo x k

x Tm 1, n y x Tm 1, n 2 x y Tm, n 1 2 x y y x Tm, n 2qo x k

Evaluating the coefficients and solving the equations simultaneously, the steady-state temperature
distribution (C), tabulated according to the node locations, is:
23.77
23.41

23.91
23.62

28.90
30.72
32.77

28.76
30.67
32.74

24.27
24.31
25.70
28.26
30.57
32.69

24.61
24.89
26.18
28.32
30.53
32.66

24.74
25.07
26.33
28.35
30.52
32.65

(b) For the prescribed conditions, the maximum allowable temperature (T24 = 40C) is reached when

q o = 1.407 105 W/m2 (14.07 W/cm2).

Options for extending this limit could include use of a copper cold plate (k 400 W/mK) and/or
increasing the convection coefficient associated with the coolant.

With k = 400 W/mK, a value of q o = 17.37 W/cm2 may be maintained.

. With k = 400 W/mK and h = 10,000 W/m2K (a practical upper limit), q o = 28.65 W/cm2.
Additional, albeit small, improvements may be realized by relocating the coolant channels closer to the
base of the cold plate.