{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "e9432603",
   "metadata": {},
   "source": [
    "# Tribonacci Numbers by Matrix Exponentiation\n",
    "\n",
    "The Tribonacci sequence used here is\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",
    "A companion matrix turns this recurrence into repeated linear\n",
    "algebra, which makes fast exponentiation available.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bb62b1af",
   "metadata": {},
   "source": [
    "## Build the transition matrix by pattern matching\n",
    "\n",
    "The first row sums the three current state entries. The\n",
    "subdiagonal shifts the older values down one position:\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",
    "The generator describes those two structural patterns rather than\n",
    "listing nine unrelated entries.\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": [
    "## Initial state\n",
    "\n",
    "The vector $B=(1,0,0)^\\mathsf T$ stores\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": [
    "## Advance the recurrence\n",
    "\n",
    "For every $n\\ge0$,\n",
    "\n",
    "$$\n",
    "A^nB=(T_{n+2},T_{n+1},T_n)^\\mathsf T.\n",
    "$$\n",
    "\n",
    "The first few full states make both the recurrence and the indexing\n",
    "convention visible.\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": [
    "## Jump directly to a large index\n",
    "\n",
    "Matrix exponentiation uses repeated squaring, so computing\n",
    "$A^{100}B$ does not require one hundred explicit recurrence steps.\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": [
    "## Takeaway\n",
    "\n",
    "Pattern matching gives a compact structural definition of the\n",
    "companion matrix, while tensor contraction performs the\n",
    "matrix-vector product. The first component of the final vector is\n",
    "$T_{102}$ under the stated indexing convention.\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
}