SRC/zgttrf.f

   1 *> \brief \b ZGTTRF
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download ZGTTRF + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgttrf.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgttrf.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgttrf.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE ZGTTRF( N, DL, D, DU, DU2, IPIV, INFO )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       INTEGER            INFO, N
  25 *       ..
  26 *       .. Array Arguments ..
  27 *       INTEGER            IPIV( * )
  28 *       COMPLEX*16         D( * ), DL( * ), DU( * ), DU2( * )
  29 *       ..
  30 *
  31 *
  32 *> \par Purpose:
  33 *  =============
  34 *>
  35 *> \verbatim
  36 *>
  37 *> ZGTTRF computes an LU factorization of a complex tridiagonal matrix A
  38 *> using elimination with partial pivoting and row interchanges.
  39 *>
  40 *> The factorization has the form
  41 *>    A = L * U
  42 *> where L is a product of permutation and unit lower bidiagonal
  43 *> matrices and U is upper triangular with nonzeros in only the main
  44 *> diagonal and first two superdiagonals.
  45 *> \endverbatim
  46 *
  47 *  Arguments:
  48 *  ==========
  49 *
  50 *> \param[in] N
  51 *> \verbatim
  52 *>          N is INTEGER
  53 *>          The order of the matrix A.
  54 *> \endverbatim
  55 *>
  56 *> \param[in,out] DL
  57 *> \verbatim
  58 *>          DL is COMPLEX*16 array, dimension (N-1)
  59 *>          On entry, DL must contain the (n-1) sub-diagonal elements of
  60 *>          A.
  61 *>
  62 *>          On exit, DL is overwritten by the (n-1) multipliers that
  63 *>          define the matrix L from the LU factorization of A.
  64 *> \endverbatim
  65 *>
  66 *> \param[in,out] D
  67 *> \verbatim
  68 *>          D is COMPLEX*16 array, dimension (N)
  69 *>          On entry, D must contain the diagonal elements of A.
  70 *>
  71 *>          On exit, D is overwritten by the n diagonal elements of the
  72 *>          upper triangular matrix U from the LU factorization of A.
  73 *> \endverbatim
  74 *>
  75 *> \param[in,out] DU
  76 *> \verbatim
  77 *>          DU is COMPLEX*16 array, dimension (N-1)
  78 *>          On entry, DU must contain the (n-1) super-diagonal elements
  79 *>          of A.
  80 *>
  81 *>          On exit, DU is overwritten by the (n-1) elements of the first
  82 *>          super-diagonal of U.
  83 *> \endverbatim
  84 *>
  85 *> \param[out] DU2
  86 *> \verbatim
  87 *>          DU2 is COMPLEX*16 array, dimension (N-2)
  88 *>          On exit, DU2 is overwritten by the (n-2) elements of the
  89 *>          second super-diagonal of U.
  90 *> \endverbatim
  91 *>
  92 *> \param[out] IPIV
  93 *> \verbatim
  94 *>          IPIV is INTEGER array, dimension (N)
  95 *>          The pivot indices; for 1 <= i <= n, row i of the matrix was
  96 *>          interchanged with row IPIV(i).  IPIV(i) will always be either
  97 *>          i or i+1; IPIV(i) = i indicates a row interchange was not
  98 *>          required.
  99 *> \endverbatim
 100 *>
 101 *> \param[out] INFO
 102 *> \verbatim
 103 *>          INFO is INTEGER
 104 *>          = 0:  successful exit
 105 *>          < 0:  if INFO = -k, the k-th argument had an illegal value
 106 *>          > 0:  if INFO = k, U(k,k) is exactly zero. The factorization
 107 *>                has been completed, but the factor U is exactly
 108 *>                singular, and division by zero will occur if it is used
 109 *>                to solve a system of equations.
 110 *> \endverbatim
 111 *
 112 *  Authors:
 113 *  ========
 114 *
 115 *> \author Univ. of Tennessee
 116 *> \author Univ. of California Berkeley
 117 *> \author Univ. of Colorado Denver
 118 *> \author NAG Ltd.
 119 *
 120 *> \date September 2012
 121 *
 122 *> \ingroup complex16GTcomputational
 123 *
 124 *  =====================================================================
 125       SUBROUTINE ZGTTRF( N, DL, D, DU, DU2, IPIV, INFO )
 126 *
 127 *  -- LAPACK computational routine (version 3.4.2) --
 128 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 129 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 130 *     September 2012
 131 *
 132 *     .. Scalar Arguments ..
 133       INTEGER            INFO, N
 134 *     ..
 135 *     .. Array Arguments ..
 136       INTEGER            IPIV( * )
 137       COMPLEX*16         D( * ), DL( * ), DU( * ), DU2( * )
 138 *     ..
 139 *
 140 *  =====================================================================
 141 *
 142 *     .. Parameters ..
 143       DOUBLE PRECISION   ZERO
 144       PARAMETER          ( ZERO = 0.0D+0 )
 145 *     ..
 146 *     .. Local Scalars ..
 147       INTEGER            I
 148       COMPLEX*16         FACT, TEMP, ZDUM
 149 *     ..
 150 *     .. External Subroutines ..
 151       EXTERNAL           XERBLA
 152 *     ..
 153 *     .. Intrinsic Functions ..
 154       INTRINSIC          ABS, DBLE, DIMAG
 155 *     ..
 156 *     .. Statement Functions ..
 157       DOUBLE PRECISION   CABS1
 158 *     ..
 159 *     .. Statement Function definitions ..
 160       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
 161 *     ..
 162 *     .. Executable Statements ..
 163 *
 164       INFO = 0
 165       IF( N.LT.0 ) THEN
 166          INFO = -1
 167          CALL XERBLA( 'ZGTTRF', -INFO )
 168          RETURN
 169       END IF
 170 *
 171 *     Quick return if possible
 172 *
 173       IF( N.EQ.0 )
 174      $   RETURN
 175 *
 176 *     Initialize IPIV(i) = i and DU2(i) = 0
 177 *
 178       DO 10 I = 1, N
 179          IPIV( I ) = I
 180    10 CONTINUE
 181       DO 20 I = 1, N - 2
 182          DU2( I ) = ZERO
 183    20 CONTINUE
 184 *
 185       DO 30 I = 1, N - 2
 186          IF( CABS1( D( I ) ).GE.CABS1( DL( I ) ) ) THEN
 187 *
 188 *           No row interchange required, eliminate DL(I)
 189 *
 190             IF( CABS1( D( I ) ).NE.ZERO ) THEN
 191                FACT = DL( I ) / D( I )
 192                DL( I ) = FACT
 193                D( I+1 ) = D( I+1 ) - FACT*DU( I )
 194             END IF
 195          ELSE
 196 *
 197 *           Interchange rows I and I+1, eliminate DL(I)
 198 *
 199             FACT = D( I ) / DL( I )
 200             D( I ) = DL( I )
 201             DL( I ) = FACT
 202             TEMP = DU( I )
 203             DU( I ) = D( I+1 )
 204             D( I+1 ) = TEMP - FACT*D( I+1 )
 205             DU2( I ) = DU( I+1 )
 206             DU( I+1 ) = -FACT*DU( I+1 )
 207             IPIV( I ) = I + 1
 208          END IF
 209    30 CONTINUE
 210       IF( N.GT.1 ) THEN
 211          I = N - 1
 212          IF( CABS1( D( I ) ).GE.CABS1( DL( I ) ) ) THEN
 213             IF( CABS1( D( I ) ).NE.ZERO ) THEN
 214                FACT = DL( I ) / D( I )
 215                DL( I ) = FACT
 216                D( I+1 ) = D( I+1 ) - FACT*DU( I )
 217             END IF
 218          ELSE
 219             FACT = D( I ) / DL( I )
 220             D( I ) = DL( I )
 221             DL( I ) = FACT
 222             TEMP = DU( I )
 223             DU( I ) = D( I+1 )
 224             D( I+1 ) = TEMP - FACT*D( I+1 )
 225             IPIV( I ) = I + 1
 226          END IF
 227       END IF
 228 *
 229 *     Check for a zero on the diagonal of U.
 230 *
 231       DO 40 I = 1, N
 232          IF( CABS1( D( I ) ).EQ.ZERO ) THEN
 233             INFO = I
 234             GO TO 50
 235          END IF
 236    40 CONTINUE
 237    50 CONTINUE
 238 *
 239       RETURN
 240 *
 241 *     End of ZGTTRF
 242 *
 243       END