毒瘤之神的考验
一、题目
二、解法
\(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(n\log n)\) 预处理出来,不多说。
按照套路这时候应该上整除分块了,但是发现后面的 \(\varphi\) 和 \(d\) 有关所以不能分块,那么仔细观察一下有没有优化的余地,发现后面枚举是到了 \(n/d\) 的,我们枚举 \(d\) 计算后面的复杂度就相当于调和级数求和,那么可以预处理。
设 \(g(d,n)=\sum_{i=1}^n \varphi(id)\) ,由于 \(dn\leq 100000\) ,所以可以 \(O(n\log n)\) 预处理,现在写出柿子:
\(O(n)\) 求显然是不行的,后面还有参数又不能分块。下面的就是奇技淫巧了,我们把这整个结构预处理成一个数组,让这个数组和 \(d\) 无关但和 \(n/d\) 相关,然后用整除分块,再平衡预处理和在线的复杂度,听起来很玄学是不是,这简直就是硬套整除分块!下面来详细讲一讲这种骚操作:
设 \(s(a,b,n)=\sum_{i=1}^nf(i)\times g(i,a)\times g(i,b)\) ,那么 \(d\in[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(T\frac{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\n",(ans+MOD)%MOD);
}
signed main()
{
init(100000);
T=read();
while(T--)
{
n=read();m=read();
if(n>m) swap(n,m);
work();
}
}