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

[ Submit][ Status][ Discuss] 

分析：
设f(n)为n所含质因子的最大幂指数
求：

现在我们的问题就在于预处理：

然而g(n)显然不是积性函数，我们就需要用一些奇技淫巧来完成：
在线性筛的过程中，我们是按照质因数分解的思路筛掉所有的非质数

我们要注意：这个函数由两部分组成

质因数分解后：

pi前面的系数一定非零即一，而p的个数会影响mu值

• 如果存在ki≠kj（i≠j），那答案只取决于最大的幂，对于其余部分通过组合数的性质，取得的数字个数为奇数和偶数的数量相等，因此和都是0，故如果存在ki≠kj（i≠j），则g(n)=0
• 如果所有的k值都相等，那么对于任意选取方案，f值都不变，由于选取奇数个元素和偶数个元素的方案数相等，和仍然为0
  但是以上情况中，有一种选取方案的f(n)=k-1
  即d=p1*p2*p3*…pk时，因此要把那个减去，根据莫比乌斯函数性质，减去的是(-1)^k，因为是减去所以最终结果为(-1)^(k+1)，故如果不存在ai≠aj，则g(n)=(-1)^(k+1)
• g(1)=0

具体计算方法：
在线性筛的时候记录一个当前最小素数的次数（t）和去掉最小素数后余下的数（last）

#include<cstdio>
#include<cstring>
#include<iostream>
#define ll long long

using namespace std;

const int N=1e7+2;
int sshu[N],tot=0,n,m,last[N],t[N];
ll g[N];
bool no[N];

void make()
{
    g[1]=0;
    for (int i=2;i<N;i++)
    {
        if (!no[i]) 
        {
            sshu[++tot]=i; 
            last[i]=t[i]=1;
            g[i]=1;
        }
        for (int j=1;j<=tot&&sshu[j]*i<N;j++)
        {
            int x=sshu[j]*i;
            no[x]=1;
            if (i%sshu[j]==0)
            {
                last[x]=last[i];
                t[x]=t[i]+1;
                if (last[x]==1) g[x]=1;
                else g[x]=(t[last[x]]==t[x] ? -g[last[x]]:0);  
                break;
            }
            last[x]=i;
            t[x]=1;
            g[x]=(t[i]==1 ? -g[i]:0);
        }
    }

    for (int i=2;i<N;i++) g[i]+=g[i-1];
}

int main()
{
    make();
    int T;
    scanf("%d",&T);
    while (T--)
    {
        scanf("%d%d",&n,&m);
        int lst;
        ll ans=0;
        for (int i=1;i<=min(n,m);i=lst+1)
        {
            lst=min(n/(n/i),m/(m/i));
            ans=ans+(ll)(g[lst]-g[i-1])*(n/i)*(m/i);
        }
        printf("%lld\n",ans);
    }
    return 0;
}