74 SUBROUTINE cqrt04(M,N,NB,RESULT)
91 COMPLEX,
ALLOCATABLE :: AF(:,:), Q(:,:),
92 $ r(:,:), work( : ), t(:,:),
93 $ cf(:,:), df(:,:), a(:,:), c(:,:), d(:,:)
94 REAL,
ALLOCATABLE :: RWORK(:)
99 parameter( zero = 0.0, one = (1.0,0.0), czero=(0.0,0.0) )
102 INTEGER INFO, J, K, L, LWORK
103 REAL ANORM, EPS, RESID, CNORM, DNORM
112 EXTERNAL slamch, clange, clansy, lsame
118 DATA iseed / 1988, 1989, 1990, 1991 /
120 eps = slamch(
'Epsilon' )
123 lwork = max(2,l)*max(2,l)*nb
127 ALLOCATE ( a(m,n), af(m,n), q(m,m), r(m,l), rwork(l),
128 $ work(lwork), t(nb,n), c(m,n), cf(m,n),
135 CALL clarnv( 2, iseed, m, a( 1, j ) )
137 CALL clacpy(
'Full', m, n, a, m, af, m )
141 CALL cgeqrt( m, n, nb, af, m, t, ldt, work, info )
145 CALL claset(
'Full', m, m, czero, one, q, m )
146 CALL cgemqrt(
'R',
'N', m, m, k, nb, af, m, t, ldt, q, m,
151 CALL claset(
'Full', m, n, czero, czero, r, m )
152 CALL clacpy(
'Upper', m, n, af, m, r, m )
156 CALL cgemm(
'C',
'N', m, n, m, -one, q, m, a, m, one, r, m )
157 anorm = clange(
'1', m, n, a, m, rwork )
158 resid = clange(
'1', m, n, r, m, rwork )
159 IF( anorm.GT.zero )
THEN
160 result( 1 ) = resid / (eps*max(1,m)*anorm)
167 CALL claset(
'Full', m, m, czero, one, r, m )
168 CALL cherk(
'U',
'C', m, m,
REAL(-ONE), Q, M,
REAL(ONE), R, M )
169 resid = clansy(
'1',
'Upper', m, r, m, rwork )
170 result( 2 ) = resid / (eps*max(1,m))
175 CALL clarnv( 2, iseed, m, c( 1, j ) )
177 cnorm = clange(
'1', m, n, c, m, rwork)
178 CALL clacpy(
'Full', m, n, c, m, cf, m )
182 CALL cgemqrt(
'L',
'N', m, n, k, nb, af, m, t, nb, cf, m,
187 CALL cgemm(
'N',
'N', m, n, m, -one, q, m, c, m, one, cf, m )
188 resid = clange(
'1', m, n, cf, m, rwork )
189 IF( cnorm.GT.zero )
THEN
190 result( 3 ) = resid / (eps*max(1,m)*cnorm)
197 CALL clacpy(
'Full', m, n, c, m, cf, m )
201 CALL cgemqrt(
'L',
'C', m, n, k, nb, af, m, t, nb, cf, m,
206 CALL cgemm(
'C',
'N', m, n, m, -one, q, m, c, m, one, cf, m )
207 resid = clange(
'1', m, n, cf, m, rwork )
208 IF( cnorm.GT.zero )
THEN
209 result( 4 ) = resid / (eps*max(1,m)*cnorm)
217 CALL clarnv( 2, iseed, n, d( 1, j ) )
219 dnorm = clange(
'1', n, m, d, n, rwork)
220 CALL clacpy(
'Full', n, m, d, n, df, n )
224 CALL cgemqrt(
'R',
'N', n, m, k, nb, af, m, t, nb, df, n,
229 CALL cgemm(
'N',
'N', n, m, m, -one, d, n, q, m, one, df, n )
230 resid = clange(
'1', n, m, df, n, rwork )
231 IF( cnorm.GT.zero )
THEN
232 result( 5 ) = resid / (eps*max(1,m)*dnorm)
239 CALL clacpy(
'Full', n, m, d, n, df, n )
243 CALL cgemqrt(
'R',
'C', n, m, k, nb, af, m, t, nb, df, n,
248 CALL cgemm(
'N',
'C', n, m, m, -one, d, n, q, m, one, df, n )
249 resid = clange(
'1', n, m, df, n, rwork )
250 IF( cnorm.GT.zero )
THEN
251 result( 6 ) = resid / (eps*max(1,m)*dnorm)
258 DEALLOCATE ( a, af, q, r, rwork, work, t, c, d, cf, df)
subroutine cqrt04(M, N, NB, RESULT)
CQRT04
subroutine cherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
CHERK
subroutine clarnv(IDIST, ISEED, N, X)
CLARNV returns a vector of random numbers from a uniform or normal distribution.
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
subroutine cgemqrt(SIDE, TRANS, M, N, K, NB, V, LDV, T, LDT, C, LDC, WORK, INFO)
CGEMQRT
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
subroutine cgeqrt(M, N, NB, A, LDA, T, LDT, WORK, INFO)
CGEQRT