【loj

目录

• description
• solution
• code
• details

---

description

已知函数 (f(k, x)(k, xin mathbb{N}_+)) 满足：

[f(k, x) = egin{cases} 1 & x = 1\ sum_{i=1}^{x-1}f(k, i) + x^k & x > 1 end{cases} ]

现在给定 (n, k)，求 (f(k, n))。

由于答案很大，你只需计算答案对 (10^9 + 7) 取模后的值。

problem link。

solution

考虑每个 (x^k) 的贡献，可以得到 (f(k, n) = n^k + sum_{i=0}^{n-1}2^{n-i-1}i^k=n^k + 2^{n-1}sum_{i=0}^{n-1}frac{i^k}{2^i})。

记 (a = frac{1}{2})，尝试求 (F_k(n) = sum_{i=0}^{n-1}a^ii^k)。

通过转下降幂构造错位相减。不妨记 (S_k(n) =sum_{i=0}^{n-1}a^{i}i^{underline k})，则 (F_k(n) = sum_{j=0}^{k}{krace j}S_{j}(n))。

[egin{aligned} S_k(n) &=sum_{i=0}^{n-1}a^{i}i^{underline k}\ aS_k(n) &= sum_{i=0}^{n-1}a^{i+1}i^{underline k} = sum_{i=1}^{n}a^i(i-1)^{underline k} \ (1-a)S_k(n) &= sum_{i=1}^{n}a^i(i^{underline k} - (i-1)^{underline k})-a^nn^{underline k} + a^{0}0^{underline k}\ &= aksum_{i=0}^{n-1}a^ii^{underline {k-1}}-a^nn^{underline k}\ S_k(n) &= frac{ak}{1-a}S_{k-1}(n)-frac{a^nn^{underline k}}{1-a}\ end{aligned} ]

当 (k = 0) 时，有 (S_0(n) = sum_{i=0}^{n-1} a^i = frac{1-a^n}{1-a})。

---

如果硬上任意模数 fft 求斯特林数是过不了的，考虑进一步优化。

考虑 (S_{k}(n)) 的形式，发现存在最高次数为 (k) 的多项式 (G_k(x)=sum g_ix^i) 满足 (S_k(n)=G_k(n) imes a^n-G_k(0))。

进一步地，存在最高次数为 (k) 的多项式 (G'_k(x)) 满足 (F_k(n)=G'_k(n) imes a^n-G'_k(0))。

可以处理出 (G'_k(i) = alpha_i imes G_{k}'(0)+eta_i) 中的 (alpha_i, eta_i)。再利用如下的 (k+1) 阶差分公式：

[Delta^{k+1}G'(0)=sum_{i=0}^{k+1}inom{k+1}{i}(-1)^{k+1-i}G'(i)=0 ]

解出 (G'_k(0))，再拉格朗日插值就可以求 (G'_k(n))。

code

#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std;

const int MAXN = 1000000;
const int MOD = int(1E9) + 7;
const int INV2 = (MOD + 1) >> 1;

inline int add(int x, int y) {x += y; return x >= MOD ? x - MOD : x;}
inline int sub(int x, int y) {x -= y; return x < 0 ? x + MOD : x;}
inline int mul(int x, int y) {return (int) (1LL * x * y % MOD);}
int mpow(int b, int p) {
	int ret;
	for(ret=1;p;p>>=1,b=mul(b,b))
		if( p & 1 ) ret = mul(ret, b);
	return ret;
}

char s[MAXN + 5]; int len, n, p, k;
void read() {
	scanf("%s%d", s, &k), len = strlen(s);
	for(int i=0;i<len;i++) {
		n = (10LL*n + s[i] - '0') % (MOD);
		p = (10LL*p + s[i] - '0') % (MOD - 1);
	}
}

int prm[MAXN + 5], pcnt; bool nprm[MAXN + 5];
int fct[MAXN + 5], ifct[MAXN + 5], pwk[MAXN + 5];
void init() {
	fct[0] = 1; for(int i=1;i<=k+3;i++) fct[i] = mul(fct[i - 1], i);
	ifct[k+3] = mpow(fct[k+3], MOD - 2);
	for(int i=k+2;i>=0;i--) ifct[i] = mul(ifct[i + 1], i + 1);
	
	pwk[1] = 1;
	for(int i=2;i<=k+3;i++) {
		if( !nprm[i] ) prm[++pcnt] = i, pwk[i] = mpow(i, k);
		for(int j=1;i*prm[j]<=k+3;j++) {
			nprm[i*prm[j]] = false;
			pwk[i*prm[j]] = mul(pwk[i], pwk[prm[j]]);
			if( i % prm[j] == 0 ) break;
		}
	}
}

int lf[MAXN + 5], rf[MAXN + 5];
int gety(int x, int k, int *y) {
	lf[0] = 1; for(int i=1;i<=k;i++) lf[i] = mul(lf[i - 1], sub(x, i));
	rf[k + 1] = 1; for(int i=k;i>=1;i--) rf[i] = mul(rf[i + 1], sub(x, i));
	
	int ret = 0;
	for(int i=1;i<=k;i++) {
		int coef = mul(mul(lf[i - 1], rf[i + 1]), mul(ifct[i - 1], ifct[k - i]));
		if( (k - i) & 1 ) ret = sub(ret, mul(coef, y[i]));
		else ret = add(ret, mul(coef, y[i]));
	}
	return ret;
}

int a[MAXN + 5], b[MAXN + 5], g[MAXN + 5];
int main() {
	read(), init();
	
	a[0] = 1;
	for(int i=1;i<=k+1;i++)
		a[i] = mul(2, a[i - 1]), b[i] = mul(2, add(b[i - 1], pwk[i - 1]));
	
	int A = 0, B = 0;
	for(int i=0;i<=k+1;i++) {
		int del = mul(ifct[i], ifct[k + 1 - i]);
		if( (k + 1 - i) & 1 )
			A = sub(A, mul(del, a[i])), B = sub(B, mul(del, b[i]));
		else A = add(A, mul(del, a[i])), B = add(B, mul(del, b[i]));
	}
		
	g[0] = mul(mpow(A, MOD - 2), MOD - B);
	for(int i=1;i<=k+1;i++) g[i] = add(mul(a[i], g[0]), b[i]);
	
	int ans = sub(mul(mpow(INV2, p), gety(n, k + 1, g)), g[0]);
	printf("%d
", add(mul(mpow(2, p+(MOD-1)-1), ans), mpow(n, k)));
}

details

lca给的solution 是根据递推式直接构造了一个多项式，然后转成了常系数线性递推（虽然最后算法的形式差不多）。

不过感觉 (sum_{i=0}^{n-1}a_ii^k) 这个形式更普遍一些（不清楚我乱说的）？