2017-06-07

最小費用流の負辺除去

Algorithm 競プロ

つい最近まで最小費用流の負辺除去、「なんか上手いことやれば出来ることもあるらしい」程度の認識だったんですが、ちゃんと考えてみたら自明やんってなったので書いておきます。
この記事を読めば多分、自明かどうかはともかくとして、かなり見通しがよくなるのでは思います。
（2019-06-27に若干書き換えました）

追記 (2021-08-15)：
ポテンシャルを求めるのが簡単な場合にはそれを使って負辺除去する方が楽です。
この記事の方法はグラフが複雑だったり、思考停止したいとき用かなぁ。
例題

最小費用流への認識

最小費用流は「始点 S から終点 T に水を流す」という問題に見えがちですが、それは少し本質とずれています。
水の例えを用いつつ言うならば「頂点間で水をやりとりして過不足を補い合う」という感じでしょうか。
実装の際に始点と終点があった方が分かりやすいことがあるだけで、定義上は必要なものではないんですね。

具体的に変数を使って書くとこんな感じです。（←最小費用流とかのLPの定式化とか見るとさらに分かりやすい）

有向グラフがある
最初、頂点には余っている水の量 $d_v$ が定まっている
- $d_v$ が負なら不足していることを表している
- 全体でみると過不足はない ( $\sum d_v = 0$ )
辺には容量 $cap_e$ とコスト $cost_e$ が定まっている
辺を使うと、からに水を移すことができる
- 移す量を $x$ とすると、 $0 \leq x \leq cap_e$ を満たさなければならず、コストは $x \cdot cost_e$ 掛かる
過不足を解消するための最小コストは？

最初から過不足がない場合は最小費用循環流問題ってやつになります。
負のコストの辺がある場合、わざわざサイクルに水を流してコストを減らすこともある点に注意。

負辺除去

本題です。
負辺がある場合は、最初にmaxまで流しておいて、コストが正の逆向きの辺があったとみなすだけで良いです。
負辺 $u \to v$ が来たときに具体的にやる操作を書くと、こんな感じです。

$d_u$ から $cap_e$ を引き、 $d_v$ に $cap_e$ を足す
答えに $cost_e \cdot cap_e$ を足しておく
容量 $cap_e$ コスト $-cost_e$ の辺 $v→u$ を追加する

かなり自然じゃない？

実装

蟻本方式の実装でやるなら、超頂点 $S', T'$ を用意して、

$d_v$ が正： $S'$ から $v$ へ容量 $d_v$ コスト 0 の辺を張る
$d_v$ が負： $v$ から $T'$ へ容量 $-d_v$ コスト 0 の辺を張る

で正の $d_v$ の和だけ流せば良いです。
実装をいじれば超頂点を作らなくてもできると思います。
あと、cost scalingってやつもおすすめらしいのでそのうち書くかも。（流量が大きくなりがちなこの手法と相性が良さそう）

おまけ

上の定式化そのままではないケースへの対処法あれこれ

SからTに流すんだけど流量が自由な場合
- T→Sに容量INFコスト0の辺を張るだけ
SからTに流すんだけど流せるだけ流さないといけない場合
- T→Sに容量INFコスト-INFの辺を張っておいて後で答えからこの辺のコストを引く
最小流量付きの辺がある場合
- 負辺除去と同じノリであらかじめ最小流量だけ流しておけば良い

補足

負辺除去には他の方法もあるので、それらを比較したwata先生のツイートを引用させていただきます。（ありがとうございます）

Bellman-Fordでポテンシャル初期化(負閉路無し) O(F m log n + nm)
適当に定数足して非負に (大きさFの最大重み二部マッチングなど) O(F m log n)
（右側の頂点たちにXのポテンシャルを付けているとも見なせる）
予めmaxまで流して負辺除去 O((F+F') m log n), F'は負の辺の容量の和

ちなみにこれらは３つとも蟻本に載っています。
シチュエーションに合わせて別の手法を選んだ方が良いこともあるかもしれません。
（特に２つ目は実装楽だし早いしよく出てくるし便利）

2017-06-03

ARC075 F 「Mirrored」速い解法

ARC075F

とりあえず桁数を1~18まで全部試します。
例えば6桁なら、

D
= fedcba - abcdef
= 99999*(f-a)
  +9990*(e-b)
  + 900*(d-c)

となるようなa~f（0≦a~f≦9, a≠0）の個数を数えれば良いです。
(f-a),(e-b),(d-c)みたいに端から決めていくんですが、
解説ではDとの差を縮めていくようなものを探索していますが、
下の位のmodに注目すると、実は(f-a),(e-b),(d-c)の組み合わせは一意に決めていけることが分かります。
例えば(f-a)を決める時は、99999*(f-a)%10 == D%10 とならないといけなくて、「0≦abs(f-a)≦9」と「999999*2>(9990+900)*9（つまり2ずれると収拾がつかなくなる）」に注意するとそういう(f-a)は一意に定まります。
(e-b)を決める時は (9990*abs(e-b)%100 == abs(D-99999*(f-a))%100) みたいな感じでどんどんやっていきます。符号はどっちかだけが正しいです。
計算量は O(log^2 D) とかです。

ソースコード

2017-05-16

ICPCのライブラリPDFの生成

ICPCのWFまたは多くのアジア地区予選では持ち込めるライブラリのページ数が決まっています。
また、PDFの右上にページ番号、左上に大学名を入れろという指定がついてきます。

PDF化するコマンド

sublimeに印刷機能がなかったし、ページ番号と大学名を入れる方法も分からなかったので調べたらenscriptっていうコマンドがあるらしい。
なければbrewとかでインストールすれば使える。
enscriptではpost scriptのファイルが吐かれるので、pstopdfコマンドを使うとPDFに出来る。

ライブラリをlibrary.cppっていう１つのファイルにまとめておけば、下みたいな感じでうまいことPDFに出来ます。

enscript --highlight=cpp --line-numbers --color --header='HOGE University (team: HOGE)                                                                                      Page $% of $=' --landscape -2 -o library.ps library.cpp
pstopdf library.ps

１つのファイルに纏める方法ですが、以下のpython スクリプトを走らせたりすればいいです。

import os

for f in os.listdir("."):
  if f == "library.cpp":
    continue
  if ".cpp" not in f:
    continue
  print("#"*50)
  l = (50-len(f)-2)//2
  r = (50-len(f)-2)-l
  print("#"*l+" "+f+" "+"#"*r)
  print("#"*50)
  print("")
  with open(f) as f:
    print(f.read())
  print("")

ページ番号とかを挿入できるツール

すでにPDFがあって、ページ番号とかを挿入したいって場合はこのサイトを使えば良さそう。
カスタム設定とかを使えば左上に大学名っていう指定も出来そう。
オンライン操作でPDFファイルにページ番号を追加www.ilovepdf.com

2017-04-18

ARC070 F「HonestOrUnkind」別解

HonestOrUnkind解いた。
面白かった、と思ってたけど想定解と違ったのでメモ。
想定解の方がシンプルで頭良い。

正直者を１人特定するところが違った。

まず人の組をN/2組作る。Nが奇数なら１人余る。
各組(x,y)について? x yを質問する。
'N'なら正嘘、嘘正、嘘嘘のどれか。（消費される人数は正≦嘘）
'Y'なら正正、嘘正、嘘嘘のどれか。（yを嘘にするためには嘘を２人消費する）
'Y'だった組のyを全部取り出した集合をTとする。
元々のN人では正＞嘘である点に注意すると、

Nが偶数：Tの中でも正＞嘘になる
Nが奇数：Tの中では正≧嘘になる

奇数のときも正＞嘘にできれば、再帰的にやることにより高々2(N/2)=N回の質問で正直者を１人見つけられる。
奇数のときは「Tのサイズが偶数なら余っていた１人をTに追加する」とすることにより、正＞嘘にできる。
いちおう証明↓
Tのサイズが奇数の場合：偶奇の性質上、正=嘘は成り立たないので正＞嘘
余りが正だった場合：元が正≧嘘なんだから正を足せば正＞嘘になる
余りが嘘だった場合：Nが偶数だった場合を思い出すと、正＞嘘になっていてつまり正が嘘より２人以上多いので、嘘を１人足しても正＞嘘
どのケースでも正＞嘘になっていることが分かる。

2017-04-17

AGC013F「Two Faced Cards」別解（TLE）

Two Faced Cardsのコンテスト中に書いてた別解です。
O(N√N logN) なのでTLEして悲しいけど紹介しておきます。
（想定解も理解したけどadhocで面白かった）

下ごしらえ

とりあえずC（とA,B）を0~Nに変換しておきます。
Ai＜Biのカードは裏で使うメリットがないのでBi=Aiとしておきます。
これでBi≦Aiが保証できます。
ここで、Bi＞Nのカードがあるとimpossibleです。
Ai＞Nのカードは後々面倒なので、Ai=Biとして後で答えから-1します。

f:id:snuke:20170417220211j:plain
入力例２の図示

クエリが１個

クエリが１個のときは以下のような貪欲法で解けます。

貪欲法

Ci=0のカードから順に相方を決めていきます。
相方は、以下のようなmultiset<int> Sを使って決めます。

Ci=xのカードを見ている時、Bi=xのカードがあればSにAiを追加する
Ci=xのカードを見ている時、Sに入っているxを全てINFに置き換える
Sに入ってる数のうちの最大値を、今見ているカードの相方として取り出す

表として取れるものがあればそれを取り、なければ表も取れるようになるのが最も遅いものを取るというイメージです。
途中でSが空なのに取り出そうとしたら失敗です。

クエリ複数

次に、クエリが複数ある場合に発展させます。

dp[k]=「Ci=kのカードには相方が要らないときの答え」というのをk=0~Nについて全て前計算しておけば、各クエリに高速に答えることができます。
dp[k]の値はどうなっているでしょうか。

まず、k=Nのときに上のような貪欲法を行うと、どこかの時点で（少なくともk=Nでは）失敗します。
失敗したiをtとすると、k＞tは全て-1です。
k=tについてはCi=tのカードをスキップしてそのまま貪欲法を続けて答えdp[t]を求めます。

あとはk＜tの場合です。
あるkについて、上の貪欲法でSから取り出した値がINFだった場合は、kをスキップしてもk+1をスキップしても答えは変わらないのでdp[k]=dp[k+1]です。（どっちにしろどうせINFを取るだけで、Sの他の部分は変わらないため）
Sから取り出した値がx(≠INF)だった場合はどうでしょうか？
xを取る代わりにスキップをして貪欲法を続けたときにxが取られるようなiをpとします。
pは必ずx以下になります。（i=xになるとxがINFに変化し、S内の最大値になるため）
このとき、kをスキップした時のSの状況はpをスキップした時とほぼ同じです。異なるのはxが取られるタイミングだけです。
したがって、p＜xのときはdp[k]=dp[p]、s=xのときはdp[k]=dp[p]+1となります。

pを全てのkについて求めておけば後ろからdpテーブルを埋めていくことができます。

ーーここまでが本質ーー

データ構造パート

pはどうやって求めればいいでしょうか？
xが取られるよりも先に取られうるものというのは以下の２通りです。（kの時点ではINFはSに入っていない点に注意）

Sにまだ残っているもの（値をyとすると、k=yでINFに変わって取れるようになる。y＜xが成り立つ点に注意）
aiがx以上のもの（bi以降なら取れる）

starry sky tree（区間加算と区間minができるsegtree）を用意します。このsegtreeでは「ある区間内で値がhoge以下のものが初めて現れる場所」をO(log N)で求めることもできます。
上記の２通りについて[取れるようになるi, N+1)に+1をしておきます。
さらに、i=0~Nについて[i,N+1)に-1をしておきます。
kの値をvとすると、[k,x)のうちv-1以下の値が初めて現れる場所がpです。（なければp=x）

上の２通りの+1を適切な順で追加したり消したりしながらpを全部求めていきたいのですが、残念ながらMoをするしかないです。
ただ、普通にMoをすると右端を1つずらしたときに増えるものが多かったり少なかったりして壊れるので、右端を「何番目に小さいbiであるか」を基準として√N個ずつのブロックに分けます。しかし、それだけだとbiが全くないゾーンが広過ぎになってしまうことがあるので、biの集合に0~N追加した集合での何番目かを基準にするなどしなければいけません。
左端については1つずらしたときには高々１つしか増減しないのでOKです。

まとめ

というわけで最後の方わりと雑でしたが、O(N√N logN)になりました。
こんな感じで６ケースだけTLEしましたΩ＼ζ°)

2017-02-03

根付き木のハッシュ

Algorithm 競プロ

りんごさん方式の記事
正当性を証明するにはユニークな多項式の形にしてmodを素数にしておけば、Schwartz–Zippel lemmaを使えて良いっぽい？
multisetのハッシュ、0になりやすいなぁと思ったけど0になる確率がat most n/modなのか。
記事での根付き木のハッシュの取り方を日本語にまとめておきます。

根付き木のハッシュ

まず、深さ i に対応する乱数 $R_i$ を [0,mod) から取っておきます。
深さ i の頂点 v について、v 以下の部分木のハッシュは「( $R_i$ + 子のハッシュ)の積」として計算します。
特に、葉のハッシュ値は 1 です。

詳細は実際に記事を読んで下さい。

僕は根付き木のハッシュをどのように取っていたかの話と、それは良くなかったという話と、ならどうすれば良かったかの話をします。

計算法

まず素数p,modを用意する。
ある頂点以下の部分木ハッシュ値は「p + 子のハッシュ値の二乗和」として下から計算していく。

（適当なハッシュ関数 h を用意して、h(子の値)の和とかにしてもよいと思う。h(x)=x^2がうまく行きそうかつ簡単だったので二乗和にしたというだけ。）
（h(x)=hoge*xとかだと、深さの分布が同じで形が違う木を同一視してしまうのでダメ）
↑冷静に考えて、二乗和もあんまりよくないので（おまけ参照）、適当な良質かつ高速なハッシュ関数を使いましょう。
適当なハッシュ関数の例としては、Murmurハッシュとかを聞きました。
CF等のhackありのコンテストで使う場合、p にあたるものとか（あるなら）ハッシュ関数のseedとかをランダムに生成すると良さそう。

妥当性について

結局、根付き木のハッシュというのは、multisetのハッシュをどう取るかという話で、
「適当なハッシュ関数」が本当に良質なら、正当なハッシュの取り方だと思います。

備考

こっちの方式の売りなのですが、いわゆる全方位木DPが簡単に書けて任意の根に対する任意の部分木のハッシュが計算できます。
（全方位木DPについてはそのうち書くかも？）（←書いた。関連記事見て。）

応用編としては、頂点にラベルがある場合は、ラベルに対する乱数的なハッシュをとり、それをpの代わりに使うと良いでしょう。

練習問題：http://codeforces.com/contest/763/problem/D https://atcoder.jp/contests/nikkei2019-2-final/tasks/nikkei2019_2_final_d
mod=1e9+7くらいでハッシュ取ると衝突するので、注意。この記事も参考に。

おまけ

~~二乗和のハッシュの正当性の証明（または嘘である証明）お待ちしてます。~~
~~ちなみに４回くらい使いましたが全部ACしてます。~~
11頂点の木どうしで撃墜できました。 A - tree hash

撃墜例（ネタバレ？）

おまけ２

りんごさん式の方の全方位化ですが、

深さとハッシュ値のペアを計算していく
$R_i$ を外に出す（下記ツイート参照）

をすると、そこまで大変でもないかもしれません。
積なので左右からの累積積を持つか、logかけて逆元（0除算は確率的に無視して良い）でやるかしないといけない点はマイナスだけど、ハッシュ関数を用意しないで良い点はプラスという感じか。

Π(x_d + c_i) ではなく x_d + Πc_i にすると定数倍良好になると見込んでいる
多項式として異なる証明は x_(d-1) に対する因数分解から再帰的に従う
— 熨斗袋 (@noshi91) 2019年12月19日

※ この記事は2019/12/19にわりと加筆修正されました。

2017-02-03

ハッシュの衝突

競プロ雑談

ローリングハッシュや木の同型判定など、競プロでハッシュを使う機会はあるけど、衝突確率とか必要な値域についてあまり考えたことがなかったので計算してみた。

ハッシュを使う場合、以下の２通りではそれぞれ衝突確率がかなり異なる。（少し考えれば当たり前の話だけど、「ハッシュ」という言葉で思考停止しがち）

２つが同じであるかを比較するだけ

２つのハッシュが衝突する確率さえ低ければ良い。
比較を何度かする場合でも、それほど大したことはない。
１回の衝突確率をp (大抵の場合 1/mod とか) として比較回数をNとすると、1回でも衝突する確率は1-(1-p)^N。
N=1e5, mod = 1e9 としたときに１回でも衝突する確率は 0.01% ほど。まぁ十分低いと言えそう？100個テストケースがあっても衝突確率は1%くらい。

種類数を求めたりする

いわゆるお誕生日。かなり衝突しやすい。
少なくとも N=1e5, mod = 1e9 としたときに答えが変わる確率は 99.3% ほど。むしろ衝突しない方がおかしいくらい。
まぁ、２つのハッシュの比較を二乗回やってると考えればわざわざ上の項と別にして書くことも無い気はするが、でもやっぱり衝突のしやすさには有意差がある。
ちなみに、１組でも衝突が起きる確率が50%を超えるようなNはだいたいO(√mod)くらいらしい。
とりあえず各N,modごとの衝突確率の表は以下の通り。（値は適当に丸めてある）
O(√mod)くらいというのも納得できる。

_mod\^N	10³	10⁴	10⁵	10⁶	10⁷
10⁹	0.05%	4.9%	99.3%	100%	100%
10¹⁰	0.005%	0.5%	39.3%	100%	100%
10¹¹	5e-4%	0.05%	4.9%	99.3%	100%
10¹²	5e-5%	0.005%	0.5%	39.3%	100%
10¹³	5e-6%	5e-4%	0.05%	4.9%	99.3%
10¹⁴	5e-7%	5e-5%	0.005%	0.5%	39.3%
10¹⁵	5e-8%	5e-6%	5e-4%	0.05%	4.9%
10¹⁶	5e-9%	5e-7%	5e-5%	0.005%	0.5%
10¹⁷	eps%	5e-8%	5e-6%	5e-4%	0.05%
10¹⁸	eps%	5e-9%	5e-7%	5e-5%	0.005%

まぁmodを1e18くらいにしとけばそうそう落ちることはないって感じかな。
ただ、それだと掛け算がオーバーフローするからバイナリ法とかで計算しないといけなくなって重そうだしなぁ。
とりあえず僕は1e9くらいのmodのpairにすることにしてライブラリを作ってみました。

1e9くらいの素数：
999999797
999999883
999999893
999999929
999999937
1000000007
1000000009
1000000021
1000000033
1000000087
...
1e18くらいの素数：
999999999999999829
999999999999999863
999999999999999877
999999999999999967
999999999999999989
1000000000000000003
1000000000000000009
1000000000000000031
1000000000000000079
1000000000000000177
...
もっと違うのが欲しければこれで生成して下さい。

注意・参考リンク

上記はハッシュが良質なものだと仮定した時の話なので、いい加減なハッシュだとこの限りではないでしょう。
いい加減ハッシュといえば、これは知っておくと良いと思います。
あとは、これの解説とか詳しいです。
あと、根付き木のハッシュについての記事を書きました。