この問題は LP 双対 、特に 費用流双対 と呼ばれるテーマを扱った問題です。

予備知識

費用流双対について簡単に解説します。

• 最小費用流問題が高速に計算できること、およびその簡単な仕組み(Bellman-Ford アルゴリズムを用いた最短路反復法)は前提知識としているので、知らない方は hos さんによる資料 などをあらかじめ読んでおくと理解が進みやすいと思います。

線形計画問題とその双対問題

線形計画問題 (LP) とは、目的関数(最大化/最小化の対象となる関数)が線形関数で、かつ変数の制約条件も線形関数の等式と不等式で記述できる問題のことを言います。
具体的に行列とベクトルを用いて表すと、以下の形式で表せる問題を言います(表現形にはいくつかバリエーションがあります)。ここで \(A\) は定数行列、\(b,c\) は定数ベクトル、\(x\) は変数ベクトルです。

\[ \begin{aligned} &\text{maximize} & &c^{T} x \\ &\text{subject to} & & Ax \leq b\\ & & & x \geq 0 \end{aligned} \]

上の問題の 双対問題 とは、以下の問題のことを言います。

\[ \begin{aligned} &\text{minimize} & &b^{T} y \\ &\text{subject to} & & A^{T}y \geq c\\ & & & y \geq 0 \end{aligned} \]

双対問題に対応する語として、元の問題を 主問題 と呼びます。

主問題と双対問題の間には次の 強双対定理 が成り立ちます。(このような性質を 強双対性 と呼びます。)

強双対定理

主問題と双対問題の両方に制約を満たす解 (実行可能解 と呼ぶ) が存在するとする。この時、両方の問題には最適解 \(x^\ast, y^\ast\) が存在して、両者の最適値は一致する(つまり、\(c^T x^\ast = b^T y^\ast\) が成り立つ)。

双対性の詳しい説明や LP のバリエーションに対する双対の取り方などについては、以下の参考資料を参考にしてください。

• wata さんによる資料
• tokoharu さんによる資料
• adamant さんの記事

最小費用流問題の双対問題

最小費用流問題は線形計画問題として記述することが出来ます。
まず定式化を行います。ここでは通常の始点と終点が存在するフロー (s-t flow) ではなく、一般化した b-flow と呼ばれる形で定式化します。すなわち次の問題です。

各辺 \(e\) に容量 \(c_e\) \((\geq 0)\) がある有向グラフ \((V, E)\) がある。
各頂点 \(v\) には \(b_v\) の供給がある(\(b_v < 0\) の場合は \(-b_v\) の需要がある)。
辺 \(e\) にフローを \(1\) 流すコストは \(w_e\) である。
辺 \(e\) に流れるフローの流量を \(f_e\) とする。\(\sum_e f_e w_e\) を最小化せよ。

この問題を線形計画問題として定式化すると次のようになります。ここで \(u\to v\) 辺を \((u,v)\) のように表します。

\[ \begin{aligned} &\min. & & \sum_{e \in E} w_e f_e \\ &\mathrm{s.t.} & & \sum_{e=(v,u)\in E} f_e - \sum_{e=(u,v)\in E} f_e = b_v & & (\forall v \in V) \\ & & & f_e \leq c_e & & (\forall e \in E) \\ & & & f_e \geq 0 & & (\forall e \in E) \end{aligned} \]

LP 双対を取ってみます。上の式から \(p_v, z_e\) を割り当てます。

\[ \begin{aligned} &\max. & & \sum_{v \in v} b_v p_v + \sum_{e \in E} c_e z_e \\ &\mathrm{s.t.} & & p_u - p_v + z_e \leq w_e & & (\forall e=(u,v) \in E)\\ & & & z_e \leq 0 & & (\forall e \in E) \end{aligned} \]

\(p_u, p_v, z_e\) を \(-1\) 倍したものに置き換えます。

\[ \begin{aligned} &\max. & & -\sum_{v \in v} b_v p_v - \sum_{e \in E} c_e z_e \\ &\mathrm{s.t.} & & p_v - p_u - z_e \leq w_e & & (\forall e=(u,v) \in E)\\ & & & z_e \geq 0 & & (\forall e \in E) \end{aligned} \]

式に注目すると、\(p_u, p_v\) を固定した時 \(z_e\) は \(\max(p_v - p_u - w_e, 0)\) を取る場合が最適です。よって \(z_e\) を消去した以下の問題に置き換えられます。

\[ \begin{aligned} &\max. & & -\sum_{v \in v} b_v p_v - \sum_{e =(u,v)\in E} c_e \max(p_v - p_u - w_e, 0) \\ \end{aligned} \]

目的関数を \(-1\) 倍して最小化問題とします。

\[ \begin{aligned} &\min. & & \sum_{v \in v} b_v p_v + \sum_{e =(u,v)\in E} c_e \max(p_v - p_u - w_e, 0) \\ \end{aligned} \]

この式が 最小費用流問題の双対問題(費用流双対) として知られている問題の定式化です。(定式化は少し異なりますが学術的には 最小費用テンション問題 と呼ばれています。) また、\(p_v\) のことを ポテンシャル と呼びます。

費用流双対の LP がどのような問題か解釈するとおおよそ次の問題になります。

各頂点 \(v\) に対して値 \(p_v\) を割り当てる。
各頂点 \(v\) に対して \(b_v p_v\) のコストがかかる。
各辺 \(e=(u,v)\) に対して \(c_e \max(0, p_v - p_u - w_e)\) の和のコストがかかる。
コストを最小化せよ。

ところで、任意の \(x\) に対して凸な折れ線(区分線形凸関数 と呼びます)は \(\max(0, x - w)\) と \(\max(0, w - x)\) の線形結合として記述することが出来ます。
よって以下の問題は費用流双対として LP で記述できて、目的関数の最適値を最小費用流を利用して計算できます。

各頂点 \(v\) に対して値 \(p_v\) を割り当てる。
各頂点 \(v\) に対して \(b_v p_v\) のコストがかかる。
各辺 \(e=(u,v)\) に対して \(F_e(x)\) を区分線形凸関数とした時 \(F_e(p_v-p_u)\) のコストがかかる。
コストを最小化せよ。

ここまでの話をまとめると次の事実がわかります。

• \(p_v-p_u\) の差分に応じて凸なコストがかかる時にコストを最小化する問題は費用流双対と呼ばれている。
• 費用流双対は LP 双対を取ると最小費用流問題になる。
• 主問題と双対問題の最適値は一致するので、費用流双対の最適値は最小費用流のアルゴリズムを用いて計算できる。

ここで次の問題が浮かびます。

• 費用流双対の最適解となるポテンシャルを復元するにはどうすればよいか？

この問題は容易ではないです。なぜならば、最小費用流のアルゴリズムが計算するのは費用流問題の最適値および最適解となる流量 \(f_e\) であってポテンシャル \(p_v\) ではないからです。(なお、実は最小費用流を計算するアルゴリズムの内部では \(p_v\) も同時に計算されていることが多いです。詳細は後述)

この問題は 相補性条件 と呼ばれる条件を用いることで解決できます。

相補性条件

この項は Taihei Oki 氏による記事 を参考に記述しました。

LP の相補性条件と呼ばれる条件について説明します。まず、LP の主問題と双対問題を再掲します。

• 主問題

\[ \begin{aligned} &\text{maximize} & &c^{T} x \\ &\text{subject to} & & Ax \leq b\\ & & & x \geq 0 \end{aligned} \]

• 双対問題

\[ \begin{aligned} &\text{minimize} & &b^{T} y \\ &\text{subject to} & & A^{T}y \geq c\\ & & & y \geq 0 \end{aligned} \]

ここで、

• \(c \leq A^T y\) かつ \(x \geq 0\) より \(c^{T} x \leq (A^T y)^T x = y^T A x\)
• \(Ax \leq b\) かつ \(y \geq 0\) より \(y^T A x \leq y^T b = b^T y\)

が成り立つので、

\[c^T x \leq y^T A x \leq b^T y\]

、すなわち (主問題の目的関数の値) \(\leq\) (双対問題の目的関数の値) が成り立ちます。これを 弱双対定理 や 弱双対性 と呼びます。
ところで、先に弱双対定理よりも強い定理である強双対定理を説明しました。すなわち最適解 \(x, y\) において \(c^T x = b^T y\) が成り立つという定理です。よって \(x, y\) が最適解の時に

\[c^T x = y^T A x = b^T y\]

が成り立つことになります。

以降では \(x\) を長さ \(n\) のベクトル、\(y\) を長さ \(m\) のベクトルとします。 \(c^T x = y^T A x\) を変形すると次式になります。

\[\sum_{i=1}^n((A^T y)_i - c_i) x_i = 0\]

ここで \(A^T y \geq c\) かつ \(x \geq 0\) より明らかに左辺は非負です。よって

• \(i=1,2,\dots,n\) に対して \(x_i = 0\) または \((A^T y)_i = c_i\)

という条件が成り立ちます。同様にして \(y^T A x = b^T y\) を変形すると

• \(j=1,2,\dots,m\) に対して \(y_j = 0\) または \((Ax)_j = b_j\)

が成り立ちます。これを 相補性条件 と呼びます。
また、解 \(x, y\) が相補性条件を満たしているならば \(c^T x = y^T A x = b^T y\) が成り立ち \(x,y\) が最適解となることも明らかです。よって \(x,y\) が最適解であることと相補性条件は同値です。これを 相補性定理 と呼びます。

相補性定理を利用して主問題の最適解集合を求めることを考えます。今、双対問題の最適解 \(y\) が 1 つわかっているとします。この時、相補性定理から、最適解 \(x\) が満たす必要十分条件は次のように記述できます。

• \(i=1,2,\dots,n\) に対して \((A^T y)_i \gt c_i\) ならば \(x_i = 0\)
• \(j=1,2,\dots,m\) に対して \(y_j \gt 0\) ならば \((Ax)_j = b_j\)

すなわち、\(x\) の最適解が満たす条件は不等式条件のいくつかを等式に変えたものになります。より詳しく言うと、双対問題においてある不等号条件が成立すれば、対応する主問題の条件が等号条件になります。このようにして、双対問題の最適解から主問題の解が満たす条件を導出することが出来ます。

費用流双対の解の復元

さて、相補性条件を利用して費用流双対の解となるポテンシャルを復元しましょう。
費用流双対を主問題とした時の LP を再掲します。(細かい -1 倍の違いについては本質的でないので無視します)

• 主問題

\[ \begin{aligned} &\min. & & \sum_{v \in v} b_v p_v + \sum_{e \in E} c_e z_e \\ &\mathrm{s.t.} & & p_v - p_u - z_e \leq w_e & & (\forall e=(u,v) \in E)\\ & & & z_e \geq 0 & & (\forall e \in E) \end{aligned} \]

• 双対問題

\[ \begin{aligned} &\min. & & \sum_{e \in E} w_e f_e \\ &\mathrm{s.t.} & & \sum_{e=(v,u)\in E} f_e - \sum_{e=(u,v)\in E} f_e = b_v & & (\forall v \in V) \\ & & & f_e \leq c_e & & (\forall e \in E) \\ & & & f_e \geq 0 & & (\forall e \in E) \end{aligned} \]

相補性条件を適用すると、\(p_v, z_e\) に関して以下の条件を得ます。

• \(f_e \lt c_e\) である \(e\) に対して \(z_e = 0\)
• \(f_e \gt 0\) である \(e=(u,v)\) に対して \(p_v - p_u - z_e = w_e\)

式変形により \(z_e\) を条件から消して整理すると次のようになります。

• \(f_e \lt c_e\) である \(e=(u,v)\) に対して \(p_v - p_u \leq w_e\)
• \(f_e \gt 0\) である \(e=(u,v)\) に対して \(p_u - p_v \leq -w_e\)

この条件を満たすポテンシャル \(p_v\) を求められればよいです。このような \(p_v\) を 1 つ求める問題は、競技プログラミングの世界では 牛ゲー(最短路問題の双対問題) としてよく知られています。すなわち、\(p_v\) は以下のアルゴリズムで計算できます。

• \(p_a - p_b \leq w_c\) という条件に対して、\(b \to a\) に重み \(w_c\) の辺を貼ったグラフを考える。
• 適当な点を始点として Bellman-Ford アルゴリズムにより各頂点の始点からの距離を計算する。
• グラフに負閉路が存在しなければ 解となる \(p_v\) を得られる。

ただしここで「グラフに負閉路が存在しなければ」という点が問題となります。上記の手順で構築したグラフに負閉路は存在しないのでしょうか？

ここで、構築したグラフは次のようなグラフです。

• \(0 \leq f_e \lt c_e\) であるとき \(u \to v\) に重み \(w_e\) の辺
• \(0 \lt f_e \leq c_e\) であるとき \(v \to u\) に重み \(-w_e\) の辺

このグラフはよく見ると最小費用流のアルゴリズム内部で構築する 残余グラフ で貼られている辺と完全に一致しています！
ここで最小費用流の性質として、費用流 \(f_e\) が最適な費用流であることと、残余グラフに負閉路が存在しないことは同値であることが知られています。よって、\(f_e\) が最適解であれば残余グラフに負閉路は存在せず、Bellman-Ford アルゴリズムによってポテンシャルを復元するアルゴリズムが正当であることが証明できました。

再度ここまでの話をまとめると次のようになります。

• \(p_v-p_u\) の差分に応じて凸なコストがかかる時にコストを最小化する問題は費用流双対と呼ばれている。
• 費用流双対は LP 双対を取ると最小費用流問題になる。
• 主問題と双対問題の最適値は一致するので、費用流双対の最適値は最小費用流のアルゴリズムを用いて計算できる。
• 費用流双対の解 \(p_v\) を求めたい場合は、相補性定理により最小費用流の解 \(f_e\) から Bellman-Ford アルゴリズムを用いて復元できる。

ここまでの説明では「最小費用流を計算するアルゴリズム」の中身に深く立ち入らずに説明していました。実際 AtCoder 環境では最小費用流はアルゴリズムの内容を知らずとも AtCoder Library(ACL) を利用すれば計算できるものです。
しかし実は、さきほど少し話した通り、最小費用流を計算するアルゴリズムによっては、最小費用流の解 \(f_e\) だけでなく費用流双対の解 \(p_v\) も同時に計算することができます。そのような最小費用流のアルゴリズムはいくつかありますが、ここでは ACL で採用されている Primal-Dual 法 について説明します。

Primal-Dual 法

以降では簡単のため、負辺(コストが負の辺)が存在しない s-t flow の最小費用流を求めるアルゴリズムを説明します。

• 負辺が存在する場合は、例えば snuke さんの記事 に書かれている手法で対処できます。
• b-flow を求めたい場合は超頂点 \(s,t\) を用意して適切に辺を貼ることで s-t flow に帰着できます。

Primal-Dual 法は 最短路反復法 (逐次最短路法とも。英語では Successive Shortest Path algorithm, SSP) の一種です。

• 資料によっては Primal-Dual 法は最短路反復法とは異なるアルゴリズムとして紹介されている点に注意してください。ここでは Dijkstra 法を用いた最短路反復法を Primal-Dual 法として定義します。

まずはじめに、 Bellman-Ford 法を用いた最短路反復法を簡単に振り返ると次のアルゴリズムになります。

• はじめ、全ての辺の流量 \(f_e\) を \(0\) とする。
• フローの流量が \(F\) 未満の間、以下の操作を繰り返す：
  • 残余グラフ上で \(s \to t\) 間の(費用に関する)最短路を求める。
  • 発見した最短路の 1 つに従ってフローを流す。

このアルゴリズムでは (初めの 1 回以外は) 残余グラフに負辺が発生することから、最短路を求めるのに Bellman-Ford アルゴリズムが必要になります。よって全体の計算量は頂点数を \(N\)、辺数を \(M\) として \(\mathrm{O}(F N M)\) になります。

最短路を求める部分の計算量を改善するためにポテンシャル \(p_v\) を利用します。 先に説明した相補性定理から、残余グラフに対して以下の性質が成り立ちます。

• 以下の条件を満たすポテンシャル \(p_v\) が存在する。
  • 残余グラフ上にコスト \(w_{uv}\) の辺 \(u\to v\) が存在する時、\(w_{uv}+p_u-p_v \geq 0\) が成り立つ。

今、残余グラフのポテンシャル \(p_v\) がわかっているとします。最短路反復法において求めたいものは頂点 \(s\) から頂点 \(t\) への最短路です。
ここで、残余グラフの辺のコストを \(w_{uv}\) から \(w_{uv}+p_u-p_v\) (これを 簡約費用 と呼びます) に置き換えたグラフを考えます。

この時、パス \(v_0 \to v_1 \to v_2 \to \dots \to v_k\) の簡約費用は

\[ \sum_{i=0}^{k-1} w_{v_i v_{i+1}} + p_{v_i} - p_{v_{i+1}} = \left(\sum_{i=0}^{k-1} w_{v_i v_{i+1}}\right) + p_{v_0} - p_{v_k} \]

となり、そのパス自体の長さおよび始点と終点のポテンシャルにのみ依存します。よって、\(s-t\) 最短路は簡約費用を用いて計算しても正しい最短路が得られることが分かります。
ここで簡約費用の非負性より、s-t 最短路の計算に簡約費用におけるグラフ上の最短路を計算することにすると Bellman-Ford 法ではなく Dijkstra 法を利用できる、という点がポイントです。

なお、残余グラフに存在する辺はフローを流すごとに変化するため、ポテンシャル \(p_v\) もフローを流すごとに更新していく必要があります。更新の方法は、実はフローを流すごとにポテンシャル \(p_v\) に簡約費用における \(s-v\) 最短距離を加算すればよいです。(証明は 岡本吉央氏の資料 に載っているのでそちらに譲ります。)

以上の内容をまとめると、次のアルゴリズムが成立します。

• はじめ、全ての辺の流量 \(f_e\) および頂点のポテンシャル \(p_v\) を \(0\) とする。
• フローの流量が \(F\) 未満の間、以下の操作を繰り返す：
  • 残余グラフ上で \(s \to t\) 間の(簡約費用に関する)最短路を求める。
  • 発見した最短路の 1 つに従ってフローを流す。(\(f_e\) の更新)
  • 各頂点のポテンシャルに簡約費用における最短距離を加算する。(\(p_v\) の更新)

このアルゴリズムでは、「主問題の解を改善する」\(\to\) 「双対問題の解を改善する」\(\to\) 「主問題の解を改善する」\(\to\) \(\dots\) というように主問題 (Primal) と双対問題 (Dual) の解の改善を交互に行っています。この様子からこのアルゴリズムは Primal-Dual 法 と呼ばれています。
計算量を考えると、最短路の計算に Dijkstra 法を用いてよいことから \(\mathrm{O}(F(N+M)\log(N+M))\) となり、Bellman-Ford 法のみを用いた最短路反復法よりも高速に動作します。

ここから先は余談ですが、Primal-Dual 法は計算量に \(F\) が登場する関係上、辺容量が大きいグラフでは高速に動作する計算量保証ができません。そこで Primal-Dual 法をさらに工夫した 容量スケーリング と呼ばれる手法を用いることで、Dijkstra 法を行う回数を \(\mathrm{O}(F)\) 回から \(\mathrm{O}(M \log F)\) 回に減らすことが出来ます。詳細は みさわさんの記事 に詳しいです。
また、ライブラリの整備が充実している競技プログラマのコードを読むと 費用スケーリング あるいは ネットワーク単体法 と呼ばれるアルゴリズムが用いられていて、実際この 2 種類のアルゴリズムが実用上抜きんでて高速なようです。(arxiv の参考文献) 最小費用流のアルゴリズムは奥が深いので色々調べてみると面白そうです。

さて、ここまでの内容をまとめると次の通りです。

• \(p_v-p_u\) の差分に応じて凸なコストがかかる時にコストを最小化する問題は費用流双対と呼ばれている。
• 費用流双対は LP 双対を取ると最小費用流問題になる。
• 主問題と双対問題の最適値は一致するので、費用流双対の最適値は最小費用流のアルゴリズムを用いて計算できる。
• 費用流双対の解 \(p_v\) を求めたい場合は以下の方法が考えられる。
  • 相補性定理により最小費用流の解 \(f_e\) から Bellman-Ford アルゴリズムを用いて復元する。
  • 最小費用流の計算に Primal-Dual 法などのアルゴリズムを利用することで \(f_e\) と \(p_v\) を同時に計算する。

以上の内容を踏まえたうえで、今回の問題の解説に入ります。

解法

解法を説明します。概要を説明すると、今回の問題は費用流双対の類型として定式化できるので、最小費用流問題に帰着してポテンシャルを復元する方針で解を求められます。

\(K = \frac{P}{Q}\) とします。グリッドを \(N^2\) 頂点 \(2N(N-1)\) 辺のグラフ \((V, E)\) とみなします。操作後に頂点 \(v\) に書かれている数を \(p_v\) とすると、答えは次の LP の最適値およびその時の最適解 \(p_v\) となります。

\[ \begin{aligned} &\min. & &\sum_{e=(u,v)\in E} \vert p_v - p_u \vert\\ &\mathrm{s.t.}& &\sum_{v \in V} \vert p_v - A_v\vert \leq K \end{aligned} \]

辺を無向辺ではなく有向辺が両方の向きに貼られたものとみなします。(すなわち、辺を有向にする代わりに倍加します。)そして式を変形すると次のようになります。ここで \(u\) から \(v\) へ向かう辺を \((u,v)\) と表します。

\[ \begin{aligned} &\min. & &\sum_{e=(u,v)\in E} \max(p_v - p_u, 0) \\ &\mathrm{s.t.}& &\sum_{v \in V} \left(\max(p_v - A_v,0)+\max(A_v - p_v, 0)\right)\leq K \end{aligned} \]

上式の \(\max\) の部分をそれぞれ \(a_e, x_v, y_v\) と置いて式を整理します。

\[ \begin{aligned} &\min. & & \sum_e a_e \\ &\mathrm{s.t.}& & p_v - p_u - a_e \leq 0& & (\forall e=(u,v) \in E)\\ & & & p_v - x_v \leq A_v & & (\forall v \in V)\\ & & & -p_v - y_v \leq -A_v& & (\forall v \in V)\\ && &\sum_v \left( x_v + y_v \right)\leq K\\ & & & a_e, x_v, y_v \geq 0 \\ \end{aligned} \]

双対を取ります。上の式から順に \(f_e, z_v, w_v,\lambda\) を割り当てると次のようになります。

\[ \begin{aligned} &\max. & & \sum_{v \in V} A_v(z_v - w_v) + K \lambda \\ &\mathrm{s.t.}& & \sum_{e=(u,v)\in E}f_e - \sum_{e=(v,u)\in E}f_e + z_v - w_v = 0 & & (\forall v \in V)\\ & & &-f_e \leq 1 & & (\forall e \in E)\\ & & &-z_v + \lambda \leq 0 & & (\forall v \in V)\\ & & &-w_v + \lambda \leq 0 & & (\forall v \in V)\\ & & &f_e, z_v, w_v, \lambda \leq 0 \end{aligned} \]

\(-1\) 倍して最小化問題にします。さらに \(f_e,z_v,w_v,\lambda\) も \(-1\) 倍したものに置き換えます。

\[ \begin{aligned} &\min. & & \sum_v A_v(z_v - w_v) + K \lambda\\ &\mathrm{s.t.}& & \sum_{e=(u,v)\in E}f_e - \sum_{e=(v,u)\in E}f_e + z_v - w_v = 0 & & (\forall v \in V)\\ & & &0 \leq f_e \leq 1 & &(\forall e \in E)\\ & & &0 \leq z_v \leq \lambda & & (\forall v \in V)\\ & & &0 \leq w_v \leq \lambda & & (\forall v \in V)\\ \end{aligned} \]

よって最小費用流問題に帰着します。すなわち、\(f(\lambda)\) を以下の最小費用循環流のコストとすると、\(U = f(\lambda) + K \lambda\) を最小化する問題になります。

• 頂点 \(S\) および \(v \in V\) を用意する。
• 頂点 \(S\) から頂点 \(v\) に容量 \(\lambda\) 、コスト \(A_i\) の辺を貼る。
• 頂点 \(v\) から頂点 \(S\) に容量 \(\lambda\) 、コスト \(-A_i\) の辺を貼る。
• 頂点 \(u\) から頂点 \(v\) に容量 \(1\) 、コスト \(0\) の辺を貼る。
• 最小費用循環流を流した時の最小費用は？

ここで \(\lambda\) は実数を取る点に注意してください。すなわち、最適解における \(\lambda\) の整数性を仮定して三分探索を行う方針は不正解となります。しかし実は、 \(f(\lambda)\) の傾きが変化する点は有限個で、しかもそのような \(\lambda\) の分母は十分小さいことが言えます。以降ではこの事実を証明します。

ここからは、\(\lambda\) を固定した時の \(f(\lambda)\) およびその時の双対問題の解(最小費用流)と主問題の解(ポテンシャル)に注目します。

まず、\(f(\lambda)\) は凸関数になります。

(証明) \(\lambda=\lambda_1, \lambda_2\) と固定した時の双対 LP の最適解を \(x_1, x_2\) と置きます。また \(0 \leq \alpha \leq 1\) とします。
\(x_3 = \alpha x_1 + (1-\alpha) x_2\) と置きます。 この時、制約条件は全て線形関数であることから \(x_3\) は実行可能解であることが確認できます。 また、この時の \(\lambda\) は \(\alpha \lambda_1 + (1-\alpha) \lambda a_2\) です。さらに、双対 LP に \(x\) を代入した時の値を \(\mathrm{DLP}(x)\) のように表すと、

\[\mathrm{DLP}(x_3) = \alpha \mathrm{DLP}(x_1) + (1-\alpha) \mathrm{DLP}(x_2)\]

が成り立ちます。よって、

\[ \begin{aligned} f(\alpha \lambda_1 + (1-\alpha) \lambda a_2) &\leq \mathrm{DLP}(x_3)\\ &= \alpha \mathrm{DLP}(x_1) + (1-\alpha) \mathrm{DLP}(x_2) \\ &= \alpha f(\lambda_1) + (1-\alpha) f(\lambda_2) \end{aligned} \]

が成り立ちます。この式から \(f(\lambda)\) が凸であることが言えるので、示されました。(証明終わり)

さらに、より強い主張として、\(f(\lambda)\) は傾きが整数性を持つ区分線形凸関数であることが証明できます。　

(証明) \(\lambda\) を固定した時に \(\displaystyle \lim_{\Delta\to +0} \frac{f(\lambda+\Delta)-f(\lambda)}{\Delta}\) が整数になることを示せば、\(f(\lambda)\) は凸関数であることと合わせて命題が示せます。よってこれを示します。
最小費用流の残余グラフを考えます。\(\lambda\) を \(\lambda + \Delta\) (\(\Delta\) は残余グラフに存在する辺の容量よりも十分小さい) に変化させるということは、\(\lambda\) における残余グラフ上に以下の辺を追加することを意味します。

• \(S\) から \(v\) に容量 \(\Delta\) 、コスト \(A_i\) の辺
• \(v\) から \(S\) に容量 \(\Delta\) 、コスト \(-A_i\) の辺

この時に新たに発生する負閉路たちの性質を考えると、\(S \to u \to \dots \to v \to S\) という形をしていて、この時に \(S\to u\) 辺と \(v\to S\) 辺の少なくとも一方は容量 \(\Delta\) である必要があります。そして、追加した辺以外の辺の容量は \(\Delta\) に比べて十分大きいので、負閉路に流れるフローの流量は \(\Delta\) です。よって、残余グラフに負閉路が存在しなくなるまで新たに負閉路にフローを流すことを繰り返すことで発生するコストは \((整数) \times \Delta\) という形になり、\(\displaystyle \lim_{\Delta\to +0} \frac{f(\lambda+\Delta)-f(\lambda)}{\Delta}\) が整数となることが示されました。(証明終わり)

さらに、費用流の形から \(f(\lambda)\) の傾きは \(-N^2 \max(A)\) 以上 \(0\) 以下であることが確認できます。よって \(f(\lambda)\) は有限個の線分からなる関数であることも言えます。

以上の議論により、\(\lambda\) の定義域 \([0,\infty)\) を \(f(\lambda)\) の傾きに応じていくつかの区間に分割することが出来ます。すなわち、 \(\lambda_0, \lambda_1, \dots, \lambda_n\) および \(K_0, K_1, \dots, K_n\) を以下の条件を満たすように取ることが出来ます。

• \((0 = \lambda_0, \lambda_1)\) における \(f(\lambda)\) の傾きは \(-K_0\) である。
• \((\lambda_1, \lambda_2)\) における \(f(\lambda)\) の傾きは \(-K_1\) である。
• \((\lambda_2, \lambda_3)\) における \(f(\lambda)\) の傾きは \(-K_2\) である。
• \(\vdots\)
• \((\lambda_n, \infty)\) における \(f(\lambda)\) の傾きは \(-K_n = 0\) である。

さて、元の LP は \(U = f(\lambda)+K \lambda\) を最小化する問題でした。(元の問題の \(U\) から符号が -1 倍されていますが同一視します。) \(f(\lambda)\) の概形から以下の事実が成り立ちます。

• \(K \gt K_0\) のとき \(\lambda = \lambda_0\) において最小値を取る。
• \(K = K_0\) のとき \(\lambda_0 \leq K \leq \lambda_1\) において最小値を取る。
• \(K_0 \gt K \gt K_1\) のとき \(\lambda = \lambda_1\) において最小値を取る。
• \(\vdots\)
• \(K_{n-1} \gt K \gt K_n\) のとき \(\lambda = \lambda_n\) において最小値を取る。
• \(K = K_n = 0\) のとき \(\lambda_n \leq \lambda\) において最小値を取る。

費用流を流したときに得られるポテンシャルに関して、主問題の解をポテンシャルと一致させるのに必要な操作のコスト (以降では単にコストと呼びます) を考えます。ある \(K\) について \(U\) を最小化する \(\lambda\) を 1 つ取った時の費用流のポテンシャルが主問題の解になっているという事実より、\(\lambda\) を固定した時の \(\lambda\) とポテンシャルのコスト \(C\) の関係は「\(K=C\) とした時において \(U\) は \(\lambda\) において最小値を取る」という関係になる必要があります。よって、ポテンシャルのコスト \(C\) について以下の事実が成り立ちます。

• \(\lambda = \lambda_0 = 0\) では \(C \geq K_0\) となる。
• \((\lambda_0, \lambda_1)\) では \(C=K_0\) となる。
• \(\lambda = \lambda_1\) では \(K_0 \geq C \geq K_1\) となる。
• \((\lambda_1, \lambda_2)\) では \(C=K_1\) となる。
• \(\vdots\)
• \((\lambda_n, \infty)\) では \(C=K_n=0\) となる。

ところで、費用流の辺重みには整数性があるので、残余グラフの最短路問題を解くことで得られるポテンシャルは常に整数値を取るため、この時の目的関数 \(U\) も整数値を取ります。よって\(K=K_i\) における目的関数 \(U\) の最小値 \(f(\lambda_i) + K_i \lambda_i\) は整数です。\((\lambda_i, \lambda_{i+1})\) における \(f(\lambda)\) は \(y=-K_i\lambda+f(\lambda_i) + K_i \lambda_i\) という式で表せる線分であることと合わせると、\(f(\lambda)\) に現れる線分は全て \(y=(整数)\lambda+(整数)\) という形の式をしていることが言えました。

よって、傾きの変化点である \(\lambda_i\) は \(y=-K_i\lambda+f(\lambda_i) + K_i \lambda_i\) と \(y=-K_{i+1}\lambda+f(\lambda_{i+1}) + K_{i+1} \lambda_{i+1}\) の交点であることから、この式を解くことで分母は \(K_i\) の差、すなわち \(N^2 \max(A)\) で抑えられることが証明できます。また、最小費用流の形から \(\lambda \geq 4\) において \(f(\lambda)\) の傾きが \(0\) になることが証明できるため、分子もまた \(4N^2 \max(A)\) で抑えられることが証明できます。以上の議論により、\(\lambda_i\) の分母・分子は十分小さいことが証明できました。

• なお、実験に基づく予想として、上手く議論すると分母・分子の大きさがより厳しい評価である \(\mathrm{O}(N^2)\) で抑えられることが示せると予想していますが、証明できていません…

ここまでの議論を元に解法を構成します。\(\lambda\) に対して Stern-Brocot 木上の二分探索を行うことを考えます。

• Stern-Brocot 木上の二分探索については ABC333-G で出題されているのでここでは詳しく説明しません。MMNMM さんによる ABC294-F のユーザ解説 なども参考になります。
• なお、今回の問題の制約では Stern-Brocot 木上を二分探索せずに一歩ずつ降りていく探索でもおそらく制約下において高速に動作して、この方が実装が簡便です。

二分探索の判定条件としては、\(\lambda\) を固定して費用流を流した時のポテンシャルのコストが \(K\) 未満かどうかを条件とします。

• 費用流の容量の整数性については、\(\lambda\) が有理数の場合でも容量全体を分母倍すれば整数流として処理できるので問題になりません。

分母を \(\mathrm{O}(N^2 \max(A))\) までで打ち切る二分探索を行うことで、

• \(\lambda_L\) : ポテンシャルのコストが \(K\) 以上である最大の分数
• \(\lambda_R\) : ポテンシャルのコストが \(K\) 未満である最小の分数

およびこの時のフローの費用 \(f(\lambda_L), f(\lambda_R)\)、ポテンシャル \(B(\lambda_L), B(\lambda_R)\) を得ることが出来ます。

• コーナーケースとして、ポテンシャルのコストが常に \(K\) 未満である場合があります。この場合は \(U=0,B=B(0)\) を出力すればよいです。

ここから \(U, B\) を求めることが出来ます。\(U\) については

\[\min(f(\lambda_L) + K \lambda_L, f(\lambda_R) + K \lambda_R)\]

から求まることは比較的明らかです。\(B\) については、実はコストが \(K\) となるように \(B(\lambda_L)\) と \(B(\lambda_R)\) の重み付き平均を取ればそれが答えとなります。(証明は単純ではないですが場合分けを丁寧に行えば証明できるので省略します。) 以上の内容を実装すればこの問題を解くことが出来ます。

計算量について考えます。費用流の流量は \(\mathrm{O}(\lambda N^2)\)、頂点数および辺数が \(\mathrm{O}(N^2)\)、 \(\lambda\) の分母分子が \(\mathrm{O}(N^2 \max(A))\) で分母倍して整数流にする点を踏まえると、ACL の mincostflow + Bellman-Ford を用いれば 1 回あたり \(\mathrm{O}(N^6 \max(A) \log N)\) で費用流およびポテンシャルを求められます。(負辺の消去を適切に行う必要があります。) これを二分探索で \(\log (N^2\max(A))\) 回行うので計算量は全体で \(\mathrm{O}(N^6 \max(A) \log N \log (N\max(A)))\) となり、十分高速に動作します。(使用するフローアルゴリズムや二分探索の方法によって \(N\) と \(\max(A)\) につく指数の和が \(4\) 乗～ \(10\) 乗になると思いますが、定数倍が十分軽いのでいずれも通ると考えられます。)
\(U\) の誤差についても簡単に評価します。\(U\) は有理数として正確な値を求めることにすれば、答えの出力に必要な誤差が発生する演算は整数から小数へのキャスト \(2\) 回と小数除算 \(1\) 回だけです。 よって \(\varepsilon\) をマシンイプシロン(double 型で \(2^{-53}\)、long double 型で \(2^{-65}\) ) とすると相対誤差は \(3\varepsilon\) で抑えられて、許容誤差に収まります。(とはいえ double 型だと際どいので C++ であれば long double 型の使用を勧めます。また、\(K\) をあらかじめ小数にして解く方針を取っている場合 \(U\) の計算方法によっては異符号の加算が発生してしまい桁落ち誤差が起きて WA の原因となります。ご注意ください。)
これは余談ですが、誤差が尋常でないほど厳しかったり \(K\) を有理数として入力で与えたりしているのは scipy.optimize.linprog や ortools.linear_solver などの LP ソルバーを用いた解法を誤差で落とすためです。(一般的な LP ソルバーは倍精度浮動小数点数型で動いているのでおそらく誤差で落ちると思います。__float128 のような四倍精度型に対応可能な LP ソルバーを持っていれば通せるかもしれません。)

練習問題

過去の ABC の問題で LP 双対や費用流双対で解ける問題を想定解・非想定解問わずまとめました。練習問題としてお使いください。

• LP 双対
  • ABC275 G
  • ABC356 G
• 費用流双対
  • ABC224 H
  • ABC275 Ex
  • ABC347 G