一、题目
二、解法
(i,j)合起来特别难受啊，看能不能把 (i,j) 拆开，由于 (varphi) 是积性函数，所以先拆成 (varphi(i)varphi(j)) 的形式，但是这么拆显然会错，错误出在对于 (i,j) 共同的公因数 (p_i) ，计入了两次 (p_i-1) ，所以我们除以 (varphi(gcd(ij))) 之后再乘上 (gcd(i,j)) 就解决了这个问题，并转化成了莫比乌斯反演的题，关系式是这样的：
这里再多说几句，我认为把上面的柿子理解成一个神奇的结论是不妥帖的，我更觉得他是一种构造，把陌生的问题转化成了莫比乌斯反演的题，推这个柿子很简单，我直接把我自己推到的草稿给出来了哦：
然后令 (f(n)=sum_{d|n}mu(frac{n}{d})frac{d}{varphi(d)}) ，这个东西可以 (O(nlog n)) 预处理出来，不多说。
按照套路这时候应该上整除分块了，但是发现后面的 (varphi) 和 (d) 有关所以不能分块，那么仔细观察一下有没有优化的余地，发现后面枚举是到了 (n/d) 的，我们枚举 (d) 计算后面的复杂度就相当于调和级数求和，那么可以预处理。
设 (g(d,n)=sum_{i=1}^n varphi(id)) ，由于 (dnleq 100000) ，所以可以 (O(nlog n)) 预处理，现在写出柿子：
(O(n)) 求显然是不行的，后面还有参数又不能分块。下面的就是奇技淫巧了，我们把这整个结构预处理成一个数组，让这个数组和 (d) 无关但和 (n/d) 相关，然后用整除分块，再平衡预处理和在线的复杂度，听起来很玄学是不是，这简直就是硬套整除分块！下面来详细讲一讲这种骚操作：
设 (s(a,b,n)=sum_{i=1}^nf(i) imes g(i,a) imes g(i,b)) ，那么 (din[l,r]) 的答案可以表示为：
直接预处理 (s) 升天，现在我们来平衡一下询问和预处理的复杂度，发现在 (l) 比较小的时候 (n/l) 比较大，但是在整除分块中那一段就很少，这就血亏，(l) 大的时候 (n/l) 比较小，这时候血赚。所以我们设置一个阈值 (B) ，用类似分块分析复杂度的方法，当 (frac{n}{l}leq B) 的时候我们预处理 (s) ，时间复杂度 (O(nB^2))
当 (frac{n}{l}>B) 的时候即 (frac{n}{B}>l) 的情况下，直接枚举 (d) 暴力求即可，时间复杂度 (O(Tfrac{n}{B}))
有两个细节在代码里注释的建议一定要看一看，实际情况下 (B=50) 就行了
#include <cstdio>
#include <vector>
#include <cstring>
#include <iostream>
using namespace std;
const int N = 100;
const int M = 100005;
const int MOD = 998244353;
#define int long long
int read()
{
int x=0,f=1;char c;
while((c=getchar())<'0' || c>'9') {if(c=='-') f=-1;}
while(c>='0' && c<='9') {x=(x<<3)+(x<<1)+(c^48);c=getchar();}
return x*f;
}
int T,n,m,ans,cnt,p[M],mu[M],phi[M],inv[M],f[M];
int *g[M],*s[N+5][N+5];
void init(int n)
{
mu[1]=phi[1]=1;
for(int i=2;i<=n;i++)
{
if(!phi[i])
{
mu[i]=-1;
phi[i]=i-1;
p[++cnt]=i;
}
for(int j=1;j<=cnt && i*p[j]<=n;j++)
{
if(i%p[j]==0)
{
phi[i*p[j]]=p[j]*phi[i];
break;
}
phi[i*p[j]]=phi[i]*phi[p[j]];
mu[i*p[j]]=-mu[i];
}
}
inv[0]=inv[1]=1;
for(int i=2;i<=n;i++) inv[i]=(MOD-MOD/i)*inv[MOD%i]%MOD;
for(int i=1;i<=n;i++)
for(int j=i;j<=n;j+=i)
f[j]=(f[j]+i*mu[j/i]*inv[phi[i]])%MOD;
for(int i=1;i<=n;i++)
{
g[i]=new int [n/i+3];g[i][0]=0;
for(int j=1;j<=n/i;j++)
g[i][j]=(g[i][j-1]+phi[i*j])%MOD;
}
for(int i=1;i<=N;i++)
for(int j=i;j<=N;j++)
{
int len=n/j;
//这是用来减少空间的，你看下面我们用到 r
//r=n/(n/l)=n/j 明白了吧 不用处理那么多
s[i][j]=new int [len+3];s[i][j][0]=0;
for(int k=1;k<=len;k++)
s[i][j][k]=(s[i][j][k-1]+f[k]*g[k][i]%MOD*g[k][j])%MOD;
}
}
void work()
{
ans=0;
for(int i=1;i<=m/N;i++) ans=(ans+f[i]*g[i][n/i]%MOD*g[i][m/i])%MOD;
//这里打n要出问题,因为后面有 m/l 所以要用 m/N
for(int l=m/N+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])%MOD;
}
printf("%lld
",(ans+MOD)%MOD);
}
signed main()
{
init(100000);
T=read();
while(T--)
{
n=read();m=read();
if(n>m) swap(n,m);
work();
}
}