Annotation of rpl/lapack/lapack/dgetrf2.f, revision 1.1
1.1 ! bertrand 1: *> \brief \b DGETRF2
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: * Definition:
! 9: * ===========
! 10: *
! 11: * RECURSIVE SUBROUTINE DGETRF2( M, N, A, LDA, IPIV, INFO )
! 12: *
! 13: * .. Scalar Arguments ..
! 14: * INTEGER INFO, LDA, M, N
! 15: * ..
! 16: * .. Array Arguments ..
! 17: * INTEGER IPIV( * )
! 18: * DOUBLE PRECISION A( LDA, * )
! 19: * ..
! 20: *
! 21: *
! 22: *> \par Purpose:
! 23: * =============
! 24: *>
! 25: *> \verbatim
! 26: *>
! 27: *> DGETRF2 computes an LU factorization of a general M-by-N matrix A
! 28: *> using partial pivoting with row interchanges.
! 29: *>
! 30: *> The factorization has the form
! 31: *> A = P * L * U
! 32: *> where P is a permutation matrix, L is lower triangular with unit
! 33: *> diagonal elements (lower trapezoidal if m > n), and U is upper
! 34: *> triangular (upper trapezoidal if m < n).
! 35: *>
! 36: *> This is the recursive version of the algorithm. It divides
! 37: *> the matrix into four submatrices:
! 38: *>
! 39: *> [ A11 | A12 ] where A11 is n1 by n1 and A22 is n2 by n2
! 40: *> A = [ -----|----- ] with n1 = min(m,n)
! 41: *> [ A21 | A22 ] n2 = n-n1
! 42: *>
! 43: *> [ A11 ]
! 44: *> The subroutine calls itself to factor [ --- ],
! 45: *> [ A12 ]
! 46: *> [ A12 ]
! 47: *> do the swaps on [ --- ], solve A12, update A22,
! 48: *> [ A22 ]
! 49: *>
! 50: *> then calls itself to factor A22 and do the swaps on A21.
! 51: *>
! 52: *> \endverbatim
! 53: *
! 54: * Arguments:
! 55: * ==========
! 56: *
! 57: *> \param[in] M
! 58: *> \verbatim
! 59: *> M is INTEGER
! 60: *> The number of rows of the matrix A. M >= 0.
! 61: *> \endverbatim
! 62: *>
! 63: *> \param[in] N
! 64: *> \verbatim
! 65: *> N is INTEGER
! 66: *> The number of columns of the matrix A. N >= 0.
! 67: *> \endverbatim
! 68: *>
! 69: *> \param[in,out] A
! 70: *> \verbatim
! 71: *> A is DOUBLE PRECISION array, dimension (LDA,N)
! 72: *> On entry, the M-by-N matrix to be factored.
! 73: *> On exit, the factors L and U from the factorization
! 74: *> A = P*L*U; the unit diagonal elements of L are not stored.
! 75: *> \endverbatim
! 76: *>
! 77: *> \param[in] LDA
! 78: *> \verbatim
! 79: *> LDA is INTEGER
! 80: *> The leading dimension of the array A. LDA >= max(1,M).
! 81: *> \endverbatim
! 82: *>
! 83: *> \param[out] IPIV
! 84: *> \verbatim
! 85: *> IPIV is INTEGER array, dimension (min(M,N))
! 86: *> The pivot indices; for 1 <= i <= min(M,N), row i of the
! 87: *> matrix was interchanged with row IPIV(i).
! 88: *> \endverbatim
! 89: *>
! 90: *> \param[out] INFO
! 91: *> \verbatim
! 92: *> INFO is INTEGER
! 93: *> = 0: successful exit
! 94: *> < 0: if INFO = -i, the i-th argument had an illegal value
! 95: *> > 0: if INFO = i, U(i,i) is exactly zero. The factorization
! 96: *> has been completed, but the factor U is exactly
! 97: *> singular, and division by zero will occur if it is used
! 98: *> to solve a system of equations.
! 99: *> \endverbatim
! 100: *
! 101: * Authors:
! 102: * ========
! 103: *
! 104: *> \author Univ. of Tennessee
! 105: *> \author Univ. of California Berkeley
! 106: *> \author Univ. of Colorado Denver
! 107: *> \author NAG Ltd.
! 108: *
! 109: *> \date November 2015
! 110: *
! 111: *> \ingroup doubleGEcomputational
! 112: *
! 113: * =====================================================================
! 114: RECURSIVE SUBROUTINE DGETRF2( M, N, A, LDA, IPIV, INFO )
! 115: *
! 116: * -- LAPACK computational routine (version 3.6.0) --
! 117: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 118: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 119: * November 2015
! 120: *
! 121: * .. Scalar Arguments ..
! 122: INTEGER INFO, LDA, M, N
! 123: * ..
! 124: * .. Array Arguments ..
! 125: INTEGER IPIV( * )
! 126: DOUBLE PRECISION A( LDA, * )
! 127: * ..
! 128: *
! 129: * =====================================================================
! 130: *
! 131: * .. Parameters ..
! 132: DOUBLE PRECISION ONE, ZERO
! 133: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
! 134: * ..
! 135: * .. Local Scalars ..
! 136: DOUBLE PRECISION SFMIN, TEMP
! 137: INTEGER I, IINFO, N1, N2
! 138: * ..
! 139: * .. External Functions ..
! 140: DOUBLE PRECISION DLAMCH
! 141: INTEGER IDAMAX
! 142: EXTERNAL DLAMCH, IDAMAX
! 143: * ..
! 144: * .. External Subroutines ..
! 145: EXTERNAL DGEMM, DSCAL, DLASWP, DTRSM, XERBLA
! 146: * ..
! 147: * .. Intrinsic Functions ..
! 148: INTRINSIC MAX, MIN
! 149: * ..
! 150: * .. Executable Statements ..
! 151: *
! 152: * Test the input parameters
! 153: *
! 154: INFO = 0
! 155: IF( M.LT.0 ) THEN
! 156: INFO = -1
! 157: ELSE IF( N.LT.0 ) THEN
! 158: INFO = -2
! 159: ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
! 160: INFO = -4
! 161: END IF
! 162: IF( INFO.NE.0 ) THEN
! 163: CALL XERBLA( 'DGETRF2', -INFO )
! 164: RETURN
! 165: END IF
! 166: *
! 167: * Quick return if possible
! 168: *
! 169: IF( M.EQ.0 .OR. N.EQ.0 )
! 170: $ RETURN
! 171:
! 172: IF ( M.EQ.1 ) THEN
! 173: *
! 174: * Use unblocked code for one row case
! 175: * Just need to handle IPIV and INFO
! 176: *
! 177: IPIV( 1 ) = 1
! 178: IF ( A(1,1).EQ.ZERO )
! 179: $ INFO = 1
! 180: *
! 181: ELSE IF( N.EQ.1 ) THEN
! 182: *
! 183: * Use unblocked code for one column case
! 184: *
! 185: *
! 186: * Compute machine safe minimum
! 187: *
! 188: SFMIN = DLAMCH('S')
! 189: *
! 190: * Find pivot and test for singularity
! 191: *
! 192: I = IDAMAX( M, A( 1, 1 ), 1 )
! 193: IPIV( 1 ) = I
! 194: IF( A( I, 1 ).NE.ZERO ) THEN
! 195: *
! 196: * Apply the interchange
! 197: *
! 198: IF( I.NE.1 ) THEN
! 199: TEMP = A( 1, 1 )
! 200: A( 1, 1 ) = A( I, 1 )
! 201: A( I, 1 ) = TEMP
! 202: END IF
! 203: *
! 204: * Compute elements 2:M of the column
! 205: *
! 206: IF( ABS(A( 1, 1 )) .GE. SFMIN ) THEN
! 207: CALL DSCAL( M-1, ONE / A( 1, 1 ), A( 2, 1 ), 1 )
! 208: ELSE
! 209: DO 10 I = 1, M-1
! 210: A( 1+I, 1 ) = A( 1+I, 1 ) / A( 1, 1 )
! 211: 10 CONTINUE
! 212: END IF
! 213: *
! 214: ELSE
! 215: INFO = 1
! 216: END IF
! 217: *
! 218: ELSE
! 219: *
! 220: * Use recursive code
! 221: *
! 222: N1 = MIN( M, N ) / 2
! 223: N2 = N-N1
! 224: *
! 225: * [ A11 ]
! 226: * Factor [ --- ]
! 227: * [ A21 ]
! 228: *
! 229: CALL DGETRF2( M, N1, A, LDA, IPIV, IINFO )
! 230:
! 231: IF ( INFO.EQ.0 .AND. IINFO.GT.0 )
! 232: $ INFO = IINFO
! 233: *
! 234: * [ A12 ]
! 235: * Apply interchanges to [ --- ]
! 236: * [ A22 ]
! 237: *
! 238: CALL DLASWP( N2, A( 1, N1+1 ), LDA, 1, N1, IPIV, 1 )
! 239: *
! 240: * Solve A12
! 241: *
! 242: CALL DTRSM( 'L', 'L', 'N', 'U', N1, N2, ONE, A, LDA,
! 243: $ A( 1, N1+1 ), LDA )
! 244: *
! 245: * Update A22
! 246: *
! 247: CALL DGEMM( 'N', 'N', M-N1, N2, N1, -ONE, A( N1+1, 1 ), LDA,
! 248: $ A( 1, N1+1 ), LDA, ONE, A( N1+1, N1+1 ), LDA )
! 249: *
! 250: * Factor A22
! 251: *
! 252: CALL DGETRF2( M-N1, N2, A( N1+1, N1+1 ), LDA, IPIV( N1+1 ),
! 253: $ IINFO )
! 254: *
! 255: * Adjust INFO and the pivot indices
! 256: *
! 257: IF ( INFO.EQ.0 .AND. IINFO.GT.0 )
! 258: $ INFO = IINFO + N1
! 259: DO 20 I = N1+1, MIN( M, N )
! 260: IPIV( I ) = IPIV( I ) + N1
! 261: 20 CONTINUE
! 262: *
! 263: * Apply interchanges to A21
! 264: *
! 265: CALL DLASWP( N1, A( 1, 1 ), LDA, N1+1, MIN( M, N), IPIV, 1 )
! 266: *
! 267: END IF
! 268: RETURN
! 269: *
! 270: * End of DGETRF2
! 271: *
! 272: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dgetrf2.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack