首页 > 代码库 > 【GDOI 2011 DAY2 T3】零什么的最讨厌了 （快速求阶乘、中国剩余定理）
【GDOI 2011 DAY2 T3】零什么的最讨厌了 （快速求阶乘、中国剩余定理）
问题描述：
林记在做数学习题的时候，经常遇到这种情况：苦思冥想了很久终于把问题解出来，结果发现答案是0，久而久之林记在得到习题答案是0的时候就没有了做出一道难题的成就感。于是林记决定：以后出题，答案一定不能是0，例如求n！最低位非零数这样的习题就很不错了。
现在林记提出了一个更难一点的问题：求n！在K进制下的最低位非零数。其中K符合一些特殊的条件：K是由若干个互不相同的质数相乘得出来的，例如K=2,3,5,6,7,10……
输入格式：
首先输入的第一行是一个整数Q，表示询问的个数。
接下来是Q个询问，每个询问有两行组成，第一行首先是一个整数nk，然后紧跟着nk个正整数：m1，m2,……mnk，则K为所有mi的乘积，输入保证mi是有若干个互不相同的质数相乘得出来的。接下来一行给出一个整数n。
输出格式：
对于每个询问，如果K不符合题目限制，则输出“I hate zero.”（不用加双引号）。否则输出K进制下n！的最低位非零数。
输入样例：
输出样例：
3
1 17
9
1 6
49
1 93
16
15
2
78
数据范围：
对于10%的数据满足：n≤10，K≤100
对于50%的数据满足：n≤100000，K≤100
对于100%的数据满足：n≤1015，K≤1015，mi≤106，Q≤20
【分析】
快速求阶乘大法！！
下面某步还用到 威尔逊定理 （p-1）！%p=p-1...... （也可以预处理出来的）
把m分解质因数，对于质数p单独求n！=p^a+b的a的值和b%p的值（logn递归搞定），然后用中国剩余定理合并起来。
long long的乘积会爆，所以...要用那个不断乘2的方法：
LL mul(LL x,LL y,LL p){ if(x==0||y==0) return 0; LL now=x,yy=y,add=0; while(yy>1) { if(yy%2!=0) add=(add+now)%p; now=(now*2)%p; yy/=2; } return (now+add)%p;}
代码如下：
1 #include<cstdio> 2 #include<cstdlib> 3 #include<cstring> 4 #include<iostream> 5 #include<algorithm> 6 #include<queue> 7 #include<cmath> 8 using namespace std; 9 #define Maxn 1000010 10 #define INF 0xfffffff 11 #define LL long long 12 13 LL nk,m[Maxn],ml,K; 14 LL n; 15 LL prime[Maxn],pl; 16 bool q[Maxn]; 17 LL mn; 18 19 LL mymin(LL x,LL y) {return x<y?x:y;} 20 21 void get_pri(LL mx) 22 { 23 pl=0; 24 memset(q,1,sizeof(q)); 25 for(LL i=2;i<=mx;i++) 26 { 27 if(q[i]) prime[++pl]=i; 28 for(int j=1;j<=pl;j++) 29 { 30 if(i*prime[j]>mx) break; 31 q[i*prime[j]]=0; 32 if(i%prime[j]==0) break; 33 } 34 } 35 } 36 37 bool div(LL mm) 38 { 39 for(LL i=1;prime[i]*prime[i]<=mm;i++) if(mm%prime[i]==0) 40 { 41 m[++ml]=prime[i]; 42 mm/=prime[i]; 43 while(mm%prime[i]==0) return 0; 44 } 45 if(mm!=1) m[++ml]=mm; 46 return 1; 47 } 48 49 bool init() 50 { 51 ml=0; 52 scanf("%lld",&nk); 53 K=1; 54 bool ok=1; 55 for(int i=1;i<=nk;i++) 56 { 57 LL mm; 58 scanf("%lld",&mm);K*=mm; 59 if(!div(mm)) ok=0; 60 } 61 sort(m+1,m+1+ml); 62 for(LL i=2;i<=ml;i++) if(m[i]==m[i-1]) {ok=0;break;} 63 scanf("%lld",&n); 64 if(!ok) return 0; 65 return 1; 66 } 67 68 LL f[Maxn]; 69 70 struct node{LL x,y;}; 71 node c[Maxn]; 72 73 node ffind(LL x,LL p) 74 { 75 node ans; 76 if(x<p) 77 { 78 ans.x=0; 79 ans.y=f[x]; 80 return ans; 81 } 82 ans.x=0;ans.y=1; 83 node t=ffind(x/p,p); 84 ans.x+=t.x; ans.y*=t.y; 85 ans.y=(ans.y*f[x%p])%p; 86 ans.x+=x/p; 87 88 LL y=(((x/p)%2)==0)?1:p-1; 89 ans.y=(ans.y*y)%p; 90 91 return ans; 92 } 93 94 LL mul(LL x,LL y,LL p) 95 { 96 if(x==0||y==0) return 0; 97 LL now=x,yy=y,add=0; 98 while(yy>1) 99 {100 if(yy%2!=0) add=(add+now)%p;101 now=(now*2)%p;102 yy/=2;103 }104 return (now+add)%p;105 }106 107 LL qpow(LL x,LL b,LL p)108 {109 LL ans=1;110 while(b)111 {112 if(b&1) ans=mul(ans,x,p);113 x=mul(x,x,p);114 b>>=1;115 }116 return ans;117 }118 119 LL d[Maxn];120 121 LL get_ans()122 {123 mn=c[1].x;124 for(LL i=2;i<=ml;i++) mn=mymin(mn,c[i].x);125 for(LL i=1;i<=ml;i++)126 {127 if(c[i].x>mn) {d[i]=0;continue;}128 LL bm=qpow(K/m[i],mn,m[i]);129 bm=qpow(bm,m[i]-2,m[i]);130 d[i]=mul(bm,c[i].y,m[i]);//(bm*c[i].y)%m[i];131 }132 LL ans=0;133 for(LL i=1;i<=ml;i++)134 {135 LL y=qpow(K/m[i],m[i]-2,m[i]);136 //ans=(ans+((K/m[i]*y)%K)*d[i])%K;137 ans=(ans+mul(mul(K/m[i],y,K),d[i],K))%K;138 }139 return ans;140 }141 142 int main()143 {144 get_pri(1000000);145 int T;146 scanf("%d",&T);147 while(T--)148 {149 if(!init()) {printf("I hate zero.\n");continue;}150 for(LL i=1;i<=ml;i++)151 {152 f[0]=1;153 for(LL j=1;j<m[i];j++) f[j]=(f[j-1]*j)%m[i];154 c[i]=ffind(n,(LL)m[i]);155 }156 printf("%lld\n",get_ans());157 }158 return 0;159 }
2016-09-04 19:37:26
【GDOI 2011 DAY2 T3】零什么的最讨厌了 （快速求阶乘、中国剩余定理）