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 ) ) /