Mod_Method_to_Share

x=var('x') K.<a>=NumberField(x^2-x-1) 

%cython include 'sage/libs/pari/decl.pxi' from sage.libs.pari.gen import pari from libc.stdint cimport uint8_t, uint_fast8_t, uint32_t, uint_fast32_t, uint_fast64_t cdef extern from "pari/pari.h": cdef void NEXT_PRIME_VIADIFF(uint32_t, uint_fast8_t *) cdef uint_fast32_t[32] shiftTab shiftTab[ 0] = 0x00000001u shiftTab[ 1] = 0x00000002u shiftTab[ 2] = 0x00000004u shiftTab[ 3] = 0x00000008u shiftTab[ 4] = 0x00000010u shiftTab[ 5] = 0x00000020u shiftTab[ 6] = 0x00000040u shiftTab[ 7] = 0x00000080u shiftTab[ 8] = 0x00000100u shiftTab[ 9] = 0x00000200u shiftTab[10] = 0x00000400u shiftTab[11] = 0x00000800u shiftTab[12] = 0x00001000u shiftTab[13] = 0x00002000u shiftTab[14] = 0x00004000u shiftTab[15] = 0x00008000u shiftTab[16] = 0x00010000u shiftTab[17] = 0x00020000u shiftTab[18] = 0x00040000u shiftTab[19] = 0x00080000u shiftTab[20] = 0x00100000u shiftTab[21] = 0x00200000u shiftTab[22] = 0x00400000u shiftTab[23] = 0x00800000u shiftTab[24] = 0x01000000u shiftTab[25] = 0x02000000u shiftTab[26] = 0x04000000u shiftTab[27] = 0x08000000u shiftTab[28] = 0x10000000u shiftTab[29] = 0x20000000u shiftTab[30] = 0x40000000u shiftTab[31] = 0x80000000u cdef uint_fast8_t[256] twoDiv cdef uint_fast32_t tempItr, tempVar twoDiv[0] = 8u for tempItr in range(1,255u): tempVar = tempItr while not tempVar&1u: twoDiv[tempItr] += 1u tempVar >>= 1u cdef uint_fast32_t exp_mod(uint_fast64_t b, uint_fast32_t e, uint_fast32_t p): cdef uint_fast64_t q if e&1u: q = b else: q = 1ull e >>= 1u while e: b *= b if b > 4294967295ull: b %= p if e&1u: q *= b if q > 4294967295ull: q %= p e >>= 1u if q > 4294967295ull: q %= p return q cdef uint_fast32_t non_residue(uint_fast32_t p): cdef uint8_t *pariPrimePtr = <uint8_t *>diffptr cdef uint32_t pariP = 0u NEXT_PRIME_VIADIFF(pariP, pariPrimePtr) while True: NEXT_PRIME_VIADIFF(pariP, pariPrimePtr) if exp_mod(p,(pariP-1u)>>1u,pariP)%pariP > 1u: return pariP cdef uint_fast32_t sqrt5_mod(uint_fast32_t p): if p&3u == 3u: return exp_mod(5ull, (p+1u)>>2u, p)%p cdef uint_fast64_t q, z if p&7u == 5u: q = exp_mod(10ull, (p-5u)>>3u, p) z = q*q if z > 1844674407370955161ull: z %= p z *= 10ull z -= 1ull if z > 4294967295ull: z %= p q *= 5ull if q > 4294967295ull: q %= p return q*z%p cdef uint_fast32_t d, dp if p&15u == 9u: q = exp_mod(10ull, (p-9u)>>4u, p) if q > 4294967295ull: q %= p z = q*q if z > 1844674407370955161ull: z %= p z *= 10ull if z > 4294967295ull: z %= p q *= 5ull if q > 4294967295ull: q %= p if z*z%p == 1ull: d = non_residue(p) dp = exp_mod(d, (p-1u)>>3u, p) q *= dp if q > 4294967295ull: q %= p dp *= dp if dp > 4294967295ull: dp %= p z *= dp if z > 4294967295ull: z %= p z -= 1ull return q*z%p p -= 1u cdef uint_fast8_t r = 0u q = 8ull while q == 8ull: q = twoDiv[(p>>r)&0xFFu] r += q q = p>>r p += 1u r -= 2u cdef uint_fast64_t v = exp_mod(non_residue(p), q, p) d = exp_mod(5ull, (q-1u)>>1u, p) cdef uint_fast64_t res = 5ull*d if res > 4294967295ull: res %= p d *= d if d > 3689348814741910323ull: d %= p d *= 5ull d %= p cdef uint_fast8_t m while not d == 1u: m = r dp = d*d%p while not dp == 1u: dp *= dp dp %= p m -= 1u z = exp_mod(v, shiftTab[m], p) res *= z if res > 4294967295ull: res %= p z *= z if z > 4294967295ull: z %= p d *= z d %= p return res%p def get_primes(p): q = sqrt5_mod(p)-1u w = ((p+1u)>>1u)*q%p q = p-1u-w if w < q: return (w,q) return (q,w) 

def split_prime_combine(p, firstCurves, secondCurves): ret = [] r1,r2 = get_primes(p) if len(firstCurves) < len(secondCurves): firstCurves, secondCurves = secondCurves, firstCurves r1,r2 = r2,r1 for E in firstCurves: liftE = [E[0],0,E[1],0] for i in range(p): for j in range(p): if ((liftE[0]-r2*liftE[1])%p,(liftE[2]-r2*liftE[3])%p) in secondCurves: ret.append(list(liftE)) liftE[2] = (liftE[2]+r1)%p liftE[3] += 1 liftE[3] = 0 liftE[0] = (liftE[0]+r1)%p liftE[1] += 1 return ret def integral_combine(n, first, second, firstCurves, secondCurves): ret = [] for E in firstCurves: for i in range(second): for j in range(second): for k in range(second): for l in range(second): if [A%second for A in E] in secondCurves: ret.append(list(E)) E[3] = (E[3]+first)%n E[2] = (E[2]+first)%n E[1] = (E[1]+first)%n E[0] = (E[0]+first)%n return ret def curves_in_GF(p): ret = {} K = GF(p) for i in range(p): for j in range(p): try: E = EllipticCurve(K,[i,j]) ap = p-len(E.points())+1 if ret.has_key(ap): ret[ap].append(E) else: ret[ap] = [E] except ArithmeticError: pass return ret def various_lifts(A,B,p): q=(p+1)/2 if p%3==1: v = (2*p+1)/3 if p%3==2: v = (p+1)/3 ret = [] for s in [0,q,q*a,q+q*a]: for r in [v*s^2,v*(s^2+1),v*(s^2-1),v*(s^2+a),v*(s^2-a),v*(s^2+a+1),v*(s^2+a-1),v*(s^2-a+1),v*(s^2-a-1)]: for t in [0,q,q*a,q+q*a]: a1 = K(2*s)[0]%p+K(2*s)[1]%p*a a2 = K(3*r-s^2)[0]%p+K(3*r-s^2)[1]%p*a if a2[0]>1: comp0=a2[0]-p else: comp0 = a2[0] if a2[1]>1: comp1 = a2[1]-p else: comp1 = a2[1] a2 = comp0+comp1*a a3 = K(2*t)[0]%p+K(2*t)[1]%p*a a4 = K(A+3*r^2-2*s*t)[0]%p+K(A+3*r^2-2*s*t)[1]%p*a a6 = K(B+r*A+r^3-t^2)[0]%p+K(B+r*A+r^3-t^2)[1]%p*a for c in [a4,a4-p,a4-p*a,a4-p-p*a]: for d in [a6,a6-p,a6-p*a,a6-p-p*a]: #print (a1,a2,a3,c,d) ret.append((a1,a2,a3,c,d)) return ret def quick_disc(ainv): a1 = ainv[0] a2 = ainv[1] a3 = ainv[2] a4 = ainv[3] a6 = ainv[4] return -(a1^2+4*a2)^2*((a1^2+4*a2)*a6-a1*a3*a4+a2*a3^2-a4^2)-8*(a1*a3+2*a4)^3-27*(a3^2+4*a6)^2+9*(a1^2+4*a2)*(a1*a3+2*a4)*(a3^2+4*a6) def ap(E,p): return E.change_ring(p.residue_field()).trace_of_frobenius() def check_aps_against(E,aps): ret = True for i,p in enumerate([11,19,29,31,41,59,61,71,79,89]): r1 = K.primes_above(p)[0] r2 = K.primes_above(p)[1] if type(aps[2*i]) != str: if ap(E,r1) != aps[2*i] or ap(E,r2) != aps[2*i+1]: ret = False return ret 

aps = ['?','?',-4,8,6,-6,-4,8,6,-6,0,-12,2,-10,0,0,8,-4,6,-18] di={K(3*a-2): aps[0], K(3*a-1): aps[1], K(-4*a+1): aps[2], K(-4*a+3): aps[3], K(-a+6): aps[4], K(a+5): aps[5], K(5*a-2): aps[6], K(5*a-3): aps[7], K(a-7): aps[8], K(a+6): aps[9], K(7*a-2): aps[10], K(7*a-5): aps[11], K(7*a-3): aps[12], K(7*a-4): aps[13], K(a-9): aps[14], K(a+8): aps[15], K(-8*a+5): aps[16], K(-8*a+3): aps[17], K(a-10): aps[18], K(a+9): aps[19]} N=275 cond = 15*a-10 found = False for p in (11,19,29,31,41,59,61,71,79,89): print 'CURRENT PRIME: ',p r1 = K.primes_above(p)[0] r2 = K.primes_above(p)[1] ap1 = di[r1.gens_reduced()[0]] ap2 = di[r2.gens_reduced()[0]] if type(ap1) == str: continue c = r1.integral_basis()[1][0] d = r2.integral_basis()[1][0] if c>d: r1,r2 = r2,r1 ap1,ap2 = ap2,ap1 modp = curves_in_GF(p) firstp = [] for E in modp[ap1]: firstp.append((int(E.a4()),int(E.a6()))) secondp = [] for E in modp[ap2]: secondp.append((int(E.a4()),int(E.a6()))) big_list_newp = split_prime_combine(p,firstp,secondp) print 'Length of big_list_new%s:'%p,len(big_list_newp) current = 1 for e in big_list_newp: print current current += 1 A = e[0]+e[1]*a B = e[2]+e[3]*a for ainv in various_lifts(A,B,p): D = quick_disc(ainv).norm() if D%N==0: E = EllipticCurve(K,list(ainv)) if E.conductor() == cond: print E, 'tentative' if check_aps_against(E,aps): print E, 'a_ps match' found = True break if found: break if found: break 

CURRENT PRIME:  11
CURRENT PRIME:  19
Length of big_list_new19: 432
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113

CURRENT PRIME:  11
CURRENT PRIME:  19
Length of big_list_new19: 432
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113