LAPACK  3.10.1
LAPACK: Linear Algebra PACKage
ddrvgex.f
Go to the documentation of this file.
1 *> \brief \b DDRVGEX
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 DDRVGE( 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 * DOUBLE PRECISION THRESH
19 * ..
20 * .. Array Arguments ..
21 * LOGICAL DOTYPE( * )
22 * INTEGER IWORK( * ), NVAL( * )
23 * DOUBLE PRECISION A( * ), AFAC( * ), ASAV( * ), B( * ),
24 * $ BSAV( * ), RWORK( * ), S( * ), WORK( * ),
25 * $ X( * ), XACT( * )
26 * ..
27 *
28 *
29 *> \par Purpose:
30 * =============
31 *>
32 *> \verbatim
33 *>
34 *> DDRVGE tests the driver routines DGESV, -SVX, and -SVXX.
35 *>
36 *> Note that this file is used only when the XBLAS are available,
37 *> otherwise ddrvge.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 DOUBLE PRECISION
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 DOUBLE PRECISION array, dimension (NMAX*NMAX)
94 *> \endverbatim
95 *>
96 *> \param[out] AFAC
97 *> \verbatim
98 *> AFAC is DOUBLE PRECISION array, dimension (NMAX*NMAX)
99 *> \endverbatim
100 *>
101 *> \param[out] ASAV
102 *> \verbatim
103 *> ASAV is DOUBLE PRECISION array, dimension (NMAX*NMAX)
104 *> \endverbatim
105 *>
106 *> \param[out] B
107 *> \verbatim
108 *> B is DOUBLE PRECISION array, dimension (NMAX*NRHS)
109 *> \endverbatim
110 *>
111 *> \param[out] BSAV
112 *> \verbatim
113 *> BSAV is DOUBLE PRECISION array, dimension (NMAX*NRHS)
114 *> \endverbatim
115 *>
116 *> \param[out] X
117 *> \verbatim
118 *> X is DOUBLE PRECISION array, dimension (NMAX*NRHS)
119 *> \endverbatim
120 *>
121 *> \param[out] XACT
122 *> \verbatim
123 *> XACT is DOUBLE PRECISION array, dimension (NMAX*NRHS)
124 *> \endverbatim
125 *>
126 *> \param[out] S
127 *> \verbatim
128 *> S is DOUBLE PRECISION array, dimension (2*NMAX)
129 *> \endverbatim
130 *>
131 *> \param[out] WORK
132 *> \verbatim
133 *> WORK is DOUBLE PRECISION array, dimension
134 *> (NMAX*max(3,NRHS))
135 *> \endverbatim
136 *>
137 *> \param[out] RWORK
138 *> \verbatim
139 *> RWORK is DOUBLE PRECISION array, dimension (2*NRHS+NMAX)
140 *> \endverbatim
141 *>
142 *> \param[out] IWORK
143 *> \verbatim
144 *> IWORK is INTEGER array, dimension (2*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 double_lin
162 *
163 * =====================================================================
164  SUBROUTINE ddrvge( 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  DOUBLE PRECISION THRESH
176 * ..
177 * .. Array Arguments ..
178  LOGICAL DOTYPE( * )
179  INTEGER IWORK( * ), NVAL( * )
180  DOUBLE PRECISION A( * ), AFAC( * ), ASAV( * ), B( * ),
181  $ bsav( * ), rwork( * ), s( * ), work( * ),
182  $ x( * ), xact( * )
183 * ..
184 *
185 * =====================================================================
186 *
187 * .. Parameters ..
188  DOUBLE PRECISION ONE, ZERO
189  PARAMETER ( ONE = 1.0d+0, zero = 0.0d+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  DOUBLE PRECISION 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  DOUBLE PRECISION RESULT( NTESTS ), BERR( NRHS ),
213  $ errbnds_n( nrhs, 3 ), errbnds_c( nrhs, 3 )
214 * ..
215 * .. External Functions ..
216  LOGICAL LSAME
217  DOUBLE PRECISION DGET06, DLAMCH, DLANGE, DLANTR, DLA_GERPVGRW
218  EXTERNAL lsame, dget06, dlamch, dlange, dlantr,
219  $ dla_gerpvgrw
220 * ..
221 * .. External Subroutines ..
222  EXTERNAL aladhd, alaerh, alasvm, derrvx, dgeequ, dgesv,
225  $ dlatms, xlaenv, dgesvxx
226 * ..
227 * .. Intrinsic Functions ..
228  INTRINSIC abs, 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 ) = 'Double 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 derrvx( 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 DLATB4 and generate a test matrix
295 * with DLATMS.
296 *
297  CALL dlatb4( path, imat, n, n, TYPE, KL, KU, ANORM, MODE,
298  $ CNDNUM, DIST )
299  rcondc = one / cndnum
300 *
301  srnamt = 'DLATMS'
302  CALL dlatms( 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 DLATMS.
307 *
308  IF( info.NE.0 ) THEN
309  CALL alaerh( path, 'DLATMS', 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 dlaset( 'Full', n, n-izero+1, zero, zero,
332  $ 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 dlacpy( '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 DGESVX (FACT = 'N' reuses
366 * the condition number from the previous iteration
367 * with FACT = 'F').
368 *
369  CALL dlacpy( '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 dgeequ( 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 dlaqge( 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 DGET04.
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 = dlange( '1', n, n, afac, lda, rwork )
407  anormi = dlange( 'I', n, n, afac, lda, rwork )
408 *
409 * Factor the matrix A.
410 *
411  CALL dgetrf( n, n, afac, lda, iwork, info )
412 *
413 * Form the inverse of A.
414 *
415  CALL dlacpy( 'Full', n, n, afac, lda, a, lda )
416  lwork = nmax*max( 3, nrhs )
417  CALL dgetri( n, a, lda, iwork, work, lwork, info )
418 *
419 * Compute the 1-norm condition number of A.
420 *
421  ainvnm = dlange( '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 = dlange( '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 dlacpy( 'Full', n, n, asav, lda, a, lda )
452 *
453 * Form an exact solution and set the right hand side.
454 *
455  srnamt = 'DLARHS'
456  CALL dlarhs( path, xtype, 'Full', trans, n, n, kl,
457  $ ku, nrhs, a, lda, xact, lda, b, lda,
458  $ iseed, info )
459  xtype = 'C'
460  CALL dlacpy( 'Full', n, nrhs, b, lda, bsav, lda )
461 *
462  IF( nofact .AND. itran.EQ.1 ) THEN
463 *
464 * --- Test DGESV ---
465 *
466 * Compute the LU factorization of the matrix and
467 * solve the system.
468 *
469  CALL dlacpy( 'Full', n, n, a, lda, afac, lda )
470  CALL dlacpy( 'Full', n, nrhs, b, lda, x, lda )
471 *
472  srnamt = 'DGESV '
473  CALL dgesv( n, nrhs, afac, lda, iwork, x, lda,
474  $ info )
475 *
476 * Check error code from DGESV .
477 *
478  IF( info.NE.izero )
479  $ CALL alaerh( path, 'DGESV ', 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 dget01( 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 dlacpy( 'Full', n, nrhs, b, lda, work,
494  $ lda )
495  CALL dget02( '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 dget04( 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 )'DGESV ', 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 DGESVX ---
522 *
523  IF( .NOT.prefac )
524  $ CALL dlaset( 'Full', n, n, zero, zero, afac,
525  $ lda )
526  CALL dlaset( 'Full', n, nrhs, zero, zero, x, lda )
527  IF( iequed.GT.1 .AND. n.GT.0 ) THEN
528 *
529 * Equilibrate the matrix if FACT = 'F' and
530 * EQUED = 'R', 'C', or 'B'.
531 *
532  CALL dlaqge( n, n, a, lda, s, s( n+1 ), rowcnd,
533  $ colcnd, amax, equed )
534  END IF
535 *
536 * Solve the system and compute the condition number
537 * and error bounds using DGESVX.
538 *
539  srnamt = 'DGESVX'
540  CALL dgesvx( fact, trans, n, nrhs, a, lda, afac,
541  $ lda, iwork, equed, s, s( n+1 ), b,
542  $ lda, x, lda, rcond, rwork,
543  $ rwork( nrhs+1 ), work, iwork( n+1 ),
544  $ info )
545 *
546 * Check the error code from DGESVX.
547 *
548  IF( info.NE.izero )
549  $ CALL alaerh( path, 'DGESVX', info, izero,
550  $ fact // trans, n, n, -1, -1, nrhs,
551  $ imat, nfail, nerrs, nout )
552 *
553 * Compare WORK(1) from DGESVX with the computed
554 * reciprocal pivot growth factor RPVGRW
555 *
556  IF( info.NE.0 ) THEN
557  rpvgrw = dlantr( 'M', 'U', 'N', info, info,
558  $ afac, lda, work )
559  IF( rpvgrw.EQ.zero ) THEN
560  rpvgrw = one
561  ELSE
562  rpvgrw = dlange( 'M', n, info, a, lda,
563  $ work ) / rpvgrw
564  END IF
565  ELSE
566  rpvgrw = dlantr( 'M', 'U', 'N', n, n, afac, lda,
567  $ work )
568  IF( rpvgrw.EQ.zero ) THEN
569  rpvgrw = one
570  ELSE
571  rpvgrw = dlange( 'M', n, n, a, lda, work ) /
572  $ rpvgrw
573  END IF
574  END IF
575  result( 7 ) = abs( rpvgrw-work( 1 ) ) /
576  $ max( work( 1 ), rpvgrw ) /
577  $ dlamch( 'E' )
578 *
579  IF( .NOT.prefac ) THEN
580 *
581 * Reconstruct matrix from factors and compute
582 * residual.
583 *
584  CALL dget01( n, n, a, lda, afac, lda, iwork,
585  $ rwork( 2*nrhs+1 ), result( 1 ) )
586  k1 = 1
587  ELSE
588  k1 = 2
589  END IF
590 *
591  IF( info.EQ.0 ) THEN
592  trfcon = .false.
593 *
594 * Compute residual of the computed solution.
595 *
596  CALL dlacpy( 'Full', n, nrhs, bsav, lda, work,
597  $ lda )
598  CALL dget02( trans, n, n, nrhs, asav, lda, x,
599  $ lda, work, lda, rwork( 2*nrhs+1 ),
600  $ result( 2 ) )
601 *
602 * Check solution from generated exact solution.
603 *
604  IF( nofact .OR. ( prefac .AND. lsame( equed,
605  $ 'N' ) ) ) THEN
606  CALL dget04( n, nrhs, x, lda, xact, lda,
607  $ rcondc, result( 3 ) )
608  ELSE
609  IF( itran.EQ.1 ) THEN
610  roldc = roldo
611  ELSE
612  roldc = roldi
613  END IF
614  CALL dget04( n, nrhs, x, lda, xact, lda,
615  $ roldc, result( 3 ) )
616  END IF
617 *
618 * Check the error bounds from iterative
619 * refinement.
620 *
621  CALL dget07( trans, n, nrhs, asav, lda, b, lda,
622  $ x, lda, xact, lda, rwork, .true.,
623  $ rwork( nrhs+1 ), result( 4 ) )
624  ELSE
625  trfcon = .true.
626  END IF
627 *
628 * Compare RCOND from DGESVX with the computed value
629 * in RCONDC.
630 *
631  result( 6 ) = dget06( rcond, rcondc )
632 *
633 * Print information about the tests that did not pass
634 * the threshold.
635 *
636  IF( .NOT.trfcon ) THEN
637  DO 40 k = k1, ntests
638  IF( result( k ).GE.thresh ) THEN
639  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
640  $ CALL aladhd( nout, path )
641  IF( prefac ) THEN
642  WRITE( nout, fmt = 9997 )'DGESVX',
643  $ fact, trans, n, equed, imat, k,
644  $ result( k )
645  ELSE
646  WRITE( nout, fmt = 9998 )'DGESVX',
647  $ fact, trans, n, imat, k, result( k )
648  END IF
649  nfail = nfail + 1
650  END IF
651  40 CONTINUE
652  nrun = nrun + 7 - k1
653  ELSE
654  IF( result( 1 ).GE.thresh .AND. .NOT.prefac )
655  $ THEN
656  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
657  $ CALL aladhd( nout, path )
658  IF( prefac ) THEN
659  WRITE( nout, fmt = 9997 )'DGESVX', fact,
660  $ trans, n, equed, imat, 1, result( 1 )
661  ELSE
662  WRITE( nout, fmt = 9998 )'DGESVX', fact,
663  $ trans, n, imat, 1, result( 1 )
664  END IF
665  nfail = nfail + 1
666  nrun = nrun + 1
667  END IF
668  IF( result( 6 ).GE.thresh ) THEN
669  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
670  $ CALL aladhd( nout, path )
671  IF( prefac ) THEN
672  WRITE( nout, fmt = 9997 )'DGESVX', fact,
673  $ trans, n, equed, imat, 6, result( 6 )
674  ELSE
675  WRITE( nout, fmt = 9998 )'DGESVX', fact,
676  $ trans, n, imat, 6, result( 6 )
677  END IF
678  nfail = nfail + 1
679  nrun = nrun + 1
680  END IF
681  IF( result( 7 ).GE.thresh ) THEN
682  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
683  $ CALL aladhd( nout, path )
684  IF( prefac ) THEN
685  WRITE( nout, fmt = 9997 )'DGESVX', fact,
686  $ trans, n, equed, imat, 7, result( 7 )
687  ELSE
688  WRITE( nout, fmt = 9998 )'DGESVX', fact,
689  $ trans, n, imat, 7, result( 7 )
690  END IF
691  nfail = nfail + 1
692  nrun = nrun + 1
693  END IF
694 *
695  END IF
696 *
697 * --- Test DGESVXX ---
698 *
699 * Restore the matrices A and B.
700 *
701  CALL dlacpy( 'Full', n, n, asav, lda, a, lda )
702  CALL dlacpy( 'Full', n, nrhs, bsav, lda, b, lda )
703 
704  IF( .NOT.prefac )
705  $ CALL dlaset( 'Full', n, n, zero, zero, afac,
706  $ lda )
707  CALL dlaset( 'Full', n, nrhs, zero, zero, x, lda )
708  IF( iequed.GT.1 .AND. n.GT.0 ) THEN
709 *
710 * Equilibrate the matrix if FACT = 'F' and
711 * EQUED = 'R', 'C', or 'B'.
712 *
713  CALL dlaqge( n, n, a, lda, s, s( n+1 ), rowcnd,
714  $ colcnd, amax, equed )
715  END IF
716 *
717 * Solve the system and compute the condition number
718 * and error bounds using DGESVXX.
719 *
720  srnamt = 'DGESVXX'
721  n_err_bnds = 3
722  CALL dgesvxx( fact, trans, n, nrhs, a, lda, afac,
723  $ lda, iwork, equed, s, s( n+1 ), b, lda, x,
724  $ lda, rcond, rpvgrw_svxx, berr, n_err_bnds,
725  $ errbnds_n, errbnds_c, 0, zero, work,
726  $ iwork( n+1 ), info )
727 *
728 * Check the error code from DGESVXX.
729 *
730  IF( info.EQ.n+1 ) GOTO 50
731  IF( info.NE.izero ) THEN
732  CALL alaerh( path, 'DGESVXX', info, izero,
733  $ fact // trans, n, n, -1, -1, nrhs,
734  $ imat, nfail, nerrs, nout )
735  GOTO 50
736  END IF
737 *
738 * Compare rpvgrw_svxx from DGESVXX with the computed
739 * reciprocal pivot growth factor RPVGRW
740 *
741 
742  IF ( info .GT. 0 .AND. info .LT. n+1 ) THEN
743  rpvgrw = dla_gerpvgrw
744  $ (n, info, a, lda, afac, lda)
745  ELSE
746  rpvgrw = dla_gerpvgrw
747  $ (n, n, a, lda, afac, lda)
748  ENDIF
749 
750  result( 7 ) = abs( rpvgrw-rpvgrw_svxx ) /
751  $ max( rpvgrw_svxx, rpvgrw ) /
752  $ dlamch( 'E' )
753 *
754  IF( .NOT.prefac ) THEN
755 *
756 * Reconstruct matrix from factors and compute
757 * residual.
758 *
759  CALL dget01( n, n, a, lda, afac, lda, iwork,
760  $ rwork( 2*nrhs+1 ), result( 1 ) )
761  k1 = 1
762  ELSE
763  k1 = 2
764  END IF
765 *
766  IF( info.EQ.0 ) THEN
767  trfcon = .false.
768 *
769 * Compute residual of the computed solution.
770 *
771  CALL dlacpy( 'Full', n, nrhs, bsav, lda, work,
772  $ lda )
773  CALL dget02( trans, n, n, nrhs, asav, lda, x,
774  $ lda, work, lda, rwork( 2*nrhs+1 ),
775  $ result( 2 ) )
776 *
777 * Check solution from generated exact solution.
778 *
779  IF( nofact .OR. ( prefac .AND. lsame( equed,
780  $ 'N' ) ) ) THEN
781  CALL dget04( n, nrhs, x, lda, xact, lda,
782  $ rcondc, result( 3 ) )
783  ELSE
784  IF( itran.EQ.1 ) THEN
785  roldc = roldo
786  ELSE
787  roldc = roldi
788  END IF
789  CALL dget04( n, nrhs, x, lda, xact, lda,
790  $ roldc, result( 3 ) )
791  END IF
792  ELSE
793  trfcon = .true.
794  END IF
795 *
796 * Compare RCOND from DGESVXX with the computed value
797 * in RCONDC.
798 *
799  result( 6 ) = dget06( rcond, rcondc )
800 *
801 * Print information about the tests that did not pass
802 * the threshold.
803 *
804  IF( .NOT.trfcon ) THEN
805  DO 45 k = k1, ntests
806  IF( result( k ).GE.thresh ) THEN
807  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
808  $ CALL aladhd( nout, path )
809  IF( prefac ) THEN
810  WRITE( nout, fmt = 9997 )'DGESVXX',
811  $ fact, trans, n, equed, imat, k,
812  $ result( k )
813  ELSE
814  WRITE( nout, fmt = 9998 )'DGESVXX',
815  $ fact, trans, n, imat, k, result( k )
816  END IF
817  nfail = nfail + 1
818  END IF
819  45 CONTINUE
820  nrun = nrun + 7 - k1
821  ELSE
822  IF( result( 1 ).GE.thresh .AND. .NOT.prefac )
823  $ THEN
824  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
825  $ CALL aladhd( nout, path )
826  IF( prefac ) THEN
827  WRITE( nout, fmt = 9997 )'DGESVXX', fact,
828  $ trans, n, equed, imat, 1, result( 1 )
829  ELSE
830  WRITE( nout, fmt = 9998 )'DGESVXX', fact,
831  $ trans, n, imat, 1, result( 1 )
832  END IF
833  nfail = nfail + 1
834  nrun = nrun + 1
835  END IF
836  IF( result( 6 ).GE.thresh ) THEN
837  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
838  $ CALL aladhd( nout, path )
839  IF( prefac ) THEN
840  WRITE( nout, fmt = 9997 )'DGESVXX', fact,
841  $ trans, n, equed, imat, 6, result( 6 )
842  ELSE
843  WRITE( nout, fmt = 9998 )'DGESVXX', fact,
844  $ trans, n, imat, 6, result( 6 )
845  END IF
846  nfail = nfail + 1
847  nrun = nrun + 1
848  END IF
849  IF( result( 7 ).GE.thresh ) THEN
850  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
851  $ CALL aladhd( nout, path )
852  IF( prefac ) THEN
853  WRITE( nout, fmt = 9997 )'DGESVXX', fact,
854  $ trans, n, equed, imat, 7, result( 7 )
855  ELSE
856  WRITE( nout, fmt = 9998 )'DGESVXX', fact,
857  $ trans, n, imat, 7, result( 7 )
858  END IF
859  nfail = nfail + 1
860  nrun = nrun + 1
861  END IF
862 *
863  END IF
864 *
865  50 CONTINUE
866  60 CONTINUE
867  70 CONTINUE
868  80 CONTINUE
869  90 CONTINUE
870 *
871 * Print a summary of the results.
872 *
873  CALL alasvm( path, nout, nfail, nrun, nerrs )
874 *
875 
876 * Test Error Bounds from DGESVXX
877 
878  CALL debchvxx( thresh, path )
879 
880  9999 FORMAT( 1x, a, ', N =', i5, ', type ', i2, ', test(', i2, ') =',
881  $ g12.5 )
882  9998 FORMAT( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
883  $ ', type ', i2, ', test(', i1, ')=', g12.5 )
884  9997 FORMAT( 1x, a, ', FACT=''', a1, ''', TRANS=''', a1, ''', N=', i5,
885  $ ', EQUED=''', a1, ''', type ', i2, ', test(', i1, ')=',
886  $ g12.5 )
887  RETURN
888 *
889 * End of DDRVGEX
890 *
891  END
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: dlaset.f:110
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:73
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 dlarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
DLARHS
Definition: dlarhs.f:205
subroutine dget02(TRANS, M, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
DGET02
Definition: dget02.f:135
subroutine dget07(TRANS, N, NRHS, A, LDA, B, LDB, X, LDX, XACT, LDXACT, FERR, CHKFERR, BERR, RESLTS)
DGET07
Definition: dget07.f:165
subroutine dget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
DGET04
Definition: dget04.f:102
subroutine dget01(M, N, A, LDA, AFAC, LDAFAC, IPIV, RWORK, RESID)
DGET01
Definition: dget01.f:107
subroutine dlatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
DLATB4
Definition: dlatb4.f:120
subroutine derrvx(PATH, NUNIT)
DERRVX
Definition: derrvx.f:55
subroutine debchvxx(THRESH, PATH)
DEBCHVXX
Definition: debchvxx.f:96
subroutine ddrvge(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK, RWORK, IWORK, NOUT)
DDRVGE
Definition: ddrvge.f:164
subroutine dlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
DLATMS
Definition: dlatms.f:321
subroutine dlaqge(M, N, A, LDA, R, C, ROWCND, COLCND, AMAX, EQUED)
DLAQGE scales a general rectangular matrix, using row and column scaling factors computed by sgeequ.
Definition: dlaqge.f:142
subroutine dgetrf(M, N, A, LDA, IPIV, INFO)
DGETRF
Definition: dgetrf.f:108
subroutine dgeequ(M, N, A, LDA, R, C, ROWCND, COLCND, AMAX, INFO)
DGEEQU
Definition: dgeequ.f:139
subroutine dgetri(N, A, LDA, IPIV, WORK, LWORK, INFO)
DGETRI
Definition: dgetri.f:114
subroutine dgesvxx(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, IWORK, INFO)
DGESVXX computes the solution to system of linear equations A * X = B for GE matrices
Definition: dgesvxx.f:540
subroutine dgesv(N, NRHS, A, LDA, IPIV, B, LDB, INFO)
DGESV computes the solution to system of linear equations A * X = B for GE matrices
Definition: dgesv.f:122
subroutine dgesvx(FACT, TRANS, N, NRHS, A, LDA, AF, LDAF, IPIV, EQUED, R, C, B, LDB, X, LDX, RCOND, FERR, BERR, WORK, IWORK, INFO)
DGESVX computes the solution to system of linear equations A * X = B for GE matrices
Definition: dgesvx.f:349