2.4 Analysis of Finite Difference Methods

2.4.4 Finite Difference Methods in Matrix Form

Measurable Outcome 2.3, Measurable Outcome 2.8

Consider the one-dimensional convection-diffusion equation,

  \(\displaystyle \dfrac {\partial U}{\partial t} + u \dfrac {\partial U}{\partial x} -\mu \dfrac {\partial ^2 U}{\partial x^2}\) \(\displaystyle =\) \(\displaystyle 0.\)   (2.73)

Approximating the spatial derivative using the central difference operators gives the following approximation at node \(i\),

  \(\displaystyle \dfrac {d U_ i}{d t} + u_ i \delta _{2x} U_ i - \mu \delta _ x^2 U_ i\) \(\displaystyle =\) \(\displaystyle 0\)   (2.74)

Assembling all of the nodal states into a single vector

  \(\displaystyle U = \left( U_0,\, U_1,\, \ldots ,\, U_{i-1},\, U_ i,\, U_{i+1},\, \ldots ,\, U_{N_ x-1},\, U_{N_ x}\right)^ T,\)   (2.75)

gives a system of coupled ODEs of the form:

  \(\displaystyle \dfrac {d U}{d t}\) \(\displaystyle =\) \(\displaystyle A U + b \label{equ:linSys}\)   (2.76)

where \(b\) will contain the boundary condition related data. The matrix \(A\) has the form:

  \(\displaystyle A\) \(\displaystyle =\) \(\displaystyle \left[\begin{array}{cccc} A_{0,0} & A_{0,1} & \ldots & A_{0,N_ x} \\ A_{1,0} & A_{1,1} & \ldots & A_{1,N_ x} \\ \vdots & \vdots & \ddots & \vdots \\ A_{N_ x,0} & A_{N_ x,1} & \ldots & A_{N_ x,N_ x} \end{array} \right].\)   (2.77)

Note that row \(i\) of this matrix contains the coefficients of the nodal values for the ODE governing node \(i\). Except for rows requiring boundary condition data, the values of \(A_{i,j}\) are related to the coefficients \(\alpha _ j\) of the finite difference approximation. Specifically, for our central difference approximation

  \(\displaystyle A_{i,i-1}\) \(\displaystyle =\) \(\displaystyle \frac{u_ i}{2 \Delta x} + \frac{\mu }{\Delta x^2},\)   (2.78)
  \(\displaystyle A_{i,i}\) \(\displaystyle =\) \(\displaystyle - 2\frac{\mu }{\Delta x^2},\)   (2.79)
  \(\displaystyle A_{i,i+1}\) \(\displaystyle =\) \(\displaystyle -\frac{u_ i}{2 \Delta x} + \frac{\mu }{\Delta x^2},\)   (2.80)

and all other entries of row \(i\) are zero. In general, the number of non-zero entries of row \(i\) will correspond to the size of the stencil of the finite difference approximations used. We refer to Equation 2.76 as being semi-discrete, since we have discretized the PDE in space but not in time. To make this a fully discrete approximation, we can apply any of the ODE integration methods that we discussed previously.

Implementation Notes When explicit methods are used to represent a FD method, then the \(A\) matrix does not need to be directly formed. Implementing a finite difference approximation at each node will suffice as we will only be interested in matrix-vector multiplication, i.e. we only need to apply \(A\).

On the other hand, implicit methods require the solution of a linear system. The solutions methods often require the representation of \(A\) as a matrix, and thus it must be constructed.