LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
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 *> maximumm 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 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 alasum(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASUM
Definition: alasum.f:73
subroutine alahd(IOUNIT, PATH)
ALAHD
Definition: alahd.f:107
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:147
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 stptrs(UPLO, TRANS, DIAG, N, NRHS, AP, B, LDB, INFO)
STPTRS
Definition: stptrs.f:130
subroutine stptri(UPLO, DIAG, N, AP, INFO)
STPTRI
Definition: stptri.f:117
subroutine stprfs(UPLO, TRANS, DIAG, N, NRHS, AP, B, LDB, X, LDX, FERR, BERR, WORK, IWORK, INFO)
STPRFS
Definition: stprfs.f:175
subroutine stpcon(NORM, UPLO, DIAG, N, AP, RCOND, WORK, IWORK, INFO)
STPCON
Definition: stpcon.f:130
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:82
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 schktp(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, AP, AINVP, B, X, XACT, WORK, RWORK, IWORK, NOUT)
SCHKTP
Definition: schktp.f:157
subroutine stpt02(UPLO, TRANS, DIAG, N, NRHS, AP, X, LDX, B, LDB, WORK, RESID)
STPT02
Definition: stpt02.f:142
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 stpt06(RCOND, RCONDC, UPLO, DIAG, N, AP, WORK, RAT)
STPT06
Definition: stpt06.f:111
subroutine stpt03(UPLO, TRANS, DIAG, N, NRHS, AP, SCALE, CNORM, TSCAL, X, LDX, B, LDB, WORK, RESID)
STPT03
Definition: stpt03.f:161
subroutine stpt01(UPLO, DIAG, N, AP, AINVP, RCOND, WORK, RESID)
STPT01
Definition: stpt01.f:108
subroutine stpt05(UPLO, TRANS, DIAG, N, NRHS, AP, B, LDB, X, LDX, XACT, LDXACT, FERR, BERR, RESLTS)
STPT05
Definition: stpt05.f:174