LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
cdrvgex.f
Go to the documentation of this file.
1*> \brief \b CDRVGEX
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 CDRVGE( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
12* A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK,
13* RWORK, IWORK, NOUT )
14*
15* .. Scalar Arguments ..
16* LOGICAL TSTERR
17* INTEGER NMAX, NN, NOUT, NRHS
18* REAL THRESH
19* ..
20* .. Array Arguments ..
21* LOGICAL DOTYPE( * )
22* INTEGER IWORK( * ), NVAL( * )
23* REAL RWORK( * ), S( * )
24* COMPLEX A( * ), AFAC( * ), ASAV( * ), B( * ),
25* $ BSAV( * ), WORK( * ), X( * ), XACT( * )
26* ..
27*
28*
29*> \par Purpose:
30* =============
31*>
32*> \verbatim
33*>
34*> CDRVGE tests the driver routines CGESV, -SVX, and -SVXX.
35*>
36*> Note that this file is used only when the XBLAS are available,
37*> otherwise cdrvge.f defines this subroutine.
38*> \endverbatim
39*
40* Arguments:
41* ==========
42*
43*> \param[in] DOTYPE
44*> \verbatim
45*> DOTYPE is LOGICAL array, dimension (NTYPES)
46*> The matrix types to be used for testing. Matrices of type j
47*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
48*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
49*> \endverbatim
50*>
51*> \param[in] NN
52*> \verbatim
53*> NN is INTEGER
54*> The number of values of N contained in the vector NVAL.
55*> \endverbatim
56*>
57*> \param[in] NVAL
58*> \verbatim
59*> NVAL is INTEGER array, dimension (NN)
60*> The values of the matrix column dimension N.
61*> \endverbatim
62*>
63*> \param[in] NRHS
64*> \verbatim
65*> NRHS is INTEGER
66*> The number of right hand side vectors to be generated for
67*> each linear system.
68*> \endverbatim
69*>
70*> \param[in] THRESH
71*> \verbatim
72*> THRESH is REAL
73*> The threshold value for the test ratios. A result is
74*> included in the output file if RESULT >= THRESH. To have
75*> every test ratio printed, use THRESH = 0.
76*> \endverbatim
77*>
78*> \param[in] TSTERR
79*> \verbatim
80*> TSTERR is LOGICAL
81*> Flag that indicates whether error exits are to be tested.
82*> \endverbatim
83*>
84*> \param[in] NMAX
85*> \verbatim
86*> NMAX is INTEGER
87*> The maximum value permitted for N, used in dimensioning the
88*> work arrays.
89*> \endverbatim
90*>
91*> \param[out] A
92*> \verbatim
93*> A is COMPLEX array, dimension (NMAX*NMAX)
94*> \endverbatim
95*>
96*> \param[out] AFAC
97*> \verbatim
98*> AFAC is COMPLEX array, dimension (NMAX*NMAX)
99*> \endverbatim
100*>
101*> \param[out] ASAV
102*> \verbatim
103*> ASAV is COMPLEX array, dimension (NMAX*NMAX)
104*> \endverbatim
105*>
106*> \param[out] B
107*> \verbatim
108*> B is COMPLEX array, dimension (NMAX*NRHS)
109*> \endverbatim
110*>
111*> \param[out] BSAV
112*> \verbatim
113*> BSAV is COMPLEX array, dimension (NMAX*NRHS)
114*> \endverbatim
115*>
116*> \param[out] X
117*> \verbatim
118*> X is COMPLEX array, dimension (NMAX*NRHS)
119*> \endverbatim
120*>
121*> \param[out] XACT
122*> \verbatim
123*> XACT is COMPLEX array, dimension (NMAX*NRHS)
124*> \endverbatim
125*>
126*> \param[out] S
127*> \verbatim
128*> S is REAL array, dimension (2*NMAX)
129*> \endverbatim
130*>
131*> \param[out] WORK
132*> \verbatim
133*> WORK is COMPLEX array, dimension
134*> (NMAX*max(3,NRHS))
135*> \endverbatim
136*>
137*> \param[out] RWORK
138*> \verbatim
139*> RWORK is REAL array, dimension (2*NRHS+NMAX)
140*> \endverbatim
141*>
142*> \param[out] IWORK
143*> \verbatim
144*> IWORK is INTEGER array, dimension (NMAX)
145*> \endverbatim
146*>
147*> \param[in] NOUT
148*> \verbatim
149*> NOUT is INTEGER
150*> The unit number for output.
151*> \endverbatim
152*
153* Authors:
154* ========
155*
156*> \author Univ. of Tennessee
157*> \author Univ. of California Berkeley
158*> \author Univ. of Colorado Denver
159*> \author NAG Ltd.
160*
161*> \ingroup complex_lin
162*
163* =====================================================================
164 SUBROUTINE cdrvge( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
165 $ A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK,
166 $ RWORK, IWORK, NOUT )
167*
168* -- LAPACK test routine --
169* -- LAPACK is a software package provided by Univ. of Tennessee, --
170* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
171*
172* .. Scalar Arguments ..
173 LOGICAL TSTERR
174 INTEGER NMAX, NN, NOUT, NRHS
175 REAL THRESH
176* ..
177* .. Array Arguments ..
178 LOGICAL DOTYPE( * )
179 INTEGER IWORK( * ), NVAL( * )
180 REAL RWORK( * ), S( * )
181 COMPLEX A( * ), AFAC( * ), ASAV( * ), B( * ),
182 $ bsav( * ), work( * ), x( * ), xact( * )
183* ..
184*
185* =====================================================================
186*
187* .. Parameters ..
188 REAL ONE, ZERO
189 PARAMETER ( ONE = 1.0e+0, zero = 0.0e+0 )
190 INTEGER NTYPES
191 parameter( ntypes = 11 )
192 INTEGER NTESTS
193 parameter( ntests = 7 )
194 INTEGER NTRAN
195 parameter( ntran = 3 )
196* ..
197* .. Local Scalars ..
198 LOGICAL EQUIL, NOFACT, PREFAC, TRFCON, ZEROT
199 CHARACTER DIST, EQUED, FACT, TRANS, TYPE, XTYPE
200 CHARACTER*3 PATH
201 INTEGER I, IEQUED, IFACT, IMAT, IN, INFO, IOFF, ITRAN,
202 $ izero, k, k1, kl, ku, lda, lwork, mode, n, nb,
203 $ nbmin, nerrs, nfact, nfail, nimat, nrun, nt,
204 $ n_err_bnds
205 REAL AINVNM, AMAX, ANORM, ANORMI, ANORMO, CNDNUM,
206 $ COLCND, RCOND, RCONDC, RCONDI, RCONDO, ROLDC,
207 $ roldi, roldo, rowcnd, rpvgrw, rpvgrw_svxx
208* ..
209* .. Local Arrays ..
210 CHARACTER EQUEDS( 4 ), FACTS( 3 ), TRANSS( NTRAN )
211 INTEGER ISEED( 4 ), ISEEDY( 4 )
212 REAL RDUM( 1 ), RESULT( NTESTS ), BERR( NRHS ),
213 $ errbnds_n( nrhs, 3 ), errbnds_c( nrhs, 3 )
214* ..
215* .. External Functions ..
216 LOGICAL LSAME
217 REAL CLANGE, CLANTR, SGET06, SLAMCH, CLA_GERPVGRW
218 EXTERNAL lsame, clange, clantr, sget06, slamch,
219 $ cla_gerpvgrw
220* ..
221* .. External Subroutines ..
222 EXTERNAL aladhd, alaerh, alasvm, cerrvx, cgeequ, cgesv,
226* ..
227* .. Intrinsic Functions ..
228 INTRINSIC abs, cmplx, max
229* ..
230* .. Scalars in Common ..
231 LOGICAL LERR, OK
232 CHARACTER*32 SRNAMT
233 INTEGER INFOT, NUNIT
234* ..
235* .. Common blocks ..
236 COMMON / infoc / infot, nunit, ok, lerr
237 COMMON / srnamc / srnamt
238* ..
239* .. Data statements ..
240 DATA iseedy / 1988, 1989, 1990, 1991 /
241 DATA transs / 'N', 'T', 'C' /
242 DATA facts / 'F', 'N', 'E' /
243 DATA equeds / 'N', 'R', 'C', 'B' /
244* ..
245* .. Executable Statements ..
246*
247* Initialize constants and the random number seed.
248*
249 path( 1: 1 ) = 'Complex precision'
250 path( 2: 3 ) = 'GE'
251 nrun = 0
252 nfail = 0
253 nerrs = 0
254 DO 10 i = 1, 4
255 iseed( i ) = iseedy( i )
256 10 CONTINUE
257*
258* Test the error exits
259*
260 IF( tsterr )
261 $ CALL cerrvx( path, nout )
262 infot = 0
263*
264* Set the block size and minimum block size for testing.
265*
266 nb = 1
267 nbmin = 2
268 CALL xlaenv( 1, nb )
269 CALL xlaenv( 2, nbmin )
270*
271* Do for each value of N in NVAL
272*
273 DO 90 in = 1, nn
274 n = nval( in )
275 lda = max( n, 1 )
276 xtype = 'N'
277 nimat = ntypes
278 IF( n.LE.0 )
279 $ nimat = 1
280*
281 DO 80 imat = 1, nimat
282*
283* Do the tests only if DOTYPE( IMAT ) is true.
284*
285 IF( .NOT.dotype( imat ) )
286 $ GO TO 80
287*
288* Skip types 5, 6, or 7 if the matrix size is too small.
289*
290 zerot = imat.GE.5 .AND. imat.LE.7
291 IF( zerot .AND. n.LT.imat-4 )
292 $ GO TO 80
293*
294* Set up parameters with CLATB4 and generate a test matrix
295* with CLATMS.
296*
297 CALL clatb4( path, imat, n, n, TYPE, kl, ku, anorm, mode,
298 $ cndnum, dist )
299 rcondc = one / cndnum
300*
301 srnamt = 'CLATMS'
302 CALL clatms( n, n, dist, iseed, TYPE, rwork, mode, cndnum,
303 $ anorm, kl, ku, 'No packing', a, lda, work,
304 $ info )
305*
306* Check error code from CLATMS.
307*
308 IF( info.NE.0 ) THEN
309 CALL alaerh( path, 'CLATMS', info, 0, ' ', n, n, -1, -1,
310 $ -1, imat, nfail, nerrs, nout )
311 GO TO 80
312 END IF
313*
314* For types 5-7, zero one or more columns of the matrix to
315* test that INFO is returned correctly.
316*
317 IF( zerot ) THEN
318 IF( imat.EQ.5 ) THEN
319 izero = 1
320 ELSE IF( imat.EQ.6 ) THEN
321 izero = n
322 ELSE
323 izero = n / 2 + 1
324 END IF
325 ioff = ( izero-1 )*lda
326 IF( imat.LT.7 ) THEN
327 DO 20 i = 1, n
328 a( ioff+i ) = zero
329 20 CONTINUE
330 ELSE
331 CALL claset( 'Full', n, n-izero+1, cmplx( zero ),
332 $ cmplx( zero ), a( ioff+1 ), lda )
333 END IF
334 ELSE
335 izero = 0
336 END IF
337*
338* Save a copy of the matrix A in ASAV.
339*
340 CALL clacpy( 'Full', n, n, a, lda, asav, lda )
341*
342 DO 70 iequed = 1, 4
343 equed = equeds( iequed )
344 IF( iequed.EQ.1 ) THEN
345 nfact = 3
346 ELSE
347 nfact = 1
348 END IF
349*
350 DO 60 ifact = 1, nfact
351 fact = facts( ifact )
352 prefac = lsame( fact, 'F' )
353 nofact = lsame( fact, 'N' )
354 equil = lsame( fact, 'E' )
355*
356 IF( zerot ) THEN
357 IF( prefac )
358 $ GO TO 60
359 rcondo = zero
360 rcondi = zero
361*
362 ELSE IF( .NOT.nofact ) THEN
363*
364* Compute the condition number for comparison with
365* the value returned by CGESVX (FACT = 'N' reuses
366* the condition number from the previous iteration
367* with FACT = 'F').
368*
369 CALL clacpy( 'Full', n, n, asav, lda, afac, lda )
370 IF( equil .OR. iequed.GT.1 ) THEN
371*
372* Compute row and column scale factors to
373* equilibrate the matrix A.
374*
375 CALL cgeequ( n, n, afac, lda, s, s( n+1 ),
376 $ rowcnd, colcnd, amax, info )
377 IF( info.EQ.0 .AND. n.GT.0 ) THEN
378 IF( lsame( equed, 'R' ) ) THEN
379 rowcnd = zero
380 colcnd = one
381 ELSE IF( lsame( equed, 'C' ) ) THEN
382 rowcnd = one
383 colcnd = zero
384 ELSE IF( lsame( equed, 'B' ) ) THEN
385 rowcnd = zero
386 colcnd = zero
387 END IF
388*
389* Equilibrate the matrix.
390*
391 CALL claqge( n, n, afac, lda, s, s( n+1 ),
392 $ rowcnd, colcnd, amax, equed )
393 END IF
394 END IF
395*
396* Save the condition number of the non-equilibrated
397* system for use in CGET04.
398*
399 IF( equil ) THEN
400 roldo = rcondo
401 roldi = rcondi
402 END IF
403*
404* Compute the 1-norm and infinity-norm of A.
405*
406 anormo = clange( '1', n, n, afac, lda, rwork )
407 anormi = clange( 'I', n, n, afac, lda, rwork )
408*
409* Factor the matrix A.
410*
411 CALL cgetrf( n, n, afac, lda, iwork, info )
412*
413* Form the inverse of A.
414*
415 CALL clacpy( 'Full', n, n, afac, lda, a, lda )
416 lwork = nmax*max( 3, nrhs )
417 CALL cgetri( n, a, lda, iwork, work, lwork, info )
418*
419* Compute the 1-norm condition number of A.
420*
421 ainvnm = clange( '1', n, n, a, lda, rwork )
422 IF( anormo.LE.zero .OR. ainvnm.LE.zero ) THEN
423 rcondo = one
424 ELSE
425 rcondo = ( one / anormo ) / ainvnm
426 END IF
427*
428* Compute the infinity-norm condition number of A.
429*
430 ainvnm = clange( 'I', n, n, a, lda, rwork )
431 IF( anormi.LE.zero .OR. ainvnm.LE.zero ) THEN
432 rcondi = one
433 ELSE
434 rcondi = ( one / anormi ) / ainvnm
435 END IF
436 END IF
437*
438 DO 50 itran = 1, ntran
439*
440* Do for each value of TRANS.
441*
442 trans = transs( itran )
443 IF( itran.EQ.1 ) THEN
444 rcondc = rcondo
445 ELSE
446 rcondc = rcondi
447 END IF
448*
449* Restore the matrix A.
450*
451 CALL clacpy( 'Full', n, n, asav, lda, a, lda )
452*
453* Form an exact solution and set the right hand side.
454*
455 srnamt = 'CLARHS'
456 CALL clarhs( path, xtype, 'Full', trans, n, n, kl,
457 $ ku, nrhs, a, lda, xact, lda, b, lda,
458 $ iseed, info )
459 xtype = 'C'
460 CALL clacpy( 'Full', n, nrhs, b, lda, bsav, lda )
461*
462 IF( nofact .AND. itran.EQ.1 ) THEN
463*
464* --- Test CGESV ---
465*
466* Compute the LU factorization of the matrix and
467* solve the system.
468*
469 CALL clacpy( 'Full', n, n, a, lda, afac, lda )
470 CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
471*
472 srnamt = 'CGESV '
473 CALL cgesv( n, nrhs, afac, lda, iwork, x, lda,
474 $ info )
475*
476* Check error code from CGESV .
477*
478 IF( info.NE.izero )
479 $ CALL alaerh( path, 'CGESV ', info, izero,
480 $ ' ', n, n, -1, -1, nrhs, imat,
481 $ nfail, nerrs, nout )
482*
483* Reconstruct matrix from factors and compute
484* residual.
485*
486 CALL cget01( n, n, a, lda, afac, lda, iwork,
487 $ rwork, result( 1 ) )
488 nt = 1
489 IF( izero.EQ.0 ) THEN
490*
491* Compute residual of the computed solution.
492*
493 CALL clacpy( 'Full', n, nrhs, b, lda, work,
494 $ lda )
495 CALL cget02( 'No transpose', n, n, nrhs, a,
496 $ lda, x, lda, work, lda, rwork,
497 $ result( 2 ) )
498*
499* Check solution from generated exact solution.
500*
501 CALL cget04( n, nrhs, x, lda, xact, lda,
502 $ rcondc, result( 3 ) )
503 nt = 3
504 END IF
505*
506* Print information about the tests that did not
507* pass the threshold.
508*
509 DO 30 k = 1, nt
510 IF( result( k ).GE.thresh ) THEN
511 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
512 $ CALL aladhd( nout, path )
513 WRITE( nout, fmt = 9999 )'CGESV ', n,
514 $ imat, k, result( k )
515 nfail = nfail + 1
516 END IF
517 30 CONTINUE
518 nrun = nrun + nt
519 END IF
520*
521* --- Test CGESVX ---
522*
523 IF( .NOT.prefac )
524 $ CALL claset( 'Full', n, n, cmplx( zero ),
525 $ cmplx( zero ), afac, lda )
526 CALL claset( 'Full', n, nrhs, cmplx( zero ),
527 $ cmplx( zero ), x, lda )
528 IF( iequed.GT.1 .AND. n.GT.0 ) THEN
529*
530* Equilibrate the matrix if FACT = 'F' and
531* EQUED = 'R', 'C', or 'B'.
532*
533 CALL claqge( n, n, a, lda, s, s( n+1 ), rowcnd,
534 $ colcnd, amax, equed )
535 END IF
536*
537* Solve the system and compute the condition number
538* and error bounds using CGESVX.
539*
540 srnamt = 'CGESVX'
541 CALL cgesvx( fact, trans, n, nrhs, a, lda, afac,
542 $ lda, iwork, equed, s, s( n+1 ), b,
543 $ lda, x, lda, rcond, rwork,
544 $ rwork( nrhs+1 ), work,
545 $ rwork( 2*nrhs+1 ), info )
546*
547* Check the error code from CGESVX.
548*
549 IF( info.NE.izero )
550 $ CALL alaerh( path, 'CGESVX', info, izero,
551 $ fact // trans, n, n, -1, -1, nrhs,
552 $ imat, nfail, nerrs, nout )
553*
554* Compare RWORK(2*NRHS+1) from CGESVX with the
555* computed reciprocal pivot growth factor RPVGRW
556*
557 IF( info.NE.0 ) THEN
558 rpvgrw = clantr( 'M', 'U', 'N', info, info,
559 $ afac, lda, rdum )
560 IF( rpvgrw.EQ.zero ) THEN
561 rpvgrw = one
562 ELSE
563 rpvgrw = clange( 'M', n, info, a, lda,
564 $ rdum ) / rpvgrw
565 END IF
566 ELSE
567 rpvgrw = clantr( 'M', 'U', 'N', n, n, afac, lda,
568 $ rdum )
569 IF( rpvgrw.EQ.zero ) THEN
570 rpvgrw = one
571 ELSE
572 rpvgrw = clange( 'M', n, n, a, lda, rdum ) /
573 $ rpvgrw
574 END IF
575 END IF
576 result( 7 ) = abs( rpvgrw-rwork( 2*nrhs+1 ) ) /
577 $ max( rwork( 2*nrhs+1 ), rpvgrw ) /
578 $ slamch( 'E' )
579*
580 IF( .NOT.prefac ) THEN
581*
582* Reconstruct matrix from factors and compute
583* residual.
584*
585 CALL cget01( n, n, a, lda, afac, lda, iwork,
586 $ rwork( 2*nrhs+1 ), result( 1 ) )
587 k1 = 1
588 ELSE
589 k1 = 2
590 END IF
591*
592 IF( info.EQ.0 ) THEN
593 trfcon = .false.
594*
595* Compute residual of the computed solution.
596*
597 CALL clacpy( 'Full', n, nrhs, bsav, lda, work,
598 $ lda )
599 CALL cget02( trans, n, n, nrhs, asav, lda, x,
600 $ lda, work, lda, rwork( 2*nrhs+1 ),
601 $ result( 2 ) )
602*
603* Check solution from generated exact solution.
604*
605 IF( nofact .OR. ( prefac .AND. lsame( equed,
606 $ 'N' ) ) ) THEN
607 CALL cget04( n, nrhs, x, lda, xact, lda,
608 $ rcondc, result( 3 ) )
609 ELSE
610 IF( itran.EQ.1 ) THEN
611 roldc = roldo
612 ELSE
613 roldc = roldi
614 END IF
615 CALL cget04( n, nrhs, x, lda, xact, lda,
616 $ roldc, result( 3 ) )
617 END IF
618*
619* Check the error bounds from iterative
620* refinement.
621*
622 CALL cget07( trans, n, nrhs, asav, lda, b, lda,
623 $ x, lda, xact, lda, rwork, .true.,
624 $ rwork( nrhs+1 ), result( 4 ) )
625 ELSE
626 trfcon = .true.
627 END IF
628*
629* Compare RCOND from CGESVX with the computed value
630* in RCONDC.
631*
632 result( 6 ) = sget06( rcond, rcondc )
633*
634* Print information about the tests that did not pass
635* the threshold.
636*
637 IF( .NOT.trfcon ) THEN
638 DO 40 k = k1, ntests
639 IF( result( k ).GE.thresh ) THEN
640 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
641 $ CALL aladhd( nout, path )
642 IF( prefac ) THEN
643 WRITE( nout, fmt = 9997 )'CGESVX',
644 $ fact, trans, n, equed, imat, k,
645 $ result( k )
646 ELSE
647 WRITE( nout, fmt = 9998 )'CGESVX',
648 $ fact, trans, n, imat, k, result( k )
649 END IF
650 nfail = nfail + 1
651 END IF
652 40 CONTINUE
653 nrun = nrun + 7 - k1
654 ELSE
655 IF( result( 1 ).GE.thresh .AND. .NOT.prefac )
656 $ THEN
657 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
658 $ CALL aladhd( nout, path )
659 IF( prefac ) THEN
660 WRITE( nout, fmt = 9997 )'CGESVX', fact,
661 $ trans, n, equed, imat, 1, result( 1 )
662 ELSE
663 WRITE( nout, fmt = 9998 )'CGESVX', fact,
664 $ trans, n, imat, 1, result( 1 )
665 END IF
666 nfail = nfail + 1
667 nrun = nrun + 1
668 END IF
669 IF( result( 6 ).GE.thresh ) THEN
670 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
671 $ CALL aladhd( nout, path )
672 IF( prefac ) THEN
673 WRITE( nout, fmt = 9997 )'CGESVX', fact,
674 $ trans, n, equed, imat, 6, result( 6 )
675 ELSE
676 WRITE( nout, fmt = 9998 )'CGESVX', fact,
677 $ trans, n, imat, 6, result( 6 )
678 END IF
679 nfail = nfail + 1
680 nrun = nrun + 1
681 END IF
682 IF( result( 7 ).GE.thresh ) THEN
683 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
684 $ CALL aladhd( nout, path )
685 IF( prefac ) THEN
686 WRITE( nout, fmt = 9997 )'CGESVX', fact,
687 $ trans, n, equed, imat, 7, result( 7 )
688 ELSE
689 WRITE( nout, fmt = 9998 )'CGESVX', fact,
690 $ trans, n, imat, 7, result( 7 )
691 END IF
692 nfail = nfail + 1
693 nrun = nrun + 1
694 END IF
695*
696 END IF
697*
698* --- Test CGESVXX ---
699*
700* Restore the matrices A and B.
701*
702
703 CALL clacpy( 'Full', n, n, asav, lda, a, lda )
704 CALL clacpy( 'Full', n, nrhs, bsav, lda, b, lda )
705
706 IF( .NOT.prefac )
707 $ CALL claset( 'Full', n, n, cmplx( zero ),
708 $ cmplx( zero ), afac, lda )
709 CALL claset( 'Full', n, nrhs, cmplx( zero ),
710 $ cmplx( zero ), x, lda )
711 IF( iequed.GT.1 .AND. n.GT.0 ) THEN
712*
713* Equilibrate the matrix if FACT = 'F' and
714* EQUED = 'R', 'C', or 'B'.
715*
716 CALL claqge( n, n, a, lda, s, s( n+1 ), rowcnd,
717 $ colcnd, amax, equed )
718 END IF
719*
720* Solve the system and compute the condition number
721* and error bounds using CGESVXX.
722*
723 srnamt = 'CGESVXX'
724 n_err_bnds = 3
725 CALL cgesvxx( fact, trans, n, nrhs, a, lda, afac,
726 $ lda, iwork, equed, s, s( n+1 ), b, lda, x,
727 $ lda, rcond, rpvgrw_svxx, berr, n_err_bnds,
728 $ errbnds_n, errbnds_c, 0, zero, work,
729 $ rwork, info )
730*
731* Check the error code from CGESVXX.
732*
733 IF( info.EQ.n+1 ) GOTO 50
734 IF( info.NE.izero ) THEN
735 CALL alaerh( path, 'CGESVXX', info, izero,
736 $ fact // trans, n, n, -1, -1, nrhs,
737 $ imat, nfail, nerrs, nout )
738 GOTO 50
739 END IF
740*
741* Compare rpvgrw_svxx from CGESVXX with the computed
742* reciprocal pivot growth factor RPVGRW
743*
744
745 IF ( info .GT. 0 .AND. info .LT. n+1 ) THEN
746 rpvgrw = cla_gerpvgrw
747 $ (n, info, a, lda, afac, lda)
748 ELSE
749 rpvgrw = cla_gerpvgrw
750 $ (n, n, a, lda, afac, lda)
751 ENDIF
752
753 result( 7 ) = abs( rpvgrw-rpvgrw_svxx ) /
754 $ max( rpvgrw_svxx, rpvgrw ) /
755 $ slamch( 'E' )
756*
757 IF( .NOT.prefac ) THEN
758*
759* Reconstruct matrix from factors and compute
760* residual.
761*
762 CALL cget01( n, n, a, lda, afac, lda, iwork,
763 $ rwork( 2*nrhs+1 ), result( 1 ) )
764 k1 = 1
765 ELSE
766 k1 = 2
767 END IF
768*
769 IF( info.EQ.0 ) THEN
770 trfcon = .false.
771*
772* Compute residual of the computed solution.
773*
774 CALL clacpy( 'Full', n, nrhs, bsav, lda, work,
775 $ lda )
776 CALL cget02( trans, n, n, nrhs, asav, lda, x,
777 $ lda, work, lda, rwork( 2*nrhs+1 ),
778 $ result( 2 ) )
779*
780* Check solution from generated exact solution.
781*
782 IF( nofact .OR. ( prefac .AND. lsame( equed,
783 $ 'N' ) ) ) THEN
784 CALL cget04( n, nrhs, x, lda, xact, lda,
785 $ rcondc, result( 3 ) )
786 ELSE
787 IF( itran.EQ.1 ) THEN
788 roldc = roldo
789 ELSE
790 roldc = roldi
791 END IF
792 CALL cget04( n, nrhs, x, lda, xact, lda,
793 $ roldc, result( 3 ) )
794 END IF
795 ELSE
796 trfcon = .true.
797 END IF
798*
799* Compare RCOND from CGESVXX with the computed value
800* in RCONDC.
801*
802 result( 6 ) = sget06( rcond, rcondc )
803*
804* Print information about the tests that did not pass
805* the threshold.
806*
807 IF( .NOT.trfcon ) THEN
808 DO 45 k = k1, ntests
809 IF( result( k ).GE.thresh ) THEN
810 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
811 $ CALL aladhd( nout, path )
812 IF( prefac ) THEN
813 WRITE( nout, fmt = 9997 )'CGESVXX',
814 $ fact, trans, n, equed, imat, k,
815 $ result( k )
816 ELSE
817 WRITE( nout, fmt = 9998 )'CGESVXX',
818 $ fact, trans, n, imat, k, result( k )
819 END IF
820 nfail = nfail + 1
821 END IF
822 45 CONTINUE
823 nrun = nrun + 7 - k1
824 ELSE
825 IF( result( 1 ).GE.thresh .AND. .NOT.prefac )
826 $ THEN
827 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
828 $ CALL aladhd( nout, path )
829 IF( prefac ) THEN
830 WRITE( nout, fmt = 9997 )'CGESVXX', fact,
831 $ trans, n, equed, imat, 1, result( 1 )
832 ELSE
833 WRITE( nout, fmt = 9998 )'CGESVXX', fact,
834 $ trans, n, imat, 1, result( 1 )
835 END IF
836 nfail = nfail + 1
837 nrun = nrun + 1
838 END IF
839 IF( result( 6 ).GE.thresh ) THEN
840 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
841 $ CALL aladhd( nout, path )
842 IF( prefac ) THEN
843 WRITE( nout, fmt = 9997 )'CGESVXX', fact,
844 $ trans, n, equed, imat, 6, result( 6 )
845 ELSE
846 WRITE( nout, fmt = 9998 )'CGESVXX', fact,
847 $ trans, n, imat, 6, result( 6 )
848 END IF
849 nfail = nfail + 1
850 nrun = nrun + 1
851 END IF
852 IF( result( 7 ).GE.thresh ) THEN
853 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
854 $ CALL aladhd( nout, path )
855 IF( prefac ) THEN
856 WRITE( nout, fmt = 9997 )'CGESVXX', fact,
857 $ trans, n, equed, imat, 7, result( 7 )
858 ELSE
859 WRITE( nout, fmt = 9998 )'CGESVXX', fact,
860 $ trans, n, imat, 7, result( 7 )
861 END IF
862 nfail = nfail + 1
863 nrun = nrun + 1
864 END IF
865*
866 END IF
867*
868 50 CONTINUE
869 60 CONTINUE
870 70 CONTINUE
871 80 CONTINUE
872 90 CONTINUE
873*
874* Print a summary of the results.
875*
876 CALL alasvm( path, nout, nfail, nrun, nerrs )
877*
878
879* Test Error Bounds for CGESVXX
880
881 CALL cebchvxx(thresh, path)
882
883 9999 FORMAT( 1x, a, ', N =', i5, ', type ', i2, ', test(', i2, ') =',
884 $ g12.5 )
885 9998 FORMAT( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
886 $ ', type ', i2, ', test(', i1, ')=', g12.5 )
887 9997 FORMAT( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
888 $ ', EQUED=''', a1, ''', type ', i2, ', test(', i1, ')=',
889 $ g12.5 )
890 RETURN
891*
892* End of CDRVGEX
893*
894 END
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
Definition alasvm.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 aladhd(iounit, path)
ALADHD
Definition aladhd.f:90
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
Definition alaerh.f:147
subroutine cdrvge(dotype, nn, nval, nrhs, thresh, tsterr, nmax, a, afac, asav, b, bsav, x, xact, s, work, rwork, iwork, nout)
CDRVGE
Definition cdrvge.f:164
subroutine cebchvxx(thresh, path)
CEBCHVXX
Definition cebchvxx.f:96
subroutine cerrvx(path, nunit)
CERRVX
Definition cerrvx.f:55
subroutine cget01(m, n, a, lda, afac, ldafac, ipiv, rwork, resid)
CGET01
Definition cget01.f:108
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 cgeequ(m, n, a, lda, r, c, rowcnd, colcnd, amax, info)
CGEEQU
Definition cgeequ.f:140
subroutine cgesv(n, nrhs, a, lda, ipiv, b, ldb, info)
Download CGESV + dependencies <a href="http://www.netlib.org/cgi-bin/netlibfiles....
Definition cgesv.f:124
subroutine cgesvx(fact, trans, n, nrhs, a, lda, af, ldaf, ipiv, equed, r, c, b, ldb, x, ldx, rcond, ferr, berr, work, rwork, info)
CGESVX computes the solution to system of linear equations A * X = B for GE matrices
Definition cgesvx.f:350
subroutine cgesvxx(fact, trans, n, nrhs, a, lda, af, ldaf, ipiv, equed, r, c, b, ldb, x, ldx, rcond, rpvgrw, berr, n_err_bnds, err_bnds_norm, err_bnds_comp, nparams, params, work, rwork, info)
CGESVXX computes the solution to system of linear equations A * X = B for GE matrices
Definition cgesvxx.f:543
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 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 claqge(m, n, a, lda, r, c, rowcnd, colcnd, amax, equed)
CLAQGE scales a general rectangular matrix, using row and column scaling factors computed by sgeequ.
Definition claqge.f:143
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