「BZOJ 1951」「SDOI2010」古代猪文

Description

猪王国的文明源远流长，博大精深。 iPig在大肥猪学校图书馆中查阅资料，得知远古时期猪文文字总个数为 (N)。当然，一种语言如果字数很多，字典也相应会很大。当时的猪王国国王考虑到如果修一本字典，规模有可能远远超过康熙字典，花费的猪力、物力将难以估量。故考虑再三没有进行这一项劳猪伤财之举。当然，猪王国的文字后来随着历史变迁逐渐进行了简化，去掉了一些不常用的字。 iPig打算研究古时某个朝代的猪文文字。根据相关文献记载，那个朝代流传的猪文文字恰好为远古时期的k分之一，其中 (k) 是 (N) 的一个正约数（可以是 (1) 和 (N)）。不过具体是哪 (k) 分之一，以及 (k) 是多少，由于历史过于久远，已经无从考证了。 iPig觉得只要符合文献，每一种能整除 (N) 的 (k) 都是有可能的。他打算考虑到所有可能的 (k) 。显然当 (k) 等于某个定值时，该朝的猪文文字个数为 (frac{N}{k}) 。然而从 (N) 个文字中保留下 (N / k) 个的情况也是相当多的。iPig预计，如果所有可能的 (k) 的所有情况数加起来为 (P) 的话，那么他研究古代文字的代价将会是 (G) 的 (P) 次方。 现在他想知道猪王国研究古代文字的代价是多少。由于iPig觉得这个数字可能是天文数字，所以你只需要告诉他答案除以 (999911659) 的余数就可以了。

Hint

( le N,G le 10^9)

Solution

简化题意：给定 (N,G)，求：

[large{ G^{ sum_{k|n} C^k_n } mod 999911659} ]

的值。

---

其实主要是求 (sum_{k|n} C^k_n) 的值。

我们发现 (999911659) 是质数，那么根据欧拉定理，可得 (G^{ sum_{k|n} C^k_n } equiv G^{ sum_{k|n} C^k_n mod varphi(999911659)}equiv G^{ sum_{k|n} C^k_n mod 999911658} pmod {999911659})

众所周知，利用 Lucas 定理可以在模数较小的情况下快速求组合数的值，而且要求模数为素数。

因此考虑将模数缩小。

首先将 (999911659) 分解质因数，得 (999911659 = 2 imes 3 imes 4679 imes 35617)，

然后用 Lucas 定理计算出 (a_1 = sum_{k|n} C^k_n mod 2,a_2 = sum_{k|n} C^k_n mod 3,a_3 = sum_{k|n} C^k_n mod 4679,a_4 = sum_{k|n} C^k_n mod 35617)

那么可以构造出方程组： (egin{cases} {x equiv a_1 pmod 2} \{x equiv a_2 pmod 3}\{x equiv a_3 pmod{4679}}\{x equiv a_4 pmod{235617}}end{cases})

发现了什么？直接可以使用 CRT 进行求解 (x) 就是 (sum_{k|n} C^k_n mod 999911658) 的值。

最后就直接快速幂。

Code

/*
 * Author : _Wallace_
 * Source : https://www.cnblogs.com/-Wallace-/
 * Problem : BZOJ 1951 SDOI2010 古代猪文
 */
#include<iostream>

using std::cin;
using std::cout;
using std::endl;
typedef long long LL;

LL fastPow(LL a,LL b,LL p) {
	if(!b) return 1ll;
	LL t=fastPow(a,b>>1,p);
	if(b&1) return t*t%p*a%p;
	else return t*t%p;
}
LL inv(const LL &b,const LL &p) {
	return fastPow(b,p-2,p);
}
LL CRT(int n,LL *a,LL *m) {
	LL M=1ll,ans=0ll;
	for(register int i=0;i<n;i++)
		M*=m[i];
	for(register int i=0;i<n;i++) {
		LL t=inv(M/m[i],m[i]);
		LL Mi=M/m[i];
		(ans+=a[i]*t%m[i]*Mi)%=M;
	}
	return ans;
}

LL fact[int(1e5+5)];
inline void initFact(int x) {
	fact[0]=fact[1]=1ll;
	for(register int i=2;i<=x;i++)
		fact[i]=(fact[i-1]*i)%x;
}
LL comb(LL n,LL m,LL mod) {
	if(n<m) return 0ll;
	return fact[n]*inv(fact[m],mod)%mod*inv(fact[n-m],mod)%mod;
}
LL Lucas(LL n,LL m,LL mod) {
	if(n<m) return 0ll;
	if(!n) return 1ll;
	return comb(n%mod,m%mod,mod)*Lucas(n/mod,m/mod,mod)%mod;
}

signed main() {
	LL n,G;
	cin>>n>>G;
	if(G%999911659==0) return cout<<0<<endl,0;
	LL m[4]={2,3,4679,35617};
	LL a[4]={0,0,0,0};
	for(register LL i=0;i<4;i++) {
		initFact(m[i]);
		for(register LL k=1;k*k<=n;k++) {
			if(n%k!=0) continue;
			(a[i]+=Lucas(n,k,m[i]))%=m[i];
			if(k*k!=n) (a[i]+=Lucas(n,n/k,m[i]))%=m[i];
		}	
	}
	LL P=CRT(4,a,m);
	cout<<fastPow(G,P,999911659)<<endl;
	return 0;
}