SRC/zgelqf.f

   1 *> \brief \b ZGELQF
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download ZGELQF + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgelqf.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgelqf.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgelqf.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE ZGELQF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       INTEGER            INFO, LDA, LWORK, M, N
  25 *       ..
  26 *       .. Array Arguments ..
  27 *       COMPLEX*16         A( LDA, * ), TAU( * ), WORK( * )
  28 *       ..
  29 *
  30 *
  31 *> \par Purpose:
  32 *  =============
  33 *>
  34 *> \verbatim
  35 *>
  36 *> ZGELQF computes an LQ factorization of a complex M-by-N matrix A:
  37 *> A = L * 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,out] A
  56 *> \verbatim
  57 *>          A is COMPLEX*16 array, dimension (LDA,N)
  58 *>          On entry, the M-by-N matrix A.
  59 *>          On exit, the elements on and below the diagonal of the array
  60 *>          contain the m-by-min(m,n) lower trapezoidal matrix L (L is
  61 *>          lower triangular if m <= n); the elements above the diagonal,
  62 *>          with the array TAU, represent the unitary matrix Q as a
  63 *>          product of elementary reflectors (see Further Details).
  64 *> \endverbatim
  65 *>
  66 *> \param[in] LDA
  67 *> \verbatim
  68 *>          LDA is INTEGER
  69 *>          The leading dimension of the array A.  LDA >= max(1,M).
  70 *> \endverbatim
  71 *>
  72 *> \param[out] TAU
  73 *> \verbatim
  74 *>          TAU is COMPLEX*16 array, dimension (min(M,N))
  75 *>          The scalar factors of the elementary reflectors (see Further
  76 *>          Details).
  77 *> \endverbatim
  78 *>
  79 *> \param[out] WORK
  80 *> \verbatim
  81 *>          WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
  82 *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
  83 *> \endverbatim
  84 *>
  85 *> \param[in] LWORK
  86 *> \verbatim
  87 *>          LWORK is INTEGER
  88 *>          The dimension of the array WORK.  LWORK >= max(1,M).
  89 *>          For optimum performance LWORK >= M*NB, where NB is the
  90 *>          optimal blocksize.
  91 *>
  92 *>          If LWORK = -1, then a workspace query is assumed; the routine
  93 *>          only calculates the optimal size of the WORK array, returns
  94 *>          this value as the first entry of the WORK array, and no error
  95 *>          message related to LWORK is issued by XERBLA.
  96 *> \endverbatim
  97 *>
  98 *> \param[out] INFO
  99 *> \verbatim
 100 *>          INFO is INTEGER
 101 *>          = 0:  successful exit
 102 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
 103 *> \endverbatim
 104 *
 105 *  Authors:
 106 *  ========
 107 *
 108 *> \author Univ. of Tennessee
 109 *> \author Univ. of California Berkeley
 110 *> \author Univ. of Colorado Denver
 111 *> \author NAG Ltd.
 112 *
 113 *> \date November 2011
 114 *
 115 *> \ingroup complex16GEcomputational
 116 *
 117 *> \par Further Details:
 118 *  =====================
 119 *>
 120 *> \verbatim
 121 *>
 122 *>  The matrix Q is represented as a product of elementary reflectors
 123 *>
 124 *>     Q = H(k)**H . . . H(2)**H H(1)**H, where k = min(m,n).
 125 *>
 126 *>  Each H(i) has the form
 127 *>
 128 *>     H(i) = I - tau * v * v**H
 129 *>
 130 *>  where tau is a complex scalar, and v is a complex vector with
 131 *>  v(1:i-1) = 0 and v(i) = 1; conjg(v(i+1:n)) is stored on exit in
 132 *>  A(i,i+1:n), and tau in TAU(i).
 133 *> \endverbatim
 134 *>
 135 *  =====================================================================
 136       SUBROUTINE ZGELQF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
 137 *
 138 *  -- LAPACK computational routine (version 3.4.0) --
 139 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 140 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 141 *     November 2011
 142 *
 143 *     .. Scalar Arguments ..
 144       INTEGER            INFO, LDA, LWORK, M, N
 145 *     ..
 146 *     .. Array Arguments ..
 147       COMPLEX*16         A( LDA, * ), TAU( * ), WORK( * )
 148 *     ..
 149 *
 150 *  =====================================================================
 151 *
 152 *     .. Local Scalars ..
 153       LOGICAL            LQUERY
 154       INTEGER            I, IB, IINFO, IWS, K, LDWORK, LWKOPT, NB,
 155      $                   NBMIN, NX
 156 *     ..
 157 *     .. External Subroutines ..
 158       EXTERNAL           XERBLA, ZGELQ2, ZLARFB, ZLARFT
 159 *     ..
 160 *     .. Intrinsic Functions ..
 161       INTRINSIC          MAX, MIN
 162 *     ..
 163 *     .. External Functions ..
 164       INTEGER            ILAENV
 165       EXTERNAL           ILAENV
 166 *     ..
 167 *     .. Executable Statements ..
 168 *
 169 *     Test the input arguments
 170 *
 171       INFO = 0
 172       NB = ILAENV( 1, 'ZGELQF', ' ', M, N, -1, -1 )
 173       LWKOPT = M*NB
 174       WORK( 1 ) = LWKOPT
 175       LQUERY = ( LWORK.EQ.-1 )
 176       IF( M.LT.0 ) THEN
 177          INFO = -1
 178       ELSE IF( N.LT.0 ) THEN
 179          INFO = -2
 180       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
 181          INFO = -4
 182       ELSE IF( LWORK.LT.MAX( 1, M ) .AND. .NOT.LQUERY ) THEN
 183          INFO = -7
 184       END IF
 185       IF( INFO.NE.0 ) THEN
 186          CALL XERBLA( 'ZGELQF', -INFO )
 187          RETURN
 188       ELSE IF( LQUERY ) THEN
 189          RETURN
 190       END IF
 191 *
 192 *     Quick return if possible
 193 *
 194       K = MIN( M, N )
 195       IF( K.EQ.0 ) THEN
 196          WORK( 1 ) = 1
 197          RETURN
 198       END IF
 199 *
 200       NBMIN = 2
 201       NX = 0
 202       IWS = M
 203       IF( NB.GT.1 .AND. NB.LT.K ) THEN
 204 *
 205 *        Determine when to cross over from blocked to unblocked code.
 206 *
 207          NX = MAX( 0, ILAENV( 3, 'ZGELQF', ' ', M, N, -1, -1 ) )
 208          IF( NX.LT.K ) THEN
 209 *
 210 *           Determine if workspace is large enough for blocked code.
 211 *
 212             LDWORK = M
 213             IWS = LDWORK*NB
 214             IF( LWORK.LT.IWS ) THEN
 215 *
 216 *              Not enough workspace to use optimal NB:  reduce NB and
 217 *              determine the minimum value of NB.
 218 *
 219                NB = LWORK / LDWORK
 220                NBMIN = MAX( 2, ILAENV( 2, 'ZGELQF', ' ', M, N, -1,
 221      $                 -1 ) )
 222             END IF
 223          END IF
 224       END IF
 225 *
 226       IF( NB.GE.NBMIN .AND. NB.LT.K .AND. NX.LT.K ) THEN
 227 *
 228 *        Use blocked code initially
 229 *
 230          DO 10 I = 1, K - NX, NB
 231             IB = MIN( K-I+1, NB )
 232 *
 233 *           Compute the LQ factorization of the current block
 234 *           A(i:i+ib-1,i:n)
 235 *
 236             CALL ZGELQ2( IB, N-I+1, A( I, I ), LDA, TAU( I ), WORK,
 237      $                   IINFO )
 238             IF( I+IB.LE.M ) THEN
 239 *
 240 *              Form the triangular factor of the block reflector
 241 *              H = H(i) H(i+1) . . . H(i+ib-1)
 242 *
 243                CALL ZLARFT( 'Forward', 'Rowwise', N-I+1, IB, A( I, I ),
 244      $                      LDA, TAU( I ), WORK, LDWORK )
 245 *
 246 *              Apply H to A(i+ib:m,i:n) from the right
 247 *
 248                CALL ZLARFB( 'Right', 'No transpose', 'Forward',
 249      $                      'Rowwise', M-I-IB+1, N-I+1, IB, A( I, I ),
 250      $                      LDA, WORK, LDWORK, A( I+IB, I ), LDA,
 251      $                      WORK( IB+1 ), LDWORK )
 252             END IF
 253    10    CONTINUE
 254       ELSE
 255          I = 1
 256       END IF
 257 *
 258 *     Use unblocked code to factor the last or only block.
 259 *
 260       IF( I.LE.K )
 261      $   CALL ZGELQ2( M-I+1, N-I+1, A( I, I ), LDA, TAU( I ), WORK,
 262      $                IINFO )
 263 *
 264       WORK( 1 ) = IWS
 265       RETURN
 266 *
 267 *     End of ZGELQF
 268 *
 269       END