LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
dchkbb.f
Go to the documentation of this file.
1 *> \brief \b DCHKBB
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE DCHKBB( NSIZES, MVAL, NVAL, NWDTHS, KK, NTYPES, DOTYPE,
12 * NRHS, ISEED, THRESH, NOUNIT, A, LDA, AB, LDAB,
13 * BD, BE, Q, LDQ, P, LDP, C, LDC, CC, WORK,
14 * LWORK, RESULT, INFO )
15 *
16 * .. Scalar Arguments ..
17 * INTEGER INFO, LDA, LDAB, LDC, LDP, LDQ, LWORK, NOUNIT,
18 * $ NRHS, NSIZES, NTYPES, NWDTHS
19 * DOUBLE PRECISION THRESH
20 * ..
21 * .. Array Arguments ..
22 * LOGICAL DOTYPE( * )
23 * INTEGER ISEED( 4 ), KK( * ), MVAL( * ), NVAL( * )
24 * DOUBLE PRECISION A( LDA, * ), AB( LDAB, * ), BD( * ), BE( * ),
25 * $ C( LDC, * ), CC( LDC, * ), P( LDP, * ),
26 * $ Q( LDQ, * ), RESULT( * ), WORK( * )
27 * ..
28 *
29 *
30 *> \par Purpose:
31 * =============
32 *>
33 *> \verbatim
34 *>
35 *> DCHKBB tests the reduction of a general real rectangular band
36 *> matrix to bidiagonal form.
37 *>
38 *> DGBBRD factors a general band matrix A as Q B P* , where * means
39 *> transpose, B is upper bidiagonal, and Q and P are orthogonal;
40 *> DGBBRD can also overwrite a given matrix C with Q* C .
41 *>
42 *> For each pair of matrix dimensions (M,N) and each selected matrix
43 *> type, an M by N matrix A and an M by NRHS matrix C are generated.
44 *> The problem dimensions are as follows
45 *> A: M x N
46 *> Q: M x M
47 *> P: N x N
48 *> B: min(M,N) x min(M,N)
49 *> C: M x NRHS
50 *>
51 *> For each generated matrix, 4 tests are performed:
52 *>
53 *> (1) | A - Q B PT | / ( |A| max(M,N) ulp ), PT = P'
54 *>
55 *> (2) | I - Q' Q | / ( M ulp )
56 *>
57 *> (3) | I - PT PT' | / ( N ulp )
58 *>
59 *> (4) | Y - Q' C | / ( |Y| max(M,NRHS) ulp ), where Y = Q' C.
60 *>
61 *> The "types" are specified by a logical array DOTYPE( 1:NTYPES );
62 *> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
63 *> Currently, the list of possible types is:
64 *>
65 *> The possible matrix types are
66 *>
67 *> (1) The zero matrix.
68 *> (2) The identity matrix.
69 *>
70 *> (3) A diagonal matrix with evenly spaced entries
71 *> 1, ..., ULP and random signs.
72 *> (ULP = (first number larger than 1) - 1 )
73 *> (4) A diagonal matrix with geometrically spaced entries
74 *> 1, ..., ULP and random signs.
75 *> (5) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
76 *> and random signs.
77 *>
78 *> (6) Same as (3), but multiplied by SQRT( overflow threshold )
79 *> (7) Same as (3), but multiplied by SQRT( underflow threshold )
80 *>
81 *> (8) A matrix of the form U D V, where U and V are orthogonal and
82 *> D has evenly spaced entries 1, ..., ULP with random signs
83 *> on the diagonal.
84 *>
85 *> (9) A matrix of the form U D V, where U and V are orthogonal and
86 *> D has geometrically spaced entries 1, ..., ULP with random
87 *> signs on the diagonal.
88 *>
89 *> (10) A matrix of the form U D V, where U and V are orthogonal and
90 *> D has "clustered" entries 1, ULP,..., ULP with random
91 *> signs on the diagonal.
92 *>
93 *> (11) Same as (8), but multiplied by SQRT( overflow threshold )
94 *> (12) Same as (8), but multiplied by SQRT( underflow threshold )
95 *>
96 *> (13) Rectangular matrix with random entries chosen from (-1,1).
97 *> (14) Same as (13), but multiplied by SQRT( overflow threshold )
98 *> (15) Same as (13), but multiplied by SQRT( underflow threshold )
99 *> \endverbatim
100 *
101 * Arguments:
102 * ==========
103 *
104 *> \param[in] NSIZES
105 *> \verbatim
106 *> NSIZES is INTEGER
107 *> The number of values of M and N contained in the vectors
108 *> MVAL and NVAL. The matrix sizes are used in pairs (M,N).
109 *> If NSIZES is zero, DCHKBB does nothing. NSIZES must be at
110 *> least zero.
111 *> \endverbatim
112 *>
113 *> \param[in] MVAL
114 *> \verbatim
115 *> MVAL is INTEGER array, dimension (NSIZES)
116 *> The values of the matrix row dimension M.
117 *> \endverbatim
118 *>
119 *> \param[in] NVAL
120 *> \verbatim
121 *> NVAL is INTEGER array, dimension (NSIZES)
122 *> The values of the matrix column dimension N.
123 *> \endverbatim
124 *>
125 *> \param[in] NWDTHS
126 *> \verbatim
127 *> NWDTHS is INTEGER
128 *> The number of bandwidths to use. If it is zero,
129 *> DCHKBB does nothing. It must be at least zero.
130 *> \endverbatim
131 *>
132 *> \param[in] KK
133 *> \verbatim
134 *> KK is INTEGER array, dimension (NWDTHS)
135 *> An array containing the bandwidths to be used for the band
136 *> matrices. The values must be at least zero.
137 *> \endverbatim
138 *>
139 *> \param[in] NTYPES
140 *> \verbatim
141 *> NTYPES is INTEGER
142 *> The number of elements in DOTYPE. If it is zero, DCHKBB
143 *> does nothing. It must be at least zero. If it is MAXTYP+1
144 *> and NSIZES is 1, then an additional type, MAXTYP+1 is
145 *> defined, which is to use whatever matrix is in A. This
146 *> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
147 *> DOTYPE(MAXTYP+1) is .TRUE. .
148 *> \endverbatim
149 *>
150 *> \param[in] DOTYPE
151 *> \verbatim
152 *> DOTYPE is LOGICAL array, dimension (NTYPES)
153 *> If DOTYPE(j) is .TRUE., then for each size in NN a
154 *> matrix of that size and of type j will be generated.
155 *> If NTYPES is smaller than the maximum number of types
156 *> defined (PARAMETER MAXTYP), then types NTYPES+1 through
157 *> MAXTYP will not be generated. If NTYPES is larger
158 *> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
159 *> will be ignored.
160 *> \endverbatim
161 *>
162 *> \param[in] NRHS
163 *> \verbatim
164 *> NRHS is INTEGER
165 *> The number of columns in the "right-hand side" matrix C.
166 *> If NRHS = 0, then the operations on the right-hand side will
167 *> not be tested. NRHS must be at least 0.
168 *> \endverbatim
169 *>
170 *> \param[in,out] ISEED
171 *> \verbatim
172 *> ISEED is INTEGER array, dimension (4)
173 *> On entry ISEED specifies the seed of the random number
174 *> generator. The array elements should be between 0 and 4095;
175 *> if not they will be reduced mod 4096. Also, ISEED(4) must
176 *> be odd. The random number generator uses a linear
177 *> congruential sequence limited to small integers, and so
178 *> should produce machine independent random numbers. The
179 *> values of ISEED are changed on exit, and can be used in the
180 *> next call to DCHKBB to continue the same random number
181 *> sequence.
182 *> \endverbatim
183 *>
184 *> \param[in] THRESH
185 *> \verbatim
186 *> THRESH is DOUBLE PRECISION
187 *> A test will count as "failed" if the "error", computed as
188 *> described above, exceeds THRESH. Note that the error
189 *> is scaled to be O(1), so THRESH should be a reasonably
190 *> small multiple of 1, e.g., 10 or 100. In particular,
191 *> it should not depend on the precision (single vs. double)
192 *> or the size of the matrix. It must be at least zero.
193 *> \endverbatim
194 *>
195 *> \param[in] NOUNIT
196 *> \verbatim
197 *> NOUNIT is INTEGER
198 *> The FORTRAN unit number for printing out error messages
199 *> (e.g., if a routine returns IINFO not equal to 0.)
200 *> \endverbatim
201 *>
202 *> \param[in,out] A
203 *> \verbatim
204 *> A is DOUBLE PRECISION array, dimension
205 *> (LDA, max(NN))
206 *> Used to hold the matrix A.
207 *> \endverbatim
208 *>
209 *> \param[in] LDA
210 *> \verbatim
211 *> LDA is INTEGER
212 *> The leading dimension of A. It must be at least 1
213 *> and at least max( NN ).
214 *> \endverbatim
215 *>
216 *> \param[out] AB
217 *> \verbatim
218 *> AB is DOUBLE PRECISION array, dimension (LDAB, max(NN))
219 *> Used to hold A in band storage format.
220 *> \endverbatim
221 *>
222 *> \param[in] LDAB
223 *> \verbatim
224 *> LDAB is INTEGER
225 *> The leading dimension of AB. It must be at least 2 (not 1!)
226 *> and at least max( KK )+1.
227 *> \endverbatim
228 *>
229 *> \param[out] BD
230 *> \verbatim
231 *> BD is DOUBLE PRECISION array, dimension (max(NN))
232 *> Used to hold the diagonal of the bidiagonal matrix computed
233 *> by DGBBRD.
234 *> \endverbatim
235 *>
236 *> \param[out] BE
237 *> \verbatim
238 *> BE is DOUBLE PRECISION array, dimension (max(NN))
239 *> Used to hold the off-diagonal of the bidiagonal matrix
240 *> computed by DGBBRD.
241 *> \endverbatim
242 *>
243 *> \param[out] Q
244 *> \verbatim
245 *> Q is DOUBLE PRECISION array, dimension (LDQ, max(NN))
246 *> Used to hold the orthogonal matrix Q computed by DGBBRD.
247 *> \endverbatim
248 *>
249 *> \param[in] LDQ
250 *> \verbatim
251 *> LDQ is INTEGER
252 *> The leading dimension of Q. It must be at least 1
253 *> and at least max( NN ).
254 *> \endverbatim
255 *>
256 *> \param[out] P
257 *> \verbatim
258 *> P is DOUBLE PRECISION array, dimension (LDP, max(NN))
259 *> Used to hold the orthogonal matrix P computed by DGBBRD.
260 *> \endverbatim
261 *>
262 *> \param[in] LDP
263 *> \verbatim
264 *> LDP is INTEGER
265 *> The leading dimension of P. It must be at least 1
266 *> and at least max( NN ).
267 *> \endverbatim
268 *>
269 *> \param[out] C
270 *> \verbatim
271 *> C is DOUBLE PRECISION array, dimension (LDC, max(NN))
272 *> Used to hold the matrix C updated by DGBBRD.
273 *> \endverbatim
274 *>
275 *> \param[in] LDC
276 *> \verbatim
277 *> LDC is INTEGER
278 *> The leading dimension of U. It must be at least 1
279 *> and at least max( NN ).
280 *> \endverbatim
281 *>
282 *> \param[out] CC
283 *> \verbatim
284 *> CC is DOUBLE PRECISION array, dimension (LDC, max(NN))
285 *> Used to hold a copy of the matrix C.
286 *> \endverbatim
287 *>
288 *> \param[out] WORK
289 *> \verbatim
290 *> WORK is DOUBLE PRECISION array, dimension (LWORK)
291 *> \endverbatim
292 *>
293 *> \param[in] LWORK
294 *> \verbatim
295 *> LWORK is INTEGER
296 *> The number of entries in WORK. This must be at least
297 *> max( LDA+1, max(NN)+1 )*max(NN).
298 *> \endverbatim
299 *>
300 *> \param[out] RESULT
301 *> \verbatim
302 *> RESULT is DOUBLE PRECISION array, dimension (4)
303 *> The values computed by the tests described above.
304 *> The values are currently limited to 1/ulp, to avoid
305 *> overflow.
306 *> \endverbatim
307 *>
308 *> \param[out] INFO
309 *> \verbatim
310 *> INFO is INTEGER
311 *> If 0, then everything ran OK.
312 *>
313 *>-----------------------------------------------------------------------
314 *>
315 *> Some Local Variables and Parameters:
316 *> ---- ----- --------- --- ----------
317 *> ZERO, ONE Real 0 and 1.
318 *> MAXTYP The number of types defined.
319 *> NTEST The number of tests performed, or which can
320 *> be performed so far, for the current matrix.
321 *> NTESTT The total number of tests performed so far.
322 *> NMAX Largest value in NN.
323 *> NMATS The number of matrices generated so far.
324 *> NERRS The number of tests which have exceeded THRESH
325 *> so far.
326 *> COND, IMODE Values to be passed to the matrix generators.
327 *> ANORM Norm of A; passed to matrix generators.
328 *>
329 *> OVFL, UNFL Overflow and underflow thresholds.
330 *> ULP, ULPINV Finest relative precision and its inverse.
331 *> RTOVFL, RTUNFL Square roots of the previous 2 values.
332 *> The following four arrays decode JTYPE:
333 *> KTYPE(j) The general type (1-10) for type "j".
334 *> KMODE(j) The MODE value to be passed to the matrix
335 *> generator for type "j".
336 *> KMAGN(j) The order of magnitude ( O(1),
337 *> O(overflow^(1/2) ), O(underflow^(1/2) )
338 *> \endverbatim
339 *
340 * Authors:
341 * ========
342 *
343 *> \author Univ. of Tennessee
344 *> \author Univ. of California Berkeley
345 *> \author Univ. of Colorado Denver
346 *> \author NAG Ltd.
347 *
348 *> \date November 2011
349 *
350 *> \ingroup double_eig
351 *
352 * =====================================================================
353  SUBROUTINE dchkbb( NSIZES, MVAL, NVAL, NWDTHS, KK, NTYPES, DOTYPE,
354  $ nrhs, iseed, thresh, nounit, a, lda, ab, ldab,
355  $ bd, be, q, ldq, p, ldp, c, ldc, cc, work,
356  $ lwork, result, info )
357 *
358 * -- LAPACK test routine (input) --
359 * -- LAPACK is a software package provided by Univ. of Tennessee, --
360 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
361 * November 2011
362 *
363 * .. Scalar Arguments ..
364  INTEGER info, lda, ldab, ldc, ldp, ldq, lwork, nounit,
365  $ nrhs, nsizes, ntypes, nwdths
366  DOUBLE PRECISION thresh
367 * ..
368 * .. Array Arguments ..
369  LOGICAL dotype( * )
370  INTEGER iseed( 4 ), kk( * ), mval( * ), nval( * )
371  DOUBLE PRECISION a( lda, * ), ab( ldab, * ), bd( * ), be( * ),
372  $ c( ldc, * ), cc( ldc, * ), p( ldp, * ),
373  $ q( ldq, * ), result( * ), work( * )
374 * ..
375 *
376 * =====================================================================
377 *
378 * .. Parameters ..
379  DOUBLE PRECISION zero, one
380  parameter( zero = 0.0d0, one = 1.0d0 )
381  INTEGER maxtyp
382  parameter( maxtyp = 15 )
383 * ..
384 * .. Local Scalars ..
385  LOGICAL badmm, badnn, badnnb
386  INTEGER i, iinfo, imode, itype, j, jcol, jr, jsize,
387  $ jtype, jwidth, k, kl, kmax, ku, m, mmax, mnmax,
388  $ mnmin, mtypes, n, nerrs, nmats, nmax, ntest,
389  $ ntestt
390  DOUBLE PRECISION amninv, anorm, cond, ovfl, rtovfl, rtunfl, ulp,
391  $ ulpinv, unfl
392 * ..
393 * .. Local Arrays ..
394  INTEGER idumma( 1 ), ioldsd( 4 ), kmagn( maxtyp ),
395  $ kmode( maxtyp ), ktype( maxtyp )
396 * ..
397 * .. External Functions ..
398  DOUBLE PRECISION dlamch
399  EXTERNAL dlamch
400 * ..
401 * .. External Subroutines ..
402  EXTERNAL dbdt01, dbdt02, dgbbrd, dlacpy, dlahd2, dlaset,
403  $ dlasum, dlatmr, dlatms, dort01, xerbla
404 * ..
405 * .. Intrinsic Functions ..
406  INTRINSIC abs, dble, max, min, sqrt
407 * ..
408 * .. Data statements ..
409  DATA ktype / 1, 2, 5*4, 5*6, 3*9 /
410  DATA kmagn / 2*1, 3*1, 2, 3, 3*1, 2, 3, 1, 2, 3 /
411  DATA kmode / 2*0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
412  $ 0, 0 /
413 * ..
414 * .. Executable Statements ..
415 *
416 * Check for errors
417 *
418  ntestt = 0
419  info = 0
420 *
421 * Important constants
422 *
423  badmm = .false.
424  badnn = .false.
425  mmax = 1
426  nmax = 1
427  mnmax = 1
428  DO 10 j = 1, nsizes
429  mmax = max( mmax, mval( j ) )
430  IF( mval( j ).LT.0 )
431  $ badmm = .true.
432  nmax = max( nmax, nval( j ) )
433  IF( nval( j ).LT.0 )
434  $ badnn = .true.
435  mnmax = max( mnmax, min( mval( j ), nval( j ) ) )
436  10 continue
437 *
438  badnnb = .false.
439  kmax = 0
440  DO 20 j = 1, nwdths
441  kmax = max( kmax, kk( j ) )
442  IF( kk( j ).LT.0 )
443  $ badnnb = .true.
444  20 continue
445 *
446 * Check for errors
447 *
448  IF( nsizes.LT.0 ) THEN
449  info = -1
450  ELSE IF( badmm ) THEN
451  info = -2
452  ELSE IF( badnn ) THEN
453  info = -3
454  ELSE IF( nwdths.LT.0 ) THEN
455  info = -4
456  ELSE IF( badnnb ) THEN
457  info = -5
458  ELSE IF( ntypes.LT.0 ) THEN
459  info = -6
460  ELSE IF( nrhs.LT.0 ) THEN
461  info = -8
462  ELSE IF( lda.LT.nmax ) THEN
463  info = -13
464  ELSE IF( ldab.LT.2*kmax+1 ) THEN
465  info = -15
466  ELSE IF( ldq.LT.nmax ) THEN
467  info = -19
468  ELSE IF( ldp.LT.nmax ) THEN
469  info = -21
470  ELSE IF( ldc.LT.nmax ) THEN
471  info = -23
472  ELSE IF( ( max( lda, nmax )+1 )*nmax.GT.lwork ) THEN
473  info = -26
474  END IF
475 *
476  IF( info.NE.0 ) THEN
477  CALL xerbla( 'DCHKBB', -info )
478  return
479  END IF
480 *
481 * Quick return if possible
482 *
483  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 .OR. nwdths.EQ.0 )
484  $ return
485 *
486 * More Important constants
487 *
488  unfl = dlamch( 'Safe minimum' )
489  ovfl = one / unfl
490  ulp = dlamch( 'Epsilon' )*dlamch( 'Base' )
491  ulpinv = one / ulp
492  rtunfl = sqrt( unfl )
493  rtovfl = sqrt( ovfl )
494 *
495 * Loop over sizes, widths, types
496 *
497  nerrs = 0
498  nmats = 0
499 *
500  DO 160 jsize = 1, nsizes
501  m = mval( jsize )
502  n = nval( jsize )
503  mnmin = min( m, n )
504  amninv = one / dble( max( 1, m, n ) )
505 *
506  DO 150 jwidth = 1, nwdths
507  k = kk( jwidth )
508  IF( k.GE.m .AND. k.GE.n )
509  $ go to 150
510  kl = max( 0, min( m-1, k ) )
511  ku = max( 0, min( n-1, k ) )
512 *
513  IF( nsizes.NE.1 ) THEN
514  mtypes = min( maxtyp, ntypes )
515  ELSE
516  mtypes = min( maxtyp+1, ntypes )
517  END IF
518 *
519  DO 140 jtype = 1, mtypes
520  IF( .NOT.dotype( jtype ) )
521  $ go to 140
522  nmats = nmats + 1
523  ntest = 0
524 *
525  DO 30 j = 1, 4
526  ioldsd( j ) = iseed( j )
527  30 continue
528 *
529 * Compute "A".
530 *
531 * Control parameters:
532 *
533 * KMAGN KMODE KTYPE
534 * =1 O(1) clustered 1 zero
535 * =2 large clustered 2 identity
536 * =3 small exponential (none)
537 * =4 arithmetic diagonal, (w/ singular values)
538 * =5 random log (none)
539 * =6 random nonhermitian, w/ singular values
540 * =7 (none)
541 * =8 (none)
542 * =9 random nonhermitian
543 *
544  IF( mtypes.GT.maxtyp )
545  $ go to 90
546 *
547  itype = ktype( jtype )
548  imode = kmode( jtype )
549 *
550 * Compute norm
551 *
552  go to( 40, 50, 60 )kmagn( jtype )
553 *
554  40 continue
555  anorm = one
556  go to 70
557 *
558  50 continue
559  anorm = ( rtovfl*ulp )*amninv
560  go to 70
561 *
562  60 continue
563  anorm = rtunfl*max( m, n )*ulpinv
564  go to 70
565 *
566  70 continue
567 *
568  CALL dlaset( 'Full', lda, n, zero, zero, a, lda )
569  CALL dlaset( 'Full', ldab, n, zero, zero, ab, ldab )
570  iinfo = 0
571  cond = ulpinv
572 *
573 * Special Matrices -- Identity & Jordan block
574 *
575 * Zero
576 *
577  IF( itype.EQ.1 ) THEN
578  iinfo = 0
579 *
580  ELSE IF( itype.EQ.2 ) THEN
581 *
582 * Identity
583 *
584  DO 80 jcol = 1, n
585  a( jcol, jcol ) = anorm
586  80 continue
587 *
588  ELSE IF( itype.EQ.4 ) THEN
589 *
590 * Diagonal Matrix, singular values specified
591 *
592  CALL dlatms( m, n, 'S', iseed, 'N', work, imode, cond,
593  $ anorm, 0, 0, 'N', a, lda, work( m+1 ),
594  $ iinfo )
595 *
596  ELSE IF( itype.EQ.6 ) THEN
597 *
598 * Nonhermitian, singular values specified
599 *
600  CALL dlatms( m, n, 'S', iseed, 'N', work, imode, cond,
601  $ anorm, kl, ku, 'N', a, lda, work( m+1 ),
602  $ iinfo )
603 *
604  ELSE IF( itype.EQ.9 ) THEN
605 *
606 * Nonhermitian, random entries
607 *
608  CALL dlatmr( m, n, 'S', iseed, 'N', work, 6, one, one,
609  $ 'T', 'N', work( n+1 ), 1, one,
610  $ work( 2*n+1 ), 1, one, 'N', idumma, kl,
611  $ ku, zero, anorm, 'N', a, lda, idumma,
612  $ iinfo )
613 *
614  ELSE
615 *
616  iinfo = 1
617  END IF
618 *
619 * Generate Right-Hand Side
620 *
621  CALL dlatmr( m, nrhs, 'S', iseed, 'N', work, 6, one, one,
622  $ 'T', 'N', work( m+1 ), 1, one,
623  $ work( 2*m+1 ), 1, one, 'N', idumma, m, nrhs,
624  $ zero, one, 'NO', c, ldc, idumma, iinfo )
625 *
626  IF( iinfo.NE.0 ) THEN
627  WRITE( nounit, fmt = 9999 )'Generator', iinfo, n,
628  $ jtype, ioldsd
629  info = abs( iinfo )
630  return
631  END IF
632 *
633  90 continue
634 *
635 * Copy A to band storage.
636 *
637  DO 110 j = 1, n
638  DO 100 i = max( 1, j-ku ), min( m, j+kl )
639  ab( ku+1+i-j, j ) = a( i, j )
640  100 continue
641  110 continue
642 *
643 * Copy C
644 *
645  CALL dlacpy( 'Full', m, nrhs, c, ldc, cc, ldc )
646 *
647 * Call DGBBRD to compute B, Q and P, and to update C.
648 *
649  CALL dgbbrd( 'B', m, n, nrhs, kl, ku, ab, ldab, bd, be,
650  $ q, ldq, p, ldp, cc, ldc, work, iinfo )
651 *
652  IF( iinfo.NE.0 ) THEN
653  WRITE( nounit, fmt = 9999 )'DGBBRD', iinfo, n, jtype,
654  $ ioldsd
655  info = abs( iinfo )
656  IF( iinfo.LT.0 ) THEN
657  return
658  ELSE
659  result( 1 ) = ulpinv
660  go to 120
661  END IF
662  END IF
663 *
664 * Test 1: Check the decomposition A := Q * B * P'
665 * 2: Check the orthogonality of Q
666 * 3: Check the orthogonality of P
667 * 4: Check the computation of Q' * C
668 *
669  CALL dbdt01( m, n, -1, a, lda, q, ldq, bd, be, p, ldp,
670  $ work, result( 1 ) )
671  CALL dort01( 'Columns', m, m, q, ldq, work, lwork,
672  $ result( 2 ) )
673  CALL dort01( 'Rows', n, n, p, ldp, work, lwork,
674  $ result( 3 ) )
675  CALL dbdt02( m, nrhs, c, ldc, cc, ldc, q, ldq, work,
676  $ result( 4 ) )
677 *
678 * End of Loop -- Check for RESULT(j) > THRESH
679 *
680  ntest = 4
681  120 continue
682  ntestt = ntestt + ntest
683 *
684 * Print out tests which fail.
685 *
686  DO 130 jr = 1, ntest
687  IF( result( jr ).GE.thresh ) THEN
688  IF( nerrs.EQ.0 )
689  $ CALL dlahd2( nounit, 'DBB' )
690  nerrs = nerrs + 1
691  WRITE( nounit, fmt = 9998 )m, n, k, ioldsd, jtype,
692  $ jr, result( jr )
693  END IF
694  130 continue
695 *
696  140 continue
697  150 continue
698  160 continue
699 *
700 * Summary
701 *
702  CALL dlasum( 'DBB', nounit, nerrs, ntestt )
703  return
704 *
705  9999 format( ' DCHKBB: ', a, ' returned INFO=', i5, '.', / 9x, 'M=',
706  $ i5, ' N=', i5, ' K=', i5, ', JTYPE=', i5, ', ISEED=(',
707  $ 3( i5, ',' ), i5, ')' )
708  9998 format( ' M =', i4, ' N=', i4, ', K=', i3, ', seed=',
709  $ 4( i4, ',' ), ' type ', i2, ', test(', i2, ')=', g10.3 )
710 *
711 * End of DCHKBB
712 *
713  END