linalg API

linalg パッケージは、Luna-Flow/linear-algebra のイミュータブルなベクトル上の関数の勾配とヤコビ行列を計算します。各ドライバーは、入力座標を 1 つずつ、双対数のベクトル上で関数を評価します。手法は linalg の設計 で説明しています。

ソース:src/linalg/linalg.mbt。

インポート

import {
  "Luna-Flow/autodiff/linalg",
  "Luna-Flow/linear-algebra/immut" @la,
}

例では、このパッケージに @linalg を、linear-algebra/immut に @la を使います。@linalg.Vector と @linalg.Matrix は @la.Vector と @la.Matrix と同じ型です。

関数の形状と事前条件

ドライバーは第 1 引数に関数を、第 2 引数に点 x∈Tnx \in T^n を取ります:

ドライバーf結果
gradient(Vector[Dual[T]]) -> Dual[T]∇f(x)∈Tn\nabla f(x) \in T^n
value_and_gradient(Vector[Dual[T]]) -> Dual[T](f(x),∇f(x))(f(x), \nabla f(x))
jacobian(Vector[Dual[T]]) -> Vector[Dual[T]]Jf(x)∈Tm×nJ_f(x) \in T^{m \times n}
value_and_jacobian(Vector[Dual[T]]) -> Vector[Dual[T]](f(x),Jf(x))(f(x), J_f(x))

4 つすべてが次を要求します:

  • f は長さ nn = x.length() のベクトルを受け取り、nn 未満のインデックスのみを読む。
  • ベクトル値の f は、どの呼び出しでも同じ長さ mm を返す。
  • f はキャプチャしたすべての値を定数(Dual::constant)として扱う。

ドライバーは形状を検証しません。違反するとインデックスエラーで中断するか、後の呼び出しがより長いベクトルを返した場合は余分な要素が黙って無視されます。

勾配

gradient

スカラー関数の勾配 ∇f(x)=(∂f/∂x1,…,∂f/∂xn)\nabla f(x) = \big(\partial f/\partial x_1, \dots, \partial f/\partial x_n\big) を計算します。

pub fn[T : @luna-generic.One + @luna-generic.Zero] gradient((@immut.Vector[@dual.Dual[T]]) -> @dual.Dual[T], @immut.Vector[T]) -> @immut.Vector[T]

各 ii について、座標 ii を Dual::new(x[i], 1)、その他の各座標を Dual::constant(x[j]) としてシードして f を 1 回呼び出し、結果の接成分を第 ii 成分として格納します。コスト:双対数上での f の評価 nn 回。n=0n = 0 のとき結果は空で、f は呼び出されません。

value_and_gradient

f(x)f(x) と ∇f(x)\nabla f(x) を返します。

pub fn[T : @luna-generic.One + @luna-generic.Zero] value_and_gradient((@immut.Vector[@dual.Dual[T]]) -> @dual.Dual[T], @immut.Vector[T]) -> (T, @immut.Vector[T])

値はすべての接成分を 0 にした追加の 1 回の呼び出しから得られるため、コストは n+1n + 1 回の評価です。

test "gradient of x0^2 + x0 x1" {
  let x = @la.Vector::from_array([2.0, 3.0])
  let f = (v : @la.Vector[@autodiff.Dual[Double]]) => v[0] * v[0] + v[0] * v[1]
  assert_eq(@linalg.gradient(f, x), @la.Vector::from_array([7.0, 2.0]))
  let (value, grad) = @linalg.value_and_gradient(f, x)
  assert_eq(value, 10.0)
  assert_eq(grad, @la.Vector::from_array([7.0, 2.0]))
}

ヤコビ行列

jacobian

ベクトル値関数のヤコビ行列 Jf(x)J_f(x) を計算します。行は出力ごと、列は入力ごとに 1 つです:

Jf(x)ij=∂fi∂xj(x),0≤i<m,0≤j<n.J_f(x)_{ij} = \frac{\partial f_i}{\partial x_j}(x), \qquad 0 \le i < m, \quad 0 \le j < n .
pub fn[T : @luna-generic.One + @luna-generic.Zero] jacobian((@immut.Vector[@dual.Dual[T]]) -> @immut.Vector[@dual.Dual[T]], @immut.Vector[T]) -> @immut.Matrix[T]

まずすべての接成分を 0 にして f を呼び出し mm を求め、次に入力 jj ごとに座標 jj をシードして 1 回ずつ呼び出します。呼び出し jj の出力の接成分が第 jj 列になります。コスト:n+1n + 1 回の評価。結果は row() =m= m、col() =n= n です。

value_and_jacobian

f(x)f(x) と Jf(x)J_f(x) を返します。

pub fn[T : @luna-generic.One + @luna-generic.Zero] value_and_jacobian((@immut.Vector[@dual.Dual[T]]) -> @immut.Vector[@dual.Dual[T]], @immut.Vector[T]) -> (@immut.Vector[T], @immut.Matrix[T])

値は接成分 0 での独自の呼び出しから、ヤコビ行列は jacobian から得られるため、コストは n+2n + 2 回の評価です。

test "jacobian of (x0 + x1, x0 x1, x0^2)" {
  let x = @la.Vector::from_array([2.0, 3.0])
  let f = (v : @la.Vector[@autodiff.Dual[Double]]) => {
    @la.Vector::from_array([v[0] + v[1], v[0] * v[1], v[0] * v[0]])
  }
  let (value, j) = @linalg.value_and_jacobian(f, x)
  assert_eq(value, @la.Vector::from_array([5.0, 6.0, 4.0]))
  assert_eq(j.row(), 3)
  assert_eq(j.col(), 2)
  assert_eq(j[1][0], 3.0) // d(x0 x1)/dx0 = x1
  assert_eq(j[2][1], 0.0) // d(x0^2)/dx1
}

再エクスポートされる型

Dual

双対数型です。dual API を参照してください。

pub using @dual {type Dual}

Vector

linear-algebra のイミュータブルなベクトルです。

pub using @immut {type Vector}

Matrix

linear-algebra のイミュータブルな行列です。

pub using @immut {type Matrix}