LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
derrgex.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine derrge (PATH, NUNIT)
 DERRGEX More...
 

Function/Subroutine Documentation

subroutine derrge ( character*3  PATH,
integer  NUNIT 
)

DERRGEX

Purpose:
 DERRGE tests the error exits for the DOUBLE PRECISION routines
 for general matrices.

 Note that this file is used only when the XBLAS are available,
 otherwise derrge.f defines this subroutine.
Parameters
[in]PATH
          PATH is CHARACTER*3
          The LAPACK path name for the routines to be tested.
[in]NUNIT
          NUNIT is INTEGER
          The unit number for output.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 60 of file derrgex.f.

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

Here is the call graph for this function: