BZOJ 3529 [Sdoi2014]数表 【莫比乌斯反演】

Description

    有一张N×m的数表，其第i行第j列（1 < =i < =n，1 < =j < =m）的数值为
能同时整除i和j的所有自然数之和。给定a，计算数表中不大于a的数之和。

HINT

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

Solution

先假设a的限制不存在

正面把答案强行写出来是这样的

其中F(i)为i的约数和，可以线性筛处理

然后慢慢化简其中的各个部分

令g(i)为1<=x<=n,1<=y<=m,gcd(x,y)=i的数对(x,y)的个数

则有

用莫比乌斯反演代换

此时只需要一个 的前缀和就可以继续愉快地分块了

然而此时还有a的限制，怎么继续做前缀和？

如果我们可以每次求前缀和的时候都只选择F( i )比a小的就好了

此时想到其实还有一种前缀和的方法树状数组，并且可以以加入的顺序控制前缀和

于是自然地离线询问，按a值排序， 按F( i )排序后依次加入树状数组

取模用了PoPoQQQ的自然溢出

#include<map>
#include<cmath>
#include<ctime>
#include<queue>
#include<stack>
#include<cstdio>
#include<climits>
#include<iomanip>
#include<cstring>
#include<cstdlib>
#include<iostream>
#include<algorithm>

#define maxp 100000
#define maxn 100000+5
#define set(a,b) memset(a,(b),sizeof(a))
#define fr(i,a,b) for(ll i=(a),_end_=(b);i<=_end_;i++)
#define rf(i,b,a) for(ll i=(a),_end_=(b);i>=_end_;i--)
#define fe(i,a,b) for(int i=first[(b)],_end_=(a);i!=_end_;i=s[i].next)
#define fec(i,a,b) for(int &i=cur[(b)],_end_=(a);i!=_end_;i=s[i].next)

using namespace std;

typedef long long ll;

struct sigma{
int number,F;
}s[maxn];

struct query{
int n,m,a;
int number;
}q[maxn];

int bit[maxn];
int prime[maxn],pri[maxn],miu[maxn];
int ans[maxn],tot=0,st=1;
int T;

void read()
{
#ifndef ONLINE_JUDGE
freopen("3529.in","r",stdin);
freopen("3529.out","w",stdout);
#endif
//cin >> T ;
scanf("%d",&T);
}

void write()
{
fr(i,1,T)
//cout << (ans[i]&0x7fffffff) << endl ;
printf("%d\n",ans[i]&0x7fffffff);
}

void get()
{
miu[1]=1;
fr(i,2,maxp){
if( !prime[i] ) pri[++tot]=i,miu[i]=-1;
int j=1;
while( j<=tot && pri[j]*i<=maxp ){
prime[i*pri[j]]=1;
if( i%pri[j]==0 ){
miu[i*pri[j]]=0;
break;
}
miu[i*pri[j]]=-miu[i];
j++;
}
}
fr(i,1,maxp){
fr(j,1,maxp/i)
s[i*j].F+=i;
s[i].number=i;
}
}

bool comp(sigma a,sigma b)
{
return a.F<b.F;
}

bool cmp(query a,query b)
{
return a.a<b.a;
}

void add(int x,int w)
{
while( x<=maxp ){
bit[x]+=w;
x+=x&-x;
}
}

int sum(int x)
{
int res=0;
while( x>0 ){
res+=bit[x];
x-=x&-x;
}
return res;
}

int calc(int x,int y)
{
if( x>y ) swap(x,y);
int res=0,pos;
int i=1;
while( i<=x ){
pos=min((x/(x/i)),(y/(y/i)));
res+=(sum(pos)-sum(i-1))*(x/i)*(y/i);
i=pos+1;
}
return res;
}

void work()
{
fr(i,1,T){
//cin >> q[i].n >> q[i].m >> q[i].a ;
scanf("%d%d%d",&q[i].n,&q[i].m,&q[i].a);
q[i].number=i;
}
get();
sort(q+1,q+T+1,cmp);
sort(s+1,s+maxp+1,comp);
fr(i,1,T){
while( st<=maxp && s[st].F<=q[i].a ){
fr(j,1,maxp/s[st].number)
add(s[st].number*j,s[st].F*miu[j]);
st++;
}
ans[q[i].number]=calc(q[i].n,q[i].m);
}
}

int main()
{
read();
work();
write();
return 0;
}