LAPACK: TESTING/LIN/ddrvac.f Source File

Go to the documentation of this file.
00001       SUBROUTINE DDRVAC( DOTYPE, NM, MVAL, NNS, NSVAL, THRESH, NMAX,
00002      $                   A, AFAC, B, X, WORK,
00003      $                   RWORK, SWORK, NOUT )
00004 *
00005 *  -- LAPACK test routine (version 3.1.2) --
00006 *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
00007 *     April 2007
00008 *
00009 *     .. Scalar Arguments ..
00010       INTEGER            NMAX, NM, NNS, NOUT
00011       DOUBLE PRECISION   THRESH
00012 *     ..
00013 *     .. Array Arguments ..
00014       LOGICAL            DOTYPE( * )
00015       INTEGER            MVAL( * ), NSVAL( * )
00016       REAL               SWORK(*)
00017       DOUBLE PRECISION   A( * ), AFAC( * ), B( * ),
00018      $                   RWORK( * ), WORK( * ), X( * )
00019 *     ..
00020 *
00021 *  Purpose
00022 *  =======
00023 *
00024 *  DDRVAC tests DSPOSV.
00025 *
00026 *  Arguments
00027 *  =========
00028 *
00029 *  DOTYPE  (input) LOGICAL array, dimension (NTYPES)
00030 *          The matrix types to be used for testing.  Matrices of type j
00031 *          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
00032 *          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
00033 *
00034 *  NM      (input) INTEGER
00035 *          The number of values of N contained in the vector MVAL.
00036 *
00037 *  MVAL    (input) INTEGER array, dimension (NM)
00038 *          The values of the matrix dimension N.
00039 *
00040 *  NNS    (input) INTEGER
00041 *          The number of values of NRHS contained in the vector NSVAL.
00042 *
00043 *  NSVAL   (input) INTEGER array, dimension (NNS)
00044 *          The values of the number of right hand sides NRHS.
00045 *
00046 *  THRESH  (input) DOUBLE PRECISION
00047 *          The threshold value for the test ratios.  A result is
00048 *          included in the output file if RESULT >= THRESH.  To have
00049 *          every test ratio printed, use THRESH = 0.
00050 *
00051 *  NMAX    (input) INTEGER
00052 *          The maximum value permitted for N, used in dimensioning the
00053 *          work arrays.
00054 *
00055 *  A       (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
00056 *
00057 *  AFAC    (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
00058 *
00059 *  B       (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
00060 *
00061 *  X       (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
00062 *
00063 *  WORK    (workspace) DOUBLE PRECISION array, dimension
00064 *                      (NMAX*max(3,NSMAX))
00065 *
00066 *  RWORK   (workspace) DOUBLE PRECISION array, dimension
00067 *                      (max(2*NMAX,2*NSMAX+NWORK))
00068 *
00069 *  SWORK   (workspace) REAL array, dimension
00070 *                      (NMAX*(NSMAX+NMAX))
00071 *
00072 *  NOUT    (input) INTEGER
00073 *          The unit number for output.
00074 *
00075 *  =====================================================================
00076 *
00077 *     .. Parameters ..
00078       DOUBLE PRECISION   ONE, ZERO
00079       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
00080       INTEGER            NTYPES
00081       PARAMETER          ( NTYPES = 9 )
00082       INTEGER            NTESTS
00083       PARAMETER          ( NTESTS = 1 )
00084 *     ..
00085 *     .. Local Scalars ..
00086       LOGICAL            ZEROT
00087       CHARACTER          DIST, TYPE, UPLO, XTYPE
00088       CHARACTER*3        PATH
00089       INTEGER            I, IM, IMAT, INFO, IOFF, IRHS, IUPLO,
00090      $                   IZERO, KL, KU, LDA, MODE, N, 
00091      $                   NERRS, NFAIL, NIMAT, NRHS, NRUN
00092       DOUBLE PRECISION   ANORM, CNDNUM
00093 *     ..
00094 *     .. Local Arrays ..
00095       CHARACTER          UPLOS( 2 )
00096       INTEGER            ISEED( 4 ), ISEEDY( 4 )
00097       DOUBLE PRECISION   RESULT( NTESTS )
00098 *     ..
00099 *     .. Local Variables ..
00100       INTEGER            ITER, KASE
00101 *     ..
00102 *     .. External Functions ..
00103       LOGICAL            LSAME
00104       EXTERNAL           LSAME
00105 *     ..
00106 *     .. External Subroutines ..
00107       EXTERNAL           ALAERH, DLACPY,
00108      $                   DLARHS, DLASET, DLATB4, DLATMS, 
00109      $                   DPOT06, DSPOSV
00110 *     ..
00111 *     .. Intrinsic Functions ..
00112       INTRINSIC          DBLE, MAX, SQRT
00113 *     ..
00114 *     .. Scalars in Common ..
00115       LOGICAL            LERR, OK
00116       CHARACTER*32       SRNAMT
00117       INTEGER            INFOT, NUNIT
00118 *     ..
00119 *     .. Common blocks ..
00120       COMMON             / INFOC / INFOT, NUNIT, OK, LERR
00121       COMMON             / SRNAMC / SRNAMT
00122 *     ..
00123 *     .. Data statements ..
00124       DATA               ISEEDY / 1988, 1989, 1990, 1991 /
00125       DATA               UPLOS / 'U', 'L' /
00126 *     ..
00127 *     .. Executable Statements ..
00128 *
00129 *     Initialize constants and the random number seed.
00130 *
00131       KASE = 0
00132       PATH( 1: 1 ) = 'Double precision'
00133       PATH( 2: 3 ) = 'PO'
00134       NRUN = 0
00135       NFAIL = 0
00136       NERRS = 0
00137       DO 10 I = 1, 4
00138          ISEED( I ) = ISEEDY( I )
00139    10 CONTINUE
00140 *
00141       INFOT = 0
00142 *
00143 *     Do for each value of N in MVAL
00144 *
00145       DO 120 IM = 1, NM
00146          N = MVAL( IM )
00147          LDA = MAX( N, 1 )
00148          NIMAT = NTYPES
00149          IF( N.LE.0 )
00150      $      NIMAT = 1
00151 *
00152          DO 110 IMAT = 1, NIMAT
00153 *
00154 *           Do the tests only if DOTYPE( IMAT ) is true.
00155 *
00156             IF( .NOT.DOTYPE( IMAT ) )
00157      $         GO TO 110
00158 *
00159 *           Skip types 3, 4, or 5 if the matrix size is too small.
00160 *
00161             ZEROT = IMAT.GE.3 .AND. IMAT.LE.5
00162             IF( ZEROT .AND. N.LT.IMAT-2 )
00163      $         GO TO 110
00164 *
00165 *           Do first for UPLO = 'U', then for UPLO = 'L'
00166 *
00167             DO 100 IUPLO = 1, 2
00168                UPLO = UPLOS( IUPLO )
00169 *
00170 *              Set up parameters with DLATB4 and generate a test matrix
00171 *              with DLATMS.
00172 *
00173                CALL DLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE,
00174      $                      CNDNUM, DIST )
00175 *
00176                SRNAMT = 'DLATMS'
00177                CALL DLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE,
00178      $                      CNDNUM, ANORM, KL, KU, UPLO, A, LDA, WORK,
00179      $                      INFO )
00180 *
00181 *              Check error code from DLATMS.
00182 *
00183                IF( INFO.NE.0 ) THEN
00184                   CALL ALAERH( PATH, 'DLATMS', INFO, 0, UPLO, N, N, -1,
00185      $                         -1, -1, IMAT, NFAIL, NERRS, NOUT )
00186                   GO TO 100
00187                END IF
00188 *
00189 *              For types 3-5, zero one row and column of the matrix to
00190 *              test that INFO is returned correctly.
00191 *
00192                IF( ZEROT ) THEN
00193                   IF( IMAT.EQ.3 ) THEN
00194                      IZERO = 1
00195                   ELSE IF( IMAT.EQ.4 ) THEN
00196                      IZERO = N
00197                   ELSE
00198                      IZERO = N / 2 + 1
00199                   END IF
00200                   IOFF = ( IZERO-1 )*LDA
00201 *
00202 *                 Set row and column IZERO of A to 0.
00203 *
00204                   IF( IUPLO.EQ.1 ) THEN
00205                      DO 20 I = 1, IZERO - 1
00206                         A( IOFF+I ) = ZERO
00207    20                CONTINUE
00208                      IOFF = IOFF + IZERO
00209                      DO 30 I = IZERO, N
00210                         A( IOFF ) = ZERO
00211                         IOFF = IOFF + LDA
00212    30                CONTINUE
00213                   ELSE
00214                      IOFF = IZERO
00215                      DO 40 I = 1, IZERO - 1
00216                         A( IOFF ) = ZERO
00217                         IOFF = IOFF + LDA
00218    40                CONTINUE
00219                      IOFF = IOFF - IZERO
00220                      DO 50 I = IZERO, N
00221                         A( IOFF+I ) = ZERO
00222    50                CONTINUE
00223                   END IF
00224                ELSE
00225                   IZERO = 0
00226                END IF
00227 *
00228                DO 60 IRHS = 1, NNS
00229                   NRHS = NSVAL( IRHS )
00230                   XTYPE = 'N'
00231 *
00232 *                 Form an exact solution and set the right hand side.
00233 *
00234                   SRNAMT = 'DLARHS'
00235                   CALL DLARHS( PATH, XTYPE, UPLO, ' ', N, N, KL, KU,
00236      $                         NRHS, A, LDA, X, LDA, B, LDA,
00237      $                         ISEED, INFO )
00238 *
00239 *                 Compute the L*L' or U'*U factorization of the
00240 *                 matrix and solve the system.
00241 *
00242                   SRNAMT = 'DSPOSV '
00243                   KASE = KASE + 1
00244 *
00245                   CALL DLACPY( 'All', N, N, A, LDA, AFAC, LDA) 
00246 *
00247                   CALL DSPOSV( UPLO, N, NRHS, AFAC, LDA, B, LDA, X, LDA,
00248      $                         WORK, SWORK, ITER, INFO )
00249 
00250                   IF (ITER.LT.0) THEN
00251                      CALL DLACPY( 'All', N, N, A, LDA, AFAC, LDA )
00252                   ENDIF
00253 *
00254 *                 Check error code from DSPOSV .
00255 *
00256                   IF( INFO.NE.IZERO ) THEN
00257 *
00258                      IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
00259      $                  CALL ALAHD( NOUT, PATH )
00260                      NERRS = NERRS + 1
00261 *
00262                      IF( INFO.NE.IZERO .AND. IZERO.NE.0 ) THEN
00263                         WRITE( NOUT, FMT = 9988 )'DSPOSV',INFO,IZERO,N,
00264      $                     IMAT
00265                      ELSE
00266                         WRITE( NOUT, FMT = 9975 )'DSPOSV',INFO,N,IMAT
00267                      END IF
00268                   END IF
00269 *
00270 *                 Skip the remaining test if the matrix is singular.
00271 *
00272                   IF( INFO.NE.0 )
00273      $               GO TO 110
00274 *
00275 *                 Check the quality of the solution
00276 *
00277                   CALL DLACPY( 'All', N, NRHS, B, LDA, WORK, LDA )
00278 *
00279                   CALL DPOT06( UPLO, N, NRHS, A, LDA, X, LDA, WORK,
00280      $               LDA, RWORK, RESULT( 1 ) )
00281 *
00282 *                 Check if the test passes the tesing.
00283 *                 Print information about the tests that did not
00284 *                 pass the testing.
00285 *
00286 *                 If iterative refinement has been used and claimed to 
00287 *                 be successful (ITER>0), we want
00288 *                 NORM1(B - A*X)/(NORM1(A)*NORM1(X)*EPS*SRQT(N)) < 1
00289 *
00290 *                 If double precision has been used (ITER<0), we want
00291 *                 NORM1(B - A*X)/(NORM1(A)*NORM1(X)*EPS) < THRES
00292 *                 (Cf. the linear solver testing routines)
00293 *
00294                   IF ((THRESH.LE.0.0E+00)
00295      $               .OR.((ITER.GE.0).AND.(N.GT.0)
00296      $               .AND.(RESULT(1).GE.SQRT(DBLE(N))))
00297      $               .OR.((ITER.LT.0).AND.(RESULT(1).GE.THRESH))) THEN
00298 *
00299                      IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) THEN
00300                         WRITE( NOUT, FMT = 8999 )'DPO'
00301                         WRITE( NOUT, FMT = '( '' Matrix types:'' )' )
00302                         WRITE( NOUT, FMT = 8979 )
00303                         WRITE( NOUT, FMT = '( '' Test ratios:'' )' )
00304                         WRITE( NOUT, FMT = 8960 )1
00305                         WRITE( NOUT, FMT = '( '' Messages:'' )' )
00306                      END IF
00307 *
00308                      WRITE( NOUT, FMT = 9998 )UPLO, N, NRHS, IMAT, 1,
00309      $                  RESULT( 1 )
00310 *
00311                      NFAIL = NFAIL + 1
00312 *
00313                   END IF
00314 *
00315                   NRUN = NRUN + 1
00316 *
00317    60          CONTINUE
00318   100       CONTINUE
00319   110    CONTINUE
00320   120 CONTINUE
00321 *
00322   130 CONTINUE
00323 *
00324 *     Print a summary of the results.
00325 *
00326       IF( NFAIL.GT.0 ) THEN
00327          WRITE( NOUT, FMT = 9996 )'DSPOSV', NFAIL, NRUN
00328       ELSE
00329          WRITE( NOUT, FMT = 9995 )'DSPOSV', NRUN
00330       END IF
00331       IF( NERRS.GT.0 ) THEN
00332          WRITE( NOUT, FMT = 9994 )NERRS
00333       END IF
00334 *
00335  9998 FORMAT( ' UPLO=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ',
00336      $      I2, ', test(', I2, ') =', G12.5 )
00337  9996 FORMAT( 1X, A6, ': ', I6, ' out of ', I6,
00338      $      ' tests failed to pass the threshold' )
00339  9995 FORMAT( /1X, 'All tests for ', A6,
00340      $      ' routines passed the threshold (', I6, ' tests run)' )
00341  9994 FORMAT( 6X, I6, ' error messages recorded' )
00342 *
00343 *     SUBNAM, INFO, INFOE, N, IMAT
00344 *
00345  9988 FORMAT( ' *** ', A6, ' returned with INFO =', I5, ' instead of ',
00346      $      I5, / ' ==> N =', I5, ', type ',
00347      $      I2 )
00348 *
00349 *     SUBNAM, INFO, N, IMAT
00350 *
00351  9975 FORMAT( ' *** Error code from ', A6, '=', I5, ' for M=', I5,
00352      $      ', type ', I2 )
00353  8999 FORMAT( / 1X, A3, ':  positive definite dense matrices' )
00354  8979 FORMAT( 4X, '1. Diagonal', 24X, '7. Last n/2 columns zero', / 4X,
00355      $      '2. Upper triangular', 16X,
00356      $      '8. Random, CNDNUM = sqrt(0.1/EPS)', / 4X,
00357      $      '3. Lower triangular', 16X, '9. Random, CNDNUM = 0.1/EPS',
00358      $      / 4X, '4. Random, CNDNUM = 2', 13X,
00359      $      '10. Scaled near underflow', / 4X, '5. First column zero',
00360      $      14X, '11. Scaled near overflow', / 4X,
00361      $      '6. Last column zero' )
00362  8960 FORMAT( 3X, I2, ': norm_1( B - A * X )  / ',
00363      $      '( norm_1(A) * norm_1(X) * EPS * SQRT(N) ) > 1 if ITERREF',
00364      $      / 4x, 'or norm_1( B - A * X )  / ',
00365      $      '( norm_1(A) * norm_1(X) * EPS ) > THRES if DPOTRF' )
00366       
00367       RETURN
00368 *
00369 *     End of DDRVAC
00370 *
00371       END