组合数
c(a,b)
根据题目给的所求组合数的数目 n ，和a，b的范围来选择正确的方法求组合数
方法一（杨辉三角）
- 数据范围
1≤n≤10000,
1≤b≤a≤2000。
mod = 1e9 + 7
由于发现a,b的范围非常小，是2000，，所以可以根据杨辉三角，预处理组合数，然后查表输出。题目是对1e9+7取模，由于是使用杨辉三角，都是加法，所以不会超出int范围。
- 时间复杂度 $O(2000 * 2000 + n)$ 所以不会tle
参考代码
#include <algorithm>
#include <iostream>
using namespace std;
const int N = 2100 ,mod = 1e9 + 7 ;
int c[N][N] ;
void init(){
for(int i = 0 ; i < N ; i++ ){
for(int j = 0 ; j <= i ; j++){
if(!j) c[i][j] = 1;
else c[i][j] = (c[i-1][j] + c[i-1][j-1]) % mod; //根据杨辉三角，记得取mod
}
}
}
int main(){
init() ; //预处理组合数
int t ;
cin >> t ;
while(t--){
int a,b;
cin >> a >> b ;
cout << c[a][b] << endl; //查表输出
}
return 0 ;
}
方法二(逆元)
- 数据范围
1≤n≤10000,
1≤b≤a≤105
mod = 1e9 + 7
我门发现a和b的范围扩大到了1e5，所以不能使用杨辉三角(会t，空间也不允许),，然而发现题目中要求模为一个素数，
那么就可以使用费马小定理快速求得逆元，，所以可以根据组合数定义，求逆元来解决。
费马小定理 如果p是一个质数，而整数a不是p的倍数，则有$a^{(p-1)}$ ≡ 1 (mod p)
那么 转化一下 ， a * $a^{p-2}$ ≡ 1 (mod p)
那么$a^{p-2}$ 就是a(mod p)的逆元
组合数 c(a,b) = (a!)/(b! * (a - b)!)
可以先预处理出所有阶乘，和阶乘的逆元
- 时间复杂度 ：预处理阶乘逆元是 5000 * log(p-2) ,所以最后时间复杂度为$O(n + 5000 * log(p - 2)$
参考代码
#include <algorithm>
#include <cstring>
#include <cstdio>
#include <iostream>
using namespace std;
typedef long long ll ;
const int N = 1e5 + 100 , mod = 1e9 + 7 ;
int n ;
int a,b;
int fact[N],infact[N] ; //fact存阶乘 ，infact存阶乘的逆元
int qmi(int a,int n,int mod){ //快速幂 求得 a^n % mod
int res = 1 ;
while(n){
if(n & 1) res = (ll)res * a % mod ;
a = (ll)a * a % mod ;
n >>= 1 ;
}
return res ;
}
int main(){
fact[0] = infact[0] = 1;
for(int i = 1 ; i < N ; i++){
fact[i] = (ll)fact[i-1] * i % mod ; //预处理前5000的阶乘和逆元
infact[i] = (ll)infact[i-1] * qmi(i,mod - 2,mod) % mod ;
}
scanf("%d",&n) ;
while(n--){
scanf("%d%d",&a,&b) ;
printf("%d\n",(ll)fact[a] * infact[b] % mod * infact[a - b] % mod ) ; //组合数定义
}
return 0;
}
方法三(卢卡斯定理)
- 数据范围
1≤n≤20,
1≤b≤a≤1018,
1≤p≤105,且p均为素数
发现，题目给的a，b很大，肯定不能预处理，所以只能直接求，直接按照定义是线性的1e18也是不行的，发现结果需要%p，p很小，所以可以使用卢卡斯定理，将求的组合数转化为（a，b）小于p组合数。
lucas定理， 递归表示为 lucas(a,b) = lucas(a/p,b/p) * C(a%p,b%p) ,如果a<p && b < p 那么直接返回 C(a,b) ;
证明自行跳转到 y总提高课视频
- 时间复杂度：lucas定理最多递归log(a)次，然后求一次(a和b都会小于p)组合数使用逆元来求是 1e5 * log(1e5)
,最后再乘n ，$O(n * 1e5 * log(le5) * log(a))$ ,其实当p很大时，lucas的log底数很大，所以不会超时
参考代码
#include <algorithm>
#include <cstring>
#include <cstdio>
#include <iostream>
using namespace std;
typedef long long ll ;
const int N = 1e5 + 100 ;
ll n ;
ll a,b,p ;
ll qmi(ll a ,ll n){ //快速幂
ll res = 1;
while(n){
if(n & 1) res = res * a % p ;
a = a * a % p ;
n >>= 1 ;
}
return res ;
}
ll C(ll a,ll b){ //使用组合数定义，结合逆元（因为p是素数）求得组合数
ll res = 1;
for(int i = a ,j = 1; j <= b ; i-- , j++){
res = res * i % p ;
res = res * qmi(j,p-2) % p ;
}
return res;
}
ll lucas(ll a,ll b){ // lucas定理递归处理
if(a < p && b < p) return C(a,b) ;
return lucas(a/p,b/p) * C(a%p,b%p) % p ;
}
int main(){
cin >> n ;
while(n--){
cin >> a >> b >> p ;
cout << lucas(a,b) << endl ;
}
return 0;
}
附 ：如果给的p不是素数，可以将组合数分解质因数，然后求得组合数（分解方法可以参考下面的高精度）
方法四(分解质因数)
- 数据范围
1≤b≤a≤5000 (无取模，需要高精度)
没什么好说的，根据组合数的定义，分解出所有的质因数
参考代码
#include <algorithm>
#include <cstring>
#include <cstdio>
#include <iostream>
using namespace std;
const int N = 5100 ;
int a,b;
int r[N] ;
int primes[N],cnt;
bool st[N] ;
void init(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] = 1 ;
if(i % primes[j] == 0) break ;
}
}
}
int get(int a,int p){ //获得从1到a中p因子的次方数，就是最多能能分解出多少p
int cnt = 0 ;
for(int i = a ; i ; i /= p) cnt += i / p ;
return cnt ;
}
void mul(int r[],int &len,int x){ //高精度数组乘以一个数(非高精度),同时更新高精度数组长度
int t = 0 ;
for(int i = 0 ; i < len ; i ++){
t += r[i] * x ;
r[i] = t % 10 ;
t /= 10 ;
}
while(t){
r[len++] = t % 10 ;
t /= 10 ;
}
}
int C(int a,int b,int r[]){
int len = 1;
r[0] = 1 ;
for(int i = 0 ; i < cnt; i++){ //枚举所有的素数（小于a）
int p = primes[i] ;
int s = get(a,p) - get(b,p) - get(a - b ,p) ; //根据定义C（a,b） = (a!) / (b ! (a - b )! )
while(s--){ //如果不是高精度可以使用快速幂
mul(r,len,p) ; //高精度就一个一个乘
}
}
return len ;
}
int main(){
init(N-1) ;
cin >> a >> b ;
int len = C(a,b,r) ;
int k = len - 1 ;
while(!r[k] && k > 0) k-- ; //输出高精度的结果
while(k >= 0) cout << r[k--] ;
cout << endl ;
return 0;
}
参考资料 y总基础课
