{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Weekly exercise 5: Improving Sieves of Eratosthenes algorithm\n",
    "\n",
    "Recall the following classic algorithm for finding all prime numbers\n",
    "up to a given upper boundary.\n",
    "\n",
    "Algorithm:\n",
    "1. initialize the list with all integers (1 to the given upper boundary)\n",
    "2. go through a list of divisors\n",
    "3. cross all numbers on the first list which are divisible by the considered divisor and not equal to it\n",
    "4. stop when all divisors (up to the given boundary) are considered, and so all divisible numbers are crossed out, leaving only primes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 1\n",
    "\n",
    "Start by implementing the following version of the algorithm, using the starter code below."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "def eratosthenes(upper, verbosity = @@@):\n",
    "    \"\"\"Return the list of primes up the given upper bound using the sieve of Eratosthenes algorithm\"\"\"\n",
    "    primes = list(range(1,@@@)) # all numbers between 1 and upper\n",
    "    divisor = @@@ # initial divisor\n",
    "    while True: # loop until break is invoked\n",
    "        divisor += 1 # next divisor\n",
    "        if divisor > upper: # checked all divisors up to upper\n",
    "            @@@\n",
    "        i=0 # check for divisibility starting from the first element of remaining list\n",
    "        while i < len(primes): # until the end of the list (which changes length inside the loop)\n",
    "            if primes[i] @@@ divisor and primes[i] @@@ divisor == 0:\n",
    "                # remove divisible, except the divisor itself\n",
    "                primes.remove(primes[i]) # remove, next is with the same index\n",
    "            else:\n",
    "                i += 1 # skip to go to next element\n",
    "        if @@@: print(\"divisor %d:\"%divisor,primes)\n",
    "    if @@@: print(\"Primes up to %d are:\"%upper,primes)\n",
    "    return @@@"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To make sure the function is running as it is supposed to, run the\n",
    "tests below and confirm that the output is as expected."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "eratosthenes(25,verbosity=True) #should print the steps of the algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "x = eratosthenes(5)\n",
    "print(x) #should print [2,3,5]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "x = eratosthenes(2500)\n",
    "print('Number of primes up to 2500 is ',len(x)) #should print \"Number of primes up to 2500 is 368\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 2\n",
    "\n",
    "Use the code snippet from lecture video 9 to collect data on the runtime of the *eratosthenes()* as\n",
    "a function of the upper bound.\n",
    "Compute 25 runtime data points using the following sequence for the upper bound: 5, 10, 15, etc,\n",
    "and using 1000 functional calls to compute the average.\n",
    "\n",
    "Make a conjecture on the computational complexity of this algorithm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "plt.rcParams['figure.figsize'] = [12, 8]\n",
    "\n",
    "N = @@@\n",
    "n,x = [0,] @@@ N, [0,] @@@ N  # lists to accumulate the data\n",
    "for i,k in [(i,@@@) for i in range(N)]:\n",
    "    n[i] = k\n",
    "    t = @@@\n",
    "    x[i] = @@@\n",
    "\n",
    "plt.plot(n,x)\n",
    "plt.xlabel('upper bound')\n",
    "plt.ylabel('average run time, sec')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 3\n",
    "\n",
    "Copy your code from above and make the following improvements to the algorithm:\n",
    "\n",
    "1. Force divisors only be prime numbers themselves  \n",
    "1. When checking for divisibility, start from the square of the divisor instead of checking each one (why can we do it this way?)  \n",
    "1. Complete the algorithm when the starting point in previous step is above the upper bound  \n",
    "\n",
    "\n",
    "Run the same tests and make sure it the results are identical to the first implementation of the algorithm."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "# Write your code here\n",
    "#"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "eratosthenes(25,verbosity=True) #should print the steps of the algorithm\n",
    "eratosthenes_better(25,verbosity=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "x = eratosthenes(5)\n",
    "print(x) #should print [2,3,5]\n",
    "x = eratosthenes_better(5)\n",
    "print(x) #should print [2,3,5]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "x = eratosthenes(2500)\n",
    "print('Number of primes up to 2500 is ',len(x)) #should print \"Number of primes up to 2500 is 368\"\n",
    "x = eratosthenes_better(2500)\n",
    "print('Number of primes up to 2500 is ',len(x)) #should print \"Number of primes up to 2500 is 368\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Task 4\n",
    "\n",
    "Repeat the run time analysis for the improved version of the algorithm, and plot the graph for the two implementations side by side."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "hide-output": false
   },
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "plt.rcParams['figure.figsize'] = [12, 8]\n",
    "\n",
    "N = @@@\n",
    "n,x,y = @@@  # lists to accumulate the data\n",
    "for i,k in [(i,@@@) for i in range(N)]:\n",
    "    n[i] = @@@\n",
    "    t = @@@\n",
    "    x[i] = @@@  # run time of original algorithm\n",
    "    t = @@@\n",
    "    y[i] = @@@  # run time of improved algorithm\n",
    "\n",
    "plt.plot(n,x,label='Original algorithm')\n",
    "plt.plot(n,y,label='Improved algorithm')\n",
    "plt.xlabel('upper bound')\n",
    "plt.ylabel('average run time, sec')\n",
    "plt.legend()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note, however, that 1 is not a prime!\n",
    "[https://www.youtube.com/watch?v=IQofiPqhJ_s](https://www.youtube.com/watch?v=IQofiPqhJ_s)"
   ]
  }
 ],
 "metadata": {
  "date": 1627475014.509599,
  "filename": "exercise05.rst",
  "kernelspec": {
   "display_name": "Python",
   "language": "python3",
   "name": "python3"
  },
  "title": "Weekly exercise 5: Improving Sieves of Eratosthenes algorithm"
 },
 "nbformat": 4,
 "nbformat_minor": 4
}