LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
cchkge.f
Go to the documentation of this file.
1*> \brief \b CCHKGE
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 CCHKGE( DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NNS,
12* NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B,
13* X, XACT, WORK, RWORK, IWORK, NOUT )
14*
15* .. Scalar Arguments ..
16* LOGICAL TSTERR
17* INTEGER NM, NMAX, NN, NNB, NNS, NOUT
18* REAL THRESH
19* ..
20* .. Array Arguments ..
21* LOGICAL DOTYPE( * )
22* INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
23* $ NVAL( * )
24* REAL RWORK( * )
25* COMPLEX A( * ), AFAC( * ), AINV( * ), B( * ),
26* $ WORK( * ), X( * ), XACT( * )
27* ..
28*
29*
30*> \par Purpose:
31* =============
32*>
33*> \verbatim
34*>
35*> CCHKGE tests CGETRF, -TRI, -TRS, -RFS, and -CON.
36*> \endverbatim
37*
38* Arguments:
39* ==========
40*
41*> \param[in] DOTYPE
42*> \verbatim
43*> DOTYPE is LOGICAL array, dimension (NTYPES)
44*> The matrix types to be used for testing. Matrices of type j
45*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
46*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
47*> \endverbatim
48*>
49*> \param[in] NM
50*> \verbatim
51*> NM is INTEGER
52*> The number of values of M contained in the vector MVAL.
53*> \endverbatim
54*>
55*> \param[in] MVAL
56*> \verbatim
57*> MVAL is INTEGER array, dimension (NM)
58*> The values of the matrix row dimension M.
59*> \endverbatim
60*>
61*> \param[in] NN
62*> \verbatim
63*> NN is INTEGER
64*> The number of values of N contained in the vector NVAL.
65*> \endverbatim
66*>
67*> \param[in] NVAL
68*> \verbatim
69*> NVAL is INTEGER array, dimension (NN)
70*> The values of the matrix column dimension N.
71*> \endverbatim
72*>
73*> \param[in] NNB
74*> \verbatim
75*> NNB is INTEGER
76*> The number of values of NB contained in the vector NBVAL.
77*> \endverbatim
78*>
79*> \param[in] NBVAL
80*> \verbatim
81*> NBVAL is INTEGER array, dimension (NNB)
82*> The values of the blocksize NB.
83*> \endverbatim
84*>
85*> \param[in] NNS
86*> \verbatim
87*> NNS is INTEGER
88*> The number of values of NRHS contained in the vector NSVAL.
89*> \endverbatim
90*>
91*> \param[in] NSVAL
92*> \verbatim
93*> NSVAL is INTEGER array, dimension (NNS)
94*> The values of the number of right hand sides NRHS.
95*> \endverbatim
96*>
97*> \param[in] THRESH
98*> \verbatim
99*> THRESH is REAL
100*> The threshold value for the test ratios. A result is
101*> included in the output file if RESULT >= THRESH. To have
102*> every test ratio printed, use THRESH = 0.
103*> \endverbatim
104*>
105*> \param[in] TSTERR
106*> \verbatim
107*> TSTERR is LOGICAL
108*> Flag that indicates whether error exits are to be tested.
109*> \endverbatim
110*>
111*> \param[in] NMAX
112*> \verbatim
113*> NMAX is INTEGER
114*> The maximum value permitted for M or N, used in dimensioning
115*> the work arrays.
116*> \endverbatim
117*>
118*> \param[out] A
119*> \verbatim
120*> A is COMPLEX array, dimension (NMAX*NMAX)
121*> \endverbatim
122*>
123*> \param[out] AFAC
124*> \verbatim
125*> AFAC is COMPLEX array, dimension (NMAX*NMAX)
126*> \endverbatim
127*>
128*> \param[out] AINV
129*> \verbatim
130*> AINV is COMPLEX array, dimension (NMAX*NMAX)
131*> \endverbatim
132*>
133*> \param[out] B
134*> \verbatim
135*> B is COMPLEX array, dimension (NMAX*NSMAX)
136*> where NSMAX is the largest entry in NSVAL.
137*> \endverbatim
138*>
139*> \param[out] X
140*> \verbatim
141*> X is COMPLEX array, dimension (NMAX*NSMAX)
142*> \endverbatim
143*>
144*> \param[out] XACT
145*> \verbatim
146*> XACT is COMPLEX array, dimension (NMAX*NSMAX)
147*> \endverbatim
148*>
149*> \param[out] WORK
150*> \verbatim
151*> WORK is COMPLEX array, dimension
152*> (NMAX*max(3,NSMAX))
153*> \endverbatim
154*>
155*> \param[out] RWORK
156*> \verbatim
157*> RWORK is REAL array, dimension
158*> (max(2*NMAX,2*NSMAX+NWORK))
159*> \endverbatim
160*>
161*> \param[out] IWORK
162*> \verbatim
163*> IWORK is INTEGER array, dimension (NMAX)
164*> \endverbatim
165*>
166*> \param[in] NOUT
167*> \verbatim
168*> NOUT is INTEGER
169*> The unit number for output.
170*> \endverbatim
171*
172* Authors:
173* ========
174*
175*> \author Univ. of Tennessee
176*> \author Univ. of California Berkeley
177*> \author Univ. of Colorado Denver
178*> \author NAG Ltd.
179*
180*> \ingroup complex_lin
181*
182* =====================================================================
183 SUBROUTINE cchkge( DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NNS,
184 $ NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B,
185 $ X, XACT, WORK, RWORK, IWORK, NOUT )
186*
187* -- LAPACK test routine --
188* -- LAPACK is a software package provided by Univ. of Tennessee, --
189* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
190*
191* .. Scalar Arguments ..
192 LOGICAL TSTERR
193 INTEGER NM, NMAX, NN, NNB, NNS, NOUT
194 REAL THRESH
195* ..
196* .. Array Arguments ..
197 LOGICAL DOTYPE( * )
198 INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
199 $ nval( * )
200 REAL RWORK( * )
201 COMPLEX A( * ), AFAC( * ), AINV( * ), B( * ),
202 $ work( * ), x( * ), xact( * )
203* ..
204*
205* =====================================================================
206*
207* .. Parameters ..
208 REAL ONE, ZERO
209 PARAMETER ( ONE = 1.0e+0, zero = 0.0e+0 )
210 INTEGER NTYPES
211 parameter( ntypes = 11 )
212 INTEGER NTESTS
213 parameter( ntests = 8 )
214 INTEGER NTRAN
215 parameter( ntran = 3 )
216* ..
217* .. Local Scalars ..
218 LOGICAL TRFCON, ZEROT
219 CHARACTER DIST, NORM, TRANS, TYPE, XTYPE
220 CHARACTER*3 PATH
221 INTEGER I, IM, IMAT, IN, INB, INFO, IOFF, IRHS, ITRAN,
222 $ izero, k, kl, ku, lda, lwork, m, mode, n, nb,
223 $ nerrs, nfail, nimat, nrhs, nrun, nt
224 REAL AINVNM, ANORM, ANORMI, ANORMO, CNDNUM, DUMMY,
225 $ RCOND, RCONDC, RCONDI, RCONDO
226* ..
227* .. Local Arrays ..
228 CHARACTER TRANSS( NTRAN )
229 INTEGER ISEED( 4 ), ISEEDY( 4 )
230 REAL RESULT( NTESTS )
231* ..
232* .. External Functions ..
233 REAL CLANGE, SGET06
234 EXTERNAL CLANGE, SGET06
235* ..
236* .. External Subroutines ..
237 EXTERNAL alaerh, alahd, alasum, cerrge, cgecon, cgerfs,
240 $ clatms, xlaenv
241* ..
242* .. Intrinsic Functions ..
243 INTRINSIC cmplx, max, min
244* ..
245* .. Scalars in Common ..
246 LOGICAL LERR, OK
247 CHARACTER*32 SRNAMT
248 INTEGER INFOT, NUNIT
249* ..
250* .. Common blocks ..
251 COMMON / infoc / infot, nunit, ok, lerr
252 COMMON / srnamc / srnamt
253* ..
254* .. Data statements ..
255 DATA iseedy / 1988, 1989, 1990, 1991 / ,
256 $ transs / 'N', 'T', 'C' /
257* ..
258* .. Executable Statements ..
259*
260* Initialize constants and the random number seed.
261*
262 path( 1: 1 ) = 'Complex precision'
263 path( 2: 3 ) = 'GE'
264 nrun = 0
265 nfail = 0
266 nerrs = 0
267 DO 10 i = 1, 4
268 iseed( i ) = iseedy( i )
269 10 CONTINUE
270*
271* Test the error exits
272*
273 CALL xlaenv( 1, 1 )
274 IF( tsterr )
275 $ CALL cerrge( path, nout )
276 infot = 0
277 CALL xlaenv( 2, 2 )
278*
279* Do for each value of M in MVAL
280*
281 DO 120 im = 1, nm
282 m = mval( im )
283 lda = max( 1, m )
284*
285* Do for each value of N in NVAL
286*
287 DO 110 in = 1, nn
288 n = nval( in )
289 xtype = 'N'
290 nimat = ntypes
291 IF( m.LE.0 .OR. n.LE.0 )
292 $ nimat = 1
293*
294 DO 100 imat = 1, nimat
295*
296* Do the tests only if DOTYPE( IMAT ) is true.
297*
298 IF( .NOT.dotype( imat ) )
299 $ GO TO 100
300*
301* Skip types 5, 6, or 7 if the matrix size is too small.
302*
303 zerot = imat.GE.5 .AND. imat.LE.7
304 IF( zerot .AND. n.LT.imat-4 )
305 $ GO TO 100
306*
307* Set up parameters with CLATB4 and generate a test matrix
308* with CLATMS.
309*
310 CALL clatb4( path, imat, m, n, TYPE, kl, ku, anorm, mode,
311 $ cndnum, dist )
312*
313 srnamt = 'CLATMS'
314 CALL clatms( m, n, dist, iseed, TYPE, rwork, mode,
315 $ cndnum, anorm, kl, ku, 'No packing', a, lda,
316 $ work, info )
317*
318* Check error code from CLATMS.
319*
320 IF( info.NE.0 ) THEN
321 CALL alaerh( path, 'CLATMS', info, 0, ' ', m, n, -1,
322 $ -1, -1, imat, nfail, nerrs, nout )
323 GO TO 100
324 END IF
325*
326* For types 5-7, zero one or more columns of the matrix to
327* test that INFO is returned correctly.
328*
329 IF( zerot ) THEN
330 IF( imat.EQ.5 ) THEN
331 izero = 1
332 ELSE IF( imat.EQ.6 ) THEN
333 izero = min( m, n )
334 ELSE
335 izero = min( m, n ) / 2 + 1
336 END IF
337 ioff = ( izero-1 )*lda
338 IF( imat.LT.7 ) THEN
339 DO 20 i = 1, m
340 a( ioff+i ) = zero
341 20 CONTINUE
342 ELSE
343 CALL claset( 'Full', m, n-izero+1, cmplx( zero ),
344 $ cmplx( zero ), a( ioff+1 ), lda )
345 END IF
346 ELSE
347 izero = 0
348 END IF
349*
350* These lines, if used in place of the calls in the DO 60
351* loop, cause the code to bomb on a Sun SPARCstation.
352*
353* ANORMO = CLANGE( 'O', M, N, A, LDA, RWORK )
354* ANORMI = CLANGE( 'I', M, N, A, LDA, RWORK )
355*
356* Do for each blocksize in NBVAL
357*
358 DO 90 inb = 1, nnb
359 nb = nbval( inb )
360 CALL xlaenv( 1, nb )
361*
362* Compute the LU factorization of the matrix.
363*
364 CALL clacpy( 'Full', m, n, a, lda, afac, lda )
365 srnamt = 'CGETRF'
366 CALL cgetrf( m, n, afac, lda, iwork, info )
367*
368* Check error code from CGETRF.
369*
370 IF( info.NE.izero )
371 $ CALL alaerh( path, 'CGETRF', info, izero, ' ', m,
372 $ n, -1, -1, nb, imat, nfail, nerrs,
373 $ nout )
374 trfcon = .false.
375*
376*+ TEST 1
377* Reconstruct matrix from factors and compute residual.
378*
379 CALL clacpy( 'Full', m, n, afac, lda, ainv, lda )
380 CALL cget01( m, n, a, lda, ainv, lda, iwork, rwork,
381 $ result( 1 ) )
382 nt = 1
383*
384*+ TEST 2
385* Form the inverse if the factorization was successful
386* and compute the residual.
387*
388 IF( m.EQ.n .AND. info.EQ.0 ) THEN
389 CALL clacpy( 'Full', n, n, afac, lda, ainv, lda )
390 srnamt = 'CGETRI'
391 nrhs = nsval( 1 )
392 lwork = nmax*max( 3, nrhs )
393 CALL cgetri( n, ainv, lda, iwork, work, lwork,
394 $ info )
395*
396* Check error code from CGETRI.
397*
398 IF( info.NE.0 )
399 $ CALL alaerh( path, 'CGETRI', info, 0, ' ', n, n,
400 $ -1, -1, nb, imat, nfail, nerrs,
401 $ nout )
402*
403* Compute the residual for the matrix times its
404* inverse. Also compute the 1-norm condition number
405* of A.
406*
407 CALL cget03( n, a, lda, ainv, lda, work, lda,
408 $ rwork, rcondo, result( 2 ) )
409 anormo = clange( 'O', m, n, a, lda, rwork )
410*
411* Compute the infinity-norm condition number of A.
412*
413 anormi = clange( 'I', m, n, a, lda, rwork )
414 ainvnm = clange( 'I', n, n, ainv, lda, rwork )
415 IF( anormi.LE.zero .OR. ainvnm.LE.zero ) THEN
416 rcondi = one
417 ELSE
418 rcondi = ( one / anormi ) / ainvnm
419 END IF
420 nt = 2
421 ELSE
422*
423* Do only the condition estimate if INFO > 0.
424*
425 trfcon = .true.
426 anormo = clange( 'O', m, n, a, lda, rwork )
427 anormi = clange( 'I', m, n, a, lda, rwork )
428 rcondo = zero
429 rcondi = zero
430 END IF
431*
432* Print information about the tests so far that did not
433* pass the threshold.
434*
435 DO 30 k = 1, nt
436 IF( result( k ).GE.thresh ) THEN
437 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
438 $ CALL alahd( nout, path )
439 WRITE( nout, fmt = 9999 )m, n, nb, imat, k,
440 $ result( k )
441 nfail = nfail + 1
442 END IF
443 30 CONTINUE
444 nrun = nrun + nt
445*
446* Skip the remaining tests if this is not the first
447* block size or if M .ne. N. Skip the solve tests if
448* the matrix is singular.
449*
450 IF( inb.GT.1 .OR. m.NE.n )
451 $ GO TO 90
452 IF( trfcon )
453 $ GO TO 70
454*
455 DO 60 irhs = 1, nns
456 nrhs = nsval( irhs )
457 xtype = 'N'
458*
459 DO 50 itran = 1, ntran
460 trans = transs( itran )
461 IF( itran.EQ.1 ) THEN
462 rcondc = rcondo
463 ELSE
464 rcondc = rcondi
465 END IF
466*
467*+ TEST 3
468* Solve and compute residual for A * X = B.
469*
470 srnamt = 'CLARHS'
471 CALL clarhs( path, xtype, ' ', trans, n, n, kl,
472 $ ku, nrhs, a, lda, xact, lda, b,
473 $ lda, iseed, info )
474 xtype = 'C'
475*
476 CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
477 srnamt = 'CGETRS'
478 CALL cgetrs( trans, n, nrhs, afac, lda, iwork,
479 $ x, lda, info )
480*
481* Check error code from CGETRS.
482*
483 IF( info.NE.0 )
484 $ CALL alaerh( path, 'CGETRS', info, 0, trans,
485 $ n, n, -1, -1, nrhs, imat, nfail,
486 $ nerrs, nout )
487*
488 CALL clacpy( 'Full', n, nrhs, b, lda, work,
489 $ lda )
490 CALL cget02( trans, n, n, nrhs, a, lda, x, lda,
491 $ work, lda, rwork, result( 3 ) )
492*
493*+ TEST 4
494* Check solution from generated exact solution.
495*
496 CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
497 $ result( 4 ) )
498*
499*+ TESTS 5, 6, and 7
500* Use iterative refinement to improve the
501* solution.
502*
503 srnamt = 'CGERFS'
504 CALL cgerfs( trans, n, nrhs, a, lda, afac, lda,
505 $ iwork, b, lda, x, lda, rwork,
506 $ rwork( nrhs+1 ), work,
507 $ rwork( 2*nrhs+1 ), info )
508*
509* Check error code from CGERFS.
510*
511 IF( info.NE.0 )
512 $ CALL alaerh( path, 'CGERFS', info, 0, trans,
513 $ n, n, -1, -1, nrhs, imat, nfail,
514 $ nerrs, nout )
515*
516 CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
517 $ result( 5 ) )
518 CALL cget07( trans, n, nrhs, a, lda, b, lda, x,
519 $ lda, xact, lda, rwork, .true.,
520 $ rwork( nrhs+1 ), result( 6 ) )
521*
522* Print information about the tests that did not
523* pass the threshold.
524*
525 DO 40 k = 3, 7
526 IF( result( k ).GE.thresh ) THEN
527 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
528 $ CALL alahd( nout, path )
529 WRITE( nout, fmt = 9998 )trans, n, nrhs,
530 $ imat, k, result( k )
531 nfail = nfail + 1
532 END IF
533 40 CONTINUE
534 nrun = nrun + 5
535 50 CONTINUE
536 60 CONTINUE
537*
538*+ TEST 8
539* Get an estimate of RCOND = 1/CNDNUM.
540*
541 70 CONTINUE
542 DO 80 itran = 1, 2
543 IF( itran.EQ.1 ) THEN
544 anorm = anormo
545 rcondc = rcondo
546 norm = 'O'
547 ELSE
548 anorm = anormi
549 rcondc = rcondi
550 norm = 'I'
551 END IF
552 srnamt = 'CGECON'
553 CALL cgecon( norm, n, afac, lda, anorm, rcond,
554 $ work, rwork, info )
555*
556* Check error code from CGECON.
557*
558 IF( info.NE.0 )
559 $ CALL alaerh( path, 'CGECON', info, 0, norm, n,
560 $ n, -1, -1, -1, imat, nfail, nerrs,
561 $ nout )
562*
563* This line is needed on a Sun SPARCstation.
564*
565 dummy = rcond
566*
567 result( 8 ) = sget06( rcond, rcondc )
568*
569* Print information about the tests that did not pass
570* the threshold.
571*
572 IF( result( 8 ).GE.thresh ) THEN
573 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
574 $ CALL alahd( nout, path )
575 WRITE( nout, fmt = 9997 )norm, n, imat, 8,
576 $ result( 8 )
577 nfail = nfail + 1
578 END IF
579 nrun = nrun + 1
580 80 CONTINUE
581 90 CONTINUE
582 100 CONTINUE
583*
584 110 CONTINUE
585 120 CONTINUE
586*
587* Print a summary of the results.
588*
589 CALL alasum( path, nout, nfail, nrun, nerrs )
590*
591 9999 FORMAT( ' M = ', i5, ', N =', i5, ', NB =', i4, ', type ', i2,
592 $ ', test(', i2, ') =', g12.5 )
593 9998 FORMAT( ' TRANS=''', a1, ''', N =', i5, ', NRHS=', i3, ', type ',
594 $ i2, ', test(', i2, ') =', g12.5 )
595 9997 FORMAT( ' NORM =''', a1, ''', N =', i5, ',', 10x, ' type ', i2,
596 $ ', test(', i2, ') =', g12.5 )
597 RETURN
598*
599* End of CCHKGE
600*
601 END
subroutine alasum(type, nout, nfail, nrun, nerrs)
ALASUM
Definition alasum.f:73
subroutine cget02(trans, m, n, nrhs, a, lda, x, ldx, b, ldb, rwork, resid)
CGET02
Definition cget02.f:134
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 xlaenv(ispec, nvalue)
XLAENV
Definition xlaenv.f:81
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 cchkge(dotype, nm, mval, nn, nval, nnb, nbval, nns, nsval, thresh, tsterr, nmax, a, afac, ainv, b, x, xact, work, rwork, iwork, nout)
CCHKGE
Definition cchkge.f:186
subroutine cerrge(path, nunit)
CERRGE
Definition cerrge.f:55
subroutine cget01(m, n, a, lda, afac, ldafac, ipiv, rwork, resid)
CGET01
Definition cget01.f:108
subroutine cget03(n, a, lda, ainv, ldainv, work, ldwork, rwork, rcond, resid)
CGET03
Definition cget03.f:110
subroutine cget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
CGET04
Definition cget04.f:102
subroutine cget07(trans, n, nrhs, a, lda, b, ldb, x, ldx, xact, ldxact, ferr, chkferr, berr, reslts)
CGET07
Definition cget07.f:166
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 cgecon(norm, n, a, lda, anorm, rcond, work, rwork, info)
CGECON
Definition cgecon.f:132
subroutine cgerfs(trans, n, nrhs, a, lda, af, ldaf, ipiv, b, ldb, x, ldx, ferr, berr, work, rwork, info)
CGERFS
Definition cgerfs.f:186
subroutine cgetrf(m, n, a, lda, ipiv, info)
CGETRF
Definition cgetrf.f:108
subroutine cgetri(n, a, lda, ipiv, work, lwork, info)
CGETRI
Definition cgetri.f:114
subroutine cgetrs(trans, n, nrhs, a, lda, ipiv, b, ldb, info)
CGETRS
Definition cgetrs.f:121
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 claset(uplo, m, n, alpha, beta, a, lda)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition claset.f:106