【lgj】【矩阵快速幂】 k 次幂

问题:求 i=1nik(k50,n1×109)\sum\limits_{i=1}^{n}i^k(k\le50,n\leq1 \times 10^9)

柿子是推不出来滴... nn 那么大,考虑矩阵。

那就要把 iki^k 转化为 (i+1)k(i+1)^k。用二项式定理展开一下吧:

(i+1)k=Ck0i0+Ck1i1++Ckkik(i+1)^k = C_k^0 i^0 + C_k^1 i^{1}+ \dots +C_k^ki^k

那我们就把 i0iki^0\sim i^k 全部丢到答案矩阵,然后再用上式填转移的矩阵。再开一个记录 kk 次幂的和就好了。

/*
	- 别摆了
	- By yfz
*/
#include<bits/stdc++.h>
using namespace std;
#define int long long
const int N = 52, mod = 1e9+7;
int n,k;
struct Matrix{
	int a[N][N];
	Matrix() {memset(a,0,sizeof a);}
	void build() {
		for(int i=0;i<=n;i++) {// 需要注意:i的取值范围为实际矩阵的大小,即 [0,k] 
			a[i][i]=1;
		}
	}
	friend Matrix operator * (const Matrix &A, const Matrix &B) {
		Matrix p;
		for(int i=0;i<=k;i++) {
			for(int j=0;j<=k;j++) {
				for(int kk=0;kk<=k;kk++) {
					p.a[i][j]=(p.a[i][j]+A.a[i][kk]*B.a[kk][j]%mod)%mod;
				}
			}
		}
		return p;
	}
};
Matrix mi(Matrix x,long long k) {
	Matrix p; p.build();
	if(!k) return p;
	p=mi(x,k>>1); p=p*p;
	return k&1?p*x:p;
}
int jc[N],inv[N];
int mii(int x,int k) {
	int p; return k?p=mii(x,k>>1),p*p%mod*(k&1?x:1)%mod:1;
}
int C(int n,int m) {
	return jc[n] * inv[m] % mod * inv[n-m] % mod;
}
signed main() {
	cin>>n>>k; k++;
	jc[0] = inv[0] = 1;
	for(int i=1;i<=51;i++) {
		jc[i] = jc[i-1] * i % mod;
		inv[i] = mii(jc[i], mod-2);
	}
	Matrix x; for(int i=0;i<=k-1;i++) { x.a[0][i] = 1;}
	Matrix p;
	for(int i=0;i<=k-1;i++) {
		for(int j=i;j<=k-1;j++) {
//			cout<<j<<" "<<i<<" "<<C(j,i)<<endl;
//			cout<<(
			p.a[i][j] = C(j,i);//)<<" ";
		}
//		cout<<endl;
	}
	p.a[k][k] = 1, p.a[k-1][k] = 1;
	p = mi(p,n); x = x*p;
	cout<<x.a[0][k];
	return 0;
}
posted @ 2024-06-19 13:40  cjrqwq  阅读(8)  评论(0编辑  收藏  举报  来源