0からはじめる計算幾何学第04回過去のプログラムの有効活用

昔書いたプログラムをもう一度使えないかと考えることは、多くのプログラマにとって（そしてエントリのネタを探す私にとっても）重要なことです。今日は、過去の遺産を有効活用して楽に問題を解いた事例を紹介しましょう。UVa Online Judge の 10012番、How Big Is It? です。

まずは問題の概要を説明しましょう。

いくつかの円が与えられます。これらの円を全て含むことが出来る最小の長方形の横幅を求めるプログラムを作成しなさい。ただし、どの円も長方形の底辺に接するように含めなければなりません。

http://icpcres.ecs.baylor.edu/onlinejudge/index.php?option=com_onlinejudge&page=show_problem&problem=953

以前のエントリで紹介した問題とほとんど同じ問題です。大きく違うところは以下の2つ。

「長方形の横幅の最小値が閾値以下か？」が「長方形の横幅の最小値はいくつか？」になった
「二つの円の間にもう一つの小さな円が入り込むことはない」という仮定がなくなった

それでは、以前のコードを書き換えてこれらの変更に対応させましょう。

まずは前者ですが、これは特筆すべきこともありません。ただ書き換えるのみです。

/*
円の半径の並び、円の個数を引数として受け取る。
円を詰められる最小の長方形の横幅を返す。
( next_permutationは順列を生成する。algorithmヘッダの中で定義されている)
*/
double execute(int *radii, int rNum){
double result;
sort(radii, radii+rNum);
result = calculate(radii, rNum);
while (next_permutation(radii, radii+rNum))
result = min(result, calculate(radii, rNum))
return result;
}

これに比べて、後者の変更は厄介です。アルゴリズムから考え直さなければいけません。

与えられた半径とその並べる順番に対してそれを詰められる長方形の最小の幅を返す、次のようなアルゴリズムを考えてみました。

1つ目の円は、長方形の左端に置く。
2つ目以降の円は、今までに長方形に詰めた円の中で中心が一番右側のものに接するように置く。ただし、それが不適切な置き方になってしまったときは、右から2番目の円に接するように置く。それでも駄目なら右から3番目、4番目・・・と試していく。
円の「不適切な置き方」というのは、その円が別の円と重なってしまう、あるいは長方形の左端のさらに左側にはみ出してしまうような置き方である。
全ての円を詰め終わったら、全ての円の右端の位置を調べ、その中で最も右側にあるものの左右方向の位置を返す。

さて、ここで不安になるのが計算量です。上記のアルゴリズムはO(n³)のオーダになっています。さらにこれを可能な円の並べ方の数だけ繰り返しますので、結局O(n³ × n!)の時間が必要になります。
一見絶望的なオーダですが、問題文によると円の数は最大で8個と書かれていますので、これでも十分間に合います。
それでは、実装してみましょう。

円の個数の最大値
const int MAX_R = 8;
/*
円の半径の配列と円の個数を引数として受け取る。
配列の通りに円を詰めたときの横幅を返す。
*/
double calculate(double *radii, int rNum){
// 円の中心の横方向の座標を格納する配列。
// 長方形の左端を0.0とする。
double position[MAX_R];
position[0]=radii[0];
for (int i=1; i<rNum; i++){
int adjacence = i;
bool invalid_place;
do {
adjacence--;
// 1番目の円の隣に(i)番目の円を置く場合
if (adjacence==0){
position[i] = max(position[i], radii[i]);
break;
}
// adjacence 番目の円の隣に i 番目の円を置く場合
position[i] = position[adjacence] + width(radii[adjacence], radius[i]);
// 不適切な置き方かどうかチェック
invalid_place = false;
// 他の円との重なりをチェック
for (int j=0; j<adjacence && !invalid_place; j++)
invalid_place = lap(radii, position, j, i);
// 長方形からのはみ出しをチェック
invalid_place = invalid_place ||  position[i]<=radii[i];
} while (invalid_place);
}
// 円の右端を計算
double result[MAX_R];
for (int i=0;i<rNum;i++)
result[i]=position[i]+radii[i];
// 最大のものを返す
return *max_element(result, result+rNum);
}
半径 r1 の円と半径 r2 の円を並べたとき、それらの中心の
横方向の距離を返す。
double width(double r1,double r2){
return sqrt((r1+r2)*(r1+r2)-(r1-r2)*(r1-r2));
}
r1_index番目の円とr2_index番目の円が重なっているかどうかを返す。
2つの円の中心同士の距離と、半径の和を比較することで求める。
bool lap(double *radii, double *position, int r1_index, int r2_index){
double r1_x = position[r1_index];
double r2_x = position[r2_index];
double r1_y = radii[r1_index];
double r2_y = radii[r2_index];
return (Molar_Morrison(r2_x-r1_x, r2_y-r1_y) < radii[r1_index]+radii[r2_index]);
}

分かりづらい点といえば、「1番目の円の隣に i 番目の円を置く場合」の部分でしょうか。
このケースでは、「置いた円が長方形をはみ出しているか」のチェックを行い、はみ出していれば長方形の左端と接するように置きなおすことになります。ここで、はみ出しているかどうかは position[i] < radii[i] を使って判定できますので、結局 position[i] = max(position[i], radii[i]); と書けば、はみ出しているかのチェックとその場合の値の書き換えが同時に出来ることになります。

接している2つの円の横幅を求める処理は、入力値がdouble型になってテーブル引きが出来なくなった関係で、毎回計算するように書き換えました。これをwidth関数としてまとめています。
また、2つの円の重なりを判定するlap関数にMolar_Morrison関数を使ってみました。これは(x, y)ベクトルの長さを計算する関数ですが、詳細はこちらをご覧ください。

このようにして、無事プログラムは完成しました。
プログラムを再利用することのメリットは、一度打ったコードをもう一度打ち直さないという作業の省略のみにとどまりません。アルゴリズムだけを読み取って再利用すれば、一度行った試行錯誤をもう一度繰り返さないという、思考の省略もできるのです。
もちろん、プログラムからアルゴリズムを読み取るためには、そのプログラムが綺麗に描かれていることが重要になってきます。「この関数をコピーして使うことはもうないから」などと言わず、常に読みやすいコードを書くことを心がけたいものです。

担当:田山