// Copyright (C) 2026 Kiyotsugu Arai
// SPDX-License-Identifier: LGPL-3.0-or-later

// IntCombinatorics.hpp
// 組合せ論的関数（階乗、二項係数など）for Int

#ifndef CALX_INT_COMBINATORICS_HPP
#define CALX_INT_COMBINATORICS_HPP

#include "IntBase.hpp"

namespace calx {

/**
 * @brief 組合せ論的関数
 *
 * 階乗、二項係数などの組合せ論で使用される関数を提供します。
 */
class IntCombinatorics {
public:
    /**
     * @brief 階乗 n!
     *
     * n! = n × (n-1) × ... × 2 × 1 を計算します。
     * ベンチマーク (2026-02-14) に基づき最適なアルゴリズムを選択:
     *   - N <= 100: テーブル参照 (O(1))
     *   - 100 < N <= 228: 単純ループ (テーブルをベースに順次乗算)
     *   - N > 228: 篩法 (素因数分解 + トーナメント乗算)
     * 注: トーナメント法は全範囲で篩法に劣るため自動選択では使用しない。
     * 詳細: TODO/NOTES_FACTORIAL_BENCHMARK.md
     *
     * @param n 非負整数
     * @return n! の値
     *
     * 特殊ケース:
     *   factorial(0) = 1（定義により）
     *   factorial(1) = 1
     *   factorial(n < 0) = NaN（負の数の階乗は未定義）
     *   factorial(NaN) = NaN
     *   factorial(∞) = ∞
     */
    static Int factorial(const Int& n);

    /**
     * @brief 二重階乗 n!!
     *
     * n!! = n × (n-2) × (n-4) × ... を計算します。
     * - n が偶数: n × (n-2) × ... × 4 × 2
     * - n が奇数: n × (n-2) × ... × 3 × 1
     *
     * @param n 非負整数
     * @return n!! の値
     *
     * 特殊ケース:
     *   doubleFactorial(0) = 1（定義により）
     *   doubleFactorial(1) = 1
     *   doubleFactorial(-1) = 1（定義により）
     *   doubleFactorial(n < -1) = NaN
     *   doubleFactorial(NaN) = NaN
     *   doubleFactorial(∞) = ∞
     *
     * 例:
     *   6!! = 6 × 4 × 2 = 48
     *   7!! = 7 × 5 × 3 × 1 = 105
     */
    static Int doubleFactorial(const Int& n);

    /**
     * @brief 二項係数 C(n, k) = n! / (k! × (n-k)!)
     *
     * 組合せの数 "n 個から k 個選ぶ場合の数" を計算します。
     *
     * @param n 非負整数
     * @param k 非負整数（0 <= k <= n）
     * @return C(n, k) の値
     *
     * 特殊ケース:
     *   binomial(n, 0) = 1（何も選ばない場合は1通り）
     *   binomial(n, n) = 1（全て選ぶ場合は1通り）
     *   binomial(n, k) = 0（k > n の場合）
     *   binomial(n, k < 0) = NaN
     *   binomial(n < 0, k) = NaN
     *   binomial(NaN, k) = NaN
     *
     * 最適化:
     *   - C(n, k) = C(n, n-k) を利用して計算量を削減
     *   - オーバーフロー回避のため、乗算と除算を交互に実行
     *
     * 例:
     *   binomial(5, 2) = 5!/(2!×3!) = 10
     *   binomial(10, 3) = 10!/(3!×7!) = 120
     */
    static Int binomial(const Int& n, const Int& k);

    /**
     * @brief フィボナッチ数 F(n)
     *
     * 高速倍加法 (fast doubling) で O(log n) 回の乗算で計算:
     *   F(2k)   = F(k) * (2*F(k+1) - F(k))
     *   F(2k+1) = F(k)^2 + F(k+1)^2
     *
     * @param n 非負整数
     * @return F(n)
     *
     * 特殊ケース:
     *   fibonacci(0) = 0, fibonacci(1) = 1
     *   fibonacci(n < 0) = NaN
     */
    static Int fibonacci(const Int& n);

    /**
     * @brief リュカ数 L(n)
     *
     * L(n) = F(n-1) + F(n+1) = 2*F(n+1) - F(n) の関係を利用。
     *
     * @param n 非負整数
     * @return L(n)
     *
     * 特殊ケース:
     *   lucas(0) = 2, lucas(1) = 1
     *   lucas(n < 0) = NaN
     */
    static Int lucas(const Int& n);

    /**
     * @brief プリモリアル primorial(n) = n 以下の素数の積
     *
     * GMP mpz_primorial_ui 相当。
     * エラトステネス篩で n 以下の素数を列挙し、バイナリ分割積で計算。
     *
     * @param n 非負整数 (int)
     * @return n 以下の全素数の積。n < 2 → 1。
     *
     * 例: primorial(10) = 2×3×5×7 = 210
     */
    static Int primorial(int n);

    /**
     * @brief m-階乗 mfac(n, m) = n · (n-m) · (n-2m) · ...
     *
     * GMP mpz_mfac_uiui 相当。
     * m=1 → n!, m=2 → n!! (doubleFactorial) に帰着。
     *
     * @param n 非負整数
     * @param m 正整数 (ステップ幅)
     * @return n · (n-m) · (n-2m) · ... (各因子が 1 以上の間)
     *
     * 例: mfac(10, 3) = 10×7×4×1 = 280
     */
    static Int mfac(int n, int m);

    // 個別アルゴリズム (テスト・ベンチマーク用に公開)

    /** テーブル参照: N <= 100 */
    static const Int& factorialTable(int n);

    /** 単純ループ: テーブル(100!)をベースに順次乗算 */
    static Int factorialSimple(int n);

    /** トーナメント法: 2の冪を分離 + 端×端バランス乗算 (ベンチマーク用、自動選択では未使用) */
    static Int factorialTournament(int n);

    /** 篩法: エラトステネス風素因数分解 + トーナメント乗算 */
    static Int factorialSieve(int n);

private:
    static constexpr int TABLE_MAX = 100;
    static constexpr int SIMPLE_MAX = 228;
};

} // namespace calx

#endif // CALX_INT_COMBINATORICS_HPP