87 INTEGER LWORK, M, N, L, NB, LDT
89 DOUBLE PRECISION RESULT(6)
95 COMPLEX*16,
ALLOCATABLE :: AF(:,:), Q(:,:),
96 $ R(:,:), RWORK(:), WORK( : ), T(:,:),
97 $ CF(:,:), DF(:,:), A(:,:), C(:,:), D(:,:)
100 DOUBLE PRECISION ZERO
101 COMPLEX*16 ONE, CZERO
102 parameter( zero = 0.0, one = (1.0,0.0), czero=(0.0,0.0) )
105 INTEGER INFO, J, K, N2, NP1,i
106 DOUBLE PRECISION ANORM, EPS, RESID, CNORM, DNORM
112 DOUBLE PRECISION DLAMCH
113 DOUBLE PRECISION ZLANGE, ZLANSY
115 EXTERNAL dlamch, zlange, zlansy, lsame
118 DATA iseed / 1988, 1989, 1990, 1991 /
120 eps = dlamch(
'Epsilon' )
132 ALLOCATE(a(m,n2),af(m,n2),q(n2,n2),r(n2,n2),rwork(n2),
133 $ work(lwork),t(nb,m),c(n2,m),cf(n2,m),
139 CALL zlaset(
'Full', m, n2, czero, czero, a, m )
140 CALL zlaset(
'Full', nb, m, czero, czero, t, nb )
142 CALL zlarnv( 2, iseed, m-j+1, a( j, j ) )
146 CALL zlarnv( 2, iseed, m, a( 1, min(n+m,m+1) + j - 1 ) )
151 CALL zlarnv( 2, iseed, m-j+1, a( j, min(n+m,n+m-l+1)
158 CALL zlacpy(
'Full', m, n2, a, m, af, m )
162 CALL ztplqt( m,n,l,nb,af,m,af(1,np1),m,t,ldt,work,info)
166 CALL zlaset(
'Full', n2, n2, czero, one, q, n2 )
167 CALL zgemlqt(
'L',
'N', n2, n2, k, nb, af, m, t, ldt, q, n2,
172 CALL zlaset(
'Full', n2, n2, czero, czero, r, n2 )
173 CALL zlacpy(
'Lower', m, n2, af, m, r, n2 )
177 CALL zgemm(
'N',
'C', m, n2, n2, -one, a, m, q, n2, one, r, n2)
178 anorm = zlange(
'1', m, n2, a, m, rwork )
179 resid = zlange(
'1', m, n2, r, n2, rwork )
180 IF( anorm.GT.zero )
THEN
181 result( 1 ) = resid / (eps*anorm*max(1,n2))
188 CALL zlaset(
'Full', n2, n2, czero, one, r, n2 )
189 CALL zherk(
'U',
'N', n2, n2, dreal(-one), q, n2, dreal(one),
191 resid = zlansy(
'1',
'Upper', n2, r, n2, rwork )
192 result( 2 ) = resid / (eps*max(1,n2))
196 CALL zlaset(
'Full', n2, m, czero, one, c, n2 )
198 CALL zlarnv( 2, iseed, n2, c( 1, j ) )
200 cnorm = zlange(
'1', n2, m, c, n2, rwork)
201 CALL zlacpy(
'Full', n2, m, c, n2, cf, n2 )
205 CALL ztpmlqt(
'L',
'N', n,m,k,l,nb,af(1, np1),m,t,ldt,cf,n2,
206 $ cf(np1,1),n2,work,info)
210 CALL zgemm(
'N',
'N', n2, m, n2, -one, q, n2, c, n2, one, cf, n2 )
211 resid = zlange(
'1', n2, m, cf, n2, rwork )
212 IF( cnorm.GT.zero )
THEN
213 result( 3 ) = resid / (eps*max(1,n2)*cnorm)
221 CALL zlacpy(
'Full', n2, m, c, n2, cf, n2 )
225 CALL ztpmlqt(
'L',
'C',n,m,k,l,nb,af(1,np1),m,t,ldt,cf,n2,
226 $ cf(np1,1),n2,work,info)
230 CALL zgemm(
'C',
'N',n2,m,n2,-one,q,n2,c,n2,one,cf,n2)
231 resid = zlange(
'1', n2, m, cf, n2, rwork )
233 IF( cnorm.GT.zero )
THEN
234 result( 4 ) = resid / (eps*max(1,n2)*cnorm)
242 CALL zlarnv( 2, iseed, m, d( 1, j ) )
244 dnorm = zlange(
'1', m, n2, d, m, rwork)
245 CALL zlacpy(
'Full', m, n2, d, m, df, m )
249 CALL ztpmlqt(
'R',
'N',m,n,k,l,nb,af(1,np1),m,t,ldt,df,m,
250 $ df(1,np1),m,work,info)
254 CALL zgemm(
'N',
'N',m,n2,n2,-one,d,m,q,n2,one,df,m)
255 resid = zlange(
'1',m, n2,df,m,rwork )
256 IF( cnorm.GT.zero )
THEN
257 result( 5 ) = resid / (eps*max(1,n2)*dnorm)
264 CALL zlacpy(
'Full',m,n2,d,m,df,m )
268 CALL ztpmlqt(
'R',
'C',m,n,k,l,nb,af(1,np1),m,t,ldt,df,m,
269 $ df(1,np1),m,work,info)
274 CALL zgemm(
'N',
'C', m, n2, n2, -one, d, m, q, n2, one, df, m )
275 resid = zlange(
'1', m, n2, df, m, rwork )
276 IF( cnorm.GT.zero )
THEN
277 result( 6 ) = resid / (eps*max(1,n2)*dnorm)
284 DEALLOCATE ( a, af, q, r, rwork, work, t, c, d, cf, df)
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
subroutine zherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
ZHERK
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine zlarnv(IDIST, ISEED, N, X)
ZLARNV returns a vector of random numbers from a uniform or normal distribution.
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 zlqt05(M, N, L, NB, RESULT)
ZLQT05
subroutine zgemlqt(SIDE, TRANS, M, N, K, MB, V, LDV, T, LDT, C, LDC, WORK, INFO)
ZGEMLQT
subroutine ztplqt(M, N, L, MB, A, LDA, B, LDB, T, LDT, WORK, INFO)
ZTPLQT
subroutine ztpmlqt(SIDE, TRANS, M, N, K, L, MB, V, LDV, T, LDT, A, LDA, B, LDB, WORK, INFO)
ZTPMLQT