I've found two good libraries for numerical integration in c++ that are easy enough for me to use. I have to integrate fairly complex functions with the possibility of singularities at endpoints and it appears that an adaptive integration method using the Wynn's epsilon algorithm is the best choice. Both libraries have an implementation for this (gsl_integration_qags and boost::numeric::quadrature::adaptive().accelerator(wynn_epsilon_algorithm)) but I can't decide which one is faster. I wonder if there are any experts who have experience with this and can help me decide. The computational cost of numerical integration is critical for my current project. (I've done some tests and found GSL to perform faster here but I'm not sure if my testing is conclusive. In general I've found Boost faster for mathematical methods.)