题解: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\) 个。那么自然会先有一个多重集排列数,使元素带上标号:

\[\frac{m!}{\prod_{i=1}^m(i!)^{a_i}a_i!} \]

多重集排列数的同时要注意 \(a_i\) 个大小相同的连通块是互不区分的。除此之外还要乘上

\[\prod_{i=1}^mf_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)!\)

所以系数就是

\[\frac{1}{m!}\frac{m!\prod_{i=1}^m((-1)^{i-1}(i-1)!)^{a_i}}{\prod_{i=1}^m(i!)^{a_i}a_i!} =\frac{1}{\prod_{i=1}^m((-1)^{i-1}i)^{a_i}a_i!} \]

这样就够了。

贡献系数 - 子任务 2

这里我们使用 Burnside 引理。对于一个方案,我们只取单调不降的,其它和它重排相同的方案我们不要。看到这个重排相同我们考虑使用 Burnside 引理:

\[|X/G|=\frac{1}{|G|}\sum_{g\in G}|X^g| \]

在这里,我们枚举这个置换 \(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!}\frac{m!\prod_{i=1}^m((i-1)!)^{a_i}}{\prod_{i=1}^m(i!)^{a_i}a_i!} =\frac{1}{\prod_{i=1}^mi^{a_i}a_i!} \]

注意,\(\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\) 的约数,其容斥系数为:

\[\prod_{i=1}^m\binom{a_i}{c_i}(-1)^{c_i} \]

这时所有的质因子就独立了。假设这个质因子是 \(j\),用生成函数进行计数:

\[[x^{A_j}]\prod_{i=1}^m\left(\frac{1-x^{i(B_j+1)}}{1-x^i}\right)^{c_i}\left(\frac{1}{1-x^i}\right)^{a_i-c_i} \]

改写成

\[[x^{A_j}]\left(\prod_{i=1}^m\left(1-x^{i(B_j+1)}\right)^{c_i}\right)\left(\prod_{i=1}^m\left(\frac{1}{1-x^i}\right)^{a_i}\right) \]

前面是关于 \(x^{B_j+1}\) 的多项式,把它算出来之后,我们枚举它的所有 \(t\) 次项系数,计算

\[[x^{A_j-t(B_j+1)}]\prod_{i=1}^m\left(\frac{1}{1-x^i}\right)^{a_i} \]

就是要计算

\[\prod_{i=1}^{|\lambda|}\frac{1}{1-x^{\lambda_i}} \]

的某项系数,结论是,记 \(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;
}
posted @ 2025-06-20 19:09  caijianhong  阅读(81)  评论(2)    收藏  举报