File:  [local] / rpl / lapack / lapack / dgelq.f
Revision 1.4: download - view: text, annotated - select for diffs - revision graph
Thu May 21 21:45:56 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 DGELQ
    2: *
    3: *  Definition:
    4: *  ===========
    5: *
    6: *       SUBROUTINE DGELQ( M, N, A, LDA, T, TSIZE, WORK, LWORK,
    7: *                         INFO )
    8: *
    9: *       .. Scalar Arguments ..
   10: *       INTEGER           INFO, LDA, M, N, TSIZE, LWORK
   11: *       ..
   12: *       .. Array Arguments ..
   13: *       DOUBLE PRECISION  A( LDA, * ), T( * ), WORK( * )
   14: *       ..
   15: *
   16: *
   17: *> \par Purpose:
   18: *  =============
   19: *>
   20: *> \verbatim
   21: *>
   22: *> DGELQ computes an LQ factorization of a real M-by-N matrix A:
   23: *>
   24: *>    A = ( L 0 ) *  Q
   25: *>
   26: *> where:
   27: *>
   28: *>    Q is a N-by-N orthogonal matrix;
   29: *>    L is an lower-triangular M-by-M matrix;
   30: *>    0 is a M-by-(N-M) zero matrix, if M < N.
   31: *>
   32: *> \endverbatim
   33: *
   34: *  Arguments:
   35: *  ==========
   36: *
   37: *> \param[in] M
   38: *> \verbatim
   39: *>          M is INTEGER
   40: *>          The number of rows of the matrix A.  M >= 0.
   41: *> \endverbatim
   42: *>
   43: *> \param[in] N
   44: *> \verbatim
   45: *>          N is INTEGER
   46: *>          The number of columns of the matrix A.  N >= 0.
   47: *> \endverbatim
   48: *>
   49: *> \param[in,out] A
   50: *> \verbatim
   51: *>          A is DOUBLE PRECISION array, dimension (LDA,N)
   52: *>          On entry, the M-by-N matrix A.
   53: *>          On exit, the elements on and below the diagonal of the array
   54: *>          contain the M-by-min(M,N) lower trapezoidal matrix L
   55: *>          (L is lower triangular if M <= N);
   56: *>          the elements above the diagonal are used to store part of the 
   57: *>          data structure to represent Q.
   58: *> \endverbatim
   59: *>
   60: *> \param[in] LDA
   61: *> \verbatim
   62: *>          LDA is INTEGER
   63: *>          The leading dimension of the array A.  LDA >= max(1,M).
   64: *> \endverbatim
   65: *>
   66: *> \param[out] T
   67: *> \verbatim
   68: *>          T is DOUBLE PRECISION array, dimension (MAX(5,TSIZE))
   69: *>          On exit, if INFO = 0, T(1) returns optimal (or either minimal 
   70: *>          or optimal, if query is assumed) TSIZE. See TSIZE for details.
   71: *>          Remaining T contains part of the data structure used to represent Q.
   72: *>          If one wants to apply or construct Q, then one needs to keep T 
   73: *>          (in addition to A) and pass it to further subroutines.
   74: *> \endverbatim
   75: *>
   76: *> \param[in] TSIZE
   77: *> \verbatim
   78: *>          TSIZE is INTEGER
   79: *>          If TSIZE >= 5, the dimension of the array T.
   80: *>          If TSIZE = -1 or -2, then a workspace query is assumed. The routine
   81: *>          only calculates the sizes of the T and WORK arrays, returns these
   82: *>          values as the first entries of the T and WORK arrays, and no error
   83: *>          message related to T or WORK is issued by XERBLA.
   84: *>          If TSIZE = -1, the routine calculates optimal size of T for the 
   85: *>          optimum performance and returns this value in T(1).
   86: *>          If TSIZE = -2, the routine calculates minimal size of T and 
   87: *>          returns this value in T(1).
   88: *> \endverbatim
   89: *>
   90: *> \param[out] WORK
   91: *> \verbatim
   92: *>          (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK))
   93: *>          On exit, if INFO = 0, WORK(1) contains optimal (or either minimal
   94: *>          or optimal, if query was assumed) LWORK.
   95: *>          See LWORK for details.
   96: *> \endverbatim
   97: *>
   98: *> \param[in] LWORK
   99: *> \verbatim
  100: *>          LWORK is INTEGER
  101: *>          The dimension of the array WORK.
  102: *>          If LWORK = -1 or -2, then a workspace query is assumed. The routine
  103: *>          only calculates the sizes of the T and WORK arrays, returns these
  104: *>          values as the first entries of the T and WORK arrays, and no error
  105: *>          message related to T or WORK is issued by XERBLA.
  106: *>          If LWORK = -1, the routine calculates optimal size of WORK for the
  107: *>          optimal performance and returns this value in WORK(1).
  108: *>          If LWORK = -2, the routine calculates minimal size of WORK and 
  109: *>          returns this value in WORK(1).
  110: *> \endverbatim
  111: *>
  112: *> \param[out] INFO
  113: *> \verbatim
  114: *>          INFO is INTEGER
  115: *>          = 0:  successful exit
  116: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
  117: *> \endverbatim
  118: *
  119: *  Authors:
  120: *  ========
  121: *
  122: *> \author Univ. of Tennessee
  123: *> \author Univ. of California Berkeley
  124: *> \author Univ. of Colorado Denver
  125: *> \author NAG Ltd.
  126: *
  127: *> \par Further Details
  128: *  ====================
  129: *>
  130: *> \verbatim
  131: *>
  132: *> The goal of the interface is to give maximum freedom to the developers for
  133: *> creating any LQ factorization algorithm they wish. The triangular 
  134: *> (trapezoidal) L has to be stored in the lower part of A. The lower part of A
  135: *> and the array T can be used to store any relevant information for applying or
  136: *> constructing the Q factor. The WORK array can safely be discarded after exit.
  137: *>
  138: *> Caution: One should not expect the sizes of T and WORK to be the same from one 
  139: *> LAPACK implementation to the other, or even from one execution to the other.
  140: *> A workspace query (for T and WORK) is needed at each execution. However, 
  141: *> for a given execution, the size of T and WORK are fixed and will not change 
  142: *> from one query to the next.
  143: *>
  144: *> \endverbatim
  145: *>
  146: *> \par Further Details particular to this LAPACK implementation:
  147: *  ==============================================================
  148: *>
  149: *> \verbatim
  150: *>
  151: *> These details are particular for this LAPACK implementation. Users should not 
  152: *> take them for granted. These details may change in the future, and are not likely
  153: *> true for another LAPACK implementation. These details are relevant if one wants
  154: *> to try to understand the code. They are not part of the interface.
  155: *>
  156: *> In this version,
  157: *>
  158: *>          T(2): row block size (MB)
  159: *>          T(3): column block size (NB)
  160: *>          T(6:TSIZE): data structure needed for Q, computed by
  161: *>                           DLASWLQ or DGELQT
  162: *>
  163: *>  Depending on the matrix dimensions M and N, and row and column
  164: *>  block sizes MB and NB returned by ILAENV, DGELQ will use either
  165: *>  DLASWLQ (if the matrix is short-and-wide) or DGELQT to compute
  166: *>  the LQ factorization.
  167: *> \endverbatim
  168: *>
  169: *  =====================================================================
  170:       SUBROUTINE DGELQ( M, N, A, LDA, T, TSIZE, WORK, LWORK,
  171:      $                  INFO )
  172: *
  173: *  -- LAPACK computational routine (version 3.9.0) --
  174: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  175: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. --
  176: *     November 2019
  177: *
  178: *     .. Scalar Arguments ..
  179:       INTEGER            INFO, LDA, M, N, TSIZE, LWORK
  180: *     ..
  181: *     .. Array Arguments ..
  182:       DOUBLE PRECISION   A( LDA, * ), T( * ), WORK( * )
  183: *     ..
  184: *
  185: *  =====================================================================
  186: *
  187: *     ..
  188: *     .. Local Scalars ..
  189:       LOGICAL            LQUERY, LMINWS, MINT, MINW
  190:       INTEGER            MB, NB, MINTSZ, NBLCKS
  191: *     ..
  192: *     .. External Functions ..
  193:       LOGICAL            LSAME
  194:       EXTERNAL           LSAME
  195: *     ..
  196: *     .. External Subroutines ..
  197:       EXTERNAL           DGELQT, DLASWLQ, XERBLA
  198: *     ..
  199: *     .. Intrinsic Functions ..
  200:       INTRINSIC          MAX, MIN, MOD
  201: *     ..
  202: *     .. External Functions ..
  203:       INTEGER            ILAENV
  204:       EXTERNAL           ILAENV
  205: *     ..
  206: *     .. Executable Statements ..
  207: *
  208: *     Test the input arguments
  209: *
  210:       INFO = 0
  211: *
  212:       LQUERY = ( TSIZE.EQ.-1 .OR. TSIZE.EQ.-2 .OR.
  213:      $           LWORK.EQ.-1 .OR. LWORK.EQ.-2 )
  214: *
  215:       MINT = .FALSE.
  216:       MINW = .FALSE.
  217:       IF( TSIZE.EQ.-2 .OR. LWORK.EQ.-2 ) THEN
  218:         IF( TSIZE.NE.-1 ) MINT = .TRUE.
  219:         IF( LWORK.NE.-1 ) MINW = .TRUE.
  220:       END IF
  221: *
  222: *     Determine the block size
  223: *
  224:       IF( MIN( M, N ).GT.0 ) THEN
  225:         MB = ILAENV( 1, 'DGELQ ', ' ', M, N, 1, -1 )
  226:         NB = ILAENV( 1, 'DGELQ ', ' ', M, N, 2, -1 )
  227:       ELSE
  228:         MB = 1
  229:         NB = N
  230:       END IF
  231:       IF( MB.GT.MIN( M, N ) .OR. MB.LT.1 ) MB = 1
  232:       IF( NB.GT.N .OR. NB.LE.M ) NB = N
  233:       MINTSZ = M + 5
  234:       IF ( NB.GT.M .AND. N.GT.M ) THEN
  235:         IF( MOD( N - M, NB - M ).EQ.0 ) THEN
  236:           NBLCKS = ( N - M ) / ( NB - M )
  237:         ELSE
  238:           NBLCKS = ( N - M ) / ( NB - M ) + 1
  239:         END IF
  240:       ELSE
  241:         NBLCKS = 1
  242:       END IF
  243: *
  244: *     Determine if the workspace size satisfies minimal size
  245: *
  246:       LMINWS = .FALSE.
  247:       IF( ( TSIZE.LT.MAX( 1, MB*M*NBLCKS + 5 ) .OR. LWORK.LT.MB*M )
  248:      $    .AND. ( LWORK.GE.M ) .AND. ( TSIZE.GE.MINTSZ )
  249:      $    .AND. ( .NOT.LQUERY ) ) THEN
  250:         IF( TSIZE.LT.MAX( 1, MB*M*NBLCKS + 5 ) ) THEN
  251:             LMINWS = .TRUE.
  252:             MB = 1
  253:             NB = N
  254:         END IF
  255:         IF( LWORK.LT.MB*M ) THEN
  256:             LMINWS = .TRUE.
  257:             MB = 1
  258:         END IF
  259:       END IF
  260: *
  261:       IF( M.LT.0 ) THEN
  262:         INFO = -1
  263:       ELSE IF( N.LT.0 ) THEN
  264:         INFO = -2
  265:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
  266:         INFO = -4
  267:       ELSE IF( TSIZE.LT.MAX( 1, MB*M*NBLCKS + 5 )
  268:      $   .AND. ( .NOT.LQUERY ) .AND. ( .NOT.LMINWS ) ) THEN
  269:         INFO = -6
  270:       ELSE IF( ( LWORK.LT.MAX( 1, M*MB ) ) .AND .( .NOT.LQUERY )
  271:      $   .AND. ( .NOT.LMINWS ) ) THEN
  272:         INFO = -8
  273:       END IF
  274: *
  275:       IF( INFO.EQ.0 ) THEN
  276:         IF( MINT ) THEN
  277:           T( 1 ) = MINTSZ
  278:         ELSE
  279:           T( 1 ) = MB*M*NBLCKS + 5
  280:         END IF
  281:         T( 2 ) = MB
  282:         T( 3 ) = NB
  283:         IF( MINW ) THEN
  284:           WORK( 1 ) = MAX( 1, N )
  285:         ELSE
  286:           WORK( 1 ) = MAX( 1, MB*M )
  287:         END IF
  288:       END IF
  289:       IF( INFO.NE.0 ) THEN
  290:         CALL XERBLA( 'DGELQ', -INFO )
  291:         RETURN
  292:       ELSE IF( LQUERY ) THEN
  293:         RETURN
  294:       END IF
  295: *
  296: *     Quick return if possible
  297: *
  298:       IF( MIN( M, N ).EQ.0 ) THEN
  299:         RETURN
  300:       END IF
  301: *
  302: *     The LQ Decomposition
  303: *
  304:       IF( ( N.LE.M ) .OR. ( NB.LE.M ) .OR. ( NB.GE.N ) ) THEN
  305:         CALL DGELQT( M, N, MB, A, LDA, T( 6 ), MB, WORK, INFO )
  306:       ELSE
  307:         CALL DLASWLQ( M, N, MB, NB, A, LDA, T( 6 ), MB, WORK,
  308:      $                LWORK, INFO )
  309:       END IF
  310: *
  311:       WORK( 1 ) = MAX( 1, MB*M )
  312: *
  313:       RETURN
  314: *
  315: *     End of DGELQ
  316: *
  317:       END

CVSweb interface <joel.bertrand@systella.fr>