Day 172-Q4 — 2026-10-08 | 容量スケーリング法の最小費用 b-フロー(Δ を半分にしながら Δ 単位で送る。負の費用・負の閉路も Δ ごとの飽和で処理)

2026-10-08 赤色 Master (Phase 8+) ★★★★★★★★★ フロー(焦点帯 M2)。Day170-Q4 は「負の辺を最初に一度だけ飽和」させてから逐次最短路(SSP)で解いた。今回は容量スケーリング: 増加回数を容量の対数で抑える方法。ただし本セッションの計測では、ランダムなグラフで SSP より速くはならなかった。その結果も含めて扱う

問題

$N$ 頂点 $M$ 辺の有向グラフがある。辺 $i$ は頂点 $u_i$ から $v_i$ へ、容量 $c_i$、単位あたり費用 $w_i$ を持つ($w_i$ は負でもよい)。頂点 $v$ には供給量 $b_v$ がある(正なら湧き出し、負なら吸い込み、$\sum b_v=0$)。

各辺の流量 $0\le f_i\le c_i$(整数)を決め、すべての頂点で「出る量 − 入る量 $=b_v$」を満たすとき、総費用 $\sum w_if_i$ の最小値を出力せよ。条件を満たす流し方が無ければ infeasible と出力せよ。負の費用の閉路があってもよい(容量が有限なので最小値は存在する)。

入力形式

N M
b_1 b_2 ... b_N
u_1 v_1 c_1 w_1
...
u_M v_M c_M w_M

制約

#条件
1$2\le N\le300$、$1\le M\le3000$
2$1\le u_i,v_i\le N$、$u_i\ne v_i$(多重辺あり)
3$0\le c_i\le10^9$、$|w_i|\le10^9$
4$|b_v|\le10^{12}$、$\sum b_v=0$

入出力例

入力例1

3 3
5 0 -5
1 2 4 1
2 3 4 1
1 3 10 3

出力例1

11

(1→2→3 に 4 単位(費用 2×4)、1→3 に 1 単位(費用 3))

入力例2

3 4
0 0 0
1 2 2 -3
2 3 2 1
3 1 1 1
2 1 5 2

出力例2

-2

(供給は 0 でも、負の閉路 1→2→1(費用 −1)に 2 単位流すと −2)

入力例3

2 1
3 -3
1 2 2 0

出力例3

infeasible

(3 例とも流量の全列挙と一致)

背景と出典

容量スケーリングは、最小費用流の「増加の回数」を容量の大きさ $U$ の対数に抑えるための方法である。Edmonds と Karp(1972)がスケーリングの考え方を最小費用流に使い、Ahuja・Magnanti・Orlin の教科書 Network Flows(1993)に「capacity scaling algorithm」として $O(M\log U\cdot(M+N\log N))$ の形で載っている、と私は記憶している(本セッションでは原典を参照していない)。

実務の実装としては、C++ のグラフライブラリ LEMON が CapacityScaling を、NetworkSimplex・CostScaling・CycleCanceling と並べて提供しており、Kovács(2015, Optimization Methods and Software)の比較実験ではネットワーク単体法と費用スケーリングが多くの場合に速かった、と記憶している(これも未確認)。本セッションの計測(後述)でも、ランダムなグラフでは容量スケーリングは SSP と同程度の速さだった。

本問は自作。容量スケーリングの価値は「平均の速さ」ではなく「容量が $10^9$ でも増加回数が容量に比例しない保証」にある、という位置づけで扱う。

概念図

入力例 1 を Δ = 8 → 4 → 2 → 1 の段階で解く(Δ 残余グラフ = 残余容量 Δ 以上の辺だけ) Δ = 8残余 ≥8 は 1→3 (10) だけ余剰 5 < 8: 送らない → Δ = 41→2, 2→3 (4), 1→3 (10)最短 1→2→3(費用 2)に 4 送る → Δ = 21→2, 2→3 は残余 0余剰 1 < 2: 送らない → Δ = 11→3 (10) が残る1→3(費用 3)に 1 送る 総費用 = 4×2 + 1×3 = 11。段階の初めに「残余 ≥Δ・被約費用が負」の辺を飽和させるので、Δ 残余グラフでは Dijkstra が使える 実測(N=300, M=3000, 容量・費用 10^9): 容量スケーリング 1.01〜1.07 s(Dijkstra 3214〜3292 回)                     SSP(ボトルネック量)0.86〜1.07 s(増加 990〜1086 回) ランダムなグラフでは速くならない。価値は「増加回数が容量の大きさに比例しない」保証

ヒント(クリックで段階的に開示)

ヒント1(方向性)

SSP(最短路に沿って流せるだけ流す)は正しいが、増加 1 回で流れる量に下限が無い。「Δ 単位ずつしか流さない」と決め、Δ を $2^{30}$ から 1 まで半分にしていく。各段階では「残余容量 Δ 以上の辺」だけのグラフで最短路を探す。

ヒント2(アプローチ)

Δ 段階の初めに、残余容量 Δ 以上で被約費用($w+p_u-p_v$)が負の辺を全部飽和させる。これで Δ 残余グラフの被約費用がすべて非負になり、Dijkstra が使える。負の費用の辺・負の閉路も、この飽和で各段階ごとに片づく(Day170-Q4 の「最初に一度だけ飽和」を、Δ ごとに繰り返す形)。

ヒント3(誘導)

余剰 $\ge\Delta$ の頂点すべてを始点に多始点 Dijkstra、最初に取り出された不足 $\le-\Delta$ の頂点で打ち切る。ポテンシャルは $p_v\mathrel{+}=\min(d_v,d_t)$(届かなかった頂点も $d_t$ で打ち切る)。経路に Δ 流す。送れなくなったら Δ を半分に。最後に余剰が残っていれば infeasible。

模範解答 (Python)

模範解答を表示する(まず自力で解いてから)
import sys
from heapq import heappush, heappop


def solve(n, b, edges):
    # 残余グラフは配列で持つ: 辺 i と i^1 が互いに逆辺
    head = [[] for _ in range(n)]
    to, cap, cost = [], [], []
    for u, v, c, w in edges:
        head[u].append(len(to)); to.append(v); cap.append(c); cost.append(w)
        head[v].append(len(to)); to.append(u); cap.append(0); cost.append(-w)
    ex = b[:]                                   # 各頂点の余剰(正 = まだ送り出す量)
    pot = [0] * n                               # ポテンシャル p。被約費用 = cost + p[u] - p[v]
    U = max([c for _, _, c, _ in edges] + [abs(x) for x in b] + [1])
    delta = 1 << (U.bit_length() - 1)
    INF = float('inf')
    total = 0
    while delta >= 1:
        # (1) 残余 delta 以上で被約費用が負の辺を飽和させる(delta 残余グラフの被約費用を非負にする)
        for u in range(n):
            for e in head[u]:
                r = cap[e]
                if r >= delta and cost[e] + pot[u] - pot[to[e]] < 0:
                    cap[e] = 0
                    cap[e ^ 1] += r
                    ex[u] -= r
                    ex[to[e]] += r
                    total += r * cost[e]
        # (2) 余剰 >= delta の頂点から不足 <= -delta の頂点へ、delta ずつ最短路で送る
        while True:
            srcs = [v for v in range(n) if ex[v] >= delta]
            if not srcs or not any(ex[v] <= -delta for v in range(n)):
                break
            dist = [INF] * n
            prev = [-1] * n
            hq = []
            for s in srcs:
                dist[s] = 0
                hq.append((0, s))
            t = -1
            while hq:
                d, u = heappop(hq)
                if d > dist[u]:
                    continue
                if ex[u] <= -delta:             # 最初に取り出された不足頂点で打ち切る
                    t = u
                    break
                pu = pot[u]
                for e in head[u]:
                    if cap[e] >= delta:
                        v = to[e]
                        nd = d + cost[e] + pu - pot[v]
                        if nd < dist[v]:
                            dist[v] = nd
                            prev[v] = e
                            heappush(hq, (nd, v))
            if t < 0:
                break                           # この delta ではもう送れない
            dt = dist[t]
            for v in range(n):                  # p += min(d, d_t): 届かなかった頂点も d_t で打ち切る
                pot[v] += dist[v] if dist[v] < dt else dt
            v = t
            while prev[v] >= 0:
                e = prev[v]
                cap[e] -= delta; cap[e ^ 1] += delta
                total += delta * cost[e]
                v = to[e ^ 1]
            ex[v] -= delta; ex[t] += delta
        delta >>= 1
    if any(ex):
        return None
    return total


def main():
    data = sys.stdin.buffer.read().split()
    n, m = int(data[0]), int(data[1])
    b = [int(x) for x in data[2:2 + n]]
    edges = []
    k = 2 + n
    for _ in range(m):
        u, v, c, w = (int(x) for x in data[k:k + 4])
        k += 4
        edges.append((u - 1, v - 1, c, w))
    ans = solve(n, b, edges)
    print("infeasible" if ans is None else ans)


if __name__ == "__main__":
    main()

計算量: 段階は $\lceil\log_2U\rceil+1$ 個($U\le10^{12}$ で 40 程度)。各段階の増加は $O(N+M)$ 回で抑えられる、と教科書の解析は述べている(記憶。未確認)。増加 1 回が Dijkstra $O(M\log N)$。$N=300$、$M=3000$、容量・費用 $10^9$ のランダムなグラフで 1.01〜1.07 s、Dijkstra 3214〜3292 回(Python 3.14.7、M1 MacBook Pro、実測)。

Step-by-Step 解説

1Δ 残余グラフ

残余容量が Δ 以上の辺だけを残したグラフを「Δ 残余グラフ」と呼ぶ。この中の経路なら必ず Δ 流せる。Δ 段階では Δ ずつしか流さないので、1 回の増加で「少なくとも Δ」進む。

2段階の初めの飽和

Δ を半分にすると、これまで使えなかった「残余 $\ge\Delta$ になった辺」が Δ 残余グラフに加わる。その中には被約費用が負の辺があり得る(前の段階で Dijkstra の対象外だったので、ポテンシャルがその辺に合っていない)。そこで、被約費用が負で残余 $\ge\Delta$ の辺は残余全部を流して飽和させる。飽和した辺は Δ 残余グラフから消え、逆辺(被約費用は正)だけが残る。供給の帳尻は余剰 $ex$ に吸収させる。負の閉路も、閉路のどこかの辺がこの飽和で消える。

3多始点 Dijkstra と打ち切り

余剰 $\ge\Delta$ の頂点すべてを距離 0 で入れる。不足 $\le-\Delta$ の頂点が最初に取り出されたらそこで止める。止めた時点で未確定の頂点 $v$ の真の距離は $d_t$ 以上なので、$p_v\mathrel{+}=\min(d_v,d_t)$ としても Δ 残余グラフの被約費用は非負のまま(辺 $(u,v)$ で、$d_u\le d_t$ なら $w_p+\min(d_u,d_t)-\min(d_v,d_t)\ge w_p+d_u-d_v\ge0$、$d_u>d_t$ なら $w_p+d_t-\min(d_v,d_t)\ge0$)。

4段階の終わりと判定

余剰 $\ge\Delta$ の頂点か不足 $\le-\Delta$ の頂点のどちらかが無くなるか、到達できなくなったら Δ を半分に。Δ=1 の段階が終わっても余剰が 0 でない頂点があれば infeasible。

5実測での位置づけ

同じ入力で、ボトルネック量を流す SSP(Day170-Q4 の方法を Dijkstra にしたもの)と比べた。$N=300$、$M=3000$ の 3 ケースで、容量スケーリングは Dijkstra 3214〜3292 回・1.01〜1.07 s、SSP は増加 990〜1086 回・0.86〜1.07 s。ランダムなグラフでは SSP の方が増加回数が約 1/3 で、時間は同程度だった。容量スケーリングを選ぶ理由は「増加回数が容量に依存しない保証」であって、平均の速さではない。

迷った他案と採らなかった理由

案計算量採否理由
流量の全列挙$\prod(c_i+1)$検算用のみ$N\le4$、$M\le5$、$c\le3$ の 400 ケース
負の辺を先に飽和+SSP(Bellman–Ford、ボトルネック量)増加回数に容量依存の保証なし検算用のみ中規模 700 ケースの参照解
負の辺を先に飽和+SSP(Dijkstra、ボトルネック量)同上比較用上の実測で同程度の速さ。理論上は増加回数が指数的になる例がある(Zadeh 1973 と記憶、未確認)
単位量ずつ流す SSP増加 $\sumb/2$ 回不採用容量 $10^9$ では $10^9$ 回以上(導出)
容量スケーリング(本解)$O(M\log U\cdot M\log N)$採用1.01〜1.07 s、容量に依存しない保証
費用スケーリング・ネットワーク単体法実用で速いとされる本問では扱わない実装が重い。LEMON の比較(記憶)では速い側

前提が変わったら:

  • 最悪ケースの保証が要らず平均で速ければよいなら、SSP(Dijkstra+ボトルネック量)で十分。本セッションのランダムなグラフでは同程度だった
  • $N,M$ が $10^5$ 規模なら Python では重い。C++ で費用スケーリングかネットワーク単体法(LEMON や OR-Tools の実装)を使う
  • 費用がすべて非負で負の閉路が無いなら、段階の初めの飽和は Δ を半分にしたときに現れる辺だけが対象になり、SSP からの書き換えがより小さい

よくあるミス

ミス原因正しい書き方
段階の初めの飽和を省く「最初に負の辺を飽和させれば十分」と思うΔ を半分にすると、ポテンシャルが合っていない辺が Δ 残余グラフに入る。全列挙 400 ケース中 52、SSP との照合 400 ケース中 272 で誤答し、そのうち 9・189 ケースは Dijkstra か経路の復元が終わらなかった(負の閉路を回り続ける。打ち切りを入れて計測)
届かなかった頂点のポテンシャルを更新しない$d_v=\infty$ を足せない$p_v\mathrel{+}=\min(d_v,d_t)$。$d_t$ で打ち切る
Dijkstra を最後まで回す全頂点の距離が要ると思う最初の不足頂点で止めてよい(上の min の扱いで正しさは保たれる)
Δ の初期値を容量だけから決める供給量 $b$ が容量より大きい場合を忘れる$U=\max(c_i,b_v)$。足りなければ infeasible になるだけだが、Δ が小さすぎると段階数が足りず遅くなる
総費用を浮動小数で持つ容量 $10^9$ × 費用 $10^9$ × 辺 3000Python の整数で持つ(出力は最大 $10^{21}$ 程度、計測用の最大ケースで $-3.05\times10^{20}$)

次のステップ

発展1: 費用スケーリング($\varepsilon$-最適性と push/relabel)を実装し、同じ入力で比べる

  • 発展2: SSP の増加回数が容量に比例して増える入力を構成する(本セッションでは Ford–Fulkerson の悪い例を費用 0 で試したが、SSP は 2 回の増加で終わり、作れなかった)
  • 発展3: 原典(Edmonds–Karp 1972、Network Flows の capacity scaling の章、Kovács 2015)を確認する(本セッションでは未確認)

検証ログ(本セッションで実行した数値)

全列挙との照合(実測): $N\le4$、$M\le5$、$c\le3$、$|w|\le5$ の 400 ケース(負の費用を含む 249、infeasible 37)で不一致 0

  • SSP(Bellman–Ford)との照合(実測): $N\le12$、$M\le30$、容量 5〜$10^9$、費用 $\pm10^9$ の 400 ケース(負の費用を含む 373、infeasible 54)と、別の 300 ケースで不一致 0
  • 誤った方法(実測): 段階の初めの飽和を省くと、全列挙 400 中 52・SSP 400 中 272 で誤答(うち停止しないもの 9・189)
  • 入力例(実測): 例 1 11、例 2 -2、例 3 infeasible(全列挙と一致)
  • 実行時間(実測、負荷平均 2〜8): $N=300$、$M=3000$、容量・費用 $10^9$(負の費用あり)の 3 ケースで 1.01〜1.07 s(Dijkstra 3214〜3292 回)、同じ入力の SSP は 0.86〜1.07 s(増加 990〜1086 回)。$N=200$、$M=2000$ では 0.42〜0.64 s(SSP 0.42〜0.75 s)
  • 扱っていない範囲: 費用スケーリング・ネットワーク単体法、SSP が遅くなる入力の構成、$10^5$ 規模、原典の確認

自己評価

自分の回答

気づき・メモ