「LOJ6482」LJJ爱数数

「LOJ6482」LJJ爱数数

解题思路 :

打表发现两个数 (a, b) 合法的充要条件是(我不管,我就是打表过的):

[a + b = ext{gcd}(a, b)^2 ]

(g = ext{gcd(a, b)}) ,那么相当于是要求:

[sum_{g=1}^{sqrt{2n}}sum_{i}[ ext{gcd}(g^2-ig, ig)=g] ]

化简一波:

[sum_{g=1}^{sqrt{2n}}sum_{i}[ ext{gcd}(g-i, i)=1] ]

根据辗转相除:

[sum_{g=1}^{sqrt{2n}}sum_{i}[ ext{gcd}(g, i)=1] ]

考虑 (i) 的上界和下界

[1 leq ig leq n \ 1 leq g^2 -ig leq n ]

解一下这两个不等式:

[ ext{max}_i =min(lfloorfrac{n}{g} floor,g - 1) \ ext{min}_i =max(g-lfloorfrac{n}{g} floor,1) ]

原来的式子相当于求:

[sum_{g=1}^{sqrt{2n}}sum_{i= ext{min}}^{max}[ ext{gcd}(g, i)=1] ]

(f(n)) 表示 ([1, n]) 之间与 (g) 互质的数的个数,反演一波可以得到:

[f(n)= sum_{d|g} lfloor frac{n}{d} floor mu(d) ]

再化简一波式子:

[sum_{g=1}^{sqrt{2n}}f(max) -f(min-1) ]

总复杂度 (O(sqrt{n}logn))

code

/*program by mangoyang*/
#pragma GCC optimize("Ofast", "inline")
#include<bits/stdc++.h>
#define inf ((int)(1e9))
#define Max(a, b) ((a) > (b) ? (a) : (b))
#define Min(a, b) ((a) < (b) ? (a) : (b))
typedef long long ll;
using namespace std;
template <class T>
inline void read(T &x){
    int f = 0, ch = 0; x = 0;
    for(; !isdigit(ch); ch = getchar()) if(ch == '-') f = 1;
    for(; isdigit(ch); ch = getchar()) x = x * 10 + ch - 48;
    if(f) x = -x;
}
const int N = 1500005;
int prime[N], d[30000005], tmp[N], mu[N], len[N], b[N], tot;
ll n, ans, m;
inline int solve(int x, int n){
	int ans = 0;
	for(register int i = len[x-1] + 1; i <= len[x]; i++) ans += mu[d[i]] * (n / d[i]);
	return ans;
}
int main(){
	read(n), mu[1] = 1;
	for(int i = 2; i < N; i++){
		if(!b[i]) prime[++tot] = i, mu[i] = -1;
		for(int j = 1; j <= tot && i * prime[j] < N; j++){
			b[i*prime[j]] = 1;
			if(i % prime[j] == 0){ mu[i*prime[j]] = 0; break; }
			mu[i*prime[j]] = -mu[i];
		}
	}
	m = (int) sqrt(2ll * n);
	for(int i = 1; i <= m; i++) if(mu[i])
		for(int j = i; j <= m; j += i) len[j]++;
	for(int i = 1; i <= m; i++) len[i] += len[i-1];
	for(int i = 1; i <= m; i++) if(mu[i])
		for(int j = i; j <= m; j += i) d[(++tmp[j])+len[j-1]] = i;
	for(int g = 1; g <= m; g++)
		ans += solve(g, Min(n / g, g - 1)) - solve(g, Max(1, g - n / g) - 1);
	cout << ans << endl;
	return 0;
}
原文地址:https://www.cnblogs.com/mangoyang/p/10173265.html