
I am trying to compute a 2 Dimensional numerical integration with the help of Boost. I've been reading the documentation for a couple of days and my CPP language knowledge is just up to legacy CPP level so I'm a noob. I saw this page <https://www.boost.org/doc/libs/1_73_0/libs/math/doc/html/math_toolkit/gauss_kronrod.html> in documentation which calculates the quadrature using* Gauss-Kronrod *rule. I was trying to build upon it and calculate a 2-d integral using Gauss-Kronrod rule. My understanding is *a 2D integral is basically 2 1D integrals with the Cartesian product.* So I was thinking If I can nest two integration function then I'm done. I've started from the aforementioned example and extended it a bit so the code for a 2d integration looks like, #include <boost/math/quadrature/gauss_kronrod.hpp> #include <iostream> int main(int argc, char *argv[]) { using namespace boost::math::quadrature; auto f1 = [](double t, double s) { return std::exp(-(t*t+s*s+t*s) / 2); }; double error; double Q = gauss_kronrod<double, 15>::integrate(gauss_kronrod<double, 61>::integrate(f1, 0, std::numeric_limits<double>::infinity(), 5), 0, std::numeric_limits<double>::infinity(), 5, 1e-9, &error); std::cout << Q << " " << error<<std::endl; return 0; } . I've just added a lambda function (which is the integrand) of two variables and nested 2 1D Gauss-Kronrod integration functions. I'm compiling this with g++ -Wall -I /path/to/boost_1_73_0 2d-gauss.cpp This gives error and compilation fails. I don't know what wrong I'm doing. How do you suggest to do such an integral? Thanks in advance.

Ritajit Kundu via Boost-users wrote:
I am trying to compute a 2 Dimensional numerical integration with the help of Boost. I've been reading the documentation for a couple of days and my CPP language knowledge is just up to legacy CPP level so I'm a noob.
I saw this page <https://www.boost.org/doc/libs/1_73_0/libs/math/doc/html/math_toolkit/gauss_kronrod.html> in documentation which calculates the quadrature using*Gauss-Kronrod *rule. I was trying to build upon it and calculate a 2-d integral using Gauss-Kronrod rule. My understanding is *a 2D integral is basically 2 1D integrals with the Cartesian product.* So I was thinking If I can nest two integration function then I'm done. I've started from the aforementioned example and extended it a bit so the code for a 2d integration looks like,
Hi Ritajit, I'm not up on combining the integrals as you would like. A math pro may come along and point it out. What you can't do is make up new functionality that the library does not support. What you will probably want to do is run with two independent gauss_kronrod calculations and combine them for your Cartesian product; code you would write on top of that. For instance, the signature,(what the integrate(..) call want to see), of the lambda is a function that has one parameter of a double and returns a double. Just as in the example: double (lambda){ return a double;} integrate(..) can not be made to work in two dimensions by changing the way you call it. It will only work on a 1 D integral.
double Q = gauss_kronrod<double, 15>::integrate(gauss_kronrod<double, 61>::integrate(f1, 0, std::numeric_limits<double>::infinity(), 5), 0, std::numeric_limits<double>::infinity(), 5, 1e-9, &error);
Look into C++ function signatures. Integrate looks like. static auto integrate(F f, Real a, Real b, unsigned max_depth = 15, Real tol = tools::root_epsilon<Real>(), Real* error = nullptr, Real* pL1 = nullptr).... I've talked about the signature of F, above you are passing a POD, double, the return value of integrate. This is what 'error: no matching function...' means Like I said, I don't know the 'math' answer here but at the least, this does compile. Best, Dan ................ #include <boost/math/quadrature/gauss_kronrod.hpp> #include <iostream> int main(int argc, char* argv[]) { using namespace boost::math::quadrature; auto f1 = [](double t) { return std::exp(-(t * t + 5 * 5 + t * 5) / 2); }; auto f2 = [](double s) { return std::exp(-(3 * 3 + s * s + 3 * s) / 2); }; double error; //for each across your matrix??? double X = gauss_kronrod<double, 15>::integrate(f1, 0, std::numeric_limits<double>::infinity(), 5, 1e-14, &error); double Y = gauss_kronrod<double, 15>::integrate(f1, 0, std::numeric_limits<double>::infinity(), 5, 1e-14, &error); std::cout << X << " " << Y << " " << error << std::endl; double error2; double Q = gauss_kronrod<double, 15>::integrate(f1, gauss_kronrod<double, 15>::integrate(f2, 0, std::numeric_limits<double>::infinity(), 5, 1e-14, &error2), std::numeric_limits<double>::infinity(), 5, 1e-14, &error); return 0; }
participants (2)
-
Dan Bloomquist
-
Ritajit Kundu