LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
zdrvls.f
Go to the documentation of this file.
00001 *> \brief \b ZDRVLS
00002 *
00003 *  =========== DOCUMENTATION ===========
00004 *
00005 * Online html documentation available at 
00006 *            http://www.netlib.org/lapack/explore-html/ 
00007 *
00008 *  Definition:
00009 *  ===========
00010 *
00011 *       SUBROUTINE ZDRVLS( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB,
00012 *                          NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B,
00013 *                          COPYB, C, S, COPYS, WORK, RWORK, IWORK, NOUT )
00014 * 
00015 *       .. Scalar Arguments ..
00016 *       LOGICAL            TSTERR
00017 *       INTEGER            NM, NN, NNB, NNS, NOUT
00018 *       DOUBLE PRECISION   THRESH
00019 *       ..
00020 *       .. Array Arguments ..
00021 *       LOGICAL            DOTYPE( * )
00022 *       INTEGER            IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
00023 *      $                   NVAL( * ), NXVAL( * )
00024 *       DOUBLE PRECISION   COPYS( * ), RWORK( * ), S( * )
00025 *       COMPLEX*16         A( * ), B( * ), C( * ), COPYA( * ), COPYB( * ),
00026 *      $                   WORK( * )
00027 *       ..
00028 *  
00029 *
00030 *> \par Purpose:
00031 *  =============
00032 *>
00033 *> \verbatim
00034 *>
00035 *> ZDRVLS tests the least squares driver routines ZGELS, CGELSX, CGELSS,
00036 *> ZGELSY and CGELSD.
00037 *> \endverbatim
00038 *
00039 *  Arguments:
00040 *  ==========
00041 *
00042 *> \param[in] DOTYPE
00043 *> \verbatim
00044 *>          DOTYPE is LOGICAL array, dimension (NTYPES)
00045 *>          The matrix types to be used for testing.  Matrices of type j
00046 *>          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
00047 *>          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
00048 *>          The matrix of type j is generated as follows:
00049 *>          j=1: A = U*D*V where U and V are random unitary matrices
00050 *>               and D has random entries (> 0.1) taken from a uniform
00051 *>               distribution (0,1). A is full rank.
00052 *>          j=2: The same of 1, but A is scaled up.
00053 *>          j=3: The same of 1, but A is scaled down.
00054 *>          j=4: A = U*D*V where U and V are random unitary matrices
00055 *>               and D has 3*min(M,N)/4 random entries (> 0.1) taken
00056 *>               from a uniform distribution (0,1) and the remaining
00057 *>               entries set to 0. A is rank-deficient.
00058 *>          j=5: The same of 4, but A is scaled up.
00059 *>          j=6: The same of 5, but A is scaled down.
00060 *> \endverbatim
00061 *>
00062 *> \param[in] NM
00063 *> \verbatim
00064 *>          NM is INTEGER
00065 *>          The number of values of M contained in the vector MVAL.
00066 *> \endverbatim
00067 *>
00068 *> \param[in] MVAL
00069 *> \verbatim
00070 *>          MVAL is INTEGER array, dimension (NM)
00071 *>          The values of the matrix row dimension M.
00072 *> \endverbatim
00073 *>
00074 *> \param[in] NN
00075 *> \verbatim
00076 *>          NN is INTEGER
00077 *>          The number of values of N contained in the vector NVAL.
00078 *> \endverbatim
00079 *>
00080 *> \param[in] NVAL
00081 *> \verbatim
00082 *>          NVAL is INTEGER array, dimension (NN)
00083 *>          The values of the matrix column dimension N.
00084 *> \endverbatim
00085 *>
00086 *> \param[in] NNB
00087 *> \verbatim
00088 *>          NNB is INTEGER
00089 *>          The number of values of NB and NX contained in the
00090 *>          vectors NBVAL and NXVAL.  The blocking parameters are used
00091 *>          in pairs (NB,NX).
00092 *> \endverbatim
00093 *>
00094 *> \param[in] NBVAL
00095 *> \verbatim
00096 *>          NBVAL is INTEGER array, dimension (NNB)
00097 *>          The values of the blocksize NB.
00098 *> \endverbatim
00099 *>
00100 *> \param[in] NXVAL
00101 *> \verbatim
00102 *>          NXVAL is INTEGER array, dimension (NNB)
00103 *>          The values of the crossover point NX.
00104 *> \endverbatim
00105 *>
00106 *> \param[in] NNS
00107 *> \verbatim
00108 *>          NNS is INTEGER
00109 *>          The number of values of NRHS contained in the vector NSVAL.
00110 *> \endverbatim
00111 *>
00112 *> \param[in] NSVAL
00113 *> \verbatim
00114 *>          NSVAL is INTEGER array, dimension (NNS)
00115 *>          The values of the number of right hand sides NRHS.
00116 *> \endverbatim
00117 *>
00118 *> \param[in] THRESH
00119 *> \verbatim
00120 *>          THRESH is DOUBLE PRECISION
00121 *>          The threshold value for the test ratios.  A result is
00122 *>          included in the output file if RESULT >= THRESH.  To have
00123 *>          every test ratio printed, use THRESH = 0.
00124 *> \endverbatim
00125 *>
00126 *> \param[in] TSTERR
00127 *> \verbatim
00128 *>          TSTERR is LOGICAL
00129 *>          Flag that indicates whether error exits are to be tested.
00130 *> \endverbatim
00131 *>
00132 *> \param[out] A
00133 *> \verbatim
00134 *>          A is COMPLEX*16 array, dimension (MMAX*NMAX)
00135 *>          where MMAX is the maximum value of M in MVAL and NMAX is the
00136 *>          maximum value of N in NVAL.
00137 *> \endverbatim
00138 *>
00139 *> \param[out] COPYA
00140 *> \verbatim
00141 *>          COPYA is COMPLEX*16 array, dimension (MMAX*NMAX)
00142 *> \endverbatim
00143 *>
00144 *> \param[out] B
00145 *> \verbatim
00146 *>          B is COMPLEX*16 array, dimension (MMAX*NSMAX)
00147 *>          where MMAX is the maximum value of M in MVAL and NSMAX is the
00148 *>          maximum value of NRHS in NSVAL.
00149 *> \endverbatim
00150 *>
00151 *> \param[out] COPYB
00152 *> \verbatim
00153 *>          COPYB is COMPLEX*16 array, dimension (MMAX*NSMAX)
00154 *> \endverbatim
00155 *>
00156 *> \param[out] C
00157 *> \verbatim
00158 *>          C is COMPLEX*16 array, dimension (MMAX*NSMAX)
00159 *> \endverbatim
00160 *>
00161 *> \param[out] S
00162 *> \verbatim
00163 *>          S is DOUBLE PRECISION array, dimension
00164 *>                      (min(MMAX,NMAX))
00165 *> \endverbatim
00166 *>
00167 *> \param[out] COPYS
00168 *> \verbatim
00169 *>          COPYS is DOUBLE PRECISION array, dimension
00170 *>                      (min(MMAX,NMAX))
00171 *> \endverbatim
00172 *>
00173 *> \param[out] WORK
00174 *> \verbatim
00175 *>          WORK is COMPLEX*16 array, dimension
00176 *>                      (MMAX*NMAX + 4*NMAX + MMAX).
00177 *> \endverbatim
00178 *>
00179 *> \param[out] RWORK
00180 *> \verbatim
00181 *>          RWORK is DOUBLE PRECISION array, dimension (5*NMAX-1)
00182 *> \endverbatim
00183 *>
00184 *> \param[out] IWORK
00185 *> \verbatim
00186 *>          IWORK is INTEGER array, dimension (15*NMAX)
00187 *> \endverbatim
00188 *>
00189 *> \param[in] NOUT
00190 *> \verbatim
00191 *>          NOUT is INTEGER
00192 *>          The unit number for output.
00193 *> \endverbatim
00194 *
00195 *  Authors:
00196 *  ========
00197 *
00198 *> \author Univ. of Tennessee 
00199 *> \author Univ. of California Berkeley 
00200 *> \author Univ. of Colorado Denver 
00201 *> \author NAG Ltd. 
00202 *
00203 *> \date November 2011
00204 *
00205 *> \ingroup complex16_lin
00206 *
00207 *  =====================================================================
00208       SUBROUTINE ZDRVLS( DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB,
00209      $                   NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B,
00210      $                   COPYB, C, S, COPYS, WORK, RWORK, IWORK, NOUT )
00211 *
00212 *  -- LAPACK test routine (version 3.4.0) --
00213 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00214 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00215 *     November 2011
00216 *
00217 *     .. Scalar Arguments ..
00218       LOGICAL            TSTERR
00219       INTEGER            NM, NN, NNB, NNS, NOUT
00220       DOUBLE PRECISION   THRESH
00221 *     ..
00222 *     .. Array Arguments ..
00223       LOGICAL            DOTYPE( * )
00224       INTEGER            IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
00225      $                   NVAL( * ), NXVAL( * )
00226       DOUBLE PRECISION   COPYS( * ), RWORK( * ), S( * )
00227       COMPLEX*16         A( * ), B( * ), C( * ), COPYA( * ), COPYB( * ),
00228      $                   WORK( * )
00229 *     ..
00230 *
00231 *  =====================================================================
00232 *
00233 *     .. Parameters ..
00234       INTEGER            NTESTS
00235       PARAMETER          ( NTESTS = 18 )
00236       INTEGER            SMLSIZ
00237       PARAMETER          ( SMLSIZ = 25 )
00238       DOUBLE PRECISION   ONE, ZERO
00239       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
00240       COMPLEX*16         CONE, CZERO
00241       PARAMETER          ( CONE = ( 1.0D+0, 0.0D+0 ),
00242      $                   CZERO = ( 0.0D+0, 0.0D+0 ) )
00243 *     ..
00244 *     .. Local Scalars ..
00245       CHARACTER          TRANS
00246       CHARACTER*3        PATH
00247       INTEGER            CRANK, I, IM, IN, INB, INFO, INS, IRANK,
00248      $                   ISCALE, ITRAN, ITYPE, J, K, LDA, LDB, LDWORK,
00249      $                   LWLSY, LWORK, M, MNMIN, N, NB, NCOLS, NERRS,
00250      $                   NFAIL, NRHS, NROWS, NRUN, RANK
00251       DOUBLE PRECISION   EPS, NORMA, NORMB, RCOND
00252 *     ..
00253 *     .. Local Arrays ..
00254       INTEGER            ISEED( 4 ), ISEEDY( 4 )
00255       DOUBLE PRECISION   RESULT( NTESTS )
00256 *     ..
00257 *     .. External Functions ..
00258       DOUBLE PRECISION   DASUM, DLAMCH, ZQRT12, ZQRT14, ZQRT17
00259       EXTERNAL           DASUM, DLAMCH, ZQRT12, ZQRT14, ZQRT17
00260 *     ..
00261 *     .. External Subroutines ..
00262       EXTERNAL           ALAERH, ALAHD, ALASVM, DAXPY, DLASRT, XLAENV,
00263      $                   ZDSCAL, ZERRLS, ZGELS, ZGELSD, ZGELSS, ZGELSX,
00264      $                   ZGELSY, ZGEMM, ZLACPY, ZLARNV, ZQRT13, ZQRT15,
00265      $                   ZQRT16
00266 *     ..
00267 *     .. Intrinsic Functions ..
00268       INTRINSIC          DBLE, MAX, MIN, SQRT
00269 *     ..
00270 *     .. Scalars in Common ..
00271       LOGICAL            LERR, OK
00272       CHARACTER*32       SRNAMT
00273       INTEGER            INFOT, IOUNIT
00274 *     ..
00275 *     .. Common blocks ..
00276       COMMON             / INFOC / INFOT, IOUNIT, OK, LERR
00277       COMMON             / SRNAMC / SRNAMT
00278 *     ..
00279 *     .. Data statements ..
00280       DATA               ISEEDY / 1988, 1989, 1990, 1991 /
00281 *     ..
00282 *     .. Executable Statements ..
00283 *
00284 *     Initialize constants and the random number seed.
00285 *
00286       PATH( 1: 1 ) = 'Zomplex precision'
00287       PATH( 2: 3 ) = 'LS'
00288       NRUN = 0
00289       NFAIL = 0
00290       NERRS = 0
00291       DO 10 I = 1, 4
00292          ISEED( I ) = ISEEDY( I )
00293    10 CONTINUE
00294       EPS = DLAMCH( 'Epsilon' )
00295 *
00296 *     Threshold for rank estimation
00297 *
00298       RCOND = SQRT( EPS ) - ( SQRT( EPS )-EPS ) / 2
00299 *
00300 *     Test the error exits
00301 *
00302       CALL XLAENV( 9, SMLSIZ )
00303       IF( TSTERR )
00304      $   CALL ZERRLS( PATH, NOUT )
00305 *
00306 *     Print the header if NM = 0 or NN = 0 and THRESH = 0.
00307 *
00308       IF( ( NM.EQ.0 .OR. NN.EQ.0 ) .AND. THRESH.EQ.ZERO )
00309      $   CALL ALAHD( NOUT, PATH )
00310       INFOT = 0
00311 *
00312       DO 140 IM = 1, NM
00313          M = MVAL( IM )
00314          LDA = MAX( 1, M )
00315 *
00316          DO 130 IN = 1, NN
00317             N = NVAL( IN )
00318             MNMIN = MIN( M, N )
00319             LDB = MAX( 1, M, N )
00320 *
00321             DO 120 INS = 1, NNS
00322                NRHS = NSVAL( INS )
00323                LWORK = MAX( 1, ( M+NRHS )*( N+2 ), ( N+NRHS )*( M+2 ),
00324      $                 M*N+4*MNMIN+MAX( M, N ), 2*N+M )
00325 *
00326                DO 110 IRANK = 1, 2
00327                   DO 100 ISCALE = 1, 3
00328                      ITYPE = ( IRANK-1 )*3 + ISCALE
00329                      IF( .NOT.DOTYPE( ITYPE ) )
00330      $                  GO TO 100
00331 *
00332                      IF( IRANK.EQ.1 ) THEN
00333 *
00334 *                       Test ZGELS
00335 *
00336 *                       Generate a matrix of scaling type ISCALE
00337 *
00338                         CALL ZQRT13( ISCALE, M, N, COPYA, LDA, NORMA,
00339      $                               ISEED )
00340                         DO 40 INB = 1, NNB
00341                            NB = NBVAL( INB )
00342                            CALL XLAENV( 1, NB )
00343                            CALL XLAENV( 3, NXVAL( INB ) )
00344 *
00345                            DO 30 ITRAN = 1, 2
00346                               IF( ITRAN.EQ.1 ) THEN
00347                                  TRANS = 'N'
00348                                  NROWS = M
00349                                  NCOLS = N
00350                               ELSE
00351                                  TRANS = 'C'
00352                                  NROWS = N
00353                                  NCOLS = M
00354                               END IF
00355                               LDWORK = MAX( 1, NCOLS )
00356 *
00357 *                             Set up a consistent rhs
00358 *
00359                               IF( NCOLS.GT.0 ) THEN
00360                                  CALL ZLARNV( 2, ISEED, NCOLS*NRHS,
00361      $                                        WORK )
00362                                  CALL ZDSCAL( NCOLS*NRHS,
00363      $                                        ONE / DBLE( NCOLS ), WORK,
00364      $                                        1 )
00365                               END IF
00366                               CALL ZGEMM( TRANS, 'No transpose', NROWS,
00367      $                                    NRHS, NCOLS, CONE, COPYA, LDA,
00368      $                                    WORK, LDWORK, CZERO, B, LDB )
00369                               CALL ZLACPY( 'Full', NROWS, NRHS, B, LDB,
00370      $                                     COPYB, LDB )
00371 *
00372 *                             Solve LS or overdetermined system
00373 *
00374                               IF( M.GT.0 .AND. N.GT.0 ) THEN
00375                                  CALL ZLACPY( 'Full', M, N, COPYA, LDA,
00376      $                                        A, LDA )
00377                                  CALL ZLACPY( 'Full', NROWS, NRHS,
00378      $                                        COPYB, LDB, B, LDB )
00379                               END IF
00380                               SRNAMT = 'ZGELS '
00381                               CALL ZGELS( TRANS, M, N, NRHS, A, LDA, B,
00382      $                                    LDB, WORK, LWORK, INFO )
00383 *
00384                               IF( INFO.NE.0 )
00385      $                           CALL ALAERH( PATH, 'ZGELS ', INFO, 0,
00386      $                                        TRANS, M, N, NRHS, -1, NB,
00387      $                                        ITYPE, NFAIL, NERRS,
00388      $                                        NOUT )
00389 *
00390 *                             Check correctness of results
00391 *
00392                               LDWORK = MAX( 1, NROWS )
00393                               IF( NROWS.GT.0 .AND. NRHS.GT.0 )
00394      $                           CALL ZLACPY( 'Full', NROWS, NRHS,
00395      $                                        COPYB, LDB, C, LDB )
00396                               CALL ZQRT16( TRANS, M, N, NRHS, COPYA,
00397      $                                     LDA, B, LDB, C, LDB, RWORK,
00398      $                                     RESULT( 1 ) )
00399 *
00400                               IF( ( ITRAN.EQ.1 .AND. M.GE.N ) .OR.
00401      $                            ( ITRAN.EQ.2 .AND. M.LT.N ) ) THEN
00402 *
00403 *                                Solving LS system
00404 *
00405                                  RESULT( 2 ) = ZQRT17( TRANS, 1, M, N,
00406      $                                         NRHS, COPYA, LDA, B, LDB,
00407      $                                         COPYB, LDB, C, WORK,
00408      $                                         LWORK )
00409                               ELSE
00410 *
00411 *                                Solving overdetermined system
00412 *
00413                                  RESULT( 2 ) = ZQRT14( TRANS, M, N,
00414      $                                         NRHS, COPYA, LDA, B, LDB,
00415      $                                         WORK, LWORK )
00416                               END IF
00417 *
00418 *                             Print information about the tests that
00419 *                             did not pass the threshold.
00420 *
00421                               DO 20 K = 1, 2
00422                                  IF( RESULT( K ).GE.THRESH ) THEN
00423                                     IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
00424      $                                 CALL ALAHD( NOUT, PATH )
00425                                     WRITE( NOUT, FMT = 9999 )TRANS, M,
00426      $                                 N, NRHS, NB, ITYPE, K,
00427      $                                 RESULT( K )
00428                                     NFAIL = NFAIL + 1
00429                                  END IF
00430    20                         CONTINUE
00431                               NRUN = NRUN + 2
00432    30                      CONTINUE
00433    40                   CONTINUE
00434                      END IF
00435 *
00436 *                    Generate a matrix of scaling type ISCALE and rank
00437 *                    type IRANK.
00438 *
00439                      CALL ZQRT15( ISCALE, IRANK, M, N, NRHS, COPYA, LDA,
00440      $                            COPYB, LDB, COPYS, RANK, NORMA, NORMB,
00441      $                            ISEED, WORK, LWORK )
00442 *
00443 *                    workspace used: MAX(M+MIN(M,N),NRHS*MIN(M,N),2*N+M)
00444 *
00445                      DO 50 J = 1, N
00446                         IWORK( J ) = 0
00447    50                CONTINUE
00448                      LDWORK = MAX( 1, M )
00449 *
00450 *                    Test ZGELSX
00451 *
00452 *                    ZGELSX:  Compute the minimum-norm solution X
00453 *                    to min( norm( A * X - B ) )
00454 *                    using a complete orthogonal factorization.
00455 *
00456                      CALL ZLACPY( 'Full', M, N, COPYA, LDA, A, LDA )
00457                      CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, B, LDB )
00458 *
00459                      SRNAMT = 'ZGELSX'
00460                      CALL ZGELSX( M, N, NRHS, A, LDA, B, LDB, IWORK,
00461      $                            RCOND, CRANK, WORK, RWORK, INFO )
00462 *
00463                      IF( INFO.NE.0 )
00464      $                  CALL ALAERH( PATH, 'ZGELSX', INFO, 0, ' ', M, N,
00465      $                               NRHS, -1, NB, ITYPE, NFAIL, NERRS,
00466      $                               NOUT )
00467 *
00468 *                    workspace used: MAX( MNMIN+3*N, 2*MNMIN+NRHS )
00469 *
00470 *                    Test 3:  Compute relative error in svd
00471 *                             workspace: M*N + 4*MIN(M,N) + MAX(M,N)
00472 *
00473                      RESULT( 3 ) = ZQRT12( CRANK, CRANK, A, LDA, COPYS,
00474      $                             WORK, LWORK, RWORK )
00475 *
00476 *                    Test 4:  Compute error in solution
00477 *                             workspace:  M*NRHS + M
00478 *
00479                      CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, WORK,
00480      $                            LDWORK )
00481                      CALL ZQRT16( 'No transpose', M, N, NRHS, COPYA,
00482      $                            LDA, B, LDB, WORK, LDWORK, RWORK,
00483      $                            RESULT( 4 ) )
00484 *
00485 *                    Test 5:  Check norm of r'*A
00486 *                             workspace: NRHS*(M+N)
00487 *
00488                      RESULT( 5 ) = ZERO
00489                      IF( M.GT.CRANK )
00490      $                  RESULT( 5 ) = ZQRT17( 'No transpose', 1, M, N,
00491      $                                NRHS, COPYA, LDA, B, LDB, COPYB,
00492      $                                LDB, C, WORK, LWORK )
00493 *
00494 *                    Test 6:  Check if x is in the rowspace of A
00495 *                             workspace: (M+NRHS)*(N+2)
00496 *
00497                      RESULT( 6 ) = ZERO
00498 *
00499                      IF( N.GT.CRANK )
00500      $                  RESULT( 6 ) = ZQRT14( 'No transpose', M, N,
00501      $                                NRHS, COPYA, LDA, B, LDB, WORK,
00502      $                                LWORK )
00503 *
00504 *                    Print information about the tests that did not
00505 *                    pass the threshold.
00506 *
00507                      DO 60 K = 3, 6
00508                         IF( RESULT( K ).GE.THRESH ) THEN
00509                            IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
00510      $                        CALL ALAHD( NOUT, PATH )
00511                            WRITE( NOUT, FMT = 9998 )M, N, NRHS, 0,
00512      $                        ITYPE, K, RESULT( K )
00513                            NFAIL = NFAIL + 1
00514                         END IF
00515    60                CONTINUE
00516                      NRUN = NRUN + 4
00517 *
00518 *                    Loop for testing different block sizes.
00519 *
00520                      DO 90 INB = 1, NNB
00521                         NB = NBVAL( INB )
00522                         CALL XLAENV( 1, NB )
00523                         CALL XLAENV( 3, NXVAL( INB ) )
00524 *
00525 *                       Test ZGELSY
00526 *
00527 *                       ZGELSY:  Compute the minimum-norm solution
00528 *                       X to min( norm( A * X - B ) )
00529 *                       using the rank-revealing orthogonal
00530 *                       factorization.
00531 *
00532                         CALL ZLACPY( 'Full', M, N, COPYA, LDA, A, LDA )
00533                         CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, B,
00534      $                               LDB )
00535 *
00536 *                       Initialize vector IWORK.
00537 *
00538                         DO 70 J = 1, N
00539                            IWORK( J ) = 0
00540    70                   CONTINUE
00541 *
00542 *                       Set LWLSY to the adequate value.
00543 *
00544                         LWLSY = MNMIN + MAX( 2*MNMIN, NB*( N+1 ),
00545      $                          MNMIN+NB*NRHS )
00546                         LWLSY = MAX( 1, LWLSY )
00547 *
00548                         SRNAMT = 'ZGELSY'
00549                         CALL ZGELSY( M, N, NRHS, A, LDA, B, LDB, IWORK,
00550      $                               RCOND, CRANK, WORK, LWLSY, RWORK,
00551      $                               INFO )
00552                         IF( INFO.NE.0 )
00553      $                     CALL ALAERH( PATH, 'ZGELSY', INFO, 0, ' ', M,
00554      $                                  N, NRHS, -1, NB, ITYPE, NFAIL,
00555      $                                  NERRS, NOUT )
00556 *
00557 *                       workspace used: 2*MNMIN+NB*NB+NB*MAX(N,NRHS)
00558 *
00559 *                       Test 7:  Compute relative error in svd
00560 *                                workspace: M*N + 4*MIN(M,N) + MAX(M,N)
00561 *
00562                         RESULT( 7 ) = ZQRT12( CRANK, CRANK, A, LDA,
00563      $                                COPYS, WORK, LWORK, RWORK )
00564 *
00565 *                       Test 8:  Compute error in solution
00566 *                                workspace:  M*NRHS + M
00567 *
00568                         CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, WORK,
00569      $                               LDWORK )
00570                         CALL ZQRT16( 'No transpose', M, N, NRHS, COPYA,
00571      $                               LDA, B, LDB, WORK, LDWORK, RWORK,
00572      $                               RESULT( 8 ) )
00573 *
00574 *                       Test 9:  Check norm of r'*A
00575 *                                workspace: NRHS*(M+N)
00576 *
00577                         RESULT( 9 ) = ZERO
00578                         IF( M.GT.CRANK )
00579      $                     RESULT( 9 ) = ZQRT17( 'No transpose', 1, M,
00580      $                                   N, NRHS, COPYA, LDA, B, LDB,
00581      $                                   COPYB, LDB, C, WORK, LWORK )
00582 *
00583 *                       Test 10:  Check if x is in the rowspace of A
00584 *                                workspace: (M+NRHS)*(N+2)
00585 *
00586                         RESULT( 10 ) = ZERO
00587 *
00588                         IF( N.GT.CRANK )
00589      $                     RESULT( 10 ) = ZQRT14( 'No transpose', M, N,
00590      $                                    NRHS, COPYA, LDA, B, LDB,
00591      $                                    WORK, LWORK )
00592 *
00593 *                       Test ZGELSS
00594 *
00595 *                       ZGELSS:  Compute the minimum-norm solution
00596 *                       X to min( norm( A * X - B ) )
00597 *                       using the SVD.
00598 *
00599                         CALL ZLACPY( 'Full', M, N, COPYA, LDA, A, LDA )
00600                         CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, B,
00601      $                               LDB )
00602                         SRNAMT = 'ZGELSS'
00603                         CALL ZGELSS( M, N, NRHS, A, LDA, B, LDB, S,
00604      $                               RCOND, CRANK, WORK, LWORK, RWORK,
00605      $                               INFO )
00606 *
00607                         IF( INFO.NE.0 )
00608      $                     CALL ALAERH( PATH, 'ZGELSS', INFO, 0, ' ', M,
00609      $                                  N, NRHS, -1, NB, ITYPE, NFAIL,
00610      $                                  NERRS, NOUT )
00611 *
00612 *                       workspace used: 3*min(m,n) +
00613 *                                       max(2*min(m,n),nrhs,max(m,n))
00614 *
00615 *                       Test 11:  Compute relative error in svd
00616 *
00617                         IF( RANK.GT.0 ) THEN
00618                            CALL DAXPY( MNMIN, -ONE, COPYS, 1, S, 1 )
00619                            RESULT( 11 ) = DASUM( MNMIN, S, 1 ) /
00620      $                                    DASUM( MNMIN, COPYS, 1 ) /
00621      $                                    ( EPS*DBLE( MNMIN ) )
00622                         ELSE
00623                            RESULT( 11 ) = ZERO
00624                         END IF
00625 *
00626 *                       Test 12:  Compute error in solution
00627 *
00628                         CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, WORK,
00629      $                               LDWORK )
00630                         CALL ZQRT16( 'No transpose', M, N, NRHS, COPYA,
00631      $                               LDA, B, LDB, WORK, LDWORK, RWORK,
00632      $                               RESULT( 12 ) )
00633 *
00634 *                       Test 13:  Check norm of r'*A
00635 *
00636                         RESULT( 13 ) = ZERO
00637                         IF( M.GT.CRANK )
00638      $                     RESULT( 13 ) = ZQRT17( 'No transpose', 1, M,
00639      $                                    N, NRHS, COPYA, LDA, B, LDB,
00640      $                                    COPYB, LDB, C, WORK, LWORK )
00641 *
00642 *                       Test 14:  Check if x is in the rowspace of A
00643 *
00644                         RESULT( 14 ) = ZERO
00645                         IF( N.GT.CRANK )
00646      $                     RESULT( 14 ) = ZQRT14( 'No transpose', M, N,
00647      $                                    NRHS, COPYA, LDA, B, LDB,
00648      $                                    WORK, LWORK )
00649 *
00650 *                       Test ZGELSD
00651 *
00652 *                       ZGELSD:  Compute the minimum-norm solution X
00653 *                       to min( norm( A * X - B ) ) using a
00654 *                       divide and conquer SVD.
00655 *
00656                         CALL XLAENV( 9, 25 )
00657 *
00658                         CALL ZLACPY( 'Full', M, N, COPYA, LDA, A, LDA )
00659                         CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, B,
00660      $                               LDB )
00661 *
00662                         SRNAMT = 'ZGELSD'
00663                         CALL ZGELSD( M, N, NRHS, A, LDA, B, LDB, S,
00664      $                               RCOND, CRANK, WORK, LWORK, RWORK,
00665      $                               IWORK, INFO )
00666                         IF( INFO.NE.0 )
00667      $                     CALL ALAERH( PATH, 'ZGELSD', INFO, 0, ' ', M,
00668      $                                  N, NRHS, -1, NB, ITYPE, NFAIL,
00669      $                                  NERRS, NOUT )
00670 *
00671 *                       Test 15:  Compute relative error in svd
00672 *
00673                         IF( RANK.GT.0 ) THEN
00674                            CALL DAXPY( MNMIN, -ONE, COPYS, 1, S, 1 )
00675                            RESULT( 15 ) = DASUM( MNMIN, S, 1 ) /
00676      $                                    DASUM( MNMIN, COPYS, 1 ) /
00677      $                                    ( EPS*DBLE( MNMIN ) )
00678                         ELSE
00679                            RESULT( 15 ) = ZERO
00680                         END IF
00681 *
00682 *                       Test 16:  Compute error in solution
00683 *
00684                         CALL ZLACPY( 'Full', M, NRHS, COPYB, LDB, WORK,
00685      $                               LDWORK )
00686                         CALL ZQRT16( 'No transpose', M, N, NRHS, COPYA,
00687      $                               LDA, B, LDB, WORK, LDWORK, RWORK,
00688      $                               RESULT( 16 ) )
00689 *
00690 *                       Test 17:  Check norm of r'*A
00691 *
00692                         RESULT( 17 ) = ZERO
00693                         IF( M.GT.CRANK )
00694      $                     RESULT( 17 ) = ZQRT17( 'No transpose', 1, M,
00695      $                                    N, NRHS, COPYA, LDA, B, LDB,
00696      $                                    COPYB, LDB, C, WORK, LWORK )
00697 *
00698 *                       Test 18:  Check if x is in the rowspace of A
00699 *
00700                         RESULT( 18 ) = ZERO
00701                         IF( N.GT.CRANK )
00702      $                     RESULT( 18 ) = ZQRT14( 'No transpose', M, N,
00703      $                                    NRHS, COPYA, LDA, B, LDB,
00704      $                                    WORK, LWORK )
00705 *
00706 *                       Print information about the tests that did not
00707 *                       pass the threshold.
00708 *
00709                         DO 80 K = 7, NTESTS
00710                            IF( RESULT( K ).GE.THRESH ) THEN
00711                               IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
00712      $                           CALL ALAHD( NOUT, PATH )
00713                               WRITE( NOUT, FMT = 9998 )M, N, NRHS, NB,
00714      $                           ITYPE, K, RESULT( K )
00715                               NFAIL = NFAIL + 1
00716                            END IF
00717    80                   CONTINUE
00718                         NRUN = NRUN + 12
00719 *
00720    90                CONTINUE
00721   100             CONTINUE
00722   110          CONTINUE
00723   120       CONTINUE
00724   130    CONTINUE
00725   140 CONTINUE
00726 *
00727 *     Print a summary of the results.
00728 *
00729       CALL ALASVM( PATH, NOUT, NFAIL, NRUN, NERRS )
00730 *
00731  9999 FORMAT( ' TRANS=''', A1, ''', M=', I5, ', N=', I5, ', NRHS=', I4,
00732      $      ', NB=', I4, ', type', I2, ', test(', I2, ')=', G12.5 )
00733  9998 FORMAT( ' M=', I5, ', N=', I5, ', NRHS=', I4, ', NB=', I4,
00734      $      ', type', I2, ', test(', I2, ')=', G12.5 )
00735       RETURN
00736 *
00737 *     End of ZDRVLS
00738 *
00739       END
 All Files Functions