AG_sage.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "cf. David Cox, John Little, Donal O’Shea. **Ideals, Varieties, and Algorithms: An Introduction to Computational Algebraic Geometry and Commutative Algebra**, Fourth Edition, Springer, 2015"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Same ways to construct Polynomial Rings"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "PolyQQ_3.<x,y,z>=PolynomialRing(QQ,3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "PolyQQ_3_b = PolynomialRing(QQ,'x,y,z')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "PolyQQ_3==PolyQQ_3_b"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Implementing \n",
    "\n",
    "$$\n",
    "f=2x^3y^2z + \\frac{3}{2}y^3z^3 - 3xyz + y^2\n",
    "$$  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 246,
   "metadata": {},
   "outputs": [],
   "source": [
    "f = 2*x**3*y**2*z+3/2*y**3*z**3-3*x*y*z+y**2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 247,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(x, y, z)\n",
      "6\n",
      "(3, 3, 3)\n",
      "{(3, 2, 1): 2, (1, 1, 1): -3, (0, 2, 0): 1, (0, 3, 3): 3/2}\n",
      "[(3, 2, 1), (0, 3, 3), (1, 1, 1), (0, 2, 0)]\n",
      "[x^3*y^2*z, y^3*z^3, x*y*z, y^2]\n",
      "(x, y, z)\n"
     ]
    }
   ],
   "source": [
    "# many (rich amount!) Python methods to try out\n",
    "print( f.args())\n",
    "print( f.degree());print(f.degrees());print(f.dict());print(f.exponents());print(f.monomials());print(f.variables())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 249,
   "metadata": {},
   "outputs": [
    {
     "ename": "NotImplementedError",
     "evalue": "",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mNotImplementedError\u001b[0m                       Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-249-f289d34bb7b8>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m()\u001b[0m\n\u001b[1;32m----> 1\u001b[1;33m \u001b[0mPolyQQ_3\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0morder\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m/home/topolo/Public/sage/src/sage/rings/ring.pyx\u001b[0m in \u001b[0;36msage.rings.ring.Ring.order (/home/topolo/Public/sage/src/build/cythonized/sage/rings/ring.c:9125)\u001b[1;34m()\u001b[0m\n\u001b[0;32m   1033\u001b[0m         \u001b[1;32mif\u001b[0m \u001b[0mself\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mis_zero\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m   1034\u001b[0m             \u001b[1;32mreturn\u001b[0m \u001b[1;36m1\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m-> 1035\u001b[1;33m         \u001b[1;32mraise\u001b[0m \u001b[0mNotImplementedError\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m   1036\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m   1037\u001b[0m     \u001b[1;32mdef\u001b[0m \u001b[0mzeta\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mn\u001b[0m\u001b[1;33m=\u001b[0m\u001b[1;36m2\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mall\u001b[0m\u001b[1;33m=\u001b[0m\u001b[0mFalse\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mNotImplementedError\u001b[0m: "
     ]
    }
   ],
   "source": [
    "PolyQQ_3.order()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Properties of Grobner Bases"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 324,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Poly_RR__2.<x,y>=PolynomialRing(RR,2,order=\"lex\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 327,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(x^3 + x*y, x*y^3)"
      ]
     },
     "execution_count": 327,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(x**5*y).quo_rem( x**2*y-y**2 )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Buchberger's Algorithm\n",
    "cf. Sec. 7 of Cox, Little, O'Shea (2015), pp. 90"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Example 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 311,
   "metadata": {},
   "outputs": [],
   "source": [
    "Poly_QQ__2.<x,y> = PolynomialRing(QQ,2,order=\"deglex\")\n",
    "f1=x**3-2*x*y\n",
    "f2=x**2*y-2*y**2+x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 312,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-x^2"
      ]
     },
     "execution_count": 312,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "make_Spoly(f1,f2).lt()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 314,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 314,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "make_Spoly(f1,f2).lt() in Ideal(f1.lt(),f2.lt())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 321,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(0, -x^2)\n",
      "(0, -x^2)\n"
     ]
    }
   ],
   "source": [
    "print( make_Spoly(f1,f2).lt().quo_rem( f1.lt() ) )\n",
    "print( make_Spoly(f1,f2).lt().quo_rem( f2.lt() ) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 332,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "F=[f1,f2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 333,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[x^3 - 2*x*y, x^2*y - 2*y^2 + x, -x^2]\n"
     ]
    }
   ],
   "source": [
    "F.append( make_Spoly(f1,f2).lt().quo_rem(f1.lt())[1])\n",
    "print(F)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 334,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from itertools import combinations\n",
    "M = [_ for _ in combinations(F,2)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 336,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-x^2"
      ]
     },
     "execution_count": 336,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "make_Spoly(M[0][0],M[0][1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 337,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-2*x*y"
      ]
     },
     "execution_count": 337,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "make_Spoly(M[1][0],M[1][1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# This is the original implementation from James Carlson, slightly modified\n",
    "# cf. https://www.math.utah.edu/~carlson/cimat/lecture1.pdf\n",
    "def div(f,g, polyring):\n",
    "    zeroelement = polyring(0)\n",
    "    n = len(g)\n",
    "    p, r, a = f,0, [0 for x in range(0,n)]\n",
    "    while p != zeroelement:\n",
    "        i, divisionoccured = 0, False\n",
    "        while i < n and divisionoccured == False:\n",
    "            fp_remainder = (p.lt()).quo_rem(g[i].lt())[0]\n",
    "            if fp_remainder == zeroelement:\n",
    "                a[i] += p.lt()//g[i].lt()\n",
    "                p -= p.lt()//g[i].lt()*g[i]\n",
    "                divisionoccured == True\n",
    "            else:\n",
    "                i += 1\n",
    "        if divisionoccured == False:\n",
    "            r += p.lt()\n",
    "            p -= p.lt()\n",
    "    return a, r"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 377,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# I understand my version here slightly better because I stepped through each step with an explicit example below; \n",
    "# I'll copy this algorithm here:\n",
    "def DivisionAlgo( p, fs, ring):\n",
    "    \"\"\" DivisionAlgo - Division algorithm\n",
    "    \"\"\"\n",
    "    s=len(fs)\n",
    "    zeroelement = ring(0)    \n",
    "    a=[zeroelement for i in range(s)]\n",
    "    r=zeroelement\n",
    "    while p != zeroelement:\n",
    "        division_FLAG=False\n",
    "        for i in range(s):\n",
    "            p_dividedby_fs = p.lm().quo_rem(fs[i].lm())\n",
    "            if p_dividedby_fs[1] == zeroelement:\n",
    "                a[i] += p.lt().quo_rem(fs[i].lt())[0]  \n",
    "#                p -= p.lt().quo_rem(fs[i].lt())[0]*fs[i]\n",
    "                p -= p.lt()//fs[i].lt()*fs[i]\n",
    "                division_FLAG=True\n",
    "                break\n",
    "        if division_FLAG is False:\n",
    "            r += p.lt()\n",
    "            p -= p.lt()\n",
    "    return r,a"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is *exactly*  \n",
    "$$  \n",
    "\\overline{S(f_1,f_3)^F}  =-2xy \\neq 0\n",
    "$$  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 341,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-2*x*y, [0, 0, 0])"
      ]
     },
     "execution_count": 341,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo( make_Spoly(M[1][0],M[1][1]), F,Poly_QQ__2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 342,
   "metadata": {},
   "outputs": [],
   "source": [
    "r,a = DivisionAlgo( make_Spoly(M[1][0],M[1][1]), F,Poly_QQ__2)\n",
    "F.append(r)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 350,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[x^3 - 2*x*y, x^2*y - 2*y^2 + x, -x^2, -2*x*y]\n",
      "(0, [0, 0, 1, 0])\n",
      "(0, [0, 0, 0, 1])\n",
      "(0, [0, 0, 0, y])\n",
      "(-2*y^2 + x, [0, 0, 0, 0])\n"
     ]
    }
   ],
   "source": [
    "print(F)\n",
    "M = [_ for _ in combinations(F,2)]\n",
    "print(DivisionAlgo(make_Spoly(M[0][0],M[0][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[1][0],M[1][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[2][0],M[2][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[3][0],M[3][1]), F, Poly_QQ__2) )\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\"Hence, we must also add $f_5 = -2y^2+x$ to our generating set.\"  pp. 91, Sec.7 Buchberger's Algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 351,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-2*y^2 + x\n"
     ]
    }
   ],
   "source": [
    "print( make_Spoly(M[3][0],M[3][1]) )\n",
    "F.append( make_Spoly(M[3][0],M[3][1]) )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 378,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[x^3 - 2*x*y, x^2*y - 2*y^2 + x, -x^2, -2*x*y, -2*y^2 + x]\n",
      "(0, [0, 0, 1, 0, 0])\n",
      "(0, [0, 0, 0, 1, 0])\n",
      "(0, [0, 0, 0, y, 0])\n",
      "(0, [1/2*x, 1, 0, y^2, -1])\n",
      "(0, [0, 0, 0, 0, 1])\n",
      "(0, [0, 0, 0, 0, 1])\n",
      "(0, [1/2, 0, 0, -1/2, y])\n",
      "(0, [0, 0, 0, 0, 0])\n",
      "(0, [1/2, 0, 0, -1/2, 0])\n",
      "(0, [0, 0, -1/2, 0, 0])\n"
     ]
    }
   ],
   "source": [
    "print(F)\n",
    "M = [_ for _ in combinations(F,2)]\n",
    "print(DivisionAlgo(make_Spoly(M[0][0],M[0][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[1][0],M[1][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[2][0],M[2][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[3][0],M[3][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[4][0],M[4][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[5][0],M[5][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[6][0],M[6][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[7][0],M[7][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[8][0],M[8][1]), F, Poly_QQ__2) )\n",
    "print(DivisionAlgo(make_Spoly(M[9][0],M[9][1]), F, Poly_QQ__2) )\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "EY:20170311 In the next 2 lines, I had to debug my DivisionAlgo, Division Algorithm implementation.  I wanted to keep this to show that you can use the very useful Python interactive interpreter in jupyter to effectively debug your program/function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 366,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(x, [1, 0, 0, 0, 0])\n",
      "(x, [1, 0, 0, 0, 0])\n",
      "-2*y^2 + x\n",
      "-2*y^2 + x\n"
     ]
    }
   ],
   "source": [
    "print(DivisionAlgo( make_Spoly( M[4][0],M[4][1] ), [_ for _ in reversed(F)],Poly_QQ__2) )\n",
    "print(DivisionAlgo( make_Spoly( M[5][0],M[5][1] ), [_ for _ in reversed(F)],Poly_QQ__2) )\n",
    "print(make_Spoly( M[4][0],M[4][1]))\n",
    "print(make_Spoly( M[5][0],M[5][1]))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 376,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-2*y^2 + x\n",
      "(0, y^2)\n",
      "False\n",
      "x\n",
      "(0, x)\n",
      "False\n",
      "0\n",
      "(0, 0)\n",
      "True\n",
      "0\n",
      "0\n",
      "(0, 0)\n",
      "True\n"
     ]
    }
   ],
   "source": [
    "p=make_Spoly(M[4][0],M[4][1])\n",
    "r=Poly_QQ__2(0)\n",
    "print(p)\n",
    "print( p.lm().quo_rem(F[0].lm()) ); print( p.lm().quo_rem(F[0].lm())[1] == Poly_QQ__2(0) );r+=p.lt(); p -= p.lt();\n",
    "print(p);print( p.lm().quo_rem(F[1].lm()) ); print( p.lm().quo_rem(F[1].lm())[1] == Poly_QQ__2(0) );r+=p.lt(); p-=p.lt()\n",
    "print(p);print( p.lm().quo_rem(F[2].lm()) ); print( p.lm().quo_rem(F[2].lm())[1] == Poly_QQ__2(0) );\n",
    "print(p.lt()//F[2].lt()*F[2])\n",
    "p-=p.lt()//F[2].lt()*F[2]\n",
    "# r+=p.lt(); p-=p.lt()\n",
    "print(p);print( p.lm().quo_rem(F[3].lm()) ); print( p.lm().quo_rem(F[3].lm())[1] == Poly_QQ__2(0) );r+=p.lt(); p-=p.lt()\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Theorem 2 (Buchberger's Algorithm) (cf. pp. 91 of Cox, Little, O'Shea (2015))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 379,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def BuchbergerAlgo(F,polyring):\n",
    "    \"\"\" BuchbergerAlgo : Buchberger's Algorithm\n",
    "    INPUT\n",
    "    =====\n",
    "    @type F : Python list of Sage Math polynomials\n",
    "    @param F : F = (f_1, \\dots, f_s)\n",
    "    \n",
    "    @type polyring : Sage Math PolynomialRing\n",
    "    \"\"\"\n",
    "    zeroelement = polyring(0)\n",
    "    pairs=[_ for _ in combinations(F,2)]\n",
    "#    NO_PAIRS = len(pairs)\n",
    "    G=F\n",
    "    G_2 = [] # denoted G' or Gprime in Cox, Little, O'Shea (2015)\n",
    "\n",
    "    # first case guaranteed to run \n",
    "    while G_2 != G:\n",
    "        G_2 = G\n",
    "        for pair in pairs:\n",
    "            p=pair[0]\n",
    "            q=pair[1]\n",
    "            r,a_list = DivisionAlgo(make_Spoly(p,q), G, polyring)\n",
    "            if r != zeroelement:\n",
    "                G.append(r)\n",
    "        # END of while, break condition being when we don't add,append, anymore to G\n",
    "    return G"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 381,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^3 - 2*x*y, x^2*y - 2*y^2 + x, -x^2, -2*x*y, -2*y^2 + x]"
      ]
     },
     "execution_count": 381,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "BuchbergerAlgo(F,Poly_QQ__2)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# From Martin R. Albrecht *DTU Crypto Group*, **Gröbner Bases**\n",
    "cf. https://martinralbrecht.files.wordpress.com/2010/07/20131022_buchberger_dtu.pdf [Grobner Bases, Albrecht](https://martinralbrecht.files.wordpress.com/2010/07/20131022_buchberger_dtu.pdf)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(7))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "f=x*y+3*y*z-3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "g=x^2-2*y^2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "x^3*y - 2*x*y^3 + 3*x^2*y*z + y^3*z - 3*x^2 - y^2"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f*g"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "x^2 + x*y - 2*y^2 + 3*y*z - 3"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f+g"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y>=PolynomialRing(FiniteField(7))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "f=2*x*y+3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x*y, 1]"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.monomials()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[2, 3]"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.coefficients()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Monomial Orderings"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(7),order=\"lex\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x^2*y*z^2 > x*y^3*z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(7), order=\"deglex\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 14,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x^2*y*z^2 > x*y^3*z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(7), order=\"degrevlex\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x^2*y*z^2 > x*y^3*z"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercises"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Poly_F7 = PolynomialRing(FiniteField(7),'x,y,z',order='deglex')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2*x^2*y^8 - 3*x^5*y*z^4 - x*y^4 + x*y*z^3 - 3*y*z + 3*x"
      ]
     },
     "execution_count": 41,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f=4*x**5*y*z**4 + 2*x**2*y**8 - x*y**4 + x*y*z**3+4*y*z+3*x\n",
    "LatexExpr(f)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 42,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2*x^2*y^8 - 3*x^5*y*z^4 - x*y^4 + x*y*z^3 - 3*y*z + 3*x"
      ]
     },
     "execution_count": 42,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Poly_F7 = PolynomialRing(FiniteField(7),'x,y,z',order='lex')\n",
    "f=4*x**5*y*z**4 + 2*x**2*y**8 - x*y**4 + x*y*z**3+4*y*z+3*x\n",
    "LatexExpr(f)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 43,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2*x^2*y^8 - 3*x^5*y*z^4 - x*y^4 + x*y*z^3 - 3*y*z + 3*x"
      ]
     },
     "execution_count": 43,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Poly_F7 = PolynomialRing(FiniteField(7),'x,y,z',order='degrevlex')\n",
    "f=4*x**5*y*z**4 + 2*x**2*y**8 - x*y**4 + x*y*z**3+4*y*z+3*x\n",
    "LatexExpr(f)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Polynomial Division\n",
    "### Polynomial Division II \n",
    "#### Algorithm 1: Polynomial Division"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 69,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Poly_F7__2 = PolynomialRing(FiniteField(7),'x,y',order='deglex')\n",
    "f = 3*x**2*y + 2*y**2+x+1\n",
    "f1 = 2*x**2+1\n",
    "f2=y+1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 70,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x^2*y\n",
      "x^2\n",
      "y\n"
     ]
    }
   ],
   "source": [
    "print(p.lm());print(f_1.lm());print(f_2.lm())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "y"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "p.lm()/f1.lm()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x^2*y + 2*y^2 + x - y + 1\n"
     ]
    }
   ],
   "source": [
    "print(p-p.lm()/f1.lm()*f1)\n",
    "#(p-p.lm()/f_1.lm()*f_1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 73,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f==(-2*y)*f1 + (2*y)*f2 + (x+1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I will try to step through the division algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 78,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "y\n",
      "-2*y\n",
      "2*y^2 + x + 2*y + 1\n"
     ]
    }
   ],
   "source": [
    "p=f\n",
    "fs = [f1,f2]\n",
    "s = len(fs)\n",
    "a = [0 for i in range(s)]\n",
    "r=0\n",
    "division_FLAG = False\n",
    "p.lm()/fs[0].lm()\n",
    "print(p.lm()/fs[0].lm())\n",
    "a[0] = a[0] + p.lt()/fs[0].lt()\n",
    "print(a[0])\n",
    "p = p - p.lt()/fs[0].lt()*fs[0]\n",
    "print(p)\n",
    "division_FLAG=True"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 80,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "True\n"
     ]
    },
    {
     "ename": "AttributeError",
     "evalue": "'sage.rings.fraction_field_element.FractionFieldElement' object has no attribute 'lm'",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mAttributeError\u001b[0m                            Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-80-1065f44a378f>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m()\u001b[0m\n\u001b[0;32m      1\u001b[0m \u001b[0mp\u001b[0m \u001b[1;32mis\u001b[0m \u001b[1;32mnot\u001b[0m \u001b[0mInteger\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[0;32m      2\u001b[0m \u001b[1;32mprint\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mp\u001b[0m \u001b[1;32mis\u001b[0m \u001b[1;32mnot\u001b[0m \u001b[0mInteger\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m----> 3\u001b[1;33m \u001b[1;32mprint\u001b[0m\u001b[1;33m(\u001b[0m \u001b[0mp\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mlm\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m/\u001b[0m\u001b[0mfs\u001b[0m\u001b[1;33m[\u001b[0m\u001b[0mInteger\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m]\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mlm\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;33m)\u001b[0m \u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m/home/topolo/Public/sage/src/sage/structure/element.pyx\u001b[0m in \u001b[0;36msage.structure.element.Element.__getattr__ (/home/topolo/Public/sage/src/build/cythonized/sage/structure/element.c:4649)\u001b[1;34m()\u001b[0m\n\u001b[0;32m    411\u001b[0m             \u001b[0mdummy_error_message\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mname\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mname\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    412\u001b[0m             \u001b[1;32mraise\u001b[0m \u001b[0mdummy_attribute_error\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 413\u001b[1;33m         \u001b[1;32mreturn\u001b[0m \u001b[0mgetattr_from_other_class\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mP\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0m_abstract_element_class\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mname\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    414\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    415\u001b[0m     \u001b[1;32mdef\u001b[0m \u001b[0m__dir__\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;32m/home/topolo/Public/sage/src/sage/structure/misc.pyx\u001b[0m in \u001b[0;36msage.structure.misc.getattr_from_other_class (/home/topolo/Public/sage/src/build/cythonized/sage/structure/misc.c:1870)\u001b[1;34m()\u001b[0m\n\u001b[0;32m    257\u001b[0m         \u001b[0mdummy_error_message\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mcls\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mtype\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    258\u001b[0m         \u001b[0mdummy_error_message\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mname\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mname\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 259\u001b[1;33m         \u001b[1;32mraise\u001b[0m \u001b[0mdummy_attribute_error\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    260\u001b[0m     \u001b[1;32mif\u001b[0m \u001b[0misinstance\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mattribute\u001b[0m\u001b[1;33m,\u001b[0m \u001b[0mmethodwrapper\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    261\u001b[0m         \u001b[0mdummy_error_message\u001b[0m\u001b[1;33m.\u001b[0m\u001b[0mcls\u001b[0m \u001b[1;33m=\u001b[0m \u001b[0mtype\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mself\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mAttributeError\u001b[0m: 'sage.rings.fraction_field_element.FractionFieldElement' object has no attribute 'lm'"
     ]
    }
   ],
   "source": [
    "p is not 0\n",
    "print(p is not 0)\n",
    "print( p.lm()/fs[0].lm() )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "I'll try `quo_rem`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 82,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-2*y, 0)"
      ]
     },
     "execution_count": 82,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "p=f\n",
    "p.lt().quo_rem( fs[0].lt() )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 136,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "y\n",
      "-2*y\n",
      "2*y^2 + x + 1\n",
      "0\n",
      "True\n",
      "2*y\n",
      "x + 1\n",
      "(0, x)\n",
      "(0, x)\n",
      "x\n",
      "1\n",
      "False\n",
      "False\n",
      "(0, 1)\n",
      "(0, 1)\n",
      "x + 1\n",
      "0\n",
      "False\n",
      "True\n"
     ]
    }
   ],
   "source": [
    "p=f\n",
    "fs = [f1,f2]\n",
    "s=len(fs)\n",
    "a=[0 for i in range(s)]\n",
    "r=0\n",
    "division_FLAG=False\n",
    "p_dividedby_fs = p.lm().quo_rem(fs[0].lm())\n",
    "print( p_dividedby_fs[0])\n",
    "if p_dividedby_fs[0] is not Poly_F7__2(0):\n",
    "    a[0] += p.lt().quo_rem(fs[0].lt())[0]\n",
    "    print(a[0])\n",
    "    p -= p.lt().quo_rem(fs[0].lt())[0]*fs[0].lt()\n",
    "    print(p)\n",
    "    division_FLAG=True\n",
    "p_dividedby_fs = p.lm().quo_rem(fs[0].lm())\n",
    "print( p_dividedby_fs[0]  )\n",
    "print( p_dividedby_fs[0] is not Poly_F7__2(0)) # notice how we need to \"FORCE\" 0 to belong to \\mathbb{F}_7[x,y]   \n",
    "if p_dividedby_fs[1] is not Poly_F7__2(0):\n",
    "    a[1] += p.lt().quo_rem(fs[1].lt())[0]\n",
    "    print(a[1])\n",
    "    p -= p.lt().quo_rem(fs[1].lt())[0]*fs[1].lt()\n",
    "    print(p)\n",
    "    division_FLAG=True\n",
    "division_FLAG=False\n",
    "p_dividedby_fs = p.lm().quo_rem(fs[0].lm())\n",
    "print( p_dividedby_fs)\n",
    "p_dividedby_fs = p.lm().quo_rem(fs[1].lm())\n",
    "print( p_dividedby_fs)\n",
    "if division_FLAG is False:\n",
    "    r += p.lt()\n",
    "    print(r)\n",
    "    p -= p.lt()\n",
    "    print(p)\n",
    "print(p is Poly_F7__2(0)) \n",
    "print(p == Poly_F7__2(0))\n",
    "division_FLAG=False\n",
    "p_dividedby_fs = p.lm().quo_rem(fs[0].lm())\n",
    "print( p_dividedby_fs)\n",
    "p_dividedby_fs = p.lm().quo_rem(fs[1].lm())\n",
    "print( p_dividedby_fs)\n",
    "if division_FLAG is False:\n",
    "    r += p.lt()\n",
    "    print(r)\n",
    "    p -= p.lt()\n",
    "    print(p)\n",
    "print(p is Poly_F7__2(0)) # vs.\n",
    "print(p == Poly_F7__2(0))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 138,
   "metadata": {},
   "outputs": [],
   "source": [
    "def DivisionAlgo( p, fs, ring):\n",
    "    \"\"\" DivisionAlgo - Division algorithm\n",
    "    \"\"\"\n",
    "    s=len(fs)\n",
    "    zeroelement = ring(0)    \n",
    "    a=[zeroelement for i in range(s)]\n",
    "    r=zeroelement\n",
    "    while p != zeroelement:\n",
    "        division_FLAG=False\n",
    "        for i in range(s):\n",
    "            p_dividedby_fs = p.lm().quo_rem(fs[i].lm())\n",
    "            if p_dividedby_fs[1] == zeroelement:\n",
    "                a[i] += p.lt().quo_rem(fs[i].lt())[0]  \n",
    "                p -= p.lt().quo_rem(fs[i].lt())[0]*fs[i].lt()\n",
    "                division_FLAG=True\n",
    "                break\n",
    "        if division_FLAG is False:\n",
    "            r += p.lt()\n",
    "            p -= p.lt()\n",
    "    return r,a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 141,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(x + 1, [-2*y, 2*y])\n"
     ]
    }
   ],
   "source": [
    "r_eg, a_eg = DivisionAlgo( f,[f1,f2],Poly_F7__2)\n",
    "print( r_eg,a_eg)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercises"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 152,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Poly_F7__2 = PolynomialRing(FiniteField(7),'x,y',order='deglex') # \\mathbb{F}_7[x,y] so it has 2 variables\n",
    "f=x**7*y**2+x**3*y**2-y+1\n",
    "f1=x*y**2-x\n",
    "f2=x-y**3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 153,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-y + 1, [x^6 + x^2, 0])"
      ]
     },
     "execution_count": 153,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo(f,[f1,f2],Poly_F7__2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 154,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Poly_F7__2 = PolynomialRing(FiniteField(7),'x,y',order='lex') # \\mathbb{F}_7[x,y] so it has 2 variables\n",
    "f=x**7*y**2+x**3*y**2-y+1\n",
    "f1=x*y**2-x\n",
    "f2=x-y**3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 155,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-y + 1, [x^6 + x^2, 0])"
      ]
     },
     "execution_count": 155,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo(f,[f1,f2],Poly_F7__2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 156,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Poly_F7__2 = PolynomialRing(FiniteField(7),'x,y',order='lex') # \\mathbb{F}_7[x,y] so it has 2 variables\n",
    "f=x*y**2-x\n",
    "f1=x*y+1\n",
    "f2=y**2-1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 157,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-x, [y, 0])"
      ]
     },
     "execution_count": 157,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo(f,[f1,f2],Poly_F7__2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 158,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(-x, [x, 0])"
      ]
     },
     "execution_count": 158,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo(f,[f2,f1],Poly_F7__2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 159,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "Poly_F7__3 = PolynomialRing(FiniteField(7),'x,y,z',order='deglex') # \\mathbb{F}_7[x,y,z] so it has 3 variables\n",
    "f=x*y**2*z**2+x*y-y*z\n",
    "f1=x-y**2\n",
    "f2=y-z**3\n",
    "f3=z**2-1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 160,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(x*y - y*z, [-x*z^2, 0, 0])"
      ]
     },
     "execution_count": 160,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo(f,[f1,f2,f3],Poly_F7__3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 161,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(x*y - y*z, [0, x*y^2, 0])"
      ]
     },
     "execution_count": 161,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo(f,[f2,f3,f1],Poly_F7__3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 162,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(x*y - y*z, [x*y^2, 0, 0])"
      ]
     },
     "execution_count": 162,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "DivisionAlgo(f,[f3,f1,f2],Poly_F7__3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Ideals\n",
    "### Definition"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 163,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z> = PolynomialRing(FiniteField(127),order='deglex')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 164,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "I=ideal(x*y+z,y^3+1,z^2-x*5-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 165,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 165,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(x*y+z)+(y^3+1) in I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 166,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 166,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x*z*(z^2-x*5-1) in I"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Varieties \n",
    "### Varieties III\n",
    "#### Example"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 167,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[{y: 0, x: 1}, {y: 0, x: 6}, {y: 1, x: 0}, {y: 6, x: 0}]"
      ]
     },
     "execution_count": 167,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P.<x,y> = PolynomialRing(FiniteField(7))\n",
    "I = Ideal(x^2+y^2-1,x*y)\n",
    "I.variety()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Exercises"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 168,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y> = PolynomialRing(FiniteField(7))\n",
    "I = Ideal(x*y+1,y**2-1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 169,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 169,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x+y in I"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$\\lbrace 1 \\rbrace$ is **not** an ideal since it doesn't include the $0$ element.  $\\lbrace x+1,y-2 \\rbrace$ is not an ideal because it doesn't include the $0$ element"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 172,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 172,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ideal( P(1) ) == Ideal( x+1,x-1,y-2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 174,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[{y: 1, x: 6}]"
      ]
     },
     "execution_count": 174,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ideal(x+1,y-1).variety()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 176,
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "(1, x) is not a valid variable.",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[1;32m<ipython-input-176-63ec9a389e83>\u001b[0m in \u001b[0;36m<module>\u001b[1;34m()\u001b[0m\n\u001b[1;32m----> 1\u001b[1;33m \u001b[0msolve\u001b[0m\u001b[1;33m(\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[0mP\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mInteger\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m,\u001b[0m\u001b[0mx\u001b[0m\u001b[1;33m+\u001b[0m\u001b[0mP\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mInteger\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[0mP\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mInteger\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;36m0\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m]\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[1;33m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[1;32m/home/topolo/Public/sage/local/lib/python2.7/site-packages/sage/symbolic/relation.pyc\u001b[0m in \u001b[0;36msolve\u001b[1;34m(f, *args, **kwds)\u001b[0m\n\u001b[0;32m    826\u001b[0m     \u001b[1;32mfor\u001b[0m \u001b[0mv\u001b[0m \u001b[1;32min\u001b[0m \u001b[0mvariables\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    827\u001b[0m         \u001b[1;32mif\u001b[0m \u001b[1;32mnot\u001b[0m \u001b[0mis_SymbolicVariable\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mv\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[1;32m--> 828\u001b[1;33m             \u001b[1;32mraise\u001b[0m \u001b[0mTypeError\u001b[0m\u001b[1;33m(\u001b[0m\u001b[1;34m\"%s is not a valid variable.\"\u001b[0m\u001b[1;33m%\u001b[0m\u001b[0mrepr\u001b[0m\u001b[1;33m(\u001b[0m\u001b[0mv\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m)\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0m\u001b[0;32m    829\u001b[0m \u001b[1;33m\u001b[0m\u001b[0m\n\u001b[0;32m    830\u001b[0m     \u001b[1;32mtry\u001b[0m\u001b[1;33m:\u001b[0m\u001b[1;33m\u001b[0m\u001b[0m\n",
      "\u001b[1;31mTypeError\u001b[0m: (1, x) is not a valid variable."
     ]
    }
   ],
   "source": [
    "solve([x+y==P(0),x+P(1)==P(0)],x,y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 177,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[[x == -1, y == 1]]"
      ]
     },
     "execution_count": 177,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x,y=var('x y')\n",
    "solve([x+y==0,x+1==0],x,y)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Gröbner Bases\n",
    "### Definition II"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 179,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "R.<x> = PolynomialRing(FiniteField(7))\n",
    "f=x^2+6"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 180,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "I=Ideal(map(P,[R.random_element() * f for _ in range(5)]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 181,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^2 - 1]"
      ]
     },
     "execution_count": 181,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.groebner_basis()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "EY:20170311 What is `map` doing?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 182,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[-x^4 + x^3 - 3*x^2 - x - 3,\n",
       " x^4 + x^3 + x^2 - x - 2,\n",
       " 3*x^4 + 2*x^3 + x^2 - 2*x + 3,\n",
       " 2*x^4 - x^3 + 3*x^2 + x + 2,\n",
       " -x^4 - x^3 + x + 1]"
      ]
     },
     "execution_count": 182,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "map(P,[R.random_element() * f for _ in range(5)])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 194,
   "metadata": {},
   "outputs": [],
   "source": [
    "P.<x,y,z> = PolynomialRing(FiniteField(127),order='deglex')\n",
    "F=Sequence([-3*y,-2*x-y-3*z+2,x+y+2*z-1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 195,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x - 1, y, z]"
      ]
     },
     "execution_count": 195,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "F.groebner_basis()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 196,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "A,v=F.coefficient_matrix()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 197,
   "metadata": {},
   "outputs": [],
   "source": [
    "A.echelonize()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 198,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x - 1     y     z]"
      ]
     },
     "execution_count": 198,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(A*v).T"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "EY : 20170311 what is going on?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 200,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[  1   0   0 126]\n",
      "[  0   1   0   0]\n",
      "[  0   0   1   0]\n",
      "[x]\n",
      "[y]\n",
      "[z]\n",
      "[1]\n"
     ]
    }
   ],
   "source": [
    "print(A); print(v)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 202,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[  0 124   0   0]\n",
      "[125 126 124   2]\n",
      "[  1   1   2 126]\n"
     ]
    }
   ],
   "source": [
    "print( F.coefficient_matrix()[0] )"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Definition IV\n",
    "You can use Sage to check if a list of polynomials is a Gröbner basis:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 203,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 203,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(7),order=\"degrevlex\")\n",
    "I=Ideal([x*y^2-z,2*x^2*z-2,x*y-z+1])\n",
    "I.basis_is_groebner()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 204,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^2 - x + z - 2, x*y - z + 1, y^2 - x - 2*y + 1, x*z - x - y + 1, y*z - y - z, z^2 - x - y - 2*z + 2]"
      ]
     },
     "execution_count": 204,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.groebner_basis()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Applications I\n",
    "#### Representation\n",
    "Reduced Gröbner bases are unique representations of an ideal with respect to a monomial ordering. So to check that two ideals are the same, you only need to compute their reduced Gröbner bases. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 205,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(7))\n",
    "I=Ideal(P.random_element() for _ in range(5))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 207,
   "metadata": {},
   "outputs": [],
   "source": [
    "J = Ideal(P.random_element() for _ in range(5))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 208,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "I=Ideal(f-f((2,3,1)) for f in I.gens())\n",
    "J=Ideal(f-f((2,3,1)) for f in J.gens())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 209,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 209,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.groebner_basis() == J.groebner_basis()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 210,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 210,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I==J"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 211,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 211,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.gens()==J.gens()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 212,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x - 2, y - 3, z - 1]"
      ]
     },
     "execution_count": 212,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.groebner_basis()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Applications II\n",
    "#### Membership\n",
    "remainder $r$ of division of any $f\\in P$ by $G$ is unique.  Hence, we can check whether $f\\in \\langle f_1,\\dots, f_s\\rangle$ by computing a Gröbner basis $G=(g_1,\\dots, g_t)$ for $\\langle f_1,\\dots , f_s\\rangle$ and check if \\overline{f}^G =0$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 220,
   "metadata": {},
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(7))\n",
    "I=Ideal(P.random_element() for _ in range(5))\n",
    "I=Ideal(f-f((2,3,1)) for f in I.gens())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 222,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "f=P.random_element()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 223,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 223,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f in I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 224,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "f = f-f((2,3,1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 225,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 225,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f in I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 226,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-3*y*z + 2*z^2 - 2*x - y"
      ]
     },
     "execution_count": 226,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.reduce(I.gens())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 227,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0"
      ]
     },
     "execution_count": 227,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "f.reduce(I.groebner_basis())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 215,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x + 1, y - 2]"
      ]
     },
     "execution_count": 215,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Applications III\n",
    "#### Solving \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 230,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z> = PolynomialRing(FiniteField(7), order=\"lex\")\n",
    "I=Ideal(P.random_element() for _ in range(5))\n",
    "I = Ideal(f-f((2,3,1)) for f in I.gens())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 231,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x + z - 3, y + 2*z + 2, z^2 + 3*z + 3]"
      ]
     },
     "execution_count": 231,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.groebner_basis()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 233,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[-y^2 + y*z + 2*z^2, 3*x*y - 2*x + y^2 + 2*z^2 - 1, -2*x^2 - 3*x*y + x*z + 3*z^2 - 3, 3*x^2 - y^2 + 2*z^2 + 2, -x*z + y^2 + 3*y*z + z^2 - 1]"
      ]
     },
     "execution_count": 233,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I=Ideal(P.random_element() for _ in range(5))\n",
    "I.gens()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercises"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 234,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x + 1, y - 2]"
      ]
     },
     "execution_count": 234,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P.<x,y>=PolynomialRing(FiniteField(7),order=\"degrevlex\")\n",
    "F=Sequence([x+1,y-2,x*y-2*x+y-2])\n",
    "F.groebner_basis()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 238,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 238,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ideal(F).basis_is_groebner()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 241,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 241,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ideal([x+1,y-2,x*y-2*x+y-2]).basis_is_groebner()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 239,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 239,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ideal( Sequence([2*x+2,3*y+1]) ).basis_is_groebner()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 240,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 240,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ideal( [2*x+2,3*y+1] ).basis_is_groebner()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Buchberger's Algorithm  \n",
    "### S-Polynomials I"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### S-Polynomials IV"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### S-Polynomials V"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### S-Polynomials VI"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 274,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y>=PolynomialRing(QQ,order='degrevlex')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 275,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "f0 = x^3-2*x*y\n",
    "f1=x**2*y-2*y**2+x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 278,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-x^2"
      ]
     },
     "execution_count": 278,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(x^3*y)//x^3*f0 - (x^3*y)//(x^2*y)*f1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that `//` is the same as `.quo_rem[0]`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 276,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "y"
      ]
     },
     "execution_count": 276,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(x**3*y)//x**3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 277,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 277,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(x**3*y)//x**3 == (x**3*y).quo_rem(x**3)[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is \n",
    "$$  \n",
    "x^{\\gamma} = \\text{LCM}{ (\\text{LM}(f), \\text{LM}(g) ) }\n",
    "$$ "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 282,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x^3*y\n"
     ]
    }
   ],
   "source": [
    "xgamma = lcm(f0.lm(), f1.lm() ); print(xgamma)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 283,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-x^2"
      ]
     },
     "execution_count": 283,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "xgamma//f0.lt()*f0 - xgamma//f1.lt()*f1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Buchberger's Criterion II"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 284,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "f1=x**3-2*x*y\n",
    "f2=x**2*y-2*y**2+x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 285,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def make_Spoly(f1,f2):\n",
    "    \"\"\" make_Spoly : make the S-polynomial \"\"\"\n",
    "    xgamma = lcm(f1.lm(),f2.lm())\n",
    "    Spoly  = xgamma//f1.lt()*f1 - xgamma//f2.lt()*f2\n",
    "    return Spoly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 286,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-x^2"
      ]
     },
     "execution_count": 286,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "make_Spoly(f1,f2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 290,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "0\n"
     ]
    }
   ],
   "source": [
    "print( make_Spoly(f1,f2)//f1.lm() )\n",
    "print( make_Spoly(f1,f2)//f2.lm() )\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Buchberger's Algorithm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 382,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def BuchbergerAlgo(F,polyring):\n",
    "    \"\"\" BuchbergerAlgo : Buchberger's Algorithm\n",
    "    INPUT\n",
    "    =====\n",
    "    @type F : Python list of Sage Math polynomials\n",
    "    @param F : F = (f_1, \\dots, f_s)\n",
    "    \n",
    "    @type polyring : Sage Math PolynomialRing\n",
    "    \"\"\"\n",
    "    zeroelement = polyring(0)\n",
    "    pairs=[_ for _ in combinations(F,2)]\n",
    "#    NO_PAIRS = len(pairs)\n",
    "    G=F\n",
    "    G_2 = [] # denoted G' or Gprime in Cox, Little, O'Shea (2015)\n",
    "\n",
    "    # first case guaranteed to run \n",
    "    while G_2 != G:\n",
    "        G_2 = G\n",
    "        for pair in pairs:\n",
    "            p=pair[0]\n",
    "            q=pair[1]\n",
    "            r,a_list = DivisionAlgo(make_Spoly(p,q), G, polyring)\n",
    "            if r != zeroelement:\n",
    "                G.append(r)\n",
    "        # END of while, break condition being when we don't add,append, anymore to G\n",
    "    return G"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Optional Improvements II  \n",
    "\n",
    "#### Definition (Buchberger's First Criterion)  \n",
    "Let $f,g\\in \\mathbb{F}[x_1,\\dots, x_n]$ with disjoint leading terms, i.e.  \n",
    "$$  \n",
    "\\text{LCM}(\\text{LM}(f), \\text{LM}(g))  = \\text{LM}(f) \\cdot \\text{LM}(g)  \n",
    "$$  \n",
    "\n",
    "Then $\\overline{S(f,g)}^G =0$  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 401,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def BuchbergerAlgo_w_1stCond(F,polyring):\n",
    "    \"\"\" BuchbergerAlgo_w_1stCond : Buchberger's Algorithm, with Buchberger's First Criterion\n",
    "    INPUT\n",
    "    =====\n",
    "    @type F : Python list of Sage Math polynomials\n",
    "    @param F : F = (f_1, \\dots, f_s)\n",
    "    \n",
    "    @type polyring : Sage Math PolynomialRing\n",
    "    \"\"\"\n",
    "    zeroelement = polyring(0)\n",
    "    pairs=[_ for _ in combinations(F,2)]\n",
    "#    NO_PAIRS = len(pairs)\n",
    "    G=F\n",
    "    G_2 = [] # denoted G' or Gprime in Cox, Little, O'Shea (2015)\n",
    "\n",
    "    counter = 0 # to count how many reductions to 0 of \\overline{S(f,g)}^G; you can comment this line out\n",
    "    \n",
    "    # first case guaranteed to run \n",
    "    while G_2 != G:\n",
    "        G_2 = G\n",
    "        for pair in pairs:\n",
    "            p=pair[0]\n",
    "            q=pair[1]\n",
    "            # Buchberger's First Criterion, check, before doing the Division Algorithm, DivisionAlgo\n",
    "            Criterion1_LHS = lcm(p.lm(), q.lm())\n",
    "            Criterion1_RHS = p.lm() * q.lm()\n",
    "\n",
    "            # you can comment out the next 2 lines\n",
    "            if Criterion1_LHS == Criterion1_RHS:\n",
    "                counter +=1\n",
    "            \n",
    "            if Criterion1_LHS != Criterion1_RHS:\n",
    "                r,a_list = DivisionAlgo(make_Spoly(p,q), G, polyring)\n",
    "            #if r != zeroelement:\n",
    "                G.append(r)\n",
    "        # END of while, break condition being when we don't add,append, anymore to G\n",
    "    print counter\n",
    "    return G"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercises"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "1. $\\mathbb{F}_7[x,y]$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 391,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[x^2*y - 1, x*y^2 - x, x^2 - y]"
      ]
     },
     "execution_count": 391,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Poly_F7__2.<x,y>=PolynomialRing(FiniteField(7),2,order=\"deglex\")\n",
    "BuchbergerAlgo_w_1stCond([x**2*y+6,x*y**2-x], Poly_F7__2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 392,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^2 - y, y^2 - 1]"
      ]
     },
     "execution_count": 392,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Sage Math sanity check\n",
    "Ideal([x**2*y+6,x*y**2-x]).groebner_basis()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 395,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^2 - y, y^2 - 1]"
      ]
     },
     "execution_count": 395,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Poly_F7__2.ideal([x**2*y+6,x*y**2-x]).groebner_basis()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 393,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^2*y - 1, x*y^2 - x, x^2 - y]"
      ]
     },
     "execution_count": 393,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "BuchbergerAlgo([x**2*y+6,x*y**2-x], Poly_F7__2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "2 . $\\mathbb{F}_7[x,y,z]$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 396,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[x^2 + y, x^4 + 2*x^2*y + y^2 + 3, -3]"
      ]
     },
     "execution_count": 396,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Poly_F7__3.<x,y,z>=PolynomialRing(FiniteField(7),3,order=\"deglex\")\n",
    "BuchbergerAlgo_w_1stCond([x**2+y,x**4+2*x**2*y+y**2+3], Poly_F7__3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 397,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^2 + y, x^4 + 2*x^2*y + y^2 + 3, -3]"
      ]
     },
     "execution_count": 397,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "BuchbergerAlgo([x**2+y,x**4+2*x**2*y+y**2+3], Poly_F7__3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 398,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[1]"
      ]
     },
     "execution_count": 398,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Poly_F7__3.ideal([x**2+y,x**4+2*x**2*y+y**2+3]).groebner_basis()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "EY : 20170312 - No solution (???)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "3 . $\\langle x-z^4,y-z^5\\rangle$ "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 402,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "[-z^4 + x, -z^5 + y, -x*z + y]"
      ]
     },
     "execution_count": 402,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "BuchbergerAlgo_w_1stCond([x-z**4,y-z**5], Poly_F7__3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 404,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[-z^4 + x, -z^5 + y, -x*z + y]"
      ]
     },
     "execution_count": 404,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "BuchbergerAlgo([x-z**4,y-z**5], Poly_F7__3)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 405,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x^4 - y^3*z, y^2*z^2 - x^3, y*z^3 - x^2, z^4 - x, x*z - y]"
      ]
     },
     "execution_count": 405,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Poly_F7__3.ideal([x-z**4,y-z**5]).groebner_basis()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Quotient Rings I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 407,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Quotient of Multivariate Polynomial Ring in x, y, z over Finite Field of size 2 by the ideal (x^2 + x, y^2 + y, z^2 + z)"
      ]
     },
     "execution_count": 407,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P.<x,y,z>=PolynomialRing(FiniteField(2))\n",
    "I=sage.rings.ideal.FieldIdeal(P)\n",
    "Q=P.quotient_ring(I); Q"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 411,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[x^2 + x, y^2 + y, z^2 + z]\n",
      "[x^2 + x, y^2 + y, z^2 + z]\n"
     ]
    }
   ],
   "source": [
    "print(I.basis);print(I.groebner_basis())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 412,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Ideal (x^2 + x, y^2 + y, z^2 + z) of Multivariate Polynomial Ring in x, y, z over Finite Field of size 2"
      ]
     },
     "execution_count": 412,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "P.<x,y,z> = BooleanPolynomialRing() # much faster!\n",
    "P.defining_ideal()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Solving Polynomial Systems with Grobner Bases  \n",
    "### Elimination Ideals I"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Elimination Theorem II  \n",
    "##### Example"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 413,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z> = PolynomialRing(FiniteField(127),order='lex')\n",
    "I = sage.rings.ideal.Cyclic(P)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 414,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Ideal (x + y + z, x*y + x*z + y*z, x*y*z - 1) of Multivariate Polynomial Ring in x, y, z over Finite Field of size 127"
      ]
     },
     "execution_count": 414,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 415,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "J=I+sage.rings.ideal.FieldIdeal(P)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These are the so-called field polynomials, that we augment our ideal with"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 416,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Ideal (x^127 - x, y^127 - y, z^127 - z) of Multivariate Polynomial Ring in x, y, z over Finite Field of size 127"
      ]
     },
     "execution_count": 416,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sage.rings.ideal.FieldIdeal(P)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Observe how this Groebner basis is each only in terms of 1 less variable (univariate variable); these are \n",
    "$$  \n",
    "G_l=G\\cap \\mathbb{F}[x_l,\\dots,x_n]  \n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 417,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(x + y + z, y^2 + y*z + z^2, z^3 - 1)"
      ]
     },
     "execution_count": 417,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "g0,g1,g2=J.groebner_basis(); g0,g1,g2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 418,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(z - 19) * (z - 1) * (z + 20)"
      ]
     },
     "execution_count": 418,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "factor(g2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 419,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(y - 1) * (y + 20)"
      ]
     },
     "execution_count": 419,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "factor(g1(x,y,19))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 420,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "x + 20"
      ]
     },
     "execution_count": 420,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "factor(g0(x,1,19))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 421,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 421,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "all(f(107,1,19)==0 for f in I.gens())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 422,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[{y: 19, z: 1, x: 107},\n",
       " {y: 107, z: 1, x: 19},\n",
       " {y: 1, z: 19, x: 107},\n",
       " {y: 107, z: 19, x: 1},\n",
       " {y: 1, z: 107, x: 19},\n",
       " {y: 19, z: 107, x: 1}]"
      ]
     },
     "execution_count": 422,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "J.variety()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### An Example II"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 431,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "x,y=var(\"x,y\")\n",
    "p1 = (x-1)**2+(y-2)**2-3**2==0\n",
    "p2 = (x+1)**2+(y-1)**2-2**2==0\n",
    "sols=solve([p1,p2],x,y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 432,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[[x == -2/5*sqrt(5) - 1, y == 4/5*sqrt(5) + 1], [x == 2/5*sqrt(5) - 1, y == -4/5*sqrt(5) + 1]]"
      ]
     },
     "execution_count": 432,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sols"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 439,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-1.89442719099992\n",
      "2.78885438199983\n",
      "-0.105572809000084\n",
      "-0.788854381999832\n"
     ]
    }
   ],
   "source": [
    "print((-2/5*sqrt(5) - 1).N());print((4/5*sqrt(5) + 1).N());print((2/5*sqrt(5) - 1).N());print((-4/5*sqrt(5) + 1).N());"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 426,
   "metadata": {},
   "outputs": [],
   "source": [
    "P.<x,y>=PolynomialRing(QQ,order='lex')\n",
    "f1=(x-1)**2+(y-2)**2-3**2\n",
    "f2=(x+1)**2+(y-1)**2-2**2\n",
    "I=Ideal(f1,f2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 427,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Ideal (x^2 - 2*x + y^2 - 4*y - 4, x^2 + 2*x + y^2 - 2*y - 2) of Multivariate Polynomial Ring in x, y over Rational Field"
      ]
     },
     "execution_count": 427,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 428,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[x + 1/2*y + 1/2, y^2 - 2*y - 11/5]"
      ]
     },
     "execution_count": 428,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.groebner_basis()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 429,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[]"
      ]
     },
     "execution_count": 429,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.variety()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 430,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[{y: -0.788854381999832, x: -0.105572809000084},\n",
       " {y: 2.78885438199983, x: -1.89442719099992}]"
      ]
     },
     "execution_count": 430,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.variety(CC)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 440,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<z>=QQ[]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 441,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[{y: -4/5*a + 1, x: 2/5*a - 1}, {y: 4/5*a + 1, x: -2/5*a - 1}]"
      ]
     },
     "execution_count": 441,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I.variety(NumberField(z^2-5,'a'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## The $\\mathbb{F}_4$ Algorithm"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Warm Up"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As a warm-up, consider a linear system of equations over $\\mathbb{F}_{127}[x,y,z]$  \n",
    "\n",
    "$$  \n",
    "f = 26y+52z + 62 = 0  \\\\\n",
    "g = 54y+119z+55=0  \\\\\n",
    "h =41x+91z + 13=0  \\\\\n",
    "$$  \n",
    "\n",
    "$\\Longrightarrow $   \n",
    "$$  \n",
    "\\left( \\begin{matrix} 0 & 26 & 52 & 62 \\\\ \n",
    "0 & 54 & 119 & 55 \\\\\n",
    "41 & 0 & 91 & 13 \\end{matrix} \\right) $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For these explicit steps, I looked up [Lecture 10](http://www.calvin.edu/~scofield/courses/m231/S11/calendar/lects/lect10.pdf).  Also, see [How to show the steps of Gauss' method](https://ask.sagemath.org/question/8840/how-to-show-the-steps-of-gauss-method/) and toy implementation by [jh](https://ask.sagemath.org/users/824/jh/)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Naive Gaussian reduction\n",
    "def gauss_method(M,rescale_leading_entry=False):\n",
    "    \"\"\"Describe the reduction to echelon form of the given matrix of rationals.\n",
    "\n",
    "    M  matrix of rationals   e.g., M = matrix(QQ, [[..], [..], ..])\n",
    "    rescale_leading_entry=False  boolean  make the leading entries to 1's\n",
    "\n",
    "    Returns: None.  Side effect: M is reduced.  Note: this is echelon form, \n",
    "    not reduced echelon form; this routine does not end the same way as does \n",
    "    M.echelon_form().\n",
    "\n",
    "    \"\"\"\n",
    "    num_rows=M.nrows()\n",
    "    num_cols=M.ncols()\n",
    "    print M    \n",
    "\n",
    "    col = 0   # all cols before this are already done\n",
    "    for row in range(0,num_rows): \n",
    "        # ?Need to swap in a nonzero entry from below\n",
    "        while (col < num_cols\n",
    "               and M[row][col] == 0): \n",
    "            for i in M.nonzero_positions_in_column(col):\n",
    "                if i > row:\n",
    "                    print \" swap row\",row+1,\"with row\",i+1\n",
    "                    M.swap_rows(row,i)\n",
    "                    print M\n",
    "                    break     \n",
    "            else:\n",
    "                col += 1\n",
    "\n",
    "        if col >= num_cols:\n",
    "            break\n",
    "\n",
    "        # Now guaranteed M[row][col] != 0\n",
    "        if (rescale_leading_entry\n",
    "           and M[row][col] != 1):\n",
    "            print \" take\",1/M[row][col],\"times row\",row+1\n",
    "            M.rescale_row(row,1/M[row][col])\n",
    "            print M\n",
    "        change_flag=False\n",
    "        for changed_row in range(row+1,num_rows):\n",
    "            if M[changed_row][col] != 0:\n",
    "                change_flag=True\n",
    "                factor=-1*M[changed_row][col]/M[row][col]\n",
    "                print \" take\",factor,\"times row\",row+1,\"plus row\",changed_row+1 \n",
    "                M.add_multiple_of_row(changed_row,row,factor)\n",
    "        if change_flag:\n",
    "            print M\n",
    "        col +=1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 489,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[  0  26  52  62]\n",
      "[  0  54 119  55]\n",
      "[ 41   0  91  13]\n",
      "[  0  26  52  62]\n",
      "[  0  54 119  55]\n",
      "[  1   0  27  22]\n",
      "[  1   0  27  22]\n",
      "[  0  54 119  55]\n",
      "[  0  26  52  62]\n"
     ]
    }
   ],
   "source": [
    "A=Matrix( FiniteField(127),[[0,26,52,62],[0,54,119,55],[41,0,91,13]])\n",
    "cA = copy(A)\n",
    "cA.rescale_row(2,1/41)\n",
    "print(A);print(cA);\n",
    "cA.swap_rows(2,0); print(cA)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 490,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1"
      ]
     },
     "execution_count": 490,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lcm(FiniteField(127)(54),FiniteField(127)(26))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 491,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[  1   0  27  22]\n",
      "[  0  54 119  55]\n",
      "[  0   1   2  61]\n",
      "[ 1  0 27 22]\n",
      "[ 0  1 61 41]\n",
      "[ 0  1  2 61]\n",
      "[ 1  0 27 22]\n",
      "[ 0  1 61 41]\n",
      "[ 0  0 68 20]\n"
     ]
    }
   ],
   "source": [
    "cA.rescale_row(2,1/26);print(cA);cA.rescale_row(1,1/54);print(cA);\n",
    "cA.add_multiple_of_row(2,1,-1);print(cA)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 492,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 1  0  0 29]\n",
      "[ 0  1  0 38]\n",
      "[ 0  0  1 75]\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "\\left(\\begin{array}{rrrr}\n",
       "1 & 0 & 0 & 29 \\\\\n",
       "0 & 1 & 0 & 38 \\\\\n",
       "0 & 0 & 1 & 75\n",
       "\\end{array}\\right)"
      ]
     },
     "execution_count": 492,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cA.rescale_row(2,1/68);cA.add_multiple_of_row(1,2,66);cA.add_multiple_of_row(0,2,100);print(cA); latex(cA)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(So) After Gaussian elimination:  \n",
    "$$\n",
    "f'=x+29=0  \\\\\n",
    "g' = y + 38  \\\\\n",
    "h' = z+75 =0 \n",
    "$$  \n",
    "$\\Longrightarrow$ \n",
    "$$  \n",
    "\\left(\\begin{array}{rrrr}\n",
    "1 & 0 & 0 & 29 \\\\\n",
    "0 & 1 & 0 & 38 \\\\\n",
    "0 & 0 & 1 & 75\n",
    "\\end{array}\\right)\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Thus, $x=-27$, $y=-38$, and $z=-75$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Generalizing to Non-Linear Polynomials"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now consider 2 polynomials in $\\mathbb{F}_{127}[x,y,z]$ with term ordering deglex"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$  \n",
    "f=x^2+2xy-2y^2+14z^2+22z  \\\\\n",
    "g = x^2 + xy + y^2 + z^2+x+2z  $$\n",
    "$\\Longrightarrow $ \n",
    "$$\n",
    "\\left( \\begin{matrix} 1 & 2 & -2 & 14 & 0 & 22 \\\\\n",
    "1 & 1 & 1 & 1 & 1 & 2 \\end{matrix} \\right)$$  \n",
    "\n",
    "$$  \n",
    "f=x^2 +4y^2 - 12z^2 + 2x-18z  \\\\\n",
    "g'=xy+-3y^2 + 13z^2 - x + 20z   \n",
    "$$  \n",
    "$$\n",
    "\\left( \\begin{matrix} 1 & 0 & 4 & -12 & 2 & - 18 \\\\\n",
    "& 1  & - 3 & 13 & -1 & 20 \\end{matrix} \\right)$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### $\\mathbb{F}_4$ Algorithm I"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 504,
   "metadata": {},
   "outputs": [],
   "source": [
    "Poly_F127__3.<x,y,z>=PolynomialRing(FiniteField(127),3,order='deglex')\n",
    "f=x**2-2*x*y-2*y**2+14*z**2\n",
    "g=x+y+2*z\n",
    "pairs=[_ for _ in combinations([f,g],2)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "EY:20170313, I'm not sure if the so-called slimgb algorithm is a version of $F_4$ Algorithm - someone please work it out and confirm this (cf. [Faugère's F4 Algorithm](https://ask.sagemath.org/question/8638/faugeres-f4-algorithm/) vs. http://doc.sagemath.org/html/en/reference/polynomial_rings/sage/rings/polynomial/multi_polynomial_ideal.html, "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 510,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[y^2 + 8*y*z + 18*z^2, x + y + 2*z]\n",
      "[y^2 + 8*y*z + 18*z^2, x + y + 2*z]\n"
     ]
    }
   ],
   "source": [
    "print( Ideal([f,g]).groebner_basis(\"singular:slimgb\"))\n",
    "print( Ideal([f,g]).groebner_basis(\"libsingular:slimgb\"))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 503,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x^2\n",
      "x\n"
     ]
    }
   ],
   "source": [
    "print( pairs[0][0].lm());print( pairs[0][1].lm()) # lm for leading monomial, it doesn't include the coefficient\n",
    "print( )\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 515,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-3*x*y - 2*x*z - 2*y^2 + 14*z^2"
      ]
     },
     "execution_count": 515,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "make_Spoly(f,g)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These are the so-called field polynomials, that we augment our ideal with"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "cf. [Term orders, Polynomials, Multivariate Polynomials and Polynomial Rings](http://doc.sagemath.org/html/en/reference/polynomial_rings/sage/rings/polynomial/term_order.html)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 299,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(QQ,3,order='lex')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 251,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 251,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x>y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 252,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 252,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x>y^2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 253,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 253,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x>1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 254,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 254,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x^1*y^2 >y^3*z^4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 255,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 255,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x^3*y^2*z^4<x^3*y^2*z^1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 265,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{(1, 0, 0): 1}\n",
      "{(0, 1, 0): 1}\n",
      "{(0, 2, 0): 1}\n",
      "{(1, 2, 0): 1}\n",
      "{(0, 3, 4): 1}\n",
      "{(3, 2, 4): 1}\n",
      "{(3, 2, 1): 1}\n"
     ]
    }
   ],
   "source": [
    "print(x.dict());print(y.dict());print((y**2).dict());print((x^1*y^2).dict());print( (y**3*z**4).dict());\n",
    "print((x^3*y^2*z^4).dict());print((x^3*y^2*z^1).dict())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 266,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "P.<x,y,z>=PolynomialRing(QQ,3,order='deglex')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 267,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 267,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x>y"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 268,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "False"
      ]
     },
     "execution_count": 268,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x>y^2*z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 269,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 269,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x>1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 270,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 270,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x^1*y^2*z^3>x^3*y^2*z^0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 271,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "True"
      ]
     },
     "execution_count": 271,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "x^2*y*z^2>x*y^3*z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 273,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{(1, 0, 0): 1}\n",
      "{(0, 1, 0): 1}\n",
      "{(0, 2, 1): 1}\n",
      "{(0, 0, 0): 1}\n",
      "{(1, 2, 3): 1}\n",
      "{(3, 2, 0): 1}\n",
      "{(2, 1, 2): 1}\n",
      "{(1, 3, 1): 1}\n"
     ]
    }
   ],
   "source": [
    "print(x.dict());print(y.dict());print((y**2*z).dict());print(P(1).dict());\n",
    "print((x^1*y^2*z^3).dict());print( (x**3*y**2*z**0).dict());\n",
    "print((x^2*y^1*z^2).dict());print((x^1*y^3*z^1).dict())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "cf. http://icm.mcs.kent.edu/reports/1995/gb.pdf"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 300,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "g1=4*x**2*z-7*y**2\n",
    "g2=x*y*z**2+3*x*z**4\n",
    "G=[g1,g2]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 292,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from itertools import combinations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 301,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Algorithm Buchberger\n",
    "M = [_ for _ in combinations(G,2)]\n",
    "pair = M.pop()\n",
    "S = make_Spoly(pair[0],pair[1])\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 302,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-3*x^2*z^4 - 7/4*y^3*z"
      ]
     },
     "execution_count": 302,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "S"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 383,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "4*x^3*y*z^3 + 12*x^3*z^5 - 7*x*y^3*z^2 - 21*x*y^2*z^4"
      ]
     },
     "execution_count": 383,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lcm(g1,g2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 8.6",
   "language": "",
   "name": "sagemath"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.15"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}