LAPACK  3.10.1
LAPACK: Linear Algebra PACKage
derrgex.f
Go to the documentation of this file.
1 *> \brief \b DERRGEX
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 DERRGE( PATH, NUNIT )
12 *
13 * .. Scalar Arguments ..
14 * CHARACTER*3 PATH
15 * INTEGER NUNIT
16 * ..
17 *
18 *
19 *> \par Purpose:
20 * =============
21 *>
22 *> \verbatim
23 *>
24 *> DERRGE tests the error exits for the DOUBLE PRECISION routines
25 *> for general matrices.
26 *>
27 *> Note that this file is used only when the XBLAS are available,
28 *> otherwise derrge.f defines this subroutine.
29 *> \endverbatim
30 *
31 * Arguments:
32 * ==========
33 *
34 *> \param[in] PATH
35 *> \verbatim
36 *> PATH is CHARACTER*3
37 *> The LAPACK path name for the routines to be tested.
38 *> \endverbatim
39 *>
40 *> \param[in] NUNIT
41 *> \verbatim
42 *> NUNIT is INTEGER
43 *> The unit number for output.
44 *> \endverbatim
45 *
46 * Authors:
47 * ========
48 *
49 *> \author Univ. of Tennessee
50 *> \author Univ. of California Berkeley
51 *> \author Univ. of Colorado Denver
52 *> \author NAG Ltd.
53 *
54 *> \ingroup double_lin
55 *
56 * =====================================================================
57  SUBROUTINE derrge( PATH, NUNIT )
58 *
59 * -- LAPACK test routine --
60 * -- LAPACK is a software package provided by Univ. of Tennessee, --
61 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
62 *
63 * .. Scalar Arguments ..
64  CHARACTER*3 PATH
65  INTEGER NUNIT
66 * ..
67 *
68 * =====================================================================
69 *
70 * .. Parameters ..
71  INTEGER NMAX, LW
72  parameter( nmax = 4, lw = 3*nmax )
73 * ..
74 * .. Local Scalars ..
75  CHARACTER EQ
76  CHARACTER*2 C2
77  INTEGER I, INFO, J, N_ERR_BNDS, NPARAMS
78  DOUBLE PRECISION ANRM, CCOND, RCOND, BERR
79 * ..
80 * .. Local Arrays ..
81  INTEGER IP( NMAX ), IW( NMAX )
82  DOUBLE PRECISION A( NMAX, NMAX ), AF( NMAX, NMAX ), B( NMAX ),
83  $ C( NMAX ), R( NMAX ), R1( NMAX ), R2( NMAX ),
84  $ W( LW ), X( NMAX ), ERR_BNDS_N( NMAX, 3 ),
85  $ ERR_BNDS_C( NMAX, 3 ), PARAMS( 1 )
86 * ..
87 * .. External Functions ..
88  LOGICAL LSAMEN
89  EXTERNAL lsamen
90 * ..
91 * .. External Subroutines ..
92  EXTERNAL alaesm, chkxer, dgbcon, dgbequ, dgbrfs, dgbtf2,
95  $ dgbequb, dgbrfsx
96 * ..
97 * .. Scalars in Common ..
98  LOGICAL LERR, OK
99  CHARACTER*32 SRNAMT
100  INTEGER INFOT, NOUT
101 * ..
102 * .. Common blocks ..
103  COMMON / infoc / infot, nout, ok, lerr
104  COMMON / srnamc / srnamt
105 * ..
106 * .. Intrinsic Functions ..
107  INTRINSIC dble
108 * ..
109 * .. Executable Statements ..
110 *
111  nout = nunit
112  WRITE( nout, fmt = * )
113  c2 = path( 2: 3 )
114 *
115 * Set the variables to innocuous values.
116 *
117  DO 20 j = 1, nmax
118  DO 10 i = 1, nmax
119  a( i, j ) = 1.d0 / dble( i+j )
120  af( i, j ) = 1.d0 / dble( i+j )
121  10 CONTINUE
122  b( j ) = 0.d0
123  r1( j ) = 0.d0
124  r2( j ) = 0.d0
125  w( j ) = 0.d0
126  x( j ) = 0.d0
127  c( j ) = 0.d0
128  r( j ) = 0.d0
129  ip( j ) = j
130  iw( j ) = j
131  20 CONTINUE
132  ok = .true.
133 *
134  IF( lsamen( 2, c2, 'GE' ) ) THEN
135 *
136 * Test error exits of the routines that use the LU decomposition
137 * of a general matrix.
138 *
139 * DGETRF
140 *
141  srnamt = 'DGETRF'
142  infot = 1
143  CALL dgetrf( -1, 0, a, 1, ip, info )
144  CALL chkxer( 'DGETRF', infot, nout, lerr, ok )
145  infot = 2
146  CALL dgetrf( 0, -1, a, 1, ip, info )
147  CALL chkxer( 'DGETRF', infot, nout, lerr, ok )
148  infot = 4
149  CALL dgetrf( 2, 1, a, 1, ip, info )
150  CALL chkxer( 'DGETRF', infot, nout, lerr, ok )
151 *
152 * DGETF2
153 *
154  srnamt = 'DGETF2'
155  infot = 1
156  CALL dgetf2( -1, 0, a, 1, ip, info )
157  CALL chkxer( 'DGETF2', infot, nout, lerr, ok )
158  infot = 2
159  CALL dgetf2( 0, -1, a, 1, ip, info )
160  CALL chkxer( 'DGETF2', infot, nout, lerr, ok )
161  infot = 4
162  CALL dgetf2( 2, 1, a, 1, ip, info )
163  CALL chkxer( 'DGETF2', infot, nout, lerr, ok )
164 *
165 * DGETRI
166 *
167  srnamt = 'DGETRI'
168  infot = 1
169  CALL dgetri( -1, a, 1, ip, w, lw, info )
170  CALL chkxer( 'DGETRI', infot, nout, lerr, ok )
171  infot = 3
172  CALL dgetri( 2, a, 1, ip, w, lw, info )
173  CALL chkxer( 'DGETRI', infot, nout, lerr, ok )
174 *
175 * DGETRS
176 *
177  srnamt = 'DGETRS'
178  infot = 1
179  CALL dgetrs( '/', 0, 0, a, 1, ip, b, 1, info )
180  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
181  infot = 2
182  CALL dgetrs( 'N', -1, 0, a, 1, ip, b, 1, info )
183  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
184  infot = 3
185  CALL dgetrs( 'N', 0, -1, a, 1, ip, b, 1, info )
186  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
187  infot = 5
188  CALL dgetrs( 'N', 2, 1, a, 1, ip, b, 2, info )
189  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
190  infot = 8
191  CALL dgetrs( 'N', 2, 1, a, 2, ip, b, 1, info )
192  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
193 *
194 * DGERFS
195 *
196  srnamt = 'DGERFS'
197  infot = 1
198  CALL dgerfs( '/', 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1, r2, w,
199  $ iw, info )
200  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
201  infot = 2
202  CALL dgerfs( 'N', -1, 0, a, 1, af, 1, ip, b, 1, x, 1, r1, r2,
203  $ w, iw, info )
204  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
205  infot = 3
206  CALL dgerfs( 'N', 0, -1, a, 1, af, 1, ip, b, 1, x, 1, r1, r2,
207  $ w, iw, info )
208  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
209  infot = 5
210  CALL dgerfs( 'N', 2, 1, a, 1, af, 2, ip, b, 2, x, 2, r1, r2, w,
211  $ iw, info )
212  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
213  infot = 7
214  CALL dgerfs( 'N', 2, 1, a, 2, af, 1, ip, b, 2, x, 2, r1, r2, w,
215  $ iw, info )
216  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
217  infot = 10
218  CALL dgerfs( 'N', 2, 1, a, 2, af, 2, ip, b, 1, x, 2, r1, r2, w,
219  $ iw, info )
220  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
221  infot = 12
222  CALL dgerfs( 'N', 2, 1, a, 2, af, 2, ip, b, 2, x, 1, r1, r2, w,
223  $ iw, info )
224  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
225 *
226 * DGERFSX
227 *
228  n_err_bnds = 3
229  nparams = 0
230  srnamt = 'DGERFSX'
231  infot = 1
232  CALL dgerfsx( '/', eq, 0, 0, a, 1, af, 1, ip, r, c, b, 1, x,
233  $ 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
234  $ nparams, params, w, iw, info )
235  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
236  infot = 2
237  eq = '/'
238  CALL dgerfsx( 'N', eq, 2, 1, a, 1, af, 2, ip, r, c, b, 2, x,
239  $ 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
240  $ nparams, params, w, iw, info )
241  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
242  infot = 3
243  eq = 'R'
244  CALL dgerfsx( 'N', eq, -1, 0, a, 1, af, 1, ip, r, c, b, 1, x,
245  $ 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
246  $ nparams, params, w, iw, info )
247  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
248  infot = 4
249  CALL dgerfsx( 'N', eq, 0, -1, a, 1, af, 1, ip, r, c, b, 1, x,
250  $ 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
251  $ nparams, params, w, iw, info )
252  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
253  infot = 6
254  CALL dgerfsx( 'N', eq, 2, 1, a, 1, af, 2, ip, r, c, b, 2, x,
255  $ 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
256  $ nparams, params, w, iw, info )
257  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
258  infot = 8
259  CALL dgerfsx( 'N', eq, 2, 1, a, 2, af, 1, ip, r, c, b, 2, x,
260  $ 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
261  $ nparams, params, w, iw, info )
262  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
263  infot = 13
264  eq = 'C'
265  CALL dgerfsx( 'N', eq, 2, 1, a, 2, af, 2, ip, r, c, b, 1, x,
266  $ 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
267  $ nparams, params, w, iw, info )
268  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
269  infot = 15
270  CALL dgerfsx( 'N', eq, 2, 1, a, 2, af, 2, ip, r, c, b, 2, x,
271  $ 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
272  $ nparams, params, w, iw, info )
273  CALL chkxer( 'DGERFSX', infot, nout, lerr, ok )
274 *
275 * DGECON
276 *
277  srnamt = 'DGECON'
278  infot = 1
279  CALL dgecon( '/', 0, a, 1, anrm, rcond, w, iw, info )
280  CALL chkxer( 'DGECON', infot, nout, lerr, ok )
281  infot = 2
282  CALL dgecon( '1', -1, a, 1, anrm, rcond, w, iw, info )
283  CALL chkxer( 'DGECON', infot, nout, lerr, ok )
284  infot = 4
285  CALL dgecon( '1', 2, a, 1, anrm, rcond, w, iw, info )
286  CALL chkxer( 'DGECON', infot, nout, lerr, ok )
287 *
288 * DGEEQU
289 *
290  srnamt = 'DGEEQU'
291  infot = 1
292  CALL dgeequ( -1, 0, a, 1, r1, r2, rcond, ccond, anrm, info )
293  CALL chkxer( 'DGEEQU', infot, nout, lerr, ok )
294  infot = 2
295  CALL dgeequ( 0, -1, a, 1, r1, r2, rcond, ccond, anrm, info )
296  CALL chkxer( 'DGEEQU', infot, nout, lerr, ok )
297  infot = 4
298  CALL dgeequ( 2, 2, a, 1, r1, r2, rcond, ccond, anrm, info )
299  CALL chkxer( 'DGEEQU', infot, nout, lerr, ok )
300 *
301 * DGEEQUB
302 *
303  srnamt = 'DGEEQUB'
304  infot = 1
305  CALL dgeequb( -1, 0, a, 1, r1, r2, rcond, ccond, anrm, info )
306  CALL chkxer( 'DGEEQUB', infot, nout, lerr, ok )
307  infot = 2
308  CALL dgeequb( 0, -1, a, 1, r1, r2, rcond, ccond, anrm, info )
309  CALL chkxer( 'DGEEQUB', infot, nout, lerr, ok )
310  infot = 4
311  CALL dgeequb( 2, 2, a, 1, r1, r2, rcond, ccond, anrm, info )
312  CALL chkxer( 'DGEEQUB', infot, nout, lerr, ok )
313 *
314  ELSE IF( lsamen( 2, c2, 'GB' ) ) THEN
315 *
316 * Test error exits of the routines that use the LU decomposition
317 * of a general band matrix.
318 *
319 * DGBTRF
320 *
321  srnamt = 'DGBTRF'
322  infot = 1
323  CALL dgbtrf( -1, 0, 0, 0, a, 1, ip, info )
324  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
325  infot = 2
326  CALL dgbtrf( 0, -1, 0, 0, a, 1, ip, info )
327  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
328  infot = 3
329  CALL dgbtrf( 1, 1, -1, 0, a, 1, ip, info )
330  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
331  infot = 4
332  CALL dgbtrf( 1, 1, 0, -1, a, 1, ip, info )
333  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
334  infot = 6
335  CALL dgbtrf( 2, 2, 1, 1, a, 3, ip, info )
336  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
337 *
338 * DGBTF2
339 *
340  srnamt = 'DGBTF2'
341  infot = 1
342  CALL dgbtf2( -1, 0, 0, 0, a, 1, ip, info )
343  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
344  infot = 2
345  CALL dgbtf2( 0, -1, 0, 0, a, 1, ip, info )
346  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
347  infot = 3
348  CALL dgbtf2( 1, 1, -1, 0, a, 1, ip, info )
349  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
350  infot = 4
351  CALL dgbtf2( 1, 1, 0, -1, a, 1, ip, info )
352  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
353  infot = 6
354  CALL dgbtf2( 2, 2, 1, 1, a, 3, ip, info )
355  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
356 *
357 * DGBTRS
358 *
359  srnamt = 'DGBTRS'
360  infot = 1
361  CALL dgbtrs( '/', 0, 0, 0, 1, a, 1, ip, b, 1, info )
362  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
363  infot = 2
364  CALL dgbtrs( 'N', -1, 0, 0, 1, a, 1, ip, b, 1, info )
365  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
366  infot = 3
367  CALL dgbtrs( 'N', 1, -1, 0, 1, a, 1, ip, b, 1, info )
368  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
369  infot = 4
370  CALL dgbtrs( 'N', 1, 0, -1, 1, a, 1, ip, b, 1, info )
371  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
372  infot = 5
373  CALL dgbtrs( 'N', 1, 0, 0, -1, a, 1, ip, b, 1, info )
374  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
375  infot = 7
376  CALL dgbtrs( 'N', 2, 1, 1, 1, a, 3, ip, b, 2, info )
377  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
378  infot = 10
379  CALL dgbtrs( 'N', 2, 0, 0, 1, a, 1, ip, b, 1, info )
380  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
381 *
382 * DGBRFS
383 *
384  srnamt = 'DGBRFS'
385  infot = 1
386  CALL dgbrfs( '/', 0, 0, 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
387  $ r2, w, iw, info )
388  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
389  infot = 2
390  CALL dgbrfs( 'N', -1, 0, 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
391  $ r2, w, iw, info )
392  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
393  infot = 3
394  CALL dgbrfs( 'N', 1, -1, 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
395  $ r2, w, iw, info )
396  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
397  infot = 4
398  CALL dgbrfs( 'N', 1, 0, -1, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
399  $ r2, w, iw, info )
400  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
401  infot = 5
402  CALL dgbrfs( 'N', 1, 0, 0, -1, a, 1, af, 1, ip, b, 1, x, 1, r1,
403  $ r2, w, iw, info )
404  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
405  infot = 7
406  CALL dgbrfs( 'N', 2, 1, 1, 1, a, 2, af, 4, ip, b, 2, x, 2, r1,
407  $ r2, w, iw, info )
408  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
409  infot = 9
410  CALL dgbrfs( 'N', 2, 1, 1, 1, a, 3, af, 3, ip, b, 2, x, 2, r1,
411  $ r2, w, iw, info )
412  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
413  infot = 12
414  CALL dgbrfs( 'N', 2, 0, 0, 1, a, 1, af, 1, ip, b, 1, x, 2, r1,
415  $ r2, w, iw, info )
416  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
417  infot = 14
418  CALL dgbrfs( 'N', 2, 0, 0, 1, a, 1, af, 1, ip, b, 2, x, 1, r1,
419  $ r2, w, iw, info )
420  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
421 *
422 * DGBRFSX
423 *
424  n_err_bnds = 3
425  nparams = 0
426  srnamt = 'DGBRFSX'
427  infot = 1
428  CALL dgbrfsx( '/', eq, 0, 0, 0, 0, a, 1, af, 1, ip, r, c, b, 1,
429  $ x, 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
430  $ nparams, params, w, iw, info )
431  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
432  infot = 2
433  eq = '/'
434  CALL dgbrfsx( 'N', eq, 2, 1, 1, 1, a, 1, af, 2, ip, r, c, b, 2,
435  $ x, 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
436  $ nparams, params, w, iw, info )
437  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
438  infot = 3
439  eq = 'R'
440  CALL dgbrfsx( 'N', eq, -1, 1, 1, 0, a, 1, af, 1, ip, r, c, b,
441  $ 1, x, 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
442  $ nparams, params, w, iw, info )
443  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
444  infot = 4
445  eq = 'R'
446  CALL dgbrfsx( 'N', eq, 2, -1, 1, 1, a, 3, af, 4, ip, r, c, b,
447  $ 1, x, 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
448  $ nparams, params, w, iw, info )
449  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
450  infot = 5
451  eq = 'R'
452  CALL dgbrfsx( 'N', eq, 2, 1, -1, 1, a, 3, af, 4, ip, r, c, b,
453  $ 1, x, 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
454  $ nparams, params, w, iw, info )
455  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
456  infot = 6
457  CALL dgbrfsx( 'N', eq, 0, 0, 0, -1, a, 1, af, 1, ip, r, c, b,
458  $ 1, x, 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
459  $ nparams, params, w, iw, info )
460  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
461  infot = 8
462  CALL dgbrfsx( 'N', eq, 2, 1, 1, 1, a, 1, af, 2, ip, r, c, b,
463  $ 2, x, 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
464  $ nparams, params, w, iw, info )
465  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
466  infot = 10
467  CALL dgbrfsx( 'N', eq, 2, 1, 1, 1, a, 3, af, 3, ip, r, c, b, 2,
468  $ x, 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
469  $ nparams, params, w, iw, info )
470  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
471  infot = 13
472  eq = 'C'
473  CALL dgbrfsx( 'N', eq, 2, 1, 1, 1, a, 3, af, 5, ip, r, c, b,
474  $ 1, x, 2, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
475  $ nparams, params, w, iw, info )
476  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
477  infot = 15
478  CALL dgbrfsx( 'N', eq, 2, 1, 1, 1, a, 3, af, 5, ip, r, c, b, 2,
479  $ x, 1, rcond, berr, n_err_bnds, err_bnds_n, err_bnds_c,
480  $ nparams, params, w, iw, info )
481  CALL chkxer( 'DGBRFSX', infot, nout, lerr, ok )
482 *
483 * DGBCON
484 *
485  srnamt = 'DGBCON'
486  infot = 1
487  CALL dgbcon( '/', 0, 0, 0, a, 1, ip, anrm, rcond, w, iw, info )
488  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
489  infot = 2
490  CALL dgbcon( '1', -1, 0, 0, a, 1, ip, anrm, rcond, w, iw,
491  $ info )
492  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
493  infot = 3
494  CALL dgbcon( '1', 1, -1, 0, a, 1, ip, anrm, rcond, w, iw,
495  $ info )
496  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
497  infot = 4
498  CALL dgbcon( '1', 1, 0, -1, a, 1, ip, anrm, rcond, w, iw,
499  $ info )
500  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
501  infot = 6
502  CALL dgbcon( '1', 2, 1, 1, a, 3, ip, anrm, rcond, w, iw, info )
503  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
504 *
505 * DGBEQU
506 *
507  srnamt = 'DGBEQU'
508  infot = 1
509  CALL dgbequ( -1, 0, 0, 0, a, 1, r1, r2, rcond, ccond, anrm,
510  $ info )
511  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
512  infot = 2
513  CALL dgbequ( 0, -1, 0, 0, a, 1, r1, r2, rcond, ccond, anrm,
514  $ info )
515  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
516  infot = 3
517  CALL dgbequ( 1, 1, -1, 0, a, 1, r1, r2, rcond, ccond, anrm,
518  $ info )
519  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
520  infot = 4
521  CALL dgbequ( 1, 1, 0, -1, a, 1, r1, r2, rcond, ccond, anrm,
522  $ info )
523  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
524  infot = 6
525  CALL dgbequ( 2, 2, 1, 1, a, 2, r1, r2, rcond, ccond, anrm,
526  $ info )
527  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
528 *
529 * DGBEQUB
530 *
531  srnamt = 'DGBEQUB'
532  infot = 1
533  CALL dgbequb( -1, 0, 0, 0, a, 1, r1, r2, rcond, ccond, anrm,
534  $ info )
535  CALL chkxer( 'DGBEQUB', infot, nout, lerr, ok )
536  infot = 2
537  CALL dgbequb( 0, -1, 0, 0, a, 1, r1, r2, rcond, ccond, anrm,
538  $ info )
539  CALL chkxer( 'DGBEQUB', infot, nout, lerr, ok )
540  infot = 3
541  CALL dgbequb( 1, 1, -1, 0, a, 1, r1, r2, rcond, ccond, anrm,
542  $ info )
543  CALL chkxer( 'DGBEQUB', infot, nout, lerr, ok )
544  infot = 4
545  CALL dgbequb( 1, 1, 0, -1, a, 1, r1, r2, rcond, ccond, anrm,
546  $ info )
547  CALL chkxer( 'DGBEQUB', infot, nout, lerr, ok )
548  infot = 6
549  CALL dgbequb( 2, 2, 1, 1, a, 2, r1, r2, rcond, ccond, anrm,
550  $ info )
551  CALL chkxer( 'DGBEQUB', infot, nout, lerr, ok )
552  END IF
553 *
554 * Print a summary line.
555 *
556  CALL alaesm( path, ok, nout )
557 *
558  RETURN
559 *
560 * End of DERRGEX
561 *
562  END
subroutine chkxer(SRNAMT, INFOT, NOUT, LERR, OK)
Definition: cblat2.f:3196
subroutine alaesm(PATH, OK, NOUT)
ALAESM
Definition: alaesm.f:63
subroutine derrge(PATH, NUNIT)
DERRGE
Definition: derrge.f:55
subroutine dgbtrs(TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB, INFO)
DGBTRS
Definition: dgbtrs.f:138
subroutine dgbequb(M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND, AMAX, INFO)
DGBEQUB
Definition: dgbequb.f:160
subroutine dgbtrf(M, N, KL, KU, AB, LDAB, IPIV, INFO)
DGBTRF
Definition: dgbtrf.f:144
subroutine dgbequ(M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND, AMAX, INFO)
DGBEQU
Definition: dgbequ.f:153
subroutine dgbrfsx(TRANS, EQUED, N, KL, KU, NRHS, AB, LDAB, AFB, LDAFB, IPIV, R, C, B, LDB, X, LDX, RCOND, BERR, N_ERR_BNDS, ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, IWORK, INFO)
DGBRFSX
Definition: dgbrfsx.f:440
subroutine dgbcon(NORM, N, KL, KU, AB, LDAB, IPIV, ANORM, RCOND, WORK, IWORK, INFO)
DGBCON
Definition: dgbcon.f:146
subroutine dgbrfs(TRANS, N, KL, KU, NRHS, AB, LDAB, AFB, LDAFB, IPIV, B, LDB, X, LDX, FERR, BERR, WORK, IWORK, INFO)
DGBRFS
Definition: dgbrfs.f:205
subroutine dgbtf2(M, N, KL, KU, AB, LDAB, IPIV, INFO)
DGBTF2 computes the LU factorization of a general band matrix using the unblocked version of the algo...
Definition: dgbtf2.f:145
subroutine dgetrf(M, N, A, LDA, IPIV, INFO)
DGETRF
Definition: dgetrf.f:108
subroutine dgecon(NORM, N, A, LDA, ANORM, RCOND, WORK, IWORK, INFO)
DGECON
Definition: dgecon.f:124
subroutine dgeequb(M, N, A, LDA, R, C, ROWCND, COLCND, AMAX, INFO)
DGEEQUB
Definition: dgeequb.f:146
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 dgetrs(TRANS, N, NRHS, A, LDA, IPIV, B, LDB, INFO)
DGETRS
Definition: dgetrs.f:121
subroutine dgetf2(M, N, A, LDA, IPIV, INFO)
DGETF2 computes the LU factorization of a general m-by-n matrix using partial pivoting with row inter...
Definition: dgetf2.f:108
subroutine dgerfs(TRANS, N, NRHS, A, LDA, AF, LDAF, IPIV, B, LDB, X, LDX, FERR, BERR, WORK, IWORK, INFO)
DGERFS
Definition: dgerfs.f:185
subroutine dgerfsx(TRANS, EQUED, N, NRHS, A, LDA, AF, LDAF, IPIV, R, C, B, LDB, X, LDX, RCOND, BERR, N_ERR_BNDS, ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, IWORK, INFO)
DGERFSX
Definition: dgerfsx.f:414