178 SUBROUTINE clasyf( UPLO, N, NB, KB, A, LDA, IPIV, W, LDW, INFO )
187 INTEGER INFO, KB, LDA, LDW, N, NB
191 COMPLEX A( lda, * ), W( ldw, * )
198 parameter( zero = 0.0e+0, one = 1.0e+0 )
200 parameter( eight = 8.0e+0, sevten = 17.0e+0 )
202 parameter( cone = ( 1.0e+0, 0.0e+0 ) )
205 INTEGER IMAX, J, JB, JJ, JMAX, JP, K, KK, KKW, KP,
207 REAL ABSAKK, ALPHA, COLMAX, ROWMAX
208 COMPLEX D11, D21, D22, R1, T, Z
213 EXTERNAL lsame, icamax
219 INTRINSIC abs, aimag, max, min,
REAL, SQRT
225 cabs1( z ) = abs(
REAL( Z ) ) + abs( AIMAG( 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 )
256 CALL
ccopy( k, a( 1, k ), 1, w( 1, kw ), 1 )
258 $ CALL
cgemv(
'No transpose', k, n-k, -cone, a( 1, k+1 ), lda,
259 $ w( k, kw+1 ), ldw, cone, w( 1, kw ), 1 )
266 absakk = cabs1( w( k, kw ) )
273 imax = icamax( k-1, w( 1, kw ), 1 )
274 colmax = cabs1( w( imax, kw ) )
279 IF( max( absakk, colmax ).EQ.zero )
THEN
287 IF( absakk.GE.alpha*colmax )
THEN
296 CALL
ccopy( imax, a( 1, imax ), 1, w( 1, kw-1 ), 1 )
297 CALL
ccopy( k-imax, a( imax, imax+1 ), lda,
298 $ w( imax+1, kw-1 ), 1 )
300 $ CALL
cgemv(
'No transpose', k, n-k, -cone,
301 $ a( 1, k+1 ), lda, w( imax, kw+1 ), ldw,
302 $ cone, w( 1, kw-1 ), 1 )
307 jmax = imax + icamax( k-imax, w( imax+1, kw-1 ), 1 )
308 rowmax = cabs1( w( jmax, kw-1 ) )
310 jmax = icamax( imax-1, w( 1, kw-1 ), 1 )
311 rowmax = max( rowmax, cabs1( w( jmax, kw-1 ) ) )
314 IF( absakk.GE.alpha*colmax*( colmax / rowmax ) )
THEN
319 ELSE IF( cabs1( w( imax, kw-1 ) ).GE.alpha*rowmax )
THEN
328 CALL
ccopy( k, w( 1, kw-1 ), 1, w( 1, kw ), 1 )
359 a( kp, kp ) = a( kk, kk )
360 CALL
ccopy( kk-1-kp, a( kp+1, kk ), 1, a( kp, kp+1 ),
363 $ CALL
ccopy( kp-1, a( 1, kk ), 1, a( 1, kp ), 1 )
371 $ CALL
cswap( n-k, a( kk, k+1 ), lda, a( kp, k+1 ),
373 CALL
cswap( n-kk+1, w( kk, kkw ), ldw, w( kp, kkw ),
377 IF( kstep.EQ.1 )
THEN
392 CALL
ccopy( k, w( 1, kw ), 1, a( 1, k ), 1 )
393 r1 = cone / a( k, k )
394 CALL
cscal( k-1, r1, a( 1, k ), 1 )
441 d11 = w( k, kw ) / d21
442 d22 = w( k-1, kw-1 ) / d21
443 t = cone / ( d11*d22-cone )
451 a( j, k-1 ) = d21*( d11*w( j, kw-1 )-w( j, kw ) )
452 a( j, k ) = d21*( d22*w( j, kw )-w( j, kw-1 ) )
458 a( k-1, k-1 ) = w( k-1, kw-1 )
459 a( k-1, k ) = w( k-1, kw )
460 a( k, k ) = w( k, kw )
468 IF( kstep.EQ.1 )
THEN
488 DO 50 j = ( ( k-1 ) / nb )*nb + 1, 1, -nb
489 jb = min( nb, k-j+1 )
493 DO 40 jj = j, j + jb - 1
494 CALL
cgemv(
'No transpose', jj-j+1, n-k, -cone,
495 $ a( j, k+1 ), lda, w( jj, kw+1 ), ldw, cone,
501 CALL
cgemm(
'No transpose',
'Transpose', j-1, jb, n-k,
502 $ -cone, a( 1, k+1 ), lda, w( j, kw+1 ), ldw,
503 $ cone, a( 1, j ), lda )
526 IF( jp.NE.jj .AND. j.LE.n )
527 $ CALL
cswap( n-j+1, a( jp, j ), lda, a( jj, j ), lda )
548 IF( ( k.GE.nb .AND. nb.LT.n ) .OR. k.GT.n )
553 CALL
ccopy( n-k+1, a( k, k ), 1, w( k, k ), 1 )
554 CALL
cgemv(
'No transpose', n-k+1, k-1, -cone, a( k, 1 ), lda,
555 $ w( k, 1 ), ldw, cone, w( k, k ), 1 )
562 absakk = cabs1( w( k, k ) )
569 imax = k + icamax( n-k, w( k+1, k ), 1 )
570 colmax = cabs1( w( imax, k ) )
575 IF( max( absakk, colmax ).EQ.zero )
THEN
583 IF( absakk.GE.alpha*colmax )
THEN
592 CALL
ccopy( imax-k, a( imax, k ), lda, w( k, k+1 ), 1 )
593 CALL
ccopy( n-imax+1, a( imax, imax ), 1, w( imax, k+1 ),
595 CALL
cgemv(
'No transpose', n-k+1, k-1, -cone, a( k, 1 ),
596 $ lda, w( imax, 1 ), ldw, cone, w( k, k+1 ),
602 jmax = k - 1 + icamax( imax-k, w( k, k+1 ), 1 )
603 rowmax = cabs1( w( jmax, k+1 ) )
605 jmax = imax + icamax( n-imax, w( imax+1, k+1 ), 1 )
606 rowmax = max( rowmax, cabs1( w( jmax, k+1 ) ) )
609 IF( absakk.GE.alpha*colmax*( colmax / rowmax ) )
THEN
614 ELSE IF( cabs1( w( imax, k+1 ) ).GE.alpha*rowmax )
THEN
623 CALL
ccopy( n-k+1, w( k, k+1 ), 1, w( k, k ), 1 )
650 a( kp, kp ) = a( kk, kk )
651 CALL
ccopy( kp-kk-1, a( kk+1, kk ), 1, a( kp, kk+1 ),
654 $ CALL
ccopy( n-kp, a( kp+1, kk ), 1, a( kp+1, kp ), 1 )
662 $ CALL
cswap( k-1, a( kk, 1 ), lda, a( kp, 1 ), lda )
663 CALL
cswap( kk, w( kk, 1 ), ldw, w( kp, 1 ), ldw )
666 IF( kstep.EQ.1 )
THEN
681 CALL
ccopy( n-k+1, w( k, k ), 1, a( k, k ), 1 )
683 r1 = cone / a( k, k )
684 CALL
cscal( n-k, r1, a( k+1, k ), 1 )
732 d11 = w( k+1, k+1 ) / d21
733 d22 = w( k, k ) / d21
734 t = cone / ( d11*d22-cone )
742 a( j, k ) = d21*( d11*w( j, k )-w( j, k+1 ) )
743 a( j, k+1 ) = d21*( d22*w( j, k+1 )-w( j, k ) )
749 a( k, k ) = w( k, k )
750 a( k+1, k ) = w( k+1, k )
751 a( k+1, k+1 ) = w( k+1, k+1 )
759 IF( kstep.EQ.1 )
THEN
780 jb = min( nb, n-j+1 )
784 DO 100 jj = j, j + jb - 1
785 CALL
cgemv(
'No transpose', j+jb-jj, k-1, -cone,
786 $ a( jj, 1 ), lda, w( jj, 1 ), ldw, cone,
793 $ CALL
cgemm(
'No transpose',
'Transpose', n-j-jb+1, jb,
794 $ k-1, -cone, a( j+jb, 1 ), lda, w( j, 1 ),
795 $ ldw, cone, a( j+jb, j ), lda )
818 IF( jp.NE.jj .AND. j.GE.1 )
819 $ CALL
cswap( j, a( jp, 1 ), lda, a( jj, 1 ), lda )
subroutine cscal(N, CA, CX, INCX)
CSCAL
subroutine clasyf(UPLO, N, NB, KB, A, LDA, IPIV, W, LDW, INFO)
CLASYF computes a partial factorization of a complex symmetric matrix using the Bunch-Kaufman diagona...
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
subroutine cgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CGEMV
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM