矩阵中不重复的元素
题意
一共有三个子任务。
已知一个\(m\times n\)的矩阵形如： \[ \begin{bmatrix}\\ a^b&(a+1)^b&(a+2)^b&\dots&(a+n-1)^b\\\\ a^{b+1}&(a+1)^{b+1}&(a+2)^{b+1}&\dots&(a+n-1)^{b+1}\\\\ a^{b+2}&(a+1)^{b+2}&(a+2)^{b+2}&\dots&(a+n-1)^{b+2}\\\\ \vdots&\vdots&\vdots&\ddots&\vdots\\\\ a^{b+m-1}&(a+1)^{b+m-1}&(a+2)^{b+m-1}&\dots&(a+n-1)^{b+m-1}\\\\ \end{bmatrix} \] 其中a,b,n,m给定，求这个矩阵中不同的元素个数
对于第一个子任务，\(2\le n,m,a,b\le 100\)
对于第二个子任务，\(2\le n,m,a,b\le 500000\)
对于第三个子任务，\(2\le n,m\le 5\times10^{15},2\le a,b\le 10^{15}\)
子任务1：
由于n,m很小，所以有两种方法：
方法1
由于\(\ln a^b=b\ln a\)（其实对于底数是其它的也一样），所以我们可以直接求出对数值然后用对数值去重。
方法2
我们可以把每一个底数分解质因数。
假设\(a=\prod p_{i}^{k_i}\),则\(a^b=a=\prod p_{i}^{bk_i}\)
然后去重。
子任务2：
由于n,m大了一些，所以我们要考虑其它解决方案。
我们发现如果某两行i,j的元素出现重复，那么i+a-1和j+a-1都可以表示成同一个数s的幂。
我们把每一行的底数用\(s^t\)来表示，其中s是一个正整数。
有很多种方法所以我们只考虑s最小的情况。
由于\(s\ge 2\)，所以在s相同的时候，t只有\(log_sn\)种可能。
然后就可以容斥。\(O(n)\)
子任务3：
由于n,m更大了，\(O(n)\)都跑不过。
我们可以反着考虑：从n*m中减去重复的个数。
同上，把每一行底数用\(s^t\)来表示。
那么在所有底数中，能用同一个s来表示的数的指数s一定有一个范围，并且在范围内的正整数都能取到。
如果要去重，显然这个范围内的整数个数>1。
这样我们只要考虑\(s\le \sqrt{a+n}\)的所有s。
对于每一个\(s\)，指数t的范围设为\([l,r]\)
所以\(1\le l<r\le log_2(a+n)\)，范围不大可以记忆化。
在容斥中，如果某一个数对最大最小值和lcm没有贡献的话，那这个数选不选贡献就会抵消。
如果某一个数选了后前面的某一个数就没有贡献的话，那这个数就不能选，不需要考虑。
然后我们就可以枚举最小最大值\(l\le l_1<r_1\le r\)，（当然这个也可以记忆化）
显然\(l_1\)和\(r_1\)必选。
接着对于每个大于l1小于r1的数，如果它的最大非本身约数\(\le l_1\)则可以考虑选，否则不考虑。
接着dfs，如果dfs到某个数时，这个数能被前面的lcm整除，则退出。
由于\(lcm(a,b)=\frac{ab}{\gcd(a,b)}\)，b很小，以及\(\gcd(a,b)=\gcd(b,a\mod b)\)，所以我们可以预处理gcd。
复杂度O(可过)
/*
Author: CNYALI_LK
LANG: C++
PROG: 1221.cpp
Mail: cnyalilk@vip.qq.com
*/
#include<bits/stdc++.h>
#define debug(...) fprintf(stderr,__VA_ARGS__)
#define DEBUG printf("Passing [%s] in LINE %lld\n",__FUNCTION__,__LINE__)
#define Debug debug("Passing [%s] in LINE %lld\n",__FUNCTION__,__LINE__)
#define all(x) x.begin(),x.end()
using namespace std;
const double eps=1e-8;
const double pi=acos(-1.0);
typedef long long ll;
typedef pair<ll,ll> pii;
template<class T>ll chkmin(T &a,T b){return a>b?a=b,1:0;}
template<class T>ll chkmax(T &a,T b){return a<b?a=b,1:0;}
template<class T>T sqr(T a){return a*a;}
template<class T>T mmin(T a,T b){return a<b?a:b;}
template<class T>T mmax(T a,T b){return a>b?a:b;}
template<class T>T aabs(T a){return a<0?-a:a;}
#define min mmin
#define max mmax
#define abs aabs
ll read(){
ll s=0,base=1;
char c;
while(!isdigit(c=getchar()))if(c=='-')base=-base;
while(isdigit(c)){s=s*10+(c^48);c=getchar();}
return s*base;
}
char WritellBuffer[1024];
template<class T>void write(T a,char end){
ll cnt=0,fu=1;
if(a<0){putchar('-');fu=-1;}
do{WritellBuffer[++cnt]=fu*(a%10)+'0';a/=10;}while(a);
while(cnt){putchar(WritellBuffer[cnt]);--cnt;}
putchar(end);
}
const ll p=1000000007;
ll lx,rx,ly,ry;
bool isn[80000000];
ll sel[80],cnt,lf,rf;
ll g[80][80],asa[80],top,mxd[80];//预处理gcd
ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}
ll dfs(ll nw,ll lm,ll flag){
if(lf/lm>=rf/lm)return 0;
if(nw>top){
ll fa=(rf/lm-lf/lm)*flag%p;
if(fa<0)return fa+p;
return fa;
}
if(!(lm%asa[nw]))return 0;//如果一个数对lcm没有贡献那么退出
ll ans=dfs(nw+1,lm,flag);
ll as=dfs(nw+1,lm*(asa[nw]/g[lm%asa[nw]][asa[nw]]),-flag);
return (ans+as)%p;
}
ll f[80][80],vis[80][80];
ll calc(ll l,ll r){
ll ans=0;
for(ll lfx=l;lfx<=r;++lfx){
rf=lfx*ry;
for(ll rfx=lfx+1;rfx<=r;++rfx){
if(vis[lfx][rfx])ans=(ans+f[lfx][rfx]) %p;//记忆化
else{
lf=rfx*ly-1;
if(lf>rf)continue;//剪枝
top=0;
for(int i=lfx+1;i<rfx;++i){//去掉某些不可能选的元素
if(mxd[i]<=lfx){
asa[++top]=i;
}
}
ll a=dfs(1,lfx*rfx/g[lfx][rfx],1);
ans=ans+a;
if(ans>=p)ans-=p;
vis[lfx][rfx]=1;
f[lfx][rfx]=a;
}
}
}
return ans;
}
ll CNT[100][100];
int main(){
#ifdef cnyali_lk
freopen("1221.in","r",stdin);
freopen("1221.out","w",stdout);
#endif
ll n,m,a,b;
m=read();n=read();a=read();b=read();
for(ll i=0;i<=70;++i)for(ll j=0;j<=70;++j){
g[i][j]=gcd(i,j);
if(i&&i<j&&!(j%i))mxd[j]=i;
}
lx=a,rx=n+a-1,ly=b,ry=b+m-1;
ll Block=(ll)(ceil(sqrt(rx)));
ll ans=(n%p)*(m%p)%p;
for(ll i=2;i<=Block;++i){
if(!isn[i]){
ll Lj=0,Rj=0;
for(unsigned long long j=1,k=i;k<=rx;k*=i,++j){
if(k<=Block)isn[k]=1;
if(lx<=k){
Rj=j;
if(!Lj)Lj=j;
}
unsigned long long ab=k;
ab*=i;
if((ab/i)!=k)break;
}
if(Lj<Rj){
++CNT[Lj][Rj];
}
}
// debug("End%lld\n",i);
}
for(ll i=1;i<=70;++i)for(ll j=i+1;j<=70;++j)if(CNT[i][j])ans=(ans-calc(i,j)*CNT[i][j]%p+p)%p;
printf("%lld\n",ans);
return 0;
}