LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine cchksy_rook ( logical, dimension( * )  DOTYPE,
integer  NN,
integer, dimension( * )  NVAL,
integer  NNB,
integer, dimension( * )  NBVAL,
integer  NNS,
integer, dimension( * )  NSVAL,
real  THRESH,
logical  TSTERR,
integer  NMAX,
complex, dimension( * )  A,
complex, dimension( * )  AFAC,
complex, dimension( * )  AINV,
complex, dimension( * )  B,
complex, dimension( * )  X,
complex, dimension( * )  XACT,
complex, dimension( * )  WORK,
real, dimension( * )  RWORK,
integer, dimension( * )  IWORK,
integer  NOUT 
)

CCHKSY_ROOK

Purpose:
 CCHKSY_ROOK tests CSYTRF_ROOK, -TRI_ROOK, -TRS_ROOK,
 and -CON_ROOK.
Parameters
[in]DOTYPE
          DOTYPE is LOGICAL array, dimension (NTYPES)
          The matrix types to be used for testing.  Matrices of type j
          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
[in]NN
          NN is INTEGER
          The number of values of N contained in the vector NVAL.
[in]NVAL
          NVAL is INTEGER array, dimension (NN)
          The values of the matrix dimension N.
[in]NNB
          NNB is INTEGER
          The number of values of NB contained in the vector NBVAL.
[in]NBVAL
          NBVAL is INTEGER array, dimension (NBVAL)
          The values of the blocksize NB.
[in]NNS
          NNS is INTEGER
          The number of values of NRHS contained in the vector NSVAL.
[in]NSVAL
          NSVAL is INTEGER array, dimension (NNS)
          The values of the number of right hand sides NRHS.
[in]THRESH
          THRESH is REAL
          The threshold value for the test ratios.  A result is
          included in the output file if RESULT >= THRESH.  To have
          every test ratio printed, use THRESH = 0.
[in]TSTERR
          TSTERR is LOGICAL
          Flag that indicates whether error exits are to be tested.
[in]NMAX
          NMAX is INTEGER
          The maximum value permitted for N, used in dimensioning the
          work arrays.
[out]A
          A is COMPLEX array, dimension (NMAX*NMAX)
[out]AFAC
          AFAC is COMPLEX array, dimension (NMAX*NMAX)
[out]AINV
          AINV is COMPLEX array, dimension (NMAX*NMAX)
[out]B
          B is COMPLEX array, dimension (NMAX*NSMAX)
          where NSMAX is the largest entry in NSVAL.
[out]X
          X is COMPLEX array, dimension (NMAX*NSMAX)
[out]XACT
          XACT is COMPLEX array, dimension (NMAX*NSMAX)
[out]WORK
          WORK is COMPLEX array, dimension (NMAX*max(3,NSMAX))
[out]RWORK
          RWORK is REAL array, dimension (max(NMAX,2*NSMAX))
[out]IWORK
          IWORK is INTEGER array, dimension (2*NMAX)
[in]NOUT
          NOUT is INTEGER
          The unit number for output.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015

Definition at line 174 of file cchksy_rook.f.

174 *
175 * -- LAPACK test routine (version 3.6.0) --
176 * -- LAPACK is a software package provided by Univ. of Tennessee, --
177 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
178 * November 2015
179 *
180 * .. Scalar Arguments ..
181  LOGICAL tsterr
182  INTEGER nmax, nn, nnb, nns, nout
183  REAL thresh
184 * ..
185 * .. Array Arguments ..
186  LOGICAL dotype( * )
187  INTEGER iwork( * ), nbval( * ), nsval( * ), nval( * )
188  REAL rwork( * )
189  COMPLEX a( * ), afac( * ), ainv( * ), b( * ),
190  $ work( * ), x( * ), xact( * )
191 * ..
192 *
193 * =====================================================================
194 *
195 * .. Parameters ..
196  REAL zero, one
197  parameter ( zero = 0.0e+0, one = 1.0e+0 )
198  REAL onehalf
199  parameter ( onehalf = 0.5e+0 )
200  REAL eight, sevten
201  parameter ( eight = 8.0e+0, sevten = 17.0e+0 )
202  COMPLEX czero
203  parameter ( czero = ( 0.0e+0, 0.0e+0 ) )
204  INTEGER ntypes
205  parameter ( ntypes = 11 )
206  INTEGER ntests
207  parameter ( ntests = 7 )
208 * ..
209 * .. Local Scalars ..
210  LOGICAL trfcon, zerot
211  CHARACTER dist, TYPE, uplo, xtype
212  CHARACTER*3 path, matpath
213  INTEGER i, i1, i2, imat, in, inb, info, ioff, irhs,
214  $ itemp, itemp2, iuplo, izero, j, k, kl, ku, lda,
215  $ lwork, mode, n, nb, nerrs, nfail, nimat, nrhs,
216  $ nrun, nt
217  REAL alpha, anorm, cndnum, const, sing_max,
218  $ sing_min, rcond, rcondc, stemp
219 * ..
220 * .. Local Arrays ..
221  CHARACTER uplos( 2 )
222  INTEGER iseed( 4 ), iseedy( 4 )
223  REAL result( ntests )
224  COMPLEX block( 2, 2 ), cdummy( 1 )
225 * ..
226 * .. External Functions ..
227  REAL clange, clansy, sget06
228  EXTERNAL clange, clansy, sget06
229 * ..
230 * .. External Subroutines ..
231  EXTERNAL alaerh, alahd, alasum, cerrsy, cgesvd, cget04,
235 * ..
236 * .. Intrinsic Functions ..
237  INTRINSIC max, min, sqrt
238 * ..
239 * .. Scalars in Common ..
240  LOGICAL lerr, ok
241  CHARACTER*32 srnamt
242  INTEGER infot, nunit
243 * ..
244 * .. Common blocks ..
245  COMMON / infoc / infot, nunit, ok, lerr
246  COMMON / srnamc / srnamt
247 * ..
248 * .. Data statements ..
249  DATA iseedy / 1988, 1989, 1990, 1991 /
250  DATA uplos / 'U', 'L' /
251 * ..
252 * .. Executable Statements ..
253 *
254 * Initialize constants and the random number seed.
255 *
256  alpha = ( one+sqrt( sevten ) ) / eight
257 *
258 * Test path
259 *
260  path( 1: 1 ) = 'Complex precision'
261  path( 2: 3 ) = 'SR'
262 *
263 * Path to generate matrices
264 *
265  matpath( 1: 1 ) = 'Complex precision'
266  matpath( 2: 3 ) = 'SY'
267 *
268  nrun = 0
269  nfail = 0
270  nerrs = 0
271  DO 10 i = 1, 4
272  iseed( i ) = iseedy( i )
273  10 CONTINUE
274 *
275 * Test the error exits
276 *
277  IF( tsterr )
278  $ CALL cerrsy( path, nout )
279  infot = 0
280 *
281 * Set the minimum block size for which the block routine should
282 * be used, which will be later returned by ILAENV
283 *
284  CALL xlaenv( 2, 2 )
285 *
286 * Do for each value of N in NVAL
287 *
288  DO 270 in = 1, nn
289  n = nval( in )
290  lda = max( n, 1 )
291  xtype = 'N'
292  nimat = ntypes
293  IF( n.LE.0 )
294  $ nimat = 1
295 *
296  izero = 0
297 *
298 * Do for each value of matrix type IMAT
299 *
300  DO 260 imat = 1, nimat
301 *
302 * Do the tests only if DOTYPE( IMAT ) is true.
303 *
304  IF( .NOT.dotype( imat ) )
305  $ GO TO 260
306 *
307 * Skip types 3, 4, 5, or 6 if the matrix size is too small.
308 *
309  zerot = imat.GE.3 .AND. imat.LE.6
310  IF( zerot .AND. n.LT.imat-2 )
311  $ GO TO 260
312 *
313 * Do first for UPLO = 'U', then for UPLO = 'L'
314 *
315  DO 250 iuplo = 1, 2
316  uplo = uplos( iuplo )
317 *
318 * Begin generate test matrix A.
319 *
320  IF( imat.NE.ntypes ) THEN
321 *
322 * Set up parameters with CLATB4 for the matrix generator
323 * based on the type of matrix to be generated.
324 *
325  CALL clatb4( matpath, imat, n, n, TYPE, kl, ku, anorm,
326  $ mode, cndnum, dist )
327 *
328 * Generate a matrix with CLATMS.
329 *
330  srnamt = 'CLATMS'
331  CALL clatms( n, n, dist, iseed, TYPE, rwork, mode,
332  $ cndnum, anorm, kl, ku, uplo, a, lda,
333  $ work, info )
334 *
335 * Check error code from CLATMS and handle error.
336 *
337  IF( info.NE.0 ) THEN
338  CALL alaerh( path, 'CLATMS', info, 0, uplo, n, n,
339  $ -1, -1, -1, imat, nfail, nerrs, nout )
340 *
341 * Skip all tests for this generated matrix
342 *
343  GO TO 250
344  END IF
345 *
346 * For matrix types 3-6, zero one or more rows and
347 * columns of the matrix to test that INFO is returned
348 * correctly.
349 *
350  IF( zerot ) THEN
351  IF( imat.EQ.3 ) THEN
352  izero = 1
353  ELSE IF( imat.EQ.4 ) THEN
354  izero = n
355  ELSE
356  izero = n / 2 + 1
357  END IF
358 *
359  IF( imat.LT.6 ) THEN
360 *
361 * Set row and column IZERO to zero.
362 *
363  IF( iuplo.EQ.1 ) THEN
364  ioff = ( izero-1 )*lda
365  DO 20 i = 1, izero - 1
366  a( ioff+i ) = czero
367  20 CONTINUE
368  ioff = ioff + izero
369  DO 30 i = izero, n
370  a( ioff ) = czero
371  ioff = ioff + lda
372  30 CONTINUE
373  ELSE
374  ioff = izero
375  DO 40 i = 1, izero - 1
376  a( ioff ) = czero
377  ioff = ioff + lda
378  40 CONTINUE
379  ioff = ioff - izero
380  DO 50 i = izero, n
381  a( ioff+i ) = czero
382  50 CONTINUE
383  END IF
384  ELSE
385  IF( iuplo.EQ.1 ) THEN
386 *
387 * Set the first IZERO rows and columns to zero.
388 *
389  ioff = 0
390  DO 70 j = 1, n
391  i2 = min( j, izero )
392  DO 60 i = 1, i2
393  a( ioff+i ) = czero
394  60 CONTINUE
395  ioff = ioff + lda
396  70 CONTINUE
397  ELSE
398 *
399 * Set the last IZERO rows and columns to zero.
400 *
401  ioff = 0
402  DO 90 j = 1, n
403  i1 = max( j, izero )
404  DO 80 i = i1, n
405  a( ioff+i ) = czero
406  80 CONTINUE
407  ioff = ioff + lda
408  90 CONTINUE
409  END IF
410  END IF
411  ELSE
412  izero = 0
413  END IF
414 *
415  ELSE
416 *
417 * For matrix kind IMAT = 11, generate special block
418 * diagonal matrix to test alternate code
419 * for the 2 x 2 blocks.
420 *
421  CALL clatsy( uplo, n, a, lda, iseed )
422 *
423  END IF
424 *
425 * End generate test matrix A.
426 *
427 *
428 * Do for each value of NB in NBVAL
429 *
430  DO 240 inb = 1, nnb
431 *
432 * Set the optimal blocksize, which will be later
433 * returned by ILAENV.
434 *
435  nb = nbval( inb )
436  CALL xlaenv( 1, nb )
437 *
438 * Copy the test matrix A into matrix AFAC which
439 * will be factorized in place. This is needed to
440 * preserve the test matrix A for subsequent tests.
441 *
442  CALL clacpy( uplo, n, n, a, lda, afac, lda )
443 *
444 * Compute the L*D*L**T or U*D*U**T factorization of the
445 * matrix. IWORK stores details of the interchanges and
446 * the block structure of D. AINV is a work array for
447 * block factorization, LWORK is the length of AINV.
448 *
449  lwork = max( 2, nb )*lda
450  srnamt = 'CSYTRF_ROOK'
451  CALL csytrf_rook( uplo, n, afac, lda, iwork, ainv,
452  $ lwork, info )
453 *
454 * Adjust the expected value of INFO to account for
455 * pivoting.
456 *
457  k = izero
458  IF( k.GT.0 ) THEN
459  100 CONTINUE
460  IF( iwork( k ).LT.0 ) THEN
461  IF( iwork( k ).NE.-k ) THEN
462  k = -iwork( k )
463  GO TO 100
464  END IF
465  ELSE IF( iwork( k ).NE.k ) THEN
466  k = iwork( k )
467  GO TO 100
468  END IF
469  END IF
470 *
471 * Check error code from CSYTRF_ROOK and handle error.
472 *
473  IF( info.NE.k)
474  $ CALL alaerh( path, 'CSYTRF_ROOK', info, k,
475  $ uplo, n, n, -1, -1, nb, imat,
476  $ nfail, nerrs, nout )
477 *
478 * Set the condition estimate flag if the INFO is not 0.
479 *
480  IF( info.NE.0 ) THEN
481  trfcon = .true.
482  ELSE
483  trfcon = .false.
484  END IF
485 *
486 *+ TEST 1
487 * Reconstruct matrix from factors and compute residual.
488 *
489  CALL csyt01_rook( uplo, n, a, lda, afac, lda, iwork,
490  $ ainv, lda, rwork, result( 1 ) )
491  nt = 1
492 *
493 *+ TEST 2
494 * Form the inverse and compute the residual,
495 * if the factorization was competed without INFO > 0
496 * (i.e. there is no zero rows and columns).
497 * Do it only for the first block size.
498 *
499  IF( inb.EQ.1 .AND. .NOT.trfcon ) THEN
500  CALL clacpy( uplo, n, n, afac, lda, ainv, lda )
501  srnamt = 'CSYTRI_ROOK'
502  CALL csytri_rook( uplo, n, ainv, lda, iwork, work,
503  $ info )
504 *
505 * Check error code from CSYTRI_ROOK and handle error.
506 *
507  IF( info.NE.0 )
508  $ CALL alaerh( path, 'CSYTRI_ROOK', info, -1,
509  $ uplo, n, n, -1, -1, -1, imat,
510  $ nfail, nerrs, nout )
511 *
512 * Compute the residual for a symmetric matrix times
513 * its inverse.
514 *
515  CALL csyt03( uplo, n, a, lda, ainv, lda, work, lda,
516  $ rwork, rcondc, result( 2 ) )
517  nt = 2
518  END IF
519 *
520 * Print information about the tests that did not pass
521 * the threshold.
522 *
523  DO 110 k = 1, nt
524  IF( result( k ).GE.thresh ) THEN
525  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
526  $ CALL alahd( nout, path )
527  WRITE( nout, fmt = 9999 )uplo, n, nb, imat, k,
528  $ result( k )
529  nfail = nfail + 1
530  END IF
531  110 CONTINUE
532  nrun = nrun + nt
533 *
534 *+ TEST 3
535 * Compute largest element in U or L
536 *
537  result( 3 ) = zero
538  stemp = zero
539 *
540  const = ( ( alpha**2-one ) / ( alpha**2-onehalf ) ) /
541  $ ( one-alpha )
542 *
543  IF( iuplo.EQ.1 ) THEN
544 *
545 * Compute largest element in U
546 *
547  k = n
548  120 CONTINUE
549  IF( k.LE.1 )
550  $ GO TO 130
551 *
552  IF( iwork( k ).GT.zero ) THEN
553 *
554 * Get max absolute value from elements
555 * in column k in in U
556 *
557  stemp = clange( 'M', k-1, 1,
558  $ afac( ( k-1 )*lda+1 ), lda, rwork )
559  ELSE
560 *
561 * Get max absolute value from elements
562 * in columns k and k-1 in U
563 *
564  stemp = clange( 'M', k-2, 2,
565  $ afac( ( k-2 )*lda+1 ), lda, rwork )
566  k = k - 1
567 *
568  END IF
569 *
570 * STEMP should be bounded by CONST
571 *
572  stemp = stemp - const + thresh
573  IF( stemp.GT.result( 3 ) )
574  $ result( 3 ) = stemp
575 *
576  k = k - 1
577 *
578  GO TO 120
579  130 CONTINUE
580 *
581  ELSE
582 *
583 * Compute largest element in L
584 *
585  k = 1
586  140 CONTINUE
587  IF( k.GE.n )
588  $ GO TO 150
589 *
590  IF( iwork( k ).GT.zero ) THEN
591 *
592 * Get max absolute value from elements
593 * in column k in in L
594 *
595  stemp = clange( 'M', n-k, 1,
596  $ afac( ( k-1 )*lda+k+1 ), lda, rwork )
597  ELSE
598 *
599 * Get max absolute value from elements
600 * in columns k and k+1 in L
601 *
602  stemp = clange( 'M', n-k-1, 2,
603  $ afac( ( k-1 )*lda+k+2 ), lda, rwork )
604  k = k + 1
605 *
606  END IF
607 *
608 * STEMP should be bounded by CONST
609 *
610  stemp = stemp - const + thresh
611  IF( stemp.GT.result( 3 ) )
612  $ result( 3 ) = stemp
613 *
614  k = k + 1
615 *
616  GO TO 140
617  150 CONTINUE
618  END IF
619 *
620 *
621 *+ TEST 4
622 * Compute largest 2-Norm (condition number)
623 * of 2-by-2 diag blocks
624 *
625  result( 4 ) = zero
626  stemp = zero
627 *
628  const = ( ( alpha**2-one ) / ( alpha**2-onehalf ) )*
629  $ ( ( one + alpha ) / ( one - alpha ) )
630 *
631  IF( iuplo.EQ.1 ) THEN
632 *
633 * Loop backward for UPLO = 'U'
634 *
635  k = n
636  160 CONTINUE
637  IF( k.LE.1 )
638  $ GO TO 170
639 *
640  IF( iwork( k ).LT.zero ) THEN
641 *
642 * Get the two singular values
643 * (real and non-negative) of a 2-by-2 block,
644 * store them in RWORK array
645 *
646  block( 1, 1 ) = afac( ( k-2 )*lda+k-1 )
647  block( 1, 2 ) = afac( (k-1)*lda+k-1 )
648  block( 2, 1 ) = block( 1, 2 )
649  block( 2, 2 ) = afac( (k-1)*lda+k )
650 *
651  CALL cgesvd( 'N', 'N', 2, 2, block, 2, rwork,
652  $ cdummy, 1, cdummy, 1,
653  $ work, 6, rwork( 3 ), info )
654 *
655 *
656  sing_max = rwork( 1 )
657  sing_min = rwork( 2 )
658 *
659  stemp = sing_max / sing_min
660 *
661 * STEMP should be bounded by CONST
662 *
663  stemp = stemp - const + thresh
664  IF( stemp.GT.result( 4 ) )
665  $ result( 4 ) = stemp
666  k = k - 1
667 *
668  END IF
669 *
670  k = k - 1
671 *
672  GO TO 160
673  170 CONTINUE
674 *
675  ELSE
676 *
677 * Loop forward for UPLO = 'L'
678 *
679  k = 1
680  180 CONTINUE
681  IF( k.GE.n )
682  $ GO TO 190
683 *
684  IF( iwork( k ).LT.zero ) THEN
685 *
686 * Get the two singular values
687 * (real and non-negative) of a 2-by-2 block,
688 * store them in RWORK array
689 *
690  block( 1, 1 ) = afac( ( k-1 )*lda+k )
691  block( 2, 1 ) = afac( ( k-1 )*lda+k+1 )
692  block( 1, 2 ) = block( 2, 1 )
693  block( 2, 2 ) = afac( k*lda+k+1 )
694 *
695  CALL cgesvd( 'N', 'N', 2, 2, block, 2, rwork,
696  $ cdummy, 1, cdummy, 1,
697  $ work, 6, rwork(3), info )
698 *
699  sing_max = rwork( 1 )
700  sing_min = rwork( 2 )
701 *
702  stemp = sing_max / sing_min
703 *
704 * STEMP should be bounded by CONST
705 *
706  stemp = stemp - const + thresh
707  IF( stemp.GT.result( 4 ) )
708  $ result( 4 ) = stemp
709  k = k + 1
710 *
711  END IF
712 *
713  k = k + 1
714 *
715  GO TO 180
716  190 CONTINUE
717  END IF
718 *
719 * Print information about the tests that did not pass
720 * the threshold.
721 *
722  DO 200 k = 3, 4
723  IF( result( k ).GE.thresh ) THEN
724  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
725  $ CALL alahd( nout, path )
726  WRITE( nout, fmt = 9999 )uplo, n, nb, imat, k,
727  $ result( k )
728  nfail = nfail + 1
729  END IF
730  200 CONTINUE
731  nrun = nrun + 2
732 *
733 * Skip the other tests if this is not the first block
734 * size.
735 *
736  IF( inb.GT.1 )
737  $ GO TO 240
738 *
739 * Do only the condition estimate if INFO is not 0.
740 *
741  IF( trfcon ) THEN
742  rcondc = zero
743  GO TO 230
744  END IF
745 *
746 * Do for each value of NRHS in NSVAL.
747 *
748  DO 220 irhs = 1, nns
749  nrhs = nsval( irhs )
750 *
751 *+ TEST 5 ( Using TRS_ROOK)
752 * Solve and compute residual for A * X = B.
753 *
754 * Choose a set of NRHS random solution vectors
755 * stored in XACT and set up the right hand side B
756 *
757  srnamt = 'CLARHS'
758  CALL clarhs( matpath, xtype, uplo, ' ', n, n,
759  $ kl, ku, nrhs, a, lda, xact, lda,
760  $ b, lda, iseed, info )
761  CALL clacpy( 'Full', n, nrhs, b, lda, x, lda )
762 *
763  srnamt = 'CSYTRS_ROOK'
764  CALL csytrs_rook( uplo, n, nrhs, afac, lda, iwork,
765  $ x, lda, info )
766 *
767 * Check error code from CSYTRS_ROOK and handle error.
768 *
769  IF( info.NE.0 )
770  $ CALL alaerh( path, 'CSYTRS_ROOK', info, 0,
771  $ uplo, n, n, -1, -1, nrhs, imat,
772  $ nfail, nerrs, nout )
773 *
774  CALL clacpy( 'Full', n, nrhs, b, lda, work, lda )
775 *
776 * Compute the residual for the solution
777 *
778  CALL csyt02( uplo, n, nrhs, a, lda, x, lda, work,
779  $ lda, rwork, result( 5 ) )
780 *
781 *+ TEST 6
782 * Check solution from generated exact solution.
783 *
784  CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
785  $ result( 6 ) )
786 *
787 * Print information about the tests that did not pass
788 * the threshold.
789 *
790  DO 210 k = 5, 6
791  IF( result( k ).GE.thresh ) THEN
792  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
793  $ CALL alahd( nout, path )
794  WRITE( nout, fmt = 9998 )uplo, n, nrhs,
795  $ imat, k, result( k )
796  nfail = nfail + 1
797  END IF
798  210 CONTINUE
799  nrun = nrun + 2
800 *
801 * End do for each value of NRHS in NSVAL.
802 *
803  220 CONTINUE
804 *
805 *+ TEST 7
806 * Get an estimate of RCOND = 1/CNDNUM.
807 *
808  230 CONTINUE
809  anorm = clansy( '1', uplo, n, a, lda, rwork )
810  srnamt = 'CSYCON_ROOK'
811  CALL csycon_rook( uplo, n, afac, lda, iwork, anorm,
812  $ rcond, work, info )
813 *
814 * Check error code from CSYCON_ROOK and handle error.
815 *
816  IF( info.NE.0 )
817  $ CALL alaerh( path, 'CSYCON_ROOK', info, 0,
818  $ uplo, n, n, -1, -1, -1, imat,
819  $ nfail, nerrs, nout )
820 *
821 * Compute the test ratio to compare values of RCOND
822 *
823  result( 7 ) = sget06( rcond, rcondc )
824 *
825 * Print information about the tests that did not pass
826 * the threshold.
827 *
828  IF( result( 7 ).GE.thresh ) THEN
829  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
830  $ CALL alahd( nout, path )
831  WRITE( nout, fmt = 9997 )uplo, n, imat, 7,
832  $ result( 7 )
833  nfail = nfail + 1
834  END IF
835  nrun = nrun + 1
836  240 CONTINUE
837 *
838  250 CONTINUE
839  260 CONTINUE
840  270 CONTINUE
841 *
842 * Print a summary of the results.
843 *
844  CALL alasum( path, nout, nfail, nrun, nerrs )
845 *
846  9999 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ', NB =', i4, ', type ',
847  $ i2, ', test ', i2, ', ratio =', g12.5 )
848  9998 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ', NRHS=', i3, ', type ',
849  $ i2, ', test(', i2, ') =', g12.5 )
850  9997 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ',', 10x, ' type ', i2,
851  $ ', test(', i2, ') =', g12.5 )
852  RETURN
853 *
854 * End of CCHKSY_ROOK
855 *
subroutine csytri_rook(UPLO, N, A, LDA, IPIV, WORK, INFO)
CSYTRI_ROOK
Definition: csytri_rook.f:131
subroutine alahd(IOUNIT, PATH)
ALAHD
Definition: alahd.f:95
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:149
subroutine clatsy(UPLO, N, X, LDX, ISEED)
CLATSY
Definition: clatsy.f:91
subroutine clarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
CLARHS
Definition: clarhs.f:211
subroutine csyt02(UPLO, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
CSYT02
Definition: csyt02.f:129
subroutine cgesvd(JOBU, JOBVT, M, N, A, LDA, S, U, LDU, VT, LDVT, WORK, LWORK, RWORK, INFO)
CGESVD computes the singular value decomposition (SVD) for GE matrices
Definition: cgesvd.f:216
subroutine cerrsy(PATH, NUNIT)
CERRSY
Definition: cerrsy.f:57
real function sget06(RCOND, RCONDC)
SGET06
Definition: sget06.f:57
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:83
subroutine csytrs_rook(UPLO, N, NRHS, A, LDA, IPIV, B, LDB, INFO)
CSYTRS_ROOK
Definition: csytrs_rook.f:138
subroutine csyt03(UPLO, N, A, LDA, AINV, LDAINV, WORK, LDWORK, RWORK, RCOND, RESID)
CSYT03
Definition: csyt03.f:128
subroutine clatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
CLATMS
Definition: clatms.f:334
subroutine csytrf_rook(UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
CSYTRF_ROOK
Definition: csytrf_rook.f:210
real function clansy(NORM, UPLO, N, A, LDA, WORK)
CLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex symmetric matrix.
Definition: clansy.f:125
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine csycon_rook(UPLO, N, A, LDA, IPIV, ANORM, RCOND, WORK, INFO)
CSYCON_ROOK
Definition: csycon_rook.f:141
subroutine clatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
CLATB4
Definition: clatb4.f:123
subroutine cget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
CGET04
Definition: cget04.f:104
subroutine csyt01_rook(UPLO, N, A, LDA, AFAC, LDAFAC, IPIV, C, LDC, RWORK, RESID)
CSYT01_ROOK
Definition: csyt01_rook.f:127
subroutine alasum(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASUM
Definition: alasum.f:75

Here is the call graph for this function:

Here is the caller graph for this function: