1*> \brief \b ZLAVHE 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* SUBROUTINE ZLAVHE( UPLO, TRANS, DIAG, N, NRHS, A, LDA, IPIV, B, 12* LDB, INFO ) 13* 14* .. Scalar Arguments .. 15* CHARACTER DIAG, TRANS, UPLO 16* INTEGER INFO, LDA, LDB, N, NRHS 17* .. 18* .. Array Arguments .. 19* INTEGER IPIV( * ) 20* COMPLEX*16 A( LDA, * ), B( LDB, * ) 21* .. 22* 23* 24*> \par Purpose: 25* ============= 26*> 27*> \verbatim 28*> 29*> ZLAVHE performs one of the matrix-vector operations 30*> x := A*x or x := A^H*x, 31*> where x is an N element vector and A is one of the factors 32*> from the symmetric factorization computed by ZHETRF. 33*> ZHETRF produces a factorization of the form 34*> U * D * U^H or L * D * L^H, 35*> where U (or L) is a product of permutation and unit upper (lower) 36*> triangular matrices, U^H (or L^H) is the conjugate transpose of 37*> U (or L), and D is Hermitian and block diagonal with 1 x 1 and 38*> 2 x 2 diagonal blocks. The multipliers for the transformations 39*> and the upper or lower triangular parts of the diagonal blocks 40*> are stored in the leading upper or lower triangle of the 2-D 41*> array A. 42*> 43*> If TRANS = 'N' or 'n', ZLAVHE multiplies either by U or U * D 44*> (or L or L * D). 45*> If TRANS = 'C' or 'c', ZLAVHE multiplies either by U^H or D * U^H 46*> (or L^H or D * L^H ). 47*> \endverbatim 48* 49* Arguments: 50* ========== 51* 52*> \verbatim 53*> UPLO - CHARACTER*1 54*> On entry, UPLO specifies whether the triangular matrix 55*> stored in A is upper or lower triangular. 56*> UPLO = 'U' or 'u' The matrix is upper triangular. 57*> UPLO = 'L' or 'l' The matrix is lower triangular. 58*> Unchanged on exit. 59*> 60*> TRANS - CHARACTER*1 61*> On entry, TRANS specifies the operation to be performed as 62*> follows: 63*> TRANS = 'N' or 'n' x := A*x. 64*> TRANS = 'C' or 'c' x := A^H*x. 65*> Unchanged on exit. 66*> 67*> DIAG - CHARACTER*1 68*> On entry, DIAG specifies whether the diagonal blocks are 69*> assumed to be unit matrices: 70*> DIAG = 'U' or 'u' Diagonal blocks are unit matrices. 71*> DIAG = 'N' or 'n' Diagonal blocks are non-unit. 72*> Unchanged on exit. 73*> 74*> N - INTEGER 75*> On entry, N specifies the order of the matrix A. 76*> N must be at least zero. 77*> Unchanged on exit. 78*> 79*> NRHS - INTEGER 80*> On entry, NRHS specifies the number of right hand sides, 81*> i.e., the number of vectors x to be multiplied by A. 82*> NRHS must be at least zero. 83*> Unchanged on exit. 84*> 85*> A - COMPLEX*16 array, dimension( LDA, N ) 86*> On entry, A contains a block diagonal matrix and the 87*> multipliers of the transformations used to obtain it, 88*> stored as a 2-D triangular matrix. 89*> Unchanged on exit. 90*> 91*> LDA - INTEGER 92*> On entry, LDA specifies the first dimension of A as declared 93*> in the calling ( sub ) program. LDA must be at least 94*> max( 1, N ). 95*> Unchanged on exit. 96*> 97*> IPIV - INTEGER array, dimension( N ) 98*> On entry, IPIV contains the vector of pivot indices as 99*> determined by ZSYTRF or ZHETRF. 100*> If IPIV( K ) = K, no interchange was done. 101*> If IPIV( K ) <> K but IPIV( K ) > 0, then row K was inter- 102*> changed with row IPIV( K ) and a 1 x 1 pivot block was used. 103*> If IPIV( K ) < 0 and UPLO = 'U', then row K-1 was exchanged 104*> with row | IPIV( K ) | and a 2 x 2 pivot block was used. 105*> If IPIV( K ) < 0 and UPLO = 'L', then row K+1 was exchanged 106*> with row | IPIV( K ) | and a 2 x 2 pivot block was used. 107*> 108*> B - COMPLEX*16 array, dimension( LDB, NRHS ) 109*> On entry, B contains NRHS vectors of length N. 110*> On exit, B is overwritten with the product A * B. 111*> 112*> LDB - INTEGER 113*> On entry, LDB contains the leading dimension of B as 114*> declared in the calling program. LDB must be at least 115*> max( 1, N ). 116*> Unchanged on exit. 117*> 118*> INFO - INTEGER 119*> INFO is the error flag. 120*> On exit, a value of 0 indicates a successful exit. 121*> A negative value, say -K, indicates that the K-th argument 122*> has an illegal value. 123*> \endverbatim 124* 125* Authors: 126* ======== 127* 128*> \author Univ. of Tennessee 129*> \author Univ. of California Berkeley 130*> \author Univ. of Colorado Denver 131*> \author NAG Ltd. 132* 133*> \date November 2011 134* 135*> \ingroup complex16_lin 136* 137* ===================================================================== 138 SUBROUTINE ZLAVHE( UPLO, TRANS, DIAG, N, NRHS, A, LDA, IPIV, B, 139 $ LDB, INFO ) 140* 141* -- LAPACK test routine (version 3.4.0) -- 142* -- LAPACK is a software package provided by Univ. of Tennessee, -- 143* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 144* November 2011 145* 146* .. Scalar Arguments .. 147 CHARACTER DIAG, TRANS, UPLO 148 INTEGER INFO, LDA, LDB, N, NRHS 149* .. 150* .. Array Arguments .. 151 INTEGER IPIV( * ) 152 COMPLEX*16 A( LDA, * ), B( LDB, * ) 153* .. 154* 155* ===================================================================== 156* 157* .. Parameters .. 158 COMPLEX*16 ONE 159 PARAMETER ( ONE = ( 1.0D+0, 0.0D+0 ) ) 160* .. 161* .. Local Scalars .. 162 LOGICAL NOUNIT 163 INTEGER J, K, KP 164 COMPLEX*16 D11, D12, D21, D22, T1, T2 165* .. 166* .. External Functions .. 167 LOGICAL LSAME 168 EXTERNAL LSAME 169* .. 170* .. External Subroutines .. 171 EXTERNAL XERBLA, ZGEMV, ZGERU, ZLACGV, ZSCAL, ZSWAP 172* .. 173* .. Intrinsic Functions .. 174 INTRINSIC ABS, DCONJG, MAX 175* .. 176* .. Executable Statements .. 177* 178* Test the input parameters. 179* 180 INFO = 0 181 IF( .NOT.LSAME( UPLO, 'U' ) .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN 182 INFO = -1 183 ELSE IF( .NOT.LSAME( TRANS, 'N' ) .AND. .NOT.LSAME( TRANS, 'C' ) ) 184 $ THEN 185 INFO = -2 186 ELSE IF( .NOT.LSAME( DIAG, 'U' ) .AND. .NOT.LSAME( DIAG, 'N' ) ) 187 $ THEN 188 INFO = -3 189 ELSE IF( N.LT.0 ) THEN 190 INFO = -4 191 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 192 INFO = -6 193 ELSE IF( LDB.LT.MAX( 1, N ) ) THEN 194 INFO = -9 195 END IF 196 IF( INFO.NE.0 ) THEN 197 CALL XERBLA( 'ZLAVHE ', -INFO ) 198 RETURN 199 END IF 200* 201* Quick return if possible. 202* 203 IF( N.EQ.0 ) 204 $ RETURN 205* 206 NOUNIT = LSAME( DIAG, 'N' ) 207*------------------------------------------ 208* 209* Compute B := A * B (No transpose) 210* 211*------------------------------------------ 212 IF( LSAME( TRANS, 'N' ) ) THEN 213* 214* Compute B := U*B 215* where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1)) 216* 217 IF( LSAME( UPLO, 'U' ) ) THEN 218* 219* Loop forward applying the transformations. 220* 221 K = 1 222 10 CONTINUE 223 IF( K.GT.N ) 224 $ GO TO 30 225 IF( IPIV( K ).GT.0 ) THEN 226* 227* 1 x 1 pivot block 228* 229* Multiply by the diagonal element if forming U * D. 230* 231 IF( NOUNIT ) 232 $ CALL ZSCAL( NRHS, A( K, K ), B( K, 1 ), LDB ) 233* 234* Multiply by P(K) * inv(U(K)) if K > 1. 235* 236 IF( K.GT.1 ) THEN 237* 238* Apply the transformation. 239* 240 CALL ZGERU( K-1, NRHS, ONE, A( 1, K ), 1, B( K, 1 ), 241 $ LDB, B( 1, 1 ), LDB ) 242* 243* Interchange if P(K) != I. 244* 245 KP = IPIV( K ) 246 IF( KP.NE.K ) 247 $ CALL ZSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) 248 END IF 249 K = K + 1 250 ELSE 251* 252* 2 x 2 pivot block 253* 254* Multiply by the diagonal block if forming U * D. 255* 256 IF( NOUNIT ) THEN 257 D11 = A( K, K ) 258 D22 = A( K+1, K+1 ) 259 D12 = A( K, K+1 ) 260 D21 = DCONJG( D12 ) 261 DO 20 J = 1, NRHS 262 T1 = B( K, J ) 263 T2 = B( K+1, J ) 264 B( K, J ) = D11*T1 + D12*T2 265 B( K+1, J ) = D21*T1 + D22*T2 266 20 CONTINUE 267 END IF 268* 269* Multiply by P(K) * inv(U(K)) if K > 1. 270* 271 IF( K.GT.1 ) THEN 272* 273* Apply the transformations. 274* 275 CALL ZGERU( K-1, NRHS, ONE, A( 1, K ), 1, B( K, 1 ), 276 $ LDB, B( 1, 1 ), LDB ) 277 CALL ZGERU( K-1, NRHS, ONE, A( 1, K+1 ), 1, 278 $ B( K+1, 1 ), LDB, B( 1, 1 ), LDB ) 279* 280* Interchange if P(K) != I. 281* 282 KP = ABS( IPIV( K ) ) 283 IF( KP.NE.K ) 284 $ CALL ZSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) 285 END IF 286 K = K + 2 287 END IF 288 GO TO 10 289 30 CONTINUE 290* 291* Compute B := L*B 292* where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m)) . 293* 294 ELSE 295* 296* Loop backward applying the transformations to B. 297* 298 K = N 299 40 CONTINUE 300 IF( K.LT.1 ) 301 $ GO TO 60 302* 303* Test the pivot index. If greater than zero, a 1 x 1 304* pivot was used, otherwise a 2 x 2 pivot was used. 305* 306 IF( IPIV( K ).GT.0 ) THEN 307* 308* 1 x 1 pivot block: 309* 310* Multiply by the diagonal element if forming L * D. 311* 312 IF( NOUNIT ) 313 $ CALL ZSCAL( NRHS, A( K, K ), B( K, 1 ), LDB ) 314* 315* Multiply by P(K) * inv(L(K)) if K < N. 316* 317 IF( K.NE.N ) THEN 318 KP = IPIV( K ) 319* 320* Apply the transformation. 321* 322 CALL ZGERU( N-K, NRHS, ONE, A( K+1, K ), 1, B( K, 1 ), 323 $ LDB, B( K+1, 1 ), LDB ) 324* 325* Interchange if a permutation was applied at the 326* K-th step of the factorization. 327* 328 IF( KP.NE.K ) 329 $ CALL ZSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) 330 END IF 331 K = K - 1 332* 333 ELSE 334* 335* 2 x 2 pivot block: 336* 337* Multiply by the diagonal block if forming L * D. 338* 339 IF( NOUNIT ) THEN 340 D11 = A( K-1, K-1 ) 341 D22 = A( K, K ) 342 D21 = A( K, K-1 ) 343 D12 = DCONJG( D21 ) 344 DO 50 J = 1, NRHS 345 T1 = B( K-1, J ) 346 T2 = B( K, J ) 347 B( K-1, J ) = D11*T1 + D12*T2 348 B( K, J ) = D21*T1 + D22*T2 349 50 CONTINUE 350 END IF 351* 352* Multiply by P(K) * inv(L(K)) if K < N. 353* 354 IF( K.NE.N ) THEN 355* 356* Apply the transformation. 357* 358 CALL ZGERU( N-K, NRHS, ONE, A( K+1, K ), 1, B( K, 1 ), 359 $ LDB, B( K+1, 1 ), LDB ) 360 CALL ZGERU( N-K, NRHS, ONE, A( K+1, K-1 ), 1, 361 $ B( K-1, 1 ), LDB, B( K+1, 1 ), LDB ) 362* 363* Interchange if a permutation was applied at the 364* K-th step of the factorization. 365* 366 KP = ABS( IPIV( K ) ) 367 IF( KP.NE.K ) 368 $ CALL ZSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) 369 END IF 370 K = K - 2 371 END IF 372 GO TO 40 373 60 CONTINUE 374 END IF 375*-------------------------------------------------- 376* 377* Compute B := A^H * B (conjugate transpose) 378* 379*-------------------------------------------------- 380 ELSE 381* 382* Form B := U^H*B 383* where U = P(m)*inv(U(m))* ... *P(1)*inv(U(1)) 384* and U^H = inv(U^H(1))*P(1)* ... *inv(U^H(m))*P(m) 385* 386 IF( LSAME( UPLO, 'U' ) ) THEN 387* 388* Loop backward applying the transformations. 389* 390 K = N 391 70 CONTINUE 392 IF( K.LT.1 ) 393 $ GO TO 90 394* 395* 1 x 1 pivot block. 396* 397 IF( IPIV( K ).GT.0 ) THEN 398 IF( K.GT.1 ) THEN 399* 400* Interchange if P(K) != I. 401* 402 KP = IPIV( K ) 403 IF( KP.NE.K ) 404 $ CALL ZSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) 405* 406* Apply the transformation 407* y = y - B' conjg(x), 408* where x is a column of A and y is a row of B. 409* 410 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 411 CALL ZGEMV( 'Conjugate', K-1, NRHS, ONE, B, LDB, 412 $ A( 1, K ), 1, ONE, B( K, 1 ), LDB ) 413 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 414 END IF 415 IF( NOUNIT ) 416 $ CALL ZSCAL( NRHS, A( K, K ), B( K, 1 ), LDB ) 417 K = K - 1 418* 419* 2 x 2 pivot block. 420* 421 ELSE 422 IF( K.GT.2 ) THEN 423* 424* Interchange if P(K) != I. 425* 426 KP = ABS( IPIV( K ) ) 427 IF( KP.NE.K-1 ) 428 $ CALL ZSWAP( NRHS, B( K-1, 1 ), LDB, B( KP, 1 ), 429 $ LDB ) 430* 431* Apply the transformations 432* y = y - B' conjg(x), 433* where x is a block column of A and y is a block 434* row of B. 435* 436 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 437 CALL ZGEMV( 'Conjugate', K-2, NRHS, ONE, B, LDB, 438 $ A( 1, K ), 1, ONE, B( K, 1 ), LDB ) 439 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 440* 441 CALL ZLACGV( NRHS, B( K-1, 1 ), LDB ) 442 CALL ZGEMV( 'Conjugate', K-2, NRHS, ONE, B, LDB, 443 $ A( 1, K-1 ), 1, ONE, B( K-1, 1 ), LDB ) 444 CALL ZLACGV( NRHS, B( K-1, 1 ), LDB ) 445 END IF 446* 447* Multiply by the diagonal block if non-unit. 448* 449 IF( NOUNIT ) THEN 450 D11 = A( K-1, K-1 ) 451 D22 = A( K, K ) 452 D12 = A( K-1, K ) 453 D21 = DCONJG( D12 ) 454 DO 80 J = 1, NRHS 455 T1 = B( K-1, J ) 456 T2 = B( K, J ) 457 B( K-1, J ) = D11*T1 + D12*T2 458 B( K, J ) = D21*T1 + D22*T2 459 80 CONTINUE 460 END IF 461 K = K - 2 462 END IF 463 GO TO 70 464 90 CONTINUE 465* 466* Form B := L^H*B 467* where L = P(1)*inv(L(1))* ... *P(m)*inv(L(m)) 468* and L^H = inv(L^H(m))*P(m)* ... *inv(L^H(1))*P(1) 469* 470 ELSE 471* 472* Loop forward applying the L-transformations. 473* 474 K = 1 475 100 CONTINUE 476 IF( K.GT.N ) 477 $ GO TO 120 478* 479* 1 x 1 pivot block 480* 481 IF( IPIV( K ).GT.0 ) THEN 482 IF( K.LT.N ) THEN 483* 484* Interchange if P(K) != I. 485* 486 KP = IPIV( K ) 487 IF( KP.NE.K ) 488 $ CALL ZSWAP( NRHS, B( K, 1 ), LDB, B( KP, 1 ), LDB ) 489* 490* Apply the transformation 491* 492 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 493 CALL ZGEMV( 'Conjugate', N-K, NRHS, ONE, B( K+1, 1 ), 494 $ LDB, A( K+1, K ), 1, ONE, B( K, 1 ), LDB ) 495 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 496 END IF 497 IF( NOUNIT ) 498 $ CALL ZSCAL( NRHS, A( K, K ), B( K, 1 ), LDB ) 499 K = K + 1 500* 501* 2 x 2 pivot block. 502* 503 ELSE 504 IF( K.LT.N-1 ) THEN 505* 506* Interchange if P(K) != I. 507* 508 KP = ABS( IPIV( K ) ) 509 IF( KP.NE.K+1 ) 510 $ CALL ZSWAP( NRHS, B( K+1, 1 ), LDB, B( KP, 1 ), 511 $ LDB ) 512* 513* Apply the transformation 514* 515 CALL ZLACGV( NRHS, B( K+1, 1 ), LDB ) 516 CALL ZGEMV( 'Conjugate', N-K-1, NRHS, ONE, 517 $ B( K+2, 1 ), LDB, A( K+2, K+1 ), 1, ONE, 518 $ B( K+1, 1 ), LDB ) 519 CALL ZLACGV( NRHS, B( K+1, 1 ), LDB ) 520* 521 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 522 CALL ZGEMV( 'Conjugate', N-K-1, NRHS, ONE, 523 $ B( K+2, 1 ), LDB, A( K+2, K ), 1, ONE, 524 $ B( K, 1 ), LDB ) 525 CALL ZLACGV( NRHS, B( K, 1 ), LDB ) 526 END IF 527* 528* Multiply by the diagonal block if non-unit. 529* 530 IF( NOUNIT ) THEN 531 D11 = A( K, K ) 532 D22 = A( K+1, K+1 ) 533 D21 = A( K+1, K ) 534 D12 = DCONJG( D21 ) 535 DO 110 J = 1, NRHS 536 T1 = B( K, J ) 537 T2 = B( K+1, J ) 538 B( K, J ) = D11*T1 + D12*T2 539 B( K+1, J ) = D21*T1 + D22*T2 540 110 CONTINUE 541 END IF 542 K = K + 2 543 END IF 544 GO TO 100 545 120 CONTINUE 546 END IF 547* 548 END IF 549 RETURN 550* 551* End of ZLAVHE 552* 553 END 554