（三）组合计数

（三）组合计数

$$\begin{gathered} 从a个元素中选择b个，有多少种取法C_a^b=\frac{a\times (a-1)\times \cdots \times (a-b+1)}{1\times 2\times 3\times \cdots \times b} \\ =\frac{a!}{b!\times (a-b)!} \\ =C_{a-1}^b+C_{a-1}^{b-1} \end{gathered}$$

组合计数 I

AcWing 885. 求组合数 I

实现思路：询问10000次，数字范围2000，直接使用递推式

• 递推式： $$C_a^b=C_{a-1}^b+C_{a-1}^{b-1}$$ 其实就是动态规划DP，O(N^2)

证明：

$$\begin{gathered} 从a个元素中存在某个元素，在选择b个元素过程中，共两种情况，包含这个元素和不包含这个元素 \\ 如果要选的b个元素包含这个元素，那么就是从剩下的a-1个元素中选b-1个，即C_{a-1}^{b-1}， \\ 如果不包含这个元素，就是从a-1个元素中选择b个，即C_{a-1}^b，综合两种情况就是C_a^b=C_{a-1}^b+C_{a-1}^{b-1} \end{gathered}$$

• 因为数字范围就2000，直接先预处理出2000范围内的所有组合数，使用二维数组存储c[][],后续只需调用就行
• 注意结果要取模，防止结果超出int范围

#include <iostream>
using namespace std;
const int N=2010,mod=1e9+7;
int c[N][N];
 
//预处理出这2000范围内的所有组合数
void init(){
    for(int i=0;i<N;i++)\
        for(int j=0;j<=i;j++)
            if(!j) c[i][j]=0;//若j为0
    		else c[i][j]=(c[i-1][j]+c[i-1][j-1])%mod;//注意取模，防止结果超出int范围
}
 
int mian(){
    int n;
    cin>>n;
    init();
    while(n--){
        int a,b;
        cin>>a>>b;
        cout<<c[a][b]<<endl;;
    }
    return 0;
}
 

组合计数 II

AcWing 886. 求组合数 II

实现思路：相比组合数I，这里数据范围增大到10^5，在预处理的过程中不能直接像组合数I一样求，否则时间复杂度过高10^10，这里使用公式+乘法逆元的方式预处理得到阶乘

$$公式C_a^b=\frac{a!}{b!\times (a-b)!}$$

• 考虑到计算过程存在取模且公式存在分式，则需要将除法转化为乘法，因为除法取模会得到错误答案，所以转化为乘法再做取模运算，这时就想到乘法逆元（因为本题模数1e9+7是质数，所以可以使用快速幂求逆元，费马定理b^(m-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 但对于除法却不成立，即(a / b) % p 不等于 (a % p / b % p) % p 。

  AcWing 886. 关于求组合数 II的一些疑问 - AcWing

• 设置一个数组fact[]，表示分子阶乘；数组infact[]，表示分母阶乘的逆元，最后通过fact[a]*infact[b]*infact[a-b]得到组合数结果

#include <iostream>
#include <algorithm>
using namespace std;
 
typedef long long LL;
const int N=1e4+10,mod=1e9+7;
 
LL fact[N],infact[N];//分子阶乘，和分母阶乘的逆元
 
//求快速幂
LL qmi(LL a,LL k,LL p){
    LL res=1;
    while(k){
        if(k & 1) res=(LL)res*a%p;
        k=>>1;
        a=(LL)a*a%p;
    }
    return res;
}
 
int main(){
    fact[0]=infact[0]=1;
    for(int i=1;i<N;i++){
        fact[i]=fact[i-1]*i%mod;//阶乘结果
        infact[i]=qmi(fact[i],mod-2,mod);//阶乘的逆元
    }
    int n;
    cin>>n;
    while(n--)
    {
        int a,b;
        cin>>a>>b;
        cout<<(LL)fact[a]*infact[b]%mod*infact[a-b]%mod<<endl;
    }
    return 0;
}

组合计数 III

AcWing 887. 求组合数 III

实现思路：询问很少，但数据范围巨大1e18，使用卢卡斯定理，时间复杂度为O(p*logN*logp)，1<p<1e5

• 卢卡斯（lucas）定理：

  $$C_a^b\equiv C_{amodp}^{bmodp}\times C_{a/p}^{b/p}(\mathrm{mod}p)$$

  $$\begin{gathered} 简要证明： \\ 把a和b转换为p进制表示 \\ a=a_k{p}^k+a_{k-1}{p}^{k-1}+\cdots +a_0{p}^0 \\ b=b_k{p}^k+b_{k-1}{p}^{k-1}+\cdots +b_0{p}^0 \\ 生成函数(1+x{)}^p=C_p^0\ast 1+C_p^1\ast x^1+C_p^2\ast x^2+\cdots +C_p^p\ast x^p\equiv 1+x^p(\mathrm{mod}p) \\ 所以有(1+x{)}^a=((1+x{)}^{p^0}{)}^{a_0}\times ((1+x{)}^{p^1}{)}^{a_1}\times ((1+x{)}^{p^2}{)}^{a_2}\times dots\times ((1+x{)}^{p^k}{)}^{a_k}=(1+x{)}^{a_0}\times (1+x^{p^1}{)}^{a_1}\times (1+x^{p^2}{)}^{a_2}\times \cdots \times (1+x^{p^k}{)}^{a_k} \\ 对比等式(1+x{)}^a=(1+x{)}^{a_0}\times (1+x^{p^1}{)}^{a_1}\times (1+x^{p^2}{)}^{a_2}\times \cdots \times (1+x^{p^k}{)}^{a_k}左右两边x^b项的系数，C_a^b\equiv C_{a_k}^{b_K}\times C_{a_{k-1}}^{b_{k-1}}\times \cdots \times C_{a_0}^{b_0}(\mathrm{mod}p) \end{gathered}$$

  lucas函数：函数的参数为a,b，若a<p且b<p则说明a,b足够小可以通过组合数函数C(a,b)直接计算，否则返回C（a%p,b%p）*lucas(a/p,b/p)（因为a,b除p可能依然很大，当足够小时就像上述情况直接返回组合数函数）

• 注意这里求组合数C，与上题不同，无需枚举处理所有组合数，直接用基本公式：

  $$C_a^b=\frac{a\ast (a-1)...\ast (a-b+1)}{b!}$$

​ 采用类似双指针算法，指针i由后向前遍历到a-b+1，j由1向前遍历到b

​ 然后再利用乘法逆元（快速幂求逆元），将除法转化为乘法，取模

#include <iostream>
#include <algorithm>
using namespace std;
typedef long long LL;
int p;//全局变量 后面函数的参数无需再传入p
 
//求快速幂
LL qmi(int a,int k){
    LL res=1;
    while(k){
        if(k & 1) res=(LL)res*a%p;
        k=>>1;
        a=(LL)a*a%p;
    }
    return res;
}
 
//基本公式求组合数（利用乘法逆元）
LL C (int a,int b){
    if(a<b) return 0;
    LL x=1,y=1;//x是分子，y是分母
    for(int i=a,j=1;j<=b;i--,j++){
        x=(LL)x*i%p;
        y=(LL)y*j%p;
    }
    return x*(LL)qmi(y,p-2,p)%p;//快速幂乘法逆元
}
 
//卢卡斯定理
LL lucas(int a,int b){
    if(a<p && b<p) return C(a,b);//a和b都小于p，直接用C求就可以
    return (LL)C(a%p,b%p)*C(a/p,b/p)%p;
}
 
 
int main()
{
    int n;
    cin>>n;
    while(n--)
    {
        LL a,b;
        cin>>a>>b>>p;
        cout<<lucas(a,b)<<endl;
    }
}

组合计数 IV

当我们需要求出组合数的真实值，而非对某个数的余数时，分解质因数的方式比较好用，同时需要用到高精度乘法

AcWing 888. 求组合数 IV

实现思路：

• 使用组合数公式

  $$C_a^b=\frac{a!}{b!\times (a-b)!}$$

• 用线筛法求出对于公式中最大的数a之前的所有的质数，之后在得到在a!的中所包含的每个质因子p对应的幂次

  $$a!中质因子p的幂次=\frac{a}{p}+\frac{a}{p^2}+.....\frac{a}{p^n}，向下取整$$

• 同理得到b!和(a-b)!中对应p的质因子的幂次，然后a!中p的幂次-b!中p幂次-(a-b)!中p的幂次=结果中p的幂次，以此类推得到结果中每个质因子pi对应的幂次

  $$最终结果C_a^b=p_1^{\alpha 1}\ast p_2^{\alpha 2}\ast p_3^{\alpha 3}\ast ....\ast p_k^{\alpha k}$$

• 再对结果进行一个高精度的乘法得到最终结果

#include <iostream>
#include <algorithm>
#include <vector>
using namespace std;
const int N=5010;
int primes[N],sum[N];//存储质数,结果中对应质因子的幂次
bool st[N];
int cnt;
 
//得到1~x中的质数   线筛法
void get_prime(int x){
    for(int i=2;i<=x;i++){
        if(!st[i]) primes[cnt++]=i;
        for(int j=0;primes[j]<=x/i;j++){
            st[primes[j]*i]=true;
            if(i%primes[j]==0) break;
        }
    }
}
 
//得到n的阶乘中质因子p的幂次
int get(int n,int p){
    int res=0;
    while(n){
        res+=n/p;
        n/=p;
    }
    return res;
}
 
//高精度乘法（第一章）
vector<int> mul(vector<int> &res,int b){
    int t=0;//表示进位
    vector<int> c;
    for(int i=0;i<res.size() || t;i++){
        if(i<res.size()) t+=res[i]*b;
        c.push_back(t%10);
        t/=10;
    }
    return c;
}
 
 
int mian(){
    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(b,p)-get(a-b,p);//得到结果对应质因子的幂次
    }
    vector<int> res;
    res.push_back(1);
    
    //循环高精度乘法得到结果
    for(int i=0;i<cnt;i++)
        for(int j=0;j<sum[i];j++)//对应幂次循环
            res=mul(res,primes[i]);
    
    //从数组末尾开始输出结果高位
    for(int i=res.size()-1;i>=0;i--) cou<<res[i];
}

卡特兰数

AcWing 889. 满足条件的01序列

实现思路：转化为求图中从原点到目标点的路径中符合某种条件的路径的方案数

假设有6个0，6个1，设0表示向右走一步，1表示向上走一步，则路径是从原点（0,0）到点(6,6)，需要12步。根据题目要求前缀中0的个数大于1的个数，转化为所走路径要满足的条件就是：路径上的每一点的坐标都必须满足横坐标x>=纵坐标y，即路径必须在图中红色斜线(y=x+1)的下方，不能与红色斜线有交叉。目标点（6,6）关于红色斜线对称的点为（5,7），则只要是从原点（0,0）走到（5,7）的路径必然会与红线有交叉，即不符合路径的要求。

（0,0）⇒ （6，6）的路径方案数：C_{12}^{6} （总方案）

（0,0）⇒ （5，7）的路径方案数：C_{12}^{5} (不符合条件的方案)

故符合条件的方案数C_{12}^{6} - C_{12}^{5}

$$\begin{gathered} 由以上推广到n个0和n个1，得到公式: \\ {C}_{2n}^n-C_{2n}^{n-1} \\ =\frac{C_{2n}^n}{n+1}，即卡特兰数 \end{gathered}$$

注：这里结果依旧要对一个质数取模，所以可以用费马定理，快速幂求逆元

#include <iostream>
#include <algorithm>
using namespace std;
const int mod=1e9+7;
typedef long long LL;
 
LL qmi(int a,int k,int p){
    LL res=1;
    while(k){
        if(k & 1) res=(LL)res*a%p;
        k=>>1;
        a=(LL)a*a%p;
    }
    return res;
}
 
//求卡特兰数   求组合数用基本公式法 和组合计数III一样 写法略有不同
LL Katelan(int n){
    int a=2*n,b=n;
    LL res=1;
    for(int i=a,j=1;j<=b;j++,i--){
        res=(LL)res*i%mod;
        res=(LL)res*qmi(j,mod-2,mod)%mod;
    }
    return (LL)res*qmi(n+1,mod-2,mod)%mod;//注意这里还有*(1/n+1)的逆元
}
 
int main(){
    int n;
    cin>>n;
    LL res=Katelan(n);
    cout<<res;
    return 0;
}