BigDigits multiple-precision arithmetic source code
BigDigits is a free library of multiple-precision arithmetic routines written in ANSI C to carry out
large natural number calculations as required in cryptography calculations.
This code has been built using the algorithms in Knuth Vol 2
and Menezes as the
primary references.
New version 2.6 released 31 March 2016: see New in Latest Version below.
2020-04-08: If using on a 64-bit non-Windows platform, check the Errata below.
The BigDigits library includes the classical multiple-precision arithmetic algorithms from Knuth: add, subtract, multiply and divide. It also includes modular multiplication, exponentiation and inversion; bit manipulation; number theory functions such as the greatest common divisor and Jacobi symbol; the Rabin-Miller Probabilistic Primality Test procedure from FIPS-186 and ANSI X9.42 to show that a large integer is probably prime; and other utilities to manipulate and handle multiple-precision numbers.



Contents
Recommended reading
- Art of Computer Programming, Volume 2: Seminumerical Algorithms by Donald E. Knuth
- Handbook of Applied Cryptography by Alfred J. Menezes, Paul C. van Oorschot and Scott A. Vanstone
- A Course in Computational Algebraic Number Theory by Henri Cohen
- C Interfaces and Implementations by David R. Hanson
- Elliptic Curves in Cryptography by Ian F. Blake, Gadiel Seroussi and Nigel P. Smart
Affiliate disclosure: we get a small commission for purchases made through the above links
- Interfaces to the library
- Documentation
- New in the latest version
- Changes in earlier versions
- Copyright and commercial use
- Using the BigDigits Software
- Download
- Files in the Download
- How to Compile
- Compiling for the "bd" library
- Creating a static library of BIGD functions
- Using the "mp" functions
- Preprocessor definitions
- Summary of files required for compiling
- Print format specifiers
- Compiling with C++
- Compiling in Linux and other Unix environments
- Compiling on a Mac
- Compiling on a 64-bit machine
- Using with an 8-bit microcontroller
- Compatibility
- Test code
- Overlapping Variables and other "Gotchas"
- A Caution About Negative Numbers
- Random Number Functions
- History
- Programmer's Comments
- Acknowledgements
- References
- Bibliography
- Errata
- Contact
Interfaces to the library
The library has two interfaces:
a complete set of functions called using the BIGD
handle (the "bd" library) which handles memory allocation automatically,
and the set of core functions (the "mp" library) which requires a bit more care but is faster
and has an option (NO_ALLOCS
) to avoid all memory allocation calls completely.
We recommend you use the bd library, which has all the functions you should need. Memory allocation is handled automatically, except the initial creation and final release of resources. There is an option to use Intel Assembler to speed up the critical single-digit multiply and divide operations, and another option to make use of native 64-bit integers.
Check for any errata.
Documentation
- Documention in CHM format (zipped, 117 kB)
- Documention in HTML format (follow the "Files" link):
- bigd.h - the "bd" functions
- bigdigits.h - the "mp" functions
- bigdRand.h - "bd" random number functions
- bigdigitsRand.h - "mp" random number functions
This documentation has been produced using the excellent Doxygen tool.
New in the latest version
Version 2.6 (Released 31 March 2016)- Added new functions to perform modular arithmetic operations, particularly useful for doing computations in a prime field for elliptic curves:
mpModSquare
,bdModSquare
to computea = x^2 (mod m)
mpModSqrt
,bdModSqrt
to compute a square root modulo a prime pmpModHalve
,bdModHalve
to computea = x / 2 (mod p)
-
mpModAdd
,bdModAdd
andmpModSub
,bdModSub
to computew = u ± v (mod m)
quickly for u and v in the range[0, m-1]
. - The function
mpModSpecial
to computeu = v (mod m)
in the special case where v is known to be in the range[0, km]
for k a small integer. This is faster thanmpModulo
.
- Added functions
mpToShort
,bdToShort
to convert a multiprecision integer to a single digit. - Added functions
mpShortIsEqual
,bdShortIsEqual
to check if a multiprecision integer is equal to a single digit. - Added the explicit constant-time comparison functions.
mpEqual_ct()
,bdIsEqual_ct()
mpIsZero_ct()
,bdIsZero_ct()
mpCompare_ct()
,bdCompare_ct()
- Reversed the change in v2.5 to make the comparison functions constant-time by default.
The original comparison functions are now not constant-time, as they were before v2.5. The user must explicitly select
the constant-time variant (ending _ct).
To summarize: the functions
mpEqual
,bdIsEqual
,mpIsZero
,bdIsZero
,mpCompare
andbdCompare
are not constant-time in v2.6.Comment: We realised after we made the change in v2.5 that the majority of comparison operations in our own code were simple checks for special cases and were not comparing secret values, so they did not need the more expensive constant-time operations. There were only a few occasions where we did comparisons on secret values and these were easily found and explicitly changed to do a constant-time comparison. So, sorry, we shouldn't have changed the default behaviour in v2.5. It's back now to how it was before.
- Added the functions
mpCompileTime
,bdCompileTime
to show the time and date of compilation. - Added the function
mpPrintBits
to print the value in 0/1 bit format. - Added the option to specify the global preprocessor definition, say,
MAX_FIXED_BIT_LENGTH=640
when using theNO_ALLOCS
option for the "mp" functions. This allows you to change the maximum fixed length of arrays without editing the core source code. - Increased the size of the fixed-length buffers in the
mpPrint*()
functions to avoid buffer overflows - Fixed an issue with the
mpPrint*()
functions when displaying zero values. - Added the
bdNewVars
andbdFreeVars
functions to allow "bulk" creation and destruction of "bd" variables in a single line of code.
Changes in earlier versions
Version 2.5 (Released 22 October 2015)- Changed license conditions to Mozilla Public License, v. 2.0.
- Fixed
bdShortMult
to catch empty u parameter. - Fixed
bdPrintBits
to print "0" if bit length is zero. - Changed comparison functions to be constant-time by default. Retained the old, quicker functions as
_q
variants. (Changed back in version 2.6). - Added constant-time modular exponentiation functions
bdModExp_ct
andmpModExp_ct
. - Added
mpPrintDecimalSigned
function to display negative numbers ("mp" only).
- Added markup for Doxygen documentation to produce a manual in CHM and HTML formats.
- Updated the "pretty-good" internal RNG functions in bigdRand and bigdigitsRand to improve the seeding of the ANSI X9.31 PRNG algorithm.
- Fixed the
NO_ALLOCS
option so it does not use malloc with thempConv
functions (thanks to Kevin Kramb for pointing this out). - Added the new functions
bdRandomOctets
,bdRandomNumber
andbdQuickRandBits
, which do various contortions you might want to do with random numbers. - Updated the test source code modules.
- Added new functions to compute the integer cube root:
bdCubeRoot
andmpCubeRoot
. - Updated and improved the integer square root functions:
bdSqrt
andmpSqrt
. - Added the more useful print functions with an optional prefix and suffix:
bdPrintHex
prints a bigdigit number in hex format.bdPrintDecimal
prints a bigdigit number in decimal format.
mpPrintHex
andmpPrintDecimal
. - Added the new function
bdPower
to compute the n-th power of a numbery = gn
. Use this with caution as you can quickly run out of memory. It's meant for small exponents like n=3. - Improved the efficiency of the greatest common divisor functions,
bdGcd
andmpGcd
. - Fixed the bug in
mpChs
(thanks to Valery Blazhnov). - Fixed the bug in
mpAlloc
which would fail if you calledbdSetDigit(b, 0)
(thanks to "Radistao").
- Added sliding-window exponentiation to speed up the ModExp functions.
- Added new bd functions:
bdJacobi
computes the Jacobi (Legendre) symbol.bdGetBit
returns the value of a bit at a given position in a number.bdSetBit
sets the bit at a given position in a number.bdNotBits
computes bitwise a = NOT a.
- Added the new functions
mpIsNegative
,mpChs
andmpAbs
in anticipation of adding full signed integer functionality in the next a later version (but see A Caution About Negative Numbers). - Added the
NO_ALLOCS
option to compile the "mp" library without using any memory allocation. - Added the
USE_64WITH32
option to use the 64-bit integer types (long long
) if available on a 32-bit machine. - Improved the zeroisation and destroy macros.
- Moved the type declarations for the exact-width types to a separate include file.
- Re-organised and re-named the ancillary source and include files (yet again!). See How to Compile.
- Deprecated the
bd*Ex
functions in favour of re-named "safe"bd*_s
functions.
- Added support for 64-bit compilers.
- Added fudge so opaque pointers will work with C++ compilers.
- Added new functions:
bdSqrt
computes the `integer' square root, s = floor(sqrt(x)).bdModPowerOf2
computes a = a (mod 2L).bdXorBits
computes bitwise a = b XOR c.bdOrBits
computes bitwise a = b OR c.bdAndBits
computes bitwise a = b AND c.bdVersion
returns version number as an integer = major*1000+minor*100+release*10+uses_asm(0|1).bdRandomBits
generates a random BIGD number up to n bits long (i.e. r = [0, 2n]) using the internal RNG.
- Removed the shift limit on
bdShiftLeft
andbdShiftRight
. This was previously limited to a maximum of 32 bits. It will now shift by any value - well, any value up toSIZE_MAX
. - Re-organised the source and include files for the internal random number generator functions to avoid having to
include the
spRandom
modules when they are not needed. See How to Compile. - Fixed up inconsistencies in the documentation for the
bdConv*
functions. - Made minor change to
mpIsPrime
to comply with ANSI X.42-2003 Annex F.1.1 Range of bases in Miller-Rabin test. - Made minor improvements to the
spBetterRand
function. - Fixed bugs in
bdConvertHex
,bdGetBit
andbdSetBit
. - Made various minor changes in the code to avoid warning messages when compiling with the
-Wall -pedantic
options.
Some speed enhancements have been introduced, including an option to use assembler calls for the critical multiply and divide operations, which makes an appreciable difference in speed on Intel processors. For more details see History. We've also added a more robust random number generator for doing tests.
Copyright and Commercial Use
As of version 2.5 (October 2015) we have re-issued this source code under the Mozilla Public License, v. 2.0.
We have no objection whatsoever to developers using this code in commercial software or as part of an open source project provided the MPL licence remains unchanged in our source code modules.
[For reference, here is the old copyright and licence notice for version 2.4 (April 2013) and earlier.]
If you use it, we'd really appreciate a link back to this page:
<a href="https://www.di-mgt.com.au/bigdigits.html">BigDigits multiple-precision arithmetic source code</a>
Please forward any comments or bug reports to <contact>. The latest version of the source code can be downloaded from www.di-mgt.com.au/bigdigits.html or here.
Using the BigDigits Software
A simple example that multiplies the 48-bit number 0xdeadbeefface (244837814106830) by 2 (we keep this example small so we can reproduce it easily)
#include "bigd.h" int main(void) { BIGD u, v, w; /* Create new BIGD objects */ u = bdNew(); v = bdNew(); w = bdNew(); /* Compute 2 * 0xdeadbeefface */ bdSetShort(u, 2); bdConvFromHex(v, "deadbeefface"); bdMultiply(w, u, v); /* Display the result in hex and decimal */ bdPrintHex("0x", w, "\n"); bdPrintDecimal("", w, "\n"); /* Free all objects we made */ bdFree(&u); bdFree(&v); bdFree(&w); return 0; }
This should display
0x1bd5b7ddff59c 489675628213660which you can verify using the Hex mode of the Windows Calculator in Scientific mode.
All multiple-precision variables need to be declared with a
BIGD v;statement and created with
v = bdNew();When finished, free the resources with
bdFree(&v);
Between calling bdNew()
and bdFree()
,
you can use the variable to carry out any multiple-precision operation using
any function in the bigd.h
interface.
Note the bdFree()
function is passed the address of the variable,
which is also set to NULL at the same time to avoid accidental reuse.
A more compact version
Use the new bdNewVars
and bdFreeVars
functions added in version 2.6 to create and free BIGD objects in bulk
(just don't forget the NULL at the end of each list like, er, we did).
#include "bigd.h" int main(void) { BIGD u, v, w; /* Create new BIGD objects */ bdNewVars(&u, &v, &w, NULL); /* Compute 2 * 0xdeadbeefface */ bdSetShort(u, 2); bdConvFromHex(v, "deadbeefface"); bdMultiply(w, u, v); /* Display the result in hex and decimal */ bdPrintHex("0x", w, "\n"); bdPrintDecimal("", w, "\n"); /* Free all objects we made */ bdFreeVars(&u, &v, &w, NULL); return 0; }
Remarks
In the bd functions, all internal memory allocation is handled automatically. The mp functions require the user to allocate fixed-length arrays of sufficient length: this is faster but less safe.
This code is all written in ANSI C (except the optional ASM and Win32-specific bits). C++ programmers should be able to incorporate this code into their programs with ease - see Compiling with C++. There is a complete list of functions in the documentation. For more details on compiling, see How to Compile.
CAUTION: there are a few "gotchas" - some variables get trashed before use so be careful using the same variable twice in a function call, see Overlapping Variables and other "Gotchas".
Download
- Latest version 2.6: BigDigits-2.6.1.zip (117 kB), released under MPLv2 on 31 March 2016.
- Manual in HTML Help (CHM) format (zipped, 135 kB).
The download contains the full source code files for the multiple-precision library functions, test programs, a Makefile for Linux/cygwin gcc users, and a list of the md5sums for all the files in the library. Check the Errata below for changes.
Files in the Download
- bigd.h
- Include file for BIGD functions, i.e. the "bd" functions with full memory management.
The only external interface you should need. Exports the
BIGD
handle and thebdigit_t
typedef for a single digit. - bigd.c
- Source code for bd functions. This is a wrapper around all the mp
functions in bigdigits.c. It completely hides all references to the mp functions and
the
DIGIT_T
typedef. - bigdigits.h
- Include file for the mp functions. Exports the
DIGIT_T
typedef. - bigdigits.c
- Source code for the mp functions.
- bigdtypes.h
- A common set of macros for exact-width types required by both
bigdigits.h
andbigd.h
. Does not need to be explicity included. - bigdRand.h
- A separate interface for the internal random number functions that rely on spRandom.
Include this if your code makes use of the
bdRandomBits
orbdRandDigit
functions. - bigdRand.c
- Separate source code for the internal RNG functions
bdRandomBits
orbdRandDigit
. - bigdigitsRand.c
- Source code for the internal RNG
spBetterRand
andmpRandomBits
functions. Called by the functions inbigdRand.c()
. We keep this separate so you can substitute your own "best" random function to taste, should you wish. - spExtra.c/.h
- Single-digit versions of some of the multiple-precision functions. Not required for either the BIGD or BIGDIGITS libraries but are included for interest. Useful as an introduction to understand how the mp functions work.
- t_*.c
- Source code for various test functions. See Test Functions.
How to Compile
Compiling for the "bd" library
To compile a program that uses the bigd functions without using the internal RNG functions, create a source file, say,bdcode.c
:
/* bdcode.c */ #include "bigd.h" int main() { BIGD b; b = bdNew(); /* ... */ bdFree(&b); return 0; }and compile your source file
bdcode.c
along with the following files
bigd.c bigd.h bigdigits.c bigdigits.h bigdtypes.h
If your program requires the internal RNG functions (bdRandDigit
or bdRandomBits
) then you will need to add
bigdRand.c bigdRand.h bigdigitsRand.c bigdigitsRand.h
Creating a static library of BIGD functions
We prefer to create a static library of BIGD functions to use in our applications.
The only interface then required is the bigd.h
include file and, if required, bigdRand.h
.
In MSVC, compile the following files to create the static library bigd.lib
:
bigd.c bigdigits.c bigdRand.c bigdigitsRand.c
To use the faster ASM alternative, define the pre-processor constant USE_SPASM or, alternatively, define USE_64WITH32 to use the native 64-bit integers.
To use the resulting library, include the lines
#include "bigd.h #include "bigdRand.h"
in your source code and link to the bigd.lib
library.
All the bd functions can then be called from your program.
Using the "mp" library
To create a program that just uses the mp functions, make a source filempcode.c
:
/* mpcode.c */ #include "bigdigits.h" int main() { DIGIT_T u[10], v[10], w[20]; /* set u and v ... */ mpMultiply(w, u, v, 10); /* ... */ return 0; }and compile
mpcode.c
along with the following files
bigdigits.c bigdigits.h bigdtypes.h
If you use the "mp" RNG functions spBetterRand
or mpRandomBits
, then you need to add
bigdigitsRand.c bigdigitsRand.h
The "mp" functions require you to allocate fixed arrays of digits instead of relying on the automatic allocation routines that the "bd" library uses. Note that some "mp" functions still use `malloc()` to create temporary arrays.
To avoid any use of memory allocation functions,
define the NO_ALLOCS
preprocessor macro. The "mp" functions will then use fixed arrays and
not use any memory allocation functions at all.
You can set the value of MAX_FIXED_DIGITS
in the
bigdigits.h
file. It is currently set to cope with integers of up to 8192 bits.
If you don't define NO_ALLOCS
, then this maximum is ignored and temporary arrays will be allocated
as required.
Preprocessor definitions
In both the "bd" and "mp" libraries you can define one of the preprocessor macros
- USE_SPASM to use the faster ASM routines (on Intel processors where the __asm option is available with your compiler); or
- USE_64WITH32 to use the 64-bit integers if available (e.g.
long long
)
In the "mp" library only, you can define NO_ALLOCS to avoid using the memory allocation functions completely. If this is set, then you can set the MAX_FIXED_DIGITS definition for the maximum integer size you want to handle (the current default is 8192 bits).
You can detect which of these preprocessor definitions has been applied by using the bdVersion
or mpVersion
functions (they return the same results).
2300 = none used 2301 = USE_SPASM 2302 = USE_64WITH32 2305 = NO_ALLOCS 2306 = NO_ALLOCS+USE_SPASM 2307 = NO_ALLOCS+USE_64WITH32
The USE_SPASM option takes precedence over USE_64WITH32. The option NO_ALLOCS cannot be used with the "bd" library, only with "mp".
Summary of files required for compiling
The internal RNG functions are bdRandomBits
or bdRandDigit
in the "bd" library; and
mpRandomBits
or spBetterRand
in the "mp" library.
Library | No internal RNG functions | Uses an internal RNG function |
---|---|---|
"mp" only |
bigdigits.c bigdigits.h bigdtypes.h |
bigdigits.c bigdigitRand.c bigdigits.h bigdtypes.h bigdigitRand.h |
Full "bd" |
bigd.c bigdigits.c bigd.h bigdigits.h bigdtypes.h |
bigd.c bigdigits.c bigdRand.c bigdigitRand.c bigd.h bigdigits.h bigdtypes.h bigdRand.h bigdigitRand.h |
Print format specifiers
We use BIGD-compatible versions of the C99 format specifiers, so you can replace print expressions of the form
bdigit_t b; /* ... */ printf("Value=%lx\n", b);with the more generic
printf("Value=%" PRIxBIGD "\n"), b);where
PRIxBIGD
is our own version of the C99 PRIx32
.
Compiling with C++
This code is written in ANSI C which is (almost, but not quite) a sub-set of C++.
The opaque pointer techniques in bigd.h
and bigd.c
based on
[HANSON97]
rely on the separate namespaces for
structures and type names,
which fall over when compiled with C++.
There is a `fudge' in the code thanks to Ian Tree which fixes this.
See the example in t_bdCPP.cpp
.
Do not try to change the extensions of the core source code files
like bigd.c
to .cpp
. This will not work.
Compiling in Linux and other Unix environments
See the Makefile in the download for Linux/Cygwin environments using gcc. These are not necessarily examples of the most efficient Makefiles possible.
Compiling on a Mac
2012-01-23: Arno Hautala reports success compiling and running the included test suite on Mac OS X.
I did have to make one addition to line 36 of bigdtypes.h to add "
|| defined(__APPLE__)
" so that stdint.h would be properly included.
Update: This change has been incorporated in version 2.4.
Compiling on a 64-bit machine
The code has been changed to use the C99 <stdint.h> types uint32_t and uint16_t instead of "unsigned long" and "unsigned short". These are implemented via some convoluted preprocessor instructions that seem to work on most test systems we've tried. You may need to tweak these. Thanks to Terry R. Friedrichsen and Allen Zhao for their help on this.
We'd be interested to hear from anyone who has changed the code to use 64-bit digits directly.
That is, changed all the uint32_t's to uint64_t, and uint16_t's to uint32_t; changed
BITS_PER_DIGIT to 64, MAX_DIGIT to UINT64_MAX, PRIu32 to PRIu64,
HIBITMASK to 0x8000 0000 0000 0000 UL
,
and so forth.
Using with an 8-bit microcontroller
Alfredo Ortega has modified the BigDigits code to run on an 8-bit microprocessor. See Implementation of the RSA algorithm in a 8-bit microcontroller.
Alternative: try compiling the "mp" library with the NO_ALLOCS
preprocessor macro defined. This removes all
calls to memory allocation functions and just uses fixed arrays.
You can define the maximum length of array you want to use.
Adding PKCS#1 and elliptic curve cryptography on binary and prime fields
2009-12-10: Dang Nguyen Duc has extended an earlier version of the BigDigits library to a more complete library for implementing public key cryptography. In particular, adding PKCS#1, elliptic curve on binary and prime fields. He has posted his code to libmcrypto: Math Library for Crypto Students.
“ Thank you very much for writing very clean and easy-to-follow code. That was much more useful for my study than using GNU MP. Along the way, I have added a few functions to your bigdigits library... I still attempt to maintain the cleanness and easy-to-follow property of your original bigdigit library without optimization or inline assembly. A study tool is my top priority. ”
Compatibility
The source code and tests have been compiled and verified on Windows 7/XP/Win98 operating systems using MSVC++ versions 5, 6, 7, 8 (Express 2005) and 9 (MSVC++ 2008); Borland C++ 5.5.1 for Win32; Cygwin using gcc version 3.3.3 on W2K; and Linux (Fedora Core 1) using gcc version 3.3.2 (Red Hat Linux 3.3.2-1).Allen Zhao reports success compiling version 2.1 in both 32- and 64-bit modes on Solaris 8, HPUX 11i/PA-RISC; AIX 4.3.3 (32bit) and AIX 5.1(32/64); IRIX 6.5 (32bit) and IRIX64 6.5 (32/64); RH Linux 7.2 / RHEL 3 / RHEL 4 with both gcc and Intel icc (v9.1) for x86/x86_64/ia64 (any combination of 32-bit/64-bit compilation). Thanks Allen. Gustavo del Dago reports success (21 Nov 2006) compiling on the OS/400 platform. Arno Hautala reports success (23 January 2012) compiling on Mac OS X.
Test Code
Test source file names start with "t_". They test various aspects of the code with perhaps somewhat boring tests.- t_bdSimple.c
- A simple example that computes 2 * 0xdeadbeefface (244837814106830) and displays the result.
- t_bdTest.c
- A selection of tests that tests most functions in the bd library at least once. Uses assert checking to catch errors.
- t_bdRSA.c
- Generates a new 1025-bit RSA key pair and tests both the RSA encryption and signing operations on random data blocks. It carries out decryption using both the inversion and the CRT methods and compares the time taken. (We use 1025 bits to test the generator for odd numbers.)
- t_bdRSA1.c
- A variant on t_bdRSA.c. This example takes a short message from the command line and encrypts it in the RSA encryption block using PKCS1-v1.5 encoding and a freshly-generated 1024-bit RSA key. Again, this just tests the functions. You could re-write this to use a fixed RSA key pair.
- t_bdRsaCrack.c
- Code to "crack" RSA if the user has done encryption in a certain (poor) way.
- t_bdRsaFactorN.c
- How to factor the RSA modulus given the secret exponent d.
- t_bdRSA_blinded.c
- Carry out RSA calculations using constant-time exponentiation and blinding to provide resistance against simple power attacks (SPA) and differential power analysis (DPA).
- t_bdDSA.c
- Reproduces the test vector of the Digital Signature Algorithm (DSA) from Appendix 5 FIPS PUB 186-2 "Digital Signature Standard (DSS)", 27 January 2000.
- t_bdRDSA.c
- Reproduces some test vectors of the rDSA algorithm from ANSI X9.31-1998 "Digital Signatures Using Reversible Public Key Cryptography for the Financial Services Industry (rDSA)", September 9, 1998. Appendix D.1 Odd e = 3 with 1024-bit n and D.3 Odd e =3 with 2048-bit n.
- t_bdRandomOctets.c
-
Examples using the
bdRandomOctets
andbdRandomNumber
functions. - t_bdQuickRandBits.c
-
Examples comparing the
bdQuickRandBits
andbdRandomBits
functions. - t_mpTest.c
- Carries out some tests on the mp library functions.
- t_mpRSA508.c
- Uses the mp library functions to reproduce the example of RSA encryption and signature using 508-bit test vectors from "Some Examples of the PKCS Standards", An RSA Laboratories Technical Note, Burton S. Kaliski Jr., 1 November 1993.
- t_mpJacobi.c
- Carries out some tests of the
mpJacobi
function using an example taken from ANSI X9.31 Appendix D.5.2. - t_spExtras.c
- Carries out some tests on the functions in
spExtras.c
. - t_bdCPP.cpp
- A simple example in C++. See also Compiling with C++.
Overlapping Variables and other "Gotchas"
The mp functions are rather raw and user-unfriendly, but faster. They assume that the user has
allocated memory of the correct length and understands the quirks of the various functions.
For example, mpMultiply
requires the two multiplicands
x and y
to be of the same ndigits
length
but the product p must be twice that length.
The user is meant to know and understand this or face the consequences of a GPF error when there is overflow.
For speed and efficiency reasons, some functions require that certain source parameters do not overlap with the destination parameter. The user is also meant to remember this and know the safe combinations. There are assert checks included in the mp source code to trap such errors.
The bd functions are designed to handle the memory allocation issues behind the scenes
but the overlap issues are still present for some functions for speed reasons.
For example, bdAdd(x, x, y)
is OK, but bdAdd(y, x, y)
is not,
and
bdMultiply(y, x, y)
is OK, but bdMultiply(x, x, y)
is not
(we took inspiration in consistency from certain C library and Windows API functions here).
If in doubt, use separate variables or use the bdAdd_s()
and bdMultiply_s()
functions
which are safe but slower (or read the instructions more carefully!).
Think of this as similar to the difference between the memcpy
and memmove
functions in the standard C string library.
A Caution About Negative Numbers
The BigDigits library is designed to work with the natural numbers ℕ; that is, the non-negative integers 0,1,2,...
.
Negative numbers will "work" with the fixed-array "mp" functions if you are happy using twos-complement arithmetic.
There are a few "experimental" functions included
(mpIsNegative()
, mpChs()
and mpAbs()
) intended for signed integers, but so far we have not developed
that very far. We added the function mpPrintDecimalSigned()
in v2.5, but only for the fixed-length "mp" functions.
But be careful with the memory-managed "bd" functions. If you go below zero with the "bd" functions the results are undefined.
You should
either make sure you always do cut-off subtraction that yields a non-negative result
(i.e. to compute a - b, first check that b ≤ a; if not return zero or flag an error)
or check the "carry/overflow" result from the bdSubtract
or bdDecrement
operation and flag an error if this is not zero.
t_bdTests.c
/* Cope, sort of, with `negative' numbers */
printf("Go ``negative''...\n");
bdSetZero(u);
overflow = bdShortSub(w, u, 2);
pr_msg("w=0-2=", w);
printf("overflow=%" PRIuBIGD "\n", overflow);
overflow = bdIncrement(w);
pr_msg("w+1=", w);
printf("overflow=%" PRIuBIGD "\n", overflow);
overflow = bdIncrement(w);
pr_msg("(Note error!) w+1=", w);
printf("overflow=%" PRIuBIGD "\n", overflow);
which gives the output
Go ``negative''... w=0-2=fffffffe overflow=1 w+1=ffffffff overflow=0 (Note error!) w+1=00000001 00000000 overflow=1Thanks to Taylor Francis for reminding us about this issue.
How you could do cut-off subtraction so the result is never negative. Check the return value to see if you underflowed.
int cutoff_subtraction(BIGD w, BIGD u, BIGD v) { int carry = 0; if (bdCompare(u, v) < 0) { bdSetZero(w); carry = 1; } else { bdSubtract(w, u, v); } return carry; }
Random Number Functions
There are three occasions when you might need to use a random number generator in connection with the BIGDIGITs package.
- To generate a few pseudo-random numbers in a limited range to carry out a few tests.
- To generate a large set of random digits to use in a large number of tests.
- To generate a cryptographically-secure random number in a security application.
- A few pseudo-random numbers for simple tests
-
In the first case, you can use the standard ANSI
rand()
function from thestdlib
library and seed it withsrand(time(NULL))
or better,srand(time(NULL) ^ clock())
You need to be careful if you want a number greater thanRAND_MAX
, which can be as small as 0x7fff on some systems. We have used this method in thebdSetRandTest
function and thespSimpleRand(DIGIT_T lower, DIGIT_T upper)
function which usesrand()
to generate a random number byte-by-byte in the range between `lower' and `upper'. This is adequate for the odd few random numbers in tests that don't need to be secure. - To fill up many arrays of digits with random values for tests
-
In the second case, using
rand()
is not suitable because it starts repeating itself too soon. For this purpose we have added thespBetterRand()
function, which is based on not-quite-cryptographically-secure techniques. This is called from the BIGD library via the synonymbdRandDigit()
and is used bybdRandomBits
, both inbigdRand.h
. This function is perfectly adequate for generating large numbers of random digits to use in tests with an infinitesimally small chance of repetition or non-random behaviour. The output should pass all FIPS-140-2 tests for randomness. The function does, however, rely on static variables and so should be used with care in a multi-threaded environment. We have kept the code for this function in a separate module so it can easily be replaced if you wish. - Cryptographically-secure RNGs
-
In the third case, you would need to call a cryptographically-secure function. In the
bdRandomSeeded()
andbdGeneratePrime()
functions - which are intended to be used in secure cryptography packages - you pass a reference to a separate RNG call-back function which can be as cryptographically-secure as you wish. Using a separate function enables you to handle multi-threaded issues and use the secure RNG of your choice. We have included a couple of example functions (which aren't secure themselves) in the tests for generating RSA keys (remember we're demonstrating a multiple-precision package, not a cryptographic module). The function has a required formatBD_RANDFUNC
with parameters to pass the output buffer, required size and seed, which your own function would need to match.
You are welcome to use the included random number functions at your own pleasure and risk. We don't want to enter protracted discussions about how "random" these functions are, but we will respond to obvious errors. Please carry out your own tests to verify the functions provided or substitute your own ones that you trust.
History
By David IrelandI know there are several packages out there that do this already. I wanted to write my own set partly to prove to myself I could do it, and partly so we could have a set we could use in our own commercial applications without restriction. You are welcome to use it in your applications, including commercial ones, provided you adhere to the copyright conditions.
I wrote the original BigDigits Version 1 back in 2001. I did most of the research and drafting while sitting overnight with my son who was 8 years old and in hospital at the time (he's fine now - look!). The original design goals were:-
- Use only "pure" ANSI C
- Keep it independent of machine endianness
- Make sure it's reasonably understandable as stand-alone code
- Avoid having to declare temporary storage where possible
The original requirement to avoid allocating memory was not practical, though. Even worse, having
a fixed maximum size constant hidden away in the include file was a cause of problems for some users.
So we've added dynamic memory allocation for those "mp" functions that needed temporary variables
of varying length, and then for the "bd" functions we added full internal memory handling.
We've tested these functions in a variety of situations for memory leakage. Provided you use the
bdFree()
before you exit, there should be no memory leaks.
We experimented with using Montgomery Exponentiation for an improved ModExp
function, but
we ended up up with little or no improvement and some very hairy code, so we dropped it.
However, adding sliding-window exponentiation in version 2.2 made some significant improvements.
Adding an mpSquare
function gives about a 40% speed improvement against using
multiplication, but again at the cost of some complicated code. Using this
improves the ModExp
function overall by about 10% on the average.
The one experiment that did produce a significant increase in speed was to use Assembler
for the frequently-called single-digit multiply and divide functions. This resulted in
up to a 40% speed improvement in tests. The original ANSI C version is still the default.
In version 2.2 you can alternatively use the USE_64WITH32
preprocessor macro to make use of
64-bit types like long long
if available.
Programmer's Comments
Given that this code has been in the public domain for years without any significant bugs being reported, we are very reluctant to make big changes to the core "mp" code. Some points of interest on Version 2:-- malloc
- Using
malloc
to allocate memory for temporary variables makes for an easier life, but causes a performance hit for those functions that call them repeatedly. Note how we have used our own internal functions likempModuloTemp
to re-use the same memory for temporary variables for the memory-intensivempModExp
function. - Assertion checking
-
The "bd" functions do assertion checks that a NULL
BIGD
parameter has not been passed. The "mp" functions do assertion checks to ensure that overlap problems will not occur. We recommend that you leave assertion checking turned on, even in your production code. - mpDivide
- Writing the Division algorithm (Algorithm D, Knuth, p272) in elegant C code was tricky.
This is probably the most tested part of the library. A suitable example to test the
"Add Back" at step D6 of the algorithm for b=232 is
u = (
7fffffff 80000001 00000000 00000000
)32, v = (80000000 80000002 00000005
)32. See the solution to 4.3.1 Example 22 in Knuth, p625. Fortunately the values given by DEK for b=216 readily convert upwards for 32-bits. - spDivide
- We are pretty certain that the "Add Back" step will never happen in our single-precision
ANSI C
spDivide
function but we are not brave enough to remove the code that handles it. We think this requires at least 3 digits in the divisor, v, and we only have two. At worst there are 15 unnecessary machine instructions; at best it catches something that 'theoretically' shouldn't happen. Comments - and proofs - would be welcome. - ASM
- Using assembler to handle the spDivide and spMultiply functions on a 32-bit machine produced a
significant improvement in performance (as you would expect when you get the processor
to use its own instruction instead of your own arcane version in ANSI C).
You still need to be careful when doing unsigned division
in assembler as the DIV operation will throw a #DE exception if there is overflow
(i.e. when the resulting quotient is larger than 32-bits). See how we handle it in the code in
spASM.c
. This exception crops up very rarely in normal operations - typical tests could go for hours without running into it. Why doesn't Intel just set an overflow flag instead of throwing an exception?
Acknowledgements
Thank you in particular to the following people for their suggestions and advice that have contributed to the development of this code:- Jesse Chisholm
- Sandy Dunlop
- Terry R. Friedrichsen
- Peter Hebert
- James Ku
- Ian Tree
- Rickard Westman
- Allen Zhao
- Valery Blazhnov
- "Radistao"
- Kevin Kramb
References
Sources of the algorithms used in BigDigits.- [KNUTH] Donald E. Knuth, Art of Computer Programming, Volume 2: Seminumerical Algorithms, 3rd ed, Addison-Wesley, 1998.
- [MENEZES] Alfred J. Menezes, Paul C. van Oorschot, Scott A Vanstone, Handbook of Applied Cryptography (Discrete Mathematics and Its Applications), CRC Press, 1997, <http://www.cacr.math.uwaterloo.ca/hac/>.
- [COHEN] Henri Cohen, A Course in Computational Algebraic Number Theory, Springer, 1996.
- [SCHNEIER] Bruce Schneier, Applied Cryptography: Protocols, Algorithms, and Source Code in C, 2nd ed, Wiley, 1996.
- [CLARK] Andrew Clark, Implementation Issues in Cryptography, <http://sky.fit/qut.edu.au/~aclark/itn536/>, (accessed March 2001).
- [GOURDON] Xavier Gourdon and Pascal Sebah, Arbitrary precision computation, <http://numbers.computation.free.fr/Constants/constants.html>, (accessed March 2001).
- [FIPS-186-2] Digital Signature Standard (DSS), FIPS PUB 186-2, U.S. Department of Commerce/National Institute of Standards and Technology, 2000 (latest PDF).
- [ANSIX942] ANSI X9.42-2003 Public Key Cryptography for the Financial Services Industry: Agreement of Symmetric Keys Using Discrete Logarithm Cryptography, American National Standards Institute, 2003.
- [CORON] Jean-Sebastian Coron, Resistance Against Differential Power Analysis for Elliptic Curve Cryptosystems, Cryptographic Hardware and Embedded Systems, Lecture Notes in Computer Science Volume 1717, 1999, pp 292-302. (PDF)
- [BLAKE] Blake, Ian F., Seroussi, G., N.P. Smart, Elliptic Curves in Cryptography, Cambridge University Press, 1999.
Bibliography
Other sources we used in devising tests and other techniques used in the code.- [ANSIX931] ANSI X9.31-1998 Digital Signatures using Reversible Public Key Cryptography for the Financial Services Industry (rDSA), American National Standards Institute, 1998.
- [FERGUSON03] Niels Ferguson and Bruce Schneier, Practical Cryptography, John Wiley, 2003. This is virtually identical to the later edition of Cryptography Engineering
- [HANSON97] David R. Hanson, C Interfaces and Implementations: Techniques for Creating Reusable Software, Addison-Wesley, 1997, <http://www.cs.princeton.edu/software/cii/>.
- [INTEL2] The Intel Architecture Software Developer's Manual, Volume 2: Instruction Set Reference, Intel Corporation, 1999.
- [ISO99] ISO/IEC 9899:1999 Programming languages - C, 2nd ed., ISO/IEC, 1999.
- [PKCS1] PKCS#1, RSA Cryptography Standard, RSA Laboratories, Version 2.1, June 2002 (republished as RFC3447).
- [PKCS-EX] Burton S Kaliski, Jr, Some Examples of the PKCS Standards, RSA Laboratories Technical Note, Nov 1993.
Errata
Last updated: 8 April 20201. Typo in bigdigits.c for uint64_t
HAVE_C99INCLUDES
to HAVE_C99_INCLUDES
.
190 /* Make sure we have a uint64_t available */ 191 #if defined (_WIN32) || defined(WIN32) 192 typedef unsigned __int64 uint64_t; 193 #elif !defined(HAVE_C99INCLUDES) && !defined(HAVE_SYS_TYPES) 194 typedef unsigned long long int uint64_t; 195 #endifIn fact, the lines 190 to 195 are redundant and can be removed.
2. Error in t_bdTest.c
NULL
to the end of the lists in lines 167 and 1052.
Before:
167 bdNewVars(&u, &v, &w, &p, &q, &r, &x, &y); // WRONG! ... 1052 bdFreeVars(&u, &v, &w, &p, &q, &r, &x, &y); // WRONG!After:
167 bdNewVars(&u, &v, &w, &p, &q, &r, &x, &y, NULL); ... 1052 bdFreeVars(&u, &v, &w, &p, &q, &r, &x, &y, NULL);
Contact
Any comments, feedback, questions: please send us a message.
This page last updated 8 April 2020