sexp: Improve argument checking of sexp parser.
[libgcrypt.git] / mpi / ec.c
1 /* ec.c -  Elliptic Curve functions
2  * Copyright (C) 2007 Free Software Foundation, Inc.
3  * Copyright (C) 2013 g10 Code GmbH
4  *
5  * This file is part of Libgcrypt.
6  *
7  * Libgcrypt is free software; you can redistribute it and/or modify
8  * it under the terms of the GNU Lesser General Public License as
9  * published by the Free Software Foundation; either version 2.1 of
10  * the License, or (at your option) any later version.
11  *
12  * Libgcrypt is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15  * GNU Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with this program; if not, see <http://www.gnu.org/licenses/>.
19  */
20
21 #include <config.h>
22 #include <stdio.h>
23 #include <stdlib.h>
24 #include <errno.h>
25
26 #include "mpi-internal.h"
27 #include "longlong.h"
28 #include "g10lib.h"
29 #include "context.h"
30 #include "ec-context.h"
31 #include "ec-internal.h"
32
33
34 #define point_init(a)  _gcry_mpi_point_init ((a))
35 #define point_free(a)  _gcry_mpi_point_free_parts ((a))
36
37
38 /* Print a point using the log functions.  If CTX is not NULL affine
39    coordinates will be printed.  */
40 void
41 _gcry_mpi_point_log (const char *name, mpi_point_t point, mpi_ec_t ctx)
42 {
43   gcry_mpi_t x, y;
44   char buf[100];
45
46   if (!point)
47     {
48       snprintf (buf, sizeof buf - 1, "%s.*", name);
49       log_mpidump (buf, NULL);
50       return;
51     }
52   snprintf (buf, sizeof buf - 1, "%s.X", name);
53
54   if (ctx)
55     {
56       x = mpi_new (0);
57       y = mpi_new (0);
58     }
59   if (!ctx || _gcry_mpi_ec_get_affine (x, y, point, ctx))
60     {
61       log_mpidump (buf, point->x);
62       buf[strlen(buf)-1] = 'Y';
63       log_mpidump (buf, point->y);
64       buf[strlen(buf)-1] = 'Z';
65       log_mpidump (buf, point->z);
66     }
67   else
68     {
69       buf[strlen(buf)-1] = 'x';
70       log_mpidump (buf, x);
71       buf[strlen(buf)-1] = 'y';
72       log_mpidump (buf, y);
73
74     }
75   if (ctx)
76     {
77       _gcry_mpi_release (x);
78       _gcry_mpi_release (y);
79     }
80 }
81
82
83 /* Create a new point option.  NBITS gives the size in bits of one
84    coordinate; it is only used to pre-allocate some resources and
85    might also be passed as 0 to use a default value.  */
86 mpi_point_t
87 _gcry_mpi_point_new (unsigned int nbits)
88 {
89   mpi_point_t p;
90
91   (void)nbits;  /* Currently not used.  */
92
93   p = xmalloc (sizeof *p);
94   _gcry_mpi_point_init (p);
95   return p;
96 }
97
98
99 /* Release the point object P.  P may be NULL. */
100 void
101 _gcry_mpi_point_release (mpi_point_t p)
102 {
103   if (p)
104     {
105       _gcry_mpi_point_free_parts (p);
106       xfree (p);
107     }
108 }
109
110
111 /* Initialize the fields of a point object.  gcry_mpi_point_free_parts
112    may be used to release the fields.  */
113 void
114 _gcry_mpi_point_init (mpi_point_t p)
115 {
116   p->x = mpi_new (0);
117   p->y = mpi_new (0);
118   p->z = mpi_new (0);
119 }
120
121
122 /* Release the parts of a point object. */
123 void
124 _gcry_mpi_point_free_parts (mpi_point_t p)
125 {
126   mpi_free (p->x); p->x = NULL;
127   mpi_free (p->y); p->y = NULL;
128   mpi_free (p->z); p->z = NULL;
129 }
130
131
132 /* Set the value from S into D.  */
133 static void
134 point_set (mpi_point_t d, mpi_point_t s)
135 {
136   mpi_set (d->x, s->x);
137   mpi_set (d->y, s->y);
138   mpi_set (d->z, s->z);
139 }
140
141
142 /* Return a copy of POINT. */
143 gcry_mpi_point_t
144 _gcry_mpi_point_copy (gcry_mpi_point_t point)
145 {
146   mpi_point_t newpoint;
147
148   newpoint = _gcry_mpi_point_new (0);
149   if (point)
150     point_set (newpoint, point);
151
152   return newpoint;
153 }
154
155
156 static void
157 point_resize (mpi_point_t p, mpi_ec_t ctx)
158 {
159   size_t nlimbs = ctx->p->nlimbs;
160
161   mpi_resize (p->x, nlimbs);
162   p->x->nlimbs = nlimbs;
163   mpi_resize (p->z, nlimbs);
164   p->z->nlimbs = nlimbs;
165
166   if (ctx->model != MPI_EC_MONTGOMERY)
167     {
168       mpi_resize (p->y, nlimbs);
169       p->y->nlimbs = nlimbs;
170     }
171 }
172
173
174 static void
175 point_swap_cond (mpi_point_t d, mpi_point_t s, unsigned long swap,
176                  mpi_ec_t ctx)
177 {
178   mpi_swap_cond (d->x, s->x, swap);
179   if (ctx->model != MPI_EC_MONTGOMERY)
180     mpi_swap_cond (d->y, s->y, swap);
181   mpi_swap_cond (d->z, s->z, swap);
182 }
183
184
185 /* Set the projective coordinates from POINT into X, Y, and Z.  If a
186    coordinate is not required, X, Y, or Z may be passed as NULL.  */
187 void
188 _gcry_mpi_point_get (gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z,
189                      mpi_point_t point)
190 {
191   if (x)
192     mpi_set (x, point->x);
193   if (y)
194     mpi_set (y, point->y);
195   if (z)
196     mpi_set (z, point->z);
197 }
198
199
200 /* Set the projective coordinates from POINT into X, Y, and Z and
201    release POINT.  If a coordinate is not required, X, Y, or Z may be
202    passed as NULL.  */
203 void
204 _gcry_mpi_point_snatch_get (gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z,
205                             mpi_point_t point)
206 {
207   mpi_snatch (x, point->x);
208   mpi_snatch (y, point->y);
209   mpi_snatch (z, point->z);
210   xfree (point);
211 }
212
213
214 /* Set the projective coordinates from X, Y, and Z into POINT.  If a
215    coordinate is given as NULL, the value 0 is stored into point.  If
216    POINT is given as NULL a new point object is allocated.  Returns
217    POINT or the newly allocated point object. */
218 mpi_point_t
219 _gcry_mpi_point_set (mpi_point_t point,
220                      gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z)
221 {
222   if (!point)
223     point = mpi_point_new (0);
224
225   if (x)
226     mpi_set (point->x, x);
227   else
228     mpi_clear (point->x);
229   if (y)
230     mpi_set (point->y, y);
231   else
232     mpi_clear (point->y);
233   if (z)
234     mpi_set (point->z, z);
235   else
236     mpi_clear (point->z);
237
238   return point;
239 }
240
241
242 /* Set the projective coordinates from X, Y, and Z into POINT.  If a
243    coordinate is given as NULL, the value 0 is stored into point.  If
244    POINT is given as NULL a new point object is allocated.  The
245    coordinates X, Y, and Z are released.  Returns POINT or the newly
246    allocated point object. */
247 mpi_point_t
248 _gcry_mpi_point_snatch_set (mpi_point_t point,
249                             gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z)
250 {
251   if (!point)
252     point = mpi_point_new (0);
253
254   if (x)
255     mpi_snatch (point->x, x);
256   else
257     mpi_clear (point->x);
258   if (y)
259     mpi_snatch (point->y, y);
260   else
261     mpi_clear (point->y);
262   if (z)
263     mpi_snatch (point->z, z);
264   else
265     mpi_clear (point->z);
266
267   return point;
268 }
269
270
271 /* W = W mod P.  */
272 static void
273 ec_mod (gcry_mpi_t w, mpi_ec_t ec)
274 {
275   if (0 && ec->dialect == ECC_DIALECT_ED25519)
276     _gcry_mpi_ec_ed25519_mod (w);
277   else if (ec->t.p_barrett)
278     _gcry_mpi_mod_barrett (w, w, ec->t.p_barrett);
279   else
280     _gcry_mpi_mod (w, w, ec->p);
281 }
282
283 static void
284 ec_addm (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
285 {
286   mpi_add (w, u, v);
287   ec_mod (w, ctx);
288 }
289
290 static void
291 ec_subm (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ec)
292 {
293   mpi_sub (w, u, v);
294   while (w->sign)
295     mpi_add (w, w, ec->p);
296   /*ec_mod (w, ec);*/
297 }
298
299 static void
300 ec_mulm (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
301 {
302   mpi_mul (w, u, v);
303   ec_mod (w, ctx);
304 }
305
306 /* W = 2 * U mod P.  */
307 static void
308 ec_mul2 (gcry_mpi_t w, gcry_mpi_t u, mpi_ec_t ctx)
309 {
310   mpi_lshift (w, u, 1);
311   ec_mod (w, ctx);
312 }
313
314 static void
315 ec_powm (gcry_mpi_t w, const gcry_mpi_t b, const gcry_mpi_t e,
316          mpi_ec_t ctx)
317 {
318   mpi_powm (w, b, e, ctx->p);
319   /* _gcry_mpi_abs (w); */
320 }
321
322
323 /* Shortcut for
324      ec_powm (B, B, mpi_const (MPI_C_TWO), ctx);
325    for easier optimization.  */
326 static void
327 ec_pow2 (gcry_mpi_t w, const gcry_mpi_t b, mpi_ec_t ctx)
328 {
329   /* Using mpi_mul is slightly faster (at least on amd64).  */
330   /* mpi_powm (w, b, mpi_const (MPI_C_TWO), ctx->p); */
331   ec_mulm (w, b, b, ctx);
332 }
333
334
335 /* Shortcut for
336      ec_powm (B, B, mpi_const (MPI_C_THREE), ctx);
337    for easier optimization.  */
338 static void
339 ec_pow3 (gcry_mpi_t w, const gcry_mpi_t b, mpi_ec_t ctx)
340 {
341   mpi_powm (w, b, mpi_const (MPI_C_THREE), ctx->p);
342 }
343
344
345 static void
346 ec_invm (gcry_mpi_t x, gcry_mpi_t a, mpi_ec_t ctx)
347 {
348   if (!mpi_invm (x, a, ctx->p))
349     {
350       log_error ("ec_invm: inverse does not exist:\n");
351       log_mpidump ("  a", a);
352       log_mpidump ("  p", ctx->p);
353     }
354 }
355 \f
356 static void
357 mpih_set_cond (mpi_ptr_t wp, mpi_ptr_t up, mpi_size_t usize, unsigned long set)
358 {
359   mpi_size_t i;
360   mpi_limb_t mask = ((mpi_limb_t)0) - set;
361   mpi_limb_t x;
362
363   for (i = 0; i < usize; i++)
364     {
365       x = mask & (wp[i] ^ up[i]);
366       wp[i] = wp[i] ^ x;
367     }
368 }
369
370 /* Routines for 2^255 - 19.  */
371
372 #define LIMB_SIZE_25519 ((256+BITS_PER_MPI_LIMB-1)/BITS_PER_MPI_LIMB)
373
374 static void
375 ec_addm_25519 (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
376 {
377   mpi_ptr_t wp, up, vp;
378   mpi_size_t wsize = LIMB_SIZE_25519;
379   mpi_limb_t n[LIMB_SIZE_25519];
380   mpi_limb_t borrow;
381
382   if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
383     log_bug ("addm_25519: different sizes\n");
384
385   memset (n, 0, sizeof n);
386   up = u->d;
387   vp = v->d;
388   wp = w->d;
389
390   _gcry_mpih_add_n (wp, up, vp, wsize);
391   borrow = _gcry_mpih_sub_n (wp, wp, ctx->p->d, wsize);
392   mpih_set_cond (n, ctx->p->d, wsize, (borrow != 0UL));
393   _gcry_mpih_add_n (wp, wp, n, wsize);
394   wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
395 }
396
397 static void
398 ec_subm_25519 (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
399 {
400   mpi_ptr_t wp, up, vp;
401   mpi_size_t wsize = LIMB_SIZE_25519;
402   mpi_limb_t n[LIMB_SIZE_25519];
403   mpi_limb_t borrow;
404
405   if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
406     log_bug ("subm_25519: different sizes\n");
407
408   memset (n, 0, sizeof n);
409   up = u->d;
410   vp = v->d;
411   wp = w->d;
412
413   borrow = _gcry_mpih_sub_n (wp, up, vp, wsize);
414   mpih_set_cond (n, ctx->p->d, wsize, (borrow != 0UL));
415   _gcry_mpih_add_n (wp, wp, n, wsize);
416   wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
417 }
418
419 static void
420 ec_mulm_25519 (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
421 {
422   mpi_ptr_t wp, up, vp;
423   mpi_size_t wsize = LIMB_SIZE_25519;
424   mpi_limb_t n[LIMB_SIZE_25519*2];
425   mpi_limb_t m[LIMB_SIZE_25519+1];
426   mpi_limb_t cy;
427   int msb;
428
429   (void)ctx;
430   if (w->nlimbs != wsize || u->nlimbs != wsize || v->nlimbs != wsize)
431     log_bug ("mulm_25519: different sizes\n");
432
433   up = u->d;
434   vp = v->d;
435   wp = w->d;
436
437   _gcry_mpih_mul_n (n, up, vp, wsize);
438   memcpy (wp, n, wsize * BYTES_PER_MPI_LIMB);
439   wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
440
441   memcpy (m, n+LIMB_SIZE_25519-1, (wsize+1) * BYTES_PER_MPI_LIMB);
442   _gcry_mpih_rshift (m, m, LIMB_SIZE_25519+1, (255 % BITS_PER_MPI_LIMB));
443
444   memcpy (n, m, wsize * BYTES_PER_MPI_LIMB);
445   cy = _gcry_mpih_lshift (m, m, LIMB_SIZE_25519, 4);
446   m[LIMB_SIZE_25519] = cy;
447   cy = _gcry_mpih_add_n (m, m, n, wsize);
448   m[LIMB_SIZE_25519] += cy;
449   cy = _gcry_mpih_add_n (m, m, n, wsize);
450   m[LIMB_SIZE_25519] += cy;
451   cy = _gcry_mpih_add_n (m, m, n, wsize);
452   m[LIMB_SIZE_25519] += cy;
453
454   cy = _gcry_mpih_add_n (wp, wp, m, wsize);
455   m[LIMB_SIZE_25519] += cy;
456
457   memset (m, 0, wsize * BYTES_PER_MPI_LIMB);
458   msb = (wp[LIMB_SIZE_25519-1] >> (255 % BITS_PER_MPI_LIMB));
459   m[0] = (m[LIMB_SIZE_25519] * 2 + msb) * 19;
460   wp[LIMB_SIZE_25519-1] &= ~((mpi_limb_t)1 << (255 % BITS_PER_MPI_LIMB));
461   _gcry_mpih_add_n (wp, wp, m, wsize);
462
463   m[0] = 0;
464   cy = _gcry_mpih_sub_n (wp, wp, ctx->p->d, wsize);
465   mpih_set_cond (m, ctx->p->d, wsize, (cy != 0UL));
466   _gcry_mpih_add_n (wp, wp, m, wsize);
467 }
468
469 static void
470 ec_mul2_25519 (gcry_mpi_t w, gcry_mpi_t u, mpi_ec_t ctx)
471 {
472   ec_addm_25519 (w, u, u, ctx);
473 }
474
475 static void
476 ec_pow2_25519 (gcry_mpi_t w, const gcry_mpi_t b, mpi_ec_t ctx)
477 {
478   ec_mulm_25519 (w, b, b, ctx);
479 }
480
481 struct field_table {
482   const char *p;
483
484   /* computation routines for the field.  */
485   void (* addm) (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx);
486   void (* subm) (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx);
487   void (* mulm) (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx);
488   void (* mul2) (gcry_mpi_t w, gcry_mpi_t u, mpi_ec_t ctx);
489   void (* pow2) (gcry_mpi_t w, const gcry_mpi_t b, mpi_ec_t ctx);
490 };
491
492 static const struct field_table field_table[] = {
493   {
494     "0x7FFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFFED",
495     ec_addm_25519,
496     ec_subm_25519,
497     ec_mulm_25519,
498     ec_mul2_25519,
499     ec_pow2_25519
500   },
501   { NULL, NULL, NULL, NULL, NULL, NULL },
502 };
503 \f
504 /* Force recomputation of all helper variables.  */
505 void
506 _gcry_mpi_ec_get_reset (mpi_ec_t ec)
507 {
508   ec->t.valid.a_is_pminus3 = 0;
509   ec->t.valid.two_inv_p = 0;
510 }
511
512
513 /* Accessor for helper variable.  */
514 static int
515 ec_get_a_is_pminus3 (mpi_ec_t ec)
516 {
517   gcry_mpi_t tmp;
518
519   if (!ec->t.valid.a_is_pminus3)
520     {
521       ec->t.valid.a_is_pminus3 = 1;
522       tmp = mpi_alloc_like (ec->p);
523       mpi_sub_ui (tmp, ec->p, 3);
524       ec->t.a_is_pminus3 = !mpi_cmp (ec->a, tmp);
525       mpi_free (tmp);
526     }
527
528   return ec->t.a_is_pminus3;
529 }
530
531
532 /* Accessor for helper variable.  */
533 static gcry_mpi_t
534 ec_get_two_inv_p (mpi_ec_t ec)
535 {
536   if (!ec->t.valid.two_inv_p)
537     {
538       ec->t.valid.two_inv_p = 1;
539       if (!ec->t.two_inv_p)
540         ec->t.two_inv_p = mpi_alloc (0);
541       ec_invm (ec->t.two_inv_p, mpi_const (MPI_C_TWO), ec);
542     }
543   return ec->t.two_inv_p;
544 }
545
546
547 static const char *curve25519_bad_points[] = {
548   "0x0000000000000000000000000000000000000000000000000000000000000000",
549   "0x0000000000000000000000000000000000000000000000000000000000000001",
550   "0x00b8495f16056286fdb1329ceb8d09da6ac49ff1fae35616aeb8413b7c7aebe0",
551   "0x57119fd0dd4e22d8868e1c58c45c44045bef839c55b1d0b1248c50a3bc959c5f",
552   "0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffec",
553   "0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffed",
554   "0x7fffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffee",
555   NULL
556 };
557
558 static gcry_mpi_t
559 scanval (const char *string)
560 {
561   gpg_err_code_t rc;
562   gcry_mpi_t val;
563
564   rc = _gcry_mpi_scan (&val, GCRYMPI_FMT_HEX, string, 0, NULL);
565   if (rc)
566     log_fatal ("scanning ECC parameter failed: %s\n", gpg_strerror (rc));
567   return val;
568 }
569
570
571 /* This function initialized a context for elliptic curve based on the
572    field GF(p).  P is the prime specifying this field, A is the first
573    coefficient.  CTX is expected to be zeroized.  */
574 static void
575 ec_p_init (mpi_ec_t ctx, enum gcry_mpi_ec_models model,
576            enum ecc_dialects dialect,
577            int flags,
578            gcry_mpi_t p, gcry_mpi_t a, gcry_mpi_t b)
579 {
580   int i;
581   static int use_barrett;
582
583   if (!use_barrett)
584     {
585       if (getenv ("GCRYPT_BARRETT"))
586         use_barrett = 1;
587       else
588         use_barrett = -1;
589     }
590
591   /* Fixme: Do we want to check some constraints? e.g.  a < p  */
592
593   ctx->model = model;
594   ctx->dialect = dialect;
595   ctx->flags = flags;
596   if (dialect == ECC_DIALECT_ED25519)
597     ctx->nbits = 256;
598   else
599     ctx->nbits = mpi_get_nbits (p);
600   ctx->p = mpi_copy (p);
601   ctx->a = mpi_copy (a);
602   ctx->b = mpi_copy (b);
603
604   ctx->t.p_barrett = use_barrett > 0? _gcry_mpi_barrett_init (ctx->p, 0):NULL;
605
606   _gcry_mpi_ec_get_reset (ctx);
607
608   if (model == MPI_EC_MONTGOMERY)
609     {
610       for (i=0; i< DIM(ctx->t.scratch) && curve25519_bad_points[i]; i++)
611         ctx->t.scratch[i] = scanval (curve25519_bad_points[i]);
612     }
613   else
614     {
615       /* Allocate scratch variables.  */
616       for (i=0; i< DIM(ctx->t.scratch); i++)
617         ctx->t.scratch[i] = mpi_alloc_like (ctx->p);
618     }
619
620   ctx->addm = ec_addm;
621   ctx->subm = ec_subm;
622   ctx->mulm = ec_mulm;
623   ctx->mul2 = ec_mul2;
624   ctx->pow2 = ec_pow2;
625
626   for (i=0; field_table[i].p; i++)
627     {
628       gcry_mpi_t f_p;
629       gpg_err_code_t rc;
630
631       rc = _gcry_mpi_scan (&f_p, GCRYMPI_FMT_HEX, field_table[i].p, 0, NULL);
632       if (rc)
633         log_fatal ("scanning ECC parameter failed: %s\n", gpg_strerror (rc));
634
635       if (!mpi_cmp (p, f_p))
636         {
637           ctx->addm = field_table[i].addm;
638           ctx->subm = field_table[i].subm;
639           ctx->mulm = field_table[i].mulm;
640           ctx->mul2 = field_table[i].mul2;
641           ctx->pow2 = field_table[i].pow2;
642           _gcry_mpi_release (f_p);
643
644           mpi_resize (ctx->a, ctx->p->nlimbs);
645           ctx->a->nlimbs = ctx->p->nlimbs;
646
647           mpi_resize (ctx->b, ctx->p->nlimbs);
648           ctx->b->nlimbs = ctx->p->nlimbs;
649
650           for (i=0; i< DIM(ctx->t.scratch) && ctx->t.scratch[i]; i++)
651             ctx->t.scratch[i]->nlimbs = ctx->p->nlimbs;
652
653           break;
654         }
655
656       _gcry_mpi_release (f_p);
657     }
658
659   /* Prepare for fast reduction.  */
660   /* FIXME: need a test for NIST values.  However it does not gain us
661      any real advantage, for 384 bits it is actually slower than using
662      mpi_mulm.  */
663 /*   ctx->nist_nbits = mpi_get_nbits (ctx->p); */
664 /*   if (ctx->nist_nbits == 192) */
665 /*     { */
666 /*       for (i=0; i < 4; i++) */
667 /*         ctx->s[i] = mpi_new (192); */
668 /*       ctx->c    = mpi_new (192*2); */
669 /*     } */
670 /*   else if (ctx->nist_nbits == 384) */
671 /*     { */
672 /*       for (i=0; i < 10; i++) */
673 /*         ctx->s[i] = mpi_new (384); */
674 /*       ctx->c    = mpi_new (384*2); */
675 /*     } */
676 }
677
678
679 static void
680 ec_deinit (void *opaque)
681 {
682   mpi_ec_t ctx = opaque;
683   int i;
684
685   _gcry_mpi_barrett_free (ctx->t.p_barrett);
686
687   /* Domain parameter.  */
688   mpi_free (ctx->p);
689   mpi_free (ctx->a);
690   mpi_free (ctx->b);
691   _gcry_mpi_point_release (ctx->G);
692   mpi_free (ctx->n);
693   mpi_free (ctx->h);
694
695   /* The key.  */
696   _gcry_mpi_point_release (ctx->Q);
697   mpi_free (ctx->d);
698
699   /* Private data of ec.c.  */
700   mpi_free (ctx->t.two_inv_p);
701
702   for (i=0; i< DIM(ctx->t.scratch); i++)
703     mpi_free (ctx->t.scratch[i]);
704
705 /*   if (ctx->nist_nbits == 192) */
706 /*     { */
707 /*       for (i=0; i < 4; i++) */
708 /*         mpi_free (ctx->s[i]); */
709 /*       mpi_free (ctx->c); */
710 /*     } */
711 /*   else if (ctx->nist_nbits == 384) */
712 /*     { */
713 /*       for (i=0; i < 10; i++) */
714 /*         mpi_free (ctx->s[i]); */
715 /*       mpi_free (ctx->c); */
716 /*     } */
717 }
718
719
720 /* This function returns a new context for elliptic curve based on the
721    field GF(p).  P is the prime specifying this field, A is the first
722    coefficient, B is the second coefficient, and MODEL is the model
723    for the curve.  This function is only used within Libgcrypt and not
724    part of the public API.
725
726    This context needs to be released using _gcry_mpi_ec_free.  */
727 mpi_ec_t
728 _gcry_mpi_ec_p_internal_new (enum gcry_mpi_ec_models model,
729                              enum ecc_dialects dialect,
730                              int flags,
731                              gcry_mpi_t p, gcry_mpi_t a, gcry_mpi_t b)
732 {
733   mpi_ec_t ctx;
734
735   ctx = xcalloc (1, sizeof *ctx);
736   ec_p_init (ctx, model, dialect, flags, p, a, b);
737
738   return ctx;
739 }
740
741
742 /* This is a variant of _gcry_mpi_ec_p_internal_new which returns an
743    public context and does some error checking on the supplied
744    arguments.  On success the new context is stored at R_CTX and 0 is
745    returned; on error NULL is stored at R_CTX and an error code is
746    returned.
747
748    The context needs to be released using gcry_ctx_release.  */
749 gpg_err_code_t
750 _gcry_mpi_ec_p_new (gcry_ctx_t *r_ctx,
751                     enum gcry_mpi_ec_models model,
752                     enum ecc_dialects dialect,
753                     int flags,
754                     gcry_mpi_t p, gcry_mpi_t a, gcry_mpi_t b)
755 {
756   gcry_ctx_t ctx;
757   mpi_ec_t ec;
758
759   *r_ctx = NULL;
760   if (!p || !a)
761     return GPG_ERR_EINVAL;
762
763   ctx = _gcry_ctx_alloc (CONTEXT_TYPE_EC, sizeof *ec, ec_deinit);
764   if (!ctx)
765     return gpg_err_code_from_syserror ();
766   ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
767   ec_p_init (ec, model, dialect, flags, p, a, b);
768
769   *r_ctx = ctx;
770   return 0;
771 }
772
773
774 void
775 _gcry_mpi_ec_free (mpi_ec_t ctx)
776 {
777   if (ctx)
778     {
779       ec_deinit (ctx);
780       xfree (ctx);
781     }
782 }
783
784
785 gcry_mpi_t
786 _gcry_mpi_ec_get_mpi (const char *name, gcry_ctx_t ctx, int copy)
787 {
788   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
789
790   return _gcry_ecc_get_mpi (name, ec, copy);
791 }
792
793
794 gcry_mpi_point_t
795 _gcry_mpi_ec_get_point (const char *name, gcry_ctx_t ctx, int copy)
796 {
797   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
798
799   (void)copy;  /* Not used.  */
800
801   return _gcry_ecc_get_point (name, ec);
802 }
803
804
805 gpg_err_code_t
806 _gcry_mpi_ec_set_mpi (const char *name, gcry_mpi_t newvalue,
807                       gcry_ctx_t ctx)
808 {
809   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
810
811   return _gcry_ecc_set_mpi (name, newvalue, ec);
812 }
813
814
815 gpg_err_code_t
816 _gcry_mpi_ec_set_point (const char *name, gcry_mpi_point_t newvalue,
817                         gcry_ctx_t ctx)
818 {
819   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
820
821   return _gcry_ecc_set_point (name, newvalue, ec);
822 }
823
824
825 /* Given an encoded point in the MPI VALUE and a context EC, decode
826  * the point according to the context and store it in RESULT.  On
827  * error an error code is return but RESULT might have been changed.
828  * If no context is given the function tries to decode VALUE by
829  * assuming a 0x04 prefixed uncompressed encoding.  */
830 gpg_err_code_t
831 _gcry_mpi_ec_decode_point (mpi_point_t result, gcry_mpi_t value, mpi_ec_t ec)
832 {
833   gcry_err_code_t rc;
834
835   if (ec && ec->dialect == ECC_DIALECT_ED25519)
836     rc = _gcry_ecc_eddsa_decodepoint (value, ec, result, NULL, NULL);
837   else if (ec && ec->model == MPI_EC_MONTGOMERY)
838     rc = _gcry_ecc_mont_decodepoint (value, ec, result);
839   else
840     rc = _gcry_ecc_os2ec (result, value);
841
842   return rc;
843 }
844
845
846 /* Compute the affine coordinates from the projective coordinates in
847    POINT.  Set them into X and Y.  If one coordinate is not required,
848    X or Y may be passed as NULL.  CTX is the usual context. Returns: 0
849    on success or !0 if POINT is at infinity.  */
850 int
851 _gcry_mpi_ec_get_affine (gcry_mpi_t x, gcry_mpi_t y, mpi_point_t point,
852                          mpi_ec_t ctx)
853 {
854   if (!mpi_cmp_ui (point->z, 0))
855     return -1;
856
857   switch (ctx->model)
858     {
859     case MPI_EC_WEIERSTRASS: /* Using Jacobian coordinates.  */
860       {
861         gcry_mpi_t z1, z2, z3;
862
863         z1 = mpi_new (0);
864         z2 = mpi_new (0);
865         ec_invm (z1, point->z, ctx);  /* z1 = z^(-1) mod p  */
866         ec_mulm (z2, z1, z1, ctx);    /* z2 = z^(-2) mod p  */
867
868         if (x)
869           ec_mulm (x, point->x, z2, ctx);
870
871         if (y)
872           {
873             z3 = mpi_new (0);
874             ec_mulm (z3, z2, z1, ctx);      /* z3 = z^(-3) mod p  */
875             ec_mulm (y, point->y, z3, ctx);
876             mpi_free (z3);
877           }
878
879         mpi_free (z2);
880         mpi_free (z1);
881       }
882       return 0;
883
884     case MPI_EC_MONTGOMERY:
885       {
886         if (x)
887           mpi_set (x, point->x);
888
889         if (y)
890           {
891             log_fatal ("%s: Getting Y-coordinate on %s is not supported\n",
892                        "_gcry_mpi_ec_get_affine", "Montgomery");
893             return -1;
894           }
895       }
896       return 0;
897
898     case MPI_EC_EDWARDS:
899       {
900         gcry_mpi_t z;
901
902         z = mpi_new (0);
903         ec_invm (z, point->z, ctx);
904
905         mpi_resize (z, ctx->p->nlimbs);
906         z->nlimbs = ctx->p->nlimbs;
907
908         if (x)
909           {
910             mpi_resize (x, ctx->p->nlimbs);
911             x->nlimbs = ctx->p->nlimbs;
912             ctx->mulm (x, point->x, z, ctx);
913           }
914         if (y)
915           {
916             mpi_resize (y, ctx->p->nlimbs);
917             y->nlimbs = ctx->p->nlimbs;
918             ctx->mulm (y, point->y, z, ctx);
919           }
920
921         _gcry_mpi_release (z);
922       }
923       return 0;
924
925     default:
926       return -1;
927     }
928 }
929
930
931 \f
932 /*  RESULT = 2 * POINT  (Weierstrass version). */
933 static void
934 dup_point_weierstrass (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
935 {
936 #define x3 (result->x)
937 #define y3 (result->y)
938 #define z3 (result->z)
939 #define t1 (ctx->t.scratch[0])
940 #define t2 (ctx->t.scratch[1])
941 #define t3 (ctx->t.scratch[2])
942 #define l1 (ctx->t.scratch[3])
943 #define l2 (ctx->t.scratch[4])
944 #define l3 (ctx->t.scratch[5])
945
946   if (!mpi_cmp_ui (point->y, 0) || !mpi_cmp_ui (point->z, 0))
947     {
948       /* P_y == 0 || P_z == 0 => [1:1:0] */
949       mpi_set_ui (x3, 1);
950       mpi_set_ui (y3, 1);
951       mpi_set_ui (z3, 0);
952     }
953   else
954     {
955       if (ec_get_a_is_pminus3 (ctx))  /* Use the faster case.  */
956         {
957           /* L1 = 3(X - Z^2)(X + Z^2) */
958           /*                          T1: used for Z^2. */
959           /*                          T2: used for the right term.  */
960           ec_pow2 (t1, point->z, ctx);
961           ec_subm (l1, point->x, t1, ctx);
962           ec_mulm (l1, l1, mpi_const (MPI_C_THREE), ctx);
963           ec_addm (t2, point->x, t1, ctx);
964           ec_mulm (l1, l1, t2, ctx);
965         }
966       else /* Standard case. */
967         {
968           /* L1 = 3X^2 + aZ^4 */
969           /*                          T1: used for aZ^4. */
970           ec_pow2 (l1, point->x, ctx);
971           ec_mulm (l1, l1, mpi_const (MPI_C_THREE), ctx);
972           ec_powm (t1, point->z, mpi_const (MPI_C_FOUR), ctx);
973           ec_mulm (t1, t1, ctx->a, ctx);
974           ec_addm (l1, l1, t1, ctx);
975         }
976       /* Z3 = 2YZ */
977       ec_mulm (z3, point->y, point->z, ctx);
978       ec_mul2 (z3, z3, ctx);
979
980       /* L2 = 4XY^2 */
981       /*                              T2: used for Y2; required later. */
982       ec_pow2 (t2, point->y, ctx);
983       ec_mulm (l2, t2, point->x, ctx);
984       ec_mulm (l2, l2, mpi_const (MPI_C_FOUR), ctx);
985
986       /* X3 = L1^2 - 2L2 */
987       /*                              T1: used for L2^2. */
988       ec_pow2 (x3, l1, ctx);
989       ec_mul2 (t1, l2, ctx);
990       ec_subm (x3, x3, t1, ctx);
991
992       /* L3 = 8Y^4 */
993       /*                              T2: taken from above. */
994       ec_pow2 (t2, t2, ctx);
995       ec_mulm (l3, t2, mpi_const (MPI_C_EIGHT), ctx);
996
997       /* Y3 = L1(L2 - X3) - L3 */
998       ec_subm (y3, l2, x3, ctx);
999       ec_mulm (y3, y3, l1, ctx);
1000       ec_subm (y3, y3, l3, ctx);
1001     }
1002
1003 #undef x3
1004 #undef y3
1005 #undef z3
1006 #undef t1
1007 #undef t2
1008 #undef t3
1009 #undef l1
1010 #undef l2
1011 #undef l3
1012 }
1013
1014
1015 /*  RESULT = 2 * POINT  (Montgomery version). */
1016 static void
1017 dup_point_montgomery (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
1018 {
1019   (void)result;
1020   (void)point;
1021   (void)ctx;
1022   log_fatal ("%s: %s not yet supported\n",
1023              "_gcry_mpi_ec_dup_point", "Montgomery");
1024 }
1025
1026
1027 /*  RESULT = 2 * POINT  (Twisted Edwards version). */
1028 static void
1029 dup_point_edwards (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
1030 {
1031 #define X1 (point->x)
1032 #define Y1 (point->y)
1033 #define Z1 (point->z)
1034 #define X3 (result->x)
1035 #define Y3 (result->y)
1036 #define Z3 (result->z)
1037 #define B (ctx->t.scratch[0])
1038 #define C (ctx->t.scratch[1])
1039 #define D (ctx->t.scratch[2])
1040 #define E (ctx->t.scratch[3])
1041 #define F (ctx->t.scratch[4])
1042 #define H (ctx->t.scratch[5])
1043 #define J (ctx->t.scratch[6])
1044
1045   /* Compute: (X_3 : Y_3 : Z_3) = 2( X_1 : Y_1 : Z_1 ) */
1046
1047   /* B = (X_1 + Y_1)^2  */
1048   ctx->addm (B, X1, Y1, ctx);
1049   ctx->pow2 (B, B, ctx);
1050
1051   /* C = X_1^2 */
1052   /* D = Y_1^2 */
1053   ctx->pow2 (C, X1, ctx);
1054   ctx->pow2 (D, Y1, ctx);
1055
1056   /* E = aC */
1057   if (ctx->dialect == ECC_DIALECT_ED25519)
1058     ctx->subm (E, ctx->p, C, ctx);
1059   else
1060     ctx->mulm (E, ctx->a, C, ctx);
1061
1062   /* F = E + D */
1063   ctx->addm (F, E, D, ctx);
1064
1065   /* H = Z_1^2 */
1066   ctx->pow2 (H, Z1, ctx);
1067
1068   /* J = F - 2H */
1069   ctx->mul2 (J, H, ctx);
1070   ctx->subm (J, F, J, ctx);
1071
1072   /* X_3 = (B - C - D) · J */
1073   ctx->subm (X3, B, C, ctx);
1074   ctx->subm (X3, X3, D, ctx);
1075   ctx->mulm (X3, X3, J, ctx);
1076
1077   /* Y_3 = F · (E - D) */
1078   ctx->subm (Y3, E, D, ctx);
1079   ctx->mulm (Y3, Y3, F, ctx);
1080
1081   /* Z_3 = F · J */
1082   ctx->mulm (Z3, F, J, ctx);
1083
1084 #undef X1
1085 #undef Y1
1086 #undef Z1
1087 #undef X3
1088 #undef Y3
1089 #undef Z3
1090 #undef B
1091 #undef C
1092 #undef D
1093 #undef E
1094 #undef F
1095 #undef H
1096 #undef J
1097 }
1098
1099
1100 /*  RESULT = 2 * POINT  */
1101 void
1102 _gcry_mpi_ec_dup_point (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
1103 {
1104   switch (ctx->model)
1105     {
1106     case MPI_EC_WEIERSTRASS:
1107       dup_point_weierstrass (result, point, ctx);
1108       break;
1109     case MPI_EC_MONTGOMERY:
1110       dup_point_montgomery (result, point, ctx);
1111       break;
1112     case MPI_EC_EDWARDS:
1113       dup_point_edwards (result, point, ctx);
1114       break;
1115     }
1116 }
1117
1118
1119 /* RESULT = P1 + P2  (Weierstrass version).*/
1120 static void
1121 add_points_weierstrass (mpi_point_t result,
1122                         mpi_point_t p1, mpi_point_t p2,
1123                         mpi_ec_t ctx)
1124 {
1125 #define x1 (p1->x    )
1126 #define y1 (p1->y    )
1127 #define z1 (p1->z    )
1128 #define x2 (p2->x    )
1129 #define y2 (p2->y    )
1130 #define z2 (p2->z    )
1131 #define x3 (result->x)
1132 #define y3 (result->y)
1133 #define z3 (result->z)
1134 #define l1 (ctx->t.scratch[0])
1135 #define l2 (ctx->t.scratch[1])
1136 #define l3 (ctx->t.scratch[2])
1137 #define l4 (ctx->t.scratch[3])
1138 #define l5 (ctx->t.scratch[4])
1139 #define l6 (ctx->t.scratch[5])
1140 #define l7 (ctx->t.scratch[6])
1141 #define l8 (ctx->t.scratch[7])
1142 #define l9 (ctx->t.scratch[8])
1143 #define t1 (ctx->t.scratch[9])
1144 #define t2 (ctx->t.scratch[10])
1145
1146   if ( (!mpi_cmp (x1, x2)) && (!mpi_cmp (y1, y2)) && (!mpi_cmp (z1, z2)) )
1147     {
1148       /* Same point; need to call the duplicate function.  */
1149       _gcry_mpi_ec_dup_point (result, p1, ctx);
1150     }
1151   else if (!mpi_cmp_ui (z1, 0))
1152     {
1153       /* P1 is at infinity.  */
1154       mpi_set (x3, p2->x);
1155       mpi_set (y3, p2->y);
1156       mpi_set (z3, p2->z);
1157     }
1158   else if (!mpi_cmp_ui (z2, 0))
1159     {
1160       /* P2 is at infinity.  */
1161       mpi_set (x3, p1->x);
1162       mpi_set (y3, p1->y);
1163       mpi_set (z3, p1->z);
1164     }
1165   else
1166     {
1167       int z1_is_one = !mpi_cmp_ui (z1, 1);
1168       int z2_is_one = !mpi_cmp_ui (z2, 1);
1169
1170       /* l1 = x1 z2^2  */
1171       /* l2 = x2 z1^2  */
1172       if (z2_is_one)
1173         mpi_set (l1, x1);
1174       else
1175         {
1176           ec_pow2 (l1, z2, ctx);
1177           ec_mulm (l1, l1, x1, ctx);
1178         }
1179       if (z1_is_one)
1180         mpi_set (l2, x2);
1181       else
1182         {
1183           ec_pow2 (l2, z1, ctx);
1184           ec_mulm (l2, l2, x2, ctx);
1185         }
1186       /* l3 = l1 - l2 */
1187       ec_subm (l3, l1, l2, ctx);
1188       /* l4 = y1 z2^3  */
1189       ec_powm (l4, z2, mpi_const (MPI_C_THREE), ctx);
1190       ec_mulm (l4, l4, y1, ctx);
1191       /* l5 = y2 z1^3  */
1192       ec_powm (l5, z1, mpi_const (MPI_C_THREE), ctx);
1193       ec_mulm (l5, l5, y2, ctx);
1194       /* l6 = l4 - l5  */
1195       ec_subm (l6, l4, l5, ctx);
1196
1197       if (!mpi_cmp_ui (l3, 0))
1198         {
1199           if (!mpi_cmp_ui (l6, 0))
1200             {
1201               /* P1 and P2 are the same - use duplicate function.  */
1202               _gcry_mpi_ec_dup_point (result, p1, ctx);
1203             }
1204           else
1205             {
1206               /* P1 is the inverse of P2.  */
1207               mpi_set_ui (x3, 1);
1208               mpi_set_ui (y3, 1);
1209               mpi_set_ui (z3, 0);
1210             }
1211         }
1212       else
1213         {
1214           /* l7 = l1 + l2  */
1215           ec_addm (l7, l1, l2, ctx);
1216           /* l8 = l4 + l5  */
1217           ec_addm (l8, l4, l5, ctx);
1218           /* z3 = z1 z2 l3  */
1219           ec_mulm (z3, z1, z2, ctx);
1220           ec_mulm (z3, z3, l3, ctx);
1221           /* x3 = l6^2 - l7 l3^2  */
1222           ec_pow2 (t1, l6, ctx);
1223           ec_pow2 (t2, l3, ctx);
1224           ec_mulm (t2, t2, l7, ctx);
1225           ec_subm (x3, t1, t2, ctx);
1226           /* l9 = l7 l3^2 - 2 x3  */
1227           ec_mul2 (t1, x3, ctx);
1228           ec_subm (l9, t2, t1, ctx);
1229           /* y3 = (l9 l6 - l8 l3^3)/2  */
1230           ec_mulm (l9, l9, l6, ctx);
1231           ec_powm (t1, l3, mpi_const (MPI_C_THREE), ctx); /* fixme: Use saved value*/
1232           ec_mulm (t1, t1, l8, ctx);
1233           ec_subm (y3, l9, t1, ctx);
1234           ec_mulm (y3, y3, ec_get_two_inv_p (ctx), ctx);
1235         }
1236     }
1237
1238 #undef x1
1239 #undef y1
1240 #undef z1
1241 #undef x2
1242 #undef y2
1243 #undef z2
1244 #undef x3
1245 #undef y3
1246 #undef z3
1247 #undef l1
1248 #undef l2
1249 #undef l3
1250 #undef l4
1251 #undef l5
1252 #undef l6
1253 #undef l7
1254 #undef l8
1255 #undef l9
1256 #undef t1
1257 #undef t2
1258 }
1259
1260
1261 /* RESULT = P1 + P2  (Montgomery version).*/
1262 static void
1263 add_points_montgomery (mpi_point_t result,
1264                        mpi_point_t p1, mpi_point_t p2,
1265                        mpi_ec_t ctx)
1266 {
1267   (void)result;
1268   (void)p1;
1269   (void)p2;
1270   (void)ctx;
1271   log_fatal ("%s: %s not yet supported\n",
1272              "_gcry_mpi_ec_add_points", "Montgomery");
1273 }
1274
1275
1276 /* RESULT = P1 + P2  (Twisted Edwards version).*/
1277 static void
1278 add_points_edwards (mpi_point_t result,
1279                     mpi_point_t p1, mpi_point_t p2,
1280                     mpi_ec_t ctx)
1281 {
1282 #define X1 (p1->x)
1283 #define Y1 (p1->y)
1284 #define Z1 (p1->z)
1285 #define X2 (p2->x)
1286 #define Y2 (p2->y)
1287 #define Z2 (p2->z)
1288 #define X3 (result->x)
1289 #define Y3 (result->y)
1290 #define Z3 (result->z)
1291 #define A (ctx->t.scratch[0])
1292 #define B (ctx->t.scratch[1])
1293 #define C (ctx->t.scratch[2])
1294 #define D (ctx->t.scratch[3])
1295 #define E (ctx->t.scratch[4])
1296 #define F (ctx->t.scratch[5])
1297 #define G (ctx->t.scratch[6])
1298 #define tmp (ctx->t.scratch[7])
1299
1300   point_resize (result, ctx);
1301
1302   /* Compute: (X_3 : Y_3 : Z_3) = (X_1 : Y_1 : Z_1) + (X_2 : Y_2 : Z_3)  */
1303
1304   /* A = Z1 · Z2 */
1305   ctx->mulm (A, Z1, Z2, ctx);
1306
1307   /* B = A^2 */
1308   ctx->pow2 (B, A, ctx);
1309
1310   /* C = X1 · X2 */
1311   ctx->mulm (C, X1, X2, ctx);
1312
1313   /* D = Y1 · Y2 */
1314   ctx->mulm (D, Y1, Y2, ctx);
1315
1316   /* E = d · C · D */
1317   ctx->mulm (E, ctx->b, C, ctx);
1318   ctx->mulm (E, E, D, ctx);
1319
1320   /* F = B - E */
1321   ctx->subm (F, B, E, ctx);
1322
1323   /* G = B + E */
1324   ctx->addm (G, B, E, ctx);
1325
1326   /* X_3 = A · F · ((X_1 + Y_1) · (X_2 + Y_2) - C - D) */
1327   ctx->addm (tmp, X1, Y1, ctx);
1328   ctx->addm (X3, X2, Y2, ctx);
1329   ctx->mulm (X3, X3, tmp, ctx);
1330   ctx->subm (X3, X3, C, ctx);
1331   ctx->subm (X3, X3, D, ctx);
1332   ctx->mulm (X3, X3, F, ctx);
1333   ctx->mulm (X3, X3, A, ctx);
1334
1335   /* Y_3 = A · G · (D - aC) */
1336   if (ctx->dialect == ECC_DIALECT_ED25519)
1337     {
1338       ctx->addm (Y3, D, C, ctx);
1339     }
1340   else
1341     {
1342       ctx->mulm (Y3, ctx->a, C, ctx);
1343       ctx->subm (Y3, D, Y3, ctx);
1344     }
1345   ctx->mulm (Y3, Y3, G, ctx);
1346   ctx->mulm (Y3, Y3, A, ctx);
1347
1348   /* Z_3 = F · G */
1349   ctx->mulm (Z3, F, G, ctx);
1350
1351
1352 #undef X1
1353 #undef Y1
1354 #undef Z1
1355 #undef X2
1356 #undef Y2
1357 #undef Z2
1358 #undef X3
1359 #undef Y3
1360 #undef Z3
1361 #undef A
1362 #undef B
1363 #undef C
1364 #undef D
1365 #undef E
1366 #undef F
1367 #undef G
1368 #undef tmp
1369 }
1370
1371
1372 /* Compute a step of Montgomery Ladder (only use X and Z in the point).
1373    Inputs:  P1, P2, and x-coordinate of DIF = P1 - P1.
1374    Outputs: PRD = 2 * P1 and  SUM = P1 + P2. */
1375 static void
1376 montgomery_ladder (mpi_point_t prd, mpi_point_t sum,
1377                    mpi_point_t p1, mpi_point_t p2, gcry_mpi_t dif_x,
1378                    mpi_ec_t ctx)
1379 {
1380   ctx->addm (sum->x, p2->x, p2->z, ctx);
1381   ctx->subm (p2->z, p2->x, p2->z, ctx);
1382   ctx->addm (prd->x, p1->x, p1->z, ctx);
1383   ctx->subm (p1->z, p1->x, p1->z, ctx);
1384   ctx->mulm (p2->x, p1->z, sum->x, ctx);
1385   ctx->mulm (p2->z, prd->x, p2->z, ctx);
1386   ctx->pow2 (p1->x, prd->x, ctx);
1387   ctx->pow2 (p1->z, p1->z, ctx);
1388   ctx->addm (sum->x, p2->x, p2->z, ctx);
1389   ctx->subm (p2->z, p2->x, p2->z, ctx);
1390   ctx->mulm (prd->x, p1->x, p1->z, ctx);
1391   ctx->subm (p1->z, p1->x, p1->z, ctx);
1392   ctx->pow2 (sum->x, sum->x, ctx);
1393   ctx->pow2 (sum->z, p2->z, ctx);
1394   ctx->mulm (prd->z, p1->z, ctx->a, ctx); /* CTX->A: (a-2)/4 */
1395   ctx->mulm (sum->z, sum->z, dif_x, ctx);
1396   ctx->addm (prd->z, p1->x, prd->z, ctx);
1397   ctx->mulm (prd->z, prd->z, p1->z, ctx);
1398 }
1399
1400
1401 /* RESULT = P1 + P2 */
1402 void
1403 _gcry_mpi_ec_add_points (mpi_point_t result,
1404                          mpi_point_t p1, mpi_point_t p2,
1405                          mpi_ec_t ctx)
1406 {
1407   switch (ctx->model)
1408     {
1409     case MPI_EC_WEIERSTRASS:
1410       add_points_weierstrass (result, p1, p2, ctx);
1411       break;
1412     case MPI_EC_MONTGOMERY:
1413       add_points_montgomery (result, p1, p2, ctx);
1414       break;
1415     case MPI_EC_EDWARDS:
1416       add_points_edwards (result, p1, p2, ctx);
1417       break;
1418     }
1419 }
1420
1421
1422 /* RESULT = P1 - P2  (Weierstrass version).*/
1423 static void
1424 sub_points_weierstrass (mpi_point_t result,
1425                         mpi_point_t p1, mpi_point_t p2,
1426                         mpi_ec_t ctx)
1427 {
1428   (void)result;
1429   (void)p1;
1430   (void)p2;
1431   (void)ctx;
1432   log_fatal ("%s: %s not yet supported\n",
1433              "_gcry_mpi_ec_sub_points", "Weierstrass");
1434 }
1435
1436
1437 /* RESULT = P1 - P2  (Montgomery version).*/
1438 static void
1439 sub_points_montgomery (mpi_point_t result,
1440                        mpi_point_t p1, mpi_point_t p2,
1441                        mpi_ec_t ctx)
1442 {
1443   (void)result;
1444   (void)p1;
1445   (void)p2;
1446   (void)ctx;
1447   log_fatal ("%s: %s not yet supported\n",
1448              "_gcry_mpi_ec_sub_points", "Montgomery");
1449 }
1450
1451
1452 /* RESULT = P1 - P2  (Twisted Edwards version).*/
1453 static void
1454 sub_points_edwards (mpi_point_t result,
1455                     mpi_point_t p1, mpi_point_t p2,
1456                     mpi_ec_t ctx)
1457 {
1458   mpi_point_t p2i = _gcry_mpi_point_new (0);
1459   point_set (p2i, p2);
1460   ctx->subm (p2i->x, ctx->p, p2i->x, ctx);
1461   add_points_edwards (result, p1, p2i, ctx);
1462   _gcry_mpi_point_release (p2i);
1463 }
1464
1465
1466 /* RESULT = P1 - P2 */
1467 void
1468 _gcry_mpi_ec_sub_points (mpi_point_t result,
1469                          mpi_point_t p1, mpi_point_t p2,
1470                          mpi_ec_t ctx)
1471 {
1472   switch (ctx->model)
1473     {
1474     case MPI_EC_WEIERSTRASS:
1475       sub_points_weierstrass (result, p1, p2, ctx);
1476       break;
1477     case MPI_EC_MONTGOMERY:
1478       sub_points_montgomery (result, p1, p2, ctx);
1479       break;
1480     case MPI_EC_EDWARDS:
1481       sub_points_edwards (result, p1, p2, ctx);
1482       break;
1483     }
1484 }
1485
1486
1487 /* Scalar point multiplication - the main function for ECC.  If takes
1488    an integer SCALAR and a POINT as well as the usual context CTX.
1489    RESULT will be set to the resulting point. */
1490 void
1491 _gcry_mpi_ec_mul_point (mpi_point_t result,
1492                         gcry_mpi_t scalar, mpi_point_t point,
1493                         mpi_ec_t ctx)
1494 {
1495   gcry_mpi_t x1, y1, z1, k, h, yy;
1496   unsigned int i, loops;
1497   mpi_point_struct p1, p2, p1inv;
1498
1499   if (ctx->model == MPI_EC_EDWARDS
1500       || (ctx->model == MPI_EC_WEIERSTRASS
1501           && mpi_is_secure (scalar)))
1502     {
1503       /* Simple left to right binary method.  Algorithm 3.27 from
1504        * {author={Hankerson, Darrel and Menezes, Alfred J. and Vanstone, Scott},
1505        *  title = {Guide to Elliptic Curve Cryptography},
1506        *  year = {2003}, isbn = {038795273X},
1507        *  url = {http://www.cacr.math.uwaterloo.ca/ecc/},
1508        *  publisher = {Springer-Verlag New York, Inc.}} */
1509       unsigned int nbits;
1510       int j;
1511
1512       nbits = mpi_get_nbits (scalar);
1513       if (ctx->model == MPI_EC_WEIERSTRASS)
1514         {
1515           mpi_set_ui (result->x, 1);
1516           mpi_set_ui (result->y, 1);
1517           mpi_set_ui (result->z, 0);
1518         }
1519       else
1520         {
1521           mpi_set_ui (result->x, 0);
1522           mpi_set_ui (result->y, 1);
1523           mpi_set_ui (result->z, 1);
1524           point_resize (point, ctx);
1525         }
1526
1527       if (mpi_is_secure (scalar))
1528         {
1529           /* If SCALAR is in secure memory we assume that it is the
1530              secret key we use constant time operation.  */
1531           mpi_point_struct tmppnt;
1532
1533           point_init (&tmppnt);
1534           point_resize (result, ctx);
1535           point_resize (&tmppnt, ctx);
1536           for (j=nbits-1; j >= 0; j--)
1537             {
1538               _gcry_mpi_ec_dup_point (result, result, ctx);
1539               _gcry_mpi_ec_add_points (&tmppnt, result, point, ctx);
1540               point_swap_cond (result, &tmppnt, mpi_test_bit (scalar, j), ctx);
1541             }
1542           point_free (&tmppnt);
1543         }
1544       else
1545         {
1546           if (ctx->model == MPI_EC_EDWARDS)
1547             {
1548               point_resize (result, ctx);
1549               point_resize (point, ctx);
1550             }
1551
1552           for (j=nbits-1; j >= 0; j--)
1553             {
1554               _gcry_mpi_ec_dup_point (result, result, ctx);
1555               if (mpi_test_bit (scalar, j))
1556                 _gcry_mpi_ec_add_points (result, result, point, ctx);
1557             }
1558         }
1559       return;
1560     }
1561   else if (ctx->model == MPI_EC_MONTGOMERY)
1562     {
1563       unsigned int nbits;
1564       int j;
1565       mpi_point_struct p1_, p2_;
1566       mpi_point_t q1, q2, prd, sum;
1567       unsigned long sw;
1568       mpi_size_t rsize;
1569
1570       /* Compute scalar point multiplication with Montgomery Ladder.
1571          Note that we don't use Y-coordinate in the points at all.
1572          RESULT->Y will be filled by zero.  */
1573
1574       nbits = mpi_get_nbits (scalar);
1575       point_init (&p1);
1576       point_init (&p2);
1577       point_init (&p1_);
1578       point_init (&p2_);
1579       mpi_set_ui (p1.x, 1);
1580       mpi_free (p2.x);
1581       p2.x  = mpi_copy (point->x);
1582       mpi_set_ui (p2.z, 1);
1583
1584       point_resize (&p1, ctx);
1585       point_resize (&p2, ctx);
1586       point_resize (&p1_, ctx);
1587       point_resize (&p2_, ctx);
1588
1589       mpi_resize (point->x, ctx->p->nlimbs);
1590       point->x->nlimbs = ctx->p->nlimbs;
1591
1592       q1 = &p1;
1593       q2 = &p2;
1594       prd = &p1_;
1595       sum = &p2_;
1596
1597       for (j=nbits-1; j >= 0; j--)
1598         {
1599           mpi_point_t t;
1600
1601           sw = mpi_test_bit (scalar, j);
1602           point_swap_cond (q1, q2, sw, ctx);
1603           montgomery_ladder (prd, sum, q1, q2, point->x, ctx);
1604           point_swap_cond (prd, sum, sw, ctx);
1605           t = q1;  q1 = prd;  prd = t;
1606           t = q2;  q2 = sum;  sum = t;
1607         }
1608
1609       mpi_clear (result->y);
1610       sw = (nbits & 1);
1611       point_swap_cond (&p1, &p1_, sw, ctx);
1612
1613       rsize = p1.z->nlimbs;
1614       MPN_NORMALIZE (p1.z->d, rsize);
1615       if (rsize == 0)
1616         {
1617           mpi_set_ui (result->x, 1);
1618           mpi_set_ui (result->z, 0);
1619         }
1620       else
1621         {
1622           z1 = mpi_new (0);
1623           ec_invm (z1, p1.z, ctx);
1624           ec_mulm (result->x, p1.x, z1, ctx);
1625           mpi_set_ui (result->z, 1);
1626           mpi_free (z1);
1627         }
1628
1629       point_free (&p1);
1630       point_free (&p2);
1631       point_free (&p1_);
1632       point_free (&p2_);
1633       return;
1634     }
1635
1636   x1 = mpi_alloc_like (ctx->p);
1637   y1 = mpi_alloc_like (ctx->p);
1638   h  = mpi_alloc_like (ctx->p);
1639   k  = mpi_copy (scalar);
1640   yy = mpi_copy (point->y);
1641
1642   if ( mpi_has_sign (k) )
1643     {
1644       k->sign = 0;
1645       ec_invm (yy, yy, ctx);
1646     }
1647
1648   if (!mpi_cmp_ui (point->z, 1))
1649     {
1650       mpi_set (x1, point->x);
1651       mpi_set (y1, yy);
1652     }
1653   else
1654     {
1655       gcry_mpi_t z2, z3;
1656
1657       z2 = mpi_alloc_like (ctx->p);
1658       z3 = mpi_alloc_like (ctx->p);
1659       ec_mulm (z2, point->z, point->z, ctx);
1660       ec_mulm (z3, point->z, z2, ctx);
1661       ec_invm (z2, z2, ctx);
1662       ec_mulm (x1, point->x, z2, ctx);
1663       ec_invm (z3, z3, ctx);
1664       ec_mulm (y1, yy, z3, ctx);
1665       mpi_free (z2);
1666       mpi_free (z3);
1667     }
1668   z1 = mpi_copy (mpi_const (MPI_C_ONE));
1669
1670   mpi_mul (h, k, mpi_const (MPI_C_THREE)); /* h = 3k */
1671   loops = mpi_get_nbits (h);
1672   if (loops < 2)
1673     {
1674       /* If SCALAR is zero, the above mpi_mul sets H to zero and thus
1675          LOOPs will be zero.  To avoid an underflow of I in the main
1676          loop we set LOOP to 2 and the result to (0,0,0).  */
1677       loops = 2;
1678       mpi_clear (result->x);
1679       mpi_clear (result->y);
1680       mpi_clear (result->z);
1681     }
1682   else
1683     {
1684       mpi_set (result->x, point->x);
1685       mpi_set (result->y, yy);
1686       mpi_set (result->z, point->z);
1687     }
1688   mpi_free (yy); yy = NULL;
1689
1690   p1.x = x1; x1 = NULL;
1691   p1.y = y1; y1 = NULL;
1692   p1.z = z1; z1 = NULL;
1693   point_init (&p2);
1694   point_init (&p1inv);
1695
1696   /* Invert point: y = p - y mod p  */
1697   point_set (&p1inv, &p1);
1698   ec_subm (p1inv.y, ctx->p, p1inv.y, ctx);
1699
1700   for (i=loops-2; i > 0; i--)
1701     {
1702       _gcry_mpi_ec_dup_point (result, result, ctx);
1703       if (mpi_test_bit (h, i) == 1 && mpi_test_bit (k, i) == 0)
1704         {
1705           point_set (&p2, result);
1706           _gcry_mpi_ec_add_points (result, &p2, &p1, ctx);
1707         }
1708       if (mpi_test_bit (h, i) == 0 && mpi_test_bit (k, i) == 1)
1709         {
1710           point_set (&p2, result);
1711           _gcry_mpi_ec_add_points (result, &p2, &p1inv, ctx);
1712         }
1713     }
1714
1715   point_free (&p1);
1716   point_free (&p2);
1717   point_free (&p1inv);
1718   mpi_free (h);
1719   mpi_free (k);
1720 }
1721
1722
1723 /* Return true if POINT is on the curve described by CTX.  */
1724 int
1725 _gcry_mpi_ec_curve_point (gcry_mpi_point_t point, mpi_ec_t ctx)
1726 {
1727   int res = 0;
1728   gcry_mpi_t x, y, w;
1729
1730   x = mpi_new (0);
1731   y = mpi_new (0);
1732   w = mpi_new (0);
1733
1734   /* Check that the point is in range.  This needs to be done here and
1735    * not after conversion to affine coordinates.  */
1736   if (mpi_cmpabs (point->x, ctx->p) >= 0)
1737     goto leave;
1738   if (mpi_cmpabs (point->y, ctx->p) >= 0)
1739     goto leave;
1740   if (mpi_cmpabs (point->z, ctx->p) >= 0)
1741     goto leave;
1742
1743   switch (ctx->model)
1744     {
1745     case MPI_EC_WEIERSTRASS:
1746       {
1747         gcry_mpi_t xxx;
1748
1749         if (_gcry_mpi_ec_get_affine (x, y, point, ctx))
1750           goto leave;
1751
1752         xxx = mpi_new (0);
1753
1754         /* y^2 == x^3 + a·x + b */
1755         ec_pow2 (y, y, ctx);
1756
1757         ec_pow3 (xxx, x, ctx);
1758         ec_mulm (w, ctx->a, x, ctx);
1759         ec_addm (w, w, ctx->b, ctx);
1760         ec_addm (w, w, xxx, ctx);
1761
1762         if (!mpi_cmp (y, w))
1763           res = 1;
1764
1765         _gcry_mpi_release (xxx);
1766       }
1767       break;
1768     case MPI_EC_MONTGOMERY:
1769       {
1770 #define xx y
1771         /* With Montgomery curve, only X-coordinate is valid.  */
1772         if (_gcry_mpi_ec_get_affine (x, NULL, point, ctx))
1773           goto leave;
1774
1775         /* The equation is: b * y^2 == x^3 + a · x^2 + x */
1776         /* We check if right hand is quadratic residue or not by
1777            Euler's criterion.  */
1778         /* CTX->A has (a-2)/4 and CTX->B has b^-1 */
1779         ec_mulm (w, ctx->a, mpi_const (MPI_C_FOUR), ctx);
1780         ec_addm (w, w, mpi_const (MPI_C_TWO), ctx);
1781         ec_mulm (w, w, x, ctx);
1782         ec_pow2 (xx, x, ctx);
1783         ec_addm (w, w, xx, ctx);
1784         ec_addm (w, w, mpi_const (MPI_C_ONE), ctx);
1785         ec_mulm (w, w, x, ctx);
1786         ec_mulm (w, w, ctx->b, ctx);
1787 #undef xx
1788         /* Compute Euler's criterion: w^(p-1)/2 */
1789 #define p_minus1 y
1790         ec_subm (p_minus1, ctx->p, mpi_const (MPI_C_ONE), ctx);
1791         mpi_rshift (p_minus1, p_minus1, 1);
1792         ec_powm (w, w, p_minus1, ctx);
1793
1794         res = !mpi_cmp_ui (w, 1);
1795 #undef p_minus1
1796       }
1797       break;
1798     case MPI_EC_EDWARDS:
1799       {
1800         if (_gcry_mpi_ec_get_affine (x, y, point, ctx))
1801           goto leave;
1802
1803         mpi_resize (w, ctx->p->nlimbs);
1804         w->nlimbs = ctx->p->nlimbs;
1805
1806         /* a · x^2 + y^2 - 1 - b · x^2 · y^2 == 0 */
1807         ctx->pow2 (x, x, ctx);
1808         ctx->pow2 (y, y, ctx);
1809         if (ctx->dialect == ECC_DIALECT_ED25519)
1810           ctx->subm (w, ctx->p, x, ctx);
1811         else
1812           ctx->mulm (w, ctx->a, x, ctx);
1813         ctx->addm (w, w, y, ctx);
1814         ctx->mulm (x, x, y, ctx);
1815         ctx->mulm (x, x, ctx->b, ctx);
1816         ctx->subm (w, w, x, ctx);
1817         if (!mpi_cmp_ui (w, 1))
1818           res = 1;
1819       }
1820       break;
1821     }
1822
1823  leave:
1824   _gcry_mpi_release (w);
1825   _gcry_mpi_release (x);
1826   _gcry_mpi_release (y);
1827
1828   return res;
1829 }
1830
1831
1832 int
1833 _gcry_mpi_ec_bad_point (gcry_mpi_point_t point, mpi_ec_t ctx)
1834 {
1835   int i;
1836   gcry_mpi_t x_bad;
1837
1838   for (i = 0; (x_bad = ctx->t.scratch[i]); i++)
1839     if (!mpi_cmp (point->x, x_bad))
1840       return 1;
1841
1842   return 0;
1843 }