【洛谷】P5348 密码解锁
很显然我们可以推导出这个式子
设(a(m))为(m)位置的值
[mu(m) = sum_{m | d} a(d) \
a(m) = sum_{m|d}mu(frac{d}{m})mu(d) \
a(m) = sum_{i = 1}^{lfloor frac{n}{m}
floor} mu(i)mu(im) \
a(m) = mu(m) sum_{i = 1}^{lfloor frac{n}{m}
floor} mu(i)^{2}[gcd(m,i) == 1]
]
而(mu(i)^{2})的本质是无平方因子数，这个可以容斥
容斥的方法是（若没有其他限制）
[ans = sum_{i = 1}^{sqrt{N}} mu(i)lfloor frac{N}{i^{2}}
floor
]
那么这里的就是
[ans = sum_{i = 1}^{sqrt{N / M}} [gcd(i,m) == 1]mu(i)sum_{j = 1}^{lfloor frac{N}{i ^ 2}
floor} [gcd(j,m) == 1]
]
前面的互质可以直接枚举
后面的互质可以通过莫比乌斯反演外加预处理M中莫比乌斯值不为0的数算出来
#include <bits/stdc++.h>
#define fi first
#define se second
#define pii pair<int,int>
#define mp make_pair
#define pb push_back
#define space putchar(' ')
#define enter putchar('
')
#define eps 1e-10
#define ba 47
#define MAXN 100005
//#define ivorysi
using namespace std;
typedef long long int64;
typedef unsigned int u32;
typedef double db;
template<class T>
void read(T &res) {
res = 0;T f = 1;char c = getchar();
while(c < '0' || c > '9') {
if(c == '-') f = -1;
c = getchar();
}
while(c >= '0' && c <= '9') {
res = res * 10 +c - '0';
c = getchar();
}
res *= f;
}
template<class T>
void out(T x) {
if(x < 0) {x = -x;putchar('-');}
if(x >= 10) {
out(x / 10);
}
putchar('0' + x % 10);
}
int64 N,MAXV;
int M,mu[1000005];
int prime[1000005],tot;
bool nonprime[1000005];
vector<pii > division;
int gcd(int a,int b) {
return b == 0 ? a : gcd(b,a % b);
}
int Mu(int x) {
if(x <= 1000000) return mu[x];
int res = 1;
for(int i = 2 ; i <= x / i ; ++i) {
if(x % i == 0) {
int c = 0;
while(x % i == 0) {x /= i;++c;}
if(c >= 2) return 0;
res = -res;
}
}
if(x != 1) res = -res;
return res;
}
int64 calc(int64 n) {
int64 res = 0;
for(auto t : division) {
if(n < t.fi) break;
res += 1LL * t.se * (n / t.fi);
}
return res;
}
void Init() {
mu[1] = 1;
for(int i = 2 ; i <= 1000000 ; ++i) {
if(!nonprime[i]) {
prime[++tot] = i;
mu[i] = -1;
}
for(int j = 1 ; j <= tot ; ++j) {
if(prime[j] > 1000000 / i) break;
nonprime[i * prime[j]] = 1;
if(i % prime[j] == 0) break;
else mu[i * prime[j]] = -mu[i];
}
}
}
void Solve() {
read(N);read(M);
if(Mu(M) == 0) {puts("0");return;}
division.clear();
for(int i = 1 ; i <= M / i ; ++i) {
if(M % i == 0) {
int j = M / i;
int x = Mu(i),y = Mu(j);
if(x) division.pb(mp(i,x));
if(i != j && y) division.pb(mp(j,y));
}
}
sort(division.begin(),division.end());
int64 T = N / M,res = 0;
for(int i = 1 ; i <= T / i ; ++i) {
if(gcd(i,M) == 1) {
res += Mu(i) * calc(T / (i * i));
}
}
res = res * Mu(M);
out(res);enter;
}
int main(){
#ifdef ivorysi
freopen("f1.in","r",stdin);
#endif
Init();
int T;
read(T);
for(int i = 1 ; i <= T ; ++i) Solve();
}