1       SUBROUTINE ZCHKGE( DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NNS,
  2      $                   NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B,
  3      $                   X, XACT, WORK, RWORK, IWORK, NOUT )
  4 *
  5 *  -- LAPACK test routine (version 3.1.1) --
  6 *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
  7 *     January 2007
  8 *
  9 *     .. Scalar Arguments ..
 10       LOGICAL            TSTERR
 11       INTEGER            NM, NMAX, NN, NNB, NNS, NOUT
 12       DOUBLE PRECISION   THRESH
 13 *     ..
 14 *     .. Array Arguments ..
 15       LOGICAL            DOTYPE( * )
 16       INTEGER            IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
 17      $                   NVAL( * )
 18       DOUBLE PRECISION   RWORK( * )
 19       COMPLEX*16         A( * ), AFAC( * ), AINV( * ), B( * ),
 20      $                   WORK( * ), X( * ), XACT( * )
 21 *     ..
 22 *
 23 *  Purpose
 24 *  =======
 25 *
 26 *  ZCHKGE tests ZGETRF, -TRI, -TRS, -RFS, and -CON.
 27 *
 28 *  Arguments
 29 *  =========
 30 *
 31 *  DOTYPE  (input) LOGICAL array, dimension (NTYPES)
 32 *          The matrix types to be used for testing.  Matrices of type j
 33 *          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
 34 *          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
 35 *
 36 *  NM      (input) INTEGER
 37 *          The number of values of M contained in the vector MVAL.
 38 *
 39 *  MVAL    (input) INTEGER array, dimension (NM)
 40 *          The values of the matrix row dimension M.
 41 *
 42 *  NN      (input) INTEGER
 43 *          The number of values of N contained in the vector NVAL.
 44 *
 45 *  NVAL    (input) INTEGER array, dimension (NN)
 46 *          The values of the matrix column dimension N.
 47 *
 48 *  NNB     (input) INTEGER
 49 *          The number of values of NB contained in the vector NBVAL.
 50 *
 51 *  NBVAL   (input) INTEGER array, dimension (NBVAL)
 52 *          The values of the blocksize NB.
 53 *
 54 *  NNS     (input) INTEGER
 55 *          The number of values of NRHS contained in the vector NSVAL.
 56 *
 57 *  NSVAL   (input) INTEGER array, dimension (NNS)
 58 *          The values of the number of right hand sides NRHS.
 59 *
 60 *  NRHS    (input) INTEGER
 61 *          The number of right hand side vectors to be generated for
 62 *          each linear system.
 63 *
 64 *  THRESH  (input) DOUBLE PRECISION
 65 *          The threshold value for the test ratios.  A result is
 66 *          included in the output file if RESULT >= THRESH.  To have
 67 *          every test ratio printed, use THRESH = 0.
 68 *
 69 *  TSTERR  (input) LOGICAL
 70 *          Flag that indicates whether error exits are to be tested.
 71 *
 72 *  NMAX    (input) INTEGER
 73 *          The maximum value permitted for M or N, used in dimensioning
 74 *          the work arrays.
 75 *
 76 *  A       (workspace) COMPLEX*16 array, dimension (NMAX*NMAX)
 77 *
 78 *  AFAC    (workspace) COMPLEX*16 array, dimension (NMAX*NMAX)
 79 *
 80 *  AINV    (workspace) COMPLEX*16 array, dimension (NMAX*NMAX)
 81 *
 82 *  B       (workspace) COMPLEX*16 array, dimension (NMAX*NSMAX)
 83 *          where NSMAX is the largest entry in NSVAL.
 84 *
 85 *  X       (workspace) COMPLEX*16 array, dimension (NMAX*NSMAX)
 86 *
 87 *  XACT    (workspace) COMPLEX*16 array, dimension (NMAX*NSMAX)
 88 *
 89 *  WORK    (workspace) COMPLEX*16 array, dimension
 90 *                      (NMAX*max(3,NSMAX))
 91 *
 92 *  RWORK   (workspace) DOUBLE PRECISION array, dimension
 93 *                      (max(2*NMAX,2*NSMAX+NWORK))
 94 *
 95 *  IWORK   (workspace) INTEGER array, dimension (NMAX)
 96 *
 97 *  NOUT    (input) INTEGER
 98 *          The unit number for output.
 99 *
100 *  =====================================================================
101 *
102 *     .. Parameters ..
103       DOUBLE PRECISION   ONE, ZERO
104       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
105       INTEGER            NTYPES
106       PARAMETER          ( NTYPES = 11 )
107       INTEGER            NTESTS
108       PARAMETER          ( NTESTS = 8 )
109       INTEGER            NTRAN
110       PARAMETER          ( NTRAN = 3 )
111 *     ..
112 *     .. Local Scalars ..
113       LOGICAL            TRFCON, ZEROT
114       CHARACTER          DIST, NORM, TRANS, TYPE, XTYPE
115       CHARACTER*3        PATH
116       INTEGER            I, IM, IMAT, IN, INB, INFO, IOFF, IRHS, ITRAN,
117      $                   IZERO, K, KL, KU, LDA, LWORK, M, MODE, N, NB,
118      $                   NERRS, NFAIL, NIMAT, NRHS, NRUN, NT
119       DOUBLE PRECISION   AINVNM, ANORM, ANORMI, ANORMO, CNDNUM, DUMMY,
120      $                   RCOND, RCONDC, RCONDI, RCONDO
121 *     ..
122 *     .. Local Arrays ..
123       CHARACTER          TRANSS( NTRAN )
124       INTEGER            ISEED( 4 ), ISEEDY( 4 )
125       DOUBLE PRECISION   RESULT( NTESTS )
126 *     ..
127 *     .. External Functions ..
128       DOUBLE PRECISION   DGET06, ZLANGE
129       EXTERNAL           DGET06, ZLANGE
130 *     ..
131 *     .. External Subroutines ..
132       EXTERNAL           ALAERH, ALAHD, ALASUM, XLAENV, ZERRGE, ZGECON,
133      $                   ZGERFS, ZGET01, ZGET02, ZGET03, ZGET04, ZGET07,
134      $                   ZGETRF, ZGETRI, ZGETRS, ZLACPY, ZLARHS, ZLASET,
135      $                   ZLATB4, ZLATMS
136 *     ..
137 *     .. Intrinsic Functions ..
138       INTRINSIC          DCMPLXMAXMIN
139 *     ..
140 *     .. Scalars in Common ..
141       LOGICAL            LERR, OK
142       CHARACTER*32       SRNAMT
143       INTEGER            INFOT, NUNIT
144 *     ..
145 *     .. Common blocks ..
146       COMMON             / INFOC / INFOT, NUNIT, OK, LERR
147       COMMON             / SRNAMC / SRNAMT
148 *     ..
149 *     .. Data statements ..
150       DATA               ISEEDY / 1988198919901991 / ,
151      $                   TRANSS / 'N''T''C' /
152 *     ..
153 *     .. Executable Statements ..
154 *
155 *     Initialize constants and the random number seed.
156 *
157       PATH( 11 ) = 'Zomplex precision'
158       PATH( 23 ) = 'GE'
159       NRUN = 0
160       NFAIL = 0
161       NERRS = 0
162       DO 10 I = 14
163          ISEED( I ) = ISEEDY( I )
164    10 CONTINUE
165 *
166 *     Test the error exits
167 *
168       CALL XLAENV( 11 )
169       IF( TSTERR )
170      $   CALL ZERRGE( PATH, NOUT )
171       INFOT = 0
172       CALL XLAENV( 22 )
173 *
174 *     Do for each value of M in MVAL
175 *
176       DO 120 IM = 1, NM
177          M = MVAL( IM )
178          LDA = MAX1, M )
179 *
180 *        Do for each value of N in NVAL
181 *
182          DO 110 IN = 1, NN
183             N = NVAL( IN )
184             XTYPE = 'N'
185             NIMAT = NTYPES
186             IF( M.LE.0 .OR. N.LE.0 )
187      $         NIMAT = 1
188 *
189             DO 100 IMAT = 1, NIMAT
190 *
191 *              Do the tests only if DOTYPE( IMAT ) is true.
192 *
193                IF.NOT.DOTYPE( IMAT ) )
194      $            GO TO 100
195 *
196 *              Skip types 5, 6, or 7 if the matrix size is too small.
197 *
198                ZEROT = IMAT.GE.5 .AND. IMAT.LE.7
199                IF( ZEROT .AND. N.LT.IMAT-4 )
200      $            GO TO 100
201 *
202 *              Set up parameters with ZLATB4 and generate a test matrix
203 *              with ZLATMS.
204 *
205                CALL ZLATB4( PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE,
206      $                      CNDNUM, DIST )
207 *
208                SRNAMT = 'ZLATMS'
209                CALL ZLATMS( M, N, DIST, ISEED, TYPE, RWORK, MODE,
210      $                      CNDNUM, ANORM, KL, KU, 'No packing', A, LDA,
211      $                      WORK, INFO )
212 *
213 *              Check error code from ZLATMS.
214 *
215                IF( INFO.NE.0 ) THEN
216                   CALL ALAERH( PATH, 'ZLATMS', INFO, 0' ', M, N, -1,
217      $                         -1-1, IMAT, NFAIL, NERRS, NOUT )
218                   GO TO 100
219                END IF
220 *
221 *              For types 5-7, zero one or more columns of the matrix to
222 *              test that INFO is returned correctly.
223 *
224                IF( ZEROT ) THEN
225                   IF( IMAT.EQ.5 ) THEN
226                      IZERO = 1
227                   ELSE IF( IMAT.EQ.6 ) THEN
228                      IZERO = MIN( M, N )
229                   ELSE
230                      IZERO = MIN( M, N ) / 2 + 1
231                   END IF
232                   IOFF = ( IZERO-1 )*LDA
233                   IF( IMAT.LT.7 ) THEN
234                      DO 20 I = 1, M
235                         A( IOFF+I ) = ZERO
236    20                CONTINUE
237                   ELSE
238                      CALL ZLASET( 'Full', M, N-IZERO+1DCMPLX( ZERO ),
239      $                            DCMPLX( ZERO ), A( IOFF+1 ), LDA )
240                   END IF
241                ELSE
242                   IZERO = 0
243                END IF
244 *
245 *              These lines, if used in place of the calls in the DO 60
246 *              loop, cause the code to bomb on a Sun SPARCstation.
247 *
248 *               ANORMO = ZLANGE( 'O', M, N, A, LDA, RWORK )
249 *               ANORMI = ZLANGE( 'I', M, N, A, LDA, RWORK )
250 *
251 *              Do for each blocksize in NBVAL
252 *
253                DO 90 INB = 1, NNB
254                   NB = NBVAL( INB )
255                   CALL XLAENV( 1, NB )
256 *
257 *                 Compute the LU factorization of the matrix.
258 *
259                   CALL ZLACPY( 'Full', M, N, A, LDA, AFAC, LDA )
260                   SRNAMT = 'ZGETRF'
261                   CALL ZGETRF( M, N, AFAC, LDA, IWORK, INFO )
262 *
263 *                 Check error code from ZGETRF.
264 *
265                   IF( INFO.NE.IZERO )
266      $               CALL ALAERH( PATH, 'ZGETRF', INFO, IZERO, ' ', M,
267      $                            N, -1-1, NB, IMAT, NFAIL, NERRS,
268      $                            NOUT )
269                   TRFCON = .FALSE.
270 *
271 *+    TEST 1
272 *                 Reconstruct matrix from factors and compute residual.
273 *
274                   CALL ZLACPY( 'Full', M, N, AFAC, LDA, AINV, LDA )
275                   CALL ZGET01( M, N, A, LDA, AINV, LDA, IWORK, RWORK,
276      $                         RESULT1 ) )
277                   NT = 1
278 *
279 *+    TEST 2
280 *                 Form the inverse if the factorization was successful
281 *                 and compute the residual.
282 *
283                   IF( M.EQ..AND. INFO.EQ.0 ) THEN
284                      CALL ZLACPY( 'Full', N, N, AFAC, LDA, AINV, LDA )
285                      SRNAMT = 'ZGETRI'
286                      NRHS = NSVAL( 1 )
287                      LWORK = NMAX*MAX3, NRHS )
288                      CALL ZGETRI( N, AINV, LDA, IWORK, WORK, LWORK,
289      $                            INFO )
290 *
291 *                    Check error code from ZGETRI.
292 *
293                      IF( INFO.NE.0 )
294      $                  CALL ALAERH( PATH, 'ZGETRI', INFO, 0' ', N, N,
295      $                               -1-1, NB, IMAT, NFAIL, NERRS,
296      $                               NOUT )
297 *
298 *                    Compute the residual for the matrix times its
299 *                    inverse.  Also compute the 1-norm condition number
300 *                    of A.
301 *
302                      CALL ZGET03( N, A, LDA, AINV, LDA, WORK, LDA,
303      $                            RWORK, RCONDO, RESULT2 ) )
304                      ANORMO = ZLANGE( 'O', M, N, A, LDA, RWORK )
305 *
306 *                    Compute the infinity-norm condition number of A.
307 *
308                      ANORMI = ZLANGE( 'I', M, N, A, LDA, RWORK )
309                      AINVNM = ZLANGE( 'I', N, N, AINV, LDA, RWORK )
310                      IF( ANORMI.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN
311                         RCONDI = ONE
312                      ELSE
313                         RCONDI = ( ONE / ANORMI ) / AINVNM
314                      END IF
315                      NT = 2
316                   ELSE
317 *
318 *                    Do only the condition estimate if INFO > 0.
319 *
320                      TRFCON = .TRUE.
321                      ANORMO = ZLANGE( 'O', M, N, A, LDA, RWORK )
322                      ANORMI = ZLANGE( 'I', M, N, A, LDA, RWORK )
323                      RCONDO = ZERO
324                      RCONDI = ZERO
325                   END IF
326 *
327 *                 Print information about the tests so far that did not
328 *                 pass the threshold.
329 *
330                   DO 30 K = 1, NT
331                      IFRESULT( K ).GE.THRESH ) THEN
332                         IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
333      $                     CALL ALAHD( NOUT, PATH )
334                         WRITE( NOUT, FMT = 9999 )M, N, NB, IMAT, K,
335      $                     RESULT( K )
336                         NFAIL = NFAIL + 1
337                      END IF
338    30             CONTINUE
339                   NRUN = NRUN + NT
340 *
341 *                 Skip the remaining tests if this is not the first
342 *                 block size or if M .ne. N.  Skip the solve tests if
343 *                 the matrix is singular.
344 *
345                   IF( INB.GT.1 .OR. M.NE.N )
346      $               GO TO 90
347                   IF( TRFCON )
348      $               GO TO 70
349 *
350                   DO 60 IRHS = 1, NNS
351                      NRHS = NSVAL( IRHS )
352                      XTYPE = 'N'
353 *
354                      DO 50 ITRAN = 1, NTRAN
355                         TRANS = TRANSS( ITRAN )
356                         IF( ITRAN.EQ.1 ) THEN
357                            RCONDC = RCONDO
358                         ELSE
359                            RCONDC = RCONDI
360                         END IF
361 *
362 *+    TEST 3
363 *                       Solve and compute residual for A * X = B.
364 *
365                         SRNAMT = 'ZLARHS'
366                         CALL ZLARHS( PATH, XTYPE, ' ', TRANS, N, N, KL,
367      $                               KU, NRHS, A, LDA, XACT, LDA, B,
368      $                               LDA, ISEED, INFO )
369                         XTYPE = 'C'
370 *
371                         CALL ZLACPY( 'Full', N, NRHS, B, LDA, X, LDA )
372                         SRNAMT = 'ZGETRS'
373                         CALL ZGETRS( TRANS, N, NRHS, AFAC, LDA, IWORK,
374      $                               X, LDA, INFO )
375 *
376 *                       Check error code from ZGETRS.
377 *
378                         IF( INFO.NE.0 )
379      $                     CALL ALAERH( PATH, 'ZGETRS', INFO, 0, TRANS,
380      $                                  N, N, -1-1, NRHS, IMAT, NFAIL,
381      $                                  NERRS, NOUT )
382 *
383                         CALL ZLACPY( 'Full', N, NRHS, B, LDA, WORK,
384      $                               LDA )
385                         CALL ZGET02( TRANS, N, N, NRHS, A, LDA, X, LDA,
386      $                               WORK, LDA, RWORK, RESULT3 ) )
387 *
388 *+    TEST 4
389 *                       Check solution from generated exact solution.
390 *
391                         CALL ZGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
392      $                               RESULT4 ) )
393 *
394 *+    TESTS 5, 6, and 7
395 *                       Use iterative refinement to improve the
396 *                       solution.
397 *
398                         SRNAMT = 'ZGERFS'
399                         CALL ZGERFS( TRANS, N, NRHS, A, LDA, AFAC, LDA,
400      $                               IWORK, B, LDA, X, LDA, RWORK,
401      $                               RWORK( NRHS+1 ), WORK,
402      $                               RWORK( 2*NRHS+1 ), INFO )
403 *
404 *                       Check error code from ZGERFS.
405 *
406                         IF( INFO.NE.0 )
407      $                     CALL ALAERH( PATH, 'ZGERFS', INFO, 0, TRANS,
408      $                                  N, N, -1-1, NRHS, IMAT, NFAIL,
409      $                                  NERRS, NOUT )
410 *
411                         CALL ZGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC,
412      $                               RESULT5 ) )
413                         CALL ZGET07( TRANS, N, NRHS, A, LDA, B, LDA, X,
414      $                               LDA, XACT, LDA, RWORK, .TRUE.,
415      $                               RWORK( NRHS+1 ), RESULT6 ) )
416 *
417 *                       Print information about the tests that did not
418 *                       pass the threshold.
419 *
420                         DO 40 K = 37
421                            IFRESULT( K ).GE.THRESH ) THEN
422                               IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
423      $                           CALL ALAHD( NOUT, PATH )
424                               WRITE( NOUT, FMT = 9998 )TRANS, N, NRHS,
425      $                           IMAT, K, RESULT( K )
426                               NFAIL = NFAIL + 1
427                            END IF
428    40                   CONTINUE
429                         NRUN = NRUN + 5
430    50                CONTINUE
431    60             CONTINUE
432 *
433 *+    TEST 8
434 *                    Get an estimate of RCOND = 1/CNDNUM.
435 *
436    70             CONTINUE
437                   DO 80 ITRAN = 12
438                      IF( ITRAN.EQ.1 ) THEN
439                         ANORM = ANORMO
440                         RCONDC = RCONDO
441                         NORM = 'O'
442                      ELSE
443                         ANORM = ANORMI
444                         RCONDC = RCONDI
445                         NORM = 'I'
446                      END IF
447                      SRNAMT = 'ZGECON'
448                      CALL ZGECON( NORM, N, AFAC, LDA, ANORM, RCOND,
449      $                            WORK, RWORK, INFO )
450 *
451 *                       Check error code from ZGECON.
452 *
453                      IF( INFO.NE.0 )
454      $                  CALL ALAERH( PATH, 'ZGECON', INFO, 0, NORM, N,
455      $                               N, -1-1-1, IMAT, NFAIL, NERRS,
456      $                               NOUT )
457 *
458 *                       This line is needed on a Sun SPARCstation.
459 *
460                      DUMMY = RCOND
461 *
462                      RESULT8 ) = DGET06( RCOND, RCONDC )
463 *
464 *                    Print information about the tests that did not pass
465 *                    the threshold.
466 *
467                      IFRESULT8 ).GE.THRESH ) THEN
468                         IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
469      $                     CALL ALAHD( NOUT, PATH )
470                         WRITE( NOUT, FMT = 9997 )NORM, N, IMAT, 8,
471      $                     RESULT8 )
472                         NFAIL = NFAIL + 1
473                      END IF
474                      NRUN = NRUN + 1
475    80             CONTINUE
476    90          CONTINUE
477   100       CONTINUE
478 *
479   110    CONTINUE
480   120 CONTINUE
481 *
482 *     Print a summary of the results.
483 *
484       CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS )
485 *
486  9999 FORMAT' M = ', I5, ', N =', I5, ', NB =', I4, ', type ', I2,
487      $      ', test(', I2, ') ='G12.5 )
488  9998 FORMAT' TRANS=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ',
489      $      I2, ', test(', I2, ') ='G12.5 )
490  9997 FORMAT' NORM =''', A1, ''', N =', I5, ','10X' type ', I2,
491      $      ', test(', I2, ') ='G12.5 )
492       RETURN
493 *
494 *     End of ZCHKGE
495 *
496       END