牛客挑战赛46 E 反演
反演
https://ac.nowcoder.com/acm/contest/9510/E
考虑使用min25。
令f(p^c)是当p为质数时的答案,根据积性函数的性质,可以进行min25筛。
但是选择的质数小于等于m时,f(p^c) != c + 1。
处理这些奇异值即可进行min25筛。
#include <cstdio>
#include <algorithm>
#include <cstring>
using namespace std;
typedef long long ll;
const int N = 1e6 + 100;
const int MOD = 998244353;
int Sqr, K, m, tot, id1[N], id2[N];
bool is_prime[N];
ll n, pri[N], w[N], g[N], spf[N];
ll tt[110], cc[110], rev[110], R;
void Sieve(int n) {
memset(is_prime, true, sizeof(is_prime));
is_prime[1] = false;
for (int i = 2; i <= n; i++) {
if (is_prime[i]) {
pri[++tot] = i;
spf[tot] = spf[tot - 1] + 1;
}
for (int j = 1; i * pri[j] <= n; j++) {
is_prime[i * pri[j]] = false;
if (i % pri[j] == 0) break;
}
}
}
ll qpow(ll x, ll n) {
ll res = 1;
while (n > 0) {
if (n & 1) res = res * x % MOD;
x = x * x % MOD;
n /= 2;
}
return res;
}
ll F(ll p, ll e) {
if (p <= K) return rev[p] * ((e + 1) * (tt[p] + 1) % MOD + cc[p]) % MOD; //奇异值
return e + 1;
}
ll S(ll x, int y) {
if (x <= 1 || pri[y] > x) return 0;
int k = (x <= Sqr) ? id1[x] : id2[n / x];
ll res = g[k] - spf[y - 1];
for (int i = y; i <= tot && 1LL * pri[i] * pri[i] <= x; i++) {
ll p1 = pri[i], p2 = pri[i] * pri[i];
for (int e = 1; p2 <= x; e++, p1 = p2, p2 *= pri[i]) {
res = (res + S(x / p1, i + 1) * F(pri[i], e) + R * F(pri[i], e + 1)) % MOD;
}
}
return res;
}
int solve() {
m = 0; tot = 0;
Sqr = sqrt(n); Sieve(Sqr);
for (ll i = 1, j; i <= n; i = j + 1) {
j = n / (n / i); w[++m] = n / i;
if (w[m] <= Sqr) id1[w[m]] = m;
else id2[n / w[m]] = m;
g[m] = (w[m] - 1) % MOD;
}
for (int j = 1; j <= tot; j++) {
for (int i = 1; i <= m && 1LL * pri[j] * pri[j] <= w[i]; i++) {
int k = (w[i] / pri[j] <= Sqr) ? id1[w[i] / pri[j]] : id2[n / (w[i] / pri[j])];
g[i] = (g[i] - (g[k] - spf[j - 1])) % MOD;
}
}
//上半部分是min25模板,下半部分处理奇异值。
R = 1;
for (int j = 2; j <= K; j++) {
if (!is_prime[j]) continue;
for (int i = j * 2; i <= K; i += j) is_prime[i] = false;
for (int i = 2; i <= K; i++) {
int x = i;
while (x % j == 0) {
x /= j;
tt[j]++;
cc[j] += tt[j];
}
}
R = R * (tt[j] + 1 + cc[j]) % MOD;
rev[j] = qpow(tt[j] + 1 + cc[j], MOD - 2);
}
ll sum = 0, cnt = 0;
for (int i = m, j = 1; i >= 1; i--) {
while (j + 1 <= K && j + 1 <= w[i]) {
j++;
if (!is_prime[j]) continue;
sum = (sum + R * rev[j] % MOD
* (2 * (tt[j] + 1) + cc[j])) % MOD;
cnt++;
}
g[i] = (2 * (g[i] - cnt) * R + sum) % MOD;
}
for (int j = 1; j <= tot; j++) {
if (pri[j] <= K) spf[j] = spf[j - 1] + R * rev[pri[j]] % MOD
* (2 * (tt[pri[j]] + 1) + cc[pri[j]]) % MOD;
else spf[j] = (spf[j - 1] + R * 2) % MOD;
}
return S(n, 1) + R;
}
int main() {
//freopen("0.txt", "r", stdin);
scanf("%lld%d", &n, &K);
ll ans = solve();
ans = (ans % MOD + MOD) % MOD;
printf("%lld\n", ans);
return 0;
}
查看10道真题和解析
美的集团公司福利 717人发布