f(x)=1^x+2^x+3^x+.....+n^x;

阮星火
2023-12-01

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
const double eps=1e-4;
const int MAX=1e6+10;
const ll mod=1e9+7;
const int D=1000105;
ll f[D],g[D],p1[D],p2[D],b[D];

namespace polysum
{
#define rep(i,a,n) for (int i=a;i<n;i++)
#define per(i,a,n) for (int i=n-1;i>=a;i--)
    ll powmod(ll a,ll b)
    {
        ll res=1;
        a%=mod;
        assert(b>=0);
        for(; b; b>>=1)
        {
            if(b&1)res=res*a%mod;
            a=a*a%mod;
        }
        return res;
    }
    ll calcn(int d,ll *a,ll n)   // a[0].. a[d] ?a[n]?
    {
        if (n<=d) return a[n];
        p1[0]=p2[0]=1;
        rep(i,0,d+1)
        {
            ll t=(n-i+mod)%mod;
            p1[i+1]=p1[i]*t%mod;
        }
        rep(i,0,d+1)
        {
            ll t=(n-d+i+mod)%mod;
            p2[i+1]=p2[i]*t%mod;
        }
        ll ans=0;
        rep(i,0,d+1)
        {
            ll t=g[i]*g[d-i]%mod*p1[i]%mod*p2[d-i]%mod*a[i]%mod;
            if ((d-i)&1) ans=(ans-t+mod)%mod;
            else ans=(ans+t)%mod;
        }
        return ans;
    }
    void init(int M)
    {
        f[0]=f[1]=g[0]=g[1]=1;
        rep(i,2,M+5) f[i]=f[i-1]*i%mod;
        g[M+4]=powmod(f[M+4],mod-2);
        per(i,1,M+4) g[i]=g[i+1]*(i+1)%mod;
    }
    ll polysum(ll m,ll *a,ll n)   // a[0].. a[m] \sum_{i=0}^{n-1} a[i]
    {
        for(int i=0; i<=m; i++) b[i]=a[i];
        b[m+1]=calcn(m,b,m+1);
        rep(i,1,m+2) b[i]=(b[i-1]+b[i])%mod;
        return calcn(m+1,b,n-1);
    }
} // polysum::init();
ll quick_mod(ll a,ll b)
{
    ll ans=1;
    while(b)
    {
        if(b&1) ans=(ans*a)%mod;
        a=(a*a)%mod;
        b/=2;
    }
    return ans%mod;
}
ll a[1000005];
int main()
{
    ll n,k;
    polysum::init(1000005);
	while(cin>>n>>k)
    {
        if(k==0) {
	cout<<n<<endl;return 0;}
	a[0]=0;
	for(int i=1;i<=k+2;i++)
	{
		a[i]=(quick_mod(i,k))%mod;
	}
	ll ans=polysum::polysum(k+2,a,n+1)%mod;
	cout<<ans%mod<<endl;
    }

    return 0;
}

 

 类似资料: