SRC/dgelq.f

   1 *
   2 *  Definition:
   3 *  ===========
   4 *
   5 *       SUBROUTINE DGELQ( M, N, A, LDA, WORK1, LWORK1, WORK2, LWORK2,
   6 *                          INFO)
   7 *
   8 *       .. Scalar Arguments ..
   9 *       INTEGER           INFO, LDA, M, N, LWORK1, LWORK2
  10 *       ..
  11 *       .. Array Arguments ..
  12 *       DOUBLE PRECISION  A( LDA, * ), WORK1( * ), WORK2( * )
  13 *       ..
  14 *
  15 *
  16 *> \par Purpose:
  17 *  =============
  18 *>
  19 *> \verbatim
  20 *>
  21 *> DGELQ computes an LQ factorization of an M-by-N matrix A,
  22 *> using DLASWLQ when A is short and wide
  23 *> (N sufficiently greater than M), and otherwise DGELQT:
  24 *> A = L * Q .
  25 *> \endverbatim
  26 *
  27 *  Arguments:
  28 *  ==========
  29 *
  30 *> \param[in] M
  31 *> \verbatim
  32 *>          M is INTEGER
  33 *>          The number of rows of the matrix A.  M >= 0.
  34 *> \endverbatim
  35 *>
  36 *> \param[in] N
  37 *> \verbatim
  38 *>          N is INTEGER
  39 *>          The number of columns of the matrix A.  N >= 0.
  40 *> \endverbatim
  41 *>
  42 *> \param[in,out] A
  43 *> \verbatim
  44 *>          A is DOUBLE PRECISION array, dimension (LDA,N)
  45 *>          On entry, the M-by-N matrix A.
  46 *>          On exit, the elements on and below the diagonal of the array
  47 *>          contain the M-by-min(M,N) lower trapezoidal matrix L
  48 *>          (L is lower triangular if M <= N);
  49 *>          the elements above the diagonal are the rows of
  50 *>          blocked V representing Q (see Further Details).
  51 *> \endverbatim
  52 *>
  53 *> \param[in] LDA
  54 *> \verbatim
  55 *>          LDA is INTEGER
  56 *>          The leading dimension of the array A.  LDA >= max(1,M).
  57 *> \endverbatim
  58 *>
  59 *> \param[out] WORK1
  60 *> \verbatim
  61 *>          WORK1 is DOUBLE PRECISION array, dimension (MAX(1,LWORK1))
  62 *>          WORK1 contains part of the data structure used to store Q.
  63 *>          WORK1(1): algorithm type = 1, to indicate output from
  64 *>                    DLASWLQ or DGELQT
  65 *>          WORK1(2): optimum size of WORK1
  66 *>          WORK1(3): minimum size of WORK1
  67 *>          WORK1(4): horizontal block size
  68 *>          WORK1(5): vertical block size
  69 *>          WORK1(6:LWORK1): data structure needed for Q, computed by
  70 *>                           DLASWLQ or DGELQT
  71 *> \endverbatim
  72 *>
  73 *> \param[in] LWORK1
  74 *> \verbatim
  75 *>          LWORK1 is INTEGER
  76 *>          The dimension of the array WORK1.
  77 *>          If LWORK1 = -1, then a query is assumed. In this case the
  78 *>          routine calculates the optimal size of WORK1 and
  79 *>          returns this value in WORK1(2),  and calculates the minimum
  80 *>          size of WORK1 and returns this value in WORK1(3).
  81 *>          No error message related to LWORK1 is issued by XERBLA when
  82 *>          LWORK1 = -1.
  83 *> \endverbatim
  84 *>
  85 *> \param[out] WORK2
  86 *> \verbatim
  87 *>         (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK2))
  88 *>
  89 *> \endverbatim
  90 *> \param[in] LWORK2
  91 *> \verbatim
  92 *>          LWORK2 is INTEGER
  93 *>          The dimension of the array WORK2.
  94 *>          If LWORK2 = -1, then a query is assumed. In this case the
  95 *>          routine calculates the optimal size of WORK2 and
  96 *>          returns this value in WORK2(1), and calculates the minimum
  97 *>          size of WORK2 and returns this value in WORK2(2).
  98 *>          No error message related to LWORK2 is issued by XERBLA when
  99 *>          LWORK2 = -1.
 100 *> \endverbatim
 101 *>
 102 *> \param[out] INFO
 103 *> \verbatim
 104 *>          INFO is INTEGER
 105 *>          = 0:  successful exit
 106 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
 107 *> \endverbatim
 108 *
 109 *  Authors:
 110 *  ========
 111 *
 112 *> \author Univ. of Tennessee
 113 *> \author Univ. of California Berkeley
 114 *> \author Univ. of Colorado Denver
 115 *> \author NAG Ltd.
 116 *
 117 *> \par Further Details:
 118 *  =====================
 119 *>
 120 *> \verbatim
 121 *>  Depending on the matrix dimensions M and N, and row and column
 122 *>  block sizes MB and NB returned by ILAENV, GELQ will use either
 123 *>  LASWLQ(if the matrix is short-and-wide) or GELQT to compute
 124 *>  the LQ decomposition.
 125 *>  The output of LASWLQ or GELQT representing Q is stored in A and in
 126 *>  array WORK1(6:LWORK1) for later use.
 127 *>  WORK1(2:5) contains the matrix dimensions M,N and block sizes MB, NB
 128 *>  which are needed to interpret A and WORK1(6:LWORK1) for later use.
 129 *>  WORK1(1)=1 indicates that the code needed to take WORK1(2:5) and
 130 *>  decide whether LASWLQ or GELQT was used is the same as used below in
 131 *>  GELQ. For a detailed description of A and WORK1(6:LWORK1), see
 132 *>  Further Details in LASWLQ or GELQT.
 133 *> \endverbatim
 134 *>
 135 *>
 136 *  =====================================================================
 137       SUBROUTINE DGELQ( M, N, A, LDA, WORK1, LWORK1, WORK2, LWORK2,
 138      $   INFO)
 139 *
 140 *  -- LAPACK computational routine (version 3.5.0) --
 141 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 142 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. --
 143 *     November 2013
 144 *
 145 *     .. Scalar Arguments ..
 146       INTEGER           INFO, LDA, M, N, LWORK1, LWORK2
 147 *     ..
 148 *     .. Array Arguments ..
 149       DOUBLE PRECISION  A( LDA, * ), WORK1( * ), WORK2( * )
 150 *     ..
 151 *
 152 *  =====================================================================
 153 *
 154 *     ..
 155 *     .. Local Scalars ..
 156       LOGICAL    LQUERY, LMINWS
 157       INTEGER    MB, NB, I, II, KK, MINLW1, NBLCKS
 158 *     ..
 159 *     .. EXTERNAL FUNCTIONS ..
 160       LOGICAL            LSAME
 161       EXTERNAL           LSAME
 162 *     .. EXTERNAL SUBROUTINES ..
 163       EXTERNAL           DGELQT, DLASWLQ, XERBLA
 164 *     .. INTRINSIC FUNCTIONS ..
 165       INTRINSIC          MAX, MIN, MOD
 166 *     ..
 167 *     .. EXTERNAL FUNCTIONS ..
 168       INTEGER            ILAENV
 169       EXTERNAL           ILAENV
 170 *     ..
 171 *     .. EXECUTABLE STATEMENTS ..
 172 *
 173 *     TEST THE INPUT ARGUMENTS
 174 *
 175       INFO = 0
 176 *
 177       LQUERY = ( LWORK1.EQ.-1 .OR. LWORK2.EQ.-1 )
 178 *
 179 *     Determine the block size
 180 *
 181       IF ( MIN(M,N).GT.0 ) THEN
 182         MB = ILAENV( 1, 'DGELQ ', ' ', M, N, 1, -1)
 183         NB = ILAENV( 1, 'DGELQ ', ' ', M, N, 2, -1)
 184       ELSE
 185         MB = 1
 186         NB = N
 187       END IF
 188       IF( MB.GT.MIN(M,N).OR.MB.LT.1) MB = 1
 189       IF( NB.GT.N.OR.NB.LE.M) NB = N
 190       MINLW1 = M + 5
 191       IF ((NB.GT.M).AND.(N.GT.M)) THEN
 192         IF(MOD(N-M, NB-M).EQ.0) THEN
 193           NBLCKS = (N-M)/(NB-M)
 194         ELSE
 195           NBLCKS = (N-M)/(NB-M) + 1
 196         END IF
 197       ELSE
 198         NBLCKS = 1
 199       END IF
 200 *
 201 *     Determine if the workspace size satisfies minimum size
 202 *
 203       LMINWS = .FALSE.
 204       IF((LWORK1.LT.MAX(1,MB*M*NBLCKS+5)
 205      $    .OR.(LWORK2.LT.MB*M)).AND.(LWORK2.GE.M).AND.(LWORK1.GE.M+5)
 206      $    .AND.(.NOT.LQUERY)) THEN
 207         IF (LWORK1.LT.MAX(1,MB*M*NBLCKS+5)) THEN
 208             LMINWS = .TRUE.
 209             MB = 1
 210         END IF
 211         IF (LWORK1.LT.MAX(1,M*NBLCKS+5)) THEN
 212             LMINWS = .TRUE.
 213             NB = N
 214         END IF
 215         IF (LWORK2.LT.MB*M) THEN
 216             LMINWS = .TRUE.
 217             MB = 1
 218         END IF
 219       END IF
 220 *
 221       IF( M.LT.0 ) THEN
 222         INFO = -1
 223       ELSE IF( N.LT.0 ) THEN
 224         INFO = -2
 225       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
 226         INFO = -4
 227       ELSE IF( LWORK1.LT.MAX( 1, MB*M*NBLCKS+5 )
 228      $   .AND.(.NOT.LQUERY).AND. (.NOT.LMINWS)) THEN
 229         INFO = -6
 230       ELSE IF( (LWORK2.LT.MAX(1,M*MB)).AND.(.NOT.LQUERY)
 231      $   .AND.(.NOT.LMINWS) ) THEN
 232         INFO = -8
 233       END IF
 234 *
 235       IF( INFO.EQ.0)  THEN
 236         WORK1(1) = 1
 237         WORK1(2) = MB*M*NBLCKS+5
 238         WORK1(3) = MINLW1
 239         WORK1(4) = MB
 240         WORK1(5) = NB
 241         WORK2(1) = MB * M
 242         WORK2(2) = M
 243       END IF
 244       IF( INFO.NE.0 ) THEN
 245         CALL XERBLA( 'DGELQ', -INFO )
 246         RETURN
 247       ELSE IF (LQUERY) THEN
 248        RETURN
 249       END IF
 250 *
 251 *     Quick return if possible
 252 *
 253       IF( MIN(M,N).EQ.0 ) THEN
 254           RETURN
 255       END IF
 256 *
 257 *     The LQ Decomposition
 258 *
 259       IF((N.LE.M).OR.(NB.LE.M).OR.(NB.GE.N)) THEN
 260         CALL DGELQT( M, N, MB, A, LDA, WORK1(6), MB, WORK2, INFO)
 261       ELSE
 262         CALL DLASWLQ( M, N, MB, NB, A, LDA, WORK1(6), MB, WORK2,
 263      $                    LWORK2, INFO)
 264       END IF
 265       RETURN
 266 *
 267 *     End of DGELQ
 268 *
 269       END