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