AOJ1133 Water Tank

方針1

いわゆるやるだけ問題ですが,頑張って通して満足するだけで終わると得られるものがなさそうな問題.
こういうので実装方針を詰める練習をしなければいけないと思う(あくまで僕の意見ですが).

まず,仕切りを外しても状態が変わらないような水槽は統合してしまってよいことがわかる.仕切りを外してよいのは,どちらか一方の水が溢れており,かつ水の高さが等しい時.

次に,シミュレーションするなら1秒毎ではなくて,どこからか水が溢れ始める時点で区切れば良さそうに見える.

溢れた時にどうするかだが,これは蛇口の場所を変えるのと同値であると考えれば見通しが良さそうにみえる.つまり,

  • 水がどちらにもあふれる場合,今見ている水槽にある蛇口を2分割して左右の水槽に分配する
  • 水が左にあふれる場合,今見ている水槽にある蛇口を左に移動させる.右にあふれるときも同様.

蛇口の移動は隣に移動させるだけじゃなくて,連続して移動しうる(移動先がすでに溢れている場合)ので,隣に移す処理を水槽の数だけループを回せば,いずれ蛇口は正しい位置に移動できる.

あとは L 回のクエリを時系列順にソートしておいて,そのタイミングが来たら答えを計算する.

入力サイズはかなり小さいのでよっぽどひどい実装でない限り計算量は気にならない.

ソースコード1

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

constexpr double eps = 1e-6;

struct tank {
    double wall_l, wall_r;
    double in;
    double w, h;
    tank() : wall_l(0), wall_r(0), in(0), w(0), h(0) {}
};

tank merge(tank const& left, tank const& right) {
    tank res;
    res.wall_l = left.wall_l;
    res.wall_r = right.wall_r;
    res.in = left.in + right.in;
    res.w = left.w + right.w;
    res.h = left.h; // (= right.h)
    assert(fabs(left.h - right.h) < eps);
    return res;
}

int main() {
    constexpr double width = 100, height = 50, depth = 30;
    int D; cin >> D;
    while(D--) {
        int N; cin >> N;
        vector<tank> ts(N + 1);
        ts[0].wall_l = ts.back().wall_r = height;
        for(int i = 0; i < N; ++i) {
            int B, H; cin >> B >> H;
            ts[i].w = B;
            ts[i].wall_r = ts[i + 1].wall_l = H;
        }
        ts[N].w = width;
        for(int i = N; i >= 1; --i) {
            ts[i].w -= ts[i - 1].w;
        }
        int M; cin >> M;
        double total_in = 0;
        for(int i = 0; i < M; ++i) {
            int F, A; cin >> F >> A;
            total_in += A;
            int x = 0;
            for(int j = 0; j <= N; ++j) {
                x += ts[j].w;
                if(x > F) {
                    ts[j].in += A / depth;
                    break;
                }
            }
        }

        int L = 0;
        cin >> L;
        vector<tuple<int, int, int>> qs;
        for(int i = 0; i < L; ++i) {
            int P, T; cin >> P >> T;
            qs.emplace_back(T, P, i);
        }
        sort(begin(qs), end(qs));

        vector<double> ans(L);
        double now_t = 0;
        int qi = 0;
        while(qi < L) {
            int T, P, idx;
            tie(T, P, idx) = qs[qi];
            if(T >= width * height * depth / total_in) {
                ans[idx] = height;
                qi++;
                continue;
            }
            double dt = T - now_t;
            for(auto& tk : ts) {
                if(tk.in > 0) {
                    double wall_h = min(tk.wall_l, tk.wall_r);
                    dt = min(dt, (wall_h - tk.h) * tk.w / tk.in);
                }
            }
            for(auto& tk : ts) {
                tk.h += dt * tk.in / tk.w;
            }
            vector<tank> nxt_ts = {ts[0]};
            for(int i = 1; i < (int)ts.size(); ++i) {
                if(fabs(ts[i].h - nxt_ts.back().h) < eps && fabs(ts[i].h - ts[i].wall_l) < eps) {
                    nxt_ts.back() = merge(nxt_ts.back(), ts[i]);
                } else {
                    nxt_ts.push_back(ts[i]);
                }
            }
            ts = move(nxt_ts);
            now_t += dt;
            for(int loop = 0; loop < 10; ++loop) {
                for(int i = 0; i < (int)ts.size(); ++i) {
                    bool to_l = fabs(ts[i].wall_l - ts[i].h) < eps,
                         to_r = fabs(ts[i].wall_r - ts[i].h) < eps;
                    if(to_l && to_r) {
                        ts[i + 1].in += ts[i].in / 2;
                        ts[i - 1].in += ts[i].in / 2;
                        ts[i].in = 0;
                    } else if(to_l) {
                        ts[i - 1].in += ts[i].in;
                        ts[i].in = 0;
                    } else if(to_r) {
                        ts[i + 1].in += ts[i].in;
                        ts[i].in = 0;
                    }
                }
            }
            if(fabs(now_t - T) < eps) {
                double x = 0;
                for(auto& tk : ts) {
                    x += tk.w;
                    if(x > P) {
                        ans[idx] = tk.h;
                        break;
                    }
                }
                qi++;
            }
        }

        for(auto h : ans) {
            cout << fixed << setprecision(6) << h << endl;
        }
    }
}

方針2

さっきのもまあ素直な実装だから辛いと言うほどではないが,もう少し上手く処理できる.
ある蛇口に注目して,初期状態から時間 T だけ水を入れた場合にどうなるか観察する.
小さい視点(直接注いでいる水槽)から視野を広げていって,最終的に全体を眺める感じなので,再帰構造が見えてくる.

pour(pos, l, r, a) := pos 番目の水槽に水量 a だけ流した時,水槽 [l..r) はどうなっているか?また,[l..r) について流しきった時に,a はどれだけ残っているか?

さて,仕切り y の位置で水があふれるということを広い視点で眺めてみる.ここでは左からあふれると仮定する.「水が仕切り y の位置で溢れてくる」ことは「仕切り y から左側で,仕切り y よりはじめて高くなる仕切りまでの水槽を考えた時,それらの水位が仕切り y とおなじになっている」ことと言い換えることができる.
つまり,水槽 [l..r) を眺める時,部分構造として,端を除いた [l..r) の仕切りの中で最も高いものを mid として,[l..mid) と [mid..r) を考えればよいということにほかならない.
そのとき,pos が [l..mid) の方にあるか [mid..r) の方にあるかで,先にどちらを処理するかを決めれば良い.
[l..r) に水を流した後に流すべき水が残っているなら,その水については一段上の再帰段階で処理するものである.

ソースコード2

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

void pour(int pos, int l, int r, double& a, vector<double> const& x, vector<double> const& h, vector<double>& ans) {
    if(r - l > 1) {
        int mid = l + 1;
        for(int i = l + 1; i < r; ++i) {
            if(h[mid] < h[i]) {
                mid = i;
            }
        }
        if(pos < mid) {
            pour(pos, l, mid, a, x, h, ans);
            pour(pos, mid, r, a, x, h, ans);
        } else {
            pour(pos, mid, r, a, x, h, ans);
            pour(pos, l, mid, a, x, h, ans);
        }
    }

    double side_h = min(h[l], h[r]);
    double add = max(0.0, (side_h - ans[l]) * (x[r] - x[l]));
    add = min(add, a);
    a -= add;
    for(int i = l; i < r; ++i) {
        ans[i] += add / (x[r] - x[l]);
    }
}

int main() {
    int D; cin >> D;
    while(D--) {
        int N; cin >> N;
        vector<double> x(N + 2), h(N + 2);
        for(int i = 1; i <= N; ++i) {
            cin >> x[i] >> h[i];
        }
        x[N + 1] = 100;
        h[0] = h[N + 1] = 50;
        int M; cin >> M;
        vector<double> in(N + 1);
        for(int i = 0; i < M; ++i) {
            int F, A; cin >> F >> A;
            int j = upper_bound(begin(x), end(x), F) - begin(x) - 1;
            in[j] += A / 30.0;
        }

        int L; cin >> L;
        while(L--) {
            int P, T; cin >> P >> T;
            vector<double> ans(N + 1);
            for(int i = 0; i < N + 1; ++i) {
                double a = in[i] * T;
                pour(i, 0, N + 1, a, x, h, ans);
            }
            int j = upper_bound(begin(x), end(x), P) - begin(x) - 1;
            cout << fixed << setprecision(6) << ans[j] << endl;
        }
    }
}

感想

本番だったら,実装方針を詰める時間に実装したほうが早いかもしれないので難しいところですが,練習で解くなら実装方針まで考えたい.