PLASMA  2.4.5
PLASMA - Parallel Linear Algebra for Scalable Multi-core Architectures
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups
cchkpo.f
Go to the documentation of this file.
1  SUBROUTINE cchkpo( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL,
2  $ thresh, tsterr, nmax, a, afac, ainv, b, x,
3  $ xact, work, rwork, nout )
4 *
5  include 'plasmaf.h'
6 *
7 * -- LAPACK test routine (version 3.1) --
8 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
9 * November 2006
10 *
11 * .. Scalar Arguments ..
12  LOGICAL tsterr
13  INTEGER nmax, nn, nnb, nns, nout
14  REAL thresh
15 * ..
16 * .. Array Arguments ..
17  LOGICAL dotype( * )
18  INTEGER nbval( * ), nsval( * ), nval( * )
19  REAL rwork( * )
20  COMPLEX a( * ), afac( * ), ainv( * ), b( * ),
21  $ work( * ), x( * ), xact( * )
22 * ..
23 *
24 * Purpose
25 * =======
26 *
27 * CCHKPO tests CPOTRF, -TRI, -TRS, -RFS, and -CON
28 *
29 * Arguments
30 * =========
31 *
32 * DOTYPE (input) LOGICAL array, dimension (NTYPES)
33 * The matrix types to be used for testing. Matrices of type j
34 * (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
35 * .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
36 *
37 * NN (input) INTEGER
38 * The number of values of N contained in the vector NVAL.
39 *
40 * NVAL (input) INTEGER array, dimension (NN)
41 * The values of the matrix dimension N.
42 *
43 * NNB (input) INTEGER
44 * The number of values of NB contained in the vector NBVAL.
45 *
46 * NBVAL (input) INTEGER array, dimension (NBVAL)
47 * The values of the blocksize NB.
48 *
49 * NNS (input) INTEGER
50 * The number of values of NRHS contained in the vector NSVAL.
51 *
52 * NSVAL (input) INTEGER array, dimension (NNS)
53 * The values of the number of right hand sides NRHS.
54 *
55 * THRESH (input) REAL
56 * The threshold value for the test ratios. A result is
57 * included in the output file if RESULT >= THRESH. To have
58 * every test ratio printed, use THRESH = 0.
59 *
60 * TSTERR (input) LOGICAL
61 * Flag that indicates whether error exits are to be tested.
62 *
63 * NMAX (input) INTEGER
64 * The maximum value permitted for N, used in dimensioning the
65 * work arrays.
66 *
67 * A (workspace) COMPLEX array, dimension (NMAX*NMAX)
68 *
69 * AFAC (workspace) COMPLEX array, dimension (NMAX*NMAX)
70 *
71 * AINV (workspace) COMPLEX array, dimension (NMAX*NMAX)
72 *
73 * B (workspace) COMPLEX array, dimension (NMAX*NSMAX)
74 * where NSMAX is the largest entry in NSVAL.
75 *
76 * X (workspace) COMPLEX array, dimension (NMAX*NSMAX)
77 *
78 * XACT (workspace) COMPLEX array, dimension (NMAX*NSMAX)
79 *
80 * WORK (workspace) COMPLEX array, dimension
81 * (NMAX*max(3,NSMAX))
82 *
83 * RWORK (workspace) REAL array, dimension
84 * (NMAX+2*NSMAX)
85 *
86 * NOUT (input) INTEGER
87 * The unit number for output.
88 *
89 * =====================================================================
90 *
91 * .. Parameters ..
92  COMPLEX czero
93  parameter( czero = ( 0.0e+0, 0.0e+0 ) )
94  INTEGER ntypes
95  parameter( ntypes = 9 )
96  INTEGER ntests
97  parameter( ntests = 8 )
98 * ..
99 * .. Local Scalars ..
100  LOGICAL zerot
101  CHARACTER dist, type, uplo, xtype
102  CHARACTER*3 path
103  INTEGER i, imat, in, inb, info, ioff, irhs, iuplo,
104  $ izero, k, kl, ku, lda, mode, n, nb, nerrs,
105  $ nfail, nimat, nrhs, nrun, plasma_uplo
106  REAL anorm, cndnum, rcond, rcondc
107 * ..
108 * .. Local Arrays ..
109  CHARACTER uplos( 2 )
110  INTEGER plasma_uplos( 2 )
111  INTEGER iseed( 4 ), iseedy( 4 )
112  REAL result( ntests )
113 * ..
114 * .. External Functions ..
115  REAL clanhe, sget06
116  EXTERNAL clanhe, sget06
117 * ..
118 * .. External Subroutines ..
119  EXTERNAL alaerh, alahd, alasum, cerrpo, cget04, clacpy,
121  $ cpot01, cpot02, cpot03, cpot05, cpotrf, cpotri,
122  $ cpotrs, xlaenv
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 * .. Intrinsic Functions ..
134  INTRINSIC max
135 * ..
136 * .. Data statements ..
137  DATA iseedy / 1988, 1989, 1990, 1991 /
138  DATA uplos / 'U', 'L' /
139  DATA plasma_uplos / plasmaupper, plasmalower /
140 * ..
141 * .. Executable Statements ..
142 *
143 * Initialize constants and the random number seed.
144 *
145  path( 1: 1 ) = 'Complex precision'
146  path( 2: 3 ) = 'PO'
147  nrun = 0
148  nfail = 0
149  nerrs = 0
150  DO 10 i = 1, 4
151  iseed( i ) = iseedy( i )
152  10 continue
153 *
154 * Test the error exits
155 *
156  IF( tsterr )
157  $ CALL cerrpo( path, nout )
158  infot = 0
159 *
160 * Do for each value of N in NVAL
161 *
162  DO 120 in = 1, nn
163  n = nval( in )
164  lda = max( n, 1 )
165  xtype = 'N'
166  nimat = ntypes
167  IF( n.LE.0 )
168  $ nimat = 1
169 *
170  izero = 0
171  DO 110 imat = 1, nimat
172 *
173 * Do the tests only if DOTYPE( IMAT ) is true.
174 *
175  IF( .NOT.dotype( imat ) )
176  $ go to 110
177 *
178 * Skip types 3, 4, or 5 if the matrix size is too small.
179 *
180  zerot = imat.GE.3 .AND. imat.LE.5
181  IF( zerot .AND. n.LT.imat-2 )
182  $ go to 110
183 *
184 * Do first for UPLO = 'U', then for UPLO = 'L'
185 *
186  DO 100 iuplo = 1, 2
187  uplo = uplos( iuplo )
188  plasma_uplo = plasma_uplos( iuplo )
189 *
190 * Set up parameters with CLATB4 and generate a test matrix
191 * with CLATMS.
192 *
193  CALL clatb4( path, imat, n, n, type, kl, ku, anorm, mode,
194  $ cndnum, dist )
195 *
196  srnamt = 'CLATMS'
197  CALL clatms( n, n, dist, iseed, type, rwork, mode,
198  $ cndnum, anorm, kl, ku, uplo, a, lda, work,
199  $ info )
200 *
201 * Check error code from CLATMS.
202 *
203  IF( info.NE.0 ) THEN
204  CALL alaerh( path, 'CLATMS', info, 0, uplo, n, n, -1,
205  $ -1, -1, imat, nfail, nerrs, nout )
206  go to 100
207  END IF
208 *
209 * For types 3-5, zero one row and column of the matrix to
210 * test that INFO is returned correctly.
211 *
212  IF( zerot ) THEN
213  IF( imat.EQ.3 ) THEN
214  izero = 1
215  ELSE IF( imat.EQ.4 ) THEN
216  izero = n
217  ELSE
218  izero = n / 2 + 1
219  END IF
220  ioff = ( izero-1 )*lda
221 *
222 * Set row and column IZERO of A to 0.
223 *
224  IF( iuplo.EQ.1 ) THEN
225  DO 20 i = 1, izero - 1
226  a( ioff+i ) = czero
227  20 continue
228  ioff = ioff + izero
229  DO 30 i = izero, n
230  a( ioff ) = czero
231  ioff = ioff + lda
232  30 continue
233  ELSE
234  ioff = izero
235  DO 40 i = 1, izero - 1
236  a( ioff ) = czero
237  ioff = ioff + lda
238  40 continue
239  ioff = ioff - izero
240  DO 50 i = izero, n
241  a( ioff+i ) = czero
242  50 continue
243  END IF
244  ELSE
245  izero = 0
246  END IF
247 *
248 * Set the imaginary part of the diagonals.
249 *
250  CALL claipd( n, a, lda+1, 0 )
251 *
252 * Do for each value of NB in NBVAL
253 *
254  DO 90 inb = 1, nnb
255  nb = nbval( inb )
256  CALL xlaenv( 1, nb )
257  IF ( (n / 25) .GT. nb ) THEN
258  goto 90
259  END IF
260  CALL plasma_set( plasma_tile_size, nb, info )
261 *
262 * Compute the L*L' or U'*U factorization of the matrix.
263 *
264  CALL clacpy( uplo, n, n, a, lda, afac, lda )
265  srnamt = 'CPOTRF'
266  CALL plasma_cpotrf( plasma_uplo, n, afac, lda, info )
267 *
268 * Check error code from CPOTRF.
269 *
270  IF( info.NE.izero ) THEN
271  CALL alaerh( path, 'CPOTRF', info, izero, uplo, n,
272  $ n, -1, -1, nb, imat, nfail, nerrs,
273  $ nout )
274  go to 90
275  END IF
276 *
277 * Skip the tests if INFO is not 0.
278 *
279  IF( info.NE.0 )
280  $ go to 90
281 *
282 *+ TEST 1
283 * Reconstruct matrix from factors and compute residual.
284 *
285  CALL clacpy( uplo, n, n, afac, lda, ainv, lda )
286  CALL cpot01( uplo, n, a, lda, ainv, lda, rwork,
287  $ result( 1 ) )
288 *
289 *+ TEST 2
290 * Form the inverse and compute the residual.
291 *
292  CALL clacpy( uplo, n, n, afac, lda, ainv, lda )
293  srnamt = 'CPOTRI'
294  CALL plasma_cpotri( plasma_uplo, n, ainv, lda,
295  $ info )
296 *
297 * Check error code from CPOTRI.
298 *
299  IF( info.NE.0 )
300  $ CALL alaerh( path, 'CPOTRI', info, 0, uplo, n, n,
301  $ -1, -1, -1, imat, nfail, nerrs, nout )
302 *
303  CALL cpot03( uplo, n, a, lda, ainv, lda, work, lda,
304  $ rwork, rcondc, result( 2 ) )
305 *
306 * Print information about the tests that did not pass
307 * the threshold.
308 *
309  DO 60 k = 1, 2
310  IF( result( k ).GE.thresh ) THEN
311  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
312  $ CALL alahd( nout, path )
313  WRITE( nout, fmt = 9999 )uplo, n, nb, imat, k,
314  $ result( k )
315  nfail = nfail + 1
316  END IF
317  60 continue
318  nrun = nrun + 2
319 *
320 * Skip the rest of the tests unless this is the first
321 * blocksize.
322 *
323  IF( inb.NE.1 )
324  $ go to 90
325 *
326  DO 80 irhs = 1, nns
327  nrhs = nsval( irhs )
328 *
329 *+ TEST 3
330 * Solve and compute residual for A * X = B .
331 *
332  srnamt = 'CLARHS'
333  CALL clarhs( path, xtype, uplo, ' ', n, n, kl, ku,
334  $ nrhs, a, lda, xact, lda, b, lda,
335  $ iseed, info )
336  CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
337 *
338  srnamt = 'CPOTRS'
339  CALL plasma_cpotrs( plasma_uplo, n, nrhs, afac,
340  $ lda, x, lda, info )
341 *
342 * Check error code from CPOTRS.
343 *
344  IF( info.NE.0 )
345  $ CALL alaerh( path, 'CPOTRS', info, 0, uplo, n,
346  $ n, -1, -1, nrhs, imat, nfail,
347  $ nerrs, nout )
348 *
349  CALL clacpy( 'Full', n, nrhs, b, lda, work, lda )
350  CALL cpot02( uplo, n, nrhs, a, lda, x, lda, work,
351  $ lda, rwork, result( 3 ) )
352 *
353 *+ TEST 4
354 * Check solution from generated exact solution.
355 *
356  CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
357  $ result( 4 ) )
358 *
359 *+ TESTS 5, 6, and 7
360 * Use iterative refinement to improve the solution.
361 *
362  srnamt = 'CPORFS'
363  CALL cporfs( uplo, n, nrhs, a, lda, afac, lda, b,
364  $ lda, x, lda, rwork, rwork( nrhs+1 ),
365  $ work, rwork( 2*nrhs+1 ), info )
366 *
367 * Check error code from CPORFS.
368 *
369  IF( info.NE.0 )
370  $ CALL alaerh( path, 'CPORFS', info, 0, uplo, n,
371  $ n, -1, -1, nrhs, imat, nfail,
372  $ nerrs, nout )
373 *
374  CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
375  $ result( 5 ) )
376  CALL cpot05( uplo, n, nrhs, a, lda, b, lda, x, lda,
377  $ xact, lda, rwork, rwork( nrhs+1 ),
378  $ result( 6 ) )
379 *
380 * Print information about the tests that did not pass
381 * the threshold.
382 *
383  DO 70 k = 3, 7
384  IF( result( k ).GE.thresh ) THEN
385  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
386  $ CALL alahd( nout, path )
387  WRITE( nout, fmt = 9998 )uplo, n, nrhs,
388  $ imat, k, result( k )
389  nfail = nfail + 1
390  END IF
391  70 continue
392  nrun = nrun + 5
393  80 continue
394 *
395 *+ TEST 8
396 * Get an estimate of RCOND = 1/CNDNUM.
397 *
398  anorm = clanhe( '1', uplo, n, a, lda, rwork )
399  srnamt = 'CPOCON'
400  CALL cpocon( uplo, n, afac, lda, anorm, rcond, work,
401  $ rwork, info )
402 *
403 * Check error code from CPOCON.
404 *
405  IF( info.NE.0 )
406  $ CALL alaerh( path, 'CPOCON', info, 0, uplo, n, n,
407  $ -1, -1, -1, imat, nfail, nerrs, nout )
408 *
409  result( 8 ) = sget06( rcond, rcondc )
410 *
411 * Print the test ratio if it is .GE. THRESH.
412 *
413  IF( result( 8 ).GE.thresh ) THEN
414  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
415  $ CALL alahd( nout, path )
416  WRITE( nout, fmt = 9997 )uplo, n, imat, 8,
417  $ result( 8 )
418  nfail = nfail + 1
419  END IF
420  nrun = nrun + 1
421  90 continue
422  100 continue
423  110 continue
424  120 continue
425 *
426 * Print a summary of the results.
427 *
428  CALL alasum( path, nout, nfail, nrun, nerrs )
429 *
430  9999 format( ' UPLO = ''', a1, ''', N =', i5, ', NB =', i4, ', type ',
431  $ i2, ', test ', i2, ', ratio =', g12.5 )
432  9998 format( ' UPLO = ''', a1, ''', N =', i5, ', NRHS=', i3, ', type ',
433  $ i2, ', test(', i2, ') =', g12.5 )
434  9997 format( ' UPLO = ''', a1, ''', N =', i5, ',', 10x, ' type ', i2,
435  $ ', test(', i2, ') =', g12.5 )
436  return
437 *
438 * End of CCHKPO
439 *
440  END