File:  [local] / rpl / lapack / lapack / dlaorhr_col_getrfnp.f
Revision 1.1: download - view: text, annotated - select for diffs - revision graph
Thu May 21 21:45:59 2020 UTC (3 years, 11 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_33, rpl-4_1_32, HEAD
Mise à jour de Lapack.

    1: *> \brief \b DLAORHR_COL_GETRFNP
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at
    6: *            http://www.netlib.org/lapack/explore-html/
    7: *
    8: *> \htmlonly
    9: *> Download DLAORHR_COL_GETRFNP + dependencies
   10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaorhr_col_getrfnp.f">
   11: *> [TGZ]</a>
   12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaorhr_col_getrfnp.f">
   13: *> [ZIP]</a>
   14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaorhr_col_getrfnp.f">
   15: *> [TXT]</a>
   16: *> \endhtmlonly
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE DLAORHR_COL_GETRFNP( M, N, A, LDA, D, INFO )
   22: *
   23: *       .. Scalar Arguments ..
   24: *       INTEGER            INFO, LDA, M, N
   25: *       ..
   26: *       .. Array Arguments ..
   27: *       DOUBLE PRECISION   A( LDA, * ), D( * )
   28: *       ..
   29: *
   30: *
   31: *> \par Purpose:
   32: *  =============
   33: *>
   34: *> \verbatim
   35: *>
   36: *> DLAORHR_COL_GETRFNP computes the modified LU factorization without
   37: *> pivoting of a real general M-by-N matrix A. The factorization has
   38: *> the form:
   39: *>
   40: *>     A - S = L * U,
   41: *>
   42: *> where:
   43: *>    S is a m-by-n diagonal sign matrix with the diagonal D, so that
   44: *>    D(i) = S(i,i), 1 <= i <= min(M,N). The diagonal D is constructed
   45: *>    as D(i)=-SIGN(A(i,i)), where A(i,i) is the value after performing
   46: *>    i-1 steps of Gaussian elimination. This means that the diagonal
   47: *>    element at each step of "modified" Gaussian elimination is
   48: *>    at least one in absolute value (so that division-by-zero not
   49: *>    not possible during the division by the diagonal element);
   50: *>
   51: *>    L is a M-by-N lower triangular matrix with unit diagonal elements
   52: *>    (lower trapezoidal if M > N);
   53: *>
   54: *>    and U is a M-by-N upper triangular matrix
   55: *>    (upper trapezoidal if M < N).
   56: *>
   57: *> This routine is an auxiliary routine used in the Householder
   58: *> reconstruction routine DORHR_COL. In DORHR_COL, this routine is
   59: *> applied to an M-by-N matrix A with orthonormal columns, where each
   60: *> element is bounded by one in absolute value. With the choice of
   61: *> the matrix S above, one can show that the diagonal element at each
   62: *> step of Gaussian elimination is the largest (in absolute value) in
   63: *> the column on or below the diagonal, so that no pivoting is required
   64: *> for numerical stability [1].
   65: *>
   66: *> For more details on the Householder reconstruction algorithm,
   67: *> including the modified LU factorization, see [1].
   68: *>
   69: *> This is the blocked right-looking version of the algorithm,
   70: *> calling Level 3 BLAS to update the submatrix. To factorize a block,
   71: *> this routine calls the recursive routine DLAORHR_COL_GETRFNP2.
   72: *>
   73: *> [1] "Reconstructing Householder vectors from tall-skinny QR",
   74: *>     G. Ballard, J. Demmel, L. Grigori, M. Jacquelin, H.D. Nguyen,
   75: *>     E. Solomonik, J. Parallel Distrib. Comput.,
   76: *>     vol. 85, pp. 3-31, 2015.
   77: *> \endverbatim
   78: *
   79: *  Arguments:
   80: *  ==========
   81: *
   82: *> \param[in] M
   83: *> \verbatim
   84: *>          M is INTEGER
   85: *>          The number of rows of the matrix A.  M >= 0.
   86: *> \endverbatim
   87: *>
   88: *> \param[in] N
   89: *> \verbatim
   90: *>          N is INTEGER
   91: *>          The number of columns of the matrix A.  N >= 0.
   92: *> \endverbatim
   93: *>
   94: *> \param[in,out] A
   95: *> \verbatim
   96: *>          A is DOUBLE PRECISION array, dimension (LDA,N)
   97: *>          On entry, the M-by-N matrix to be factored.
   98: *>          On exit, the factors L and U from the factorization
   99: *>          A-S=L*U; the unit diagonal elements of L are not stored.
  100: *> \endverbatim
  101: *>
  102: *> \param[in] LDA
  103: *> \verbatim
  104: *>          LDA is INTEGER
  105: *>          The leading dimension of the array A.  LDA >= max(1,M).
  106: *> \endverbatim
  107: *>
  108: *> \param[out] D
  109: *> \verbatim
  110: *>          D is DOUBLE PRECISION array, dimension min(M,N)
  111: *>          The diagonal elements of the diagonal M-by-N sign matrix S,
  112: *>          D(i) = S(i,i), where 1 <= i <= min(M,N). The elements can
  113: *>          be only plus or minus one.
  114: *> \endverbatim
  115: *>
  116: *> \param[out] INFO
  117: *> \verbatim
  118: *>          INFO is INTEGER
  119: *>          = 0:  successful exit
  120: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
  121: *> \endverbatim
  122: *>
  123: *  Authors:
  124: *  ========
  125: *
  126: *> \author Univ. of Tennessee
  127: *> \author Univ. of California Berkeley
  128: *> \author Univ. of Colorado Denver
  129: *> \author NAG Ltd.
  130: *
  131: *> \date November 2019
  132: *
  133: *> \ingroup doubleGEcomputational
  134: *
  135: *> \par Contributors:
  136: *  ==================
  137: *>
  138: *> \verbatim
  139: *>
  140: *> November 2019, Igor Kozachenko,
  141: *>                Computer Science Division,
  142: *>                University of California, Berkeley
  143: *>
  144: *> \endverbatim
  145: *
  146: *  =====================================================================
  147:       SUBROUTINE DLAORHR_COL_GETRFNP( M, N, A, LDA, D, INFO )
  148:       IMPLICIT NONE
  149: *
  150: *  -- LAPACK computational routine (version 3.9.0) --
  151: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  152: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  153: *     November 2019
  154: *
  155: *     .. Scalar Arguments ..
  156:       INTEGER            INFO, LDA, M, N
  157: *     ..
  158: *     .. Array Arguments ..
  159:       DOUBLE PRECISION   A( LDA, * ), D( * )
  160: *     ..
  161: *
  162: *  =====================================================================
  163: *
  164: *     .. Parameters ..
  165:       DOUBLE PRECISION   ONE
  166:       PARAMETER          ( ONE = 1.0D+0 )
  167: *     ..
  168: *     .. Local Scalars ..
  169:       INTEGER            IINFO, J, JB, NB
  170: *     ..
  171: *     .. External Subroutines ..
  172:       EXTERNAL           DGEMM, DLAORHR_COL_GETRFNP2, DTRSM, XERBLA
  173: *     ..
  174: *     .. External Functions ..
  175:       INTEGER            ILAENV
  176:       EXTERNAL           ILAENV
  177: *     ..
  178: *     .. Intrinsic Functions ..
  179:       INTRINSIC          MAX, MIN
  180: *     ..
  181: *     .. Executable Statements ..
  182: *
  183: *     Test the input parameters.
  184: *
  185:       INFO = 0
  186:       IF( M.LT.0 ) THEN
  187:          INFO = -1
  188:       ELSE IF( N.LT.0 ) THEN
  189:          INFO = -2
  190:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
  191:          INFO = -4
  192:       END IF
  193:       IF( INFO.NE.0 ) THEN
  194:          CALL XERBLA( 'DLAORHR_COL_GETRFNP', -INFO )
  195:          RETURN
  196:       END IF
  197: *
  198: *     Quick return if possible
  199: *
  200:       IF( MIN( M, N ).EQ.0 )
  201:      $   RETURN
  202: *
  203: *     Determine the block size for this environment.
  204: *
  205: 
  206:       NB = ILAENV( 1, 'DLAORHR_COL_GETRFNP', ' ', M, N, -1, -1 )
  207: 
  208:       IF( NB.LE.1 .OR. NB.GE.MIN( M, N ) ) THEN
  209: *
  210: *        Use unblocked code.
  211: *
  212:          CALL DLAORHR_COL_GETRFNP2( M, N, A, LDA, D, INFO )
  213:       ELSE
  214: *
  215: *        Use blocked code.
  216: *
  217:          DO J = 1, MIN( M, N ), NB
  218:             JB = MIN( MIN( M, N )-J+1, NB )
  219: *
  220: *           Factor diagonal and subdiagonal blocks.
  221: *
  222:             CALL DLAORHR_COL_GETRFNP2( M-J+1, JB, A( J, J ), LDA,
  223:      $                                 D( J ), IINFO )
  224: *
  225:             IF( J+JB.LE.N ) THEN
  226: *
  227: *              Compute block row of U.
  228: *
  229:                CALL DTRSM( 'Left', 'Lower', 'No transpose', 'Unit', JB,
  230:      $                     N-J-JB+1, ONE, A( J, J ), LDA, A( J, J+JB ),
  231:      $                     LDA )
  232:                IF( J+JB.LE.M ) THEN
  233: *
  234: *                 Update trailing submatrix.
  235: *
  236:                   CALL DGEMM( 'No transpose', 'No transpose', M-J-JB+1,
  237:      $                        N-J-JB+1, JB, -ONE, A( J+JB, J ), LDA,
  238:      $                        A( J, J+JB ), LDA, ONE, A( J+JB, J+JB ),
  239:      $                        LDA )
  240:                END IF
  241:             END IF
  242:          END DO
  243:       END IF
  244:       RETURN
  245: *
  246: *     End of DLAORHR_COL_GETRFNP
  247: *
  248:       END

CVSweb interface <joel.bertrand@systella.fr>