51nod 1822 序列求和 V5

题面

组合数没了，怎么说嘛

题解

当 \(r = 0\) 或 \(r = 1\) 时是 trivial 的，下面假设 \(r \neq 0, 1\)。

设 \(S_k(n) = \sum_{i=1}^n i^k r^i\)。

那么有

\[\begin{aligned} S_k(n) &= \sum_{i=1}^n i^k r^i\\ &= \frac 1{r - 1} \left( n^kr^{n + 1} + \sum_{i=1}^n r^i((i - 1)^k - i^k) \right)\\ &= \frac 1{r - 1} \left( n^kr^{n + 1} + \sum_{i=1}^n r^i \sum_{j=0}^{k - 1} (-1)^{k-j} i^j \binom kj \right)\\ &= \frac 1{r - 1} \left( n^kr^{n + 1} + \sum_{j=0}^{k - 1} (-1)^{k-j} \binom kj \sum_{i=1}^n r^i i^j \right)\\ &= \frac 1{r - 1} \left( n^kr^{n + 1} + \sum_{j=0}^{k - 1} (-1)^{k-j} \binom kj S_j(n) \right) \end{aligned} \]

不过好像只能 \(\mathcal O(k^2)\) 算。

但是如果将 \(S_k(n)\) 的封闭形式展开，也许就可以发现下面这个规律：

引理 1 存在一个 \(k\) 次多项式 \(F_k(x)\) 满足 \(S_k(n) = r^{n} F_k(n) - F_k(0)\)。

当 \(k = 0\) 时显然成立。

当 \(k > 0\) 时：

\[\begin{aligned} S_k(n) &= \frac 1{r - 1} \left( n^kr^{n + 1} + \sum_{j=0}^{k - 1} (-1)^{k-j} \binom kj S_j(n) \right)\\ &= \frac 1{r - 1} \left( n^kr^{n + 1} + \sum_{j=0}^{k - 1} (-1)^{k-j} \binom kj (r^n F_j(n) - F_j(0)) \right)\\ &= r^n \left(\frac{r n^k + \sum_{j=0}^{k - 1} (-1)^{k-j} \binom kj F_j(n)}{r-1} \right) - \left(\frac{ \sum_{j=0}^{k - 1} (-1)^{k-j} \binom kj F_j(0)}{r-1} \right)\\ \end{aligned} \]

也就是说，只要 \(F_k(n) = \frac{r n^k + \sum_{j=0}^{k - 1} (-1)^{k-j} \binom kj F_j(n)}{r-1}\) 就是合法的，这显然是一个 \(k\) 次多项式。

考虑如何计算 \(F_k(n)\)，首先由定义立即可以得出：

\[F_k(n+1)=\frac{F_k(n)}r + (n + 1)^k \]

我们只需要知道 \(F_k(0)\) 就可以插值出答案了。

引理 2 对任意一个 \(n\) 次多项式 \(f(x)\)，满足如下恒等式：

\[\sum_{i=0}^{n+1} (-1)^i \binom {n + 1} i f(i) = 0 \]

定义平移算子 \(\mathrm E f(x) = f(x + 1)\)，恒等算子 \(1 f(x) = f(x)\)，差分算子 \(\Delta f(x) = f(x + 1) - f(x)\)。

那么有 \(\Delta^n = (\mathrm E - 1)^n = \sum_{i=0}^n (-1)^i \binom ni \mathrm E^i 1^{n-i}\)。

所以 \(\Delta^{n+1} f(0) = \sum_{i=0}^{n+1} (-1)^i \binom {n+1}i f(i)\)，由于 \(f\) 是 \(n\) 次多项式，所以 \(\Delta^{n+1} f(0) = 0\)，得证。

由于所有的 \(F_k(n)\) 都可以用 \(F_k(0)\) 表示，所以就可以得出一个和 \(F_k(0)\) 有关的方程，解出来即可。

代码

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <vector>
#define file(x) freopen(#x".in", "r", stdin), freopen(#x".out", "w", stdout)

template <typename T> inline T Read()
{
	T 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;
}

#define read Read<int>
#define readll Read<ll>

using ll = long long;
const int N(2e5 + 10), Mod(985661441);
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 T, K, R, prime[N], not_prime[N], cnt, Pow[N], f[N], fac[N], inv[N];
long long n;

inline int C(int n, int m) { return 1ll * fac[n] * inv[m] % Mod * inv[n - m] % Mod; }
int Lagrange(long long x, int m)
{
	if (x <= m) return f[x];
	static int g[N], h[N], y, ans; g[0] = h[m + 1] = 1, y = x % Mod, ans = 0;
	for (int i = 1; i <= m; i++) g[i] = 1ll * g[i - 1] * upd(y - i) % Mod;
	for (int i = m; i; i--) h[i] = 1ll * h[i + 1] * upd(y - i) % Mod;
	for (int i = 1, op = (m & 1) ? 1 : Mod - 1; i <= m; op = Mod - op, i++)
		ans = (ans + 1ll * op * f[i] % Mod * inv[i - 1] % Mod * inv[m - i] % Mod * g[i - 1] % Mod * h[i + 1]) % Mod;
	return ans;
}

int main()
{
#ifndef ONLINE_JUDGE
	file(cpp);
#endif
	K = 2e5 + 9, fac[0] = 1;
	for (int i = 1; i <= K; i++) fac[i] = 1ll * fac[i - 1] * i % Mod;
	inv[K] = fastpow(fac[K], Mod - 2);
	for (int i = K; i; i--) inv[i - 1] = 1ll * inv[i] * i % Mod;
	for (T = read(); T--; )
	{
		K = read(), n = readll(), R = readll() % Mod;
		Pow[1] = not_prime[1] = 1, cnt = f[0] = 0;
		for (int i = 2; i <= K + 2; i++)
		{
			if (!not_prime[i]) prime[++cnt] = i, Pow[i] = fastpow(i, K);
			for (int j = 1; j <= cnt && i * prime[j] <= K + 2; j++)
			{
				not_prime[i * prime[j]] = 1, Pow[i * prime[j]] = 1ll * Pow[i] * Pow[prime[j]] % Mod;
				if (i % prime[j] == 0) break;
			}
		}
		if (R == 0) { puts("0"); continue; }
		if (R == 1)
		{
			for (int i = 1; i <= K + 2; i++) f[i] = upd(f[i - 1] + Pow[i] - Mod);
			printf("%d\n", Lagrange(n, K + 2)); continue;
		}
		int k = 1, b = 0, sk = 1, sb = 0, invR = fastpow(R, Mod - 2);
		for (int i = 1; i <= K + 1; i++)
		{
			k = 1ll * k * invR % Mod, b = (1ll * b * invR + Pow[i]) % Mod;
			if (i & 1) sk = (sk - 1ll * C(K + 1, i) * k) % Mod, sb = (sb - 1ll * C(K + 1, i) * b) % Mod;
			else sk = (sk + 1ll * C(K + 1, i) * k) % Mod, sb = (sb + 1ll * C(K + 1, i) * b) % Mod;
		}
		f[0] = 1ll * sb * fastpow(sk, Mod - 2) % Mod, f[0] = upd(-f[0]);
		for (int i = 1; i <= K + 1; i++) f[i] = (1ll * f[i - 1] * invR + Pow[i]) % Mod;
		printf("%d\n", upd(1ll * Lagrange(n, K + 1) * fastpow(R, n % (Mod - 1)) % Mod - f[0]));
	}
	return 0;
}