Boost logo

Boost :

From: Kirill Kudashkin (kirill.kudashkin_at_[hidden])
Date: 2021-08-05 08:21:27


Dear all,

I have encountered a problem while using Boost C++ libraries (
multiprecision/mpc.hpp, math/quadrature) which I cannot resolve.

This problems appears when running the following piece of code

||

|#include <boost/math/constants/constants.hpp> #include
<boost/math/quadrature/exp_sinh.hpp> #include
<boost/multiprecision/mpc.hpp> #include <boost/multiprecision/mpfr.hpp>
#include <boost/math/special_functions/fpclassify.hpp> namespace mpns =
boost::multiprecision; typedef
mpns::number<boost::multiprecision::mpc_complex_backend <50> > mpc_type
; typedef mpc_type::value_type mp_type ; int main() { using
boost::math::constants::root_pi ; mpc_type z{mp_type(1),mp_type(1)} ;
auto f = [&](mp_type x) { //actually I did not test if all these
conditions are needed... if (boost::math::fpclassify<mp_type> (x) ==
FP_ZERO) { return mpc_type(mp_type(1)) ; }; if
(boost::math::fpclassify<mp_type> (x) == FP_INFINITE) { return
mpc_type(mp_type(0)) ; }; mp_type x2 = mpns::pow(x,2U) ; if
(boost::math::fpclassify<mp_type> (x2) == FP_ZERO) { return
mpc_type(mp_type(1)) ; }; if (boost::math::fpclassify<mp_type> (x2) ==
FP_INFINITE) { return mpc_type(mp_type(0)) ; }; mpc_type ex2 =
mpns::exp(-z*x2) ; if (boost::math::fpclassify<mpc_type> (ex2) ==
FP_ZERO) { return mpc_type(mp_type(0)) ; }; return ex2 ; } ; mp_type
termination = boost::math::tools::root_epsilon <mp_type> () ; mp_type
error; mp_type L1; size_t max_halvings = 12;
boost::math::quadrature::exp_sinh<mp_type> integrator(max_halvings) ;
mpc_type res = integrator.integrate(f,termination,&error,&L1) ; mpc_type
div = 2U*mpns::sqrt(z) ; mpc_type result = (mpc_type(root_pi<mp_type>
())/div) - res ; return 0; } |

||Namely, the problem appears when evaluating the exponent.

The code is compiled without any problems, but hangs on computing

|    mpc_type ex2 = mpns::exp(-z*x2) ;|

It was possible to go into the function while debugging.

It crashed eventually, when the compiler reached "unable to resolve
non-existing file.../src/init2.c".

I let it run for a while checking the call stack. It appears that the
programs is stuck atcalling ___gmpn_sub_n_coreisbr().

Or better, in the call stuck it appeared frequently.

I did several tests.

The real-valued version (with mpfr_backend <2000>) works without any
problems.

cpp_complex_backend <100> works fine, but I could not test higher
precision yet, since it is much slower.

I found a work-around which is not ideal. Namely,

mpc_type ex2 = mpc_type(mpns::exp(x2)) ;

    mpc_type ezx2 = mpns::pow(ex2,-z) ;

works fine with mpc_complex_backend <2000>.

The lambda-function as it is declared in the code above works and yields
correct numbers.

The problem appears when it is called from the integrator.

My system.

I use openSUSE Tumbeleweed.

I installed boost libraries as provided by
|boost_1_76_0-gnu-mpich-hpc-devel.|

Complier -> g++, cppStandard gnu+17.

I check that I have all libs, e.g. libgmp, libmpc, libmpfr and etc.

All headers are present (the programs is complied without any warnings).

Sorry for a lengthy question.

Thank you in advance and

Best wishes,

Kirill

||


Boost list run by bdawes at acm.org, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk