156 SUBROUTINE zget54( N, A, LDA, B, LDB, S, LDS, T, LDT, U, LDU, V,
157 $ ldv, work, result )
165 INTEGER LDA, LDB, LDS, LDT, LDU, LDV, N
166 DOUBLE PRECISION RESULT
169 COMPLEX*16 A( lda, * ), B( ldb, * ), S( lds, * ),
170 $ t( ldt, * ), u( ldu, * ), v( ldv, * ),
177 DOUBLE PRECISION ZERO, ONE
178 parameter( zero = 0.0d+0, one = 1.0d+0 )
179 COMPLEX*16 CZERO, CONE
180 parameter( czero = ( 0.0d+0, 0.0d+0 ),
181 $ cone = ( 1.0d+0, 0.0d+0 ) )
184 DOUBLE PRECISION ABNORM, ULP, UNFL, WNORM
187 DOUBLE PRECISION DUM( 1 )
190 DOUBLE PRECISION DLAMCH, ZLANGE
191 EXTERNAL dlamch, zlange
197 INTRINSIC dble, max, min
207 unfl = dlamch(
'Safe minimum' )
208 ulp = dlamch(
'Epsilon' )*dlamch(
'Base' )
212 CALL
zlacpy(
'Full', n, n, a, lda, work, n )
213 CALL
zlacpy(
'Full', n, n, b, ldb, work( n*n+1 ), n )
214 abnorm = max( zlange(
'1', n, 2*n, work, n, dum ), unfl )
218 CALL
zlacpy(
' ', n, n, a, lda, work, n )
219 CALL
zgemm(
'N',
'N', n, n, n, cone, u, ldu, s, lds, czero,
222 CALL
zgemm(
'N',
'C', n, n, n, -cone, work( n*n+1 ), n, v, ldv,
227 CALL
zlacpy(
' ', n, n, b, ldb, work( n*n+1 ), n )
228 CALL
zgemm(
'N',
'N', n, n, n, cone, u, ldu, t, ldt, czero,
229 $ work( 2*n*n+1 ), n )
231 CALL
zgemm(
'N',
'C', n, n, n, -cone, work( 2*n*n+1 ), n, v, ldv,
232 $ cone, work( n*n+1 ), n )
236 wnorm = zlange(
'1', n, 2*n, work, n, dum )
238 IF( abnorm.GT.wnorm )
THEN
239 result = ( wnorm / abnorm ) / ( 2*n*ulp )
241 IF( abnorm.LT.one )
THEN
242 result = ( min( wnorm, 2*n*abnorm ) / abnorm ) / ( 2*n*ulp )
244 result = min( wnorm / abnorm, dble( 2*n ) ) / ( 2*n*ulp )
subroutine zget54(N, A, LDA, B, LDB, S, LDS, T, LDT, U, LDU, V, LDV, WORK, RESULT)
ZGET54
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM