149 SUBROUTINE sget51( ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK,
158 INTEGER itype, lda, ldb, ldu, ldv, n
162 REAL a( lda, * ),
b( ldb, * ), u( ldu, * ),
163 $ v( ldv, * ), work( * )
170 parameter( zero = 0.0, one = 1.0e0, ten = 10.0e0 )
173 INTEGER jcol, jdiag, jrow
174 REAL anorm, ulp, unfl, wnorm
184 INTRINSIC max, min, real
194 unfl =
slamch(
'Safe minimum' )
199 IF( itype.LT.1 .OR. itype.GT.3 )
THEN
204 IF( itype.LE.2 )
THEN
208 anorm = max(
slange(
'1', n, n, a, lda, work ), unfl )
210 IF( itype.EQ.1 )
THEN
214 CALL
slacpy(
' ', n, n, a, lda, work, n )
215 CALL
sgemm(
'N',
'N', n, n, n, one, u, ldu,
b, ldb, zero,
216 $ work( n**2+1 ), n )
218 CALL
sgemm(
'N',
'C', n, n, n, -one, work( n**2+1 ), n, v,
219 $ ldv, one, work, n )
225 CALL
slacpy(
' ', n, n,
b, ldb, work, n )
229 work( jrow+n*( jcol-1 ) ) = work( jrow+n*( jcol-1 ) )
237 wnorm =
slange(
'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,
REAL( N ) ) / ( n*ulp )
255 CALL
sgemm(
'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(
slange(
'1', n, n, work, n, work( n**2+1 ) ),
264 $
REAL( N ) ) / ( n*ulp )
REAL function slamch(CMACH)
SLAMCH
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real b(3) integer i
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
REAL function slange(NORM, M, N, A, LDA, WORK)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
subroutine sgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
SGEMM
subroutine sget51(ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK, RESULT)
SGET51