洛谷P5548 [BJ United Round #3] 押韵
题面
题解
容易发现答案就是 \(\displaystyle [x^{nd}] \left(\sum_{i \geq 0} [d|i] \frac {x^i}{i!}\right)^k\)，对其进行单位根反演就是 \(\displaystyle [x^{nd}] \left(\frac 1d \sum_{j=0}^{d-1} \exp \omega_d^i x\right)^k\)。
\(d = 1, 2, 3\)
同复读机。
\(d = 4\)
通过观察可以发现其实是对于所有 \(\exp (a + \omega b)\) 求出它的系数，相当于你现在在 \((0, 0)\)，可以上下左右走，走 \(k\) 步恰好到 \((a, b)\) 的方案数。将坐标系旋转一下可以发现两维就独立了，方案数就是 \(\displaystyle \binom k {\frac {k + |a + b|} 2} \binom k {\frac {k + |a - b|} 2}\)。
当然也可以像 \(d = 6\) 一样推。
\(d = 6\)
此时可以发现任意 \(\omega^k\) 都可以用 \(1\) 和 \(\omega\) 表示出来。
那么问题也就变成了从 \((0, 0)\) 开始走 \(k\) 步，可以走 \((1, 0), (0, 1), (-1, 1), (-1, 0), (0, -1), (1, -1)\)，走到 \((a, b)\) 的方案数。
考虑到负数的下标不好处理，于是整体加一个向量 \((1, 1)\)。
令 \(F(x, y) = x + x^2 + x^2y + y + y^2 + y^2x\)
这样现在的问题就是求出 \(G(x, y) = F(x, y)^k\) 的系数。
对 \(x\) 求导，有 \(\displaystyle \frac {\partial G} {\partial x} = k F^{k - 1} \frac {\partial F} {\partial x}\)，令两边系数对应相等，有：
将 \(g_{n + 1, m - 1}\) 移到一边然后令 \(n \gets n - 1, m \gets m + 1\) 可得：
递推的边界是 \(g_{i, k - i} = \binom ki\)，\(g_{i, 0} = g_{0, i} = \binom k{i - k}\)，这样就做完了。
#include <cstdio>
#include <cstring>
#include <algorithm>
#include <vector>
#define file(x) freopen(#x".in", "r", stdin), freopen(#x".out", "w", stdout)
inline int read()
{
int data = 0, w = 1; char ch = getchar();
while (ch != '-' && (ch < '0' || ch > '9')) ch = getchar();
if (ch == '-') w = -1, ch = getchar();
while (ch >= '0' && ch <= '9') data = data * 10 + (ch ^ 48), ch = getchar();
return data * w;
}
const int N(2010), Mod(1049874433);
inline int upd(const int &x) { return x + (x >> 31 & Mod); }
int fastpow(int x, int y)
{
int ans = 1;
for (; y; y >>= 1, x = 1ll * x * x % Mod)
if (y & 1) ans = 1ll * ans * x % Mod;
return ans;
}
int n, K, D, w, pw[10], fac[N], inv[N], dK;
int C(int n, int m) { return 1ll * fac[n] * inv[m] % Mod * inv[n - m] % Mod; }
namespace Subtask1
{
int a[10], ans, res;
void calc()
{
for (int i = res = 0; i < D; i++) res = (res + 1ll * a[i] * pw[i]) % Mod;
res = fastpow(res, n);
for (int i = 0; i < D; i++) res = 1ll * res * inv[a[i]] % Mod;
ans = (ans + 1ll * res * fac[K]) % Mod;
}
void dfs(int x, int k)
{
if (x == D - 1) return a[x] = k, calc();
for (int i = 0; i <= k; i++) a[x] = i, dfs(x + 1, k - i), a[x] = 0;
}
void main()
{
dfs(0, K);
printf("%lld\n", 1ll * ans * dK % Mod);
}
}
namespace Subtask2
{
void main()
{
int ans = 0;
for (int i = 1; i <= K; i++)
for (int j = 0; i + j <= K; j++)
{
if ((i + j + K) & 1) continue;
int cnt = 1ll * C(K, (K + i + j) >> 1) * C(K, (K + std::abs(i - j)) >> 1) % Mod;
ans = (ans + 1ll * fastpow((i + 1ll * j * w) % Mod, n) * cnt) % Mod;
}
printf("%lld\n", 4ll * ans * dK % Mod);
}
}
namespace Subtask3
{
int g[N << 1][N << 1], m, ans, rev[N << 1];
void main()
{
m = K << 1, rev[1] = 1;
for (int i = 2; i <= m; i++) rev[i] = 1ll * (Mod - Mod / i) * rev[Mod % i] % Mod;
for (int i = 0; i <= K; i++) g[i][K - i] = C(K, i);
for (int i = K + 1; i <= m; i++) g[i][0] = g[0][i] = C(K, i - K);
for (int i = 1; i <= m; i++)
for (int j = std::max(K - i + 1, 1); j <= m; j++)
{
g[i][j] = 1ll * (K - i + 1) * (g[i - 1][j + 1] + g[i - 1][j - 1]) % Mod;
if (i > 1) g[i][j] = (g[i][j] + 1ll * (((K + 1) << 1) - i) * (g[i - 2][j + 1] + g[i - 2][j])) % Mod;
g[i][j] = (g[i][j] - 1ll * i * g[i][j - 1]) % Mod;
g[i][j] = 1ll * g[i][j] * rev[i] % Mod;
}
for (int i = 0; i <= m; i++) for (int j = 0; j <= m; j++) if (g[i][j])
ans = (ans + 1ll * g[i][j] * fastpow((i - K + 1ll * (j - K) * w) % Mod, n)) % Mod;
printf("%lld\n", 1ll * upd(ans) * dK % Mod);
}
}
int main()
{
#ifndef ONLINE_JUDGE
file(cpp);
#endif
n = read(), K = read(), D = read(), n = 1ll * n * D % (Mod - 1);
w = fastpow(7, (Mod - 1) / D), fac[0] = pw[0] = 1;
for (int i = 1; i < D; i++) pw[i] = 1ll * pw[i - 1] * w % Mod;
for (int i = 1; i <= K; i++) fac[i] = 1ll * fac[i - 1] * i % Mod;
inv[K] = fastpow(fac[K], Mod - 2), dK = fastpow(fastpow(D, K), Mod - 2);
for (int i = K; i; i--) inv[i - 1] = 1ll * inv[i] * i % Mod;
if (D <= 3) Subtask1::main();
if (D == 4) Subtask2::main();
if (D == 6) Subtask3::main();
return 0;
}
