Annotation of rpl/lapack/lapack/dgeqr.f, revision 1.3

1.1       bertrand    1: *
                      2: *  Definition:
                      3: *  ===========
                      4: *
                      5: *       SUBROUTINE DGEQR( 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: *> DGEQR computes a QR 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 above the diagonal of the array
                     43: *>          contain the min(M,N)-by-N upper trapezoidal matrix R
                     44: *>          (R is upper triangular if M >= N);
                     45: *>          the elements below 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 QR factorization algorithm they wish. The triangular 
                    123: *> (trapezoidal) R has to be stored in the upper 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: *>                           DLATSQR or DGEQRT
                    151: *>
                    152: *>  Depending on the matrix dimensions M and N, and row and column
                    153: *>  block sizes MB and NB returned by ILAENV, DGEQR will use either
                    154: *>  DLATSQR (if the matrix is tall-and-skinny) or DGEQRT to compute
                    155: *>  the QR factorization.
                    156: *>
                    157: *> \endverbatim
                    158: *>
                    159: *  =====================================================================
                    160:       SUBROUTINE DGEQR( M, N, A, LDA, T, TSIZE, WORK, LWORK,
                    161:      $                  INFO )
                    162: *
                    163: *  -- LAPACK computational routine (version 3.7.0) --
                    164: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                    165: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. --
                    166: *     December 2016
                    167: *
                    168: *     .. Scalar Arguments ..
                    169:       INTEGER            INFO, LDA, M, N, TSIZE, LWORK
                    170: *     ..
                    171: *     .. Array Arguments ..
                    172:       DOUBLE PRECISION   A( LDA, * ), T( * ), WORK( * )
                    173: *     ..
                    174: *
                    175: *  =====================================================================
                    176: *
                    177: *     ..
                    178: *     .. Local Scalars ..
                    179:       LOGICAL            LQUERY, LMINWS, MINT, MINW
                    180:       INTEGER            MB, NB, MINTSZ, NBLCKS
                    181: *     ..
                    182: *     .. External Functions ..
                    183:       LOGICAL            LSAME
                    184:       EXTERNAL           LSAME
                    185: *     ..
                    186: *     .. External Subroutines ..
                    187:       EXTERNAL           DLATSQR, DGEQRT, XERBLA
                    188: *     ..
                    189: *     .. Intrinsic Functions ..
                    190:       INTRINSIC          MAX, MIN, MOD
                    191: *     ..
                    192: *     .. External Functions ..
                    193:       INTEGER            ILAENV
                    194:       EXTERNAL           ILAENV
                    195: *     ..
                    196: *     .. Executable Statements ..
                    197: *
                    198: *     Test the input arguments
                    199: *
                    200:       INFO = 0
                    201: *
                    202:       LQUERY = ( TSIZE.EQ.-1 .OR. TSIZE.EQ.-2 .OR.
                    203:      $           LWORK.EQ.-1 .OR. LWORK.EQ.-2 )
                    204: *
                    205:       MINT = .FALSE.
                    206:       MINW = .FALSE.
                    207:       IF( TSIZE.EQ.-2 .OR. LWORK.EQ.-2 ) THEN
                    208:         IF( TSIZE.NE.-1 ) MINT = .TRUE.
                    209:         IF( LWORK.NE.-1 ) MINW = .TRUE.
                    210:       END IF
                    211: *
                    212: *     Determine the block size
                    213: *
                    214:       IF( MIN( M, N ).GT.0 ) THEN
                    215:         MB = ILAENV( 1, 'DGEQR ', ' ', M, N, 1, -1 )
                    216:         NB = ILAENV( 1, 'DGEQR ', ' ', M, N, 2, -1 )
                    217:       ELSE
                    218:         MB = M
                    219:         NB = 1
                    220:       END IF
                    221:       IF( MB.GT.M .OR. MB.LE.N ) MB = M
                    222:       IF( NB.GT.MIN( M, N ) .OR. NB.LT.1 ) NB = 1
                    223:       MINTSZ = N + 5
                    224:       IF( MB.GT.N .AND. M.GT.N ) THEN
                    225:         IF( MOD( M - N, MB - N ).EQ.0 ) THEN
                    226:           NBLCKS = ( M - N ) / ( MB - N )
                    227:         ELSE
                    228:           NBLCKS = ( M - N ) / ( MB - N ) + 1
                    229:         END IF
                    230:       ELSE
                    231:         NBLCKS = 1
                    232:       END IF
                    233: *
                    234: *     Determine if the workspace size satisfies minimal size
                    235: *
                    236:       LMINWS = .FALSE.
                    237:       IF( ( TSIZE.LT.MAX( 1, NB*N*NBLCKS + 5 ) .OR. LWORK.LT.NB*N )
                    238:      $    .AND. ( LWORK.GE.N ) .AND. ( TSIZE.GE.MINTSZ )
                    239:      $    .AND. ( .NOT.LQUERY ) ) THEN
                    240:         IF( TSIZE.LT.MAX( 1, NB*N*NBLCKS + 5 ) ) THEN
                    241:           LMINWS = .TRUE.
                    242:           NB = 1
                    243:           MB = M
                    244:         END IF
                    245:         IF( LWORK.LT.NB*N ) THEN
                    246:           LMINWS = .TRUE.
                    247:           NB = 1
                    248:         END IF
                    249:       END IF
                    250: *
                    251:       IF( M.LT.0 ) THEN
                    252:         INFO = -1
                    253:       ELSE IF( N.LT.0 ) THEN
                    254:         INFO = -2
                    255:       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
                    256:         INFO = -4
                    257:       ELSE IF( TSIZE.LT.MAX( 1, NB*N*NBLCKS + 5 )
                    258:      $   .AND. ( .NOT.LQUERY ) .AND. ( .NOT.LMINWS ) ) THEN
                    259:         INFO = -6
                    260:       ELSE IF( ( LWORK.LT.MAX( 1, N*NB ) ) .AND. ( .NOT.LQUERY )
                    261:      $   .AND. ( .NOT.LMINWS ) ) THEN
                    262:         INFO = -8
                    263:       END IF
                    264: *
                    265:       IF( INFO.EQ.0 ) THEN
                    266:         IF( MINT ) THEN
                    267:           T( 1 ) = MINTSZ
                    268:         ELSE
                    269:           T( 1 ) = NB*N*NBLCKS + 5
                    270:         END IF
                    271:         T( 2 ) = MB
                    272:         T( 3 ) = NB
                    273:         IF( MINW ) THEN
                    274:           WORK( 1 ) = MAX( 1, N )
                    275:         ELSE
                    276:           WORK( 1 ) = MAX( 1, NB*N )
                    277:         END IF
                    278:       END IF
                    279:       IF( INFO.NE.0 ) THEN
                    280:         CALL XERBLA( 'DGEQR', -INFO )
                    281:         RETURN
                    282:       ELSE IF( LQUERY ) THEN
                    283:         RETURN
                    284:       END IF
                    285: *
                    286: *     Quick return if possible
                    287: *
                    288:       IF( MIN( M, N ).EQ.0 ) THEN
                    289:         RETURN
                    290:       END IF
                    291: *
                    292: *     The QR Decomposition
                    293: *
                    294:       IF( ( M.LE.N ) .OR. ( MB.LE.N ) .OR. ( MB.GE.M ) ) THEN
                    295:         CALL DGEQRT( M, N, NB, A, LDA, T( 6 ), NB, WORK, INFO )
                    296:       ELSE
                    297:         CALL DLATSQR( M, N, MB, NB, A, LDA, T( 6 ), NB, WORK,
                    298:      $                LWORK, INFO )
                    299:       END IF
                    300: *
                    301:       WORK( 1 ) = MAX( 1, NB*N )
                    302: *
                    303:       RETURN
                    304: *
                    305: *     End of DGEQR
                    306: *
                    307:       END

CVSweb interface <joel.bertrand@systella.fr>