SRC/cgeqrt.f

   1 *> \brief \b CGEQRT
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download CGEQRT + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cgeqrt.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cgeqrt.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cgeqrt.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE CGEQRT( M, N, NB, A, LDA, T, LDT, WORK, INFO )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       INTEGER INFO, LDA, LDT, M, N, NB
  25 *       ..
  26 *       .. Array Arguments ..
  27 *       COMPLEX A( LDA, * ), T( LDT, * ), WORK( * )
  28 *       ..
  29 *
  30 *
  31 *> \par Purpose:
  32 *  =============
  33 *>
  34 *> \verbatim
  35 *>
  36 *> CGEQRT computes a blocked QR factorization of a complex M-by-N matrix A
  37 *> using the compact WY representation of Q.
  38 *> \endverbatim
  39 *
  40 *  Arguments:
  41 *  ==========
  42 *
  43 *> \param[in] M
  44 *> \verbatim
  45 *>          M is INTEGER
  46 *>          The number of rows of the matrix A.  M >= 0.
  47 *> \endverbatim
  48 *>
  49 *> \param[in] N
  50 *> \verbatim
  51 *>          N is INTEGER
  52 *>          The number of columns of the matrix A.  N >= 0.
  53 *> \endverbatim
  54 *>
  55 *> \param[in] NB
  56 *> \verbatim
  57 *>          NB is INTEGER
  58 *>          The block size to be used in the blocked QR.  MIN(M,N) >= NB >= 1.
  59 *> \endverbatim
  60 *>
  61 *> \param[in,out] A
  62 *> \verbatim
  63 *>          A is COMPLEX array, dimension (LDA,N)
  64 *>          On entry, the M-by-N matrix A.
  65 *>          On exit, the elements on and above the diagonal of the array
  66 *>          contain the min(M,N)-by-N upper trapezoidal matrix R (R is
  67 *>          upper triangular if M >= N); the elements below the diagonal
  68 *>          are the columns of V.
  69 *> \endverbatim
  70 *>
  71 *> \param[in] LDA
  72 *> \verbatim
  73 *>          LDA is INTEGER
  74 *>          The leading dimension of the array A.  LDA >= max(1,M).
  75 *> \endverbatim
  76 *>
  77 *> \param[out] T
  78 *> \verbatim
  79 *>          T is COMPLEX array, dimension (LDT,MIN(M,N))
  80 *>          The upper triangular block reflectors stored in compact form
  81 *>          as a sequence of upper triangular blocks.  See below
  82 *>          for further details.
  83 *> \endverbatim
  84 *>
  85 *> \param[in] LDT
  86 *> \verbatim
  87 *>          LDT is INTEGER
  88 *>          The leading dimension of the array T.  LDT >= NB.
  89 *> \endverbatim
  90 *>
  91 *> \param[out] WORK
  92 *> \verbatim
  93 *>          WORK is COMPLEX array, dimension (NB*N)
  94 *> \endverbatim
  95 *>
  96 *> \param[out] INFO
  97 *> \verbatim
  98 *>          INFO is INTEGER
  99 *>          = 0:  successful exit
 100 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
 101 *> \endverbatim
 102 *
 103 *  Authors:
 104 *  ========
 105 *
 106 *> \author Univ. of Tennessee
 107 *> \author Univ. of California Berkeley
 108 *> \author Univ. of Colorado Denver
 109 *> \author NAG Ltd.
 110 *
 111 *> \date November 2013
 112 *
 113 *> \ingroup complexGEcomputational
 114 *
 115 *> \par Further Details:
 116 *  =====================
 117 *>
 118 *> \verbatim
 119 *>
 120 *>  The matrix V stores the elementary reflectors H(i) in the i-th column
 121 *>  below the diagonal. For example, if M=5 and N=3, the matrix V is
 122 *>
 123 *>               V = (  1       )
 124 *>                   ( v1  1    )
 125 *>                   ( v1 v2  1 )
 126 *>                   ( v1 v2 v3 )
 127 *>                   ( v1 v2 v3 )
 128 *>
 129 *>  where the vi's represent the vectors which define H(i), which are returned
 130 *>  in the matrix A.  The 1's along the diagonal of V are not stored in A.
 131 *>
 132 *>  Let K=MIN(M,N).  The number of blocks is B = ceiling(K/NB), where each
 133 *>  block is of order NB except for the last block, which is of order
 134 *>  IB = K - (B-1)*NB.  For each of the B blocks, a upper triangular block
 135 *>  reflector factor is computed: T1, T2, ..., TB.  The NB-by-NB (and IB-by-IB
 136 *>  for the last block) T's are stored in the NB-by-N matrix T as
 137 *>
 138 *>               T = (T1 T2 ... TB).
 139 *> \endverbatim
 140 *>
 141 *  =====================================================================
 142       SUBROUTINE CGEQRT( M, N, NB, A, LDA, T, LDT, WORK, INFO )
 143 *
 144 *  -- LAPACK computational routine (version 3.5.0) --
 145 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 146 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 147 *     November 2013
 148 *
 149 *     .. Scalar Arguments ..
 150       INTEGER INFO, LDA, LDT, M, N, NB
 151 *     ..
 152 *     .. Array Arguments ..
 153       COMPLEX A( LDA, * ), T( LDT, * ), WORK( * )
 154 *     ..
 155 *
 156 * =====================================================================
 157 *
 158 *     ..
 159 *     .. Local Scalars ..
 160       INTEGER    I, IB, IINFO, K
 161       LOGICAL    USE_RECURSIVE_QR
 162       PARAMETER( USE_RECURSIVE_QR=.TRUE. )
 163 *     ..
 164 *     .. External Subroutines ..
 165       EXTERNAL   CGEQRT2, CGEQRT3, CLARFB, XERBLA
 166 *     ..
 167 *     .. Executable Statements ..
 168 *
 169 *     Test the input arguments
 170 *
 171       INFO = 0
 172       IF( M.LT.0 ) THEN
 173          INFO = -1
 174       ELSE IF( N.LT.0 ) THEN
 175          INFO = -2
 176       ELSE IF( NB.LT.1 .OR. ( NB.GT.MIN(M,N) .AND. MIN(M,N).GT.0 ) )THEN
 177          INFO = -3
 178       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
 179          INFO = -5
 180       ELSE IF( LDT.LT.NB ) THEN
 181          INFO = -7
 182       END IF
 183       IF( INFO.NE.0 ) THEN
 184          CALL XERBLA( 'CGEQRT', -INFO )
 185          RETURN
 186       END IF
 187 *
 188 *     Quick return if possible
 189 *
 190       K = MIN( M, N )
 191       IF( K.EQ.0 ) RETURN
 192 *
 193 *     Blocked loop of length K
 194 *
 195       DO I = 1, K,  NB
 196          IB = MIN( K-I+1, NB )
 197 *
 198 *     Compute the QR factorization of the current block A(I:M,I:I+IB-1)
 199 *
 200          IF( USE_RECURSIVE_QR ) THEN
 201             CALL CGEQRT3( M-I+1, IB, A(I,I), LDA, T(1,I), LDT, IINFO )
 202          ELSE
 203             CALL CGEQRT2( M-I+1, IB, A(I,I), LDA, T(1,I), LDT, IINFO )
 204          END IF
 205          IF( I+IB.LE.N ) THEN
 206 *
 207 *     Update by applying H**H to A(I:M,I+IB:N) from the left
 208 *
 209             CALL CLARFB( 'L', 'C', 'F', 'C', M-I+1, N-I-IB+1, IB,
 210      $                   A( I, I ), LDA, T( 1, I ), LDT,
 211      $                   A( I, I+IB ), LDA, WORK , N-I-IB+1 )
 212          END IF
 213       END DO
 214       RETURN
 215 *
 216 *     End of CGEQRT
 217 *
 218       END