AOJ 2344 Multi Ending Story

解法

木DPで解きます.
まず部分木に注目しましょう.
冷静に考えると,この部分木のすべての葉を回収するのには,以下の3通りだけしかありません.この部分木の根を $v$,左右の子をそれぞれ $l, r$ とします.

  • $v$ で quick save を使わず,それぞれの子の部分木でいい感じに quick save を用いて最適に回収する
  • $v$ で quick save した後,まず $l$ の部分木の葉をすべて回収(この時 quick save は用いない)して,その後 $r$ の部分木で quick save を用いつつ最適に回収する
  • 2つ目のパターンで左右を入れ替えたパターン

そこで以下のように定めます.
$dp[v] := v$ が根の部分木にあるすべての葉を quick save をいい感じに用いて最適に回収したときのコスト
$dp_2[v] := v$が根の部分木にあるすべての葉を, $v$ で quick save をしたとき(子の部分木では一切 quick save しない)の回収のコスト.ただし真の根から $v$ までのコストは除くものとする.
$cnt[v] := v$ の部分木に含まれる葉の個数

1つ目のパターンは,$v$ で quick save しないので真の根から $v$ までのコスト $depth[v]$ が余計にかかります.なので,このときのコストは
$dp[l] + dp[r] + depth[v] + 2$
となります.

2つ目(3つ目も同様)のパターンでは,$l$ の部分木の回収に $dp_2[l]$ と,
$v$ から $l$ に $cnt[l]$ 回移動することになるのでその分だけのコストがまずかかります.その後は,$r$ にうつって最適に回収するので,そのコストは $dp[r] + 1$ です.よって総コストは
$dp_2[l] + cnt[l] + dp[r] + 1$
となります.

あとはこれらのうちで最小のものが $dp[v]$ の値となります.

ソースコード

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

using ll = long long;

// スタック拡張します

int N;
ll dp[500001], dp2[500001];
int leaf_cnt[500001];
int g[500001][2];

ll dfs(int v, int depth) {
    for(int i = 0; i < 2; ++i) {
        int to = g[v][i];
        if(to == N) {
            leaf_cnt[v]++;
        } else {
            dfs(to, depth + 1);
            leaf_cnt[v] += leaf_cnt[to];
            dp2[v] += dp2[to];
        }
    }
    dp2[v] += leaf_cnt[v];

    ll& res = dp[v];
    int l = g[v][0], r = g[v][1];
    dp[v] = min({dp2[l] + (leaf_cnt[l] == 0 ? 1 : leaf_cnt[l]) + 1 + dp[r],
                 dp[l] + dp2[r] + (leaf_cnt[r] == 0 ? 1 : leaf_cnt[r]) + 1,
                 dp[l] + depth + dp[r] + 2});
    return dp[v];
}

int main() {
    cin >> N;
    for(int i = 0; i < N - 1; ++i) {
        int a, b;
        cin >> a >> b;
        a -= (a != N);
        b -= (b != N);
        g[i][0] = a;
        g[i][1] = b;
    }
    cout << dfs(0, 0) << endl;    
}

感想

メモリの制約が厳しい上に,普通にスタックオーバーフローするので†スタック拡張†します.
DP部分は普通に難しかった…….

AOJ 2598 Website Tour

解法

移動に時間がかからないので,各強連結成分においては,好きな広告を(制限以内なら)好きな回数だけ見ることが出来ます.これは典型的な個数制限付きナップザック問題です.

あとは強連結成分分解をしたグラフで,トポロジカル順序でDP(DAGの上のDP)をします.
各強連結成分では,先に述べた個数制限付きナップザック問題を解きますが,自己辺が無いようなそれ自身だけからなる強連結成分は,自分の広告を何回も見ることはどうやっても出来ないので,個数制限を 1 と改めておきます.
あとは普通のナップザック問題で,
$dp[i][j] := $ 頂点 $i$ まで見たときに,時間の総和が $j$ 以下であるような広告の見方をしたときの最大値
として,各頂点で答えが求まったら,子に対して $dp[child][j] = dp[i][j]$ とかするだけでよいです.

ソースコード

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

using pii = pair<int, int>;

struct edge {
    int from, to;
};

using edges = std::vector<edge>;
using graph = std::vector<edges>;

void add_edge(graph& g, int from, int to) {
    g[from].push_back((edge){from, to});
}

int scc(graph& G, std::vector<int>& cmp) {
    int V = G.size();
    std::vector<std::vector<int>> g(V), rg(V);
    std::vector<bool> used(V, false);
    std::vector<int> vs;
    cmp.resize(V);
    for(int i=0; i<V; ++i) {
        for(auto e : G[i]) {
            g[i].push_back(e.to);
            rg[e.to].push_back(i);
        }
    }
    std::function<void(int)> dfs = [&g, &vs, &used, &dfs](int v) {
        used[v] = true;
        for(auto i : g[v]) {
            if(!used[i]) {
                dfs(i);
            }
        }
        vs.push_back(v);
    };
    std::function<void(int, int)> rdfs = [&rg, &cmp, &used, &rdfs](int v, int k) {
        used[v] = true;
        cmp[v] = k;
        for(int i : rg[v]) {
            if(!used[i]) {
                rdfs(i, k);
            }
        }
    };
    for(int v=0; v<V; ++v) {
        if(!used[v]) {
            dfs(v);
        }
    }
    std::fill(used.begin(), used.end(), false);
    int k = 0;
    for(int i=vs.size()-1; i>=0; --i) {
        if(!used[vs[i]]) {
            rdfs(vs[i], k++);
        }
    }
    return k;
}

std::vector<std::vector<int>> build_graph(graph const& g, std::vector<int> const& cmp, int K) {
    int V = g.size();
    std::vector<std::set<int>> s(K);
    std::vector<std::vector<int>> res(K);
    for(int i=0; i<V; ++i) {
        for(auto e : g[i]) {
            s[cmp[i]].insert(cmp[e.to]);
        }
    }
    for(int i=0; i<K; ++i) {
        for(auto j : s[i]) {
            if(i != j) {
                res[i].push_back(j);
            }
        }
    }
    return res;
}


int main() {
    int N, M, T;
    while(cin >> N >> M >> T, N) {
        graph g(N);
        vector<int> p(N), t(N), k(N);
        vector<bool> self(N);
        for(int i = 0; i < N; ++i) {
            cin >> p[i] >> t[i] >> k[i];
        }
        for(int i = 0; i < M; ++i) {
            int a, b;
            cin >> a >> b;
            a--; b--;
            add_edge(g, a, b);
            if(a == b) {
                self[a] = true;
            }
        }
        
        vector<int> cmp;
        int const K = scc(g, cmp);
        auto g2 = build_graph(g, cmp, K);
        vector<vector<int>> C(K);
        for(int i = 0; i < N; ++i) {
            C[cmp[i]].push_back(i);
        }

        vector<vector<int>> dp(K, vector<int>(T + 1));
        for(int i = 0; i < K; ++i) {
            if(C[i].size() == 1 && !self[C[i][0]]) {
                k[C[i][0]] = 1;
            }
            for(auto v : C[i]) {
                for(int j = 0; j < t[v]; ++j) {
                    deque<pii> deq;
                    for(int l = 0; l * t[v] + j <= T; ++l) {
                        int val = dp[i][l * t[v] + j] - l * p[v];
                        while(deq.size() > 0 && deq.back().second <= val) {
                            deq.pop_back();
                        }
                        deq.emplace_back(l, val);
                        dp[i][l * t[v] + j] = deq.front().second + l * p[v];
                        if(deq.front().first == l - k[v]) {
                            deq.pop_front();
                        }
                    }
                }
                for(auto to : g2[i]) {
                    for(int j = 0; j <= T; ++j) {
                        dp[to][j] = max(dp[to][j], dp[i][j]);
                    }
                }
            }
        }

        int res = 0;
        for(int i = 0; i < K; ++i) {
            res = max(res, dp[i][T]);
        }
        cout << res << endl;
    }
}

感想

最初何を血迷ったか,DFSしながら各連結成分ごとに dp2[i][j] を計算し,dp[i][j] = max(dp[ch][j - t], dp2[i][t]) みたいに実装しようとしていて O(nT^2) でだめだなーとかアホな事を考えていた.

AOJ 1338 Clock Hands

解法

まず,それぞれの針がなす角を求めましょう.短針が1周するのに$H$時間かかるとします.
今時刻が $h:m:s$ である時,針が上を向いている状態を 0 度とすると,長針,短針,秒針はそれぞれ
$\displaystyle ang_h = \frac{360}{H} \times h + \frac{6}{H} \times m + \frac{s}{10H}$
$\displaystyle ang_m = 6m + \frac{s}{10}$
$ang_s = 6s$
の角度にあります.
簡単にするために,$ang_h = 0$ となるように時計を回転させます.
$ang_h = 0$
$\displaystyle ang_m = 6m + \frac{s}{10} - \left(\frac{360}{H} \times h + \frac{6}{H} \times m + \frac{s}{10H}\right)$
$\displaystyle ang_s = 6s- \left(\frac{360}{H} \times h + \frac{6}{H} \times m + \frac{s}{10H}\right)$
これらの値が負になっていれば,360 を加えて辻褄を合わせます.

秒針と長針の位置関係は二通りありますが,とりあえず $ang_s < ang_m $の場合を考えましょう.
このとき,秒針との短針および長針の間の角はそれぞれ
$\displaystyle ang_s - ang_h = 6s- \left(\frac{360}{H} \times h + \frac{6}{H} \times m + \frac{s}{10H}\right)$
$\displaystyle ang_m - ang_s = 6m - \frac{59s}{10}$
これらが仮に等しいとして等式を立て,式を整理すると,
$\displaystyle s = \frac{n}{d} = \frac{60m(H + 1) + 3600h}{119H - 1}$
が成り立ちます.秒針と長針の位置関係が異なる場合も,似た式になります.

つまり,目的の時刻において,$d$ は $119H - 1$ の約数であるということが分かります.
ここから,$d = 119H - 1$ とおいてしまって,その上で n を少しずつ増やしていく全探索が考えられます.実際,秒針は高々2周程度する間に,目的の時刻は必ず存在します.(最後に $n$ と $d$ を gcd で割れば良い)

浮動小数点を扱ったり,分数をきちんと持つのは面倒なので,すべての式に $10 \times H \times d$ をかけてしまって,一周が $360 \times 10 \times H \times d$ 度であるということにすれば実装が楽です.

計算量も,$H <= 100$ という制約から,十分間に合います.

ソースコード

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

using ll = long long;

int main() {
    int H, h0, m0, s0;
    while(cin >> H >> h0 >> m0 >> s0, H) {
        const ll d = 119 * H - 1;
        ll n = s0 * d;
        ll ho = h0, mo = m0;
        while(true) {
            ll ang_h = 3600 * ho * d + 60 * mo * d + n;
            ll ang_m = 60 * H * mo * d + n * H;
            ll ang_s = 60 * H * n;
            ang_m = (ang_m - ang_h) % (360 * 10 * H * d);
            ang_s = (ang_s - ang_h) % (360 * 10 * H * d);
            if(ang_m < 0) {
                ang_m += 360 * 10 * H * d;
            }
            if(ang_s < 0) {
                ang_s += 360 * 10 * H * d;
            }
            ang_h = 0;

            if(ang_h != ang_m) {
                ll ang1, ang2;
                if(ang_m < ang_s) {
                    ang1 = 360 * 10 * H * d - ang_s;
                    ang2 = ang_s - ang_m;
                } else {
                    ang1 = ang_s;
                    ang2 = ang_m - ang_s;
                }
                if(ang1 == ang2) {
                    break;
                }
            }

            n++;
            if(n == 60 * d) {
                n = 0;
                mo++;
            }
            if(mo == 60) {
                mo = 0;
                ho++;
            }
            if(ho == H) {
                ho = 0;
            }
        }
        ll g = __gcd(n, d);
        cout << ho << ' ' << mo << ' ' << n / g << ' ' << d / g << endl;
    }
}

感想

THE ICPC って感じの問題だった.

ARC057 D - 全域木

解法

メモ化再帰で解きます.
$rec(C) := $連結成分の状態が $C$ であるような状態から始めたときに,条件を満たすような辺の張り方の総数

$A[i]$ を小さい方から張って順に張っていくとき,$A[i]$に含まれない辺は,それまでにクラスカル法の要領で作られていた連結成分同士を結ばないように辺をはる必要があります.
逆に,$A[i]$ に含まれる辺は,連結成分同士を結ぶように辺を張ることになります.

今 $A[i]$ の辺を張ろうとした時,仮に連結成分の集合が $C = \{C_1, C_2, \ldots, C_k\}$ であったとしましょう.この時,重みが $A[i - 1] + 1$ から $A[i] - 1$ である辺は,$\displaystyle \sum_{i=1}^{k} \frac{|C_i|(|C_i| - 1)}{2} - (N - |C|)$ 本のまだ辺が張られていないところから適当に選ぶ必要があるので,これは単純に順列の総数を求めれば良いです.($N - |C|$ はMSTにすでに使われている辺の数)
次に $A[i]$ を張るときは,連結成分から適当に2個選んでそれらを併合して,再帰的に解きます.例えば $C_i$ と $C_j$ を併合するときは,どの頂点同士を結ぶかで $|C_i| \times |C_j|$ 通りあります.
最後にすべての併合を試したものの総和に,$A[i - 1]$ から $A[i]$ までの辺の張り方をかけると答えになります.

$C$ の状態数は $N = 30$ の時高々5000程度(分割数という)なので,十分間に合います.

実装面では,番兵として $A[i]$ に 0 と $\frac{N(N - 1)}{2}$ を最初と最後に入れておくと楽かも?

ソースコード

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

using pii = pair<int, int>;
using ll = long long;

constexpr ll M = 1e9 + 7;

map<vector<int>, ll> memo;

ll rec(vector<int> C, vector<int> const& A) {
    sort(C.begin(), C.end());
    if(memo.count(C) == 1) {
        return memo[C];
    }

    int const N = A.size() - 1;

    ll t = 0; // 連結成分内で,使うと閉路が出来てしまう辺
    for(auto sz : C) {
        t += sz * (sz - 1) / 2;
    }
    int now = N - C.size(); // 次に加える辺
    t -= A[now];
    ll unuse = 1; // MSTで使わない辺の張り方
    for(int i = 0; i < A[now + 1] - A[now] - 1; ++i, --t) {
        unuse = unuse * t % M;
    }
    if(unuse == 0) { // どうしても A[i] 以外の辺を使うハメになる
        return memo[C] = 0;
    }

    // 連結成分同士の結び方
    if(C.size() == 1) {
        return memo[C] = unuse;
    } else {
        ll res = 0;
        int const m = C.size();
        for(int i = 0; i < m; ++i) {
            for(int j = i + 1; j < m; ++j) {
                vector<int> nC;
                nC.push_back(C[i] + C[j]);
                for(int k = 0; k < m; ++k) {
                    if(k != i && k != j) {
                        nC.push_back(C[k]);
                    }
                }
                res = (res + rec(nC, A) * C[i] * C[j]) % M;
            }
        }
        return memo[C] = (res * unuse) % M;
    }
}

int main() {
    int N;
    cin >> N;
    vector<int> A(N + 1);
    for(int i = 1; i < N; ++i) {
        cin >> A[i];
    }
    A[N] = N * (N - 1) / 2 + 1; // 番兵
    vector<int> C(N, 1);
    cout << rec(C, A) << endl;
}

AOJ 2436 Card

解法

カードの数字同士を足すわけではなくそのまま並べるので,桁に注目したいなーという気持ちになります.
しかも $a[i]$ は高々4桁なので嬉しいです.

こういう問題は,すべての数字をまとめてやるのは大変な気がするので,ある $a[i]$に着目してやってみましょう.
つまり,$a[i]$ を $j$ 桁目に並べるような並べ方が何通りあるかさえわかれば,それを $C$ 通りだとすると,$a[i] \times 10^j \times C$ を答えに加えるということをすればよいことになります.

さて,例えば $X$ が $j$ 桁めにあるようなものをカウントしましょう.
それは以下のような感じになります.
.......($X$ の前に数字が並ぶ) $X$ ($X$ のあとに $j$ 桁並ぶ)........
当たり前ですが,前と後ろで独立に考えられます.まず後ろを考えましょう.
後ろは以下のようなdpで求めることが出来ます.
$dp[i][j] := X$ 以外の $i$ 枚のカードを並べて $j$ 桁にする場合の数(ゼロリーディングを含める)
遷移の式は簡単ですね.
しかしこのままでは,注目する数 $X$,今追加して並べようとする数 $a[j]$,すでに並べている枚数 $k$,すでに出来ている桁数 $l$,の4つでループが回るので間に合いません.
しかしよく考えると,同じ桁数である $a[i]$ と $a[j]$ について dp の値は全く同じであることがわかります.なので,この計算はループの外に逃がすことが出来ます.以降,dp のキーに桁数を追加して説明することとします.

さて,dp テーブルが求まったら,前に並ぶ数字も考えて最終的な答えを求めにいきます.
今注目している数字が $a[i]$ で,$a[i]$ が $d[i]$ 桁であり,後ろで $j$ 枚のカードを使っていて,それが $k$ 桁であり,前に $l$ 枚のカードを並べるとき,$\displaystyle \sum_{l = 0}^{n - 1 - j}(_{n-1-j}C_{l} \times l! \times a[i] \times 10^j \times dp[d[i]][j][k])$
残念ながら,これも4重ループになっていますが,先と同様,同じ桁である $a[i]$ をまとめて処理してやるとオーダーが落とせます.$sum[x] := x$桁である$a[i]$の総和とすれば,$a[i]$ の部分を $sum[i]$ にしてやればOKです.

最後にゼロリーディングの処理ですが,これは最初に0である$a[i]$を1つだけ取り除いたカード列について同じ計算をして,最後に0を先頭にくっつけたものを引く,という処理で実現できます.先頭の0をどの0にするか,を書け忘れないように.

計算量は O(4 * N^3) になります.

ソースコード

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

using pii = pair<int, int>;
using ll = long long;

constexpr ll M = 1e9 + 7;

ll ten_pow[1024];
ll fact[256];
ll comb[256][256];

void init() {
    ten_pow[0] = 1;
    for(int i = 1; i < 1024; ++i) {
        ten_pow[i] = (ten_pow[i - 1] * 10) % M;
    }
    fact[0] = 1;
    for(int i = 1; i < 256; ++i) {
        fact[i] = (fact[i - 1] * i) % M;
    }
    for(int i = 0; i < 256; ++i) {
        for(int j = 0; j < 256; ++j) {
            if(i == j || j == 0) {
                comb[i][j] = 1;
            } else {
                comb[i][j] = (comb[i - 1][j - 1] + comb[i - 1][j]) % M;
            }
        }
    }
}

ll solve(vector<int> const& a, vector<int> const& d) {
    int const n = a.size();

    vector<vector<vector<ll>>> dp(5, vector<vector<ll>>(n, vector<ll>(n * 4 + 1)));
    ll res = 0;
    vector<ll> sum(5);
    for(int i = 0; i < n; ++i) {
        sum[d[i]] = (sum[d[i]] + a[i]) % M;
    }
    for(int i = 0; i < n; ++i) {
        if(dp[d[i]][0][0] == 0) {
            dp[d[i]][0][0] = 1;
            for(int j = 0; j < n; ++j) {
                if(i == j) {
                    continue;
                }
                for(int k = n - 1; k >= 0; --k) {
                    for(int l = 0; l < n * 4 + 1; ++l) {
                        if(dp[d[i]][k][l] == 0) {
                            continue;
                        }
                        dp[d[i]][k + 1][l + d[j]] = (dp[d[i]][k + 1][l + d[j]] + dp[d[i]][k][l] * (k + 1)) % M;
                    }
                }
            }
        }
    }
    for(int i = 0; i < 5; ++i) {
        for(int j = 0; j < n; ++j) {
            for(int k = 0; k < 4 * n + 1; ++k) {
                if(dp[i][j][k] == 0) {
                    continue;
                }
                for(int l = 0; l <= n - 1 - j; ++l) {
                    ll t = (sum[i] * dp[i][j][k]) % M;
                    t = (t * ten_pow[k]) % M;
                    t = (t * comb[n - 1 - j][l]) % M;
                    t = (t * fact[l]) % M;
                    res = (res + t) % M;
                }
            }
        }
    }

    return res;
}

int main() {
    init();

    int n;
    cin >> n;
    int zero_cnt = 0;
    vector<int> a(n);
    vector<int> d(n);
    for(int i = 0; i < n; ++i) {
        cin >> a[i];
        d[i] = to_string(a[i]).size();
        zero_cnt += (a[i] == 0);
    }

    ll res = solve(a, d);
    if(zero_cnt != 0) {
        int idx = find(begin(a), end(a), 0) - begin(a);
        a.erase(a.begin() + idx);
        d.erase(d.begin() + idx);
        ll leading_zero = (zero_cnt * solve(a, d)) % M;
        res = (res - leading_zero + M) % M;
    }

    cout << res << endl;
}

感想

自明なDPから考えて,そこから少しずつオーダー落としていくと楽だった.
でもまあ結構めんどくさかった.

AOJ 2725 Live Programming

問題概要

催し物を時間 $T$ だけ開催することになった.
そこで歌う歌の候補が $N$ 個あり,それぞれの歌はパラメーター $t[i], p[i], f[i]$ を持つ.$t[i]$ は歌の長さ,$p[i]$ は歌の満足度,$f[i]$ は歌の特徴量を表します.これらの歌からいくつか選んで,歌うことにする.
この時,以下の制約を満たすような歌の列 $\{s_1, ..., s_k\}$ で,満足度が最大を取るときの満足度を求めよ.

  • 歌う歌の時間の総和が $T$ 以下である.
  • 最初に歌う歌の満足度は $p[s_0]$ とする.
  • 以降に歌う歌の満足度は,$p[s_i] - (f[s_{i-1}] - f[s_i])^2$ である.$s_{i-1}$ は,直前に歌った歌である.

・制約
1 <= N <= 4000
1 <= T <= 4000

解法

DPっぽい感じがするのでDPでやるんですが,その前に考えるべきことがあります.

今回問題になるのは,「どの歌の集合を選ぶか」と「歌をどの順番で並べるか」の2つあり,このままではしんどいです.
では,仮に歌の集合を決め打ちした場合,それらをどう並べるのがいいか考えてみましょう.これは明らかに,f[i] の大きさでソートされた順番に歌うのが最適です.

つまり,あらかじめ f[i] の値で歌をソートしてしまえば,以下のDPでできることは明白です.
$dp[i][j] := i$ 番目の歌を最後に歌うとして,総時間が $j$ であるときの最大満足度
このDPの遷移は以下のとおりです.
$\displaystyle dp[i][j] = \max_{k < i}\{dp[k][j - t[i]] + p[i] - (f[k] - f[i])^2\}$
式を整理して,さらに(別に置き換えなくてもできますが)max を min に置き換えると以下の式になります.
$\displaystyle dp[i][j] = -\min_{k < i}\{-dp[k][j - t[i]] + (f[k])^2 - 2f[k] \times f[i]\} + p[i] - (f[i])^2$
これはConvexHullTrickチャンスなので高速化でき,計算量 O(NT) になります.

注意すべき点は j のループは大きい方から小さい方に回すことです.これは普通のナップザックと同じです.
こういう細かいミスをなくすには,先に自明なDPを書いてから,あとで高速化するとよいですね.(慣れているなら別)

ソースコード

#include <bits/stdc++.h>
using namespace std;
 
using ll = long long;
 
class convex_hull_trick {
    using T = long long;
 
    bool check(int f1, int f2, T aa, T bb) {
        return (a[f2] - a[f1]) * (bb - b[f2]) >= (b[f2] - b[f1]) * (aa - a[f2]);
    }
    bool check(int f1, int f2, int f3) {
        return (a[f2] - a[f1]) * (b[f3] - b[f2]) >= (b[f2] - b[f1]) * (a[f3] - a[f2]);
    }
 
    T f(int idx, T x) {
        return a[idx] * x + b[idx];
    }
 
public:
    convex_hull_trick() = default;
    // a : sorted sequence
    convex_hull_trick(std::vector<T> const& a_, std::vector<T> const& b_) {
        assert(a_.size() == b_.size());
        for(int i = 0; i < a_.size(); ++i) {
            add_f(a_[i], b_[i]);
        }
    }
 
    void add_f(T aa, T bb) {
        while(a.size() >= 2 && check(a.size() - 2, a.size() - 1, aa, bb)) {
            a.pop_back();
            b.pop_back();
        }
        a.push_back(aa);
        b.push_back(bb);
    }
    T min(T x) {
        assert(a.size() > 0);
        while(a.size() >= 2 && f(0, x) >= f(1, x)) {
            a.pop_front();
            b.pop_front();
        }
        return a[0] * x + b[0];
    }
 
    bool empty() const {
        return a.size() == 0;
    }
 
private:
    std::deque<T> a, b; // functions
};
 
int main() {
    int N, T;
    cin >> N >> T;
    // (f, t, p)
    vector<tuple<ll, ll, ll>> songs;
    for(int i = 0; i < N; ++i) {
        ll t, p, f;
        cin >> t >> p >> f;
        if(t > T) {
            continue;
        }
        songs.emplace_back(f, t, p);
    }
    N = songs.size();
    sort(songs.begin(), songs.end());
 
    vector<vector<ll>> dp(N + 1, vector<ll>(T + 1));
    vector<convex_hull_trick> cht(T + 1);
    ll res = 0;
    for(int i = 0; i < N; ++i) {
        ll fi, ti, pi;
        tie(fi, ti, pi) = songs[i];
        for(int j = T; j >= ti; --j) {
            if(cht[j - ti].empty()) {
                continue;
            }
            dp[i][j] = -cht[j - ti].min(fi) - fi * fi + pi;
            res = max(res, dp[i][j]);
            cht[j].add_f(-2 * fi, -dp[i][j] + fi * fi);
        }
        dp[i][ti] = max(dp[i][ti], pi);
        res = max(res, dp[i][ti]);
        cht[ti].add_f(-2 * fi, -dp[i][ti] + fi * fi);
    }
    cout << res << endl;
}

Codeforces Round #185 (Div. 1) B. Cats Transport

問題概要

N個の地点があり,順に 0 から N-1 とする.地点 i と i + 1 の間の距離は d[i] である.
また,猫が m 匹いて,猫 j は時刻 t[j] に地点 h[j] にあらわれるものとする.それぞれの猫は,餌をもらえるまで出現してからその地点で待機し続ける.
今,猫に餌を上げる人が p 人いて,それぞれ地点 0 から N-1 の方向に単位時間あたり1だけ移動する.(ただし初期位置は自由.つまり負の時刻に地点 0 を出発するような人がいても良い.これは制約に書いてないが…)
移動中に餌を待っている猫と遭遇した場合,その猫に必ず餌をやる.ただし,餌をやる時間は無視できるものとする.
では,全ての猫に餌をやるように人を配置した時に,それぞれの猫が餌を待つ時間の総和を最小化せよ.

・制約
2 <= N <= 2 * 10^5
1 <= M <= 10^5
1 <= p <= 100

解法

まず,今のままでは位置と時間という軸が2つあり面倒なので,これをひとつにまとめます.人は単位時間に1だけすすむので,時刻 t[j] に地点 h[j] に現れるという情報を a[j] := t[j] - pos[j] とまとめてしまいます(pos[j] は地点 h[j] の座標).a[j] は,時刻 a[j] に人を地点 0 に配置すれば,その人は猫 j が現れるよりあとに通過する(つまり餌をあげることができる)ことを意味します.
こうして得られた a[j] を昇順にソートしておきます.
時刻 a[j] に人を地点 0 に配置した場合,それ以前の直近で人が配置された時刻を a[k] をすると,k + 1 <= j となるような全ての猫について餌をやることになります.

ここからDPパートです.
p が小さいので,以下のような dp が自然に思い浮かびます.
dp[i][j] := i 人目まで使い,かつ i 人目を時刻 a[j] に出発させた時の,j 以前の猫の待ち時間の総和の最小値
このdpの更新は以下のようになります.
$\displaystyle dp[i][j] = \min_{k < j} \{dp[i - 1][k] + \sum_{l=k+1}^{j}(a[j] - a[l])\}$
式を変形すると,以下のようにできます.
$\displaystyle dp[i][j] = \min_{k < j} \{dp[i - 1][k] - k \times a[j] - \sum_{l=k+1}^{j}a[l]\} + j \times a[j]$
さらに累積和 $\displaystyle S[j] = \sum_{l=0}^{j-1} a[l]$を用いて整理すると
$\displaystyle dp[i][j] = \min_{k < j} \{dp[i - 1][k] - k \times a[j] + S[k+1] \} + j \times a[j] - S[j+1]$
このままだと O(pM^2) ですが,この形のDPは ConvexHullTrick という方法で高速化できます.(この方法自体の解説は省略します.)
よって O(pM) でこの問題を解くことができます.

ソースコード

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

class convex_hull_trick {
    using T = long long;

    bool check(int f1, int f2, T aa, T bb) {
        return (a[f2] - a[f1]) * (bb - b[f2]) >= (b[f2] - b[f1]) * (aa - a[f2]);
    }
    bool check(int f1, int f2, int f3) {
        return (a[f2] - a[f1]) * (b[f3] - b[f2]) >= (b[f2] - b[f1]) * (a[f3] - a[f2]);
    }

    T f(int idx, T x) {
        return a[idx] * x + b[idx];
    }

public:
    convex_hull_trick() = default;
    // a : sorted sequence
    convex_hull_trick(std::vector<T> const& a_, std::vector<T> const& b_) {
        assert(a_.size() == b_.size());
        for(int i = 0; i < a_.size(); ++i) {
            add_f(a_[i], b_[i]);
        }
    }

    void add_f(T aa, T bb) {
        while(a.size() >= 2 && check(a.size() - 2, a.size() - 1, aa, bb)) {
            a.pop_back();
            b.pop_back();
        }
        a.push_back(aa);
        b.push_back(bb);
    }
    T min(T x) {
        assert(a.size() > 0);
        while(a.size() >= 2 && f(0, x) >= f(1, x)) {
            a.pop_front();
            b.pop_front();
        }
        return a[0] * x + b[0];
    }

private:
    std::deque<T> a, b; // functions
};

using namespace std;

using ll = long long;

constexpr ll INF = 2e18;

int main() {
    ll n, m, p;
    cin >> n >> m >> p;
    vector<ll> d(n - 1);
    for(int i = 0; i < n - 1; ++i) {
        cin >> d[i];
    }
    vector<ll> pos(n);
    for(int i = 0; i < n - 1; ++i) {
        pos[i + 1] = pos[i] + d[i];
    }
    vector<ll> x(m);
    for(int i = 0; i < m; ++i) {
        ll h, t;
        cin >> h >> t;
        x[i] = t - pos[h - 1];
    }
    sort(x.begin(), x.end());
    vector<ll> sum(m + 1);
    for(int i = 0; i < m; ++i) {
        sum[i + 1] = sum[i] + x[i];
    }

    vector<vector<ll>> dp(p + 1, vector<ll>(m, INF));
    vector<convex_hull_trick> cht(p + 1);
    for(int i = 0; i < m; ++i) {
        dp[1][i] = i * x[i] - sum[i];
        cht[1].add_f(-i, dp[1][i] + sum[i + 1]);
    }
    for(int i = 2; i <= p; ++i) {
        for(int j = 0; j < m; ++j) {
            if(j == 0) {
                dp[i][j] = 0;
            } else {
                dp[i][j] = cht[i - 1].min(x[j]) + j * x[j] - sum[j + 1];
            }
            cht[i].add_f(-j, dp[i][j] + sum[j + 1]);
        }
    }
    cout << dp[p][m - 1] << endl;
}

感想

ConvexHullTrickを初めて書いたので大変だった.かなり典型問題なので勉強になりました.