[CF932E] Team Work

[CF932E] Team Work - 第二类斯特林数

Description

求 (sum_{i=1}^n C(n,i)i^k,) (n le 10^9, k le 5000)

Solution

定理：(m^n=sumlimits_{i=1}^mP_m^iS_n^i)

证明：将n个不同的球装进m个不同的箱子里，允许箱子为空的方案数。用第二类斯特林数算（枚举非空的箱子的个数，由于是有标号的，所以还要乘一个排列数）和直接表示为 (m^n)

第二类斯特林数表示把 n 个元素分成 k 个非空集合的方案数，集合内是无序的，集合之间也是无序的。

[sumlimits_{i=0}^nC_n^ii^k\=sumlimits_{i=0}^nC_n^isumlimits_{j=0}^iP_i^jS_k^j\=sumlimits_{i=0}^nC_n^isumlimits_{j=0}^iC_i^jj!S_k^j\=sumlimits_{j=0}^nj!S_k^jsumlimits_{i=j}^nC_n^iC_i^j\=sumlimits_{j=0}^kS_k^jj!C_n^j2^{n-j} ]

因为斯特林数只用刀 (k le 5000)，所以用递推式暴力预处理就可以了 (S(n,k)=S(n-1,k-1)+kS(n-1,k))

#include <bits/stdc++.h>
using namespace std;

#define int long long

const int mod = 1e9 + 7;

namespace math_mod
{
    signed c__[5005][5005];
    int fac__[1000005];

    int qpow(int p, int q)
    {
        return (q & 1 ? p : 1) * (q ? qpow(p * p % mod, q / 2) : 1) % mod;
    }

    int inv(int p)
    {
        return qpow(p, mod - 2);
    }

    int fac(int p)
    {
        if (p <= 1000000)
            return fac__[p];
        if (p == 0)
            return 1;
        return p * fac(p - 1) % mod;
    }

    int __fac(int p)
    {
        return fac(p);
    }

    int ncr(int n, int r)
    {
        if (r < 0 || r > n)
            return 0;
        return fac(n) * inv(fac(r)) % mod * inv(fac(n - r)) % mod;
    }

    void math_presolve()
    {
        fac__[0] = 1;
        for (int i = 1; i <= 1000000; i++)
        {
            fac__[i] = fac__[i - 1] * i % mod;
        }
        for (int i = 0; i <= 5000; 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], c__[i][j] %= mod;
        }
    }

    int __c(int n, int r)
    {
        if (r < 0 || r > n)
            return 0;
        if (n > 5000)
            return ncr(n, r);
        return c__[n][r];
    }
}

using namespace math_mod;

signed s[5005][5005];

signed main()
{
    math_presolve();
    s[1][1] = 1;
    for (int i = 2; i <= 5000; i++)
        for (int j = 1; j <= i; j++)
            s[i][j] = (s[i - 1][j - 1] + 1ll * j * s[i - 1][j]) % mod;
    int ans = 0;
    int n, k;
    cin >> n >> k;
    int C = 1;
    for (int j = 0; j <= k; j++)
    {
        ans += fac(j) * s[k][j] % mod * C % mod * qpow(2, n - j) % mod;
        ans %= mod;
        C = (C * (n - j ) % mod * inv(j + 1) % mod) % mod;
    }
    cout << ans << endl;
}