adaptint <- function ( f, a, b, n = 10, tol = 1.0e-6 )

#*****************************************************************************80
#
## adaptint estimates an integral using an adaptive midpoint rule.
#
#  Licensing:
#
#    Copyright 2016 James P. Howard, II
#
#    The computer code and data files on this web page are distributed under
#    https://opensource.org/licenses/BSD-2-Clause, the BSD-2-Clause license.
#
#  Modified:
#
#    07 March 2020
#
#  Author:
#
#    Original R code by James Howard;
#    Modifications by John Burkardt.
#
#  Reference:
#
#    James Howard,
#    Computational Methods for Numerical Analysis with R,
#    CRC Press, 2017
#    ISBN13: 978-1-4987-2363-3.
#
{
  source ( "/home/burkardt/public_html/r_src/midpt/midpt.R" )

  if ( n == 1 )
  {
    area <- midpt ( f, a, b, m = 2 )
  }
  else
  {
    q1 <- midpt ( f, a, b, m = 1 )
    q2 <- midpt ( f, a, b, m = 2 )
    if ( 3.0 * tol < abs ( q1 - q2 ) )
    {
      n <- n - 1
      tol <- tol / 2.0
      c <- ( a + b ) / 2.0
      lt <- adaptint ( f, a, c, n = n, tol = tol )
      rt <- adaptint ( f, c, b, n = n, tol = tol )
      area <- lt + rt
    }
    else
    {
      area <- q2
    }
  }

  return ( area )
}