173 SUBROUTINE slasyf( UPLO, N, NB, KB, A, LDA, IPIV, W, LDW,
182 INTEGER INFO, KB, LDA, LDW, N, NB
186 REAL A( LDA, * ), W( LDW, * )
193 parameter( zero = 0.0e+0, one = 1.0e+0 )
195 parameter( eight = 8.0e+0, sevten = 17.0e+0 )
198 INTEGER IMAX, J, JJ, JMAX, JP, K, KK, KKW, KP,
200 REAL ABSAKK, ALPHA, COLMAX, D11, D21, D22, R1,
206 EXTERNAL lsame, isamax
212 INTRINSIC abs, max, min, sqrt
220 alpha = ( one+sqrt( sevten ) ) / eight
222 IF( lsame( uplo,
'U' ) )
THEN
238 IF( ( k.LE.n-nb+1 .AND. nb.LT.n ) .OR. k.LT.1 )
243 CALL scopy( k, a( 1, k ), 1, w( 1, kw ), 1 )
245 $
CALL sgemv(
'No transpose', k, n-k, -one, a( 1, k+1 ),
247 $ w( k, kw+1 ), ldw, one, w( 1, kw ), 1 )
254 absakk = abs( w( k, kw ) )
261 imax = isamax( k-1, w( 1, kw ), 1 )
262 colmax = abs( w( imax, kw ) )
267 IF( max( absakk, colmax ).EQ.zero )
THEN
275 IF( absakk.GE.alpha*colmax )
THEN
284 CALL scopy( imax, a( 1, imax ), 1, w( 1, kw-1 ), 1 )
285 CALL scopy( k-imax, a( imax, imax+1 ), lda,
286 $ w( imax+1, kw-1 ), 1 )
288 $
CALL sgemv(
'No transpose', k, n-k, -one, a( 1,
290 $ lda, w( imax, kw+1 ), ldw, one,
296 jmax = imax + isamax( k-imax, w( imax+1, kw-1 ), 1 )
297 rowmax = abs( w( jmax, kw-1 ) )
299 jmax = isamax( imax-1, w( 1, kw-1 ), 1 )
300 rowmax = max( rowmax, abs( w( jmax, kw-1 ) ) )
303 IF( absakk.GE.alpha*colmax*( colmax / rowmax ) )
THEN
308 ELSE IF( abs( w( imax, kw-1 ) ).GE.alpha*rowmax )
THEN
317 CALL scopy( k, w( 1, kw-1 ), 1, w( 1, kw ), 1 )
348 a( kp, kp ) = a( kk, kk )
349 CALL scopy( kk-1-kp, a( kp+1, kk ), 1, a( kp, kp+1 ),
352 $
CALL scopy( kp-1, a( 1, kk ), 1, a( 1, kp ), 1 )
360 $
CALL sswap( n-k, a( kk, k+1 ), lda, a( kp, k+1 ),
362 CALL sswap( n-kk+1, w( kk, kkw ), ldw, w( kp, kkw ),
366 IF( kstep.EQ.1 )
THEN
381 CALL scopy( k, w( 1, kw ), 1, a( 1, k ), 1 )
383 CALL sscal( k-1, r1, a( 1, k ), 1 )
430 d11 = w( k, kw ) / d21
431 d22 = w( k-1, kw-1 ) / d21
432 t = one / ( d11*d22-one )
440 a( j, k-1 ) = d21*( d11*w( j, kw-1 )-w( j, kw ) )
441 a( j, k ) = d21*( d22*w( j, kw )-w( j, kw-1 ) )
447 a( k-1, k-1 ) = w( k-1, kw-1 )
448 a( k-1, k ) = w( k-1, kw )
449 a( k, k ) = w( k, kw )
457 IF( kstep.EQ.1 )
THEN
475 CALL sgemmtr(
'Upper',
'No transpose',
'Transpose', k, n-k,
476 $ -one, a( 1, k+1 ), lda, w( 1, kw+1 ), ldw,
477 $ one, a( 1, 1 ), lda )
499 IF( jp.NE.jj .AND. j.LE.n )
500 $
CALL sswap( n-j+1, a( jp, j ), lda, a( jj, j ), lda )
521 IF( ( k.GE.nb .AND. nb.LT.n ) .OR. k.GT.n )
526 CALL scopy( n-k+1, a( k, k ), 1, w( k, k ), 1 )
527 CALL sgemv(
'No transpose', n-k+1, k-1, -one, a( k, 1 ),
529 $ w( k, 1 ), ldw, one, w( k, k ), 1 )
536 absakk = abs( w( k, k ) )
543 imax = k + isamax( n-k, w( k+1, k ), 1 )
544 colmax = abs( w( imax, k ) )
549 IF( max( absakk, colmax ).EQ.zero )
THEN
557 IF( absakk.GE.alpha*colmax )
THEN
566 CALL scopy( imax-k, a( imax, k ), lda, w( k, k+1 ),
568 CALL scopy( n-imax+1, a( imax, imax ), 1, w( imax,
571 CALL sgemv(
'No transpose', n-k+1, k-1, -one, a( k,
573 $ lda, w( imax, 1 ), ldw, one, w( k, k+1 ), 1 )
578 jmax = k - 1 + isamax( imax-k, w( k, k+1 ), 1 )
579 rowmax = abs( w( jmax, k+1 ) )
581 jmax = imax + isamax( n-imax, w( imax+1, k+1 ), 1 )
582 rowmax = max( rowmax, abs( w( jmax, k+1 ) ) )
585 IF( absakk.GE.alpha*colmax*( colmax / rowmax ) )
THEN
590 ELSE IF( abs( w( imax, k+1 ) ).GE.alpha*rowmax )
THEN
599 CALL scopy( n-k+1, w( k, k+1 ), 1, w( k, k ), 1 )
626 a( kp, kp ) = a( kk, kk )
627 CALL scopy( kp-kk-1, a( kk+1, kk ), 1, a( kp, kk+1 ),
630 $
CALL scopy( n-kp, a( kp+1, kk ), 1, a( kp+1, kp ),
639 $
CALL sswap( k-1, a( kk, 1 ), lda, a( kp, 1 ), lda )
640 CALL sswap( kk, w( kk, 1 ), ldw, w( kp, 1 ), ldw )
643 IF( kstep.EQ.1 )
THEN
658 CALL scopy( n-k+1, w( k, k ), 1, a( k, k ), 1 )
661 CALL sscal( n-k, r1, a( k+1, k ), 1 )
709 d11 = w( k+1, k+1 ) / d21
710 d22 = w( k, k ) / d21
711 t = one / ( d11*d22-one )
719 a( j, k ) = d21*( d11*w( j, k )-w( j, k+1 ) )
720 a( j, k+1 ) = d21*( d22*w( j, k+1 )-w( j, k ) )
726 a( k, k ) = w( k, k )
727 a( k+1, k ) = w( k+1, k )
728 a( k+1, k+1 ) = w( k+1, k+1 )
736 IF( kstep.EQ.1 )
THEN
754 CALL sgemmtr(
'Lower',
'No transpose',
'Transpose', n-k+1,
755 $ k-1, -one, a( k, 1 ), lda, w( k, 1 ), ldw,
756 $ one, a( k, k ), lda )
778 IF( jp.NE.jj .AND. j.GE.1 )
779 $
CALL sswap( j, a( jp, 1 ), lda, a( jj, 1 ), lda )