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
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 )
LOGICAL function lsame(CA, CB)
LSAME
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
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
INTEGER function icamax(N, CX, INCX)
ICAMAX