crypto: special case ∞+a, a+∞ and a+a in p224.

crypto/p224.cc

Index: crypto/p224.cc
diff --git a/crypto/p224.cc b/crypto/p224.cc
index 575b51f4ce946cf06b4443ddf68cda8667240752..b277e75ce4d24f82074439908eb999a23a0b2652 100644
--- a/crypto/p224.cc
+++ b/crypto/p224.cc
@@ -32,6 +32,46 @@ using base::NetToHost32;
 
 using crypto::p224::FieldElement;
 
+// kP is the P224 prime.
+const FieldElement kP = {
+  1, 0, 0, 268431360,
+  268435455, 268435455, 268435455, 268435455,
+};
+
+void Contract(FieldElement* inout);
+
+// IsZero returns 0xffffffff if a == 0 mod p and 0 otherwise.
+uint32_t IsZero(const FieldElement& a) {
+  FieldElement minimal;
+  memcpy(&minimal, &a, sizeof(minimal));
+  Contract(&minimal);
+
+  uint32 is_zero = 0, is_p = 0;
+  for (unsigned i = 0; i < 8; i++) {
+    is_zero |= minimal[i];
+    is_p |= minimal[i] - kP[i];
+  }
+
+  // If either is_zero or is_p is 0, then we should return 1.
+  is_zero |= is_zero >> 16;
+  is_zero |= is_zero >> 8;
+  is_zero |= is_zero >> 4;
+  is_zero |= is_zero >> 2;
+  is_zero |= is_zero >> 1;
+
+  is_p |= is_p >> 16;
+  is_p |= is_p >> 8;
+  is_p |= is_p >> 4;
+  is_p |= is_p >> 2;
+  is_p |= is_p >> 1;
+
+  // For is_zero and is_p, the LSB is 0 iff all the bits are zero.
+  is_zero &= is_p & 1;
+  is_zero = (~is_zero) << 31;
+  is_zero = static_cast<int32>(is_zero) >> 31;
+  return is_zero;
+}
+
 // Add computes *out = a+b
 //
 // a[i] + b[i] < 2**32
@@ -318,7 +358,7 @@ void Contract(FieldElement* inout) {
   // true.
   uint32 top4AllOnes = 0xffffffffu;

[willchan no longer on Chromium, 2012/07/25 23:35:26]

I know you're going to hate me but...
nit: Top4AllOnes is proper camel casing naming style in Chromium.

[agl, 2012/07/27 17:21:06]

No, that's my screw up. I debugged this code in Go and translated across. I've
changed to top_4_all_ones as its a local variable.

   for (int i = 4; i < 8; i++) {
-    top4AllOnes &= (out[i] & kBottom28Bits) - 1;
+    top4AllOnes &= out[i];
   }
   top4AllOnes |= 0xf0000000;
   // Now we replicate any zero bits to all the bits in top4AllOnes.
@@ -338,7 +378,7 @@ void Contract(FieldElement* inout) {
   bottom3NonZero |= bottom3NonZero >> 2;
   bottom3NonZero |= bottom3NonZero >> 1;
   bottom3NonZero =
-      static_cast<uint32>(static_cast<int32>(bottom3NonZero << 31) >> 31);
+      static_cast<uint32>(static_cast<int32>(bottom3NonZero) >> 31);
 
   // Everything depends on the value of out[3].
   //    If it's > 0xffff000 and top4AllOnes != 0 then the whole value is >= p
@@ -375,18 +415,15 @@ void Contract(FieldElement* inout) {
 
 using crypto::p224::Point;
 
-// kP is the P224 prime.
-const FieldElement kP = {
-  1, 0, 0, 268431360,
-  268435455, 268435455, 268435455, 268435455,
-};
-
 // kB is parameter of the elliptic curve.
 const FieldElement kB = {
   55967668, 11768882, 265861671, 185302395,
   39211076, 180311059, 84673715, 188764328,
 };
 
+void CopyConditional(Point* out, const Point& a, uint32 mask);

[willchan no longer on Chromium, 2012/07/25 23:35:26]

Google style is to put output parameters last. That said, you seem to have
flouted this Google style rule in the whole file, so it's best that you stay
consistently inconsistent :)

[agl, 2012/07/27 17:21:06]

Well, I'm just being consistent with all the memcpys and every other system
function :)

+void DoubleJacobian(Point* out, const Point& a);
+
 // AddJacobian computes *out = a+b where a != b.
 void AddJacobian(Point *out,
                  const Point& a,
@@ -394,6 +431,9 @@ void AddJacobian(Point *out,
   // See http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#addition-add-2007-bl
   FieldElement z1z1, z2z2, u1, u2, s1, s2, h, i, j, r, v;
 
+  uint32 z1_is_zero = IsZero(a.z);
+  uint32 z2_is_zero = IsZero(b.z);
+
   // Z1Z1 = Z1²
   Square(&z1z1, a.z);
 
@@ -417,6 +457,7 @@ void AddJacobian(Point *out,
   // H = U2-U1
   Subtract(&h, u2, u1);
   Reduce(&h);
+  uint32 x_equal = IsZero(h);
 
   // I = (2*H)²
   for (int j = 0; j < 8; j++) {
@@ -430,6 +471,15 @@ void AddJacobian(Point *out,
   // r = 2*(S2-S1)
   Subtract(&r, s2, s1);
   Reduce(&r);
+  uint32 y_equal = IsZero(r);
+
+  if (x_equal && y_equal && !z1_is_zero && !z2_is_zero) {
+    // The two input points are the same therefore we must use the dedicated
+    // doubling function as the slope of the line is undefined.
+    DoubleJacobian(out, a);
+    return;
+  }
+
   for (int i = 0; i < 8; i++) {
     r[i] <<= 1;
   }
@@ -467,6 +517,9 @@ void AddJacobian(Point *out,
   Mul(&z1z1, z1z1, r);
   Subtract(&out->y, z1z1, s1);
   Reduce(&out->y);
+
+  CopyConditional(out, a, z2_is_zero);
+  CopyConditional(out, b, z1_is_zero);
 }
 
 // DoubleJacobian computes *out = a+a.
@@ -542,16 +595,13 @@ void ScalarMult(Point* out, const Point& a,
   memset(out, 0, sizeof(*out));
   Point tmp;
 
-  uint32 first_bit = 0xffffffff;
   for (size_t i = 0; i < scalar_len; i++) {
     for (unsigned int bit_num = 0; bit_num < 8; bit_num++) {
       DoubleJacobian(out, *out);
       uint32 bit = static_cast<uint32>(static_cast<int32>(
           (((scalar[i] >> (7 - bit_num)) & 1) << 31) >> 31));
       AddJacobian(&tmp, a, *out);
-      CopyConditional(out, a, first_bit & bit);
-      CopyConditional(out, tmp, ~first_bit & bit);
-      first_bit = first_bit & ~bit;
+      CopyConditional(out, tmp, bit);
     }
   }
 }
@@ -628,6 +678,13 @@ bool Point::SetFromString(const base::StringPiece& in) {
 std::string Point::ToString() const {
   FieldElement zinv, zinv_sq, x, y;
 
+  // If this is the point at infinity we return a string of all zeros.
+  if (IsZero(this->z)) {
+    char zeros[56];

[willchan no longer on Chromium, 2012/07/25 23:35:26]

Why isn't this just a static const? It'd POD, so the linker should initialize it
to 0. That said, I guess this is only used by tests, so it's not perf critical.
So whatever.

[agl, 2012/07/27 17:21:06]

Good point. Done.

+    memset(zeros, 0, sizeof(zeros));
+    return std::string(zeros, sizeof(zeros));
+  }
+
   Invert(&zinv, this->z);
   Square(&zinv_sq, zinv);
   Mul(&x, this->x, zinv_sq);