bzoj3309: DZY Loves Math ，莫比乌斯反演

bzoj3309: DZY Loves Math

Description

对于正整数n，定义f(n)为n所含质因子的最大幂指数。例如f(1960)=f(2^3 * 5^1 * 7^2)=3, f(10007)=1, f(1)=0。
给定正整数a,b，求sigma(sigma(f(gcd(i,j)))) (i=1..a, j=1..b)。

Input

第一行一个数T，表示询问数。
接下来T行，每行两个数a,b，表示一个询问。

Output

对于每一个询问，输出一行一个非负整数作为回答。

Sample Input

4
7558588 9653114
6514903 4451211
7425644 1189442
6335198 4957

Sample Output

35793453939901
14225956593420
4332838845846
15400094813

HINT

【数据规模】
T<=10000
1<=a,b<=10^7

分析

$$ans=\sum\limits_i^n \sum\limits_j^m f(gcd(i,j))$$
$$=\sum\limits_d f(d) \sum\limits_i^n \sum\limits_j^m [gcd(i,j)==d]$$
$$=\sum\limits_d f(d) \sum\limits_i^{\lfloor\frac{n}{d}\rfloor} \sum\limits_j^{\lfloor\frac{m}{d}\rfloor} [gcd(i,j)==1]$$
$$=\sum\limits_d f(d) \sum\limits_i^{\lfloor\frac{n}{d}\rfloor} \sum\limits_j^{\lfloor\frac{m}{d}\rfloor} \sum_{d|i,d|j} \mu(d)$$
$$=\sum_d f(d) \sum_k \mu(k) \lfloor\frac{m}{kd}\rfloor \lfloor\frac{n}{kd}\rfloor$$
$$=\sum_d f(d) \sum_k \mu(k) \lfloor\frac{m}{kd}\rfloor \lfloor\frac{n}{kd}\rfloor$$
$$=\sum_D \lfloor\frac{m}{D}\rfloor \lfloor\frac{n}{D}\rfloor \sum_{d|D} f(d) \mu(\frac{D}{d})$$
一波推导之后，我们成功找到了一个 $O(nlogn+\sqrt n +\sqrt m)$的算法。
发现会超时，瓶颈在于 $O(logn)$的复杂度。
发现 $g(D)=\sum_{d|D} f(d) \mu(\frac{D}{d})$是一个狄利克雷卷积，根据套路肯定是线性筛时乱搞。
我们发掘 $g(D)$的性质，打表找规律。
0 1 1 1 1 -1 1 1 1 -1 1 0 1 -1 -1 1 1 0
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18
突然发现0的项很少，然后找0的项
1 12 18 20 24 28 40 44 45 48 50
不难发现，这些数各个素因子的次幂不同。
其次我们发现 $g$这个数组只有1和-1，凭直觉就知道和素因子个数有关。
发现奇数个不同的素因子为1，偶数个不同的素因子为0
假设这是ioi赛制，你交上去过了，就不管它了
考虑证明。
将D标准分解，$D=p_1^{a_1}p_2^{a_2} \cdots p_k^{a_k}$
我们其实把D中的素数分成两半，一半给f，一半给 $\mu$
由于如果 $\mu$函数如果里面塞超过两个以上相同素数，那么 $\mu$就是零，没有意义。所以每次最多往 $\mu$里塞同种素数中的一个。
假设 $f(D)=a$，那么我们选数最多导致 $f(\frac{D}{d})=a或a-1$
我们把D的素因子分为两个集合A，B，其中 $A=\{p_i|a_i=a\},B=\{p_i|a_i \neq a\}$
非常显然地，对于每个从A集合中选若干个数的方案，从B集合选若干个数得到的所有方案的选数的个数的奇偶方案相同。由于 $\mu$的性质会抵消掉。
只有一种特殊情况，就是 $B=\emptyset$，这个时候考虑A集合的选取，假设我们不管f值的变化，那么A集合的奇偶方案仍然是相同的，会两两被 $\mu$抵消，可是当A集合全部选的时候，f的值会减一，所以最后的答案会是1或-1，由其因数个数决定。
证明完以上性质后，我们在线性筛的时候递推就好了，具体看代码。

代码

/**************************************************************
    Problem: 3309
    User: 2014lvzelong
    Language: C++
    Result: Accepted
    Time:11360 ms
    Memory:206368 kb
****************************************************************/

#include<iostream>
#include<cstdlib>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<map>
#include<cmath>
using namespace std;
const int N = 1e7 + 10;
int read() {
    char ch = getchar(); int x = 0, f = 1;
    for(;ch < '0' || ch > '9'; ch = getchar()) if(ch == '-') f  = -1;
    for(;ch >= '0' && ch <= '9'; ch = getchar()) x = (x << 1) + (x << 3) - '0' + ch;
    return x * f;
}
bool vis[N]; int pri[N], cnt[N], mn[N], top, n, m;
long long g[N]; 
void Init(int N) {
    for(int i = 2;i <= N; ++i) {
        if(!vis[i]) {cnt[i] = 1; mn[i] = i; g[i] = 1; pri[++top] = i;}
        for(int j = 1, c = (i << 1);j <= top && c <= N; c = i * pri[++j]) {
            vis[c] = true;
            if(i % pri[j]) {
                cnt[c] = 1; mn[c] = pri[j];
                g[c] = (cnt[i] == 1 ? -g[i] : 0);
            }
            else {
                cnt[c] = cnt[i] + 1;
                mn[c] = mn[i] * pri[j]; int d = i / mn[i];
                if(d == 1) g[c] = 1;
                else g[c] = (cnt[c] == cnt[d] ? -g[d] : 0);
                break;
            }
        }
    }
    for(int i = 1;i <= N; ++i) g[i] += g[i - 1];
}
int main() {
    Init(1e7);
    for(int T = read(); T--; ) {
        n = read(); m = read(); long long ans = 0;
        for(int i = 1, pos;i <= min(n, m); i = pos + 1) {
            pos = min(n / (n / i), m / (m / i));
            ans += 1LL * (n / i) * (m / i) * (g[pos] - g[i - 1]) ;
        }
        printf("%lld\n", ans);
    }
    return 0;
}