LAPACK  3.4.1
LAPACK: Linear Algebra PACKage
dchkbb.f
Go to the documentation of this file.
00001 *> \brief \b DCHKBB
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 DCHKBB( NSIZES, MVAL, NVAL, NWDTHS, KK, NTYPES, DOTYPE,
00012 *                          NRHS, ISEED, THRESH, NOUNIT, A, LDA, AB, LDAB,
00013 *                          BD, BE, Q, LDQ, P, LDP, C, LDC, CC, WORK,
00014 *                          LWORK, RESULT, INFO )
00015 * 
00016 *       .. Scalar Arguments ..
00017 *       INTEGER            INFO, LDA, LDAB, LDC, LDP, LDQ, LWORK, NOUNIT,
00018 *      $                   NRHS, NSIZES, NTYPES, NWDTHS
00019 *       DOUBLE PRECISION   THRESH
00020 *       ..
00021 *       .. Array Arguments ..
00022 *       LOGICAL            DOTYPE( * )
00023 *       INTEGER            ISEED( 4 ), KK( * ), MVAL( * ), NVAL( * )
00024 *       DOUBLE PRECISION   A( LDA, * ), AB( LDAB, * ), BD( * ), BE( * ),
00025 *      $                   C( LDC, * ), CC( LDC, * ), P( LDP, * ),
00026 *      $                   Q( LDQ, * ), RESULT( * ), WORK( * )
00027 *       ..
00028 *  
00029 *
00030 *> \par Purpose:
00031 *  =============
00032 *>
00033 *> \verbatim
00034 *>
00035 *> DCHKBB tests the reduction of a general real rectangular band
00036 *> matrix to bidiagonal form.
00037 *>
00038 *> DGBBRD factors a general band matrix A as  Q B P* , where * means
00039 *> transpose, B is upper bidiagonal, and Q and P are orthogonal;
00040 *> DGBBRD can also overwrite a given matrix C with Q* C .
00041 *>
00042 *> For each pair of matrix dimensions (M,N) and each selected matrix
00043 *> type, an M by N matrix A and an M by NRHS matrix C are generated.
00044 *> The problem dimensions are as follows
00045 *>    A:          M x N
00046 *>    Q:          M x M
00047 *>    P:          N x N
00048 *>    B:          min(M,N) x min(M,N)
00049 *>    C:          M x NRHS
00050 *>
00051 *> For each generated matrix, 4 tests are performed:
00052 *>
00053 *> (1)   | A - Q B PT | / ( |A| max(M,N) ulp ), PT = P'
00054 *>
00055 *> (2)   | I - Q' Q | / ( M ulp )
00056 *>
00057 *> (3)   | I - PT PT' | / ( N ulp )
00058 *>
00059 *> (4)   | Y - Q' C | / ( |Y| max(M,NRHS) ulp ), where Y = Q' C.
00060 *>
00061 *> The "types" are specified by a logical array DOTYPE( 1:NTYPES );
00062 *> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
00063 *> Currently, the list of possible types is:
00064 *>
00065 *> The possible matrix types are
00066 *>
00067 *> (1)  The zero matrix.
00068 *> (2)  The identity matrix.
00069 *>
00070 *> (3)  A diagonal matrix with evenly spaced entries
00071 *>      1, ..., ULP  and random signs.
00072 *>      (ULP = (first number larger than 1) - 1 )
00073 *> (4)  A diagonal matrix with geometrically spaced entries
00074 *>      1, ..., ULP  and random signs.
00075 *> (5)  A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
00076 *>      and random signs.
00077 *>
00078 *> (6)  Same as (3), but multiplied by SQRT( overflow threshold )
00079 *> (7)  Same as (3), but multiplied by SQRT( underflow threshold )
00080 *>
00081 *> (8)  A matrix of the form  U D V, where U and V are orthogonal and
00082 *>      D has evenly spaced entries 1, ..., ULP with random signs
00083 *>      on the diagonal.
00084 *>
00085 *> (9)  A matrix of the form  U D V, where U and V are orthogonal and
00086 *>      D has geometrically spaced entries 1, ..., ULP with random
00087 *>      signs on the diagonal.
00088 *>
00089 *> (10) A matrix of the form  U D V, where U and V are orthogonal and
00090 *>      D has "clustered" entries 1, ULP,..., ULP with random
00091 *>      signs on the diagonal.
00092 *>
00093 *> (11) Same as (8), but multiplied by SQRT( overflow threshold )
00094 *> (12) Same as (8), but multiplied by SQRT( underflow threshold )
00095 *>
00096 *> (13) Rectangular matrix with random entries chosen from (-1,1).
00097 *> (14) Same as (13), but multiplied by SQRT( overflow threshold )
00098 *> (15) Same as (13), but multiplied by SQRT( underflow threshold )
00099 *> \endverbatim
00100 *
00101 *  Arguments:
00102 *  ==========
00103 *
00104 *> \param[in] NSIZES
00105 *> \verbatim
00106 *>          NSIZES is INTEGER
00107 *>          The number of values of M and N contained in the vectors
00108 *>          MVAL and NVAL.  The matrix sizes are used in pairs (M,N).
00109 *>          If NSIZES is zero, DCHKBB does nothing.  NSIZES must be at
00110 *>          least zero.
00111 *> \endverbatim
00112 *>
00113 *> \param[in] MVAL
00114 *> \verbatim
00115 *>          MVAL is INTEGER array, dimension (NSIZES)
00116 *>          The values of the matrix row dimension M.
00117 *> \endverbatim
00118 *>
00119 *> \param[in] NVAL
00120 *> \verbatim
00121 *>          NVAL is INTEGER array, dimension (NSIZES)
00122 *>          The values of the matrix column dimension N.
00123 *> \endverbatim
00124 *>
00125 *> \param[in] NWDTHS
00126 *> \verbatim
00127 *>          NWDTHS is INTEGER
00128 *>          The number of bandwidths to use.  If it is zero,
00129 *>          DCHKBB does nothing.  It must be at least zero.
00130 *> \endverbatim
00131 *>
00132 *> \param[in] KK
00133 *> \verbatim
00134 *>          KK is INTEGER array, dimension (NWDTHS)
00135 *>          An array containing the bandwidths to be used for the band
00136 *>          matrices.  The values must be at least zero.
00137 *> \endverbatim
00138 *>
00139 *> \param[in] NTYPES
00140 *> \verbatim
00141 *>          NTYPES is INTEGER
00142 *>          The number of elements in DOTYPE.   If it is zero, DCHKBB
00143 *>          does nothing.  It must be at least zero.  If it is MAXTYP+1
00144 *>          and NSIZES is 1, then an additional type, MAXTYP+1 is
00145 *>          defined, which is to use whatever matrix is in A.  This
00146 *>          is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
00147 *>          DOTYPE(MAXTYP+1) is .TRUE. .
00148 *> \endverbatim
00149 *>
00150 *> \param[in] DOTYPE
00151 *> \verbatim
00152 *>          DOTYPE is LOGICAL array, dimension (NTYPES)
00153 *>          If DOTYPE(j) is .TRUE., then for each size in NN a
00154 *>          matrix of that size and of type j will be generated.
00155 *>          If NTYPES is smaller than the maximum number of types
00156 *>          defined (PARAMETER MAXTYP), then types NTYPES+1 through
00157 *>          MAXTYP will not be generated.  If NTYPES is larger
00158 *>          than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
00159 *>          will be ignored.
00160 *> \endverbatim
00161 *>
00162 *> \param[in] NRHS
00163 *> \verbatim
00164 *>          NRHS is INTEGER
00165 *>          The number of columns in the "right-hand side" matrix C.
00166 *>          If NRHS = 0, then the operations on the right-hand side will
00167 *>          not be tested. NRHS must be at least 0.
00168 *> \endverbatim
00169 *>
00170 *> \param[in,out] ISEED
00171 *> \verbatim
00172 *>          ISEED is INTEGER array, dimension (4)
00173 *>          On entry ISEED specifies the seed of the random number
00174 *>          generator. The array elements should be between 0 and 4095;
00175 *>          if not they will be reduced mod 4096.  Also, ISEED(4) must
00176 *>          be odd.  The random number generator uses a linear
00177 *>          congruential sequence limited to small integers, and so
00178 *>          should produce machine independent random numbers. The
00179 *>          values of ISEED are changed on exit, and can be used in the
00180 *>          next call to DCHKBB to continue the same random number
00181 *>          sequence.
00182 *> \endverbatim
00183 *>
00184 *> \param[in] THRESH
00185 *> \verbatim
00186 *>          THRESH is DOUBLE PRECISION
00187 *>          A test will count as "failed" if the "error", computed as
00188 *>          described above, exceeds THRESH.  Note that the error
00189 *>          is scaled to be O(1), so THRESH should be a reasonably
00190 *>          small multiple of 1, e.g., 10 or 100.  In particular,
00191 *>          it should not depend on the precision (single vs. double)
00192 *>          or the size of the matrix.  It must be at least zero.
00193 *> \endverbatim
00194 *>
00195 *> \param[in] NOUNIT
00196 *> \verbatim
00197 *>          NOUNIT is INTEGER
00198 *>          The FORTRAN unit number for printing out error messages
00199 *>          (e.g., if a routine returns IINFO not equal to 0.)
00200 *> \endverbatim
00201 *>
00202 *> \param[in,out] A
00203 *> \verbatim
00204 *>          A is DOUBLE PRECISION array, dimension
00205 *>                            (LDA, max(NN))
00206 *>          Used to hold the matrix A.
00207 *> \endverbatim
00208 *>
00209 *> \param[in] LDA
00210 *> \verbatim
00211 *>          LDA is INTEGER
00212 *>          The leading dimension of A.  It must be at least 1
00213 *>          and at least max( NN ).
00214 *> \endverbatim
00215 *>
00216 *> \param[out] AB
00217 *> \verbatim
00218 *>          AB is DOUBLE PRECISION array, dimension (LDAB, max(NN))
00219 *>          Used to hold A in band storage format.
00220 *> \endverbatim
00221 *>
00222 *> \param[in] LDAB
00223 *> \verbatim
00224 *>          LDAB is INTEGER
00225 *>          The leading dimension of AB.  It must be at least 2 (not 1!)
00226 *>          and at least max( KK )+1.
00227 *> \endverbatim
00228 *>
00229 *> \param[out] BD
00230 *> \verbatim
00231 *>          BD is DOUBLE PRECISION array, dimension (max(NN))
00232 *>          Used to hold the diagonal of the bidiagonal matrix computed
00233 *>          by DGBBRD.
00234 *> \endverbatim
00235 *>
00236 *> \param[out] BE
00237 *> \verbatim
00238 *>          BE is DOUBLE PRECISION array, dimension (max(NN))
00239 *>          Used to hold the off-diagonal of the bidiagonal matrix
00240 *>          computed by DGBBRD.
00241 *> \endverbatim
00242 *>
00243 *> \param[out] Q
00244 *> \verbatim
00245 *>          Q is DOUBLE PRECISION array, dimension (LDQ, max(NN))
00246 *>          Used to hold the orthogonal matrix Q computed by DGBBRD.
00247 *> \endverbatim
00248 *>
00249 *> \param[in] LDQ
00250 *> \verbatim
00251 *>          LDQ is INTEGER
00252 *>          The leading dimension of Q.  It must be at least 1
00253 *>          and at least max( NN ).
00254 *> \endverbatim
00255 *>
00256 *> \param[out] P
00257 *> \verbatim
00258 *>          P is DOUBLE PRECISION array, dimension (LDP, max(NN))
00259 *>          Used to hold the orthogonal matrix P computed by DGBBRD.
00260 *> \endverbatim
00261 *>
00262 *> \param[in] LDP
00263 *> \verbatim
00264 *>          LDP is INTEGER
00265 *>          The leading dimension of P.  It must be at least 1
00266 *>          and at least max( NN ).
00267 *> \endverbatim
00268 *>
00269 *> \param[out] C
00270 *> \verbatim
00271 *>          C is DOUBLE PRECISION array, dimension (LDC, max(NN))
00272 *>          Used to hold the matrix C updated by DGBBRD.
00273 *> \endverbatim
00274 *>
00275 *> \param[in] LDC
00276 *> \verbatim
00277 *>          LDC is INTEGER
00278 *>          The leading dimension of U.  It must be at least 1
00279 *>          and at least max( NN ).
00280 *> \endverbatim
00281 *>
00282 *> \param[out] CC
00283 *> \verbatim
00284 *>          CC is DOUBLE PRECISION array, dimension (LDC, max(NN))
00285 *>          Used to hold a copy of the matrix C.
00286 *> \endverbatim
00287 *>
00288 *> \param[out] WORK
00289 *> \verbatim
00290 *>          WORK is DOUBLE PRECISION array, dimension (LWORK)
00291 *> \endverbatim
00292 *>
00293 *> \param[in] LWORK
00294 *> \verbatim
00295 *>          LWORK is INTEGER
00296 *>          The number of entries in WORK.  This must be at least
00297 *>          max( LDA+1, max(NN)+1 )*max(NN).
00298 *> \endverbatim
00299 *>
00300 *> \param[out] RESULT
00301 *> \verbatim
00302 *>          RESULT is DOUBLE PRECISION array, dimension (4)
00303 *>          The values computed by the tests described above.
00304 *>          The values are currently limited to 1/ulp, to avoid
00305 *>          overflow.
00306 *> \endverbatim
00307 *>
00308 *> \param[out] INFO
00309 *> \verbatim
00310 *>          INFO is INTEGER
00311 *>          If 0, then everything ran OK.
00312 *>
00313 *>-----------------------------------------------------------------------
00314 *>
00315 *>       Some Local Variables and Parameters:
00316 *>       ---- ----- --------- --- ----------
00317 *>       ZERO, ONE       Real 0 and 1.
00318 *>       MAXTYP          The number of types defined.
00319 *>       NTEST           The number of tests performed, or which can
00320 *>                       be performed so far, for the current matrix.
00321 *>       NTESTT          The total number of tests performed so far.
00322 *>       NMAX            Largest value in NN.
00323 *>       NMATS           The number of matrices generated so far.
00324 *>       NERRS           The number of tests which have exceeded THRESH
00325 *>                       so far.
00326 *>       COND, IMODE     Values to be passed to the matrix generators.
00327 *>       ANORM           Norm of A; passed to matrix generators.
00328 *>
00329 *>       OVFL, UNFL      Overflow and underflow thresholds.
00330 *>       ULP, ULPINV     Finest relative precision and its inverse.
00331 *>       RTOVFL, RTUNFL  Square roots of the previous 2 values.
00332 *>               The following four arrays decode JTYPE:
00333 *>       KTYPE(j)        The general type (1-10) for type "j".
00334 *>       KMODE(j)        The MODE value to be passed to the matrix
00335 *>                       generator for type "j".
00336 *>       KMAGN(j)        The order of magnitude ( O(1),
00337 *>                       O(overflow^(1/2) ), O(underflow^(1/2) )
00338 *> \endverbatim
00339 *
00340 *  Authors:
00341 *  ========
00342 *
00343 *> \author Univ. of Tennessee 
00344 *> \author Univ. of California Berkeley 
00345 *> \author Univ. of Colorado Denver 
00346 *> \author NAG Ltd. 
00347 *
00348 *> \date November 2011
00349 *
00350 *> \ingroup double_eig
00351 *
00352 *  =====================================================================
00353       SUBROUTINE DCHKBB( NSIZES, MVAL, NVAL, NWDTHS, KK, NTYPES, DOTYPE,
00354      $                   NRHS, ISEED, THRESH, NOUNIT, A, LDA, AB, LDAB,
00355      $                   BD, BE, Q, LDQ, P, LDP, C, LDC, CC, WORK,
00356      $                   LWORK, RESULT, INFO )
00357 *
00358 *  -- LAPACK test routine (input) --
00359 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00360 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00361 *     November 2011
00362 *
00363 *     .. Scalar Arguments ..
00364       INTEGER            INFO, LDA, LDAB, LDC, LDP, LDQ, LWORK, NOUNIT,
00365      $                   NRHS, NSIZES, NTYPES, NWDTHS
00366       DOUBLE PRECISION   THRESH
00367 *     ..
00368 *     .. Array Arguments ..
00369       LOGICAL            DOTYPE( * )
00370       INTEGER            ISEED( 4 ), KK( * ), MVAL( * ), NVAL( * )
00371       DOUBLE PRECISION   A( LDA, * ), AB( LDAB, * ), BD( * ), BE( * ),
00372      $                   C( LDC, * ), CC( LDC, * ), P( LDP, * ),
00373      $                   Q( LDQ, * ), RESULT( * ), WORK( * )
00374 *     ..
00375 *
00376 *  =====================================================================
00377 *
00378 *     .. Parameters ..
00379       DOUBLE PRECISION   ZERO, ONE
00380       PARAMETER          ( ZERO = 0.0D0, ONE = 1.0D0 )
00381       INTEGER            MAXTYP
00382       PARAMETER          ( MAXTYP = 15 )
00383 *     ..
00384 *     .. Local Scalars ..
00385       LOGICAL            BADMM, BADNN, BADNNB
00386       INTEGER            I, IINFO, IMODE, ITYPE, J, JCOL, JR, JSIZE,
00387      $                   JTYPE, JWIDTH, K, KL, KMAX, KU, M, MMAX, MNMAX,
00388      $                   MNMIN, MTYPES, N, NERRS, NMATS, NMAX, NTEST,
00389      $                   NTESTT
00390       DOUBLE PRECISION   AMNINV, ANORM, COND, OVFL, RTOVFL, RTUNFL, ULP,
00391      $                   ULPINV, UNFL
00392 *     ..
00393 *     .. Local Arrays ..
00394       INTEGER            IDUMMA( 1 ), IOLDSD( 4 ), KMAGN( MAXTYP ),
00395      $                   KMODE( MAXTYP ), KTYPE( MAXTYP )
00396 *     ..
00397 *     .. External Functions ..
00398       DOUBLE PRECISION   DLAMCH
00399       EXTERNAL           DLAMCH
00400 *     ..
00401 *     .. External Subroutines ..
00402       EXTERNAL           DBDT01, DBDT02, DGBBRD, DLACPY, DLAHD2, DLASET,
00403      $                   DLASUM, DLATMR, DLATMS, DORT01, XERBLA
00404 *     ..
00405 *     .. Intrinsic Functions ..
00406       INTRINSIC          ABS, DBLE, MAX, MIN, SQRT
00407 *     ..
00408 *     .. Data statements ..
00409       DATA               KTYPE / 1, 2, 5*4, 5*6, 3*9 /
00410       DATA               KMAGN / 2*1, 3*1, 2, 3, 3*1, 2, 3, 1, 2, 3 /
00411       DATA               KMODE / 2*0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
00412      $                   0, 0 /
00413 *     ..
00414 *     .. Executable Statements ..
00415 *
00416 *     Check for errors
00417 *
00418       NTESTT = 0
00419       INFO = 0
00420 *
00421 *     Important constants
00422 *
00423       BADMM = .FALSE.
00424       BADNN = .FALSE.
00425       MMAX = 1
00426       NMAX = 1
00427       MNMAX = 1
00428       DO 10 J = 1, NSIZES
00429          MMAX = MAX( MMAX, MVAL( J ) )
00430          IF( MVAL( J ).LT.0 )
00431      $      BADMM = .TRUE.
00432          NMAX = MAX( NMAX, NVAL( J ) )
00433          IF( NVAL( J ).LT.0 )
00434      $      BADNN = .TRUE.
00435          MNMAX = MAX( MNMAX, MIN( MVAL( J ), NVAL( J ) ) )
00436    10 CONTINUE
00437 *
00438       BADNNB = .FALSE.
00439       KMAX = 0
00440       DO 20 J = 1, NWDTHS
00441          KMAX = MAX( KMAX, KK( J ) )
00442          IF( KK( J ).LT.0 )
00443      $      BADNNB = .TRUE.
00444    20 CONTINUE
00445 *
00446 *     Check for errors
00447 *
00448       IF( NSIZES.LT.0 ) THEN
00449          INFO = -1
00450       ELSE IF( BADMM ) THEN
00451          INFO = -2
00452       ELSE IF( BADNN ) THEN
00453          INFO = -3
00454       ELSE IF( NWDTHS.LT.0 ) THEN
00455          INFO = -4
00456       ELSE IF( BADNNB ) THEN
00457          INFO = -5
00458       ELSE IF( NTYPES.LT.0 ) THEN
00459          INFO = -6
00460       ELSE IF( NRHS.LT.0 ) THEN
00461          INFO = -8
00462       ELSE IF( LDA.LT.NMAX ) THEN
00463          INFO = -13
00464       ELSE IF( LDAB.LT.2*KMAX+1 ) THEN
00465          INFO = -15
00466       ELSE IF( LDQ.LT.NMAX ) THEN
00467          INFO = -19
00468       ELSE IF( LDP.LT.NMAX ) THEN
00469          INFO = -21
00470       ELSE IF( LDC.LT.NMAX ) THEN
00471          INFO = -23
00472       ELSE IF( ( MAX( LDA, NMAX )+1 )*NMAX.GT.LWORK ) THEN
00473          INFO = -26
00474       END IF
00475 *
00476       IF( INFO.NE.0 ) THEN
00477          CALL XERBLA( 'DCHKBB', -INFO )
00478          RETURN
00479       END IF
00480 *
00481 *     Quick return if possible
00482 *
00483       IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 .OR. NWDTHS.EQ.0 )
00484      $   RETURN
00485 *
00486 *     More Important constants
00487 *
00488       UNFL = DLAMCH( 'Safe minimum' )
00489       OVFL = ONE / UNFL
00490       ULP = DLAMCH( 'Epsilon' )*DLAMCH( 'Base' )
00491       ULPINV = ONE / ULP
00492       RTUNFL = SQRT( UNFL )
00493       RTOVFL = SQRT( OVFL )
00494 *
00495 *     Loop over sizes, widths, types
00496 *
00497       NERRS = 0
00498       NMATS = 0
00499 *
00500       DO 160 JSIZE = 1, NSIZES
00501          M = MVAL( JSIZE )
00502          N = NVAL( JSIZE )
00503          MNMIN = MIN( M, N )
00504          AMNINV = ONE / DBLE( MAX( 1, M, N ) )
00505 *
00506          DO 150 JWIDTH = 1, NWDTHS
00507             K = KK( JWIDTH )
00508             IF( K.GE.M .AND. K.GE.N )
00509      $         GO TO 150
00510             KL = MAX( 0, MIN( M-1, K ) )
00511             KU = MAX( 0, MIN( N-1, K ) )
00512 *
00513             IF( NSIZES.NE.1 ) THEN
00514                MTYPES = MIN( MAXTYP, NTYPES )
00515             ELSE
00516                MTYPES = MIN( MAXTYP+1, NTYPES )
00517             END IF
00518 *
00519             DO 140 JTYPE = 1, MTYPES
00520                IF( .NOT.DOTYPE( JTYPE ) )
00521      $            GO TO 140
00522                NMATS = NMATS + 1
00523                NTEST = 0
00524 *
00525                DO 30 J = 1, 4
00526                   IOLDSD( J ) = ISEED( J )
00527    30          CONTINUE
00528 *
00529 *              Compute "A".
00530 *
00531 *              Control parameters:
00532 *
00533 *                  KMAGN  KMODE        KTYPE
00534 *              =1  O(1)   clustered 1  zero
00535 *              =2  large  clustered 2  identity
00536 *              =3  small  exponential  (none)
00537 *              =4         arithmetic   diagonal, (w/ singular values)
00538 *              =5         random log   (none)
00539 *              =6         random       nonhermitian, w/ singular values
00540 *              =7                      (none)
00541 *              =8                      (none)
00542 *              =9                      random nonhermitian
00543 *
00544                IF( MTYPES.GT.MAXTYP )
00545      $            GO TO 90
00546 *
00547                ITYPE = KTYPE( JTYPE )
00548                IMODE = KMODE( JTYPE )
00549 *
00550 *              Compute norm
00551 *
00552                GO TO ( 40, 50, 60 )KMAGN( JTYPE )
00553 *
00554    40          CONTINUE
00555                ANORM = ONE
00556                GO TO 70
00557 *
00558    50          CONTINUE
00559                ANORM = ( RTOVFL*ULP )*AMNINV
00560                GO TO 70
00561 *
00562    60          CONTINUE
00563                ANORM = RTUNFL*MAX( M, N )*ULPINV
00564                GO TO 70
00565 *
00566    70          CONTINUE
00567 *
00568                CALL DLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA )
00569                CALL DLASET( 'Full', LDAB, N, ZERO, ZERO, AB, LDAB )
00570                IINFO = 0
00571                COND = ULPINV
00572 *
00573 *              Special Matrices -- Identity & Jordan block
00574 *
00575 *                 Zero
00576 *
00577                IF( ITYPE.EQ.1 ) THEN
00578                   IINFO = 0
00579 *
00580                ELSE IF( ITYPE.EQ.2 ) THEN
00581 *
00582 *                 Identity
00583 *
00584                   DO 80 JCOL = 1, N
00585                      A( JCOL, JCOL ) = ANORM
00586    80             CONTINUE
00587 *
00588                ELSE IF( ITYPE.EQ.4 ) THEN
00589 *
00590 *                 Diagonal Matrix, singular values specified
00591 *
00592                   CALL DLATMS( M, N, 'S', ISEED, 'N', WORK, IMODE, COND,
00593      $                         ANORM, 0, 0, 'N', A, LDA, WORK( M+1 ),
00594      $                         IINFO )
00595 *
00596                ELSE IF( ITYPE.EQ.6 ) THEN
00597 *
00598 *                 Nonhermitian, singular values specified
00599 *
00600                   CALL DLATMS( M, N, 'S', ISEED, 'N', WORK, IMODE, COND,
00601      $                         ANORM, KL, KU, 'N', A, LDA, WORK( M+1 ),
00602      $                         IINFO )
00603 *
00604                ELSE IF( ITYPE.EQ.9 ) THEN
00605 *
00606 *                 Nonhermitian, random entries
00607 *
00608                   CALL DLATMR( M, N, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
00609      $                         'T', 'N', WORK( N+1 ), 1, ONE,
00610      $                         WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, KL,
00611      $                         KU, ZERO, ANORM, 'N', A, LDA, IDUMMA,
00612      $                         IINFO )
00613 *
00614                ELSE
00615 *
00616                   IINFO = 1
00617                END IF
00618 *
00619 *              Generate Right-Hand Side
00620 *
00621                CALL DLATMR( M, NRHS, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
00622      $                      'T', 'N', WORK( M+1 ), 1, ONE,
00623      $                      WORK( 2*M+1 ), 1, ONE, 'N', IDUMMA, M, NRHS,
00624      $                      ZERO, ONE, 'NO', C, LDC, IDUMMA, IINFO )
00625 *
00626                IF( IINFO.NE.0 ) THEN
00627                   WRITE( NOUNIT, FMT = 9999 )'Generator', IINFO, N,
00628      $               JTYPE, IOLDSD
00629                   INFO = ABS( IINFO )
00630                   RETURN
00631                END IF
00632 *
00633    90          CONTINUE
00634 *
00635 *              Copy A to band storage.
00636 *
00637                DO 110 J = 1, N
00638                   DO 100 I = MAX( 1, J-KU ), MIN( M, J+KL )
00639                      AB( KU+1+I-J, J ) = A( I, J )
00640   100             CONTINUE
00641   110          CONTINUE
00642 *
00643 *              Copy C
00644 *
00645                CALL DLACPY( 'Full', M, NRHS, C, LDC, CC, LDC )
00646 *
00647 *              Call DGBBRD to compute B, Q and P, and to update C.
00648 *
00649                CALL DGBBRD( 'B', M, N, NRHS, KL, KU, AB, LDAB, BD, BE,
00650      $                      Q, LDQ, P, LDP, CC, LDC, WORK, IINFO )
00651 *
00652                IF( IINFO.NE.0 ) THEN
00653                   WRITE( NOUNIT, FMT = 9999 )'DGBBRD', IINFO, N, JTYPE,
00654      $               IOLDSD
00655                   INFO = ABS( IINFO )
00656                   IF( IINFO.LT.0 ) THEN
00657                      RETURN
00658                   ELSE
00659                      RESULT( 1 ) = ULPINV
00660                      GO TO 120
00661                   END IF
00662                END IF
00663 *
00664 *              Test 1:  Check the decomposition A := Q * B * P'
00665 *                   2:  Check the orthogonality of Q
00666 *                   3:  Check the orthogonality of P
00667 *                   4:  Check the computation of Q' * C
00668 *
00669                CALL DBDT01( M, N, -1, A, LDA, Q, LDQ, BD, BE, P, LDP,
00670      $                      WORK, RESULT( 1 ) )
00671                CALL DORT01( 'Columns', M, M, Q, LDQ, WORK, LWORK,
00672      $                      RESULT( 2 ) )
00673                CALL DORT01( 'Rows', N, N, P, LDP, WORK, LWORK,
00674      $                      RESULT( 3 ) )
00675                CALL DBDT02( M, NRHS, C, LDC, CC, LDC, Q, LDQ, WORK,
00676      $                      RESULT( 4 ) )
00677 *
00678 *              End of Loop -- Check for RESULT(j) > THRESH
00679 *
00680                NTEST = 4
00681   120          CONTINUE
00682                NTESTT = NTESTT + NTEST
00683 *
00684 *              Print out tests which fail.
00685 *
00686                DO 130 JR = 1, NTEST
00687                   IF( RESULT( JR ).GE.THRESH ) THEN
00688                      IF( NERRS.EQ.0 )
00689      $                  CALL DLAHD2( NOUNIT, 'DBB' )
00690                      NERRS = NERRS + 1
00691                      WRITE( NOUNIT, FMT = 9998 )M, N, K, IOLDSD, JTYPE,
00692      $                  JR, RESULT( JR )
00693                   END IF
00694   130          CONTINUE
00695 *
00696   140       CONTINUE
00697   150    CONTINUE
00698   160 CONTINUE
00699 *
00700 *     Summary
00701 *
00702       CALL DLASUM( 'DBB', NOUNIT, NERRS, NTESTT )
00703       RETURN
00704 *
00705  9999 FORMAT( ' DCHKBB: ', A, ' returned INFO=', I5, '.', / 9X, 'M=',
00706      $      I5, ' N=', I5, ' K=', I5, ', JTYPE=', I5, ', ISEED=(',
00707      $      3( I5, ',' ), I5, ')' )
00708  9998 FORMAT( ' M =', I4, ' N=', I4, ', K=', I3, ', seed=',
00709      $      4( I4, ',' ), ' type ', I2, ', test(', I2, ')=', G10.3 )
00710 *
00711 *     End of DCHKBB
00712 *
00713       END
 All Files Functions