Module 5    : NUMERICAL METHODS IN HEAT CONDUCTION
Lecture 15 : FINITE DIFFERENCE FORMULATION OF DIFFERENTIAL EQUATIONS

The numerical methods for solving differential equations are based on replacing the differential equations by algebraic equations. In the case of the popular finite difference method, this is done by replacing the derivatives by differences. Below we will demonstrate this with both first and second order derivatives. But first we give a motivational example.

Derivates are the building blocks of differential equations, and thus we first give a brief review of derivatives. Consider a function f that depends on x , as shown in Figure. 5.1.

Figure 5.1 The derivative of the function at a point represents the slope of the function at that point

The first derivative of f(x) at a point is equivalent to the slope of a line tangent to the curve at that point and is defined as

(5.4)

which is the ratio of the increment of the function to the temperature of the independent variable as . If we don't take the indicated limit, we will have the following approximate relation for the derivative:

(5.5)

This approximate expression of the derivative in terms of differences is the finite difference form of the first derivative. The equation above can also be obtained by writing the Taylor series expansion of the function f about the point x,

(5.6)

and neglecting all the terms in the expansion except the first two. The first term neglected is proportional to , and thus the error involved in each step of this approximation is also proportional to . However, the commutative error involved after M steps in the direction of length L is proportional to since . Therefore, the smaller the , the smaller the error, and thus the more accurate the approximation.

Now consider steady one dimensional heat transfer in a plane wall of thickness L with heat generation. The wall is subdivided into M equal sections of thickness in the x -direction, separated by M+1 points 0,1,2,…..m-1,m,m+1,…..,M called nodes or nodal points, as shown in Figure 5.2. The x -coordinate of any point m is simply , and the temperature at that point is simply T(x m )=Tm.

Figure 5.2 Schematic of the nodes and the nodal temperatures used in the development of the finite difference formulation of heat transfer in a plane wall

The heat conduction equation involves the second derivatives of temperature with respect to the space variables, such as , and the finite difference formulation is based on replacing the second derivatives by appropriate differences. But we need to start the process with first derivatives. Using Equation 5.5, the first derivative of temperature at the midpoints and of the sections surrounding the node m can be expressed as

and

(5.7)

Noting that the second derivative is simply the derivative of the first derivative, the second derivative of temperature at node m can be expressed as

=

(5.8)

which is the finite difference representation of the second derivative at a general internal node m . Note that the second derivative of temperature at a node m is expressed in terms of the temperatures at node m and its two neighboring nodes. Then the differential equation

(5.9)

which is the governing equation for steady one dimensional heat transfer in a plane wall with heat generation and constant thermal conductivity, can be expressed in the finite difference form as (Figure. 5.3)

(5.10)

where is the rate of heat generation at node m . If the surface temperatures Toand TM are specified, the application of this equation to each of the M-1 interior nodes results in M-1 equations for the determination of M-1 unknown temperatures at the interior nodes. Solving these equations simultaneously gives the temperature values at the nodes. If the temperatures at the outer surfaces are not known, then we need to obtain two more equations in a similar manner using the specified boundary conditions. Then the unknown temperatures at M+1 nodes are determined by solving the resulting system of M+1 equations in M+1 unknowns simultaneously.
Figure 5.3 The differential equation is valid at every point of a medium, whereas the finite difference equation is valid at discrete points (the nodes) only.

Note that the boundary conditions have no effect on the finite difference formulation of interior nodes of the medium. This is not surprising since the control volume used in the development of the formulation does not involve any part of the boundary. You may recall that the boundary conditions had no effect ton the differential equation of heat conduction in the medium either.

The finite difference formulation above can easily be extended to two-or-three-dimensional heat transfer problems by replacing each second derivative by a difference equation in that direction. For example, the finite difference formulation for steady two dimensional heat conduction in a region with heat generation and constant thermal conductivity can be expressed in rectangular coordinates as (Figure. 5.4)

Figure. 5.4 Finite Difference Mesh for Two-dimensional conduction in rectangular co-ordinates

(5.11)

For m = 1,2,3,….M-1 and n = 1,2,3,…N-1 at any interior node (m,n) . Note that a rectangular region that is divided into M equal subregions in the x -direction and N equal subregions in the y -direction has a total of (M+1)(N+1) nodes, and Equation 5.10 can be used to obtain the finite difference equations at (M-1)(N-1) of these nodes (i.e., all nodes except those at the boundaries).

The finite difference formulation is given above to demonstrate how difference equations are obtained from differential equations. However, we will use the energy balance approach in the following sections to obtain the numerical formulation because it is more intuitive and can handle boundary conditions more easily. Besides, the energy balance approach does not require having the differential equation before the analysis.