BZOJ 1951 [Sdoi2010]古代猪文 欧拉定理+(扩展)lucas定理+逆元+快速幂

Description

“在那山的那边海的那边有一群小肥猪。他们活泼又聪明,他们调皮又灵敏。他们自由自在生活在那绿色的大草坪,他们善良勇敢相互都关心……” ——选自猪王国民歌 很久很久以前,在山的那边海的那边的某片风水宝地曾经存在过一个猪王国。猪王国地理位置偏僻,实施的是适应当时社会的自给自足的庄园经济,很少与外界联系,商贸活动就更少了。因此也很少有其他动物知道这样一个王国。 猪王国虽然不大,但是土地肥沃,屋舍俨然。如果一定要拿什么与之相比的话,那就只能是东晋陶渊明笔下的大家想象中的桃花源了。猪王勤政爱民,猪民安居乐业,邻里和睦相处,国家秩序井然,经济欣欣向荣,社会和谐稳定。和谐的社会带给猪民们对工作火红的热情和对未来的粉色的憧憬。 小猪iPig是猪王国的一个很普通的公民。小猪今年10岁了,在大肥猪学校上小学三年级。和大多数猪一样,他不是很聪明,因此经常遇到很多或者稀奇古怪或者旁人看来轻而易举的事情令他大伤脑筋。小猪后来参加了全猪信息学奥林匹克竞赛(Pig Olympiad in Informatics, POI),取得了不错的名次,最终保送进入了猪王国大学(Pig Kingdom University, PKU)深造。 现在的小猪已经能用计算机解决简单的问题了,比如能用P++语言编写程序计算出A + B的值。这个“成就”已经成为了他津津乐道的话题。当然,不明真相的同学们也开始对他刮目相看啦~ 小猪的故事就将从此展开,伴随大家两天时间,希望大家能够喜欢小猪。 题目描述 猪王国的文明源远流长,博大精深。 iPig在大肥猪学校图书馆中查阅资料,得知远古时期猪文文字总个数为N。当然,一种语言如果字数很多,字典也相应会很大。当时的猪王国国王考虑到如果修一本字典,规模有可能远远超过康熙字典,花费的猪力、物力将难以估量。故考虑再三没有进行这一项劳猪伤财之举。当然,猪王国的文字后来随着历史变迁逐渐进行了简化,去掉了一些不常用的字。 iPig打算研究古时某个朝代的猪文文字。根据相关文献记载,那个朝代流传的猪文文字恰好为远古时期的k分之一,其中k是N的一个正约数(可以是1和N)。不过具体是哪k分之一,以及k是多少,由于历史过于久远,已经无从考证了。 iPig觉得只要符合文献,每一种能整除N的k都是有可能的。他打算考虑到所有可能的k。显然当k等于某个定值时,该朝的猪文文字个数为N / k。然而从N个文字中保留下N / k个的情况也是相当多的。iPig预计,如果所有可能的k的所有情况数加起来为P的话,那么他研究古代文字的代价将会是G的P次方。 现在他想知道猪王国研究古代文字的代价是多少。由于iPig觉得这个数字可能是天文数字,所以你只需要告诉他答案除以999911659的余数就可以了。

Input

有且仅有一行:两个数N、G,用一个空格分开。

Output

有且仅有一行:一个数,表示答案除以999911659的余数。

Sample Input

4 2

Sample Output

2048

HINT

10%的数据中,1 <= N <= 50;
20%的数据中,1 <= N <= 1000;
40%的数据中,1 <= N <= 100000;
100%的数据中,1 <= G <= 1000000000,1 <= N <= 1000000000。

















数论大集合。。
题意就是对于n,以及它的所有约数d,求
G^sigma(C(n,n/d))  %999911659
经测试,999911659是个质数。。
那么由于指数不能直接对999911659取模。。所以得用一下欧拉定理,
也就是sigma(C(n,n/d))对phi(999911659)取模,
phi(999911659)=999911658。。
然后999911658当然不是个质数。。于是求组合数的时候,得用到扩展lucas定理。
扩展lucas定理的大致内容,,
就是把模数分解质因数后得到p1^k1*p2^k2*……*ps^ks,
那么分别求出组合数%p1^k1,%p2^k2……因为p1……ps是互质的,
直接用中国剩余定理合并就好了。
还有一个问题就是求p1^k1之类的,由于是质数的x次幂,可以考虑从阶乘入手,
把阶乘中的p1因子都晒去,然后可以发现是以p1^k1为循环节的,
而且不符合循环节部分也只有p1^k1长度。。暴力就好了。
具体的在代码里面(代码好乱啊->A->)。

吐槽一下。。由于999911658=2*3*4679*35617。。
所以ki都是1,分别求的时候分别用lucas定理就好了= =
一开始硬生生地用扩展lucas,结果一直TLE,一个点要1.5s……
换成lucas连0.5s都不到了= c =





直接分别lucas:
#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int 
	mod=999911659;
//999911658 = 2*3*4679*35617
int fac[4][35618],inv[4][35618];
void Pre(int num,int mod){
	inv[num][0]=inv[num][1]=fac[num][0]=1;
	for (int i=1;i<mod;i++) fac[num][i]=fac[num][i-1]*i%mod;
	for (int i=2;i<mod;i++)
		inv[num][i]=(ll)(mod-mod/i)*inv[num][mod%i]%mod;
	for (int i=2;i<mod;i++) inv[num][i]=(ll)inv[num][i-1]*inv[num][i]%mod;
}
int ksm(int x,int y,int mod){
	int z=1;
	while (y){
		if (y&1) z=(ll)z*x%mod;
		y>>=1;
		x=(ll)x*x%mod;
	}
	return z;
}
int exgcd(int a,int b,int &x,int &y){
	if (!b){
		x=1,y=0;
		return a;
	}
	int tt=exgcd(b,a%b,x,y),tmp=x;
	x=y,y=tmp-a/b*y;
	return tt;
}
int INV(int a,int b){
	int t1,t2;
	exgcd(a,b,t1,t2);
	t1=(t1+b)%b;
	return t1;
}
int Lucas(int num,int n,int m,int mod){
	if (n<m) return 0;
	if (n<mod && m<mod)
		return (ll)fac[num][n]*inv[num][m]%mod*inv[num][n-m]%mod;
	return (ll)Lucas(num,n/mod,m/mod,mod)*Lucas(num,n%mod,m%mod,mod)%mod;
}
int C(int nn,int mm){
	int t=0,t1=1,x=mod-1;
	for (int n=nn,m=mm;m;n/=2,m/=2)
		t1=(ll)t1*Lucas(0,n%2,m%2,2)*(x/2)%x*INV(x/2,2)%x;
	t=((ll)t1+t)%x,t1=1;
	for (int n=nn,m=mm;m;n/=3,m/=3)
		t1=(ll)t1*Lucas(1,n%3,m%3,3)*(x/3)%x*INV(x/3,3)%x;
	t=((ll)t1+t)%x,t1=1;
	for (int n=nn,m=mm;m;n/=4679,m/=4679)
		t1=(ll)t1*Lucas(2,n%4679,m%4679,4679)*(x/4679)%x*INV(x/4679,4679)%x;
	t=((ll)t1+t)%x,t1=1;
	for (int n=nn,m=mm;m;n/=35617,m/=35617)
		t1=(ll)t1*Lucas(3,n%35617,m%35617,35617)*(x/35617)%x*INV(x/35617,35617)%x;
	return ((ll)t1+t)%x;
}
int main(){
	Pre(0,2),Pre(1,3),Pre(2,4679),Pre(3,35617);
	int n,G;
	scanf("%d%d",&n,&G);
	if (G==mod) return puts("0"),0;
	int t=sqrt(n),ans=0;
	for (int i=1;i<=t;i++)
		if (n%i==0){
			ans=((ll)ans+C(n,i))%(mod-1);
			if (n/i!=i) ans=((ll)ans+C(n,n/i))%(mod-1);
		}
	printf("%d\n",ksm(G,ans,mod));
	return 0;
}


扩展lucas做法:
#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int 
    mod=999911659;
//999911658 = 2*3*4679*35617
int ksm(int x,int y,int mod){
    int z=1;
    while (y){
        if (y&1) z=(ll)z*x%mod;
        y>>=1;
        x=(ll)x*x%mod;
    }
    return z;
}
int inv(int a,int b){return ksm(a,b-2,b);}
int Fac(int n,int pi,int pk){
    if (!n) return 1;
    int t=1;
    if (n/pk){
        for (int i=2;i<=pk;i++)
            if (i%pi) t=(ll)t*i%pk;
        t=ksm(t,n/pk,pk);
    }
    for (int i=2;i<=n%pk;i++)
        if (i%pi) t=(ll)t*i%pk;
    return (ll)t*Fac(n/pi,pi,pk)%pk;
}
int Lucas(int n,int m,int pi,int pk){
    if (n<m) return 0;
    int x=Fac(n,pi,pk),y=Fac(m,pi,pk),z=Fac(n-m,pi,pk);
    int p1=0;
    for (int p=n;p;p/=pi) p1+=p/pi;
    for (int p=m;p;p/=pi) p1-=p/pi;
    for (int p=n-m;p;p/=pi) p1-=p/pi;
    return (ll)x*ksm(pi,p1,pk)%pk*inv(y,pk)%pk*inv(z,pk)%pk;
}
int C(int n,int m){
    int t=0;
    t=(t+(ll)Lucas(n,m,2,2)*(mod-1)/2*inv((mod-1)/2,2))%(mod-1);
    t=(t+(ll)Lucas(n,m,3,3)*(mod-1)/3*inv((mod-1)/3,3))%(mod-1);
    t=(t+(ll)Lucas(n,m,4679,4679)*(mod-1)/4679*inv((mod-1)/4679,4679))%(mod-1);
    t=(t+(ll)Lucas(n,m,35617,35617)*(mod-1)/35617*inv((mod-1)/35617,35617))%(mod-1);
    return t;
}
int main(){
    int n,G;
    scanf("%d%d",&n,&G);
    if (G==mod) return puts("0"),0;
    int t=sqrt(n),ans=0;
    for (int i=1;i<=t;i++)
        if (n%i==0){
            ans=((ll)ans+C(n,i))%(mod-1);
            if (n/i!=i) ans=((ll)ans+C(n,n/i))%(mod-1);
        }
    printf("%d\n",ksm(G,ans,mod));
    return 0;
}


  • 0
    点赞
  • 2
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值