CF439E Devu and Birthday Celebration题解
做法：DP + 暴力容斥
-
分析
观察题面给出的两个要求，要求二感觉不太好计算，那么就先计算满足要求一的数列个数。
尝试 DP，设
f_{i,j} 表示和为j 的i 位数列个数，那么f_{i,j}=\sum_{k=i-1}^{j-1}f_{i-1,k} 。因为
f_{i,j-1}=\sum_{k=i-1}^{j-2}f_{i-1,k} ，所以f_{i,j}=f_{i,j-1}+f_{i-1,j-1} 。我们发现这个式子非常熟悉，因为杨辉三角的递推式子也是这个，于是我们就可以直接套用杨辉三角的通项公式，所以f_{i,j}=\binom{j-1}{i-1} 。结合题面给出的变量，符合第一个条件的数列个数为
\binom{n-1}{f-1} 。接下来考虑条件二，然后我们可以发现，一个
\gcd 不为1 的数列之和显然为\gcd 的倍数且数列里每一位上的数a_i 均可以表示为a_i \times \gcd 。那么如果我们设g_{i,j,k} 表示一个长度为k 的数列的和为j 且\gcd 为i 的倍数时的数列数，那么g_{i,j,k}=[i\;|\;j] \times [\frac{j}{i} \geq k] \times \binom{\frac{j}{i}-1}{k-1} 。那么我们是否可以枚举能够整除数列和的质数，然后减去它们的
g 呢？答案显然是否定的，因为当\gcd 为合数时，会出现被多个g 统计然后被减去多次的情况。于是考虑容斥。然后我们发现
2\times3\times5\times7\times11\times13\times17=510510>10^5 ，所以对于一个数而言最多存在6 个互不相同的质因子，于是预处理出所有因子只有1 \sim 6 个互不相同的质因子的数，然后进行容斥（减一个质因子的g ，加两个质因子的g ，再减去三个质因子的g\dots ）。如何预处理？暴力六层循环即可，接下来分析暴力六层循环的时间复杂度。
我们可以在当前质数之积大于
10^5 的时候跳出当前循环，于是我们的总操作次数就等于积合法处理数字数和判断数+ 积非法时判断数。积非法的时候判断数小于等于处理数字数，因为每一个积要么是从一个处理数字乘一个新质数转移而来的，要么是从一个处理数字增大了一个质数转移而来的。因为一个积非法后会直接跳出循环所以不可能从非法数字转移来，同时由于不一定每个处理数字增大之后一定非法，所以积非法的时候判断数小于等于处理数字数。
积合法处理数字数和判断数一定等于处理数字乘以一个常数，所以我们只要分析出处理数字的规模即可。因为算数基本定理，每个数的质因数分解是唯一的，所以不同质因数分解对应的数是不同的。可以发现我们枚举的几个质数就是积的质因数分解，而我们枚举时可以使质数序列严格递增保证不会重复出现两次枚举的质数均相同（因为就算前面和后面均相同，当前位上每个质数也只能出现一次，而且因为严格递增所以当前位上的质数也不可能在前面或者后面再次出现导致交换后实际上是重复的），所以这些质数之积，即处理数字是互不相同且均
\leq 10^5 的，所以处理数字的数量最多为10^5 个。积非法的时候判断数小于等于处理数字数，积合法处理数字数和判断数为处理数字数乘以一个常数，所以总操作次数的规模和处理数字的规模相同，均为
10^5 级别，即\mathcal{O(w)} 的时间复杂度（这里w 为值域，大小为10^5 ，后文中的w 为同一定义）。所以最后的时间复杂度为
\mathcal{O(w+q\sqrt{n}\log mod)} ，可以通过本题。 -
代码
#include <iostream>
#define int long long
using namespace std;
constexpr int MAXN(1000007);
constexpr int mod(1000000007);
int fct[MAXN], vis[MAXN], pr[MAXN], num[MAXN];
int q, cnt;
inline void read(int &temp) { cin >> temp; }
inline void DealFact() { fct[0] = 1; for (int i(1); i <= 100000; ++i) fct[i] = fct[i - 1] * i % mod; }
inline void DealPrime() {
for (int i(2); i <= 100000; ++i) {
if (!vis[i]) pr[++cnt] = i;
for (int j(1); j <= cnt; ++j) {
if (i * pr[j] > 100000) break;
vis[i * pr[j]] = 1;
if (i % pr[j] == 0) break;
}
}
}
inline void DealRc() {
for (int a(1); a <= cnt; ++a) {
num[pr[a]] = 1;
for (int b(a + 1); b <= cnt; ++b) {
if (pr[a] * pr[b] > 100000) break;
num[pr[a] * pr[b]] = 2;
for (int c(b + 1); c <= cnt; ++c) {
if (pr[a] * pr[b] * pr[c] > 100000) break;
num[pr[a] * pr[b] * pr[c]] = 3;
for (int d(c + 1); d <= cnt; ++d) {
if (pr[a] * pr[b] * pr[c] * pr[d] > 100000) break;
num[pr[a] * pr[b] * pr[c] * pr[d]] = 4;
for (int e(d + 1); e <= cnt; ++e) {
if (pr[a] * pr[b] * pr[c] * pr[d] * pr[e] > 100000) break;
num[pr[a] * pr[b] * pr[c] * pr[d] * pr[e]] = 5;
for (int f(e + 1); f <= cnt; ++f) {
if (pr[a] * pr[b] * pr[c] * pr[d] * pr[e] * pr[f] > 100000) break;
num[pr[a] * pr[b] * pr[c] * pr[d] * pr[e] * pr[f]] = 6;
}}}}}}
}
inline int ksm(int base, int k) {
int res(1);
while (k) {
if (k & 1) res = res * base % mod;
base = base * base % mod, k >>= 1;
}
return res;
}
inline int C(int n, int m) { return fct[n] * ksm(fct[n - m] * fct[m] % mod, mod - 2) % mod; }
inline int calc(int n, int k) {
int res = C(n - 1, k - 1);
for (int i(1); i * i <= n; ++i) {
if (n % i == 0 && num[i] && n / i >= k) {
if (num[i] % 2 == 1) res = ((res - C(n / i - 1, k - 1)) % mod + mod) % mod;
if (num[i] % 2 == 0) res = (res + C(n / i - 1, k - 1)) % mod;
}
if (n % i == 0 && num[n / i] && i * i != n && i >= k) {
if (num[n / i] % 2 == 1) res = ((res - C(i - 1, k - 1)) % mod + mod) % mod;
if (num[n / i] % 2 == 0) res = (res + C(i - 1, k - 1)) % mod;
}
}
return res;
}
signed main() {
ios::sync_with_stdio(false), cin.tie(nullptr), cout.tie(nullptr);
read(q);
DealFact(), DealPrime(), DealRc();
for (int i(1), x, y; i <= q; ++i) {
read(x), read(y);
cout << calc(x, y) << '\n';
}
return 0;
}