Arnoldi法(Lanczos法)についてのメモ

動機

久しぶりのブログエントリになります。
休んでる間何もしてなかったかというと、勿論そうではなくて、比重がインプットに偏っていたのが原因でした。
偏りなくできるのが理想ですが、あまり時間がない時はひとつのことに集中するのが、経験上効率的だと思っているので、今度書く時もまた半年ぶりとかになる可能性もありますが、まぁ自分用のメモみたいなものなので、気楽にやっていきたいと思います。

前置きはこの辺にして、本題に入ります。
最近、個人的な興味があって固有値の求め方を調べていました。
固有値は、特異値分解において重要な役割を持ち、そして特異値分解は主成分分析、潜在意味インデキシングといった応用の基本となるアルゴリズムです。特異値分解自体は、それ以外にも非常に多数の応用が考えられますが、詳細については省略します。

実際に固有値を求める方法は、学部時代に固有方程式を手計算で解き、また数値解析の授業でべき乗法とかヤコビ法を授業で聞いたなぁという程度でしたが、自然言語処理等でよく使われる大規模疎行列でも高速に計算可能と言われる方法について、一度はちゃんと勉強しておかなければいけないと思ったのが、今回調査したきっかけとなります。

なお、実用上は大規模疎行列の特異値分解や固有値を求めることが必要な状況であれば、ARPACKやredsvdを利用すれば良いと思いますので、別に興味がなければ知らないままでいいと思います。

事前準備

そもそも固有値は行列に対して、下記の式で成り立つような及びベクトルが存在する時にを固有値と呼びます。

固有値を考える上では複素ベクトル空間で考えなければならないので、の随伴行列(転置行列の複素共役行列)をと表すとします。
ここで、 (は単位ベクトル)が成り立つ行列をユニタリ行列と言います。
ユニタリ行列は、各列ベクトル(または各行ベクトル)をとすると、内積が以下のように記述できるので、各列ベクトルは正規直交基底であると言えます。

正規直交基底の定義
(内積に下記の関係が成り立つベクトル列であること)

こういった正規直交基底は、グラムシュミットの正規直交化法というアルゴリズムで任意の独立なベクトル列から生成することが可能です。

グラムシュミットの正規直交化法

1. Start: Compute  := ||. If  = 0 stop, else  := /

2. Loop: For j = 2, . . . , p do:

    Compute  := ( , ) for i = 1, 2, . . . , j − 1

     :=  - 

     := ||

    If  = 0 then stop, else  := /

固有値を求めるフレームワーク

大規模行列の固有値を効率的に求める際に一般的に、下記のようなフレームワークがあります。

1. を部分空間に射影し、新しい行列を得る
2. の固有値を求める

ここで、の固有値がの固有値の良い近似であることを証明すれば、このフレームワークは十分に機能することとなります。
一般に1の過程はReductionと呼ばれ、よく名前の上がるArnoldi法やLanczos法は、1の過程でKrylov部分空間に射影する手法のことであり、の固有値を求める際はArnoldi法であればQR分解などまた別の手法を用います。

Krylov部分空間

Krylov部分空間は、べき乗法のアイディアを自然に発展させたものと考えることができます。
べき乗法は、初期値を適当なベクトルを取ると、

という計算を繰り返すことによりが絶対値最大の固有値の固有ベクトルに収束する性質を利用した固有値の求め方です。(証明はこちら)

べき乗法の場合、途中の計算結果は捨ててしまいますが、Krylov部分空間を用いた手法は途中の計算結果を利用して下記のようなベクトル列を構成します。

このベクトル列を基底とする空間をKrylov部分空間と言います。
Krylov部分空間は、べき乗法と同様、絶対値最大の固有値をよく近似することが可能であり、またそのままべき乗法を利用する場合は最大固有値しか得られないという問題があるのですが、一度Krylov部分空間に射影してから固有値を求めるという方法を取ることにより、より高速に計算できるようになります。

Arnoldi法

それでは実際に非対称行列をArnoldi法を用いてKrylov部分空間に射影する方法について考えます。
対称行列の場合は、より計算が簡便なLanczos法を用いますが、これはArnoldi法の特殊ケースであると考えることができます。

をKrylov部分空間に射影するアルゴリズムは下記のように計算ができます。

1. Start: Choose  where 

2. Loop: For j = 1, . . .  do:

    Compute  := ( , ) for i = 1, 2, . . . , j

     :=  - 

     := ||

    If  = 0 then stop, else  :=  / 

このアルゴリズムが具体的にどういうことを行っているのかということを説明すると、まずn回目のイテレーションで、Krylov部分空間を構成するベクトル列からなる行列をと置き、AをKrylov部分空間に射影した行列をと置くと、は各ベクトルをグラムシュミットの正規直交化法により正規化し、下記のようにを求めていることとなります。

求められたは、ヘッセンベルグ標準形と呼ばれる形式となります。

ヘッセンベルグ標準形

ヘッセンベルグ標準形はギブンス回転を用いたQR分解を実施することで固有値を高速に計算可能です。

の固有値 の固有値の証明

では、次にの固有値がの固有値の近似となるかを考えます。
n回目のを考えると、誤差項をと置くと

が成り立ちます。
ここでの固有ベクトルを、固有値をと置いて、上式の両辺にをかけます。

ただし、

Arnoldi法の収束条件はこの誤差項が十分に0に近づく時ですので、

が成り立ち、これはの固有値を求める式と一致します。

まとめ

• Lanczos法やArnoldi法は、固有値を直接求めるアルゴリズムではなく、Krylov部分空間に射影するアルゴリズム
• Krylov部分空間はべき乗法のアイディアを発展させたもの
• 実用上は既存パッケージを使えばいい

Arnoldi法、Lanczos法はいくつも派生があり、非対称行列にも適用可能なLanczos法なども存在します。
部分空間に射影する方法は、redsvdが実装している方法も同様で、ランダム行列を使って射影を行い、特異値分解を行っているだけですのでアルゴリズム自体は実装も簡単です。ただ、理論的に近似誤差を求めているのですが、丁寧に説明が書かれているにも関わらず、まだよく理解できていません。特異値分解や固有値の求め方は、どちらかというと工学よりも理学(数学)よりですので、理論が難しくなりますね・・・。

また説明の際は簡単のため、各行列, ベクトルの次元を省略しています。
厳密な説明は参考文献を参照してください。

参考文献

数値解析(森 正武)
数値計算の基礎
Numerical Methods for Large Eigenvalue Problems
英語版wikipedia(Arnoldi Iteration)
固有値問題ノートの補足
固有値計算アルゴリズムの最近の進展

Giraphの解剖記録

導入

Giraphのソースコードを読んで見たので、そのメモとして残しておきます。

GiraphというのはHadoop上で動作する大規模グラフ処理システムのことです。 オリジナルはPregelといわれるGoogleが論文を出したBulk synchronous parallel(BSP)を用いたグラフ処理システムのようです。

もう少し詳しい導入が必要な方は、smly先生のTokyoWebMiningでの発表を見てください。

さて、あまりに適当な導入を飛ばしてまずはコード例を見ていきます。

コード例を見てみる

Giraph自体は、グラフ処理を行うためのフレームワークですので、具体的なアプリケーションコードは用途に合わせて書いてあげる必要があります。

具体例として、Giraphのソースコードにはたくさんのexampleがあるので、その中で最短経路を求めるプログラムを見てみます。 細かい部分は見ず、今は薄目でざーっと眺めてみます。

public class SimpleShortestPathsVertex extends
    EdgeListVertex {
  public static final String SOURCE_ID = "SimpleShortestPathsVertex.sourceId";
  public static final long SOURCE_ID_DEFAULT = 1;
  private static final Logger LOG =
      Logger.getLogger(SimpleShortestPathsVertex.class);

  private boolean isSource() {
    return getVertexId().get() ==
        getContext().getConfiguration().getLong(SOURCE_ID,
            SOURCE_ID_DEFAULT);
  }

  @Override
  public void compute(Iterator msgIterator) {
    if (getSuperstep() == 0) {
      setVertexValue(new DoubleWritable(Double.MAX_VALUE));
    }
    double minDist = isSource() ? 0d : Double.MAX_VALUE;
    while (msgIterator.hasNext()) {
      minDist = Math.min(minDist, msgIterator.next().get());
    }
    if (LOG.isDebugEnabled()) {
      LOG.debug("Vertex " + getVertexId() + " got minDist = " + minDist +
          " vertex value = " + getVertexValue());
    }
    if (minDist < getVertexValue().get()) {
      setVertexValue(new DoubleWritable(minDist));
      for (LongWritable targetVertexId : this) {
        FloatWritable edgeValue = getEdgeValue(targetVertexId);
        if (LOG.isDebugEnabled()) {
          LOG.debug("Vertex " + getVertexId() + " sent to " +
              targetVertexId + " = " +
              (minDist + edgeValue.get()));
        }
        sendMsg(targetVertexId,
            new DoubleWritable(minDist + edgeValue.get()));
      }
    }
    voteToHalt();
  }
}

どうやらDijkstraを実装しており、グラフの各頂点ごとにcomputeというメソッドが実行されるようです。 処理の主体が頂点であり、各頂点からは、隣接する頂点に対してメッセージを送信したり、受け取ったりすることで最短経路を求めることができるようです。 おそらくメッセージのやりとりがなくなった段階で収束したと見なして処理が終了するのでしょう。

それではこのコードは実際どのように呼ばれるかをジョブの投入から見ていきたいと思います。

ジョブの投入からの処理の流れ

実のところGiraphはHadoopのMapReduceを利用しています。 ですが、一般的なMapReduceアプリケーションのように、MapReduceを何度も実行して上記のcompute関数を呼び出している、というわけではありません。

MapReduceもReducerは使わずMapperだけが使用され、たくさんのサーバ上でプログラムを実行するためのディスパッチャとして利用されているだけです。

メッセージのやり取りについても、独自にhadoop.ipc.RPCを利用してネットワーク通信方法を構築しており、また各サーバ間の同期にはzookeeperを利用しています。

それではcomputeがどのように呼ばれるか、処理の流れを追っていきます。

Flow-1
()内は処理が書かれているファイル名及び関数です。
1. main関数からMapReduceジョブを実行(GiraphRunner.java main)
  -> 2. MapReduceジョブの各種設定及びドライバーコードを実行(GiraphJob.java run)
    -> 3. Mapperタスクの並列起動(GraphMapper.java run)
      -> 4. 各Mapperタスクの役割(後述)の決定と初期化(GraphMapper.java prepare)
        -> 5. スーパーステップ(後述)の開始(GraphMapper.java map)
        -> 6. 自分に割り当てられたグラフの頂点取得(GraphMapper.java map)
        -> 7. 各頂点ごとにcomputeメソッドを実行(GraphMapper.java map)
        -> 8. 全処理が終了していなければ次のスーパーステップを実行するため5に戻る(GraphMapper.java map)

Mapperタスクの役割や、スーパーステップという聞きなれない言葉が出てきました。 スーパーステップというのは何でしょうか？

Bulk synchronous parallel(BSP)

最初にも紹介しましたがGiraphはPregelのオープンソース実装です。PregelはBSPを利用しているという話もありました。 BSPというのは並列アルゴリズムを実行するための計算モデルです。

次の画像はBSPの処理を表現したものです。

実際のところそう難しくはなさそうです。下記の処理をスーパーステップという1単位として、何度もスーパーステップを繰り返すことにより、上記の最短経路のような問題を解くことができます。

1. Concurrent computation: 各サーバがそれぞれローカルの情報を元に計算
2. Communication: 各サーバがデータを交換
3. Barrier synchronisation: 全ての計算及びデータ交換が終了するまで待つ

上記のFlow-1では、5-8が1つのスーパーステップになっていることがわかります。「Concurrent Computation」及び「Communication」は最短経路プログラムでのcomputeメソッドが行なっていますね。「Barrier Synchronisation」はどこに相当するかというと、Flow-1 8の時に準備ができるまで待つという処理が含まれており、zookeeperに今回のスーパーステップ終了だよと書きこまれれば各サーバーはその情報を参照して次のステップへ進んでいきます。

それではその「スーパーステップ終了の合図」は誰が書き込むのでしょうか？　次の「各Mapperの役割」で説明します。

各Mapperの役割

スーパーステップを実行をするためには、各頂点における計算が終了したかどうかといった同期情報を管理・監視する処理も必要であることに気づきます。 Giraphはこのような処理を実行するために、下記のような役割があります。

これらの役割は重複して持つことも可能ですし、各サーバが別々に持つこともできます。

• Master: スーパーステップの管理・監視を行う。1~複数台(パラメータで設定可能)で起動される。(BspServiceMaster.java, MasterThread.java)
• Zookeeper: Masterノードと同じサーバで実行されることが多い. Zookeeperを起動して分散協調管理を行う。1~複数台(パラメータで設定可能)で起動される。
• Worker: 「compute」メソッドを実行するMapperタスク(Graph Mapper.java, BspService.java, BspServiceWorker.java)

注意点ですが、どの役割もMapperタスクで実行されるので、HadoopのいわゆるMasterノードと混同しないように気をつけてください。

ネットワーク通信

最後に各サーバ間でどのような通信が行われるかを説明します。

スーパーステップの同期や、グラフの頂点の割り当て情報などはZookeeperを通じてMaster, Worker間でやり取りをします。

各ワーカー間で送られるメッセージについては、毎回メッセージが送信されるわけではありません。 ある程度大きな塊となった段階で送るか、それともスーパーステップの終了段階で送るかの２通りにわかれます。 こちらの処理は、hadoop.ipc.RPCを通じて、各サーバに対して直接メッセージの送受信を行います。

例えば5台、ワーカが起動しているサーバがあれば、ネットワーク通信は最大5x5=25回、特にスーパーステップが完了するあたりで実行されるということになります。 (ソースコード的には、BasicRPCCommunications.java のPeerFlushExecutorやLargeMessageFlushExecutorでproxy.putMsgと自然にリモートメソッド呼び出し(RPC)しているあたりがネットワーク通信に該当します。)

感想

ソースコードを熟読したわけではなくて、さらっと見ただけですので間違えているかもしれません。そもそもSVNのTRUNKですので、公開した段階で既に修正されている部分があるかもしれません。 必ずこちらのメモを鵜呑みにはせずに自分でも調べるようにしてください。

それから、コード読んでてブログにする段階で見つけてしまったんですが、Giraphの内部の仕組みを解説するパワポはこのスライドが詳しいですね。。。なので多分続きは書きません。

あと、どこかで「キリンさんが好きです。でも象さんの方がもっと好きです」というベタベタなネタを入れたかったんですが、入れる場所がありませんでした。

参考URL

Giraph Introducing Apache Giraph for Large Scale Graph Processing グラフ問題とバルク同期並列の常識をGiraphで体得 TokyoWebMiningでの発表 Wikipedia - Bulk synchronous parallel HadoopのIPC/RPC

Hadoop開発者/管理者トレーニングに参加してきました。

前置き

少し時間が立ってしまいましたが、Hadoop開発者/管理者向けトレーニングに参加してきましたのでその感想をつらつらと書きたいと思います。

なお、Hadoop開発者/管理者向けトレーニングというのは、Clouderaという会社が実施しているセミナーです。 Clouderaは、社員には相当数のHadoopコミッタが在籍していて、Hadoopのサポート・コンサルティングを行なっている会社です。

正直、自分がお金を出したわけでもありませんし、仕事の中で勉強できるというありえないぐらい恵まれた時間だったので、受けることができて良かったなぁとしか思っていません。 (その分、会社にはちゃんと還元しないとと思っていますし、それだけのものを得られるようにと、可能な限り真剣に講義を受けたつもりですが……) そのため、かなり好意的な感想になっていると自覚していますので、冷静に受講するかどうかを検討したい人はその分を割り引いて読んでいただければと思います。

トレーニングで何ができるようになるか

トレーニングの中身で実際どのような内容を勉強するかについては、下記のリンクから確認できますので詳細はそちらをご覧頂きたいと思いますが、さらっと説明します。

開発者向けトレーニング 管理者向けトレーニング

開発者向けトレーニングでは、Mapper, Reducer, Combiner, Partitionerといった基本プログラムの書き方から、自作Writable, 自作InputFormatの書き方まで学びます。 変わったところでは、MapReduceを使った幅優先探索の手法や、セカンダリソート、2種類のデータの結合(SQLでいうJOIN)をMapReduceで実装する方法なんかまでカリキュラムに入っています。 (結局JOINするならJavaで書くと大変だから、Hive使いましょう…だったりしますが…)

管理者向けトレーニングでは、推奨ハードウェア構成や、チューニングが必要なパラメータであったり、Hadoopで利用可能なジョブスケジューラの説明とその利用シーンなど、主に運用面での注意点・ノウハウを学ぶことができます。

個人的な意見ですが、「Hadoopを使えば何ができるか？」を知りたい人は開発者向けトレーニングを、「Hadoopを実運用で使っていくにあたって知っておかなければならないこと」を知りたい人は管理者向けトレーニングを受講するのが良いのかなと思います。

ですので、まだ検討段階で、会社で使うかも？というぐらいの人は、現在システム管理者やプロジェクトマネージャの方でも、まず開発者向けトレーニングを受講して、Hadoopが得意なこと、得意じゃないことを把握するのがいいんじゃないかなと思いました。

トレーニングの流れ

トレーニングは、上記の内容について、2時間程度の説明を聞いてから、30分〜1時間程度の演習実施(実際にプログラムを書いたり、パッケージをインストールして設定したり)を3~4日間ひたすら繰り返す、という流れです。

開発者トレーニングの前提条件にはプログラミング経験者が対象と書かれていますが、プログラミングの演習問題は解答も一緒に配られるので、プログラミングに自信がなくても大丈夫です。例えば新人研修の中で取り入れたとしても、問題なくついていけるんじゃないかなーと思います。また、Hadoopの経験もゼロでOKです。実際、受講した人の半分ぐらいは初めてHadoopを触ったという方だったように思います。

受講してよかったこと、残念だったこと

個人的にこのトレーニングを受けてよかったと思うところを箇条書きにしておきます。

1. 講師の説明がわかりやすく、理解できるまで教えてくれる。 Clouderaはつい先日(4月26日)に日本法人ができたばかりの会社ですが、講師の方は他社でも元々トレーナーをされていた方で、トレーニング経験豊富なせいか、予備校講師の説明を受けているような感じでわかりやすかったです。 (2日目まで気付かなかったのですが、自分はこの方の講義を受けるのは2回目でした。)

また、ちゃんと理解できるまで先に進まず丁寧に教えてくれます。その時間も含めてカリキュラムが計算されているようで、時間がないので説明を飛ばしたりというようなこともありませんでした。

1. 少しでも疑問に思ったことは、どんどん質問ができる。 自分が受けた開発者向けトレーニングの時は特にそうだったかと思うのですが、受講者はちょっと気になったことでも、どんどん質問していました。例えばnamenodeがeditsファイルに書き込む際は1エントリごとにディスクとシンクしているの？とか、とかdatanode間のプロトコルやポートの一覧は？とか、SqoopのコネクタはMySQLだとJDBC接続しかなさそうだけど遅いの？(いや、MySQL用にちゃんとチューニングされてます)とか。 講師の方がその場ですぐ答えるのが難しい場合は、トレーニング終了までに調べて教えてもらえます。

1. まとまった時間をとって一つの技術について勉強するのは効率が良い。 なかなか業務内で一つの技術を勉強するまとまった時間というのは取れないものなので、トレーニングなどでじっくり勉強できる時間が持てるというのは良いと思います。

1. ハンズオン演習で触りだけでも触れておくと、後で導入する時がスムーズ 例えばFlumeのセットアップをその後仕事場で行う必要があったのですが、Getting Startedを終えているので、とっかかりの1,2時間はかかる調査を省略できました。

一応、残念なところも書いておきます。

1. 演習の時間も実際はそれほどないため、トライアンドエラーやトラブルシューティングによるノウハウはあまりたまらないかなーと思います。そのあたりは、結局トレーニングでは解消が難しく、業務を通して実際に使っていかないとわからないところかもしれないですね。

1. トレーニング費用は高いので、個人での受講は難しいところ。。。実際4日間朝から晩までなので、他のトレーニングと同様妥当な額かもしれませんが、会社で受ける場合は4間も研修に出してくれるような時間的余裕があるかという問題もありますね。。

まとめ

最後にまとめです。 講義内容のレベル的には、オライリーから出ているTom Whiteさんの象本を熟読してバリバリHadoopを使っているという人には少し物足りないかなと思います。

ですが、既にHadoopを使い始めている人には基礎固めとして、これからHadoopを使い始める人にはとっかかりとして、金銭的、時間的な問題をクリアできるようであれば積極的に活用するといいんじゃないかな、と思いました。

大川をわたりました。

どうでもいいことですが、山本一力の作品に「大川わたり」という作品があります。 この作品は借金が原因で大川(現在の隅田川)をわたれなくなり、新しい生活を送ることになった大工のお話です。

そして、これまたどうでもいいことですが、１ヶ月程前、ヤフー株式会社を退職しました。 練馬から江東区の方に引越しをして、今は毎日隅田川を渡って新しい会社に通っています。

写真は橋の上から写した朝焼けの中の隅田川です。綺麗ですね。

Hopscotch Hashingを実装してみた

動機

数年前Cuckoo Hashingも知らないの？と友人から信じられないという目で見られた経験から、いつかCuckoo Hashingを超えるアルゴリズムを実装してやるという(考案してやるという程は高くない)野望を持っていました。 そこで今回はHopscotch Hashingを実装してみました。

Herlihy, M., N. Shavit, and M. Tzafrir, “Hopscotch Hashing,” Proceedings of the 22nd International Symposium on Distributed Computing (DISC), pp. 350-364, Arcachon, France, September 2008.

Hopscotch HashingはCuckoo Hashingと同様、定数時間での検索を可能にするハッシュテーブルのデータ構造です。 Cuckoo Hashingについては、ここや「日本語入力を支える技術」でも説明されている(らしい)ので詳細な説明は省略します。

Hopscotch HashingのCuckoo Hashingに勝るメリットは、高い充填率(ハッシュテーブル上の空き要素の数が少ないこと)でも高速に動作することらしいです。 (ただ、実装した限りでは、ハッシュ関数が悪いせいかあまり充填率を上げることができず、あまりメリットを感じられませんでした。) Hopscotch Hashingのアルゴリズムは、下記の通りです。

アルゴリズム

検索方法

Hopscotch Hashingでは、ハッシュテーブルを構成する各バケットにHop Informationと言われるビットマップの索引を持ちます。 Hop Informationは各バケットからH-1以内(Hは論文では32を指定)に、そのバケットのハッシュ番号を持つエントリが格納されて るかどうかを0/1のビットで保持します。

検索対象のキーは、キーから計算されたハッシュ値が指すバケットのH-1以内にあることが保証されているため、Hop Informationを元にH-1以内のエントリを調べるだけで必ず検索対象があるかないかを確認できます。 これにより、最悪計算量がO(H)という定数時間での検索が可能となります。

1. キーからハッシュ値を計算して、ハッシュ値が指すバケットのHop Infromationを調べる
2. Hop Informationでフラグが立っているバケットに検索対象のキーが格納されているはずなので、そのバケットを検索する

Hを4とした場合に、検索アルゴリズムを図示した例は下記です。 (図1)

新規キーの挿入

「キーから計算されたハッシュ値が指すバケットのH-1以内に、検索対象のキーが格納されたバケットがある」という制約を保証するために、挿入する場合は少し複雑な処理が必要となります。 Cuckoo Hashing

1. キーからハッシュ値を計算して、ハッシュ値の指すバケットから、線形探索で未使用バケットを探索
2. 未使用バケットがハッシュ値の指すバケットからH-1以内にあるならばそのバケットに挿入する
3. H-1より離れているならば、ハッシュ値に近い場所に空きバケットを作るため、未使用バケットに移動可能なバケットを調べ移動を繰り返す。
   (移動可能なバケットを探す際に、Hop Informationがあるため、ハッシュ値を再計算する必要は ない)
   もしそのようなバケットが見つからなければ、ハッシュテーブルのサイズを変更して再度ハッシュを計算しなおす。

検索と同様Hを4とした場合の挿入アルゴリズムを図示した例を示します。 (図2)

(図3)

(図4)

結果

リサイズ部分の実装が適当だったためか、unordered_mapになかなか勝てず微妙な感じでしたが、ソースコードはこちらです。 最初は実装も簡単で速度が早いが衝突しやすいハッシュ関数を使っていたのですが、あまりに衝突が多すぎたため結局unordered_mapで使われているハッシュ関数を利用することにしました。ハッシュ関数ってやはり大事ですね。。。 ハッシュ関数さえあれば、C++03のunordered_mapが利用できない環境でもお手軽にそこそこ高速なハッシュテーブル実装を書けるというところがメリットであるように思います。

また、Hopscotch Hashingは並行処理可能なHashMapとしても考案されたものであるため、マルチスレッド環境下で本来の実力を発揮できるかと思いますが、今回はそこまで実験していません。今後の課題としたいと思います。

冒頭で書いたCuckoo Hashingについては、あまり詳しくないのですが、今でも研究されているようでたくさんの亜種が発生しているようです。(d-ary cuckoo hashingや、bucketize version Cuckoo Hashingなど) Hopscotch Hashingもその亜種の一つであると言えるでしょう。 今回の実装を通じて、古典的アルゴリズムであるハッシュテーブルでさえ、今も研究され、日々進歩しているということを忘れないよう心に留めておく必要があると実感しました。

全探索の解き方

動機

ふと、自分が全探索な問題を解く時にどう解いているかをまとめて見ようと思った。 もちろん自分もまだまだなわけなので、解説というものにはならないし、あくまで自分の考えをまとめただけである。

ちなみに自分が少しでも全探索についてちゃんと理解できたかなと思えたのは、競技プログラミングを始めてからだと思う。 それまでも、もちろんスタック・キュー・再帰なんかを使って、深さ優先探索や幅優先探索を書くことがあったけど、今みたいに書く前にどう書けばいいかを整理できていたわけではなくて、消したり書いたりして長時間考え、いろいろ苦労していた気がする。

基本テンプレ

基本的に自分の場合、最小手数を求めるであったり、最短経路を求めるようなものでない限り全部を探索する必要があるものは、深さ優先探索で解くことが多い。再帰でかけるので、コード量が少なく感じるからだと思う。 まず自分は何も考えず以下のように書く。 (まぁ、lispとかで再帰を書く時のイディオムと同じ。)

  //関数名はrecもしくはdfsに決めてる。
  int rec(状態変数){
     if(<終了条件>){

     }
     <選択肢の列挙(recの再帰呼び出し)>
  }

関数名とか、これだけでも書き方を決めておくと、時間少ない中悩む必要がないので楽だったりする。 幅優先探索の場合も下みたいなテンプレートを脳内に用意してたりする。 (Stateはいちいち書くのが面倒なのでpairを使ったりする。)

struct State {
  <状態変数>
};

  queue Q;
  Q.push(<初期状態>)
  while(!Q.empty()){
    State s = Q.front(); Q.pop();

    <選択肢の列挙(Q.pushの呼び出し)>
  }

問題パターン

列挙パターン 例) 1,2,3,4の４数字からなる10桁の数字列を列挙せよ。

  vector > result; //列挙した結果を格納するコンテナ
  char cand[] = {1,2,3,4};

  //状態変数 i は、今現在何桁目までを着目してるかを確認
  //状態変数 v は、現在までに選択した各桁の数字を格納
  void rec(int i, vector &v){
     if(i == 10){
        //列挙対象に何らかの条件がある場合は、最終的に条件に当てはまるかどうかをチェックする。
        //例 末尾が4の数字列の場合
        //if(v.back != 4) return;
        result.push_back(v);
        return;
     }
     for(int j = 0; j < 4; j++){
        //枝刈りする場合はここでもチェックする。
        v.push_back(cand[j]);
        rec(i+1,v);
        v.pop_back(); //後始末を忘れずに
     }
  }

二次元探索パターン i,jの２つの状態を持つが基本同じ。 例) 数独パズルを解け

//左上から各行ごとに、スキャンしていく
//状態変数は現在着目しているセルの(i行,j列)と、数字が埋められた盤面をvとして持つ。
void rec(int i, int j, vector > &v){
  //再下端までいったら終了。条件に合うものが出力される。
  if(i == 9){
    cout << v << endl;
    return;
  }
  //右端までいったら1行下に。
  if(j == 9){
    rec(i+1,0,v);
    return;
  }
  //既に数字が入力されている場合は次へ
  if(v[i][j] != -1){
    rec(i,j+1,v);
    return;
  }
  for(int k = 1; k <= 9; k++){
    //その数字を入力しても問題ないかどうかをチェックして再帰。
    //数独なので、上下左右及び9マス内に同じ数字があるかどうか。
    if(check(i,j,k,v)){
      v[i][j] = k;
      rec(i,j+1,v);
      v[i][j] = -1; //初期値に一応戻す(戻さなくても動くはず)
    }
  }
}

ループパターン ループしてしまう場合は、既にチェックしたかどうかを管理する必要がある。

  //訪れたことがあったらこれ以上探索しない。
  vector > visited;
  void rec(int i,,...){
     if(<終了条件>){
        
     }
     int new_i = i + to[hoge];
     //遷移先が既に訪問済みの場合はこれ以上探索しない
     if(visited[new_i]) return;
     visited[new_i] = true;
     //訪問済みじゃない場合は探索を続ける。
     rec(new_i ...);
  }

まとめ

• 問題パターンごとに書き方を考えておくと楽ですよ。

個人的には、Project Euler100問解いた時に、全探索系の問題を考えまくったので、その経験が結構生きてる気がする。 あと、夢の中でDPを解いてたほどDPについて考えまくった週があったりとか、結局遠回りながらも時間をかけて考えないと成長できないんだよなぁと思ったり。 (だからって時間をかけるのがいいかというと、社会人的にこればかりに時間をかけられないので、いろいろ難しいところ・・・)

とまぁ、最近微妙なエントリしか書いてないなと思いつつも、このエントリも微妙なままこれで終わる。

兎と亀のアルゴリズムをふと考えてみた

リストの循環チェックアルゴリズムといえば、兎と亀のアルゴリズムが有名だけど、ふと疑問がわいたのでさらっと検証してみた話。 時間も労力もかけたエントリではないけど、まぁいいかと。

兎と亀のアルゴリズムというのは、リスト構造を２つずつ辿るポインタ(兎)と、１つずつ辿るポインタ(亀)を用意して、兎が亀に追いつけば、そのリストは循環していると判断できるというアルゴリズム。

その話から最初に考えた擬似コードが下。兎が２つ進む間にちゃんと亀に追いついてるかどうかをチェックしている。

bool check(List *listp){
  List *turtlep = listp;
  List *rabitp = listp;

  while(1){
    turtlep = turtlep->next();
    rabitp = rabitp->next();
    if(rabitp == turtlep) return true;
    rabitp = rabitp->next();
    if(rabitp == turtlep) return true;
    if(rabitp == listp.end()) return false;
  }
}

しかし、このコードはあまり綺麗じゃない。 で、よく書かれるコードは下の方なのかなーと思う。

bool check(List *listp){
  List *turtlep = listp;
  List *rabitp = listp;

  while(1){
    rabitp = rabitp->next()->next();
    turtlep = turtlep->next();
    if(rabitp == listp.end()) return false;
    if(rabitp == turtlep) return true;
  }
}

ただ、こちらの場合だと、兎が亀を追い越し続ける場合があるんじゃね？と思ってしまった。

そこで実際に何ターン(ループの実行回数)で兎が亀と同じ地点に立つかを計算してみる。 亀がリストの循環ループ部分に入った時に、兎が循環ループのk番目の地点にいるとする。

すると、それからtターン経過した時の兎と亀の位置はそれぞれ

従って、兎と亀の距離は、となるから、循環ループを構成するノード数をmとすると、がmの倍数になれば追いつく計算である。 つまり、循環部分で兎が亀に追いつくまでに必要なターン数はたかだかmのため、循環ループにたどり着くまでのノード数も考えると、最終的な計算量としては、O(n)で良いことがわかる。(nはノード数)

周期の一致とかを考えると最小公倍数かなぁとか一瞬考えてしまい、結構時間がかかることもある？とか思ってしまったけど、ちょっと考えてみるとなんか普通の答えで拍子抜けした話でした。