LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zdrvsyx.f
Go to the documentation of this file.
1*> \brief \b ZDRVSYX
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 ZDRVSY( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
12* A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK,
13* NOUT )
14*
15* .. Scalar Arguments ..
16* LOGICAL TSTERR
17* INTEGER NMAX, NN, NOUT, NRHS
18* DOUBLE PRECISION THRESH
19* ..
20* .. Array Arguments ..
21* LOGICAL DOTYPE( * )
22* INTEGER IWORK( * ), 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*> ZDRVSY tests the driver routines ZSYSV, -SVX, and -SVXX.
35*>
36*> Note that this file is used only when the XBLAS are available,
37*> otherwise zdrvsy.f defines this subroutine.
38*> \endverbatim
39*
40* Arguments:
41* ==========
42*
43*> \param[in] DOTYPE
44*> \verbatim
45*> DOTYPE is LOGICAL array, dimension (NTYPES)
46*> The matrix types to be used for testing. Matrices of type j
47*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
48*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
49*> \endverbatim
50*>
51*> \param[in] NN
52*> \verbatim
53*> NN is INTEGER
54*> The number of values of N contained in the vector NVAL.
55*> \endverbatim
56*>
57*> \param[in] NVAL
58*> \verbatim
59*> NVAL is INTEGER array, dimension (NN)
60*> The values of the matrix dimension N.
61*> \endverbatim
62*>
63*> \param[in] NRHS
64*> \verbatim
65*> NRHS is INTEGER
66*> The number of right hand side vectors to be generated for
67*> each linear system.
68*> \endverbatim
69*>
70*> \param[in] THRESH
71*> \verbatim
72*> THRESH is DOUBLE PRECISION
73*> The threshold value for the test ratios. A result is
74*> included in the output file if RESULT >= THRESH. To have
75*> every test ratio printed, use THRESH = 0.
76*> \endverbatim
77*>
78*> \param[in] TSTERR
79*> \verbatim
80*> TSTERR is LOGICAL
81*> Flag that indicates whether error exits are to be tested.
82*> \endverbatim
83*>
84*> \param[in] NMAX
85*> \verbatim
86*> NMAX is INTEGER
87*> The maximum value permitted for N, used in dimensioning the
88*> work arrays.
89*> \endverbatim
90*>
91*> \param[out] A
92*> \verbatim
93*> A is COMPLEX*16 array, dimension (NMAX*NMAX)
94*> \endverbatim
95*>
96*> \param[out] AFAC
97*> \verbatim
98*> AFAC is COMPLEX*16 array, dimension (NMAX*NMAX)
99*> \endverbatim
100*>
101*> \param[out] AINV
102*> \verbatim
103*> AINV is COMPLEX*16 array, dimension (NMAX*NMAX)
104*> \endverbatim
105*>
106*> \param[out] B
107*> \verbatim
108*> B is COMPLEX*16 array, dimension (NMAX*NRHS)
109*> \endverbatim
110*>
111*> \param[out] X
112*> \verbatim
113*> X is COMPLEX*16 array, dimension (NMAX*NRHS)
114*> \endverbatim
115*>
116*> \param[out] XACT
117*> \verbatim
118*> XACT is COMPLEX*16 array, dimension (NMAX*NRHS)
119*> \endverbatim
120*>
121*> \param[out] WORK
122*> \verbatim
123*> WORK is COMPLEX*16 array, dimension
124*> (NMAX*max(2,NRHS))
125*> \endverbatim
126*>
127*> \param[out] RWORK
128*> \verbatim
129*> RWORK is DOUBLE PRECISION array, dimension (2*NMAX+2*NRHS)
130*> \endverbatim
131*>
132*> \param[out] IWORK
133*> \verbatim
134*> IWORK is INTEGER array, dimension (NMAX)
135*> \endverbatim
136*>
137*> \param[in] NOUT
138*> \verbatim
139*> NOUT is INTEGER
140*> The unit number for output.
141*> \endverbatim
142*
143* Authors:
144* ========
145*
146*> \author Univ. of Tennessee
147*> \author Univ. of California Berkeley
148*> \author Univ. of Colorado Denver
149*> \author NAG Ltd.
150*
151*> \ingroup complex16_lin
152*
153* =====================================================================
154 SUBROUTINE zdrvsy( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX,
155 $ A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK,
156 $ NOUT )
157*
158* -- LAPACK test routine --
159* -- LAPACK is a software package provided by Univ. of Tennessee, --
160* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
161*
162* .. Scalar Arguments ..
163 LOGICAL TSTERR
164 INTEGER NMAX, NN, NOUT, NRHS
165 DOUBLE PRECISION THRESH
166* ..
167* .. Array Arguments ..
168 LOGICAL DOTYPE( * )
169 INTEGER IWORK( * ), NVAL( * )
170 DOUBLE PRECISION RWORK( * )
171 COMPLEX*16 A( * ), AFAC( * ), AINV( * ), B( * ),
172 $ work( * ), x( * ), xact( * )
173* ..
174*
175* =====================================================================
176*
177* .. Parameters ..
178 DOUBLE PRECISION ONE, ZERO
179 PARAMETER ( ONE = 1.0d+0, zero = 0.0d+0 )
180 INTEGER NTYPES, NTESTS
181 parameter( ntypes = 11, ntests = 6 )
182 INTEGER NFACT
183 parameter( nfact = 2 )
184* ..
185* .. Local Scalars ..
186 LOGICAL ZEROT
187 CHARACTER DIST, EQUED, FACT, TYPE, UPLO, XTYPE
188 CHARACTER*3 PATH
189 INTEGER I, I1, I2, IFACT, IMAT, IN, INFO, IOFF, IUPLO,
190 $ izero, j, k, k1, kl, ku, lda, lwork, mode, n,
191 $ nb, nbmin, nerrs, nfail, nimat, nrun, nt,
192 $ n_err_bnds
193 DOUBLE PRECISION AINVNM, ANORM, CNDNUM, RCOND, RCONDC,
194 $ RPVGRW_SVXX
195* ..
196* .. Local Arrays ..
197 CHARACTER FACTS( NFACT ), UPLOS( 2 )
198 INTEGER ISEED( 4 ), ISEEDY( 4 )
199 DOUBLE PRECISION RESULT( NTESTS ), BERR( NRHS ),
200 $ errbnds_n( nrhs, 3 ), errbnds_c( nrhs, 3 )
201* ..
202* .. External Functions ..
203 DOUBLE PRECISION DGET06, ZLANSY
204 EXTERNAL DGET06, ZLANSY
205* ..
206* .. External Subroutines ..
207 EXTERNAL aladhd, alaerh, alasvm, xlaenv, zerrvx, zget04,
211* ..
212* .. Scalars in Common ..
213 LOGICAL LERR, OK
214 CHARACTER*32 SRNAMT
215 INTEGER INFOT, NUNIT
216* ..
217* .. Common blocks ..
218 COMMON / infoc / infot, nunit, ok, lerr
219 COMMON / srnamc / srnamt
220* ..
221* .. Intrinsic Functions ..
222 INTRINSIC dcmplx, max, min
223* ..
224* .. Data statements ..
225 DATA iseedy / 1988, 1989, 1990, 1991 /
226 DATA uplos / 'U', 'L' / , facts / 'F', 'N' /
227* ..
228* .. Executable Statements ..
229*
230* Initialize constants and the random number seed.
231*
232 path( 1: 1 ) = 'Zomplex precision'
233 path( 2: 3 ) = 'SY'
234 nrun = 0
235 nfail = 0
236 nerrs = 0
237 DO 10 i = 1, 4
238 iseed( i ) = iseedy( i )
239 10 CONTINUE
240 lwork = max( 2*nmax, nmax*nrhs )
241*
242* Test the error exits
243*
244 IF( tsterr )
245 $ CALL zerrvx( path, nout )
246 infot = 0
247*
248* Set the block size and minimum block size for testing.
249*
250 nb = 1
251 nbmin = 2
252 CALL xlaenv( 1, nb )
253 CALL xlaenv( 2, nbmin )
254*
255* Do for each value of N in NVAL
256*
257 DO 180 in = 1, nn
258 n = nval( in )
259 lda = max( n, 1 )
260 xtype = 'N'
261 nimat = ntypes
262 IF( n.LE.0 )
263 $ nimat = 1
264*
265 DO 170 imat = 1, nimat
266*
267* Do the tests only if DOTYPE( IMAT ) is true.
268*
269 IF( .NOT.dotype( imat ) )
270 $ GO TO 170
271*
272* Skip types 3, 4, 5, or 6 if the matrix size is too small.
273*
274 zerot = imat.GE.3 .AND. imat.LE.6
275 IF( zerot .AND. n.LT.imat-2 )
276 $ GO TO 170
277*
278* Do first for UPLO = 'U', then for UPLO = 'L'
279*
280 DO 160 iuplo = 1, 2
281 uplo = uplos( iuplo )
282*
283 IF( imat.NE.ntypes ) THEN
284*
285* Set up parameters with ZLATB4 and generate a test
286* matrix with ZLATMS.
287*
288 CALL zlatb4( path, imat, n, n, TYPE, kl, ku, anorm,
289 $ mode, cndnum, dist )
290*
291 srnamt = 'ZLATMS'
292 CALL zlatms( n, n, dist, iseed, TYPE, rwork, mode,
293 $ cndnum, anorm, kl, ku, uplo, a, lda,
294 $ work, info )
295*
296* Check error code from ZLATMS.
297*
298 IF( info.NE.0 ) THEN
299 CALL alaerh( path, 'ZLATMS', info, 0, uplo, n, n,
300 $ -1, -1, -1, imat, nfail, nerrs, nout )
301 GO TO 160
302 END IF
303*
304* For types 3-6, zero one or more rows and columns of
305* the matrix to test that INFO is returned correctly.
306*
307 IF( zerot ) THEN
308 IF( imat.EQ.3 ) THEN
309 izero = 1
310 ELSE IF( imat.EQ.4 ) THEN
311 izero = n
312 ELSE
313 izero = n / 2 + 1
314 END IF
315*
316 IF( imat.LT.6 ) THEN
317*
318* Set row and column IZERO to zero.
319*
320 IF( iuplo.EQ.1 ) THEN
321 ioff = ( izero-1 )*lda
322 DO 20 i = 1, izero - 1
323 a( ioff+i ) = zero
324 20 CONTINUE
325 ioff = ioff + izero
326 DO 30 i = izero, n
327 a( ioff ) = zero
328 ioff = ioff + lda
329 30 CONTINUE
330 ELSE
331 ioff = izero
332 DO 40 i = 1, izero - 1
333 a( ioff ) = zero
334 ioff = ioff + lda
335 40 CONTINUE
336 ioff = ioff - izero
337 DO 50 i = izero, n
338 a( ioff+i ) = zero
339 50 CONTINUE
340 END IF
341 ELSE
342 IF( iuplo.EQ.1 ) THEN
343*
344* Set the first IZERO rows to zero.
345*
346 ioff = 0
347 DO 70 j = 1, n
348 i2 = min( j, izero )
349 DO 60 i = 1, i2
350 a( ioff+i ) = zero
351 60 CONTINUE
352 ioff = ioff + lda
353 70 CONTINUE
354 ELSE
355*
356* Set the last IZERO rows to zero.
357*
358 ioff = 0
359 DO 90 j = 1, n
360 i1 = max( j, izero )
361 DO 80 i = i1, n
362 a( ioff+i ) = zero
363 80 CONTINUE
364 ioff = ioff + lda
365 90 CONTINUE
366 END IF
367 END IF
368 ELSE
369 izero = 0
370 END IF
371 ELSE
372*
373* IMAT = NTYPES: Use a special block diagonal matrix to
374* test alternate code for the 2-by-2 blocks.
375*
376 CALL zlatsy( uplo, n, a, lda, iseed )
377 END IF
378*
379 DO 150 ifact = 1, nfact
380*
381* Do first for FACT = 'F', then for other values.
382*
383 fact = facts( ifact )
384*
385* Compute the condition number for comparison with
386* the value returned by ZSYSVX.
387*
388 IF( zerot ) THEN
389 IF( ifact.EQ.1 )
390 $ GO TO 150
391 rcondc = zero
392*
393 ELSE IF( ifact.EQ.1 ) THEN
394*
395* Compute the 1-norm of A.
396*
397 anorm = zlansy( '1', uplo, n, a, lda, rwork )
398*
399* Factor the matrix A.
400*
401 CALL zlacpy( uplo, n, n, a, lda, afac, lda )
402 CALL zsytrf( uplo, n, afac, lda, iwork, work,
403 $ lwork, info )
404*
405* Compute inv(A) and take its norm.
406*
407 CALL zlacpy( uplo, n, n, afac, lda, ainv, lda )
408 lwork = (n+nb+1)*(nb+3)
409 CALL zsytri2( uplo, n, ainv, lda, iwork, work,
410 $ lwork, info )
411 ainvnm = zlansy( '1', uplo, n, ainv, lda, rwork )
412*
413* Compute the 1-norm condition number of A.
414*
415 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
416 rcondc = one
417 ELSE
418 rcondc = ( one / anorm ) / ainvnm
419 END IF
420 END IF
421*
422* Form an exact solution and set the right hand side.
423*
424 srnamt = 'ZLARHS'
425 CALL zlarhs( path, xtype, uplo, ' ', n, n, kl, ku,
426 $ nrhs, a, lda, xact, lda, b, lda, iseed,
427 $ info )
428 xtype = 'C'
429*
430* --- Test ZSYSV ---
431*
432 IF( ifact.EQ.2 ) THEN
433 CALL zlacpy( uplo, n, n, a, lda, afac, lda )
434 CALL zlacpy( 'Full', n, nrhs, b, lda, x, lda )
435*
436* Factor the matrix and solve the system using ZSYSV.
437*
438 srnamt = 'ZSYSV '
439 CALL zsysv( uplo, n, nrhs, afac, lda, iwork, x,
440 $ lda, work, lwork, info )
441*
442* Adjust the expected value of INFO to account for
443* pivoting.
444*
445 k = izero
446 IF( k.GT.0 ) THEN
447 100 CONTINUE
448 IF( iwork( k ).LT.0 ) THEN
449 IF( iwork( k ).NE.-k ) THEN
450 k = -iwork( k )
451 GO TO 100
452 END IF
453 ELSE IF( iwork( k ).NE.k ) THEN
454 k = iwork( k )
455 GO TO 100
456 END IF
457 END IF
458*
459* Check error code from ZSYSV .
460*
461 IF( info.NE.k ) THEN
462 CALL alaerh( path, 'ZSYSV ', info, k, uplo, n,
463 $ n, -1, -1, nrhs, imat, nfail,
464 $ nerrs, nout )
465 GO TO 120
466 ELSE IF( info.NE.0 ) THEN
467 GO TO 120
468 END IF
469*
470* Reconstruct matrix from factors and compute
471* residual.
472*
473 CALL zsyt01( uplo, n, a, lda, afac, lda, iwork,
474 $ ainv, lda, rwork, result( 1 ) )
475*
476* Compute residual of the computed solution.
477*
478 CALL zlacpy( 'Full', n, nrhs, b, lda, work, lda )
479 CALL zsyt02( uplo, n, nrhs, a, lda, x, lda, work,
480 $ lda, rwork, result( 2 ) )
481*
482* Check solution from generated exact solution.
483*
484 CALL zget04( n, nrhs, x, lda, xact, lda, rcondc,
485 $ result( 3 ) )
486 nt = 3
487*
488* Print information about the tests that did not pass
489* the threshold.
490*
491 DO 110 k = 1, nt
492 IF( result( k ).GE.thresh ) THEN
493 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
494 $ CALL aladhd( nout, path )
495 WRITE( nout, fmt = 9999 )'ZSYSV ', uplo, n,
496 $ imat, k, result( k )
497 nfail = nfail + 1
498 END IF
499 110 CONTINUE
500 nrun = nrun + nt
501 120 CONTINUE
502 END IF
503*
504* --- Test ZSYSVX ---
505*
506 IF( ifact.EQ.2 )
507 $ CALL zlaset( uplo, n, n, dcmplx( zero ),
508 $ dcmplx( zero ), afac, lda )
509 CALL zlaset( 'Full', n, nrhs, dcmplx( zero ),
510 $ dcmplx( zero ), x, lda )
511*
512* Solve the system and compute the condition number and
513* error bounds using ZSYSVX.
514*
515 srnamt = 'ZSYSVX'
516 CALL zsysvx( fact, uplo, n, nrhs, a, lda, afac, lda,
517 $ iwork, b, lda, x, lda, rcond, rwork,
518 $ rwork( nrhs+1 ), work, lwork,
519 $ rwork( 2*nrhs+1 ), info )
520*
521* Adjust the expected value of INFO to account for
522* pivoting.
523*
524 k = izero
525 IF( k.GT.0 ) THEN
526 130 CONTINUE
527 IF( iwork( k ).LT.0 ) THEN
528 IF( iwork( k ).NE.-k ) THEN
529 k = -iwork( k )
530 GO TO 130
531 END IF
532 ELSE IF( iwork( k ).NE.k ) THEN
533 k = iwork( k )
534 GO TO 130
535 END IF
536 END IF
537*
538* Check the error code from ZSYSVX.
539*
540 IF( info.NE.k ) THEN
541 CALL alaerh( path, 'ZSYSVX', info, k, fact // uplo,
542 $ n, n, -1, -1, nrhs, imat, nfail,
543 $ nerrs, nout )
544 GO TO 150
545 END IF
546*
547 IF( info.EQ.0 ) THEN
548 IF( ifact.GE.2 ) THEN
549*
550* Reconstruct matrix from factors and compute
551* residual.
552*
553 CALL zsyt01( uplo, n, a, lda, afac, lda, iwork,
554 $ ainv, lda, rwork( 2*nrhs+1 ),
555 $ result( 1 ) )
556 k1 = 1
557 ELSE
558 k1 = 2
559 END IF
560*
561* Compute residual of the computed solution.
562*
563 CALL zlacpy( 'Full', n, nrhs, b, lda, work, lda )
564 CALL zsyt02( uplo, n, nrhs, a, lda, x, lda, work,
565 $ lda, rwork( 2*nrhs+1 ), result( 2 ) )
566*
567* Check solution from generated exact solution.
568*
569 CALL zget04( n, nrhs, x, lda, xact, lda, rcondc,
570 $ result( 3 ) )
571*
572* Check the error bounds from iterative refinement.
573*
574 CALL zpot05( uplo, n, nrhs, a, lda, b, lda, x, lda,
575 $ xact, lda, rwork, rwork( nrhs+1 ),
576 $ result( 4 ) )
577 ELSE
578 k1 = 6
579 END IF
580*
581* Compare RCOND from ZSYSVX with the computed value
582* in RCONDC.
583*
584 result( 6 ) = dget06( rcond, rcondc )
585*
586* Print information about the tests that did not pass
587* the threshold.
588*
589 DO 140 k = k1, 6
590 IF( result( k ).GE.thresh ) THEN
591 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
592 $ CALL aladhd( nout, path )
593 WRITE( nout, fmt = 9998 )'ZSYSVX', fact, uplo,
594 $ n, imat, k, result( k )
595 nfail = nfail + 1
596 END IF
597 140 CONTINUE
598 nrun = nrun + 7 - k1
599*
600* --- Test ZSYSVXX ---
601*
602* Restore the matrices A and B.
603*
604 IF( ifact.EQ.2 )
605 $ CALL zlaset( uplo, n, n, dcmplx( zero ),
606 $ dcmplx( zero ), afac, lda )
607 CALL zlaset( 'Full', n, nrhs, dcmplx( zero ),
608 $ dcmplx( zero ), x, lda )
609*
610* Solve the system and compute the condition number
611* and error bounds using ZSYSVXX.
612*
613 srnamt = 'ZSYSVXX'
614 n_err_bnds = 3
615 equed = 'N'
616 CALL zsysvxx( fact, uplo, n, nrhs, a, lda, afac,
617 $ lda, iwork, equed, work( n+1 ), b, lda, x,
618 $ lda, rcond, rpvgrw_svxx, berr, n_err_bnds,
619 $ errbnds_n, errbnds_c, 0, zero, work,
620 $ rwork, info )
621*
622* Adjust the expected value of INFO to account for
623* pivoting.
624*
625 k = izero
626 IF( k.GT.0 ) THEN
627 135 CONTINUE
628 IF( iwork( k ).LT.0 ) THEN
629 IF( iwork( k ).NE.-k ) THEN
630 k = -iwork( k )
631 GO TO 135
632 END IF
633 ELSE IF( iwork( k ).NE.k ) THEN
634 k = iwork( k )
635 GO TO 135
636 END IF
637 END IF
638*
639* Check the error code from ZSYSVXX.
640*
641 IF( info.NE.k .AND. info.LE.n ) THEN
642 CALL alaerh( path, 'ZSYSVXX', info, k,
643 $ fact // uplo, n, n, -1, -1, nrhs, imat, nfail,
644 $ nerrs, nout )
645 GO TO 150
646 END IF
647*
648 IF( info.EQ.0 ) THEN
649 IF( ifact.GE.2 ) THEN
650*
651* Reconstruct matrix from factors and compute
652* residual.
653*
654 CALL zsyt01( uplo, n, a, lda, afac, lda, iwork,
655 $ ainv, lda, rwork(2*nrhs+1),
656 $ result( 1 ) )
657 k1 = 1
658 ELSE
659 k1 = 2
660 END IF
661*
662* Compute residual of the computed solution.
663*
664 CALL zlacpy( 'Full', n, nrhs, b, lda, work, lda )
665 CALL zsyt02( uplo, n, nrhs, a, lda, x, lda, work,
666 $ lda, rwork( 2*nrhs+1 ), result( 2 ) )
667 result( 2 ) = 0.0
668*
669* Check solution from generated exact solution.
670*
671 CALL zget04( n, nrhs, x, lda, xact, lda, rcondc,
672 $ result( 3 ) )
673*
674* Check the error bounds from iterative refinement.
675*
676 CALL zpot05( uplo, n, nrhs, a, lda, b, lda, x, lda,
677 $ xact, lda, rwork, rwork( nrhs+1 ),
678 $ result( 4 ) )
679 ELSE
680 k1 = 6
681 END IF
682*
683* Compare RCOND from ZSYSVXX with the computed value
684* in RCONDC.
685*
686 result( 6 ) = dget06( rcond, rcondc )
687*
688* Print information about the tests that did not pass
689* the threshold.
690*
691 DO 85 k = k1, 6
692 IF( result( k ).GE.thresh ) THEN
693 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
694 $ CALL aladhd( nout, path )
695 WRITE( nout, fmt = 9998 )'ZSYSVXX',
696 $ fact, uplo, n, imat, k,
697 $ result( k )
698 nfail = nfail + 1
699 END IF
700 85 CONTINUE
701 nrun = nrun + 7 - k1
702*
703 150 CONTINUE
704*
705 160 CONTINUE
706 170 CONTINUE
707 180 CONTINUE
708*
709* Print a summary of the results.
710*
711 CALL alasvm( path, nout, nfail, nrun, nerrs )
712*
713
714* Test Error Bounds from ZSYSVXX
715
716 CALL zebchvxx(thresh, path)
717
718 9999 FORMAT( 1x, a, ', UPLO=''', a1, ''', N =', i5, ', type ', i2,
719 $ ', test ', i2, ', ratio =', g12.5 )
720 9998 FORMAT( 1x, a, ', FACT=''', a1, ''', UPLO=''', a1, ''', N =', i5,
721 $ ', type ', i2, ', test ', i2, ', ratio =', g12.5 )
722 RETURN
723*
724* End of ZDRVSYX
725*
726 END
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
Definition alasvm.f:73
subroutine xlaenv(ispec, nvalue)
XLAENV
Definition xlaenv.f:81
subroutine zlarhs(path, xtype, uplo, trans, m, n, kl, ku, nrhs, a, lda, x, ldx, b, ldb, iseed, info)
ZLARHS
Definition zlarhs.f:208
subroutine aladhd(iounit, path)
ALADHD
Definition aladhd.f:90
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
Definition alaerh.f:147
subroutine zsysv(uplo, n, nrhs, a, lda, ipiv, b, ldb, work, lwork, info)
ZSYSV computes the solution to system of linear equations A * X = B for SY matrices
Definition zsysv.f:171
subroutine zsysvx(fact, uplo, n, nrhs, a, lda, af, ldaf, ipiv, b, ldb, x, ldx, rcond, ferr, berr, work, lwork, rwork, info)
ZSYSVX computes the solution to system of linear equations A * X = B for SY matrices
Definition zsysvx.f:285
subroutine zsysvxx(fact, uplo, n, nrhs, a, lda, af, ldaf, ipiv, equed, s, b, ldb, x, ldx, rcond, rpvgrw, berr, n_err_bnds, err_bnds_norm, err_bnds_comp, nparams, params, work, rwork, info)
ZSYSVXX computes the solution to system of linear equations A * X = B for SY matrices
Definition zsysvxx.f:506
subroutine zsytrf(uplo, n, a, lda, ipiv, work, lwork, info)
ZSYTRF
Definition zsytrf.f:182
subroutine zsytri2(uplo, n, a, lda, ipiv, work, lwork, info)
ZSYTRI2
Definition zsytri2.f:127
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
subroutine zlaset(uplo, m, n, alpha, beta, a, lda)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition zlaset.f:106
subroutine zdrvsy(dotype, nn, nval, nrhs, thresh, tsterr, nmax, a, afac, ainv, b, x, xact, work, rwork, iwork, nout)
ZDRVSY
Definition zdrvsy.f:153
subroutine zebchvxx(thresh, path)
ZEBCHVXX
Definition zebchvxx.f:96
subroutine zerrvx(path, nunit)
ZERRVX
Definition zerrvx.f:55
subroutine zget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
ZGET04
Definition zget04.f:102
subroutine zlatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
ZLATB4
Definition zlatb4.f:121
subroutine zlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
ZLATMS
Definition zlatms.f:332
subroutine zlatsy(uplo, n, x, ldx, iseed)
ZLATSY
Definition zlatsy.f:89
subroutine zpot05(uplo, n, nrhs, a, lda, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
ZPOT05
Definition zpot05.f:165
subroutine zsyt01(uplo, n, a, lda, afac, ldafac, ipiv, c, ldc, rwork, resid)
ZSYT01
Definition zsyt01.f:125
subroutine zsyt02(uplo, n, nrhs, a, lda, x, ldx, b, ldb, rwork, resid)
ZSYT02
Definition zsyt02.f:127