ElligatorSwift + integrated x-only DH

[sipa]

Builds on top of #979, #1118. Replaces #982.

This implements encoding of curve points using the ElligatorSwift algorithm, using 4 new API calls:

• secp256k1_ellswift_encode, which converts a public key to a 64-byte pseudorandom encoding.
• secp256k1_ellswift_decode, the reverse operation to convert back to normal public keys.
• secp256k1_ellswift_create, which can be seen as a combination of secp256k1_ec_pubkey_create + secp256k1_ellswift_encode, but is somewhat safer.
• secp256k1_ellswift_xdh, which implements x-only Diffie-Hellman directly on top of 64-byte encoded public keys, and more efficiently than decoding + invoking normal ECDH.

The scheme matches that of the SwiftEC paper (https://eprint.iacr.org/2022/759), with two changes (remapping undefined inputs, and encoding the Y parity in the u/t values themselves rather than in a separate bit). To decode an ElligatorSwift 64-byte encoded public key:

• Interpret the encoding as two field elements $(u, t)$ in big-endian 32-byte encoding each, reducing modulo $p$.
• If $u=0$, let $u=1$ instead.
• If $t=0$, let $t=1$ instead.
• If $u^3 + 7 + t^2 = 0$, let $t = 2t$ instead.
• Let $X = \dfrac{u^3 + 7 - t^2}{2t}$.
• Let $Y = \dfrac{X + t}{\sqrt{-3}u}$.
• Let $x$ be the first of $(u + 4Y^2, \dfrac{-X}{2Y} - \dfrac{u}{2}, \dfrac{X}{2Y} - \dfrac{u}{2})$ which is on the curve (at least one will be).
• Let $y$ be the corresponding Y coordinate, with the same parity as (the original) $t$.
• Return $(x,y)$.

This is significantly faster than the Elligator Squared code in #982.

Relevant benchmark (AMD Ryzen 5950X, GCC 12.2.0, default config options; frequency fixed at 2.80 GHz):

Benchmark                     ,    Min(us)    ,    Avg(us)    ,    Max(us)    

ec_keygen                     ,    27.9       ,    28.0       ,    28.0    
ecdh                          ,    53.9       ,    53.9       ,    54.0    
ellswift_encode               ,    24.5       ,    24.5       ,    24.6    
ellswift_decode               ,    11.1       ,    11.1       ,    11.1    
ellswift_keygen               ,    52.5       ,    52.5       ,    52.6    
ellswift_ecdh                 ,    57.9       ,    58.0       ,    58.0    

[real-or-random]

(just skimming)

Do you think there's a nice way to avoid code duplication between divsteps and posdivsteps?

[sipa]

@real-or-random Re posdivsteps see #979 which this PR is based on.

[paulmillr]

Looks awesome! Wanted to clarify this bit, for implementors of EllSwift in other languages / for standardization purpose.

Inputs (u,t) where u=0, t=0, or u^2+t^3+7=0, are remapped to other finite points, rather than outputting infinity.

To which points exactly? It's not immediately obvious from the code.

[sipa]

@paulmillr They are remapped as follows:

• If $t=0$, set $t=1$ instead
• If $u=0$, set $u=1$ instead
• If $u^2+t^3+B=0$, set $t=2t$ instead
• Run the normal algorithm

(for applicable odd-ordered curves with A=0, this covers all otherwise unmapped points, but this remapping doesn't work for every ellswift-compatible curve).

[sipa]

Added test vectors.

[sipa]

I've made a number of improvements:

• Addressed comments so far
• Added test vectors that are verified against the paper author's Sage code (https://github.com/Jchavezsaab/SwiftEC).
• Documented the algorithms line-by-line
• Renamed functions to better match the paper's names.
• Split up the commits
• A few small performance optimizations
• Factored out helper functions for checking if an X coordinate is on the curve (directly, and when given as a fraction).

I think it's ready for more review.

[sipa]

Rebased on updated #979.

[sipa]

Added an explanation of the algorithm and its relation to the paper, in doc/ellswift.md. I've also changed some of the code to better match the algorithm as described in that document (I found some simplifications while writing it).

[theStack]

Verified that the code in the decoding/encoding functions secp256k1_ellswift_xswiftec_frac_var and secp256k1_ellswift_xswiftec_inv_var match the algorithm written in comments each (https://github.com/sipa/secp256k1/blob/55efb1cb525c1b264c57c35203d1b83109948d5e/src/modules/ellswift/main_impl.h#L74-L84 and https://github.com/sipa/secp256k1/blob/55efb1cb525c1b264c57c35203d1b83109948d5e/src/modules/ellswift/main_impl.h#L170-L189). Still planning to do a pen-and-paper review for the algorithm simplifications steps / substitutions pointed out in secp256k1_ellswift_xswiftec_frac_var.

Left two nits below (unrelated to the functions above), regarding API doc and the ellswift benchmark.

[real-or-random]

ACK 55efb1c

Here's a branch with some fixup comments, with commits prepared for git rebase --interactive --autosquash: https://github.com/real-or-random/secp256k1/tree/202206_ellswift_dh
You are welcome to apply or reject; ACK anyway. My initial plan was to have only comment nits, but I found two further refinements on the way, see the commits there.

*edit: It's probably clear to you, but I recommend first taking whatever fixup commits you want from my branch, and then addressing other review commits from @theStack here in order to minimize conflicts.

[sipa]

Squashed in @real-or-random's fixups:

diff --git a/include/secp256k1_ellswift.h b/include/secp256k1_ellswift.h
index b0284d2b..d20fc0cc 100644
--- a/include/secp256k1_ellswift.h
+++ b/include/secp256k1_ellswift.h
@@ -60,11 +60,11 @@ extern "C" {
  *           data:       arbitrary data pointer that is passed through
  */
 typedef int (*secp256k1_ellswift_xdh_hash_function)(
-  unsigned char *output,
-  const unsigned char *x32,
-  const unsigned char *ell_a64,
-  const unsigned char *ell_b64,
-  void *data
+    unsigned char *output,
+    const unsigned char *x32,
+    const unsigned char *ell_a64,
+    const unsigned char *ell_b64,
+    void *data
 );
 
 /** An implementation of an secp256k1_ellswift_xdh_hash_function which uses
diff --git a/src/group_impl.h b/src/group_impl.h
index 8f6e05e0..dcd171f5 100644
--- a/src/group_impl.h
+++ b/src/group_impl.h
@@ -823,8 +823,7 @@ static int secp256k1_ge_is_in_correct_subgroup(const secp256k1_ge* ge) {
 #endif
 }
 
-static int secp256k1_ge_x_on_curve_var(const secp256k1_fe *x)
-{
+static int secp256k1_ge_x_on_curve_var(const secp256k1_fe *x) {
     secp256k1_fe c;
     secp256k1_fe_sqr(&c, x);
     secp256k1_fe_mul(&c, &c, x);
diff --git a/src/modules/ellswift/main_impl.h b/src/modules/ellswift/main_impl.h
index b021be3f..103ab40b 100644
--- a/src/modules/ellswift/main_impl.h
+++ b/src/modules/ellswift/main_impl.h
@@ -24,112 +24,112 @@ static const secp256k1_fe secp256k1_ellswift_c4 = SECP256K1_FE_CONST(0x851695d4,
 static void secp256k1_ellswift_xswiftec_frac_var(secp256k1_fe *xn, secp256k1_fe *xd, const secp256k1_fe *u, const secp256k1_fe *t) {
     /* The implemented algorithm is the following (all operations in GF(p)):
      *
-     * - c0 = sqrt(-3) = 0xa2d2ba93507f1df233770c2a797962cc61f6d15da14ecd47d8d27ae1cd5f852
-     * - If u=0, set u=1.
-     * - If t=0, set t=1.
-     * - If u^3+7+t^2 = 0, set t=2*t.
-     * - Let X=(u^3+7-t^2)/(2*t)
-     * - Let Y=(X+t)/(c0*u)
-     * - If x3=u+4*Y^2 is a valid x coordinate, return x3.
-     * - If x2=(-X/Y-u)/2 is a valid x coordinate, return x2.
-     * - Return x1=(X/Y-u)/2 (which is now guaranteed to be a valid x coordinate).
+     * - Let c0 = sqrt(-3) = 0xa2d2ba93507f1df233770c2a797962cc61f6d15da14ecd47d8d27ae1cd5f852.
+     * - If u = 0, set u = 1.
+     * - If t = 0, set t = 1.
+     * - If u^3+7+t^2 = 0, set t = 2*t.
+     * - Let X = (u^3+7-t^2)/(2*t).
+     * - Let Y = (X+t)/(c0*u).
+     * - If x3 = u+4*Y^2 is a valid x coordinate, return it.
+     * - If x2 = (-X/Y-u)/2 is a valid x coordinate, return it.
+     * - Return x1 = (X/Y-u)/2 (which is now guaranteed to be a valid x coordinate).
      *
      * Introducing s=t^2, g=u^3+7, and simplifying x1=-(x2+u) we get:
      *
      * - Let c0 = ...
-     * - If u=0, set u=1.
-     * - If t=0, set t=1.
-     * - Let s=t^2
-     * - Let g=u^3+7
-     * - If g+s=0, set t=2*t, s=4*s
-     * - Let X=(g-s)/(2*t)
-     * - Let Y=(X+t)/(c0*u) = (g+s)/(2*c0*t*u)
-     * - If x3=u+4*Y^2 is a valid x coordinate, return x3.
-     * - If x2=(-X/Y-u)/2 is a valid x coordinate, return it.
-     * - Return x1=-(x2+u).
+     * - If u = 0, set u = 1.
+     * - If t = 0, set t = 1.
+     * - Let s = t^2
+     * - Let g = u^3+7
+     * - If g+s = 0, set t = 2*t, s = 4*s
+     * - Let X = (g-s)/(2*t).
+     * - Let Y = (X+t)/(c0*u) = (g+s)/(2*c0*t*u).
+     * - If x3 = u+4*Y^2 is a valid x coordinate, return it.
+     * - If x2 = (-X/Y-u)/2 is a valid x coordinate, return it.
+     * - Return x1 = -(x2+u).
      *
      * Now substitute Y^2 = -(g+s)^2/(12*s*u^2) and X/Y = c0*u*(g-s)/(g+s). This
      * means X and Y do not need to be evaluated explicitly anymore.
      *
      * - ...
-     * - If g+s=0, set s=4*s
-     * - If x3=u-(g+s)^2/(3*s*u^2) is a valid x coordinate, return it.
-     * - If x2=(-c0*u*(g-s)/(g+s)-u)/2 is a valid x coordinate, return it.
-     * - Return x1=-(x2+u).
+     * - If g+s = 0, set s = 4*s.
+     * - If x3 = u-(g+s)^2/(3*s*u^2) is a valid x coordinate, return it.
+     * - If x2 = (-c0*u*(g-s)/(g+s)-u)/2 is a valid x coordinate, return it.
+     * - Return x1 = -(x2+u).
      *
      * Simplifying x2 using 2 additional constants:
      *
-     * - c1 = (c0-1)/2 = 0x851695d49a83f8ef919bb86153cbcb16630fb68aed0a766a3ec693d68e6afa40
-     * - c2 = (-c0-1)/2 = 0x7ae96a2b657c07106e64479eac3434e99cf0497512f58995c1396c28719501ee
+     * - Let c1 = (c0-1)/2 = 0x851695d49a83f8ef919bb86153cbcb16630fb68aed0a766a3ec693d68e6afa40.
+     * - Let c2 = (-c0-1)/2 = 0x7ae96a2b657c07106e64479eac3434e99cf0497512f58995c1396c28719501ee.
      * - ...
-     * - If x2=u*(c1*s+c2*g)/(g+s) is a valid x coordinate, return it.
+     * - If x2 = u*(c1*s+c2*g)/(g+s) is a valid x coordinate, return it.
      * - ...
      *
      * Writing x3 as a fraction:
      *
      * - ...
-     * - If x3=(3*s*u^3-(g+s)^2)/(3*s*u^2)
+     * - If x3 = (3*s*u^3-(g+s)^2)/(3*s*u^2) ...
      * - ...
 
      * Overall, we get:
      *
-     * - c1 = 0x851695d49a83f8ef919bb86153cbcb16630fb68aed0a766a3ec693d68e6afa40
-     * - c2 = 0x7ae96a2b657c07106e64479eac3434e99cf0497512f58995c1396c28719501ee
-     * - If u=0, set u=1.
-     * - If t=0, set s=1, else set s=t^2
-     * - Let g=u^3+7
-     * - If g+s=0, set s=4*s
-     * - If x3=(3*s*u^3-(g+s)^2)/(3*s*u^2) is a valid x coordinate, return it.
-     * - If x2=u*(c1*s+c2*g)/(g+s) is a valid x coordinate, return it.
-     * - Return x1=-(x2+u)
+     * - Let c1 = 0x851695d49a83f8ef919bb86153cbcb16630fb68aed0a766a3ec693d68e6afa40.
+     * - Let c2 = 0x7ae96a2b657c07106e64479eac3434e99cf0497512f58995c1396c28719501ee.
+     * - If u = 0, set u = 1.
+     * - If t = 0, set s = 1, else set s = t^2.
+     * - Let g = u^3+7.
+     * - If g+s = 0, set s = 4*s.
+     * - If x3 = (3*s*u^3-(g+s)^2)/(3*s*u^2) is a valid x coordinate, return it.
+     * - If x2 = u*(c1*s+c2*g)/(g+s) is a valid x coordinate, return it.
+     * - Return x1 = -(x2+u).
      */
     secp256k1_fe u1, s, g, p, d, n, l;
     u1 = *u;
     if (EXPECT(secp256k1_fe_normalizes_to_zero_var(&u1), 0)) u1 = secp256k1_fe_one;
     secp256k1_fe_sqr(&s, t);
     if (EXPECT(secp256k1_fe_normalizes_to_zero_var(t), 0)) s = secp256k1_fe_one;
-    secp256k1_fe_sqr(&l, &u1); /* l = u^2 */
-    secp256k1_fe_mul(&g, &l, &u1); /* g = u^3 */
-    secp256k1_fe_add_int(&g, SECP256K1_B); /* g = u^3 + 7 */
-    p = g; /* p = g */
-    secp256k1_fe_add(&p, &s); /* p = g+s */
+    secp256k1_fe_sqr(&l, &u1);                                   /* l = u^2 */
+    secp256k1_fe_mul(&g, &l, &u1);                               /* g = u^3 */
+    secp256k1_fe_add_int(&g, SECP256K1_B);                       /* g = u^3 + 7 */
+    p = g;                                                       /* p = g */
+    secp256k1_fe_add(&p, &s);                                    /* p = g+s */
     if (EXPECT(secp256k1_fe_normalizes_to_zero_var(&p), 0)) {
-        secp256k1_fe_mul_int(&s, 4); /* s = 4*s */
-        /* recompute p = g+s */
-        p = g; /* p = g */
-        secp256k1_fe_add(&p, &s); /* p = g+s */
-    }
-    secp256k1_fe_mul(&d, &s, &l); /* d = s*u^2 */
-    secp256k1_fe_mul_int(&d, 3); /* d = 3*s*u^2 */
-    secp256k1_fe_sqr(&l, &p); /* l = (g+s)^2 */
-    secp256k1_fe_negate(&l, &l, 1); /* l = -(g+s)^2 */
-    secp256k1_fe_mul(&n, &d, &u1); /* n = 3*s*u^3 */
-    secp256k1_fe_add(&n, &l); /* n = 3*s*u^3-(g+s)^2 */
-    if (secp256k1_ge_x_frac_on_curve_var(&n, &d)) {
-        /* Return n/d = (3*s*u^3-(g+s)^2)/(3*s*u^2) */
+        secp256k1_fe_mul_int(&s, 4);
+        /* Recompute p = g+s */
+        p = g;                                                   /* p = g */
+        secp256k1_fe_add(&p, &s);                                /* p = g+s */
+    }                                              
+    secp256k1_fe_mul(&d, &s, &l);                                /* d = s*u^2 */
+    secp256k1_fe_mul_int(&d, 3);                                 /* d = 3*s*u^2 */
+    secp256k1_fe_sqr(&l, &p);                                    /* l = (g+s)^2 */
+    secp256k1_fe_negate(&l, &l, 1);                              /* l = -(g+s)^2 */
+    secp256k1_fe_mul(&n, &d, &u1);                               /* n = 3*s*u^3 */
+    secp256k1_fe_add(&n, &l);                                    /* n = 3*s*u^3-(g+s)^2 */
+    if (secp256k1_ge_x_frac_on_curve_var(&n, &d)) {              
+        /* Return x3 = n/d = (3*s*u^3-(g+s)^2)/(3*s*u^2) */
         *xn = n;
         *xd = d;
         return;
     }
     *xd = p;
-    secp256k1_fe_mul(&l, &secp256k1_ellswift_c1, &s); /* l = c1*s */
-    secp256k1_fe_mul(&n, &secp256k1_ellswift_c2, &g); /* n = c2*g */
-    secp256k1_fe_add(&n, &l); /* n = c1*s+c2*g */
-    secp256k1_fe_mul(&n, &n, &u1); /* n = u*(c1*s+c2*g) */
-    /* Possible optimization: in the invocation below, d^2 = (g+s)^2 is computed,
+    secp256k1_fe_mul(&l, &secp256k1_ellswift_c1, &s);            /* l = c1*s */
+    secp256k1_fe_mul(&n, &secp256k1_ellswift_c2, &g);            /* n = c2*g */
+    secp256k1_fe_add(&n, &l);                                    /* n = c1*s+c2*g */
+    secp256k1_fe_mul(&n, &n, &u1);                               /* n = u*(c1*s+c2*g) */
+    /* Possible optimization: in the invocation below, p^2 = (g+s)^2 is computed,
      * which we already have computed above. This could be deduplicated. */
     if (secp256k1_ge_x_frac_on_curve_var(&n, &p)) {
-        /* Return n/p = u*(c1*s+c2*g)/(g+s) */
+        /* Return x2 = n/p = u*(c1*s+c2*g)/(g+s) */
         *xn = n;
         return;
     }
-    secp256k1_fe_mul(&l, &p, &u1); /* l = u*(g+s) */
-    secp256k1_fe_add(&n, &l); /* n = u*(c1*s+c2*g)+u*(g+s) */
-    secp256k1_fe_negate(xn, &n, 2); /* n = -u*(c1*s+c2*g)-u*(g+s) */
+    secp256k1_fe_mul(&l, &p, &u1);                               /* l = u*(g+s) */
+    secp256k1_fe_add(&n, &l);                                    /* n = u*(c1*s+c2*g)+u*(g+s) */
+    secp256k1_fe_negate(xn, &n, 2);                              /* n = -u*(c1*s+c2*g)-u*(g+s) */
 #ifdef VERIFY
     VERIFY_CHECK(secp256k1_ge_x_frac_on_curve_var(xn, &p));
 #endif
-    /* Return n/p = -(u*(c1*s+c2*g)/(g+s)+u) */
+    /* Return x3 = n/p = -(u*(c1*s+c2*g)/(g+s)+u) */
 }
 
 /** Decode ElligatorSwift encoding (u, t) to X coordinate. */
@@ -182,12 +182,12 @@ static int secp256k1_ellswift_xswiftec_inv_var(secp256k1_fe *t, const secp256k1_
      *   - If s=0, fail.
      *   - Let v=(r/s-u)/2.
      * - Let w=sqrt(s).
-     * - If (c & 5) = 0: return -w*(c3*u + v)
-     * - If (c & 5) = 1: return  w*(c4*u + v)
-     * - If (c & 5) = 4: return  w*(c3*u + v)
-     * - If (c & 5) = 5: return -w*(c4*u + v)
+     * - If (c & 5) = 0: return -w*(c3*u + v).
+     * - If (c & 5) = 1: return  w*(c4*u + v).
+     * - If (c & 5) = 4: return  w*(c3*u + v).
+     * - If (c & 5) = 5: return -w*(c4*u + v).
      */
-    secp256k1_fe x = *x_in, u = *u_in, u2, g, v, s, m, r, q;
+    secp256k1_fe x = *x_in, u = *u_in, g, v, s, m, r, q;
     int ret;
 
     secp256k1_fe_normalize_weak(&x);
@@ -205,28 +205,28 @@ static int secp256k1_ellswift_xswiftec_inv_var(secp256k1_fe *t, const secp256k1_
         /* If -u-x is a valid X coordinate, fail. This would yield an encoding that roundtrips
          * back under the x3 formula instead (which has priority over x1 and x2, so the decoding
          * would not match x). */
-        m = x; /* m = x */
-        secp256k1_fe_add(&m, &u); /* m = u+x */
-        secp256k1_fe_negate(&m, &m, 2); /* m = -u-x */
-        /* test if (-u-x) is a valid X coordinate. If so, fail. */
+        m = x;                                          /* m = x */
+        secp256k1_fe_add(&m, &u);                       /* m = u+x */
+        secp256k1_fe_negate(&m, &m, 2);                 /* m = -u-x */
+        /* Test if (-u-x) is a valid X coordinate. If so, fail. */
         if (secp256k1_ge_x_on_curve_var(&m)) return 0;
 
         /* Let s = -(u^3 + 7)/(u^2 + u*x + x^2) [first part] */
-        secp256k1_fe_sqr(&s, &m); /* s = (u+x)^2 */
-        secp256k1_fe_negate(&s, &s, 1); /* s= -(u+x)^2 */
-        secp256k1_fe_mul(&m, &u, &x); /* m = u*x */
-        secp256k1_fe_add(&s, &m); /* s = -(u^2 + u*x + x^2) */
-
-        /* Note that at this point, s=0 is impossible. If it were the case:
-         *    s = -(u^2 + u*x + x^2) = 0.
-         * => u^2 + u*x + x^2 = 0
-         * => (u + 2*x) * (u^2 + u*x + x^2) = 0
+        secp256k1_fe_sqr(&s, &m);                       /* s = (u+x)^2 */
+        secp256k1_fe_negate(&s, &s, 1);                 /* s = -(u+x)^2 */
+        secp256k1_fe_mul(&m, &u, &x);                   /* m = u*x */
+        secp256k1_fe_add(&s, &m);                       /* s = -(u^2 + u*x + x^2) */
+
+        /* Note that at this point, s = 0 is impossible. If it were the case:
+         *             s = -(u^2 + u*x + x^2) = 0
+         * =>                 u^2 + u*x + x^2 = 0
+         * =>   (u + 2*x) * (u^2 + u*x + x^2) = 0
          * => 2*x^3 + 3*x^2*u + 3*x*u^2 + u^3 = 0
-         * => (x + u)^3 + x^3 = 0
-         * => x^3 = -(x + u)^3
-         * => x^3 + B = (-u - x)^3 + B
+         * =>                 (x + u)^3 + x^3 = 0
+         * =>                             x^3 = -(x + u)^3
+         * =>                         x^3 + B = (-u - x)^3 + B
          *
-         * However, We know x^3 + B is square (because x is on the curve) and
+         * However, we know x^3 + B is square (because x is on the curve) and
          * that (-u-x)^3 + B is not square (the secp256k1_ge_x_on_curve_var(&m)
          * test above would have failed). This is a contradiction, and thus the
          * assumption s=0 is false. */
@@ -237,15 +237,15 @@ static int secp256k1_ellswift_xswiftec_inv_var(secp256k1_fe *t, const secp256k1_
         /* If s is not square, fail. We have not fully computed s yet, but s is square iff
          * -(u^3+7)*(u^2+u*x+x^2) is square (because a/b is square iff a*b is square and b is
          * nonzero). */
-        secp256k1_fe_sqr(&g, &u); /* g = u^2 */
-        secp256k1_fe_mul(&g, &g, &u); /* g = u^3 */
-        secp256k1_fe_add_int(&g, SECP256K1_B); /* g = u^3+7 */
-        secp256k1_fe_mul(&m, &s, &g); /* m = -(u^3 + 7)*(u^2 + u*x + x^2) */
+        secp256k1_fe_sqr(&g, &u);                       /* g = u^2 */
+        secp256k1_fe_mul(&g, &g, &u);                   /* g = u^3 */
+        secp256k1_fe_add_int(&g, SECP256K1_B);          /* g = u^3+7 */
+        secp256k1_fe_mul(&m, &s, &g);                   /* m = -(u^3 + 7)*(u^2 + u*x + x^2) */
         if (!secp256k1_fe_is_square_var(&m)) return 0;
 
         /* Let s = -(u^3 + 7)/(u^2 + u*x + x^2) [second part] */
-        secp256k1_fe_inv_var(&s, &s); /* s = -1/(u^2 + u*x + x^2) [cannot be div by 0] */
-        secp256k1_fe_mul(&s, &s, &g); /* s = -(u^3 + 7)/(u^2 + u*x + x^2) */
+        secp256k1_fe_inv_var(&s, &s);                   /* s = -1/(u^2 + u*x + x^2) [no div by 0] */
+        secp256k1_fe_mul(&s, &s, &g);                   /* s = -(u^3 + 7)/(u^2 + u*x + x^2) */
 
         /* Let v = x. */
         v = x;
@@ -253,61 +253,60 @@ static int secp256k1_ellswift_xswiftec_inv_var(secp256k1_fe *t, const secp256k1_
         /* c is in {2, 3, 6, 7}. In this case we look for an inverse under the x3 formula. */
 
         /* Let s = x-u. */
-        secp256k1_fe_negate(&m, &u, 1); /* m = -u */
-        s = m; /* s = -u */
-        secp256k1_fe_add(&s, &x); /* s = x-u */
+        secp256k1_fe_negate(&m, &u, 1);                 /* m = -u */
+        s = m;                                          /* s = -u */
+        secp256k1_fe_add(&s, &x);                       /* s = x-u */
 
         /* If s is not square, fail. */
         if (!secp256k1_fe_is_square_var(&s)) return 0;
 
         /* Let r = sqrt(-s*(4*(u^3+7)+3*u^2*s)); fail if it doesn't exist. */
-        secp256k1_fe_sqr(&u2, &u); /* u2 = u^2 */
-        secp256k1_fe_mul(&g, &u2, &u); /* g = u^3 */
-        secp256k1_fe_add_int(&g, SECP256K1_B); /* g = u^3+7 */
-        secp256k1_fe_normalize_weak(&g);
-        secp256k1_fe_mul_int(&g, 4); /* g = 4*(u^3+7) */
-        secp256k1_fe_mul_int(&u2, 3); /* u2 = 3*u^2 */
-        secp256k1_fe_mul(&q, &s, &u2); /* q = 3*s*u^2 */
-        secp256k1_fe_add(&q, &g); /* q = 4*(u^3+7)+3*s*u^2 */
-        secp256k1_fe_mul(&q, &q, &s); /* q = s*(4*(u^3+7)+3*u^2*s) */
-        secp256k1_fe_negate(&q, &q, 1); /* q = -s*(4*(u^3+7)+3*u^2*s) */
+        secp256k1_fe_sqr(&g, &u);                       /* g = u^2 */
+        secp256k1_fe_mul(&q, &s, &g);                   /* q = s*u^2 */
+        secp256k1_fe_mul_int(&q, 3);                    /* q = 3*s*u^2 */
+        secp256k1_fe_mul(&g, &g, &u);                   /* g = u^3 */
+        secp256k1_fe_mul_int(&g, 4);                    /* g = 4*u^3 */
+        secp256k1_fe_add_int(&g, 4 * SECP256K1_B);      /* g = 4*(u^3+7) */
+        secp256k1_fe_add(&q, &g);                       /* q = 4*(u^3+7)+3*s*u^2 */
+        secp256k1_fe_mul(&q, &q, &s);                   /* q = s*(4*(u^3+7)+3*u^2*s) */
+        secp256k1_fe_negate(&q, &q, 1);                 /* q = -s*(4*(u^3+7)+3*u^2*s) */
         if (!secp256k1_fe_is_square_var(&q)) return 0;
-        ret = secp256k1_fe_sqrt(&r, &q); /* r = sqrt(-s*(4*(u^3+7)+3*u^2*s)) */
+        ret = secp256k1_fe_sqrt(&r, &q);                /* r = sqrt(-s*(4*(u^3+7)+3*u^2*s)) */
         VERIFY_CHECK(ret);
 
         /* If (c & 1) = 1 and r = 0, fail. */
         if (EXPECT((c & 1) && secp256k1_fe_normalizes_to_zero_var(&r), 0)) return 0;
 
-        /* If s=0, fail. */
+        /* If s = 0, fail. */
         if (EXPECT(secp256k1_fe_normalizes_to_zero_var(&s), 0)) return 0;
 
-        /* Let v=(r/s-u)/2. */
-        secp256k1_fe_inv_var(&v, &s); /* v=1/s [cannot be div by 0] */
-        secp256k1_fe_mul(&v, &v, &r); /* v=r/s */
-        secp256k1_fe_add(&v, &m); /* v=r/s-u */
-        secp256k1_fe_half(&v); /* v=(r/s-u)/2 */
+        /* Let v = (r/s-u)/2. */
+        secp256k1_fe_inv_var(&v, &s);                   /* v = 1/s [no div by 0] */
+        secp256k1_fe_mul(&v, &v, &r);                   /* v = r/s */
+        secp256k1_fe_add(&v, &m);                       /* v = r/s-u */
+        secp256k1_fe_half(&v);                          /* v = (r/s-u)/2 */
     }
 
-    /* Let w=sqrt(s). */
-    ret = secp256k1_fe_sqrt(&m, &s); /* m = sqrt(s) = w */
+    /* Let w = sqrt(s). */
+    ret = secp256k1_fe_sqrt(&m, &s);                    /* m = sqrt(s) = w */
     VERIFY_CHECK(ret);
 
     /* Return logic. */
     if ((c & 5) == 0 || (c & 5) == 5) {
-        secp256k1_fe_negate(&m, &m, 1); /* m = -w */
+        secp256k1_fe_negate(&m, &m, 1);                 /* m = -w */
     }
     /* Now m = {-w if c&5=0 or c&5=5; w otherwise}. */
     secp256k1_fe_mul(&u, &u, c&1 ? &secp256k1_ellswift_c4 : &secp256k1_ellswift_c3);
     /* u = {c4 if c&1=1; c3 otherwise}*u */
-    secp256k1_fe_add(&u, &v); /* u = {c4 if c&1=1; c3 otherwise}*u + v */
+    secp256k1_fe_add(&u, &v);                           /* u = {c4 if c&1=1; c3 otherwise}*u + v */
     secp256k1_fe_mul(t, &m, &u);
     return 1;
 }
 
 /** Use SHA256 as a PRNG, returning SHA256(hasher || cnt).
  *
- * hasher is a SHA256 object which a incrementing 4-byte counter is added to generate randomness.
- * Adding 13 bytes (4 bytes for counter, plus 9 bytes for the SHA256 padding) cannot cross a
+ * hasher is a SHA256 object to which an incrementing 4-byte counter is written to generate randomness.
+ * Writing 13 bytes (4 bytes for counter, plus 9 bytes for the SHA256 padding) cannot cross a
  * 64-byte block size boundary (to make sure it only triggers a single SHA256 compression). */
 static void secp256k1_ellswift_prng(unsigned char* out32, const secp256k1_sha256 *hasher, uint32_t cnt) {
     secp256k1_sha256 hash = *hasher;
diff --git a/src/tests.c b/src/tests.c
index 3ce87bfe..8ada3f86 100644
--- a/src/tests.c
+++ b/src/tests.c
@@ -3945,7 +3945,7 @@ static void test_ge(void) {
         free(ge_set_all);
     }
 
-    /* Test all elements have X coordinates on the curve. */
+    /* Test that all elements have X coordinates on the curve. */
     for (i = 1; i < 4 * runs + 1; i++) {
         secp256k1_fe n;
         CHECK(secp256k1_ge_x_on_curve_var(&ge[i].x));
@@ -3954,7 +3954,7 @@ static void test_ge(void) {
         CHECK(secp256k1_ge_x_frac_on_curve_var(&n, &zf));
     }
 
-    /* Test correspondence secp256k1_ge_x{,_frac}_on_curve_var with ge_set_xo. */
+    /* Test correspondence of secp256k1_ge_x{,_frac}_on_curve_var with ge_set_xo. */
     {
         secp256k1_fe n;
         secp256k1_ge q;

[sipa]

Addressed comments by @theStack and @real-or-random:

diff --git a/include/secp256k1_ellswift.h b/include/secp256k1_ellswift.h
index d20fc0cc..3851f930 100644
--- a/include/secp256k1_ellswift.h
+++ b/include/secp256k1_ellswift.h
@@ -145,7 +145,7 @@ SECP256K1_API int secp256k1_ellswift_decode(
  * followed by secp256k1_ellswift_encode. It is safer, as it uses the secret
  * key as entropy for the encoding (supplemented with auxrnd32, if provided).
  *
- * Like secp256k1_ellswift_encode, this function does not guaranteed that the
+ * Like secp256k1_ellswift_encode, this function does not guarantee that the
  * computed encoding is stable across versions of the library, even if all
  * arguments (including auxrnd32) are the same.
  */
diff --git a/src/modules/ellswift/bench_impl.h b/src/modules/ellswift/bench_impl.h
index 1b0ce6a7..b16a3a36 100644
--- a/src/modules/ellswift/bench_impl.h
+++ b/src/modules/ellswift/bench_impl.h
@@ -94,7 +94,6 @@ void run_ellswift_bench(int iters, int argc, char **argv) {
 
     /* create a context with signing capabilities */
     data.ctx = secp256k1_context_create(SECP256K1_CONTEXT_NONE);
-    memset(data.rnd64, 11, sizeof(data.rnd64));
 
     if (d || have_flag(argc, argv, "ellswift") || have_flag(argc, argv, "encode") || have_flag(argc, argv, "ellswift_encode")) run_benchmark("ellswift_encode", bench_ellswift_encode, bench_ellswift_setup, NULL, &data, 10, iters);
     if (d || have_flag(argc, argv, "ellswift") || have_flag(argc, argv, "decode") || have_flag(argc, argv, "ellswift_decode")) run_benchmark("ellswift_decode", bench_ellswift_decode, bench_ellswift_setup, NULL, &data, 10, iters);
diff --git a/src/modules/ellswift/main_impl.h b/src/modules/ellswift/main_impl.h
index 103ab40b..374c70cd 100644
--- a/src/modules/ellswift/main_impl.h
+++ b/src/modules/ellswift/main_impl.h
@@ -427,6 +427,7 @@ int secp256k1_ellswift_encode(const secp256k1_context *ctx, unsigned char *ell64
         return 1;
     }
     /* Only reached in case the provided pubkey is invalid. */
+    memset(ell64, 0, 64);
     return 0;
 }

[sipa]

And dropped some trailing spaces:

diff --git a/src/modules/ellswift/main_impl.h b/src/modules/ellswift/main_impl.h
index 374c70cd..00bb8a3d 100644
--- a/src/modules/ellswift/main_impl.h
+++ b/src/modules/ellswift/main_impl.h
@@ -98,14 +98,14 @@ static void secp256k1_ellswift_xswiftec_frac_var(secp256k1_fe *xn, secp256k1_fe
         /* Recompute p = g+s */
         p = g;                                                   /* p = g */
         secp256k1_fe_add(&p, &s);                                /* p = g+s */
-    }                                              
+    }
     secp256k1_fe_mul(&d, &s, &l);                                /* d = s*u^2 */
     secp256k1_fe_mul_int(&d, 3);                                 /* d = 3*s*u^2 */
     secp256k1_fe_sqr(&l, &p);                                    /* l = (g+s)^2 */
     secp256k1_fe_negate(&l, &l, 1);                              /* l = -(g+s)^2 */
     secp256k1_fe_mul(&n, &d, &u1);                               /* n = 3*s*u^3 */
     secp256k1_fe_add(&n, &l);                                    /* n = 3*s*u^3-(g+s)^2 */
-    if (secp256k1_ge_x_frac_on_curve_var(&n, &d)) {              
+    if (secp256k1_ge_x_frac_on_curve_var(&n, &d)) {
         /* Return x3 = n/d = (3*s*u^3-(g+s)^2)/(3*s*u^2) */
         *xn = n;
         *xd = d;

[Davidson-Souza]

tACK 90e360a. Full testing backlog:

• Iterated through the implementation, checking all the math (the verbose comments helped a lot!)
• Indepently generated all constants using python/sage (Python 3.10.10 and SageMath version 9.8, Release Date: 2023-02-11)
• A "manual mutation test" for the valgrind ctime tests, adding some dummy, secret dependent branching, it indeed causes memcheck to complain about uninitialized values. Using VALGRIND_MAKE_MEM_DEFINED silences those warnings, as expected.
• Executed all tests, including the memcheck ones (OS and hardware below). Compilation using both gcc and clang, configure script bellow.
• The documentation is sound and very helpful

Test environments:

Config script: [CC=clang] ./configure --with-valgrind=yes --enable-module-ecdh --enable-module-ellswift --enable-experimental --no-recursion

Machine 1: Arch Linux on a AMD Ryze 5 (Linux 6.1.21-hardened1-1-hardened #1 SMP PREEMPT_DYNAMIC x86_64 GNU/Linux )
gcc version 12.2.1 20230201
clang version 14.0.7
GNU libc stable release version 2.37.
valgrind-3.20.0

Machine 2: Raspbian on a Raspbery Pi 4b (Linux raspberrypi 6.1.21-v8+ #1642 SMP PREEMPT aarch64 GNU/Linux)
gcc version 10.2.1 20210110 (Debian 10.2.1-6)
Debian clang version 11.0.1-2
valgrind-3.16.1
(Debian GLIBC 2.31-13+rpt2+rpi1+deb11u5) stable release version 2.31

Benchmarks (only for the new module)

Machine 1, clang

ellswift_encode               ,    23.8       ,    24.1       ,    24.9    
ellswift_decode               ,    10.9       ,    11.4       ,    11.6    
ellswift_keygen               ,    55.3       ,    58.8       ,    61.7    
ellswift_ecdh                 ,    63.3       ,    64.1       ,    65.3 

Machine 1, gcc

ellswift_encode               ,    22.6       ,    22.7       ,    22.9    
ellswift_decode               ,    10.7       ,    10.8       ,    11.0    
ellswift_keygen               ,    54.3       ,    54.7       ,    55.9    
ellswift_ecdh                 ,    62.7       ,    63.2       ,    64.1

Machine 2, clang

ellswift_encode               ,    63.7       ,    63.7       ,    63.8    
ellswift_decode               ,    33.7       ,    33.8       ,    34.0    
ellswift_keygen               ,   177.0       ,   177.0       ,   178.0    
ellswift_ecdh                 ,   233.0       ,   234.0       ,   234.0

Machine 2, gcc

ellswift_encode               ,    66.4       ,    66.4       ,    66.5    
ellswift_decode               ,    34.9       ,    34.9       ,    34.9    
ellswift_keygen               ,   182.0       ,   182.0       ,   183.0    
ellswift_ecdh                 ,   230.0       ,   231.0       ,   231.0  

[real-or-random]

ACK 90e360a

[jonasnick]

ACK 90e360a

[jonasnick]

This has been merged as it received adequate reviews and enough ACKs (as discussed in IRC). Post-merge ACKs are still welcome. If any issues still arise, they can be addressed in follow-up PRs.

[stratospher]

Suggested change

	* Returns: 1: shared secret was succesfully computed
	* Returns: 1: shared secret was successfully computed

[stratospher]

swiftEC paper shows this:
$\dfrac{g(u) + h(u)(Y_0(u) - X_0(u)t)^2}{X_0(u)(1 + h(u)t^2)}$

[stratospher]

tried using WolframAlpha(my first time) for this condition - link.

EDIT: sorry for the spam! got it when i plugged in the right equation.

[stratospher]

(from L332 actually) how does this happen? (when Y = 0 in x1, x2 computation or something else?)

EDIT: L332 needs to be updated to $w(u+2v) = 2X_0(u)$. (this happens when X = X0)

[stratospher]

why not $sign(w)$?

[stratospher]

post-merge ACK 90e360a.
rereviewed the code and skimmed through the swiftEC paper.