{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "# Foundations of Computational Economics #12\n",
    "\n",
    "by Fedor Iskhakov, ANU\n",
    "\n",
    "<img src=\"_static/img/dag3logo.png\" style=\"width:256px;\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "fragment"
    }
   },
   "source": [
    "## Enumeration of discrete compositions\n",
    "\n",
    "<img src=\"_static/img/lab.png\" style=\"width:64px;\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "subslide"
    }
   },
   "source": [
    "<img src=\"_static/img/youtube.png\" style=\"width:65px;\">\n",
    "\n",
    "[https://youtu.be/eU2WRHBTFBw](https://youtu.be/eU2WRHBTFBw)\n",
    "\n",
    "Description: Combinatorial enumeration. Python generators."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### Generators in Python\n",
    "\n",
    "**yield** keyword is a special kind of **return** from a function\n",
    "\n",
    "- “pauses” the function execution  \n",
    "- waits fro the next iteration request from the caller  \n",
    "- saves its state for the next call  \n",
    "\n",
    "\n",
    "Use generators instead of lists when only one element of the list is needed at any time\n",
    "\n",
    "- **saves memory**: alternative is a list of objects to be generated  \n",
    "- overhead in keeping the state of the generator function  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "#### Iterator vs Iterable vs Generator\n",
    "\n",
    "Iterator:\n",
    "\n",
    "- object that can return its members one at a time  \n",
    "- has __iter__ and __next__ methods  \n",
    "- support StopIteration error  \n",
    "- can be used in for loops  \n",
    "\n",
    "\n",
    "Iterable:\n",
    "\n",
    "- object that can be converted to an iterator with iter() function  \n",
    "- can be used in for loops directly  \n",
    "- examples: lists, tuples, strings  \n",
    "\n",
    "\n",
    "Generator:\n",
    "\n",
    "- a particular kind of iterator  \n",
    "- can be implemented as a function with yield returns  \n",
    "- or as *comprehension expression* similar to list comprehension  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "#### Note on range()\n",
    "\n",
    "Range object is an itarable, therefore\n",
    "\n",
    "- iter() function converts it to iterator  \n",
    "- it can be used directly in for loops  \n",
    "- but it is more than a generator: can be used in many other tasks besides looping  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "hide-output": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "initial call returned 1\n",
      "loop call returned 5\n",
      "loop call returned 15\n",
      "loop call returned 25\n",
      "loop finished gracefully\n"
     ]
    },
    {
     "ename": "StopIteration",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mStopIteration\u001b[0m                             Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-4-cb7c6099be2a>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[1;32m     13\u001b[0m     \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'loop call returned %d'\u001b[0m \u001b[0;34m%\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0mi\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m     14\u001b[0m \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'loop finished gracefully'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m---> 15\u001b[0;31m \u001b[0mnext\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mg\u001b[0m\u001b[0;34m)\u001b[0m  \u001b[0;31m# impossible to get any output once generator is done\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[0;31mStopIteration\u001b[0m: "
     ]
    }
   ],
   "source": [
    "# simple examples of generators\n",
    "def test_generator():\n",
    "    '''Testing generator\n",
    "    '''\n",
    "    yield 1\n",
    "    yield 5\n",
    "    yield 15\n",
    "    yield 25\n",
    "\n",
    "g = test_generator()\n",
    "print('initial call returned %d' % (next(g)))\n",
    "for i in g:\n",
    "    print('loop call returned %d' % (i))\n",
    "print('loop finished gracefully')\n",
    "next(g)  # impossible to get any output once generator is done"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "hide-output": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "generator at step 0\n",
      "main code at step 0\n",
      "generator at step 1\n",
      "main code at step 1\n",
      "generator at step 2\n",
      "main code at step 2\n",
      "generator at step 3\n",
      "main code at step 3\n",
      "generator at step 4\n",
      "main code at step 4\n",
      "generator at step 5\n",
      "main code at step 5\n",
      "generator at step 6\n",
      "main code at step 6\n",
      "generator at step 7\n",
      "main code at step 7\n",
      "generator at step 8\n",
      "main code at step 8\n",
      "generator at step 9\n",
      "main code at step 9\n"
     ]
    }
   ],
   "source": [
    "def test_generator_steps():\n",
    "    '''Testing generator\n",
    "    '''\n",
    "    for i in range(10):\n",
    "        print('generator at step %d'%(i))\n",
    "        yield i\n",
    "\n",
    "g = test_generator_steps()\n",
    "for i in g:\n",
    "    print('main code at step %d' % (i))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### Discrete compositions\n",
    "\n",
    "$$\n",
    "(p_1, p_2, \\dots, p_n) \\text{ such that } p_i \\in \\mathbb{Z}, 0 \\le p_i \\le M, \\sum_{i=1}^{n}p_i = M\n",
    "$$\n",
    "\n",
    "$ (p_1,p_2,\\dots,p_n) $ is **composition** of number $ M $ into\n",
    "$ n $ parts."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "#### Number of compositions\n",
    "\n",
    "- composition corresponds to cutting the interval of length $ M $ into $ n $ parts  \n",
    "- but cut points can be at the same place to allow for zeros  \n",
    "- instead let interval of length $ M+n $ be cut  \n",
    "- discrete $ \\rightarrow $ have $ M+n-1 $ points for $ n-1 $ cuts  \n",
    "- no overlaps in the latter scheme  \n",
    "\n",
    "\n",
    "$$\n",
    "\\text{Number of compositions} = {M+n-1 \\choose n-1} = \\frac{(M+n-1)!}{n!(M-1)!}\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "hide-output": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n=  2, M=  2 --> NC=3\n",
      "n=  2, M= 10 --> NC=11\n",
      "n=  2, M=100 --> NC=101\n",
      "n=  5, M= 10 --> NC=1001\n",
      "n=  5, M=100 --> NC=4598126\n",
      "n= 10, M=100 --> NC=4263421511270\n",
      "n= 50, M=100 --> NC=6709553636577312758557068362648666505216\n"
     ]
    }
   ],
   "source": [
    "import scipy.special\n",
    "def number_compositions(n,M):\n",
    "    '''Returns the number of discrete compositions for given parameters\n",
    "    '''\n",
    "    return int(scipy.special.comb(M+n-1,n-1)) # M+n-1 choose n-1\n",
    "\n",
    "print('n=%3d, M=%3d --> NC=%d'%(2,2,number_compositions(2,2)))\n",
    "print('n=%3d, M=%3d --> NC=%d'%(2,10,number_compositions(2,10)))\n",
    "print('n=%3d, M=%3d --> NC=%d'%(2,100,number_compositions(2,100)))\n",
    "print('n=%3d, M=%3d --> NC=%d'%(5,10,number_compositions(5,10)))\n",
    "print('n=%3d, M=%3d --> NC=%d'%(5,100,number_compositions(5,100)))\n",
    "print('n=%3d, M=%3d --> NC=%d'%(10,100,number_compositions(10,100)))\n",
    "print('n=%3d, M=%3d --> NC=%d'%(50,100,number_compositions(50,100)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "#### Lexicographical order\n",
    "\n",
    "Composition $ p = (p_1, p_2, \\dots, p_n) $ is greater than\n",
    "composition $ p' = (p'_1, p'_2, \\dots, p'_n) $ in lexicographical sense\n",
    "iff for some $ J \\in \\{1,\\dots,n\\} $ $ p_j=p'_j $ for all $ j<J $\n",
    "and $ p_J>p'_J $.\n",
    "\n",
    "$ j=n $ is the *lowest digit*, $ j=1 $ is the *highest digit*\n",
    "\n",
    "Composition $ p $ is greater that $ p' $ iff the *highest different digit*\n",
    "of $ p $ is greater than that of $ p' $."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "#### Examples of lexicographical order\n",
    "\n",
    "- numbers in any base system  \n",
    "- words in a dictionary  \n",
    "- compositions to be generated  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "hide-output": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "def compositions(N,m):\n",
    "    '''Iterable on compositions of N with m parts\n",
    "    Returns the generator (to be used in for loops)\n",
    "    '''\n",
    "    pass"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "hide-output": false,
    "scrolled": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  0: [0, 0, 0, 8]\n",
      "  1: [0, 0, 1, 7]\n",
      "  2: [0, 0, 2, 6]\n",
      "  3: [0, 0, 3, 5]\n",
      "  4: [0, 0, 4, 4]\n",
      "  5: [0, 0, 5, 3]\n",
      "  6: [0, 0, 6, 2]\n",
      "  7: [0, 0, 7, 1]\n",
      "  8: [0, 0, 8, 0]\n",
      "  9: [0, 1, 0, 7]\n",
      " 10: [0, 1, 1, 6]\n",
      " 11: [0, 1, 2, 5]\n",
      " 12: [0, 1, 3, 4]\n",
      " 13: [0, 1, 4, 3]\n",
      " 14: [0, 1, 5, 2]\n",
      " 15: [0, 1, 6, 1]\n",
      " 16: [0, 1, 7, 0]\n",
      " 17: [0, 2, 0, 6]\n",
      " 18: [0, 2, 1, 5]\n",
      " 19: [0, 2, 2, 4]\n",
      " 20: [0, 2, 3, 3]\n",
      " 21: [0, 2, 4, 2]\n",
      " 22: [0, 2, 5, 1]\n",
      " 23: [0, 2, 6, 0]\n",
      " 24: [0, 3, 0, 5]\n",
      " 25: [0, 3, 1, 4]\n",
      " 26: [0, 3, 2, 3]\n",
      " 27: [0, 3, 3, 2]\n",
      " 28: [0, 3, 4, 1]\n",
      " 29: [0, 3, 5, 0]\n",
      " 30: [0, 4, 0, 4]\n",
      " 31: [0, 4, 1, 3]\n",
      " 32: [0, 4, 2, 2]\n",
      " 33: [0, 4, 3, 1]\n",
      " 34: [0, 4, 4, 0]\n",
      " 35: [0, 5, 0, 3]\n",
      " 36: [0, 5, 1, 2]\n",
      " 37: [0, 5, 2, 1]\n",
      " 38: [0, 5, 3, 0]\n",
      " 39: [0, 6, 0, 2]\n",
      " 40: [0, 6, 1, 1]\n",
      " 41: [0, 6, 2, 0]\n",
      " 42: [0, 7, 0, 1]\n",
      " 43: [0, 7, 1, 0]\n",
      " 44: [0, 8, 0, 0]\n",
      " 45: [1, 0, 0, 7]\n",
      " 46: [1, 0, 1, 6]\n",
      " 47: [1, 0, 2, 5]\n",
      " 48: [1, 0, 3, 4]\n",
      " 49: [1, 0, 4, 3]\n",
      " 50: [1, 0, 5, 2]\n",
      " 51: [1, 0, 6, 1]\n",
      " 52: [1, 0, 7, 0]\n",
      " 53: [1, 1, 0, 6]\n",
      " 54: [1, 1, 1, 5]\n",
      " 55: [1, 1, 2, 4]\n",
      " 56: [1, 1, 3, 3]\n",
      " 57: [1, 1, 4, 2]\n",
      " 58: [1, 1, 5, 1]\n",
      " 59: [1, 1, 6, 0]\n",
      " 60: [1, 2, 0, 5]\n",
      " 61: [1, 2, 1, 4]\n",
      " 62: [1, 2, 2, 3]\n",
      " 63: [1, 2, 3, 2]\n",
      " 64: [1, 2, 4, 1]\n",
      " 65: [1, 2, 5, 0]\n",
      " 66: [1, 3, 0, 4]\n",
      " 67: [1, 3, 1, 3]\n",
      " 68: [1, 3, 2, 2]\n",
      " 69: [1, 3, 3, 1]\n",
      " 70: [1, 3, 4, 0]\n",
      " 71: [1, 4, 0, 3]\n",
      " 72: [1, 4, 1, 2]\n",
      " 73: [1, 4, 2, 1]\n",
      " 74: [1, 4, 3, 0]\n",
      " 75: [1, 5, 0, 2]\n",
      " 76: [1, 5, 1, 1]\n",
      " 77: [1, 5, 2, 0]\n",
      " 78: [1, 6, 0, 1]\n",
      " 79: [1, 6, 1, 0]\n",
      " 80: [1, 7, 0, 0]\n",
      " 81: [2, 0, 0, 6]\n",
      " 82: [2, 0, 1, 5]\n",
      " 83: [2, 0, 2, 4]\n",
      " 84: [2, 0, 3, 3]\n",
      " 85: [2, 0, 4, 2]\n",
      " 86: [2, 0, 5, 1]\n",
      " 87: [2, 0, 6, 0]\n",
      " 88: [2, 1, 0, 5]\n",
      " 89: [2, 1, 1, 4]\n",
      " 90: [2, 1, 2, 3]\n",
      " 91: [2, 1, 3, 2]\n",
      " 92: [2, 1, 4, 1]\n",
      " 93: [2, 1, 5, 0]\n",
      " 94: [2, 2, 0, 4]\n",
      " 95: [2, 2, 1, 3]\n",
      " 96: [2, 2, 2, 2]\n",
      " 97: [2, 2, 3, 1]\n",
      " 98: [2, 2, 4, 0]\n",
      " 99: [2, 3, 0, 3]\n",
      "100: [2, 3, 1, 2]\n",
      "101: [2, 3, 2, 1]\n",
      "102: [2, 3, 3, 0]\n",
      "103: [2, 4, 0, 2]\n",
      "104: [2, 4, 1, 1]\n",
      "105: [2, 4, 2, 0]\n",
      "106: [2, 5, 0, 1]\n",
      "107: [2, 5, 1, 0]\n",
      "108: [2, 6, 0, 0]\n",
      "109: [3, 0, 0, 5]\n",
      "110: [3, 0, 1, 4]\n",
      "111: [3, 0, 2, 3]\n",
      "112: [3, 0, 3, 2]\n",
      "113: [3, 0, 4, 1]\n",
      "114: [3, 0, 5, 0]\n",
      "115: [3, 1, 0, 4]\n",
      "116: [3, 1, 1, 3]\n",
      "117: [3, 1, 2, 2]\n",
      "118: [3, 1, 3, 1]\n",
      "119: [3, 1, 4, 0]\n",
      "120: [3, 2, 0, 3]\n",
      "121: [3, 2, 1, 2]\n",
      "122: [3, 2, 2, 1]\n",
      "123: [3, 2, 3, 0]\n",
      "124: [3, 3, 0, 2]\n",
      "125: [3, 3, 1, 1]\n",
      "126: [3, 3, 2, 0]\n",
      "127: [3, 4, 0, 1]\n",
      "128: [3, 4, 1, 0]\n",
      "129: [3, 5, 0, 0]\n",
      "130: [4, 0, 0, 4]\n",
      "131: [4, 0, 1, 3]\n",
      "132: [4, 0, 2, 2]\n",
      "133: [4, 0, 3, 1]\n",
      "134: [4, 0, 4, 0]\n",
      "135: [4, 1, 0, 3]\n",
      "136: [4, 1, 1, 2]\n",
      "137: [4, 1, 2, 1]\n",
      "138: [4, 1, 3, 0]\n",
      "139: [4, 2, 0, 2]\n",
      "140: [4, 2, 1, 1]\n",
      "141: [4, 2, 2, 0]\n",
      "142: [4, 3, 0, 1]\n",
      "143: [4, 3, 1, 0]\n",
      "144: [4, 4, 0, 0]\n",
      "145: [5, 0, 0, 3]\n",
      "146: [5, 0, 1, 2]\n",
      "147: [5, 0, 2, 1]\n",
      "148: [5, 0, 3, 0]\n",
      "149: [5, 1, 0, 2]\n",
      "150: [5, 1, 1, 1]\n",
      "151: [5, 1, 2, 0]\n",
      "152: [5, 2, 0, 1]\n",
      "153: [5, 2, 1, 0]\n",
      "154: [5, 3, 0, 0]\n",
      "155: [6, 0, 0, 2]\n",
      "156: [6, 0, 1, 1]\n",
      "157: [6, 0, 2, 0]\n",
      "158: [6, 1, 0, 1]\n",
      "159: [6, 1, 1, 0]\n",
      "160: [6, 2, 0, 0]\n",
      "161: [7, 0, 0, 1]\n",
      "162: [7, 0, 1, 0]\n",
      "163: [7, 1, 0, 0]\n",
      "164: [8, 0, 0, 0]\n"
     ]
    }
   ],
   "source": [
    "n, M = 4, 8\n",
    "for i,k in enumerate(compositions(M,n)):\n",
    "    print('%3d'%i,end=\": \")\n",
    "    print(k)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "hide-output": false,
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "outputs": [],
   "source": [
    "def compositions(N,m):\n",
    "    '''Iterable on compositions of N with m parts\n",
    "    Returns the generator (to be used in for loops)\n",
    "    '''\n",
    "    cmp=[0,]*m\n",
    "    cmp[m-1]=N  # initial composition is all to the last\n",
    "    yield cmp\n",
    "    while cmp[0]!=N:\n",
    "        i=m-1\n",
    "        while cmp[i]==0: i-=1  # find lowest non-zero digit\n",
    "        cmp[i-1] = cmp[i-1]+1  # increment next digit\n",
    "        cmp[m-1] = cmp[i]-1    # the rest to the lowest\n",
    "        if i!=m-1: cmp[i] = 0  # maintain cost sum\n",
    "        yield cmp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "slideshow": {
     "slide_type": "slide"
    }
   },
   "source": [
    "### Further learning resources\n",
    "\n",
    "- Generators vs. lists in Python\n",
    "  [https://www.youtube.com/watch?v=bD05uGo_sVI](https://www.youtube.com/watch?v=bD05uGo_sVI)  \n",
    "- Iterators and generators in Trey Hunner blog\n",
    "  [https://treyhunner.com/2018/06/how-to-make-an-iterator-in-python/](https://treyhunner.com/2018/06/how-to-make-an-iterator-in-python/)  \n",
    "- Some examples of combinatorial optimization\n",
    "  [https://neos-guide.org/content/combinatorial-optimization](https://neos-guide.org/content/combinatorial-optimization)  "
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Slideshow",
  "date": 1612589584.68963,
  "download_nb": false,
  "filename": "12_alg_combinatorics.rst",
  "filename_with_path": "12_alg_combinatorics",
  "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.7.6"
  },
  "title": "Foundations of Computational Economics #12"
 },
 "nbformat": 4,
 "nbformat_minor": 4
}