51NOD 1258 自然数幂和

题目链接

description

\(T\) 次询问，每次给定 \(n,k\)，求 \(\sum\limits_{i=1}^n i^k\) 模 1e9+7.

\(n\leq 10^{18},k\leq 5\times 10^4\)

solution

可以拉插用什么多项式

考虑将 \(n\) 带入 \(f(x)=\sum\limits_{i=1}^x i^k\)。可以证明，\(f(x)\) 是一个 \(k+1\) 次多项式。

一种感性的理解方式是 \(i^k\) 积分后是个 \(k+1\) 次函数

于是我们可以通过拉格朗日插值法构造出这个多项式在 \(n\) 处的取值。

我们容易线性筛求出 \(f(x)\) 在 \(1,2,\dots,k,k+1,k+2\) 处的点值，根据 \(n+1\) 个点确定一个 \(n\) 次多项式，我们可以通过这些点值构造出这个多项式。

设 \(t=k+2\)，\(s_p=\sum\limits_{i=1}^p i^k\)，对于 \(1\leq p\leq t,p\in \mathbb{Z}\)，构造函数 \(g_p(x)=s_p\prod\limits_{i=1,i\neq p}^t \dfrac{(x-i)}{(p-i)}\)。不难验证，\(\forall x\in[1,t]\cap\mathbb{Z},x\neq p\) 有 \(g_p(x)=0\) 且 \(g_p(p)=s_p\)。

于是可以得出 \(f(x)=\sum\limits_{p=1}^t g_p(x)=\sum\limits_{p=1}^t s_p \prod\limits_{i=1,i\neq p}^t \dfrac{(x-i)}{(p-i)}\)。

于是 \(f(n)=\sum\limits_{p=1}^t s_p \prod\limits_{i=1,i\neq p}^t \dfrac{(n-i)}{(p-i)}\)

线性筛出 \(s_p\)，预处理阶乘以及 \(n-i\) 的前缀后缀积，即可做到严格线性。

时间复杂度 \(O(Tk)\)

code

#include<bits/stdc++.h>

using namespace std;

const int N=(1<<18)+10,mod=1e9+7;

int f[N],n,k,m;
bitset<N> st;

int ksm(int a,int b){
  int ret=1;
  while(b){
    if(b&1) ret=ret*1ll*a%mod;
    a=a*1ll*a%mod;
    b>>=1;
  }
  return ret;
}

vector<int> primes;

int pre[N],sfx[N],jc[N];

void solve(){
  long long nn;
  cin>>nn;
  nn%=mod;
  cin>>k;
  n=nn;

  primes.clear();
  st.reset();
  m=k+2;
  f[1]=1;
  for(int i=2; i<=m; i++){
    if(!st[i]){
      f[i]=ksm(i,k);
      primes.emplace_back(i);
    }
    for(auto j:primes){
      if(i*j>m) break;
      st[i*j]=1;
      f[i*j]=f[i]*1ll*f[j]%mod;
      if(i%j==0){
        break;
      }
    }
  }
  for(int i=2; i<=m; i++){
    f[i]=(f[i-1]+f[i])%mod;
  }

  jc[0]=pre[0]=sfx[m+1]=1;
  for(int i=1; i<=m; i++){
    jc[i]=jc[i-1]*1ll*i%mod;
    pre[i]=pre[i-1]*1ll*(n-i)%mod;
  }

  for(int i=m; i; i--){
    sfx[i]=sfx[i+1]*1ll*(n-i)%mod;
    jc[i]=i==m?ksm(jc[i],mod-2):jc[i+1]*1ll*(i+1)%mod;
  }

  int ans=0;
  for(int i=1; i<=m; i++){
    int t=sfx[i+1]*1ll*pre[i-1]%mod*jc[i-1]%mod*jc[m-i]%mod;
    if((m-i)&1) t=mod-t;
    ans=(ans+t*1ll*f[i]%mod)%mod;
  }

  cout<<ans<<'\n';
}

int main(){

  int T;
  cin>>T;
  while(T--) solve();

  return 0;
}