Atto Round 1 (Codeforces Round 1041, Div. 1 + Div. 2)

Atto Round 1 (Codeforces Round 1041, Div. 1 + Div. 2)

F

首先，我们需要理解 $f(a)$ 的本质。

$f(a) = (\text{所有最大值之和}) - (\text{所有最大值之后紧邻的元素之和}) - a_1$

我们可以分类算贡献。

我们枚举数组的最大值可能为 $x$（其中 $1 \le x \le m$），然后计算所有最大值为 $x$ 的数组对这三项总和的贡献。

一个关键的子问题是：

有多少个长度为 $n$ 的数组 $a$，满足 $\sum_{i=1}^n a_i = m$，且对所有 $i$ 都有 $0 \le a_i \le x$

我们定义函数 $g(n, m, x)$ 为这个问题的答案。这个问题可以用隔板法和容斥原理解决。

1. 无上限（隔板法）：如果不考虑 $a_i \le x$ 的上限，问题就变成了将 $m$ 个无差别小球放入 $n$ 个有差别盒子，允许空盒。方案数是 $\binom{m+n-1}{n-1}$。

2. 有上限（容斥原理）：我们从无上限的总方案数中，减去至少有一个 $a_i > x$ 的方案数，加上至少有两个 $a_j, a_k > x$ 的方案数，以此类推。

   • 强制 $t$ 个变量大于 $x$（即 $\ge x+1$），我们先给这 $t$ 个变量各分配 $x+1$，剩下 $m-t(x+1)$ 再用隔板法分配。
   • 选择哪 $t$ 个变量有 $\binom{n}{t}$ 种方法。

最终得到的公式为：

$g(n, m, x) = \sum_{t=0}^{\lfloor \frac{m}{x+1} \rfloor} (-1)^t \binom{n}{t} \binom{m - t(x+1) + n - 1}{n - 1}$

这个函数可以在 $\mathcal{O}(\frac{m}{x})$ 的时间内计算。我们可以预处理阶乘和逆元来快速计算组合数。

我们固定一个最大值 $x$，并计算所有最大值恰好为 $x$ 的数组，它们对总和的贡献。

我们对每个 $x \in [0, m]$ 进行迭代。在每次迭代中，我们计算所有以 $x$ 为某个最大值的数组集合的贡献。

令 $g_1 = g(n-1, m-x, x)$ 和 $g_2 = g(n-2, m-2x, x)$。
对于这个数组集（固定 $a_n=x$, $a_i \le x$）：

1. 最大值之和：
   • $a_n$ 贡献：$x \cdot g_1$
   • $a_i=x, i<n$ 贡献：有 $n-1$ 个这样的位置，每个位置成为 $x$ 的数组有 $g_2$ 个。总贡献 $(n-1) \cdot x \cdot g_2$。
   • 总和：$x \cdot g_1 + (n-1)x \cdot g_2$
2. a_1 之和：
   • 在 $a_n=x$ 的前提下，剩下 $n-1$ 个元素和为 $m-x$。由对称性，$a_1$ 的总和是 $\frac{m-x}{n-1} \cdot g_1$。
3. 最大值之后元素之和：
   • 如果 $a_{i-1}=x, i<n$，我们需要 $a_i$ 的和。剩下 $n-2$ 个元素和为 $m-2x$。由对称性，$a_i$ 的和是 $\frac{m-2x}{n-2} \cdot g_2$。有 $n-2$ 个这样的 $i \in [2, n-1]$。总和是 $(m-2x)g_2$。
   • 如果 $a_{n-1}=x$，我们需要 $a_n$ 的和。$a_n$ 固定为 $x$。数组数量为 $g_2$。总和是 $x \cdot g_2$。

将这些项针对每个 $x$ 累加起来，就得到了 $f(a)$ 的总和。

总复杂度为 $\sum_{x=1}^m \mathcal{O}(m/x) = \mathcal{O}(m \log m)$。预处理阶乘需要 $\mathcal{O}(m+n)$。

#include <bits/stdc++.h>
using i64 = long long;
using u64 = unsigned long long;
using i128 = __int128;
#define int i64
#define pb push_back
#define ep emplace
#define eb emplace_back
using namespace std;
int read(int x = 0, int f = 0, char ch = getchar()) {
    while (ch < 48 or 57 < ch) f = ch == 45, ch = getchar();
    while (48 <= ch and ch <= 57) x = x * 10 + ch - 48, ch = getchar();
    return f ? -x : x;
}
#define debug(x) cout << #x << " = " << x << "\n";
#define vdebug(a) cout << #a << " = "; for (auto x : a) cout << x << " "; cout << "\n";
template <class T1, class T2> ostream &operator<<(ostream &os, const std::pair<T1, T2> &a) { return os << "(" << a.first << ", " << a.second << ")"; };
template <class T> ostream &operator<<(ostream &os, const vector<T> &as) { const int sz = as.size(); os << "["; for (int i = 0; i < sz; i++) { if (i >= 256) { os << ", ..."; break; } if (i > 0) { os << ", "; } os << as[i]; } return os << "]"; }
template <class T> void pv(T a, T b) { for (T i = a; i != b; i++) cerr << *i << " "; cerr << '\n'; }
using pii = pair<int, int>;
const int inf = 1e18;

template <class T>
constexpr T power(T a, i64 b) { T res {1}; for (; b; b /= 2, a *= a) if (b % 2) res *= a; return res; }
constexpr i64 mul(i64 a, i64 b, i64 p) { i64 res = a * b - (i64)(1.L * a * b / p) * p; res %= p; if (res < 0) res += p; return res; }
template <i64 P>
struct MInt {
    i64 x;
    constexpr MInt() : x {0} {}
    constexpr MInt(i64 x) : x {norm(x % getMod())} {}
    static i64 Mod;
    constexpr static i64 getMod() { return P > 0 ? P : Mod; }
    constexpr static void setMod(i64 Mod_) { Mod = Mod_; }
    constexpr i64 norm(i64 x) const { if (x < 0) x += getMod(); if (x >= getMod()) x -= getMod(); return x; }
    constexpr i64 val() const { return x; }
    constexpr MInt operator-() const { MInt res; res.x = norm(getMod() - x); return res; }
    constexpr MInt inv() const { return power(*this, getMod() - 2); }
    constexpr MInt &operator*=(MInt rhs) & { if (getMod() < (1ULL << 31)) x = x * rhs.x % int(getMod()); else x = mul(x, rhs.x, getMod()); return *this; }
    constexpr MInt &operator+=(MInt rhs) & { x = norm(x + rhs.x); return *this; }
    constexpr MInt &operator-=(MInt rhs) & { x = norm(x - rhs.x); return *this; }
    constexpr MInt &operator/=(MInt rhs) & { return *this *= rhs.inv(); }
    friend constexpr MInt operator*(MInt lhs, MInt rhs) { MInt res = lhs; res *= rhs; return res; }
    friend constexpr MInt operator+(MInt lhs, MInt rhs) { MInt res = lhs; res += rhs; return res; }
    friend constexpr MInt operator-(MInt lhs, MInt rhs) { MInt res = lhs; res -= rhs; return res; }
    friend constexpr MInt operator/(MInt lhs, MInt rhs) { MInt res = lhs; res /= rhs; return res; }
    friend constexpr std::istream &operator>>(std::istream &is, MInt &a) { i64 v; is >> v; a = MInt(v); return is; }
    friend constexpr std::ostream &operator<<(std::ostream &os, const MInt &a) { return os << a.val(); }
    friend constexpr bool operator==(MInt lhs, MInt rhs) { return lhs.val() == rhs.val(); }
    friend constexpr bool operator!=(MInt lhs, MInt rhs) { return lhs.val() != rhs.val(); }
    friend constexpr bool operator<(MInt lhs, MInt rhs) { return lhs.val() < rhs.val(); }
};
template <>
i64 MInt<0>::Mod = 1e9 + 7;
constexpr int P = 1e9 + 7;
using Z = MInt<P>;

struct Comb {
    int n;
    std::vector<Z> _fac, _invfac, _inv;
    Comb() : n{0}, _fac{1}, _invfac{1}, _inv{0} {}
    Comb(int n) : Comb() { init(n); }
    void init(int m) {
        m = std::min<i64>(m, Z::getMod() - 1);
        if (m <= n) return;
        _fac.resize(m + 1);
        _invfac.resize(m + 1);
        _inv.resize(m + 1);
        for (int i = n + 1; i <= m; i++) _fac[i] = _fac[i - 1] * i;
        _invfac[m] = _fac[m].inv();
        for (int i = m; i > n; i--) _invfac[i - 1] = _invfac[i] * i, _inv[i] = _invfac[i] * _fac[i - 1];
        n = m;
    }
    Z fac(int m) { if (m > n) init(2 * m); return _fac[m]; }
    Z invfac(int m) { if (m > n) init(2 * m); return _invfac[m]; }
    Z inv(int m) { if (m > n) init(2 * m); return _inv[m]; }
    Z binom(int n, int m) { if (n < m || m < 0) return 0; return fac(n) * invfac(m) * invfac(n - m); }
} comb;

Z g(int n, int m, int x) {
    Z res = 0;
    int upper = m / (x + 1);
    for (int t = 0; t <= upper; t++) {
        int coef = t & 1 ? -1 : 1;
        res += coef * comb.binom(n, t) * comb.binom(m + n - 1 - t * (x + 1), n - 1);
    }
    return res;
}

void solve() {
    int n = read(), m = read();
    Z ans = 0;

    for (int x = 0; x <= m; x++) {
        Z g1 = g(n - 1, m - x, x);
        Z g2 = g(n - 2, m - x - x, x);

        // sum of maximums
        ans += x * g1;
        ans += Z(x) * (n - 1) * g2;

        // sum of a_1
        ans -= Z(m - x) / Z(n - 1) * g1;

        // sum of right after maximums
        ans -= (m - 2 * x) * g2;
        ans -= x * g2;
    }

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

signed main() {
    for (int T = read(); T--; solve());
    // solve();
    return 0;
}