Loops like this one:
[openssl.git] / crypto / ec / ec_point.c
1 /*
2  *
3  *      ec_point.c
4  *
5  *      Elliptic Curve Arithmetic Functions
6  *
7  *      Copyright (C) Lenka Fibikova 2000
8  *
9  *
10  */
11
12 #include <stdio.h>
13 #include <stdlib.h>
14 #include <assert.h>
15 #include <memory.h>
16
17 #include <openssl/bn.h>
18
19 #include "../bn/bn_mont2.h" /* XXX */
20 #include "ec.h"
21
22
23 EC_POINT *ECP_new()
24         {
25         EC_POINT *ret;
26
27         ret=(EC_POINT *)malloc(sizeof(EC_POINT));
28         if (ret == NULL) return NULL;
29         ret->X = BN_new();
30         ret->Y = BN_new();
31         ret->Z = BN_new();
32         ret->is_in_mont = 0;
33
34         if (ret->X == NULL || ret->Y == NULL || ret->Z == NULL)
35                 {
36                 if (ret->X != NULL) BN_free(ret->X);
37                 if (ret->Y != NULL) BN_free(ret->Y);
38                 if (ret->Z != NULL) BN_free(ret->Z);
39                 free(ret);
40                 return(NULL);
41                 }
42         return(ret);
43         }
44
45
46 void ECP_clear_free(EC_POINT *P)
47         {
48         if (P == NULL) return;
49         
50         P->is_in_mont = 0;
51         if (P->X != NULL) BN_clear_free(P->X);
52         if (P->Y != NULL) BN_clear_free(P->Y);
53         if (P->Z != NULL) BN_clear_free(P->Z);
54         free(P);
55         }
56
57
58 void ECP_clear_free_precompute(ECP_PRECOMPUTE *prec)
59         {
60         int i;
61         int max;
62
63         if (prec == NULL) return;
64         if (prec->Pi != NULL)
65                 {
66                 max = 1;
67                 max <<= (prec->r - 1);
68
69                 for (i = 0; i < max; i++)
70                         {
71                         if (prec->Pi[i] != NULL) ECP_clear_free(prec->Pi[i]);
72                         }
73                 }
74         free(prec);
75         }
76
77
78 int ECP_is_on_ec(EC_POINT *P, EC *E, BN_CTX *ctx)
79         {
80         BIGNUM *n0, *n1, *n2, *p;
81         int Pnorm;
82         int ret = -1;
83
84         assert(P != NULL);
85         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
86
87         assert(E != NULL);
88         assert(E->A != NULL && E->B != NULL && E->p != NULL);
89
90         assert(ctx != NULL);
91
92         assert(!P->is_in_mont);
93
94         if (ECP_is_infty(P)) return 1;
95         
96         BN_CTX_start(ctx);
97         n0 = BN_CTX_get(ctx);
98         n1 = BN_CTX_get(ctx);
99         n2 = BN_CTX_get(ctx);
100         if (n2 == NULL)
101                 goto err;
102
103         p = E->p;
104
105         Pnorm = (ECP_is_norm(P));
106
107         if (!Pnorm)
108                 {
109                 if (!BN_mod_mul(n0, P->Z, P->Z, p, ctx)) goto err;
110                 if (!BN_mod_mul(n1, n0, n0, p, ctx)) goto err;
111                 if (!BN_mod_mul(n2, n0, n1, p, ctx)) goto err;
112                 }
113
114         if (!BN_mod_mul(n0, P->X, P->X, p, ctx)) goto err;
115         if (!BN_mod_mul(n0, n0, P->X, p, ctx)) goto err;
116
117         if (Pnorm)
118                 {
119                 if (!BN_mod_mul(n1, P->X, E->A, p, ctx)) goto err;
120                 }
121         else
122                 {
123                 if (!BN_mod_mul(n1, n1, P->X, p, ctx)) goto err;
124                 if (!BN_mod_mul(n1, n1, E->A, p, ctx)) goto err;
125                 }
126         if (!BN_mod_add(n0, n0, n1, p, ctx)) goto err;
127
128         if (Pnorm)
129                 {
130                 if (!BN_mod_add(n0, n0, E->B, p, ctx)) goto err;
131                 }
132         else
133                 {
134                 if (!BN_mod_mul(n2, n2, E->B,  p, ctx)) goto err;
135                 if (!BN_mod_add(n0, n0, n2, p, ctx)) goto err;
136                 }
137         
138         if (!BN_mod_mul(n1, P->Y, P->Y, p, ctx)) goto err;
139
140         if (BN_cmp(n0, n1))
141                 ret = 0;
142         else
143                 ret = 1;
144
145 err:
146         BN_CTX_end(ctx);
147         return ret;
148         }
149
150
151 EC_POINT *ECP_generate(BIGNUM *x, BIGNUM *z,EC *E, BN_CTX *ctx)
152 /* x == NULL || z = 0  -> point of infinity     */
153 /* z == NULL || z = 1  -> normalized            */
154         {
155         BIGNUM *n0, *n1;
156         EC_POINT *ret = NULL;
157         int Pnorm, Pinfty, X0, A0;
158
159         assert(E != NULL);
160         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
161
162         assert(ctx != NULL);
163
164         Pinfty = (x == NULL);
165         Pnorm = (z == NULL);
166         if (!Pnorm)
167                 {
168                 Pnorm = BN_is_one(z);
169                 Pinfty = (Pinfty || BN_is_zero(z));
170                 }
171
172         if (Pinfty)
173                 {
174                 if ((ret = ECP_new()) == NULL) return NULL;
175                 if (!BN_zero(ret->Z))
176                         {
177                         ECP_clear_free(ret);
178                         return NULL;
179                         }
180                 return ret;
181                 }
182
183         X0 = BN_is_zero(x);
184         A0 = BN_is_zero(E->A);
185
186         if ((ret = ECP_new()) == NULL) return NULL;
187
188         ret->is_in_mont = 0;
189
190         BN_CTX_start(ctx);
191         n0 = BN_CTX_get(ctx);
192         n1 = BN_CTX_get(ctx);
193         if (n1 == NULL) goto err;
194
195         if (!BN_zero(n0)) goto err;
196         if (!BN_zero(n1)) goto err;
197
198         if (!X0)
199                 {
200                 if (!BN_mod_sqr(n0, x, E->p, ctx)) goto err;
201                 if (!BN_mod_mul(n0, n0, x, E->p, ctx)) goto err;        /* x^3 */
202                 }
203
204         if (!X0 && !A0)
205                 {
206                 if (!BN_mod_mul(n1, E->A, x, E->p, ctx)) goto err;      /* Ax */
207                 if (!BN_mod_add(n0, n0, n1, E->p, ctx)) goto err;       /* x^3 + Ax */
208                 }
209
210         if (!BN_is_zero(E->B))
211                 if (!BN_mod_add(n0, n0, E->B, E->p, ctx)) goto err;     /* x^3 + Ax +B */
212
213         if (!BN_mod_sqrt(ret->Y, n0, E->p, ctx)) goto err;
214         if (BN_copy(ret->X, x) == NULL) goto err;
215         
216         if (Pnorm)
217                 {
218                 if (!BN_one(ret->Z)) goto err;
219                 }
220         else
221                 {
222                 if (BN_copy(ret->Z, z) == NULL) goto err;
223                 if (!BN_mod_sqr(n0, z, E->p, ctx)) goto err;
224                 if (!BN_mod_mul(ret->X, ret->X, n0, E->p, ctx)) goto err;
225                 if (!BN_mod_mul(n0, n0, z, E->p, ctx)) goto err;
226                 if (!BN_mod_mul(ret->Y, ret->Y, n0, E->p, ctx)) goto err;
227                 }
228
229 #ifdef TEST
230         if (!ECP_is_on_ec(ret, E, ctx)) goto err;
231 #endif
232         
233         BN_CTX_end(ctx);
234         return ret;
235
236 err:
237         if (ret != NULL) ECP_clear_free(ret);
238         BN_CTX_end(ctx);
239         return NULL;
240         }
241
242
243 int ECP_ecp2bin(EC_POINT *P, unsigned char *to, int form)
244 /* form =       1 ... compressed
245                 2 ... uncompressed
246                 3 ... hybrid */
247         {
248         int bytes, bx, by;
249
250         assert (P != NULL);
251         assert (P->X != NULL && P->Y != NULL && P->Z != NULL);
252         assert (!P->is_in_mont);
253         assert (ECP_is_norm(P) || ECP_is_infty(P));
254         assert (to != NULL);
255         assert (form > 0 && form < 4);
256
257         if (BN_is_zero(P->Z))
258                 {
259                 to[0] = 0;
260                 return 1;
261                 }
262
263         bx = BN_num_bytes(P->X);
264         if (form == 1 ) bytes = bx + 1;
265         else
266                 {
267                 by = BN_num_bytes(P->Y);
268                 bytes = (bx > by ? bx : by);
269                 bytes = bytes * 2 + 1;
270                 }
271         memset(to, 0, bytes);
272
273         switch (form)
274                 {
275         case 1: to[0] = 2;      break;
276         case 2: to[0] = 4;      break;
277         case 3: to[0] = 6;      break;
278                 }
279         if (form != 2) to[0] += BN_is_bit_set(P->Y, 0);
280
281         
282         if ((BN_bn2bin(P->X, to + 1)) != bx) return 0;
283         if (form != 1)
284                 {
285                 if ((BN_bn2bin(P->Y, to + bx + 1)) != by) return 0;
286                 }
287
288         return bytes;
289         }
290
291
292 int ECP_bin2ecp(unsigned char *from, int len, EC_POINT *P, EC *E, BN_CTX *ctx)
293         {
294         int y;
295         BIGNUM *x;
296         EC_POINT *pp;
297
298         assert (E != NULL);
299         assert (E->A != NULL && E->B != NULL && E->p != NULL);
300         assert (!E->is_in_mont);
301
302         assert (ctx != NULL);
303         assert (from != NULL);
304         assert (P != NULL);
305         assert (P->X != NULL && P->Y != NULL && P->Z != NULL);
306
307         if (len == 1 && from[0] != 0) return 0;
308
309         if (len == 0 || len == 1)
310                 {
311                 if (!BN_zero(P->Z)) return 0;
312                 return 1;
313                 }
314
315         switch (from[0])
316                 {
317         case 2:
318         case 3:
319                 y = from[0] - 2;
320                 if ((x = BN_new()) == NULL) return 0;
321                 if (BN_bin2bn(from + 1, len - 1, x) == NULL) return 0;
322
323                 pp = ECP_generate(x, NULL, E, ctx);
324                 BN_clear_free(x);
325                 if (pp == NULL) return 0;
326
327                 ECP_copy(P, pp);
328                 ECP_clear_free(pp);
329
330                 if (BN_is_bit_set(P->Y, 0) != y)
331                         if (!BN_sub(P->Y, E->p, P->Y)) return 0;
332                 break;
333
334         case 4:
335         case 6:
336         case 7:
337                 y = (len - 1)/2;
338                 if (BN_bin2bn(from + 1, y, P->X) == NULL) return 0;
339                 if (BN_bin2bn(from + y + 1, y, P->Y) == NULL) return 0;
340                 if (!BN_set_word(P->Z, 1)) return 0;
341                 break;
342
343         default:
344                 assert(0);
345
346                 }
347
348         if (!ECP_is_on_ec(P, E, ctx)) return 0;
349         return 1;
350         }
351
352
353 int ECP_normalize(EC_POINT *P, EC *E, BN_CTX *ctx)
354         {
355         BIGNUM *z, *zm;
356
357         assert (P != NULL);
358         assert (P->X != NULL && P->Y != NULL && P->Z != NULL);
359
360         assert (E != NULL);
361         assert (E->A != NULL && E->B != NULL && E->p != NULL);
362
363         assert (ctx != NULL);
364
365         if (ECP_is_norm(P)) return 1;
366         if (ECP_is_infty(P)) return 0;
367
368         if ((zm = BN_mod_inverse(P->Z, P->Z, E->p, ctx)) == NULL) return 0;
369
370         assert(!P->is_in_mont);
371
372
373         BN_CTX_start(ctx);
374         z = BN_CTX_get(ctx);
375         if (z == NULL) goto err;
376
377         if (!BN_mod_mul(z, zm, zm, E->p, ctx)) goto err;
378         if (!BN_mod_mul(P->X, P->X, z, E->p, ctx)) goto err;
379
380         if (!BN_mod_mul(z, z, zm, E->p, ctx)) goto err;
381         if (!BN_mod_mul(P->Y, P->Y, z, E->p, ctx)) goto err;
382
383         if (!BN_one(P->Z)) goto err;
384
385         if (zm != NULL) BN_clear_free(zm);
386
387         BN_CTX_end(ctx);
388         return 1;
389
390 err:
391         if (zm != NULL) BN_clear_free(zm);
392         BN_CTX_end(ctx);
393         return 0;
394         }
395
396
397 int ECP_copy(EC_POINT *R, EC_POINT *P)
398         {
399         assert(P != NULL);
400         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
401
402         assert(R != NULL);
403         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
404
405         if (BN_copy(R->X, P->X) == NULL) return 0;
406         if (BN_copy(R->Y, P->Y) == NULL) return 0;
407         if (BN_copy(R->Z, P->Z) == NULL) return 0;
408         R->is_in_mont = P->is_in_mont;
409
410         return 1;
411         }
412
413
414 EC_POINT *ECP_dup(EC_POINT *P)
415         {
416         EC_POINT *ret;
417
418         ret = ECP_new();
419         if (ret == NULL) return NULL;
420
421         if (!ECP_copy(ret, P))
422                 {
423                 ECP_clear_free(ret);
424                 return(NULL);
425                 }
426
427         return(ret);
428         }
429
430
431 EC_POINT *ECP_minus(EC_POINT *P, BIGNUM *p) /* mont || non-mont */
432         {
433         EC_POINT *ret;
434
435         assert(P != NULL);
436         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
437
438         assert(p != NULL);
439
440         assert(BN_cmp(P->Y, p) < 0);
441
442         ret = ECP_dup(P);
443         if (ret == NULL) return NULL;
444
445         if (BN_is_zero(ret->Y)) return ret;
446
447         if (!BN_sub(ret->Y, p, ret->Y))
448                 {
449                 ECP_clear_free(ret);
450                 return NULL;
451                 }
452
453         return ret;
454         }
455
456
457 #ifdef SIMPLE
458 int ECP_cmp(EC_POINT *P, EC_POINT *Q, BIGNUM *p, BN_CTX *ctx)
459 /* return values:
460         -2 ... error
461          0 ... P = Q
462         -1 ... P = -Q
463          1 ... else
464 */
465         {
466         BIGNUM *n0, *n1, *n2, *n3, *n4;
467         int Pnorm, Qnorm;
468
469         assert(P != NULL);
470         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
471
472         assert(Q != NULL);
473         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
474
475         assert(p != NULL);
476         assert(ctx != NULL);
477
478         assert(!P->is_in_mont);
479         assert(!Q->is_in_mont);
480
481         if (ECP_is_infty(P) && ECP_is_infty(Q)) return 0;
482         if (ECP_is_infty(P) || ECP_is_infty(Q)) return 1;
483
484         
485         Pnorm = (ECP_is_norm(P));
486         Qnorm = (ECP_is_norm(Q));
487         
488         BN_CTX_start(ctx);
489         n0 = BN_CTX_get(ctx);
490         n1 = BN_CTX_get(ctx);
491         n2 = BN_CTX_get(ctx);
492         n3 = BN_CTX_get(ctx);
493         n4 = BN_CTX_get(ctx);
494         if (n4 == NULL) goto err;
495         
496         if (Qnorm)
497                 {
498                 if (BN_copy(n1, P->X) == NULL) goto err;                        /* L1 = x_p */
499                 if (BN_copy(n2, P->Y) == NULL) goto err;                        /* L2 = y_p */
500                 }
501         else
502                 {
503                 if (!BN_sqr(n0, Q->Z, ctx)) goto err;
504                 if (!BN_mod_mul(n1, P->X, n0, p, ctx)) goto err;        /* L1 = x_p * z_q^2 */
505
506                 if (!BN_mod_mul(n0, n0, Q->Z, p, ctx)) goto err;
507                 if (!BN_mod_mul(n2, P->Y, n0, p, ctx)) goto err;        /* L2 = y_p * z_q^3 */
508                 }
509
510         if (Pnorm)
511                 {
512                 if (BN_copy(n3, Q->X) == NULL) goto err;                        /* L3 = x_q */
513                 if (BN_copy(n4, Q->Y) == NULL) goto err;                        /* L4 = y_q */
514                 }
515         else
516                 {
517                 if (!BN_sqr(n0, P->Z, ctx)) goto err;
518                 if (!BN_mod_mul(n3, Q->X, n0, p, ctx)) goto err;        /* L3 = x_q * z_p^2 */
519
520                 if (!BN_mod_mul(n0, n0, P->Z, p, ctx)) goto err;
521                 if (!BN_mod_mul(n4, Q->Y, n0, p, ctx)) goto err;        /* L4 = y_q * z_p^3 */
522                 }
523         
524         if (!BN_mod_sub(n0, n1, n3, p, ctx)) goto err;                  /* L5 = L1 - L3 */
525
526         if (!BN_is_zero(n0))
527                 {
528                 BN_CTX_end(ctx);
529                 return 1;
530                 }
531         
532         if (!BN_mod_sub(n0, n2, n4, p, ctx)) goto err;                  /* L6 = L2 - L4 */
533
534         if (!BN_is_zero(n0))
535                 {
536                 BN_CTX_end(ctx);
537                 return -1;
538                 }
539
540         BN_CTX_end(ctx);
541         return 0;
542
543 err:
544         BN_CTX_end(ctx);
545         return -2;
546         }
547
548
549 int ECP_double(EC_POINT *R, EC_POINT *P, EC *E, BN_CTX *ctx)
550 /* R <- 2P (on E) */
551         {
552         BIGNUM *n0, *n1, *n2, *n3, *p;
553         int Pnorm, A0;
554
555         assert(P != NULL);
556         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
557
558         assert(R != NULL);
559         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
560
561         assert(E != NULL);
562         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
563
564         assert(ctx != NULL);
565
566         assert(!P->is_in_mont);
567
568         if (ECP_is_infty(P))
569                 {
570                 if (!BN_zero(R->Z)) return 0;
571                 return 1;
572                 }
573
574         Pnorm = (ECP_is_norm(P));
575         A0 = (BN_is_zero(E->A));
576
577         BN_CTX_start(ctx);
578         n0 = BN_CTX_get(ctx);
579         n1 = BN_CTX_get(ctx);
580         n2 = BN_CTX_get(ctx);
581         n3 = BN_CTX_get(ctx);
582         if (n3 == NULL) goto err;
583
584         p = E->p;
585
586         /* L1 */
587         if (Pnorm || A0)
588                 {
589                 if (!BN_mod_sqr(n1, P->X, p, ctx)) goto err;
590                 if (!BN_mul_word(n1, 3)) goto err;
591                 if (!A0)                                   /* if A = 0: L1 = 3 * x^2 + a * z^4 = 3 * x ^2    */
592                         if (!BN_mod_add(n1, n1, E->A, p, ctx)) goto err; /* L1 = 3 * x^2 + a * z^4 = 3 * x^2 + a */
593                 }
594         else
595                 {
596                 if (!BN_mod_sqr(n0, P->Z, p, ctx)) goto err;
597                 if (!BN_mod_mul(n0, n0, n0, p, ctx)) goto err;
598                 if (!BN_mod_mul(n0, n0, E->A, p, ctx)) goto err;
599                 if (!BN_mod_sqr(n1, P->X, p, ctx)) goto err;
600                 if (!BN_mul_word(n1, 3)) goto err;
601                 if (!BN_mod_add(n1, n1, n0, p, ctx)) goto err;          /* L1 = 3 * x^2 + a * z^4 */
602                 }
603
604         /* Z */
605         if (Pnorm)
606                 {
607                 if (BN_copy(n0, P->Y) == NULL) goto err;
608                 }
609         else
610                 {
611                 if (!BN_mod_mul(n0, P->Y, P->Z, p, ctx)) goto err;
612                 }
613         if (!BN_lshift1(n0, n0)) goto err;
614         if (!BN_smod(R->Z, n0, p, ctx)) goto err;                               /* Z = 2 * y * z */
615
616         /* L2 */
617         if (!BN_mod_sqr(n3, P->Y, p, ctx)) goto err;
618         if (!BN_mod_mul(n2, P->X, n3, p, ctx)) goto err;
619         if (!BN_lshift(n2, n2, 2)) goto err;
620         if (!BN_smod(n2, n2, p, ctx)) goto err;                                 /* L2 = 4 * x * y^2 */
621
622         /* X */
623         if (!BN_lshift1(n0, n2)) goto err;
624         if (!BN_mod_sqr(R->X, n1, p, ctx)) goto err;
625         if (!BN_mod_sub(R->X, R->X, n0, p, ctx)) goto err;              /* X = L1^2 - 2 * L2 */
626         
627         /* L3 */
628         if (!BN_mod_sqr(n0, n3, p, ctx)) goto err;
629         if (!BN_lshift(n3, n0, 3)) goto err;
630         if (!BN_smod(n3, n3, p, ctx)) goto err;                                 /* L3 = 8 * y^4 */
631         
632         /* Y */
633         if (!BN_mod_sub(n0, n2, R->X, p, ctx)) goto err;
634         if (!BN_mod_mul(n0, n1, n0, p, ctx)) goto err;
635         if (!BN_mod_sub(R->Y, n0, n3, p, ctx)) goto err;                /* Y = L1 * (L2 - X) - L3 */
636
637
638 #ifdef TEST
639         if (!ECP_is_on_ec(R, E, ctx)) return 0;
640 #endif
641
642         BN_CTX_end(ctx);
643         return 1;
644
645 err:
646         BN_CTX_end(ctx);
647         return 0;
648         }
649
650
651 int ECP_add(EC_POINT *R, EC_POINT *P, EC_POINT *Q, EC *E, BN_CTX *ctx)
652 /* R <- P + Q (on E) */
653         {
654         BIGNUM *n0, *n1, *n2, *n3, *n4, *n5, *n6, *p;
655         int Pnorm, Qnorm;
656
657         assert(P != NULL);
658         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
659
660         assert(Q != NULL);
661         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
662
663         assert(R != NULL);
664         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
665
666         assert(E != NULL);
667         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
668         assert(!BN_is_zero(E->h));;
669
670         assert(ctx != NULL);
671
672         assert(!P->is_in_mont);
673         assert(!Q->is_in_mont);
674
675         if (P == Q) return ECP_double(R, P, E, ctx);
676
677         if (ECP_is_infty(P)) return ECP_copy(R, Q);
678         if (ECP_is_infty(Q)) return ECP_copy(R, P);
679         
680         Pnorm = (ECP_is_norm(P));
681         Qnorm = (ECP_is_norm(Q));
682         
683         BN_CTX_start(ctx);
684         n0 = BN_CTX_get(ctx);
685         n1 = BN_CTX_get(ctx);
686         n2 = BN_CTX_get(ctx);
687         n3 = BN_CTX_get(ctx);
688         n4 = BN_CTX_get(ctx);
689         n5 = BN_CTX_get(ctx);
690         n6 = BN_CTX_get(ctx);
691         if (n6 == NULL) goto err;
692
693         p = E->p;
694         
695         /* L1; L2 */
696         if (Qnorm)
697                 {
698                 if (BN_copy(n1, P->X) == NULL) goto err;         /* L1 = x_p */
699                 if (BN_copy(n2, P->Y) == NULL) goto err;         /* L2 = y_p */
700                 }
701         else
702                 {
703                 if (!BN_sqr(n0, Q->Z, ctx)) goto err;
704                 if (!BN_mod_mul(n1, P->X, n0, p, ctx)) goto err; /* L1 = x_p * z_q^2 */
705
706                 if (!BN_mod_mul(n0, n0, Q->Z, p, ctx)) goto err;
707                 if (!BN_mod_mul(n2, P->Y, n0, p, ctx)) goto err; /* L2 = y_p * z_q^3 */
708                 }
709
710         /* L3; L4 */
711         if (Pnorm)
712                 {
713                 if (BN_copy(n3, Q->X) == NULL) goto err;         /* L3 = x_q */
714                 if (BN_copy(n4, Q->Y) == NULL) goto err;         /* L4 = y_q */
715                 }
716         else
717                 {
718                 if (!BN_sqr(n0, P->Z, ctx)) goto err;
719                 if (!BN_mod_mul(n3, Q->X, n0, p, ctx)) goto err; /* L3 = x_q * z_p^2 */
720
721                 if (!BN_mod_mul(n0, n0, P->Z, p, ctx)) goto err;
722                 if (!BN_mod_mul(n4, Q->Y, n0, p, ctx)) goto err; /* L4 = y_q * z_p^3 */
723                 }
724
725         /* L5; L6 */
726         if (!BN_mod_sub(n5, n1, n3, p, ctx)) goto err;          /* L5 = L1 - L3 */
727         if (!BN_mod_sub(n6, n2, n4, p, ctx)) goto err;          /* L6 = L2 - L4 */
728
729         /* pata */
730         if (BN_is_zero(n5))
731                 {
732                 if (BN_is_zero(n6))     /* P = Q => P + Q = 2P */
733                         {
734                         BN_CTX_end(ctx);
735                         return ECP_double(R, P, E, ctx);
736                         }
737                 else                             /* P = -Q => P + Q = \infty */
738                         {
739                         BN_CTX_end(ctx);
740                         if (!BN_zero(R->Z)) return 0;
741                         return 1;
742                         }
743                 }
744
745         /* L7; L8 */
746         if (!BN_mod_add(n1, n1, n3, p, ctx)) goto err;          /* L7 = L1 + L3 */
747         if (!BN_mod_add(n2, n2, n4, p, ctx)) goto err;          /* L8 = L2 + L4 */
748
749         /* Z */
750         if (Pnorm)
751                 {
752                 if (BN_copy(n0, Q->Z) == NULL) goto err;
753                 }
754         else
755                 {
756                 if (!BN_mod_mul(n0, P->Z, Q->Z, p, ctx)) goto err;
757                 }
758         if (!BN_mod_mul(R->Z, n0, n5, p, ctx)) goto err;        /* Z = z_p * z_q * L_5 */
759
760         /* X */
761         if (!BN_mod_sqr(n0, n6, p, ctx)) goto err;
762         if (!BN_mod_sqr(n4, n5, p, ctx)) goto err;
763         if (!BN_mod_mul(n3, n1, n4, p, ctx)) goto err;
764         if (!BN_mod_sub(R->X, n0, n3, p, ctx)) goto err;        /* X = L6^2 - L5^2 * L7 */
765         
766         /* L9 */
767         if (!BN_lshift1(n0, R->X)) goto err;
768         if (!BN_mod_sub(n0, n3, n0, p, ctx)) goto err;          /* L9 = L5^2 * L7 - 2X */
769
770         /* Y */
771         if (!BN_mod_mul(n0, n0, n6, p, ctx)) goto err;
772         if (!BN_mod_mul(n5, n4, n5, p, ctx)) goto err;
773         if (!BN_mod_mul(n1, n2, n5, p, ctx)) goto err;
774         if (!BN_mod_sub(n0, n0, n1, p, ctx)) goto err;
775         if (!BN_mod_mul(R->Y, n0, E->h, p, ctx)) goto err;      /* Y = (L6 * L9 - L8 * L5^3) / 2 */
776
777
778
779 #ifdef TEST
780         if (!ECP_is_on_ec(R, E, ctx)) return 0;
781 #endif
782
783         BN_CTX_end(ctx);
784         return 1;
785
786 err:
787         BN_CTX_end(cxt);
788         return 0;
789         }
790
791
792 ECP_PRECOMPUTE *ECP_precompute(int r, EC_POINT *P, EC *E, BN_CTX *ctx)
793         {
794         ECP_PRECOMPUTE *ret;
795         EC_POINT *P2;
796         int i, max;
797
798         assert(r > 2);
799         assert(!P->is_in_mont);
800         assert(!E->is_in_mont);
801
802         ret=(ECP_PRECOMPUTE *)malloc(sizeof(ECP_PRECOMPUTE));
803         if (ret == NULL) return NULL;
804
805         max = 1;
806         max <<= (r - 1);
807
808         ret->r = 0;
809
810         ret->Pi=(EC_POINT **)malloc(sizeof(EC_POINT *) * max);
811         if (ret->Pi == NULL) goto err;
812
813         
814         /* P2 = [2]P */
815         if ((P2 = ECP_new()) == NULL) goto err;
816         if (!ECP_double(P2, P, E, ctx)) goto err;
817
818         /* P_0 = P */
819         if((ret->Pi[0] = ECP_dup(P)) == NULL) goto err;
820
821
822         /* P_i = P_(i-1) + P2 */
823         for (i = 1; i < max; i++)
824                 {
825                 if ((ret->Pi[i] = ECP_new()) == NULL) goto err;
826                 
827                 if (!ECP_add(ret->Pi[i], P2, ret->Pi[i - 1], E, ctx)) goto err;
828                 }
829
830         ret->r = r;
831         ECP_clear_free(P2);
832
833         return ret;
834
835 err:
836         ECP_clear_free(P2);
837         ECP_clear_free_precompute(ret);
838         return NULL;
839         }
840
841
842 int ECP_multiply(EC_POINT *R, BIGNUM *k, ECP_PRECOMPUTE *prec, EC *E, BN_CTX *ctx)
843 /* R = [k]P */
844         {
845         int j;
846         int t, nextw, h, r;
847
848         assert(R != NULL);
849         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
850
851         assert(E != NULL);
852         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
853
854         assert(k != NULL);
855         assert(!k->neg);
856
857         assert(ctx != NULL);
858         assert(prec != NULL);
859
860         assert(!E->is_in_mont);
861
862         if (BN_is_zero(k))
863                 {
864                 if (!BN_zero(R->Z)) return 0;
865                 R->is_in_mont = 0;
866                 return 1;
867                 }
868
869
870         j = BN_num_bits(k);
871         j--;
872
873         r = prec->r;
874
875         if (!BN_zero(R->Z)) return 0;
876         R->is_in_mont = 0;
877
878         while(j >= 0)
879                 {
880                 if (!BN_is_bit_set(k, j))
881                         {
882                         if (!ECP_double(R, R, E, ctx)) return 0;
883                         j--;
884                         }
885                 else
886                         {
887                         nextw = j - r;
888                         if (nextw < -1) nextw = -1;
889                         t = nextw + 1;                  
890                         while(!BN_is_bit_set(k, t))
891                                 t++;
892
893                         if (!ECP_double(R, R, E, ctx)) return 0;
894
895                         j--;
896                         if (j < t) h = 0;
897                         else
898                                 {
899                                 h = 1;
900                                 for(; j > t; j--)
901                                         {
902                                         h <<= 1;
903                                         if (BN_is_bit_set(k, j)) h++;
904                                         if (!ECP_double(R, R, E, ctx)) return 0;
905                                         }
906                                 if (!ECP_double(R, R, E, ctx)) return 0;
907                                 j--;
908                                 }
909
910                         if (!ECP_add(R, R, prec->Pi[h], E, ctx)) return 0;
911
912                         for (; j > nextw; j--)
913                                 {
914                                 if (!ECP_double(R, R, E, ctx)) return 0;
915                                 }
916
917                         }
918                 }
919         
920         return 1;
921         }
922
923 #endif /* SIMPLE */
924
925
926 #ifdef MONTGOMERY
927
928 int ECP_to_montgomery(EC_POINT *P, BN_MONTGOMERY *mont, BN_CTX *ctx)
929         {
930         assert(P != NULL);
931         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
932
933         assert(mont != NULL);
934         assert(mont->p != NULL);
935
936         assert(ctx != NULL);
937
938         if (P->is_in_mont) return 1;
939
940         if (!BN_lshift(P->X, P->X, mont->R_num_bits)) return 0;
941         if (!BN_mod(P->X, P->X, mont->p, ctx)) return 0;
942
943         if (!BN_lshift(P->Y, P->Y, mont->R_num_bits)) return 0;
944         if (!BN_mod(P->Y, P->Y, mont->p, ctx)) return 0;
945
946         if (!BN_lshift(P->Z, P->Z, mont->R_num_bits)) return 0;
947         if (!BN_mod(P->Z, P->Z, mont->p, ctx)) return 0;
948
949         P->is_in_mont = 1;
950         return 1;
951         }
952
953
954 int ECP_from_montgomery(EC_POINT *P, BN_MONTGOMERY *mont, BN_CTX *ctx)
955         {
956
957         assert(P != NULL);
958         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
959
960         assert(mont != NULL);
961         assert(mont->p != NULL);
962
963         assert(ctx != NULL);
964
965         if (!P->is_in_mont) return 1;
966
967         if (!BN_mont_red(P->X, mont)) return 0;
968         if (!BN_mont_red(P->Y, mont)) return 0;
969         if (!BN_mont_red(P->Z, mont)) return 0;
970
971         P->is_in_mont = 0;
972         return 1;
973         }
974
975
976 int ECP_mont_cmp(EC_POINT *P, EC_POINT *Q, BN_MONTGOMERY *mont, BN_CTX *ctx)
977 /* return values:
978         -2 ... error
979          0 ... P = Q
980         -1 ... P = -Q
981          1 ... else
982 */
983         {
984         BIGNUM *n0, *n1, *n2, *n3, *n4, *n5, *p;
985
986         assert(P != NULL);
987         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
988
989         assert(Q != NULL);
990         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
991
992         assert(mont != NULL);
993         assert(mont->p != NULL);
994
995         assert(ctx != NULL);
996
997         if (!P->is_in_mont)
998                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
999
1000         if (!Q->is_in_mont)
1001                 if (!ECP_to_montgomery(Q, mont, ctx)) return 0;
1002
1003
1004         if (ECP_is_infty(P) && ECP_is_infty(Q)) return 0;
1005         if (ECP_is_infty(P) || ECP_is_infty(Q)) return 1;
1006
1007         
1008         BN_CTX_start(ctx);
1009         n0 = BN_CTX_get(ctx);
1010         n1 = BN_CTX_get(ctx);
1011         n2 = BN_CTX_get(ctx);
1012         n3 = BN_CTX_get(ctx);
1013         n4 = BN_CTX_get(ctx);
1014         n5 = BN_CTX_get(ctx);
1015         if (n5 == 0) goto err;
1016
1017
1018         p = mont->p;
1019         
1020
1021         if (!BN_mont_mod_mul(n5, Q->Z, Q->Z, mont)) goto err;
1022         if (!BN_mont_mod_mul(n1, P->X, n5, mont)) goto err;     /* L1 = x_p * z_q^2 */
1023
1024         if (!BN_mont_mod_mul(n0, n5, Q->Z, mont)) goto err;
1025         if (!BN_mont_mod_mul(n2, P->Y, n0, mont)) goto err;     /* L2 = y_p * z_q^3 */
1026
1027         if (!BN_mont_mod_mul(n5, P->Z, P->Z, mont)) goto err;
1028         if (!BN_mont_mod_mul(n3, Q->X, n5, mont)) goto err;     /* L3 = x_q * z_p^2 */
1029
1030         if (!BN_mont_mod_mul(n0, n5, P->Z, mont)) goto err;
1031         if (!BN_mont_mod_mul(n4, Q->Y, n0, mont)) goto err;     /* L4 = y_q * z_p^3 */
1032
1033
1034         if (!BN_mod_sub_quick(n0, n1, n3, p)) goto err;                 /* L5 = L1 - L3 */
1035
1036         if (!BN_is_zero(n0))
1037                 {
1038                 BN_CTX_end(ctx);
1039                 return 1;
1040                 }
1041         
1042         if (!BN_mod_sub_quick(n0, n2, n4, p)) goto err;                 /* L6 = L2 - L4 */
1043
1044         if (!BN_is_zero(n0))
1045                 {
1046                 BN_CTX_end(ctx);
1047                 return -1;
1048                 }
1049
1050         BN_CTX_end(ctx);
1051         return 0;
1052
1053 err:
1054         BN_CTX_end(ctx);
1055         return -2;
1056         }
1057
1058
1059 int ECP_mont_double(EC_POINT *R, EC_POINT *P, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1060 /* R <- 2P (on E) */
1061         {
1062         BIGNUM *n0, *n1, *n2, *n3, *p;
1063
1064         assert(P != NULL);
1065         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
1066
1067         assert(R != NULL);
1068         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1069
1070         assert(E != NULL);
1071         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1072
1073         assert(ctx != NULL);
1074         
1075         if (!P->is_in_mont)
1076                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1077
1078         if (!E->is_in_mont)
1079                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1080
1081         R->is_in_mont = 1;
1082
1083         if (ECP_is_infty(P))
1084                 {
1085                 if (!BN_zero(R->Z)) return 0;
1086                 return 1;
1087                 }
1088
1089
1090         BN_CTX_start(ctx);
1091         n0 = BN_CTX_get(ctx);
1092         n1 = BN_CTX_get(ctx);
1093         n2 = BN_CTX_get(ctx);
1094         n3 = BN_CTX_get(ctx);
1095         if (n3 == 0) goto err;
1096
1097         p = E->p;
1098
1099         /* L1 */
1100         if (!BN_mont_mod_mul(n0, P->Z, P->Z, mont)) goto err;
1101         if (!BN_mont_mod_mul(n2, n0, n0, mont)) goto err;
1102         if (!BN_mont_mod_mul(n0, n2, E->A, mont)) goto err;
1103         if (!BN_mont_mod_mul(n1, P->X, P->X, mont)) goto err;
1104         if (!BN_mod_lshift1_quick(n2, n1, p)) goto err;
1105         if (!BN_mod_add_quick(n1, n1, n2, p)) goto err;
1106         if (!BN_mod_add_quick(n1, n1, n0, p)) goto err;         /* L1 = 3 * x^2 + a * z^4 */
1107
1108         /* Z */
1109         if (!BN_mont_mod_mul(n0, P->Y, P->Z, mont)) goto err;
1110         if (!BN_mod_lshift1_quick(R->Z, n0, p)) goto err;               /* Z = 2 * y * z */
1111
1112         /* L2 */
1113         if (!BN_mont_mod_mul(n3, P->Y, P->Y, mont)) goto err;
1114         if (!BN_mont_mod_mul(n2, P->X, n3, mont)) goto err;
1115         if (!BN_mod_lshift_quick(n2, n2, 2, p)) goto err;               /* L2 = 4 * x * y^2 */
1116
1117         /* X */
1118         if (!BN_mod_lshift1_quick(n0, n2, p)) goto err;
1119         if (!BN_mont_mod_mul(R->X, n1, n1, mont)) goto err;
1120         if (!BN_mod_sub_quick(R->X, R->X, n0, p)) goto err;     /* X = L1^2 - 2 * L2 */
1121         
1122         /* L3 */
1123         if (!BN_mont_mod_mul(n0, n3, n3, mont)) goto err;
1124         if (!BN_mod_lshift_quick(n3, n0, 3, p)) goto err;               /* L3 = 8 * y^4 */
1125
1126         
1127         /* Y */
1128         if (!BN_mod_sub_quick(n2, n2, R->X, p)) goto err;
1129         if (!BN_mont_mod_mul(n0, n1, n2, mont)) goto err;
1130         if (!BN_mod_sub_quick(R->Y, n0, n3, p)) goto err;               /* Y = L1 * (L2 - X) - L3 */
1131
1132         BN_CTX_end(ctx);
1133         return 1;
1134
1135 err:
1136         BN_CTX_end(ctx);
1137         return 0;
1138         }
1139
1140
1141 int ECP_mont_add(EC_POINT *R, EC_POINT *P, EC_POINT *Q, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1142 /* R <- P + Q (on E) */
1143         {
1144         BIGNUM *n0, *n1, *n2, *n3, *n4, *n5, *n6, *p;
1145
1146         assert(P != NULL);
1147         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
1148
1149         assert(Q != NULL);
1150         assert(Q->X != NULL && Q->Y != NULL && Q->Z != NULL);
1151
1152         assert(R != NULL);
1153         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1154
1155         assert(E != NULL);
1156         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1157         assert(!BN_is_zero(E->h));;
1158
1159         assert(ctx != NULL);
1160
1161         if (!Q->is_in_mont)
1162                 if (!ECP_to_montgomery(Q, mont, ctx)) return 0;
1163
1164         if (!P->is_in_mont)
1165                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1166
1167         if (!E->is_in_mont)
1168                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1169
1170         if (P == Q) return ECP_mont_double(R, P, E, mont, ctx);
1171
1172         if (ECP_is_infty(P)) return ECP_copy(R, Q);
1173         if (ECP_is_infty(Q)) return ECP_copy(R, P);
1174         
1175
1176         BN_CTX_start(ctx);
1177         n0 = BN_CTX_get(ctx);
1178         n1 = BN_CTX_get(ctx);
1179         n2 = BN_CTX_get(ctx);
1180         n3 = BN_CTX_get(ctx);
1181         n4 = BN_CTX_get(ctx);
1182         n5 = BN_CTX_get(ctx);
1183         n6 = BN_CTX_get(ctx);
1184         if (n6 == NULL) goto err;
1185
1186
1187         p = E->p;
1188
1189         R->is_in_mont = 1;
1190         
1191         /* L1; L2 */
1192         if (!BN_mont_mod_mul(n6, Q->Z, Q->Z, mont)) goto err;
1193         if (!BN_mont_mod_mul(n1, P->X, n6, mont)) goto err;     /* L1 = x_p * z_q^2 */
1194
1195         if (!BN_mont_mod_mul(n0, n6, Q->Z, mont)) goto err;
1196         if (!BN_mont_mod_mul(n2, P->Y, n0, mont)) goto err;     /* L2 = y_p * z_q^3 */
1197
1198
1199         /* L3; L4 */
1200         if (!BN_mont_mod_mul(n6, P->Z, P->Z, mont)) goto err;
1201         if (!BN_mont_mod_mul(n3, Q->X, n6, mont)) goto err;     /* L3 = x_q * z_p^2 */
1202
1203         if (!BN_mont_mod_mul(n0, n6, P->Z, mont)) goto err;
1204         if (!BN_mont_mod_mul(n4, Q->Y, n0, mont)) goto err;     /* L4 = y_q * z_p^3 */
1205
1206
1207         /* L5; L6 */
1208         if (!BN_mod_sub_quick(n5, n1, n3, p)) goto err;                 /* L5 = L1 - L3 */
1209         if (!BN_mod_sub_quick(n6, n2, n4, p)) goto err;                 /*L6 = L2 - L4 */
1210
1211
1212         /* pata */
1213         if (BN_is_zero(n5))
1214                 {
1215                 if (BN_is_zero(n6))  /* P = Q => P + Q = 2P */
1216                         {
1217                         BN_CTX_end(ctx);
1218                         return ECP_mont_double(R, P, E, mont, ctx);
1219                         }
1220                 else                             /* P = -Q => P + Q = \infty */
1221                         {
1222                         BN_CTX_end(ctx);
1223                         if (!BN_zero(R->Z)) return 0;
1224                         return 1;
1225                         }
1226                 }
1227
1228         /* L7; L8 */
1229         if (!BN_mod_add_quick(n1, n1, n3, p)) goto err;                 /* L7 = L1 + L3 */
1230         if (!BN_mod_add_quick(n2, n2, n4, p)) goto err;                 /* L8 = L2 + L4 */
1231
1232
1233         /* Z */
1234         if (!BN_mont_mod_mul(n0, P->Z, Q->Z, mont)) goto err;
1235         if (!BN_mont_mod_mul(R->Z, n0, n5, mont)) goto err;     /* Z = z_p * z_q * L_5 */
1236
1237
1238         /* X */
1239         if (!BN_mont_mod_mul(n0, n6, n6, mont)) goto err;
1240         if (!BN_mont_mod_mul(n4, n5, n5, mont)) goto err;
1241         if (!BN_mont_mod_mul(n3, n1, n4, mont)) goto err;
1242         if (!BN_mod_sub_quick(R->X, n0, n3, p)) goto err;                       /* X = L6^2 - L5^2 * L7 */
1243
1244         
1245         /* L9 */
1246         if (!BN_mod_lshift1_quick(n0, R->X, p)) goto err;
1247         if (!BN_mod_sub_quick(n3, n3, n0, p)) goto err;                 /* L9 = L5^2 * L7 - 2X */
1248
1249
1250         /* Y */
1251         if (!BN_mont_mod_mul(n0, n3, n6, mont)) goto err;
1252         if (!BN_mont_mod_mul(n6, n4, n5, mont)) goto err;
1253         if (!BN_mont_mod_mul(n1, n2, n6, mont)) goto err;
1254         if (!BN_mod_sub_quick(n0, n0, n1, p)) goto err;
1255         if (!BN_mont_mod_mul(R->Y, n0, E->h, mont)) goto err;   /* Y = (L6 * L9 - L8 * L5^3) / 2 */
1256
1257
1258         BN_CTX_end(ctx);
1259         return 1;
1260
1261 err:
1262         BN_CTX_end(ctx);
1263         return 0;
1264         }
1265
1266
1267 ECP_PRECOMPUTE *ECP_mont_precompute(int r, EC_POINT *P, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1268         {
1269         ECP_PRECOMPUTE *ret;
1270         EC_POINT *P2;
1271         int i, max;
1272
1273         assert(r > 2);
1274         assert(r < sizeof(unsigned int) * 8 - 1);
1275
1276         assert(mont != NULL);
1277         assert(mont->p != NULL);
1278         
1279         if (!P->is_in_mont)
1280                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1281
1282         if (!E->is_in_mont)
1283                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1284
1285         ret=(ECP_PRECOMPUTE *)malloc(sizeof(ECP_PRECOMPUTE));
1286         if (ret == NULL) return NULL;
1287
1288         max = 1;
1289         max <<= (r - 1);
1290
1291         ret->r = 0;
1292
1293         ret->Pi=(EC_POINT **)malloc(sizeof(EC_POINT *) * max);
1294         if (ret->Pi == NULL) goto err;
1295
1296         
1297         /* P2 = [2]P */
1298         if ((P2 = ECP_new()) == NULL) goto err;
1299         if (!ECP_mont_double(P2, P, E, mont, ctx)) goto err;
1300
1301         /* P_0 = P */
1302         if((ret->Pi[0] = ECP_dup(P)) == NULL) goto err;
1303
1304
1305         /* P_i = P_(i-1) + P2 */
1306         for (i = 1; i < max; i++)
1307                 {
1308                 if ((ret->Pi[i] = ECP_new()) == NULL) goto err;
1309                 if (!ECP_mont_add(ret->Pi[i], P2, ret->Pi[i - 1], E, mont, ctx)) goto err;
1310                 }
1311
1312         ret->r = r;
1313         ECP_clear_free(P2);
1314
1315         return ret;
1316
1317 err:
1318         ECP_clear_free(P2);
1319         ECP_clear_free_precompute(ret);
1320         return NULL;
1321         }
1322
1323
1324 int ECP_mont_multiply(EC_POINT *R, BIGNUM *k, ECP_PRECOMPUTE *prec, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1325 /* R = [k]P   P = prec->Pi[0]*/
1326         {
1327         int j;
1328         int t, nextw, h, r;
1329
1330         assert(R != NULL);
1331         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1332
1333         assert(E != NULL);
1334         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1335
1336         assert(k != NULL);
1337         assert(!k->neg);
1338
1339         assert(ctx != NULL);
1340         assert(prec != NULL);
1341
1342         assert(mont != NULL);
1343         assert(mont->p != NULL);
1344
1345         if (!E->is_in_mont)
1346                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1347
1348
1349         if (BN_is_zero(k))
1350                 {
1351                 if (!BN_zero(R->Z)) return 0;
1352                 R->is_in_mont = 1;
1353                 return 1;
1354                 }
1355
1356         j = BN_num_bits(k);
1357         j--;
1358
1359         r = prec->r;
1360
1361         if (!BN_zero(R->Z)) return 0;
1362         R->is_in_mont = 1;
1363
1364         while(j >= 0)
1365                 {
1366                 if (!BN_is_bit_set(k, j))
1367                         {
1368                         if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1369                         j--;
1370                         }
1371                 else
1372                         {
1373                         nextw = j - r;
1374                         if (nextw < -1) nextw = -1;
1375                         t = nextw + 1;                  
1376                         while(!BN_is_bit_set(k, t))
1377                                 t++;
1378
1379                         if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1380
1381                         j--;
1382                         if (j < t) h = 0;
1383                         else
1384                                 {
1385                                 h = 1;
1386                                 for(; j > t; j--)
1387                                         {
1388                                         h <<= 1;
1389                                         if (BN_is_bit_set(k, j)) h++;
1390                                         if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1391                                         }
1392                                 if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1393                                 j--;
1394                                 }
1395
1396                         if (!ECP_mont_add(R, R, prec->Pi[h], E, mont, ctx)) return 0;
1397
1398                         for (; j > nextw; j--)
1399                                 {
1400                                 if (!ECP_mont_double(R, R, E, mont, ctx)) return 0;
1401                                 }
1402
1403                         }
1404                 }
1405
1406         return 1;
1407         }
1408
1409
1410 int ECP_mont_multiply2(EC_POINT *R, BIGNUM *k, EC_POINT *P, EC *E, BN_MONTGOMERY *mont, BN_CTX *ctx)
1411 /* R = [k]P */
1412         {
1413         int j, hj, kj;
1414         BIGNUM *h;
1415         EC_POINT *mP;
1416
1417         assert(R != NULL);
1418         assert(R->X != NULL && R->Y != NULL && R->Z != NULL);
1419
1420         assert(P != NULL);
1421         assert(P->X != NULL && P->Y != NULL && P->Z != NULL);
1422
1423         assert(E != NULL);
1424         assert(E->A != NULL && E->B != NULL && E->p != NULL && E->h != NULL);
1425
1426         assert(k != NULL);
1427         assert(!k->neg);
1428
1429         assert(ctx != NULL);
1430
1431         assert(mont != NULL);
1432         assert(mont->p != NULL);
1433
1434         if (!E->is_in_mont)
1435                 if (!EC_to_montgomery(E, mont, ctx)) return 0;
1436
1437         if (!P->is_in_mont)
1438                 if (!ECP_to_montgomery(P, mont, ctx)) return 0;
1439
1440
1441         if (BN_is_zero(k))
1442                 {
1443                 if (!BN_zero(R->Z)) return 0;
1444                 R->is_in_mont = 1;
1445                 return 1;
1446                 }
1447
1448         if ((h = BN_dup(k)) == NULL) return 0;
1449         
1450         if (!BN_lshift1(h, h)) goto err;
1451         if (!BN_add(h, h, k)) goto err;
1452
1453         if (!ECP_copy(R, P)) goto err;
1454         if ((mP = ECP_mont_minus(P, mont)) == NULL) goto err;
1455
1456         for(j = BN_num_bits(h) - 2; j > 0; j--)
1457                 {
1458                 if (!ECP_mont_double(R, R, E, mont, ctx)) goto err;
1459                 kj = BN_is_bit_set(k, j);
1460                 hj = BN_is_bit_set(h, j);
1461                 if (hj == 1 && kj == 0)
1462                         if (!ECP_mont_add(R, R, P, E, mont, ctx)) goto err;
1463                 if (hj == 0 && kj == 1)
1464                         if (!ECP_mont_add(R, R, mP, E, mont, ctx)) goto err;
1465                 }
1466
1467         if (h != NULL) BN_free(h);
1468         if (mP != NULL) ECP_clear_free(mP);
1469         return 1;
1470
1471 err:
1472         if (h != NULL) BN_free(h);
1473         if (mP != NULL) ECP_clear_free(mP);
1474         return 0;
1475         }
1476
1477 #endif /* MONTGOMERY */