【HDU4676】Sum Of Gcd-莫队算法+欧拉函数

测试地址：Sum Of Gcd
题目大意：给定一个1$1$~n$n$的全排列A$A$，若干个询问，每次询问给出一个区间[l,r]$[l,r]$，要求得出∑l≤i<j≤rgcd(Ai,Aj)$\sum _{l\le i<j\le r}gcd(A_i,A_j)$的值。
做法：本题需要用到莫队算法+欧拉函数。
我们不好直接统计各个数对的最大公约数出现的次数，但是我们可以比较简单地统计各个数对的公约数出现的次数。由于数对的公约数都是数对最大公约数的约数（有点绕……），我们很自然地思考一个数和它的约数之间有什么容易统计的关系。
我们可以用一个结论将一个数字和它的约数建立统计关系：
n=∑d|nφ(d)$n=\sum _{d|n}\phi (d)$
那么有：
∑l≤i<j≤rgcd(Ai,Aj)$\sum _{l\le i<j\le r}gcd(A_i,A_j)$
=∑l≤i<j≤r∑d|gcd(Ai,Aj)φ(d)$=\sum _{l\le i<j\le r}\sum _{d|gcd(A_i,A_j)}\phi (d)$
=∑nd=1φ(d)∑l≤i<j≤r[d|Ai][d|Aj]$=\sum _{d=1}^n\phi (d)\sum _{l\le i<j\le r}[d|A_i][d|A_j]$
显然式子中的∑l≤i<j≤r[d|Ai][d|Aj]$\sum _{l\le i<j\le r}[d|A_i][d|A_j]$就是指d$d$作为区间内数对的公约数出现的次数。而统计每个数作为公约数出现的次数的方法就很简单了，你只需要统计每个数i$i$是区间内多少个数的约数，令这个数值为cnti$cn{t}_i$，那么它作为公约数出现的次数就是cnti(cnti−1)2$\frac{cn{t}_i(cn{t}_i-1)}{2}$，这样就将问题转化为了对区间内所有数的约数出现的次数的统计。
接下来我们就可以用莫队算法解决了。我们预先处理出所有数的所有约数以及欧拉函数，然后每新增或删除一个元素，对它的所有约数更新cntd$cn{t}_d$。很显然，对于一个约数d$d$，在新增一个元素时，答案增加(cntd−1)φ(d)$(cn{t}_d-1)\phi (d)$，在删除一个元素时，答案减少cntdφ(d)$cn{t}_d\phi (d)$，注意这里的cntd$cn{t}_d$是更新后的值。每新增或删除一个元素的复杂度大概在O(logn)$O(\log n)$左右，所以总的时间复杂度是O(nn−−√logn)$O(n\sqrt{n}\log n)$，可以通过此题。
（最近在学习分块，据说莫队算法是分块的一种，所以也做一道题吧……）
以下是本人代码：

#include <bits/stdc++.h>
#define ll long long
using namespace std;
int T,n,m,blocklen,block[20010],a[20010],fac[20010][210];
ll phi[20010],p[20010],cnt[20010],ans[20010],sum;
bool prime[20010]={0};
struct Query
{
    int id,l,r;
}q[20010];

void calc_phi()
{
    prime[1]=1;p[0]=0;
    for(int i=1;i<=20000;i++) phi[i]=1;
    for(int i=2;i<=20000;i++)
    {
        if (!prime[i]) phi[i]=i-1,p[++p[0]]=i;
        for(int j=1;j<=p[0]&&i*p[j]<=20000;j++)
        {
            prime[i*p[j]]=1;
            if (i%p[j]==0) {phi[i*p[j]]=phi[i]*p[j];break;}
            phi[i*p[j]]=phi[i]*(p[j]-1);
        }
    }
}

bool cmp(Query a,Query b)
{
    if (block[a.l]!=block[b.l]) return block[a.l]<block[b.l];
    else return a.r<b.r;
}

void expand(int x,ll add)
{
    for(int i=1;i<=fac[a[x]][0];i++)
    {
        int j=fac[a[x]][i];
        if (add==1) sum+=cnt[j]*phi[j];
        else sum-=(cnt[j]-1)*phi[j];
        cnt[j]+=add;
    }
}

void Mo()
{
    int l=1,r=0;
    memset(cnt,0,sizeof(cnt));
    sum=0;
    for(int i=1;i<=m;i++)
    {
        while (q[i].l<l) expand(--l,1);
        while (q[i].r>r) expand(++r,1);
        while (q[i].l>l) expand(l++,-1);
        while (q[i].r<r) expand(r--,-1);
        ans[q[i].id]=sum;
    }
}

int main()
{
    scanf("%d",&T);
    calc_phi();
    for(int i=1;i<=20000;i++)
    {
        fac[i][0]=0;
        for(int j=1;j*j<=i;j++)
            if (i%j==0)
            {
                fac[i][++fac[i][0]]=j;
                if (j*j!=i) fac[i][++fac[i][0]]=i/j;
            }
    }

    int t=1;
    while(T--)
    {
        scanf("%d",&n);
        for(int i=1;i<=n;i++)
            scanf("%d",&a[i]);
        scanf("%d",&m);
        blocklen=(int)sqrt(n);
        for(int i=1;i<=n;i++)
            block[i]=i/blocklen;
        for(int i=1;i<=m;i++)
        {
            scanf("%d%d",&q[i].l,&q[i].r);
            q[i].id=i;
        }

        sort(q+1,q+m+1,cmp);
        Mo();

        printf("Case #%d:\n",t);
        for(int i=1;i<=m;i++)
            printf("%lld\n",ans[i]); 
        t++;
    }

    return 0;
}