001:       SUBROUTINE ZGTTRF( N, DL, D, DU, DU2, IPIV, INFO )
002: *
003: *  -- LAPACK routine (version 3.2) --
004: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
005: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
006: *     November 2006
007: *
008: *     .. Scalar Arguments ..
009:       INTEGER            INFO, N
010: *     ..
011: *     .. Array Arguments ..
012:       INTEGER            IPIV( * )
013:       COMPLEX*16         D( * ), DL( * ), DU( * ), DU2( * )
014: *     ..
015: *
016: *  Purpose
017: *  =======
018: *
019: *  ZGTTRF computes an LU factorization of a complex tridiagonal matrix A
020: *  using elimination with partial pivoting and row interchanges.
021: *
022: *  The factorization has the form
023: *     A = L * U
024: *  where L is a product of permutation and unit lower bidiagonal
025: *  matrices and U is upper triangular with nonzeros in only the main
026: *  diagonal and first two superdiagonals.
027: *
028: *  Arguments
029: *  =========
030: *
031: *  N       (input) INTEGER
032: *          The order of the matrix A.
033: *
034: *  DL      (input/output) COMPLEX*16 array, dimension (N-1)
035: *          On entry, DL must contain the (n-1) sub-diagonal elements of
036: *          A.
037: *
038: *          On exit, DL is overwritten by the (n-1) multipliers that
039: *          define the matrix L from the LU factorization of A.
040: *
041: *  D       (input/output) COMPLEX*16 array, dimension (N)
042: *          On entry, D must contain the diagonal elements of A.
043: *
044: *          On exit, D is overwritten by the n diagonal elements of the
045: *          upper triangular matrix U from the LU factorization of A.
046: *
047: *  DU      (input/output) COMPLEX*16 array, dimension (N-1)
048: *          On entry, DU must contain the (n-1) super-diagonal elements
049: *          of A.
050: *
051: *          On exit, DU is overwritten by the (n-1) elements of the first
052: *          super-diagonal of U.
053: *
054: *  DU2     (output) COMPLEX*16 array, dimension (N-2)
055: *          On exit, DU2 is overwritten by the (n-2) elements of the
056: *          second super-diagonal of U.
057: *
058: *  IPIV    (output) INTEGER array, dimension (N)
059: *          The pivot indices; for 1 <= i <= n, row i of the matrix was
060: *          interchanged with row IPIV(i).  IPIV(i) will always be either
061: *          i or i+1; IPIV(i) = i indicates a row interchange was not
062: *          required.
063: *
064: *  INFO    (output) INTEGER
065: *          = 0:  successful exit
066: *          < 0:  if INFO = -k, the k-th argument had an illegal value
067: *          > 0:  if INFO = k, U(k,k) is exactly zero. The factorization
068: *                has been completed, but the factor U is exactly
069: *                singular, and division by zero will occur if it is used
070: *                to solve a system of equations.
071: *
072: *  =====================================================================
073: *
074: *     .. Parameters ..
075:       DOUBLE PRECISION   ZERO
076:       PARAMETER          ( ZERO = 0.0D+0 )
077: *     ..
078: *     .. Local Scalars ..
079:       INTEGER            I
080:       COMPLEX*16         FACT, TEMP, ZDUM
081: *     ..
082: *     .. External Subroutines ..
083:       EXTERNAL           XERBLA
084: *     ..
085: *     .. Intrinsic Functions ..
086:       INTRINSIC          ABS, DBLE, DIMAG
087: *     ..
088: *     .. Statement Functions ..
089:       DOUBLE PRECISION   CABS1
090: *     ..
091: *     .. Statement Function definitions ..
092:       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
093: *     ..
094: *     .. Executable Statements ..
095: *
096:       INFO = 0
097:       IF( N.LT.0 ) THEN
098:          INFO = -1
099:          CALL XERBLA( 'ZGTTRF', -INFO )
100:          RETURN
101:       END IF
102: *
103: *     Quick return if possible
104: *
105:       IF( N.EQ.0 )
106:      $   RETURN
107: *
108: *     Initialize IPIV(i) = i and DU2(i) = 0
109: *
110:       DO 10 I = 1, N
111:          IPIV( I ) = I
112:    10 CONTINUE
113:       DO 20 I = 1, N - 2
114:          DU2( I ) = ZERO
115:    20 CONTINUE
116: *
117:       DO 30 I = 1, N - 2
118:          IF( CABS1( D( I ) ).GE.CABS1( DL( I ) ) ) THEN
119: *
120: *           No row interchange required, eliminate DL(I)
121: *
122:             IF( CABS1( D( I ) ).NE.ZERO ) THEN
123:                FACT = DL( I ) / D( I )
124:                DL( I ) = FACT
125:                D( I+1 ) = D( I+1 ) - FACT*DU( I )
126:             END IF
127:          ELSE
128: *
129: *           Interchange rows I and I+1, eliminate DL(I)
130: *
131:             FACT = D( I ) / DL( I )
132:             D( I ) = DL( I )
133:             DL( I ) = FACT
134:             TEMP = DU( I )
135:             DU( I ) = D( I+1 )
136:             D( I+1 ) = TEMP - FACT*D( I+1 )
137:             DU2( I ) = DU( I+1 )
138:             DU( I+1 ) = -FACT*DU( I+1 )
139:             IPIV( I ) = I + 1
140:          END IF
141:    30 CONTINUE
142:       IF( N.GT.1 ) THEN
143:          I = N - 1
144:          IF( CABS1( D( I ) ).GE.CABS1( DL( I ) ) ) THEN
145:             IF( CABS1( D( I ) ).NE.ZERO ) THEN
146:                FACT = DL( I ) / D( I )
147:                DL( I ) = FACT
148:                D( I+1 ) = D( I+1 ) - FACT*DU( I )
149:             END IF
150:          ELSE
151:             FACT = D( I ) / DL( I )
152:             D( I ) = DL( I )
153:             DL( I ) = FACT
154:             TEMP = DU( I )
155:             DU( I ) = D( I+1 )
156:             D( I+1 ) = TEMP - FACT*D( I+1 )
157:             IPIV( I ) = I + 1
158:          END IF
159:       END IF
160: *
161: *     Check for a zero on the diagonal of U.
162: *
163:       DO 40 I = 1, N
164:          IF( CABS1( D( I ) ).EQ.ZERO ) THEN
165:             INFO = I
166:             GO TO 50
167:          END IF
168:    40 CONTINUE
169:    50 CONTINUE
170: *
171:       RETURN
172: *
173: *     End of ZGTTRF
174: *
175:       END
176: