cleanup ordp.sage code

script/research/ec/ordp.sage

@@ -1,68 +1,74 @@
 from tabulate import tabulate
 # This is a more usable version of valuate.sage, less instructional
 
-K.<x, y> = GF(11)[]
-Px, Py = K(2), K(4)
-S = K.quotient(y^2 - x^3 - 4*x).fraction_field()
-X, Y = S(x), S(y)
-
-EC_A = 4
-EC_B = 0
-EC = y^2 - x^3 - EC_A*x - EC_B
-
-original_f = (y - 2*x)^2
-b0, b1, b2 = [(x - Px), (y - Py), 1]
-
 # Return components for basis
-def decomp(f):
+def decomp(f, basis):
+    b0, b1, b2 = basis
     a0, r = f.quo_rem(b0)
     a1, r = r.quo_rem(b1)
     a2, r = r.quo_rem(b2)
     assert r == 0
     return [a0, a1, a2]
 
-def comp(comps):
-    return sum(a*b for a, b in zip(comps, (b0, b1, b2)))
+def comp(comps, basis):
+    return sum(a*b for a, b in zip(comps, basis))
 
-assert comp(decomp(original_f)) == original_f
-
-# so we can replace (y - Py) with this
-Ef = b0^2 + binomial(3,2)*Px*b0^1 + (3*Px^2 + EC_A)
-Eg = (y + Py)
-assert EC == b1*Eg - b0*Ef
-
-def apply_reduction(a, g):
+def apply_reduction(a, g, Ef, Eg):
     assert a[2] == 0
     a[0] = a[0]*Eg + a[1]*Ef
     a[1] = 0
     g[0] *= Eg
 
-k = 0
-a = [original_f, 0, 0]
-g = [1]
+# EC_A, EC_B must be defined before calling this function
+def ordp(P, original_f, debug=True):
+    EC = y^2 - x^3 - EC_A*x - EC_B
+
+    Px, Py = P
+    b0, b1, b2 = basis = [(x - Px), (y - Py), 1]
+
+    # so we can replace (y - Py) with this
+    Ef = b0^2 + binomial(3,2)*Px*b0^1 + (3*Px^2 + EC_A)
+    Eg = (y + Py)
+    assert EC == b1*Eg - b0*Ef
+
+    k = 0
+    a = [original_f, 0, 0]
+    g = [1]
+
+    table = []
+    table.append(("", "a", "g", "k"))
+
+    def log(step_name, a, g, k):
+        table.append((step_name, str(a), str(g), k))
+
+    log("start", a, g, k)
 
-table = []
-table.append(("", "a", "g", "k"))
+    while True:
+        f = a[0]
+        a = decomp(f, basis)
+        log("decomp", a, g, k)
 
-def log(step_name, a, g, k):
-    table.append((step_name, str(a), str(g), k))
+        # Check remainder
+        if a[2] != 0:
+            break
 
-log("start", a, g, k)
+        # We can apply a reduction
+        k += 1
 
-while True:
-    f = a[0]
-    a = decomp(f)
-    log("decomp", a, g, k)
+        apply_reduction(a, g, Ef, Eg)
+        log("reduce", a, g, k)
 
-    # Check remainder
-    if a[2] != 0:
-        break
+    if debug:
+        print(tabulate(table))
 
-    # We can apply a reduction
-    k += 1
+    return k
 
-    apply_reduction(a, g)
-    log("reduce", a, g, k)
+if __name__ == "__main__":
+    K.<x, y> = GF(11)[]
+    EC_A = 4
+    EC_B = 0
+    P = (2, 4)
+    f = y - 2*x
+    k = ordp(P, f)
+    print(f"k = {k}")
 
-print(tabulate(table))
-print(f"k = {k}")