include/apfel/integrator.h

0001 //
0002 // APFEL++ 2017
0003 //
0004 // Author: Valerio Bertone: valerio.bertone@cern.ch
0005 //
0006
0007 #pragma once
0008
0009 #include <functional>
0010 #include <vector>
0011
0012 namespace apfel
0013 {
0014   /**
0015    * @brief The Integrator class performs unidimensional numerical
0016    * integrations using the Guassian quadrature.
0017    */
0018   class Integrator
0019   {
0020   public:
0021     /**
0022      * @name Enumerator for the currently available integration
0023      * methods.
0024      */
0025     enum IntegrationMethod: int {GAUSS_LEGENDRE, GAUSS_KRONROD};
0026
0027     /**
0028      * @name Constructors
0029      * List of constructors.
0030      */
0031     ///@{
0032     /**
0033      * @brief The Integrator constructor.
0034      * @param func: The function of one variable to be integrated
0035      * @param method: The integration method to be used (default: GAUSS_KRONROD)
0036      */
0037     Integrator(std::function<double(double const&)> const& func, IntegrationMethod const& method = GAUSS_KRONROD);
0038
0039     /**
0040      * @brief Function that integrates the integrand with a given
0041      * relative accuracy using the method defined in the constructor.
0042      * @param xmin: the lower bound integration bound
0043      * @param xmax: the upper bound integration bound
0044      * @param eps: the required relative accuracy
0045      * @return the value of the integral
0046      */
0047     double integrate(double const& xmin, double const& xmax, double const& eps) const;
0048
0049     /**
0050      * @brief Function that integrates the integrand with a given
0051      * relative accuracy using a set of fixed point on the integration
0052      * range.
0053      * @param xmin: the lower bound integration bound
0054      * @param xmax: the upper bound integration bound
0055      * @param FixPts: the vector of fixed points of the integration
0056      * @param eps: the required relative accuracy
0057      * @return the value of the integral
0058      */
0059     double integrate(double const& xmin, double const& xmax, std::vector<double> const& FixPts, double const& eps) const;
0060
0061     /**
0062      * @brief Function that integrates the integrand with a given
0063      * relative accuracy using a set of fixed point on the integration
0064      * range. In this case, the vector also contains the lower and
0065      * upper bounds.
0066      * @param FixPts: the vector of fixed points of the integration
0067      * @param eps: the required relative accuracy
0068      * integration: 0 for 8/7 points for Gauss-Legendre/Kronrod, 1 for
0069      * 16/15 points for Gauss-Legendre/Kronrod (default: 1)
0070      * @return the value of the integral
0071      * @note In order to make this function as fast as possible, no
0072      * manipulation of the vector "FixPts" is made. To be used with
0073      * care.
0074      */
0075     double integrate(std::vector<double> const& FixPts, double const& eps) const;
0076
0077     /**
0078      * @brief Function that integrates the integrand with a given
0079      * relative accuracy using the method defined in the constructor.
0080      * @param xmin: the lower bound integration bound
0081      * @param xmax: the upper bound integration bound
0082      * @param n: index associated to the number of points used for the
0083      * integration: 0 for 8/7 points for Gauss-Legendre/Kronrod, 1 for
0084      * 16/15 points for Gauss-Legendre/Kronrod (default: 1)
0085      * @return the value of the integral
0086      */
0087     double integrate(double const& xmin, double const& xmax, int const& n = 1) const;
0088
0089     /**
0090      * @brief Function that integrates the integrand with a given
0091      * relative accuracy using a set of fixed point on the integration
0092      * range.
0093      * @param xmin: the lower bound integration bound
0094      * @param xmax: the upper bound integration bound
0095      * @param FixPts: the vector of fixed points of the integration
0096      * @param n: index associated to the number of points used for the
0097      * integration: 0 for 8/7 points for Gauss-Legendre/Kronrod, 1 for
0098      * 16/15 points for Gauss-Legendre/Kronrod (default: 1)
0099      * @return the value of the integral
0100      */
0101     double integrate(double const& xmin, double const& xmax, std::vector<double> const& FixPts, int const& n = 1) const;
0102
0103     /**
0104      * @brief Function that integrates the integrand with a given
0105      * relative accuracy using a set of fixed point on the integration
0106      * range. In this case, the vector also contains the lower and
0107      * upper bounds.
0108      * @param FixPts: the vector of fixed points of the integration
0109      * @param n: index associated to the number of points used for the
0110      * integration: 0 for 8/7 points for Gauss-Legendre/Kronrod, 1 for
0111      * 16/15 points for Gauss-Legendre/Kronrod (default: 1)
0112      * @return the value of the integral
0113      * @note In order to make this function as fast as possible, no
0114      * manipulation of the vector "FixPts" is made. To be used with
0115      * care.
0116      */
0117     double integrate(std::vector<double> const& FixPts, int const& n = 1) const;
0118
0119     /**
0120      * @brief Function for the integrand.
0121      * @param x: the integration variable
0122      * @return the integrand evaluated at x
0123      */
0124     double integrand(double const& x) const { return _func(x); };
0125
0126     /**
0127      * @brief Function that returns the integration method.
0128      */
0129     IntegrationMethod Method() const { return _method; };
0130
0131   private:
0132     /**
0133      * @brief Function that integrates the integrand using the Gauss-Legendre method.
0134      * @param xmin: the lower bound integration bound
0135      * @param xmax: the upper bound integration bound
0136      * @return a pair containing the value of the integral computed
0137      * with the 16-point method and the relative difference w.r.t. the
0138      * 8-point one.
0139      */
0140     std::pair<double, double> integrateGL(double const& xmin, double const& xmax) const;
0141
0142     /**
0143      * @brief Function that integrates the integrand using the Gauss-Kronrod method.
0144      * @param xmin: the lower bound integration bound
0145      * @param xmax: the upper bound integration bound
0146      * @return a pair containing the value of the integral computed
0147      * with the 15-point method and the relative difference w.r.t. the
0148      * 7-point one.
0149      */
0150     std::pair<double, double> integrateGK(double const& xmin, double const& xmax) const;
0151
0152     /**
0153      * @brief Function that integrates the integrand using the
0154      * Gauss-Legendre method (no accuracy estimate).
0155      * @param xmin: the lower bound integration bound
0156      * @param xmax: the upper bound integration bound
0157      * @param n: index associated to the number of points used for the
0158      * integration: 0 for 8 points for Gauss-Legendr, 1 for 16 points
0159      * @return the value of the integral
0160      */
0161     double integrateGL(double const& xmin, double const& xmax, int const& n) const;
0162
0163     /**
0164      * @brief Function that integrates the integrand using the
0165      * Gauss-Kronrod method (no accuracy estimate).
0166      * @param xmin: the lower bound integration bound
0167      * @param xmax: the upper bound integration bound
0168      * @param n: index associated to the number of points used for the
0169      * integration: 0 for 7 points for Gauss-Legendr, 1 for 15 points
0170      * @return the value of the integral
0171      */
0172     double integrateGK(double const& xmin, double const& xmax, int const& n) const;
0173
0174   private:
0175     std::function<double(double const&)> const _func;   //!< The integrand function
0176     IntegrationMethod                    const _method; //!< The integration method
0177   };
0178 }