315 SUBROUTINE dchksb2stg( NSIZES, NN, NWDTHS, KK, NTYPES, DOTYPE,
316 $ iseed, thresh, nounit, a, lda, sd, se, d1,
317 $ d2, d3, u, ldu, work, lwork, result, info )
325 INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES,
327 DOUBLE PRECISION THRESH
331 INTEGER ISEED( 4 ), KK( * ), NN( * )
332 DOUBLE PRECISION A( lda, * ), RESULT( * ), SD( * ), SE( * ),
333 $ d1( * ), d2( * ), d3( * ),
334 $ u( ldu, * ), work( * )
340 DOUBLE PRECISION ZERO, ONE, TWO, TEN
341 parameter ( zero = 0.0d0, one = 1.0d0, two = 2.0d0,
343 DOUBLE PRECISION HALF
344 parameter ( half = one / two )
346 parameter ( maxtyp = 15 )
349 LOGICAL BADNN, BADNNB
350 INTEGER I, IINFO, IMODE, ITYPE, J, JC, JCOL, JR, JSIZE,
351 $ jtype, jwidth, k, kmax, lh, lw, mtypes, n,
352 $ nerrs, nmats, nmax, ntest, ntestt
353 DOUBLE PRECISION ANINV, ANORM, COND, OVFL, RTOVFL, RTUNFL,
354 $ temp1, temp2, temp3, temp4, ulp, ulpinv, unfl
357 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KMAGN( maxtyp ),
358 $ kmode( maxtyp ), ktype( maxtyp )
361 DOUBLE PRECISION DLAMCH
369 INTRINSIC abs, dble, max, min, sqrt
372 DATA ktype / 1, 2, 5*4, 5*5, 3*8 /
373 DATA kmagn / 2*1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1,
375 DATA kmode / 2*0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
390 nmax = max( nmax, nn( j ) )
398 kmax = max( kmax, kk( j ) )
402 kmax = min( nmax-1, kmax )
406 IF( nsizes.LT.0 )
THEN
408 ELSE IF( badnn )
THEN
410 ELSE IF( nwdths.LT.0 )
THEN
412 ELSE IF( badnnb )
THEN
414 ELSE IF( ntypes.LT.0 )
THEN
416 ELSE IF( lda.LT.kmax+1 )
THEN
418 ELSE IF( ldu.LT.nmax )
THEN
420 ELSE IF( ( max( lda, nmax )+1 )*nmax.GT.lwork )
THEN
425 CALL xerbla(
'DCHKSBSTG', -info )
431 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 .OR. nwdths.EQ.0 )
436 unfl = dlamch(
'Safe minimum' )
438 ulp = dlamch(
'Epsilon' )*dlamch(
'Base' )
440 rtunfl = sqrt( unfl )
441 rtovfl = sqrt( ovfl )
448 DO 190 jsize = 1, nsizes
450 aninv = one / dble( max( 1, n ) )
452 DO 180 jwidth = 1, nwdths
456 k = max( 0, min( n-1, k ) )
458 IF( nsizes.NE.1 )
THEN
459 mtypes = min( maxtyp, ntypes )
461 mtypes = min( maxtyp+1, ntypes )
464 DO 170 jtype = 1, mtypes
465 IF( .NOT.dotype( jtype ) )
471 ioldsd( j ) = iseed( j )
491 IF( mtypes.GT.maxtyp )
494 itype = ktype( jtype )
495 imode = kmode( jtype )
499 GO TO ( 40, 50, 60 )kmagn( jtype )
506 anorm = ( rtovfl*ulp )*aninv
510 anorm = rtunfl*n*ulpinv
515 CALL dlaset(
'Full', lda, n, zero, zero, a, lda )
517 IF( jtype.LE.15 )
THEN
520 cond = ulpinv*aninv / ten
527 IF( itype.EQ.1 )
THEN
530 ELSE IF( itype.EQ.2 )
THEN
535 a( k+1, jcol ) = anorm
538 ELSE IF( itype.EQ.4 )
THEN
542 CALL dlatms( n, n,
'S', iseed,
'S', work, imode, cond,
543 $ anorm, 0, 0,
'Q', a( k+1, 1 ), lda,
544 $ work( n+1 ), iinfo )
546 ELSE IF( itype.EQ.5 )
THEN
550 CALL dlatms( n, n,
'S', iseed,
'S', work, imode, cond,
551 $ anorm, k, k,
'Q', a, lda, work( n+1 ),
554 ELSE IF( itype.EQ.7 )
THEN
558 CALL dlatmr( n, n,
'S', iseed,
'S', work, 6, one, one,
559 $
'T',
'N', work( n+1 ), 1, one,
560 $ work( 2*n+1 ), 1, one,
'N', idumma, 0, 0,
561 $ zero, anorm,
'Q', a( k+1, 1 ), lda,
564 ELSE IF( itype.EQ.8 )
THEN
568 CALL dlatmr( n, n,
'S', iseed,
'S', work, 6, one, one,
569 $
'T',
'N', work( n+1 ), 1, one,
570 $ work( 2*n+1 ), 1, one,
'N', idumma, k, k,
571 $ zero, anorm,
'Q', a, lda, idumma, iinfo )
573 ELSE IF( itype.EQ.9 )
THEN
577 CALL dlatms( n, n,
'S', iseed,
'P', work, imode, cond,
578 $ anorm, k, k,
'Q', a, lda, work( n+1 ),
581 ELSE IF( itype.EQ.10 )
THEN
587 CALL dlatms( n, n,
'S', iseed,
'P', work, imode, cond,
588 $ anorm, 1, 1,
'Q', a( k, 1 ), lda,
589 $ work( n+1 ), iinfo )
591 temp1 = abs( a( k, i ) ) /
592 $ sqrt( abs( a( k+1, i-1 )*a( k+1, i ) ) )
593 IF( temp1.GT.half )
THEN
594 a( k, i ) = half*sqrt( abs( a( k+1,
595 $ i-1 )*a( k+1, i ) ) )
604 IF( iinfo.NE.0 )
THEN
605 WRITE( nounit, fmt = 9999 )
'Generator', iinfo, n,
615 CALL dlacpy(
' ', k+1, n, a, lda, work, lda )
618 CALL dsbtrd(
'V',
'U', n, k, work, lda, sd, se, u, ldu,
619 $ work( lda*n+1 ), iinfo )
621 IF( iinfo.NE.0 )
THEN
622 WRITE( nounit, fmt = 9999 )
'DSBTRD(U)', iinfo, n,
625 IF( iinfo.LT.0 )
THEN
635 CALL dsbt21(
'Upper', n, k, 1, a, lda, sd, se, u, ldu,
636 $ work, result( 1 ) )
650 CALL dcopy( n, sd, 1, d1, 1 )
652 $
CALL dcopy( n-1, se, 1, work, 1 )
654 CALL dsteqr(
'N', n, d1, work, work( n+1 ), ldu,
655 $ work( n+1 ), iinfo )
656 IF( iinfo.NE.0 )
THEN
657 WRITE( nounit, fmt = 9999 )
'DSTEQR(N)', iinfo, n,
660 IF( iinfo.LT.0 )
THEN
673 CALL dlaset(
'Full', n, 1, zero, zero, sd, 1 )
674 CALL dlaset(
'Full', n, 1, zero, zero, se, 1 )
675 CALL dlacpy(
' ', k+1, n, a, lda, u, ldu )
678 CALL dsytrd_sb2st(
'N',
'N',
"U", n, k, u, ldu, sd, se,
679 $ work, lh, work( lh+1 ), lw, iinfo )
683 CALL dcopy( n, sd, 1, d2, 1 )
685 $
CALL dcopy( n-1, se, 1, work, 1 )
687 CALL dsteqr(
'N', n, d2, work, work( n+1 ), ldu,
688 $ work( n+1 ), iinfo )
689 IF( iinfo.NE.0 )
THEN
690 WRITE( nounit, fmt = 9999 )
'DSTEQR(N)', iinfo, n,
693 IF( iinfo.LT.0 )
THEN
705 DO 110 jr = 0, min( k, n-jc )
706 a( jr+1, jc ) = a( k+1-jr, jc+jr )
709 DO 140 jc = n + 1 - k, n
710 DO 130 jr = min( k, n-jc ) + 1, k
717 CALL dlacpy(
' ', k+1, n, a, lda, work, lda )
720 CALL dsbtrd(
'V',
'L', n, k, work, lda, sd, se, u, ldu,
721 $ work( lda*n+1 ), iinfo )
723 IF( iinfo.NE.0 )
THEN
724 WRITE( nounit, fmt = 9999 )
'DSBTRD(L)', iinfo, n,
727 IF( iinfo.LT.0 )
THEN
738 CALL dsbt21(
'Lower', n, k, 1, a, lda, sd, se, u, ldu,
739 $ work, result( 3 ) )
746 CALL dlaset(
'Full', n, 1, zero, zero, sd, 1 )
747 CALL dlaset(
'Full', n, 1, zero, zero, se, 1 )
748 CALL dlacpy(
' ', k+1, n, a, lda, u, ldu )
751 CALL dsytrd_sb2st(
'N',
'N',
"L", n, k, u, ldu, sd, se,
752 $ work, lh, work( lh+1 ), lw, iinfo )
756 CALL dcopy( n, sd, 1, d3, 1 )
758 $
CALL dcopy( n-1, se, 1, work, 1 )
760 CALL dsteqr(
'N', n, d3, work, work( n+1 ), ldu,
761 $ work( n+1 ), iinfo )
762 IF( iinfo.NE.0 )
THEN
763 WRITE( nounit, fmt = 9999 )
'DSTEQR(N)', iinfo, n,
766 IF( iinfo.LT.0 )
THEN
785 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
786 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
787 temp3 = max( temp3, abs( d1( j ) ), abs( d3( j ) ) )
788 temp4 = max( temp4, abs( d1( j )-d3( j ) ) )
791 result(5) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
792 result(6) = temp4 / max( unfl, ulp*max( temp3, temp4 ) )
797 ntestt = ntestt + ntest
802 IF( result( jr ).GE.thresh )
THEN
807 IF( nerrs.EQ.0 )
THEN
808 WRITE( nounit, fmt = 9998 )
'DSB'
809 WRITE( nounit, fmt = 9997 )
810 WRITE( nounit, fmt = 9996 )
811 WRITE( nounit, fmt = 9995 )
'Symmetric'
812 WRITE( nounit, fmt = 9994 )
'orthogonal',
'''',
813 $
'transpose', (
'''', j = 1, 6 )
816 WRITE( nounit, fmt = 9993 )n, k, ioldsd, jtype,
827 CALL dlasum(
'DSB', nounit, nerrs, ntestt )
830 9999
FORMAT(
' DCHKSBSTG: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
831 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
833 9998
FORMAT( / 1x, a3,
834 $
' -- Real Symmetric Banded Tridiagonal Reduction Routines' )
835 9997
FORMAT(
' Matrix types (see DCHKSBSTG for details): ' )
837 9996
FORMAT( /
' Special Matrices:',
838 $ /
' 1=Zero matrix. ',
839 $
' 5=Diagonal: clustered entries.',
840 $ /
' 2=Identity matrix. ',
841 $
' 6=Diagonal: large, evenly spaced.',
842 $ /
' 3=Diagonal: evenly spaced entries. ',
843 $
' 7=Diagonal: small, evenly spaced.',
844 $ /
' 4=Diagonal: geometr. spaced entries.' )
845 9995
FORMAT(
' Dense ', a,
' Banded Matrices:',
846 $ /
' 8=Evenly spaced eigenvals. ',
847 $
' 12=Small, evenly spaced eigenvals.',
848 $ /
' 9=Geometrically spaced eigenvals. ',
849 $
' 13=Matrix with random O(1) entries.',
850 $ /
' 10=Clustered eigenvalues. ',
851 $
' 14=Matrix with large random entries.',
852 $ /
' 11=Large, evenly spaced eigenvals. ',
853 $
' 15=Matrix with small random entries.' )
855 9994
FORMAT( /
' Tests performed: (S is Tridiag, U is ', a,
',',
856 $ / 20x, a,
' means ', a,
'.', /
' UPLO=''U'':',
857 $ /
' 1= | A - U S U', a1,
' | / ( |A| n ulp ) ',
858 $
' 2= | I - U U', a1,
' | / ( n ulp )', /
' UPLO=''L'':',
859 $ /
' 3= | A - U S U', a1,
' | / ( |A| n ulp ) ',
860 $
' 4= | I - U U', a1,
' | / ( n ulp )' /
' Eig check:',
861 $ /
' 5= | D1 - D2',
'',
' | / ( |D1| ulp ) ',
862 $
' 6= | D1 - D3',
'',
' | / ( |D1| ulp ) ' )
863 9993
FORMAT(
' N=', i5,
', K=', i4,
', seed=', 4( i4,
',' ),
' type ',
864 $ i2,
', test(', i2,
')=', g10.3 )
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
DSTEQR
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
subroutine dlatmr(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, RSIGN, GRADE, DL, MODEL, CONDL, DR, MODER, CONDR, PIVTNG, IPIVOT, KL, KU, SPARSE, ANORM, PACK, A, LDA, IWORK, INFO)
DLATMR
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 dsbt21(UPLO, N, KA, KS, A, LDA, D, E, U, LDU, WORK, RESULT)
DSBT21
subroutine dsbtrd(VECT, UPLO, N, KD, AB, LDAB, D, E, Q, LDQ, WORK, INFO)
DSBTRD
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine dlasum(TYPE, IOUNIT, IE, NRUN)
DLASUM
subroutine dlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
DLATMS
subroutine dchksb2stg(NSIZES, NN, NWDTHS, KK, NTYPES, DOTYPE, ISEED, THRESH, NOUNIT, A, LDA, SD, SE, D1, D2, D3, U, LDU, WORK, LWORK, RESULT, INFO)
DCHKSBSTG