【算法】计算组合数的四种常用方法

【算法】计算组合数的四种常用方法

算法一：$C_a^b=C_{a-1}^{b-1}+C_{a-1}^b$

解析：

$C_a^b$可以理解为从CSDN网站a位大佬中选出b位大佬的总方案数，a位大佬中选取b位大佬的总方案数有两种情况：

• 情况一：选取“你”大佬；那么此时只需要从剩下的a-1位大佬中选出b-1位大佬
• 情况二：不选取“你”大佬，那么此时从剩下的a-1位大佬中选出b位大佬

根据加法计数原理有：$$C_a^b=C_{a-1}^{b-1}+C_{a-1}^b$$

时间复杂度为$O(n^2)$

代码如下：

#include<iostream>
using namespace std;
const int mod = 1e9+7;
long long f[2010][2010];
int main()
{
    //预处理
    for(int i=0;i<=2000;i++)
    {
        for(int j=0;j<=i;j++)
        {
            if(!j) f[i][j]=1;
            else f[i][j]=(f[i-1][j-1]+f[i-1][j])%mod;
        }
    }
    int n;
    cin>>n;
    while(n--)
    {
        int a,b;
        cin>>a>>b;
        cout<<f[a][b]<<endl;
    }
}

算法二：$C_a^b=\frac{a!}{b!(a-b)!}=a!\ast b{!}^{-1}\ast (a-b){!}^{-1}$

解析：

$b{!}^{-1}$和$(a-b){!}^{-1}$表示逆元，可以用快速幂算法求解

• $b{!}^{-1}=((b-1)!b{)}^{-1}=(b-1){!}^{-1}\ast b^{-1}$
• $b^{-1}=b^{1e9+5}(根据费马小定理)$

由上两式得：$$b{!}^{-1}=(b-1){!}^{-1}\ast b^{1e9+5}=(b-1){!}^{-1}\ast b^{mod-2}$$

时间复杂度为 $O(a\ast log(mod))$

代码如下：

#include<iostream>
using namespace std;
const int mod=1e9+7,N=1e5+10;
typedef long long LL;
long long fac[N];//存储数
long long infac[N];//存储逆元
//快速幂算法
int quick_pow(int a, int k, int p)
{
    int res = 1;
    while (k)
    {
        if (k & 1) res = (LL)res * a % p;
        a = (LL)a * a % p;
        k >>= 1;
    }
    return res;
}
int main()
{
    int n;
    // 预处理
    fac[0]=infac[0]=1;
    for(int i=1;i<=1e5;i++)
    {
        fac[i]=fac[i-1]*i%mod;
        infac[i]=(LL)infac[i - 1] * quick_pow(i,mod-2,mod)%mod;
    }
    cin>>n;
    while(n--)
    {
        int a,b;
        cin>>a>>b;
        cout<<(LL)fac[a] * infac[b] % mod * infac[a - b] % mod<<endl;
    }
}

算法三：卢卡斯定理$C_a^b(lucas)\equiv C_{\frac{a}{p}}^{\frac{b}{p}}{C}_{amodp}^{bmodp}(modp)$

代码如下：

#include<iostream>
using namespace std;
typedef long long LL;
int quickPow(int a,int  b,int p)
{
    int t=1;
    while(b)
    {
        if(b&1) t=(LL)t*a%p;
        a=(LL)a*a%p;
        b>>=1;
    }
    return t;
}
int C(int a,int b,int p)
{
    if(b>a) return 0;
    int res=1;
    for(int i=a,j=1;j<=b;i--,j++)
    {
        res=(LL)res*i%p;//求a * (a-1) * .... * (a-b+1)
        cout<<i<<endl;
        res=(LL)res*quickPow(j,p-2,p)%p;//b的逆元
    }
    return res;
}
int lucas(LL a, LL b, int p)
{
    if (a < p && b < p) return C(a, b, p);
    return (LL)C(a % p, b % p, p) * lucas(a / p, b / p, p) % p;
}

int main()
{
    int n;
    cin>>n;
    while(n--)
    {
        long long a,b,p;
        cin>>a>>b>>p;
        cout<<lucas(a,b,p)<<endl;
    }
}

算法四：高精度-素数组合

解析：

具体步骤：

• 筛素数（1～5000）
• 求每个质数的次数
• 用高精度把所有质因子乘上

代码如下：

#include<iostream>
#include<algorithm>
#include<vector>

using namespace std;

const int N=5010;

int primes[N],cnt;
int sum[N];
bool st[N];

void get_primes(int n)
{
    for(int i=2;i<=n;i++)
    {
        if(!st[i])primes[cnt++]=i;
        for(int j=0;primes[j]*i<=n;j++)
        {
            st[primes[j]*i]=true;
            if(i%primes[j]==0)break;//==0每次漏
        }
    }
}
// 对p的各个<=a的次数算整除下取整倍数
int get(int n,int p)
{
    int res =0;
    while(n)
    {
        res+=n/p;
        n/=p;
    }
    return res;
}
//高精度乘
vector<int> mul(vector<int> a, int b)
{
    vector<int> c;
    int t = 0;
    for (int i = 0; i < a.size(); i ++ )
    {
        t += a[i] * b;
        c.push_back(t % 10);
        t /= 10;
    }
    while (t)
    {
        c.push_back(t % 10);
        t /= 10;
    }
    // while(C.size()>1 && C.back()==0) C.pop_back();//考虑b==0时才有pop多余的0 b!=0不需要这行
    return c;
}

int main()
{
    int a,b;
    cin >> a >> b;
    get_primes(a);

    for(int i=0;i<cnt;i++)
    {
        int p = primes[i];
        sum[i] = get(a,p)-get(a-b,p)-get(b,p);//是a-b不是b-a
    }

    vector<int> res;
    res.push_back(1);

    for (int i = 0; i < cnt; i ++ )
        for (int j = 0; j < sum[i]; j ++ )//primes[i]的次数
            res = mul(res, primes[i]);

    for (int i = res.size() - 1; i >= 0; i -- )cout<<res[i]<<" ";
    return 0;
}