add edwards curve and point

src/edwards.rs

@@ -0,0 +1,270 @@
+use crate::traits::{Curve, Point};
+use crypto_bigint::{AddMod, ConstZero, Integer, InvMod, MulMod, NonZero, SubMod, U512, Zero};
+use std::ops::Div;
+
+#[derive(Clone, Copy)]
+struct EdwardsCurve<N>
+where
+    N: AddMod<Output = N> + SubMod<Output = N> + MulMod<Output = N> + InvMod<Output = N> + Copy,
+{
+    pub d: N,
+    pub p: N,
+    generator_x: N,
+    generator_y: N,
+}
+
+/// The edwards curve Ed448-Goldilocks
+const ED_MOD: U512 = U512::ONE
+    .shl(448)
+    .sub_mod(&U512::ONE.shl(224), &U512::MAX)
+    .sub_mod(&U512::ONE, &U512::MAX);
+const ED448_GOLDILOCKS: EdwardsCurve<U512> = EdwardsCurve::new(
+    U512::ZERO.sub_mod(&U512::from_u16(39081), &ED_MOD),
+    ED_MOD,
+    U512::from_be_hex(
+        "0000000000000000aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa955555555555555555555555555555555555555555555555555555555",
+    ),
+    U512::from_be_hex(
+        "0000000000000000ae05e9634ad7048db359d6205086c2b0036ed7a035884dd7b7e36d728ad8c4b80d6565833a2a3098bbbcb2bed1cda06bdaeafbcdea9386ed",
+    ),
+);
+
+impl<N> EdwardsCurve<N>
+where
+    N: AddMod<Output = N>
+        + SubMod<Output = N>
+        + MulMod<Output = N>
+        + InvMod<Output = N>
+        + Copy
+        + ConstZero
+        + crypto_bigint::subtle::ConstantTimeEq
+        + From<u8>,
+{
+    pub const fn new(d: N, p: N, x: N, y: N) -> EdwardsCurve<N> {
+        EdwardsCurve {
+            d,
+            p,
+            generator_x: x,
+            generator_y: y,
+        }
+    }
+}
+
+impl<N> Curve<EdwardsPoint<N>, N> for EdwardsCurve<N>
+where
+    N: AddMod<Output = N>
+        + SubMod<Output = N>
+        + MulMod<Output = N>
+        + InvMod<Output = N>
+        + Integer
+        + Copy
+        + ConstZero
+        + Div<Output = N>,
+{
+    fn point(self, x: N) -> EdwardsPoint<N> {
+        EdwardsPoint {
+            x,
+            y: N::one(),
+            z: N::one(),
+            curve: self,
+        }
+    }
+
+    fn generator(self) -> EdwardsPoint<N> {
+        EdwardsPoint {
+            x: self.generator_x,
+            y: self.generator_y,
+            z: N::one(),
+            curve: self,
+        }
+    }
+}
+
+#[derive(Clone, Copy)]
+pub struct EdwardsPoint<N>
+where
+    N: AddMod<Output = N> + SubMod<Output = N> + MulMod<Output = N> + InvMod<Output = N> + Copy,
+{
+    x: N,
+    y: N,
+    z: N,
+    curve: EdwardsCurve<N>,
+}
+
+impl<N> Point<N> for EdwardsPoint<N>
+where
+    N: AddMod<Output = N>
+        + SubMod<Output = N>
+        + MulMod<Output = N>
+        + InvMod<Output = N>
+        + Copy
+        + crypto_bigint::subtle::ConstantTimeEq
+        + From<u8>
+        + Integer
+        + Div<Output = N>,
+{
+    fn get_x(self) -> N {
+        self.x.mul_mod(
+            &self.z.inv_mod(&self.curve.p).unwrap(),
+            &NonZero::new(self.curve.p).unwrap(),
+        )
+    }
+
+    fn add(&self, rhs: &Self, _: &Self) -> Self {
+        let a = self.z.mul_mod(&rhs.z, &self.curve.p);
+        let b = a.mul_mod(&a, &self.curve.p);
+        let c = self.x.mul_mod(&rhs.x, &self.curve.p);
+        let d = self.y.mul_mod(&rhs.y, &self.curve.p);
+        let e = self
+            .curve
+            .d
+            .mul_mod(&c.mul_mod(&d, &self.curve.p), &self.curve.p);
+        let f = b.sub_mod(&e, &self.curve.p);
+        let g = b.add_mod(&e, &self.curve.p);
+
+        let x = a.mul_mod(&f, &self.curve.p).mul_mod(
+            &self
+                .x
+                .add_mod(&self.y, &self.curve.p)
+                .mul_mod(&rhs.x.add_mod(&rhs.y, &self.curve.p), &self.curve.p)
+                .sub_mod(&c, &self.curve.p)
+                .sub_mod(&d, &self.curve.p),
+            &self.curve.p,
+        );
+        let y = a
+            .mul_mod(&g, &self.curve.p)
+            .mul_mod(&d.sub_mod(&c, &self.curve.p), &self.curve.p);
+        let z = g.mul_mod(&f, &self.curve.p);
+
+        EdwardsPoint {
+            x,
+            y,
+            z,
+            curve: self.curve,
+        }
+    }
+
+    fn double(&self) -> Self {
+        /*
+        let a = self.x.add_mod(&self.y, &self.curve.p);
+        let b = a.mul_mod(&a, &self.curve.p);
+        let c = self.x.mul_mod(&self.x, &self.curve.p);
+        let d = self.y.mul_mod(&self.y, &self.curve.p);
+        let e = c.add_mod(&d, &self.curve.p);
+        let h = self.z.mul_mod(&self.z, &self.curve.p);
+        let j = e.sub_mod(&h.mul_mod(&N::from(2u8), &self.curve.p), &self.curve.p);
+
+        let x = b.sub_mod(&e, &self.curve.p);
+        let y = e.mul_mod(&c.sub_mod(&d, &self.curve.p), &self.curve.p);
+        let z = e.mul_mod(&j, &self.curve.p);
+
+        EdwardsPoint {
+            x,
+            y,
+            z,
+            curve: self.curve,
+        }
+         */
+
+        self.add(self, self)
+    }
+
+    fn mul(&self, x: &N) -> Self {
+        self.ct_ladder(x)
+    }
+
+    fn neg(&self) -> Self {
+        EdwardsPoint {
+            x: N::zero().sub_mod(&self.x, &self.curve.p),
+            y: self.y,
+            z: self.z,
+            curve: self.curve,
+        }
+    }
+}
+
+impl<N> EdwardsPoint<N>
+where
+    N: AddMod<Output = N>
+        + SubMod<Output = N>
+        + MulMod<Output = N>
+        + InvMod<Output = N>
+        + Copy
+        + ConstZero
+        + crypto_bigint::subtle::ConstantTimeEq,
+{
+    fn get_y(self) -> N {
+        self.y
+            .mul_mod(&self.z.inv_mod(&self.curve.p).unwrap(), &self.curve.p)
+    }
+}
+
+ladder_impl!(EdwardsPoint, swap_edwards);
+
+fn swap_edwards<N>(a: &mut EdwardsPoint<N>, b: &mut EdwardsPoint<N>, c: u8, max: N)
+where
+    N: AddMod + SubMod + MulMod + InvMod<Output = N> + Copy + Integer + std::ops::Div<Output = N>,
+{
+    let m = max * N::from(c);
+
+    let x = m & (a.x ^ b.x);
+    let y = m & (a.y ^ b.y);
+    let z = m & (a.z ^ b.z);
+
+    a.x ^= x;
+    b.x ^= x;
+
+    a.y ^= y;
+    b.y ^= y;
+
+    a.z ^= z;
+    b.z ^= z;
+}
+
+#[cfg(test)]
+mod test {
+    use super::*;
+
+    #[test]
+    fn double() {
+        let x = U512::from_be_hex(
+            "000000000000000049dcbc5c6c0cce2c1419a17226f929ea255a09cf4e0891c693fda4be70c74cc301b7bdf1515dd8ba21aee1798949e120e2ce42ac48ba7f30",
+        );
+        let y = U512::from_be_hex(
+            "0000000000000000d49077e4accde527164b33a5de021b979cb7c02f0457d845c90dc3227b8a5bc1c0d8f97ea1ca9472b5d444285d0d4f5b32e236f86de51839",
+        );
+        let p = ED448_GOLDILOCKS.generator().double();
+
+        assert_eq!(p.get_x(), x);
+        assert_eq!(p.get_y(), y);
+    }
+
+    #[test]
+    fn explicit_double() {
+        let x = U512::from_be_hex(
+            "000000000000000049dcbc5c6c0cce2c1419a17226f929ea255a09cf4e0891c693fda4be70c74cc301b7bdf1515dd8ba21aee1798949e120e2ce42ac48ba7f30",
+        );
+        let y = U512::from_be_hex(
+            "0000000000000000d49077e4accde527164b33a5de021b979cb7c02f0457d845c90dc3227b8a5bc1c0d8f97ea1ca9472b5d444285d0d4f5b32e236f86de51839",
+        );
+        let p = ED448_GOLDILOCKS.generator();
+        let p = p.add(&p, &p);
+
+        assert_eq!(p.get_x(), x);
+        assert_eq!(p.get_y(), y);
+    }
+
+    #[test]
+    fn add() {
+        let p = ED448_GOLDILOCKS.generator();
+        let q = p.neg();
+
+        let x = U512::ZERO;
+        let y = U512::ONE;
+
+        let sum = p.add(&q, &p);
+
+        assert_eq!(sum.get_x(), x);
+        assert_eq!(sum.get_y(), y);
+    }
+}

src/lib.rs

@@ -1,231 +1,18 @@
+#[macro_use]
+mod macros;
+
+mod edwards;
+mod montgomery;
+mod traits;
+
 pub use crypto_bigint::{ConstChoice, Encoding, NonZero, U256};
 use std::fmt::{Display, Formatter};
-use std::ops::{Div, Mul};
 
-#[cfg(feature = "rand")] use crypto_bigint::{rand_core::OsRng, Random};
+#[cfg(feature = "rand")]
-
+use crypto_bigint::{Random, rand_core::OsRng};
-
-/// By^2 = x^3 + Ax^2 + x
-#[derive(Clone, Copy, Debug)]
-pub struct MontgomeryCurve {
-    pub a: U256,
-    pub b: U256,
-    pub p: U256,
-}
-
-impl MontgomeryCurve {
-    pub const fn new(a: U256, b: U256, p: U256) -> MontgomeryCurve {
-        MontgomeryCurve { a, b, p }
-    }
-
-    pub fn y(&self, x: U256) -> U256 {
-        ((x*x*x + self.a*x*x + x)/self.b).sqrt()
-    }
-
-    pub fn point(self, x: U256) -> MontgomeryPoint {
-        MontgomeryPoint{
-            x,
-            z: U256::ONE,
-            curve: self,
-        }
-    }
-}
-
-pub fn clamp_u256(x: U256) -> U256 {
-    let mut bytes: [u8; 32] = x.to_le_bytes();
-
-    // clamp value
-    bytes[0] &= 248;
-    bytes[31] &= 127;
-    bytes[31] |= 64;
-
-    U256::from_le_bytes(bytes)
-}
-
-#[derive(Clone, Copy, Debug)]
-pub struct MontgomeryPoint {
-    x: U256,
-    z: U256,
-
-    curve: MontgomeryCurve,
-}
-
-impl MontgomeryPoint {
-    pub fn add(&self, rhs: MontgomeryPoint, neg: MontgomeryPoint) -> MontgomeryPoint {
-        let p = &NonZero::new(self.curve.p).unwrap();
-
-        let mut v_0 = self.x.add_mod(&self.z, p); // 1: V_0 = X_P + Z_P
-        let mut v_1 = rhs.x.sub_mod(&rhs.z, p);   // 2: V_1 = X_Q - Z_Q
-        v_1 = v_1.mul_mod(&v_0, p);               // 3: V_1 = V_1 * V_0
-
-        v_0 = self.x.sub_mod(&self.z, p);       // 4: V_0 = X_P - Z_P
-        let mut v_2 = rhs.x.add_mod(&rhs.z, p); // 5: V_2 = X_Q + Z_Q
-        v_2 = v_2.mul_mod(&v_0, p);             // 6: V_2 = V_2 * V_0
-
-        let mut v_3 = v_1.add_mod(&v_2, p); // 7: V_3 = V_1 + V_2
-        v_3 = v_3.mul_mod(&v_3, p);         // 8: V_3 = V_3^2
-
-        let mut v_4 = v_1.sub_mod(&v_2, p); // 9: V_4 = V_1 - V_2
-        v_4 = v_4.mul_mod(&v_4, p);         // 10: v_4 = v_4^2
-
-        let x = neg.z.mul_mod(&v_3, p); // 11: X_⨁ = Z_⊝ * V_3
-        let z = neg.x.mul_mod(&v_4, p); // 12: Z_⨁ = X_⊝ * V_4
-
-        MontgomeryPoint {
-            x,
-            z,
-            curve: self.curve,
-        }
-    }
-
-    pub fn double(&self) -> MontgomeryPoint {
-        let p = &NonZero::new(self.curve.p).unwrap();
-
-        let mut v_1 = self.x.add_mod(&self.z, p);
-        v_1 = v_1.mul_mod(&v_1, p);
-
-        let mut v_2 = self.x.sub_mod(&self.z, p);
-        v_2 = v_2.mul_mod(&v_2, p);
-
-        let x_2p = v_1.mul_mod(&v_2, p);
-
-        v_1 = v_1.sub_mod(&v_2, p);
-        let mut v_3 = self.curve.a.add_mod(&U256::from_u8(2), p).div(U256::from_u8(4)).mul_mod(&v_1, p);
-
-        v_3 = v_3.add_mod(&v_2, p);
-
-        let z_2p = v_1.mul_mod(&v_3, p);
-
-        MontgomeryPoint {
-            x: x_2p,
-            z: z_2p,
-            curve: self.curve,
-        }
-    }
-
-    pub fn get_x(&self) -> U256 {
-        self.x.mul_mod(
-            &self.z.inv_mod(&self.curve.p).unwrap(),
-            &NonZero::new(self.curve.p).unwrap()
-        )
-    }
-
-    pub fn ladder(self, rhs: U256) -> MontgomeryPoint {
-        let mut x_0 = self;
-        let mut x_1 = x_0.double();
-
-        for i in (0..rhs.bits() - 1).rev() {
-            if rhs.bit(i).eq(&ConstChoice::FALSE) {
-                x_1 = x_0.add(x_1, self);
-                x_0 = x_0.double();
-            } else {
-                x_0 = x_0.add(x_1, self);
-                x_1 = x_1.double();
-            }
-        }
-
-        x_0
-    }
-
-    pub fn ct_ladder(self, rhs: U256) -> MontgomeryPoint {
-        let mut x_0 = self;
-        let mut x_1 = x_0.double();
-
-        let mut prev: u8 = 0;
-        for i in (0..rhs.bits() - 1).rev() {
-            let curr: u8 = bool::from(rhs.bit(i)).into();
-            swap(&mut x_0, &mut x_1, curr ^ prev);
-            prev = curr;
-
-            x_1 = x_0.add(x_1, self);
-            x_0 = x_0.double();
-        }
-        swap(&mut x_0, &mut x_1, prev);
-
-        x_0
-    }
-}
-
-const SWAP_MASKS: [U256; 2] = [U256::ZERO, U256::from_be_hex("FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFF")];
-fn swap(a: &mut MontgomeryPoint, b: &mut MontgomeryPoint, c: u8) {
-    let m = SWAP_MASKS[c as usize];
-
-    let x = m & (a.x ^ b.x);
-    let z = m & (a.z ^ b.z);
-
-    *a = MontgomeryPoint{
-        x: a.x ^ x,
-        z: a.z ^ z,
-        curve: a.curve,
-    };
-
-    *b = MontgomeryPoint{
-        x: b.x ^ x,
-        z: b.z ^ z,
-        curve: b.curve,
-    }
-}
-
-impl Display for MontgomeryPoint {
-    fn fmt(&self, f: &mut Formatter<'_>) -> std::fmt::Result {
-        f.write_fmt(format_args!("({} : {})", self.x, self.z))
-    }
-}
-
-impl Mul<U256> for MontgomeryPoint {
-    type Output = MontgomeryPoint;
-
-    /// (not)constant-time montgomery ladder
-    fn mul(self, rhs: U256) -> Self::Output {
-        self.ct_ladder(rhs)
-    }
-}
-
-/*
-fn montgomery_ladder(x: U256, n: U256) -> U256 {
-    let mut x_1 = x;
-    let mut x_2 = U256::ONE;
-    let mut z_2 = U256::ZERO;
-    let mut x_3 = x;
-    let mut z_3 = U256::ONE;
-
-    let mut prevbit = 0u8;
-
-    for i in (0..(n.bits()-1)).rev() {
-        let bit = bool::from(n.bit(i)) as u8;
-        let b = bit ^ prevbit;
-        prevbit = bit;
-
-        // CSwap
-        match b & 1 == 0 {
-            true => {
-                swap(&mut x_2, &mut x_3);
-                swap(&mut z_2, &mut z_3);
-            }
-            false => {
-                swap(&mut x_2, &mut x_3);
-                swap(&mut x_2, &mut x_3);
-            }
-        }
-
-        //ladder_step(&mut x_2, &mut z_2, &mut x_3, &mut z_3, x_1);
-    }
-
-    U256::ZERO
-}
- */
-
-
-/// # D. J. Bernstein's Curve25519
-/// from the curve: \
-/// $y^2 = x*(x^2 + 486662x + 1)$ \
-/// where A = 486662, B = 1, p = 2^255 - 19
-pub const CURVE_25519: MontgomeryCurve = MontgomeryCurve::new(
-    U256::from_u32(486662u32),
-    U256::ONE,
-    U256::ONE.shl(255).sub_mod(&U256::from_u8(19), &U256::MAX),
-);
 
+use crate::traits::{Curve, Point};
+pub(crate) use montgomery::*;
 
 #[derive(Debug, Clone, Copy, Eq, PartialEq)]
 pub struct Public(U256);
@@ -287,15 +74,17 @@ impl From<U256> for Secret {
 #[cfg(test)]
 mod tests {
     use super::*;
-
+    const BIG_VALID_VALUE: U256 = U256::ONE
-    const BIG_VALID_VALUE: U256 = U256::ONE.shl(251).add_mod(&U256::from_u64(1232039487203*8), &CURVE_25519.p);
+        .shl(251)
+        .add_mod(&U256::from_u64(1232039487203 * 8), &CURVE_25519.p);
 
     #[test]
     fn public_key_derivation() {
         let a = Secret(U256::from_u8(3));
         let p = Public::from(&a);
 
-        let x_k = U256::from_be_hex("1c12bc1a6d57abe645534d91c21bba64f8824e67621c0859c00a03affb713c12");
+        let x_k =
+            U256::from_be_hex("1c12bc1a6d57abe645534d91c21bba64f8824e67621c0859c00a03affb713c12");
         assert_eq!(p.0, x_k);
     }
 
@@ -306,7 +95,8 @@ mod tests {
 
         let x = s.diffie_hellman(&p);
 
-        let x_k = U256::from_be_hex("0933dc6ed7122bdf2a5f1b0516e218b743868769778787cfa1ac0ae2dfa89891");
+        let x_k =
+            U256::from_be_hex("0933dc6ed7122bdf2a5f1b0516e218b743868769778787cfa1ac0ae2dfa89891");
         assert_eq!(x.0, x_k);
     }
 
@@ -325,160 +115,4 @@ mod tests {
 
         assert_eq!(sa, sb);
     }
-
-    #[test]
-    fn direct_diffie_hellman() {
-        let o = CURVE_25519.point(U256::from_u8(9));
-
-        // multiples of 8
-        let a = clamp_u256(U256::from_u32(1 << 9));
-        let b = clamp_u256(U256::from_u32(1 << 8));
-
-        let pa = o.ladder(a).get_x();
-        let pb = o.ladder(b).get_x();
-
-        let pap = CURVE_25519.point(pa);
-        let sa = pap.ladder(b).get_x();
-
-        let pbp = CURVE_25519.point(pb);
-        let sb = pbp.ladder(a).get_x();
-
-        assert_eq!(sa, sb);
-    }
-
-    #[test]
-    fn montgomery_point_double() {
-        let point = CURVE_25519.point(U256::from_u8(2));
-
-        let dbl = point.double();
-
-        let x_k = U256::from_be_hex("3c0c7855eae09525fb741a3649b05a2da70d0e6b0ff35660a3feba622566fc9c");
-
-        assert_eq!(dbl.get_x(), x_k);
-    }
-
-    #[test]
-    fn montgomery_point_big_double() {
-        let point = CURVE_25519.point(BIG_VALID_VALUE);
-
-        let dbl = point.double();
-
-        let x_k = U256::from_be_hex("3f4d5612f4df4042680d41c7f5d9673fef3590c0c40bc04959848a88a117f311");
-
-        assert_eq!(dbl.get_x(), x_k);
-    }
-
-    #[test]
-    fn montgomery_point_get_x() {
-        let point = MontgomeryPoint{
-            x: U256::from_u8(2),
-            z: U256::from_u8(3),
-            curve: CURVE_25519,
-        };
-
-        let k_x = U256::from_be_hex("2aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa5");
-        assert_eq!(point.get_x(), k_x);
-    }
-
-    #[test]
-    fn montgomery_point_add() {
-        let point_p = CURVE_25519.point(U256::from_u8(2));
-        let point_q = point_p.double();
-
-        let sum = point_p.add(point_q, point_p);
-
-        let x_k = U256::from_be_hex("56f0b3e7e53fed658bf39e9f8691055ba8a58a935482b00b21b5678c16261fec");
-        assert_eq!(sum.get_x(), x_k);
-    }
-
-    #[test]
-    fn montgomery_point_big_add() {
-        let point_p = CURVE_25519.point(BIG_VALID_VALUE);
-        let point_q = point_p.double();
-
-        let sum = point_p.add(point_q, point_p);
-
-        let x_k = U256::from_be_hex("1431b13075c5e0c75f2aad0f604f7286fbeac04cdc6099a63c54cdd96cbd3333");
-        assert_eq!(sum.get_x(), x_k);
-    }
-
-    #[test]
-    fn montgomery_ladder_double() {
-        let p = MontgomeryPoint{
-            x: U256::from_u8(2),
-            z: U256::from_u8(3),
-            curve: CURVE_25519,
-        };
-
-        let n = U256::from_u8(2);
-        let a = p * n;
-        let b = p.double();
-
-        assert_eq!(a.get_x(), b.get_x());
-    }
-
-    #[test]
-    fn montgomery_ladder() {
-        let p = CURVE_25519.point(U256::from_u8(2));
-
-        let n = U256::from_u8(4);
-
-        let a = p.ladder(n);
-
-        let x_k = U256::from_be_hex("2eceef1936e6df00c49e7aedac94446cc3b156165b50f247a15fdcee5e065582");
-        assert_eq!(a.get_x(), x_k);
-    }
-
-    #[test]
-    fn large_montgomery_ladder() {
-        let p = CURVE_25519.point(U256::from_u8(9));
-
-        let n = BIG_VALID_VALUE;
-
-        let a = p.ladder(n);
-
-        let x_k = U256::from_be_hex("0b379a815f36aca005dd1b19b1d483cd73fad06225c6ad927fde9316214a68c1");
-        assert_eq!(a.get_x(), x_k);
-    }
-
-    #[test]
-    fn constant_time_montgomery_ladder() {
-        let p = CURVE_25519.point(U256::from_u8(9));
-
-        let n = BIG_VALID_VALUE;
-
-        let a = p.ct_ladder(n);
-        let b = p.ladder(n);
-
-        assert_eq!(b.get_x(), a.get_x());
-    }
-
-    #[test]
-    fn montgomery_ladder_commutative() {
-        let p = CURVE_25519.point(U256::from_u8(2));
-
-        let n_1 = U256::from_u8(3);
-        let n_2 = U256::from_u8(4);
-
-        let a = p * n_1 * n_2;
-        let b = p * n_2 * n_1;
-
-        assert_eq!(a.get_x(), b.get_x());
-    }
-
-    #[test]
-    fn montgomery_point_conversion() {
-        let p = CURVE_25519.point(U256::from_u8(9));
-
-        let n_1 = U256::from_u8(3);
-        let n_2 = U256::from_u8(4);
-
-        let a = p * n_1;
-
-        let x_a = a.get_x();
-        let pa = CURVE_25519.point(x_a) * n_2;
-
-
-        assert_eq!(pa.get_x(), (a * n_2).get_x());
-    }
 }

src/montgomery.rs

@@ -0,0 +1,394 @@
+use crate::traits::{Curve, Point};
+use crypto_bigint::subtle::{Choice, ConstantTimeEq};
+use crypto_bigint::{AddMod, ConstChoice, Encoding, NonZero, SubMod, U256};
+use crypto_bigint::{Integer, InvMod, MulMod};
+use std::fmt::{Display, Formatter};
+use std::ops::{Add, Div, Mul};
+
+/// A curve of the form `By^2 = x^3 + Ax^2 + x`
+#[derive(Clone, Copy, Debug)]
+pub struct MontgomeryCurve<N>
+where
+    N: AddMod + SubMod + MulMod + InvMod + Integer + Copy,
+{
+    pub a: N,
+    pub b: N,
+    pub p: N,
+}
+
+impl<N> MontgomeryCurve<N>
+where
+    N: AddMod + SubMod + MulMod + InvMod + Copy + Integer,
+{
+    pub const fn new(a: N, b: N, p: N) -> MontgomeryCurve<N> {
+        MontgomeryCurve { a, b, p }
+    }
+
+    /*/// Returns the absolute y-value of a point on the curve at a given x-coordinate.
+    pub fn y(&self, x: N) -> N
+    where    <N as Mul>::Output: Mul<N>,
+            <<N as Mul>::Output as Mul<N>>::Output: Add,
+           <<<N as Mul>::Output as Mul<N>>::Output as Add>::Output: Add<N>,
+          <<<<N as Mul>::Output as Mul<N>>::Output as Add>::Output as Add<N>>::Output: Div<N>
+    {
+        ((x*x*x + self.a*x*x + x)/self.b).sqrt()
+    }*/
+}
+
+impl<N> Curve<MontgomeryPoint<N>, N> for MontgomeryCurve<N>
+where
+    N: AddMod + SubMod + MulMod + InvMod<Output = N> + Copy + Integer + Div<Output = N>,
+{
+    /// Get the point at an x-coordinate with parameters from this curve.
+    fn point(self, x: N) -> MontgomeryPoint<N> {
+        MontgomeryPoint {
+            x,
+            z: N::one(),
+            curve: self,
+        }
+    }
+
+    fn generator(self) -> MontgomeryPoint<N> {
+        self.point(N::from(9u8))
+    }
+}
+
+pub fn clamp_u256(x: U256) -> U256 {
+    let mut bytes: [u8; 32] = x.to_le_bytes();
+
+    // clamp value
+    bytes[0] &= 248;
+    bytes[31] &= 127;
+    bytes[31] |= 64;
+
+    U256::from_le_bytes(bytes)
+}
+
+#[derive(Copy, Debug, Clone)]
+pub struct MontgomeryPoint<N>
+where
+    N: AddMod + SubMod + Div<Output = N> + MulMod<Output = N> + InvMod<Output = N> + Copy + Integer,
+{
+    pub x: N,
+    pub z: N,
+
+    pub curve: MontgomeryCurve<N>,
+}
+
+impl<N> Point<N> for MontgomeryPoint<N>
+where
+    N: AddMod
+        + SubMod
+        + Div<Output = N>
+        + MulMod<Output = N>
+        + InvMod<Output = N>
+        + Copy
+        + Integer
+        + From<u8>,
+{
+    fn get_x(self) -> N {
+        self.x.mul_mod(
+            &self.z.inv_mod(&self.curve.p).unwrap(),
+            &NonZero::new(self.curve.p).unwrap(),
+        )
+    }
+
+    fn add(&self, rhs: &Self, neg: &Self) -> Self {
+        let p = &NonZero::new(self.curve.p).unwrap();
+
+        let mut v_0 = self.x.add_mod(&self.z, p); // 1: V_0 = X_P + Z_P
+        let mut v_1 = rhs.x.sub_mod(&rhs.z, p); // 2: V_1 = X_Q - Z_Q
+        v_1 = v_1.mul_mod(&v_0, p); // 3: V_1 = V_1 * V_0
+
+        v_0 = self.x.sub_mod(&self.z, p); // 4: V_0 = X_P - Z_P
+        let mut v_2 = rhs.x.add_mod(&rhs.z, p); // 5: V_2 = X_Q + Z_Q
+        v_2 = v_2.mul_mod(&v_0, p); // 6: V_2 = V_2 * V_0
+
+        let mut v_3 = v_1.add_mod(&v_2, p); // 7: V_3 = V_1 + V_2
+        v_3 = v_3.mul_mod(&v_3, p); // 8: V_3 = V_3^2
+
+        let mut v_4 = v_1.sub_mod(&v_2, p); // 9: V_4 = V_1 - V_2
+        v_4 = v_4.mul_mod(&v_4, p); // 10: v_4 = v_4^2
+
+        let x = neg.z.mul_mod(&v_3, p); // 11: X_⨁ = Z_⊝ * V_3
+        let z = neg.x.mul_mod(&v_4, p); // 12: Z_⨁ = X_⊝ * V_4
+
+        MontgomeryPoint {
+            x,
+            z,
+            curve: self.curve,
+        }
+    }
+
+    fn double(&self) -> MontgomeryPoint<N> {
+        let p = &NonZero::new(self.curve.p).unwrap();
+
+        let mut v_1 = self.x.add_mod(&self.z, p);
+        v_1 = v_1.mul_mod(&v_1, p);
+
+        let mut v_2 = self.x.sub_mod(&self.z, p);
+        v_2 = v_2.mul_mod(&v_2, p);
+
+        let x_2p = v_1.mul_mod(&v_2, p);
+
+        v_1 = v_1.sub_mod(&v_2, p);
+        let mut v_3 = self
+            .curve
+            .a
+            .add_mod(&N::from(2u8), p)
+            .div(N::from(4u8))
+            .mul_mod(&v_1, p);
+
+        v_3 = v_3.add_mod(&v_2, p);
+
+        let z_2p = v_1.mul_mod(&v_3, p);
+
+        MontgomeryPoint {
+            x: x_2p,
+            z: z_2p,
+            curve: self.curve,
+        }
+    }
+
+    fn mul(&self, other: &N) -> Self {
+        self.ct_ladder(other)
+    }
+
+    fn neg(&self) -> Self {
+        *self
+    }
+}
+
+ladder_impl!(MontgomeryPoint, swap_montgomery);
+
+fn swap_montgomery<N>(a: &mut MontgomeryPoint<N>, b: &mut MontgomeryPoint<N>, c: u8, max: N)
+where
+    N: AddMod + SubMod + MulMod + InvMod<Output = N> + Copy + Integer + std::ops::Div<Output = N>,
+{
+    let m = max * N::from(c);
+
+    let x = m & (a.x ^ b.x);
+    let z = m & (a.z ^ b.z);
+
+    *a = MontgomeryPoint {
+        x: a.x ^ x,
+        z: a.z ^ z,
+        curve: a.curve,
+    };
+
+    *b = MontgomeryPoint {
+        x: b.x ^ x,
+        z: b.z ^ z,
+        curve: b.curve,
+    }
+}
+
+impl<N> Display for MontgomeryPoint<N>
+where
+    N: AddMod
+        + SubMod
+        + MulMod
+        + InvMod<Output = N>
+        + Copy
+        + Integer
+        + Display
+        + std::ops::Div<Output = N>,
+{
+    fn fmt(&self, f: &mut Formatter<'_>) -> std::fmt::Result {
+        f.write_fmt(format_args!("({} : {})", self.x, self.z))
+    }
+}
+
+impl<N> Mul<N> for MontgomeryPoint<N>
+where
+    N: AddMod + SubMod + MulMod + InvMod<Output = N> + Copy + Integer + std::ops::Div<Output = N>,
+{
+    type Output = MontgomeryPoint<N>;
+
+    /// constant-time montgomery ladder
+    fn mul(self, rhs: N) -> Self::Output {
+        self.ct_ladder(&rhs)
+    }
+}
+
+/// # D. J. Bernstein's Curve25519
+/// from the curve: \
+/// $y^2 = x*(x^2 + 486662x + 1)$ \
+/// where A = 486662, B = 1, p = 2^255 - 19
+pub const CURVE_25519: MontgomeryCurve<U256> = MontgomeryCurve::new(
+    U256::from_u32(486662u32),
+    U256::ONE,
+    U256::ONE.shl(255).sub_mod(&U256::from_u8(19), &U256::MAX),
+);
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+    use crate::traits::{Curve, Point};
+
+    const BIG_VALID_VALUE: U256 = U256::ONE
+        .shl(251)
+        .add_mod(&U256::from_u64(1232039487203 * 8), &CURVE_25519.p);
+
+    #[test]
+    fn direct_diffie_hellman() {
+        let o = CURVE_25519.point(U256::from_u8(9));
+
+        // multiples of 8
+        let a = clamp_u256(U256::from_u32(1 << 9));
+        let b = clamp_u256(U256::from_u32(1 << 8));
+
+        let pa = o.ladder(&a).get_x();
+        let pb = o.ladder(&b).get_x();
+
+        let pap = CURVE_25519.point(pa);
+        let sa = pap.ladder(&b).get_x();
+
+        let pbp = CURVE_25519.point(pb);
+        let sb = pbp.ladder(&a).get_x();
+
+        assert_eq!(sa, sb);
+    }
+
+    #[test]
+    fn montgomery_point_double() {
+        let point = CURVE_25519.point(U256::from_u8(2));
+
+        let dbl = point.double();
+
+        let x_k =
+            U256::from_be_hex("3c0c7855eae09525fb741a3649b05a2da70d0e6b0ff35660a3feba622566fc9c");
+
+        assert_eq!(dbl.get_x(), x_k);
+    }
+
+    #[test]
+    fn montgomery_point_big_double() {
+        let point = CURVE_25519.point(BIG_VALID_VALUE);
+
+        let dbl = point.double();
+
+        let x_k =
+            U256::from_be_hex("3f4d5612f4df4042680d41c7f5d9673fef3590c0c40bc04959848a88a117f311");
+
+        assert_eq!(dbl.get_x(), x_k);
+    }
+
+    #[test]
+    fn montgomery_point_get_x() {
+        let point = MontgomeryPoint {
+            x: U256::from_u8(2),
+            z: U256::from_u8(3),
+            curve: CURVE_25519,
+        };
+
+        let k_x =
+            U256::from_be_hex("2aaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaaa5");
+        assert_eq!(point.get_x(), k_x);
+    }
+
+    #[test]
+    fn montgomery_point_add() {
+        let point_p = CURVE_25519.point(U256::from_u8(2));
+        let point_q = point_p.double();
+
+        let sum = point_p.add(&point_q, &point_p);
+
+        let x_k =
+            U256::from_be_hex("56f0b3e7e53fed658bf39e9f8691055ba8a58a935482b00b21b5678c16261fec");
+        assert_eq!(sum.get_x(), x_k);
+    }
+
+    #[test]
+    fn montgomery_point_big_add() {
+        let point_p = CURVE_25519.point(BIG_VALID_VALUE);
+        let point_q = point_p.double();
+
+        let sum = point_p.add(&point_q, &point_p);
+
+        let x_k =
+            U256::from_be_hex("1431b13075c5e0c75f2aad0f604f7286fbeac04cdc6099a63c54cdd96cbd3333");
+        assert_eq!(sum.get_x(), x_k);
+    }
+
+    #[test]
+    fn montgomery_ladder_double() {
+        let p = MontgomeryPoint {
+            x: U256::from_u8(2),
+            z: U256::from_u8(3),
+            curve: CURVE_25519,
+        };
+
+        let n = U256::from_u8(2);
+        let a = p * n;
+        let b = p.double();
+
+        assert_eq!(a.get_x(), b.get_x());
+    }
+
+    #[test]
+    fn montgomery_ladder() {
+        let p = CURVE_25519.point(U256::from_u8(2));
+
+        let n = U256::from_u8(4);
+
+        let a = p.ladder(&n);
+
+        let x_k =
+            U256::from_be_hex("2eceef1936e6df00c49e7aedac94446cc3b156165b50f247a15fdcee5e065582");
+        assert_eq!(a.get_x(), x_k);
+    }
+
+    #[test]
+    fn large_montgomery_ladder() {
+        let p = CURVE_25519.point(U256::from_u8(9));
+
+        let n = BIG_VALID_VALUE;
+
+        let a = p.ladder(&n);
+
+        let x_k =
+            U256::from_be_hex("0b379a815f36aca005dd1b19b1d483cd73fad06225c6ad927fde9316214a68c1");
+        assert_eq!(a.get_x(), x_k);
+    }
+
+    #[test]
+    fn constant_time_montgomery_ladder() {
+        let p = CURVE_25519.point(U256::from_u8(9));
+
+        let n = BIG_VALID_VALUE;
+
+        let a = p.ct_ladder(&n);
+        let b = p.ladder(&n);
+
+        assert_eq!(b.get_x(), a.get_x());
+    }
+
+    #[test]
+    fn montgomery_ladder_commutative() {
+        let p = CURVE_25519.point(U256::from_u8(2));
+
+        let n_1 = U256::from_u8(3);
+        let n_2 = U256::from_u8(4);
+
+        let a = p * n_1 * n_2;
+        let b = p * n_2 * n_1;
+
+        assert_eq!(a.get_x(), b.get_x());
+    }
+
+    #[test]
+    fn montgomery_point_conversion() {
+        let p = CURVE_25519.point(U256::from_u8(9));
+
+        let n_1 = U256::from_u8(3);
+        let n_2 = U256::from_u8(4);
+
+        let a = p * n_1;
+
+        let x_a = a.get_x();
+        let pa = CURVE_25519.point(x_a) * n_2;
+
+        assert_eq!(pa.get_x(), (a * n_2).get_x());
+    }
+}

src/traits.rs

@@ -0,0 +1,22 @@
+use crypto_bigint::{AddMod, InvMod, MulMod, SubMod, U256};
+
+pub trait Curve<T, N>
+where
+    T: Point<N>,
+    N: AddMod<Output = N> + SubMod<Output = N> + MulMod<Output = N> + InvMod<Output = N> + Copy,
+{
+    fn point(self, x: N) -> T;
+    fn generator(self) -> T;
+}
+
+pub trait Point<N>
+where
+    N: AddMod + SubMod + MulMod + InvMod + Copy,
+{
+    fn get_x(self) -> N;
+
+    fn add(&self, rhs: &Self, neg: &Self) -> Self;
+    fn double(&self) -> Self;
+    fn mul(&self, x: &N) -> Self;
+    fn neg(&self) -> Self;
+}