LAPACK  3.9.0
LAPACK: Linear Algebra PACKage
cdrvbd.f
Go to the documentation of this file.
1 *> \brief \b CDRVBD
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 CDRVBD( NSIZES, MM, NN, NTYPES, DOTYPE, ISEED, THRESH,
12 * A, LDA, U, LDU, VT, LDVT, ASAV, USAV, VTSAV, S,
13 * SSAV, E, WORK, LWORK, RWORK, IWORK, NOUNIT,
14 * INFO )
15 *
16 * .. Scalar Arguments ..
17 * INTEGER INFO, LDA, LDU, LDVT, LWORK, NOUNIT, NSIZES,
18 * $ NTYPES
19 * REAL THRESH
20 * ..
21 * .. Array Arguments ..
22 * LOGICAL DOTYPE( * )
23 * INTEGER ISEED( 4 ), IWORK( * ), MM( * ), NN( * )
24 * REAL E( * ), RWORK( * ), S( * ), SSAV( * )
25 * COMPLEX A( LDA, * ), ASAV( LDA, * ), U( LDU, * ),
26 * $ USAV( LDU, * ), VT( LDVT, * ),
27 * $ VTSAV( LDVT, * ), WORK( * )
28 * ..
29 *
30 *
31 *> \par Purpose:
32 * =============
33 *>
34 *> \verbatim
35 *>
36 *> CDRVBD checks the singular value decomposition (SVD) driver CGESVD,
37 *> CGESDD, CGESVJ, CGEJSV, CGESVDX, and CGESVDQ.
38 *>
39 *> CGESVD and CGESDD factors A = U diag(S) VT, where U and VT are
40 *> unitary and diag(S) is diagonal with the entries of the array S on
41 *> its diagonal. The entries of S are the singular values, nonnegative
42 *> and stored in decreasing order. U and VT can be optionally not
43 *> computed, overwritten on A, or computed partially.
44 *>
45 *> A is M by N. Let MNMIN = min( M, N ). S has dimension MNMIN.
46 *> U can be M by M or M by MNMIN. VT can be N by N or MNMIN by N.
47 *>
48 *> When CDRVBD is called, a number of matrix "sizes" (M's and N's)
49 *> and a number of matrix "types" are specified. For each size (M,N)
50 *> and each type of matrix, and for the minimal workspace as well as
51 *> workspace adequate to permit blocking, an M x N matrix "A" will be
52 *> generated and used to test the SVD routines. For each matrix, A will
53 *> be factored as A = U diag(S) VT and the following 12 tests computed:
54 *>
55 *> Test for CGESVD:
56 *>
57 *> (1) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
58 *>
59 *> (2) | I - U'U | / ( M ulp )
60 *>
61 *> (3) | I - VT VT' | / ( N ulp )
62 *>
63 *> (4) S contains MNMIN nonnegative values in decreasing order.
64 *> (Return 0 if true, 1/ULP if false.)
65 *>
66 *> (5) | U - Upartial | / ( M ulp ) where Upartial is a partially
67 *> computed U.
68 *>
69 *> (6) | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
70 *> computed VT.
71 *>
72 *> (7) | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
73 *> vector of singular values from the partial SVD
74 *>
75 *> Test for CGESDD:
76 *>
77 *> (8) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
78 *>
79 *> (9) | I - U'U | / ( M ulp )
80 *>
81 *> (10) | I - VT VT' | / ( N ulp )
82 *>
83 *> (11) S contains MNMIN nonnegative values in decreasing order.
84 *> (Return 0 if true, 1/ULP if false.)
85 *>
86 *> (12) | U - Upartial | / ( M ulp ) where Upartial is a partially
87 *> computed U.
88 *>
89 *> (13) | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
90 *> computed VT.
91 *>
92 *> (14) | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
93 *> vector of singular values from the partial SVD
94 *>
95 *> Test for CGESVDQ:
96 *>
97 *> (36) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
98 *>
99 *> (37) | I - U'U | / ( M ulp )
100 *>
101 *> (38) | I - VT VT' | / ( N ulp )
102 *>
103 *> (39) S contains MNMIN nonnegative values in decreasing order.
104 *> (Return 0 if true, 1/ULP if false.)
105 *>
106 *> Test for CGESVJ:
107 *>
108 *> (15) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
109 *>
110 *> (16) | I - U'U | / ( M ulp )
111 *>
112 *> (17) | I - VT VT' | / ( N ulp )
113 *>
114 *> (18) S contains MNMIN nonnegative values in decreasing order.
115 *> (Return 0 if true, 1/ULP if false.)
116 *>
117 *> Test for CGEJSV:
118 *>
119 *> (19) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
120 *>
121 *> (20) | I - U'U | / ( M ulp )
122 *>
123 *> (21) | I - VT VT' | / ( N ulp )
124 *>
125 *> (22) S contains MNMIN nonnegative values in decreasing order.
126 *> (Return 0 if true, 1/ULP if false.)
127 *>
128 *> Test for CGESVDX( 'V', 'V', 'A' )/CGESVDX( 'N', 'N', 'A' )
129 *>
130 *> (23) | A - U diag(S) VT | / ( |A| max(M,N) ulp )
131 *>
132 *> (24) | I - U'U | / ( M ulp )
133 *>
134 *> (25) | I - VT VT' | / ( N ulp )
135 *>
136 *> (26) S contains MNMIN nonnegative values in decreasing order.
137 *> (Return 0 if true, 1/ULP if false.)
138 *>
139 *> (27) | U - Upartial | / ( M ulp ) where Upartial is a partially
140 *> computed U.
141 *>
142 *> (28) | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
143 *> computed VT.
144 *>
145 *> (29) | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
146 *> vector of singular values from the partial SVD
147 *>
148 *> Test for CGESVDX( 'V', 'V', 'I' )
149 *>
150 *> (30) | U' A VT''' - diag(S) | / ( |A| max(M,N) ulp )
151 *>
152 *> (31) | I - U'U | / ( M ulp )
153 *>
154 *> (32) | I - VT VT' | / ( N ulp )
155 *>
156 *> Test for CGESVDX( 'V', 'V', 'V' )
157 *>
158 *> (33) | U' A VT''' - diag(S) | / ( |A| max(M,N) ulp )
159 *>
160 *> (34) | I - U'U | / ( M ulp )
161 *>
162 *> (35) | I - VT VT' | / ( N ulp )
163 *>
164 *> The "sizes" are specified by the arrays MM(1:NSIZES) and
165 *> NN(1:NSIZES); the value of each element pair (MM(j),NN(j))
166 *> specifies one size. The "types" are specified by a logical array
167 *> DOTYPE( 1:NTYPES ); if DOTYPE(j) is .TRUE., then matrix type "j"
168 *> will be generated.
169 *> Currently, the list of possible types is:
170 *>
171 *> (1) The zero matrix.
172 *> (2) The identity matrix.
173 *> (3) A matrix of the form U D V, where U and V are unitary and
174 *> D has evenly spaced entries 1, ..., ULP with random signs
175 *> on the diagonal.
176 *> (4) Same as (3), but multiplied by the underflow-threshold / ULP.
177 *> (5) Same as (3), but multiplied by the overflow-threshold * ULP.
178 *> \endverbatim
179 *
180 * Arguments:
181 * ==========
182 *
183 *> \param[in] NSIZES
184 *> \verbatim
185 *> NSIZES is INTEGER
186 *> The number of sizes of matrices to use. If it is zero,
187 *> CDRVBD does nothing. It must be at least zero.
188 *> \endverbatim
189 *>
190 *> \param[in] MM
191 *> \verbatim
192 *> MM is INTEGER array, dimension (NSIZES)
193 *> An array containing the matrix "heights" to be used. For
194 *> each j=1,...,NSIZES, if MM(j) is zero, then MM(j) and NN(j)
195 *> will be ignored. The MM(j) values must be at least zero.
196 *> \endverbatim
197 *>
198 *> \param[in] NN
199 *> \verbatim
200 *> NN is INTEGER array, dimension (NSIZES)
201 *> An array containing the matrix "widths" to be used. For
202 *> each j=1,...,NSIZES, if NN(j) is zero, then MM(j) and NN(j)
203 *> will be ignored. The NN(j) values must be at least zero.
204 *> \endverbatim
205 *>
206 *> \param[in] NTYPES
207 *> \verbatim
208 *> NTYPES is INTEGER
209 *> The number of elements in DOTYPE. If it is zero, CDRVBD
210 *> does nothing. It must be at least zero. If it is MAXTYP+1
211 *> and NSIZES is 1, then an additional type, MAXTYP+1 is
212 *> defined, which is to use whatever matrices are in A and B.
213 *> This is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
214 *> DOTYPE(MAXTYP+1) is .TRUE. .
215 *> \endverbatim
216 *>
217 *> \param[in] DOTYPE
218 *> \verbatim
219 *> DOTYPE is LOGICAL array, dimension (NTYPES)
220 *> If DOTYPE(j) is .TRUE., then for each size (m,n), a matrix
221 *> of type j will be generated. If NTYPES is smaller than the
222 *> maximum number of types defined (PARAMETER MAXTYP), then
223 *> types NTYPES+1 through MAXTYP will not be generated. If
224 *> NTYPES is larger than MAXTYP, DOTYPE(MAXTYP+1) through
225 *> DOTYPE(NTYPES) will be ignored.
226 *> \endverbatim
227 *>
228 *> \param[in,out] ISEED
229 *> \verbatim
230 *> ISEED is INTEGER array, dimension (4)
231 *> On entry ISEED specifies the seed of the random number
232 *> generator. The array elements should be between 0 and 4095;
233 *> if not they will be reduced mod 4096. Also, ISEED(4) must
234 *> be odd. The random number generator uses a linear
235 *> congruential sequence limited to small integers, and so
236 *> should produce machine independent random numbers. The
237 *> values of ISEED are changed on exit, and can be used in the
238 *> next call to CDRVBD to continue the same random number
239 *> sequence.
240 *> \endverbatim
241 *>
242 *> \param[in] THRESH
243 *> \verbatim
244 *> THRESH is REAL
245 *> A test will count as "failed" if the "error", computed as
246 *> described above, exceeds THRESH. Note that the error
247 *> is scaled to be O(1), so THRESH should be a reasonably
248 *> small multiple of 1, e.g., 10 or 100. In particular,
249 *> it should not depend on the precision (single vs. double)
250 *> or the size of the matrix. It must be at least zero.
251 *> \endverbatim
252 *>
253 *> \param[out] A
254 *> \verbatim
255 *> A is COMPLEX array, dimension (LDA,max(NN))
256 *> Used to hold the matrix whose singular values are to be
257 *> computed. On exit, A contains the last matrix actually
258 *> used.
259 *> \endverbatim
260 *>
261 *> \param[in] LDA
262 *> \verbatim
263 *> LDA is INTEGER
264 *> The leading dimension of A. It must be at
265 *> least 1 and at least max( MM ).
266 *> \endverbatim
267 *>
268 *> \param[out] U
269 *> \verbatim
270 *> U is COMPLEX array, dimension (LDU,max(MM))
271 *> Used to hold the computed matrix of right singular vectors.
272 *> On exit, U contains the last such vectors actually computed.
273 *> \endverbatim
274 *>
275 *> \param[in] LDU
276 *> \verbatim
277 *> LDU is INTEGER
278 *> The leading dimension of U. It must be at
279 *> least 1 and at least max( MM ).
280 *> \endverbatim
281 *>
282 *> \param[out] VT
283 *> \verbatim
284 *> VT is COMPLEX array, dimension (LDVT,max(NN))
285 *> Used to hold the computed matrix of left singular vectors.
286 *> On exit, VT contains the last such vectors actually computed.
287 *> \endverbatim
288 *>
289 *> \param[in] LDVT
290 *> \verbatim
291 *> LDVT is INTEGER
292 *> The leading dimension of VT. It must be at
293 *> least 1 and at least max( NN ).
294 *> \endverbatim
295 *>
296 *> \param[out] ASAV
297 *> \verbatim
298 *> ASAV is COMPLEX array, dimension (LDA,max(NN))
299 *> Used to hold a different copy of the matrix whose singular
300 *> values are to be computed. On exit, A contains the last
301 *> matrix actually used.
302 *> \endverbatim
303 *>
304 *> \param[out] USAV
305 *> \verbatim
306 *> USAV is COMPLEX array, dimension (LDU,max(MM))
307 *> Used to hold a different copy of the computed matrix of
308 *> right singular vectors. On exit, USAV contains the last such
309 *> vectors actually computed.
310 *> \endverbatim
311 *>
312 *> \param[out] VTSAV
313 *> \verbatim
314 *> VTSAV is COMPLEX array, dimension (LDVT,max(NN))
315 *> Used to hold a different copy of the computed matrix of
316 *> left singular vectors. On exit, VTSAV contains the last such
317 *> vectors actually computed.
318 *> \endverbatim
319 *>
320 *> \param[out] S
321 *> \verbatim
322 *> S is REAL array, dimension (max(min(MM,NN)))
323 *> Contains the computed singular values.
324 *> \endverbatim
325 *>
326 *> \param[out] SSAV
327 *> \verbatim
328 *> SSAV is REAL array, dimension (max(min(MM,NN)))
329 *> Contains another copy of the computed singular values.
330 *> \endverbatim
331 *>
332 *> \param[out] E
333 *> \verbatim
334 *> E is REAL array, dimension (max(min(MM,NN)))
335 *> Workspace for CGESVD.
336 *> \endverbatim
337 *>
338 *> \param[out] WORK
339 *> \verbatim
340 *> WORK is COMPLEX array, dimension (LWORK)
341 *> \endverbatim
342 *>
343 *> \param[in] LWORK
344 *> \verbatim
345 *> LWORK is INTEGER
346 *> The number of entries in WORK. This must be at least
347 *> MAX(3*MIN(M,N)+MAX(M,N)**2,5*MIN(M,N),3*MAX(M,N)) for all
348 *> pairs (M,N)=(MM(j),NN(j))
349 *> \endverbatim
350 *>
351 *> \param[out] RWORK
352 *> \verbatim
353 *> RWORK is REAL array,
354 *> dimension ( 5*max(max(MM,NN)) )
355 *> \endverbatim
356 *>
357 *> \param[out] IWORK
358 *> \verbatim
359 *> IWORK is INTEGER array, dimension at least 8*min(M,N)
360 *> \endverbatim
361 *>
362 *> \param[in] NOUNIT
363 *> \verbatim
364 *> NOUNIT is INTEGER
365 *> The FORTRAN unit number for printing out error messages
366 *> (e.g., if a routine returns IINFO not equal to 0.)
367 *> \endverbatim
368 *>
369 *> \param[out] INFO
370 *> \verbatim
371 *> INFO is INTEGER
372 *> If 0, then everything ran OK.
373 *> -1: NSIZES < 0
374 *> -2: Some MM(j) < 0
375 *> -3: Some NN(j) < 0
376 *> -4: NTYPES < 0
377 *> -7: THRESH < 0
378 *> -10: LDA < 1 or LDA < MMAX, where MMAX is max( MM(j) ).
379 *> -12: LDU < 1 or LDU < MMAX.
380 *> -14: LDVT < 1 or LDVT < NMAX, where NMAX is max( NN(j) ).
381 *> -29: LWORK too small.
382 *> If CLATMS, or CGESVD returns an error code, the
383 *> absolute value of it is returned.
384 *> \endverbatim
385 *
386 * Authors:
387 * ========
388 *
389 *> \author Univ. of Tennessee
390 *> \author Univ. of California Berkeley
391 *> \author Univ. of Colorado Denver
392 *> \author NAG Ltd.
393 *
394 *> \date June 2016
395 *
396 *> \ingroup complex_eig
397 *
398 * =====================================================================
399  SUBROUTINE cdrvbd( NSIZES, MM, NN, NTYPES, DOTYPE, ISEED, THRESH,
400  $ A, LDA, U, LDU, VT, LDVT, ASAV, USAV, VTSAV, S,
401  $ SSAV, E, WORK, LWORK, RWORK, IWORK, NOUNIT,
402  $ INFO )
403 *
404 * -- LAPACK test routine (version 3.7.0) --
405 * -- LAPACK is a software package provided by Univ. of Tennessee, --
406 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
407 * June 2016
408 *
409  IMPLICIT NONE
410 *
411 * .. Scalar Arguments ..
412  INTEGER INFO, LDA, LDU, LDVT, LWORK, NOUNIT, NSIZES,
413  $ NTYPES
414  REAL THRESH
415 * ..
416 * .. Array Arguments ..
417  LOGICAL DOTYPE( * )
418  INTEGER ISEED( 4 ), IWORK( * ), MM( * ), NN( * )
419  REAL E( * ), RWORK( * ), S( * ), SSAV( * )
420  COMPLEX A( LDA, * ), ASAV( LDA, * ), U( LDU, * ),
421  $ usav( ldu, * ), vt( ldvt, * ),
422  $ vtsav( ldvt, * ), work( * )
423 * ..
424 *
425 * =====================================================================
426 *
427 * .. Parameters ..
428  REAL ZERO, ONE, TWO, HALF
429  PARAMETER ( ZERO = 0.0e0, one = 1.0e0, two = 2.0e0,
430  $ half = 0.5e0 )
431  COMPLEX CZERO, CONE
432  parameter( czero = ( 0.0e+0, 0.0e+0 ),
433  $ cone = ( 1.0e+0, 0.0e+0 ) )
434  INTEGER MAXTYP
435  parameter( maxtyp = 5 )
436 * ..
437 * .. Local Scalars ..
438  LOGICAL BADMM, BADNN
439  CHARACTER JOBQ, JOBU, JOBVT, RANGE
440  INTEGER I, IINFO, IJQ, IJU, IJVT, IL, IU, ITEMP,
441  $ iwspc, iwtmp, j, jsize, jtype, lswork, m,
442  $ minwrk, mmax, mnmax, mnmin, mtypes, n,
443  $ nerrs, nfail, nmax, ns, nsi, nsv, ntest,
444  $ ntestf, ntestt, lrwork
445  REAL ANORM, DIF, DIV, OVFL, RTUNFL, ULP, ULPINV,
446  $ UNFL, VL, VU
447 * ..
448 * .. Local Scalars for CGESVDQ ..
449  INTEGER LIWORK, NUMRANK
450 * ..
451 * .. Local Arrays ..
452  CHARACTER CJOB( 4 ), CJOBR( 3 ), CJOBV( 2 )
453  INTEGER IOLDSD( 4 ), ISEED2( 4 )
454  REAL RESULT( 39 )
455 * ..
456 * .. External Functions ..
457  REAL SLAMCH, SLARND
458  EXTERNAL SLAMCH, SLARND
459 * ..
460 * .. External Subroutines ..
461  EXTERNAL alasvm, xerbla, cbdt01, cbdt05, cgesdd,
464 * ..
465 * .. Intrinsic Functions ..
466  INTRINSIC abs, real, max, min
467 * ..
468 * .. Scalars in Common ..
469  CHARACTER*32 SRNAMT
470 * ..
471 * .. Common blocks ..
472  COMMON / srnamc / srnamt
473 * ..
474 * .. Data statements ..
475  DATA cjob / 'N', 'O', 'S', 'A' /
476  DATA cjobr / 'A', 'V', 'I' /
477  DATA cjobv / 'N', 'V' /
478 * ..
479 * .. Executable Statements ..
480 *
481 * Check for errors
482 *
483  info = 0
484 *
485 * Important constants
486 *
487  nerrs = 0
488  ntestt = 0
489  ntestf = 0
490  badmm = .false.
491  badnn = .false.
492  mmax = 1
493  nmax = 1
494  mnmax = 1
495  minwrk = 1
496  DO 10 j = 1, nsizes
497  mmax = max( mmax, mm( j ) )
498  IF( mm( j ).LT.0 )
499  $ badmm = .true.
500  nmax = max( nmax, nn( j ) )
501  IF( nn( j ).LT.0 )
502  $ badnn = .true.
503  mnmax = max( mnmax, min( mm( j ), nn( j ) ) )
504  minwrk = max( minwrk, max( 3*min( mm( j ),
505  $ nn( j ) )+max( mm( j ), nn( j ) )**2, 5*min( mm( j ),
506  $ nn( j ) ), 3*max( mm( j ), nn( j ) ) ) )
507  10 CONTINUE
508 *
509 * Check for errors
510 *
511  IF( nsizes.LT.0 ) THEN
512  info = -1
513  ELSE IF( badmm ) THEN
514  info = -2
515  ELSE IF( badnn ) THEN
516  info = -3
517  ELSE IF( ntypes.LT.0 ) THEN
518  info = -4
519  ELSE IF( lda.LT.max( 1, mmax ) ) THEN
520  info = -10
521  ELSE IF( ldu.LT.max( 1, mmax ) ) THEN
522  info = -12
523  ELSE IF( ldvt.LT.max( 1, nmax ) ) THEN
524  info = -14
525  ELSE IF( minwrk.GT.lwork ) THEN
526  info = -21
527  END IF
528 *
529  IF( info.NE.0 ) THEN
530  CALL xerbla( 'CDRVBD', -info )
531  RETURN
532  END IF
533 *
534 * Quick return if nothing to do
535 *
536  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
537  $ RETURN
538 *
539 * More Important constants
540 *
541  unfl = slamch( 'S' )
542  ovfl = one / unfl
543  ulp = slamch( 'E' )
544  ulpinv = one / ulp
545  rtunfl = sqrt( unfl )
546 *
547 * Loop over sizes, types
548 *
549  nerrs = 0
550 *
551  DO 310 jsize = 1, nsizes
552  m = mm( jsize )
553  n = nn( jsize )
554  mnmin = min( m, n )
555 *
556  IF( nsizes.NE.1 ) THEN
557  mtypes = min( maxtyp, ntypes )
558  ELSE
559  mtypes = min( maxtyp+1, ntypes )
560  END IF
561 *
562  DO 300 jtype = 1, mtypes
563  IF( .NOT.dotype( jtype ) )
564  $ GO TO 300
565  ntest = 0
566 *
567  DO 20 j = 1, 4
568  ioldsd( j ) = iseed( j )
569  20 CONTINUE
570 *
571 * Compute "A"
572 *
573  IF( mtypes.GT.maxtyp )
574  $ GO TO 50
575 *
576  IF( jtype.EQ.1 ) THEN
577 *
578 * Zero matrix
579 *
580  CALL claset( 'Full', m, n, czero, czero, a, lda )
581  DO 30 i = 1, min( m, n )
582  s( i ) = zero
583  30 CONTINUE
584 *
585  ELSE IF( jtype.EQ.2 ) THEN
586 *
587 * Identity matrix
588 *
589  CALL claset( 'Full', m, n, czero, cone, a, lda )
590  DO 40 i = 1, min( m, n )
591  s( i ) = one
592  40 CONTINUE
593 *
594  ELSE
595 *
596 * (Scaled) random matrix
597 *
598  IF( jtype.EQ.3 )
599  $ anorm = one
600  IF( jtype.EQ.4 )
601  $ anorm = unfl / ulp
602  IF( jtype.EQ.5 )
603  $ anorm = ovfl*ulp
604  CALL clatms( m, n, 'U', iseed, 'N', s, 4, real( mnmin ),
605  $ anorm, m-1, n-1, 'N', a, lda, work, iinfo )
606  IF( iinfo.NE.0 ) THEN
607  WRITE( nounit, fmt = 9996 )'Generator', iinfo, m, n,
608  $ jtype, ioldsd
609  info = abs( iinfo )
610  RETURN
611  END IF
612  END IF
613 *
614  50 CONTINUE
615  CALL clacpy( 'F', m, n, a, lda, asav, lda )
616 *
617 * Do for minimal and adequate (for blocking) workspace
618 *
619  DO 290 iwspc = 1, 4
620 *
621 * Test for CGESVD
622 *
623  iwtmp = 2*min( m, n )+max( m, n )
624  lswork = iwtmp + ( iwspc-1 )*( lwork-iwtmp ) / 3
625  lswork = min( lswork, lwork )
626  lswork = max( lswork, 1 )
627  IF( iwspc.EQ.4 )
628  $ lswork = lwork
629 *
630  DO 60 j = 1, 35
631  result( j ) = -one
632  60 CONTINUE
633 *
634 * Factorize A
635 *
636  IF( iwspc.GT.1 )
637  $ CALL clacpy( 'F', m, n, asav, lda, a, lda )
638  srnamt = 'CGESVD'
639  CALL cgesvd( 'A', 'A', m, n, a, lda, ssav, usav, ldu,
640  $ vtsav, ldvt, work, lswork, rwork, iinfo )
641  IF( iinfo.NE.0 ) THEN
642  WRITE( nounit, fmt = 9995 )'GESVD', iinfo, m, n,
643  $ jtype, lswork, ioldsd
644  info = abs( iinfo )
645  RETURN
646  END IF
647 *
648 * Do tests 1--4
649 *
650  CALL cbdt01( m, n, 0, asav, lda, usav, ldu, ssav, e,
651  $ vtsav, ldvt, work, rwork, result( 1 ) )
652  IF( m.NE.0 .AND. n.NE.0 ) THEN
653  CALL cunt01( 'Columns', mnmin, m, usav, ldu, work,
654  $ lwork, rwork, result( 2 ) )
655  CALL cunt01( 'Rows', mnmin, n, vtsav, ldvt, work,
656  $ lwork, rwork, result( 3 ) )
657  END IF
658  result( 4 ) = 0
659  DO 70 i = 1, mnmin - 1
660  IF( ssav( i ).LT.ssav( i+1 ) )
661  $ result( 4 ) = ulpinv
662  IF( ssav( i ).LT.zero )
663  $ result( 4 ) = ulpinv
664  70 CONTINUE
665  IF( mnmin.GE.1 ) THEN
666  IF( ssav( mnmin ).LT.zero )
667  $ result( 4 ) = ulpinv
668  END IF
669 *
670 * Do partial SVDs, comparing to SSAV, USAV, and VTSAV
671 *
672  result( 5 ) = zero
673  result( 6 ) = zero
674  result( 7 ) = zero
675  DO 100 iju = 0, 3
676  DO 90 ijvt = 0, 3
677  IF( ( iju.EQ.3 .AND. ijvt.EQ.3 ) .OR.
678  $ ( iju.EQ.1 .AND. ijvt.EQ.1 ) )GO TO 90
679  jobu = cjob( iju+1 )
680  jobvt = cjob( ijvt+1 )
681  CALL clacpy( 'F', m, n, asav, lda, a, lda )
682  srnamt = 'CGESVD'
683  CALL cgesvd( jobu, jobvt, m, n, a, lda, s, u, ldu,
684  $ vt, ldvt, work, lswork, rwork, iinfo )
685 *
686 * Compare U
687 *
688  dif = zero
689  IF( m.GT.0 .AND. n.GT.0 ) THEN
690  IF( iju.EQ.1 ) THEN
691  CALL cunt03( 'C', m, mnmin, m, mnmin, usav,
692  $ ldu, a, lda, work, lwork, rwork,
693  $ dif, iinfo )
694  ELSE IF( iju.EQ.2 ) THEN
695  CALL cunt03( 'C', m, mnmin, m, mnmin, usav,
696  $ ldu, u, ldu, work, lwork, rwork,
697  $ dif, iinfo )
698  ELSE IF( iju.EQ.3 ) THEN
699  CALL cunt03( 'C', m, m, m, mnmin, usav, ldu,
700  $ u, ldu, work, lwork, rwork, dif,
701  $ iinfo )
702  END IF
703  END IF
704  result( 5 ) = max( result( 5 ), dif )
705 *
706 * Compare VT
707 *
708  dif = zero
709  IF( m.GT.0 .AND. n.GT.0 ) THEN
710  IF( ijvt.EQ.1 ) THEN
711  CALL cunt03( 'R', n, mnmin, n, mnmin, vtsav,
712  $ ldvt, a, lda, work, lwork,
713  $ rwork, dif, iinfo )
714  ELSE IF( ijvt.EQ.2 ) THEN
715  CALL cunt03( 'R', n, mnmin, n, mnmin, vtsav,
716  $ ldvt, vt, ldvt, work, lwork,
717  $ rwork, dif, iinfo )
718  ELSE IF( ijvt.EQ.3 ) THEN
719  CALL cunt03( 'R', n, n, n, mnmin, vtsav,
720  $ ldvt, vt, ldvt, work, lwork,
721  $ rwork, dif, iinfo )
722  END IF
723  END IF
724  result( 6 ) = max( result( 6 ), dif )
725 *
726 * Compare S
727 *
728  dif = zero
729  div = max( real( mnmin )*ulp*s( 1 ),
730  $ slamch( 'Safe minimum' ) )
731  DO 80 i = 1, mnmin - 1
732  IF( ssav( i ).LT.ssav( i+1 ) )
733  $ dif = ulpinv
734  IF( ssav( i ).LT.zero )
735  $ dif = ulpinv
736  dif = max( dif, abs( ssav( i )-s( i ) ) / div )
737  80 CONTINUE
738  result( 7 ) = max( result( 7 ), dif )
739  90 CONTINUE
740  100 CONTINUE
741 *
742 * Test for CGESDD
743 *
744  iwtmp = 2*mnmin*mnmin + 2*mnmin + max( m, n )
745  lswork = iwtmp + ( iwspc-1 )*( lwork-iwtmp ) / 3
746  lswork = min( lswork, lwork )
747  lswork = max( lswork, 1 )
748  IF( iwspc.EQ.4 )
749  $ lswork = lwork
750 *
751 * Factorize A
752 *
753  CALL clacpy( 'F', m, n, asav, lda, a, lda )
754  srnamt = 'CGESDD'
755  CALL cgesdd( 'A', m, n, a, lda, ssav, usav, ldu, vtsav,
756  $ ldvt, work, lswork, rwork, iwork, iinfo )
757  IF( iinfo.NE.0 ) THEN
758  WRITE( nounit, fmt = 9995 )'GESDD', iinfo, m, n,
759  $ jtype, lswork, ioldsd
760  info = abs( iinfo )
761  RETURN
762  END IF
763 *
764 * Do tests 1--4
765 *
766  CALL cbdt01( m, n, 0, asav, lda, usav, ldu, ssav, e,
767  $ vtsav, ldvt, work, rwork, result( 8 ) )
768  IF( m.NE.0 .AND. n.NE.0 ) THEN
769  CALL cunt01( 'Columns', mnmin, m, usav, ldu, work,
770  $ lwork, rwork, result( 9 ) )
771  CALL cunt01( 'Rows', mnmin, n, vtsav, ldvt, work,
772  $ lwork, rwork, result( 10 ) )
773  END IF
774  result( 11 ) = 0
775  DO 110 i = 1, mnmin - 1
776  IF( ssav( i ).LT.ssav( i+1 ) )
777  $ result( 11 ) = ulpinv
778  IF( ssav( i ).LT.zero )
779  $ result( 11 ) = ulpinv
780  110 CONTINUE
781  IF( mnmin.GE.1 ) THEN
782  IF( ssav( mnmin ).LT.zero )
783  $ result( 11 ) = ulpinv
784  END IF
785 *
786 * Do partial SVDs, comparing to SSAV, USAV, and VTSAV
787 *
788  result( 12 ) = zero
789  result( 13 ) = zero
790  result( 14 ) = zero
791  DO 130 ijq = 0, 2
792  jobq = cjob( ijq+1 )
793  CALL clacpy( 'F', m, n, asav, lda, a, lda )
794  srnamt = 'CGESDD'
795  CALL cgesdd( jobq, m, n, a, lda, s, u, ldu, vt, ldvt,
796  $ work, lswork, rwork, iwork, iinfo )
797 *
798 * Compare U
799 *
800  dif = zero
801  IF( m.GT.0 .AND. n.GT.0 ) THEN
802  IF( ijq.EQ.1 ) THEN
803  IF( m.GE.n ) THEN
804  CALL cunt03( 'C', m, mnmin, m, mnmin, usav,
805  $ ldu, a, lda, work, lwork, rwork,
806  $ dif, iinfo )
807  ELSE
808  CALL cunt03( 'C', m, mnmin, m, mnmin, usav,
809  $ ldu, u, ldu, work, lwork, rwork,
810  $ dif, iinfo )
811  END IF
812  ELSE IF( ijq.EQ.2 ) THEN
813  CALL cunt03( 'C', m, mnmin, m, mnmin, usav, ldu,
814  $ u, ldu, work, lwork, rwork, dif,
815  $ iinfo )
816  END IF
817  END IF
818  result( 12 ) = max( result( 12 ), dif )
819 *
820 * Compare VT
821 *
822  dif = zero
823  IF( m.GT.0 .AND. n.GT.0 ) THEN
824  IF( ijq.EQ.1 ) THEN
825  IF( m.GE.n ) THEN
826  CALL cunt03( 'R', n, mnmin, n, mnmin, vtsav,
827  $ ldvt, vt, ldvt, work, lwork,
828  $ rwork, dif, iinfo )
829  ELSE
830  CALL cunt03( 'R', n, mnmin, n, mnmin, vtsav,
831  $ ldvt, a, lda, work, lwork,
832  $ rwork, dif, iinfo )
833  END IF
834  ELSE IF( ijq.EQ.2 ) THEN
835  CALL cunt03( 'R', n, mnmin, n, mnmin, vtsav,
836  $ ldvt, vt, ldvt, work, lwork, rwork,
837  $ dif, iinfo )
838  END IF
839  END IF
840  result( 13 ) = max( result( 13 ), dif )
841 *
842 * Compare S
843 *
844  dif = zero
845  div = max( real( mnmin )*ulp*s( 1 ),
846  $ slamch( 'Safe minimum' ) )
847  DO 120 i = 1, mnmin - 1
848  IF( ssav( i ).LT.ssav( i+1 ) )
849  $ dif = ulpinv
850  IF( ssav( i ).LT.zero )
851  $ dif = ulpinv
852  dif = max( dif, abs( ssav( i )-s( i ) ) / div )
853  120 CONTINUE
854  result( 14 ) = max( result( 14 ), dif )
855  130 CONTINUE
856 
857 *
858 * Test CGESVDQ
859 * Note: CGESVDQ only works for M >= N
860 *
861  result( 36 ) = zero
862  result( 37 ) = zero
863  result( 38 ) = zero
864  result( 39 ) = zero
865 *
866  IF( m.GE.n ) THEN
867  iwtmp = 2*mnmin*mnmin + 2*mnmin + max( m, n )
868  lswork = iwtmp + ( iwspc-1 )*( lwork-iwtmp ) / 3
869  lswork = min( lswork, lwork )
870  lswork = max( lswork, 1 )
871  IF( iwspc.EQ.4 )
872  $ lswork = lwork
873 *
874  CALL clacpy( 'F', m, n, asav, lda, a, lda )
875  srnamt = 'CGESVDQ'
876 *
877  lrwork = max(2, m, 5*n)
878  liwork = max( n, 1 )
879  CALL cgesvdq( 'H', 'N', 'N', 'A', 'A',
880  $ m, n, a, lda, ssav, usav, ldu,
881  $ vtsav, ldvt, numrank, iwork, liwork,
882  $ work, lwork, rwork, lrwork, iinfo )
883 *
884  IF( iinfo.NE.0 ) THEN
885  WRITE( nounit, fmt = 9995 )'CGESVDQ', iinfo, m, n,
886  $ jtype, lswork, ioldsd
887  info = abs( iinfo )
888  RETURN
889  END IF
890 *
891 * Do tests 36--39
892 *
893  CALL cbdt01( m, n, 0, asav, lda, usav, ldu, ssav, e,
894  $ vtsav, ldvt, work, rwork, result( 36 ) )
895  IF( m.NE.0 .AND. n.NE.0 ) THEN
896  CALL cunt01( 'Columns', m, m, usav, ldu, work,
897  $ lwork, rwork, result( 37 ) )
898  CALL cunt01( 'Rows', n, n, vtsav, ldvt, work,
899  $ lwork, rwork, result( 38 ) )
900  END IF
901  result( 39 ) = zero
902  DO 199 i = 1, mnmin - 1
903  IF( ssav( i ).LT.ssav( i+1 ) )
904  $ result( 39 ) = ulpinv
905  IF( ssav( i ).LT.zero )
906  $ result( 39 ) = ulpinv
907  199 CONTINUE
908  IF( mnmin.GE.1 ) THEN
909  IF( ssav( mnmin ).LT.zero )
910  $ result( 39 ) = ulpinv
911  END IF
912  END IF
913 *
914 * Test CGESVJ
915 * Note: CGESVJ only works for M >= N
916 *
917  result( 15 ) = zero
918  result( 16 ) = zero
919  result( 17 ) = zero
920  result( 18 ) = zero
921 *
922  IF( m.GE.n ) THEN
923  iwtmp = 2*mnmin*mnmin + 2*mnmin + max( m, n )
924  lswork = iwtmp + ( iwspc-1 )*( lwork-iwtmp ) / 3
925  lswork = min( lswork, lwork )
926  lswork = max( lswork, 1 )
927  lrwork = max(6,n)
928  IF( iwspc.EQ.4 )
929  $ lswork = lwork
930 *
931  CALL clacpy( 'F', m, n, asav, lda, usav, lda )
932  srnamt = 'CGESVJ'
933  CALL cgesvj( 'G', 'U', 'V', m, n, usav, lda, ssav,
934  & 0, a, ldvt, work, lwork, rwork,
935  & lrwork, iinfo )
936 *
937 * CGESVJ returns V not VH
938 *
939  DO j=1,n
940  DO i=1,n
941  vtsav(j,i) = conjg(a(i,j))
942  END DO
943  END DO
944 *
945  IF( iinfo.NE.0 ) THEN
946  WRITE( nounit, fmt = 9995 )'GESVJ', iinfo, m, n,
947  $ jtype, lswork, ioldsd
948  info = abs( iinfo )
949  RETURN
950  END IF
951 *
952 * Do tests 15--18
953 *
954  CALL cbdt01( m, n, 0, asav, lda, usav, ldu, ssav, e,
955  $ vtsav, ldvt, work, rwork, result( 15 ) )
956  IF( m.NE.0 .AND. n.NE.0 ) THEN
957  CALL cunt01( 'Columns', m, m, usav, ldu, work,
958  $ lwork, rwork, result( 16 ) )
959  CALL cunt01( 'Rows', n, n, vtsav, ldvt, work,
960  $ lwork, rwork, result( 17 ) )
961  END IF
962  result( 18 ) = zero
963  DO 131 i = 1, mnmin - 1
964  IF( ssav( i ).LT.ssav( i+1 ) )
965  $ result( 18 ) = ulpinv
966  IF( ssav( i ).LT.zero )
967  $ result( 18 ) = ulpinv
968  131 CONTINUE
969  IF( mnmin.GE.1 ) THEN
970  IF( ssav( mnmin ).LT.zero )
971  $ result( 18 ) = ulpinv
972  END IF
973  END IF
974 *
975 * Test CGEJSV
976 * Note: CGEJSV only works for M >= N
977 *
978  result( 19 ) = zero
979  result( 20 ) = zero
980  result( 21 ) = zero
981  result( 22 ) = zero
982  IF( m.GE.n ) THEN
983  iwtmp = 2*mnmin*mnmin + 2*mnmin + max( m, n )
984  lswork = iwtmp + ( iwspc-1 )*( lwork-iwtmp ) / 3
985  lswork = min( lswork, lwork )
986  lswork = max( lswork, 1 )
987  IF( iwspc.EQ.4 )
988  $ lswork = lwork
989  lrwork = max( 7, n + 2*m)
990 *
991  CALL clacpy( 'F', m, n, asav, lda, vtsav, lda )
992  srnamt = 'CGEJSV'
993  CALL cgejsv( 'G', 'U', 'V', 'R', 'N', 'N',
994  & m, n, vtsav, lda, ssav, usav, ldu, a, ldvt,
995  & work, lwork, rwork,
996  & lrwork, iwork, iinfo )
997 *
998 * CGEJSV returns V not VH
999 *
1000  DO 133 j=1,n
1001  DO 132 i=1,n
1002  vtsav(j,i) = conjg(a(i,j))
1003  132 END DO
1004  133 END DO
1005 *
1006  IF( iinfo.NE.0 ) THEN
1007  WRITE( nounit, fmt = 9995 )'GEJSV', iinfo, m, n,
1008  $ jtype, lswork, ioldsd
1009  info = abs( iinfo )
1010  RETURN
1011  END IF
1012 *
1013 * Do tests 19--22
1014 *
1015  CALL cbdt01( m, n, 0, asav, lda, usav, ldu, ssav, e,
1016  $ vtsav, ldvt, work, rwork, result( 19 ) )
1017  IF( m.NE.0 .AND. n.NE.0 ) THEN
1018  CALL cunt01( 'Columns', m, m, usav, ldu, work,
1019  $ lwork, rwork, result( 20 ) )
1020  CALL cunt01( 'Rows', n, n, vtsav, ldvt, work,
1021  $ lwork, rwork, result( 21 ) )
1022  END IF
1023  result( 22 ) = zero
1024  DO 134 i = 1, mnmin - 1
1025  IF( ssav( i ).LT.ssav( i+1 ) )
1026  $ result( 22 ) = ulpinv
1027  IF( ssav( i ).LT.zero )
1028  $ result( 22 ) = ulpinv
1029  134 CONTINUE
1030  IF( mnmin.GE.1 ) THEN
1031  IF( ssav( mnmin ).LT.zero )
1032  $ result( 22 ) = ulpinv
1033  END IF
1034  END IF
1035 *
1036 * Test CGESVDX
1037 *
1038 * Factorize A
1039 *
1040  CALL clacpy( 'F', m, n, asav, lda, a, lda )
1041  srnamt = 'CGESVDX'
1042  CALL cgesvdx( 'V', 'V', 'A', m, n, a, lda,
1043  $ vl, vu, il, iu, ns, ssav, usav, ldu,
1044  $ vtsav, ldvt, work, lwork, rwork,
1045  $ iwork, iinfo )
1046  IF( iinfo.NE.0 ) THEN
1047  WRITE( nounit, fmt = 9995 )'GESVDX', iinfo, m, n,
1048  $ jtype, lswork, ioldsd
1049  info = abs( iinfo )
1050  RETURN
1051  END IF
1052 *
1053 * Do tests 1--4
1054 *
1055  result( 23 ) = zero
1056  result( 24 ) = zero
1057  result( 25 ) = zero
1058  CALL cbdt01( m, n, 0, asav, lda, usav, ldu, ssav, e,
1059  $ vtsav, ldvt, work, rwork, result( 23 ) )
1060  IF( m.NE.0 .AND. n.NE.0 ) THEN
1061  CALL cunt01( 'Columns', mnmin, m, usav, ldu, work,
1062  $ lwork, rwork, result( 24 ) )
1063  CALL cunt01( 'Rows', mnmin, n, vtsav, ldvt, work,
1064  $ lwork, rwork, result( 25 ) )
1065  END IF
1066  result( 26 ) = zero
1067  DO 140 i = 1, mnmin - 1
1068  IF( ssav( i ).LT.ssav( i+1 ) )
1069  $ result( 26 ) = ulpinv
1070  IF( ssav( i ).LT.zero )
1071  $ result( 26 ) = ulpinv
1072  140 CONTINUE
1073  IF( mnmin.GE.1 ) THEN
1074  IF( ssav( mnmin ).LT.zero )
1075  $ result( 26 ) = ulpinv
1076  END IF
1077 *
1078 * Do partial SVDs, comparing to SSAV, USAV, and VTSAV
1079 *
1080  result( 27 ) = zero
1081  result( 28 ) = zero
1082  result( 29 ) = zero
1083  DO 170 iju = 0, 1
1084  DO 160 ijvt = 0, 1
1085  IF( ( iju.EQ.0 .AND. ijvt.EQ.0 ) .OR.
1086  $ ( iju.EQ.1 .AND. ijvt.EQ.1 ) ) GO TO 160
1087  jobu = cjobv( iju+1 )
1088  jobvt = cjobv( ijvt+1 )
1089  range = cjobr( 1 )
1090  CALL clacpy( 'F', m, n, asav, lda, a, lda )
1091  srnamt = 'CGESVDX'
1092  CALL cgesvdx( jobu, jobvt, 'A', m, n, a, lda,
1093  $ vl, vu, il, iu, ns, ssav, u, ldu,
1094  $ vt, ldvt, work, lwork, rwork,
1095  $ iwork, iinfo )
1096 *
1097 * Compare U
1098 *
1099  dif = zero
1100  IF( m.GT.0 .AND. n.GT.0 ) THEN
1101  IF( iju.EQ.1 ) THEN
1102  CALL cunt03( 'C', m, mnmin, m, mnmin, usav,
1103  $ ldu, u, ldu, work, lwork, rwork,
1104  $ dif, iinfo )
1105  END IF
1106  END IF
1107  result( 27 ) = max( result( 27 ), dif )
1108 *
1109 * Compare VT
1110 *
1111  dif = zero
1112  IF( m.GT.0 .AND. n.GT.0 ) THEN
1113  IF( ijvt.EQ.1 ) THEN
1114  CALL cunt03( 'R', n, mnmin, n, mnmin, vtsav,
1115  $ ldvt, vt, ldvt, work, lwork,
1116  $ rwork, dif, iinfo )
1117  END IF
1118  END IF
1119  result( 28 ) = max( result( 28 ), dif )
1120 *
1121 * Compare S
1122 *
1123  dif = zero
1124  div = max( real( mnmin )*ulp*s( 1 ),
1125  $ slamch( 'Safe minimum' ) )
1126  DO 150 i = 1, mnmin - 1
1127  IF( ssav( i ).LT.ssav( i+1 ) )
1128  $ dif = ulpinv
1129  IF( ssav( i ).LT.zero )
1130  $ dif = ulpinv
1131  dif = max( dif, abs( ssav( i )-s( i ) ) / div )
1132  150 CONTINUE
1133  result( 29) = max( result( 29 ), dif )
1134  160 CONTINUE
1135  170 CONTINUE
1136 *
1137 * Do tests 8--10
1138 *
1139  DO 180 i = 1, 4
1140  iseed2( i ) = iseed( i )
1141  180 CONTINUE
1142  IF( mnmin.LE.1 ) THEN
1143  il = 1
1144  iu = max( 1, mnmin )
1145  ELSE
1146  il = 1 + int( ( mnmin-1 )*slarnd( 1, iseed2 ) )
1147  iu = 1 + int( ( mnmin-1 )*slarnd( 1, iseed2 ) )
1148  IF( iu.LT.il ) THEN
1149  itemp = iu
1150  iu = il
1151  il = itemp
1152  END IF
1153  END IF
1154  CALL clacpy( 'F', m, n, asav, lda, a, lda )
1155  srnamt = 'CGESVDX'
1156  CALL cgesvdx( 'V', 'V', 'I', m, n, a, lda,
1157  $ vl, vu, il, iu, nsi, s, u, ldu,
1158  $ vt, ldvt, work, lwork, rwork,
1159  $ iwork, iinfo )
1160  IF( iinfo.NE.0 ) THEN
1161  WRITE( nounit, fmt = 9995 )'GESVDX', iinfo, m, n,
1162  $ jtype, lswork, ioldsd
1163  info = abs( iinfo )
1164  RETURN
1165  END IF
1166 *
1167  result( 30 ) = zero
1168  result( 31 ) = zero
1169  result( 32 ) = zero
1170  CALL cbdt05( m, n, asav, lda, s, nsi, u, ldu,
1171  $ vt, ldvt, work, result( 30 ) )
1172  IF( m.NE.0 .AND. n.NE.0 ) THEN
1173  CALL cunt01( 'Columns', m, nsi, u, ldu, work,
1174  $ lwork, rwork, result( 31 ) )
1175  CALL cunt01( 'Rows', nsi, n, vt, ldvt, work,
1176  $ lwork, rwork, result( 32 ) )
1177  END IF
1178 *
1179 * Do tests 11--13
1180 *
1181  IF( mnmin.GT.0 .AND. nsi.GT.1 ) THEN
1182  IF( il.NE.1 ) THEN
1183  vu = ssav( il ) +
1184  $ max( half*abs( ssav( il )-ssav( il-1 ) ),
1185  $ ulp*anorm, two*rtunfl )
1186  ELSE
1187  vu = ssav( 1 ) +
1188  $ max( half*abs( ssav( ns )-ssav( 1 ) ),
1189  $ ulp*anorm, two*rtunfl )
1190  END IF
1191  IF( iu.NE.ns ) THEN
1192  vl = ssav( iu ) - max( ulp*anorm, two*rtunfl,
1193  $ half*abs( ssav( iu+1 )-ssav( iu ) ) )
1194  ELSE
1195  vl = ssav( ns ) - max( ulp*anorm, two*rtunfl,
1196  $ half*abs( ssav( ns )-ssav( 1 ) ) )
1197  END IF
1198  vl = max( vl,zero )
1199  vu = max( vu,zero )
1200  IF( vl.GE.vu ) vu = max( vu*2, vu+vl+half )
1201  ELSE
1202  vl = zero
1203  vu = one
1204  END IF
1205  CALL clacpy( 'F', m, n, asav, lda, a, lda )
1206  srnamt = 'CGESVDX'
1207  CALL cgesvdx( 'V', 'V', 'V', m, n, a, lda,
1208  $ vl, vu, il, iu, nsv, s, u, ldu,
1209  $ vt, ldvt, work, lwork, rwork,
1210  $ iwork, iinfo )
1211  IF( iinfo.NE.0 ) THEN
1212  WRITE( nounit, fmt = 9995 )'GESVDX', iinfo, m, n,
1213  $ jtype, lswork, ioldsd
1214  info = abs( iinfo )
1215  RETURN
1216  END IF
1217 *
1218  result( 33 ) = zero
1219  result( 34 ) = zero
1220  result( 35 ) = zero
1221  CALL cbdt05( m, n, asav, lda, s, nsv, u, ldu,
1222  $ vt, ldvt, work, result( 33 ) )
1223  IF( m.NE.0 .AND. n.NE.0 ) THEN
1224  CALL cunt01( 'Columns', m, nsv, u, ldu, work,
1225  $ lwork, rwork, result( 34 ) )
1226  CALL cunt01( 'Rows', nsv, n, vt, ldvt, work,
1227  $ lwork, rwork, result( 35 ) )
1228  END IF
1229 *
1230 * End of Loop -- Check for RESULT(j) > THRESH
1231 *
1232  ntest = 0
1233  nfail = 0
1234  DO 190 j = 1, 39
1235  IF( result( j ).GE.zero )
1236  $ ntest = ntest + 1
1237  IF( result( j ).GE.thresh )
1238  $ nfail = nfail + 1
1239  190 CONTINUE
1240 *
1241  IF( nfail.GT.0 )
1242  $ ntestf = ntestf + 1
1243  IF( ntestf.EQ.1 ) THEN
1244  WRITE( nounit, fmt = 9999 )
1245  WRITE( nounit, fmt = 9998 )thresh
1246  ntestf = 2
1247  END IF
1248 *
1249  DO 200 j = 1, 39
1250  IF( result( j ).GE.thresh ) THEN
1251  WRITE( nounit, fmt = 9997 )m, n, jtype, iwspc,
1252  $ ioldsd, j, result( j )
1253  END IF
1254  200 CONTINUE
1255 *
1256  nerrs = nerrs + nfail
1257  ntestt = ntestt + ntest
1258 *
1259  290 CONTINUE
1260 *
1261  300 CONTINUE
1262  310 CONTINUE
1263 *
1264 * Summary
1265 *
1266  CALL alasvm( 'CBD', nounit, nerrs, ntestt, 0 )
1267 *
1268  9999 FORMAT( ' SVD -- Complex Singular Value Decomposition Driver ',
1269  $ / ' Matrix types (see CDRVBD for details):',
1270  $ / / ' 1 = Zero matrix', / ' 2 = Identity matrix',
1271  $ / ' 3 = Evenly spaced singular values near 1',
1272  $ / ' 4 = Evenly spaced singular values near underflow',
1273  $ / ' 5 = Evenly spaced singular values near overflow',
1274  $ / / ' Tests performed: ( A is dense, U and V are unitary,',
1275  $ / 19x, ' S is an array, and Upartial, VTpartial, and',
1276  $ / 19x, ' Spartial are partially computed U, VT and S),', / )
1277  9998 FORMAT( ' Tests performed with Test Threshold = ', f8.2,
1278  $ / ' CGESVD: ', /
1279  $ ' 1 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
1280  $ / ' 2 = | I - U**T U | / ( M ulp ) ',
1281  $ / ' 3 = | I - VT VT**T | / ( N ulp ) ',
1282  $ / ' 4 = 0 if S contains min(M,N) nonnegative values in',
1283  $ ' decreasing order, else 1/ulp',
1284  $ / ' 5 = | U - Upartial | / ( M ulp )',
1285  $ / ' 6 = | VT - VTpartial | / ( N ulp )',
1286  $ / ' 7 = | S - Spartial | / ( min(M,N) ulp |S| )',
1287  $ / ' CGESDD: ', /
1288  $ ' 8 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
1289  $ / ' 9 = | I - U**T U | / ( M ulp ) ',
1290  $ / '10 = | I - VT VT**T | / ( N ulp ) ',
1291  $ / '11 = 0 if S contains min(M,N) nonnegative values in',
1292  $ ' decreasing order, else 1/ulp',
1293  $ / '12 = | U - Upartial | / ( M ulp )',
1294  $ / '13 = | VT - VTpartial | / ( N ulp )',
1295  $ / '14 = | S - Spartial | / ( min(M,N) ulp |S| )',
1296  $ / ' CGESVJ: ', /
1297  $ / '15 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
1298  $ / '16 = | I - U**T U | / ( M ulp ) ',
1299  $ / '17 = | I - VT VT**T | / ( N ulp ) ',
1300  $ / '18 = 0 if S contains min(M,N) nonnegative values in',
1301  $ ' decreasing order, else 1/ulp',
1302  $ / ' CGESJV: ', /
1303  $ / '19 = | A - U diag(S) VT | / ( |A| max(M,N) ulp )',
1304  $ / '20 = | I - U**T U | / ( M ulp ) ',
1305  $ / '21 = | I - VT VT**T | / ( N ulp ) ',
1306  $ / '22 = 0 if S contains min(M,N) nonnegative values in',
1307  $ ' decreasing order, else 1/ulp',
1308  $ / ' CGESVDX(V,V,A): ', /
1309  $ '23 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
1310  $ / '24 = | I - U**T U | / ( M ulp ) ',
1311  $ / '25 = | I - VT VT**T | / ( N ulp ) ',
1312  $ / '26 = 0 if S contains min(M,N) nonnegative values in',
1313  $ ' decreasing order, else 1/ulp',
1314  $ / '27 = | U - Upartial | / ( M ulp )',
1315  $ / '28 = | VT - VTpartial | / ( N ulp )',
1316  $ / '29 = | S - Spartial | / ( min(M,N) ulp |S| )',
1317  $ / ' CGESVDX(V,V,I): ',
1318  $ / '30 = | U**T A VT**T - diag(S) | / ( |A| max(M,N) ulp )',
1319  $ / '31 = | I - U**T U | / ( M ulp ) ',
1320  $ / '32 = | I - VT VT**T | / ( N ulp ) ',
1321  $ / ' CGESVDX(V,V,V) ',
1322  $ / '33 = | U**T A VT**T - diag(S) | / ( |A| max(M,N) ulp )',
1323  $ / '34 = | I - U**T U | / ( M ulp ) ',
1324  $ / '35 = | I - VT VT**T | / ( N ulp ) ',
1325  $ ' CGESVDQ(H,N,N,A,A',
1326  $ / '36 = | A - U diag(S) VT | / ( |A| max(M,N) ulp ) ',
1327  $ / '37 = | I - U**T U | / ( M ulp ) ',
1328  $ / '38 = | I - VT VT**T | / ( N ulp ) ',
1329  $ / '39 = 0 if S contains min(M,N) nonnegative values in',
1330  $ ' decreasing order, else 1/ulp',
1331  $ / / )
1332  9997 FORMAT( ' M=', i5, ', N=', i5, ', type ', i1, ', IWS=', i1,
1333  $ ', seed=', 4( i4, ',' ), ' test(', i2, ')=', g11.4 )
1334  9996 FORMAT( ' CDRVBD: ', a, ' returned INFO=', i6, '.', / 9x, 'M=',
1335  $ i6, ', N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ),
1336  $ i5, ')' )
1337  9995 FORMAT( ' CDRVBD: ', a, ' returned INFO=', i6, '.', / 9x, 'M=',
1338  $ i6, ', N=', i6, ', JTYPE=', i6, ', LSWORK=', i6, / 9x,
1339  $ 'ISEED=(', 3( i5, ',' ), i5, ')' )
1340 *
1341  RETURN
1342 *
1343 * End of CDRVBD
1344 *
1345  END
alasvm
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:75
cunt03
subroutine cunt03(RC, MU, MV, N, K, U, LDU, V, LDV, WORK, LWORK, RWORK, RESULT, INFO)
CUNT03
Definition: cunt03.f:164
cdrvbd
subroutine cdrvbd(NSIZES, MM, NN, NTYPES, DOTYPE, ISEED, THRESH, A, LDA, U, LDU, VT, LDVT, ASAV, USAV, VTSAV, S, SSAV, E, WORK, LWORK, RWORK, IWORK, NOUNIT, INFO)
CDRVBD
Definition: cdrvbd.f:403
clacpy
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
cbdt01
subroutine cbdt01(M, N, KD, A, LDA, Q, LDQ, D, E, PT, LDPT, WORK, RWORK, RESID)
CBDT01
Definition: cbdt01.f:148
cbdt05
subroutine cbdt05(M, N, A, LDA, S, NS, U, LDU, VT, LDVT, WORK, RESID)
CBDT05
Definition: cbdt05.f:127
cgesdd
subroutine cgesdd(JOBZ, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, RWORK, IWORK, INFO)
CGESDD
Definition: cgesdd.f:228
cgesvdx
subroutine cgesvdx(JOBU, JOBVT, RANGE, M, N, A, LDA, VL, VU, IL, IU, NS, S, U, LDU, VT, LDVT, WORK, LWORK, RWORK, IWORK, INFO)
CGESVDX computes the singular value decomposition (SVD) for GE matrices
Definition: cgesvdx.f:272
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
cgesvd
subroutine cgesvd(JOBU, JOBVT, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, RWORK, INFO)
CGESVD computes the singular value decomposition (SVD) for GE matrices
Definition: cgesvd.f:216
clatms
subroutine clatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
CLATMS
Definition: clatms.f:334
cunt01
subroutine cunt01(ROWCOL, M, N, U, LDU, WORK, LWORK, RWORK, RESID)
CUNT01
Definition: cunt01.f:128
claset
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.
Definition: claset.f:108
cgesvdq
subroutine cgesvdq(JOBA, JOBP, JOBR, JOBU, JOBV, M, N, A, LDA, S, U, LDU, V, LDV, NUMRANK, IWORK, LIWORK, CWORK, LCWORK, RWORK, LRWORK, INFO)
CGESVDQ computes the singular value decomposition (SVD) with a QR-Preconditioned QR SVD Method for GE...
Definition: cgesvdq.f:415
cgesvj
subroutine cgesvj(JOBA, JOBU, JOBV, M, N, A, LDA, SVA, MV, V, LDV, CWORK, LWORK, RWORK, LRWORK, INFO)
CGESVJ
Definition: cgesvj.f:353
cgejsv
subroutine cgejsv(JOBA, JOBU, JOBV, JOBR, JOBT, JOBP, M, N, A, LDA, SVA, U, LDU, V, LDV, CWORK, LWORK, RWORK, LRWORK, IWORK, INFO)
CGEJSV
Definition: cgejsv.f:570