236 SUBROUTINE slatrs( UPLO, TRANS, DIAG, NORMIN, N, A, LDA, X, SCALE,
244 CHARACTER DIAG, NORMIN, TRANS, UPLO
249 REAL A( LDA, * ), CNORM( * ), X( * )
256 parameter( zero = 0.0e+0, half = 0.5e+0, one = 1.0e+0 )
259 LOGICAL NOTRAN, NOUNIT, UPPER
260 INTEGER I, IMAX, J, JFIRST, JINC, JLAST
261 REAL BIGNUM, GROW, REC, SMLNUM, SUMJ, TJJ, TJJS,
262 $ tmax, tscal, uscal, xbnd, xj, xmax
267 REAL SASUM, SDOT, SLAMCH, SLANGE
268 EXTERNAL lsame, isamax, sasum, sdot, slamch, slange
274 INTRINSIC abs, max, min
279 upper = lsame( uplo,
'U' )
280 notran = lsame( trans,
'N' )
281 nounit = lsame( diag,
'N' )
285 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
287 ELSE IF( .NOT.notran .AND. .NOT.lsame( trans,
'T' ) .AND. .NOT.
288 $ lsame( trans,
'C' ) )
THEN
290 ELSE IF( .NOT.nounit .AND. .NOT.lsame( diag,
'U' ) )
THEN
292 ELSE IF( .NOT.lsame( normin,
'Y' ) .AND. .NOT.
293 $ lsame( normin,
'N' ) )
THEN
295 ELSE IF( n.LT.0 )
THEN
297 ELSE IF( lda.LT.max( 1, n ) )
THEN
301 CALL xerbla(
'SLATRS', -info )
313 smlnum = slamch(
'Safe minimum' ) / slamch(
'Precision' )
314 bignum = one / smlnum
316 IF( lsame( normin,
'N' ) )
THEN
325 cnorm( j ) = sasum( j-1, a( 1, j ), 1 )
332 cnorm( j ) = sasum( n-j, a( j+1, j ), 1 )
341 imax = isamax( n, cnorm, 1 )
343 IF( tmax.LE.bignum )
THEN
350 IF ( tmax.LE.slamch(
'Overflow') )
THEN
352 tscal = one / ( smlnum*tmax )
353 CALL sscal( n, tscal, cnorm, 1 )
365 tmax = max( slange(
'M', j-1, 1, a( 1, j ), 1, sumj ),
373 tmax = max( slange(
'M', n-j, 1, a( j+1, j ), 1,
378 IF( tmax.LE.slamch(
'Overflow') )
THEN
379 tscal = one / ( smlnum*tmax )
381 IF( cnorm( j ).LE.slamch(
'Overflow') )
THEN
382 cnorm( j ) = cnorm( j )*tscal
389 cnorm( j ) = cnorm( j ) +
390 $ tscal * abs( a( i, j ) )
394 cnorm( j ) = cnorm( j ) +
395 $ tscal * abs( a( i, j ) )
403 CALL strsv( uplo, trans, diag, n, a, lda, x, 1 )
412 j = isamax( n, x, 1 )
429 IF( tscal.NE.one )
THEN
441 grow = one / max( xbnd, smlnum )
443 DO 30 j = jfirst, jlast, jinc
452 tjj = abs( a( j, j ) )
453 xbnd = min( xbnd, min( one, tjj )*grow )
454 IF( tjj+cnorm( j ).GE.smlnum )
THEN
458 grow = grow*( tjj / ( tjj+cnorm( j ) ) )
473 grow = min( one, one / max( xbnd, smlnum ) )
474 DO 40 j = jfirst, jlast, jinc
483 grow = grow*( one / ( one+cnorm( j ) ) )
502 IF( tscal.NE.one )
THEN
514 grow = one / max( xbnd, smlnum )
516 DO 60 j = jfirst, jlast, jinc
525 xj = one + cnorm( j )
526 grow = min( grow, xbnd / xj )
530 tjj = abs( a( j, j ) )
532 $ xbnd = xbnd*( tjj / xj )
534 grow = min( grow, xbnd )
541 grow = min( one, one / max( xbnd, smlnum ) )
542 DO 70 j = jfirst, jlast, jinc
551 xj = one + cnorm( j )
558 IF( ( grow*tscal ).GT.smlnum )
THEN
563 CALL strsv( uplo, trans, diag, n, a, lda, x, 1 )
568 IF( xmax.GT.bignum )
THEN
573 scale = bignum / xmax
574 CALL sscal( n, scale, x, 1 )
582 DO 100 j = jfirst, jlast, jinc
588 tjjs = a( j, j )*tscal
595 IF( tjj.GT.smlnum )
THEN
599 IF( tjj.LT.one )
THEN
600 IF( xj.GT.tjj*bignum )
THEN
605 CALL sscal( n, rec, x, 1 )
610 x( j ) = x( j ) / tjjs
612 ELSE IF( tjj.GT.zero )
THEN
616 IF( xj.GT.tjj*bignum )
THEN
621 rec = ( tjj*bignum ) / xj
622 IF( cnorm( j ).GT.one )
THEN
627 rec = rec / cnorm( j )
629 CALL sscal( n, rec, x, 1 )
633 x( j ) = x( j ) / tjjs
655 IF( cnorm( j ).GT.( bignum-xmax )*rec )
THEN
660 CALL sscal( n, rec, x, 1 )
663 ELSE IF( xj*cnorm( j ).GT.( bignum-xmax ) )
THEN
667 CALL sscal( n, half, x, 1 )
677 CALL saxpy( j-1, -x( j )*tscal, a( 1, j ), 1, x,
679 i = isamax( j-1, x, 1 )
688 CALL saxpy( n-j, -x( j )*tscal, a( j+1, j ), 1,
690 i = j + isamax( n-j, x( j+1 ), 1 )
700 DO 140 j = jfirst, jlast, jinc
707 rec = one / max( xmax, one )
708 IF( cnorm( j ).GT.( bignum-xj )*rec )
THEN
714 tjjs = a( j, j )*tscal
719 IF( tjj.GT.one )
THEN
723 rec = min( one, rec*tjj )
726 IF( rec.LT.one )
THEN
727 CALL sscal( n, rec, x, 1 )
734 IF( uscal.EQ.one )
THEN
740 sumj = sdot( j-1, a( 1, j ), 1, x, 1 )
741 ELSE IF( j.LT.n )
THEN
742 sumj = sdot( n-j, a( j+1, j ), 1, x( j+1 ), 1 )
750 sumj = sumj + ( a( i, j )*uscal )*x( i )
752 ELSE IF( j.LT.n )
THEN
754 sumj = sumj + ( a( i, j )*uscal )*x( i )
759 IF( uscal.EQ.tscal )
THEN
764 x( j ) = x( j ) - sumj
767 tjjs = a( j, j )*tscal
777 IF( tjj.GT.smlnum )
THEN
781 IF( tjj.LT.one )
THEN
782 IF( xj.GT.tjj*bignum )
THEN
787 CALL sscal( n, rec, x, 1 )
792 x( j ) = x( j ) / tjjs
793 ELSE IF( tjj.GT.zero )
THEN
797 IF( xj.GT.tjj*bignum )
THEN
801 rec = ( tjj*bignum ) / xj
802 CALL sscal( n, rec, x, 1 )
806 x( j ) = x( j ) / tjjs
825 x( j ) = x( j ) / tjjs - sumj
827 xmax = max( xmax, abs( x( j ) ) )
830 scale = scale / tscal
835 IF( tscal.NE.one )
THEN
836 CALL sscal( n, one / tscal, cnorm, 1 )
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine slatrs(UPLO, TRANS, DIAG, NORMIN, N, A, LDA, X, SCALE, CNORM, INFO)
SLATRS solves a triangular system of equations with the scale factor set to prevent overflow.
subroutine sscal(N, SA, SX, INCX)
SSCAL
subroutine saxpy(N, SA, SX, INCX, SY, INCY)
SAXPY
subroutine strsv(UPLO, TRANS, DIAG, N, A, LDA, X, INCX)
STRSV