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
117 DATA iseed / 1988, 1989, 1990, 1991 /
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
double precision function zlansy(NORM, UPLO, N, A, LDA, WORK)
ZLANSY 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 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
logical function lsame(CA, CB)
LSAME
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.
double precision function zlange(NORM, M, N, A, LDA, WORK)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
double precision function dlamch(CMACH)
DLAMCH
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine zgeqrt(M, N, NB, A, LDA, T, LDT, WORK, INFO)
ZGEQRT