行列の累乗によるトリボナッチ数の計算

ここでは、次のトリボナッチ数列を扱います。

$$ T_0=0,\quad T_1=0,\quad T_2=1,\qquad T_{n+3}=T_{n+2}+T_{n+1}+T_n. $$

コンパニオン行列を使うと、この漸化式を線形代数の演算の反復として 表せるため、高速な累乗計算を利用できます。

パターンマッチで遷移行列を組み立てる

第1行は現在の状態に含まれる三つの成分を足し合わせます。下副対角線は、 過去の値を一つ下の位置へ移します。

$$ A= \begin{pmatrix} 1&1&1\\ 1&0&0\\ 0&1&0 \end{pmatrix}. $$

ジェネレーターでは、互いに無関係な九つの成分を列挙する代わりに、 この二つの構造的なパターンを記述します。

def m : Integer := 3

def A : Matrix Integer :=
  generateTensor
    (\match as list integer with
      | [#1, _] -> 1
      | [$x, #(x - 1)] -> 1
      | _ -> 0)
    [m, m]
A
$\begin{pmatrix} 1 & 1 & 1 \\ 1 & 0 & 0 \\ 0 & 1 & 0 \\ \end{pmatrix}$

初期状態

ベクトル $B=(1,0,0)^\mathsf T$ には $(T_2,T_1,T_0)^\mathsf T$ が格納されています。

def B : Vector Integer :=
  generateTensor
    (\[x] -> if x = 1 then 1 else 0)
    [m]

def tribonacciState (n : Integer) : Vector Integer :=
  MV.* (M.power A n) B
B
$\begin{pmatrix} 1 \\ 0 \\ 0\\ \end{pmatrix}$

漸化式を進める

任意の $n\ge0$ に対して、

$$ A^nB=(T_{n+2},T_{n+1},T_n)^\mathsf T. $$

最初のいくつかの状態ベクトルをすべて表示することで、漸化式と添字の 取り方を同時に確認できます。

tribonacciState 1
$\begin{pmatrix} 1 \\ 1 \\ 0\\ \end{pmatrix}$
tribonacciState 3
$\begin{pmatrix} 4 \\ 2 \\ 1\\ \end{pmatrix}$
tribonacciState 5
$\begin{pmatrix} 13 \\ 7 \\ 4\\ \end{pmatrix}$

大きな添字へ一気に進む

行列の累乗には二乗を繰り返す方法を使うため、$A^{100}B$ の計算に 漸化式を明示的に100回適用する必要はありません。

tribonacciState 100
$\begin{pmatrix} 180396380815100901214157639 \\ 98079530178586034536500564 \\ 53324762928098149064722658\\ \end{pmatrix}$

まとめ

パターンマッチによりコンパニオン行列の構造を簡潔に定義し、テンソル縮約に より行列とベクトルの積を計算できます。ここで採用した添字の規約では、 最終ベクトルの第1成分が $T_{102}$ です。

リンク

Egison 数学ノート目次に戻る