http://codeforces.com/problemset/problem/392/C
Description
Everyone knows what the Fibonacci sequence is. This sequence can be defined by the recurrence relation:
F1 = 1, F2 = 2, Fi = Fi - 1 + Fi - 2 (i > 2).
We'll define a new number sequence Ai(k) by the formula:
Ai(k) = Fi × ik (i ≥ 1).
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).
Sample Input
Input
1 1
Output
1
Input
4 1
Output
34
Input
5 2
Output
316
Input
7 4
Output
73825 设: 为了方便讨论,先不考虑mod的问题:读题发现,n很大,k很小,且和斐波那契相关,所以考虑用矩阵把n不断变小。现在找关系:由牛顿二项式定理:
有:
为了方便,进一步换元:
完成了从n+1变成n,关系矩阵:
当最右边的矩阵n=1时,那一列全是1.
哈哈,接下来就是矩阵快速幂。(第一次用Latex写这么大的矩阵,我的手~~~)
#include <iostream>
#include <cstdio>
#include <cstring>
using namespace std;
typedef long long LL;
const int mod=1e9+7;
int len;
LL n,k;
LL c[45][45];
struct matrix{
LL m[100][100];
}mc,I;
matrix multi(matrix a,matrix b){
matrix ans;
for(int i=0;i<len;i++){
for(int j=0;j<len;j++){
ans.m[i][j]=0;
for(int k=0;k<len;k++){
ans.m[i][j]=(ans.m[i][j]+a.m[i][k]*b.m[k][j]%mod)%mod;
}
}
}
return ans;
}
matrix quick(matrix a, LL n){
matrix ans=I,t=a;
while(n>0){
if(n&1) ans=multi(ans,t);
t=multi(t,t);
n>>=1;
}
return ans;
}
int main()
{
c[0][0]=1;
c[1][0]=c[1][1]=1;
for(int i=2;i<=40;i++){
for(int j=0;j<=i;j++){
if(j==0 || j==i) c[i][j]=1;
else c[i][j]=c[i-1][j-1]+c[i-1][j];
}
}
for(int i=0;i<=40;i++){
I.m[i][i]=1;
}
while(cin>>n>>k){
memset(mc.m,0,sizeof(mc.m));
len=2*k+3;
mc.m[0][0]=1;
for(int j=0;j<=k;j++){
mc.m[0][j+1]=c[k][j];
mc.m[0][j+k+2]=c[k][j];
}
mc.m[1][1]=mc.m[1][k+2]=1;
for(int i=1;i<=k;i++){
for(int j=0;j<=i;j++){
mc.m[i+1][j+1]=c[i][j];
mc.m[i+1][j+k+2]=c[i][j];
}
}
mc.m[k+2][1]=1;
for(int i=1;i<=k;i++){
for(int j=0;j<=i;j++){
mc.m[i+k+2][j+1]=c[i][j];
}
}
mc=quick(mc,n-1);
LL ans=0;
for(int i=0;i<len;i++){
ans=(ans+mc.m[0][i]%mod)%mod;
}
printf("%I64d\n",ans);
}
return 0;
}