440
Appendix A. Numerical Miscellany
When it does not lead to ambiguity, the grid-point indices are omitted from
finite-difference formulae expressed in operator notation. As an example,
expandsto
-
J
J
q,'! - q,'!
+ J+I
J -I=O.
2tlt
2tlx
Simple finite-difference formulae are generally written out in expanded form.
A.2 Tridiagonal Solvers
Implicit finite-difference schemes for one-dimensional problems often lead to
tridiagonal systems of linear algebraic equations . Tridiagonal systems can be
solved very efficiently using the Thomas tridiagonal algorithm (lsaacson and Keller
1966). The algorithm can be extended to the periodic case as discussed in Strikwerda (1989). FORTRAN programs for both the standard and periodic cases are
given below. A well-behaved solution will be obtained using the following algorithms whenever the tridiagonal systems are diagonally dominant, which will be
the case iffor all j , Ibjl > lajl + ICjl, where a j, b], and Cj are, respectively, the
subdiagonal, diagonal, and superdiagonal entries in row j.
A.2.J Code for a Tridiagonal Solver
subroutine tridiag(jmx,a ,b,c,f,q)
c
Solves a standard tr idiagonal system
c
c
Definition of the variables:
c
c
c
c
c
c
jmx
a
b
c
f
q
dimens ion of all the follow ing arrays
sub (lower) diagonalcenter diagonal
super (upper) diagonalright hand side
work array provided by calling program
c
c
c
a(l) and c(jmx) need not be initialized
The output is i n f ; a, b, and c are unchanged
real a ( * ) ,b (*) , c ( *) , f ( * ) ,q ( * ) , p
integer j, jmx
c(jmx)=O .