UOJ#54 BZOJ3434 [WC2014]时空穿梭

题目描述

小 X 驾驶着他的飞船准备穿梭过一个 \(n\) 维空间，这个空间里每个点的坐标可以用 \(n\) 个实数表示，即 \((x_1,x_2,\dots,x_n)\)。

为了穿过这个空间，小 X 需要在这个空间中选取 \(c\)（\(c\geq 2\)）个点作为飞船停留的地方，而这些点需要满足以下三个条件：

每个点的每一维坐标均为正整数，且第 \(i\) 维坐标不超过 \(m_i\)。

第 \(i+1\)（\(1\leq i<c\)）个点的第 \(j\)（\(1\leq j\leq n\)）维坐标必须严格大于第 \(i\) 个点的第 \(j\) 维坐标。

存在一条直线经过所选的所有点。在这个 \(n\) 维空间里，一条直线可以用 \(2n\) 个实数 \(p_1,p_2,\dots,p_n,v_1,v_2,\dots,v_n\) 表示。直线经过点 \((x_1,x_2,\dots,x_n)\)，当且仅当存在实数 \(t\)，使得对 \(i=1\dots n\) 均满足 \(x_i=p_i+v_it\)。

小 X 还没有确定他的最终方案，请你帮他计算一下一共有多少种不同的方案满足他的要求。由于答案可能会很大，你只需要输出答案 \(mod10007\) 后的值。

输入格式

第一行包含一个正整数 \(T\)，表示有 \(T\) 组数据需要求解。

每组数据包含两行，第一行包含两个正整数 \(n,c\)（\(c\geq 2\)），分别表示空间的维数和需要选择的暂停点的个数。

第二行包含 \(n\) 个正整数，依次表示 \(m_1,m_2,…,m_n\)。

输出格式

共 \(T\) 行，每行包含一个非负整数，依次对应每组数据的答案。

样例

input

output

explanation

对于第一组数据，共有两种可行的方案：一种选择 \((1,1),(2,2),(3,3)\)，另一种选择 \((1,2),(2,3),(3,4)\)。

限制与约定

\[T \leq 1000, n \leq 11, c \leq 20, m_i \leq 10^5
\]

时间限制: 1s, 空间限制: 512MB.

题解

首先，我们发现只要满足\(x_i < y_i, i = 1,2,\dots,n\) ，（其中\((x_1, x_2,\dots,x_n),(y_1, y_2,\dots,y_n)\)分别是第一个点和最后一个点的坐标）且点的不同排列算同一种即可。

那么，我们枚举这两个点之后，以其为端点的线段上共有（除去两端）\(\gcd(y_1-x_1,y_2-x_2,\dots.y_n-x_n)\)个点。其中选出\(c-2\)个点的方式有\(C_{\gcd(y_1-x_1,y_2-x_2,\dots.y_n-x_n)}^{c-2}\)种，所以总方案数是

\[ans=\sum_{\substack{1\leq x_i < y_i\leq m_i \\ i = 1,2,\dots,n}}C_{\gcd(y_1-x_1,y_2-x_2,\dots,y_n-x_n)}^{c-2}
\]

我们发现，所有\(p_i=y_i-x_i\)固定之后，方案数恰有

\[C_{\gcd(p_1, p_2,\dots,p_n)}^{c-2}\prod_{i=1}^n(m_i-p_i)
\]

种，这是因为\(p_i\)固定后\(x_i\)只能取\(1,2,\dots,m_i-p_i\).

于是我们有

\[ans=\sum_{\substack{1\leq q_i < m_i \\ i = 1,2,\dots,n}}C_{\gcd(p_1, p_2,\dots,p_n)}^{c-2}\prod_{i=1}^n(m_i-p_i)
\]

优先枚举\(d=\gcd(p_1, p_2,\dots,p_n)\)（下式中\(m_{min}\)表示\(\min(m_1,m_2,\dots,m_n)\):

\[ans=\sum_{d=1}^{m_{min}-1}C_d^{c-2}\sum_{\substack{1\leq q_i < m_i \\ i = 1,2,\dots,n\\ \gcd(p_1, p_2,\dots,p_n)=d}}\prod_{i=1}^n(m_i-p_i)
\]

令

\[f(d)=\sum_{\substack{1\leq q_i < m_i \\ i = 1,2,\dots,n\\ \gcd(p_1, p_2,\dots,p_n)=d}}\prod_{i=1}^n(m_i-p_i)
\]

则有

\[\begin{aligned}
\sum_{d\mid l}f(l)&=\sum_{\substack{1\leq q_i < m_i \\ i = 1,2,\dots,n\\d\mid\gcd(p_1, p_2,\dots,p_n)}}\prod_{i=1}^n(m_i-p_i)\\
&=\sum_{\substack{1\leq q_i'd < m_i \\ i = 1,2,\dots,n}}\prod_{i=1}^n(m_i-p_i'd)\\
&=\prod_{i=1}^n\sum_{p'=1}^{\left\lfloor\frac{m_i-1}d\right\rfloor}(m_i-p'd)\\
&=\prod_{i=1}^n\left[m_i\left\lfloor\frac{m_i-1}d\right\rfloor-\frac{d\left\lfloor\frac{m_i-1}d\right\rfloor\left(1+\left\lfloor\frac{m_i-1}d\right\rfloor\right)}2\right]
\end{aligned}\]

最后一步是等差数列求和。

由莫比乌斯反演可知：

\[f(d)=\sum_{d\mid l}\mu\left(\frac ld \right)\prod_{i=1}^n\left[m_i\left\lfloor\frac{m_i-1}d\right\rfloor-\frac{d\left\lfloor\frac{m_i-1}d\right\rfloor\left(1+\left\lfloor\frac{m_i-1}d\right\rfloor\right)}2\right]
\]

于是有

\[\begin{aligned}
ans&=\sum_{d=1}^{m_{min}-1}C_d^{c-2}\sum_{d\mid l}\mu\left(\frac ld \right)\prod_{i=1}^n\left[m_i\left\lfloor\frac{m_i-1}l\right\rfloor-\frac{l\left\lfloor\frac{m_i-1}l\right\rfloor\left(1+\left\lfloor\frac{m_i-1}l\right\rfloor\right)}2\right]\\
&=\sum_{l=1}^{m_{min}-1}\prod_{i=1}^n\left[m_i\left\lfloor\frac{m_i-1}l\right\rfloor-\frac{l\left\lfloor\frac{m_i-1}l\right\rfloor\left(1+\left\lfloor\frac{m_i-1}l\right\rfloor\right)}2\right]\sum_{d\mid l}\mu\left(\frac ld \right)C_d^{c-2}
\end{aligned}\]

易知\(\prod_{i=1}^n\left[m_i\left\lfloor\frac{m_i-1}l\right\rfloor-\frac{l\left\lfloor\frac{m_i-1}l\right\rfloor\left(1+\left\lfloor\frac{m_i-1}l\right\rfloor\right)}2\right]\)在将所有\(\left\lfloor\frac{m_i-1}l\right\rfloor\)作为常数看待的情况下是一个关于\(l\)的不超过\(n\)次（是因为\(mod10007\)之后有可能次数更低）的多项式，而且引起至少一个\(\left\lfloor\frac{m_i-1}l\right\rfloor\)改变的\(l\)有\(O(n\sqrt m)\)个，于是可以对每一个所有\(\left\lfloor\frac{m_i-1}l\right\rfloor\)都不改变的段分别求和。由于此时它是关于\(l\)的多项式，我们只需预处理出所有的\(S_{c,p,t}=\sum_{l=1}^tl^p\sum_{d\mid l}\mu\left(\frac ld \right)C_d^{c-2}\)即可。至于如何在只改变一个因式时快速维护成绩中各项的系数，可以用线段树实现（代码中为zkw线段树）。

完。

Code

#include <algorithm>

#include <cstdio>

const int N = 12;

const int C = 21;

const int M = 100001;

const int mod = 10007;

int fac[mod], inv[mod];

int mu[M], f[C][M], S[C][M][N];

int p[N * 1000], m[N];

inline int CC(int a, int b) {

  if (!b) return 1;

  if (a % mod < b % mod) return 0;

  return CC(a / mod, b / mod) * fac[a % mod] % mod * inv[fac[b % mod]] % mod * inv[fac[(a - b) % mod]] % mod;

}

int a[64][N], po[64];

int mdp[N];

int n, mm, c;

void upd(int x) {

  int l = x * 2, r = x * 2 + 1;

  po[x] = po[l] + po[r];

  for (int i = 0; i <= po[x]; ++i) a[x][i] = 0;

  for (int i = 0; i <= po[l]; ++i)

    for (int j = 0; j <= po[r]; ++j)

      a[x][i + j] = (a[x][i + j] + a[l][i] * a[r][j] % mod) % mod;

}

int sum(int l, int r) {

  int ans = 0;

  for (int i = 0; i <= po[1]; ++i)

    ans = (ans + a[1][i] * (S[c][r][i] - S[c][l - 1][i]) % mod) % mod;

  return ans;

}

int main() {

  fac[0] = 1;

  for (int i = 1; i < mod; ++i) fac[i] = fac[i - 1] * i % mod;

  inv[1] = 1;

  for (int i = 2; i < mod; ++i) inv[i] = mod - (mod / i) * inv[mod % i] % mod;

  mu[1] = 1;

  for (int i = 1; i < M; ++i)

    for (int j = i * 2; j < M; j += i)

      mu[j] -= mu[i];

  for (int c = 2; c < C; ++c) {

    for (int i = 1; i < M; ++i)

      f[c][i] = 0;

    for (int d = c - 1; d < M; ++d) {

      int cc = CC(d - 1, c - 2);

      for (int i = 1; i * d < M; ++i)

        f[c][i * d] = (f[c][i * d] + cc * mu[i]) % mod;

    }

    for (int j = 0; j < N; ++j)

      S[c][0][j] = 0;

    for (int i = 1; i < M; ++i)

      for (int j = 0, ij = 1; j < N; ++j, ij = ij * (i % mod) % mod)

        S[c][i][j] = (S[c][i - 1][j] + f[c][i] * ij) % mod;

  }

  int T;

  scanf("%d", &T);

  while (T--) {

    mm = 10000000;

    scanf("%d%d", &n, &c);

    for (int i = 0; i < n; ++i)

      scanf("%d", &m[i]), mm = std::min(mm, m[i]);

    if (mm < c) {

      printf("0\n");

      continue;

    } if (n == 1) {

      printf("%d\n", CC(m[0], c));

      continue;

    }

    int t = 0;

    p[t++] = 1;

    for (int i = 0; i < n; ++i) {

      int last = 1;

      while (last < m[i] - 1) {

        last = (m[i] - 1) / ((m[i] - 1) / last) + 1;

        p[t++] = last;

      }

    }

    std::sort(p, p + t);

    t = (int)(std::unique(p, p + t) - p);

    p[t] = mm;

    int d = 1;

    while (d < n) d <<= 1;

    for (int i = 0; i < n; ++i) mdp[i] = m[i];

    for (int i = 1; i <= d * 2; ++i) a[i][po[i] = 0] = 1;

    int ans = 0;

    for (int i = 0; p[i] < mm && i < t; ++i) {

      int l = p[i], r = p[i + 1] - 1; //[l, r]

      for (int j = 0; j < n; ++j)

        if ((m[j] - 1) / l != mdp[j]) {

          mdp[j] = (m[j] - 1) / l;

          int k = d + j;

          po[k] = 1;

          a[k][0] = m[j] % mod * (mdp[j] % mod) % mod;

          a[k][1] = -mdp[j] % mod * ((mdp[j] + 1) % mod) % mod * inv[2] % mod;

          while (k /= 2)

            upd(k);

        }

      ans = (ans + sum(l, r)) % mod;

    }

    printf("%d\n", (ans + mod) % mod);

  }

  return 0;

}

巴特西