LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
ctzrzf.f
Go to the documentation of this file.
00001 *> \brief \b CTZRZF
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *> \htmlonly
00009 *> Download CTZRZF + dependencies 
00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ctzrzf.f"> 
00011 *> [TGZ]</a> 
00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ctzrzf.f"> 
00013 *> [ZIP]</a> 
00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ctzrzf.f"> 
00015 *> [TXT]</a>
00016 *> \endhtmlonly 
00017 *
00018 *  Definition:
00019 *  ===========
00020 *
00021 *       SUBROUTINE CTZRZF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
00022 * 
00023 *       .. Scalar Arguments ..
00024 *       INTEGER            INFO, LDA, LWORK, M, N
00025 *       ..
00026 *       .. Array Arguments ..
00027 *       COMPLEX            A( LDA, * ), TAU( * ), WORK( * )
00028 *       ..
00029 *  
00030 *
00031 *> \par Purpose:
00032 *  =============
00033 *>
00034 *> \verbatim
00035 *>
00036 *> CTZRZF reduces the M-by-N ( M<=N ) complex upper trapezoidal matrix A
00037 *> to upper triangular form by means of unitary transformations.
00038 *>
00039 *> The upper trapezoidal matrix A is factored as
00040 *>
00041 *>    A = ( R  0 ) * Z,
00042 *>
00043 *> where Z is an N-by-N unitary matrix and R is an M-by-M upper
00044 *> triangular matrix.
00045 *> \endverbatim
00046 *
00047 *  Arguments:
00048 *  ==========
00049 *
00050 *> \param[in] M
00051 *> \verbatim
00052 *>          M is INTEGER
00053 *>          The number of rows of the matrix A.  M >= 0.
00054 *> \endverbatim
00055 *>
00056 *> \param[in] N
00057 *> \verbatim
00058 *>          N is INTEGER
00059 *>          The number of columns of the matrix A.  N >= M.
00060 *> \endverbatim
00061 *>
00062 *> \param[in,out] A
00063 *> \verbatim
00064 *>          A is COMPLEX array, dimension (LDA,N)
00065 *>          On entry, the leading M-by-N upper trapezoidal part of the
00066 *>          array A must contain the matrix to be factorized.
00067 *>          On exit, the leading M-by-M upper triangular part of A
00068 *>          contains the upper triangular matrix R, and elements M+1 to
00069 *>          N of the first M rows of A, with the array TAU, represent the
00070 *>          unitary matrix Z as a product of M elementary reflectors.
00071 *> \endverbatim
00072 *>
00073 *> \param[in] LDA
00074 *> \verbatim
00075 *>          LDA is INTEGER
00076 *>          The leading dimension of the array A.  LDA >= max(1,M).
00077 *> \endverbatim
00078 *>
00079 *> \param[out] TAU
00080 *> \verbatim
00081 *>          TAU is COMPLEX array, dimension (M)
00082 *>          The scalar factors of the elementary reflectors.
00083 *> \endverbatim
00084 *>
00085 *> \param[out] WORK
00086 *> \verbatim
00087 *>          WORK is COMPLEX array, dimension (MAX(1,LWORK))
00088 *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
00089 *> \endverbatim
00090 *>
00091 *> \param[in] LWORK
00092 *> \verbatim
00093 *>          LWORK is INTEGER
00094 *>          The dimension of the array WORK.  LWORK >= max(1,M).
00095 *>          For optimum performance LWORK >= M*NB, where NB is
00096 *>          the optimal blocksize.
00097 *>
00098 *>          If LWORK = -1, then a workspace query is assumed; the routine
00099 *>          only calculates the optimal size of the WORK array, returns
00100 *>          this value as the first entry of the WORK array, and no error
00101 *>          message related to LWORK is issued by XERBLA.
00102 *> \endverbatim
00103 *>
00104 *> \param[out] INFO
00105 *> \verbatim
00106 *>          INFO is INTEGER
00107 *>          = 0:  successful exit
00108 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
00109 *> \endverbatim
00110 *
00111 *  Authors:
00112 *  ========
00113 *
00114 *> \author Univ. of Tennessee 
00115 *> \author Univ. of California Berkeley 
00116 *> \author Univ. of Colorado Denver 
00117 *> \author NAG Ltd. 
00118 *
00119 *> \date April 2012
00120 *
00121 *> \ingroup complexOTHERcomputational
00122 *
00123 *> \par Contributors:
00124 *  ==================
00125 *>
00126 *>    A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA
00127 *
00128 *> \par Further Details:
00129 *  =====================
00130 *>
00131 *> \verbatim
00132 *>
00133 *>  The N-by-N matrix Z can be computed by
00134 *>
00135 *>     Z =  Z(1)*Z(2)* ... *Z(M)
00136 *>
00137 *>  where each N-by-N Z(k) is given by
00138 *>
00139 *>     Z(k) = I - tau(k)*v(k)*v(k)**H
00140 *>
00141 *>  with v(k) is the kth row vector of the M-by-N matrix
00142 *>
00143 *>     V = ( I   A(:,M+1:N) )
00144 *>
00145 *>  I is the M-by-M identity matrix, A(:,M+1:N) 
00146 *>  is the output stored in A on exit from DTZRZF,
00147 *>  and tau(k) is the kth element of the array TAU.
00148 *>
00149 *> \endverbatim
00150 *>
00151 *  =====================================================================
00152       SUBROUTINE CTZRZF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
00153 *
00154 *  -- LAPACK computational routine (version 3.4.1) --
00155 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00156 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00157 *     April 2012
00158 *
00159 *     .. Scalar Arguments ..
00160       INTEGER            INFO, LDA, LWORK, M, N
00161 *     ..
00162 *     .. Array Arguments ..
00163       COMPLEX            A( LDA, * ), TAU( * ), WORK( * )
00164 *     ..
00165 *
00166 *  =====================================================================
00167 *
00168 *     .. Parameters ..
00169       COMPLEX            ZERO
00170       PARAMETER          ( ZERO = ( 0.0E+0, 0.0E+0 ) )
00171 *     ..
00172 *     .. Local Scalars ..
00173       LOGICAL            LQUERY
00174       INTEGER            I, IB, IWS, KI, KK, LDWORK, LWKMIN, LWKOPT,
00175      $                   M1, MU, NB, NBMIN, NX
00176 *     ..
00177 *     .. External Subroutines ..
00178       EXTERNAL           XERBLA, CLARZB, CLARZT, CLATRZ
00179 *     ..
00180 *     .. Intrinsic Functions ..
00181       INTRINSIC          MAX, MIN
00182 *     ..
00183 *     .. External Functions ..
00184       INTEGER            ILAENV
00185       EXTERNAL           ILAENV
00186 *     ..
00187 *     .. Executable Statements ..
00188 *
00189 *     Test the input arguments
00190 *
00191       INFO = 0
00192       LQUERY = ( LWORK.EQ.-1 )
00193       IF( M.LT.0 ) THEN
00194          INFO = -1
00195       ELSE IF( N.LT.M ) THEN
00196          INFO = -2
00197       ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
00198          INFO = -4
00199       END IF
00200 *
00201       IF( INFO.EQ.0 ) THEN
00202          IF( M.EQ.0 .OR. M.EQ.N ) THEN
00203             LWKOPT = 1
00204             LWKMIN = 1
00205          ELSE
00206 *
00207 *           Determine the block size.
00208 *
00209             NB = ILAENV( 1, 'CGERQF', ' ', M, N, -1, -1 )
00210             LWKOPT = M*NB
00211             LWKMIN = MAX( 1, M )
00212          END IF
00213          WORK( 1 ) = LWKOPT
00214 *
00215          IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) THEN
00216             INFO = -7
00217          END IF
00218       END IF
00219 *
00220       IF( INFO.NE.0 ) THEN
00221          CALL XERBLA( 'CTZRZF', -INFO )
00222          RETURN
00223       ELSE IF( LQUERY ) THEN
00224          RETURN
00225       END IF
00226 *
00227 *     Quick return if possible
00228 *
00229       IF( M.EQ.0 ) THEN
00230          RETURN
00231       ELSE IF( M.EQ.N ) THEN
00232          DO 10 I = 1, N
00233             TAU( I ) = ZERO
00234    10    CONTINUE
00235          RETURN
00236       END IF
00237 *
00238       NBMIN = 2
00239       NX = 1
00240       IWS = M
00241       IF( NB.GT.1 .AND. NB.LT.M ) THEN
00242 *
00243 *        Determine when to cross over from blocked to unblocked code.
00244 *
00245          NX = MAX( 0, ILAENV( 3, 'CGERQF', ' ', M, N, -1, -1 ) )
00246          IF( NX.LT.M ) THEN
00247 *
00248 *           Determine if workspace is large enough for blocked code.
00249 *
00250             LDWORK = M
00251             IWS = LDWORK*NB
00252             IF( LWORK.LT.IWS ) THEN
00253 *
00254 *              Not enough workspace to use optimal NB:  reduce NB and
00255 *              determine the minimum value of NB.
00256 *
00257                NB = LWORK / LDWORK
00258                NBMIN = MAX( 2, ILAENV( 2, 'CGERQF', ' ', M, N, -1,
00259      $                 -1 ) )
00260             END IF
00261          END IF
00262       END IF
00263 *
00264       IF( NB.GE.NBMIN .AND. NB.LT.M .AND. NX.LT.M ) THEN
00265 *
00266 *        Use blocked code initially.
00267 *        The last kk rows are handled by the block method.
00268 *
00269          M1 = MIN( M+1, N )
00270          KI = ( ( M-NX-1 ) / NB )*NB
00271          KK = MIN( M, KI+NB )
00272 *
00273          DO 20 I = M - KK + KI + 1, M - KK + 1, -NB
00274             IB = MIN( M-I+1, NB )
00275 *
00276 *           Compute the TZ factorization of the current block
00277 *           A(i:i+ib-1,i:n)
00278 *
00279             CALL CLATRZ( IB, N-I+1, N-M, A( I, I ), LDA, TAU( I ),
00280      $                   WORK )
00281             IF( I.GT.1 ) THEN
00282 *
00283 *              Form the triangular factor of the block reflector
00284 *              H = H(i+ib-1) . . . H(i+1) H(i)
00285 *
00286                CALL CLARZT( 'Backward', 'Rowwise', N-M, IB, A( I, M1 ),
00287      $                      LDA, TAU( I ), WORK, LDWORK )
00288 *
00289 *              Apply H to A(1:i-1,i:n) from the right
00290 *
00291                CALL CLARZB( 'Right', 'No transpose', 'Backward',
00292      $                      'Rowwise', I-1, N-I+1, IB, N-M, A( I, M1 ),
00293      $                      LDA, WORK, LDWORK, A( 1, I ), LDA,
00294      $                      WORK( IB+1 ), LDWORK )
00295             END IF
00296    20    CONTINUE
00297          MU = I + NB - 1
00298       ELSE
00299          MU = M
00300       END IF
00301 *
00302 *     Use unblocked code to factor the last or only block
00303 *
00304       IF( MU.GT.0 )
00305      $   CALL CLATRZ( MU, N, N-M, A, LDA, TAU, WORK )
00306 *
00307       WORK( 1 ) = LWKOPT
00308 *
00309       RETURN
00310 *
00311 *     End of CTZRZF
00312 *
00313       END
 All Files Functions