1 SUBROUTINE DCHKGE( DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NNS,
2 $ NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B,
3 $ X, XACT, WORK, RWORK, IWORK, NOUT )
4 *
5 * -- LAPACK test routine (version 3.2.1) --
6 *
7 * -- April 2009 --
8 *
9 * -- LAPACK is a software package provided by Univ. of Tennessee, --
10 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
11 *
12 * .. Scalar Arguments ..
13 LOGICAL TSTERR
14 INTEGER NM, NMAX, NN, NNB, NNS, NOUT
15 DOUBLE PRECISION THRESH
16 * ..
17 * .. Array Arguments ..
18 LOGICAL DOTYPE( * )
19 INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
20 $ NVAL( * )
21 DOUBLE PRECISION A( * ), AFAC( * ), AINV( * ), B( * ),
22 $ RWORK( * ), WORK( * ), X( * ), XACT( * )
23 * ..
24 *
25 * Purpose
26 * =======
27 *
28 * DCHKGE tests DGETRF, -TRI, -TRS, -RFS, and -CON.
29 *
30 * Arguments
31 * =========
32 *
33 * DOTYPE (input) LOGICAL array, dimension (NTYPES)
34 * The matrix types to be used for testing. Matrices of type j
35 * (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
36 * .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
37 *
38 * NM (input) INTEGER
39 * The number of values of M contained in the vector MVAL.
40 *
41 * MVAL (input) INTEGER array, dimension (NM)
42 * The values of the matrix row dimension M.
43 *
44 * NN (input) INTEGER
45 * The number of values of N contained in the vector NVAL.
46 *
47 * NVAL (input) INTEGER array, dimension (NN)
48 * The values of the matrix column dimension N.
49 *
50 * NNB (input) INTEGER
51 * The number of values of NB contained in the vector NBVAL.
52 *
53 * NBVAL (input) INTEGER array, dimension (NBVAL)
54 * The values of the blocksize NB.
55 *
56 * NNS (input) INTEGER
57 * The number of values of NRHS contained in the vector NSVAL.
58 *
59 * NSVAL (input) INTEGER array, dimension (NNS)
60 * The values of the number of right hand sides NRHS.
61 *
62 * THRESH (input) DOUBLE PRECISION
63 * The threshold value for the test ratios. A result is
64 * included in the output file if RESULT >= THRESH. To have
65 * every test ratio printed, use THRESH = 0.
66 *
67 * TSTERR (input) LOGICAL
68 * Flag that indicates whether error exits are to be tested.
69 *
70 * NMAX (input) INTEGER
71 * The maximum value permitted for M or N, used in dimensioning
72 * the work arrays.
73 *
74 * A (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
75 *
76 * AFAC (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
77 *
78 * AINV (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
79 *
80 * B (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
81 * where NSMAX is the largest entry in NSVAL.
82 *
83 * X (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
84 *
85 * XACT (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
86 *
87 * WORK (workspace) DOUBLE PRECISION array, dimension
88 * (NMAX*max(3,NSMAX))
89 *
90 * RWORK (workspace) DOUBLE PRECISION array, dimension
91 * (max(2*NMAX,2*NSMAX+NWORK))
92 *
93 * IWORK (workspace) INTEGER array, dimension (2*NMAX)
94 *
95 * NOUT (input) INTEGER
96 * The unit number for output.
97 *
98 * =====================================================================
99 *
100 * .. Parameters ..
101 DOUBLE PRECISION ONE, ZERO
102 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
103 INTEGER NTYPES
104 PARAMETER ( NTYPES = 11 )
105 INTEGER NTESTS
106 PARAMETER ( NTESTS = 8 )
107 INTEGER NTRAN
108 PARAMETER ( NTRAN = 3 )
109 * ..
110 * .. Local Scalars ..
111 LOGICAL TRFCON, ZEROT
112 CHARACTER DIST, NORM, TRANS, TYPE, XTYPE
113 CHARACTER*3 PATH
114 INTEGER I, IM, IMAT, IN, INB, INFO, IOFF, IRHS, ITRAN,
115 $ IZERO, K, KL, KU, LDA, LWORK, M, MODE, N, NB,
116 $ NERRS, NFAIL, NIMAT, NRHS, NRUN, NT
117 DOUBLE PRECISION AINVNM, ANORM, ANORMI, ANORMO, CNDNUM, DUMMY,
118 $ RCOND, RCONDC, RCONDI, RCONDO
119 * ..
120 * .. Local Arrays ..
121 CHARACTER TRANSS( NTRAN )
122 INTEGER ISEED( 4 ), ISEEDY( 4 )
123 DOUBLE PRECISION RESULT( NTESTS )
124 * ..
125 * .. External Functions ..
126 DOUBLE PRECISION DGET06, DLANGE
127 EXTERNAL DGET06, DLANGE
128 * ..
129 * .. External Subroutines ..
130 EXTERNAL ALAERH, ALAHD, ALASUM, DERRGE, DGECON, DGERFS,
131 $ DGET01, DGET02, DGET03, DGET04, DGET07, DGETRF,
132 $ DGETRI, DGETRS, DLACPY, DLARHS, DLASET, DLATB4,
133 $ DLATMS, XLAENV
134 * ..
135 * .. Intrinsic Functions ..
136 INTRINSIC MAX, MIN
137 * ..
138 * .. Scalars in Common ..
139 LOGICAL LERR, OK
140 CHARACTER*32 SRNAMT
141 INTEGER INFOT, NUNIT
142 * ..
143 * .. Common blocks ..
144 COMMON / INFOC / INFOT, NUNIT, OK, LERR
145 COMMON / SRNAMC / SRNAMT
146 * ..
147 * .. Data statements ..
148 DATA ISEEDY / 1988, 1989, 1990, 1991 / ,
149 $ TRANSS / 'N', 'T', 'C' /
150 * ..
151 * .. Executable Statements ..
152 *
153 * Initialize constants and the random number seed.
154 *
155 PATH( 1: 1 ) = 'Double precision'
156 PATH( 2: 3 ) = 'GE'
157 NRUN = 0
158 NFAIL = 0
159 NERRS = 0
160 DO 10 I = 1, 4
161 ISEED( I ) = ISEEDY( I )
162 10 CONTINUE
163 *
164 * Test the error exits
165 *
166 CALL XLAENV( 1, 1 )
167 IF( TSTERR )
168 $ CALL DERRGE( PATH, NOUT )
169 INFOT = 0
170 CALL XLAENV( 2, 2 )
171 *
172 * Do for each value of M in MVAL
173 *
174 DO 120 IM = 1, NM
175 M = MVAL( IM )
176 LDA = MAX( 1, M )
177 *
178 * Do for each value of N in NVAL
179 *
180 DO 110 IN = 1, NN
181 N = NVAL( IN )
182 XTYPE = 'N'
183 NIMAT = NTYPES
184 IF( M.LE.0 .OR. N.LE.0 )
185 $ NIMAT = 1
186 *
187 DO 100 IMAT = 1, NIMAT
188 *
189 * Do the tests only if DOTYPE( IMAT ) is true.
190 *
191 IF( .NOT.DOTYPE( IMAT ) )
192 $ GO TO 100
193 *
194 * Skip types 5, 6, or 7 if the matrix size is too small.
195 *
196 ZEROT = IMAT.GE.5 .AND. IMAT.LE.7
197 IF( ZEROT .AND. N.LT.IMAT-4 )
198 $ GO TO 100
199 *
200 * Set up parameters with DLATB4 and generate a test matrix
201 * with DLATMS.
202 *
203 CALL DLATB4( PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE,
204 $ CNDNUM, DIST )
205 *
206 SRNAMT = 'DLATMS'
207 CALL DLATMS( M, N, DIST, ISEED, TYPE, RWORK, MODE,
208 $ CNDNUM, ANORM, KL, KU, 'No packing', A, LDA,
209 $ WORK, INFO )
210 *
211 * Check error code from DLATMS.
212 *
213 IF( INFO.NE.0 ) THEN
214 CALL ALAERH( PATH, 'DLATMS', INFO, 0, ' ', M, N, -1,
215 $ -1, -1, IMAT, NFAIL, NERRS, NOUT )
216 GO TO 100
217 END IF
218 *
219 * For types 5-7, zero one or more columns of the matrix to
220 * test that INFO is returned correctly.
221 *
222 IF( ZEROT ) THEN
223 IF( IMAT.EQ.5 ) THEN
224 IZERO = 1
225 ELSE IF( IMAT.EQ.6 ) THEN
226 IZERO = MIN( M, N )
227 ELSE
228 IZERO = MIN( M, N ) / 2 + 1
229 END IF
230 IOFF = ( IZERO-1 )*LDA
231 IF( IMAT.LT.7 ) THEN
232 DO 20 I = 1, M
233 A( IOFF+I ) = ZERO
234 20 CONTINUE
235 ELSE
236 CALL DLASET( 'Full', M, N-IZERO+1, ZERO, ZERO,
237 $ A( IOFF+1 ), LDA )
238 END IF
239 ELSE
240 IZERO = 0
241 END IF
242 *
243 * These lines, if used in place of the calls in the DO 60
244 * loop, cause the code to bomb on a Sun SPARCstation.
245 *
246 * ANORMO = DLANGE( 'O', M, N, A, LDA, RWORK )
247 * ANORMI = DLANGE( 'I', M, N, A, LDA, RWORK )
248 *
249 * Do for each blocksize in NBVAL
250 *
251 DO 90 INB = 1, NNB
252 NB = NBVAL( INB )
253 CALL XLAENV( 1, NB )
254 *
255 * Compute the LU factorization of the matrix.
256 *
257 CALL DLACPY( 'Full', M, N, A, LDA, AFAC, LDA )
258 SRNAMT = 'DGETRF'
259 CALL DGETRF( M, N, AFAC, LDA, IWORK, INFO )
260 *
261 * Check error code from DGETRF.
262 *
263 IF( INFO.NE.IZERO )
264 $ CALL ALAERH( PATH, 'DGETRF', INFO, IZERO, ' ', M,
265 $ N, -1, -1, NB, IMAT, NFAIL, NERRS,
266 $ NOUT )
267 TRFCON = .FALSE.
268 *
269 *+ TEST 1
270 * Reconstruct matrix from factors and compute residual.
271 *
272 CALL DLACPY( 'Full', M, N, AFAC, LDA, AINV, LDA )
273 CALL DGET01( M, N, A, LDA, AINV, LDA, IWORK, RWORK,
274 $ RESULT( 1 ) )
275 NT = 1
276 *
277 *+ TEST 2
278 * Form the inverse if the factorization was successful
279 * and compute the residual.
280 *
281 IF( M.EQ.N .AND. INFO.EQ.0 ) THEN
282 CALL DLACPY( 'Full', N, N, AFAC, LDA, AINV, LDA )
283 SRNAMT = 'DGETRI'
284 NRHS = NSVAL( 1 )
285 LWORK = NMAX*MAX( 3, NRHS )
286 CALL DGETRI( N, AINV, LDA, IWORK, WORK, LWORK,
287 $ INFO )
288 *
289 * Check error code from DGETRI.
290 *
291 IF( INFO.NE.0 )
292 $ CALL ALAERH( PATH, 'DGETRI', INFO, 0, ' ', N, N,
293 $ -1, -1, NB, IMAT, NFAIL, NERRS,
294 $ NOUT )
295 *
296 * Compute the residual for the matrix times its
297 * inverse. Also compute the 1-norm condition number
298 * of A.
299 *
300 CALL DGET03( N, A, LDA, AINV, LDA, WORK, LDA,
301 $ RWORK, RCONDO, RESULT( 2 ) )
302 ANORMO = DLANGE( 'O', M, N, A, LDA, RWORK )
303 *
304 * Compute the infinity-norm condition number of A.
305 *
306 ANORMI = DLANGE( 'I', M, N, A, LDA, RWORK )
307 AINVNM = DLANGE( 'I', N, N, AINV, LDA, RWORK )
308 IF( ANORMI.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN
309 RCONDI = ONE
310 ELSE
311 RCONDI = ( ONE / ANORMI ) / AINVNM
312 END IF
313 NT = 2
314 ELSE
315 *
316 * Do only the condition estimate if INFO > 0.
317 *
318 TRFCON = .TRUE.
319 ANORMO = DLANGE( 'O', M, N, A, LDA, RWORK )
320 ANORMI = DLANGE( 'I', M, N, A, LDA, RWORK )
321 RCONDO = ZERO
322 RCONDI = ZERO
323 END IF
324 *
325 * Print information about the tests so far that did not
326 * pass the threshold.
327 *
328 DO 30 K = 1, NT
329 IF( RESULT( K ).GE.THRESH ) THEN
330 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
331 $ CALL ALAHD( NOUT, PATH )
332 WRITE( NOUT, FMT = 9999 )M, N, NB, IMAT, K,
333 $ RESULT( K )
334 NFAIL = NFAIL + 1
335 END IF
336 30 CONTINUE
337 NRUN = NRUN + NT
338 *
339 * Skip the remaining tests if this is not the first
340 * block size or if M .ne. N. Skip the solve tests if
341 * the matrix is singular.
342 *
343 IF( INB.GT.1 .OR. M.NE.N )
344 $ GO TO 90
345 IF( TRFCON )
346 $ GO TO 70
347 *
348 DO 60 IRHS = 1, NNS
349 NRHS = NSVAL( IRHS )
350 XTYPE = 'N'
351 *
352 DO 50 ITRAN = 1, NTRAN
353 TRANS = TRANSS( ITRAN )
354 IF( ITRAN.EQ.1 ) THEN
355 RCONDC = RCONDO
356 ELSE
357 RCONDC = RCONDI
358 END IF
359 *
360 *+ TEST 3
361 * Solve and compute residual for A * X = B.
362 *
363 SRNAMT = 'DLARHS'
364 CALL DLARHS( PATH, XTYPE, ' ', TRANS, N, N, KL,
365 $ KU, NRHS, A, LDA, XACT, LDA, B,
366 $ LDA, ISEED, INFO )
367 XTYPE = 'C'
368 *
369 CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA )
370 SRNAMT = 'DGETRS'
371 CALL DGETRS( TRANS, N, NRHS, AFAC, LDA, IWORK,
372 $ X, LDA, INFO )
373 *
374 * Check error code from DGETRS.
375 *
376 IF( INFO.NE.0 )
377 $ CALL ALAERH( PATH, 'DGETRS', INFO, 0, TRANS,
378 $ N, N, -1, -1, NRHS, IMAT, NFAIL,
379 $ NERRS, NOUT )
380 *
381 CALL DLACPY( 'Full', N, NRHS, B, LDA, WORK,
382 $ LDA )
383 CALL DGET02( TRANS, N, N, NRHS, A, LDA, X, LDA,
384 $ WORK, LDA, RWORK, RESULT( 3 ) )
385 *
386 *+ TEST 4
387 * Check solution from generated exact solution.
388 *
389 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
390 $ RESULT( 4 ) )
391 *
392 *+ TESTS 5, 6, and 7
393 * Use iterative refinement to improve the
394 * solution.
395 *
396 SRNAMT = 'DGERFS'
397 CALL DGERFS( TRANS, N, NRHS, A, LDA, AFAC, LDA,
398 $ IWORK, B, LDA, X, LDA, RWORK,
399 $ RWORK( NRHS+1 ), WORK,
400 $ IWORK( N+1 ), INFO )
401 *
402 * Check error code from DGERFS.
403 *
404 IF( INFO.NE.0 )
405 $ CALL ALAERH( PATH, 'DGERFS', INFO, 0, TRANS,
406 $ N, N, -1, -1, NRHS, IMAT, NFAIL,
407 $ NERRS, NOUT )
408 *
409 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
410 $ RESULT( 5 ) )
411 CALL DGET07( TRANS, N, NRHS, A, LDA, B, LDA, X,
412 $ LDA, XACT, LDA, RWORK, .TRUE.,
413 $ RWORK( NRHS+1 ), RESULT( 6 ) )
414 *
415 * Print information about the tests that did not
416 * pass the threshold.
417 *
418 DO 40 K = 3, 7
419 IF( RESULT( K ).GE.THRESH ) THEN
420 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
421 $ CALL ALAHD( NOUT, PATH )
422 WRITE( NOUT, FMT = 9998 )TRANS, N, NRHS,
423 $ IMAT, K, RESULT( K )
424 NFAIL = NFAIL + 1
425 END IF
426 40 CONTINUE
427 NRUN = NRUN + 5
428 50 CONTINUE
429 60 CONTINUE
430 *
431 *+ TEST 8
432 * Get an estimate of RCOND = 1/CNDNUM.
433 *
434 70 CONTINUE
435 DO 80 ITRAN = 1, 2
436 IF( ITRAN.EQ.1 ) THEN
437 ANORM = ANORMO
438 RCONDC = RCONDO
439 NORM = 'O'
440 ELSE
441 ANORM = ANORMI
442 RCONDC = RCONDI
443 NORM = 'I'
444 END IF
445 SRNAMT = 'DGECON'
446 CALL DGECON( NORM, N, AFAC, LDA, ANORM, RCOND,
447 $ WORK, IWORK( N+1 ), INFO )
448 *
449 * Check error code from DGECON.
450 *
451 IF( INFO.NE.0 )
452 $ CALL ALAERH( PATH, 'DGECON', INFO, 0, NORM, N,
453 $ N, -1, -1, -1, IMAT, NFAIL, NERRS,
454 $ NOUT )
455 *
456 * This line is needed on a Sun SPARCstation.
457 *
458 DUMMY = RCOND
459 *
460 RESULT( 8 ) = DGET06( RCOND, RCONDC )
461 *
462 * Print information about the tests that did not pass
463 * the threshold.
464 *
465 IF( RESULT( 8 ).GE.THRESH ) THEN
466 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
467 $ CALL ALAHD( NOUT, PATH )
468 WRITE( NOUT, FMT = 9997 )NORM, N, IMAT, 8,
469 $ RESULT( 8 )
470 NFAIL = NFAIL + 1
471 END IF
472 NRUN = NRUN + 1
473 80 CONTINUE
474 90 CONTINUE
475 100 CONTINUE
476 110 CONTINUE
477 120 CONTINUE
478 *
479 * Print a summary of the results.
480 *
481 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS )
482 *
483 9999 FORMAT( ' M = ', I5, ', N =', I5, ', NB =', I4, ', type ', I2,
484 $ ', test(', I2, ') =', G12.5 )
485 9998 FORMAT( ' TRANS=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ',
486 $ I2, ', test(', I2, ') =', G12.5 )
487 9997 FORMAT( ' NORM =''', A1, ''', N =', I5, ',', 10X, ' type ', I2,
488 $ ', test(', I2, ') =', G12.5 )
489 RETURN
490 *
491 * End of DCHKGE
492 *
493 END
2 $ NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B,
3 $ X, XACT, WORK, RWORK, IWORK, NOUT )
4 *
5 * -- LAPACK test routine (version 3.2.1) --
6 *
7 * -- April 2009 --
8 *
9 * -- LAPACK is a software package provided by Univ. of Tennessee, --
10 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
11 *
12 * .. Scalar Arguments ..
13 LOGICAL TSTERR
14 INTEGER NM, NMAX, NN, NNB, NNS, NOUT
15 DOUBLE PRECISION THRESH
16 * ..
17 * .. Array Arguments ..
18 LOGICAL DOTYPE( * )
19 INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
20 $ NVAL( * )
21 DOUBLE PRECISION A( * ), AFAC( * ), AINV( * ), B( * ),
22 $ RWORK( * ), WORK( * ), X( * ), XACT( * )
23 * ..
24 *
25 * Purpose
26 * =======
27 *
28 * DCHKGE tests DGETRF, -TRI, -TRS, -RFS, and -CON.
29 *
30 * Arguments
31 * =========
32 *
33 * DOTYPE (input) LOGICAL array, dimension (NTYPES)
34 * The matrix types to be used for testing. Matrices of type j
35 * (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
36 * .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
37 *
38 * NM (input) INTEGER
39 * The number of values of M contained in the vector MVAL.
40 *
41 * MVAL (input) INTEGER array, dimension (NM)
42 * The values of the matrix row dimension M.
43 *
44 * NN (input) INTEGER
45 * The number of values of N contained in the vector NVAL.
46 *
47 * NVAL (input) INTEGER array, dimension (NN)
48 * The values of the matrix column dimension N.
49 *
50 * NNB (input) INTEGER
51 * The number of values of NB contained in the vector NBVAL.
52 *
53 * NBVAL (input) INTEGER array, dimension (NBVAL)
54 * The values of the blocksize NB.
55 *
56 * NNS (input) INTEGER
57 * The number of values of NRHS contained in the vector NSVAL.
58 *
59 * NSVAL (input) INTEGER array, dimension (NNS)
60 * The values of the number of right hand sides NRHS.
61 *
62 * THRESH (input) DOUBLE PRECISION
63 * The threshold value for the test ratios. A result is
64 * included in the output file if RESULT >= THRESH. To have
65 * every test ratio printed, use THRESH = 0.
66 *
67 * TSTERR (input) LOGICAL
68 * Flag that indicates whether error exits are to be tested.
69 *
70 * NMAX (input) INTEGER
71 * The maximum value permitted for M or N, used in dimensioning
72 * the work arrays.
73 *
74 * A (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
75 *
76 * AFAC (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
77 *
78 * AINV (workspace) DOUBLE PRECISION array, dimension (NMAX*NMAX)
79 *
80 * B (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
81 * where NSMAX is the largest entry in NSVAL.
82 *
83 * X (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
84 *
85 * XACT (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX)
86 *
87 * WORK (workspace) DOUBLE PRECISION array, dimension
88 * (NMAX*max(3,NSMAX))
89 *
90 * RWORK (workspace) DOUBLE PRECISION array, dimension
91 * (max(2*NMAX,2*NSMAX+NWORK))
92 *
93 * IWORK (workspace) INTEGER array, dimension (2*NMAX)
94 *
95 * NOUT (input) INTEGER
96 * The unit number for output.
97 *
98 * =====================================================================
99 *
100 * .. Parameters ..
101 DOUBLE PRECISION ONE, ZERO
102 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
103 INTEGER NTYPES
104 PARAMETER ( NTYPES = 11 )
105 INTEGER NTESTS
106 PARAMETER ( NTESTS = 8 )
107 INTEGER NTRAN
108 PARAMETER ( NTRAN = 3 )
109 * ..
110 * .. Local Scalars ..
111 LOGICAL TRFCON, ZEROT
112 CHARACTER DIST, NORM, TRANS, TYPE, XTYPE
113 CHARACTER*3 PATH
114 INTEGER I, IM, IMAT, IN, INB, INFO, IOFF, IRHS, ITRAN,
115 $ IZERO, K, KL, KU, LDA, LWORK, M, MODE, N, NB,
116 $ NERRS, NFAIL, NIMAT, NRHS, NRUN, NT
117 DOUBLE PRECISION AINVNM, ANORM, ANORMI, ANORMO, CNDNUM, DUMMY,
118 $ RCOND, RCONDC, RCONDI, RCONDO
119 * ..
120 * .. Local Arrays ..
121 CHARACTER TRANSS( NTRAN )
122 INTEGER ISEED( 4 ), ISEEDY( 4 )
123 DOUBLE PRECISION RESULT( NTESTS )
124 * ..
125 * .. External Functions ..
126 DOUBLE PRECISION DGET06, DLANGE
127 EXTERNAL DGET06, DLANGE
128 * ..
129 * .. External Subroutines ..
130 EXTERNAL ALAERH, ALAHD, ALASUM, DERRGE, DGECON, DGERFS,
131 $ DGET01, DGET02, DGET03, DGET04, DGET07, DGETRF,
132 $ DGETRI, DGETRS, DLACPY, DLARHS, DLASET, DLATB4,
133 $ DLATMS, XLAENV
134 * ..
135 * .. Intrinsic Functions ..
136 INTRINSIC MAX, MIN
137 * ..
138 * .. Scalars in Common ..
139 LOGICAL LERR, OK
140 CHARACTER*32 SRNAMT
141 INTEGER INFOT, NUNIT
142 * ..
143 * .. Common blocks ..
144 COMMON / INFOC / INFOT, NUNIT, OK, LERR
145 COMMON / SRNAMC / SRNAMT
146 * ..
147 * .. Data statements ..
148 DATA ISEEDY / 1988, 1989, 1990, 1991 / ,
149 $ TRANSS / 'N', 'T', 'C' /
150 * ..
151 * .. Executable Statements ..
152 *
153 * Initialize constants and the random number seed.
154 *
155 PATH( 1: 1 ) = 'Double precision'
156 PATH( 2: 3 ) = 'GE'
157 NRUN = 0
158 NFAIL = 0
159 NERRS = 0
160 DO 10 I = 1, 4
161 ISEED( I ) = ISEEDY( I )
162 10 CONTINUE
163 *
164 * Test the error exits
165 *
166 CALL XLAENV( 1, 1 )
167 IF( TSTERR )
168 $ CALL DERRGE( PATH, NOUT )
169 INFOT = 0
170 CALL XLAENV( 2, 2 )
171 *
172 * Do for each value of M in MVAL
173 *
174 DO 120 IM = 1, NM
175 M = MVAL( IM )
176 LDA = MAX( 1, M )
177 *
178 * Do for each value of N in NVAL
179 *
180 DO 110 IN = 1, NN
181 N = NVAL( IN )
182 XTYPE = 'N'
183 NIMAT = NTYPES
184 IF( M.LE.0 .OR. N.LE.0 )
185 $ NIMAT = 1
186 *
187 DO 100 IMAT = 1, NIMAT
188 *
189 * Do the tests only if DOTYPE( IMAT ) is true.
190 *
191 IF( .NOT.DOTYPE( IMAT ) )
192 $ GO TO 100
193 *
194 * Skip types 5, 6, or 7 if the matrix size is too small.
195 *
196 ZEROT = IMAT.GE.5 .AND. IMAT.LE.7
197 IF( ZEROT .AND. N.LT.IMAT-4 )
198 $ GO TO 100
199 *
200 * Set up parameters with DLATB4 and generate a test matrix
201 * with DLATMS.
202 *
203 CALL DLATB4( PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE,
204 $ CNDNUM, DIST )
205 *
206 SRNAMT = 'DLATMS'
207 CALL DLATMS( M, N, DIST, ISEED, TYPE, RWORK, MODE,
208 $ CNDNUM, ANORM, KL, KU, 'No packing', A, LDA,
209 $ WORK, INFO )
210 *
211 * Check error code from DLATMS.
212 *
213 IF( INFO.NE.0 ) THEN
214 CALL ALAERH( PATH, 'DLATMS', INFO, 0, ' ', M, N, -1,
215 $ -1, -1, IMAT, NFAIL, NERRS, NOUT )
216 GO TO 100
217 END IF
218 *
219 * For types 5-7, zero one or more columns of the matrix to
220 * test that INFO is returned correctly.
221 *
222 IF( ZEROT ) THEN
223 IF( IMAT.EQ.5 ) THEN
224 IZERO = 1
225 ELSE IF( IMAT.EQ.6 ) THEN
226 IZERO = MIN( M, N )
227 ELSE
228 IZERO = MIN( M, N ) / 2 + 1
229 END IF
230 IOFF = ( IZERO-1 )*LDA
231 IF( IMAT.LT.7 ) THEN
232 DO 20 I = 1, M
233 A( IOFF+I ) = ZERO
234 20 CONTINUE
235 ELSE
236 CALL DLASET( 'Full', M, N-IZERO+1, ZERO, ZERO,
237 $ A( IOFF+1 ), LDA )
238 END IF
239 ELSE
240 IZERO = 0
241 END IF
242 *
243 * These lines, if used in place of the calls in the DO 60
244 * loop, cause the code to bomb on a Sun SPARCstation.
245 *
246 * ANORMO = DLANGE( 'O', M, N, A, LDA, RWORK )
247 * ANORMI = DLANGE( 'I', M, N, A, LDA, RWORK )
248 *
249 * Do for each blocksize in NBVAL
250 *
251 DO 90 INB = 1, NNB
252 NB = NBVAL( INB )
253 CALL XLAENV( 1, NB )
254 *
255 * Compute the LU factorization of the matrix.
256 *
257 CALL DLACPY( 'Full', M, N, A, LDA, AFAC, LDA )
258 SRNAMT = 'DGETRF'
259 CALL DGETRF( M, N, AFAC, LDA, IWORK, INFO )
260 *
261 * Check error code from DGETRF.
262 *
263 IF( INFO.NE.IZERO )
264 $ CALL ALAERH( PATH, 'DGETRF', INFO, IZERO, ' ', M,
265 $ N, -1, -1, NB, IMAT, NFAIL, NERRS,
266 $ NOUT )
267 TRFCON = .FALSE.
268 *
269 *+ TEST 1
270 * Reconstruct matrix from factors and compute residual.
271 *
272 CALL DLACPY( 'Full', M, N, AFAC, LDA, AINV, LDA )
273 CALL DGET01( M, N, A, LDA, AINV, LDA, IWORK, RWORK,
274 $ RESULT( 1 ) )
275 NT = 1
276 *
277 *+ TEST 2
278 * Form the inverse if the factorization was successful
279 * and compute the residual.
280 *
281 IF( M.EQ.N .AND. INFO.EQ.0 ) THEN
282 CALL DLACPY( 'Full', N, N, AFAC, LDA, AINV, LDA )
283 SRNAMT = 'DGETRI'
284 NRHS = NSVAL( 1 )
285 LWORK = NMAX*MAX( 3, NRHS )
286 CALL DGETRI( N, AINV, LDA, IWORK, WORK, LWORK,
287 $ INFO )
288 *
289 * Check error code from DGETRI.
290 *
291 IF( INFO.NE.0 )
292 $ CALL ALAERH( PATH, 'DGETRI', INFO, 0, ' ', N, N,
293 $ -1, -1, NB, IMAT, NFAIL, NERRS,
294 $ NOUT )
295 *
296 * Compute the residual for the matrix times its
297 * inverse. Also compute the 1-norm condition number
298 * of A.
299 *
300 CALL DGET03( N, A, LDA, AINV, LDA, WORK, LDA,
301 $ RWORK, RCONDO, RESULT( 2 ) )
302 ANORMO = DLANGE( 'O', M, N, A, LDA, RWORK )
303 *
304 * Compute the infinity-norm condition number of A.
305 *
306 ANORMI = DLANGE( 'I', M, N, A, LDA, RWORK )
307 AINVNM = DLANGE( 'I', N, N, AINV, LDA, RWORK )
308 IF( ANORMI.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN
309 RCONDI = ONE
310 ELSE
311 RCONDI = ( ONE / ANORMI ) / AINVNM
312 END IF
313 NT = 2
314 ELSE
315 *
316 * Do only the condition estimate if INFO > 0.
317 *
318 TRFCON = .TRUE.
319 ANORMO = DLANGE( 'O', M, N, A, LDA, RWORK )
320 ANORMI = DLANGE( 'I', M, N, A, LDA, RWORK )
321 RCONDO = ZERO
322 RCONDI = ZERO
323 END IF
324 *
325 * Print information about the tests so far that did not
326 * pass the threshold.
327 *
328 DO 30 K = 1, NT
329 IF( RESULT( K ).GE.THRESH ) THEN
330 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
331 $ CALL ALAHD( NOUT, PATH )
332 WRITE( NOUT, FMT = 9999 )M, N, NB, IMAT, K,
333 $ RESULT( K )
334 NFAIL = NFAIL + 1
335 END IF
336 30 CONTINUE
337 NRUN = NRUN + NT
338 *
339 * Skip the remaining tests if this is not the first
340 * block size or if M .ne. N. Skip the solve tests if
341 * the matrix is singular.
342 *
343 IF( INB.GT.1 .OR. M.NE.N )
344 $ GO TO 90
345 IF( TRFCON )
346 $ GO TO 70
347 *
348 DO 60 IRHS = 1, NNS
349 NRHS = NSVAL( IRHS )
350 XTYPE = 'N'
351 *
352 DO 50 ITRAN = 1, NTRAN
353 TRANS = TRANSS( ITRAN )
354 IF( ITRAN.EQ.1 ) THEN
355 RCONDC = RCONDO
356 ELSE
357 RCONDC = RCONDI
358 END IF
359 *
360 *+ TEST 3
361 * Solve and compute residual for A * X = B.
362 *
363 SRNAMT = 'DLARHS'
364 CALL DLARHS( PATH, XTYPE, ' ', TRANS, N, N, KL,
365 $ KU, NRHS, A, LDA, XACT, LDA, B,
366 $ LDA, ISEED, INFO )
367 XTYPE = 'C'
368 *
369 CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA )
370 SRNAMT = 'DGETRS'
371 CALL DGETRS( TRANS, N, NRHS, AFAC, LDA, IWORK,
372 $ X, LDA, INFO )
373 *
374 * Check error code from DGETRS.
375 *
376 IF( INFO.NE.0 )
377 $ CALL ALAERH( PATH, 'DGETRS', INFO, 0, TRANS,
378 $ N, N, -1, -1, NRHS, IMAT, NFAIL,
379 $ NERRS, NOUT )
380 *
381 CALL DLACPY( 'Full', N, NRHS, B, LDA, WORK,
382 $ LDA )
383 CALL DGET02( TRANS, N, N, NRHS, A, LDA, X, LDA,
384 $ WORK, LDA, RWORK, RESULT( 3 ) )
385 *
386 *+ TEST 4
387 * Check solution from generated exact solution.
388 *
389 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
390 $ RESULT( 4 ) )
391 *
392 *+ TESTS 5, 6, and 7
393 * Use iterative refinement to improve the
394 * solution.
395 *
396 SRNAMT = 'DGERFS'
397 CALL DGERFS( TRANS, N, NRHS, A, LDA, AFAC, LDA,
398 $ IWORK, B, LDA, X, LDA, RWORK,
399 $ RWORK( NRHS+1 ), WORK,
400 $ IWORK( N+1 ), INFO )
401 *
402 * Check error code from DGERFS.
403 *
404 IF( INFO.NE.0 )
405 $ CALL ALAERH( PATH, 'DGERFS', INFO, 0, TRANS,
406 $ N, N, -1, -1, NRHS, IMAT, NFAIL,
407 $ NERRS, NOUT )
408 *
409 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
410 $ RESULT( 5 ) )
411 CALL DGET07( TRANS, N, NRHS, A, LDA, B, LDA, X,
412 $ LDA, XACT, LDA, RWORK, .TRUE.,
413 $ RWORK( NRHS+1 ), RESULT( 6 ) )
414 *
415 * Print information about the tests that did not
416 * pass the threshold.
417 *
418 DO 40 K = 3, 7
419 IF( RESULT( K ).GE.THRESH ) THEN
420 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
421 $ CALL ALAHD( NOUT, PATH )
422 WRITE( NOUT, FMT = 9998 )TRANS, N, NRHS,
423 $ IMAT, K, RESULT( K )
424 NFAIL = NFAIL + 1
425 END IF
426 40 CONTINUE
427 NRUN = NRUN + 5
428 50 CONTINUE
429 60 CONTINUE
430 *
431 *+ TEST 8
432 * Get an estimate of RCOND = 1/CNDNUM.
433 *
434 70 CONTINUE
435 DO 80 ITRAN = 1, 2
436 IF( ITRAN.EQ.1 ) THEN
437 ANORM = ANORMO
438 RCONDC = RCONDO
439 NORM = 'O'
440 ELSE
441 ANORM = ANORMI
442 RCONDC = RCONDI
443 NORM = 'I'
444 END IF
445 SRNAMT = 'DGECON'
446 CALL DGECON( NORM, N, AFAC, LDA, ANORM, RCOND,
447 $ WORK, IWORK( N+1 ), INFO )
448 *
449 * Check error code from DGECON.
450 *
451 IF( INFO.NE.0 )
452 $ CALL ALAERH( PATH, 'DGECON', INFO, 0, NORM, N,
453 $ N, -1, -1, -1, IMAT, NFAIL, NERRS,
454 $ NOUT )
455 *
456 * This line is needed on a Sun SPARCstation.
457 *
458 DUMMY = RCOND
459 *
460 RESULT( 8 ) = DGET06( RCOND, RCONDC )
461 *
462 * Print information about the tests that did not pass
463 * the threshold.
464 *
465 IF( RESULT( 8 ).GE.THRESH ) THEN
466 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
467 $ CALL ALAHD( NOUT, PATH )
468 WRITE( NOUT, FMT = 9997 )NORM, N, IMAT, 8,
469 $ RESULT( 8 )
470 NFAIL = NFAIL + 1
471 END IF
472 NRUN = NRUN + 1
473 80 CONTINUE
474 90 CONTINUE
475 100 CONTINUE
476 110 CONTINUE
477 120 CONTINUE
478 *
479 * Print a summary of the results.
480 *
481 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS )
482 *
483 9999 FORMAT( ' M = ', I5, ', N =', I5, ', NB =', I4, ', type ', I2,
484 $ ', test(', I2, ') =', G12.5 )
485 9998 FORMAT( ' TRANS=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ',
486 $ I2, ', test(', I2, ') =', G12.5 )
487 9997 FORMAT( ' NORM =''', A1, ''', N =', I5, ',', 10X, ' type ', I2,
488 $ ', test(', I2, ') =', G12.5 )
489 RETURN
490 *
491 * End of DCHKGE
492 *
493 END