LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches
cchksp.f
Go to the documentation of this file.
1*> \brief \b CCHKSP
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 CCHKSP( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
12* NMAX, A, AFAC, AINV, 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 RWORK( * )
24* COMPLEX A( * ), AFAC( * ), AINV( * ), B( * ),
25* \$ WORK( * ), X( * ), XACT( * )
26* ..
27*
28*
29*> \par Purpose:
30* =============
31*>
32*> \verbatim
33*>
34*> CCHKSP tests CSPTRF, -TRI, -TRS, -RFS, and -CON
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 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 maximum value permitted for N, used in dimensioning the
90*> work arrays.
91*> \endverbatim
92*>
93*> \param[out] A
94*> \verbatim
95*> A is COMPLEX array, dimension
96*> (NMAX*(NMAX+1)/2)
97*> \endverbatim
98*>
99*> \param[out] AFAC
100*> \verbatim
101*> AFAC is COMPLEX array, dimension
102*> (NMAX*(NMAX+1)/2)
103*> \endverbatim
104*>
105*> \param[out] AINV
106*> \verbatim
107*> AINV is COMPLEX array, dimension
108*> (NMAX*(NMAX+1)/2)
109*> \endverbatim
110*>
111*> \param[out] B
112*> \verbatim
113*> B is COMPLEX array, dimension (NMAX*NSMAX)
114*> where NSMAX is the largest entry in NSVAL.
115*> \endverbatim
116*>
117*> \param[out] X
118*> \verbatim
119*> X is COMPLEX array, dimension (NMAX*NSMAX)
120*> \endverbatim
121*>
122*> \param[out] XACT
123*> \verbatim
124*> XACT is COMPLEX array, dimension (NMAX*NSMAX)
125*> \endverbatim
126*>
127*> \param[out] WORK
128*> \verbatim
129*> WORK is COMPLEX array, dimension
130*> (NMAX*max(2,NSMAX))
131*> \endverbatim
132*>
133*> \param[out] RWORK
134*> \verbatim
135*> RWORK is REAL array,
136*> dimension (NMAX+2*NSMAX)
137*> \endverbatim
138*>
139*> \param[out] IWORK
140*> \verbatim
141*> IWORK is INTEGER array, dimension (NMAX)
142*> \endverbatim
143*>
144*> \param[in] NOUT
145*> \verbatim
146*> NOUT is INTEGER
147*> The unit number for output.
148*> \endverbatim
149*
150* Authors:
151* ========
152*
153*> \author Univ. of Tennessee
154*> \author Univ. of California Berkeley
155*> \author Univ. of Colorado Denver
156*> \author NAG Ltd.
157*
158*> \ingroup complex_lin
159*
160* =====================================================================
161 SUBROUTINE cchksp( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR,
162 \$ NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK,
163 \$ IWORK, NOUT )
164*
165* -- LAPACK test routine --
166* -- LAPACK is a software package provided by Univ. of Tennessee, --
167* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
168*
169* .. Scalar Arguments ..
170 LOGICAL TSTERR
171 INTEGER NMAX, NN, NNS, NOUT
172 REAL THRESH
173* ..
174* .. Array Arguments ..
175 LOGICAL DOTYPE( * )
176 INTEGER IWORK( * ), NSVAL( * ), NVAL( * )
177 REAL RWORK( * )
178 COMPLEX A( * ), AFAC( * ), AINV( * ), B( * ),
179 \$ work( * ), x( * ), xact( * )
180* ..
181*
182* =====================================================================
183*
184* .. Parameters ..
185 REAL ZERO
186 PARAMETER ( ZERO = 0.0e+0 )
187 INTEGER NTYPES
188 parameter( ntypes = 11 )
189 INTEGER NTESTS
190 parameter( ntests = 8 )
191* ..
192* .. Local Scalars ..
193 LOGICAL TRFCON, ZEROT
194 CHARACTER DIST, PACKIT, TYPE, UPLO, XTYPE
195 CHARACTER*3 PATH
196 INTEGER I, I1, I2, IMAT, IN, INFO, IOFF, IRHS, IUPLO,
197 \$ izero, j, k, kl, ku, lda, mode, n, nerrs,
198 \$ nfail, nimat, npp, nrhs, nrun, nt
199 REAL ANORM, CNDNUM, RCOND, RCONDC
200* ..
201* .. Local Arrays ..
202 CHARACTER UPLOS( 2 )
203 INTEGER ISEED( 4 ), ISEEDY( 4 )
204 REAL RESULT( NTESTS )
205* ..
206* .. External Functions ..
207 LOGICAL LSAME
208 REAL CLANSP, SGET06
209 EXTERNAL lsame, clansp, sget06
210* ..
211* .. External Subroutines ..
212 EXTERNAL alaerh, alahd, alasum, ccopy, cerrsy, cget04,
215 \$ csptri, csptrs
216* ..
217* .. Intrinsic Functions ..
218 INTRINSIC max, min
219* ..
220* .. Scalars in Common ..
221 LOGICAL LERR, OK
222 CHARACTER*32 SRNAMT
223 INTEGER INFOT, NUNIT
224* ..
225* .. Common blocks ..
226 COMMON / infoc / infot, nunit, ok, lerr
227 COMMON / srnamc / srnamt
228* ..
229* .. Data statements ..
230 DATA iseedy / 1988, 1989, 1990, 1991 /
231 DATA uplos / 'U', 'L' /
232* ..
233* .. Executable Statements ..
234*
235* Initialize constants and the random number seed.
236*
237 path( 1: 1 ) = 'Complex precision'
238 path( 2: 3 ) = 'SP'
239 nrun = 0
240 nfail = 0
241 nerrs = 0
242 DO 10 i = 1, 4
243 iseed( i ) = iseedy( i )
244 10 CONTINUE
245*
246* Test the error exits
247*
248 IF( tsterr )
249 \$ CALL cerrsy( path, nout )
250 infot = 0
251*
252* Do for each value of N in NVAL
253*
254 DO 170 in = 1, nn
255 n = nval( in )
256 lda = max( n, 1 )
257 xtype = 'N'
258 nimat = ntypes
259 IF( n.LE.0 )
260 \$ nimat = 1
261*
262 DO 160 imat = 1, nimat
263*
264* Do the tests only if DOTYPE( IMAT ) is true.
265*
266 IF( .NOT.dotype( imat ) )
267 \$ GO TO 160
268*
269* Skip types 3, 4, 5, or 6 if the matrix size is too small.
270*
271 zerot = imat.GE.3 .AND. imat.LE.6
272 IF( zerot .AND. n.LT.imat-2 )
273 \$ GO TO 160
274*
275* Do first for UPLO = 'U', then for UPLO = 'L'
276*
277 DO 150 iuplo = 1, 2
278 uplo = uplos( iuplo )
279 IF( lsame( uplo, 'U' ) ) THEN
280 packit = 'C'
281 ELSE
282 packit = 'R'
283 END IF
284*
285 IF( imat.NE.ntypes ) THEN
286*
287* Set up parameters with CLATB4 and generate a test
288* matrix with CLATMS.
289*
290 CALL clatb4( path, imat, n, n, TYPE, kl, ku, anorm,
291 \$ mode, cndnum, dist )
292*
293 srnamt = 'CLATMS'
294 CALL clatms( n, n, dist, iseed, TYPE, rwork, mode,
295 \$ cndnum, anorm, kl, ku, packit, a, lda,
296 \$ work, info )
297*
298* Check error code from CLATMS.
299*
300 IF( info.NE.0 ) THEN
301 CALL alaerh( path, 'CLATMS', info, 0, uplo, n, n,
302 \$ -1, -1, -1, imat, nfail, nerrs, nout )
303 GO TO 150
304 END IF
305*
306* For types 3-6, zero one or more rows and columns of
307* the matrix to test that INFO is returned correctly.
308*
309 IF( zerot ) THEN
310 IF( imat.EQ.3 ) THEN
311 izero = 1
312 ELSE IF( imat.EQ.4 ) THEN
313 izero = n
314 ELSE
315 izero = n / 2 + 1
316 END IF
317*
318 IF( imat.LT.6 ) THEN
319*
320* Set row and column IZERO to zero.
321*
322 IF( iuplo.EQ.1 ) THEN
323 ioff = ( izero-1 )*izero / 2
324 DO 20 i = 1, izero - 1
325 a( ioff+i ) = zero
326 20 CONTINUE
327 ioff = ioff + izero
328 DO 30 i = izero, n
329 a( ioff ) = zero
330 ioff = ioff + i
331 30 CONTINUE
332 ELSE
333 ioff = izero
334 DO 40 i = 1, izero - 1
335 a( ioff ) = zero
336 ioff = ioff + n - i
337 40 CONTINUE
338 ioff = ioff - izero
339 DO 50 i = izero, n
340 a( ioff+i ) = zero
341 50 CONTINUE
342 END IF
343 ELSE
344 IF( iuplo.EQ.1 ) THEN
345*
346* Set the first IZERO rows and columns to zero.
347*
348 ioff = 0
349 DO 70 j = 1, n
350 i2 = min( j, izero )
351 DO 60 i = 1, i2
352 a( ioff+i ) = zero
353 60 CONTINUE
354 ioff = ioff + j
355 70 CONTINUE
356 ELSE
357*
358* Set the last IZERO rows and columns to zero.
359*
360 ioff = 0
361 DO 90 j = 1, n
362 i1 = max( j, izero )
363 DO 80 i = i1, n
364 a( ioff+i ) = zero
365 80 CONTINUE
366 ioff = ioff + n - j
367 90 CONTINUE
368 END IF
369 END IF
370 ELSE
371 izero = 0
372 END IF
373 ELSE
374*
375* Use a special block diagonal matrix to test alternate
376* code for the 2 x 2 blocks.
377*
378 CALL clatsp( uplo, n, a, iseed )
379 END IF
380*
381* Compute the L*D*L' or U*D*U' factorization of the matrix.
382*
383 npp = n*( n+1 ) / 2
384 CALL ccopy( npp, a, 1, afac, 1 )
385 srnamt = 'CSPTRF'
386 CALL csptrf( uplo, n, afac, iwork, info )
387*
388* Adjust the expected value of INFO to account for
389* pivoting.
390*
391 k = izero
392 IF( k.GT.0 ) THEN
393 100 CONTINUE
394 IF( iwork( k ).LT.0 ) THEN
395 IF( iwork( k ).NE.-k ) THEN
396 k = -iwork( k )
397 GO TO 100
398 END IF
399 ELSE IF( iwork( k ).NE.k ) THEN
400 k = iwork( k )
401 GO TO 100
402 END IF
403 END IF
404*
405* Check error code from CSPTRF.
406*
407 IF( info.NE.k )
408 \$ CALL alaerh( path, 'CSPTRF', info, k, uplo, n, n, -1,
409 \$ -1, -1, imat, nfail, nerrs, nout )
410 IF( info.NE.0 ) THEN
411 trfcon = .true.
412 ELSE
413 trfcon = .false.
414 END IF
415*
416*+ TEST 1
417* Reconstruct matrix from factors and compute residual.
418*
419 CALL cspt01( uplo, n, a, afac, iwork, ainv, lda, rwork,
420 \$ result( 1 ) )
421 nt = 1
422*
423*+ TEST 2
424* Form the inverse and compute the residual.
425*
426 IF( .NOT.trfcon ) THEN
427 CALL ccopy( npp, afac, 1, ainv, 1 )
428 srnamt = 'CSPTRI'
429 CALL csptri( uplo, n, ainv, iwork, work, info )
430*
431* Check error code from CSPTRI.
432*
433 IF( info.NE.0 )
434 \$ CALL alaerh( path, 'CSPTRI', info, 0, uplo, n, n,
435 \$ -1, -1, -1, imat, nfail, nerrs, nout )
436*
437 CALL cspt03( uplo, n, a, ainv, work, lda, rwork,
438 \$ rcondc, result( 2 ) )
439 nt = 2
440 END IF
441*
442* Print information about the tests that did not pass
443* the threshold.
444*
445 DO 110 k = 1, nt
446 IF( result( k ).GE.thresh ) THEN
447 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
448 \$ CALL alahd( nout, path )
449 WRITE( nout, fmt = 9999 )uplo, n, imat, k,
450 \$ result( k )
451 nfail = nfail + 1
452 END IF
453 110 CONTINUE
454 nrun = nrun + nt
455*
456* Do only the condition estimate if INFO is not 0.
457*
458 IF( trfcon ) THEN
459 rcondc = zero
460 GO TO 140
461 END IF
462*
463 DO 130 irhs = 1, nns
464 nrhs = nsval( irhs )
465*
466*+ TEST 3
467* Solve and compute residual for A * X = B.
468*
469 srnamt = 'CLARHS'
470 CALL clarhs( path, xtype, uplo, ' ', n, n, kl, ku,
471 \$ nrhs, a, lda, xact, lda, b, lda, iseed,
472 \$ info )
473 CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
474*
475 srnamt = 'CSPTRS'
476 CALL csptrs( uplo, n, nrhs, afac, iwork, x, lda,
477 \$ info )
478*
479* Check error code from CSPTRS.
480*
481 IF( info.NE.0 )
482 \$ CALL alaerh( path, 'CSPTRS', info, 0, uplo, n, n,
483 \$ -1, -1, nrhs, imat, nfail, nerrs,
484 \$ nout )
485*
486 CALL clacpy( 'Full', n, nrhs, b, lda, work, lda )
487 CALL cspt02( uplo, n, nrhs, a, x, lda, work, lda,
488 \$ rwork, result( 3 ) )
489*
490*+ TEST 4
491* Check solution from generated exact solution.
492*
493 CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
494 \$ result( 4 ) )
495*
496*+ TESTS 5, 6, and 7
497* Use iterative refinement to improve the solution.
498*
499 srnamt = 'CSPRFS'
500 CALL csprfs( uplo, n, nrhs, a, afac, iwork, b, lda, x,
501 \$ lda, rwork, rwork( nrhs+1 ), work,
502 \$ rwork( 2*nrhs+1 ), info )
503*
504* Check error code from CSPRFS.
505*
506 IF( info.NE.0 )
507 \$ CALL alaerh( path, 'CSPRFS', info, 0, uplo, n, n,
508 \$ -1, -1, nrhs, imat, nfail, nerrs,
509 \$ nout )
510*
511 CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
512 \$ result( 5 ) )
513 CALL cppt05( uplo, n, nrhs, a, b, lda, x, lda, xact,
514 \$ lda, rwork, rwork( nrhs+1 ),
515 \$ result( 6 ) )
516*
517* Print information about the tests that did not pass
518* the threshold.
519*
520 DO 120 k = 3, 7
521 IF( result( k ).GE.thresh ) THEN
522 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
523 \$ CALL alahd( nout, path )
524 WRITE( nout, fmt = 9998 )uplo, n, nrhs, imat,
525 \$ k, result( k )
526 nfail = nfail + 1
527 END IF
528 120 CONTINUE
529 nrun = nrun + 5
530 130 CONTINUE
531*
532*+ TEST 8
533* Get an estimate of RCOND = 1/CNDNUM.
534*
535 140 CONTINUE
536 anorm = clansp( '1', uplo, n, a, rwork )
537 srnamt = 'CSPCON'
538 CALL cspcon( uplo, n, afac, iwork, anorm, rcond, work,
539 \$ info )
540*
541* Check error code from CSPCON.
542*
543 IF( info.NE.0 )
544 \$ CALL alaerh( path, 'CSPCON', info, 0, uplo, n, n, -1,
545 \$ -1, -1, imat, nfail, nerrs, nout )
546*
547 result( 8 ) = sget06( rcond, rcondc )
548*
549* Print the test ratio if it is .GE. THRESH.
550*
551 IF( result( 8 ).GE.thresh ) THEN
552 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
553 \$ CALL alahd( nout, path )
554 WRITE( nout, fmt = 9999 )uplo, n, imat, 8,
555 \$ result( 8 )
556 nfail = nfail + 1
557 END IF
558 nrun = nrun + 1
559 150 CONTINUE
560 160 CONTINUE
561 170 CONTINUE
562*
563* Print a summary of the results.
564*
565 CALL alasum( path, nout, nfail, nrun, nerrs )
566*
567 9999 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ', type ', i2, ', test ',
568 \$ i2, ', ratio =', g12.5 )
569 9998 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ', NRHS=', i3, ', type ',
570 \$ i2, ', test(', i2, ') =', g12.5 )
571 RETURN
572*
573* End of CCHKSP
574*
575 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 cchksp(dotype, nn, nval, nns, nsval, thresh, tsterr, nmax, a, afac, ainv, b, x, xact, work, rwork, iwork, nout)
CCHKSP
Definition cchksp.f:164
subroutine cerrsy(path, nunit)
CERRSY
Definition cerrsy.f:55
subroutine cget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
CGET04
Definition cget04.f:102
subroutine clatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
CLATB4
Definition clatb4.f:121
subroutine clatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
CLATMS
Definition clatms.f:332
subroutine clatsp(uplo, n, x, iseed)
CLATSP
Definition clatsp.f:84
subroutine cppt05(uplo, n, nrhs, ap, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
CPPT05
Definition cppt05.f:157
subroutine cspt01(uplo, n, a, afac, ipiv, c, ldc, rwork, resid)
CSPT01
Definition cspt01.f:112
subroutine cspt02(uplo, n, nrhs, a, x, ldx, b, ldb, rwork, resid)
CSPT02
Definition cspt02.f:123
subroutine cspt03(uplo, n, a, ainv, work, ldw, rwork, rcond, resid)
CSPT03
Definition cspt03.f:110
subroutine ccopy(n, cx, incx, cy, incy)
CCOPY
Definition ccopy.f:81
subroutine cspcon(uplo, n, ap, ipiv, anorm, rcond, work, info)
CSPCON
Definition cspcon.f:118
subroutine csprfs(uplo, n, nrhs, ap, afp, ipiv, b, ldb, x, ldx, ferr, berr, work, rwork, info)
CSPRFS
Definition csprfs.f:180
subroutine csptrf(uplo, n, ap, ipiv, info)
CSPTRF
Definition csptrf.f:158
subroutine csptri(uplo, n, ap, ipiv, work, info)
CSPTRI
Definition csptri.f:109
subroutine csptrs(uplo, n, nrhs, ap, ipiv, b, ldb, info)
CSPTRS
Definition csptrs.f:115
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