問題描述

熱情好客的小猴子請森林中的朋友們吃飯，他的朋友被編號為 \(1\sim N\) 每個到來的朋友都會帶給他一些禮物：大香蕉。其中，第一個朋友會帶給他 1 個大香蕉，之後，每一個朋友到來以後，都會帶給他之前所有人帶來的禮物個數再加他的編號的 K 次方那麼多個。所以，假設 K=2，前幾位朋友帶來的禮物個數分別是：

\(1,5,15,37,83,\ldots\)

假設 K=3，前幾位朋友帶來的禮物個數分別是：

\(1,9,37,111,\ldots\)

現在，小猴子好奇自己到底能收到第 N 個朋友多少禮物，因此拜託於你了。

已知 N,K，請輸出第 N 個朋友送的禮物個數對 \(10^9+7\)

輸入格式

第一行，兩個整數 N,K。

輸出格式

一個整數，表示第 N 個朋友送的禮物個數對 \(10^9+7\) 取模的結果。

樣例輸入

1234567890000 3

樣例輸出

891659731

資料範圍

100% 的資料：\(N \le 10^{18}\)，\(K \le 10\)。

解析

設 \(f_i\)表示第 \(i\) 個朋友送的禮物數量，\(sum_i\) 表示前 \(i\) 個朋友送的禮物個數之和。不難得到 \(f_i=sum_{i-1}+i^k\)。由此，我們可以得到遞推式：

\[sum_i=2sum_{i-1}+i^k \]

由於 \(n\le 10^{18}\)，我們顯然要用矩陣快速冪解決這個問題。但 \(i^k\)

\[\left[ \begin{matrix} 2 & C_k^0 & C_k^1 & C_k^2 & ... &C_k^k\\ 0 & C_k^0 & C_k^1 & C_k^2 & ... &C_k^k \\ 0 & 0 & C_{k-1}^0 & c_{k-1}^1 & ... &C_{k-1}^{k-1} \\ .&.&.&.&.&.\\ 0&0&0&0&...&C_0^0 \end{matrix} \right]\times \left[ \begin{matrix} s_i\\i^k\\i^{k-1}\\...\\i^0 \end{matrix} \right]= \left[ \begin{matrix} s_{i+1}\\{(i+1)}^k\\{(i+1)}^{k-1}\\...\\{(i+1)}^0 \end{matrix} \right] \]

矩陣快速冪即可，利用字首和的性質得到答案。注意n=1或2時需特判。

程式碼

#include <iostream>
#include <cstdio>
#include <cstring>
#define int long long
#define N 20
using namespace std;
const int mod=1000000007;
struct Matrix{
	int n,m,a[N][N];
}ans,s;
Matrix operator * (Matrix a,Matrix b)
{
	Matrix c;
	c.n=a.n;c.m=b.m;
	memset(c.a,0,sizeof(c.a));
	for(int i=1;i<=c.n;i++){
		for(int j=1;j<=c.m;j++){
			for(int k=1;k<=a.m;k++) c.a[i][j]=(c.a[i][j]+a.a[i][k]*b.a[k][j]%mod)%mod;
		}
	}
	return c;
}
int n,k,i,j,c[N][N],sum1,sum2;
Matrix poww(Matrix a,int b)
{
	Matrix ans=a,base=a;
	b--;
	while(b){
		if(b&1) ans=ans*base;
		base=base*base;
		b>>=1;
	}
	return ans;
}
signed main()
{
	scanf("%lld%lld",&n,&k);
	if(n==1){
		puts("1");
		return 0;
	}
	if(n==2){
		int ans=1;
		for(i=1;i<=k;i++) ans=ans*2%mod;
		printf("%lld\n",ans+1);
		return 0;
	}
	c[0][0]=1;
	for(i=1;i<=k;i++){
		c[i][0]=1;
		for(j=1;j<=i;j++) c[i][j]=(c[i-1][j-1]+c[i-1][j])%mod;
	}
	s.a[1][1]=2;s.n=s.m=k+2;
	for(i=2;i<=k+2;i++) s.a[1][i]=c[k][i-2];
	for(i=2;i<=k+2;i++){
		for(j=i;j<=k+2;j++) s.a[i][j]=c[k-i+2][j-i];
	}
	ans.n=k+2;ans.m=1;
	for(i=1;i<=k+2;i++) ans.a[i][1]=1;
	ans=poww(s,n-2)*ans;
	sum1=ans.a[1][1];
	ans=s*ans;
	sum2=ans.a[1][1];
	printf("%lld\n",(sum2-sum1+mod)%mod);
	return 0;
}