 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE ZLARFB( SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV, 00002 $ T, LDT, C, LDC, WORK, LDWORK ) 00003 IMPLICIT NONE 00004 * 00005 * -- LAPACK auxiliary routine (version 3.3.1) -- 00006 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00007 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00008 * -- April 2011 -- 00009 * 00010 * .. Scalar Arguments .. 00011 CHARACTER DIRECT, SIDE, STOREV, TRANS 00012 INTEGER K, LDC, LDT, LDV, LDWORK, M, N 00013 * .. 00014 * .. Array Arguments .. 00015 COMPLEX*16 C( LDC, * ), T( LDT, * ), V( LDV, * ), 00016 $ WORK( LDWORK, * ) 00017 * .. 00018 * 00019 * Purpose 00020 * ======= 00021 * 00022 * ZLARFB applies a complex block reflector H or its transpose H**H to a 00023 * complex M-by-N matrix C, from either the left or the right. 00024 * 00025 * Arguments 00026 * ========= 00027 * 00028 * SIDE (input) CHARACTER*1 00029 * = 'L': apply H or H**H from the Left 00030 * = 'R': apply H or H**H from the Right 00031 * 00032 * TRANS (input) CHARACTER*1 00033 * = 'N': apply H (No transpose) 00034 * = 'C': apply H**H (Conjugate transpose) 00035 * 00036 * DIRECT (input) CHARACTER*1 00037 * Indicates how H is formed from a product of elementary 00038 * reflectors 00039 * = 'F': H = H(1) H(2) . . . H(k) (Forward) 00040 * = 'B': H = H(k) . . . H(2) H(1) (Backward) 00041 * 00042 * STOREV (input) CHARACTER*1 00043 * Indicates how the vectors which define the elementary 00044 * reflectors are stored: 00045 * = 'C': Columnwise 00046 * = 'R': Rowwise 00047 * 00048 * M (input) INTEGER 00049 * The number of rows of the matrix C. 00050 * 00051 * N (input) INTEGER 00052 * The number of columns of the matrix C. 00053 * 00054 * K (input) INTEGER 00055 * The order of the matrix T (= the number of elementary 00056 * reflectors whose product defines the block reflector). 00057 * 00058 * V (input) COMPLEX*16 array, dimension 00059 * (LDV,K) if STOREV = 'C' 00060 * (LDV,M) if STOREV = 'R' and SIDE = 'L' 00061 * (LDV,N) if STOREV = 'R' and SIDE = 'R' 00062 * The matrix V. See Further Details. 00063 * 00064 * LDV (input) INTEGER 00065 * The leading dimension of the array V. 00066 * If STOREV = 'C' and SIDE = 'L', LDV >= max(1,M); 00067 * if STOREV = 'C' and SIDE = 'R', LDV >= max(1,N); 00068 * if STOREV = 'R', LDV >= K. 00069 * 00070 * T (input) COMPLEX*16 array, dimension (LDT,K) 00071 * The triangular K-by-K matrix T in the representation of the 00072 * block reflector. 00073 * 00074 * LDT (input) INTEGER 00075 * The leading dimension of the array T. LDT >= K. 00076 * 00077 * C (input/output) COMPLEX*16 array, dimension (LDC,N) 00078 * On entry, the M-by-N matrix C. 00079 * On exit, C is overwritten by H*C or H**H*C or C*H or C*H**H. 00080 * 00081 * LDC (input) INTEGER 00082 * The leading dimension of the array C. LDC >= max(1,M). 00083 * 00084 * WORK (workspace) COMPLEX*16 array, dimension (LDWORK,K) 00085 * 00086 * LDWORK (input) INTEGER 00087 * The leading dimension of the array WORK. 00088 * If SIDE = 'L', LDWORK >= max(1,N); 00089 * if SIDE = 'R', LDWORK >= max(1,M). 00090 * 00091 * Further Details 00092 * =============== 00093 * 00094 * The shape of the matrix V and the storage of the vectors which define 00095 * the H(i) is best illustrated by the following example with n = 5 and 00096 * k = 3. The elements equal to 1 are not stored; the corresponding 00097 * array elements are modified but restored on exit. The rest of the 00098 * array is not used. 00099 * 00100 * DIRECT = 'F' and STOREV = 'C': DIRECT = 'F' and STOREV = 'R': 00101 * 00102 * V = ( 1 ) V = ( 1 v1 v1 v1 v1 ) 00103 * ( v1 1 ) ( 1 v2 v2 v2 ) 00104 * ( v1 v2 1 ) ( 1 v3 v3 ) 00105 * ( v1 v2 v3 ) 00106 * ( v1 v2 v3 ) 00107 * 00108 * DIRECT = 'B' and STOREV = 'C': DIRECT = 'B' and STOREV = 'R': 00109 * 00110 * V = ( v1 v2 v3 ) V = ( v1 v1 1 ) 00111 * ( v1 v2 v3 ) ( v2 v2 v2 1 ) 00112 * ( 1 v2 v3 ) ( v3 v3 v3 v3 1 ) 00113 * ( 1 v3 ) 00114 * ( 1 ) 00115 * 00116 * ===================================================================== 00117 * 00118 * .. Parameters .. 00119 COMPLEX*16 ONE 00120 PARAMETER ( ONE = ( 1.0D+0, 0.0D+0 ) ) 00121 * .. 00122 * .. Local Scalars .. 00123 CHARACTER TRANST 00124 INTEGER I, J, LASTV, LASTC 00125 * .. 00126 * .. External Functions .. 00127 LOGICAL LSAME 00128 INTEGER ILAZLR, ILAZLC 00129 EXTERNAL LSAME, ILAZLR, ILAZLC 00130 * .. 00131 * .. External Subroutines .. 00132 EXTERNAL ZCOPY, ZGEMM, ZLACGV, ZTRMM 00133 * .. 00134 * .. Intrinsic Functions .. 00135 INTRINSIC DCONJG 00136 * .. 00137 * .. Executable Statements .. 00138 * 00139 * Quick return if possible 00140 * 00141 IF( M.LE.0 .OR. N.LE.0 ) 00142 $ RETURN 00143 * 00144 IF( LSAME( TRANS, 'N' ) ) THEN 00145 TRANST = 'C' 00146 ELSE 00147 TRANST = 'N' 00148 END IF 00149 * 00150 IF( LSAME( STOREV, 'C' ) ) THEN 00151 * 00152 IF( LSAME( DIRECT, 'F' ) ) THEN 00153 * 00154 * Let V = ( V1 ) (first K rows) 00155 * ( V2 ) 00156 * where V1 is unit lower triangular. 00157 * 00158 IF( LSAME( SIDE, 'L' ) ) THEN 00159 * 00160 * Form H * C or H**H * C where C = ( C1 ) 00161 * ( C2 ) 00162 * 00163 LASTV = MAX( K, ILAZLR( M, K, V, LDV ) ) 00164 LASTC = ILAZLC( LASTV, N, C, LDC ) 00165 * 00166 * W := C**H * V = (C1**H * V1 + C2**H * V2) (stored in WORK) 00167 * 00168 * W := C1**H 00169 * 00170 DO 10 J = 1, K 00171 CALL ZCOPY( LASTC, C( J, 1 ), LDC, WORK( 1, J ), 1 ) 00172 CALL ZLACGV( LASTC, WORK( 1, J ), 1 ) 00173 10 CONTINUE 00174 * 00175 * W := W * V1 00176 * 00177 CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', 00178 $ LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00179 IF( LASTV.GT.K ) THEN 00180 * 00181 * W := W + C2**H *V2 00182 * 00183 CALL ZGEMM( 'Conjugate transpose', 'No transpose', 00184 $ LASTC, K, LASTV-K, ONE, C( K+1, 1 ), LDC, 00185 $ V( K+1, 1 ), LDV, ONE, WORK, LDWORK ) 00186 END IF 00187 * 00188 * W := W * T**H or W * T 00189 * 00190 CALL ZTRMM( 'Right', 'Upper', TRANST, 'Non-unit', 00191 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00192 * 00193 * C := C - V * W**H 00194 * 00195 IF( M.GT.K ) THEN 00196 * 00197 * C2 := C2 - V2 * W**H 00198 * 00199 CALL ZGEMM( 'No transpose', 'Conjugate transpose', 00200 $ LASTV-K, LASTC, K, 00201 $ -ONE, V( K+1, 1 ), LDV, WORK, LDWORK, 00202 $ ONE, C( K+1, 1 ), LDC ) 00203 END IF 00204 * 00205 * W := W * V1**H 00206 * 00207 CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose', 00208 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00209 * 00210 * C1 := C1 - W**H 00211 * 00212 DO 30 J = 1, K 00213 DO 20 I = 1, LASTC 00214 C( J, I ) = C( J, I ) - DCONJG( WORK( I, J ) ) 00215 20 CONTINUE 00216 30 CONTINUE 00217 * 00218 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00219 * 00220 * Form C * H or C * H**H where C = ( C1 C2 ) 00221 * 00222 LASTV = MAX( K, ILAZLR( N, K, V, LDV ) ) 00223 LASTC = ILAZLR( M, LASTV, C, LDC ) 00224 * 00225 * W := C * V = (C1*V1 + C2*V2) (stored in WORK) 00226 * 00227 * W := C1 00228 * 00229 DO 40 J = 1, K 00230 CALL ZCOPY( LASTC, C( 1, J ), 1, WORK( 1, J ), 1 ) 00231 40 CONTINUE 00232 * 00233 * W := W * V1 00234 * 00235 CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', 00236 $ LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00237 IF( LASTV.GT.K ) THEN 00238 * 00239 * W := W + C2 * V2 00240 * 00241 CALL ZGEMM( 'No transpose', 'No transpose', 00242 $ LASTC, K, LASTV-K, 00243 $ ONE, C( 1, K+1 ), LDC, V( K+1, 1 ), LDV, 00244 $ ONE, WORK, LDWORK ) 00245 END IF 00246 * 00247 * W := W * T or W * T**H 00248 * 00249 CALL ZTRMM( 'Right', 'Upper', TRANS, 'Non-unit', 00250 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00251 * 00252 * C := C - W * V**H 00253 * 00254 IF( LASTV.GT.K ) THEN 00255 * 00256 * C2 := C2 - W * V2**H 00257 * 00258 CALL ZGEMM( 'No transpose', 'Conjugate transpose', 00259 $ LASTC, LASTV-K, K, 00260 $ -ONE, WORK, LDWORK, V( K+1, 1 ), LDV, 00261 $ ONE, C( 1, K+1 ), LDC ) 00262 END IF 00263 * 00264 * W := W * V1**H 00265 * 00266 CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose', 00267 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00268 * 00269 * C1 := C1 - W 00270 * 00271 DO 60 J = 1, K 00272 DO 50 I = 1, LASTC 00273 C( I, J ) = C( I, J ) - WORK( I, J ) 00274 50 CONTINUE 00275 60 CONTINUE 00276 END IF 00277 * 00278 ELSE 00279 * 00280 * Let V = ( V1 ) 00281 * ( V2 ) (last K rows) 00282 * where V2 is unit upper triangular. 00283 * 00284 IF( LSAME( SIDE, 'L' ) ) THEN 00285 * 00286 * Form H * C or H**H * C where C = ( C1 ) 00287 * ( C2 ) 00288 * 00289 LASTV = MAX( K, ILAZLR( M, K, V, LDV ) ) 00290 LASTC = ILAZLC( LASTV, N, C, LDC ) 00291 * 00292 * W := C**H * V = (C1**H * V1 + C2**H * V2) (stored in WORK) 00293 * 00294 * W := C2**H 00295 * 00296 DO 70 J = 1, K 00297 CALL ZCOPY( LASTC, C( LASTV-K+J, 1 ), LDC, 00298 $ WORK( 1, J ), 1 ) 00299 CALL ZLACGV( LASTC, WORK( 1, J ), 1 ) 00300 70 CONTINUE 00301 * 00302 * W := W * V2 00303 * 00304 CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00305 $ LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00306 $ WORK, LDWORK ) 00307 IF( LASTV.GT.K ) THEN 00308 * 00309 * W := W + C1**H*V1 00310 * 00311 CALL ZGEMM( 'Conjugate transpose', 'No transpose', 00312 $ LASTC, K, LASTV-K, 00313 $ ONE, C, LDC, V, LDV, 00314 $ ONE, WORK, LDWORK ) 00315 END IF 00316 * 00317 * W := W * T**H or W * T 00318 * 00319 CALL ZTRMM( 'Right', 'Lower', TRANST, 'Non-unit', 00320 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00321 * 00322 * C := C - V * W**H 00323 * 00324 IF( LASTV.GT.K ) THEN 00325 * 00326 * C1 := C1 - V1 * W**H 00327 * 00328 CALL ZGEMM( 'No transpose', 'Conjugate transpose', 00329 $ LASTV-K, LASTC, K, 00330 $ -ONE, V, LDV, WORK, LDWORK, 00331 $ ONE, C, LDC ) 00332 END IF 00333 * 00334 * W := W * V2**H 00335 * 00336 CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose', 00337 $ 'Unit', LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00338 $ WORK, LDWORK ) 00339 * 00340 * C2 := C2 - W**H 00341 * 00342 DO 90 J = 1, K 00343 DO 80 I = 1, LASTC 00344 C( LASTV-K+J, I ) = C( LASTV-K+J, I ) - 00345 $ DCONJG( WORK( I, J ) ) 00346 80 CONTINUE 00347 90 CONTINUE 00348 * 00349 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00350 * 00351 * Form C * H or C * H**H where C = ( C1 C2 ) 00352 * 00353 LASTV = MAX( K, ILAZLR( N, K, V, LDV ) ) 00354 LASTC = ILAZLR( M, LASTV, C, LDC ) 00355 * 00356 * W := C * V = (C1*V1 + C2*V2) (stored in WORK) 00357 * 00358 * W := C2 00359 * 00360 DO 100 J = 1, K 00361 CALL ZCOPY( LASTC, C( 1, LASTV-K+J ), 1, 00362 $ WORK( 1, J ), 1 ) 00363 100 CONTINUE 00364 * 00365 * W := W * V2 00366 * 00367 CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00368 $ LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00369 $ WORK, LDWORK ) 00370 IF( LASTV.GT.K ) THEN 00371 * 00372 * W := W + C1 * V1 00373 * 00374 CALL ZGEMM( 'No transpose', 'No transpose', 00375 $ LASTC, K, LASTV-K, 00376 $ ONE, C, LDC, V, LDV, ONE, WORK, LDWORK ) 00377 END IF 00378 * 00379 * W := W * T or W * T**H 00380 * 00381 CALL ZTRMM( 'Right', 'Lower', TRANS, 'Non-unit', 00382 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00383 * 00384 * C := C - W * V**H 00385 * 00386 IF( LASTV.GT.K ) THEN 00387 * 00388 * C1 := C1 - W * V1**H 00389 * 00390 CALL ZGEMM( 'No transpose', 'Conjugate transpose', 00391 $ LASTC, LASTV-K, K, -ONE, WORK, LDWORK, V, LDV, 00392 $ ONE, C, LDC ) 00393 END IF 00394 * 00395 * W := W * V2**H 00396 * 00397 CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose', 00398 $ 'Unit', LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00399 $ WORK, LDWORK ) 00400 * 00401 * C2 := C2 - W 00402 * 00403 DO 120 J = 1, K 00404 DO 110 I = 1, LASTC 00405 C( I, LASTV-K+J ) = C( I, LASTV-K+J ) 00406 $ - WORK( I, J ) 00407 110 CONTINUE 00408 120 CONTINUE 00409 END IF 00410 END IF 00411 * 00412 ELSE IF( LSAME( STOREV, 'R' ) ) THEN 00413 * 00414 IF( LSAME( DIRECT, 'F' ) ) THEN 00415 * 00416 * Let V = ( V1 V2 ) (V1: first K columns) 00417 * where V1 is unit upper triangular. 00418 * 00419 IF( LSAME( SIDE, 'L' ) ) THEN 00420 * 00421 * Form H * C or H**H * C where C = ( C1 ) 00422 * ( C2 ) 00423 * 00424 LASTV = MAX( K, ILAZLC( K, M, V, LDV ) ) 00425 LASTC = ILAZLC( LASTV, N, C, LDC ) 00426 * 00427 * W := C**H * V**H = (C1**H * V1**H + C2**H * V2**H) (stored in WORK) 00428 * 00429 * W := C1**H 00430 * 00431 DO 130 J = 1, K 00432 CALL ZCOPY( LASTC, C( J, 1 ), LDC, WORK( 1, J ), 1 ) 00433 CALL ZLACGV( LASTC, WORK( 1, J ), 1 ) 00434 130 CONTINUE 00435 * 00436 * W := W * V1**H 00437 * 00438 CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose', 00439 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00440 IF( LASTV.GT.K ) THEN 00441 * 00442 * W := W + C2**H*V2**H 00443 * 00444 CALL ZGEMM( 'Conjugate transpose', 00445 $ 'Conjugate transpose', LASTC, K, LASTV-K, 00446 $ ONE, C( K+1, 1 ), LDC, V( 1, K+1 ), LDV, 00447 $ ONE, WORK, LDWORK ) 00448 END IF 00449 * 00450 * W := W * T**H or W * T 00451 * 00452 CALL ZTRMM( 'Right', 'Upper', TRANST, 'Non-unit', 00453 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00454 * 00455 * C := C - V**H * W**H 00456 * 00457 IF( LASTV.GT.K ) THEN 00458 * 00459 * C2 := C2 - V2**H * W**H 00460 * 00461 CALL ZGEMM( 'Conjugate transpose', 00462 $ 'Conjugate transpose', LASTV-K, LASTC, K, 00463 $ -ONE, V( 1, K+1 ), LDV, WORK, LDWORK, 00464 $ ONE, C( K+1, 1 ), LDC ) 00465 END IF 00466 * 00467 * W := W * V1 00468 * 00469 CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00470 $ LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00471 * 00472 * C1 := C1 - W**H 00473 * 00474 DO 150 J = 1, K 00475 DO 140 I = 1, LASTC 00476 C( J, I ) = C( J, I ) - DCONJG( WORK( I, J ) ) 00477 140 CONTINUE 00478 150 CONTINUE 00479 * 00480 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00481 * 00482 * Form C * H or C * H**H where C = ( C1 C2 ) 00483 * 00484 LASTV = MAX( K, ILAZLC( K, N, V, LDV ) ) 00485 LASTC = ILAZLR( M, LASTV, C, LDC ) 00486 * 00487 * W := C * V**H = (C1*V1**H + C2*V2**H) (stored in WORK) 00488 * 00489 * W := C1 00490 * 00491 DO 160 J = 1, K 00492 CALL ZCOPY( LASTC, C( 1, J ), 1, WORK( 1, J ), 1 ) 00493 160 CONTINUE 00494 * 00495 * W := W * V1**H 00496 * 00497 CALL ZTRMM( 'Right', 'Upper', 'Conjugate transpose', 00498 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00499 IF( LASTV.GT.K ) THEN 00500 * 00501 * W := W + C2 * V2**H 00502 * 00503 CALL ZGEMM( 'No transpose', 'Conjugate transpose', 00504 $ LASTC, K, LASTV-K, ONE, C( 1, K+1 ), LDC, 00505 $ V( 1, K+1 ), LDV, ONE, WORK, LDWORK ) 00506 END IF 00507 * 00508 * W := W * T or W * T**H 00509 * 00510 CALL ZTRMM( 'Right', 'Upper', TRANS, 'Non-unit', 00511 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00512 * 00513 * C := C - W * V 00514 * 00515 IF( LASTV.GT.K ) THEN 00516 * 00517 * C2 := C2 - W * V2 00518 * 00519 CALL ZGEMM( 'No transpose', 'No transpose', 00520 $ LASTC, LASTV-K, K, 00521 $ -ONE, WORK, LDWORK, V( 1, K+1 ), LDV, 00522 $ ONE, C( 1, K+1 ), LDC ) 00523 END IF 00524 * 00525 * W := W * V1 00526 * 00527 CALL ZTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00528 $ LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00529 * 00530 * C1 := C1 - W 00531 * 00532 DO 180 J = 1, K 00533 DO 170 I = 1, LASTC 00534 C( I, J ) = C( I, J ) - WORK( I, J ) 00535 170 CONTINUE 00536 180 CONTINUE 00537 * 00538 END IF 00539 * 00540 ELSE 00541 * 00542 * Let V = ( V1 V2 ) (V2: last K columns) 00543 * where V2 is unit lower triangular. 00544 * 00545 IF( LSAME( SIDE, 'L' ) ) THEN 00546 * 00547 * Form H * C or H**H * C where C = ( C1 ) 00548 * ( C2 ) 00549 * 00550 LASTV = MAX( K, ILAZLC( K, M, V, LDV ) ) 00551 LASTC = ILAZLC( LASTV, N, C, LDC ) 00552 * 00553 * W := C**H * V**H = (C1**H * V1**H + C2**H * V2**H) (stored in WORK) 00554 * 00555 * W := C2**H 00556 * 00557 DO 190 J = 1, K 00558 CALL ZCOPY( LASTC, C( LASTV-K+J, 1 ), LDC, 00559 $ WORK( 1, J ), 1 ) 00560 CALL ZLACGV( LASTC, WORK( 1, J ), 1 ) 00561 190 CONTINUE 00562 * 00563 * W := W * V2**H 00564 * 00565 CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose', 00566 $ 'Unit', LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00567 $ WORK, LDWORK ) 00568 IF( LASTV.GT.K ) THEN 00569 * 00570 * W := W + C1**H * V1**H 00571 * 00572 CALL ZGEMM( 'Conjugate transpose', 00573 $ 'Conjugate transpose', LASTC, K, LASTV-K, 00574 $ ONE, C, LDC, V, LDV, ONE, WORK, LDWORK ) 00575 END IF 00576 * 00577 * W := W * T**H or W * T 00578 * 00579 CALL ZTRMM( 'Right', 'Lower', TRANST, 'Non-unit', 00580 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00581 * 00582 * C := C - V**H * W**H 00583 * 00584 IF( LASTV.GT.K ) THEN 00585 * 00586 * C1 := C1 - V1**H * W**H 00587 * 00588 CALL ZGEMM( 'Conjugate transpose', 00589 $ 'Conjugate transpose', LASTV-K, LASTC, K, 00590 $ -ONE, V, LDV, WORK, LDWORK, ONE, C, LDC ) 00591 END IF 00592 * 00593 * W := W * V2 00594 * 00595 CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', 00596 $ LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00597 $ WORK, LDWORK ) 00598 * 00599 * C2 := C2 - W**H 00600 * 00601 DO 210 J = 1, K 00602 DO 200 I = 1, LASTC 00603 C( LASTV-K+J, I ) = C( LASTV-K+J, I ) - 00604 $ DCONJG( WORK( I, J ) ) 00605 200 CONTINUE 00606 210 CONTINUE 00607 * 00608 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00609 * 00610 * Form C * H or C * H**H where C = ( C1 C2 ) 00611 * 00612 LASTV = MAX( K, ILAZLC( K, N, V, LDV ) ) 00613 LASTC = ILAZLR( M, LASTV, C, LDC ) 00614 * 00615 * W := C * V**H = (C1*V1**H + C2*V2**H) (stored in WORK) 00616 * 00617 * W := C2 00618 * 00619 DO 220 J = 1, K 00620 CALL ZCOPY( LASTC, C( 1, LASTV-K+J ), 1, 00621 $ WORK( 1, J ), 1 ) 00622 220 CONTINUE 00623 * 00624 * W := W * V2**H 00625 * 00626 CALL ZTRMM( 'Right', 'Lower', 'Conjugate transpose', 00627 $ 'Unit', LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00628 $ WORK, LDWORK ) 00629 IF( LASTV.GT.K ) THEN 00630 * 00631 * W := W + C1 * V1**H 00632 * 00633 CALL ZGEMM( 'No transpose', 'Conjugate transpose', 00634 $ LASTC, K, LASTV-K, ONE, C, LDC, V, LDV, ONE, 00635 $ WORK, LDWORK ) 00636 END IF 00637 * 00638 * W := W * T or W * T**H 00639 * 00640 CALL ZTRMM( 'Right', 'Lower', TRANS, 'Non-unit', 00641 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00642 * 00643 * C := C - W * V 00644 * 00645 IF( LASTV.GT.K ) THEN 00646 * 00647 * C1 := C1 - W * V1 00648 * 00649 CALL ZGEMM( 'No transpose', 'No transpose', 00650 $ LASTC, LASTV-K, K, -ONE, WORK, LDWORK, V, LDV, 00651 $ ONE, C, LDC ) 00652 END IF 00653 * 00654 * W := W * V2 00655 * 00656 CALL ZTRMM( 'Right', 'Lower', 'No transpose', 'Unit', 00657 $ LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00658 $ WORK, LDWORK ) 00659 * 00660 * C1 := C1 - W 00661 * 00662 DO 240 J = 1, K 00663 DO 230 I = 1, LASTC 00664 C( I, LASTV-K+J ) = C( I, LASTV-K+J ) 00665 $ - WORK( I, J ) 00666 230 CONTINUE 00667 240 CONTINUE 00668 * 00669 END IF 00670 * 00671 END IF 00672 END IF 00673 * 00674 RETURN 00675 * 00676 * End of ZLARFB 00677 * 00678 END
1.7.3 
