87 INTEGER LWORK, M, N, L, NB, LDT
95 COMPLEX,
ALLOCATABLE :: AF(:,:), Q(:,:),
96 $ R(:,:), RWORK(:), WORK( : ), T(:,:),
97 $ CF(:,:), DF(:,:), A(:,:), C(:,:), D(:,:)
102 parameter( zero = 0.0, one = (1.0,0.0), czero=(0.0,0.0) )
105 INTEGER INFO, J, K, N2, NP1,i
106 REAL ANORM, EPS, RESID, CNORM, DNORM
115 EXTERNAL slamch, clange, clansy, lsame
118 DATA iseed / 1988, 1989, 1990, 1991 /
120 eps = slamch(
'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 claset(
'Full', m, n2, czero, czero, a, m )
140 CALL claset(
'Full', nb, m, czero, czero, t, nb )
142 CALL clarnv( 2, iseed, m-j+1, a( j, j ) )
146 CALL clarnv( 2, iseed, m, a( 1, min(n+m,m+1) + j - 1 ) )
151 CALL clarnv( 2, iseed, m-j+1, a( j, min(n+m,n+m-l+1)
158 CALL clacpy(
'Full', m, n2, a, m, af, m )
162 CALL ctplqt( m,n,l,nb,af,m,af(1,np1),m,t,ldt,work,info)
166 CALL claset(
'Full', n2, n2, czero, one, q, n2 )
167 CALL cgemlqt(
'L',
'N', n2, n2, k, nb, af, m, t, ldt, q, n2,
172 CALL claset(
'Full', n2, n2, czero, czero, r, n2 )
173 CALL clacpy(
'Lower', m, n2, af, m, r, n2 )
177 CALL cgemm(
'N',
'C', m, n2, n2, -one, a, m, q, n2, one, r, n2)
178 anorm = clange(
'1', m, n2, a, m, rwork )
179 resid = clange(
'1', m, n2, r, n2, rwork )
180 IF( anorm.GT.zero )
THEN
181 result( 1 ) = resid / (eps*anorm*max(1,n2))
188 CALL claset(
'Full', n2, n2, czero, one, r, n2 )
189 CALL cherk(
'U',
'N', n2, n2, real(-one), q, n2, real(one),
191 resid = clansy(
'1',
'Upper', n2, r, n2, rwork )
192 result( 2 ) = resid / (eps*max(1,n2))
196 CALL claset(
'Full', n2, m, czero, one, c, n2 )
198 CALL clarnv( 2, iseed, n2, c( 1, j ) )
200 cnorm = clange(
'1', n2, m, c, n2, rwork)
201 CALL clacpy(
'Full', n2, m, c, n2, cf, n2 )
205 CALL ctpmlqt(
'L',
'N', n,m,k,l,nb,af(1, np1),m,t,ldt,cf,n2,
206 $ cf(np1,1),n2,work,info)
210 CALL cgemm(
'N',
'N', n2, m, n2, -one, q, n2, c, n2, one, cf, n2 )
211 resid = clange(
'1', n2, m, cf, n2, rwork )
212 IF( cnorm.GT.zero )
THEN
213 result( 3 ) = resid / (eps*max(1,n2)*cnorm)
221 CALL clacpy(
'Full', n2, m, c, n2, cf, n2 )
225 CALL ctpmlqt(
'L',
'C',n,m,k,l,nb,af(1,np1),m,t,ldt,cf,n2,
226 $ cf(np1,1),n2,work,info)
230 CALL cgemm(
'C',
'N',n2,m,n2,-one,q,n2,c,n2,one,cf,n2)
231 resid = clange(
'1', n2, m, cf, n2, rwork )
233 IF( cnorm.GT.zero )
THEN
234 result( 4 ) = resid / (eps*max(1,n2)*cnorm)
242 CALL clarnv( 2, iseed, m, d( 1, j ) )
244 dnorm = clange(
'1', m, n2, d, m, rwork)
245 CALL clacpy(
'Full', m, n2, d, m, df, m )
249 CALL ctpmlqt(
'R',
'N',m,n,k,l,nb,af(1,np1),m,t,ldt,df,m,
250 $ df(1,np1),m,work,info)
254 CALL cgemm(
'N',
'N',m,n2,n2,-one,d,m,q,n2,one,df,m)
255 resid = clange(
'1',m, n2,df,m,rwork )
256 IF( cnorm.GT.zero )
THEN
257 result( 5 ) = resid / (eps*max(1,n2)*dnorm)
264 CALL clacpy(
'Full',m,n2,d,m,df,m )
268 CALL ctpmlqt(
'R',
'C',m,n,k,l,nb,af(1,np1),m,t,ldt,df,m,
269 $ df(1,np1),m,work,info)
274 CALL cgemm(
'N',
'C', m, n2, n2, -one, d, m, q, n2, one, df, m )
275 resid = clange(
'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 cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
subroutine cherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
CHERK
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 clarnv(IDIST, ISEED, N, X)
CLARNV returns a vector of random numbers from a uniform or normal distribution.
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine clqt05(M, N, L, NB, RESULT)
CLQT05
subroutine cgemlqt(SIDE, TRANS, M, N, K, MB, V, LDV, T, LDT, C, LDC, WORK, INFO)
CGEMLQT
subroutine ctplqt(M, N, L, MB, A, LDA, B, LDB, T, LDT, WORK, INFO)
CTPLQT
subroutine ctpmlqt(SIDE, TRANS, M, N, K, L, MB, V, LDV, T, LDT, A, LDA, B, LDB, WORK, INFO)
CTPMLQT