Загадочное исключение stackoverflow в boost::multiprecision::gmp_float::operator=?

Я пытаюсь рассчитать число пи для одного из своих университетских проектов, используя Рамануджана. формула для произвольного количества цифр после точки с плавающей запятой. Для работы я использую библиотеку boost::multiprecision, которая представляет собой просто оболочку вокруг mpfr и mpir, которые я уже установил на свой компьютер.

Пока все хорошо, но я либо что-то упускаю, либо делаю что-то не так, так как моя функция вычисления (та, которая вычисляет каждую итерацию суммы в формуле) время от времени выдает StackOverflow exception. (когда он начинается, он последователен)

Вот как это выглядит

boost::multiprecision::mpf_float calculatePi(int start, int end)
{
    boost::multiprecision::mpf_float partition = 0;
    for (; start < end; ++start)
    {
        boost::multiprecision::mpf_float n = 
            factorial(4 * start) * boost::multiprecision::mpf_float(1103 + 26390 * start);
        boost::multiprecision::mpf_float d = 
            boost::multiprecision::pow(factorial(start), 4) * pow((boost::multiprecision::mpf_float)396, 4 * start); <-- this is where stackoverflow exception is thrown

        partition += (n / d);
    }

    return (2 * boost::multiprecision::sqrt((boost::multiprecision::mpf_float)2) / 9801) * partition;
}

Я признаю, что я слишком зелен во всех вычислениях/преобразованиях произвольной точности и библиотеках, поэтому я могу что-то упустить.

Я не могу показать вам полный стек вызовов, потому что он слишком длинный, но после

ParallelPi.exe!boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>::**do_multiplies**<boost::multiprecision::detail::expression<boost::multiprecision::detail::function,boost::multiprecision::detail::number_kind_floating_pointpow_funct<boost::multiprecision::backends::gmp_float<0> >,boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>,int,void>,boost::multiprecision::detail::function>(const boost::multiprecision::detail::expression<boost::multiprecision::detail::function,boost::multiprecision::detail::number_kind_floating_pointpow_funct<boost::multiprecision::backends::gmp_float<0> >,boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>,int,void> & e, const boost::multiprecision::detail::function & __formal) Line 1754    C++

Затем следуют 3603 последовательных вызова

ParallelPi.exe!boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>::**operator=**<boost::multiprecision::detail::function,boost::multiprecision::detail::number_kind_floating_pointpow_funct<boost::multiprecision::backends::gmp_float<0> >,boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>,int,void>(const boost::multiprecision::detail::expression<boost::multiprecision::detail::function,boost::multiprecision::detail::number_kind_floating_pointpow_funct<boost::multiprecision::backends::gmp_float<0> >,boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>,int,void> & e) Line 216  C++
ParallelPi.exe!boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>::number<boost::multiprecision::backends::gmp_float<0>,1><boost::multiprecision::detail::function,boost::multiprecision::detail::number_kind_floating_pointpow_funct<boost::multiprecision::backends::gmp_float<0> >,boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>,int,void>(const boost::multiprecision::detail::expression<boost::multiprecision::detail::function,boost::multiprecision::detail::number_kind_floating_pointpow_funct<boost::multiprecision::backends::gmp_float<0> >,boost::multiprecision::number<boost::multiprecision::backends::gmp_float<0>,1>,int,void> & e, void * __formal) Line 324   C++

только для того, чтобы привести к исключениям stackoverflow через несколько вызовов позже в

ParallelPi.exe!boost::multiprecision::backends::gmp_float<0>::precision() Line 630  C++

Как ни странно, я не припомню, чтобы функция calculatePi менялась до того, как она «заработала».

Не могли бы вы, ребята, помочь мне расшифровать, что происходит, потому что я, кажется, потерялся?

Моя точность по умолчанию для boost::multiprecision::mpf_float составляет 150 цифр после с плавающей запятой, что раньше не было проблемой (я помню, что вычислил некоторое число с точностью 10 000 с плавающей запятой, а не SO или что-то еще)


Запросы

1.Код моего факториала

boost::multiprecision::mpf_float factorial(int n)
{
    boost::multiprecision::mpf_float fact = 1;
    for (int i = 1; i <= n; ++i)
        fact *= i;

    return fact;
}

person kuskmen    schedule 21.05.2019    source источник
comment
У меня нет ответа, но вы можете значительно улучшить производительность, повторно используя некоторые термины. [(n+1)!]^4 = (n+1)^4 * (n!)^4, 396^[4(n+1)] = 394^4 * 396^(4n). 150!^4 имеет более 1000 цифр, 394^(4n) более 1500. Возможно, это скроет проблему.   -  person Quimby    schedule 21.05.2019
comment
В настоящее время я работаю над тем, чтобы заставить работать самое тупое решение, а затем приступаю к оптимизации, тем не менее, спасибо за предложение, которое я буду иметь в виду, когда дойду до этого. Что касается переполнения цифр, раньше у меня не было этой проблемы, плюс это происходит на самой первой итерации цикла, когда start равен 0.   -  person kuskmen    schedule 21.05.2019
comment
@kuskmen: Где твое определение факториала?   -  person user14717    schedule 22.05.2019
comment
Вероятно, это не связано с вашей проблемой, но я предлагаю вам использовать mpfr_float вместо mpf_float. У него также есть то преимущество, что он уже предоставляет pi, поэтому вы можете сравнить его со своим результатом.   -  person Marc Glisse    schedule 22.05.2019
comment
@MarcGlisse, меняющий mpf_float на mpfr_float, устраняет проблему, и теперь я не получаю исключение ТАК, но я полагаю, что вопрос остается открытым, почему ТАК с mpf_float   -  person kuskmen    schedule 22.05.2019


Ответы (1)


Это ошибка последних версий Boost.Multiprecision. Сейчас я пытаюсь разобраться, но см. https://github.com/boostorg/multiprecision/issues/164.

Обратите внимание, что использование mpfr вместо mpf устранило бы проблему.

person John Maddock    schedule 07.11.2019