wiedemann_tower.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "97bc9dd4-1dbe-44e2-b2de-8e5ce2312c52",
   "metadata": {
    "editable": true,
    "slideshow": {
     "slide_type": ""
    },
    "tags": []
   },
   "source": [
    "# Wiedemann Tower: Recursive Construction\n",
    "\n",
    "## References\n",
    "\n",
    "- [1] Wiedemann, Doug. \"An iterated quadratic extension of GF (2).\" Fibonacci Quart 26.4 (1988): 290-295."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "b112f166-44f1-44c3-b498-0c2eea0394d3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Univariate Polynomial Ring in X0 over Finite Field of size 2 (using GF2X)"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "F2 = GF(2)\n",
    "F2\n",
    "\n",
    "# Define the polynomial ring F2[X]\n",
    "\n",
    "F2_PolyR.<X0> = F2[]; F2_PolyR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "681f4468-4337-408a-bb28-2a719d242510",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Finite Field in a of size 2^2"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Define F4: the degree-2 extension over F2\n",
    "\n",
    "F4.<a> = F2.extension(X0^2 + X0 + 1, 'a'); F4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "c73ee646-1108-42e1-bfd6-0e8484557cfc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Univariate Polynomial Ring in X1 over Finite Field in a of size 2^2"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Define the polynomial ring F4[X]\n",
    "\n",
    "F4_PolyR.<X1> = F4[]; F4_PolyR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "cbdccef0-a49d-483e-b5c7-e62d56089ac5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Univariate Quotient Polynomial Ring in b over Finite Field in a of size 2^2 with modulus b^2 + a*b + 1"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Define F16: the degree-2 extension over F4\n",
    "\n",
    "F16.<b> = F4.extension(X1^2 + a * X1 + 1, 'b'); F16"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "c66007b6-17b2-4d9d-a5ca-8e93a968fcc4",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define F16[X]\n",
    "\n",
    "F16_PolyR.<X2> = F16[]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "900475ec-8718-427d-982c-17abe48d2e9c",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Univariate Quotient Polynomial Ring in c over Univariate Quotient Polynomial Ring in b over Finite Field in a of size 2^2 with modulus b^2 + a*b + 1 with modulus c^2 + b*c + 1"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Define F256: the degree-2 extension over F16\n",
    "\n",
    "F256.<c> = F16.extension(X2^2 + b * X2 + 1, 'c'); F256"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "c0dd1b1b-e9fa-4775-b3f2-40c8a6c6395b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Univariate Polynomial Ring in X3 over Univariate Quotient Polynomial Ring in c over Univariate Quotient Polynomial Ring in b over Finite Field in a of size 2^2 with modulus b^2 + a*b + 1 with modulus c^2 + b*c + 1"
      ]
     },
     "execution_count": 8,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Define F256[X]\n",
    "\n",
    "F256_PolyR.<X3> = F256[]; F256_PolyR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "967393b4-1cc6-471b-8cc5-07298f9b44f6",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Univariate Quotient Polynomial Ring in d over Univariate Quotient Polynomial Ring in c over Univariate Quotient Polynomial Ring in b over Finite Field in a of size 2^2 with modulus b^2 + a*b + 1 with modulus c^2 + b*c + 1 with modulus d^2 + c*d + 1"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Define F65536: the degree-2 extension over F256\n",
    "\n",
    "F65536.<d> = F256.extension(X3^2 + c * X3 + 1, 'd'); F65536"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "2ad0c56e-0c8d-4dac-83c3-33ac86657d1e",
   "metadata": {},
   "outputs": [],
   "source": [
    "# def intFromF4(t):\n",
    "#     if t == 0:\n",
    "#         return '00'\n",
    "#     elif t == 1:\n",
    "#         return '01'\n",
    "#     elif t == a:\n",
    "#         return '10'\n",
    "#     else: \n",
    "#         return '11'\n",
    "\n",
    "# def intFromF16(t):\n",
    "#     t0=t.list()[0]\n",
    "#     t1=t.list()[1]\n",
    "#     return intFromF4(t1)+\":\"+intFromF4(t0)\n",
    "\n",
    "def log(s, t):\n",
    "    for i in range(1,16):\n",
    "        if t^i == s:\n",
    "            return i\n",
    "\n",
    "def conjugates(t):\n",
    "    r = []\n",
    "    for i in [1,2,4,8,16,32,64,128,256]:\n",
    "        t2 = t^i\n",
    "        if t2 in r:\n",
    "            continue\n",
    "        r.append(t^i)\n",
    "    return r\n",
    "\n",
    "def ord(t):\n",
    "    card = t.parent().order()\n",
    "    \n",
    "    for i in range(1,card):\n",
    "        if t**i == 1:\n",
    "            return i\n",
    "    return \"False\"\n",
    "    \n",
    "def logtable(prim_elem):\n",
    "    for i in range(1, 16):  # Iterate from 0 to p-1\n",
    "        elem = prim_elem**i\n",
    "        print(\"{:3},ord={:3} t={}\".format(i, ord(elem), elem))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "61246ab3-1a86-4eae-b16a-d21348cbd6c7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  1,ord= 15 t=b + 1\n",
      "  2,ord= 15 t=a*b\n",
      "  3,ord=  5 t=b + a\n",
      "  4,ord= 15 t=b + a + 1\n",
      "  5,ord=  3 t=a\n",
      "  6,ord=  5 t=a*b + a\n",
      "  7,ord= 15 t=(a + 1)*b\n",
      "  8,ord= 15 t=a*b + a + 1\n",
      "  9,ord=  5 t=a*b + 1\n",
      " 10,ord=  3 t=a + 1\n",
      " 11,ord= 15 t=(a + 1)*b + a + 1\n",
      " 12,ord=  5 t=b\n",
      " 13,ord= 15 t=(a + 1)*b + 1\n",
      " 14,ord= 15 t=(a + 1)*b + a\n",
      " 15,ord=  1 t=1\n"
     ]
    }
   ],
   "source": [
    "logtable(b+1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "26e47ab3-0ff0-4d78-b19e-ff5dfa97c901",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Vitalik's Puzzle\n",
    "\n",
    "# 2 + 2 = 0\n",
    "\n",
    "(a + 1) + (a + 1) == 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "a7d08f5d-258c-4999-b704-cc941d9ee65d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Vitalik's Puzzle\n",
    "# Basis is [1, a, b, a * b]\n",
    "\n",
    "# 2 * 6 ?= 11\n",
    "#\n",
    "#\n",
    "# 2 = (0010), |--> < (0, 0, 1, 0), (`a*b`, `b`, `a`, `1`) > = (a)\n",
    "# 6 = (0110), |--> < (0, 1, 1, 0), (`a*b`, `b`, `a`, `1`) > = (a + b)\n",
    "# 11= (1011), |--> < (1, 0, 1, 1), (`a*b`, `b`, `a`, `1`) > = a*b + a + 1\n",
    "#\n",
    "\n",
    "(a) * (a + b) == (a*b + a + 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "599d5d54-9627-49f7-ba6f-cbcc6301812f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Vitalik's Puzzle\n",
    "#\n",
    "# Basis is [1, a, b, a * b]\n",
    "\n",
    "# 3 / 5 ?= 9  ===> 5 * 9 ?= 3\n",
    "#\n",
    "#\n",
    "# 5 = (0101), |--> < (0, 1, 0, 1), (`a*b`, `b`, `a`, `1`) > = (b + 1)\n",
    "# 9 = (1001), |--> < (1, 0, 0, 1), (`a*b`, `b`, `a`, `1`) > = (a*b + 1)\n",
    "# 3 = (0011), |--> < (0, 0, 1, 1), (`a*b`, `b`, `a`, `1`) > = (a + 1)\n",
    "#\n",
    "\n",
    "(a*b + 1) * (b + 1) == (a + 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "72af4830-ee9c-41a9-ad5f-78ed14e43ccc",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Vitalik's Puzzle\n",
    "#\n",
    "# Basis is [1, a, b, a * b]\n",
    "\n",
    "# 3^2 + 8^2 ?= 11^2\n",
    "#\n",
    "# 3 = (0011), |--> < (0, 0, 1, 1), (`a*b`, `b`, `a`, `1`) > = (a + 1)\n",
    "# 8 = (1000), |--> < (1, 0, 0, 0), (`a*b`, `b`, `a`, `1`) > = (a*b)\n",
    "# 11= (1011), |--> < (1, 0, 1, 1), (`a*b`, `b`, `a`, `1`) > = a*b + a + 1\n",
    "\n",
    "\n",
    "(a + 1)^2 + (a * b)^2 == (a*b + a + 1)^2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "55da4a0a-1ed9-4fae-a901-4392b89343f3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(True, True, True, True)"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Fun Fact: The relations between (a, b, c, d)\n",
    "#\n",
    "#   x_{i+1} + x_{i+1}^{-1} == x_i\n",
    "#\n",
    "\n",
    "a + 1/a == 1, b + 1/b == a, c + 1/c == b, d + 1/d == c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "b953605e-556d-4195-8714-7aad4b6a1209",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(3, 5, 17, 257)"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Check the orders of the roots\n",
    "#\n",
    "\n",
    "ord(a), ord(b), ord(c), ord(d)\n",
    "\n",
    "# Output:\n",
    "#   (3, 5, 17, 257)\n",
    "# \n",
    "# Clearly, `b` is in order 5, not a primitive element of F16\n",
    "# Neither `c` nor `d` is primitive in F256 or F65536"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "49882517-cdd0-42e5-8f07-364a3eea8842",
   "metadata": {},
   "outputs": [],
   "source": [
    "# `3, 5, 17, 257` are so called Fermat Numbers FN0, FN1, FN2, FN3, \n",
    "#\n",
    "#        FNk = 2^{2^k} + 1\n",
    "#\n",
    "# Fact: FN0, FN1, FN2, FN3, FN4 are prime \n",
    "#   However, FN5 is composite, FN5 = 2^{32} + 1 = 4,294,967,297 = 641 * 6,700,417\n",
    "#   (see more: https://en.wikipedia.org/wiki/Fermat_number)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "ea40d29b-ea3c-4bba-b668-88bb478c2172",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(True, True, True)"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Unproven Fact: a primitive element of F_{2^{2^{n+1}}} is `x0 * x1 * ... * xn`\n",
    "#\n",
    "#  F4   : x0              = `a`\n",
    "#  F16  : x1 * x0         = `b * a`\n",
    "#  F256 : x2 * x1 * x0    = `c * b * a`\n",
    "\n",
    "ord(a) == (2^(2^1) - 1), ord(b*a) == (2^4 - 1), ord(c*b*a) == (2^8 - 1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "060c4b0e-0117-4121-8bad-c2067a59efc5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# How to compute the minimal polynomial of `xi` over F2? \n",
    "#\n",
    "#    f0(a) = 0,  f0 = X0^2 + X0 + 1 is the minimal polynomial of `a` over F2\n",
    "\n",
    "f0 = X0^2 + X0 + 1 \n",
    "f0.subs(X0=a)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "1d5a1d08-1806-4051-9e04-cee3d718d7de",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#   f1(b) = 0  ==> f1 = X1^2 + a * X1 + 1 is the minimal polynomial of `b` over F4\n",
    "#   f0(a) = f0(b + 1/b) = 0\n",
    "#\n",
    "# Therefore, \n",
    "#\n",
    "#     b^2 + b^(-2) + b + b^(-1) + 1 = 0\n",
    "# \n",
    "# Multiply (b^2) to the both sides:\n",
    "#\n",
    "#     b^4 + 1 + b^3 + b + b^2 = 0\n",
    "#\n",
    "f1 = X0^4 + X0^3 + X0^2 + X0 + 1; f1.subs(X0 = b) == 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "id": "70473924-95bc-42f8-90fd-9b85b9533ed2",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[b, a*b + 1, b + a, a*b + a]"
      ]
     },
     "execution_count": 22,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# The alternative way: compute the minimal polynomial of `b` from conjugates\n",
    "\n",
    "conjugates(b)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "id": "041928bb-1818-41f5-ad2a-26e890392fb7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "X^4 + X^3 + X^2 + X + 1"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# The result will be (X - b) * (X - (a*b +1)) * (X - (b+a)) * (X - (a*b + a))\n",
    "# All terms in F16\\F2 will be canceled out, and the coefficients of the polynomial are all in F2\n",
    "#\n",
    "Mb = ((X2 - b) * (X2 - a*b - 1) * (X2 - b - a) * (X2 - a*b -a)).change_variable_name(var('X')); Mb\n",
    "\n",
    "# Output: \n",
    "#   X^4 + X^3 + X^2 + X + 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "id": "c422c564-0e85-4df4-8a5d-430f6cb687c2",
   "metadata": {},
   "outputs": [],
   "source": [
    "# An iterated alg to compute the minimal polynomial of 'xi' over F2 \n",
    "#\n",
    "#  A0 = 1 + X^2,  B0 = X\n",
    "#  A{k+1} = Ak^2 + Bk^2,  B{k+1} = An*Bn\n",
    "#  \n",
    "#  minimal_poly(xn) = An + Bn\n",
    "#\n",
    "\n",
    "def minimal_poly_iter(i):\n",
    "    A = 1 + X0^2\n",
    "    B = X0\n",
    "    for i in range(0,i-1):\n",
    "        Ak = A^2 + B^2\n",
    "        Bk = A * B\n",
    "        A  = Ak\n",
    "        B  = Bk\n",
    "    return (A + B).change_variable_name(var('X'))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "id": "2c750290-ca2e-4657-b47b-344f0b8e5625",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 25,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Compute minimal_poly(b)\n",
    "\n",
    "minimal_poly_iter(2) == Mb"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "id": "48a3ab3b-416d-46b7-9e92-cea1822a9981",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "X^8 + X^7 + X^6 + X^4 + X^2 + X + 1"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Compute minimal_poly(c)\n",
    "\n",
    "Mc = minimal_poly_iter(3); Mc"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "id": "a042241e-6e2b-455a-928f-5400c9db71a5",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Mc.subs(X=c) == 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "id": "775a29ad-6b19-465d-b0a1-930bce42f93a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[c,\n",
       " b*c + 1,\n",
       " (a*b + a)*c + a*b,\n",
       " (a*b + 1)*c + a + 1,\n",
       " c + b,\n",
       " b*c + a*b,\n",
       " (a*b + a)*c + (a + 1)*b + a,\n",
       " (a*b + 1)*c + a*b + 1]"
      ]
     },
     "execution_count": 28,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Double check Mc: the minimal polynomial of `c` over F2\n",
    "\n",
    "conjugates(c)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "id": "bb2415c1-b6b9-4e41-bbc8-e3fc8fb03ffd",
   "metadata": {},
   "outputs": [],
   "source": [
    "Mc_alt = F256.one()\n",
    "for t in conjugates(c): Mc_alt = Mc_alt * (X3 + t)\n",
    "Mc_alt = Mc_alt.change_variable_name(var('X'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "id": "e048122b-aa6b-4790-9d1f-05745b0085f4",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Mc_alt == Mc"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0f4cd344-5a8a-4c0c-8181-50ee85c0c597",
   "metadata": {},
   "source": [
    "## TODO: trace and basis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "045f6f94-f982-4961-b283-c1fc15424ac8",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 10.3",
   "language": "sage",
   "name": "sagemath-10.3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}