Description
猪王国的文明源远流长,博大精深。
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的余数。
Range
10%的数据中,$1 <= N <= 50$;
20%的数据中,$1 <= N <= 1000$;
40%的数据中,$1 <= N <= 100000$;
100%的数据中,$1 <= G <= 1000000000,1 <= N <= 1000000000$。
Solution
先一句话题意:求 $G^{\Sigma _{d \mid n} C_n^d} \; \% 999911659$ 的值。
若 G =999911659,则上式结果为 0。否则,因为999911659是质数,所以 G,n 互质。由欧拉定理推论得:
$$G^{\Sigma _{d \mid n} C_n^d} \equiv G^{\Sigma _{d \mid n} C_n^d \; \% \; 999911658} \quad (mod\;999911659)$$
尝试分解质因数,可以发现 $999911658 = 2 \times 3 \times 4679 \times 35617$。因为所以质因子的指数都为 1,所以我们枚举 $n$ 的约数 $d$,然后运用 $Lucas$ 求出组合数 $C_n^d$ ,分别计算出 $\Sigma _{d \mid n} C_n^d$ 对 2,3,4679,35617 四个质数取模的结果。然后中国剩余定理合并即可。
快速幂非递归版本 1A,递归版本无限 wa。玄学错误=.=
Code
1 #include<cstdio> 2 #define mod 999911659 3 #define int long long 4 5 int n,g; 6 int m[6],r[6]; 7 int factorcnt; 8 int inv[100005]; 9 int fac[100005]; 10 int factor[10005]; 11 12 void init_factor(){ 13 for(int i=1;i*i<=n;i++){ 14 if(n%i) continue; 15 factor[++factorcnt]=i; 16 if(i*i!=n) factor[++factorcnt]=n/i; 17 } 18 } 19 20 int lucas(int x,int y,int p){ 21 if(x<y) return 0; 22 if(x>=p) return (lucas(x%p,y%p,p)*lucas(x/p,y/p,p))%p; 23 return (fac[x]*inv[fac[y]]*inv[fac[x-y]])%p; 24 } 25 26 int calc(int x){ 27 fac[0]=fac[1]=1; 28 inv[0]=inv[1]=1; 29 for(int i=2;i<=m[x];i++) fac[i]=(fac[i-1]*i)%m[x]; 30 for(int i=2;i<=m[x];i++) inv[i]=((m[x]-m[x]/i)*inv[m[x]%i])%m[x]; 31 for(int i=1;i<=factorcnt;i++) 32 (r[x]+=lucas(n,factor[i],m[x]))%=m[x]; 33 } 34 35 int exgcd(int a,int b,int &x,int &y){ 36 if(!b){ 37 x=1; y=0; 38 return a; 39 } 40 int c=exgcd(b,a%b,x,y); 41 int t=x; 42 x=y; 43 y=t-a/b*y; 44 return c; 45 } 46 47 int innv(int a,int b){ 48 int x,y; 49 exgcd(a,b,x,y); 50 return (x%b+b)%b; 51 } 52 53 int CRT(){ 54 int M=1,ans=0; 55 for(int i=1;i<=4;i++) M*=m[i]; 56 for(int i=1;i<=4;i++) 57 (ans+=(M/m[i])*innv(M/m[i],m[i])*r[i])%=M; 58 return (ans%M+M)%M; 59 } 60 61 int ksm(int x,int y){ 62 /*if(y==0) return 1; 63 if(y==1) return x%mod; 64 int c=ksm(x,y>>1); 65 if(y&1) return (x*c*c)%mod; 66 return (c*c)%mod;*/ 67 int r=1; 68 for(;y;y>>=1){ 69 if(y&1) (r*=x)%=mod; 70 (x*=x)%=mod; 71 } 72 return r; 73 } 74 75 signed main(){ 76 scanf("%lld%lld",&n,&g); 77 if(g==mod){ 78 puts("0"); 79 return 0; 80 } 81 m[1]=2; m[2]=3; m[3]=4679; m[4]=35617; 82 init_factor(); 83 for(int i=1;i<=4;i++) calc(i); 84 printf("%lld\n",ksm(g,CRT())%mod); 85 return 0; 86 }
原文地址:https://www.cnblogs.com/YoungNeal/p/8955639.html