LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches
cchkhp.f
Go to the documentation of this file.
1*> \brief \b CCHKHP
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 CCHKHP( 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*> CCHKHP tests CHPTRF, -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 cchkhp( 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 = 10 )
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 CLANHP, SGET06
209 EXTERNAL lsame, clanhp, sget06
210* ..
211* .. External Subroutines ..
212 EXTERNAL alaerh, alahd, alasum, ccopy, cerrsy, cget04,
215 \$ cppt03, cppt05
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 ) = 'HP'
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 izero = 0
263 DO 160 imat = 1, nimat
264*
265* Do the tests only if DOTYPE( IMAT ) is true.
266*
267 IF( .NOT.dotype( imat ) )
268 \$ GO TO 160
269*
270* Skip types 3, 4, 5, or 6 if the matrix size is too small.
271*
272 zerot = imat.GE.3 .AND. imat.LE.6
273 IF( zerot .AND. n.LT.imat-2 )
274 \$ GO TO 160
275*
276* Do first for UPLO = 'U', then for UPLO = 'L'
277*
278 DO 150 iuplo = 1, 2
279 uplo = uplos( iuplo )
280 IF( lsame( uplo, 'U' ) ) THEN
281 packit = 'C'
282 ELSE
283 packit = 'R'
284 END IF
285*
286* Set up parameters with CLATB4 and generate a test matrix
287* with CLATMS.
288*
289 CALL clatb4( path, imat, n, n, TYPE, kl, ku, anorm, mode,
290 \$ cndnum, dist )
291*
292 srnamt = 'CLATMS'
293 CALL clatms( n, n, dist, iseed, TYPE, rwork, mode,
294 \$ cndnum, anorm, kl, ku, packit, a, lda, work,
295 \$ info )
296*
297* Check error code from CLATMS.
298*
299 IF( info.NE.0 ) THEN
300 CALL alaerh( path, 'CLATMS', info, 0, uplo, n, n, -1,
301 \$ -1, -1, imat, nfail, nerrs, nout )
302 GO TO 150
303 END IF
304*
305* For types 3-6, zero one or more rows and columns of
306* the matrix to test that INFO is returned correctly.
307*
308 IF( zerot ) THEN
309 IF( imat.EQ.3 ) THEN
310 izero = 1
311 ELSE IF( imat.EQ.4 ) THEN
312 izero = n
313 ELSE
314 izero = n / 2 + 1
315 END IF
316*
317 IF( imat.LT.6 ) THEN
318*
319* Set row and column IZERO to zero.
320*
321 IF( iuplo.EQ.1 ) THEN
322 ioff = ( izero-1 )*izero / 2
323 DO 20 i = 1, izero - 1
324 a( ioff+i ) = zero
325 20 CONTINUE
326 ioff = ioff + izero
327 DO 30 i = izero, n
328 a( ioff ) = zero
329 ioff = ioff + i
330 30 CONTINUE
331 ELSE
332 ioff = izero
333 DO 40 i = 1, izero - 1
334 a( ioff ) = zero
335 ioff = ioff + n - i
336 40 CONTINUE
337 ioff = ioff - izero
338 DO 50 i = izero, n
339 a( ioff+i ) = zero
340 50 CONTINUE
341 END IF
342 ELSE
343 ioff = 0
344 IF( iuplo.EQ.1 ) THEN
345*
346* Set the first IZERO rows and columns to zero.
347*
348 DO 70 j = 1, n
349 i2 = min( j, izero )
350 DO 60 i = 1, i2
351 a( ioff+i ) = zero
352 60 CONTINUE
353 ioff = ioff + j
354 70 CONTINUE
355 ELSE
356*
357* Set the last IZERO rows and columns to zero.
358*
359 DO 90 j = 1, n
360 i1 = max( j, izero )
361 DO 80 i = i1, n
362 a( ioff+i ) = zero
363 80 CONTINUE
364 ioff = ioff + n - j
365 90 CONTINUE
366 END IF
367 END IF
368 ELSE
369 izero = 0
370 END IF
371*
372* Set the imaginary part of the diagonals.
373*
374 IF( iuplo.EQ.1 ) THEN
375 CALL claipd( n, a, 2, 1 )
376 ELSE
377 CALL claipd( n, a, n, -1 )
378 END IF
379*
380* Compute the L*D*L' or U*D*U' factorization of the matrix.
381*
382 npp = n*( n+1 ) / 2
383 CALL ccopy( npp, a, 1, afac, 1 )
384 srnamt = 'CHPTRF'
385 CALL chptrf( uplo, n, afac, iwork, info )
386*
387* Adjust the expected value of INFO to account for
388* pivoting.
389*
390 k = izero
391 IF( k.GT.0 ) THEN
392 100 CONTINUE
393 IF( iwork( k ).LT.0 ) THEN
394 IF( iwork( k ).NE.-k ) THEN
395 k = -iwork( k )
396 GO TO 100
397 END IF
398 ELSE IF( iwork( k ).NE.k ) THEN
399 k = iwork( k )
400 GO TO 100
401 END IF
402 END IF
403*
404* Check error code from CHPTRF.
405*
406 IF( info.NE.k )
407 \$ CALL alaerh( path, 'CHPTRF', info, k, uplo, n, n, -1,
408 \$ -1, -1, imat, nfail, nerrs, nout )
409 IF( info.NE.0 ) THEN
410 trfcon = .true.
411 ELSE
412 trfcon = .false.
413 END IF
414*
415*+ TEST 1
416* Reconstruct matrix from factors and compute residual.
417*
418 CALL chpt01( uplo, n, a, afac, iwork, ainv, lda, rwork,
419 \$ result( 1 ) )
420 nt = 1
421*
422*+ TEST 2
423* Form the inverse and compute the residual.
424*
425 IF( .NOT.trfcon ) THEN
426 CALL ccopy( npp, afac, 1, ainv, 1 )
427 srnamt = 'CHPTRI'
428 CALL chptri( uplo, n, ainv, iwork, work, info )
429*
430* Check error code from CHPTRI.
431*
432 IF( info.NE.0 )
433 \$ CALL alaerh( path, 'CHPTRI', info, 0, uplo, n, n,
434 \$ -1, -1, -1, imat, nfail, nerrs, nout )
435*
436 CALL cppt03( uplo, n, a, ainv, work, lda, rwork,
437 \$ rcondc, result( 2 ) )
438 nt = 2
439 END IF
440*
441* Print information about the tests that did not pass
442* the threshold.
443*
444 DO 110 k = 1, nt
445 IF( result( k ).GE.thresh ) THEN
446 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
447 \$ CALL alahd( nout, path )
448 WRITE( nout, fmt = 9999 )uplo, n, imat, k,
449 \$ result( k )
450 nfail = nfail + 1
451 END IF
452 110 CONTINUE
453 nrun = nrun + nt
454*
455* Do only the condition estimate if INFO is not 0.
456*
457 IF( trfcon ) THEN
458 rcondc = zero
459 GO TO 140
460 END IF
461*
462 DO 130 irhs = 1, nns
463 nrhs = nsval( irhs )
464*
465*+ TEST 3
466* Solve and compute residual for A * X = B.
467*
468 srnamt = 'CLARHS'
469 CALL clarhs( path, xtype, uplo, ' ', n, n, kl, ku,
470 \$ nrhs, a, lda, xact, lda, b, lda, iseed,
471 \$ info )
472 xtype = 'C'
473 CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
474*
475 srnamt = 'CHPTRS'
476 CALL chptrs( uplo, n, nrhs, afac, iwork, x, lda,
477 \$ info )
478*
479* Check error code from CHPTRS.
480*
481 IF( info.NE.0 )
482 \$ CALL alaerh( path, 'CHPTRS', 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 cppt02( 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 = 'CHPRFS'
500 CALL chprfs( 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 CHPRFS.
505*
506 IF( info.NE.0 )
507 \$ CALL alaerh( path, 'CHPRFS', 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 = clanhp( '1', uplo, n, a, rwork )
537 srnamt = 'CHPCON'
538 CALL chpcon( uplo, n, afac, iwork, anorm, rcond, work,
539 \$ info )
540*
541* Check error code from CHPCON.
542*
543 IF( info.NE.0 )
544 \$ CALL alaerh( path, 'CHPCON', 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 CCHKHP
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 cchkhp(dotype, nn, nval, nns, nsval, thresh, tsterr, nmax, a, afac, ainv, b, x, xact, work, rwork, iwork, nout)
CCHKHP
Definition cchkhp.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 chpt01(uplo, n, a, afac, ipiv, c, ldc, rwork, resid)
CHPT01
Definition chpt01.f:113
subroutine claipd(n, a, inda, vinda)
CLAIPD
Definition claipd.f:83
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 cppt02(uplo, n, nrhs, a, x, ldx, b, ldb, rwork, resid)
CPPT02
Definition cppt02.f:123
subroutine cppt03(uplo, n, a, ainv, work, ldwork, rwork, rcond, resid)
CPPT03
Definition cppt03.f:110
subroutine cppt05(uplo, n, nrhs, ap, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
CPPT05
Definition cppt05.f:157
subroutine ccopy(n, cx, incx, cy, incy)
CCOPY
Definition ccopy.f:81
subroutine chpcon(uplo, n, ap, ipiv, anorm, rcond, work, info)
CHPCON
Definition chpcon.f:118
subroutine chprfs(uplo, n, nrhs, ap, afp, ipiv, b, ldb, x, ldx, ferr, berr, work, rwork, info)
CHPRFS
Definition chprfs.f:180
subroutine chptrf(uplo, n, ap, ipiv, info)
CHPTRF
Definition chptrf.f:159
subroutine chptri(uplo, n, ap, ipiv, work, info)
CHPTRI
Definition chptri.f:109
subroutine chptrs(uplo, n, nrhs, ap, ipiv, b, ldb, info)
CHPTRS
Definition chptrs.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