LAPACK  3.7.0
LAPACK: Linear Algebra PACKage
zdrvls.f
Go to the documentation of this file.
1 *> \brief \b ZDRVLS
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 ZDRVLS( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB,
12 * NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B,
13 * COPYB, C, S, COPYS, NOUT )
14 *
15 * .. Scalar Arguments ..
16 * LOGICAL TSTERR
17 * INTEGER NM, NN, NNB, NNS, NOUT
18 * DOUBLE PRECISION THRESH
19 * ..
20 * .. Array Arguments ..
21 * LOGICAL DOTYPE( * )
22 * INTEGER MVAL( * ), NBVAL( * ), NSVAL( * ),
23 * $ NVAL( * ), NXVAL( * )
24 * DOUBLE PRECISION COPYS( * ), S( * )
25 * COMPLEX*16 A( * ), B( * ), C( * ), COPYA( * ), COPYB( * )
26 * ..
27 *
28 *
29 *> \par Purpose:
30 * =============
31 *>
32 *> \verbatim
33 *>
34 *> ZDRVLS tests the least squares driver routines ZGELS, ZGETSLS, ZGELSS, ZGELSY
35 *> and ZGELSD.
36 *> \endverbatim
37 *
38 * Arguments:
39 * ==========
40 *
41 *> \param[in] DOTYPE
42 *> \verbatim
43 *> DOTYPE is LOGICAL array, dimension (NTYPES)
44 *> The matrix types to be used for testing. Matrices of type j
45 *> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
46 *> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
47 *> The matrix of type j is generated as follows:
48 *> j=1: A = U*D*V where U and V are random unitary matrices
49 *> and D has random entries (> 0.1) taken from a uniform
50 *> distribution (0,1). A is full rank.
51 *> j=2: The same of 1, but A is scaled up.
52 *> j=3: The same of 1, but A is scaled down.
53 *> j=4: A = U*D*V where U and V are random unitary matrices
54 *> and D has 3*min(M,N)/4 random entries (> 0.1) taken
55 *> from a uniform distribution (0,1) and the remaining
56 *> entries set to 0. A is rank-deficient.
57 *> j=5: The same of 4, but A is scaled up.
58 *> j=6: The same of 5, but A is scaled down.
59 *> \endverbatim
60 *>
61 *> \param[in] NM
62 *> \verbatim
63 *> NM is INTEGER
64 *> The number of values of M contained in the vector MVAL.
65 *> \endverbatim
66 *>
67 *> \param[in] MVAL
68 *> \verbatim
69 *> MVAL is INTEGER array, dimension (NM)
70 *> The values of the matrix row dimension M.
71 *> \endverbatim
72 *>
73 *> \param[in] NN
74 *> \verbatim
75 *> NN is INTEGER
76 *> The number of values of N contained in the vector NVAL.
77 *> \endverbatim
78 *>
79 *> \param[in] NVAL
80 *> \verbatim
81 *> NVAL is INTEGER array, dimension (NN)
82 *> The values of the matrix column dimension N.
83 *> \endverbatim
84 *>
85 *> \param[in] NNB
86 *> \verbatim
87 *> NNB is INTEGER
88 *> The number of values of NB and NX contained in the
89 *> vectors NBVAL and NXVAL. The blocking parameters are used
90 *> in pairs (NB,NX).
91 *> \endverbatim
92 *>
93 *> \param[in] NBVAL
94 *> \verbatim
95 *> NBVAL is INTEGER array, dimension (NNB)
96 *> The values of the blocksize NB.
97 *> \endverbatim
98 *>
99 *> \param[in] NXVAL
100 *> \verbatim
101 *> NXVAL is INTEGER array, dimension (NNB)
102 *> The values of the crossover point NX.
103 *> \endverbatim
104 *>
105 *> \param[in] NNS
106 *> \verbatim
107 *> NNS is INTEGER
108 *> The number of values of NRHS contained in the vector NSVAL.
109 *> \endverbatim
110 *>
111 *> \param[in] NSVAL
112 *> \verbatim
113 *> NSVAL is INTEGER array, dimension (NNS)
114 *> The values of the number of right hand sides NRHS.
115 *> \endverbatim
116 *>
117 *> \param[in] THRESH
118 *> \verbatim
119 *> THRESH is DOUBLE PRECISION
120 *> The threshold value for the test ratios. A result is
121 *> included in the output file if RESULT >= THRESH. To have
122 *> every test ratio printed, use THRESH = 0.
123 *> \endverbatim
124 *>
125 *> \param[in] TSTERR
126 *> \verbatim
127 *> TSTERR is LOGICAL
128 *> Flag that indicates whether error exits are to be tested.
129 *> \endverbatim
130 *>
131 *> \param[out] A
132 *> \verbatim
133 *> A is COMPLEX*16 array, dimension (MMAX*NMAX)
134 *> where MMAX is the maximum value of M in MVAL and NMAX is the
135 *> maximum value of N in NVAL.
136 *> \endverbatim
137 *>
138 *> \param[out] COPYA
139 *> \verbatim
140 *> COPYA is COMPLEX*16 array, dimension (MMAX*NMAX)
141 *> \endverbatim
142 *>
143 *> \param[out] B
144 *> \verbatim
145 *> B is COMPLEX*16 array, dimension (MMAX*NSMAX)
146 *> where MMAX is the maximum value of M in MVAL and NSMAX is the
147 *> maximum value of NRHS in NSVAL.
148 *> \endverbatim
149 *>
150 *> \param[out] COPYB
151 *> \verbatim
152 *> COPYB is COMPLEX*16 array, dimension (MMAX*NSMAX)
153 *> \endverbatim
154 *>
155 *> \param[out] C
156 *> \verbatim
157 *> C is COMPLEX*16 array, dimension (MMAX*NSMAX)
158 *> \endverbatim
159 *>
160 *> \param[out] S
161 *> \verbatim
162 *> S is DOUBLE PRECISION array, dimension
163 *> (min(MMAX,NMAX))
164 *> \endverbatim
165 *>
166 *> \param[out] COPYS
167 *> \verbatim
168 *> COPYS is DOUBLE PRECISION array, dimension
169 *> (min(MMAX,NMAX))
170 *> \endverbatim
171 *>
172 *> \param[in] NOUT
173 *> \verbatim
174 *> NOUT is INTEGER
175 *> The unit number for output.
176 *> \endverbatim
177 *
178 * Authors:
179 * ========
180 *
181 *> \author Univ. of Tennessee
182 *> \author Univ. of California Berkeley
183 *> \author Univ. of Colorado Denver
184 *> \author NAG Ltd.
185 *
186 *> \date December 2016
187 *
188 *> \ingroup complex16_lin
189 *
190 * =====================================================================
191  SUBROUTINE zdrvls( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB,
192  $ nbval, nxval, thresh, tsterr, a, copya, b,
193  $ copyb, c, s, copys, nout )
194 *
195 * -- LAPACK test routine (version 3.7.0) --
196 * -- LAPACK is a software package provided by Univ. of Tennessee, --
197 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
198 * December 2016
199 *
200 * .. Scalar Arguments ..
201  LOGICAL TSTERR
202  INTEGER NM, NN, NNB, NNS, NOUT
203  DOUBLE PRECISION THRESH
204 * ..
205 * .. Array Arguments ..
206  LOGICAL DOTYPE( * )
207  INTEGER MVAL( * ), NBVAL( * ), NSVAL( * ),
208  $ nval( * ), nxval( * )
209  DOUBLE PRECISION COPYS( * ), S( * )
210  COMPLEX*16 A( * ), B( * ), C( * ), COPYA( * ), COPYB( * )
211 * ..
212 *
213 * =====================================================================
214 *
215 * .. Parameters ..
216  INTEGER NTESTS
217  parameter ( ntests = 16 )
218  INTEGER SMLSIZ
219  parameter ( smlsiz = 25 )
220  DOUBLE PRECISION ONE, ZERO
221  parameter ( one = 1.0d+0, zero = 0.0d+0 )
222  COMPLEX*16 CONE, CZERO
223  parameter ( cone = ( 1.0d+0, 0.0d+0 ),
224  $ czero = ( 0.0d+0, 0.0d+0 ) )
225 * ..
226 * .. Local Scalars ..
227  CHARACTER TRANS
228  CHARACTER*3 PATH
229  INTEGER CRANK, I, IM, IMB, IN, INB, INFO, INS, IRANK,
230  $ iscale, itran, itype, j, k, lda, ldb, ldwork,
231  $ lwlsy, lwork, m, mnmin, n, nb, ncols, nerrs,
232  $ nfail, nrhs, nrows, nrun, rank, mb,
233  $ mmax, nmax, nsmax, liwork, lrwork,
234  $ lwork_zgels, lwork_zgetsls, lwork_zgelss,
235  $ lwork_zgelsy, lwork_zgelsd,
236  $ lrwork_zgelsy, lrwork_zgelss, lrwork_zgelsd
237  DOUBLE PRECISION EPS, NORMA, NORMB, RCOND
238 * ..
239 * .. Local Arrays ..
240  INTEGER ISEED( 4 ), ISEEDY( 4 ), IWORKQUERY
241  DOUBLE PRECISION RESULT( ntests ), RWORKQUERY
242  COMPLEX*16 WORKQUERY
243 * ..
244 * .. Allocatable Arrays ..
245  COMPLEX*16, ALLOCATABLE :: WORK (:)
246  DOUBLE PRECISION, ALLOCATABLE :: RWORK (:)
247  INTEGER, ALLOCATABLE :: IWORK (:)
248 * ..
249 * .. External Functions ..
250  DOUBLE PRECISION DASUM, DLAMCH, ZQRT12, ZQRT14, ZQRT17
251  EXTERNAL dasum, dlamch, zqrt12, zqrt14, zqrt17
252 * ..
253 * .. External Subroutines ..
254  EXTERNAL alaerh, alahd, alasvm, daxpy, dlasrt, xlaenv,
255  $ zdscal, zerrls, zgels, zgelsd, zgelss,
256  $ zgelsy, zgemm, zlacpy, zlarnv, zqrt13, zqrt15,
257  $ zqrt16, zgetsls
258 * ..
259 * .. Intrinsic Functions ..
260  INTRINSIC dble, max, min, int, sqrt
261 * ..
262 * .. Scalars in Common ..
263  LOGICAL LERR, OK
264  CHARACTER*32 SRNAMT
265  INTEGER INFOT, IOUNIT
266 * ..
267 * .. Common blocks ..
268  COMMON / infoc / infot, iounit, ok, lerr
269  COMMON / srnamc / srnamt
270 * ..
271 * .. Data statements ..
272  DATA iseedy / 1988, 1989, 1990, 1991 /
273 * ..
274 * .. Executable Statements ..
275 *
276 * Initialize constants and the random number seed.
277 *
278  path( 1: 1 ) = 'Zomplex precision'
279  path( 2: 3 ) = 'LS'
280  nrun = 0
281  nfail = 0
282  nerrs = 0
283  DO 10 i = 1, 4
284  iseed( i ) = iseedy( i )
285  10 CONTINUE
286  eps = dlamch( 'Epsilon' )
287 *
288 * Threshold for rank estimation
289 *
290  rcond = sqrt( eps ) - ( sqrt( eps )-eps ) / 2
291 *
292 * Test the error exits
293 *
294  CALL xlaenv( 9, smlsiz )
295  IF( tsterr )
296  $ CALL zerrls( path, nout )
297 *
298 * Print the header if NM = 0 or NN = 0 and THRESH = 0.
299 *
300  IF( ( nm.EQ.0 .OR. nn.EQ.0 ) .AND. thresh.EQ.zero )
301  $ CALL alahd( nout, path )
302  infot = 0
303 *
304 * Compute maximal workspace needed for all routines
305 *
306  nmax = 0
307  mmax = 0
308  nsmax = 0
309  DO i = 1, nm
310  IF ( mval( i ).GT.mmax ) THEN
311  mmax = mval( i )
312  END IF
313  ENDDO
314  DO i = 1, nn
315  IF ( nval( i ).GT.nmax ) THEN
316  nmax = nval( i )
317  END IF
318  ENDDO
319  DO i = 1, nns
320  IF ( nsval( i ).GT.nsmax ) THEN
321  nsmax = nsval( i )
322  END IF
323  ENDDO
324  m = mmax
325  n = nmax
326  nrhs = nsmax
327  lda = max( 1, m )
328  ldb = max( 1, m, n )
329  mnmin = max( min( m, n ), 1 )
330 *
331 * Compute workspace needed for routines
332 * ZQRT14, ZQRT17 (two side cases), ZQRT15 and ZQRT12
333 *
334  lwork = max( ( m+n )*nrhs,
335  $ ( n+nrhs )*( m+2 ), ( m+nrhs )*( n+2 ),
336  $ max( m+mnmin, nrhs*mnmin,2*n+m ),
337  $ max( m*n+4*mnmin+max(m,n), m*n+2*mnmin+4*n ) )
338 *
339 * Compute workspace needed for ZGELS
340  CALL zgels( 'N', m, n, nrhs, a, lda, b, ldb,
341  $ workquery, -1, info )
342  lwork_zgels = int( workquery )
343 * Compute workspace needed for ZGETSLS
344  CALL zgetsls( 'N', m, n, nrhs, a, lda, b, ldb,
345  $ workquery, -1, info )
346  lwork_zgetsls = int( workquery )
347 * Compute workspace needed for ZGELSY
348  CALL zgelsy( m, n, nrhs, a, lda, b, ldb, iworkquery,
349  $ rcond, crank, workquery, -1, rwork, info )
350  lwork_zgelsy = int( workquery )
351  lrwork_zgelsy = 2*n
352 * Compute workspace needed for ZGELSS
353  CALL zgelss( m, n, nrhs, a, lda, b, ldb, s,
354  $ rcond, crank, workquery, -1 , rwork, info )
355  lwork_zgelss = int( workquery )
356  lrwork_zgelss = 5*mnmin
357 * Compute workspace needed for ZGELSD
358  CALL zgelsd( m, n, nrhs, a, lda, b, ldb, s, rcond, crank,
359  $ workquery, -1, rworkquery, iworkquery, info )
360  lwork_zgelsd = int( workquery )
361  lrwork_zgelsd = int( rworkquery )
362 * Compute LIWORK workspace needed for ZGELSY and ZGELSD
363  liwork = max( 1, n, iworkquery )
364 * Compute LRWORK workspace needed for ZGELSY, ZGELSS and ZGELSD
365  lrwork = max( 1, lrwork_zgelsy, lrwork_zgelss, lrwork_zgelsd )
366 * Compute LWORK workspace needed for all functions
367  lwork = max( 1, lwork, lwork_zgels, lwork_zgetsls, lwork_zgelsy,
368  $ lwork_zgelss, lwork_zgelsd )
369  lwlsy = lwork
370 *
371  ALLOCATE( work( lwork ) )
372  ALLOCATE( iwork( liwork ) )
373  ALLOCATE( rwork( lrwork ) )
374 *
375  DO 140 im = 1, nm
376  m = mval( im )
377  lda = max( 1, m )
378 *
379  DO 130 in = 1, nn
380  n = nval( in )
381  mnmin = max(min( m, n ),1)
382  ldb = max( 1, m, n )
383  mb = (mnmin+1)
384 *
385  DO 120 ins = 1, nns
386  nrhs = nsval( ins )
387 *
388  DO 110 irank = 1, 2
389  DO 100 iscale = 1, 3
390  itype = ( irank-1 )*3 + iscale
391  IF( .NOT.dotype( itype ) )
392  $ GO TO 100
393 *
394  IF( irank.EQ.1 ) THEN
395 *
396 * Test ZGELS
397 *
398 * Generate a matrix of scaling type ISCALE
399 *
400  CALL zqrt13( iscale, m, n, copya, lda, norma,
401  $ iseed )
402  DO 40 inb = 1, nnb
403  nb = nbval( inb )
404  CALL xlaenv( 1, nb )
405  CALL xlaenv( 3, nxval( inb ) )
406 *
407  DO 30 itran = 1, 2
408  IF( itran.EQ.1 ) THEN
409  trans = 'N'
410  nrows = m
411  ncols = n
412  ELSE
413  trans = 'C'
414  nrows = n
415  ncols = m
416  END IF
417  ldwork = max( 1, ncols )
418 *
419 * Set up a consistent rhs
420 *
421  IF( ncols.GT.0 ) THEN
422  CALL zlarnv( 2, iseed, ncols*nrhs,
423  $ work )
424  CALL zdscal( ncols*nrhs,
425  $ one / dble( ncols ), work,
426  $ 1 )
427  END IF
428  CALL zgemm( trans, 'No transpose', nrows,
429  $ nrhs, ncols, cone, copya, lda,
430  $ work, ldwork, czero, b, ldb )
431  CALL zlacpy( 'Full', nrows, nrhs, b, ldb,
432  $ copyb, ldb )
433 *
434 * Solve LS or overdetermined system
435 *
436  IF( m.GT.0 .AND. n.GT.0 ) THEN
437  CALL zlacpy( 'Full', m, n, copya, lda,
438  $ a, lda )
439  CALL zlacpy( 'Full', nrows, nrhs,
440  $ copyb, ldb, b, ldb )
441  END IF
442  srnamt = 'ZGELS '
443  CALL zgels( trans, m, n, nrhs, a, lda, b,
444  $ ldb, work, lwork, info )
445 *
446  IF( info.NE.0 )
447  $ CALL alaerh( path, 'ZGELS ', info, 0,
448  $ trans, m, n, nrhs, -1, nb,
449  $ itype, nfail, nerrs,
450  $ nout )
451 *
452 * Check correctness of results
453 *
454  ldwork = max( 1, nrows )
455  IF( nrows.GT.0 .AND. nrhs.GT.0 )
456  $ CALL zlacpy( 'Full', nrows, nrhs,
457  $ copyb, ldb, c, ldb )
458  CALL zqrt16( trans, m, n, nrhs, copya,
459  $ lda, b, ldb, c, ldb, rwork,
460  $ result( 1 ) )
461 *
462  IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
463  $ ( itran.EQ.2 .AND. m.LT.n ) ) THEN
464 *
465 * Solving LS system
466 *
467  result( 2 ) = zqrt17( trans, 1, m, n,
468  $ nrhs, copya, lda, b, ldb,
469  $ copyb, ldb, c, work,
470  $ lwork )
471  ELSE
472 *
473 * Solving overdetermined system
474 *
475  result( 2 ) = zqrt14( trans, m, n,
476  $ nrhs, copya, lda, b, ldb,
477  $ work, lwork )
478  END IF
479 *
480 * Print information about the tests that
481 * did not pass the threshold.
482 *
483  DO 20 k = 1, 2
484  IF( result( k ).GE.thresh ) THEN
485  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
486  $ CALL alahd( nout, path )
487  WRITE( nout, fmt = 9999 )trans, m,
488  $ n, nrhs, nb, itype, k,
489  $ result( k )
490  nfail = nfail + 1
491  END IF
492  20 CONTINUE
493  nrun = nrun + 2
494  30 CONTINUE
495  40 CONTINUE
496 *
497 *
498 * Test ZGETSLS
499 *
500 * Generate a matrix of scaling type ISCALE
501 *
502  CALL zqrt13( iscale, m, n, copya, lda, norma,
503  $ iseed )
504  DO 65 inb = 1, nnb
505  mb = nbval( inb )
506  CALL xlaenv( 1, mb )
507  DO 62 imb = 1, nnb
508  nb = nbval( imb )
509  CALL xlaenv( 2, nb )
510 *
511  DO 60 itran = 1, 2
512  IF( itran.EQ.1 ) THEN
513  trans = 'N'
514  nrows = m
515  ncols = n
516  ELSE
517  trans = 'C'
518  nrows = n
519  ncols = m
520  END IF
521  ldwork = max( 1, ncols )
522 *
523 * Set up a consistent rhs
524 *
525  IF( ncols.GT.0 ) THEN
526  CALL zlarnv( 2, iseed, ncols*nrhs,
527  $ work )
528  CALL zscal( ncols*nrhs,
529  $ one / dble( ncols ), work,
530  $ 1 )
531  END IF
532  CALL zgemm( trans, 'No transpose', nrows,
533  $ nrhs, ncols, cone, copya, lda,
534  $ work, ldwork, czero, b, ldb )
535  CALL zlacpy( 'Full', nrows, nrhs, b, ldb,
536  $ copyb, ldb )
537 *
538 * Solve LS or overdetermined system
539 *
540  IF( m.GT.0 .AND. n.GT.0 ) THEN
541  CALL zlacpy( 'Full', m, n, copya, lda,
542  $ a, lda )
543  CALL zlacpy( 'Full', nrows, nrhs,
544  $ copyb, ldb, b, ldb )
545  END IF
546  srnamt = 'ZGETSLS '
547  CALL zgetsls( trans, m, n, nrhs, a,
548  $ lda, b, ldb, work, lwork, info )
549  IF( info.NE.0 )
550  $ CALL alaerh( path, 'ZGETSLS ', info, 0,
551  $ trans, m, n, nrhs, -1, nb,
552  $ itype, nfail, nerrs,
553  $ nout )
554 *
555 * Check correctness of results
556 *
557  ldwork = max( 1, nrows )
558  IF( nrows.GT.0 .AND. nrhs.GT.0 )
559  $ CALL zlacpy( 'Full', nrows, nrhs,
560  $ copyb, ldb, c, ldb )
561  CALL zqrt16( trans, m, n, nrhs, copya,
562  $ lda, b, ldb, c, ldb, work,
563  $ result( 15 ) )
564 *
565  IF( ( itran.EQ.1 .AND. m.GE.n ) .OR.
566  $ ( itran.EQ.2 .AND. m.LT.n ) ) THEN
567 *
568 * Solving LS system
569 *
570  result( 16 ) = zqrt17( trans, 1, m, n,
571  $ nrhs, copya, lda, b, ldb,
572  $ copyb, ldb, c, work,
573  $ lwork )
574  ELSE
575 *
576 * Solving overdetermined system
577 *
578  result( 16 ) = zqrt14( trans, m, n,
579  $ nrhs, copya, lda, b, ldb,
580  $ work, lwork )
581  END IF
582 *
583 * Print information about the tests that
584 * did not pass the threshold.
585 *
586  DO 50 k = 15, 16
587  IF( result( k ).GE.thresh ) THEN
588  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
589  $ CALL alahd( nout, path )
590  WRITE( nout, fmt = 9997 )trans, m,
591  $ n, nrhs, mb, nb, itype, k,
592  $ result( k )
593  nfail = nfail + 1
594  END IF
595  50 CONTINUE
596  nrun = nrun + 2
597  60 CONTINUE
598  62 CONTINUE
599  65 CONTINUE
600  END IF
601 *
602 * Generate a matrix of scaling type ISCALE and rank
603 * type IRANK.
604 *
605  CALL zqrt15( iscale, irank, m, n, nrhs, copya, lda,
606  $ copyb, ldb, copys, rank, norma, normb,
607  $ iseed, work, lwork )
608 *
609 * workspace used: MAX(M+MIN(M,N),NRHS*MIN(M,N),2*N+M)
610 *
611  ldwork = max( 1, m )
612 *
613 * Loop for testing different block sizes.
614 *
615  DO 90 inb = 1, nnb
616  nb = nbval( inb )
617  CALL xlaenv( 1, nb )
618  CALL xlaenv( 3, nxval( inb ) )
619 *
620 * Test ZGELSY
621 *
622 * ZGELSY: Compute the minimum-norm solution
623 * X to min( norm( A * X - B ) )
624 * using the rank-revealing orthogonal
625 * factorization.
626 *
627  CALL zlacpy( 'Full', m, n, copya, lda, a, lda )
628  CALL zlacpy( 'Full', m, nrhs, copyb, ldb, b,
629  $ ldb )
630 *
631 * Initialize vector IWORK.
632 *
633  DO 70 j = 1, n
634  iwork( j ) = 0
635  70 CONTINUE
636 *
637  srnamt = 'ZGELSY'
638  CALL zgelsy( m, n, nrhs, a, lda, b, ldb, iwork,
639  $ rcond, crank, work, lwlsy, rwork,
640  $ info )
641  IF( info.NE.0 )
642  $ CALL alaerh( path, 'ZGELSY', info, 0, ' ', m,
643  $ n, nrhs, -1, nb, itype, nfail,
644  $ nerrs, nout )
645 *
646 * workspace used: 2*MNMIN+NB*NB+NB*MAX(N,NRHS)
647 *
648 * Test 3: Compute relative error in svd
649 * workspace: M*N + 4*MIN(M,N) + MAX(M,N)
650 *
651  result( 3 ) = zqrt12( crank, crank, a, lda,
652  $ copys, work, lwork, rwork )
653 *
654 * Test 4: Compute error in solution
655 * workspace: M*NRHS + M
656 *
657  CALL zlacpy( 'Full', m, nrhs, copyb, ldb, work,
658  $ ldwork )
659  CALL zqrt16( 'No transpose', m, n, nrhs, copya,
660  $ lda, b, ldb, work, ldwork, rwork,
661  $ result( 4 ) )
662 *
663 * Test 5: Check norm of r'*A
664 * workspace: NRHS*(M+N)
665 *
666  result( 5 ) = zero
667  IF( m.GT.crank )
668  $ result( 5 ) = zqrt17( 'No transpose', 1, m,
669  $ n, nrhs, copya, lda, b, ldb,
670  $ copyb, ldb, c, work, lwork )
671 *
672 * Test 6: Check if x is in the rowspace of A
673 * workspace: (M+NRHS)*(N+2)
674 *
675  result( 6 ) = zero
676 *
677  IF( n.GT.crank )
678  $ result( 6 ) = zqrt14( 'No transpose', m, n,
679  $ nrhs, copya, lda, b, ldb,
680  $ work, lwork )
681 *
682 * Test ZGELSS
683 *
684 * ZGELSS: Compute the minimum-norm solution
685 * X to min( norm( A * X - B ) )
686 * using the SVD.
687 *
688  CALL zlacpy( 'Full', m, n, copya, lda, a, lda )
689  CALL zlacpy( 'Full', m, nrhs, copyb, ldb, b,
690  $ ldb )
691  srnamt = 'ZGELSS'
692  CALL zgelss( m, n, nrhs, a, lda, b, ldb, s,
693  $ rcond, crank, work, lwork, rwork,
694  $ info )
695 *
696  IF( info.NE.0 )
697  $ CALL alaerh( path, 'ZGELSS', info, 0, ' ', m,
698  $ n, nrhs, -1, nb, itype, nfail,
699  $ nerrs, nout )
700 *
701 * workspace used: 3*min(m,n) +
702 * max(2*min(m,n),nrhs,max(m,n))
703 *
704 * Test 7: Compute relative error in svd
705 *
706  IF( rank.GT.0 ) THEN
707  CALL daxpy( mnmin, -one, copys, 1, s, 1 )
708  result( 7 ) = dasum( mnmin, s, 1 ) /
709  $ dasum( mnmin, copys, 1 ) /
710  $ ( eps*dble( mnmin ) )
711  ELSE
712  result( 7 ) = zero
713  END IF
714 *
715 * Test 8: Compute error in solution
716 *
717  CALL zlacpy( 'Full', m, nrhs, copyb, ldb, work,
718  $ ldwork )
719  CALL zqrt16( 'No transpose', m, n, nrhs, copya,
720  $ lda, b, ldb, work, ldwork, rwork,
721  $ result( 8 ) )
722 *
723 * Test 9: Check norm of r'*A
724 *
725  result( 9 ) = zero
726  IF( m.GT.crank )
727  $ result( 9 ) = zqrt17( 'No transpose', 1, m,
728  $ n, nrhs, copya, lda, b, ldb,
729  $ copyb, ldb, c, work, lwork )
730 *
731 * Test 10: Check if x is in the rowspace of A
732 *
733  result( 10 ) = zero
734  IF( n.GT.crank )
735  $ result( 10 ) = zqrt14( 'No transpose', m, n,
736  $ nrhs, copya, lda, b, ldb,
737  $ work, lwork )
738 *
739 * Test ZGELSD
740 *
741 * ZGELSD: Compute the minimum-norm solution X
742 * to min( norm( A * X - B ) ) using a
743 * divide and conquer SVD.
744 *
745  CALL xlaenv( 9, 25 )
746 *
747  CALL zlacpy( 'Full', m, n, copya, lda, a, lda )
748  CALL zlacpy( 'Full', m, nrhs, copyb, ldb, b,
749  $ ldb )
750 *
751  srnamt = 'ZGELSD'
752  CALL zgelsd( m, n, nrhs, a, lda, b, ldb, s,
753  $ rcond, crank, work, lwork, rwork,
754  $ iwork, info )
755  IF( info.NE.0 )
756  $ CALL alaerh( path, 'ZGELSD', info, 0, ' ', m,
757  $ n, nrhs, -1, nb, itype, nfail,
758  $ nerrs, nout )
759 *
760 * Test 11: Compute relative error in svd
761 *
762  IF( rank.GT.0 ) THEN
763  CALL daxpy( mnmin, -one, copys, 1, s, 1 )
764  result( 11 ) = dasum( mnmin, s, 1 ) /
765  $ dasum( mnmin, copys, 1 ) /
766  $ ( eps*dble( mnmin ) )
767  ELSE
768  result( 11 ) = zero
769  END IF
770 *
771 * Test 12: Compute error in solution
772 *
773  CALL zlacpy( 'Full', m, nrhs, copyb, ldb, work,
774  $ ldwork )
775  CALL zqrt16( 'No transpose', m, n, nrhs, copya,
776  $ lda, b, ldb, work, ldwork, rwork,
777  $ result( 12 ) )
778 *
779 * Test 13: Check norm of r'*A
780 *
781  result( 13 ) = zero
782  IF( m.GT.crank )
783  $ result( 13 ) = zqrt17( 'No transpose', 1, m,
784  $ n, nrhs, copya, lda, b, ldb,
785  $ copyb, ldb, c, work, lwork )
786 *
787 * Test 14: Check if x is in the rowspace of A
788 *
789  result( 14 ) = zero
790  IF( n.GT.crank )
791  $ result( 14 ) = zqrt14( 'No transpose', m, n,
792  $ nrhs, copya, lda, b, ldb,
793  $ work, lwork )
794 *
795 * Print information about the tests that did not
796 * pass the threshold.
797 *
798  DO 80 k = 3, 14
799  IF( result( k ).GE.thresh ) THEN
800  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
801  $ CALL alahd( nout, path )
802  WRITE( nout, fmt = 9998 )m, n, nrhs, nb,
803  $ itype, k, result( k )
804  nfail = nfail + 1
805  END IF
806  80 CONTINUE
807  nrun = nrun + 12
808 *
809  90 CONTINUE
810  100 CONTINUE
811  110 CONTINUE
812  120 CONTINUE
813  130 CONTINUE
814  140 CONTINUE
815 *
816 * Print a summary of the results.
817 *
818  CALL alasvm( path, nout, nfail, nrun, nerrs )
819 *
820  9999 FORMAT( ' TRANS=''', a1, ''', M=', i5, ', N=', i5, ', NRHS=', i4,
821  $ ', NB=', i4, ', type', i2, ', test(', i2, ')=', g12.5 )
822  9998 FORMAT( ' M=', i5, ', N=', i5, ', NRHS=', i4, ', NB=', i4,
823  $ ', type', i2, ', test(', i2, ')=', g12.5 )
824  9997 FORMAT( ' TRANS=''', a1,' M=', i5, ', N=', i5, ', NRHS=', i4,
825  $ ', MB=', i4,', NB=', i4,', type', i2,
826  $ ', test(', i2, ')=', g12.5 )
827 *
828  DEALLOCATE( work )
829  DEALLOCATE( iwork )
830  DEALLOCATE( rwork )
831  RETURN
832 *
833 * End of ZDRVLS
834 *
835  END
alasvm
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:75
alahd
subroutine alahd(IOUNIT, PATH)
ALAHD
Definition: alahd.f:107
zqrt16
subroutine zqrt16(TRANS, M, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
ZQRT16
Definition: zqrt16.f:135
zgels
subroutine zgels(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO)
ZGELS solves overdetermined or underdetermined systems for GE matrices
Definition: zgels.f:184
alaerh
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:149
zgetsls
subroutine zgetsls(TRANS, M, N, NRHS, A, LDA, B, LDB, WORK, LWORK, INFO)
Definition: zgetsls.f:162
zgelsy
subroutine zgelsy(M, N, NRHS, A, LDA, B, LDB, JPVT, RCOND, RANK, WORK, LWORK, RWORK, INFO)
ZGELSY solves overdetermined or underdetermined systems for GE matrices
Definition: zgelsy.f:212
zgelss
subroutine zgelss(M, N, NRHS, A, LDA, B, LDB, S, RCOND, RANK, WORK, LWORK, RWORK, INFO)
ZGELSS solves overdetermined or underdetermined systems for GE matrices
Definition: zgelss.f:180
zgelsd
subroutine zgelsd(M, N, NRHS, A, LDA, B, LDB, S, RCOND, RANK, WORK, LWORK, RWORK, IWORK, INFO)
ZGELSD computes the minimum-norm solution to a linear least squares problem for GE matrices ...
Definition: zgelsd.f:227
zlacpy
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:105
zgemm
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:189
daxpy
subroutine daxpy(N, DA, DX, INCX, DY, INCY)
DAXPY
Definition: daxpy.f:54
zdrvls
subroutine zdrvls(DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB, NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B, COPYB, C, S, COPYS, NOUT)
ZDRVLS
Definition: zdrvls.f:194
xlaenv
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:83
zqrt13
subroutine zqrt13(SCALE, M, N, A, LDA, NORMA, ISEED)
ZQRT13
Definition: zqrt13.f:93
zerrls
subroutine zerrls(PATH, NUNIT)
ZERRLS
Definition: zerrls.f:57
dlasrt
subroutine dlasrt(ID, N, D, INFO)
DLASRT sorts numbers in increasing or decreasing order.
Definition: dlasrt.f:90
zdscal
subroutine zdscal(N, DA, ZX, INCX)
ZDSCAL
Definition: zdscal.f:54
zlarnv
subroutine zlarnv(IDIST, ISEED, N, X)
ZLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition: zlarnv.f:101
zscal
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL
Definition: zscal.f:54
zqrt15
subroutine zqrt15(SCALE, RKSEL, M, N, NRHS, A, LDA, B, LDB, S, RANK, NORMA, NORMB, ISEED, WORK, LWORK)
ZQRT15
Definition: zqrt15.f:151