Description
Solution
首先想到要把 /(/varphi(ij)/) 拆开，这里有个公式
/[/varphi(ij)=/dfrac{/varphi(i)/varphi(j)/gcd(i,j)}{/varphi(/gcd(i,j))}
/]
考虑证明，有
/[/begin{aligned} /varphi(i)/varphi(j) &= i/prod/limits_{p|i,p/in /mathtt{prime}}/frac{p-1}{p}/cdot j /prod/limits_{p|j,p/in/mathtt{prime}}/frac{p-1}{p}// &= ij/prod/limits_{p|ij,p/in/mathtt{prime}}/frac{p-1}{p}/cdot /prod/limits_{p|/gcd(i,j),p/in/mathtt{prime}}/frac{p-1}{p} /end{aligned}
/]
也就是 /(i,j/) 的并加上 /(i,j/) 的交。 接下来进入推柿子环节
/[/begin{aligned} /sum/limits_{i=1}^{n}/sum/limits_{j=1}^{m} /varphi(ij) &= /sum/limits_{i=1}^{n}/sum/limits_{j=1}^{m} /frac{/varphi(i)/varphi(j)/gcd(i,j)}{/varphi(/gcd(i,j))} //&= /sum/limits_{d=1}^{n}/sum/limits_{i=1}^{n}/sum/limits_{j=1}^{m} /frac{/varphi(i)/varphi(j)d[/gcd(i,j)=d]}{/varphi(d)} //&= /sum/limits_{d=1}^{n}/frac{d}{/varphi(d)}/sum/limits_{i=1}^{n}/sum/limits_{j=1}^{m }/varphi(i)/varphi(j)[/gcd(i,j)=d] //&= /sum/limits_{d=1}^{n}/frac{d}{/varphi(d)}/sum_{t=1}^{n/d}/mu(t)/cdot /sum/limits_{i=1}^{n/(dt)}/sum/limits_{j=1}^{m/(dt)}/varphi(idt)/varphi(jdt) //&= /sum/limits_{k=1}^{n}/sum/limits_{d|k} /frac{d/cdot/mu(k/d)}{/varphi(d)} /sum/limits_{i=1}^{n/k }/varphi(ik)/sum/limits_{j=1}^{m/k }/varphi(jk) /end{aligned}
/]
记 /(/displaystyle f(k)=/sum/limits_{d|k} /frac{d/cdot/mu(k/d)}{/varphi(d)},g(k,n)=/sum/limits_{i=1}^{n}/varphi(ik)/)，容易发现这两个函数都容易在调和级数的复杂度内预处理出来。然而 /(/mathcal O(nT)/) 仍然是不现实的。如果是 CF 我就直接上了。
考虑怎么优化，还是将优化重点放在整除分块上
/[/sum_{i=1}^n/sum_{j=1}^m/varphi(ij)=/sum_{k=1}^nf(k)/cdot g(k,n/k)/cdot g(k,m/k)
/]
记 /(/displaystyle h(a,b,n)=/sum_{k=1}^nf(k)/cdot g(k,a)/cdot g(k,b)/)，于是可以整除分块
/[/sum_{i=1}^n/sum_{j=1}^m/varphi(ij)=/sum_{n/l=n/r,/ m/l=m/r}h(n/r,m/r,r)-h(n/r,m/r,l-1)
/]
问题是，预处理 /(h/) 的复杂度是很高的，这里可以考虑 根号分治 —— 设定阈值 /(B/)，将 /(a,b/leqslant B/) 的 /(h/) 预处理出来，不难发现这是 /(/mathcal O(nB^2)/) 的。当查询的时候，若 /(n/r/leqslant B/) 就直接 /(/mathcal O(1)/) 查询，否则有 /(r/leqslant n/B/)，干脆不差分直接暴力算就是 /(/mathcal O(n/B)/) 的。
复杂度 /(/mathcal O/left(n/ln n+nB^2+/left(n^{1/2}+n/B/right)T/right)/).
块长可以用均值不等式来算，总复杂度算出来大概是 /(/mathcal O(6/cdot 10^7)/).
Code
# include <cstdio>
# include <cctype>
# define print(x,y) write(x), putchar(y)
template <class T>
inline T read(const T sample) {
T x=0; char s; bool f=0;
while(!isdigit(s=getchar())) f|=(s=='-');
for(; isdigit(s); s=getchar()) x=(x<<1)+(x<<3)+(s^48);
return f? -x: x;
}
template <class T>
inline void write(T x) {
static int writ[50], w_tp=0;
if(x<0) putchar('-'), x=-x;
do writ[++w_tp]=x-x/10*10, x/=10; while(x);
while(putchar(writ[w_tp--]^48), w_tp);
}
# include <vector>
# include <iostream>
using namespace std;
const int B = 20;
const int maxn = 1e5+5;
const int mod = 998244353;
int inv(int x,int y=mod-2,int r=1) {
for(; y; y>>=1, x=1ll*x*x%mod)
if(y&1) r=1ll*r*x%mod; return r;
}
bool is[maxn];
vector <int> g[maxn]; int h[B+1][B+1][maxn];
int phi[maxn], mu[maxn], pc, p[maxn], f[maxn];
int beelzebul(int n,int m) {
int ans=0; if(n>m) swap(n,m);
for(int i=1;i<=m/B+1;++i)
ans = (ans+1ll*f[i]*g[i][n/i]%mod*g[i][m/i]%mod)%mod;
for(int l=m/B+2, r; l<=n; l=r+1) {
r = min(n, min(n/(n/l),m/(m/l)));
ans = (0ll+ans+h[n/r][m/r][r]-h[n/r][m/r][l-1])%mod;
}
return (ans+mod)%mod;
}
int func(int i,int a,int b) {
if(a>=g[i].size() || b>=g[i].size()) return 0;
return 1ll*f[i]*g[i][a]%mod*g[i][b]%mod;
}
void sieve() {
phi[1]=mu[1]=1;
for(int i=2;i<=maxn-5;++i) {
if(!is[i]) p[++pc]=i,
mu[i]=-1, phi[i]=i-1;
for(int j=1; j<=pc && i*p[j]<=maxn-5; ++j) {
is[i*p[j]] = true, mu[i*p[j]]=-mu[i];
if(i%p[j]==0) {
phi[i*p[j]] = phi[i]*p[j];
mu[i*p[j]]=0; break;
} phi[i*p[j]] = phi[i]*(p[j]-1);
}
}
for(int i=1;i<=maxn-5;++i) {
const int coe = 1ll*i*inv(phi[i])%mod;
g[i].emplace_back(0);
for(int j=i;j<=maxn-5;j+=i)
f[j] = (1ll*coe*mu[j/i]+f[j])%mod,
g[i].emplace_back((g[i][j/i-1]+phi[j])%mod);
}
for(int a=1;a<=B;++a) for(int b=1;b<=B;++b)
for(int i=1;i<=maxn-5;++i)
h[a][b][i] = (h[a][b][i-1]+func(i,a,b))%mod;
}
int main() {
sieve();
for(int T=read(9); T; --T) {
int n=read(9), m=read(9);
print(beelzebul(n,m),'/n');
}
return 0;
}
原创文章，作者：kirin，如若转载，请注明出处：https://blog.ytso.com/tech/pnotes/271656.html