


function result = p00_simpson ( prob, int_num )

#*****************************************************************************80
#
## p00_simpson() applies the composite Simpson rule to integrate a function.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license.
#
#  Modified:
#
#    30 December 2011
#
#  Author:
#
#    John Burkardt
#
#  Input:
#
#    integer PROB, the problem index.
#
#    integer INT_NUM, the number of subintervals.
#
#  Output:
#
#    real RESULT, the approximate integral.
#
  a, b = p00_limits ( prob )

  h = ( b - a ) / ( int_num )
#
#  This section of the code was rewritten to vectorize.
#
  x(1:2*int_num+1) = ( ( 2*int_num - ( 0 : 2*int_num ) ) * a   ...
                   + (               ( 0 : 2*int_num ) ) * b ) ...
                     / ( 2*int_num                         )
   
  if ( 1 == int_num )
    result = ...
                    p00_fun ( prob, 1,           x(1) ) ...
      + 4.0 *       p00_fun ( prob, int_num,     x(2) ) ...
      +             p00_fun ( prob, 1,           x(3) )
  else
    result = ...
                    p00_fun ( prob, 1,           x(1)               ) ...
      + 2.0 * sum ( p00_fun ( prob, int_num - 1, x(3:2:2*int_num-1) ) ) ...
      + 4.0 * sum ( p00_fun ( prob, int_num,     x(2:2:2*int_num)   ) )...
      +             p00_fun ( prob, 1,           x(2*int_num+1)     )
   end

  result = result * h / 6.0

  return
end




function test_int_test02 ( )

#*****************************************************************************80
#
## test_int_test02() applies a composite Simpson rule to finite interval 1D problems.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license. 
#
#  Modified:
#
#    12 November 2009
#
#  Author:
#
#    John Burkardt
#
  fprintf ( 1, '\n' )
  fprintf ( 1, 'test_int_test02():\n' )
  fprintf ( 1, '  Composite Simpson rule,\n' )
  fprintf ( 1, '  for 1D finite interval problems.\n' )

  prob_num = p00_prob_num ( )

  fprintf ( 1, '\n' )
  fprintf ( 1, '  Problem       Exact\n' )
  fprintf ( 1, '         Ints   Approx       Error\n' )
#
#  Pick a problem.
#
  for prob = 1 : prob_num
#
#  Some problems have singularities that kill the calculation.
#
    exact = p00_exact ( prob )

    fprintf ( 1, '\n' )
    fprintf ( 1, '  #4d        #14f\n', prob, exact )
#
#  Pick a number of subintervals.
#
    for int_log = 0 : 7

      int_num = 2^int_log

      result = p00_simpson ( prob, int_num )

      error = abs ( exact - result )

      fprintf ( 1, '        #4d  #14f  #14e\n', int_num, result, error )

    end

  end

  return
end
function test_int_test03 ( )

#*****************************************************************************80
#
## test_int_test03() applies a Monte Carlo rule to finite interval 1D problems.
#
#  Licensing:
#
#    This code is distributed under the GNU LGPL license. 
#
#  Modified:
#
#    12 November 2009
#
#  Author:
#
#    John Burkardt
#
  fprintf ( 1, '\n' )
  fprintf ( 1, 'test_int_test03():\n' )
  fprintf ( 1, '  Monte Carlo rule,\n' )
  fprintf ( 1, '  for 1D finite interval problems.\n' )

  prob_num = p00_prob_num ( )

  fprintf ( 1, '\n' )
  fprintf ( 1, '  Problem       Exact\n' )
  fprintf ( 1, '         Pts    Approx       Error\n' )
#
#  Pick a problem.
#
  for prob = 1 : prob_num

    exact = p00_exact ( prob )

    fprintf ( 1, '\n' )
    fprintf ( 1, '  #4d        #14f\n', prob, exact )
#
#  Pick a number of points.
#
    for int_log = 0 : 7

      int_num = 2^int_log

      result = p00_montecarlo ( prob, int_num )

      error = abs ( exact - result )

      fprintf ( 1, '        #4d  #14f  #14e\n', int_num, result, error )

    end

  end

  return
end

