# Einsum: numpy 風のテンソル縮約
:::{container} prog-cpp
Hi-QUBO は **`qbpp::einsum(subscript, arrays...)`** を提供します.
これは numpy の
[アインシュタイン縮約](https://en.wikipedia.org/wiki/Einstein_notation)
と同じ記法で,整数・変数・項・式の多次元配列を 1 行で縮約できる関数です.
テンプレート引数 `OutDim` は出力配列の次元数を指定し,
subscript の出力ラベル数と実行時に照合されます.
## C++ でなぜ `OutDim` が必要か
Hi-QUBO の多次元配列は **`Array`** という型で表現されます.
ここで次元 `Dim` は**テンプレート引数(コンパイル時定数)**であり,
`Array<1, Expr>` と `Array<2, Expr>` は **まったく別の型** です.
`einsum` のコードを生成する時点で,コンパイラは戻り値の型を確定する
必要があります.
ところが,出力の次元は subscript 文字列(例: `"ij,jk->ik"` なら 2 次元,
`"i,i->"` ならスカラー)から決まります.subscript は `const char*` 引数として
**実行時にしか中身を見ることができない** ため,**C++ コンパイラには出力次元を
推論する手段がありません**.
そのため,呼び出し側がテンプレート引数として明示的に `OutDim` を指定する
必要があります:
```{include} /../programFiles/markDown/advanced/Einstein-sum.md
:start-after:
:end-before:
```
`OutDim` と subscript の実際の出力ラベル数が一致しているかは実行時に検査され,
不一致ならエラー終了します.誤ったテンプレート引数が「形の違う配列」として
黙って通ってしまうことはありません.
Python 版でこの引数が不要なのは,Python のオブジェクトが次元情報を実行時に
保持しているためです.Python バインディング側で subscript を解析して,
正しい次元の出力配列を自動的に組み立てています.
## subscript の文法
```{include} /../programFiles/markDown/advanced/Einstein-sum.md
:start-after:
:end-before:
```
- 各 **label** は ASCII 1 文字(`,`・`-`・`>`・空白を除く)です.
- 各入力配列は次元数とちょうど同じ数のラベルを持つ必要があります.
- 入力に現れて出力に現れないラベルは **縮約(総和)** されます.
- 入力と出力の両方に現れるラベルは自由軸として残ります.
- **同一入力内に同じラベルが 2 度現れる** と,その 2 つの軸が結合されます
(trace や対角抽出に使います).
- 暗黙形式 `"ij,jk"`(`->` を省略)では,全入力中にちょうど 1 度だけ
現れるラベルをアルファベット順に並べたものが出力になります(numpy と同じ仕様).
- 右辺が空(`"i,i->"`)の場合は **スカラー出力**(`OutDim == 0`)になります.
## 出力型
- **すべての入力が整数配列**(`Array`)の場合,結果も
整数配列 `Array` になります.`OutDim == 0` のときは
`coeff_t` のスカラーが返ります.
- それ以外(`Var`, `Term`, `Expr` を 1 つでも含む)場合は
`Array` が返ります.`OutDim == 0` のときは `Expr` のスカラーです.
## 使用例
以下のプログラムは,`einsum` の代表的な使い方を示します.
```{literalinclude} /../programFiles/cppPrograms/advanced/Einstein-sum-program1.cpp
:language: cpp
:caption: Einstein-sum-program1.cpp
```
このプログラムは以下を出力します:
```{include} /../programFiles/markDown/advanced/Einstein-sum.md
:start-after:
:end-before:
```
### double フロントエンド
[double フロントエンド](../basic/variables-and-expressions.md#実数double係数)(`DOUBLE_TYPE*`)では `coeff_t` は `double` であり,
係数配列は実数値を保持します.`qbpp::array()` にも `double` 値を直接渡せます
(1次元・2次元の初期化子リストと `std::vector`).係数行列が数値計算の結果として
得られる場合 — 例えば FMQA 型のブラックボックス最適化で学習した代理モデルの二次形式を
ソルバーに渡す場合 — に便利です:
```{literalinclude} /../programFiles/cppPrograms/advanced/Einstein-sum-program2.cpp
:language: cpp
:caption: Einstein-sum-program2.cpp
```
同じ二次形式は二重 `for` ループでも構築できます.C++ では手書きループと `einsum` の
性能はほぼ同等ですが,`einsum` は縮約全体を 1 行で表現でき,コードの見通しが良くなります.
## 3 つ以上の入力
`einsum` は任意個数の入力配列を受け取れます.組合せ最適化での代表例は
**二次割当問題(QAP)** 形の目的関数
$\sum_{a,k,l} f_a\, d_{kl}\, x_{a,k}\, x_{a,l}$ です:
```{include} /../programFiles/markDown/advanced/Einstein-sum.md
:start-after:
:end-before:
```
この 1 行が 4 重の for ループを置き換えます.大規模な場合は内部で複数の
CPU スレッドにより並列に計算されます.
## どのような場面で使うか
目的関数や制約が「テンソル添字でインデックスされた積の総和」として書ける
場合は,`einsum` を使うのが最も簡潔です.明示的な多重ループに比べて,
- 数式構造を直接表現でき,
- 添字計算のミスを避けられ,
- 大規模配列では内部でマルチスレッド化されて高速です.
単純な全要素総和や軸ごとの総和には **`qbpp::sum()`** や
**`qbpp::vector_sum()`**([Sum 関数](../advanced/sum-functions-for-arrays.md)を参照)の方が直接的です.
複数の配列の積を取る,あるいはインデックス間の関係が複雑になった時点で
`einsum` への切り替えを検討してください.
:::
:::{container} prog-python
PyQBPP は **`qbpp.einsum(subscript, *arrays)`** を提供します.
これは numpy の
[アインシュタイン縮約](https://en.wikipedia.org/wiki/Einstein_notation)
と同じ記法で,整数・変数・項・式の多次元配列を 1 行で縮約できる関数です.
出力配列の次元は subscript から自動的に推論されます.
C++ 版の `qbpp::einsum(...)` のようにテンプレート引数で次元を
指定する必要はなく,Python 版は subscript と入力配列のみを渡します.
## subscript の文法
```{include} /../programFiles/markDown/advanced/Einstein-sum.md
:start-after:
:end-before:
```
- 各 **label** は ASCII 1 文字(`,`・`-`・`>`・空白を除く)です.
- 各入力配列は次元数とちょうど同じ数のラベルを持つ必要があります.
- 入力に現れて出力に現れないラベルは **縮約(総和)** されます.
- 入力と出力の両方に現れるラベルは自由軸として残ります.
- **同一入力内に同じラベルが 2 度現れる** と,その 2 つの軸が結合されます
(trace や対角抽出に使います).
- 暗黙形式 `"ij,jk"`(`->` を省略)では,全入力中にちょうど 1 度だけ
現れるラベルをアルファベット順に並べたものが出力になります(numpy と同じ仕様).
- 右辺が空(`"i,i->"`)の場合は **スカラー出力** になります.
## 出力型
- **すべての入力が整数配列**の場合,結果も整数配列になります.
出力次元が 0 のときは `int` のスカラーが返ります.
- それ以外(`Var`, `Term`, `Expr` を 1 つでも含む)場合は
`Expr` の配列が返ります.出力次元が 0 のときは `Expr` のスカラーです.
## 使用例
以下のプログラムは,`einsum` の代表的な使い方を示します.
```{literalinclude} /../programFiles/pythonPrograms/advanced/Einstein-sum-program1.py
:language: python
:caption: Einstein-sum-program1.py
```
このプログラムは以下を出力します:
```{include} /../programFiles/markDown/advanced/Einstein-sum.md
:start-after:
:end-before:
```
### double フロントエンドと numpy 入力
[double フロントエンド](../basic/variables-and-expressions.md#実数double係数)のモジュール(`pyqbpp.d`,`pyqbpp.dc64e64`,
`pyqbpp.dc128e128`)では,係数配列は Python の `float` 値を保持します.`qbpp.array()` には
float のリストに加えて **numpy の ndarray** を直接渡すこともでき,ネイティブコードで変換されます.
また,ndarray は(通常のリストと同様に)**`einsum` に直接渡す**こともでき,自動的に変換されます.
数値計算で得られた密な係数行列 — 例えば FMQA 型のブラックボックス最適化で学習した
代理モデルの二次形式 — を `einsum` に渡す最速の方法です:
```{literalinclude} /../programFiles/pythonPrograms/advanced/Einstein-sum-program2.py
:language: python
:caption: Einstein-sum-program2.py
```
同じ二次形式を二重 Python ループで構築すると項ごとにライブラリを呼び出すことになりますが,
`einsum` は縮約全体をネイティブコードで実行するため,密な行列では通常数倍高速です.
自動変換は**呼び出しのたびに**行われます(numpy 配列は書き換え可能なため,変換結果は
キャッシュされません).同じ係数行列を多数の `einsum` 呼び出しで使う場合は,
`qbpp.array(W)` で一度変換し,変換済み配列を再利用してください.
## 3 つ以上の入力
`einsum` は任意個数の入力配列を受け取れます.組合せ最適化での代表例は
**二次割当問題(QAP)** 形の目的関数
$\sum_{a,k,l} f_a\, d_{kl}\, x_{a,k}\, x_{a,l}$ です:
```{include} /../programFiles/markDown/advanced/Einstein-sum.md
:start-after:
:end-before:
```
この 1 行が 4 重の for ループを置き換えます.大規模な場合は内部で複数の
CPU スレッドにより並列に計算されます.
## どのような場面で使うか
目的関数や制約が「テンソル添字でインデックスされた積の総和」として書ける
場合は,`einsum` を使うのが最も簡潔です.明示的な多重ループに比べて,
- 数式構造を直接表現でき,
- 添字計算のミスを避けられ,
- 大規模配列では内部でマルチスレッド化されて高速です.
単純な全要素総和や軸ごとの総和には **`qbpp.sum()`** や
**`qbpp.vector_sum()`**([Sum 関数](../advanced/sum-functions-for-arrays.md)を参照)の方が直接的です.
複数の配列の積を取る,あるいはインデックス間の関係が複雑になった時点で
`einsum` への切り替えを検討してください.
:::