174 SUBROUTINE ddrvgb( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, A, LA,
175 $ afb, lafb, asav, b, bsav, x, xact, s, work,
176 $ rwork, iwork, nout )
185 INTEGER LA, LAFB, NN, NOUT, NRHS
186 DOUBLE PRECISION THRESH
190 INTEGER IWORK( * ), NVAL( * )
191 DOUBLE PRECISION A( * ), AFB( * ), ASAV( * ), B( * ), BSAV( * ),
192 $ rwork( * ), s( * ), work( * ), x( * ),
199 DOUBLE PRECISION ONE, ZERO
200 parameter( one = 1.0d+0, zero = 0.0d+0 )
202 parameter( ntypes = 8 )
204 parameter( ntests = 7 )
206 parameter( ntran = 3 )
209 LOGICAL EQUIL, NOFACT, PREFAC, TRFCON, ZEROT
210 CHARACTER DIST, EQUED, FACT, TRANS,
TYPE, XTYPE
212 INTEGER I, I1, I2, IEQUED, IFACT, IKL, IKU, IMAT, IN,
213 $ info, ioff, itran, izero, j, k, k1, kl, ku,
214 $ lda, ldafb, ldb, mode, n, nb, nbmin, nerrs,
215 $ nfact, nfail, nimat, nkl, nku, nrun, nt,
217 DOUBLE PRECISION AINVNM, AMAX, ANORM, ANORMI, ANORMO, ANRMPV,
218 $ cndnum, colcnd, rcond, rcondc, rcondi, rcondo,
219 $ roldc, roldi, roldo, rowcnd, rpvgrw,
223 CHARACTER EQUEDS( 4 ), FACTS( 3 ), TRANSS( ntran )
224 INTEGER ISEED( 4 ), ISEEDY( 4 )
225 DOUBLE PRECISION RESULT( ntests ), BERR( nrhs ),
226 $ errbnds_n( nrhs, 3 ), errbnds_c( nrhs, 3 )
230 DOUBLE PRECISION DGET06, DLAMCH, DLANGB, DLANGE, DLANTB,
232 EXTERNAL lsame, dget06, dlamch, dlangb, dlange, dlantb,
242 INTRINSIC abs, max, min
250 COMMON / infoc / infot, nunit, ok, lerr
251 COMMON / srnamc / srnamt
254 DATA iseedy / 1988, 1989, 1990, 1991 /
255 DATA transs /
'N',
'T',
'C' /
256 DATA facts /
'F',
'N',
'E' /
257 DATA equeds /
'N',
'R',
'C',
'B' /
263 path( 1: 1 ) =
'Double precision'
269 iseed( i ) = iseedy( i )
275 $ CALL
derrvx( path, nout )
294 nkl = max( 1, min( n, 4 ) )
309 ELSE IF( ikl.EQ.2 )
THEN
311 ELSE IF( ikl.EQ.3 )
THEN
313 ELSE IF( ikl.EQ.4 )
THEN
324 ELSE IF( iku.EQ.2 )
THEN
326 ELSE IF( iku.EQ.3 )
THEN
328 ELSE IF( iku.EQ.4 )
THEN
336 ldafb = 2*kl + ku + 1
337 IF( lda*n.GT.la .OR. ldafb*n.GT.lafb )
THEN
338 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
339 $ CALL
aladhd( nout, path )
340 IF( lda*n.GT.la )
THEN
341 WRITE( nout, fmt = 9999 )la, n, kl, ku,
345 IF( ldafb*n.GT.lafb )
THEN
346 WRITE( nout, fmt = 9998 )lafb, n, kl, ku,
353 DO 120 imat = 1, nimat
357 IF( .NOT.dotype( imat ) )
362 zerot = imat.GE.2 .AND. imat.LE.4
363 IF( zerot .AND. n.LT.imat-1 )
369 CALL
dlatb4( path, imat, n, n,
TYPE, KL, KU, ANORM,
370 $ mode, cndnum, dist )
371 rcondc = one / cndnum
374 CALL
dlatms( n, n, dist, iseed,
TYPE, RWORK, MODE,
375 $ cndnum, anorm, kl, ku,
'Z', a, lda, work,
381 CALL
alaerh( path,
'DLATMS', info, 0,
' ', n, n,
382 $ kl, ku, -1, imat, nfail, nerrs, nout )
393 ELSE IF( imat.EQ.3 )
THEN
398 ioff = ( izero-1 )*lda
400 i1 = max( 1, ku+2-izero )
401 i2 = min( kl+ku+1, ku+1+( n-izero ) )
407 DO 30 i = max( 1, ku+2-j ),
408 $ min( kl+ku+1, ku+1+( n-j ) )
418 CALL
dlacpy(
'Full', kl+ku+1, n, a, lda, asav, lda )
421 equed = equeds( iequed )
422 IF( iequed.EQ.1 )
THEN
428 DO 100 ifact = 1, nfact
429 fact = facts( ifact )
430 prefac = lsame( fact,
'F' )
431 nofact = lsame( fact,
'N' )
432 equil = lsame( fact,
'E' )
440 ELSE IF( .NOT.nofact )
THEN
447 CALL
dlacpy(
'Full', kl+ku+1, n, asav, lda,
448 $ afb( kl+1 ), ldafb )
449 IF( equil .OR. iequed.GT.1 )
THEN
454 CALL
dgbequ( n, n, kl, ku, afb( kl+1 ),
455 $ ldafb, s, s( n+1 ), rowcnd,
456 $ colcnd, amax, info )
457 IF( info.EQ.0 .AND. n.GT.0 )
THEN
458 IF( lsame( equed,
'R' ) )
THEN
461 ELSE IF( lsame( equed,
'C' ) )
THEN
464 ELSE IF( lsame( equed,
'B' ) )
THEN
471 CALL
dlaqgb( n, n, kl, ku, afb( kl+1 ),
472 $ ldafb, s, s( n+1 ),
473 $ rowcnd, colcnd, amax,
488 anormo = dlangb(
'1', n, kl, ku, afb( kl+1 ),
490 anormi = dlangb(
'I', n, kl, ku, afb( kl+1 ),
495 CALL
dgbtrf( n, n, kl, ku, afb, ldafb, iwork,
500 CALL
dlaset(
'Full', n, n, zero, one, work,
503 CALL
dgbtrs(
'No transpose', n, kl, ku, n,
504 $ afb, ldafb, iwork, work, ldb,
509 ainvnm = dlange(
'1', n, n, work, ldb,
511 IF( anormo.LE.zero .OR. ainvnm.LE.zero )
THEN
514 rcondo = ( one / anormo ) / ainvnm
520 ainvnm = dlange(
'I', n, n, work, ldb,
522 IF( anormi.LE.zero .OR. ainvnm.LE.zero )
THEN
525 rcondi = ( one / anormi ) / ainvnm
529 DO 90 itran = 1, ntran
533 trans = transs( itran )
534 IF( itran.EQ.1 )
THEN
542 CALL
dlacpy(
'Full', kl+ku+1, n, asav, lda,
549 CALL
dlarhs( path, xtype,
'Full', trans, n,
550 $ n, kl, ku, nrhs, a, lda, xact,
551 $ ldb, b, ldb, iseed, info )
553 CALL
dlacpy(
'Full', n, nrhs, b, ldb, bsav,
556 IF( nofact .AND. itran.EQ.1 )
THEN
563 CALL
dlacpy(
'Full', kl+ku+1, n, a, lda,
564 $ afb( kl+1 ), ldafb )
565 CALL
dlacpy(
'Full', n, nrhs, b, ldb, x,
569 CALL
dgbsv( n, kl, ku, nrhs, afb, ldafb,
570 $ iwork, x, ldb, info )
575 $ CALL
alaerh( path,
'DGBSV ', info,
576 $ izero,
' ', n, n, kl, ku,
577 $ nrhs, imat, nfail, nerrs,
583 CALL
dgbt01( n, n, kl, ku, a, lda, afb,
584 $ ldafb, iwork, work,
587 IF( izero.EQ.0 )
THEN
592 CALL
dlacpy(
'Full', n, nrhs, b, ldb,
594 CALL
dgbt02(
'No transpose', n, n, kl,
595 $ ku, nrhs, a, lda, x, ldb,
596 $ work, ldb, result( 2 ) )
601 CALL
dget04( n, nrhs, x, ldb, xact,
602 $ ldb, rcondc, result( 3 ) )
610 IF( result( k ).GE.thresh )
THEN
611 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
612 $ CALL
aladhd( nout, path )
613 WRITE( nout, fmt = 9997 )
'DGBSV ',
614 $ n, kl, ku, imat, k, result( k )
624 $ CALL
dlaset(
'Full', 2*kl+ku+1, n, zero,
626 CALL
dlaset(
'Full', n, nrhs, zero, zero, x,
628 IF( iequed.GT.1 .AND. n.GT.0 )
THEN
633 CALL
dlaqgb( n, n, kl, ku, a, lda, s,
634 $ s( n+1 ), rowcnd, colcnd,
642 CALL
dgbsvx( fact, trans, n, kl, ku, nrhs, a,
643 $ lda, afb, ldafb, iwork, equed,
644 $ s, s( n+1 ), b, ldb, x, ldb,
645 $ rcond, rwork, rwork( nrhs+1 ),
646 $ work, iwork( n+1 ), info )
651 $ CALL
alaerh( path,
'DGBSVX', info, izero,
652 $ fact // trans, n, n, kl, ku,
653 $ nrhs, imat, nfail, nerrs,
662 DO 60 i = max( ku+2-j, 1 ),
663 $ min( n+ku+1-j, kl+ku+1 )
664 anrmpv = max( anrmpv,
665 $ abs( a( i+( j-1 )*lda ) ) )
668 rpvgrw = dlantb(
'M',
'U',
'N', info,
669 $ min( info-1, kl+ku ),
670 $ afb( max( 1, kl+ku+2-info ) ),
672 IF( rpvgrw.EQ.zero )
THEN
675 rpvgrw = anrmpv / rpvgrw
678 rpvgrw = dlantb(
'M',
'U',
'N', n, kl+ku,
680 IF( rpvgrw.EQ.zero )
THEN
683 rpvgrw = dlangb(
'M', n, kl, ku, a,
684 $ lda, work ) / rpvgrw
687 result( 7 ) = abs( rpvgrw-work( 1 ) ) /
688 $ max( work( 1 ), rpvgrw ) /
691 IF( .NOT.prefac )
THEN
696 CALL
dgbt01( n, n, kl, ku, a, lda, afb,
697 $ ldafb, iwork, work,
709 CALL
dlacpy(
'Full', n, nrhs, bsav, ldb,
711 CALL
dgbt02( trans, n, n, kl, ku, nrhs,
712 $ asav, lda, x, ldb, work, ldb,
718 IF( nofact .OR. ( prefac .AND.
719 $ lsame( equed,
'N' ) ) )
THEN
720 CALL
dget04( n, nrhs, x, ldb, xact,
721 $ ldb, rcondc, result( 3 ) )
723 IF( itran.EQ.1 )
THEN
728 CALL
dget04( n, nrhs, x, ldb, xact,
729 $ ldb, roldc, result( 3 ) )
735 CALL
dgbt05( trans, n, kl, ku, nrhs, asav,
736 $ lda, b, ldb, x, ldb, xact,
737 $ ldb, rwork, rwork( nrhs+1 ),
746 result( 6 ) = dget06( rcond, rcondc )
751 IF( .NOT.trfcon )
THEN
753 IF( result( k ).GE.thresh )
THEN
754 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
755 $ CALL
aladhd( nout, path )
757 WRITE( nout, fmt = 9995 )
758 $
'DGBSVX', fact, trans, n, kl,
759 $ ku, equed, imat, k,
762 WRITE( nout, fmt = 9996 )
763 $
'DGBSVX', fact, trans, n, kl,
764 $ ku, imat, k, result( k )
771 IF( result( 1 ).GE.thresh .AND. .NOT.
773 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
774 $ CALL
aladhd( nout, path )
776 WRITE( nout, fmt = 9995 )
'DGBSVX',
777 $ fact, trans, n, kl, ku, equed,
778 $ imat, 1, result( 1 )
780 WRITE( nout, fmt = 9996 )
'DGBSVX',
781 $ fact, trans, n, kl, ku, imat, 1,
787 IF( result( 6 ).GE.thresh )
THEN
788 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
789 $ CALL
aladhd( nout, path )
791 WRITE( nout, fmt = 9995 )
'DGBSVX',
792 $ fact, trans, n, kl, ku, equed,
793 $ imat, 6, result( 6 )
795 WRITE( nout, fmt = 9996 )
'DGBSVX',
796 $ fact, trans, n, kl, ku, imat, 6,
802 IF( result( 7 ).GE.thresh )
THEN
803 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
804 $ CALL
aladhd( nout, path )
806 WRITE( nout, fmt = 9995 )
'DGBSVX',
807 $ fact, trans, n, kl, ku, equed,
808 $ imat, 7, result( 7 )
810 WRITE( nout, fmt = 9996 )
'DGBSVX',
811 $ fact, trans, n, kl, ku, imat, 7,
824 CALL
dlacpy(
'Full', kl+ku+1, n, asav, lda, a,
826 CALL
dlacpy(
'Full', n, nrhs, bsav, ldb, b, ldb )
829 $ CALL
dlaset(
'Full', 2*kl+ku+1, n, zero, zero,
831 CALL
dlaset(
'Full', n, nrhs, zero, zero, x, ldb )
832 IF( iequed.GT.1 .AND. n.GT.0 )
THEN
837 CALL
dlaqgb( n, n, kl, ku, a, lda, s, s( n+1 ),
838 $ rowcnd, colcnd, amax, equed )
846 CALL
dgbsvxx( fact, trans, n, kl, ku, nrhs, a, lda,
847 $ afb, ldafb, iwork, equed, s, s( n+1 ), b, ldb,
848 $ x, ldb, rcond, rpvgrw_svxx, berr, n_err_bnds,
849 $ errbnds_n, errbnds_c, 0, zero, work,
850 $ iwork( n+1 ), info )
854 IF( info.EQ.n+1 ) goto 90
855 IF( info.NE.izero )
THEN
856 CALL
alaerh( path,
'DGBSVXX', info, izero,
857 $ fact // trans, n, n, -1, -1, nrhs,
858 $ imat, nfail, nerrs, nout )
866 IF ( info .GT. 0 .AND. info .LT. n+1 )
THEN
867 rpvgrw = dla_gbrpvgrw(n, kl, ku, info, a, lda,
870 rpvgrw = dla_gbrpvgrw(n, kl, ku, n, a, lda,
874 result( 7 ) = abs( rpvgrw-rpvgrw_svxx ) /
875 $ max( rpvgrw_svxx, rpvgrw ) /
878 IF( .NOT.prefac )
THEN
883 CALL
dgbt01( n, n, kl, ku, a, lda, afb, ldafb,
884 $ iwork, work, result( 1 ) )
895 CALL
dlacpy(
'Full', n, nrhs, bsav, ldb, work,
897 CALL
dgbt02( trans, n, n, kl, ku, nrhs, asav,
898 $ lda, x, ldb, work, ldb,
903 IF( nofact .OR. ( prefac .AND. lsame( equed,
905 CALL
dget04( n, nrhs, x, ldb, xact, ldb,
906 $ rcondc, result( 3 ) )
908 IF( itran.EQ.1 )
THEN
913 CALL
dget04( n, nrhs, x, ldb, xact, ldb,
914 $ roldc, result( 3 ) )
923 result( 6 ) = dget06( rcond, rcondc )
928 IF( .NOT.trfcon )
THEN
930 IF( result( k ).GE.thresh )
THEN
931 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
932 $ CALL
aladhd( nout, path )
934 WRITE( nout, fmt = 9995 )
'DGBSVXX',
935 $ fact, trans, n, kl, ku, equed,
936 $ imat, k, result( k )
938 WRITE( nout, fmt = 9996 )
'DGBSVXX',
939 $ fact, trans, n, kl, ku, imat, k,
947 IF( result( 1 ).GE.thresh .AND. .NOT.prefac )
949 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
950 $ CALL
aladhd( nout, path )
952 WRITE( nout, fmt = 9995 )
'DGBSVXX', fact,
953 $ trans, n, kl, ku, equed, imat, 1,
956 WRITE( nout, fmt = 9996 )
'DGBSVXX', fact,
957 $ trans, n, kl, ku, imat, 1,
963 IF( result( 6 ).GE.thresh )
THEN
964 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
965 $ CALL
aladhd( nout, path )
967 WRITE( nout, fmt = 9995 )
'DGBSVXX', fact,
968 $ trans, n, kl, ku, equed, imat, 6,
971 WRITE( nout, fmt = 9996 )
'DGBSVXX', fact,
972 $ trans, n, kl, ku, imat, 6,
978 IF( result( 7 ).GE.thresh )
THEN
979 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
980 $ CALL
aladhd( nout, path )
982 WRITE( nout, fmt = 9995 )
'DGBSVXX', fact,
983 $ trans, n, kl, ku, equed, imat, 7,
986 WRITE( nout, fmt = 9996 )
'DGBSVXX', fact,
987 $ trans, n, kl, ku, imat, 7,
1005 CALL
alasvm( path, nout, nfail, nrun, nerrs )
1011 9999
FORMAT(
' *** In DDRVGB, LA=', i5,
' is too small for N=', i5,
1012 $
', KU=', i5,
', KL=', i5, /
' ==> Increase LA to at least ',
1014 9998
FORMAT(
' *** In DDRVGB, LAFB=', i5,
' is too small for N=', i5,
1015 $
', KU=', i5,
', KL=', i5, /
1016 $
' ==> Increase LAFB to at least ', i5 )
1017 9997
FORMAT( 1x, a,
', N=', i5,
', KL=', i5,
', KU=', i5,
', type ',
1018 $ i1,
', test(', i1,
')=', g12.5 )
1019 9996
FORMAT( 1x, a,
'( ''', a1,
''',''', a1,
''',', i5,
',', i5,
',',
1020 $ i5,
',...), type ', i1,
', test(', i1,
')=', g12.5 )
1021 9995
FORMAT( 1x, a,
'( ''', a1,
''',''', a1,
''',', i5,
',', i5,
',',
1022 $ i5,
',...), EQUED=''', a1,
''', type ', i1,
', test(', i1,
subroutine dgbt01(M, N, KL, KU, A, LDA, AFAC, LDAFAC, IPIV, WORK, RESID)
DGBT01
subroutine dgbsv(N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB, INFO)
DGBSV computes the solution to system of linear equations A * X = B for GB matrices (simple driver) ...
subroutine dlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
DLATMS
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
subroutine derrvx(PATH, NUNIT)
DERRVX
subroutine alaerh(PATH, SUBNAM, INFO, INFOE, OPTS, M, N, KL, KU, N5, IMAT, NFAIL, NERRS, NOUT)
ALAERH
subroutine dlarhs(PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, ISEED, INFO)
DLARHS
subroutine dlatb4(PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, CNDNUM, DIST)
DLATB4
subroutine dgbtrf(M, N, KL, KU, AB, LDAB, IPIV, INFO)
DGBTRF
subroutine ddrvgb(DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, A, LA, AFB, LAFB, ASAV, B, BSAV, X, XACT, S, WORK, RWORK, IWORK, NOUT)
DDRVGB
subroutine dgbequb(M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND, AMAX, INFO)
DGBEQUB
subroutine dgbsvx(FACT, TRANS, N, KL, KU, NRHS, AB, LDAB, AFB, LDAFB, IPIV, EQUED, R, C, B, LDB, X, LDX, RCOND, FERR, BERR, WORK, IWORK, INFO)
DGBSVX computes the solution to system of linear equations A * X = B for GB matrices ...
subroutine dlaqgb(M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND, AMAX, EQUED)
DLAQGB scales a general band matrix, using row and column scaling factors computed by sgbequ...
subroutine xlaenv(ISPEC, NVALUE)
XLAENV
subroutine dgbtrs(TRANS, N, KL, KU, NRHS, AB, LDAB, IPIV, B, LDB, INFO)
DGBTRS
subroutine dgbequ(M, N, KL, KU, AB, LDAB, R, C, ROWCND, COLCND, AMAX, INFO)
DGBEQU
subroutine debchvxx(THRESH, PATH)
DEBCHVXX
subroutine aladhd(IOUNIT, PATH)
ALADHD
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dget04(N, NRHS, X, LDX, XACT, LDXACT, RCOND, RESID)
DGET04
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...
subroutine dgbsvxx(FACT, TRANS, N, KL, KU, NRHS, AB, LDAB, AFB, LDAFB, IPIV, EQUED, R, C, B, LDB, X, LDX, RCOND, RPVGRW, BERR, N_ERR_BNDS, ERR_BNDS_NORM, ERR_BNDS_COMP, NPARAMS, PARAMS, WORK, IWORK, INFO)
DGBSVXX computes the solution to system of linear equations A * X = B for GB matrices ...
subroutine dgbt05(TRANS, N, KL, KU, NRHS, AB, LDAB, B, LDB, X, LDX, XACT, LDXACT, FERR, BERR, RESLTS)
DGBT05
subroutine dgbt02(TRANS, M, N, KL, KU, NRHS, A, LDA, X, LDX, B, LDB, RESID)
DGBT02