题目大意
定义斐波那契数列：
[F_0=0,F_1=1\
F_n=F_{n-1}+F_{n-2} (n>1)
]
给定整数(N,C,K)。请计算下式的值：
[(F_{0})^K+(F_{C})^K+(F_{2C})^K+dots +(F_{NC})^{K}
]
因为答案太大，请输出其对(1000000009) ((10^9+9)) 取模的值。
数据范围：(1leq N,Cleq 10^{18})，(1leq Kleq 10^5)。
本题题解
首先需要了解斐波那契数列的通项公式：
[F_n=frac{1}{sqrt{5}}left(left(frac{1+sqrt{5}}{2}
ight)^{n}-left(frac{1-sqrt{5}}{2}
ight)^{n}
ight)
]
通过暴力，可以找到(mod10^9+9)意义下(sqrt{5})的值是(383008016)。
const int MOD=1e9+9;
for(int i=0;i<MOD;++i){
if((long long)i*i%MOD == 5){
cout<<i<<endl;//383008016
break;
}
}
考虑求答案。先把(frac{1}{sqrt{5}})提取出来，即：
[egin{align}
ext{answer}&=sum_{i=0}^{N}(F_{iC})^{K}\
&=sum_{i=0}^{N}left(frac{1}{sqrt{5}}left(left(frac{1+sqrt{5}}{2}
ight)^{iC}-left(frac{1-sqrt{5}}{2}
ight)^{iC}
ight)
ight)^K\
&=left(frac{1}{sqrt{5}}
ight)^Kcdot sum_{i=0}^{N}left(left(frac{1+sqrt{5}}{2}
ight)^{iC}-left(frac{1-sqrt{5}}{2}
ight)^{iC}
ight)^{K}
end{align}
]
设(A=left(frac{1+sqrt{5}}{2} ight)^C), (B=left(frac{1-sqrt{5}}{2} ight)^C)。则：
[ ext{answer}=left(frac{1}{sqrt{5}}
ight)^Ksum_{i=0}^{N}(A^i-B^i)^K
]
问题转化为求：
[sum_{i=0}^{N}(A^i-B^i)^K
]
把后面的东西用二项式定理展开，得到：
[=sum_{i=0}^{N}sum_{j=0}^{K}{Kchoose j}A^{ij}(-B^i)^{K-j}
]
交换和式，得：
[=sum_{j=0}^{K}{Kchoose j}(-1)^{K-j}sum_{i=0}^{N}(A^j)^i(B^{K-j})^i
]
后面的(sum_{i=0}^{N}(A^j)^i(B^{K-j})^i)，是一个等比数列求和，可以(O(log N))求出来。（即(sum_{i=0}^{n}x^i=frac{x^{n+1}-1}{x-1})，(x eq 1)）。复杂度里的(log)来自求(x^{n+1})时使用的快速幂。
于是总时间复杂度就是(O(Klog N))的了。
参考代码：
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int MOD=1e9+9;
const int MAXK=1e5;
int pow_mod(int x,int i){
int y=1;
while(i){
if(i&1)y=(ll)y*x%MOD;
x=(ll)x*x%MOD;
i>>=1;
}
return y;
}
ll N,C;
int K,fac[MAXK+5],ifac[MAXK+5];
int comb(int n,int k){
if(n<k)return 0;
return (ll)fac[n]*ifac[k]%MOD*ifac[n-k]%MOD;
}
int main(){
ios::sync_with_stdio(0);
fac[0]=1;
for(int i=1;i<=MAXK;++i)fac[i]=(ll)fac[i-1]*i%MOD;
ifac[MAXK]=pow_mod(fac[MAXK],MOD-2);
for(int i=MAXK-1;i>=0;--i)ifac[i]=(ll)(i+1)*ifac[i+1]%MOD;
int T;cin>>T;while(T--){
cin>>N>>C>>K;
int A=691504013,B=308495997;
A=pow_mod(A,C%(MOD-1));
B=pow_mod(B,C%(MOD-1));
int a=1,b=pow_mod(B,K);
int ib=pow_mod(B,MOD-2);
int ans=0;
for(int j=0;j<=K;++j){
int x=(ll)a*b%MOD;
if(x==1)
x=(N+1)%MOD;
else
x=(ll)(pow_mod(x,(N+1)%(MOD-1))-1+MOD)%MOD * pow_mod((x-1+MOD)%MOD,MOD-2) % MOD;
if((K-j)&1)
x=(x==0?x:MOD-x);
ans=((ll)ans+(ll)comb(K,j)*x)%MOD;
a=(ll)a*A%MOD;
b=(ll)b*ib%MOD;
}
// for(int i=0;i<=N;++i){
// ans=(ans+pow_mod((pow_mod(A,i)-pow_mod(B,i)+MOD)%MOD,K))%MOD;
// }
int mul=276601605;
mul=pow_mod(mul,K);
ans=(ll)ans*mul%MOD;
cout<<ans<<endl;
}
return 0;
}