SP19997 MOON2

题目描述:求

[sum_{i=1}^ni^kr^i ]

对某个质数取模。(T)组数据。

数据范围:(n,rle 10^{18},sum kle 2.56 imes 10^6)


数论神题。。。

我们设(S_k(n)=sum_{i=1}^{n-1}i^kr^i),通过推通项公式,你会发现这个结论:

[S_k(n)=r^nF_k(n)-F_k(0) ]

其中(F_k(n))是一个关于(n)且次数不高于(k)的多项式。

如何证明呢?首先你可以通过错位相减法得到一个递推公式,然后就可以数学归纳法了。

然后我们就是要求出(F_k(0),F_k(1),F_k(2),ldots,F_k(k+1)),这样就可以通过拉格朗日插值求出(F_k(n))

我们注意到(S_k(n)-S_k(n-1)=r^nF_k(n)-r^{n-1}F_k(n-1)=r^{n-1}(n-1)^k),所以(F_k(n)=dfrac{F_k(n-1)+(n-1)^k}{r})

然后我们知道,(k)次多项式的(k+1)阶差分为(0),即

[sum_{i=0}^{k+1}(-1)^iinom{k+1}{i}F_k(k+1-i)=0 ]

于是可以设(F_k(0)=x),然后不停代入,然后用这个式子解方程。具体实现可以写一个二项式binom来做。

注意特判一些特殊情况。时间复杂度(O(k))

code
#include<bits/stdc++.h>
#define Rint register int
using namespace std;
typedef long long LL;
const int N = 1000012, mod = 1e9 + 7;
inline int kasumi(int a, int b){
	int res = 1;
	while(b){
		if(b & 1) res = (LL) res * a % mod;
		a = (LL) a * a % mod; b >>= 1;
	}
	return res;
}
inline void upd(int &a, int b){a += b; if(a >= mod) a -= mod;}
inline int add(int a, int b){return (a + b >= mod) ? (a + b - mod) : (a + b);}
inline int sub(int a, int b){return (a < b) ? (a + mod - b) : (a - b);}
int T, k, fac[N], inv[N];
LL n, r;
struct binom {
	int a, b;
	inline binom(int _a = 0, int _b = 0): a(_a), b(_b){}
	inline binom operator + (const binom &o) const {return binom(add(a, o.a), add(b, o.b));}
	inline binom operator - (const binom &o) const {return binom(sub(a, o.a), sub(b, o.b));}
	inline binom operator * (int v) const {return binom((LL) a * v % mod, (LL) b * v % mod);}
} S[N];
inline int lagrange(int *A, int k, LL n){
	static int pre[N], suf[N];
	pre[0] = suf[k + 1] = 1;
	for(Rint i = 1;i <= k;i ++) pre[i] = (n - i + 1) % mod * pre[i - 1] % mod;
	for(Rint i = k;i;i --) suf[i] = (n - i) % mod * suf[i + 1] % mod;
	int ans = 0;
	for(Rint i = 0;i <= k;i ++){
		int tmp = (LL) inv[i] * inv[k - i] % mod * pre[i] % mod * suf[i + 1] % mod;
		if((k - i) & 1) tmp = mod - tmp;
		upd(ans, (LL) A[i] * tmp % mod);
	}
	return ans;
}
int pri[N], tot, po[N], minp[N];
bool notp[N];
inline void init(int n){
	notp[0] = notp[1] = 1;
	for(Rint i = 2;i <= n;i ++){
		if(!notp[i]) pri[++ tot] = i;
		for(Rint j = 1;j <= tot && i * pri[j] <= n;j ++){
			notp[i * pri[j]] = 1; minp[i * pri[j]] = pri[j];
			if(!(i % pri[j])) break;
		}
	}
	fac[0] = 1;
	for(Rint i = 1;i <= n;i ++) fac[i] = (LL) fac[i - 1] * i % mod;
	inv[n] = kasumi(fac[n], mod - 2);
	for(Rint i = n;i;i --) inv[i - 1] = (LL) inv[i] * i % mod;
}
inline void calc(int k){
	po[0] = 0; po[1] = 1;
	for(Rint i = 2;i <= k + 1;i ++)
		if(notp[i]) po[i] = (LL) po[minp[i]] * po[i / minp[i]] % mod;
		else po[i] = kasumi(i, k); 
}
inline int C(int n, int m){if(n < 0 || m < 0 || n < m) return 0; return (LL) fac[n] * inv[m] % mod * inv[n - m] % mod;}
inline int solve(LL n, int k, int r){
	static int F[N];
	int invr = kasumi(r, mod - 2);
	if(!r) return 0;
	memset(F, 0, (k + 2) << 2);
	calc(k);
	if(r == 1){
		for(Rint i = 1;i <= k + 1;i ++) F[i] = add(F[i - 1], po[i]);
		return n <= k + 1 ? F[n] : lagrange(F, k + 1, n);
	}
	if(!n) return 0;
	if(!k) return (LL) sub(kasumi(r, n + 1), r) * kasumi(r - 1, mod - 2) % mod;
	S[0] = binom(1, 0);
	for(Rint i = 1;i <= k + 1;i ++)
		S[i] = (S[i - 1] + binom(0, po[i - 1])) * invr;
	binom qwq;
	for(Rint i = 0;i <= k + 1;i ++)
		if(i & 1) qwq = qwq + S[k + 1 - i] * C(k + 1, i);
		else qwq = qwq - S[k + 1 - i] * C(k + 1, i);
	F[0] = (LL) (mod - qwq.b) * kasumi(qwq.a, mod - 2) % mod;
	for(Rint i = 1;i <= k + 1;i ++) F[i] = (LL) add(F[i - 1], po[i - 1]) * invr % mod;
	int ans = lagrange(F, k + 1, n + 1);
	return (LL) sub((LL) ans * kasumi(r, (n + 1) % (mod - 1)) % mod, F[0]);
}
int main(){
	scanf("%d", &T); init(1000010);
	while(T --){
		scanf("%lld%lld%d", &n, &r, &k);
		printf("%d
", solve(n, k, r % mod));
	}
}
原文地址:https://www.cnblogs.com/AThousandMoons/p/11779232.html