# 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}

$$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+1)=\frac{F_k(n)}r + (n + 1)^k$

$\sum_{i=0}^{n+1} (-1)^i \binom {n + 1} i f(i) = 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;
}

posted @ 2021-02-22 14:42  xgzc  阅读(152)  评论(0编辑  收藏  举报