Path: blob/master/sage/schemes/hyperelliptic_curves/hyperelliptic_finite_field.py
4156 views
r"""1Hyperelliptic curves over a finite field23EXAMPLES::45sage: K.<a> = GF(9, 'a')6sage: x = polygen(K)7sage: C = HyperellipticCurve(x^7 - x^5 - 2, x^2 + a)8sage: C._points_fast_sqrt()9[(0 : 1 : 0), (2*a : 2*a + 2 : 1), (2*a : 2*a : 1), (a + 1 : a : 1), (a + 1 : a + 1 : 1), (2 : a + 1 : 1), (1 : a + 1 : 1)]10"""1112#*****************************************************************************13# Copyright (C) 2006 David Kohel <[email protected]>14# Copyright (C) 2007 Robert Bradshaw <[email protected]>15# Copyright (C) 2010 Alyson Deines <[email protected]>, Marina Gresham16# <[email protected]>, Gagan Sekhon <[email protected]>17# Distributed under the terms of the GNU General Public License (GPL)18# http://www.gnu.org/licenses/19#*****************************************************************************2021from sage.rings.all import ZZ, RR, binomial22import hyperelliptic_generic23from sage.schemes.hyperelliptic_curves.hypellfrob import hypellfrob24from sage.misc.cachefunc import cached_method25from sage.matrix.constructor import identity_matrix, matrix26from sage.misc.functional import rank272829class HyperellipticCurve_finite_field(hyperelliptic_generic.HyperellipticCurve_generic):3031def _frobenius_coefficient_bound(self):32"""33Computes bound on number of p-adic digits needed to recover34frobenius polynomial, i.e. returns B so that knowledge of35a_1, ..., a_g modulo p^B determine frobenius polynomial uniquely.3637TESTS::3839sage: R.<t> = PolynomialRing(GF(37))40sage: HyperellipticCurve(t^3 + t + 1)._frobenius_coefficient_bound()41142sage: HyperellipticCurve(t^5 + t + 1)._frobenius_coefficient_bound()43244sage: HyperellipticCurve(t^7 + t + 1)._frobenius_coefficient_bound()4534647sage: R.<t> = PolynomialRing(GF(next_prime(10^9)))48sage: HyperellipticCurve(t^3 + t + 1)._frobenius_coefficient_bound()49150sage: HyperellipticCurve(t^5 + t + 1)._frobenius_coefficient_bound()51252sage: HyperellipticCurve(t^7 + t + 1)._frobenius_coefficient_bound()53254sage: HyperellipticCurve(t^9 + t + 1)._frobenius_coefficient_bound()55356sage: HyperellipticCurve(t^11 + t + 1)._frobenius_coefficient_bound()57358sage: HyperellipticCurve(t^13 + t + 1)._frobenius_coefficient_bound()59460"""61assert self.base_ring().is_finite()62p = self.base_ring().characteristic()63q = self.base_ring().order()64sqrtq = RR(q).sqrt()65g = self.genus()6667# note: this bound is from Kedlaya's paper, but he tells me it's not68# the best possible69M = 2 * binomial(2*g, g) * sqrtq**g70B = ZZ(M.ceil()).exact_log(p)71if p**B < M:72B += 173return B747576def _frobenius_matrix(self, N=None):77"""78Compute p-adic frobenius matrix to precision p^N. If N not supplied,79a default value is selected, which is the minimum needed to recover80the charpoly unambiguously.8182Currently only implemented using hypellfrob, which means only works83over GF(p^1), and must have p > (2g+1)(2N-1).8485TESTS::8687sage: R.<t> = PolynomialRing(GF(37))88sage: H = HyperellipticCurve(t^5 + t + 2)89sage: H._frobenius_matrix()90[1258 + O(37^2) 925 + O(37^2) 132 + O(37^2) 587 + O(37^2)]91[1147 + O(37^2) 814 + O(37^2) 241 + O(37^2) 1011 + O(37^2)]92[1258 + O(37^2) 1184 + O(37^2) 1105 + O(37^2) 482 + O(37^2)]93[1073 + O(37^2) 999 + O(37^2) 772 + O(37^2) 929 + O(37^2)]9495"""96p = self.base_ring().characteristic()97f, h = self.hyperelliptic_polynomials()98if h != 0:99# need y^2 = f(x)100raise NotImplementedError, "only implemented for curves y^2 = f(x)"101102sign = 1103if not f.is_monic():104# at this time we need a monic f105c = f.leading_coefficient()106f = f / c107if c.is_square():108# solutions of $y^2 = c * f(x)$ correspond naturally to109# solutions of $(sqrt(c) y)^2 = f(x)$110pass111else:112# we'll count points on a twist and then correct by untwisting...113sign = -1114assert f.is_monic()115116if N is None:117N = self._frobenius_coefficient_bound()118119matrix_of_frobenius = hypellfrob(p, N, f)120matrix_of_frobenius = sign * matrix_of_frobenius121return matrix_of_frobenius122123124def frobenius_polynomial(self):125"""126Charpoly of frobenius, as an element of ZZ[x].127128TESTS::129130sage: R.<t> = PolynomialRing(GF(37))131sage: H = HyperellipticCurve(t^5 + t + 2)132sage: H.frobenius_polynomial()133x^4 + x^3 - 52*x^2 + 37*x + 1369134135A quadratic twist:136137::138139sage: H = HyperellipticCurve(2*t^5 + 2*t + 4)140sage: H.frobenius_polynomial()141x^4 - x^3 - 52*x^2 - 37*x + 1369142143"""144p = self.base_ring().characteristic()145g = self.genus()146N = self._frobenius_coefficient_bound()147# compute chapoly over ZZ and then reduce back148# (because charpoly of p-adic matrices sometimes loses precision)149M = self._frobenius_matrix(N=N).change_ring(ZZ)150151# get a_g, ..., a_0 in ZZ (i.e. with correct signs)152f = M.charpoly().list()[g:2*g+1]153ppow = p**N154f = [x % ppow for x in f]155f = [x if 2*x < ppow else x - ppow for x in f]156157# get a_{2g}, ..., a_{g+1}158f = [f[g-i] * p**(g-i) for i in range(g)] + f159160return ZZ['x'](f)161162163164def _points_fast_sqrt(self):165"""166Count points by enumerating over x and solving the resulting167quadratic for y.168169EXAMPLES::170171sage: K.<a> = GF(9, 'a')172sage: x = polygen(K)173sage: C = HyperellipticCurve(x^7 - 1, x^2 + a)174sage: C._points_fast_sqrt()175[(0 : 1 : 0), (2 : a + 1 : 1), (a : 2*a + 1 : 1), (2*a + 2 : 2*a : 1), (2*a + 2 : 1 : 1), (1 : 2*a + 2 : 1), (1 : 0 : 1)]176sage: K.<a> = GF(49, 'a')177sage: x = polygen(K)178sage: C = HyperellipticCurve(x^5 - x^2 - 1, x^2 + a)179sage: len(C._points_fast_sqrt())18031181182TESTS::183184sage: x = polygen(GF(16, 'a'))185sage: C = HyperellipticCurve(x^5 - x + 1, x^2 + x + 1)186sage: set(C._points_fast_sqrt()) == set(C._points_cache_sqrt())187True188sage: x = polygen(GF(19))189sage: C = HyperellipticCurve(x^5 + 5*x^2 + 1, x + 1)190sage: set(C._points_fast_sqrt()) == set(C._points_cache_sqrt())191True192sage: x = polygen(GF(13))193sage: C = HyperellipticCurve(x^3 + x^2 - 1)194sage: C._points_fast_sqrt()195[(0 : 1 : 0), (0 : 5 : 1), (0 : 8 : 1), (1 : 1 : 1), (1 : 12 : 1), (3 : 3 : 1), (3 : 10 : 1), (4 : 1 : 1), (4 : 12 : 1), (6 : 2 : 1), (6 : 11 : 1), (7 : 1 : 1), (7 : 12 : 1), (8 : 4 : 1), (8 : 9 : 1), (9 : 4 : 1), (9 : 9 : 1), (12 : 5 : 1), (12 : 8 : 1)]196sage: set(C._points_fast_sqrt()) == set(C._points_cache_sqrt())197True198"""199# For givaro finite fields, taking square roots is very fast200# so no need to cache as in prime case201K = self.base_ring()202f, h = self.hyperelliptic_polynomials()203one = K(1)204205# start with the points at infinity206P = self.defining_polynomial()207if not P(K(0), K(1), K(0)):208# (0:1:0) is a point on the curve209points = [self.point([K(0), K(1), K(0)], check=True)]210else:211points=[]212if P.degree() > 2:213# P(1, y, 0) = r*y + s214s = P(K(1), K(0), K(0))215r = P(K(1), K(1), K(0)) - s216if r: # r not zero217points.append(self.point([K(1), -s/r, K(0)], check=True))218# the case r = 0 need not be considered219elif K.characteristic() == 2: # deg(P) = 2 and char(K) = 2220# quadratic equation doesn't work in characteristic 2 so use brute221# force222points += [self.point([K(1), y, K(0)], check=True) for y in K \223if not P(K(1), y, K(0))]224else: # deg(P) = 2 and char(K) not 2225# P(1, y, 0) = y^2 + r*y + s226s = -f[2]227r = h[1]228d = r**2/4 - s229if not d: # d = 0230points.append(self.point([K(1), -r/2, K(0)], check=True))231elif d.is_square():232sqrtd = d.sqrt()233points.append(self.point([K(1), -r/2+sqrtd, K(0)], check=True))234points.append(self.point([K(1), -r/2-sqrtd, K(0)], check=True))235236if K.characteristic() == 2:237# quadratic equation doesn't work in characteristic 2238if h.is_zero():239for x in K:240points.append(self.point([x, f(x).sqrt(), one], check=True))241else:242R = f.base_ring()243a_sqrts = { } # Artin-Schreier 2-roots244for x in K:245a_sqrts[x**2 + x] = x # char 2 => x^2 - x == x^2 + x246for x in K:247b = h(x)248c = f(x)249if b:250try:251r = a_sqrts[c / b**2]252points.append(self.point([x, r*b, one], check=True))253points.append(self.point([x, r*b+b, one], check=True))254except KeyError:255# y^2 + by + c irreducible, so yields no points256pass257else: # b == 0258points.append(self.point([x, c.sqrt(), one], check=True))259elif h.is_zero():260# special case to save work if we are of the form y^2 = f(x)261for x in K:262y2 = f(x)263if not y2: # y = 0264points.append(self.point([x, y2, one], check=True))265elif y2.is_square():266y = y2.sqrt()267points.append(self.point([x, y, one], check=True))268points.append(self.point([x, -y, one], check=True))269else:270b = -h/2271D = b*b + f272for x in K:273Dval = D(x)274if not Dval: # D(x) = 0275points.append(self.point([x, b(x), one], check=True))276elif Dval.is_square():277sqrtD = Dval.sqrt()278v = b(x)279points.append(self.point([x, v+sqrtD, one], check=True))280points.append(self.point([x, v-sqrtD, one], check=True))281return points282283def _points_cache_sqrt(self, brute_force=False):284"""285Count points by enumerating over x and solving the resulting286quadratic for y.287288Caches all square roots ahead of time by squaring every element of289the field. Elements must have an __index__ method.290291EXAMPLES::292293sage: x = polygen(GF(7))294sage: C = HyperellipticCurve(x^3 + x^2 - 1)295sage: C._points_cache_sqrt()296[(0 : 1 : 0), (1 : 6 : 1), (1 : 1 : 1), (2 : 5 : 1), (2 : 2 : 1), (3 : 0 : 1), (4 : 4 : 1), (4 : 3 : 1), (5 : 4 : 1), (5 : 3 : 1)]297sage: set(C._points_cache_sqrt()) == set(C._points_cache_sqrt(brute_force=True))298True299"""300K = self.base_ring()301if K.characteristic() != 2:302# cache the squares (faster than O(p) sqrts)303square_roots = [None] * len(K)304for x in K:305square_roots[x*x] = x306f, h = self.hyperelliptic_polynomials()307one = K(1)308309# start with the points at infinity310P = self.defining_polynomial()311if not P(K(0), K(1), K(0)):312# (0:1:0) is a point on the curve313points = [self.point([K(0), K(1), K(0)], check=True)]314else:315points=[]316if P.degree() > 2:317# P(1, y, 0) = r*y + s318s = P(K(1), K(0), K(0))319r = P(K(1), K(1), K(0)) - s320if r: # r not zero321points.append(self.point([K(1), -s/r, K(0)], check=True))322# the case r = 0 need not be considered323elif K.characteristic() == 2: # deg(P) = 2 and char(K) = 2324# quadratic equation doesn't work in characteristic 2 so use brute325# force326points += [self.point([K(1), y, K(0)], check=True) for y in K \327if not P(K(1), y, K(0))]328else: # deg(P) = 2 and char(K) not 2329# P(1, y, 0) = y^2 + r*y + s330s = -f[2]331r = h[1]332d = r**2/4 - s333sqrtd = square_roots[d]334if not d: # d = 0335points.append(self.point([K(1), -r/2, K(0)], check=True))336elif sqrtd is not None:337points.append(self.point([K(1), -r/2+sqrtd, K(0)], check=True))338points.append(self.point([K(1), -r/2-sqrtd, K(0)], check=True))339340if K.characteristic() == 2 or brute_force:341# quadratic equation doesn't work in characteristic 2342# but there are only 4 affine points, so just test them343f = self.defining_polynomial()344points += [self.point([x, y, one], check=True) for x in K for y in K if not f(x, y, one)]345elif h.is_zero():346# special case to save work if we are of the form y^2 = f(x)347for x in K:348y2 = f(x)349y = square_roots[y2]350if not y2: # y = 0351points.append(self.point([x, y2, one], check=True))352elif y is not None:353points.append(self.point([x, y, one], check=True))354points.append(self.point([x, -y, one], check=True))355else:356b = -h/2357D = b*b + f # this is really disc/4358for x in K:359Dval = D(x)360sqrtD = square_roots[Dval]361if not Dval: # D(x) = 0362points.append(self.point([x, b(x), one], check=True))363elif sqrtD is not None:364v = b(x)365points.append(self.point([x, v+sqrtD, one], check=True))366points.append(self.point([x, v-sqrtD, one], check=True))367return points368369370def points(self):371r"""372All the points on this hyperelliptic curve.373374EXAMPLES::375376sage: x = polygen(GF(7))377sage: C = HyperellipticCurve(x^7 - x^2 - 1)378sage: C.points()379[(0 : 1 : 0), (2 : 5 : 1), (2 : 2 : 1), (3 : 0 : 1), (4 : 6 : 1), (4 : 1 : 1), (5 : 0 : 1), (6 : 5 : 1), (6 : 2 : 1)]380381::382383sage: x = polygen(GF(121, 'a'))384sage: C = HyperellipticCurve(x^5 + x - 1, x^2 + 2)385sage: len(C.points())386122387388Conics are allowed (the issue reported at #11800 has been resolved)::389390sage: R.<x> = GF(7)[]391sage: H = HyperellipticCurve(3*x^2 + 5*x + 1)392sage: H.points()393[(0 : 6 : 1), (0 : 1 : 1), (1 : 4 : 1), (1 : 3 : 1), (2 : 4 : 1), (2 : 3 : 1), (3 : 6 : 1), (3 : 1 : 1)]394395The method currently lists points on the plane projective model, that396is the closure in $\mathbb{P}^2$ of the curve defined by $y^2+hy=f$.397This means that one point $(0:1:0)$ at infinity is returned if the398degree of the curve is at least 4 and $\deg(f)>\deg(h)+1$. This point399is a singular point of the plane model. Later implementations may400consider a smooth model instead since that would be a more relevant401object. Then, for a curve whose only singularity is at $(0:1:0)$, the402point at infinity would be replaced by a number of rational points of403the smooth model. We illustrate this with an example of a genus 2404hyperelliptic curve::405406sage: R.<x>=GF(11)[]407sage: H = HyperellipticCurve(x*(x+1)*(x+2)*(x+3)*(x+4)*(x+5))408sage: H.points()409[(0 : 1 : 0), (0 : 0 : 1), (1 : 7 : 1), (1 : 4 : 1), (5 : 7 : 1), (5 : 4 : 1), (6 : 0 : 1), (7 : 0 : 1), (8 : 0 : 1), (9 : 0 : 1), (10 : 0 : 1)]410411The plane model of the genus 2 hyperelliptic curve in the above example412is the curve in $\mathbb{P}^2$ defined by $y^2z^4=g(x,z)$ where413$g(x,z)=x(x+z)(x+2z)(x+3z)(x+4z)(x+5z).$ This model has one point at414infinity $(0:1:0)$ which is also the only singular point of the plane415model. In contrast, the hyperelliptic curve is smooth and imbeds via416the equation $y^2=g(x,z)$ into weighted projected space417$\mathbb{P}(1,3,1)$. The latter model has two points at infinity:418$(1:1:0)$ and $(1:-1:0)$.419"""420from sage.rings.finite_rings.constructor import zech_log_bound421try:422return self.__points423except AttributeError: pass424425if self.base_ring().is_prime_field():426self.__points = self._points_cache_sqrt()427else:428if self.base_ring().order() < zech_log_bound:429self.__points = self._points_fast_sqrt() # this is fast using Zech logarithms430else:431self.__points = self._points_cache_sqrt()432433return self.__points434435436#This where Cartier Matrix is actually computed. This is either called by E.Cartier_matrix, E.a_number, or E.Hasse_Witt437@cached_method438def _Cartier_matrix_cached(self):439r"""440INPUT:441442- 'E' - Hyperelliptic Curve of the form `y^2 = f(x)` over a443finite field, `\mathbb{F}_q`444445OUTPUT:446447- matrix(Fq,M)' The matrix $M = (c_(pi-j)), f(x)^((p-1)/2) = \sum c_i x^i$448- 'Coeff' List of Coeffs of F, this is needed for Hasse-Witt function.449- 'g' genus of the curve self, this is needed by a-number.450- 'Fq' is the base field of self, and it is needed for Hasse-Witt451- 'p' is the char(Fq), this is needed for Hasse-Witt.452- 'E' The initial elliptic curve to check some caching conditions.453454EXAMPLES::455456sage: K.<x>=GF(9,'x')[]457sage: C=HyperellipticCurve(x^7-1,0)458sage: C._Cartier_matrix_cached()459(460[0 0 2]461[0 0 0]462[0 1 0], [2, 0, 0, 0, 0, 0, 0, 1, 0], 3, Finite Field in x of size 3^2, 3, Hyperelliptic Curve over Finite Field in x of size 3^2 defined by y^2 = x^7 + 2463)464sage: K.<x>=GF(49,'x')[]465sage: C=HyperellipticCurve(x^5+1,0)466sage: C._Cartier_matrix_cached()467(468[0 3]469[0 0], [1, 0, 0, 0, 0, 3, 0, 0, 0, 0, 3, 0, 0, 0, 0, 1], 2, Finite Field in x of size 7^2, 7, Hyperelliptic Curve over Finite Field in x of size 7^2 defined by y^2 = x^5 + 1470)471472473sage: P.<x>=GF(9,'a')[]474sage: C=HyperellipticCurve(x^29+1,0)475sage: C._Cartier_matrix_cached()476(477[0 0 1 0 0 0 0 0 0 0 0 0 0 0]478[0 0 0 0 0 1 0 0 0 0 0 0 0 0]479[0 0 0 0 0 0 0 0 1 0 0 0 0 0]480[0 0 0 0 0 0 0 0 0 0 0 1 0 0]481[0 0 0 0 0 0 0 0 0 0 0 0 0 0]482[0 0 0 0 0 0 0 0 0 0 0 0 0 0]483[0 0 0 0 0 0 0 0 0 0 0 0 0 0]484[0 0 0 0 0 0 0 0 0 0 0 0 0 0]485[0 0 0 0 0 0 0 0 0 0 0 0 0 0]486[1 0 0 0 0 0 0 0 0 0 0 0 0 0]487[0 0 0 1 0 0 0 0 0 0 0 0 0 0]488[0 0 0 0 0 0 1 0 0 0 0 0 0 0]489[0 0 0 0 0 0 0 0 0 1 0 0 0 0]490[0 0 0 0 0 0 0 0 0 0 0 0 1 0], [1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], 14, Finite Field in a of size 3^2, 3, Hyperelliptic Curve over Finite Field in a of size 3^2 defined by y^2 = x^29 + 1491)492493TESTS::494495sage: K.<x>=GF(2,'x')[]496sage: C=HyperellipticCurve(x^7-1,x)497sage: C._Cartier_matrix_cached()498Traceback (most recent call last):499...500ValueError: p must be odd501502503sage: K.<x>=GF(5,'x')[]504sage: C=HyperellipticCurve(x^7-1,4)505sage: C._Cartier_matrix_cached()506Traceback (most recent call last):507...508ValueError: E must be of the form y^2 = f(x)509510511sage: K.<x>=GF(5,'x')[]512sage: C=HyperellipticCurve(x^8-1,0)513sage: C._Cartier_matrix_cached()514Traceback (most recent call last):515...516ValueError: In this implementation the degree of f must be odd517518519sage: K.<x>=GF(5,'x')[]520sage: C=HyperellipticCurve(x^5+1,0,check_squarefree=False)521sage: C._Cartier_matrix_cached()522Traceback (most recent call last):523...524ValueError: curve is not smooth525526"""527528#Compute the finite field and prime p.529Fq=self.base_ring();530p=Fq.characteristic()531#checks532533if p == 2:534raise ValueError, "p must be odd";535536537g = self.genus()538539#retrieve the function f(x) ,where y^2=f(x)540f,h = self.hyperelliptic_polynomials()541#This implementation only deals with h=0542if h!=0:543raise ValueError, "E must be of the form y^2 = f(x)"544545d = f.degree()546#this implementation is for odd degree only, even degree will be handled later.547if d%2 == 0:548raise ValueError, "In this implementation the degree of f must be odd"549#Compute resultant to make sure no repeated roots550df=f.derivative()551R=df.resultant(f)552if R == 0:553raise ValueError, "curve is not smooth"554555#computing F, since the entries of the matrix are c_i where F= \sum c_i x^i556557F = f**((p-1)/2)558559#coefficients returns a_0, ... , a_n where f(x) = a_n x^n + ... + a_0560561Coeff = F.list()562563#inserting zeros when necessary-- that is, when deg(F) < p*g-1, (simplified if p <2g-1)564#which is the highest powered coefficient needed for our matrix565#So we don't have enough coefficients we add extra zeros to have the same poly,566#but enough coeff.567568zeros = [0 for i in range(p*g-len(Coeff))]569Coeff = Coeff + zeros570571# compute each row of matrix as list and then M=list of lists(rows)572573M=[];574for j in range(1,g+1):575H=[Coeff[i] for i in range((p*j-1), (p*j-g-1),-1)]576M.append(H);577return matrix(Fq,M), Coeff, g, Fq,p, self578579580#This is what is called from command line581def Cartier_matrix(self):582r"""583INPUT:584585- ``E`` : Hyperelliptic Curve of the form `y^2 = f(x)` over a finite field, `\mathbb{F}_q`586587OUTPUT:588589590- ``M``: The matrix `M = (c_{pi-j})`, where `c_i` are the coeffients of `f(x)^{(p-1)/2} = \sum c_i x^i`591592Reference-N. Yui. On the Jacobian varieties of hyperelliptic curves over fields of characteristic `p > 2`.593594EXAMPLES::595596sage: K.<x>=GF(9,'x')[]597sage: C=HyperellipticCurve(x^7-1,0)598sage: C.Cartier_matrix()599[0 0 2]600[0 0 0]601[0 1 0]602603sage: K.<x>=GF(49,'x')[]604sage: C=HyperellipticCurve(x^5+1,0)605sage: C.Cartier_matrix()606[0 3]607[0 0]608609sage: P.<x>=GF(9,'a')[]610sage: E=HyperellipticCurve(x^29+1,0)611sage: E.Cartier_matrix()612[0 0 1 0 0 0 0 0 0 0 0 0 0 0]613[0 0 0 0 0 1 0 0 0 0 0 0 0 0]614[0 0 0 0 0 0 0 0 1 0 0 0 0 0]615[0 0 0 0 0 0 0 0 0 0 0 1 0 0]616[0 0 0 0 0 0 0 0 0 0 0 0 0 0]617[0 0 0 0 0 0 0 0 0 0 0 0 0 0]618[0 0 0 0 0 0 0 0 0 0 0 0 0 0]619[0 0 0 0 0 0 0 0 0 0 0 0 0 0]620[0 0 0 0 0 0 0 0 0 0 0 0 0 0]621[1 0 0 0 0 0 0 0 0 0 0 0 0 0]622[0 0 0 1 0 0 0 0 0 0 0 0 0 0]623[0 0 0 0 0 0 1 0 0 0 0 0 0 0]624[0 0 0 0 0 0 0 0 0 1 0 0 0 0]625[0 0 0 0 0 0 0 0 0 0 0 0 1 0]626627TESTS::628629sage: K.<x>=GF(2,'x')[]630sage: C=HyperellipticCurve(x^7-1,x)631sage: C.Cartier_matrix()632Traceback (most recent call last):633...634ValueError: p must be odd635636637sage: K.<x>=GF(5,'x')[]638sage: C=HyperellipticCurve(x^7-1,4)639sage: C.Cartier_matrix()640Traceback (most recent call last):641...642ValueError: E must be of the form y^2 = f(x)643644645sage: K.<x>=GF(5,'x')[]646sage: C=HyperellipticCurve(x^8-1,0)647sage: C.Cartier_matrix()648Traceback (most recent call last):649...650ValueError: In this implementation the degree of f must be odd651652653sage: K.<x>=GF(5,'x')[]654sage: C=HyperellipticCurve(x^5+1,0,check_squarefree=False)655sage: C.Cartier_matrix()656Traceback (most recent call last):657...658ValueError: curve is not smooth659660"""661#checking first that Cartier matrix is not already cached. Since662#it can be called by either Hasse_Witt or a_number.663#This way it does not matter which function is called first664#in the code.665# Trac Ticket #11115: Why shall we waste time by studying666# the cache manually? We only need to check whether the cached667# data belong to self.668M, Coeffs,g, Fq, p, E= self._Cartier_matrix_cached()669if E!=self:670self._Cartier_matrix_cached.clear_cache()671M, Coeffs,g, Fq, p, E= self._Cartier_matrix_cached()672return M673674#This where Hasse_Witt is actually computed. This is either called by E.Hasse_Witt or p_rank675@cached_method676def _Hasse_Witt_cached(self):677r"""678INPUT:679680- ``E`` : Hyperelliptic Curve of the form `y^2 = f(x)` over a finite field, `\mathbb{F}_q`681682OUTPUT:683684- ``N`` : The matrix `N = M M^p \dots M^{p^{g-1}}` where `M = c_{pi-j}, f(x)^{(p-1)/2} = \sum c_i x^i`685- ``E`` : The initial curve to check some caching conditions.686687EXAMPLES::688689sage: K.<x>=GF(9,'x')[]690sage: C=HyperellipticCurve(x^7-1,0)691sage: C._Hasse_Witt_cached()692(693[0 0 0]694[0 0 0]695[0 0 0], Hyperelliptic Curve over Finite Field in x of size 3^2 defined by y^2 = x^7 + 2696)697698sage: K.<x>=GF(49,'x')[]699sage: C=HyperellipticCurve(x^5+1,0)700sage: C._Hasse_Witt_cached()701(702[0 0]703[0 0], Hyperelliptic Curve over Finite Field in x of size 7^2 defined by y^2 = x^5 + 1704)705706sage: P.<x>=GF(9,'a')[]707sage: C=HyperellipticCurve(x^29+1,0)708sage: C._Hasse_Witt_cached()709(710[0 0 0 0 0 0 0 0 0 0 0 0 0 0]711[0 0 0 0 0 0 0 0 0 0 0 0 0 0]712[0 0 0 0 0 0 0 0 0 0 0 0 0 0]713[0 0 0 0 0 0 0 0 0 0 0 0 0 0]714[0 0 0 0 0 0 0 0 0 0 0 0 0 0]715[0 0 0 0 0 0 0 0 0 0 0 0 0 0]716[0 0 0 0 0 0 0 0 0 0 0 0 0 0]717[0 0 0 0 0 0 0 0 0 0 0 0 0 0]718[0 0 0 0 0 0 0 0 0 0 0 0 0 0]719[0 0 0 0 0 0 0 0 0 0 0 0 0 0]720[0 0 0 0 0 0 0 0 0 0 0 0 0 0]721[0 0 0 0 0 0 0 0 0 0 0 0 0 0]722[0 0 0 0 0 0 0 0 0 0 0 0 0 0]723[0 0 0 0 0 0 0 0 0 0 0 0 0 0], Hyperelliptic Curve over Finite Field in a of size 3^2 defined by y^2 = x^29 + 1724)725726"""727728729730# If Cartier Matrix is already cached for this curve, use that or evaluate it to get M,731#Coeffs, genus, Fq=base field of self, p=char(Fq). This is so we have one less matrix to732#compute.733734#We use caching here since Cartier matrix is needed to compute Hasse Witt. So if the Cartier735#is already computed it is stored in list A. If it was not cached (i.e. A is empty), we simply736#compute it. If it is cached then we need to make sure that we have the correct one. So check737#which curve does the matrix correspond to. Since caching stores a lot of stuff, we only check738#the last entry in A. If it does not match, clear A and compute Cartier.739#740#Since Trac Ticket #11115, there is a different cache for methods741#that don't accept arguments. Anyway, the easiest is to call742#the cached method and simply see whether the data belong to self.743M, Coeffs,g, Fq, p, E= self._Cartier_matrix_cached()744if E!=self:745self._Cartier_matrix_cached.clear_cache()746M, Coeffs,g, Fq, p, E= self._Cartier_matrix_cached()747748#This compute the action of p^kth Frobenius on list of coefficients749def frob_mat(Coeffs, k):750a = p**k751mat = []752Coeffs_pow = [c**a for c in Coeffs]753for i in range(1,g+1):754H=[(Coeffs[j]) for j in range((p*i-1), (p*i - g-1), -1)]755mat.append(H);756return matrix(Fq,mat)757758#Computes all the different possible action of frobenius on matrix M and stores in list Mall759Mall = [M] + [frob_mat(Coeffs,k) for k in range(1,g)]760761#initial N=I, so we can go through Mall and multiply all matrices with I and762#get the Hasse-Witt matrix.763N = identity_matrix(Fq,g)764for l in Mall:765N = N*l;766return N, E767768#This is the function which is actually called by command line769def Hasse_Witt(self):770r"""771INPUT:772773- ``E`` : Hyperelliptic Curve of the form `y^2 = f(x)` over a finite field, `\mathbb{F}_q`774775OUTPUT:776777- ``N`` : The matrix `N = M M^p \dots M^{p^{g-1}}` where `M = c_{pi-j}`, and `f(x)^{(p-1)/2} = \sum c_i x^i`778779780781Reference-N. Yui. On the Jacobian varieties of hyperelliptic curves over fields of characteristic `p > 2`.782783EXAMPLES::784785sage: K.<x>=GF(9,'x')[]786sage: C=HyperellipticCurve(x^7-1,0)787sage: C.Hasse_Witt()788[0 0 0]789[0 0 0]790[0 0 0]791792sage: K.<x>=GF(49,'x')[]793sage: C=HyperellipticCurve(x^5+1,0)794sage: C.Hasse_Witt()795[0 0]796[0 0]797798sage: P.<x>=GF(9,'a')[]799sage: E=HyperellipticCurve(x^29+1,0)800sage: E.Hasse_Witt()801[0 0 0 0 0 0 0 0 0 0 0 0 0 0]802[0 0 0 0 0 0 0 0 0 0 0 0 0 0]803[0 0 0 0 0 0 0 0 0 0 0 0 0 0]804[0 0 0 0 0 0 0 0 0 0 0 0 0 0]805[0 0 0 0 0 0 0 0 0 0 0 0 0 0]806[0 0 0 0 0 0 0 0 0 0 0 0 0 0]807[0 0 0 0 0 0 0 0 0 0 0 0 0 0]808[0 0 0 0 0 0 0 0 0 0 0 0 0 0]809[0 0 0 0 0 0 0 0 0 0 0 0 0 0]810[0 0 0 0 0 0 0 0 0 0 0 0 0 0]811[0 0 0 0 0 0 0 0 0 0 0 0 0 0]812[0 0 0 0 0 0 0 0 0 0 0 0 0 0]813[0 0 0 0 0 0 0 0 0 0 0 0 0 0]814[0 0 0 0 0 0 0 0 0 0 0 0 0 0]815816817"""818# Since Trac Ticket #11115, there is a special819# type of cached for those methods that don't820# accept arguments. We want to get Hasse-Witt821# from the cache - but apparently it could be822# that the cached value does not belong to self.823# So, the easiest is:824N, E= self._Hasse_Witt_cached()825if E!=self:826self._Hasse_Witt_cached.clear_cache()827N, E= self._Hasse_Witt_cached()828return N829830def a_number(self):831r"""832INPUT:833834- ``E``: Hyperelliptic Curve of the form `y^2 = f(x)` over a finite field, `\mathbb{F}_q`835836OUTPUT:837838- ``a`` : a-number839840841EXAMPLES::842843sage: K.<x>=GF(49,'x')[]844sage: C=HyperellipticCurve(x^5+1,0)845sage: C.a_number()8461847848sage: K.<x>=GF(9,'x')[]849sage: C=HyperellipticCurve(x^7-1,0)850sage: C.a_number()8511852853sage: P.<x>=GF(9,'a')[]854sage: E=HyperellipticCurve(x^29+1,0)855sage: E.a_number()8565857858859860"""861#We use caching here since Cartier matrix is needed to compute a_number. So if the Cartier862#is already computed it is stored in list A. If it was not cached (i.e. A is empty), we simply863#compute it. If it is cached then we need to make sure that we have the correct one. So check864#which curve does the matrix correspond to. Since caching stores a lot of stuff, we only check865#the last entry in A. If it does not match, clear A and compute Cartier.866# Since Trac Ticket #11115, there is a special cache for methods867# that don't accept arguments. The easiest is: Call the cached868# method, and test whether the last entry is self.869M,Coeffs,g, Fq, p,E= self._Cartier_matrix_cached()870if E != self:871self._Cartier_matrix_cached.clear_cache()872M,Coeffs,g, Fq, p,E= self._Cartier_matrix_cached()873a=g-rank(M);874return a;875876def p_rank(self):877r"""878INPUT:879880- ``E`` : Hyperelliptic Curve of the form `y^2 = f(x)` over a finite field, `\mathbb{F}_q`881882OUTPUT:883884- ``pr`` :p-rank885886887EXAMPLES::888889sage: K.<x>=GF(49,'x')[]890sage: C=HyperellipticCurve(x^5+1,0)891sage: C.p_rank()8920893894sage: K.<x>=GF(9,'x')[]895sage: C=HyperellipticCurve(x^7-1,0)896sage: C.p_rank()8970898899sage: P.<x>=GF(9,'a')[]900sage: E=HyperellipticCurve(x^29+1,0)901sage: E.p_rank()9020903"""904#We use caching here since Hasse Witt is needed to compute p_rank. So if the Hasse Witt905#is already computed it is stored in list A. If it was not cached (i.e. A is empty), we simply906#compute it. If it is cached then we need to make sure that we have the correct one. So check907#which curve does the matrix correspond to. Since caching stores a lot of stuff, we only check908#the last entry in A. If it does not match, clear A and compute Hasse Witt.909# However, it seems a waste of time to manually analyse the cache910# -- See Trac Ticket #11115911N,E= self._Hasse_Witt_cached()912if E!=self:913self._Hasse_Witt_cached.clear_cache()914N,E= self._Hasse_Witt_cached()915pr=rank(N);916return pr917918919920