10

5. Numerical Analysis for Heat Conduction

1. Why Numerical Methods

›  applicable to complex geometries and boundary conditions, and variable thermal properties
›  the advanced mathematical skills are not required
›  easily fit to practical “what-if” optimization problems
›  the advance in computing technologies

2. Steady-State Heat Conduction

­  Numerical Formulation for Interior Nodes

›  1D heat conduction

where:

\ or

m = 1,2,…,M-1

›  2D and 3D heat conduction

Consider the discretized element (Dx´Dy´1) around the node (m,n):

where:

Assuming Dx=Dy, we can obtain after substituting into the above energy balance equation:

m = 1,2,…,M-1 n = 1,2,…,N-1

­  Boundary Conditions

Take an internal corner with surface convection as a sample:

where:

Assuming Dx=Dy, we can obtain after substituting into the above energy balance equation:

­  Solution Methods

›  Matrix Inversion Method/Gaussian Method

Noting nodal equations are linear algebraic equations, they can be written as:

or

where:

= = =

By inverting the matrix , the solution matrix is:

where: =

›  Gauss-Seidel Method
  1. Re-order the nodal equations as much as possible so that |aii| > |aij|, where
  2. Express the in explicit form with the rest of nodal temp.
  1. Assume the initial for each node reasonably

4.  Calculate the new by substituting the assumed temperatures

5.  Repeat the above procedure to calculate the new with the most recent estimates of until the convergence criterion is satisfied

3. Transient Heat Conduction

­  FD Discretization

Space discretization:

Time discretization:

= incremental time step

= nodal temperature at the beginning of the time step

= nodal temperature at the end of the time step

­  Explicit FD Method

Consider a 2D inner node (m,n) at a time step (tp+1 –tp):

where:

Substituting them into the energy equation leads to:

If and , it can be further simplified as:

m = 1, 2, …, M-1 n = 1, 2, …, N-1

Note: see Table 5.2 for other transient FD nodal equations including the nodes subject to the boundary conditions

›  Solution Technique

Noting that is explicitly expressed as ’s, the solution technique is:

1) write the FD equation at each node,

2) starting at , we know all from initial conditions,

3) at , the temperature change commences,

4) solve the FD equations for at all the nodes,

5) renew the with the most current temperature at the node,

6) keep track of the time

7) repeat the above process in next time steps until the final time is reached.

›  Stability Criterion

Note that if is too large or if is too small, the coefficient would become negative, which may cause the solution oscillating and being physically impossible.

\ must be so chosen that the coefficients of !

For the above 2D interior nodes, the stability criterion is:

­  Implicit FD Method

Instead of FD formulation at tp, we replace the spatial derivatives at tp+1,

In similar way, we can write the energy balance FD equation:

If and , it can be further simplified as:

a) See Table 5.2 for other implicit FD equations

b) The above equations are implicit, thus, must be solved simultaneously, e.g., using Gauss/Seidal method

c)  The solution of the implicit FD method is unconditionally stable for any chosen Dx and Dt.

4. Some Remarks

­  Irregular Boundaries

Approximate it with a series of simple volume elements.

­  Controlling the Numerical Errors

›  Types of numerical errors

Discretization error Round-off error

›  Controlling the numerical errors

1)  Discretization error: check the solution accuracy so that it is grid-independent (e.g, Dx® Dx/2, Dx/4, ¼,).

2)  Round-off error: check the solution accuracy with more significant digits.

3)  Reduce the grid size in regions of importance.

­  In summary

›  Overlay the system with a rectangular grid
›  Place nodes at the intersections of the gridlines
›  Associate an incremental area (volume) with each node

›  Write energy balance on each area in terms of nodal T

›  Solve the resultant simultaneous equations for each of the nodal temperatures