libs/math/example/ooura_fourier_integrals_example.cpp

   1 // Copyright Paul A. Bristow, 2019
   2 // Copyright Nick Thompson, 2019
   3
   4 // Use, modification and distribution are subject to the
   5 // Boost Software License, Version 1.0.
   6 // (See accompanying file LICENSE_1_0.txt
   7 // or copy at http://www.boost.org/LICENSE_1_0.txt)
   8
   9 #ifdef BOOST_NO_CXX11_LAMBDAS
  10 #  error "This example requires a C++11 compiler that supports lambdas.  Try C++11 or later."
  11 #endif
  12
  13 //#define BOOST_MATH_INSTRUMENT_OOURA // or -DBOOST_MATH_INSTRUMENT_OOURA etc for diagnostics.
  14
  15 #include <boost/math/quadrature/ooura_fourier_integrals.hpp>
  16 #include <boost/math/constants/constants.hpp>  // For pi (including for multiprecision types, if used.)
  17
  18 #include <cmath>
  19 #include <iostream>
  20 #include <limits>
  21 #include <iostream>
  22
  23 int main()
  24 {
  25   try
  26   {
  27   std::cout.precision(std::numeric_limits<double>::max_digits10); // Show all potentially significant digits.
  28
  29   using boost::math::quadrature::ooura_fourier_sin;
  30   using boost::math::constants::half_pi;
  31
  32 //[ooura_fourier_integrals_example_1
  33   ooura_fourier_sin<double>integrator = ooura_fourier_sin<double>();
  34   // Use the default tolerance root_epsilon and eight levels for type double.
  35
  36   auto f = [](double x)
  37   { // Simple reciprocal function for sinc.
  38     return 1 / x;
  39   };
  40
  41   double omega = 1;
  42   std::pair<double, double> result = integrator.integrate(f, omega);
  43   std::cout << "Integral = " << result.first << ", relative error estimate " << result.second << std::endl;
  44
  45 //] [/ooura_fourier_integrals_example_1]
  46
  47 //[ooura_fourier_integrals_example_2
  48
  49   constexpr double expected = half_pi<double>();
  50   std::cout << "pi/2 =     " << expected << ", difference " << result.first - expected << std::endl;
  51 //] [/ooura_fourier_integrals_example_2]
  52   }
  53   catch (std::exception const & ex)
  54   {
  55     // Lacking try&catch blocks, the program will abort after any throw, whereas the
  56     // message below from the thrown exception will give some helpful clues as to the cause of the problem.
  57     std::cout << "\n""Message from thrown exception was:\n   " << ex.what() << std::endl;
  58   }
  59 } // int main()
  60
  61 /*
  62
  63 //[ooura_fourier_integrals_example_output_1
  64
  65 integral = 1.5707963267948966, relative error estimate 1.2655356398390254e-11
  66 pi/2 =     1.5707963267948966, difference 0
  67
  68 //] [/ooura_fourier_integrals_example_output_1]
  69
  70
  71 //[ooura_fourier_integrals_example_diagnostic_output_1
  72
  73 ooura_fourier_sin with relative error goal 1.4901161193847656e-08 & 8 levels.
  74 h = 1.000000000000000, I_h = 1.571890732004545 = 0x1.92676e56d853500p+0, absolute error estimate = nan
  75 h = 0.500000000000000, I_h = 1.570793292491940 = 0x1.921f825c076f600p+0, absolute error estimate = 1.097439512605325e-03
  76 h = 0.250000000000000, I_h = 1.570796326814776 = 0x1.921fb54458acf00p+0, absolute error estimate = 3.034322835882008e-06
  77 h = 0.125000000000000, I_h = 1.570796326794897 = 0x1.921fb54442d1800p+0, absolute error estimate = 1.987898734512328e-11
  78 Integral = 1.570796326794897e+00, relative error estimate 1.265535639839025e-11
  79 pi/2 =     1.570796326794897e+00, difference 0.000000000000000e+00
  80
  81 //] [/ooura_fourier_integrals_example_diagnostic_output_1]
  82
  83 */