麻雀の配牌は運が悪いだけ？三麻・四麻のシャンテン数と連続確率を検証

// 麻雀の初期配牌：シャンテン数別の組合せ数を厳密に全数計算
// 対象：三麻108枚／四麻136枚 × 子13枚／親14枚（第一打前）
// 方法：合法な完成形を列挙し、距離変換と牌種別の組合せ重みで集計
// 作成：AIによる作成支援 

#include <algorithm>
#include <array>
#include <cassert>
#include <chrono>
#include <cstdint>
#include <fstream>
#include <functional>
#include <iostream>
#include <limits>
#include <map>
#include <random>
#include <string>
#include <unordered_map>
#include <vector>

using U = __uint128_t;  // 巨大な組合せ数を扱う128ビット整数
using D = std::array<unsigned char, 10>;
constexpr int INF = 99;

// 128ビット整数を10進数文字列に変換
std::string dec(U x) {
    if (!x) return "0";

    std::string s;
    while (x) {
        s += char('0' + x % 10);
        x /= 10;
    }
    std::reverse(s.begin(), s.end());
    return s;
}

// 組合せ C(n, k)
U choose(int n, int k) {
    if (k < 0 || k > n) return 0;

    U r = 1;
    for (int j = 1; j <= k; j++) {
        r = r * (n - j + 1) / j;
    }
    return r;
}

// 牌グループの状態をまとめるためのシグネチャ
struct Sig {
    D d{};
    int kinds = 0, pairs = 0, term = 0, tp = 0;

    bool operator==(Sig const& b) const {
        return d == b.d && kinds == b.kinds && pairs == b.pairs &&
               term == b.term && tp == b.tp;
    }
};

struct Hash {
    size_t operator()(Sig const& s) const {
        size_t r = 0;
        for (auto v : s.d) r = r * 131 + v;
        for (int v : {s.kinds, s.pairs, s.term, s.tp}) r = r * 131 + v;
        return r;
    }
};

struct Entry {
    Sig sig;
    uint64_t w;
};

// ─────────────────────────────────────────────
// 1. 牌種別の合法な完成形と距離テーブルを作る
//    数牌9種／字牌7種／三麻萬子2種を別々に前計算
// ─────────────────────────────────────────────
struct Table {
    int n, size;
    bool runs, allterm;
    std::vector<int> pow5;
    std::array<std::vector<unsigned char>, 10> dist;
    std::array<std::vector<Entry>, 15> groups;

    Table(int n_, bool runs_, bool term_)
        : n(n_), size(1), runs(runs_), allterm(term_) {

        // 各牌の所持枚数0〜4を5進数で表す
        for (int i = 0; i < n; i++) {
            pow5.push_back(size);
            size *= 5;
        }
        for (auto& a : dist) a.assign(size, INF);

        // 刻子と、許可されている牌種では順子を列挙
        std::vector<std::vector<int>> melds;
        for (int i = 0; i < n; i++) {
            std::vector<int> a(n);
            a[i] = 3;
            melds.push_back(a);
        }
        if (runs) {
            for (int i = 0; i + 2 < n; i++) {
                std::vector<int> a(n);
                a[i] = a[i + 1] = a[i + 2] = 1;
                melds.push_back(a);
            }
        }

        // 4面子＋雀頭までの合法な目標形を列挙する
        std::vector<int> a(n);
        std::function<void(int, int)> rec = [&](int start, int m) {
            int code = 0;
            for (int i = 0; i < n; i++) code += a[i] * pow5[i];

            dist[m][code] = 0;
            for (int i = 0; i < n; i++) {
                if (a[i] <= 2) dist[5 + m][code + 2 * pow5[i]] = 0;
            }
            if (m == 4) return;

            for (int j = start; j < (int)melds.size(); j++) {
                bool ok = true;
                for (int i = 0; i < n; i++) {
                    if (a[i] + melds[j][i] > 4) ok = false;
                }
                if (!ok) continue;

                for (int i = 0; i < n; i++) a[i] += melds[j][i];
                rec(j, m + 1);
                for (int i = 0; i < n; i++) a[i] -= melds[j][i];
            }
        };
        rec(0, 0);

        // 足りない牌の枚数を、1牌種ずつ距離変換で求める
        for (int z = 0; z < 10; z++) {
            size_t targets = std::count(dist[z].begin(), dist[z].end(), 0);
            std::cerr << "TARGETS n=" << n << " h=" << z / 5
                      << " m=" << z % 5 << " count=" << targets << "\n";

            // min_y f(y) + Σ max(y_i - c_i, 0) を座標ごとに計算
            for (int step : pow5) {
                for (int base = 0; base < size; base += 5 * step) {
                    for (int off = 0; off < step; off++) {
                        int p = base + off;
                        for (int c = 1; c < 5; c++) {
                            dist[z][p + c * step] = std::min(
                                dist[z][p + c * step],
                                dist[z][p + (c - 1) * step]
                            );
                        }
                        for (int c = 3; c >= 0; c--) {
                            dist[z][p + c * step] = std::min<int>(
                                dist[z][p + c * step],
                                dist[z][p + (c + 1) * step] + 1
                            );
                        }
                    }
                }
            }
        }

        // 物理牌4枚からの選び方 C(4, c) で重み付けする
        std::array<std::unordered_map<Sig, uint64_t, Hash>, 15> agg;
        constexpr uint64_t w[5] = {1, 4, 6, 4, 1};
        size_t vectors = 0;

        for (int code = 0; code < size; code++) {
            Sig s;
            int v = code, total = 0;
            uint64_t weight = 1;

            for (int i = 0; i < n; i++) {
                int c = v % 5;
                v /= 5;
                total += c;
                weight *= w[c];
                s.kinds += c > 0;
                s.pairs += c >= 2;
                if (allterm || i == 0 || i == n - 1) {
                    s.term += c > 0;
                    s.tp += c >= 2;
                }
            }
            if (total > 14) continue;

            for (int z = 0; z < 10; z++) s.d[z] = dist[z][code];
            agg[total][s] += weight;
            vectors++;
        }

        size_t signatures = 0;
        for (int k = 0; k <= 14; k++) {
            U sum = 0;
            for (auto const& [s, wg] : agg[k]) {
                groups[k].push_back({s, wg});
                sum += wg;
            }
            assert(sum == choose(4 * n, k));
            signatures += groups[k].size();
        }
        std::cerr << "LOCAL n=" << n << " valid_vectors=" << vectors
                  << " signatures=" << signatures << "\n";
    }
};

// ─────────────────────────────────────────────
// 2. 牌グループどうしの最小距離を合成する
//    数牌・字牌の対称性による倍率は使用しない
// ─────────────────────────────────────────────
D combine(D const& a, D const& b) {
    D r;
    r.fill(INF);

    for (int ha = 0; ha < 2; ha++) {
        for (int hb = 0; ha + hb < 2; hb++) {
            for (int ma = 0; ma < 5; ma++) {
                for (int mb = 0; ma + mb < 5; mb++) {
                    r[5 * (ha + hb) + ma + mb] = std::min<int>(
                        r[5 * (ha + hb) + ma + mb],
                        a[5 * ha + ma] + b[5 * hb + mb]
                    );
                }
            }
        }
    }
    return r;
}

// ─────────────────────────────────────────────
// 3. 三麻・四麻の13枚／14枚を全数集計
//    distance = シャンテン数 + 1
// ─────────────────────────────────────────────
void count(std::string mode, int h, std::array<Table const*, 4> ts) {
    auto start = std::chrono::steady_clock::now();
    std::array<U, 8> counts{};
    std::array<std::array<U, 15>, 3> shapes{};
    uint64_t products = 0;

    for (int n0 = 0; n0 <= h; n0++) {
        for (int n1 = 0; n0 + n1 <= h; n1++) {
            for (int n2 = 0; n0 + n1 + n2 <= h; n2++) {
                int n3 = h - n0 - n1 - n2;

                for (auto const& a : ts[0]->groups[n0]) {
                    for (auto const& b : ts[1]->groups[n1]) {
                        auto ab = combine(a.sig.d, b.sig.d);

                        for (auto const& c : ts[2]->groups[n2]) {
                            auto abc = combine(ab, c.sig.d);

                            for (auto const& d : ts[3]->groups[n3]) {
                                int normal = INF;
                                for (int head = 0; head <= 1; head++) {
                                    for (int m = 0; m <= 4; m++) {
                                        normal = std::min(
                                            normal,
                                            abc[head * 5 + m] +
                                            d.sig.d[(1 - head) * 5 + 4 - m]
                                        );
                                    }
                                }

                                int kinds = a.sig.kinds + b.sig.kinds +
                                            c.sig.kinds + d.sig.kinds;
                                int pairs = a.sig.pairs + b.sig.pairs +
                                            c.sig.pairs + d.sig.pairs;
                                int term = a.sig.term + b.sig.term +
                                           c.sig.term + d.sig.term;
                                int tp = a.sig.tp + b.sig.tp +
                                         c.sig.tp + d.sig.tp;

                                // 通常形・七対子・国士無双の最小距離
                                int distance = std::min({
                                    normal,
                                    7 - pairs + std::max(0, 7 - kinds),
                                    14 - term - (tp > 0)
                                });
                                assert(distance >= 0 && distance < 8);

                                U weight = U(a.w) * b.w * c.w * d.w;
                                counts[distance] += weight;
                                shapes[0][normal] += weight;
                                shapes[1][7 - pairs + std::max(0, 7 - kinds)] += weight;
                                shapes[2][14 - term - (tp > 0)] += weight;
                                products++;
                            }
                        }
                    }
                }
            }
        }
    }

    // 標準出力：mode, tiles, shanten, count
    U total = 0;
    for (int i = 0; i < 8; i++) {
        total += counts[i];
        std::cout << mode << "," << h << "," << i - 1 << ","
                  << dec(counts[i]) << "\n";
    }
    assert(total == choose(mode == "three" ? 108 : 136, h));

    // 補助出力：役種別の距離分布をshapes.csvに保存
    std::ofstream sf("shapes.csv", std::ios::app);
    for (int kind = 0; kind < 3; kind++) {
        for (int d = 0; d < 15; d++) {
            if (shapes[kind][d]) {
                sf << mode << "," << h << "," << kind << "," << d - 1
                   << "," << dec(shapes[kind][d]) << "\n";
            }
        }
    }
    std::cerr << "COUNT " << mode << " " << h << " total=" << dec(total)
              << " ordered_products=" << products << " seconds="
              << std::chrono::duration<double>(
                     std::chrono::steady_clock::now() - start).count()
              << "\n";
}

// ─────────────────────────────────────────────
// 4. 実行：三麻／四麻 × 13枚／14枚の4ケース
// ─────────────────────────────────────────────
int main(int argc, char** argv) {
    auto start = std::chrono::steady_clock::now();
    std::ofstream("shapes.csv")
        << "mode,tiles,shape_0normal_1seven_2orphans,shanten,count\n";

    Table suit(9, true, false);
    Table honor(7, false, true);
    Table man(2, false, true);

    // 引数に出力先を指定すると、局所距離テーブルをバイナリ保存
    if (argc > 1) {
        std::string path = argv[1];
        for (auto t : {&suit, &honor, &man}) {
            std::ofstream f(
                path + "/table" + std::to_string(t->n) + ".bin",
                std::ios::binary
            );
            for (auto const& a : t->dist) {
                f.write((char const*)a.data(), a.size());
            }
        }
    }

    std::cout << "mode,tiles,shanten,count\n";
    for (int h : {13, 14}) {
        count("three", h, {&man, &suit, &suit, &honor});
        count("four", h, {&suit, &suit, &suit, &honor});
    }

    std::cerr << "ALL seconds="
              << std::chrono::duration<double>(
                     std::chrono::steady_clock::now() - start).count()
              << "\n";
}