description
给你一个数组({a_i})以及(X,Y)，问你有多少对((i,j))满足存在一个(vin mathbb{N}^+)使得 (gcd(a_i,v)=X,mbox{lcm}(a_j,v)=Y)。(i,j)交换顺序被视为不同的数对，(i,j)可以相等。
(nle2 imes10^5,a_i,X,Yle10^{18})
sol
首先一定要有(X|Y)，不然肯定无解。
假设(Y=P_1^{mx_1}P_2^{mx_2}...P_k^{mx_k},X=P_1^{mn_1}P_2^{mn_2}...P_k^{mn_k})。
（(mn_ile mx_i)，(mn_i)可以为(0)）
考虑如果(gcd(a_i,v)=X)，那么首先有(X|a_i)，然后可以把(a_i)分解成(a_i=P_1^{c_1}P_2^{c_2}...P_k^{c_k})。如果有某个(c_i>mn_i)，那么(v)中(P_i)的次数就一定要是(mn_i)（不然(gcd)中(P_i)的次数就大于(mn_i)了），否则可以任意。
我们状压记录对于每个(a_i)，与它可以满足条件的(v)对于(iin[1,k])的(P_i)的次数是否一定要是(mn_i)。考虑到(kle15)，所以状态量也只有(2^{15})。
然后我们考虑一个(j)可以和哪些(i)匹配。同理对于(mbox{lcm}(a_j,v)=Y)，我们也可以用一个(2^{15})的状态表示(P_i)的次数是否一定要取(mx_i)。
剩下是工作就只有一个高维前缀和了。注意特判(mn_i=mx_i)的情况。
至于质因数分解。原题题解中有一种很神奇的方法。。。然而直接上(mbox{Pollard_Rho})不好吗。。。
code
#include<cstdio>
#include<algorithm>
#include<vector>
#include<ctime>
using namespace std;
#define ll long long
ll gi(){
ll x=0,w=1;char ch=getchar();
while ((ch<'0'||ch>'9')&&ch!='-') ch=getchar();
if (ch=='-') w=0,ch=getchar();
while (ch>='0'&&ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar();
return w?x:-x;
}
const int N = 2e5+5;
const int M = 15;
int n,tot,all,mx[M],mn[M],c[1<<M];
ll G,L,a[N],P[M*5],ans;
ll mul(ll x,ll y,ll mod){
x%=mod;y%=mod;ll res=0;
while (y){
if (y&1) {res+=x;if (res>=mod) res-=mod;}
x+=x;if (x>=mod) x-=mod;y>>=1;
}
return res;
}
ll fastpow(ll x,ll y,ll mod){
ll res=1;
while (y){
if (y&1) res=mul(res,x,mod);
x=mul(x,x,mod);y>>=1;
}
return res;
}
ll f[]={2,3,5,7,11,13,17,19,23,29};
bool MR(ll p){
for (int i=0;i<10;++i){
if (p<=f[i]) break;
if (fastpow(f[i],p-1,p)!=1) return false;
ll pp=p-1;
while (~pp&1){
pp>>=1;ll y=fastpow(f[i],pp,p);
if (mul(y,y,p)==1&&y!=1&&y!=p-1) return false;
}
}
return true;
}
ll PR(ll n,ll c){
ll i=0,k=2,x,y;x=y=1+rand()%(n-1);
while (1){
x=(mul(x,x,n)+c)%n;
ll d=__gcd((y-x+n)%n,n);
if (d!=1&&d!=n) return d;
if (x==y) return n;
if (++i==k) y=x,k<<=1;
}
}
void fact(ll n,ll c){
if (n==1) return;
if (MR(n)) {P[tot++]=n;return;}
ll p=n;while (p==n) p=PR(n,c--);
fact(p,233);fact(n/p,233);
}
int main(){
srand(19260817);n=gi();G=gi();L=gi();
for (int i=1;i<=n;++i) a[i]=gi();
if (L%G) return puts("0"),0;
fact(L,233);sort(P,P+tot);tot=unique(P,P+tot)-P;
for (int i=0;i<tot;++i){
ll x=G;while (x%P[i]==0) ++mn[i],x/=P[i];
x=L;while (x%P[i]==0) ++mx[i],x/=P[i];
}
all=(1<<tot)-1;
for (int i=1;i<=n;++i)
if (a[i]%G==0){
int zt=0;
for (int j=0;j<tot;++j){
ll x=a[i];int d=0;
while (x%P[j]==0) ++d,x/=P[j];
if (d>mn[j]) zt|=1<<j;
}
++c[zt];
}
for (int j=0;j<tot;++j)
for (int i=0;i<=all;++i)
if (i&(1<<j)) c[i]+=c[i^(1<<j)];// sum of subset
for (int i=1;i<=n;++i)
if (L%a[i]==0){
int zt=0;
for (int j=0;j<tot;++j){
ll x=a[i];int d=0;
while (x%P[j]==0) ++d,x/=P[j];
if (d==mx[j]||mx[j]==mn[j]) zt|=1<<j;
}
ans+=c[zt];
}
printf("%I64d
",ans);return 0;
}