      PROGRAM AREA
************************************************************************
* Program to approximate the integral of a function over the interval  *
* [A,B] using the rectangle method with altitudes chosen at the        *
* midpoints of the subintervals.  Identifiers used are:                *
*     A, B   : the endpoints of the interval of integration            *
*     N      : the number of subintervals used                         *
*     I      : counter                                                 *
*     DELX   : the length of the subintervals                          *
*     MID    : the midpoint of one of the subintervals                 *
*     Y      : the value of the function at X                          *
*     SUM    : the approximating sum                                   *
*     F      : the integrand                                           *
*                                                                      *
* Input:  A, B, and N                                                  *
* Output: Approximation to integral of F on [A, B]                     *
************************************************************************

      REAL F, A, B, DELX, MID, Y, SUM
      INTEGER N, I

      PRINT *, 'ENTER THE INTERVAL ENDPOINTS AND THE # OF SUBINTERVALS'
      READ *, A, B, N
      DELX = (B - A) / REAL(N)

* Initialize the approximating SUM and set MID equal to
* the midpoint of the first subinterval

      SUM = 0
      MID  = A + DELX / 2.0

* Now calculate and display the sum

      DO 10 I = 1, N

*         Calculate the value of the function at MID
          Y = F(MID)

          SUM = SUM + Y
          MID = MID + DELX
10    CONTINUE
      SUM = DELX * SUM
      PRINT *, 'APPROXIMATE VALUE USING ', N, ' SUBINTERVALS IS ', SUM
      END

** F(X) *********
* The integrand *
*****************

      REAL FUNCTION F(X)
      REAL X

      F = X**2 + 1

      END

