      SUBROUTINE Q1DA(A,B,EPS,R,E,KF,IFLAG)
C***BEGIN PROLOGUE  Q1DA
C***DATE WRITTEN   821018   (YYMMDD)
C***REVISION DATE  870525   (YYMMDD)
C***CATEGORY NO.  H2A1A1
C***KEYWORDS   ADAPTIVE  QUADRATURE, AUTOMATIC  QUADRATURE
C***AUTHOR  KAHANER, DAVID K., SCIENTIFIC COMPUTING DIVISION, NBS.
C***PURPOSE  Approximates one dimensional integrals of user defined
C            functions, easy to use.
C
C***DESCRIPTION
C       Q1DA IS A SUBROUTINE FOR THE AUTOMATIC EVALUATION
C            OF THE DEFINITE INTEGRAL OF A USER DEFINED FUNCTION
C            OF ONE VARIABLE.
C
C       From the book "Numerical Methods and Software"
C            by  D. Kahaner, C. Moler, S. Nash
C               Prentice Hall 1988
C
C       A R G U M E N T S   I N   T H E   C A L L   S E Q U E N C E
C
C       A
C       B     (INPUT) THE ENDPOINTS OF THE INTEGRATION INTERVAL
C       EPS   (INPUT) THE ACCURACY TO WHICH YOU WANT THE INTEGRAL
C                COMPUTED.  IF YOU WANT 2 DIGITS OF ACCURACY SET
C                EPS=.01, FOR 3 DIGITS SET EPS=.001, ETC.
C                EPS MUST BE POSITIVE.
C       R     (OUTPUT) Q1DA'S BEST ESTIMATE OF YOUR INTEGRAL
C       E     (OUTPUT) AN ESTIMATE OF ABS(INTEGRAL-R)
C       KF    (OUTPUT) THE COST OF THE INTEGRATION, MEASURED IN
C                   NUMBER OF EVALUATIONS OF YOUR INTEGRAND.
C                   KF WILL ALWAYS BE AT LEAST 30.
C       IFLAG (OUTPUT) TERMINATION FLAG...POSSIBLE VALUES ARE
C               0   NORMAL COMPLETION, E SATISFIES
C                        E<EPS  AND  E<EPS*ABS(R)
C               1   NORMAL COMPLETION, E SATISFIES
C                        E<EPS, BUT E>EPS*ABS(R)
C               2   NORMAL COMPLETION, E SATISFIES
C                        E<EPS*ABS(R), BUT E>EPS
C               3   NORMAL COMPLETION BUT EPS WAS TOO SMALL TO
C                     SATISFY ABSOLUTE OR RELATIVE ERROR REQUEST.
C
C               4   ABORTED CALCULATION BECAUSE OF SERIOUS ROUNDING
C                     ERROR.  PROBABLY E AND R ARE CONSISTENT.
C               5   ABORTED CALCULATION BECAUSE OF INSUFFICIENT STORAGE.
C                     R AND E ARE CONSISTENT.
C               6   ABORTED CALCULATION BECAUSE OF SERIOUS DIFFICULTIES
C                     MEETING YOUR ERROR REQUEST.
C               7   ABORTED CALCULATION BECAUSE EPS WAS SET <= 0.0
C
C            NOTE...IF IFLAG=3, 4, 5 OR 6 CONSIDER USING Q1DAX INSTEAD.
C
C    W H E R E   I S   Y O U R   I N T E G R A N D ?
C
C        YOU MUST WRITE A FORTRAN FUNCTION, CALLED F, TO EVALUATE
C        THE INTEGRAND.  USUALLY THIS LOOKS LIKE...
C                 FUNCTION F(X)
C                    F=(EVALUATE THE INTEGRAND AT THE POINT X)
C                    RETURN
C                 END
C
C
C    T Y P I C A L   P R O B L E M   S E T U P
C
C          A=0.0
C          B=1.0          (SET INTERVAL ENDPOINTS TO [0,1])
C          EPS=0.001       (SET ACCURACY REQUEST FOR 3 DIGITS)
C          CALL Q1DA(A,B,EPS,R,E,KF,IFLAG)
C          END
C          FUNCTION F(X)
C            F=SIN(2.*X)-SQRT(X)     (FOR EXAMPLE)
C            RETURN
C          END
C      FOR THIS SAMPLE PROBLEM, THE OUTPUT IS
C  0.0    1.0     .001    .041406750    .69077E-07    30    0
C
C    R E M A R K   I.
C
C           A SMALL AMOUT OF RANDOMIZATION IS BUILT INTO THIS PROGRAM.
C           CALLING Q1DA A FEW TIMES IN SUCCESSION WILL GIVE DIFFERENT
C           BUT HOPEFULLY CONSISTENT RESULTS.
C
C   R E M A R K   II.
C
C           THIS ROUTINE IS DESIGNED FOR INTEGRATION OVER A FINITE
C           INTERVAL.  THUS THE INPUT ARGUMENTS A AND B MUST BE
C           VALID REAL NUMBERS ON YOUR COMPUTER.  IF YOU WANT TO DO
C           AN INTEGRAL OVER AN INFINITE INTERVAL SET A OR B OR BOTH
C           LARGE ENOUGH SO THAT THE INTERVAL [A,B] CONTAINS MOST OF
C           THE INTEGRAND.  CARE IS NECESSARY, HOWEVER.  FOR EXAMPLE,
C           TO INTEGRATE EXP(-X*X) ON THE ENTIRE REAL LINE ONE COULD
C           TAKE A=-20., B=20. OR SIMILAR VALUES TO GET GOOD RESULTS.
C           IF YOU TOOK A=-1.E10 AND B=+1.E10 TWO BAD THINGS WOULD
C           OCCUR. FIRST, YOU WILL CERTAINLY GET AN ERROR MESSAGE FROM
C           THE EXP ROUTINE, AS ITS ARGUMENT IS TOO SMALL.  OTHER
C           THINGS COULD HAPPEN TOO, FOR EXAMPLE AN UNDERFLOW.
C           SECOND, EVEN IF THE ARITHMETIC WORKED PROPERLY Q1DA WILL
C           SURELY GIVE AN INCORRECT ANSWER, BECAUSE ITS FIRST TRY
C           AT SAMPLING THE INTEGRAND IS BASED ON YOUR SCALING AND
C           IT IS VERY UNLIKELY TO SELECT EVALUATION POINTS IN THE
C           INFINITESMALLY SMALL INTERVAL [-20,20] WHERE ALL THE
C           INTEGRAND IS CONCENTRATED, WHEN A, B ARE SO LARGE.
C
C    M O R E   F L E X I B I L I T Y
C
C           Q1DA IS AN EASY TO USE DRIVER FOR ANOTHER PROGRAM, Q1DAX.
C           Q1DAX PROVIDES SEVERAL OPTIONS WHICH ARE NOT AVAILABLE
C                WITH Q1DA.
C
C***REFERENCES  (NONE)
C***ROUTINES CALLED  Q1DAX
C***END PROLOGUE  Q1DA
C
      REAL A,B,E,EPS,F,FMAX,FMIN,R,W(50,6)
      LOGICAL RST
      EXTERNAL F
C
C***FIRST EXECUTUTABLE STATEMENT  Q1DA
      NINT=1
      RST = .FALSE.
      NMAX=50
      CALL Q1DAX(F,A,B,EPS,R,E,NINT,RST,W,NMAX,FMIN,FMAX,KF,IFLAG)
      RETURN
      END
