LAPACK  3.10.1
LAPACK: Linear Algebra PACKage

◆ dchkge()

subroutine dchkge ( logical, dimension( * )  DOTYPE,
integer  NM,
integer, dimension( * )  MVAL,
integer  NN,
integer, dimension( * )  NVAL,
integer  NNB,
integer, dimension( * )  NBVAL,
integer  NNS,
integer, dimension( * )  NSVAL,
double precision  THRESH,
logical  TSTERR,
integer  NMAX,
double precision, dimension( * )  A,
double precision, dimension( * )  AFAC,
double precision, dimension( * )  AINV,
double precision, dimension( * )  B,
double precision, dimension( * )  X,
double precision, dimension( * )  XACT,
double precision, dimension( * )  WORK,
double precision, dimension( * )  RWORK,
integer, dimension( * )  IWORK,
integer  NOUT 
)

DCHKGE

Purpose:
 DCHKGE tests DGETRF, -TRI, -TRS, -RFS, and -CON.
Parameters
[in]DOTYPE
          DOTYPE is LOGICAL array, dimension (NTYPES)
          The matrix types to be used for testing.  Matrices of type j
          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
[in]NM
          NM is INTEGER
          The number of values of M contained in the vector MVAL.
[in]MVAL
          MVAL is INTEGER array, dimension (NM)
          The values of the matrix row dimension M.
[in]NN
          NN is INTEGER
          The number of values of N contained in the vector NVAL.
[in]NVAL
          NVAL is INTEGER array, dimension (NN)
          The values of the matrix column dimension N.
[in]NNB
          NNB is INTEGER
          The number of values of NB contained in the vector NBVAL.
[in]NBVAL
          NBVAL is INTEGER array, dimension (NNB)
          The values of the blocksize NB.
[in]NNS
          NNS is INTEGER
          The number of values of NRHS contained in the vector NSVAL.
[in]NSVAL
          NSVAL is INTEGER array, dimension (NNS)
          The values of the number of right hand sides NRHS.
[in]THRESH
          THRESH is DOUBLE PRECISION
          The threshold value for the test ratios.  A result is
          included in the output file if RESULT >= THRESH.  To have
          every test ratio printed, use THRESH = 0.
[in]TSTERR
          TSTERR is LOGICAL
          Flag that indicates whether error exits are to be tested.
[in]NMAX
          NMAX is INTEGER
          The maximum value permitted for M or N, used in dimensioning
          the work arrays.
[out]A
          A is DOUBLE PRECISION array, dimension (NMAX*NMAX)
[out]AFAC
          AFAC is DOUBLE PRECISION array, dimension (NMAX*NMAX)
[out]AINV
          AINV is DOUBLE PRECISION array, dimension (NMAX*NMAX)
[out]B
          B is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
          where NSMAX is the largest entry in NSVAL.
[out]X
          X is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
[out]XACT
          XACT is DOUBLE PRECISION array, dimension (NMAX*NSMAX)
[out]WORK
          WORK is DOUBLE PRECISION array, dimension
                      (NMAX*max(3,NSMAX))
[out]RWORK
          RWORK is DOUBLE PRECISION array, dimension
                      (max(2*NMAX,2*NSMAX+NWORK))
[out]IWORK
          IWORK is INTEGER array, dimension (2*NMAX)
[in]NOUT
          NOUT is INTEGER
          The unit number for output.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 182 of file dchkge.f.

185 *
186 * -- LAPACK test routine --
187 * -- LAPACK is a software package provided by Univ. of Tennessee, --
188 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
189 *
190 * .. Scalar Arguments ..
191  LOGICAL TSTERR
192  INTEGER NM, NMAX, NN, NNB, NNS, NOUT
193  DOUBLE PRECISION THRESH
194 * ..
195 * .. Array Arguments ..
196  LOGICAL DOTYPE( * )
197  INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ),
198  $ NVAL( * )
199  DOUBLE PRECISION A( * ), AFAC( * ), AINV( * ), B( * ),
200  $ RWORK( * ), WORK( * ), X( * ), XACT( * )
201 * ..
202 *
203 * =====================================================================
204 *
205 * .. Parameters ..
206  DOUBLE PRECISION ONE, ZERO
207  parameter( one = 1.0d+0, zero = 0.0d+0 )
208  INTEGER NTYPES
209  parameter( ntypes = 11 )
210  INTEGER NTESTS
211  parameter( ntests = 8 )
212  INTEGER NTRAN
213  parameter( ntran = 3 )
214 * ..
215 * .. Local Scalars ..
216  LOGICAL TRFCON, ZEROT
217  CHARACTER DIST, NORM, TRANS, TYPE, XTYPE
218  CHARACTER*3 PATH
219  INTEGER I, IM, IMAT, IN, INB, INFO, IOFF, IRHS, ITRAN,
220  $ IZERO, K, KL, KU, LDA, LWORK, M, MODE, N, NB,
221  $ NERRS, NFAIL, NIMAT, NRHS, NRUN, NT
222  DOUBLE PRECISION AINVNM, ANORM, ANORMI, ANORMO, CNDNUM, DUMMY,
223  $ RCOND, RCONDC, RCONDI, RCONDO
224 * ..
225 * .. Local Arrays ..
226  CHARACTER TRANSS( NTRAN )
227  INTEGER ISEED( 4 ), ISEEDY( 4 )
228  DOUBLE PRECISION RESULT( NTESTS )
229 * ..
230 * .. External Functions ..
231  DOUBLE PRECISION DGET06, DLANGE
232  EXTERNAL dget06, dlange
233 * ..
234 * .. External Subroutines ..
235  EXTERNAL alaerh, alahd, alasum, derrge, dgecon, dgerfs,
238  $ dlatms, xlaenv
239 * ..
240 * .. Intrinsic Functions ..
241  INTRINSIC max, min
242 * ..
243 * .. Scalars in Common ..
244  LOGICAL LERR, OK
245  CHARACTER*32 SRNAMT
246  INTEGER INFOT, NUNIT
247 * ..
248 * .. Common blocks ..
249  COMMON / infoc / infot, nunit, ok, lerr
250  COMMON / srnamc / srnamt
251 * ..
252 * .. Data statements ..
253  DATA iseedy / 1988, 1989, 1990, 1991 / ,
254  $ transs / 'N', 'T', 'C' /
255 * ..
256 * .. Executable Statements ..
257 *
258 * Initialize constants and the random number seed.
259 *
260  path( 1: 1 ) = 'Double precision'
261  path( 2: 3 ) = 'GE'
262  nrun = 0
263  nfail = 0
264  nerrs = 0
265  DO 10 i = 1, 4
266  iseed( i ) = iseedy( i )
267  10 CONTINUE
268 *
269 * Test the error exits
270 *
271  CALL xlaenv( 1, 1 )
272  IF( tsterr )
273  $ CALL derrge( path, nout )
274  infot = 0
275  CALL xlaenv( 2, 2 )
276 *
277 * Do for each value of M in MVAL
278 *
279  DO 120 im = 1, nm
280  m = mval( im )
281  lda = max( 1, m )
282 *
283 * Do for each value of N in NVAL
284 *
285  DO 110 in = 1, nn
286  n = nval( in )
287  xtype = 'N'
288  nimat = ntypes
289  IF( m.LE.0 .OR. n.LE.0 )
290  $ nimat = 1
291 *
292  DO 100 imat = 1, nimat
293 *
294 * Do the tests only if DOTYPE( IMAT ) is true.
295 *
296  IF( .NOT.dotype( imat ) )
297  $ GO TO 100
298 *
299 * Skip types 5, 6, or 7 if the matrix size is too small.
300 *
301  zerot = imat.GE.5 .AND. imat.LE.7
302  IF( zerot .AND. n.LT.imat-4 )
303  $ GO TO 100
304 *
305 * Set up parameters with DLATB4 and generate a test matrix
306 * with DLATMS.
307 *
308  CALL dlatb4( path, imat, m, n, TYPE, KL, KU, ANORM, MODE,
309  $ CNDNUM, DIST )
310 *
311  srnamt = 'DLATMS'
312  CALL dlatms( m, n, dist, iseed, TYPE, RWORK, MODE,
313  $ CNDNUM, ANORM, KL, KU, 'No packing', A, LDA,
314  $ WORK, INFO )
315 *
316 * Check error code from DLATMS.
317 *
318  IF( info.NE.0 ) THEN
319  CALL alaerh( path, 'DLATMS', info, 0, ' ', m, n, -1,
320  $ -1, -1, imat, nfail, nerrs, nout )
321  GO TO 100
322  END IF
323 *
324 * For types 5-7, zero one or more columns of the matrix to
325 * test that INFO is returned correctly.
326 *
327  IF( zerot ) THEN
328  IF( imat.EQ.5 ) THEN
329  izero = 1
330  ELSE IF( imat.EQ.6 ) THEN
331  izero = min( m, n )
332  ELSE
333  izero = min( m, n ) / 2 + 1
334  END IF
335  ioff = ( izero-1 )*lda
336  IF( imat.LT.7 ) THEN
337  DO 20 i = 1, m
338  a( ioff+i ) = zero
339  20 CONTINUE
340  ELSE
341  CALL dlaset( 'Full', m, n-izero+1, zero, zero,
342  $ a( ioff+1 ), lda )
343  END IF
344  ELSE
345  izero = 0
346  END IF
347 *
348 * These lines, if used in place of the calls in the DO 60
349 * loop, cause the code to bomb on a Sun SPARCstation.
350 *
351 * ANORMO = DLANGE( 'O', M, N, A, LDA, RWORK )
352 * ANORMI = DLANGE( 'I', M, N, A, LDA, RWORK )
353 *
354 * Do for each blocksize in NBVAL
355 *
356  DO 90 inb = 1, nnb
357  nb = nbval( inb )
358  CALL xlaenv( 1, nb )
359 *
360 * Compute the LU factorization of the matrix.
361 *
362  CALL dlacpy( 'Full', m, n, a, lda, afac, lda )
363  srnamt = 'DGETRF'
364  CALL dgetrf( m, n, afac, lda, iwork, info )
365 *
366 * Check error code from DGETRF.
367 *
368  IF( info.NE.izero )
369  $ CALL alaerh( path, 'DGETRF', info, izero, ' ', m,
370  $ n, -1, -1, nb, imat, nfail, nerrs,
371  $ nout )
372  trfcon = .false.
373 *
374 *+ TEST 1
375 * Reconstruct matrix from factors and compute residual.
376 *
377  CALL dlacpy( 'Full', m, n, afac, lda, ainv, lda )
378  CALL dget01( m, n, a, lda, ainv, lda, iwork, rwork,
379  $ result( 1 ) )
380  nt = 1
381 *
382 *+ TEST 2
383 * Form the inverse if the factorization was successful
384 * and compute the residual.
385 *
386  IF( m.EQ.n .AND. info.EQ.0 ) THEN
387  CALL dlacpy( 'Full', n, n, afac, lda, ainv, lda )
388  srnamt = 'DGETRI'
389  nrhs = nsval( 1 )
390  lwork = nmax*max( 3, nrhs )
391  CALL dgetri( n, ainv, lda, iwork, work, lwork,
392  $ info )
393 *
394 * Check error code from DGETRI.
395 *
396  IF( info.NE.0 )
397  $ CALL alaerh( path, 'DGETRI', info, 0, ' ', n, n,
398  $ -1, -1, nb, imat, nfail, nerrs,
399  $ nout )
400 *
401 * Compute the residual for the matrix times its
402 * inverse. Also compute the 1-norm condition number
403 * of A.
404 *
405  CALL dget03( n, a, lda, ainv, lda, work, lda,
406  $ rwork, rcondo, result( 2 ) )
407  anormo = dlange( 'O', m, n, a, lda, rwork )
408 *
409 * Compute the infinity-norm condition number of A.
410 *
411  anormi = dlange( 'I', m, n, a, lda, rwork )
412  ainvnm = dlange( 'I', n, n, ainv, lda, rwork )
413  IF( anormi.LE.zero .OR. ainvnm.LE.zero ) THEN
414  rcondi = one
415  ELSE
416  rcondi = ( one / anormi ) / ainvnm
417  END IF
418  nt = 2
419  ELSE
420 *
421 * Do only the condition estimate if INFO > 0.
422 *
423  trfcon = .true.
424  anormo = dlange( 'O', m, n, a, lda, rwork )
425  anormi = dlange( 'I', m, n, a, lda, rwork )
426  rcondo = zero
427  rcondi = zero
428  END IF
429 *
430 * Print information about the tests so far that did not
431 * pass the threshold.
432 *
433  DO 30 k = 1, nt
434  IF( result( k ).GE.thresh ) THEN
435  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
436  $ CALL alahd( nout, path )
437  WRITE( nout, fmt = 9999 )m, n, nb, imat, k,
438  $ result( k )
439  nfail = nfail + 1
440  END IF
441  30 CONTINUE
442  nrun = nrun + nt
443 *
444 * Skip the remaining tests if this is not the first
445 * block size or if M .ne. N. Skip the solve tests if
446 * the matrix is singular.
447 *
448  IF( inb.GT.1 .OR. m.NE.n )
449  $ GO TO 90
450  IF( trfcon )
451  $ GO TO 70
452 *
453  DO 60 irhs = 1, nns
454  nrhs = nsval( irhs )
455  xtype = 'N'
456 *
457  DO 50 itran = 1, ntran
458  trans = transs( itran )
459  IF( itran.EQ.1 ) THEN
460  rcondc = rcondo
461  ELSE
462  rcondc = rcondi
463  END IF
464 *
465 *+ TEST 3
466 * Solve and compute residual for A * X = B.
467 *
468  srnamt = 'DLARHS'
469  CALL dlarhs( path, xtype, ' ', trans, n, n, kl,
470  $ ku, nrhs, a, lda, xact, lda, b,
471  $ lda, iseed, info )
472  xtype = 'C'
473 *
474  CALL dlacpy( 'Full', n, nrhs, b, lda, x, lda )
475  srnamt = 'DGETRS'
476  CALL dgetrs( trans, n, nrhs, afac, lda, iwork,
477  $ x, lda, info )
478 *
479 * Check error code from DGETRS.
480 *
481  IF( info.NE.0 )
482  $ CALL alaerh( path, 'DGETRS', info, 0, trans,
483  $ n, n, -1, -1, nrhs, imat, nfail,
484  $ nerrs, nout )
485 *
486  CALL dlacpy( 'Full', n, nrhs, b, lda, work,
487  $ lda )
488  CALL dget02( trans, n, n, nrhs, a, lda, x, lda,
489  $ work, lda, rwork, result( 3 ) )
490 *
491 *+ TEST 4
492 * Check solution from generated exact solution.
493 *
494  CALL dget04( n, nrhs, x, lda, xact, lda, rcondc,
495  $ result( 4 ) )
496 *
497 *+ TESTS 5, 6, and 7
498 * Use iterative refinement to improve the
499 * solution.
500 *
501  srnamt = 'DGERFS'
502  CALL dgerfs( trans, n, nrhs, a, lda, afac, lda,
503  $ iwork, b, lda, x, lda, rwork,
504  $ rwork( nrhs+1 ), work,
505  $ iwork( n+1 ), info )
506 *
507 * Check error code from DGERFS.
508 *
509  IF( info.NE.0 )
510  $ CALL alaerh( path, 'DGERFS', info, 0, trans,
511  $ n, n, -1, -1, nrhs, imat, nfail,
512  $ nerrs, nout )
513 *
514  CALL dget04( n, nrhs, x, lda, xact, lda, rcondc,
515  $ result( 5 ) )
516  CALL dget07( trans, n, nrhs, a, lda, b, lda, x,
517  $ lda, xact, lda, rwork, .true.,
518  $ rwork( nrhs+1 ), result( 6 ) )
519 *
520 * Print information about the tests that did not
521 * pass the threshold.
522 *
523  DO 40 k = 3, 7
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 = 9998 )trans, n, nrhs,
528  $ imat, k, result( k )
529  nfail = nfail + 1
530  END IF
531  40 CONTINUE
532  nrun = nrun + 5
533  50 CONTINUE
534  60 CONTINUE
535 *
536 *+ TEST 8
537 * Get an estimate of RCOND = 1/CNDNUM.
538 *
539  70 CONTINUE
540  DO 80 itran = 1, 2
541  IF( itran.EQ.1 ) THEN
542  anorm = anormo
543  rcondc = rcondo
544  norm = 'O'
545  ELSE
546  anorm = anormi
547  rcondc = rcondi
548  norm = 'I'
549  END IF
550  srnamt = 'DGECON'
551  CALL dgecon( norm, n, afac, lda, anorm, rcond,
552  $ work, iwork( n+1 ), info )
553 *
554 * Check error code from DGECON.
555 *
556  IF( info.NE.0 )
557  $ CALL alaerh( path, 'DGECON', info, 0, norm, n,
558  $ n, -1, -1, -1, imat, nfail, nerrs,
559  $ nout )
560 *
561 * This line is needed on a Sun SPARCstation.
562 *
563  dummy = rcond
564 *
565  result( 8 ) = dget06( rcond, rcondc )
566 *
567 * Print information about the tests that did not pass
568 * the threshold.
569 *
570  IF( result( 8 ).GE.thresh ) THEN
571  IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
572  $ CALL alahd( nout, path )
573  WRITE( nout, fmt = 9997 )norm, n, imat, 8,
574  $ result( 8 )
575  nfail = nfail + 1
576  END IF
577  nrun = nrun + 1
578  80 CONTINUE
579  90 CONTINUE
580  100 CONTINUE
581  110 CONTINUE
582  120 CONTINUE
583 *
584 * Print a summary of the results.
585 *
586  CALL alasum( path, nout, nfail, nrun, nerrs )
587 *
588  9999 FORMAT( ' M = ', i5, ', N =', i5, ', NB =', i4, ', type ', i2,
589  $ ', test(', i2, ') =', g12.5 )
590  9998 FORMAT( ' TRANS=''', a1, ''', N =', i5, ', NRHS=', i3, ', type ',
591  $ i2, ', test(', i2, ') =', g12.5 )
592  9997 FORMAT( ' NORM =''', a1, ''', N =', i5, ',', 10x, ' type ', i2,
593  $ ', test(', i2, ') =', g12.5 )
594  RETURN
595 *
596 * End of DCHKGE
597 *
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: dlaset.f:110
subroutine alasum(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASUM
Definition: alasum.f:73
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
Definition: xlaenv.f:81
subroutine alahd(IOUNIT, PATH)
ALAHD
Definition: alahd.f:107
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
Definition: alaerh.f:147
subroutine dlarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
DLARHS
Definition: dlarhs.f:205
subroutine dget02(TRANS, M, N, NRHS, A, LDA, X, LDX, B, LDB, RWORK, RESID)
DGET02
Definition: dget02.f:135
subroutine dget07(TRANS, N, NRHS, A, LDA, B, LDB, X, LDX, XACT, LDXACT, FERR, CHKFERR, BERR, RESLTS)
DGET07
Definition: dget07.f:165
subroutine dget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
DGET04
Definition: dget04.f:102
subroutine dget01(M, N, A, LDA, AFAC, LDAFAC, IPIV, RWORK, RESID)
DGET01
Definition: dget01.f:107
subroutine dget03(N, A, LDA, AINV, LDAINV, WORK, LDWORK, RWORK, RCOND, RESID)
DGET03
Definition: dget03.f:109
subroutine dlatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
DLATB4
Definition: dlatb4.f:120
subroutine derrge(PATH, NUNIT)
DERRGE
Definition: derrge.f:55
double precision function dget06(RCOND, RCONDC)
DGET06
Definition: dget06.f:55
subroutine dlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
DLATMS
Definition: dlatms.f:321
double precision function dlange(NORM, M, N, A, LDA, WORK)
DLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: dlange.f:114
subroutine dgetrf(M, N, A, LDA, IPIV, INFO)
DGETRF
Definition: dgetrf.f:108
subroutine dgecon(NORM, N, A, LDA, ANORM, RCOND, WORK, IWORK, INFO)
DGECON
Definition: dgecon.f:124
subroutine dgetri(N, A, LDA, IPIV, WORK, LWORK, INFO)
DGETRI
Definition: dgetri.f:114
subroutine dgetrs(TRANS, N, NRHS, A, LDA, IPIV, B, LDB, INFO)
DGETRS
Definition: dgetrs.f:121
subroutine dgerfs(TRANS, N, NRHS, A, LDA, AF, LDAF, IPIV, B, LDB, X, LDX, FERR, BERR, WORK, IWORK, INFO)
DGERFS
Definition: dgerfs.f:185
Here is the call graph for this function:
Here is the caller graph for this function: