Little Pony and Elements of Harmony

题目描述

The Elements of Harmony are six supernatural artifacts representing subjective aspects of harmony. They are arguably the most powerful force in Equestria. The inside of Elements of Harmony can be seen as a complete graph with $n$ vertices labeled from $0$ to $n-1$, where $n$ is a power of two, equal to $2^m$.
The energy in Elements of Harmony is in constant movement. According to the ancient book, the energy of vertex $u$ in time $i$ ($e_i[u]$) equals to:
$$
e_i[u]=\sum_v e_{i-1}[v] \cdot b[f(u, v)]
$$
Here $b[]$ is the transformation coefficient – an array of $m+1$ integers and $f(u, v)$ is the number of ones in the binary representation of number ($x \oplus y$).
Given the transformation coefficient and the energy distribution at time $0$ ($e_0[]$). Help Twilight Sparkle predict the energy distribution at time $t$ ($e_t[]$). The answer can be quite large, so output it modulo $p$.

题意概述

给定$m$。令$n=2^m$,$f(x, y)$表示$x \oplus y$(异或)在二进制下$1$的个数。给定一个长为$n$的数组$e_0$和一个长为$(m+1)$的数组$b$。有递推关系$e_i[u]=\sum_v e_{i-1}[v] \cdot b[f(u, v)]$。求$e_t$在模$p$意义下的值。
数据范围:$1 \le m \le 20, \; 0 \le t \le 10^{18}, \; 2 \le p \le 10^9, \; 1 \le e_0[i] \le 10^9, \; 0 \le b[i] \le 10^9$。

算法分析

令$g(x)$表示$x$在二进制下$1$的个数,$c[i]=b[g(i)]$。那么递推式可以表示为
$$
\begin{align}
e_i[u]&=\sum_v e_{i-1}[v] \cdot c[u \oplus v] \\
e_i[u \oplus v]&=\sum_{u, v} e_{i-1}[u] \cdot c[v]
\end{align}
$$
这是一个异或卷积的式子,可以用FWT快速求值。由于逆FWT时需要除以$2^m$,而在模$p$意义下$2^m$可能不存在逆元,所以需要在模$2^mp$意义下计算,最后将结果除以$2^m$。

代码

/*
 * The light at the end of the tunnel is the headlight of an approaching train.
 */

#include <cstdio>

#define int long long

void read(int &n) {
  char c;
  for (; (c = getchar()) < '0' || c > '9';)
    ;
  for (n = c - '0'; (c = getchar()) >= '0' && c <= '9'; (n *= 10) += c - '0')
    ;
}

static const int N = 1100000;
int p, b[25], c[N], e[N], f[N];

int mul(int a, int b) {
  return (a * b - (int)(((double)a * b + 0.5) / p) * p) % p;
}

int power(int a, int b) {
  int ret = 1;
  for (; b; b >>= 1)
    b & 1 && (ret = mul(ret, a)), a = mul(a, a);
  return ret;
}

void fwt(int *a, int n) {
  for (int i = 1; i < n; i <<= 1)
    for (int j = 0; j < n; j += i << 1)
      for (int k = 0; k < i; ++k) {
        int x = a[j + k], y = a[j + k + i];
        a[j + k] = (x + y) % p, a[j + k + i] = (p + x - y) % p;
      }
}

signed main() {
  int n, m, t;
  read(m), read(t), read(p), n = 1 << m, p <<= m;
  for (int i = 0; i < n; ++i)
    read(e[i]), f[i] = f[i >> 1] + (i & 1);
  for (int i = 0; i <= m; ++i)
    read(b[i]);
  for (int i = 0; i < n; ++i)
    c[i] = b[f[i]];
  fwt(e, n), fwt(c, n);
  for (int i = 0; i < n; ++i)
    e[i] = mul(e[i], power(c[i], t));
  fwt(e, n);
  for (int i = 0; i < n; ++i)
    printf("%lld\n", e[i] >> m);
  return 0;
}

The Sum of Unitary Divisors

题目描述

A natural number $d$ is a unitary divisor of $n$ if $d$ is a divisor of $n$ and if $d$ and ${n \over d}$ are coprime.
Let $\sigma^{\ast}(n)$ be the sum of the unitary divisors of $n$. For example, $\sigma^{\ast}(1)=1, \sigma^{\ast}(2)=3$ and $\sigma^{\ast}(6)=12$.
Let $S(n)=\sum_{i=1}^n \sigma^{\ast}(i)$.
Given $n$, find $S(n) \bmod 2^{64}$.

题意概述

定义$d$为$n$的单位因子当且仅当$d \mid n$且$(d, {n \over d})=1$。
令$\sigma^{\ast}(n)$为$n$的所有单位因子之和,$S(n)=\sum_{i=1}^n \sigma^{\ast}(i)$。给定$n$,求$S(n) \bmod 2^{64}$。
数据范围:$1 \le n \le 5 \times 10^{13}$。

算法分析

推式子
$$
\begin{align}
S(n)&=\sum_{i=1}^n \sigma^{\ast}(i) \\
&=\sum_{i=1}^n i \sum_{j=1}^{\lfloor n/i \rfloor} [(i, j)=1]
\end{align}
$$
利用莫比乌斯反演
$$
\begin{align}
S(n)&=\sum_{d=1}^n \mu(d) \sum_{i=1}^{\lfloor n/d \rfloor} di \sum_{j=1}^{\lfloor n/d^2i \rfloor} 1 \\
&=\sum_{d=1}^n d \cdot \mu(d) \sum_{i=1}^{\lfloor n/d \rfloor} i \cdot \left\lfloor {n \over d^2i} \right\rfloor
\end{align}
$$
显然,当$d^2 \gt n$时不产生贡献。那么
$$
\begin{align}
S(n)=\sum_{d=1}^{\lfloor \sqrt{n} \rfloor} d \cdot \mu(d) \sum_{i=1}^{\lfloor n/d^2 \rfloor} i \cdot \left\lfloor {n \over d^2i} \right\rfloor
\end{align}
$$
问题变成了求$F(n)=\sum_{i=1}^n i \cdot \mu(i)$和$G(n)=\sum_{i=1}^n i \cdot \left\lfloor {n \over i} \right\rfloor$。
对于$F$,我们有$(I \cdot \mu) \ast I=\varepsilon$。所以
$$
\begin{align}
F(n)=1-\sum_{i=2}^n i \cdot F\left(\left\lfloor {n \over i} \right\rfloor\right)
\end{align}
$$
这样$F$和$G$都可以用分块的技巧求出。而求$S$只要枚举$d$,将$\left\lfloor {n \over d^2} \right\rfloor$相等的$d$一起计算。

代码

#include <cmath>
#include <cstdio>

#define int long long

static const int N = 10000005;
static const int M = 50005;
int prime[N], mui[N], g[M];
bool vis[N];

void init() {
  int top = 0; mui[1] = 1;
  for (int i = 2; i < N; ++ i) {
    if (! vis[i]) prime[top ++] = i, mui[i] = -1;
    for (int j = 0; j < top; ++ j) {
      int k = i * prime[j]; if (k >= N) break; vis[k] = true;
      if (i % prime[j]) mui[k] = -mui[i]; else { mui[k] = 0; break; }
    }
  }
  for (int i = 2; i < N; ++ i) (mui[i] *= i) += mui[i - 1];
}

int get_sum(int n) {
  return n & 1 ? ((n + 1) >> 1) * n : (n >> 1) * (n + 1);
}

static const int HASH = 23456789;
int numr, h[HASH];
struct Record {
  int n, v, nxt;
} r[HASH];

void add_rec(int n, int v) {
  r[++ numr] = (Record) { n, v, h[n % HASH] }, h[n % HASH] = numr;
}

int get_f(int n) {
  if (n < N) return mui[n];
  for (int i = h[n % HASH]; i; i = r[i].nxt)
    if (r[i].n == n) return r[i].v;
  int ret = 1;
  for (int i = 2, j; i <= n; i = j + 1)
    j = n / (n / i), ret -= (get_sum(j) - get_sum(i - 1)) * get_f(n / i);
  return add_rec(n, ret), ret;
}

int get_g(int n) {
  if (n < M && g[n]) return g[n];
  int ret = 0;
  for (int i = 1, j; i <= n; i = j + 1)
    j = n / (n / i), ret += (get_sum(j) - get_sum(i - 1)) * (n / i);
  return ret;
}

int calc(int n) {
  int ret = 0;
  for (int i = 1, j; i * i <= n; i = j + 1)
    j = sqrt(n / (n / i / i) + 0.01), ret += (get_f(j) - get_f(i - 1)) * get_g(n / i / i);
  return ret;
}

signed main() {
  int T; scanf("%lld", &T), init();
  for (int i = 1; i < M; ++ i) g[i] = get_g(i);
  for (int n; T --;) scanf("%lld", &n), printf("%llu\n", calc(n));
  return 0;
}

Counting Divisors (Square)

题目描述

Let $\sigma_0(n)$ be the number of positive divisors of $n$.
For example, $\sigma_0(1)=1, \sigma_0(2)=2$ and $\sigma_0(6)=4$.
Let $S_2(n)=\sum_{i=1}^n \sigma_0(i^2)$.
Given $N$, find $S_2(N)$.

题意概述

令$\sigma_0(n)$表示$n$的正因数个数,$S_2(n)=\sum_{i=1}^n \sigma_0(i^2)$。给定$N$,求$S_2(N)$。
数据范围:$1 \le N \le 10^{12}$。

算法分析

为了方便,我们令$f(n)=\sigma_0(n^2)$,$\omega(n)$表示$n$的不同质因子个数。
设$n=\prod_{k=1}^{\omega(n)} p_k^{a_k}$,其中$p_k$表示不同质数。
先来推$f(n)$。对于$n$的每个因数$d$,将$d$的某些质因子$p_k$的次数加上$a_k$得到$d’$,有$d’ \nmid n, d’ \mid n^2$。这样的方案数有$(2^{\omega(d)}-1)$种,再算上$d$本身,可得
$$
\begin{align}
f(n)=\sum_{d \mid n} 2^{\omega(d)}
\end{align}
$$
化简$2^{\omega(d)}$。它等价于$d$的每个质因数要么选一个,要么不选。那么它就等于$d$的因数中square-free的数的个数。因此
$$
\begin{align}
2^{\omega(d)}=\sum_{i \mid d} \mu^2(i)
\end{align}
$$
回到原式
$$
\begin{align}
S_2(n)=\sum_{i=1}^n \sum_{j \mid i} \sum_{k \mid j} \mu^2(k)
\end{align}
$$
还是毫无头绪。我们继续化简$f(n)$。考虑它的卷积形式$f=\mu^2 \ast 1 \ast 1=\mu^2 \ast \sigma_0$。所以
$$
\begin{align}
S_2(n)&=\sum_{i=1}^n \sum_{j \mid i} \mu^2(j) \sigma_0\left({i \over j}\right) \\
&=\sum_{ij \le n} \mu^2(i) \sigma_0(j)
\end{align}
$$
那么只要用反演求出$\mu^2$的前缀和、用分块求出$\sigma_0$的前缀和,就可以利用分块的技巧了。具体如下
$$
\begin{align}
\sum_{i=1}^n \mu^2(i)&=\sum_{i=1}^n [i\text{最大平方因子}=1] \\
&=\sum_{i=1}^n \sum_{j^2 \mid i} \mu(j) \\
&=\sum_{i=1}^{\lfloor \sqrt{n} \rfloor} \mu(i) \left\lfloor {n \over i^2} \right\rfloor \\
\sum_{i=1}^n \sigma_0(i)&=\sum_{i=1}^n \left\lfloor {n \over i} \right\rfloor
\end{align}
$$

代码

#include <cmath>
#include <cstdio>
#include <cstring>

void read(long long &t) {
  char c; while ((c = getchar()) < '0' || c > '9') ; t = c - '0';
  while ((c = getchar()) >= '0' && c <= '9') (t *= 10) += c - '0';
}

int N;
char mui[100000000];
int d[100000000], sum[100000000], prime[10000000];
bool vis[100000000];
long long n, sigma[100000000];

void init() {
  int top = 0; mui[1] = 1, sigma[1] = 1;
  for (int i = 2; i < N; ++ i) {
    if (! vis[i]) prime[top ++] = i, mui[i] = -1, sigma[i] = 2, d[i] = 1;
    for (int j = 0; j < top; ++ j) {
      int k = i * prime[j]; if (k >= N) break;
      vis[k] = true;
      if (i % prime[j]) mui[k] = -mui[i], sigma[k] = sigma[i] << 1, d[k] = 1;
      else { mui[k] = 0, sigma[k] = sigma[i] / (d[i] + 1) * (d[i] + 2), d[k] = d[i] + 1; break; }
    }
  }
  for (int i = 1; i < N; ++ i) sum[i] = sum[i - 1] + mui[i] * mui[i], sigma[i] += sigma[i - 1];
}

inline long long get_mui(long long n) {
  if (n < N) return sum[n];
  int q = sqrt(n); long long ret = 0;
  for (register int i = 1; i <= q; ++ i) if (mui[i]) ret += mui[i] * (n / i / i);
  return ret;
}

inline long long get_sigma(long long n) {
  if (n < N) return sigma[n];
  long long ret = 0;
  for (register long long i = 1, j; i <= n; i = j + 1)
    j = n / (n / i), ret += (j - i + 1) * (n / i);
  return ret;
}

long long calc(long long n) {
  int q = sqrt(n); long long ret = 0;
  for (register int i = 1; i <= q; ++ i) if (mui[i]) ret += get_sigma(n / i);
  for (register long long i = q + 1, j, cur, last = sum[q]; i <= n; i = j + 1, last = cur)
    j = n / (n / i), cur = get_mui(j), ret += (cur - last) * get_sigma(n / i);
  return ret;
}

signed main() {
  long long T; read(T);
  if (T > 800) N = 20000; else N = 100000000;
  init();
  while (T --) read(n), printf("%lld\n", calc(n));
  return 0;
}

GCD Extreme (Hard)

题目描述

Let $G(n)=\sum_{i=1}^n \sum_{j=i+1}^n (i, j)$.
For example, $G(1)=0, G(2)=(1, 2)=1, G(3)=(1, 2)+(1, 3)+(2, 3)=3$.
Given $N$, find $G(N)$ modulo $2^{64}$.

题意概述

令$G(n)=\sum_{i=1}^n \sum_{j=i+1}^n (i, j)$。给定$N$,求$G(N) \bmod 2^{64}$。
数据范围:$1 \le N \le 235711131719$。

算法分析

推式子
$$
\begin{align}
G(n)&=\sum_{i=1}^n \sum_{j=i+1}^n (i, j)=\sum_{i=2}^n \sum_{j=1}^{i-1} (i, j) \\
&=\sum_{d=1}^n d \sum_{i=2}^{\lfloor n/d \rfloor} \sum_{j=1}^{i-1} [(i, j)=1] \\
&=\sum_{d=1}^n d \sum_{i=2}^{\lfloor n/d \rfloor} \varphi(i)
\end{align}
$$
令$S(n)=\sum_{i=1}^n \varphi(i)$。由于$\varphi \ast 1=I$,因此
$$
\begin{align}
\sum_{i=1}^n (\varphi \ast 1)(i)&={n(n+1) \over 2} \\
&=\sum_{i=1}^n \sum_{j \mid i} \varphi(j) \\
&=\sum_{ij \le n} \varphi(j) \\
&=\sum_{i=1}^n \sum_{j=1}^{\lfloor n/i \rfloor} \varphi(j) \\
&=\sum_{i=1}^n S\left(\left\lfloor {n \over i} \right\rfloor\right)
\end{align}
$$
那么
$$
\begin{align}
S(n)&={n(n+1) \over 2}-\sum_{i=2}^n S\left(\left\lfloor {n \over i} \right\rfloor\right) \\
G(n)&=\sum_{d=1}^n d \cdot \left(S\left(\left\lfloor {n \over d} \right\rfloor\right)-1\right)
\end{align}
$$
接下来就可以利用分块的技巧计算。

代码

#include <map>
#include <cstdio>

#define int long long

static const int N = 10000000;
int prime[N], phi[N];
bool vis[N];

void init() {
  int top = 0; phi[1] = 1;
  for (int i = 2; i < N; ++ i) {
    if (! vis[i]) prime[top ++] = i, phi[i] = i - 1;
    for (int j = 0; j < top; ++ j) {
      int k = i * prime[j]; if (k >= N) break; vis[k] = true;
      if (i % prime[j]) phi[k] = phi[i] * (prime[j] - 1); else { phi[k] = phi[i] * prime[j]; break; }
    }
  }
  for (int i = 2; i < N; ++ i) phi[i] += phi[i - 1];
}

int get_sum(int n) { return n & 1 ? ((n + 1) >> 1) * n : (n >> 1) * (n + 1); }

std :: map <int, int> mp;

int get_phi(int n) {
  if (n < N) return phi[n];
  if (mp.count(n)) return mp[n];
  int ret = get_sum(n);
  for (int i = 2, j; i <= n; i = j + 1)
    j = n / (n / i), ret -= (j - i + 1) * get_phi(n / i);
  return mp[n] = ret;
}

int calc(int n) {
  int ret = 0;
  for (int i = 1, j; i <= n; i = j + 1)
    j = n / (n / i), ret += (get_sum(j) - get_sum(i - 1)) * (get_phi(n / i) - 1);
  return ret;
}

signed main() {
  int T; scanf("%lld", &T), init();
  for (int n; T --;) scanf("%lld", &n), printf("%llu\n", calc(n));
  return 0;
}

The Child and Binary Tree

题目描述

Our child likes computer science very much, especially he likes binary trees.
Consider the sequence of $n$ distinct positive integers: $c_1, c_2, \ldots, c_n$. The child calls a vertex-weighted rooted binary tree good if and only if for every vertex $v$, the weight of $v$ is in the set ${c_1, c_2, \ldots, c_n}$. Also our child thinks that the weight of a vertex-weighted tree is the sum of all vertices’ weights.
Given an integer $m$, can you for all $s \; (1 \le s \le m)$ calculate the number of good vertex-weighted rooted binary trees with weight $s$? Please, check the samples for better understanding what trees are considered different.
We only want to know the answer modulo $998244353$ ($7 \times 17 \times 2^{23}+1$, a prime number).

题意概述

有$n$种点,每种点有无限个,第$i$种点的权值为$c_i$。定义一棵二叉树的权值等于它所有点的权值之和。求对于所有$s \in [1, m]$,权值为$s$的二叉树有几棵。两棵二叉树不同当且仅当它们左子树或右子树不同,或者根节点权值不同。
数据范围:$1 \le n, m, c_i \le 10^5$。

算法分析

令$f(x)$表示权值为$x$的二叉树个数,$F(x)$为其生成函数($F(x)=\sum_{i \ge 0} f(i)x^i$)。
令$C(x)$为给定$c$的集合的生成函数($C(x)=\sum_{i=1}^n x^{c_i}$)。
根据DP转移方程,易知
$$
f(x)=\sum_{w \in {c_1, c_2, \ldots, c_n}} \sum_{i=0}^{x-w} f(i)f(x-w-i)
$$

$$
F(x)=C(x)F(x)^2+1
$$
解得
$$
F(x)={1 \pm \sqrt{1-4C(x)} \over 2C(x)}={2 \over 1 \pm \sqrt{1-4C(x)}}
$$
显然,若取减号,则当$x$趋近$0$时分母为$0$,因此只能取加号。接着就是多项式开根和多项式求逆了。

  • 多项式求逆:
    求$GF \equiv 1 \pmod {x^n}$。

    假设已知$G_0F \equiv 1 \pmod {x^{\lceil n/2 \rceil}}$
    $G-G_0 \equiv 0 \pmod {x^{\lceil n/2 \rceil}}$
    $G^2-2GG_0+G_0^2 \equiv 0 \pmod {x^n}$
    $G-2G_0+G_0^2F \equiv 0 \pmod {x^n}$
    $G \equiv 2G_0-G_0^2F \pmod {x^n}$

  • 多项式开根:
    求$G^2 \equiv F \pmod {x^n}$。

    假设已知$G_0^2 \equiv F \pmod {x^{\lceil n/2 \rceil}}$
    $(G_0^2-F)^2 \equiv 0 \pmod {x^n}$
    $(G_0^2+F)^2 \equiv 4G_0^2F \pmod {x^n}$
    $\left({G_0^2+F \over 2G_0}\right)^2 \equiv F \pmod {x^n}$
    $G \equiv {G_0+G_0^{-1}F \over 2} \pmod {x^n}$

代码

#include <algorithm>
#include <cstdio>
#include <cstring>

static const int N = 500000;
static const int MOD = 998244353;
static const int G = 3;
static const int INV2 = 499122177;
int n, m, c[N], C[N], rev[N], wn[N], tmp[N], tmp2[N], tmp3[N];

int power(int a, int b) {
  int ret = 1;
  for (a %= MOD, b %= MOD - 1; b; b >>= 1)
    b & 1 && (ret = 1ll * ret * a % MOD), a = 1ll * a * a % MOD;
  return ret;
}

void init(int &n) {
  int m = n << 1, l = 0;
  for (n = 1; n < m; n <<= 1, ++l)
    ;
  for (int i = 1; i < n; ++i)
    rev[i] = rev[i >> 1] >> 1 | (i & 1) << (l - 1);
}

void ntt(int *a, int n, bool inv) {
  for (int i = 0; i < n; ++i)
    if (i < rev[i])
      std::swap(a[i], a[rev[i]]);
  wn[0] = 1, wn[1] = power(G, (MOD - 1) / n);
  for (int i = 2; i < n >> 1; ++i)
    wn[i] = 1ll * wn[i - 1] * wn[1] % MOD;
  for (int i = 1; i < n; i <<= 1)
    for (int j = 0; j < n; j += i << 1)
      for (int k = 0; k < i; ++k) {
        int x = a[j + k], y = 1ll * wn[n / (i << 1) * k] * a[j + k + i] % MOD;
        a[j + k] = (x + y) % MOD, a[j + k + i] = (MOD + x - y) % MOD;
      }
  if (inv) {
    for (int i = 1; i < n >> 1; ++i)
      std::swap(a[i], a[n - i]);
    int rec = power(n, MOD - 2);
    for (int i = 0; i < n; ++i)
      a[i] = 1ll * a[i] * rec % MOD;
  }
}

void get_inv(int *f, int *g, int n) {
  if (n == 1)
    return void(g[0] = power(f[0], MOD - 2));
  int rec = n;
  get_inv(f, g, (n + 1) >> 1), init(n);
  for (int i = (rec + 1) >> 1; i < n; ++i)
    g[i] = 0;
  for (int i = 0; i < rec; ++i)
    tmp[i] = f[i];
  for (int i = rec; i < n; ++i)
    tmp[i] = 0;
  ntt(g, n, 0), ntt(tmp, n, 0);
  for (int i = 0; i < n; ++i)
    g[i] = 1ll * g[i] * (MOD + 2 - 1ll * g[i] * tmp[i] % MOD) % MOD;
  ntt(g, n, 1);
  for (int i = rec; i < n; ++i)
    g[i] = 0;
}

void get_sqrt(int *f, int *g, int n) {
  if (n == 1)
    return void(g[0] = 1);
  int rec = n;
  get_sqrt(f, g, (n + 1) >> 1);
  for (int i = 0; i<(n + 1)>> 1; ++i)
    tmp2[i] = g[i];
  for (int i = (n + 1) >> 1; i < n; ++i)
    tmp2[i] = 0;
  get_inv(tmp2, tmp3, n), init(n);
  for (int i = (rec + 1) >> 1; i < n; ++i)
    g[i] = 0;
  for (int i = 0; i < rec; ++i)
    tmp2[i] = f[i];
  for (int i = rec; i < n; ++i)
    tmp2[i] = 0;
  ntt(tmp2, n, 0), ntt(tmp3, n, 0);
  for (int i = 0; i < n; ++i)
    tmp3[i] = 1ll * tmp3[i] * tmp2[i] % MOD;
  ntt(tmp3, n, 1);
  for (int i = 0; i < rec; ++i)
    g[i] = 1ll * (g[i] + tmp3[i]) * INV2 % MOD;
  for (int i = rec; i < n; ++i)
    g[i] = 0;
}

int main() {
  scanf("%d%d", &n, &m);
  for (int i = 0; i < n; ++i)
    scanf("%d", &c[i]);
  for (int i = 0; i < n; ++i)
    if (c[i] <= m)
      ++C[c[i]];
  for (int i = 1; i <= m; ++i)
    C[i] = (MOD - (C[i] << 2)) % MOD;
  C[0] = 1;
  get_sqrt(C, c, m + 1), ++c[0], get_inv(c, C, m + 1);
  for (int i = 1; i <= m; ++i)
    printf("%d\n", (C[i] << 1) % MOD);
  return 0;
}

序列求和 V4

题意概述

给定$n$和$k$,求$\sum_{i=1}^n i^k \bmod 10^9+7$。有$T$组数据。
数据范围:$1 \le T \le 500, \; 1 \le n \le 10^{18}, \; 1 \le k \le 50000$。

算法分析1

设答案为$f(n)$。这是一个$(k+1)$次多项式,只要确定$(k+2)$个点就可以确定$f$。
我们可以令$x=1\ldots k+2$,在$O(k\log k)$时间内计算出它们对应的$y$。
利用拉格朗日插值法。假设有$n$个点$(x_i, y_i)$,那么
$$f(x)=\sum_{i=1}^n {\prod_{j \in [1, n] \land j \neq i} (x-x_j) \over \prod_{j \in [1, n] \land j \neq i} (x_i-x_j)}y_i$$
考虑分子部分。在给定$x$的情况下,这一部分可以$O(k)$预处理$O(1)$求值。
考虑分母部分。它形如$\prod_{j \in [i-k-2, i-1] \land j \neq 0} j$,可以在插值时$O(1)$维护。
总时间复杂度为$O(Tk\log k)$。

代码1

#include <cstdio>

#define int long long

void read(int &n) {
  char c; while ((c = getchar()) < '0' || c > '9') ; n = c - '0';
  while ((c = getchar()) >= '0' && c <= '9') (n *= 10) += c - '0';
}

static const int K = 50005;
static const int MOD = 1000000007;
int n, k, a[K], p[K], q[K], inv[K], base[K];

int power(int a, int b) {
  int ret = 1; a %= MOD, b %= MOD - 1;
  while (b) b & 1 && ((ret *= a) %= MOD), (a *= a) %= MOD, b >>= 1;
  return ret;
}

int calc(int n) {
  if (n <= k + 2) return a[n]; n %= MOD;
  int w = power(base[k + 2], MOD - 2), ans = 0;
  p[0] = q[k + 3] = 1;
  for (int i = 1; i <= k + 2; ++ i) p[i] = p[i - 1] * (n - i) % MOD;
  for (int i = k + 2; i; -- i) q[i] = q[i + 1] * (n - i) % MOD;
  for (int i = 1; i <= k + 2; ++ i) (ans += a[i] * w % MOD * p[i - 1] % MOD * q[i + 1]) %= MOD, w = w * (i - k - 2) % MOD * inv[i] % MOD;
  return (ans + MOD) % MOD;
}

signed main() {
  int T; read(T), base[1] = 1;
  for (int i = 1; i < K; ++ i) inv[i] = power(i, MOD - 2);
  for (int i = 2; i < K; ++ i) base[i] = base[i - 1] * (MOD + 1 - i) % MOD;
  while (T --) {
    read(n), read(k);
    for (int i = 1; i < k + 3; ++ i) a[i] = (a[i - 1] + power(i, k)) % MOD;
    printf("%lld\n", calc(n));
  }
  return 0;
}

算法分析2

此题也可以用伯努利数解决,因为
$$\sum_{i=1}^n i^k={1 \over k+1} \sum_{i=0}^k (-1)^i {k+1 \choose i} B_in^{k+1-i}$$
其中伯努利数由其指数型生成函数${x \over e^x-1}$定义。易知
$$\sum_{i=0}^{\infty} B_i{x^i \over i!}={x \over e^x-1}={x \over \sum_{i=1}^{\infty} {x^i \over i!}}={1 \over \sum_{i=0}^{\infty} {x^i \over (i+1)!}}$$
只要对分母求多项式逆元,即可得到伯努利数,时间复杂度$O(k\log k)$。其它部分均可$O(k)$预处理,答案也可以$O(k)$计算。因此总时间复杂度是$O(k\log k+Tk)$。

代码2

#include <algorithm>
#include <cmath>
#include <cstdio>
#include <cstring>

#define int long long

static int const N = 200000;
static int const MOD = 1000000007;
static int const M = 31622;
static double const PI = acos(-1);
int b[N], rev[N], fac[N], inv[N];
class complex {
private:
  double x, y;

public:
  complex(double _x = 0, double _y = 0) : x(_x), y(_y) {}
  double real() { return x; }
  friend complex operator+(complex const &a, complex const &b) {
    return complex(a.x + b.x, a.y + b.y);
  }
  friend complex operator-(complex const &a, complex const &b) {
    return complex(a.x - b.x, a.y - b.y);
  }
  friend complex operator*(complex const &a, complex const &b) {
    return complex(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x);
  }
} wn[N], A[N], B[N], C[N], D[N], E[N], F[N], G[N];

int power(int a, int b) {
  int ret = 1;
  for (; b; b >>= 1)
    b & 1 && ((ret *= a) %= MOD), (a *= a) %= MOD;
  return ret;
}

int init(int m) {
  int n = 1, l = 0;
  for (; n <= m; n <<= 1, ++l)
    ;
  for (int i = 1; i < n; ++i)
    rev[i] = rev[i >> 1] >> 1 | (i & 1) << l - 1;
  for (int i = 0; i < n >> 1; ++i)
    wn[i] = complex(cos(2 * PI / n * i), sin(2 * PI / n * i));
  return n;
}

void fft(complex *a, int n, bool inv = 0) {
  for (int i = 0; i < n; ++i)
    if (i < rev[i])
      std::swap(a[i], a[rev[i]]);
  for (int i = 1; i < n; i <<= 1)
    for (int j = 0; j < n; j += i << 1)
      for (int k = 0; k < i; ++k) {
        complex x = a[j + k], y = wn[n / (i << 1) * k] * a[j + k + i];
        a[j + k] = x + y, a[j + k + i] = x - y;
      }
  if (inv) {
    std::reverse(a + 1, a + n);
    for (int i = 0; i < n; ++i)
      a[i] = a[i].real() / n;
  }
}

int get_mod(double t) {
  return (int)(fmod(t, MOD) + MOD + 0.5) % MOD;
}

void get_inv(int *f, int *g, int n) {
  if (n == 1)
    return void(g[0] = power(f[0], MOD - 2));
  get_inv(f, g, n + 1 >> 1);
  int len = init(n << 1);
  for (int i = 0; i < n + 1 >> 1; ++i)
    A[i] = g[i] / M, B[i] = g[i] % M;
  for (int i = n + 1 >> 1; i < len; ++i)
    A[i] = B[i] = 0;
  for (int i = 0; i < n; ++i)
    C[i] = f[i] / M, D[i] = f[i] % M;
  for (int i = n; i < len; ++i)
    C[i] = D[i] = 0;
  fft(A, len), fft(B, len), fft(C, len), fft(D, len);
  for (int i = 0; i < len; ++i) {
    E[i] = 0 - A[i] * C[i];
    F[i] = 0 - A[i] * D[i] - B[i] * C[i];
    G[i] = 2 - B[i] * D[i];
  }
  fft(E, len, 1), fft(F, len, 1), fft(G, len, 1);
  for (int i = 0; i < n; ++i) {
    int x = get_mod(E[i].real()) * M % MOD * M % MOD;
    int y = get_mod(F[i].real()) * M % MOD;
    int z = get_mod(G[i].real());
    int w = (x + y + z) % MOD;
    C[i] = w / M, D[i] = w % M;
  }
  for (int i = n; i < len; ++i)
    C[i] = D[i] = 0;
  fft(C, len), fft(D, len);
  for (int i = 0; i < len; ++i) {
    E[i] = A[i] * C[i];
    F[i] = A[i] * D[i] + B[i] * C[i];
    G[i] = B[i] * D[i];
  }
  fft(E, len, 1), fft(F, len, 1), fft(G, len, 1);
  for (int i = 0; i < n; ++i) {
    int x = get_mod(E[i].real()) * M % MOD * M % MOD;
    int y = get_mod(F[i].real()) * M % MOD;
    int z = get_mod(G[i].real());
    g[i] = (x + y + z) % MOD;
  }
}

int get_c(int n, int m) {
  return fac[n] * inv[m] % MOD * inv[n - m] % MOD;
}

signed main() {
  int T;
  fac[0] = 1;
  for (int i = 1; i < N; ++i)
    fac[i] = fac[i - 1] * i % MOD;
  inv[N - 1] = power(fac[N - 1], MOD - 2);
  for (int i = N - 1; i; --i)
    inv[i - 1] = inv[i] * i % MOD;
  get_inv(inv + 1, b, 50001);
  for (int i = 0; i < 50001; ++i)
    (b[i] *= fac[i]) %= MOD;
  for (scanf("%lld", &T); T--;) {
    int n, k, ans = 0;
    scanf("%lld%lld", &n, &k), n %= MOD;
    for (int i = k, f = k & 1 ? -1 : 1, p = n; ~i; --i, f = -f, (p *= n) %= MOD)
      ans += (MOD + f) * get_c(k + 1, i) % MOD * b[i] % MOD * p % MOD;
    printf("%lld\n", ans % MOD * power(k + 1, MOD - 2) % MOD);
  }
  return 0;
}

Panda’s Birthday Present

题目描述

On Panda’s Birthday party, he received a strange present from Jason. The present is a black box with $4$ dices in it which is used to play a game. The dice in the box is unusual. Instead of the digits, only red or blue is painted on each side of the dice. Before the first round of the game, the box can repaint every side of every dice red or blue with equal probability. Then for each round of the game, the box will roll the $4$ dices and tell the player the number of red side facing up, which is the point the player get. Now, Panda has play it for two rounds and he tell you the point he has got for each round. Can you tell him the expected point he can get for next round.

题意概述

一个盒子里有四个骰子,初始时它们每一面都被等概率染成红色或蓝色(你并不知道具体的染色方案)。每次你会摇动盒子,然后观察有几个骰子朝上的面是红色的。给定前两次观察到的值$p, q$,求第三次观察到的值的期望。

算法分析

为了方便起见,我们假设$0^0=1$。首先可以知道的结论:

  1. 对于一个骰子,它有$t$个面被染成红色的概率$P(x=t)={{6 \choose t} \over 2^6}={{6 \choose t} \over 64}$;
  2. 扔一个骰子$n$次,其中有$k$次朝上的面是红色的概率为$\sum_{i=0}^6 P(x=i) \cdot {i^k(6-i)^{n-k} \over 6^n}$。

考虑不管是一次操作中的多个骰子还是多次操作,它们其实是独立重复实验,与顺序无关。因此,只要$p+q$是确定的,那么第三次的期望也是确定的。令$1$表示红色面朝上、$0$表示蓝色面朝上,那么我们所要求的就是$P(??1|??)$,其中$??$表示前两次骰子朝上的面分别是什么。
根据贝叶斯定理,$P(H|E)={P(H)P(E|H) \over P(E)}$,所以$P(111|11)={P(111)P(11|111) \over P(11)}$。经计算可得
$$
\begin{align}
P(111|11)&={9 \over 14} \\
P(101|10)&={1 \over 2} \\
P(001|00)&={5 \over 14}
\end{align}
$$
再次考虑独立重复实验。因为它们可以任意交换顺序,所以可以直接将前两次朝上的面是红色的骰子两两并在一起算,朝上的面是蓝色的骰子两两并在一起算,若有剩下的一对红色蓝色则将它们并在一起。
具体来说,若$2 \mid p+q$,则期望为
$${p+q \over 2} \cdot {9 \over 14}+\left(4-{p+q \over 2}\right) \cdot {5 \over 14}={p+q+10 \over 7}$$
否则,期望为
$${p+q-1 \over 2} \cdot {9 \over 14}+\left(4-{p+q+1 \over 2}\right) \cdot {5 \over 14}+{1 \over 2}={p+q+10 \over 7}$$

Max and Min

题目描述

Two kittens, Max and Min, play with a pair of non-negative integers $x$ and $y$. As you can guess from their names, kitten Max loves to maximize and kitten Min loves to minimize. As part of this game Min wants to make sure that both numbers, $x$ and $y$ became negative at the same time, and kitten Max tries to prevent him from doing so.
Each kitten has a set of pairs of integers available to it. Kitten Max has $n$ pairs of non-negative integers $(a_i, b_i)$, and kitten Min has $m$ pairs of non-negative integers $(c_j, d_j)$. As kitten Max makes a move, it can take any available pair $(a_i, b_i)$ and add $a_i$ to $x$ and $b_i$ to $y$, and kitten Min can take any available pair $(c_j, d_j)$ and subtract $c_j$ from $x$ and $d_j$ from $y$. Each kitten can use each pair multiple times during distinct moves.
Max moves first. Kitten Min is winning if at some moment both numbers $a, b$ are negative simultaneously. Otherwise, the winner of the game is kitten Max. Determine which kitten wins if both of them play optimally.

题意概述

有两个人Max和Min,初始时有两个数$x$和$y$。Max有$n$个数对$(a_i, b_i)$,Min有$m$个数对$(c_j, d_j)$。两个人轮流操作,Max先手。Max可以选择一个数对$(a_i, b_i)$,将$x$加上$a_i$、$y$加上$b_i$;Min可以选择一个数对$(c_j, d_j)$,将$x$减去$c_j$、$y$减去$d_j$。若在某一时刻$x, y$同时为负,则Min获胜,否则Max获胜。求谁必胜。
数据范围:$1 \le n, m \le 10^5, \; 1 \le x, y, a_i, b_i, c_j, d_j \le 10^9$。

算法分析

对于两个向量$v_1, v_2$,我们用$(v_1, v_2)$表示$v_1 \cdot v_2$。
我们把所有数对看成从原点出发的向量。假设有一个新的向量$(a, b)$。若
$$
\begin{align}
&\exists (a, b), \; a, b \ge 0 \land a+b \gt 0, \\
&\exists (a_i, b_i), \; \forall (c_j, d_j), \; ((a_i, b_i), (a, b)) \ge ((c_j, d_j), (a, b))
\end{align}
$$
那么Max必胜。证明如下:

$\because$初始时$((x, y), (a, b))=ax+by \gt 0$
$((a_i-c_j, b_i-d_j), (a, b))=((a_i, b_i), (a, b))-((c_j, d_j), (a, b)) \ge 0$
$\therefore$在任意时刻$((x’, y’), (a, b)) \gt 0$
$\therefore x’, y’$不全为负数

由此得证。
若$\exists (c_{j_0}, d_{j_0}), (c_{j_1}, d_{j_1})$,点$(a_i, b_i)$严格在点$(c_{j_0}, d_{j_0})$、$(c_{j_1}, d_{j_1})$和原点所组成的三角形内部,那么根据上述证明,$(a_i, b_i)$对Max没有贡献。特别的,若
$$\exists (c_j, d_j), \; (a_i, b_i) \parallel (c_j, d_j) \land |(a_i, b_i)| \lt |(c_j, d_j)|$$
那么$(a_i, b_i)$也没有贡献。
算法已经很清晰了——对Min的向量求凸包,判断是否所有Max的向量都严格在凸包内部。假设Min的向量中$c_j$的最大值为$mx_c$,$d_j$的最大值为$mx_d$,凸包中还需加入三个点$(0, 0), (mx_c, 0), (0, mx_d)$。
具体实现时,可以将Max的向量拿来一起求凸包。若求出来的凸包上有Max的向量,则Max必胜,否则Min必胜。

代码

#include <vector>
#include <cstdio>
#include <algorithm>

#define int long long

typedef const int ci;

inline int max(ci &a, ci &b) {
  return a > b ? a : b;
}

static ci N = 100005;
int st[N << 1];
std :: vector <int> rec;
struct Point {
  int type, x, y;
  Point(ci &_x = 0, ci &_y = 0) : type(0), x(_x), y(_y) {}
  bool operator == (const Point &p) const {
    return type == p.type && x == p.x && y == p.y;
  }
  bool operator < (const Point &p) const {
    return x == p.x ? y < p.y : x < p.x;
  }
  Point operator - (const Point &p) const {
    return Point(x - p.x, y - p.y);
  }
  int operator * (const Point &p) const {
    return x * p.y - y * p.x;
  }
} p[N << 1];

signed main() {
  int n, m, tot, mxa = 0, mxb = 0; scanf("%lld%lld%*d%*d", &n, &m), tot = n + m + 3;
  for (int i = 0; i < n; ++ i) scanf("%lld%lld", &p[i].x, &p[i].y), p[i].type = 1;
  for (int i = n; i < n + m; ++ i) scanf("%lld%lld", &p[i].x, &p[i].y), p[i].type = 2, mxa = max(mxa, p[i].x), mxb = max(mxb, p[i].y);
  p[n + m] = Point(0, 0), p[n + m + 1] = Point(mxa, 0), p[n + m + 2] = Point(0, mxb), std :: sort(p, p + tot);
  for (int i = 0, j; i < tot; i = j) {
    j = i; int rec = 0;
    while (j < tot && p[i].x == p[j].x && p[i].y == p[j].y) ++ j;
    for (int k = i; k < j; ++ k) rec |= p[k].type;
    if (rec & 1) for (int k = i; k < j; ++ k) p[k].type = 1;
  }
  tot = std :: unique(p, p + tot) - p;
  int sz = 0;
  for (int i = 0; i < tot; ++ i) {
    while (sz > 1 && (p[st[sz - 1]] - p[st[sz - 2]]) * (p[i] - p[st[sz - 1]]) < 0) -- sz;
    st[sz ++] = i;
  }
  for (int i = 0; i < sz; ++ i) rec.push_back(st[i]); sz = 0;
  for (int i = tot - 1; ~ i; -- i) {
    while (sz > 1 && (p[st[sz - 1]] - p[st[sz - 2]]) * (p[i] - p[st[sz - 1]]) < 0) -- sz;
    st[sz ++] = i;
  }
  for (int i = 1; i < sz - 1; ++ i) rec.push_back(st[i]);
  for (int i = rec.size() - 1; ~ i; -- i) if (p[rec[i]].type == 1) return printf("Max\n"), 0;
  return printf("Min\n"), 0;
}

Matching Names

题目描述

Teachers of one programming summer school decided to make a surprise for the students by giving them names in the style of the “Hobbit” movie. Each student must get a pseudonym maximally similar to his own name. The pseudonym must be a name of some character of the popular saga and now the teachers are busy matching pseudonyms to student names.
There are $n$ students in a summer school. Teachers chose exactly $n$ pseudonyms for them. Each student must get exactly one pseudonym corresponding to him. Let us determine the relevance of a pseudonym $b$ to a student with name $a$ as the length of the largest common prefix $a$ and $b$. We will represent such value as $LCP(a, b)$. Then we can determine the quality of matching of the pseudonyms to students as a sum of relevances of all pseudonyms to the corresponding students.
Find the matching between students and pseudonyms with the maximum quality.

题意概述

给定$2n$个字符串$s_i$,前$n$个为一组,后$n$个为一组。要求将它们建立匹配(不能匹配同一组的字符串),使得两两匹配的字符串的$LCP$的长度之和最大。求最大值及方案。
数据范围:$1 \le n \le 10^5, \; 1 \le \sum |s_i| \le 8 \times 10^5$。

算法分析

首先对所有串建立Trie树,并记录下每个节点是哪些串的前缀。
考虑贪心——在Trie树上DFS(后序遍历),若在一个节点上有未匹配的不同组的串,则直接将它们匹配。这样做的正确性是显而易见的:若它们在之后匹配或与其他串匹配,则一定没有直接匹配更优。
具体实现时可以用vector启发式合并来维护每个节点上的值。

代码

#include <cstdio>
#include <vector>
#include <cstring>
#include <algorithm>

int min(int a, int b) {
  return a < b ? a : b;
}

static const int N = 800005;
int numn = 1, tot;
char s[N];
struct Trie {
  int child[26];
  std :: vector <int> rec[2];
} node[N];
std :: vector <std :: pair <int, int> > ans;

void insert(char s[], int id, int val) {
  int len = strlen(s), root = 0;
  for (int i = 0; i < len; ++ i) {
    if (! node[root].child[s[i] - 'a']) node[root].child[s[i] - 'a'] = numn ++;
    root = node[root].child[s[i] - 'a'];
  }
  node[root].rec[val].push_back(id);
}

void merge(std :: vector <int> &a, std :: vector <int> &b) {
  if (a.size() < b.size()) swap(a, b);
  a.insert(a.end(), b.begin(), b.end());
}

void dfs(int t, int dep) {
  for (int i = 0; i < 26; ++ i)
    if (node[t].child[i]) {
      dfs(node[t].child[i], dep + 1);
      merge(node[t].rec[0], node[node[t].child[i]].rec[0]);
      merge(node[t].rec[1], node[node[t].child[i]].rec[1]);
    }
  while (! node[t].rec[0].empty() && ! node[t].rec[1].empty()) {
    tot += dep;
    ans.push_back(std :: make_pair(node[t].rec[0][node[t].rec[0].size() - 1], node[t].rec[1][node[t].rec[1].size() - 1]));
    node[t].rec[0].pop_back(), node[t].rec[1].pop_back();
  }
}

int main() {
  int n; scanf("%d", &n);
  for (int i = 0; i < n; ++ i) scanf(" %s", s), insert(s, i + 1, 0);
  for (int i = 0; i < n; ++ i) scanf(" %s", s), insert(s, i + 1, 1);
  dfs(0, 0);
  printf("%d\n", tot);
  for (int i = 0; i < n; ++ i) printf("%d %d\n", ans[i].first, ans[i].second);
  return 0;
}

Similarity of Necklaces 2

题意概述

给定四个数组$M, P, L, R$,要求构造一个数组$T$,使得
$$\sum_{i=1}^N M_iT_i=0 \; (L_i \le T_i \le R_i)$$
求$\sum_{i=1}^N P_iT_i$的最大值。
数据范围:$1 \le N \le 200, \; 1 \le M_i \le 20, \; 0 \le P_i \le 10^5, \; -25 \le L_i \lt R_i \le 25$。

算法分析

令$D_i=T_i-L_i$,那么限制条件变为$\sum_{i=1}^N M_iD_i=-\sum_{i=1}^N M_iL_i \; (0 \le D_i \le R_i-L_i)$,要求的值变为$\sum_{i=1}^N P_iD_i+\sum_{i=1}^N P_iL_i$。这就相当于有$N$种物品,第$i$种物品有$(R_i-L_i)$个,每一个的体积为$M_i$,价值为$P_i$,背包的体积为$-\sum_{i=1}^N M_iL_i$,求物品恰好装满背包时的最大价值。这是分组背包问题,用单调队列优化DP即可。

代码

#include <cstdio>
#include <cstring>

int min(int a, int b) {
  return a < b ? a : b;
}

static const int N = 205;
int p[N], m[N], u[N], l[N];
int f[N << 10], q1[N << 10], q2[N << 10];

int main() {
  int n;
  while (~ scanf("%d", &n)) {
    int v = 0, w = 0;
    for (int i = 0; i < n; ++ i) {
      scanf("%d%d%d%d", &p[i], &m[i], &l[i], &u[i]);
      if (l[i]) u[i] -= l[i], v -= l[i] * m[i], w -= l[i] * p[i];
    }
    memset(f, -0x3f, sizeof f), f[0] = 0;
    for (int i = 0; i < n; ++ i) {
      u[i] = min(u[i], v / m[i]);
      for (int d = 0; d < m[i]; ++ d) {
        int l = 1, r = 0;
        for (int j = 0; j <= (v - d) / m[i]; ++ j) {
          int val = f[j * m[i] + d] - j * p[i];
          while (l <= r && q1[r] < val) -- r;
          q2[++ r] = j, q1[r] = val;
          if (j - q2[l] > u[i]) ++ l;
          f[j * m[i] + d] = q1[l] + j * p[i];
        }
      }
    }
    printf("%d\n", f[v] - w);
  }
  return 0;
}