BSOJ4129【BZOJ2693】jzptab
问题
结论
结论一:
结论二:(进一步优化)
分析
先将式子化简:
这时我们单独计算后半部分,这里没有使用莫比乌斯反演而是直接用 \(\large \mu*I=\varepsilon\) 的性质:(这部分其实就是上一个问题)
其中,\(\large G(x)=\frac {x(x+1)}{2}\) ,表示等差数列求和。
很容易发现,对于 \(\large S\) 函数,可以通过整除分块来计算,按照 \(\large\lfloor\frac nd\rfloor\) 以及 \(\large \lfloor\frac md\rfloor\) 来分块即可,前面的部分就是预处理一下前缀和即可。
现在代回原式子:
然后我们可以发现,原式还需要一次除法分块,而第二次除法分块将借助第一次的作为上界进行分块,这就是二次除法分块,这样做的复杂度可以证明是 \(\large O(n^{\frac 34})\) 的。
接下来考虑进一步优化:
考虑枚举 \(\large dk\) :
交换求和符号(先枚举一个数再枚举其倍数等价于先枚举其倍数再枚举倍数的因数):
考虑先枚举 \(\large \frac xd\) :
直接把 \(\large x\) 提出来:
于是后面那部分就是 \(\large F(x)=\sum\limits_{d\vert x}id(x)\mu(x)\) ,因为 \(id\) 和 \(\mu\) 都是积性函数,所以 \(F\) 也是积性函数。
所以可以考虑线性筛出来,前面部分就是一个整除分块,直接做即可。
时间复杂度是 \(\large O(V+m\sqrt{V})\) ,其中 \(V\) 是值域,\(m\) 是询问次数,如果使用杜教筛可以优化到 \(\large O(V^{\frac 23}+m\sqrt{V})\) 。
代码
#include<bits/stdc++.h>
using namespace std;
template <typename T>
inline void read(T &x){
x=0;bool f=false;char ch=getchar();
while(!isdigit(ch)){f|=ch=='-';ch=getchar();}
while(isdigit(ch)){x=x*10+(ch^48);ch=getchar();}
x=f?-x:x;
return ;
}
template <typename T>
inline void write(T x){
if(x<0) x=-x,putchar('-');
if(x>9) write(x/10);
putchar(x%10^48);
return ;
}
#define int long long
#define ll long long
#define ull unsigned long long
#define rep(i,x,y) for(int i=(x);i<=(y);i++)
#define dep(i,y,x) for(int i=(y);i>=(x);i--)
const int N=1e7+5,M=2e5+5,MOD=1e8+9;
int n,m,k;
int prime[N],f[N],cnt,pre[N],low[N];
ll val[N],Ans;
bool isprime[N];
inline int inc(int x,int y){x+=y;return x>=MOD?x-MOD:x;}
inline void incc(int &x,int y){x+=y;if(x>=MOD) x-=MOD;}
bool vis[N];
inline void GetPrimes(int d){
vis[1]=true;f[1]=low[1]=pre[1]=1;//对1进行定义
for(ll i=2;i<=d;i++){
if(!vis[i]) prime[++cnt]=low[i]=i,f[i]=1-i;//对质数进行定义
pre[i]=inc(pre[i-1],i*f[i]%MOD);
for(ll j=1;j<=cnt&&i*prime[j]<=d;j++){
vis[i*prime[j]]=true;
if(i%prime[j]==0){
low[i*prime[j]]=low[i]*prime[j];
if (low[i]==i) f[i*prime[j]]=1-prime[j];//对质数的若干次幂进行定义(一般由f[i]递推);
else f[i*prime[j]]=f[i/low[i]]*f[low[i]*prime[j]];//把当前质因子的所有项剃掉,然后把新的加回来
break;
}
low[i*prime[j]]=prime[j];
f[i*prime[j]]=f[i]*f[prime[j]];//通过 f(n),f(m) 计算 f(nm)
}
}
return ;
}
inline int GetSum(int x){return (1ll*x*(x+1)/2)%MOD;}
inline int Sum(int x,int y){
int res=0;
if(x>y) swap(x,y);
for(int l=1,r;l<=x;l=r+1){
r=min(x/(x/l),y/(y/l));
incc(res,1ll*(inc((pre[r]-pre[l-1])%MOD,MOD))*GetSum(x/l)%MOD*GetSum(y/l)%MOD);
}
return res;
}
inline int Calc(int x,int y){
int res=0;
for(int l=1,r;l<=x;l=r+1){
r=min(x/(x/l),y/(y/l));
incc(res,1ll*inc((pre[r]-pre[l-1])%MOD,MOD)%MOD*GetSum(x/l)%MOD*GetSum(y/l)%MOD);
}
return res;
}
signed main(){
GetPrimes(1e7);
int t;read(t);
while(t--){
read(n),read(m);
if(n>m) swap(n,m);
write(Calc(n,m)),putchar('\n');
}
return 0;
}