题解 P5320 【[BJOI2019]勘破神机】
_rqy
2019-04-21 15:34:36
先看 $m=2$。众所周知的是 $2\times n$ 填骨牌的方案是斐波那契数列 $f_{n+1}$，其中 $f_0=0,f_1=1, f_n=f_{n-1}+f_{n-2}$，而要求的东西就是 $\sum_{n=l}^r{f_n\choose k}$ （为了方便我们假设方案是 $f_n$ 而不是 $f_{n+1}$，反正只需要把 $l,r$ 都加上 $1$ 就可以了）
然而斐波那契数列的通项公式为
$$f_n=\frac1{\sqrt5}\left(\frac{1+\sqrt5}{2}\right)^n-\frac1{\sqrt5}\left(\frac{1-\sqrt5}{2}\right)^n$$
设 $A=\frac1{\sqrt5},B=-\frac1{\sqrt5},x=\frac{1+\sqrt5}{2},y=\frac{1-\sqrt5}{2}$，那么有 $f_n=Ax^n+Bx^n$
再加上下降幂与阶乘幂之间的转换
$$x^{\underline k}=\sum_{i=0}^k(-1)^{k-i}s(k,i)x^i$$
（其中 $s(k,i)$ 为第一类斯特林数）
于是
$$\begin{aligned}&\sum_{n=l}^r{f_n\choose k}\\=&\frac1{k!}\sum_{n=l}^r\sum_{i=0}^k(-1)^{k-i}s(k,i)f_n^i\\=&\frac1{k!}\sum_{n=l}^r\sum_{i=0}^k(-1)^{k-i}s(k,i)(Ax^n+By^n)^i\\=&\frac1{k!}\sum_{n=l}^r\sum_{i=0}^k(-1)^{k-i}s(k,i)\sum_{j=0}^i{i\choose j}A^jB^{i-j}(x^jy^{i-j})^n\\=&\frac1{k!}\sum_{i=0}^k(-1)^{k-i}s(k,i)\sum_{j=0}^i{i\choose j}A^jB^{i-j}\sum_{n=l}^r(x^jy^{i-j})^n\\\end{aligned}$$
只需要枚举 $i,j$，然后后面那个 $\sum$ 就是一个等比数列。
不过还有一个问题就是， $\sqrt5$ 这个东西在 $\bmod998244353$ 意义下不存在。
我们对剩余系进行扩域，把每个数表示成 $a+b\sqrt5$ 的形式，就可以做了。
对于 $m=3$，可以发现 $n$ 是奇数的时候方案数是 $0$，而偶数的时候（令 $g_n$ 表示 $2n$ 的答案）可以搞出递推式 $g_n=4g_{n-1}-g_{n-2}$，从而解出通项公式 $g_n=\frac{3+\sqrt3}6(2+\sqrt3)^n+\frac{3-\sqrt3}6(2-\sqrt3)^n$，仍然用上述方法计算即可。
（不过有一个本题中没什么用处的就是，实际上模素数意义下当 $\sqrt{x}$ 不存在的时候 $\sqrt{-x}$ 一定是存在的，因此可以将 $a+b\sqrt x$ 写成 $a+b'i$，其中 $b'=b\sqrt{-x}$ 而 $i^2=1$ （模意义复数！））
附代码：
```cpp
#include <algorithm>
#include <cstdio>
#include <cstring>
typedef long long LL;
const int mod = 998244353;
LL pow_mod(LL a, LL b) {
LL ans = 1;
for (a %= mod; b; b >>= 1, a = a * a % mod)
if (b & 1) ans = ans * a % mod;
return ans;
}
template <int a>
struct Sq { // x + y*sqrt(a) (mod 998244353)
typedef Sq<a> sq;
LL x, y;
Sq(LL x = 0, LL y = 0) : x(x % mod), y(y % mod) {}
inline sq conj() const { return sq(x, -y); }
// assume a is small (a * mod * mod < 2^63)
friend inline sq operator+(const sq &l, const sq &r) { return sq(l.x + r.x, l.y + r.y); }
friend inline sq operator-(const sq &l, const sq &r) { return sq(l.x - r.x, l.y - r.y); }
friend inline sq operator*(const sq &l, const sq &r) { return sq(l.x * r.x + a * l.y * r.y, l.x * r.y + l.y * r.x); }
friend inline sq operator/(const sq &l, const sq &r) { return sq(l.x * r.x - a * l.y * r.y, l.y * r.x - l.x * r.y) * pow_mod(r.x * r.x - a * r.y * r.y, mod - 2); }
inline sq& operator+=(const sq &t) { x = (x + t.x) % mod; y = (y + t.y) % mod; return *this; }
inline sq& operator-=(const sq &t) { x = (x - t.x) % mod; y = (y - t.y) % mod; return *this; }
inline sq& operator*=(const sq &t) { return *this = *this * t; }
sq operator^(LL t) const {
sq ans = 1, x = *this;
for (; t; t >>= 1, x *= x) if (t & 1) ans *= x;
return ans;
}
inline sq sum_pow(LL t) { return (1 + this->pow(t + 1)) / (1 + *this); } // (x^(t+1)-1) / (x-1)
inline void put() { printf("%lld+%llds%d", (x + mod) % mod, (y + mod) % mod, a); }
};
const int K = 550;
LL C[K][K], ss[K][K], fac[K];
template <int a>
Sq<a> Solve(Sq<a> A, Sq<a> p, Sq<a> B, Sq<a> q, int k, LL l, LL r) {
++r;
Sq<a> ans = 0, qr = q^r, ql = q^l, pr = p^r, pl = p^l;
Sq<a> ppl = 1, ppr = 1, pp = 1, pA = 1;
for (int u = 0; u <= k; ++u) {
Sq<a> pql = 1, pqr = 1, pq = 1, pB = 1;
for (int v = 0; u + v <= k; ++v) {
Sq<a> quq = pp * pq - 1, qvq = ppr * pqr - ppl * pql;
Sq<a> qaq = !quq.x && !quq.y ? (r - l) % mod : qvq / quq;
Sq<a> qwq = ss[k][u + v] * C[u + v][u] % mod * pA * pB * qaq;
ans += qwq;
pB *= B; pql *= ql; pqr *= qr; pq *= q;
}
pA *= A; ppl *= pl; ppr *= pr; pp *= p;
}
return ans * pow_mod(fac[k], mod - 2);
}
void Init() {
for (int i = 0; i < K; ++i)
for (int j = C[i][0] = 1; j <= i; ++j)
C[i][j] = (C[i - 1][j] + C[i - 1][j - 1]) % mod;
fac[0] = 1;
for (int i = 1; i < K; ++i) fac[i] = i * fac[i - 1] % mod;
ss[0][0] = 1;
for (int i = 0; i + 1 < K; ++i) {
ss[i + 1][0] = - i * ss[i][0];
for (int j = 1; j <= i + 1; ++j)
ss[i + 1][j] = (ss[i][j - 1] - i * ss[i][j]) % mod;
}
}
int main() {
Init();
int T, m;
scanf("%d%d", &T, &m);
while (T--) {
LL l, r; int k;
scanf("%lld%lld%d", &l, &r, &k);
if (m == 2) {
Sq<5> p = Sq<5>(1, 1) / 2, t = Sq<5>(0, 1) / 5;
printf("%lld\n", (Solve(t, p, t.conj(), p.conj(), k, l + 1, r + 1).x * pow_mod(r - l + 1, mod - 2) % mod + mod) % mod);
} else {
Sq<3> p = Sq<3>(2, 1), t = Sq<3>(3, 1) / 6;
printf("%lld\n", (Solve(t, p, t.conj(), p.conj(), k, (l + 1) / 2, r / 2).x * pow_mod(r - l + 1, mod - 2) % mod + mod) % mod);
}
}
}
```