不妨设(nle m)
[egin{aligned}
ans&=sum_{i=1}^nsum_{j=1}^mvarphi(ij)\
&=sum_{i=1}^nsum_{j=1}^mfrac{varphi(i)varphi(j)gcd(i,j)}{varphi(gcd(i,j))}\
&=sum_{d=1}^{n}frac{d}{varphi(d)}sum_{i=1}^nsum_{j=1}^mvarphi(i)varphi(j)[d=gcd(i,j)]
\
f(d)&=sum_{i=1}^nsum_{j=1}^mvarphi(i)varphi(j)[d=gcd(i,j)]
\
g(d)&=sum_{d|x}f(x)\
&=sum_{i=1}^nsum_{j=1}^mvarphi(i)varphi(j)[d|gcd(i,j)]\
&=sum_{i=1}^{n/d}varphi(id)sum_{j=1}^{m/d}varphi(jd)
\
f(d)&=sum_{d|x}mu(frac{x}d)g(x)
\
ans&=sum_{d=1}^{n}frac{d}{varphi(d)}f(d)\
&=sum_{d=1}^{n}frac{d}{varphi(d)}sum_{d|x}mu(frac{x}d)g(x)\
\
end{aligned}
]
预处理(g(d)),然后暴力做(ans),复杂度都是(nlog n)的。
是的这并不能过……
[egin{aligned}
ans&=sum_{d=1}^nfrac{d}{varphi(d)}sum_{d|x}mu(frac{x}d)sum_{i=1}^{n/x}varphi(ix)sum_{j=1}^{m/x}varphi(jx)\
&=sum_{x=1}^nsum_{i=1}^{n/x}varphi(ix)sum_{j=1}^{m/x}varphi(jx)sum_{d|x}frac{d}{varphi(d)}mu(frac{x}d)\
G(x,t)&=sum_{i=1}^tvarphi(ix)\
&=G(x,t-1)+varphi(tx)\
F(x)&=sum_{d|x}frac{d}{varphi(d)}mu(frac{x}d)\
S(n,m,t)&=sum_{x=1}^tG(x,n/x)G(x,m/x)F(x)\
&=S(n,m,t-1)+G(t,n/t)G(t,m/t)F(t)
end{aligned}
]
埃氏筛求出(G,F)以及一部分的(S),剩下的暴力算。如果设块大小为(B),则预处理、计算的(S)复杂度为(nB^2+T(sqrt n-sqrt{n/B} +n/B)),然后实践出真知地求出较优的(B)。
#include <bits/stdc++.h>
using namespace std;
const int N=1e5;
const int B=40;
const int P=998244353;
bool vis[N+1];
int pri[N+1],phi[N+1],mu[N+1],tot;
int inv[N+1];
int F[N+1],*G[N+1],*S[B+1][B+1];
void initial() {
phi[1]=mu[1]=1;
for(int i=2; i<=N; ++i) {
if(!vis[i]) pri[++tot]=i,phi[i]=i-1,mu[i]=-1;
for(int j=1; j<=tot&&i*pri[j]<=N; ++j) {
vis[i*pri[j]]=1;
if(i%pri[j]==0) {phi[i*pri[j]]=phi[i]*pri[j]; break;}
phi[i*pri[j]]=phi[i]*(pri[j]-1);
mu[i*pri[j]]=-mu[i];
}
}
inv[1]=1;
for(int i=2; i<=N; ++i) inv[i]=1LL*inv[P%i]*(P-P/i)%P;
for(int i=1; i<=N; ++i)
for(int j=1; j<=N/i; ++j)
F[i*j]=(F[i*j]+1LL*i*inv[phi[i]]%P*mu[j]%P+P)%P;
for(int i=1; i<=N; ++i) {
G[i]=new int[N/i+1];
G[i][0]=0;
for(int j=1; j<=N/i; ++j)
G[i][j]=(G[i][j-1]+phi[i*j])%P;
}
for(int x=1; x<=B; ++x)
for(int y=1; y<=B; ++y) {
int T=N/max(x,y);
S[x][y]=new int[T+1];
S[x][y][0]=0;
for(int t=1; t<=T; ++t)
S[x][y][t]=(S[x][y][t-1]+1LL*F[t]*G[t][x]%P*G[t][y]%P)%P;
}
}
int solve(int n,int m) {
if(m<n) swap(n,m);
int ans=0;
for(int i=1; i<=m/B; ++i)
ans=(ans+1LL*F[i]*G[i][n/i]%P*G[i][m/i]%P)%P;
for(int l=m/B+1,r; l<=n; l=r+1) {
r=min(n/(n/l),m/(m/l));
ans=(ans+(S[n/l][m/l][r]-S[n/l][m/l][l-1]+P)%P)%P;
}
return ans;
}
int main() {
initial();
int T,n,m;
scanf("%d",&T);
while(T--) {
scanf("%d%d",&n,&m);
printf("%d
",solve(n,m));
}
return 0;
}