Week 6 (2021/11/10)¶
今日やること
- 構造体
- ビット演算
- ライブラリツアー
構造体¶
今日は構造体をまず見ていきましょう。これまでは変数や配列は単一の型に関する話でした。 これに対し、様々な型の変数をまとめて一つのまとまりにしたものを構造体と呼びます。
構造体の基本¶
構造体の基本のコードを見てみましょう。ここでは二次元空間の一点を表すデータ表現を考え、それを構造体で表します。
#include <stdio.h>
struct point {
int x;
int y;
};
int main () {
struct point pt;
pt.x = 10;
pt.y = 30;
printf("%d %d\n", pt.x, pt.y); // 10 30
}
struct pointという書き方で、新たに構造体を宣言できます。
ここでpointという部分は構造体タグと呼ばれます。構造体タグは構造体の名前です。
構造体の中で命名された変数はメンバと呼ばれます。
ここではint xとint yをメンバとして持ちます。
この宣言の部分はあくまで宣言だけなので、メモリを使いません。
構造体は、struct point ptといった表記で、一つ生成することができます。
ここでは、
新たな型 struct pointの変数ptを一個作る、という意味です。
ptのメンバ変数は、ドットで読み書きすることができます。
コラム
このように、構造体宣言によって定義された変数がptとして作られたとき、ptのことをインスタンス（実体）と呼びます。
インスタンスは、実際にメモリを消費します。これは構造体によらない話で、型によって作られる「もの」をインスタンスと言います。
さて、いくつかの表記を勉強しましょう。 以下のように、構造体は宣言と同時にインスタンスを確保することもできます。しかし、この表記は使う機会はあまりないかもしれません。
// 以下の表記で、(1) 構造体の宣言と(2) インスタンスptの作成 を同時に行う
struct point {
int x;
int y;
} pt;
インスタンスの初期化は、初期化の演算子を使って次のようにも書けます。変数の順番は、宣言のときと同じです。
struct point pt = {10, 30};
// 上は、以下と同じ
// struct point pt;
// pt.x = 10;
// pt.y = 30;
また、C99からは次のようにメンバ変数を指定して初期化する表記が可能になりました。
struct point pt = {.x = 10, .y = 30};
struct point pt = {.y = 30, .x = 10};
次のように、初期化以外の部分では初期化構文を使えません。その場合、型名を最初にカッコで書くと 可能になります。この機能を複合リテラルと言います。
struct point pt;
// pt = {.x = 10, .y = 30}; // これはダメ
pt = (struct point){.x = 10, .y = 30}; // これはOK
いちいちstruct pointと二単語書くのがめんどくさいときは、typedefという機能をつかって、
言い換えることができます。次の例では構造体の宣言時にstruct pointをPointと言い換えています。
これにより、Pointを新たな型だと思って、Point ptのようにより直感的にインスタンスを作れます。
struct point {
int x;
int y;
};
// このようにtypedefすることで、
// 「struct point」と書く代わりに「Point」と書けるようになる
typedef struct point Point;
int main() {
// 「Point」と書ける
Point pt = {.x = 10, .y = 30};
}
また、typedefに関して、以下のようにまとめて書くこともできます。
// (1) pointの定義と(2) typedef を一つにまとめた書き方
typedef struct point {
int x;
int y;
} Point;
int main() {
// 「Point」と書ける
Point pt = {.x = 10, .y = 30};
}
コラム
typedef 旧型名 新型名 の表記により、構造体に限らず型名を新しく作ることができます。例えば
typedef int JPY; // 今後、JPYはintと同じ意味
JPY salary = 100;
printf("my salary: %d\n", salary); // 100
構造体は別の構造体に代入できます。この場合、メンバ変数は全てコピーされます。
また、==や!=を使うことはできません。
Point pt = {.x = 10, .y = 30};
Point pt2 = pt; // pt2.x = pt.x および pt2.y = pt.y が行われる
printf("%d %d\n", pt.x, pt.y); // 10 30
printf("%d %d\n", pt2.x, pt2.y); // 10 30
// 次のようなことはできない
// if (pt == pt2) { ... }
構造体の中に配列をもったり、構造体自身を持ったりすることもできます。Pointの構造体が定義されているとして、
Pointを二個、および名前をもった構造体Rectは以下のように作れます。
#define MAX_NAME_LEN 20
typedef struct rect {
Point pt1;
Point pt2;
char name[MAX_NAME_LEN + 1];
} Rect;
int main () {
Rect r = {.pt1 = {.x = 10, .y = 30},
.pt2 = {.x = 25, .y = 12},
.name = "sample_rect"};
printf("r.pt1 = %d, %d\n", r.pt1.x, r.pt1.y); // 10 30
printf("r.pt2 = %d, %d\n", r.pt2.x, r.pt2.y); // 25 12
printf("r.name = %s\n", r.name); // sample_rect
}
r.pt1.xのようにドットを連続させることで、メンバの構造体のメンバにアクセスできます。
ここで、注意として、メンバの配列はchar name[MAX_NAME_LEN + 1]のように明示的に配列の長さを決定する必要があります。
char name[] のようにしてはダメです。
コラム
ここで構造体の非常に奇妙な特性を紹介します：構造体がメンバに配列をもつとする。このとき、構造体をコピー すると、メンバ配列の中身も全てコピーされます。 まず通常の配列を思い出しましょう。配列は、単純な代入でコピーすることができません。
int arr[3] = {1, 2, 3};
int arr2[3] = arr; // ダメ
一方で、構造体はコピー可能だと上で述べました。その場合、各メンバ変数がコピーされます。 ここで、構造体が配列をもつとします。ここで構造体をコピーすると、なんと 配列の中身が全て自動でコピーされるのです。 これはポインタが関係するのではなく、全てコピーです。 なので、次の例が示すように、別の構造体にコピーしたあとは、内部の配列はそれぞれ別々です。
typedef struct a {
int arr[3];
} A;
int main () {
A a1 = {.arr = {1, 2, 3}};
A a2 = a1; // a1.arrの中身全てがa2.arrにコピーされる！
for (int i = 0; i < 3; ++i) {
printf("%d, %d\n", a1.arr[i], a2.arr[i]);
// 1, 1
// 2, 2
// 3, 3
}
a2.arr[0] = 100; // 値を編集
for (int i = 0; i < 3; ++i) {
printf("%d, %d\n", a1.arr[i], a2.arr[i]);
// 1, 100 編集がa2にだけ反映されている
// 2, 2
// 3, 3
}
}
構造体と関数¶
さて、構造体を関数で扱う方式を見てみましょう。構造体は通常の変数と同じように関数に渡すことができます。
上で述べたPointが定義されているとして、
void print_point(Point pt) {
printf("pt: (%d, %d)\n", pt.x, pt.y);
}
Point pt = {.x = 10, .y = 3};
print_point(pt); // pt: (10, 3)
また、構造体を返すこともできます。次の例は、二つのPointを受け取り、x, yのそれぞれで和をとり、結果を返すものです。
Point add_point(Point pt1, Point pt2) {
pt1.x += pt2.x;
pt1.y += pt2.y;
return pt1;
}
Point pt_a = {.x = 10, .y = 2};
Point pt_b = {.x = 3, .y = 5};
Point out = add_point(pt_a, pt_b);
print_point(pt_a); // pt: (10, 2)
print_point(pt_b); // pt: (3, 5)
print_point(out); // pt: (13, 7)
print_point(pt_a); // pt: (10, 2) // pt_aは変更されない
pt1をいくら変更しても、呼び出し元のpt_aは変更されません。
この「コピーで渡す」という特性には注意が必要です。例えば巨大な配列をメンバとしても持つ
BigArrayという構造体を考えます。
そしてそれを引数として受け取るdo_nothing関数を考えます。
typedef struct big_array {
double val[1000];
} BigArray;
void do_nothing(BigArray ba) {
printf("do_nothing\n");
}
BigArray tmp;
do_nothing(tmp); // 構造体のコピー、すなわち巨大配列のコピー
returnするときも同じです。構造体が巨大である場合、それを関数内で作ってreturnすると、
大がかりなコピーが発生します。
やってみよう（20分）
上記を写経してみましょう
構造体の配列¶
構造体を要素としてもつ配列を作ることも出来ます。これはintなどの配列を作るときと同様に出来るので、わかりやすいです。
例えばPointを3つ作りたいときは次のようにします。
Point points[3];
Point points[3] = {
{.x = 3, .y = 4},
{.x = 5, .y = 6},
{.x = 7, .y = 8}
};
for (int i = 0; i < 3; ++i) {
print_point(points[i]);
}
// pt: (3, 4)
// pt: (5, 6)
// pt: (7, 8)
コラム
これまで明示的に触れませんでしたが、もともとCでは配列の長さはコンパイル時に大きさが決定される必要がありました。 ですが、C99では可変長配列（Variable-Length Array; VLA）という機能が導入され、配列の長さを変数で与えたり出来るようになりました。 例えば次のコードです。
int a[3]; // 通常の配列
int N = 3;
int b[N]; // 可変長配列
int a[3] = {1, 2, 3}; // OK
int a2[] = {1, 2, 3}; // OK。上の単なる省略表記
int N = 3;
int b[N] = {1, 2, 3}; // ダメ
そして、このことは構造体配列でも同様で、次の表記は出来ません。
int N;
Point points[N] = { // ダメ。可変長配列なので初期化表記を使えない
{.x = 3, .y = 4},
{.x = 5, .y = 6},
{.x = 7, .y = 8}
};
構造体とポインタ¶
さて、構造体へのポインタも、通常の変数のポインタと同様に扱うことができます。
Point pt = {.x = 3, .y = 10};
Point *q;
q = &pt; // q はPoint型変数ptへのポインタ
(*q).xのようにします。カッコが必要なことに注意してください。
*q.xだと*(q.x)の意味になってしまいます。この「ポインタが指す構造体のメンバ変数へのアクセス」はよく使うため、特殊な記法が
用意されており、q->xと書けます。これをまとめると次のようになります。
printf("%d\n", (*q).x); // 3. メンバへのアクセス。
printf("%d\n", q->x); // 3. 上と全く同じ処理。「->」という特殊な記法
さて、構造体はポインタをメンバ変数として持つこともできます。 その場合、構造体をコピーしたときにコピーされるのはポインタつまりアドレスなので、注意が必要です。 まず下の例を見てみましょう。
#define MAX_LEN 6
typedef struct str_by_arr {
char s[MAX_LEN + 1];
} StrByArr;
StrByArrを考えます。その実体は配列になっています。
これに対し次の処理を考えます。
StrByArr x1 = {.s = "hoge"};
StrByArr x2 = x1;
strcpy(x1.s, "fuga");
printf("x1.s: %s\n", x1.s); // x1.s: fuga
printf("x2.s: %s\n", x2.s); // x2.s: hoge
x1を作り、その内容をx2にコピーします。
先ほど勉強した通り、ここでは配列内全てのコピーが発生するということに注意しましょう。
その後、x1の中身をhogeからfugaに変更します。
ここではx1とx2は全く違うものなので、当然、x2の中身は変更されません。
これを図示すると次のようになります。
一方で、この自作string構造体を次のようにポインタで実現するとします。
typedef struct str_by_ptr {
char *s;
} StrByPtr;
buffという配列に文字列を確保し、そこを指すy1を考えます。
char buff[] = "hoge";
StrByPtr y1 = {.s = buff};
StrByPtr y2 = y1;
strcpy(buff, "fuga");
printf("y1.s: %s\n", y1.s); // y1.s: fuga
printf("y2.s: %s\n", y2.s); // y2.s: fuga
StrByArrのときと同様にy1をy2にコピーします。しかし、
このときにコピーされているのは「ポインタ変数そのもの」です。これにより、y1もy2も
同じbuffを指すようになります。そのため、もしbuffを変更してしまうと、
y1もy2も変更することになります。これを図示すると次のようになります。
ここで実際にy1, y2のメンバのs（指し示しているアドレス）
を確認すると、同じであることがわかります。
printf("%p %p\n", y1.s, y2.s); // 0x7ffece6aaa03 0x7ffece6aaa03
やってみよう（20分）
上記を写経してみましょう
クイズ
- 上記で考えた
PointおよびRectが定義されているとし、次のコードを考えますこのとき、次の4つはそれぞれ何を意味するでしょうか？Rect r = {.pt1 = {.x = 10, .y = 30}, .pt2 = {.x = 25, .y = 12}, .name = "sample test"}; Rect *q = &r;r.pt1.x,q->pt1.x,(r.pt1).x,(q->pt1).x - 以下のように
Point ptとRect rが与えられたとき、ptがrの中に含まれれば1、そうでなければ0となるような処理を書いてください。 ここで、r.pt1は必ず左上（xもyも値が小さい）、r.pt2は必ず右下（xもyも値が大きい）とします。 例えば次の場合は含まれるので1が返ります。Rect r = {.pt1 = {1, 1}, .pt2 = {3, 4}}; Point pt = {.x = 2, .y = 2}; - 上記について、
r.pt1とr.pt2の条件が外れた、次の場合も対応してください。どちらの場合でも、結果は1になります。Rect r = {.pt1 = {3, 4}, .pt2 = {1, 1}};Rect r = {.pt1 = {1, 4}, .pt2 = {3, 1}};
答え
- 全て同じ。
10 - 例えば次
int ptinrect(Point p, Rect r) { // K&R, p158 // pがrに含まれれば1を、そうでなければ0を返す。 // rは標準形式である必要がある、 return p.x >= r.pt1.x && p.x < r.pt2.x && p.y >= r.pt1.y && p.y < r.pt2.y; } - 例えば、次のような標準形式変換関数を作る
そのうえで、次のように呼ぶ
int min(a, b) { if (a < b) { return a; } else { return b; } } int max(a, b) { if (a < b) { return b; } else { return a; } } Rect canonrect(Rect r) { // K&R, p169 // rを標準形式にして返す。標準形式とは、pt1は左上を指し、pt2は右下を指す形式 Rect tmp; tmp.pt1.x = min(r.pt1.x, r.pt2.x); tmp.pt1.y = min(r.pt1.y, r.pt2.y); tmp.pt2.x = max(r.pt1.x, r.pt2.x); tmp.pt2.y = max(r.pt1.y, r.pt2.y); return tmp; }printf("%d\n", ptinrect(pt, canonrect(r)));
ビット演算¶
さて、次はビット演算を考えましょう。プログラミングにおいて、全ての変数の
実体は0か1という数字が並んでいる列です。この0, 1に対して処理を行うことをビット演算と言います。
ビット演算はローレベルな処理です。普段のプログラミングで出会う機会は少ないかもしれませんが、速度やメモリ効率を追求するときに切っても切れないです。
ある種、C/C++らしい処理だと言えます。
また、ビット演算を知っていると、それを使わなければ書きづらいような処理を簡単に書けるときもあります。
以下では、簡単のためunsigned char (非負の8ビット整数。8桁の二進数)をたくさん使っていきます。
なぜunsignedなものにするかというと、signedな整数は最上位桁を負の数の表現に使うために、ビット演算には適さない場面があるからです。その点についても説明します。
基本的な演算¶
まずは準備として、unsighed charを8ビットの二進数として表示してくれる以下の関数を考えましょう。
void print_bit_uchar(unsigned char c) {
for (int i = 0; i < 8; ++i) {
if (c & 1 << (7 - i)) {
printf("1");
} else {
printf("0");
}
}
printf("\n");
}
unsigned charの変数を渡すと、その二進数表記を表示してくれます。
unsigned char a = 3;
printf("%u\n", a); // 3
print_bit_uchar(a); // 00000011
さて、ビット演算のための演算子を見ていきましょう。
unsigned char a = 3; // 00000011
unsigned char b = 6; // 00000110
unsigned char x;
// ビットごとのAND
x = a & b;
print_bit_uchar(x); // 00000010
printf("%u\n", x); // 2
// ビットごとのOR
x = a | b;
print_bit_uchar(x); // 00000111
printf("%u\n", x); // 7
// ビットごとのXOR
x = a ^ b;
print_bit_uchar(x); // 00000101
printf("%u\n", x); // 5
printfすると整数です。これを先ほどのビット表示関数に入れると、結果を2進数で可視化してくれます。
次に、ビット演算独特のシフト演算子、およびビット反転を見ていきましょう。
unsigned char a = 35; // 00100011
unsigned char x;
// 左シフト
x = a << 2; // 右側は0詰めされる
print_bit_uchar(x); // 10001100
printf("%u\n", x); // 140
// 右シフト
x = a >> 2; // unsignedの場合、左側は0詰めされる
print_bit_uchar(x); // 00001000
printf("%u\n", x); // 8
// 1の補数（ビット反転）
x = ~a; // 単項演算子
print_bit_uchar(x); // 11011100
printf("%u\n", x); // 220
また、「&」「|」「^」「<<」「>>」は全て「=」と組み合わせて代入演算子とすることができます。
a <<= 2; // a = a << 2 の略
print_bit_uchar(a); // 10001100
printf("%u\n", a); // 140
コラム
実はgccとclangではバイナリ定数表記が準備されています。
すなわち、先頭に0bをつけ、後ろに0, 1を並べることで、
二進数の整数を直接表記することができます。
例として、以下のaとbは全くおなじものです
unsigned char a = 12;
unsigned char b = 0b00001100;
print_bit_uchar(a); // 00001100
print_bit_uchar(b); // 00001100
2の補数¶
さて、ここではビット演算そのものからは少し話がそれるのですが、重要な概念である「2の補数」を勉強しましょう。
unsignedな整数は、二進数が直接表現されているため、ビット列と数字の対応は自明です。
ところが、singedな整数の場合、負の数はどのようなビット列として表されるのでしょうか？
現代の多くのコンピュータでは、「2の補数」という表現形式が採用されています。以下を見てみましょう。
unsigned char uch = 193;
char ch = -63;
printf("uch: %d, ch: %d\n", uch, ch); // uch: 193, ch: -63
print_bit_uchar(uch); // 11000001
print_bit_uchar(ch); // 11000001
ここで、chとuchはどちらも同じビット列になります。uchのほうは自明ですが、chのほうはなぜ-63になるのでしょうか？
ここで、2の補数の考え方においては、最上位のケタをマイナスにします。
11000001がunsigned charの場合- \(2^7 * 1 + 2^6 * 1 + 2^5 * 0 + 2^4* 0 + 2^3 * 0 + 2^2 * 0 + 2^1* 0 + 2^0 * 1 = 193\)
11000001がsigned charの場合- \(-2^7 * 1 + 2^6 * 1 + 2^5 * 0 + 2^4* 0 + 2^3 * 0 + 2^2 * 0 + 2^1* 0 + 2^0 * 1 = -63\)
すなわち、この場合一番最初の1は、unsignedの場合は\(2^7\)を表しますが、signedの場合は\(-2^7\)を表します。
このように、signedな変数では、2の補数という独特な形式で負の数を表します。
ちなみに、上の定義から自明なように、「全ての桁が1」の場合、unsignedにおいてはそれは最大の整数になります。一方、2の補数の場合、それは必ず-1になります。
11111111がunsigned charの場合- \(2^7 * 1 + 2^6 * 1 + 2^5 * 1 + 2^4* 1 + 2^3 * 1 + 2^2 * 1 + 2^1* 1 + 2^0 * 1 = 255\)
11111111がsigned charの場合- \(-2^7 * 1 + 2^6 * 1 + 2^5 * 1 + 2^4* 1 + 2^3 * 1 + 2^2 * 1 + 2^1* 1 + 2^0 * 1 = -1\)
また、上の定義から明らかなように、2の補数においては、
- 最小値：一番上の桁が1で、それ以外が0
0b10000000 = -128
- 最大値：一番上の桁が0で、それ以外が1
0b01111111 = 127
になります。
さて、なぜこのようなややこしい表現になっているのでしょうか？例えば、もっとシンプルに、一番上の桁が0ならば正、1ならば負、としてはダメだったのでしょうか？
実はその場合、0b00000000（+0）と0b10000000（-0）はどちらもゼロです。同じゼロに二つのビット列を割り振ってしまうため、効率が悪いのです。
やってみよう（20分）
ここまでをを写経してみましょう
クイズ
次の例を考えましょう。
int a = 345;
int z = a + ~a;
zは何になりますか？また、aが違う値のとき、zはどうなりますか？
上記の例はMIT Open Courseware, 6.172, Fall 2018, Performance Engineering of Software System, Prof. Charles Leiserson and Prof. Julian Shun, #3, Bit Hacksで紹介されていたものです。この講義はめちゃくちゃオススメです。
答え
aがなんであれ、zは必ず-1になります。なぜなら、
aと~aは、二進数表現すると、最上位桁も含めた全ての桁について、0,1が反転しています。- なので、
a + ~aの算術足し算は、2進数表現において桁上げが発生しません。 - よって、かならず全ビットが1になります。
- 上で説明した通り、2の補数においては全ビットが1の場合、その符号付整数は-1です。
複雑な操作¶
さて、それではビット演算を組み合わせることで、様々な操作を実現していきましょうKing, p512。まずは「i番目の位置にビットをセットする」方法を学びましょう。
// 00000001 を左側にiビットずらす。つまり右からiビット目を1にする
// 注意：一番右は0ビット目
int i = 3;
unsigned char mask = 1 << i;
print_bit_uchar(mask); // 00001000
unsigned char a = 0b00000101;
a = a | mask; // aのうち、maskの位置にビットを立てる
print_bit_uchar(a); // 00001101
aに対し、右からi番目のビットを立てることができます。
これは次のようにまとめて表記できます。
a |= 1 << i; // a のうち右からiビット目を1にする
次に、「右からiビット目を0にする」方式を学びましょう。
// 右からiビット目を1にしたものを反転。つまり右からiビット目が0でそれ以外1
int i = 3;
unsigned char mask = ~(1 << i);
print_bit_uchar(mask); // 11110111
unsigned char a = 0b00001110;
a = a & mask; // aのうち、0のところは必ず0。それ以外はmaskに等しい
print_bit_uchar(a); // 00000110
aの右からi番目をクリア（0にする）ことができます。これをまとめると下記です。
a &= ~(1 << i); // a のうち右からiビット目を0にする
次に、「i番目のビットが立っているかどうかの判定」を学びましょう。
// 右からiビット目を1
int i = 3;
unsigned char mask = 1 << i;
print_bit_uchar(mask); // 00001000
unsigned char a = 0b00001110;
// aに対し、maskの位置（右からiビット目）に1がたっている場合：
// そこは1になる。それ以外は0になる（maskが0なので）
// 結局、a & mask はmaskと同じ値、すなわち非0になるので、その真偽値は真になる
if (a & mask) {
printf("a's %d-bit is on\n", i);
}
unsigned char b = 0b00000110;
// bに対し、maskの位置（右からiビット目）に1がたっていない場合：
// そこは0になる。それ以外も0になる（maskが0なので）
// 結局、b & mask は0（偽）になる
if (b & mask) {
printf("a's %d-bit is on\n", i); // ここには到達しない
}
// aの右からiビット目に1が立っている場合、という意味
if (a & 1 << i) { ... }
コラム
あらためてprint_bit_ucharを見てみましょう。
for (int i = 0; i < 8; ++i) {
if (c & 1 << (7 - i)) {
printf("1");
} else {
printf("0");
}
}
1 << (7 - i)の部分は、00000001を左に7-iビット分だけシフトさせています。
つまり、i = 0のとき10000000、i=1のとき01000000、といった具合です。
つまり、左からi番目に1が立った数を作っています。
ちなみに、演算子の優先順位としてはc & 1 << 7 - iでもよいのですが、わかりやすさのためにカッコをつけています。
これとターゲットのcのビット毎ANDをとるということは、ターゲットの整数の
左からi番目が1なのかどうか調べている、ということです。
ビット演算の例¶
さて、ここまで色々学んできましたが、果たしてビット演算というのは具体的にどういうときに使うのでしょうか？ 以下にいくつかの例を紹介します。
フラグ管理¶
古典的な例は、フラグの管理です。ここでは、お店の麺類メニューをどう表現するかを考えます。 以下のように、各麺類をある数字として宣言します。
#define RAMEN 0b000000001
#define SOBA 0b000000010
#define UDON 0b000000100
#define PASTA 0b000001000
さて、下記を見てみましょう。
unsigned char menu = 0;
// SOBA bitを立てる
menu |= SOBA; // 00000010
// UDON bitを立てる
menu |= UDON; // 00000110
if (menu & SOBA) {
printf("We're selling SOBA\n");
}
&演算子で簡単にチェックできます。
この表記の利点はなんでしょうか？
- メモリ効率の良さ。上記に対しラーメンを売っているときは
int ramen_flag = 1のように別々に変数を割り当てることもできるのですが、 多くのメモリ領域を使ってしまいます。これを上記のようにフラグで表現すればメモリの節約になります。 - 状態をまとめて表現。上記のようなフラグ管理は状態をまとめて一つの整数で表現出来るため、便利です。
このように、フラグ管理のためにビット演算は使われてきました。 特に、メモリ量が十分でなかった時代によく使われていたと思います。 現代ではメモリ量を切り詰める目的でフラグ管理をすることは少ないかもしれませんが、覚えておくとよいです。
冪集合のイテレーション¶
さて、別の応用例も考えてみましょう。ある集合\(\mathcal{X}\)に対し、その冪集合（部分集合全てを集めた集合）\(2^\mathcal{X}\)を考えます。 例えば、\(\mathcal{X} = \{A, B, C\}\) としたとき、\(2^\mathcal{X} = \{ \emptyset, \{ A \}, \{ B \}, \{ C \}, \{ A, B \}, \{ A, C \}, \{ B, C \}, \{ A, B, C \} \}\)です。 さて、ここで、\(2^\mathcal{X}\)の要素全てに対して何か処理をしたいときを考えましょう。 それはプログラミングとしてどう記述すればよいでしょうか？実はこのように部分集合の列挙を簡潔に書くのは簡単ではありません。 ここでビット演算の登場です。整数の各桁を集合の要素だと考えます。ここで「整数に1を足す」という処理をよく考えると、それは 可能な組み合わせを列挙していることと同じだと気付きます。 よって、「冪集合の要素の列挙」は、「要素数分のビットに対し、インクリメントを行う」と同じです。これを元に列挙を行ったのが次です。
int n = 3; // 要素は3つ
// 0b00000000 <= b < 0b00001000 の区間をイテレーションする
for (unsigned char b = 0; b < (1 << n); ++b) {
printf("b: ");
print_bit_uchar(b);
printf("{");
if (b & 0b00000001) { // 右から0番目にビットが立っていたら「A」
printf("A, ");
}
if (b & 0b00000010) { // 右から1番目にビットが立っていたら「B」
printf("B, ");
}
if (b & 0b00000100) { // 右から2番目にビットが立っていたら「C」
printf("C, ");
}
printf("}\n");
}
b: 00000000
{}
b: 00000001
{A, }
b: 00000010
{B, }
b: 00000011
{A, B, }
b: 00000100
{C, }
b: 00000101
{A, C, }
b: 00000110
{B, C, }
b: 00000111
{A, B, C, }
XORスワップ¶
次に、XORスワップを体験してみましょう。 通常、変数をスワップするためには一次変数を作る必要があります。
unsigned int x = 3;
unsigned int y = 5;
unsigned int tmp = x;
x = y;
y = tmp;
printf("x: %u, y: %u\n", x, y); // x: 5, y: 3
ところが、実は以下のようにXOR演算を３回行うとxとyをスワップできます。
unsigned int x = 3;
unsigned int y = 5;
x = x ^ y;
y = x ^ y;
x = x ^ y;
printf("x: %u, y: %u\n", x, y); // x: 5, y: 3
これはなぜでしょうか？全て書き下してみるとわかります。xor演算は各桁に対して独立に適用されるので、ある桁についてのみ考えます。 ある桁は当然ビットなので1か0しか値を持たないので、全ての可能性を書き下せます。
x |
y |
x^y |
(x^y)^y |
(x^y)^((x^y)^y) |
|---|---|---|---|---|
| １行目 | ２行目。xに等しい |
３行目。(x^y)^xに等しく、それはyに等しい |
||
| 0 | 0 | 0 | 0 | 0 |
| 1 | 0 | 1 | 1 | 0 |
| 0 | 1 | 1 | 0 | 1 |
| 1 | 1 | 0 | 1 | 1 |
上の表を見ているとわかる通り、xor演算には以下の性質があります。
- xor演算は交換法則が成り立つ:
x^y = y^x - xor演算は結合法則が成り立つ（カッコをつけたり外してもいい）:
(x^y)^z = x^(y^z) - 同じxor演算を施すと元に戻る:
(x^y)^y = x
ちなみに、昔は上記のようなビットトリックを用いることにメリットがあったそうなのですが、現在のプロセッサではこのようにビット演算を用いて書くよりも、 普通に書いてコンパイラに最適化させるほうが性能が良い場合が多いです。
こういった、ビットを扱う処理をビットトリックと呼んだりします。 黒魔術とも言える様々な難解ビットトリックが世の中には存在します。 こちらなどを見ると色々と載っています。日本語の記事ではこちらが非常にわかりやすいです。
やってみよう（30分）
上記を写経してみましょう
クイズ
- ビット単位の和を計算する
&と、論理演算子の&&を間違えないようにしましょう。同様に、|と||を間違えないようにしましょう。 これらは別者ですが、時として同じ値になることがあるので、注意が必要です。以下の場合を考えましょう。unsigned char i = 1; unsigned char j = 0;とき、上記の4演算を試すとどうなるでしょうか。unsigned char i = 3; unsigned char j = 3;とき、上記の4演算を試すとどうなるでしょうか。
- ハミング距離のコードを書いてみましょう。ハミング距離とは、二つのビット列が与えられたときに、各桁について、値が同じなら0、違えば1、として、それをまとめたものです。
ここでは、ビット列の長さは8とします。
まずはビット演算ではなく、
intの配列として考えてみましょう。以下のような関数hamming_dist_array(int x1[], x2[])を考えてくださいint v1[8] = {0, 0, 1, 0, 0, 1, 1, 0}; int v2[8] = {0, 1, 0, 0, 1, 0, 1, 0}; int d = hamming_dist_array(v1, v2); // 4になる - 次に、8桁のビット列を
unsigned charで表しましょう。これに対するハミング距離計算の関数hamming_dist(unsigned char c1, unsigned char c2)を書いてみましょう。ヒント：まず各桁における値の差異を求め、それを最後に集計すればよいです。unsigned char v1 = 0b00100110; unsigned char v2 = 0b01001010; int d = hamming_dist(v1, v2); // 4になる - 最後に、
hamming_distを高速化するアイデアを考えてみてください。
答え
- 以下のようになります。
- 最初の例では、ビット演算と論理演算では同じになります
unsigned char i = 1; unsigned char j = 0; printf("%d\n", i & j); // 0 printf("%d\n", i && j); // 0 printf("%d\n", i | j); // 1 printf("%d\n", i || j); // 1 - 二つ目の場合は、ビット演算と論理演算で結果が違います。
unsigned char i = 3; unsigned char j = 3; printf("%d\n", i & j); // 3 printf("%d\n", i && j); // 1 printf("%d\n", i | j); // 3 printf("%d\n", i || j); // 1
- 最初の例では、ビット演算と論理演算では同じになります
- 例えば以下のようになります。
int hamming_dist_array(int x1[], int x2[]) { int d = 0; for (int i = 0; i < 8; ++i) { if (x1[i] != x2[i]) { ++d; } } return d; } - 例えば以下のようになります。
int popcount(unsigned char a) { // a中に出現する1の回数をカウントする // K&R, p62 int cnt; for(cnt = 0; a != 0; a >>= 1) { if (a & 1) { ++cnt; } } return cnt; } int hamming_dist(unsigned char c1, unsigned char c2) { unsigned char diff = c1 ^ c2; // まずはXORで差分を一気にもとめる return popcount(diff); // そのうち、1が出現している回数を数える } - 一つの例は、
popcountを高速化することです。popcountは入力としてunsigned charを取ります。unsigned charは256通りの整数しか表現できないので、あらかじめ全ての結果のパターンを計算しておき、 表として持っておくことができます。計算時にはその表を見に行くことで、高速にpopcountを実行できます。 まず次のような表を事前に用意しておきます。そして、次のように// この表は「ビットを数える・探すアルゴリズム」から引用 // http://www.nminoru.jp/~nminoru/programming/bitcount.html int count_table[256] = { 0, 1, 1, 2, 1, 2, 2, 3, 1, 2, 2, 3, 2, 3, 3, 4, 1, 2, 2, 3, 2, 3, 3, 4, 2, 3, 3, 4, 3, 4, 4, 5, 1, 2, 2, 3, 2, 3, 3, 4, 2, 3, 3, 4, 3, 4, 4, 5, 2, 3, 3, 4, 3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 1, 2, 2, 3, 2, 3, 3, 4, 2, 3, 3, 4, 3, 4, 4, 5, 2, 3, 3, 4, 3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 2, 3, 3, 4, 3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 1, 2, 2, 3, 2, 3, 3, 4, 2, 3, 3, 4, 3, 4, 4, 5, 2, 3, 3, 4, 3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 2, 3, 3, 4, 3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 2, 3, 3, 4, 3, 4, 4, 5, 3, 4, 4, 5, 4, 5, 5, 6, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 3, 4, 4, 5, 4, 5, 5, 6, 4, 5, 5, 6, 5, 6, 6, 7, 4, 5, 5, 6, 5, 6, 6, 7, 5, 6, 6, 7, 6, 7, 7, 8, };popcountを変更すればOKですさて、このpopcountという処理は重要な処理なので、これまでに多くの計算方式が考えられてきました。 例えば先ほどのビットトリック集のCounting bits setというところを見ると色々書いてありますので、眺めてみると面白いです。 そして、このpopcountの機能はハードウェア上で実装されている場合もあります。その場合、例えばGCCがそのハードウェアを直接扱う関数を 準備していたりします（int popcount_table(unsigned char a) { return count_table[a]; }int __builtin_popcount (unsigned int x)）。これを使うと、上記のように アルゴリズムで高速化するよりも、速いです。
ライブラリツアー¶
Cでは講義で扱わなかった様々な便利関数があります。それらは適切なファイルをincludeすると使えるようになります。以下に見てみましょう。
これらのファイルの詳細はたとえばこういうところを見ると書いてあります。
数学関数：math.h¶
math.hをincludeすると、様々な数学関数が使えます。
math.hを使う場合、コンパイル時に、gcc -lm main.cのように-lmオプションが必要かもしれません。
printf("PI: %f\n", M_PI); // 3.141593 円周率が定義されています
printf("cos(PI): %f\n", cos(M_PI)); // -1 コサイン
printf("exp(3): %f\n", exp(3.0)); // 20.085537 exp
printf("pow(3.0, 2.0: %f\n", pow(3.0, 2.0)); // 9.000000 累乗
printf("sqrt(9.0): %f\n", sqrt(9.0)); // 3.000000 ルート
printf("fabs(-2.5): %f\n", fabs(-2.5)); // 2.500000 絶対値
整数値の限界：limits.h¶
limits.hをインクルードすると、整数値の限界に関する情報を入手することが出来ます。
たとえばINT_MINという表記で、int型がとりうる最小値を知ることが出来ます。
printf("INT_MIN = %d\n", INT_MIN); // -2147483648
printf("INT_MAX = %d\n", INT_MAX); // 2147483647
小数に関する情報：float.h¶
float.hをインクルードすると、小数に関する情報を入手することが出来ます。
たとえばFLT_MINという表記で、float型がとりうる最小値を知ることが出来ます。
printf("FLT_MIN = %e\n", FLT_MIN); // 1.175494e-38
printf("FLT_MAX = %e\n", FLT_MAX); // 3.402823e+38
文字判定：ctype.h¶
ctype.hをインクルードすると、文字に関する様々な関数
が使えるようになります。
// アルファベットかどうかの判定
printf("isalpha('a'): %d\n", isalpha('a')); // 1024 非0の値になる
printf("isalpha('3'): %d\n", isalpha('3')); // 0
// 数字かどうかの判定
printf("isdigit('a'): %d\n", isdigit('a')); // 0
printf("isdigit('3'): %d\n", isdigit('3')); // 2048 非0の値になる
ローレベルなメモリコピー：string.h¶
string.hをincludeすると文字列に対する様々な関数が使えることは既に述べました。
それ以外にも、ローレベルで直接配列をコピーするものがあります。
char src[] = {'a', 'b', '\0', 'c', 'd'};
char dst1[10];
char dst2[10];
// 文字列のコピー。'\0'が出るところまでをコピー
strcpy(dst1, src); // 'a', 'b', '\0'
// 文字列に限らない、任意のメモリ分量のコピー
// 三番目の引数はコピーする`char`の個数
memcpy(dst2, src, 5); // 'a', 'b', '\0', 'c', 'd'
警告：assert.h¶
assert.hには、警告に関する関数が含まれます。
これは、次のように、「必ずある条件を満たしてほしい」ということを示すために使います。
int N = 10;
int a[N];
int n = getchar();
n -= '0';
// assertの中身が真であれば何もしない。偽であればエラーになる。
// この場合、キーボードから入力された文字がN未満の数字であればOK
assert(0 <= n && n < N);
a[n] = 123;
printf("a[%d]: %d\n", n, a[n]);
乱数、system：stdlib.h¶
stdlib.hには様々な機能があることを既に勉強してきましたが、
乱数も扱うことができます。
乱数とは、その名の通りランダムな数です。
ここではrand()関数によりランダムなintが返ります。
また、取りうる値の最大値としてRAND_MAXも準備されているため、
それで割ることで疑似的に0 - 1のレンジの乱数も作れます。
ここでsrand(seed)というのは、乱数の種（seed）を与えるものです。
疑似乱数生成器は同じseedに対しては同じ乱数列を与えます。
このseedを同じにすれば、乱数であっても結果が再現できます。
srand(123);
for (int n = 0; n < 4; ++n) {
printf("%d\n", rand());
}
// 128959393
// 1692901013
// 436085873
// 748533630
for (int n = 0; n < 4; ++n) {
double r = (double) rand() / RAND_MAX;
printf("%f\n", r);
}
// 0.361609
// 0.134639
// 0.375968
// 0.259322
また、システムに関する関数も準備されています。
例えばsystem(cmd)関数はcmdをターミナルから実行したものと同じことをする関数です。
これは危険なのであまり使わないほうがいいです。
system("ls"); // ターミナルでlsと打った時とおなじように、現在のディレクトリのファイルを表示
やってみよう（15分）
上記を写経してみましょう
宿題¶
それでは今週の宿題です。
- 締切は次の授業の前日の深夜23:59までです（今回の場合、2021/11/23, 23:59）
- 宿題リンクをslackで配付します。このリンクは公開しないでください。
week6_1¶
二次元点群のICP (Iterative Closest Point)を計算しましょう。 点群とは、物体の形状を点の集合であらわすというものです。三次元データに対して、その表面をサンプリングしたようなものです。 二つの同じ物体の点群が与えられたとき、その位置をピッタリ合わせることを点群位置合わせといいます。 例えば3Dスキャナで部屋を撮影しそれを点群で表現し、そのあと部屋にある椅子を撮影し点群で表現したとき、元の部屋点群のなかのどこに椅子があったか探すときに使ったりします。 そのための最も基本的なアルゴリズムがICPです。今回は二次元の場合で、かつ回転運動がないという非常に簡略化された場合を考えてみましょう。
以下の図を見て下さい。
ここでは、赤い飛行機の概形が、\(N\)個の点による点群\(\mathcal{P}_\mathrm{airplane} = \{\mathbf{p}_n\}_{n=1}^N\)として表されています。ここで一つの点は単純なx,y座標です：\(\mathbf{p}_n \in \mathbb{R}^2\)
ここで、飛行機の先端部分の点群である\(\mathcal{Q}_\mathrm{head} = \{\mathbf{q}_m\}_{m=1}^M\)を考えます。ここでは点の個数は\(M\)であり、\(\mathbf{q}_m \in \mathbb{R}^2\)です。 この緑の先端部分は、少し赤い飛行機からずれています。 このとき、この二つの点群がピッタリ合うように、先端部部分を移動させるというのが、今回のタスクです。
それでは宿題リポジトリを見てみましょう。
ここにはairplane.txtやhead1.txtといった点群データがあります。これらは、
- 一行目が点の個数
- 二行目以降が、点の座標。一行が「x y」というフォーマット
となっています。
ここで、補助スクリプトのplot.pyを実行してみましょう。このpythonファイルは、二つの点群ファイルを引数にとります。
python plot.py airplane.txt head1.txt
out.pngというファイルが生成されます。vscode上でそれをクリックすると、中身が見えます。
これは、airplane.txtの点群データとhead1.txtのデータを重ねて可視化してくれる便利スクリプトです。これを使って、自分の結果も確認できます。
一つ目の引数の点群が赤で、二つ目の引数の点群が緑で表示されます。
そして、距離計算補助スクリプトのpc_dist_1_to_2.pyも実行してみましょう。
python pc_dist_1_to_2.py head1.txt airplane.txt
これは、二つの点群がどの程度離れているかを計算してくれます。以下のようになります。
Distance from head1.txt to airplane.txt is 10.8
さて、ここでmain.cをコンパイルして実行してみましょう。
./a.out head1.txt airplane.txt
head1.txtを位置合わせすることが目標です。
こうすると、緑の先端を移動させたout.txtが生成されます。これをチェックしてみましょう。
python plot.py airplane.txt out.txt
コードの中では、点を表すPointや、点群を表すPointCloudといった構造体が定義されています。これらがどういう意味なのか、コードを読んで確認してみてください。
さて、今回の問題を考えるうえで、まず点の対応というものを考えてみましょう。先端部分の各点から、最も近い（ユークリッド距離、あるいはその二乗が小さい）点を飛行機の中から探してみましょう。 それは次のように書けます。
これが、「\(\mathbf{q}\)から見て、飛行機の中で対応する点の添え字」ということになります。
ここで、先端側からみた飛行機までの距離（エラー。小さくしたいもの）は、対応点までの距離の平均になります。
これを最小化してください。
この簡単な方法は、残差ベクトルの平均である\(\mathbf{t} \in \mathbb{R}^2\)を求めるというものです。
この\(\mathbf{t}\)を全ての\(\mathbf{q}_m\)に足すと、位置が更新されます。
よって、最終的にこの問題は、イケている併進ベクトル\(\mathbf{t}\)をどう求めるかということに帰着されます。
そして、上記のプロセスを何度も繰り返すと、だんだん正解位置に近づいて行ってくれます。
上記は今回のように併進だけの場合で、かつ初期位置がとても近くないとうまくいかないナイーブな方法です。 もしより高度な方法を知りたければ、関連書籍などを見てみるといいでしょう。
いつ戻り、main関数はいじらないでください。関数は自由に追加していいです。cの標準のライブラリは使っていいです（math.hとか）。
head1.txt¶
まずは簡単なhead1.txtをやってみましょう。これは実際に飛行機の先頭部分を切り出して少し移動させたものなので、かなりうまく行きます。
自動採点は以下です。ここでthreの部分は合格の閾値をあらわしています。今回は1以下ならOKです。
入力：
./a.out head1.txt airplane.txt
python pc_dist_1_to_2.py out.txt airplane.txt --thre 1
出力（最終行）：
Less than 1.0. Ok!
head2.txt¶
次は少し難しいhead2.txtです。これは同じ飛行機図形から別のサンプリングを行っているので、すこし難しいです。厳密対応にはなりません。
入力：
./a.out head2.txt airplane.txt
python pc_dist_1_to_2.py out.txt airplane.txt --thre 3
出力（最終行）：
Less than 3.0. Ok!
head3.txt¶
最後に、移動幅を大きくしたhead3.txtです。これは緩い閾値になっています。
入力：
./a.out head3.txt airplane.txt
python pc_dist_1_to_2.py out.txt airplane.txt --thre 5
出力（最終行）：
Less than 5.0. Ok!
松井の手元ではもっと難しい例（初期位置がかなり離れている、など。） 難しい例でもうまくいくアルゴリズムになっていると、加点するかもしれません。
答えの例¶
答え
week6_1
#include <stdio.h>
#include <stdlib.h>
#define N_MAX 2000
// 関数は自由に追加していいです。
// Point: 点を表す構造体
typedef struct point {
double x;
double y;
} Point;
// PointCloud: 点群を表す構造体。
// N: 点群の個数
// pts: 点群を表す配列。N_MAX個だけ確保されるが、実際に使っている点は上位N個だけ。
typedef struct point_cloud {
int N;
Point pts[N_MAX];
} PointCloud;
double L2sqr(Point p1, Point p2) {
// 二点のユークリッド距離の二乗を返す
return (p1.x - p2.x) * (p1.x - p2.x) + (p1.y - p2.y) * (p1.y - p2.y);
}
int NN(Point p, PointCloud pc) {
// 点pとポイントクラウドpcが与えられたとき、pc中でpにもっとも近いもののインデックスを返す
// min/maxを探す定石は、「ベスト候補」を用意しておき、ループ中でそれを更新するもの。
int min_n = 0; // ベスト候補のID
double min_dist = L2sqr(p, pc.pts[0]); // ベスト候補の値。初期値は一個目の要素にしておけば安全
// 全点を舐める
for (int n = 0; n < pc.N; ++n) {
double dist = L2sqr(p, pc.pts[n]);
// ベスト候補よりも良い場合、更新
if (dist < min_dist) {
min_dist = dist;
min_n = n;
}
}
return min_n;
}
Point ICP(PointCloud pc1, PointCloud pc2) {
// 最終残差
Point res_final = {.x = 0, .y = 0};
// 適当に10回行う
for (int t = 0; t < 10; ++t) {
// あるイテレーション中で考える残差。
Point res = {.x = 0, .y = 0};
// pc1の全点に対し
for (int n = 0; n < pc1.N; ++n) {
// pc2中で一番近いものを求める
int nn = NN(pc1.pts[n], pc2);
// 残差を集めていく
res.x += pc2.pts[nn].x - pc1.pts[n].x;
res.y += pc2.pts[nn].y - pc1.pts[n].y;
}
// 残差は平均値にしておく
res.x /= pc1.N;
res.y /= pc1.N;
// pc1を更新
for (int n = 0; n < pc1.N; ++n) {
pc1.pts[n].x += res.x;
pc1.pts[n].y += res.y;
}
// 最終残差の更新
res_final.x += res.x;
res_final.y += res.y;
}
return res_final;
}
// 点群を表示する関数
void ShowPointCloud(PointCloud pc) {
printf("N: %d\n", pc.N);
for (int n = 0; n < pc.N; ++n) {
printf("%3.2f %3.2f\n", pc.pts[n].x, pc.pts[n].y);
}
}
// 点群を読みこむ関数
PointCloud ReadPointCloud(char filename[]) {
FILE *fp;
fp = fopen(filename, "r");
if (fp == NULL) {
printf("Cannot open %s\n", filename);
exit(1);
}
PointCloud pc;
fscanf(fp, "%d", &(pc.N));
for (int n = 0; n < pc.N; ++n) {
fscanf(fp, "%lf %lf", &(pc.pts[n].x), &(pc.pts[n].y));
}
fclose(fp);
return pc;
}
// 点群を書き込む関数
void WritePointCloud(char filename[], PointCloud pc) {
FILE *fp;
fp = fopen(filename, "w");
if (fp == NULL) {
printf("Cannot open %s\n", filename);
exit(1);
}
fprintf(fp, "%d\n", pc.N);
for (int n = 0; n < pc.N; ++n) {
fprintf(fp, "%.1f %.1f\n", pc.pts[n].x, pc.pts[n].y);
}
fclose(fp);
}
int main(int argc, char *argv[]) {
// ==================== main関数は触らない ==================
if (argc != 3) {
printf("Error. Set two argment. For example: $ ./a.out head1.txt airplane.txt\n");
return 1;
}
PointCloud pc1 = ReadPointCloud(argv[1]);
PointCloud pc2 = ReadPointCloud(argv[2]);
// 残差ベクトルを作る
Point t = ICP(pc1, pc2);
// 残差ベクトルをpc1の全ての点に足す
for (int n = 0; n < pc1.N; ++n) {
pc1.pts[n].x += t.x;
pc1.pts[n].y += t.y;
}
// それを保存
WritePointCloud("out.txt", pc1);
return 0;
}