bzoj 3529: [Sdoi2014]数表 （反演）

题目描述

传送门

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

题解

先把题目化成式子。

$$\sum\limits_{i=1}^n\sum\limits_{j=1}^m h(gcd(i,j))[h(gcd(i,j))<=a]$$ 其中 $h(x)$表示x的约数和，可以用线性筛求解，也可以用 $O(nlogn)$的时间求解。
我们先考虑没有a的限制怎么做。
$$\sum\limits_{d=1}^n \sum\limits_{i=1}^n\sum\limits_{j=1}^m h(d) [gcd(i,j)=d]$$
$$\sum\limits_{d=1}^n h(d)\sum\limits_{i=1}^n \mu(i)*{\lfloor {n\over i*d}\rfloor} *{\lfloor{m\over i*d}\rfloor}$$
$$\sum\limits_{T=1}^n{\lfloor {n\over T}\rfloor} *{\lfloor{m\over T}\rfloor} \sum\limits_{d|n} \mu({T \over d}) *h(d)$$
那么我们可以考虑将所有的询问离线，按照a排序。然后将所有的数按照h的值排序。每次计算某个询问之前，就将h值小于等于ai的数全部加入，然后用树状数组维护 $f(n)=\sum\limits_{d|n} \mu({T \over d}) *h(d)$ 的前缀和。

代码

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
#define N 100000
#define LL long long
#define inf 1000000000
using namespace std;
int pd[N+3],prime[N+3],mu[N+3],tr[N+3],n,m,T,ans[N];
struct data {
    int x,y,val,id;
}a[N],num[N];
int cmp(data a,data b)
{
    return a.val<b.val;
}
void init()
{
    mu[1]=1;
    for (int i=2;i<=N;i++) {
        if (!pd[i]) {
            prime[++prime[0]]=i;
            mu[i]=-1;
        }
        for (int j=1;j<=prime[0];j++){
            if (prime[j]*i>N) break;
            pd[prime[j]*i]=1;
            if (i%prime[j]==0) {
                mu[i*prime[j]]=0;
                break;
            } 
            mu[i*prime[j]]=-mu[i];
        }
    }
    for (int i=1;i<=N;i++){
     num[i].x=i;
     for (int j=i;j<=N;j+=i)
      num[j].val+=i;
    }
    sort(num+1,num+N+1,cmp);
    //cout<<num[1].x<<" "<<num[1].val<<endl;
}
int lowbit(int i)
{
    return i&(-i);
}
int query(int x)
{
    int ans=0;
    for (int i=x;i>=1;i-=lowbit(i))
     ans+=tr[i];
    return ans;
}
void change(int x,int val)
{
    for (int i=x;i<=N;i+=lowbit(i))
     tr[i]+=val;
}
void solve(int i)
{
    int x=num[i].x; int val=num[i].val;
    //cout<<x<<" "<<val<<endl;
    if (val<0) return;
    for (int j=1;j*x<=N;j++)
     change(j*x,mu[j]*val);
}
int main()
{
    freopen("a.in","r",stdin);
    freopen("my.out","w",stdout);
    scanf("%d",&T);
    for (int i=1;i<=T;i++) {
     scanf("%d%d%d",&a[i].x,&a[i].y,&a[i].val),a[i].id=i;
     if (a[i].x>a[i].y) swap(a[i].x,a[i].y);
    }
    sort(a+1,a+T+1,cmp);
    int l=1; init();
    for (int k=1;k<=T;k++) {
        while (num[l].val<=a[k].val&&l<=N) 
         solve(l++);
        int sum=0; n=a[k].x; m=a[k].y;
        for (int i=1,j;i<=n;i=j+1) {
            j=min(n/(n/i),m/(m/i)); 
            sum=sum+(n/i)*(m/i)*(query(j)-query(i-1));
        }
        ans[a[k].id]=sum;
    }
    LL p=(LL)(1<<30)*2;
    for (int i=1;i<=T;i++) {
        LL t=ans[i];
        printf("%lld\n",(t%p+p)%p);
    }
}