SP14175 题解

UPD:增加了利用斯特林数性质的正宗推式子做法。

本题首A。

首先有个显然的式子:

\[\begin{aligned}\sum_{i}\sum_{j|i}j^k&=\sum_{j}\sum_{j|i}j^k\\&=\sum_{j}j^k\left[\frac{n}{j}\right]\end{aligned} \]

然后你兴高采烈地准备拿拉格朗日插值去水题的时候,你发现一个悲惨的事实:这屑模数不是质数……

出题人插头地爬!

于是我们考虑找到

\[\sum_{i=1}^{n}i^k \]

的通项公式。

首先我们有个引理:

\(d_{(k)}f(x)\) 表示 \(f(x)\)\(k\) 阶差分,有:

\[\sum_{i=1}^{n}f(i)=\sum_{i=1}^{n}\binom{n}{i}d_{(i-1)}f(1) \]

写个简单的证明:

显然有

\[d_{(i-1)}f(1)=\sum_{k=1}^{i}(-1)^{i-k}f(k)\binom{i-1}{k-1} \]

所以

\[\begin{aligned}\sum_{i=1}^{n}\binom{n}{i}d_{(i-1)}f(1)&=\sum_{i=1}^{n}\binom{n}{i}\sum_{k=1}^{i}(-1)^{i-k}f(k)\binom{i-1}{k-1}\\&=\sum_{k=1}^{n}f(k)\sum_{i=k}^{n}\binom{i-1}{k-1}(-1)^{i-k}\\&=\sum_{k=1}^{n}f(k)\end{aligned} \]

得证!

所以我们就可以计算啦!

\[F(k)=\sum_{i=1}^{n}i^k \]

考虑把每个\(F(k)\)都写成形如

\[\sum_{i=1}^{k}\binom{n}{i}a_{k,i} \]

的形式。

先讲一下怎么计算\(\binom{n}{i}\)

考虑\(O(k\pi(k))\)预处理出\(i!\)的质因子,然后暴力除,复杂度\(O(k\log k)\)

下面考虑怎么计算\(a_{n,k}\)

第一种方法:由定义可知:

\[a_{n,k}=d_{(i-1)}f(1)=\sum_{k=1}^{i}(-1)^{i-k}f(k)\binom{i-1}{k-1} \]

暴力计算。复杂度\(O(k^3+Tk^2\log k\sqrt n)\)

当然如果你不把这玩意展开你可以做到预处理\(O(k^2)\)

第二种方法:考虑手算出几组小数据,然后试图扔到oeis里面。

\[F(0)=\binom{n}{1} \]

\[F(1)=\binom{n}{1}+\binom{n}{2} \]

\[F(2)=\binom{n}{1}+3\binom{n}{2}+2\binom{n}{3} \]

\[F(3)=\binom{n}{1}+7\binom{n}{2}+12\binom{n}{3}+6\binom{n}{4} \]

取出系数1 1 1 1 3 2 1 7 12 6,扔进 oeis,发现它就是 A028246

然而他提供的东西太少了……所以我们还是要写个小程序(雾

翻到 formula 板块,我们惊奇的发现:

\[a_{n,k}=\begin{Bmatrix}n\\k \end{Bmatrix}(k-1)! \]

于是我们考虑在\(O(k^2)\)的时间内预处理出\(a_{n,k}\),然后就可以做了。

总复杂度:\(O(k^2+Tk^2\log k\sqrt n)\)

上面那个式子的推导过程太复杂了……我们需要一种新的推式子的做法。

我们关注一下第二类斯特林数的定义,可以得到一个显然但是很重要的式子:

\[m^n=\sum_{k=0}^{m}\binom{n}{k}\begin{Bmatrix}m\\k \end{Bmatrix}k! \]

\(m\) 求和:

\[\sum_{i=1}^{n}i^k=\sum_{i=1}^{n}\sum_{j=0}^{i}\binom{k}{j}\begin{Bmatrix}i\\j \end{Bmatrix}j! \]

代码:

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll N=11;
ll T,n,m,k;
ll a[N+5][N+5],fac[N+5],pr[5]={2,3,5,7,11},mi[N+5][5];
char ss[1<<17],*A=ss,*B=ss;
inline char gc(){return A==B&&(B=(A=ss)+fread(ss,1,1<<17,stdin),A==B)?-1:*A++;}
inline void read(ll&x){
    x=0;char s=gc();
    while(!isdigit(s))s=gc();
    while(isdigit(s))x=x*10+(s^'0'),s=gc();
}
void print(ll x){
    if(x<0)putchar('-'),x=-x;
    if(x>9)print(x/10);
    putchar(x%10+'0');
}
namespace Solution{
	void init(){
		a[1][1]=1;
		for(ll i=2;i<=N;i++)
			for(ll j=1;j<=i;j++)
				a[i][j]=j*a[i-1][j]+a[i-1][j-1];
		fac[0]=1;
		for(ll i=1;i<=N;i++)fac[i]=fac[i-1]*i;
		for(ll i=1;i<=N;i++){
			ll t=i;
			for(ll j=0;j<5;j++){
				mi[i][j]=mi[i-1][j];
				while(t%pr[j]==0)mi[i][j]++,t/=pr[j];
			}
		}
		for(ll i=1;i<=N;i++)
			for(ll j=1;j<=i;j++)
				a[i][j]*=fac[j-1];
	}
	ll t[10],ans;
	ll C(ll n,ll k){
		ll ret=1;
		for(ll i=0;i<k;i++)t[i]=n-i;
		for(ll i=0;i<5;i++){
			ll cnt=mi[k][i];
			for(ll j=n-n/pr[i]*pr[i];j<k&&cnt;j+=pr[i]){
				while(cnt&&t[j]%pr[i]==0)t[j]/=pr[i],cnt--;
			}
		}
		for(ll i=0;i<k;i++)ret=(__int128)ret*t[i]%m;
		return ret;
	}
	ll calc(ll n,ll k){
		ll ret=0;
		for(ll i=1;i<=k+1;i++)ret=(ret+(__int128)C(n,i)*a[k+1][i]%m)%m;
		return ret;
	}
	void solve(){
		read(n),read(k),read(m);ans=0;
		for(ll l=1,r;l<=n;l=r+1){
			r=n/(n/l);
			ans=(ans+(__int128)n/l*((calc(r,k)-calc(l-1,k)+m)%m)%m)%m;
		}
		print(ans),puts("");
	}
} 
int main(){
	Solution::init();
	for(read(T);T--;){
		Solution::solve();
	}
	return 0;
}
posted @ 2020-08-12 19:51  happydef  阅读(225)  评论(5编辑  收藏  举报