PUBLIC INTERFACE / ROUTINES / CHANGES / ERRORS / REFERENCES / NOTES

Last updated on .


Module Tridiagonal

     Version: v2.0
     Date:    October 5, 1999
     Contact: Isaac Held or Bruce Wyman

OVERVIEW


    Solves the tridiagonal system of equations.


    The following schematic represents the the system of equations solved,
    where X is the solution.

    | B(1)  A(1)   0     0                .......            0    |  |X(1)|   |D(1)| 
    | C(2)  B(2)  A(2)   0                .......            0    |  |X(2)|   |D(2)|
    |  0    C(3)  B(3)  A(3)  0           .......            0    |  | .. |   | .. |
    |  ..........................................                 |  | .. | = | .. |
    |  ..........................................                 |  | .. |   | .. |
    |                                  C(N-2) B(N-2) A(N-2)  0    |  | .. |   | .. |
    |                                    0    C(N-1) B(N-1) A(N-1)|  | .. |   | .. |
    |                                    0      0    C(N)   B(N)  |  |X(N)|   |D(N)|


OTHER MODULES USED


     None.


PUBLIC INTERFACE


   use tridiagonal_mod [, only: tri_invert, close_tridiagonal ]

      tri_invert:        Sets up and solves the tridiagonal 
                         system of equations.

      close_tridiagonal: Releases memory used by the solver.


PUBLIC ROUTINES


call tri_invert ( X, D, A, B, C )

  Input

      D    The right-hand side term, see the schematic above.
             [real, dimension(:,:,:)]

  Output

      X    The solution to the tridiagonal system of equations.
             [real, dimension(:,:,:)]

  Optional Input/Output

      A,B,C  The left-hand-side terms (matrix), see the schematic above.
             If A is not present, it is assumed that the matrix (A,B.C)
             has not been changed since the last call to tri_invert.
                [real, dimension(:,:,:)]

  Notes:

      For simplicity, A and C are assumed to be dimensioned the same size 
      as B, D, and X, although any input values for A(N) and C(1) are ignored.
      There are no checks to make sure the sizes agree.

      The value of A(N) is modified on output, and B and C are unchanged.

--------------------------------------------------------------

call close_tridiagonal

   There are no arguments to this routine.


CHANGE HISTORY


     No recent changes.


ERROR MESSAGES


     None.


REFERENCES


     None.


KNOWN BUGS


   * Optional arguments A,B,C have no intent declaration,
     so the default intent is inout. The value of A(N) is modified
     on output, and B and C are unchanged.


NOTES


  The following private allocatable arrays save the relevant information
  if one recalls tri_invert without changing (A,B,C):

        allocate ( e  (size(x,1), size(x,2), size(x,3)) )
        allocate ( g  (size(x,1), size(x,2), size(x,3)) )
        allocate ( cc (size(x,1), size(x,2), size(x,3)) )
        allocate ( bb (size(x,1), size(x,2)) )

  This storage is deallocated when close_tridiagonal is called.


FUTURE PLANS


     Maybe a cleaner version?