参考博客:https://www.90yang.com/2019ccpc-xiangtan-f-neko-and-function/
思维比较巧妙,首先从贡献的角度计算 f ( n , k ) f(n,k) f(n,k),枚举 n 的每一个质因子,设其幂次为 e e e,问题等价于将这 e e e 个 p p p 放到 k k k 个盒子里。
由于盒子不允许为空,这个贡献很难计算,考虑盒子可以为空,用容斥得到最后的答案:令 g ( n , k ) g(n,k) g(n,k) 为 f ( n , k ) f(n,k) f(n,k) 去掉 a i > 1 a_i > 1 ai>1 的限制。
没有这个限制,每一个素数的贡献独立,即 g ( n , k ) g(n,k) g(n,k) 是一个积性函数
考虑容斥得到
f
(
n
,
k
)
f(n,k)
f(n,k):枚举有几个盒子为空,容斥一下即可:
f
(
n
,
k
)
=
∑
x
=
0
k
(
−
1
)
k
−
x
∗
C
(
k
,
x
)
∗
g
(
n
,
x
)
\displaystyle f(n,k) = \sum_{x = 0}^k(-1)^{k - x}*C(k,x)*g(n,x)
f(n,k)=x=0∑k(−1)k−x∗C(k,x)∗g(n,x)
∑
i
=
1
n
f
(
i
,
k
)
=
∑
x
=
0
k
(
−
1
)
k
−
x
∗
C
(
k
,
x
)
∗
∑
i
=
1
n
g
(
i
,
x
)
\displaystyle \sum_{i = 1}^nf(i,k) = \sum_{x = 0}^k(-1)^{k - x}*C(k,x)*\sum_{i = 1}^ng(i,x)
i=1∑nf(i,k)=x=0∑k(−1)k−x∗C(k,x)∗i=1∑ng(i,x)
由于
g
(
n
,
k
)
g(n,k)
g(n,k) 是积性函数,考虑用min_25筛得到
∑
i
=
1
n
g
(
i
,
x
)
\displaystyle\sum_{i = 1}^ng(i,x)
i=1∑ng(i,x)
由于这题卡常,我最终也没有卡过去,待update
#include<bits/stdc++.h>
using namespace std;
const int maxn = 1e6 + 10;
typedef long long ll;
const int mod = 1e9 + 7;
int num,tot,k;
bool ispri[maxn];
int sum[maxn],pri[maxn],w[maxn],g[maxn],id1[maxn],id2[maxn];
int n,sqr,fact[maxn],ifact[maxn],ans[maxn];
inline void sieve(int x) {
num = 0; ispri[0] = ispri[1] = true;
for (int i = 2; i <= x; i++) {
if (!ispri[i]) {
pri[++num] = i;
sum[num] = (sum[num - 1] + 1) % mod;
}
for (int j = 1; j <= num && i * pri[j] <= x; j++) {
ispri[i * pri[j]] = true;
if (i % pri[j] == 0) break;
}
}
}
int fpow(int a,int b) {
int r = 1;
while (b) {
if (b & 1) r = 1ll * r * a % mod;
b >>= 1;
a = 1ll * a * a % mod;
}
return r;
}
inline int C(int x,int y) {
if (y > x || y < 0) return 0;
return 1ll * fact[x] * ifact[y] % mod * ifact[x - y] % mod;
}
inline int S(ll x,int y,int z) {
if (pri[y] >= x) return 0;
int k = x <= sqr ? id1[x] : id2[n / x];
int ans = 1ll * z * (g[k] - sum[y] + mod) % mod;
for (int i = y + 1; i <= num && pri[i] <= x / pri[i]; i++)
for (ll e = 1, pe = pri[i]; pe <= x; e++, pe = pe * pri[i])
ans = (ans + 1ll * (S(x / pe,i,z) + (e != 1)) % mod * C(e + z - 1,z - 1) % mod) % mod;
return ans;
}
int main() {
int N = maxn - 10;
sieve(N);
fact[0] = 1;
for (int i = 1; i <= N; i++)
fact[i] = 1ll * fact[i - 1] * i % mod;
ifact[N] = fpow(fact[N],mod - 2);
for (int i = N - 1; i >= 0; i--)
ifact[i] = 1ll * ifact[i + 1] * (i + 1) % mod;
while(scanf("%d%d",&n,&k) != EOF) {
if (k == 1) {
printf("%d\n",(n - 1) % mod);
continue;
}
else if (n < (1<<k)) {
puts("0");
continue;
}
else if (n < (1 << k) + (1 << (k - 1))) {
puts("1");
continue;
}
tot = 0;
sqr = sqrt(n);
for (int i = 1,j; i <= n; i = j + 1) {
j = n / (n / i);
w[++tot] = n / i;
int p = n / i % mod;
g[tot] = (p - 1) % mod;
if (n / i <= sqr) id1[n / i] = tot;
else id2[j] = tot;
}
for (int i = 1; i <= num; i++) {
for (int j = 1; j <= tot && pri[i] <= w[j] / pri[i]; j++) {
int x = w[j] / pri[i];
x = x <= sqr ? id1[x] : id2[n / x];
g[j] -= (g[x] - sum[i - 1] + mod) % mod;
if (g[j] < 0) g[j] += mod;
}
}
int res = 0, t = 1;
for (int i = k; i >= 0; i--) {
if (k == 0) res = (res + C(k,i) * t % mod) % mod;
else if (k == 1) res = (res + 1ll * n * C(k,i) % mod * t % mod) % mod;
else res = (res + 1ll * S(n,0,i) * C(k,i) % mod * t % mod) % mod;
if (res < 0) res += mod;
t *= -1;
}
printf("%d\n",res);
}
return 0;
}