* F07HFF Example Program Text * Mark 21 Release. NAG Copyright 2004. * .. Parameters .. INTEGER NIN, NOUT PARAMETER (NIN=5,NOUT=6) INTEGER KDMAX, NMAX PARAMETER (KDMAX=4,NMAX=8) INTEGER LDAB PARAMETER (LDAB=KDMAX+1) CHARACTER UPLO PARAMETER (UPLO='U') * .. Local Scalars .. DOUBLE PRECISION AMAX, BIG, SCOND, SJ, SMALL INTEGER I, IFAIL, INFO, J, JJ, KD, N * .. Local Arrays .. DOUBLE PRECISION AB(LDAB,NMAX), S(NMAX) * .. External Functions .. DOUBLE PRECISION X02AJF, X02AMF INTEGER X02BHF EXTERNAL X02AJF, X02AMF, X02BHF * .. External Subroutines .. EXTERNAL DPBEQU, X04CEF * .. Intrinsic Functions .. INTRINSIC MAX, MIN * .. Executable Statements .. WRITE (NOUT,*) 'F07HFF Example Program Results' WRITE (NOUT,*) * Skip heading in data file READ (NIN,*) READ (NIN,*) N, KD IF (N.LE.NMAX .AND. KD.LE.KDMAX) THEN * * Read the upper or lower triangular part of the band matrix A * from data file * IF (UPLO.EQ.'U') THEN DO 20 I = 1, N READ (NIN,*) (AB(KD+1+I-J,J),J=I,MIN(N,I+KD)) 20 CONTINUE ELSE IF (UPLO.EQ.'L') THEN DO 40 I = 1, N READ (NIN,*) (AB(1+I-J,J),J=MAX(1,I-KD),I) 40 CONTINUE END IF * * Print the matrix A * IFAIL = 0 IF (UPLO.EQ.'U') THEN CALL X04CEF(N,N,0,KD,AB,LDAB,'Matrix A',IFAIL) ELSE IF (UPLO.EQ.'L') THEN CALL X04CEF(N,N,KD,0,AB,LDAB,'Matrix A',IFAIL) END IF WRITE (NOUT,*) * * Compute diagonal scaling factors * CALL DPBEQU(UPLO,N,KD,AB,LDAB,S,SCOND,AMAX,INFO) * IF (INFO.GT.0) THEN WRITE (NOUT,99999) 'Diagonal element', INFO, + ' of A is non positive' ELSE * * Print SCOND, AMAX and the scale factors * WRITE (NOUT,99998) 'SCOND = ', SCOND, ', AMAX = ', AMAX WRITE (NOUT,*) WRITE (NOUT,*) 'Diagonal scaling factors' WRITE (NOUT,99997) (S(I),I=1,N) WRITE (NOUT,*) * * Compute values close to underflow and overflow * SMALL = X02AMF()/(X02AJF()*X02BHF()) BIG = 1.0D0/SMALL IF ((SCOND.LT.0.1D0) .OR. (AMAX.LT.SMALL) .OR. (AMAX.GT.BIG) + ) THEN * * Scale A * IF (UPLO.EQ.'U') THEN DO 80 J = 1, N SJ = S(J) JJ = KD + 1 - J DO 60 I = MAX(1,J-KD), J AB(I+JJ,J) = S(I)*AB(I+JJ,J)*SJ 60 CONTINUE 80 CONTINUE ELSE IF (UPLO.EQ.'L') THEN DO 120 J = 1, N SJ = S(J) JJ = 1 - J DO 100 I = J, MIN(N,J+KD) AB(I+JJ,J) = S(I)*AB(I+JJ,J)*SJ 100 CONTINUE 120 CONTINUE END IF * * Print the scaled matrix * IFAIL = 0 IF (UPLO.EQ.'U') THEN CALL X04CEF(N,N,0,KD,AB,LDAB,'Scaled matrix',IFAIL) ELSE IF (UPLO.EQ.'L') THEN CALL X04CEF(N,N,KD,0,AB,LDAB,'Scaled matrix',IFAIL) END IF END IF END IF ELSE WRITE (NOUT,*) 'NMAX and/or KDMAX too small' END IF * 99999 FORMAT (1X,A,I4,A) 99998 FORMAT (1X,2(A,1P,E7.1)) 99997 FORMAT ((1X,1P,7E11.1)) END