Return to zgetc2.f CVS log

Up to [local] / rpl / lapack / lapack

Cohérence Lapack.

    1: *> \brief \b ZGETC2 computes the LU factorization with complete pivoting of the general n-by-n matrix.
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at 
    6: *            http://www.netlib.org/lapack/explore-html/ 
    7: *
    8: *> \htmlonly
    9: *> Download ZGETC2 + dependencies 
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgetc2.f"> 
   11: *> [TGZ]</a> 
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgetc2.f"> 
   13: *> [ZIP]</a> 
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgetc2.f"> 
   15: *> [TXT]</a>
   16: *> \endhtmlonly 
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE ZGETC2( N, A, LDA, IPIV, JPIV, INFO )
   22: * 
   23: *       .. Scalar Arguments ..
   24: *       INTEGER            INFO, LDA, N
   25: *       ..
   26: *       .. Array Arguments ..
   27: *       INTEGER            IPIV( * ), JPIV( * )
   28: *       COMPLEX*16         A( LDA, * )
   29: *       ..
   30: *  
   31: *
   32: *> \par Purpose:
   33: *  =============
   34: *>
   35: *> \verbatim
   36: *>
   37: *> ZGETC2 computes an LU factorization, using complete pivoting, of the
   38: *> n-by-n matrix A. The factorization has the form A = P * L * U * Q,
   39: *> where P and Q are permutation matrices, L is lower triangular with
   40: *> unit diagonal elements and U is upper triangular.
   41: *>
   42: *> This is a level 1 BLAS version of the algorithm.
   43: *> \endverbatim
   44: *
   45: *  Arguments:
   46: *  ==========
   47: *
   48: *> \param[in] N
   49: *> \verbatim
   50: *>          N is INTEGER
   51: *>          The order of the matrix A. N >= 0.
   52: *> \endverbatim
   53: *>
   54: *> \param[in,out] A
   55: *> \verbatim
   56: *>          A is COMPLEX*16 array, dimension (LDA, N)
   57: *>          On entry, the n-by-n matrix to be factored.
   58: *>          On exit, the factors L and U from the factorization
   59: *>          A = P*L*U*Q; the unit diagonal elements of L are not stored.
   60: *>          If U(k, k) appears to be less than SMIN, U(k, k) is given the
   61: *>          value of SMIN, giving a nonsingular perturbed system.
   62: *> \endverbatim
   63: *>
   64: *> \param[in] LDA
   65: *> \verbatim
   66: *>          LDA is INTEGER
   67: *>          The leading dimension of the array A.  LDA >= max(1, N).
   68: *> \endverbatim
   69: *>
   70: *> \param[out] IPIV
   71: *> \verbatim
   72: *>          IPIV is INTEGER array, dimension (N).
   73: *>          The pivot indices; for 1 <= i <= N, row i of the
   74: *>          matrix has been interchanged with row IPIV(i).
   75: *> \endverbatim
   76: *>
   77: *> \param[out] JPIV
   78: *> \verbatim
   79: *>          JPIV is INTEGER array, dimension (N).
   80: *>          The pivot indices; for 1 <= j <= N, column j of the
   81: *>          matrix has been interchanged with column JPIV(j).
   82: *> \endverbatim
   83: *>
   84: *> \param[out] INFO
   85: *> \verbatim
   86: *>          INFO is INTEGER
   87: *>           = 0: successful exit
   88: *>           > 0: if INFO = k, U(k, k) is likely to produce overflow if
   89: *>                one tries to solve for x in Ax = b. So U is perturbed
   90: *>                to avoid the overflow.
   91: *> \endverbatim
   92: *
   93: *  Authors:
   94: *  ========
   95: *
   96: *> \author Univ. of Tennessee 
   97: *> \author Univ. of California Berkeley 
   98: *> \author Univ. of Colorado Denver 
   99: *> \author NAG Ltd. 
  100: *
  101: *> \date June 2016
  102: *
  103: *> \ingroup complex16GEauxiliary
  104: *
  105: *> \par Contributors:
  106: *  ==================
  107: *>
  108: *>     Bo Kagstrom and Peter Poromaa, Department of Computing Science,
  109: *>     Umea University, S-901 87 Umea, Sweden.
  110: *
  111: *  =====================================================================
  112:       SUBROUTINE ZGETC2( N, A, LDA, IPIV, JPIV, INFO )
  113: *
  114: *  -- LAPACK auxiliary routine (version 3.6.1) --
  115: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  116: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  117: *     June 2016
  118: *
  119: *     .. Scalar Arguments ..
  120:       INTEGER            INFO, LDA, N
  121: *     ..
  122: *     .. Array Arguments ..
  123:       INTEGER            IPIV( * ), JPIV( * )
  124:       COMPLEX*16         A( LDA, * )
  125: *     ..
  126: *
  127: *  =====================================================================
  128: *
  129: *     .. Parameters ..
  130:       DOUBLE PRECISION   ZERO, ONE
  131:       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
  132: *     ..
  133: *     .. Local Scalars ..
  134:       INTEGER            I, IP, IPV, J, JP, JPV
  135:       DOUBLE PRECISION   BIGNUM, EPS, SMIN, SMLNUM, XMAX
  136: *     ..
  137: *     .. External Subroutines ..
  138:       EXTERNAL           ZGERU, ZSWAP
  139: *     ..
  140: *     .. External Functions ..
  141:       DOUBLE PRECISION   DLAMCH
  142:       EXTERNAL           DLAMCH
  143: *     ..
  144: *     .. Intrinsic Functions ..
  145:       INTRINSIC          ABS, DCMPLX, MAX
  146: *     ..
  147: *     .. Executable Statements ..
  148: *
  149:       INFO = 0
  150: *
  151: *     Quick return if possible
  152: *
  153:       IF( N.EQ.0 )
  154:      $   RETURN
  155: *
  156: *     Set constants to control overflow
  157: *
  158:       EPS = DLAMCH( 'P' )
  159:       SMLNUM = DLAMCH( 'S' ) / EPS
  160:       BIGNUM = ONE / SMLNUM
  161:       CALL DLABAD( SMLNUM, BIGNUM )
  162: *
  163: *     Handle the case N=1 by itself
  164: *
  165:       IF( N.EQ.1 ) THEN
  166:          IPIV( 1 ) = 1
  167:          JPIV( 1 ) = 1
  168:          IF( ABS( A( 1, 1 ) ).LT.SMLNUM ) THEN
  169:             INFO = 1
  170:             A( 1, 1 ) = DCMPLX( SMLNUM, ZERO )
  171:          END IF
  172:          RETURN
  173:       END IF
  174: *
  175: *     Factorize A using complete pivoting.
  176: *     Set pivots less than SMIN to SMIN
  177: *
  178:       DO 40 I = 1, N - 1
  179: *
  180: *        Find max element in matrix A
  181: *
  182:          XMAX = ZERO
  183:          DO 20 IP = I, N
  184:             DO 10 JP = I, N
  185:                IF( ABS( A( IP, JP ) ).GE.XMAX ) THEN
  186:                   XMAX = ABS( A( IP, JP ) )
  187:                   IPV = IP
  188:                   JPV = JP
  189:                END IF
  190:    10       CONTINUE
  191:    20    CONTINUE
  192:          IF( I.EQ.1 )
  193:      $      SMIN = MAX( EPS*XMAX, SMLNUM )
  194: *
  195: *        Swap rows
  196: *
  197:          IF( IPV.NE.I )
  198:      $      CALL ZSWAP( N, A( IPV, 1 ), LDA, A( I, 1 ), LDA )
  199:          IPIV( I ) = IPV
  200: *
  201: *        Swap columns
  202: *
  203:          IF( JPV.NE.I )
  204:      $      CALL ZSWAP( N, A( 1, JPV ), 1, A( 1, I ), 1 )
  205:          JPIV( I ) = JPV
  206: *
  207: *        Check for singularity
  208: *
  209:          IF( ABS( A( I, I ) ).LT.SMIN ) THEN
  210:             INFO = I
  211:             A( I, I ) = DCMPLX( SMIN, ZERO )
  212:          END IF
  213:          DO 30 J = I + 1, N
  214:             A( J, I ) = A( J, I ) / A( I, I )
  215:    30    CONTINUE
  216:          CALL ZGERU( N-I, N-I, -DCMPLX( ONE ), A( I+1, I ), 1,
  217:      $               A( I, I+1 ), LDA, A( I+1, I+1 ), LDA )
  218:    40 CONTINUE
  219: *
  220:       IF( ABS( A( N, N ) ).LT.SMIN ) THEN
  221:          INFO = N
  222:          A( N, N ) = DCMPLX( SMIN, ZERO )
  223:       END IF
  224: *
  225: *     Set last pivots to N
  226: *
  227:       IPIV( N ) = N
  228:       JPIV( N ) = N
  229: *
  230:       RETURN
  231: *
  232: *     End of ZGETC2
  233: *
  234:       END