189 SUBROUTINE dpbrfs( UPLO, N, KD, NRHS, AB, LDAB, AFB, LDAFB, B,
190 $ ldb, x, ldx, ferr, berr, work, iwork, info )
199 INTEGER INFO, KD, LDAB, LDAFB, LDB, LDX, N, NRHS
203 DOUBLE PRECISION AB( ldab, * ), AFB( ldafb, * ), B( ldb, * ),
204 $ berr( * ), ferr( * ), work( * ), x( ldx, * )
211 parameter ( itmax = 5 )
212 DOUBLE PRECISION ZERO
213 parameter ( zero = 0.0d+0 )
215 parameter ( one = 1.0d+0 )
217 parameter ( two = 2.0d+0 )
218 DOUBLE PRECISION THREE
219 parameter ( three = 3.0d+0 )
223 INTEGER COUNT, I, J, K, KASE, L, NZ
224 DOUBLE PRECISION EPS, LSTRES, S, SAFE1, SAFE2, SAFMIN, XK
233 INTRINSIC abs, max, min
237 DOUBLE PRECISION DLAMCH
238 EXTERNAL lsame, dlamch
245 upper = lsame( uplo,
'U' )
246 IF( .NOT.upper .AND. .NOT.lsame( uplo,
'L' ) )
THEN
248 ELSE IF( n.LT.0 )
THEN
250 ELSE IF( kd.LT.0 )
THEN
252 ELSE IF( nrhs.LT.0 )
THEN
254 ELSE IF( ldab.LT.kd+1 )
THEN
256 ELSE IF( ldafb.LT.kd+1 )
THEN
258 ELSE IF( ldb.LT.max( 1, n ) )
THEN
260 ELSE IF( ldx.LT.max( 1, n ) )
THEN
264 CALL xerbla(
'DPBRFS', -info )
270 IF( n.EQ.0 .OR. nrhs.EQ.0 )
THEN
280 nz = min( n+1, 2*kd+2 )
281 eps = dlamch(
'Epsilon' )
282 safmin = dlamch(
'Safe minimum' )
298 CALL dcopy( n, b( 1, j ), 1, work( n+1 ), 1 )
299 CALL dsbmv( uplo, n, kd, -one, ab, ldab, x( 1, j ), 1, one,
312 work( i ) = abs( b( i, j ) )
320 xk = abs( x( k, j ) )
322 DO 40 i = max( 1, k-kd ), k - 1
323 work( i ) = work( i ) + abs( ab( l+i, k ) )*xk
324 s = s + abs( ab( l+i, k ) )*abs( x( i, j ) )
326 work( k ) = work( k ) + abs( ab( kd+1, k ) )*xk + s
331 xk = abs( x( k, j ) )
332 work( k ) = work( k ) + abs( ab( 1, k ) )*xk
334 DO 60 i = k + 1, min( n, k+kd )
335 work( i ) = work( i ) + abs( ab( l+i, k ) )*xk
336 s = s + abs( ab( l+i, k ) )*abs( x( i, j ) )
338 work( k ) = work( k ) + s
343 IF( work( i ).GT.safe2 )
THEN
344 s = max( s, abs( work( n+i ) ) / work( i ) )
346 s = max( s, ( abs( work( n+i ) )+safe1 ) /
347 $ ( work( i )+safe1 ) )
358 IF( berr( j ).GT.eps .AND. two*berr( j ).LE.lstres .AND.
359 $ count.LE.itmax )
THEN
363 CALL dpbtrs( uplo, n, kd, 1, afb, ldafb, work( n+1 ), n,
365 CALL daxpy( n, one, work( n+1 ), 1, x( 1, j ), 1 )
394 IF( work( i ).GT.safe2 )
THEN
395 work( i ) = abs( work( n+i ) ) + nz*eps*work( i )
397 work( i ) = abs( work( n+i ) ) + nz*eps*work( i ) + safe1
403 CALL dlacn2( n, work( 2*n+1 ), work( n+1 ), iwork, ferr( j ),
410 CALL dpbtrs( uplo, n, kd, 1, afb, ldafb, work( n+1 ), n,
413 work( n+i ) = work( n+i )*work( i )
415 ELSE IF( kase.EQ.2 )
THEN
420 work( n+i ) = work( n+i )*work( i )
422 CALL dpbtrs( uplo, n, kd, 1, afb, ldafb, work( n+1 ), n,
432 lstres = max( lstres, abs( x( i, j ) ) )
435 $ ferr( j ) = ferr( j ) / lstres
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
subroutine dpbrfs(UPLO, N, KD, NRHS, AB, LDAB, AFB, LDAFB, B, LDB, X, LDX, FERR, BERR, WORK, IWORK, INFO)
DPBRFS
subroutine daxpy(N, DA, DX, INCX, DY, INCY)
DAXPY
subroutine dlacn2(N, V, X, ISGN, EST, KASE, ISAVE)
DLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vec...
subroutine dsbmv(UPLO, N, K, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DSBMV
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine dpbtrs(UPLO, N, KD, NRHS, AB, LDAB, B, LDB, INFO)
DPBTRS