2010年1月22日金曜日

C++0x: 構造体を格納したSTLコンテナに対してソート・探索・削除などのアルゴリズムを適用する

このブログ記事をはてなブックマークに追加

先日、構造体を格納したコンテナに対して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 };

今回のコードでは関数オブジェクトを多用していたのでラムダ式が非常に有効だった。ラムダ式を使えば関数オブジェクトを用意しなくてもその場ですぐに記述できるからだ。

以下のコードは関数オブジェクトを用意してsort関数を使用した例だ。

class LessAn { public: bool operator()(const A& a, const A& b) { return a.n < b.n; } }; sort(A_list.begin(), A_list.end(), LessAn());

これに対しラムダ式を使えば、以下のように一行で済ませることができる。

sort(A_list.begin(), A_list.end(), [](const A& a, const A& b) { return a.n < b.n; });

このように、C++0xは機能の一部を使うだけでもとても便利だ。規格が新しくなるたびに新たに覚えなくてはならないことが増えるので、それを面倒と見る向きもあるかもしれないが、そのコストに見合うだけの価値がC++0xにはあると思う。

最後にC++0xで書き換えたソースコードを示す。クラスによる関数オブジェクトがごっそり無くなっていることに気付くだろう。前回からの変更部分は赤で示している。変更前のコードや実行時の出力結果については前回の記事を参考にして欲しい。

#include <iostream> #include <iomanip> #include <vector> #include <algorithm> #include <numeric> #include <cstdlib> using namespace std; const int num_A = 10; // A_listコンテナ内の構造体Aの数. const int max_n = 30; // 構造体A内の配列の最大数. // 構造体A : 要素数n, 要素の配列pを持つ. struct A { int n; int *p; A(int n_, int *p_) : n(n_), p(p_) { for (int i = 0; i < n; i++) p[i] = rand() % max_n; } }; // 構造体Aの表示. void print(const A& a) { cout << setw(3) << a.n << " :"; for (int i = 0; i < a.n; i++) cout << " " << a.p[i]; cout << endl; } // A_listコンテナの表示. void print(const vector<A>& a_list) { cout << " n : p[]" << endl; for (auto p = a_list.begin(); p != a_list.end(); ++p) print(*p); cout << endl; } // 任意の数値を素数判定する関数オブジェクト. class IsPrime { private: static vector<int> primes; int i, j; public: bool operator()(int n) { if (n < 2) return false; if (binary_search(primes.begin(), primes.end(), n)) return true; for (auto p = primes.begin(); p != primes.end(); ++p) { if (n < (*p) * (*p)) return true; if (n % *p == 0) return false; } for (i = *(primes.end() - 1) + 2; i * i <= n; i += 2) { for (j = 0; i > primes[j] * primes[j] && i % primes[j] != 0; j++); if (i < primes[j] * primes[j]) { primes.push_back(i); if (n % i == 0) return false; } } return true; } }; // 素数を格納する静的メンバコンテナ. vector<int> IsPrime::primes = { 2, 3, 5, 7 }; int main() { vector<A> A_list; // 構造体Aのvectorコンテナ A_list. // A_listを乱数で初期化. for (int i = 0; i < num_A; i++) { int n = rand() % max_n + 1; A_list.push_back(A(n, new int[n])); } cout << "初期化した構造体Aのコンテナ(A_list):" << endl; print(A_list); cout << "構造体Aの配列(A.p)をソートして重複データを削除:" << endl; for_each(A_list.begin(), A_list.end(), [](A& a) { sort(a.p, a.p + a.n); a.n = unique(a.p, a.p + a.n) - a.p; }); print(A_list); cout << "構造体Aの配列(A.p)の要素数(A.n)でA_listをソート:" << endl; sort(A_list.begin(), A_list.end(), [](const A& a, const A& b) { return a.n < b.n; }); print(A_list); cout << "構造体Aの配列(A.p)を素数のみにする:" << endl; for_each(A_list.begin(), A_list.end(), [](A& a) { a.n = stable_partition(a.p, a.p + a.n, IsPrime()) - a.p; }); print(A_list); cout << "構造体Aの配列要素の和が偶数ならA_listから削除:" << endl; A_list.erase(remove_if(A_list.begin(), A_list.end(), [](const A& a) { return accumulate(a.p, a.p + a.n, 0) % 2 == 0; }), A_list.end()); print(A_list); cout << "再び構造体Aの配列(A.p)の要素数(A.n)でA_listをソート:" << endl; sort(A_list.begin(), A_list.end(), [](const A& a, const A& b) { return a.n < b.n; }); print(A_list); cout << "構造体Aの配列(A.p)内で指定した数値と一致する要素を持つ構造体Aを探す:" << endl; for (int i = 0; i < 10; i++) { auto p = find_if(A_list.begin(), A_list.end(), [i](const A& a) { return find(a.p, a.p + a.n, i) != a.p + a.n; }); cout << setw(3) << i << " : "; if (p != A_list.end()) print(*p); else cout << "指定の数値を要素として持つ構造体Aは見つかりませんでした." << endl; } return 0; }

続きを読む...

2010年1月19日火曜日

整数の分割

このブログ記事をはてなブックマークに追加

整数の分割(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; } for (int i = min(n, max_n); i > 0; i--) { v.push_back(i); partitions(n - i, i, v, vv); v.pop_back(); } } int main(int argc, char* argv[]) { if (argc < 2) { cerr << "Usage: " << argv[0] << " number [max_number_in_partitions]" << endl; exit(1); } int n = atoi(argv[1]); int max_n = argc > 2 ? atoi(argv[2]) : n; vector<vector<int> > vv; partitions(n, max_n, vector<int>(), vv); cout << n << " (max " << max_n << "): " << vv.size() << endl; for (vector<vector<int> >::const_iterator p = vv.begin(); p != vv.end(); ++p) { for (vector<int>::const_iterator q = p->begin(); q != p->end(); ++q) cout << " " << *q; cout << endl; } return 0; }

余談だが、ActiveState CodeでPythonのジェネレータを使った以下のようなコードを見つけた。さらに同じ作者がその2年後に書いた記事から、別アルゴリズムのソースコード(tarball)を読むことができる。

def partitions(n): # base case of recursion: zero is the sum of the empty list if n == 0: yield [] return # modify partitions of n-1 to form partitions of n for p in partitions(n-1): yield [1] + p if p and (len(p) < 2 or p[1] > p[0]): yield [p[0] + 1] + p[1:]

続きを読む...

2010年1月12日火曜日

C++: 騎士の巡歴と周遊

このブログ記事をはてなブックマークに追加

騎士の巡歴(Knight's Tour)」を解くプログラムをC++で書いてみた。騎士の巡歴とはチェスボード上の駒「ナイト」を移動させてすべてのマスを一回ずつ通過させる問題だ。特に最初の駒の位置に戻ってくる解を「騎士の周遊(Closed Knight's Tour)」と呼ぶ。

一応、C++らしいコードを心掛けてみたけど、思ったより長くなってしまった。騎士の巡歴ではWarnsdorffのアルゴリズム、騎士の周遊においてはそれに加えてSchwenkの定理を用いている。また、Puzzle DE Programmingで書かれている方法も参考にさせてもらった。

まず、ボード上のマスをNode構造体のポインタで表現し、*nodeとする。そのnodeからナイトが動けるマスをvector<Node*> next;として保持する。また、一度通ったnodenode->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 16 c

出力結果:

1 254 31 232 205 196 29 194 203 106 27 104 187 108 25 102 32 229 256 253 30 235 204 197 28 193 202 107 26 103 180 109 255 2 231 206 233 252 195 238 201 198 105 186 181 188 101 24 230 33 228 219 236 213 234 251 192 247 190 199 112 179 110 177 3 224 209 214 207 220 237 246 239 200 185 182 189 176 23 100 34 227 218 223 210 241 212 167 250 191 248 113 184 111 178 67 225 4 215 208 217 166 221 240 245 160 183 172 175 68 99 22 128 35 226 165 222 211 242 161 168 249 170 149 114 173 66 69 5 164 129 216 143 162 157 244 159 148 115 174 171 98 21 64 36 127 136 163 156 243 144 141 152 169 150 93 116 65 70 55 133 6 155 130 137 142 153 158 147 92 117 80 97 56 63 20 126 37 132 135 154 145 140 87 118 151 94 73 62 71 54 57 7 134 125 138 131 86 119 146 91 74 81 96 79 58 19 50 38 123 40 85 120 139 88 45 82 95 76 61 72 51 16 53 41 8 121 124 43 10 83 90 75 12 47 78 59 14 49 18 122 39 42 9 84 89 44 11 46 77 60 13 48 17 52 15



24x12マスで(10, 5)の位置から騎士の巡歴を解く場合は以下のようにする。

knight 24 12 10 5

出力結果:

24 27 106 225 272 29 108 231 250 31 110 113 105 226 25 28 107 230 267 30 109 112 249 32 26 23 224 271 268 273 232 251 266 247 114 111 211 104 227 288 229 270 259 276 233 252 33 248 22 223 210 269 280 285 274 265 258 263 246 115 103 212 281 228 287 260 277 262 275 234 253 34 214 21 222 209 284 279 286 257 264 245 116 235 193 102 213 282 219 256 261 278 241 254 35 244 20 215 194 221 208 283 240 255 204 243 236 117 101 192 153 218 239 220 207 242 237 200 205 36 152 19 216 195 190 1 238 203 206 175 118 201 159 100 191 154 217 196 189 178 199 202 37 124 18 151 158 187 172 179 2 197 174 125 176 119 99 160 155 184 157 188 173 168 177 198 123 38 150 17 148 161 186 171 180 3 126 169 120 129 95 98 185 156 183 146 167 170 165 128 39 122 16 149 96 147 162 181 164 127 4 121 130 69 97 94 137 182 135 142 145 166 131 70 5 40 138 15 92 141 144 163 134 73 66 57 68 59 93 88 139 136 83 74 143 132 71 60 41 6 14 81 84 91 140 133 72 65 56 67 58 49 87 78 89 82 75 64 55 46 61 50 7 42 80 13 76 85 90 11 62 53 44 9 48 51 77 86 79 12 63 54 45 10 47 52 43 8



示した画像は出力した解をPython+PILで画像に変換したものだ。以下のように使用する。

show_knight.py 解のデータファイル 出力画像ファイル [幅=400 高さ=400]

上記の例では以下のようにする。

show_knight.py knight_16x16.txt knight_16x16.png
show_knight.py knight_24x12.txt knight_24x12.png 300 600

また、今回のプログラムでは200x200マスの騎士の周遊でもほぼ瞬時に解くことができた。先頭の画像はその解だ。多くの場合は問題なく解を見つけることができるが、場合によっては困難なこともある。そのようなときはスタート位置を変更することで解けることがある。

以下にプログラムのソースコードを示す。

knight.cpp

// Solver of Knight's Tour using Warnsdorff's algorithm and Schwenk's theorem #include <iostream> #include <iomanip> #include <vector> #include <algorithm> #include <utility> #include <cmath> #include <cstdlib> using namespace std; struct Node { int row, col; bool visited; vector<Node*> next; Node(int r, int c) : row(r), col(c), visited(false) { } }; class IsUnvisited { public: bool operator()(const Node* a) { return !a->visited; } }; class IsVisited { public: bool operator()(const Node* a) { return a->visited; } }; class NotEqualUnvisited { private: const int cnt; public: NotEqualUnvisited(const Node* a) : cnt(count_if(a->next.begin(), a->next.end(), IsUnvisited())) { } bool operator()(const Node* a) { return count_if(a->next.begin(), a->next.end(), IsUnvisited()) != cnt; } }; class LessMovable { public: bool operator()(const Node* a, const Node* b) { int cnt_a = count_if(a->next.begin(), a->next.end(), IsUnvisited()); int cnt_b = count_if(b->next.begin(), b->next.end(), IsUnvisited()); return cnt_a < cnt_b; } }; class KnightTour { private: int nrow, ncol; vector<pair<int, int> > moves; vector<Node*> nodes, tour, best_tour; bool is_closed; public: KnightTour(int r, int c, bool closed) : nrow(r), ncol(c), is_closed(closed) { // Knight can move two horizontally and one vertically, // or one horizontally and two vertically. moves.push_back(make_pair( 2, 1)); moves.push_back(make_pair( 1, 2)); moves.push_back(make_pair( 2, -1)); moves.push_back(make_pair( 1, -2)); moves.push_back(make_pair(-2, 1)); moves.push_back(make_pair(-1, 2)); moves.push_back(make_pair(-2, -1)); moves.push_back(make_pair(-1, -2)); } void search(Node* node); void run(int start_pos); void print(); }; void KnightTour::search(Node* node) { if (node->visited) return; if (best_tour.size() == nrow * ncol) { // for Closed Knight's Tour if (is_closed) { if (find(best_tour.back()->next.begin(), best_tour.back()->next.end(), best_tour.front()) != best_tour.back()->next.end()) return; for (vector<Node*>::iterator p = best_tour.back()->next.begin(); p != best_tour.back()->next.end(); ++p) { vector<Node*>::iterator q = find(best_tour.begin(), best_tour.end(), *p) + 1; vector<Node*>::iterator r = find(best_tour.front()->next.begin(), best_tour.front()->next.end(), *q); if (r != best_tour.front()->next.end()) { reverse(q, best_tour.end()); return; } } best_tour.clear(); } return; } node->visited = true; tour.push_back(node); if (best_tour.size() < tour.size()) best_tour = tour; // Warnsdorff's algorithm sort(node->next.begin(), node->next.end(), LessMovable()); vector<Node*> next(node->next); next.erase(remove_if(next.begin(), next.end(), IsVisited()), next.end()); if (!next.empty()) next.erase(remove_if(next.begin(), next.end(), NotEqualUnvisited(next.front())), next.end()); for (vector<Node*>::iterator p = next.begin(); p != next.end(); ++p) search(*p); node->visited = false; tour.pop_back(); } void KnightTour::run(int start_pos) { // initialization for nodes for (int i = 0; i < nrow; i++) for (int j = 0; j < ncol; j++) nodes.push_back(new Node(i, j)); for (vector<Node*>::iterator p = nodes.begin(); p != nodes.end(); ++p) { for (vector<pair<int, int> >::iterator q = moves.begin(); q != moves.end(); ++q) { int r = (*p)->row + q->first; int c = (*p)->col + q->second; if (r >= 0 && r < nrow && c >= 0 && c < ncol) (*p)->next.push_back(nodes[r*ncol+c]); } } // Schwenk's theorem if (is_closed && (nrow * ncol % 2 == 1 || (min(nrow, ncol) == 2 || min(nrow, ncol) == 4) || (min(nrow, ncol) == 3 && (max(nrow, ncol) == 4 || max(nrow, ncol) == 6 || max(nrow, ncol) == 8)))) return; if (!is_closed && nrow * ncol % 2 == 1 && start_pos % 2 == 1) return; search(nodes[start_pos]); } void KnightTour::print() { if (best_tour.size() < nrow * ncol) { cout << "No solution." << endl; return; } vector<vector<int> > board(nrow, vector<int>(ncol)); int cnt = 1; for (vector<Node*>::iterator p = best_tour.begin(); p != best_tour.end(); ++p) board[(*p)->row][(*p)->col] = cnt++; int width = static_cast<int>(log10(static_cast<double>(nrow * ncol)) + 2); for (int i = 0; i < nrow; i++) { for (int j = 0; j < ncol; j++) cout << setw(width) << board[i][j]; cout << endl; } } int main(int argc, char* argv[]) { if (argc <= 2) { cerr << "Usage: " << argv[0] << " nrow ncol [start_row=0 start_col=0] [closed]" << endl; exit(1); } int nrow = atoi(argv[1]); int ncol = atoi(argv[2]); int r = 0, c = 0; if (argc > 4) { r = atoi(argv[3]); c = atoi(argv[4]); } bool is_closed = false; if (argv[argc-1][0] == 'c') is_closed = true; KnightTour kt(nrow, ncol, is_closed); kt.run(r * ncol + c); kt.print(); return 0; }

show_knight.py

#!/usr/bin/env python import sys, os, Image, ImageDraw def draw_board(board, nrow, ncol, size, offset_x, offset_y, out_image): im = Image.new("RGB", size) im.paste((255, 255, 255)) draw = ImageDraw.Draw(im) for i in range(0, size[0], offset_x): draw.line((i, 0, i, size[1]), fill=0) for i in range(0, size[1], offset_y): draw.line((0, i, size[0], i), fill=0) for i in range(0, size[0], offset_x): for j in range(0, size[1], offset_y): if (i // offset_x + j // offset_y) % 2 == 1: im.paste((190, 190, 190), (i + 1, j + 1, i + offset_x, j + offset_y)) for i in range(nrow * ncol - 1): draw.line((board[i] % ncol * offset_x + offset_x // 2, board[i] // ncol * offset_y + offset_y // 2, board[i+1] % ncol * offset_x + offset_x // 2, board[i+1] // ncol * offset_y + offset_y // 2), fill=(255, 0, 0)) if sorted([abs(board[nrow*ncol-1] % ncol - board[0] % ncol), abs(board[nrow*ncol-1] // ncol - board[0] // ncol)]) == [1, 2]: draw.line((board[nrow*ncol-1] % ncol * offset_x + offset_x // 2, board[nrow*ncol-1] // ncol * offset_y + offset_y // 2, board[0] % ncol * offset_x + offset_x // 2, board[0] // ncol * offset_y + offset_y // 2), fill=(255, 0, 0)) im.save(out_image) def main(args): if len(args) < 3: print >>sys.stderr, "Usage: %s in_datafile out_imagefile [size_x=400 size_y=400]" % os.path.basename(args[0]) sys.exit(1) in_file, out_image = args[1:3] if len(args) < 5: size = (400, 400) else: size = map(int, args[3:5]) data = [map(int, l.strip().split()) for l in file(in_file)] nrow, ncol = len(data), len(data[0]) board = range(nrow * ncol) cnt = 0 for r in data: for c in r: board[c-1] = cnt cnt += 1 offset_x, offset_y = size[0] // ncol, size[1] // nrow size = (offset_x * ncol + 1, offset_y * nrow + 1) draw_board(board, nrow, ncol, size, offset_x, offset_y, out_image) if __name__ == "__main__": main(sys.argv)

続きを読む...

2010年1月1日金曜日

2010年を迎えて

このブログ記事をはてなブックマークに追加

明けましておめでとうございます。

今年でこのブログも4年目になるなぁ。早いものだ。自分が初めてインターネットに触れたのが1994年で、それ以前はNIFTY-ServeやPC-VAN、草の根BBSなど利用していた。インターネットを利用し始めたころはメールアドレスを持っていても周りが誰も使っていなかったので全くの無用の長物だった。その頃、みんなが電子メールを持つようになれば素晴らしい世界になるのにと考えたものだが、それから数年でその「素晴らしい世界」になったわけだ。

自分のウェブサイトを作ったのが1996年から、ブログを書き始めたのは2000年2月14日からだけど、最初の頃はHTMLを直接編集していた。その後、Movable Typeに移行してしばらく使っていた。それからいくつかのサイトやブログなどに手を出し、科学情報ばかりを扱った「科学随想録」や以前ハマったMMOのEverQuestのギルド用フォーラムを作成したりしつつ、今の「良いもの。悪いもの。」に落ち着いた。

このブログはプログラミング関連情報をメインとして扱っているけど、昨年は、Google Wave、Android、Twitter、mixiアプリなどが興味深かったな。今年はどんな技術が出てくるんだろう。今からワクワクする。ブログ記事は自分が楽しめるものを書く、書きたい時が書き時というスタンスだけど、他の人にも楽しんでもらえるとしたら嬉しい。因みに自分のプログラミングスタイルの変遷は「プログラミングができるということは、人生を楽にできるということ」に書いてある。

これからもどうぞよろしくお願いします。

続きを読む...