CSA Round #82 RBubbleSort

誰もブログを書いてなくて困ったのでメモを書きのこしておきます.

https://csacademy.com/contest/round-82/task/rbubblesort/

1~Nの順列が与えられる.ランダムシャッフルする操作と任意の隣り合う2要素をスワップする操作がK回まで出来る.最適な操作をしたときの反転数の期待値を求めよ.という問題.

K≦10^9という制約だが,K≧反転数ならスワップで反転数を0に出来るので,K<反転数の場合をかんがえる.

f(k+1)=1回シャッフルしたあとk回まで操作できるときの期待値

とすると,シャッフルしたとき反転数がfloor(f(k))+k以下なら,のこりのk回をスワップに使うのが良い.一方,反転数がfloor(f(k))+kより大きいならばまたシャッフルした方がいい.

そうするとt=floor(f(k))+kとしたとき,f(k+1)=(反転数がtより大きい確率)×f(k)+(反転数がk+1以上t以下でない場合は0と見なしたときの期待値)-k×(反転数がk+1以上t以下である確率)と計算できる.ただし答えるのは,答えをa/bとしたときa×b^(-1) mod 10^9+7なので同時に計算しておく.

f(K)>(元の反転数)-Kなら(元の反転数)-K,そうでないならf(k)が答えとなる.

#include "bits/stdc++.h"
using namespace std;
typedef long long ll;
const ll MOD = 1000000007LL;
ll mod_pow(ll a, ll b) {
    if (b == 0) return 1;
    if (b % 2 == 0) {
        ll x = mod_pow(a, b / 2);
        return (x*x) % MOD;
    }
    return (a*mod_pow(a, b - 1)) % MOD;
}
int P[300];
ll inv[301];
double p[301][50100];
double s[301][50100];
double e[301][50100];
ll pp[301][50100];
ll ss[301][50100];
ll ee[301][50100];

void solve() {
    int N, K;
    cin >> N >> K;
    for (int i = 0; i < N; i++) cin >> P[i];
    int inv_n = 0;
    for (int i = 0; i < N; i++) {
        for (int j = i + 1; j < N; j++) if (P[i] > P[j]) inv_n++;
    }
    if (inv_n <= K) {
        cout << 0 << endl;
        return;
    }
    if (K == 0) {
        cout << inv_n << endl;
        return;
    }
    double E = e[N][N*(N - 1) / 2];
    ll EE = ee[N][N*(N - 1) / 2];
    for (int i = 1; i < K; i++) {
        int t = min((int)E + i, N*(N - 1) / 2);
        E = (1 - s[N][t])*E + (e[N][t] - e[N][i]) - (s[N][t] - s[N][i])*i;
        ll tmp = (((1 + MOD - ss[N][t]) % MOD)*EE) % MOD;
        (tmp += ee[N][t] + MOD - ee[N][i]) %= MOD;
        (tmp += MOD - (((ss[N][t] + MOD - ss[N][i]) % MOD)*i) % MOD) %= MOD;
        EE = tmp;
    }
    if (E > inv_n - K) {
        cout << inv_n - K << endl;
    }
    else {
        cout << EE << endl;
    }
}
int main() {
    for (int i = 1; i <= 300; i++) inv[i] = mod_pow(i, MOD - 2);
    p[1][0] = s[1][0] = 1;
    pp[1][0] = ss[1][0] = 1;
    for (int i = 1; i < 300; i++) {
        for (int j = 0; j <= i*(i + 1) / 2; j++) {
            p[i + 1][j] = s[i][min(j, i*(i - 1) / 2)];
            if (j >= i + 1) p[i + 1][j] -= s[i][j - (i + 1)];
            p[i + 1][j] /= (i + 1);
            s[i + 1][j] = p[i + 1][j];
            e[i + 1][j] = p[i + 1][j] * j;
            if (j > 0) {
                s[i + 1][j] += s[i + 1][j - 1];
                e[i + 1][j] += e[i + 1][j - 1];
            }
            pp[i + 1][j] = ss[i][min(j, i*(i - 1) / 2)];
            if (j >= i + 1) (pp[i + 1][j] += MOD - ss[i][j - (i + 1)]) %= MOD;
            (pp[i + 1][j] *= inv[i + 1]) %= MOD;
            ss[i + 1][j] = pp[i + 1][j];
            ee[i + 1][j] = (pp[i + 1][j] * j) % MOD;
            if (j > 0) {
                (ss[i + 1][j] += ss[i + 1][j - 1]) %= MOD;
                (ee[i + 1][j] += ee[i + 1][j - 1]) %= MOD;
            }
        }
    }
    int T;
    cin >> T;
    for (int t = 0; t < T; t++) {
        solve();
    }
}