luoguP3312 [SDOI2014]数表

题意

默认(nleqslant m)。

设(f(i))表示(i)的约数和，因为是积性函数，可以用线性筛求。

先不考虑(a)的限制，我们推下式子：
(sumlimits_{i=1}^{n}sumlimits_{j=1}^{m}f(gcd(i,j)))
枚举(gcd(i,j))
(sumlimits_{d=1}^{n}f(d)sumlimits_{i=1}^{n}sumlimits_{j=1}^{m}[gcd(i,j)=d])
之后就莫反：
(sumlimits_{d=1}^{n}f(d)sumlimits_{i=1}^{frac{n}{d}}sumlimits_{j=1}^{frac{m}{d}}[gcd(i,j)=1])
(sumlimits_{d=1}^{n}f(d)sumlimits_{i=1}^{frac{n}{d}}sumlimits_{j=1}^{frac{m}{d}}sumlimits_{x|gcd(i,j)}mu(x))
(sumlimits_{d=1}^{n}f(d)sumlimits_{x=1}^{frac{n}{d}}mu(x)sumlimits_{i=1}^{frac{n}{d}}sumlimits_{j=1}^{frac{m}{d}}[x|gcd(i,j)])
(sumlimits_{d=1}^{n}f(d)sumlimits_{x=1}^{frac{n}{d}}mu(x)*frac{n}{d*x}*frac{m}{d*x})
设(T=d*x)：
(sumlimits_{T=1}^{n}frac{n}{T}*frac{m}{T}sumlimits_{d|T}f(d)*mu(frac{T}{d}))
除法分块加预处理(sumlimits_{d|T}f(d)*mu(frac{T}{d}))即可。

对于(a)的限制，我们离线按(a)排序，同时用树状数组维护即可。

code:

#include<bits/stdc++.h>
using namespace std;
#define pli pair<ll,int>
#define mkp make_pair
#define fir first
#define sec second
typedef long long ll;
const int maxn=1e5+10;
const int maxq=2*1e4+10;
const ll mod=2147483648;
int Q;
int mu[maxn];
ll g[maxn],ans[maxq];
bool vis[maxn];
pli f[maxn];
vector<int>prime;
struct Query{int n,m,lim,id;}qr[maxq];
struct Tree_arry
{
	#define lowbit(x) (x&-x)
	ll a[maxn];
	inline void add(int x,ll k){for(int i=x;i<=100000;i+=lowbit(i))a[i]=(a[i]+k)%mod;}
	inline ll query(int x){ll res=0;for(int i=x;i;i-=lowbit(i))res=(res+a[i])%mod;return res;}
}tr;
inline bool cmp(Query x,Query y){return x.lim<y.lim;}
inline void pre_work(int n)
{
	mu[1]=1;f[1]=mkp(1,1);vis[1]=1;
	for(int i=2;i<=n;i++)
	{
		if(!vis[i])prime.push_back(i),mu[i]=-1,f[i]=mkp(i+1,i),g[i]=i+1;
		for(unsigned int j=0;j<prime.size()&&i*prime[j]<=n;j++)
		{
			vis[i*prime[j]]=1;
			if(i%prime[j]==0)
			{
				mu[i*prime[j]]=0;
				g[i*prime[j]]=g[i]*prime[j]+1;
				f[i*prime[j]]=mkp(f[i].fir/g[i]*g[i*prime[j]],i*prime[j]);
				break;
			}
			mu[i*prime[j]]=-mu[i];
			f[i*prime[j]]=mkp(f[i].fir*f[prime[j]].fir,i*prime[j]);
			g[i*prime[j]]=prime[j]+1;
		}
	}
}
inline void work(int x)
{
	for(int i=1;i*f[x].sec<=100000;i++)tr.add(i*f[x].sec,(f[x].fir*mu[i]%mod+mod)%mod);
}
inline ll calc(int n,int m)
{
	ll res=0;
	if(n>m)swap(n,m);
	for(int l=1,r;l<=n;l=r+1)
	{
		r=min(n/(n/l),m/(m/l));
		res=(res+((tr.query(r)-tr.query(l-1))%mod+mod)%mod*(n/l)%mod*(m/l)%mod)%mod;
	}
	return res;
}
int main()
{
	pre_work(100000);
	sort(f+1,f+100000+1);
	scanf("%d",&Q);
	for(int i=1;i<=Q;i++)scanf("%d%d%d",&qr[i].n,&qr[i].m,&qr[i].lim),qr[i].id=i;
	sort(qr+1,qr+Q+1,cmp);
	for(int i=1,j=1;i<=Q;i++)
	{
		while(f[j].fir<=qr[i].lim&&j<=100000)work(j),j++;
		ans[qr[i].id]=calc(qr[i].n,qr[i].m);
	}
	for(int i=1;i<=Q;i++)printf("%lld
",ans[i]);
	return 0;
}