![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b ZLARZB 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download ZLARZB + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlarzb.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlarzb.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlarzb.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE ZLARZB( SIDE, TRANS, DIRECT, STOREV, M, N, K, L, V, 00022 * LDV, T, LDT, C, LDC, WORK, LDWORK ) 00023 * 00024 * .. Scalar Arguments .. 00025 * CHARACTER DIRECT, SIDE, STOREV, TRANS 00026 * INTEGER K, L, LDC, LDT, LDV, LDWORK, M, N 00027 * .. 00028 * .. Array Arguments .. 00029 * COMPLEX*16 C( LDC, * ), T( LDT, * ), V( LDV, * ), 00030 * $ WORK( LDWORK, * ) 00031 * .. 00032 * 00033 * 00034 *> \par Purpose: 00035 * ============= 00036 *> 00037 *> \verbatim 00038 *> 00039 *> ZLARZB applies a complex block reflector H or its transpose H**H 00040 *> to a complex distributed M-by-N C from the left or the right. 00041 *> 00042 *> Currently, only STOREV = 'R' and DIRECT = 'B' are supported. 00043 *> \endverbatim 00044 * 00045 * Arguments: 00046 * ========== 00047 * 00048 *> \param[in] SIDE 00049 *> \verbatim 00050 *> SIDE is CHARACTER*1 00051 *> = 'L': apply H or H**H from the Left 00052 *> = 'R': apply H or H**H from the Right 00053 *> \endverbatim 00054 *> 00055 *> \param[in] TRANS 00056 *> \verbatim 00057 *> TRANS is CHARACTER*1 00058 *> = 'N': apply H (No transpose) 00059 *> = 'C': apply H**H (Conjugate transpose) 00060 *> \endverbatim 00061 *> 00062 *> \param[in] DIRECT 00063 *> \verbatim 00064 *> DIRECT is CHARACTER*1 00065 *> Indicates how H is formed from a product of elementary 00066 *> reflectors 00067 *> = 'F': H = H(1) H(2) . . . H(k) (Forward, not supported yet) 00068 *> = 'B': H = H(k) . . . H(2) H(1) (Backward) 00069 *> \endverbatim 00070 *> 00071 *> \param[in] STOREV 00072 *> \verbatim 00073 *> STOREV is CHARACTER*1 00074 *> Indicates how the vectors which define the elementary 00075 *> reflectors are stored: 00076 *> = 'C': Columnwise (not supported yet) 00077 *> = 'R': Rowwise 00078 *> \endverbatim 00079 *> 00080 *> \param[in] M 00081 *> \verbatim 00082 *> M is INTEGER 00083 *> The number of rows of the matrix C. 00084 *> \endverbatim 00085 *> 00086 *> \param[in] N 00087 *> \verbatim 00088 *> N is INTEGER 00089 *> The number of columns of the matrix C. 00090 *> \endverbatim 00091 *> 00092 *> \param[in] K 00093 *> \verbatim 00094 *> K is INTEGER 00095 *> The order of the matrix T (= the number of elementary 00096 *> reflectors whose product defines the block reflector). 00097 *> \endverbatim 00098 *> 00099 *> \param[in] L 00100 *> \verbatim 00101 *> L is INTEGER 00102 *> The number of columns of the matrix V containing the 00103 *> meaningful part of the Householder reflectors. 00104 *> If SIDE = 'L', M >= L >= 0, if SIDE = 'R', N >= L >= 0. 00105 *> \endverbatim 00106 *> 00107 *> \param[in] V 00108 *> \verbatim 00109 *> V is COMPLEX*16 array, dimension (LDV,NV). 00110 *> If STOREV = 'C', NV = K; if STOREV = 'R', NV = L. 00111 *> \endverbatim 00112 *> 00113 *> \param[in] LDV 00114 *> \verbatim 00115 *> LDV is INTEGER 00116 *> The leading dimension of the array V. 00117 *> If STOREV = 'C', LDV >= L; if STOREV = 'R', LDV >= K. 00118 *> \endverbatim 00119 *> 00120 *> \param[in] T 00121 *> \verbatim 00122 *> T is COMPLEX*16 array, dimension (LDT,K) 00123 *> The triangular K-by-K matrix T in the representation of the 00124 *> block reflector. 00125 *> \endverbatim 00126 *> 00127 *> \param[in] LDT 00128 *> \verbatim 00129 *> LDT is INTEGER 00130 *> The leading dimension of the array T. LDT >= K. 00131 *> \endverbatim 00132 *> 00133 *> \param[in,out] C 00134 *> \verbatim 00135 *> C is COMPLEX*16 array, dimension (LDC,N) 00136 *> On entry, the M-by-N matrix C. 00137 *> On exit, C is overwritten by H*C or H**H*C or C*H or C*H**H. 00138 *> \endverbatim 00139 *> 00140 *> \param[in] LDC 00141 *> \verbatim 00142 *> LDC is INTEGER 00143 *> The leading dimension of the array C. LDC >= max(1,M). 00144 *> \endverbatim 00145 *> 00146 *> \param[out] WORK 00147 *> \verbatim 00148 *> WORK is COMPLEX*16 array, dimension (LDWORK,K) 00149 *> \endverbatim 00150 *> 00151 *> \param[in] LDWORK 00152 *> \verbatim 00153 *> LDWORK is INTEGER 00154 *> The leading dimension of the array WORK. 00155 *> If SIDE = 'L', LDWORK >= max(1,N); 00156 *> if SIDE = 'R', LDWORK >= max(1,M). 00157 *> \endverbatim 00158 * 00159 * Authors: 00160 * ======== 00161 * 00162 *> \author Univ. of Tennessee 00163 *> \author Univ. of California Berkeley 00164 *> \author Univ. of Colorado Denver 00165 *> \author NAG Ltd. 00166 * 00167 *> \date November 2011 00168 * 00169 *> \ingroup complex16OTHERcomputational 00170 * 00171 *> \par Contributors: 00172 * ================== 00173 *> 00174 *> A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA 00175 * 00176 *> \par Further Details: 00177 * ===================== 00178 *> 00179 *> \verbatim 00180 *> \endverbatim 00181 *> 00182 * ===================================================================== 00183 SUBROUTINE ZLARZB( SIDE, TRANS, DIRECT, STOREV, M, N, K, L, V, 00184 $ LDV, T, LDT, C, LDC, WORK, LDWORK ) 00185 * 00186 * -- LAPACK computational routine (version 3.4.0) -- 00187 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00188 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00189 * November 2011 00190 * 00191 * .. Scalar Arguments .. 00192 CHARACTER DIRECT, SIDE, STOREV, TRANS 00193 INTEGER K, L, LDC, LDT, LDV, LDWORK, M, N 00194 * .. 00195 * .. Array Arguments .. 00196 COMPLEX*16 C( LDC, * ), T( LDT, * ), V( LDV, * ), 00197 $ WORK( LDWORK, * ) 00198 * .. 00199 * 00200 * ===================================================================== 00201 * 00202 * .. Parameters .. 00203 COMPLEX*16 ONE 00204 PARAMETER ( ONE = ( 1.0D+0, 0.0D+0 ) ) 00205 * .. 00206 * .. Local Scalars .. 00207 CHARACTER TRANST 00208 INTEGER I, INFO, J 00209 * .. 00210 * .. External Functions .. 00211 LOGICAL LSAME 00212 EXTERNAL LSAME 00213 * .. 00214 * .. External Subroutines .. 00215 EXTERNAL XERBLA, ZCOPY, ZGEMM, ZLACGV, ZTRMM 00216 * .. 00217 * .. Executable Statements .. 00218 * 00219 * Quick return if possible 00220 * 00221 IF( M.LE.0 .OR. N.LE.0 ) 00222 $ RETURN 00223 * 00224 * Check for currently supported options 00225 * 00226 INFO = 0 00227 IF( .NOT.LSAME( DIRECT, 'B' ) ) THEN 00228 INFO = -3 00229 ELSE IF( .NOT.LSAME( STOREV, 'R' ) ) THEN 00230 INFO = -4 00231 END IF 00232 IF( INFO.NE.0 ) THEN 00233 CALL XERBLA( 'ZLARZB', -INFO ) 00234 RETURN 00235 END IF 00236 * 00237 IF( LSAME( TRANS, 'N' ) ) THEN 00238 TRANST = 'C' 00239 ELSE 00240 TRANST = 'N' 00241 END IF 00242 * 00243 IF( LSAME( SIDE, 'L' ) ) THEN 00244 * 00245 * Form H * C or H**H * C 00246 * 00247 * W( 1:n, 1:k ) = C( 1:k, 1:n )**H 00248 * 00249 DO 10 J = 1, K 00250 CALL ZCOPY( N, C( J, 1 ), LDC, WORK( 1, J ), 1 ) 00251 10 CONTINUE 00252 * 00253 * W( 1:n, 1:k ) = W( 1:n, 1:k ) + ... 00254 * C( m-l+1:m, 1:n )**H * V( 1:k, 1:l )**T 00255 * 00256 IF( L.GT.0 ) 00257 $ CALL ZGEMM( 'Transpose', 'Conjugate transpose', N, K, L, 00258 $ ONE, C( M-L+1, 1 ), LDC, V, LDV, ONE, WORK, 00259 $ LDWORK ) 00260 * 00261 * W( 1:n, 1:k ) = W( 1:n, 1:k ) * T**T or W( 1:m, 1:k ) * T 00262 * 00263 CALL ZTRMM( 'Right', 'Lower', TRANST, 'Non-unit', N, K, ONE, T, 00264 $ LDT, WORK, LDWORK ) 00265 * 00266 * C( 1:k, 1:n ) = C( 1:k, 1:n ) - W( 1:n, 1:k )**H 00267 * 00268 DO 30 J = 1, N 00269 DO 20 I = 1, K 00270 C( I, J ) = C( I, J ) - WORK( J, I ) 00271 20 CONTINUE 00272 30 CONTINUE 00273 * 00274 * C( m-l+1:m, 1:n ) = C( m-l+1:m, 1:n ) - ... 00275 * V( 1:k, 1:l )**H * W( 1:n, 1:k )**H 00276 * 00277 IF( L.GT.0 ) 00278 $ CALL ZGEMM( 'Transpose', 'Transpose', L, N, K, -ONE, V, LDV, 00279 $ WORK, LDWORK, ONE, C( M-L+1, 1 ), LDC ) 00280 * 00281 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00282 * 00283 * Form C * H or C * H**H 00284 * 00285 * W( 1:m, 1:k ) = C( 1:m, 1:k ) 00286 * 00287 DO 40 J = 1, K 00288 CALL ZCOPY( M, C( 1, J ), 1, WORK( 1, J ), 1 ) 00289 40 CONTINUE 00290 * 00291 * W( 1:m, 1:k ) = W( 1:m, 1:k ) + ... 00292 * C( 1:m, n-l+1:n ) * V( 1:k, 1:l )**H 00293 * 00294 IF( L.GT.0 ) 00295 $ CALL ZGEMM( 'No transpose', 'Transpose', M, K, L, ONE, 00296 $ C( 1, N-L+1 ), LDC, V, LDV, ONE, WORK, LDWORK ) 00297 * 00298 * W( 1:m, 1:k ) = W( 1:m, 1:k ) * conjg( T ) or 00299 * W( 1:m, 1:k ) * T**H 00300 * 00301 DO 50 J = 1, K 00302 CALL ZLACGV( K-J+1, T( J, J ), 1 ) 00303 50 CONTINUE 00304 CALL ZTRMM( 'Right', 'Lower', TRANS, 'Non-unit', M, K, ONE, T, 00305 $ LDT, WORK, LDWORK ) 00306 DO 60 J = 1, K 00307 CALL ZLACGV( K-J+1, T( J, J ), 1 ) 00308 60 CONTINUE 00309 * 00310 * C( 1:m, 1:k ) = C( 1:m, 1:k ) - W( 1:m, 1:k ) 00311 * 00312 DO 80 J = 1, K 00313 DO 70 I = 1, M 00314 C( I, J ) = C( I, J ) - WORK( I, J ) 00315 70 CONTINUE 00316 80 CONTINUE 00317 * 00318 * C( 1:m, n-l+1:n ) = C( 1:m, n-l+1:n ) - ... 00319 * W( 1:m, 1:k ) * conjg( V( 1:k, 1:l ) ) 00320 * 00321 DO 90 J = 1, L 00322 CALL ZLACGV( K, V( 1, J ), 1 ) 00323 90 CONTINUE 00324 IF( L.GT.0 ) 00325 $ CALL ZGEMM( 'No transpose', 'No transpose', M, L, K, -ONE, 00326 $ WORK, LDWORK, V, LDV, ONE, C( 1, N-L+1 ), LDC ) 00327 DO 100 J = 1, L 00328 CALL ZLACGV( K, V( 1, J ), 1 ) 00329 100 CONTINUE 00330 * 00331 END IF 00332 * 00333 RETURN 00334 * 00335 * End of ZLARZB 00336 * 00337 END