「SDOI2010」古代猪文 Lucas+中国剩余定理

题意:

猪王国的文明源远流长,博大精深。

iPig在大肥猪学校图书馆中查阅资料,得知远古时期猪文文字总个数为N。当然,一种语言如果字数很多,字典也相应会很大。当时的猪王国国王考虑到如果修一本字典,规模有可能远远超过康熙字典,花费的猪力、物力将难以估量。故考虑再三没有进行这一项劳猪伤财之举。当然,猪王国的文字后来随着历史变迁逐渐进行了简化,去掉了一些不常用的字。

iPig打算研究古时某个朝代的猪文文字。根据相关文献记载,那个朝代流传的猪文文字恰好为远古时期的k分之一,其中k是N的一个正约数(可以是1和N)。不过具体是哪k分之一,以及k是多少,由于历史过于久远,已经无从考证了。

iPig觉得只要符合文献,每一种能整除N的k都是有可能的。他打算考虑到所有可能的k。显然当k等于某个定值时,该朝的猪文文字个数为(frac{N}{k}) 。然而从N个文字中保留下(frac{N}{k})
个的情况也是相当多的。iPig预计,如果所有可能的k的所有情况数加起来为P的话,那么他研究古代文字的代价将会是G的P次方。

现在他想知道猪王国研究古代文字的代价是多少。由于iPig觉得这个数字可能是天文数字,所以你只需要告诉他答案除以999911659的余数就可以了。

输入:

输入文件ancient.in有且仅有一行:两个数N,G,用一个空格分开。

输出:

输出文件ancient.out有且仅有一行:一个数,表示答案除以999911659的余数。

思路:

(表示考试的时候是完全不会)
1、首先(P)会很大,所以需要有一个模数(x)满足 (G^{P mod x} equiv G^P(mod\,Mod))

由费马小定理 (G^{Mod-1}equiv1(mod\,Mod))可知 (x=Mod-1)
因为可以将(P)拆成(frac{P}{x})(x)相乘再乘上(P%x),前面的数在对于(Mod)取模的情况下都为(1),不会对答案产生影响

2、题目中要求计算的是(G^{sum_{x|n}C_{n}^{x}}) 因此我们需要处理出来组合数。

但是题目中的(n)的范围最大到(1000000000),用递推和逆元都无法直接运算,不过可以发现(Mod-1)可以拆成四个质数的乘积,即(999911658=2*3*4679*35617)这里最大的也只有(35617),此时就可以选择采用(Lucas) 现学的不会证明

也就是这个样子:

int lucas(int x,int y,int p) {
	if(y==0)return 1;
	return 1ll*calc(x%a[p],y%a[p],p)*lucas(x/a[p],y/a[p],p)%a[p];
}

(上面的p是指这四个模数)

3、对于每个模数我们都求出了一个组合数,那么问题就转换成了找到一个数,使得

对于(i∈[1,4])都满足,(xequiv b_i(mod \,a_i))

这条式子也就可以用中国剩余定理求解√

上述提到的几个知识点应该挺重要的,虽然我考试的时候一个不会一个不熟

注意点:如果 (G=Mod)答案应为(0),但根据Lucas求解会得出(1),需要特别判一下

唔如果卡时间的话,可以选择线性推逆元。

代码如下:

#include<bits/stdc++.h>
#define Mod 999911659
using namespace std;
int n,G,a[]= {2,3,4679,35617},b[4],pr[4][35620],ni[4][35620];
int mul(int x,int y,int mod) {
	int res=1;
	while(y) {
		if(y&1)res=1ll*res*x%mod;
		x=1ll*x*x%mod,y>>=1;
	}
	return res;
}
int calc(int x,int y,int p) {
	return 1ll*pr[p][x]*ni[p][y]%a[p]*ni[p][x-y]%a[p];
}
int lucas(int x,int y,int p) {
	if(y==0)return 1;
	return 1ll*calc(x%a[p],y%a[p],p)*lucas(x/a[p],y/a[p],p)%a[p];
}
int ans=0;
void exgcd(int A,int B,int &x,int &y) {
	if(B==0)x=1,y=0;
	else exgcd(B,A%B,y,x),y-=(A/B)*x;
}
void Calc(int x,int y) {
	int M=1;
	for(int i=0; i<4; i++)b[i]=lucas(x,y,i),M*=a[i];
	int tot=0;
	for(int i=0; i<4; i++) {
		int A=0,B=0,p=M/a[i];
		exgcd(p,a[i],A,B);//(A*p)%a[i]=1
		A=(A%a[i]+a[i])%a[i];
		tot=(tot+1ll*A*b[i]%M*p%M)%M;
	}
	ans=(ans+tot)%(Mod-1);
}
int main() {
	scanf("%d%d",&n,&G);
	int k=sqrt(n);
	for(int i=0; i<4; i++) {
		pr[i][0]=1,ni[i][0]=1;
		for(int j=1; j<=a[i]; j++)pr[i][j]=1ll*pr[i][j-1]*j%a[i],ni[i][j]=mul(pr[i][j],a[i]-2,a[i]);
	}
	if(Mod==G)printf("0
");
	else {
		for(int i=1; i<=k; i++) {
			if(n%i)continue;
			int a=i,b=n/i;//c(n,a)+c(n,b)
			Calc(n,a);
			if(a!=b)Calc(n,b);
		}
		printf("%d
",mul(G,ans,Mod));
	}
	return 0;
}
原文地址:https://www.cnblogs.com/cly1231/p/11355743.html