a" />

【BZOJ3529】數表(SDOI2014)-莫比烏斯反演+樹狀數組


測試地址:數表
做法:本題需要用到莫比烏斯反演+樹狀數組。
首先忽略 a 的限制,我們要求的是:
a n s = i = 1 n j = 1 m D ( gcd ( i , j ) )
其中 D ( i ) i 的約數和。
我們把上式改成枚舉公因數 d ,不妨設 n m ,則:
a n s = d = 1 n D ( d ) i = 1 n j = 1 m [ gcd ( i , j ) = d ]
f ( d ) = i = 1 n j = 1 m [ gcd ( i , j ) = d ] ,我們發現這就是求 gcd 等於 d 的數對個數(這簡直是句廢話)。
又令 F ( d ) = d | i f ( i ) ,顯然 F ( d ) gcd 等於 d 倍數的數對個數,那么也有 F ( d ) = n d m d 。根據莫比烏斯反演定理的第二種形式,有:
f ( d ) = d | i μ ( i d ) F ( i ) = d | i μ ( i d ) n i m i
把這個式子帶進 a n s 那個式子,有:
a n s = d = 1 n D ( d ) d | i μ ( i d ) n i m i
交換 d , i 的位置,有:
a n s = i = 1 n n i m i d | i D ( d ) μ ( i d )
顯然我們如果預處理出 g ( i ) = d | i D ( d ) μ ( i d ) 的前綴和,就可以用數論分塊處理每個詢問了。
那么現在我們考慮 a 的限制,實際上 a 是在限制只有 D ( i ) a D ( i ) g ( i ) 有貢獻,因此我們把所有詢問按 a 從小到大排序,對新產生的貢獻暴力在 g 中單點修改,因為 D ( i ) g ( j ) 產生貢獻當且僅當 i | j ,那么修改次數顯然是 n log n 級別的,然后要在數論分塊中維護 g 的前綴和查詢,這個顯然可以用常數小又好寫的樹狀數組解決。
那么我們就解決了此題,時間復雜度為 O ( Q n log n + n log 2 n ) 。注意到這題模數很特殊,直接用int自然溢出,最后輸出時對 2 31 1 取個按位與就行了。
以下是本人代碼:

#include <bits/stdc++.h>
using namespace std;
int T,p[100010];
int maxn,D[100010],mu[100010],s[100010],ans[20010];
int prime[100010];
bool vis[100010]={0};
struct query
{
    int id,n,m,a;
}q[20010];

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

bool cmpp(int a,int b)
{
    return D[a]<D[b];
}

int lowbit(int x)
{
    return x&(-x);
}

void add(int x,int c)
{
    for(int i=x;i<=maxn;i+=lowbit(i))
        s[i]+=c;
}

int sum(int x)
{
    int ans=0;
    for(int i=x;i;i-=lowbit(i))
        ans+=s[i];
    return ans;
}

void init()
{
    scanf("%d",&T);
    maxn=0;
    for(int i=1;i<=T;i++)
    {
        scanf("%lld%lld%lld",&q[i].n,&q[i].m,&q[i].a);
        if (q[i].n>q[i].m) swap(q[i].n,q[i].m);
        maxn=max(maxn,q[i].n);
        q[i].id=i;
    }
    sort(q+1,q+T+1,cmp);

    for(int i=1;i<=maxn;i++)
        for(int j=1;i*j<=maxn;j++)
            D[i*j]=D[i*j]+i;
    for(int i=1;i<=maxn;i++)
        p[i]=i;
    sort(p+1,p+maxn+1,cmpp);
}

void calc_mu()
{
    mu[1]=1;
    prime[0]=0;
    for(int i=2;i<=maxn;i++)
    {
        if (!vis[i])
        {
            prime[++prime[0]]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=prime[0]&&i*prime[j]<=maxn;j++)
        {
            vis[i*prime[j]]=1;
            if (i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            mu[i*prime[j]]=-mu[i];
        }
    }
}

int main()
{
    init();
    calc_mu();

    int nowp=1;
    for(int i=1;i<=T;i++)
    {
        while(D[p[nowp]]<=q[i].a)
        {
            for(int j=1;j*p[nowp]<=maxn;j++)
                add(j*p[nowp],D[p[nowp]]*mu[j]);
            nowp++;
        }
        int n=q[i].n,m=q[i].m;
        ans[q[i].id]=0;
        for(int j=n;j>=1;j=max(n/(n/j+1),m/(m/j+1)))
        {
            int l=max(n/(n/j+1)+1,m/(m/j+1)+1),r=j;
            ans[q[i].id]=ans[q[i].id]+(sum(r)-sum(l-1))*(n/j)*(m/j);
        }
    }

    for(int i=1;i<=T;i++)
        printf("%d\n",ans[i]&2147483647);

    return 0;
}

注意!

本站转载的文章为个人学习借鉴使用,本站对版权不负任何法律责任。如果侵犯了您的隐私权益,请联系我们删除。



 
粤ICP备14056181号  © 2014-2020 ITdaan.com