146 SUBROUTINE dglmts( N, M, P, A, AF, LDA, B, BF, LDB, D, DF, X, U,
147 $ work, lwork, rwork, result )
155 INTEGER LDA, LDB, LWORK, M, N, P
156 DOUBLE PRECISION RESULT
162 DOUBLE PRECISION A( lda, * ), AF( lda, * ), B( ldb, * ),
163 $ bf( ldb, * ), d( * ), df( * ), rwork( * ),
164 $ u( * ), work( lwork ), x( * )
167 DOUBLE PRECISION ZERO, ONE
168 parameter ( zero = 0.0d+0, one = 1.0d+0 )
172 DOUBLE PRECISION ANORM, BNORM, DNORM, EPS, UNFL, XNORM, YNORM
175 DOUBLE PRECISION DASUM, DLAMCH, DLANGE
176 EXTERNAL dasum, dlamch, dlange
187 eps = dlamch(
'Epsilon' )
188 unfl = dlamch(
'Safe minimum' )
189 anorm = max( dlange(
'1', n, m, a, lda, rwork ), unfl )
190 bnorm = max( dlange(
'1', n, p, b, ldb, rwork ), unfl )
195 CALL dlacpy(
'Full', n, m, a, lda, af, lda )
196 CALL dlacpy(
'Full', n, p, b, ldb, bf, ldb )
197 CALL dcopy( n, d, 1, df, 1 )
201 CALL dggglm( n, m, p, af, lda, bf, ldb, df, x, u, work, lwork,
210 CALL dcopy( n, d, 1, df, 1 )
211 CALL dgemv(
'No transpose', n, m, -one, a, lda, x, 1, one, df, 1 )
213 CALL dgemv(
'No transpose', n, p, -one, b, ldb, u, 1, one, df, 1 )
215 dnorm = dasum( n, df, 1 )
216 xnorm = dasum( m, x, 1 ) + dasum( p, u, 1 )
217 ynorm = anorm + bnorm
219 IF( xnorm.LE.zero )
THEN
222 result = ( ( dnorm / ynorm ) / xnorm ) / eps
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
subroutine dgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
DGEMV
subroutine dglmts(N, M, P, A, AF, LDA, B, BF, LDB, D, DF, X, U, WORK, LWORK, RWORK, RESULT)
DGLMTS
subroutine dggglm(N, M, P, A, LDA, B, LDB, D, X, Y, WORK, LWORK, INFO)
DGGGLM