【分割統治法】マスター法で常勝！Strassen、マージソートなどの計算量証明【アルゴリズムイントロダクション】


# 導入

## 本記事で目指すこと
本記事は『**アルゴリズムイントロダクション第4版**』**第4章**「**分割統治**」の輪読会用資料です。
40ページある章を60分の輪読会で説明するため、的を絞った紹介となります。
今回は「**Strassenのアルゴリズムを理解し、その計算量をマスター法で示す**」ことを一旦のゴールとします。

## 構成

まず、マージソートを例に**分割統治とアルゴリズム的漸化式の復習**をします。
次に、**正方行列乗算**の3つのアルゴリズムを概観し、**Strassenのアルゴリズム**を理解します。
その後、分割統治型の漸化式を解く強力な手法である**マスター法**をマスターします。
最後に、マスター法を使ってマージソートとStrassenのアルゴリズムの計算量を証明します。

Appendixとして、より汎用的な手法である**Akra-Bazzi法**の紹介や、アルゴリズムイントロダクションの章末問題の解を掲載しています。

## 本記事で目指さないこと

『アルゴリズムイントロダクション第4版』第4章「分割統治」の全ては説明しません。特に、
- 置き換え法を使う際の上手な推定や証明のコツ
- 再帰木を用いて計算量を見積るテクニック
- 連続マスター定理の証明

など、一部重要な箇所をバッサリ省略しています。

## 前提
『アルゴリズムイントロダクション第4版』第3章までの知識を説明なく用います。
特に、計算量の漸近記法と、マージソートの理解を仮定します。

https://rsk0315.hatenablog.com/entry/2021/10/13/235627

https://algo-logic.info/merge-sort/




# 分割統治法の復習

分割統治法は、**大きな問題をより小さい問題に分割し**、分解したそれぞれの問題を（直接高速に解けるサイズまで）再帰的に更に分割し、**小さい問題の解を用いて最初の大きな問題を解く手法**です。

以下の用語が用いられます。

**基底段階**（base case）：サイズが十分に小さい時。再帰なしに直接問題を解く
**再帰段階**（recursive case）：再帰の各レベルで以下の3つの段階を適用して問題を解く

再帰段階では、再帰の各レベルで以下の処理を行います。

- **分割**（divide）：問題をいくつかの部分問題に分ける
    - これらの部分問題は、元の問題と同じであり、かつサイズが元の問題より小さい
- **統治**（conqure）：部分問題を再帰的に解く
- **結合**（combine）：部分問題の解を組み合わせて元の問題の解を得る


```分割統治法の例：マージソート.rs
fn merge_sort(arr: &mut [i32]) {
    /* 基底段階 */
    if arr.len() <= 1 {
        return;
    }
    /* 再帰段階 */
    // 分割 Θ(1)
    let mid = arr.len() / 2;
    let (left, right) = arr.split_at_mut(mid);
    // 統治
    merge_sort(left);
    merge_sort(right);
    // 結合 Θ(n)
    merge(arr, mid);
}
```
[マージソートのコード全体はこちら](https://github.com/BinomialSheep/rust-algo-intro/blob/main/ch4/marge_sort/src/main.rs)

「分割」の実行時間は$\Theta(1)$、「結合」の実行時間は$\Theta(n)$なので、マージソートの実行時間は下記のように表すことができます。

```math
T(n) =
\begin{cases}
  \Theta(1) & n \leq 1 \\
  T(\lfloor n/2 \rfloor) + T(\lceil n/2 \rceil) + \Theta(n) & n > 1
\end{cases}
```

今後は基底段階では$\Theta(1)$であることを前提とし、省略します。また、$\lfloor n/2 \rfloor$と$\lceil n/2 \rceil$の差はたかだか1であり、基本的に無視できるため、常に$n/2$と書くことにします（Akra-Bazzi 漸化式の節で再度扱います）
つまり、今後は単に$$T(n) = 2T(n/2) + \Theta(n)$$と書くことにします（$2T(n/2)$の係数$2$は無視できません）
この解が$T(n) = \Theta(n \log n)$であることは有名で、数学的帰納法や再帰木を使っても証明できますが、本記事で紹介するマスター法で瞬殺できます。


# 正方行列乗算

行数と列数が等しい行列を正方行列といい、正方行列同士の掛け算を正方行列乗算と言います。
定義通りの正方行列乗算の実行時間は$\Theta(n^3)$ですが、Strassenのアルゴリズムというより高速な手法が知られています。

## 定義通りの乗算

$n\times n$型の正方行列$A = \lbrace a_{ij} \rbrace, B = \lbrace a_{ij} \rbrace$の積を$C=A・B$とすると、すべての$i, j = 1,2, ..., n$に対して、その要素$c_{ij}$は、次のように定義されます。

```math
c_{ij} = \sum_{k=1}^n a_{ik} b_{kj}
```

プログラム上は行列を2次元配列で表現することにします。3重ループで定義通りに計算することができ、その計算量は$\Theta(n^3)$です。

```
fn naive_matrix_multiply(matrix_a: &SquareMatrix, matrix_b: &SquareMatrix) -> SquareMatrix {
    let size = matrix_a.size;
    // Θ(n^2)
    let mut result = vec![vec![0; size]; size];
    // Θ(n^3)
    for i in 0..size {
        for j in 0..size {
            for k in 0..size {
                result[i][j] += matrix_a.data[i][k] * matrix_b.data[k][j];
            }
        }
    }
    SquareMatrix {
        size: size,
        data: result,
    }
}
```

## 単純な分割統治法

分割統治法を用いて行列積$A・B$を計算する方法を考えます。
1つの$n\times n$型行列を（n>1である限り）4つの$n/2 \times n/2$行列に分割して、分割した各行列について再帰的に計算してから結合することで、求めたかった$n\times n$型行列を求めることができます。簡単のため、ここではnは2の冪乗であると仮定します。


```math
A =
\begin{pmatrix}
A_{11} & A_{12} \\
A_{21} & A_{22}
\end{pmatrix}, \quad
B =
\begin{pmatrix}
B_{11} & B_{12} \\
B_{21} & B_{22}
\end{pmatrix}, \quad
C =
\begin{pmatrix}
C_{11} & C_{12} \\
C_{21} & C_{22}
\end{pmatrix}
```

に分割することができ、行列積$A×B$は
```math
\begin{align}
\begin{pmatrix}
C_{11} & C_{12} \\
C_{21} & C_{22}
\end{pmatrix}
&=
\begin{pmatrix}
A_{11} & A_{12} \\
A_{21} & A_{22}
\end{pmatrix}
\begin{pmatrix}
B_{11} & B_{12} \\
B_{21} & B_{22}
\end{pmatrix} \notag \\
&=
\begin{pmatrix}
A_{11} \cdot B_{11} + A_{12} \cdot B_{21} & A_{11} \cdot B_{12} + A_{12} \cdot B_{22} \\
A_{21} \cdot B_{11} + A_{22} \cdot B_{21} & A_{21} \cdot B_{12} + A_{22} \cdot B_{22}
\end{pmatrix}
\end{align}
```

と計算できます。
つまり、$n/2 \times n/2$型の8つの乗算と、4つの加算です。
乗算については再帰的に分割統治法を適用し、加算については愚直に$\Theta(n^2)$で計算すると、このアルゴリズムの漸化式は
$$T(n) = 8T(n/2) + \Theta(n^2)$$
となります。これを後述のマスター定理で解くと、$T(n)=\Theta(n^3)$であり、この分割統治法では定義通りの計算と漸近計算量は変わりません。

```
// 再帰を使ったΘ(n^3)の行列乗算
fn matrix_multiply_recursive(
    matrix_a: &SquareMatrix,
    matrix_b: &SquareMatrix,
    result: &mut SquareMatrix,
    n: usize,
) {
    // 基底段階
    if n == 1 {
        result.data[0][0] += matrix_a.data[0][0] * matrix_b.data[0][0];
        return;
    }

    // 分割
    let (a11, a12, a21, a22) = matrix_a.split();
    let (b11, b12, b21, b22) = matrix_b.split();

    let half = n / 2;
    let mut c11 = SquareMatrix::zero(half);
    let mut c12 = SquareMatrix::zero(half);
    let mut c21 = SquareMatrix::zero(half);
    let mut c22 = SquareMatrix::zero(half);

    // 統治
    matrix_multiply_recursive(&a11, &b11, &mut c11, half);
    matrix_multiply_recursive(&a12, &b21, &mut c11, half);
    matrix_multiply_recursive(&a11, &b12, &mut c12, half);
    matrix_multiply_recursive(&a12, &b22, &mut c12, half);
    matrix_multiply_recursive(&a21, &b11, &mut c21, half);
    matrix_multiply_recursive(&a22, &b21, &mut c21, half);
    matrix_multiply_recursive(&a21, &b12, &mut c22, half);
    matrix_multiply_recursive(&a22, &b22, &mut c22, half);

    // 結合
    result.merge(&c11, &c12, &c21, &c22);
}
```

## 分割統治法（Strassenのアルゴリズム）

さて、いよいよStrassenのアルゴリズムを紹介します。
Strassenのアルゴリズムは先ほどの$\Theta(n^3)$の分割統治法から更に乗算を減らすことで、$\Theta(n^3)$よりも高速に（つまり小文字のオーを使って$o(n^3)$で）行列積を求めることができる驚くべきアルゴリズムです。

Strassenのアルゴリズムのアイディアは、先ほどの分割統治が「8回の乗算と4回の加算」だったのに対して、数式を弄って「7回の乗算と18回の加算」に変形してしまうというものです。
行列乗算は行列加算よりも重い計算なので、乗算を1回減らすために14回加算を増やしてもお釣りが来るんですね。


```math
\begin{array}{ll}
S_1 = B_{12} - B_{22} & S_6 = B_{11} + B_{22} \\
S_2 = A_{11} + A_{12} & S_7 = A_{12} - A_{22} \\
S_3 = A_{21} + A_{22} & S_8 = B_{21} + B_{22} \\
S_4 = B_{21} - B_{11} & S_9 = A_{11} - A_{21} \\
S_5 = A_{11} + A_{22} & S_{10} = B_{11} + B_{12} \\
\\
P_1 = A_{11} \cdot S_1 & P_5 = S_5 \cdot S_6 \\
P_2 = S_2 \cdot B_{22} & P_6 = S_7 \cdot S_8 \\
P_3 = S_3 \cdot B_{11} & P_7 = S_9 \cdot S_{10} \\
P_4 = A_{22} \cdot S_4\\
\\
C_1 = P_5 + P_4 - P_2 + P_6 \\
C_2 = P_1 + P_2 \\
C_3 = P_3 + P_4 \\
C_4 = P_5 + P_1 - P_3 - P_7
\end{array}
```

こんなものをどうやって思いつくかはさておき、丁寧に$C_i$の式を追えば正しいことは容易に確かめられます。

このアルゴリズムの漸化式は
$$T(n) = 7T(n/2) + \Theta(n^2)$$
となります。マスター定理で解くと、$T(n)=\Theta(n^{\log7})=o(n^{2.81})$であり、これは定義通りの　$\Theta(n^3)$よりも高速なアルゴリズムです。

```strassen.rs
// StrassenのアルゴリズムΘ(n^log7)
fn strassen(matrix_a: &SquareMatrix, matrix_b: &SquareMatrix, result: &mut SquareMatrix, n: usize) {
    // 基底段階
    if n == 1 {
        result.data[0][0] += matrix_a.data[0][0] * matrix_b.data[0][0];
        return;
    }
    // 分割
    let (a11, a12, a21, a22) = matrix_a.split();
    let (b11, b12, b21, b22) = matrix_b.split();

    let half = n / 2;
    // 統治
    // 加算パート
    let s1 = b12.clone() - b22.clone();
    let s2 = a11.clone() + a12.clone();
    let s3 = a21.clone() + a22.clone();
    let s4 = b21.clone() - b11.clone();
    let s5 = a11.clone() + a22.clone();
    let s6 = b11.clone() + b22.clone();
    let s7 = a12.clone() - a22.clone();
    let s8 = b21.clone() + b22.clone();
    let s9 = a11.clone() - a21.clone();
    let s10 = b11.clone() + b12.clone();
    // 乗算パート
    let mut p1 = SquareMatrix::zero(half);
    let mut p2 = SquareMatrix::zero(half);
    let mut p3 = SquareMatrix::zero(half);
    let mut p4 = SquareMatrix::zero(half);
    let mut p5 = SquareMatrix::zero(half);
    let mut p6 = SquareMatrix::zero(half);
    let mut p7 = SquareMatrix::zero(half);
    strassen(&a11, &s1, &mut p1, half);
    strassen(&s2, &b22, &mut p2, half);
    strassen(&s3, &b11, &mut p3, half);
    strassen(&a22, &s4, &mut p4, half);
    strassen(&s5, &s6, &mut p5, half);
    strassen(&s7, &s8, &mut p6, half);
    strassen(&s9, &s10, &mut p7, half);
    // 加算パート
    let c11 = p5.clone() + p4.clone() - p2.clone() + p6.clone();
    let c12 = p1.clone() + p2.clone();
    let c21 = p3.clone() + p4.clone();
    let c22 = p5.clone() + p1.clone() - p3.clone() - p7.clone();

    // 結合
    result.merge(&c11, &c12, &c21, &c22);
}
```
[3種の正方行列乗算のコード全体](https://github.com/BinomialSheep/rust-algo-intro/blob/main/ch4/matrix_multiply/src/main.rs)


# アルゴリズム的漸化式を解く4つの方法

『アルゴリズムイントロダクション』では、アルゴリズム的漸化式を解く方法として以下の4つが紹介されています。
1. 置き換え法（substitution method）
    - 推測し、その推測が正しいことを数学的帰納法を用いて証明する
    - 最も一般的に使えるが、良い推測と数学的な証明が要求される
2. 再帰木法（recursion-tree method）
    - 漸化式を再帰木として可視化する
    - 厳密さに欠けるが、直感的な議論や、置き換え法のための推定を得るのに役立つ
3. マスター法
    - $T(n)=aT(n/b)+f(n)$型の漸化式に対して適用できる
    - 適用可能な場合に一番簡単な手法かつ、実際頻繁に適用できる
4. Akra-Bazzi法
    - マスター法よりも複雑な漸化式に適用できる
    - マスター法よりも難しい

1, 2はよくある話なので省略します（書籍にはちょっとしたテクニックや落とし穴の回避も解説されていてオススメです）
本記事では主にマスター法を紹介し、AppendixでAkra-Bazzi法を扱います。

# マスター法

マスター法は以下の形の漸化式（マスター漸化式）を、マスター定理を適用して解く手法です。

$$T(n) = aT(n/b) + f(n)$$

$a>0$と$b>1$は定数であり、「サイズ$n$の問題を、サイズ$n/b$を持つ$a$個の部分問題に分割して解く分割統治法」を表しています。
$f(n)$は駆動関数と呼び、分割統治法の分割や結合に掛かるコストと対応しています。


「Strassenのアルゴリズム」の実行時間の漸化式は

$$T(n) = 7T(n/2) + \Theta(n^2)$$

なので、$a=7, b=2, f(n)=\Theta(n^2)$のマスター漸化式です。

## マスター定理
マスター漸化式
$$T(n) = aT(n/b) + f(n)$$
は、以下の3ケースに当てはまる場合、$T(n)$の漸近解が直ちに求まります。

##### ケース1
$f(n)=O(n^{\log_b a - \varepsilon})$を満たす定数$\varepsilon > 0$が存在するならば、$T(n)=\Theta(n^{\log_b a})$
##### ケース2
$f(n)=\Theta (n^{\log_b a \log_2^k n})$を満たす定数$k≧0$が存在するならば、$T(n)=\Theta(n^{\log_b a} \log_2^{k+1}n)$

##### ケース3
$f(n)=\Omega(n^{\log_b a + \varepsilon}) $を満たす定数$\varepsilon > 0$が存在し、かつ$f(n)$が次の正則条件を満たすならば、$T(n)=\Theta(f(n))$

**正則条件**：ある定数$c<1$と十分大きな全ての$n$に対して、$af(n/b)≦cf(n)$

マスター定理において何度も出てくる$$n^{\log_b a}$$は分水界関数（water-shed function）と呼ばれます。


### 直感的なイメージ

マスター漸化式
$$T(n) = aT(n/b) + f(n)$$
は、統治の$aT(n/b)$に対応する分水界関数$n^{\log_b a}$と、分割や結合に対応する駆動関数の$f(n)$の比較で3パターンに分類できます。
**パターンA**：分水界関数が駆動関数よりも漸近的に速く増加するパターン
**パターンB**：分水界関数と駆動関数が同じくらいの速度で増加するパターン
**パターンC**：分水界関数よりも駆動関数が漸近的に速く増加するパターン

パターンAの場合、分水界関数がボトルネックなので、$T(n)=\Theta(n^{log_b a})$になりそうです。
パターンCの場合、駆動関数がボトルネックなので、$T(n)=\Theta(f(n))$になりそうです。
パターンBの場合は同じくらいなので、例えば$T(n)=\Theta(n^{\log_b a} \log_2 n)$になったりしそうです。（※ケース2でk=0の場合に相当）

ただしこれはあくまで直感的なイメージで、不正確です。
たとえばケース1は実際には、分水界関数が駆動関数よりも**多項式的に**速く漸近的に増加しなければならないことを要請しています（「マスター法が適用できない例」の節で説明）

### マスター法の適用例
#### １．Strassenのアルゴリズム
$$T(n) = 7T(n/2) + \Theta(n^2)$$

分水界関数は$n^{\log_2 7}$、駆動関数は$\Theta(n^2)$です。$\log_2 7$は約2.81なので、適当に$\varepsilon = 0.1$などとすると、$\Theta(n^2)=O(n^{\log_2 7 - 0.1})$より、ケース1が適用できます。

> **ケース1（再掲）**
> $f(n)=O(n^{\log_b a - \varepsilon})$を満たす定数$\varepsilon > 0$が存在するならば、$T(n)=\Theta(n^{\log_b a})$

よって、漸化式の解は$$T(n)=\Theta(n^{\log_b a})=\Theta(n^{\log_2 7})$$




#### ２．マージソート
$$T(n) = 2T(n/2) + \Theta(n)$$

分水界関数は$n^{\log_2 2} = n$、駆動関数は$\Theta(n)$です。マスター漸化式において$a=2, b=2$であり、$k=0$とするとマスター定理のケース2における$f(n)=\Theta (n^{\log_b a \log_2^k n})$を満足するので、ケース2が適用できます。

> **ケース2（再掲）**
> $f(n)=\Theta (n^{\log_b a \log_2^k n})$を満たす定数$k≧0$が存在するならば、$T(n)=\Theta(n^{\log_b a} \log_2^{k+1}n)$


よって、漸化式の解は$$T(n)=\Theta(n^{\log_b a} \log_2^{k+1}n)=\Theta(n \log_2 n)$$

#### ３. ケース3になるやつ

$$T(n) = 3T(n/4) + n \log_2 n$$
について考えます。

分水界関数は$n^{\log_4 3}$、駆動関数は$n \log_2 n$です。$\log_4 3 < 0.793$なので、$\varepsilon = 0.2$などとすれば$n \log_2 n=\Omega(n^{\log_4 3 + 0.2})$を満たします。
また、適当に$c=3/4$などを代入することで、正則条件も満たすので、ケース3が適用できます。
<details><summary>正則条件を満たすことの確認</summary>

$a=3, b=4, f(n)=n \log_2 n$なので、$af(n/b) = (3n/4)\log_2 (n/4)$.
また、$c=3/4$とすると$cf(n)=(3n/4)\log_2 n$.
</details>

> **ケース3（再掲）**
> $f(n)=\Omega(n^{\log_b a + \varepsilon}) $を満たす定数$\varepsilon > 0$が存在し、かつ$f(n)$が次の正則条件を満たすならば、$T(n)=\Theta(f(n))$
> 正則条件（再掲）：ある定数$c<1$と十分大きな全ての$n$に対して、$af(n/b)≦cf(n)$

よって、漸化式の解は
$$T(n) = \Theta(n \log_2 n)$$


より多くの例については練習問題4.5-1や章末問題を参照してください。

### マスター法が適用できない例

#### マスター漸化式の形をしていない場合
マスター法は、「サイズ$n/b$を持つ$a$個の部分問題に分割して解く」タイプの問題にのみ適用できます。$a,b$が定数でない場合や、分割サイズが均等でない場合には適用できません。
たとえばクイックソートは、ピボットをランダムに選んだ場合、分割のバランスが常に$n/2$になるとは限りません。分割割合が$cn$と$(1-c)n$とすると、漸化式は次のような形になります。

$$T(n)=T(cn)+T((1-c)n)+O(n)$$

これはマスター法では解けず、Akra-Bazzi法が必要となる典型的な例です。

https://note.com/cute_orchid39/n/nf1490085101f

#### 分水界関数と駆動関数の漸近的な比較ができない場合

無限個の$n$の値に対して$f(n)\gg n^{\log_b a}$が成り立ち、また異なる無限個の$n$の値に対して$f(n)\ll n^{\log_b a}$が成り立つような場合です。実際のアルゴリズムで出現する大部分の駆動関数は、分水界関数と比較可能です。

#### ケース1とケース2の"隙間"の場合
> **ケース1（再掲）**
> $f(n)=O(n^{\log_b a - \varepsilon})$を満たす定数$\varepsilon > 0$が存在するならば、$T(n)=\Theta(n^{\log_b a})$
> **ケース2（再掲）**
> $f(n)=\Theta (n^{\log_b a \log_2^k n})$を満たす定数$k≧0$が存在するならば、$T(n)=\Theta(n^{\log_b a} \log_2^{k+1}n)$


ケース1は、分水界関数が駆動関数よりも多項式的に速く増加する場合に適用できます。**対数多項式的**に速く増加する場合にはケース1は適用できず、ケース2も適用できません。

$$T(n)=2T(n/2)+n/\log_2 n$$

がその例です。分水界関数$n^{\log_2 2} = n$は駆動関数よりも速く増加しますが、任意の$\varepsilon > 0$に対して$n/\log_2 n=O(n^{1 - \varepsilon})$が成り立たず、ケース1は適用できません。
また、f(n)=$\Theta (n /\log_2 n)$であり、これはケース2で$k=-1$の場合に相当しますが、$k$は非負が条件なので、ケース2も適用できません。


同様に、ケース2とケース3の間にも隙間があります。

※ちなみにこの解は$\Theta(n\log_2
\log_2 n)$です（アルゴリズムイントロダクション練習問題4.6-3）


# Appendix
## Akra-Bazzi 法

Akra-Bazzi法はマスター法よりも複雑な分割統治型漸化式の求解に使われる手法です。
以下のアルゴリズム的漸化式のクラスをAkra-Bazzi漸化式と言います。

$$T(n) = f(n) + \sum_{i=1}^{k} a_i T\left(\frac{n}{b_i}\right)$$

ここで、$k$は正の定数であり、すべての定数$a_1, a_2,...,a_k \in \mathbb{R} $は正の実数、すべての定数$b_1, b_2,...,b_k \in \mathbb{R} $は1より大きい実数です。駆動関数$f(n)$は十分大きい非負の実数上で定義され、関数値も非負です。
難しいですが、分割統治法において分割後のサイズが異なるようなアルゴリズムの実行時間の漸化式も記述できるよということです。

### フロア関数とシーリング関数を無視していい条件

Akra-Bazzi法ではフロア関数、シーリング関数を一般には無視できませんが、駆動関数$f(n)$が下記の多項式的増加条件（polynomial-growth condition）を満たす場合は、無視することができます。

十分大きなすべての正の実数上で定義された関数$f(n)$において、次の命題を満足するような定数$\hat{n}$が存在するとき、「関数$f(n)$は多項式的増加条件を満たす」といいいます。

命題：どの定数$\phi ≧1$に対しても、すべての$1≦\psi≦\phi$および$n≧\hat{n}$に対して、$f(n)/d≦f(\psi n)≦df(n)$を満足する定数d>1が存在する

なお、これは必要十分条件ではなく、多項式的増加条件を満たさない場合でも、フロア関数やシーリング関数を無視しても解に影響しないことはあります。

### Akra-Bazzi法
#### ステップ０
必要なら、多項式的増加条件を満たすことを確認します。
#### ステップ１
$\sum_{i=1}^{k}\frac{a_i}{b_i^p}=1$
を満足するような実数$p$を見つけます。左辺の値は、$p\to -\infty$で$\infty$であり、$p$が増加すれば単調減少し、$p\to \infty$で$0$となるので、等式を満たす$p$は常に存在します。

※補足
$b_i>1$なので、各$a_i(\frac{1}{b_i})^p$が上記の性質を満たし、式全体としても常に$p$は存在する。

#### ステップ２
このとき、漸化式の解は
$$T(n) = \Theta \left( n^p \left( 1 + \int_{1}^{n} \frac{f(x)}{x^{p+1}} dx \right) \right)$$

になります（ほよ～？）

### 適用例

$$T(n)=T(n/5)+T(7n/10)+n$$
について考えます（これは[線形最悪時間選択アルゴリズム](https://ja.m.wikipedia.org/wiki/中央値の中央値)の計算量で、『アルゴリズムイントロダクション』では9章で扱います）

> Akra-Bazzi漸化式（再掲）
> $$T(n) = f(n) + \sum_{i=1}^{k} a_i T\left(\frac{n}{b_i}\right)$$

解きたい漸化式は、$a_1=a_2=1, b_1=5, b_2=7/10, f(n)=n$に場合のAkra-Bazzi漸化式です。

#### ステップ０
フロア関数、シーリング関数はないので省略します。

#### ステップ１
$$\left(\frac{1}{5}\right)^p+\left(\frac{7}{10}\right)^p=1$$
を満たすような$p$を決定します。
[二分探索で求まり](https://chatgpt.com/share/679a3bb8-ffbc-800e-83d2-a31f2b2ac629)、$p=0.83978..$なのですが、Akra-Bazzi法で解を得る上では、$p$に$0$と$1$を代入して　$0<p<1$が分かれば十分です。

#### ステップ２

```math
\begin{align}
T(n) &= \Theta \left( n^p \left( 1 + \int_{1}^{n} \frac{f(x)}{x^{p+1}} dx \right) \right) \\
&= 
\Theta \left( n^p \left( 1 + \int_{1}^{n} \ {x^{-p}} dx \right) \right)　f(x)=xを代入した \\
&= 
\Theta \left( n^p \left( 1 + \left[ \frac{x^{1-p}}{1-p} \right]_1^n \right) \right)　-p \ne -1より\\
&= 
\Theta \left( n^p \left( 1 + \left(\frac{n^{1-p}}{1-p} - \frac{1}{1-p}\right)\right)\right)\\
&= 
\Theta (n^p \cdot\Theta(n^{1-p}))\\
&= 
\Theta (n)
\end{align}
```

お疲れ様でした。



## Strassen周りの補足

### なぜこんなものが思いつくのか……

https://joisino.hatenablog.com/entry/2017/09/03/200000
### より漸近的に高速な正方行列乗算アルゴリズム

Strassenが1969年に発表したあとも研究は続けられており、Gallのアルゴリズム（$O(n^{2.372863})$）などが知られています。Gall先生は東大で博士号を取得したあとも日本で研究を続けており、現在は名古屋大の教授です。

[京都大学大学院 情報学研究科通信情報システム専攻 François Le Gall 特定准教授インタビュー－量子計算の考え方を用いて世界最速の“行列のかけ算”を実現－](https://www.nistep.go.jp/wp/wp-content/uploads/NISTEP-STIH4-1-00112.pdf)

### で、実際速いのか
速くないです。実際、NumPyなどのライブラリの行列乗算で用いられていません。実用的に高速な行列乗算については、「BLASのGEMM (General Matrix-Matrix Multiply) ルーチン」などでGoogleやChatGPTに訊いてみてください。


## 分割統治で高速に解ける有名問題

### 最近点対問題
https://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_5_A&lang=ja

### 高速フーリエ変換（FFT）

https://ja.tech.jar.jp/cellular/fft3.html

### カラツバ法
多倍長乗算。式変形して分割乗算の回数を減らすと高速になるというメカニズムはStrassenに似ています。

https://www3.cuc.ac.jp/~miyata/mathprg/poly/karatsuba.html

### 最大部分配列問題
分割統治法を用いて$\Theta(n\log_2 n)$で解けますが、Kadane's Algorithmという$\Theta(n)$で解けるアルゴリズムの方が有名です。
『アルゴリズムイントロダクション 第3版』には載っていたようで、
[公式サイト](https://mitpress.mit.edu/9780262046305/introduction-to-algorithms/)のResources>Additional Material>Materials Removed from 3eの「The maximum-subarray problem.pdf」から英語で読めます。日本語では『珠玉のプログラミング』に載っていた記憶です。

https://ark4rk.hatenablog.com/entry/2018/01/08/002508

### 木上の分割統治
重心分解して木のサイズを半分以下にして解くやつです。

https://qiita.com/drken/items/4b4c3f1824339b090202

## 問題回答集

『アルゴリズムイントロダクション 第4版』の問題回答。


### 練習問題4.1
#### 4.1-1
ちょうど2べきになるように行列のサイズを拡張し、0を埋める。
この操作を行ってもサイズはたかだか2倍にしかならないので、漸近計算量は変わらない。

### 練習問題4.2

#### 4.2-3
漸化式は
$$T(n) = kT(n/3) + \Theta(n^2)$$
となる。$o(n^{\log_2 7})$となる最大の$k$について考えるので、マスター法のケース1を考えればよい。分水界関数は$n^{\log_3 k}$。
$\log_3 21<\log_2 7<\log_3 22$なので、k=21が最大であり、その実行時間は$\Theta(n^{\log_3 21})=\Theta(n^{1+\log_3 7})$

#### 4.2-4
マスター法で、それぞれの実行時間は$\Theta(n^{\log_{68} 143640}), \Theta(n^{\log_{70} 132464}), \Theta(n^{\log_{72} 155424})$となります。

| 分割数 | 乗算回数 | 指数|
|--------------|---------------|--------------|
| 68           | 132,464       | 2.795128     |
| 70           | 143,640       | **2.795123** |
| 72           | 155,424       | 2.795147     |
| 2 (Strassen) | 7             | 2.807355     |

#### 4.2-5
実
$P_1=(a+b)(c+d), P_2=ac, P_3=bd$を求めておく。
実数部：$P_2-P_3 = ac - bd$
虚数部：$P_1-P_2-P_3 = (ac+ad+bc+bd)-ac-bd=ad+bc$
より求まった。

#### 4.2-6
行列$A,B$の積を求めたいとする。
行列サイズを2倍にして、以下のような平方を計算してから右下の領域を取り出せばよい。

```math
\begin{pmatrix}
O & B \\
A & O
\end{pmatrix}
\begin{pmatrix}
O & B \\
A & O
\end{pmatrix}
=
\begin{pmatrix}
BA & O \\
O & AB
\end{pmatrix}
```

### 練習問題4.5
#### 4.5_1
マスター法を用いて以下の漸化式の解を求める。

分水界関数はいずれも
$$n^{\log_4 2}=\sqrt{n}$$

であり、駆動関数と比較して適切なケースを適用するだけ。

##### a. T(n)=2T(n/4)+1
$$T(n)=\Theta(\sqrt{n})$$

##### b. T(n)=2T(n/4)+√n
$$T(n)=\Theta(\sqrt{n}\log_2 n)$$

##### c. T(n)=2T(n/4)+√n (log_2 n)^2
$$T(n)=\Theta(\sqrt{n}\log_2^3 n)$$

##### d T(n)=2T(n/4)+n
$$T(n)=\Theta(n)$$

##### e T(n)=2T(n/4)+n^2
$$T(n)=\Theta(n^2)$$



#### 4.5_2
4.2-3と同様。
$\log_4 49=\log_2 7$なので、$a=48$が最大。


#### 4.7-5
Akra-Bazzi法を用いて次の漸化式の解を求めよ。

> Akra-Bazzi法（再掲）
$\sum_{i=1}^{k}\frac{a_i}{b_i^p}=1$
を満足するような実数$p$を見つけます。
このとき、漸化式の解は
$$T(n) = \Theta \left( n^p \left( 1 + \int_{1}^{n} \frac{f(x)}{x^{p+1}} dx \right) \right)$$
になります

##### a.

$$T(n)=T(n/2)+T(n/3)+T(n/6)+n \log_2 n$$

$p=1$です。

```math
\begin{align}
T(n) &= \Theta \left( n^p \left( 1 + \int_{1}^{n} \frac{f(x)}{x^{p+1}} dx \right) \right) \\
&= 
\Theta \left( n \left( 1 + \int_{1}^{n} \frac{x \log_2 x}{x^2} dx \right) \right)　f(x)=x\log_2 x, p=1を代入した \\
&= 
\Theta \left( n \left( 1 + \int_{1}^{n} \frac{\log_2 x}{x} dx \right) \right) \\
&= 
\Theta \left( n \left( 1 + \left[ \frac{(\log_2 x)^2}{2} \right]_1^n \right) \right)　置換積分（t=\log_2 x）\\
&= 
\Theta \left( n \left( 1 + \frac{(\log_2 n)^2}{2}\right)\right)　\log_2 1は0なので \\
&= 
\Theta (n \log_2^2 n)
\end{align}
```




##### b.

$$T(n)=3T(n/3)+8T(n/4)+n^2/\log_2 n$$

##### c.

$$T(n)=(2/3)T(n/3)+(1/3)T(2n/3)+\log_2 n$$

##### d.

$$T(n)=(1/3)T(n/3)+1/n$$

##### e.

$$T(n)=3T(n/3)+3T(2n/3)+n^2$$




### 章末問題4-3 変数変換
#### a
$$S(m)=2S(m/2)+\Theta(m)$$
#### b
マスター法ケース2より
$$S(m)=\Theta(m \log_2 m)$$

#### c
$$T(n)=S(\log_2 n)=\Theta(\log_2 n \log_2\log_2 n)$$

### 章末問題4-5 フィボナッチ数
有名問題かつ本章とあまり関係がないので略。

#### a
$$z+zF(z)+z^2F(z)=(1+F_0)z+(F_1+F_0)z^2+(F_1+F_0)z^3+...$$

フィボナッチ数の定義より$F_{n+2}=F_{n+1}+F_n$なので、$z^2$以降はよい。さらに$F_0=0, F_1=1$を考えれば、$$F(z)=0+z+F_2z^2+F_3z^3$$
なので、一致している。

#### b
1つ目はaから従う。2つ目は解の公式。3つ目は部分分数分解。

#### c

(b)の部分分数分解した各項をべき級数に展開する。

#### d
$|\hat{\phi}|<1$なので、iが大きくなると0に近づく。

#### e
帰納的に示せる

### 章末問題4-6 LSIチップのテスト

#### a
「少なくとも半分のチップが不良なことは分かるが、不良のチップの具体的な枚数は分からない」という問題の意図だと仮定する。
不良チップが常にウソを付けば$\binom{n}{2}$通り試しても「お互いにお互いを正常と報告するチップの集合」が2つに分かれるだけで、どちらの集合が正常なチップの集合か分からない。
（※単純に常にウソを付くだけだと、「集合のサイズが大きい方が不良チップの集合！」と言えるハックがあるが、不良チップは共謀して上手く集合の数やサイズを調整できることに留意。例えば正常チップの集合とサイズを合わせればよい）

#### b
$n$が偶数の場合を考える。テキトーに好きなもの同士で二人組を作らせて、$n/2$回テストを行う。
1. 互いに正常と報告した場合、両方正常か両方不良
2. 一方を不良と報告した場合、少なくとも一方は不良
3. 互いに異常と報告した場合、少なくとも一方は不良

なので、結果が2か3のペアを取り除いたとしても、「正常枚数>不良枚数」は維持される。
つまり1のペアの集合だけで「正常枚数>不良枚数」は維持される。
ここで1のペアを別々の集合に入れると、それぞれの集合は「正常枚数>不良枚数」であり、集合のサイズは高々$n/2$になるので、示された。

奇数の場合も少し考えると、正常と不良どちらが余ったとしても、偶数の場合と同戦略を取ったあとで片方のグループに足せば、そのグループは条件を満たすことが分かる。

#### c

bで得た元のサイズの半分以下の集合のうち1つに対して再帰的に適用すればよい。
$O(\log n)$回で「正常枚数>不良枚数」のまま3枚とか4枚になるので、最後は全て試して多数決をする。

$T(n) = T(n/2) + \Theta(n)$を解いて$T(n) = \Theta(n)$で1つの正常なチップが求まる。

#### d

cで得た1枚の正常なチップで$n-1$回テストをすれば全ての正常なチップを特定できるので、「1つの特定」「1つ特定した後の全ての特定」がそれぞれ$\Theta(n)$で、全体も$\Theta(n)$で求まる。



### 章末問題4-7 モンジュ行列
Monge行列と、各行の最左の最小要素を求める分割統治アルゴリズム。有名問題。
区間DP高速化など、暖色レベルの高度典型に応用例がいろいろあるらしい。

https://lorent-kyopro.hatenablog.com/entry/2021/04/04/133958

https://hackmd.io/@tatyam-prime/monge1