another code 221 for FM16

sage/tnfs/curve/fotiadis_martindale.py

@@ -31,7 +31,7 @@ Curves from ePrint 2019/555
 Also from 2018/969
 """
 
-allowed_code = [1, 10, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29]
+allowed_code = [1, 10, 17, 18, 19, 20, 21, 22, 221, 23, 24, 25, 26, 27, 28, 29]
 
 def polynomial_params_from_code(code):
     
@@ -82,7 +82,7 @@ def polynomial_params_from_code(code):
         #yx = (tx-2)* ((K(-1/D).sqrt()).polynomial()) % rx
         #px = (tx**2 + D*yx**2)/4
         
-        rx = 36*x**4 + 36*x**3 +18*x**2 + 6*x + 1 
+        rx = 36*x**4 + 36*x**3 +18*x**2 + 6*x + 1
         tx = (-6*x**2) + 1
         yx = 48*x**3 + 30*x**2 + 12*x + 3
         px = 1728*x**6 + 2160*x**5 + 1548*x**4 + 756*x**3 + 240*x**2 + 54*x + 7
@@ -113,7 +113,7 @@ def polynomial_params_from_code(code):
         D=3
         m=30
         u_mod_m = [8,23]
-        rx = (x**4 -2*x**3 -3*x**2 + 4*x + 13)/225 # because when 
+        rx = (x**4 -2*x**3 -3*x**2 + 4*x + 13)/225
         tx = (-x**3+4*x**2+5*x+6)/15
         yx = (x**3 - 4*x**2 + 5*x + 4)/15
         px = (x**6-8*x**5+21*x**4-17*x**3+13*x**2+45*x+21)/225
@@ -128,7 +128,7 @@ def polynomial_params_from_code(code):
         D=3
         m=285
         u_mod_m = [209,266]
-        rx = (x**4 - 37*x**2 + 361)/9025 # because for u=209, 266 mod 285, p(u)  in ZZ and 9025 | r(u)
+        rx = (x**4 - 37*x**2 + 361)/9025 # because for u=209, 266 mod 285, p(u) in ZZ and 9025 | r(u)
         tx = (-2*x**3 + 17*x + 95)/95
         yx = (8*x**3 + 38*x**2 - 163*x - 703)/285
         px = (x**6 + 8*x**5 - 18*x**4 - 326*x**3 - 342*x**2 + 3143*x + 6859)/1425
@@ -144,35 +144,59 @@ def polynomial_params_from_code(code):
         D=3
         m=3
         u_mod_m = [0]
+        # by design, p = 3,4 mod 5 but p = 1 mod 5 is not possible
         rx = x**8 - x**7 + x**5 - x**4 + x**3 - x + 1
+        # rx == cyclotomic_polynomial(15)
         tx = x**8 - x**7 + x**5 + x**3 - x + 2
-        # where (tx-2)/sqrt(-3) mod rx is
+        assert tx == x**4 + 1 + rx
+        # where yx0 = (tx-2)/sqrt(-3) mod rx is
         yx0 = (2*x**7 - 2*x**6 - 2*x**5 + x**4 - 2*x**3 + 2*x**2 - 3)/3
         yx = (3*x**8 - 5*x**7 + 2*x**6 + 5*x**5 - 4*x**4 + 5*x**3 - 2*x**2 - 3*x + 6)/3
         assert yx == -yx0 + rx
         px = (3*x**16 - 9*x**15 + 10*x**14 + 4*x**13 - 19*x**12 + 24*x**11 - 14*x**10 - 6*x**9 + 37*x**8 - 36*x**7 + 8*x**6 + 19*x**5 - 20*x**4 + 21*x**3 - 3*x**2 - 12*x + 12)/3
         cx = (3*x**8 - 6*x**7 + 4*x**6 + 5*x**5 - 5*x**4 + 6*x**3 - 3*x**2 + 9)/3
-        zeta_k_mod_rx = tx-rx-1
-        exp_tr = 1
+        zeta_k_mod_rx = x
+        exp_tr = 4
         betax = (1245*x**15 - 2397*x**14 + 940*x**13 + 3840*x**12 - 4777*x**11 + 4318*x**10 - 294*x**9 - 6426*x**8 + 12547*x**7 - 3262*x**6 - 4656*x**5 + 6061*x**4 - 5922*x**3 + 3855*x**2 + 4269*x - 5514)/1044
         lambx = x**5
 
+    elif code == 221:
+        # another family with rho = 2 r = Phi_15(x) and this time tr = (x^2+1)-r
+        k=15
+        D=3
+        m=15
+        # u = 1,2 mod 3 so that px is an integer and u = 1,2,3,4 mod 5 so that p = 1 mod 5
+        u_mod_m = [1, 2, 4, 7, 8, 11, 13, 14]
+        rho = 2
+        px = (3*x**16 - 9*x**15 + 7*x**14 + 9*x**13 - 20*x**12 + 15*x**11 - 5*x**10 - 11*x**9 + 27*x**8 - 20*x**7 - 3*x**6 + 13*x**5 - 11*x**4 + 9*x**3 + x**2 - 6*x + 4)/3
+        rx = x**8 - x**7 + x**5 - x**4 + x**3 - x + 1
+        zeta_k_mod_rx = x**2
+        exp_tr = 2
+        tx = -x**8 + x**7 - x**5 + x**4 - x**3 + x**2 + x
+        assert tx == - rx + x**2 + 1
+        yx = (-3*x**8 + 5*x**7 - 5*x**5 + 3*x**4 - 3*x**3 + x**2 + 3*x - 4)/3
+        cx = (3*x**8 - 6*x**7 + x**6 + 7*x**5 - 4*x**4 + x**3 - 4*x**2 - 2*x + 7)/3
+        betax = (18*x**15 - 60*x**14 + 51*x**13 + 70*x**12 - 180*x**11 + 128*x**10 + 41*x**9 - 208*x**8 + 224*x**7 - 59*x**6 - 112*x**5 + 163*x**4 - 113*x**3 + 83*x - 60)/22
+        lambx = x**5
+
     elif code == 22:
         # table 2 p. 9 of eprint 2019/555 with x -> x/3
         k=15
         D=3
-        m=3
-        u_mod_m = [0]
+        m=15
+        u_mod_m = [0, 3] # x = 0 mod 3 so that p is an integer, u = 0,3 mod 15 so that p = 1 mod 5
         rx = x**8 - x**7 + x**5 - x**4 + x**3 - x + 1
+        # rx == cyclotomic_polynomial(15)
         tx = x**8 - x**7 + x**5 + x**3 - x + 2
-        # where (tx-2)/sqrt(-3) mod rx is
+        assert tx == x**4 + 1 + rx
+        # where yx0 = (tx-2)/sqrt(-3) mod rx is
         yx0 = (2*x**7 - 2*x**6 - 2*x**5 + x**4 - 2*x**3 + 2*x**2 - 3)/3
         yx = x * (3*x**7 - x**6 - 2*x**5 + x**4 - 2*x**3 + x**2 + 2*x - 3)/3
         assert yx == yx0 + rx
         px = (3*x**16 - 3*x**15 - 2*x**14 + 4*x**13 - 4*x**12 + 3*x**11 + 4*x**10 - 9*x**9 + 7*x**8 - 4*x**6 + 7*x**5 - 2*x**4 + 3*x**2 - 3*x + 3)/3
         cx = x**2 * (3*x**6 - 2*x**4 - x**3 - 2*x**2 + 3)/3
-        zeta_k_mod_rx = tx-rx-1
-        exp_tr = 1
+        zeta_k_mod_rx = x
+        exp_tr = 4
         betax = (588*x**15 - 564*x**14 - 323*x**13 + 96*x**12 - 412*x**11 + 745*x**10 + 477*x**9 - 354*x**8 + 376*x**7 - 547*x**6 + 165*x**5 + 70*x**4 + 255*x**3 + 174*x**2 - 366*x - 378)/501
         lambx = x**5
 
@@ -227,7 +251,7 @@ def polynomial_params_from_code(code):
         yx = (3*x**6 + x**4 - x**3 - 2*x + 2)/3
         zeta_k_mod_rx = x
         exp_tr = 13
-        lambx = x**6
+        lambx = x**3-1
         betax = (-93*x**11 + 48*x**10 + 240*x**9 + 267*x**8 - 235*x**7 + 84*x**6 - 139*x**5 + 191*x**4 + 131*x**3 + 266*x**2 - 411*x - 112)/342
 
     elif code==26: