题解:QOJ4264 / UOJ99 / P11739 [集训队互测 2015] 普罗达科特
题解:QOJ4264 [集训队互测 2015] 普罗达科特
题目描述
令 \(N = \prod_{i = 1}^{n} p_i^{A_i}\),\(M = \prod_{i = 1}^{n} p_i^{B_i}\),其中 \(p_i\) 是两两不同的素数。
求将 \(N\) 表示成 \(m\) 个正整数(原题写的是 \(k\) 个)的乘积的方案数,也就是 \(N = \prod_{i = 1}^k x_i\) 的解的个数,答案对 \(10^9 + 21\) 取模。
对于子问题 \(1\),要求对于所有整数 \(i\) 满足 \(1 \leq i < m\),都有 \(x_i < x_{i + 1}\)。
对于子问题 \(2\),要求对于所有整数 \(i\) 满足 \(1 \leq i < m\),都有 \(x_i \leq x_{i + 1}\)。
对于两个子问题都要求对于所有整数 \(i\) 满足 \(1 \leq i \leq m\) 都有 \(x_i \nmid M\),即 \(x_i\) 不是 \(M\) 的约数。
\(n\leq 50, A_i, B_i\leq 10^{18}, m\leq 25\)。
题解
分为两步:以枚举 \(m\) 的所有拆分 \(\lambda\vdash m\) 也就是 \(m=\lambda_1+\lambda_2+\cdots+\lambda_{|\lambda|}\) 为桥梁,分别计算这个拆分的贡献系数和方案数。拆分表达的意思是我们将所有 \(x_i\) 划分成一些等价类,等价类中的 \(x_i\) 是相等的。等价类大小排序之后就和这个拆分一样。不同等价类的数也可能相等。元素是没有标号的,只有等价类的划分带来的相等关系。然后计算这样的限制下有多少方案。计算完了以后再考虑乘上一定的系数作为两个问题的答案。
贡献系数 - 子问题 1
对于子问题 1,首先改成求使所有 \(x_i\) 互不相同的方案数,最后除去 \(m!\)。有一个容斥就是考虑 \(m\) 个数之间的相等关系,如果两个 \(x_i=x_j\) 我们就在 \((i,j)\) 之间连一条边,我们枚举所有完全图 \(\mathbb K_m\) 的子图,根据连通性区分等价类(一个等价类里的数连通),容斥系数就是 \((-1)^{|E|}\)(\(|E|\) 是边数)。由于最后要枚举的是拆分,考虑算一下一个拆分具体对应的系数。假设这个拆分中值为 \(i\) 的有 \(a_i\) 个。那么自然会先有一个多重集排列数,使元素带上标号:
多重集排列数的同时要注意 \(a_i\) 个大小相同的连通块是互不区分的。除此之外还要乘上
这里 \(f_i\) 表示大小为 \(i\) 的连通图的 \((-1)^{|E|}\) 之和是多少。这个是一个简单的东西:
为了计算 \(f_n\),容斥,首先计算大小为 \(n\) 的所有图的 \((-1)^{|E|}\) 之和,这是 \([n=1]\),然后枚举 \(1\) 号点所在连通块的大小,也就是:\(f_n=[n=1]-\sum_{i=1}^n\binom{n-1}{i-1}f_i[n-i=1]\)。
后面这个系数就很搞笑,可以删掉改成:\(f_n=[n=1]-(n-1)f_{n-1}\),又有 \(f_1=1\),所以 \(f_n=(-1)^{n-1}(n-1)!\)。
所以系数就是
这样就够了。
贡献系数 - 子任务 2
这里我们使用 Burnside 引理。对于一个方案,我们只取单调不降的,其它和它重排相同的方案我们不要。看到这个重排相同我们考虑使用 Burnside 引理:
在这里,我们枚举这个置换 \(g\in G\) 的置换环拆分得到的大小序列,也就是一个 \(m\) 的拆分 \(\lambda\vdash m\),它对应 \(\frac{m!\prod_{i=1}^m(a_i-1)!}{\prod_{i=1}^m(i!)^{a_i}a_i!}\) 个置换(和上面一样的多重集排列数,再考虑每个置换环内部的填法,也就是圆排列数),在这个 \(\lambda\vdash m\) 的拆分下,我们计算每个置换环里的元素都相同的方案数,乘上 \(\lambda\) 对应的置换个数,这个数就是这个拆分的贡献系数。
注意,\(\frac{1}{m!}\) 是 Burnside 引理带来的。
方案数
以枚举 \(m\) 的所有拆分 \(\lambda\vdash m\) 也就是 \(m=\lambda_1+\lambda_2+\cdots+\lambda_{|\lambda|}\) 为桥梁,分别计算这个拆分的贡献系数和方案数。拆分表达的意思是我们将所有 \(x_i\) 划分成一些等价类,等价类中的 \(x_i\) 是相等的。等价类大小排序之后就和这个拆分一样。不同等价类的数也可能相等。元素是没有标号的,只有等价类的划分带来的相等关系。然后计算这样的限制下有多少方案。
考虑对一个拆分 \(\lambda\vdash m\) 计算方案数。假设这个拆分中值为 \(i\) 的有 \(a_i\) 个。首先解决 \(M\) 的限制,枚举一个数组 \(c_i\)(暴力枚举)表示大小为 \(i\) 的等价类里面有多少个是 \(M\) 的约数,其容斥系数为:
这时所有的质因子就独立了。假设这个质因子是 \(j\),用生成函数进行计数:
改写成
前面是关于 \(x^{B_j+1}\) 的多项式,把它算出来之后,我们枚举它的所有 \(t\) 次项系数,计算
就是要计算
的某项系数,结论是,记 \(L=\operatorname{LCM}(\lambda_1, \lambda_2, \cdots, \lambda_{|\lambda|})\),则它的第 \(kL+r\) 项系数是关于 \(k\) 的 \(|\lambda|-1\) 次多项式。
这个可以感性理解一下,假如 \(\lambda_i=1\),这时就是求长度为 \(|\lambda|\) 的和为 \(kL+r\) 的非负整数数组的个数,这是一个插板法组合数 \(\binom{kL+r+|\lambda|-1}{|\lambda|-1}\),如果将 \(k\) 视作变量,就是关于 \(k\) 的 \(|\lambda|-1\) 次多项式(因为是这么多次的下降幂多项式)。如果 \(\lambda\) 有别的取值,经过一些无聊的转化后也能写成若干个组合数的线性组合,依然是 \(|\lambda|-1\) 次多项式。
使用拉格朗日插值计算即可。枚举 \(r\),计算 \(r,k+r,2k+r,\cdots,(|\lambda|-1)k\) 这些位置的系数,然后拿出来当作点值,需要计算下一项系数的时候做拉格朗日插值即可。
计算系数可以算出分母后求逆。
拉格朗日插值部分需要记忆化(?)。
代码
uoj 上过不去。懒得改。
#include <bits/stdc++.h>
using namespace std;
#ifdef LOCAL
#define debug(...) fprintf(stderr, ##__VA_ARGS__)
#else
#define endl "\n"
#define debug(...) void(0)
#endif
using LL = long long;
#define gcd __gcd
template <unsigned umod>
struct modint {/*{{{*/
static constexpr int mod = umod;
unsigned v;
modint() = default;
template <class T, enable_if_t<is_integral<T>::value, int> = 0>
modint(const T& y) : v((unsigned)(y % mod + (is_signed<T>() && y < 0 ? mod : 0))) {}
modint& operator+=(const modint& rhs) { v += rhs.v; if (v >= umod) v -= umod; return *this; }
modint& operator-=(const modint& rhs) { v -= rhs.v; if (v >= umod) v += umod; return *this; }
modint& operator*=(const modint& rhs) { v = (unsigned)(1ull * v * rhs.v % umod); return *this; }
modint& operator/=(const modint& rhs) { assert(rhs.v); return *this *= qpow(rhs, mod - 2); }
friend modint operator+(modint lhs, const modint& rhs) { return lhs += rhs; }
friend modint operator-(modint lhs, const modint& rhs) { return lhs -= rhs; }
friend modint operator*(modint lhs, const modint& rhs) { return lhs *= rhs; }
friend modint operator/(modint lhs, const modint& rhs) { return lhs /= rhs; }
template <class T> friend modint qpow(modint a, T b) {
modint r = 1;
for (assert(b >= 0); b; b >>= 1, a *= a) if (b & 1) r *= a;
return r;
}
friend int raw(const modint& self) { return self.v; }
friend ostream& operator<<(ostream& os, const modint& self) { return os << raw(self); }
explicit operator bool() const { return v != 0; }
modint operator-() const { return modint(0) - *this; }
bool operator==(const modint& rhs) const { return v == rhs.v; }
bool operator!=(const modint& rhs) const { return v != rhs.v; }
};/*}}}*/
using mint = modint<1000000021>;
template <int N>
struct C_prime {/*{{{*/
mint fac[N + 10], ifac[N + 10];
C_prime() {
fac[0] = 1;
for (int i = 1; i <= N; i++) fac[i] = fac[i - 1] * i;
ifac[N] = 1 / fac[N];
for (int i = N; i >= 1; i--) ifac[i - 1] = ifac[i] * i;
}
mint operator()(int n, int m) const { return n >= m ? fac[n] * ifac[m] * ifac[n - m] : 0; }
};/*}}}*/
C_prime<10000010> binom;
int n, m;
LL A[60], B[60];
int p[30], a[30], c[30];
mint calc(int typ) {
mint res = 1;
for (int i = 1; i <= m; i++) if (a[i]) {
auto coe = typ == +1 || i % 2 ? +1 : -1;
res *= qpow(mint(coe * i), a[i]) * binom.fac[a[i]];
}
return 1 / res;
}
vector<mint> inverse(const vector<mint> &F, int deg) {
vector<mint> G(deg);
mint inv_F0 = G[0] = 1 / F[0];
for (int i = 1; i < deg; i++) {
for (int j = 1; j <= min(i, (int)F.size() - 1); j++) G[i] -= G[i - j] * F[j];
G[i] *= inv_F0;
}
return G;
}
struct machine {
vector<mint> vec;
mint query(mint x) {
int m = (int)vec.size();
if (raw(x) < m) return vec[raw(x)];
mint res = 0, pre = 1;
for (int i = 0; i < m; i++) {
mint val = pre * binom.ifac[i] * binom.ifac[m - 1 - i] * (i & 1 ? -1 : +1);
pre *= i - x;
res = res * (i - x) + val * vec[i];
}
return res;
}
} mac[2010];
int L;
int vis[60][60], tim;
void init_macs(int cnt) {
++tim;
L = 1;
for (int i = 1; i <= cnt; i++) L = L / gcd(L, p[i]) * p[i];
vector<mint> F(m + 1);
F[0] = 1;
for (int i = 1; i <= cnt; i++) {
for (int j = m; j >= p[i]; j--) F[j] -= F[j - p[i]];
}
auto inv = inverse(F, L * m);
for (int r = 0; r < L; r++) {
mac[r].vec.resize(m);
for (int i = 0; i < m; i++) mac[r].vec[i] = inv[i * L + r];
}
}
mint query(LL pos) {
return pos < 0 ? 0 : mac[pos % L].query(pos / L);
}
mint mem[60][60];
mint query2(int j, int t) {
if (vis[j][t] < tim) vis[j][t] = tim, mem[j][t] = query(A[j] - t * (B[j] + 1));
return mem[j][t];
}
mint solve2() {
mint res = 1;
for (int i = 1; i <= m; i++) res *= binom(a[i], c[i]) * (c[i] & 1 ? -1 : +1);
vector<mint> F(m + 1);
F[0] = 1;
for (int i = 1; i <= m; i++) {
for (int t = 1; t <= c[i]; t++) {
for (int j = m; j >= i; j--) F[j] -= F[j - i];
}
}
for (int j = 1; j <= n; j++) {
mint ret = 0;
for (int t = 0; t <= min((LL)m, A[j] / (B[j] + 1)); t++)
//ret += F[t] * query(A[j] - t * (B[j] + 1));
ret += F[t] * query2(j, t);
res *= ret;
}
return res;
}
mint dfs2(int now) {
if (now > m) return solve2();
mint res = 0;
for (int &i = c[now] = 0; i <= a[now]; i++) res += dfs2(now + 1);
return res;
}
mint solve(int cnt) {
init_macs(cnt);
return dfs2(1);
}
mint ans1 = 0, ans2 = 0;
void dfs(int rest, int lst, int cnt) {
if (!rest) {
memset(a, 0, sizeof a);
for (int i = 1; i <= cnt; i++) a[p[i]] += 1;
auto res = solve(cnt);
ans1 += calc(-1) * res;
ans2 += calc(+1) * res;
return ;
}
++cnt;
for (int &i = p[cnt] = lst; i <= rest; i++) {
dfs(rest - i, i, cnt);
}
}
int main() {
#ifndef LOCAL
cin.tie(nullptr)->sync_with_stdio(false);
#endif
cin >> n >> m;
for (int i = 1; i <= n; i++) cin >> A[i];
for (int i = 1; i <= n; i++) cin >> B[i];
dfs(m, 1, 0);
cout << ans1 << " " << ans2 << endl;
return 0;
}
本文来自博客园,作者:caijianhong,转载请注明原文链接:https://www.cnblogs.com/caijianhong/p/18938859/solution-qoj4264