【题目链接】
【思路要点】
- 给出结论:
定义函数 χ ( x ) ( x ∈ N + ) \chi(x)\ (x\in\mathbb{N^{+}}) χ(x) (x∈N+) ,满足 χ ( x ) = { 1 x ≡ 1 ( m o d 4 ) − 1 x ≡ 3 ( m o d 4 ) 0 x ≡ 0 ( m o d 2 ) \chi(x)=\left\{\begin{array}{rcl}1&&{x\equiv 1\ (mod\ 4)}\\-1&&{x\equiv 3\ (mod\ 4)}\\0&&{x\equiv 0\ (mod\ 2)}\end{array} \right. χ(x)=⎩⎨⎧1−10x≡1 (mod 4)x≡3 (mod 4)x≡0 (mod 2)
令 x 2 = p 1 k 1 p 2 k 2 p 3 k 3 . . . p s k s x^2=p_1^{k_1}p_2^{k_2}p_3^{k_3}...p_s^{k_s} x2=p1k1p2k2p3k3...psks ,其中 p 1 , p 2 , p 3 , . . . , p s p_1,p_2,p_3,...,p_s p1,p2,p3,...,ps 均为质数, k 1 , k 2 , k 3 , . . . , k s k_1,k_2,k_3,...,k_s k1,k2,k3,...,ks 均为正整数。那么有 f ( x ) = 4 × ∏ i = 1 s ∑ j = 0 k i χ ( p i j ) = 4 × ∑ i ∣ x χ ( i ) f(x)=4\times\prod_{i=1}^{s}\sum_{j=0}^{k_i}\chi(p_i^j)=4\times\sum_{i\mid x}\chi(i) f(x)=4×i=1∏sj=0∑kiχ(pij)=4×i∣x∑χ(i)- 因此 f ( x ) f(x) f(x) 除去 4 4 4 后为积性函数,使用传统的积性函数求和即可。
- 关于结论的证明可以参考笔者 2019 2019 2019 年的国家集训队论文。
- 时间复杂度 O ( M i n 25 ( N ) ) O(Min25(N)) O(Min25(N)) 。
【代码】
#include<bits/stdc++.h> using namespace std; const int MAXN = 1e6 + 5; const int MAXLOG = 256; typedef long long ll; typedef long double ld; typedef unsigned long long ull; template <typename T> void chkmax(T &x, T y) {x = max(x, y); } template <typename T> void chkmin(T &x, T y) {x = min(x, y); } template <typename T> void read(T &x) { x = 0; int f = 1; char c = getchar(); for (; !isdigit(c); c = getchar()) if (c == '-') f = -f; for (; isdigit(c); c = getchar()) x = x * 10 + c - '0'; x *= f; } template <typename T> void write(T x) { if (x < 0) x = -x, putchar('-'); if (x > 9) write(x / 10); putchar(x % 10 + '0'); } template <typename T> void writeln(T x) { write(x); puts(""); } ll n, k, val[MAXN], cnt[MAXN][2]; int m, P, k3, home1[MAXN], home2[MAXN]; int limit, tot, tpow[MAXN], prime[MAXN], pre[MAXN][2]; void update(int &x, int y) { x += y; if (x >= P) x -= P; } int power(int x, ll y) { if (y == 0) return 1; int tmp = power(x, y / 2); if (y % 2 == 0) return 1ll * tmp * tmp % P; else return 1ll * tmp * tmp % P * x % P; } int solve(ll x, int y) { if (x <= 1 || prime[y] > x) return 0; int ans = 0, pos = 0; if (x <= limit) pos = home1[x]; else pos = home2[n / x]; ans = ((cnt[pos][1] - pre[y - 1][1]) + (cnt[pos][0] - pre[y - 1][0]) % P * tpow[3]) % P; if (prime[y] == 2) update(ans, 1); for (int i = y; i <= tot && 1ll * prime[i] * prime[i] <= x; i++) { ll now = prime[i], nxt = 1ll * prime[i] * prime[i]; for (int j = 1; nxt <= x; j++, now *= prime[i], nxt *= prime[i]) { if (prime[i] % 4 == 1) update(ans, (1ll * tpow[2 * j + 1] * solve(x / now, i + 1) + tpow[2 * j + 3]) % P); else update(ans, (solve(x / now, i + 1) + 1) % P); } } return ans; } void sieve(int n) { static int f[MAXN]; for (int i = 2; i <= n; i++) { if (f[i] == 0) prime[++tot] = f[i] = i; for (int j = 1; j <= tot && prime[j] <= f[i]; j++) { int tmp = prime[j] * i; if (tmp > n) break; f[tmp] = prime[j]; } } for (int i = 1; i <= tot; i++) { pre[i][0] = pre[i - 1][0]; pre[i][1] = pre[i - 1][1]; if (prime[i] % 4 == 1) pre[i][0]++; if (prime[i] % 4 == 3) pre[i][1]++; } } int main() { read(n), read(k), read(P); sieve(limit = sqrt(n) + 5); for (ll i = 1, nxt; i <= n; i = nxt + 1) { ll tmp = n / i; val[++m] = tmp; if (tmp <= limit) home1[tmp] = m; else home2[n / tmp] = m; nxt = n / val[m]; cnt[m][0] = (tmp - 1) / 4; cnt[m][1] = (tmp + 1) / 4; } for (int i = 2; i <= tot; i++) for (int j = 1; 1ll * prime[i] * prime[i] <= val[j]; j++) { ll tmp = val[j] / prime[i]; int pos = 0; if (tmp <= limit) pos = home1[tmp]; else pos = home2[n / tmp]; if (prime[i] % 4 == 1) { cnt[j][0] -= cnt[pos][0] - (pre[i][0] - 1); cnt[j][1] -= cnt[pos][1] - pre[i][1]; } else { cnt[j][0] -= cnt[pos][1] - (pre[i][1] - 1); cnt[j][1] -= cnt[pos][0] - pre[i][0]; } } for (int i = 1; i <= MAXLOG; i++) tpow[i] = power(i, k); writeln(1ll * tpow[4] * (solve(n, 1) + 1) % P); return 0; }