Improve Edwards448 decompression checks

Author: daxpedda

Cargo.lock

@@ -39,6 +39,7 @@ serde = ["dep:serdect", "ed448?/serde_bytes"]
 [dev-dependencies]
 hex-literal = "1"
 hex = "0.4"
+proptest = "1"
 rand_core = { version = "0.9", features = ["os_rng"] }
 rand_chacha = "0.9"
 serde_bare = "0.5"

ed448-goldilocks/Cargo.toml

@@ -5,7 +5,7 @@ use core::ops::{Add, AddAssign, Mul, MulAssign, Neg, Sub, SubAssign};
 
 use crate::curve::scalar_mul::variable_base;
 use crate::curve::twedwards::extended::ExtendedPoint as TwistedExtendedPoint;
-use crate::field::FieldElement;
+use crate::field::{ConstMontyType, FieldElement};
 use crate::*;
 use elliptic_curve::{
     CurveGroup, Error,
@@ -240,7 +240,7 @@ impl CompressedEdwardsY {
     /// - if the input point has nonzero torsion component.
     pub fn decompress(&self) -> CtOption<EdwardsPoint> {
         self.decompress_unchecked()
-            .and_then(|pt| CtOption::new(pt, pt.is_on_curve() & pt.is_torsion_free()))
+            .and_then(|pt| CtOption::new(pt, pt.is_in_subgroup()))
     }
 
     /// View this `CompressedEdwardsY` as an array of bytes.
@@ -661,6 +661,59 @@ impl EdwardsPoint {
         AffinePoint { x, y }
     }
 
+    // See https://eprint.iacr.org/2022/1164.
+    fn is_in_subgroup(self) -> Choice {
+        const A: FieldElement = FieldElement(ConstMontyType::new(&U448::from_be_hex(
+            "fffffffffffffffffffffffffffffffffffffffffffffffffffffffefffffffffffffffffffffffffffffffffffffffffffffffffffeceaf",
+        )));
+        const A1: FieldElement = FieldElement(ConstMontyType::new(&U448::from_u64(156320)));
+        const MINUS_SQRT_B1: FieldElement = FieldElement(ConstMontyType::new(&U448::from_be_hex(
+            "749a7410536c225f1025ca374176557d7839611d691caad26d74a1fca5cfad15f196642c0a4484b67f321025577cc6b5a6f443c2eaa36327",
+        )));
+
+        let mut e = self.X * (self.Z - self.Y);
+        let ee = e.square();
+        let mut u = FieldElement::A_PLUS_TWO_OVER_FOUR * (self.Z + self.Y) * e * self.X;
+        let w = self.Z.double() * (self.Z - self.Y);
+
+        let u2 = u.double().double();
+        let w2 = w.double();
+
+        let mut w1 = u2.sqrt();
+        let mut ok = w1.square().ct_eq(&u2);
+        let u1 = (u2 - A1 * ee - w1 * w2).half();
+
+        // If `u1` happens not to be a square, then `sqrt(u1)` returns `sqrt(-u1)`
+        // in that case (since we are in a finite field GF(q) with q = 3 mod 4,
+        // if `u1` is not a square then `-u1` must be a square). In such a case, we
+        // should replace `(u1,w1)` with `((B1*e^4)/u1, -w1)`. To avoid the division,
+        // we instead switch to an isomorphic curve; namely:
+        //   u2 = B1*(e^4)*u1
+        //   w2 = -w1*u1
+        //   e2 = e*u1
+        // Then:
+        //   w = sqrt(u2) = sqrt(-B1)*(e^2)*sqrt(-u1)
+        //   u = (w^2 - A*e^2 - w*w1)/2
+        let mut w = u1.sqrt();
+        let u1_is_square = w.square().ct_eq(&u1);
+        w1.conditional_assign(&-(w1 * u1), !u1_is_square);
+        e.conditional_assign(&(e * u1), !u1_is_square);
+        w.conditional_assign(&(MINUS_SQRT_B1 * ee * w), !u1_is_square);
+        u = (w.square() - A * e.square() - w * w1).half();
+
+        ok &= u.is_square();
+
+        // If the source point was a low-order point, then the computations
+        // above are incorrect. We handle this case here; among the
+        // low-order points, only the neutral point is in the prime-order
+        // subgroup.
+        let is_low_order = self.X.is_zero() | self.Y.is_zero();
+        let is_neutral = self.Y.ct_eq(&self.Z);
+        ok ^= is_low_order & (ok ^ is_neutral);
+
+        ok
+    }
+
     /// Edwards_Isogeny is derived from the doubling formula
     /// XXX: There is a duplicate method in the twisted edwards module to compute the dual isogeny
     /// XXX: Not much point trying to make it generic I think. So what we can do is optimise each respective isogeny method for a=1 or a = -1 (currently, I just made it really slow and simple)
@@ -972,6 +1025,8 @@ mod tests {
     use super::*;
     use elliptic_curve::Field;
     use hex_literal::hex;
+    use proptest::prelude::any;
+    use proptest::proptest;
     use rand_core::TryRngCore;
 
     fn hex_to_field(hex: &'static str) -> FieldElement {
@@ -1267,4 +1322,15 @@ mod tests {
 
         assert_eq!(computed_commitment, expected_commitment);
     }
+
+    proptest! {
+        #[test]
+        fn is_in_subgroup(
+           bytes in any::<[u8; 57]>()
+        ) {
+            let scalar = EdwardsScalar::from_bytes_mod_order(&bytes.into());
+            let point = EdwardsPoint::mul_by_generator(&scalar);
+            assert_eq!(point.is_in_subgroup().unwrap_u8(), 1);
+        }
+    }
 }

ed448-goldilocks/src/edwards/extended.rs

@@ -9,7 +9,7 @@ use crate::{
 use elliptic_curve::{
     array::Array,
     bigint::{
-        Integer, NonZero, U448, U704,
+        Integer, NonZero, U448, U704, Zero,
         consts::{U56, U84, U88},
     },
     group::cofactor::CofactorGroup,
@@ -261,6 +261,10 @@ impl FieldElement {
         self.0.retrieve().is_odd()
     }
 
+    pub fn is_zero(&self) -> Choice {
+        self.0.is_zero()
+    }
+
     /// Inverts a field element
     /// Previous chain length: 462, new length 460
     pub fn invert(&self) -> Self {
@@ -368,6 +372,10 @@ impl FieldElement {
         (inv_sqrt_x * u, zero_u | is_res)
     }
 
+    pub(crate) fn half(&self) -> FieldElement {
+        Self(self.0.div_by_2())
+    }
+
     pub(crate) fn map_to_curve_elligator2(&self) -> AffinePoint {
         let mut t1 = self.square(); // 1.   t1 = u^2
         t1 *= Self::Z; // 2.   t1 = Z * t1              // Z * u^2