LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
schktb.f
Go to the documentation of this file.
1*> \brief \b SCHKTB
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 SCHKTB( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
12* NMAX, AB, AINV, B, X, XACT, WORK, RWORK, IWORK,
13* NOUT )
14*
15* .. Scalar Arguments ..
16* LOGICAL TSTERR
17* INTEGER NMAX, NN, NNS, NOUT
18* REAL THRESH
19* ..
20* .. Array Arguments ..
21* LOGICAL DOTYPE( * )
22* INTEGER IWORK( * ), NSVAL( * ), NVAL( * )
23* REAL AB( * ), AINV( * ), B( * ), RWORK( * ),
24* $ WORK( * ), X( * ), XACT( * )
25* ..
26*
27*
28*> \par Purpose:
29* =============
30*>
31*> \verbatim
32*>
33*> SCHKTB tests STBTRS, -RFS, and -CON, and SLATBS.
34*> \endverbatim
35*
36* Arguments:
37* ==========
38*
39*> \param[in] DOTYPE
40*> \verbatim
41*> DOTYPE is LOGICAL array, dimension (NTYPES)
42*> The matrix types to be used for testing. Matrices of type j
43*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
44*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
45*> \endverbatim
46*>
47*> \param[in] NN
48*> \verbatim
49*> NN is INTEGER
50*> The number of values of N contained in the vector NVAL.
51*> \endverbatim
52*>
53*> \param[in] NVAL
54*> \verbatim
55*> NVAL is INTEGER array, dimension (NN)
56*> The values of the matrix column dimension N.
57*> \endverbatim
58*>
59*> \param[in] NNS
60*> \verbatim
61*> NNS is INTEGER
62*> The number of values of NRHS contained in the vector NSVAL.
63*> \endverbatim
64*>
65*> \param[in] NSVAL
66*> \verbatim
67*> NSVAL is INTEGER array, dimension (NNS)
68*> The values of the number of right hand sides NRHS.
69*> \endverbatim
70*>
71*> \param[in] THRESH
72*> \verbatim
73*> THRESH is REAL
74*> The threshold value for the test ratios. A result is
75*> included in the output file if RESULT >= THRESH. To have
76*> every test ratio printed, use THRESH = 0.
77*> \endverbatim
78*>
79*> \param[in] TSTERR
80*> \verbatim
81*> TSTERR is LOGICAL
82*> Flag that indicates whether error exits are to be tested.
83*> \endverbatim
84*>
85*> \param[in] NMAX
86*> \verbatim
87*> NMAX is INTEGER
88*> The leading dimension of the work arrays.
89*> NMAX >= the maximum value of N in NVAL.
90*> \endverbatim
91*>
92*> \param[out] AB
93*> \verbatim
94*> AB is REAL array, dimension (NMAX*NMAX)
95*> \endverbatim
96*>
97*> \param[out] AINV
98*> \verbatim
99*> AINV is REAL array, dimension (NMAX*NMAX)
100*> \endverbatim
101*>
102*> \param[out] B
103*> \verbatim
104*> B is REAL array, dimension (NMAX*NSMAX)
105*> where NSMAX is the largest entry in NSVAL.
106*> \endverbatim
107*>
108*> \param[out] X
109*> \verbatim
110*> X is REAL array, dimension (NMAX*NSMAX)
111*> \endverbatim
112*>
113*> \param[out] XACT
114*> \verbatim
115*> XACT is REAL array, dimension (NMAX*NSMAX)
116*> \endverbatim
117*>
118*> \param[out] WORK
119*> \verbatim
120*> WORK is REAL array, dimension
121*> (NMAX*max(3,NSMAX))
122*> \endverbatim
123*>
124*> \param[out] RWORK
125*> \verbatim
126*> RWORK is REAL array, dimension
127*> (max(NMAX,2*NSMAX))
128*> \endverbatim
129*>
130*> \param[out] IWORK
131*> \verbatim
132*> IWORK is INTEGER array, dimension (NMAX)
133*> \endverbatim
134*>
135*> \param[in] NOUT
136*> \verbatim
137*> NOUT is INTEGER
138*> The unit number for output.
139*> \endverbatim
140*
141* Authors:
142* ========
143*
144*> \author Univ. of Tennessee
145*> \author Univ. of California Berkeley
146*> \author Univ. of Colorado Denver
147*> \author NAG Ltd.
148*
149*> \ingroup single_lin
150*
151* =====================================================================
152 SUBROUTINE schktb( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
153 $ NMAX, AB, AINV, B, X, XACT, WORK, RWORK, IWORK,
154 $ NOUT )
155*
156* -- LAPACK test routine --
157* -- LAPACK is a software package provided by Univ. of Tennessee, --
158* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
159*
160* .. Scalar Arguments ..
161 LOGICAL TSTERR
162 INTEGER NMAX, NN, NNS, NOUT
163 REAL THRESH
164* ..
165* .. Array Arguments ..
166 LOGICAL DOTYPE( * )
167 INTEGER IWORK( * ), NSVAL( * ), NVAL( * )
168 REAL AB( * ), AINV( * ), B( * ), RWORK( * ),
169 $ work( * ), x( * ), xact( * )
170* ..
171*
172* =====================================================================
173*
174* .. Parameters ..
175 INTEGER NTYPE1, NTYPES
176 PARAMETER ( NTYPE1 = 9, ntypes = 17 )
177 INTEGER NTESTS
178 parameter( ntests = 8 )
179 INTEGER NTRAN
180 parameter( ntran = 3 )
181 REAL ONE, ZERO
182 parameter( one = 1.0e+0, zero = 0.0e+0 )
183* ..
184* .. Local Scalars ..
185 CHARACTER DIAG, NORM, TRANS, UPLO, XTYPE
186 CHARACTER*3 PATH
187 INTEGER I, IDIAG, IK, IMAT, IN, INFO, IRHS, ITRAN,
188 $ iuplo, j, k, kd, lda, ldab, n, nerrs, nfail,
189 $ nimat, nimat2, nk, nrhs, nrun
190 REAL AINVNM, ANORM, RCOND, RCONDC, RCONDI, RCONDO,
191 $ SCALE
192* ..
193* .. Local Arrays ..
194 CHARACTER TRANSS( NTRAN ), UPLOS( 2 )
195 INTEGER ISEED( 4 ), ISEEDY( 4 )
196 REAL RESULT( NTESTS )
197* ..
198* .. External Functions ..
199 LOGICAL LSAME
200 REAL SLANTB, SLANTR
201 EXTERNAL lsame, slantb, slantr
202* ..
203* .. External Subroutines ..
204 EXTERNAL alaerh, alahd, alasum, scopy, serrtr, sget04,
207 $ stbtrs
208* ..
209* .. Scalars in Common ..
210 LOGICAL LERR, OK
211 CHARACTER*32 SRNAMT
212 INTEGER INFOT, IOUNIT
213* ..
214* .. Common blocks ..
215 COMMON / infoc / infot, iounit, ok, lerr
216 COMMON / srnamc / srnamt
217* ..
218* .. Intrinsic Functions ..
219 INTRINSIC max, min
220* ..
221* .. Data statements ..
222 DATA iseedy / 1988, 1989, 1990, 1991 /
223 DATA uplos / 'U', 'L' / , transs / 'N', 'T', 'C' /
224* ..
225* .. Executable Statements ..
226*
227* Initialize constants and the random number seed.
228*
229 path( 1: 1 ) = 'Single precision'
230 path( 2: 3 ) = 'TB'
231 nrun = 0
232 nfail = 0
233 nerrs = 0
234 DO 10 i = 1, 4
235 iseed( i ) = iseedy( i )
236 10 CONTINUE
237*
238* Test the error exits
239*
240 IF( tsterr )
241 $ CALL serrtr( path, nout )
242 infot = 0
243*
244 DO 140 in = 1, nn
245*
246* Do for each value of N in NVAL
247*
248 n = nval( in )
249 lda = max( 1, n )
250 xtype = 'N'
251 nimat = ntype1
252 nimat2 = ntypes
253 IF( n.LE.0 ) THEN
254 nimat = 1
255 nimat2 = ntype1 + 1
256 END IF
257*
258 nk = min( n+1, 4 )
259 DO 130 ik = 1, nk
260*
261* Do for KD = 0, N, (3N-1)/4, and (N+1)/4. This order makes
262* it easier to skip redundant values for small values of N.
263*
264 IF( ik.EQ.1 ) THEN
265 kd = 0
266 ELSE IF( ik.EQ.2 ) THEN
267 kd = max( n, 0 )
268 ELSE IF( ik.EQ.3 ) THEN
269 kd = ( 3*n-1 ) / 4
270 ELSE IF( ik.EQ.4 ) THEN
271 kd = ( n+1 ) / 4
272 END IF
273 ldab = kd + 1
274*
275 DO 90 imat = 1, nimat
276*
277* Do the tests only if DOTYPE( IMAT ) is true.
278*
279 IF( .NOT.dotype( imat ) )
280 $ GO TO 90
281*
282 DO 80 iuplo = 1, 2
283*
284* Do first for UPLO = 'U', then for UPLO = 'L'
285*
286 uplo = uplos( iuplo )
287*
288* Call SLATTB to generate a triangular test matrix.
289*
290 srnamt = 'SLATTB'
291 CALL slattb( imat, uplo, 'No transpose', diag, iseed,
292 $ n, kd, ab, ldab, x, work, info )
293*
294* Set IDIAG = 1 for non-unit matrices, 2 for unit.
295*
296 IF( lsame( diag, 'N' ) ) THEN
297 idiag = 1
298 ELSE
299 idiag = 2
300 END IF
301*
302* Form the inverse of A so we can get a good estimate
303* of RCONDC = 1/(norm(A) * norm(inv(A))).
304*
305 CALL slaset( 'Full', n, n, zero, one, ainv, lda )
306 IF( lsame( uplo, 'U' ) ) THEN
307 DO 20 j = 1, n
308 CALL stbsv( uplo, 'No transpose', diag, j, kd,
309 $ ab, ldab, ainv( ( j-1 )*lda+1 ), 1 )
310 20 CONTINUE
311 ELSE
312 DO 30 j = 1, n
313 CALL stbsv( uplo, 'No transpose', diag, n-j+1,
314 $ kd, ab( ( j-1 )*ldab+1 ), ldab,
315 $ ainv( ( j-1 )*lda+j ), 1 )
316 30 CONTINUE
317 END IF
318*
319* Compute the 1-norm condition number of A.
320*
321 anorm = slantb( '1', uplo, diag, n, kd, ab, ldab,
322 $ rwork )
323 ainvnm = slantr( '1', uplo, diag, n, n, ainv, lda,
324 $ rwork )
325 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
326 rcondo = one
327 ELSE
328 rcondo = ( one / anorm ) / ainvnm
329 END IF
330*
331* Compute the infinity-norm condition number of A.
332*
333 anorm = slantb( 'I', uplo, diag, n, kd, ab, ldab,
334 $ rwork )
335 ainvnm = slantr( 'I', uplo, diag, n, n, ainv, lda,
336 $ rwork )
337 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
338 rcondi = one
339 ELSE
340 rcondi = ( one / anorm ) / ainvnm
341 END IF
342*
343 DO 60 irhs = 1, nns
344 nrhs = nsval( irhs )
345 xtype = 'N'
346*
347 DO 50 itran = 1, ntran
348*
349* Do for op(A) = A, A**T, or A**H.
350*
351 trans = transs( itran )
352 IF( itran.EQ.1 ) THEN
353 norm = 'O'
354 rcondc = rcondo
355 ELSE
356 norm = 'I'
357 rcondc = rcondi
358 END IF
359*
360*+ TEST 1
361* Solve and compute residual for op(A)*x = b.
362*
363 srnamt = 'SLARHS'
364 CALL slarhs( path, xtype, uplo, trans, n, n, kd,
365 $ idiag, nrhs, ab, ldab, xact, lda,
366 $ b, lda, iseed, info )
367 xtype = 'C'
368 CALL slacpy( 'Full', n, nrhs, b, lda, x, lda )
369*
370 srnamt = 'STBTRS'
371 CALL stbtrs( uplo, trans, diag, n, kd, nrhs, ab,
372 $ ldab, x, lda, info )
373*
374* Check error code from STBTRS.
375*
376 IF( info.NE.0 )
377 $ CALL alaerh( path, 'STBTRS', info, 0,
378 $ uplo // trans // diag, n, n, kd,
379 $ kd, nrhs, imat, nfail, nerrs,
380 $ nout )
381*
382 CALL stbt02( uplo, trans, diag, n, kd, nrhs, ab,
383 $ ldab, x, lda, b, lda, work,
384 $ result( 1 ) )
385*
386*+ TEST 2
387* Check solution from generated exact solution.
388*
389 CALL sget04( n, nrhs, x, lda, xact, lda, rcondc,
390 $ result( 2 ) )
391*
392*+ TESTS 3, 4, and 5
393* Use iterative refinement to improve the solution
394* and compute error bounds.
395*
396 srnamt = 'STBRFS'
397 CALL stbrfs( uplo, trans, diag, n, kd, nrhs, ab,
398 $ ldab, b, lda, x, lda, rwork,
399 $ rwork( nrhs+1 ), work, iwork,
400 $ info )
401*
402* Check error code from STBRFS.
403*
404 IF( info.NE.0 )
405 $ CALL alaerh( path, 'STBRFS', info, 0,
406 $ uplo // trans // diag, n, n, kd,
407 $ kd, nrhs, imat, nfail, nerrs,
408 $ nout )
409*
410 CALL sget04( n, nrhs, x, lda, xact, lda, rcondc,
411 $ result( 3 ) )
412 CALL stbt05( uplo, trans, diag, n, kd, nrhs, ab,
413 $ ldab, b, lda, x, lda, xact, lda,
414 $ rwork, rwork( nrhs+1 ),
415 $ result( 4 ) )
416*
417* Print information about the tests that did not
418* pass the threshold.
419*
420 DO 40 k = 1, 5
421 IF( result( k ).GE.thresh ) THEN
422 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
423 $ CALL alahd( nout, path )
424 WRITE( nout, fmt = 9999 )uplo, trans,
425 $ diag, n, kd, nrhs, imat, k, result( k )
426 nfail = nfail + 1
427 END IF
428 40 CONTINUE
429 nrun = nrun + 5
430 50 CONTINUE
431 60 CONTINUE
432*
433*+ TEST 6
434* Get an estimate of RCOND = 1/CNDNUM.
435*
436 DO 70 itran = 1, 2
437 IF( itran.EQ.1 ) THEN
438 norm = 'O'
439 rcondc = rcondo
440 ELSE
441 norm = 'I'
442 rcondc = rcondi
443 END IF
444 srnamt = 'STBCON'
445 CALL stbcon( norm, uplo, diag, n, kd, ab, ldab,
446 $ rcond, work, iwork, info )
447*
448* Check error code from STBCON.
449*
450 IF( info.NE.0 )
451 $ CALL alaerh( path, 'STBCON', info, 0,
452 $ norm // uplo // diag, n, n, kd, kd,
453 $ -1, imat, nfail, nerrs, nout )
454*
455 CALL stbt06( rcond, rcondc, uplo, diag, n, kd, ab,
456 $ ldab, rwork, result( 6 ) )
457*
458* Print information about the tests that did not pass
459* the threshold.
460*
461 IF( result( 6 ).GE.thresh ) THEN
462 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
463 $ CALL alahd( nout, path )
464 WRITE( nout, fmt = 9998 ) 'STBCON', norm, uplo,
465 $ diag, n, kd, imat, 6, result( 6 )
466 nfail = nfail + 1
467 END IF
468 nrun = nrun + 1
469 70 CONTINUE
470 80 CONTINUE
471 90 CONTINUE
472*
473* Use pathological test matrices to test SLATBS.
474*
475 DO 120 imat = ntype1 + 1, nimat2
476*
477* Do the tests only if DOTYPE( IMAT ) is true.
478*
479 IF( .NOT.dotype( imat ) )
480 $ GO TO 120
481*
482 DO 110 iuplo = 1, 2
483*
484* Do first for UPLO = 'U', then for UPLO = 'L'
485*
486 uplo = uplos( iuplo )
487 DO 100 itran = 1, ntran
488*
489* Do for op(A) = A, A**T, and A**H.
490*
491 trans = transs( itran )
492*
493* Call SLATTB to generate a triangular test matrix.
494*
495 srnamt = 'SLATTB'
496 CALL slattb( imat, uplo, trans, diag, iseed, n, kd,
497 $ ab, ldab, x, work, info )
498*
499*+ TEST 7
500* Solve the system op(A)*x = b
501*
502 srnamt = 'SLATBS'
503 CALL scopy( n, x, 1, b, 1 )
504 CALL slatbs( uplo, trans, diag, 'N', n, kd, ab,
505 $ ldab, b, scale, rwork, info )
506*
507* Check error code from SLATBS.
508*
509 IF( info.NE.0 )
510 $ CALL alaerh( path, 'SLATBS', info, 0,
511 $ uplo // trans // diag // 'N', n, n,
512 $ kd, kd, -1, imat, nfail, nerrs,
513 $ nout )
514*
515 CALL stbt03( uplo, trans, diag, n, kd, 1, ab, ldab,
516 $ scale, rwork, one, b, lda, x, lda,
517 $ work, result( 7 ) )
518*
519*+ TEST 8
520* Solve op(A)*x = b again with NORMIN = 'Y'.
521*
522 CALL scopy( n, x, 1, b, 1 )
523 CALL slatbs( uplo, trans, diag, 'Y', n, kd, ab,
524 $ ldab, b, scale, rwork, info )
525*
526* Check error code from SLATBS.
527*
528 IF( info.NE.0 )
529 $ CALL alaerh( path, 'SLATBS', info, 0,
530 $ uplo // trans // diag // 'Y', n, n,
531 $ kd, kd, -1, imat, nfail, nerrs,
532 $ nout )
533*
534 CALL stbt03( uplo, trans, diag, n, kd, 1, ab, ldab,
535 $ scale, rwork, one, b, lda, x, lda,
536 $ work, result( 8 ) )
537*
538* Print information about the tests that did not pass
539* the threshold.
540*
541 IF( result( 7 ).GE.thresh ) THEN
542 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
543 $ CALL alahd( nout, path )
544 WRITE( nout, fmt = 9997 )'SLATBS', uplo, trans,
545 $ diag, 'N', n, kd, imat, 7, result( 7 )
546 nfail = nfail + 1
547 END IF
548 IF( result( 8 ).GE.thresh ) THEN
549 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
550 $ CALL alahd( nout, path )
551 WRITE( nout, fmt = 9997 )'SLATBS', uplo, trans,
552 $ diag, 'Y', n, kd, imat, 8, result( 8 )
553 nfail = nfail + 1
554 END IF
555 nrun = nrun + 2
556 100 CONTINUE
557 110 CONTINUE
558 120 CONTINUE
559 130 CONTINUE
560 140 CONTINUE
561*
562* Print a summary of the results.
563*
564 CALL alasum( path, nout, nfail, nrun, nerrs )
565*
566 9999 FORMAT( ' UPLO=''', a1, ''', TRANS=''', a1, ''',
567 $ DIAG=''', a1, ''', N=', i5, ', KD=', i5, ', NRHS=', i5,
568 $ ', type ', i2, ', test(', i2, ')=', g12.5 )
569 9998 FORMAT( 1x, a, '( ''', a1, ''', ''', a1, ''', ''', a1, ''',',
570 $ i5, ',', i5, ', ... ), type ', i2, ', test(', i2, ')=',
571 $ g12.5 )
572 9997 FORMAT( 1x, a, '( ''', a1, ''', ''', a1, ''', ''', a1, ''', ''',
573 $ a1, ''',', i5, ',', i5, ', ... ), type ', i2, ', test(',
574 $ i1, ')=', g12.5 )
575 RETURN
576*
577* End of SCHKTB
578*
579 END
subroutine alasum(type, nout, nfail, nrun, nerrs)
ALASUM
Definition alasum.f:73
subroutine slarhs(path, xtype, uplo, trans, m, n, kl, ku, nrhs, a, lda, x, ldx, b, ldb, iseed, info)
SLARHS
Definition slarhs.f:205
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
Definition alaerh.f:147
subroutine alahd(iounit, path)
ALAHD
Definition alahd.f:107
subroutine scopy(n, sx, incx, sy, incy)
SCOPY
Definition scopy.f:82
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
Definition slacpy.f:103
subroutine slaset(uplo, m, n, alpha, beta, a, lda)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition slaset.f:110
subroutine slatbs(uplo, trans, diag, normin, n, kd, ab, ldab, x, scale, cnorm, info)
SLATBS solves a triangular banded system of equations.
Definition slatbs.f:242
subroutine stbcon(norm, uplo, diag, n, kd, ab, ldab, rcond, work, iwork, info)
STBCON
Definition stbcon.f:143
subroutine stbrfs(uplo, trans, diag, n, kd, nrhs, ab, ldab, b, ldb, x, ldx, ferr, berr, work, iwork, info)
STBRFS
Definition stbrfs.f:188
subroutine stbsv(uplo, trans, diag, n, k, a, lda, x, incx)
STBSV
Definition stbsv.f:189
subroutine stbtrs(uplo, trans, diag, n, kd, nrhs, ab, ldab, b, ldb, info)
STBTRS
Definition stbtrs.f:146
subroutine schktb(dotype, nn, nval, nns, nsval, thresh, tsterr, nmax, ab, ainv, b, x, xact, work, rwork, iwork, nout)
SCHKTB
Definition schktb.f:155
subroutine serrtr(path, nunit)
SERRTR
Definition serrtr.f:55
subroutine sget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
SGET04
Definition sget04.f:102
subroutine slattb(imat, uplo, trans, diag, iseed, n, kd, ab, ldab, b, work, info)
SLATTB
Definition slattb.f:135
subroutine stbt02(uplo, trans, diag, n, kd, nrhs, ab, ldab, x, ldx, b, ldb, work, resid)
STBT02
Definition stbt02.f:154
subroutine stbt03(uplo, trans, diag, n, kd, nrhs, ab, ldab, scale, cnorm, tscal, x, ldx, b, ldb, work, resid)
STBT03
Definition stbt03.f:175
subroutine stbt05(uplo, trans, diag, n, kd, nrhs, ab, ldab, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
STBT05
Definition stbt05.f:189
subroutine stbt06(rcond, rcondc, uplo, diag, n, kd, ab, ldab, work, rat)
STBT06
Definition stbt06.f:125