LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ sdrvbd()

subroutine sdrvbd ( integer  NSIZES,
integer, dimension( * )  MM,
integer, dimension( * )  NN,
integer  NTYPES,
logical, dimension( * )  DOTYPE,
integer, dimension( 4 )  ISEED,
real  THRESH,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( ldu, * )  U,
integer  LDU,
real, dimension( ldvt, * )  VT,
integer  LDVT,
real, dimension( lda, * )  ASAV,
real, dimension( ldu, * )  USAV,
real, dimension( ldvt, * )  VTSAV,
real, dimension( * )  S,
real, dimension( * )  SSAV,
real, dimension( * )  E,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  NOUT,
integer  INFO 
)

SDRVBD

Purpose:
 SDRVBD checks the singular value decomposition (SVD) drivers
 SGESVD, SGESDD, SGESVDQ, SGESVJ, SGEJSV, and DGESVDX.

 Both SGESVD and SGESDD factor A = U diag(S) VT, where U and VT are
 orthogonal and diag(S) is diagonal with the entries of the array S
 on its diagonal. The entries of S are the singular values,
 nonnegative and stored in decreasing order.  U and VT can be
 optionally not computed, overwritten on A, or computed partially.

 A is M by N. Let MNMIN = min( M, N ). S has dimension MNMIN.
 U can be M by M or M by MNMIN. VT can be N by N or MNMIN by N.

 When SDRVBD is called, a number of matrix "sizes" (M's and N's)
 and a number of matrix "types" are specified.  For each size (M,N)
 and each type of matrix, and for the minimal workspace as well as
 workspace adequate to permit blocking, an  M x N  matrix "A" will be
 generated and used to test the SVD routines.  For each matrix, A will
 be factored as A = U diag(S) VT and the following 12 tests computed:

 Test for SGESVD:

 (1)    | A - U diag(S) VT | / ( |A| max(M,N) ulp )

 (2)    | I - U'U | / ( M ulp )

 (3)    | I - VT VT' | / ( N ulp )

 (4)    S contains MNMIN nonnegative values in decreasing order.
        (Return 0 if true, 1/ULP if false.)

 (5)    | U - Upartial | / ( M ulp ) where Upartial is a partially
        computed U.

 (6)    | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
        computed VT.

 (7)    | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
        vector of singular values from the partial SVD

 Test for SGESDD:

 (8)    | A - U diag(S) VT | / ( |A| max(M,N) ulp )

 (9)    | I - U'U | / ( M ulp )

 (10)   | I - VT VT' | / ( N ulp )

 (11)   S contains MNMIN nonnegative values in decreasing order.
        (Return 0 if true, 1/ULP if false.)

 (12)   | U - Upartial | / ( M ulp ) where Upartial is a partially
        computed U.

 (13)   | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
        computed VT.

 (14)   | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
        vector of singular values from the partial SVD

 Test for SGESVDQ:

 (36)   | A - U diag(S) VT | / ( |A| max(M,N) ulp )

 (37)   | I - U'U | / ( M ulp )

 (38)   | I - VT VT' | / ( N ulp )

 (39)   S contains MNMIN nonnegative values in decreasing order.
        (Return 0 if true, 1/ULP if false.)

 Test for SGESVJ:

 (15)   | A - U diag(S) VT | / ( |A| max(M,N) ulp )

 (16)   | I - U'U | / ( M ulp )

 (17)   | I - VT VT' | / ( N ulp )

 (18)   S contains MNMIN nonnegative values in decreasing order.
        (Return 0 if true, 1/ULP if false.)

 Test for SGEJSV:

 (19)   | A - U diag(S) VT | / ( |A| max(M,N) ulp )

 (20)   | I - U'U | / ( M ulp )

 (21)   | I - VT VT' | / ( N ulp )

 (22)   S contains MNMIN nonnegative values in decreasing order.
        (Return 0 if true, 1/ULP if false.)

 Test for SGESVDX( 'V', 'V', 'A' )/SGESVDX( 'N', 'N', 'A' )

 (23)   | A - U diag(S) VT | / ( |A| max(M,N) ulp )

 (24)   | I - U'U | / ( M ulp )

 (25)   | I - VT VT' | / ( N ulp )

 (26)   S contains MNMIN nonnegative values in decreasing order.
        (Return 0 if true, 1/ULP if false.)

 (27)   | U - Upartial | / ( M ulp ) where Upartial is a partially
        computed U.

 (28)   | VT - VTpartial | / ( N ulp ) where VTpartial is a partially
        computed VT.

 (29)   | S - Spartial | / ( MNMIN ulp |S| ) where Spartial is the
        vector of singular values from the partial SVD

 Test for SGESVDX( 'V', 'V', 'I' )

 (30)   | U' A VT''' - diag(S) | / ( |A| max(M,N) ulp )

 (31)   | I - U'U | / ( M ulp )

 (32)   | I - VT VT' | / ( N ulp )

 Test for SGESVDX( 'V', 'V', 'V' )

 (33)   | U' A VT''' - diag(S) | / ( |A| max(M,N) ulp )

 (34)   | I - U'U | / ( M ulp )

 (35)   | I - VT VT' | / ( N ulp )

 The "sizes" are specified by the arrays MM(1:NSIZES) and
 NN(1:NSIZES); the value of each element pair (MM(j),NN(j))
 specifies one size.  The "types" are specified by a logical array
 DOTYPE( 1:NTYPES ); if DOTYPE(j) is .TRUE., then matrix type "j"
 will be generated.
 Currently, the list of possible types is:

 (1)  The zero matrix.
 (2)  The identity matrix.
 (3)  A matrix of the form  U D V, where U and V are orthogonal and
      D has evenly spaced entries 1, ..., ULP with random signs
      on the diagonal.
 (4)  Same as (3), but multiplied by the underflow-threshold / ULP.
 (5)  Same as (3), but multiplied by the overflow-threshold * ULP.
Parameters
[in]NSIZES
          NSIZES is INTEGER
          The number of matrix sizes (M,N) contained in the vectors
          MM and NN.
[in]MM
          MM is INTEGER array, dimension (NSIZES)
          The values of the matrix row dimension M.
[in]NN
          NN is INTEGER array, dimension (NSIZES)
          The values of the matrix column dimension N.
[in]NTYPES
          NTYPES is INTEGER
          The number of elements in DOTYPE.   If it is zero, SDRVBD
          does nothing.  It must be at least zero.  If it is MAXTYP+1
          and NSIZES is 1, then an additional type, MAXTYP+1 is
          defined, which is to use whatever matrices are in A and B.
          This is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
          DOTYPE(MAXTYP+1) is .TRUE. .
[in]DOTYPE
          DOTYPE is LOGICAL array, dimension (NTYPES)
          If DOTYPE(j) is .TRUE., then for each size (m,n), a matrix
          of type j will be generated.  If NTYPES is smaller than the
          maximum number of types defined (PARAMETER MAXTYP), then
          types NTYPES+1 through MAXTYP will not be generated.  If
          NTYPES is larger than MAXTYP, DOTYPE(MAXTYP+1) through
          DOTYPE(NTYPES) will be ignored.
[in,out]ISEED
          ISEED is INTEGER array, dimension (4)
          On entry, the seed of the random number generator.  The array
          elements should be between 0 and 4095; if not they will be
          reduced mod 4096.  Also, ISEED(4) must be odd.
          On exit, ISEED is changed and can be used in the next call to
          SDRVBD to continue the same random number sequence.
[in]THRESH
          THRESH is REAL
          The threshold value for the test ratios.  A result is
          included in the output file if RESULT >= THRESH.  The test
          ratios are scaled to be O(1), so THRESH should be a small
          multiple of 1, e.g., 10 or 100.  To have every test ratio
          printed, use THRESH = 0.
[out]A
          A is REAL array, dimension (LDA,NMAX)
          where NMAX is the maximum value of N in NN.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,MMAX),
          where MMAX is the maximum value of M in MM.
[out]U
          U is REAL array, dimension (LDU,MMAX)
[in]LDU
          LDU is INTEGER
          The leading dimension of the array U.  LDU >= max(1,MMAX).
[out]VT
          VT is REAL array, dimension (LDVT,NMAX)
[in]LDVT
          LDVT is INTEGER
          The leading dimension of the array VT.  LDVT >= max(1,NMAX).
[out]ASAV
          ASAV is REAL array, dimension (LDA,NMAX)
[out]USAV
          USAV is REAL array, dimension (LDU,MMAX)
[out]VTSAV
          VTSAV is REAL array, dimension (LDVT,NMAX)
[out]S
          S is REAL array, dimension
                      (max(min(MM,NN)))
[out]SSAV
          SSAV is REAL array, dimension
                      (max(min(MM,NN)))
[out]E
          E is REAL array, dimension
                      (max(min(MM,NN)))
[out]WORK
          WORK is REAL array, dimension (LWORK)
[in]LWORK
          LWORK is INTEGER
          The number of entries in WORK.  This must be at least
          max(3*MN+MX,5*MN-4)+2*MN**2 for all pairs
          pairs  (MN,MX)=( min(MM(j),NN(j), max(MM(j),NN(j)) )
[out]IWORK
          IWORK is INTEGER array, dimension at least 8*min(M,N)
[in]NOUT
          NOUT is INTEGER
          The FORTRAN unit number for printing out error messages
          (e.g., if a routine returns IINFO not equal to 0.)
[out]INFO
          INFO is INTEGER
          If 0, then everything ran OK.
           -1: NSIZES < 0
           -2: Some MM(j) < 0
           -3: Some NN(j) < 0
           -4: NTYPES < 0
           -7: THRESH < 0
          -10: LDA < 1 or LDA < MMAX, where MMAX is max( MM(j) ).
          -12: LDU < 1 or LDU < MMAX.
          -14: LDVT < 1 or LDVT < NMAX, where NMAX is max( NN(j) ).
          -21: LWORK too small.
          If  SLATMS, or SGESVD returns an error code, the
              absolute value of it is returned.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 363 of file sdrvbd.f.

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 *
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
real function slarnd(IDIST, ISEED)
SLARND
Definition: slarnd.f:73
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 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
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:68
Here is the call graph for this function:
Here is the caller graph for this function: