LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
cchkaa.f
Go to the documentation of this file.
1 *> \brief \b CCHKAA
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * PROGRAM CCHKAA
12 *
13 *
14 *> \par Purpose:
15 * =============
16 *>
17 *> \verbatim
18 *>
19 *> CCHKAA is the main test program for the COMPLEX linear equation
20 *> routines.
21 *>
22 *> The program must be driven by a short data file. The first 15 records
23 *> (not including the first comment line) specify problem dimensions
24 *> and program options using list-directed input. The remaining lines
25 *> specify the LAPACK test paths and the number of matrix types to use
26 *> in testing. An annotated example of a data file can be obtained by
27 *> deleting the first 3 characters from the following 42 lines:
28 *> Data file for testing COMPLEX LAPACK linear equation routines
29 *> 7 Number of values of M
30 *> 0 1 2 3 5 10 16 Values of M (row dimension)
31 *> 7 Number of values of N
32 *> 0 1 2 3 5 10 16 Values of N (column dimension)
33 *> 1 Number of values of NRHS
34 *> 2 Values of NRHS (number of right hand sides)
35 *> 5 Number of values of NB
36 *> 1 3 3 3 20 Values of NB (the blocksize)
37 *> 1 0 5 9 1 Values of NX (crossover point)
38 *> 3 Number of values of RANK
39 *> 30 50 90 Values of rank (as a % of N)
40 *> 30.0 Threshold value of test ratio
41 *> T Put T to test the LAPACK routines
42 *> T Put T to test the driver routines
43 *> T Put T to test the error exits
44 *> CGE 11 List types on next line if 0 < NTYPES < 11
45 *> CGB 8 List types on next line if 0 < NTYPES < 8
46 *> CGT 12 List types on next line if 0 < NTYPES < 12
47 *> CPO 9 List types on next line if 0 < NTYPES < 9
48 *> CPO 9 List types on next line if 0 < NTYPES < 9
49 *> CPP 9 List types on next line if 0 < NTYPES < 9
50 *> CPB 8 List types on next line if 0 < NTYPES < 8
51 *> CPT 12 List types on next line if 0 < NTYPES < 12
52 *> CHE 10 List types on next line if 0 < NTYPES < 10
53 *> CHR 10 List types on next line if 0 < NTYPES < 10
54 *> CHP 10 List types on next line if 0 < NTYPES < 10
55 *> CSY 11 List types on next line if 0 < NTYPES < 11
56 *> CSR 11 List types on next line if 0 < NTYPES < 11
57 *> CSP 11 List types on next line if 0 < NTYPES < 11
58 *> CTR 18 List types on next line if 0 < NTYPES < 18
59 *> CTP 18 List types on next line if 0 < NTYPES < 18
60 *> CTB 17 List types on next line if 0 < NTYPES < 17
61 *> CQR 8 List types on next line if 0 < NTYPES < 8
62 *> CRQ 8 List types on next line if 0 < NTYPES < 8
63 *> CLQ 8 List types on next line if 0 < NTYPES < 8
64 *> CQL 8 List types on next line if 0 < NTYPES < 8
65 *> CQP 6 List types on next line if 0 < NTYPES < 6
66 *> CTZ 3 List types on next line if 0 < NTYPES < 3
67 *> CLS 6 List types on next line if 0 < NTYPES < 6
68 *> CEQ
69 *> CQT
70 *> CQX
71 *> \endverbatim
72 *
73 * Parameters:
74 * ==========
75 *
76 *> \verbatim
77 *> NMAX INTEGER
78 *> The maximum allowable value for M and N.
79 *>
80 *> MAXIN INTEGER
81 *> The number of different values that can be used for each of
82 *> M, N, NRHS, NB, NX and RANK
83 *>
84 *> MAXRHS INTEGER
85 *> The maximum number of right hand sides
86 *>
87 *> MATMAX INTEGER
88 *> The maximum number of matrix types to use for testing
89 *>
90 *> NIN INTEGER
91 *> The unit number for input
92 *>
93 *> NOUT INTEGER
94 *> The unit number for output
95 *> \endverbatim
96 *
97 * Authors:
98 * ========
99 *
100 *> \author Univ. of Tennessee
101 *> \author Univ. of California Berkeley
102 *> \author Univ. of Colorado Denver
103 *> \author NAG Ltd.
104 *
105 *> \date November 2015
106 *
107 *> \ingroup complex_lin
108 *
109 * =====================================================================
110  PROGRAM cchkaa
111 *
112 * -- LAPACK test routine (version 3.6.0) --
113 * -- LAPACK is a software package provided by Univ. of Tennessee, --
114 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
115 * November 2015
116 *
117 * =====================================================================
118 *
119 * .. Parameters ..
120  INTEGER NMAX
121  parameter ( nmax = 132 )
122  INTEGER MAXIN
123  parameter ( maxin = 12 )
124  INTEGER MAXRHS
125  parameter ( maxrhs = 16 )
126  INTEGER MATMAX
127  parameter ( matmax = 30 )
128  INTEGER NIN, NOUT
129  parameter ( nin = 5, nout = 6 )
130  INTEGER KDMAX
131  parameter ( kdmax = nmax+( nmax+1 ) / 4 )
132 * ..
133 * .. Local Scalars ..
134  LOGICAL FATAL, TSTCHK, TSTDRV, TSTERR
135  CHARACTER C1
136  CHARACTER*2 C2
137  CHARACTER*3 PATH
138  CHARACTER*10 INTSTR
139  CHARACTER*72 ALINE
140  INTEGER I, IC, J, K, LA, LAFAC, LDA, NB, NM, NMATS, NN,
141  $ nnb, nnb2, nns, nrhs, ntypes, nrank,
142  $ vers_major, vers_minor, vers_patch
143  REAL EPS, S1, S2, THREQ, THRESH
144 * ..
145 * .. Local Arrays ..
146  LOGICAL DOTYPE( matmax )
147  INTEGER IWORK( 25*nmax ), MVAL( maxin ),
148  $ nbval( maxin ), nbval2( maxin ),
149  $ nsval( maxin ), nval( maxin ), nxval( maxin ),
150  $ rankval( maxin ), piv( nmax )
151  REAL RWORK( 150*nmax+2*maxrhs ), S( 2*nmax )
152  COMPLEX A( ( kdmax+1 )*nmax, 7 ), B( nmax*maxrhs, 4 ),
153  $ work( nmax, nmax+maxrhs+10 )
154 * ..
155 * .. External Functions ..
156  LOGICAL LSAME, LSAMEN
157  REAL SECOND, SLAMCH
158  EXTERNAL lsame, lsamen, second, slamch
159 * ..
160 * .. External Subroutines ..
161  EXTERNAL alareq, cchkeq, cchkgb, cchkge, cchkgt, cchkhe,
169 
170 * ..
171 * .. Scalars in Common ..
172  LOGICAL LERR, OK
173  CHARACTER*32 SRNAMT
174  INTEGER INFOT, NUNIT
175 * ..
176 * .. Arrays in Common ..
177  INTEGER IPARMS( 100 )
178 * ..
179 * .. Common blocks ..
180  COMMON / claenv / iparms
181  COMMON / infoc / infot, nunit, ok, lerr
182  COMMON / srnamc / srnamt
183 * ..
184 * .. Data statements ..
185  DATA threq / 2.0 / , intstr / '0123456789' /
186 * ..
187 * .. Executable Statements ..
188 *
189  s1 = second( )
190  lda = nmax
191  fatal = .false.
192 *
193 * Read a dummy line.
194 *
195  READ( nin, fmt = * )
196 *
197 * Report values of parameters.
198 *
199  CALL ilaver( vers_major, vers_minor, vers_patch )
200  WRITE( nout, fmt = 9994 ) vers_major, vers_minor, vers_patch
201 *
202 * Read the values of M
203 *
204  READ( nin, fmt = * )nm
205  IF( nm.LT.1 ) THEN
206  WRITE( nout, fmt = 9996 )' NM ', nm, 1
207  nm = 0
208  fatal = .true.
209  ELSE IF( nm.GT.maxin ) THEN
210  WRITE( nout, fmt = 9995 )' NM ', nm, maxin
211  nm = 0
212  fatal = .true.
213  END IF
214  READ( nin, fmt = * )( mval( i ), i = 1, nm )
215  DO 10 i = 1, nm
216  IF( mval( i ).LT.0 ) THEN
217  WRITE( nout, fmt = 9996 )' M ', mval( i ), 0
218  fatal = .true.
219  ELSE IF( mval( i ).GT.nmax ) THEN
220  WRITE( nout, fmt = 9995 )' M ', mval( i ), nmax
221  fatal = .true.
222  END IF
223  10 CONTINUE
224  IF( nm.GT.0 )
225  $ WRITE( nout, fmt = 9993 )'M ', ( mval( i ), i = 1, nm )
226 *
227 * Read the values of N
228 *
229  READ( nin, fmt = * )nn
230  IF( nn.LT.1 ) THEN
231  WRITE( nout, fmt = 9996 )' NN ', nn, 1
232  nn = 0
233  fatal = .true.
234  ELSE IF( nn.GT.maxin ) THEN
235  WRITE( nout, fmt = 9995 )' NN ', nn, maxin
236  nn = 0
237  fatal = .true.
238  END IF
239  READ( nin, fmt = * )( nval( i ), i = 1, nn )
240  DO 20 i = 1, nn
241  IF( nval( i ).LT.0 ) THEN
242  WRITE( nout, fmt = 9996 )' N ', nval( i ), 0
243  fatal = .true.
244  ELSE IF( nval( i ).GT.nmax ) THEN
245  WRITE( nout, fmt = 9995 )' N ', nval( i ), nmax
246  fatal = .true.
247  END IF
248  20 CONTINUE
249  IF( nn.GT.0 )
250  $ WRITE( nout, fmt = 9993 )'N ', ( nval( i ), i = 1, nn )
251 *
252 * Read the values of NRHS
253 *
254  READ( nin, fmt = * )nns
255  IF( nns.LT.1 ) THEN
256  WRITE( nout, fmt = 9996 )' NNS', nns, 1
257  nns = 0
258  fatal = .true.
259  ELSE IF( nns.GT.maxin ) THEN
260  WRITE( nout, fmt = 9995 )' NNS', nns, maxin
261  nns = 0
262  fatal = .true.
263  END IF
264  READ( nin, fmt = * )( nsval( i ), i = 1, nns )
265  DO 30 i = 1, nns
266  IF( nsval( i ).LT.0 ) THEN
267  WRITE( nout, fmt = 9996 )'NRHS', nsval( i ), 0
268  fatal = .true.
269  ELSE IF( nsval( i ).GT.maxrhs ) THEN
270  WRITE( nout, fmt = 9995 )'NRHS', nsval( i ), maxrhs
271  fatal = .true.
272  END IF
273  30 CONTINUE
274  IF( nns.GT.0 )
275  $ WRITE( nout, fmt = 9993 )'NRHS', ( nsval( i ), i = 1, nns )
276 *
277 * Read the values of NB
278 *
279  READ( nin, fmt = * )nnb
280  IF( nnb.LT.1 ) THEN
281  WRITE( nout, fmt = 9996 )'NNB ', nnb, 1
282  nnb = 0
283  fatal = .true.
284  ELSE IF( nnb.GT.maxin ) THEN
285  WRITE( nout, fmt = 9995 )'NNB ', nnb, maxin
286  nnb = 0
287  fatal = .true.
288  END IF
289  READ( nin, fmt = * )( nbval( i ), i = 1, nnb )
290  DO 40 i = 1, nnb
291  IF( nbval( i ).LT.0 ) THEN
292  WRITE( nout, fmt = 9996 )' NB ', nbval( i ), 0
293  fatal = .true.
294  END IF
295  40 CONTINUE
296  IF( nnb.GT.0 )
297  $ WRITE( nout, fmt = 9993 )'NB ', ( nbval( i ), i = 1, nnb )
298 *
299 * Set NBVAL2 to be the set of unique values of NB
300 *
301  nnb2 = 0
302  DO 60 i = 1, nnb
303  nb = nbval( i )
304  DO 50 j = 1, nnb2
305  IF( nb.EQ.nbval2( j ) )
306  $ GO TO 60
307  50 CONTINUE
308  nnb2 = nnb2 + 1
309  nbval2( nnb2 ) = nb
310  60 CONTINUE
311 *
312 * Read the values of NX
313 *
314  READ( nin, fmt = * )( nxval( i ), i = 1, nnb )
315  DO 70 i = 1, nnb
316  IF( nxval( i ).LT.0 ) THEN
317  WRITE( nout, fmt = 9996 )' NX ', nxval( i ), 0
318  fatal = .true.
319  END IF
320  70 CONTINUE
321  IF( nnb.GT.0 )
322  $ WRITE( nout, fmt = 9993 )'NX ', ( nxval( i ), i = 1, nnb )
323 *
324 * Read the values of RANKVAL
325 *
326  READ( nin, fmt = * )nrank
327  IF( nn.LT.1 ) THEN
328  WRITE( nout, fmt = 9996 )' NRANK ', nrank, 1
329  nrank = 0
330  fatal = .true.
331  ELSE IF( nn.GT.maxin ) THEN
332  WRITE( nout, fmt = 9995 )' NRANK ', nrank, maxin
333  nrank = 0
334  fatal = .true.
335  END IF
336  READ( nin, fmt = * )( rankval( i ), i = 1, nrank )
337  DO i = 1, nrank
338  IF( rankval( i ).LT.0 ) THEN
339  WRITE( nout, fmt = 9996 )' RANK ', rankval( i ), 0
340  fatal = .true.
341  ELSE IF( rankval( i ).GT.100 ) THEN
342  WRITE( nout, fmt = 9995 )' RANK ', rankval( i ), 100
343  fatal = .true.
344  END IF
345  END DO
346  IF( nrank.GT.0 )
347  $ WRITE( nout, fmt = 9993 )'RANK % OF N',
348  $ ( rankval( i ), i = 1, nrank )
349 *
350 * Read the threshold value for the test ratios.
351 *
352  READ( nin, fmt = * )thresh
353  WRITE( nout, fmt = 9992 )thresh
354 *
355 * Read the flag that indicates whether to test the LAPACK routines.
356 *
357  READ( nin, fmt = * )tstchk
358 *
359 * Read the flag that indicates whether to test the driver routines.
360 *
361  READ( nin, fmt = * )tstdrv
362 *
363 * Read the flag that indicates whether to test the error exits.
364 *
365  READ( nin, fmt = * )tsterr
366 *
367  IF( fatal ) THEN
368  WRITE( nout, fmt = 9999 )
369  stop
370  END IF
371 *
372 * Calculate and print the machine dependent constants.
373 *
374  eps = slamch( 'Underflow threshold' )
375  WRITE( nout, fmt = 9991 )'underflow', eps
376  eps = slamch( 'Overflow threshold' )
377  WRITE( nout, fmt = 9991 )'overflow ', eps
378  eps = slamch( 'Epsilon' )
379  WRITE( nout, fmt = 9991 )'precision', eps
380  WRITE( nout, fmt = * )
381  nrhs = nsval( 1 )
382 *
383  80 CONTINUE
384 *
385 * Read a test path and the number of matrix types to use.
386 *
387  READ( nin, fmt = '(A72)', end = 140 )aline
388  path = aline( 1: 3 )
389  nmats = matmax
390  i = 3
391  90 CONTINUE
392  i = i + 1
393  IF( i.GT.72 )
394  $ GO TO 130
395  IF( aline( i: i ).EQ.' ' )
396  $ GO TO 90
397  nmats = 0
398  100 CONTINUE
399  c1 = aline( i: i )
400  DO 110 k = 1, 10
401  IF( c1.EQ.intstr( k: k ) ) THEN
402  ic = k - 1
403  GO TO 120
404  END IF
405  110 CONTINUE
406  GO TO 130
407  120 CONTINUE
408  nmats = nmats*10 + ic
409  i = i + 1
410  IF( i.GT.72 )
411  $ GO TO 130
412  GO TO 100
413  130 CONTINUE
414  c1 = path( 1: 1 )
415  c2 = path( 2: 3 )
416 *
417 * Check first character for correct precision.
418 *
419  IF( .NOT.lsame( c1, 'Complex precision' ) ) THEN
420  WRITE( nout, fmt = 9990 )path
421 *
422  ELSE IF( nmats.LE.0 ) THEN
423 *
424 * Check for a positive number of tests requested.
425 *
426  WRITE( nout, fmt = 9989 )path
427 *
428  ELSE IF( lsamen( 2, c2, 'GE' ) ) THEN
429 *
430 * GE: general matrices
431 *
432  ntypes = 11
433  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
434 *
435  IF( tstchk ) THEN
436  CALL cchkge( dotype, nm, mval, nn, nval, nnb2, nbval2, nns,
437  $ nsval, thresh, tsterr, lda, a( 1, 1 ),
438  $ a( 1, 2 ), a( 1, 3 ), b( 1, 1 ), b( 1, 2 ),
439  $ b( 1, 3 ), work, rwork, iwork, nout )
440  ELSE
441  WRITE( nout, fmt = 9989 )path
442  END IF
443 *
444  IF( tstdrv ) THEN
445  CALL cdrvge( dotype, nn, nval, nrhs, thresh, tsterr, lda,
446  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
447  $ b( 1, 2 ), b( 1, 3 ), b( 1, 4 ), s, work,
448  $ rwork, iwork, nout )
449  ELSE
450  WRITE( nout, fmt = 9988 )path
451  END IF
452 *
453  ELSE IF( lsamen( 2, c2, 'GB' ) ) THEN
454 *
455 * GB: general banded matrices
456 *
457  la = ( 2*kdmax+1 )*nmax
458  lafac = ( 3*kdmax+1 )*nmax
459  ntypes = 8
460  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
461 *
462  IF( tstchk ) THEN
463  CALL cchkgb( dotype, nm, mval, nn, nval, nnb2, nbval2, nns,
464  $ nsval, thresh, tsterr, a( 1, 1 ), la,
465  $ a( 1, 3 ), lafac, b( 1, 1 ), b( 1, 2 ),
466  $ b( 1, 3 ), work, rwork, iwork, nout )
467  ELSE
468  WRITE( nout, fmt = 9989 )path
469  END IF
470 *
471  IF( tstdrv ) THEN
472  CALL cdrvgb( dotype, nn, nval, nrhs, thresh, tsterr,
473  $ a( 1, 1 ), la, a( 1, 3 ), lafac, a( 1, 6 ),
474  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), b( 1, 4 ), s,
475  $ work, rwork, iwork, nout )
476  ELSE
477  WRITE( nout, fmt = 9988 )path
478  END IF
479 *
480  ELSE IF( lsamen( 2, c2, 'GT' ) ) THEN
481 *
482 * GT: general tridiagonal matrices
483 *
484  ntypes = 12
485  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
486 *
487  IF( tstchk ) THEN
488  CALL cchkgt( dotype, nn, nval, nns, nsval, thresh, tsterr,
489  $ a( 1, 1 ), a( 1, 2 ), b( 1, 1 ), b( 1, 2 ),
490  $ b( 1, 3 ), work, rwork, iwork, nout )
491  ELSE
492  WRITE( nout, fmt = 9989 )path
493  END IF
494 *
495  IF( tstdrv ) THEN
496  CALL cdrvgt( dotype, nn, nval, nrhs, thresh, tsterr,
497  $ a( 1, 1 ), a( 1, 2 ), b( 1, 1 ), b( 1, 2 ),
498  $ b( 1, 3 ), work, rwork, iwork, nout )
499  ELSE
500  WRITE( nout, fmt = 9988 )path
501  END IF
502 *
503  ELSE IF( lsamen( 2, c2, 'PO' ) ) THEN
504 *
505 * PO: positive definite matrices
506 *
507  ntypes = 9
508  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
509 *
510  IF( tstchk ) THEN
511  CALL cchkpo( dotype, nn, nval, nnb2, nbval2, nns, nsval,
512  $ thresh, tsterr, lda, a( 1, 1 ), a( 1, 2 ),
513  $ a( 1, 3 ), b( 1, 1 ), b( 1, 2 ), b( 1, 3 ),
514  $ work, rwork, nout )
515  ELSE
516  WRITE( nout, fmt = 9989 )path
517  END IF
518 *
519  IF( tstdrv ) THEN
520  CALL cdrvpo( dotype, nn, nval, nrhs, thresh, tsterr, lda,
521  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
522  $ b( 1, 2 ), b( 1, 3 ), b( 1, 4 ), s, work,
523  $ rwork, nout )
524  ELSE
525  WRITE( nout, fmt = 9988 )path
526  END IF
527 *
528  ELSE IF( lsamen( 2, c2, 'PS' ) ) THEN
529 *
530 * PS: positive semi-definite matrices
531 *
532  ntypes = 9
533 *
534  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
535 *
536  IF( tstchk ) THEN
537  CALL cchkps( dotype, nn, nval, nnb2, nbval2, nrank,
538  $ rankval, thresh, tsterr, lda, a( 1, 1 ),
539  $ a( 1, 2 ), a( 1, 3 ), piv, work, rwork,
540  $ nout )
541  ELSE
542  WRITE( nout, fmt = 9989 )path
543  END IF
544 *
545  ELSE IF( lsamen( 2, c2, 'PP' ) ) THEN
546 *
547 * PP: positive definite packed matrices
548 *
549  ntypes = 9
550  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
551 *
552  IF( tstchk ) THEN
553  CALL cchkpp( dotype, nn, nval, nns, nsval, thresh, tsterr,
554  $ lda, a( 1, 1 ), a( 1, 2 ), a( 1, 3 ),
555  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), work, rwork,
556  $ nout )
557  ELSE
558  WRITE( nout, fmt = 9989 )path
559  END IF
560 *
561  IF( tstdrv ) THEN
562  CALL cdrvpp( dotype, nn, nval, nrhs, thresh, tsterr, lda,
563  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
564  $ b( 1, 2 ), b( 1, 3 ), b( 1, 4 ), s, work,
565  $ rwork, nout )
566  ELSE
567  WRITE( nout, fmt = 9988 )path
568  END IF
569 *
570  ELSE IF( lsamen( 2, c2, 'PB' ) ) THEN
571 *
572 * PB: positive definite banded matrices
573 *
574  ntypes = 8
575  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
576 *
577  IF( tstchk ) THEN
578  CALL cchkpb( dotype, nn, nval, nnb2, nbval2, nns, nsval,
579  $ thresh, tsterr, lda, a( 1, 1 ), a( 1, 2 ),
580  $ a( 1, 3 ), b( 1, 1 ), b( 1, 2 ), b( 1, 3 ),
581  $ work, rwork, nout )
582  ELSE
583  WRITE( nout, fmt = 9989 )path
584  END IF
585 *
586  IF( tstdrv ) THEN
587  CALL cdrvpb( dotype, nn, nval, nrhs, thresh, tsterr, lda,
588  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
589  $ b( 1, 2 ), b( 1, 3 ), b( 1, 4 ), s, work,
590  $ rwork, nout )
591  ELSE
592  WRITE( nout, fmt = 9988 )path
593  END IF
594 *
595  ELSE IF( lsamen( 2, c2, 'PT' ) ) THEN
596 *
597 * PT: positive definite tridiagonal matrices
598 *
599  ntypes = 12
600  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
601 *
602  IF( tstchk ) THEN
603  CALL cchkpt( dotype, nn, nval, nns, nsval, thresh, tsterr,
604  $ a( 1, 1 ), s, a( 1, 2 ), b( 1, 1 ), b( 1, 2 ),
605  $ b( 1, 3 ), work, rwork, nout )
606  ELSE
607  WRITE( nout, fmt = 9989 )path
608  END IF
609 *
610  IF( tstdrv ) THEN
611  CALL cdrvpt( dotype, nn, nval, nrhs, thresh, tsterr,
612  $ a( 1, 1 ), s, a( 1, 2 ), b( 1, 1 ), b( 1, 2 ),
613  $ b( 1, 3 ), work, rwork, nout )
614  ELSE
615  WRITE( nout, fmt = 9988 )path
616  END IF
617 *
618  ELSE IF( lsamen( 2, c2, 'HE' ) ) THEN
619 *
620 * HE: Hermitian indefinite matrices,
621 * with partial (Bunch-Kaufman) pivoting algorithm
622 *
623  ntypes = 10
624  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
625 *
626  IF( tstchk ) THEN
627  CALL cchkhe( dotype, nn, nval, nnb2, nbval2, nns, nsval,
628  $ thresh, tsterr, lda, a( 1, 1 ), a( 1, 2 ),
629  $ a( 1, 3 ), b( 1, 1 ), b( 1, 2 ), b( 1, 3 ),
630  $ work, rwork, iwork, nout )
631  ELSE
632  WRITE( nout, fmt = 9989 )path
633  END IF
634 *
635  IF( tstdrv ) THEN
636  CALL cdrvhe( dotype, nn, nval, nrhs, thresh, tsterr, lda,
637  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
638  $ b( 1, 2 ), b( 1, 3 ), work, rwork, iwork,
639  $ nout )
640  ELSE
641  WRITE( nout, fmt = 9988 )path
642  END IF
643 *
644  ELSE IF( lsamen( 2, c2, 'HR' ) ) THEN
645 *
646 * HR: Hermitian indefinite matrices,
647 * with "rook" (bounded Bunch-Kaufman) pivoting algorithm
648 *
649  ntypes = 10
650  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
651 *
652  IF( tstchk ) THEN
653  CALL cchkhe_rook(dotype, nn, nval, nnb2, nbval2, nns, nsval,
654  $ thresh, tsterr, lda, a( 1, 1 ), a( 1, 2 ),
655  $ a( 1, 3 ), b( 1, 1 ), b( 1, 2 ), b( 1, 3 ),
656  $ work, rwork, iwork, nout )
657  ELSE
658  WRITE( nout, fmt = 9989 )path
659  END IF
660 *
661  IF( tstdrv ) THEN
662  CALL cdrvhe_rook( dotype, nn, nval, nrhs, thresh, tsterr,
663  $ lda, a( 1, 1 ), a( 1, 2 ), a( 1, 3 ),
664  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), work,
665  $ rwork, iwork, nout )
666  ELSE
667  WRITE( nout, fmt = 9988 )path
668  END IF
669 *
670  ELSE IF( lsamen( 2, c2, 'HP' ) ) THEN
671 *
672 * HP: Hermitian indefinite packed matrices,
673 * with partial (Bunch-Kaufman) pivoting algorithm
674 *
675  ntypes = 10
676  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
677 *
678  IF( tstchk ) THEN
679  CALL cchkhp( dotype, nn, nval, nns, nsval, thresh, tsterr,
680  $ lda, a( 1, 1 ), a( 1, 2 ), a( 1, 3 ),
681  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), work, rwork,
682  $ iwork, nout )
683  ELSE
684  WRITE( nout, fmt = 9989 )path
685  END IF
686 *
687  IF( tstdrv ) THEN
688  CALL cdrvhp( dotype, nn, nval, nrhs, thresh, tsterr, lda,
689  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
690  $ b( 1, 2 ), b( 1, 3 ), work, rwork, iwork,
691  $ nout )
692  ELSE
693  WRITE( nout, fmt = 9988 )path
694  END IF
695 *
696  ELSE IF( lsamen( 2, c2, 'SY' ) ) THEN
697 *
698 * SY: symmetric indefinite matrices,
699 * with partial (Bunch-Kaufman) pivoting algorithm
700 *
701  ntypes = 11
702  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
703 *
704  IF( tstchk ) THEN
705  CALL cchksy( dotype, nn, nval, nnb2, nbval2, nns, nsval,
706  $ thresh, tsterr, lda, a( 1, 1 ), a( 1, 2 ),
707  $ a( 1, 3 ), b( 1, 1 ), b( 1, 2 ), b( 1, 3 ),
708  $ work, rwork, iwork, nout )
709  ELSE
710  WRITE( nout, fmt = 9989 )path
711  END IF
712 *
713  IF( tstdrv ) THEN
714  CALL cdrvsy( dotype, nn, nval, nrhs, thresh, tsterr, lda,
715  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
716  $ b( 1, 2 ), b( 1, 3 ), work, rwork, iwork,
717  $ nout )
718  ELSE
719  WRITE( nout, fmt = 9988 )path
720  END IF
721 *
722  ELSE IF( lsamen( 2, c2, 'SR' ) ) THEN
723 *
724 * SR: symmetric indefinite matrices,
725 * with "rook" (bounded Bunch-Kaufman) pivoting algorithm
726 *
727  ntypes = 11
728  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
729 *
730  IF( tstchk ) THEN
731  CALL cchksy_rook(dotype, nn, nval, nnb2, nbval2, nns, nsval,
732  $ thresh, tsterr, lda, a( 1, 1 ), a( 1, 2 ),
733  $ a( 1, 3 ), b( 1, 1 ), b( 1, 2 ), b( 1, 3 ),
734  $ work, rwork, iwork, nout )
735  ELSE
736  WRITE( nout, fmt = 9989 )path
737  END IF
738 *
739  IF( tstdrv ) THEN
740  CALL cdrvsy_rook( dotype, nn, nval, nrhs, thresh, tsterr,
741  $ lda, a( 1, 1 ), a( 1, 2 ), a( 1, 3 ),
742  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), work,
743  $ rwork, iwork, nout )
744  ELSE
745  WRITE( nout, fmt = 9988 )path
746  END IF
747 *
748  ELSE IF( lsamen( 2, c2, 'SP' ) ) THEN
749 *
750 * SP: symmetric indefinite packed matrices,
751 * with partial (Bunch-Kaufman) pivoting algorithm
752 *
753  ntypes = 11
754  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
755 *
756  IF( tstchk ) THEN
757  CALL cchksp( dotype, nn, nval, nns, nsval, thresh, tsterr,
758  $ lda, a( 1, 1 ), a( 1, 2 ), a( 1, 3 ),
759  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), work, rwork,
760  $ iwork, nout )
761  ELSE
762  WRITE( nout, fmt = 9989 )path
763  END IF
764 *
765  IF( tstdrv ) THEN
766  CALL cdrvsp( dotype, nn, nval, nrhs, thresh, tsterr, lda,
767  $ a( 1, 1 ), a( 1, 2 ), a( 1, 3 ), b( 1, 1 ),
768  $ b( 1, 2 ), b( 1, 3 ), work, rwork, iwork,
769  $ nout )
770  ELSE
771  WRITE( nout, fmt = 9988 )path
772  END IF
773 *
774  ELSE IF( lsamen( 2, c2, 'TR' ) ) THEN
775 *
776 * TR: triangular matrices
777 *
778  ntypes = 18
779  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
780 *
781  IF( tstchk ) THEN
782  CALL cchktr( dotype, nn, nval, nnb2, nbval2, nns, nsval,
783  $ thresh, tsterr, lda, a( 1, 1 ), a( 1, 2 ),
784  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), work, rwork,
785  $ nout )
786  ELSE
787  WRITE( nout, fmt = 9989 )path
788  END IF
789 *
790  ELSE IF( lsamen( 2, c2, 'TP' ) ) THEN
791 *
792 * TP: triangular packed matrices
793 *
794  ntypes = 18
795  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
796 *
797  IF( tstchk ) THEN
798  CALL cchktp( dotype, nn, nval, nns, nsval, thresh, tsterr,
799  $ lda, a( 1, 1 ), a( 1, 2 ), b( 1, 1 ),
800  $ b( 1, 2 ), b( 1, 3 ), work, rwork, nout )
801  ELSE
802  WRITE( nout, fmt = 9989 )path
803  END IF
804 *
805  ELSE IF( lsamen( 2, c2, 'TB' ) ) THEN
806 *
807 * TB: triangular banded matrices
808 *
809  ntypes = 17
810  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
811 *
812  IF( tstchk ) THEN
813  CALL cchktb( dotype, nn, nval, nns, nsval, thresh, tsterr,
814  $ lda, a( 1, 1 ), a( 1, 2 ), b( 1, 1 ),
815  $ b( 1, 2 ), b( 1, 3 ), work, rwork, nout )
816  ELSE
817  WRITE( nout, fmt = 9989 )path
818  END IF
819 *
820  ELSE IF( lsamen( 2, c2, 'QR' ) ) THEN
821 *
822 * QR: QR factorization
823 *
824  ntypes = 8
825  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
826 *
827  IF( tstchk ) THEN
828  CALL cchkqr( dotype, nm, mval, nn, nval, nnb, nbval, nxval,
829  $ nrhs, thresh, tsterr, nmax, a( 1, 1 ),
830  $ a( 1, 2 ), a( 1, 3 ), a( 1, 4 ), a( 1, 5 ),
831  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), b( 1, 4 ),
832  $ work, rwork, iwork, nout )
833  ELSE
834  WRITE( nout, fmt = 9989 )path
835  END IF
836 *
837  ELSE IF( lsamen( 2, c2, 'LQ' ) ) THEN
838 *
839 * LQ: LQ factorization
840 *
841  ntypes = 8
842  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
843 *
844  IF( tstchk ) THEN
845  CALL cchklq( dotype, nm, mval, nn, nval, nnb, nbval, nxval,
846  $ nrhs, thresh, tsterr, nmax, a( 1, 1 ),
847  $ a( 1, 2 ), a( 1, 3 ), a( 1, 4 ), a( 1, 5 ),
848  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), b( 1, 4 ),
849  $ work, rwork, nout )
850  ELSE
851  WRITE( nout, fmt = 9989 )path
852  END IF
853 *
854  ELSE IF( lsamen( 2, c2, 'QL' ) ) THEN
855 *
856 * QL: QL factorization
857 *
858  ntypes = 8
859  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
860 *
861  IF( tstchk ) THEN
862  CALL cchkql( dotype, nm, mval, nn, nval, nnb, nbval, nxval,
863  $ nrhs, thresh, tsterr, nmax, a( 1, 1 ),
864  $ a( 1, 2 ), a( 1, 3 ), a( 1, 4 ), a( 1, 5 ),
865  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), b( 1, 4 ),
866  $ work, rwork, nout )
867  ELSE
868  WRITE( nout, fmt = 9989 )path
869  END IF
870 
871 *
872  ELSE IF( lsamen( 2, c2, 'RQ' ) ) THEN
873 *
874 * RQ: RQ factorization
875 *
876  ntypes = 8
877  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
878 *
879  IF( tstchk ) THEN
880  CALL cchkrq( dotype, nm, mval, nn, nval, nnb, nbval, nxval,
881  $ nrhs, thresh, tsterr, nmax, a( 1, 1 ),
882  $ a( 1, 2 ), a( 1, 3 ), a( 1, 4 ), a( 1, 5 ),
883  $ b( 1, 1 ), b( 1, 2 ), b( 1, 3 ), b( 1, 4 ),
884  $ work, rwork, iwork, nout )
885  ELSE
886  WRITE( nout, fmt = 9989 )path
887  END IF
888 *
889  ELSE IF( lsamen( 2, c2, 'EQ' ) ) THEN
890 *
891 * EQ: Equilibration routines for general and positive definite
892 * matrices (THREQ should be between 2 and 10)
893 *
894  IF( tstchk ) THEN
895  CALL cchkeq( threq, nout )
896  ELSE
897  WRITE( nout, fmt = 9989 )path
898  END IF
899 *
900  ELSE IF( lsamen( 2, c2, 'TZ' ) ) THEN
901 *
902 * TZ: Trapezoidal matrix
903 *
904  ntypes = 3
905  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
906 *
907  IF( tstchk ) THEN
908  CALL cchktz( dotype, nm, mval, nn, nval, thresh, tsterr,
909  $ a( 1, 1 ), a( 1, 2 ), s( 1 ),
910  $ b( 1, 1 ), work, rwork, nout )
911  ELSE
912  WRITE( nout, fmt = 9989 )path
913  END IF
914 *
915  ELSE IF( lsamen( 2, c2, 'QP' ) ) THEN
916 *
917 * QP: QR factorization with pivoting
918 *
919  ntypes = 6
920  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
921 *
922  IF( tstchk ) THEN
923  CALL cchkq3( dotype, nm, mval, nn, nval, nnb, nbval, nxval,
924  $ thresh, a( 1, 1 ), a( 1, 2 ), s( 1 ),
925  $ b( 1, 1 ), work, rwork, iwork, nout )
926  ELSE
927  WRITE( nout, fmt = 9989 )path
928  END IF
929 
930 *
931  ELSE IF( lsamen( 2, c2, 'LS' ) ) THEN
932 *
933 * LS: Least squares drivers
934 *
935  ntypes = 6
936  CALL alareq( path, nmats, dotype, ntypes, nin, nout )
937 *
938  IF( tstdrv ) THEN
939  CALL cdrvls( dotype, nm, mval, nn, nval, nns, nsval, nnb,
940  $ nbval, nxval, thresh, tsterr, a( 1, 1 ),
941  $ a( 1, 2 ), a( 1, 3 ), a( 1, 4 ), a( 1, 5 ),
942  $ s( 1 ), s( nmax+1 ), work, rwork, iwork,
943  $ nout )
944  ELSE
945  WRITE( nout, fmt = 9989 )path
946  END IF
947 *
948  ELSE IF( lsamen( 2, c2, 'QT' ) ) THEN
949 *
950 * QT: QRT routines for general matrices
951 *
952  IF( tstchk ) THEN
953  CALL cchkqrt( thresh, tsterr, nm, mval, nn, nval, nnb,
954  $ nbval, nout )
955  ELSE
956  WRITE( nout, fmt = 9989 )path
957  END IF
958 *
959  ELSE IF( lsamen( 2, c2, 'QX' ) ) THEN
960 *
961 * QX: QRT routines for triangular-pentagonal matrices
962 *
963  IF( tstchk ) THEN
964  CALL cchkqrtp( thresh, tsterr, nm, mval, nn, nval, nnb,
965  $ nbval, nout )
966  ELSE
967  WRITE( nout, fmt = 9989 )path
968  END IF
969 *
970  ELSE
971 *
972  WRITE( nout, fmt = 9990 )path
973  END IF
974 *
975 * Go back to get another input line.
976 *
977  GO TO 80
978 *
979 * Branch to this line when the last record is read.
980 *
981  140 CONTINUE
982  CLOSE ( nin )
983  s2 = second( )
984  WRITE( nout, fmt = 9998 )
985  WRITE( nout, fmt = 9997 )s2 - s1
986 *
987  9999 FORMAT( / ' Execution not attempted due to input errors' )
988  9998 FORMAT( / ' End of tests' )
989  9997 FORMAT( ' Total time used = ', f12.2, ' seconds', / )
990  9996 FORMAT( ' Invalid input value: ', a4, '=', i6, '; must be >=',
991  $ i6 )
992  9995 FORMAT( ' Invalid input value: ', a4, '=', i6, '; must be <=',
993  $ i6 )
994  9994 FORMAT( ' Tests of the COMPLEX LAPACK routines ',
995  $ / ' LAPACK VERSION ', i1, '.', i1, '.', i1,
996  $ / / ' The following parameter values will be used:' )
997  9993 FORMAT( 4x, a4, ': ', 10i6, / 11x, 10i6 )
998  9992 FORMAT( / ' Routines pass computational tests if test ratio is ',
999  $ 'less than', f8.2, / )
1000  9991 FORMAT( ' Relative machine ', a, ' is taken to be', e16.6 )
1001  9990 FORMAT( / 1x, a3, ': Unrecognized path name' )
1002  9989 FORMAT( / 1x, a3, ' routines were not tested' )
1003  9988 FORMAT( / 1x, a3, ' driver routines were not tested' )
1004 *
1005 * End of CCHKAA
1006 *
1007  END
subroutine cdrvgb(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, A, LA, AFB, LAFB, ASAV, B, BSAV, X, XACT, S, WORK, RWORK, IWORK, NOUT)
CDRVGB
Definition: cdrvgb.f:174
subroutine cdrvhp(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CDRVHP
Definition: cdrvhp.f:159
subroutine cchksy(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKSY
Definition: cchksy.f:173
subroutine cchkqrt(THRESH, TSTERR, NM, MVAL, NN, NVAL, NNB, NBVAL, NOUT)
CCHKQRT
Definition: cchkqrt.f:104
subroutine cdrvpb(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK, RWORK, NOUT)
CDRVPB
Definition: cdrvpb.f:161
subroutine cdrvsy(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CDRVSY
Definition: cdrvsy.f:155
subroutine cchkqrtp(THRESH, TSTERR, NM, MVAL, NN, NVAL, NNB, NBVAL, NOUT)
CCHKQRTP
Definition: cchkqrtp.f:104
subroutine cdrvhe(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CDRVHE
Definition: cdrvhe.f:155
subroutine cchkpp(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, NOUT)
CCHKPP
Definition: cchkpp.f:161
subroutine cchkeq(THRESH, NOUT)
CCHKEQ
Definition: cchkeq.f:56
subroutine cchkql(DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NXVAL, NRHS, THRESH, TSTERR, NMAX, A, AF, AQ, AL, AC, B, X, XACT, TAU, WORK, RWORK, NOUT)
CCHKQL
Definition: cchkql.f:198
subroutine cchkhe_rook(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKHE_ROOK
Definition: cchkhe_rook.f:174
subroutine alareq(PATH, NMATS, DOTYPE, NTYPES, NIN, NOUT)
ALAREQ
Definition: alareq.f:92
subroutine cchktr(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AINV, B, X, XACT, WORK, RWORK, NOUT)
CCHKTR
Definition: cchktr.f:165
subroutine cdrvhe_rook(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CDRVHE_ROOK
Definition: cdrvhe_rook.f:155
subroutine cdrvls(DOTYPE, NM, MVAL, NN, NVAL, NNS, NSVAL, NNB, NBVAL, NXVAL, THRESH, TSTERR, A, COPYA, B, COPYB, C, S, COPYS, WORK, RWORK, IWORK, NOUT)
CDRVLS
Definition: cdrvls.f:213
program cchkaa
CCHKAA
Definition: cchkaa.f:110
subroutine cdrvsy_rook(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CDRVSY_ROOK
Definition: cdrvsy_rook.f:154
subroutine cchktp(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, AP, AINVP, B, X, XACT, WORK, RWORK, NOUT)
CCHKTP
Definition: cchktp.f:153
subroutine cdrvge(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK, RWORK, IWORK, NOUT)
CDRVGE
Definition: cdrvge.f:166
subroutine cchklq(DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NXVAL, NRHS, THRESH, TSTERR, NMAX, A, AF, AQ, AL, AC, B, X, XACT, TAU, WORK, RWORK, NOUT)
CCHKLQ
Definition: cchklq.f:198
subroutine cchkgb(DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, A, LA, AFAC, LAFAC, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKGB
Definition: cchkgb.f:193
subroutine cdrvsp(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CDRVSP
Definition: cdrvsp.f:159
subroutine cchkgt(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, A, AF, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKGT
Definition: cchkgt.f:149
subroutine cchkrq(DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NXVAL, NRHS, THRESH, TSTERR, NMAX, A, AF, AQ, AR, AC, B, X, XACT, TAU, WORK, RWORK, IWORK, NOUT)
CCHKRQ
Definition: cchkrq.f:203
subroutine cchkge(DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKGE
Definition: cchkge.f:188
subroutine cdrvgt(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, A, AF, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CDRVGT
Definition: cdrvgt.f:141
subroutine cchktb(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, AB, AINV, B, X, XACT, WORK, RWORK, NOUT)
CCHKTB
Definition: cchktb.f:151
subroutine cdrvpt(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, A, D, E, B, X, XACT, WORK, RWORK, NOUT)
CDRVPT
Definition: cdrvpt.f:142
subroutine cchkpt(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, A, D, E, B, X, XACT, WORK, RWORK, NOUT)
CCHKPT
Definition: cchkpt.f:149
subroutine ilaver(VERS_MAJOR, VERS_MINOR, VERS_PATCH)
ILAVER returns the LAPACK version.
Definition: ilaver.f:50
subroutine cchktz(DOTYPE, NM, MVAL, NN, NVAL, THRESH, TSTERR, A, COPYA, S, TAU, WORK, RWORK, NOUT)
CCHKTZ
Definition: cchktz.f:139
subroutine cchkhe(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKHE
Definition: cchkhe.f:173
subroutine cchkqr(DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NXVAL, NRHS, THRESH, TSTERR, NMAX, A, AF, AQ, AR, AC, B, X, XACT, TAU, WORK, RWORK, IWORK, NOUT)
CCHKQR
Definition: cchkqr.f:203
subroutine cchkps(DOTYPE, NN, NVAL, NNB, NBVAL, NRANK, RANKVAL, THRESH, TSTERR, NMAX, A, AFAC, PERM, PIV, WORK, RWORK, NOUT)
CCHKPS
Definition: cchkps.f:156
subroutine cdrvpo(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK, RWORK, NOUT)
CDRVPO
Definition: cdrvpo.f:161
subroutine cdrvpp(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, A, AFAC, ASAV, B, BSAV, X, XACT, S, WORK, RWORK, NOUT)
CDRVPP
Definition: cdrvpp.f:161
subroutine cchkpb(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, NOUT)
CCHKPB
Definition: cchkpb.f:170
subroutine cchkq3(DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NXVAL, THRESH, A, COPYA, S, TAU, WORK, RWORK, IWORK, NOUT)
CCHKQ3
Definition: cchkq3.f:160
subroutine cchksy_rook(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKSY_ROOK
Definition: cchksy_rook.f:174
subroutine cchksp(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKSP
Definition: cchksp.f:166
subroutine cchkpo(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, NOUT)
CCHKPO
Definition: cchkpo.f:170
subroutine cchkhp(DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
CCHKHP
Definition: cchkhp.f:166