SPOJ 20173 DIVCNT2 - Counting Divisors (square)（数论+莫比乌斯反演）

65 阅读 0 评论 43 点赞

我是靠谱客的博主危机猎豹，这篇文章主要介绍SPOJ 20173 DIVCNT2 - Counting Divisors (square)（数论+莫比乌斯反演），现在分享给大家，希望可以做个参考。

Description

给出一整数 n n ，求 $s (n) = \sum_{i = 1}^{n} σ_{0} (i^{2})$ ，其中 σ0(x) σ 0 ( x ) 为 x x 的因子数

Input

第一行一整数 $t$ 表示用例组数，每组用例输入一整数 n(1≤T≤104,1≤n≤1012) n ( 1 ≤ T ≤ 10 4 , 1 ≤ n ≤ 10 12 )

Output

输出 S(n) S ( n )

Sample Input

5
1
2
3
10
100

Sample Output

1
4
7
48
1194

Solution

对 x x 质因子分解有 $x = p_{1}^{a_{1}} . . . p_{k}^{a_{k}}$ ，进而

σ 0 (x 2) = \prod i = 1 k (2 a i + 1) = \sum S \subset {1, . . ., k} 2 | S | \prod i \in S a i = \sum d | x 2 ω (d) σ 0 ( x 2 ) = \prod i = 1 k ( 2 a i + 1 ) = \sum S \subset { 1 , . . . , k } 2 | S | \prod i \in S a i = \sum d | x 2 ω ( d )

其中

ω(d) ω ( d ) 表示

d d 的素因子个数

（累乘 $\prod_{i \in s} a_{i}$ 即为素因子集合为 S S 且每个素因子的指数不超过 $a_{i}$ 的数字个数，也即 x x 的满足素因子集合为 $s$ 的因子个数，而这些数字每个对答案贡献均为 2|S| 2 | S | ，也即为这些数字的素因子个数）

注意到 2ω(d) 2 ω ( d ) 即为 d d 的为无平方因子数的因子个数，故有

2^{ω (d)} = \sum_{j | d} μ^{2} (j)

进而

S (n) = \sum i = 1 n \sum d | i \sum j | d μ 2 (j) = \sum d = 1 n ⌊ n d ⌋ \sum j | d μ 2 (j) = \sum j = 1 μ 2 (j) \sum i = 1 ⌊ n j ⌋ ⌊ n i j ⌋ S ( n ) = \sum i = 1 n \sum d | i \sum j | d μ 2 ( j ) = \sum d = 1 n ⌊ n d ⌋ \sum j | d μ 2 ( j ) = \sum j = 1 μ 2 ( j ) \sum i = 1 ⌊ n j ⌋ ⌊ n i j ⌋

由于

⌊nj⌋ ⌊ n j ⌋ 的取值只有

O(n−−√) O ( n ) 种，对每种取值可以快速计算

μ2(j) μ 2 ( j ) 的区间和以及

∑i=1n⌊ni⌋ ∑ i = 1 n ⌊ n i ⌋ 即可

设 F(n)=∑i=1nμ2(i) F ( n ) = ∑ i = 1 n μ 2 ( i ) ，则有

F (n) = \sum i = 1 n μ 2 (i) = \sum i = 1 n \sum d 2 | i μ (d) = \sum d = 1 n \sqrt μ (d) ⌊ n d 2 ⌋ F ( n ) = \sum i = 1 n μ 2 ( i ) = \sum i = 1 n \sum d 2 | i μ ( d ) = \sum d = 1 n μ ( d ) ⌊ n d 2 ⌋

故可以

O(n−−√) O ( n ) 求出

F(n) F ( n )

设 G(n)=∑i=1n⌊ni⌋=∑i=1nσ0(i) G ( n ) = ∑ i = 1 n ⌊ n i ⌋ = ∑ i = 1 n σ 0 ( i ) ，则 G(n) G ( n ) 为积性函数，满足 G(1)=1,G(p)=2 G ( 1 ) = 1 , G ( p ) = 2 ，且分块后可以 O(n−−√) O ( n ) 求出 G(n) G ( n )

预处理 n23 n 2 3 范围内的 F,G F , G ，对于较大的 n n 可以 $o (\sqrt{n})$ 求出 F(n),G(n) F ( n ) , G ( n )

Code

复制代码

#include<cstdio>
#include<algorithm>
#include<vector>
#include<cmath>
using namespace std;
typedef long long ll;
#define maxn 100000005
int m=100000000;
int mu[maxn],num[maxn],f[maxn],p[maxn/10],res=0;
ll g[maxn];
bool vis[maxn];
void get_prime(int n)
{
mu[1]=1,f[1]=1,g[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i])p[res++]=i,mu[i]=-1,g[i]=2,num[i]=1;
f[i]=f[i-1]+(mu[i]?1:0);
for(int j=0;j<res&&i*p[j]<=n;j++)
{
vis[i*p[j]]=1;
if(i%p[j])
{
mu[i*p[j]]=-mu[i];
g[i*p[j]]=g[i]*2;
num[i*p[j]]=1;
}
else
{
mu[i*p[j]]=0;
num[i*p[j]]=num[i]+1;
g[i*p[j]]=g[i]/(num[i]+1)*(num[i]+2);
break;
}
}
}
for(int i=2;i<=n;i++)g[i]+=g[i-1];
}
ll F(ll n)
{
if(n<=m)return f[n];
int nn=(int)sqrt(n+0.5);
ll ans=0;
for(int d=1;d<=nn;d++)
if(mu[d])ans+=mu[d]*(n/d/d);
return ans;
}
ll G(ll n)
{
if(n<=m)return g[n];
ll ans=0;
for(ll i=1,pre;i<=n;i=pre+1)
{
pre=n/(n/i);
ans+=(pre-i+1)*(n/i);
}
return ans;
}
ll a[10005],mx;
int main()
{
int T;
scanf("%d",&T);
for(int i=1;i<=T;i++)scanf("%lld",&a[i]),mx=max(mx,a[i]);
mx=min(mx,(ll)m);
get_prime(mx);
for(int j=1;j<=T;j++)
{
ll ans=0,last=0,now,n=a[j];
for(ll i=1,pre;i<=n;i=pre+1)
{
pre=n/(n/i);
now=F(pre);
ans+=(now-last)*G(n/i);
last=now;
}
printf("%lldn",ans);
}
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
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
#include<cstdio>
#include<algorithm>
#include<vector>
#include<cmath>
using namespace std;
typedef long long ll;
#define maxn 100000005
int m=100000000;
int mu[maxn],num[maxn],f[maxn],p[maxn/10],res=0;
ll g[maxn];
bool vis[maxn];
void get_prime(int n)
{
mu[1]=1,f[1]=1,g[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i])p[res++]=i,mu[i]=-1,g[i]=2,num[i]=1;
f[i]=f[i-1]+(mu[i]?1:0);
for(int j=0;j<res&&i*p[j]<=n;j++)
{
vis[i*p[j]]=1;
if(i%p[j])
{
mu[i*p[j]]=-mu[i];
g[i*p[j]]=g[i]*2;
num[i*p[j]]=1;
}
else
{
mu[i*p[j]]=0;
num[i*p[j]]=num[i]+1;
g[i*p[j]]=g[i]/(num[i]+1)*(num[i]+2);
break;
}
}
}
for(int i=2;i<=n;i++)g[i]+=g[i-1];
}
ll F(ll n)
{
if(n<=m)return f[n];
int nn=(int)sqrt(n+0.5);
ll ans=0;
for(int d=1;d<=nn;d++)
if(mu[d])ans+=mu[d]*(n/d/d);
return ans;
}
ll G(ll n)
{
if(n<=m)return g[n];
ll ans=0;
for(ll i=1,pre;i<=n;i=pre+1)
{
pre=n/(n/i);
ans+=(pre-i+1)*(n/i);
}
return ans;
}
ll a[10005],mx;
int main()
{
int T;
scanf("%d",&T);
for(int i=1;i<=T;i++)scanf("%lld",&a[i]),mx=max(mx,a[i]);
mx=min(mx,(ll)m);
get_prime(mx);
for(int j=1;j<=T;j++)
{
ll ans=0,last=0,now,n=a[j];
for(ll i=1,pre;i<=n;i=pre+1)
{
pre=n/(n/i);
now=F(pre);
ans+=(now-last)*G(n/i);
last=now;
}
printf("%lldn",ans);
}
return 0;
}