LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
dchkgt.f
Go to the documentation of this file.
00001 *> \brief \b DCHKGT
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 DCHKGT( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
00012 *                          A, AF, B, X, XACT, WORK, RWORK, IWORK, NOUT )
00013 * 
00014 *       .. Scalar Arguments ..
00015 *       LOGICAL            TSTERR
00016 *       INTEGER            NN, NNS, NOUT
00017 *       DOUBLE PRECISION   THRESH
00018 *       ..
00019 *       .. Array Arguments ..
00020 *       LOGICAL            DOTYPE( * )
00021 *       INTEGER            IWORK( * ), NSVAL( * ), NVAL( * )
00022 *       DOUBLE PRECISION   A( * ), AF( * ), B( * ), RWORK( * ), WORK( * ),
00023 *      $                   X( * ), XACT( * )
00024 *       ..
00025 *  
00026 *
00027 *> \par Purpose:
00028 *  =============
00029 *>
00030 *> \verbatim
00031 *>
00032 *> DCHKGT tests DGTTRF, -TRS, -RFS, and -CON
00033 *> \endverbatim
00034 *
00035 *  Arguments:
00036 *  ==========
00037 *
00038 *> \param[in] DOTYPE
00039 *> \verbatim
00040 *>          DOTYPE is LOGICAL array, dimension (NTYPES)
00041 *>          The matrix types to be used for testing.  Matrices of type j
00042 *>          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
00043 *>          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
00044 *> \endverbatim
00045 *>
00046 *> \param[in] NN
00047 *> \verbatim
00048 *>          NN is INTEGER
00049 *>          The number of values of N contained in the vector NVAL.
00050 *> \endverbatim
00051 *>
00052 *> \param[in] NVAL
00053 *> \verbatim
00054 *>          NVAL is INTEGER array, dimension (NN)
00055 *>          The values of the matrix dimension N.
00056 *> \endverbatim
00057 *>
00058 *> \param[in] NNS
00059 *> \verbatim
00060 *>          NNS is INTEGER
00061 *>          The number of values of NRHS contained in the vector NSVAL.
00062 *> \endverbatim
00063 *>
00064 *> \param[in] NSVAL
00065 *> \verbatim
00066 *>          NSVAL is INTEGER array, dimension (NNS)
00067 *>          The values of the number of right hand sides NRHS.
00068 *> \endverbatim
00069 *>
00070 *> \param[in] THRESH
00071 *> \verbatim
00072 *>          THRESH is DOUBLE PRECISION
00073 *>          The threshold value for the test ratios.  A result is
00074 *>          included in the output file if RESULT >= THRESH.  To have
00075 *>          every test ratio printed, use THRESH = 0.
00076 *> \endverbatim
00077 *>
00078 *> \param[in] TSTERR
00079 *> \verbatim
00080 *>          TSTERR is LOGICAL
00081 *>          Flag that indicates whether error exits are to be tested.
00082 *> \endverbatim
00083 *>
00084 *> \param[out] A
00085 *> \verbatim
00086 *>          A is DOUBLE PRECISION array, dimension (NMAX*4)
00087 *> \endverbatim
00088 *>
00089 *> \param[out] AF
00090 *> \verbatim
00091 *>          AF is DOUBLE PRECISION array, dimension (NMAX*4)
00092 *> \endverbatim
00093 *>
00094 *> \param[out] B
00095 *> \verbatim
00096 *>          B is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
00097 *>          where NSMAX is the largest entry in NSVAL.
00098 *> \endverbatim
00099 *>
00100 *> \param[out] X
00101 *> \verbatim
00102 *>          X is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
00103 *> \endverbatim
00104 *>
00105 *> \param[out] XACT
00106 *> \verbatim
00107 *>          XACT is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
00108 *> \endverbatim
00109 *>
00110 *> \param[out] WORK
00111 *> \verbatim
00112 *>          WORK is DOUBLE PRECISION array, dimension
00113 *>                      (NMAX*max(3,NSMAX))
00114 *> \endverbatim
00115 *>
00116 *> \param[out] RWORK
00117 *> \verbatim
00118 *>          RWORK is DOUBLE PRECISION array, dimension
00119 *>                      (max(NMAX,2*NSMAX))
00120 *> \endverbatim
00121 *>
00122 *> \param[out] IWORK
00123 *> \verbatim
00124 *>          IWORK is INTEGER array, dimension (2*NMAX)
00125 *> \endverbatim
00126 *>
00127 *> \param[in] NOUT
00128 *> \verbatim
00129 *>          NOUT is INTEGER
00130 *>          The unit number for output.
00131 *> \endverbatim
00132 *
00133 *  Authors:
00134 *  ========
00135 *
00136 *> \author Univ. of Tennessee 
00137 *> \author Univ. of California Berkeley 
00138 *> \author Univ. of Colorado Denver 
00139 *> \author NAG Ltd. 
00140 *
00141 *> \date November 2011
00142 *
00143 *> \ingroup double_lin
00144 *
00145 *  =====================================================================
00146       SUBROUTINE DCHKGT( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
00147      $                   A, AF, B, X, XACT, WORK, RWORK, IWORK, NOUT )
00148 *
00149 *  -- LAPACK test routine (version 3.4.0) --
00150 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00151 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00152 *     November 2011
00153 *
00154 *     .. Scalar Arguments ..
00155       LOGICAL            TSTERR
00156       INTEGER            NN, NNS, NOUT
00157       DOUBLE PRECISION   THRESH
00158 *     ..
00159 *     .. Array Arguments ..
00160       LOGICAL            DOTYPE( * )
00161       INTEGER            IWORK( * ), NSVAL( * ), NVAL( * )
00162       DOUBLE PRECISION   A( * ), AF( * ), B( * ), RWORK( * ), WORK( * ),
00163      $                   X( * ), XACT( * )
00164 *     ..
00165 *
00166 *  =====================================================================
00167 *
00168 *     .. Parameters ..
00169       DOUBLE PRECISION   ONE, ZERO
00170       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
00171       INTEGER            NTYPES
00172       PARAMETER          ( NTYPES = 12 )
00173       INTEGER            NTESTS
00174       PARAMETER          ( NTESTS = 7 )
00175 *     ..
00176 *     .. Local Scalars ..
00177       LOGICAL            TRFCON, ZEROT
00178       CHARACTER          DIST, NORM, TRANS, TYPE
00179       CHARACTER*3        PATH
00180       INTEGER            I, IMAT, IN, INFO, IRHS, ITRAN, IX, IZERO, J,
00181      $                   K, KL, KOFF, KU, LDA, M, MODE, N, NERRS, NFAIL,
00182      $                   NIMAT, NRHS, NRUN
00183       DOUBLE PRECISION   AINVNM, ANORM, COND, RCOND, RCONDC, RCONDI,
00184      $                   RCONDO
00185 *     ..
00186 *     .. Local Arrays ..
00187       CHARACTER          TRANSS( 3 )
00188       INTEGER            ISEED( 4 ), ISEEDY( 4 )
00189       DOUBLE PRECISION   RESULT( NTESTS ), Z( 3 )
00190 *     ..
00191 *     .. External Functions ..
00192       DOUBLE PRECISION   DASUM, DGET06, DLANGT
00193       EXTERNAL           DASUM, DGET06, DLANGT
00194 *     ..
00195 *     .. External Subroutines ..
00196       EXTERNAL           ALAERH, ALAHD, ALASUM, DCOPY, DERRGE, DGET04,
00197      $                   DGTCON, DGTRFS, DGTT01, DGTT02, DGTT05, DGTTRF,
00198      $                   DGTTRS, DLACPY, DLAGTM, DLARNV, DLATB4, DLATMS,
00199      $                   DSCAL
00200 *     ..
00201 *     .. Intrinsic Functions ..
00202       INTRINSIC          MAX
00203 *     ..
00204 *     .. Scalars in Common ..
00205       LOGICAL            LERR, OK
00206       CHARACTER*32       SRNAMT
00207       INTEGER            INFOT, NUNIT
00208 *     ..
00209 *     .. Common blocks ..
00210       COMMON             / INFOC / INFOT, NUNIT, OK, LERR
00211       COMMON             / SRNAMC / SRNAMT
00212 *     ..
00213 *     .. Data statements ..
00214       DATA               ISEEDY / 0, 0, 0, 1 / , TRANSS / 'N', 'T',
00215      $                   'C' /
00216 *     ..
00217 *     .. Executable Statements ..
00218 *
00219       PATH( 1: 1 ) = 'Double precision'
00220       PATH( 2: 3 ) = 'GT'
00221       NRUN = 0
00222       NFAIL = 0
00223       NERRS = 0
00224       DO 10 I = 1, 4
00225          ISEED( I ) = ISEEDY( I )
00226    10 CONTINUE
00227 *
00228 *     Test the error exits
00229 *
00230       IF( TSTERR )
00231      $   CALL DERRGE( PATH, NOUT )
00232       INFOT = 0
00233 *
00234       DO 110 IN = 1, NN
00235 *
00236 *        Do for each value of N in NVAL.
00237 *
00238          N = NVAL( IN )
00239          M = MAX( N-1, 0 )
00240          LDA = MAX( 1, N )
00241          NIMAT = NTYPES
00242          IF( N.LE.0 )
00243      $      NIMAT = 1
00244 *
00245          DO 100 IMAT = 1, NIMAT
00246 *
00247 *           Do the tests only if DOTYPE( IMAT ) is true.
00248 *
00249             IF( .NOT.DOTYPE( IMAT ) )
00250      $         GO TO 100
00251 *
00252 *           Set up parameters with DLATB4.
00253 *
00254             CALL DLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE,
00255      $                   COND, DIST )
00256 *
00257             ZEROT = IMAT.GE.8 .AND. IMAT.LE.10
00258             IF( IMAT.LE.6 ) THEN
00259 *
00260 *              Types 1-6:  generate matrices of known condition number.
00261 *
00262                KOFF = MAX( 2-KU, 3-MAX( 1, N ) )
00263                SRNAMT = 'DLATMS'
00264                CALL DLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE, COND,
00265      $                      ANORM, KL, KU, 'Z', AF( KOFF ), 3, WORK,
00266      $                      INFO )
00267 *
00268 *              Check the error code from DLATMS.
00269 *
00270                IF( INFO.NE.0 ) THEN
00271                   CALL ALAERH( PATH, 'DLATMS', INFO, 0, ' ', N, N, KL,
00272      $                         KU, -1, IMAT, NFAIL, NERRS, NOUT )
00273                   GO TO 100
00274                END IF
00275                IZERO = 0
00276 *
00277                IF( N.GT.1 ) THEN
00278                   CALL DCOPY( N-1, AF( 4 ), 3, A, 1 )
00279                   CALL DCOPY( N-1, AF( 3 ), 3, A( N+M+1 ), 1 )
00280                END IF
00281                CALL DCOPY( N, AF( 2 ), 3, A( M+1 ), 1 )
00282             ELSE
00283 *
00284 *              Types 7-12:  generate tridiagonal matrices with
00285 *              unknown condition numbers.
00286 *
00287                IF( .NOT.ZEROT .OR. .NOT.DOTYPE( 7 ) ) THEN
00288 *
00289 *                 Generate a matrix with elements from [-1,1].
00290 *
00291                   CALL DLARNV( 2, ISEED, N+2*M, A )
00292                   IF( ANORM.NE.ONE )
00293      $               CALL DSCAL( N+2*M, ANORM, A, 1 )
00294                ELSE IF( IZERO.GT.0 ) THEN
00295 *
00296 *                 Reuse the last matrix by copying back the zeroed out
00297 *                 elements.
00298 *
00299                   IF( IZERO.EQ.1 ) THEN
00300                      A( N ) = Z( 2 )
00301                      IF( N.GT.1 )
00302      $                  A( 1 ) = Z( 3 )
00303                   ELSE IF( IZERO.EQ.N ) THEN
00304                      A( 3*N-2 ) = Z( 1 )
00305                      A( 2*N-1 ) = Z( 2 )
00306                   ELSE
00307                      A( 2*N-2+IZERO ) = Z( 1 )
00308                      A( N-1+IZERO ) = Z( 2 )
00309                      A( IZERO ) = Z( 3 )
00310                   END IF
00311                END IF
00312 *
00313 *              If IMAT > 7, set one column of the matrix to 0.
00314 *
00315                IF( .NOT.ZEROT ) THEN
00316                   IZERO = 0
00317                ELSE IF( IMAT.EQ.8 ) THEN
00318                   IZERO = 1
00319                   Z( 2 ) = A( N )
00320                   A( N ) = ZERO
00321                   IF( N.GT.1 ) THEN
00322                      Z( 3 ) = A( 1 )
00323                      A( 1 ) = ZERO
00324                   END IF
00325                ELSE IF( IMAT.EQ.9 ) THEN
00326                   IZERO = N
00327                   Z( 1 ) = A( 3*N-2 )
00328                   Z( 2 ) = A( 2*N-1 )
00329                   A( 3*N-2 ) = ZERO
00330                   A( 2*N-1 ) = ZERO
00331                ELSE
00332                   IZERO = ( N+1 ) / 2
00333                   DO 20 I = IZERO, N - 1
00334                      A( 2*N-2+I ) = ZERO
00335                      A( N-1+I ) = ZERO
00336                      A( I ) = ZERO
00337    20             CONTINUE
00338                   A( 3*N-2 ) = ZERO
00339                   A( 2*N-1 ) = ZERO
00340                END IF
00341             END IF
00342 *
00343 *+    TEST 1
00344 *           Factor A as L*U and compute the ratio
00345 *              norm(L*U - A) / (n * norm(A) * EPS )
00346 *
00347             CALL DCOPY( N+2*M, A, 1, AF, 1 )
00348             SRNAMT = 'DGTTRF'
00349             CALL DGTTRF( N, AF, AF( M+1 ), AF( N+M+1 ), AF( N+2*M+1 ),
00350      $                   IWORK, INFO )
00351 *
00352 *           Check error code from DGTTRF.
00353 *
00354             IF( INFO.NE.IZERO )
00355      $         CALL ALAERH( PATH, 'DGTTRF', INFO, IZERO, ' ', N, N, 1,
00356      $                      1, -1, IMAT, NFAIL, NERRS, NOUT )
00357             TRFCON = INFO.NE.0
00358 *
00359             CALL DGTT01( N, A, A( M+1 ), A( N+M+1 ), AF, AF( M+1 ),
00360      $                   AF( N+M+1 ), AF( N+2*M+1 ), IWORK, WORK, LDA,
00361      $                   RWORK, RESULT( 1 ) )
00362 *
00363 *           Print the test ratio if it is .GE. THRESH.
00364 *
00365             IF( RESULT( 1 ).GE.THRESH ) THEN
00366                IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
00367      $            CALL ALAHD( NOUT, PATH )
00368                WRITE( NOUT, FMT = 9999 )N, IMAT, 1, RESULT( 1 )
00369                NFAIL = NFAIL + 1
00370             END IF
00371             NRUN = NRUN + 1
00372 *
00373             DO 50 ITRAN = 1, 2
00374                TRANS = TRANSS( ITRAN )
00375                IF( ITRAN.EQ.1 ) THEN
00376                   NORM = 'O'
00377                ELSE
00378                   NORM = 'I'
00379                END IF
00380                ANORM = DLANGT( NORM, N, A, A( M+1 ), A( N+M+1 ) )
00381 *
00382                IF( .NOT.TRFCON ) THEN
00383 *
00384 *                 Use DGTTRS to solve for one column at a time of inv(A)
00385 *                 or inv(A^T), computing the maximum column sum as we
00386 *                 go.
00387 *
00388                   AINVNM = ZERO
00389                   DO 40 I = 1, N
00390                      DO 30 J = 1, N
00391                         X( J ) = ZERO
00392    30                CONTINUE
00393                      X( I ) = ONE
00394                      CALL DGTTRS( TRANS, N, 1, AF, AF( M+1 ),
00395      $                            AF( N+M+1 ), AF( N+2*M+1 ), IWORK, X,
00396      $                            LDA, INFO )
00397                      AINVNM = MAX( AINVNM, DASUM( N, X, 1 ) )
00398    40             CONTINUE
00399 *
00400 *                 Compute RCONDC = 1 / (norm(A) * norm(inv(A))
00401 *
00402                   IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN
00403                      RCONDC = ONE
00404                   ELSE
00405                      RCONDC = ( ONE / ANORM ) / AINVNM
00406                   END IF
00407                   IF( ITRAN.EQ.1 ) THEN
00408                      RCONDO = RCONDC
00409                   ELSE
00410                      RCONDI = RCONDC
00411                   END IF
00412                ELSE
00413                   RCONDC = ZERO
00414                END IF
00415 *
00416 *+    TEST 7
00417 *              Estimate the reciprocal of the condition number of the
00418 *              matrix.
00419 *
00420                SRNAMT = 'DGTCON'
00421                CALL DGTCON( NORM, N, AF, AF( M+1 ), AF( N+M+1 ),
00422      $                      AF( N+2*M+1 ), IWORK, ANORM, RCOND, WORK,
00423      $                      IWORK( N+1 ), INFO )
00424 *
00425 *              Check error code from DGTCON.
00426 *
00427                IF( INFO.NE.0 )
00428      $            CALL ALAERH( PATH, 'DGTCON', INFO, 0, NORM, N, N, -1,
00429      $                         -1, -1, IMAT, NFAIL, NERRS, NOUT )
00430 *
00431                RESULT( 7 ) = DGET06( RCOND, RCONDC )
00432 *
00433 *              Print the test ratio if it is .GE. THRESH.
00434 *
00435                IF( RESULT( 7 ).GE.THRESH ) THEN
00436                   IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
00437      $               CALL ALAHD( NOUT, PATH )
00438                   WRITE( NOUT, FMT = 9997 )NORM, N, IMAT, 7,
00439      $               RESULT( 7 )
00440                   NFAIL = NFAIL + 1
00441                END IF
00442                NRUN = NRUN + 1
00443    50       CONTINUE
00444 *
00445 *           Skip the remaining tests if the matrix is singular.
00446 *
00447             IF( TRFCON )
00448      $         GO TO 100
00449 *
00450             DO 90 IRHS = 1, NNS
00451                NRHS = NSVAL( IRHS )
00452 *
00453 *              Generate NRHS random solution vectors.
00454 *
00455                IX = 1
00456                DO 60 J = 1, NRHS
00457                   CALL DLARNV( 2, ISEED, N, XACT( IX ) )
00458                   IX = IX + LDA
00459    60          CONTINUE
00460 *
00461                DO 80 ITRAN = 1, 3
00462                   TRANS = TRANSS( ITRAN )
00463                   IF( ITRAN.EQ.1 ) THEN
00464                      RCONDC = RCONDO
00465                   ELSE
00466                      RCONDC = RCONDI
00467                   END IF
00468 *
00469 *                 Set the right hand side.
00470 *
00471                   CALL DLAGTM( TRANS, N, NRHS, ONE, A, A( M+1 ),
00472      $                         A( N+M+1 ), XACT, LDA, ZERO, B, LDA )
00473 *
00474 *+    TEST 2
00475 *                 Solve op(A) * X = B and compute the residual.
00476 *
00477                   CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA )
00478                   SRNAMT = 'DGTTRS'
00479                   CALL DGTTRS( TRANS, N, NRHS, AF, AF( M+1 ),
00480      $                         AF( N+M+1 ), AF( N+2*M+1 ), IWORK, X,
00481      $                         LDA, INFO )
00482 *
00483 *                 Check error code from DGTTRS.
00484 *
00485                   IF( INFO.NE.0 )
00486      $               CALL ALAERH( PATH, 'DGTTRS', INFO, 0, TRANS, N, N,
00487      $                            -1, -1, NRHS, IMAT, NFAIL, NERRS,
00488      $                            NOUT )
00489 *
00490                   CALL DLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA )
00491                   CALL DGTT02( TRANS, N, NRHS, A, A( M+1 ), A( N+M+1 ),
00492      $                         X, LDA, WORK, LDA, RESULT( 2 ) )
00493 *
00494 *+    TEST 3
00495 *                 Check solution from generated exact solution.
00496 *
00497                   CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
00498      $                         RESULT( 3 ) )
00499 *
00500 *+    TESTS 4, 5, and 6
00501 *                 Use iterative refinement to improve the solution.
00502 *
00503                   SRNAMT = 'DGTRFS'
00504                   CALL DGTRFS( TRANS, N, NRHS, A, A( M+1 ), A( N+M+1 ),
00505      $                         AF, AF( M+1 ), AF( N+M+1 ),
00506      $                         AF( N+2*M+1 ), IWORK, B, LDA, X, LDA,
00507      $                         RWORK, RWORK( NRHS+1 ), WORK,
00508      $                         IWORK( N+1 ), INFO )
00509 *
00510 *                 Check error code from DGTRFS.
00511 *
00512                   IF( INFO.NE.0 )
00513      $               CALL ALAERH( PATH, 'DGTRFS', INFO, 0, TRANS, N, N,
00514      $                            -1, -1, NRHS, IMAT, NFAIL, NERRS,
00515      $                            NOUT )
00516 *
00517                   CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
00518      $                         RESULT( 4 ) )
00519                   CALL DGTT05( TRANS, N, NRHS, A, A( M+1 ), A( N+M+1 ),
00520      $                         B, LDA, X, LDA, XACT, LDA, RWORK,
00521      $                         RWORK( NRHS+1 ), RESULT( 5 ) )
00522 *
00523 *                 Print information about the tests that did not pass
00524 *                 the threshold.
00525 *
00526                   DO 70 K = 2, 6
00527                      IF( RESULT( K ).GE.THRESH ) THEN
00528                         IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
00529      $                     CALL ALAHD( NOUT, PATH )
00530                         WRITE( NOUT, FMT = 9998 )TRANS, N, NRHS, IMAT,
00531      $                     K, RESULT( K )
00532                         NFAIL = NFAIL + 1
00533                      END IF
00534    70             CONTINUE
00535                   NRUN = NRUN + 5
00536    80          CONTINUE
00537    90       CONTINUE
00538 *
00539   100    CONTINUE
00540   110 CONTINUE
00541 *
00542 *     Print a summary of the results.
00543 *
00544       CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS )
00545 *
00546  9999 FORMAT( 12X, 'N =', I5, ',', 10X, ' type ', I2, ', test(', I2,
00547      $      ') = ', G12.5 )
00548  9998 FORMAT( ' TRANS=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ',
00549      $      I2, ', test(', I2, ') = ', G12.5 )
00550  9997 FORMAT( ' NORM =''', A1, ''', N =', I5, ',', 10X, ' type ', I2,
00551      $      ', test(', I2, ') = ', G12.5 )
00552       RETURN
00553 *
00554 *     End of DCHKGT
00555 *
00556       END
 All Files Functions