莫反做题记

一些（两个）常用结论：

$$\sum _{d|n}\mu (d)=[n=1]$$

$$\sum _{d|n}\phi (d)=n$$

反演公式:

给定函数 $f(x)$ ，定义函数 $g(n)=\sum _{d|n}f(d)$

则有：

$$g(x)=\sum _{d|n}\mu (\frac{x}{d})f(d)$$

[POI2007] ZAP-Queries

即：

$$\begin{array}{rl}ans& =\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=d]\\ & =\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }[gcd(i,j)=1]\\ & =\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\sum _{k|gcd(i,j)}\mu (k)\\ & =\sum _{k=1}^{\lfloor \frac{n}{d}\rfloor }\mu (k)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }[k|i][k|j]\\ & =\sum _{k=1}^{\lfloor \frac{n}{d}\rfloor }\mu (k)\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor \end{array}$$

然后再除法分块即可。
或者令 $T=kd$ ，则：

$$\begin{array}{rl}ans& =\sum _{d|T}^n\mu (\frac{T}{d})\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \end{array}$$

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
const int N=5e4+100;
int miu[N],prime[N],tot,t,n,m,d;
bool vis[N];
void get_miu(){
	miu[1]=1;vis[1]=1;
	for(int i=2;i<N;i++){
		if(!vis[i])prime[++tot]=i,miu[i]=-1;
		for(int j=1;j<=tot&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j])miu[i*prime[j]]=-miu[i];
			else{
				miu[i*prime[j]]=yhl;break;
			}
		}
	}
	for(int i=1;i<N;i++)miu[i]+=miu[i-1];
}
int ghfz(int n,int m){
	int ans=yhl,r=yhl;
	for(int i=1;i<=n;i=r+1){
		r=min(n/(n/i),m/(m/i));
		ans+=(miu[r]-miu[i-1])*(n/i)*(m/i);
	}
	return ans;
}
int main(){
	get_miu();
	scanf("%d",&t);
	while(t--){
		scanf("%d%d%d",&n,&m,&d);
		n/=d,m/=d;if(m<n)swap(n,m);
		printf("%d\n",ghfz(n,m));
	}
	return yhl;
}

$\mathrm{三}\mathrm{倍}\mathrm{经}\mathrm{验}$ 双亲数 [HAOI2011] Problem b，后面这个拿容斥搞一下即可。

YY的GCD

$$\begin{array}{rl}ans& =\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=p,p\in prime]\\ & =\sum _{p\in prime}\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=p]\\ & =\sum _{p\in prime}\sum _{i=1}^{\lfloor \frac{n}{p}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{p}\rfloor }[gcd(i,j)=1]\\ & =\sum _{p\in prime}\sum _{i=1}^{\lfloor \frac{p}{p}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{p}\rfloor }\sum _{k|gcd(i,j}\mu (k)\\ & =\sum _{p\in prime}\sum _{k=1}^{\lfloor \frac{n}{p}\rfloor }\mu (k)\sum _{i=1}^{\lfloor \frac{n}{p}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{p}\rfloor }[k|i][k|j]\\ & =\sum _{p\in prime}\sum _{k=1}^{\lfloor \frac{n}{p}\rfloor }\mu (k)\lfloor \frac{n}{pk}\rfloor \lfloor \frac{m}{pk}\rfloor \end{array}$$

令 $T=pk$ ，再变为枚举 $T$ 则：

$$\begin{array}{rl}ans& =\sum _{T=1}^n\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{k|T,k\in prime}\mu (\frac{T}{k})\end{array}$$
后面部分预处理出来，然后拿前缀和维护一下即可，再用根号分治就好了。

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=1e7+100;
int miu[N],prime[N],tot,t,n,m,d,f[N],s[N];
bool vis[N];
void get_miu(){
	miu[1]=1;vis[1]=1;
	for(int i=2;i<N;i++){
		if(!vis[i])prime[++tot]=i,miu[i]=-1;
		for(int j=1;j<=tot&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j])miu[i*prime[j]]=-miu[i];
			else{miu[i*prime[j]]=yhl;break;}
		}
	}
	for(int i=1;i<=tot;i++)
		for(int j=1;j*prime[i]<N;j++)
			f[j*prime[i]]+=miu[j];
	for(int i=1;i<N;i++)s[i]=s[i-1]+f[i];
}
int ghfz(int n,int m){
	int ans=yhl,r=yhl;
	for(int i=1;i<=n;i=r+1){
		r=min(n/(n/i),m/(m/i));
		ans+=(s[r]-s[i-1])*(n/i)*(m/i);
	}
	return ans;
}
signed main(){
	get_miu();
	scanf("%lld",&t);
	while(t--){
		scanf("%lld%lld",&n,&m);
		if(m<n)swap(n,m);
		printf("%lld\n",ghfz(n,m));
	}
	return yhl;
}

[SDOI2015] 约数个数和

即 $ans=\sum _{i=1}^n\sum _{j=1}^{m}d(ij)$ ，其中 $d(x)$ 为 $x$ 的约数个数。

一个结论：$d(ij)=\sum _{a|i}\sum _{b|j}[gcd(a,b)=1]$ ，证明我也不知道😕。

那么：

$$\begin{array}{rl}ans& =\sum _{i=1}^n\sum _{j=1}^m\sum _{a|i}\sum _{b|j}[gcd(a,b)=1]\\ & =\sum _{a=1}^n\sum _{b=1}^m\sum _{i=1}^n\sum _{j=1}^m[a|i][b|j][gcd(a,b)=1]\\ & =\sum _{a=1}^n\sum _{b=1}^m\lfloor \frac{n}{a}\rfloor \lfloor \frac{m}{b}\rfloor [gcd(a,b)=1]\\ & =\sum _{a=1}^n\sum _{b=1}^m\lfloor \frac{n}{a}\rfloor \lfloor \frac{m}{b}\rfloor \sum _{d=1}^n\mu (d)[d|a][d|b]\\ & =\sum _{d=1}^n\mu (d)\sum _{a=1}^{\lfloor \frac{n}{d}\rfloor }\lfloor \frac{n}{ad}\rfloor \sum _{b=1}^{\lfloor \frac{m}{d}\rfloor }\lfloor \frac{m}{bd}\rfloor \\ & =\sum _{d=1}^n\mu (d)\sum _{a=1}^{\lfloor \frac{n}{d}\rfloor }\lfloor \frac{\lfloor \frac{n}{d}\rfloor }{a}\rfloor \sum _{b=1}^{\lfloor \frac{m}{d}\rfloor }\lfloor \frac{\lfloor \frac{m}{d}\rfloor }{b}\rfloor \end{array}$$

可以观察到，后面两部分形式相同，可以 $O(n)$ 预处理出来，再维护一个 $\mu$ 的前缀和，然后根号分治即可。时间复杂度 $O(T\sqrt{n})$ 。

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=5e4+100;
int miu[N],prime[N],tot,t,n,m,d,f[N];
bool vis[N];
void get_miu(){
	miu[1]=1;vis[1]=1;
	for(int i=2;i<N;i++){
		if(!vis[i])prime[++tot]=i,miu[i]=-1;
		for(int j=1;j<=tot&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j])miu[i*prime[j]]=-miu[i];
			else{miu[i*prime[j]]=yhl;break;}
		}
	}
	for(int i=1;i<N;i++){
		int ans=0,r=0;
		for(int j=1;j<=i;j=r+1){
			r=(i/(i/j));
			ans+=(r-j+1)*(i/j);
		}
		f[i]=ans;
	}
	for(int i=1;i<N;i++)miu[i]+=miu[i-1];
}
int ghfz(int n,int m){
	int ans=yhl,r=yhl;
	for(int i=1;i<=n;i=r+1){
		r=min(n/(n/i),m/(m/i));
		ans+=(miu[r]-miu[i-1])*f[n/i]*f[m/i];
	}
	return ans;
}
signed main(){
	get_miu();
	scanf("%lld",&t);
	while(t--){
		scanf("%lld%lld",&n,&m);
		if(m<n)swap(n,m);
		printf("%lld\n",ghfz(n,m));
	}
	return yhl;
}

[SDOI2017] 数字表格

首先，斐波那契数列有个很神奇的性质（下文用$F$表示斐波那契数列）：若 $F_i|F_j$ ，则 $i|j$ 。（其实反过来也成立）。
那么：

$$\begin{array}{rl}ans& =\prod _{i=1}^n\prod _{j=1}^m{F}_{gcd(i,j)}\\ & =\prod _{i=1}^n\prod _{j=1}^{m}gcd(F_i,F_j)\\ & =\prod _{d=1}^n{F}_d\prod _{i=1}^n\prod _{j=1}^m[gcd(i,j)=d]\\ & =\prod _{d=1}^n{F}_d^{\sum _{i=1}^n\sum _{j=1}^{m}gcd(i,j)=d}\\ & =\prod _{d=1}^n{F}_d^{\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }[gcd(i,j)=1]}\\ & =\prod _{d=1}^n{F}_d^{\sum _{k=1}^{\lfloor \frac{n}{d}\rfloor }\mu (k)\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor }\end{array}$$

与上上题相似，我们仍令 $T=kd$ ，则：

$$\begin{array}{rl}ans& =\prod _{d=1}^n{F}_d^{\sum _{d|T}^n\mu (\frac{T}{d})\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor }\\ & =\prod _{T=1}^n\prod _{d|T}(F_d^{\mu (\frac{T}{d})}{)}^{\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor }\end{array}$$

前半部分用根号分治解决，后半部分用前缀积解决即可。

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
#define int long long 
const int N=1e6+10,p=1e9+7;
int f[N],g[N],n,t,m,prime[N],miu[N],F[N];
bool vis[N];
inline int qum(int a,int b){
	int ans=1;
	while(b){
		if(b&1)ans=1ll*ans*a%p;
		a=1ll*a*a%p;
		b>>=1;
	}
	return ans;
}
void init(){
	f[1]=g[1]=F[yhl]=F[1]=1;
	for(int i=2;i<N;i++)f[i]=(f[i-1]+f[i-2])%p;
	vis[1]=1;miu[1]=1;
	for(int i=2;i<N;i++){
		g[i]=qum(f[i],p-2);F[i]=1;
		if(!vis[i])prime[++prime[yhl]]=i,miu[i]=-1;
		for(int j=1;j<=prime[yhl]&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j])miu[i*prime[j]]=-miu[i];
			else{miu[i*prime[j]]=yhl;break;}
		}
	}
	for(int i=1;i<N;i++){
		if(!miu[i])continue;
		for(int j=i;j<N;j+=i){
			F[j]=1ll*F[j]*(miu[i]==1?f[j/i]:g[j/i])%p;
		}
	}
	for(int i=2;i<N;i++)F[i]=1ll*F[i]*F[i-1]%p;
}
int ghfz(int n,int m){
	int ans=1,r=yhl;
	for(int i=1;i<=n;i=r+1){
		r=min(n/(n/i),m/(m/i));
		int ny=1ll*F[r]*qum(F[i-1],p-2)%p;
		//注意用逆元来实现前缀积。 
		ans=1ll*ans*qum(ny,1ll*(n/i)*(m/i)%(p-1))%p;
		//这里用费马小定理省去一些不必要的计算。 
	}
	return ans;
}
signed main(){
	init();
	scanf("%lld",&t);
	while(t--){
		scanf("%lld%lld",&n,&m);
		if(n>m)swap(n,m);
		printf("%lld\n",ghfz(n,m));
	}
	return yhl;
}

数表

即：

$$\begin{array}{rl}ans& =\sum _{i=1}^n\sum _{j=1}^m[{\sigma }_{gcd(i,j)<=a}]{\sigma }_{gcd(i,j)}\\ & =\sum _{d=1}^n[{\sigma }_d<=a]{\sigma }_d\sum _{i=1}^n\sum _{j=1}^m[gcd(i,j)=d]\\ & =\sum _{d=1}^n[{\sigma }_d<=a]{\sigma }_d\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }\sum _{k|gcd(i,j)}\mu (k)\\ & =\sum _{d=1}^n[{\sigma }_d<=a]{\sigma }_d\sum _{k=1}^{\lfloor \frac{n}{d}\rfloor }\mu (k)\lfloor \frac{n}{kd}\rfloor \lfloor \frac{m}{kd}\rfloor \end{array}$$

《典》。仍令 $T=kd$ ，则：

$$\begin{array}{rl}ans& =\sum _{d=1}^n[{\sigma }_d<=a]{\sigma }_d\sum _{T=1}^n\mu (\lfloor \frac{T}{d}\rfloor )\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \\ & =\sum _{T=1}^n\lfloor \frac{n}{T}\rfloor \lfloor \frac{m}{T}\rfloor \sum _{d|T}[{\sigma }_d<=a]{\sigma }_d\mu (\frac{T}{d})\end{array}$$

维护一下后半部分即可，但是要注意并不是所有的 $d$ 都会对这部分产生贡献。所以我们可以考虑离线，按照 $a$ 的值进行排序，使其递增，这样就可以用树状数组不断地进行修改，从而维护出后半部分。

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=1e5+10,M=1e4+10,p=(1ll<<31);
struct stu{
	int n,m,a,id;
	bool operator<(const stu stu1)const{
		return a<stu1.a;
	}
}s[N];
struct sta{
	int m,id;
	bool operator<(const sta stu1)const{
		return m<stu1.m;
	}
}di[N];
inline int lb(int x){return x&-x;}
int miu[N],t,prime[N],dd[N],ans[N],c[N];
bool vis[N];
void init(){
	vis[1]=1;miu[1]=1;
	for(int i=2;i<N;i++){
		if(!vis[i]){
			prime[++prime[yhl]]=i;
			miu[i]=-1;
		}
		for(int j=1;j<=prime[yhl]&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(i%prime[j])miu[i*prime[j]]=-miu[i];
			else break;
		}
	}
	for(int i=1;i<N;i++)
		for(int j=i;j<N;j+=i)dd[j]+=i;
	for(int i=1;i<N;i++)di[i].id=i,di[i].m=dd[i];
	sort(di+1,di+N);
}
void add(int x,int val){
	while(x<N){
		c[x]=(c[x]+val)%p;
		x+=lb(x);
	}
}
int get_sum(int x){
	int ans=yhl;
	while(x){
		ans=(ans+c[x])%p;
		x-=lb(x);
	}
	return ans;
}
void Add(int x){
	for(int i=1;i*x<N;i++)add(i*x,miu[i]*dd[x]%p);
}
int ghfz(int n,int m){
	if(n>m)swap(n,m);
	int ans=yhl,r=yhl;
	for(int i=1;i<=n;i=r+1){
		r=min(n/(n/i),m/(m/i));
		ans=(1ll*(1ll*(n/i)*(m/i)%p)*(get_sum(r)-get_sum(i-1)+p)%p+ans)%p;
	}
	return ans;
}
signed main(){
	scanf("%lld",&t);
	for(int i=1;i<=t;i++)scanf("%lld%lld%lld",&s[i].n,&s[i].m,&s[i].a),s[i].id=i;
	sort(s+1,s+1+t);
	int r=yhl;init();
	for(int i=1;i<=t;i++){
		while(di[r+1].m<=s[i].a&&r<N-1)++r,Add(di[r].id);
		ans[s[i].id]=ghfz(s[i].n,s[i].m)%p;
	}
	for(int i=1;i<=t;i++)printf("%lld\n",ans[i]);
	return yhl;
}

P6156 简单题

即：

$$\begin{array}{rl}ans& =\sum _{i=1}^n\sum _{j=1}^n{\mu }^2(gcd(i,j))gcd(i,j)(i+j{)}^k\\ & =\sum _{d=1}^n{\mu }^2(d)d^{k+1}\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }(i+j{)}^k[gcd(i,j)=1]\\ & =\sum _{d=1}^n{\mu }^2(d)d^{k+1}\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{d}\rfloor }(i+j{)}^k\sum _{t|gcd(i,j)}\mu (t)\\ & =\sum _{d=1}^n{\mu }^2(d)d^{k+1}\sum _{t=1}^{\lfloor \frac{n}{d}\rfloor }\mu (t)t^k\sum _{i=1}^{\lfloor \frac{n}{td}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{td}\rfloor }(i+j{)}^k\end{array}$$

到这里我觉得就可以做了。前两个部分分别 $O(n)$ 预处理，注意充分利用我们所筛出来的 $i^k$ 。
令$f(x)=x^k$ ，观察可知，它是一个积性函数，直接套在线性筛里可以做到近似 $O(n)$ 的复杂度。还可以用费马小定理将 $k$ 缩小到 $[0,p-1]$ 的范围。
在计算柿子时，我们套两个数论分块，这样就是近似 $O(n)$ 的时间复杂度（我认为🌚）。像我这种喜欢 define int long long 的，可以考虑重复利用数组，比如线性筛之后你的 $prime$ 数组就没用了，就可以拿来预处理东西。
upd:学长好像已经算出我的时间复杂度为 $O(n^{\frac{3}{4}})$ 了🌚。

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
#define int long long
#define it int
const int N=5e6+100,p=998244353;
int n,k,g[N],s[N<<1],op[N];
it prime[4900000],miu[N];
bool vis[N<<1];
inline int qum(register int a,register int b){
	register int ans=1;
	while(b){
		(b&1)&&(ans=1ll*ans*a%p);
		a=1ll*a*a%p;
		b>>=1;
	}
	return ans;
}
inline void init(){
	vis[1]=1;miu[1]=1;s[1]=1;
	for(register int i=2;i<=n*2;i++){
		if(!vis[i]){
			if(i<=n)miu[i]=-1;
			prime[++prime[yhl]]=i,s[i]=qum(i,k);
		}
		for(register int j=1;j<=prime[yhl]&&i*prime[j]<=2*n;j++){
			vis[i*prime[j]]=1;s[i*prime[j]]=1ll*s[i]*s[prime[j]]%p;
			if(i%prime[j])if(i*prime[j]<=n)miu[i*prime[j]]=-miu[i];
			else break;
		}
	}
	for( int i=1;i<=n*2;i++)s[i]=(s[i]+s[i-1])%p;
	for(register int i=1;i<=n;i++){
		g[i]=miu[i]?(g[i-1]+1ll*(s[i]-s[i-1]+p)*i%p)%p:g[i-1];
	}
	for(register int i=1;i<=n;i++)op[i]=(op[i-1]+(s[i*2]+s[i*2-1]-2*s[i])+p)%p;
	register int last=yhl,llst=yhl;
	for(register int i=1;i<=n;i++){
		llst=last;
		last=s[i];
		s[i]=(s[i-1]+miu[i]*(s[i]-llst))%p;
		if(s[i]<0)s[i]+=p;
	}
}
inline int Ghfz(register int n){
	register int ans=yhl,r=yhl;
	for(register int i=1;i<=n;i=r+1){
		r=n/(n/i);
		ans+=1ll*(s[r]-s[i-1]+p)*op[n/i]%p;
	}
	return ans%p;
}
inline int ghfz(register int n){
	register int ans=yhl,r=yhl;
	for(register int i=1;i<=n;i=r+1){
		r=n/(n/i);
		register int ret=Ghfz(n/i);
		ans+=1ll*(g[r]-g[i-1]+p)*ret%p;
	}
	return ans%p;
}
signed main(){
	scanf("%lld%lld",&n,&k);k%=(p-1);
	init();
	printf("%lld",ghfz(n));
	return yhl;
}

正解改出来了，这是另一道题令 $T=dt$ ，则：

$$\begin{array}{rl}ans& =\sum _{T=1}^n\sum _{d|T}{\mu }^2(d)d^{k+1}\mu (\frac{T}{d})(\frac{T}{d}{)}^k\sum _{i=1}^{\lfloor \frac{n}{T}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{T}\rfloor }(i+j{)}^k\\ & =\sum _{T=1}^n{T}^k\sum _{i=1}^{\lfloor \frac{n}{T}\rfloor }\sum _{j=1}^{\lfloor \frac{n}{T}\rfloor }(i+j{)}^k\sum _{d|T}{\mu }^2(d)\mu (\frac{T}{d})d\end{array}$$

中间那部分还是那么维护，后面的维护《见仁见智了》。

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
const int N=1e7+10;
#define int unsigned int
int prime[N<<1],f[N],n,t,k,g[N<<1],op[N];
bool vis[N<<1];
int qum(int a,int b){
	int ans=1;
	while(b){
		if(b&1)ans=ans*a;
		a=a*a;
		b>>=1;
	}
	return ans;
}
void init(){
	vis[1]=1;g[1]=1;f[1]=1;
	for(int i=2;i<=2*n;i++){
		if(!vis[i]){
			prime[++prime[yhl]]=i;g[i]=qum(i,k);
			if(i<=n)f[i]=i-1;
		}
		for(int j=1;j<=prime[yhl]&&i*prime[j]<=2*n;j++){
			vis[i*prime[j]]=1;g[i*prime[j]]=g[i]*g[prime[j]];
			if(i%prime[j]==yhl){
				if(i*prime[j]<=n)
				if((i/prime[j])%prime[j])f[i*prime[j]]=-prime[j]*f[i/prime[j]];
				break;
			}
			if(i*prime[j]<=n)f[i*prime[j]]=f[i]*(prime[j]-1);
		}
	}
	for(int i=1;i<=n;i++)f[i]=(f[i-1]+f[i]*g[i]);
	for(int i=1;i<=n*2;i++)g[i]+=g[i-1];
	for(int i=1;i<=n;i++)
		op[i]=(op[i-1]+2*(g[i*2]-g[i])-(g[i*2]-g[i*2-1]));
}
int read(){
	int ans=0;char ch=getchar();
	while(ch<'0'||ch>'9')ch=getchar();
	while(ch>='0'&&ch<='9')ans=(ans<<3)+(ans<<1)+(ch^48),ch=getchar();
	return ans;
}
void w(int x){
	if(x>9)w(x/10);
	putchar(x%10+48);
}
signed main(){
	t=read(),n=read(),k=read();
	init();
	while(t--){
		n=read();
		int ans=yhl,r;
		for(int i=1;i<=n;i=r+1){
			r=n/(n/i);
			ans=(ans+1ll*op[n/i]*(f[r]-f[i-1]));
		}
		w(ans);puts("");
	}
	return yhl;
}

Crash的数字表格 / JZPTAB

即：

$$\begin{array}{rl}ans& =\sum _{i=1}^n\sum _{j=1}^{m}lcm(i,j)\\ & =\sum _{i=1}^n\sum _{j=1}^m\frac{ij}{gcd(i,j)}\\ & =\sum _{d=1}^n\sum _{i=1}^n\sum _{j=1}^m\frac{ij}{d}[gcd(i,j)=d]\\ & =\sum _{d=1}^n\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{d}\rfloor }ijd\sum _{k|gcd(i,j)}\mu (k)\\ & =\sum _{d=1}^{n}d\sum _{k=1}^{\lfloor \frac{n}{d}\rfloor }\mu (k)k^2\sum _{i=1}^{\lfloor \frac{n}{dk}\rfloor }\sum _{j=1}^{\lfloor \frac{m}{dk}\rfloor }ij\end{array}$$

令 $T=kd$ ，$f(x)=\sum _{i=1}^{x}i=\frac{x\times (x+1)}{2}$ 则：

$$\begin{array}{rl}ans& =\sum _{T=1}^{n}f(\lfloor \frac{n}{T}\rfloor )f(\lfloor \frac{m}{T}\rfloor )\sum _{d|T}\mu (d)d^2\frac{T}{d}\\ & =\sum _{T=1}^{n}f(\lfloor \frac{n}{T}\rfloor )f(\lfloor \frac{m}{T}\rfloor )\sum _{d|T}\mu (d)dT\end{array}$$

后面单独预处理即可。

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=1e7+10,p=1e8+9;
int prime[N],miu[N],n,t,m,f[N],calc[N];
bool vis[N];
inline int read(){
	register int ans=0;char ch=getchar();
	while(ch<'0'||ch>'9')ch=getchar();
	while(ch>='0'&&ch<='9')ans=(ans<<3)+(ans<<1)+(ch^48),ch=getchar();
	return ans;
}
inline void w(register int x){
	if(x>9)w(x/10);
	putchar(x%10+48);
}
inline void init(){
	vis[1]=1;miu[1]=1;
	for(register int i=2;i<N;i++){
		if(!vis[i]){
			prime[++prime[yhl]]=i;
			miu[i]=-1;
		}
		for(register int j=1;j<=prime[yhl]&&i*prime[j]<N;j++){
			vis[i*prime[j]]=1;
			if(!(i%prime[j]))break;
			miu[i*prime[j]]=-miu[i];
		}
	}
	for(register int d=1;d<N;d++){
		for(register int j=1;j*d<N;j++){
			if(miu[d])f[d*j]=(f[d*j]+1ll*(1ll*((1ll*miu[d]*d%p)+p)*d%p)*j%p)%p;
			
		}
	}
	for(register int i=1;i<N;i++)f[i]=(f[i]+f[i-1])%p;
	for(register int i=1;i<N;i++)calc[i]=(calc[i-1]+i)%p;
	//这个其实没必要预处理。 
}
signed main(){
	init();
	scanf("%lld",&t);
	while(t--){
		scanf("%lld%lld",&n,&m);
		if(n>m)swap(n,m);
		int ans=yhl,r=yhl;
		for(register int i=1;i<=n;i=r+1){
			r=min(m/(m/i),n/(n/i));
			ans=(ans+1ll*(1ll*(f[r]-f[i-1]+p)*calc[n/i]%p)*calc[m/i]%p)%p;
		}
		w(ans);puts("");
	}
	return yhl;
}

看到这你就会发现，虽然这叫莫比乌斯反演做题记，但是我并没用反演来解决问题，而更多的是用 $\mu$ 函数的性质。

《请无视上面那句话》，终于是用上开头的反演柿子了。

镜中的野兽

首先可以发现：若干个数的最小公倍数一定是最大公因数的倍数，那么我们可以在 $m$ 的因数里枚举这个最大公因数，记这个最大公因数是 $d$ 。

那么，问题就变成了：求元素个数为 $n$ 、这些元素的 $gcd$ 与 $lcm$ 分别为 $d$ 和 $m-d$ 的集合个数。

《看起来和之前那个问题一样》，但是这个时候我们发现对于每个合法的集合，我们把其中每个元素都除以一个 $d$ ，问题就可以转化为：

求大小为 $n$ ，所有元素的 $gcd$ 和 $lcm$ 分别是 $1$ 和 $\frac{m}{d}-1$ 的集合个数。

就可以将其抽象成数学柿子：

$$f(x)=\sum \sum \sum \cdots \sum [gcd=1][lcm=x]$$

那么

$$ans=\sum _{d|m}f(\frac{m}{d}-1)$$

现在再定义两个函数 $g、h$ 。

$$g(x)=\sum \sum \sum \cdots \sum [gcd|x][lcm=x]$$

$$h(x)=\sum \sum \sum \cdots \sum [gcd|x][lcm|x]$$

思考 $h$ 函数的意义，其实就是在 $x$ 的因数里随便选择 $n$ 个，那答案就是 $C_{d(x)}^n$ ，这里 $d(x)$ 是 $x$ 的约数个数。
然后还有两个柿子：

$$g(x)=\sum _{d|x}f(d)$$

$$h(x)=\sum _{d|x}g(d)$$

证明

然后就有：

$$f(x)=\sum _{d|x}\mu (\frac{x}{d})g(d)$$

$$g(x)=\sum _{d|x}\mu (\frac{x}{d})h(d)$$

整理可得：

$$f(x)=\sum _{d|x}\mu (\frac{x}{d})\sum _{k|d}\mu (\frac{d}{k})C_{d(k)}^n$$

所以：

$$ans=\sum _{x|m}\sum _{d|(\frac{m}{x}-1)}\mu (\frac{\frac{m}{x}-1}{d})\sum _{k|d}\mu (\frac{d}{k})C_{d(k)}^n$$

点击查看代码

#define yhl 0
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N=1e5+10;
int n,m,p,ans;
int ny[N],fact[N];
unordered_map<int,int>d,miu;
inline int C(int n,int m){if(n<m)return 0;return fact[n]*ny[m]%p*ny[n-m]%p;}
int qum(int a,int b){
    int ans=1;
    while(b){
        (b&1)&&(ans=ans*a%p);
        a=a*a%p;
        b>>=1;
    }
    return ans;
}
void init(){
    fact[yhl]=1;ny[yhl]=1;
    for(int i=1;i<N;i++){
        fact[i]=fact[i-1]*i%p;
    }
    ny[N-1]=qum(fact[N-1],p-2);
    for(int i=N-2;i;i--){
        ny[i]=ny[i+1]*(i+1)%p;
    }
}
int get_d(int x){
    if(d[x])return d[x];
    int tot=yhl;
    for(int i=1;i*i<=x;i++){
        if(!(x%i))
        tot+=1+(i*i!=x);
    }
    return d[x]=tot; 
}
int get_miu(int x){
    if(miu[x])return miu[x]-2;
    int tot=(x!=1),op=x;
    for(int i=2;i*i<=x;i++){
        if(!(x%i)){
            tot++;
            x/=i;
            if(x%i==yhl){
                miu[op]=2;
                return yhl;
            }
        }
    }
    miu[op]=(tot&1)?-1:1;
    miu[op]+=2;
    return (tot&1)?-1:1;
}
int get_g(int x){
    int tot=yhl;
    for(int i=1;i*i<=x;i++){
        if(x%i==yhl){
            int j=x/i;
            if(get_miu(j))tot=(tot+get_miu(j)*C(get_d(i),n)+p)%p;
            if(j!=i&&(get_miu(i)))tot=(tot+get_miu(i)*C(get_d(j),n)+p)%p;
        }
    }
    return tot;
}
int get_f(int x){
    int tot=yhl;
    for(int i=1;i*i<=x;i++){
        if(x%i==yhl){
            int j=x/i;
            if(get_miu(j))tot=(tot+get_miu(j)*get_g(i)+p)%p;
            if(j!=i&&get_miu(i))tot=(tot+get_miu(i)*get_g(j)+p)%p;
        }
    }
    return tot;
}
signed main(){
    // freopen("1.in","r",stdin);
    scanf("%lld%lld%lld",&n,&m,&p);
    init();
    for(int i=1;i*i<=m;i++){
        if(m%i==yhl){
            int j=m/i;
            ans=(ans+get_f(j-1))%p;
            if(j!=i)ans=(ans+get_f(i-1))%p;
        }
    }
    printf("%lld",ans);
    return yhl;
}

《感谢ppllxx_9G的指正🌚》

$ToBeContinued\dots \dots$