莫比乌斯反演练习

HDU - 1695 : GCD

CodeForces - 235E : Number Challenge

BZOJ 2154: Crash的数字表格

BZOJ 3529: [Sdoi2014]数表

BZOJ 2440: [中山市选2011]完全平方数

BZOJ 2301: [HAOI2011]Problem b

HDU - 4746 : Mophues

HDU - 5212 : Code

HDU - 5321 : Beautiful Set

HDU - 4947 : GCD Array

HDU - 5663 : Hillan and the girl

UESTC - 811 : GCD

HDU - 5656 : CA Loves GCD

ACdream - 1114 : Number theory

HDU - 5468 : Puzzled Elena

HDU - 2841 : Visible Trees

POJ - 3904 : Sky Code

莫比乌斯反演式:

$\sum_{d|n} u(d) = [n=1]$ ,其中$d = \prod_{i=1}^{k} p_i^{e_i}$ ,$u(d)$ 为莫比乌斯函数,见下式
$$
u(d)=\left{
\begin{array}{rcl}
1 &{d=1}\\
0 &{\exists e_i > 1}\\
(-1)^k &{\forall e_i=1}\\
\end{array} \right.
$$

HDU - 1695 : GCD

题目大意:$T$ 次询问,求$1 \leqslant x \leqslant b,1 \leqslant y \leqslant d$ 的所有$gcd(x,y)=k$ 的个数,$(x_0,y_0)$ 与$(y_0,x_0)$ 被视作一样。

莫比乌斯反演

设$g(k)=\sum_{i=1}^b\sum_{j=1}^d[gcd(i,j)=k],f(n)=\sum_{n|k}g(k)$

$f(n)$ 的意义为$gcd(i,j)$ 能被$n$ 整除的且符合条件的点对个数

设$n\leqslant m$ ,则$f(n)= \frac{ \lfloor \frac{b}{n} \rfloor (\lfloor \frac{b}{n} \rfloor+1)}{2}+ \lfloor \frac{b}{n} \rfloor(\lfloor \frac{d}{n} \rfloor-\lfloor \frac{b}{n} \rfloor)$ 。

故$g(n)=\sum_{n|k}u(\frac{k}{n})f(k)$ 代入即可求得。

复杂度$O(Tn)$ 。

代码如下:

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
39
#include <iostream>
using namespace std;
typedef long long ll;
const int N=100005;
int p[N],vis[N]={1,1},tot,u[N]={0,1};
void init(int n){
for(int i=2;i<n;++i){
if(!vis[i]){
p[tot++]=i;
u[i]=-1;
}
for(int j=0;j<tot&&i*p[j]<n;++j){
vis[i*p[j]]=1;
if(i%p[j]==0){
u[i*p[j]]=0;
break;
}else u[i*p[j]]=-u[i];
}
}
}
ll T,a,b,c,d,k;
int main(void){
init(N);
cin>>T;
for(int t=1;t<=T;++t){
cin>>a>>b>>c>>d>>k;
ll ans=0;
if(!k){
cout<<"Case "<<t<<": "<<ans<<"\n";
continue;
}
if(b>d)swap(b,d);
for(int i=k;i<=b;i+=k){
ll x=b/i,y=d/i;
ans+=u[i/k]*(x*(x+1)/2+x*(y-x));
}
cout<<"Case "<<t<<": "<<ans<<"\n";
}
}

CodeForces - 235E : Number Challenge

待更

BZOJ 2154: Crash的数字表格

题目大意:求$\sum_{i=1}^n\sum_{j=1}^mlcm(i,j)$ 。

莫比乌斯反演

$$
f(n,m)=\sum_{i=1}^n\sum_{j=1}^m[i,j]
=\sum_{i=1}^n\sum_{j=1}^m\frac{i\times j}{(i,j)}\\
=\sum_{d=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m\frac{i\times j}{(i,j)}[(i,j)=d]\\
=\sum_{d=1}^{min(n,m)}\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j\times d[(i,j)=1]\\
=\sum_{d=1}^{min(n,m)}d\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}i\times j[(i,j)=1]\\
$$

设$g(n,m)=\sum_{i=1}^n\sum_{j=1}^mi\times j[(i,j)=1]$ ,则有
$$
g(n,m)=\sum_{d=1}^{min(n,m)}u(d)\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}i\times j\times d^2\\
=\sum_{d=1}^{min(n,m)}u(d)d^2\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}i\times j\\
=\sum_{d=1}^{min(n,m)}u(d)d^2\times \frac{(\lfloor \frac{n}{d} \rfloor+1)\lfloor \frac{n}{d} \rfloor}{2}\times \frac{(\lfloor \frac{m}{d} \rfloor+1)\lfloor \frac{m}{d} \rfloor}{2}\\
$$

$$
\therefore f(n,m)=\sum_{d=1}^{min(n,m)}d\times g(\lfloor \frac{n}{d}\rfloor,\lfloor \frac{m}{d}\rfloor)
$$

求$f(n,m)$ 和$g(n,m)$ 的复杂度均为$O(\sqrt{n})$ ,总复杂度为$O(n)$ 。

代码如下:

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
39
40
41
42
43
44
45
46
47
48
49
50
51
#include <cstdio>
using namespace std;
typedef long long ll;
const int N=1e7+5;
const int P=20101009;
int p[N],k,pre[N];
char mu[N]={0,1};
bool vis[N]={1,1};
void init(int n){
for(int i=2;i<n;++i){
if(!vis[i]){
mu[i]=-1;
p[k++]=i;
}
for(int j=0;j<k&&p[j]*i<n;++j){
vis[i*p[j]]=1;
if(i%p[j]==0){
mu[i*p[j]]=0;
break;
}else mu[i*p[j]]=-mu[i];
}
}
for(int i=1;i<n;++i)
pre[i]=((pre[i-1]+(ll)mu[i]*i*i%P)%P+P)%P;
}
int n,m;
int Min(int a,int b){return a<b?a:b;}
int Max(int a,int b){return a>b?a:b;}
int g(int n,int m){
int ans=0,d,td,len=Min(n,m);
for(d=1;d<=len;d=td+1){
td=Min(n/(n/d),m/(m/d));
ll temp=(((pre[td]-pre[d-1])%P*((1LL*(n/d+1)*(n/d)/2)%P))%P*((1LL*(m/d+1)*(m/d)/2)%P))%P;
ans=((ans+temp)%P+P)%P;
}
return ans;
}
int f(int n,int m){
int ans=0,d,td,len=Min(n,m);
for(d=1;d<=len;d=td+1){
td=Min(n/(n/d),m/(m/d));
ll temp=1LL*(td+d)*(td-d+1)/2%P*g(n/d,m/d)%P;
ans=((ans+temp)%P+P)%P;
}
return ans;
}
int main(void){
init(N);
scanf("%d%d",&n,&m);
printf("%d\n",f(n,m));
}

BZOJ 3529: [Sdoi2014]数表

待更

BZOJ 2440: [中山市选2011]完全平方数

题目大意:求第$K$ 个非平方数的倍数的数。

二分+容斥+莫比乌斯函数

做这题前学习了下容斥原理:

$|\bar{A_1} \cap \bar{A_2}\cap … \cap \bar{A_n}|=|S|-\sum|A_i|+\sum|A_i\cap A_j|+…+(-1)^n\sum|A_1\cap A_2 \cap … \cap A_n|$

设$f(n)$ 是$n$ 个素数平方的倍数的数的个数,那么$f(0)=|S|-f(1)+f(2)+…+(-1)^nf(n)$

不难看出$f(x)$ 前的系数即为莫比乌斯函数$u(x)$ 。

于是二分第$K$ 个数$x$ ,检验$f(x)$ 是否等于$K$ 即可。

代码如下:

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
39
40
#include <cstdio>
typedef long long ll;
const int N=100005;
ll p[N],vis[N]={1,1},k,mu[N]={0,1};
void init(ll n){
for(ll i=2;i<n;++i){
if(!vis[i]){
p[k++]=i;
mu[i]=-1;
}
for(ll j=0;j<k&&i*p[j]<n;++j){
vis[i*p[j]]=1;
if(i%p[j]==0){
mu[i*p[j]]=0;
break;
}mu[i*p[j]]=-mu[i];
}
}
}
ll check(ll n){
ll ans=n;
for(ll i=2;i*i<=n;++i)
ans+=mu[i]*(n/(i*i));
return ans;
}
int main(void){
init(N);
ll T,n,l,r,mid;
scanf("%lld",&T);
while(T--){
scanf("%lld",&n);
l=1,r=10000000000;
while(l<r){
mid=(l+r)>>1;
if(check(mid)<n)l=mid+1;
else r=mid;
}
printf("%lld\n",l);
}
}

BZOJ 2301: [HAOI2011]Problem b

题目大意:$T$ 次询问,每次询问有多少数对$(x,y)$ 满足$a\leqslant x \leqslant b,c\leqslant y \leqslant d$ 且$gcd(x,y)=k$

莫比乌斯反演

和上面HDU 1695类似,可求$1 \leqslant x \leqslant n,1 \leqslant y \leqslant m$ 且$gcd(x,y)=k$ 的点对有$F(n,m,k)$ 个。

那么简单地容斥下,$a\leqslant x \leqslant b,c\leqslant y \leqslant d$ 且$gcd(x,y)=k$ 的点对就有$F(b,d,k)-F(a-1,d,k)-F(b,c-1,k)+F(a-1,c-1,k)$ 个。

然而这道题的数据范围更大,用之前的方法直接暴力求不可行,采用分块优化。

$\lfloor \frac{N}{d} \rfloor$ 的取值并不是连续的,其最多只有$2\lfloor \sqrt{N} \rfloor$ 种,同理$\lfloor \frac{M}{d} \rfloor$ 的取值最多$2\lfloor \sqrt{M} \rfloor$ 种。故$\lfloor \frac{N}{d} \rfloor$ 和$\lfloor \frac{M}{d} \rfloor$ 均相同的区间只有$2\lfloor \sqrt{N} \rfloor+2\lfloor \sqrt{M} \rfloor$ 个,所以$d$ 的枚举量其实是$O(\sqrt{N}+\sqrt{M})$ 的。

复杂度$O(T\sqrt{n})$ 。

代码如下:

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
39
40
41
42
43
44
#include <cstdio>
using namespace std;
typedef long long ll;
const int N=50005;
ll p[N],vis[N]={1,1},tot,u[N]={0,1},pre[N];
void init(int n){
for(int i=2;i<n;++i){
if(!vis[i]){
p[tot++]=i;
u[i]=-1;
}
for(int j=0;j<tot&&i*p[j]<n;++j){
vis[i*p[j]]=1;
if(i%p[j]==0){
u[i*p[j]]=0;
break;
}else u[i*p[j]]=-u[i];
}
}
for(int i=1;i<n;++i)
pre[i]=pre[i-1]+u[i];
}
ll Min(ll a,ll b){
return a<b?a:b;
}
ll solve(ll a,ll b,ll k){
if(!k)return 0;
ll ans=0,l;
for(ll i=1;i<=Min(a,b);i=l+1){
l=Min(a/(a/i),b/(b/i));
ans+=(pre[l]-pre[i-1])*(a/i)*(b/i);
}
return ans;
}
ll T,a,b,c,d,k;
int main(void){
init(N);
scanf("%lld",&T);
while(T--){
scanf("%lld%lld%lld%lld%lld",&a,&b,&c,&d,&k);
ll ans=solve(b/k,d/k,k)-solve((a-1)/k,d/k,k)-solve(b/k,(c-1)/k,k)+solve((a-1)/k,(c-1)/k,k);
printf("%lld\n",ans);
}
}

HDU - 4746 : Mophues

待更

HDU - 5212 : Code

题目大意:求$f(n)=\sum_{i=1}^n\sum_{j=1}^n(a_i,a_j)[(a_i,a_j)-1]$

莫比乌斯反演

$$
f(n)=\sum_{i=1}^n\sum_{j=1}^n(a_i,a_j)[(a_i,a_j)-1]\\
=\sum_{d=1}^nd(d-1)\sum_{i=1}^n\sum_{j=1}^n[(a_i,a_j)=d]\\
设g(n,m)=\sum_{i=1}^n\sum_{j=1}^n[(a_i,a_j)=m]\\
=\sum_{i=1}^n\sum_{j=1}^n[\frac{(a_i,a_j)}{m}=1]\\
=\sum_{i=1}^n\sum_{j=1}^n\sum_{d|\frac{(a_i,a_j)}{m}}u(d)
$$

设$h(d)$ 表示能被$d$ 整除的数的个数,则
$$
g(n,m)=\sum_{d=1}^{\frac{max(a_i)}{m}}u(d)h^2(md)\\
f(n)=\sum_{d=1}^nd(d-1)\times g(n,d)
$$
复杂度$O(nlogn)$ 。

代码如下:

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
39
40
41
42
43
44
45
46
47
48
49
50
#include <cstring>
#include <cstdio>
using namespace std;
typedef long long ll;
const int N=10005;
const int P=10007;
ll p[N],vis[N]={1,1},k,mu[N]={0,1};
ll n,a[N],h[N],t;
void init(ll n){
for(ll i=2;i<n;++i){
if(!vis[i]){
mu[i]=-1;
p[k++]=i;
}
for(ll j=0;j<k&&p[j]*i<n;++j){
vis[i*p[j]]=1;
if(i%p[j]==0){
mu[i*p[j]]=0;
break;
}else mu[i*p[j]]=-mu[i];
}
}
}
ll g(ll n,ll m){
ll ans=0;
for(ll d=1;d<=10000/m;++d)
ans=(ans+mu[d]*h[m*d]*h[m*d])%P;
return ans;
}
ll f(ll n){
ll ans=0;
for(ll d=1;d<=10000;++d)
ans=(ans+d*(d-1)*g(n,d))%P;
return ans;
}
int main(void){
init(N);
while(~scanf("%I64d",&n)){
memset(a,0,sizeof(a));
memset(h,0,sizeof(h));
for(ll i=0;i<n;++i){
scanf("%I64d",&t);
a[t]++;
}
for(ll i=1;i<=10000;++i)
for(ll j=i;j<=10000;j+=i)
h[i]=(h[i]+a[j])%P;
printf("%I64d\n",(f(n)+P)%P);
}
}

HDU - 5321 : Beautiful Set

待更

HDU - 4947 : GCD Array

待更

HDU - 5663 : Hillan and the girl

待更

UESTC - 811 : GCD

待更

HDU - 5656 : CA Loves GCD

待更

ACdream - 1114 : Number theory

题目大意:给出$n$ 个数,问满足$gcd(a_i,a_j)=1$ 且$i<j$ 的点对有多少个。

莫比乌斯反演

$$
f(n)=\sum_{i=1}^n\sum_{j=1}^n[(a_i,a_j)=1]\\
f(n)=\sum_{i=1}^n\sum_{j=1}^n\sum_{d|(a_i,a_j)}u(d)
$$

设$g(d)$ 表示能被$d$ 整除的数的个数,则有
$$
f(n)=\sum_{d=1}^{max(a_i)}u(d)g^2(d)
$$
复杂度$O(nlogn)$ 。

代码如下:

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
39
40
41
42
43
44
45
46
47
48
49
50
#include <cstdio>
#include <cstring>
using namespace std;
typedef long long ll;
const int N=222222+5;
ll p[N],vis[N]={1,1},k,mu[N]={0,1};
void init(ll n){
for(ll i=2;i<n;++i){
if(!vis[i]){
p[k++]=i;
mu[i]=-1;
}
for(ll j=0;j<k&&p[j]*i<n;++j){
vis[i*p[j]]=1;
if(i%p[j]==0){
mu[i*p[j]]=0;
break;
}else mu[i*p[j]]=-mu[i];
}
}
}
ll n,a[N],d[N],t;
ll Max(ll a,ll b){
return a>b?a:b;
}
ll solve(ll n){
ll ans=0;
for(ll i=1;i<=n;++i)
ans+=mu[i]*d[i];
return ans;
}
int main(void){
init(N);
while(~scanf("%lld",&n)){
memset(a,0,sizeof(a));
memset(d,0,sizeof(d));
ll mx=0;
for(ll i=0;i<n;++i){
scanf("%lld",&t);
mx=Max(mx,t);
a[t]++;
}
for(ll i=1;i<=mx;++i)
for(ll j=i;j<=mx;j+=i)
d[i]+=a[j];
for(ll i=1;i<=mx;++i)
d[i]=d[i]*(d[i]-1)/2;
printf("%lld\n",solve(mx));
}
}

HDU - 5468 : Puzzled Elena

待更

HDU - 2841 : Visible Trees

题目大意:问从$(0,0)$ 点看第一象限$n \times m$ 范围,最多能看到几个点。

莫比乌斯反演

$$
f(n,m)=\sum_{i=1}^n\sum_{j=1}^m[(i,j=1)]\\
f(n,m)=\sum_{i=1}^n\sum_{j=1}^m\sum_{d|(i,j)}u(d)\\
f(n,m)=\sum_{d=1}^{min(n,m)}u(d)\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}1\\
$$

分块优化后复杂度为$O(\sqrt{n})$ 。

代码如下:

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
39
40
#include <iostream>
#include <cstdio>
using namespace std;
typedef long long ll;
const int N=100005;
ll T,n,m,p[N],vis[N]={1,1},k,mu[N]={0,1},pu[N];
void init(int n){
for(int i=2;i<n;++i){
if(!vis[i]){
p[k++]=i;
mu[i]=-1;
}
for(int j=0;j<k&&i*p[j]<n;++j){
vis[i*p[j]]=1;
if(i%p[j]==0){
mu[i*p[j]]=0;
break;
}else mu[i*p[j]]=-mu[i];
}
}
for(int i=1;i<n;++i)
pu[i]=pu[i-1]+mu[i];
}
ll solve(ll n,ll m){
ll ans=0,len=min(n,m);
for(ll d=1;d<=len;){
ll td=min(n/(n/d),m/(m/d));
ans+=(pu[td]-pu[d-1])*(n/d)*(m/d);
d=td+1;
}
return ans;
}
int main(void){
init(N);
scanf("%I64d",&T);
while(T--){
scanf("%I64d%I64d",&n,&m);
printf("%I64d\n",solve(n,m));
}
}

POJ - 3904 : Sky Code

待更

Comments

Your browser is out-of-date!

Update your browser to view this website correctly. Update my browser now

×