P3327 [SDOI2015]约数个数和

题目描述

d(x)x 的约数个数,给定 n,m。求
\sum_{i=1}^n\sum_{j=1}^md(ij)
多组数据 1\leq T,n,m \leq 50000

Solution

d(ij)=\sum\limits_{x|i}\sum\limits_{y|j}[(x,y)=1]

原题即求
\sum_{i=1}^n\sum_{j=1}^m\sum\limits_{x|i}\sum\limits_{y|j}[(x,y)=1]
交换求和顺序,有
\sum_{x=1}^n\sum_{y=1}^m\sum_{x|i}^{i\leq n}\sum_{y|j}^{j\leq m}[(x,y)=1]
=\sum_{x=1}^n\sum_{y=1}^m\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor[(x,y)=1]

考虑这个式子还是不太好求,我们利用莫比乌斯函数来求解
f(d)=\sum_{x=1}^n\sum_{y=1}^m\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor[(x,y)=d]

F(s)=\sum_{s|d}f(d)=\sum_{x=1}^n\sum_{y=1}^m\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor[(x,y)|d]

说明后面式子如果有贡献,一定有 x|dy|d
F(s)=\sum_{s|d}f(d)=\sum_{x=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{y=1}^{\lfloor\frac{m}{d}\rfloor}\lfloor\frac{n}{x}\rfloor\lfloor\frac{m}{y}\rfloor

于是有
f(s)=\sum_{s|d}\mu(\frac{d}{s})F(d)

f(1)=\sum_{d=1}^{\min(n,m)}\mu(d)F(d)

利用整除分块优化一下即可

#include<bits/stdc++.h>
#define ll long long
const int N = 51000;
bool vis[N];
int prime[N];
int Mob[N];
int g[N],pre[N];
ll h[N];
void Mobius_sieve(){
     int cnt = 0;
     vis[1] = 1;
     Mob[1] = 1;
     for(int i = 2; i < N; i++){
        if(!vis[i])
            prime[cnt++] = i, Mob[i] = - 1 , g[i] = 1;
        for(int j = 0; j < cnt && 1LL * prime[j] * i < N; j++){
            vis[prime[j] * i] = 1;
            g[prime[j]*i] = g[i] + 1;
            Mob[i * prime[j]] = (i % prime[j] ? -Mob[i]: 0);
            if(i % prime[j] == 0)
                break;
        }
    }
}
ll H(int x){
    ll ans = 0 ;
    for(int l=1,r=0;l<=x;l=r+1){
        r = std::min(x,x/(x/l));
        ans += 1ll*(x/l)*(r-l+1);
    }
    return ans;
}
int main(){
    Mobius_sieve();
    int t = 0 ;
    for(int i=1;i<N;++i) h[i] = H(i);
    for(int i=1;i<N;++i) pre[i] = pre[i-1]+Mob[i];
    scanf("%d",&t);
    while(t--){
        ll f = 0 ;
        int n,m;
        scanf("%d%d",&n,&m);
        if(n>m) std::swap(n,m);
        for(int l=1,r=0;l<=n;l=r+1){
            r = std::min(n,std::min(n/(n/l),m/(m/l)));
//          f += 1ll*Mob[i]*h[n/i]*h[m/i];
            f += 1ll*(pre[r]-pre[l-1])*h[n/l]*h[m/l];
        }
        printf("%lld\n",f);
    }
}

0 条评论

发表评论

Avatar placeholder

您的电子邮箱地址不会被公开。 必填项已用*标注