Факториальная функция через рекурсию с использованием R с Rcpp

Мой основной вопрос: почему результаты различаются для этих четырех реализаций факториальной функции и, более конкретно, почему функции начинают различаться при n = 13?

library(Rcpp)
cppFunction('   int facCpp(int n)
                {
                    if (n==0) return 1;
                    if (n==1) return 1;
                    return n*facCpp(n-1);
                }
            ')



cppFunction('   double fac2Cpp(int n)
                {
                    if (n==0) return 1;
                    if (n==1) return 1;
                    return n*fac2Cpp(n-1);
                }
            ')

cppFunction('   long int fac3Cpp(long int n)
                {
                    if (n==0) return 1;
                    if (n==1) return 1;
                    return n*fac3Cpp(n-1);
                }
            ')

c(factorial(12),prod(1:12),facCpp(12),fac2Cpp(12),fac3Cpp(12))
c(factorial(13),prod(1:13),facCpp(13),fac2Cpp(13),fac3Cpp(13))
c(factorial(20),prod(1:20),facCpp(20),fac2Cpp(20),fac3Cpp(20))
c(factorial(40),prod(1:40),facCpp(40),fac2Cpp(40),fac3Cpp(40))

Я понимаю, что вопрос, возможно, является дубликатом, поскольку ответы, вероятно, предлагаются здесь Rcpp, создание кадра данных с вектором long long, который также показывает, почему функции начинают различаться для f(13)

2^31-1>facCpp(12)
#> [1] TRUE
2^31-1>13*facCpp(12)
#> [1] FALSE


c(factorial(12),prod(1:12),facCpp(12),fac2Cpp(12),fac3Cpp(12))
#>[1] 479001600 479001600 479001600 479001600 479001600
c(factorial(13),prod(1:13),facCpp(13),fac2Cpp(13),fac3Cpp(13))
#> [1] 6227020800 6227020800 1932053504 6227020800 1932053504
c(factorial(20),prod(1:20),facCpp(20),fac2Cpp(20),fac3Cpp(20))
#> [1]  2.432902e+18  2.432902e+18 -2.102133e+09  2.432902e+18 -2.102133e+09

r recursion rcpp factorial
person Jesper Hybel Pedersen    schedule 17.04.2014    source источник
comment
12! не переполняет ни int, ни long int. 13! переполняет int, но не long int. 40! переполняет оба. То, что вы видите, является переполнением в вашем коде C++.   -  person josliber♦    schedule 17.04.2014
comment
@josilber спасибо за ответ, который проясняет часть, связанную с C++. Но если вы запустите код, вы заметите, что fac3Cpp(13) также возвращает неправильный результат. Итак, я предполагаю, что это связано с тем, что пакет Rcpp не поддерживает long int, поэтому fac3Cpp становится facCpp.   -  person Jesper Hybel Pedersen    schedule 17.04.2014
comment
@user2055639 fac3Cpp(13) возвращает правильный результат на моей машине. Однако стандарт С++ допускает эту разницу.   -  person Dason    schedule 17.04.2014
comment
Как и @Dason, fac3Cpp(13) дает правильный результат на моем компьютере. Пожалуйста, обновите свой пост, указав результат, который вы видите, чтобы получить помощь по устранению несоответствия, которое вы видите.   -  person josliber♦    schedule 17.04.2014
comment
@Dason ad fac3Cpp(13) возвращает правильный результат на моей машине, это меня озадачивает? Однако стандарт С++ допускает эту разницу. не уверен, что я понимаю, что вы стремитесь здесь ?? Думаю, в конце концов я ищу несколько общее рабочее правило, когда не использовать int. Так что я думаю, что в целом нельзя быть уверенным, что long int решит проблему. ... Я обновлю свой пост выводом.   -  person Jesper Hybel Pedersen    schedule 17.04.2014

Ответы (1)


arrow_upward
3
arrow_downward

Вы по сути делаете это неправильно. Смотрите страницу справки R для факториала:

‘factorial(x)’ (x! для неотрицательного целого числа ‘x’) определяется как ‘gamma(x+1)’, а ‘lfactorial’ как ‘lgamma(x+1)’.

Вы не должны вычислять это таким образом. Почему? Ну посмотри на это:

R> evalCpp("INT_MAX")
[1] 2147483647
R>

Вы столкнетесь с числовым переполнением. Отсюда и другой алгоритм, реализованный, например, в функции factorial() R, которая просто выполняет gamma(x+1). И вы можете сделать это и на C++:

R> cppFunction('double myFac(int x) { return(R::gammafn(x+1.0)); }')
R> myFac(4)
[1] 24
R> myFac(12)
[1] 479001600
R> myFac(13)
[1] 6227020800
R> myFac(20)
[1] 2.4329e+18
R> myFac(40)
[1] 8.15915e+47
R>
person Dirk Eddelbuettel    schedule 17.04.2014
comment
@Eddelbuettel, спасибо за ответ, важно видеть, как вы должны это делать. На самом деле я не стремился к правильной реализации факториала, а просто дурачился, чтобы начать использовать Rcpp. Если у вас есть какие-либо мысли о том, почему использование long int в fac3Cpp возвращает неправильный ответ, пожалуйста, не сдерживайтесь? - person Jesper Hybel Pedersen; 18.04.2014
comment
У меня не было времени подробно разбираться, и я бы заподозрил некоторую зависимость от версии ОС и компилятора, поскольку мы все знаем, каким должен быть правильный ответ. Я бы настроил алгоритмы таким образом, чтобы вы не сталкивались с недоливом или переполнением. [Кстати: это Eddelbuettel с двумя D] - person Dirk Eddelbuettel; 18.04.2014
comment
@user2055639 user2055639 сравните evalCpp("LONG_MAX") и factorial(21) - вы все равно будете переполняться после 20. double может хранить гораздо больший диапазон значений. - person Kevin Ushey; 18.04.2014