HDU3398 String

11 May 2019 | duanyll | Tags: OI 题解 数论

有毒。。。谁想得到要这样做。。。

题意

一个字符串只由0，1组成，且0有m个，1有n个，要求该字符串中任意的前缀中1的个数不能小于0的个数，问这样的字符串一共有多少个。结果对20100501取膜。

分析

首先总方案数是$C^{n}_{n+m}$，考虑不合法方案数。

此处有一个很有意思的抽象（谁想得到啊），就是在平面直角坐标系中，从$(0,0)$开始，在字符串中遇到0就向右走1，遇到1就向上走1，最后要走到$(m,n)$。那么，合法的路径一定不会与直线$y<x$即$y=x-1$相交，不合法的路径就一定会与它相交。

但是这样还是没法求，因为起点和终点都在直线的同侧。我们可以找到$(0,0)$关于$y=x-1$的对称点$(1,-1)$，这样从$(1,-1)$出发到$(m,n)$的所有路线（取一个对称就是原来的路径）就会与直线相交了，可以直接用排列组合的公式（相当于问从$(0,0)$走到$(m-1,n+1)$的路线数）求得不合法路径数为$C^{m-1}_{n+m}$。

故最终答案是

\[C^{n}_{n+m} - C^{m-1}_{n+m}\]

膜质数意义下的组合数可用卢卡斯定理快速求出

\[C^{n}_{m}\%p=(C^{n/p}_{m/p}*C^{n\%p}_{m\%p})\%p\]

可是20100501不是质数。展开上式得

\[\frac{(n+1-m)(n+m)!}{m!(n+1)!}\]

有除法，不能用逆元。应该用唯一分解定理对上下约分再把素因子乘起来。

阶乘的质因数分解

对于n的阶乘中数x的次数，1到n中含有x的一次项的有$x,2x,3x,\cdots,\lfloor \frac{n}{x} \rfloor x$共计$\lfloor \frac{n}{x} \rfloor$个，以此类推，含有二次项的有$\lfloor \frac{n}{x^2} \rfloor$个，含有三次项的有$\lfloor \frac{n}{x^3} \rfloor$个……最终想的次数就是以上数列求和，可以递推解决。

int count_fact(int n, int x) {
    int ret = 0;
    while (n > 0) {
        n /= x;
        ret += n;
    }
    return ret;
}

代码

#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>
using namespace std;

typedef long long int64;

const int MAXN = 1e6+10;
const int MAXA = 2e6+10;
const int64 MOD = 20100501;
const int AMXA = 2e6;

int64 prime[MAXN];
bool isntp[MAXA];
int pcnt;

void make_prime() {
    pcnt = 0;
    memset(isntp, false, sizeof(isntp));
    for (int i = 2; i < AMXA; ++i) {
        if (!isntp[i]) {
            prime[++pcnt] = i;
        }
        for (int j = 1; j <= pcnt; ++j) {
            if (i * prime[j] > AMXA) {
                break;
            }
            isntp[i * prime[j]] = true;
            if (i % prime[j] == 0) {
                break;
            }
        }
    }
}

template <typename T>
T pow_mod(T a, T b, T MOD) {
    T res = 1;
    while (b) {
        if (b & 1) res = res * a % MOD;
        a = a * a % MOD;
        b /= 2;
    }
    return res;
}

int64 count_fact(int64 n, int64 x) {
    int64 ret = 0;
    while (n > 0) {
        n /= x;
        ret += n;
    }
    return ret;
}

int main(){
	make_prime();
	int tcnt;
	cin >> tcnt;
	for(int T = 1;T<=tcnt;T++){
		int n,m;
		cin >> n >> m;
		int64 ans = 1;
		int64 res = n-m+1;	//原式中n-m+1项 
		for(int i = 1;i<=pcnt && prime[i] <= (m+n);i++){
			int64 cnt = 0;
			while(res % prime[i] == 0){
				cnt++;
				res /= prime[i];
			}
			cnt += count_fact(n+m,prime[i]);
			cnt -= count_fact(m,prime[i]);
			cnt -= count_fact(n+1,prime[i]);
			ans *= pow_mod(prime[i],cnt,MOD);
			ans %= MOD;
		}
		cout << ans << endl;
	}
}