258 SUBROUTINE clasyf_rk( UPLO, N, NB, KB, A, LDA, E, IPIV, W, LDW,
267 INTEGER INFO, KB, LDA, LDW, N, NB
271 COMPLEX A( LDA, * ), E( * ), W( LDW, * )
278 parameter( zero = 0.0e+0, one = 1.0e+0 )
280 parameter( eight = 8.0e+0, sevten = 17.0e+0 )
282 parameter( cone = ( 1.0e+0, 0.0e+0 ),
283 $ czero = ( 0.0e+0, 0.0e+0 ) )
287 INTEGER IMAX, ITEMP, J, JB, JJ, JMAX, K, KK, KW, KKW,
289 REAL ABSAKK, ALPHA, COLMAX, ROWMAX, SFMIN, STEMP
290 COMPLEX D11, D12, D21, D22, R1, T, Z
296 EXTERNAL lsame, icamax, slamch
302 INTRINSIC abs, aimag, max, min, real, sqrt
308 cabs1( z ) = abs( real( z ) ) + abs( aimag( z ) )
316 alpha = ( one+sqrt( sevten ) ) / eight
320 sfmin = slamch(
'S' )
322 IF( lsame( uplo,
'U' ) )
THEN
344 IF( ( k.LE.n-nb+1 .AND. nb.LT.n ) .OR. k.LT.1 )
352 CALL ccopy( k, a( 1, k ), 1, w( 1, kw ), 1 )
354 $
CALL cgemv(
'No transpose', k, n-k, -cone, a( 1, k+1 ),
355 $ lda, w( k, kw+1 ), ldw, cone, w( 1, kw ), 1 )
360 absakk = cabs1( w( k, kw ) )
367 imax = icamax( k-1, w( 1, kw ), 1 )
368 colmax = cabs1( w( imax, kw ) )
373 IF( max( absakk, colmax ).EQ.zero )
THEN
380 CALL ccopy( k, w( 1, kw ), 1, a( 1, k ), 1 )
396 IF( .NOT.( absakk.LT.alpha*colmax ) )
THEN
415 CALL ccopy( imax, a( 1, imax ), 1, w( 1, kw-1 ),
417 CALL ccopy( k-imax, a( imax, imax+1 ), lda,
418 $ w( imax+1, kw-1 ), 1 )
421 $
CALL cgemv(
'No transpose', k, n-k, -cone,
422 $ a( 1, k+1 ), lda, w( imax, kw+1 ), ldw,
423 $ cone, w( 1, kw-1 ), 1 )
430 jmax = imax + icamax( k-imax, w( imax+1, kw-1 ),
432 rowmax = cabs1( w( jmax, kw-1 ) )
438 itemp = icamax( imax-1, w( 1, kw-1 ), 1 )
439 stemp = cabs1( w( itemp, kw-1 ) )
440 IF( stemp.GT.rowmax )
THEN
450 IF( .NOT.(cabs1( w( imax, kw-1 ) ).LT.alpha*rowmax ) )
460 CALL ccopy( k, w( 1, kw-1 ), 1, w( 1, kw ), 1 )
467 ELSE IF( ( p.EQ.jmax ) .OR. ( rowmax.LE.colmax ) )
486 CALL ccopy( k, w( 1, kw-1 ), 1, w( 1, kw ), 1 )
492 IF( .NOT. done )
GOTO 12
504 IF( ( kstep.EQ.2 ) .AND. ( p.NE.k ) )
THEN
508 CALL ccopy( k-p, a( p+1, k ), 1, a( p, p+1 ), lda )
509 CALL ccopy( p, a( 1, k ), 1, a( 1, p ), 1 )
514 CALL cswap( n-k+1, a( k, k ), lda, a( p, k ), lda )
515 CALL cswap( n-kk+1, w( k, kkw ), ldw, w( p, kkw ),
525 a( kp, k ) = a( kk, k )
526 CALL ccopy( k-1-kp, a( kp+1, kk ), 1, a( kp, kp+1 ),
528 CALL ccopy( kp, a( 1, kk ), 1, a( 1, kp ), 1 )
533 CALL cswap( n-kk+1, a( kk, kk ), lda, a( kp, kk ),
535 CALL cswap( n-kk+1, w( kk, kkw ), ldw, w( kp, kkw ),
539 IF( kstep.EQ.1 )
THEN
549 CALL ccopy( k, w( 1, kw ), 1, a( 1, k ), 1 )
551 IF( cabs1( a( k, k ) ).GE.sfmin )
THEN
552 r1 = cone / a( k, k )
553 CALL cscal( k-1, r1, a( 1, k ), 1 )
554 ELSE IF( a( k, k ).NE.czero )
THEN
556 a( ii, k ) = a( ii, k ) / a( k, k )
581 d11 = w( k, kw ) / d12
582 d22 = w( k-1, kw-1 ) / d12
583 t = cone / ( d11*d22-cone )
585 a( j, k-1 ) = t*( (d11*w( j, kw-1 )-w( j, kw ) ) /
587 a( j, k ) = t*( ( d22*w( j, kw )-w( j, kw-1 ) ) /
596 a( k-1, k-1 ) = w( k-1, kw-1 )
598 a( k, k ) = w( k, kw )
599 e( k ) = w( k-1, kw )
610 IF( kstep.EQ.1 )
THEN
628 CALL cgemmtr(
'Upper',
'No transpose',
'Transpose', k, n-k,
629 $ -cone, a( 1, k+1 ), lda, w( 1, kw+1 ), ldw,
630 $ cone, a( 1, 1 ), lda )
653 IF( ( k.GE.nb .AND. nb.LT.n ) .OR. k.GT.n )
661 CALL ccopy( n-k+1, a( k, k ), 1, w( k, k ), 1 )
663 $
CALL cgemv(
'No transpose', n-k+1, k-1, -cone, a( k, 1 ),
664 $ lda, w( k, 1 ), ldw, cone, w( k, k ), 1 )
669 absakk = cabs1( w( k, k ) )
676 imax = k + icamax( n-k, w( k+1, k ), 1 )
677 colmax = cabs1( w( imax, k ) )
682 IF( max( absakk, colmax ).EQ.zero )
THEN
689 CALL ccopy( n-k+1, w( k, k ), 1, a( k, k ), 1 )
705 IF( .NOT.( absakk.LT.alpha*colmax ) )
THEN
724 CALL ccopy( imax-k, a( imax, k ), lda, w( k, k+1 ),
726 CALL ccopy( n-imax+1, a( imax, imax ), 1,
727 $ w( imax, k+1 ), 1 )
729 $
CALL cgemv(
'No transpose', n-k+1, k-1, -cone,
730 $ a( k, 1 ), lda, w( imax, 1 ), ldw,
731 $ cone, w( k, k+1 ), 1 )
738 jmax = k - 1 + icamax( imax-k, w( k, k+1 ), 1 )
739 rowmax = cabs1( w( jmax, k+1 ) )
745 itemp = imax + icamax( n-imax, w( imax+1, k+1 ),
747 stemp = cabs1( w( itemp, k+1 ) )
748 IF( stemp.GT.rowmax )
THEN
758 IF( .NOT.( cabs1( w( imax, k+1 ) ).LT.alpha*rowmax ) )
768 CALL ccopy( n-k+1, w( k, k+1 ), 1, w( k, k ),
776 ELSE IF( ( p.EQ.jmax ) .OR. ( rowmax.LE.colmax ) )
795 CALL ccopy( n-k+1, w( k, k+1 ), 1, w( k, k ),
802 IF( .NOT. done )
GOTO 72
810 IF( ( kstep.EQ.2 ) .AND. ( p.NE.k ) )
THEN
814 CALL ccopy( p-k, a( k, k ), 1, a( p, k ), lda )
815 CALL ccopy( n-p+1, a( p, k ), 1, a( p, p ), 1 )
820 CALL cswap( k, a( k, 1 ), lda, a( p, 1 ), lda )
821 CALL cswap( kk, w( k, 1 ), ldw, w( p, 1 ), ldw )
830 a( kp, k ) = a( kk, k )
831 CALL ccopy( kp-k-1, a( k+1, kk ), 1, a( kp, k+1 ),
833 CALL ccopy( n-kp+1, a( kp, kk ), 1, a( kp, kp ), 1 )
837 CALL cswap( kk, a( kk, 1 ), lda, a( kp, 1 ), lda )
838 CALL cswap( kk, w( kk, 1 ), ldw, w( kp, 1 ), ldw )
841 IF( kstep.EQ.1 )
THEN
851 CALL ccopy( n-k+1, w( k, k ), 1, a( k, k ), 1 )
853 IF( cabs1( a( k, k ) ).GE.sfmin )
THEN
854 r1 = cone / a( k, k )
855 CALL cscal( n-k, r1, a( k+1, k ), 1 )
856 ELSE IF( a( k, k ).NE.czero )
THEN
858 a( ii, k ) = a( ii, k ) / a( k, k )
882 d11 = w( k+1, k+1 ) / d21
883 d22 = w( k, k ) / d21
884 t = cone / ( d11*d22-cone )
886 a( j, k ) = t*( ( d11*w( j, k )-w( j, k+1 ) ) /
888 a( j, k+1 ) = t*( ( d22*w( j, k+1 )-w( j, k ) ) /
897 a( k, k ) = w( k, k )
899 a( k+1, k+1 ) = w( k+1, k+1 )
911 IF( kstep.EQ.1 )
THEN
929 CALL cgemmtr(
'Lower',
'No transpose',
'Transpose', n-k+1,
930 $ k-1, -cone, a( k, 1 ), lda, w( k, 1 ), ldw,
931 $ cone, a( k, k ), lda )