Luogu3327 [SDOI2015]约数个数和

Luogu3327 [SDOI2015]约数个数和

蒟蒻的数论真的没救了。。。

结论：

\[d(xy)=\sum_{i|x} \sum_{j|y} [\gcd(i,j)=1] \]

考虑一个数的约数个数，显然有（\(p_i\)为质因数）：

\[n=p_1^{a_1}p_2^{a_2}\cdots p_m^{a_m}\\ d(n)=\prod_{i=1}^m a_i+1 \]

假如因数\(p\)在\(x\)中的最大幂次因子为\(p^a\)，在\(y\)中的最大幂次因子为\(p^b\)。

那么其贡献为\(a+b+1\)

观察原式\(d(xy)=\sum_{i|x} \sum_{j|y} [\gcd(i,j)=1]\)，当仅取\(x\)中\(p\)因子时，有\(a\)种方案；当仅取\(y\)中\(p\)因子时，有\(b\)种方案，不取有\(1\)种方案，总共的贡献为\(a+b+1\)。由于不同因子间相互独立，根据乘法原理，总贡献的统计是正确的。

\[\sum_{x=1}^n \sum_{y=1}^m \sum_{i|x} \sum_{j|y} [\gcd(i,j)=1]\\ =\sum_{i=1}^n \sum_{j=1}^m \lfloor \frac{n}{i} \rfloor \lfloor \frac{m}{j} \rfloor [\gcd(i,j)=1]\\ =\sum_{i=1}^n \sum_{j=1}^m \sum_{d|i,d|j} \mu(d) \lfloor \frac{n}{i} \rfloor \lfloor \frac{m}{j} \rfloor \\ =\sum_{d=1}^n \mu(d) \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^{\lfloor \frac{m}{d} \rfloor} \lfloor \frac{n}{id} \rfloor \lfloor \frac{m}{jd} \rfloor\\ =\sum_{d=1}^n \mu(d) \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \lfloor \frac{n}{id} \rfloor \sum_{j=1}^{\lfloor \frac{m}{d} \rfloor} \lfloor \frac{m}{jd} \rfloor \]

预处理：

\[f(n)=\sum_{i=1}^n \lfloor \frac{n}{i} \rfloor \]

总时间复杂度：\(O(n \sqrt n+T \sqrt n )\)。

\(Code:\)

#include<iostream>
#include<cstdio>
#include<algorithm>
#define N 50000
#define ll long long
using namespace std;
bool pri[N];
int T,n,m;
int cnt=0,prime[N],mu[N],smu[N];
ll ans=0,s[N];
void check(int n)
{
    int l,r;
    for (l=1;l<=n;l=r+1)
    {
        r=n/(n/l);
        s[n]+=(ll)(r-l+1)*(n/l);
    }
}
int main()
{
    mu[1]=1;
    for (int i=2;i<=N;++i)
    {
        if (!pri[i])
        {
            prime[++cnt]=i;
            mu[i]=-1;
        }
        for (int j=1;j<=cnt;++j)
        {
            if ((ll)i*prime[j]>N)
                break;
            int g=i*prime[j];
            pri[g]=true;
            if (i % prime[j] == 0)
            {
                mu[g]=0;
                break;
            }
            mu[g]=-mu[i];
        }
    }
    for (int i=1;i<=N;++i)
        check(i),smu[i]=smu[i-1]+mu[i];
    scanf("%d",&T);
    while (T--)
    {
        scanf("%d%d",&n,&m);
        if (n<m)
            swap(n,m);
        ans=0;
        int l,r;
        for (l=1;l<=m;l=r+1)
        {
            r=min(n/(n/l),m/(m/l));
            ans+=s[n/l]*s[m/l]*(smu[r]-smu[l-1]);
        }
        printf("%lld\n",ans);
    }
    return 0;
}