💡

再帰やループを使わずに組合せを求める

に公開
組合せ論
combinatorics
gosper
tech

再帰やループを使わずに組合せを求める

はじめに

世の中の言語や環境の中には，順列や組合せの生成といった高度なコレクション操作を標準で提供してくれないものが少なからずあります．かかる折には自前でジェネレータを定義するものですが，再帰だのループだのと書いているうちに目が回って疲れてしまいます．しかし，こと組合せの生成に関しては心配無用です．繰り返しを用いずにビット演算だけでこれを実現する方法があります．

組合せと二進数

組合せというのは，集団の中から何個か要素を取り出す時，どんなバリエーションが有り得るのかを考える分野です．例えば，7 色のボールが入った箱から 3 個取り出す場合，赤青緑だったり青黄橙といったふうに複数の組合せを考えることができます．この例だと組合せは全部で {}_{7}C_{3} = 35 通りあります．

中高の数学で習う単純な問題ですが，実際に計算しようとすると「何個」だの「取り出す」だのと余計な概念が思考の邪魔をしてきます．そこで，発想を転じてシンプルに考えることを試みるのであります．各要素は「箱の中にある」と「選ばれている」の 2 通りの状態しか持ちません．故に，全要素数を n とすれば，要素それぞれに 1 bit のフラグを与えることで，組合せを n 桁の二進数で表すことができます．7 色のボールの例で言えば，3 個取り出す場合の状態は 0001110 とか 0101010 といった具合に表現されます．

組合せを列挙するビット演算

このようにして組合せと二進数を同一視すると，n 個の中から k 個選ぶ組合せを列挙するとは，すなわち「1 を k 個含む n 桁の二進数」を数え上げる話に他ならぬことが判ります．そのような数 (以降 x とする) を小さい順に求めていくのが今回の目標です．難しそうですが，これは補数と排他的論理和の性質を利用することで案外あっさりと片付いてしまうのであります．

では，早速 JavaScript で実装してみます．x を与えると次の x を返す関数 nextCombination を定義します．JavaScript ではビット演算子のオペランドは 32 bit 整数として扱われますので，31 個以上 [1] を扱う場合は BigInt を用いましょう．

nextCombination.js
const nextCombination = x => {
  const a = x & -x;
  const b = a + x;

  return (((b ^ x) / a) >> 2) + b;
};

まず，x とその符号を反転した -x の論理積を求めます．これは，最も右にある 1 だけを残して残りの桁を全部 0 にする操作です．例えば，視野を 8 bit にして x を 0011 1000 と置くとします．符号を反転させた -x は，2 の補数表現により 1100 1000 となります．この二つの論理積を求めると，最も右側の 1 だけが都合よく残って 0000 1000 が得られます．これを a とします．

次に，aと x の和を求めて b とします．先程と同様に 8 bit の例で考えてみます．aと x，つまり 0000 1000 と 0011 1000 を足し合わせると 0100 0000 になります．この加算によって，右端の 1 で繰り上がりが生じていることが判ります．この b の値が「次の x」の種になります．

ただし，繰り上がりの影響により，このままでは b は「1 を k 個含む」という条件を満たしません．ここで，最後の (((b ^ x) / a) >> 2) + b という式が登場し，足りない 1 を埋め合わせてくれます．一つずつ分解して見ていきましょう．

  1. b ^ x

    これは，x と比較してどのビットが変化したのかを表すマスクです．二つの値の排他的論理和 (XOR) をとることで，変化した箇所が 1，変化していない箇所が 0 であるような値が得られます．先程の例では 0100 0000 と 0011 1000 の XOR で 0111 1000 となります．

  2. / a

    x から b への変化を最下位ビットから詰め直す操作です．a は最も右側の 1 だけを取り出した値，すなわち 1 を 1 つだけ含む数ですから，これは 2 の冪です．2 の冪で割るということは，右から数えた 1 の位置だけ，言い換えれば繰り上がりが起こった桁の分だけ右へシフトすることに相当します．例で求めた b ^ x (0111 1000) をa (0000 1000) で割ると，0000 1111 が得られます．a に含まれる 1 は右 (0 番目) から数えて 3 番目にあるので，右へ 3 回シフトした結果と考えることが出来ます．

  3. >> 2

    右端の 2 bit を捨てます．(b ^ x) / a の値には，繰り上がりの副作用として生じた余分な 2 bit が付いているからです．この演算の結果に含まれる 1 の個数は，x から b に変化する際に消えてしまった 1 の個数を表します．先程得られた 0000 1111 を 2 回右シフトしてみると 0000 0011 が得られます．この値には 1 が 2 個含まれることから，b には 1 が 2 個だけ足りないことが判ります．

  4. + b

    不足している 1 を b に補う操作です．右シフトによって得られた不足分をそのまま b に足してやると，これこそが「次の x」になります．例の値で計算すると，0000 0011 と 0100 0000 を足した 0100 0011 が答えです．

実際に計算する

nextCombination を用いて組合せを列挙します．7 個の中から 3 個選ぶ組合せを考えます．

const n = 7;
const k = 3;

const first = (1 << k) - 1;
const second = nextCombination(first);

// 文字列表現のための関数
const show = x => x.toString(2).padStart(n, "0");

console.log(show(first));         // "0000111"
console.log(show(second));  // "0001011"

1 が 3 個含まれるような 7 桁の二進整数が，きちんと得られています．first が (1 << k) - 1 で初期化されているのは気になりますね．まず 1 を左へ k 回シフトすると 2 の冪が作られます．そこから 1 を引けば，1 が k 個連続する最小の数が得られます．十進数において 10 の冪から 1 を引くと 9 のゾロ目が得られますが，これと同じことです．

任意回数で列挙するには，当然ながら再帰かループを用いる必要があります．

// 組合せ n C k を列挙するジェネレータ
function* combinations(n, k) {
  let x = (1 << k) - 1;

  const until = 1 << n;

  while (x < until) {
    yield x;

    x = nextCombination(x);
  }
}

// 7 C 3 の組合せ (35 通り)
const all = [...combinations(n, k)];

console.log(all.length);  // 35

// 全部表示
all.forEach(e => console.log(show(e)));

ここで得られた値について，各桁の 1 を対応するオブジェクトなり値なりに置換してやれば，様々な組合せ問題に応用できましょう．

Gosper’s Hack

この手法は Gosper’s Hack [2] と呼ばれるもので，米国の数学者 Ralph William Gosper Jr. さんに由来するとされています．ただ，これが果たして当人の考案であるのかどうかは不明であり，この手法を解説した論文も今のところ見当たりません．そもそもインターネットのない時代に生まれたアルゴリズムですので，一次資料は米国の大学かどこかに眠っているのでしょう．この手法は今やインターネット上における伝承的な知識と化しており，StackOverflow をはじめとするオンラインコミュニティにて様々な亜種が引用され続けています．

補足

  • 本記事の実装では k \ge 1 を想定しており，nextCombination は k = 0 の場合に正しい値を生成しません．数学的な厳密さを要求するのであれば {}_{n}C_{0} = 1 になるよう場合分けが必要です．n \le 0, k < n のような条件下でも同様．

  • JavaScript の除算 / は浮動小数点演算です．今回は b ^ x が a で割り切れるため結果的に整数になっています．

  • 算術右シフト演算 >> は符号を伝播します．(b ^ x) / a は非負であるため，明示的に論理右シフト >>> を用いても結果は同じになります．

脚注

  1. 先頭ビットは符号ビットなので，二進 31 桁の値を Number で扱うと左シフトの結果が負数になってしまいます． ↩︎

  2. Gosper’s Algorithm とは異なります． ↩︎

Discussion