bzoj3529: [Sdoi2014]数表 莫比乌斯反演

题意:求\(\sum_{i=1}^n\sum_{j=1}^nf(gcd(i,j))(gcd(i,j)<=a),f(x)是x的因子和函数\)
先考虑没有限制的情况,考虑枚举gcd为x,那么有\(\sum_{x=1}^{min(n,m)}f(x)\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)==x]\)
可以发现右侧就是最裸的莫比乌斯反演,那么\(\sum_{x=1}^{min(n,m)}f(x)\sum_{d=1}^{min(\lfloor \frac{n}{x} \rfloor,\lfloor \frac{m}{x} \rfloor)}\mu(d)*{\lfloor \frac{n}{x*d} \rfloor}*{\lfloor \frac{m}{x*d} \rfloor}\)
考虑枚举q=x*d,那么\(\sum_{q=1}^{min(n,m)}{\lfloor \frac{n}{q} \rfloor}*{\lfloor \frac{m}{q} \rfloor}\sum_{x|q}f(x)*\mu(\frac{q}{x})\)
发现后面是个狄利克雷卷积\(g=f*\mu\),f(x)是积性函数可以线性预处理,对于a的限制我们可以离线询问按a排序,然后维护一个树状数组每次扫到询问就把g更新到树状数组中,然后整除分块更新答案
取模可以直接爆int,然后&0x7fffffff,即(1<<32)-1

/**************************************************************
    Problem: 3529
    User: walfy
    Language: C++
    Result: Accepted
    Time:3404 ms
    Memory:6076 kb
****************************************************************/
 
//#pragma GCC optimize(2)
//#pragma GCC optimize(3)
//#pragma GCC optimize(4)
//#pragma GCC optimize("unroll-loops")
//#pragma comment(linker, "/stack:200000000")
//#pragma GCC optimize("Ofast,no-stack-protector")
//#pragma GCC target("sse,sse2,sse3,ssse3,sse4,popcnt,abm,mmx,avx,tune=native")
#include<bits/stdc++.h>
#define fi first
#define se second
#define db double
#define mp make_pair
#define pb push_back
#define pi acos(-1.0)
#define ll long long
#define vi vector<int>
#define mod 1000000007
#define ld long double
#define C 0.5772156649
#define ls l,m,rt<<1
#define rs m+1,r,rt<<1|1
#define pll pair<ll,ll>
#define pil pair<int,ll>
#define pli pair<ll,int>
#define pii pair<int,int>
//#define cd complex<double>
#define ull unsigned long long
#define base 1000000000000000000
#define Max(a,b) ((a)>(b)?(a):(b))
#define Min(a,b) ((a)<(b)?(a):(b))
#define fin freopen("a.txt","r",stdin)
#define fout freopen("a.txt","w",stdout)
#define fio ios::sync_with_stdio(false);cin.tie(0)
template<typename T>
inline T const& MAX(T const &a,T const &b){return a>b?a:b;}
template<typename T>
inline T const& MIN(T const &a,T const &b){return a<b?a:b;}
inline void add(ll &a,ll b){a+=b;if(a>=mod)a-=mod;}
inline void sub(ll &a,ll b){a-=b;if(a<0)a+=mod;}
inline ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}
inline ll qp(ll a,ll b){ll ans=1;while(b){if(b&1)ans=ans*a%mod;a=a*a%mod,b>>=1;}return ans;}
inline ll qp(ll a,ll b,ll c){ll ans=1;while(b){if(b&1)ans=ans*a%c;a=a*a%c,b>>=1;}return ans;}
 
using namespace std;
 
const double eps=1e-8;
const ll INF=0x3f3f3f3f3f3f3f3f;
const int N=100000+10,maxn=400000+10,inf=0x3f3f3f3f;
 
int prime[N],cnt,sum[N],mi[N],mu[N];
bool mark[N];
void init()
{
    sum[1]=mu[1]=1;
    for(int i=2;i<N;i++)
    {
        if(!mark[i])prime[++cnt]=i,mu[i]=-1,sum[i]=i+1,mi[i]=i+1;
        for(int j=1;j<=cnt&&i*prime[j]<N;j++)
        {
            mark[i*prime[j]]=1;
            if(i%prime[j]==0)
            {
                mi[i*prime[j]]=mi[i]*prime[j]+1;
                sum[i*prime[j]]=sum[i]/mi[i]*mi[i*prime[j]];
                mu[i*prime[j]]=0;
                break;
            }
            sum[i*prime[j]]=sum[i]*(prime[j]+1);
            mu[i*prime[j]]=-mu[i];
            mi[i*prime[j]]=1+prime[j];
        }
    }
}
struct bit{
    int sum[N];
    void update(int i,int v)
    {
        for(;i<N;i+=i&(-i))sum[i]+=v;
    }
    int query(int i)
    {
        int ans=0;
        for(;i;i-=i&(-i))ans+=sum[i];
        return ans;
    }
}b;
int ans[N];
struct node{
    int n,m,a,id;
    bool operator <(const node&rhs)const{
        return a<rhs.a;
    }
}p[N];
struct point{
    int d,id;
    bool operator <(const point &rhs)const{
        return d<rhs.d;
    }
}f[N];
int main()
{
    init();
    int q,ma=0;scanf("%d",&q);
    for(int i=1;i<=q;i++)
    {
        scanf("%d%d%d",&p[i].n,&p[i].m,&p[i].a),p[i].id=i;
        if(p[i].n>p[i].m)swap(p[i].n,p[i].m);
        ma=MAX(ma,p[i].n);
    }
    for(int i=1;i<=ma;i++)f[i].d=sum[i],f[i].id=i;
    sort(p+1,p+1+q);
    sort(f+1,f+1+ma);
    for(int now=0,i=1;i<=q;i++)
    {
        while(now+1<=ma&&f[now+1].d<=p[i].a)
        {
            now++;
            for(int j=f[now].id;j<=ma;j+=f[now].id)
                b.update(j,f[now].d*mu[j/f[now].id]);
        }
        for(int j=1,k;j<=p[i].n;j=k+1)
        {
            k=MIN(p[i].n/(p[i].n/j),p[i].m/(p[i].m/j));
            ans[p[i].id]+=(p[i].n/j)*(p[i].m/j)*(b.query(k)-b.query(j-1));
        }
        ans[p[i].id]&=0x7fffffff;
    }
    for(int i=1;i<=q;i++)
        printf("%d\n",ans[i]);
    return 0;
}
/********************
 
********************/