File:  [local] / rpl / lapack / lapack / dgelq.f
Revision 1.2: download - view: text, annotated - select for diffs - revision graph
Sat Jun 17 11:06:15 2017 UTC (6 years, 10 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_27, rpl-4_1_26, HEAD
Cohérence.

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

CVSweb interface <joel.bertrand@systella.fr>