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