for q in range(40):
  q = Integer(q)
  if q < 2: continue
  if not q.is_prime_power(): continue
  if q.is_power_of(2): continue

  Fq.<alpha> = GF(q)
  def chi(a): return a^((q-1)//2)

  for A in Fq:
    for B in Fq:
      if A*B*(A^2-4*B) == 0: continue
      for u in Fq:
        if u.is_square(): continue
        R = [r for r in Fq if 1+u*r^2!=0 and A^2*u*r^2!=B*(1+u*r^2)^2]

        def psi(r):
          if r == 0: return (0,0)
          v = -A/(1+u*r^2)
          eps = chi(v^3+A*v^2+B*v)
          x = eps*v-(1-eps)*A/2
          y = -eps*sqrt(x^3+A*x^2+B*x)
          return (x,y)

        psigraph = [(r,psi(r)) for r in R]

        for (r,p) in psigraph:
          preimages = [r2 for (r2,p2) in psigraph if p2 == p]
          if r == 0:
            print preimages == [0]
          else:
            print preimages == [r,-r] or preimages == [-r,r]

        images = set()
        for x in Fq:
          for y in Fq:
            if y^2 == x^3+A*x^2+B*x:
              if x == -A: continue
              if y == 0 and x != 0: continue
              if not (-u*x*(x+A)).is_square(): continue
              images.add((x,y))

        psiR = set()
        for r in R:
          psiR.add(psi(r))

        print psiR == images

        for x,y in psiR:
          if y == sqrt(y*y):
            print (-x/((x+A)*u)).is_square()
            r = sqrt(-x/((x+A)*u))
            print psi(r) == (x,y)
          else:
            print (-(x+A)/(u*x)).is_square()
            r = sqrt(-(x+A)/(u*x))
            print psi(r) == (x,y)