组合数的计算（+线性筛，快速幂的运用）

题目链接https://codeforces.com/contest/1475/problem/E

题意:

在n个博主中选择k个博主使得k个博主的粉丝和最大, 输出最大方案数.

思路:

选择k个最大粉丝数的博主就可以了. 直接$sort$.
但是为什么会出现多种选法的情况呢. 那就取决于选出k个博主中,粉丝数量最少的那个博主了. 如果选出的最低粉丝数是$f$, 那么如果拥有$f$粉丝数的博主就一个, 么情况必然只有一种; 但是, 如果有多个呢, 情况就会出现多个.

我们对会出现多个$f$个粉丝的情况进行讨论. 假设有$f$个粉丝的博主总数为$B$, 被选中的博主中,有$A$个博主有$f$个粉丝的话就可能出现以下情况:
粉丝数列 F a n s [ n ] = { z , y , x , . . . h , g , f . . . f , f , f , f , f } Fans[n] = \{z, y, x, ...h, g, f... f, f, f, f, f\} Fans[n]={z,y,x,...h,g,f...f,f,f,f,f}其中, $f$有B个.

当被选中的$f$粉丝数的博主数量少于等于$f$粉丝数的博主总数时(即$A<=B$), 就是从B个人中随意选择A个人,就是组合数$C\left(\genfrac{}{}{0ex}{}{A}{B}\right)$. 而且$A<=B$是一种必然情况.

所以在实现时, 只需要先算出在没有加上粉丝数为$f$的博主时的人数, 用选人的总数$k$中需要的粉丝数为$f$的博主时的人数$A$, 再计算出粉丝数为$f$的博主时的人数的总数$B$, 调用组合数模板就过了.

在交题时我用的网上找到的组合数模板, 模板如下:

ll fac[maxn]={1,1},inv[maxn]={1,1},f[maxn]={1,1};
ll cmn(ll a,ll b){
    if(b>a)return 0;
    return fac[a]*inv[b]%M*inv[a-b]%M;
}
void init(){//快速计算阶乘的逆元
    for(int i=2;i<maxn;i++){
        fac[i]=fac[i-1]*i%M;
        f[i]=(M-M/i)*f[M%i]%M;
        inv[i]=inv[i-1]*f[i]%M;
    }
}

此模板运用组合数公式 $C\left(\genfrac{}{}{0ex}{}{n}{m}\right)$ = $\frac{n!}{m!(n-m)!}$.计算出范内所有数的阶乘,然后运用此公式,调用cnm函数计算出答案即可.

在之后的补题中, 我又发现一种有趣的算法求组合数.

对于任何一个数, 都可以看成若干素数的乘积. 在范围内, 素数的值, 个数是确定的, 那么就需要求出在阶乘内, 出现的所有素数和所有素数的总数, 将式子简化为素数阶数相加减之后再相乘，就能快速且数据不溢出的条件下求出组合数答案.

比如: $17!={\prod }_{i=1}^{17}i=2^{15}\ast 3^6\ast 5^3\ast 7^2\ast 11^1\ast 13^1\ast 17^1$。

这个式子怎么计算的呢？拿2举例，对于整数$[1,17]$来说，含有因子2的数有 $17/2=8个$，被舍弃的$17-8=9个$数中，是不能被2整除的，所以不含因子2.那么对于这8个数字，还有$8/2=4个$数字还含有因子2，以此类推直到这个范围的整数的因子2的个数全部被找出为止，代码就为：

while(num){
    temp += num / 2; //temp记录素数在范围内出现的词次数。
    num /= 2;
}

由于范围内所有整数是乘法的形式出现的，那么素数因子2出现的次数在式子中就转换为指数形式即： $\ast 2^{temp}$.其他素数因子出现的次数以此类推。这个素数因子的范围呢当然是$[2,n]$啦。

这个就有一点点类似整数除法的思想. 对于任意一个因子$x$, 在$[1,y]$中, 含有这个因子的数一共有$y/x$个, 这是整数除法的部分定义,

将组合数的式子成功转换之后，由于式子结果需要取模，如果式子的运算结果可能太大以至于直接溢出，所以这种情况不能直接调用pow函数，需要运用快速幂，具体的思想和实现见这篇博文.

代码：

#include <bits/stdc++.h>

using namespace std;

typedef long long ll;
const int maxn = 1010;
const ll M = 1e9 + 7;
int b[maxn];
int mPow[maxn], nPow[maxn], n_mPow[maxn];

ll ksm(int x, int n){
    ll res = 1;
    while(n) {
        if(n & 1) res = res * x % M;
        x = x * x % M;
        n >>= 1;
    }
    return res;
}

ll prime[maxn];
ll visit[maxn];
void isPrime(){ //线性筛，O(n)
    memset(visit, 0, maxn);
    memset(prime, 0, maxn);
    for(int i = 2; i <= maxn; i ++ ){
        if(!visit[i]) prime[++ prime[0]] = i;
        for(int j = 1; j <= prime[0] && i * prime[j] <= maxn; j ++ ){
            visit[i * prime[j]] = 1;
            if(i % prime[j] == 0) break;
        }
    }
}

void powNum(int n, int POW[]){
    for(int i = 1; prime[i] <= n; i ++ ){
        int temp = 0, now = n;
        while(now){
            temp += now / prime[i];
            now /= prime[i];
        }
        POW[i] = temp;
    }
}

ll cnm(int n, int m) {
    memset(nPow, 0, maxn);
    memset(mPow, 0, maxn);
    memset(n_mPow, 0, maxn);
    powNum(n, nPow);
    powNum(m, mPow);
    powNum(n - m, n_mPow);
    ll res = 1;
    int MaxPrime = sqrt(n);
    for(int i = 1; prime[i] <= n; i ++ )
        res = res * ksm(prime[i], nPow[i] - mPow[i] - n_mPow[i]) % M;
    return res;
}

int main()
{
    ios::sync_with_stdio(0);
    cin.tie(0); cout.tie(0);

    #ifndef ONLINE_JUDGE
        freopen("test.txt", "r", stdin);
    #endif // ONLINE_JUDGE
    int t; cin >> t;
    isPrime();
    while(t --){
        int n, k; cin >> n >> k;
        for(int i = 0; i < n; i ++ ) cin >> b[i];
        sort(b, b + n);
        int f = b[n - k], B = 1;
        for(int i = n - k + 1; i < n; i ++ ){
            if(f != b[i]) break;
            else B ++;
        }
        int A = B;
        for(int i = n - k - 1; i >= 0; i -- ){
            if(f != b[i]) break;
            else B ++;
        }
        cout << cnm(B, A) << "\n";
    }
    return 0;
}