BZOJ4652 [Noi2016]循环之美 数论 + 莫比乌斯反演 + 杜教筛

Posted Mychael

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了BZOJ4652 [Noi2016]循环之美 数论 + 莫比乌斯反演 + 杜教筛相关的知识,希望对你有一定的参考价值。

题目链接

BZOJ

题解

orz
此题太优美了

我们令\\(\\frac{x}{y}\\)为最简分数,则\\(x \\perp y\\)即,\\(gcd(x,y) = 1\\)

先不管\\(k\\)进制,我们知道\\(10\\)进制下如果\\(\\frac{x}{y}\\)是纯循环的,只要\\(2 \\perp y\\)\\(5 \\perp y\\)
可以猜想在\\(k\\)进制下同样成立
证明:
\\(\\frac{x}{y}\\)为纯循环小数,设其循环节长度为\\(l\\),那么一定满足

\\[\\{ \\frac{xk^{l}}{y} \\} = \\{ \\frac{x}{y}\\} \\]

其中\\(\\{x\\}\\)指小数部分
可以写成

\\[xk^{l}- \\lfloor \\frac{xk^{l}}{y} \\rfloor * y= x - \\lfloor \\frac{x}{y} \\rfloor * y \\]

可以看出就是同余的形式

\\[xk^{l} \\equiv x \\pmod y \\]

由于\\(x \\perp y\\)

\\[k^{l} \\equiv 1 \\pmod y \\]

要使存在\\(l\\),使得式子成立,那么一定\\(k \\perp y\\)

所以我们有

\\[\\begin{aligned} ans &= \\sum\\limits_{i = 1}^{n} \\sum\\limits_{j = 1}^{m} [i \\perp j] [j \\perp k] \\\\ &= \\sum\\limits_{j = 1}^{m} [j \\perp k] \\sum\\limits_{i = 1}^{n} [(i,j) == 1] \\\\ &= \\sum\\limits_{j = 1}^{m} [j \\perp k] \\sum\\limits_{i = 1}^{n} \\epsilon((i,j)) \\\\ &= \\sum\\limits_{j = 1}^{m} [j \\perp k] \\sum\\limits_{i = 1}^{n} \\sum\\limits_{d|(i,j)} \\mu(d) \\\\ &= \\sum\\limits_{d = 1}^{n} \\mu(d) \\sum\\limits_{d|i}^{n} \\sum\\limits_{d|j}^{m} [j \\perp k] \\\\ &= \\sum\\limits_{d = 1}^{n} [d \\perp k] \\mu(d) \\lfloor \\frac{n}{d} \\rfloor \\sum\\limits_{j}^{\\lfloor \\frac{m}{d} \\rfloor} [j \\perp k] \\\\ \\end{aligned} \\]

我们令

\\[f(n) = \\sum\\limits_{i}^{n} [i \\perp k] \\]

根据\\(gcd(i + k,k) = gcd(i,k)\\),我们有

\\[f(n) = \\lfloor \\frac{n}{k} \\rfloor f(k) + f(n \\mod k) \\]

所以我们只需要\\(O(klogk)\\)暴力计算\\(f(1....k)\\)就可以\\(O(1)\\)计算后面的式子了

现在考虑前面的

\\[\\sum\\limits_{d = 1}^{n} [d \\perp k] \\mu(d) \\]

为了使能整除分块,我们必须算出其前缀和

\\[g(n,k) = \\sum\\limits_{i = 1}^{n} [i \\perp k] \\mu(d) \\]

用类似上面同样的方法:

\\[\\begin{aligned} g(n,k) &= \\sum\\limits_{i = 1}^{n} [i \\perp k] \\mu(i) \\\\ &= \\sum\\limits_{i = 1}^{n} \\mu(i) \\sum\\limits_{d|i,d|k} \\mu(d) \\\\ &= \\sum\\limits_{d | k} \\mu(d) \\sum\\limits_{d|i} \\mu(i) \\\\ &= \\sum\\limits_{d | k} \\mu(d) \\sum\\limits_{i = 1}^{\\lfloor \\frac{n}{d} \\rfloor} \\mu(id) \\\\ &= \\sum\\limits_{d | k} \\mu(d) \\sum\\limits_{i = 1}^{\\lfloor \\frac{n}{d} \\rfloor} [i \\perp d]\\mu(i) * \\mu(d) \\\\ &= \\sum\\limits_{d | k} \\mu(d)^2 \\sum\\limits_{i = 1}^{\\lfloor \\frac{n}{d} \\rfloor} [i \\perp d]\\mu(i) \\\\ &= \\sum\\limits_{d | k} \\mu(d)^2 g(\\lfloor \\frac{n}{d} \\rfloor,d) \\\\ \\end{aligned} \\]

就可以递归求解
\\(n = 0\\)时,\\(g(0,k) = 0\\)
\\(k = 1\\)时,\\(g(n,1) = \\sum\\limits_{i = 1}^{n} \\mu(i)\\),上杜教筛即可
因为\\(\\lfloor \\frac{n}{d} \\rfloor\\)只有\\(\\sqrt{n}\\)种取值,\\(d\\)只能取\\(k\\)的因子,记其数量为\\(p\\)
那么求\\(g(n,k)\\)总的复杂度为\\(O(p\\sqrt{n} + n^{\\frac{2}{3}})\\)

于是乎我们就解决这道题了
总复杂度\\(O(p\\sqrt{n} + n^{\\frac{2}{3}} + \\sqrt{n} + klogk)\\)

#include<iostream>
#include<cstdio>
#include<cmath>
#include<map>
#include<cstring>
#include<algorithm>
#define LL long long int
#define Redge(u) for (int k = h[u],to; k; k = ed[k].nxt)
#define REP(i,n) for (int i = 1; i <= (n); i++)
#define BUG(s,n) for (int i = 1; i <= (n); i++) cout<<s[i]<<\' \'; puts("");
#define mp(a,b) make_pair<LL,LL>(a,b)
using namespace std;
const int maxn = 1000005,maxm = 100005,INF = 1000000000;
map<LL,LL> _mu;
map<LL,LL>::iterator it;
map<pair<LL,LL>,LL> G;
map<pair<LL,LL>,LL>::iterator IT;
int N;
int p[maxn],pi,isn[maxn];
LL Smu[maxn],f[maxn],mu[maxn];
LL n,m,K;
int gcd(int a,int b){return b ? gcd(b,a % b) : a;}
void init(){
	N = 1000000;
	mu[1] = 1;
	for (int i = 2; i <= N; i++){
		if (!isn[i]) p[++pi] = i,mu[i] = -1;
		for (int j = 1; j <= pi && i * p[j] <= N; j++){
			isn[i * p[j]] = true;
			if (i % p[j] == 0){
				mu[i * p[j]] = 0;
				break;
			}
			mu[i * p[j]] = -mu[i];
		}
	}
	for (int i = 1; i <= N; i++) Smu[i] = Smu[i - 1] + mu[i];
	for (int i = 1; i <= K; i++) f[i] = f[i - 1] + (gcd(i,K) == 1);
}
LL S(LL n){
	if (n <= N) return Smu[n];
	if ((it = _mu.find(n)) != _mu.end()) return it->second;
	LL ans = 1;
	for (LL i = 2,nxt; i <= n; i = nxt + 1){
		nxt = n / (n / i);
		ans -= (nxt - i + 1) * S(n / i);
	}
	return _mu[n] = ans;
}
LL F(LL n){
	return (n / K) * f[K] + f[n % K];
}
LL g(LL n,LL k){
	if ((IT = G.find(mp(n,k))) != G.end())
		return IT->second;
	if (n == 0) return 0;
	if (k == 1) return S(n);
	LL ans = 0;
	for (LL i = 1; i * i <= k; i++){
		if (k % i == 0){
			if (mu[i]) ans += g(n / i,i);
			if (i * i != k && mu[k / i])
				ans += g(n / (k / i),k / i);
		}
	}
	return G[mp(n,k)] = ans;
}
int main(){
	cin >> n >> m >> K;
	init();
	LL ans = 0,now,last = 0;
	for (LL i = 1,nxt; i <= min(n,m); i = nxt + 1){
		nxt = min(n / (n / i),m / (m / i));
		now = g(nxt,K);
		ans += 1ll * (now - last) * (n / i) * F(m / i);
		last = now;
	}
	cout << ans << endl;
	return 0;
}

以上是关于BZOJ4652 [Noi2016]循环之美 数论 + 莫比乌斯反演 + 杜教筛的主要内容,如果未能解决你的问题,请参考以下文章

[BZOJ4652][NOI2016]循环之美

Bzoj4652: [Noi2016]循环之美

BZOJ4652 [Noi2016]循环之美 数论 + 莫比乌斯反演 + 杜教筛

NOI 2016 循环之美

UOJ #221 NOI2016 循环之美

[NOI2016]循环之美