|
|
## $C_a^b$的多种场景下的求法
|
|
|
|
|
|
### [一、$AcWing$ $885$. 求组合数 $I$](https://www.acwing.com/problem/content/887/)
|
|
|
|
|
|
理论依据：$\large C_a^b=C_{a-1}^b+C_{a-1}^{b-1}$
|
|
|
|
|
|
> <font color='red'><b>限制条件：
|
|
|
> $\large a<=2000,b<=2000$</b></font>
|
|
|
>
|
|
|
> **受限原因**：需要二维数组，$2000*2000=4000000=4e6$,占据太大的空间。
|
|
|
|
|
|
|
|
|
|
|
|
**证明:**
|
|
|
有$a$个苹果，现在需要选出$b$个苹果,一共有多少种选法呢？
|
|
|
**答**：走到第一个苹果面前，面临两个选择：
|
|
|
|
|
|
- 选择它
|
|
|
- 放弃它
|
|
|
|
|
|
如果选择它，将要面对$a-1$个苹果中选$b-1$个苹果的问题
|
|
|
|
|
|
如果放弃它，将要面对$a-1$个苹果中选$b$个苹果的问题
|
|
|
|
|
|
根据加法原理，两种选择方法加在一起，就是选法总数，也就是上面的递推式。
|
|
|
|
|
|
```c++
|
|
|
#include <bits/stdc++.h>
|
|
|
|
|
|
using namespace std;
|
|
|
const int N = 2010;
|
|
|
const int p = 1e9 + 7;
|
|
|
int c[N][N];
|
|
|
|
|
|
void init() {
|
|
|
/*
|
|
|
功能：第0列初始化为1，现实意义：C(n,0)=1
|
|
|
理解：从n个元素中,每次取出0个元素,也就是每次都不取出元素,有几种取法?
|
|
|
很明显只有1种取法,就是什么都不取这样一种情况。类似:C(n,0)=C(n,n)=1
|
|
|
这是递归的base case,如果没有这个1打底，以后再怎么加法原理都是白费
|
|
|
*/
|
|
|
for (int i = 1; i < N; i++) c[i][0] = 1, c[i][i] = 1; // base case
|
|
|
|
|
|
for (int i = 1; i < N; i++)
|
|
|
for (int j = 1; j < i; j++)
|
|
|
// 递归式最重要！这是根据数学公式得到的递推式
|
|
|
// 从递推式由右向左看过去，需要由小到大填充i,j,并且，i>j
|
|
|
c[i][j] = (c[i - 1][j] + c[i - 1][j - 1]) % p;
|
|
|
}
|
|
|
|
|
|
int n;
|
|
|
|
|
|
int main() {
|
|
|
init();
|
|
|
cin >> n;
|
|
|
while (n--) {
|
|
|
int a, b;
|
|
|
cin >> a >> b;
|
|
|
cout << c[a][b] << endl;
|
|
|
}
|
|
|
return 0;
|
|
|
}
|
|
|
```
|
|
|
|
|
|
|
|
|
|
|
|
### [二、$AcWing$ $886$. 求组合数 $II$](https://www.acwing.com/problem/content/888/)
|
|
|
|
|
|
上个办法把$C_a^b$的值预处理出来,用的是$c[N][N]$,$N$最大$2010$.
|
|
|
本题的$a$,$b$都是上限$10^5$,如果按上题来，就是$c[10^5][10^5]$,
|
|
|
直接报 $Memory\ Limit\ Exceeded$. 所以不能直接递推求出所有解。
|
|
|
|
|
|
> <font color='red'><b>限制条件：</b>
|
|
|
> $a<=1e5,b<=1e5$</font>
|
|
|
|
|
|
理论依据：$\large C_a^b=\frac{a!}{(a-b)! * b!}$
|
|
|
|
|
|
举栗子：
|
|
|
$\large C_3^2=\frac{3 \times 2}{2 \times 1}=3$
|
|
|
|
|
|
|
|
|
|
|
|
#### 前导知识
|
|
|
|
|
|
##### 模逆元
|
|
|
|
|
|
定义与用途
|
|
|
|
|
|
若$b,p$互质，$b \times b^{-1} \equiv 1 \ (mod \ p)$,称$b^{-1}$是$b$的模$p$逆元。
|
|
|
|
|
|
有 $$a \equiv a \ (mod \ p) \ \equiv a \times 1 \ (mod \ p) \equiv a \times (b \times b^{-1}) \ (mod \ p)$$
|
|
|
|
|
|
若$b \mid a$,则$\frac{a}{b} \equiv a \times b^{-1} \ (mod \ p)$
|
|
|
|
|
|
|
|
|
|
|
|
##### 阶乘模逆元
|
|
|
|
|
|
若
|
|
|
|
|
|
$a \times a^{-1} \equiv 1 \ (mod \ p) \\ b \times b^{-1} \equiv 1 \ (mod \ p)$
|
|
|
|
|
|
则
|
|
|
|
|
|
$$ab \times a^{-1} b^{-1} \equiv 1 \ (mod \ p)$$
|
|
|
|
|
|
根据逆元定义，则
|
|
|
|
|
|
|
|
|
|
|
|
$$(ab)^{-1} \equiv a^{-1}b^{-1} \ (mod \ p)$$
|
|
|
|
|
|
|
|
|
|
|
|
因此
|
|
|
|
|
|
$\large \displaystyle (n!)^{-1}= \prod_{i=1}^{n}i^{-1} \ (mod \ p)$
|
|
|
|
|
|
|
|
|
|
|
|
本题要求结果 $mod\ (1e9+7)$, **模是质数(~~一般取模的都是质数~~)**，如果有除法，可以直接使用 **费马小定理+快速幂** 求逆元转化为乘法求解。
|
|
|
|
|
|
现在要求计算出 $a!$, ${b!}^{-1}$, ${(a-b)!}^{-1}$，可以使用两个数组来递推得到：
|
|
|
|
|
|
$\large fact[i] = i!\ \% \ p$
|
|
|
|
|
|
$\large infact[i] = (i!)^{-1}\ \% \ p = ((i-1)!^{-1}*i^{-1}\ ) \% \ p$
|
|
|
|
|
|
> **解释**：上面推导的阶乘模逆元
|
|
|
|
|
|
所以： $C_a^b=\frac{a!}{(a-b)! * b!} = (fact[a] * infact[a-b] * infact[b] )\ \% \ p$
|
|
|
|
|
|
|
|
|
#### C++ 代码
|
|
|
```cpp
|
|
|
#include <bits/stdc++.h>
|
|
|
|
|
|
using namespace std;
|
|
|
#define int long long
|
|
|
#define endl "\n"
|
|
|
const int N = 100010;
|
|
|
const int p = 1e9 + 7;
|
|
|
|
|
|
int fact[N]; // 阶乘结果数组
|
|
|
int infact[N]; // 阶乘逆元结果数组
|
|
|
|
|
|
// 快速幂
|
|
|
int qmi(int a, int k) {
|
|
|
int res = 1;
|
|
|
while (k) {
|
|
|
if (k & 1) res = res * a % p;
|
|
|
a = a * a % p;
|
|
|
k >>= 1;
|
|
|
}
|
|
|
return res;
|
|
|
}
|
|
|
|
|
|
signed main() {
|
|
|
// base case
|
|
|
fact[0] = infact[0] = 1;
|
|
|
|
|
|
for (int i = 1; i < N; i++) {
|
|
|
fact[i] = fact[i - 1] * i % p;
|
|
|
infact[i] = infact[i - 1] * qmi(i, p - 2) % p;
|
|
|
}
|
|
|
|
|
|
int n;
|
|
|
cin >> n;
|
|
|
while (n--) {
|
|
|
int a, b;
|
|
|
cin >> a >> b;
|
|
|
printf("%lld\n", fact[a] * infact[b] % p * infact[a - b] % p);
|
|
|
}
|
|
|
}
|
|
|
```
|
|
|
|
|
|
### [$AcWing$ $887$. 求组合数 $III$](https://www.acwing.com/problem/content/889/)
|
|
|
|
|
|
<font color='red'><b>
|
|
|
数据范围：</b>$b<=a<=1e18$,$p<=1e5$,$p$是质数 </font>
|
|
|
|
|
|
本题中$a,b$都太大了，如果像上一题一样预处理出$fact[],infact[]$,直接就炸了，肯定行不通！
|
|
|
|
|
|
还好，它的$p$比较小，而且是质数，我们有一种武器叫$Lucas$公式：
|
|
|
|
|
|
#### 1.$Lucas$公式
|
|
|
$\large C_a^b \equiv C_{a\%p}^{b\%p} * C_{a \div p}^{b \div p} (mod \ p)$
|
|
|
|
|
|
> <font color='red'><b>
|
|
|
> 注意：
|
|
|
> 这个$p$是小于$1e5$的，看到要求模$1e9+7$可不能用$Lucas$定理啊！不是干那个的啊！适合$a,b$很大，$p$很小的场景啊！！</b></font>
|
|
|
|
|
|
#### 2.此种情况下如何求解$C_a^b$
|
|
|
|
|
|
$C_a^b=\frac{a!}{(a-b)! \times b!}=\frac{a \times (a-1) \times (a-2) \times ...\times(a-b+1) \times (a-b) \times ... \times 1}{(a-b) \times (a-b-1) \times ... \times 1 \times b!}= \frac{a \times (a-1) \times (a-2) \times ...(a-b+1)}{b!}=\frac{a \times (a-1) \times (a-2) \times ...(a-b+1)}{b \times (b-1) \times (b-2) \times ... \times 2 \times 1}$
|
|
|
|
|
|
根据这个变形:
|
|
|
- 对于$a$来讲，需要变量$j$从$a$一直遍历到$a-b+1$
|
|
|
- 对于$b$来讲，需要变量$i$从$1$一起遍历到$b$
|
|
|
|
|
|
从$a$倒着走到$a-b+1$是$b$次($a>=b$)，比如$10$遍历到$8$，就是$10,9,8$三次，即$b=3$,也就是$10-3+1=8$。
|
|
|
而$i$也是遍历$b$次，太巧了！它们两个可以在一个$b$次的循环中一起变动 (代码能省则省啊~)！
|
|
|
|
|
|
``` c++
|
|
|
int C(int a, int b) {
|
|
|
if (a < b) return 0;
|
|
|
int down = 1, up = 1;
|
|
|
for (int i = a, j = 1; j <= b; i--, j++) {
|
|
|
up = up * i % p;
|
|
|
down = down * j % p;
|
|
|
}
|
|
|
return up * qmi(down, p - 2, p) % p;
|
|
|
}
|
|
|
|
|
|
```
|
|
|
|
|
|
#### $Code$
|
|
|
```cpp
|
|
|
#include <bits/stdc++.h>
|
|
|
|
|
|
using namespace std;
|
|
|
#define int long long
|
|
|
#define endl "\n"
|
|
|
|
|
|
// 快速幂
|
|
|
int qmi(int a, int k, int p) {
|
|
|
int res = 1;
|
|
|
while (k) {
|
|
|
if (k & 1) res = res * a % p;
|
|
|
a = a * a % p;
|
|
|
k >>= 1;
|
|
|
}
|
|
|
return res;
|
|
|
}
|
|
|
|
|
|
// 功能：组合数模板
|
|
|
int C(int a, int b, int p) {
|
|
|
if (a < b) return 0;
|
|
|
int down = 1, up = 1;
|
|
|
|
|
|
for (int i = a, j = 1; j <= b; i--, j++) {
|
|
|
up = up * i % p;
|
|
|
down = down * j % p;
|
|
|
}
|
|
|
return up * qmi(down, p - 2, p) % p;
|
|
|
}
|
|
|
|
|
|
// Lucas公式
|
|
|
int lucas(int a, int b, int p) {
|
|
|
if (a < p && b < p) return C(a, b, p);
|
|
|
return C(a % p, b % p, p) * lucas(a / p, b / p, p) % p; // 递归套用公式
|
|
|
}
|
|
|
|
|
|
int n, p;
|
|
|
|
|
|
signed main() {
|
|
|
cin >> n;
|
|
|
while (n--) { // n组询问
|
|
|
int a, b;
|
|
|
cin >> a >> b >> p;
|
|
|
cout << lucas(a, b, p) << endl; // 利用lucas公式计算组合数
|
|
|
}
|
|
|
return 0;
|
|
|
}
|
|
|
```
|
|
|
|
|
|
|
|
|
### [四、$AcWing$ $888$. 求组合数 $IV$](https://www.acwing.com/problem/content/890/)
|
|
|
|
|
|
**不让取模怎么办?**
|
|
|
|
|
|
#### 1、与前三道题的区别
|
|
|
不再是数据上限的问题了，而是一直不$mod$,有多大都保留，这无疑可能会爆$long\ long$,需要使用高精度。
|
|
|
|
|
|
#### 2、公式一【从定义出发的组合数公式】：
|
|
|
(1) $\large C_a^b=\frac{a \times (a-1) \times ... \times(a-b+1)}{b\times (b-1) \times ... \times 1}=\frac{a \times (a-1) \times ... \times (a-b+1) \times (a-b)!}{ b \times (b-1) \times ... \times 1 \times (a-b)!} =\frac{a!}{b! \times (a-b)!}$
|
|
|
|
|
|
#### 3、是不是我们利用高精度+上面的组合数公式直接算就行了？
|
|
|
不是的，因为$yxc$(大雪菜老师)说，这样的效率太低，不可以，需要再想一个好办法。
|
|
|
|
|
|
#### 4、公式二 【算术基本定理】
|
|
|
算术基本定理:$\large C_a^b=p_1^{\alpha_1} \times p_2^{\alpha_2} \times p_3^{\alpha_3} ... \times p_k^{\alpha_k}$
|
|
|
其中$p_1$,$p_2$...是质数，$\alpha 1$,$\alpha 2$,...是指质数$p_1$,$p_2$,...的个数。 如果能分解质因数成功的话，那么就可以通过 **高精度乘法** 解决掉这个问题。
|
|
|
|
|
|
(1) 那么$p_1$和$p_2$这些东东怎么求？
|
|
|
**思路**：因为$a,b$都是在$[1,5000]$之间的数，所以可以通过筛质数的方法，提前获取到一个质数数组，然后逐个去看，是不是含有这个质数。就能知道有哪些$p_1,p_2,...$了。
|
|
|
|
|
|
(2) $\alpha_1$,$\alpha_2,...,\alpha_k$又该怎么求呢？
|
|
|
|
|
|
#### 5、公式三【求质数$p$在$a!$中出现的次数】
|
|
|
|
|
|
$\large cnt=\lfloor \frac{a}{p} \rfloor + \lfloor \frac{a}{p^2} \rfloor + \lfloor \frac{a}{p^3} \rfloor + ...+ \lfloor \frac{a}{p^k} \rfloor$
|
|
|
|
|
|
[参考资料](https://blog.csdn.net/spidy_harker/article/details/88414504)
|
|
|
|
|
|
计算$n$的阶乘中质数$p$的个数：
|
|
|
```cpp {.line-numbers}
|
|
|
int get(int n, int p) {
|
|
|
int res = 0;
|
|
|
while (n) {
|
|
|
res += n / p;
|
|
|
n /= p;
|
|
|
}
|
|
|
return res;
|
|
|
}
|
|
|
```
|
|
|
|
|
|
**举个栗子**：
|
|
|
|
|
|
$ 5!=1 \times 2 \times 3 \times 4 \times 5$,求$2$的个数，
|
|
|
则$1 \sim 5$之间含$2$的个数就是 $\lfloor \frac{5}{2} \rfloor$,就是$2$个。
|
|
|
则$1 \sim 5$之间含$2^2$的个数就是 $\lfloor \frac{5}{2^2} \rfloor$,就是$1$个。
|
|
|
$2^3$就大于$5$，就不再继续。那么一共的个数就是$2+1=3$个。表示在 $5!$这个数字中，存在$3$个$2$，就是$2$中有$1$个$2$，$4$中有两个$2$。
|
|
|
|
|
|
#### 6、算法流程
|
|
|
(1)、筛素数，把$1-5000$之间的素数筛出来。
|
|
|
(2)、计算$C_a^b$中每个已求得素数的个数。
|
|
|
(3)、利用高精度，计算$C_a^b$$=p_1^{\alpha1} \times p_2^{\alpha2} \times p_3^{\alpha3} ... \times p_k^{\alpha k}$的值。
|
|
|
|
|
|
#### $Code$
|
|
|
|
|
|
```cpp
|
|
|
#include <bits/stdc++.h>
|
|
|
|
|
|
using namespace std;
|
|
|
const int N = 5010;
|
|
|
|
|
|
int primes[N], cnt;
|
|
|
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] <= n / i; j++) {
|
|
|
st[primes[j] * i] = true;
|
|
|
if (i % primes[j] == 0) break;
|
|
|
}
|
|
|
}
|
|
|
}
|
|
|
|
|
|
// 高精乘低精
|
|
|
void mul(int a[], int &al, int b) {
|
|
|
int t = 0;
|
|
|
for (int i = 1; i <= al; i++) {
|
|
|
t += a[i] * b;
|
|
|
a[i] = t % 10;
|
|
|
t /= 10;
|
|
|
}
|
|
|
while (t) {
|
|
|
a[++al] = t % 10;
|
|
|
t /= 10;
|
|
|
}
|
|
|
}
|
|
|
|
|
|
/**
|
|
|
* 功能：n的阶乘中包含的质因子p的个数
|
|
|
* @param n n的阶乘
|
|
|
* @param p 质因子p
|
|
|
* @return 有多少个
|
|
|
*/
|
|
|
int get(int n, int p) {
|
|
|
int cnt = 0;
|
|
|
while (n) { // p^1的个数，p^2的个数,p^3的个数...
|
|
|
cnt += n / p;
|
|
|
n /= p;
|
|
|
}
|
|
|
return cnt;
|
|
|
}
|
|
|
|
|
|
// C(a,b)的结果，高精度保存到c数组，同时，返回c数组的长度len
|
|
|
void C(int a, int b, int c[], int &cl) {
|
|
|
// 乘法的基数是1
|
|
|
c[1] = 1, cl = 1; // 由于高精度数组中只有一位，是1，所以长度也是1
|
|
|
|
|
|
for (int i = 0; i < cnt; i++) { // 枚举区间内所有质数
|
|
|
int p = primes[i];
|
|
|
/*
|
|
|
C(a,b)=a!/(b! * (a-b)!)
|
|
|
a!中有多少个质数因子p
|
|
|
减去(a-b)!的多少个质数因子p,
|
|
|
再减去b!的质数因子p的个数，就是总个数
|
|
|
s记录了p这个质数因子出现的次数
|
|
|
*/
|
|
|
int s = get(a, p) - get(b, p) - get(a - b, p);
|
|
|
while (s--) mul(c, cl, p); // 不断的乘p,结果保存到数组c中。len将带回c的有效长度
|
|
|
}
|
|
|
}
|
|
|
|
|
|
int a, b;
|
|
|
int c[N], cl;
|
|
|
|
|
|
int main() {
|
|
|
cin >> a >> b;
|
|
|
// 筛质数
|
|
|
get_primes(N - 1);
|
|
|
|
|
|
// 计算组合数
|
|
|
C(a, b, c, cl);
|
|
|
|
|
|
// 输出高精度结果
|
|
|
for (int i = cl; i >= 1; i--) printf("%d", c[i]);
|
|
|
return 0;
|
|
|
}
|
|
|
```
|
|
|
#### 7、这段代码如何理解？
|
|
|
```
|
|
|
s = get(a, p) - get(a - b, p) - get(b, p);
|
|
|
```
|
|
|
|
|
|
看公式一：
|
|
|
$C_a^b=\frac{a \times (a-1) \times ... \times(a-b+1)}{b\times (b-1) \times ... \times 1}$ $=\frac{a \times (a-1) \times ... \times (a-b+1) \times (a-b)!}{ b \times (b-1) \times ... \times 1 \times (a-b)!}$ $=\frac{a!}{b! \times (a-b)!}$
|
|
|
$a!$中质数$p$的个数，减去$b!$中质数$p$的个数，再减去$(a-b)!$中质数$p$的个数，就是公因子消掉的意思。 |
