178 SUBROUTINE zlahef( UPLO, N, NB, KB, A, LDA, IPIV, W, LDW, INFO )
187 INTEGER info, kb, lda, ldw, n, nb
191 COMPLEX*16 a( lda, * ), w( ldw, * )
197 DOUBLE PRECISION zero, one
198 parameter( zero = 0.0d+0, one = 1.0d+0 )
200 parameter( cone = ( 1.0d+0, 0.0d+0 ) )
201 DOUBLE PRECISION eight, sevten
202 parameter( eight = 8.0d+0, sevten = 17.0d+0 )
205 INTEGER imax,
j, jb, jj, jmax, jp, k, kk, kkw, kp,
207 DOUBLE PRECISION absakk, alpha, colmax, r1, rowmax, t
208 COMPLEX*16 d11, d21, d22, z
219 INTRINSIC abs, dble, dconjg, dimag, max, min, sqrt
222 DOUBLE PRECISION cabs1
225 cabs1( z ) = abs( dble( z ) ) + abs( dimag( z ) )
233 alpha = ( one+sqrt( sevten ) ) / eight
235 IF(
lsame( uplo,
'U' ) )
THEN
251 IF( ( k.LE.n-nb+1 .AND. nb.LT.n ) .OR. k.LT.1 )
258 CALL
zcopy( k-1, a( 1, k ), 1, w( 1, kw ), 1 )
259 w( k, kw ) = dble( a( k, k ) )
261 CALL
zgemv(
'No transpose', k, n-k, -cone, a( 1, k+1 ), lda,
262 $ w( k, kw+1 ), ldw, cone, w( 1, kw ), 1 )
263 w( k, kw ) = dble( w( k, kw ) )
269 absakk = abs( dble( w( k, kw ) ) )
276 imax =
izamax( k-1, w( 1, kw ), 1 )
277 colmax = cabs1( w( imax, kw ) )
282 IF( max( absakk, colmax ).EQ.zero )
THEN
289 a( k, k ) = dble( a( k, k ) )
297 IF( absakk.GE.alpha*colmax )
THEN
309 CALL
zcopy( imax-1, a( 1, imax ), 1, w( 1, kw-1 ), 1 )
310 w( imax, kw-1 ) = dble( a( imax, imax ) )
311 CALL
zcopy( k-imax, a( imax, imax+1 ), lda,
312 $ w( imax+1, kw-1 ), 1 )
313 CALL
zlacgv( k-imax, w( imax+1, kw-1 ), 1 )
315 CALL
zgemv(
'No transpose', k, n-k, -cone,
316 $ a( 1, k+1 ), lda, w( imax, kw+1 ), ldw,
317 $ cone, w( 1, kw-1 ), 1 )
318 w( imax, kw-1 ) = dble( w( imax, kw-1 ) )
325 jmax = imax +
izamax( k-imax, w( imax+1, kw-1 ), 1 )
326 rowmax = cabs1( w( jmax, kw-1 ) )
328 jmax =
izamax( imax-1, w( 1, kw-1 ), 1 )
329 rowmax = max( rowmax, cabs1( w( jmax, kw-1 ) ) )
333 IF( absakk.GE.alpha*colmax*( colmax / rowmax ) )
THEN
340 ELSE IF( abs( dble( w( imax, kw-1 ) ) ).GE.alpha*rowmax )
350 CALL
zcopy( k, w( 1, kw-1 ), 1, w( 1, kw ), 1 )
389 a( kp, kp ) = dble( a( kk, kk ) )
390 CALL
zcopy( kk-1-kp, a( kp+1, kk ), 1, a( kp, kp+1 ),
392 CALL
zlacgv( kk-1-kp, a( kp, kp+1 ), lda )
394 $ CALL
zcopy( kp-1, a( 1, kk ), 1, a( 1, kp ), 1 )
402 $ CALL
zswap( n-k, a( kk, k+1 ), lda, a( kp, k+1 ),
404 CALL
zswap( n-kk+1, w( kk, kkw ), ldw, w( kp, kkw ),
408 IF( kstep.EQ.1 )
THEN
426 CALL
zcopy( k, w( 1, kw ), 1, a( 1, k ), 1 )
433 r1 = one / dble( a( k, k ) )
434 CALL
zdscal( k-1, r1, a( 1, k ), 1 )
438 CALL
zlacgv( k-1, w( 1, kw ), 1 )
503 d11 = w( k, kw ) / dconjg( d21 )
504 d22 = w( k-1, kw-1 ) / d21
505 t = one / ( dble( d11*d22 )-one )
513 a(
j, k-1 ) = d21*( d11*w(
j, kw-1 )-w(
j, kw ) )
514 a(
j, k ) = dconjg( d21 )*
515 $ ( d22*w(
j, kw )-w(
j, kw-1 ) )
521 a( k-1, k-1 ) = w( k-1, kw-1 )
522 a( k-1, k ) = w( k-1, kw )
523 a( k, k ) = w( k, kw )
527 CALL
zlacgv( k-1, w( 1, kw ), 1 )
528 CALL
zlacgv( k-2, w( 1, kw-1 ), 1 )
536 IF( kstep.EQ.1 )
THEN
557 DO 50
j = ( ( k-1 ) / nb )*nb + 1, 1, -nb
558 jb = min( nb, k-
j+1 )
562 DO 40 jj =
j,
j + jb - 1
563 a( jj, jj ) = dble( a( jj, jj ) )
564 CALL
zgemv(
'No transpose', jj-
j+1, n-k, -cone,
565 $ a(
j, k+1 ), lda, w( jj, kw+1 ), ldw, cone,
567 a( jj, jj ) = dble( a( jj, jj ) )
572 CALL
zgemm(
'No transpose',
'Transpose',
j-1, jb, n-k,
573 $ -cone, a( 1, k+1 ), lda, w(
j, kw+1 ), ldw,
574 $ cone, a( 1,
j ), lda )
597 IF( jp.NE.jj .AND.
j.LE.n )
598 $ CALL
zswap( n-
j+1, a( jp,
j ), lda, a( jj,
j ), lda )
619 IF( ( k.GE.nb .AND. nb.LT.n ) .OR. k.GT.n )
626 w( k, k ) = dble( a( k, k ) )
628 $ CALL
zcopy( n-k, a( k+1, k ), 1, w( k+1, k ), 1 )
629 CALL
zgemv(
'No transpose', n-k+1, k-1, -cone, a( k, 1 ), lda,
630 $ w( k, 1 ), ldw, cone, w( k, k ), 1 )
631 w( k, k ) = dble( w( k, k ) )
636 absakk = abs( dble( w( k, k ) ) )
643 imax = k +
izamax( n-k, w( k+1, k ), 1 )
644 colmax = cabs1( w( imax, k ) )
649 IF( max( absakk, colmax ).EQ.zero )
THEN
656 a( k, k ) = dble( a( k, k ) )
664 IF( absakk.GE.alpha*colmax )
THEN
676 CALL
zcopy( imax-k, a( imax, k ), lda, w( k, k+1 ), 1 )
677 CALL
zlacgv( imax-k, w( k, k+1 ), 1 )
678 w( imax, k+1 ) = dble( a( imax, imax ) )
680 $ CALL
zcopy( n-imax, a( imax+1, imax ), 1,
681 $ w( imax+1, k+1 ), 1 )
682 CALL
zgemv(
'No transpose', n-k+1, k-1, -cone, a( k, 1 ),
683 $ lda, w( imax, 1 ), ldw, cone, w( k, k+1 ),
685 w( imax, k+1 ) = dble( w( imax, k+1 ) )
691 jmax = k - 1 +
izamax( imax-k, w( k, k+1 ), 1 )
692 rowmax = cabs1( w( jmax, k+1 ) )
694 jmax = imax +
izamax( n-imax, w( imax+1, k+1 ), 1 )
695 rowmax = max( rowmax, cabs1( w( jmax, k+1 ) ) )
699 IF( absakk.GE.alpha*colmax*( colmax / rowmax ) )
THEN
706 ELSE IF( abs( dble( w( imax, k+1 ) ) ).GE.alpha*rowmax )
716 CALL
zcopy( n-k+1, w( k, k+1 ), 1, w( k, k ), 1 )
751 a( kp, kp ) = dble( a( kk, kk ) )
752 CALL
zcopy( kp-kk-1, a( kk+1, kk ), 1, a( kp, kk+1 ),
754 CALL
zlacgv( kp-kk-1, a( kp, kk+1 ), lda )
756 $ CALL
zcopy( n-kp, a( kp+1, kk ), 1, a( kp+1, kp ), 1 )
764 $ CALL
zswap( k-1, a( kk, 1 ), lda, a( kp, 1 ), lda )
765 CALL
zswap( kk, w( kk, 1 ), ldw, w( kp, 1 ), ldw )
768 IF( kstep.EQ.1 )
THEN
786 CALL
zcopy( n-k+1, w( k, k ), 1, a( k, k ), 1 )
793 r1 = one / dble( a( k, k ) )
794 CALL
zdscal( n-k, r1, a( k+1, k ), 1 )
798 CALL
zlacgv( n-k, w( k+1, k ), 1 )
863 d11 = w( k+1, k+1 ) / d21
864 d22 = w( k, k ) / dconjg( d21 )
865 t = one / ( dble( d11*d22 )-one )
873 a(
j, k ) = dconjg( d21 )*
874 $ ( d11*w(
j, k )-w(
j, k+1 ) )
875 a(
j, k+1 ) = d21*( d22*w(
j, k+1 )-w(
j, k ) )
881 a( k, k ) = w( k, k )
882 a( k+1, k ) = w( k+1, k )
883 a( k+1, k+1 ) = w( k+1, k+1 )
887 CALL
zlacgv( n-k, w( k+1, k ), 1 )
888 CALL
zlacgv( n-k-1, w( k+2, k+1 ), 1 )
896 IF( kstep.EQ.1 )
THEN
918 jb = min( nb, n-
j+1 )
922 DO 100 jj =
j,
j + jb - 1
923 a( jj, jj ) = dble( a( jj, jj ) )
924 CALL
zgemv(
'No transpose',
j+jb-jj, k-1, -cone,
925 $ a( jj, 1 ), lda, w( jj, 1 ), ldw, cone,
927 a( jj, jj ) = dble( a( jj, jj ) )
933 $ CALL
zgemm(
'No transpose',
'Transpose', n-
j-jb+1, jb,
934 $ k-1, -cone, a(
j+jb, 1 ), lda, w(
j, 1 ),
935 $ ldw, cone, a(
j+jb,
j ), lda )
958 IF( jp.NE.jj .AND.
j.GE.1 )
959 $ CALL
zswap(
j, a( jp, 1 ), lda, a( jj, 1 ), lda )
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
LOGICAL function lsame(CA, CB)
LSAME
subroutine zgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZGEMV
subroutine zlahef(UPLO, N, NB, KB, A, LDA, IPIV, W, LDW, INFO)
ZLAHEF computes a partial factorization of a complex Hermitian indefinite matrix using the Bunch-Kauf...
subroutine zdscal(N, DA, ZX, INCX)
ZDSCAL
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
subroutine zlacgv(N, X, INCX)
ZLACGV conjugates a complex vector.
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
INTEGER function izamax(N, ZX, INCX)
IZAMAX