mpi.c
Go to the documentation of this file.
1 /**
2  * @file mpi.c
3  * @brief MPI (Multiple Precision Integer Arithmetic)
4  *
5  * @section License
6  *
7  * SPDX-License-Identifier: GPL-2.0-or-later
8  *
9  * Copyright (C) 2010-2023 Oryx Embedded SARL. All rights reserved.
10  *
11  * This file is part of CycloneCRYPTO Open.
12  *
13  * This program is free software; you can redistribute it and/or
14  * modify it under the terms of the GNU General Public License
15  * as published by the Free Software Foundation; either version 2
16  * of the License, or (at your option) any later version.
17  *
18  * This program is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with this program; if not, write to the Free Software Foundation,
25  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
26  *
27  * @author Oryx Embedded SARL (www.oryx-embedded.com)
28  * @version 2.2.4
29  **/
30 
31 //Switch to the appropriate trace level
32 #define TRACE_LEVEL CRYPTO_TRACE_LEVEL
33 
34 //Dependencies
35 #include "core/crypto.h"
36 #include "mpi/mpi.h"
37 #include "debug.h"
38 
39 //Check crypto library configuration
40 #if (MPI_SUPPORT == ENABLED)
41 
42 
43 /**
44  * @brief Initialize a multiple precision integer
45  * @param[in,out] r Pointer to the multiple precision integer to be initialized
46  **/
47 
48 void mpiInit(Mpi *r)
49 {
50  //Initialize structure
51  r->sign = 1;
52  r->size = 0;
53  r->data = NULL;
54 }
55 
56 
57 /**
58  * @brief Release a multiple precision integer
59  * @param[in,out] r Pointer to the multiple precision integer to be freed
60  **/
61 
62 void mpiFree(Mpi *r)
63 {
64  //Any memory previously allocated?
65  if(r->data != NULL)
66  {
67  //Erase contents before releasing memory
68  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
69  cryptoFreeMem(r->data);
70  }
71 
72  //Set size to zero
73  r->size = 0;
74  r->data = NULL;
75 }
76 
77 
78 /**
79  * @brief Adjust the size of multiple precision integer
80  * @param[in,out] r A multiple precision integer whose size is to be increased
81  * @param[in] size Desired size in words
82  * @return Error code
83  **/
84 
86 {
87  uint_t *data;
88 
89  //Ensure the parameter is valid
90  size = MAX(size, 1);
91 
92  //Check the current size
93  if(r->size >= size)
94  return NO_ERROR;
95 
96  //Allocate a memory buffer
98  //Failed to allocate memory?
99  if(data == NULL)
100  return ERROR_OUT_OF_MEMORY;
101 
102  //Clear buffer contents
103  osMemset(data, 0, size * MPI_INT_SIZE);
104 
105  //Any data to copy?
106  if(r->size > 0)
107  {
108  //Copy original data
109  osMemcpy(data, r->data, r->size * MPI_INT_SIZE);
110  //Free previously allocated memory
111  cryptoFreeMem(r->data);
112  }
113 
114  //Update the size of the multiple precision integer
115  r->size = size;
116  r->data = data;
117 
118  //Successful operation
119  return NO_ERROR;
120 }
121 
122 
123 /**
124  * @brief Get the actual length in words
125  * @param[in] a Pointer to a multiple precision integer
126  * @return The actual length in words
127  **/
128 
130 {
131  int_t i;
132 
133  //Check whether the specified multiple precision integer is empty
134  if(a->size == 0)
135  return 0;
136 
137  //Start from the most significant word
138  for(i = a->size - 1; i >= 0; i--)
139  {
140  //Loop as long as the current word is zero
141  if(a->data[i] != 0)
142  break;
143  }
144 
145  //Return the actual length
146  return i + 1;
147 }
148 
149 
150 /**
151  * @brief Get the actual length in bytes
152  * @param[in] a Pointer to a multiple precision integer
153  * @return The actual byte count
154  **/
155 
157 {
158  uint_t n;
159  uint32_t m;
160 
161  //Check whether the specified multiple precision integer is empty
162  if(a->size == 0)
163  return 0;
164 
165  //Start from the most significant word
166  for(n = a->size - 1; n > 0; n--)
167  {
168  //Loop as long as the current word is zero
169  if(a->data[n] != 0)
170  break;
171  }
172 
173  //Get the current word
174  m = a->data[n];
175  //Convert the length to a byte count
176  n *= MPI_INT_SIZE;
177 
178  //Adjust the byte count
179  for(; m != 0; m >>= 8)
180  {
181  n++;
182  }
183 
184  //Return the actual length in bytes
185  return n;
186 }
187 
188 
189 /**
190  * @brief Get the actual length in bits
191  * @param[in] a Pointer to a multiple precision integer
192  * @return The actual bit count
193  **/
194 
196 {
197  uint_t n;
198  uint32_t m;
199 
200  //Check whether the specified multiple precision integer is empty
201  if(a->size == 0)
202  return 0;
203 
204  //Start from the most significant word
205  for(n = a->size - 1; n > 0; n--)
206  {
207  //Loop as long as the current word is zero
208  if(a->data[n] != 0)
209  break;
210  }
211 
212  //Get the current word
213  m = a->data[n];
214  //Convert the length to a bit count
215  n *= MPI_INT_SIZE * 8;
216 
217  //Adjust the bit count
218  for(; m != 0; m >>= 1)
219  {
220  n++;
221  }
222 
223  //Return the actual length in bits
224  return n;
225 }
226 
227 
228 /**
229  * @brief Set the bit value at the specified index
230  * @param[in] r Pointer to a multiple precision integer
231  * @param[in] index Position of the bit to be written
232  * @param[in] value Bit value
233  * @return Error code
234  **/
235 
237 {
238  error_t error;
239  uint_t n1;
240  uint_t n2;
241 
242  //Retrieve the position of the bit to be written
243  n1 = index / (MPI_INT_SIZE * 8);
244  n2 = index % (MPI_INT_SIZE * 8);
245 
246  //Ajust the size of the multiple precision integer if necessary
247  error = mpiGrow(r, n1 + 1);
248  //Failed to adjust the size?
249  if(error)
250  return error;
251 
252  //Set bit value
253  if(value)
254  r->data[n1] |= (1 << n2);
255  else
256  r->data[n1] &= ~(1 << n2);
257 
258  //No error to report
259  return NO_ERROR;
260 }
261 
262 
263 /**
264  * @brief Get the bit value at the specified index
265  * @param[in] a Pointer to a multiple precision integer
266  * @param[in] index Position where to read the bit
267  * @return The actual bit value
268  **/
269 
271 {
272  uint_t n1;
273  uint_t n2;
274 
275  //Retrieve the position of the bit to be read
276  n1 = index / (MPI_INT_SIZE * 8);
277  n2 = index % (MPI_INT_SIZE * 8);
278 
279  //Index out of range?
280  if(n1 >= a->size)
281  return 0;
282 
283  //Return the actual bit value
284  return (a->data[n1] >> n2) & 0x01;
285 }
286 
287 
288 /**
289  * @brief Compare two multiple precision integers
290  * @param[in] a The first multiple precision integer to be compared
291  * @param[in] b The second multiple precision integer to be compared
292  * @return Comparison result
293  **/
294 
295 int_t mpiComp(const Mpi *a, const Mpi *b)
296 {
297  uint_t m;
298  uint_t n;
299 
300  //Determine the actual length of A and B
301  m = mpiGetLength(a);
302  n = mpiGetLength(b);
303 
304  //Compare lengths
305  if(!m && !n)
306  return 0;
307  else if(m > n)
308  return a->sign;
309  else if(m < n)
310  return -b->sign;
311 
312  //Compare signs
313  if(a->sign > 0 && b->sign < 0)
314  return 1;
315  else if(a->sign < 0 && b->sign > 0)
316  return -1;
317 
318  //Then compare values
319  while(n--)
320  {
321  if(a->data[n] > b->data[n])
322  return a->sign;
323  else if(a->data[n] < b->data[n])
324  return -a->sign;
325  }
326 
327  //Multiple precision integers are equals
328  return 0;
329 }
330 
331 
332 /**
333  * @brief Compare a multiple precision integer with an integer
334  * @param[in] a Multiple precision integer to be compared
335  * @param[in] b Integer to be compared
336  * @return Comparison result
337  **/
338 
340 {
341  uint_t value;
342  Mpi t;
343 
344  //Initialize a temporary multiple precision integer
345  value = (b >= 0) ? b : -b;
346  t.sign = (b >= 0) ? 1 : -1;
347  t.size = 1;
348  t.data = &value;
349 
350  //Return comparison result
351  return mpiComp(a, &t);
352 }
353 
354 
355 /**
356  * @brief Compare the absolute value of two multiple precision integers
357  * @param[in] a The first multiple precision integer to be compared
358  * @param[in] b The second multiple precision integer to be compared
359  * @return Comparison result
360  **/
361 
362 int_t mpiCompAbs(const Mpi *a, const Mpi *b)
363 {
364  uint_t m;
365  uint_t n;
366 
367  //Determine the actual length of A and B
368  m = mpiGetLength(a);
369  n = mpiGetLength(b);
370 
371  //Compare lengths
372  if(!m && !n)
373  return 0;
374  else if(m > n)
375  return 1;
376  else if(m < n)
377  return -1;
378 
379  //Then compare values
380  while(n--)
381  {
382  if(a->data[n] > b->data[n])
383  return 1;
384  else if(a->data[n] < b->data[n])
385  return -1;
386  }
387 
388  //Operands are equals
389  return 0;
390 }
391 
392 
393 /**
394  * @brief Copy a multiple precision integer
395  * @param[out] r Pointer to a multiple precision integer (destination)
396  * @param[in] a Pointer to a multiple precision integer (source)
397  * @return Error code
398  **/
399 
401 {
402  error_t error;
403  uint_t n;
404 
405  //R and A are the same instance?
406  if(r == a)
407  return NO_ERROR;
408 
409  //Determine the actual length of A
410  n = mpiGetLength(a);
411 
412  //Ajust the size of the destination operand
413  error = mpiGrow(r, n);
414  //Any error to report?
415  if(error)
416  return error;
417 
418  //Clear the contents of the multiple precision integer
419  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
420  //Let R = A
421  osMemcpy(r->data, a->data, n * MPI_INT_SIZE);
422  //Set the sign of R
423  r->sign = a->sign;
424 
425  //Successful operation
426  return NO_ERROR;
427 }
428 
429 
430 /**
431  * @brief Set the value of a multiple precision integer
432  * @param[out] r Pointer to a multiple precision integer
433  * @param[in] a Value to be assigned to the multiple precision integer
434  * @return Error code
435  **/
436 
438 {
439  error_t error;
440 
441  //Ajust the size of the destination operand
442  error = mpiGrow(r, 1);
443  //Failed to adjust the size?
444  if(error)
445  return error;
446 
447  //Clear the contents of the multiple precision integer
448  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
449  //Set the value or R
450  r->data[0] = (a >= 0) ? a : -a;
451  //Set the sign of R
452  r->sign = (a >= 0) ? 1 : -1;
453 
454  //Successful operation
455  return NO_ERROR;
456 }
457 
458 
459 /**
460  * @brief Generate a random value
461  * @param[out] r Pointer to a multiple precision integer
462  * @param[in] length Desired length in bits
463  * @param[in] prngAlgo PRNG algorithm
464  * @param[in] prngContext Pointer to the PRNG context
465  * @return Error code
466  **/
467 
468 error_t mpiRand(Mpi *r, uint_t length, const PrngAlgo *prngAlgo,
469  void *prngContext)
470 {
471  error_t error;
472  uint_t m;
473  uint_t n;
474 
475  //Compute the required length, in words
476  n = (length + (MPI_INT_SIZE * 8) - 1) / (MPI_INT_SIZE * 8);
477  //Number of bits in the most significant word
478  m = length % (MPI_INT_SIZE * 8);
479 
480  //Ajust the size of the multiple precision integer if necessary
481  error = mpiGrow(r, n);
482  //Failed to adjust the size?
483  if(error)
484  return error;
485 
486  //Clear the contents of the multiple precision integer
487  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
488  //Set the sign of R
489  r->sign = 1;
490 
491  //Generate a random pattern
492  error = prngAlgo->read(prngContext, (uint8_t *) r->data, n * MPI_INT_SIZE);
493  //Any error to report?
494  if(error)
495  return error;
496 
497  //Remove the meaningless bits in the most significant word
498  if(n > 0 && m > 0)
499  {
500  r->data[n - 1] &= (1 << m) - 1;
501  }
502 
503  //Successful operation
504  return NO_ERROR;
505 }
506 
507 
508 /**
509  * @brief Generate a random value in the range 1 to p-1
510  * @param[out] r Pointer to a multiple precision integer
511  * @param[in] p The upper bound of the range
512  * @param[in] prngAlgo PRNG algorithm
513  * @param[in] prngContext Pointer to the PRNG context
514  * @return Error code
515  **/
516 
517 error_t mpiRandRange(Mpi *r, const Mpi *p, const PrngAlgo *prngAlgo,
518  void *prngContext)
519 {
520  error_t error;
521  uint_t n;
522  Mpi a;
523 
524  //Make sure p is greater than 1
525  if(mpiCompInt(p, 1) <= 0)
527 
528  //Initialize multiple precision integer
529  mpiInit(&a);
530 
531  //Get the actual length of p
532  n = mpiGetBitLength(p);
533 
534  //Generate extra random bits so that the bias produced by the modular
535  //reduction is negligible
536  MPI_CHECK(mpiRand(r, n + 64, prngAlgo, prngContext));
537 
538  //Compute r = (r mod (p - 1)) + 1
539  MPI_CHECK(mpiSubInt(&a, p, 1));
540  MPI_CHECK(mpiMod(r, r, &a));
541  MPI_CHECK(mpiAddInt(r, r, 1));
542 
543 end:
544  //Release previously allocated memory
545  mpiFree(&a);
546 
547  //Return status code
548  return error;
549 }
550 
551 
552 /**
553  * @brief Test whether a number is probable prime
554  * @param[in] a Pointer to a multiple precision integer
555  * @return Error code
556  **/
557 
558 __weak_func error_t mpiCheckProbablePrime(const Mpi *a)
559 {
560  //The algorithm is implemented by hardware
561  return ERROR_NOT_IMPLEMENTED;
562 }
563 
564 
565 /**
566  * @brief Octet string to integer conversion
567  *
568  * Converts an octet string to a non-negative integer
569  *
570  * @param[out] r Non-negative integer resulting from the conversion
571  * @param[in] data Octet string to be converted
572  * @param[in] length Length of the octet string
573  * @param[in] format Input format
574  * @return Error code
575  **/
576 
577 error_t mpiImport(Mpi *r, const uint8_t *data, uint_t length, MpiFormat format)
578 {
579  error_t error;
580  uint_t i;
581 
582  //Check input format
583  if(format == MPI_FORMAT_LITTLE_ENDIAN)
584  {
585  //Skip trailing zeroes
586  while(length > 0 && data[length - 1] == 0)
587  {
588  length--;
589  }
590 
591  //Ajust the size of the multiple precision integer
592  error = mpiGrow(r, (length + MPI_INT_SIZE - 1) / MPI_INT_SIZE);
593 
594  //Check status code
595  if(!error)
596  {
597  //Clear the contents of the multiple precision integer
598  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
599  //Set sign
600  r->sign = 1;
601 
602  //Import data
603  for(i = 0; i < length; i++, data++)
604  {
605  r->data[i / MPI_INT_SIZE] |= *data << ((i % MPI_INT_SIZE) * 8);
606  }
607  }
608  }
609  else if(format == MPI_FORMAT_BIG_ENDIAN)
610  {
611  //Skip leading zeroes
612  while(length > 1 && *data == 0)
613  {
614  data++;
615  length--;
616  }
617 
618  //Ajust the size of the multiple precision integer
619  error = mpiGrow(r, (length + MPI_INT_SIZE - 1) / MPI_INT_SIZE);
620 
621  //Check status code
622  if(!error)
623  {
624  //Clear the contents of the multiple precision integer
625  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
626  //Set sign
627  r->sign = 1;
628 
629  //Start from the least significant byte
630  data += length - 1;
631 
632  //Import data
633  for(i = 0; i < length; i++, data--)
634  {
635  r->data[i / MPI_INT_SIZE] |= *data << ((i % MPI_INT_SIZE) * 8);
636  }
637  }
638  }
639  else
640  {
641  //Report an error
642  error = ERROR_INVALID_PARAMETER;
643  }
644 
645  //Return status code
646  return error;
647 }
648 
649 
650 /**
651  * @brief Integer to octet string conversion
652  *
653  * Converts an integer to an octet string of a specified length
654  *
655  * @param[in] a Non-negative integer to be converted
656  * @param[out] data Octet string resulting from the conversion
657  * @param[in] length Intended length of the resulting octet string
658  * @param[in] format Output format
659  * @return Error code
660  **/
661 
662 error_t mpiExport(const Mpi *a, uint8_t *data, uint_t length, MpiFormat format)
663 {
664  uint_t i;
665  uint_t n;
666  error_t error;
667 
668  //Initialize status code
669  error = NO_ERROR;
670 
671  //Check input format
672  if(format == MPI_FORMAT_LITTLE_ENDIAN)
673  {
674  //Get the actual length in bytes
675  n = mpiGetByteLength(a);
676 
677  //Make sure the output buffer is large enough
678  if(n <= length)
679  {
680  //Clear output buffer
681  osMemset(data, 0, length);
682 
683  //Export data
684  for(i = 0; i < n; i++, data++)
685  {
686  *data = a->data[i / MPI_INT_SIZE] >> ((i % MPI_INT_SIZE) * 8);
687  }
688  }
689  else
690  {
691  //Report an error
692  error = ERROR_INVALID_LENGTH;
693  }
694  }
695  else if(format == MPI_FORMAT_BIG_ENDIAN)
696  {
697  //Get the actual length in bytes
698  n = mpiGetByteLength(a);
699 
700  //Make sure the output buffer is large enough
701  if(n <= length)
702  {
703  //Clear output buffer
704  osMemset(data, 0, length);
705 
706  //Point to the least significant word
707  data += length - 1;
708 
709  //Export data
710  for(i = 0; i < n; i++, data--)
711  {
712  *data = a->data[i / MPI_INT_SIZE] >> ((i % MPI_INT_SIZE) * 8);
713  }
714  }
715  else
716  {
717  //Report an error
718  error = ERROR_INVALID_LENGTH;
719  }
720  }
721  else
722  {
723  //Report an error
724  error = ERROR_INVALID_PARAMETER;
725  }
726 
727  //Return status code
728  return error;
729 }
730 
731 
732 /**
733  * @brief Multiple precision addition
734  * @param[out] r Resulting integer R = A + B
735  * @param[in] a First operand A
736  * @param[in] b Second operand B
737  * @return Error code
738  **/
739 
740 error_t mpiAdd(Mpi *r, const Mpi *a, const Mpi *b)
741 {
742  error_t error;
743  int_t sign;
744 
745  //Retrieve the sign of A
746  sign = a->sign;
747 
748  //Both operands have the same sign?
749  if(a->sign == b->sign)
750  {
751  //Perform addition
752  error = mpiAddAbs(r, a, b);
753  //Set the sign of the resulting number
754  r->sign = sign;
755  }
756  //Operands have different signs?
757  else
758  {
759  //Compare the absolute value of A and B
760  if(mpiCompAbs(a, b) >= 0)
761  {
762  //Perform subtraction
763  error = mpiSubAbs(r, a, b);
764  //Set the sign of the resulting number
765  r->sign = sign;
766  }
767  else
768  {
769  //Perform subtraction
770  error = mpiSubAbs(r, b, a);
771  //Set the sign of the resulting number
772  r->sign = -sign;
773  }
774  }
775 
776  //Return status code
777  return error;
778 }
779 
780 
781 /**
782  * @brief Add an integer to a multiple precision integer
783  * @param[out] r Resulting integer R = A + B
784  * @param[in] a First operand A
785  * @param[in] b Second operand B
786  * @return Error code
787  **/
788 
790 {
791  uint_t value;
792  Mpi t;
793 
794  //Convert the second operand to a multiple precision integer
795  value = (b >= 0) ? b : -b;
796  t.sign = (b >= 0) ? 1 : -1;
797  t.size = 1;
798  t.data = &value;
799 
800  //Perform addition
801  return mpiAdd(r, a, &t);
802 }
803 
804 
805 /**
806  * @brief Multiple precision subtraction
807  * @param[out] r Resulting integer R = A - B
808  * @param[in] a First operand A
809  * @param[in] b Second operand B
810  * @return Error code
811  **/
812 
813 error_t mpiSub(Mpi *r, const Mpi *a, const Mpi *b)
814 {
815  error_t error;
816  int_t sign;
817 
818  //Retrieve the sign of A
819  sign = a->sign;
820 
821  //Both operands have the same sign?
822  if(a->sign == b->sign)
823  {
824  //Compare the absolute value of A and B
825  if(mpiCompAbs(a, b) >= 0)
826  {
827  //Perform subtraction
828  error = mpiSubAbs(r, a, b);
829  //Set the sign of the resulting number
830  r->sign = sign;
831  }
832  else
833  {
834  //Perform subtraction
835  error = mpiSubAbs(r, b, a);
836  //Set the sign of the resulting number
837  r->sign = -sign;
838  }
839  }
840  //Operands have different signs?
841  else
842  {
843  //Perform addition
844  error = mpiAddAbs(r, a, b);
845  //Set the sign of the resulting number
846  r->sign = sign;
847  }
848 
849  //Return status code
850  return error;
851 }
852 
853 
854 /**
855  * @brief Subtract an integer from a multiple precision integer
856  * @param[out] r Resulting integer R = A - B
857  * @param[in] a First operand A
858  * @param[in] b Second operand B
859  * @return Error code
860  **/
861 
863 {
864  uint_t value;
865  Mpi t;
866 
867  //Convert the second operand to a multiple precision integer
868  value = (b >= 0) ? b : -b;
869  t.sign = (b >= 0) ? 1 : -1;
870  t.size = 1;
871  t.data = &value;
872 
873  //Perform subtraction
874  return mpiSub(r, a, &t);
875 }
876 
877 
878 /**
879  * @brief Helper routine for multiple precision addition
880  * @param[out] r Resulting integer R = |A + B|
881  * @param[in] a First operand A
882  * @param[in] b Second operand B
883  * @return Error code
884  **/
885 
886 error_t mpiAddAbs(Mpi *r, const Mpi *a, const Mpi *b)
887 {
888  error_t error;
889  uint_t i;
890  uint_t n;
891  uint_t c;
892  uint_t d;
893 
894  //R and B are the same instance?
895  if(r == b)
896  {
897  //Swap A and B
898  const Mpi *t = a;
899  a = b;
900  b = t;
901  }
902  //R is neither A nor B?
903  else if(r != a)
904  {
905  //Copy the first operand to R
906  MPI_CHECK(mpiCopy(r, a));
907  }
908 
909  //Determine the actual length of B
910  n = mpiGetLength(b);
911  //Extend the size of the destination register as needed
912  MPI_CHECK(mpiGrow(r, n));
913 
914  //The result is always positive
915  r->sign = 1;
916  //Clear carry bit
917  c = 0;
918 
919  //Add operands
920  for(i = 0; i < n; i++)
921  {
922  //Add carry bit
923  d = r->data[i] + c;
924  //Update carry bit
925  if(d != 0) c = 0;
926  //Perform addition
927  d += b->data[i];
928  //Update carry bit
929  if(d < b->data[i]) c = 1;
930  //Save result
931  r->data[i] = d;
932  }
933 
934  //Loop as long as the carry bit is set
935  for(i = n; c && i < r->size; i++)
936  {
937  //Add carry bit
938  r->data[i] += c;
939  //Update carry bit
940  if(r->data[i] != 0) c = 0;
941  }
942 
943  //Check the final carry bit
944  if(c && n >= r->size)
945  {
946  //Extend the size of the destination register
947  MPI_CHECK(mpiGrow(r, n + 1));
948  //Add carry bit
949  r->data[n] = 1;
950  }
951 
952 end:
953  //Return status code
954  return error;
955 }
956 
957 
958 /**
959  * @brief Helper routine for multiple precision subtraction
960  * @param[out] r Resulting integer R = |A - B|
961  * @param[in] a First operand A
962  * @param[in] b Second operand B
963  * @return Error code
964  **/
965 
966 error_t mpiSubAbs(Mpi *r, const Mpi *a, const Mpi *b)
967 {
968  error_t error;
969  uint_t c;
970  uint_t d;
971  uint_t i;
972  uint_t m;
973  uint_t n;
974 
975  //Check input parameters
976  if(mpiCompAbs(a, b) < 0)
977  {
978  //Swap A and B if necessary
979  const Mpi *t = b;
980  a = b;
981  b = t;
982  }
983 
984  //Determine the actual length of A
985  m = mpiGetLength(a);
986  //Determine the actual length of B
987  n = mpiGetLength(b);
988 
989  //Extend the size of the destination register as needed
990  MPI_CHECK(mpiGrow(r, m));
991 
992  //The result is always positive
993  r->sign = 1;
994  //Clear carry bit
995  c = 0;
996 
997  //Subtract operands
998  for(i = 0; i < n; i++)
999  {
1000  //Read first operand
1001  d = a->data[i];
1002 
1003  //Check the carry bit
1004  if(c)
1005  {
1006  //Update carry bit
1007  if(d != 0) c = 0;
1008  //Propagate carry bit
1009  d -= 1;
1010  }
1011 
1012  //Update carry bit
1013  if(d < b->data[i]) c = 1;
1014  //Perform subtraction
1015  r->data[i] = d - b->data[i];
1016  }
1017 
1018  //Loop as long as the carry bit is set
1019  for(i = n; c && i < m; i++)
1020  {
1021  //Update carry bit
1022  if(a->data[i] != 0) c = 0;
1023  //Propagate carry bit
1024  r->data[i] = a->data[i] - 1;
1025  }
1026 
1027  //R and A are not the same instance?
1028  if(r != a)
1029  {
1030  //Copy the remaining words
1031  for(; i < m; i++)
1032  {
1033  r->data[i] = a->data[i];
1034  }
1035 
1036  //Zero the upper part of R
1037  for(; i < r->size; i++)
1038  {
1039  r->data[i] = 0;
1040  }
1041  }
1042 
1043 end:
1044  //Return status code
1045  return error;
1046 }
1047 
1048 
1049 /**
1050  * @brief Left shift operation
1051  * @param[in,out] r The multiple precision integer to be shifted to the left
1052  * @param[in] n The number of bits to shift
1053  * @return Error code
1054  **/
1055 
1057 {
1058  error_t error;
1059  uint_t i;
1060 
1061  //Number of 32-bit words to shift
1062  uint_t n1 = n / (MPI_INT_SIZE * 8);
1063  //Number of bits to shift
1064  uint_t n2 = n % (MPI_INT_SIZE * 8);
1065 
1066  //Check parameters
1067  if(!r->size || !n)
1068  return NO_ERROR;
1069 
1070  //Increase the size of the multiple-precision number
1071  error = mpiGrow(r, r->size + (n + 31) / 32);
1072  //Check return code
1073  if(error)
1074  return error;
1075 
1076  //First, shift words
1077  if(n1 > 0)
1078  {
1079  //Process the most significant words
1080  for(i = r->size - 1; i >= n1; i--)
1081  {
1082  r->data[i] = r->data[i - n1];
1083  }
1084 
1085  //Fill the rest with zeroes
1086  for(i = 0; i < n1; i++)
1087  {
1088  r->data[i] = 0;
1089  }
1090  }
1091 
1092  //Then shift bits
1093  if(n2 > 0)
1094  {
1095  //Process the most significant words
1096  for(i = r->size - 1; i >= 1; i--)
1097  {
1098  r->data[i] = (r->data[i] << n2) | (r->data[i - 1] >> (32 - n2));
1099  }
1100 
1101  //The least significant word requires a special handling
1102  r->data[0] <<= n2;
1103  }
1104 
1105  //Shift operation is complete
1106  return NO_ERROR;
1107 }
1108 
1109 
1110 /**
1111  * @brief Right shift operation
1112  * @param[in,out] r The multiple precision integer to be shifted to the right
1113  * @param[in] n The number of bits to shift
1114  * @return Error code
1115  **/
1116 
1118 {
1119  uint_t i;
1120  uint_t m;
1121 
1122  //Number of 32-bit words to shift
1123  uint_t n1 = n / (MPI_INT_SIZE * 8);
1124  //Number of bits to shift
1125  uint_t n2 = n % (MPI_INT_SIZE * 8);
1126 
1127  //Check parameters
1128  if(n1 >= r->size)
1129  {
1130  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
1131  return NO_ERROR;
1132  }
1133 
1134  //First, shift words
1135  if(n1 > 0)
1136  {
1137  //Process the least significant words
1138  for(m = r->size - n1, i = 0; i < m; i++)
1139  {
1140  r->data[i] = r->data[i + n1];
1141  }
1142 
1143  //Fill the rest with zeroes
1144  for(i = m; i < r->size; i++)
1145  {
1146  r->data[i] = 0;
1147  }
1148  }
1149 
1150  //Then shift bits
1151  if(n2 > 0)
1152  {
1153  //Process the least significant words
1154  for(m = r->size - n1 - 1, i = 0; i < m; i++)
1155  {
1156  r->data[i] = (r->data[i] >> n2) | (r->data[i + 1] << (32 - n2));
1157  }
1158 
1159  //The most significant word requires a special handling
1160  r->data[m] >>= n2;
1161  }
1162 
1163  //Shift operation is complete
1164  return NO_ERROR;
1165 }
1166 
1167 
1168 /**
1169  * @brief Multiple precision multiplication
1170  * @param[out] r Resulting integer R = A * B
1171  * @param[in] a First operand A
1172  * @param[in] b Second operand B
1173  * @return Error code
1174  **/
1175 
1176 __weak_func error_t mpiMul(Mpi *r, const Mpi *a, const Mpi *b)
1177 {
1178  error_t error;
1179  int_t i;
1180  int_t m;
1181  int_t n;
1182  Mpi ta;
1183  Mpi tb;
1184 
1185  //Initialize multiple precision integers
1186  mpiInit(&ta);
1187  mpiInit(&tb);
1188 
1189  //R and A are the same instance?
1190  if(r == a)
1191  {
1192  //Copy A to TA
1193  MPI_CHECK(mpiCopy(&ta, a));
1194  //Use TA instead of A
1195  a = &ta;
1196  }
1197 
1198  //R and B are the same instance?
1199  if(r == b)
1200  {
1201  //Copy B to TB
1202  MPI_CHECK(mpiCopy(&tb, b));
1203  //Use TB instead of B
1204  b = &tb;
1205  }
1206 
1207  //Determine the actual length of A and B
1208  m = mpiGetLength(a);
1209  n = mpiGetLength(b);
1210 
1211  //Adjust the size of R
1212  MPI_CHECK(mpiGrow(r, m + n));
1213  //Set the sign of R
1214  r->sign = (a->sign == b->sign) ? 1 : -1;
1215 
1216  //Clear the contents of the destination integer
1217  osMemset(r->data, 0, r->size * MPI_INT_SIZE);
1218 
1219  //Perform multiplication
1220  if(m < n)
1221  {
1222  for(i = 0; i < m; i++)
1223  {
1224  mpiMulAccCore(&r->data[i], b->data, n, a->data[i]);
1225  }
1226  }
1227  else
1228  {
1229  for(i = 0; i < n; i++)
1230  {
1231  mpiMulAccCore(&r->data[i], a->data, m, b->data[i]);
1232  }
1233  }
1234 
1235 end:
1236  //Release multiple precision integers
1237  mpiFree(&ta);
1238  mpiFree(&tb);
1239 
1240  //Return status code
1241  return error;
1242 }
1243 
1244 
1245 /**
1246  * @brief Multiply a multiple precision integer by an integer
1247  * @param[out] r Resulting integer R = A * B
1248  * @param[in] a First operand A
1249  * @param[in] b Second operand B
1250  * @return Error code
1251  **/
1252 
1254 {
1255  uint_t value;
1256  Mpi t;
1257 
1258  //Convert the second operand to a multiple precision integer
1259  value = (b >= 0) ? b : -b;
1260  t.sign = (b >= 0) ? 1 : -1;
1261  t.size = 1;
1262  t.data = &value;
1263 
1264  //Perform multiplication
1265  return mpiMul(r, a, &t);
1266 }
1267 
1268 
1269 /**
1270  * @brief Multiple precision division
1271  * @param[out] q The quotient Q = A / B
1272  * @param[out] r The remainder R = A mod B
1273  * @param[in] a The dividend A
1274  * @param[in] b The divisor B
1275  * @return Error code
1276  **/
1277 
1278 error_t mpiDiv(Mpi *q, Mpi *r, const Mpi *a, const Mpi *b)
1279 {
1280  error_t error;
1281  uint_t m;
1282  uint_t n;
1283  Mpi c;
1284  Mpi d;
1285  Mpi e;
1286 
1287  //Check whether the divisor is equal to zero
1288  if(!mpiCompInt(b, 0))
1289  return ERROR_INVALID_PARAMETER;
1290 
1291  //Initialize multiple precision integers
1292  mpiInit(&c);
1293  mpiInit(&d);
1294  mpiInit(&e);
1295 
1296  MPI_CHECK(mpiCopy(&c, a));
1297  MPI_CHECK(mpiCopy(&d, b));
1298  MPI_CHECK(mpiSetValue(&e, 0));
1299 
1300  m = mpiGetBitLength(&c);
1301  n = mpiGetBitLength(&d);
1302 
1303  if(m > n)
1304  MPI_CHECK(mpiShiftLeft(&d, m - n));
1305 
1306  while(n++ <= m)
1307  {
1308  MPI_CHECK(mpiShiftLeft(&e, 1));
1309 
1310  if(mpiComp(&c, &d) >= 0)
1311  {
1312  MPI_CHECK(mpiSetBitValue(&e, 0, 1));
1313  MPI_CHECK(mpiSub(&c, &c, &d));
1314  }
1315 
1316  MPI_CHECK(mpiShiftRight(&d, 1));
1317  }
1318 
1319  if(q != NULL)
1320  MPI_CHECK(mpiCopy(q, &e));
1321 
1322  if(r != NULL)
1323  MPI_CHECK(mpiCopy(r, &c));
1324 
1325 end:
1326  //Release previously allocated memory
1327  mpiFree(&c);
1328  mpiFree(&d);
1329  mpiFree(&e);
1330 
1331  //Return status code
1332  return error;
1333 }
1334 
1335 
1336 /**
1337  * @brief Divide a multiple precision integer by an integer
1338  * @param[out] q The quotient Q = A / B
1339  * @param[out] r The remainder R = A mod B
1340  * @param[in] a The dividend A
1341  * @param[in] b The divisor B
1342  * @return Error code
1343  **/
1344 
1346 {
1347  uint_t value;
1348  Mpi t;
1349 
1350  //Convert the divisor to a multiple precision integer
1351  value = (b >= 0) ? b : -b;
1352  t.sign = (b >= 0) ? 1 : -1;
1353  t.size = 1;
1354  t.data = &value;
1355 
1356  //Perform division
1357  return mpiDiv(q, r, a, &t);
1358 }
1359 
1360 
1361 /**
1362  * @brief Modulo operation
1363  * @param[out] r Resulting integer R = A mod P
1364  * @param[in] a The multiple precision integer to be reduced
1365  * @param[in] p The modulus P
1366  * @return Error code
1367  **/
1368 
1369 error_t mpiMod(Mpi *r, const Mpi *a, const Mpi *p)
1370 {
1371  error_t error;
1372  int_t sign;
1373  uint_t m;
1374  uint_t n;
1375  Mpi c;
1376 
1377  //Make sure the modulus is positive
1378  if(mpiCompInt(p, 0) <= 0)
1379  return ERROR_INVALID_PARAMETER;
1380 
1381  //Initialize multiple precision integer
1382  mpiInit(&c);
1383 
1384  //Save the sign of A
1385  sign = a->sign;
1386  //Determine the actual length of A
1387  m = mpiGetBitLength(a);
1388  //Determine the actual length of P
1389  n = mpiGetBitLength(p);
1390 
1391  //Let R = A
1392  MPI_CHECK(mpiCopy(r, a));
1393 
1394  if(m >= n)
1395  {
1396  MPI_CHECK(mpiCopy(&c, p));
1397  MPI_CHECK(mpiShiftLeft(&c, m - n));
1398 
1399  while(mpiCompAbs(r, p) >= 0)
1400  {
1401  if(mpiCompAbs(r, &c) >= 0)
1402  {
1403  MPI_CHECK(mpiSubAbs(r, r, &c));
1404  }
1405 
1406  MPI_CHECK(mpiShiftRight(&c, 1));
1407  }
1408  }
1409 
1410  if(sign < 0)
1411  {
1412  MPI_CHECK(mpiSubAbs(r, p, r));
1413  }
1414 
1415 end:
1416  //Release previously allocated memory
1417  mpiFree(&c);
1418 
1419  //Return status code
1420  return error;
1421 }
1422 
1423 
1424 
1425 /**
1426  * @brief Modular addition
1427  * @param[out] r Resulting integer R = A + B mod P
1428  * @param[in] a The first operand A
1429  * @param[in] b The second operand B
1430  * @param[in] p The modulus P
1431  * @return Error code
1432  **/
1433 
1434 error_t mpiAddMod(Mpi *r, const Mpi *a, const Mpi *b, const Mpi *p)
1435 {
1436  error_t error;
1437 
1438  //Perform modular addition
1439  MPI_CHECK(mpiAdd(r, a, b));
1440  MPI_CHECK(mpiMod(r, r, p));
1441 
1442 end:
1443  //Return status code
1444  return error;
1445 }
1446 
1447 
1448 /**
1449  * @brief Modular subtraction
1450  * @param[out] r Resulting integer R = A - B mod P
1451  * @param[in] a The first operand A
1452  * @param[in] b The second operand B
1453  * @param[in] p The modulus P
1454  * @return Error code
1455  **/
1456 
1457 error_t mpiSubMod(Mpi *r, const Mpi *a, const Mpi *b, const Mpi *p)
1458 {
1459  error_t error;
1460 
1461  //Perform modular subtraction
1462  MPI_CHECK(mpiSub(r, a, b));
1463  MPI_CHECK(mpiMod(r, r, p));
1464 
1465 end:
1466  //Return status code
1467  return error;
1468 }
1469 
1470 
1471 /**
1472  * @brief Modular multiplication
1473  * @param[out] r Resulting integer R = A * B mod P
1474  * @param[in] a The first operand A
1475  * @param[in] b The second operand B
1476  * @param[in] p The modulus P
1477  * @return Error code
1478  **/
1479 
1480 __weak_func error_t mpiMulMod(Mpi *r, const Mpi *a, const Mpi *b, const Mpi *p)
1481 {
1482  error_t error;
1483 
1484  //Perform modular multiplication
1485  MPI_CHECK(mpiMul(r, a, b));
1486  MPI_CHECK(mpiMod(r, r, p));
1487 
1488 end:
1489  //Return status code
1490  return error;
1491 }
1492 
1493 
1494 /**
1495  * @brief Modular inverse
1496  * @param[out] r Resulting integer R = A^-1 mod P
1497  * @param[in] a The multiple precision integer A
1498  * @param[in] p The modulus P
1499  * @return Error code
1500  **/
1501 
1502 __weak_func error_t mpiInvMod(Mpi *r, const Mpi *a, const Mpi *p)
1503 {
1504  error_t error;
1505  Mpi b;
1506  Mpi c;
1507  Mpi q0;
1508  Mpi r0;
1509  Mpi t;
1510  Mpi u;
1511  Mpi v;
1512 
1513  //Initialize multiple precision integers
1514  mpiInit(&b);
1515  mpiInit(&c);
1516  mpiInit(&q0);
1517  mpiInit(&r0);
1518  mpiInit(&t);
1519  mpiInit(&u);
1520  mpiInit(&v);
1521 
1522  MPI_CHECK(mpiCopy(&b, p));
1523  MPI_CHECK(mpiCopy(&c, a));
1524  MPI_CHECK(mpiSetValue(&u, 0));
1525  MPI_CHECK(mpiSetValue(&v, 1));
1526 
1527  while(mpiCompInt(&c, 0) > 0)
1528  {
1529  MPI_CHECK(mpiDiv(&q0, &r0, &b, &c));
1530 
1531  MPI_CHECK(mpiCopy(&b, &c));
1532  MPI_CHECK(mpiCopy(&c, &r0));
1533 
1534  MPI_CHECK(mpiCopy(&t, &v));
1535  MPI_CHECK(mpiMul(&q0, &q0, &v));
1536  MPI_CHECK(mpiSub(&v, &u, &q0));
1537  MPI_CHECK(mpiCopy(&u, &t));
1538  }
1539 
1540  if(mpiCompInt(&b, 1))
1541  {
1543  }
1544 
1545  if(mpiCompInt(&u, 0) > 0)
1546  {
1547  MPI_CHECK(mpiCopy(r, &u));
1548  }
1549  else
1550  {
1551  MPI_CHECK(mpiAdd(r, &u, p));
1552  }
1553 
1554 end:
1555  //Release previously allocated memory
1556  mpiFree(&b);
1557  mpiFree(&c);
1558  mpiFree(&q0);
1559  mpiFree(&r0);
1560  mpiFree(&t);
1561  mpiFree(&u);
1562  mpiFree(&v);
1563 
1564  //Return status code
1565  return error;
1566 }
1567 
1568 
1569 /**
1570  * @brief Modular exponentiation
1571  * @param[out] r Resulting integer R = A ^ E mod P
1572  * @param[in] a Pointer to a multiple precision integer
1573  * @param[in] e Exponent
1574  * @param[in] p Modulus
1575  * @return Error code
1576  **/
1577 
1578 __weak_func error_t mpiExpMod(Mpi *r, const Mpi *a, const Mpi *e, const Mpi *p)
1579 {
1580  error_t error;
1581  int_t i;
1582  int_t j;
1583  int_t n;
1584  uint_t d;
1585  uint_t k;
1586  uint_t u;
1587  Mpi b;
1588  Mpi c2;
1589  Mpi t;
1590  Mpi s[8];
1591 
1592  //Initialize multiple precision integers
1593  mpiInit(&b);
1594  mpiInit(&c2);
1595  mpiInit(&t);
1596 
1597  //Initialize precomputed values
1598  for(i = 0; (uint_t) i < arraysize(s); i++)
1599  {
1600  mpiInit(&s[i]);
1601  }
1602 
1603  //Very small exponents are often selected with low Hamming weight.
1604  //The sliding window mechanism should be disabled in that case
1605  d = (mpiGetBitLength(e) <= 32) ? 1 : 4;
1606 
1607  //Even modulus?
1608  if(mpiIsEven(p))
1609  {
1610  //Let B = A^2
1611  MPI_CHECK(mpiMulMod(&b, a, a, p));
1612  //Let S[0] = A
1613  MPI_CHECK(mpiCopy(&s[0], a));
1614 
1615  //Precompute S[i] = A^(2 * i + 1)
1616  for(i = 1; i < (1 << (d - 1)); i++)
1617  {
1618  MPI_CHECK(mpiMulMod(&s[i], &s[i - 1], &b, p));
1619  }
1620 
1621  //Let R = 1
1622  MPI_CHECK(mpiSetValue(r, 1));
1623 
1624  //The exponent is processed in a left-to-right fashion
1625  i = mpiGetBitLength(e) - 1;
1626 
1627  //Perform sliding window exponentiation
1628  while(i >= 0)
1629  {
1630  //The sliding window exponentiation algorithm decomposes E
1631  //into zero and nonzero windows
1632  if(!mpiGetBitValue(e, i))
1633  {
1634  //Compute R = R^2
1635  MPI_CHECK(mpiMulMod(r, r, r, p));
1636  //Next bit to be processed
1637  i--;
1638  }
1639  else
1640  {
1641  //Find the longest window
1642  n = MAX(i - d + 1, 0);
1643 
1644  //The least significant bit of the window must be equal to 1
1645  while(!mpiGetBitValue(e, n)) n++;
1646 
1647  //The algorithm processes more than one bit per iteration
1648  for(u = 0, j = i; j >= n; j--)
1649  {
1650  //Compute R = R^2
1651  MPI_CHECK(mpiMulMod(r, r, r, p));
1652  //Compute the relevant index to be used in the precomputed table
1653  u = (u << 1) | mpiGetBitValue(e, j);
1654  }
1655 
1656  //Perform a single multiplication per iteration
1657  MPI_CHECK(mpiMulMod(r, r, &s[u >> 1], p));
1658  //Next bit to be processed
1659  i = n - 1;
1660  }
1661  }
1662  }
1663  else
1664  {
1665  //Compute the smaller C = (2^32)^k such as C > P
1666  k = mpiGetLength(p);
1667 
1668  //Compute C^2 mod P
1669  MPI_CHECK(mpiSetValue(&c2, 1));
1670  MPI_CHECK(mpiShiftLeft(&c2, 2 * k * (MPI_INT_SIZE * 8)));
1671  MPI_CHECK(mpiMod(&c2, &c2, p));
1672 
1673  //Let B = A * C mod P
1674  if(mpiComp(a, p) >= 0)
1675  {
1676  MPI_CHECK(mpiMod(&b, a, p));
1677  MPI_CHECK(mpiMontgomeryMul(&b, &b, &c2, k, p, &t));
1678  }
1679  else
1680  {
1681  MPI_CHECK(mpiMontgomeryMul(&b, a, &c2, k, p, &t));
1682  }
1683 
1684  //Let R = B^2 * C^-1 mod P
1685  MPI_CHECK(mpiMontgomeryMul(r, &b, &b, k, p, &t));
1686  //Let S[0] = B
1687  MPI_CHECK(mpiCopy(&s[0], &b));
1688 
1689  //Precompute S[i] = B^(2 * i + 1) * C^-1 mod P
1690  for(i = 1; i < (1 << (d - 1)); i++)
1691  {
1692  MPI_CHECK(mpiMontgomeryMul(&s[i], &s[i - 1], r, k, p, &t));
1693  }
1694 
1695  //Let R = C mod P
1696  MPI_CHECK(mpiCopy(r, &c2));
1697  MPI_CHECK(mpiMontgomeryRed(r, r, k, p, &t));
1698 
1699  //The exponent is processed in a left-to-right fashion
1700  i = mpiGetBitLength(e) - 1;
1701 
1702  //Perform sliding window exponentiation
1703  while(i >= 0)
1704  {
1705  //The sliding window exponentiation algorithm decomposes E
1706  //into zero and nonzero windows
1707  if(!mpiGetBitValue(e, i))
1708  {
1709  //Compute R = R^2 * C^-1 mod P
1710  MPI_CHECK(mpiMontgomeryMul(r, r, r, k, p, &t));
1711  //Next bit to be processed
1712  i--;
1713  }
1714  else
1715  {
1716  //Find the longest window
1717  n = MAX(i - d + 1, 0);
1718 
1719  //The least significant bit of the window must be equal to 1
1720  while(!mpiGetBitValue(e, n)) n++;
1721 
1722  //The algorithm processes more than one bit per iteration
1723  for(u = 0, j = i; j >= n; j--)
1724  {
1725  //Compute R = R^2 * C^-1 mod P
1726  MPI_CHECK(mpiMontgomeryMul(r, r, r, k, p, &t));
1727  //Compute the relevant index to be used in the precomputed table
1728  u = (u << 1) | mpiGetBitValue(e, j);
1729  }
1730 
1731  //Compute R = R * T[u/2] * C^-1 mod P
1732  MPI_CHECK(mpiMontgomeryMul(r, r, &s[u >> 1], k, p, &t));
1733  //Next bit to be processed
1734  i = n - 1;
1735  }
1736  }
1737 
1738  //Compute R = R * C^-1 mod P
1739  MPI_CHECK(mpiMontgomeryRed(r, r, k, p, &t));
1740  }
1741 
1742 end:
1743  //Release multiple precision integers
1744  mpiFree(&b);
1745  mpiFree(&c2);
1746  mpiFree(&t);
1747 
1748  //Release precomputed values
1749  for(i = 0; (uint_t) i < arraysize(s); i++)
1750  {
1751  mpiFree(&s[i]);
1752  }
1753 
1754  //Return status code
1755  return error;
1756 }
1757 
1758 
1759 /**
1760  * @brief Modular exponentiation (fast calculation)
1761  * @param[out] r Resulting integer R = A ^ E mod P
1762  * @param[in] a Pointer to a multiple precision integer
1763  * @param[in] e Exponent
1764  * @param[in] p Modulus
1765  * @return Error code
1766  **/
1767 
1768 __weak_func error_t mpiExpModFast(Mpi *r, const Mpi *a, const Mpi *e, const Mpi *p)
1769 {
1770  //Perform modular exponentiation
1771  return mpiExpMod(r, a, e, p);
1772 }
1773 
1774 
1775 /**
1776  * @brief Modular exponentiation (regular calculation)
1777  * @param[out] r Resulting integer R = A ^ E mod P
1778  * @param[in] a Pointer to a multiple precision integer
1779  * @param[in] e Exponent
1780  * @param[in] p Modulus
1781  * @return Error code
1782  **/
1783 
1784 __weak_func error_t mpiExpModRegular(Mpi *r, const Mpi *a, const Mpi *e, const Mpi *p)
1785 {
1786  //Perform modular exponentiation
1787  return mpiExpMod(r, a, e, p);
1788 }
1789 
1790 
1791 /**
1792  * @brief Montgomery multiplication
1793  * @param[out] r Resulting integer R = A * B / 2^k mod P
1794  * @param[in] a An integer A such as 0 <= A < 2^k
1795  * @param[in] b An integer B such as 0 <= B < 2^k
1796  * @param[in] k An integer k such as P < 2^k
1797  * @param[in] p Modulus P
1798  * @param[in] t An preallocated integer T (for internal operation)
1799  * @return Error code
1800  **/
1801 
1802 error_t mpiMontgomeryMul(Mpi *r, const Mpi *a, const Mpi *b, uint_t k,
1803  const Mpi *p, Mpi *t)
1804 {
1805  error_t error;
1806  uint_t i;
1807  uint_t m;
1808  uint_t n;
1809  uint_t q;
1810 
1811  //Use Newton's method to compute the inverse of P[0] mod 2^32
1812  for(m = 2 - p->data[0], i = 0; i < 4; i++)
1813  {
1814  m = m * (2 - m * p->data[0]);
1815  }
1816 
1817  //Precompute -1/P[0] mod 2^32;
1818  m = ~m + 1;
1819 
1820  //We assume that B is always less than 2^k
1821  n = MIN(b->size, k);
1822 
1823  //Make sure T is large enough
1824  MPI_CHECK(mpiGrow(t, 2 * k + 1));
1825  //Let T = 0
1826  MPI_CHECK(mpiSetValue(t, 0));
1827 
1828  //Perform Montgomery multiplication
1829  for(i = 0; i < k; i++)
1830  {
1831  //Check current index
1832  if(i < a->size)
1833  {
1834  //Compute q = ((T[i] + A[i] * B[0]) * m) mod 2^32
1835  q = (t->data[i] + a->data[i] * b->data[0]) * m;
1836  //Compute T = T + A[i] * B
1837  mpiMulAccCore(t->data + i, b->data, n, a->data[i]);
1838  }
1839  else
1840  {
1841  //Compute q = (T[i] * m) mod 2^32
1842  q = t->data[i] * m;
1843  }
1844 
1845  //Compute T = T + q * P
1846  mpiMulAccCore(t->data + i, p->data, k, q);
1847  }
1848 
1849  //Compute R = T / 2^(32 * k)
1850  MPI_CHECK(mpiShiftRight(t, k * (MPI_INT_SIZE * 8)));
1851  MPI_CHECK(mpiCopy(r, t));
1852 
1853  //A final subtraction is required
1854  if(mpiComp(r, p) >= 0)
1855  {
1856  MPI_CHECK(mpiSub(r, r, p));
1857  }
1858 
1859 end:
1860  //Return status code
1861  return error;
1862 }
1863 
1864 
1865 /**
1866  * @brief Montgomery reduction
1867  * @param[out] r Resulting integer R = A / 2^k mod P
1868  * @param[in] a An integer A such as 0 <= A < 2^k
1869  * @param[in] k An integer k such as P < 2^k
1870  * @param[in] p Modulus P
1871  * @param[in] t An preallocated integer T (for internal operation)
1872  * @return Error code
1873  **/
1874 
1875 error_t mpiMontgomeryRed(Mpi *r, const Mpi *a, uint_t k, const Mpi *p, Mpi *t)
1876 {
1877  uint_t value;
1878  Mpi b;
1879 
1880  //Let B = 1
1881  value = 1;
1882  b.sign = 1;
1883  b.size = 1;
1884  b.data = &value;
1885 
1886  //Compute R = A / 2^k mod P
1887  return mpiMontgomeryMul(r, a, &b, k, p, t);
1888 }
1889 
1890 
1891 #if (MPI_ASM_SUPPORT == DISABLED)
1892 
1893 /**
1894  * @brief Multiply-accumulate operation
1895  * @param[out] r Resulting integer
1896  * @param[in] a First operand A
1897  * @param[in] m Size of A in words
1898  * @param[in] b Second operand B
1899  **/
1900 
1901 void mpiMulAccCore(uint_t *r, const uint_t *a, int_t m, const uint_t b)
1902 {
1903  int_t i;
1904  uint32_t c;
1905  uint32_t u;
1906  uint32_t v;
1907  uint64_t p;
1908 
1909  //Clear variables
1910  c = 0;
1911  u = 0;
1912  v = 0;
1913 
1914  //Perform multiplication
1915  for(i = 0; i < m; i++)
1916  {
1917  p = (uint64_t) a[i] * b;
1918  u = (uint32_t) p;
1919  v = (uint32_t) (p >> 32);
1920 
1921  u += c;
1922  if(u < c) v++;
1923 
1924  u += r[i];
1925  if(u < r[i]) v++;
1926 
1927  r[i] = u;
1928  c = v;
1929  }
1930 
1931  //Propagate carry
1932  for(; c != 0; i++)
1933  {
1934  r[i] += c;
1935  c = (r[i] < c);
1936  }
1937 }
1938 
1939 #endif
1940 
1941 
1942 /**
1943  * @brief Display the contents of a multiple precision integer
1944  * @param[in] stream Pointer to a FILE object that identifies an output stream
1945  * @param[in] prepend String to prepend to the left of each line
1946  * @param[in] a Pointer to a multiple precision integer
1947  **/
1948 
1949 void mpiDump(FILE *stream, const char_t *prepend, const Mpi *a)
1950 {
1951  uint_t i;
1952 
1953  //Process each word
1954  for(i = 0; i < a->size; i++)
1955  {
1956  //Beginning of a new line?
1957  if(i == 0 || ((a->size - i - 1) % 8) == 7)
1958  fprintf(stream, "%s", prepend);
1959 
1960  //Display current data
1961  fprintf(stream, "%08X ", a->data[a->size - 1 - i]);
1962 
1963  //End of current line?
1964  if(((a->size - i - 1) % 8) == 0 || i == (a->size - 1))
1965  fprintf(stream, "\r\n");
1966  }
1967 }
1968 
1969 #endif
uint8_t length
Definition: coap_common.h:193
error_t mpiShiftLeft(Mpi *r, uint_t n)
Left shift operation.
Definition: mpi.c:1056
error_t mpiSubInt(Mpi *r, const Mpi *a, int_t b)
Subtract an integer from a multiple precision integer.
Definition: mpi.c:862
__weak_func error_t mpiExpMod(Mpi *r, const Mpi *a, const Mpi *e, const Mpi *p)
Modular exponentiation.
Definition: mpi.c:1578
uint8_t a
Definition: ndp.h:409
uint8_t data[]
Definition: ethernet.h:220
Arbitrary precision integer.
Definition: mpi.h:70
signed int int_t
Definition: compiler_port.h:49
error_t mpiMontgomeryRed(Mpi *r, const Mpi *a, uint_t k, const Mpi *p, Mpi *t)
Montgomery reduction.
Definition: mpi.c:1875
@ ERROR_NOT_IMPLEMENTED
Definition: error.h:66
#define PrngAlgo
Definition: crypto.h:861
error_t mpiSubAbs(Mpi *r, const Mpi *a, const Mpi *b)
Helper routine for multiple precision subtraction.
Definition: mpi.c:966
error_t mpiDivInt(Mpi *q, Mpi *r, const Mpi *a, int_t b)
Divide a multiple precision integer by an integer.
Definition: mpi.c:1345
uint8_t p
Definition: ndp.h:298
uint8_t t
Definition: lldp_ext_med.h:210
__weak_func error_t mpiMulMod(Mpi *r, const Mpi *a, const Mpi *b, const Mpi *p)
Modular multiplication.
Definition: mpi.c:1480
#define mpiIsEven(a)
Definition: mpi.h:45
error_t mpiSetValue(Mpi *r, int_t a)
Set the value of a multiple precision integer.
Definition: mpi.c:437
@ ERROR_OUT_OF_MEMORY
Definition: error.h:63
error_t mpiSetBitValue(Mpi *r, uint_t index, uint_t value)
Set the bit value at the specified index.
Definition: mpi.c:236
error_t mpiAddAbs(Mpi *r, const Mpi *a, const Mpi *b)
Helper routine for multiple precision addition.
Definition: mpi.c:886
void mpiDump(FILE *stream, const char_t *prepend, const Mpi *a)
Display the contents of a multiple precision integer.
Definition: mpi.c:1949
error_t mpiRandRange(Mpi *r, const Mpi *p, const PrngAlgo *prngAlgo, void *prngContext)
Generate a random value in the range 1 to p-1.
Definition: mpi.c:517
error_t mpiShiftRight(Mpi *r, uint_t n)
Right shift operation.
Definition: mpi.c:1117
error_t mpiRand(Mpi *r, uint_t length, const PrngAlgo *prngAlgo, void *prngContext)
Generate a random value.
Definition: mpi.c:468
error_t mpiImport(Mpi *r, const uint8_t *data, uint_t length, MpiFormat format)
Octet string to integer conversion.
Definition: mpi.c:577
int_t mpiCompAbs(const Mpi *a, const Mpi *b)
Compare the absolute value of two multiple precision integers.
Definition: mpi.c:362
void mpiInit(Mpi *r)
Initialize a multiple precision integer.
Definition: mpi.c:48
uint8_t r
Definition: ndp.h:344
error_t mpiMod(Mpi *r, const Mpi *a, const Mpi *p)
Modulo operation.
Definition: mpi.c:1369
__weak_func error_t mpiExpModRegular(Mpi *r, const Mpi *a, const Mpi *e, const Mpi *p)
Modular exponentiation (regular calculation)
Definition: mpi.c:1784
error_t mpiExport(const Mpi *a, uint8_t *data, uint_t length, MpiFormat format)
Integer to octet string conversion.
Definition: mpi.c:662
@ ERROR_INVALID_PARAMETER
Invalid parameter.
Definition: error.h:47
#define osMemcpy(dest, src, length)
Definition: os_port.h:140
error_t mpiSubMod(Mpi *r, const Mpi *a, const Mpi *b, const Mpi *p)
Modular subtraction.
Definition: mpi.c:1457
error_t mpiSub(Mpi *r, const Mpi *a, const Mpi *b)
Multiple precision subtraction.
Definition: mpi.c:813
__weak_func error_t mpiInvMod(Mpi *r, const Mpi *a, const Mpi *p)
Modular inverse.
Definition: mpi.c:1502
@ MPI_FORMAT_LITTLE_ENDIAN
Definition: mpi.h:60
error_t
Error codes.
Definition: error.h:43
#define MPI_CHECK(f)
Definition: mpi.h:42
error_t mpiAdd(Mpi *r, const Mpi *a, const Mpi *b)
Multiple precision addition.
Definition: mpi.c:740
uint8_t value[]
Definition: tcp.h:367
@ ERROR_FAILURE
Generic error code.
Definition: error.h:45
MPI (Multiple Precision Integer Arithmetic)
@ ERROR_INVALID_LENGTH
Definition: error.h:111
General definitions for cryptographic algorithms.
uint8_t u
Definition: lldp_ext_med.h:211
#define MIN(a, b)
Definition: os_port.h:65
error_t mpiDiv(Mpi *q, Mpi *r, const Mpi *a, const Mpi *b)
Multiple precision division.
Definition: mpi.c:1278
uint_t mpiGetBitLength(const Mpi *a)
Get the actual length in bits.
Definition: mpi.c:195
error_t mpiCopy(Mpi *r, const Mpi *a)
Copy a multiple precision integer.
Definition: mpi.c:400
MpiFormat
MPI import/export format.
Definition: mpi.h:59
error_t mpiAddInt(Mpi *r, const Mpi *a, int_t b)
Add an integer to a multiple precision integer.
Definition: mpi.c:789
uint8_t b[6]
Definition: ethernet.h:190
#define MAX(a, b)
Definition: os_port.h:69
uint_t mpiGetLength(const Mpi *a)
Get the actual length in words.
Definition: mpi.c:129
char char_t
Definition: compiler_port.h:48
error_t mpiMontgomeryMul(Mpi *r, const Mpi *a, const Mpi *b, uint_t k, const Mpi *p, Mpi *t)
Montgomery multiplication.
Definition: mpi.c:1802
uint8_t m
Definition: ndp.h:302
uint8_t n
#define cryptoFreeMem(p)
Definition: crypto.h:714
uint8_t s
@ MPI_FORMAT_BIG_ENDIAN
Definition: mpi.h:61
#define cryptoAllocMem(size)
Definition: crypto.h:709
#define MPI_INT_SIZE
Definition: mpi.h:39
int_t mpiComp(const Mpi *a, const Mpi *b)
Compare two multiple precision integers.
Definition: mpi.c:295
unsigned int uint_t
Definition: compiler_port.h:50
__weak_func error_t mpiMul(Mpi *r, const Mpi *a, const Mpi *b)
Multiple precision multiplication.
Definition: mpi.c:1176
#define osMemset(p, value, length)
Definition: os_port.h:134
uint_t mpiGetBitValue(const Mpi *a, uint_t index)
Get the bit value at the specified index.
Definition: mpi.c:270
error_t mpiAddMod(Mpi *r, const Mpi *a, const Mpi *b, const Mpi *p)
Modular addition.
Definition: mpi.c:1434
error_t mpiMulInt(Mpi *r, const Mpi *a, int_t b)
Multiply a multiple precision integer by an integer.
Definition: mpi.c:1253
int_t mpiCompInt(const Mpi *a, int_t b)
Compare a multiple precision integer with an integer.
Definition: mpi.c:339
error_t mpiGrow(Mpi *r, uint_t size)
Adjust the size of multiple precision integer.
Definition: mpi.c:85
void mpiMulAccCore(uint_t *r, const uint_t *a, int_t m, const uint_t b)
Multiply-accumulate operation.
Definition: mpi.c:1901
__weak_func error_t mpiExpModFast(Mpi *r, const Mpi *a, const Mpi *e, const Mpi *p)
Modular exponentiation (fast calculation)
Definition: mpi.c:1768
@ NO_ERROR
Success.
Definition: error.h:44
uint8_t c
Definition: ndp.h:512
Debugging facilities.
#define arraysize(a)
Definition: os_port.h:73
uint_t mpiGetByteLength(const Mpi *a)
Get the actual length in bytes.
Definition: mpi.c:156
void mpiFree(Mpi *r)
Release a multiple precision integer.
Definition: mpi.c:62
__weak_func error_t mpiCheckProbablePrime(const Mpi *a)
Test whether a number is probable prime.
Definition: mpi.c:558