洛谷P3327 [SDOI2015]约数个数和 题解

洛谷P3327 [SDOI2015]约数个数和 题解

题目链接:P3327 [SDOI2015]约数个数和

题意:设 d ( x ) d(x) d(x) x x x 的约数个数,给定 n , m n,m n,m ,求
∑ i = 1 n ∑ j = 1 m d ( i j ) \sum_{i=1}^{n}\sum_{j=1}^{m}d(ij) i=1nj=1md(ij)

感觉我前几篇题解的无解释推导有点不可读

所以这篇还是写的清楚一点叭

这题需要一个引理(感觉没必要记住证明的方法)

引理:设 d ( x ) d(x) d(x) x x x 的约数个数,则有
d ( i j ) = ∑ x ∣ i ∑ y ∣ j [ gcd ⁡ ( x , y ) = 1 ] d(ij)=\sum_{x \mid i}\sum_{y \mid j}[\gcd(x,y)=1] d(ij)=xiyj[gcd(x,y)=1]

证明:摘自 https://siyuan.blog.luogu.org/solution-p3327

我们考虑把每个因子一一映射。

如果 i j ij ij 的因子 k k k 中有一个因子 p c p^c pc i i i 中有因子 p a p^a pa j j j 中有因子 p b p^b pb。我们规定:

  • 如果 c ≤ a c\le a ca,那么在 i i i 中选择。
  • 如果 c > a c>a c>a,那么我们把 c c c 减去 a a a,在 j j j 中选择 p c − a p^{c-a} pca(在 jj 中选择 p e p^e pe 表示的是 p a + e p^{a+e} pa+e

对于 i j ij ij 的因子 k k k 的其他因子同理。于是对于任何一个 k k k 有一个唯一的映射,且每一个选择对应着唯一的 k k k

通过如上过程,我们发现:对于 i j ij ij的因子 k = ∏ p i c i k=\prod {p_i}^{c_i} k=pici,我们不可能同时在 i i i j j j 中选择 p i p_i pi(优先在 i i i 中选择,如果不够就只在 j j j 中选择不够的指数),故 x x x y y y 必须互质。

等式得证。

有了这个引理,就可以开始搞了

不妨假设 n ≤ m n \le m nm
∑ i = 1 n ∑ j = 1 m d ( i j ) \sum_{i=1}^{n}\sum_{j=1}^{m}d(ij) i=1nj=1md(ij)
代入
∑ i = 1 n ∑ j = 1 m ∑ x ∣ i ∑ y ∣ j [ gcd ⁡ ( x , y ) = 1 ] \sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{x\mid i}\sum_{y \mid j}[\gcd(x,y)=1] i=1nj=1mxiyj[gcd(x,y)=1]
莫比乌斯反演
∑ i = 1 n ∑ j = 1 m ∑ x ∣ i ∑ y ∣ j ∑ d ∣ gcd ⁡ ( x , y ) μ ( d ) \sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{x\mid i}\sum_{y \mid j}\sum_{d\mid \gcd(x,y)}\mu(d) i=1nj=1mxiyjdgcd(x,y)μ(d)
这个 d ∣ gcd ⁡ ( x , y ) d \mid \gcd(x,y) dgcd(x,y) 看着很难受,把它换掉
∑ i = 1 n ∑ j = 1 m ∑ x ∣ i ∑ y ∣ j ∑ d = 1 n μ ( d ) [ d ∣ gcd ⁡ ( x , y ) ] \sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{x\mid i}\sum_{y \mid j}\sum_{d=1}^{n}\mu(d)[d\mid \gcd(x,y)] i=1nj=1mxiyjd=1nμ(d)[dgcd(x,y)]
啊它怎么又回来了,不急,继续推

移一下
∑ d = 1 n μ ( d ) ∑ i = 1 n ∑ j = 1 m ∑ x ∣ i ∑ y ∣ j [ d ∣ gcd ⁡ ( x , y ) ] \sum_{d=1}^{n}\mu(d)\sum_{i=1}^{n}\sum_{j=1}^{m}\sum_{x\mid i}\sum_{y \mid j}[d\mid \gcd(x,y)] d=1nμ(d)i=1nj=1mxiyj[dgcd(x,y)]
可以发现这个正整数 x x x 的取值范围在 [ 1 , n ] [1,n] [1,n]

而它的出现次数其实就是 [ 1 , n ] [1,n] [1,n] 中有多少个数有因数 x x x

也就是 ∑ x = 1 n [ x ∣ n ] = ∑ x = 1 n ⌊ n x ⌋ \sum_{x=1}^{n}[x \mid n] = \sum_{x=1}^{n}\left\lfloor\frac{n}{x}\right\rfloor x=1n[xn]=x=1nxn

y y y 同理,则可以得到
∑ d = 1 n μ ( d ) ∑ x = 1 n ∑ y = 1 m ⌊ n x ⌋ ⌊ m y ⌋ [ d ∣ gcd ⁡ ( x , y ) ] \sum_{d=1}^{n}\mu(d)\sum_{x=1}^{n}\sum_{y=1}^{m}\left\lfloor{\dfrac{n}{x}}\right\rfloor\left\lfloor{\dfrac{m}{y}}\right\rfloor[d\mid \gcd(x,y)] d=1nμ(d)x=1ny=1mxnym[dgcd(x,y)]
考虑消掉后面那个东西

注意到 d ∣ gcd ⁡ ( d x , d y ) d \mid \gcd(dx,dy) dgcd(dx,dy) 恒成立,又因为 ∑ k ∈ K a k = ∑ p ( k ) ∈ K a p ( k ) \sum_{k\in K}a_k = \sum_{p(k) \in K}a_{p(k)} kKak=p(k)Kap(k)


∑ d = 1 n μ ( d ) ∑ d x = 1 n ∑ d y = 1 m ⌊ n d x ⌋ ⌊ m d y ⌋ \sum_{d=1}^{n}\mu(d)\sum_{dx=1}^{n}\sum_{dy=1}^{m}\left\lfloor{\dfrac{n}{dx}}\right\rfloor\left\lfloor{\dfrac{m}{dy}}\right\rfloor d=1nμ(d)dx=1ndy=1mdxndym
整理一下就是
∑ d = 1 n μ ( d ) ∑ x = 1 ⌊ n d ⌋ ⌊ n d x ⌋ ∑ y = 1 ⌊ m d ⌋ ⌊ m d y ⌋ \sum_{d=1}^{n}\mu(d)\sum_{x=1}^{\left\lfloor{\frac{n}{d}}\right\rfloor}\left\lfloor{\dfrac{n}{dx}}\right\rfloor\sum_{y=1}^{\left\lfloor{\frac{m}{d}}\right\rfloor}\left\lfloor{\dfrac{m}{dy}}\right\rfloor d=1nμ(d)x=1dndxny=1dmdym
然后就可以数论分块啦!

代码:

#include <bits/stdc++.h>
using namespace std;
#define int long long
typedef long long ll;
#define INF 0x3f3f3f3f3f3f3f3f
namespace FastIO
{
    #define gc() readchar()
    #define pc(a) putchar(a)
    #define SIZ (int)(1e6+15)
    char buf1[SIZ],*p1,*p2;
    char readchar()
    {
        if(p1==p2)p1=buf1,p2=buf1+fread(buf1,1,SIZ,stdin);
        return p1==p2?EOF:*p1++;
    }
    template<typename T>void read(T &k)
    {
        char ch=gc();T x=0,f=1;
        while(!isdigit(ch)){if(ch=='-')f=-1;ch=gc();}
        while(isdigit(ch)){x=(x<<1)+(x<<3)+(ch^48);ch=gc();}
        k=x*f;
    }
    template<typename T>void write(T k)
    {
        if(k<0){k=-k;pc('-');}
        static T stk[66];T top=0;
        do{stk[top++]=k%10,k/=10;}while(k);
        while(top){pc(stk[--top]+'0');}
    }
}using namespace FastIO;
#define N (int)(5e4+15)
int prime[N],sum[N],pcnt,mu[N],g[N];
bool ck[N];
int solve(int n)
{
    int res=0;
    for(int l=1,r; l<=n; l=r+1)
    {
        r=(n/(n/l));
        res+=(r-l+1)*(n/l);
    }
    return res;
}
void Mobius()
{
    mu[1]=1;
    for(int i=2; i<=N-5; i++)
    {
        if(!ck[i])
        {
            prime[++pcnt]=i;
            mu[i]=-1;
        }
        for(int j=1; j<=pcnt&&i*prime[j]<=N-5; j++)
        {
            int pos=i*prime[j];
            ck[pos]=1;
            if(i%prime[j])
            {
                mu[pos]=-mu[i];
            }else
            {
                mu[pos]=0;
                break;
            }
        }
    }
    for(int i=1; i<=N-5; i++)
        sum[i]+=sum[i-1]+mu[i];
    for(int i=1; i<=N-5; i++)
        g[i]=solve(i);
}
int Q,n,m;
signed main()
{
    Mobius();
    read(Q);
    while(Q--)
    {
        read(n);read(m);
        if(n>m)swap(n,m);
        int res=0;
        for(int l=1,r; l<=n; l=r+1)
        {
            r=min(n/(n/l),m/(m/l));
            res+=(sum[r]-sum[l-1])*g[n/l]*g[m/l];
        }
        write(res);pc('\n');
    }
    return 0;
}

转载请说明出处

  • 2
    点赞
  • 3
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值