core 設計

このページでは Complex[T] の代数を導出し、各 luna-generic インスタンスが法則を満たす条件を述べ、ジェネリック型の背後にある選択(可変フィールド、スケーリングしない除算の式、浮動小数点解析との分離)を説明します。

設計目標

エコシステムのあらゆるスカラーに対する単一のジェネリックな複素数型であり、そのインスタンスは構成が持つ代数構造をちょうど過不足なく表明し、IEEE の特殊値、分岐切断、超越関数はバックエンドパッケージに委ねます。

数学的背景

構成

可換環 $R$ に対し、 $R$ 上の複素数は多項式環を $x^2 + 1$ で生成されるイデアルで割った剰余環です:

R[i] = R[x]/(x^2 + 1), \qquad i = x + (x^2 + 1), \qquad i^2 = -1 .

モニック多項式 $x^2 + 1$ で割ると次数 2 未満の剰余が一意に残るため、すべての元は一意な $a, b \in R$ によって $a + bi$ と表されます。これがペア (re, im) です。多項式として計算し $i^2 = -1$ で簡約すると、このパッケージの演算が得られます:

\begin{aligned} (a + bi) \pm (c + di) &= (a \pm c) + (b \pm d)i, \\ (a + bi)(c + di) &= ac + (ad + bc)i + bd\,i^2 = (ac - bd) + (ad + bc)i . \end{aligned}

可換環の剰余環として、 $R[i]$ は零元 $0 + 0i$ 、単位元 $1 + 0i$ を持つ可換環であり、 $R$ は $a \mapsto a + 0i$ によって埋め込まれます。

共役とノルム

$\overline{a + bi} = a - bi$ は $x \mapsto -x$ から誘導される写像で、イデアル $(x^2 + 1)$ を保つため、位数 2 の環自己同型です:

\overline{z + w} = \bar z + \bar w, \qquad \overline{zw} = \bar z\,\bar w, \qquad \bar{\bar z} = z .

ノルム $N(z) = z\bar z = (a + bi)(a - bi) = a^2 + b^2$ は $R$ に属し、乗法的です: $N(zw) = zw\,\bar z\bar w = N(z)N(w)$ 。

逆元と除算

$N(z)$ が $R$ で可逆ならば $z \cdot \bar z\,N(z)^{-1} = 1$ なので、

z^{-1} = \frac{\bar z}{N(z)} = \frac{a}{a^2 + b^2} - \frac{b}{a^2 + b^2}\,i, \qquad \frac{a + bi}{c + di} = \frac{(ac + bd) + (bc - ad)i}{c^2 + d^2} .

逆に、 $z$ が可逆ならば $N(z)N(z^{-1}) = N(1) = 1$ なので、 $N(z)$ は可逆です。したがって $z$ が単元であるのは、 $a^2 + b^2$ が単元であるときに限ります。

構成が体になる場合

$K$ を体とします。 $K[x]/(x^2 + 1)$ が体になるのは、 $x^2 + 1$ が $K$ 上既約であるとき、すなわち $K$ において $-1$ が平方数でないときに限ります。直接確かめると、 $N(z) = a^2 + b^2 = 0$ となるゼロでない $z$ は $b \ne 0$ でなければならず、そのとき $(a/b)^2 = -1$ です。逆に $s^2 = -1$ ならば、 $(s + i)(s - i) = s^2 + 1 = 0$ が零因子を与えます。

$K = \mathbb R$ (Float、Double)では $-1$ は平方数ではなく、 $\mathbb R[i] = \mathbb C$ は体です。
$K = \mathbb C$ (Complex[Double])では $-1 = i^2$ は平方数なので、Complex[Complex[Double]] $\cong \mathbb C[j]/(j^2 + 1) \cong \mathbb C \times \mathbb C$ は零因子を持つ環です。例えば $(1 + ij)(1 - ij) = 1 - i^2 j^2 = 0$ です。

設計上の決定

インスタンスは `T` の構造に従う

問題。 Complex[T] はどの luna-generic トレイトを実装してよいか?

選択。 Zero、One、AddMonoid、AddGroup は成分ごとに作用するため、T の対応する演算だけを必要とします。MulMonoid、Semiring、Ring は T : Ring を必要とします。積が減算( $ac - bd$ )を使い、 $R[i]$ の環の法則が $R$ の環の法則を必要とするためです。Inverse、MulGroup、Field、Div 演算子は $N(z)^{-1}$ のために T : Field を必要とします。Conjugate は Neg だけを必要とします。

Field インスタンスはすべての T : Field に対して宣言されています。前節により、これが法則を満たすのは T において $-1$ が平方数でない場合で、実スカラー型はこれに該当します。T = Complex[Double] では該当せず、inv はゼロでない零因子に対して中断します。トレイト制約では「−1 は平方数でない」を表現できないため、このインスタンスはここで呼び出し側を信頼します。¹1 これは Luna Flow の「法則を満たすインスタンスだけを実装する」というルールにおける既知の欠落です。テストスイートでは、入れ子の複素数はノルムが可逆な値に対してのみ使っています。

4 回の乗算

選択肢。 4 回の乗算による教科書どおりの積 $(ac - bd) + (ad + bc)i$ 。あるいは Gauss の 3 回乗算形式 $k_1 = c(a + b)$ 、 $k_2 = a(d - c)$ 、 $k_3 = b(c + d)$ による積 $(k_1 - k_3) + (k_1 + k_2)i$ 。

選択。 4 回の乗算。浮動小数点では、教科書どおりの形式のノルム相対誤差は高々 $\sqrt5\,u$ です。²2 R. Brent, C. Percival and P. Zimmermann, “Error bounds on complex floating-point multiplication”, Mathematics of Computation 76 (2007). 一方、3 回乗算形式は $a + b$ や $d - c$ のような桁落ちしやすい和を加えるため、成分ごとの精度が劣ります。また、ジェネリックな T では乗算が加算より高コストとは限りません。

成分ごとに見ると、 $|\theta_2| \le \gamma_2 = 2u/(1 - 2u)$ として $\mathrm{fl}(ac - bd) = ac(1 + \theta_2) - bd(1 + \theta_2')$ なので、

|\mathrm{fl}(ac - bd) - (ac - bd)| \le \gamma_2\,(|ac| + |bd|),

これは $ac$ と $bd$ が打ち消し合わない限り、結果に比べて小さい値です。

ジェネリックなコアにおける教科書どおりの除算

問題。 除算には $N(w)^{-1}$ が必要です。浮動小数点では、 $c + di$ よりはるかに早く $c^2 + d^2$ がオーバーフローまたはアンダーフローします。

選択肢。 スケーリングを行うアルゴリズム(Smith の方法や 2 のべきによるスケーリング)は、ジェネリックな体にはない比較と絶対値を必要とします。教科書どおりの式は体の演算だけで済みます。

選択。 ジェネリックな Div と Inverse は、T の Inverse::inv を使った教科書どおりの式を使います。厳密な体では正確です。Double では $c^2 + d^2$ が範囲内にある限り正しく、Double の Inverse はゼロ(アンダーフローしたノルムを含む)で中断します。スケーリングを行い特殊値を考慮した除算は float_backend にあります。

可変フィールド

問題。 複素数の値はループ内で更新されることが多く(累積器、漸化式など)、ステップごとに新しい値を割り当てるのは無駄です。

選択。 Complex[T] は mut re と mut im を持つ pub(all) 構造体で、set、set_re、set_im を備えます。すべての演算は新しい値を返し、入力を変更しないため、セッターを呼ばないコードは複素数を値として扱えます。変更を行う呼び出し側は、構造体が参照で共有されることに注意する必要があります。

コアに解析関数を置かない

$\sqrt z$ 、 $\log z$ 、 $\sin z$ には、順序、絶対値、実スカラーの超越関数、分枝の選択、IEEE の特殊値が必要です。ジェネリックな体にはこれらのどれもありません。そのため、ルートパッケージは代数までにとどめ、float_backend パッケージが Double 向けの解析を提供します。

テキスト形式とデバッグ形式

Show は T のテキストを使って re + imi を出力します( $1 - 2i$ なら 1 + -2i)。これが型の正規のテキストであり、to_string が昇格されたメソッドとして残っている理由です。Debug はテストと診断のためにレコード形式を出力します。その他のトレイトメソッドの MoonBit 0.10 における昇格は src/extends.mbt で明示的に行っています。

正しさと不変条件

$R[i]$ は剰余環なので、T : Ring に対して環の法則が成り立ちます。テストスイートでは整数上で結合律、単位元、 $\overline{zw} = \bar z\bar w$ を確認しています。
厳密な T では、z * z.conjugate() は $N(z) + 0i$ にちょうど等しくなります。
z * w / w == z と z.inv() * z == 1 は、 $N(w)$ 、 $N(z)$ が単元であれば厳密な体で成り立ち、Double では範囲内であれば丸め誤差の範囲で成り立ちます。
どの演算も引数を変更しません。変更するのは 3 つのセッターだけです。

採用しなかった代替案

Double 専用の複素数型。 Float、整数(ガウス整数)、厳密な有理数のために型を重複させることになります。
不変フィールド。 値のセマンティクスはすっきりしますが、re と im をインプレースで更新している既存の呼び出し側にとっては破壊的変更になります。
コアでのスケーリング除算。 ジェネリックな体にはない、順序と大きさの演算が必要になります。

範囲外

解析関数や超越関数、特殊値の処理、分岐切断は扱いません。float_backend を参照してください。
コアには極形式の表現はありません。
T によって $R[i]$ が体になるかどうかは確認しません。Field インスタンスはそれを信頼します。
浮動小数点の T に対するオーバーフローに安全な除算はありません。

これは Luna Flow の「法則を満たすインスタンスだけを実装する」というルールにおける既知の欠落です。テストスイートでは、入れ子の複素数はノルムが可逆な値に対してのみ使っています。 ↩
R. Brent, C. Percival and P. Zimmermann, “Error bounds on complex floating-point multiplication”, Mathematics of Computation 76 (2007). ↩