146 INTEGER INFO, LDA, LDB, N, NRHS
150 DOUBLE PRECISION A( lda, * ), B( ldb, * )
157 parameter( one = 1.0d+0 )
162 DOUBLE PRECISION AK, AKM1, AKM1K, BK, BKM1, DENOM
177 upper = lsame( uplo,
'U' )
178 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
180 ELSE IF( n.LT.0 )
THEN
182 ELSE IF( nrhs.LT.0 )
THEN
184 ELSE IF( lda.LT.max( 1, n ) )
THEN
186 ELSE IF( ldb.LT.max( 1, n ) )
THEN
190 CALL
xerbla(
'DSYTRS_ROOK', -info )
196 IF( n.EQ.0 .OR. nrhs.EQ.0 )
216 IF( ipiv( k ).GT.0 )
THEN
224 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
229 CALL
dger( k-1, nrhs, -one, a( 1, k ), 1, b( k, 1 ), ldb,
234 CALL
dscal( nrhs, one / a( k, k ), b( k, 1 ), ldb )
244 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
248 $ CALL
dswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ), ldb )
254 CALL
dger( k-2, nrhs, -one, a( 1, k ), 1, b( k, 1 ),
255 $ ldb, b( 1, 1 ), ldb )
256 CALL
dger( k-2, nrhs, -one, a( 1, k-1 ), 1, b( k-1, 1 ),
257 $ ldb, b( 1, 1 ), ldb )
263 akm1 = a( k-1, k-1 ) / akm1k
264 ak = a( k, k ) / akm1k
265 denom = akm1*ak - one
267 bkm1 = b( k-1, j ) / akm1k
268 bk = b( k, j ) / akm1k
269 b( k-1, j ) = ( ak*bkm1-bk ) / denom
270 b( k, j ) = ( akm1*bk-bkm1 ) / denom
291 IF( ipiv( k ).GT.0 )
THEN
299 $ CALL
dgemv(
'Transpose', k-1, nrhs, -one, b,
300 $ ldb, a( 1, k ), 1, one, b( k, 1 ), ldb )
306 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
316 CALL
dgemv(
'Transpose', k-1, nrhs, -one, b,
317 $ ldb, a( 1, k ), 1, one, b( k, 1 ), ldb )
318 CALL
dgemv(
'Transpose', k-1, nrhs, -one, b,
319 $ ldb, a( 1, k+1 ), 1, one, b( k+1, 1 ), ldb )
326 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
330 $ CALL
dswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ), ldb )
355 IF( ipiv( k ).GT.0 )
THEN
363 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
369 $ CALL
dger( n-k, nrhs, -one, a( k+1, k ), 1, b( k, 1 ),
370 $ ldb, b( k+1, 1 ), ldb )
374 CALL
dscal( nrhs, one / a( k, k ), b( k, 1 ), ldb )
384 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
388 $ CALL
dswap( nrhs, b( k+1, 1 ), ldb, b( kp, 1 ), ldb )
394 CALL
dger( n-k-1, nrhs, -one, a( k+2, k ), 1, b( k, 1 ),
395 $ ldb, b( k+2, 1 ), ldb )
396 CALL
dger( n-k-1, nrhs, -one, a( k+2, k+1 ), 1,
397 $ b( k+1, 1 ), ldb, b( k+2, 1 ), ldb )
403 akm1 = a( k, k ) / akm1k
404 ak = a( k+1, k+1 ) / akm1k
405 denom = akm1*ak - one
407 bkm1 = b( k, j ) / akm1k
408 bk = b( k+1, j ) / akm1k
409 b( k, j ) = ( ak*bkm1-bk ) / denom
410 b( k+1, j ) = ( akm1*bk-bkm1 ) / denom
431 IF( ipiv( k ).GT.0 )
THEN
439 $ CALL
dgemv(
'Transpose', n-k, nrhs, -one, b( k+1, 1 ),
440 $ ldb, a( k+1, k ), 1, one, b( k, 1 ), ldb )
446 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
456 CALL
dgemv(
'Transpose', n-k, nrhs, -one, b( k+1, 1 ),
457 $ ldb, a( k+1, k ), 1, one, b( k, 1 ), ldb )
458 CALL
dgemv(
'Transpose', n-k, nrhs, -one, b( k+1, 1 ),
459 $ ldb, a( k+1, k-1 ), 1, one, b( k-1, 1 ),
467 $ CALL
dswap( nrhs, b( k, 1 ), ldb, b( kp, 1 ), ldb )
471 $ CALL
dswap( nrhs, b( k-1, 1 ), ldb, b( kp, 1 ), ldb )
subroutine dger(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
DGER
subroutine dsytrs_rook(UPLO, N, NRHS, A, LDA, IPIV, B, LDB, INFO)
DSYTRS_ROOK
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine dscal(N, DA, DX, INCX)
DSCAL
subroutine dswap(N, DX, INCX, DY, INCY)
DSWAP
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV