莫比乌斯反演

# 前置芝士 —— 莫比乌斯函数

# 定义

μ(n)={1n=1(1)kn无平方因数,n=p1p2p3...pk0n有大于1的平方因数μ(n) = \begin{cases} 1 & n = 1 \\ (-1)^k & n\text{无平方因数,}n = p_1p_2p_3...p_k \\ 0 & n\text{有大于}1\text{的平方因数} \\ \end{cases}

可以简化为:
nn 无平方因数时: μ(n)=(1)n的不同质因子的个数μ(n) = (-1)^{n的不同质因子的个数}
其他情况: μ(n)=0\qquad\;\;\;μ(n) = 0

# 性质

正常情况下在 nnxx 数个不同质因子,mmyy 数个不同质因子时

  • xx 奇,yy 奇,n×mn \times m 的质因子个数 =x+y== x + y = 偶, μ(n)μ(m)=(1)(1)=1μ(n) * μ(m) = (-1) * (-1) = 1
  • xx 奇,yy 偶,n×mn \times m 的质因子个数 =x+y== x + y = 奇, μ(n)μ(m)=(1)1=1μ(n) * μ(m) = (-1) * \quad1 \;\;= -1
  • xx 偶,yy 奇,n×mn \times m 的质因子个数 =x+y== x + y = 奇, μ(n)μ(m)=1(1)=1μ(n) * μ(m) = \quad1 \;\;* (-1) = -1
  • xx 偶,yy 偶,n×mn \times m 的质因子个数 =x+y== x + y = 偶, μ(n)μ(m)=11=1μ(n) * μ(m) = \quad1 \;\;* \quad1 \;\;= 1
    可以看出莫比乌斯函数是个积性函数

但是特殊情况例如 n=m=2n = m = 2
μ(n)=μ(m)=1μ(n) = μ(m) = -1
μ(nm)=0!=(1)(1)=μ(n)μ(m)μ(n * m) = 0\;\;!= (-1) * (-1) = μ(n) * μ(m)
所以莫比乌斯函数不是完全积性函数

# 推论

有两个狄利克雷卷积 (opens new window)出来的推论

1.(μId)(n)=ϕ(n)1.\;(\mu * Id)(n)=\phi(n)

2.(μ1)(n)=ε(n)2.\;(\mu * 1)(n)=\varepsilon(n)

(具体证明请看上面的狄利克雷卷积传送门

# 程序

线性筛打表:

const int maxn = 2005;

bool isprime[maxn];
ll mu[maxn];//Mobius函数表
vector<ll> prime;

inline void Mobius(){//线性筛
        isprime[0] = isprime[1] = 1;
        mu[1] = 1;//特判mu[i] = 1
        for(ll i = 2; i <= maxn; i ++){
                if( !isprime[i] ) mu[i] = -1, prime.push_back(i);//质数的质因子只有自己,所以为-1
                for(ll j = 0; j < prime.size() && i * prime[j] <= maxn; j ++){
                        isprime[i * prime[j]] = 1;
                        if(i % prime[j] == 0) break;
                        mu[i * prime[j]] = -mu[i];//积性函数性质: (i * prime[j])多出来一个质数因数(prime[j]),修正为 (-1) * mu[i]
                }
        }
        //剩余的没筛到的是其他情况,为0
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19

# 概述

# 概念

莫反是一种利用莫比乌斯函数的积性性质,对方程进行计算用时简化的一种方法

# 思想

(上文中性质的利用)

# 反演定理

设有两个方程 f(x)f(x)F(x)F(x) ,有以下两种反演方式

  • F(n)=dnf(d)F(n) = \sum_{d|n}f(d) \Rightarrow f(n)=dnμ(d)F(nd)f(n) = \sum_{d|n}{}μ(d)F(\frac nd)
证明

dnμ(d)F(nd)=dnμ(d)indf(i)=inf(i)dniμ(d)=f(n)\sum_{d|n}\mu(d)F(\frac nd) = \sum_{d|n}\mu(d)\sum_{i|\frac nd}f(i) = \sum_{i|n}f(i)\sum_{d|\frac ni}\mu(d) = f(n)

  • F(n)=ndf(d)F(n) = \sum_{n|d}f(d) \Rightarrow f(n)=ndμ(dn)F(d)f(n) = \sum_{n|d}{}μ(\frac dn)F(d)
证明

ndμ(dn)F(d)=ndμ(dn)dif(n)=(d=dn)nif(i)dinμ(d)=f(n)\sum_{n|d}\mu(\frac dn)F(d) = \sum_{n|d}\mu(\frac dn)\sum_{d|i}f(n)=(d' = \frac dn)\sum_{n|i}f(i)\sum_{d'|\frac in}\mu(d') = f(n)

# 实例

# 题目

UVA10214 《Trees in a Wood.》传送门 (opens new window) 在这里插入图片描述

# 思路

[SDOI2008]仪仗队 (opens new window)很像
在一个象限内
都是让求

i=1Nj=1M[gcd(i,j)=1]\sum_{i = 1}^{N}\sum_{j=1}^{M}[gcd(i,j)=1]

所以我们设置

f(n)=i=1Nj=1M[gcd(i,j)=n],f(1)=?f(n) = \sum_{i = 1}^{N}\sum_{j=1}^{M}[gcd(i,j)=n],\quad f(1) = \;?

但是因为 f(1)f(1) 比较难求,所以我们同时要设置一个满足 F(n)=ndf(d)F(n) = \sum_{n|d}f(d)F(n)F(n)

F(n)=i=1Nj=1M[ngcd(i,j)],F(1)=i=1Nj=1M1F(n) = \sum_{i = 1}^{N}\sum_{j=1}^{M}[n | gcd(i,j)],\quad F(1) = \sum_{i = 1}^{N}\sum_{j=1}^{M}1

F(n)=ndf(d),F(1)=d=1min(N,M)f(d)\therefore F(n) = \sum_{n|d}f(d),\quad F(1) = \sum_{d = 1}^{min(N,M)}f(d)

f(n)=ndμ(dn)F(d),f(1)=d=1min(N,M)μ(d)F(d)\therefore f(n) = \sum_{n|d}\mu(\frac dn)F(d),\quad f(1) = \sum_{d=1}^{min(N,M)}\mu(d)F(d)

1dmin(N,M)\because 1 \le d \le min(N, M)

F(d)=NdMd\therefore F(d) = \left \lfloor \frac Nd \right \rfloor * \left \lfloor \frac Md \right \rfloor

f(1)=d=1min(N,M)μ(d)NdMd\therefore f(1) = \sum_{d=1}^{min(N,M)}\mu(d) * \left \lfloor \frac Nd \right \rfloor * \left \lfloor \frac Md \right \rfloor

由于四个象限 + 四个坐标轴,所以分子为 4d=1min(N,M)μ(d)ndmd+44 * \sum_{d = 1}^{min(N, M)}μ(d)*\left \lfloor \frac nd \right \rfloor *\left \lfloor \frac md \right \rfloor + 4
分母则是所有的树 (N2+1)(M2+1)1(N * 2 + 1) * (M * 2 + 1) - 1

答案则是 4d=1min(N,M)μ(d)ndmd+4(N2+1)(M2+1)1\frac {4 * \sum_{d = 1}^{min(N, M)}μ(d)*\left \lfloor \frac nd \right \rfloor *\left \lfloor \frac md \right \rfloor + 4}{(N * 2 + 1) * (M * 2 + 1) - 1} 保留 77 位小数

# 程序

const int maxn = 2005;

bool isprime[maxn];
ll mu[maxn];//Mobius函数表
ll n, m;
vector<ll> prime;

inline void Mobius(){//线性筛
        isprime[0] = isprime[1] = 1;
        mu[1] = 1;//特判mu[i] = 1
        for(ll i = 2; i <= maxn; i ++){
                if( !isprime[i] ) mu[i] = -1, prime.push_back(i);//质数的质因子只有自己,所以为-1
                for(ll j = 0; j < prime.size() && i * prime[j] <= maxn; j ++){
                        isprime[i * prime[j]] = 1;
                        if(i % prime[j] == 0) break;
                        mu[i * prime[j]] = -mu[i];//积性函数性质: (i * prime[j])多出来一个质数因数(prime[j]),修正为 (-1) * mu[i]
                }
        }
        //剩余的没筛到的是其他情况,为0
}

inline void solve(){
        ll res = 0;
        for(ll d = 1; d <= MIN(n, m); d ++){
                res += mu[d] * (n / d) * (m / d);
        }
        res = res * 4 + 4;//四个象限 + 坐标轴的四个贡献
        ll down = (n * 2 + 1) * (m * 2 + 1) - 1;//分母,矩阵所有树 - 原点
        printf("%.7f\n", res * 1.0 / down);
}

int main () {
        Mobius();
        while(scanf("%lld%lld", &n, &m) == 2, n || m){
                solve();
        }
        return 0;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
Last Updated: 10/14/2023, 7:51:49 PM