Discription
Everyone knows what the Fibonacci sequence is. This sequence can be defined by the recurrence relation:
We'll define a new number sequence Ai(k) by the formula:
In this problem, your task is to calculate the following sum: A1(k) + A2(k) + ... + An(k). The answer can be very large, so print it modulo 1000000007 (109 + 7).
Input
The first line contains two space-separated integers n, k (1 ≤ n ≤ 1017; 1 ≤ k ≤ 40).
Output
Print a single integer — the sum of the first n elements of the sequence Ai(k)modulo 1000000007 (109 + 7).
Examples
Input
1 1
Output
1
Input
4 1
Output
34
Input
5 2
Output
316
Input
7 4
Output
73825
之前做过一道 需要求 f[i] * i 生成函数形式的题,在那道题的题解里(http://www.cnblogs.com/JYYHH/p/8822572.html)已经证明过了
这个玩意的生成函数的分母是 (1-x-x^2)^2 。。
当然,更普遍的,我们可以证明 f[i] * i^k 的生成函数表示的分母是 (1-x-x^2)^(k+1) ,这个用二项式定理解一下多项式闭形式就ojbk了。
于是我们可以得到这个函数的递推式,于是就可以直接预处理出前若干项然后直接用矩阵做了。
你问我它还要求前缀和???? 这个是矩阵的常规操作啊23333
#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int maxn=93;
const int ha=1000000007;
int K,A[maxn],n,F[maxn];
inline int add(int x,int y){ x+=y; return x>=ha?x-ha:x;}
inline int ksm(int x,int y){ int an=1; for(;y;y>>=1,x=x*(ll)x%ha) if(y&1) an=an*(ll)x%ha; return an;}
struct node{
int a[maxn][maxn];
inline void clear(){ memset(a,0,sizeof(a));}
inline void BASE(){ clear(); for(int i=n+1;i;i--) a[i][i]=1;}
node operator *(const node &u)const{
node r; r.clear();
for(int k=n+1;k;k--)
for(int i=n+1;i;i--)
for(int j=n+1;j;j--) r.a[i][j]=add(r.a[i][j],a[i][k]*(ll)u.a[k][j]%ha);
return r;
}
}X,ANS;
ll N;
inline void init(){
n=2,A[0]=1,A[1]=A[2]=ha-1;
for(int i=1;i<=K;i++){
for(int j=n;j>=0;j--){
A[j+2]=add(A[j+2],ha-A[j]);
A[j+1]=add(A[j+1],ha-A[j]);
}
n+=2;
}
for(int i=1;i<=n;i++) A[i]=ha-A[i];
F[0]=F[1]=1;
for(int i=2;i<=n;i++) F[i]=add(F[i-1],F[i-2]);
for(int i=1;i<=n;i++) F[i]=F[i]*(ll)ksm(i,K)%ha;
}
inline void build(){
X.clear(),ANS.BASE();
for(int i=1;i<n;i++) X.a[i][i+1]=1;
for(int i=1;i<=n;i++) X.a[i][1]=X.a[i][n+1]=A[i];
X.a[n+1][n+1]=1;
}
inline int calc(){
int ans=0,cnt=0;
if(N<=n) for(int i=1;i<=N;i++) ans=add(ans,F[i]);
else{
N-=n,memset(A,0,sizeof(A));
for(int i=1;i<=n;i++){
A[i]=F[n-i+1];
A[n+1]=add(A[n+1],F[i]);
}
for(;N;N>>=1,X=X*X) if(N&1) ANS=ANS*X;
for(int i=n+1;i;i--) ans=add(ans,A[i]*(ll)ANS.a[i][n+1]%ha);
}
return ans;
}
inline void solve(){
cin>>N>>K;
init(),build();
printf("%d\n",calc());
}
int main(){
solve();
return 0;
}