51nod 1236 序列求和 V3

Fib(n)表示斐波那契数列的第n项，Fib(n) = Fib(n-1) + Fib(n-2)。Fib(0) = 0, Fib(1) = 1。

(1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, ...)

F(n, k) = Fib(n)^k（Fib(n)的k次幂）。

S(n, k) = F(1, k) + F(2, k) + ...... F(n, k)。

例如：S(4, 2) = 1^2 + 1^2 + 2^2 + 3^2 = 15。

给出n和k，求S(n, k) Mod 1000000009的结果。

输入

第1行：一个数T，表示后面用作输入测试的数的数量。（1 <= T <= 10)
第2 - T + 1行：每行2个数，N, K中间用空格分割。(1 <= N <= 10^18, 1 <= K <= 100000)

输出

输出共T行，对应S(n, k) Mod 1000000009的结果。

$$
f(i)=\frac{1}{\sqrt 5}\left[\left(\frac{1+\sqrt 5}{2}\right)^i-\left(\frac{1-\sqrt 5}{2}\right)^i\right]
$$

设 $a=\frac{1+\sqrt 5}{2},b=\frac{1-\sqrt 5}{2}$

$$
\begin{aligned}
&\sum_{i=1}^{n}\left\{\frac{1}{\sqrt 5}\left[\left(\frac{1+\sqrt 5}{2}\right)^i-\left(\frac{1-\sqrt 5}{2}\right)^i\right]\right\}^k \\
=&\left(\frac{1}{\sqrt 5}\right)^k\sum_{i=1}^{n} \sum_{j=0}^{k} {k \choose j} a^{ij}b^{i(k-j)}(-1)^{k-j} \\
=&\left(\frac{1}{\sqrt 5}\right)^k \sum_{j=0}^{k} {k \choose j} (-1)^{k-j} \sum_{i=1}^{n}\left(a^{j}b^{k-j}\right)^i \\
\end{aligned}
$$

设 $q=a^jb^{k-j}$，则可以化为：

$$
\left(\frac{1}{\sqrt 5}\right)^k \sum_{j=0}^{k} {k \choose j} (-1)^{k-j} \frac{q^{n+1}-q}{q-1}
$$

同时有：$\sqrt 5 \equiv 616991993 \pmod {10^9+9}$

然后似乎就做完了啊

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

const int mod = 1e9 + 9;

const int s5 = 616991993;

ll pw(ll a, ll b) {
    ll res = 1;
    for( ; b ; b >>= 1, a = a * a % mod)
        if(b & 1)
            res = res * a % mod;
    return res;
}

const int N = 1e5 + 10;

ll fac[N], invfac[N];

ll C(int n, int m) {
    return n < m ? 0 : fac[n] * invfac[m] % mod * invfac[n - m] % mod;
}

const ll a = (1 + s5) * pw(2, mod - 2) % mod;
const ll b = (1 - s5) * pw(2, mod - 2) % mod;

void sol() {
    ll n, k; cin >> n >> k;
    ll ans = 0;
    for(int j = 0 ; j <= k ; ++ j) {
        ll q = pw(a, j) * pw(b, k - j) % mod;
        if(q == 1) {
            ans += C(k, j) * pw(-1, k - j) % mod * (n % mod) % mod;
        } else {
            ans += C(k, j) * pw(-1, k - j) % mod * ((pw(q, n + 1) - q) % mod * pw(q - 1, mod - 2) % mod) % mod;
        }
        ans %= mod;
    }
    ans = ans * pw(pw(s5, mod - 2), k) % mod;
    cout << (ans % mod + mod) % mod << endl;
}

int main() {
    ios :: sync_with_stdio(0);
    const int T = 1e5;
    fac[0] = invfac[0] = 1;
    for(int i = 1 ; i <= T ; ++ i) fac[i] = fac[i - 1] * i % mod;
    invfac[T] = pw(fac[T], mod - 2);
    for(int i = T - 1 ; i >= 1 ; -- i) invfac[i] = invfac[i + 1] * (i + 1) % mod;
    int fafa; cin >> fafa;
    while(fafa --) sol();
}