242 SUBROUTINE chetf2_rk( UPLO, N, A, LDA, E, IPIV, INFO )
255 COMPLEX A( lda, * ), E( * )
262 parameter ( zero = 0.0e+0, one = 1.0e+0 )
264 parameter ( eight = 8.0e+0, sevten = 17.0e+0 )
266 parameter ( czero = ( 0.0e+0, 0.0e+0 ) )
270 INTEGER I, II, IMAX, ITEMP, J, JMAX, K, KK, KP, KSTEP,
272 REAL ABSAKK, ALPHA, COLMAX, D, D11, D22, R1, STEMP,
274 COMPLEX D12, D21, T, WK, WKM1, WKP1, Z
281 EXTERNAL lsame, icamax, slamch, slapy2
287 INTRINSIC abs, aimag, cmplx, conjg, max,
REAL, SQRT
293 cabs1( z ) = abs(
REAL( Z ) ) + abs( AIMAG( z ) )
300 upper = lsame( uplo,
'U' )
301 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
303 ELSE IF( n.LT.0 )
THEN
305 ELSE IF( lda.LT.max( 1, n ) )
THEN
309 CALL xerbla(
'CHETF2_RK', -info )
315 alpha = ( one+sqrt( sevten ) ) / eight
319 sfmin = slamch(
'S' )
346 absakk = abs(
REAL( A( K, K ) ) )
353 imax = icamax( k-1, a( 1, k ), 1 )
354 colmax = cabs1( a( imax, k ) )
359 IF( ( max( absakk, colmax ).EQ.zero ) )
THEN
366 a( k, k ) =
REAL( A( K, K ) )
383 IF( .NOT.( absakk.LT.alpha*colmax ) )
THEN
405 jmax = imax + icamax( k-imax, a( imax, imax+1 ),
407 rowmax = cabs1( a( imax, jmax ) )
413 itemp = icamax( imax-1, a( 1, imax ), 1 )
414 stemp = cabs1( a( itemp, imax ) )
415 IF( stemp.GT.rowmax )
THEN
426 IF( .NOT.( abs(
REAL( A( IMAX, IMAX ) ) )
427 $ .LT.alpha*rowmax ) )
THEN
439 ELSE IF( ( p.EQ.jmax ) .OR. ( rowmax.LE.colmax ) )
461 IF( .NOT.done )
GOTO 12
476 IF( ( kstep.EQ.2 ) .AND. ( p.NE.k ) )
THEN
479 $
CALL cswap( p-1, a( 1, k ), 1, a( 1, p ), 1 )
481 DO 14 j = p + 1, k - 1
482 t = conjg( a( j, k ) )
483 a( j, k ) = conjg( a( p, j ) )
487 a( p, k ) = conjg( a( p, k ) )
489 r1 =
REAL( A( K, K ) )
490 a( k, k ) =
REAL( A( P, P ) )
497 $
CALL cswap( n-k, a( k, k+1 ), lda, a( p, k+1 ), lda )
507 $
CALL cswap( kp-1, a( 1, kk ), 1, a( 1, kp ), 1 )
509 DO 15 j = kp + 1, kk - 1
510 t = conjg( a( j, kk ) )
511 a( j, kk ) = conjg( a( kp, j ) )
515 a( kp, kk ) = conjg( a( kp, kk ) )
517 r1 =
REAL( A( KK, KK ) )
518 a( kk, kk ) =
REAL( A( KP, KP ) )
521 IF( kstep.EQ.2 )
THEN
523 a( k, k ) =
REAL( A( K, K ) )
526 a( k-1, k ) = a( kp, k )
534 $
CALL cswap( n-k, a( kk, k+1 ), lda, a( kp, k+1 ),
539 a( k, k ) =
REAL( A( K, K ) )
541 $ a( k-1, k-1 ) =
REAL( A( K-1, K-1 ) )
546 IF( kstep.EQ.1 )
THEN
559 IF( abs(
REAL( A( K, K ) ) ).GE.sfmin ) then
565 d11 = one /
REAL( A( K, K ) )
566 CALL cher( uplo, k-1, -d11, a( 1, k ), 1, a, lda )
570 CALL csscal( k-1, d11, a( 1, k ), 1 )
575 d11 =
REAL( A( K, K ) )
577 a( ii, k ) = a( ii, k ) / d11
585 CALL cher( uplo, k-1, -d11, a( 1, k ), 1, a, lda )
612 d = slapy2(
REAL( A( K-1, K ) ),
613 $ aimag( a( k-1, k ) ) )
615 d22 = a( k-1, k-1 ) / d
616 d12 = a( k-1, k ) / d
617 tt = one / ( d11*d22-one )
619 DO 30 j = k - 2, 1, -1
623 wkm1 = tt*( d11*a( j, k-1 )-conjg( d12 )*
625 wk = tt*( d22*a( j, k )-d12*a( j, k-1 ) )
630 a( i, j ) = a( i, j ) -
631 $ ( a( i, k ) / d )*conjg( wk ) -
632 $ ( a( i, k-1 ) / d )*conjg( wkm1 )
638 a( j, k-1 ) = wkm1 / d
640 a( j, j ) = cmplx(
REAL( A( J, J ) ), ZERO )
661 IF( kstep.EQ.1 )
THEN
699 absakk = abs(
REAL( A( K, K ) ) )
706 imax = k + icamax( n-k, a( k+1, k ), 1 )
707 colmax = cabs1( a( imax, k ) )
712 IF( max( absakk, colmax ).EQ.zero )
THEN
719 a( k, k ) =
REAL( A( K, K ) )
736 IF( .NOT.( absakk.LT.alpha*colmax ) )
THEN
758 jmax = k - 1 + icamax( imax-k, a( imax, k ), lda )
759 rowmax = cabs1( a( imax, jmax ) )
765 itemp = imax + icamax( n-imax, a( imax+1, imax ),
767 stemp = cabs1( a( itemp, imax ) )
768 IF( stemp.GT.rowmax )
THEN
779 IF( .NOT.( abs(
REAL( A( IMAX, IMAX ) ) )
780 $ .LT.alpha*rowmax ) )
THEN
792 ELSE IF( ( p.EQ.jmax ) .OR. ( rowmax.LE.colmax ) )
815 IF( .NOT.done )
GOTO 42
830 IF( ( kstep.EQ.2 ) .AND. ( p.NE.k ) )
THEN
833 $
CALL cswap( n-p, a( p+1, k ), 1, a( p+1, p ), 1 )
835 DO 44 j = k + 1, p - 1
836 t = conjg( a( j, k ) )
837 a( j, k ) = conjg( a( p, j ) )
841 a( p, k ) = conjg( a( p, k ) )
843 r1 =
REAL( A( K, K ) )
844 a( k, k ) =
REAL( A( P, P ) )
851 $
CALL cswap( k-1, a( k, 1 ), lda, a( p, 1 ), lda )
861 $
CALL cswap( n-kp, a( kp+1, kk ), 1, a( kp+1, kp ), 1 )
863 DO 45 j = kk + 1, kp - 1
864 t = conjg( a( j, kk ) )
865 a( j, kk ) = conjg( a( kp, j ) )
869 a( kp, kk ) = conjg( a( kp, kk ) )
871 r1 =
REAL( A( KK, KK ) )
872 a( kk, kk ) =
REAL( A( KP, KP ) )
875 IF( kstep.EQ.2 )
THEN
877 a( k, k ) =
REAL( A( K, K ) )
880 a( k+1, k ) = a( kp, k )
888 $
CALL cswap( k-1, a( kk, 1 ), lda, a( kp, 1 ), lda )
892 a( k, k ) =
REAL( A( K, K ) )
894 $ a( k+1, k+1 ) =
REAL( A( K+1, K+1 ) )
899 IF( kstep.EQ.1 )
THEN
914 IF( abs(
REAL( A( K, K ) ) ).GE.sfmin ) then
920 d11 = one /
REAL( A( K, K ) )
921 CALL cher( uplo, n-k, -d11, a( k+1, k ), 1,
922 $ a( k+1, k+1 ), lda )
926 CALL csscal( n-k, d11, a( k+1, k ), 1 )
931 d11 =
REAL( A( K, K ) )
933 a( ii, k ) = a( ii, k ) / d11
941 CALL cher( uplo, n-k, -d11, a( k+1, k ), 1,
942 $ a( k+1, k+1 ), lda )
970 d = slapy2(
REAL( A( K+1, K ) ),
971 $ aimag( a( k+1, k ) ) )
972 d11 =
REAL( A( K+1, K+1 ) ) / D
973 d22 =
REAL( A( K, K ) ) / D
974 d21 = a( k+1, k ) / d
975 tt = one / ( d11*d22-one )
981 wk = tt*( d11*a( j, k )-d21*a( j, k+1 ) )
982 wkp1 = tt*( d22*a( j, k+1 )-conjg( d21 )*
988 a( i, j ) = a( i, j ) -
989 $ ( a( i, k ) / d )*conjg( wk ) -
990 $ ( a( i, k+1 ) / d )*conjg( wkp1 )
996 a( j, k+1 ) = wkp1 / d
998 a( j, j ) = cmplx(
REAL( A( J, J ) ), ZERO )
1007 e( k ) = a( k+1, k )
1019 IF( kstep.EQ.1 )
THEN
subroutine cher(UPLO, N, ALPHA, X, INCX, A, LDA)
CHER
subroutine chetf2_rk(UPLO, N, A, LDA, E, IPIV, INFO)
CHETF2_RK computes the factorization of a complex Hermitian indefinite matrix using the bounded Bunch...
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
subroutine csscal(N, SA, CX, INCX)
CSSCAL