LibreOJ 2185 「SDOI2015」约数个数和

考虑把 d(xy) 这个玩意拆开。
发现 d(xy)=∑i|x∑j|y[gcd(i,j)=1]。
然鹅并不会证明……这里传送到 Siyuan 的证明。
于是题目就是问的 n∑i=1m∑j=1∑x|i∑y|j[gcd(x,y)=1]。

对于上面式子，我们可以转化一下思路：枚举 x,y，计算其对给定范围内所有数的贡献。
于是有

n∑i=1m∑j=1∑x|i∑y|j[gcd(x,y)=1]=n∑i=1m∑j=1⌊ni⌋⌊mj⌋[gcd(i,j)=1]

然后设 f(x)=n∑i=1m∑j=1⌊ni⌋⌊mj⌋[gcd(i,j)=x],F(x)=∑x|if(i)=n∑i=1m∑j=1⌊ni⌋⌊mj⌋[x|gcd(i,j)]。

n∑i=1m∑j=1⌊ni⌋⌊mj⌋[x|gcd(i,j)]=∑xi≤n∑xj≤m⌊nxi⌋⌊mxj⌋

那么有 f(x)=∑x|iμ(ix)F(i)，此题求 f(1)=n∑i=1μ(i)F(i)=n∑i=1μ(i)∑xi≤n⌊⌊nx⌋i⌋∑yi≤m⌊⌊my⌋i⌋。

对于任意的 n，预处理一下 n∑i=1⌊ni⌋ 和 μ 的前缀和然后数论分块。

代码：

#include <cstdio>
#include <algorithm>
using namespace std;
const int N = 5e4;
int t;
int cnt,prime[N + 5],mu[N + 5],vis[N + 5];
long long sum[N + 5];
int main()
{
    mu[1] = 1;
    for(register int i = 2;i <= N;++i)
    {
        if(!vis[i])
            prime[++cnt] = i,mu[i] = -1;
        for(register int j = 1;j <= cnt && i * prime[j] <= N;++j)
        {
            vis[i * prime[j]] = 1;
            if(!(i % prime[j]))
                break;
            else
                mu[i * prime[j]] = -mu[i];
        }
    }
    for(register int i = 1;i <= N;++i)
    {
        mu[i] += mu[i - 1];
        for(register int l = 1,r;l <= i;l = r + 1)
        {
            r = i / (i / l);
            sum[i] += (long long)(i / l) * (r - l + 1);
        }
    }
    scanf("%d",&t);
    int n,m;
    long long ans;
    while(t--)
    {
        ans = 0;
        scanf("%d%d",&n,&m);
        if(n > m)
            swap(n,m);
        for(register int l = 1,r;l <= n;l = r + 1)
        {
            r = min(n / (n / l),m / (m / l));
            ans += (long long)(mu[r] - mu[l - 1]) * sum[n / l] * sum[m / l];
        }
        printf("%lld\n",ans);
    }
}