LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
schktp.f
Go to the documentation of this file.
1*> \brief \b SCHKTP
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 SCHKTP( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
12* NMAX, AP, AINVP, B, X, XACT, WORK, RWORK,
13* IWORK, 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 AINVP( * ), AP( * ), B( * ), RWORK( * ),
24* $ WORK( * ), X( * ), XACT( * )
25* ..
26*
27*
28*> \par Purpose:
29* =============
30*>
31*> \verbatim
32*>
33*> SCHKTP tests STPTRI, -TRS, -RFS, and -CON, and SLATPS
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. NMAX >= the
89*> maximum value of N in NVAL.
90*> \endverbatim
91*>
92*> \param[out] AP
93*> \verbatim
94*> AP is REAL array, dimension
95*> (NMAX*(NMAX+1)/2)
96*> \endverbatim
97*>
98*> \param[out] AINVP
99*> \verbatim
100*> AINVP is REAL array, dimension
101*> (NMAX*(NMAX+1)/2)
102*> \endverbatim
103*>
104*> \param[out] B
105*> \verbatim
106*> B is REAL array, dimension (NMAX*NSMAX)
107*> where NSMAX is the largest entry in NSVAL.
108*> \endverbatim
109*>
110*> \param[out] X
111*> \verbatim
112*> X is REAL array, dimension (NMAX*NSMAX)
113*> \endverbatim
114*>
115*> \param[out] XACT
116*> \verbatim
117*> XACT is REAL array, dimension (NMAX*NSMAX)
118*> \endverbatim
119*>
120*> \param[out] WORK
121*> \verbatim
122*> WORK is REAL array, dimension
123*> (NMAX*max(3,NSMAX))
124*> \endverbatim
125*>
126*> \param[out] IWORK
127*> \verbatim
128*> IWORK is INTEGER array, dimension (NMAX)
129*> \endverbatim
130*>
131*> \param[out] RWORK
132*> \verbatim
133*> RWORK is REAL array, dimension
134*> (max(NMAX,2*NSMAX))
135*> \endverbatim
136*>
137*> \param[in] NOUT
138*> \verbatim
139*> NOUT is INTEGER
140*> The unit number for output.
141*> \endverbatim
142*
143* Authors:
144* ========
145*
146*> \author Univ. of Tennessee
147*> \author Univ. of California Berkeley
148*> \author Univ. of Colorado Denver
149*> \author NAG Ltd.
150*
151*> \ingroup single_lin
152*
153* =====================================================================
154 SUBROUTINE schktp( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
155 $ NMAX, AP, AINVP, B, X, XACT, WORK, RWORK,
156 $ IWORK, NOUT )
157*
158* -- LAPACK test routine --
159* -- LAPACK is a software package provided by Univ. of Tennessee, --
160* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
161*
162* .. Scalar Arguments ..
163 LOGICAL TSTERR
164 INTEGER NMAX, NN, NNS, NOUT
165 REAL THRESH
166* ..
167* .. Array Arguments ..
168 LOGICAL DOTYPE( * )
169 INTEGER IWORK( * ), NSVAL( * ), NVAL( * )
170 REAL AINVP( * ), AP( * ), B( * ), RWORK( * ),
171 $ work( * ), x( * ), xact( * )
172* ..
173*
174* =====================================================================
175*
176* .. Parameters ..
177 INTEGER NTYPE1, NTYPES
178 PARAMETER ( NTYPE1 = 10, ntypes = 18 )
179 INTEGER NTESTS
180 parameter( ntests = 9 )
181 INTEGER NTRAN
182 parameter( ntran = 3 )
183 REAL ONE, ZERO
184 parameter( one = 1.0e+0, zero = 0.0e+0 )
185* ..
186* .. Local Scalars ..
187 CHARACTER DIAG, NORM, TRANS, UPLO, XTYPE
188 CHARACTER*3 PATH
189 INTEGER I, IDIAG, IMAT, IN, INFO, IRHS, ITRAN, IUPLO,
190 $ k, lap, lda, n, nerrs, nfail, nrhs, nrun
191 REAL AINVNM, ANORM, RCOND, RCONDC, RCONDI, RCONDO,
192 $ SCALE
193* ..
194* .. Local Arrays ..
195 CHARACTER TRANSS( NTRAN ), UPLOS( 2 )
196 INTEGER ISEED( 4 ), ISEEDY( 4 )
197 REAL RESULT( NTESTS )
198* ..
199* .. External Functions ..
200 LOGICAL LSAME
201 REAL SLANTP
202 EXTERNAL lsame, slantp
203* ..
204* .. External Subroutines ..
205 EXTERNAL alaerh, alahd, alasum, scopy, serrtr, sget04,
208 $ stptrs
209* ..
210* .. Scalars in Common ..
211 LOGICAL LERR, OK
212 CHARACTER*32 SRNAMT
213 INTEGER INFOT, IOUNIT
214* ..
215* .. Common blocks ..
216 COMMON / infoc / infot, iounit, ok, lerr
217 COMMON / srnamc / srnamt
218* ..
219* .. Intrinsic Functions ..
220 INTRINSIC max
221* ..
222* .. Data statements ..
223 DATA iseedy / 1988, 1989, 1990, 1991 /
224 DATA uplos / 'U', 'L' / , transs / 'N', 'T', 'C' /
225* ..
226* .. Executable Statements ..
227*
228* Initialize constants and the random number seed.
229*
230 path( 1: 1 ) = 'Single precision'
231 path( 2: 3 ) = 'TP'
232 nrun = 0
233 nfail = 0
234 nerrs = 0
235 DO 10 i = 1, 4
236 iseed( i ) = iseedy( i )
237 10 CONTINUE
238*
239* Test the error exits
240*
241 IF( tsterr )
242 $ CALL serrtr( path, nout )
243 infot = 0
244*
245 DO 110 in = 1, nn
246*
247* Do for each value of N in NVAL
248*
249 n = nval( in )
250 lda = max( 1, n )
251 lap = lda*( lda+1 ) / 2
252 xtype = 'N'
253*
254 DO 70 imat = 1, ntype1
255*
256* Do the tests only if DOTYPE( IMAT ) is true.
257*
258 IF( .NOT.dotype( imat ) )
259 $ GO TO 70
260*
261 DO 60 iuplo = 1, 2
262*
263* Do first for UPLO = 'U', then for UPLO = 'L'
264*
265 uplo = uplos( iuplo )
266*
267* Call SLATTP to generate a triangular test matrix.
268*
269 srnamt = 'SLATTP'
270 CALL slattp( imat, uplo, 'No transpose', diag, iseed, n,
271 $ ap, x, work, info )
272*
273* Set IDIAG = 1 for non-unit matrices, 2 for unit.
274*
275 IF( lsame( diag, 'N' ) ) THEN
276 idiag = 1
277 ELSE
278 idiag = 2
279 END IF
280*
281*+ TEST 1
282* Form the inverse of A.
283*
284 IF( n.GT.0 )
285 $ CALL scopy( lap, ap, 1, ainvp, 1 )
286 srnamt = 'STPTRI'
287 CALL stptri( uplo, diag, n, ainvp, info )
288*
289* Check error code from STPTRI.
290*
291 IF( info.NE.0 )
292 $ CALL alaerh( path, 'STPTRI', info, 0, uplo // diag, n,
293 $ n, -1, -1, -1, imat, nfail, nerrs, nout )
294*
295* Compute the infinity-norm condition number of A.
296*
297 anorm = slantp( 'I', uplo, diag, n, ap, rwork )
298 ainvnm = slantp( 'I', uplo, diag, n, ainvp, rwork )
299 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
300 rcondi = one
301 ELSE
302 rcondi = ( one / anorm ) / ainvnm
303 END IF
304*
305* Compute the residual for the triangular matrix times its
306* inverse. Also compute the 1-norm condition number of A.
307*
308 CALL stpt01( uplo, diag, n, ap, ainvp, rcondo, rwork,
309 $ result( 1 ) )
310*
311* Print the test ratio if it is .GE. THRESH.
312*
313 IF( result( 1 ).GE.thresh ) THEN
314 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
315 $ CALL alahd( nout, path )
316 WRITE( nout, fmt = 9999 )uplo, diag, n, imat, 1,
317 $ result( 1 )
318 nfail = nfail + 1
319 END IF
320 nrun = nrun + 1
321*
322 DO 40 irhs = 1, nns
323 nrhs = nsval( irhs )
324 xtype = 'N'
325*
326 DO 30 itran = 1, ntran
327*
328* Do for op(A) = A, A**T, or A**H.
329*
330 trans = transs( itran )
331 IF( itran.EQ.1 ) THEN
332 norm = 'O'
333 rcondc = rcondo
334 ELSE
335 norm = 'I'
336 rcondc = rcondi
337 END IF
338*
339*+ TEST 2
340* Solve and compute residual for op(A)*x = b.
341*
342 srnamt = 'SLARHS'
343 CALL slarhs( path, xtype, uplo, trans, n, n, 0,
344 $ idiag, nrhs, ap, lap, xact, lda, b,
345 $ lda, iseed, info )
346 xtype = 'C'
347 CALL slacpy( 'Full', n, nrhs, b, lda, x, lda )
348*
349 srnamt = 'STPTRS'
350 CALL stptrs( uplo, trans, diag, n, nrhs, ap, x,
351 $ lda, info )
352*
353* Check error code from STPTRS.
354*
355 IF( info.NE.0 )
356 $ CALL alaerh( path, 'STPTRS', info, 0,
357 $ uplo // trans // diag, n, n, -1,
358 $ -1, -1, imat, nfail, nerrs, nout )
359*
360 CALL stpt02( uplo, trans, diag, n, nrhs, ap, x,
361 $ lda, b, lda, work, result( 2 ) )
362*
363*+ TEST 3
364* Check solution from generated exact solution.
365*
366 CALL sget04( n, nrhs, x, lda, xact, lda, rcondc,
367 $ result( 3 ) )
368*
369*+ TESTS 4, 5, and 6
370* Use iterative refinement to improve the solution and
371* compute error bounds.
372*
373 srnamt = 'STPRFS'
374 CALL stprfs( uplo, trans, diag, n, nrhs, ap, b,
375 $ lda, x, lda, rwork, rwork( nrhs+1 ),
376 $ work, iwork, info )
377*
378* Check error code from STPRFS.
379*
380 IF( info.NE.0 )
381 $ CALL alaerh( path, 'STPRFS', info, 0,
382 $ uplo // trans // diag, n, n, -1,
383 $ -1, nrhs, imat, nfail, nerrs,
384 $ nout )
385*
386 CALL sget04( n, nrhs, x, lda, xact, lda, rcondc,
387 $ result( 4 ) )
388 CALL stpt05( uplo, trans, diag, n, nrhs, ap, b,
389 $ lda, x, lda, xact, lda, rwork,
390 $ rwork( nrhs+1 ), result( 5 ) )
391*
392* Print information about the tests that did not pass
393* the threshold.
394*
395 DO 20 k = 2, 6
396 IF( result( k ).GE.thresh ) THEN
397 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
398 $ CALL alahd( nout, path )
399 WRITE( nout, fmt = 9998 )uplo, trans, diag,
400 $ n, nrhs, imat, k, result( k )
401 nfail = nfail + 1
402 END IF
403 20 CONTINUE
404 nrun = nrun + 5
405 30 CONTINUE
406 40 CONTINUE
407*
408*+ TEST 7
409* Get an estimate of RCOND = 1/CNDNUM.
410*
411 DO 50 itran = 1, 2
412 IF( itran.EQ.1 ) THEN
413 norm = 'O'
414 rcondc = rcondo
415 ELSE
416 norm = 'I'
417 rcondc = rcondi
418 END IF
419*
420 srnamt = 'STPCON'
421 CALL stpcon( norm, uplo, diag, n, ap, rcond, work,
422 $ iwork, info )
423*
424* Check error code from STPCON.
425*
426 IF( info.NE.0 )
427 $ CALL alaerh( path, 'STPCON', info, 0,
428 $ norm // uplo // diag, n, n, -1, -1,
429 $ -1, imat, nfail, nerrs, nout )
430*
431 CALL stpt06( rcond, rcondc, uplo, diag, n, ap, rwork,
432 $ result( 7 ) )
433*
434* Print the test ratio if it is .GE. THRESH.
435*
436 IF( result( 7 ).GE.thresh ) THEN
437 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
438 $ CALL alahd( nout, path )
439 WRITE( nout, fmt = 9997 ) 'STPCON', norm, uplo,
440 $ diag, n, imat, 7, result( 7 )
441 nfail = nfail + 1
442 END IF
443 nrun = nrun + 1
444 50 CONTINUE
445 60 CONTINUE
446 70 CONTINUE
447*
448* Use pathological test matrices to test SLATPS.
449*
450 DO 100 imat = ntype1 + 1, ntypes
451*
452* Do the tests only if DOTYPE( IMAT ) is true.
453*
454 IF( .NOT.dotype( imat ) )
455 $ GO TO 100
456*
457 DO 90 iuplo = 1, 2
458*
459* Do first for UPLO = 'U', then for UPLO = 'L'
460*
461 uplo = uplos( iuplo )
462 DO 80 itran = 1, ntran
463*
464* Do for op(A) = A, A**T, or A**H.
465*
466 trans = transs( itran )
467*
468* Call SLATTP to generate a triangular test matrix.
469*
470 srnamt = 'SLATTP'
471 CALL slattp( imat, uplo, trans, diag, iseed, n, ap, x,
472 $ work, info )
473*
474*+ TEST 8
475* Solve the system op(A)*x = b.
476*
477 srnamt = 'SLATPS'
478 CALL scopy( n, x, 1, b, 1 )
479 CALL slatps( uplo, trans, diag, 'N', n, ap, b, scale,
480 $ rwork, info )
481*
482* Check error code from SLATPS.
483*
484 IF( info.NE.0 )
485 $ CALL alaerh( path, 'SLATPS', info, 0,
486 $ uplo // trans // diag // 'N', n, n,
487 $ -1, -1, -1, imat, nfail, nerrs, nout )
488*
489 CALL stpt03( uplo, trans, diag, n, 1, ap, scale,
490 $ rwork, one, b, lda, x, lda, work,
491 $ result( 8 ) )
492*
493*+ TEST 9
494* Solve op(A)*x = b again with NORMIN = 'Y'.
495*
496 CALL scopy( n, x, 1, b( n+1 ), 1 )
497 CALL slatps( uplo, trans, diag, 'Y', n, ap, b( n+1 ),
498 $ scale, rwork, info )
499*
500* Check error code from SLATPS.
501*
502 IF( info.NE.0 )
503 $ CALL alaerh( path, 'SLATPS', info, 0,
504 $ uplo // trans // diag // 'Y', n, n,
505 $ -1, -1, -1, imat, nfail, nerrs, nout )
506*
507 CALL stpt03( uplo, trans, diag, n, 1, ap, scale,
508 $ rwork, one, b( n+1 ), lda, x, lda, work,
509 $ result( 9 ) )
510*
511* Print information about the tests that did not pass
512* the threshold.
513*
514 IF( result( 8 ).GE.thresh ) THEN
515 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
516 $ CALL alahd( nout, path )
517 WRITE( nout, fmt = 9996 )'SLATPS', uplo, trans,
518 $ diag, 'N', n, imat, 8, result( 8 )
519 nfail = nfail + 1
520 END IF
521 IF( result( 9 ).GE.thresh ) THEN
522 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
523 $ CALL alahd( nout, path )
524 WRITE( nout, fmt = 9996 )'SLATPS', uplo, trans,
525 $ diag, 'Y', n, imat, 9, result( 9 )
526 nfail = nfail + 1
527 END IF
528 nrun = nrun + 2
529 80 CONTINUE
530 90 CONTINUE
531 100 CONTINUE
532 110 CONTINUE
533*
534* Print a summary of the results.
535*
536 CALL alasum( path, nout, nfail, nrun, nerrs )
537*
538 9999 FORMAT( ' UPLO=''', a1, ''', DIAG=''', a1, ''', N=', i5,
539 $ ', type ', i2, ', test(', i2, ')= ', g12.5 )
540 9998 FORMAT( ' UPLO=''', a1, ''', TRANS=''', a1, ''', DIAG=''', a1,
541 $ ''', N=', i5, ''', NRHS=', i5, ', type ', i2, ', test(',
542 $ i2, ')= ', g12.5 )
543 9997 FORMAT( 1x, a, '( ''', a1, ''', ''', a1, ''', ''', a1, ''',',
544 $ i5, ', ... ), type ', i2, ', test(', i2, ')=', g12.5 )
545 9996 FORMAT( 1x, a, '( ''', a1, ''', ''', a1, ''', ''', a1, ''', ''',
546 $ a1, ''',', i5, ', ... ), type ', i2, ', test(', i2, ')=',
547 $ g12.5 )
548 RETURN
549*
550* End of SCHKTP
551*
552 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 slatps(uplo, trans, diag, normin, n, ap, x, scale, cnorm, info)
SLATPS solves a triangular system of equations with the matrix held in packed storage.
Definition slatps.f:229
subroutine stpcon(norm, uplo, diag, n, ap, rcond, work, iwork, info)
STPCON
Definition stpcon.f:130
subroutine stprfs(uplo, trans, diag, n, nrhs, ap, b, ldb, x, ldx, ferr, berr, work, iwork, info)
STPRFS
Definition stprfs.f:175
subroutine stptri(uplo, diag, n, ap, info)
STPTRI
Definition stptri.f:117
subroutine stptrs(uplo, trans, diag, n, nrhs, ap, b, ldb, info)
STPTRS
Definition stptrs.f:130
subroutine schktp(dotype, nn, nval, nns, nsval, thresh, tsterr, nmax, ap, ainvp, b, x, xact, work, rwork, iwork, nout)
SCHKTP
Definition schktp.f:157
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 slattp(imat, uplo, trans, diag, iseed, n, a, b, work, info)
SLATTP
Definition slattp.f:125
subroutine stpt01(uplo, diag, n, ap, ainvp, rcond, work, resid)
STPT01
Definition stpt01.f:108
subroutine stpt02(uplo, trans, diag, n, nrhs, ap, x, ldx, b, ldb, work, resid)
STPT02
Definition stpt02.f:142
subroutine stpt03(uplo, trans, diag, n, nrhs, ap, scale, cnorm, tscal, x, ldx, b, ldb, work, resid)
STPT03
Definition stpt03.f:161
subroutine stpt05(uplo, trans, diag, n, nrhs, ap, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
STPT05
Definition stpt05.f:174
subroutine stpt06(rcond, rcondc, uplo, diag, n, ap, work, rat)
STPT06
Definition stpt06.f:111