迷いませんか?

継続しないを座右の銘に

戻すDPわからない

はてなでの数式の扱い良くわからなくて意味わからないことになっている

戻すDPわからない

随分前から名前を見てたものの学ぼうとしてなかった
コドフォに出たことでTLでめっちゃ見えたので流石にやろうと思ったら理解に苦しんだので自分のためだけにブログを書く
実際には深く理解してないので意味のわからないことしか書けないけど

yukicoderのこれを見て理解できる人は すごいと思う.これが理解できるような人になりたい人生だった.

概要

DPで数え上げをするときに特定の要素を特別扱いしたくなることがある
そのとき,要素の順番を変えてもDPの結果が同じになるような場合ではその要素がなかったときの結果が求められる
これは例えばdp[i][j] := i番目まで見てi-1番目がj個残ってるときみたいな定義になるとダメ,多分(自信がない)
この場合でも当然,N個目を取り除いた場合が知りたければN-1個目までしか 処理しなければ求められるけど,戻すDPが適用できる場合はそれが何個目の要素のときでもできるよというのが本質

戻すDPが適用できる条件は次の2つ
1. 要素の順番が数え上げに関係しない
2. DPの漸化式を式変形することでdp[i] = dp[i+1] * ... + ...みたいな雰囲気ができる

このとき,2.ができるためには情報が消えてはいけないので遷移に場合分けがあったらダメ

3問戻すDPを使う問題を解いたのでそのときの気持ちと解法を下に書く

注文の多い高橋商店

とりあえず簡単すぎると言われている方が解けないことには話にならないので解く.
僕ぐらいになるとここで詰まってしまうから悲しい
単純に考えるとDPの遷移は  {dp_{i, j} = i種類目まで見てj個選んでいるときの通り数} として一行目のようになっている
このままだと計算量大爆発( {O(NM^2)}ぐらい)なため,累積和を使うことで削減できるがめんどくさい
ここで {j-1}で計算したものが二行目のようになっていることを考えると引いて式変形することで三行目のような遷移になり単純になる
このとき,漸化式は左辺に {dp_{i+1, j}}が来るようなもらう形になってるとなんか戻しやすい気がする.そんな関係ないかも

{\displaystyle
dp_{i+1, j} = \sum_{k=0}^{a_i} dp_{i,j-k} \\
dp_{i+1, j-1} = \sum_{k=0}^{a_i} dp_{i,j-k-1} \\
dp_{i+1, j} - dp_{i+1, j-1} = \sum_{k=0}^{a_i} dp_{i,j-k} - \sum_{k=0}^{a_i} dp_{i,j-k-1} \\
\Leftrightarrow dp_{i+1, j} = dp_{i, j} - dp_{i, j-a_i-1} + dp_{i+1, j-1}
}

よって普通にDPができたので,これを戻せれば勝ち
ちなみに添字が0未満になるところは0として考える.これは場合分けになるから戻せなくなるのでは?みたいな気持ちになるけど0未満になる場合はありえなくて実際0通りなので問題ない.
僕はget(dp, i, j)みたいな感じで呼んだらどっちか0未満なら0を返すみたいなのを用意した

戻した結果として {dp2_{i, j} = i種類目を使わないでj個選ぶ通り数}となるdp2を求める
ここで漸化式を式変形して左辺に {dp_{i, j}}が来るようにすると {dp_{i, j} = dp_{i+1, j} - dp_{i+1, j-1} + dp_{i, j-a_i-1}}となる.
なので,i番目を戻したい(除いた場合を考えたい)ときに行う遷移は {dp2_{i, j} = dp_{N, j} - dp_{N, j-1} + dp2_{i, j-a_i-1}}となり,右辺で {dp2}を参照するときの添字を考えて {j}が小さい方から回すことになる.
 {dp2}が求まったので,クエリの {k_i}種類目をちょうど {x_i}個使う時というのは {k_i}種類目を使わないで {M-x_i}個選ぶのと同じであることから答えは {dp2_{k_i, M-x_i}}となる.
コードは下のやつなんですが,この問題だけModIntを使っていない

int main(void){
    for (auto& x : dp) for (auto& y : x) y = 0;
    for (auto& x : dp2) for (auto& y : x) y = 0;
    cin >> N >> M >> Q;
    a.resize(N);
    for (auto& x : a) cin >> x;
    dp[0][0] = 1;
 
    REP(i, N){
        REP(j, M+1){
            dp[i+1][j] = (get(dp, i+1, j-1) - get(dp, i, j-a[i]-1) + get(dp, i, j))%mod;
        }
    }
    REP(i, N){
        REP(j, M+1){
            dp2[i][j] = (get(dp, N, j) - get(dp, N, j-1) + get(dp2, i, j-a[i]-1))%mod;
            dp2[i][j] = (dp2[i][j] % mod + mod) %mod;
        }
    }
    while(Q--){
        int64 k, x;
        scanf("%lld %lld", &k, &x); k--;
        printf("%lld\n", dp2[k][M-x]);
    }
}

生放送とBGM

これは戻すDPなのかなあとなるけど戻すDPのところで紹介されてたので戻すDPだと思いながら考える(最悪)
特定の曲に注目して考えることになるので考えると,ある曲が再生されるときの通り数を全ての曲に対して求め足し合わせると答えになることがわかる

特定の曲が再生されるかどうかはその曲が再生され始める時間が {L-1}以下なら良いので {dp_{i} = 曲の合計時間がiとなるような通り数}にして,戻した後に {\sum_{k=0}^{L-1} dp_{k}}を求めれば良さそうとなる
しかしこの問題では全部の順列に対して求める必要があるので,時間が {k}になる組み合わせではなく順列が知りたい.更に特定の曲の後に流れる予定の部分についても順番が区別されるため,この両方を求めるために特定の曲までに何曲流したのかも欲しい.

なので {dp_{i, j} = i曲流した合計時間がj}としてDPを行う.
これを戻すと {dp2_{k, i, j} = k曲目を除いてi曲流したときの合計時間がj} という配列が手に入るため,答えは  {\sum_{k=0}^{N} \sum_{i=0}^{N-1} \sum_{j=0}^{L-1} dp2_{k, i, j} \times i! \times (N-i-1)!} を求めた後に期待値なので {N!}で割ればいい
僕は初心者なのでmodの問題でもないのにどうやってそんな大きい階乗をとれば...と思ってたけどdoubleでやれば大丈夫だった(それはそう)

肝心のDPの遷移と戻し方だが,普通のDPは {dp'_{i, j} = dp_{i, j} + dp_{i-1, j-S_i}}でできる.
このときに配るDPをやろうとして,更に時間の合計がL以上になる部分はまとめていいだろ〜みたいな気持ちになり  {dp'_{i+1, min(L, j+S_i)} = dp_{i+1, min(L, j+S_I)} + dp_{i, j}} とかやると, {min}が実質的に場合分けになっているので戻せなくなる.この場合は戻すときに {dp_{N, L}} にどこから足されたかがわからなくなってしまうため
どちらの場合でも,この問題では一回戻すため一回戻せるだけの余裕をもたせる必要がある.

一曲の時間の最大値を {S_{max}} とおくと配るDPなら {min} のとり方を {min(L+S_{max}, j+S_i)} にすれば良くて(普通のDPをする際には {j} {L-1} までしか回す必要がないので配列を十分にとっておけばminを取らなくてもいい)
もらうDPなら {j} {L+S_{max}} まで回せばいい

このとき,戻すDPの遷移は {dp2_{k, i, j} = dp_{N, j} - dp2_{k, i-1, j-S_k}} でできる
今回も {j} が小さい方から回すことになる

コード,initでfactを初期化している.
さっきのもそうだったけどdpの添字や配列名があってないのは空間計算量削減とかあるので許してください.

int main(void){
    init();
    cin >> N >> L;
    L *= 60;
    S.resize(N);
    REP(i, N){
        int64 m, s;
        scanf("%d:%d", &m, &s);
        S[i] = m*60+s;
    }
    vector<vector<int64>> dp(N+1, vector<int64>(L+60*60+1, 0));
    dp[0][0] = 1;
    REP(i, N){
        vector<vector<int64>> dp2(N+1, vector<int64>(L+60*60+1, 0));
        REP(j, N+1){
            REP(k, L+60*60+1){
                dp2[j][k] = get(dp, j-1, k-S[i]) + get(dp, j, k);
            }
        }
        swap(dp, dp2);
    }
    double res = 0;
    REP(i, N){
        vector<vector<int64>> dp2(N+1, vector<int64>(L+60*60+1, 0));
        REP(j, N+1){
            REP(k, L+60*60+1){
                dp2[j][k] = dp[j][k] - get(dp2, j-1, k-S[i]);
                if (k < L && N-j-1 >= 0)
                    res += dp2[j][k] * fact[j] * fact[N-j-1] / fact[N];
            }
        }
    }
    cout << fixed << setprecision(10) << res << endl;
}

Destroy the Colony

問題概要

僕が読解に苦しんだのでこれだけ問題概要
偶数長の英文字(大小)からなる文字列 {S}が与えられる
これの並び替えのうち,文字列を前後半に分けた時に同じ文字が同じ側にある並び替えは嬉しい(Destroy the Colonyできる)
このとき次のクエリが与えられるので処理しなさい

  •  {S_x} {S_y}が同じ側になるような並び替えは何通りあるか出力せよ

 {S_x} {S_y}が同じなら普通に全部の嬉しい並び替えを答えればいいみたいな感じ
クエリの種類が実質 {52^2}種類しかない( {S_x} {S_y}が逆でも答え一緒なので実際は半分強)ので先に計算しておけば勝ち

明らかに特別扱いしたいものがあるので戻すDPだなとなるが,2つあるので当然2回戻す

答えの求め方が僕の考えたTLEするやつと解説の方法がある.
解説のほうが明らかに明快(僕には難しくて難解だった)なのでそっちを先に書いて,僕の考えたやつは最後に書く(読みたい人は読んでって感じ)

解説解

それぞれの文字種がいくつ現れるかを {c_i}とする. この時, {|S|/2} 個ずつになるように分ける組み合わせがわかれば,そのときの組み合わせを {x} とすると,前半と後半それぞれについて重複順列を求めてかけ合わせたものに {x} をかければいいので下のようになる

 {\displaystyle
ans = \frac{x \times (|S|/2)! \times (|S|/2)!}{\prod_{i=0}^{51} c_i!}
}

 {x}以外の部分はどういう分け方でも共通なのであらかじめ計算しておく
じゃあ {x}を求めればいいねとなり,これは上の二問よりも簡単では?みたいな気持ちになる
ここで {dp_{i, j} = i番目まで見て前半にj個あるときの分け方の通り数}とすると {dp_{i+1, j} = dp_{i, j} + dp_{i, j-c_i}}となる.
これを戻す形にすると {dp_{i, j} = dp_{i+1, j} - dp_{i, j-c_i}}になる
 {j}の範囲としては,2回戻す分必要かつそれぞれの文字種の個数が決まってないので {|S|}まで回す必要がある?(ないかも)

文字種 {a}で一回戻したときを {dp1},文字種 {b}で更に戻したときを {dp2}とすると, {dp1_{i} = dp_{N, i} - dp1_{i-c_a}} {dp1}がわかり,全く同じことをして {dp2_{i} = dp1_{i} - dp2_{i-c_b}} {dp2}が求まる
 {a=b}なら {x = dp1_{|S|/2} = dp1_{|S|/2-c_a}}, そうでなければ {x = dp2_{|S|/2} = dp2_{|S|/2-c_a-c_b}}となる
 {a}について {dp1}を求め,それを用いて各 {b}について {dp2}を求めればいいので,計算量は使用されている文字の種類数を {C}とすると {O(C|S| + C^2|S|)}となり, {O(C^2|S|)}である.
戻すときの処理が全く同じなのでそこを関数化すると簡単にかけるらしい

コード,僕はそれをしていない
MintがModIntで,fというのはFactorialライブラリ,f.comb(a, b) {_aC_b}を返す.
めっちゃコードが長くて何が起こってるかわからないね
これが1970msぐらい

int main(void){
    vector<int64> c(52, 0);
    cin >> s >> q;
    n = s.size();
    REP(i, s.size()) c[ctoi(s[i])]++;

    Mint W = f[n/2]*f[n/2];
    REP(i, 52) W /= f[c[i]];

    vector<Mint> dp(n+1, 0);
    dp[0] = 1;
    REP(i, 52){
        if (c[i] == 0) continue;
        vector<Mint> dp2(n+1, 0);
        REP(j, n+1){
            dp2[j] = get(dp, j-c[i]) + get(dp, j);
        }
        swap(dp, dp2);
    }
    vector<vector<Mint>> ans(52, vector<Mint>(52, 0));
    REP(i, 52){
        if (c[i] == 0) continue;
        vector<Mint> t1(n+1, 0);
        REP(j, n+1){
            t1[j] = get(dp, j) - get(t1, j-c[i]);
        }
        REP(j, 52){
            if (c[j] == 0) continue;
            if (i == j) {
                ans[i][j] = t1[n/2];
                continue;
            }
            vector<Mint> t2(n+1, 0);
            REP(k, n+1) {
                t2[k] = get(t1, k) - get(t2, k-c[j]);
            }
            ans[i][j] = t2[n/2];
        }
    }
    while(q--){
        int32 x, y;
        scanf("%d%d", &x, &y); x--; y--;
        x = ctoi(s[x]); y = ctoi(s[y]);
        cout << ans[x][y] * W * 2 << endl;
    }
}

僕の解法

DPは万能なので,DPで組み合わせだけ求めてそれに並び替えをかけなくても最初から並び替えた後を求められるじゃ~んみたいな気持ちになっていた
つまり, {dp_{i, j} = i種類目まで見て前半にj個あるときの前半後半含めて並べ方} をした.
ここで,既に {i}個あるところにこれを並び替えず同じ種類のものを {j}個追加したときの並べ方は,ある種類のもの {i}個と別の種類のもの {j}個の並べ方と同じなので {\frac{(i+j)!}{i!j!} = _{i+j}C_{j}}となる.

 {sum_i = \sum_{k = 0}^{i-1} c_i}とする
前半に {i}種類目を追加するか後半に追加するかの二通りの遷移があるので合わせて {dp_{i+1, j} = dp_{i, j-c_i} \times _jC_{c_i} + dp_{i, j} \times _{sum_i-j+c_i}C_{c_i}} となる
これを戻す形にすると下のような形になる.
戻す文字種を先ほどと同じように {a, b}とすると,文字種を並び替えて一番後ろから {a, b}にして,分母の {sum_i-j+c_i}が一回目は {sum_a+c_a = N}であることから {N-j},二回目は {sum_b+c_b = sum_a-c_b + c_b = N-c_a}となる.
こっちが最初に書いたやつなので答えの部分がすこしごちゃごちゃしてます.

 {\displaystyle
dp_{i, j} = \frac{dp_{i+1, j} - dp_{i, j-c_i} \times _jC_{c_i}}{_{sum_i-j+c_i}C_{c_i}}
}

二回戻した後, {dp_{|S|/2}}が求まったら取り除いている分を追加したときの並び替えを答えないといけないので, {dp_{|S|/2} \times _{|S|/2-c_a}C_{c_b} \times _{|S|/2}C_{c_a}}が答えになります
後は組み合わせが0になってゼロ割してしまわないように気をつけて下のコードのようになります
さっきのが1970msぐらいなので,遷移が定数倍重くなってるこっちは当然TLEします.
お疲れ様でした,多分test9までは通ってるので解法としてはあっているはず

int main(void){
    cin.tie(0);
    ios::sync\_with\_stdio(false);
    vector<int64> c(52, 0), sum(53, 0);
    cin >> s >> q;
    n = s.size();
    REP(i, s.size()) c[ctoi(s[i])]++;
    REP(i, 52) sum[i+1] += sum[i] + c[i];

    vector<Mint> dp(n+1, 0);
    dp[0] = 1;
    REP(i, 52){
        if (c[i] == 0) continue;
        vector<Mint> dp2(n+1, 0);
        REP(j, n+1){
            dp2[j] = dp[j] * f.comb(sum[i]-j+c[i], c[i]) + get(dp, j-c[i]) * f.comb(j, c[i]);
        }
        swap(dp, dp2);
    }
    vector<vector<Mint>> res(52, vector<Mint>(52, 0));
    vector<Mint> t1(n+1, 0);
    vector<Mint> t2(n+1, 0);
    REP(i, 52){
        if (c[i] == 0) continue;
        REP(j, n+1) {
            if (n-j >= c[i]) {
                t1[j] = (get(dp, j) - get(t1, j-c[i]) * f.comb(j, c[i])) / f.comb(n-j, c[i]);
            }
        }
        FOR(j, i, 52){
            if (c[j] == 0) continue;
            if (i == j) {
                res[i][j] = t1[n/2] * f.comb(n/2, c[i]) * 2;
                continue;
            }
            REP(k, n+1){
                if (n-c[i]-k >= c[j]) {
                    t2[k] = (get(t1, k) - get(t2, k-c[j]) * f.comb(k, c[j])) / f.comb(n-c[i]-k, c[j]);
                }
            }
            if (c[i]+c[j] <= n/2)
                res[i][j] = res[j][i] = t2[n/2] * f.comb(n/2-c[i], c[j]) * f.comb(n/2, c[i]) * 2;
        }
    }
    while(q--) {
        int64 x, y;
        cin >> x >> y; x--; y--;
        x = ctoi(s[x]); y = ctoi(s[y]);
        if (x == y) {
            cout << res[x][y] << endl;
        } else if (c[x] + c[y] > n/2) {
            cout << 0 << endl;
        } else {
            cout << res[x][y]  << endl;
        }
    }
}