LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
cblat3.f
Go to the documentation of this file.
1 *> \brief \b CBLAT3
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 CBLAT3
12 *
13 *
14 *> \par Purpose:
15 * =============
16 *>
17 *> \verbatim
18 *>
19 *> Test program for the COMPLEX Level 3 Blas.
20 *>
21 *> The program must be driven by a short data file. The first 14 records
22 *> of the file are read using list-directed input, the last 9 records
23 *> are read using the format ( A6, L2 ). An annotated example of a data
24 *> file can be obtained by deleting the first 3 characters from the
25 *> following 23 lines:
26 *> 'cblat3.out' NAME OF SUMMARY OUTPUT FILE
27 *> 6 UNIT NUMBER OF SUMMARY FILE
28 *> 'CBLAT3.SNAP' NAME OF SNAPSHOT OUTPUT FILE
29 *> -1 UNIT NUMBER OF SNAPSHOT FILE (NOT USED IF .LT. 0)
30 *> F LOGICAL FLAG, T TO REWIND SNAPSHOT FILE AFTER EACH RECORD.
31 *> F LOGICAL FLAG, T TO STOP ON FAILURES.
32 *> T LOGICAL FLAG, T TO TEST ERROR EXITS.
33 *> 16.0 THRESHOLD VALUE OF TEST RATIO
34 *> 6 NUMBER OF VALUES OF N
35 *> 0 1 2 3 5 9 VALUES OF N
36 *> 3 NUMBER OF VALUES OF ALPHA
37 *> (0.0,0.0) (1.0,0.0) (0.7,-0.9) VALUES OF ALPHA
38 *> 3 NUMBER OF VALUES OF BETA
39 *> (0.0,0.0) (1.0,0.0) (1.3,-1.1) VALUES OF BETA
40 *> CGEMM T PUT F FOR NO TEST. SAME COLUMNS.
41 *> CHEMM T PUT F FOR NO TEST. SAME COLUMNS.
42 *> CSYMM T PUT F FOR NO TEST. SAME COLUMNS.
43 *> CTRMM T PUT F FOR NO TEST. SAME COLUMNS.
44 *> CTRSM T PUT F FOR NO TEST. SAME COLUMNS.
45 *> CHERK T PUT F FOR NO TEST. SAME COLUMNS.
46 *> CSYRK T PUT F FOR NO TEST. SAME COLUMNS.
47 *> CHER2K T PUT F FOR NO TEST. SAME COLUMNS.
48 *> CSYR2K T PUT F FOR NO TEST. SAME COLUMNS.
49 *>
50 *> Further Details
51 *> ===============
52 *>
53 *> See:
54 *>
55 *> Dongarra J. J., Du Croz J. J., Duff I. S. and Hammarling S.
56 *> A Set of Level 3 Basic Linear Algebra Subprograms.
57 *>
58 *> Technical Memorandum No.88 (Revision 1), Mathematics and
59 *> Computer Science Division, Argonne National Laboratory, 9700
60 *> South Cass Avenue, Argonne, Illinois 60439, US.
61 *>
62 *> -- Written on 8-February-1989.
63 *> Jack Dongarra, Argonne National Laboratory.
64 *> Iain Duff, AERE Harwell.
65 *> Jeremy Du Croz, Numerical Algorithms Group Ltd.
66 *> Sven Hammarling, Numerical Algorithms Group Ltd.
67 *>
68 *> 10-9-00: Change STATUS='NEW' to 'UNKNOWN' so that the testers
69 *> can be run multiple times without deleting generated
70 *> output files (susan)
71 *> \endverbatim
72 *
73 * Authors:
74 * ========
75 *
76 *> \author Univ. of Tennessee
77 *> \author Univ. of California Berkeley
78 *> \author Univ. of Colorado Denver
79 *> \author NAG Ltd.
80 *
81 *> \date April 2012
82 *
83 *> \ingroup complex_blas_testing
84 *
85 * =====================================================================
86  PROGRAM cblat3
87 *
88 * -- Reference BLAS test routine (version 3.4.1) --
89 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
90 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
91 * April 2012
92 *
93 * =====================================================================
94 *
95 * .. Parameters ..
96  INTEGER nin
97  parameter( nin = 5 )
98  INTEGER nsubs
99  parameter( nsubs = 9 )
100  COMPLEX zero, one
101  parameter( zero = ( 0.0, 0.0 ), one = ( 1.0, 0.0 ) )
102  REAL rzero
103  parameter( rzero = 0.0 )
104  INTEGER nmax
105  parameter( nmax = 65 )
106  INTEGER nidmax, nalmax, nbemax
107  parameter( nidmax = 9, nalmax = 7, nbemax = 7 )
108 * .. Local Scalars ..
109  REAL eps, err, thresh
110  INTEGER i, isnum, j, n, nalf, nbet, nidim, nout, ntra
111  LOGICAL fatal, ltestt, rewi, same, sfatal, trace,
112  $ tsterr
113  CHARACTER*1 transa, transb
114  CHARACTER*6 snamet
115  CHARACTER*32 snaps, summry
116 * .. Local Arrays ..
117  COMPLEX aa( nmax*nmax ), ab( nmax, 2*nmax ),
118  $ alf( nalmax ), as( nmax*nmax ),
119  $ bb( nmax*nmax ), bet( nbemax ),
120  $ bs( nmax*nmax ), c( nmax, nmax ),
121  $ cc( nmax*nmax ), cs( nmax*nmax ), ct( nmax ),
122  $ w( 2*nmax )
123  REAL g( nmax )
124  INTEGER idim( nidmax )
125  LOGICAL ltest( nsubs )
126  CHARACTER*6 snames( nsubs )
127 * .. External Functions ..
128  REAL sdiff
129  LOGICAL lce
130  EXTERNAL sdiff, lce
131 * .. External Subroutines ..
132  EXTERNAL cchk1, cchk2, cchk3, cchk4, cchk5, cchke, cmmch
133 * .. Intrinsic Functions ..
134  INTRINSIC max, min
135 * .. Scalars in Common ..
136  INTEGER infot, noutc
137  LOGICAL lerr, ok
138  CHARACTER*6 srnamt
139 * .. Common blocks ..
140  common /infoc/infot, noutc, ok, lerr
141  common /srnamc/srnamt
142 * .. Data statements ..
143  DATA snames/'CGEMM ', 'CHEMM ', 'CSYMM ', 'CTRMM ',
144  $ 'CTRSM ', 'CHERK ', 'CSYRK ', 'CHER2K',
145  $ 'CSYR2K'/
146 * .. Executable Statements ..
147 *
148 * Read name and unit number for summary output file and open file.
149 *
150  READ( nin, fmt = * )summry
151  READ( nin, fmt = * )nout
152  OPEN( nout, file = summry )
153  noutc = nout
154 *
155 * Read name and unit number for snapshot output file and open file.
156 *
157  READ( nin, fmt = * )snaps
158  READ( nin, fmt = * )ntra
159  trace = ntra.GE.0
160  IF( trace )THEN
161  OPEN( ntra, file = snaps )
162  END IF
163 * Read the flag that directs rewinding of the snapshot file.
164  READ( nin, fmt = * )rewi
165  rewi = rewi.AND.trace
166 * Read the flag that directs stopping on any failure.
167  READ( nin, fmt = * )sfatal
168 * Read the flag that indicates whether error exits are to be tested.
169  READ( nin, fmt = * )tsterr
170 * Read the threshold value of the test ratio
171  READ( nin, fmt = * )thresh
172 *
173 * Read and check the parameter values for the tests.
174 *
175 * Values of N
176  READ( nin, fmt = * )nidim
177  IF( nidim.LT.1.OR.nidim.GT.nidmax )THEN
178  WRITE( nout, fmt = 9997 )'N', nidmax
179  go to 220
180  END IF
181  READ( nin, fmt = * )( idim( i ), i = 1, nidim )
182  DO 10 i = 1, nidim
183  IF( idim( i ).LT.0.OR.idim( i ).GT.nmax )THEN
184  WRITE( nout, fmt = 9996 )nmax
185  go to 220
186  END IF
187  10 continue
188 * Values of ALPHA
189  READ( nin, fmt = * )nalf
190  IF( nalf.LT.1.OR.nalf.GT.nalmax )THEN
191  WRITE( nout, fmt = 9997 )'ALPHA', nalmax
192  go to 220
193  END IF
194  READ( nin, fmt = * )( alf( i ), i = 1, nalf )
195 * Values of BETA
196  READ( nin, fmt = * )nbet
197  IF( nbet.LT.1.OR.nbet.GT.nbemax )THEN
198  WRITE( nout, fmt = 9997 )'BETA', nbemax
199  go to 220
200  END IF
201  READ( nin, fmt = * )( bet( i ), i = 1, nbet )
202 *
203 * Report values of parameters.
204 *
205  WRITE( nout, fmt = 9995 )
206  WRITE( nout, fmt = 9994 )( idim( i ), i = 1, nidim )
207  WRITE( nout, fmt = 9993 )( alf( i ), i = 1, nalf )
208  WRITE( nout, fmt = 9992 )( bet( i ), i = 1, nbet )
209  IF( .NOT.tsterr )THEN
210  WRITE( nout, fmt = * )
211  WRITE( nout, fmt = 9984 )
212  END IF
213  WRITE( nout, fmt = * )
214  WRITE( nout, fmt = 9999 )thresh
215  WRITE( nout, fmt = * )
216 *
217 * Read names of subroutines and flags which indicate
218 * whether they are to be tested.
219 *
220  DO 20 i = 1, nsubs
221  ltest( i ) = .false.
222  20 continue
223  30 READ( nin, fmt = 9988, END = 60 )snamet, ltestt
224  DO 40 i = 1, nsubs
225  IF( snamet.EQ.snames( i ) )
226  $ go to 50
227  40 continue
228  WRITE( nout, fmt = 9990 )snamet
229  stop
230  50 ltest( i ) = ltestt
231  go to 30
232 *
233  60 continue
234  CLOSE ( nin )
235 *
236 * Compute EPS (the machine precision).
237 *
238  eps = epsilon(rzero)
239  WRITE( nout, fmt = 9998 )eps
240 *
241 * Check the reliability of CMMCH using exact data.
242 *
243  n = min( 32, nmax )
244  DO 100 j = 1, n
245  DO 90 i = 1, n
246  ab( i, j ) = max( i - j + 1, 0 )
247  90 continue
248  ab( j, nmax + 1 ) = j
249  ab( 1, nmax + j ) = j
250  c( j, 1 ) = zero
251  100 continue
252  DO 110 j = 1, n
253  cc( j ) = j*( ( j + 1 )*j )/2 - ( ( j + 1 )*j*( j - 1 ) )/3
254  110 continue
255 * CC holds the exact result. On exit from CMMCH CT holds
256 * the result computed by CMMCH.
257  transa = 'N'
258  transb = 'N'
259  CALL cmmch( transa, transb, n, 1, n, one, ab, nmax,
260  $ ab( 1, nmax + 1 ), nmax, zero, c, nmax, ct, g, cc,
261  $ nmax, eps, err, fatal, nout, .true. )
262  same = lce( cc, ct, n )
263  IF( .NOT.same.OR.err.NE.rzero )THEN
264  WRITE( nout, fmt = 9989 )transa, transb, same, err
265  stop
266  END IF
267  transb = 'C'
268  CALL cmmch( transa, transb, n, 1, n, one, ab, nmax,
269  $ ab( 1, nmax + 1 ), nmax, zero, c, nmax, ct, g, cc,
270  $ nmax, eps, err, fatal, nout, .true. )
271  same = lce( cc, ct, n )
272  IF( .NOT.same.OR.err.NE.rzero )THEN
273  WRITE( nout, fmt = 9989 )transa, transb, same, err
274  stop
275  END IF
276  DO 120 j = 1, n
277  ab( j, nmax + 1 ) = n - j + 1
278  ab( 1, nmax + j ) = n - j + 1
279  120 continue
280  DO 130 j = 1, n
281  cc( n - j + 1 ) = j*( ( j + 1 )*j )/2 -
282  $ ( ( j + 1 )*j*( j - 1 ) )/3
283  130 continue
284  transa = 'C'
285  transb = 'N'
286  CALL cmmch( transa, transb, n, 1, n, one, ab, nmax,
287  $ ab( 1, nmax + 1 ), nmax, zero, c, nmax, ct, g, cc,
288  $ nmax, eps, err, fatal, nout, .true. )
289  same = lce( cc, ct, n )
290  IF( .NOT.same.OR.err.NE.rzero )THEN
291  WRITE( nout, fmt = 9989 )transa, transb, same, err
292  stop
293  END IF
294  transb = 'C'
295  CALL cmmch( transa, transb, n, 1, n, one, ab, nmax,
296  $ ab( 1, nmax + 1 ), nmax, zero, c, nmax, ct, g, cc,
297  $ nmax, eps, err, fatal, nout, .true. )
298  same = lce( cc, ct, n )
299  IF( .NOT.same.OR.err.NE.rzero )THEN
300  WRITE( nout, fmt = 9989 )transa, transb, same, err
301  stop
302  END IF
303 *
304 * Test each subroutine in turn.
305 *
306  DO 200 isnum = 1, nsubs
307  WRITE( nout, fmt = * )
308  IF( .NOT.ltest( isnum ) )THEN
309 * Subprogram is not to be tested.
310  WRITE( nout, fmt = 9987 )snames( isnum )
311  ELSE
312  srnamt = snames( isnum )
313 * Test error exits.
314  IF( tsterr )THEN
315  CALL cchke( isnum, snames( isnum ), nout )
316  WRITE( nout, fmt = * )
317  END IF
318 * Test computations.
319  infot = 0
320  ok = .true.
321  fatal = .false.
322  go to( 140, 150, 150, 160, 160, 170, 170,
323  $ 180, 180 )isnum
324 * Test CGEMM, 01.
325  140 CALL cchk1( snames( isnum ), eps, thresh, nout, ntra, trace,
326  $ rewi, fatal, nidim, idim, nalf, alf, nbet, bet,
327  $ nmax, ab, aa, as, ab( 1, nmax + 1 ), bb, bs, c,
328  $ cc, cs, ct, g )
329  go to 190
330 * Test CHEMM, 02, CSYMM, 03.
331  150 CALL cchk2( snames( isnum ), eps, thresh, nout, ntra, trace,
332  $ rewi, fatal, nidim, idim, nalf, alf, nbet, bet,
333  $ nmax, ab, aa, as, ab( 1, nmax + 1 ), bb, bs, c,
334  $ cc, cs, ct, g )
335  go to 190
336 * Test CTRMM, 04, CTRSM, 05.
337  160 CALL cchk3( snames( isnum ), eps, thresh, nout, ntra, trace,
338  $ rewi, fatal, nidim, idim, nalf, alf, nmax, ab,
339  $ aa, as, ab( 1, nmax + 1 ), bb, bs, ct, g, c )
340  go to 190
341 * Test CHERK, 06, CSYRK, 07.
342  170 CALL cchk4( snames( isnum ), eps, thresh, nout, ntra, trace,
343  $ rewi, fatal, nidim, idim, nalf, alf, nbet, bet,
344  $ nmax, ab, aa, as, ab( 1, nmax + 1 ), bb, bs, c,
345  $ cc, cs, ct, g )
346  go to 190
347 * Test CHER2K, 08, CSYR2K, 09.
348  180 CALL cchk5( snames( isnum ), eps, thresh, nout, ntra, trace,
349  $ rewi, fatal, nidim, idim, nalf, alf, nbet, bet,
350  $ nmax, ab, aa, as, bb, bs, c, cc, cs, ct, g, w )
351  go to 190
352 *
353  190 IF( fatal.AND.sfatal )
354  $ go to 210
355  END IF
356  200 continue
357  WRITE( nout, fmt = 9986 )
358  go to 230
359 *
360  210 continue
361  WRITE( nout, fmt = 9985 )
362  go to 230
363 *
364  220 continue
365  WRITE( nout, fmt = 9991 )
366 *
367  230 continue
368  IF( trace )
369  $ CLOSE ( ntra )
370  CLOSE ( nout )
371  stop
372 *
373  9999 format( ' ROUTINES PASS COMPUTATIONAL TESTS IF TEST RATIO IS LES',
374  $ 'S THAN', f8.2 )
375  9998 format( ' RELATIVE MACHINE PRECISION IS TAKEN TO BE', 1p, e9.1 )
376  9997 format( ' NUMBER OF VALUES OF ', a, ' IS LESS THAN 1 OR GREATER ',
377  $ 'THAN ', i2 )
378  9996 format( ' VALUE OF N IS LESS THAN 0 OR GREATER THAN ', i2 )
379  9995 format( ' TESTS OF THE COMPLEX LEVEL 3 BLAS', //' THE F',
380  $ 'OLLOWING PARAMETER VALUES WILL BE USED:' )
381  9994 format( ' FOR N ', 9i6 )
382  9993 format( ' FOR ALPHA ',
383  $ 7( '(', f4.1, ',', f4.1, ') ', : ) )
384  9992 format( ' FOR BETA ',
385  $ 7( '(', f4.1, ',', f4.1, ') ', : ) )
386  9991 format( ' AMEND DATA FILE OR INCREASE ARRAY SIZES IN PROGRAM',
387  $ /' ******* TESTS ABANDONED *******' )
388  9990 format( ' SUBPROGRAM NAME ', a6, ' NOT RECOGNIZED', /' ******* T',
389  $ 'ESTS ABANDONED *******' )
390  9989 format( ' ERROR IN CMMCH - IN-LINE DOT PRODUCTS ARE BEING EVALU',
391  $ 'ATED WRONGLY.', /' CMMCH WAS CALLED WITH TRANSA = ', a1,
392  $ ' AND TRANSB = ', a1, /' AND RETURNED SAME = ', l1, ' AND ',
393  $ 'ERR = ', f12.3, '.', /' THIS MAY BE DUE TO FAULTS IN THE ',
394  $ 'ARITHMETIC OR THE COMPILER.', /' ******* TESTS ABANDONED ',
395  $ '*******' )
396  9988 format( a6, l2 )
397  9987 format( 1x, a6, ' WAS NOT TESTED' )
398  9986 format( /' END OF TESTS' )
399  9985 format( /' ******* FATAL ERROR - TESTS ABANDONED *******' )
400  9984 format( ' ERROR-EXITS WILL NOT BE TESTED' )
401 *
402 * End of CBLAT3.
403 *
404  END
405  SUBROUTINE cchk1( SNAME, EPS, THRESH, NOUT, NTRA, TRACE, REWI,
406  $ fatal, nidim, idim, nalf, alf, nbet, bet, nmax,
407  $ a, aa, as, b, bb, bs, c, cc, cs, ct, g )
408 *
409 * Tests CGEMM.
410 *
411 * Auxiliary routine for test program for Level 3 Blas.
412 *
413 * -- Written on 8-February-1989.
414 * Jack Dongarra, Argonne National Laboratory.
415 * Iain Duff, AERE Harwell.
416 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
417 * Sven Hammarling, Numerical Algorithms Group Ltd.
418 *
419 * .. Parameters ..
420  COMPLEX zero
421  parameter( zero = ( 0.0, 0.0 ) )
422  REAL rzero
423  parameter( rzero = 0.0 )
424 * .. Scalar Arguments ..
425  REAL eps, thresh
426  INTEGER nalf, nbet, nidim, nmax, nout, ntra
427  LOGICAL fatal, rewi, trace
428  CHARACTER*6 sname
429 * .. Array Arguments ..
430  COMPLEX a( nmax, nmax ), aa( nmax*nmax ), alf( nalf ),
431  $ as( nmax*nmax ), b( nmax, nmax ),
432  $ bb( nmax*nmax ), bet( nbet ), bs( nmax*nmax ),
433  $ c( nmax, nmax ), cc( nmax*nmax ),
434  $ cs( nmax*nmax ), ct( nmax )
435  REAL g( nmax )
436  INTEGER idim( nidim )
437 * .. Local Scalars ..
438  COMPLEX alpha, als, beta, bls
439  REAL err, errmax
440  INTEGER i, ia, ib, ica, icb, ik, im, in, k, ks, laa,
441  $ lbb, lcc, lda, ldas, ldb, ldbs, ldc, ldcs, m,
442  $ ma, mb, ms, n, na, nargs, nb, nc, ns
443  LOGICAL null, reset, same, trana, tranb
444  CHARACTER*1 tranas, tranbs, transa, transb
445  CHARACTER*3 ich
446 * .. Local Arrays ..
447  LOGICAL isame( 13 )
448 * .. External Functions ..
449  LOGICAL lce, lceres
450  EXTERNAL lce, lceres
451 * .. External Subroutines ..
452  EXTERNAL cgemm, cmake, cmmch
453 * .. Intrinsic Functions ..
454  INTRINSIC max
455 * .. Scalars in Common ..
456  INTEGER infot, noutc
457  LOGICAL lerr, ok
458 * .. Common blocks ..
459  common /infoc/infot, noutc, ok, lerr
460 * .. Data statements ..
461  DATA ich/'NTC'/
462 * .. Executable Statements ..
463 *
464  nargs = 13
465  nc = 0
466  reset = .true.
467  errmax = rzero
468 *
469  DO 110 im = 1, nidim
470  m = idim( im )
471 *
472  DO 100 in = 1, nidim
473  n = idim( in )
474 * Set LDC to 1 more than minimum value if room.
475  ldc = m
476  IF( ldc.LT.nmax )
477  $ ldc = ldc + 1
478 * Skip tests if not enough room.
479  IF( ldc.GT.nmax )
480  $ go to 100
481  lcc = ldc*n
482  null = n.LE.0.OR.m.LE.0
483 *
484  DO 90 ik = 1, nidim
485  k = idim( ik )
486 *
487  DO 80 ica = 1, 3
488  transa = ich( ica: ica )
489  trana = transa.EQ.'T'.OR.transa.EQ.'C'
490 *
491  IF( trana )THEN
492  ma = k
493  na = m
494  ELSE
495  ma = m
496  na = k
497  END IF
498 * Set LDA to 1 more than minimum value if room.
499  lda = ma
500  IF( lda.LT.nmax )
501  $ lda = lda + 1
502 * Skip tests if not enough room.
503  IF( lda.GT.nmax )
504  $ go to 80
505  laa = lda*na
506 *
507 * Generate the matrix A.
508 *
509  CALL cmake( 'GE', ' ', ' ', ma, na, a, nmax, aa, lda,
510  $ reset, zero )
511 *
512  DO 70 icb = 1, 3
513  transb = ich( icb: icb )
514  tranb = transb.EQ.'T'.OR.transb.EQ.'C'
515 *
516  IF( tranb )THEN
517  mb = n
518  nb = k
519  ELSE
520  mb = k
521  nb = n
522  END IF
523 * Set LDB to 1 more than minimum value if room.
524  ldb = mb
525  IF( ldb.LT.nmax )
526  $ ldb = ldb + 1
527 * Skip tests if not enough room.
528  IF( ldb.GT.nmax )
529  $ go to 70
530  lbb = ldb*nb
531 *
532 * Generate the matrix B.
533 *
534  CALL cmake( 'GE', ' ', ' ', mb, nb, b, nmax, bb,
535  $ ldb, reset, zero )
536 *
537  DO 60 ia = 1, nalf
538  alpha = alf( ia )
539 *
540  DO 50 ib = 1, nbet
541  beta = bet( ib )
542 *
543 * Generate the matrix C.
544 *
545  CALL cmake( 'GE', ' ', ' ', m, n, c, nmax,
546  $ cc, ldc, reset, zero )
547 *
548  nc = nc + 1
549 *
550 * Save every datum before calling the
551 * subroutine.
552 *
553  tranas = transa
554  tranbs = transb
555  ms = m
556  ns = n
557  ks = k
558  als = alpha
559  DO 10 i = 1, laa
560  as( i ) = aa( i )
561  10 continue
562  ldas = lda
563  DO 20 i = 1, lbb
564  bs( i ) = bb( i )
565  20 continue
566  ldbs = ldb
567  bls = beta
568  DO 30 i = 1, lcc
569  cs( i ) = cc( i )
570  30 continue
571  ldcs = ldc
572 *
573 * Call the subroutine.
574 *
575  IF( trace )
576  $ WRITE( ntra, fmt = 9995 )nc, sname,
577  $ transa, transb, m, n, k, alpha, lda, ldb,
578  $ beta, ldc
579  IF( rewi )
580  $ rewind ntra
581  CALL cgemm( transa, transb, m, n, k, alpha,
582  $ aa, lda, bb, ldb, beta, cc, ldc )
583 *
584 * Check if error-exit was taken incorrectly.
585 *
586  IF( .NOT.ok )THEN
587  WRITE( nout, fmt = 9994 )
588  fatal = .true.
589  go to 120
590  END IF
591 *
592 * See what data changed inside subroutines.
593 *
594  isame( 1 ) = transa.EQ.tranas
595  isame( 2 ) = transb.EQ.tranbs
596  isame( 3 ) = ms.EQ.m
597  isame( 4 ) = ns.EQ.n
598  isame( 5 ) = ks.EQ.k
599  isame( 6 ) = als.EQ.alpha
600  isame( 7 ) = lce( as, aa, laa )
601  isame( 8 ) = ldas.EQ.lda
602  isame( 9 ) = lce( bs, bb, lbb )
603  isame( 10 ) = ldbs.EQ.ldb
604  isame( 11 ) = bls.EQ.beta
605  IF( null )THEN
606  isame( 12 ) = lce( cs, cc, lcc )
607  ELSE
608  isame( 12 ) = lceres( 'GE', ' ', m, n, cs,
609  $ cc, ldc )
610  END IF
611  isame( 13 ) = ldcs.EQ.ldc
612 *
613 * If data was incorrectly changed, report
614 * and return.
615 *
616  same = .true.
617  DO 40 i = 1, nargs
618  same = same.AND.isame( i )
619  IF( .NOT.isame( i ) )
620  $ WRITE( nout, fmt = 9998 )i
621  40 continue
622  IF( .NOT.same )THEN
623  fatal = .true.
624  go to 120
625  END IF
626 *
627  IF( .NOT.null )THEN
628 *
629 * Check the result.
630 *
631  CALL cmmch( transa, transb, m, n, k,
632  $ alpha, a, nmax, b, nmax, beta,
633  $ c, nmax, ct, g, cc, ldc, eps,
634  $ err, fatal, nout, .true. )
635  errmax = max( errmax, err )
636 * If got really bad answer, report and
637 * return.
638  IF( fatal )
639  $ go to 120
640  END IF
641 *
642  50 continue
643 *
644  60 continue
645 *
646  70 continue
647 *
648  80 continue
649 *
650  90 continue
651 *
652  100 continue
653 *
654  110 continue
655 *
656 * Report result.
657 *
658  IF( errmax.LT.thresh )THEN
659  WRITE( nout, fmt = 9999 )sname, nc
660  ELSE
661  WRITE( nout, fmt = 9997 )sname, nc, errmax
662  END IF
663  go to 130
664 *
665  120 continue
666  WRITE( nout, fmt = 9996 )sname
667  WRITE( nout, fmt = 9995 )nc, sname, transa, transb, m, n, k,
668  $ alpha, lda, ldb, beta, ldc
669 *
670  130 continue
671  return
672 *
673  9999 format( ' ', a6, ' PASSED THE COMPUTATIONAL TESTS (', i6, ' CALL',
674  $ 'S)' )
675  9998 format( ' ******* FATAL ERROR - PARAMETER NUMBER ', i2, ' WAS CH',
676  $ 'ANGED INCORRECTLY *******' )
677  9997 format( ' ', a6, ' COMPLETED THE COMPUTATIONAL TESTS (', i6, ' C',
678  $ 'ALLS)', /' ******* BUT WITH MAXIMUM TEST RATIO', f8.2,
679  $ ' - SUSPECT *******' )
680  9996 format( ' ******* ', a6, ' FAILED ON CALL NUMBER:' )
681  9995 format( 1x, i6, ': ', a6, '(''', a1, ''',''', a1, ''',',
682  $ 3( i3, ',' ), '(', f4.1, ',', f4.1, '), A,', i3, ', B,', i3,
683  $ ',(', f4.1, ',', f4.1, '), C,', i3, ').' )
684  9994 format( ' ******* FATAL ERROR - ERROR-EXIT TAKEN ON VALID CALL *',
685  $ '******' )
686 *
687 * End of CCHK1.
688 *
689  END
690  SUBROUTINE cchk2( SNAME, EPS, THRESH, NOUT, NTRA, TRACE, REWI,
691  $ fatal, nidim, idim, nalf, alf, nbet, bet, nmax,
692  $ a, aa, as, b, bb, bs, c, cc, cs, ct, g )
693 *
694 * Tests CHEMM and CSYMM.
695 *
696 * Auxiliary routine for test program for Level 3 Blas.
697 *
698 * -- Written on 8-February-1989.
699 * Jack Dongarra, Argonne National Laboratory.
700 * Iain Duff, AERE Harwell.
701 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
702 * Sven Hammarling, Numerical Algorithms Group Ltd.
703 *
704 * .. Parameters ..
705  COMPLEX zero
706  parameter( zero = ( 0.0, 0.0 ) )
707  REAL rzero
708  parameter( rzero = 0.0 )
709 * .. Scalar Arguments ..
710  REAL eps, thresh
711  INTEGER nalf, nbet, nidim, nmax, nout, ntra
712  LOGICAL fatal, rewi, trace
713  CHARACTER*6 sname
714 * .. Array Arguments ..
715  COMPLEX a( nmax, nmax ), aa( nmax*nmax ), alf( nalf ),
716  $ as( nmax*nmax ), b( nmax, nmax ),
717  $ bb( nmax*nmax ), bet( nbet ), bs( nmax*nmax ),
718  $ c( nmax, nmax ), cc( nmax*nmax ),
719  $ cs( nmax*nmax ), ct( nmax )
720  REAL g( nmax )
721  INTEGER idim( nidim )
722 * .. Local Scalars ..
723  COMPLEX alpha, als, beta, bls
724  REAL err, errmax
725  INTEGER i, ia, ib, ics, icu, im, in, laa, lbb, lcc,
726  $ lda, ldas, ldb, ldbs, ldc, ldcs, m, ms, n, na,
727  $ nargs, nc, ns
728  LOGICAL conj, left, null, reset, same
729  CHARACTER*1 side, sides, uplo, uplos
730  CHARACTER*2 ichs, ichu
731 * .. Local Arrays ..
732  LOGICAL isame( 13 )
733 * .. External Functions ..
734  LOGICAL lce, lceres
735  EXTERNAL lce, lceres
736 * .. External Subroutines ..
737  EXTERNAL chemm, cmake, cmmch, csymm
738 * .. Intrinsic Functions ..
739  INTRINSIC max
740 * .. Scalars in Common ..
741  INTEGER infot, noutc
742  LOGICAL lerr, ok
743 * .. Common blocks ..
744  common /infoc/infot, noutc, ok, lerr
745 * .. Data statements ..
746  DATA ichs/'LR'/, ichu/'UL'/
747 * .. Executable Statements ..
748  conj = sname( 2: 3 ).EQ.'HE'
749 *
750  nargs = 12
751  nc = 0
752  reset = .true.
753  errmax = rzero
754 *
755  DO 100 im = 1, nidim
756  m = idim( im )
757 *
758  DO 90 in = 1, nidim
759  n = idim( in )
760 * Set LDC to 1 more than minimum value if room.
761  ldc = m
762  IF( ldc.LT.nmax )
763  $ ldc = ldc + 1
764 * Skip tests if not enough room.
765  IF( ldc.GT.nmax )
766  $ go to 90
767  lcc = ldc*n
768  null = n.LE.0.OR.m.LE.0
769 * Set LDB to 1 more than minimum value if room.
770  ldb = m
771  IF( ldb.LT.nmax )
772  $ ldb = ldb + 1
773 * Skip tests if not enough room.
774  IF( ldb.GT.nmax )
775  $ go to 90
776  lbb = ldb*n
777 *
778 * Generate the matrix B.
779 *
780  CALL cmake( 'GE', ' ', ' ', m, n, b, nmax, bb, ldb, reset,
781  $ zero )
782 *
783  DO 80 ics = 1, 2
784  side = ichs( ics: ics )
785  left = side.EQ.'L'
786 *
787  IF( left )THEN
788  na = m
789  ELSE
790  na = n
791  END IF
792 * Set LDA to 1 more than minimum value if room.
793  lda = na
794  IF( lda.LT.nmax )
795  $ lda = lda + 1
796 * Skip tests if not enough room.
797  IF( lda.GT.nmax )
798  $ go to 80
799  laa = lda*na
800 *
801  DO 70 icu = 1, 2
802  uplo = ichu( icu: icu )
803 *
804 * Generate the hermitian or symmetric matrix A.
805 *
806  CALL cmake( sname( 2: 3 ), uplo, ' ', na, na, a, nmax,
807  $ aa, lda, reset, zero )
808 *
809  DO 60 ia = 1, nalf
810  alpha = alf( ia )
811 *
812  DO 50 ib = 1, nbet
813  beta = bet( ib )
814 *
815 * Generate the matrix C.
816 *
817  CALL cmake( 'GE', ' ', ' ', m, n, c, nmax, cc,
818  $ ldc, reset, zero )
819 *
820  nc = nc + 1
821 *
822 * Save every datum before calling the
823 * subroutine.
824 *
825  sides = side
826  uplos = uplo
827  ms = m
828  ns = n
829  als = alpha
830  DO 10 i = 1, laa
831  as( i ) = aa( i )
832  10 continue
833  ldas = lda
834  DO 20 i = 1, lbb
835  bs( i ) = bb( i )
836  20 continue
837  ldbs = ldb
838  bls = beta
839  DO 30 i = 1, lcc
840  cs( i ) = cc( i )
841  30 continue
842  ldcs = ldc
843 *
844 * Call the subroutine.
845 *
846  IF( trace )
847  $ WRITE( ntra, fmt = 9995 )nc, sname, side,
848  $ uplo, m, n, alpha, lda, ldb, beta, ldc
849  IF( rewi )
850  $ rewind ntra
851  IF( conj )THEN
852  CALL chemm( side, uplo, m, n, alpha, aa, lda,
853  $ bb, ldb, beta, cc, ldc )
854  ELSE
855  CALL csymm( side, uplo, m, n, alpha, aa, lda,
856  $ bb, ldb, beta, cc, ldc )
857  END IF
858 *
859 * Check if error-exit was taken incorrectly.
860 *
861  IF( .NOT.ok )THEN
862  WRITE( nout, fmt = 9994 )
863  fatal = .true.
864  go to 110
865  END IF
866 *
867 * See what data changed inside subroutines.
868 *
869  isame( 1 ) = sides.EQ.side
870  isame( 2 ) = uplos.EQ.uplo
871  isame( 3 ) = ms.EQ.m
872  isame( 4 ) = ns.EQ.n
873  isame( 5 ) = als.EQ.alpha
874  isame( 6 ) = lce( as, aa, laa )
875  isame( 7 ) = ldas.EQ.lda
876  isame( 8 ) = lce( bs, bb, lbb )
877  isame( 9 ) = ldbs.EQ.ldb
878  isame( 10 ) = bls.EQ.beta
879  IF( null )THEN
880  isame( 11 ) = lce( cs, cc, lcc )
881  ELSE
882  isame( 11 ) = lceres( 'GE', ' ', m, n, cs,
883  $ cc, ldc )
884  END IF
885  isame( 12 ) = ldcs.EQ.ldc
886 *
887 * If data was incorrectly changed, report and
888 * return.
889 *
890  same = .true.
891  DO 40 i = 1, nargs
892  same = same.AND.isame( i )
893  IF( .NOT.isame( i ) )
894  $ WRITE( nout, fmt = 9998 )i
895  40 continue
896  IF( .NOT.same )THEN
897  fatal = .true.
898  go to 110
899  END IF
900 *
901  IF( .NOT.null )THEN
902 *
903 * Check the result.
904 *
905  IF( left )THEN
906  CALL cmmch( 'N', 'N', m, n, m, alpha, a,
907  $ nmax, b, nmax, beta, c, nmax,
908  $ ct, g, cc, ldc, eps, err,
909  $ fatal, nout, .true. )
910  ELSE
911  CALL cmmch( 'N', 'N', m, n, n, alpha, b,
912  $ nmax, a, nmax, beta, c, nmax,
913  $ ct, g, cc, ldc, eps, err,
914  $ fatal, nout, .true. )
915  END IF
916  errmax = max( errmax, err )
917 * If got really bad answer, report and
918 * return.
919  IF( fatal )
920  $ go to 110
921  END IF
922 *
923  50 continue
924 *
925  60 continue
926 *
927  70 continue
928 *
929  80 continue
930 *
931  90 continue
932 *
933  100 continue
934 *
935 * Report result.
936 *
937  IF( errmax.LT.thresh )THEN
938  WRITE( nout, fmt = 9999 )sname, nc
939  ELSE
940  WRITE( nout, fmt = 9997 )sname, nc, errmax
941  END IF
942  go to 120
943 *
944  110 continue
945  WRITE( nout, fmt = 9996 )sname
946  WRITE( nout, fmt = 9995 )nc, sname, side, uplo, m, n, alpha, lda,
947  $ ldb, beta, ldc
948 *
949  120 continue
950  return
951 *
952  9999 format( ' ', a6, ' PASSED THE COMPUTATIONAL TESTS (', i6, ' CALL',
953  $ 'S)' )
954  9998 format( ' ******* FATAL ERROR - PARAMETER NUMBER ', i2, ' WAS CH',
955  $ 'ANGED INCORRECTLY *******' )
956  9997 format( ' ', a6, ' COMPLETED THE COMPUTATIONAL TESTS (', i6, ' C',
957  $ 'ALLS)', /' ******* BUT WITH MAXIMUM TEST RATIO', f8.2,
958  $ ' - SUSPECT *******' )
959  9996 format( ' ******* ', a6, ' FAILED ON CALL NUMBER:' )
960  9995 format( 1x, i6, ': ', a6, '(', 2( '''', a1, ''',' ), 2( i3, ',' ),
961  $ '(', f4.1, ',', f4.1, '), A,', i3, ', B,', i3, ',(', f4.1,
962  $ ',', f4.1, '), C,', i3, ') .' )
963  9994 format( ' ******* FATAL ERROR - ERROR-EXIT TAKEN ON VALID CALL *',
964  $ '******' )
965 *
966 * End of CCHK2.
967 *
968  END
969  SUBROUTINE cchk3( SNAME, EPS, THRESH, NOUT, NTRA, TRACE, REWI,
970  $ fatal, nidim, idim, nalf, alf, nmax, a, aa, as,
971  $ b, bb, bs, ct, g, c )
972 *
973 * Tests CTRMM and CTRSM.
974 *
975 * Auxiliary routine for test program for Level 3 Blas.
976 *
977 * -- Written on 8-February-1989.
978 * Jack Dongarra, Argonne National Laboratory.
979 * Iain Duff, AERE Harwell.
980 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
981 * Sven Hammarling, Numerical Algorithms Group Ltd.
982 *
983 * .. Parameters ..
984  COMPLEX zero, one
985  parameter( zero = ( 0.0, 0.0 ), one = ( 1.0, 0.0 ) )
986  REAL rzero
987  parameter( rzero = 0.0 )
988 * .. Scalar Arguments ..
989  REAL eps, thresh
990  INTEGER nalf, nidim, nmax, nout, ntra
991  LOGICAL fatal, rewi, trace
992  CHARACTER*6 sname
993 * .. Array Arguments ..
994  COMPLEX a( nmax, nmax ), aa( nmax*nmax ), alf( nalf ),
995  $ as( nmax*nmax ), b( nmax, nmax ),
996  $ bb( nmax*nmax ), bs( nmax*nmax ),
997  $ c( nmax, nmax ), ct( nmax )
998  REAL g( nmax )
999  INTEGER idim( nidim )
1000 * .. Local Scalars ..
1001  COMPLEX alpha, als
1002  REAL err, errmax
1003  INTEGER i, ia, icd, ics, ict, icu, im, in, j, laa, lbb,
1004  $ lda, ldas, ldb, ldbs, m, ms, n, na, nargs, nc,
1005  $ ns
1006  LOGICAL left, null, reset, same
1007  CHARACTER*1 diag, diags, side, sides, tranas, transa, uplo,
1008  $ uplos
1009  CHARACTER*2 ichd, ichs, ichu
1010  CHARACTER*3 icht
1011 * .. Local Arrays ..
1012  LOGICAL isame( 13 )
1013 * .. External Functions ..
1014  LOGICAL lce, lceres
1015  EXTERNAL lce, lceres
1016 * .. External Subroutines ..
1017  EXTERNAL cmake, cmmch, ctrmm, ctrsm
1018 * .. Intrinsic Functions ..
1019  INTRINSIC max
1020 * .. Scalars in Common ..
1021  INTEGER infot, noutc
1022  LOGICAL lerr, ok
1023 * .. Common blocks ..
1024  common /infoc/infot, noutc, ok, lerr
1025 * .. Data statements ..
1026  DATA ichu/'UL'/, icht/'NTC'/, ichd/'UN'/, ichs/'LR'/
1027 * .. Executable Statements ..
1028 *
1029  nargs = 11
1030  nc = 0
1031  reset = .true.
1032  errmax = rzero
1033 * Set up zero matrix for CMMCH.
1034  DO 20 j = 1, nmax
1035  DO 10 i = 1, nmax
1036  c( i, j ) = zero
1037  10 continue
1038  20 continue
1039 *
1040  DO 140 im = 1, nidim
1041  m = idim( im )
1042 *
1043  DO 130 in = 1, nidim
1044  n = idim( in )
1045 * Set LDB to 1 more than minimum value if room.
1046  ldb = m
1047  IF( ldb.LT.nmax )
1048  $ ldb = ldb + 1
1049 * Skip tests if not enough room.
1050  IF( ldb.GT.nmax )
1051  $ go to 130
1052  lbb = ldb*n
1053  null = m.LE.0.OR.n.LE.0
1054 *
1055  DO 120 ics = 1, 2
1056  side = ichs( ics: ics )
1057  left = side.EQ.'L'
1058  IF( left )THEN
1059  na = m
1060  ELSE
1061  na = n
1062  END IF
1063 * Set LDA to 1 more than minimum value if room.
1064  lda = na
1065  IF( lda.LT.nmax )
1066  $ lda = lda + 1
1067 * Skip tests if not enough room.
1068  IF( lda.GT.nmax )
1069  $ go to 130
1070  laa = lda*na
1071 *
1072  DO 110 icu = 1, 2
1073  uplo = ichu( icu: icu )
1074 *
1075  DO 100 ict = 1, 3
1076  transa = icht( ict: ict )
1077 *
1078  DO 90 icd = 1, 2
1079  diag = ichd( icd: icd )
1080 *
1081  DO 80 ia = 1, nalf
1082  alpha = alf( ia )
1083 *
1084 * Generate the matrix A.
1085 *
1086  CALL cmake( 'TR', uplo, diag, na, na, a,
1087  $ nmax, aa, lda, reset, zero )
1088 *
1089 * Generate the matrix B.
1090 *
1091  CALL cmake( 'GE', ' ', ' ', m, n, b, nmax,
1092  $ bb, ldb, reset, zero )
1093 *
1094  nc = nc + 1
1095 *
1096 * Save every datum before calling the
1097 * subroutine.
1098 *
1099  sides = side
1100  uplos = uplo
1101  tranas = transa
1102  diags = diag
1103  ms = m
1104  ns = n
1105  als = alpha
1106  DO 30 i = 1, laa
1107  as( i ) = aa( i )
1108  30 continue
1109  ldas = lda
1110  DO 40 i = 1, lbb
1111  bs( i ) = bb( i )
1112  40 continue
1113  ldbs = ldb
1114 *
1115 * Call the subroutine.
1116 *
1117  IF( sname( 4: 5 ).EQ.'MM' )THEN
1118  IF( trace )
1119  $ WRITE( ntra, fmt = 9995 )nc, sname,
1120  $ side, uplo, transa, diag, m, n, alpha,
1121  $ lda, ldb
1122  IF( rewi )
1123  $ rewind ntra
1124  CALL ctrmm( side, uplo, transa, diag, m,
1125  $ n, alpha, aa, lda, bb, ldb )
1126  ELSE IF( sname( 4: 5 ).EQ.'SM' )THEN
1127  IF( trace )
1128  $ WRITE( ntra, fmt = 9995 )nc, sname,
1129  $ side, uplo, transa, diag, m, n, alpha,
1130  $ lda, ldb
1131  IF( rewi )
1132  $ rewind ntra
1133  CALL ctrsm( side, uplo, transa, diag, m,
1134  $ n, alpha, aa, lda, bb, ldb )
1135  END IF
1136 *
1137 * Check if error-exit was taken incorrectly.
1138 *
1139  IF( .NOT.ok )THEN
1140  WRITE( nout, fmt = 9994 )
1141  fatal = .true.
1142  go to 150
1143  END IF
1144 *
1145 * See what data changed inside subroutines.
1146 *
1147  isame( 1 ) = sides.EQ.side
1148  isame( 2 ) = uplos.EQ.uplo
1149  isame( 3 ) = tranas.EQ.transa
1150  isame( 4 ) = diags.EQ.diag
1151  isame( 5 ) = ms.EQ.m
1152  isame( 6 ) = ns.EQ.n
1153  isame( 7 ) = als.EQ.alpha
1154  isame( 8 ) = lce( as, aa, laa )
1155  isame( 9 ) = ldas.EQ.lda
1156  IF( null )THEN
1157  isame( 10 ) = lce( bs, bb, lbb )
1158  ELSE
1159  isame( 10 ) = lceres( 'GE', ' ', m, n, bs,
1160  $ bb, ldb )
1161  END IF
1162  isame( 11 ) = ldbs.EQ.ldb
1163 *
1164 * If data was incorrectly changed, report and
1165 * return.
1166 *
1167  same = .true.
1168  DO 50 i = 1, nargs
1169  same = same.AND.isame( i )
1170  IF( .NOT.isame( i ) )
1171  $ WRITE( nout, fmt = 9998 )i
1172  50 continue
1173  IF( .NOT.same )THEN
1174  fatal = .true.
1175  go to 150
1176  END IF
1177 *
1178  IF( .NOT.null )THEN
1179  IF( sname( 4: 5 ).EQ.'MM' )THEN
1180 *
1181 * Check the result.
1182 *
1183  IF( left )THEN
1184  CALL cmmch( transa, 'N', m, n, m,
1185  $ alpha, a, nmax, b, nmax,
1186  $ zero, c, nmax, ct, g,
1187  $ bb, ldb, eps, err,
1188  $ fatal, nout, .true. )
1189  ELSE
1190  CALL cmmch( 'N', transa, m, n, n,
1191  $ alpha, b, nmax, a, nmax,
1192  $ zero, c, nmax, ct, g,
1193  $ bb, ldb, eps, err,
1194  $ fatal, nout, .true. )
1195  END IF
1196  ELSE IF( sname( 4: 5 ).EQ.'SM' )THEN
1197 *
1198 * Compute approximation to original
1199 * matrix.
1200 *
1201  DO 70 j = 1, n
1202  DO 60 i = 1, m
1203  c( i, j ) = bb( i + ( j - 1 )*
1204  $ ldb )
1205  bb( i + ( j - 1 )*ldb ) = alpha*
1206  $ b( i, j )
1207  60 continue
1208  70 continue
1209 *
1210  IF( left )THEN
1211  CALL cmmch( transa, 'N', m, n, m,
1212  $ one, a, nmax, c, nmax,
1213  $ zero, b, nmax, ct, g,
1214  $ bb, ldb, eps, err,
1215  $ fatal, nout, .false. )
1216  ELSE
1217  CALL cmmch( 'N', transa, m, n, n,
1218  $ one, c, nmax, a, nmax,
1219  $ zero, b, nmax, ct, g,
1220  $ bb, ldb, eps, err,
1221  $ fatal, nout, .false. )
1222  END IF
1223  END IF
1224  errmax = max( errmax, err )
1225 * If got really bad answer, report and
1226 * return.
1227  IF( fatal )
1228  $ go to 150
1229  END IF
1230 *
1231  80 continue
1232 *
1233  90 continue
1234 *
1235  100 continue
1236 *
1237  110 continue
1238 *
1239  120 continue
1240 *
1241  130 continue
1242 *
1243  140 continue
1244 *
1245 * Report result.
1246 *
1247  IF( errmax.LT.thresh )THEN
1248  WRITE( nout, fmt = 9999 )sname, nc
1249  ELSE
1250  WRITE( nout, fmt = 9997 )sname, nc, errmax
1251  END IF
1252  go to 160
1253 *
1254  150 continue
1255  WRITE( nout, fmt = 9996 )sname
1256  WRITE( nout, fmt = 9995 )nc, sname, side, uplo, transa, diag, m,
1257  $ n, alpha, lda, ldb
1258 *
1259  160 continue
1260  return
1261 *
1262  9999 format( ' ', a6, ' PASSED THE COMPUTATIONAL TESTS (', i6, ' CALL',
1263  $ 'S)' )
1264  9998 format( ' ******* FATAL ERROR - PARAMETER NUMBER ', i2, ' WAS CH',
1265  $ 'ANGED INCORRECTLY *******' )
1266  9997 format( ' ', a6, ' COMPLETED THE COMPUTATIONAL TESTS (', i6, ' C',
1267  $ 'ALLS)', /' ******* BUT WITH MAXIMUM TEST RATIO', f8.2,
1268  $ ' - SUSPECT *******' )
1269  9996 format( ' ******* ', a6, ' FAILED ON CALL NUMBER:' )
1270  9995 format( 1x, i6, ': ', a6, '(', 4( '''', a1, ''',' ), 2( i3, ',' ),
1271  $ '(', f4.1, ',', f4.1, '), A,', i3, ', B,', i3, ') ',
1272  $ ' .' )
1273  9994 format( ' ******* FATAL ERROR - ERROR-EXIT TAKEN ON VALID CALL *',
1274  $ '******' )
1275 *
1276 * End of CCHK3.
1277 *
1278  END
1279  SUBROUTINE cchk4( SNAME, EPS, THRESH, NOUT, NTRA, TRACE, REWI,
1280  $ fatal, nidim, idim, nalf, alf, nbet, bet, nmax,
1281  $ a, aa, as, b, bb, bs, c, cc, cs, ct, g )
1282 *
1283 * Tests CHERK and CSYRK.
1284 *
1285 * Auxiliary routine for test program for Level 3 Blas.
1286 *
1287 * -- Written on 8-February-1989.
1288 * Jack Dongarra, Argonne National Laboratory.
1289 * Iain Duff, AERE Harwell.
1290 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
1291 * Sven Hammarling, Numerical Algorithms Group Ltd.
1292 *
1293 * .. Parameters ..
1294  COMPLEX zero
1295  parameter( zero = ( 0.0, 0.0 ) )
1296  REAL rone, rzero
1297  parameter( rone = 1.0, rzero = 0.0 )
1298 * .. Scalar Arguments ..
1299  REAL eps, thresh
1300  INTEGER nalf, nbet, nidim, nmax, nout, ntra
1301  LOGICAL fatal, rewi, trace
1302  CHARACTER*6 sname
1303 * .. Array Arguments ..
1304  COMPLEX a( nmax, nmax ), aa( nmax*nmax ), alf( nalf ),
1305  $ as( nmax*nmax ), b( nmax, nmax ),
1306  $ bb( nmax*nmax ), bet( nbet ), bs( nmax*nmax ),
1307  $ c( nmax, nmax ), cc( nmax*nmax ),
1308  $ cs( nmax*nmax ), ct( nmax )
1309  REAL g( nmax )
1310  INTEGER idim( nidim )
1311 * .. Local Scalars ..
1312  COMPLEX alpha, als, beta, bets
1313  REAL err, errmax, ralpha, rals, rbeta, rbets
1314  INTEGER i, ia, ib, ict, icu, ik, in, j, jc, jj, k, ks,
1315  $ laa, lcc, lda, ldas, ldc, ldcs, lj, ma, n, na,
1316  $ nargs, nc, ns
1317  LOGICAL conj, null, reset, same, tran, upper
1318  CHARACTER*1 trans, transs, transt, uplo, uplos
1319  CHARACTER*2 icht, ichu
1320 * .. Local Arrays ..
1321  LOGICAL isame( 13 )
1322 * .. External Functions ..
1323  LOGICAL lce, lceres
1324  EXTERNAL lce, lceres
1325 * .. External Subroutines ..
1326  EXTERNAL cherk, cmake, cmmch, csyrk
1327 * .. Intrinsic Functions ..
1328  INTRINSIC cmplx, max, real
1329 * .. Scalars in Common ..
1330  INTEGER infot, noutc
1331  LOGICAL lerr, ok
1332 * .. Common blocks ..
1333  common /infoc/infot, noutc, ok, lerr
1334 * .. Data statements ..
1335  DATA icht/'NC'/, ichu/'UL'/
1336 * .. Executable Statements ..
1337  conj = sname( 2: 3 ).EQ.'HE'
1338 *
1339  nargs = 10
1340  nc = 0
1341  reset = .true.
1342  errmax = rzero
1343 *
1344  DO 100 in = 1, nidim
1345  n = idim( in )
1346 * Set LDC to 1 more than minimum value if room.
1347  ldc = n
1348  IF( ldc.LT.nmax )
1349  $ ldc = ldc + 1
1350 * Skip tests if not enough room.
1351  IF( ldc.GT.nmax )
1352  $ go to 100
1353  lcc = ldc*n
1354 *
1355  DO 90 ik = 1, nidim
1356  k = idim( ik )
1357 *
1358  DO 80 ict = 1, 2
1359  trans = icht( ict: ict )
1360  tran = trans.EQ.'C'
1361  IF( tran.AND..NOT.conj )
1362  $ trans = 'T'
1363  IF( tran )THEN
1364  ma = k
1365  na = n
1366  ELSE
1367  ma = n
1368  na = k
1369  END IF
1370 * Set LDA to 1 more than minimum value if room.
1371  lda = ma
1372  IF( lda.LT.nmax )
1373  $ lda = lda + 1
1374 * Skip tests if not enough room.
1375  IF( lda.GT.nmax )
1376  $ go to 80
1377  laa = lda*na
1378 *
1379 * Generate the matrix A.
1380 *
1381  CALL cmake( 'GE', ' ', ' ', ma, na, a, nmax, aa, lda,
1382  $ reset, zero )
1383 *
1384  DO 70 icu = 1, 2
1385  uplo = ichu( icu: icu )
1386  upper = uplo.EQ.'U'
1387 *
1388  DO 60 ia = 1, nalf
1389  alpha = alf( ia )
1390  IF( conj )THEN
1391  ralpha = REAL( alpha )
1392  alpha = cmplx( ralpha, rzero )
1393  END IF
1394 *
1395  DO 50 ib = 1, nbet
1396  beta = bet( ib )
1397  IF( conj )THEN
1398  rbeta = REAL( beta )
1399  beta = cmplx( rbeta, rzero )
1400  END IF
1401  null = n.LE.0
1402  IF( conj )
1403  $ null = null.OR.( ( k.LE.0.OR.ralpha.EQ.
1404  $ rzero ).AND.rbeta.EQ.rone )
1405 *
1406 * Generate the matrix C.
1407 *
1408  CALL cmake( sname( 2: 3 ), uplo, ' ', n, n, c,
1409  $ nmax, cc, ldc, reset, zero )
1410 *
1411  nc = nc + 1
1412 *
1413 * Save every datum before calling the subroutine.
1414 *
1415  uplos = uplo
1416  transs = trans
1417  ns = n
1418  ks = k
1419  IF( conj )THEN
1420  rals = ralpha
1421  ELSE
1422  als = alpha
1423  END IF
1424  DO 10 i = 1, laa
1425  as( i ) = aa( i )
1426  10 continue
1427  ldas = lda
1428  IF( conj )THEN
1429  rbets = rbeta
1430  ELSE
1431  bets = beta
1432  END IF
1433  DO 20 i = 1, lcc
1434  cs( i ) = cc( i )
1435  20 continue
1436  ldcs = ldc
1437 *
1438 * Call the subroutine.
1439 *
1440  IF( conj )THEN
1441  IF( trace )
1442  $ WRITE( ntra, fmt = 9994 )nc, sname, uplo,
1443  $ trans, n, k, ralpha, lda, rbeta, ldc
1444  IF( rewi )
1445  $ rewind ntra
1446  CALL cherk( uplo, trans, n, k, ralpha, aa,
1447  $ lda, rbeta, cc, ldc )
1448  ELSE
1449  IF( trace )
1450  $ WRITE( ntra, fmt = 9993 )nc, sname, uplo,
1451  $ trans, n, k, alpha, lda, beta, ldc
1452  IF( rewi )
1453  $ rewind ntra
1454  CALL csyrk( uplo, trans, n, k, alpha, aa,
1455  $ lda, beta, cc, ldc )
1456  END IF
1457 *
1458 * Check if error-exit was taken incorrectly.
1459 *
1460  IF( .NOT.ok )THEN
1461  WRITE( nout, fmt = 9992 )
1462  fatal = .true.
1463  go to 120
1464  END IF
1465 *
1466 * See what data changed inside subroutines.
1467 *
1468  isame( 1 ) = uplos.EQ.uplo
1469  isame( 2 ) = transs.EQ.trans
1470  isame( 3 ) = ns.EQ.n
1471  isame( 4 ) = ks.EQ.k
1472  IF( conj )THEN
1473  isame( 5 ) = rals.EQ.ralpha
1474  ELSE
1475  isame( 5 ) = als.EQ.alpha
1476  END IF
1477  isame( 6 ) = lce( as, aa, laa )
1478  isame( 7 ) = ldas.EQ.lda
1479  IF( conj )THEN
1480  isame( 8 ) = rbets.EQ.rbeta
1481  ELSE
1482  isame( 8 ) = bets.EQ.beta
1483  END IF
1484  IF( null )THEN
1485  isame( 9 ) = lce( cs, cc, lcc )
1486  ELSE
1487  isame( 9 ) = lceres( sname( 2: 3 ), uplo, n,
1488  $ n, cs, cc, ldc )
1489  END IF
1490  isame( 10 ) = ldcs.EQ.ldc
1491 *
1492 * If data was incorrectly changed, report and
1493 * return.
1494 *
1495  same = .true.
1496  DO 30 i = 1, nargs
1497  same = same.AND.isame( i )
1498  IF( .NOT.isame( i ) )
1499  $ WRITE( nout, fmt = 9998 )i
1500  30 continue
1501  IF( .NOT.same )THEN
1502  fatal = .true.
1503  go to 120
1504  END IF
1505 *
1506  IF( .NOT.null )THEN
1507 *
1508 * Check the result column by column.
1509 *
1510  IF( conj )THEN
1511  transt = 'C'
1512  ELSE
1513  transt = 'T'
1514  END IF
1515  jc = 1
1516  DO 40 j = 1, n
1517  IF( upper )THEN
1518  jj = 1
1519  lj = j
1520  ELSE
1521  jj = j
1522  lj = n - j + 1
1523  END IF
1524  IF( tran )THEN
1525  CALL cmmch( transt, 'N', lj, 1, k,
1526  $ alpha, a( 1, jj ), nmax,
1527  $ a( 1, j ), nmax, beta,
1528  $ c( jj, j ), nmax, ct, g,
1529  $ cc( jc ), ldc, eps, err,
1530  $ fatal, nout, .true. )
1531  ELSE
1532  CALL cmmch( 'N', transt, lj, 1, k,
1533  $ alpha, a( jj, 1 ), nmax,
1534  $ a( j, 1 ), nmax, beta,
1535  $ c( jj, j ), nmax, ct, g,
1536  $ cc( jc ), ldc, eps, err,
1537  $ fatal, nout, .true. )
1538  END IF
1539  IF( upper )THEN
1540  jc = jc + ldc
1541  ELSE
1542  jc = jc + ldc + 1
1543  END IF
1544  errmax = max( errmax, err )
1545 * If got really bad answer, report and
1546 * return.
1547  IF( fatal )
1548  $ go to 110
1549  40 continue
1550  END IF
1551 *
1552  50 continue
1553 *
1554  60 continue
1555 *
1556  70 continue
1557 *
1558  80 continue
1559 *
1560  90 continue
1561 *
1562  100 continue
1563 *
1564 * Report result.
1565 *
1566  IF( errmax.LT.thresh )THEN
1567  WRITE( nout, fmt = 9999 )sname, nc
1568  ELSE
1569  WRITE( nout, fmt = 9997 )sname, nc, errmax
1570  END IF
1571  go to 130
1572 *
1573  110 continue
1574  IF( n.GT.1 )
1575  $ WRITE( nout, fmt = 9995 )j
1576 *
1577  120 continue
1578  WRITE( nout, fmt = 9996 )sname
1579  IF( conj )THEN
1580  WRITE( nout, fmt = 9994 )nc, sname, uplo, trans, n, k, ralpha,
1581  $ lda, rbeta, ldc
1582  ELSE
1583  WRITE( nout, fmt = 9993 )nc, sname, uplo, trans, n, k, alpha,
1584  $ lda, beta, ldc
1585  END IF
1586 *
1587  130 continue
1588  return
1589 *
1590  9999 format( ' ', a6, ' PASSED THE COMPUTATIONAL TESTS (', i6, ' CALL',
1591  $ 'S)' )
1592  9998 format( ' ******* FATAL ERROR - PARAMETER NUMBER ', i2, ' WAS CH',
1593  $ 'ANGED INCORRECTLY *******' )
1594  9997 format( ' ', a6, ' COMPLETED THE COMPUTATIONAL TESTS (', i6, ' C',
1595  $ 'ALLS)', /' ******* BUT WITH MAXIMUM TEST RATIO', f8.2,
1596  $ ' - SUSPECT *******' )
1597  9996 format( ' ******* ', a6, ' FAILED ON CALL NUMBER:' )
1598  9995 format( ' THESE ARE THE RESULTS FOR COLUMN ', i3 )
1599  9994 format( 1x, i6, ': ', a6, '(', 2( '''', a1, ''',' ), 2( i3, ',' ),
1600  $ f4.1, ', A,', i3, ',', f4.1, ', C,', i3, ') ',
1601  $ ' .' )
1602  9993 format( 1x, i6, ': ', a6, '(', 2( '''', a1, ''',' ), 2( i3, ',' ),
1603  $ '(', f4.1, ',', f4.1, ') , A,', i3, ',(', f4.1, ',', f4.1,
1604  $ '), C,', i3, ') .' )
1605  9992 format( ' ******* FATAL ERROR - ERROR-EXIT TAKEN ON VALID CALL *',
1606  $ '******' )
1607 *
1608 * End of CCHK4.
1609 *
1610  END
1611  SUBROUTINE cchk5( SNAME, EPS, THRESH, NOUT, NTRA, TRACE, REWI,
1612  $ fatal, nidim, idim, nalf, alf, nbet, bet, nmax,
1613  $ ab, aa, as, bb, bs, c, cc, cs, ct, g, w )
1614 *
1615 * Tests CHER2K and CSYR2K.
1616 *
1617 * Auxiliary routine for test program for Level 3 Blas.
1618 *
1619 * -- Written on 8-February-1989.
1620 * Jack Dongarra, Argonne National Laboratory.
1621 * Iain Duff, AERE Harwell.
1622 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
1623 * Sven Hammarling, Numerical Algorithms Group Ltd.
1624 *
1625 * .. Parameters ..
1626  COMPLEX zero, one
1627  parameter( zero = ( 0.0, 0.0 ), one = ( 1.0, 0.0 ) )
1628  REAL rone, rzero
1629  parameter( rone = 1.0, rzero = 0.0 )
1630 * .. Scalar Arguments ..
1631  REAL eps, thresh
1632  INTEGER nalf, nbet, nidim, nmax, nout, ntra
1633  LOGICAL fatal, rewi, trace
1634  CHARACTER*6 sname
1635 * .. Array Arguments ..
1636  COMPLEX aa( nmax*nmax ), ab( 2*nmax*nmax ),
1637  $ alf( nalf ), as( nmax*nmax ), bb( nmax*nmax ),
1638  $ bet( nbet ), bs( nmax*nmax ), c( nmax, nmax ),
1639  $ cc( nmax*nmax ), cs( nmax*nmax ), ct( nmax ),
1640  $ w( 2*nmax )
1641  REAL g( nmax )
1642  INTEGER idim( nidim )
1643 * .. Local Scalars ..
1644  COMPLEX alpha, als, beta, bets
1645  REAL err, errmax, rbeta, rbets
1646  INTEGER i, ia, ib, ict, icu, ik, in, j, jc, jj, jjab,
1647  $ k, ks, laa, lbb, lcc, lda, ldas, ldb, ldbs,
1648  $ ldc, ldcs, lj, ma, n, na, nargs, nc, ns
1649  LOGICAL conj, null, reset, same, tran, upper
1650  CHARACTER*1 trans, transs, transt, uplo, uplos
1651  CHARACTER*2 icht, ichu
1652 * .. Local Arrays ..
1653  LOGICAL isame( 13 )
1654 * .. External Functions ..
1655  LOGICAL lce, lceres
1656  EXTERNAL lce, lceres
1657 * .. External Subroutines ..
1658  EXTERNAL cher2k, cmake, cmmch, csyr2k
1659 * .. Intrinsic Functions ..
1660  INTRINSIC cmplx, conjg, max, real
1661 * .. Scalars in Common ..
1662  INTEGER infot, noutc
1663  LOGICAL lerr, ok
1664 * .. Common blocks ..
1665  common /infoc/infot, noutc, ok, lerr
1666 * .. Data statements ..
1667  DATA icht/'NC'/, ichu/'UL'/
1668 * .. Executable Statements ..
1669  conj = sname( 2: 3 ).EQ.'HE'
1670 *
1671  nargs = 12
1672  nc = 0
1673  reset = .true.
1674  errmax = rzero
1675 *
1676  DO 130 in = 1, nidim
1677  n = idim( in )
1678 * Set LDC to 1 more than minimum value if room.
1679  ldc = n
1680  IF( ldc.LT.nmax )
1681  $ ldc = ldc + 1
1682 * Skip tests if not enough room.
1683  IF( ldc.GT.nmax )
1684  $ go to 130
1685  lcc = ldc*n
1686 *
1687  DO 120 ik = 1, nidim
1688  k = idim( ik )
1689 *
1690  DO 110 ict = 1, 2
1691  trans = icht( ict: ict )
1692  tran = trans.EQ.'C'
1693  IF( tran.AND..NOT.conj )
1694  $ trans = 'T'
1695  IF( tran )THEN
1696  ma = k
1697  na = n
1698  ELSE
1699  ma = n
1700  na = k
1701  END IF
1702 * Set LDA to 1 more than minimum value if room.
1703  lda = ma
1704  IF( lda.LT.nmax )
1705  $ lda = lda + 1
1706 * Skip tests if not enough room.
1707  IF( lda.GT.nmax )
1708  $ go to 110
1709  laa = lda*na
1710 *
1711 * Generate the matrix A.
1712 *
1713  IF( tran )THEN
1714  CALL cmake( 'GE', ' ', ' ', ma, na, ab, 2*nmax, aa,
1715  $ lda, reset, zero )
1716  ELSE
1717  CALL cmake( 'GE', ' ', ' ', ma, na, ab, nmax, aa, lda,
1718  $ reset, zero )
1719  END IF
1720 *
1721 * Generate the matrix B.
1722 *
1723  ldb = lda
1724  lbb = laa
1725  IF( tran )THEN
1726  CALL cmake( 'GE', ' ', ' ', ma, na, ab( k + 1 ),
1727  $ 2*nmax, bb, ldb, reset, zero )
1728  ELSE
1729  CALL cmake( 'GE', ' ', ' ', ma, na, ab( k*nmax + 1 ),
1730  $ nmax, bb, ldb, reset, zero )
1731  END IF
1732 *
1733  DO 100 icu = 1, 2
1734  uplo = ichu( icu: icu )
1735  upper = uplo.EQ.'U'
1736 *
1737  DO 90 ia = 1, nalf
1738  alpha = alf( ia )
1739 *
1740  DO 80 ib = 1, nbet
1741  beta = bet( ib )
1742  IF( conj )THEN
1743  rbeta = REAL( beta )
1744  beta = cmplx( rbeta, rzero )
1745  END IF
1746  null = n.LE.0
1747  IF( conj )
1748  $ null = null.OR.( ( k.LE.0.OR.alpha.EQ.
1749  $ zero ).AND.rbeta.EQ.rone )
1750 *
1751 * Generate the matrix C.
1752 *
1753  CALL cmake( sname( 2: 3 ), uplo, ' ', n, n, c,
1754  $ nmax, cc, ldc, reset, zero )
1755 *
1756  nc = nc + 1
1757 *
1758 * Save every datum before calling the subroutine.
1759 *
1760  uplos = uplo
1761  transs = trans
1762  ns = n
1763  ks = k
1764  als = alpha
1765  DO 10 i = 1, laa
1766  as( i ) = aa( i )
1767  10 continue
1768  ldas = lda
1769  DO 20 i = 1, lbb
1770  bs( i ) = bb( i )
1771  20 continue
1772  ldbs = ldb
1773  IF( conj )THEN
1774  rbets = rbeta
1775  ELSE
1776  bets = beta
1777  END IF
1778  DO 30 i = 1, lcc
1779  cs( i ) = cc( i )
1780  30 continue
1781  ldcs = ldc
1782 *
1783 * Call the subroutine.
1784 *
1785  IF( conj )THEN
1786  IF( trace )
1787  $ WRITE( ntra, fmt = 9994 )nc, sname, uplo,
1788  $ trans, n, k, alpha, lda, ldb, rbeta, ldc
1789  IF( rewi )
1790  $ rewind ntra
1791  CALL cher2k( uplo, trans, n, k, alpha, aa,
1792  $ lda, bb, ldb, rbeta, cc, ldc )
1793  ELSE
1794  IF( trace )
1795  $ WRITE( ntra, fmt = 9993 )nc, sname, uplo,
1796  $ trans, n, k, alpha, lda, ldb, beta, ldc
1797  IF( rewi )
1798  $ rewind ntra
1799  CALL csyr2k( uplo, trans, n, k, alpha, aa,
1800  $ lda, bb, ldb, beta, cc, ldc )
1801  END IF
1802 *
1803 * Check if error-exit was taken incorrectly.
1804 *
1805  IF( .NOT.ok )THEN
1806  WRITE( nout, fmt = 9992 )
1807  fatal = .true.
1808  go to 150
1809  END IF
1810 *
1811 * See what data changed inside subroutines.
1812 *
1813  isame( 1 ) = uplos.EQ.uplo
1814  isame( 2 ) = transs.EQ.trans
1815  isame( 3 ) = ns.EQ.n
1816  isame( 4 ) = ks.EQ.k
1817  isame( 5 ) = als.EQ.alpha
1818  isame( 6 ) = lce( as, aa, laa )
1819  isame( 7 ) = ldas.EQ.lda
1820  isame( 8 ) = lce( bs, bb, lbb )
1821  isame( 9 ) = ldbs.EQ.ldb
1822  IF( conj )THEN
1823  isame( 10 ) = rbets.EQ.rbeta
1824  ELSE
1825  isame( 10 ) = bets.EQ.beta
1826  END IF
1827  IF( null )THEN
1828  isame( 11 ) = lce( cs, cc, lcc )
1829  ELSE
1830  isame( 11 ) = lceres( 'HE', uplo, n, n, cs,
1831  $ cc, ldc )
1832  END IF
1833  isame( 12 ) = ldcs.EQ.ldc
1834 *
1835 * If data was incorrectly changed, report and
1836 * return.
1837 *
1838  same = .true.
1839  DO 40 i = 1, nargs
1840  same = same.AND.isame( i )
1841  IF( .NOT.isame( i ) )
1842  $ WRITE( nout, fmt = 9998 )i
1843  40 continue
1844  IF( .NOT.same )THEN
1845  fatal = .true.
1846  go to 150
1847  END IF
1848 *
1849  IF( .NOT.null )THEN
1850 *
1851 * Check the result column by column.
1852 *
1853  IF( conj )THEN
1854  transt = 'C'
1855  ELSE
1856  transt = 'T'
1857  END IF
1858  jjab = 1
1859  jc = 1
1860  DO 70 j = 1, n
1861  IF( upper )THEN
1862  jj = 1
1863  lj = j
1864  ELSE
1865  jj = j
1866  lj = n - j + 1
1867  END IF
1868  IF( tran )THEN
1869  DO 50 i = 1, k
1870  w( i ) = alpha*ab( ( j - 1 )*2*
1871  $ nmax + k + i )
1872  IF( conj )THEN
1873  w( k + i ) = conjg( alpha )*
1874  $ ab( ( j - 1 )*2*
1875  $ nmax + i )
1876  ELSE
1877  w( k + i ) = alpha*
1878  $ ab( ( j - 1 )*2*
1879  $ nmax + i )
1880  END IF
1881  50 continue
1882  CALL cmmch( transt, 'N', lj, 1, 2*k,
1883  $ one, ab( jjab ), 2*nmax, w,
1884  $ 2*nmax, beta, c( jj, j ),
1885  $ nmax, ct, g, cc( jc ), ldc,
1886  $ eps, err, fatal, nout,
1887  $ .true. )
1888  ELSE
1889  DO 60 i = 1, k
1890  IF( conj )THEN
1891  w( i ) = alpha*conjg( ab( ( k +
1892  $ i - 1 )*nmax + j ) )
1893  w( k + i ) = conjg( alpha*
1894  $ ab( ( i - 1 )*nmax +
1895  $ j ) )
1896  ELSE
1897  w( i ) = alpha*ab( ( k + i - 1 )*
1898  $ nmax + j )
1899  w( k + i ) = alpha*
1900  $ ab( ( i - 1 )*nmax +
1901  $ j )
1902  END IF
1903  60 continue
1904  CALL cmmch( 'N', 'N', lj, 1, 2*k, one,
1905  $ ab( jj ), nmax, w, 2*nmax,
1906  $ beta, c( jj, j ), nmax, ct,
1907  $ g, cc( jc ), ldc, eps, err,
1908  $ fatal, nout, .true. )
1909  END IF
1910  IF( upper )THEN
1911  jc = jc + ldc
1912  ELSE
1913  jc = jc + ldc + 1
1914  IF( tran )
1915  $ jjab = jjab + 2*nmax
1916  END IF
1917  errmax = max( errmax, err )
1918 * If got really bad answer, report and
1919 * return.
1920  IF( fatal )
1921  $ go to 140
1922  70 continue
1923  END IF
1924 *
1925  80 continue
1926 *
1927  90 continue
1928 *
1929  100 continue
1930 *
1931  110 continue
1932 *
1933  120 continue
1934 *
1935  130 continue
1936 *
1937 * Report result.
1938 *
1939  IF( errmax.LT.thresh )THEN
1940  WRITE( nout, fmt = 9999 )sname, nc
1941  ELSE
1942  WRITE( nout, fmt = 9997 )sname, nc, errmax
1943  END IF
1944  go to 160
1945 *
1946  140 continue
1947  IF( n.GT.1 )
1948  $ WRITE( nout, fmt = 9995 )j
1949 *
1950  150 continue
1951  WRITE( nout, fmt = 9996 )sname
1952  IF( conj )THEN
1953  WRITE( nout, fmt = 9994 )nc, sname, uplo, trans, n, k, alpha,
1954  $ lda, ldb, rbeta, ldc
1955  ELSE
1956  WRITE( nout, fmt = 9993 )nc, sname, uplo, trans, n, k, alpha,
1957  $ lda, ldb, beta, ldc
1958  END IF
1959 *
1960  160 continue
1961  return
1962 *
1963  9999 format( ' ', a6, ' PASSED THE COMPUTATIONAL TESTS (', i6, ' CALL',
1964  $ 'S)' )
1965  9998 format( ' ******* FATAL ERROR - PARAMETER NUMBER ', i2, ' WAS CH',
1966  $ 'ANGED INCORRECTLY *******' )
1967  9997 format( ' ', a6, ' COMPLETED THE COMPUTATIONAL TESTS (', i6, ' C',
1968  $ 'ALLS)', /' ******* BUT WITH MAXIMUM TEST RATIO', f8.2,
1969  $ ' - SUSPECT *******' )
1970  9996 format( ' ******* ', a6, ' FAILED ON CALL NUMBER:' )
1971  9995 format( ' THESE ARE THE RESULTS FOR COLUMN ', i3 )
1972  9994 format( 1x, i6, ': ', a6, '(', 2( '''', a1, ''',' ), 2( i3, ',' ),
1973  $ '(', f4.1, ',', f4.1, '), A,', i3, ', B,', i3, ',', f4.1,
1974  $ ', C,', i3, ') .' )
1975  9993 format( 1x, i6, ': ', a6, '(', 2( '''', a1, ''',' ), 2( i3, ',' ),
1976  $ '(', f4.1, ',', f4.1, '), A,', i3, ', B,', i3, ',(', f4.1,
1977  $ ',', f4.1, '), C,', i3, ') .' )
1978  9992 format( ' ******* FATAL ERROR - ERROR-EXIT TAKEN ON VALID CALL *',
1979  $ '******' )
1980 *
1981 * End of CCHK5.
1982 *
1983  END
1984  SUBROUTINE cchke( ISNUM, SRNAMT, NOUT )
1985 *
1986 * Tests the error exits from the Level 3 Blas.
1987 * Requires a special version of the error-handling routine XERBLA.
1988 * A, B and C should not need to be defined.
1989 *
1990 * Auxiliary routine for test program for Level 3 Blas.
1991 *
1992 * -- Written on 8-February-1989.
1993 * Jack Dongarra, Argonne National Laboratory.
1994 * Iain Duff, AERE Harwell.
1995 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
1996 * Sven Hammarling, Numerical Algorithms Group Ltd.
1997 *
1998 * 3-19-92: Initialize ALPHA, BETA, RALPHA, and RBETA (eca)
1999 * 3-19-92: Fix argument 12 in calls to CSYMM and CHEMM
2000 * with INFOT = 9 (eca)
2001 *
2002 * .. Scalar Arguments ..
2003  INTEGER isnum, nout
2004  CHARACTER*6 srnamt
2005 * .. Scalars in Common ..
2006  INTEGER infot, noutc
2007  LOGICAL lerr, ok
2008 * .. Parameters ..
2009  REAL one, two
2010  parameter( one = 1.0e0, two = 2.0e0 )
2011 * .. Local Scalars ..
2012  COMPLEX alpha, beta
2013  REAL ralpha, rbeta
2014 * .. Local Arrays ..
2015  COMPLEX a( 2, 1 ), b( 2, 1 ), c( 2, 1 )
2016 * .. External Subroutines ..
2017  EXTERNAL cgemm, chemm, cher2k, cherk, chkxer, csymm,
2018  $ csyr2k, csyrk, ctrmm, ctrsm
2019 * .. Common blocks ..
2020  common /infoc/infot, noutc, ok, lerr
2021 * .. Executable Statements ..
2022 * OK is set to .FALSE. by the special version of XERBLA or by CHKXER
2023 * if anything is wrong.
2024  ok = .true.
2025 * LERR is set to .TRUE. by the special version of XERBLA each time
2026 * it is called, and is then tested and re-set by CHKXER.
2027  lerr = .false.
2028 *
2029 * Initialize ALPHA, BETA, RALPHA, and RBETA.
2030 *
2031  alpha = cmplx( one, -one )
2032  beta = cmplx( two, -two )
2033  ralpha = one
2034  rbeta = two
2035 *
2036  go to( 10, 20, 30, 40, 50, 60, 70, 80,
2037  $ 90 )isnum
2038  10 infot = 1
2039  CALL cgemm( '/', 'N', 0, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2040  CALL chkxer( srnamt, infot, nout, lerr, ok )
2041  infot = 1
2042  CALL cgemm( '/', 'C', 0, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2043  CALL chkxer( srnamt, infot, nout, lerr, ok )
2044  infot = 1
2045  CALL cgemm( '/', 'T', 0, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2046  CALL chkxer( srnamt, infot, nout, lerr, ok )
2047  infot = 2
2048  CALL cgemm( 'N', '/', 0, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2049  CALL chkxer( srnamt, infot, nout, lerr, ok )
2050  infot = 2
2051  CALL cgemm( 'C', '/', 0, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2052  CALL chkxer( srnamt, infot, nout, lerr, ok )
2053  infot = 2
2054  CALL cgemm( 'T', '/', 0, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2055  CALL chkxer( srnamt, infot, nout, lerr, ok )
2056  infot = 3
2057  CALL cgemm( 'N', 'N', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2058  CALL chkxer( srnamt, infot, nout, lerr, ok )
2059  infot = 3
2060  CALL cgemm( 'N', 'C', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2061  CALL chkxer( srnamt, infot, nout, lerr, ok )
2062  infot = 3
2063  CALL cgemm( 'N', 'T', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2064  CALL chkxer( srnamt, infot, nout, lerr, ok )
2065  infot = 3
2066  CALL cgemm( 'C', 'N', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2067  CALL chkxer( srnamt, infot, nout, lerr, ok )
2068  infot = 3
2069  CALL cgemm( 'C', 'C', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2070  CALL chkxer( srnamt, infot, nout, lerr, ok )
2071  infot = 3
2072  CALL cgemm( 'C', 'T', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2073  CALL chkxer( srnamt, infot, nout, lerr, ok )
2074  infot = 3
2075  CALL cgemm( 'T', 'N', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2076  CALL chkxer( srnamt, infot, nout, lerr, ok )
2077  infot = 3
2078  CALL cgemm( 'T', 'C', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2079  CALL chkxer( srnamt, infot, nout, lerr, ok )
2080  infot = 3
2081  CALL cgemm( 'T', 'T', -1, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2082  CALL chkxer( srnamt, infot, nout, lerr, ok )
2083  infot = 4
2084  CALL cgemm( 'N', 'N', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2085  CALL chkxer( srnamt, infot, nout, lerr, ok )
2086  infot = 4
2087  CALL cgemm( 'N', 'C', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2088  CALL chkxer( srnamt, infot, nout, lerr, ok )
2089  infot = 4
2090  CALL cgemm( 'N', 'T', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2091  CALL chkxer( srnamt, infot, nout, lerr, ok )
2092  infot = 4
2093  CALL cgemm( 'C', 'N', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2094  CALL chkxer( srnamt, infot, nout, lerr, ok )
2095  infot = 4
2096  CALL cgemm( 'C', 'C', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2097  CALL chkxer( srnamt, infot, nout, lerr, ok )
2098  infot = 4
2099  CALL cgemm( 'C', 'T', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2100  CALL chkxer( srnamt, infot, nout, lerr, ok )
2101  infot = 4
2102  CALL cgemm( 'T', 'N', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2103  CALL chkxer( srnamt, infot, nout, lerr, ok )
2104  infot = 4
2105  CALL cgemm( 'T', 'C', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2106  CALL chkxer( srnamt, infot, nout, lerr, ok )
2107  infot = 4
2108  CALL cgemm( 'T', 'T', 0, -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2109  CALL chkxer( srnamt, infot, nout, lerr, ok )
2110  infot = 5
2111  CALL cgemm( 'N', 'N', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2112  CALL chkxer( srnamt, infot, nout, lerr, ok )
2113  infot = 5
2114  CALL cgemm( 'N', 'C', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2115  CALL chkxer( srnamt, infot, nout, lerr, ok )
2116  infot = 5
2117  CALL cgemm( 'N', 'T', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2118  CALL chkxer( srnamt, infot, nout, lerr, ok )
2119  infot = 5
2120  CALL cgemm( 'C', 'N', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2121  CALL chkxer( srnamt, infot, nout, lerr, ok )
2122  infot = 5
2123  CALL cgemm( 'C', 'C', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2124  CALL chkxer( srnamt, infot, nout, lerr, ok )
2125  infot = 5
2126  CALL cgemm( 'C', 'T', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2127  CALL chkxer( srnamt, infot, nout, lerr, ok )
2128  infot = 5
2129  CALL cgemm( 'T', 'N', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2130  CALL chkxer( srnamt, infot, nout, lerr, ok )
2131  infot = 5
2132  CALL cgemm( 'T', 'C', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2133  CALL chkxer( srnamt, infot, nout, lerr, ok )
2134  infot = 5
2135  CALL cgemm( 'T', 'T', 0, 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2136  CALL chkxer( srnamt, infot, nout, lerr, ok )
2137  infot = 8
2138  CALL cgemm( 'N', 'N', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 2 )
2139  CALL chkxer( srnamt, infot, nout, lerr, ok )
2140  infot = 8
2141  CALL cgemm( 'N', 'C', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 2 )
2142  CALL chkxer( srnamt, infot, nout, lerr, ok )
2143  infot = 8
2144  CALL cgemm( 'N', 'T', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 2 )
2145  CALL chkxer( srnamt, infot, nout, lerr, ok )
2146  infot = 8
2147  CALL cgemm( 'C', 'N', 0, 0, 2, alpha, a, 1, b, 2, beta, c, 1 )
2148  CALL chkxer( srnamt, infot, nout, lerr, ok )
2149  infot = 8
2150  CALL cgemm( 'C', 'C', 0, 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2151  CALL chkxer( srnamt, infot, nout, lerr, ok )
2152  infot = 8
2153  CALL cgemm( 'C', 'T', 0, 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2154  CALL chkxer( srnamt, infot, nout, lerr, ok )
2155  infot = 8
2156  CALL cgemm( 'T', 'N', 0, 0, 2, alpha, a, 1, b, 2, beta, c, 1 )
2157  CALL chkxer( srnamt, infot, nout, lerr, ok )
2158  infot = 8
2159  CALL cgemm( 'T', 'C', 0, 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2160  CALL chkxer( srnamt, infot, nout, lerr, ok )
2161  infot = 8
2162  CALL cgemm( 'T', 'T', 0, 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2163  CALL chkxer( srnamt, infot, nout, lerr, ok )
2164  infot = 10
2165  CALL cgemm( 'N', 'N', 0, 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2166  CALL chkxer( srnamt, infot, nout, lerr, ok )
2167  infot = 10
2168  CALL cgemm( 'C', 'N', 0, 0, 2, alpha, a, 2, b, 1, beta, c, 1 )
2169  CALL chkxer( srnamt, infot, nout, lerr, ok )
2170  infot = 10
2171  CALL cgemm( 'T', 'N', 0, 0, 2, alpha, a, 2, b, 1, beta, c, 1 )
2172  CALL chkxer( srnamt, infot, nout, lerr, ok )
2173  infot = 10
2174  CALL cgemm( 'N', 'C', 0, 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2175  CALL chkxer( srnamt, infot, nout, lerr, ok )
2176  infot = 10
2177  CALL cgemm( 'C', 'C', 0, 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2178  CALL chkxer( srnamt, infot, nout, lerr, ok )
2179  infot = 10
2180  CALL cgemm( 'T', 'C', 0, 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2181  CALL chkxer( srnamt, infot, nout, lerr, ok )
2182  infot = 10
2183  CALL cgemm( 'N', 'T', 0, 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2184  CALL chkxer( srnamt, infot, nout, lerr, ok )
2185  infot = 10
2186  CALL cgemm( 'C', 'T', 0, 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2187  CALL chkxer( srnamt, infot, nout, lerr, ok )
2188  infot = 10
2189  CALL cgemm( 'T', 'T', 0, 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2190  CALL chkxer( srnamt, infot, nout, lerr, ok )
2191  infot = 13
2192  CALL cgemm( 'N', 'N', 2, 0, 0, alpha, a, 2, b, 1, beta, c, 1 )
2193  CALL chkxer( srnamt, infot, nout, lerr, ok )
2194  infot = 13
2195  CALL cgemm( 'N', 'C', 2, 0, 0, alpha, a, 2, b, 1, beta, c, 1 )
2196  CALL chkxer( srnamt, infot, nout, lerr, ok )
2197  infot = 13
2198  CALL cgemm( 'N', 'T', 2, 0, 0, alpha, a, 2, b, 1, beta, c, 1 )
2199  CALL chkxer( srnamt, infot, nout, lerr, ok )
2200  infot = 13
2201  CALL cgemm( 'C', 'N', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2202  CALL chkxer( srnamt, infot, nout, lerr, ok )
2203  infot = 13
2204  CALL cgemm( 'C', 'C', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2205  CALL chkxer( srnamt, infot, nout, lerr, ok )
2206  infot = 13
2207  CALL cgemm( 'C', 'T', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2208  CALL chkxer( srnamt, infot, nout, lerr, ok )
2209  infot = 13
2210  CALL cgemm( 'T', 'N', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2211  CALL chkxer( srnamt, infot, nout, lerr, ok )
2212  infot = 13
2213  CALL cgemm( 'T', 'C', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2214  CALL chkxer( srnamt, infot, nout, lerr, ok )
2215  infot = 13
2216  CALL cgemm( 'T', 'T', 2, 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2217  CALL chkxer( srnamt, infot, nout, lerr, ok )
2218  go to 100
2219  20 infot = 1
2220  CALL chemm( '/', 'U', 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2221  CALL chkxer( srnamt, infot, nout, lerr, ok )
2222  infot = 2
2223  CALL chemm( 'L', '/', 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2224  CALL chkxer( srnamt, infot, nout, lerr, ok )
2225  infot = 3
2226  CALL chemm( 'L', 'U', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2227  CALL chkxer( srnamt, infot, nout, lerr, ok )
2228  infot = 3
2229  CALL chemm( 'R', 'U', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2230  CALL chkxer( srnamt, infot, nout, lerr, ok )
2231  infot = 3
2232  CALL chemm( 'L', 'L', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2233  CALL chkxer( srnamt, infot, nout, lerr, ok )
2234  infot = 3
2235  CALL chemm( 'R', 'L', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2236  CALL chkxer( srnamt, infot, nout, lerr, ok )
2237  infot = 4
2238  CALL chemm( 'L', 'U', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2239  CALL chkxer( srnamt, infot, nout, lerr, ok )
2240  infot = 4
2241  CALL chemm( 'R', 'U', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2242  CALL chkxer( srnamt, infot, nout, lerr, ok )
2243  infot = 4
2244  CALL chemm( 'L', 'L', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2245  CALL chkxer( srnamt, infot, nout, lerr, ok )
2246  infot = 4
2247  CALL chemm( 'R', 'L', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2248  CALL chkxer( srnamt, infot, nout, lerr, ok )
2249  infot = 7
2250  CALL chemm( 'L', 'U', 2, 0, alpha, a, 1, b, 2, beta, c, 2 )
2251  CALL chkxer( srnamt, infot, nout, lerr, ok )
2252  infot = 7
2253  CALL chemm( 'R', 'U', 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2254  CALL chkxer( srnamt, infot, nout, lerr, ok )
2255  infot = 7
2256  CALL chemm( 'L', 'L', 2, 0, alpha, a, 1, b, 2, beta, c, 2 )
2257  CALL chkxer( srnamt, infot, nout, lerr, ok )
2258  infot = 7
2259  CALL chemm( 'R', 'L', 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2260  CALL chkxer( srnamt, infot, nout, lerr, ok )
2261  infot = 9
2262  CALL chemm( 'L', 'U', 2, 0, alpha, a, 2, b, 1, beta, c, 2 )
2263  CALL chkxer( srnamt, infot, nout, lerr, ok )
2264  infot = 9
2265  CALL chemm( 'R', 'U', 2, 0, alpha, a, 1, b, 1, beta, c, 2 )
2266  CALL chkxer( srnamt, infot, nout, lerr, ok )
2267  infot = 9
2268  CALL chemm( 'L', 'L', 2, 0, alpha, a, 2, b, 1, beta, c, 2 )
2269  CALL chkxer( srnamt, infot, nout, lerr, ok )
2270  infot = 9
2271  CALL chemm( 'R', 'L', 2, 0, alpha, a, 1, b, 1, beta, c, 2 )
2272  CALL chkxer( srnamt, infot, nout, lerr, ok )
2273  infot = 12
2274  CALL chemm( 'L', 'U', 2, 0, alpha, a, 2, b, 2, beta, c, 1 )
2275  CALL chkxer( srnamt, infot, nout, lerr, ok )
2276  infot = 12
2277  CALL chemm( 'R', 'U', 2, 0, alpha, a, 1, b, 2, beta, c, 1 )
2278  CALL chkxer( srnamt, infot, nout, lerr, ok )
2279  infot = 12
2280  CALL chemm( 'L', 'L', 2, 0, alpha, a, 2, b, 2, beta, c, 1 )
2281  CALL chkxer( srnamt, infot, nout, lerr, ok )
2282  infot = 12
2283  CALL chemm( 'R', 'L', 2, 0, alpha, a, 1, b, 2, beta, c, 1 )
2284  CALL chkxer( srnamt, infot, nout, lerr, ok )
2285  go to 100
2286  30 infot = 1
2287  CALL csymm( '/', 'U', 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2288  CALL chkxer( srnamt, infot, nout, lerr, ok )
2289  infot = 2
2290  CALL csymm( 'L', '/', 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2291  CALL chkxer( srnamt, infot, nout, lerr, ok )
2292  infot = 3
2293  CALL csymm( 'L', 'U', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2294  CALL chkxer( srnamt, infot, nout, lerr, ok )
2295  infot = 3
2296  CALL csymm( 'R', 'U', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2297  CALL chkxer( srnamt, infot, nout, lerr, ok )
2298  infot = 3
2299  CALL csymm( 'L', 'L', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2300  CALL chkxer( srnamt, infot, nout, lerr, ok )
2301  infot = 3
2302  CALL csymm( 'R', 'L', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2303  CALL chkxer( srnamt, infot, nout, lerr, ok )
2304  infot = 4
2305  CALL csymm( 'L', 'U', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2306  CALL chkxer( srnamt, infot, nout, lerr, ok )
2307  infot = 4
2308  CALL csymm( 'R', 'U', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2309  CALL chkxer( srnamt, infot, nout, lerr, ok )
2310  infot = 4
2311  CALL csymm( 'L', 'L', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2312  CALL chkxer( srnamt, infot, nout, lerr, ok )
2313  infot = 4
2314  CALL csymm( 'R', 'L', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2315  CALL chkxer( srnamt, infot, nout, lerr, ok )
2316  infot = 7
2317  CALL csymm( 'L', 'U', 2, 0, alpha, a, 1, b, 2, beta, c, 2 )
2318  CALL chkxer( srnamt, infot, nout, lerr, ok )
2319  infot = 7
2320  CALL csymm( 'R', 'U', 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2321  CALL chkxer( srnamt, infot, nout, lerr, ok )
2322  infot = 7
2323  CALL csymm( 'L', 'L', 2, 0, alpha, a, 1, b, 2, beta, c, 2 )
2324  CALL chkxer( srnamt, infot, nout, lerr, ok )
2325  infot = 7
2326  CALL csymm( 'R', 'L', 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2327  CALL chkxer( srnamt, infot, nout, lerr, ok )
2328  infot = 9
2329  CALL csymm( 'L', 'U', 2, 0, alpha, a, 2, b, 1, beta, c, 2 )
2330  CALL chkxer( srnamt, infot, nout, lerr, ok )
2331  infot = 9
2332  CALL csymm( 'R', 'U', 2, 0, alpha, a, 1, b, 1, beta, c, 2 )
2333  CALL chkxer( srnamt, infot, nout, lerr, ok )
2334  infot = 9
2335  CALL csymm( 'L', 'L', 2, 0, alpha, a, 2, b, 1, beta, c, 2 )
2336  CALL chkxer( srnamt, infot, nout, lerr, ok )
2337  infot = 9
2338  CALL csymm( 'R', 'L', 2, 0, alpha, a, 1, b, 1, beta, c, 2 )
2339  CALL chkxer( srnamt, infot, nout, lerr, ok )
2340  infot = 12
2341  CALL csymm( 'L', 'U', 2, 0, alpha, a, 2, b, 2, beta, c, 1 )
2342  CALL chkxer( srnamt, infot, nout, lerr, ok )
2343  infot = 12
2344  CALL csymm( 'R', 'U', 2, 0, alpha, a, 1, b, 2, beta, c, 1 )
2345  CALL chkxer( srnamt, infot, nout, lerr, ok )
2346  infot = 12
2347  CALL csymm( 'L', 'L', 2, 0, alpha, a, 2, b, 2, beta, c, 1 )
2348  CALL chkxer( srnamt, infot, nout, lerr, ok )
2349  infot = 12
2350  CALL csymm( 'R', 'L', 2, 0, alpha, a, 1, b, 2, beta, c, 1 )
2351  CALL chkxer( srnamt, infot, nout, lerr, ok )
2352  go to 100
2353  40 infot = 1
2354  CALL ctrmm( '/', 'U', 'N', 'N', 0, 0, alpha, a, 1, b, 1 )
2355  CALL chkxer( srnamt, infot, nout, lerr, ok )
2356  infot = 2
2357  CALL ctrmm( 'L', '/', 'N', 'N', 0, 0, alpha, a, 1, b, 1 )
2358  CALL chkxer( srnamt, infot, nout, lerr, ok )
2359  infot = 3
2360  CALL ctrmm( 'L', 'U', '/', 'N', 0, 0, alpha, a, 1, b, 1 )
2361  CALL chkxer( srnamt, infot, nout, lerr, ok )
2362  infot = 4
2363  CALL ctrmm( 'L', 'U', 'N', '/', 0, 0, alpha, a, 1, b, 1 )
2364  CALL chkxer( srnamt, infot, nout, lerr, ok )
2365  infot = 5
2366  CALL ctrmm( 'L', 'U', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2367  CALL chkxer( srnamt, infot, nout, lerr, ok )
2368  infot = 5
2369  CALL ctrmm( 'L', 'U', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2370  CALL chkxer( srnamt, infot, nout, lerr, ok )
2371  infot = 5
2372  CALL ctrmm( 'L', 'U', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2373  CALL chkxer( srnamt, infot, nout, lerr, ok )
2374  infot = 5
2375  CALL ctrmm( 'R', 'U', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2376  CALL chkxer( srnamt, infot, nout, lerr, ok )
2377  infot = 5
2378  CALL ctrmm( 'R', 'U', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2379  CALL chkxer( srnamt, infot, nout, lerr, ok )
2380  infot = 5
2381  CALL ctrmm( 'R', 'U', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2382  CALL chkxer( srnamt, infot, nout, lerr, ok )
2383  infot = 5
2384  CALL ctrmm( 'L', 'L', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2385  CALL chkxer( srnamt, infot, nout, lerr, ok )
2386  infot = 5
2387  CALL ctrmm( 'L', 'L', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2388  CALL chkxer( srnamt, infot, nout, lerr, ok )
2389  infot = 5
2390  CALL ctrmm( 'L', 'L', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2391  CALL chkxer( srnamt, infot, nout, lerr, ok )
2392  infot = 5
2393  CALL ctrmm( 'R', 'L', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2394  CALL chkxer( srnamt, infot, nout, lerr, ok )
2395  infot = 5
2396  CALL ctrmm( 'R', 'L', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2397  CALL chkxer( srnamt, infot, nout, lerr, ok )
2398  infot = 5
2399  CALL ctrmm( 'R', 'L', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2400  CALL chkxer( srnamt, infot, nout, lerr, ok )
2401  infot = 6
2402  CALL ctrmm( 'L', 'U', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2403  CALL chkxer( srnamt, infot, nout, lerr, ok )
2404  infot = 6
2405  CALL ctrmm( 'L', 'U', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2406  CALL chkxer( srnamt, infot, nout, lerr, ok )
2407  infot = 6
2408  CALL ctrmm( 'L', 'U', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2409  CALL chkxer( srnamt, infot, nout, lerr, ok )
2410  infot = 6
2411  CALL ctrmm( 'R', 'U', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2412  CALL chkxer( srnamt, infot, nout, lerr, ok )
2413  infot = 6
2414  CALL ctrmm( 'R', 'U', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2415  CALL chkxer( srnamt, infot, nout, lerr, ok )
2416  infot = 6
2417  CALL ctrmm( 'R', 'U', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2418  CALL chkxer( srnamt, infot, nout, lerr, ok )
2419  infot = 6
2420  CALL ctrmm( 'L', 'L', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2421  CALL chkxer( srnamt, infot, nout, lerr, ok )
2422  infot = 6
2423  CALL ctrmm( 'L', 'L', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2424  CALL chkxer( srnamt, infot, nout, lerr, ok )
2425  infot = 6
2426  CALL ctrmm( 'L', 'L', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2427  CALL chkxer( srnamt, infot, nout, lerr, ok )
2428  infot = 6
2429  CALL ctrmm( 'R', 'L', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2430  CALL chkxer( srnamt, infot, nout, lerr, ok )
2431  infot = 6
2432  CALL ctrmm( 'R', 'L', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2433  CALL chkxer( srnamt, infot, nout, lerr, ok )
2434  infot = 6
2435  CALL ctrmm( 'R', 'L', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2436  CALL chkxer( srnamt, infot, nout, lerr, ok )
2437  infot = 9
2438  CALL ctrmm( 'L', 'U', 'N', 'N', 2, 0, alpha, a, 1, b, 2 )
2439  CALL chkxer( srnamt, infot, nout, lerr, ok )
2440  infot = 9
2441  CALL ctrmm( 'L', 'U', 'C', 'N', 2, 0, alpha, a, 1, b, 2 )
2442  CALL chkxer( srnamt, infot, nout, lerr, ok )
2443  infot = 9
2444  CALL ctrmm( 'L', 'U', 'T', 'N', 2, 0, alpha, a, 1, b, 2 )
2445  CALL chkxer( srnamt, infot, nout, lerr, ok )
2446  infot = 9
2447  CALL ctrmm( 'R', 'U', 'N', 'N', 0, 2, alpha, a, 1, b, 1 )
2448  CALL chkxer( srnamt, infot, nout, lerr, ok )
2449  infot = 9
2450  CALL ctrmm( 'R', 'U', 'C', 'N', 0, 2, alpha, a, 1, b, 1 )
2451  CALL chkxer( srnamt, infot, nout, lerr, ok )
2452  infot = 9
2453  CALL ctrmm( 'R', 'U', 'T', 'N', 0, 2, alpha, a, 1, b, 1 )
2454  CALL chkxer( srnamt, infot, nout, lerr, ok )
2455  infot = 9
2456  CALL ctrmm( 'L', 'L', 'N', 'N', 2, 0, alpha, a, 1, b, 2 )
2457  CALL chkxer( srnamt, infot, nout, lerr, ok )
2458  infot = 9
2459  CALL ctrmm( 'L', 'L', 'C', 'N', 2, 0, alpha, a, 1, b, 2 )
2460  CALL chkxer( srnamt, infot, nout, lerr, ok )
2461  infot = 9
2462  CALL ctrmm( 'L', 'L', 'T', 'N', 2, 0, alpha, a, 1, b, 2 )
2463  CALL chkxer( srnamt, infot, nout, lerr, ok )
2464  infot = 9
2465  CALL ctrmm( 'R', 'L', 'N', 'N', 0, 2, alpha, a, 1, b, 1 )
2466  CALL chkxer( srnamt, infot, nout, lerr, ok )
2467  infot = 9
2468  CALL ctrmm( 'R', 'L', 'C', 'N', 0, 2, alpha, a, 1, b, 1 )
2469  CALL chkxer( srnamt, infot, nout, lerr, ok )
2470  infot = 9
2471  CALL ctrmm( 'R', 'L', 'T', 'N', 0, 2, alpha, a, 1, b, 1 )
2472  CALL chkxer( srnamt, infot, nout, lerr, ok )
2473  infot = 11
2474  CALL ctrmm( 'L', 'U', 'N', 'N', 2, 0, alpha, a, 2, b, 1 )
2475  CALL chkxer( srnamt, infot, nout, lerr, ok )
2476  infot = 11
2477  CALL ctrmm( 'L', 'U', 'C', 'N', 2, 0, alpha, a, 2, b, 1 )
2478  CALL chkxer( srnamt, infot, nout, lerr, ok )
2479  infot = 11
2480  CALL ctrmm( 'L', 'U', 'T', 'N', 2, 0, alpha, a, 2, b, 1 )
2481  CALL chkxer( srnamt, infot, nout, lerr, ok )
2482  infot = 11
2483  CALL ctrmm( 'R', 'U', 'N', 'N', 2, 0, alpha, a, 1, b, 1 )
2484  CALL chkxer( srnamt, infot, nout, lerr, ok )
2485  infot = 11
2486  CALL ctrmm( 'R', 'U', 'C', 'N', 2, 0, alpha, a, 1, b, 1 )
2487  CALL chkxer( srnamt, infot, nout, lerr, ok )
2488  infot = 11
2489  CALL ctrmm( 'R', 'U', 'T', 'N', 2, 0, alpha, a, 1, b, 1 )
2490  CALL chkxer( srnamt, infot, nout, lerr, ok )
2491  infot = 11
2492  CALL ctrmm( 'L', 'L', 'N', 'N', 2, 0, alpha, a, 2, b, 1 )
2493  CALL chkxer( srnamt, infot, nout, lerr, ok )
2494  infot = 11
2495  CALL ctrmm( 'L', 'L', 'C', 'N', 2, 0, alpha, a, 2, b, 1 )
2496  CALL chkxer( srnamt, infot, nout, lerr, ok )
2497  infot = 11
2498  CALL ctrmm( 'L', 'L', 'T', 'N', 2, 0, alpha, a, 2, b, 1 )
2499  CALL chkxer( srnamt, infot, nout, lerr, ok )
2500  infot = 11
2501  CALL ctrmm( 'R', 'L', 'N', 'N', 2, 0, alpha, a, 1, b, 1 )
2502  CALL chkxer( srnamt, infot, nout, lerr, ok )
2503  infot = 11
2504  CALL ctrmm( 'R', 'L', 'C', 'N', 2, 0, alpha, a, 1, b, 1 )
2505  CALL chkxer( srnamt, infot, nout, lerr, ok )
2506  infot = 11
2507  CALL ctrmm( 'R', 'L', 'T', 'N', 2, 0, alpha, a, 1, b, 1 )
2508  CALL chkxer( srnamt, infot, nout, lerr, ok )
2509  go to 100
2510  50 infot = 1
2511  CALL ctrsm( '/', 'U', 'N', 'N', 0, 0, alpha, a, 1, b, 1 )
2512  CALL chkxer( srnamt, infot, nout, lerr, ok )
2513  infot = 2
2514  CALL ctrsm( 'L', '/', 'N', 'N', 0, 0, alpha, a, 1, b, 1 )
2515  CALL chkxer( srnamt, infot, nout, lerr, ok )
2516  infot = 3
2517  CALL ctrsm( 'L', 'U', '/', 'N', 0, 0, alpha, a, 1, b, 1 )
2518  CALL chkxer( srnamt, infot, nout, lerr, ok )
2519  infot = 4
2520  CALL ctrsm( 'L', 'U', 'N', '/', 0, 0, alpha, a, 1, b, 1 )
2521  CALL chkxer( srnamt, infot, nout, lerr, ok )
2522  infot = 5
2523  CALL ctrsm( 'L', 'U', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2524  CALL chkxer( srnamt, infot, nout, lerr, ok )
2525  infot = 5
2526  CALL ctrsm( 'L', 'U', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2527  CALL chkxer( srnamt, infot, nout, lerr, ok )
2528  infot = 5
2529  CALL ctrsm( 'L', 'U', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2530  CALL chkxer( srnamt, infot, nout, lerr, ok )
2531  infot = 5
2532  CALL ctrsm( 'R', 'U', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2533  CALL chkxer( srnamt, infot, nout, lerr, ok )
2534  infot = 5
2535  CALL ctrsm( 'R', 'U', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2536  CALL chkxer( srnamt, infot, nout, lerr, ok )
2537  infot = 5
2538  CALL ctrsm( 'R', 'U', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2539  CALL chkxer( srnamt, infot, nout, lerr, ok )
2540  infot = 5
2541  CALL ctrsm( 'L', 'L', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2542  CALL chkxer( srnamt, infot, nout, lerr, ok )
2543  infot = 5
2544  CALL ctrsm( 'L', 'L', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2545  CALL chkxer( srnamt, infot, nout, lerr, ok )
2546  infot = 5
2547  CALL ctrsm( 'L', 'L', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2548  CALL chkxer( srnamt, infot, nout, lerr, ok )
2549  infot = 5
2550  CALL ctrsm( 'R', 'L', 'N', 'N', -1, 0, alpha, a, 1, b, 1 )
2551  CALL chkxer( srnamt, infot, nout, lerr, ok )
2552  infot = 5
2553  CALL ctrsm( 'R', 'L', 'C', 'N', -1, 0, alpha, a, 1, b, 1 )
2554  CALL chkxer( srnamt, infot, nout, lerr, ok )
2555  infot = 5
2556  CALL ctrsm( 'R', 'L', 'T', 'N', -1, 0, alpha, a, 1, b, 1 )
2557  CALL chkxer( srnamt, infot, nout, lerr, ok )
2558  infot = 6
2559  CALL ctrsm( 'L', 'U', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2560  CALL chkxer( srnamt, infot, nout, lerr, ok )
2561  infot = 6
2562  CALL ctrsm( 'L', 'U', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2563  CALL chkxer( srnamt, infot, nout, lerr, ok )
2564  infot = 6
2565  CALL ctrsm( 'L', 'U', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2566  CALL chkxer( srnamt, infot, nout, lerr, ok )
2567  infot = 6
2568  CALL ctrsm( 'R', 'U', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2569  CALL chkxer( srnamt, infot, nout, lerr, ok )
2570  infot = 6
2571  CALL ctrsm( 'R', 'U', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2572  CALL chkxer( srnamt, infot, nout, lerr, ok )
2573  infot = 6
2574  CALL ctrsm( 'R', 'U', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2575  CALL chkxer( srnamt, infot, nout, lerr, ok )
2576  infot = 6
2577  CALL ctrsm( 'L', 'L', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2578  CALL chkxer( srnamt, infot, nout, lerr, ok )
2579  infot = 6
2580  CALL ctrsm( 'L', 'L', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2581  CALL chkxer( srnamt, infot, nout, lerr, ok )
2582  infot = 6
2583  CALL ctrsm( 'L', 'L', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2584  CALL chkxer( srnamt, infot, nout, lerr, ok )
2585  infot = 6
2586  CALL ctrsm( 'R', 'L', 'N', 'N', 0, -1, alpha, a, 1, b, 1 )
2587  CALL chkxer( srnamt, infot, nout, lerr, ok )
2588  infot = 6
2589  CALL ctrsm( 'R', 'L', 'C', 'N', 0, -1, alpha, a, 1, b, 1 )
2590  CALL chkxer( srnamt, infot, nout, lerr, ok )
2591  infot = 6
2592  CALL ctrsm( 'R', 'L', 'T', 'N', 0, -1, alpha, a, 1, b, 1 )
2593  CALL chkxer( srnamt, infot, nout, lerr, ok )
2594  infot = 9
2595  CALL ctrsm( 'L', 'U', 'N', 'N', 2, 0, alpha, a, 1, b, 2 )
2596  CALL chkxer( srnamt, infot, nout, lerr, ok )
2597  infot = 9
2598  CALL ctrsm( 'L', 'U', 'C', 'N', 2, 0, alpha, a, 1, b, 2 )
2599  CALL chkxer( srnamt, infot, nout, lerr, ok )
2600  infot = 9
2601  CALL ctrsm( 'L', 'U', 'T', 'N', 2, 0, alpha, a, 1, b, 2 )
2602  CALL chkxer( srnamt, infot, nout, lerr, ok )
2603  infot = 9
2604  CALL ctrsm( 'R', 'U', 'N', 'N', 0, 2, alpha, a, 1, b, 1 )
2605  CALL chkxer( srnamt, infot, nout, lerr, ok )
2606  infot = 9
2607  CALL ctrsm( 'R', 'U', 'C', 'N', 0, 2, alpha, a, 1, b, 1 )
2608  CALL chkxer( srnamt, infot, nout, lerr, ok )
2609  infot = 9
2610  CALL ctrsm( 'R', 'U', 'T', 'N', 0, 2, alpha, a, 1, b, 1 )
2611  CALL chkxer( srnamt, infot, nout, lerr, ok )
2612  infot = 9
2613  CALL ctrsm( 'L', 'L', 'N', 'N', 2, 0, alpha, a, 1, b, 2 )
2614  CALL chkxer( srnamt, infot, nout, lerr, ok )
2615  infot = 9
2616  CALL ctrsm( 'L', 'L', 'C', 'N', 2, 0, alpha, a, 1, b, 2 )
2617  CALL chkxer( srnamt, infot, nout, lerr, ok )
2618  infot = 9
2619  CALL ctrsm( 'L', 'L', 'T', 'N', 2, 0, alpha, a, 1, b, 2 )
2620  CALL chkxer( srnamt, infot, nout, lerr, ok )
2621  infot = 9
2622  CALL ctrsm( 'R', 'L', 'N', 'N', 0, 2, alpha, a, 1, b, 1 )
2623  CALL chkxer( srnamt, infot, nout, lerr, ok )
2624  infot = 9
2625  CALL ctrsm( 'R', 'L', 'C', 'N', 0, 2, alpha, a, 1, b, 1 )
2626  CALL chkxer( srnamt, infot, nout, lerr, ok )
2627  infot = 9
2628  CALL ctrsm( 'R', 'L', 'T', 'N', 0, 2, alpha, a, 1, b, 1 )
2629  CALL chkxer( srnamt, infot, nout, lerr, ok )
2630  infot = 11
2631  CALL ctrsm( 'L', 'U', 'N', 'N', 2, 0, alpha, a, 2, b, 1 )
2632  CALL chkxer( srnamt, infot, nout, lerr, ok )
2633  infot = 11
2634  CALL ctrsm( 'L', 'U', 'C', 'N', 2, 0, alpha, a, 2, b, 1 )
2635  CALL chkxer( srnamt, infot, nout, lerr, ok )
2636  infot = 11
2637  CALL ctrsm( 'L', 'U', 'T', 'N', 2, 0, alpha, a, 2, b, 1 )
2638  CALL chkxer( srnamt, infot, nout, lerr, ok )
2639  infot = 11
2640  CALL ctrsm( 'R', 'U', 'N', 'N', 2, 0, alpha, a, 1, b, 1 )
2641  CALL chkxer( srnamt, infot, nout, lerr, ok )
2642  infot = 11
2643  CALL ctrsm( 'R', 'U', 'C', 'N', 2, 0, alpha, a, 1, b, 1 )
2644  CALL chkxer( srnamt, infot, nout, lerr, ok )
2645  infot = 11
2646  CALL ctrsm( 'R', 'U', 'T', 'N', 2, 0, alpha, a, 1, b, 1 )
2647  CALL chkxer( srnamt, infot, nout, lerr, ok )
2648  infot = 11
2649  CALL ctrsm( 'L', 'L', 'N', 'N', 2, 0, alpha, a, 2, b, 1 )
2650  CALL chkxer( srnamt, infot, nout, lerr, ok )
2651  infot = 11
2652  CALL ctrsm( 'L', 'L', 'C', 'N', 2, 0, alpha, a, 2, b, 1 )
2653  CALL chkxer( srnamt, infot, nout, lerr, ok )
2654  infot = 11
2655  CALL ctrsm( 'L', 'L', 'T', 'N', 2, 0, alpha, a, 2, b, 1 )
2656  CALL chkxer( srnamt, infot, nout, lerr, ok )
2657  infot = 11
2658  CALL ctrsm( 'R', 'L', 'N', 'N', 2, 0, alpha, a, 1, b, 1 )
2659  CALL chkxer( srnamt, infot, nout, lerr, ok )
2660  infot = 11
2661  CALL ctrsm( 'R', 'L', 'C', 'N', 2, 0, alpha, a, 1, b, 1 )
2662  CALL chkxer( srnamt, infot, nout, lerr, ok )
2663  infot = 11
2664  CALL ctrsm( 'R', 'L', 'T', 'N', 2, 0, alpha, a, 1, b, 1 )
2665  CALL chkxer( srnamt, infot, nout, lerr, ok )
2666  go to 100
2667  60 infot = 1
2668  CALL cherk( '/', 'N', 0, 0, ralpha, a, 1, rbeta, c, 1 )
2669  CALL chkxer( srnamt, infot, nout, lerr, ok )
2670  infot = 2
2671  CALL cherk( 'U', 'T', 0, 0, ralpha, a, 1, rbeta, c, 1 )
2672  CALL chkxer( srnamt, infot, nout, lerr, ok )
2673  infot = 3
2674  CALL cherk( 'U', 'N', -1, 0, ralpha, a, 1, rbeta, c, 1 )
2675  CALL chkxer( srnamt, infot, nout, lerr, ok )
2676  infot = 3
2677  CALL cherk( 'U', 'C', -1, 0, ralpha, a, 1, rbeta, c, 1 )
2678  CALL chkxer( srnamt, infot, nout, lerr, ok )
2679  infot = 3
2680  CALL cherk( 'L', 'N', -1, 0, ralpha, a, 1, rbeta, c, 1 )
2681  CALL chkxer( srnamt, infot, nout, lerr, ok )
2682  infot = 3
2683  CALL cherk( 'L', 'C', -1, 0, ralpha, a, 1, rbeta, c, 1 )
2684  CALL chkxer( srnamt, infot, nout, lerr, ok )
2685  infot = 4
2686  CALL cherk( 'U', 'N', 0, -1, ralpha, a, 1, rbeta, c, 1 )
2687  CALL chkxer( srnamt, infot, nout, lerr, ok )
2688  infot = 4
2689  CALL cherk( 'U', 'C', 0, -1, ralpha, a, 1, rbeta, c, 1 )
2690  CALL chkxer( srnamt, infot, nout, lerr, ok )
2691  infot = 4
2692  CALL cherk( 'L', 'N', 0, -1, ralpha, a, 1, rbeta, c, 1 )
2693  CALL chkxer( srnamt, infot, nout, lerr, ok )
2694  infot = 4
2695  CALL cherk( 'L', 'C', 0, -1, ralpha, a, 1, rbeta, c, 1 )
2696  CALL chkxer( srnamt, infot, nout, lerr, ok )
2697  infot = 7
2698  CALL cherk( 'U', 'N', 2, 0, ralpha, a, 1, rbeta, c, 2 )
2699  CALL chkxer( srnamt, infot, nout, lerr, ok )
2700  infot = 7
2701  CALL cherk( 'U', 'C', 0, 2, ralpha, a, 1, rbeta, c, 1 )
2702  CALL chkxer( srnamt, infot, nout, lerr, ok )
2703  infot = 7
2704  CALL cherk( 'L', 'N', 2, 0, ralpha, a, 1, rbeta, c, 2 )
2705  CALL chkxer( srnamt, infot, nout, lerr, ok )
2706  infot = 7
2707  CALL cherk( 'L', 'C', 0, 2, ralpha, a, 1, rbeta, c, 1 )
2708  CALL chkxer( srnamt, infot, nout, lerr, ok )
2709  infot = 10
2710  CALL cherk( 'U', 'N', 2, 0, ralpha, a, 2, rbeta, c, 1 )
2711  CALL chkxer( srnamt, infot, nout, lerr, ok )
2712  infot = 10
2713  CALL cherk( 'U', 'C', 2, 0, ralpha, a, 1, rbeta, c, 1 )
2714  CALL chkxer( srnamt, infot, nout, lerr, ok )
2715  infot = 10
2716  CALL cherk( 'L', 'N', 2, 0, ralpha, a, 2, rbeta, c, 1 )
2717  CALL chkxer( srnamt, infot, nout, lerr, ok )
2718  infot = 10
2719  CALL cherk( 'L', 'C', 2, 0, ralpha, a, 1, rbeta, c, 1 )
2720  CALL chkxer( srnamt, infot, nout, lerr, ok )
2721  go to 100
2722  70 infot = 1
2723  CALL csyrk( '/', 'N', 0, 0, alpha, a, 1, beta, c, 1 )
2724  CALL chkxer( srnamt, infot, nout, lerr, ok )
2725  infot = 2
2726  CALL csyrk( 'U', 'C', 0, 0, alpha, a, 1, beta, c, 1 )
2727  CALL chkxer( srnamt, infot, nout, lerr, ok )
2728  infot = 3
2729  CALL csyrk( 'U', 'N', -1, 0, alpha, a, 1, beta, c, 1 )
2730  CALL chkxer( srnamt, infot, nout, lerr, ok )
2731  infot = 3
2732  CALL csyrk( 'U', 'T', -1, 0, alpha, a, 1, beta, c, 1 )
2733  CALL chkxer( srnamt, infot, nout, lerr, ok )
2734  infot = 3
2735  CALL csyrk( 'L', 'N', -1, 0, alpha, a, 1, beta, c, 1 )
2736  CALL chkxer( srnamt, infot, nout, lerr, ok )
2737  infot = 3
2738  CALL csyrk( 'L', 'T', -1, 0, alpha, a, 1, beta, c, 1 )
2739  CALL chkxer( srnamt, infot, nout, lerr, ok )
2740  infot = 4
2741  CALL csyrk( 'U', 'N', 0, -1, alpha, a, 1, beta, c, 1 )
2742  CALL chkxer( srnamt, infot, nout, lerr, ok )
2743  infot = 4
2744  CALL csyrk( 'U', 'T', 0, -1, alpha, a, 1, beta, c, 1 )
2745  CALL chkxer( srnamt, infot, nout, lerr, ok )
2746  infot = 4
2747  CALL csyrk( 'L', 'N', 0, -1, alpha, a, 1, beta, c, 1 )
2748  CALL chkxer( srnamt, infot, nout, lerr, ok )
2749  infot = 4
2750  CALL csyrk( 'L', 'T', 0, -1, alpha, a, 1, beta, c, 1 )
2751  CALL chkxer( srnamt, infot, nout, lerr, ok )
2752  infot = 7
2753  CALL csyrk( 'U', 'N', 2, 0, alpha, a, 1, beta, c, 2 )
2754  CALL chkxer( srnamt, infot, nout, lerr, ok )
2755  infot = 7
2756  CALL csyrk( 'U', 'T', 0, 2, alpha, a, 1, beta, c, 1 )
2757  CALL chkxer( srnamt, infot, nout, lerr, ok )
2758  infot = 7
2759  CALL csyrk( 'L', 'N', 2, 0, alpha, a, 1, beta, c, 2 )
2760  CALL chkxer( srnamt, infot, nout, lerr, ok )
2761  infot = 7
2762  CALL csyrk( 'L', 'T', 0, 2, alpha, a, 1, beta, c, 1 )
2763  CALL chkxer( srnamt, infot, nout, lerr, ok )
2764  infot = 10
2765  CALL csyrk( 'U', 'N', 2, 0, alpha, a, 2, beta, c, 1 )
2766  CALL chkxer( srnamt, infot, nout, lerr, ok )
2767  infot = 10
2768  CALL csyrk( 'U', 'T', 2, 0, alpha, a, 1, beta, c, 1 )
2769  CALL chkxer( srnamt, infot, nout, lerr, ok )
2770  infot = 10
2771  CALL csyrk( 'L', 'N', 2, 0, alpha, a, 2, beta, c, 1 )
2772  CALL chkxer( srnamt, infot, nout, lerr, ok )
2773  infot = 10
2774  CALL csyrk( 'L', 'T', 2, 0, alpha, a, 1, beta, c, 1 )
2775  CALL chkxer( srnamt, infot, nout, lerr, ok )
2776  go to 100
2777  80 infot = 1
2778  CALL cher2k( '/', 'N', 0, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2779  CALL chkxer( srnamt, infot, nout, lerr, ok )
2780  infot = 2
2781  CALL cher2k( 'U', 'T', 0, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2782  CALL chkxer( srnamt, infot, nout, lerr, ok )
2783  infot = 3
2784  CALL cher2k( 'U', 'N', -1, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2785  CALL chkxer( srnamt, infot, nout, lerr, ok )
2786  infot = 3
2787  CALL cher2k( 'U', 'C', -1, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2788  CALL chkxer( srnamt, infot, nout, lerr, ok )
2789  infot = 3
2790  CALL cher2k( 'L', 'N', -1, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2791  CALL chkxer( srnamt, infot, nout, lerr, ok )
2792  infot = 3
2793  CALL cher2k( 'L', 'C', -1, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2794  CALL chkxer( srnamt, infot, nout, lerr, ok )
2795  infot = 4
2796  CALL cher2k( 'U', 'N', 0, -1, alpha, a, 1, b, 1, rbeta, c, 1 )
2797  CALL chkxer( srnamt, infot, nout, lerr, ok )
2798  infot = 4
2799  CALL cher2k( 'U', 'C', 0, -1, alpha, a, 1, b, 1, rbeta, c, 1 )
2800  CALL chkxer( srnamt, infot, nout, lerr, ok )
2801  infot = 4
2802  CALL cher2k( 'L', 'N', 0, -1, alpha, a, 1, b, 1, rbeta, c, 1 )
2803  CALL chkxer( srnamt, infot, nout, lerr, ok )
2804  infot = 4
2805  CALL cher2k( 'L', 'C', 0, -1, alpha, a, 1, b, 1, rbeta, c, 1 )
2806  CALL chkxer( srnamt, infot, nout, lerr, ok )
2807  infot = 7
2808  CALL cher2k( 'U', 'N', 2, 0, alpha, a, 1, b, 1, rbeta, c, 2 )
2809  CALL chkxer( srnamt, infot, nout, lerr, ok )
2810  infot = 7
2811  CALL cher2k( 'U', 'C', 0, 2, alpha, a, 1, b, 1, rbeta, c, 1 )
2812  CALL chkxer( srnamt, infot, nout, lerr, ok )
2813  infot = 7
2814  CALL cher2k( 'L', 'N', 2, 0, alpha, a, 1, b, 1, rbeta, c, 2 )
2815  CALL chkxer( srnamt, infot, nout, lerr, ok )
2816  infot = 7
2817  CALL cher2k( 'L', 'C', 0, 2, alpha, a, 1, b, 1, rbeta, c, 1 )
2818  CALL chkxer( srnamt, infot, nout, lerr, ok )
2819  infot = 9
2820  CALL cher2k( 'U', 'N', 2, 0, alpha, a, 2, b, 1, rbeta, c, 2 )
2821  CALL chkxer( srnamt, infot, nout, lerr, ok )
2822  infot = 9
2823  CALL cher2k( 'U', 'C', 0, 2, alpha, a, 2, b, 1, rbeta, c, 1 )
2824  CALL chkxer( srnamt, infot, nout, lerr, ok )
2825  infot = 9
2826  CALL cher2k( 'L', 'N', 2, 0, alpha, a, 2, b, 1, rbeta, c, 2 )
2827  CALL chkxer( srnamt, infot, nout, lerr, ok )
2828  infot = 9
2829  CALL cher2k( 'L', 'C', 0, 2, alpha, a, 2, b, 1, rbeta, c, 1 )
2830  CALL chkxer( srnamt, infot, nout, lerr, ok )
2831  infot = 12
2832  CALL cher2k( 'U', 'N', 2, 0, alpha, a, 2, b, 2, rbeta, c, 1 )
2833  CALL chkxer( srnamt, infot, nout, lerr, ok )
2834  infot = 12
2835  CALL cher2k( 'U', 'C', 2, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2836  CALL chkxer( srnamt, infot, nout, lerr, ok )
2837  infot = 12
2838  CALL cher2k( 'L', 'N', 2, 0, alpha, a, 2, b, 2, rbeta, c, 1 )
2839  CALL chkxer( srnamt, infot, nout, lerr, ok )
2840  infot = 12
2841  CALL cher2k( 'L', 'C', 2, 0, alpha, a, 1, b, 1, rbeta, c, 1 )
2842  CALL chkxer( srnamt, infot, nout, lerr, ok )
2843  go to 100
2844  90 infot = 1
2845  CALL csyr2k( '/', 'N', 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2846  CALL chkxer( srnamt, infot, nout, lerr, ok )
2847  infot = 2
2848  CALL csyr2k( 'U', 'C', 0, 0, alpha, a, 1, b, 1, beta, c, 1 )
2849  CALL chkxer( srnamt, infot, nout, lerr, ok )
2850  infot = 3
2851  CALL csyr2k( 'U', 'N', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2852  CALL chkxer( srnamt, infot, nout, lerr, ok )
2853  infot = 3
2854  CALL csyr2k( 'U', 'T', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2855  CALL chkxer( srnamt, infot, nout, lerr, ok )
2856  infot = 3
2857  CALL csyr2k( 'L', 'N', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2858  CALL chkxer( srnamt, infot, nout, lerr, ok )
2859  infot = 3
2860  CALL csyr2k( 'L', 'T', -1, 0, alpha, a, 1, b, 1, beta, c, 1 )
2861  CALL chkxer( srnamt, infot, nout, lerr, ok )
2862  infot = 4
2863  CALL csyr2k( 'U', 'N', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2864  CALL chkxer( srnamt, infot, nout, lerr, ok )
2865  infot = 4
2866  CALL csyr2k( 'U', 'T', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2867  CALL chkxer( srnamt, infot, nout, lerr, ok )
2868  infot = 4
2869  CALL csyr2k( 'L', 'N', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2870  CALL chkxer( srnamt, infot, nout, lerr, ok )
2871  infot = 4
2872  CALL csyr2k( 'L', 'T', 0, -1, alpha, a, 1, b, 1, beta, c, 1 )
2873  CALL chkxer( srnamt, infot, nout, lerr, ok )
2874  infot = 7
2875  CALL csyr2k( 'U', 'N', 2, 0, alpha, a, 1, b, 1, beta, c, 2 )
2876  CALL chkxer( srnamt, infot, nout, lerr, ok )
2877  infot = 7
2878  CALL csyr2k( 'U', 'T', 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2879  CALL chkxer( srnamt, infot, nout, lerr, ok )
2880  infot = 7
2881  CALL csyr2k( 'L', 'N', 2, 0, alpha, a, 1, b, 1, beta, c, 2 )
2882  CALL chkxer( srnamt, infot, nout, lerr, ok )
2883  infot = 7
2884  CALL csyr2k( 'L', 'T', 0, 2, alpha, a, 1, b, 1, beta, c, 1 )
2885  CALL chkxer( srnamt, infot, nout, lerr, ok )
2886  infot = 9
2887  CALL csyr2k( 'U', 'N', 2, 0, alpha, a, 2, b, 1, beta, c, 2 )
2888  CALL chkxer( srnamt, infot, nout, lerr, ok )
2889  infot = 9
2890  CALL csyr2k( 'U', 'T', 0, 2, alpha, a, 2, b, 1, beta, c, 1 )
2891  CALL chkxer( srnamt, infot, nout, lerr, ok )
2892  infot = 9
2893  CALL csyr2k( 'L', 'N', 2, 0, alpha, a, 2, b, 1, beta, c, 2 )
2894  CALL chkxer( srnamt, infot, nout, lerr, ok )
2895  infot = 9
2896  CALL csyr2k( 'L', 'T', 0, 2, alpha, a, 2, b, 1, beta, c, 1 )
2897  CALL chkxer( srnamt, infot, nout, lerr, ok )
2898  infot = 12
2899  CALL csyr2k( 'U', 'N', 2, 0, alpha, a, 2, b, 2, beta, c, 1 )
2900  CALL chkxer( srnamt, infot, nout, lerr, ok )
2901  infot = 12
2902  CALL csyr2k( 'U', 'T', 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2903  CALL chkxer( srnamt, infot, nout, lerr, ok )
2904  infot = 12
2905  CALL csyr2k( 'L', 'N', 2, 0, alpha, a, 2, b, 2, beta, c, 1 )
2906  CALL chkxer( srnamt, infot, nout, lerr, ok )
2907  infot = 12
2908  CALL csyr2k( 'L', 'T', 2, 0, alpha, a, 1, b, 1, beta, c, 1 )
2909  CALL chkxer( srnamt, infot, nout, lerr, ok )
2910 *
2911  100 IF( ok )THEN
2912  WRITE( nout, fmt = 9999 )srnamt
2913  ELSE
2914  WRITE( nout, fmt = 9998 )srnamt
2915  END IF
2916  return
2917 *
2918  9999 format( ' ', a6, ' PASSED THE TESTS OF ERROR-EXITS' )
2919  9998 format( ' ******* ', a6, ' FAILED THE TESTS OF ERROR-EXITS *****',
2920  $ '**' )
2921 *
2922 * End of CCHKE.
2923 *
2924  END
2925  SUBROUTINE cmake( TYPE, UPLO, DIAG, M, N, A, NMAX, AA, LDA, RESET,
2926  $ transl )
2927 *
2928 * Generates values for an M by N matrix A.
2929 * Stores the values in the array AA in the data structure required
2930 * by the routine, with unwanted elements set to rogue value.
2931 *
2932 * TYPE is 'GE', 'HE', 'SY' or 'TR'.
2933 *
2934 * Auxiliary routine for test program for Level 3 Blas.
2935 *
2936 * -- Written on 8-February-1989.
2937 * Jack Dongarra, Argonne National Laboratory.
2938 * Iain Duff, AERE Harwell.
2939 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
2940 * Sven Hammarling, Numerical Algorithms Group Ltd.
2941 *
2942 * .. Parameters ..
2943  COMPLEX zero, one
2944  parameter( zero = ( 0.0, 0.0 ), one = ( 1.0, 0.0 ) )
2945  COMPLEX rogue
2946  parameter( rogue = ( -1.0e10, 1.0e10 ) )
2947  REAL rzero
2948  parameter( rzero = 0.0 )
2949  REAL rrogue
2950  parameter( rrogue = -1.0e10 )
2951 * .. Scalar Arguments ..
2952  COMPLEX transl
2953  INTEGER lda, m, n, nmax
2954  LOGICAL reset
2955  CHARACTER*1 diag, uplo
2956  CHARACTER*2 type
2957 * .. Array Arguments ..
2958  COMPLEX a( nmax, * ), aa( * )
2959 * .. Local Scalars ..
2960  INTEGER i, ibeg, iend, j, jj
2961  LOGICAL gen, her, lower, sym, tri, unit, upper
2962 * .. External Functions ..
2963  COMPLEX cbeg
2964  EXTERNAL cbeg
2965 * .. Intrinsic Functions ..
2966  INTRINSIC cmplx, conjg, real
2967 * .. Executable Statements ..
2968  gen = type.EQ.'GE'
2969  her = type.EQ.'HE'
2970  sym = type.EQ.'SY'
2971  tri = type.EQ.'TR'
2972  upper = ( her.OR.sym.OR.tri ).AND.uplo.EQ.'U'
2973  lower = ( her.OR.sym.OR.tri ).AND.uplo.EQ.'L'
2974  unit = tri.AND.diag.EQ.'U'
2975 *
2976 * Generate data in array A.
2977 *
2978  DO 20 j = 1, n
2979  DO 10 i = 1, m
2980  IF( gen.OR.( upper.AND.i.LE.j ).OR.( lower.AND.i.GE.j ) )
2981  $ THEN
2982  a( i, j ) = cbeg( reset ) + transl
2983  IF( i.NE.j )THEN
2984 * Set some elements to zero
2985  IF( n.GT.3.AND.j.EQ.n/2 )
2986  $ a( i, j ) = zero
2987  IF( her )THEN
2988  a( j, i ) = conjg( a( i, j ) )
2989  ELSE IF( sym )THEN
2990  a( j, i ) = a( i, j )
2991  ELSE IF( tri )THEN
2992  a( j, i ) = zero
2993  END IF
2994  END IF
2995  END IF
2996  10 continue
2997  IF( her )
2998  $ a( j, j ) = cmplx( REAL( A( J, J ) ), rzero )
2999  IF( tri )
3000  $ a( j, j ) = a( j, j ) + one
3001  IF( unit )
3002  $ a( j, j ) = one
3003  20 continue
3004 *
3005 * Store elements in array AS in data structure required by routine.
3006 *
3007  IF( type.EQ.'GE' )THEN
3008  DO 50 j = 1, n
3009  DO 30 i = 1, m
3010  aa( i + ( j - 1 )*lda ) = a( i, j )
3011  30 continue
3012  DO 40 i = m + 1, lda
3013  aa( i + ( j - 1 )*lda ) = rogue
3014  40 continue
3015  50 continue
3016  ELSE IF( type.EQ.'HE'.OR.type.EQ.'SY'.OR.type.EQ.'TR' )THEN
3017  DO 90 j = 1, n
3018  IF( upper )THEN
3019  ibeg = 1
3020  IF( unit )THEN
3021  iend = j - 1
3022  ELSE
3023  iend = j
3024  END IF
3025  ELSE
3026  IF( unit )THEN
3027  ibeg = j + 1
3028  ELSE
3029  ibeg = j
3030  END IF
3031  iend = n
3032  END IF
3033  DO 60 i = 1, ibeg - 1
3034  aa( i + ( j - 1 )*lda ) = rogue
3035  60 continue
3036  DO 70 i = ibeg, iend
3037  aa( i + ( j - 1 )*lda ) = a( i, j )
3038  70 continue
3039  DO 80 i = iend + 1, lda
3040  aa( i + ( j - 1 )*lda ) = rogue
3041  80 continue
3042  IF( her )THEN
3043  jj = j + ( j - 1 )*lda
3044  aa( jj ) = cmplx( REAL( AA( JJ ) ), rrogue )
3045  END IF
3046  90 continue
3047  END IF
3048  return
3049 *
3050 * End of CMAKE.
3051 *
3052  END
3053  SUBROUTINE cmmch( TRANSA, TRANSB, M, N, KK, ALPHA, A, LDA, B, LDB,
3054  $ beta, c, ldc, ct, g, cc, ldcc, eps, err, fatal,
3055  $ nout, mv )
3056 *
3057 * Checks the results of the computational tests.
3058 *
3059 * Auxiliary routine for test program for Level 3 Blas.
3060 *
3061 * -- Written on 8-February-1989.
3062 * Jack Dongarra, Argonne National Laboratory.
3063 * Iain Duff, AERE Harwell.
3064 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
3065 * Sven Hammarling, Numerical Algorithms Group Ltd.
3066 *
3067 * .. Parameters ..
3068  COMPLEX zero
3069  parameter( zero = ( 0.0, 0.0 ) )
3070  REAL rzero, rone
3071  parameter( rzero = 0.0, rone = 1.0 )
3072 * .. Scalar Arguments ..
3073  COMPLEX alpha, beta
3074  REAL eps, err
3075  INTEGER kk, lda, ldb, ldc, ldcc, m, n, nout
3076  LOGICAL fatal, mv
3077  CHARACTER*1 transa, transb
3078 * .. Array Arguments ..
3079  COMPLEX a( lda, * ), b( ldb, * ), c( ldc, * ),
3080  $ cc( ldcc, * ), ct( * )
3081  REAL g( * )
3082 * .. Local Scalars ..
3083  COMPLEX cl
3084  REAL erri
3085  INTEGER i, j, k
3086  LOGICAL ctrana, ctranb, trana, tranb
3087 * .. Intrinsic Functions ..
3088  INTRINSIC abs, aimag, conjg, max, REAL, sqrt
3089 * .. Statement Functions ..
3090  REAL abs1
3091 * .. Statement Function definitions ..
3092  abs1( cl ) = abs( REAL( CL ) ) + abs( aimag( cl ) )
3093 * .. Executable Statements ..
3094  trana = transa.EQ.'T'.OR.transa.EQ.'C'
3095  tranb = transb.EQ.'T'.OR.transb.EQ.'C'
3096  ctrana = transa.EQ.'C'
3097  ctranb = transb.EQ.'C'
3098 *
3099 * Compute expected result, one column at a time, in CT using data
3100 * in A, B and C.
3101 * Compute gauges in G.
3102 *
3103  DO 220 j = 1, n
3104 *
3105  DO 10 i = 1, m
3106  ct( i ) = zero
3107  g( i ) = rzero
3108  10 continue
3109  IF( .NOT.trana.AND..NOT.tranb )THEN
3110  DO 30 k = 1, kk
3111  DO 20 i = 1, m
3112  ct( i ) = ct( i ) + a( i, k )*b( k, j )
3113  g( i ) = g( i ) + abs1( a( i, k ) )*abs1( b( k, j ) )
3114  20 continue
3115  30 continue
3116  ELSE IF( trana.AND..NOT.tranb )THEN
3117  IF( ctrana )THEN
3118  DO 50 k = 1, kk
3119  DO 40 i = 1, m
3120  ct( i ) = ct( i ) + conjg( a( k, i ) )*b( k, j )
3121  g( i ) = g( i ) + abs1( a( k, i ) )*
3122  $ abs1( b( k, j ) )
3123  40 continue
3124  50 continue
3125  ELSE
3126  DO 70 k = 1, kk
3127  DO 60 i = 1, m
3128  ct( i ) = ct( i ) + a( k, i )*b( k, j )
3129  g( i ) = g( i ) + abs1( a( k, i ) )*
3130  $ abs1( b( k, j ) )
3131  60 continue
3132  70 continue
3133  END IF
3134  ELSE IF( .NOT.trana.AND.tranb )THEN
3135  IF( ctranb )THEN
3136  DO 90 k = 1, kk
3137  DO 80 i = 1, m
3138  ct( i ) = ct( i ) + a( i, k )*conjg( b( j, k ) )
3139  g( i ) = g( i ) + abs1( a( i, k ) )*
3140  $ abs1( b( j, k ) )
3141  80 continue
3142  90 continue
3143  ELSE
3144  DO 110 k = 1, kk
3145  DO 100 i = 1, m
3146  ct( i ) = ct( i ) + a( i, k )*b( j, k )
3147  g( i ) = g( i ) + abs1( a( i, k ) )*
3148  $ abs1( b( j, k ) )
3149  100 continue
3150  110 continue
3151  END IF
3152  ELSE IF( trana.AND.tranb )THEN
3153  IF( ctrana )THEN
3154  IF( ctranb )THEN
3155  DO 130 k = 1, kk
3156  DO 120 i = 1, m
3157  ct( i ) = ct( i ) + conjg( a( k, i ) )*
3158  $ conjg( b( j, k ) )
3159  g( i ) = g( i ) + abs1( a( k, i ) )*
3160  $ abs1( b( j, k ) )
3161  120 continue
3162  130 continue
3163  ELSE
3164  DO 150 k = 1, kk
3165  DO 140 i = 1, m
3166  ct( i ) = ct( i ) + conjg( a( k, i ) )*b( j, k )
3167  g( i ) = g( i ) + abs1( a( k, i ) )*
3168  $ abs1( b( j, k ) )
3169  140 continue
3170  150 continue
3171  END IF
3172  ELSE
3173  IF( ctranb )THEN
3174  DO 170 k = 1, kk
3175  DO 160 i = 1, m
3176  ct( i ) = ct( i ) + a( k, i )*conjg( b( j, k ) )
3177  g( i ) = g( i ) + abs1( a( k, i ) )*
3178  $ abs1( b( j, k ) )
3179  160 continue
3180  170 continue
3181  ELSE
3182  DO 190 k = 1, kk
3183  DO 180 i = 1, m
3184  ct( i ) = ct( i ) + a( k, i )*b( j, k )
3185  g( i ) = g( i ) + abs1( a( k, i ) )*
3186  $ abs1( b( j, k ) )
3187  180 continue
3188  190 continue
3189  END IF
3190  END IF
3191  END IF
3192  DO 200 i = 1, m
3193  ct( i ) = alpha*ct( i ) + beta*c( i, j )
3194  g( i ) = abs1( alpha )*g( i ) +
3195  $ abs1( beta )*abs1( c( i, j ) )
3196  200 continue
3197 *
3198 * Compute the error ratio for this result.
3199 *
3200  err = zero
3201  DO 210 i = 1, m
3202  erri = abs1( ct( i ) - cc( i, j ) )/eps
3203  IF( g( i ).NE.rzero )
3204  $ erri = erri/g( i )
3205  err = max( err, erri )
3206  IF( err*sqrt( eps ).GE.rone )
3207  $ go to 230
3208  210 continue
3209 *
3210  220 continue
3211 *
3212 * If the loop completes, all results are at least half accurate.
3213  go to 250
3214 *
3215 * Report fatal error.
3216 *
3217  230 fatal = .true.
3218  WRITE( nout, fmt = 9999 )
3219  DO 240 i = 1, m
3220  IF( mv )THEN
3221  WRITE( nout, fmt = 9998 )i, ct( i ), cc( i, j )
3222  ELSE
3223  WRITE( nout, fmt = 9998 )i, cc( i, j ), ct( i )
3224  END IF
3225  240 continue
3226  IF( n.GT.1 )
3227  $ WRITE( nout, fmt = 9997 )j
3228 *
3229  250 continue
3230  return
3231 *
3232  9999 format( ' ******* FATAL ERROR - COMPUTED RESULT IS LESS THAN HAL',
3233  $ 'F ACCURATE *******', /' EXPECTED RE',
3234  $ 'SULT COMPUTED RESULT' )
3235  9998 format( 1x, i7, 2( ' (', g15.6, ',', g15.6, ')' ) )
3236  9997 format( ' THESE ARE THE RESULTS FOR COLUMN ', i3 )
3237 *
3238 * End of CMMCH.
3239 *
3240  END
3241  LOGICAL FUNCTION lce( RI, RJ, LR )
3242 *
3243 * Tests if two arrays are identical.
3244 *
3245 * Auxiliary routine for test program for Level 3 Blas.
3246 *
3247 * -- Written on 8-February-1989.
3248 * Jack Dongarra, Argonne National Laboratory.
3249 * Iain Duff, AERE Harwell.
3250 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
3251 * Sven Hammarling, Numerical Algorithms Group Ltd.
3252 *
3253 * .. Scalar Arguments ..
3254  INTEGER lr
3255 * .. Array Arguments ..
3256  COMPLEX ri( * ), rj( * )
3257 * .. Local Scalars ..
3258  INTEGER i
3259 * .. Executable Statements ..
3260  DO 10 i = 1, lr
3261  IF( ri( i ).NE.rj( i ) )
3262  $ go to 20
3263  10 continue
3264  lce = .true.
3265  go to 30
3266  20 continue
3267  lce = .false.
3268  30 return
3269 *
3270 * End of LCE.
3271 *
3272  END
3273  LOGICAL FUNCTION lceres( TYPE, UPLO, M, N, AA, AS, LDA )
3274 *
3275 * Tests if selected elements in two arrays are equal.
3276 *
3277 * TYPE is 'GE' or 'HE' or 'SY'.
3278 *
3279 * Auxiliary routine for test program for Level 3 Blas.
3280 *
3281 * -- Written on 8-February-1989.
3282 * Jack Dongarra, Argonne National Laboratory.
3283 * Iain Duff, AERE Harwell.
3284 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
3285 * Sven Hammarling, Numerical Algorithms Group Ltd.
3286 *
3287 * .. Scalar Arguments ..
3288  INTEGER lda, m, n
3289  CHARACTER*1 uplo
3290  CHARACTER*2 type
3291 * .. Array Arguments ..
3292  COMPLEX aa( lda, * ), as( lda, * )
3293 * .. Local Scalars ..
3294  INTEGER i, ibeg, iend, j
3295  LOGICAL upper
3296 * .. Executable Statements ..
3297  upper = uplo.EQ.'U'
3298  IF( type.EQ.'GE' )THEN
3299  DO 20 j = 1, n
3300  DO 10 i = m + 1, lda
3301  IF( aa( i, j ).NE.as( i, j ) )
3302  $ go to 70
3303  10 continue
3304  20 continue
3305  ELSE IF( type.EQ.'HE'.OR.type.EQ.'SY' )THEN
3306  DO 50 j = 1, n
3307  IF( upper )THEN
3308  ibeg = 1
3309  iend = j
3310  ELSE
3311  ibeg = j
3312  iend = n
3313  END IF
3314  DO 30 i = 1, ibeg - 1
3315  IF( aa( i, j ).NE.as( i, j ) )
3316  $ go to 70
3317  30 continue
3318  DO 40 i = iend + 1, lda
3319  IF( aa( i, j ).NE.as( i, j ) )
3320  $ go to 70
3321  40 continue
3322  50 continue
3323  END IF
3324 *
3325  lceres = .true.
3326  go to 80
3327  70 continue
3328  lceres = .false.
3329  80 return
3330 *
3331 * End of LCERES.
3332 *
3333  END
3334  COMPLEX FUNCTION cbeg( RESET )
3335 *
3336 * Generates complex numbers as pairs of random numbers uniformly
3337 * distributed between -0.5 and 0.5.
3338 *
3339 * Auxiliary routine for test program for Level 3 Blas.
3340 *
3341 * -- Written on 8-February-1989.
3342 * Jack Dongarra, Argonne National Laboratory.
3343 * Iain Duff, AERE Harwell.
3344 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
3345 * Sven Hammarling, Numerical Algorithms Group Ltd.
3346 *
3347 * .. Scalar Arguments ..
3348  LOGICAL reset
3349 * .. Local Scalars ..
3350  INTEGER i, ic, j, mi, mj
3351 * .. Save statement ..
3352  SAVE i, ic, j, mi, mj
3353 * .. Intrinsic Functions ..
3354  INTRINSIC cmplx
3355 * .. Executable Statements ..
3356  IF( reset )THEN
3357 * Initialize local variables.
3358  mi = 891
3359  mj = 457
3360  i = 7
3361  j = 7
3362  ic = 0
3363  reset = .false.
3364  END IF
3365 *
3366 * The sequence of values of I or J is bounded between 1 and 999.
3367 * If initial I or J = 1,2,3,6,7 or 9, the period will be 50.
3368 * If initial I or J = 4 or 8, the period will be 25.
3369 * If initial I or J = 5, the period will be 10.
3370 * IC is used to break up the period by skipping 1 value of I or J
3371 * in 6.
3372 *
3373  ic = ic + 1
3374  10 i = i*mi
3375  j = j*mj
3376  i = i - 1000*( i/1000 )
3377  j = j - 1000*( j/1000 )
3378  IF( ic.GE.5 )THEN
3379  ic = 0
3380  go to 10
3381  END IF
3382  cbeg = cmplx( ( i - 500 )/1001.0, ( j - 500 )/1001.0 )
3383  return
3384 *
3385 * End of CBEG.
3386 *
3387  END
3388  REAL FUNCTION sdiff( X, Y )
3389 *
3390 * Auxiliary routine for test program for Level 3 Blas.
3391 *
3392 * -- Written on 8-February-1989.
3393 * Jack Dongarra, Argonne National Laboratory.
3394 * Iain Duff, AERE Harwell.
3395 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
3396 * Sven Hammarling, Numerical Algorithms Group Ltd.
3397 *
3398 * .. Scalar Arguments ..
3399  REAL x, y
3400 * .. Executable Statements ..
3401  sdiff = x - y
3402  return
3403 *
3404 * End of SDIFF.
3405 *
3406  END
3407  SUBROUTINE chkxer( SRNAMT, INFOT, NOUT, LERR, OK )
3408 *
3409 * Tests whether XERBLA has detected an error when it should.
3410 *
3411 * Auxiliary routine for test program for Level 3 Blas.
3412 *
3413 * -- Written on 8-February-1989.
3414 * Jack Dongarra, Argonne National Laboratory.
3415 * Iain Duff, AERE Harwell.
3416 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
3417 * Sven Hammarling, Numerical Algorithms Group Ltd.
3418 *
3419 * .. Scalar Arguments ..
3420  INTEGER infot, nout
3421  LOGICAL lerr, ok
3422  CHARACTER*6 srnamt
3423 * .. Executable Statements ..
3424  IF( .NOT.lerr )THEN
3425  WRITE( nout, fmt = 9999 )infot, srnamt
3426  ok = .false.
3427  END IF
3428  lerr = .false.
3429  return
3430 *
3431  9999 format( ' ***** ILLEGAL VALUE OF PARAMETER NUMBER ', i2, ' NOT D',
3432  $ 'ETECTED BY ', a6, ' *****' )
3433 *
3434 * End of CHKXER.
3435 *
3436  END
3437  SUBROUTINE xerbla( SRNAME, INFO )
3438 *
3439 * This is a special version of XERBLA to be used only as part of
3440 * the test program for testing error exits from the Level 3 BLAS
3441 * routines.
3442 *
3443 * XERBLA is an error handler for the Level 3 BLAS routines.
3444 *
3445 * It is called by the Level 3 BLAS routines if an input parameter is
3446 * invalid.
3447 *
3448 * Auxiliary routine for test program for Level 3 Blas.
3449 *
3450 * -- Written on 8-February-1989.
3451 * Jack Dongarra, Argonne National Laboratory.
3452 * Iain Duff, AERE Harwell.
3453 * Jeremy Du Croz, Numerical Algorithms Group Ltd.
3454 * Sven Hammarling, Numerical Algorithms Group Ltd.
3455 *
3456 * .. Scalar Arguments ..
3457  INTEGER info
3458  CHARACTER*6 srname
3459 * .. Scalars in Common ..
3460  INTEGER infot, nout
3461  LOGICAL lerr, ok
3462  CHARACTER*6 srnamt
3463 * .. Common blocks ..
3464  common /infoc/infot, nout, ok, lerr
3465  common /srnamc/srnamt
3466 * .. Executable Statements ..
3467  lerr = .true.
3468  IF( info.NE.infot )THEN
3469  IF( infot.NE.0 )THEN
3470  WRITE( nout, fmt = 9999 )info, infot
3471  ELSE
3472  WRITE( nout, fmt = 9997 )info
3473  END IF
3474  ok = .false.
3475  END IF
3476  IF( srname.NE.srnamt )THEN
3477  WRITE( nout, fmt = 9998 )srname, srnamt
3478  ok = .false.
3479  END IF
3480  return
3481 *
3482  9999 format( ' ******* XERBLA WAS CALLED WITH INFO = ', i6, ' INSTEAD',
3483  $ ' OF ', i2, ' *******' )
3484  9998 format( ' ******* XERBLA WAS CALLED WITH SRNAME = ', a6, ' INSTE',
3485  $ 'AD OF ', a6, ' *******' )
3486  9997 format( ' ******* XERBLA WAS CALLED WITH INFO = ', i6,
3487  $ ' *******' )
3488 *
3489 * End of XERBLA
3490 *
3491  END
3492