[Sdoi2014]数表

题意

3529: [Sdoi2014]数表

Time Limit: 10 Sec Memory Limit: 512 MB
Submit: 3009 Solved: 1559
[Submit][Status][Discuss]

Description

有一张 n×m 的数表，其第 i 行第 j 列（1 <= i <= n, 1 <= j <= m）的数值为

能同时整除 i 和 j 的所有自然数之和。给定 a , 计算数表中不大于 a 的数之和。

Input

输入包含多组数据。

输入的第一行一个整数Q表示测试点内的数据组数

接下来Q行，每行三个整数n，m，a(|a| < =10^9)描述一组数据。

1 < =N．m < =10^5 ， 1 < =Q < =2×10^4

Output

对每组数据，输出一行一个整数，表示答案模2^31的值。

Sample Input

2
4 4 3
10 10 5

Sample Output

20
148

HINT

Source

Round 1 Day 1

[Submit][Status][Discuss]
HOME Back

分析

参照MashiroSky的题解。

题意

多组询问，每次给出({n,m,a})。求$${sum_{i=1}^nsum_{j=1}^m[sum_{d|i,d|j}d<=a]sum_{d|i,d|j}d}$$

Solution

有({a})的限制很不好做，我们不妨先来看看如果没有({a})的限制，是个什么情况。

[egin{aligned} & sum_{i=1}^nsum_{j=1}^msum_{d|i,d|j}d \ =&sum_{i=1}^nsum_{j=1}^msum_{d=1}^{min(n,m)}f(d)[gcd(i,j)=d] end{aligned} ]

这里({f(d)})表示({d})的约数和，由约数和定理：$${d=p_1^{a_1}p_2^{a_2}p_3^{a_3}······p_k^{a_k}}$$

[{f(d)=(p_1^0+p_1^1+······+p_1^{a_1})(p_2^0+p_2^1+······+p_2^{a_2})······(p_k^0+p_k^1+······+p_k^{a_k})} ]

我们发现当({p,q})互质时，({f(pq)=f(p)f(q)})，({f})是个积性函数。也就是说({f})可以线性筛求解，具体求法可以看代码。回到上面的问题，我们继续变化这个式子。

[egin{aligned} &sum_{i=1}^nsum_{j=1}^msum_{d=1}^{min(n,m)}f(d)[gcd(i,j)=d] \ =&sum_{d=1}^{min(n,m)}f(d)sum_{i=1}^{lfloor{n/d} floor}sum_{j=1}^{lfloor{m/d} floor}[gcd(i,j)=1] \ =&sum_{d=1}^{min(n,m)}f(d)sum_{t=1}^{min(lfloor{n/d} floor,lfloor{m/d} floor)}μ(t)lfloorfrac{n}{dt} floorlfloorfrac{m}{dt} floor end{aligned} ]

我们令({Q=dt})，先枚举({Q})，得到：$${sum_{Q=1}^{min(n,m)}lfloorfrac{n}{Q} floorlfloorfrac{m}{Q} floorsum_{d|Q}f(d)μ(frac{Q}{d})}$$

我们对({f})和({μ})函数的狄利克雷卷积({g(Q)=sum_{d|Q}f(d)μ(frac{Q}{d})})预处理前缀和，前面的({lfloorfrac{n}{Q} floorlfloorfrac{m}{Q} floor})分块，那么这个式子的计算就达到了({sqrt{n}})。

考虑如何加入限制({a})。

我们可以离线，将询问按照({a})从小到大排序，({f})函数按照从小到大排序。询问的时候，我们就动态维护一个树状数组({c})。树状数组的每一位({c[x]})存的是卷积({g(x)})中({f(d)})不超过({a})的数之和。数据组数为({T})，询问复杂度({O(Tsqrt{n}logn)})

细节

这题复杂度有点卡，对于取模运算我们可以对({int})自然溢，最后再({&2147483647})就可以了。

代码

#include<bits/stdc++.h>
#define rg register
#define il inline
#define co const
template<class T>il T read(){
    rg T data=0,w=1;
    rg char ch=getchar();
    while(!isdigit(ch)){
        if(ch=='-') w=-1;
        ch=getchar();
    }
    while(isdigit(ch))
        data=data*10+ch-'0',ch=getchar();
    return data*w;
}
template<class T>il T read(rg T&x){
    return x=read<T>();
}
typedef long long ll;

co int maxn=1e5+1,INF=0x7fffffff;
int mu[maxn],ans[maxn],c[maxn],p[maxn],t[maxn],g[maxn];
struct F{int d,num;}f[maxn];
struct Q{int n,m,a,id;}q[maxn];
bool cmpT(co Q&a,co Q&b) {return a.a<b.a;}
bool cmpt(co F&a,co F&b) {return a.d<b.d;}
#define lowbit(x) (x&-x)
int power(int a,int b){
	int res=1;
	while(b){
		if(b&1) res*=a;
		a*=a,b>>=1;
	}
	return res;
}
void add(int x,int val){
	for(int i=x;i<maxn;i+=lowbit(i)) c[i]+=val;
}
int query(int x){
	int s=0;
	for(int i=x;i>=1;i-=lowbit(i)) s+=c[i];
	return s;
}
int main(){
//	freopen(".in","r",stdin),freopen(".out","w",stdout);
	mu[1]=1,f[1].d=f[1].num=1;
	for(int i=2;i<maxn;++i){
		f[i].num=i;
		if(!p[i]) mu[i]=-1,f[i].d=t[i]=1+i,g[i]=1,p[++p[0]]=i;
		for(int j=1;j<=p[0]&&i*p[j]<maxn;++j){
			p[i*p[j]]=1;
			if(i%p[j]==0){
				mu[i*p[j]]=0;
				g[i*p[j]]=g[i]+1;
				t[i*p[j]]=t[i]+power(p[j],g[i]+1);
				f[i*p[j]].d=f[i].d/t[i]*t[i*p[j]];
			}
			else{
				mu[i*p[j]]=-mu[i];
				f[i*p[j]].d=f[i].d*f[p[j]].d;
				g[i*p[j]]=1;t[i*p[j]]=p[j]+1;
			}
		}
	}
	int T=read<int>();
	for(int i=1;i<=T;++i) read(q[i].n),read(q[i].m),read(q[i].a),q[i].id=i;
	std::sort(q+1,q+T+1,cmpT);
	std::sort(f+1,f+maxn,cmpt);
	for(int now=0,i=1;i<=T;++i){
		while(now+1<maxn&&f[now+1].d<=q[i].a){
			++now;
			for(int j=1;j*f[now].num<maxn;++j)
				add(j*f[now].num,mu[j]*f[now].d);
		}
		int n=q[i].n,m=q[i].m;
		if(n>m) std::swap(n,m);
		for(int j=1,k;j<=n;j=k+1){
			k=std::min(n/(n/j),m/(m/j));
			ans[q[i].id]+=(n/j)*(m/j)*(query(k)-query(j-1));
		}
		ans[q[i].id]&=INF;
	}
	for(int i=1;i<=T;++i) printf("%d
",ans[i]);
	return 0;
}