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