{ "cells": [ { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "# For an elliptic curve, with notation as in the paper,\n", "# this block of code computes the preimage of a divisor \n", "# under the map Frac(\\overline{R}) \\rightarrow \\Div^0(\\overline}{E})\n", "# if such a preimage exists.\n", "\n", "def lineBetweenPoints(P, Q):\n", " a1, b1 = P.xy()\n", " a2, b2 = Q.xy()\n", " if (a1 == a2):\n", " return (x-a1)\n", " return ((b2-b1)/(a2-a1))*(x-a1) - (y-b1)\n", "\n", "def lineTangentTo(P, E):\n", " a, b = P.xy()\n", " if (P.order() == 2):\n", " return (x-a)\n", " defining_array = E.a_invariants()\n", " a1 = defining_array[0]\n", " a2 = defining_array[1]\n", " a3 = defining_array[2]\n", " a4 = defining_array[3]\n", " a6 = defining_array[4]\n", " slope = (3*a^2 + 2*a2*a - a1*b + a4)/(2*b + a1*a+a3)\n", " return slope*(x-a) - (y-b)\n", "\n", "class Divisor:\n", " def __init__(self, E, dict, base_point):\n", " self.E = E\n", " self.base_point = base_point\n", " self.dict = {}\n", " for key, value in dict.items():\n", " if value != 0:\n", " self.dict[key] = value\n", "\n", " def __add__(self, another_divisor):\n", " new_dict = {}\n", " for key, value in self.dict.items():\n", " if value is not 0:\n", " new_dict[key] = value\n", "\n", " for key, value in another_divisor.dict.items():\n", " if key in new_dict.keys() and new_dict[key] + value is not 0:\n", " new_dict[key] = new_dict[key] + value\n", " elif value is not 0:\n", " new_dict[key] = value\n", " return Divisor(self.E, new_dict, base_point)\n", "\n", " def __sub__(self, another_divisor):\n", " new_dict = {}\n", " for key, value in self.dict.items():\n", " if value is not 0:\n", " new_dict[key] = value\n", "\n", " for key, value in another_divisor.dict.items():\n", " if key in new_dict.keys() and new_dict[key] - value is not 0:\n", " new_dict[key] = new_dict[key] - value\n", " elif value is not 0:\n", " new_dict[key] = -value\n", " return Divisor(self.E, new_dict, base_point)\n", "\n", " def __eq__(self, another_divisor):\n", " if (self.E is not another_divisor.E):\n", " return False\n", " elif self.dict == another_divisor.dict:\n", " return True\n", " else:\n", " return False\n", "\n", " def __ne__(self, another_divisor):\n", " return not (self == another_divisor)\n", "\n", " def coeff(self, point):\n", " return self.dict[point]\n", "\n", " def support(self):\n", " return self.dict.keys()\n", "\n", " def getMNPQ(self):\n", " m = 0\n", " n = 0\n", " P = None\n", " Q = None\n", " for point in self.support():\n", " if (self.coeff(point) > 0 and point != base_point):\n", " P = point\n", " m = self.coeff(point)\n", " elif (self.coeff(point) < 0 and point != base_point):\n", " Q = point\n", " n = -self.coeff(point)\n", " return m, n, P, Q\n", "\n", " def getPair(self, comparisonOperator, equalityOperator):\n", " points = set()\n", " for point in self.support():\n", " if point != base_point and (comparisonOperator(self.coeff(point), 0)):\n", " points.add(point)\n", " for P in points:\n", " for Q in points:\n", " if (P != Q and equalityOperator(P,-Q)):\n", " return (P,Q)\n", "\n", " def existsPair(self, comparisonOperator, equalityOperator):\n", " points = set()\n", " for point in self.support():\n", " if point != base_point and (comparisonOperator(self.coeff(point), 0)):\n", " points.add(point)\n", " for P in points:\n", " for Q in points:\n", " if (P != Q and equalityOperator(P,-Q)):\n", " return True\n", " return False\n", "\n", "def interpolate(D, base_point, E):\n", " # start with rational function 1\n", " # while D is not 0\n", " rational_function = 1\n", " E = D.E\n", " zero_divisor = Divisor(E,{},base_point)\n", " while (D != zero_divisor):\n", " #print \"loop\"\n", " #print D.dict\n", " #print rational_function\n", " # if P, Q have coefficients > 0 and P != -Q, subtract (

+ + <-(P+Q)> - 3<0>) and multiply rational function by line between P and Q\n", " if D.existsPair(operator.gt, operator.ne):\n", " #print \"case 1\"\n", " P, Q = D.getPair(operator.gt, operator.ne)\n", " sub_divisor = Divisor(E, {P : 1, Q : 1, -(P+Q) : 1, base_point: -3}, base_point)\n", " D = D - sub_divisor\n", " rational_function = rational_function*lineBetweenPoints(P,Q)\n", " # else if P, Q have coefficients > 0 and P=-Q subtract (

+ - 2<0>) and multiply rational function by line between P and Q\n", " elif D.existsPair(operator.gt, operator.eq):\n", " #print \"case 2\"\n", " P, Q = D.getPair(operator.gt, operator.eq)\n", " sub_divisor = Divisor(E, {P:1,Q:1, base_point: -2}, base_point)\n", " D = D - sub_divisor\n", " rational_function = rational_function*lineBetweenPoints(P,Q)\n", " # else if P, Q have coefficients < 0 and P != -Q, add (

+ + <-(P+Q)> - 3<0>) and divide rational function by line between P and Q\n", " elif D.existsPair(operator.lt, operator.ne):\n", " #print \"case 3\"\n", " P, Q = D.getPair(operator.lt, operator.ne)\n", " add_divisor = Divisor(E, {P:1,Q:1,-(P+Q):1,base_point:-3}, base_point)\n", " D = D+add_divisor\n", " rational_function = rational_function/lineBetweenPoints(P,Q)\n", " # else if P, Q have coefficients < 0 and P=-Q add (

+ - 2<0>) and divide rational function by line between P and Q\n", " elif D.existsPair(operator.lt, operator.eq):\n", " #print \"case 4\"\n", " P, Q = D.getPair(operator.lt, operator.eq)\n", " add_divisor = Divisor(E, {P:1,Q:1,base_point:-2}, base_point)\n", " D = D+add_divisor\n", " rational_function = rational_function/lineBetweenPoints(P,Q)\n", " # else if D is of form m

- n + o<0> (o may be positive or negative, m >= 0, n>= 0)\n", " else:\n", " m, n, P, Q = D.getMNPQ()\n", " # if m >= 2 and P not order 2, then subtract (2

+ <-2P> - 3<0>) and multiply rational function by tangent line at P\n", " if m >= 2 and P.order() != 2:\n", " #print \"case 5\"\n", " if P.order() == 3:\n", " # -2P = P\n", " sub_divisor = Divisor(E, {P:3, base_point:-3}, base_point)\n", " else:\n", " sub_divisor = Divisor(E, {P:2, (-2*P):1, base_point:-3}, base_point)\n", " D = D - sub_divisor\n", " rational_function = rational_function*lineTangentTo(P, E)\n", " # else if m >= 2 and P of order 2, then subtract (2

- 2<0>) and multiply rational function by tangent line at P\n", " elif m >= 2 and P.order() == 2:\n", " #print \"case 6\"\n", " sub_divisor = Divisor(E, {P:2,base_point:-2}, base_point)\n", " D = D - sub_divisor\n", " rational_function = rational_function*lineTangentTo(P, E)\n", " # else if n >=2 and Q not of order 2, then add (2 + <-2Q> - 3<0>) and divide rational function by tangent line at Q\n", " elif n >=2 and Q.order() != 2:\n", " #print \"case 7\"\n", " if Q.order() == 3:\n", " # -2Q = Q\n", " add_divisor = Divisor(E, {Q:3, base_point: -3}, base_point)\n", " else:\n", " add_divisor = Divisor(E, {Q:2, (-2*Q):1, base_point: -3}, base_point)\n", " D = D + add_divisor\n", " rational_function = rational_function/lineTangentTo(Q, E)\n", " # else if n >= 2 and Q of order 2, then add (2 - 2<0>) and divide rational function by tangent line at Q\n", " elif n >=2 and Q.order() == 2:\n", " #print \"case 8\"\n", " add_divisor = Divisor(E, {Q:2, base_point:-2}, base_point)\n", " D = D + add_divisor\n", " rational_function = rational_function/lineTangentTo(Q, E)\n", " # else if m=n=1, add ( + <-Q> - 2<0>) and divide rational function by line through Q and -Q\n", " elif m==1 and n==1:\n", " #print \"case 9\"\n", " add_divisor = Divisor(E, {Q:1, -Q:1, base_point:-2}, base_point)\n", " D = D + add_divisor\n", " rational_function = rational_function/lineBetweenPoints(Q,-Q)\n", " else:\n", " print \"cannot interpolate the following divisor\"\n", " print D.dict\n", " print rational_function\n", " return None\n", " return rational_function" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [ "# Example 1: computing the units of the Fermat 3-curve\n", "from sage.modules.free_module_integer import IntegerLattice\n", "from itertools import product\n", "\n", "# construct the curve\n", "a = PolynomialRing(RationalField(), 'a').gen()\n", "K. = NumberField(a^2 + a + 1)\n", "R. = PolynomialRing(K,3)\n", "mor = EllipticCurve_from_cubic(x^3+y^3-z^3, [1,-1,0], morphism=True)\n", "E = EllipticCurve_from_cubic(x^3+y^3-z^3, [1,-1,0], morphism=False)\n", "\n", "# construct the boundary points\n", "\n", "T1 = E(mor([0,1,1]))\n", "T2 = E(mor([0,b,1]))\n", "T3 = E(mor([0,b^2,1]))\n", "T4 = E(mor([1,0,1]))\n", "T5 = E(mor([b,0,1]))\n", "T6 = E(mor([b^2,0,1]))\n", "T7 = E(mor([-1,1,0]))\n", "T8 = E(mor([-b,1,0]))\n", "T9 = E(mor([-b^2,1,0]))\n", "\n", "# ommitted -- we compute the lattice of relations, and get \n", "# the following matrix using the algorithm described\n", "# in the paper\n", "\n", "#[ 1 0 0 0 0 1 0 0 -2] 1\n", "#[ 0 1 0 0 0 1 0 1 -3] 2\n", "#[ 0 0 1 0 0 1 0 2 -4] 3\n", "#[ 0 0 0 1 0 2 0 2 -5] 4\n", "#[ 0 0 0 0 1 2 0 1 -4] 5\n", "#[ 0 0 0 0 0 3 0 0 -3] 6\n", "#[ 0 0 0 0 0 0 1 1 -2] 7\n", "#[ 0 0 0 0 0 0 0 3 -3] 8\n", "\n", "base_point = E([0,1,0])\n", "\n", "#interpolate each divisor\n", "print interpolate(Divisor(E, {T1:1, T6:1, T9:-2}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T2:1, T6:1, T8:1, T9:-3}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T3:1, T6:1, T8:2, T9:-4}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T4:1, T6:2, T8:2, T9:-5}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T5:1, T6:2, T8:1, T9:-4}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T6:3, T9:-3}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T7:1, T8:1, T9:-2}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T8:3, T9:-3}, base_point), base_point, E)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ 0.570850711738291 -0.570850711738291]\n", "[-0.570850711738291 0.570850711738291]\n" ] } ], "source": [ "# Example 2: computing the units of an elliptic curve\n", "# with defining equation y^2 = x^3 - x + root(2)\n", "from sage.modules.free_module_integer import IntegerLattice\n", "from itertools import product\n", "\n", "# construct an appropriate number field\n", "a = PolynomialRing(RationalField(), 'a').gen()\n", "M. = NumberField(a^4 - 2)\n", "L. = M.extension(a^2 + 1)\n", "K. = L.absolute_field()\n", "\n", "# note that taking the absolute_field creates an injection of fields\n", "# so we are now working in K\n", "fm = K.structure()[1]\n", "i = fm(prei)\n", "root4of2 = fm(preroot4of2)\n", "root2 = root4of2^2\n", "\n", "# construct the elliptic curve\n", "R. = PolynomialRing(K,3)\n", "E = EllipticCurve([0,0,0,-1, root2])\n", "\n", "# construct the points\n", "T1 = E([0,1,0])\n", "T2 = E([-root2,0,1])\n", "T3 = E([(1+i)/root2,0,1])\n", "T4 = E([(1-i)/root2,0,1])\n", "\n", "Q1 = E([0,root4of2,1])\n", "Q2 = E([0,-root4of2,1])\n", "\n", "# we compute that\n", "# T1.order() = 1\n", "# T2.order() = 2\n", "# T3.order() = 2\n", "# T4.order() = 2\n", "# Q1.order() = infinity\n", "# Q2.order() = infinity\n", "\n", "print E.height_pairing_matrix([Q1,Q2])\n", "\n", "# the matrix will be\n", "# [ 0.570850711738291 -0.570850711738291]\n", "# [-0.570850711738291 0.570850711738291]\n", "# which indicates the relation Q1 + Q2 = 0" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Free module of degree 6 and rank 4 over Integer Ring\n", "Echelon basis matrix:\n", "[ 1 1 1 1 -2 -2]\n", "[ 0 2 0 0 -1 -1]\n", "[ 0 0 2 0 -1 -1]\n", "[ 0 0 0 2 -1 -1]\n", "[\n", "(1, 1, 1, 1, -2, -2),\n", "(0, 2, 0, 0, -1, -1),\n", "(0, 0, 2, 0, -1, -1),\n", "(0, 0, 0, 2, -1, -1)\n", "]\n" ] } ], "source": [ "# Example 2 continued \n", "\n", "# from before we have Q1 + Q2 = 0 \n", "# the only nontrivial relation between the torsion points are:\n", "# T2 + T3 + T4 = 0\n", "# T1 = 0\n", "# 2*T2 = 0\n", "# 2*T3 = 0\n", "# 2*T4 = 0\n", "\n", "\n", "L = IntegerLattice([[1,0,0,0,0,0],\n", " [0,2,0,0,0,0],\n", " [0,0,2,0,0,0],\n", " [0,0,0,2,0,0],\n", " [0,1,1,1,0,0],\n", " [0,0,0,0,1,1]])\n", "\n", "H = IntegerLattice([[1,-1,0,0,0,0],\n", " [0,1,-1,0,0,0],\n", " [0,0,1,-1,0,0],\n", " [0,0,0,1,-1,0],\n", " [0,0,0,0,1,-1]])\n", "\n", "# print our basis of relations\n", "print H.intersection(L).basis()\n", "\n" ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "(-y)/x^2\n", "(x + (-1/12*d^6 - 5/12*d^4 - 1/12*d^2 - 17/12))/x\n", "(x + (-1/8*d^7 + 1/24*d^6 - 1/2*d^5 + 5/24*d^4 - 3/8*d^3 + 1/24*d^2 - 15/4*d + 17/24))/x\n", "(x + (1/8*d^7 + 1/24*d^6 + 1/2*d^5 + 5/24*d^4 + 3/8*d^3 + 1/24*d^2 + 15/4*d + 17/24))/x\n" ] } ], "source": [ "# Example 2 continued\n", "# interpolate\n", "# [ 1 1 1 1 -2 -2]\n", "# [ 0 2 0 0 -1 -1]\n", "# [ 0 0 2 0 -1 -1]\n", "# [ 0 0 0 2 -1 -1]\n", "\n", "base_point = E([0,1,0])\n", "print interpolate(Divisor(E, {T1:1, T2:1, T3:1, T4:1, Q1:-2, Q2:-2}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T2:2, Q1:-1, Q2:-1}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T3:2, Q1:-1, Q2:-1}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T4:2, Q1:-1, Q2:-1}, base_point), base_point, E)" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "preroot4of2^2\n", "-1/2*preroot4of2^2*prei - 1/2*preroot4of2^2\n", "1/2*preroot4of2^2*prei - 1/2*preroot4of2^2\n" ] } ], "source": [ "# Example 2 continued\n", "# translate the results back from coefficients in K \n", "# to coefficients in L to get our final unit group generators\n", "revfm = K.structure()[0]\n", "print revfm(-1/12*d^6 - 5/12*d^4 - 1/12*d^2 - 17/12)\n", "print revfm(-1/8*d^7 + 1/24*d^6 - 1/2*d^5 + 5/24*d^4 - 3/8*d^3 + 1/24*d^2 - 15/4*d + 17/24)\n", "print revfm(1/8*d^7 + 1/24*d^6 + 1/2*d^5 + 5/24*d^4 + 3/8*d^3 + 1/24*d^2 + 15/4*d + 17/24)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "438976/225\n", "\n", "torsion orders:\n", "1\n", "2\n", "2\n", "2\n", "+Infinity\n", "+Infinity\n", "\n", "the final units:\n", "(-y)/x^2\n", "(x - 1)/x\n", "(x + 1)/x\n", "(x - 4)/x\n" ] } ], "source": [ "# Example 3: computing units of an elliptic curve\n", "# defined by equation y^2 = (x-1)(x+1)(x-4)\n", "\n", "# construct the curve\n", "E = EllipticCurve([0,-4,0,-1,4])\n", "# note that the j-invariant is negative under the 5-adic valuation and\n", "# thus the intrinsic tropicalization will be interesting\n", "print E.j_invariant()\n", "# construct the poitns\n", "T1 = E([0,1,0])\n", "T2 = E([1,0,1])\n", "T3 = E([-1,0,1])\n", "T4 = E([4,0,1])\n", "\n", "Q1 = E([0,2,1])\n", "Q2 = E([0,-2,1])\n", "\n", "# compute their orders\n", "print \"\\ntorsion orders:\"\n", "print T1.order()\n", "print T2.order()\n", "print T3.order()\n", "print T4.order()\n", "print Q1.order()\n", "print Q2.order()\n", "\n", "# the lattice of relations is\n", "# interpolate\n", "# [ 1 1 1 1 -2 -2]\n", "# [ 0 2 0 0 -1 -1]\n", "# [ 0 0 2 0 -1 -1]\n", "# [ 0 0 0 2 -1 -1]\n", "\n", "base_point = E([0,1,0])\n", "R. = PolynomialRing(QQ,3)\n", "\n", "# interpolate our divisors\n", "print \"\\nthe final units:\"\n", "print interpolate(Divisor(E, {T1:1, T2:1, T3:1, T4:1, Q1:-2, Q2:-2}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T2:2, Q1:-1, Q2:-1}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T3:2, Q1:-1, Q2:-1}, base_point), base_point, E)\n", "print interpolate(Divisor(E, {T4:2, Q1:-1, Q2:-1}, base_point), base_point, E)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 8.2", "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.14" } }, "nbformat": 4, "nbformat_minor": 2 }