LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dchkqp3rk.f
Go to the documentation of this file.
1*> \brief \b DCHKQP3RK
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8* Definition:
9* ===========
10*
11* SUBROUTINE DCHKQP3RK( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL,
12* $ NNB, NBVAL, NXVAL, THRESH, A, COPYA,
13* $ B, COPYB, S, TAU,
14* $ WORK, IWORK, NOUT )
15* IMPLICIT NONE
16*
17* -- LAPACK test routine --
18* -- LAPACK is a software package provided by Univ. of Tennessee, --
19* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
20*
21* .. Scalar Arguments ..
22* INTEGER NM, NN, NNS, NNB, NOUT
23* DOUBLE PRECISION THRESH
24* ..
25* .. Array Arguments ..
26* LOGICAL DOTYPE( * )
27* INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
28* $ NVAL( * ), NXVAL( * )
29* DOUBLE PRECISION A( * ), COPYA( * ), B( * ), COPYB( * ),
30* $ S( * ), TAU( * ), WORK( * )
31* ..
32*
33*
34*> \par Purpose:
35* =============
36*>
37*> \verbatim
38*>
39*> DCHKQP3RK tests DGEQP3RK.
40*> \endverbatim
41*
42* Arguments:
43* ==========
44*
45*> \param[in] DOTYPE
46*> \verbatim
47*> DOTYPE is LOGICAL array, dimension (NTYPES)
48*> The matrix types to be used for testing. Matrices of type j
49*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
50*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
51*> \endverbatim
52*>
53*> \param[in] NM
54*> \verbatim
55*> NM is INTEGER
56*> The number of values of M contained in the vector MVAL.
57*> \endverbatim
58*>
59*> \param[in] MVAL
60*> \verbatim
61*> MVAL is INTEGER array, dimension (NM)
62*> The values of the matrix row dimension M.
63*> \endverbatim
64*>
65*> \param[in] NN
66*> \verbatim
67*> NN is INTEGER
68*> The number of values of N contained in the vector NVAL.
69*> \endverbatim
70*>
71*> \param[in] NVAL
72*> \verbatim
73*> NVAL is INTEGER array, dimension (NN)
74*> The values of the matrix column dimension N.
75*> \endverbatim
76*>
77*> \param[in] NNS
78*> \verbatim
79*> NNS is INTEGER
80*> The number of values of NRHS contained in the vector NSVAL.
81*> \endverbatim
82*>
83*> \param[in] NSVAL
84*> \verbatim
85*> NSVAL is INTEGER array, dimension (NNS)
86*> The values of the number of right hand sides NRHS.
87*> \endverbatim
88*>
89*> \param[in] NNB
90*> \verbatim
91*> NNB is INTEGER
92*> The number of values of NB and NX contained in the
93*> vectors NBVAL and NXVAL. The blocking parameters are used
94*> in pairs (NB,NX).
95*> \endverbatim
96*>
97*> \param[in] NBVAL
98*> \verbatim
99*> NBVAL is INTEGER array, dimension (NNB)
100*> The values of the blocksize NB.
101*> \endverbatim
102*>
103*> \param[in] NXVAL
104*> \verbatim
105*> NXVAL is INTEGER array, dimension (NNB)
106*> The values of the crossover point NX.
107*> \endverbatim
108*>
109*> \param[in] THRESH
110*> \verbatim
111*> THRESH is DOUBLE PRECISION
112*> The threshold value for the test ratios. A result is
113*> included in the output file if RESULT >= THRESH. To have
114*> every test ratio printed, use THRESH = 0.
115*> \endverbatim
116*>
117*> \param[out] A
118*> \verbatim
119*> A is DOUBLE PRECISION array, dimension (MMAX*NMAX)
120*> where MMAX is the maximum value of M in MVAL and NMAX is the
121*> maximum value of N in NVAL.
122*> \endverbatim
123*>
124*> \param[out] COPYA
125*> \verbatim
126*> COPYA is DOUBLE PRECISION array, dimension (MMAX*NMAX)
127*> \endverbatim
128*>
129*> \param[out] B
130*> \verbatim
131*> B is DOUBLE PRECISION array, dimension (MMAX*NSMAX)
132*> where MMAX is the maximum value of M in MVAL and NSMAX is the
133*> maximum value of NRHS in NSVAL.
134*> \endverbatim
135*>
136*> \param[out] COPYB
137*> \verbatim
138*> COPYB is DOUBLE PRECISION array, dimension (MMAX*NSMAX)
139*> \endverbatim
140*>
141*> \param[out] S
142*> \verbatim
143*> S is DOUBLE PRECISION array, dimension
144*> (min(MMAX,NMAX))
145*> \endverbatim
146*>
147*> \param[out] TAU
148*> \verbatim
149*> TAU is DOUBLE PRECISION array, dimension (MMAX)
150*> \endverbatim
151*>
152*> \param[out] WORK
153*> \verbatim
154*> WORK is DOUBLE PRECISION array, dimension
155*> (MMAX*NMAX + 4*NMAX + MMAX)
156*> \endverbatim
157*>
158*> \param[out] IWORK
159*> \verbatim
160*> IWORK is INTEGER array, dimension (2*NMAX)
161*> \endverbatim
162*>
163*> \param[in] NOUT
164*> \verbatim
165*> NOUT is INTEGER
166*> The unit number for output.
167*> \endverbatim
168*
169* Authors:
170* ========
171*
172*> \author Univ. of Tennessee
173*> \author Univ. of California Berkeley
174*> \author Univ. of Colorado Denver
175*> \author NAG Ltd.
176*
177*> \ingroup double_lin
178*
179* =====================================================================
180 SUBROUTINE dchkqp3rk( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL,
181 $ NNB, NBVAL, NXVAL, THRESH, A, COPYA,
182 $ B, COPYB, S, TAU,
183 $ WORK, IWORK, NOUT )
184 IMPLICIT NONE
185*
186* -- LAPACK test routine --
187* -- LAPACK is a software package provided by Univ. of Tennessee, --
188* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
189*
190* .. Scalar Arguments ..
191 INTEGER NM, NN, NNB, NNS, NOUT
192 DOUBLE PRECISION THRESH
193* ..
194* .. Array Arguments ..
195 LOGICAL DOTYPE( * )
196 INTEGER IWORK( * ), NBVAL( * ), MVAL( * ), NVAL( * ),
197 $ NSVAL( * ), NXVAL( * )
198 DOUBLE PRECISION A( * ), COPYA( * ), B( * ), COPYB( * ),
199 $ s( * ), tau( * ), work( * )
200* ..
201*
202* =====================================================================
203*
204* .. Parameters ..
205 INTEGER NTYPES
206 PARAMETER ( NTYPES = 19 )
207 integer ntests
208 parameter( ntests = 5 )
209 DOUBLE PRECISION ONE, ZERO, BIGNUM
210 parameter( one = 1.0d+0, zero = 0.0d+0,
211 $ bignum = 1.0d+38 )
212* ..
213* .. Local Scalars ..
214 CHARACTER DIST, TYPE
215 CHARACTER*3 PATH
216 INTEGER I, IHIGH, ILOW, IM, IMAT, IN, INC_ZERO,
217 $ inb, ind_offset_gen,
218 $ ind_in, ind_out, ins, info,
219 $ istep, j, j_inc, j_first_nz, jb_zero,
220 $ kfact, kl, kmax, ku, lda, lw, lwork,
221 $ lwork_mqr, m, minmn, minmnb_gen, mode, n,
222 $ nb, nb_zero, nerrs, nfail, nb_gen, nrhs,
223 $ nrun, nx, t
224 DOUBLE PRECISION ANORM, CNDNUM, EPS, ABSTOL, RELTOL,
225 $ DTEMP, MAXC2NRMK, RELMAXC2NRMK
226* ..
227* .. Local Arrays ..
228 INTEGER ISEED( 4 ), ISEEDY( 4 )
229 DOUBLE PRECISION RESULT( NTESTS ), RDUMMY( 1 )
230* ..
231* .. External Functions ..
232 DOUBLE PRECISION DLAMCH, DQPT01, DQRT11, DQRT12, DLANGE,
233 $ DLAPY2
234 EXTERNAL DLAMCH, DQPT01, DQRT11, DQRT12, DLANGE
235* ..
236* .. External Subroutines ..
237 EXTERNAL alaerh, alahd, alasum, daxpy, dgeqp3rk,
240* ..
241* .. Intrinsic Functions ..
242 INTRINSIC abs, dble, max, min, mod
243* ..
244* .. Scalars in Common ..
245 LOGICAL LERR, OK
246 CHARACTER*32 SRNAMT
247 INTEGER INFOT, IOUNIT
248* ..
249* .. Common blocks ..
250 COMMON / infoc / infot, iounit, ok, lerr
251 COMMON / srnamc / srnamt
252* ..
253* .. Data statements ..
254 DATA iseedy / 1988, 1989, 1990, 1991 /
255* ..
256* .. Executable Statements ..
257*
258* Initialize constants and the random number seed.
259*
260 path( 1: 1 ) = 'Double precision'
261 path( 2: 3 ) = 'QK'
262 nrun = 0
263 nfail = 0
264 nerrs = 0
265 DO i = 1, 4
266 iseed( i ) = iseedy( i )
267 END DO
268 eps = dlamch( 'Epsilon' )
269 infot = 0
270*
271 DO im = 1, nm
272*
273* Do for each value of M in MVAL.
274*
275 m = mval( im )
276 lda = max( 1, m )
277*
278 DO in = 1, nn
279*
280* Do for each value of N in NVAL.
281*
282 n = nval( in )
283 minmn = min( m, n )
284 lwork = max( 1, m*max( m, n )+4*minmn+max( m, n ),
285 $ m*n + 2*minmn + 4*n )
286*
287 DO ins = 1, nns
288 nrhs = nsval( ins )
289*
290* Set up parameters with DLATB4 and generate
291* M-by-NRHS B matrix with DLATMS.
292* IMAT = 14:
293* Random matrix, CNDNUM = 2, NORM = ONE,
294* MODE = 3 (geometric distribution of singular values).
295*
296 CALL dlatb4( path, 14, m, nrhs, TYPE, kl, ku, anorm,
297 $ mode, cndnum, dist )
298*
299 srnamt = 'DLATMS'
300 CALL dlatms( m, nrhs, dist, iseed, TYPE, s, mode,
301 $ cndnum, anorm, kl, ku, 'No packing',
302 $ copyb, lda, work, info )
303
304
305*
306* Check error code from DLATMS.
307*
308 IF( info.NE.0 ) THEN
309 CALL alaerh( path, 'DLATMS', info, 0, ' ', m,
310 $ nrhs, -1, -1, -1, 6, nfail, nerrs,
311 $ nout )
312 cycle
313 END IF
314*
315 DO imat = 1, ntypes
316*
317* Do the tests only if DOTYPE( IMAT ) is true.
318*
319 IF( .NOT.dotype( imat ) )
320 $ cycle
321*
322* The type of distribution used to generate the random
323* eigen-/singular values:
324* ( 'S' for symmetric distribution ) => UNIFORM( -1, 1 )
325*
326* Do for each type of NON-SYMMETRIC matrix: CNDNUM NORM MODE
327* 1. Zero matrix
328* 2. Random, Diagonal, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
329* 3. Random, Upper triangular, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
330* 4. Random, Lower triangular, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
331* 5. Random, First column is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
332* 6. Random, Last MINMN column is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
333* 7. Random, Last N column is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
334* 8. Random, Middle column in MINMN is zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
335* 9. Random, First half of MINMN columns are zero, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
336* 10. Random, Last columns are zero starting from MINMN/2+1, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
337* 11. Random, Half MINMN columns in the middle are zero starting
338* from MINMN/2-(MINMN/2)/2+1, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
339* 12. Random, Odd columns are ZERO, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
340* 13. Random, Even columns are ZERO, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
341* 14. Random, CNDNUM = 2 CNDNUM = 2 ONE 3 ( geometric distribution of singular values )
342* 15. Random, CNDNUM = sqrt(0.1/EPS) CNDNUM = BADC1 = sqrt(0.1/EPS) ONE 3 ( geometric distribution of singular values )
343* 16. Random, CNDNUM = 0.1/EPS CNDNUM = BADC2 = 0.1/EPS ONE 3 ( geometric distribution of singular values )
344* 17. Random, CNDNUM = 0.1/EPS, CNDNUM = BADC2 = 0.1/EPS ONE 2 ( one small singular value, S(N)=1/CNDNUM )
345* one small singular value S(N)=1/CNDNUM
346* 18. Random, CNDNUM = 2, scaled near underflow CNDNUM = 2 SMALL = SAFMIN
347* 19. Random, CNDNUM = 2, scaled near overflow CNDNUM = 2 LARGE = 1.0/( 0.25 * ( SAFMIN / EPS ) ) 3 ( geometric distribution of singular values )
348*
349 IF( imat.EQ.1 ) THEN
350*
351* Matrix 1: Zero matrix
352*
353 CALL dlaset( 'Full', m, n, zero, zero, copya, lda )
354 DO i = 1, minmn
355 s( i ) = zero
356 END DO
357*
358 ELSE IF( (imat.GE.2 .AND. imat.LE.4 )
359 $ .OR. (imat.GE.14 .AND. imat.LE.19 ) ) THEN
360*
361* Matrices 2-5.
362*
363* Set up parameters with DLATB4 and generate a test
364* matrix with DLATMS.
365*
366 CALL dlatb4( path, imat, m, n, TYPE, kl, ku, anorm,
367 $ mode, cndnum, dist )
368*
369 srnamt = 'DLATMS'
370 CALL dlatms( m, n, dist, iseed, TYPE, s, mode,
371 $ cndnum, anorm, kl, ku, 'No packing',
372 $ copya, lda, work, info )
373*
374* Check error code from DLATMS.
375*
376 IF( info.NE.0 ) THEN
377 CALL alaerh( path, 'DLATMS', info, 0, ' ', m, n,
378 $ -1, -1, -1, imat, nfail, nerrs,
379 $ nout )
380 cycle
381 END IF
382*
383 CALL dlaord( 'Decreasing', minmn, s, 1 )
384*
385 ELSE IF( minmn.GE.2
386 $ .AND. imat.GE.5 .AND. imat.LE.13 ) THEN
387*
388* Rectangular matrices 5-13 that contain zero columns,
389* only for matrices MINMN >=2.
390*
391* JB_ZERO is the column index of ZERO block.
392* NB_ZERO is the column block size of ZERO block.
393* NB_GEN is the column blcok size of the
394* generated block.
395* J_INC in the non_zero column index increment
396* for matrix 12 and 13.
397* J_FIRS_NZ is the index of the first non-zero
398* column.
399*
400 IF( imat.EQ.5 ) THEN
401*
402* First column is zero.
403*
404 jb_zero = 1
405 nb_zero = 1
406 nb_gen = n - nb_zero
407*
408 ELSE IF( imat.EQ.6 ) THEN
409*
410* Last column MINMN is zero.
411*
412 jb_zero = minmn
413 nb_zero = 1
414 nb_gen = n - nb_zero
415*
416 ELSE IF( imat.EQ.7 ) THEN
417*
418* Last column N is zero.
419*
420 jb_zero = n
421 nb_zero = 1
422 nb_gen = n - nb_zero
423*
424 ELSE IF( imat.EQ.8 ) THEN
425*
426* Middle column in MINMN is zero.
427*
428 jb_zero = minmn / 2 + 1
429 nb_zero = 1
430 nb_gen = n - nb_zero
431*
432 ELSE IF( imat.EQ.9 ) THEN
433*
434* First half of MINMN columns is zero.
435*
436 jb_zero = 1
437 nb_zero = minmn / 2
438 nb_gen = n - nb_zero
439*
440 ELSE IF( imat.EQ.10 ) THEN
441*
442* Last columns are zero columns,
443* starting from (MINMN / 2 + 1) column.
444*
445 jb_zero = minmn / 2 + 1
446 nb_zero = n - jb_zero + 1
447 nb_gen = n - nb_zero
448*
449 ELSE IF( imat.EQ.11 ) THEN
450*
451* Half of the columns in the middle of MINMN
452* columns is zero, starting from
453* MINMN/2 - (MINMN/2)/2 + 1 column.
454*
455 jb_zero = minmn / 2 - (minmn / 2) / 2 + 1
456 nb_zero = minmn / 2
457 nb_gen = n - nb_zero
458*
459 ELSE IF( imat.EQ.12 ) THEN
460*
461* Odd-numbered columns are zero,
462*
463 nb_gen = n / 2
464 nb_zero = n - nb_gen
465 j_inc = 2
466 j_first_nz = 2
467*
468 ELSE IF( imat.EQ.13 ) THEN
469*
470* Even-numbered columns are zero.
471*
472 nb_zero = n / 2
473 nb_gen = n - nb_zero
474 j_inc = 2
475 j_first_nz = 1
476*
477 END IF
478*
479*
480* 1) Set the first NB_ZERO columns in COPYA(1:M,1:N)
481* to zero.
482*
483 CALL dlaset( 'Full', m, nb_zero, zero, zero,
484 $ copya, lda )
485*
486* 2) Generate an M-by-(N-NB_ZERO) matrix with the
487* chosen singular value distribution
488* in COPYA(1:M,NB_ZERO+1:N).
489*
490 CALL dlatb4( path, imat, m, nb_gen, TYPE, kl, ku,
491 $ anorm, mode, cndnum, dist )
492*
493 srnamt = 'DLATMS'
494*
495 ind_offset_gen = nb_zero * lda
496*
497 CALL dlatms( m, nb_gen, dist, iseed, TYPE, s, mode,
498 $ cndnum, anorm, kl, ku, 'No packing',
499 $ copya( ind_offset_gen + 1 ), lda,
500 $ work, info )
501*
502* Check error code from DLATMS.
503*
504 IF( info.NE.0 ) THEN
505 CALL alaerh( path, 'DLATMS', info, 0, ' ', m,
506 $ nb_gen, -1, -1, -1, imat, nfail,
507 $ nerrs, nout )
508 cycle
509 END IF
510*
511* 3) Swap the gererated colums from the right side
512* NB_GEN-size block in COPYA into correct column
513* positions.
514*
515 IF( imat.EQ.6
516 $ .OR. imat.EQ.7
517 $ .OR. imat.EQ.8
518 $ .OR. imat.EQ.10
519 $ .OR. imat.EQ.11 ) THEN
520*
521* Move by swapping the generated columns
522* from the right NB_GEN-size block from
523* (NB_ZERO+1:NB_ZERO+JB_ZERO)
524* into columns (1:JB_ZERO-1).
525*
526 DO j = 1, jb_zero-1, 1
527 CALL dswap( m,
528 $ copya( ( nb_zero+j-1)*lda+1), 1,
529 $ copya( (j-1)*lda + 1 ), 1 )
530 END DO
531*
532 ELSE IF( imat.EQ.12 .OR. imat.EQ.13 ) THEN
533*
534* ( IMAT = 12, Odd-numbered ZERO columns. )
535* Swap the generated columns from the right
536* NB_GEN-size block into the even zero colums in the
537* left NB_ZERO-size block.
538*
539* ( IMAT = 13, Even-numbered ZERO columns. )
540* Swap the generated columns from the right
541* NB_GEN-size block into the odd zero colums in the
542* left NB_ZERO-size block.
543*
544 DO j = 1, nb_gen, 1
545 ind_out = ( nb_zero+j-1 )*lda + 1
546 ind_in = ( j_inc*(j-1)+(j_first_nz-1) )*lda
547 $ + 1
548 CALL dswap( m,
549 $ copya( ind_out ), 1,
550 $ copya( ind_in), 1 )
551 END DO
552*
553 END IF
554*
555* 5) Order the singular values generated by
556* DLAMTS in decreasing order and add trailing zeros
557* that correspond to zero columns.
558* The total number of singular values is MINMN.
559*
560 minmnb_gen = min( m, nb_gen )
561*
562 DO i = minmnb_gen+1, minmn
563 s( i ) = zero
564 END DO
565*
566 ELSE
567*
568* IF(MINMN.LT.2) skip this size for this matrix type.
569*
570 cycle
571 END IF
572*
573* Initialize a copy array for a pivot array for DGEQP3RK.
574*
575 DO i = 1, n
576 iwork( i ) = 0
577 END DO
578*
579 DO inb = 1, nnb
580*
581* Do for each pair of values (NB,NX) in NBVAL and NXVAL.
582*
583 nb = nbval( inb )
584 CALL xlaenv( 1, nb )
585 nx = nxval( inb )
586 CALL xlaenv( 3, nx )
587*
588* We do MIN(M,N)+1 because we need a test for KMAX > N,
589* when KMAX is larger than MIN(M,N), KMAX should be
590* KMAX = MIN(M,N)
591*
592 DO kmax = 0, min(m,n)+1
593*
594* Get a working copy of COPYA into A( 1:M,1:N ).
595* Get a working copy of COPYB into A( 1:M, (N+1):NRHS ).
596* Get a working copy of COPYB into into B( 1:M, 1:NRHS ).
597* Get a working copy of IWORK(1:N) awith zeroes into
598* which is going to be used as pivot array IWORK( N+1:2N ).
599* NOTE: IWORK(2N+1:3N) is going to be used as a WORK array
600* for the routine.
601*
602 CALL dlacpy( 'All', m, n, copya, lda, a, lda )
603 CALL dlacpy( 'All', m, nrhs, copyb, lda,
604 $ a( lda*n + 1 ), lda )
605 CALL dlacpy( 'All', m, nrhs, copyb, lda,
606 $ b, lda )
607 CALL icopy( n, iwork( 1 ), 1, iwork( n+1 ), 1 )
608*
609 abstol = -1.0
610 reltol = -1.0
611*
612* Compute the QR factorization with pivoting of A
613*
614 lw = max( 1, max( 2*n + nb*( n+nrhs+1 ),
615 $ 3*n + nrhs - 1 ) )
616*
617* Compute DGEQP3RK factorization of A.
618*
619 srnamt = 'DGEQP3RK'
620 CALL dgeqp3rk( m, n, nrhs, kmax, abstol, reltol,
621 $ a, lda, kfact, maxc2nrmk,
622 $ relmaxc2nrmk, iwork( n+1 ), tau,
623 $ work, lw, iwork( 2*n+1 ), info )
624*
625* Check error code from DGEQP3RK.
626*
627 IF( info.LT.0 )
628 $ CALL alaerh( path, 'DGEQP3RK', info, 0, ' ',
629 $ m, n, nx, -1, nb, imat,
630 $ nfail, nerrs, nout )
631*
632* Compute test 1:
633*
634* This test in only for the full rank factorization of
635* the matrix A.
636*
637* Array S(1:min(M,N)) contains svd(A) the sigular values
638* of the original matrix A in decreasing absolute value
639* order. The test computes svd(R), the vector sigular
640* values of the upper trapezoid of A(1:M,1:N) that
641* contains the factor R, in decreasing order. The test
642* returns the ratio:
643*
644* 2-norm(svd(R) - svd(A)) / ( max(M,N) * 2-norm(svd(A)) * EPS )
645*
646 IF( kfact.EQ.minmn ) THEN
647*
648 result( 1 ) = dqrt12( m, n, a, lda, s, work,
649 $ lwork )
650*
651 DO t = 1, 1
652 IF( result( t ).GE.thresh ) THEN
653 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
654 $ CALL alahd( nout, path )
655 WRITE( nout, fmt = 9999 ) 'DGEQP3RK', m, n,
656 $ nrhs, kmax, abstol, reltol, nb, nx,
657 $ imat, t, result( t )
658 nfail = nfail + 1
659 END IF
660 END DO
661 nrun = nrun + 1
662*
663* End test 1
664*
665 END IF
666*
667* Compute test 2:
668*
669* The test returns the ratio:
670*
671* 1-norm( A*P - Q*R ) / ( max(M,N) * 1-norm(A) * EPS )
672*
673 result( 2 ) = dqpt01( m, n, kfact, copya, a, lda, tau,
674 $ iwork( n+1 ), work, lwork )
675*
676* Compute test 3:
677*
678* The test returns the ratio:
679*
680* 1-norm( Q**T * Q - I ) / ( M * EPS )
681*
682 result( 3 ) = dqrt11( m, kfact, a, lda, tau, work,
683 $ lwork )
684*
685* Print information about the tests that did not pass
686* the threshold.
687*
688 DO t = 2, 3
689 IF( result( t ).GE.thresh ) THEN
690 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
691 $ CALL alahd( nout, path )
692 WRITE( nout, fmt = 9999 ) 'DGEQP3RK', m, n,
693 $ nrhs, kmax, abstol, reltol,
694 $ nb, nx, imat, t, result( t )
695 nfail = nfail + 1
696 END IF
697 END DO
698 nrun = nrun + 2
699*
700* Compute test 4:
701*
702* This test is only for the factorizations with the
703* rank greater than 2.
704* The elements on the diagonal of R should be non-
705* increasing.
706*
707* The test returns the ratio:
708*
709* Returns 1.0D+100 if abs(R(K+1,K+1)) > abs(R(K,K)),
710* K=1:KFACT-1
711*
712 IF( min(kfact, minmn).GE.2 ) THEN
713*
714 DO j = 1, kfact-1, 1
715
716 dtemp = (( abs( a( (j-1)*m+j ) ) -
717 $ abs( a( (j)*m+j+1 ) ) ) /
718 $ abs( a(1) ) )
719*
720 IF( dtemp.LT.zero ) THEN
721 result( 4 ) = bignum
722 END IF
723*
724 END DO
725*
726* Print information about the tests that did not
727* pass the threshold.
728*
729 DO t = 4, 4
730 IF( result( t ).GE.thresh ) THEN
731 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
732 $ CALL alahd( nout, path )
733 WRITE( nout, fmt = 9999 ) 'DGEQP3RK',
734 $ m, n, nrhs, kmax, abstol, reltol,
735 $ nb, nx, imat, t,
736 $ result( t )
737 nfail = nfail + 1
738 END IF
739 END DO
740 nrun = nrun + 1
741*
742* End test 4.
743*
744 END IF
745*
746* Compute test 5:
747*
748* This test in only for matrix A with min(M,N) > 0.
749*
750* The test returns the ratio:
751*
752* 1-norm(Q**T * B - Q**T * B ) /
753* ( M * EPS )
754*
755* (1) Compute B:=Q**T * B in the matrix B.
756*
757 IF( minmn.GT.0 ) THEN
758*
759 lwork_mqr = max(1, nrhs)
760 CALL dormqr( 'Left', 'Transpose',
761 $ m, nrhs, kfact, a, lda, tau, b, lda,
762 $ work, lwork_mqr, info )
763*
764 DO i = 1, nrhs
765*
766* Compare N+J-th column of A and J-column of B.
767*
768 CALL daxpy( m, -one, a( ( n+i-1 )*lda+1 ), 1,
769 $ b( ( i-1 )*lda+1 ), 1 )
770 END DO
771*
772 result( 5 ) =
773 $ abs(
774 $ dlange( 'One-norm', m, nrhs, b, lda, rdummy ) /
775 $ ( dble( m )*dlamch( 'Epsilon' ) )
776 $ )
777*
778* Print information about the tests that did not pass
779* the threshold.
780*
781 DO t = 5, 5
782 IF( result( t ).GE.thresh ) THEN
783 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
784 $ CALL alahd( nout, path )
785 WRITE( nout, fmt = 9999 ) 'DGEQP3RK', m, n,
786 $ nrhs, kmax, abstol, reltol,
787 $ nb, nx, imat, t, result( t )
788 nfail = nfail + 1
789 END IF
790 END DO
791 nrun = nrun + 1
792*
793* End compute test 5.
794*
795 END IF
796*
797* END DO KMAX = 1, MIN(M,N)+1
798*
799 END DO
800*
801* END DO for INB = 1, NNB
802*
803 END DO
804*
805* END DO for IMAT = 1, NTYPES
806*
807 END DO
808*
809* END DO for INS = 1, NNS
810*
811 END DO
812*
813* END DO for IN = 1, NN
814*
815 END DO
816*
817* END DO for IM = 1, NM
818*
819 END DO
820*
821* Print a summary of the results.
822*
823 CALL alasum( path, nout, nfail, nrun, nerrs )
824*
825 9999 FORMAT( 1x, a, ' M =', i5, ', N =', i5, ', NRHS =', i5,
826 $ ', KMAX =', i5, ', ABSTOL =', g12.5,
827 $ ', RELTOL =', g12.5, ', NB =', i4, ', NX =', i4,
828 $ ', type ', i2, ', test ', i2, ', ratio =', g12.5 )
829*
830* End of DCHKQP3RK
831*
832 END
subroutine alasum(type, nout, nfail, nrun, nerrs)
ALASUM
Definition alasum.f:73
subroutine xlaenv(ispec, nvalue)
XLAENV
Definition xlaenv.f:81
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
Definition alaerh.f:147
subroutine alahd(iounit, path)
ALAHD
Definition alahd.f:107
subroutine dchkqp3rk(dotype, nm, mval, nn, nval, nns, nsval, nnb, nbval, nxval, thresh, a, copya, b, copyb, s, tau, work, iwork, nout)
DCHKQP3RK
Definition dchkqp3rk.f:184
subroutine dgeqp3rk(m, n, nrhs, kmax, abstol, reltol, a, lda, k, maxc2nrmk, relmaxc2nrmk, jpiv, tau, work, lwork, iwork, info)
DGEQP3RK computes a truncated Householder QR factorization with column pivoting of a real m-by-n matr...
Definition dgeqp3rk.f:585
subroutine dlaord(job, n, x, incx)
DLAORD
Definition dlaord.f:73
subroutine dlatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
DLATB4
Definition dlatb4.f:120
subroutine dlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
DLATMS
Definition dlatms.f:321
subroutine daxpy(n, da, dx, incx, dy, incy)
DAXPY
Definition daxpy.f:89
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:103
subroutine dlaset(uplo, m, n, alpha, beta, a, lda)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition dlaset.f:110
subroutine dswap(n, dx, incx, dy, incy)
DSWAP
Definition dswap.f:82
subroutine dormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
DORMQR
Definition dormqr.f:167
subroutine icopy(n, sx, incx, sy, incy)
ICOPY
Definition icopy.f:75