小朋友学算法（3）：求组合数

绝地灬酷狼 2022-05-23 08:28 208阅读 0赞

关于组合的介绍，可以参考[小朋友学奥数（12）：组合][12]

\#一、利用基本公式，递归

#include <iostream>
    using namespace std;
    typedef long long ll;
    
    ll combination(ll n, ll k)
    {
        if(0 == k)
        {
            return 1;
        }
    
        return combination(n, k - 1) * (n - k + 1)/k;
    }
    
    int main()
    {
        cout << combination(10, 3) << endl;
        return 0;
    }

运行结果：

120

分析：  
 C(10, 3)  
= C(10, 2） \* 8 / 3  
= C(10, 1) \* 9 \* 8 / (3 \* 2)  
= C(10, 0) \* 10 \* 9 \* 8 / (3 \* 2 \* 1)  
= 1 \* 10 \* 9 \* 8 / (3 \* 2 \* 1)  
= 120

\#二、利用基本性质，递归  
利用公式C(n, k) = C(n - 1, k) + C(n - 1, k - 1)

#include <iostream>
    using namespace std;
    typedef long long ll;
    
    ll combination(ll n, ll k)
    {
        if(0 == k || n == k)
        {
            return 1;
        }
    
        if(1 == k)
        {
            return n;
        }
    
        return combination(n - 1, k) + combination(n - 1, k - 1);
    }
    
    int main()
    {
        cout << combination(10, 3)<< endl;
        return 0;
    }

\#三、逆元+快速幂解法  
\#\#（一）基本概念  
上面两种方法都使用了递归方法，递归方法有个缺陷，就是在数据较大时效率较低。所以这里要介绍一个种新的求组合算法。在了解此算法之前，要先了解一些概念。

\#\#\#1 同余  
同余是数论中的重要概念。  
给定一个正整数m，如果两个整数a和b满足a-b能够被m整除，即(a-b)/m得到一个整数，那么就称整数a与b对模m同余，记作a≡b(mod m)

例1：4 ≡ 9 (mod 5)，即4和9对模5同余

例2：13 ≡ 23（mod 10)，即13和23对模10同余

\#\#\#2 模的加减乘除运算  
取模运算的等价变形适合加法、减法、乘法  
(a + b) % p = (a % p + b % p) % p  
(a - b) % p = (a % p - b % p) % p  
( a \* b) % p = (a % p \* b % p) % p

例3：(30 + 40) % 11 = 70 % 11 = 4  
(30% 11 + 40%11) % 11 = (8 + 7) % 11 = 15 % 11 = 4

例4：(80 - 20) % 7 = 60 % 7 = 4  
(80 % 7 - 20 % 7) % 7 = (3 - 6) % 7 = -3 % 7 = 4 （取模是让商尽可能小，所以这里有 -3 / 7 = -1 …… 4)

例5：(18 \* 20) % 7 = 360 % 7 = 3  
（18%7 \* 20%7）% 7 = （4 \* 6）% 7 = 3

但是，取模运算的等价变形不符合除法  
a/b % p ≠ (a%p / b%p) % p

例6：（100 % 20）% 11 = 5 % 11 = 5  
（100%11 / 20%11) % 11 = (1 / 9) % 11 = 0 % 11 = 0

\#\#\#3 逆元  
逆元：对于a和p，若gcd(a, p) = 1（a和p互素）且 a\*b%p≡1，则称b为a%p的逆元。

那这个逆元有什么用呢？试想一下求(a / b)%p，如果你知道b%p的逆元是c，那么就可以转变成(a/b)%p = (a/b) \* 1 % p = (a / b) \* (b\* c % p) % p = a\*c % p = (a%p) (c%p) % p，这样的话，除法运算就可以转化为乘法运算。

那怎么求逆元呢？这时候就要引入强大的费马小定理！

费马小定理：对于a和素数p，满足a^(p-1) % p ≡ 1

接着因为a^(p−1) = a^(p−2) \* a，所以有a^(p−2) \* a % p ≡ 1。  
对比逆元的定义可得，a^(p−2)就是a的逆元。  
所以问题就转换成求解a^(p−2)，即变成求快速幂的问题了。

\#\#\#4 快速幂  
这部分的内容可以参考 [小朋友学算法（6）：求幂pow函数的四种实现方式][6_pow] 中的第四种方法

\#\#（二）逆元 + 快速幂求组合思路  
现在目标是求C(n, m) %p，p为素数（经典p=1e9+7）。  
虽然有C(n, m) = n! / \[m! (n - m)!\]，但由于取模的性质对于除法不适用，则有

![1.png][]

所以需要利用逆元把“除法”转换成“乘法”，才能借助取模的性质计算组合数。

求解C(n, m)%p的步骤：

（1）通过循环，预先算好所有小于max\_number的阶乘（%p）的结果，存到fac\[max\_number\]里 (fac\[i\] = i! % p)  
（2）求m! % p的逆元（即求fac\[m\]的逆元）：根据费马小定理，x%p的逆元为x^(p−2)， 因此通过快速幂，求解fac\[m\]^(p−2) % p，记为M  
（3）求(n-m)! % p的逆元：同理就是求解fac\[n−m\]^(p−2) % p，记为NM  
（4）C(n, m) % p = ((fac\[n\] \* M) % p \* NM) % p

\#\#（三） 算法代码

using namespace std;
    
    #define MAX_NUMBER 100000
    //快速幂求x^n%mod
    long long quick_pow(long long x, long long n, long long mod)
    {
        long long res = 1;
        while (n)
        {
            if (n & 1)
            {
                res = res * x % mod;
            }
    
            x = x * x % mod;
            n >>= 1;
        }
    
        return res;
    }
    
    long long fac[MAX_NUMBER+5];
    long long n, m, p;
    
    int main()
    {
        while (~scanf("%lld %lld %lld", &n, &m, &p))
        {
            //预处理求fac，fac[i] = i!%p
            fac[0] = 1;
            for (int i = 1; i <= n; i++)
            {
                fac[i] = fac[i - 1] * i % p;
            }
            //组合数 = n!*(m!%p的逆元)*((n-m)!%p的逆元)%p
            printf("%lld\n", fac[n] * quick_pow(fac[m], p - 2, p) % p * quick_pow(fac[n - m], p - 2, p) % p);
        }
    }

运行结果：

10000 5000 100000007
    93446621

>TopCoder & Codeforces & AtCoder交流QQ群：648202993 >更多内容请关注微信公众号 !\[wechat\_public\_header.jpg\](https://img-blog.csdnimg.cn/img\_convert/640562211be0f192c77098248fb8b943.png)

[12]: http://www.jianshu.com/p/5237bae01960
[6_pow]: https://www.jianshu.com/p/1f4260ca6736
[1.png]: /images/20220523/1c094574521d4529b682177dbdd29619.png