莫比乌斯反演与数论分块

数学反演总结

1 数论分块

1.1 定理

定理 1.1.1

\[\left\lfloor\frac {\left\lfloor \frac{x}{b} \right\rfloor}{c} \right\rfloor=\left\lfloor\frac{x}{bc}\right\rfloor \]

其中 \(x\in\R,b,c\in\N\)

  • 证明:设 \(x=kbc+r\),其中 \(r\in[0,bc)\) ,我们把 \(x\) 带入左边式子可以的得到:

    \[k+\left\lfloor\frac {\left\lfloor \frac{r}{b} \right\rfloor}{c} \right\rfloor \]

    因为:

    \[0\le \left\lfloor\frac {\left\lfloor \frac{r}{b} \right\rfloor}{c} \right\rfloor\le \left\lfloor\frac {r}{bc} \right\rfloor=0 \]

    所以我们可以得到:

    \[\left\lfloor\frac {\left\lfloor \frac{r}{b} \right\rfloor}{c} \right\rfloor=k \]

    显然,右边式子也等于 \(k\),所以结论成立。

1.2 算法讲解

其实数论分块(又被称为除法分块)是可以在 \(\sqrt n\) 的复杂度内算出:

\[\sum\limits_{i=1}^n\lfloor \frac{n}i \rfloor \]

这是因为我们注意到和式求和的部分的取值只有 \(2\sqrt n\) 中,证明这个结论只需要讨论 \(i\le \sqrt n\) 和 \(i> \sqrt n\)​ 两种情况即可。

一般的,对于:

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

我们要求其值相等的左右边界,只需要让:

\[\lfloor\frac{n}{f(l)}\rfloor=\lfloor\frac{n}{f(r)}\rfloor \]

由数论分块的正确性可以得到:

\[f(r)=\left\lfloor \frac{n}{\lfloor \frac{n}{f(l)}\rfloor}\right\rfloor \]

我们求解 \(f(l)\) 并解方程解出 \(r\) 为多少即可得出数论分块区间。

1.3 例题

1.3.1 P2261 [CQOI2007]余数求和

我们直接推式子:

\[\sum\limits_{i=1}^nk\bmod i=\sum\limits_{i=1}^nk-\lfloor \frac ki\rfloor\times i=nk-\sum\limits_{i=1}^ni\lfloor\frac ik\rfloor \]

可以应用数论分块,注意后面要特判是否等于 \(0\)。

代码:

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define uint unsigned int
#define ull unsigned long long
#define N number
#define M number
using namespace std;

const int INF=0x3f3f3f3f;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

inline ll calc(int l,int r){
    return 1ll*(r+l)*(r-l+1)/2;
}

ll ans=0,n,k;

int main(){
    read(n);read(k);
    for(int i=1;i<=n;){
        int j;
        if(k/i>0) j=k/(k/i);
        else j=n;
        j=min((ll)j,n);
        ans+=(k/i)*calc(i,j);
        i=j+1;
    }
    printf("%lld\n",n*k-ans);
    return 0;
}

代码第 \(31\) 行 \(32\) 行为特判。

1.3.2 CF1485C Floor and Mod

推式子:

\[\left\lfloor \frac{a}{b} \right\rfloor=a\bmod b=a-\lfloor\frac{a}{b}\rfloor\times b\\ \Rightarrow \lfloor\frac{a}{b}\rfloor\times (b+1) =a\Rightarrow \lfloor\frac{a}{b}\rfloor=\frac{a}{b+1} \]

所以我们只需要枚举 \(b\) ,同时统计有多少 \(a\) 是 \(b+1\) 的倍数,这个数应该是 \(\frac{x}{b+1}\) 。

但是需要注意的一点是,我们需要保证 \(\lfloor\frac{a}{b}\rfloor=\frac{a}{b+1}<b\) ,所以我们有 \(a<b^2+b\) ,所以上面这个式子应该变成:\(\frac{\min(x,b^2+b-1)}{b+1}\)。

我们根号分治并分开讨论,数论分块即可。复杂度为 \(O(\sqrt n)\)。

代码:

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define int long long
#define uint unsigned int
#define ull unsigned long long
#define N number
#define M number
using namespace std;

const int INF=0x3f3f3f3f;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

int t;
ll x,y;

inline ll calc(ll x,ll y){
    ll ans=0;ll b=1;
    for(;b*b+b-1<=x&&b<=y;b++) ans+=(b*b+b-1)/(b+1);
    for(int i=b;i<=y;){
        int j=(x/(i+1)==0)?y:(x/(x/(i+1))-1);
        j=min((ll)j,y);
        ans+=(x/(i+1))*(j-i+1);
        i=j+1;
    }
    return ans;
}

signed main(){
    // freopen("my.in","r",stdin);
    // freopen("my.out","w",stdout);
    read(t);
    while(t--){
        read(x);read(y);
        printf("%lld\n",calc(x,y));
    }
    return 0;
}

2 莫比乌斯反演(进阶版)

以前博主整理过一篇莫比乌斯反演,讲的比较基础,这里补充一下那篇博客里没有的内容。

2.1 常用公式

2.1.1 最大公因数,最小公倍数

\[\sum\limits_{i=1}^n\sum\limits_{j=1}^m\text{lcm}(i,j)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\frac{ij}{\gcd(i,j)} \\ =\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{d|i,d|j} [\gcd(\frac id,\frac jd)=1]\frac{ij}d\\ =\sum\limits_{d=1}^{\min(n,m)}d\sum\limits_{i=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{m}{d} \right\rfloor}[\gcd(i,j)=1]ij\\ \]

令 \(f(n,m)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m[\gcd(i,j)=1]ij\),那么有:

\[f(n,m)=\sum\limits_{i=1}^n\sum\limits_{i=1}^m\epsilon(\gcd(i,j))ij\\ =\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{d|\gcd(i,j)}\mu(d)ij\\ =\sum\limits_{d=1}^{\min(n,m)}\mu(d)d^2\sum\limits_{i=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{m}{d} \right\rfloor}ij\\ \]

令 \(g(n,m)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m ij\),显然,我们有 \(g(n,m)=\frac{(n+1)n}{2}\frac{(m+1)m}{2}\) ,可以 \(O(1)\) 计算,所以通过处理前缀和,\(f(n,m)\) 可以在 \(O(\sqrt n+\sqrt m)\) 的时间1复杂度算出,所以原式的复杂度可以做到 \(O(n+m+\sqrt n+\sqrt m)\) 。

2.1.2 除数函数

公式 2.1.2.1
  • \(\sigma_k(n)=\sum\limits_{d|n}d^k\)​​ ,令 \(id_k(n)=n^k\)​​ ,那么我们有 \(\sigma_k=id_k*I\)​​
公式2.1.2.2
  • \(\sigma_0(xy)=\sum\limits_{i|x}\sum\limits_{j|y}[\gcd(i,j)=1]\)

    证明:我们考虑将 \(xy\)​​ 的每一个因子都与互质的 \(i,j\)​​ 意义映射。设 \(xy\)​​ 的因子为 \(k\)​​ ,设把 \(k\)​​ 唯一分解之后得到的式子为 \(\prod\limits_{i=1}^mp_i^{\alpha_i}\)​​​ 我们考虑每一个质数,设当前的质数为 \(p\)​​ ,设这个质数在 \(k\)​​ 中的指数为 \(\alpha\)​​ ,设在 \(x\)​​ 中的指数为 \(a\)​​ 设在 \(y\)​​ 中的指数为 \(b\)​​​ 。我们规定,如果 \(\alpha \le a\)​ ,那么令从 \(x\)​ 中拿取 \(p^{\alpha}\)​ ,否则就从 \(y\)​ 中拿取 \(p^{\alpha -a}\)​ 次方。这样,从 \(x\)​​ 中拿出来数就与从 \(y\)​​ 中拿出来的数互质。同时每一对互质的数,都可以映射到这么一个因子上。所以等式成立。证毕。这是一个重要的式子。

定理 2.1.2.1
  • 除数函数是积性函数。

    证明:如果 \(x,y\) 互质,那么 \(xy\) 的一个因子 \(k\) ,\(k\) 一定可以分解成两个数,这两个数互质。也就是说,\(k\) 一部分来自 \(x\) ,一部分来自 \(y\)​ 。证毕。

    上面这个积性函数的性质,揭示我们除数函数也是可以线性筛的。不过我们有计算一些特殊除数函数前缀和更快的做法。

公式 2.1.2.3
  • \(\sum\limits_{i=1}^n\sigma_0(i)=\sum\limits_{j=1}^n\left\lfloor \frac{n}{j} \right\rfloor\)​

    证明:

    \[\sum\limits_{i=1}^n\sigma_0(i)=\sum\limits_{i=1}^n\sum\limits_{j=1}^i[j|i]=\sum\limits_{j=1}^n\sum\limits_{i=1}^{\left\lfloor \frac{n}{j} \right\rfloor}1=\sum\limits_{j=1}^n \left\lfloor \frac{n}{j} \right\rfloor \]

    整除分块即可。

公式 2.1.2.4
  • \(\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sigma_0(ij)=\sum\limits_{d=1}^{\min(n,m)}\mu(d)(\sum\limits_{x=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\left\lfloor \frac{n}{xd} \right\rfloor)(\sum\limits_{y=1}^{\left\lfloor \frac{m}{d} \right\rfloor}\left\lfloor \frac{m}{yd} \right\rfloor)\)​

    证明:

    \[\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sigma_0(ij)=\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{x|i}\sum\limits_{y|j}\epsilon(\gcd(x,y))\\ =\sum\limits_{x=1}^n\sum\limits_{y=1}^m \epsilon(\gcd(x,y)) \sum\limits_{i=1}^{\left\lfloor \frac{n}{x} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{m}{y} \right\rfloor}1\\ =\sum\limits_{x=1}^n\sum\limits_{y=1}^m\sum\limits_{d|\gcd(x,y)}\mu(d)\left\lfloor \frac{n}{x} \right\rfloor\left\lfloor \frac{m}{y} \right\rfloor\\ =\sum\limits_{d=1}^{\min(n,m)}\sum\limits_{x=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sum\limits_{y=1}^{\left\lfloor \frac{m}{d} \right\rfloor}\mu(d)\left\lfloor \frac{n}{xd} \right\rfloor\left\lfloor \frac{m}{yd} \right\rfloor\\ =\sum\limits_{d=1}^{\min(n,m)}\mu(d)(\sum\limits_{x=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\left\lfloor \frac{n}{xd} \right\rfloor)(\sum\limits_{y=1}^{\left\lfloor \frac{m}{d} \right\rfloor}\left\lfloor \frac{m}{yd} \right\rfloor) \]

公式 2.1.2.5
  • \(\sigma_1(xy)=\sum\limits_{i|x}\sum\limits_{j|y}\frac{jx}{i}[\gcd(i,j)=1]\)​​

    证明:同样还是构造一一对应来证明,构造方法和上面证明 \(2.1.2.2\) 相同,设:

    \[x=\prod\limits_{i=1}^np_i^{\alpha_i},y=\prod\limits_{j=1}^np_j^{\beta_j}\\ i=\prod\limits_{j=1}^kp_{i_j}^{a_{i_j}},j=\prod_{i=1}^{n-k}p_{j_i}^{a_{j_i}} \]

    我们看 \(\frac{jx}i\)​ 的值表示出来是多少,我们发现这个值是:

    \[\prod\limits_{j=1}^k p_{i_j}^{\alpha_{i_j}-a_{i_j}}\times \prod\limits_{i=1}^{n-k}p_{j_i}^{\alpha_{j_i}+a_{j_i}} \]

    这个 \(i,j\)​ 对应的因子是 \(\prod_{j=1}^kp_{i_j}^{a_{i_j}}\prod_{p_{j_i}}^{n-k}p_{j_i}^{\alpha_{j_i}+a_{j_i}}\)​,不难发现,这个因子与上面相比,右边相同,左边相乘是 \(x\)​​​ ,所以这两个因子是一一对应的。显然这个式子是正确的。

公式 2.1.2.6
  • \(\sum\limits_{i=1}^n\sigma_1(i)=\sum\limits_{i=1}^ni\left\lfloor \frac{n}{i} \right\rfloor\)

    证明:

    \[\sum\limits_{i=1}^n\sigma_1(i)=\sum\limits_{i=1}^n\sum\limits_{j=1}^ij[j|i]=\sum\limits_{j=1}^n\sum\limits_{i=1}^{\left\lfloor \frac{n}{j} \right\rfloor}j=\sum\limits_{i=1}^nj\left\lfloor \frac{n}{j} \right\rfloor \]

公式 2.1.2.7
  • \(\sum\limits_{i=1}^n\sum\limits_{j=1}^n\sigma_1(ij)=\sum\limits_{d=1}^n\mu(d)d(\sum\limits_{i=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sigma_1(i))^2\)​

    证明:

    \[\sum\limits_{i=1}^n\sum\limits_{j=1}^n\sigma_1(ij)=\sum\limits_{i=1}^n\sum\limits_{j=1}^n\sum\limits_{x|i}\sum\limits_{y|j}\frac{jx}{y}[\gcd(x,y)=1]\\ =\sum\limits_{x=1}^{n}\sum\limits_{y=1}^n\sum\limits_{i=1}^{\left\lfloor \frac{n}{x} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{n}{y} \right\rfloor}jx\times \epsilon(\gcd(x,y))\\ =\sum\limits_{d=1}^n\mu(d)d \sum\limits_{x=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sum\limits_{y=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sum\limits_{i=1}^{\left\lfloor \frac{n}{xd} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{n}{yd} \right\rfloor}jx\\ =\sum\limits_{d=1}^{n}\mu(d)d\sum\limits_{x=1}^{\left\lfloor \frac{n}{d} \right\rfloor}x\left\lfloor \frac{n}{xd} \right\rfloor\sum\limits_{y=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{n}{yd} \right\rfloor}j\\ =\sum\limits_{d=1}^{n}\mu(d)d\sum\limits_{x=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sigma_1(x)\sum\limits_{y=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sigma_1(y)\\ =\sum\limits_{d=1}^{n}\mu(d)d(\sum\limits_{i=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sigma_1(i))^2 \]

公式 2.1.2.8
  • \(\sum\limits_{i=1}^n\mu^2(i)=\sum\limits_{i=1}^{\sqrt n}\mu(i)\times\left\lfloor \frac{n}{i^2} \right\rfloor\)

    这个公式需要构造来证明,具体证明过程请参照例题 \(2.2.4\)。

公式 2.1.2.9
  • \(\sum\limits_{i=1}^n\mu^2(i)\left\lfloor \sqrt{\frac{n}{i}} \right\rfloor\)​​

    这个东西的证明和上面类似,我们同样考虑构造。考虑 \(1\) 到 \(n\) 内所有不能被完全平方数整除的数组成的集合 \(P\) ,那么 \(1\) 到 \(n\) 内的每一个数都可以表示成 \(xy^2\) 的形式,其中 \(x\in P,y\le \sqrt{\frac nx}\),上面的式子其实等价于 \(\sum\limits_{p\in P}\left\lfloor \sqrt{\frac{n}{i}} \right\rfloor\)​,根据 \(y\) 的取值范围,不难发现这个式子的正确性。

公式 2.1.2.10
  • \(\sum\limits_{i=1}^n\sum\limits_{j=1}^n[\gcd(i,j)=1]=\sum\limits_{d=1}^n\mu(d)\left\lfloor \frac{n}{d} \right\rfloor=2S_{\phi}(\left\lfloor \frac{n}{d} \right\rfloor)-1\)

    证明:

    这里只证明后面的那个式子。

    设 \(f(n)=\sum\limits_{i=1}^n\sum\limits_{j=1}^n[\gcd(i,j)=1],g(n)=\sum\limits_{i=1}^n\sum\limits_{j=1}^i[\gcd(i,j)=1]\),那么 \(g(n)=S_{\phi}(i)\),那么容易验证有:\(f(n)=2g(n)-1\) ,减 \(1\) 是因为 \((1,1)\) 被算重了。

公式 2.1.2.11
  • \(\prod\limits_{i=1}^n\prod\limits_{j=1}^n\gcd^2(i,j)=(\prod\limits_{i=1}^nd^{S_{\phi}(\left\lfloor \frac{n}{d} \right\rfloor)})\)

    证明:

    \[\prod\limits_{i=1}^n\prod\limits_{j=1}^n\gcd{^2}(i,j)=\prod\limits_{i=1}^nd^{\sum\limits_{i=1}^n\sum\limits_{j=1}^n[\gcd(i,j)=d]}\\ =\prod\limits_{i=1}^nd^{\sum\limits_{i=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{n}{d} \right\rfloor}\gcd(i,j)=1} \]

    由公式 \(2.1.2.10\) 可以知道,该公式成立。

2.2 例题

上面的每一个公式按顺序依次标号。

2.2.1 洛谷P3935

运用上面公式 \(2.1.2.3\)​​​ 即可。

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define uint unsigned int
#define ull unsigned long long
#define N number
#define M number
using namespace std;

const int INF=0x3f3f3f3f;
const int mod=998244353;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

ll l,r;

inline int Solve(ll k){
    if(k<1) return 0;
    int ans=0;
    for(ll i=1,j;i<=k;i=j+1){
        j=k/(k/i);
        ans=1ll*(ans+1ll*k/i*(j-i+1)%mod)%mod;
    }
    return ans;
}

int main(){
    read(l);read(r);
    printf("%d\n",((Solve(r)-Solve(l-1))%mod+mod)%mod);
    return 0;
}

2.2.2 洛谷P3327

公式 \(2.1.2.4\)​

代码:

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define uint unsigned int
#define ull unsigned long long
#define N 50010
#define M number
using namespace std;

const int INF=0x3f3f3f3f;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

bool NotPrime[N];
int Prime[N],tail,Mu[N],SumMu[N];

inline void Getmu(int n){
    NotPrime[1]=1;Mu[1]=1;
    for(int i=2;i<=n;i++){
        if(!NotPrime[i]) Prime[++tail]=i,Mu[i]=-1;
        for(int j=1;j<=tail&&i*Prime[j]<=n;j++){
            int k=i*Prime[j];NotPrime[k]=1;
            if(i%Prime[j]==0){Mu[k]=0;break;}
            else Mu[k]=Mu[i]*Mu[Prime[j]];
        }
    }
    for(int i=1;i<=n;i++) SumMu[i]=SumMu[i-1]+Mu[i];
}

ll f[N];

inline ll Get(int k){
    ll ans=0;
    for(int i=1,j;i<=k;i=j+1){
        j=k/(k/i);ans+=(k/i)*(j-i+1);
    }
    return ans;
}

inline ll Solve(int n,int m){
    int minn=min(n,m);
    ll ans=0;
    for(int i=1,j;i<=minn;i=j+1){
        int j1=n/(n/i),j2=m/(m/i);j=min(j1,j2);
        ans+=(SumMu[j]-SumMu[i-1])*f[n/i]*f[m/i];
    }
    return ans;
}

int t;

int main(){
    read(t);Getmu(50000);
    for(int i=1;i<=50000;i++) f[i]=Get(i);
    while(t--){
        int n,m;read(n);read(m);
        printf("%lld\n",Solve(n,m));
    }
    return 0;
}

2.2.3 51nod1220

公式 \(2.1.2.7\)

代码:

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define int long long
#define uint unsigned int
#define ull unsigned long long
#define N 1000010
#define M number
using namespace std;

const int INF=0x3f3f3f3f;
const int mod=1e9+7;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

map<int,ll> Map,f;
const ll inv2=500000004;

inline ll CalcSumSigma_1(int n){
    ll nowans=0;
    for(int i=1,j;i<=n;i=j+1){
        j=n/(n/i);
        nowans=(nowans+1ll*(n/i)*(i+j)%mod*(j-i+1)%mod*inv2%mod)%mod;
    }
    return nowans;
}

int Prime[N],F[N],SumF[N],tail,MAXN;
bool NotPrime[N];

inline void GetSumF(){
    int n=1000000;MAXN=n;
    NotPrime[1]=1;F[1]=1;
    for(int i=2;i<=n;i++){
        if(!NotPrime[i]) Prime[++tail]=i,F[i]=-i;
        for(int j=1,k;j<=tail&&(k=i*Prime[j])<=n;j++){
            NotPrime[k]=1;
            if(i%Prime[j]==0){F[k]=0;break;}
            else F[k]=F[Prime[j]]*F[i]%mod;
        }
    }
    for(int i=1;i<=n;i++) SumF[i]=(SumF[i-1]+F[i])%mod;
}

inline ll CalcSumF(int n){
    if(n<=MAXN) return SumF[n];
    if(f[n]) return f[n];
    ll now=(n>=1);
    for(int i=2,j;i<=n;i=j+1){
        j=n/(n/i);
        now=(now-1ll*CalcSumF(n/i)*(i+j)%mod*(j-i+1)%mod*inv2%mod)%mod;
    }
    f[n]=now;
    return now;
}

int n;

inline void GetSumSigma_1(){
    for(int i=1,j;i<=n;i=j+1){
        j=n/(n/i);
        ll val=CalcSumSigma_1(n/i)%mod;
        Map[n/i]=val*val%mod;
    }
}

inline void PreWork(){
    GetSumSigma_1();GetSumF();
}

signed main(){
    // freopen("my.in","r",stdin);
    // freopen("my.out","w",stdout);
    read(n);
    PreWork();
    // cout<<"here\n";
    // for(int i=1;i<=100;i++){
    //     printf("%lld\n",SumF[i]);
    // }
    ll ans=0;
    for(int i=1,j;i<=n;i=j+1){
        j=n/(n/i);
        // printf("i:%lld\n",i);
        // cout<<Map[n/i]<<endl;
        // cout<<(CalcSumF(j)-CalcSumF(i-1))<<endl;
        ans=(ans+1ll*(CalcSumF(j)-CalcSumF(i-1))*Map[n/i]%mod)%mod;
    }
    printf("%lld\n",(ans%mod+mod)%mod);
    return 0;
}

2.2.4 洛谷P4318

不难发现,答案具有可二分性。我们规定,如果一个数不能被完全平方数整除,则称其合法,否则称其不合法。那么 \(n\) 以内的合法的数一共有多少呢?不难发现,答案应该是:

\[\sum\limits_{i=1}^n\mu^2(i) \]

即为 \(1\) 到 \(n\) 中所有 \(\mu(x)\) 不为 \(0\) 的数。

我们考虑用容斥来计算这个东西,不难发现这个东西其实也是可以这样计算的:

含 \(0\) 个质数的平方的倍数减去含 \(1\) 个质数的平方的倍数加上含 \(2\) 个质数平方的倍数......,考虑莫比乌斯函数的定义,我们发现上面这个式子其实也是:

\[\sum\limits_{i=1}^{i\times i\le n}\mu(i)\times\left\lfloor \frac{n}{i^2} \right\rfloor \]

所以我们就可以二分这个答案,然后在 \(O(\sqrt n)\) 的时间复杂度内判断合不合法。

代码:

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define uint unsigned int
#define ull unsigned long long
#define N 60010
#define M number
using namespace std;

const int INF=0x3f3f3f3f;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

int Prime[N],tail,Mu[N];
bool NotPrime[N];

inline void PreWork(int n){
    NotPrime[1]=1;Mu[1]=1;
    for(int i=2;i<=n;i++){
        if(!NotPrime[i]) Prime[++tail]=i,Mu[i]=-1;
        for(int j=1,k;j<=tail&&(k=i*Prime[j])<=n;j++){
            NotPrime[k]=1;
            if(i%Prime[j]==0) break;
            else Mu[k]=-Mu[i];
        }
    }
}

int t;

inline bool Check(int mid,int k){
    int res=0;
    for(ll i=1;i*i<=mid;i++){
        res=res+Mu[i]*(mid/(i*i));
    }
    return res>=k;
}

int main(){
    // freopen("my.in","r",stdin);
    // freopen("my.out","w",stdout);
    PreWork(60000);
    read(t);
    while(t--){
        ll k;read(k);
        ll l=k,r=2*k;
        while(l<r){
            // cout<<"l"<<" "<<l<<endl;

            int mid=(l+r)>>1;
            if(Check(mid,k)) r=mid;
            else l=mid+1;
        }
        printf("%lld\n",l);
    }
}

2.2.5 SP5971 LCMSUM - LCM Sum

要求计算 \(\sum_{i=1}^nlcm(i,n)\) ,我们尝试化简这个东西:

\[\sum\limits_{i=1}^nlcm(i,n)=\sum\limits_{i=1}^n\frac{in}{\gcd(i,n)}=\sum\limits_{i=1}^n\sum\limits_{d|i,d|n}[\gcd(\frac di,\frac nd)=1]\times \frac{in}{d}\\ =\sum\limits_{d|n}\sum\limits_{i=1}^{\left\lfloor \frac{n}{d} \right\rfloor}in[\gcd(i,\frac nd)=1] \]

我们令 \(g(x)=\sum\limits_{i=1}^x[\gcd(i,x)=1]\times i\) 即在 \(1\) 到 \(x\) 总所有与 \(x\) 互质的数的个数,因为 \(\gcd(i,n)=1\Leftrightarrow\gcd(n-i,i)\) ,所以互质的数总是两两一对,且和为 \(n\) 。

所以这个和应该是 \(\frac{\phi(x)x}{2}\)。

所以上面那个式子可以化简为:

\[n\sum\limits_{d|n}\frac{\phi(d)d}{2} \]

我们可以用 \(n\log n\) 的时间复杂度完成这个事情。

这个地方有一点不太严谨,当 \(d=1\) 的时候注意 \(g(1)=1\),这个是个特殊的。

#include<bits/stdc++.h>
#define dd double
#define ld long double
#define ll long long
#define uint unsigned int
#define ull unsigned long long
#define N 1000100
#define M number
using namespace std;

const int INF=0x3f3f3f3f;

template<typename T> inline void read(T &x) {
    x=0; int f=1;
    char c=getchar();
    for(;!isdigit(c);c=getchar()) if(c == '-') f=-f;
    for(;isdigit(c);c=getchar()) x=x*10+c-'0';
    x*=f;
}

int Prime[N],tail,Phi[N];
ll g[N],f[N];
bool NotPrime[N];

inline void PreWork(int n){
    NotPrime[1]=1;Phi[1]=1;
    for(int i=2;i<=n;i++){
        if(!NotPrime[i]) Prime[++tail]=i,Phi[i]=i-1;
        for(int j=1,k;j<=tail&&(k=i*Prime[j])<=n;j++){
            NotPrime[k]=1;
            if(i%Prime[j]==0){
                Phi[k]=Phi[i]*Prime[j];break;
            }
            else Phi[k]=Phi[i]*Phi[Prime[j]];
        }
    }
    // for(int i=1;i<=5;i++) printf("%d ",Phi[i]);puts("");
    g[1]=1;
    for(int i=2;i<=n;i++) g[i]=1ll*Phi[i]*i/2;
    // for(int i=1;i<=5;i++) printf("%lld ",g[i]);puts("");
    for(int i=1;i<=n;i++)
        for(int j=i;j<=n;j+=i) f[j]+=g[i];
    for(int i=1;i<=n;i++) f[i]*=i;
}

int t,n;

int main(){
    // freopen("my.in","r",stdin);
    // freopen("my.out","w",stdout);
    PreWork(1000000);
    read(t);
    while(t--){
        read(n);printf("%lld\n",f[n]);
    }
}//

2.2.6 P3768 简单的数学题

首先我们尝试化简式子:

\[\sum\limits_{i=1}^n\sum\limits_{j=1}^nij\gcd(i,j)=\sum\limits_{k=1}^n\sum\limits_{k|i,i\le n}\sum\limits_{k|j,j\le n}ijk[\gcd(\left\lfloor \frac{k}{i} \right\rfloor,\left\lfloor \frac{k}{j} \right\rfloor)=1]\\ =\sum\limits_{k=1}^n\sum\limits_{i=1}^{\left\lfloor \frac{n}{k} \right\rfloor}\sum\limits_{j=1}^{\left\lfloor \frac{n}{k} \right\rfloor}ijk^3[\gcd(i,j)=1] \]

令 \(g(n)=\sum_{i=1}^n\sum_{j=1}^nij[\gcd(i,j)=1]\),我们可以得到:

\[g(n)=\sum\limits_{i=1}^ni\sum\limits_{j=1}^nj[\gcd(i,j)=1]=[n\ge 1]+\sum\limits_{i=2}^n\frac{\phi(i)i}{2} \]

所以原来的式子我们可以化简为:

\[\sum\limits_{k=1}^nk^3g(\left\lfloor \frac{n}{k} \right\rfloor) \]

由于 \(n\) 是 \(1e9\) 级别的,所以我们可以用杜教筛筛一下 \(g(n)\) 。然后用整除分块。因为杜教筛也是整除分块计算的,所以杜教筛外面套一个整除分块没有问题。

莫比乌斯函数扩展

莫比乌斯函数非卷积形式。

对于数论函数和完全积性函数 \(t\) 且 \(t(1)=1\) ,有以下等价式子:

\[f(n)=\sum\limits_{i=1}^nt(i)g(\left\lfloor \frac{n}{i} \right\rfloor)\Leftrightarrow g(n)=\sum\limits_{i=1}^nt(i)\mu(i)f(\left\lfloor \frac{n}{i} \right\rfloor) \]

  • 证明:

    \[g(n)=\sum\limits_{i=1}^nt(i)\mu(i)f(\left\lfloor \frac{n}{i} \right\rfloor)\\ =\sum\limits_{i=1}^nt(i)\mu(i)\sum\limits_{j=1}^{\left\lfloor \frac{n}{i} \right\rfloor}t(j)g(\left\lfloor \frac{n}{ij} \right\rfloor)\\ =\sum\limits_{d=1}^n\sum\limits_{i|d} t(i)\mu(i)t(\frac{d}{i})g(\left\lfloor \frac{n}{d} \right\rfloor)\\ =\sum\limits_{d=1}^ng(\left\lfloor \frac{n}{d} \right\rfloor)t(d)\sum\limits_{i|d}t(i)\\ =\sum\limits_{d=1}^ng(\left\lfloor \frac{n}{d} \right\rfloor)t(d)\epsilon(d)\\ =g(n) \]

上一篇:【luogu P4213】【模板】杜教筛(Sum)(数学)(整除分块)


下一篇:海量IT资料 + 各种平台下的Oracle安装文件 + 公开课录像 + 各种视频教程资料