{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "

cs1001.py , Tel Aviv University, Fall 2018-2019

\n", "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Recitation 7" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We reviewed some properties of prime numbers and used them for primality testing. We reviewed the Diffie-Hellman protocol for finding a shared secret key and also tried to crack it. \n", "Then, we discussed OOP and special methods that allow us to define operators in Python.\n", "\n", "### Takeaways:\n", "
    \n", "
  1. The probabilistic function is_prime, that uses Fermat's primality test, can be used to detect primes quickly and efficiently, but has a (very small) probability of error. Its time complexity is $O(n^3)$, where $n$ is the number of bits of the input.
    2. \n", "
  2. The DH protocol relies on two main principles: the following equality $(g^{a}\\mod p)^b \\mod p = g^{ab} \\mod p $ and the (believed) hardness of the discrete log problem (given $g,p$, and $x = g^{a} \\mod p$ finding $a$ is hard). Make sure you understand the details of the protocol.
    3. \n", "
  3. OOP allows us to create classes at will and to define complex objects. Remember that the most important thing is to find a good inner representation of the object so that you will have an easy time working with it.
    4. \n", "
" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Code for printing several outputs in one cell (not part of the recitation):" ] }, { "cell_type": "code", "execution_count": 66, "metadata": {}, "outputs": [], "source": [ "from IPython.core.interactiveshell import InteractiveShell\n", "InteractiveShell.ast_node_interactivity = \"all\"" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Primes and Diffie-Hellman" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Primality test using Ferma's witness\n", "\n", "Fermat's little theorem: if $p$ is prime and $1 < a < p$, then $a^{p-1} (\\textrm{mod}\\ p) \\equiv 1$\n", "\n", "On the other hand: If $p$ is composite, write $p = ab$ for some $a,b>1$. Now, assume toward contradiction that $a^{ab-1} (\\textrm{mod}\\ ab) \\equiv 1$. It follows that $ ab ~|~a^{ab-1} -1$ and therefore $a~|~a^{ab-1} -1$. However, as clearly $a ~|~a^{ab-1}$ we must have $a~|~1$, a contradiction.\n", "\n", "We can use Fermat's little theorem in order test whether a given number is prime. Note that if the number has $n$ bits than testing all possible $a$-s will require $O(2^n)$ iterations (a lot!).\n", "\n", "Instead, we will try 100 random $a$-s in the range and see if one works as a witness." ] }, { "cell_type": "code", "execution_count": 67, "metadata": {}, "outputs": [], "source": [ "import random\n", "\n", "def is_prime(m, show_witness=False):\n", "\n", " \"\"\" probabilistic test for m's compositeness \"\"\"''\n", "\n", " for i in range(0, 100):\n", " a = random.randint(1, m - 1) # a is a random integer in [1..m-1]\n", " if pow(a, m - 1, m) != 1:\n", " if show_witness: # caller wishes to see a witness\n", " print(m, \"is composite\", \"\\n\", a, \"is a witness, i=\", i + 1)\n", " return False\n", "\n", " return True\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For $a,b,c$ of at most $n$ bits each, time complexity of modpower is $O(n^3)$" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "def modpower(a, b, c):\n", " \"\"\" computes a**b modulo c, using iterated squaring \"\"\"\n", " result = 1\n", " while b > 0: # while b is nonzero\n", " if b % 2 == 1: # b is odd\n", " result = (result * a) % c\n", " a = (a*a) % c\n", " b = b // 2\n", " return result" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Runtime analysis:\n", "\n", "* The main loop runs over $b$, dividing $b \\to b/2$ at each iteration, so it runs $O(n)$ times.\n", "* In each iteration we do: \n", " * One operation of $b%2$ in $O(1)$ time\n", " * One operation of $b/2$ in $O(1)$ time\n", " * At most two multiplication and two modulu operations\n", " * Multiplication of two $n$ bit numbers runs in time $O(n^2)$\n", " * Modulu can be implemented by addition, division and multiplication: $a \\textrm{ mod } b = a - (a // b) b$ and division runs in time $O(n^2)$ same as multiplication\n", " * Finally, the modulu operation keeps all numbers at most $n$ bits, thus the running time does not increase with each iteration\n", "* In total - $O(n^3)$\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### The probability of error:\n", "First, notice that if the function says that an imput number $m$ is not prime, then it is true. \n", "The function can make a mistake only is the case where a number $m$ is not prime, and is excidentally categorized by the function as prime. This can happen if all $100$ $a$'s that the function tried were not witnesses. \n", "\n", "A quick computation shows that if $m$ is **not** a Charmichael number then at least $\\frac{1}{2}$ of all possible $a$s are witnesses, so in almost all cases the probability for error is $(\\frac{1}{2})^{100}$ (this is extremely low)." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Testing the prime number theorem: For a large n, a number of n bits is prime with a prob. of O(1/n)\n", "We decide on the size of the sample (to avoid testing all possible $2^{n-1}$ numbers of $n$ bits) and test whether each number we sample is prime. Then we divide the number of primes with the size of the sample." ] }, { "cell_type": "code", "execution_count": 68, "metadata": {}, "outputs": [], "source": [ "def prob_prime(n, sample):\n", " cnt = 0\n", " for i in range(sample):\n", " m = random.randint(2**(n-1), 2**n - 1)\n", " cnt += is_prime(m)\n", " return cnt / sample\n", "\n" ] }, { "cell_type": "code", "execution_count": 69, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1.0" ] }, "execution_count": 69, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "0.5057" ] }, "execution_count": 69, "metadata": {}, "output_type": "execute_result" } ], "source": [ "prob_prime(2, 10**4)\n", "prob_prime(3, 10**4)\n" ] }, { "cell_type": "code", "execution_count": 70, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.0148" ] }, "execution_count": 70, "metadata": {}, "output_type": "execute_result" } ], "source": [ "prob_prime(100, 10**4)" ] }, { "cell_type": "code", "execution_count": 71, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "0.0081" ] }, "execution_count": 71, "metadata": {}, "output_type": "execute_result" } ], "source": [ "prob_prime(200, 10**4)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Diffie Hellman from lecture" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### The protocol as code" ] }, { "cell_type": "code", "execution_count": 72, "metadata": {}, "outputs": [], "source": [ "def DH_exchange(p):\n", " \"\"\" generates a shared DH key \"\"\"\n", " g = random.randint(1, p - 1)\n", " a = random.randint(1, p - 1)# Alice's secret\n", " b = random.randint(1, p - 1)# Bob's secret\n", " x = pow(g, a, p)\n", " y = pow(g, b, p)\n", " key_A = pow(y, a, p)\n", " key_B = pow(x, b, p)\n", " #the next line is different from lecture\n", " return g, a, b, x, y, key_A #key_A=key_B" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Find a prime number" ] }, { "cell_type": "code", "execution_count": 73, "metadata": {}, "outputs": [], "source": [ "def find_prime(n):\n", " \"\"\" find random n-bit long prime \"\"\"\n", " while(True):\n", " candidate = random.randrange(2**(n-1), 2**n)\n", " if is_prime(candidate):\n", " return candidate" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Demostration:" ] }, { "cell_type": "code", "execution_count": 98, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "557\n" ] }, { "data": { "text/plain": [ "(247, 88, 353, 367, 2, 206)" ] }, "execution_count": 98, "metadata": {}, "output_type": "execute_result" } ], "source": [ "import random\n", "p = find_prime(10)\n", "print(p)\n", "g, a, b, x, y, key = DH_exchange(p)\n", "g, a, b, x, y, key" ] }, { "cell_type": "code", "execution_count": 99, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "367\n", "206\n" ] } ], "source": [ "print(pow(g, a, p))\n", "print(pow(x, b, p))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Crack the Diffie Hellman code\n", "There is no known way to find $a$ efficiently, so we try the naive one: iterating over all $a$-s and cheking whether the equation $g^a \\mod p = x$ holds for them. \n", "\n", "If we found $a'$ that satisfies the condition but is not the original $a$, does it matter?\n", "\n", "The time complexity of crack_DH is $O(2^nn^3)$" ] }, { "cell_type": "code", "execution_count": 76, "metadata": {}, "outputs": [], "source": [ "def crack_DH(p, g, x):\n", " ''' find secret \"a\" that satisfies g**a%p == x\n", " Not feasible for large p '''\n", " for a in range(1, p - 1):\n", " if a % 100000 == 0:\n", " print(a) #just to estimate running time\n", " if pow(g, a, p) == x:\n", " return a\n", " return None #should never get here" ] }, { "cell_type": "code", "execution_count": 100, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "88\n" ] }, { "data": { "text/plain": [ "88" ] }, "execution_count": 100, "metadata": {}, "output_type": "execute_result" } ], "source": [ "#g, a, b, x, y, key = DH_exchange(p)\n", "print(a)\n", "crack_DH(p, g, x)\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Cracking DH with different values of $a$ (different private keys)\n", "\n", "The algorithm crack_DH can return a different private key ($a$) than the one chosen by Alice, i.e. - $crack\\_DH(p,g,x) = a' \\neq a$, however, in this case we have, by definition of the cracking algorithm: $g^{a'} \\textrm{ mod } p = g^{a} \\textrm{ mod } p$, thus:\n", "\n", "$$y^{a'} \\textrm{ mod } p = \\left(g^{b}\\right)^{a'} \\textrm{ mod } p = \\left(g^{a'}\\right)^{b} \\textrm{ mod } p = \\left(g^{a} \\textrm{ mod } p \\right)^{b} \\textrm{ mod } p = g^{ab} \\textrm{ mod } p$$\n", "\n", "I.e. - we can still compute the shared secret!" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Trying to crack the protocol with a 100 bit prime" ] }, { "cell_type": "code", "execution_count": 78, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "972762959355461077758092346931\n", "711085486200126306605788797672 301337617630740085784318464474 953829050690741734973058379801 673125365307331009336129574105 220291338338430064459432090221 203856719177576397001814724913\n", "100000\n", "200000\n", "300000\n", "400000\n", "500000\n", "600000\n", "700000\n", "800000\n", "900000\n", "1000000\n", "1100000\n", "1200000\n", "1300000\n", "1400000\n", "1500000\n", "1600000\n" ] }, { "ename": "KeyboardInterrupt", "evalue": "", "output_type": "error", "traceback": [ "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[1;31mKeyboardInterrupt\u001b[0m Traceback (most recent call last)", "\u001b[1;32m\u001b[0m in \u001b[0;36m\u001b[1;34m()\u001b[0m\n\u001b[0;32m 5\u001b[0m \u001b[0mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mg\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0ma\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mb\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0my\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mkey\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 6\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 7\u001b[1;33m \u001b[0mcrack_DH\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mp\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mg\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[1;32m\u001b[0m in \u001b[0;36mcrack_DH\u001b[1;34m(p, g, x)\u001b[0m\n\u001b[0;32m 3\u001b[0m Not feasible for large p '''\n\u001b[0;32m 4\u001b[0m \u001b[1;32mfor\u001b[0m \u001b[0ma\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mrange\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m1\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mp\u001b[0m \u001b[1;33m-\u001b[0m \u001b[1;36m1\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 5\u001b[1;33m \u001b[1;32mif\u001b[0m \u001b[0ma\u001b[0m \u001b[1;33m%\u001b[0m \u001b[1;36m100000\u001b[0m \u001b[1;33m==\u001b[0m \u001b[1;36m0\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m 6\u001b[0m \u001b[0mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0ma\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;31m#just to estimate running time\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 7\u001b[0m \u001b[1;32mif\u001b[0m \u001b[0mpow\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mg\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0ma\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mp\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;33m==\u001b[0m \u001b[0mx\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n", "\u001b[1;31mKeyboardInterrupt\u001b[0m: " ] } ], "source": [ "import random\n", "p = find_prime(100)\n", "print(p)\n", "g, a, b, x, y, key = DH_exchange(p)\n", "print(g, a, b, x, y, key)\n", "\n", "crack_DH(p, g, x)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Analyzing the nubmer of years it will take to crack the protocol if $a$ is found at the end (assuming iterating over 100000 $a$s takes a second)" ] }, { "cell_type": "code", "execution_count": 79, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "301337617630740085784318464474" ] }, "execution_count": 79, "metadata": {}, "output_type": "execute_result" } ], "source": [ "a" ] }, { "cell_type": "code", "execution_count": 80, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "9.55535317195396e+16" ] }, "execution_count": 80, "metadata": {}, "output_type": "execute_result" } ], "source": [ "a//100000/60/60/24/365\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## OOP" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "A link to special method names in python: http://getpython3.com/diveintopython3/special-method-names.html" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Basic polynomial class\n", "\n", "We represent a polynomial as a list of coefficients, where the $i$'th element in the list is the coefficient of $x^i$" ] }, { "cell_type": "code", "execution_count": 81, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(1*x^0) + (1*x^1)" ] }, "execution_count": 81, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "(2*x^0) + (3*x^1) + (4*x^2)" ] }, "execution_count": 81, "metadata": {}, "output_type": "execute_result" } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " \n", " \n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", "f = Polynomial([1, 1])\n", "g = Polynomial([2, 3, 4])\n", "\n", "f\n", "g" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Add degree and evaluate functions (these are special to class polynomial and are not part of Python's recognized functions)" ] }, { "cell_type": "code", "execution_count": 82, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "1" ] }, "execution_count": 82, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "6" ] }, "execution_count": 82, "metadata": {}, "output_type": "execute_result" } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", " def degree(self):\n", " return len(self.coeffs) - 1\n", "\n", " def evaluate(self, x):\n", " result = 0\n", " for i, c in enumerate(self.coeffs):\n", " result += c * pow(x, i)\n", " return result\n", "\n", "f = Polynomial([1, 1])\n", "f.degree()\n", "f.evaluate(5)" ] }, { "cell_type": "code", "execution_count": 83, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "False" ] }, "execution_count": 83, "metadata": {}, "output_type": "execute_result" } ], "source": [ "p = Polynomial([1,0,1])\n", "q = Polynomial([1,0,1])\n", "p == q" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Add \\__eq\\__ function for the '==' operator\n", "Without it, Python compares memory adderesses" ] }, { "cell_type": "code", "execution_count": 84, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "False" ] }, "execution_count": 84, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "True" ] }, "execution_count": 84, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "True" ] }, "execution_count": 84, "metadata": {}, "output_type": "execute_result" } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", " def degree(self):\n", " return len(self.coeffs) - 1\n", "\n", " def evaluate(self, x):\n", " result = 0\n", " for i, c in enumerate(self.coeffs):\n", " result += c * pow(x, i)\n", " return result\n", "\n", " def __eq__(self, other):\n", " assert isinstance(other, Polynomial) \n", " return self.coeffs == other.coeffs\n", "\n", "f = Polynomial([1, 1])\n", "f2 = Polynomial([1, 1])\n", "g = Polynomial([2, 3, 4])\n", "\n", "f == g\n", "f == f2\n", "f.__eq__(f2)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Add the \\__add\\__ function for the '+' operator" ] }, { "cell_type": "code", "execution_count": 85, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(3*x^0) + (4*x^1) + (4*x^2)" ] }, "execution_count": 85, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "(101*x^0) + (1*x^1)" ] }, "execution_count": 85, "metadata": {}, "output_type": "execute_result" }, { "ename": "TypeError", "evalue": "unsupported operand type(s) for +: 'int' and 'Polynomial'", "output_type": "error", "traceback": [ "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m", "\u001b[1;31mTypeError\u001b[0m Traceback (most recent call last)", "\u001b[1;32m\u001b[0m in \u001b[0;36m\u001b[1;34m()\u001b[0m\n\u001b[0;32m 46\u001b[0m \u001b[0mf\u001b[0m\u001b[1;33m+\u001b[0m\u001b[0mg\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m 47\u001b[0m \u001b[0mf\u001b[0m\u001b[1;33m+\u001b[0m\u001b[1;36m100\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m---> 48\u001b[1;33m \u001b[1;36m100\u001b[0m\u001b[1;33m+\u001b[0m\u001b[0mf\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m", "\u001b[1;31mTypeError\u001b[0m: unsupported operand type(s) for +: 'int' and 'Polynomial'" ] } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", " def degree(self):\n", " return len(self.coeffs) - 1\n", "\n", " def evaluate(self, x):\n", " result = 0\n", " for i, c in enumerate(self.coeffs):\n", " result += c * pow(x, i)\n", " return result\n", "\n", " def __eq__(self, other):\n", " assert isinstance(other, Polynomial) \n", " return self.coeffs == other.coeffs\n", "\n", " def __add__(self, other):\n", " assert isinstance(other, (Polynomial, int)) \n", " if isinstance(other, int):\n", " terms = [c for c in self.coeffs]\n", " terms[0] += other\n", " return Polynomial(terms)\n", " #else, other is a polynomial\n", " size = min(self.degree(), other.degree()) + 1\n", " terms = [self.coeffs[i] + other.coeffs[i] for i in range(size)]\n", " terms += self.coeffs[size : self.degree() + 1]\n", " terms += other.coeffs[size : other.degree() + 1]\n", " #remove leading zeros\n", " last = len(terms) - 1\n", " while last > 0 and terms[last] == 0:\n", " last -= 1\n", " return Polynomial(terms[:last+1])\n", "\n", "f = Polynomial([1, 1])\n", "g = Polynomial([2, 3, 4])\n", "\n", "f+g\n", "f+100\n", "100+f" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Add \\__radd\\__ function for right addition with int\n", "Since there is no support for addition with a Polynomial in the class int, python needs a definition of right addition for int in class Polynomial\n", "\n", "When executing the command $f + 3$, python will first search for an \\__add\\__ method in class int that can support a second operand of type polynomial (self = int, other = polynomial).\n", "If such method does not exist, python will then search for an \\__radd\\__ method in class polynomial that can support a second operand of type int (self= polynomial, other = int).\n", "If such method does not exist, an exception is thrown.\n", "\n", "Since + is commutative, then we set \\__radd\\__ to be identical to \\__add\\__\n" ] }, { "cell_type": "code", "execution_count": 86, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(2*x^0) + (1*x^1)" ] }, "execution_count": 86, "metadata": {}, "output_type": "execute_result" } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", " def degree(self):\n", " return len(self.coeffs) - 1\n", "\n", " def evaluate(self, x):\n", " result = 0\n", " for i, c in enumerate(self.coeffs):\n", " result += c * pow(x, i)\n", " return result\n", "\n", " def __eq__(self, other):\n", " assert isinstance(other, Polynomial) \n", " return self.coeffs == other.coeffs\n", "\n", " def __add__(self, other):\n", " assert isinstance(other, (Polynomial, int)) \n", " if isinstance(other, int):\n", " terms = [c for c in self.coeffs]\n", " terms[0] += other\n", " return Polynomial(terms)\n", " #else, other is a polynomial\n", " size = min(self.degree(), other.degree()) + 1\n", " terms = [self.coeffs[i] + other.coeffs[i] for i in range(size)]\n", " terms += self.coeffs[size : self.degree() + 1]\n", " terms += other.coeffs[size : other.degree() + 1]\n", " #remove leading zeros\n", " last = len(terms) - 1\n", " while last > 0 and terms[last] == 0:\n", " last -= 1\n", " return Polynomial(terms[:last+1])\n", "\n", " __radd__ = __add__ #to allow int+Polynomial\n", "\n", "f = Polynomial([1, 1])\n", "1 + f" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Add \\__neg\\__ function for unary minus (placing a '-' before a polynomial object)\n" ] }, { "cell_type": "code", "execution_count": 87, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(-1*x^0) + (-1*x^1)" ] }, "execution_count": 87, "metadata": {}, "output_type": "execute_result" } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", " def degree(self):\n", " return len(self.coeffs) - 1\n", "\n", " def evaluate(self, x):\n", " result = 0\n", " for i, c in enumerate(self.coeffs):\n", " result += c * pow(x, i)\n", " return result\n", "\n", " def __eq__(self, other):\n", " assert isinstance(other, Polynomial) \n", " return self.coeffs == other.coeffs\n", "\n", " def __add__(self, other):\n", " assert isinstance(other, (Polynomial, int)) \n", " if isinstance(other, int):\n", " terms = [c for c in self.coeffs]\n", " terms[0] += other\n", " return Polynomial(terms)\n", " #else, other is a polynomial\n", " size = min(self.degree(), other.degree()) + 1\n", " terms = [self.coeffs[i] + other.coeffs[i] for i in range(size)]\n", " terms += self.coeffs[size : self.degree() + 1]\n", " terms += other.coeffs[size : other.degree() + 1]\n", " #remove leading zeros\n", " last = len(terms) - 1\n", " while last > 0 and terms[last] == 0:\n", " last -= 1\n", " return Polynomial(terms[:last+1])\n", "\n", " __radd__ = __add__ #to allow int+Polynomial\n", " \n", " def __neg__(self):\n", " return Polynomial([-co for co in self.coeffs])\n", "\n", "f = Polynomial([1, 1])\n", "-f" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Add \\__sub\\__ function to subtract one poly form another or an int from a poly\n", "Using the fact we already have \\__neg\\__ and \\__add\\__ we define $f-g = f + (-g)$" ] }, { "cell_type": "code", "execution_count": 88, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(1*x^0) + (2*x^1) + (4*x^2)" ] }, "execution_count": 88, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "(1*x^0) + (3*x^1) + (4*x^2)" ] }, "execution_count": 88, "metadata": {}, "output_type": "execute_result" } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", " def degree(self):\n", " return len(self.coeffs) - 1\n", "\n", " def evaluate(self, x):\n", " result = 0\n", " for i, c in enumerate(self.coeffs):\n", " result += c * pow(x, i)\n", " return result\n", "\n", " def __eq__(self, other):\n", " assert isinstance(other, Polynomial) \n", " return self.coeffs == other.coeffs\n", "\n", " def __add__(self, other):\n", " assert isinstance(other, (Polynomial, int)) \n", " if isinstance(other, int):\n", " terms = [c for c in self.coeffs]\n", " terms[0] += other\n", " return Polynomial(terms)\n", " #else, other is a polynomial\n", " size = min(self.degree(), other.degree()) + 1\n", " terms = [self.coeffs[i] + other.coeffs[i] for i in range(size)]\n", " terms += self.coeffs[size : self.degree() + 1]\n", " terms += other.coeffs[size : other.degree() + 1]\n", " #remove leading zeros\n", " last = len(terms) - 1\n", " while last > 0 and terms[last] == 0:\n", " last -= 1\n", " return Polynomial(terms[:last+1])\n", "\n", " __radd__ = __add__ #to allow int+Polynomial\n", " \n", " def __neg__(self):\n", " return Polynomial([-co for co in self.coeffs])\n", " \n", " def __sub__(self, other):\n", " assert isinstance(other, (int, Polynomial)) \n", " return (self + (-other))\n", " \n", "f = Polynomial([1, 1])\n", "g = Polynomial([2, 3, 4])\n", "\n", "g - f\n", "g - 1" ] }, { "cell_type": "code", "execution_count": 89, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(2*x^0) + (3*x^1) + (4*x^2)" ] }, "execution_count": 89, "metadata": {}, "output_type": "execute_result" } ], "source": [ "g" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Add \\__rsub\\__ function for right subtraction\n", "Using the fact we already have \\__neg\\__ and \\__sub\\__ we define $f-g = -(g-f)$\n", "\n", "since oprator \\- is not commutative, then we cannot set \\__rsub\\__ to be identical to \\__sub\\__\n", "\n", "Note that the following implementation of \\__rsub\\__ is promlematic, because it will keep calling itself (infinitely)\n", "\n", " def __rsub__(self, other):\n", " return other-self #PROBLEMATIC!\n", "\n", "\n" ] }, { "cell_type": "code", "execution_count": 91, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(-1*x^1)" ] }, "execution_count": 91, "metadata": {}, "output_type": "execute_result" }, { "data": { "text/plain": [ "(-1*x^0) + (-3*x^1) + (-4*x^2)" ] }, "execution_count": 91, "metadata": {}, "output_type": "execute_result" } ], "source": [ "class Polynomial():\n", " def __init__(self, coeffs):\n", " self.coeffs = coeffs\n", "\n", " def __repr__(self):\n", " terms = [\" + (\" + str(self.coeffs[k]) + \"*x^\" + \\\n", " str(k) + \")\" \\\n", " for k in range(len(self.coeffs)) \\\n", " if self.coeffs[k] != 0]\n", " terms = \"\".join(terms)\n", " if terms == \"\":\n", " return \"0\"\n", " else:\n", " return terms[3:] #discard leftmost '+'\n", "\n", " def degree(self):\n", " return len(self.coeffs) - 1\n", "\n", " def evaluate(self, x):\n", " result = 0\n", " for i, c in enumerate(self.coeffs):\n", " result += c * pow(x, i)\n", " return result\n", "\n", " def __eq__(self, other):\n", " assert isinstance(other, Polynomial) \n", " return self.coeffs == other.coeffs\n", "\n", " def __add__(self, other):\n", " assert isinstance(other, (Polynomial, int)) \n", " if isinstance(other, int):\n", " terms = [c for c in self.coeffs]\n", " terms[0] += other\n", " return Polynomial(terms)\n", " #else, other is a polynomial\n", " size = min(self.degree(), other.degree()) + 1\n", " terms = [self.coeffs[i] + other.coeffs[i] for i in range(size)]\n", " terms += self.coeffs[size : self.degree() + 1]\n", " terms += other.coeffs[size : other.degree() + 1]\n", " #remove leading zeros\n", " last = len(terms) - 1\n", " while last > 0 and terms[last] == 0:\n", " last -= 1\n", " return Polynomial(terms[:last+1])\n", "\n", " __radd__ = __add__ #to allow int+Polynomial\n", "\n", " def __neg__(self):\n", " return Polynomial([-co for co in self.coeffs])\n", "\n", " def __sub__(self, other):\n", " assert isinstance(other, (int, Polynomial)) \n", " return (self + (-other))\n", "\n", " #__rsub__ = __sub__ #not what we want...\n", "\n", " def __rsub__(self, other):\n", " # other - self\n", " return -(self - other) #why not just other-self ?\n", "\n", "f = Polynomial([1, 1])\n", "g = Polynomial([2, 3, 4])\n", "\n", "1 - f\n", "1 - g" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "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.6.2" } }, "nbformat": 4, "nbformat_minor": 1 }