伯努利数

Posted 2020-12-05 acjiumeng

是用来算 (sumi^k) 的
https://blog.csdn.net/acdreamers/article/details/38929067
先预处理逆元组合数

void init()
{
    inv[1]=1;
    for(ll i=2;i<N;i++)
        inv[i]=(mod-mod/i)*inv[mod%i]%mod;
    for(int i=0;i<N;i++)
    {
        c[i][0]=c[i][i]=1;
        for(int j=1;j<i;j++)
            c[i][j]=(c[i-1][j]+c[i-1][j-1])%mod;
    }
}

然后预处理伯努利数b

ll getb(int n)
{
    if(b[n]!=-1)return b[n];
    if(n==0)return b[0]=1;
    ll ans=0;
    for(int i=0;i<n;i++)
        add(ans,c[n+1][i]*getb(i)%mod);
    ans=-ans*inv[n+1]%mod;
    ans=(ans%mod+mod)%mod;
    return b[n]=ans;
}

最后是求(sumn^k)

ll get(ll n,ll k)
{
    ll ans=0;
    for(int i=1;i<=k+1;i++)
        add(ans,c[k+1][i]*b[k+1-i]%mod*qp(n+1,i)%mod);
    return ans*inv[k+1]%mod;
}

完整代码

//#pragma comment(linker, "/stack:200000000")
//#pragma GCC optimize("Ofast,no-stack-protector")
//#pragma GCC target("sse,sse2,sse3,ssse3,sse4,popcnt,abm,mmx,avx,tune=native")
//#pragma GCC optimize("unroll-loops")
#include<bits/stdc++.h>
#define fi first
#define se second
#define mp make_pair
#define pb push_back
#define pi acos(-1.0)
#define ll long long
#define vi vector<int>
#define mod 1000000007
#define C 0.5772156649
#define ls l,m,rt<<1
#define rs m+1,r,rt<<1|1
#define pil pair<int,ll>
#define pli pair<ll,int>
#define pii pair<int,int>
#define cd complex<double>
#define ull unsigned long long
#define base 1000000000000000000
#define fio ios::sync_with_stdio(false);cin.tie(0)
inline void add(ll &a,ll b){a+=b;if(a>=mod)a-=mod;}
inline ll gcd(ll a,ll b){return b?gcd(b,a%b):a;}
inline ll qp(ll a,ll b){ll ans=1;while(b){if(b&1)ans=ans*a%mod;a=a*a%mod,b>>=1;}return ans;}
 
using namespace std;
 
const double g=10.0,eps=1e-12;
const int N=2000+10,maxn=1500+10,inf=0x3f3f3f3f,INF=0x3f3f3f3f3f3f3f3f;
 
int n;
ll a[N],c[N][N],b[N],inv[N];
void init()
{
    inv[1]=1;
    for(ll i=2;i<N;i++)
        inv[i]=(mod-mod/i)*inv[mod%i]%mod;
    for(int i=0;i<N;i++)
    {
        c[i][0]=c[i][i]=1;
        for(int j=1;j<i;j++)
            c[i][j]=(c[i-1][j]+c[i-1][j-1])%mod;
    }
}
ll getb(int n)
{
    if(b[n]!=-1)return b[n];
    if(n==0)return b[0]=1;
    ll ans=0;
    for(int i=0;i<n;i++)
        add(ans,c[n+1][i]*getb(i)%mod);
    ans=-ans*inv[n+1]%mod;
    ans=(ans%mod+mod)%mod;
    return b[n]=ans;
}
ll get(ll n,ll k)
{
    ll ans=0;
    for(int i=1;i<=k+1;i++)
        add(ans,c[k+1][i]*b[k+1-i]%mod*qp(n+1,i)%mod);
    return ans*inv[k+1]%mod;
}
int main()
{
    init();
    memset(b,-1,sizeof b);
    getb(2000);
    while(~scanf("%d",&n))
    {
        for(int i=1;i<=n;i++)scanf("%lld",&a[i]);
        sort(a+1,a+1+n);
        a[0]=0;
        ll now=1,ans=0;
        for(int i=1;i<=n;i++)
        {
            if(a[i]!=a[i-1])
            {
                ll te=(get(a[i],n-i+2)-get(a[i-1],n-i+2))%mod;
                te-=(get(a[i]-1,n-i+2)-get(a[i-1]-1,n-i+2))%mod;
                te-=(get(a[i]-1,n-i+1)-get(a[i-1]-1,n-i+1))%mod;
                te=(te%mod+mod)%mod;
                add(ans,now*te%mod);
                now=now*a[i]%mod;
            }
        }
        printf("%lld
",ans);
    }
    return 0;
}
/***********************
2
1 2
5
2 3 3 3 3
***********************/