先描述一下问题
已知m,n
求
∑
(
i
,
j
)
=
k
1
<
=
i
<
=
m
,
1
<
=
j
<
=
n
1
\sum_{(i,j)=k}^{1<=i<=m,1<=j<=n}1
(i,j)=k∑1<=i<=m,1<=j<=n1
用文字描述即为已知整数n,m,
1<=i<=m,1<=j<=n,求有多少对(i,j)满足i,j的最大公约数为k
算法1
首先最简单的算法即为暴力枚举
枚举所有符合条件的数对,判断是否满足要求即可
时间复杂度
o
(
m
n
)
o(mn)
o(mn)
算法1的优化
显然原式等价于
∑
(
i
,
j
)
=
1
1
<
=
i
<
=
m
/
k
,
1
<
=
j
<
=
n
/
k
1
\sum_{(i,j)=1}^{1<=i<=m/k,1<=j<=n/k}1
(i,j)=1∑1<=i<=m/k,1<=j<=n/k1
这样时间复杂度为
o
(
m
n
/
k
2
)
o(mn/k^2)
o(mn/k2)
算法2
显然这个时间复杂度是难以接受的
那么我们可以怎么做呢?
先考虑另一个问题,
有多少对(i,j)有公因子k
这个显然是
[
m
k
]
[
n
k
]
[ {\frac m k} ][{\frac n k }]
[km][kn]
那么这跟正确答案相差多少?
显然这个数包含了最大公约数为2k,3k,4k…的情况,
那么我们减去这些数就好了
为保证无后效性,我们不得不采用倒推,即先推出n/k*k,再依次往下,直到枚举到k
由考虑最差情况,调和级数可知 时间复杂度为O(nlogn)
例题:
【NOI2010】能量采集
https://www.luogu.org/problemnew/show/P1447
在分析清楚思路之后
代码异常简洁
#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <algorithm>
using namespace std;
#define ll long long
#define N 100007
ll n, m, ans, f[N];
int main()
{
cin >> n >> m;
if (n > m)
swap(m, n);
for (int i = n; i >= 1; i--)
{
f[i] = (m / i) * (n / i);
for (int j = i * 2; j <= n; j += i)
{
f[i] -= f[j];
}
ans += f[i] * i;
}
cout << 2 * ans - n * m;
return 0;
}
算法3
莫比乌斯反演
由算术基本定理知,任意正整数都可以拆分为多个质数的乘积
p
i
p_i
pi表示质数
则任意正整数
n
=
p
1
a
1
p
2
a
2
.
.
.
p
q
a
q
n = p_1^{a_1}p_2^{a_2}...p_q^{a_q}
n=p1a1p2a2...pqaq
引入莫比乌斯函数
μ
(
n
)
=
{
1
,
n
=
1
(
−
1
)
r
,
n
=
p
1
p
2
p
3
.
.
p
r
,
p
i
为
互
不
相
同
的
素
数
0
,
其
他
情
况
μ(n)= \begin {cases} 1,n = 1 \\(-1)^r, n=p_1p_2p_3..p_r,p_i为互不相同的素数\\0,其他情况\end{cases}
μ(n)=⎩⎪⎨⎪⎧1,n=1(−1)r,n=p1p2p3..pr,pi为互不相同的素数0,其他情况
根据mobius反转定理
如
果
F
(
n
)
=
∑
d
∣
n
f
(
n
)
,
则
有
f
(
n
)
=
∑
d
∣
n
μ
(
d
)
F
(
n
d
)
如果F(n)=\sum_{d|n}f(n),则有f(n)=\sum_{d|n}μ(d)F(\frac n d)
如果F(n)=d∣n∑f(n),则有f(n)=d∣n∑μ(d)F(dn)
一种等价形式
如
果
F
(
n
)
=
∑
n
∣
d
f
(
d
)
,
则
有
f
(
n
)
=
∑
n
∣
d
μ
(
d
)
F
(
d
n
)
如果F(n)=\sum_{n|d}f(d),则有f(n)=\sum_{n|d}μ(d)F(\frac d n)
如果F(n)=n∣d∑f(d),则有f(n)=n∣d∑μ(d)F(nd)
使用后一种形式
设F(x)为k|(i,j)的对数
那么
F
(
d
)
=
[
n
k
]
[
m
k
]
=
∑
d
∣
i
f
(
i
)
F(d)=[\frac n k][\frac m k]=\sum_{d|i}f(i)
F(d)=[kn][km]=∑d∣if(i)
反演后得
f
(
d
)
=
∑
d
∣
i
μ
(
i
)
F
(
i
d
)
f(d)=\sum_{d|i}μ(i)F(\frac i d)
f(d)=∑d∣iμ(i)F(di)
于是时间复杂度就变成O(n)了
算法4
仔细观察发现(可以打表)
[
n
k
]
[\frac n k]
[kn]只有
n
\sqrt n
n种取值
再利用μ(i)的前缀和
可以优化时间复杂度至
O
(
n
)
O(\sqrt n)
O(n)
(预处理O(n))
例题
[HAOI2011]Problem b
https://www.luogu.org/problemnew/solution/P2522
题目描述
对于给出的n个询问,每次求有多少个数对(x,y),满足a≤x≤b,c≤y≤d,且gcd(x,y) = k,gcd(x,y)函数为x和y的最大公约数。
100%的数据满足:1≤n≤50000,1≤a≤b≤50000,1≤c≤d≤50000,1≤k≤50000
套用算法4,再利用容斥原理
#include <iostream>
#include <cstdio>
#include <algorithm>
#include <cstdio>
#include <cstdlib>
using namespace std;
#define N 50000
#define ll long long
int t,k,a,b,c,d,mu[N+6],sumu[N+7],prime[N+6],tot;
bool notprime[N+6];
void init()
{
sumu[1]=mu[1]=1;
notprime[1]=1;
for(int i = 2; i <= N; i ++ )
{
if(!notprime[i])
{
prime[++tot]=i;
mu[i]=-1;
}
sumu[i]=sumu[i-1]+mu[i];
for(int j = 1;i*prime[j]<=N&&j <= tot; j ++)
{
notprime[i*prime[j]]=1;
if(i%prime[j]!=0)
{
mu[i*prime[j]]=-mu[i];
}
else break;
}
}
}
ll f(int m,int n)
{
ll ans = 0;
int last ;
if(m>n)swap(m,n);
for(int i = 1 ; i <= m; i =last +1)
{
last = min(n/(n/i),m/(m/i));
ans += ((ll)sumu[last]-sumu[i-1])*(n/i)*(m/i);
}
return ans;
}
int main()
{
init();
scanf("%d",&t);
while(t--)
{
scanf("%d%d%d%d%d",&a,&b,&c,&d,&k);
a--;c--;
a/=k;b/=k;c/=k;d/=k;
printf("%lld\n",f(b,d)-f(a,d)-f(b,c)+f(a,c));
}
}