寻找大数的组合值答案

【问题标题】：To find combination value of large numbers寻找大数的组合值
【发布时间】：2014-06-27 06:00:08
【问题描述】：

我想为大整数找到（n 选择 r），我还必须找出那个数字的模。

long long int choose(int a,int b)
{
    if (b > a)
        return (-1);
    if(b==0 || a==1 || b==a)
        return(1);
    else
    {
        long long int r = ((choose(a-1,b))%10000007+(choose(a-1,b-  1))%10000007)%10000007;
        return r;
    }
}

我正在使用这段代码，但我得到了 TLE。如果有其他方法可以做到这一点，请告诉我。

【问题讨论】：

这是为了一些代码挑战吗？常量10000007 看起来有些眼熟：stackoverflow.com/questions/14189713/…
@YuHao Time Length Exceeded
计算 binom.coeff 的效率非常低。你说你需要计算大数——你有没有试过估计你的choose()必须被调用多少次来计算大n的choose(n,n/2)？
你确定你有正确的模数吗？ 100000007 是一个更常见的模数编程竞赛，因为您经常希望利用它是素数的事实。使用递归的记忆版本，您的运行时间将是O(a*b)。如果模数是素数 p，则可以在O((a+b)*log p) 中求解。
告诉a和b的约束

标签： c++ combinations

【解决方案1】：

我还没有评论的声誉，但我想指出，rock321987 的回答效果很好：

快速且正确直到并包括 C(62, 31)

但无法处理所有输出适合 uint64_t 的输入。作为证据，请尝试：

C(67, 33) = 14,226,520,737,620,288,370 (verify correctness and size)

不幸的是，另一个实现吐出不正确的 8,829,174,638,479,413。还有其他计算 nCr 的方法不会像这样中断，但这里真正的问题是没有尝试利用模数。

请注意 p = 10000007 是素数，这使我们能够利用所有整数都具有逆模 p 的事实，并且该逆是唯一的。此外，我们可以很快找到逆。另一个问题有一个关于如何做到这一点的答案here，我在下面复制了。

这很方便，因为：

x/y mod p == x*(y inverse) mod p;和
xy mod p == (x mod p)(y mod p)

稍微修改一下其他代码，将问题概括如下：

#include <iostream>
#include <assert.h>

// p MUST be prime and less than 2^63
uint64_t inverseModp(uint64_t a, uint64_t p) {
    assert(p < (1ull << 63));
    assert(a < p);
    assert(a != 0);
    uint64_t ex = p-2, result = 1;
    while (ex > 0) {
        if (ex % 2 == 1) {
            result = (result*a) % p;
        }
        a = (a*a) % p;
        ex /= 2;
    }
    return result;
}

// p MUST be prime
uint32_t nCrModp(uint32_t n, uint32_t r, uint32_t p)
{
    assert(r <= n);
    if (r > n-r) r = n-r;
    if (r == 0) return 1;
    if(n/p - (n-r)/p > r/p) return 0;

    uint64_t result = 1; //intermediary results may overflow 32 bits

    for (uint32_t i = n, x = 1; i > r; --i, ++x) {
        if( i % p != 0) {
            result *= i % p;
            result %= p;
        }
        if( x % p != 0) {
            result *= inverseModp(x % p, p);
            result %= p;
        }
    }
    return result;
}

int main() {
    uint32_t smallPrime = 17;
    uint32_t medNum = 3001;
    uint32_t halfMedNum = medNum >> 1;
    std::cout << nCrModp(medNum, halfMedNum, smallPrime) << std::endl;

    uint32_t bigPrime = 4294967291ul; // 2^32-5 is largest prime < 2^32
    uint32_t bigNum = 1ul << 24;
    uint32_t halfBigNum = bigNum >> 1;
    std::cout << nCrModp(bigNum, halfBigNum, bigPrime) << std::endl;
}

如果您愿意等待，这应该会为任何一组 32 位输入产生结果。为了证明这一点，我已经包含了 24 位 n 和最大 32 位素数的计算。我的普通电脑花了大约 13 秒来计算这个。对照wolfram alpha检查答案，但要注意它可能会超过那里的“标准计算时间”。

如果 p 远小于 (n-r) 其中 r

【讨论】：

【解决方案2】：

nCr = n! / (r! * (n-r)!) {! = 阶乘}

现在选择r or n - r，使它们中的任何一个都是最小值

#include <cstdio>
#include <cmath>

#define MOD 10000007

int main()
{
    int n, r, i, x = 1;
    long long int res = 1;
    scanf("%d%d", &n, &r);
    int mini = fmin(r, (n - r));//minimum of r,n-r

    for (i = n;i > mini;i--) {
        res = (res * i) / x;
        x++;
    }
    printf("%lld\n", res % MOD);
    return 0;
}

如果n 和r 的值不太高

，它将适用于编程比赛所需的大多数情况

时间复杂度：- O(min(r, n - r))

限制 :- 对于 C/C++ 等语言，如果

n > 60（大约）

因为没有数据类型可以存储最终值..

【讨论】：

中心二项式系数增长exponentially。所以，res 会很快溢出。你不能简单地在最后做res%MOD 并期望它起作用。
这就是为什么我提到n 和k 的约束不是太高的原因

【解决方案3】：

nCr 的展开式总是可以简化为整数的乘积。这是通过取消分母中的项来完成的。这种方法应用在下面给出的函数中。

这个函数的时间复杂度为 O(n^2 * log(n))。这将在 1 秒内计算 n

#include <numeric>
#include <algorithm>
int M=1e7+7;
int ncr(int n, int r)
{
    r=min(r,n-r);
    int A[r],i,j,B[r];
    iota(A,A+r,n-r+1);  //initializing A starting from n-r+1 to n
    iota(B,B+r,1);      //initializing B starting from 1 to r

    int g;
    for(i=0;i<r;i++)
    for(j=0;j<r;j++)
    {
        if(B[i]==1)
            break;
        g=__gcd(B[i], A[j] );
        A[j]/=g;
        B[i]/=g;
    }
    long long ans=1;
    for(i=0;i<r;i++)
        ans=(ans*A[i])%M;
    return ans;
}

【讨论】：