$N \le 6$ ゲーは確かに面白いかもしれない。

問題概要

長さ $N$ の整数列 $A_{1}, \dots, A_{N}$ が与えらる。

同じく長さ $N$ の整数列 $X$ は、各 $i$ について独立に、 $1 \le X_{i} \le A_{i}$ を満たす整数の一様分布からランダムに選ばれる。

このとき、 $X$ の最長増加部分列の長さの期待値を mod 1000000007 で計算せよ。

制約

$1 \le N \le 6$
$1 \le A_{i} \le 10^{9}$

考えたこと

期待値と言われると一瞬「期待値の線形性」を活用した解法を考えたくなる。つまり、各要素に対して「それが LIS にどれだけ寄与するか」を考察して、その要素が LIS に寄与する確率を計算して...といった具合だ。でも、「LIS に寄与する確率」とかよくわからない。

そこで方針転換。 $N \le 6$ という制約がすべてで、 $X$ の大小関係を全パターン網羅できる！たとえば $N = 3$ であれば、 $X$ の各項の大小関係は本質的に以下のパターンしかない。

一般の $N$ に対しても、

$N$ 個の対象をグループ分けして
それらを順序づける

方法の個数は順序付きベル数と呼ばれる。 $N = 6$ でも、4683 通りと少ない。そこまでしなくても、上界を $N^{N}$ 通りと見積もることができて、これでも十分少ない。

各大小関係に対しては、

0, 1, 2, ... の順に
その値をもつような添字 i についての $A_{i}$ の最小値をとってきてできる数列を $R$ とする

という風にしたとき、次のような問題に帰着される。

数列 $R_{1}, \dots, R_{K}$ が与えられる。

$1 \le X_{i} \le R_{i}$
$X_{1} \lt \dots \lt X_{K}$

を満たすような数列 $X$ を数え上げよ。

後半パート

後半パートについては、実はほぼまったく同じ問題がこどふぉで出たことがあった！！！違いは

こどふぉ問題では $X$ は広義単調増加
今回の問題では $X$ は狭義単調増加

という部分のみ。

drken1215.hatenablog.com

ここから引っ張ってきた！下の実装でも、計算量は $O(N^{5})$ で十分間に合う (それが 4683 通り)。ただ

二項係数計算を効率化
累積和で高速化

をやれば $O(N^{3})$ まで高速化できる模様。

コード

#include <bits/stdc++.h>
using namespace std;
template<class T> inline bool chmax(T& a, T b) { if (a < b) { a = b; return 1; } return 0; }
template<class T> inline bool chmin(T& a, T b) { if (a > b) { a = b; return 1; } return 0; }

// modint
template<int MOD> struct Fp {
    long long val;
    constexpr Fp(long long v = 0) noexcept : val(v % MOD) {
        if (val < 0) val += MOD;
    }
    constexpr int getmod() const { return MOD; }
    constexpr Fp operator - () const noexcept {
        return val ? MOD - val : 0;
    }
    constexpr Fp operator + (const Fp& r) const noexcept { return Fp(*this) += r; }
    constexpr Fp operator - (const Fp& r) const noexcept { return Fp(*this) -= r; }
    constexpr Fp operator * (const Fp& r) const noexcept { return Fp(*this) *= r; }
    constexpr Fp operator / (const Fp& r) const noexcept { return Fp(*this) /= r; }
    constexpr Fp& operator += (const Fp& r) noexcept {
        val += r.val;
        if (val >= MOD) val -= MOD;
        return *this;
    }
    constexpr Fp& operator -= (const Fp& r) noexcept {
        val -= r.val;
        if (val < 0) val += MOD;
        return *this;
    }
    constexpr Fp& operator *= (const Fp& r) noexcept {
        val = val * r.val % MOD;
        return *this;
    }
    constexpr Fp& operator /= (const Fp& r) noexcept {
        long long a = r.val, b = MOD, u = 1, v = 0;
        while (b) {
            long long t = a / b;
            a -= t * b, swap(a, b);
            u -= t * v, swap(u, v);
        }
        val = val * u % MOD;
        if (val < 0) val += MOD;
        return *this;
    }
    constexpr bool operator == (const Fp& r) const noexcept {
        return this->val == r.val;
    }
    constexpr bool operator != (const Fp& r) const noexcept {
        return this->val != r.val;
    }
    friend constexpr istream& operator >> (istream &is, Fp<MOD>& x) noexcept {
        is >> x.val;
        x.val %= MOD;
        if (x.val < 0) x.val += MOD;
        return is;
    }
    friend constexpr ostream& operator << (ostream &os, const Fp<MOD>& x) noexcept {
        return os << x.val;
    }
    friend constexpr Fp<MOD> modpow(const Fp<MOD> &a, long long n) noexcept {
        if (n == 0) return 1;
        auto t = modpow(a, n / 2);
        t = t * t;
        if (n & 1) t = t * a;
        return t;
    }
};

const int MOD = 1000000007;
using mint = Fp<MOD>;

// Binomial coefficient
template<class T> struct BiCoef {
    vector<T> fact_, inv_, finv_;
    constexpr BiCoef() {}
    constexpr BiCoef(int n) noexcept : fact_(n, 1), inv_(n, 1), finv_(n, 1) {
        init(n);
    }
    constexpr void init(int n) noexcept {
        fact_.assign(n, 1), inv_.assign(n, 1), finv_.assign(n, 1);
        int MOD = fact_[0].getmod();
        for(int i = 2; i < n; i++){
            fact_[i] = fact_[i-1] * i;
            inv_[i] = -inv_[MOD%i] * (MOD/i);
            finv_[i] = finv_[i-1] * inv_[i];
        }
    }
    constexpr T com(int n, int k) const noexcept {
        if (n < k || n < 0 || k < 0) return 0;
        return fact_[n] * finv_[k] * finv_[n-k];
    }
    constexpr T fact(int n) const noexcept {
        if (n < 0) return 0;
        return fact_[n];
    }
    constexpr T inv(int n) const noexcept {
        if (n < 0) return 0;
        return inv_[n];
    }
    constexpr T finv(int n) const noexcept {
        if (n < 0) return 0;
        return finv_[n];
    }
};

// n is big
mint com(long long n, long long r, const BiCoef<mint> &bc) {
    mint res = 1;
    for (long long i = n; i >= n-r+1; --i) res *= i;
    res *= bc.finv(r);
    return res;
}

// #x (x[1] < x[2] < ... < x[N], L[i] <= x[i] < R[i])
mint calc_increasing(vector<long long> &L, vector<long long> &R, const BiCoef<mint> &bc) {
    int N = (int)L.size();
    vector<long long> alt = {-1, 1LL<<60};
    for (auto l : L) alt.push_back(l);
    for (auto r : R) alt.push_back(r);
    sort(alt.begin(), alt.end());
    alt.erase(unique(alt.begin(), alt.end()), alt.end());
    int M = alt.size();
    vector<vector<mint>> dp(N+1, vector<mint>(M+1, 0));
    dp[0][0] = 1;
    for (int i = 1; i <= N; ++i) {
        for (int j = 0; j < M; ++j) {
            for (int pj = 0; pj < j; ++pj) {
                for (int pi = i-1; pi >= 0; --pi) {
                    if (L[pi] >= alt[j+1] || R[pi] <= alt[j]) break;
                    long long num = i - pi;
                    long long width = alt[j+1] - alt[j];
                    mint fac = com(width, num, bc); // com(width + num - 1, num, bc);
                    dp[i][j] += dp[pi][pj] * fac;
                }
            }
        }
    }
    mint res = 0;
    for (int j = 0; j <= M; ++j) res += dp[N][j];
    return res;
}

template<class T> int LIS(vector<T> a,  bool is_strong = true) {
    const T INF = 1<<30; // to be set appropriately
    int n = (int)a.size();
    vector<T> dp(n, INF);
    for (int i = 0; i < n; ++i) {
        if (is_strong) *lower_bound(dp.begin(), dp.end(), a[i]) = a[i];
        else *upper_bound(dp.begin(), dp.end(), a[i]) = a[i];
    }
    return lower_bound(dp.begin(), dp.end(), INF) - dp.begin();
}

void enumerate(int N, vector<int> &vec, vector<vector<int>> &res) {
    if (vec.size() == N) {
        set<int> se;
        for (auto v : vec) se.insert(v);
        int ma = *se.rbegin();
        for (int i = 0; i <= ma; ++i) if (!se.count(i)) return;
        res.push_back(vec);
        return;
    }
    for (int i = 0; i < N; ++i) {
        vec.push_back(i);
        enumerate(N, vec, res);
        vec.pop_back();
    }
}

int main() {
    BiCoef<mint> bc(210000);
    int N;
    cin >> N;
    vector<long long> A(N);
    for (int i = 0; i < N; ++i) cin >> A[i];
    vector<vector<int>> pats;
    vector<int> vec;
    enumerate(N, vec, pats);
    
    mint res = 0;
    for (auto pat : pats) {
        long long lis = LIS(pat);
        map<int, long long> ma;
        for (int i = 0; i < pat.size(); ++i) {
            if (ma.count(pat[i])) chmin(ma[pat[i]], A[i]);
            else ma[pat[i]] = A[i];
        }
        vector<long long> L(ma.size(), 0), R;
        for (auto it : ma) R.push_back(it.second);
        mint num = calc_increasing(L, R, bc);
        res += num * lis;
    }
    mint all = 1;
    for (int i = 0; i < N; ++i) all *= A[i];
    cout << res / all << endl;
}