# 機能設計書 14-組み合わせ数学

## 概要

本ドキュメントは、Julia言語における組み合わせ数学機能の設計を記述する。factorial, binomial, gcd, lcm, gcdx, invmod, powermod等の組み合わせ論・整数論関数を提供する。

### 本機能の処理概要

組み合わせ数学機能は、整数論と組み合わせ論の基本的な関数群を提供する。階乗・二項係数・最大公約数・最小公倍数・拡張ユークリッド互除法・モジュラ逆元・冪剰余など、暗号学や数学計算で頻繁に使用される関数を効率的に実装している。

**業務上の目的・背景**：暗号アルゴリズムの実装（RSA、楕円曲線暗号等）、組み合わせ最適化、統計計算（二項係数）、数論研究などで整数論関数は不可欠である。Juliaの型システムを活用し、固定長整数（Int64/Int128）向けの高速実装とBigInt向けのGMP実装を透過的に提供する。

**機能の利用シーン**：（1）暗号鍵生成でのgcd/invmod/powermod、（2）確率・統計計算でのfactorial/binomial、（3）有理数の約分でのgcd、（4）信号処理のFFTサイズ最適化でのnextprod、（5）順列・置換操作でのisperm/invperm/permute!。

**主要な処理内容**：
1. 階乗（factorial）：ルックアップテーブルベースの高速計算（Int64: 20!, Int128: 34!まで）
2. 最大公約数（gcd）：二進GCD（Stein's algorithm）による高速計算
3. 最小公倍数（lcm）：gcdベースの計算（checked_abs/checked_mulでオーバーフロー検出）
4. 拡張ユークリッド互除法（gcdx）：ベズー係数の計算
5. モジュラ逆元（invmod）：gcdxベースの計算
6. 冪剰余（powermod）：効率的なmodular exponentiation
7. 順列操作（isperm, invperm, permute!, invpermute!）
8. 次の積因数（nextprod）

**関連システム・外部連携**：BigInt向けの実装ではGMPライブラリ（gmp.jl）の関数を使用。

**権限による制御**：特になし。

## 関連画面

| 画面No | 画面名 | 関連種別 | 関連する操作・処理 |
|--------|--------|----------|------------------|
| 2 | Juliaプロンプト（julia>） | 参照画面 | REPL上で整数論関数を対話的に実行 |

## 機能種別

計算処理（整数論・組み合わせ論関数）

## 入力仕様

### 入力パラメータ

| パラメータ名 | 型 | 必須 | 説明 | バリデーション |
|-------------|-----|-----|------|---------------|
| n | Integer | Yes | factorial/binomialの引数 | factorial: 0 <= n <= 20(Int64)/34(Int128)。BigIntは制限なし |
| a, b | Integer/Rational | Yes | gcd/lcm/gcdxの引数 | - |
| n, m | Integer | Yes | invmodの引数 | gcd(n,m)==1かつm!=0 |
| x, p, m | Integer | Yes | powermodの引数 | - |
| v | AbstractVector | Yes | isperm/invperm/permute!の引数 | permute!は1-basedインデックス |

### 入力データソース

Julia式からの直接入力。

## 出力仕様

### 出力データ

| 項目名 | 型 | 説明 |
|--------|-----|------|
| factorial_result | Integer | 階乗の計算結果 |
| gcd_result | Integer/Rational | 最大公約数 |
| lcm_result | Integer/Rational | 最小公倍数 |
| gcdx_result | Tuple{Integer,Integer,Integer} | (gcd, u, v) ベズー係数 |
| invmod_result | Integer | モジュラ逆元 |
| powermod_result | Integer | 冪剰余結果 |

### 出力先

演算結果はメモリ上に返却。

## 処理フロー

### 処理シーケンス

```
1. gcd(a, b) - 二進GCDアルゴリズム
   └─ trailing_zeros計算 → abs → 交互剰余+シフト → 結果にkビットシフト
2. factorial(n) - テーブルルックアップ
   └─ 範囲チェック → _fact_table64/128から直接取得
3. gcdx(a, b) - 拡張ユークリッド互除法
   └─ divrem反復 → ベズー係数(s0,t0)更新 → 符号調整
4. invmod(n, m) - gcdxベースの逆元計算
   └─ gcdx(n, m) → gcd==1チェック → 符号調整
```

### フローチャート

```mermaid
flowchart TD
    A[gcd入力 a, b] --> B{a==0 or b==0?}
    B -->|Yes| C[abs of 非ゼロ側を返却]
    B -->|No| D{typemin チェック}
    D -->|a==typemin かつ a==b| E[OverflowError]
    D -->|OK| F[_gcd: 二進GCDアルゴリズム]
    F --> G[trailing_zeros計算]
    G --> H[abs取得, b右シフト]
    H --> I{a != 0?}
    I -->|Yes| J[a右シフト → absdiff → min/max更新]
    J --> I
    I -->|No| K[b << k を返却]
```

## ビジネスルール

### 業務ルール

| ルールNo | ルール名 | 内容 | 適用条件 |
|---------|---------|------|---------|
| BR-01 | 階乗テーブル上限 | Int64: 20!, UInt64: 20!, Int128: 33!, UInt128: 34! | factorial呼び出し時 |
| BR-02 | gcd正値保証 | gcdは常に非負値を返す（checked_abs使用） | 全gcd呼び出し |
| BR-03 | lcmゼロルール | 引数にゼロがあればlcmは0 | lcm(0, x) |
| BR-04 | gcdxベズー | u*a + v*b = gcd(a,b) を満たすu,vを返す | gcdx呼び出し時 |
| BR-05 | invmod前提条件 | gcd(n,m)==1 かつ m!=0 が必要 | invmod呼び出し時 |
| BR-06 | 二進GCD性能 | BitInteger向けにStein's algorithmを使用（約1.7x高速） | gcd(Int64,Int64)等 |

### 計算ロジック

- 二進GCD: trailing_zeros(a), trailing_zeros(b)を計算し、k=min(za,zb)。a,bを交互にabsdiffで更新し、最後にb<<kを返す
- 階乗テーブル: _fact_table64はTuple{20,Int64}、_fact_table128はTuple{34,UInt128}として事前計算
- 拡張GCD: divrem反復でs0*a0 + t0*b0 == a のループ不変条件を維持
- @assume_effects :terminates_locally でコンパイラの定数畳み込みを有効化（Rationalコンストラクタ最適化用）

## データベース操作仕様

該当なし

## エラー処理

### エラーケース一覧

| エラーコード | エラー種別 | 発生条件 | 対処方法 |
|------------|----------|---------|---------|
| - | DomainError | factorial(n)でn < 0 | 非負整数を指定 |
| - | OverflowError | factorial(n)でnがテーブル上限超過 | BigIntを使用 |
| - | OverflowError | gcd(typemin(T), typemin(T)) | BigIntを使用 |
| - | DomainError | invmod(n, 0) | m != 0を指定 |
| - | DomainError | invmod(n, m)でgcd(n,m)!=1 | 互いに素な値を指定 |
| - | ArgumentError | isperm/invpermで不正な置換 | 有効な置換を指定 |

### リトライ仕様

該当なし

## トランザクション仕様

該当なし

## パフォーマンス要件

- factorialはテーブルルックアップで O(1)
- 二進GCDはランダムInt64で約1.7x、Int128で約2.1x高速（従来のユークリッド互除法比）
- @assume_effects :terminates_locally による定数畳み込み最適化

## セキュリティ考慮事項

- checked_abs/checked_mulでオーバーフローを検出
- typemin(Signed)の取り扱いに特別な注意（abs(typemin(Int64))はオーバーフロー）

## 備考

- gcd/lcmはRational型引数もサポート（Julia 1.4以降）
- gcdxは3引数以上もサポート（Julia 1.12以降、再帰的に2引数版を呼び出し）
- nextprodはBase.DSPで使用されるためcombinatorics.jlに残置

---

## コードリーディングガイド

本機能を理解するために参照すべきファイルと、推奨する読み解き順序を以下に示す。

### 推奨読解順序

#### Step 1: データ構造を理解する

階乗テーブルと型制約を把握する。

| 順序 | ファイル | パス | 読解ポイント |
|-----|---------|------|-------------|
| 1-1 | combinatorics.jl | `base/combinatorics.jl` | _fact_table64/128（5-19行目）：事前計算された階乗テーブル |

**読解のコツ**: テーブルはVector→Tupleに変換して格納。getfieldでインデックスアクセスする。

#### Step 2: 整数論関数を理解する

gcd, lcm, gcdx, invmodの実装。

| 順序 | ファイル | パス | 読解ポイント |
|-----|---------|------|-------------|
| 2-1 | intfuncs.jl | `base/intfuncs.jl` | gcd(a, b)二進GCD（80-95行目）：Stein's algorithm |
| 2-2 | intfuncs.jl | `base/intfuncs.jl` | lcm(a, b)（136-144行目）：gcdベース計算 |
| 2-3 | intfuncs.jl | `base/intfuncs.jl` | gcdx(a, b)拡張ユークリッド互除法（232-257行目） |
| 2-4 | intfuncs.jl | `base/intfuncs.jl` | invmod(n, m)（286行目以降） |

**主要処理フロー**:
- **80-95行目**: _gcd: trailing_zeros→absシフト→交互absdiff→b<<k
- **232-257行目**: gcdx: 特殊ケース判定→divremループ→符号調整

#### Step 3: 階乗と組み合わせ関数を理解する

factorial, binomialの実装。

| 順序 | ファイル | パス | 読解ポイント |
|-----|---------|------|-------------|
| 3-1 | combinatorics.jl | `base/combinatorics.jl` | factorial_lookup（21-30行目）：テーブル参照 |
| 3-2 | combinatorics.jl | `base/combinatorics.jl` | factorial型別ディスパッチ（32-41行目） |

#### Step 4: 順列操作を理解する

isperm, invperm, permute!, invpermute!の実装。

| 順序 | ファイル | パス | 読解ポイント |
|-----|---------|------|-------------|
| 4-1 | combinatorics.jl | `base/combinatorics.jl` | isperm（76-101行目）：置換検証 |
| 4-2 | combinatorics.jl | `base/combinatorics.jl` | invperm（286-316行目）：逆置換計算 |
| 4-3 | combinatorics.jl | `base/combinatorics.jl` | permute!/invpermute!（214, 243行目） |
| 4-4 | combinatorics.jl | `base/combinatorics.jl` | nextprod（337-373行目）：次の積因数 |

### プログラム呼び出し階層図

```
gcd(a, b) [BitInteger版]
    |
    +-- checked_abs(a/b)       [ゼロチェック]
    +-- _gcd(a, b)             [二進GCDアルゴリズム]
            +-- trailing_zeros(a/b)
            +-- abs(a/b)
            +-- absdiff(a, b)  [差の絶対値]

factorial(n) [Int64版]
    |
    +-- factorial_lookup(n, _fact_table64, 20)
            +-- getfield(table, idx)

gcdx(a, b) [Integer版]
    |
    +-- divrem(a, b)           [反復]
    +-- isnegative(a) → abs/符号反転

invmod(n, m)
    |
    +-- gcdx(n, m)
    +-- mod(result, m)
```

### データフロー図

```
[入力]                [処理]                         [出力]

Int64 a, b    ─────> _gcd (二進GCD)         ─────> Int64 gcd
Int64 n       ─────> _fact_table64[n]       ─────> Int64 n!
Int64 a, b    ─────> gcdx (拡張ユークリッド) ─────> (gcd, u, v)
Int64 n, m    ─────> invmod (gcdxベース)    ─────> Int64 逆元
BigInt a, b   ─────> MPZ.gcd (GMP)          ─────> BigInt gcd
BigInt n      ─────> MPZ.fac_ui (GMP)       ─────> BigInt n!
```

### 関連ファイル一覧

| ファイル | パス | 種別 | 役割 |
|---------|------|------|------|
| combinatorics.jl | `base/combinatorics.jl` | ソース | factorial, isperm, invperm, permute!, nextprod |
| intfuncs.jl | `base/intfuncs.jl` | ソース | gcd, lcm, gcdx, invmod, powermod |
| gmp.jl | `base/gmp.jl` | ソース | BigInt向けgcd/factorial/binomial/powermod/gcdx/invmod |
| checked.jl | `base/checked.jl` | ソース | checked_abs, checked_mul（オーバーフロー検出） |
