Add some special case handling for ec_twin_mult.

Previously we were computing incorrect answers for `ec_twin_mult j S k T`
in the following 4 special cases: `S == 0`, `T == 0`, `S+T == 0` and
`S-T == 0`. In these cases, we must be more careful about performing
normalization in the twin multiply procedure.

lib/PrimeEC.cry

@@ -67,6 +67,15 @@ ec_is_nonsingular b = (fromInteger 4) * a^^3 + (fromInteger 27) * b^^2 != 0
 ec_is_identity : {p} (prime p, p > 3) => ProjectivePoint p -> Bit
 ec_is_identity S = S.z == 0 /\ ~(S.x == 0 /\ S.y == 0)
 
+/**
+ * Test two projective points for equality, up to the equivalence relation
+ * on projective points.
+ */
+ec_equal : {p} (prime p, p > 3) => ProjectivePoint p -> ProjectivePoint p -> Bit
+ec_equal S T =
+  (S.z == 0 /\ T.z == 0) \/
+  (S.z != 0 /\ T.z != 0 /\ ec_affinify S == ec_affinify T)
+
 /**
  * Compute a projective representative for the given affine point.
  */

src/Cryptol/PrimeEC.hs

@@ -263,6 +263,11 @@ ec_sub p s t = ec_add p s u
   where u = t{ py = Integer.minusBigNat (primeMod p) (py t) }
 {-# INLINE ec_sub #-}
 
+
+ec_negate :: PrimeModulus -> ProjectivePoint -> ProjectivePoint
+ec_negate p s = s{ py = Integer.minusBigNat (primeMod p) (py s) }
+{-# INLINE ec_negate #-}
+
 -- | Compute the elliptic curve group addition operation
 --   for values known not to be the identity.
 --   In other words, if @S@ and @T@ are projective points on a curve,
@@ -381,10 +386,41 @@ ec_mult p d s
 --   and normalize all four values.  This can be done jointly
 --   in a more efficient way than computing the necessary
 --   field inverses separately.
+--   When given points S and T, the returned tuple contains
+--   normalized representations for (S, T, S+T, S-T).
+--
+--   Note there are some special cases that must be handled separately.
 normalizeForTwinMult ::
   PrimeModulus -> ProjectivePoint -> ProjectivePoint ->
   (ProjectivePoint, ProjectivePoint, ProjectivePoint, ProjectivePoint)
-normalizeForTwinMult p s t = (s',t',spt',smt')
+normalizeForTwinMult p s t
+     -- S == 0 && T == 0
+   | Integer.isZeroBigNat a && Integer.isZeroBigNat b =
+        (zro, zro, zro, zro)
+
+     -- S == 0 && T != 0
+   | Integer.isZeroBigNat a =
+        let tnorm = ec_normalize p t
+         in (zro, tnorm, tnorm, ec_negate p tnorm)
+
+     -- T == 0 && S != 0
+   | Integer.isZeroBigNat b =
+        let snorm = ec_normalize p s
+         in (snorm, zro, snorm, snorm)
+
+     -- S+T == 0, both != 0
+   | Integer.isZeroBigNat c =
+        let snorm = ec_normalize p s
+         in (snorm, ec_negate p snorm, zro, ec_double p snorm)
+
+     -- S-T == 0, both != 0
+   | Integer.isZeroBigNat d =
+        let snorm = ec_normalize p s
+         in (snorm, snorm, ec_double p snorm, zro)
+
+     -- S, T, S+T and S-T all != 0
+   | otherwise = (s',t',spt',smt')
+
   where
   spt = ec_add p s t
   smt = ec_sub p s t

tests/regression/twinmult.cry

@@ -0,0 +1,19 @@
+import PrimeEC
+
+zro : ProjectivePoint 7
+zro = { x = 1, y = 1, z = 0 }
+
+property twin_mult_zro j k =
+  ec_equal`{7} (ec_twin_mult j zro k zro) zro
+
+property twin_mult_zro1 j k S =
+  ec_equal`{7} (ec_twin_mult j zro k S) (ec_mult k S)
+
+property twin_mult_zro2 j k S =
+  ec_equal`{7} (ec_twin_mult j S k zro) (ec_mult j S)
+
+property twin_mult_same j k S =
+  ec_equal`{7} (ec_twin_mult j S k S) (ec_add (ec_mult j S) (ec_mult k S))
+
+property twin_mult_neg j k S =
+  ec_equal`{7} (ec_twin_mult j S k (ec_negate S)) (ec_sub (ec_mult j S) (ec_mult k S))

tests/regression/twinmult.icry

@@ -0,0 +1,3 @@
+:l twinmult.cry
+
+:exhaust