LAPACK  3.10.1
LAPACK: Linear Algebra PACKage

◆ derrge()

subroutine derrge ( character*3  PATH,
integer  NUNIT 
)

DERRGE

DERRGEX

Purpose:
 DERRGE tests the error exits for the DOUBLE PRECISION routines
 for general matrices.
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.
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.

Definition at line 54 of file derrge.f.

55 *
56 * -- LAPACK test routine --
57 * -- LAPACK is a software package provided by Univ. of Tennessee, --
58 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
59 *
60 * .. Scalar Arguments ..
61  CHARACTER*3 PATH
62  INTEGER NUNIT
63 * ..
64 *
65 * =====================================================================
66 *
67 * .. Parameters ..
68  INTEGER NMAX, LW
69  parameter( nmax = 4, lw = 3*nmax )
70 * ..
71 * .. Local Scalars ..
72  CHARACTER*2 C2
73  INTEGER I, INFO, J
74  DOUBLE PRECISION ANRM, CCOND, RCOND
75 * ..
76 * .. Local Arrays ..
77  INTEGER IP( NMAX ), IW( NMAX )
78  DOUBLE PRECISION A( NMAX, NMAX ), AF( NMAX, NMAX ), B( NMAX ),
79  $ R1( NMAX ), R2( NMAX ), W( LW ), X( NMAX )
80 * ..
81 * .. External Functions ..
82  LOGICAL LSAMEN
83  EXTERNAL lsamen
84 * ..
85 * .. External Subroutines ..
86  EXTERNAL alaesm, chkxer, dgbcon, dgbequ, dgbrfs, dgbtf2,
88  $ dgetrf, dgetri, dgetrs
89 * ..
90 * .. Scalars in Common ..
91  LOGICAL LERR, OK
92  CHARACTER*32 SRNAMT
93  INTEGER INFOT, NOUT
94 * ..
95 * .. Common blocks ..
96  COMMON / infoc / infot, nout, ok, lerr
97  COMMON / srnamc / srnamt
98 * ..
99 * .. Intrinsic Functions ..
100  INTRINSIC dble
101 * ..
102 * .. Executable Statements ..
103 *
104  nout = nunit
105  WRITE( nout, fmt = * )
106  c2 = path( 2: 3 )
107 *
108 * Set the variables to innocuous values.
109 *
110  DO 20 j = 1, nmax
111  DO 10 i = 1, nmax
112  a( i, j ) = 1.d0 / dble( i+j )
113  af( i, j ) = 1.d0 / dble( i+j )
114  10 CONTINUE
115  b( j ) = 0.d0
116  r1( j ) = 0.d0
117  r2( j ) = 0.d0
118  w( j ) = 0.d0
119  x( j ) = 0.d0
120  ip( j ) = j
121  iw( j ) = j
122  20 CONTINUE
123  ok = .true.
124 *
125  IF( lsamen( 2, c2, 'GE' ) ) THEN
126 *
127 * Test error exits of the routines that use the LU decomposition
128 * of a general matrix.
129 *
130 * DGETRF
131 *
132  srnamt = 'DGETRF'
133  infot = 1
134  CALL dgetrf( -1, 0, a, 1, ip, info )
135  CALL chkxer( 'DGETRF', infot, nout, lerr, ok )
136  infot = 2
137  CALL dgetrf( 0, -1, a, 1, ip, info )
138  CALL chkxer( 'DGETRF', infot, nout, lerr, ok )
139  infot = 4
140  CALL dgetrf( 2, 1, a, 1, ip, info )
141  CALL chkxer( 'DGETRF', infot, nout, lerr, ok )
142 *
143 * DGETF2
144 *
145  srnamt = 'DGETF2'
146  infot = 1
147  CALL dgetf2( -1, 0, a, 1, ip, info )
148  CALL chkxer( 'DGETF2', infot, nout, lerr, ok )
149  infot = 2
150  CALL dgetf2( 0, -1, a, 1, ip, info )
151  CALL chkxer( 'DGETF2', infot, nout, lerr, ok )
152  infot = 4
153  CALL dgetf2( 2, 1, a, 1, ip, info )
154  CALL chkxer( 'DGETF2', infot, nout, lerr, ok )
155 *
156 * DGETRI
157 *
158  srnamt = 'DGETRI'
159  infot = 1
160  CALL dgetri( -1, a, 1, ip, w, lw, info )
161  CALL chkxer( 'DGETRI', infot, nout, lerr, ok )
162  infot = 3
163  CALL dgetri( 2, a, 1, ip, w, lw, info )
164  CALL chkxer( 'DGETRI', infot, nout, lerr, ok )
165 *
166 * DGETRS
167 *
168  srnamt = 'DGETRS'
169  infot = 1
170  CALL dgetrs( '/', 0, 0, a, 1, ip, b, 1, info )
171  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
172  infot = 2
173  CALL dgetrs( 'N', -1, 0, a, 1, ip, b, 1, info )
174  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
175  infot = 3
176  CALL dgetrs( 'N', 0, -1, a, 1, ip, b, 1, info )
177  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
178  infot = 5
179  CALL dgetrs( 'N', 2, 1, a, 1, ip, b, 2, info )
180  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
181  infot = 8
182  CALL dgetrs( 'N', 2, 1, a, 2, ip, b, 1, info )
183  CALL chkxer( 'DGETRS', infot, nout, lerr, ok )
184 *
185 * DGERFS
186 *
187  srnamt = 'DGERFS'
188  infot = 1
189  CALL dgerfs( '/', 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1, r2, w,
190  $ iw, info )
191  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
192  infot = 2
193  CALL dgerfs( 'N', -1, 0, a, 1, af, 1, ip, b, 1, x, 1, r1, r2,
194  $ w, iw, info )
195  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
196  infot = 3
197  CALL dgerfs( 'N', 0, -1, a, 1, af, 1, ip, b, 1, x, 1, r1, r2,
198  $ w, iw, info )
199  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
200  infot = 5
201  CALL dgerfs( 'N', 2, 1, a, 1, af, 2, ip, b, 2, x, 2, r1, r2, w,
202  $ iw, info )
203  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
204  infot = 7
205  CALL dgerfs( 'N', 2, 1, a, 2, af, 1, ip, b, 2, x, 2, r1, r2, w,
206  $ iw, info )
207  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
208  infot = 10
209  CALL dgerfs( 'N', 2, 1, a, 2, af, 2, ip, b, 1, x, 2, r1, r2, w,
210  $ iw, info )
211  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
212  infot = 12
213  CALL dgerfs( 'N', 2, 1, a, 2, af, 2, ip, b, 2, x, 1, r1, r2, w,
214  $ iw, info )
215  CALL chkxer( 'DGERFS', infot, nout, lerr, ok )
216 *
217 * DGECON
218 *
219  srnamt = 'DGECON'
220  infot = 1
221  CALL dgecon( '/', 0, a, 1, anrm, rcond, w, iw, info )
222  CALL chkxer( 'DGECON', infot, nout, lerr, ok )
223  infot = 2
224  CALL dgecon( '1', -1, a, 1, anrm, rcond, w, iw, info )
225  CALL chkxer( 'DGECON', infot, nout, lerr, ok )
226  infot = 4
227  CALL dgecon( '1', 2, a, 1, anrm, rcond, w, iw, info )
228  CALL chkxer( 'DGECON', infot, nout, lerr, ok )
229 *
230 * DGEEQU
231 *
232  srnamt = 'DGEEQU'
233  infot = 1
234  CALL dgeequ( -1, 0, a, 1, r1, r2, rcond, ccond, anrm, info )
235  CALL chkxer( 'DGEEQU', infot, nout, lerr, ok )
236  infot = 2
237  CALL dgeequ( 0, -1, a, 1, r1, r2, rcond, ccond, anrm, info )
238  CALL chkxer( 'DGEEQU', infot, nout, lerr, ok )
239  infot = 4
240  CALL dgeequ( 2, 2, a, 1, r1, r2, rcond, ccond, anrm, info )
241  CALL chkxer( 'DGEEQU', infot, nout, lerr, ok )
242 *
243  ELSE IF( lsamen( 2, c2, 'GB' ) ) THEN
244 *
245 * Test error exits of the routines that use the LU decomposition
246 * of a general band matrix.
247 *
248 * DGBTRF
249 *
250  srnamt = 'DGBTRF'
251  infot = 1
252  CALL dgbtrf( -1, 0, 0, 0, a, 1, ip, info )
253  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
254  infot = 2
255  CALL dgbtrf( 0, -1, 0, 0, a, 1, ip, info )
256  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
257  infot = 3
258  CALL dgbtrf( 1, 1, -1, 0, a, 1, ip, info )
259  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
260  infot = 4
261  CALL dgbtrf( 1, 1, 0, -1, a, 1, ip, info )
262  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
263  infot = 6
264  CALL dgbtrf( 2, 2, 1, 1, a, 3, ip, info )
265  CALL chkxer( 'DGBTRF', infot, nout, lerr, ok )
266 *
267 * DGBTF2
268 *
269  srnamt = 'DGBTF2'
270  infot = 1
271  CALL dgbtf2( -1, 0, 0, 0, a, 1, ip, info )
272  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
273  infot = 2
274  CALL dgbtf2( 0, -1, 0, 0, a, 1, ip, info )
275  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
276  infot = 3
277  CALL dgbtf2( 1, 1, -1, 0, a, 1, ip, info )
278  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
279  infot = 4
280  CALL dgbtf2( 1, 1, 0, -1, a, 1, ip, info )
281  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
282  infot = 6
283  CALL dgbtf2( 2, 2, 1, 1, a, 3, ip, info )
284  CALL chkxer( 'DGBTF2', infot, nout, lerr, ok )
285 *
286 * DGBTRS
287 *
288  srnamt = 'DGBTRS'
289  infot = 1
290  CALL dgbtrs( '/', 0, 0, 0, 1, a, 1, ip, b, 1, info )
291  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
292  infot = 2
293  CALL dgbtrs( 'N', -1, 0, 0, 1, a, 1, ip, b, 1, info )
294  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
295  infot = 3
296  CALL dgbtrs( 'N', 1, -1, 0, 1, a, 1, ip, b, 1, info )
297  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
298  infot = 4
299  CALL dgbtrs( 'N', 1, 0, -1, 1, a, 1, ip, b, 1, info )
300  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
301  infot = 5
302  CALL dgbtrs( 'N', 1, 0, 0, -1, a, 1, ip, b, 1, info )
303  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
304  infot = 7
305  CALL dgbtrs( 'N', 2, 1, 1, 1, a, 3, ip, b, 2, info )
306  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
307  infot = 10
308  CALL dgbtrs( 'N', 2, 0, 0, 1, a, 1, ip, b, 1, info )
309  CALL chkxer( 'DGBTRS', infot, nout, lerr, ok )
310 *
311 * DGBRFS
312 *
313  srnamt = 'DGBRFS'
314  infot = 1
315  CALL dgbrfs( '/', 0, 0, 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
316  $ r2, w, iw, info )
317  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
318  infot = 2
319  CALL dgbrfs( 'N', -1, 0, 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
320  $ r2, w, iw, info )
321  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
322  infot = 3
323  CALL dgbrfs( 'N', 1, -1, 0, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
324  $ r2, w, iw, info )
325  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
326  infot = 4
327  CALL dgbrfs( 'N', 1, 0, -1, 0, a, 1, af, 1, ip, b, 1, x, 1, r1,
328  $ r2, w, iw, info )
329  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
330  infot = 5
331  CALL dgbrfs( 'N', 1, 0, 0, -1, a, 1, af, 1, ip, b, 1, x, 1, r1,
332  $ r2, w, iw, info )
333  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
334  infot = 7
335  CALL dgbrfs( 'N', 2, 1, 1, 1, a, 2, af, 4, ip, b, 2, x, 2, r1,
336  $ r2, w, iw, info )
337  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
338  infot = 9
339  CALL dgbrfs( 'N', 2, 1, 1, 1, a, 3, af, 3, ip, b, 2, x, 2, r1,
340  $ r2, w, iw, info )
341  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
342  infot = 12
343  CALL dgbrfs( 'N', 2, 0, 0, 1, a, 1, af, 1, ip, b, 1, x, 2, r1,
344  $ r2, w, iw, info )
345  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
346  infot = 14
347  CALL dgbrfs( 'N', 2, 0, 0, 1, a, 1, af, 1, ip, b, 2, x, 1, r1,
348  $ r2, w, iw, info )
349  CALL chkxer( 'DGBRFS', infot, nout, lerr, ok )
350 *
351 * DGBCON
352 *
353  srnamt = 'DGBCON'
354  infot = 1
355  CALL dgbcon( '/', 0, 0, 0, a, 1, ip, anrm, rcond, w, iw, info )
356  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
357  infot = 2
358  CALL dgbcon( '1', -1, 0, 0, a, 1, ip, anrm, rcond, w, iw,
359  $ info )
360  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
361  infot = 3
362  CALL dgbcon( '1', 1, -1, 0, a, 1, ip, anrm, rcond, w, iw,
363  $ info )
364  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
365  infot = 4
366  CALL dgbcon( '1', 1, 0, -1, a, 1, ip, anrm, rcond, w, iw,
367  $ info )
368  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
369  infot = 6
370  CALL dgbcon( '1', 2, 1, 1, a, 3, ip, anrm, rcond, w, iw, info )
371  CALL chkxer( 'DGBCON', infot, nout, lerr, ok )
372 *
373 * DGBEQU
374 *
375  srnamt = 'DGBEQU'
376  infot = 1
377  CALL dgbequ( -1, 0, 0, 0, a, 1, r1, r2, rcond, ccond, anrm,
378  $ info )
379  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
380  infot = 2
381  CALL dgbequ( 0, -1, 0, 0, a, 1, r1, r2, rcond, ccond, anrm,
382  $ info )
383  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
384  infot = 3
385  CALL dgbequ( 1, 1, -1, 0, a, 1, r1, r2, rcond, ccond, anrm,
386  $ info )
387  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
388  infot = 4
389  CALL dgbequ( 1, 1, 0, -1, a, 1, r1, r2, rcond, ccond, anrm,
390  $ info )
391  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
392  infot = 6
393  CALL dgbequ( 2, 2, 1, 1, a, 2, r1, r2, rcond, ccond, anrm,
394  $ info )
395  CALL chkxer( 'DGBEQU', infot, nout, lerr, ok )
396  END IF
397 *
398 * Print a summary line.
399 *
400  CALL alaesm( path, ok, nout )
401 *
402  RETURN
403 *
404 * End of DERRGE
405 *
subroutine chkxer(SRNAMT, INFOT, NOUT, LERR, OK)
Definition: cblat2.f:3196
logical function lsamen(N, CA, CB)
LSAMEN
Definition: lsamen.f:74
subroutine alaesm(PATH, OK, NOUT)
ALAESM
Definition: alaesm.f:63
subroutine dgbtrs(TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB, INFO)
DGBTRS
Definition: dgbtrs.f:138
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 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 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
Here is the call graph for this function:
Here is the caller graph for this function: