HW3_EE_4940.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 68,
   "id": "99953193-e838-4ac6-b361-f6fa997db102",
   "metadata": {},
   "outputs": [],
   "source": [
    "#Import libraries\n",
    "import sympy as sp\n",
    "import numpy as np\n",
    "np.set_printoptions(legacy='1.25')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f6d7be78",
   "metadata": {},
   "source": [
    "## Problem 1"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ae303dc9",
   "metadata": {},
   "source": [
    "Consider the matrix\n",
    "\n",
    "$ A = \\begin{bmatrix}\n",
    "3 & 1 & 6 \\\\\n",
    "3 & 0 & 4 \\\\\n",
    "1 & 3 & 2 \n",
    "\\end{bmatrix}  $\n",
    "\n",
    "in GF(7).\n",
    "\n",
    "(a) Find the determinant of $A$. Is it invertible?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "b7fa1405",
   "metadata": {},
   "source": [
    ">A matrix is invertible if it is a square $n$ by $n$ matrix, and the determinant is non-zero."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "id": "ec3ecba6",
   "metadata": {},
   "outputs": [],
   "source": [
    "def find_modulo_inverse(input, modulo) : #Finds inverse in modulo.\n",
    "    #e.g. 2 inverse in mod 3 is 2 b/c 2*2 = 4(mod 3) = 1\n",
    "    #==Inputs== \n",
    "    #input: value to find the inverse of\n",
    "    #modulo: value of the modulus\n",
    "    i=modulo\n",
    "    while i>0 :\n",
    "        if np.mod(input*i,modulo) == 1 :\n",
    "            return i\n",
    "        i-=1\n",
    "    return 0\n",
    "\n",
    "def recursive_det(matrix, p): # Recursively finds the determinant of a matrix in GF(p)\n",
    "    #base case: 2x2 matrix\n",
    "    #output: determinant of a matrix.\n",
    "\n",
    "    if len(matrix) == 2 : #base case\n",
    "        return round(np.linalg.det(matrix))\n",
    "    else :\n",
    "        i=0\n",
    "        j=0\n",
    "        while j<len(matrix) :\n",
    "            coefficient = (-1)**(i+j)\n",
    "            submatrix = np.delete((np.delete(matrix, i,axis=0)),j,axis=1)\n",
    "            value = recursive_det(submatrix,p) #recursion\n",
    "            value2 = coefficient*value*matrix[i][j]\n",
    "            if i==0 :\n",
    "                if j==0 : accum = value2 % p\n",
    "                else : accum = (accum + value2) % p\n",
    "            j+=1\n",
    "\n",
    "    return accum"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "id": "025d8e7c",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The determinant of A = 2\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[3,1,6],[3,0,4],[1,3,2]])\n",
    "\n",
    "print(f'The determinant of A = {recursive_det(A,7)}')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7cf081ee",
   "metadata": {},
   "source": [
    "> The determinant of A = 2, which is not zero, therefore A is invertible.\n",
    "\n",
    "> I used numpy's numpy.linalg.det() function which returns the determinant of a matrix. I take the modulus of the result."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7074ab8e",
   "metadata": {},
   "source": [
    "(b) Recall that in real numbers, the inverse of $A$ can be found as\n",
    "\n",
    "$ A^{-1} = \\frac{1}{det(A)} \\begin{bmatrix}\n",
    "det(A(1, 1)) & − det(A(1, 2)) & det(A(1, 3)) \\\\\n",
    "− det(A(2, 1)) & det(A(2, 2)) & − det(A(2, 3)) \\\\\n",
    "det(A(3, 1)) & − det(A(3, 2)) & det(A(3, 3)) \n",
    "\\end{bmatrix}  $\n",
    "\n",
    "where $A(i, j)$ is the matrix obtained by removing row $i$ and column $j$ of matrix $A$. For instance,\n",
    "\n",
    "$ A(1,2) = \\begin{bmatrix}\n",
    "3 & 4 \\\\\n",
    "1 & 2 \n",
    "\\end{bmatrix}  $\n",
    "\n",
    "Note that $±$ det($A(i, j)$) is still an integer, and hence, you can find the inverse of $A$ in GF(7) by replacing division by det($A$) by multiplying by in the multiplicative inverse of det($A$) in GF(7). Now, use the method above to find the inverse of $A$ in GF(7). You may use MATLAB or any other software. Describe your method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "id": "71d41856",
   "metadata": {},
   "outputs": [],
   "source": [
    "def inverse_rec(matrix,p) :\n",
    "    i=0\n",
    "    output = np.zeros([len(matrix),len(matrix)],dtype=int)\n",
    "    det = recursive_det(matrix,p)\n",
    "    inv_det = find_modulo_inverse(det, p)\n",
    "    while i<len(matrix):\n",
    "        j=0\n",
    "        while j<len(matrix) :\n",
    "            submatrix = np.delete((np.delete(matrix, i,axis=0)),j,axis=1)\n",
    "            #print(output)\n",
    "            #print(submatrix)\n",
    "            val = recursive_det(submatrix,p)*((-1)**(i+j))\n",
    "            #print(f'val={val}')\n",
    "            output[i][j] = inv_det*val % p\n",
    "            j+=1\n",
    "        i+=1\n",
    "    output = np.array(output)\n",
    "    return output.transpose()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "id": "dd5adb63",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "A^-1 = \n",
      "[[1 1 2]\n",
      " [6 0 3]\n",
      " [1 3 2]]\n",
      "\n",
      "A*A^-1 = \n",
      "[[1 0 0]\n",
      " [0 1 0]\n",
      " [0 0 1]]\n",
      "2\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[3,1,6],[3,0,4],[1,3,2]])\n",
    "A_inv = inverse_rec(A,7)\n",
    "print('A^-1 = ')\n",
    "print(A_inv)\n",
    "\n",
    "print('')\n",
    "#Testing to see if inverse produced identity matrix.\n",
    "print('A*A^-1 = ')\n",
    "print(np.matmul(A,A_inv) % 7)\n",
    "\n",
    "print(recursive_det(A,7))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "01103641",
   "metadata": {},
   "source": [
    "> $ A^{-1} = \\begin{bmatrix}\n",
    "1 & 1 & 2 \\\\\n",
    "6 & 0 & 3 \\\\\n",
    "1 & 3 & 2 \n",
    "\\end{bmatrix}  $\n",
    "\n",
    "> $A \\cdot A^{-1} =\n",
    "\\begin{bmatrix}\n",
    "1 & 0 & 0 \\\\\n",
    "0 & 1 & 0 \\\\\n",
    "0 & 0 & 1 \n",
    "\\end{bmatrix}$\n",
    "\n",
    "Using a RREF calculator, https://www.emathhelp.net/en/calculators/linear-algebra/reduced-row-echelon-form-rref-calculator/?i=%5B%5B3%2C2%2C5%5D%2C%5B3%2C5%2C6%5D%2C%5B2%2C3%2C1%5D%5D&reduced=on ,\n",
    "\n",
    "we get:\n",
    "\n",
    "$A \\cdot A^{-1} = \\begin{bmatrix}\n",
    "1 & 0 & 0 \\\\\n",
    "0 & 1 & 0 \\\\\n",
    "0 & 0 & 1 \n",
    "\\end{bmatrix}$\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0cc7cee6",
   "metadata": {},
   "source": [
    "## Problem 2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "22ab9ef3",
   "metadata": {},
   "source": [
    "Consider a general block cipher algorithm that maps $n$ bits to $n$ bits.\n",
    "\n",
    "(a) How many different such mapping (including the identity map) exist?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7b5dfe95",
   "metadata": {},
   "source": [
    "> $2^{N}!$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f8c91efb",
   "metadata": {},
   "source": [
    "(b) Let $K(n)$ be your answer to part (a). Bob and Alice want to communicate privately. In order to set the mapping, they shared a long table including all possible $K(n)$ mappings and associated a number between $0$ and $K(n) − 1$ to each. Then, they wrote the corresponding number for each mapping in Binary representation. For instance, if $K(n) = 1000$, mapping number $347$ is represented by $’0101011011’$. Then, each time Bob wants to send a private message to Alice, he chooses one of these mappings and sends her its binary representation through a secure channel. How many bits does Bob have to send for each mapping?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "fc02a722",
   "metadata": {},
   "source": [
    "> ceiling($log_{2}(2^{N}!-1)$). Because we have $2^{N}!$ mappings, and number x can be represented by ceiling($log_2(x)$) bits"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "bd7c6c7d",
   "metadata": {},
   "source": [
    "\n",
    "(c) Recall that we mentioned in class that $n$-to-$n$ bits mappings can be also represented by $n$ × $2n$ bits (similar to Table. 4.1 of the textbook). Is this number less or greater than what you found in part (b)? Explain the discrepancy."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "302f01a5",
   "metadata": {},
   "source": [
    "> This number would be greater, because for each possible mapping, there are twice as many choices to choose from."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "94f424c1",
   "metadata": {},
   "source": [
    "## Problem 3"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a327fc4b",
   "metadata": {},
   "source": [
    "Consider the GF($2^{8}$) generated by $q(x) = x^{8} + x^{7} + x^{2} + x + 1$. Evaluate the following quantities:\n",
    "\n",
    "(a) $\\{BC\\} + \\{2E\\}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "id": "d1343ec6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0x92\n"
     ]
    }
   ],
   "source": [
    "def hex_add(A,B,p) : #Adds A to B in GF(p), where A and B are hexadecimal values.\n",
    "    A_array = [x for x in bin(A)]\n",
    "    B_array = [x for x in bin(B)]\n",
    "    A_array.pop(0) #get rid of '0b' from the arrays\n",
    "    A_array.pop(0)\n",
    "    B_array.pop(0)\n",
    "    B_array.pop(0)\n",
    "    A_array = [int(x) for x in A_array] #convert each element to integer\n",
    "    B_array = [int(x) for x in B_array]\n",
    "\n",
    "\n",
    "    while len(A_array) < 8 : #add zeros to make array an 8 bit number\n",
    "        A_array.insert(0,0)\n",
    "    while len(B_array) < 8 :\n",
    "        B_array.insert(0,0)\n",
    "\n",
    "    output = [0] * 8\n",
    "    i=0\n",
    "    while i<8 :\n",
    "        output[i] = (A_array[i] + B_array[i])%p #add each bit, modulus 2\n",
    "        i+=1\n",
    "    return hex(int('0b' + ''.join(str(x) for x in output), 2))\n",
    "\n",
    "print(hex_add(0xBC,0x2E,2))\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e9a2c932",
   "metadata": {},
   "source": [
    "> $\\{BC\\} + \\{2E\\} = \\{92\\}$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1f720f23",
   "metadata": {},
   "source": [
    "(b) $\\{0A\\} \\cdot \\{A3\\}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "id": "bf9f8fed",
   "metadata": {},
   "outputs": [],
   "source": [
    "#These are mostly altered versions of my code from HW2.\n",
    "\n",
    "def find_degree(polynomial) : #Finds the degree of a polynomial vector by ignoring finding the first non-zero value.\n",
    "    shift = 0\n",
    "    for i in polynomial :\n",
    "        if i!=0 :\n",
    "            return len(polynomial)-shift-1\n",
    "        else : \n",
    "            shift +=1\n",
    "    return 0\n",
    "\n",
    "def remove_zeros(polynomial) : #Removes all zeros in a polynomial vector preceeding the first non-zero element.\n",
    "    #e.g. [0, 0, 0, 2, 0, 3, 0] --> [2, 0, 3, 0]\n",
    "\n",
    "    shift = 0\n",
    "    for i in polynomial :\n",
    "        if i!=0 :\n",
    "            while shift>0 :\n",
    "                polynomial.pop(0)\n",
    "                shift-=1\n",
    "            break\n",
    "        else : \n",
    "            shift +=1\n",
    "    return polynomial\n",
    "\n",
    "def MyGFAdd(A, B, M, p) : # Add A and B modulo p.\n",
    "    #==INPUTS==\n",
    "    #prime: prime number from p in GF(p)\n",
    "    #A: an array of polynomial coefficients e.g. 2x^2 + 1 --> [2 0 1]\n",
    "    #B: same formate as A\n",
    "    #M: do not care\n",
    "    #p: p in GF(p), prime number\n",
    "    #==OUTPUT==\n",
    "    #Returns polynomial coefficient array that is \n",
    "\n",
    "\n",
    "    degree_A = find_degree(A) #highest degree of A\n",
    "    degree_B = find_degree(B) #highest degree of B\n",
    "    degree_diff = abs(degree_A - degree_B)\n",
    "    A = remove_zeros(A) #remove preceeding zeros\n",
    "    B = remove_zeros(B)\n",
    "\n",
    "    #Add A and B differently depending on which one is longer or shorter.\n",
    "    if degree_A > degree_B : \n",
    "        #print('A > B')\n",
    "        output = [0] * len(A) # A is a longer array\n",
    "        i=0\n",
    "        while i<len(A) :\n",
    "            if i < degree_diff : output[i] = A[i] #don't start adding B until A and B are at the same degree.\n",
    "            else : output[i] = A[i]+B[i-degree_diff]\n",
    "            i+=1\n",
    "\n",
    "    else :\n",
    "        if degree_A < degree_B : \n",
    "            #print('A < B')\n",
    "            output = [0] * len(B) # B is a logner array\n",
    "            i=0\n",
    "            while i<len(B) :\n",
    "                if i < degree_diff : output[i] = B[i] #don't start adding B until A and B are at the same degree.\n",
    "                else : output[i] = B[i]+A[i-degree_diff]\n",
    "                i+=1\n",
    "        else : #degree_A = degree_B\n",
    "            #print('A=B')\n",
    "            output = [0] * len(A)\n",
    "            i=0\n",
    "            for element in A :\n",
    "                output[i] = element+B[i]\n",
    "                i+=1\n",
    "    \n",
    "    # apply modulus to output.\n",
    "    i=0\n",
    "    while i<len(output) :\n",
    "        output[i] = np.mod(output[i], p)\n",
    "        i+=1\n",
    "\n",
    "    return output\n",
    "\n",
    "def single_mult(k,degree_k,A) : #Multiplies one value into a full polynomial\n",
    "    #==hexS==\n",
    "    #k: scalar coefficent of single value.\n",
    "    #degree_k: degree of \n",
    "    #A: polynomial coefficient vector.\n",
    "    #==OUTPUT==\n",
    "    #Returns a polynomial coefficient array: A*kx^(degree)\n",
    "    degree_A = find_degree(A)\n",
    "    A = remove_zeros(A)\n",
    "    output = [0] * (degree_A+degree_k+1)\n",
    "\n",
    "    if k == 0 :\n",
    "        return output\n",
    "\n",
    "    i=0\n",
    "    while i<len(A):\n",
    "        output[i] = k*A[i]\n",
    "        i+=1\n",
    "\n",
    "    return output\n",
    "\n",
    "def raw_poly_mult_GF(A, B, GF) :\n",
    "    #Multiplies two polynomials, ignoring degree > M\n",
    "    degree_A = find_degree(A)\n",
    "    degree_B = find_degree(B)\n",
    "    \n",
    "    if degree_A > degree_B : \n",
    "        #print('A > B')\n",
    "        partial_output = [[0 for _ in range(degree_A+degree_B+1)] for _ in range(len(B))]\n",
    "        i=0\n",
    "        for element in B :\n",
    "                partial_output[i] = single_mult(element,(degree_B-i),A)\n",
    "                i+=1\n",
    "\n",
    "    else :\n",
    "        if degree_A < degree_B : \n",
    "            #print('A < B')\n",
    "            partial_output = [[0 for _ in range(degree_A+degree_B+1)] for _ in range(len(A))]\n",
    "            i=0\n",
    "            for element in A :\n",
    "                partial_output[i] = single_mult(element,(degree_A-i),B)\n",
    "                i+=1\n",
    "\n",
    "        else : #degree_A = degree_B\n",
    "            #print('A=B')\n",
    "            partial_output = [[0 for _ in range(degree_A+degree_B+1)] for _ in range(len(A))]\n",
    "            i=0\n",
    "            for element in A :\n",
    "                partial_output[i] = single_mult(element,(degree_A-i),B)\n",
    "                i+=1\n",
    "\n",
    "\n",
    "    max_degree = 0\n",
    "    for element in partial_output :\n",
    "        deg = find_degree(element)\n",
    "        if deg > max_degree :\n",
    "            max_degree = deg\n",
    "    #print(f'max_deg={max_degree}')\n",
    "\n",
    "    output = [0] * (max_degree+1)\n",
    "\n",
    "    output = partial_output[0]\n",
    "    #print(f'partial_output={partial_output}')\n",
    "    i=1\n",
    "    while i<(len(partial_output)) :\n",
    "        output = MyGFAdd(output,partial_output[i],0,GF)\n",
    "        #print(f'output={output}')\n",
    "        i+=1\n",
    "\n",
    "\n",
    "    return output\n",
    "\n",
    "def find_highest_coefficient(M,p) :\n",
    "    #Finds what the highest degree coefficient is equal to. Assumes that that coefficient = 1.\n",
    "    #e.g. \n",
    "    # hex: [1, 1, 1] output: [2, 2]\n",
    "    #      (x^2 + x + 1)   =   (2x+2)\n",
    "\n",
    "    M[0] = 0\n",
    "    M = list(np.multiply(M,(-1)))\n",
    "    M = np.mod(remove_zeros(M),p)\n",
    "    return M\n",
    "\n",
    "def MyGFMult(A, B, M, p) :\n",
    "    #Returns the multiplication of A and B in GF(p^k)\n",
    "    #==INPUTS==\n",
    "    #A: coefficient polynomial array A\n",
    "    #B: coefficient polynomial array B\n",
    "    #M: m(x) is an irreducible polynomial in GF(p^k)\n",
    "    #p: prime number p in GF(p)\n",
    "    #==OUTPUTS==\n",
    "    #Returns a polynomial coefficient array of the form Poly [c_n, c_n-1, ..., c1 , c0] --> Poly = c_n*x^n + c_(n-1)*x^(n-1) + ... c1*x + c0\n",
    "\n",
    "\n",
    "    #find the degree of each polynomial\n",
    "    A_array = [x for x in bin(A)]\n",
    "    B_array = [x for x in bin(B)]\n",
    "    A_array.pop(0) #get rid of '0b' from the arrays\n",
    "    A_array.pop(0)\n",
    "    B_array.pop(0)\n",
    "    B_array.pop(0)\n",
    "    A_array = [int(x) for x in A_array] #convert each element to integer\n",
    "    B_array = [int(x) for x in B_array]\n",
    "\n",
    "    degree_A = find_degree(A_array)\n",
    "    degree_B = find_degree(B_array)\n",
    "    degree_M = find_degree(M)\n",
    "    degree_diff = abs(degree_A - degree_B)\n",
    "\n",
    "    #get rid of any preceeding zeros\n",
    "    A = remove_zeros(A_array)\n",
    "    B = remove_zeros(B_array)\n",
    "    M = remove_zeros(M)\n",
    "\n",
    "    raw_mult = raw_poly_mult_GF(A,B,p)\n",
    "    raw_degree = find_degree(raw_mult)\n",
    "    inv_M = find_highest_coefficient(M,p) # not actually an inverse, but rather what the highest degree of M equals. e.g. M == x^2 + x + 1, inv_M = 2x+2\n",
    "    \n",
    "    #if there needs to be substitutions.\n",
    "    #e.g. for GF(3), M == x^2 + x + 1, (x^2 = 2x+2), and raw_mult == x^3 --> raw_mult = x(x^2) = x(2x+2) = 2x^2 + 2x = 4x + 4 + 2x = 6x + 4 = 1 (mod 3)\n",
    "    if raw_degree >= degree_M : \n",
    "        degree_diff = raw_degree - degree_diff\n",
    "        i=0\n",
    "        while degree_diff>=0 :\n",
    "            single_elem = [0] *(len(raw_mult)-len(M)+1)\n",
    "            single_elem[0] = raw_mult[0]\n",
    "\n",
    "            raw_mult[0] = 0 #remove highest degree of raw_mult for substitution\n",
    "\n",
    "            #(single_elem*inv_M) is the replacement for the highest degree power of raw_mult, whose result\n",
    "            #is added back to raw_mult to fully apply the subsitution\n",
    "            raw_mult = MyGFAdd(raw_poly_mult_GF(inv_M,single_elem,p),raw_mult,M,p) #\n",
    "            raw_mult = remove_zeros(raw_mult) #important for reducing len(raw_mult) to proper size.\n",
    "            \n",
    "            raw_degree = find_degree(raw_mult)\n",
    "            degree_diff = raw_degree - degree_M\n",
    "            \n",
    "            i+=1\n",
    "\n",
    "    return raw_mult"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 75,
   "id": "320b47d6",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Answer = \n",
      "0xcb\n"
     ]
    }
   ],
   "source": [
    "Mx = [1,1,0,0,0,0,1,1,1]\n",
    "p = 2\n",
    "print('Answer = ')\n",
    "print(hex(int('0b'+''.join([str(x) for x in MyGFMult(0x0A,0xA3,Mx, p)]),2)))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2c1e6414",
   "metadata": {},
   "source": [
    "> $\\{0A\\} \\cdot \\{A3\\} = \\{CB\\}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 76,
   "id": "aa5e3044",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[1, 0, 1, 1, 0, 0, 1, 0]"
      ]
     },
     "execution_count": 76,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Mx2 = [1,0,0,0,1,1,0,1,1]\n",
    "p = 2\n",
    "MyGFMult(0x03,0x6E,Mx2, p)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1b22dcf6",
   "metadata": {},
   "source": [
    "(c) $\\{84\\}^{−1} = \\{73\\}$. See problem 4d for calculation."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d6af11a7",
   "metadata": {},
   "source": [
    "## Problem 4"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d3e2283e",
   "metadata": {},
   "source": [
    "Consider GF($3^{3}$) generated by $m(x) = x^{3} + 2x^{2} + 1$, and let $\\alpha$ be a root of $m(x)$. Let $a = 2\\alpha^{2} + \\alpha + 1$. We have already seen that we can find the multiplicative inverse of $a$ as follows: Let $b = b_{2}\\alpha^{2} + b_{1}\\alpha + b_{0}$ is $a^{−1}$. Then, we there should exists some $q$ such that\n",
    "\n",
    "$a \\cdot b = m(α) \\cdot q + 1.$\n",
    "\n",
    "(a) We can argue that $q$ should be of the form of $q = q_{1}\\alpha + q_{0}$. Why $q$ does not have any term of $\\alpha^{2}$?"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7f1c4cad",
   "metadata": {},
   "source": [
    ">Because then $m(\\alpha)\\cdot q$ would be of the form: $c_{4} \\alpha ^{4} + c_{3}\\alpha ^{3} + c_{2}\\alpha ^{2} c_{1} \\alpha + c_{0}$ which is a higher degree than $a \\cdot b$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7e88c6f9",
   "metadata": {},
   "source": [
    "Now, to solve for (1), we can write\n",
    "\n",
    "$a \\cdot b = (2\\alpha^{2} + \\alpha + 1) \\cdot (b_{2}\\alpha^{2} + b_{1}\\alpha + b_{0}) = (2b_{2})\\alpha^{4} + (2b_{1} + b_{2})\\alpha^{3} + (2b_{0} + b_{1} + b_{0})\\alpha^{2} + (b_{1} + b_{0})\\alpha + b_{0}$\n",
    "\n",
    "$m(α) \\cdot q + 1 = q_{1}\\alpha^{4} + (2q_{1} + q_{0})\\alpha^{3} + (2q_{0})\\alpha^{2} + q_{1}\\alpha + (q_{0}) + 1.$\n",
    "\n",
    "So, we need to determine $(b_{2}, b_{1}, b_{0})$ and $(q_{1}, q_{0})$ such that\n",
    "\n",
    "$2b_{2} = q_{1}$\n",
    "\n",
    "$2b_{1} + b_{2} = 2q_{1} + q_{0}$\n",
    "\n",
    "$2b_{0} + b_{1} + b_{0} = 2q_{0}$\n",
    "\n",
    "$b_{1} + b_{0} = q_{1}$\n",
    "\n",
    "$b_{0} = q_{0} + 1.$\n",
    "\n",
    "(b) Write the equations in (2) as a linear system, i.e, $A \\cdot x = b$, where $x = [b_{0}, b_{1}, b_{2}, q_{0}, q_{1}]^{T}$ . Find $A$ and $b$."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6eb32bf5",
   "metadata": {},
   "source": [
    "> $ A = \\begin{bmatrix}\n",
    "0 & 0 & 2 & 0 & 2 \\\\\n",
    "0 & 2 & 1 & 2 & 1 \\\\\n",
    "2 & 1 & 1 & 1 & 0 \\\\\n",
    "1 & 1 & 0 & 0 & 2 \\\\\n",
    "1 & 0 & 0 & 2 & 0 \\\\\n",
    "\\end{bmatrix}\n",
    ";   b = \\begin{bmatrix}\n",
    "0 \\\\\n",
    "0 \\\\\n",
    "0 \\\\\n",
    "0 \\\\\n",
    "1\n",
    "\\end{bmatrix}\n",
    "$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c896edf6",
   "metadata": {},
   "source": [
    "(c) Use Problem 1 to find the inverse of $A$ and solve the system for $x$. What is $b$?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "id": "5ee49ddd",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "0.0\n"
     ]
    }
   ],
   "source": [
    "A = [[1,2,3,4],[5,6,7,8],[9,10,11,12],[16,15,14,13]]\n",
    "print(recursive_det(A,3))\n",
    "print(np.linalg.det(A)%3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "id": "5dc484df",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2\n",
      "2\n",
      "A = \n",
      "[[0 0 2 0 2]\n",
      " [0 2 1 2 1]\n",
      " [2 1 1 1 0]\n",
      " [1 1 0 0 2]\n",
      " [1 0 0 2 0]]\n",
      "\n",
      "A^-1 = \n",
      "[[2 0 2 1 2]\n",
      " [0 2 1 2 2]\n",
      " [0 1 0 1 2]\n",
      " [2 0 2 1 1]\n",
      " [2 2 0 2 1]]\n",
      "\n",
      "A*A^-1 = \n",
      "[[1 0 0 0 0]\n",
      " [0 1 0 0 0]\n",
      " [0 0 1 0 0]\n",
      " [0 0 0 1 0]\n",
      " [0 0 0 0 1]]\n",
      "\n",
      "x = \n",
      "[[2]\n",
      " [2]\n",
      " [2]\n",
      " [1]\n",
      " [1]]\n"
     ]
    }
   ],
   "source": [
    "A = [[0,0,2,0,2],[0,2,1,2,1],[2,1,1,1,0],[1,1,0,0,2],[1,0,0,2,0]] #\"A\" from 4b)\n",
    "print(recursive_det(A,3))\n",
    "print(round(np.linalg.det(A)%3))\n",
    "b = [[0],[0],[0],[0],[1]] #\"b\" from 4b)\n",
    "inv_a = inverse_rec(A,3)\n",
    "print('A = ')\n",
    "print(np.array(A))\n",
    "print('')\n",
    "print('A^-1 = ')\n",
    "print(inv_a)\n",
    "print('')\n",
    "print('A*A^-1 = ')\n",
    "print(np.matmul(inv_a,A)%3)\n",
    "print('')\n",
    "print('x = ')\n",
    "print(np.matmul(inv_a,b))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "681c2c9e",
   "metadata": {},
   "source": [
    "> $ A^{-1} = \\begin{bmatrix}\n",
    "2 & 0 & 2 & 1 & 2 \\\\\n",
    "0 & 2 & 1 & 2 & 2 \\\\\n",
    "0 & 1 & 0 & 1 & 2 \\\\\n",
    "2 & 0 & 2 & 1 & 1 \\\\\n",
    "2 & 2 & 0 & 2 & 1 \\\\\n",
    "\\end{bmatrix} $\n",
    "\n",
    "> $ A \\cdot A^{-1} = \n",
    "\\begin{bmatrix}\n",
    "1 & 0 & 0 & 0 & 0 \\\\\n",
    "0 & 1 & 0 & 0 & 0 \\\\\n",
    "0 & 0 & 1 & 0 & 0 \\\\\n",
    "0 & 0 & 0 & 1 & 0 \\\\\n",
    "0 & 0 & 0 & 0 & 1 \\\\\n",
    "\\end{bmatrix}$\n",
    "\n",
    "> $ x = \n",
    "\\begin{bmatrix}\n",
    "2 \\\\\n",
    "2 \\\\\n",
    "2 \\\\\n",
    "1 \\\\\n",
    "1 \\\\\n",
    "\\end{bmatrix}$\n",
    "\n",
    "$b_{0} = 2, b_{1} = 2, b_{2} = 2$\n",
    "\n",
    "$b = 2\\alpha^{2} + 2\\alpha + 2$"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2a0f42b3",
   "metadata": {},
   "source": [
    "(d) Write a computer program that computes the multiplicative inverse in a general finite field. Consider GF($p^{k}$) generated by $m(x) = m_{k}x^{k} + · · · + m_{1}x + m_{0}$ which is represented as $M = [m_{k}, m_{k−1}, . . . , m_{1}, m_{0}]$. For a value $a = a_{k−1}\\alpha^{k−1} + · · · + a_{1}\\alpha + a_{0}$, which is represented as $A = [a_{k−1}, . . . , a_{1}, a_{0}]$, compute $b = b_{k−1}\\alpha^{k−1} + · · · + b_{1}\\alpha + b_{0}$ such that $a \\cdot b = 1$ and return $B = [b_{k−1}, . . . , b_{1}, b_{0}]$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 79,
   "id": "8319ee27",
   "metadata": {},
   "outputs": [],
   "source": [
    "def poly_str(polynomial) :\n",
    "    #Returns a string statement to be printed. Provides a more readable version than just printing a polynomial.\n",
    "    i=0\n",
    "    output = ''\n",
    "    while i<len(polynomial) :\n",
    "        j=0\n",
    "        while j<len(polynomial[i]) :\n",
    "            if polynomial[i][j] != [0] :\n",
    "                if j==0 :\n",
    "                    output = output + '('\n",
    "                output = output + str(polynomial[i][j][0])+'*'+str(polynomial[i][j][1])\n",
    "                if j!=(len(polynomial[i])-1) :\n",
    "                    output = output + ' + '\n",
    "            j+=1\n",
    "        output = output + ')'\n",
    "        if i!=len(polynomial)-1 :\n",
    "            output = output + '*α^' + str(len(polynomial)-i-1)\n",
    "            output = output + ' + ' \n",
    "        i+=1\n",
    "\n",
    "    return output\n",
    "\n",
    "\n",
    "def mult_inverse_GF(a, M, p) :\n",
    "    #Returns the multiplicative inverse of a in a general finite field GF(p^k) generated by m(x) = mkx^k + ... m1x + m0\n",
    "    #INPUTS:\n",
    "    #a[]: polynomial coefficients in GF(p) to find inverse of\n",
    "    #M[]: polynomial coefficients of m(x)\n",
    "    #p: p in GF(p^k), integer.\n",
    "\n",
    "    degree_a = find_degree(a)\n",
    "    degree_M = find_degree(M)\n",
    "    mult = [[] for x in range(0,2*len(a)-1)] # matrix result of a*b.\n",
    "\n",
    "    b = [['b' + str(k-1), k-1] for k in range(degree_M,0,-1)] # [['bk-1',k-1],['bk-1',k-2],...,['b1',1],['b0',0]]\n",
    "    \n",
    "    if degree_a<degree_M-1 : \n",
    "        #add zeros to make \"a\" of degree k-1, where k is the degree of M.\n",
    "        diff=degree_M-degree_a\n",
    "        i=0\n",
    "        while i<diff : \n",
    "            #add zeros to make \"a\" of degree k-1, where k is the degree of M.\n",
    "            a.insert(0,0)\n",
    "            i+=1\n",
    "\n",
    "    i=0\n",
    "    while i<len(a) : #Multiply coefficents of a by coefficients of b, but don't add up like values yet.\n",
    "        #Store them in 2d array \"mult\", of the form: mult[<descending polynomial dexponent>][<numeric coefficient>,<any of {b0, b1, ..., bk-1}>]\n",
    "        # e.g. if mult = [[0],[[2,'b1'],[1,'b2']],[1,'b0']] ---> 0*alpha^2 + (2b1 + 1b2)*alpha^1 + 1b0\n",
    "        j=0\n",
    "        while j<len(b) :\n",
    "            exp_a = len(a)-i #exponent value of current element of a\n",
    "            exp_b = b[j][1] #exponent value of current element of b\n",
    "            exp_total = exp_a + exp_b\n",
    "            if a[i] != 0 : # do not multiply if coefficient of a=0\n",
    "                mult[len(mult)-exp_total].append([a[i],b[j][0]])\n",
    "            #print(mult)\n",
    "            j+=1\n",
    "        i+=1\n",
    "\n",
    "\n",
    "    i=0\n",
    "    while i<len(mult) : #replace unused elements with zero\n",
    "        if mult[i] == [] : mult[i] = [[0]]\n",
    "        i+=1\n",
    "\n",
    "    i=0\n",
    "    while i<len(mult) : #remove preceeding zeros.\n",
    "        #e.g. if mult = [[[0]],[[1,'b0],[2,'b1']], [[0]]]\n",
    "        #   then mult = [[[1,'b0],[2,'b1']], [[0]]]\n",
    "        if mult[i] != [[0]] :\n",
    "            break\n",
    "        mult.pop(0)\n",
    "        i+=1\n",
    "\n",
    "    #Replace all elements > degree M with elements in lower degree.\n",
    "    inv_M = find_highest_coefficient(M,p) # not actually an inverse, but rather what the highest degree of M equals. e.g. M == x^2 + x + 1, inv_M = 2x+2\n",
    "    degree_mult = find_degree(mult)\n",
    "    if degree_mult >= degree_M : \n",
    "        degree_diff = degree_mult - degree_M\n",
    "        while degree_diff>=0 :\n",
    "            mult = remove_zeros(mult)\n",
    "            single_elem = [0]\n",
    "            single_elem[0] = mult[0]\n",
    "            mult[0] = 0 #remove highest degree of raw_mult for substitution\n",
    "            mult = remove_zeros(mult)\n",
    "            if single_elem[0][0] != [0] :\n",
    "                #(single_elem*inv_M) is the replacement for the highest degree power of raw_mult, whose result\n",
    "                #is added back to raw_mult to fully apply the subsitution\n",
    "                i=0\n",
    "                while i<len(inv_M) :\n",
    "                    j=0\n",
    "                    while j<len(single_elem[0]) :\n",
    "                        exp_M = len(inv_M)-i-1 #exponent value of current element of a\n",
    "                        exp_single_elem = len(mult)-degree_M #exponent value of current element of b\n",
    "                        exp_total = exp_M + exp_single_elem\n",
    "                        index = len(mult) - exp_total - 1\n",
    "\n",
    "                        if inv_M[i] != 0: # do not multiply if coefficient of a=0\n",
    "                            mult[index].append([inv_M[i]*single_elem[0][j][0],single_elem[0][j][1]])\n",
    "                        j+=1\n",
    "                    i+=1\n",
    "            degree_mult = find_degree(mult)\n",
    "            degree_diff = degree_mult - degree_M # loop condition\n",
    "\n",
    "    i=0\n",
    "    while i<len(mult) :\n",
    "        j=0\n",
    "        while j<len(mult[i]) :\n",
    "            k=j\n",
    "            if mult[i][j] != [0] :\n",
    "                total = mult[i][j][0] # running total for how many times that item occurs\n",
    "                coefficient = mult[i][j][1] #e.g. \"b0\" or \"b1\"\n",
    "                while k<len(mult[i]) :\n",
    "                    if mult[i][k] != [0] :\n",
    "                        if (coefficient == mult[i][k][1]) & (k!=j) : \n",
    "                            total = total + mult[i][k][0]\n",
    "                            mult[i][k] = [0]\n",
    "                    k+=1\n",
    "                mult[i][j][0] = total\n",
    "            j+=1\n",
    "        i+=1\n",
    "\n",
    "    #remove empty values\n",
    "    i=0\n",
    "    while i<len(mult) :\n",
    "        j=0\n",
    "        while j<len(mult[i]) :\n",
    "            if mult[i][j] == [0] :\n",
    "                mult[i].pop(j)\n",
    "                j-=1\n",
    "            j+=1\n",
    "        i+=1\n",
    "\n",
    "\n",
    "    solution_matrix = [[0 for x in range(0,len(mult)+1)] for y in range(0,len(mult))]\n",
    "    i=0\n",
    "    while i<len(mult) :\n",
    "        j=0\n",
    "        while j<len(mult[i]) :\n",
    "            coefficient_index = int(mult[i][j][1][1:None])\n",
    "            solution_matrix[i][coefficient_index] = np.mod((mult[i][j][0]),p)\n",
    "            j+=1\n",
    "        if i!=(len(mult)-1) :\n",
    "            solution_matrix[i][len(solution_matrix)] = 0\n",
    "        else : solution_matrix[i][len(solution_matrix)] = 1\n",
    "        i+=1\n",
    "\n",
    "    #Create a copy of solution_matrix, called solution_matrix_2 which uses the Matrix Inverse Method\n",
    "    solution_matrix_2 = [[0 for x in range(0,len(mult))] for y in range(0,len(mult))]\n",
    "    i=0\n",
    "    while i<len(solution_matrix) :\n",
    "        j=0\n",
    "        while j<(len(solution_matrix)) :\n",
    "            solution_matrix_2[i][j] = solution_matrix[i][j]\n",
    "            j+=1\n",
    "        i+=1\n",
    "\n",
    "    #Calculate the RREF of solution_matrix to get the answer\n",
    "    solution_matrix = sp.Matrix(solution_matrix)\n",
    "    RREF = solution_matrix.rref()\n",
    "    output = [0 for x in range(0,len(RREF[0][0:None,0]))]\n",
    "    i=0\n",
    "    while i<(len(RREF[0][0:None,0])) :\n",
    "        j=0\n",
    "        while j<(len(RREF[0][i,0:None])-1) :\n",
    "            output[len(output)-i-1] = np.mod(RREF[0][i,(len(RREF[0][i,0:None])-1)],p)\n",
    "            j+=1\n",
    "\n",
    "        i+=1\n",
    "\n",
    "    #Calculate answer using Inverse Matrix\n",
    "    inv_matrix = inverse_rec(solution_matrix_2,p)\n",
    "    one_matrix = [[0] for x in range(0,len(solution_matrix_2))]\n",
    "    one_matrix[len(one_matrix)-1] = [1] # matrix of form: [[0], [0], [0], ... , [0], [1]] (one_matrix = 1 in a system of equations.)\n",
    "    output2 = np.matmul(np.array(inv_matrix),np.array(one_matrix))\n",
    "    output2 = np.flip(output2.transpose())\n",
    "\n",
    "\n",
    "    return output, output2\n",
    "\n",
    "#prints out 8 bit array into hexadecimal based on the answer format from mult_inverse_GF()\n",
    "def print_mult_inv_ans(input, answer) :\n",
    "    input_int = int(str(hex(int('0b'+''.join([str(x) for x in input]),2))),0)\n",
    "    # method1 = int(str(hex(int('0b'+''.join([str(x) for x in answer[0]]),2))),0)\n",
    "    method2 = int(str(hex(int('0b'+''.join([str(x) for x in answer[1][0]]),2))),0)\n",
    "    lb = '{'\n",
    "    rb = '}'\n",
    "    print(f'Finding {lb}{hex(input_int)}{rb}^(-1),')\n",
    "    # print(f'using RREF, answer = {lb}{hex(method1)}{rb}')\n",
    "    print(f'using Matrix Inverse, answer = {lb}{hex(method2)}{rb}')\n",
    "    return"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "id": "5da7f414",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Finding {0xb}^(-1),\n",
      "using Matrix Inverse, answer = {0x5}\n"
     ]
    }
   ],
   "source": [
    "a =   [1,0,1,1] #{0xb}\n",
    "M = [1,0,0,1,1]\n",
    "answer = mult_inverse_GF(a,M,2)\n",
    "print_mult_inv_ans(a,answer)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 81,
   "id": "439ff865",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Finding {0x5}^(-1),\n",
      "using Matrix Inverse, answer = {0xb}\n"
     ]
    }
   ],
   "source": [
    "b =   [0,1,0,1] #{0x5}\n",
    "M = [1,0,0,1,1]\n",
    "answer2 = mult_inverse_GF(b,M,2)\n",
    "print_mult_inv_ans(b,answer2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "dfded8c8",
   "metadata": {},
   "source": [
    "### Calculating 3c"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 82,
   "id": "017f8b28",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Finding {0x84}^(-1),\n",
      "using Matrix Inverse, answer = {0x73}\n"
     ]
    }
   ],
   "source": [
    "\n",
    "a = [1,0,0,0,0,1,0,0] # = {84}\n",
    "Mx = [1,1,0,0,0,0,1,1,1] # m(x) = x^8 + x^7 + x^2 + x + 1\n",
    "answer = mult_inverse_GF(a,Mx,2)\n",
    "print_mult_inv_ans(a,answer) #print it out all nice\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 83,
   "id": "41b69556",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Finding {0x73}^(-1),\n",
      "using Matrix Inverse, answer = {0x84}\n"
     ]
    }
   ],
   "source": [
    "input = [0,1,1,1,0,0,1,1] # = {73}\n",
    "Mx = [1,1,0,0,0,0,1,1,1] # m(x) = x^8 + x^7 + x^2 + x + 1\n",
    "answer2 = mult_inverse_GF(input,Mx,2)\n",
    "print_mult_inv_ans(input,answer2) #print it out all nice"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "429da908",
   "metadata": {},
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "56278cd9",
   "metadata": {},
   "source": [
    "## Problem 5"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9ccd902c",
   "metadata": {},
   "source": [
    "Recall the S-box used in the Advanced Encryption Standard (AES). It maps a pair $\\{xy\\}$ to $\\{x^{′}y^{′}\\}$ = $S(\\{xy\\})$ according to the following steps:\n",
    "\n",
    "• You first find $\\{uv\\} = \\{xy\\}^{−1}$.\n",
    "\n",
    "• Write $\\{uv\\}$ in binary $\\{uv\\} = [b_{7}, . . . , b_{0}]$.\n",
    "\n",
    "• Compute\n",
    "\n",
    "$ \\begin{bmatrix}\n",
    "b_{0}^{'} \\\\\n",
    "b_{1}^{'} \\\\\n",
    "b_{2}^{'} \\\\\n",
    "b_{3}^{'} \\\\\n",
    "b_{4}^{'} \\\\\n",
    "b_{0}^{'} \\\\\n",
    "b_{0}^{'} \\\\\n",
    "b_{0}^{'} \\\\\n",
    "\\end{bmatrix} = \\begin{bmatrix}\n",
    "1 & 0 & 0 & 0 & 1 & 1 & 1 & 1 \\\\\n",
    "1 & 1 & 0 & 0 & 0 & 1 & 1 & 1 \\\\\n",
    "1 & 1 & 1 & 0 & 0 & 0 & 1 & 1 \\\\\n",
    "1 & 1 & 1 & 1 & 0 & 0 & 0 & 1 \\\\\n",
    "1 & 1 & 1 & 1 & 1 & 0 & 0 & 0 \\\\\n",
    "0 & 1 & 1 & 1 & 1 & 1 & 0 & 0 \\\\\n",
    "0 & 0 & 1 & 1 & 1 & 1 & 1 & 0 \\\\\n",
    "0 & 0 & 0 & 1 & 1 & 1 & 1 & 1 \n",
    "\\end{bmatrix} \\cdot\n",
    "\\begin{bmatrix}\n",
    "b_{0} \\\\\n",
    "b_{1} \\\\\n",
    "b_{2} \\\\\n",
    "b_{3} \\\\\n",
    "b_{4} \\\\\n",
    "b_{5} \\\\\n",
    "b_{6} \\\\\n",
    "b_{7} \n",
    "\\end{bmatrix} +\n",
    "\\begin{bmatrix}\n",
    "1 \\\\\n",
    "1 \\\\\n",
    "0 \\\\\n",
    "0 \\\\\n",
    "0 \\\\\n",
    "1 \\\\\n",
    "1 \\\\\n",
    "0 \n",
    "\\end{bmatrix}$ (mod 2)\n",
    "\n",
    "• Compute $\\{x^{′}y^{′}\\} = [b^{′}_{7}, . . . , b^{′}_{0}]$.\n",
    "\n",
    "Write a computer program to generate the S-box. It should take $\\{xy\\}$ and return $\\{x^{′}y^{′}\\}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 196,
   "id": "89b389ef",
   "metadata": {},
   "outputs": [],
   "source": [
    "def hex_to_bin_list(hex, length) :\n",
    "    #Converts hexadecimal value(or a decimal) \n",
    "    hex = list(bin(hex))\n",
    "    hex.pop(0)\n",
    "    hex.pop(0)\n",
    "    if len(hex)<length:\n",
    "        diff = (length-len(hex))\n",
    "        for i in range(0,diff) :\n",
    "            hex.insert(0,'0')\n",
    "\n",
    "    for i in range(0,len(hex)) :\n",
    "        hex[i] = int(hex[i])\n",
    "    return hex\n",
    "\n",
    "def bin_list_to_hex(bin) :\n",
    "    #converts a binary list to a hexidecimal number.\n",
    "    return hex(int('0b'+''.join([str(x) for x in bin]),2))\n",
    "\n",
    "def S(input) :\n",
    "    M = [1,0,0,0,1,1,0,1,1]\n",
    "    l = np.array([[1],[1],[0],[0],[0],[1],[1],[0]]) #63\n",
    "    S_array = np.array([[1,0,0,0,1,1,1,1],[1,1,0,0,0,1,1,1],[1,1,1,0,0,0,1,1],[1,1,1,1,0,0,0,1],[1,1,1,1,1,0,0,0],[0,1,1,1,1,1,0,0],[0,0,1,1,1,1,1,0],[0,0,0,1,1,1,1,1]])\n",
    "\n",
    "    input = hex_to_bin_list(input, len(M))\n",
    "    inverse = mult_inverse_GF(input,M,2)\n",
    "    inverse = np.vstack(np.flip(np.array(inverse[1][0])))\n",
    "    output = np.matmul(S_array,inverse) % 2\n",
    "    output = (output + l) % 2\n",
    "    output = np.flip(np.reshape(output,8))\n",
    "    return bin_list_to_hex(output)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 262,
   "id": "40efdc81",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Input = {0x95}\n",
      "Output = {0x2a}\n"
     ]
    }
   ],
   "source": [
    "# Testing S(input) to see if it works.\n",
    "\n",
    "inp = 0x95 #{95}\n",
    "outp = S(inp) #expected answer: {2A}\n",
    "lb = '{'\n",
    "rb = '}'\n",
    "print(f'Input = {lb}{hex(inp)}{rb}')\n",
    "print(f'Output = {lb}{outp}{rb}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 254,
   "id": "63d5f69d",
   "metadata": {},
   "outputs": [],
   "source": [
    "def generate_SBox() :\n",
    "    matrix = [['00' for x in range(0,16)] for y in range(0,16)]\n",
    "    for i in range(0,16) :\n",
    "        for j in range(0,16) :\n",
    "            value = np.concatenate((np.array(hex_to_bin_list(i, 4)),np.array(hex_to_bin_list(j, 4))),axis=0)\n",
    "            value = int(bin_list_to_hex(value),16)\n",
    "            matrix[i][j] = S((value))\n",
    "    return np.array(matrix)\n",
    "            "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 255,
   "id": "bd895fe3",
   "metadata": {},
   "outputs": [],
   "source": [
    "#This takes 13.5 minutes to compute!!\n",
    "SBox_array = generate_SBox()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 261,
   "id": "062478a7",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "   | 0    1    2    3    4    5    6    7    8    9    a    b    c    d    e    f   \n",
      "---|-------------------------------------------------------------------------------\n",
      " 0 |0x63 0x7c 0x77 0x7b 0xf2 0x6b 0x6f 0xc5 0x30 0x01 0x67 0x2b 0xfe 0xd7 0xab 0x76 \n",
      " 1 |0xca 0x82 0xc9 0x7d 0xfa 0x59 0x47 0xf0 0xad 0xd4 0xa2 0xaf 0x9c 0xa4 0x72 0xc0 \n",
      " 2 |0xb7 0xfd 0x93 0x26 0x36 0x3f 0xf7 0xcc 0x34 0xa5 0xe5 0xf1 0x71 0xd8 0x31 0x15 \n",
      " 3 |0x04 0xc7 0x23 0xc3 0x18 0x96 0x05 0x9a 0x07 0x12 0x80 0xe2 0xeb 0x27 0xb2 0x75 \n",
      " 4 |0x09 0x83 0x2c 0x1a 0x1b 0x6e 0x5a 0xa0 0x52 0x3b 0xd6 0xb3 0x29 0xe3 0x2f 0x84 \n",
      " 5 |0x53 0xd1 0x00 0xed 0x20 0xfc 0xb1 0x5b 0x6a 0xcb 0xbe 0x39 0x4a 0x4c 0x58 0xcf \n",
      " 6 |0xd0 0xef 0xaa 0xfb 0x43 0x4d 0x33 0x85 0x45 0xf9 0x02 0x7f 0x50 0x3c 0x9f 0xa8 \n",
      " 7 |0x51 0xa3 0x40 0x8f 0x92 0x9d 0x38 0xf5 0xbc 0xb6 0xda 0x21 0x10 0xff 0xf3 0xd2 \n",
      " 8 |0xcd 0x0c 0x13 0xec 0x5f 0x97 0x44 0x17 0xc4 0xa7 0x7e 0x3d 0x64 0x5d 0x19 0x73 \n",
      " 9 |0x60 0x81 0x4f 0xdc 0x22 0x2a 0x90 0x88 0x46 0xee 0xb8 0x14 0xde 0x5e 0x0b 0xdb \n",
      " a |0xe0 0x32 0x3a 0x0a 0x49 0x06 0x24 0x5c 0xc2 0xd3 0xac 0x62 0x91 0x95 0xe4 0x79 \n",
      " b |0xe7 0xc8 0x37 0x6d 0x8d 0xd5 0x4e 0xa9 0x6c 0x56 0xf4 0xea 0x65 0x7a 0xae 0x08 \n",
      " c |0xba 0x78 0x25 0x2e 0x1c 0xa6 0xb4 0xc6 0xe8 0xdd 0x74 0x1f 0x4b 0xbd 0x8b 0x8a \n",
      " d |0x70 0x3e 0xb5 0x66 0x48 0x03 0xf6 0x0e 0x61 0x35 0x57 0xb9 0x86 0xc1 0x1d 0x9e \n",
      " e |0xe1 0xf8 0x98 0x11 0x69 0xd9 0x8e 0x94 0x9b 0x1e 0x87 0xe9 0xce 0x55 0x28 0xdf \n",
      " f |0x8c 0xa1 0x89 0x0d 0xbf 0xe6 0x42 0x68 0x41 0x99 0x2d 0x0f 0xb0 0x54 0xbb 0x16 \n"
     ]
    }
   ],
   "source": [
    "def print_SBox(SBox) : # Prints an S box out.\n",
    "    topline = '   |'\n",
    "    for num in range(0,16) :\n",
    "            topline = topline + ' ' + hex(num)[2] + '   '\n",
    "    print(topline)\n",
    "    divider = '---|'\n",
    "    for y in range(0,79) :\n",
    "        divider = divider + '-'\n",
    "    print(divider)\n",
    "    for i in range(0,len(SBox)) :\n",
    "        string = ' ' + hex(i)[2] + ' |'\n",
    "        for j in range(0,len(SBox[i])) :\n",
    "            value = SBox[i][j]\n",
    "            if len(value)<4 : # adds in necessary zeros e.g. 0x7 --> 0x07\n",
    "                 value = value[:2] + '0' + value[2:None]\n",
    "            string = string + value + ' '\n",
    "        print(string)\n",
    "    return\n",
    "\n",
    "print_SBox(SBox_array)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "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.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}