SRC/cgemqr.f

   1 *
   2 *  Definition:
   3 *  ===========
   4 *
   5 *      SUBROUTINE CGEMQR( SIDE, TRANS, M, N, K, A, LDA, WORK1,
   6 *     $                     LWORK1, C, LDC, WORK2, LWORK2, INFO )
   7 *
   8 *
   9 *     .. Scalar Arguments ..
  10 *      CHARACTER         SIDE, TRANS
  11 *      INTEGER           INFO, LDA, M, N, K, LDT, LWORK1, LWORK2, LDC
  12 *     ..
  13 *     .. Array Arguments ..
  14 *      COMPLEX        A( LDA, * ), WORK1( * ), C(LDC, * ),
  15 *     $                  WORK2( * )
  16 *> \par Purpose:
  17 *  =============
  18 *>
  19 *> \verbatim
  20 *>
  21 *>      CGEMQR overwrites the general real M-by-N matrix C with
  22 *>
  23 *>
  24 *>                      SIDE = 'L'     SIDE = 'R'
  25 *>      TRANS = 'N':      Q * C          C * Q
  26 *>      TRANS = 'T':      Q**T * C       C * Q**T
  27 *>      where Q is a complex orthogonal matrix defined as the product
  28 *>      of blocked elementary reflectors computed by tall skinny
  29 *>      QR factorization (CGEQR)
  30 *> \endverbatim
  31 *
  32 *  Arguments:
  33 *  ==========
  34 *
  35 *> \param[in] SIDE
  36 *>          SIDE is CHARACTER*1
  37 *>          = 'L': apply Q or Q**T from the Left;
  38 *>          = 'R': apply Q or Q**T from the Right.
  39 *>
  40 *> \param[in] TRANS
  41 *>          TRANS is CHARACTER*1
  42 *>          = 'N':  No transpose, apply Q;
  43 *>          = 'T':  Transpose, apply Q**T.
  44 *> \param[in] M
  45 *> \verbatim
  46 *>          M is INTEGER
  47 *>          The number of rows of the matrix A.  M >=0.
  48 *> \endverbatim
  49 *>
  50 *> \param[in] N
  51 *> \verbatim
  52 *>          N is INTEGER
  53 *>          The number of columns of the matrix C. M >= N >= 0.
  54 *> \endverbatim
  55 *>
  56 *> \param[in] K
  57 *> \verbatim
  58 *>          K is INTEGER
  59 *>          The number of elementary reflectors whose product defines
  60 *>          the matrix Q.
  61 *>          N >= K >= 0;
  62 *>
  63 *> \endverbatim
  64 *>
  65 *> \param[in,out] A
  66 *> \verbatim
  67 *>          A is COMPLEX array, dimension (LDA,K)
  68 *>          The i-th column must contain the vector which defines the
  69 *>          blockedelementary reflector H(i), for i = 1,2,...,k, as
  70 *>          returned by DGETSQR in the first k columns of
  71 *>          its array argument A.
  72 *> \endverbatim
  73 *>
  74 *> \param[in] LDA
  75 *> \verbatim
  76 *>          LDA is INTEGER
  77 *>          The leading dimension of the array A.
  78 *>          If SIDE = 'L', LDA >= max(1,M);
  79 *>          if SIDE = 'R', LDA >= max(1,N).
  80 *> \endverbatim
  81 *>
  82 *> \param[in] WORK1
  83 *> \verbatim
  84 *>          WORK1 is COMPLEX array, dimension (MAX(1,LWORK1)) as
  85 *>          it is returned by GEQR.
  86 *> \endverbatim
  87 *>
  88 *> \param[in] LWORK1
  89 *> \verbatim
  90 *>          LWORK1 is INTEGER
  91 *>          The dimension of the array WORK1.
  92 *> \endverbatim
  93 *>
  94 *> \param[in,out] C
  95 *>          C is COMPLEX array, dimension (LDC,N)
  96 *>          On entry, the M-by-N matrix C.
  97 *>          On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q.
  98 *>
  99 *> \param[in] LDC
 100 *>          LDC is INTEGER
 101 *>          The leading dimension of the array C. LDC >= max(1,M).
 102 *>
 103 *> \param[out] WORK2
 104 *> \verbatim
 105 *>         (workspace) COMPLEX array, dimension (MAX(1,LWORK2))
 106 *>
 107 *> \endverbatim
 108 *> \param[in] LWORK2
 109 *> \verbatim
 110 *>          LWORK2 is INTEGER
 111 *>          The dimension of the array WORK2.
 112 *>          If LWORK2 = -1, then a workspace query is assumed; the routine
 113 *>          only calculates the optimal size of the WORK2 array, returns
 114 *>          this value as the third entry of the WORK2 array (WORK2(1)),
 115 *>          and no error message related to LWORK2 is issued by XERBLA.
 116 *>
 117 *> \endverbatim
 118 *> \param[out] INFO
 119 *> \verbatim
 120 *>          INFO is INTEGER
 121 *>          = 0:  successful exit
 122 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
 123 *> \endverbatim
 124 *
 125 *  Authors:
 126 *  ========
 127 *
 128 *> \author Univ. of Tennessee
 129 *> \author Univ. of California Berkeley
 130 *> \author Univ. of Colorado Denver
 131 *> \author NAG Ltd.
 132 *
 133 *> \par Further Details:
 134 *  =====================
 135 *>
 136 *> \verbatim
 137 *>  Depending on the matrix dimensions M and N, and row and column
 138 *>  block sizes MB and NB returned by ILAENV, GEQR will use either
 139 *>  LATSQR (if the matrix is tall-and-skinny) or GEQRT to compute
 140 *>  the QR decomposition.
 141 *>  The output of LATSQR or GEQRT representing Q is stored in A and in
 142 *>  array WORK1(6:LWORK1) for later use.
 143 *>  WORK1(2:5) contains the matrix dimensions M,N and block sizes MB,NB
 144 *>  which are needed to interpret A and WORK1(6:LWORK1) for later use.
 145 *>  WORK1(1)=1 indicates that the code needed to take WORK1(2:5) and
 146 *>  decide whether LATSQR or GEQRT was used is the same as used below in
 147 *>  GEQR. For a detailed description of A and WORK1(6:LWORK1), see
 148 *>  Further Details in LATSQR or GEQRT.
 149 *> \endverbatim
 150 *>
 151 *  =====================================================================
 152       SUBROUTINE CGEMQR( SIDE, TRANS, M, N, K, A, LDA, WORK1, LWORK1,
 153      $        C, LDC, WORK2, LWORK2, INFO )
 154 *
 155 *  -- LAPACK computational routine (version 3.5.0) --
 156 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 157 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 158 *     November 2013
 159 *
 160 *     .. Scalar Arguments ..
 161       CHARACTER         SIDE, TRANS
 162       INTEGER           INFO, LDA, M, N, K, LWORK1, LWORK2, LDC
 163 *     ..
 164 *     .. Array Arguments ..
 165       COMPLEX           A( LDA, * ), WORK1( * ), C(LDC, * ),
 166      $               WORK2( * )
 167 *     ..
 168 *
 169 * =====================================================================
 170 *
 171 *     ..
 172 *     .. Local Scalars ..
 173       LOGICAL    LEFT, RIGHT, TRAN, NOTRAN, LQUERY
 174       INTEGER    MB, NB, I, II, KK, LW, NBLCKS, MN
 175 *     ..
 176 *     .. External Functions ..
 177       LOGICAL            LSAME
 178       EXTERNAL           LSAME
 179 *     .. External Subroutines ..
 180       EXTERNAL           CGEMQRT, CLAMTSQR, XERBLA
 181 *     .. Intrinsic Functions ..
 182       INTRINSIC          INT, MAX, MIN, MOD
 183 *     ..
 184 *     .. Executable Statements ..
 185 *
 186 *     Test the input arguments
 187 *
 188       LQUERY  = LWORK2.LT.0
 189       NOTRAN  = LSAME( TRANS, 'N' )
 190       TRAN    = LSAME( TRANS, 'C' )
 191       LEFT    = LSAME( SIDE, 'L' )
 192       RIGHT   = LSAME( SIDE, 'R' )
 193 *
 194       MB = INT(WORK1(4))
 195       NB = INT(WORK1(5))
 196       IF(LEFT) THEN
 197         LW = N * NB
 198         MN = M
 199       ELSE IF(RIGHT) THEN
 200         LW = MB * NB
 201         MN = N
 202       END IF
 203 *
 204       IF ((MB.GT.K).AND.(MN.GT.K)) THEN
 205           IF(MOD(MN-K, MB-K).EQ.0) THEN
 206              NBLCKS = (MN-K)/(MB-K)
 207           ELSE
 208              NBLCKS = (MN-K)/(MB-K) + 1
 209           END IF
 210       ELSE
 211         NBLCKS = 1
 212       END IF
 213 *
 214       INFO = 0
 215       IF( .NOT.LEFT .AND. .NOT.RIGHT ) THEN
 216          INFO = -1
 217       ELSE IF( .NOT.TRAN .AND. .NOT.NOTRAN ) THEN
 218          INFO = -2
 219       ELSE IF( M.LT.0 ) THEN
 220         INFO = -3
 221       ELSE IF( N.LT.0 ) THEN
 222         INFO = -4
 223       ELSE IF( K.LT.0 ) THEN
 224         INFO = -5
 225       ELSE IF( LDA.LT.MAX( 1, K ) ) THEN
 226         INFO = -7
 227       ELSE IF( LWORK1.LT.MAX( 1, NB*K*NBLCKS+5 )) THEN
 228         INFO = -9
 229       ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
 230          INFO = -11
 231       ELSE IF(( LWORK2.LT.MAX(1,LW)).AND.(.NOT.LQUERY)) THEN
 232         INFO = -13
 233       END IF
 234 *
 235 *     Determine the block size if it is tall skinny or short and wide
 236 *
 237       IF( INFO.EQ.0)  THEN
 238           WORK2(1) = LW
 239       END IF
 240 *
 241       IF( INFO.NE.0 ) THEN
 242         CALL XERBLA( 'CGEMQR', -INFO )
 243         RETURN
 244       ELSE IF (LQUERY) THEN
 245        RETURN
 246       END IF
 247 *
 248 *     Quick return if possible
 249 *
 250       IF( MIN(M,N,K).EQ.0 ) THEN
 251         RETURN
 252       END IF
 253 *
 254       IF((LEFT.AND.M.LE.K).OR.(RIGHT.AND.N.LE.K).OR.(MB.LE.K).OR.
 255      $   (MB.GE.MAX(M,N,K))) THEN
 256         CALL CGEMQRT( SIDE, TRANS, M, N, K, NB, A, LDA,
 257      $        WORK1(6), NB, C, LDC, WORK2, INFO)
 258       ELSE
 259         CALL CLAMTSQR( SIDE, TRANS, M, N, K, MB, NB, A, LDA, WORK1(6),
 260      $      NB, C, LDC, WORK2, LWORK2, INFO )
 261       END IF
 262 *
 263       WORK2(1) = LW
 264       RETURN
 265 *
 266 *     End of CGEMQR
 267 *
 268       END