Change exhaustive test groups so they have a point with X=1

This enables testing overflow is correctly encoded in the recid, and
likely triggers more edge cases.

Also introduce a Sage script to generate the parameters.

sage/gen_exhaustive_groups.sage

@@ -0,0 +1,129 @@
+# Define field size and field
+P = 2^256 - 2^32 - 977
+F = GF(P)
+BETA = F(0x7ae96a2b657c07106e64479eac3434e99cf0497512f58995c1396c28719501ee)
+
+assert(BETA != F(1) and BETA^3 == F(1))
+
+orders_done = set()
+results = {}
+first = True
+for b in range(1, P):
+    # There are only 6 curves (up to isomorphism) of the form y^2=x^3+B. Stop once we have tried all.
+    if len(orders_done) == 6:
+        break
+
+    E = EllipticCurve(F, [0, b])
+    print("Analyzing curve y^2 = x^3 + %i" % b)
+    n = E.order()
+    # Skip curves with an order we've already tried
+    if n in orders_done:
+        print("- Isomorphic to earlier curve")
+        continue
+    orders_done.add(n)
+    # Skip curves isomorphic to the real secp256k1
+    if n.is_pseudoprime():
+        print(" - Isomorphic to secp256k1")
+        continue
+
+    print("- Finding subgroups")
+
+    # Find what prime subgroups exist
+    for f, _ in n.factor():
+        print("- Analyzing subgroup of order %i" % f)
+        # Skip subgroups of order >1000
+        if f < 4 or f > 1000:
+            print("  - Bad size")
+            continue
+
+        # Iterate over X coordinates until we find one that is on the curve, has order f,
+        # and for which curve isomorphism exists that maps it to X coordinate 1.
+        for x in range(1, P):
+            # Skip X coordinates not on the curve, and construct the full point otherwise.
+            if not E.is_x_coord(x):
+                continue
+            G = E.lift_x(F(x))
+
+            print("  - Analyzing (multiples of) point with X=%i" % x)
+
+            # Skip points whose order is not a multiple of f. Project the point to have
+            # order f otherwise.
+            if (G.order() % f):
+                print("    - Bad order")
+                continue
+            G = G * (G.order() // f)
+
+            # Find lambda for endomorphism. Skip if none can be found.
+            lam = None
+            for l in Integers(f)(1).nth_root(3, all=True):
+                if int(l)*G == E(BETA*G[0], G[1]):
+                    lam = int(l)
+                    break
+            if lam is None:
+                print("    - No endomorphism for this subgroup")
+                break
+
+            # Now look for an isomorphism of the curve that gives this point an X
+            # coordinate equal to 1.
+            # If (x,y) is on y^2 = x^3 + b, then (a^2*x, a^3*y) is on y^2 = x^3 + a^6*b.
+            # So look for m=a^2=1/x.
+            m = F(1)/G[0]
+            if not m.is_square():
+                print("    - No curve isomorphism maps it to a point with X=1")
+                continue
+            a = m.sqrt()
+            rb = a^6*b
+            RE = EllipticCurve(F, [0, rb])
+
+            # Use as generator twice the image of G under the above isormorphism.
+            # This means that generator*(1/2 mod f) will have X coordinate 1.
+            RG = RE(1, a^3*G[1]) * 2
+            # And even Y coordinate.
+            if int(RG[1]) % 2:
+                RG = -RG
+            assert(RG.order() == f)
+            assert(lam*RG == RE(BETA*RG[0], RG[1]))
+
+            # We have found curve RE:y^2=x^3+rb with generator RG of order f. Remember it
+            results[f] = {"b": rb, "G": RG, "lambda": lam}
+            print("    - Found solution")
+            break
+
+    print("")
+
+print("")
+print("")
+print("/* To be put in src/group_impl.h: */")
+first = True
+for f in sorted(results.keys()):
+    b = results[f]["b"]
+    G = results[f]["G"]
+    print("#  %s EXHAUSTIVE_TEST_ORDER == %i" % ("if" if first else "elif", f))
+    first = False
+    print("static const secp256k1_ge secp256k1_ge_const_g = SECP256K1_GE_CONST(")
+    print("    0x%08x, 0x%08x, 0x%08x, 0x%08x," % tuple((int(G[0]) >> (32 * (7 - i))) & 0xffffffff for i in range(4)))
+    print("    0x%08x, 0x%08x, 0x%08x, 0x%08x," % tuple((int(G[0]) >> (32 * (7 - i))) & 0xffffffff for i in range(4, 8)))
+    print("    0x%08x, 0x%08x, 0x%08x, 0x%08x," % tuple((int(G[1]) >> (32 * (7 - i))) & 0xffffffff for i in range(4)))
+    print("    0x%08x, 0x%08x, 0x%08x, 0x%08x" % tuple((int(G[1]) >> (32 * (7 - i))) & 0xffffffff for i in range(4, 8)))
+    print(");")
+    print("static const secp256k1_fe secp256k1_fe_const_b = SECP256K1_FE_CONST(")
+    print("    0x%08x, 0x%08x, 0x%08x, 0x%08x," % tuple((int(b) >> (32 * (7 - i))) & 0xffffffff for i in range(4)))
+    print("    0x%08x, 0x%08x, 0x%08x, 0x%08x" % tuple((int(b) >> (32 * (7 - i))) & 0xffffffff for i in range(4, 8)))
+    print(");")
+print("#  else")
+print("#    error No known generator for the specified exhaustive test group order.")
+print("#  endif")
+
+print("")
+print("")
+print("/* To be put in src/scalar_impl.h: */")
+first = True
+for f in sorted(results.keys()):
+    lam = results[f]["lambda"]
+    print("#  %s EXHAUSTIVE_TEST_ORDER == %i" % ("if" if first else "elif", f))
+    first = False
+    print("#    define EXHAUSTIVE_TEST_LAMBDA %i" % lam)
+print("#  else")
+print("#    error No known lambda for the specified exhaustive test group order.")
+print("#  endif")
+print("")

src/group_impl.h

@@ -11,52 +11,38 @@
 #include "field.h"
 #include "group.h"
 
-/* These points can be generated in sage as follows:
+/* These exhaustive group test orders and generators are chosen such that:
+ * - The field size is equal to that of secp256k1, so field code is the same.
+ * - The curve equation is of the form y^2=x^3+B for some constant B.
+ * - The subgroup has a generator 2*P, where P.x=1.
+ * - The subgroup has size less than 1000 to permit exhaustive testing.
+ * - The subgroup admits an endomorphism of the form lambda*(x,y) == (beta*x,y).
  *
- * 0. Setup a worksheet with the following parameters.
+ * These parameters are generated using sage/gen_exhaustive_groups.sage.
- *   b = 4  # whatever secp256k1_fe_const_b will be set to
- *   F = FiniteField (0xFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFEFFFFFC2F)
- *   C = EllipticCurve ([F (0), F (b)])
- *
- * 1. Determine all the small orders available to you. (If there are
- *    no satisfactory ones, go back and change b.)
- *   print C.order().factor(limit=1000)
- *
- * 2. Choose an order as one of the prime factors listed in the above step.
- *    (You can also multiply some to get a composite order, though the
- *    tests will crash trying to invert scalars during signing.) We take a
- *    random point and scale it to drop its order to the desired value.
- *    There is some probability this won't work; just try again.
- *   order = 199
- *   P = C.random_point()
- *   P = (int(P.order()) / int(order)) * P
- *   assert(P.order() == order)
- *
- * 3. Print the values. You'll need to use a vim macro or something to
- *    split the hex output into 4-byte chunks.
- *   print "%x %x" % P.xy()
  */
 #if defined(EXHAUSTIVE_TEST_ORDER)
-#  if EXHAUSTIVE_TEST_ORDER == 199
+#  if EXHAUSTIVE_TEST_ORDER == 13
 static const secp256k1_ge secp256k1_ge_const_g = SECP256K1_GE_CONST(
-    0xFA7CC9A7, 0x0737F2DB, 0xA749DD39, 0x2B4FB069,
+    0xc3459c3d, 0x35326167, 0xcd86cce8, 0x07a2417f,
-    0x3B017A7D, 0xA808C2F1, 0xFB12940C, 0x9EA66C18,
+    0x5b8bd567, 0xde8538ee, 0x0d507b0c, 0xd128f5bb,
-    0x78AC123A, 0x5ED8AEF3, 0x8732BC91, 0x1F3A2868,
+    0x8e467fec, 0xcd30000a, 0x6cc1184e, 0x25d382c2,
-    0x48DF246C, 0x808DAE72, 0xCFE52572, 0x7F0501ED
+    0xa2f4494e, 0x2fbe9abc, 0x8b64abac, 0xd005fb24
 );
-
+static const secp256k1_fe secp256k1_fe_const_b = SECP256K1_FE_CONST(
-static const secp256k1_fe secp256k1_fe_const_b = SECP256K1_FE_CONST(0, 0, 0, 0, 0, 0, 0, 4);
+    0x3d3486b2, 0x159a9ca5, 0xc75638be, 0xb23a69bc,
-
+    0x946a45ab, 0x24801247, 0xb4ed2b8e, 0x26b6a417
-#  elif EXHAUSTIVE_TEST_ORDER == 13
+);
+#  elif EXHAUSTIVE_TEST_ORDER == 199
 static const secp256k1_ge secp256k1_ge_const_g = SECP256K1_GE_CONST(
-    0xedc60018, 0xa51a786b, 0x2ea91f4d, 0x4c9416c0,
+    0x226e653f, 0xc8df7744, 0x9bacbf12, 0x7d1dcbf9,
-    0x9de54c3b, 0xa1316554, 0x6cf4345c, 0x7277ef15,
+    0x87f05b2a, 0xe7edbd28, 0x1f564575, 0xc48dcf18,
-    0x54cb1b6b, 0xdc8c1273, 0x087844ea, 0x43f4603e,
+    0xa13872c2, 0xe933bb17, 0x5d9ffd5b, 0xb5b6e10c,
-    0x0eaf9a43, 0xf6effe55, 0x939f806d, 0x37adf8ac
+    0x57fe3c00, 0xbaaaa15a, 0xe003ec3e, 0x9c269bae
+);
+static const secp256k1_fe secp256k1_fe_const_b = SECP256K1_FE_CONST(
+    0x2cca28fa, 0xfc614b80, 0x2a3db42b, 0x00ba00b1,
+    0xbea8d943, 0xdace9ab2, 0x9536daea, 0x0074defb
 );
-
-static const secp256k1_fe secp256k1_fe_const_b = SECP256K1_FE_CONST(0, 0, 0, 0, 0, 0, 0, 2);
-
 #  else
 #    error No known generator for the specified exhaustive test group order.
 #  endif

src/modules/recovery/tests_exhaustive_impl.h

@@ -25,6 +25,7 @@ void test_exhaustive_recovery_sign(const secp256k1_context *ctx, const secp256k1
                 unsigned char sk32[32], msg32[32];
                 int expected_recid;
                 int recid;
+                int overflow;
                 secp256k1_scalar_set_int(&msg, i);
                 secp256k1_scalar_set_int(&sk, j);
                 secp256k1_scalar_get_b32(sk32, &sk);
@@ -34,17 +35,18 @@ void test_exhaustive_recovery_sign(const secp256k1_context *ctx, const secp256k1
 
                 /* Check directly */
                 secp256k1_ecdsa_recoverable_signature_load(ctx, &r, &s, &recid, &rsig);
-                r_from_k(&expected_r, group, k);
+                r_from_k(&expected_r, group, k, &overflow);
                 CHECK(r == expected_r);
                 CHECK((k * s) % EXHAUSTIVE_TEST_ORDER == (i + r * j) % EXHAUSTIVE_TEST_ORDER ||
                       (k * (EXHAUSTIVE_TEST_ORDER - s)) % EXHAUSTIVE_TEST_ORDER == (i + r * j) % EXHAUSTIVE_TEST_ORDER);
                 /* The recid's second bit is for conveying overflow (R.x value >= group order).
                  * In the actual secp256k1 this is an astronomically unlikely event, but in the
-                 * small group used here, it will always be the case.
+                 * small group used here, it will be the case for all points except the ones where
+                 * R.x=1 (which the group is specifically selected to have).
                  * Note that this isn't actually useful; full recovery would need to convey
                  * floor(R.x / group_order), but only one bit is used as that is sufficient
                  * in the real group. */
-                expected_recid = 2;
+                expected_recid = overflow ? 2 : 0;
                 r_dot_y_normalized = group[k].y;
                 secp256k1_fe_normalize(&r_dot_y_normalized);
                 /* Also the recovery id is flipped depending if we hit the low-s branch */
@@ -61,7 +63,7 @@ void test_exhaustive_recovery_sign(const secp256k1_context *ctx, const secp256k1
                 /* Note that we compute expected_r *after* signing -- this is important
                  * because our nonce-computing function function might change k during
                  * signing. */
-                r_from_k(&expected_r, group, k);
+                r_from_k(&expected_r, group, k, NULL);
                 CHECK(r == expected_r);
                 CHECK((k * s) % EXHAUSTIVE_TEST_ORDER == (i + r * j) % EXHAUSTIVE_TEST_ORDER ||
                       (k * (EXHAUSTIVE_TEST_ORDER - s)) % EXHAUSTIVE_TEST_ORDER == (i + r * j) % EXHAUSTIVE_TEST_ORDER);
@@ -104,7 +106,7 @@ void test_exhaustive_recovery_verify(const secp256k1_context *ctx, const secp256
                     should_verify = 0;
                     for (k = 0; k < EXHAUSTIVE_TEST_ORDER; k++) {
                         secp256k1_scalar check_x_s;
-                        r_from_k(&check_x_s, group, k);
+                        r_from_k(&check_x_s, group, k, NULL);
                         if (r_s == check_x_s) {
                             secp256k1_scalar_set_int(&s_times_k_s, k);
                             secp256k1_scalar_mul(&s_times_k_s, &s_times_k_s, &s_s);

src/scalar_impl.h

@@ -253,6 +253,7 @@ static void secp256k1_scalar_inverse_var(secp256k1_scalar *r, const secp256k1_sc
 }
 
 #ifdef USE_ENDOMORPHISM
+/* These parameters are generated using sage/gen_exhaustive_groups.sage. */
 #if defined(EXHAUSTIVE_TEST_ORDER)
 #  if EXHAUSTIVE_TEST_ORDER == 13
 #    define EXHAUSTIVE_TEST_LAMBDA 9

src/tests_exhaustive.c

@@ -215,14 +215,14 @@ void test_exhaustive_ecmult_multi(const secp256k1_context *ctx, const secp256k1_
     secp256k1_scratch_destroy(&ctx->error_callback, scratch);
 }
 
-void r_from_k(secp256k1_scalar *r, const secp256k1_ge *group, int k) {
+void r_from_k(secp256k1_scalar *r, const secp256k1_ge *group, int k, int* overflow) {
     secp256k1_fe x;
     unsigned char x_bin[32];
     k %= EXHAUSTIVE_TEST_ORDER;
     x = group[k].x;
     secp256k1_fe_normalize(&x);
     secp256k1_fe_get_b32(x_bin, &x);
-    secp256k1_scalar_set_b32(r, x_bin, NULL);
+    secp256k1_scalar_set_b32(r, x_bin, overflow);
 }
 
 void test_exhaustive_verify(const secp256k1_context *ctx, const secp256k1_ge *group) {
@@ -250,7 +250,7 @@ void test_exhaustive_verify(const secp256k1_context *ctx, const secp256k1_ge *gr
                     should_verify = 0;
                     for (k = 0; k < EXHAUSTIVE_TEST_ORDER; k++) {
                         secp256k1_scalar check_x_s;
-                        r_from_k(&check_x_s, group, k);
+                        r_from_k(&check_x_s, group, k, NULL);
                         if (r_s == check_x_s) {
                             secp256k1_scalar_set_int(&s_times_k_s, k);
                             secp256k1_scalar_mul(&s_times_k_s, &s_times_k_s, &s_s);
@@ -297,7 +297,7 @@ void test_exhaustive_sign(const secp256k1_context *ctx, const secp256k1_ge *grou
                 /* Note that we compute expected_r *after* signing -- this is important
                  * because our nonce-computing function function might change k during
                  * signing. */
-                r_from_k(&expected_r, group, k);
+                r_from_k(&expected_r, group, k, NULL);
                 CHECK(r == expected_r);
                 CHECK((k * s) % EXHAUSTIVE_TEST_ORDER == (i + r * j) % EXHAUSTIVE_TEST_ORDER ||
                       (k * (EXHAUSTIVE_TEST_ORDER - s)) % EXHAUSTIVE_TEST_ORDER == (i + r * j) % EXHAUSTIVE_TEST_ORDER);