BZOJ 3157/3516 国王奇遇记 [数学]

题意

\[\sum_{i=0}^n i ^ m \times m ^ i \]

其中 \(1 \leqslant n \leqslant 10^9\) ， \(1 \leqslant m \leqslant N\) ，其中的 \(N\) 是因为这个题有三个版本

\[\begin{align*} &bzoj3157 - N \leqslant 200 \\ &bzoj3516 - N \leqslant 1000 \\ &b bzoj4126 - N \leqslant 500000 \\ \end{align*} \]

但我只会 \(O(m^2)\) 的，所以先不提第三题（挖坑待填）

分析

发现 \(n\) 很大，不能直接枚举 \(i\) ，考虑从 \(m\) 入手 递推

记 \(f_{n,a} = \sum_{i=0}^n i^a m^i\) （这里 \(i\) 从 \(0\) 开始方便后面转换，\(i =0\) 对答案没有影响），题目要求 \(f_{n,m}\) ，考虑转移

\[\begin{align*} f_{n,a} &= \sum_{i=0}^n i^a m^i \\ &= m \sum_{i=0}^n i^a m^{i-1} \ \ \ //提出一个m \\ &= m \sum_{i=0}^{n-1} (i+1)^a m^i \ \ \ //用 \ i \ 表示 \ i+1 \\ \end{align*} \]

当 \(i = -1\) 时 \((i+1)^a =0\) 大力不理

再根据二项式定理 \((a+b)^n = \sum_{i=0}^n C_n^i a^i b^{n-i}\)

\[\begin{align*} &= m \sum_{i=0}^{n-1} m^i \ \sum_{k=0}^a C_a^k i^k \\ &= m \sum_{k=0}^a C_a^k \ \sum_{i=0}^{n-1} i^k m^i //交换枚举顺序 \\ &= m \sum_{k=0}^a C_a^k \ f_{n-1,k} \\ &= m \sum_{k=0}^a C_a^k \ (f_{n,k} - n^k m^n) \end{align*} \]

得到了\(O(m^2)\) 的转移方程，完结撒花

好像 \(O(m)\) 的要多项式插值和快速傅里叶变换求解，不会，不想填坑了

等等！完结不了，注意：

1、上式的 \(k\) 有取到 \(a\) 的情况，解个关于 \(f_{n,i}\) 方程就好了

2、\(f[0] = \sum_{i=0}^n m^i\) 等比数列，先一算，特判 \(m=1\) 时不能用等比数列求和公式，直接 \(ans = \sum_{i=1}^n i = \frac{n \times (n+1)}{2}\)

代码

#include<bits/stdc++.h>
using namespace std;

namespace TYC
{
	typedef long long ll;

	const ll p=1e9+7;
	const int N=1010;

	int n,m;
	ll f[N],Pow[N],C[N][N];

	inline ll qpow(ll x,ll tim)
	{
		ll ans=1;
		for(;tim;tim>>=1,x=x*x%p)
			if(tim&1) ans=ans*x%p;
		return ans;
	}

	#define Mod(x) (x)>=p?(x)-p:(x)

	void init()
	{
		C[0][0]=1;
		for(int i=1;i<=m;i++)
		{
			C[i][0]=1;
			for(int j=1;j<=m;j++)
				C[i][j]=Mod(C[i-1][j-1]+C[i-1][j]);
		}
		Pow[0]=1;
		for(int i=1;i<=m;i++) Pow[i]=Pow[i-1]*n%p;
	}

	void work()
	{
		scanf("%d%d",&n,&m);
		if(m==1)
		{
			printf("%lld\n",(ll)(n+1)*n/2%p);
			return;
		}
		init();
		ll t1=qpow(m,n),t2=t1*m%p,inv=qpow(m-1,p-2);
		f[0]=(Mod(qpow(m,n+1)-1+p))*qpow(m-1,p-2)%p;
		for(int i=1;i<=m;i++)
		{
			ll sum=0;
			for(int j=0;j<i;j++)
				sum=Mod(sum+(C[i][j]*(f[j]-Pow[j]*t1%p+p)%p));
			if(m!=1) f[i]=(Pow[i]*t2-m*sum+p)%p*inv%p;
		}
		printf("%lld\n",f[m]);
	}
}

int main()
{
	TYC::work();
	return 0;
}