SPOJ Open Contest 2007

http://www.spoj.pl/ZFUN07/


4/15〜4/27までSPOJというところでプログラミングコンテストが開かれていたので、
うわさを聞きつけて私も参加しました。
結果としてはSUD以外で満点が取れてほぼパーフェクトに優勝できました。
なんだかコード書いている最中に制限時間が来たので、
ちょっともやもや感も残っていますが、勝ったのでとりあえず良しです。
2位のoxyさんとは新旧R高校理科部対決を繰り広げていました。
若い人とムキになって争うというのもなかなか楽しいものでありました。
問題の方は、これがなかなか面白い問題が多くて良かったです。
AcceptされればそれでOKの問題だけではなくて、
最適解に近いものを出したら高得点とか、
高速に解を出力したら高得点とか、
いろいろ考えさせられるものが多くて楽しかったです。
とくに今回いくつかの問題で速度的なチューニングが要求されて、
それで何日かはひたすらチューニングをやっていたので、
かなりその辺で勉強になりました。
さて期限も来たことなので、
とりあえず、いろいろとそれなりに頑張った成果を放出してみることにします。
それもっといい方法があるよ、とかありましたら
こっそり教えていただけると幸いです。

  • SUD

数独を高速に解く。
ICPCで書いた超適当なソースをそのまま送ったので、
汚い上に長くて高速化する気が起きず…。
速度に比較して得点の変動もとても少ないので、
まあよろしいのではないかと。

  • ELC

SuperCon 2000とまったく同じ問題。
私はこのコンテストで一応優勝しているので、絶対に負けるわけにはいかんよなあというか。
仕様としては点の数が3000。
さらに速度が考慮される。
スコア=(200+かかった秒数)*総長/200
方針はまずは最小全域木を求める。
じつのところ最小全域木は点がランダムにばら撒かれている時は
最適解の3%程度の近似解になってしまうので、
実行時間がスコアに影響されることから、速度が速くないと勝てない。
とくに実行時間が6秒を越えると速度によるペナルティー
中継点を入れることによる総長の改善を上回ってしまう。
よってまず最小全域木をなるべく短い時間で求めなければならない。


最小全域木は図形的特長を用いるとO(nlogn)で作ることも出来るそうなのだが、
私が用いたのはスタンダードなO(n^2)の方法。
というか、プリム法。
プリム法はSuperCon2000でプリムと独立に発見したアルゴリズムなので、
アルゴリズム中でも最も思い入れが深いアルゴリズムです。
最初は20秒ほどかかっていたのが、最適化に次ぐ最適化により、
最終的には1.2秒程度にまで高速化できた。
用いた主な技法は、既に構成した木から最短距離の点を見つけたときに、
そこにフラグを立てて、見つけたという印をつけるという代わりに、
見つけた点を配列の前方に移動させて前のほうから確定点を
埋めていくという方法(これにより配列をなめる回数が半分になる上に
フラグが不要になる。そのかわり入れ替え前の位置をそれぞれ覚えておかなければならない)
と、座標値を整数で計算すること。
整数で計算する問題としては、
座標は小数で与えられるので、厳密な最適解は出なくなるが、
座標値の絶対値が大きいので(最終的な総長も大きい)、
これによる劣化はほとんど誤差の範囲である。


続いて最小全域木に中継点を追加していく処理。
これは、各点をなめて、折れ線を検出したら、その角度が120度以下の場合、
折れ線が構成する三角形の中の点のうち各頂点との為す角が120度になるような点を
選び、その点を中継点となるようにする。
既に追加された点は、そこから伸びている三辺の為す角が120度ずつになるように移動する。
これを収束するまで繰り返す。
計算量はこちらはO(nlogn)程度。
これも収束したと判断する誤差を甘くすると線形で速度が速くなるが、
あまり甘くしすぎると誤差のびのほうが速度の伸びを上回るので適切なポイントを選ぶ。
改善が10未満だった時に打ち切るぐらいが丁度良かった。
あとは、中継点の追加を改善距離が大きいものから順にするようにすると
かなりスコアが改善する。


最終的に、7年越しの課題としては満足のいくものが出来た。
当時のより速くてスコアが良くて、さらにコードが綺麗なものが出来たと思う。
http://fxp.infoseek.ne.jp/tmp/spoj/ELC.cpp

画像からランダムインパルスノイズを除去せよという問題。
これは辛かった。
たくさん論文をあさって、たくさんアルゴリズムを実装したけど、
ほとんどゴミになった。
この問題においては解像度の低い画像が用いられるので、
ほとんどの論文のアルゴリズムでは
画像のディテールがぼろぼろに失われて、さっぱりスコアが伸びない。
評価でそれっぽくノイズが取れている図が載っていても
あんまりこの問題には使えなかった。
というか、そもそも、ウインドウサイズ3*3のアルゴリズムとかって、
たったそんな範囲のピクセルを調べて
その画素がノイズかどうか判定するだなんて、
そんなの無理に決まってる。
人間だって、広範囲に認識してその画素がノイズかどうか分かるのだから、
9ピクセルだけ見たって分からない。
結局自分で書いたアルゴリズムの方が良く伸びる。
最終的に用いた方法は、
・ライン復元
・ノイズ粗取り
・取りにくいノイズ取り
・とりこぼし取り
の4フェーズによるもの。


ラインの復元は、もともとあったと思われるラインを復元する。
解像度が低いデータなので、幅1のラインが存在した場合、
それにノイズが載ると、ノイズでないところがノイズと判定されたりする。
さらに悪いことにノイズと判定されたらメジアンフィルタでは
このようなラインは絶対に復元できない。
さらに、ラインなので、周りと画素値が大きく違うことが多く、
これがぼろぼろになるとスコアは大幅に悪くなる。
よって、まず最初にラインを復元する。
縦横にラインを走査してきわめて近い画素値が長く並んだところを
ラインと認識して塗りつぶす。
途中、乱数ノイズが載っているので、いくつか違う画素値があっても飛ばす。
とりあえずサブミットしたやつでは2つまで無視する。


ノイズの粗取りでは、絶対にノイズと思われるピクセルの復元を行う。
4近傍全てのピクセルとある程度大きく離れている時などに誤差と判定する。
実際にこれはほとんど間違いなく誤差を判定する。
平坦部のノイズはかなり取れる。
エッジ部を含む取りこぼしがかなり多いが、それは後のフェーズで削除する。


取りこぼしノイズ取りはノイズ除去のメイン部分。
方法としては、あるピクセルの縦横の平均、上下の平均、
上下左右それぞれ2ピクセルからの連続から推定される自然な画素値、
それぞれとの実際の画素値との差の絶対値の最小値をスコアとして、
そのスコアの大きいものからノイズと判断して
メジアンで埋める。
閾値は入力で与えられる乱数度から決める。
具体的には(Q*W*H-(上2フェーズで復元した画素数))/300ほど復元したところで打ち切る。
これでそれらしき部分が復元される。
もちろん誤判定や取りこぼしもあるが、
誤判定はどうしようもないとして、
取りこぼしは次のフェーズで取り去る。


取りこぼし取りは最後のノイズ除去を行う。
ここでは誤判定は許されないが、そもそも大きなノイズはほとんど取られているので、
2つめのフェーズの閾値を下げたような処理を行う。


これがそのソースである。
ファイル名に6と付いているが、本当は9まで作った。
でも後続のはあんまり良くなかった。
http://fxp.infoseek.ne.jp/tmp/spoj/DIP6.cpp

  • SQRT2

さて、SQRT2である。
なんだか、妙に日本人率が高い。
問題は√2を小数点200万桁まで求めよというもの。
まず最初にお断りだが、FFTは必須になる。
はじめはHaskellで適当に書いて10万桁ぐらいで辛そうだなと思って、
C言語に変えた。
C++じゃなくてCなのは覚悟の表れであろう。
伝家の宝刀C言語を抜いたらからにはもう負けるわけにはいきませんね、という。
それでも最初は感覚をつかむために分割統治によるアルゴリズムを実装して、
ああ、本当に桁数2倍で速度3倍になるなあとか思いながら、
この方法は10万桁強で断念した。


http://momonga.t.u-tokyo.ac.jp/~ooura/fftman/index.html
FFTは基本的にこのページに載っているものを参考にした。
再帰+スクランブラーのタイプで、
メモリを走査する順を工夫してキャッシュ効率を上げて
ちょこちょこ高速化していると、
30万桁、60万桁と桁が伸びていった。
120万桁がなかなか通らなくて苦労していたら、
速度を一気に5倍ぐらいにするアイデアが出てきて
それであっという間に戦いは終わってしまった。
最初、イテレーション回数が問題になって、
ニュートン法イテレーションごとに有効桁数が2倍になるので、
初期に持つ桁数が多ければイテレーション回数が減らせる。
それで、最初に100桁とか200桁とか持っておけば回数が大分減らせて、
60万桁を通した時は600桁ぐらい埋め込んでいた。
それからFFTを高速化して、120万桁が1200桁ぐらい埋め込んでも
さらにもう一回イテレーションを減らさないと通らない、
しかし、コードサイズ制限が4096バイトだから、もうこれを2倍にはできない、
本当にどうしたものか、となっているところに、
ピコーンと天啓がやってきて、
え、そんなの60万桁計算した時は60万桁まで求まってるやん?
その60万桁を用いて1回イテレーションすればよろし、となって、
もっと言うと、60万桁求めるのも30万桁から1回のイテレーションでよくて、
要するに、桁数を倍倍にしていけば、
計算総時間としては200万桁のFFT2回分ぐらいで終わってしまう、
という驚愕の事実にたどりついてしまった。
ええと、もしかしたらこの辺のことは知らないのは私だけだったりするのでしょうか。
とにもかくにも桁数を倍倍にしながら計算するようにすると
本当にすんなりと通ってしまった。
高速化のアイデアはまだまだいくつもあったので、
ちょっと残念なところではあったが。
それをとりあえずここに書きなぐっておくと、
まず、多倍長の掛け算という観点から言うと、
FFTにおけるスクランブラーは完全に削除できるということ。
FFTにおける乗算では、FFT→各要素の掛け算→IFFT
でもって、FFT
スクランブル→バタフライ、
バタフライ→スクランブル、
どちらの順序でも記述できる。
その際にバタフライ演算は若干異なる必要がある。
スクランブラーは二回掛けると恒等写像になるので、
バタフライ→スクランブル→各要素の掛け算→スクランブル→バタフライ
というようにFFTとIFFTを実装すれば、
バタフライ→各要素の掛け算→バタフライ
となってスクランブラーは外せる。
スクランブルの処理は完全にキャッシュを外れる処理なので、かなり遅く、
実際のところ提出したプログラムでも実行時間の2割以上を占めていたので、
これを外せば2割以上高速になることになる。
あとは、実数のFFTをまともに実装すれば
きっとかなり速くなるはずだ。
送ったプログラムは2つの実数列を実数部と虚数部に埋める実装になっているので、
半分の桁数のFFT2回の方がおそらく高速だろう。


そして、提出したプログラムはこちら。
http://fxp.infoseek.ne.jp/tmp/spoj/SQRT2_8.c
256桁ほど答えが埋め込んであるが、これは別に速度のためではなくて、
実際のところ、ニュートン法でもきっちり有効桁数が倍になるということはなくて、
1、2桁情報が落ちてしまうので、
最初に1桁から始めると純粋に倍倍にしていったのでは、
最終的にかなり有効桁数が少なくなってしまう。
それを防ぐには桁数の少ない時は2回ずつ繰り返すなどすればいいのだが、
面倒だったので、もう最初からある程度桁数を確保してしまうことにした。

        • -


というわけで、全体的になかなか面白い問題で、勉強にもなりました。
とくにチューニングはもう任せてください。
上にあげた問題以外も、なんというか、実装が重い問題はなくて、
どれも思いつけば一瞬な感じの問題なので、
ひらめく快感というものを思い起こさせてくれました。
さて、次の大きなお祭りはICFPですかねえ。
いまからわくわくしっぱなしです。