149 SUBROUTINE dget51( ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK,
158 INTEGER itype, lda, ldb, ldu, ldv, n
159 DOUBLE PRECISION result
162 DOUBLE PRECISION a( lda, * ),
b( ldb, * ), u( ldu, * ),
163 $ v( ldv, * ), work( * )
169 DOUBLE PRECISION zero, one, ten
170 parameter( zero = 0.0d0, one = 1.0d0, ten = 10.0d0 )
173 INTEGER jcol, jdiag, jrow
174 DOUBLE PRECISION anorm, ulp, unfl, wnorm
184 INTRINSIC dble, max, min
194 unfl =
dlamch(
'Safe minimum' )
199 IF( itype.LT.1 .OR. itype.GT.3 )
THEN
204 IF( itype.LE.2 )
THEN
208 anorm = max(
dlange(
'1', n, n, a, lda, work ), unfl )
210 IF( itype.EQ.1 )
THEN
214 CALL
dlacpy(
' ', n, n, a, lda, work, n )
215 CALL
dgemm(
'N',
'N', n, n, n, one, u, ldu,
b, ldb, zero,
216 $ work( n**2+1 ), n )
218 CALL
dgemm(
'N',
'C', n, n, n, -one, work( n**2+1 ), n, v,
219 $ ldv, one, work, n )
225 CALL
dlacpy(
' ', n, n,
b, ldb, work, n )
229 work( jrow+n*( jcol-1 ) ) = work( jrow+n*( jcol-1 ) )
237 wnorm =
dlange(
'1', n, n, work, n, work( n**2+1 ) )
239 IF( anorm.GT.wnorm )
THEN
240 result = ( wnorm / anorm ) / ( n*ulp )
242 IF( anorm.LT.one )
THEN
243 result = ( min( wnorm, n*anorm ) / anorm ) / ( n*ulp )
245 result = min( wnorm / anorm, dble( n ) ) / ( n*ulp )
255 CALL
dgemm(
'N',
'C', n, n, n, one, u, ldu, u, ldu, zero, work,
259 work( ( n+1 )*( jdiag-1 )+1 ) = work( ( n+1 )*( jdiag-1 )+
263 result = min(
dlange(
'1', n, n, work, n, work( n**2+1 ) ),
264 $ dble( n ) ) / ( n*ulp )
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real b(3) integer i
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 ...
subroutine dget51(ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK, RESULT)
DGET51
double precision function dlamch(CMACH)
DLAMCH
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.