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