74 SUBROUTINE zqrt04(M,N,NB,RESULT)
85 DOUBLE PRECISION RESULT(6)
91 COMPLEX*16,
ALLOCATABLE :: AF(:,:), Q(:,:),
92 $ r(:,:), rwork(:), work( : ), t(:,:),
93 $ cf(:,:), df(:,:), a(:,:), c(:,:), d(:,:)
98 parameter( zero = 0.0, one = (1.0,0.0), czero=(0.0,0.0) )
101 INTEGER INFO, J, K, L, LWORK
102 DOUBLE PRECISION ANORM, EPS, RESID, CNORM, DNORM
108 DOUBLE PRECISION DLAMCH
109 DOUBLE PRECISION ZLANGE, ZLANSY
111 EXTERNAL dlamch, zlange, zlansy, lsame
117 DATA iseed / 1988, 1989, 1990, 1991 /
119 eps = dlamch(
'Epsilon' )
122 lwork = max(2,l)*max(2,l)*nb
126 ALLOCATE ( a(m,n), af(m,n), q(m,m), r(m,l), rwork(l),
127 $ work(lwork), t(nb,n), c(m,n), cf(m,n),
134 CALL zlarnv( 2, iseed, m, a( 1, j ) )
136 CALL zlacpy(
'Full', m, n, a, m, af, m )
140 CALL zgeqrt( m, n, nb, af, m, t, ldt, work, info )
144 CALL zlaset(
'Full', m, m, czero, one, q, m )
145 CALL zgemqrt(
'R',
'N', m, m, k, nb, af, m, t, ldt, q, m,
150 CALL zlaset(
'Full', m, n, czero, czero, r, m )
151 CALL zlacpy(
'Upper', m, n, af, m, r, m )
155 CALL zgemm(
'C',
'N', m, n, m, -one, q, m, a, m, one, r, m )
156 anorm = zlange(
'1', m, n, a, m, rwork )
157 resid = zlange(
'1', m, n, r, m, rwork )
158 IF( anorm.GT.zero )
THEN
159 result( 1 ) = resid / (eps*max(1,m)*anorm)
166 CALL zlaset(
'Full', m, m, czero, one, r, m )
167 CALL zherk(
'U',
'C', m, m, dreal(-one), q, m, dreal(one), r, m )
168 resid = zlansy(
'1',
'Upper', m, r, m, rwork )
169 result( 2 ) = resid / (eps*max(1,m))
174 CALL zlarnv( 2, iseed, m, c( 1, j ) )
176 cnorm = zlange(
'1', m, n, c, m, rwork)
177 CALL zlacpy(
'Full', m, n, c, m, cf, m )
181 CALL zgemqrt(
'L',
'N', m, n, k, nb, af, m, t, nb, cf, m,
186 CALL zgemm(
'N',
'N', m, n, m, -one, q, m, c, m, one, cf, m )
187 resid = zlange(
'1', m, n, cf, m, rwork )
188 IF( cnorm.GT.zero )
THEN
189 result( 3 ) = resid / (eps*max(1,m)*cnorm)
196 CALL zlacpy(
'Full', m, n, c, m, cf, m )
200 CALL zgemqrt(
'L',
'C', m, n, k, nb, af, m, t, nb, cf, m,
205 CALL zgemm(
'C',
'N', m, n, m, -one, q, m, c, m, one, cf, m )
206 resid = zlange(
'1', m, n, cf, m, rwork )
207 IF( cnorm.GT.zero )
THEN
208 result( 4 ) = resid / (eps*max(1,m)*cnorm)
216 CALL zlarnv( 2, iseed, n, d( 1, j ) )
218 dnorm = zlange(
'1', n, m, d, n, rwork)
219 CALL zlacpy(
'Full', n, m, d, n, df, n )
223 CALL zgemqrt(
'R',
'N', n, m, k, nb, af, m, t, nb, df, n,
228 CALL zgemm(
'N',
'N', n, m, m, -one, d, n, q, m, one, df, n )
229 resid = zlange(
'1', n, m, df, n, rwork )
230 IF( cnorm.GT.zero )
THEN
231 result( 5 ) = resid / (eps*max(1,m)*dnorm)
238 CALL zlacpy(
'Full', n, m, d, n, df, n )
242 CALL zgemqrt(
'R',
'C', n, m, k, nb, af, m, t, nb, df, n,
247 CALL zgemm(
'N',
'C', n, m, m, -one, d, n, q, m, one, df, n )
248 resid = zlange(
'1', n, m, df, n, rwork )
249 IF( cnorm.GT.zero )
THEN
250 result( 6 ) = resid / (eps*max(1,m)*dnorm)
257 DEALLOCATE ( a, af, q, r, rwork, work, t, c, d, cf, df)
subroutine zherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
ZHERK
subroutine zqrt04(M, N, NB, RESULT)
ZQRT04
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine zgemqrt(SIDE, TRANS, M, N, K, NB, V, LDV, T, LDT, C, LDC, WORK, INFO)
ZGEMQRT
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
subroutine zlarnv(IDIST, ISEED, N, X)
ZLARNV returns a vector of random numbers from a uniform or normal distribution.
subroutine zgeqrt(M, N, NB, A, LDA, T, LDT, WORK, INFO)
ZGEQRT