【数学】gcd

上个阶段太颓也没怎么写题解导致直接炸掉，所以这个阶段争取继续写吧（咕咕咕）

题目内容

有 (n) 个正整数 (x_1) ~ (x_n)，初始时状态均为未选。 有 (m) 个操作，每个操作给定一个编号 (i)，将 (x_i) 的选取状态取反。 每次操作后，你需要求出选取的数中有多少个互质的无序数对。

数据范围

(n,mleq 200000,x_ileq 500000)

思路

做数学题总是很需要磨。

首先解释一波题意，选取状态取反的是 (x_i)，互质的无序数对指的也是 (x) 间形成的数对。

首先考虑不修改怎么做，那么问题就是有 (n) 个正整数，求其中互质的数对的个数。

互质即 (gcd(a,b)=1)。那我们可以定义 (f[i]) 表示 (gcd=i) 的数对个数。考虑如何求出 (f[i])。定义 (g[i]) 表示 (gcd) 为 (i) 的倍数的数对个数。那么显然有：

[g[i]=sumlimits_{ivert d}f[d] ]

嗯，这个显然是莫比乌斯反演的第二个式子。那我们就可以得到 (f[i]) 的计算式了：

[f[i]=sumlimits_{ivert d}mu(frac{d}{i})g[d] ]

那么现在的问题就是如何求出 (g[i])。定义 (s[i]) 表示 (i) 的倍数的个数。我们可以知道，(i) 的倍数组成的数对的 (gcd) 一定也是 (i) 的倍数。那么我们将其两两配对得到的方案数即是 (g[i])：

[g[i]=cfrac{s[i] imes(s[i]-1)}{2} ]

而 (s[i]) 很好求，当我们加入或删除一个数的时候，直接枚举其因数修改就可以了。记录答案只需要删掉原来的再加上现在的。那么带修改的我们也完成了。

最后的答案就是 (f[1])，所以得到结论，(1) ~ (n) 中互质的数对的个数即是：

[sumlimits_{d=1}^nmu(d)g[d] ]

时间复杂度 (O(msqrt{max x}))，需要一点点卡常。

Code

#include <bits/stdc++.h>
using namespace std;
const int maxn=5e5+10;
int n,m;
int a[maxn];
long long ans;
long long s[maxn],g[maxn];
bool vis[maxn];

inline int read(){
    int x=0;bool fopt=1;char ch=getchar();
    for(;!isdigit(ch);ch=getchar())if(ch=='-')fopt=0;
    for(;isdigit(ch);ch=getchar())x=(x<<3)+(x<<1)+ch-48;
    return fopt?x:-x;
}

int pri[maxn],mobius[maxn];
bool notpri[maxn];
inline void EulerSieve(){
    notpri[0]=1;notpri[1]=1;mobius[1]=1;
    for(int i=2;i<=5e5;i++){
        if(!notpri[i]){
            pri[++pri[0]]=i;
            mobius[i]=-1;
        }
        for(int j=1;j<=pri[j]&&i*pri[j]<=5e5;j++){
            notpri[i*pri[j]]=1;
            if(i%pri[j]==0){
                mobius[i*pri[j]]=0;
                break;
            }else mobius[i*pri[j]]=-mobius[i];
        }
    }
}

signed main(){
#ifndef LOCAL
    freopen("gcd.in","r",stdin);
    freopen("gcd.out","w",stdout);
#endif
    n=read();m=read();
    for(int i=1;i<=n;i++)
        a[i]=read();
    EulerSieve();
    for(register int i=1;i<=m;i++){
        int x=read();vis[x]=!vis[x];
        int d=vis[x]?1:-1;
        for(register int j=1;j*j<=a[x];j++){
            if(a[x]%j)continue;
            s[j]+=d;
            ans-=mobius[j]*g[j];
            g[j]=s[j]*(s[j]-1)/2;
            ans+=mobius[j]*g[j];
            if(j*j!=a[x]){
                s[a[x]/j]+=d;
                ans-=mobius[a[x]/j]*g[a[x]/j];
                g[a[x]/j]=s[a[x]/j]*(s[a[x]/j]-1)/2;
                ans+=mobius[a[x]/j]*g[a[x]/j];
            }
        }
        printf("%lld
",ans);
    }
    return 0;
}