数论--p次方前n项和&伯努利数--Sum of Maximum

牛客网暑期ACM多校训练营（第一场）F Sum of Maximum

https://www.nowcoder.com/acm/contest/139/F

给定数组a，a[i]表示第i个位置可以在[1,a[i]]中任意取，求在每种情况下，数组最大值 的和。

即数组最大值max * max为最大值的次数，求和。

解：
排序。

a1,a2,a3, ... ,ai = x,aj = y, ... ,an

考虑到y时，mula = a1 * a2 * ... * ai

z = n - j + 1，表示y及后面还有z个数

max = x + 1，情况有mula * [(x + 1)^z - x^z]，贡献度(x + 1) * mula * [(x + 1)^z - x^z]

            mula表示前面i个数可以在[1,a[k]] (1<=k<=i) 中任意取。

            (x+1)^z表示后面z个数可以在[1,a[k]] (j<=k<=n)中任意取，但是

            -x^z表示要减去后面z个数都在[1,x]中任意取的情况，即此时没有数达到x + 1,最大值不为x +1。

max = x + 2，情况有mula * [(x + 2)^z - (x + 1)^z]，贡献度(x + 2) * mula * [(x + 2)^z - (x + 1)^z]

...

max = y，情况有mula * [y^z - (y - 1)^z]，贡献度y * mula * [y^z - (y - 1)^z]

贡献度求和

当x + 1 < y时，可以化简为

mula * {y^(z+1) - x^(z+1) - f(y-1,z) + f(x-1,z)} //f(n,z) = 1^z + 2^z + 3^z + ... + n^z

 

p次方前n项和f(n,z) = 1^z + 2^z + 3^z + ... + n^z

 

伯努利数Bj

 

[m=0]表示m为0时取1，否则取0

#include <iostream>
#include <cstdio>
#include <algorithm>
using namespace std;
const int maxn = 1000 + 5;
typedef long long ll;
const int mod = 1e9 + 7;

ll B[maxn];
ll C[maxn][maxn];//C(n,k)
ll inv[maxn];
void exgcd(ll a,ll b,ll &x,ll &y,ll &d)
{
    if(b == 0){
        d = a;x = 1;y = 0;
        return ;
    }
    exgcd(b, a % b, y, x, d);
    y -= a / b * x;
}
void init()
{
    //组合数
    C[0][0] = 1;
    for (int n = 1; n < maxn; n ++) {
        C[n][0] = C[n][n] = 1;
        for (int k = 1; k < n; k ++) {
            C[n][k] = (C[n - 1][k - 1] + C[n - 1][k]) % mod;
        }
    }
    //逆元
    ll x,y,d;
    inv[1] = 1;
    for (int i = 2; i < maxn; i ++) {
        inv[i] = (mod - mod / i) * inv[mod % i] % mod;
    }
    //伯努利数//虽然伯努利数是小数，但是这里需要取模，除法变成乘法。因此，伯努利数p/q = p * inv(q) % mod，为整数
    B[0] = 1;
    for (int n = 1; n < maxn; n ++) {//B[2*n+1] = 0
        ll ans = 0;
        for (int k = 0; k < n ; k ++) { //ans -= Com(n,k) * B[k] / (n - k + 1);
            ans = ans - (C[n][k] * B[k]) % mod * inv[n - k + 1] % mod;
            ans = (ans % mod + mod) % mod;
        }
        B[n] = ans;
    }
}

int a[maxn];
int n;
ll qk_pow(ll a,ll b)
{
    ll ans = 1;
    while (b) {
        if(b & 1) ans = ans * a % mod;
        b >>= 1;
        a = a * a % mod;
    }
    return ans;
}
ll f(int n,int p)//f(n) = 1 + 2^p + 3^p + 4^p + ..
{
   ll ans = 0;
        int sf = 1;
        for (int j = 0; j <= p; j ++) {
            ans = ans + (sf * C[p + 1][j] % mod) * B[j] % mod * qk_pow(n,p + 1 - j) ;
            ans = (ans % mod + mod) % mod;
            sf = -sf;
        }
   return ans * inv[p + 1] % mod;
}

int main()
{
    init();
    while (scanf("%d",&n) != EOF) {
        for (int i = 1; i <= n; i ++) {
            scanf("%d",&a[i]);
        }//a[0] = 0
        sort(a + 1,a + n + 1);
        ll ans = 0;
        ll mula = 1;
        int p = 1;
        for (int i = 1; i <= n; i = p) {
            while(p <= n && a[p] == a[i]) p ++;
            int z = n - i + 1;
            
            /*公式化简*//*化简后发现需要求sum(k^p),通过伯努利系数计算*/
            int x = a[i - 1],y = a[i];
            ll t = 0;
            
            if(x + 1 < y){//t = y^(z+1) - x^(z+1) - f(y-1,z) + f(x-1,z)
               t = (qk_pow(y, z + 1) - qk_pow(x, z + 1) ) % mod;
                t = ((t - f(y - 1,z) % mod) + f(x - 1,z) % mod + mod) % mod;
            }
            else if(x + 1 == y){//t = (x+1)^(z+1) - x^(z+1) - x^z
                t = qk_pow(x +1 , z + 1) - qk_pow(x, z + 1);
                t = (t % mod + mod) - qk_pow(x, z);
                t = (t % mod + mod) % mod;
            }
           ans =( mula % mod  * t % mod + ans) % mod;
            //i = p的话中间a[i]的连乘会漏掉
            for (int j = i; j < p; j ++ ) {
                mula = mula * a[j] % mod;
            }

        }
        printf("%lld\n",ans);
    }
    return 0;
}