|
Home | Libraries | People | FAQ | More |
The globally adaptive integrators are functions that can be used to integrate integrands with more difficult behaviour.
The minimum arguments required are the integrand, the lower and upper limits of integration, a result which is used to return the approximation of the integral, and an error estimate, which is used to return an estimate of the error in the integral approximation.
By default, the accuracy requested is a relative accuracy of 1e-7 (is it useful to have a default? at least one of the relative and absolute accuracies should be specified).
Requested accuracy can be set using the relative_accuracy, and absolute_accuracy method. If both are set, the more restrictive value is used.
The sample program adaptive_accuracy.cpp illustrates its use.
#include <boost/numeric/quadrature/adaptive.hpp> #include <boost/numeric/quadrature/kronrodgauss.hpp> #include <iostream> #include <cmath> namespace quadrature=boost::numeric::quadrature; struct f { double operator()(double x) const { return x/std::log(x); } }; int main() { double answer, error_estimate; // integrate x/log(x) on [0,1] quadrature::adaptive().relative_accuracy(1e-5).absolute_accuracy(1e-7) (f(), 0., 1., answer, error_estimate); std::cout << "integral(x/log(x)) on [0,1] is " << answer << " with error estimate " << error_estimate << std::endl; return EXIT_SUCCESS; }
An integration kernel can be supplied using the kernel method, and this is used to perform the integrations over the subintervals determined by the global adaptive integrator.
The sample program adaptive_kernel.cpp illustrates its use.
#include <boost/numeric/quadrature/adaptive.hpp> #include <boost/numeric/quadrature/rms.hpp> #include <iostream> #include <cmath> namespace quadrature=boost::numeric::quadrature; struct f { double operator()(double x) const { return x/std::log(x); } }; int main() { double answer, error_estimate; quadrature::rms<19,13> kernel; // integrate x/log(x) on [0,1] quadrature::adaptive().kernel(kernel)(f(), 0., 1., answer, error_estimate); std::cout << "integral(x/log(x)) on [0,1] is " << answer << " with error estimate " << error_estimate << std::endl; return EXIT_SUCCESS; }
The maximum number of sub-intervals the algorithm should use is specified using the max_intervals method, and defaults to 100.
The sample program adaptive_max_intervals.cpp illustrates its use.
#include <boost/numeric/quadrature/adaptive.hpp> #include <boost/numeric/quadrature/kronrodgauss.hpp> #include <iostream> #include <cmath> namespace quadrature=boost::numeric::quadrature; struct f { double operator()(double x) const { return x/std::log(x); } }; int main() { double answer, error_estimate; // integrate x/log(x) on [0,1] quadrature::adaptive().max_intervals(50)(f(), 0., 1., answer, error_estimate); std::cout << "integral(x/log(x)) on [0,1] is " << answer << " with error estimate " << error_estimate << std::endl; return EXIT_SUCCESS; }
If you have many integrands to integrate, and width to reduce the memory allocation overhead of allocating interval storage for each integral, then you can specify the interval storage using the intervals method. In this case the maximum number of intervals is determined by the size of the intervals storage.
The value type of the storage container is given by the kernel_interval trait.
The storage container should have low complexity for reordering, and defaults to a list.
The sample program adaptive_intervals.cpp illustrates its use.
#include <boost/numeric/quadrature/adaptive.hpp> #include <boost/numeric/quadrature/kronrodgauss.hpp> #include <iostream> #include <list> #include <cmath> namespace quadrature=boost::numeric::quadrature; struct f { double operator()(double x) const { return x/std::log(x); } }; int main() { double answer, error_estimate; typedef quadrature::kronrod_gauss<21> kernel_type; kernel_type kernel; typedef quadrature::kernel_interval<kernel_type, double, double>::type interval_type; std::list<interval_type> interval_storage(50); // integrate x/log(x) on [0,1] using 21 point kronrod gauss quadrature quadrature::adaptive().intervals(interval_storage).kernel(kernel) (f(), 0., 1., answer, error_estimate); std::cout << "integral(x/log(x)) on [0,1] is " << answer << " with error estimate " << error_estimate << std::endl; return EXIT_SUCCESS; }
A recorder can be specified with the recorder method and is passed through to the integration kernel each time an integration is performed.
The sample program adaptive_recorder.cpp illustrates its use.
#include <boost/numeric/quadrature/adaptive.hpp> #include <boost/numeric/quadrature/kronrodgauss.hpp> #include <boost/numeric/quadrature/stream_recorder.hpp> #include <iostream> #include <cmath> namespace quadrature=boost::numeric::quadrature; struct f { double operator()(double x) const { return x/std::log(x); } }; int main() { double answer, error_estimate; // declare the recorder quadrature::stream_recorder recorder(std::cout); // integrate x/log(x) on [0,1] quadrature::adaptive().recorder(recorder) (f(), 0., 1., answer, error_estimate); std::cout << "integral(x/log(x)) on [0,1] is " << answer << " with error estimate " << error_estimate << std::endl; return EXIT_SUCCESS; }
Information on the performance of the adaptive integrator can be collected by passing an information object using the info method.
Currently there is only one information scheme implemented, adaptive_info in adaptive.hpp.
The sample program adaptive_info.cpp illustrates its use.
#include <boost/numeric/quadrature/adaptive.hpp> #include <boost/numeric/quadrature/kronrodgauss.hpp> #include <iostream> #include <cmath> namespace quadrature=boost::numeric::quadrature; struct f { double operator()(double x) const { return x/std::log(x); } }; int main() { double answer, error_estimate; // declare the info data quadrature::adaptive_info info; // integrate x/log(x) on [0,1] quadrature::adaptive().info(info)(f(), 0., 1., answer, error_estimate); std::cout << "integral(x/log(x)) on [0,1] is " << answer << " with error estimate " << error_estimate << " with #kernel integrations " << info.num_kernel_evaluations() << " with #intervals " << info.num_intervals() << std::endl; return EXIT_SUCCESS; }
The repeated integration over sub-intervals leads to a series of integral approximations, that can be treated as a mathematical series. There are a number of algorithms that can be applied to a series to predict it's final value, and these can be applied to the integral approximations to reduce the number of integrand evaluations required.
A (possibly stateful) function object can be passed to the adaptive integrator using the accelerator method.
At the moment only one scheme is implemented, wynn_epsilon_algorithm in epsilon.hpp. When the adaptive integrator is used with this accelerator, the routine is equivalent to QUADPACK library's QAGS.
The sample program adaptive_wynn.cpp illustrates its use.
#include <boost/numeric/quadrature/adaptive.hpp> #include <boost/numeric/quadrature/kronrodgauss.hpp> #include <boost/numeric/quadrature/epsilon.hpp> #include <iostream> #include <list> #include <cmath> namespace quadrature=boost::numeric::quadrature; struct f { double operator()(double x) const { return x/std::log(x); } }; int main() { double answer, error_estimate; // declare Wynn's epsilon algorithm quadrature::wynn_epsilon_algorithm<double> epsilon; // declare the info data quadrature::adaptive_info info; // integrate x/log(x) on [0,1] quadrature::adaptive().accelerator(epsilon).info(info) (f(), 0., 1., answer, error_estimate); std::cout << "integral(x/log(x)) on [0,1] is " << answer << " with error estimate " << error_estimate << " with #kernel integrations " << info.num_kernel_evaluations() << " with #intervals " << info.num_intervals() << std::endl; return EXIT_SUCCESS; }
| Copyright © 2007 Hugo Duncan |