83 SUBROUTINE ctsqr01(TSSW, M, N, MB, NB, RESULT)
101 COMPLEX,
ALLOCATABLE :: AF(:,:), Q(:,:),
102 $ r(:,:), rwork(:), work( : ), t(:),
103 $ cf(:,:), df(:,:), a(:,:), c(:,:), d(:,:), lq(:,:)
108 parameter( zero = 0.0, one = (1.0,0.0), czero=(0.0,0.0) )
111 LOGICAL TESTZEROS, TS
112 INTEGER INFO, J, K, L, LWORK, TSIZE, MNB
113 REAL ANORM, EPS, RESID, CNORM, DNORM
117 COMPLEX TQUERY( 5 ), WORKQUERY
120 REAL SLAMCH, CLANGE, CLANSY
123 EXTERNAL slamch, clange, clansy, lsame, ilaenv
131 COMMON / srnamc / srnamt
134 DATA iseed / 1988, 1989, 1990, 1991 /
138 ts = lsame(tssw,
'TS')
144 eps = slamch(
'Epsilon' )
152 ALLOCATE ( a(m,n), af(m,n), q(l,l), r(m,l), rwork(l),
154 $ d(n,m), df(n,m), lq(l,n) )
159 CALL clarnv( 2, iseed, m, a( 1, j ) )
164 CALL clarnv( 2, iseed, m/2, a( m/4, j ) )
168 CALL clacpy(
'Full', m, n, a, m, af, m )
174 CALL cgeqr( m, n, af, m, tquery, -1, workquery, -1, info )
175 tsize = int( tquery( 1 ) )
176 lwork = int( workquery )
177 CALL cgemqr(
'L',
'N', m, m, k, af, m, tquery, tsize, cf, m,
178 $ workquery, -1, info)
179 lwork = max( lwork, int( workquery ) )
180 CALL cgemqr(
'L',
'N', m, n, k, af, m, tquery, tsize, cf, m,
181 $ workquery, -1, info)
182 lwork = max( lwork, int( workquery ) )
183 CALL cgemqr(
'L',
'C', m, n, k, af, m, tquery, tsize, cf, m,
184 $ workquery, -1, info)
185 lwork = max( lwork, int( workquery ) )
186 CALL cgemqr(
'R',
'N', n, m, k, af, m, tquery, tsize, df, n,
187 $ workquery, -1, info)
188 lwork = max( lwork, int( workquery ) )
189 CALL cgemqr(
'R',
'C', n, m, k, af, m, tquery, tsize, df, n,
190 $ workquery, -1, info)
191 lwork = max( lwork, int( workquery ) )
192 ALLOCATE ( t( tsize ) )
193 ALLOCATE ( work( lwork ) )
195 CALL cgeqr( m, n, af, m, t, tsize, work, lwork, info )
199 CALL claset(
'Full', m, m, czero, one, q, m )
201 CALL cgemqr(
'L',
'N', m, m, k, af, m, t, tsize, q, m,
202 $ work, lwork, info )
206 CALL claset(
'Full', m, n, czero, czero, r, m )
207 CALL clacpy(
'Upper', m, n, af, m, r, m )
211 CALL cgemm(
'C',
'N', m, n, m, -one, q, m, a, m, one, r, m )
212 anorm = clange(
'1', m, n, a, m, rwork )
213 resid = clange(
'1', m, n, r, m, rwork )
214 IF( anorm.GT.zero )
THEN
215 result( 1 ) = resid / (eps*max(1,m)*anorm)
222 CALL claset(
'Full', m, m, czero, one, r, m )
223 CALL cherk(
'U',
'C', m, m,
REAL(-ONE), Q, M,
REAL(ONE), R, M )
224 resid = clansy(
'1',
'Upper', m, r, m, rwork )
225 result( 2 ) = resid / (eps*max(1,m))
230 CALL clarnv( 2, iseed, m, c( 1, j ) )
232 cnorm = clange(
'1', m, n, c, m, rwork)
233 CALL clacpy(
'Full', m, n, c, m, cf, m )
238 CALL cgemqr(
'L',
'N', m, n, k, af, m, t, tsize, cf, m,
243 CALL cgemm(
'N',
'N', m, n, m, -one, q, m, c, m, one, cf, m )
244 resid = clange(
'1', m, n, cf, m, rwork )
245 IF( cnorm.GT.zero )
THEN
246 result( 3 ) = resid / (eps*max(1,m)*cnorm)
253 CALL clacpy(
'Full', m, n, c, m, cf, m )
258 CALL cgemqr(
'L',
'C', m, n, k, af, m, t, tsize, cf, m,
263 CALL cgemm(
'C',
'N', m, n, m, -one, q, m, c, m, one, cf, m )
264 resid = clange(
'1', m, n, cf, m, rwork )
265 IF( cnorm.GT.zero )
THEN
266 result( 4 ) = resid / (eps*max(1,m)*cnorm)
274 CALL clarnv( 2, iseed, n, d( 1, j ) )
276 dnorm = clange(
'1', n, m, d, n, rwork)
277 CALL clacpy(
'Full', n, m, d, n, df, n )
282 CALL cgemqr(
'R',
'N', n, m, k, af, m, t, tsize, df, n,
287 CALL cgemm(
'N',
'N', n, m, m, -one, d, n, q, m, one, df, n )
288 resid = clange(
'1', n, m, df, n, rwork )
289 IF( dnorm.GT.zero )
THEN
290 result( 5 ) = resid / (eps*max(1,m)*dnorm)
297 CALL clacpy(
'Full', n, m, d, n, df, n )
301 CALL cgemqr(
'R',
'C', n, m, k, af, m, t, tsize, df, n,
306 CALL cgemm(
'N',
'C', n, m, m, -one, d, n, q, m, one, df, n )
307 resid = clange(
'1', n, m, df, n, rwork )
308 IF( cnorm.GT.zero )
THEN
309 result( 6 ) = resid / (eps*max(1,m)*dnorm)
317 CALL cgelq( m, n, af, m, tquery, -1, workquery, -1, info )
318 tsize = int( tquery( 1 ) )
319 lwork = int( workquery )
320 CALL cgemlq(
'R',
'N', n, n, k, af, m, tquery, tsize, q, n,
321 $ workquery, -1, info )
322 lwork = max( lwork, int( workquery ) )
323 CALL cgemlq(
'L',
'N', n, m, k, af, m, tquery, tsize, df, n,
324 $ workquery, -1, info)
325 lwork = max( lwork, int( workquery ) )
326 CALL cgemlq(
'L',
'C', n, m, k, af, m, tquery, tsize, df, n,
327 $ workquery, -1, info)
328 lwork = max( lwork, int( workquery ) )
329 CALL cgemlq(
'R',
'N', m, n, k, af, m, tquery, tsize, cf, m,
330 $ workquery, -1, info)
331 lwork = max( lwork, int( workquery ) )
332 CALL cgemlq(
'R',
'C', m, n, k, af, m, tquery, tsize, cf, m,
333 $ workquery, -1, info)
334 lwork = max( lwork, int( workquery ) )
335 ALLOCATE ( t( tsize ) )
336 ALLOCATE ( work( lwork ) )
338 CALL cgelq( m, n, af, m, t, tsize, work, lwork, info )
343 CALL claset(
'Full', n, n, czero, one, q, n )
345 CALL cgemlq(
'R',
'N', n, n, k, af, m, t, tsize, q, n,
346 $ work, lwork, info )
350 CALL claset(
'Full', m, n, czero, czero, lq, l )
351 CALL clacpy(
'Lower', m, n, af, m, lq, l )
355 CALL cgemm(
'N',
'C', m, n, n, -one, a, m, q, n, one, lq, l )
356 anorm = clange(
'1', m, n, a, m, rwork )
357 resid = clange(
'1', m, n, lq, l, rwork )
358 IF( anorm.GT.zero )
THEN
359 result( 1 ) = resid / (eps*max(1,n)*anorm)
366 CALL claset(
'Full', n, n, czero, one, lq, l )
367 CALL cherk(
'U',
'C', n, n,
REAL(-ONE), Q, N,
REAL(ONE), LQ, L)
368 resid = clansy(
'1',
'Upper', n, lq, l, rwork )
369 result( 2 ) = resid / (eps*max(1,n))
374 CALL clarnv( 2, iseed, n, d( 1, j ) )
376 dnorm = clange(
'1', n, m, d, n, rwork)
377 CALL clacpy(
'Full', n, m, d, n, df, n )
381 CALL cgemlq(
'L',
'N', n, m, k, af, m, t, tsize, df, n,
386 CALL cgemm(
'N',
'N', n, m, n, -one, q, n, d, n, one, df, n )
387 resid = clange(
'1', n, m, df, n, rwork )
388 IF( dnorm.GT.zero )
THEN
389 result( 3 ) = resid / (eps*max(1,n)*dnorm)
396 CALL clacpy(
'Full', n, m, d, n, df, n )
400 CALL cgemlq(
'L',
'C', n, m, k, af, m, t, tsize, df, n,
405 CALL cgemm(
'C',
'N', n, m, n, -one, q, n, d, n, one, df, n )
406 resid = clange(
'1', n, m, df, n, rwork )
407 IF( dnorm.GT.zero )
THEN
408 result( 4 ) = resid / (eps*max(1,n)*dnorm)
416 CALL clarnv( 2, iseed, m, c( 1, j ) )
418 cnorm = clange(
'1', m, n, c, m, rwork)
419 CALL clacpy(
'Full', m, n, c, m, cf, m )
423 CALL cgemlq(
'R',
'N', m, n, k, af, m, t, tsize, cf, m,
428 CALL cgemm(
'N',
'N', m, n, n, -one, c, m, q, n, one, cf, m )
429 resid = clange(
'1', n, m, df, n, rwork )
430 IF( cnorm.GT.zero )
THEN
431 result( 5 ) = resid / (eps*max(1,n)*cnorm)
438 CALL clacpy(
'Full', m, n, c, m, cf, m )
442 CALL cgemlq(
'R',
'C', m, n, k, af, m, t, tsize, cf, m,
447 CALL cgemm(
'N',
'C', m, n, n, -one, c, m, q, n, one, cf, m )
448 resid = clange(
'1', m, n, cf, m, rwork )
449 IF( cnorm.GT.zero )
THEN
450 result( 6 ) = resid / (eps*max(1,n)*cnorm)
459 DEALLOCATE ( a, af, q, r, rwork, work, t, c, d, cf, df)
subroutine cgemqr(SIDE, TRANS, M, N, K, A, LDA, T, TSIZE, C, LDC, WORK, LWORK, INFO)
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 cgelq(M, N, A, LDA, T, TSIZE, WORK, LWORK, INFO)
subroutine clarnv(IDIST, ISEED, N, X)
CLARNV returns a vector of random numbers from a uniform or normal distribution.
subroutine cherk(UPLO, TRANS, N, K, ALPHA, A, LDA, BETA, C, LDC)
CHERK
subroutine cgemlq(SIDE, TRANS, M, N, K, A, LDA, T, TSIZE, C, LDC, WORK, LWORK, INFO)
subroutine cgeqr(M, N, A, LDA, T, TSIZE, WORK, LWORK, INFO)
subroutine ctsqr01(TSSW, M, N, MB, NB, RESULT)
CTSQR01
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