PROGRAM f07gffe ! F07GFF Example Program Text ! Mark 23 Release. NAG Copyright 2011. ! .. Use Statements .. USE nag_library, ONLY : dppequ, dscal, f06fcf, nag_wp, x02ajf, x02amf, & x02bhf, x04ccf ! .. Implicit None Statement .. IMPLICIT NONE ! .. Parameters .. REAL (KIND=nag_wp), PARAMETER :: one = 1.0_nag_wp REAL (KIND=nag_wp), PARAMETER :: thresh = 0.1_nag_wp INTEGER, PARAMETER :: nin = 5, nout = 6 CHARACTER (1), PARAMETER :: uplo = 'U' ! .. Local Scalars .. REAL (KIND=nag_wp) :: amax, big, scond, small INTEGER :: i, ifail, info, j, jinc, jj, n ! .. Local Arrays .. REAL (KIND=nag_wp), ALLOCATABLE :: ap(:), s(:) ! .. Intrinsic Functions .. INTRINSIC real ! .. Executable Statements .. WRITE (nout,*) 'F07GFF Example Program Results' WRITE (nout,*) FLUSH (nout) ! Skip heading in data file READ (nin,*) READ (nin,*) n ALLOCATE (ap((n*(n+1))/2),s(n)) ! Read the upper or lower triangular part of the matrix A from ! data file IF (uplo=='U') THEN READ (nin,*) ((ap(i+(j*(j-1))/2),j=i,n),i=1,n) ELSE IF (uplo=='L') THEN READ (nin,*) ((ap(i+((2*n-j)*(j-1))/2),j=1,i),i=1,n) END IF ! Print the matrix A ! ifail: behaviour on error exit ! =0 for hard exit, =1 for quiet-soft, =-1 for noisy-soft ifail = 0 CALL x04ccf(uplo,'Non-unit diagonal',n,ap,'Matrix A',ifail) WRITE (nout,*) ! Compute diagonal scaling factors ! The NAG name equivalent of dppequ is f07gff CALL dppequ(uplo,n,ap,s,scond,amax,info) IF (info>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(1:n) WRITE (nout,*) FLUSH (nout) ! Compute values close to underflow and overflow small = x02amf()/(x02ajf()*real(x02bhf(),kind=nag_wp)) big = one/small IF ((scondbig)) THEN ! Scale A IF (uplo=='U') THEN ! The NAG name equivalent of dscal is f06edf jj = 1 DO j = 1, n CALL dscal(j,s(j),ap(jj),1) CALL f06fcf(j,s,1,ap(jj),1) jj = jj + j END DO ELSE IF (uplo=='L') THEN jj = 1 jinc = n DO j = 1, n CALL dscal(jinc,s(j),ap(jj),1) CALL f06fcf(jinc,s(j),1,ap(jj),1) jj = jj + jinc jinc = jinc - 1 END DO END IF ! Print the scaled matrix ifail = 0 CALL x04ccf(uplo,'Non-unit diagonal',n,ap,'Scaled matrix',ifail) END IF END IF 99999 FORMAT (1X,A,I4,A) 99998 FORMAT (1X,2(A,1P,E8.1)) 99997 FORMAT ((1X,1P,7E11.1)) END PROGRAM f07gffe