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