immut/sparse の設計

設計目標

SparsePolynomial[A] は多変数多項式を単項式から係数への順序付きマップとして格納します。アルゴリズムが多項式全体を走査するよりも「xαx^\alpha の係数は何か」を問い合わせることが多い場合に使う表現であり、イミュータブルな値のまま、その答えを対数時間で返します。

数学的背景

多項式は有限台の関数 f:N(∞)→Rf : \mathbb{N}^{(\infty)} \to R、α↦fα\alpha \mapsto f_\alpha です (core の設計 を参照)。その台 supp⁡f={α∣fα≠0}\operatorname{supp} f = \{ \alpha \mid f_\alpha \neq 0 \} は有限であり、ff は台への制限によって決まります。疎な多項式はまさにその制限、すなわち次の有限部分写像を格納します。

supp⁡f→R∖{0},α↦fα,\operatorname{supp} f \to R \setminus \{0\}, \qquad \alpha \mapsto f_\alpha ,

したがって係数の検索は部分写像の評価であり、存在しないキーは 00 を意味します。2 つの多項式が等しいのは、これらの部分写像が等しいときに限ります。

設計上の判断

単項式順序をキーとする平衡探索木

問題。 指数ベクトルによる検索には索引が必要です。ハッシュマップなら期待 O(1)O(1) の検索が得られますが順序がなく、ソート済み配列なら O(log⁡m)O(\log m) で検索できますが挿入は O(m)O(m) です。

選択。 項は @sorted_map.SortedMap[ExponentVector, A]、すなわち ExponentVector::compare で順序付けられた AVL 木に格納されます。AVL 木は高さを 1.44log⁡2(m+2)1.44 \log_2(m + 2) 未満に保つため、get のキー比較は O(log⁡m)O(\log m) 回で、各比較はベクトル長に対して O(ℓ)O(\ell) です。11 Adelson-Velsky と Landis、1962 年。高さの上界は、高さ hh の AVL 木の最小ノード数に関する Fibonacci 型の漸化式から従います。 マップは TermPolynomial と同じ単項式順序で順序付けられるので、反復すると正規形の項リストが得られます。ただしここでは 昇順 で、定数項が最初、先頭項が最後です。

木は他の演算のコスト特性にも合っています。mm 個のソート済みの項から構築するのは O(mlog⁡m)O(m \log m) で、ソートと同じです。

正規形の内容: 零値を持たない

不変条件。 格納される値はすべて非零です。構築では項表現の正規化 (ソート、等しいキーのマージ、和が零のものの除去) を再利用し、残ったものを挿入します。neg と scale は零と等しいと比較される結果をスキップします (零因子があると起こりえます)。キーは正規形の指数ベクトルで、値は非零なので、格納されたマップは supp⁡f\operatorname{supp} f 上の部分写像 α↦fα\alpha \mapsto f_\alpha そのもの であり、

p == q  ⟺  p.to_terms() == q.to_terms()  ⟺  p=q.\texttt{p == q} \iff \texttt{p.to\_terms() == q.to\_terms()} \iff p = q .

したがって get(α) が None を返すことは、ちょうど fα=0f_\alpha = 0 を意味します。

カプセル化によるイミュータビリティ

SortedMap はミュータブルな構造ですが、それを保持するフィールドはプライベートであり、このパッケージのどの関数も、多項式を返した後にマップを変更しません。add_term を含むすべての演算は新しいマップを構築します。これにより永続木なしで値セマンティクスが得られますが、その代償として項を 1 つ加えるだけでマップが O(mlog⁡m)O(m \log m) で再構築されます。項を 1 つずつ加えるワークロードでは mutable/sparse を使うべきです。その add_term_inplace は木を O(log⁡m)O(\log m) で更新します。

項表現と共有する算術

加算と乗算は項リスト (* では mnmn 個すべての積) を集め、TermPolynomial とまったく同じように共有の正規化を通してマップを再構築します。したがって両方の表現は同じ正規形の結果を計算し、コストも木への挿入の定数倍を除いて同じです。eval は同じ項ごとの評価準同型であり、同じアリティの前提条件を持ちます。

2 つの多変数表現

2 つの型は同じ多項式を記述し、アクセス経路だけが異なります。

演算TermPolynomialSparsePolynomial
xαx^\alpha の係数to_terms() の O(m)O(m) 走査O(log⁡m)O(\log m) get
先頭項最初の要素最後の要素
反復順序降順昇順
mm 個の項からの構築O(mlog⁡m)O(m \log m)O(mlog⁡m)O(m \log m)
項を 1 つ追加 (イミュータブル)+ 経由で O(mlog⁡m)O(m \log m)O(mlog⁡m)O(m \log m) add_term
+, *O((m+n)log⁡(m+n))O((m+n)\log(m+n)), O(mnlog⁡(mn))O(mn\log(mn))同じ
メモリフラットな配列項ごとに木のノード 1 つ

変換は to_terms() を経由し、これにより実装パッケージは互いに独立に保たれます。両者が出会うのはファサードと ContextPolynomial です。

正しさ / 不変条件

  • キー は正規形の ExponentVector 値、値 は非零です。
  • 等価性 は昇順の項リストを比較し、多項式の等価性と一致します。
  • 検索。 get(α) == Some(c) であるのは fα=c≠0f_\alpha = c \neq 0 のときに限り、None であるのは fα=0f_\alpha = 0 のときに限ります。get_checked は同じ関数です。
  • 一致。 同じ入力項に対して、SparsePolynomial と TermPolynomial は同じ項集合を保持し、+、-、*、pow、scale、eval の結果が一致します。
  • 値セマンティクス。 既存の SparsePolynomial を変更する公開関数はありません。

採用しなかった代替案

  • ハッシュマップによるストレージ。 期待 O(1)O(1) で検索できますが、反復順序が任意になるため、等価性、表示、先頭項のいずれにもソートが必要になります。
  • 永続 (パスコピー) 木。 コアライブラリの @immut/sorted_map を使えば、変更なしで add_term を O(log⁡m)O(\log m) にできます。このパッケージは代わりにミュータブルな SortedMap をプライベートフィールドの背後に置き、逐次的な構築は mutable/sparse に任せています。
  • ストレージフラグ付きの単一の多変数型。 これは ContextPolynomial が内部で行っていることです。添字で変数を指定するレベルでは、2 つの明示的な型にすることで、それぞれのコストモデルがシグネチャに現れます。

境界

  • Compare や Hash のインスタンスはありません。疎な多項式はマップのキーにできません。
  • add_term は逐次的ではなく、マップを再構築します。
  • 変数は位置ベースです。名前付き変数と代入は ContextPolynomial の担当です。
  • 除算や Gröbner 基底のアルゴリズムはありません。
  • 指数は UInt でオーバーフロー時にラップアラウンドします。係数は A の意味に従います。

Footnotes

  1. Adelson-Velsky と Landis、1962 年。高さの上界は、高さ hh の AVL 木の最小ノード数に関する Fibonacci 型の漸化式から従います。 ↩