ファッションドリーマーのビンゴゲームでquadrupleが出る確率を計算する

ファッションドリーマーは、イヴのミューズのコーデをコクーンでルカットするゲームです。

このゲームにはおまけのミニゲームとしてビンゴゲームがあり、ゲーム内で色々と行動すると入手できるビンゴチケットを使って遊ぶことができます。以下のような仕様です。
- 1から75までの数字のうち、ランダムな25個を5x5のマスに配置したビンゴシートがある。(中央が最初から空いていたりはしない)
- ビンゴチケットを1枚使うと数字が1つ抽選されて、その数字がシートにある場合は、印がつけられる。
- 縦・横・斜めのラインに5つ印が揃ったら当たりで、ゲーム内ポイントが入手できる。当たったらビンゴシートや抽選状況はリセットされ、次はまた最初から始まる。

そして、このビンゴゲームの大きな特徴として、同時に複数のラインが揃った場合にボーナスポイントがもらえるというルールがあります。
最大で縦・横・斜め2つの4列が同時に揃う可能性があり、quadrupleと呼ばれています。
ただし、どこかのラインが揃った時点で当たりとなって終了してしまうので、揃ってしまわないように交差する形でリーチをかけ、交点を最後に開ける必要があります。

このビンゴゲームはビンゴチケットさえあれば繰り返し数字を抽選できるため、いずれはどこかのラインが揃って当たりになりますが、そのときquadrupleになる確率を計算してみます。


## 解くべき問題の定義

考えやすくするために、まずはできるだけ問題を単純化します。

前提として、数字の抽選は無作為で偏りがないものとします。(実際には偏っている可能性は十分にありますが、確かめるのは大変なので諦めます)

数字は75個あり、うち25個のみがビンゴシートに書かれているため、50個はハズレです。しかし、当たりになるまでひたすら回すと考えると、ハズレ数字が途中のどのタイミングでいくつ抽選されようと関係がないです。
また、どの25個の数字がどのような並びでビンゴシートに書かれていたとしても、数字の抽選が無作為であれば結果には関係がないです。
そこで、
- ビンゴシートの25個のマス目に1から25の数字が順番に並んでいる
- 抽選される数字は1から25までのいずれか
としても結果には影響がないでしょう。このとき、数字が抽選される順序は、 25! 通りとなり、それらは全て同じ確率で発生すると考えられます。

## 解き方
数字が抽選される順序は、 25! 通りであり、それらを全て試せば答えは得られるはずですが、25! = 15511210043330985984000000 と多すぎるので難しいです。

現実的な時間で解くために、動的計画法と呼ばれる方法を使います。
ビンゴシートの状態は、25箇所それぞれについて抽選済みかどうかを1か0で区別する2進値として表現することができます。
例えば、`0000000000000000000000000` を最初の状態、`0000000000001000000000000` は真ん中の1マスだけが抽選済みの状態、`1111111111111111111111111` は全てのマスが抽選済みの状態(実際にはその前にどこかのラインが揃って終了するため、この状態にはなりません)、のようになります。
このとき、取りうる状態の数は 2^25 = 33554432 通りです。
最初の状態 `0000000000000000000000000` から、1つ目の数字を抽選することを考えてみます。このとき、状態は `0000000000000000000000001` や `1000000000000000000000000` など、1箇所だけ0から1に変わることとなり、結果の状態は25種類のいずれかとなります。

その25通りの状態それぞれについて、2つ目の数字を抽選すると、それぞれ抽選結果には24通りの可能性があるので、結果は 25 * 24 通りです。これを繰り返すことで、 25 * 24 * 23 * ... * 1 となり、 25! 通りとなります。

ところが、実は2回目の抽選が終わった時点を考えてみると、抽選された順番が異なるだけで結果の状態は同じであるパターンがあります。例えば次の2つです。
`0000000000000000000000000` ---1回目の抽選--> `0000000000000000000000001` ---2回目の抽選--> `0000000000000000000000011`
`0000000000000000000000000` ---1回目の抽選--> `0000000000000000000000010` ---2回目の抽選--> `0000000000000000000000011`
この後は3回目の抽選を行うわけですが、3回目の抽選により起きる状態変化は、2回目の抽選が終わった時点での状態にのみ依存し、それまでに数字が抽選された順番には関係がありません。
なので、2回の抽選が終わった時点で状態が 0000000000000000000000011 であるものは2パターンである、というようにパターンの数だけを覚えておくことで、計算をまとめられます。
状態の種類は最大でも 2^25 = 33554432 なので、抽選の回数を掛けて 33554432 * 25 = 838860800 通りの可能性を計算すればよくなり、これは、状態でまとめずに全ての可能性を計算した場合の 25! = 15511210043330985984000000 と比べてずっと小さく、現実的に計算可能です。
なお、すでにどこかのラインが揃っている場合はビンゴゲームが終了するので、その後の抽選では状態が変化しないようにします。

## 結果
というわけで計算するプログラムをC++で書きました。大きな数を扱うため、128ビット整数を使っています。

```bingo_dreamer.cpp
#include <bits/stdc++.h>
using namespace std;

__int128 factorical(int x)
{
    __int128 f = 1;
    for (int i = 1; i <= x; ++i) {
        f *= i;
    }
    return f;
}

int pow(int x, int y)
{
    __int128 p = 1;
    for (int i = 0; i < y; ++i) {
        p *= x;
    }
    return p;
}

const int SHEET_SIZE = 5;
// 次のように数字が並んでいるとして、 縦と横がそれぞれ5つ、斜めが2つの12ラインの当たりパターンを定義
// 1  6 11 16 21
// 2  7 12 17 22
// 3  8 13 18 23
// 4  9 14 19 24
// 5 10 15 20 25
//
const int BINGO_LINES[] = { (1 << 0) + (1 << 1) + (1 << 2) + (1 << 3) + (1 << 4),
    (1 << 5) + (1 << 6) + (1 << 7) + (1 << 8) + (1 << 9),
    (1 << 10) + (1 << 11) + (1 << 12) + (1 << 13) + (1 << 14),
    (1 << 15) + (1 << 16) + (1 << 17) + (1 << 18) + (1 << 19),
    (1 << 20) + (1 << 21) + (1 << 22) + (1 << 23) + (1 << 24),
    (1 << 0) + (1 << 5) + (1 << 10) + (1 << 15) + (1 << 20),
    (1 << 1) + (1 << 6) + (1 << 11) + (1 << 16) + (1 << 21),
    (1 << 2) + (1 << 7) + (1 << 12) + (1 << 17) + (1 << 22),
    (1 << 3) + (1 << 8) + (1 << 13) + (1 << 18) + (1 << 23),
    (1 << 4) + (1 << 9) + (1 << 14) + (1 << 19) + (1 << 24),
    (1 << 0) + (1 << 6) + (1 << 12) + (1 << 18) + (1 << 24),
    (1 << 4) + (1 << 8) + (1 << 12) + (1 << 16) + (1 << 20) };

// シートの状態の種類数 2^25
const int STATE_LENGTH = pow(2, SHEET_SIZE* SHEET_SIZE);

// 抽選順序の総数 25!
const __int128 ALL_COUNT = factorical(SHEET_SIZE * SHEET_SIZE);

int main()
{
    // 状態それぞれについて、何列揃っているかを事前に計算しておく
    vector<int> hit_lines(STATE_LENGTH, 0);
    int max_hit_lines = 0;
    for (int i = 0; i < STATE_LENGTH; ++i) {
        int h = 0;
        for (int j = 0; j < sizeof(BINGO_LINES) / sizeof(int); ++j) {
            if ((i & BINGO_LINES[j]) == BINGO_LINES[j])
                ++h;
        }
        hit_lines.at(i) = h;
        if (h > max_hit_lines)
            max_hit_lines = h;
    }

    vector<__int128> dp(STATE_LENGTH, 0);
    // 最初はどこも開いていない状態が1通りだけある
    dp.at(0) = 1;

    // 状態遷移を計算する
    for (int i = 0; i < 25; ++i) {
        for (int j = STATE_LENGTH - 1; j >= 0; --j) {
            if (dp.at(j)) {
                if (hit_lines.at(j))
                    dp.at(j) = dp.at(j) * (25 - i);
                else {
                    for (int k = 0; k < 25; ++k) {
                        if ((1 << k) & j)
                            continue;
                        int nj = j + (1 << k);
                        dp.at(nj) += dp.at(j);
                    }
                    dp.at(j) = 0;
                }
            }
        }
    }

    // 同時に揃っている列数ごとに集計
    vector<__int128> hit_count(max_hit_lines + 1, 0);
    for (int i = 0; i < STATE_LENGTH; ++i) {
        if (dp.at(i))
            hit_count.at(hit_lines.at(i)) += dp.at(i);
    }

    for (int i = 1; i <= max_hit_lines; ++i) {
        if (hit_count.at(i) != 0)
            cout << i << " : " << ((double)((__float128)(hit_count.at(i)) / (__float128)ALL_COUNT)) << endl;
    }
}

```

実行結果です。同時に揃う列の数と、その確率が出力されています。
```
$ g++ -O2 -o bingo_dreamer bingo_dreamer.cpp 
$ time ./bingo_dreamer 
1 : 0.920612
2 : 0.0765579
3 : 0.00279262
4 : 3.70375e-05

real    0m2.964s
user    0m2.887s
sys     0m0.070
```

計算によると、2列が同時に揃う double は約7.7%で発生するようなので、何度もビンゴゲームをやっていればいずれ遭遇できそうです。
quadrupleは約0.0037%で、およそ3万回に1回の低確率です。自分では見れそうにないですが、プレイヤー全体では目撃している人もいそうな数字ですね。

なお、そもそもビンゴゲームよりも他のことをしたほうが圧倒的に楽しい上にポイントも多く得られる模様です……