# How do I multiply a Twisted Edwards point in Montgomery space?

EdDSA (and ed25519) signatures require a scalar multiplication. Currently, I do this directly in Twisted Edwards space. (The code can be found in my crypto library.) My research and my tests suggest it would be quite a bit faster to do that multiplication in Montgomery space instead. This would mean 4 steps:

1. Convert the Twisted Edwards point to Montgomery space.
2. Perform the Montgomery ladder.
3. Recover the Y coordinate.
4. Convert back to Twisted Edwards space.

Steps 1, 2, and 4 seem pretty simple (even straightforward). Step 3 is more complex, but is documented in this paper, which I have followed. Here is my attempt, which unfortunately gives the wrong result:

typedef i32 fe[10]
typedef struct { fe X; fe Y; fe Z; fe T; } ge;

sv ge_scalarmult(ge *p, const ge *q, const u8 scalar[32])
{
// convert q to montgomery format
fe x1, y1, z1, x2, z2, x3, z3, t1, t2, t3, t4;
fe_sub(z1, q->Z, q->Y);  fe_mul(z1, z1, q->X);  fe_invert(z1, z1);
fe_mul(x1, q->X, t1  );  fe_mul(x1, x1, z1);
fe_mul(y1, q->Z, t1  );  fe_mul(y1, y1, z1);
fe_1(z1); // coherence

// montgomery scalarmult
// The same ladder is used for x25519,
// it comes from the ref10 implementation.
// Field elements are modified in-place
x25519_ladder(x1, x2, z2, x3, z3, scalar);

// recover the y1 coordinate
fe_mul(t1, x1, z2);  // t1 = x1 * z2
fe_add(t2, x2, t1);  // t2 = x2 + t1
fe_sub(t3, x2, t1);  // t3 = x2 − t1
fe_sq (t3, t3);      // t3 = t3^2
fe_mul(t3, t3, x3);  // t3 = t3 * x3
fe_mul973324(t1, z2);// t1 = 2a * z2
fe_add(t2, t2, t1);  // t2 = t2 + t1
fe_mul(t4, x1, x2);  // t4 = x1 * x2
fe_add(t4, t4, z2);  // t4 = t4 + z2
fe_mul(t2, t2, t4);  // t2 = t2 * t4
fe_mul(t1, t1, z2);  // t1 = t1 * z2
fe_sub(t2, t2, t1);  // t2 = t2 − t1
fe_mul(t2, t2, z3);  // t2 = t2 * z3
fe_add(t1, y1, y1);  // t1 = y1 + y1
fe_mul(t1, t1, z2);  // t1 = t1 * z2
fe_mul(t1, t1, z3);  // t1 = t1 * z3
fe_mul(x1, t1, x2);  // x1 = t1 * x2
fe_sub(y1, t2, t3);  // y1 = t2 − t3
fe_mul(z1, t1, z2);  // z1 = t1 * z2

// convert back to twisted edwards
fe_sub(t1  , x1, z1);
fe_add(t2  , x1, z1);
fe_mul(p->X, x1, t2);
fe_mul(p->Y, y1, t1);
fe_mul(p->Z, y1, t2);
fe_mul(p->T, x1, t1);
}


Can someone tell me what I did wrong?

• I can't see anything immediately wrong, but I wonder if you are mixing some coventions from Okeya-Sakurai and another source. Check out Algorithm 5 here. Does that match your algorithm? – user47922 Jul 11 '17 at 16:45
• Scrap the last comment, that was me getting used to the interface… So, I think it matches, if B equals 1 in ed25519. I don't recall precisely, but I think it does. The error doesn't seem to lie there… – Loup Vaillant Jul 11 '17 at 21:27

I finally got the answer from Hacker News. My mistake came from the conversion code. I used a formula from Wikipedia, which is slightly different from the correct one. Simply put, I missed a multiplication by a constant (specifically, the square root of -486664, or 6853475219497561581579357271197624642482790079785650197046958215289687604742 —the other square root also works.)

The correct code for converting to Montgomery format is this (note the inversion, needed to ensure that z1 == 1, which allows a faster Montgomery ladder):

// sqrt(-486664), in 26/25 bit limbs.
fe K = { 54885894, 25242303, 55597453,  9067496, 51808079,
33312638, 25456129, 14121551, 54921728,  3972023 };

// Actual conversion
fe x1, y1, z1, x2, z2, x3, z3, t1, t2, t3, t4;
fe_sub(z1, q->Z, q->Y);  fe_mul(z1, z1, q->X);  fe_invert(z1, z1);
fe_mul(x1, q->X, t1  );  fe_mul(x1, x1, z1);
fe_mul(y1, q->Z, t1  );  fe_mul(y1, y1, z1);
fe_mul(y1, K, y1);  // missing multiplication
fe_1(z1);


(Even though we use projective coordinates, the Montgomery ladder is faster when it can assume z1 is equal to one. This compensates the cost of the inversion.)

The correct code for converting back to Twisted Edwards format is this:

fe_sub(t1  , x1, z1);
fe_add(t2  , x1, z1);
fe_mul(x1  , K , x1);  // missing multiplication
fe_mul(p->X, x1, t2);
fe_mul(p->Y, y1, t1);
fe_mul(p->Z, y1, t2);
fe_mul(p->T, x1, t1);