Annotation of rpl/lapack/lapack/zgtsv.f, revision 1.3

1.1       bertrand    1:       SUBROUTINE ZGTSV( N, NRHS, DL, D, DU, B, LDB, INFO )
                      2: *
                      3: *  -- LAPACK routine (version 3.2) --
                      4: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                      5: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                      6: *     November 2006
                      7: *
                      8: *     .. Scalar Arguments ..
                      9:       INTEGER            INFO, LDB, N, NRHS
                     10: *     ..
                     11: *     .. Array Arguments ..
                     12:       COMPLEX*16         B( LDB, * ), D( * ), DL( * ), DU( * )
                     13: *     ..
                     14: *
                     15: *  Purpose
                     16: *  =======
                     17: *
                     18: *  ZGTSV  solves the equation
                     19: *
                     20: *     A*X = B,
                     21: *
                     22: *  where A is an N-by-N tridiagonal matrix, by Gaussian elimination with
                     23: *  partial pivoting.
                     24: *
                     25: *  Note that the equation  A'*X = B  may be solved by interchanging the
                     26: *  order of the arguments DU and DL.
                     27: *
                     28: *  Arguments
                     29: *  =========
                     30: *
                     31: *  N       (input) INTEGER
                     32: *          The order of the matrix A.  N >= 0.
                     33: *
                     34: *  NRHS    (input) INTEGER
                     35: *          The number of right hand sides, i.e., the number of columns
                     36: *          of the matrix B.  NRHS >= 0.
                     37: *
                     38: *  DL      (input/output) COMPLEX*16 array, dimension (N-1)
                     39: *          On entry, DL must contain the (n-1) subdiagonal elements of
                     40: *          A.
                     41: *          On exit, DL is overwritten by the (n-2) elements of the
                     42: *          second superdiagonal of the upper triangular matrix U from
                     43: *          the LU factorization of A, in DL(1), ..., DL(n-2).
                     44: *
                     45: *  D       (input/output) COMPLEX*16 array, dimension (N)
                     46: *          On entry, D must contain the diagonal elements of A.
                     47: *          On exit, D is overwritten by the n diagonal elements of U.
                     48: *
                     49: *  DU      (input/output) COMPLEX*16 array, dimension (N-1)
                     50: *          On entry, DU must contain the (n-1) superdiagonal elements
                     51: *          of A.
                     52: *          On exit, DU is overwritten by the (n-1) elements of the first
                     53: *          superdiagonal of U.
                     54: *
                     55: *  B       (input/output) COMPLEX*16 array, dimension (LDB,NRHS)
                     56: *          On entry, the N-by-NRHS right hand side matrix B.
                     57: *          On exit, if INFO = 0, the N-by-NRHS solution matrix X.
                     58: *
                     59: *  LDB     (input) INTEGER
                     60: *          The leading dimension of the array B.  LDB >= max(1,N).
                     61: *
                     62: *  INFO    (output) INTEGER
                     63: *          = 0:  successful exit
                     64: *          < 0:  if INFO = -i, the i-th argument had an illegal value
                     65: *          > 0:  if INFO = i, U(i,i) is exactly zero, and the solution
                     66: *                has not been computed.  The factorization has not been
                     67: *                completed unless i = N.
                     68: *
                     69: *  =====================================================================
                     70: *
                     71: *     .. Parameters ..
                     72:       COMPLEX*16         ZERO
                     73:       PARAMETER          ( ZERO = ( 0.0D+0, 0.0D+0 ) )
                     74: *     ..
                     75: *     .. Local Scalars ..
                     76:       INTEGER            J, K
                     77:       COMPLEX*16         MULT, TEMP, ZDUM
                     78: *     ..
                     79: *     .. Intrinsic Functions ..
                     80:       INTRINSIC          ABS, DBLE, DIMAG, MAX
                     81: *     ..
                     82: *     .. External Subroutines ..
                     83:       EXTERNAL           XERBLA
                     84: *     ..
                     85: *     .. Statement Functions ..
                     86:       DOUBLE PRECISION   CABS1
                     87: *     ..
                     88: *     .. Statement Function definitions ..
                     89:       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
                     90: *     ..
                     91: *     .. Executable Statements ..
                     92: *
                     93:       INFO = 0
                     94:       IF( N.LT.0 ) THEN
                     95:          INFO = -1
                     96:       ELSE IF( NRHS.LT.0 ) THEN
                     97:          INFO = -2
                     98:       ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
                     99:          INFO = -7
                    100:       END IF
                    101:       IF( INFO.NE.0 ) THEN
                    102:          CALL XERBLA( 'ZGTSV ', -INFO )
                    103:          RETURN
                    104:       END IF
                    105: *
                    106:       IF( N.EQ.0 )
                    107:      $   RETURN
                    108: *
                    109:       DO 30 K = 1, N - 1
                    110:          IF( DL( K ).EQ.ZERO ) THEN
                    111: *
                    112: *           Subdiagonal is zero, no elimination is required.
                    113: *
                    114:             IF( D( K ).EQ.ZERO ) THEN
                    115: *
                    116: *              Diagonal is zero: set INFO = K and return; a unique
                    117: *              solution can not be found.
                    118: *
                    119:                INFO = K
                    120:                RETURN
                    121:             END IF
                    122:          ELSE IF( CABS1( D( K ) ).GE.CABS1( DL( K ) ) ) THEN
                    123: *
                    124: *           No row interchange required
                    125: *
                    126:             MULT = DL( K ) / D( K )
                    127:             D( K+1 ) = D( K+1 ) - MULT*DU( K )
                    128:             DO 10 J = 1, NRHS
                    129:                B( K+1, J ) = B( K+1, J ) - MULT*B( K, J )
                    130:    10       CONTINUE
                    131:             IF( K.LT.( N-1 ) )
                    132:      $         DL( K ) = ZERO
                    133:          ELSE
                    134: *
                    135: *           Interchange rows K and K+1
                    136: *
                    137:             MULT = D( K ) / DL( K )
                    138:             D( K ) = DL( K )
                    139:             TEMP = D( K+1 )
                    140:             D( K+1 ) = DU( K ) - MULT*TEMP
                    141:             IF( K.LT.( N-1 ) ) THEN
                    142:                DL( K ) = DU( K+1 )
                    143:                DU( K+1 ) = -MULT*DL( K )
                    144:             END IF
                    145:             DU( K ) = TEMP
                    146:             DO 20 J = 1, NRHS
                    147:                TEMP = B( K, J )
                    148:                B( K, J ) = B( K+1, J )
                    149:                B( K+1, J ) = TEMP - MULT*B( K+1, J )
                    150:    20       CONTINUE
                    151:          END IF
                    152:    30 CONTINUE
                    153:       IF( D( N ).EQ.ZERO ) THEN
                    154:          INFO = N
                    155:          RETURN
                    156:       END IF
                    157: *
                    158: *     Back solve with the matrix U from the factorization.
                    159: *
                    160:       DO 50 J = 1, NRHS
                    161:          B( N, J ) = B( N, J ) / D( N )
                    162:          IF( N.GT.1 )
                    163:      $      B( N-1, J ) = ( B( N-1, J )-DU( N-1 )*B( N, J ) ) / D( N-1 )
                    164:          DO 40 K = N - 2, 1, -1
                    165:             B( K, J ) = ( B( K, J )-DU( K )*B( K+1, J )-DL( K )*
                    166:      $                  B( K+2, J ) ) / D( K )
                    167:    40    CONTINUE
                    168:    50 CONTINUE
                    169: *
                    170:       RETURN
                    171: *
                    172: *     End of ZGTSV
                    173: *
                    174:       END

CVSweb interface <joel.bertrand@systella.fr>