莫比乌斯函数定义:
设 n = p 1 k 1 ⋅ p 2 k 2 ⋅ ⋯ ⋅ p m k m n = p_1 ^ {k_1} \cdot p_2 ^ {k_2} \cdot\cdots\cdot p_m ^ {k_m} n=p1k1⋅p2k2⋅⋯⋅pmkm,其中 p 为素数,则定义如下:
μ ( n ) = { 1 n = 1 ( − 1 ) m ∏ i = 1 m k i = 1 0 otherwise ( k i > 1 ) \mu(n) = \begin{cases} 1 & n = 1 \\ (-1) ^ m & \prod\limits_{i = 1} ^ {m} k_i = 1 \\ 0 & \textrm{otherwise}(k_i \gt 1) \end{cases} μ(n)=⎩⎪⎪⎨⎪⎪⎧1(−1)m0n=1i=1∏mki=1otherwise(ki>1)
可知有平方因子的数的莫比乌斯函数值为 0。
- 由于Mobius函数是积性函数,即
μ
(
a
⋅
b
)
=
μ
(
a
)
⋅
μ
(
b
)
,
(
a
,
b
)
=
1
\mu(a \cdot b)=\mu(a)\cdot \mu(b), (a,b)=1
μ(a⋅b)=μ(a)⋅μ(b),(a,b)=1,所以
可以在 O ( n ) O(n) O(n)的时间内筛出 [ 1 , n ] [1,n] [1,n]的函数值,一个数 x = p 1 k 1 ⋅ p 2 k 2 ⋅ ⋯ ⋅ p m k m x= p_1 ^ {k_1} \cdot p_2 ^ {k_2} \cdot\dots\cdot p_m ^ {k_m} x=p1k1⋅p2k2⋅⋯⋅pmkm,
若 k i = 1 , k j > 1 , j > i k_i=1 , k_j>1, j>i ki=1,kj>1,j>i 则x只会被 x p j m o d    p j = 0 \frac x {p_j} \mod p_j=0 pjxmodpj=0 删除
int mu[MAXN],prime[MAXN],cnt;
bool noprime[MAXN];
void mublus(){
mu[1]=1;
for(int i=2;i<MAXN;++i){
if(!noprime[i]){
prime[++cnt]=i;
mu[i]=-1;//质数为-1
}
for(int j=1;j<=cnt&&prime[j]*i<MAXN;++j){
noprime[prime[j]*i]=1;
if(i%prime[j]==0){
mu[prime[j]*i]=0;// prime[j]*i有因子prime[j]^2 ,故为0;
break;
}
else mu[prime[j]*i]=-mu[i];// i 比prime[j]*i多了一个质因子。
}
}
}
- ∑ d ∣ n μ ( d ) = [ n = 1 ] \sum_{d\mid n}\mu(d) = [n = 1] ∑d∣nμ(d)=[n=1]
证明:当 n = 1 时有 μ ( 1 ) = 1 \mu(1) = 1 μ(1)=1,显然成立。否则设 n = p 1 k 1 ⋅ p 2 k 2 ⋅ ⋯ ⋅ p m k m n = p_1 ^ {k_1} \cdot p_2 ^ {k_2}\cdot\cdots\cdot p_m ^ {k_m} n=p1k1⋅p2k2⋅⋯⋅pmkm , d = p 1 x 1 ⋅ p 2 x 2 ⋅ ⋯ ⋅ p m x m d = p_1 ^ {x_1} \cdot p_2 ^ {x_2} \cdot\cdots\cdot p_m ^ {x_m} d=p1x1⋅p2x2⋅⋯⋅pmxm
根据 μ ( n ) \mu(n) μ(n) 的定义,只需考虑 x i = 0 x_i = 0 xi=0 或 x i = 1 x_i = 1 xi=1的情况。我们设 d 中存在 r 个 x i x_i xi为 1,那么有:
∑ d ∣ n μ ( d ) = ∑ r = 0 m ( m r ) ( − 1 ) r ( n ≠ 1 ) \sum_{d\mid n}\mu(d) = \sum_{r = 0}^m\binom m r(-1)^r(n \neq 1) ∑d∣nμ(d)=∑r=0m(rm)(−1)r(n̸=1)
根据二项式定理:
(
x
+
y
)
n
=
∑
k
=
0
n
(
n
k
)
x
n
−
k
y
k
(
x
+
y
)
(x+y)^{n}=\sum _{k=0}^{n}\binom n kx^{n-k}y^{k} (x+y)
(x+y)n=∑k=0n(kn)xn−kyk(x+y)
我们令 x = 1, y = -1,即得证:
∑ r = 0 m ( m r ) ( − 1 ) r = ( 1 − 1 ) m = 0 ( n ≠ 1 ) \sum_{r = 0}^m\binom m r(-1)^r = (1 - 1)^m = 0(n \neq 1) ∑r=0m(rm)(−1)r=(1−1)m=0(n̸=1)
莫比乌斯反演
假设对于数论函数 f ( n ) f(n) f(n) 和 F ( n ) F(n) F(n),满足以下关系式: F ( n ) = ∑ d ∣ n f ( d ) F(n)=\sum _{{d|n}}f(d) F(n)=∑d∣nf(d)
则将其默比乌斯反演公式定义为:
f
(
n
)
=
∑
d
∣
n
μ
(
d
)
F
(
n
d
)
f(n)=\sum _{d|n}\mu (d) F({\frac {n}{d}})
f(n)=d∣n∑μ(d)F(dn)
或者满足:
F
(
n
)
=
∑
n
∣
d
f
(
d
)
F(n)=\sum_{n|d}f(d)
F(n)=∑n∣df(d)
有
f
(
n
)
=
∑
n
∣
d
μ
(
d
n
)
F
(
d
)
f(n)=\sum_{n|d}\mu(\frac d n)F(d)
f(n)=n∣d∑μ(nd)F(d)
1.求区间[1,a], [1,b]中gcd=k的个数
∑
i
=
1
a
∑
j
=
1
b
[
g
c
d
(
i
,
j
)
=
k
]
\sum_{i=1}^{a}\sum_{j=1}^{b}[gcd(i,j)=k]
i=1∑aj=1∑b[gcd(i,j)=k]
即求
∑
i
=
1
⌊
a
k
⌋
∑
j
=
1
⌊
b
k
⌋
[
g
c
d
(
i
,
j
)
=
1
]
\sum_{i=1}^{\lfloor\frac a k\rfloor}\sum_{j=1}^{\lfloor\frac b k\rfloor}[gcd(i,j)=1]
i=1∑⌊ka⌋j=1∑⌊kb⌋[gcd(i,j)=1]
令 f ( k ) f(k) f(k)为满足gcd=k的个数, F ( k ) F(k) F(k)为满足 g c d = x ⋅ k , x > 0 gcd=x \cdot k,x>0 gcd=x⋅k,x>0 的个数.
则有 F ( n ) = ∑ n ∣ d f ( d ) F(n)=\sum_{n|d}f(d) F(n)=∑n∣df(d) (d是n的倍数, ∑ f ( d ) \sum f(d) ∑f(d)显然就是F的定义)
所以有反演公式: f ( n ) = ∑ n ∣ d μ ( d n ) F ( d ) f(n)=\sum_{n|d}\mu(\frac d n)F(d) f(n)=∑n∣dμ(nd)F(d)
显然, F ( n ) = ⌊ a n ⌋ ⋅ ⌊ b n ⌋ F(n)=\lfloor\frac a n\rfloor\cdot\lfloor\frac b n\rfloor F(n)=⌊na⌋⋅⌊nb⌋
(a/n是第一个区间里n的倍数的个数,两者相乘则为 i,j 都是n的倍数的个数,所以gcd(i,j)是n的倍数,显然就是F的定义)
所以 f ( n ) = ∑ n ∣ d μ ( d n ) F ( d ) = ∑ n ∣ d μ ( d n ) ⋅ ⌊ a d ⌋ ⋅ ⌊ b d ⌋ f(n)=\sum_{n|d}\mu(\frac d n)F(d)=\sum_{n|d}\mu(\frac d n)\cdot\lfloor\frac a d\rfloor\cdot\lfloor\frac b d\rfloor f(n)=∑n∣dμ(nd)F(d)=∑n∣dμ(nd)⋅⌊da⌋⋅⌊db⌋
由于k=1
有 f ( 1 ) = ∑ 1 ∣ d μ ( d ) F ( d ) = ∑ 1 ∣ d μ ( d ) ⋅ ⌊ a d ⌋ ⋅ ⌊ b d ⌋ f(1)=\sum_{1|d}\mu(d)F(d)=\sum_{1|d}\mu(d)\cdot\lfloor\frac a d\rfloor\cdot\lfloor\frac b d\rfloor f(1)=∑1∣dμ(d)F(d)=∑1∣dμ(d)⋅⌊da⌋⋅⌊db⌋
求[a,b],[c,d]里gcd=k的个数,a=c=1,gcd(1,2)和gcd(2,1)算作一个。
a
n
s
=
f
(
1
)
=
∑
1
∣
i
min
(
b
,
d
)
μ
(
i
)
F
(
i
)
=
∑
1
∣
i
μ
(
i
)
⋅
⌊
b
i
⌋
⋅
⌊
d
i
⌋
ans=f(1)=\sum_{1|i}^{\min(b,d)}\mu(i)F(i)=\sum_{1|i}\mu(i)\cdot\lfloor\frac b i\rfloor\cdot\lfloor\frac d i\rfloor
ans=f(1)=∑1∣imin(b,d)μ(i)F(i)=∑1∣iμ(i)⋅⌊ib⌋⋅⌊id⌋
为何sum的上界是min(b,d),因为i大于min时,F(i)=0,且F(i)表示两个区间里各取一个数,满足gcd=i的倍数,若i大于其中一个区间的上界,则不能满足gcd=i的倍数。
最后,ans是有重复的,例如[1,5]和[1,10]两个区间,[1,3]和[3,1]算作两次,但是[1,6],[6,1]中[6,1]不可能倍计算,所以只算了一次。所以,在[1,min(b,d)]里的答案都计算了两次。
#include<bits/stdc++.h>
using namespace std;
#define Init(arr,val) memset(arr,val,sizeof(arr))
const int inf=0x3f3f3f3f,mod=1000000007,MAXN=1e5+8;
typedef long long ll;
int mu[MAXN],prime[MAXN],cnt;
bool noprime[MAXN];
void mubius(){
mu[1]=1;
for(int i=2;i<MAXN;++i){
if(!noprime[i]){
prime[++cnt]=i;
mu[i]=-1;
}
for(int j=1;j<=cnt&&prime[j]*i<MAXN;++j){
noprime[prime[j]*i]=1;
if(i%prime[j]==0){
mu[prime[j]*i]=0;
break;
}
else mu[prime[j]*i]=-mu[i];
}
}
}
int main() {
std::ios::sync_with_stdio(0);
std::cin.tie(0);
std::cout.tie(0);
mubius();
int T,a,b,c,d,k;
cin>>T;
for(int Case=1;Case<=T;++Case){
cin>>a>>b>>c>>d>>k;
if(k==0){//没有gcd=0
cout<<"Case "<<Case<<": 0\n";
continue;
}
b/=k,d/=k;//将gcd=k化成gcd=1.
ll ans1=0,ans2=0;
if(b>d)swap(b,d);//b,d里b是最小值。
for(int i=1;i<=b;++i){
ans1+=1ll*mu[i]*(b/i)*(d/i);
ans2+=1ll*mu[i]*(b/i)*(b/i);
}
cout<<"Case "<<Case<<": "<<ans1-(ans2>>1)<<endl;
}
return 0;
}
这样做是 O(n)。
例如,b=10000,d=11000,i>6000时,
⌊
b
i
⌋
=
⌊
d
i
⌋
=
1
\lfloor\frac b i\rfloor=\lfloor\frac d i\rfloor=1
⌊ib⌋=⌊id⌋=1,但是我们还是计算了4000次。
所以,不如先求 μ \mu μ 的前缀和,再求出对于b,d两数, ⌊ b i ⌋ , ⌊ d i ⌋ \lfloor\frac b i\rfloor, \lfloor\frac d i\rfloor ⌊ib⌋,⌊id⌋均不变的范围,在该范围内,
a n s = s u m ( μ ) ⋅ ⌊ b i ⌋ ⋅ ⌊ d i ⌋ ans=sum(\mu)\cdot\lfloor\frac b i\rfloor\cdot\lfloor\frac d i\rfloor ans=sum(μ)⋅⌊ib⌋⋅⌊id⌋
对于一个在区间[1 ,x]内的位置 i,设 ⌊ x i ⌋ = k \lfloor\frac x i\rfloor=k ⌊ix⌋=k ,对于向下取整都为k的最后一个位置last,
⌊ x l a s t ⌋ = k \lfloor\frac x {last}\rfloor=k ⌊lastx⌋=k,且 ⌊ x l a s t + 1 ⌋ = k − 1 , x m o d    ( l a s t + 1 ) = 0 \lfloor\frac {x} {last+1}\rfloor=k-1,x\mod ({last+1})=0 ⌊last+1x⌋=k−1,xmod(last+1)=0
得: l a s t = x ⌊ x i ⌋ − 1 − 1 last=\frac{x}{\lfloor\frac x i\rfloor -1}-1 last=⌊ix⌋−1x−1,因为x和 ⌊ x i ⌋ − 1 \lfloor\frac x i\rfloor -1 ⌊ix⌋−1是整除关系,把-1去掉,分子相当于+1,所以 x ⌊ x i ⌋ − 1 \frac{x}{\lfloor\frac x i\rfloor }-1 ⌊ix⌋x−1不是整除,
比原来少了1,再加上1就和原来一样了,所以为了方便,减1全部去掉。
最后再和另一个区间里的取最小值。
for(int i=1;i<=b;++i){
ans1+=1ll*mu[i]*(b/i)*(d/i);
ans2+=1ll*mu[i]*(b/i)*(b/i);
}
//由上面变成了下面,只需处理一下前缀和:
//for(int i=1;i<MAXN;++i)sumu[i]=sumu[i-1]+mu[i];
for(int last,i=1;i<=b;i=last+1){
last=min(b/(b/i),d/(d/i));
ans1+=1ll*(sumu[last]-sumu[i-1])*(b/i)*(d/i);
ans2+=1ll*(sumu[last]-sumu[i-1])*(b/i)*(b/i);
}
复杂度: O ( n + m ) O(\sqrt n + \sqrt m) O(n+m)
P2257 YY的GCD
题意:给定N, M,求1<=x<=N, 1<=y<=M且gcd(x, y)为质数的(x, y)有多少对。
还是一样,设 f ( d ) = ∑ i = 1 N ∑ j = 1 M [ g c d ( i , j ) = d ] f(d)=\sum_{i=1}^N \sum_{j=1}^M[gcd(i,j)=d] f(d)=∑i=1N∑j=1M[gcd(i,j)=d] 假设 N ≤ M N\le M N≤M
F ( n ) = ∑ i = 1 N ∑ j = 1 M [ n ∣ g c d ( i , j ) ] = ⌊ N n ⌋ ⌊ M n ⌋ F(n)=\sum_{i=1}^N \sum_{j=1}^M[n|gcd(i,j)]=\lfloor\frac N n\rfloor\lfloor\frac M n\rfloor F(n)=∑i=1N∑j=1M[n∣gcd(i,j)]=⌊nN⌋⌊nM⌋
由反演公式得到:
f ( n ) = ∑ n ∣ d μ ( d n ) F ( d ) = ∑ n ∣ d μ ( d n ) ⌊ N d ⌋ ⌊ M d ⌋ f(n)=\sum_{n|d}\mu(\frac d n)F(d)=\sum_{n|d}\mu(\frac d n)\lfloor\frac N d\rfloor\lfloor\frac M d\rfloor f(n)=∑n∣dμ(nd)F(d)=∑n∣dμ(nd)⌊dN⌋⌊dM⌋
a n s = ∑ p ∈ p r i m e f ( p ) = ∑ p ∈ p r i m e ∑ p ∣ d μ ( d p ) ⌊ N d ⌋ ⌊ M d ⌋ ans=\sum_{p\in prime}f(p)=\sum_{p\in prime}\sum_{p|d}\mu(\frac d p)\lfloor\frac N d\rfloor\lfloor\frac M d\rfloor ans=∑p∈primef(p)=∑p∈prime∑p∣dμ(pd)⌊dN⌋⌊dM⌋
= ∑ p ∈ p r i m e ∑ k = 1 μ ( k ) ⌊ N d p ⌋ ⌊ M k p ⌋ , d : = k p =\sum_{p\in prime}\sum_{k=1}\mu(k)\lfloor\frac N {dp}\rfloor\lfloor\frac M {kp}\rfloor, d:=kp =∑p∈prime∑k=1μ(k)⌊dpN⌋⌊kpM⌋,d:=kp
令 T = k p T=kp T=kp,有:
a n s = ∑ p ∈ p r i m e ∑ k = 1 μ ( k ) ⌊ N T ⌋ ⌊ M T ⌋ = ∑ p ∈ p r i m e ∑ p ∣ T μ ( T p ) ⌊ N T ⌋ ⌊ M T ⌋ ans=\sum_{p\in prime}\sum_{k=1}\mu(k)\lfloor\frac N {T}\rfloor\lfloor\frac M {T}\rfloor=\sum_{p\in prime}\sum_{p|T}\mu(\frac T p)\lfloor\frac N {T}\rfloor\lfloor\frac M {T}\rfloor ans=∑p∈prime∑k=1μ(k)⌊TN⌋⌊TM⌋=∑p∈prime∑p∣Tμ(pT)⌊TN⌋⌊TM⌋
重点来了,如何交换两个求和符号的顺序?画图大法好。
p
=
p
1
p=p_1
p=p1时,
T
1
=
p
1
,
2
p
1
,
3
p
1
,
…
,
n
1
p
1
T_1=p_1,2p_1,3p_1,\dots,n_1p_1
T1=p1,2p1,3p1,…,n1p1
p
=
p
2
p=p_2
p=p2时,
T
2
=
p
2
,
2
p
2
,
3
p
2
,
…
,
n
2
p
2
T_2=p_2,2p_2,3p_2,\dots,n_2p_2
T2=p2,2p2,3p2,…,n2p2
…
\dots
…
p
=
p
n
p=p_n
p=pn时,
T
n
=
p
n
,
2
p
n
,
3
p
n
,
…
,
n
p
n
T_n=p_n,2p_n,3p_n,\dots,np_n
Tn=pn,2pn,3pn,…,npn
其中
T
i
≤
N
T_i\le N
Ti≤N
所以,T从最小的质数
p
1
p_1
p1开始枚举,由于任意合数
x
x
x,都能分解为
x
=
k
i
p
i
=
k
j
p
j
x=k_ip_i=k_jp_j
x=kipi=kjpj,所以T遍历
p
1
p_1
p1开始到x的
所有数, T p \frac T p pT枚举的是 k i … k_i \dots ki…,所以有:
a n s = ∑ p ∈ p r i m e ∑ p ∣ T μ ( T p ) ⌊ N T ⌋ ⌊ M T ⌋ ans=\sum_{p\in prime}\sum_{p|T}\mu(\frac T p)\lfloor\frac N {T}\rfloor\lfloor\frac M {T}\rfloor ans=∑p∈prime∑p∣Tμ(pT)⌊TN⌋⌊TM⌋
= ∑ T = 1 x ∑ p ∣ T , p ∈ p r i m e μ ( T p ) ⌊ N T ⌋ ⌊ M T ⌋ =\sum_{T=1}^{x}\sum_{p|T, p\in prime}\mu(\frac T p)\lfloor\frac N {T}\rfloor\lfloor\frac M {T}\rfloor =∑T=1x∑p∣T,p∈primeμ(pT)⌊TN⌋⌊TM⌋
= ∑ T = 1 x ⌊ N T ⌋ ⌊ M T ⌋ ∑ p ∣ T , p ∈ p r i m e μ ( T p ) =\sum_{T=1}^{x}\lfloor\frac N {T}\rfloor\lfloor\frac M {T}\rfloor\sum_{p|T, p\in prime}\mu(\frac T p) =∑T=1x⌊TN⌋⌊TM⌋∑p∣T,p∈primeμ(pT)(T从2开始也是一样的,T=1时后面的求和等于0)
右边的求和可以通过预处理获得前缀和,对于 ∑ p ∣ T , p ∈ p r i m e μ ( T p ) \sum_{p|T, p\in prime}\mu(\frac T p) ∑p∣T,p∈primeμ(pT),我们考虑每个质数p对T=kp的贡献。
#include<bits/stdc++.h>
using namespace std;
#define Init(arr,val) memset(arr,val,sizeof(arr))
const int inf=0x3f3f3f3f,mod=10007,MAXN=1e7+8;
typedef long long ll;
int mu[MAXN],prime[MAXN],cnt;
bool noprime[MAXN];
int sumu[MAXN];
void Mobius() {
mu[1]=1;
for(int i=2; i<MAXN; ++i) {
if(!noprime[i]) {
prime[++cnt]=i;
mu[i]=-1;
}
for(int j=1; j<=cnt&&prime[j]*i<MAXN; ++j) {
noprime[prime[j]*i]=1;
if(i%prime[j]==0) {
mu[i*prime[j]]=0;
break;
} else
mu[i*prime[j]]=-mu[i];
}
}
for(int i=1; i<=cnt; ++i)
for(int j=1; j*prime[i]<MAXN; ++j)
sumu[j*prime[i]]+=mu[j];
for(int i=1; i<MAXN; ++i)
sumu[i]+=sumu[i-1];
}
int main() {
std::ios::sync_with_stdio(0);
std::cin.tie(0);
std::cout.tie(0);
Mobius();
int t;
cin>>t;
while(t--){
int x,y,last;
cin>>x>>y;
ll ans=0;
if(x>y)swap(x,y);
for(int i=1;i<=x;i=last+1){
last=min(x/(x/i),y/(y/i));
ans+=1ll*(x/last)*(y/last)*(sumu[last]-sumu[i-1]);
}
cout<<ans<<endl;
}
return 0;
}