[组合数取模][中国剩余定理]

[BZOJ2142]礼物

这里引用Po姐的例子

19!%9=(1*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18*19) %9

=(1*2*4*5*7*8*10*11*13*14*16*17*19)*3^6*(1*2*3*4*5*6) %9

可以把和3有关的东西提出来,剩下的对于9就有逆元了(此时应该Exgcd)

最后再把3乘上去

也叫扩展Lucas

 

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
#include <cmath>
#define maxn 100
using namespace std;

typedef long long ll;
ll n, m, n1, n2, M, ans[maxn];

int a[maxn];
int p[maxn], p_c[maxn], Nm;
void FJ(ll x){
	int q = sqrt(x);
	for(int i = 2; i <= q; i ++){
		if(x % i == 0){
			++ Nm;
			p[Nm] = i;
			p_c[Nm] = 1;
			while(x % i == 0){
				x /= i;
				p_c[Nm] *= i;
			}
		}
	}
	if(x != 1){
		Nm ++;
		p[Nm] = p_c[Nm] = x;
	}
}

void Exgcd(ll a, ll b, ll& d, ll& x, ll& y){
	if(b == 0){d = a, x = 1, y = 0; return;}
	Exgcd(b, a % b, d, y, x); y -= x * (a / b);
}

ll Crt(){
	ll ret = 0, d, x, y;
	for(int i = 1; i <= Nm; i ++){
		ll w = M / p_c[i] * ans[i];
		Exgcd(M / p_c[i], p_c[i], d, x, y);
		ret += w * x, ret %= M;
	}
	return (ret + M) % M;
}

ll power_mod(ll a, ll b, ll md){
	ll ret = 1;
	while(b > 0){
		if(b & 1)ret = ret * a % md;
		b >>= 1;
		a = a * a % md;
	}return ret;
}

typedef pair<ll, int> pli;
ll TMP[maxn];
pli ask(ll n, int k){
	if(n == 0)return make_pair(1, 0);
	pli nw(0, n / p[k]);
	ll ans = 1;
	ans = power_mod(TMP[k], n / p_c[k], p_c[k]);
	for(ll i = n / p_c[k] * p_c[k] + 1; i <= n; i ++)
		if(i % p[k])ans = ans * i % p_c[k];
	pli t = ask(n / p[k], k);
	nw.first = ans * t.first % p_c[k];
	nw.second += t.second;
	return nw;
}

ll inv(ll a, ll b){
	ll d, x, y;
	Exgcd(a, b, d, x, y);
	return (x % b + b) % b;
}

ll get(ll n, ll m){
	if(n < m)return 0;
	pli nw, a, b, c;
	for(int i = 1; i <= Nm; i ++){
		a = ask(n, i), b = ask(m, i), c = ask(n - m, i);
		nw = make_pair(a.first * inv(b.first, p_c[i]) % p_c[i] * inv(c.first, p_c[i]) % p_c[i], a.second - b.second - c.second);
		ans[i] = power_mod(p[i], nw.second, p_c[i]) * nw.first % p_c[i];
	}
	return Crt();
}

ll ret;

void dfs(int dep, int sum, bool t){
	if(dep == n1 + 1){
		if(t)ret -= get(m - sum - 1, n - 1);
		else ret += get(m - sum - 1, n - 1);
		ret %= M;
		return;
	}
	dfs(dep + 1, sum, t);
	dfs(dep + 1, sum + a[dep], t ^ 1);
}

ll solve(){
	for(int i = n1 + 1; i <= n1 + n2; i ++)
		m -= a[i] - 1;
	ret = 0;
	dfs(1, 0, 0);
	return (ret + M) % M;
}

int main(){
	int test;
	scanf("%d%lld", &test, &M), FJ(M);
	
	for(int i = 1; i <= Nm; i ++){
		ll nw = 1;
		for(int j = 1; j < p_c[i]; j ++)
			if(j % p[i])nw = nw * j % p_c[i];
		TMP[i] = nw;
	}
	
	while(test --){
		scanf("%lld%lld%lld%lld", &n, &n1, &n2, &m);
		for(int i = 1; i <= n1 + n2; i ++)
			scanf("%d", &a[i]);
		printf("%lld\n", solve());
	}
	return 0;
}

  

 这里贴一下扩展CRT(和上述并没有什么关系)

首先来计算两组 x ≡ r1 ( mod a1 ) ; x ≡  r2 ( mod a2 ) ; 定义变量 k1 k2 得到 x = k1 * a1 + r1 ; x = k2 * a2 + r2 ; 由上面的等式得到 k1 * a1 + r1 = k2 * a2 + r2 ; 转化为 k1*a1 = (r2 - r1) + k2 *a2 ; 对左右取模a2,因为 (k2*a2)%s2 = 0 ,所以等式转化为  k1 * a1 ≡ ( r2 - r1 ) (mod a2) ;使用扩展欧几里得可以求解到 k1的值(判断是否存在k1的值),将k1带回到 x1 = k1 * a1 + r1 ;得到同时满足于{ x = k1 * a1 + r1 ; x = k2 * a2 + r2 ; }的一个特解 , 所以 x ≡ x1 (mod lcm(a1,a2) ) ; 也就是 x ≡ ( k1*a1+r1 ) ( mod ( a1*a2/d ) );这样也就将两个同余式转化为了一个,通过不断的转化,将k个等式合并为一个 ,用扩展欧几里得求出最小的正解x 

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
#define maxn 1010
using namespace std;
typedef long long ll;
ll Mm[maxn], Rr[maxn], n;

void Exgcd(ll a, ll b, ll &d, ll &x, ll &y){
	if(b == 0){d = a, x = 1, y = 0; return;}
	Exgcd(b, a % b, d, y, x); y -= x * (a / b);
}

ll Excrt(){
	ll m1 = Mm[0], r1 = Rr[0], m2, r2;
	ll c, d, x, y;
	for(ll i = 1; i < n; i ++){
		m2 = Mm[i], r2 = Rr[i];
		Exgcd(m1, m2, d, x, y);//d = gcd(m1, m2), x * m1 + y * m2 = d
		c = r2 - r1;
		if(c % d)return -1;
		ll t = m2 / d;//x * m1 + y * m2 = r2 - r1
		x = (c / d * x % t + t) % t;
		r1 = m1 * x + r1;
		m1 = m1 * m2 / d;
	}
	if(n == 1 && r1 == 0)return m1;
	return r1;
}

int main(){
	ll test, Case = 0;
	scanf("%lld", &test);
	while(test --){
		scanf("%lld", &n);
		for(ll i = 0; i < n; i ++)scanf("%lld", &Mm[i]);
		for(ll i = 0; i < n; i ++)scanf("%lld", &Rr[i]);
		printf("Case %lld: %lld\n", ++ Case, Excrt());
	}
	
	return 0;
}

hdu1573

#include <cstdio>
#include <iostream>
using namespace std;
#define maxn 100

typedef long long ll;
int n, N;

ll Mm[maxn], Rr[maxn];

void Exgcd(ll a, ll b, ll& d, ll& x, ll& y){
	if(b == 0){d = a, x = 1, y = 0; return;}
	Exgcd(b, a % b, d, y, x); y -= x * (a / b);
}

ll Excrt(){
	ll m1 = Mm[0], r1 = Rr[0], m2, r2;
	ll c, d, x, y;
	for(int i = 1; i < n; i ++){
		m2 = Mm[i], r2 = Rr[i];
		Exgcd(m1, m2, d, x, y);//x * m1 + y * m2 = d
		c = r2 - r1;
		if(c % d)return -1;
		ll t = m2 / d;
		x = (c / d * x % t + t) % t;
		r1 = m1 * x + r1;
		m1 = m1 * m2 / d;
	}
	if(r1 == 0)r1 = m1;
	if(r1 > N)return 0;	
	return (N - r1) / m1 + 1;
}

int main(){
	int test;
	scanf("%d", &test);
	while(test --){
		scanf("%d%d", &N, &n);
		for(int i = 0; i < n; i ++)scanf("%lld", &Mm[i]);
		for(int i = 0; i < n; i ++)scanf("%lld", &Rr[i]);
		ll nw = Excrt();
		if(~nw)printf("%lld\n", nw);
		else puts("0");
	}
	return 0;
}

  

posted @ 2016-05-30 15:52  _Horizon  阅读(503)  评论(0编辑  收藏  举报