HDU - 4676 Sum Of Gcd（莫比乌斯反演+莫队）

传送门

---

题目大意

给出数组$a$是$n$的排列，给定$l,r$，求出区间$[a_l,a_r]$所有的${\sum }_{i=l}^r{\sum }_{j=l}^{r}gcd(a_i,a_j)$

解题思路

在莫比乌斯反演那里学到了这个公式：$n={\sum }_{d|n}\phi (d)$

代入后得：${\sum }_{i=l}^r{\sum }_{j=l}^r{\sum }_{d|gcd(i,j)}\phi (d)$

然后枚举$d\in [1,n]$，得到${\sum }_{d=1}^n\phi (d){\sum }_{i=l}^r{\sum }_{j=l}^r[d|a_i][d|a_j]$

化简到这里似乎没办法继续化简了。该式子的意思是对于每一个因数$d$，从区间$[l,r]$选出两个数$a_i,a_j$，当且仅当$a_i,a_j$都是它的倍数的方案数为$x$，那么对答案的贡献为$\frac{x(x-1)}{2}\ast \phi (d)$

观察到多个区间查询，考虑到莫队，但是莫队需要$O(1)$修改，而如果枚举因数的话就不是$O(1)$了，但是此题的数据范围较小，最多只有几十个因数，实际上莫队也不会$TLE$

那么具体的做法就是在修改时枚举每一个因数，然后更新答案，这里有一个$trick$，尽管答案是$\frac{x(x-1)}{2}\ast \phi (d)$，但是我们是逐个更新的，考虑$x+1$时的贡献为$\frac{(x+1)x}{2}$，而在$x$的基础上$\frac{x(x-1)}{2}+x$恰好等于$\frac{(x+1)x}{2}$，那么我们在$add()$函数中先$sum+=x\ast phi[d]$再$x++$；同理在$del()$函数中先$x--$再$sum-=x\ast phi[d]$

//
// Created by Happig on 2020/9/28
//
#include <bits/stdc++.h>
#include <unordered_map>
#include <unordered_set>

using namespace std;
#define fi first
#define se second
#define pb push_back
#define ins insert
#define Vector Point
#define ENDL "\n"
#define lowbit(x) (x&(-x))
#define mkp(x, y) make_pair(x,y)
#define mem(a, x) memset(a,x,sizeof a);
typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
typedef pair<int, int> pii;
typedef pair<ll, ll> pll;
typedef pair<double, double> pdd;
const double eps = 1e-8;
const double pi = acos(-1.0);
const int inf = 0x3f3f3f3f;
const double dinf = 1e300;
const ll INF = 1e18;
const int Mod = 1e9 + 7;
const int maxn = 2e4 + 10;

int a[maxn], cnt[maxn], phi[maxn];
ll ans[maxn], res;
int n, m;
vector<int> prime, fac[maxn];
bitset<maxn> vis;


struct node {
    int l, r, id, pos;

    bool operator<(const node &p) {
        return pos == p.pos ? r < p.r : pos < p.pos;
    }
} q[maxn];


void init() {
    vis.reset(), prime.clear();
    phi[1] = 1;
    for (int i = 2; i < maxn; i++) {
        if (!vis[i]) {
            prime.push_back(i);
            phi[i] = i - 1;
        }
        for (int j = 0; j < prime.size() && i * prime[j] < maxn; j++) {
            vis[i * prime[j]] = 1;
            if (i % prime[j]) {
                phi[i * prime[j]] = phi[i] * phi[prime[j]];
            } else {
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
        }
    }
    for (int i = 1; i < maxn; i++) {
        for (int j = i; j < maxn; j += i)
            fac[j].push_back(i);
    }
}

inline void add(int pos) {
    for (auto i: fac[a[pos]]) {
        res += cnt[i] * phi[i];
        cnt[i]++;
    }
}

inline void del(int pos) {
    for (auto i: fac[a[pos]]) {
        cnt[i]--;
        res -= cnt[i] * phi[i];
    }
}

void solve() {
    int l = 1, r = 0;
    res = 0;
    memset(cnt, 0, sizeof cnt);
    for (int i = 1; i <= m; i++) {
        int ql = q[i].l, qr = q[i].r;
        while (l < ql) del(l++);
        while (l > ql) add(--l);
        while (r < qr) add(++r);
        while (r > qr) del(r--);
        ans[q[i].id] = res;
    }
    for (int i = 1; i <= m; i++) cout << ans[i] << ENDL;
}

int main() {
    //freopen("in.txt","r",stdin);
    //freopen("out.txt","w",stdout);
    ios_base::sync_with_stdio(0), cin.tie(0), cout.tie(0);
    int t, kase = 0;
    init();
    cin >> t;
    while (t--) {
        cin >> n;
        for (int i = 1; i <= n; i++) {
            cin >> a[i];
        }
        int sz = sqrt(n);
        cin >> m;
        for (int i = 1; i <= m; i++) {
            cin >> q[i].l >> q[i].r;
            q[i].id = i, q[i].pos = (q[i].l - 1) / sz + 1;
        }
        sort(q + 1, q + 1 + m);
        cout << "Case #" << ++kase << ":" << ENDL;
        solve();
    }
    return 0;
}