 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE CLARFB( 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 C( LDC, * ), T( LDT, * ), V( LDV, * ), 00016 $ WORK( LDWORK, * ) 00017 * .. 00018 * 00019 * Purpose 00020 * ======= 00021 * 00022 * CLARFB 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 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 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 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 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 ONE 00120 PARAMETER ( ONE = ( 1.0E+0, 0.0E+0 ) ) 00121 * .. 00122 * .. Local Scalars .. 00123 CHARACTER TRANST 00124 INTEGER I, J, LASTV, LASTC 00125 * .. 00126 * .. External Functions .. 00127 LOGICAL LSAME 00128 INTEGER ILACLR, ILACLC 00129 EXTERNAL LSAME, ILACLR, ILACLC 00130 * .. 00131 * .. External Subroutines .. 00132 EXTERNAL CCOPY, CGEMM, CLACGV, CTRMM 00133 * .. 00134 * .. Intrinsic Functions .. 00135 INTRINSIC CONJG 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, ILACLR( M, K, V, LDV ) ) 00164 LASTC = ILACLC( 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 CCOPY( LASTC, C( J, 1 ), LDC, WORK( 1, J ), 1 ) 00172 CALL CLACGV( LASTC, WORK( 1, J ), 1 ) 00173 10 CONTINUE 00174 * 00175 * W := W * V1 00176 * 00177 CALL CTRMM( '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 CGEMM( '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 CTRMM( '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 CGEMM( 'No transpose', 'Conjugate transpose', 00200 $ LASTV-K, LASTC, K, -ONE, V( K+1, 1 ), LDV, 00201 $ WORK, LDWORK, ONE, C( K+1, 1 ), LDC ) 00202 END IF 00203 * 00204 * W := W * V1**H 00205 * 00206 CALL CTRMM( 'Right', 'Lower', 'Conjugate transpose', 00207 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00208 * 00209 * C1 := C1 - W**H 00210 * 00211 DO 30 J = 1, K 00212 DO 20 I = 1, LASTC 00213 C( J, I ) = C( J, I ) - CONJG( WORK( I, J ) ) 00214 20 CONTINUE 00215 30 CONTINUE 00216 * 00217 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00218 * 00219 * Form C * H or C * H**H where C = ( C1 C2 ) 00220 * 00221 LASTV = MAX( K, ILACLR( N, K, V, LDV ) ) 00222 LASTC = ILACLR( M, LASTV, C, LDC ) 00223 * 00224 * W := C * V = (C1*V1 + C2*V2) (stored in WORK) 00225 * 00226 * W := C1 00227 * 00228 DO 40 J = 1, K 00229 CALL CCOPY( LASTC, C( 1, J ), 1, WORK( 1, J ), 1 ) 00230 40 CONTINUE 00231 * 00232 * W := W * V1 00233 * 00234 CALL CTRMM( 'Right', 'Lower', 'No transpose', 'Unit', 00235 $ LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00236 IF( LASTV.GT.K ) THEN 00237 * 00238 * W := W + C2 * V2 00239 * 00240 CALL CGEMM( 'No transpose', 'No transpose', 00241 $ LASTC, K, LASTV-K, 00242 $ ONE, C( 1, K+1 ), LDC, V( K+1, 1 ), LDV, 00243 $ ONE, WORK, LDWORK ) 00244 END IF 00245 * 00246 * W := W * T or W * T**H 00247 * 00248 CALL CTRMM( 'Right', 'Upper', TRANS, 'Non-unit', 00249 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00250 * 00251 * C := C - W * V**H 00252 * 00253 IF( LASTV.GT.K ) THEN 00254 * 00255 * C2 := C2 - W * V2**H 00256 * 00257 CALL CGEMM( 'No transpose', 'Conjugate transpose', 00258 $ LASTC, LASTV-K, K, 00259 $ -ONE, WORK, LDWORK, V( K+1, 1 ), LDV, 00260 $ ONE, C( 1, K+1 ), LDC ) 00261 END IF 00262 * 00263 * W := W * V1**H 00264 * 00265 CALL CTRMM( 'Right', 'Lower', 'Conjugate transpose', 00266 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00267 * 00268 * C1 := C1 - W 00269 * 00270 DO 60 J = 1, K 00271 DO 50 I = 1, LASTC 00272 C( I, J ) = C( I, J ) - WORK( I, J ) 00273 50 CONTINUE 00274 60 CONTINUE 00275 END IF 00276 * 00277 ELSE 00278 * 00279 * Let V = ( V1 ) 00280 * ( V2 ) (last K rows) 00281 * where V2 is unit upper triangular. 00282 * 00283 IF( LSAME( SIDE, 'L' ) ) THEN 00284 * 00285 * Form H * C or H**H * C where C = ( C1 ) 00286 * ( C2 ) 00287 * 00288 LASTV = MAX( K, ILACLR( M, K, V, LDV ) ) 00289 LASTC = ILACLC( LASTV, N, C, LDC ) 00290 * 00291 * W := C**H * V = (C1**H * V1 + C2**H * V2) (stored in WORK) 00292 * 00293 * W := C2**H 00294 * 00295 DO 70 J = 1, K 00296 CALL CCOPY( LASTC, C( LASTV-K+J, 1 ), LDC, 00297 $ WORK( 1, J ), 1 ) 00298 CALL CLACGV( LASTC, WORK( 1, J ), 1 ) 00299 70 CONTINUE 00300 * 00301 * W := W * V2 00302 * 00303 CALL CTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00304 $ LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00305 $ WORK, LDWORK ) 00306 IF( LASTV.GT.K ) THEN 00307 * 00308 * W := W + C1**H*V1 00309 * 00310 CALL CGEMM( 'Conjugate transpose', 'No transpose', 00311 $ LASTC, K, LASTV-K, ONE, C, LDC, V, LDV, 00312 $ ONE, WORK, LDWORK ) 00313 END IF 00314 * 00315 * W := W * T**H or W * T 00316 * 00317 CALL CTRMM( 'Right', 'Lower', TRANST, 'Non-unit', 00318 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00319 * 00320 * C := C - V * W**H 00321 * 00322 IF( LASTV.GT.K ) THEN 00323 * 00324 * C1 := C1 - V1 * W**H 00325 * 00326 CALL CGEMM( 'No transpose', 'Conjugate transpose', 00327 $ LASTV-K, LASTC, K, -ONE, V, LDV, WORK, LDWORK, 00328 $ ONE, C, LDC ) 00329 END IF 00330 * 00331 * W := W * V2**H 00332 * 00333 CALL CTRMM( 'Right', 'Upper', 'Conjugate transpose', 00334 $ 'Unit', LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00335 $ WORK, LDWORK ) 00336 * 00337 * C2 := C2 - W**H 00338 * 00339 DO 90 J = 1, K 00340 DO 80 I = 1, LASTC 00341 C( LASTV-K+J, I ) = C( LASTV-K+J, I ) - 00342 $ CONJG( WORK( I, J ) ) 00343 80 CONTINUE 00344 90 CONTINUE 00345 * 00346 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00347 * 00348 * Form C * H or C * H**H where C = ( C1 C2 ) 00349 * 00350 LASTV = MAX( K, ILACLR( N, K, V, LDV ) ) 00351 LASTC = ILACLR( M, LASTV, C, LDC ) 00352 * 00353 * W := C * V = (C1*V1 + C2*V2) (stored in WORK) 00354 * 00355 * W := C2 00356 * 00357 DO 100 J = 1, K 00358 CALL CCOPY( LASTC, C( 1, LASTV-K+J ), 1, 00359 $ WORK( 1, J ), 1 ) 00360 100 CONTINUE 00361 * 00362 * W := W * V2 00363 * 00364 CALL CTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00365 $ LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00366 $ WORK, LDWORK ) 00367 IF( LASTV.GT.K ) THEN 00368 * 00369 * W := W + C1 * V1 00370 * 00371 CALL CGEMM( 'No transpose', 'No transpose', 00372 $ LASTC, K, LASTV-K, 00373 $ ONE, C, LDC, V, LDV, ONE, WORK, LDWORK ) 00374 END IF 00375 * 00376 * W := W * T or W * T**H 00377 * 00378 CALL CTRMM( 'Right', 'Lower', TRANS, 'Non-unit', 00379 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00380 * 00381 * C := C - W * V**H 00382 * 00383 IF( LASTV.GT.K ) THEN 00384 * 00385 * C1 := C1 - W * V1**H 00386 * 00387 CALL CGEMM( 'No transpose', 'Conjugate transpose', 00388 $ LASTC, LASTV-K, K, -ONE, WORK, LDWORK, V, LDV, 00389 $ ONE, C, LDC ) 00390 END IF 00391 * 00392 * W := W * V2**H 00393 * 00394 CALL CTRMM( 'Right', 'Upper', 'Conjugate transpose', 00395 $ 'Unit', LASTC, K, ONE, V( LASTV-K+1, 1 ), LDV, 00396 $ WORK, LDWORK ) 00397 * 00398 * C2 := C2 - W 00399 * 00400 DO 120 J = 1, K 00401 DO 110 I = 1, LASTC 00402 C( I, LASTV-K+J ) = C( I, LASTV-K+J ) 00403 $ - WORK( I, J ) 00404 110 CONTINUE 00405 120 CONTINUE 00406 END IF 00407 END IF 00408 * 00409 ELSE IF( LSAME( STOREV, 'R' ) ) THEN 00410 * 00411 IF( LSAME( DIRECT, 'F' ) ) THEN 00412 * 00413 * Let V = ( V1 V2 ) (V1: first K columns) 00414 * where V1 is unit upper triangular. 00415 * 00416 IF( LSAME( SIDE, 'L' ) ) THEN 00417 * 00418 * Form H * C or H**H * C where C = ( C1 ) 00419 * ( C2 ) 00420 * 00421 LASTV = MAX( K, ILACLC( K, M, V, LDV ) ) 00422 LASTC = ILACLC( LASTV, N, C, LDC ) 00423 * 00424 * W := C**H * V**H = (C1**H * V1**H + C2**H * V2**H) (stored in WORK) 00425 * 00426 * W := C1**H 00427 * 00428 DO 130 J = 1, K 00429 CALL CCOPY( LASTC, C( J, 1 ), LDC, WORK( 1, J ), 1 ) 00430 CALL CLACGV( LASTC, WORK( 1, J ), 1 ) 00431 130 CONTINUE 00432 * 00433 * W := W * V1**H 00434 * 00435 CALL CTRMM( 'Right', 'Upper', 'Conjugate transpose', 00436 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00437 IF( LASTV.GT.K ) THEN 00438 * 00439 * W := W + C2**H*V2**H 00440 * 00441 CALL CGEMM( 'Conjugate transpose', 00442 $ 'Conjugate transpose', LASTC, K, LASTV-K, 00443 $ ONE, C( K+1, 1 ), LDC, V( 1, K+1 ), LDV, 00444 $ ONE, WORK, LDWORK ) 00445 END IF 00446 * 00447 * W := W * T**H or W * T 00448 * 00449 CALL CTRMM( 'Right', 'Upper', TRANST, 'Non-unit', 00450 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00451 * 00452 * C := C - V**H * W**H 00453 * 00454 IF( LASTV.GT.K ) THEN 00455 * 00456 * C2 := C2 - V2**H * W**H 00457 * 00458 CALL CGEMM( 'Conjugate transpose', 00459 $ 'Conjugate transpose', LASTV-K, LASTC, K, 00460 $ -ONE, V( 1, K+1 ), LDV, WORK, LDWORK, 00461 $ ONE, C( K+1, 1 ), LDC ) 00462 END IF 00463 * 00464 * W := W * V1 00465 * 00466 CALL CTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00467 $ LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00468 * 00469 * C1 := C1 - W**H 00470 * 00471 DO 150 J = 1, K 00472 DO 140 I = 1, LASTC 00473 C( J, I ) = C( J, I ) - CONJG( WORK( I, J ) ) 00474 140 CONTINUE 00475 150 CONTINUE 00476 * 00477 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00478 * 00479 * Form C * H or C * H**H where C = ( C1 C2 ) 00480 * 00481 LASTV = MAX( K, ILACLC( K, N, V, LDV ) ) 00482 LASTC = ILACLR( M, LASTV, C, LDC ) 00483 * 00484 * W := C * V**H = (C1*V1**H + C2*V2**H) (stored in WORK) 00485 * 00486 * W := C1 00487 * 00488 DO 160 J = 1, K 00489 CALL CCOPY( LASTC, C( 1, J ), 1, WORK( 1, J ), 1 ) 00490 160 CONTINUE 00491 * 00492 * W := W * V1**H 00493 * 00494 CALL CTRMM( 'Right', 'Upper', 'Conjugate transpose', 00495 $ 'Unit', LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00496 IF( LASTV.GT.K ) THEN 00497 * 00498 * W := W + C2 * V2**H 00499 * 00500 CALL CGEMM( 'No transpose', 'Conjugate transpose', 00501 $ LASTC, K, LASTV-K, ONE, C( 1, K+1 ), LDC, 00502 $ V( 1, K+1 ), LDV, ONE, WORK, LDWORK ) 00503 END IF 00504 * 00505 * W := W * T or W * T**H 00506 * 00507 CALL CTRMM( 'Right', 'Upper', TRANS, 'Non-unit', 00508 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00509 * 00510 * C := C - W * V 00511 * 00512 IF( LASTV.GT.K ) THEN 00513 * 00514 * C2 := C2 - W * V2 00515 * 00516 CALL CGEMM( 'No transpose', 'No transpose', 00517 $ LASTC, LASTV-K, K, 00518 $ -ONE, WORK, LDWORK, V( 1, K+1 ), LDV, 00519 $ ONE, C( 1, K+1 ), LDC ) 00520 END IF 00521 * 00522 * W := W * V1 00523 * 00524 CALL CTRMM( 'Right', 'Upper', 'No transpose', 'Unit', 00525 $ LASTC, K, ONE, V, LDV, WORK, LDWORK ) 00526 * 00527 * C1 := C1 - W 00528 * 00529 DO 180 J = 1, K 00530 DO 170 I = 1, LASTC 00531 C( I, J ) = C( I, J ) - WORK( I, J ) 00532 170 CONTINUE 00533 180 CONTINUE 00534 * 00535 END IF 00536 * 00537 ELSE 00538 * 00539 * Let V = ( V1 V2 ) (V2: last K columns) 00540 * where V2 is unit lower triangular. 00541 * 00542 IF( LSAME( SIDE, 'L' ) ) THEN 00543 * 00544 * Form H * C or H**H * C where C = ( C1 ) 00545 * ( C2 ) 00546 * 00547 LASTV = MAX( K, ILACLC( K, M, V, LDV ) ) 00548 LASTC = ILACLC( LASTV, N, C, LDC ) 00549 * 00550 * W := C**H * V**H = (C1**H * V1**H + C2**H * V2**H) (stored in WORK) 00551 * 00552 * W := C2**H 00553 * 00554 DO 190 J = 1, K 00555 CALL CCOPY( LASTC, C( LASTV-K+J, 1 ), LDC, 00556 $ WORK( 1, J ), 1 ) 00557 CALL CLACGV( LASTC, WORK( 1, J ), 1 ) 00558 190 CONTINUE 00559 * 00560 * W := W * V2**H 00561 * 00562 CALL CTRMM( 'Right', 'Lower', 'Conjugate transpose', 00563 $ 'Unit', LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00564 $ WORK, LDWORK ) 00565 IF( LASTV.GT.K ) THEN 00566 * 00567 * W := W + C1**H * V1**H 00568 * 00569 CALL CGEMM( 'Conjugate transpose', 00570 $ 'Conjugate transpose', LASTC, K, LASTV-K, 00571 $ ONE, C, LDC, V, LDV, ONE, WORK, LDWORK ) 00572 END IF 00573 * 00574 * W := W * T**H or W * T 00575 * 00576 CALL CTRMM( 'Right', 'Lower', TRANST, 'Non-unit', 00577 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00578 * 00579 * C := C - V**H * W**H 00580 * 00581 IF( LASTV.GT.K ) THEN 00582 * 00583 * C1 := C1 - V1**H * W**H 00584 * 00585 CALL CGEMM( 'Conjugate transpose', 00586 $ 'Conjugate transpose', LASTV-K, LASTC, K, 00587 $ -ONE, V, LDV, WORK, LDWORK, ONE, C, LDC ) 00588 END IF 00589 * 00590 * W := W * V2 00591 * 00592 CALL CTRMM( 'Right', 'Lower', 'No transpose', 'Unit', 00593 $ LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00594 $ WORK, LDWORK ) 00595 * 00596 * C2 := C2 - W**H 00597 * 00598 DO 210 J = 1, K 00599 DO 200 I = 1, LASTC 00600 C( LASTV-K+J, I ) = C( LASTV-K+J, I ) - 00601 $ CONJG( WORK( I, J ) ) 00602 200 CONTINUE 00603 210 CONTINUE 00604 * 00605 ELSE IF( LSAME( SIDE, 'R' ) ) THEN 00606 * 00607 * Form C * H or C * H**H where C = ( C1 C2 ) 00608 * 00609 LASTV = MAX( K, ILACLC( K, N, V, LDV ) ) 00610 LASTC = ILACLR( M, LASTV, C, LDC ) 00611 * 00612 * W := C * V**H = (C1*V1**H + C2*V2**H) (stored in WORK) 00613 * 00614 * W := C2 00615 * 00616 DO 220 J = 1, K 00617 CALL CCOPY( LASTC, C( 1, LASTV-K+J ), 1, 00618 $ WORK( 1, J ), 1 ) 00619 220 CONTINUE 00620 * 00621 * W := W * V2**H 00622 * 00623 CALL CTRMM( 'Right', 'Lower', 'Conjugate transpose', 00624 $ 'Unit', LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00625 $ WORK, LDWORK ) 00626 IF( LASTV.GT.K ) THEN 00627 * 00628 * W := W + C1 * V1**H 00629 * 00630 CALL CGEMM( 'No transpose', 'Conjugate transpose', 00631 $ LASTC, K, LASTV-K, ONE, C, LDC, V, LDV, ONE, 00632 $ WORK, LDWORK ) 00633 END IF 00634 * 00635 * W := W * T or W * T**H 00636 * 00637 CALL CTRMM( 'Right', 'Lower', TRANS, 'Non-unit', 00638 $ LASTC, K, ONE, T, LDT, WORK, LDWORK ) 00639 * 00640 * C := C - W * V 00641 * 00642 IF( LASTV.GT.K ) THEN 00643 * 00644 * C1 := C1 - W * V1 00645 * 00646 CALL CGEMM( 'No transpose', 'No transpose', 00647 $ LASTC, LASTV-K, K, -ONE, WORK, LDWORK, V, LDV, 00648 $ ONE, C, LDC ) 00649 END IF 00650 * 00651 * W := W * V2 00652 * 00653 CALL CTRMM( 'Right', 'Lower', 'No transpose', 'Unit', 00654 $ LASTC, K, ONE, V( 1, LASTV-K+1 ), LDV, 00655 $ WORK, LDWORK ) 00656 * 00657 * C1 := C1 - W 00658 * 00659 DO 240 J = 1, K 00660 DO 230 I = 1, LASTC 00661 C( I, LASTV-K+J ) = C( I, LASTV-K+J ) 00662 $ - WORK( I, J ) 00663 230 CONTINUE 00664 240 CONTINUE 00665 * 00666 END IF 00667 * 00668 END IF 00669 END IF 00670 * 00671 RETURN 00672 * 00673 * End of CLARFB 00674 * 00675 END
1.7.3 
