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
118 DATA iseed / 1988, 1989, 1990, 1991 /
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 cherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
CHERK
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
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
real function clansy(NORM, UPLO, N, A, LDA, WORK)
CLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex symmetric matrix.
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
real function slamch(CMACH)
SLAMCH
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
logical function lsame(CA, CB)
LSAME
subroutine cgeqrt(M, N, NB, A, LDA, T, LDT, WORK, INFO)
CGEQRT