* F08GAF Example Program Text * Mark 21 Release. NAG Copyright 2004. * .. Parameters .. INTEGER NIN, NOUT PARAMETER (NIN=5,NOUT=6) INTEGER NMAX PARAMETER (NMAX=10) CHARACTER UPLO PARAMETER (UPLO='U') * .. Local Scalars .. DOUBLE PRECISION EERRBD, EPS INTEGER I, INFO, J, N * .. Local Arrays .. DOUBLE PRECISION AP((NMAX*(NMAX+1))/2), DUMMY(1,1), W(NMAX), + WORK(3*NMAX) * .. External Functions .. DOUBLE PRECISION X02AJF EXTERNAL X02AJF * .. External Subroutines .. EXTERNAL DSPEV * .. Intrinsic Functions .. INTRINSIC ABS, MAX * .. Executable Statements .. WRITE (NOUT,*) 'F08GAF Example Program Results' WRITE (NOUT,*) * Skip heading in data file READ (NIN,*) READ (NIN,*) N IF (N.LE.NMAX) THEN * * Read the upper or lower triangular part of the matrix A from * data file * IF (UPLO.EQ.'U') THEN READ (NIN,*) ((AP(I+(J*(J-1))/2),J=I,N),I=1,N) ELSE IF (UPLO.EQ.'L') THEN READ (NIN,*) ((AP(I+((2*N-J)*(J-1))/2),J=1,I),I=1,N) END IF * * Solve the symmetric eigenvalue problem * CALL DSPEV('No vectors',UPLO,N,AP,W,DUMMY,1,WORK,INFO) * IF (INFO.EQ.0) THEN * * Print solution * WRITE (NOUT,*) 'Eigenvalues' WRITE (NOUT,99999) (W(J),J=1,N) * * Get the machine precision, EPS and compute the approximate * error bound for the computed eigenvalues. Note that for * the 2-norm, max( abs(W(i)) ) = norm(A), and since the * eigenvalues are returned in ascending order * max( abs(W(i)) ) = max( abs(W(1)), abs(W(n))) * EPS = X02AJF() EERRBD = EPS*MAX(ABS(W(1)),ABS(W(N))) * * Print the approximate error bound for the eigenvalues * WRITE (NOUT,*) WRITE (NOUT,*) 'Error estimate for the eigenvalues' WRITE (NOUT,99998) EERRBD ELSE WRITE (NOUT,99997) 'Failure in DSPEV. INFO =', INFO END IF ELSE WRITE (NOUT,*) 'NMAX too small' END IF STOP * 99999 FORMAT (3X,(8F8.4)) 99998 FORMAT (4X,1P,6E11.1) 99997 FORMAT (1X,A,I4) END