PolyaCount.ipynb

polyacount.ipynb

{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "colab": {
      "name": "PolyaCount.ipynb",
      "provenance": [],
      "collapsed_sections": [],
      "include_colab_link": true
    },
    "kernelspec": {
      "name": "python3",
      "display_name": "Python 3"
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "view-in-github",
        "colab_type": "text"
      },
      "source": [
        "<a href=\"https://colab.research.google.com/gist/jtb/8d73053c1a4aab3f92c31ffafd172899/polyacount.ipynb\" target=\"_parent\"><img src=\"https://colab.research.google.com/assets/colab-badge.svg\" alt=\"Open In Colab\"/></a>"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "c7QxAHoQ2L7t"
      },
      "source": [
        "from sympy import *\n",
        "from sympy.combinatorics.permutations import Permutation\n",
        "from collections import Iterable\n",
        "\n",
        "# Generates the cyclic group of order n\n",
        "# from the specified rotate operation.\n",
        "def cyclicGen(n, rotate, init = None):\n",
        "  if init:\n",
        "    element = init\n",
        "  else:\n",
        "    element = rotate**0\n",
        "  for i in range(n):\n",
        "    yield Permutation(element)\n",
        "    element = element*rotate\n",
        "\n",
        "# Generates the dihedral group of order 2*n\n",
        "# from the specified rotate and flip operations.\n",
        "def dihedralGen(n, rotate, flip, init = None):\n",
        "  yield from cyclicGen(n, rotate, init)\n",
        "  if init:\n",
        "    initFlipped = init*flip\n",
        "  else:\n",
        "    initFlipped = flip\n",
        "  yield from cyclicGen(n, rotate, initFlipped)\n",
        "  \n",
        "# Generates an arbitrary group from a\n",
        "# set of generators\n",
        "def groupGen(*generators):\n",
        "  # ensure non empty\n",
        "  if not generators: raise ValueError(\"Empty list of generators\")\n",
        "  n = generators[0].size;\n",
        "  # ensure all permutations are same size\n",
        "  for item in generators:\n",
        "    if item.size != n:\n",
        "      raise ValueError(\"Size of permutations don't agree\")\n",
        "  \n",
        "  # Create a set of elements in the group, starting with the identity.\n",
        "  newElement = Permutation(size=n)\n",
        "  results = { newElement }\n",
        "  yield newElement\n",
        "  while True:\n",
        "    # create the next generation of elements\n",
        "    nextGen = set()\n",
        "    for i in results:\n",
        "      for j in generators:\n",
        "        p = i*j\n",
        "        if(p not in results and p not in nextGen):\n",
        "          nextGen.add(p)\n",
        "          yield p\n",
        "    if not nextGen: break\n",
        "    results = results.union(nextGen)\n",
        "    \n",
        "def permutation2monomial(groupElement, vars):\n",
        "  group_structure = groupElement.cycle_structure\n",
        "  ans = 1\n",
        "  for cycle_size, card in group_structure.items():\n",
        "    if isinstance(vars, str):\n",
        "      ans = ans*Symbol(vars + str(cycle_size))**card\n",
        "    elif isinstance(vars, Iterable):\n",
        "      ans = ans*sum([v**cycle_size for v in vars])**card\n",
        "    else: # number or symbol\n",
        "      ans = ans*(vars**card)\n",
        "  return simplify(ans)\n",
        "\n",
        "def polya(group, vars):\n",
        "  n = 0\n",
        "  ans = 0\n",
        "  for g in group:\n",
        "    n = n + 1\n",
        "    ans = ans + permutation2monomial(g, vars)\n",
        "  if n == 0:\n",
        "    raise ValueError(\"Group is empty. Perhaps you already used the generator.\")\n",
        "  return expand(simplify(ans / n))\n",
        "\n",
        "def containsVar(f, x):\n",
        "  return simplify(f - f.subs(x, 0))"
      ],
      "execution_count": 1,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "FWeLpslo38C6",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "f8376eda-c016-4b21-ca0c-7351a6e2f9c1"
      },
      "source": [
        "# Cyclic Group Examples \n",
        "\n",
        "# List elements of C3, starting at different elements\n",
        "rotate = Permutation([1, 2, 0])\n",
        "print([element for element in cyclicGen(3, rotate)])\n",
        "print([element for element in cyclicGen(3, rotate, rotate)])\n",
        "print([element for element in cyclicGen(3, rotate, rotate*rotate)])\n",
        "\n",
        "# Example: Rotational symmetries of verticies\n",
        "# of a square in the Euclidean plane.\n",
        "rotate = Permutation([1, 2, 3, 0])\n",
        "C_4 = cyclicGen(4, rotate)\n",
        "f = polya(C_4, \"a\")\n",
        "print(f)\n",
        "\n",
        "# Example: Rotation symmetries of sides and diagonals\n",
        "# of a square in the Euclidean plane\n",
        "rotate = Permutation([1, 2, 3, 0, 5, 4])\n",
        "C_4 = cyclicGen(4, rotate)\n",
        "f = polya(C_4, \"a\")\n",
        "print(f)\n",
        "\n",
        "# Example: Rotational symmetries of ordered pairs\n",
        "# of a square in the Euclidean plane, including (v,v)\n",
        "rotate = Permutation([1, 2, 3, 0, 5, 6, 7, 4, 9, 10, 11, 8, 13, 14, 15, 12])\n",
        "C_4 = cyclicGen(4, rotate)\n",
        "f = polya(C_4, \"a\")\n",
        "print(f)"
      ],
      "execution_count": 2,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "[Permutation(2), Permutation(0, 1, 2), Permutation(0, 2, 1)]\n",
            "[Permutation(0, 1, 2), Permutation(0, 2, 1), Permutation(2)]\n",
            "[Permutation(0, 2, 1), Permutation(2), Permutation(0, 1, 2)]\n",
            "a1**4/4 + a2**2/4 + a4/2\n",
            "a1**6/4 + a1**2*a2**2/4 + a2*a4/2\n",
            "a1**16/4 + a2**8/4 + a4**4/2\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "9LdIWsWjAqRk",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "7ae4578b-da4e-4f7b-956c-4784a7244e70"
      },
      "source": [
        "# Dihedral Group Examples\n",
        "# http://mathworld.wolfram.com/DihedralGroup.html\n",
        "\n",
        "# Z(D_1)\n",
        "rotate = Permutation([0])\n",
        "flip = Permutation([0])\n",
        "D_1 = dihedralGen(1, rotate, flip)\n",
        "print(polya(D_1, 'x'))\n",
        "\n",
        "# Z(D_2)\n",
        "rotate = Permutation([1, 0])\n",
        "flip = Permutation([0, 1])\n",
        "D_2 = dihedralGen(2, rotate, flip)\n",
        "print(polya(D_2, 'x'))\n",
        "\n",
        "# Z(D_3)\n",
        "rotate = Permutation([1, 2, 0])\n",
        "flip = Permutation([0, 2, 1])\n",
        "D_3 = dihedralGen(3, rotate, flip)\n",
        "print(polya(D_3, 'x'))\n",
        "\n",
        "# Z(D_4)\n",
        "rotate = Permutation([1, 2, 3, 0])\n",
        "flip = Permutation([1, 0, 3, 2])\n",
        "D_4 = dihedralGen(4, rotate, flip)\n",
        "print(polya(D_4, 'x'))\n",
        "\n",
        "# Z(D_5)\n",
        "rotate = Permutation([1, 2, 3, 4, 0])\n",
        "flip = Permutation([0, 4, 3, 2, 1])\n",
        "D_5 = dihedralGen(5, rotate, flip)\n",
        "print(polya(D_5, 'x'))\n"
      ],
      "execution_count": 3,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "x1\n",
            "x1**2/2 + x2/2\n",
            "x1**3/6 + x1*x2/2 + x3/3\n",
            "x1**4/8 + x1**2*x2/4 + 3*x2**2/8 + x4/4\n",
            "x1**5/10 + x1*x2**2/2 + 2*x5/5\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "AQCJvjnNHBWB",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "22b45b33-bfd4-491e-8fbc-7fe64eb1807c"
      },
      "source": [
        "# Generic Group Examples\n",
        "\n",
        "# Print permutations in index notation to\n",
        "# make it easier to see permutations.\n",
        "# Save original state so we can return back to it\n",
        "# once we are done.\n",
        "print_cycle_bool = Permutation.print_cyclic\n",
        "Permutation.print_cyclic = False\n",
        "\n",
        "# Symmetric Group Example\n",
        "swap = Permutation([1, 0, 2])\n",
        "rotate = Permutation([1, 2, 0])\n",
        "print([element for element in groupGen(swap, rotate)])\n",
        "\n",
        "# Edge Permutation Group\n",
        "# https://en.wikipedia.org/wiki/Cycle_index\n",
        "swap = Permutation([2, 1, 0])\n",
        "rotate = Permutation([1, 2, 0])\n",
        "K_3 = groupGen(rotate, swap)\n",
        "#print(polya(K_3, 'a'))\n",
        "\n",
        "x = Symbol('x')\n",
        "y = Symbol('y')\n",
        "z = Symbol('z')\n",
        "print(polya(K_3, [x,y]))\n",
        "\n",
        "swap = Permutation([0, 4, 2, 5, 1, 3])\n",
        "rotate = Permutation([1, 2, 3, 0, 5, 4])\n",
        "K_4 = groupGen(rotate, swap)\n",
        "print(polya(K_4, 'a'))\n",
        "\n",
        "# Return print_cyclic to original state\n",
        "Permutation.print_cyclic = print_cycle_bool\n",
        "\n",
        "# Choose(3, 2)\n",
        "Z = groupGen(Permutation(2))\n",
        "print(polya(Z, [x,1]).coeff(x**2))\n",
        "\n",
        "# 2x2x2 2-generator\n",
        "# http://sporadic.stanford.edu/bump/match/morepolished.pdf\n",
        "#        8     9\n",
        "#      +----+----+\n",
        "#  7  / 2  /  3 /|\n",
        "#    +----+----+ |\n",
        "# 6 / 1  / 0  /|10\n",
        "#  +----+----+ | +\n",
        "#  |    |    |11/| 14\n",
        "#  | 5  | 4  |/|12\n",
        "#  +----+----+ |/  15\n",
        "#       |    |13\n",
        "#       | 17 |/ 16\n",
        "#       +----+\n",
        "\n",
        "# Denote by R a clockwise rotation of the right face,\n",
        "# and let U denote a clockwise rotation of the top face,\n",
        "# so that G = <R, U>.\n",
        "R = Permutation([9, 1, 2, 14, 3, 5, 6, 7, 8, 15, 12, 10, 13, 11, 16, 17, 4, 0])\n",
        "U = Permutation([1, 2, 3, 0, 6, 7, 8, 9, 10, 11, 4, 5, 12, 13, 14, 15, 16, 17])\n",
        "G = groupGen(R, U)\n",
        "# Count order of group.\n",
        "print(sum(1 for g in G))\n",
        "\n"
      ],
      "execution_count": 7,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "[Permutation([0, 1, 2]), Permutation([1, 0, 2]), Permutation([1, 2, 0]), Permutation([2, 1, 0]), Permutation([0, 2, 1]), Permutation([2, 0, 1])]\n",
            "x**3 + x**2*y + x*y**2 + y**3\n",
            "a1**6/24 + 3*a1**2*a2**2/8 + a2*a4/4 + a3**2/3\n",
            "3\n",
            "29160\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Cp5-Fp9wKncp",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 51
        },
        "outputId": "49f5cf09-8874-43e7-b143-d3c99ecdfb0a"
      },
      "source": [
        "# Tic-Tac-Toe\n",
        "x = Symbol('x')\n",
        "y = Symbol('y')\n",
        "z = Symbol('z')\n",
        "rotate = Permutation([2, 5, 8, 1, 4, 7, 0, 3, 6])\n",
        "flip = Permutation([2, 1, 0, 5, 4, 3, 8, 7, 6])\n",
        "ticTacToe = dihedralGen(4, rotate, flip)\n",
        "#print(polya(ticTacToe, [x, 1]).coeff(x**5))\n",
        "\n",
        "ticTacToe2 = groupGen(rotate, flip)\n",
        "print(polya(ticTacToe2, [x, 1]).coeff(x**5))\n",
        "\n",
        "print(polya(ticTacToe, [x, y, 1]).coeff(x*y))"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "23\n",
            "12\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "tZxLqt5PnqyK",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 55
        },
        "outputId": "2160cbb0-3edd-41bf-f91b-55b0c16fb495"
      },
      "source": [
        "# Example: 2-colorings of 5x5 grid\n",
        "x, y = symbols('x y')\n",
        "rotate = Permutation(\n",
        "    [4, 9, 14, 19, 24,\n",
        "     3, 8, 13, 18, 23,\n",
        "     2, 7, 12, 17, 22,\n",
        "     1, 6, 11, 16, 21,\n",
        "     0, 5, 10, 15, 20])\n",
        "C_4 = cyclicGen(4, rotate)\n",
        "f = polya(C_4, [x, y])\n",
        "print(f)"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "x**25 + 7*x**24*y + 78*x**23*y**2 + 578*x**22*y**3 + 3182*x**21*y**4 + 13302*x**20*y**5 + 44330*x**19*y**6 + 120230*x**18*y**7 + 270525*x**17*y**8 + 510875*x**16*y**9 + 817388*x**15*y**10 + 1114548*x**14*y**11 + 1300316*x**13*y**12 + 1300316*x**12*y**13 + 1114548*x**11*y**14 + 817388*x**10*y**15 + 510875*x**9*y**16 + 270525*x**8*y**17 + 120230*x**7*y**18 + 44330*x**6*y**19 + 13302*x**5*y**20 + 3182*x**4*y**21 + 578*x**3*y**22 + 78*x**2*y**23 + 7*x*y**24 + y**25\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "eKSU6w28jAdX",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 139
        },
        "outputId": "36b91f29-6a02-4336-d097-eadf0e295794"
      },
      "source": [
        "# Cube\n",
        "one, two, three, four, five, six = symbols('one, two, three, four, five, six')\n",
        "right = Permutation([1, 2, 3, 0, 4, 5])\n",
        "up = Permutation([4, 1, 5, 3, 2, 0])\n",
        "dice = groupGen(right, up)\n",
        "\n",
        "# Number of cubes with six different colors.\n",
        "f = polya(dice, [one, two, three, four, five, six])\n",
        "print(f.as_poly().nth(1,1,1,1,1,1))\n",
        "\n",
        "# Number of cubes with up to m colors.\n",
        "# https://www.jstor.org/stable/2688630\n",
        "# https://oeis.org/A047780/internal\n",
        "# https://math.stackexchange.com/questions/460735/coloring-a-cube-with-4-colors\n",
        "#print(polya(groupGen(right, up), [1]*6)) # number of ways to color with up to 6 colors\n",
        "m = Symbol('m')\n",
        "g = polya(groupGen(right, up), m)\n",
        "# Up to 6 colors\n",
        "print(g.subs(m, 6))\n",
        "# Up to 4 colors\n",
        "print(g.subs(m, 4))\n",
        "\n",
        "# Exactly 4 colors\n",
        "g = polya(groupGen(right, up), [one, two, three, four])\n",
        "h = containsVar(containsVar(containsVar(containsVar(g, one), two), three), four)\n",
        "print(h)\n",
        "print(h.subs([(one, 1), (two, 1), (three, 1), (four, 1)]))\n",
        "\n",
        "# Print cyclic index\n",
        "print(polya(groupGen(right, up), 'a'))"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "30\n",
            "2226\n",
            "240\n",
            "four*one*three*two*(5*four**2 + 8*four*one + 8*four*three + 8*four*two + 5*one**2 + 8*one*three + 8*one*two + 5*three**2 + 8*three*two + 5*two**2)\n",
            "68\n",
            "a1**6/24 + a1**2*a2**2/8 + a1**2*a4/4 + a2**3/4 + a3**2/3\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "KFJOe5cJosoK",
        "colab": {
          "base_uri": "https://localhost:8080/"
        },
        "outputId": "6a763538-ebf5-4f2f-fd92-98aff70ac504"
      },
      "source": [
        "# Example from math.stackexchange\n",
        "# https://math.stackexchange.com/questions/416683/number-of-edge-colorings-in-a-tetrahedron-with-three-colors-is-my-solution-corr\n",
        "\n",
        "# Improper edge coloring of tetrahedron with 3 colors\n",
        "right = Permutation([1, 2, 0, 4, 5, 3])\n",
        "down = Permutation([2, 5, 3, 0, 1, 4])\n",
        "tetraV = groupGen(right, down)\n",
        "print(polya(tetraV, [1,1,1]))"
      ],
      "execution_count": null,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "87\n"
          ],
          "name": "stdout"
        }
      ]
    }
  ]
}