Mophues HDU - 4746【莫比乌斯反演+前缀和】

题意：

求当 \(1\leq a\leq n,1\leq b\leq m\)，有多少对 \((a,b)\) 满足 \(gcd(a,b)\) 的质因子个数不大于 \(p\)。

分析：

根据经验，莫比乌斯反演主要要会构造函数，改变枚举变量，变量代换，经常需要各种预处理来优化时间复杂度。
令 \(f(x)\) 表示 \(x\) 的质因子个数。
那么，问题可以转化为求（\(n\leq m\)）：

\[ans=\sum_{g=1}^{n}{[f(g)\leq p]\sum_{i=1}^{n}{\sum_{j=1}^{m}{[gcd(i,j)=g]}}} \]

下面开始化简：

\[=\sum_{g=1}^{n}{[f(g)\leq p]\sum_{i=1}^{\lfloor \frac{n}{g} \rfloor}{\mu(i)*\lfloor \frac{n}{i*g} \rfloor*\lfloor \frac{m}{i*g} \rfloor}} \]

令 \(T=i*g\)，则：

\[=\sum_{T=1}^{n}{\sum_{g|T}^{T}{[f(g)\leq p]*\mu(\lfloor \frac{T}{g}\rfloor)*\lfloor \frac{n}{T} \rfloor*\lfloor \frac{m}{T} \rfloor}} \]

\[=\sum_{T=1}^{n}{\lfloor \frac{n}{T} \rfloor*\lfloor \frac{m}{T} \rfloor \sum_{g|T}^{T}{[f(g)\leq p]*\mu(\lfloor \frac{T}{g}\rfloor)}} \]

观察上式发现，后面部分可以先预处理出来，前面的可以用数论分块优化。
考虑到数据范围，素因子的个数不会大于 \(20\) ，因此当 \(p\geq 20\)，可以直接输出结果为：\(n\times m\)；
预处理部分，因为 \(p\) 已经很小，用二维数组 \(a[i][k]\) 表示当 \(T=i,p=k\) 时，
\(\sum_{g|T}^{T}{[f(g)\leq p]*\mu(\lfloor \frac{T}{g}\rfloor)}\)的结果。然后借助前缀和，在数论分块时优化。

注意：一开始二维数组的第二维开的是 \(30\)，交了一发后发现 \(tle\)，找了好久才发现第二维开 \(20\)，就不会 \(t\) 了。为什么数组开大了，还会超时？

代码：

#include <bits/stdc++.h>
#define pb push_back
using namespace std;
typedef long long ll;
const int N=5e5+5;
const int maxn=5e5;
int mu[N],num[N],prime[N],total;
bool vis[N];
ll a[N][20],sum[N][20];//如果第二维开30就会tle
int divide(int x)
{
    int cnt=0;
    for(int i=1;i<=total&&1LL*prime[i]*prime[i]<=x;i++)
    {
        while(x%prime[i]==0)
        {
            cnt++;
            x/=prime[i];
        }
    }
    if(x>1)
        cnt++;
    return cnt;
}
void init()
{//预处理莫比乌斯函数:
    mu[1]=1;
    total=0;
    for(int i=2;i<N;i++)
    {
        if(!vis[i])
        {
            prime[++total]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=total&&1LL*i*prime[j]<N;j++)
        {
            vis[i*prime[j]]=1;
            if(i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                break;
            }
            else
                mu[i*prime[j]]=-mu[i];
        }
    }
    for(int i=1;i<N;i++)
        num[i]=divide(i);
    for(int i=1;i<N;i++)//g
    {
        for(int j=i;j<N;j+=i)
            a[j][num[i]]+=mu[j/i];
    }
    for(int i=1;i<N;i++)
    {
        for(int k=0;k<20;k++)
        {
            if(k>0)
                a[i][k]+=a[i][k-1];
            sum[i][k]=sum[i-1][k]+a[i][k];
        }
    }
}
ll solve(int n,int m,int p)
{
    if(p>=20)
        return 1LL*n*m;
    ll res=0;
    if(n>m) swap(n,m);
    for(int l=1,r=1;l<=n;l=r+1)
    {
        r=min(n/(n/l),m/(m/l));
        ll t=1LL*(n/l)*(m/l)*(sum[r][p]-sum[l-1][p]);
        res+=t;
    }
    return res;
}
int main()
{
    init();
    int q,n,m,p;
    scanf("%d",&q);
    while(q--)
    {
        scanf("%d%d%d",&n,&m,&p);
        printf("%lld\n",solve(n,m,p));
    }
    return 0;
}