题目链接
先考虑 $x^a=y^b$ 意味着什么
记 $p(x)$ 为最小的 $p$, 使得 $\exists k , p^k=x$。
那么 $x^a=p(x)^{k(x) \times a}$
则 $ x^a=y^b $ 等价于 $ p(x)^{k(x) \times a}=p(y)^{k(y) \times b}$ 等价于 $ p(x) = p(y) \land k(x) \times a = k(y) \times b $
考虑枚举 $p$ ,计算所有 $p(x)=p$ 的 $x$ 带来的重复
那么 使 $x$ 在给定区间 的 $k(x)$ 在一个区间中
于是问题转化为,求 $ card \{ i \times j \, | \, i \in [l_1,r_1] \, , j \in [l_2,r_2] \} $其中$l_1$,$r_1$为$O(log_2(n))$级别,$l_2$,$r_2$为 $O(n)$ 级别
容斥。
先加上在 $[l_1,r_1]$ 中任选一个数,并在 $[l_2,r_2]$ 任选一个数的方案数
再减去在 $[l_1,r_1]$ 任选两个数 $x,y$,并在 $[l_2,r_2]$ 任选两个数 $a,b$ ,使得$x \times a=y \times b$的方案数
再加上在 $[l_1,r_1]$ 任选三个数 $x,y,z$,并在 $[l_2,r_2]$ 任选三个数 $a,b,c$ ,使得$x \times a=y \times b=z \times c$的方案数
......
可以发现 当 $[l_1,r_1]$ 中的数确定之后, $[l_2,r_2]$ 的方案数是可以算出来的 因为乘积一定为 $lcm$ 的倍数
设为 $x$ 倍 则$l_2 \le x \times lcm/max \, \& \& \, x \times lcm/min \le r_2$ 于是 $l_2 \times max \le x \times lcm \le r_2 \times min$
所以 $x$ 的取值数=$r_2 \times min/lcm - l_2 \times max/lcm +1$
于是暴力容斥复杂度 $=O(2^{r_1-l_1+1})=O(n)$
$l_1$,$r_1$一共有 $ log_2^2(n) $种
对于所有的一样的$l_1$,我们可以从小到大枚举$r_1$计算。 总复杂度还是 $O(n)$
当然还过不了
注意到 如果一个数选了之后对于$lcm,max,min$没有影响 那么他选和不选的贡献就会抵消
所以先枚举$max,min$ 之后$dfs$ 如果碰到对$lcm$无贡献的数,就直接退出
同时,如果选了它导致前面的数成为对$lcm$无贡献的数,那么它就不选
这样优化之后复杂度就很好了 李陶冶说"比较接近于光滑数的数量(其实比这个还要低一个数量级)"
但我不知道光滑数是什么东西 总之O(跑的过)
还有个问题,就是如何求出有多少个$p(x)$会使得$k(x)$在 $[l_1,r_1]$ 中
我是筛出$ \sqrt{n} $以内所有$p(x)=x$的$x$ 然后枚举 求更优的方法?
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int D=1e9+7,Lo=52;
int cnt[Lo+2][Lo+2];
ll L1,R1,l2,r2,ans,now;
void init()
{
ll m,n,a,b;
cin>>m>>n>>a>>b;
ans=(m%D)*(n%D)%D;
L1=a;R1=a+n-1;l2=b;r2=b+m-1;
}
int Gcd[Lo+2][Lo+2];
int gcd(int x,int y)
{
while(y)swap(x%=y,y);
return x;
}
int l1,r1,mn;
ll L,R;//L(=l2*mx)<=x*Lcm<=R(=r2*mn)
int mx_p[Lo+5],st[Lo+5],*top;
#define len (R/Lcm - (L+Lcm-1)/Lcm +1)
void dfs(int *h,const ll &Lcm,bool ji)
{
if(Lcm>R)return ;
if(h>top)
{
if(ji)now-=len;
else now+=len;
return ;
}
int x=*h;
if(Lcm%x==0) return ;
//if(len<=0) return ;
dfs(h+1,Lcm,ji);
dfs(h+1,Lcm*(x/Gcd[x][Lcm%x]),ji^1);
}
const int U=8e7;
bool mark[U];
int main()
{
//freopen("1.in","r",stdin);freopen("1.out","w",stdout);
init();
int mp=sqrt(R1);
for(int p=2;p<=mp;++p)
if(!mark[p])
{
int k=0;
ll x=1;
while(x<L1&&x<=R1/p) { x*=p;if(x<=mp)mark[x]=1;++k; }
if(x<L1||x>R1)continue;
int l1=k;
while(x<=R1/p) { x*=p;if(x<=mp)mark[x]=1;++k; }
if(k>l1)
++cnt[l1][k];
}
for(int i=1;i<=Lo;++i)
for(int j=0;j<i;++j)
Gcd[i][j]=gcd(i,j);
for(int i=1;i<=Lo;++i)
for(int j=1;j<i;++j)
if(i%j==0) mx_p[i]=j;
for(l1=1;l1<=Lo;++l1)
{
int mr=0;
for(r1=Lo;r1>l1;--r1)
if(cnt[l1][r1]){mr=r1;break;}
now=0;
for(r1=l1+1;r1<=mr;++r1)
{
L=l2*r1;
for(mn=l1;mn<r1;++mn)
{
R=r2*mn;
if(L>R)continue;
int Lcm=r1*(mn/gcd(r1,mn));
int i=mn+1;
for(;i<r1;++i)
if(Lcm%i==0) break;
if(i<r1)continue;
top=st-1;
for(int i=mn+1;i<r1;++i)
if(mx_p[i]<=mn)
*++top=i;
dfs(st,Lcm,0);
}
now%=D;
ans-=cnt[l1][r1]*now%D;
}
}
cout<<(ans%D+D)%D;
}