题解 P5176 【公约数】

· · 题解

感觉标准题解里的式子推错了。。。我来写一篇详细一点的题解吧。

注意：本文所有除法均为向下取整。

公式、引理及证明

引理：gcd(ij,ik,jk)=\frac {gcd(i,j)gcd(i,k)gcd(j,k)}{gcd(i,j,k)}

证明：考虑每个数质因数的幂次，设 i,j,k 关于质因子 p 的幂次分别为 a,b,c ,不妨设 a\le b\le c 。

易证 gcd(i,j) 关于 p 的幂次为 min(a,b) 。

则原式关于 p 的幂次的式子等于

min(a+b,a+c,b+c)=min(a,b)+min(a,c)+min(b,c)-min(a,b,c)

即 a+b=a+a+b-a ，显然成立。

若两个数关于各质因子的幂次都相等，则这两个数相等。

得证。

所以

Ans=\sum_{i=1}^n \sum_{j=1}^m \sum_{k=1}^p gcd^2(i,j)+gcd^2(i,k)+gcd^2(j,k)

发现每个 gcd 中都只有两个参数，第三个参数不会影响 gcd 的值，只会影响 gcd 的值出现的次数，分别枚举，得

Ans=p\times \sum_{i=1}^n\sum_{j=1}^m gcd^2(i,j) +m\times \sum_{i=1}^n\sum_{j=1}^p gcd^2(i,j) +n\times \sum_{i=1}^m\sum_{j=1}^p gcd^2(i,j)

发现要求的 gcd 的平方和都是形如 \sum_{i=1}^n\sum_{j=1}^m gcd^2(i,j) 的，莫比乌斯反演如下：

枚举 gcd 的值

\sum_{i=1}^n \sum_{j=1}^m gcd^2(i,j)=\sum_{d=1}^{min(n,m)} d^2\sum_{i=1}^n\sum_{j=1}^m [gcd(i,j)=d]

对于方括号内的内容，如果式子的值为真，则方括号的返回值为 1 ，否则为 0 。

$=\sum_{d=1}^{min(n,m)} d^2 \sum_{i=1}^{\frac n d} \sum_{j=1}^{\frac m d} [gcd(i,j)=1]

现在题目转化为求 \sum_{i=1}^a \sum_{j=1}^b [gcd(i,j)=1] 的值了。

我们设 f(x)=\sum_{i=1}^a\sum_{j=1}^b [gcd(i,j)=x] ， F(x)=\sum_{i=1}^a\sum_{j=1}^b [x|gcd(i,j)] ，易得 F(x)=\sum_{x|d} f(d)

而因为 gcd(i,j) 是 x 的倍数，当且仅当 i,j 均为 x 的倍数，所以 F(x)=\frac n x \times \frac m x

由倍数反演的式子 f(1)=\sum_{d=1}\mu(d)F(d)

=\sum_{d=1}^{min(n,m)} \mu(d)\times \frac n d\times \frac m d \therefore \sum_{i=1}^n \sum_{j=1}^m gcd^2(i,j)=\sum_{d=1}^{min(n,m)} d^2 \sum_{i=1}^{\frac n d} \sum_{j=1}^{\frac m d} [gcd(i,j)=1] =\sum_{d=1}^{min(n,m)} d^2 \sum_{g=1}^{min(\frac n d,\frac m d)} \mu(g)\times \frac n {gd}\times \frac m {gd}

枚举 T=gd ，同时枚举 d|T ，则 g=\frac T d

=\sum_{T=1}^{min(n,m)} \frac n T \times \frac m T \times \sum_{d|T} d^2\mu(\frac T d)

观察到 \frac n T 和 \frac m T 都最多分别只有 O(\sqrt n) 种取值（默认 n,m 同阶），可以用整除分块使枚举 T 的值的时间复杂度降为 O(\sqrt n) ，现在的问题是，如何在合理的预处理时间内，在 O(1) 的时间复杂度下求出 \sum_{d|x} d^2\mu(\frac x d) 的值或前缀和。

观察到 f(x)=d^2 和 g=\mu 都是 积性函数 。所求值的形式为两个函数的卷积，所以也是积性函数。考虑到 O(n) 的预处理时间是可以接受的，借助它是积性函数的性质，我们可以利用 线性筛 计算出该函数的值。

记 h(x)=\sum_{d|x} d^2\mu(\frac x d) 。

根据定义计算得， h(1)=1,h(p)=1^2\times \mu(p)+p^2\times \mu(1)=p^2-1 （p 是质数）。

由积性函数的性质得， h(xy)=h(x)h(y)~(gcd(x,y)=1) 。

现在，我们只需考虑线性筛的最后一个问题，若 x=\prod_{i=1}^k p_i^{a_i}(p_i<p_{i+1}) ，则 h(x\times p_1) 如何表示？

引理： h(x\times p_1)=h(x)\times p_1^2

证明：观察 h 函数的形式。因为 p_1 是 x 的最小素数因子，所以 x\times p_1 一定含有 p_1 的平方因子。考虑 \mu 值对答案的贡献。\mu 值为 0 时，即 \frac x d 有平方因子时，该 d 对答案无贡献。

若一个 x 的约数 d 对 h(x) 的答案有贡献，则 d\times p_1 对 h(x\times p_1) 的答案也有贡献，否则即有平方因子，即只有这些 d\times p_1 对答案有贡献。因为乘上 p_1 后贡献的值乘了平方，所以答案会乘上 p_1^2 。

得证。

同理，该结论也能扩展到任意 x 的质因子 p_k 上。

同样的，也可以类比 \phi(x)=\sum_{d|x} d\mu(\frac x d) 的线性筛法得出 h(x) 的线性筛方法。

代码展示

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int maxn=20000010;
typedef long long ll;
const ll mod=1000000007;
ll T,n,m,p,cur,pri[maxn],f[maxn],sum[maxn];
bool tf[maxn];
ll calc(ll n,ll m)
{
    ll ret=0;
    //for(ll T=1;T<=min(n,m);T++)ret=(ret+(n/T)*(m/T)%mod*f[T]%mod)%mod;
    for(ll i=1,j;i<=min(n,m);i=j+1)
    {
        j=min(n/(n/i),m/(m/i));
        ret=(ret+(n/i)*(m/i)%mod*((sum[j]-sum[i-1]+mod)%mod)%mod)%mod;
    }
    return ret;
}
int main()
{
    f[1]=1;
    for(ll i=2;i<=20000000;i++)
    {
        if(!tf[i]){pri[++cur]=i;f[i]=(i*i%mod-1+mod)%mod;}
        for(ll j=1;j<=cur&&i<=20000000/pri[j];j++)
        {
            tf[i*pri[j]]=true;
            if(i%pri[j]==0){f[i*pri[j]]=f[i]*pri[j]%mod*pri[j]%mod;break;}
            else f[i*pri[j]]=f[i]*f[pri[j]]%mod;
        }
    }
    for(ll i=1;i<=20000000;i++)sum[i]=(sum[i-1]+f[i])%mod;
    scanf("%lld",&T);
    while(T--)
    {
        scanf("%lld%lld%lld",&n,&m,&p);
        printf("%lld\n",(p*calc(n,m)%mod+m*calc(n,p)%mod+n*calc(m,p)%mod)%mod);
    }
    return 0;
}