Fixlets to debian package build rules.
[deb_shairplay.git] / src / lib / crypto / bigint.c
1 /*
2 * Copyright (c) 2007, Cameron Rich
3 *
4 * All rights reserved.
5 *
6 * Redistribution and use in source and binary forms, with or without
7 * modification, are permitted provided that the following conditions are met:
8 *
9 * * Redistributions of source code must retain the above copyright notice,
10 * this list of conditions and the following disclaimer.
11 * * Redistributions in binary form must reproduce the above copyright notice,
12 * this list of conditions and the following disclaimer in the documentation
13 * and/or other materials provided with the distribution.
14 * * Neither the name of the axTLS project nor the names of its contributors
15 * may be used to endorse or promote products derived from this software
16 * without specific prior written permission.
17 *
18 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
19 * "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
20 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
21 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR
22 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
23 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
24 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
25 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
26 * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
27 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
28 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29 */
30
31 /**
32 * @defgroup bigint_api Big Integer API
33 * @brief The bigint implementation as used by the axTLS project.
34 *
35 * The bigint library is for RSA encryption/decryption as well as signing.
36 * This code tries to minimise use of malloc/free by maintaining a small
37 * cache. A bigint context may maintain state by being made "permanent".
38 * It be be later released with a bi_depermanent() and bi_free() call.
39 *
40 * It supports the following reduction techniques:
41 * - Classical
42 * - Barrett
43 * - Montgomery
44 *
45 * It also implements the following:
46 * - Karatsuba multiplication
47 * - Squaring
48 * - Sliding window exponentiation
49 * - Chinese Remainder Theorem (implemented in rsa.c).
50 *
51 * All the algorithms used are pretty standard, and designed for different
52 * data bus sizes. Negative numbers are not dealt with at all, so a subtraction
53 * may need to be tested for negativity.
54 *
55 * This library steals some ideas from Jef Poskanzer
56 * <http://cs.marlboro.edu/term/cs-fall02/algorithms/crypto/RSA/bigint>
57 * and GMP <http://www.swox.com/gmp>. It gets most of its implementation
58 * detail from "The Handbook of Applied Cryptography"
59 * <http://www.cacr.math.uwaterloo.ca/hac/about/chap14.pdf>
60 * @{
61 */
62
63 #include <stdlib.h>
64 #include <limits.h>
65 #include <string.h>
66 #include <stdio.h>
67 #include <time.h>
68 #include "os_port.h"
69 #include "bigint.h"
70
71 #define V1 v->comps[v->size-1] /**< v1 for division */
72 #define V2 v->comps[v->size-2] /**< v2 for division */
73 #define U(j) tmp_u->comps[tmp_u->size-j-1] /**< uj for division */
74 #define Q(j) quotient->comps[quotient->size-j-1] /**< qj for division */
75
76 static bigint *bi_int_multiply(BI_CTX *ctx, bigint *bi, comp i);
77 static bigint *bi_int_divide(BI_CTX *ctx, bigint *biR, comp denom);
78 static bigint *alloc(BI_CTX *ctx, int size);
79 static bigint *trim(bigint *bi);
80 static void more_comps(bigint *bi, int n);
81 #if defined(CONFIG_BIGINT_KARATSUBA) || defined(CONFIG_BIGINT_BARRETT) || \
82 defined(CONFIG_BIGINT_MONTGOMERY)
83 static bigint *comp_right_shift(bigint *biR, int num_shifts);
84 static bigint *comp_left_shift(bigint *biR, int num_shifts);
85 #endif
86
87 #ifdef CONFIG_BIGINT_CHECK_ON
88 static void check(const bigint *bi);
89 #else
90 #define check(A) /**< disappears in normal production mode */
91 #endif
92
93
94 /**
95 * @brief Start a new bigint context.
96 * @return A bigint context.
97 */
98 BI_CTX *bi_initialize(void)
99 {
100 /* calloc() sets everything to zero */
101 BI_CTX *ctx = (BI_CTX *)calloc(1, sizeof(BI_CTX));
102
103 /* the radix */
104 ctx->bi_radix = alloc(ctx, 2);
105 ctx->bi_radix->comps[0] = 0;
106 ctx->bi_radix->comps[1] = 1;
107 bi_permanent(ctx->bi_radix);
108 return ctx;
109 }
110
111 /**
112 * @brief Close the bigint context and free any resources.
113 *
114 * Free up any used memory - a check is done if all objects were not
115 * properly freed.
116 * @param ctx [in] The bigint session context.
117 */
118 void bi_terminate(BI_CTX *ctx)
119 {
120 bi_depermanent(ctx->bi_radix);
121 bi_free(ctx, ctx->bi_radix);
122
123 if (ctx->active_count != 0)
124 {
125 #ifdef CONFIG_SSL_FULL_MODE
126 printf("bi_terminate: there were %d un-freed bigints\n",
127 ctx->active_count);
128 #endif
129 abort();
130 }
131
132 bi_clear_cache(ctx);
133 free(ctx);
134 }
135
136 /**
137 *@brief Clear the memory cache.
138 */
139 void bi_clear_cache(BI_CTX *ctx)
140 {
141 bigint *p, *pn;
142
143 if (ctx->free_list == NULL)
144 return;
145
146 for (p = ctx->free_list; p != NULL; p = pn)
147 {
148 pn = p->next;
149 free(p->comps);
150 free(p);
151 }
152
153 ctx->free_count = 0;
154 ctx->free_list = NULL;
155 }
156
157 /**
158 * @brief Increment the number of references to this object.
159 * It does not do a full copy.
160 * @param bi [in] The bigint to copy.
161 * @return A reference to the same bigint.
162 */
163 bigint *bi_copy(bigint *bi)
164 {
165 check(bi);
166 if (bi->refs != PERMANENT)
167 bi->refs++;
168 return bi;
169 }
170
171 /**
172 * @brief Simply make a bigint object "unfreeable" if bi_free() is called on it.
173 *
174 * For this object to be freed, bi_depermanent() must be called.
175 * @param bi [in] The bigint to be made permanent.
176 */
177 void bi_permanent(bigint *bi)
178 {
179 check(bi);
180 if (bi->refs != 1)
181 {
182 #ifdef CONFIG_SSL_FULL_MODE
183 printf("bi_permanent: refs was not 1\n");
184 #endif
185 abort();
186 }
187
188 bi->refs = PERMANENT;
189 }
190
191 /**
192 * @brief Take a permanent object and make it eligible for freedom.
193 * @param bi [in] The bigint to be made back to temporary.
194 */
195 void bi_depermanent(bigint *bi)
196 {
197 check(bi);
198 if (bi->refs != PERMANENT)
199 {
200 #ifdef CONFIG_SSL_FULL_MODE
201 printf("bi_depermanent: bigint was not permanent\n");
202 #endif
203 abort();
204 }
205
206 bi->refs = 1;
207 }
208
209 /**
210 * @brief Free a bigint object so it can be used again.
211 *
212 * The memory itself it not actually freed, just tagged as being available
213 * @param ctx [in] The bigint session context.
214 * @param bi [in] The bigint to be freed.
215 */
216 void bi_free(BI_CTX *ctx, bigint *bi)
217 {
218 check(bi);
219 if (bi->refs == PERMANENT)
220 {
221 return;
222 }
223
224 if (--bi->refs > 0)
225 {
226 return;
227 }
228
229 bi->next = ctx->free_list;
230 ctx->free_list = bi;
231 ctx->free_count++;
232
233 if (--ctx->active_count < 0)
234 {
235 #ifdef CONFIG_SSL_FULL_MODE
236 printf("bi_free: active_count went negative "
237 "- double-freed bigint?\n");
238 #endif
239 abort();
240 }
241 }
242
243 /**
244 * @brief Convert an (unsigned) integer into a bigint.
245 * @param ctx [in] The bigint session context.
246 * @param i [in] The (unsigned) integer to be converted.
247 *
248 */
249 bigint *int_to_bi(BI_CTX *ctx, comp i)
250 {
251 bigint *biR = alloc(ctx, 1);
252 biR->comps[0] = i;
253 return biR;
254 }
255
256 /**
257 * @brief Do a full copy of the bigint object.
258 * @param ctx [in] The bigint session context.
259 * @param bi [in] The bigint object to be copied.
260 */
261 bigint *bi_clone(BI_CTX *ctx, const bigint *bi)
262 {
263 bigint *biR = alloc(ctx, bi->size);
264 check(bi);
265 memcpy(biR->comps, bi->comps, bi->size*COMP_BYTE_SIZE);
266 return biR;
267 }
268
269 /**
270 * @brief Perform an addition operation between two bigints.
271 * @param ctx [in] The bigint session context.
272 * @param bia [in] A bigint.
273 * @param bib [in] Another bigint.
274 * @return The result of the addition.
275 */
276 bigint *bi_add(BI_CTX *ctx, bigint *bia, bigint *bib)
277 {
278 int n;
279 comp carry = 0;
280 comp *pa, *pb;
281
282 check(bia);
283 check(bib);
284
285 n = max(bia->size, bib->size);
286 more_comps(bia, n+1);
287 more_comps(bib, n);
288 pa = bia->comps;
289 pb = bib->comps;
290
291 do
292 {
293 comp sl, rl, cy1;
294 sl = *pa + *pb++;
295 rl = sl + carry;
296 cy1 = sl < *pa;
297 carry = cy1 | (rl < sl);
298 *pa++ = rl;
299 } while (--n != 0);
300
301 *pa = carry; /* do overflow */
302 bi_free(ctx, bib);
303 return trim(bia);
304 }
305
306 /**
307 * @brief Perform a subtraction operation between two bigints.
308 * @param ctx [in] The bigint session context.
309 * @param bia [in] A bigint.
310 * @param bib [in] Another bigint.
311 * @param is_negative [out] If defined, indicates that the result was negative.
312 * is_negative may be null.
313 * @return The result of the subtraction. The result is always positive.
314 */
315 bigint *bi_subtract(BI_CTX *ctx,
316 bigint *bia, bigint *bib, int *is_negative)
317 {
318 int n = bia->size;
319 comp *pa, *pb, carry = 0;
320
321 check(bia);
322 check(bib);
323
324 more_comps(bib, n);
325 pa = bia->comps;
326 pb = bib->comps;
327
328 do
329 {
330 comp sl, rl, cy1;
331 sl = *pa - *pb++;
332 rl = sl - carry;
333 cy1 = sl > *pa;
334 carry = cy1 | (rl > sl);
335 *pa++ = rl;
336 } while (--n != 0);
337
338 if (is_negative) /* indicate a negative result */
339 {
340 *is_negative = carry;
341 }
342
343 bi_free(ctx, trim(bib)); /* put bib back to the way it was */
344 return trim(bia);
345 }
346
347 /**
348 * Perform a multiply between a bigint an an (unsigned) integer
349 */
350 static bigint *bi_int_multiply(BI_CTX *ctx, bigint *bia, comp b)
351 {
352 int j = 0, n = bia->size;
353 bigint *biR = alloc(ctx, n + 1);
354 comp carry = 0;
355 comp *r = biR->comps;
356 comp *a = bia->comps;
357
358 check(bia);
359
360 /* clear things to start with */
361 memset(r, 0, ((n+1)*COMP_BYTE_SIZE));
362
363 do
364 {
365 long_comp tmp = *r + (long_comp)a[j]*b + carry;
366 *r++ = (comp)tmp; /* downsize */
367 carry = (comp)(tmp >> COMP_BIT_SIZE);
368 } while (++j < n);
369
370 *r = carry;
371 bi_free(ctx, bia);
372 return trim(biR);
373 }
374
375 /**
376 * @brief Does both division and modulo calculations.
377 *
378 * Used extensively when doing classical reduction.
379 * @param ctx [in] The bigint session context.
380 * @param u [in] A bigint which is the numerator.
381 * @param v [in] Either the denominator or the modulus depending on the mode.
382 * @param is_mod [n] Determines if this is a normal division (0) or a reduction
383 * (1).
384 * @return The result of the division/reduction.
385 */
386 bigint *bi_divide(BI_CTX *ctx, bigint *u, bigint *v, int is_mod)
387 {
388 int n = v->size, m = u->size-n;
389 int j = 0, orig_u_size = u->size;
390 uint8_t mod_offset = ctx->mod_offset;
391 comp d;
392 bigint *quotient, *tmp_u;
393 comp q_dash;
394
395 check(u);
396 check(v);
397
398 /* if doing reduction and we are < mod, then return mod */
399 if (is_mod && bi_compare(v, u) > 0)
400 {
401 bi_free(ctx, v);
402 return u;
403 }
404
405 quotient = alloc(ctx, m+1);
406 tmp_u = alloc(ctx, n+1);
407 v = trim(v); /* make sure we have no leading 0's */
408 d = (comp)((long_comp)COMP_RADIX/(V1+1));
409
410 /* clear things to start with */
411 memset(quotient->comps, 0, ((quotient->size)*COMP_BYTE_SIZE));
412
413 /* normalise */
414 if (d > 1)
415 {
416 u = bi_int_multiply(ctx, u, d);
417
418 if (is_mod)
419 {
420 v = ctx->bi_normalised_mod[mod_offset];
421 }
422 else
423 {
424 v = bi_int_multiply(ctx, v, d);
425 }
426 }
427
428 if (orig_u_size == u->size) /* new digit position u0 */
429 {
430 more_comps(u, orig_u_size + 1);
431 }
432
433 do
434 {
435 /* get a temporary short version of u */
436 memcpy(tmp_u->comps, &u->comps[u->size-n-1-j], (n+1)*COMP_BYTE_SIZE);
437
438 /* calculate q' */
439 if (U(0) == V1)
440 {
441 q_dash = COMP_RADIX-1;
442 }
443 else
444 {
445 q_dash = (comp)(((long_comp)U(0)*COMP_RADIX + U(1))/V1);
446
447 if (v->size > 1 && V2)
448 {
449 /* we are implementing the following:
450 if (V2*q_dash > (((U(0)*COMP_RADIX + U(1) -
451 q_dash*V1)*COMP_RADIX) + U(2))) ... */
452 comp inner = (comp)((long_comp)COMP_RADIX*U(0) + U(1) -
453 (long_comp)q_dash*V1);
454 if ((long_comp)V2*q_dash > (long_comp)inner*COMP_RADIX + U(2))
455 {
456 q_dash--;
457 }
458 }
459 }
460
461 /* multiply and subtract */
462 if (q_dash)
463 {
464 int is_negative;
465 tmp_u = bi_subtract(ctx, tmp_u,
466 bi_int_multiply(ctx, bi_copy(v), q_dash), &is_negative);
467 more_comps(tmp_u, n+1);
468
469 Q(j) = q_dash;
470
471 /* add back */
472 if (is_negative)
473 {
474 Q(j)--;
475 tmp_u = bi_add(ctx, tmp_u, bi_copy(v));
476
477 /* lop off the carry */
478 tmp_u->size--;
479 v->size--;
480 }
481 }
482 else
483 {
484 Q(j) = 0;
485 }
486
487 /* copy back to u */
488 memcpy(&u->comps[u->size-n-1-j], tmp_u->comps, (n+1)*COMP_BYTE_SIZE);
489 } while (++j <= m);
490
491 bi_free(ctx, tmp_u);
492 bi_free(ctx, v);
493
494 if (is_mod) /* get the remainder */
495 {
496 bi_free(ctx, quotient);
497 return bi_int_divide(ctx, trim(u), d);
498 }
499 else /* get the quotient */
500 {
501 bi_free(ctx, u);
502 return trim(quotient);
503 }
504 }
505
506 /*
507 * Perform an integer divide on a bigint.
508 */
509 static bigint *bi_int_divide(BI_CTX *ctx, bigint *biR, comp denom)
510 {
511 int i = biR->size - 1;
512 long_comp r = 0;
513
514 check(biR);
515
516 do
517 {
518 r = (r<<COMP_BIT_SIZE) + biR->comps[i];
519 biR->comps[i] = (comp)(r / denom);
520 r %= denom;
521 } while (--i >= 0);
522
523 return trim(biR);
524 }
525
526 #ifdef CONFIG_BIGINT_MONTGOMERY
527 /**
528 * There is a need for the value of integer N' such that B^-1(B-1)-N^-1N'=1,
529 * where B^-1(B-1) mod N=1. Actually, only the least significant part of
530 * N' is needed, hence the definition N0'=N' mod b. We reproduce below the
531 * simple algorithm from an article by Dusse and Kaliski to efficiently
532 * find N0' from N0 and b */
533 static comp modular_inverse(bigint *bim)
534 {
535 int i;
536 comp t = 1;
537 comp two_2_i_minus_1 = 2; /* 2^(i-1) */
538 long_comp two_2_i = 4; /* 2^i */
539 comp N = bim->comps[0];
540
541 for (i = 2; i <= COMP_BIT_SIZE; i++)
542 {
543 if ((long_comp)N*t%two_2_i >= two_2_i_minus_1)
544 {
545 t += two_2_i_minus_1;
546 }
547
548 two_2_i_minus_1 <<= 1;
549 two_2_i <<= 1;
550 }
551
552 return (comp)(COMP_RADIX-t);
553 }
554 #endif
555
556 #if defined(CONFIG_BIGINT_KARATSUBA) || defined(CONFIG_BIGINT_BARRETT) || \
557 defined(CONFIG_BIGINT_MONTGOMERY)
558 /**
559 * Take each component and shift down (in terms of components)
560 */
561 static bigint *comp_right_shift(bigint *biR, int num_shifts)
562 {
563 int i = biR->size-num_shifts;
564 comp *x = biR->comps;
565 comp *y = &biR->comps[num_shifts];
566
567 check(biR);
568
569 if (i <= 0) /* have we completely right shifted? */
570 {
571 biR->comps[0] = 0; /* return 0 */
572 biR->size = 1;
573 return biR;
574 }
575
576 do
577 {
578 *x++ = *y++;
579 } while (--i > 0);
580
581 biR->size -= num_shifts;
582 return biR;
583 }
584
585 /**
586 * Take each component and shift it up (in terms of components)
587 */
588 static bigint *comp_left_shift(bigint *biR, int num_shifts)
589 {
590 int i = biR->size-1;
591 comp *x, *y;
592
593 check(biR);
594
595 if (num_shifts <= 0)
596 {
597 return biR;
598 }
599
600 more_comps(biR, biR->size + num_shifts);
601
602 x = &biR->comps[i+num_shifts];
603 y = &biR->comps[i];
604
605 do
606 {
607 *x-- = *y--;
608 } while (i--);
609
610 memset(biR->comps, 0, num_shifts*COMP_BYTE_SIZE); /* zero LS comps */
611 return biR;
612 }
613 #endif
614
615 /**
616 * @brief Allow a binary sequence to be imported as a bigint.
617 * @param ctx [in] The bigint session context.
618 * @param data [in] The data to be converted.
619 * @param size [in] The number of bytes of data.
620 * @return A bigint representing this data.
621 */
622 bigint *bi_import(BI_CTX *ctx, const uint8_t *data, int size)
623 {
624 bigint *biR = alloc(ctx, (size+COMP_BYTE_SIZE-1)/COMP_BYTE_SIZE);
625 int i, j = 0, offset = 0;
626
627 memset(biR->comps, 0, biR->size*COMP_BYTE_SIZE);
628
629 for (i = size-1; i >= 0; i--)
630 {
631 biR->comps[offset] += data[i] << (j*8);
632
633 if (++j == COMP_BYTE_SIZE)
634 {
635 j = 0;
636 offset ++;
637 }
638 }
639
640 return trim(biR);
641 }
642
643 #ifdef CONFIG_SSL_FULL_MODE
644 /**
645 * @brief The testharness uses this code to import text hex-streams and
646 * convert them into bigints.
647 * @param ctx [in] The bigint session context.
648 * @param data [in] A string consisting of hex characters. The characters must
649 * be in upper case.
650 * @return A bigint representing this data.
651 */
652 bigint *bi_str_import(BI_CTX *ctx, const char *data)
653 {
654 int size = strlen(data);
655 bigint *biR = alloc(ctx, (size+COMP_NUM_NIBBLES-1)/COMP_NUM_NIBBLES);
656 int i, j = 0, offset = 0;
657 memset(biR->comps, 0, biR->size*COMP_BYTE_SIZE);
658
659 for (i = size-1; i >= 0; i--)
660 {
661 int num = (data[i] <= '9') ? (data[i] - '0') : (data[i] - 'A' + 10);
662 biR->comps[offset] += num << (j*4);
663
664 if (++j == COMP_NUM_NIBBLES)
665 {
666 j = 0;
667 offset ++;
668 }
669 }
670
671 return biR;
672 }
673
674 void bi_print(const char *label, bigint *x)
675 {
676 int i, j;
677
678 if (x == NULL)
679 {
680 printf("%s: (null)\n", label);
681 return;
682 }
683
684 printf("%s: (size %d)\n", label, x->size);
685 for (i = x->size-1; i >= 0; i--)
686 {
687 for (j = COMP_NUM_NIBBLES-1; j >= 0; j--)
688 {
689 comp mask = 0x0f << (j*4);
690 comp num = (x->comps[i] & mask) >> (j*4);
691 putc((num <= 9) ? (num + '0') : (num + 'A' - 10), stdout);
692 }
693 }
694
695 printf("\n");
696 }
697 #endif
698
699 /**
700 * @brief Take a bigint and convert it into a byte sequence.
701 *
702 * This is useful after a decrypt operation.
703 * @param ctx [in] The bigint session context.
704 * @param x [in] The bigint to be converted.
705 * @param data [out] The converted data as a byte stream.
706 * @param size [in] The maximum size of the byte stream. Unused bytes will be
707 * zeroed.
708 */
709 void bi_export(BI_CTX *ctx, bigint *x, uint8_t *data, int size)
710 {
711 int i, j, k = size-1;
712
713 check(x);
714 memset(data, 0, size); /* ensure all leading 0's are cleared */
715
716 for (i = 0; i < x->size; i++)
717 {
718 for (j = 0; j < COMP_BYTE_SIZE; j++)
719 {
720 comp mask = 0xff << (j*8);
721 int num = (x->comps[i] & mask) >> (j*8);
722 data[k--] = num;
723
724 if (k < 0)
725 {
726 goto buf_done;
727 }
728 }
729 }
730 buf_done:
731
732 bi_free(ctx, x);
733 }
734
735 /**
736 * @brief Pre-calculate some of the expensive steps in reduction.
737 *
738 * This function should only be called once (normally when a session starts).
739 * When the session is over, bi_free_mod() should be called. bi_mod_power()
740 * relies on this function being called.
741 * @param ctx [in] The bigint session context.
742 * @param bim [in] The bigint modulus that will be used.
743 * @param mod_offset [in] There are three moduluii that can be stored - the
744 * standard modulus, and its two primes p and q. This offset refers to which
745 * modulus we are referring to.
746 * @see bi_free_mod(), bi_mod_power().
747 */
748 void bi_set_mod(BI_CTX *ctx, bigint *bim, int mod_offset)
749 {
750 int k = bim->size;
751 comp d = (comp)((long_comp)COMP_RADIX/(bim->comps[k-1]+1));
752 #ifdef CONFIG_BIGINT_MONTGOMERY
753 bigint *R, *R2;
754 #endif
755
756 ctx->bi_mod[mod_offset] = bim;
757 bi_permanent(ctx->bi_mod[mod_offset]);
758 ctx->bi_normalised_mod[mod_offset] = bi_int_multiply(ctx, bim, d);
759 bi_permanent(ctx->bi_normalised_mod[mod_offset]);
760
761 #if defined(CONFIG_BIGINT_MONTGOMERY)
762 /* set montgomery variables */
763 R = comp_left_shift(bi_clone(ctx, ctx->bi_radix), k-1); /* R */
764 R2 = comp_left_shift(bi_clone(ctx, ctx->bi_radix), k*2-1); /* R^2 */
765 ctx->bi_RR_mod_m[mod_offset] = bi_mod(ctx, R2); /* R^2 mod m */
766 ctx->bi_R_mod_m[mod_offset] = bi_mod(ctx, R); /* R mod m */
767
768 bi_permanent(ctx->bi_RR_mod_m[mod_offset]);
769 bi_permanent(ctx->bi_R_mod_m[mod_offset]);
770
771 ctx->N0_dash[mod_offset] = modular_inverse(ctx->bi_mod[mod_offset]);
772
773 #elif defined (CONFIG_BIGINT_BARRETT)
774 ctx->bi_mu[mod_offset] =
775 bi_divide(ctx, comp_left_shift(
776 bi_clone(ctx, ctx->bi_radix), k*2-1), ctx->bi_mod[mod_offset], 0);
777 bi_permanent(ctx->bi_mu[mod_offset]);
778 #endif
779 }
780
781 /**
782 * @brief Used when cleaning various bigints at the end of a session.
783 * @param ctx [in] The bigint session context.
784 * @param mod_offset [in] The offset to use.
785 * @see bi_set_mod().
786 */
787 void bi_free_mod(BI_CTX *ctx, int mod_offset)
788 {
789 bi_depermanent(ctx->bi_mod[mod_offset]);
790 bi_free(ctx, ctx->bi_mod[mod_offset]);
791 #if defined (CONFIG_BIGINT_MONTGOMERY)
792 bi_depermanent(ctx->bi_RR_mod_m[mod_offset]);
793 bi_depermanent(ctx->bi_R_mod_m[mod_offset]);
794 bi_free(ctx, ctx->bi_RR_mod_m[mod_offset]);
795 bi_free(ctx, ctx->bi_R_mod_m[mod_offset]);
796 #elif defined(CONFIG_BIGINT_BARRETT)
797 bi_depermanent(ctx->bi_mu[mod_offset]);
798 bi_free(ctx, ctx->bi_mu[mod_offset]);
799 #endif
800 bi_depermanent(ctx->bi_normalised_mod[mod_offset]);
801 bi_free(ctx, ctx->bi_normalised_mod[mod_offset]);
802 }
803
804 /**
805 * Perform a standard multiplication between two bigints.
806 *
807 * Barrett reduction has no need for some parts of the product, so ignore bits
808 * of the multiply. This routine gives Barrett its big performance
809 * improvements over Classical/Montgomery reduction methods.
810 */
811 static bigint *regular_multiply(BI_CTX *ctx, bigint *bia, bigint *bib,
812 int inner_partial, int outer_partial)
813 {
814 int i = 0, j;
815 int n = bia->size;
816 int t = bib->size;
817 bigint *biR = alloc(ctx, n + t);
818 comp *sr = biR->comps;
819 comp *sa = bia->comps;
820 comp *sb = bib->comps;
821
822 check(bia);
823 check(bib);
824
825 /* clear things to start with */
826 memset(biR->comps, 0, ((n+t)*COMP_BYTE_SIZE));
827
828 do
829 {
830 long_comp tmp;
831 comp carry = 0;
832 int r_index = i;
833 j = 0;
834
835 if (outer_partial && outer_partial-i > 0 && outer_partial < n)
836 {
837 r_index = outer_partial-1;
838 j = outer_partial-i-1;
839 }
840
841 do
842 {
843 if (inner_partial && r_index >= inner_partial)
844 {
845 break;
846 }
847
848 tmp = sr[r_index] + ((long_comp)sa[j])*sb[i] + carry;
849 sr[r_index++] = (comp)tmp; /* downsize */
850 carry = tmp >> COMP_BIT_SIZE;
851 } while (++j < n);
852
853 sr[r_index] = carry;
854 } while (++i < t);
855
856 bi_free(ctx, bia);
857 bi_free(ctx, bib);
858 return trim(biR);
859 }
860
861 #ifdef CONFIG_BIGINT_KARATSUBA
862 /*
863 * Karatsuba improves on regular multiplication due to only 3 multiplications
864 * being done instead of 4. The additional additions/subtractions are O(N)
865 * rather than O(N^2) and so for big numbers it saves on a few operations
866 */
867 static bigint *karatsuba(BI_CTX *ctx, bigint *bia, bigint *bib, int is_square)
868 {
869 bigint *x0, *x1;
870 bigint *p0, *p1, *p2;
871 int m;
872
873 if (is_square)
874 {
875 m = (bia->size + 1)/2;
876 }
877 else
878 {
879 m = (max(bia->size, bib->size) + 1)/2;
880 }
881
882 x0 = bi_clone(ctx, bia);
883 x0->size = m;
884 x1 = bi_clone(ctx, bia);
885 comp_right_shift(x1, m);
886 bi_free(ctx, bia);
887
888 /* work out the 3 partial products */
889 if (is_square)
890 {
891 p0 = bi_square(ctx, bi_copy(x0));
892 p2 = bi_square(ctx, bi_copy(x1));
893 p1 = bi_square(ctx, bi_add(ctx, x0, x1));
894 }
895 else /* normal multiply */
896 {
897 bigint *y0, *y1;
898 y0 = bi_clone(ctx, bib);
899 y0->size = m;
900 y1 = bi_clone(ctx, bib);
901 comp_right_shift(y1, m);
902 bi_free(ctx, bib);
903
904 p0 = bi_multiply(ctx, bi_copy(x0), bi_copy(y0));
905 p2 = bi_multiply(ctx, bi_copy(x1), bi_copy(y1));
906 p1 = bi_multiply(ctx, bi_add(ctx, x0, x1), bi_add(ctx, y0, y1));
907 }
908
909 p1 = bi_subtract(ctx,
910 bi_subtract(ctx, p1, bi_copy(p2), NULL), bi_copy(p0), NULL);
911
912 comp_left_shift(p1, m);
913 comp_left_shift(p2, 2*m);
914 return bi_add(ctx, p1, bi_add(ctx, p0, p2));
915 }
916 #endif
917
918 /**
919 * @brief Perform a multiplication operation between two bigints.
920 * @param ctx [in] The bigint session context.
921 * @param bia [in] A bigint.
922 * @param bib [in] Another bigint.
923 * @return The result of the multiplication.
924 */
925 bigint *bi_multiply(BI_CTX *ctx, bigint *bia, bigint *bib)
926 {
927 check(bia);
928 check(bib);
929
930 #ifdef CONFIG_BIGINT_KARATSUBA
931 if (min(bia->size, bib->size) < MUL_KARATSUBA_THRESH)
932 {
933 return regular_multiply(ctx, bia, bib, 0, 0);
934 }
935
936 return karatsuba(ctx, bia, bib, 0);
937 #else
938 return regular_multiply(ctx, bia, bib, 0, 0);
939 #endif
940 }
941
942 #ifdef CONFIG_BIGINT_SQUARE
943 /*
944 * Perform the actual square operion. It takes into account overflow.
945 */
946 static bigint *regular_square(BI_CTX *ctx, bigint *bi)
947 {
948 int t = bi->size;
949 int i = 0, j;
950 bigint *biR = alloc(ctx, t*2+1);
951 comp *w = biR->comps;
952 comp *x = bi->comps;
953 long_comp carry;
954 memset(w, 0, biR->size*COMP_BYTE_SIZE);
955
956 do
957 {
958 long_comp tmp = w[2*i] + (long_comp)x[i]*x[i];
959 w[2*i] = (comp)tmp;
960 carry = tmp >> COMP_BIT_SIZE;
961
962 for (j = i+1; j < t; j++)
963 {
964 uint8_t c = 0;
965 long_comp xx = (long_comp)x[i]*x[j];
966 if ((COMP_MAX-xx) < xx)
967 c = 1;
968
969 tmp = (xx<<1);
970
971 if ((COMP_MAX-tmp) < w[i+j])
972 c = 1;
973
974 tmp += w[i+j];
975
976 if ((COMP_MAX-tmp) < carry)
977 c = 1;
978
979 tmp += carry;
980 w[i+j] = (comp)tmp;
981 carry = tmp >> COMP_BIT_SIZE;
982
983 if (c)
984 carry += COMP_RADIX;
985 }
986
987 tmp = w[i+t] + carry;
988 w[i+t] = (comp)tmp;
989 w[i+t+1] = tmp >> COMP_BIT_SIZE;
990 } while (++i < t);
991
992 bi_free(ctx, bi);
993 return trim(biR);
994 }
995
996 /**
997 * @brief Perform a square operation on a bigint.
998 * @param ctx [in] The bigint session context.
999 * @param bia [in] A bigint.
1000 * @return The result of the multiplication.
1001 */
1002 bigint *bi_square(BI_CTX *ctx, bigint *bia)
1003 {
1004 check(bia);
1005
1006 #ifdef CONFIG_BIGINT_KARATSUBA
1007 if (bia->size < SQU_KARATSUBA_THRESH)
1008 {
1009 return regular_square(ctx, bia);
1010 }
1011
1012 return karatsuba(ctx, bia, NULL, 1);
1013 #else
1014 return regular_square(ctx, bia);
1015 #endif
1016 }
1017 #endif
1018
1019 /**
1020 * @brief Compare two bigints.
1021 * @param bia [in] A bigint.
1022 * @param bib [in] Another bigint.
1023 * @return -1 if smaller, 1 if larger and 0 if equal.
1024 */
1025 int bi_compare(bigint *bia, bigint *bib)
1026 {
1027 int r, i;
1028
1029 check(bia);
1030 check(bib);
1031
1032 if (bia->size > bib->size)
1033 r = 1;
1034 else if (bia->size < bib->size)
1035 r = -1;
1036 else
1037 {
1038 comp *a = bia->comps;
1039 comp *b = bib->comps;
1040
1041 /* Same number of components. Compare starting from the high end
1042 * and working down. */
1043 r = 0;
1044 i = bia->size - 1;
1045
1046 do
1047 {
1048 if (a[i] > b[i])
1049 {
1050 r = 1;
1051 break;
1052 }
1053 else if (a[i] < b[i])
1054 {
1055 r = -1;
1056 break;
1057 }
1058 } while (--i >= 0);
1059 }
1060
1061 return r;
1062 }
1063
1064 /*
1065 * Allocate and zero more components. Does not consume bi.
1066 */
1067 static void more_comps(bigint *bi, int n)
1068 {
1069 if (n > bi->max_comps)
1070 {
1071 bi->max_comps = max(bi->max_comps * 2, n);
1072 bi->comps = (comp*)realloc(bi->comps, bi->max_comps * COMP_BYTE_SIZE);
1073 }
1074
1075 if (n > bi->size)
1076 {
1077 memset(&bi->comps[bi->size], 0, (n-bi->size)*COMP_BYTE_SIZE);
1078 }
1079
1080 bi->size = n;
1081 }
1082
1083 /*
1084 * Make a new empty bigint. It may just use an old one if one is available.
1085 * Otherwise get one off the heap.
1086 */
1087 static bigint *alloc(BI_CTX *ctx, int size)
1088 {
1089 bigint *biR;
1090
1091 /* Can we recycle an old bigint? */
1092 if (ctx->free_list != NULL)
1093 {
1094 biR = ctx->free_list;
1095 ctx->free_list = biR->next;
1096 ctx->free_count--;
1097
1098 if (biR->refs != 0)
1099 {
1100 #ifdef CONFIG_SSL_FULL_MODE
1101 printf("alloc: refs was not 0\n");
1102 #endif
1103 abort(); /* create a stack trace from a core dump */
1104 }
1105
1106 more_comps(biR, size);
1107 }
1108 else
1109 {
1110 /* No free bigints available - create a new one. */
1111 biR = (bigint *)malloc(sizeof(bigint));
1112 biR->comps = (comp*)malloc(size * COMP_BYTE_SIZE);
1113 biR->max_comps = size; /* give some space to spare */
1114 }
1115
1116 biR->size = size;
1117 biR->refs = 1;
1118 biR->next = NULL;
1119 ctx->active_count++;
1120 return biR;
1121 }
1122
1123 /*
1124 * Work out the highest '1' bit in an exponent. Used when doing sliding-window
1125 * exponentiation.
1126 */
1127 static int find_max_exp_index(bigint *biexp)
1128 {
1129 int i = COMP_BIT_SIZE-1;
1130 comp shift = COMP_RADIX/2;
1131 comp test = biexp->comps[biexp->size-1]; /* assume no leading zeroes */
1132
1133 check(biexp);
1134
1135 do
1136 {
1137 if (test & shift)
1138 {
1139 return i+(biexp->size-1)*COMP_BIT_SIZE;
1140 }
1141
1142 shift >>= 1;
1143 } while (i-- != 0);
1144
1145 return -1; /* error - must have been a leading 0 */
1146 }
1147
1148 /*
1149 * Is a particular bit is an exponent 1 or 0? Used when doing sliding-window
1150 * exponentiation.
1151 */
1152 static int exp_bit_is_one(bigint *biexp, int offset)
1153 {
1154 comp test = biexp->comps[offset / COMP_BIT_SIZE];
1155 int num_shifts = offset % COMP_BIT_SIZE;
1156 comp shift = 1;
1157 int i;
1158
1159 check(biexp);
1160
1161 for (i = 0; i < num_shifts; i++)
1162 {
1163 shift <<= 1;
1164 }
1165
1166 return (test & shift) != 0;
1167 }
1168
1169 #ifdef CONFIG_BIGINT_CHECK_ON
1170 /*
1171 * Perform a sanity check on bi.
1172 */
1173 static void check(const bigint *bi)
1174 {
1175 if (bi->refs <= 0)
1176 {
1177 printf("check: zero or negative refs in bigint\n");
1178 abort();
1179 }
1180
1181 if (bi->next != NULL)
1182 {
1183 printf("check: attempt to use a bigint from "
1184 "the free list\n");
1185 abort();
1186 }
1187 }
1188 #endif
1189
1190 /*
1191 * Delete any leading 0's (and allow for 0).
1192 */
1193 static bigint *trim(bigint *bi)
1194 {
1195 check(bi);
1196
1197 while (bi->comps[bi->size-1] == 0 && bi->size > 1)
1198 {
1199 bi->size--;
1200 }
1201
1202 return bi;
1203 }
1204
1205 #if defined(CONFIG_BIGINT_MONTGOMERY)
1206 /**
1207 * @brief Perform a single montgomery reduction.
1208 * @param ctx [in] The bigint session context.
1209 * @param bixy [in] A bigint.
1210 * @return The result of the montgomery reduction.
1211 */
1212 bigint *bi_mont(BI_CTX *ctx, bigint *bixy)
1213 {
1214 int i = 0, n;
1215 uint8_t mod_offset = ctx->mod_offset;
1216 bigint *bim = ctx->bi_mod[mod_offset];
1217 comp mod_inv = ctx->N0_dash[mod_offset];
1218
1219 check(bixy);
1220
1221 if (ctx->use_classical) /* just use classical instead */
1222 {
1223 return bi_mod(ctx, bixy);
1224 }
1225
1226 n = bim->size;
1227
1228 do
1229 {
1230 bixy = bi_add(ctx, bixy, comp_left_shift(
1231 bi_int_multiply(ctx, bim, bixy->comps[i]*mod_inv), i));
1232 } while (++i < n);
1233
1234 comp_right_shift(bixy, n);
1235
1236 if (bi_compare(bixy, bim) >= 0)
1237 {
1238 bixy = bi_subtract(ctx, bixy, bim, NULL);
1239 }
1240
1241 return bixy;
1242 }
1243
1244 #elif defined(CONFIG_BIGINT_BARRETT)
1245 /*
1246 * Stomp on the most significant components to give the illusion of a "mod base
1247 * radix" operation
1248 */
1249 static bigint *comp_mod(bigint *bi, int mod)
1250 {
1251 check(bi);
1252
1253 if (bi->size > mod)
1254 {
1255 bi->size = mod;
1256 }
1257
1258 return bi;
1259 }
1260
1261 /**
1262 * @brief Perform a single Barrett reduction.
1263 * @param ctx [in] The bigint session context.
1264 * @param bi [in] A bigint.
1265 * @return The result of the Barrett reduction.
1266 */
1267 bigint *bi_barrett(BI_CTX *ctx, bigint *bi)
1268 {
1269 bigint *q1, *q2, *q3, *r1, *r2, *r;
1270 uint8_t mod_offset = ctx->mod_offset;
1271 bigint *bim = ctx->bi_mod[mod_offset];
1272 int k = bim->size;
1273
1274 check(bi);
1275 check(bim);
1276
1277 /* use Classical method instead - Barrett cannot help here */
1278 if (bi->size > k*2)
1279 {
1280 return bi_mod(ctx, bi);
1281 }
1282
1283 q1 = comp_right_shift(bi_clone(ctx, bi), k-1);
1284
1285 /* do outer partial multiply */
1286 q2 = regular_multiply(ctx, q1, ctx->bi_mu[mod_offset], 0, k-1);
1287 q3 = comp_right_shift(q2, k+1);
1288 r1 = comp_mod(bi, k+1);
1289
1290 /* do inner partial multiply */
1291 r2 = comp_mod(regular_multiply(ctx, q3, bim, k+1, 0), k+1);
1292 r = bi_subtract(ctx, r1, r2, NULL);
1293
1294 /* if (r >= m) r = r - m; */
1295 if (bi_compare(r, bim) >= 0)
1296 {
1297 r = bi_subtract(ctx, r, bim, NULL);
1298 }
1299
1300 return r;
1301 }
1302 #endif /* CONFIG_BIGINT_BARRETT */
1303
1304 #ifdef CONFIG_BIGINT_SLIDING_WINDOW
1305 /*
1306 * Work out g1, g3, g5, g7... etc for the sliding-window algorithm
1307 */
1308 static void precompute_slide_window(BI_CTX *ctx, int window, bigint *g1)
1309 {
1310 int k = 1, i;
1311 bigint *g2;
1312
1313 for (i = 0; i < window-1; i++) /* compute 2^(window-1) */
1314 {
1315 k <<= 1;
1316 }
1317
1318 ctx->g = (bigint **)malloc(k*sizeof(bigint *));
1319 ctx->g[0] = bi_clone(ctx, g1);
1320 bi_permanent(ctx->g[0]);
1321 g2 = bi_residue(ctx, bi_square(ctx, ctx->g[0])); /* g^2 */
1322
1323 for (i = 1; i < k; i++)
1324 {
1325 ctx->g[i] = bi_residue(ctx, bi_multiply(ctx, ctx->g[i-1], bi_copy(g2)));
1326 bi_permanent(ctx->g[i]);
1327 }
1328
1329 bi_free(ctx, g2);
1330 ctx->window = k;
1331 }
1332 #endif
1333
1334 /**
1335 * @brief Perform a modular exponentiation.
1336 *
1337 * This function requires bi_set_mod() to have been called previously. This is
1338 * one of the optimisations used for performance.
1339 * @param ctx [in] The bigint session context.
1340 * @param bi [in] The bigint on which to perform the mod power operation.
1341 * @param biexp [in] The bigint exponent.
1342 * @return The result of the mod exponentiation operation
1343 * @see bi_set_mod().
1344 */
1345 bigint *bi_mod_power(BI_CTX *ctx, bigint *bi, bigint *biexp)
1346 {
1347 int i = find_max_exp_index(biexp), j, window_size = 1;
1348 bigint *biR = int_to_bi(ctx, 1);
1349
1350 #if defined(CONFIG_BIGINT_MONTGOMERY)
1351 uint8_t mod_offset = ctx->mod_offset;
1352 if (!ctx->use_classical)
1353 {
1354 /* preconvert */
1355 bi = bi_mont(ctx,
1356 bi_multiply(ctx, bi, ctx->bi_RR_mod_m[mod_offset])); /* x' */
1357 bi_free(ctx, biR);
1358 biR = ctx->bi_R_mod_m[mod_offset]; /* A */
1359 }
1360 #endif
1361
1362 check(bi);
1363 check(biexp);
1364
1365 #ifdef CONFIG_BIGINT_SLIDING_WINDOW
1366 for (j = i; j > 32; j /= 5) /* work out an optimum size */
1367 window_size++;
1368
1369 /* work out the slide constants */
1370 precompute_slide_window(ctx, window_size, bi);
1371 #else /* just one constant */
1372 ctx->g = (bigint **)malloc(sizeof(bigint *));
1373 ctx->g[0] = bi_clone(ctx, bi);
1374 ctx->window = 1;
1375 bi_permanent(ctx->g[0]);
1376 #endif
1377
1378 /* if sliding-window is off, then only one bit will be done at a time and
1379 * will reduce to standard left-to-right exponentiation */
1380 do
1381 {
1382 if (exp_bit_is_one(biexp, i))
1383 {
1384 int l = i-window_size+1;
1385 int part_exp = 0;
1386
1387 if (l < 0) /* LSB of exponent will always be 1 */
1388 l = 0;
1389 else
1390 {
1391 while (exp_bit_is_one(biexp, l) == 0)
1392 l++; /* go back up */
1393 }
1394
1395 /* build up the section of the exponent */
1396 for (j = i; j >= l; j--)
1397 {
1398 biR = bi_residue(ctx, bi_square(ctx, biR));
1399 if (exp_bit_is_one(biexp, j))
1400 part_exp++;
1401
1402 if (j != l)
1403 part_exp <<= 1;
1404 }
1405
1406 part_exp = (part_exp-1)/2; /* adjust for array */
1407 biR = bi_residue(ctx, bi_multiply(ctx, biR, ctx->g[part_exp]));
1408 i = l-1;
1409 }
1410 else /* square it */
1411 {
1412 biR = bi_residue(ctx, bi_square(ctx, biR));
1413 i--;
1414 }
1415 } while (i >= 0);
1416
1417 /* cleanup */
1418 for (i = 0; i < ctx->window; i++)
1419 {
1420 bi_depermanent(ctx->g[i]);
1421 bi_free(ctx, ctx->g[i]);
1422 }
1423
1424 free(ctx->g);
1425 bi_free(ctx, bi);
1426 bi_free(ctx, biexp);
1427 #if defined CONFIG_BIGINT_MONTGOMERY
1428 return ctx->use_classical ? biR : bi_mont(ctx, biR); /* convert back */
1429 #else /* CONFIG_BIGINT_CLASSICAL or CONFIG_BIGINT_BARRETT */
1430 return biR;
1431 #endif
1432 }
1433
1434 #ifdef CONFIG_SSL_CERT_VERIFICATION
1435 /**
1436 * @brief Perform a modular exponentiation using a temporary modulus.
1437 *
1438 * We need this function to check the signatures of certificates. The modulus
1439 * of this function is temporary as it's just used for authentication.
1440 * @param ctx [in] The bigint session context.
1441 * @param bi [in] The bigint to perform the exp/mod.
1442 * @param bim [in] The temporary modulus.
1443 * @param biexp [in] The bigint exponent.
1444 * @return The result of the mod exponentiation operation
1445 * @see bi_set_mod().
1446 */
1447 bigint *bi_mod_power2(BI_CTX *ctx, bigint *bi, bigint *bim, bigint *biexp)
1448 {
1449 bigint *biR, *tmp_biR;
1450
1451 /* Set up a temporary bigint context and transfer what we need between
1452 * them. We need to do this since we want to keep the original modulus
1453 * which is already in this context. This operation is only called when
1454 * doing peer verification, and so is not expensive :-) */
1455 BI_CTX *tmp_ctx = bi_initialize();
1456 bi_set_mod(tmp_ctx, bi_clone(tmp_ctx, bim), BIGINT_M_OFFSET);
1457 tmp_biR = bi_mod_power(tmp_ctx,
1458 bi_clone(tmp_ctx, bi),
1459 bi_clone(tmp_ctx, biexp));
1460 biR = bi_clone(ctx, tmp_biR);
1461 bi_free(tmp_ctx, tmp_biR);
1462 bi_free_mod(tmp_ctx, BIGINT_M_OFFSET);
1463 bi_terminate(tmp_ctx);
1464
1465 bi_free(ctx, bi);
1466 bi_free(ctx, bim);
1467 bi_free(ctx, biexp);
1468 return biR;
1469 }
1470 #endif
1471
1472 #ifdef CONFIG_BIGINT_CRT
1473 /**
1474 * @brief Use the Chinese Remainder Theorem to quickly perform RSA decrypts.
1475 *
1476 * @param ctx [in] The bigint session context.
1477 * @param bi [in] The bigint to perform the exp/mod.
1478 * @param dP [in] CRT's dP bigint
1479 * @param dQ [in] CRT's dQ bigint
1480 * @param p [in] CRT's p bigint
1481 * @param q [in] CRT's q bigint
1482 * @param qInv [in] CRT's qInv bigint
1483 * @return The result of the CRT operation
1484 */
1485 bigint *bi_crt(BI_CTX *ctx, bigint *bi,
1486 bigint *dP, bigint *dQ,
1487 bigint *p, bigint *q, bigint *qInv)
1488 {
1489 bigint *m1, *m2, *h;
1490
1491 /* Montgomery has a condition the 0 < x, y < m and these products violate
1492 * that condition. So disable Montgomery when using CRT */
1493 #if defined(CONFIG_BIGINT_MONTGOMERY)
1494 ctx->use_classical = 1;
1495 #endif
1496 ctx->mod_offset = BIGINT_P_OFFSET;
1497 m1 = bi_mod_power(ctx, bi_copy(bi), dP);
1498
1499 ctx->mod_offset = BIGINT_Q_OFFSET;
1500 m2 = bi_mod_power(ctx, bi, dQ);
1501
1502 h = bi_subtract(ctx, bi_add(ctx, m1, p), bi_copy(m2), NULL);
1503 h = bi_multiply(ctx, h, qInv);
1504 ctx->mod_offset = BIGINT_P_OFFSET;
1505 h = bi_residue(ctx, h);
1506 #if defined(CONFIG_BIGINT_MONTGOMERY)
1507 ctx->use_classical = 0; /* reset for any further operation */
1508 #endif
1509 return bi_add(ctx, m2, bi_multiply(ctx, q, h));
1510 }
1511 #endif
1512 /** @} */