File:  [local] / rpl / lapack / lapack / dgetf2.f
Revision 1.2: download - view: text, annotated - select for diffs - revision graph
Wed Apr 21 13:45:14 2010 UTC (14 years, 1 month ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_0_17, rpl-4_0_16, rpl-4_0_15, HEAD
En route pour la 4.0.15 !

    1:       SUBROUTINE DGETF2( M, N, A, LDA, IPIV, 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, LDA, M, N
   10: *     ..
   11: *     .. Array Arguments ..
   12:       INTEGER            IPIV( * )
   13:       DOUBLE PRECISION   A( LDA, * )
   14: *     ..
   15: *
   16: *  Purpose
   17: *  =======
   18: *
   19: *  DGETF2 computes an LU factorization of a general m-by-n matrix A
   20: *  using partial pivoting with row interchanges.
   21: *
   22: *  The factorization has the form
   23: *     A = P * L * U
   24: *  where P is a permutation matrix, L is lower triangular with unit
   25: *  diagonal elements (lower trapezoidal if m > n), and U is upper
   26: *  triangular (upper trapezoidal if m < n).
   27: *
   28: *  This is the right-looking Level 2 BLAS version of the algorithm.
   29: *
   30: *  Arguments
   31: *  =========
   32: *
   33: *  M       (input) INTEGER
   34: *          The number of rows of the matrix A.  M >= 0.
   35: *
   36: *  N       (input) INTEGER
   37: *          The number of columns of the matrix A.  N >= 0.
   38: *
   39: *  A       (input/output) DOUBLE PRECISION array, dimension (LDA,N)
   40: *          On entry, the m by n matrix to be factored.
   41: *          On exit, the factors L and U from the factorization
   42: *          A = P*L*U; the unit diagonal elements of L are not stored.
   43: *
   44: *  LDA     (input) INTEGER
   45: *          The leading dimension of the array A.  LDA >= max(1,M).
   46: *
   47: *  IPIV    (output) INTEGER array, dimension (min(M,N))
   48: *          The pivot indices; for 1 <= i <= min(M,N), row i of the
   49: *          matrix was interchanged with row IPIV(i).
   50: *
   51: *  INFO    (output) INTEGER
   52: *          = 0: successful exit
   53: *          < 0: if INFO = -k, the k-th argument had an illegal value
   54: *          > 0: if INFO = k, U(k,k) is exactly zero. The factorization
   55: *               has been completed, but the factor U is exactly
   56: *               singular, and division by zero will occur if it is used
   57: *               to solve a system of equations.
   58: *
   59: *  =====================================================================
   60: *
   61: *     .. Parameters ..
   62:       DOUBLE PRECISION   ONE, ZERO
   63:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
   64: *     ..
   65: *     .. Local Scalars ..
   66:       DOUBLE PRECISION   SFMIN 
   67:       INTEGER            I, J, JP
   68: *     ..
   69: *     .. External Functions ..
   70:       DOUBLE PRECISION   DLAMCH      
   71:       INTEGER            IDAMAX
   72:       EXTERNAL           DLAMCH, IDAMAX
   73: *     ..
   74: *     .. External Subroutines ..
   75:       EXTERNAL           DGER, DSCAL, DSWAP, XERBLA
   76: *     ..
   77: *     .. Intrinsic Functions ..
   78:       INTRINSIC          MAX, MIN
   79: *     ..
   80: *     .. Executable Statements ..
   81: *
   82: *     Test the input parameters.
   83: *
   84:       INFO = 0
   85:       IF( M.LT.0 ) THEN
   86:          INFO = -1
   87:       ELSE IF( N.LT.0 ) THEN
   88:          INFO = -2
   89:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
   90:          INFO = -4
   91:       END IF
   92:       IF( INFO.NE.0 ) THEN
   93:          CALL XERBLA( 'DGETF2', -INFO )
   94:          RETURN
   95:       END IF
   96: *
   97: *     Quick return if possible
   98: *
   99:       IF( M.EQ.0 .OR. N.EQ.0 )
  100:      $   RETURN
  101: *
  102: *     Compute machine safe minimum 
  103:   104:       SFMIN = DLAMCH('S')  
  105: *
  106:       DO 10 J = 1, MIN( M, N )
  107: *
  108: *        Find pivot and test for singularity.
  109: *
  110:          JP = J - 1 + IDAMAX( M-J+1, A( J, J ), 1 )
  111:          IPIV( J ) = JP
  112:          IF( A( JP, J ).NE.ZERO ) THEN
  113: *
  114: *           Apply the interchange to columns 1:N.
  115: *
  116:             IF( JP.NE.J )
  117:      $         CALL DSWAP( N, A( J, 1 ), LDA, A( JP, 1 ), LDA )
  118: *
  119: *           Compute elements J+1:M of J-th column.
  120: *
  121:             IF( J.LT.M ) THEN 
  122:                IF( ABS(A( J, J )) .GE. SFMIN ) THEN 
  123:                   CALL DSCAL( M-J, ONE / A( J, J ), A( J+1, J ), 1 ) 
  124:                ELSE 
  125:                  DO 20 I = 1, M-J 
  126:                     A( J+I, J ) = A( J+I, J ) / A( J, J ) 
  127:    20            CONTINUE 
  128:                END IF 
  129:             END IF 
  130: *
  131:          ELSE IF( INFO.EQ.0 ) THEN
  132: *
  133:             INFO = J
  134:          END IF
  135: *
  136:          IF( J.LT.MIN( M, N ) ) THEN
  137: *
  138: *           Update trailing submatrix.
  139: *
  140:             CALL DGER( M-J, N-J, -ONE, A( J+1, J ), 1, A( J, J+1 ), LDA,
  141:      $                 A( J+1, J+1 ), LDA )
  142:          END IF
  143:    10 CONTINUE
  144:       RETURN
  145: *
  146: *     End of DGETF2
  147: *
  148:       END

CVSweb interface <joel.bertrand@systella.fr>