• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 // Copyright (c) 2012 The Chromium Authors. All rights reserved.
2 // Use of this source code is governed by a BSD-style license that can be
3 // found in the LICENSE file.
4 
5 // This is an implementation of the P224 elliptic curve group. It's written to
6 // be short and simple rather than fast, although it's still constant-time.
7 //
8 // See http://www.imperialviolet.org/2010/12/04/ecc.html ([1]) for background.
9 
10 #include "crypto/p224.h"
11 
12 #include <stddef.h>
13 #include <stdint.h>
14 #include <string.h>
15 
16 #include "base/sys_byteorder.h"
17 
18 namespace {
19 
20 using base::HostToNet32;
21 using base::NetToHost32;
22 
23 // Field element functions.
24 //
25 // The field that we're dealing with is ℤ/pℤ where p = 2**224 - 2**96 + 1.
26 //
27 // Field elements are represented by a FieldElement, which is a typedef to an
28 // array of 8 uint32_t's. The value of a FieldElement, a, is:
29 //   a[0] + 2**28·a[1] + 2**56·a[1] + ... + 2**196·a[7]
30 //
31 // Using 28-bit limbs means that there's only 4 bits of headroom, which is less
32 // than we would really like. But it has the useful feature that we hit 2**224
33 // exactly, making the reflections during a reduce much nicer.
34 
35 using crypto::p224::FieldElement;
36 
37 // kP is the P224 prime.
38 const FieldElement kP = {
39   1, 0, 0, 268431360,
40   268435455, 268435455, 268435455, 268435455,
41 };
42 
43 void Contract(FieldElement* inout);
44 
45 // IsZero returns 0xffffffff if a == 0 mod p and 0 otherwise.
IsZero(const FieldElement & a)46 uint32_t IsZero(const FieldElement& a) {
47   FieldElement minimal;
48   memcpy(&minimal, &a, sizeof(minimal));
49   Contract(&minimal);
50 
51   uint32_t is_zero = 0, is_p = 0;
52   for (unsigned i = 0; i < 8; i++) {
53     is_zero |= minimal[i];
54     is_p |= minimal[i] - kP[i];
55   }
56 
57   // If either is_zero or is_p is 0, then we should return 1.
58   is_zero |= is_zero >> 16;
59   is_zero |= is_zero >> 8;
60   is_zero |= is_zero >> 4;
61   is_zero |= is_zero >> 2;
62   is_zero |= is_zero >> 1;
63 
64   is_p |= is_p >> 16;
65   is_p |= is_p >> 8;
66   is_p |= is_p >> 4;
67   is_p |= is_p >> 2;
68   is_p |= is_p >> 1;
69 
70   // For is_zero and is_p, the LSB is 0 iff all the bits are zero.
71   is_zero &= is_p & 1;
72   is_zero = (~is_zero) << 31;
73   is_zero = static_cast<int32_t>(is_zero) >> 31;
74   return is_zero;
75 }
76 
77 // Add computes *out = a+b
78 //
79 // a[i] + b[i] < 2**32
Add(FieldElement * out,const FieldElement & a,const FieldElement & b)80 void Add(FieldElement* out, const FieldElement& a, const FieldElement& b) {
81   for (int i = 0; i < 8; i++) {
82     (*out)[i] = a[i] + b[i];
83   }
84 }
85 
86 static const uint32_t kTwo31p3 = (1u << 31) + (1u << 3);
87 static const uint32_t kTwo31m3 = (1u << 31) - (1u << 3);
88 static const uint32_t kTwo31m15m3 = (1u << 31) - (1u << 15) - (1u << 3);
89 // kZero31ModP is 0 mod p where bit 31 is set in all limbs so that we can
90 // subtract smaller amounts without underflow. See the section "Subtraction" in
91 // [1] for why.
92 static const FieldElement kZero31ModP = {
93   kTwo31p3, kTwo31m3, kTwo31m3, kTwo31m15m3,
94   kTwo31m3, kTwo31m3, kTwo31m3, kTwo31m3
95 };
96 
97 // Subtract computes *out = a-b
98 //
99 // a[i], b[i] < 2**30
100 // out[i] < 2**32
Subtract(FieldElement * out,const FieldElement & a,const FieldElement & b)101 void Subtract(FieldElement* out, const FieldElement& a, const FieldElement& b) {
102   for (int i = 0; i < 8; i++) {
103     // See the section on "Subtraction" in [1] for details.
104     (*out)[i] = a[i] + kZero31ModP[i] - b[i];
105   }
106 }
107 
108 static const uint64_t kTwo63p35 = (1ull << 63) + (1ull << 35);
109 static const uint64_t kTwo63m35 = (1ull << 63) - (1ull << 35);
110 static const uint64_t kTwo63m35m19 = (1ull << 63) - (1ull << 35) - (1ull << 19);
111 // kZero63ModP is 0 mod p where bit 63 is set in all limbs. See the section
112 // "Subtraction" in [1] for why.
113 static const uint64_t kZero63ModP[8] = {
114     kTwo63p35,    kTwo63m35, kTwo63m35, kTwo63m35,
115     kTwo63m35m19, kTwo63m35, kTwo63m35, kTwo63m35,
116 };
117 
118 static const uint32_t kBottom28Bits = 0xfffffff;
119 
120 // LargeFieldElement also represents an element of the field. The limbs are
121 // still spaced 28-bits apart and in little-endian order. So the limbs are at
122 // 0, 28, 56, ..., 392 bits, each 64-bits wide.
123 typedef uint64_t LargeFieldElement[15];
124 
125 // ReduceLarge converts a LargeFieldElement to a FieldElement.
126 //
127 // in[i] < 2**62
ReduceLarge(FieldElement * out,LargeFieldElement * inptr)128 void ReduceLarge(FieldElement* out, LargeFieldElement* inptr) {
129   LargeFieldElement& in(*inptr);
130 
131   for (int i = 0; i < 8; i++) {
132     in[i] += kZero63ModP[i];
133   }
134 
135   // Eliminate the coefficients at 2**224 and greater while maintaining the
136   // same value mod p.
137   for (int i = 14; i >= 8; i--) {
138     in[i-8] -= in[i];  // reflection off the "+1" term of p.
139     in[i-5] += (in[i] & 0xffff) << 12;  // part of the "-2**96" reflection.
140     in[i-4] += in[i] >> 16;  // the rest of the "-2**96" reflection.
141   }
142   in[8] = 0;
143   // in[0..8] < 2**64
144 
145   // As the values become small enough, we start to store them in |out| and use
146   // 32-bit operations.
147   for (int i = 1; i < 8; i++) {
148     in[i+1] += in[i] >> 28;
149     (*out)[i] = static_cast<uint32_t>(in[i] & kBottom28Bits);
150   }
151   // Eliminate the term at 2*224 that we introduced while keeping the same
152   // value mod p.
153   in[0] -= in[8];  // reflection off the "+1" term of p.
154   (*out)[3] += static_cast<uint32_t>(in[8] & 0xffff) << 12;  // "-2**96" term
155   (*out)[4] += static_cast<uint32_t>(in[8] >> 16);  // rest of "-2**96" term
156   // in[0] < 2**64
157   // out[3] < 2**29
158   // out[4] < 2**29
159   // out[1,2,5..7] < 2**28
160 
161   (*out)[0] = static_cast<uint32_t>(in[0] & kBottom28Bits);
162   (*out)[1] += static_cast<uint32_t>((in[0] >> 28) & kBottom28Bits);
163   (*out)[2] += static_cast<uint32_t>(in[0] >> 56);
164   // out[0] < 2**28
165   // out[1..4] < 2**29
166   // out[5..7] < 2**28
167 }
168 
169 // Mul computes *out = a*b
170 //
171 // a[i] < 2**29, b[i] < 2**30 (or vice versa)
172 // out[i] < 2**29
Mul(FieldElement * out,const FieldElement & a,const FieldElement & b)173 void Mul(FieldElement* out, const FieldElement& a, const FieldElement& b) {
174   LargeFieldElement tmp;
175   memset(&tmp, 0, sizeof(tmp));
176 
177   for (int i = 0; i < 8; i++) {
178     for (int j = 0; j < 8; j++) {
179       tmp[i + j] += static_cast<uint64_t>(a[i]) * static_cast<uint64_t>(b[j]);
180     }
181   }
182 
183   ReduceLarge(out, &tmp);
184 }
185 
186 // Square computes *out = a*a
187 //
188 // a[i] < 2**29
189 // out[i] < 2**29
Square(FieldElement * out,const FieldElement & a)190 void Square(FieldElement* out, const FieldElement& a) {
191   LargeFieldElement tmp;
192   memset(&tmp, 0, sizeof(tmp));
193 
194   for (int i = 0; i < 8; i++) {
195     for (int j = 0; j <= i; j++) {
196       uint64_t r = static_cast<uint64_t>(a[i]) * static_cast<uint64_t>(a[j]);
197       if (i == j) {
198         tmp[i+j] += r;
199       } else {
200         tmp[i+j] += r << 1;
201       }
202     }
203   }
204 
205   ReduceLarge(out, &tmp);
206 }
207 
208 // Reduce reduces the coefficients of in_out to smaller bounds.
209 //
210 // On entry: a[i] < 2**31 + 2**30
211 // On exit: a[i] < 2**29
Reduce(FieldElement * in_out)212 void Reduce(FieldElement* in_out) {
213   FieldElement& a = *in_out;
214 
215   for (int i = 0; i < 7; i++) {
216     a[i+1] += a[i] >> 28;
217     a[i] &= kBottom28Bits;
218   }
219   uint32_t top = a[7] >> 28;
220   a[7] &= kBottom28Bits;
221 
222   // top < 2**4
223   // Constant-time: mask = (top != 0) ? 0xffffffff : 0
224   uint32_t mask = top;
225   mask |= mask >> 2;
226   mask |= mask >> 1;
227   mask <<= 31;
228   mask = static_cast<uint32_t>(static_cast<int32_t>(mask) >> 31);
229 
230   // Eliminate top while maintaining the same value mod p.
231   a[0] -= top;
232   a[3] += top << 12;
233 
234   // We may have just made a[0] negative but, if we did, then we must
235   // have added something to a[3], thus it's > 2**12. Therefore we can
236   // carry down to a[0].
237   a[3] -= 1 & mask;
238   a[2] += mask & ((1<<28) - 1);
239   a[1] += mask & ((1<<28) - 1);
240   a[0] += mask & (1<<28);
241 }
242 
243 // Invert calcuates *out = in**-1 by computing in**(2**224 - 2**96 - 1), i.e.
244 // Fermat's little theorem.
Invert(FieldElement * out,const FieldElement & in)245 void Invert(FieldElement* out, const FieldElement& in) {
246   FieldElement f1, f2, f3, f4;
247 
248   Square(&f1, in);                        // 2
249   Mul(&f1, f1, in);                       // 2**2 - 1
250   Square(&f1, f1);                        // 2**3 - 2
251   Mul(&f1, f1, in);                       // 2**3 - 1
252   Square(&f2, f1);                        // 2**4 - 2
253   Square(&f2, f2);                        // 2**5 - 4
254   Square(&f2, f2);                        // 2**6 - 8
255   Mul(&f1, f1, f2);                       // 2**6 - 1
256   Square(&f2, f1);                        // 2**7 - 2
257   for (int i = 0; i < 5; i++) {           // 2**12 - 2**6
258     Square(&f2, f2);
259   }
260   Mul(&f2, f2, f1);                       // 2**12 - 1
261   Square(&f3, f2);                        // 2**13 - 2
262   for (int i = 0; i < 11; i++) {          // 2**24 - 2**12
263     Square(&f3, f3);
264   }
265   Mul(&f2, f3, f2);                       // 2**24 - 1
266   Square(&f3, f2);                        // 2**25 - 2
267   for (int i = 0; i < 23; i++) {          // 2**48 - 2**24
268     Square(&f3, f3);
269   }
270   Mul(&f3, f3, f2);                       // 2**48 - 1
271   Square(&f4, f3);                        // 2**49 - 2
272   for (int i = 0; i < 47; i++) {          // 2**96 - 2**48
273     Square(&f4, f4);
274   }
275   Mul(&f3, f3, f4);                       // 2**96 - 1
276   Square(&f4, f3);                        // 2**97 - 2
277   for (int i = 0; i < 23; i++) {          // 2**120 - 2**24
278     Square(&f4, f4);
279   }
280   Mul(&f2, f4, f2);                       // 2**120 - 1
281   for (int i = 0; i < 6; i++) {           // 2**126 - 2**6
282     Square(&f2, f2);
283   }
284   Mul(&f1, f1, f2);                       // 2**126 - 1
285   Square(&f1, f1);                        // 2**127 - 2
286   Mul(&f1, f1, in);                       // 2**127 - 1
287   for (int i = 0; i < 97; i++) {          // 2**224 - 2**97
288     Square(&f1, f1);
289   }
290   Mul(out, f1, f3);                       // 2**224 - 2**96 - 1
291 }
292 
293 // Contract converts a FieldElement to its minimal, distinguished form.
294 //
295 // On entry, in[i] < 2**29
296 // On exit, in[i] < 2**28
Contract(FieldElement * inout)297 void Contract(FieldElement* inout) {
298   FieldElement& out = *inout;
299 
300   // Reduce the coefficients to < 2**28.
301   for (int i = 0; i < 7; i++) {
302     out[i+1] += out[i] >> 28;
303     out[i] &= kBottom28Bits;
304   }
305   uint32_t top = out[7] >> 28;
306   out[7] &= kBottom28Bits;
307 
308   // Eliminate top while maintaining the same value mod p.
309   out[0] -= top;
310   out[3] += top << 12;
311 
312   // We may just have made out[0] negative. So we carry down. If we made
313   // out[0] negative then we know that out[3] is sufficiently positive
314   // because we just added to it.
315   for (int i = 0; i < 3; i++) {
316     uint32_t mask = static_cast<uint32_t>(static_cast<int32_t>(out[i]) >> 31);
317     out[i] += (1 << 28) & mask;
318     out[i+1] -= 1 & mask;
319   }
320 
321   // We might have pushed out[3] over 2**28 so we perform another, partial
322   // carry chain.
323   for (int i = 3; i < 7; i++) {
324     out[i+1] += out[i] >> 28;
325     out[i] &= kBottom28Bits;
326   }
327   top = out[7] >> 28;
328   out[7] &= kBottom28Bits;
329 
330   // Eliminate top while maintaining the same value mod p.
331   out[0] -= top;
332   out[3] += top << 12;
333 
334   // There are two cases to consider for out[3]:
335   //   1) The first time that we eliminated top, we didn't push out[3] over
336   //      2**28. In this case, the partial carry chain didn't change any values
337   //      and top is zero.
338   //   2) We did push out[3] over 2**28 the first time that we eliminated top.
339   //      The first value of top was in [0..16), therefore, prior to eliminating
340   //      the first top, 0xfff1000 <= out[3] <= 0xfffffff. Therefore, after
341   //      overflowing and being reduced by the second carry chain, out[3] <=
342   //      0xf000. Thus it cannot have overflowed when we eliminated top for the
343   //      second time.
344 
345   // Again, we may just have made out[0] negative, so do the same carry down.
346   // As before, if we made out[0] negative then we know that out[3] is
347   // sufficiently positive.
348   for (int i = 0; i < 3; i++) {
349     uint32_t mask = static_cast<uint32_t>(static_cast<int32_t>(out[i]) >> 31);
350     out[i] += (1 << 28) & mask;
351     out[i+1] -= 1 & mask;
352   }
353 
354   // The value is < 2**224, but maybe greater than p. In order to reduce to a
355   // unique, minimal value we see if the value is >= p and, if so, subtract p.
356 
357   // First we build a mask from the top four limbs, which must all be
358   // equal to bottom28Bits if the whole value is >= p. If top_4_all_ones
359   // ends up with any zero bits in the bottom 28 bits, then this wasn't
360   // true.
361   uint32_t top_4_all_ones = 0xffffffffu;
362   for (int i = 4; i < 8; i++) {
363     top_4_all_ones &= out[i];
364   }
365   top_4_all_ones |= 0xf0000000;
366   // Now we replicate any zero bits to all the bits in top_4_all_ones.
367   top_4_all_ones &= top_4_all_ones >> 16;
368   top_4_all_ones &= top_4_all_ones >> 8;
369   top_4_all_ones &= top_4_all_ones >> 4;
370   top_4_all_ones &= top_4_all_ones >> 2;
371   top_4_all_ones &= top_4_all_ones >> 1;
372   top_4_all_ones =
373       static_cast<uint32_t>(static_cast<int32_t>(top_4_all_ones << 31) >> 31);
374 
375   // Now we test whether the bottom three limbs are non-zero.
376   uint32_t bottom_3_non_zero = out[0] | out[1] | out[2];
377   bottom_3_non_zero |= bottom_3_non_zero >> 16;
378   bottom_3_non_zero |= bottom_3_non_zero >> 8;
379   bottom_3_non_zero |= bottom_3_non_zero >> 4;
380   bottom_3_non_zero |= bottom_3_non_zero >> 2;
381   bottom_3_non_zero |= bottom_3_non_zero >> 1;
382   bottom_3_non_zero =
383       static_cast<uint32_t>(static_cast<int32_t>(bottom_3_non_zero) >> 31);
384 
385   // Everything depends on the value of out[3].
386   //    If it's > 0xffff000 and top_4_all_ones != 0 then the whole value is >= p
387   //    If it's = 0xffff000 and top_4_all_ones != 0 and bottom_3_non_zero != 0,
388   //      then the whole value is >= p
389   //    If it's < 0xffff000, then the whole value is < p
390   uint32_t n = out[3] - 0xffff000;
391   uint32_t out_3_equal = n;
392   out_3_equal |= out_3_equal >> 16;
393   out_3_equal |= out_3_equal >> 8;
394   out_3_equal |= out_3_equal >> 4;
395   out_3_equal |= out_3_equal >> 2;
396   out_3_equal |= out_3_equal >> 1;
397   out_3_equal =
398       ~static_cast<uint32_t>(static_cast<int32_t>(out_3_equal << 31) >> 31);
399 
400   // If out[3] > 0xffff000 then n's MSB will be zero.
401   uint32_t out_3_gt =
402       ~static_cast<uint32_t>(static_cast<int32_t>(n << 31) >> 31);
403 
404   uint32_t mask =
405       top_4_all_ones & ((out_3_equal & bottom_3_non_zero) | out_3_gt);
406   out[0] -= 1 & mask;
407   out[3] -= 0xffff000 & mask;
408   out[4] -= 0xfffffff & mask;
409   out[5] -= 0xfffffff & mask;
410   out[6] -= 0xfffffff & mask;
411   out[7] -= 0xfffffff & mask;
412 }
413 
414 
415 // Group element functions.
416 //
417 // These functions deal with group elements. The group is an elliptic curve
418 // group with a = -3 defined in FIPS 186-3, section D.2.2.
419 
420 using crypto::p224::Point;
421 
422 // kB is parameter of the elliptic curve.
423 const FieldElement kB = {
424   55967668, 11768882, 265861671, 185302395,
425   39211076, 180311059, 84673715, 188764328,
426 };
427 
428 void CopyConditional(Point* out, const Point& a, uint32_t mask);
429 void DoubleJacobian(Point* out, const Point& a);
430 
431 // AddJacobian computes *out = a+b where a != b.
AddJacobian(Point * out,const Point & a,const Point & b)432 void AddJacobian(Point *out,
433                  const Point& a,
434                  const Point& b) {
435   // See http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#addition-add-2007-bl
436   FieldElement z1z1, z2z2, u1, u2, s1, s2, h, i, j, r, v;
437 
438   uint32_t z1_is_zero = IsZero(a.z);
439   uint32_t z2_is_zero = IsZero(b.z);
440 
441   // Z1Z1 = Z1²
442   Square(&z1z1, a.z);
443 
444   // Z2Z2 = Z2²
445   Square(&z2z2, b.z);
446 
447   // U1 = X1*Z2Z2
448   Mul(&u1, a.x, z2z2);
449 
450   // U2 = X2*Z1Z1
451   Mul(&u2, b.x, z1z1);
452 
453   // S1 = Y1*Z2*Z2Z2
454   Mul(&s1, b.z, z2z2);
455   Mul(&s1, a.y, s1);
456 
457   // S2 = Y2*Z1*Z1Z1
458   Mul(&s2, a.z, z1z1);
459   Mul(&s2, b.y, s2);
460 
461   // H = U2-U1
462   Subtract(&h, u2, u1);
463   Reduce(&h);
464   uint32_t x_equal = IsZero(h);
465 
466   // I = (2*H)²
467   for (int k = 0; k < 8; k++) {
468     i[k] = h[k] << 1;
469   }
470   Reduce(&i);
471   Square(&i, i);
472 
473   // J = H*I
474   Mul(&j, h, i);
475   // r = 2*(S2-S1)
476   Subtract(&r, s2, s1);
477   Reduce(&r);
478   uint32_t y_equal = IsZero(r);
479 
480   if (x_equal && y_equal && !z1_is_zero && !z2_is_zero) {
481     // The two input points are the same therefore we must use the dedicated
482     // doubling function as the slope of the line is undefined.
483     DoubleJacobian(out, a);
484     return;
485   }
486 
487   for (int k = 0; k < 8; k++) {
488     r[k] <<= 1;
489   }
490   Reduce(&r);
491 
492   // V = U1*I
493   Mul(&v, u1, i);
494 
495   // Z3 = ((Z1+Z2)²-Z1Z1-Z2Z2)*H
496   Add(&z1z1, z1z1, z2z2);
497   Add(&z2z2, a.z, b.z);
498   Reduce(&z2z2);
499   Square(&z2z2, z2z2);
500   Subtract(&out->z, z2z2, z1z1);
501   Reduce(&out->z);
502   Mul(&out->z, out->z, h);
503 
504   // X3 = r²-J-2*V
505   for (int k = 0; k < 8; k++) {
506     z1z1[k] = v[k] << 1;
507   }
508   Add(&z1z1, j, z1z1);
509   Reduce(&z1z1);
510   Square(&out->x, r);
511   Subtract(&out->x, out->x, z1z1);
512   Reduce(&out->x);
513 
514   // Y3 = r*(V-X3)-2*S1*J
515   for (int k = 0; k < 8; k++) {
516     s1[k] <<= 1;
517   }
518   Mul(&s1, s1, j);
519   Subtract(&z1z1, v, out->x);
520   Reduce(&z1z1);
521   Mul(&z1z1, z1z1, r);
522   Subtract(&out->y, z1z1, s1);
523   Reduce(&out->y);
524 
525   CopyConditional(out, a, z2_is_zero);
526   CopyConditional(out, b, z1_is_zero);
527 }
528 
529 // DoubleJacobian computes *out = a+a.
DoubleJacobian(Point * out,const Point & a)530 void DoubleJacobian(Point* out, const Point& a) {
531   // See http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#doubling-dbl-2001-b
532   FieldElement delta, gamma, beta, alpha, t;
533 
534   Square(&delta, a.z);
535   Square(&gamma, a.y);
536   Mul(&beta, a.x, gamma);
537 
538   // alpha = 3*(X1-delta)*(X1+delta)
539   Add(&t, a.x, delta);
540   for (int i = 0; i < 8; i++) {
541           t[i] += t[i] << 1;
542   }
543   Reduce(&t);
544   Subtract(&alpha, a.x, delta);
545   Reduce(&alpha);
546   Mul(&alpha, alpha, t);
547 
548   // Z3 = (Y1+Z1)²-gamma-delta
549   Add(&out->z, a.y, a.z);
550   Reduce(&out->z);
551   Square(&out->z, out->z);
552   Subtract(&out->z, out->z, gamma);
553   Reduce(&out->z);
554   Subtract(&out->z, out->z, delta);
555   Reduce(&out->z);
556 
557   // X3 = alpha²-8*beta
558   for (int i = 0; i < 8; i++) {
559           delta[i] = beta[i] << 3;
560   }
561   Reduce(&delta);
562   Square(&out->x, alpha);
563   Subtract(&out->x, out->x, delta);
564   Reduce(&out->x);
565 
566   // Y3 = alpha*(4*beta-X3)-8*gamma²
567   for (int i = 0; i < 8; i++) {
568           beta[i] <<= 2;
569   }
570   Reduce(&beta);
571   Subtract(&beta, beta, out->x);
572   Reduce(&beta);
573   Square(&gamma, gamma);
574   for (int i = 0; i < 8; i++) {
575           gamma[i] <<= 3;
576   }
577   Reduce(&gamma);
578   Mul(&out->y, alpha, beta);
579   Subtract(&out->y, out->y, gamma);
580   Reduce(&out->y);
581 }
582 
583 // CopyConditional sets *out=a if mask is 0xffffffff. mask must be either 0 of
584 // 0xffffffff.
CopyConditional(Point * out,const Point & a,uint32_t mask)585 void CopyConditional(Point* out, const Point& a, uint32_t mask) {
586   for (int i = 0; i < 8; i++) {
587     out->x[i] ^= mask & (a.x[i] ^ out->x[i]);
588     out->y[i] ^= mask & (a.y[i] ^ out->y[i]);
589     out->z[i] ^= mask & (a.z[i] ^ out->z[i]);
590   }
591 }
592 
593 // ScalarMult calculates *out = a*scalar where scalar is a big-endian number of
594 // length scalar_len and != 0.
ScalarMult(Point * out,const Point & a,const uint8_t * scalar,size_t scalar_len)595 void ScalarMult(Point* out,
596                 const Point& a,
597                 const uint8_t* scalar,
598                 size_t scalar_len) {
599   memset(out, 0, sizeof(*out));
600   Point tmp;
601 
602   for (size_t i = 0; i < scalar_len; i++) {
603     for (unsigned int bit_num = 0; bit_num < 8; bit_num++) {
604       DoubleJacobian(out, *out);
605       uint32_t bit = static_cast<uint32_t>(static_cast<int32_t>(
606           (((scalar[i] >> (7 - bit_num)) & 1) << 31) >> 31));
607       AddJacobian(&tmp, a, *out);
608       CopyConditional(out, tmp, bit);
609     }
610   }
611 }
612 
613 // Get224Bits reads 7 words from in and scatters their contents in
614 // little-endian form into 8 words at out, 28 bits per output word.
Get224Bits(uint32_t * out,const uint32_t * in)615 void Get224Bits(uint32_t* out, const uint32_t* in) {
616   out[0] = NetToHost32(in[6]) & kBottom28Bits;
617   out[1] = ((NetToHost32(in[5]) << 4) |
618             (NetToHost32(in[6]) >> 28)) & kBottom28Bits;
619   out[2] = ((NetToHost32(in[4]) << 8) |
620             (NetToHost32(in[5]) >> 24)) & kBottom28Bits;
621   out[3] = ((NetToHost32(in[3]) << 12) |
622             (NetToHost32(in[4]) >> 20)) & kBottom28Bits;
623   out[4] = ((NetToHost32(in[2]) << 16) |
624             (NetToHost32(in[3]) >> 16)) & kBottom28Bits;
625   out[5] = ((NetToHost32(in[1]) << 20) |
626             (NetToHost32(in[2]) >> 12)) & kBottom28Bits;
627   out[6] = ((NetToHost32(in[0]) << 24) |
628             (NetToHost32(in[1]) >> 8)) & kBottom28Bits;
629   out[7] = (NetToHost32(in[0]) >> 4) & kBottom28Bits;
630 }
631 
632 // Put224Bits performs the inverse operation to Get224Bits: taking 28 bits from
633 // each of 8 input words and writing them in big-endian order to 7 words at
634 // out.
Put224Bits(uint32_t * out,const uint32_t * in)635 void Put224Bits(uint32_t* out, const uint32_t* in) {
636   out[6] = HostToNet32((in[0] >> 0) | (in[1] << 28));
637   out[5] = HostToNet32((in[1] >> 4) | (in[2] << 24));
638   out[4] = HostToNet32((in[2] >> 8) | (in[3] << 20));
639   out[3] = HostToNet32((in[3] >> 12) | (in[4] << 16));
640   out[2] = HostToNet32((in[4] >> 16) | (in[5] << 12));
641   out[1] = HostToNet32((in[5] >> 20) | (in[6] << 8));
642   out[0] = HostToNet32((in[6] >> 24) | (in[7] << 4));
643 }
644 
645 }  // anonymous namespace
646 
647 namespace crypto {
648 
649 namespace p224 {
650 
SetFromString(const base::StringPiece & in)651 bool Point::SetFromString(const base::StringPiece& in) {
652   if (in.size() != 2*28)
653     return false;
654   const uint32_t* inwords = reinterpret_cast<const uint32_t*>(in.data());
655   Get224Bits(x, inwords);
656   Get224Bits(y, inwords + 7);
657   memset(&z, 0, sizeof(z));
658   z[0] = 1;
659 
660   // Check that the point is on the curve, i.e. that y² = x³ - 3x + b.
661   FieldElement lhs;
662   Square(&lhs, y);
663   Contract(&lhs);
664 
665   FieldElement rhs;
666   Square(&rhs, x);
667   Mul(&rhs, x, rhs);
668 
669   FieldElement three_x;
670   for (int i = 0; i < 8; i++) {
671     three_x[i] = x[i] * 3;
672   }
673   Reduce(&three_x);
674   Subtract(&rhs, rhs, three_x);
675   Reduce(&rhs);
676 
677   ::Add(&rhs, rhs, kB);
678   Contract(&rhs);
679   return memcmp(&lhs, &rhs, sizeof(lhs)) == 0;
680 }
681 
ToString() const682 std::string Point::ToString() const {
683   FieldElement zinv, zinv_sq, xx, yy;
684 
685   // If this is the point at infinity we return a string of all zeros.
686   if (IsZero(this->z)) {
687     static const char zeros[56] = {0};
688     return std::string(zeros, sizeof(zeros));
689   }
690 
691   Invert(&zinv, this->z);
692   Square(&zinv_sq, zinv);
693   Mul(&xx, x, zinv_sq);
694   Mul(&zinv_sq, zinv_sq, zinv);
695   Mul(&yy, y, zinv_sq);
696 
697   Contract(&xx);
698   Contract(&yy);
699 
700   uint32_t outwords[14];
701   Put224Bits(outwords, xx);
702   Put224Bits(outwords + 7, yy);
703   return std::string(reinterpret_cast<const char*>(outwords), sizeof(outwords));
704 }
705 
ScalarMult(const Point & in,const uint8_t * scalar,Point * out)706 void ScalarMult(const Point& in, const uint8_t* scalar, Point* out) {
707   ::ScalarMult(out, in, scalar, 28);
708 }
709 
710 // kBasePoint is the base point (generator) of the elliptic curve group.
711 static const Point kBasePoint = {
712   {22813985, 52956513, 34677300, 203240812,
713    12143107, 133374265, 225162431, 191946955},
714   {83918388, 223877528, 122119236, 123340192,
715    266784067, 263504429, 146143011, 198407736},
716   {1, 0, 0, 0, 0, 0, 0, 0},
717 };
718 
ScalarBaseMult(const uint8_t * scalar,Point * out)719 void ScalarBaseMult(const uint8_t* scalar, Point* out) {
720   ::ScalarMult(out, kBasePoint, scalar, 28);
721 }
722 
Add(const Point & a,const Point & b,Point * out)723 void Add(const Point& a, const Point& b, Point* out) {
724   AddJacobian(out, a, b);
725 }
726 
Negate(const Point & in,Point * out)727 void Negate(const Point& in, Point* out) {
728   // Guide to elliptic curve cryptography, page 89 suggests that (X : X+Y : Z)
729   // is the negative in Jacobian coordinates, but it doesn't actually appear to
730   // be true in testing so this performs the negation in affine coordinates.
731   FieldElement zinv, zinv_sq, y;
732   Invert(&zinv, in.z);
733   Square(&zinv_sq, zinv);
734   Mul(&out->x, in.x, zinv_sq);
735   Mul(&zinv_sq, zinv_sq, zinv);
736   Mul(&y, in.y, zinv_sq);
737 
738   Subtract(&out->y, kP, y);
739   Reduce(&out->y);
740 
741   memset(&out->z, 0, sizeof(out->z));
742   out->z[0] = 1;
743 }
744 
745 }  // namespace p224
746 
747 }  // namespace crypto
748