def SY(n):
    L=[n]
    S=n
    if n==1:
        return L
    while S!=1:
        if S%2==0:
            S=S/2
            L.append(S)
        else :
            S=3*S+1
            L.append(S)
    return L
        
#Here I define the Syracuse sequence
    

def PSY(n):
    A.<X>=PolynomialRing(CC)
    P=A(SY(n))
    return P

#Here I define the Polynomial PSY(n) associated with the Syracuse sequence starting with the number n

PSY(21).roots(RR)

def PSY_EvenRoots(n):
    L=[]
    for i in range (2,n):
        N=2*i
        if PSY(N).roots(RR)!=[]:
            for j in range (len(PSY(N).roots(RR))):
                if -0.5<=PSY(N).roots(RR)[j][0]<=0.0:
                    L.append(N)
    return L

#I look at PSY(n) for n an even number, and see wether the polynomial has a root in the place of interest

def PSY_OddRoots(n):
    L=[]
    for i in range (2,n):
        N=2*i+1
        if PSY(N).roots(RR)!=[]:
            for j in range (len(PSY(N).roots(RR))):
                if -0.5<=PSY(N).roots(RR)[j][0]<=0.0:
                    L.append(N)
    return L

#Same but for n an odd number
        

def psy(n):
    A.<X>=PolynomialRing(CC)
    L=SY(n)[::-1]
    return A(L)

#Here I define the "reverse polynomial"

def psy_goodroots(n):
    L=[]
    for i in range(len(psy(n).roots(CC))):
        if -0.5<psy(n).roots(CC)[i][0].real()<-0.35 and 0.29<psy(n).roots(CC)[i][0].imag()<0.35:
               L.append(psy(n).roots(CC)[i][0])    
    return L

psy_goodroots(17)

complex_plot((psy(1000000000000000000007)),(-1,1),(-1,1))