142 COMPLEX*16 A( lda, * ), WORK( * )
149 COMPLEX*16 CONE, CZERO
150 parameter( one = 1.0d+0, cone = ( 1.0d+0, 0.0d+0 ),
151 $ czero = ( 0.0d+0, 0.0d+0 ) )
155 INTEGER J, K, KP, KSTEP
156 DOUBLE PRECISION AK, AKP1, D, T
157 COMPLEX*16 AKKP1, TEMP
162 EXTERNAL lsame, zdotc
168 INTRINSIC abs, dconjg, max, dble
175 upper = lsame( uplo,
'U' )
176 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
178 ELSE IF( n.LT.0 )
THEN
180 ELSE IF( lda.LT.max( 1, n ) )
THEN
184 CALL
xerbla(
'ZHETRI_ROOK', -info )
199 DO 10 info = n, 1, -1
200 IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.czero )
208 IF( ipiv( info ).GT.0 .AND. a( info, info ).EQ.czero )
229 IF( ipiv( k ).GT.0 )
THEN
235 a( k, k ) = one / dble( a( k, k ) )
240 CALL
zcopy( k-1, a( 1, k ), 1, work, 1 )
241 CALL
zhemv( uplo, k-1, -cone, a, lda, work, 1, czero,
243 a( k, k ) = a( k, k ) - dble( zdotc( k-1, work, 1, a( 1,
253 t = abs( a( k, k+1 ) )
254 ak = dble( a( k, k ) ) / t
255 akp1 = dble( a( k+1, k+1 ) ) / t
256 akkp1 = a( k, k+1 ) / t
257 d = t*( ak*akp1-one )
259 a( k+1, k+1 ) = ak / d
260 a( k, k+1 ) = -akkp1 / d
265 CALL
zcopy( k-1, a( 1, k ), 1, work, 1 )
266 CALL
zhemv( uplo, k-1, -cone, a, lda, work, 1, czero,
268 a( k, k ) = a( k, k ) - dble( zdotc( k-1, work, 1, a( 1,
270 a( k, k+1 ) = a( k, k+1 ) -
271 $ zdotc( k-1, a( 1, k ), 1, a( 1, k+1 ), 1 )
272 CALL
zcopy( k-1, a( 1, k+1 ), 1, work, 1 )
273 CALL
zhemv( uplo, k-1, -cone, a, lda, work, 1, czero,
275 a( k+1, k+1 ) = a( k+1, k+1 ) -
276 $ dble( zdotc( k-1, work, 1, a( 1, k+1 ),
282 IF( kstep.EQ.1 )
THEN
291 $ CALL zswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
293 DO 40 j = kp + 1, k - 1
294 temp = dconjg( a( j, k ) )
295 a( j, k ) = dconjg( a( kp, j ) )
299 a( kp, k ) = dconjg( a( kp, k ) )
302 a( k, k ) = a( kp, kp )
316 $ CALL zswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
318 DO 50 j = kp + 1, k - 1
319 temp = dconjg( a( j, k ) )
320 a( j, k ) = dconjg( a( kp, j ) )
324 a( kp, k ) = dconjg( a( kp, k ) )
327 a( k, k ) = a( kp, kp )
331 a( k, k+1 ) = a( kp, k+1 )
342 $ CALL zswap( kp-1, a( 1, k ), 1, a( 1, kp ), 1 )
344 DO 60 j = kp + 1, k - 1
345 temp = dconjg( a( j, k ) )
346 a( j, k ) = dconjg( a( kp, j ) )
350 a( kp, k ) = dconjg( a( kp, k ) )
353 a( k, k ) = a( kp, kp )
377 IF( ipiv( k ).GT.0 )
THEN
383 a( k, k ) = one / dble( a( k, k ) )
388 CALL
zcopy( n-k, a( k+1, k ), 1, work, 1 )
389 CALL
zhemv( uplo, n-k, -cone, a( k+1, k+1 ), lda, work,
390 $ 1, czero, a( k+1, k ), 1 )
391 a( k, k ) = a( k, k ) - dble( zdotc( n-k, work, 1,
401 t = abs( a( k, k-1 ) )
402 ak = dble( a( k-1, k-1 ) ) / t
403 akp1 = dble( a( k, k ) ) / t
404 akkp1 = a( k, k-1 ) / t
405 d = t*( ak*akp1-one )
406 a( k-1, k-1 ) = akp1 / d
408 a( k, k-1 ) = -akkp1 / d
413 CALL
zcopy( n-k, a( k+1, k ), 1, work, 1 )
414 CALL
zhemv( uplo, n-k, -cone, a( k+1, k+1 ), lda, work,
415 $ 1, czero, a( k+1, k ), 1 )
416 a( k, k ) = a( k, k ) - dble( zdotc( n-k, work, 1,
418 a( k, k-1 ) = a( k, k-1 ) -
419 $ zdotc( n-k, a( k+1, k ), 1, a( k+1, k-1 ),
421 CALL
zcopy( n-k, a( k+1, k-1 ), 1, work, 1 )
422 CALL
zhemv( uplo, n-k, -cone, a( k+1, k+1 ), lda, work,
423 $ 1, czero, a( k+1, k-1 ), 1 )
424 a( k-1, k-1 ) = a( k-1, k-1 ) -
425 $ dble( zdotc( n-k, work, 1, a( k+1, k-1 ),
431 IF( kstep.EQ.1 )
THEN
440 $ CALL zswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
442 DO 90 j = k + 1, kp - 1
443 temp = dconjg( a( j, k ) )
444 a( j, k ) = dconjg( a( kp, j ) )
448 a( kp, k ) = dconjg( a( kp, k ) )
451 a( k, k ) = a( kp, kp )
465 $ CALL zswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
467 DO 100 j = k + 1, kp - 1
468 temp = dconjg( a( j, k ) )
469 a( j, k ) = dconjg( a( kp, j ) )
473 a( kp, k ) = dconjg( a( kp, k ) )
476 a( k, k ) = a( kp, kp )
480 a( k, k-1 ) = a( kp, k-1 )
491 $ CALL zswap( n-kp, a( kp+1, k ), 1, a( kp+1, kp ), 1 )
493 DO 110 j = k + 1, kp - 1
494 temp = dconjg( a( j, k ) )
495 a( j, k ) = dconjg( a( kp, j ) )
499 a( kp, k ) = dconjg( a( kp, k ) )
502 a( k, k ) = a( kp, kp )
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine zhetri_rook(UPLO, N, A, LDA, IPIV, WORK, INFO)
ZHETRI_ROOK computes the inverse of HE matrix using the factorization obtained with the bounded Bunch...
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
subroutine zhemv(UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZHEMV