140 SUBROUTINE dbdt01( M, N, KD, A, LDA, Q, LDQ, D, E, PT, LDPT, WORK,
149 INTEGER kd, lda, ldpt, ldq, m, n
150 DOUBLE PRECISION resid
153 DOUBLE PRECISION a( lda, * ), d( * ), e( * ), pt( ldpt, * ),
154 $ q( ldq, * ), work( * )
160 DOUBLE PRECISION zero, one
161 parameter( zero = 0.0d+0, one = 1.0d+0 )
165 DOUBLE PRECISION anorm, eps
175 INTRINSIC dble, max, min
181 IF( m.LE.0 .OR. n.LE.0 )
THEN
193 IF( kd.NE.0 .AND. m.GE.n )
THEN
198 CALL
dcopy( m, a( 1,
j ), 1, work, 1 )
200 work( m+i ) = d( i )*pt( i,
j ) + e( i )*pt( i+1,
j )
202 work( m+n ) = d( n )*pt( n,
j )
203 CALL
dgemv(
'No transpose', m, n, -one, q, ldq,
204 $ work( m+1 ), 1, one, work, 1 )
205 resid = max( resid,
dasum( m, work, 1 ) )
207 ELSE IF( kd.LT.0 )
THEN
212 CALL
dcopy( m, a( 1,
j ), 1, work, 1 )
214 work( m+i ) = d( i )*pt( i,
j ) + e( i )*pt( i+1,
j )
216 work( m+m ) = d( m )*pt( m,
j )
217 CALL
dgemv(
'No transpose', m, m, -one, q, ldq,
218 $ work( m+1 ), 1, one, work, 1 )
219 resid = max( resid,
dasum( m, work, 1 ) )
226 CALL
dcopy( m, a( 1,
j ), 1, work, 1 )
227 work( m+1 ) = d( 1 )*pt( 1,
j )
229 work( m+i ) = e( i-1 )*pt( i-1,
j ) +
232 CALL
dgemv(
'No transpose', m, m, -one, q, ldq,
233 $ work( m+1 ), 1, one, work, 1 )
234 resid = max( resid,
dasum( m, work, 1 ) )
243 CALL
dcopy( m, a( 1,
j ), 1, work, 1 )
245 work( m+i ) = d( i )*pt( i,
j )
247 CALL
dgemv(
'No transpose', m, n, -one, q, ldq,
248 $ work( m+1 ), 1, one, work, 1 )
249 resid = max( resid,
dasum( m, work, 1 ) )
253 CALL
dcopy( m, a( 1,
j ), 1, work, 1 )
255 work( m+i ) = d( i )*pt( i,
j )
257 CALL
dgemv(
'No transpose', m, m, -one, q, ldq,
258 $ work( m+1 ), 1, one, work, 1 )
259 resid = max( resid,
dasum( m, work, 1 ) )
266 anorm =
dlange(
'1', m, n, a, lda, work )
267 eps =
dlamch(
'Precision' )
269 IF( anorm.LE.zero )
THEN
273 IF( anorm.GE.resid )
THEN
274 resid = ( resid / anorm ) / ( dble( n )*eps )
276 IF( anorm.LT.one )
THEN
277 resid = ( min( resid, dble( n )*anorm ) / anorm ) /
280 resid = min( resid / anorm, dble( n ) ) /
subroutine dbdt01(M, N, KD, A, LDA, Q, LDQ, D, E, PT, LDPT, WORK, RESID)
DBDT01
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
double precision function dlange(NORM, M, N, A, LDA, WORK)
DLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
double precision function dlamch(CMACH)
DLAMCH
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
double precision function dasum(N, DX, INCX)
DASUM
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV