# Listings included in the paper 
# "A bootstrap for the number of F_{q^r}-rational points on a curve over F_q"
# by S. Molina, N. Sayols and S. Xambó. arXiv ...


# Function XN(q,X,k)
# An improved procedure to get the number X_k of points of a curve C/F_q
# in arbitrary extensions F_{q^k} taking X=[X_1,...,X_g] as input, g the genus of C.
# It outputs [X_1,...,X_g,...,X_k]

# We assume that Q_ is an implementation of the rational numbers and that
# x >> Q_ is the inclusion of the integer x in Q_
def XN(q,X,k):
    g = len(X)
    if k<=g: return X[:k]
    X = [0]+X           # trick so that X[j] refers to F_{q^j}
    X = [x>>Q_ for x in X]
    S = [q**(j)+1-X[j] for j in range(1,g+1)] # Newton sums
    S = [0]+S             # similar trick
    # Computation of c1,...,cg; set c0=1 
    c = [1>>Q_]           
    for j in range(1,g+1):
        cj = S[j]
        for i in range(1,j):
            cj += c[i]*S[j-i]
        c += [-cj/j]
    # Add c_{g+i}, for i=1,...,g
    for i in range(1,g+1):
        c += [q**i*c[g-i]]
    # Find Sj for j = g+1,...,k
    for j in range(g+1,k+1):
        if j>2*g:
            Sj=0
        else: 
            Sj = j*c[j]
        for i in range(1,j):
            if i>2*g: break
            Sj += c[i]*S[j-i]
        S += [-Sj]
    # Find X[i] for i = g+1,...,k
    for i in range(g+1,k+1):
        X += [q**i+1-S[i]]
    return create_vector(X[1:])
    
# This function implements the Deuring's algorithm to list the possible cardinals #E(F_q)
# of the elliptic curves E/F_q. Our main reference here has been Waterhouse-1969. 
# We have split the computation in two parts: the function Deuring_offsets(q) 
# (which computes the list of integers t in the segment [-m,m], m=\floor{2\sqrt{q}},
# such that #E(F_q)=q+1-t for some E), and Deuring_set(q), that outputs the list in question.    

# def Deuring_offsets(q):
    P = prime_factors(q)  # prime_factors(12) => [2, 2, 3]
    p = P[0]; n = len(P)
    m = int(2*sqrt(q))
    D = [t for t in range(-m,m+1) if gcd(p,t)==1]
    if n%2==0: 
        r = p**(n//2)
        D += [-2*r,2*r] 
        if p%3 != 1: 
            D += [-r,r]
    if n%2 and (p==2 or p==3):
        r = p**((n+1)//2)
        D += [-r,r]
    if n%2 or (n%2==0 and p%4!=1):
        D += [0]
    return sorted([t for t in D])
#
def Deuring_set(q):
    D =deuring_offsets(q)
    return [t+q+1 for t in D]