LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Files Functions Typedefs Macros
zchksy_rook.f
Go to the documentation of this file.
1 *> \brief \b ZCHKSY_ROOK
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE ZCHKSY_ROOK( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL,
12 * THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X,
13 * XACT, WORK, RWORK, IWORK, NOUT )
14 *
15 * .. Scalar Arguments ..
16 * LOGICAL TSTERR
17 * INTEGER NMAX, NN, NNB, NNS, NOUT
18 * DOUBLE PRECISION THRESH
19 * ..
20 * .. Array Arguments ..
21 * LOGICAL DOTYPE( * )
22 * INTEGER IWORK( * ), NBVAL( * ), NSVAL( * ), NVAL( * )
23 * DOUBLE PRECISION RWORK( * )
24 * COMPLEX*16 A( * ), AFAC( * ), AINV( * ), B( * ),
25 * $ WORK( * ), X( * ), XACT( * )
26 * ..
27 *
28 *
29 *> \par Purpose:
30 * =============
31 *>
32 *> \verbatim
33 *>
34 *> ZCHKSY_ROOK tests ZSYTRF_ROOK, -TRI_ROOK, -TRS_ROOK,
35 *> and -CON_ROOK.
36 *> \endverbatim
37 *
38 * Arguments:
39 * ==========
40 *
41 *> \param[in] DOTYPE
42 *> \verbatim
43 *> DOTYPE is LOGICAL array, dimension (NTYPES)
44 *> The matrix types to be used for testing. Matrices of type j
45 *> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
46 *> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
47 *> \endverbatim
48 *>
49 *> \param[in] NN
50 *> \verbatim
51 *> NN is INTEGER
52 *> The number of values of N contained in the vector NVAL.
53 *> \endverbatim
54 *>
55 *> \param[in] NVAL
56 *> \verbatim
57 *> NVAL is INTEGER array, dimension (NN)
58 *> The values of the matrix dimension N.
59 *> \endverbatim
60 *>
61 *> \param[in] NNB
62 *> \verbatim
63 *> NNB is INTEGER
64 *> The number of values of NB contained in the vector NBVAL.
65 *> \endverbatim
66 *>
67 *> \param[in] NBVAL
68 *> \verbatim
69 *> NBVAL is INTEGER array, dimension (NBVAL)
70 *> The values of the blocksize NB.
71 *> \endverbatim
72 *>
73 *> \param[in] NNS
74 *> \verbatim
75 *> NNS is INTEGER
76 *> The number of values of NRHS contained in the vector NSVAL.
77 *> \endverbatim
78 *>
79 *> \param[in] NSVAL
80 *> \verbatim
81 *> NSVAL is INTEGER array, dimension (NNS)
82 *> The values of the number of right hand sides NRHS.
83 *> \endverbatim
84 *>
85 *> \param[in] THRESH
86 *> \verbatim
87 *> THRESH is DOUBLE PRECISION
88 *> The threshold value for the test ratios. A result is
89 *> included in the output file if RESULT >= THRESH. To have
90 *> every test ratio printed, use THRESH = 0.
91 *> \endverbatim
92 *>
93 *> \param[in] TSTERR
94 *> \verbatim
95 *> TSTERR is LOGICAL
96 *> Flag that indicates whether error exits are to be tested.
97 *> \endverbatim
98 *>
99 *> \param[in] NMAX
100 *> \verbatim
101 *> NMAX is INTEGER
102 *> The maximum value permitted for N, used in dimensioning the
103 *> work arrays.
104 *> \endverbatim
105 *>
106 *> \param[out] A
107 *> \verbatim
108 *> A is COMPLEX*16 array, dimension (NMAX*NMAX)
109 *> \endverbatim
110 *>
111 *> \param[out] AFAC
112 *> \verbatim
113 *> AFAC is COMPLEX*16 array, dimension (NMAX*NMAX)
114 *> \endverbatim
115 *>
116 *> \param[out] AINV
117 *> \verbatim
118 *> AINV is COMPLEX*16 array, dimension (NMAX*NMAX)
119 *> \endverbatim
120 *>
121 *> \param[out] B
122 *> \verbatim
123 *> B is COMPLEX*16 array, dimension (NMAX*NSMAX)
124 *> where NSMAX is the largest entry in NSVAL.
125 *> \endverbatim
126 *>
127 *> \param[out] X
128 *> \verbatim
129 *> X is COMPLEX*16 array, dimension (NMAX*NSMAX)
130 *> \endverbatim
131 *>
132 *> \param[out] XACT
133 *> \verbatim
134 *> XACT is COMPLEX*16 array, dimension (NMAX*NSMAX)
135 *> \endverbatim
136 *>
137 *> \param[out] WORK
138 *> \verbatim
139 *> WORK is COMPLEX*16 array, dimension (NMAX*max(3,NSMAX))
140 *> \endverbatim
141 *>
142 *> \param[out] RWORK
143 *> \verbatim
144 *> RWORK is DOUBLE PRECISION array, dimension (max(NMAX,2*NSMAX))
145 *> \endverbatim
146 *>
147 *> \param[out] IWORK
148 *> \verbatim
149 *> IWORK is INTEGER array, dimension (2*NMAX)
150 *> \endverbatim
151 *>
152 *> \param[in] NOUT
153 *> \verbatim
154 *> NOUT is INTEGER
155 *> The unit number for output.
156 *> \endverbatim
157 *
158 * Authors:
159 * ========
160 *
161 *> \author Univ. of Tennessee
162 *> \author Univ. of California Berkeley
163 *> \author Univ. of Colorado Denver
164 *> \author NAG Ltd.
165 *
166 *> \date November 2013
167 *
168 *> \ingroup complex16_lin
169 *
170 * =====================================================================
171  SUBROUTINE zchksy_rook( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL,
172  $ thresh, tsterr, nmax, a, afac, ainv, b, x,
173  $ xact, work, rwork, iwork, nout )
174 *
175 * -- LAPACK test routine (version 3.5.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 2013
179 *
180 * .. Scalar Arguments ..
181  LOGICAL TSTERR
182  INTEGER NMAX, NN, NNB, NNS, NOUT
183  DOUBLE PRECISION THRESH
184 * ..
185 * .. Array Arguments ..
186  LOGICAL DOTYPE( * )
187  INTEGER IWORK( * ), NBVAL( * ), NSVAL( * ), NVAL( * )
188  DOUBLE PRECISION RWORK( * )
189  COMPLEX*16 A( * ), AFAC( * ), AINV( * ), B( * ),
190  $ work( * ), x( * ), xact( * )
191 * ..
192 *
193 * =====================================================================
194 *
195 * .. Parameters ..
196  DOUBLE PRECISION ZERO, ONE
197  parameter( zero = 0.0d+0, one = 1.0d+0 )
198  DOUBLE PRECISION ONEHALF
199  parameter( onehalf = 0.5d+0 )
200  DOUBLE PRECISION EIGHT, SEVTEN
201  parameter( eight = 8.0d+0, sevten = 17.0d+0 )
202  COMPLEX*16 CZERO
203  parameter( czero = ( 0.0d+0, 0.0d+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  DOUBLE PRECISION ALPHA, ANORM, CNDNUM, CONST, DTEMP, LAM_MAX,
218  $ lam_min, rcond, rcondc
219 * ..
220 * .. Local Arrays ..
221  CHARACTER UPLOS( 2 )
222  INTEGER ISEED( 4 ), ISEEDY( 4 )
223  DOUBLE PRECISION RESULT( ntests )
224  COMPLEX*16 BLOCK( 2, 2 ), ZDUMMY( 1 )
225 * ..
226 * .. External Functions ..
227  DOUBLE PRECISION DGET06, ZLANGE, ZLANSY
228  EXTERNAL dget06, zlange, zlansy
229 * ..
230 * .. External Subroutines ..
231  EXTERNAL alaerh, alahd, alasum, zerrsy, zgeevx, zget04,
235 * ..
236 * .. Intrinsic Functions ..
237  INTRINSIC abs, 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 ) = 'Zomplex precision'
261  path( 2: 3 ) = 'SR'
262 *
263 * Path to generate matrices
264 *
265  matpath( 1: 1 ) = 'Zomplex 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 zerrsy( 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 ZLATB4 for the matrix generator
323 * based on the type of matrix to be generated.
324 *
325  CALL zlatb4( matpath, imat, n, n, TYPE, KL, KU, ANORM,
326  $ mode, cndnum, dist )
327 *
328 * Generate a matrix with ZLATMS.
329 *
330  srnamt = 'ZLATMS'
331  CALL zlatms( n, n, dist, iseed, TYPE, RWORK, MODE,
332  $ cndnum, anorm, kl, ku, uplo, a, lda,
333  $ work, info )
334 *
335 * Check error code from ZLATMS and handle error.
336 *
337  IF( info.NE.0 ) THEN
338  CALL alaerh( path, 'ZLATMS', 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 zlatsy( 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 zlacpy( 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 = 'ZSYTRF_ROOK'
451  CALL zsytrf_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 ZSYTRF_ROOK and handle error.
472 *
473  IF( info.NE.k)
474  $ CALL alaerh( path, 'ZSYTRF_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 zsyt01_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 zlacpy( uplo, n, n, afac, lda, ainv, lda )
501  srnamt = 'ZSYTRI_ROOK'
502  CALL zsytri_rook( uplo, n, ainv, lda, iwork, work,
503  $ info )
504 *
505 * Check error code from ZSYTRI_ROOK and handle error.
506 *
507  IF( info.NE.0 )
508  $ CALL alaerh( path, 'ZSYTRI_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 zsyt03( 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  dtemp = 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  dtemp = zlange( '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  dtemp = zlange( 'M', k-2, 2,
565  $ afac( ( k-2 )*lda+1 ), lda, rwork )
566  k = k - 1
567 *
568  END IF
569 *
570 * DTEMP should be bounded by CONST
571 *
572  dtemp = dtemp - const + thresh
573  IF( dtemp.GT.result( 3 ) )
574  $ result( 3 ) = dtemp
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  dtemp = zlange( '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  dtemp = zlange( 'M', n-k-1, 2,
603  $ afac( ( k-1 )*lda+k+2 ), lda, rwork )
604  k = k + 1
605 *
606  END IF
607 *
608 * DTEMP should be bounded by CONST
609 *
610  dtemp = dtemp - const + thresh
611  IF( dtemp.GT.result( 3 ) )
612  $ result( 3 ) = dtemp
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 of 2-by-2 diag blocks
623 *
624  result( 4 ) = zero
625  dtemp = zero
626 *
627  const = ( ( alpha**2-one ) / ( alpha**2-onehalf ) )*
628  $ ( ( one + alpha ) / ( one - alpha ) )
629 *
630  IF( iuplo.EQ.1 ) THEN
631 *
632 * Loop backward for UPLO = 'U'
633 *
634  k = n
635  160 CONTINUE
636  IF( k.LE.1 )
637  $ go to 170
638 *
639  IF( iwork( k ).LT.zero ) THEN
640 *
641 * Get the two eigenvalues of a 2-by-2 block,
642 * store them in WORK array
643 *
644  block( 1, 1 ) = afac( ( k-2 )*lda+k-1 )
645  block( 2, 1 ) = afac( ( k-2 )*lda+k )
646  block( 1, 2 ) = block( 2, 1 )
647  block( 2, 2 ) = afac( (k-1)*lda+k )
648 *
649  CALL zgeevx( 'N', 'N', 'N', 'N', 2, block,
650  $ 2, work, zdummy, 1, zdummy, 1,
651  $ itemp, itemp2, rwork, dtemp,
652  $ rwork( 3 ), rwork( 5 ), work( 3 ),
653  $ 4, rwork( 7 ), info )
654 *
655  lam_max = max( abs( work( 1 ) ),
656  $ abs( work( 2 ) ) )
657  lam_min = min( abs( work( 1 ) ),
658  $ abs( work( 2 ) ) )
659 *
660  dtemp = lam_max / lam_min
661 *
662 * DTEMP should be bounded by CONST
663 *
664  dtemp = abs( dtemp ) - const + thresh
665  IF( dtemp.GT.result( 4 ) )
666  $ result( 4 ) = dtemp
667  k = k - 1
668 *
669  END IF
670 *
671  k = k - 1
672 *
673  go to 160
674  170 CONTINUE
675 *
676  ELSE
677 *
678 * Loop forward for UPLO = 'L'
679 *
680  k = 1
681  180 CONTINUE
682  IF( k.GE.n )
683  $ go to 190
684 *
685  IF( iwork( k ).LT.zero ) THEN
686 *
687 * Get the two eigenvalues of a 2-by-2 block,
688 * store them in WORK 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 zgeevx( 'N', 'N', 'N', 'N', 2, block,
696  $ 2, work, zdummy, 1, zdummy, 1,
697  $ itemp, itemp2, rwork, dtemp,
698  $ rwork( 3 ), rwork( 5 ), work( 3 ),
699  $ 4, rwork( 7 ), info )
700 *
701  lam_max = max( abs( work( 1 ) ),
702  $ abs( work( 2 ) ) )
703  lam_min = min( abs( work( 1 ) ),
704  $ abs( work( 2 ) ) )
705 *
706  dtemp = lam_max / lam_min
707 *
708 * DTEMP should be bounded by CONST
709 *
710  dtemp = abs( dtemp ) - const + thresh
711  IF( dtemp.GT.result( 4 ) )
712  $ result( 4 ) = dtemp
713  k = k + 1
714 *
715  END IF
716 *
717  k = k + 1
718 *
719  go to 180
720  190 CONTINUE
721  END IF
722 *
723 * Print information about the tests that did not pass
724 * the threshold.
725 *
726  DO 200 k = 3, 4
727  IF( result( k ).GE.thresh ) THEN
728  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
729  $ CALL alahd( nout, path )
730  WRITE( nout, fmt = 9999 )uplo, n, nb, imat, k,
731  $ result( k )
732  nfail = nfail + 1
733  END IF
734  200 CONTINUE
735  nrun = nrun + 2
736 *
737 * Skip the other tests if this is not the first block
738 * size.
739 *
740  IF( inb.GT.1 )
741  $ go to 240
742 *
743 * Do only the condition estimate if INFO is not 0.
744 *
745  IF( trfcon ) THEN
746  rcondc = zero
747  go to 230
748  END IF
749 *
750 * Do for each value of NRHS in NSVAL.
751 *
752  DO 220 irhs = 1, nns
753  nrhs = nsval( irhs )
754 *
755 *+ TEST 5 ( Using TRS_ROOK)
756 * Solve and compute residual for A * X = B.
757 *
758 * Choose a set of NRHS random solution vectors
759 * stored in XACT and set up the right hand side B
760 *
761  srnamt = 'ZLARHS'
762  CALL zlarhs( matpath, xtype, uplo, ' ', n, n,
763  $ kl, ku, nrhs, a, lda, xact, lda,
764  $ b, lda, iseed, info )
765  CALL zlacpy( 'Full', n, nrhs, b, lda, x, lda )
766 *
767  srnamt = 'ZSYTRS_ROOK'
768  CALL zsytrs_rook( uplo, n, nrhs, afac, lda, iwork,
769  $ x, lda, info )
770 *
771 * Check error code from ZSYTRS_ROOK and handle error.
772 *
773  IF( info.NE.0 )
774  $ CALL alaerh( path, 'ZSYTRS_ROOK', info, 0,
775  $ uplo, n, n, -1, -1, nrhs, imat,
776  $ nfail, nerrs, nout )
777 *
778  CALL zlacpy( 'Full', n, nrhs, b, lda, work, lda )
779 *
780 * Compute the residual for the solution
781 *
782  CALL zsyt02( uplo, n, nrhs, a, lda, x, lda, work,
783  $ lda, rwork, result( 5 ) )
784 *
785 *+ TEST 6
786 * Check solution from generated exact solution.
787 *
788  CALL zget04( n, nrhs, x, lda, xact, lda, rcondc,
789  $ result( 6 ) )
790 *
791 * Print information about the tests that did not pass
792 * the threshold.
793 *
794  DO 210 k = 5, 6
795  IF( result( k ).GE.thresh ) THEN
796  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
797  $ CALL alahd( nout, path )
798  WRITE( nout, fmt = 9998 )uplo, n, nrhs,
799  $ imat, k, result( k )
800  nfail = nfail + 1
801  END IF
802  210 CONTINUE
803  nrun = nrun + 2
804 *
805 * End do for each value of NRHS in NSVAL.
806 *
807  220 CONTINUE
808 *
809 *+ TEST 7
810 * Get an estimate of RCOND = 1/CNDNUM.
811 *
812  230 CONTINUE
813  anorm = zlansy( '1', uplo, n, a, lda, rwork )
814  srnamt = 'ZSYCON_ROOK'
815  CALL zsycon_rook( uplo, n, afac, lda, iwork, anorm,
816  $ rcond, work, info )
817 *
818 * Check error code from ZSYCON_ROOK and handle error.
819 *
820  IF( info.NE.0 )
821  $ CALL alaerh( path, 'ZSYCON_ROOK', info, 0,
822  $ uplo, n, n, -1, -1, -1, imat,
823  $ nfail, nerrs, nout )
824 *
825 * Compute the test ratio to compare values of RCOND
826 *
827  result( 7 ) = dget06( rcond, rcondc )
828 *
829 * Print information about the tests that did not pass
830 * the threshold.
831 *
832  IF( result( 7 ).GE.thresh ) THEN
833  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
834  $ CALL alahd( nout, path )
835  WRITE( nout, fmt = 9997 )uplo, n, imat, 7,
836  $ result( 7 )
837  nfail = nfail + 1
838  END IF
839  nrun = nrun + 1
840  240 CONTINUE
841 *
842  250 CONTINUE
843  260 CONTINUE
844  270 CONTINUE
845 *
846 * Print a summary of the results.
847 *
848  CALL alasum( path, nout, nfail, nrun, nerrs )
849 *
850  9999 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ', NB =', i4, ', type ',
851  $ i2, ', test ', i2, ', ratio =', g12.5 )
852  9998 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ', NRHS=', i3, ', type ',
853  $ i2, ', test(', i2, ') =', g12.5 )
854  9997 FORMAT( ' UPLO = ''', a1, ''', N =', i5, ',', 10x, ' type ', i2,
855  $ ', test(', i2, ') =', g12.5 )
856  RETURN
857 *
858 * End of ZCHKSY_ROOK
859 *
860  END
subroutine zsyt01_rook(UPLO, N, A, LDA, AFAC, LDAFAC, IPIV, C, LDC, RWORK, RESID)
ZSYT01_ROOK
Definition: zsyt01_rook.f:125
subroutine alahd(IOUNIT, PATH)
ALAHD
Definition: alahd.f:94
subroutine zget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
ZGET04
Definition: zget04.f:103
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:104
subroutine zsyt02(UPLO, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
ZSYT02
Definition: zsyt02.f:127
subroutine zerrsy(PATH, NUNIT)
ZERRSY
Definition: zerrsy.f:56
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:147
subroutine zsytrf_rook(UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO)
ZSYTRF_ROOK
Definition: zsytrf_rook.f:209
subroutine zchksy_rook(DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, NOUT)
ZCHKSY_ROOK
Definition: zchksy_rook.f:171
subroutine zlarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
ZLARHS
Definition: zlarhs.f:209
subroutine alasum(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASUM
Definition: alasum.f:74
subroutine zlatsy(UPLO, N, X, LDX, ISEED)
ZLATSY
Definition: zlatsy.f:90
subroutine zsytrs_rook(UPLO, N, NRHS, A, LDA, IPIV, B, LDB, INFO)
ZSYTRS_ROOK
Definition: zsytrs_rook.f:136
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:82
subroutine zsycon_rook(UPLO, N, A, LDA, IPIV, ANORM, RCOND, WORK, INFO)
ZSYCON_ROOK
Definition: zsycon_rook.f:139
subroutine zgeevx(BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, WORK, LWORK, RWORK, INFO)
ZGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...
Definition: zgeevx.f:284
subroutine zsyt03(UPLO, N, A, LDA, AINV, LDAINV, WORK, LDWORK, RWORK, RCOND, RESID)
ZSYT03
Definition: zsyt03.f:126
subroutine zlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
ZLATMS
Definition: zlatms.f:332
subroutine zsytri_rook(UPLO, N, A, LDA, IPIV, WORK, INFO)
ZSYTRI_ROOK
Definition: zsytri_rook.f:130
subroutine zlatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
ZLATB4
Definition: zlatb4.f:121