from sage.all import *

def multi_iter(files):
    F0 = open(files[0]).xreadlines()
    F = [open(file).xreadlines() for file in files[1:]]
    for L in F0:
        try:
            yield [L] + [X.next() for X in F]
        except:
            return
    

def iterate(file='surj.1-1000'):
    for line in open(file).xreadlines():
        N,iso,num,r,ainvs,v=line.split()
        N=int(N)
        num=int(num)
        r=int(r)
        ainvs=eval(ainvs)
        v=eval(v)
        yield N,iso,r,ainvs,v

def nontor_nonsurj(file='surj.1-1000', pmin=5):
    X = []
    for E in iterate(file):
        print E
        for p, why in E[4]:
            if p >= pmin and (not isinstance(why,str) or why.find("tor") == -1):
                X.append(E)
                break
    return X

def bad_kato(omit_isog=True):
    X = []
    for surj, isog in multi_iter(['surj.1-1000', 'allisog.1-1000']):
        N,iso,_,r,ainvs,v = surj.split()
        N = int(N); r=int(r); ainvs=eval(ainvs); v=eval(v)
        N1, iso1, _, ainvs1, _, matrix = isog.split()
        N1 = int(N1); ainvs1=eval(ainvs1); matrix=eval(matrix)
        assert N == N1, ainvs == ainvs1
        B = []
        for p in matrix[0]:
            p = int(p)
            if p != 0 and (not p in B):
                B.append(p)
        B.sort()
        if r > 0:
            continue
        for p, why in v:
            if p == 0: continue
            if N % p != 0 and p >= 5 and (not omit_isog or not p in B):
                X.append((N,iso,p,why,B))
    return X

def geom_compgrp_prod(E):
    B = 1
    for p, e in factor(E.j_invariant()):
        if e < 0:
            B *= e
    return B

def nonsurj_irred(file='surj.1-1000', r=0, max_level=1000):
    """
    List the non-surjective irreducible reps associated to rank r curves.
    """
    D = {}
    for L in open("isogdegs.1-40000").xreadlines():
        N,iso,num,P=L.split()
        N = int(N)
        if N > max_level:
            break
        D["%s%s"%(N,iso)] = int(P)

    X = []
    for E in iterate(file):
        if E[2] != r:
            continue
        if E[0] > max_level:
            break
        label ="%s%s"%(E[0],E[1])
        e = EllipticCurve(label)
        if e.has_cm():
            continue
        B = D[label] # * prod([f.torsion_order() for f in e.isogeny_class()]) #* e.tamagawa_product() #* geom_compgrp_prod(e)
        for p, why in E[4]:
            if p == 2: continue
            if p == 0:
                X.append((p,E))
            elif B % p != 0:
                X.append((p,E))
                print (p,E)
    return X

def isogeny_degs(file="allisog.1-40000"):
    O = open("isog_degrees-%s"%file,"w")
    I = open(file,"r")
    for L in I.xreadlines():
        N, iso, num, ainvs, curves, isog = L.split()
        B = 1
        G = sum(eval(isog),[])
        for n in G:
            if n:
                for p, _ in factor(n):
                    if B % p:
                        B *= p
        s = "%s %s %s \t%s"%(N,iso,num,B)                        
        O.write(s+"\n")
        O.flush()
        print s