208 COMPLEX*16 a( lda, * )
214 DOUBLE PRECISION zero, one
215 parameter( zero = 0.0d+0, one = 1.0d+0 )
216 DOUBLE PRECISION eight, sevten
217 parameter( eight = 8.0d+0, sevten = 17.0d+0 )
221 INTEGER i, ii, imax, itemp,
j, jmax, k, kk, kp, kstep,
223 DOUBLE PRECISION absakk, alpha, colmax, d, d11, d22, r1, dtemp,
225 COMPLEX*16 d12, d21, t, wk, wkm1, wkp1, z
238 INTRINSIC abs, dble, dcmplx, dconjg, dimag, max, sqrt
241 DOUBLE PRECISION cabs1
244 cabs1( z ) = abs( dble( z ) ) + abs( dimag( z ) )
251 upper =
lsame( uplo,
'U' )
252 IF( .NOT.upper .AND. .NOT.
lsame( uplo,
'L' ) )
THEN
254 ELSE IF( n.LT.0 )
THEN
256 ELSE IF( lda.LT.max( 1, n ) )
THEN
260 CALL
xerbla(
'ZHETF2_ROOK', -info )
266 alpha = ( one+sqrt( sevten ) ) / eight
292 absakk = abs( dble( a( k, k ) ) )
299 imax =
izamax( k-1, a( 1, k ), 1 )
300 colmax = cabs1( a( imax, k ) )
305 IF( ( max( absakk, colmax ).EQ.zero ) )
THEN
312 a( k, k ) = dble( a( k, k ) )
323 IF( .NOT.( absakk.LT.alpha*colmax ) )
THEN
345 jmax = imax +
izamax( k-imax, a( imax, imax+1 ),
347 rowmax = cabs1( a( imax, jmax ) )
353 itemp =
izamax( imax-1, a( 1, imax ), 1 )
354 dtemp = cabs1( a( itemp, imax ) )
355 IF( dtemp.GT.rowmax )
THEN
366 IF( .NOT.( abs( dble( a( imax, imax ) ) )
367 $ .LT.alpha*rowmax ) )
THEN
379 ELSE IF( ( p.EQ.jmax ) .OR. ( rowmax.LE.colmax ) )
401 IF( .NOT.done ) goto 12
416 IF( ( kstep.EQ.2 ) .AND. ( p.NE.k ) )
THEN
419 $ CALL
zswap( p-1, a( 1, k ), 1, a( 1, p ), 1 )
421 DO 14
j = p + 1, k - 1
422 t = dconjg( a(
j, k ) )
423 a(
j, k ) = dconjg( a( p,
j ) )
427 a( p, k ) = dconjg( a( p, k ) )
429 r1 = dble( a( k, k ) )
430 a( k, k ) = dble( a( p, p ) )
440 $ CALL
zswap( kp-1, a( 1, kk ), 1, a( 1, kp ), 1 )
442 DO 15
j = kp + 1, kk - 1
443 t = dconjg( a(
j, kk ) )
444 a(
j, kk ) = dconjg( a( kp,
j ) )
448 a( kp, kk ) = dconjg( a( kp, kk ) )
450 r1 = dble( a( kk, kk ) )
451 a( kk, kk ) = dble( a( kp, kp ) )
454 IF( kstep.EQ.2 )
THEN
456 a( k, k ) = dble( a( k, k ) )
459 a( k-1, k ) = a( kp, k )
464 a( k, k ) = dble( a( k, k ) )
466 $ a( k-1, k-1 ) = dble( a( k-1, k-1 ) )
471 IF( kstep.EQ.1 )
THEN
484 IF( abs( dble( a( k, k ) ) ).GE.sfmin )
THEN
490 d11 = one / dble( a( k, k ) )
491 CALL
zher( uplo, k-1, -d11, a( 1, k ), 1, a, lda )
495 CALL
zdscal( k-1, d11, a( 1, k ), 1 )
500 d11 = dble( a( k, k ) )
502 a( ii, k ) = a( ii, k ) / d11
510 CALL
zher( uplo, k-1, -d11, a( 1, k ), 1, a, lda )
532 d =
dlapy2( dble( a( k-1, k ) ),
533 $ dimag( a( k-1, k ) ) )
535 d22 = a( k-1, k-1 ) / d
536 d12 = a( k-1, k ) / d
537 tt = one / ( d11*d22-one )
539 DO 30
j = k - 2, 1, -1
543 wkm1 = tt*( d11*a(
j, k-1 )-dconjg( d12 )*
545 wk = tt*( d22*a(
j, k )-d12*a(
j, k-1 ) )
550 a( i,
j ) = a( i,
j ) -
551 $ ( a( i, k ) / d )*dconjg( wk ) -
552 $ ( a( i, k-1 ) / d )*dconjg( wkm1 )
558 a(
j, k-1 ) = wkm1 / d
560 a(
j,
j ) = dcmplx( dble( a(
j,
j ) ), zero )
572 IF( kstep.EQ.1 )
THEN
604 absakk = abs( dble( a( k, k ) ) )
611 imax = k +
izamax( n-k, a( k+1, k ), 1 )
612 colmax = cabs1( a( imax, k ) )
617 IF( max( absakk, colmax ).EQ.zero )
THEN
624 a( k, k ) = dble( a( k, k ) )
635 IF( .NOT.( absakk.LT.alpha*colmax ) )
THEN
657 jmax = k - 1 +
izamax( imax-k, a( imax, k ), lda )
658 rowmax = cabs1( a( imax, jmax ) )
664 itemp = imax +
izamax( n-imax, a( imax+1, imax ),
666 dtemp = cabs1( a( itemp, imax ) )
667 IF( dtemp.GT.rowmax )
THEN
678 IF( .NOT.( abs( dble( a( imax, imax ) ) )
679 $ .LT.alpha*rowmax ) )
THEN
691 ELSE IF( ( p.EQ.jmax ) .OR. ( rowmax.LE.colmax ) )
714 IF( .NOT.done ) goto 42
729 IF( ( kstep.EQ.2 ) .AND. ( p.NE.k ) )
THEN
732 $ CALL
zswap( n-p, a( p+1, k ), 1, a( p+1, p ), 1 )
734 DO 44
j = k + 1, p - 1
735 t = dconjg( a(
j, k ) )
736 a(
j, k ) = dconjg( a( p,
j ) )
740 a( p, k ) = dconjg( a( p, k ) )
742 r1 = dble( a( k, k ) )
743 a( k, k ) = dble( a( p, p ) )
753 $ CALL
zswap( n-kp, a( kp+1, kk ), 1, a( kp+1, kp ), 1 )
755 DO 45
j = kk + 1, kp - 1
756 t = dconjg( a(
j, kk ) )
757 a(
j, kk ) = dconjg( a( kp,
j ) )
761 a( kp, kk ) = dconjg( a( kp, kk ) )
763 r1 = dble( a( kk, kk ) )
764 a( kk, kk ) = dble( a( kp, kp ) )
767 IF( kstep.EQ.2 )
THEN
769 a( k, k ) = dble( a( k, k ) )
772 a( k+1, k ) = a( kp, k )
777 a( k, k ) = dble( a( k, k ) )
779 $ a( k+1, k+1 ) = dble( a( k+1, k+1 ) )
784 IF( kstep.EQ.1 )
THEN
799 IF( abs( dble( a( k, k ) ) ).GE.sfmin )
THEN
805 d11 = one / dble( a( k, k ) )
806 CALL
zher( uplo, n-k, -d11, a( k+1, k ), 1,
807 $ a( k+1, k+1 ), lda )
811 CALL
zdscal( n-k, d11, a( k+1, k ), 1 )
816 d11 = dble( a( k, k ) )
818 a( ii, k ) = a( ii, k ) / d11
826 CALL
zher( uplo, n-k, -d11, a( k+1, k ), 1,
827 $ a( k+1, k+1 ), lda )
850 d =
dlapy2( dble( a( k+1, k ) ),
851 $ dimag( a( k+1, k ) ) )
852 d11 = dble( a( k+1, k+1 ) ) / d
853 d22 = dble( a( k, k ) ) / d
854 d21 = a( k+1, k ) / d
855 tt = one / ( d11*d22-one )
861 wk = tt*( d11*a(
j, k )-d21*a(
j, k+1 ) )
862 wkp1 = tt*( d22*a(
j, k+1 )-dconjg( d21 )*
868 a( i,
j ) = a( i,
j ) -
869 $ ( a( i, k ) / d )*dconjg( wk ) -
870 $ ( a( i, k+1 ) / d )*dconjg( wkp1 )
876 a(
j, k+1 ) = wkp1 / d
878 a(
j,
j ) = dcmplx( dble( a(
j,
j ) ), zero )
890 IF( kstep.EQ.1 )
THEN
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
subroutine zhetf2_rook(UPLO, N, A, LDA, IPIV, INFO)
ZHETF2_ROOK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bun...
integer function izamax(N, ZX, INCX)
IZAMAX
subroutine zdscal(N, DA, ZX, INCX)
ZDSCAL
subroutine xerbla(SRNAME, INFO)
XERBLA
logical function lsame(CA, CB)
LSAME
double precision function dlapy2(X, Y)
DLAPY2 returns sqrt(x2+y2).
double precision function dlamch(CMACH)
DLAMCH
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine zher(UPLO, N, ALPHA, X, INCX, A, LDA)
ZHER