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