log: Try to print s-expressions in a more compact format.
[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
32
33 #define point_init(a)  _gcry_mpi_point_init ((a))
34 #define point_free(a)  _gcry_mpi_point_free_parts ((a))
35
36
37 /* Print a point using the log fucntions.  If CTX is not NULL affine
38    coordinates will be printed.  */
39 void
40 _gcry_mpi_point_log (const char *name, mpi_point_t point, mpi_ec_t ctx)
41 {
42   gcry_mpi_t x, y;
43   char buf[100];
44
45   if (!point)
46     {
47       snprintf (buf, sizeof buf - 1, "%s.*", name);
48       log_mpidump (buf, NULL);
49       return;
50     }
51   snprintf (buf, sizeof buf - 1, "%s.X", name);
52
53   if (ctx)
54     {
55       x = gcry_mpi_new (0);
56       y = gcry_mpi_new (0);
57     }
58   if (!ctx || _gcry_mpi_ec_get_affine (x, y, point, ctx))
59     {
60       log_mpidump (buf, point->x);
61       buf[strlen(buf)-1] = 'Y';
62       log_mpidump (buf, point->y);
63       buf[strlen(buf)-1] = 'Z';
64       log_mpidump (buf, point->z);
65     }
66   else
67     {
68       buf[strlen(buf)-1] = 'x';
69       log_mpidump (buf, x);
70       buf[strlen(buf)-1] = 'y';
71       log_mpidump (buf, y);
72
73     }
74   if (ctx)
75     {
76       gcry_mpi_release (x);
77       gcry_mpi_release (y);
78     }
79 }
80
81
82 /* Create a new point option.  NBITS gives the size in bits of one
83    coordinate; it is only used to pre-allocate some resources and
84    might also be passed as 0 to use a default value.  */
85 mpi_point_t
86 gcry_mpi_point_new (unsigned int nbits)
87 {
88   mpi_point_t p;
89
90   (void)nbits;  /* Currently not used.  */
91
92   p = gcry_xmalloc (sizeof *p);
93   _gcry_mpi_point_init (p);
94   return p;
95 }
96
97
98 /* Release the point object P.  P may be NULL. */
99 void
100 gcry_mpi_point_release (mpi_point_t p)
101 {
102   if (p)
103     {
104       _gcry_mpi_point_free_parts (p);
105       gcry_free (p);
106     }
107 }
108
109
110 /* Initialize the fields of a point object.  gcry_mpi_point_free_parts
111    may be used to release the fields.  */
112 void
113 _gcry_mpi_point_init (mpi_point_t p)
114 {
115   p->x = mpi_new (0);
116   p->y = mpi_new (0);
117   p->z = mpi_new (0);
118 }
119
120
121 /* Release the parts of a point object. */
122 void
123 _gcry_mpi_point_free_parts (mpi_point_t p)
124 {
125   mpi_free (p->x); p->x = NULL;
126   mpi_free (p->y); p->y = NULL;
127   mpi_free (p->z); p->z = NULL;
128 }
129
130
131 /* Set the value from S into D.  */
132 static void
133 point_set (mpi_point_t d, mpi_point_t s)
134 {
135   mpi_set (d->x, s->x);
136   mpi_set (d->y, s->y);
137   mpi_set (d->z, s->z);
138 }
139
140
141 /* Set the projective coordinates from POINT into X, Y, and Z.  If a
142    coordinate is not required, X, Y, or Z may be passed as NULL.  */
143 void
144 gcry_mpi_point_get (gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z,
145                     mpi_point_t point)
146 {
147   if (x)
148     mpi_set (x, point->x);
149   if (y)
150     mpi_set (y, point->y);
151   if (z)
152     mpi_set (z, point->z);
153 }
154
155
156 /* Set the projective coordinates from POINT into X, Y, and Z and
157    release POINT.  If a coordinate is not required, X, Y, or Z may be
158    passed as NULL.  */
159 void
160 gcry_mpi_point_snatch_get (gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z,
161                            mpi_point_t point)
162 {
163   mpi_snatch (x, point->x);
164   mpi_snatch (y, point->y);
165   mpi_snatch (z, point->z);
166   gcry_free (point);
167 }
168
169
170 /* Set the projective coordinates from X, Y, and Z into POINT.  If a
171    coordinate is given as NULL, the value 0 is stored into point.  If
172    POINT is given as NULL a new point object is allocated.  Returns
173    POINT or the newly allocated point object. */
174 mpi_point_t
175 gcry_mpi_point_set (mpi_point_t point,
176                     gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z)
177 {
178   if (!point)
179     point = gcry_mpi_point_new (0);
180
181   if (x)
182     mpi_set (point->x, x);
183   else
184     mpi_clear (point->x);
185   if (y)
186     mpi_set (point->y, y);
187   else
188     mpi_clear (point->y);
189   if (z)
190     mpi_set (point->z, z);
191   else
192     mpi_clear (point->z);
193
194   return point;
195 }
196
197
198 /* Set the projective coordinates from X, Y, and Z into POINT.  If a
199    coordinate is given as NULL, the value 0 is stored into point.  If
200    POINT is given as NULL a new point object is allocated.  The
201    coordinates X, Y, and Z are released.  Returns POINT or the newly
202    allocated point object. */
203 mpi_point_t
204 gcry_mpi_point_snatch_set (mpi_point_t point,
205                            gcry_mpi_t x, gcry_mpi_t y, gcry_mpi_t z)
206 {
207   if (!point)
208     point = gcry_mpi_point_new (0);
209
210   if (x)
211     mpi_snatch (point->x, x);
212   else
213     mpi_clear (point->x);
214   if (y)
215     mpi_snatch (point->y, y);
216   else
217     mpi_clear (point->y);
218   if (z)
219     mpi_snatch (point->z, z);
220   else
221     mpi_clear (point->z);
222
223   return point;
224 }
225
226
227 static void
228 ec_addm (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
229 {
230   mpi_addm (w, u, v, ctx->p);
231 }
232
233 static void
234 ec_subm (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
235 {
236   mpi_subm (w, u, v, ctx->p);
237 }
238
239 static void
240 ec_mulm (gcry_mpi_t w, gcry_mpi_t u, gcry_mpi_t v, mpi_ec_t ctx)
241 {
242 #if 0
243   /* NOTE: This code works only for limb sizes of 32 bit.  */
244   mpi_limb_t *wp, *sp;
245
246   if (ctx->nist_nbits == 192)
247     {
248       mpi_mul (w, u, v);
249       mpi_resize (w, 12);
250       wp = w->d;
251
252       sp = ctx->s[0]->d;
253       sp[0*2+0] = wp[0*2+0];
254       sp[0*2+1] = wp[0*2+1];
255       sp[1*2+0] = wp[1*2+0];
256       sp[1*2+1] = wp[1*2+1];
257       sp[2*2+0] = wp[2*2+0];
258       sp[2*2+1] = wp[2*2+1];
259
260       sp = ctx->s[1]->d;
261       sp[0*2+0] = wp[3*2+0];
262       sp[0*2+1] = wp[3*2+1];
263       sp[1*2+0] = wp[3*2+0];
264       sp[1*2+1] = wp[3*2+1];
265       sp[2*2+0] = 0;
266       sp[2*2+1] = 0;
267
268       sp = ctx->s[2]->d;
269       sp[0*2+0] = 0;
270       sp[0*2+1] = 0;
271       sp[1*2+0] = wp[4*2+0];
272       sp[1*2+1] = wp[4*2+1];
273       sp[2*2+0] = wp[4*2+0];
274       sp[2*2+1] = wp[4*2+1];
275
276       sp = ctx->s[3]->d;
277       sp[0*2+0] = wp[5*2+0];
278       sp[0*2+1] = wp[5*2+1];
279       sp[1*2+0] = wp[5*2+0];
280       sp[1*2+1] = wp[5*2+1];
281       sp[2*2+0] = wp[5*2+0];
282       sp[2*2+1] = wp[5*2+1];
283
284       ctx->s[0]->nlimbs = 6;
285       ctx->s[1]->nlimbs = 6;
286       ctx->s[2]->nlimbs = 6;
287       ctx->s[3]->nlimbs = 6;
288
289       mpi_add (ctx->c, ctx->s[0], ctx->s[1]);
290       mpi_add (ctx->c, ctx->c, ctx->s[2]);
291       mpi_add (ctx->c, ctx->c, ctx->s[3]);
292
293       while ( mpi_cmp (ctx->c, ctx->p ) >= 0 )
294         mpi_sub ( ctx->c, ctx->c, ctx->p );
295       mpi_set (w, ctx->c);
296     }
297   else if (ctx->nist_nbits == 384)
298     {
299       int i;
300       mpi_mul (w, u, v);
301       mpi_resize (w, 24);
302       wp = w->d;
303
304 #define NEXT(a) do { ctx->s[(a)]->nlimbs = 12; \
305                      sp = ctx->s[(a)]->d; \
306                      i = 0; } while (0)
307 #define X(a) do { sp[i++] = wp[(a)];} while (0)
308 #define X0(a) do { sp[i++] = 0; } while (0)
309       NEXT(0);
310       X(0);X(1);X(2);X(3);X(4);X(5);X(6);X(7);X(8);X(9);X(10);X(11);
311       NEXT(1);
312       X0();X0();X0();X0();X(21);X(22);X(23);X0();X0();X0();X0();X0();
313       NEXT(2);
314       X(12);X(13);X(14);X(15);X(16);X(17);X(18);X(19);X(20);X(21);X(22);X(23);
315       NEXT(3);
316       X(21);X(22);X(23);X(12);X(13);X(14);X(15);X(16);X(17);X(18);X(19);X(20);
317       NEXT(4);
318       X0();X(23);X0();X(20);X(12);X(13);X(14);X(15);X(16);X(17);X(18);X(19);
319       NEXT(5);
320       X0();X0();X0();X0();X(20);X(21);X(22);X(23);X0();X0();X0();X0();
321       NEXT(6);
322       X(20);X0();X0();X(21);X(22);X(23);X0();X0();X0();X0();X0();X0();
323       NEXT(7);
324       X(23);X(12);X(13);X(14);X(15);X(16);X(17);X(18);X(19);X(20);X(21);X(22);
325       NEXT(8);
326       X0();X(20);X(21);X(22);X(23);X0();X0();X0();X0();X0();X0();X0();
327       NEXT(9);
328       X0();X0();X0();X(23);X(23);X0();X0();X0();X0();X0();X0();X0();
329 #undef X0
330 #undef X
331 #undef NEXT
332       mpi_add (ctx->c, ctx->s[0], ctx->s[1]);
333       mpi_add (ctx->c, ctx->c, ctx->s[1]);
334       mpi_add (ctx->c, ctx->c, ctx->s[2]);
335       mpi_add (ctx->c, ctx->c, ctx->s[3]);
336       mpi_add (ctx->c, ctx->c, ctx->s[4]);
337       mpi_add (ctx->c, ctx->c, ctx->s[5]);
338       mpi_add (ctx->c, ctx->c, ctx->s[6]);
339       mpi_sub (ctx->c, ctx->c, ctx->s[7]);
340       mpi_sub (ctx->c, ctx->c, ctx->s[8]);
341       mpi_sub (ctx->c, ctx->c, ctx->s[9]);
342
343       while ( mpi_cmp (ctx->c, ctx->p ) >= 0 )
344         mpi_sub ( ctx->c, ctx->c, ctx->p );
345       while ( ctx->c->sign )
346         mpi_add ( ctx->c, ctx->c, ctx->p );
347       mpi_set (w, ctx->c);
348     }
349   else
350 #endif /*0*/
351     mpi_mulm (w, u, v, ctx->p);
352 }
353
354 static void
355 ec_powm (gcry_mpi_t w, const gcry_mpi_t b, const gcry_mpi_t e,
356          mpi_ec_t ctx)
357 {
358   mpi_powm (w, b, e, ctx->p);
359   /* _gcry_mpi_abs (w); */
360 }
361
362
363 /* Shortcut for
364      ec_powm (B, B, mpi_const (MPI_C_TWO), ctx);
365    for easier optimization.  */
366 static void
367 ec_pow2 (gcry_mpi_t w, const gcry_mpi_t b, mpi_ec_t ctx)
368 {
369   /* Using mpi_mul is slightly faster (at least on amd64).  */
370   /* mpi_powm (w, b, mpi_const (MPI_C_TWO), ctx->p); */
371   mpi_mulm (w, b, b, ctx->p);
372 }
373
374
375 static void
376 ec_invm (gcry_mpi_t x, gcry_mpi_t a, mpi_ec_t ctx)
377 {
378   if (!mpi_invm (x, a, ctx->p))
379     {
380       log_error ("ec_invm: inverse does not exist:\n");
381       log_mpidump ("  a", a);
382       log_mpidump ("  p", ctx->p);
383     }
384 }
385
386
387 /* Force recomputation of all helper variables.  */
388 void
389 _gcry_mpi_ec_get_reset (mpi_ec_t ec)
390 {
391   ec->t.valid.a_is_pminus3 = 0;
392   ec->t.valid.two_inv_p = 0;
393 }
394
395
396 /* Accessor for helper variable.  */
397 static int
398 ec_get_a_is_pminus3 (mpi_ec_t ec)
399 {
400   gcry_mpi_t tmp;
401
402   if (!ec->t.valid.a_is_pminus3)
403     {
404       ec->t.valid.a_is_pminus3 = 1;
405       tmp = mpi_alloc_like (ec->p);
406       mpi_sub_ui (tmp, ec->p, 3);
407       ec->t.a_is_pminus3 = !mpi_cmp (ec->a, tmp);
408       mpi_free (tmp);
409     }
410
411   return ec->t.a_is_pminus3;
412 }
413
414
415 /* Accessor for helper variable.  */
416 static gcry_mpi_t
417 ec_get_two_inv_p (mpi_ec_t ec)
418 {
419   if (!ec->t.valid.two_inv_p)
420     {
421       ec->t.valid.two_inv_p = 1;
422       if (!ec->t.two_inv_p)
423         ec->t.two_inv_p = mpi_alloc (0);
424       ec_invm (ec->t.two_inv_p, mpi_const (MPI_C_TWO), ec);
425     }
426   return ec->t.two_inv_p;
427 }
428
429
430
431 /* This function initialized a context for elliptic curve based on the
432    field GF(p).  P is the prime specifying this field, A is the first
433    coefficient.  CTX is expected to be zeroized.  */
434 static void
435 ec_p_init (mpi_ec_t ctx, enum gcry_mpi_ec_models model,
436            enum ecc_dialects dialect,
437            gcry_mpi_t p, gcry_mpi_t a, gcry_mpi_t b)
438 {
439   int i;
440
441   /* Fixme: Do we want to check some constraints? e.g.  a < p  */
442
443   ctx->model = model;
444   ctx->dialect = dialect;
445   if (dialect == ECC_DIALECT_ED25519)
446     ctx->nbits = 256;
447   else
448     ctx->nbits = mpi_get_nbits (p);
449   ctx->p = mpi_copy (p);
450   ctx->a = mpi_copy (a);
451   if (b && model == MPI_EC_TWISTEDEDWARDS)
452     ctx->b = mpi_copy (b);
453
454   _gcry_mpi_ec_get_reset (ctx);
455
456   /* Allocate scratch variables.  */
457   for (i=0; i< DIM(ctx->t.scratch); i++)
458     ctx->t.scratch[i] = mpi_alloc_like (ctx->p);
459
460   /* Prepare for fast reduction.  */
461   /* FIXME: need a test for NIST values.  However it does not gain us
462      any real advantage, for 384 bits it is actually slower than using
463      mpi_mulm.  */
464 /*   ctx->nist_nbits = mpi_get_nbits (ctx->p); */
465 /*   if (ctx->nist_nbits == 192) */
466 /*     { */
467 /*       for (i=0; i < 4; i++) */
468 /*         ctx->s[i] = mpi_new (192); */
469 /*       ctx->c    = mpi_new (192*2); */
470 /*     } */
471 /*   else if (ctx->nist_nbits == 384) */
472 /*     { */
473 /*       for (i=0; i < 10; i++) */
474 /*         ctx->s[i] = mpi_new (384); */
475 /*       ctx->c    = mpi_new (384*2); */
476 /*     } */
477 }
478
479
480 static void
481 ec_deinit (void *opaque)
482 {
483   mpi_ec_t ctx = opaque;
484   int i;
485
486   /* Domain parameter.  */
487   mpi_free (ctx->p);
488   mpi_free (ctx->a);
489   mpi_free (ctx->b);
490   gcry_mpi_point_release (ctx->G);
491   mpi_free (ctx->n);
492
493   /* The key.  */
494   gcry_mpi_point_release (ctx->Q);
495   mpi_free (ctx->d);
496
497   /* Private data of ec.c.  */
498   mpi_free (ctx->t.two_inv_p);
499
500   for (i=0; i< DIM(ctx->t.scratch); i++)
501     mpi_free (ctx->t.scratch[i]);
502
503 /*   if (ctx->nist_nbits == 192) */
504 /*     { */
505 /*       for (i=0; i < 4; i++) */
506 /*         mpi_free (ctx->s[i]); */
507 /*       mpi_free (ctx->c); */
508 /*     } */
509 /*   else if (ctx->nist_nbits == 384) */
510 /*     { */
511 /*       for (i=0; i < 10; i++) */
512 /*         mpi_free (ctx->s[i]); */
513 /*       mpi_free (ctx->c); */
514 /*     } */
515 }
516
517
518 /* This function returns a new context for elliptic curve based on the
519    field GF(p).  P is the prime specifying this field, A is the first
520    coefficient, B is the second coefficient, and MODEL is the model
521    for the curve.  This function is only used within Libgcrypt and not
522    part of the public API.
523
524    This context needs to be released using _gcry_mpi_ec_free.  */
525 mpi_ec_t
526 _gcry_mpi_ec_p_internal_new (enum gcry_mpi_ec_models model,
527                              enum ecc_dialects dialect,
528                              gcry_mpi_t p, gcry_mpi_t a, gcry_mpi_t b)
529 {
530   mpi_ec_t ctx;
531
532   ctx = gcry_xcalloc (1, sizeof *ctx);
533   ec_p_init (ctx, model, dialect, p, a, b);
534
535   return ctx;
536 }
537
538
539 /* This is a variant of _gcry_mpi_ec_p_internal_new which returns an
540    public contect and does some error checking on the supplied
541    arguments.  On success the new context is stored at R_CTX and 0 is
542    returned; on error NULL is stored at R_CTX and an error code is
543    returned.
544
545    The context needs to be released using gcry_ctx_release.  */
546 gpg_err_code_t
547 _gcry_mpi_ec_p_new (gcry_ctx_t *r_ctx,
548                     enum gcry_mpi_ec_models model,
549                     enum ecc_dialects dialect,
550                     gcry_mpi_t p, gcry_mpi_t a, gcry_mpi_t b)
551 {
552   gcry_ctx_t ctx;
553   mpi_ec_t ec;
554
555   *r_ctx = NULL;
556   if (!p || !a || !mpi_cmp_ui (a, 0))
557     return GPG_ERR_EINVAL;
558
559   ctx = _gcry_ctx_alloc (CONTEXT_TYPE_EC, sizeof *ec, ec_deinit);
560   if (!ctx)
561     return gpg_err_code_from_syserror ();
562   ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
563   ec_p_init (ec, model, dialect, p, a, b);
564
565   *r_ctx = ctx;
566   return 0;
567 }
568
569
570 void
571 _gcry_mpi_ec_free (mpi_ec_t ctx)
572 {
573   if (ctx)
574     {
575       ec_deinit (ctx);
576       gcry_free (ctx);
577     }
578 }
579
580
581 gcry_mpi_t
582 _gcry_mpi_ec_get_mpi (const char *name, gcry_ctx_t ctx, int copy)
583 {
584   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
585
586   return _gcry_ecc_get_mpi (name, ec, copy);
587 }
588
589
590 gcry_mpi_point_t
591 _gcry_mpi_ec_get_point (const char *name, gcry_ctx_t ctx, int copy)
592 {
593   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
594
595   (void)copy;  /* Not used.  */
596
597   return _gcry_ecc_get_point (name, ec);
598 }
599
600
601 gpg_err_code_t
602 _gcry_mpi_ec_set_mpi (const char *name, gcry_mpi_t newvalue,
603                       gcry_ctx_t ctx)
604 {
605   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
606
607   return _gcry_ecc_set_mpi (name, newvalue, ec);
608 }
609
610
611 gpg_err_code_t
612 _gcry_mpi_ec_set_point (const char *name, gcry_mpi_point_t newvalue,
613                         gcry_ctx_t ctx)
614 {
615   mpi_ec_t ec = _gcry_ctx_get_pointer (ctx, CONTEXT_TYPE_EC);
616
617   return _gcry_ecc_set_point (name, newvalue, ec);
618 }
619
620
621 /* Compute the affine coordinates from the projective coordinates in
622    POINT.  Set them into X and Y.  If one coordinate is not required,
623    X or Y may be passed as NULL.  CTX is the usual context. Returns: 0
624    on success or !0 if POINT is at infinity.  */
625 int
626 _gcry_mpi_ec_get_affine (gcry_mpi_t x, gcry_mpi_t y, mpi_point_t point,
627                          mpi_ec_t ctx)
628 {
629   if (!mpi_cmp_ui (point->z, 0))
630     return -1;
631
632   switch (ctx->model)
633     {
634     case MPI_EC_WEIERSTRASS: /* Using Jacobian coordinates.  */
635       {
636         gcry_mpi_t z1, z2, z3;
637
638         z1 = mpi_new (0);
639         z2 = mpi_new (0);
640         ec_invm (z1, point->z, ctx);  /* z1 = z^(-1) mod p  */
641         ec_mulm (z2, z1, z1, ctx);    /* z2 = z^(-2) mod p  */
642
643         if (x)
644           ec_mulm (x, point->x, z2, ctx);
645
646         if (y)
647           {
648             z3 = mpi_new (0);
649             ec_mulm (z3, z2, z1, ctx);      /* z3 = z^(-3) mod p  */
650             ec_mulm (y, point->y, z3, ctx);
651             mpi_free (z3);
652           }
653
654         mpi_free (z2);
655         mpi_free (z1);
656       }
657       return 0;
658
659     case MPI_EC_MONTGOMERY:
660       {
661         log_fatal ("%s: %s not yet supported\n",
662                    "_gcry_mpi_ec_get_affine", "Montgomery");
663       }
664       return -1;
665
666     case MPI_EC_TWISTEDEDWARDS:
667       {
668         gcry_mpi_t z;
669
670         z = mpi_new (0);
671         ec_invm (z, point->z, ctx);
672
673         if (x)
674           ec_mulm (x, point->x, z, ctx);
675         if (y)
676           ec_mulm (y, point->y, z, ctx);
677
678         gcry_mpi_release (z);
679       }
680       return 0;
681
682     default:
683       return -1;
684     }
685 }
686
687
688 \f
689 /*  RESULT = 2 * POINT  (Weierstrass version). */
690 static void
691 dup_point_weierstrass (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
692 {
693 #define x3 (result->x)
694 #define y3 (result->y)
695 #define z3 (result->z)
696 #define t1 (ctx->t.scratch[0])
697 #define t2 (ctx->t.scratch[1])
698 #define t3 (ctx->t.scratch[2])
699 #define l1 (ctx->t.scratch[3])
700 #define l2 (ctx->t.scratch[4])
701 #define l3 (ctx->t.scratch[5])
702
703   if (!mpi_cmp_ui (point->y, 0) || !mpi_cmp_ui (point->z, 0))
704     {
705       /* P_y == 0 || P_z == 0 => [1:1:0] */
706       mpi_set_ui (x3, 1);
707       mpi_set_ui (y3, 1);
708       mpi_set_ui (z3, 0);
709     }
710   else
711     {
712       if (ec_get_a_is_pminus3 (ctx))  /* Use the faster case.  */
713         {
714           /* L1 = 3(X - Z^2)(X + Z^2) */
715           /*                          T1: used for Z^2. */
716           /*                          T2: used for the right term.  */
717           ec_pow2 (t1, point->z, ctx);
718           ec_subm (l1, point->x, t1, ctx);
719           ec_mulm (l1, l1, mpi_const (MPI_C_THREE), ctx);
720           ec_addm (t2, point->x, t1, ctx);
721           ec_mulm (l1, l1, t2, ctx);
722         }
723       else /* Standard case. */
724         {
725           /* L1 = 3X^2 + aZ^4 */
726           /*                          T1: used for aZ^4. */
727           ec_pow2 (l1, point->x, ctx);
728           ec_mulm (l1, l1, mpi_const (MPI_C_THREE), ctx);
729           ec_powm (t1, point->z, mpi_const (MPI_C_FOUR), ctx);
730           ec_mulm (t1, t1, ctx->a, ctx);
731           ec_addm (l1, l1, t1, ctx);
732         }
733       /* Z3 = 2YZ */
734       ec_mulm (z3, point->y, point->z, ctx);
735       ec_mulm (z3, z3, mpi_const (MPI_C_TWO), ctx);
736
737       /* L2 = 4XY^2 */
738       /*                              T2: used for Y2; required later. */
739       ec_pow2 (t2, point->y, ctx);
740       ec_mulm (l2, t2, point->x, ctx);
741       ec_mulm (l2, l2, mpi_const (MPI_C_FOUR), ctx);
742
743       /* X3 = L1^2 - 2L2 */
744       /*                              T1: used for L2^2. */
745       ec_pow2 (x3, l1, ctx);
746       ec_mulm (t1, l2, mpi_const (MPI_C_TWO), ctx);
747       ec_subm (x3, x3, t1, ctx);
748
749       /* L3 = 8Y^4 */
750       /*                              T2: taken from above. */
751       ec_pow2 (t2, t2, ctx);
752       ec_mulm (l3, t2, mpi_const (MPI_C_EIGHT), ctx);
753
754       /* Y3 = L1(L2 - X3) - L3 */
755       ec_subm (y3, l2, x3, ctx);
756       ec_mulm (y3, y3, l1, ctx);
757       ec_subm (y3, y3, l3, ctx);
758     }
759
760 #undef x3
761 #undef y3
762 #undef z3
763 #undef t1
764 #undef t2
765 #undef t3
766 #undef l1
767 #undef l2
768 #undef l3
769 }
770
771
772 /*  RESULT = 2 * POINT  (Montgomery version). */
773 static void
774 dup_point_montgomery (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
775 {
776   (void)result;
777   (void)point;
778   (void)ctx;
779   log_fatal ("%s: %s not yet supported\n",
780              "_gcry_mpi_ec_dup_point", "Montgomery");
781 }
782
783
784 /*  RESULT = 2 * POINT  (Twisted Edwards version). */
785 static void
786 dup_point_twistededwards (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
787 {
788 #define X1 (point->x)
789 #define Y1 (point->y)
790 #define Z1 (point->z)
791 #define X3 (result->x)
792 #define Y3 (result->y)
793 #define Z3 (result->z)
794 #define B (ctx->t.scratch[0])
795 #define C (ctx->t.scratch[1])
796 #define D (ctx->t.scratch[2])
797 #define E (ctx->t.scratch[3])
798 #define F (ctx->t.scratch[4])
799 #define H (ctx->t.scratch[5])
800 #define J (ctx->t.scratch[6])
801
802   /* Compute: (X_3 : Y_3 : Z_3) = 2( X_1 : Y_1 : Z_1 ) */
803
804   /* B = (X_1 + Y_1)^2  */
805   ec_addm (B, X1, Y1, ctx);
806   ec_pow2 (B, B, ctx);
807
808   /* C = X_1^2 */
809   /* D = Y_1^2 */
810   ec_pow2 (C, X1, ctx);
811   ec_pow2 (D, Y1, ctx);
812
813   /* E = aC */
814   ec_mulm (E, ctx->a, C, ctx);
815
816   /* F = E + D */
817   ec_addm (F, E, D, ctx);
818
819   /* H = Z_1^2 */
820   ec_pow2 (H, Z1, ctx);
821
822   /* J = F - 2H */
823   ec_mulm (J, H, mpi_const (MPI_C_TWO), ctx);
824   ec_subm (J, F, J, ctx);
825
826   /* X_3 = (B - C - D) · J */
827   ec_subm (X3, B, C, ctx);
828   ec_subm (X3, X3, D, ctx);
829   ec_mulm (X3, X3, J, ctx);
830
831   /* Y_3 = F · (E - D) */
832   ec_subm (Y3, E, D, ctx);
833   ec_mulm (Y3, Y3, F, ctx);
834
835   /* Z_3 = F · J */
836   ec_mulm (Z3, F, J, ctx);
837
838 #undef X1
839 #undef Y1
840 #undef Z1
841 #undef X3
842 #undef Y3
843 #undef Z3
844 #undef B
845 #undef C
846 #undef D
847 #undef E
848 #undef F
849 #undef H
850 #undef J
851 }
852
853
854 /*  RESULT = 2 * POINT  */
855 void
856 _gcry_mpi_ec_dup_point (mpi_point_t result, mpi_point_t point, mpi_ec_t ctx)
857 {
858   switch (ctx->model)
859     {
860     case MPI_EC_WEIERSTRASS:
861       dup_point_weierstrass (result, point, ctx);
862       break;
863     case MPI_EC_MONTGOMERY:
864       dup_point_montgomery (result, point, ctx);
865       break;
866     case MPI_EC_TWISTEDEDWARDS:
867       dup_point_twistededwards (result, point, ctx);
868       break;
869     }
870 }
871
872
873 /* RESULT = P1 + P2  (Weierstrass version).*/
874 static void
875 add_points_weierstrass (mpi_point_t result,
876                         mpi_point_t p1, mpi_point_t p2,
877                         mpi_ec_t ctx)
878 {
879 #define x1 (p1->x    )
880 #define y1 (p1->y    )
881 #define z1 (p1->z    )
882 #define x2 (p2->x    )
883 #define y2 (p2->y    )
884 #define z2 (p2->z    )
885 #define x3 (result->x)
886 #define y3 (result->y)
887 #define z3 (result->z)
888 #define l1 (ctx->t.scratch[0])
889 #define l2 (ctx->t.scratch[1])
890 #define l3 (ctx->t.scratch[2])
891 #define l4 (ctx->t.scratch[3])
892 #define l5 (ctx->t.scratch[4])
893 #define l6 (ctx->t.scratch[5])
894 #define l7 (ctx->t.scratch[6])
895 #define l8 (ctx->t.scratch[7])
896 #define l9 (ctx->t.scratch[8])
897 #define t1 (ctx->t.scratch[9])
898 #define t2 (ctx->t.scratch[10])
899
900   if ( (!mpi_cmp (x1, x2)) && (!mpi_cmp (y1, y2)) && (!mpi_cmp (z1, z2)) )
901     {
902       /* Same point; need to call the duplicate function.  */
903       _gcry_mpi_ec_dup_point (result, p1, ctx);
904     }
905   else if (!mpi_cmp_ui (z1, 0))
906     {
907       /* P1 is at infinity.  */
908       mpi_set (x3, p2->x);
909       mpi_set (y3, p2->y);
910       mpi_set (z3, p2->z);
911     }
912   else if (!mpi_cmp_ui (z2, 0))
913     {
914       /* P2 is at infinity.  */
915       mpi_set (x3, p1->x);
916       mpi_set (y3, p1->y);
917       mpi_set (z3, p1->z);
918     }
919   else
920     {
921       int z1_is_one = !mpi_cmp_ui (z1, 1);
922       int z2_is_one = !mpi_cmp_ui (z2, 1);
923
924       /* l1 = x1 z2^2  */
925       /* l2 = x2 z1^2  */
926       if (z2_is_one)
927         mpi_set (l1, x1);
928       else
929         {
930           ec_pow2 (l1, z2, ctx);
931           ec_mulm (l1, l1, x1, ctx);
932         }
933       if (z1_is_one)
934         mpi_set (l2, x2);
935       else
936         {
937           ec_pow2 (l2, z1, ctx);
938           ec_mulm (l2, l2, x2, ctx);
939         }
940       /* l3 = l1 - l2 */
941       ec_subm (l3, l1, l2, ctx);
942       /* l4 = y1 z2^3  */
943       ec_powm (l4, z2, mpi_const (MPI_C_THREE), ctx);
944       ec_mulm (l4, l4, y1, ctx);
945       /* l5 = y2 z1^3  */
946       ec_powm (l5, z1, mpi_const (MPI_C_THREE), ctx);
947       ec_mulm (l5, l5, y2, ctx);
948       /* l6 = l4 - l5  */
949       ec_subm (l6, l4, l5, ctx);
950
951       if (!mpi_cmp_ui (l3, 0))
952         {
953           if (!mpi_cmp_ui (l6, 0))
954             {
955               /* P1 and P2 are the same - use duplicate function.  */
956               _gcry_mpi_ec_dup_point (result, p1, ctx);
957             }
958           else
959             {
960               /* P1 is the inverse of P2.  */
961               mpi_set_ui (x3, 1);
962               mpi_set_ui (y3, 1);
963               mpi_set_ui (z3, 0);
964             }
965         }
966       else
967         {
968           /* l7 = l1 + l2  */
969           ec_addm (l7, l1, l2, ctx);
970           /* l8 = l4 + l5  */
971           ec_addm (l8, l4, l5, ctx);
972           /* z3 = z1 z2 l3  */
973           ec_mulm (z3, z1, z2, ctx);
974           ec_mulm (z3, z3, l3, ctx);
975           /* x3 = l6^2 - l7 l3^2  */
976           ec_pow2 (t1, l6, ctx);
977           ec_pow2 (t2, l3, ctx);
978           ec_mulm (t2, t2, l7, ctx);
979           ec_subm (x3, t1, t2, ctx);
980           /* l9 = l7 l3^2 - 2 x3  */
981           ec_mulm (t1, x3, mpi_const (MPI_C_TWO), ctx);
982           ec_subm (l9, t2, t1, ctx);
983           /* y3 = (l9 l6 - l8 l3^3)/2  */
984           ec_mulm (l9, l9, l6, ctx);
985           ec_powm (t1, l3, mpi_const (MPI_C_THREE), ctx); /* fixme: Use saved value*/
986           ec_mulm (t1, t1, l8, ctx);
987           ec_subm (y3, l9, t1, ctx);
988           ec_mulm (y3, y3, ec_get_two_inv_p (ctx), ctx);
989         }
990     }
991
992 #undef x1
993 #undef y1
994 #undef z1
995 #undef x2
996 #undef y2
997 #undef z2
998 #undef x3
999 #undef y3
1000 #undef z3
1001 #undef l1
1002 #undef l2
1003 #undef l3
1004 #undef l4
1005 #undef l5
1006 #undef l6
1007 #undef l7
1008 #undef l8
1009 #undef l9
1010 #undef t1
1011 #undef t2
1012 }
1013
1014
1015 /* RESULT = P1 + P2  (Montgomery version).*/
1016 static void
1017 add_points_montgomery (mpi_point_t result,
1018                        mpi_point_t p1, mpi_point_t p2,
1019                        mpi_ec_t ctx)
1020 {
1021   (void)result;
1022   (void)p1;
1023   (void)p2;
1024   (void)ctx;
1025   log_fatal ("%s: %s not yet supported\n",
1026              "_gcry_mpi_ec_add_points", "Montgomery");
1027 }
1028
1029
1030 /* RESULT = P1 + P2  (Twisted Edwards version).*/
1031 static void
1032 add_points_twistededwards (mpi_point_t result,
1033                            mpi_point_t p1, mpi_point_t p2,
1034                            mpi_ec_t ctx)
1035 {
1036 #define X1 (p1->x)
1037 #define Y1 (p1->y)
1038 #define Z1 (p1->z)
1039 #define X2 (p2->x)
1040 #define Y2 (p2->y)
1041 #define Z2 (p2->z)
1042 #define X3 (result->x)
1043 #define Y3 (result->y)
1044 #define Z3 (result->z)
1045 #define A (ctx->t.scratch[0])
1046 #define B (ctx->t.scratch[1])
1047 #define C (ctx->t.scratch[2])
1048 #define D (ctx->t.scratch[3])
1049 #define E (ctx->t.scratch[4])
1050 #define F (ctx->t.scratch[5])
1051 #define G (ctx->t.scratch[6])
1052 #define tmp (ctx->t.scratch[7])
1053
1054   /* Compute: (X_3 : Y_3 : Z_3) = (X_1 : Y_1 : Z_1) + (X_2 : Y_2 : Z_3)  */
1055
1056   /* A = Z1 · Z2 */
1057   ec_mulm (A, Z1, Z2, ctx);
1058
1059   /* B = A^2 */
1060   ec_pow2 (B, A, ctx);
1061
1062   /* C = X1 · X2 */
1063   ec_mulm (C, X1, X2, ctx);
1064
1065   /* D = Y1 · Y2 */
1066   ec_mulm (D, Y1, Y2, ctx);
1067
1068   /* E = d · C · D */
1069   ec_mulm (E, ctx->b, C, ctx);
1070   ec_mulm (E, E, D, ctx);
1071
1072   /* F = B - E */
1073   ec_subm (F, B, E, ctx);
1074
1075   /* G = B + E */
1076   ec_addm (G, B, E, ctx);
1077
1078   /* X_3 = A · F · ((X_1 + Y_1) · (X_2 + Y_2) - C - D) */
1079   ec_addm (tmp, X1, Y1, ctx);
1080   ec_addm (X3, X2, Y2, ctx);
1081   ec_mulm (X3, X3, tmp, ctx);
1082   ec_subm (X3, X3, C, ctx);
1083   ec_subm (X3, X3, D, ctx);
1084   ec_mulm (X3, X3, F, ctx);
1085   ec_mulm (X3, X3, A, ctx);
1086
1087   /* Y_3 = A · G · (D - aC) */
1088   ec_mulm (Y3, ctx->a, C, ctx);
1089   ec_subm (Y3, D, Y3, ctx);
1090   ec_mulm (Y3, Y3, G, ctx);
1091   ec_mulm (Y3, Y3, A, ctx);
1092
1093   /* Z_3 = F · G */
1094   ec_mulm (Z3, F, G, ctx);
1095
1096
1097 #undef X1
1098 #undef Y1
1099 #undef Z1
1100 #undef X2
1101 #undef Y2
1102 #undef Z2
1103 #undef X3
1104 #undef Y3
1105 #undef Z3
1106 #undef A
1107 #undef B
1108 #undef C
1109 #undef D
1110 #undef E
1111 #undef F
1112 #undef G
1113 #undef tmp
1114 }
1115
1116
1117 /* RESULT = P1 + P2 */
1118 void
1119 _gcry_mpi_ec_add_points (mpi_point_t result,
1120                          mpi_point_t p1, mpi_point_t p2,
1121                          mpi_ec_t ctx)
1122 {
1123   switch (ctx->model)
1124     {
1125     case MPI_EC_WEIERSTRASS:
1126       add_points_weierstrass (result, p1, p2, ctx);
1127       break;
1128     case MPI_EC_MONTGOMERY:
1129       add_points_montgomery (result, p1, p2, ctx);
1130       break;
1131     case MPI_EC_TWISTEDEDWARDS:
1132       add_points_twistededwards (result, p1, p2, ctx);
1133       break;
1134     }
1135 }
1136
1137
1138 /* Scalar point multiplication - the main function for ECC.  If takes
1139    an integer SCALAR and a POINT as well as the usual context CTX.
1140    RESULT will be set to the resulting point. */
1141 void
1142 _gcry_mpi_ec_mul_point (mpi_point_t result,
1143                         gcry_mpi_t scalar, mpi_point_t point,
1144                         mpi_ec_t ctx)
1145 {
1146   gcry_mpi_t x1, y1, z1, k, h, yy;
1147   unsigned int i, loops;
1148   mpi_point_struct p1, p2, p1inv;
1149
1150   if (ctx->model == MPI_EC_TWISTEDEDWARDS)
1151     {
1152       /* Simple left to right binary method.  GECC Algorithm 3.27 */
1153       unsigned int nbits;
1154       int j;
1155
1156       nbits = mpi_get_nbits (scalar);
1157       mpi_set_ui (result->x, 0);
1158       mpi_set_ui (result->y, 1);
1159       mpi_set_ui (result->z, 1);
1160
1161       for (j=nbits-1; j >= 0; j--)
1162         {
1163           _gcry_mpi_ec_dup_point (result, result, ctx);
1164           if (mpi_test_bit (scalar, j) == 1)
1165             _gcry_mpi_ec_add_points (result, result, point, ctx);
1166         }
1167       return;
1168     }
1169
1170   x1 = mpi_alloc_like (ctx->p);
1171   y1 = mpi_alloc_like (ctx->p);
1172   h  = mpi_alloc_like (ctx->p);
1173   k  = mpi_copy (scalar);
1174   yy = mpi_copy (point->y);
1175
1176   if ( mpi_has_sign (k) )
1177     {
1178       k->sign = 0;
1179       ec_invm (yy, yy, ctx);
1180     }
1181
1182   if (!mpi_cmp_ui (point->z, 1))
1183     {
1184       mpi_set (x1, point->x);
1185       mpi_set (y1, yy);
1186     }
1187   else
1188     {
1189       gcry_mpi_t z2, z3;
1190
1191       z2 = mpi_alloc_like (ctx->p);
1192       z3 = mpi_alloc_like (ctx->p);
1193       ec_mulm (z2, point->z, point->z, ctx);
1194       ec_mulm (z3, point->z, z2, ctx);
1195       ec_invm (z2, z2, ctx);
1196       ec_mulm (x1, point->x, z2, ctx);
1197       ec_invm (z3, z3, ctx);
1198       ec_mulm (y1, yy, z3, ctx);
1199       mpi_free (z2);
1200       mpi_free (z3);
1201     }
1202   z1 = mpi_copy (mpi_const (MPI_C_ONE));
1203
1204   mpi_mul (h, k, mpi_const (MPI_C_THREE)); /* h = 3k */
1205   loops = mpi_get_nbits (h);
1206   if (loops < 2)
1207     {
1208       /* If SCALAR is zero, the above mpi_mul sets H to zero and thus
1209          LOOPs will be zero.  To avoid an underflow of I in the main
1210          loop we set LOOP to 2 and the result to (0,0,0).  */
1211       loops = 2;
1212       mpi_clear (result->x);
1213       mpi_clear (result->y);
1214       mpi_clear (result->z);
1215     }
1216   else
1217     {
1218       mpi_set (result->x, point->x);
1219       mpi_set (result->y, yy);
1220       mpi_set (result->z, point->z);
1221     }
1222   mpi_free (yy); yy = NULL;
1223
1224   p1.x = x1; x1 = NULL;
1225   p1.y = y1; y1 = NULL;
1226   p1.z = z1; z1 = NULL;
1227   point_init (&p2);
1228   point_init (&p1inv);
1229
1230   for (i=loops-2; i > 0; i--)
1231     {
1232       _gcry_mpi_ec_dup_point (result, result, ctx);
1233       if (mpi_test_bit (h, i) == 1 && mpi_test_bit (k, i) == 0)
1234         {
1235           point_set (&p2, result);
1236           _gcry_mpi_ec_add_points (result, &p2, &p1, ctx);
1237         }
1238       if (mpi_test_bit (h, i) == 0 && mpi_test_bit (k, i) == 1)
1239         {
1240           point_set (&p2, result);
1241           /* Invert point: y = p - y mod p  */
1242           point_set (&p1inv, &p1);
1243           ec_subm (p1inv.y, ctx->p, p1inv.y, ctx);
1244           _gcry_mpi_ec_add_points (result, &p2, &p1inv, ctx);
1245         }
1246     }
1247
1248   point_free (&p1);
1249   point_free (&p2);
1250   point_free (&p1inv);
1251   mpi_free (h);
1252   mpi_free (k);
1253 }
1254
1255
1256 /* Return true if POINT is on the curve described by CTX.  */
1257 int
1258 _gcry_mpi_ec_curve_point (gcry_mpi_point_t point, mpi_ec_t ctx)
1259 {
1260   int res = 0;
1261   gcry_mpi_t x, y, w;
1262
1263   x = mpi_new (0);
1264   y = mpi_new (0);
1265   w = mpi_new (0);
1266
1267   if (_gcry_mpi_ec_get_affine (x, y, point, ctx))
1268     return 0;
1269
1270   switch (ctx->model)
1271     {
1272     case MPI_EC_WEIERSTRASS:
1273       log_fatal ("%s: %s not yet supported\n",
1274                  "_gcry_mpi_ec_curve_point", "Weierstrass");
1275       break;
1276     case MPI_EC_MONTGOMERY:
1277       log_fatal ("%s: %s not yet supported\n",
1278                  "_gcry_mpi_ec_curve_point", "Montgomery");
1279       break;
1280     case MPI_EC_TWISTEDEDWARDS:
1281       {
1282         /* a · x^2 + y^2 - 1 - b · x^2 · y^2 == 0 */
1283         ec_pow2 (x, x, ctx);
1284         ec_pow2 (y, y, ctx);
1285         ec_mulm (w, ctx->a, x, ctx);
1286         ec_addm (w, w, y, ctx);
1287         ec_subm (w, w, mpi_const (MPI_C_ONE), ctx);
1288         ec_mulm (x, x, y, ctx);
1289         ec_mulm (x, x, ctx->b, ctx);
1290         ec_subm (w, w, x, ctx);
1291         if (!mpi_cmp_ui (w, 0))
1292           res = 1;
1293       }
1294       break;
1295     }
1296
1297   gcry_mpi_release (w);
1298   gcry_mpi_release (x);
1299   gcry_mpi_release (y);
1300
1301   return res;
1302 }