良いもの。悪いもの。

先日、 構造体を格納したコンテナに対してSTLのアルゴリズムを適用する方法 について記したが、今回はそれを C++0x で書き直してみた。C++0xの規格はまだきちんと定まっていないが、草案ではかなりのところまでつめられてきており、また対応するC++コンパイラも出てきていることから今のうちにC++0xに慣れておいた方が良いと思う。今回のプログラムで行う内容については前回と全く同じだが、C++0xを使うことによりかなりすっきり書くことができた。因みに現時点の最新の規格(草案)は N3000 (PDF)になる。日本語での資料は C++0xの言語拡張まとめ (Faith and Brave - C++で遊ぼう) や 本の虫 が詳しい。余談だが、2010年になっても0xから1xに変わることはないようだ。C++の生みの親であるBjarne Stroustrupによれば、 0xのxは16進数と考えて欲しい とのこと。 C++0xでは多くの便利な機能が追加されているが、今回はラムダ式、auto型指定子、初期化子リストを使用した。使用したコンパイラは Visual Studio 2010 beta2 (VS2010) 、 GCC 4.5 (snapshot 20100114) 、 Intel Compiler 11.1 だが、残念ながらVS2010とIntel Compiler 11.1では初期化子リストの機能が実装されていないので、それらを使用する場合、今回示したコードのうち素数を格納する静的メンバコンテナの初期化部分は変更しなくてはならない。また、GCCについては最新リリース版である4.4.2で ラムダ式が使えない ので、ここでは開発版の4.5を使用している。その他のコンパイラも含めてC++0xへの対応は C++0xCompilerSupport で確認できる。 初期化子リストを使用しない: static int init_primes[] = { 2, 3, 5, 7 }; vector<int> IsPrime::primes(init_primes, init_primes + 4); 初期化子リストを使用する: vector<int> IsPrime::primes = { 2, 3, 5, 7 }; 今回の...

整数の分割(integer partitions) を列挙するプログラムを書いてみた。ここでは再帰を用いた比較的単純な実装を示すが、より効率的なアルゴリズムを考えるとなると奥が深い。 整数の分割とは、ある整数を自然数の和で表したものだ。例えば分割する整数を5とすると求める自然数は以下のようになる。 5 4 1 3 2 3 1 1 2 2 1 2 1 1 1 1 1 1 1 1 今回のC++で書いたコードは単純で理解しやすいと思うが、高速化などは考えていない。もしよりよいアルゴリズムに興味があるのなら Antoine Zoghbiu and Ivan Stojmenovic. Fast Algorithms for Generating Integer Partitions. Intern. J. Computer Math. , 70 , 319-332 (1998) (PDF)あたりが参考になるかも。 使い方は以下の通り。 partitions 分割したい整数 [ 分割した要素の最大値 ] 分割したい整数が5であれば、 partitions 5 とする。出力は次のようになる。 5 (max 5): 7 5 4 1 3 2 3 1 1 2 2 1 2 1 1 1 1 1 1 1 1 また、分割したい整数が5で、分割した要素の最大値を3とする場合は、 partitions 5 3 で、以下のような出力を得る。 5 (max 3): 5 3 2 3 1 1 2 2 1 2 1 1 1 1 1 1 1 1 partitions.cpp #include <iostream> #include <vector> #include <algorithm> #include <cstdlib> using namespace std; void partitions(int n, int max_n, vector<int>& v, vector<vector<int> >& vv) { if (n == 0) { vv.push_back(v); return; ...

「 騎士の巡歴(Knight's Tour) 」を解くプログラムをC++で書いてみた。騎士の巡歴とはチェスボード上の駒「ナイト」を移動させてすべてのマスを一回ずつ通過させる問題だ。特に最初の駒の位置に戻ってくる解を「騎士の周遊(Closed Knight's Tour)」と呼ぶ。 一応、C++らしいコードを心掛けてみたけど、思ったより長くなってしまった。騎士の巡歴では Warnsdorffのアルゴリズム 、騎士の周遊においてはそれに加えて Schwenkの定理 を用いている。また、 Puzzle DE Programming で書かれている方法も参考にさせてもらった。 まず、ボード上のマスをNode構造体のポインタで表現し、 *node とする。その node からナイトが動けるマスを vector<Node*> next; として保持する。また、一度通った node は node->visited にそれを記録しておく。 next を順に巡っていくが、Warndorffのアルゴリズムにより、移動先となる next 内の node のうち、動けるマスの最も少ない node を移動先とする。すべてのマスを通れば探索終了だ。 また、騎士の周遊については、最初にSchwenkの定理から解があるかどうかを判断する。解があれば、騎士の巡歴と同様の手順でまず解を見つける。ここで、見つけた解の最初のマスからナイトが動けるマスをA、最後のマスからナイトが動けるマスをBとし、BからAの順となるマスがあるかどうか調べる。もしあれば、そのAのマスから最後のマスまでの順路を逆にすればそれが騎士の周遊の解となる。もし、BからAの順となるマスがなければ、動けるマスが見つかるまでバックトラックにより探索する。この手順については上記の Puzzle DE Programming が詳しい。 プログラムの使い方は以下の通り。 knight 行数 列数 [駒のスタート行=0 駒のスタート列=0 巡歴か周遊か?] スタート位置はデフォルトで(0, 0)で、何も指定がなければ騎士の巡歴を求める。騎士の周遊の場合は、コマンドラインの最後に c を付ける。 例えば、16x16マス上で騎士の周遊を解く場合は以下のようにする。 knight 16...

明けましておめでとうございます。 今年でこのブログも4年目になるなぁ。早いものだ。自分が初めてインターネットに触れたのが1994年で、それ以前はNIFTY-ServeやPC-VAN、草の根BBSなど利用していた。インターネットを利用し始めたころはメールアドレスを持っていても周りが誰も使っていなかったので全くの無用の長物だった。その頃、みんなが電子メールを持つようになれば素晴らしい世界になるのにと考えたものだが、それから数年でその「素晴らしい世界」になったわけだ。 自分のウェブサイトを作ったのが1996年から、ブログを書き始めたのは2000年2月14日からだけど、最初の頃はHTMLを直接編集していた。その後、Movable Typeに移行してしばらく使っていた。それからいくつかのサイトやブログなどに手を出し、科学情報ばかりを扱った「科学随想録」や以前ハマったMMOのEverQuestのギルド用フォーラムを作成したりしつつ、今の「良いもの。悪いもの。」に落ち着いた。 このブログはプログラミング関連情報をメインとして扱っているけど、昨年は、Google Wave、Android、Twitter、mixiアプリなどが興味深かったな。今年はどんな技術が出てくるんだろう。今からワクワクする。ブログ記事は自分が楽しめるものを書く、書きたい時が書き時というスタンスだけど、他の人にも楽しんでもらえるとしたら嬉しい。因みに自分のプログラミングスタイルの変遷は「 プログラミングができるということは、人生を楽にできるということ 」に書いてある。 これからもどうぞよろしくお願いします。

前回 に引き続き今回も Processing による日本全国コンビニ店舗分布地図について書いてみる。前回は大きさが固定された地図で拡大・縮小や移動ができなかったので、高解像度インタラクティブ版として、それをできるようにしてみた。 まず、日本地図がPNG画像なのでこれをベクタ画像に変更する。フォーマットはSVGだ。前回と同様に「 カビパン男と私 」で提供されている日本地図のSVGファイルを若干加工して利用させてもらった。 次にこのSVGデータをProcessingで利用する方法について述べる。「 ビジュアライジング・データ 」にもSVGの利用方法(processing.candy.*)が書いてあるのだが、実はここに書いてある情報は古くて現在では使用できない。現在では、PShapeを使ってSVGを利用する。 PShape mapShape = loadShape("japan.svg"); こんな感じだ。画面の(x, y)座標に表示するにはdraw関数内で以下のように記述すればいい。 shape(mapShape, x, y); 詳しくは Processingのリファレンス に書かれているので、下記のソースコードなどを参考にして調べてみて欲しい。 さて、これでSVG画像が利用できるようになったので、今度はこれを使って、 マウスの左ボタンのドラッグで平行移動 を行えるようにしてみる。やり方としてはmouseDragged関数を使う。移動したピクセル分だけ画像も移動する。mouseDragged関数内で移動分のオフセットを変数に入れるようにし、上述のshape関数でオフセット分だけ移動した位置に描画すればいい。 次に画像の拡大・縮小を行う。 マウスの右ボタンのドラッグで操作 する。 マウスを上に移動すれば拡大、下に移動すれば縮小 だ。SVG画像の縮小・拡大はscaleメソッドを使う。例えば、上記の例で元の画像を1.5倍にしたければ以下のようにする。 mapShape.scale(1.5); scaleメソッドによる拡大・縮小は画像の位置やそのオフセットなどもちゃんと正しい位置になるようにしなければならないのだが、それらについてはソースコードを参考にして欲しい。 上記はSVG画像による日本...

「ビジュアライジング・データ Processingによる情報視覚化手法」を読んでみたのだが、 Processing の情報処理とその視覚化がとても興味深いものだった。そこで、日本地図を使って何らかの情報を視覚化したいと考え、それならば全国に存在するコンビニ店舗の分布を調べてみようと思い、Processingで作ってみることにした。 因みに以前にも「 Processingで分子動力学計算 」や「 Processingを使ってWebカメラを監視カメラにする 」などのブログ記事を書いているように、プログラミング言語Processingはビジュアル関連で広範に応用でき、そのポテンシャルはとても大きい。 まずは、コンビニの店舗情報が必要なので、gooの コンビニ店舗検索 を利用し、ここから住所情報を得ることにした。ただ、一度に検索できる件数は5,000件までなので5,000件以上存在するコンビニに対してはPythonで都道府県ごとに検索してマージした。データはHTMLで書かれているのでPythonのre.findallのパターンマッチを使って一気に必要なデータに変換した。 次に、Google Mapsの HTTPリクエスト経由のジオコーディング を使って住所情報を経度・緯度に変換した。この辺もPythonを使えばちょちょいのちょいだ。以下のAPIを利用する。 http://maps.google.com/maps/geo?q=住所情報&output=json&sensor=false&key=APIキー ただし、変換は1日につき15,000件の制限が掛かっているので4万件以上ある住所情報を取得するには3日かかることになる。もっともIPアドレスによる制限なので別のIPを使えば問題ないようだ。 日本地図上にデータを表示するためには、経度・緯度情報からXY座標に変換する必要がある。日本地図は国土地理院によるとユニバーサル横メルカトル(UTM)図法が使われているらしい。簡単に言うと世界地図でよく使われているメルカトル図法の赤道に合わせている中心線を6°ごとの経度線に合わせて作成する方法らしい。そのための変換関数をPythonで以下のように作成した。本当は縮尺の微調整が入るようなのだが、今回はそこまで厳密にする必要はないのでその辺は省い...

以前に、 はじめてのCUDAプログラミングで分子動力学計算 というブログ記事を書いたことがある。最近このなんちゃって分子動力学計算(MD)プログラムのソースコードを読み直してみたのだが、かなりひどい。しかも、「CUDAプログラミング」でググってみると、この記事が2番目にくる。こんないい加減なコードを参考にされたら読んだ人にも迷惑が掛かるので修正することにした。それほどCUDAに慣れているわけではないが、前回のコードよりはましだと思う。 それにしてもコードの内容がひどい。意味もなく__syncthreads()が入っているし、複数のスレッドから同じグローバルメモリに書き込みしてるし、レジスタを活用してないし、計算の順序は非効率だし、ダメダメだ。 そこでまず、CUDA Visual Profilerで関数のパフォーマンス測定を行った。因みにこのプロファイラは CUDAプログラミングツールキット に含まれているもので、CUDAでプログラミングを行うには必須だと思う。使い方は簡単で、Windowsであればcudaprof.exeを実行して、FileメニューからNew...を選び、プロジェクトの名前と場所を設定して、実行するCUDAプログラムを指定するだけだ。引数が必要ならそれも指定しておく。デフォルトでは4回実行され解析される。それぞれの関数がどれだけ時間がかかったのかプロットされるので、どこがボトルネックになっているのか一目瞭然だ。 まあ、予想していた通り、Calc関数が全体の99.8%を占めていたので、ここから修正を行った。まず、iとjを使っていたループをjだけにした。iはスレッド数で分割されていたが、ブロックと合わせてそれを消した。次に、jのループ内で不必要にグローバルメモリにアクセスしないようにした。例えばchg[i]というグローバル変数はループ外でローカル変数に入れてそれを使うようにする。あとは、できるだけ除算などの演算を減らすようにした。 次に、系全体を中心に戻す関数があるのだが、中心座標を求める部分はホスト側の関数で実装し、系全体を戻す関数のみをCUDAの関数とした。もともとこの処理は毎回やらなくてもよいものなので、指定したステップ毎に1回行うようにした。今回は100ステップに1回としている。なので、全体の処理時間から見れば無視しても良く、頑張...