1*> \brief \b ZCHKHE_AA
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 ZCHKHE_AA( 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            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*> ZCHKHE_AA tests ZHETRF_AA, -TRS_AA.
35*> \endverbatim
36*
37*  Arguments:
38*  ==========
39*
40*> \param[in] DOTYPE
41*> \verbatim
42*>          DOTYPE is LOGICAL array, dimension (NTYPES)
43*>          The matrix types to be used for testing.  Matrices of type j
44*>          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
45*>          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
46*> \endverbatim
47*>
48*> \param[in] NN
49*> \verbatim
50*>          NN is INTEGER
51*>          The number of values of N contained in the vector NVAL.
52*> \endverbatim
53*>
54*> \param[in] NVAL
55*> \verbatim
56*>          NVAL is INTEGER array, dimension (NN)
57*>          The values of the matrix dimension N.
58*> \endverbatim
59*>
60*> \param[in] NNB
61*> \verbatim
62*>          NNB is INTEGER
63*>          The number of values of NB contained in the vector NBVAL.
64*> \endverbatim
65*>
66*> \param[in] NBVAL
67*> \verbatim
68*>          NBVAL is INTEGER array, dimension (NNB)
69*>          The values of the blocksize NB.
70*> \endverbatim
71*>
72*> \param[in] NNS
73*> \verbatim
74*>          NNS is INTEGER
75*>          The number of values of NRHS contained in the vector NSVAL.
76*> \endverbatim
77*>
78*> \param[in] NSVAL
79*> \verbatim
80*>          NSVAL is INTEGER array, dimension (NNS)
81*>          The values of the number of right hand sides NRHS.
82*> \endverbatim
83*>
84*> \param[in] THRESH
85*> \verbatim
86*>          THRESH is DOUBLE PRECISION
87*>          The threshold value for the test ratios.  A result is
88*>          included in the output file if RESULT >= THRESH.  To have
89*>          every test ratio printed, use THRESH = 0.
90*> \endverbatim
91*>
92*> \param[in] TSTERR
93*> \verbatim
94*>          TSTERR is LOGICAL
95*>          Flag that indicates whether error exits are to be tested.
96*> \endverbatim
97*>
98*> \param[in] NMAX
99*> \verbatim
100*>          NMAX is INTEGER
101*>          The maximum value permitted for N, used in dimensioning the
102*>          work arrays.
103*> \endverbatim
104*>
105*> \param[out] A
106*> \verbatim
107*>          A is COMPLEX*16 array, dimension (NMAX*NMAX)
108*> \endverbatim
109*>
110*> \param[out] AFAC
111*> \verbatim
112*>          AFAC is COMPLEX*16 array, dimension (NMAX*NMAX)
113*> \endverbatim
114*>
115*> \param[out] AINV
116*> \verbatim
117*>          AINV is COMPLEX*16 array, dimension (NMAX*NMAX)
118*> \endverbatim
119*>
120*> \param[out] B
121*> \verbatim
122*>          B is COMPLEX*16 array, dimension (NMAX*NSMAX)
123*>          where NSMAX is the largest entry in NSVAL.
124*> \endverbatim
125*>
126*> \param[out] X
127*> \verbatim
128*>          X is COMPLEX*16 array, dimension (NMAX*NSMAX)
129*> \endverbatim
130*>
131*> \param[out] XACT
132*> \verbatim
133*>          XACT is COMPLEX*16 array, dimension (NMAX*NSMAX)
134*> \endverbatim
135*>
136*> \param[out] WORK
137*> \verbatim
138*>          WORK is COMPLEX*16 array, dimension (NMAX*max(3,NSMAX))
139*> \endverbatim
140*>
141*> \param[out] RWORK
142*> \verbatim
143*>          RWORK is DOUBLE PRECISION array, dimension (max(NMAX,2*NSMAX))
144*> \endverbatim
145*>
146*> \param[out] IWORK
147*> \verbatim
148*>          IWORK is INTEGER array, dimension (NMAX)
149*> \endverbatim
150*>
151*> \param[in] NOUT
152*> \verbatim
153*>          NOUT is INTEGER
154*>          The unit number for output.
155*> \endverbatim
156*
157*  Authors:
158*  ========
159*
160*> \author Univ. of Tennessee
161*> \author Univ. of California Berkeley
162*> \author Univ. of Colorado Denver
163*> \author NAG Ltd.
164*
165*> \ingroup complex16_lin
166*
167*  =====================================================================
168      SUBROUTINE ZCHKHE_AA( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL,
169     $                      THRESH, TSTERR, NMAX, A, AFAC, AINV, B,
170     $                      X, XACT, WORK, RWORK, IWORK, NOUT )
171*
172*  -- LAPACK test routine --
173*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
174*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
175*
176      IMPLICIT NONE
177*
178*     .. Scalar Arguments ..
179      LOGICAL            TSTERR
180      INTEGER            NMAX, NN, NNB, NNS, NOUT
181      DOUBLE PRECISION   THRESH
182*     ..
183*     .. Array Arguments ..
184      LOGICAL            DOTYPE( * )
185      INTEGER            IWORK( * ), NBVAL( * ), NSVAL( * ), NVAL( * )
186      DOUBLE PRECISION   RWORK( * )
187      COMPLEX*16         A( * ), AFAC( * ), AINV( * ), B( * ),
188     $                   WORK( * ), X( * ), XACT( * )
189*     ..
190*
191*  =====================================================================
192*
193*     .. Parameters ..
194      DOUBLE PRECISION   ZERO
195      PARAMETER          ( ZERO = 0.0D+0 )
196      COMPLEX*16         CZERO
197      PARAMETER          ( CZERO = ( 0.0D+0, 0.0D+0 ) )
198      INTEGER            NTYPES
199      PARAMETER          ( NTYPES = 10 )
200      INTEGER            NTESTS
201      PARAMETER          ( NTESTS = 9 )
202*     ..
203*     .. Local Scalars ..
204      LOGICAL            ZEROT
205      CHARACTER          DIST, TYPE, UPLO, XTYPE
206      CHARACTER*3        PATH, MATPATH
207      INTEGER            I, I1, I2, IMAT, IN, INB, INFO, IOFF, IRHS,
208     $                   IUPLO, IZERO, J, K, KL, KU, LDA, LWORK, MODE,
209     $                   N, NB, NERRS, NFAIL, NIMAT, NRHS, NRUN, NT
210      DOUBLE PRECISION   ANORM, CNDNUM
211*     ..
212*     .. Local Arrays ..
213      CHARACTER          UPLOS( 2 )
214      INTEGER            ISEED( 4 ), ISEEDY( 4 )
215      DOUBLE PRECISION   RESULT( NTESTS )
216*     ..
217*     .. External Subroutines ..
218      EXTERNAL           ALAERH, ALAHD, ALASUM, XLAENV, ZERRHE,
219     $                   ZHET01_AA, ZHETRF_AA, ZHETRS_AA, ZLACPY,
220     $                   ZLAIPD, ZLARHS, ZLATB4, ZLATMS, ZPOT02
221*     ..
222*     .. Intrinsic Functions ..
223      INTRINSIC          MAX, MIN
224*     ..
225*     .. Scalars in Common ..
226      LOGICAL            LERR, OK
227      CHARACTER*32       SRNAMT
228      INTEGER            INFOT, NUNIT
229*     ..
230*     .. Common blocks ..
231      COMMON             / INFOC / INFOT, NUNIT, OK, LERR
232      COMMON             / SRNAMC / SRNAMT
233*     ..
234*     .. Data statements ..
235      DATA               ISEEDY / 1988, 1989, 1990, 1991 /
236      DATA               UPLOS / 'U', 'L' /
237*     ..
238*     .. Executable Statements ..
239*
240*     Initialize constants and the random number seed.
241*
242*     Test path
243*
244      PATH( 1: 1 ) = 'Zomplex precision'
245      PATH( 2: 3 ) = 'HA'
246*
247*     Path to generate matrices
248*
249      MATPATH( 1: 1 ) = 'Zomplex precision'
250      MATPATH( 2: 3 ) = 'HE'
251      NRUN = 0
252      NFAIL = 0
253      NERRS = 0
254      DO 10 I = 1, 4
255         ISEED( I ) = ISEEDY( I )
256   10 CONTINUE
257*
258*     Test the error exits
259*
260      IF( TSTERR )
261     $   CALL ZERRHE( PATH, NOUT )
262      INFOT = 0
263*
264*     Set the minimum block size for which the block routine should
265*     be used, which will be later returned by ILAENV
266*
267      CALL XLAENV( 2, 2 )
268*
269*     Do for each value of N in NVAL
270*
271      DO 180 IN = 1, NN
272         N = NVAL( IN )
273         IF( N .GT. NMAX ) THEN
274            NFAIL = NFAIL + 1
275            WRITE(NOUT, 9995) 'M ', N, NMAX
276            GO TO 180
277         END IF
278         LDA = MAX( N, 1 )
279         XTYPE = 'N'
280         NIMAT = NTYPES
281         IF( N.LE.0 )
282     $      NIMAT = 1
283*
284         IZERO = 0
285         DO 170 IMAT = 1, NIMAT
286*
287*           Do the tests only if DOTYPE( IMAT ) is true.
288*
289            IF( .NOT.DOTYPE( IMAT ) )
290     $         GO TO 170
291*
292*           Skip types 3, 4, 5, or 6 if the matrix size is too small.
293*
294            ZEROT = IMAT.GE.3 .AND. IMAT.LE.6
295            IF( ZEROT .AND. N.LT.IMAT-2 )
296     $         GO TO 170
297*
298*           Do first for UPLO = 'U', then for UPLO = 'L'
299*
300            DO 160 IUPLO = 1, 2
301               UPLO = UPLOS( IUPLO )
302*
303*              Set up parameters with ZLATB4 for the matrix generator
304*              based on the type of matrix to be generated.
305*
306               CALL ZLATB4( MATPATH, IMAT, N, N, TYPE, KL, KU,
307     $                       ANORM, MODE, CNDNUM, DIST )
308*
309*              Generate a matrix with ZLATMS.
310*
311               SRNAMT = 'ZLATMS'
312               CALL ZLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE,
313     $                      CNDNUM, ANORM, KL, KU, UPLO, A, LDA, WORK,
314     $                      INFO )
315*
316*              Check error code from ZLATMS and handle error.
317*
318               IF( INFO.NE.0 ) THEN
319                  CALL ALAERH( PATH, 'ZLATMS', INFO, 0, UPLO, N, N, -1,
320     $                         -1, -1, IMAT, NFAIL, NERRS, NOUT )
321*
322*                 Skip all tests for this generated matrix
323*
324                  GO TO 160
325               END IF
326*
327*              For types 3-6, zero one or more rows and columns of
328*              the matrix to test that INFO is returned correctly.
329*
330               IF( ZEROT ) THEN
331                  IF( IMAT.EQ.3 ) THEN
332                     IZERO = 1
333                  ELSE IF( IMAT.EQ.4 ) THEN
334                     IZERO = N
335                  ELSE
336                     IZERO = N / 2 + 1
337                  END IF
338*
339                  IF( IMAT.LT.6 ) THEN
340*
341*                    Set row and column IZERO to zero.
342*
343                     IF( IUPLO.EQ.1 ) THEN
344                        IOFF = ( IZERO-1 )*LDA
345                        DO 20 I = 1, IZERO - 1
346                           A( IOFF+I ) = CZERO
347   20                   CONTINUE
348                        IOFF = IOFF + IZERO
349                        DO 30 I = IZERO, N
350                           A( IOFF ) = CZERO
351                           IOFF = IOFF + LDA
352   30                   CONTINUE
353                     ELSE
354                        IOFF = IZERO
355                        DO 40 I = 1, IZERO - 1
356                           A( IOFF ) = CZERO
357                           IOFF = IOFF + LDA
358   40                   CONTINUE
359                        IOFF = IOFF - IZERO
360                        DO 50 I = IZERO, N
361                           A( IOFF+I ) = CZERO
362   50                   CONTINUE
363                     END IF
364                  ELSE
365                     IF( IUPLO.EQ.1 ) THEN
366*
367*                       Set the first IZERO rows and columns to zero.
368*
369                        IOFF = 0
370                        DO 70 J = 1, N
371                           I2 = MIN( J, IZERO )
372                           DO 60 I = 1, I2
373                              A( IOFF+I ) = CZERO
374   60                      CONTINUE
375                           IOFF = IOFF + LDA
376   70                   CONTINUE
377                        IZERO = 1
378                     ELSE
379*
380*                       Set the last IZERO rows and columns to zero.
381*
382                        IOFF = 0
383                        DO 90 J = 1, N
384                           I1 = MAX( J, IZERO )
385                           DO 80 I = I1, N
386                              A( IOFF+I ) = CZERO
387   80                      CONTINUE
388                           IOFF = IOFF + LDA
389   90                   CONTINUE
390                     END IF
391                  END IF
392               ELSE
393                  IZERO = 0
394               END IF
395*
396*              End generate test matrix A.
397*
398*
399*              Set the imaginary part of the diagonals.
400*
401               CALL ZLAIPD( N, A, LDA+1, 0 )
402*
403*              Do for each value of NB in NBVAL
404*
405               DO 150 INB = 1, NNB
406*
407*                 Set the optimal blocksize, which will be later
408*                 returned by ILAENV.
409*
410                  NB = NBVAL( INB )
411                  CALL XLAENV( 1, NB )
412*
413*                 Copy the test matrix A into matrix AFAC which
414*                 will be factorized in place. This is needed to
415*                 preserve the test matrix A for subsequent tests.
416*
417                  CALL ZLACPY( UPLO, N, N, A, LDA, AFAC, LDA )
418*
419*                 Compute the L*D*L**T or U*D*U**T factorization of the
420*                 matrix. IWORK stores details of the interchanges and
421*                 the block structure of D. AINV is a work array for
422*                 block factorization, LWORK is the length of AINV.
423*
424                  LWORK = MAX( 1, ( NB+1 )*LDA )
425                  SRNAMT = 'ZHETRF_AA'
426                  CALL ZHETRF_AA( UPLO, N, AFAC, LDA, IWORK, AINV,
427     $                            LWORK, INFO )
428*
429*                 Adjust the expected value of INFO to account for
430*                 pivoting.
431*
432c                  IF( IZERO.GT.0 ) THEN
433c                     J = 1
434c                     K = IZERO
435c  100                CONTINUE
436c                     IF( J.EQ.K ) THEN
437c                        K = IWORK( J )
438c                     ELSE IF( IWORK( J ).EQ.K ) THEN
439c                        K = J
440c                     END IF
441c                     IF( J.LT.K ) THEN
442c                        J = J + 1
443c                        GO TO 100
444c                     END IF
445c                  ELSE
446                     K = 0
447c                  END IF
448*
449*                 Check error code from ZHETRF and handle error.
450*
451                  IF( INFO.NE.K ) THEN
452                     CALL ALAERH( PATH, 'ZHETRF_AA', INFO, K, UPLO,
453     $                            N, N, -1, -1, NB, IMAT, NFAIL, NERRS,
454     $                            NOUT )
455                  END IF
456*
457*+    TEST 1
458*                 Reconstruct matrix from factors and compute residual.
459*
460                  CALL ZHET01_AA( UPLO, N, A, LDA, AFAC, LDA, IWORK,
461     $                            AINV, LDA, RWORK, RESULT( 1 ) )
462                  NT = 1
463*
464*
465*                 Print information about the tests that did not pass
466*                 the threshold.
467*
468                  DO 110 K = 1, NT
469                     IF( RESULT( K ).GE.THRESH ) THEN
470                        IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
471     $                     CALL ALAHD( NOUT, PATH )
472                        WRITE( NOUT, FMT = 9999 )UPLO, N, NB, IMAT, K,
473     $                     RESULT( K )
474                        NFAIL = NFAIL + 1
475                     END IF
476  110             CONTINUE
477                  NRUN = NRUN + NT
478*
479*                 Skip solver test if INFO is not 0.
480*
481                  IF( INFO.NE.0 ) THEN
482                     GO TO 140
483                  END IF
484*
485*                 Do for each value of NRHS in NSVAL.
486*
487                  DO 130 IRHS = 1, NNS
488                     NRHS = NSVAL( IRHS )
489*
490*+    TEST 2 (Using TRS)
491*                 Solve and compute residual for  A * X = B.
492*
493*                    Choose a set of NRHS random solution vectors
494*                    stored in XACT and set up the right hand side B
495*
496                     SRNAMT = 'ZLARHS'
497                     CALL ZLARHS( MATPATH, XTYPE, UPLO, ' ', N, N,
498     $                             KL, KU, NRHS, A, LDA, XACT, LDA,
499     $                            B, LDA, ISEED, INFO )
500                     CALL ZLACPY( 'Full', N, NRHS, B, LDA, X, LDA )
501*
502                     SRNAMT = 'ZHETRS_AA'
503                     LWORK = MAX( 1, 3*N-2 )
504                     CALL ZHETRS_AA( UPLO, N, NRHS, AFAC, LDA, IWORK,
505     $                               X, LDA, WORK, LWORK, INFO )
506*
507*                    Check error code from ZHETRS and handle error.
508*
509                     IF( INFO.NE.0 ) THEN
510                        IF( IZERO.EQ.0 ) THEN
511                           CALL ALAERH( PATH, 'ZHETRS_AA', INFO, 0,
512     $                                  UPLO, N, N, -1, -1, NRHS, IMAT,
513     $                                  NFAIL, NERRS, NOUT )
514                        END IF
515                     ELSE
516*
517                        CALL ZLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA
518     $                               )
519*
520*                       Compute the residual for the solution
521*
522                        CALL ZPOT02( UPLO, N, NRHS, A, LDA, X, LDA,
523     $                               WORK, LDA, RWORK, RESULT( 2 ) )
524*
525*                       Print information about the tests that did not pass
526*                       the threshold.
527*
528                        DO 120 K = 2, 2
529                           IF( RESULT( K ).GE.THRESH ) THEN
530                              IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
531     $                           CALL ALAHD( NOUT, PATH )
532                              WRITE( NOUT, FMT = 9998 )UPLO, N, NRHS,
533     $                           IMAT, K, RESULT( K )
534                              NFAIL = NFAIL + 1
535                           END IF
536  120                   CONTINUE
537                     END IF
538                     NRUN = NRUN + 1
539*
540*                 End do for each value of NRHS in NSVAL.
541*
542  130             CONTINUE
543  140             CONTINUE
544  150          CONTINUE
545  160       CONTINUE
546  170    CONTINUE
547  180 CONTINUE
548*
549*     Print a summary of the results.
550*
551      CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS )
552*
553 9999 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ', NB =', I4, ', type ',
554     $      I2, ', test ', I2, ', ratio =', G12.5 )
555 9998 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ', NRHS=', I3, ', type ',
556     $      I2, ', test(', I2, ') =', G12.5 )
557c 9997 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ',', 10X, ' type ', I2,
558c     $      ', test(', I2, ') =', G12.5 )
559 9995 FORMAT( ' Invalid input value: ', A4, '=', I6, '; must be <=',
560     $      I6 )
561      RETURN
562*
563*     End of ZCHKHE_AA
564*
565      END
566