考虑枚举x，分情况讨论
当Q=0时，无论x取什么都有对应的y合法，个数为gcd(x,P)个
否则对于x有合法的y当且仅当gcd(x,P)|Q，并且数量为gcd(x,P)，这个仔细思考还是不难理解的
ans=sigema(1~P-1)x [gcd(x,P)|Q]*gcd(x,P)
=sigema(d|gcd(P,Q))d d*(sigema(1~P/d)x [gcd(x,P/d)==1]) (如果Q等于0把d的限制去掉就完事了)
后面明显是欧拉函数，=sigema(d|gcd(P,Q))d d*phi(P/d)
然后是一个很类似狄利克雷卷积的形式，当积性函数好了
那么拆分出来，把每个质因子分开搞，设A为gcd(P,Q)在全局第k个质因子的出现次数，B为P在全局第k个质因子出现的次数
=product(1~plen)k [sigema(0~A)i pk^i*phi(pk^(B-i))]
=product(1~plen)k [sigema(0~A)i pk^i* (pk^(B-i)*(pk-1)/pk)]
=product(1~plen)k [sigema(0~A)i pk^B-1*(pk-1)]
=product(1~plen)k [sigema(0~A)i pk^B-1*(pk-1)]
=product(1~plen)k (A+1)*pk^B-1*(pk-1)
这样就可以快速算完了，但是有一点需要注意，当i==B时，phi(1)不能表示成1*(1-1)/1的形式，要特判一下A==B的情况
#include<cstdio> #include<iostream> #include<cstring> #include<cstdlib> #include<algorithm> #include<cmath> using namespace std; typedef long long LL; const LL mod=1e9+7; LL gcd(LL a,LL b) { if(a<0)a=-a; if(a==0)return b; return gcd(b%a,a); } LL quick_mul(LL A,LL p,LL mo) { LL ret=0; while(p!=0) { if(p%2==1)ret=(ret+A)%mo; A=(A+A)%mo;p/=2; } return ret; } LL quick_pow(LL A,LL p,LL mo) { LL ret=1; while(p!=0) { if(p%2==1)ret=quick_mul(ret,A,mo); A=quick_mul(A,A,mo);p/=2; } return ret; } //-----------------------------------def---------------------------------------------- bool check(LL a,LL t,LL mi,LL n) { LL d=quick_pow(a,t,n); if(d==1||d==n-1)return true; while(mi--) { d=quick_mul(d,d,n); if(d==n-1)return true; } return false; } int prr[10]={2,3,5,7,11,13,17,19,23,29}; bool miller_rabin(LL n) { if(n<=1)return false; for(int i=0;i<=9;i++) { if(n==prr[i])return true; if(n%prr[i]==0)return false; } LL t=n-1;int mi=0; while(t%2==0)t/=2,mi++; for(int i=1;i<=20;i++) { LL a=rand()%(n-1)+1; if(!check(a,t,mi,n))return false; } return true; } //~~~~~~~~~~~~~~~miller_rabin~~~~~~~~~~~~~~ LL pollard_rho(LL n) { LL a=rand()%n,b=a,c=rand()%n; LL i=1,k=2; while(1) { a=(quick_mul(a,a,n)+c)%n; LL d=gcd(b-a,n); if(d!=1&&d!=n)return d; if(a==b)return n; i++; if(i==k)b=a,k*=2; } } //~~~~~~~~~~~~~~~pollard_rho~~~~~~~~~~~~~~~ int pr;LL prime[110]; void getprime(LL n) { if(n==1)return ; if(miller_rabin(n)){prime[++pr]=n;return ;} LL d=n; while(d>=n)d=pollard_rho(d); getprime(d); getprime(n/d); } LL plen,p[1100],tt[1100];LL c[2][1100]; void merge() { int i=1,j=1,k=0; while(i<=pr&&j<=plen) { if(prime[i]<p[j])tt[++k]=prime[i++]; else if(prime[i]>p[j])tt[++k]=p[j++]; else tt[++k]=p[j],i++,j++; } while(i<=pr)tt[++k]=prime[i++]; while(j<=plen)tt[++k]=p[j++]; memcpy(p,tt,sizeof(p));plen=k; } int n;LL s[2][110]; void getci(int w) { for(int i=1;i<=n;i++) { LL u=s[w][i]; for(int j=1;j<=plen;j++) while(u%p[j]==0)u/=p[j],c[w][j]++; } if(w==1)for(int j=1;j<=plen;j++)c[w][j]=min(c[0][j],c[1][j]); } //------------------------------------------divi-------------------------------------------------- int main() { freopen("4.in","r",stdin); freopen("a.out","w",stdout); scanf("%d",&n); for(int i=1;i<=n;i++) { scanf("%lld",&s[0][i]); pr=0,getprime(s[0][i]); sort(prime+1,prime+pr+1); pr=unique(prime+1,prime+pr+1)-prime-1; merge(); } bool bk=false; for(int i=1;i<=n;i++) { scanf("%lld",&s[1][i]); if(s[1][i]==0){bk=true;break;} pr=0,getprime(s[1][i]); sort(prime+1,prime+pr+1); pr=unique(prime+1,prime+pr+1)-prime-1; merge(); } if(bk==true) { getci(0); memcpy(c[1],c[0],sizeof(c[1])); } else getci(0),getci(1); LL ans=1; for(int i=1;i<=plen;i++) { LL A=c[1][i],B=c[0][i]; if(A<B)ans=ans*(quick_pow(p[i],B-1,mod)*(p[i]-1)%mod*(A+1)%mod)%mod; else ans=ans*( (quick_pow(p[i],B-1,mod)*(p[i]-1)%mod*A%mod+quick_pow(p[i],A,mod))%mod )%mod; } printf("%lld ",ans); return 0; }