1234.sagews

def PnList(n,ell):
    """
    list of elements of P^n with entries in Z/ell*Z
    """
    V = GF(ell)^(n+1)
    V1 = V.subspaces(1)
    S = list(V1)
    return [u.basis()[0] for u in S]

def mw_gens(E):
   return E.gens()

def conj4(E,D,q):
   P, Q, M, N = mw_gens(E)
   if D == -7:
       ED = EllipticCurve([ 0, 0, 0, -5028723, -4308309810])
       R = ED(3355,128780)
   elif D == -55:
       ED = EllipticCurve([ 0, 0, 0, -310446675, -2089781471250])
       R = ED(25906, -2693296)
   elif D == -71:
       ED = EllipticCurve([0, 0, 0, -517342707, -4495601960370])
       R = ED(72891, 18576108) 
   else:
       R = E.quadratic_twist(D).gens()[0]
   p1 = PnList(3,q)
   divpols = ([(ZZ(uv[0])*P+ZZ(uv[1])*Q +ZZ(uv[2])*M + ZZ(uv[3])*N).division_points(q,poly_only=True)
for uv in p1] +
             [R.division_points(q,poly_only=True)])
   def is_in_ker(h, c):
       for p in c:
            if max([ZZ(k.poldegree()) for k in pari(h).factormod(p).list()[0]]) >2:
               return False
       return True            
       #return sum([[k.degree()]*e for k, e in h.factor_mod(p)],[])

   @cached_function
   def tau_conj(c):
       v = [g for g in divpols if is_in_ker(g, c)]
       if len(v) > 1:
           return 0
       else:
           return p1[divpols.index(v[0])]
           
   return tau_conj

def find_l(E, D, p, bound):
    """
    """
    K.<a> = QuadraticField(D)
    for q in prime_range(3,bound):
        if ( (gcd(E.ap(q),q+1))%p == 0 and len(K.factor(q)) == 1 and D%q):
            return q 

def ell_list(D, p, l, bound):
    """
    computes list of ell_i given D, p, l, and a bound on the ell_i
    """        
    K.<a> = QuadraticField(D)
    L = []
    for q in prime_range(l+1,bound):
        if len(K.factor(q))==1:
            if (gcd(E.ap(q),q+1))%(p) == 0:
                L.append(q)
    return L
                
def ell_tuples(L, r):
    """
    given list L and rank of curve r, returns list of (r-1)-tuples (distinct elements) of L without repetition
    """ 
    return Subsets(L,r-1) 

E = elliptic_curves.rank(4)[0]
D = -7
q = 3
tau = conj4(E,D,q)
r = E.rank()
l = find_l(E,D,q,20)
print "l is %s"%l
L = ell_tuples(ell_list(D, q, l, 200),r)
for i in L:
    try:
        print i, tau(i)
    except PariError:
        print "Error"

l is 17 {41, 83, 101} (1, 2, 2, 1) {41, 83, 59} (1, 2, 2, 1) {41, 83, 89} 0 {41, 59, 101} (1, 2, 2, 1) {41, 101, 89} (1, 1, 2, 2) {41, 59, 89} (1, 1, 0, 0) {83, 59, 101} 0 {89, 83, 101} (1, 2, 2, 0) {89, 83, 59} (1, 2, 2, 0) {89, 59, 101} (1, 2, 2, 0)