{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "e9432603",
   "metadata": {},
   "source": [
    "# 行列の累乗によるトリボナッチ数の計算\n",
    "\n",
    "ここでは、次のトリボナッチ数列を扱います。\n",
    "\n",
    "$$\n",
    "T_0=0,\\quad T_1=0,\\quad T_2=1,\\qquad\n",
    "T_{n+3}=T_{n+2}+T_{n+1}+T_n.\n",
    "$$\n",
    "\n",
    "コンパニオン行列を使うと、この漸化式を線形代数の演算の反復として\n",
    "表せるため、高速な累乗計算を利用できます。\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bb62b1af",
   "metadata": {},
   "source": [
    "## パターンマッチで遷移行列を組み立てる\n",
    "\n",
    "第1行は現在の状態に含まれる三つの成分を足し合わせます。下副対角線は、\n",
    "過去の値を一つ下の位置へ移します。\n",
    "\n",
    "$$\n",
    "A=\n",
    "\\begin{pmatrix}\n",
    "1&1&1\\\\\n",
    "1&0&0\\\\\n",
    "0&1&0\n",
    "\\end{pmatrix}.\n",
    "$$\n",
    "\n",
    "ジェネレーターでは、互いに無関係な九つの成分を列挙する代わりに、\n",
    "この二つの構造的なパターンを記述します。\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "f4e0074f",
   "metadata": {},
   "outputs": [],
   "source": [
    "def m : Integer := 3\n",
    "\n",
    "def A : Matrix Integer :=\n",
    "  generateTensor\n",
    "    (\\match as list integer with\n",
    "      | [#1, _] -> 1\n",
    "      | [$x, #(x - 1)] -> 1\n",
    "      | _ -> 0)\n",
    "    [m, m]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "17391d59",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "$\\begin{pmatrix} 1 & 1 & 1 \\\\ 1 & 0 & 0 \\\\ 0 & 1 & 0 \\\\ \\end{pmatrix}$"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "A\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0f9ae517",
   "metadata": {},
   "source": [
    "## 初期状態\n",
    "\n",
    "ベクトル $B=(1,0,0)^\\mathsf T$ には\n",
    "$(T_2,T_1,T_0)^\\mathsf T$ が格納されています。\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "f226502c",
   "metadata": {},
   "outputs": [],
   "source": [
    "def B : Vector Integer :=\n",
    "  generateTensor\n",
    "    (\\[x] -> if x = 1 then 1 else 0)\n",
    "    [m]\n",
    "\n",
    "def tribonacciState (n : Integer) : Vector Integer :=\n",
    "  MV.* (M.power A n) B\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "8845e48f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "$\\begin{pmatrix} 1 \\\\ 0 \\\\ 0\\\\ \\end{pmatrix}$"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "B\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "38f64217",
   "metadata": {},
   "source": [
    "## 漸化式を進める\n",
    "\n",
    "任意の $n\\ge0$ に対して、\n",
    "\n",
    "$$\n",
    "A^nB=(T_{n+2},T_{n+1},T_n)^\\mathsf T.\n",
    "$$\n",
    "\n",
    "最初のいくつかの状態ベクトルをすべて表示することで、漸化式と添字の\n",
    "取り方を同時に確認できます。\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "24cfc0dc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "$\\begin{pmatrix} 1 \\\\ 1 \\\\ 0\\\\ \\end{pmatrix}$"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "tribonacciState 1\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "b47d3f9f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "$\\begin{pmatrix} 4 \\\\ 2 \\\\ 1\\\\ \\end{pmatrix}$"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "tribonacciState 3\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "66b99ab3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "$\\begin{pmatrix} 13 \\\\ 7 \\\\ 4\\\\ \\end{pmatrix}$"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "tribonacciState 5\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "771937dd",
   "metadata": {},
   "source": [
    "## 大きな添字へ一気に進む\n",
    "\n",
    "行列の累乗には二乗を繰り返す方法を使うため、$A^{100}B$ の計算に\n",
    "漸化式を明示的に100回適用する必要はありません。\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "f8ef5481",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "$\\begin{pmatrix} 180396380815100901214157639 \\\\ 98079530178586034536500564 \\\\ 53324762928098149064722658\\\\ \\end{pmatrix}$"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "tribonacciState 100\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a23f25f0",
   "metadata": {},
   "source": [
    "## まとめ\n",
    "\n",
    "パターンマッチによりコンパニオン行列の構造を簡潔に定義し、テンソル縮約に\n",
    "より行列とベクトルの積を計算できます。ここで採用した添字の規約では、\n",
    "最終ベクトルの第1成分が $T_{102}$ です。\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Egison",
   "language": "egison",
   "name": "egison"
  },
  "language_info": {
   "codemirror_mode": "egison",
   "file_extension": ".egi",
   "mimetype": "text/x-egison",
   "name": "egison"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}