{ "cells": [ { "cell_type": "code", "execution_count": 121, "metadata": { "collapsed": false }, "outputs": [], "source": [ "class WeilObject(object):\n", " # Konstruktor: Gibt elliptische Kurve mit vorberechneter Weil-Funktion zurück.\n", " # @param n .. Natürliche Zahl aus dem Bereich {4,...,10,12}\n", " # self._E = Elliptische Kurve X_1(n) in den Variablen x, y und Parameter c\n", " # self._P = n-Torsionspunkt (0,0) auf E\n", " # self._n = natürliche Zahl n\n", " # self._f = Weil-Funktion E -> C mit Divisor div(f)=n{P} - n[OO]\n", " # self._s = Wohlformatierter String der Weil-Funktion f\n", " #\n", " def __init__(self, n):\n", " self._n = n\n", " x, y, c = var('x y c');\n", " \n", " if n == 4:\n", " self._E = EllipticCurve([1, c, c, 0, 0]); # X_1(4)\n", " self._P = self._E(0,0); # 4-Torsionspunkt (0,0)\n", " self._f = x^2 - y; # TODO: http://trac.sagemath.org/5930\n", " self._s = \"x^2 - y\";\n", " return;\n", " \n", " if n == 5:\n", " self._E = EllipticCurve([c+1, c, c, 0, 0]);\n", " self._P = self._E(0,0);\n", " self._f = -x^2 + x*y + y;\n", " self._s = \"-x^2 + x*y + y\";\n", " return;\n", "\n", " if n == 6:\n", " self._E = EllipticCurve([1-c, -c*(c+1), -c*(c+1), 0, 0]);\n", " self._P = self._E(0,0);\n", " self._f = y^2 - (c + 1)*x*y - (c + 1)^2*y + (c + 1)^2*x^2;\n", " self._s = \"y^2 - (c + 1)*x*y - (c + 1)^2*y + (c + 1)^2*x^2\";\n", " return;\n", " \n", " if n == 7:\n", " self._E = EllipticCurve([1-c-c^2, c^2*(c+1), c^2*(c+1), 0, 0]);\n", " self._P = self._E(0,0);\n", " self._f = (c - 1)*y^2 + x^2*y - c^3*x*y + c^4*y - c^4*x^2;\n", " self._s = \"(c - 1)*y^2 + x^2*y - c^3*x*y + c^4*y - c^4*x^2\";\n", " return;\n", "\n", " if n == 8:\n", " self._E = EllipticCurve([1-2*c^2, -c*(2*c+1)*(c+1)^2, -c*(2*c+1)*(c+1)^3, 0, 0]);\n", " self._P = self._E(0,0);\n", " self._f = (x*y^2 + (2*c + 3)*(c + 1)^3*y^2 - 2*(c + 1)^2*x^2*y - (2*c + 1)^2*(c + 1)^4*x*y\n", " - (2*c + 1)^2*(c + 1)^7*y + (2*c + 1)^2*(c + 1)^6*x^2);\n", " self._s = (\"x*y^2 \\n+ (2*c + 3)*(c + 1)^3*y^2 \\n- 2*(c + 1)^2*x^2*y \\n- (2*c + 1)^2*(c + 1)^4*x*y \\n\" +\n", " \"- (2*c + 1)^2*(c + 1)^7*y \\n+ (2*c + 1)^2*(c + 1)^6*x^2\");\n", " return;\n", "\n", " if n == 9:\n", " self._E = EllipticCurve([c^3 + c^2 + 1, c^2 * (c+1) * (c^2+c+1), \n", " c^2 * (c+1) * (c^2+c+1), 0, 0]);\n", " self._P = self._E(0,0);\n", " self._f = (y^3 + (c - 1)*(c^2 + c + 1)*x*y^2 + (c^2 + c + 1)^2*(c^3 + 2*c - 1)*y^2\n", " - (2*c - 1)*(c^2 + c + 1)^2*x^2*y + c^4*(c^2 + c + 1)^3*x*y\n", " + c^4*(c^2 + c + 1)^4*y - c^4*(c^2 + c + 1)^4*x^2);\n", " self._s = (\"y^3 \\n+ (c - 1)*(c^2 + c + 1)*x*y^2 \\n+ (c^2 + c + 1)^2*(c^3 + 2*c - 1)*y^2 \\n\" +\n", " \"- (2*c - 1)*(c^2 + c + 1)^2*x^2*y \\n+ c^4*(c^2 + c + 1)^3*x*y \\n\" +\n", " \"+ c^4*(c^2 + c + 1)^4*y \\n- c^4*(c^2 + c + 1)^4*x^2\");\n", " return;\n", "\n", " if n == 10:\n", " self._E = EllipticCurve([-c^3 - 2*c^2 + 4*c + 4, (c+1) * (c+2) * c^3, \n", " (c+1) * (c+2) * (c^2+6*c+4) * c^3, 0, 0]);\n", " self._P = self._E(0,0);\n", " self._f = (2*(c^2 - 2*c - 2)*y^3 + x^2*y^2 - (2*c + 1)*c^4*x*y^2\n", " + (c^3 + 16*c^2 + 22*c + 8)*c^6*y^2 - (3*c + 2)*c^6*x^2*y\n", " + (c + 1)^2*c^10*x*y - (c + 1)^2*(c^2 + 6*c + 4)*c^12*y\n", " + (c + 1)^2*c^12*x^2);\n", " self._s = (\"2*(c^2 - 2*c - 2)*y^3 \\n+ x^2*y^2 \\n- (2*c + 1)*c^4*x*y^2 \\n\" +\n", " \"+ (c^3 + 16*c^2 \\n+ 22*c + 8)*c^6*y^2 \\n- (3*c + 2)*c^6*x^2*y \\n\" +\n", " \"+ (c + 1)^2*c^10*x*y \\n- (c + 1)^2*(c^2 + 6*c + 4)*c^12*y \\n\" +\n", " \"+ (c + 1)^2*c^12*x^2\");\n", " return;\n", "\n", " if n == 12:\n", " self._E = EllipticCurve([6*c^4 - 8*c^3 + 2*c^2 + 2*c - 1, \n", " -c * (c-1)^2 * (2*c-1) * (2*c^2-2*c+1) * (3*c^2-3*c+1),\n", " -c * (c-1)^5 * (2*c-1) * (2*c^2-2*c+1) * (3*c^2-3*c+1),\n", " 0, 0]);\n", " self._P = self._E(0,0);\n", " self._f = (y^4\n", " + (6*c^2 - 8*c + 3)*(2*c^2 - 2*c + 1)*x*y^3\n", " + (c - 1)^4*(2*c^2 - 2*c + 1)^2*(36*c^4 - 90*c^3 + 96*c^2 - 49*c + 10)*y^3\n", " + 3*(c - 1)^2*(5*c^2 - 6*c + 2)*(2*c^2 - 2*c + 1)^2*x^2*y^2\n", " + (c - 1)^4*(14*c^2 - 16*c + 5)*(3*c^2 - 3*c + 1)^2*(2*c^2 - 2*c + 1)^3*x*y^2\n", " + (c - 1)^8*(3*c^2 - 3*c + 1)^2*(2*c^2 - 2*c + 1)^4*(12*c^4 - 42*c^3 + 57*c^2 - 33*c + 7)*y^2\n", " + 2*(c - 1)^6*(9*c^2 - 10*c + 3)*(3*c^2 - 3*c + 1)^2*(2*c^2 - 2*c + 1)^4*x^2*y\n", " + (2*c - 1)^2*(c - 1)^8*(3*c^2 - 3*c + 1)^4*(2*c^2 - 2*c + 1)^5*x*y\n", " - (2*c - 1)^2*(c - 1)^13*(3*c^2 - 3*c + 1)^4*(2*c^2 - 2*c + 1)^6*y\n", " + (2*c - 1)^2*(c - 1)^10*(3*c^2 - 3*c + 1)^4*(2*c^2 - 2*c + 1)^6*x^2);\n", " self._s = (\"y^4 \\n\" +\n", " \"+ (6*c^2 - 8*c + 3)*(2*c^2 - 2*c + 1)*x*y^3 \\n\" +\n", " \"+ (c - 1)^4*(2*c^2 - 2*c + 1)^2*(36*c^4 - 90*c^3 + 96*c^2 - 49*c + 10)*y^3 \\n\" +\n", " \"+ 3*(c - 1)^2*(5*c^2 - 6*c + 2)*(2*c^2 - 2*c + 1)^2*x^2*y^2 \\n\" +\n", " \"+ (c - 1)^4*(14*c^2 - 16*c + 5)*(3*c^2 - 3*c + 1)^2*(2*c^2 - 2*c + 1)^3*x*y^2 \\n\" +\n", " \"+ (c - 1)^8*(3*c^2 - 3*c + 1)^2*(2*c^2 - 2*c + 1)^4*(12*c^4 - 42*c^3 + 57*c^2 - 33*c + 7)*y^2 \\n\" +\n", " \"+ 2*(c - 1)^6*(9*c^2 - 10*c + 3)*(3*c^2 - 3*c + 1)^2*(2*c^2 - 2*c + 1)^4*x^2*y \\n\" +\n", " \"+ (2*c - 1)^2*(c - 1)^8*(3*c^2 - 3*c + 1)^4*(2*c^2 - 2*c + 1)^5*x*y \\n\" +\n", " \"- (2*c - 1)^2*(c - 1)^13*(3*c^2 - 3*c + 1)^4*(2*c^2 - 2*c + 1)^6*y \\n\" +\n", " \"+ (2*c - 1)^2*(c - 1)^10*(3*c^2 - 3*c + 1)^4*(2*c^2 - 2*c + 1)^6*x^2\");\n", " return;\n", " \n", " # Assert im Konstruktor erzeugen, damit kein Objekt für nicht unterstützte n-Werte erzeugt wird\n", " assert false;\n", " \n", " # Gibt die String-Repräsentation des Objektes zurück.\n", " def __repr__(self):\n", " str = \"Weil function info object of elliptic curve with %i-torsion point (0,0):\\n- \"%(self._n)\n", " str = str + self._E.__repr__() + \"\\n\"\n", " str = str + \"- %i-torsion point P = \"%(self._n) + self._P.__repr__() + \"\\n\"\n", " str = str + \"- Weil function f(x,y) = \" + self._f.__repr__() + \"\\n\"\n", " return str\n", " \n", " # Gibt die Weil-Funktion als wohlformatierten String auf der Konsole aus.\n", " def PrintNiceString(self):\n", " print self._s\n", " \n", " # Gibt den Wert der Weilfunktion f im Punkt (xval, yval) zurück.\n", " def f(self, xval, yval):\n", " return self._f(x = xval, y = yval)\n", " \n", " # Gibt f(-(x,y)) zurück, wobei mit -(x,y) das Inverse bzgl. der Addition auf der elliptischen Kurve gemeint ist.\n", " def GetMinusWeilFunction(self):\n", " a1, a2, a3, a4, a6 = self._E.a_invariants();\n", " return self.f(x, -a1 * x - a3 - y);\n", "\n", " # Reduziert den Ausdruck im Quotienten-Ring modulo der Gleichung der elliptischen Kurve.\n", " # @param term .. Element aus Q(x,y,c) = Quotient aus zwei Polynomen mit\n", " # Variablen x und y und Parameter c\n", " # @returns term modulo Q(x,y,c)/I, wobei Q(x,y,c) = Polynom/Polynom ist \n", " # und I das Ideal der Gleichung der elliptischen Kurve E\n", " def ReduceInQuotientRing(self, term):\n", " # Koeffizienten der elliptischen Kurve bestimmen\n", " a1, a2, a3, a4, a6 = self._E.a_invariants();\n", " # Polynomring, Ideal und Quotientenring bestimmen\n", " Q. = QQ['x','y','c'];\n", " J = Q.ideal(y^2 + a1*x*y + a3*y - x^3 - a2*x^2 - a4*x - a6);\n", " R. = QuotientRing(Q, J);\n", " # In werden die Variablen ersetzt gemäß x -> X, y -> Y, c -> C.\n", " termXYZ = sage_eval(str(term), locals={'x':X, 'y':Y, 'c':C});\n", " # Reduzieren im Ideal modulo der Gleichung der elliptischen Kurve.\n", " reduc = R(termXYZ).reduce(J.gens());\n", " # Rückbenennung X -> x, Y -> y, C -> c\n", " return sage_eval(str(reduc), locals={'X':x, 'Y':y, 'C':c}); \n", " \n", " # Statische Methode.\n", " # Ersetzt in alle Vorkommnisse von durch ().\n", " # Unterstützt werden die Variablen x, y und c.\n", " def ReplaceSubTerm(self, term, search, repl):\n", " term_str = str(term);\n", " term_repl = term_str.replace(str(search), \"(\" + str(repl) + \")\");\n", " return sage_eval(term_repl, locals={'x':x, 'y':y, 'c':c}) \n", " \n", " # Drückt x^m für m >= 3 durch kleinere x-Potenzen aus.\n", " def XPower(self, m):\n", " a1, a2, a3, a4, a6 = self._E.a_invariants();\n", " if m <= 2:\n", " return x^m\n", " if m >= 3:\n", " return expand(x^(m-3) * (y^2 + a1*x*y + a3*y - a2*x^2 - a4*x - a6))\n", " \n", " # Ersetzt in alle Vorkommnisse von durch Potenzen x^2, x, 1.\n", " # Parameter term: Polynom in x, y und c.\n", " def ReplaceHighXMonoms(self, term):\n", " # Polynomring als \"Monom-Wörterbuch\"\n", " R. = sage.rings.polynomial.multi_polynomial_ring.MPolynomialRing_polydict(QQ, 3, order='lex');\n", " \n", " # term als Polynom in x, y, c interpreteren\n", " poly = R(term)\n", " \n", " # Maximalgrad von x im Polynom bestimmen.\n", " max_x_degree = poly.degree(x)\n", " \n", " # Iteration von max_x_degree bis 3\n", " for e in range(max_x_degree, 2, -1): \n", " # Ersetze x^e durch kleinere Potenzen.\n", " term = self.ReplaceSubTerm(term, x^e, self.XPower(e)); \n", " # Ausmultiplizieren, damit im nächsten Schritt x^(e-1) ersetzt werden kann.\n", " term = expand(term); \n", " return term;\n", " \n", " # Drückt y^m für m >= 2 durch kleinere y-Potenzen aus, FALLS x = 1 gilt.\n", " def YPowerX1(self, m):\n", " a1, a2, a3, a4, a6 = self._E.a_invariants();\n", " if m < 2:\n", " return y^m\n", " if m >= 2:\n", " return expand(y^(m-2) * ((-a1-a3)* y + a2 + a4 + a6 + 1)) # Achtung: x ist 1 !!! \n", " \n", " # Ersetzt in alle Vorkommnisse von durch Potenzen y^2, y, 1,\n", " # FALLS x = 1 gilt.\n", " # Parameter term: Polynom in y und c.\n", " def ReplaceHighYMonomsX1(self, term):\n", " # Polynomring als \"Monom-Wörterbuch\"\n", " R. = sage.rings.polynomial.multi_polynomial_ring.MPolynomialRing_polydict(QQ, 3, order='lex');\n", " \n", " # term als Polynom in y, c interpreteren\n", " poly = R(term)\n", " \n", " # Maximalgrad von y im Polynom bestimmen.\n", " max_y_degree = poly.degree(y)\n", " \n", " # Iteration von max_y_degree bis 2\n", " for e in range(max_y_degree, 1, -1): \n", " # Ersetze y^e durch kleinere Potenzen.\n", " term = self.ReplaceSubTerm(term, y^e, self.YPowerX1(e)); \n", " # Ausmultiplizieren, damit im nächsten Schritt y^(e-1) ersetzt werden kann.\n", " term = expand(term); \n", " return term;\n", " \n", " # Bestimmt die komplexe Zahl alpha^n, die die folgende Eigenschaft besitzt:\n", " # f((x,y)) * f(-(x,y)) = alpha^n * x^n, wobei mit -(x,y) das Inverse bzgl. der Addition\n", " # auf der elliptischen Kurve gemeint ist.\n", " # -> Wenn alles korrekt ist, muss alpha = (-1)^n gelten.\n", " #\n", " def CalcAlpha(self):\n", " # Minus-Funktion g der Weil-Funktion berechnen, g(Q) = f(-Q) für alle Q in E\n", " g(x, y) = self.GetMinusWeilFunction()\n", " \n", " # Produkt aus Weil-Funktion und Minus-Funktion bestimmen\n", " p(x, y) = self.f(x, y) * g(x, y)\n", " \n", " # alpha^n berechnen, indem spezieller Wert x = 1 eingesetzt wird\n", " alpha = p(x = 1, y = y)\n", "\n", " # ausmultiplizieren\n", " alpha = expand(alpha)\n", " \n", " # Potenzen y^m für m >= 2 ersetzen durch Linearkombination aus 1 und y;\n", " # die Variable y muss aus dem Term komplett herausfallen, so dass eine Konstante übrig bleibt\n", " alpha = self.ReplaceHighYMonomsX1(alpha)\n", " \n", " # Den Wert (alpha^n) zurückgeben.\n", " return alpha\n", " \n", " # Überprüft die Integrität des Objektes:\n", " # - Ist P wirklich ein n-Torsionspunkt von E ?\n", " # - Gilt wirklich alpha = (-1)^n ?\n", " # - Gilt wirklich f((x,y)) * f(-(x,y)) = alpha * x^n ? \n", " def Verify(self):\n", " # 1. Ist P wirklich ein n-Torsionspunkt von E ?\n", " assert (self._n * self._P == self._E(0));\n", " \n", " # 2. Gilt wirklich alpha^n = (-1)^n ?\n", " alpha = self.CalcAlpha();\n", " assert (alpha == (-1)^self._n);\n", " \n", " # 3. Gilt wirklich f((x,y)) * f(-(x,y)) = alpha^n * x^n ?\n", " \n", " # 3a. Minus-Funktion g der Weil-Funktion berechnen, g(Q) = f(-Q) für alle Q in E\n", " g(x, y) = self.GetMinusWeilFunction()\n", " \n", " # 3b. Produkt aus Weil-Funktion und Minus-Funktion bestimmen\n", " p(x, y) = self.f(x, y) * g(x, y) \n", " \n", " #3c. Differenz p - (alpha^n * x^n) im Ideal reduzieren\n", " diff = p(x,y)-(alpha * x^self._n)\n", " diff = self.ReduceInQuotientRing(diff)\n", " \n", " #3d. Differenz muss Null sein\n", " assert (diff == 0);\n", " \n", " # 4. Alles OK.\n", " return true; \n" ] }, { "cell_type": "code", "execution_count": 122, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "Weil function info object of elliptic curve with 5-torsion point (0,0):\n", "- Elliptic Curve defined by y^2 + (c+1)*x*y + c*y = x^3 + c*x^2 over Symbolic Ring\n", "- 5-torsion point P = (0 : 0 : 1)\n", "- Weil function f(x,y) = -x^2 + x*y + y\n" ] }, "execution_count": 122, "metadata": {}, "output_type": "execute_result" } ], "source": [ "WeilObject(5)" ] }, { "cell_type": "code", "execution_count": 123, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "(c + 1)^2*x^2 - ((c + 1)*c + (c - 1)*x - y)*(c + 1)^2 - ((c + 1)*c + (c - 1)*x - y)*(c + 1)*x + ((c + 1)*c + (c - 1)*x - y)^2" ] }, "execution_count": 123, "metadata": {}, "output_type": "execute_result" } ], "source": [ "WeilObject(6).GetMinusWeilFunction()" ] }, { "cell_type": "code", "execution_count": 124, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "-1" ] }, "execution_count": 124, "metadata": {}, "output_type": "execute_result" } ], "source": [ "WeilObject(7).CalcAlpha()" ] }, { "cell_type": "code", "execution_count": 125, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "x*y^2 \n", "+ (2*c + 3)*(c + 1)^3*y^2 \n", "- 2*(c + 1)^2*x^2*y \n", "- (2*c + 1)^2*(c + 1)^4*x*y \n", "- (2*c + 1)^2*(c + 1)^7*y \n", "+ (2*c + 1)^2*(c + 1)^6*x^2\n" ] } ], "source": [ "WeilObject(8).PrintNiceString()" ] }, { "cell_type": "code", "execution_count": 126, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "True\n", "CPU times: user 9.63 s, sys: 0 ns, total: 9.63 s\n", "Wall time: 12 s\n" ] } ], "source": [ "%time result = WeilObject(9).Verify(); print result" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": true }, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 7.0", "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.10" } }, "nbformat": 4, "nbformat_minor": 0 }