1 SUBROUTINE CDRVRFP( NOUT, NN, NVAL, NNS, NSVAL, NNT, NTVAL,
2 + THRESH, A, ASAV, AFAC, AINV, B,
3 + BSAV, XACT, X, ARF, ARFINV,
4 + C_WORK_CLATMS, C_WORK_CPOT02,
5 + C_WORK_CPOT03, S_WORK_CLATMS, S_WORK_CLANHE,
6 + S_WORK_CPOT01, S_WORK_CPOT02, S_WORK_CPOT03 )
7 *
8 * -- LAPACK test routine (version 3.2.0) --
9 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
10 * November 2008
11 *
12 * .. Scalar Arguments ..
13 INTEGER NN, NNS, NNT, NOUT
14 REAL THRESH
15 * ..
16 * .. Array Arguments ..
17 INTEGER NVAL( NN ), NSVAL( NNS ), NTVAL( NNT )
18 COMPLEX A( * )
19 COMPLEX AINV( * )
20 COMPLEX ASAV( * )
21 COMPLEX B( * )
22 COMPLEX BSAV( * )
23 COMPLEX AFAC( * )
24 COMPLEX ARF( * )
25 COMPLEX ARFINV( * )
26 COMPLEX XACT( * )
27 COMPLEX X( * )
28 COMPLEX C_WORK_CLATMS( * )
29 COMPLEX C_WORK_CPOT02( * )
30 COMPLEX C_WORK_CPOT03( * )
31 REAL S_WORK_CLATMS( * )
32 REAL S_WORK_CLANHE( * )
33 REAL S_WORK_CPOT01( * )
34 REAL S_WORK_CPOT02( * )
35 REAL S_WORK_CPOT03( * )
36 * ..
37 *
38 * Purpose
39 * =======
40 *
41 * CDRVRFP tests the LAPACK RFP routines:
42 * CPFTRF, CPFTRS, and CPFTRI.
43 *
44 * This testing routine follow the same tests as CDRVPO (test for the full
45 * format Symmetric Positive Definite solver).
46 *
47 * The tests are performed in Full Format, convertion back and forth from
48 * full format to RFP format are performed using the routines CTRTTF and
49 * CTFTTR.
50 *
51 * First, a specific matrix A of size N is created. There is nine types of
52 * different matrixes possible.
53 * 1. Diagonal 6. Random, CNDNUM = sqrt(0.1/EPS)
54 * 2. Random, CNDNUM = 2 7. Random, CNDNUM = 0.1/EPS
55 * *3. First row and column zero 8. Scaled near underflow
56 * *4. Last row and column zero 9. Scaled near overflow
57 * *5. Middle row and column zero
58 * (* - tests error exits from CPFTRF, no test ratios are computed)
59 * A solution XACT of size N-by-NRHS is created and the associated right
60 * hand side B as well. Then CPFTRF is called to compute L (or U), the
61 * Cholesky factor of A. Then L (or U) is used to solve the linear system
62 * of equations AX = B. This gives X. Then L (or U) is used to compute the
63 * inverse of A, AINV. The following four tests are then performed:
64 * (1) norm( L*L' - A ) / ( N * norm(A) * EPS ) or
65 * norm( U'*U - A ) / ( N * norm(A) * EPS ),
66 * (2) norm(B - A*X) / ( norm(A) * norm(X) * EPS ),
67 * (3) norm( I - A*AINV ) / ( N * norm(A) * norm(AINV) * EPS ),
68 * (4) ( norm(X-XACT) * RCOND ) / ( norm(XACT) * EPS ),
69 * where EPS is the machine precision, RCOND the condition number of A, and
70 * norm( . ) the 1-norm for (1,2,3) and the inf-norm for (4).
71 * Errors occur when INFO parameter is not as expected. Failures occur when
72 * a test ratios is greater than THRES.
73 *
74 * Arguments
75 * =========
76 *
77 * NOUT (input) INTEGER
78 * The unit number for output.
79 *
80 * NN (input) INTEGER
81 * The number of values of N contained in the vector NVAL.
82 *
83 * NVAL (input) INTEGER array, dimension (NN)
84 * The values of the matrix dimension N.
85 *
86 * NNS (input) INTEGER
87 * The number of values of NRHS contained in the vector NSVAL.
88 *
89 * NSVAL (input) INTEGER array, dimension (NNS)
90 * The values of the number of right-hand sides NRHS.
91 *
92 * NNT (input) INTEGER
93 * The number of values of MATRIX TYPE contained in the vector NTVAL.
94 *
95 * NTVAL (input) INTEGER array, dimension (NNT)
96 * The values of matrix type (between 0 and 9 for PO/PP/PF matrices).
97 *
98 * THRESH (input) REAL
99 * The threshold value for the test ratios. A result is
100 * included in the output file if RESULT >= THRESH. To have
101 * every test ratio printed, use THRESH = 0.
102 *
103 * A (workspace) COMPLEX array, dimension (NMAX*NMAX)
104 *
105 * ASAV (workspace) COMPLEX array, dimension (NMAX*NMAX)
106 *
107 * AFAC (workspace) COMPLEX array, dimension (NMAX*NMAX)
108 *
109 * AINV (workspace) COMPLEX array, dimension (NMAX*NMAX)
110 *
111 * B (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
112 *
113 * BSAV (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
114 *
115 * XACT (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
116 *
117 * X (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
118 *
119 * ARF (workspace) COMPLEX array, dimension ((NMAX*(NMAX+1))/2)
120 *
121 * ARFINV (workspace) COMPLEX array, dimension ((NMAX*(NMAX+1))/2)
122 *
123 * C_WORK_CLATMS (workspace) COMPLEX array, dimension ( 3*NMAX )
124 *
125 * C_WORK_CPOT02 (workspace) COMPLEX array, dimension ( NMAX*MAXRHS )
126 *
127 * C_WORK_CPOT03 (workspace) COMPLEX array, dimension ( NMAX*NMAX )
128 *
129 * S_WORK_CLATMS (workspace) REAL array, dimension ( NMAX )
130 *
131 * S_WORK_CLANHE (workspace) REAL array, dimension ( NMAX )
132 *
133 * S_WORK_CPOT01 (workspace) REAL array, dimension ( NMAX )
134 *
135 * S_WORK_CPOT02 (workspace) REAL array, dimension ( NMAX )
136 *
137 * S_WORK_CPOT03 (workspace) REAL array, dimension ( NMAX )
138 *
139 * =====================================================================
140 *
141 * .. Parameters ..
142 REAL ONE, ZERO
143 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 )
144 INTEGER NTESTS
145 PARAMETER ( NTESTS = 4 )
146 * ..
147 * .. Local Scalars ..
148 LOGICAL ZEROT
149 INTEGER I, INFO, IUPLO, LDA, LDB, IMAT, NERRS, NFAIL,
150 + NRHS, NRUN, IZERO, IOFF, K, NT, N, IFORM, IIN,
151 + IIT, IIS
152 CHARACTER DIST, CTYPE, UPLO, CFORM
153 INTEGER KL, KU, MODE
154 REAL ANORM, AINVNM, CNDNUM, RCONDC
155 * ..
156 * .. Local Arrays ..
157 CHARACTER UPLOS( 2 ), FORMS( 2 )
158 INTEGER ISEED( 4 ), ISEEDY( 4 )
159 REAL RESULT( NTESTS )
160 * ..
161 * .. External Functions ..
162 REAL CLANHE
163 EXTERNAL CLANHE
164 * ..
165 * .. External Subroutines ..
166 EXTERNAL ALADHD, ALAERH, ALASVM, CGET04, CTFTTR, CLACPY,
167 + CLAIPD, CLARHS, CLATB4, CLATMS, CPFTRI, CPFTRF,
168 + CPFTRS, CPOT01, CPOT02, CPOT03, CPOTRI, CPOTRF,
169 + CTRTTF
170 * ..
171 * .. Scalars in Common ..
172 CHARACTER*32 SRNAMT
173 * ..
174 * .. Common blocks ..
175 COMMON / SRNAMC / SRNAMT
176 * ..
177 * .. Data statements ..
178 DATA ISEEDY / 1988, 1989, 1990, 1991 /
179 DATA UPLOS / 'U', 'L' /
180 DATA FORMS / 'N', 'C' /
181 * ..
182 * .. Executable Statements ..
183 *
184 * Initialize constants and the random number seed.
185 *
186 NRUN = 0
187 NFAIL = 0
188 NERRS = 0
189 DO 10 I = 1, 4
190 ISEED( I ) = ISEEDY( I )
191 10 CONTINUE
192 *
193 DO 130 IIN = 1, NN
194 *
195 N = NVAL( IIN )
196 LDA = MAX( N, 1 )
197 LDB = MAX( N, 1 )
198 *
199 DO 980 IIS = 1, NNS
200 *
201 NRHS = NSVAL( IIS )
202 *
203 DO 120 IIT = 1, NNT
204 *
205 IMAT = NTVAL( IIT )
206 *
207 * If N.EQ.0, only consider the first type
208 *
209 IF( N.EQ.0 .AND. IIT.GT.1 ) GO TO 120
210 *
211 * Skip types 3, 4, or 5 if the matrix size is too small.
212 *
213 IF( IMAT.EQ.4 .AND. N.LE.1 ) GO TO 120
214 IF( IMAT.EQ.5 .AND. N.LE.2 ) GO TO 120
215 *
216 * Do first for UPLO = 'U', then for UPLO = 'L'
217 *
218 DO 110 IUPLO = 1, 2
219 UPLO = UPLOS( IUPLO )
220 *
221 * Do first for CFORM = 'N', then for CFORM = 'C'
222 *
223 DO 100 IFORM = 1, 2
224 CFORM = FORMS( IFORM )
225 *
226 * Set up parameters with CLATB4 and generate a test
227 * matrix with CLATMS.
228 *
229 CALL CLATB4( 'CPO', IMAT, N, N, CTYPE, KL, KU,
230 + ANORM, MODE, CNDNUM, DIST )
231 *
232 SRNAMT = 'CLATMS'
233 CALL CLATMS( N, N, DIST, ISEED, CTYPE,
234 + S_WORK_CLATMS,
235 + MODE, CNDNUM, ANORM, KL, KU, UPLO, A,
236 + LDA, C_WORK_CLATMS, INFO )
237 *
238 * Check error code from CLATMS.
239 *
240 IF( INFO.NE.0 ) THEN
241 CALL ALAERH( 'CPF', 'CLATMS', INFO, 0, UPLO, N,
242 + N, -1, -1, -1, IIT, NFAIL, NERRS,
243 + NOUT )
244 GO TO 100
245 END IF
246 *
247 * For types 3-5, zero one row and column of the matrix to
248 * test that INFO is returned correctly.
249 *
250 ZEROT = IMAT.GE.3 .AND. IMAT.LE.5
251 IF( ZEROT ) THEN
252 IF( IIT.EQ.3 ) THEN
253 IZERO = 1
254 ELSE IF( IIT.EQ.4 ) THEN
255 IZERO = N
256 ELSE
257 IZERO = N / 2 + 1
258 END IF
259 IOFF = ( IZERO-1 )*LDA
260 *
261 * Set row and column IZERO of A to 0.
262 *
263 IF( IUPLO.EQ.1 ) THEN
264 DO 20 I = 1, IZERO - 1
265 A( IOFF+I ) = ZERO
266 20 CONTINUE
267 IOFF = IOFF + IZERO
268 DO 30 I = IZERO, N
269 A( IOFF ) = ZERO
270 IOFF = IOFF + LDA
271 30 CONTINUE
272 ELSE
273 IOFF = IZERO
274 DO 40 I = 1, IZERO - 1
275 A( IOFF ) = ZERO
276 IOFF = IOFF + LDA
277 40 CONTINUE
278 IOFF = IOFF - IZERO
279 DO 50 I = IZERO, N
280 A( IOFF+I ) = ZERO
281 50 CONTINUE
282 END IF
283 ELSE
284 IZERO = 0
285 END IF
286 *
287 * Set the imaginary part of the diagonals.
288 *
289 CALL CLAIPD( N, A, LDA+1, 0 )
290 *
291 * Save a copy of the matrix A in ASAV.
292 *
293 CALL CLACPY( UPLO, N, N, A, LDA, ASAV, LDA )
294 *
295 * Compute the condition number of A (RCONDC).
296 *
297 IF( ZEROT ) THEN
298 RCONDC = ZERO
299 ELSE
300 *
301 * Compute the 1-norm of A.
302 *
303 ANORM = CLANHE( '1', UPLO, N, A, LDA,
304 + S_WORK_CLANHE )
305 *
306 * Factor the matrix A.
307 *
308 CALL CPOTRF( UPLO, N, A, LDA, INFO )
309 *
310 * Form the inverse of A.
311 *
312 CALL CPOTRI( UPLO, N, A, LDA, INFO )
313 *
314 * Compute the 1-norm condition number of A.
315 *
316 AINVNM = CLANHE( '1', UPLO, N, A, LDA,
317 + S_WORK_CLANHE )
318 RCONDC = ( ONE / ANORM ) / AINVNM
319 *
320 * Restore the matrix A.
321 *
322 CALL CLACPY( UPLO, N, N, ASAV, LDA, A, LDA )
323 *
324 END IF
325 *
326 * Form an exact solution and set the right hand side.
327 *
328 SRNAMT = 'CLARHS'
329 CALL CLARHS( 'CPO', 'N', UPLO, ' ', N, N, KL, KU,
330 + NRHS, A, LDA, XACT, LDA, B, LDA,
331 + ISEED, INFO )
332 CALL CLACPY( 'Full', N, NRHS, B, LDA, BSAV, LDA )
333 *
334 * Compute the L*L' or U'*U factorization of the
335 * matrix and solve the system.
336 *
337 CALL CLACPY( UPLO, N, N, A, LDA, AFAC, LDA )
338 CALL CLACPY( 'Full', N, NRHS, B, LDB, X, LDB )
339 *
340 SRNAMT = 'CTRTTF'
341 CALL CTRTTF( CFORM, UPLO, N, AFAC, LDA, ARF, INFO )
342 SRNAMT = 'CPFTRF'
343 CALL CPFTRF( CFORM, UPLO, N, ARF, INFO )
344 *
345 * Check error code from CPFTRF.
346 *
347 IF( INFO.NE.IZERO ) THEN
348 *
349 * LANGOU: there is a small hick here: IZERO should
350 * always be INFO however if INFO is ZERO, ALAERH does not
351 * complain.
352 *
353 CALL ALAERH( 'CPF', 'CPFSV ', INFO, IZERO,
354 + UPLO, N, N, -1, -1, NRHS, IIT,
355 + NFAIL, NERRS, NOUT )
356 GO TO 100
357 END IF
358 *
359 * Skip the tests if INFO is not 0.
360 *
361 IF( INFO.NE.0 ) THEN
362 GO TO 100
363 END IF
364 *
365 SRNAMT = 'CPFTRS'
366 CALL CPFTRS( CFORM, UPLO, N, NRHS, ARF, X, LDB,
367 + INFO )
368 *
369 SRNAMT = 'CTFTTR'
370 CALL CTFTTR( CFORM, UPLO, N, ARF, AFAC, LDA, INFO )
371 *
372 * Reconstruct matrix from factors and compute
373 * residual.
374 *
375 CALL CLACPY( UPLO, N, N, AFAC, LDA, ASAV, LDA )
376 CALL CPOT01( UPLO, N, A, LDA, AFAC, LDA,
377 + S_WORK_CPOT01, RESULT( 1 ) )
378 CALL CLACPY( UPLO, N, N, ASAV, LDA, AFAC, LDA )
379 *
380 * Form the inverse and compute the residual.
381 *
382 IF(MOD(N,2).EQ.0)THEN
383 CALL CLACPY( 'A', N+1, N/2, ARF, N+1, ARFINV,
384 + N+1 )
385 ELSE
386 CALL CLACPY( 'A', N, (N+1)/2, ARF, N, ARFINV,
387 + N )
388 END IF
389 *
390 SRNAMT = 'CPFTRI'
391 CALL CPFTRI( CFORM, UPLO, N, ARFINV , INFO )
392 *
393 SRNAMT = 'CTFTTR'
394 CALL CTFTTR( CFORM, UPLO, N, ARFINV, AINV, LDA,
395 + INFO )
396 *
397 * Check error code from CPFTRI.
398 *
399 IF( INFO.NE.0 )
400 + CALL ALAERH( 'CPO', 'CPFTRI', INFO, 0, UPLO, N,
401 + N, -1, -1, -1, IMAT, NFAIL, NERRS,
402 + NOUT )
403 *
404 CALL CPOT03( UPLO, N, A, LDA, AINV, LDA,
405 + C_WORK_CPOT03, LDA, S_WORK_CPOT03,
406 + RCONDC, RESULT( 2 ) )
407 *
408 * Compute residual of the computed solution.
409 *
410 CALL CLACPY( 'Full', N, NRHS, B, LDA,
411 + C_WORK_CPOT02, LDA )
412 CALL CPOT02( UPLO, N, NRHS, A, LDA, X, LDA,
413 + C_WORK_CPOT02, LDA, S_WORK_CPOT02,
414 + RESULT( 3 ) )
415 *
416 * Check solution from generated exact solution.
417 *
418 CALL CGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
419 + RESULT( 4 ) )
420 NT = 4
421 *
422 * Print information about the tests that did not
423 * pass the threshold.
424 *
425 DO 60 K = 1, NT
426 IF( RESULT( K ).GE.THRESH ) THEN
427 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
428 + CALL ALADHD( NOUT, 'CPF' )
429 WRITE( NOUT, FMT = 9999 )'CPFSV ', UPLO,
430 + N, IIT, K, RESULT( K )
431 NFAIL = NFAIL + 1
432 END IF
433 60 CONTINUE
434 NRUN = NRUN + NT
435 100 CONTINUE
436 110 CONTINUE
437 120 CONTINUE
438 980 CONTINUE
439 130 CONTINUE
440 *
441 * Print a summary of the results.
442 *
443 CALL ALASVM( 'CPF', NOUT, NFAIL, NRUN, NERRS )
444 *
445 9999 FORMAT( 1X, A6, ', UPLO=''', A1, ''', N =', I5, ', type ', I1,
446 + ', test(', I1, ')=', G12.5 )
447 *
448 RETURN
449 *
450 * End of CDRVRFP
451 *
452 END
2 + THRESH, A, ASAV, AFAC, AINV, B,
3 + BSAV, XACT, X, ARF, ARFINV,
4 + C_WORK_CLATMS, C_WORK_CPOT02,
5 + C_WORK_CPOT03, S_WORK_CLATMS, S_WORK_CLANHE,
6 + S_WORK_CPOT01, S_WORK_CPOT02, S_WORK_CPOT03 )
7 *
8 * -- LAPACK test routine (version 3.2.0) --
9 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
10 * November 2008
11 *
12 * .. Scalar Arguments ..
13 INTEGER NN, NNS, NNT, NOUT
14 REAL THRESH
15 * ..
16 * .. Array Arguments ..
17 INTEGER NVAL( NN ), NSVAL( NNS ), NTVAL( NNT )
18 COMPLEX A( * )
19 COMPLEX AINV( * )
20 COMPLEX ASAV( * )
21 COMPLEX B( * )
22 COMPLEX BSAV( * )
23 COMPLEX AFAC( * )
24 COMPLEX ARF( * )
25 COMPLEX ARFINV( * )
26 COMPLEX XACT( * )
27 COMPLEX X( * )
28 COMPLEX C_WORK_CLATMS( * )
29 COMPLEX C_WORK_CPOT02( * )
30 COMPLEX C_WORK_CPOT03( * )
31 REAL S_WORK_CLATMS( * )
32 REAL S_WORK_CLANHE( * )
33 REAL S_WORK_CPOT01( * )
34 REAL S_WORK_CPOT02( * )
35 REAL S_WORK_CPOT03( * )
36 * ..
37 *
38 * Purpose
39 * =======
40 *
41 * CDRVRFP tests the LAPACK RFP routines:
42 * CPFTRF, CPFTRS, and CPFTRI.
43 *
44 * This testing routine follow the same tests as CDRVPO (test for the full
45 * format Symmetric Positive Definite solver).
46 *
47 * The tests are performed in Full Format, convertion back and forth from
48 * full format to RFP format are performed using the routines CTRTTF and
49 * CTFTTR.
50 *
51 * First, a specific matrix A of size N is created. There is nine types of
52 * different matrixes possible.
53 * 1. Diagonal 6. Random, CNDNUM = sqrt(0.1/EPS)
54 * 2. Random, CNDNUM = 2 7. Random, CNDNUM = 0.1/EPS
55 * *3. First row and column zero 8. Scaled near underflow
56 * *4. Last row and column zero 9. Scaled near overflow
57 * *5. Middle row and column zero
58 * (* - tests error exits from CPFTRF, no test ratios are computed)
59 * A solution XACT of size N-by-NRHS is created and the associated right
60 * hand side B as well. Then CPFTRF is called to compute L (or U), the
61 * Cholesky factor of A. Then L (or U) is used to solve the linear system
62 * of equations AX = B. This gives X. Then L (or U) is used to compute the
63 * inverse of A, AINV. The following four tests are then performed:
64 * (1) norm( L*L' - A ) / ( N * norm(A) * EPS ) or
65 * norm( U'*U - A ) / ( N * norm(A) * EPS ),
66 * (2) norm(B - A*X) / ( norm(A) * norm(X) * EPS ),
67 * (3) norm( I - A*AINV ) / ( N * norm(A) * norm(AINV) * EPS ),
68 * (4) ( norm(X-XACT) * RCOND ) / ( norm(XACT) * EPS ),
69 * where EPS is the machine precision, RCOND the condition number of A, and
70 * norm( . ) the 1-norm for (1,2,3) and the inf-norm for (4).
71 * Errors occur when INFO parameter is not as expected. Failures occur when
72 * a test ratios is greater than THRES.
73 *
74 * Arguments
75 * =========
76 *
77 * NOUT (input) INTEGER
78 * The unit number for output.
79 *
80 * NN (input) INTEGER
81 * The number of values of N contained in the vector NVAL.
82 *
83 * NVAL (input) INTEGER array, dimension (NN)
84 * The values of the matrix dimension N.
85 *
86 * NNS (input) INTEGER
87 * The number of values of NRHS contained in the vector NSVAL.
88 *
89 * NSVAL (input) INTEGER array, dimension (NNS)
90 * The values of the number of right-hand sides NRHS.
91 *
92 * NNT (input) INTEGER
93 * The number of values of MATRIX TYPE contained in the vector NTVAL.
94 *
95 * NTVAL (input) INTEGER array, dimension (NNT)
96 * The values of matrix type (between 0 and 9 for PO/PP/PF matrices).
97 *
98 * THRESH (input) REAL
99 * The threshold value for the test ratios. A result is
100 * included in the output file if RESULT >= THRESH. To have
101 * every test ratio printed, use THRESH = 0.
102 *
103 * A (workspace) COMPLEX array, dimension (NMAX*NMAX)
104 *
105 * ASAV (workspace) COMPLEX array, dimension (NMAX*NMAX)
106 *
107 * AFAC (workspace) COMPLEX array, dimension (NMAX*NMAX)
108 *
109 * AINV (workspace) COMPLEX array, dimension (NMAX*NMAX)
110 *
111 * B (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
112 *
113 * BSAV (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
114 *
115 * XACT (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
116 *
117 * X (workspace) COMPLEX array, dimension (NMAX*MAXRHS)
118 *
119 * ARF (workspace) COMPLEX array, dimension ((NMAX*(NMAX+1))/2)
120 *
121 * ARFINV (workspace) COMPLEX array, dimension ((NMAX*(NMAX+1))/2)
122 *
123 * C_WORK_CLATMS (workspace) COMPLEX array, dimension ( 3*NMAX )
124 *
125 * C_WORK_CPOT02 (workspace) COMPLEX array, dimension ( NMAX*MAXRHS )
126 *
127 * C_WORK_CPOT03 (workspace) COMPLEX array, dimension ( NMAX*NMAX )
128 *
129 * S_WORK_CLATMS (workspace) REAL array, dimension ( NMAX )
130 *
131 * S_WORK_CLANHE (workspace) REAL array, dimension ( NMAX )
132 *
133 * S_WORK_CPOT01 (workspace) REAL array, dimension ( NMAX )
134 *
135 * S_WORK_CPOT02 (workspace) REAL array, dimension ( NMAX )
136 *
137 * S_WORK_CPOT03 (workspace) REAL array, dimension ( NMAX )
138 *
139 * =====================================================================
140 *
141 * .. Parameters ..
142 REAL ONE, ZERO
143 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 )
144 INTEGER NTESTS
145 PARAMETER ( NTESTS = 4 )
146 * ..
147 * .. Local Scalars ..
148 LOGICAL ZEROT
149 INTEGER I, INFO, IUPLO, LDA, LDB, IMAT, NERRS, NFAIL,
150 + NRHS, NRUN, IZERO, IOFF, K, NT, N, IFORM, IIN,
151 + IIT, IIS
152 CHARACTER DIST, CTYPE, UPLO, CFORM
153 INTEGER KL, KU, MODE
154 REAL ANORM, AINVNM, CNDNUM, RCONDC
155 * ..
156 * .. Local Arrays ..
157 CHARACTER UPLOS( 2 ), FORMS( 2 )
158 INTEGER ISEED( 4 ), ISEEDY( 4 )
159 REAL RESULT( NTESTS )
160 * ..
161 * .. External Functions ..
162 REAL CLANHE
163 EXTERNAL CLANHE
164 * ..
165 * .. External Subroutines ..
166 EXTERNAL ALADHD, ALAERH, ALASVM, CGET04, CTFTTR, CLACPY,
167 + CLAIPD, CLARHS, CLATB4, CLATMS, CPFTRI, CPFTRF,
168 + CPFTRS, CPOT01, CPOT02, CPOT03, CPOTRI, CPOTRF,
169 + CTRTTF
170 * ..
171 * .. Scalars in Common ..
172 CHARACTER*32 SRNAMT
173 * ..
174 * .. Common blocks ..
175 COMMON / SRNAMC / SRNAMT
176 * ..
177 * .. Data statements ..
178 DATA ISEEDY / 1988, 1989, 1990, 1991 /
179 DATA UPLOS / 'U', 'L' /
180 DATA FORMS / 'N', 'C' /
181 * ..
182 * .. Executable Statements ..
183 *
184 * Initialize constants and the random number seed.
185 *
186 NRUN = 0
187 NFAIL = 0
188 NERRS = 0
189 DO 10 I = 1, 4
190 ISEED( I ) = ISEEDY( I )
191 10 CONTINUE
192 *
193 DO 130 IIN = 1, NN
194 *
195 N = NVAL( IIN )
196 LDA = MAX( N, 1 )
197 LDB = MAX( N, 1 )
198 *
199 DO 980 IIS = 1, NNS
200 *
201 NRHS = NSVAL( IIS )
202 *
203 DO 120 IIT = 1, NNT
204 *
205 IMAT = NTVAL( IIT )
206 *
207 * If N.EQ.0, only consider the first type
208 *
209 IF( N.EQ.0 .AND. IIT.GT.1 ) GO TO 120
210 *
211 * Skip types 3, 4, or 5 if the matrix size is too small.
212 *
213 IF( IMAT.EQ.4 .AND. N.LE.1 ) GO TO 120
214 IF( IMAT.EQ.5 .AND. N.LE.2 ) GO TO 120
215 *
216 * Do first for UPLO = 'U', then for UPLO = 'L'
217 *
218 DO 110 IUPLO = 1, 2
219 UPLO = UPLOS( IUPLO )
220 *
221 * Do first for CFORM = 'N', then for CFORM = 'C'
222 *
223 DO 100 IFORM = 1, 2
224 CFORM = FORMS( IFORM )
225 *
226 * Set up parameters with CLATB4 and generate a test
227 * matrix with CLATMS.
228 *
229 CALL CLATB4( 'CPO', IMAT, N, N, CTYPE, KL, KU,
230 + ANORM, MODE, CNDNUM, DIST )
231 *
232 SRNAMT = 'CLATMS'
233 CALL CLATMS( N, N, DIST, ISEED, CTYPE,
234 + S_WORK_CLATMS,
235 + MODE, CNDNUM, ANORM, KL, KU, UPLO, A,
236 + LDA, C_WORK_CLATMS, INFO )
237 *
238 * Check error code from CLATMS.
239 *
240 IF( INFO.NE.0 ) THEN
241 CALL ALAERH( 'CPF', 'CLATMS', INFO, 0, UPLO, N,
242 + N, -1, -1, -1, IIT, NFAIL, NERRS,
243 + NOUT )
244 GO TO 100
245 END IF
246 *
247 * For types 3-5, zero one row and column of the matrix to
248 * test that INFO is returned correctly.
249 *
250 ZEROT = IMAT.GE.3 .AND. IMAT.LE.5
251 IF( ZEROT ) THEN
252 IF( IIT.EQ.3 ) THEN
253 IZERO = 1
254 ELSE IF( IIT.EQ.4 ) THEN
255 IZERO = N
256 ELSE
257 IZERO = N / 2 + 1
258 END IF
259 IOFF = ( IZERO-1 )*LDA
260 *
261 * Set row and column IZERO of A to 0.
262 *
263 IF( IUPLO.EQ.1 ) THEN
264 DO 20 I = 1, IZERO - 1
265 A( IOFF+I ) = ZERO
266 20 CONTINUE
267 IOFF = IOFF + IZERO
268 DO 30 I = IZERO, N
269 A( IOFF ) = ZERO
270 IOFF = IOFF + LDA
271 30 CONTINUE
272 ELSE
273 IOFF = IZERO
274 DO 40 I = 1, IZERO - 1
275 A( IOFF ) = ZERO
276 IOFF = IOFF + LDA
277 40 CONTINUE
278 IOFF = IOFF - IZERO
279 DO 50 I = IZERO, N
280 A( IOFF+I ) = ZERO
281 50 CONTINUE
282 END IF
283 ELSE
284 IZERO = 0
285 END IF
286 *
287 * Set the imaginary part of the diagonals.
288 *
289 CALL CLAIPD( N, A, LDA+1, 0 )
290 *
291 * Save a copy of the matrix A in ASAV.
292 *
293 CALL CLACPY( UPLO, N, N, A, LDA, ASAV, LDA )
294 *
295 * Compute the condition number of A (RCONDC).
296 *
297 IF( ZEROT ) THEN
298 RCONDC = ZERO
299 ELSE
300 *
301 * Compute the 1-norm of A.
302 *
303 ANORM = CLANHE( '1', UPLO, N, A, LDA,
304 + S_WORK_CLANHE )
305 *
306 * Factor the matrix A.
307 *
308 CALL CPOTRF( UPLO, N, A, LDA, INFO )
309 *
310 * Form the inverse of A.
311 *
312 CALL CPOTRI( UPLO, N, A, LDA, INFO )
313 *
314 * Compute the 1-norm condition number of A.
315 *
316 AINVNM = CLANHE( '1', UPLO, N, A, LDA,
317 + S_WORK_CLANHE )
318 RCONDC = ( ONE / ANORM ) / AINVNM
319 *
320 * Restore the matrix A.
321 *
322 CALL CLACPY( UPLO, N, N, ASAV, LDA, A, LDA )
323 *
324 END IF
325 *
326 * Form an exact solution and set the right hand side.
327 *
328 SRNAMT = 'CLARHS'
329 CALL CLARHS( 'CPO', 'N', UPLO, ' ', N, N, KL, KU,
330 + NRHS, A, LDA, XACT, LDA, B, LDA,
331 + ISEED, INFO )
332 CALL CLACPY( 'Full', N, NRHS, B, LDA, BSAV, LDA )
333 *
334 * Compute the L*L' or U'*U factorization of the
335 * matrix and solve the system.
336 *
337 CALL CLACPY( UPLO, N, N, A, LDA, AFAC, LDA )
338 CALL CLACPY( 'Full', N, NRHS, B, LDB, X, LDB )
339 *
340 SRNAMT = 'CTRTTF'
341 CALL CTRTTF( CFORM, UPLO, N, AFAC, LDA, ARF, INFO )
342 SRNAMT = 'CPFTRF'
343 CALL CPFTRF( CFORM, UPLO, N, ARF, INFO )
344 *
345 * Check error code from CPFTRF.
346 *
347 IF( INFO.NE.IZERO ) THEN
348 *
349 * LANGOU: there is a small hick here: IZERO should
350 * always be INFO however if INFO is ZERO, ALAERH does not
351 * complain.
352 *
353 CALL ALAERH( 'CPF', 'CPFSV ', INFO, IZERO,
354 + UPLO, N, N, -1, -1, NRHS, IIT,
355 + NFAIL, NERRS, NOUT )
356 GO TO 100
357 END IF
358 *
359 * Skip the tests if INFO is not 0.
360 *
361 IF( INFO.NE.0 ) THEN
362 GO TO 100
363 END IF
364 *
365 SRNAMT = 'CPFTRS'
366 CALL CPFTRS( CFORM, UPLO, N, NRHS, ARF, X, LDB,
367 + INFO )
368 *
369 SRNAMT = 'CTFTTR'
370 CALL CTFTTR( CFORM, UPLO, N, ARF, AFAC, LDA, INFO )
371 *
372 * Reconstruct matrix from factors and compute
373 * residual.
374 *
375 CALL CLACPY( UPLO, N, N, AFAC, LDA, ASAV, LDA )
376 CALL CPOT01( UPLO, N, A, LDA, AFAC, LDA,
377 + S_WORK_CPOT01, RESULT( 1 ) )
378 CALL CLACPY( UPLO, N, N, ASAV, LDA, AFAC, LDA )
379 *
380 * Form the inverse and compute the residual.
381 *
382 IF(MOD(N,2).EQ.0)THEN
383 CALL CLACPY( 'A', N+1, N/2, ARF, N+1, ARFINV,
384 + N+1 )
385 ELSE
386 CALL CLACPY( 'A', N, (N+1)/2, ARF, N, ARFINV,
387 + N )
388 END IF
389 *
390 SRNAMT = 'CPFTRI'
391 CALL CPFTRI( CFORM, UPLO, N, ARFINV , INFO )
392 *
393 SRNAMT = 'CTFTTR'
394 CALL CTFTTR( CFORM, UPLO, N, ARFINV, AINV, LDA,
395 + INFO )
396 *
397 * Check error code from CPFTRI.
398 *
399 IF( INFO.NE.0 )
400 + CALL ALAERH( 'CPO', 'CPFTRI', INFO, 0, UPLO, N,
401 + N, -1, -1, -1, IMAT, NFAIL, NERRS,
402 + NOUT )
403 *
404 CALL CPOT03( UPLO, N, A, LDA, AINV, LDA,
405 + C_WORK_CPOT03, LDA, S_WORK_CPOT03,
406 + RCONDC, RESULT( 2 ) )
407 *
408 * Compute residual of the computed solution.
409 *
410 CALL CLACPY( 'Full', N, NRHS, B, LDA,
411 + C_WORK_CPOT02, LDA )
412 CALL CPOT02( UPLO, N, NRHS, A, LDA, X, LDA,
413 + C_WORK_CPOT02, LDA, S_WORK_CPOT02,
414 + RESULT( 3 ) )
415 *
416 * Check solution from generated exact solution.
417 *
418 CALL CGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
419 + RESULT( 4 ) )
420 NT = 4
421 *
422 * Print information about the tests that did not
423 * pass the threshold.
424 *
425 DO 60 K = 1, NT
426 IF( RESULT( K ).GE.THRESH ) THEN
427 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
428 + CALL ALADHD( NOUT, 'CPF' )
429 WRITE( NOUT, FMT = 9999 )'CPFSV ', UPLO,
430 + N, IIT, K, RESULT( K )
431 NFAIL = NFAIL + 1
432 END IF
433 60 CONTINUE
434 NRUN = NRUN + NT
435 100 CONTINUE
436 110 CONTINUE
437 120 CONTINUE
438 980 CONTINUE
439 130 CONTINUE
440 *
441 * Print a summary of the results.
442 *
443 CALL ALASVM( 'CPF', NOUT, NFAIL, NRUN, NERRS )
444 *
445 9999 FORMAT( 1X, A6, ', UPLO=''', A1, ''', N =', I5, ', type ', I1,
446 + ', test(', I1, ')=', G12.5 )
447 *
448 RETURN
449 *
450 * End of CDRVRFP
451 *
452 END