[e6f27aa56f] Add files I missed in my previous commit

10 files changed, 755 insertions, 0 deletions

diff --git a/libtommath/bn_mp_export.c b/libtommath/bn_mp_export.c
new file mode 100644
index 0000000..ac4c2f9
--- /dev/null
+++ b/libtommath/bn_mp_export.c
@@ -0,0 +1,88 @@
+#include <tommath_private.h>
+#ifdef BN_MP_EXPORT_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* based on gmp's mpz_export.
+ * see http://gmplib.org/manual/Integer-Import-and-Export.html
+ */
+int mp_export(void* rop, size_t* countp, int order, size_t size, 
+                                int endian, size_t nails, mp_int* op) {
+	int result;
+	size_t odd_nails, nail_bytes, i, j, bits, count;
+	unsigned char odd_nail_mask;
+
+	mp_int t;
+
+	if ((result = mp_init_copy(&t, op)) != MP_OKAY) {
+		return result;
+	}
+
+	if (endian == 0) {
+		union {
+			unsigned int i;
+			char c[4];
+		} lint;
+		lint.i = 0x01020304;
+
+		endian = (lint.c[0] == 4) ? -1 : 1;
+	}
+
+	odd_nails = (nails % 8);
+	odd_nail_mask = 0xff;
+	for (i = 0; i < odd_nails; ++i) {
+		odd_nail_mask ^= (1 << (7 - i));
+	}
+	nail_bytes = nails / 8;
+
+	bits = mp_count_bits(&t);
+	count = (bits / ((size * 8) - nails)) + (((bits % ((size * 8) - nails)) != 0) ? 1 : 0);
+
+	for (i = 0; i < count; ++i) {
+		for (j = 0; j < size; ++j) {
+			unsigned char* byte = (
+				(unsigned char*)rop + 
+				(((order == -1) ? i : ((count - 1) - i)) * size) +
+				((endian == -1) ? j : ((size - 1) - j))
+			);
+
+			if (j >= (size - nail_bytes)) {
+				*byte = 0;
+				continue;
+			}
+
+			*byte = (unsigned char)((j == ((size - nail_bytes) - 1)) ? (t.dp[0] & odd_nail_mask) : (t.dp[0] & 0xFF));
+
+			if ((result = mp_div_2d(&t, ((j == ((size - nail_bytes) - 1)) ? (8 - odd_nails) : 8), &t, NULL)) != MP_OKAY) {
+				mp_clear(&t);
+				return result;
+			}
+		}
+	}
+
+	mp_clear(&t);
+
+	if (countp != NULL) {
+		*countp = count;
+	}
+
+	return MP_OKAY;
+}
+
+#endif
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */
diff --git a/libtommath/bn_mp_expt_d_ex.c b/libtommath/bn_mp_expt_d_ex.c
new file mode 100644
index 0000000..649d224
--- /dev/null
+++ b/libtommath/bn_mp_expt_d_ex.c
@@ -0,0 +1,83 @@
+#include <tommath_private.h>
+#ifdef BN_MP_EXPT_D_EX_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* calculate c = a**b  using a square-multiply algorithm */
+int mp_expt_d_ex (mp_int * a, mp_digit b, mp_int * c, int fast)
+{
+  int     res;
+  unsigned int x;
+
+  mp_int  g;
+
+  if ((res = mp_init_copy (&g, a)) != MP_OKAY) {
+    return res;
+  }
+
+  /* set initial result */
+  mp_set (c, 1);
+
+  if (fast != 0) {
+    while (b > 0) {
+      /* if the bit is set multiply */
+      if ((b & 1) != 0) {
+        if ((res = mp_mul (c, &g, c)) != MP_OKAY) {
+          mp_clear (&g);
+          return res;
+        }
+      }
+
+      /* square */
+      if (b > 1) {
+        if ((res = mp_sqr (&g, &g)) != MP_OKAY) {
+          mp_clear (&g);
+          return res;
+        }
+      }
+
+      /* shift to next bit */
+      b >>= 1;
+    }
+  }
+  else {
+    for (x = 0; x < DIGIT_BIT; x++) {
+      /* square */
+      if ((res = mp_sqr (c, c)) != MP_OKAY) {
+        mp_clear (&g);
+        return res;
+      }
+
+      /* if the bit is set multiply */
+      if ((b & (mp_digit) (((mp_digit)1) << (DIGIT_BIT - 1))) != 0) {
+        if ((res = mp_mul (c, &g, c)) != MP_OKAY) {
+           mp_clear (&g);
+           return res;
+        }
+      }
+
+      /* shift to next bit */
+      b <<= 1;
+    }
+  } /* if ... else */
+
+  mp_clear (&g);
+  return MP_OKAY;
+}
+#endif
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */
diff --git a/libtommath/bn_mp_get_long.c b/libtommath/bn_mp_get_long.c
new file mode 100644
index 0000000..7c3d0fe
--- /dev/null
+++ b/libtommath/bn_mp_get_long.c
@@ -0,0 +1,41 @@
+#include <tommath_private.h>
+#ifdef BN_MP_GET_LONG_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* get the lower unsigned long of an mp_int, platform dependent */
+unsigned long mp_get_long(mp_int * a)
+{
+  int i;
+  unsigned long res;
+
+  if (a->used == 0) {
+     return 0;
+  }
+
+  /* get number of digits of the lsb we have to read */
+  i = MIN(a->used,(int)(((sizeof(unsigned long) * CHAR_BIT) + DIGIT_BIT - 1) / DIGIT_BIT)) - 1;
+
+  /* get most significant digit of result */
+  res = DIGIT(a,i);
+
+#if (ULONG_MAX != 0xffffffffuL) || (DIGIT_BIT < 32)
+  while (--i >= 0) {
+    res = (res << DIGIT_BIT) | DIGIT(a,i);
+  }
+#endif
+  return res;
+}
+#endif
diff --git a/libtommath/bn_mp_get_long_long.c b/libtommath/bn_mp_get_long_long.c
new file mode 100644
index 0000000..4b959e6
--- /dev/null
+++ b/libtommath/bn_mp_get_long_long.c
@@ -0,0 +1,41 @@
+#include <tommath_private.h>
+#ifdef BN_MP_GET_LONG_LONG_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* get the lower unsigned long long of an mp_int, platform dependent */
+unsigned long long mp_get_long_long (mp_int * a)
+{
+  int i;
+  unsigned long long res;
+
+  if (a->used == 0) {
+     return 0;
+  }
+
+  /* get number of digits of the lsb we have to read */
+  i = MIN(a->used,(int)(((sizeof(unsigned long long) * CHAR_BIT) + DIGIT_BIT - 1) / DIGIT_BIT)) - 1;
+
+  /* get most significant digit of result */
+  res = DIGIT(a,i);
+
+#if DIGIT_BIT < 64
+  while (--i >= 0) {
+    res = (res << DIGIT_BIT) | DIGIT(a,i);
+  }
+#endif
+  return res;
+}
+#endif
diff --git a/libtommath/bn_mp_import.c b/libtommath/bn_mp_import.c
new file mode 100644
index 0000000..dd4b8e6
--- /dev/null
+++ b/libtommath/bn_mp_import.c
@@ -0,0 +1,73 @@
+#include <tommath_private.h>
+#ifdef BN_MP_IMPORT_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* based on gmp's mpz_import.
+ * see http://gmplib.org/manual/Integer-Import-and-Export.html
+ */
+int mp_import(mp_int* rop, size_t count, int order, size_t size, 
+                            int endian, size_t nails, const void* op) {
+	int result;
+	size_t odd_nails, nail_bytes, i, j;
+	unsigned char odd_nail_mask;
+
+	mp_zero(rop);
+
+	if (endian == 0) {
+		union {
+			unsigned int i;
+			char c[4];
+		} lint;
+		lint.i = 0x01020304;
+
+		endian = (lint.c[0] == 4) ? -1 : 1;
+	}
+
+	odd_nails = (nails % 8);
+	odd_nail_mask = 0xff;
+	for (i = 0; i < odd_nails; ++i) {
+		odd_nail_mask ^= (1 << (7 - i));
+	}
+	nail_bytes = nails / 8;
+
+	for (i = 0; i < count; ++i) {
+		for (j = 0; j < (size - nail_bytes); ++j) {
+			unsigned char byte = *(
+					(unsigned char*)op + 
+					(((order == 1) ? i : ((count - 1) - i)) * size) +
+					((endian == 1) ? (j + nail_bytes) : (((size - 1) - j) - nail_bytes))
+				);
+
+			if (
+				(result = mp_mul_2d(rop, ((j == 0) ? (8 - odd_nails) : 8), rop)) != MP_OKAY) {
+				return result;
+			}
+
+			rop->dp[0] |= (j == 0) ? (byte & odd_nail_mask) : byte;
+			rop->used  += 1;
+		}
+	}
+
+	mp_clamp(rop);
+
+	return MP_OKAY;
+}
+
+#endif
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */
diff --git a/libtommath/bn_mp_n_root_ex.c b/libtommath/bn_mp_n_root_ex.c
new file mode 100644
index 0000000..79d1dfb
--- /dev/null
+++ b/libtommath/bn_mp_n_root_ex.c
@@ -0,0 +1,132 @@
+#include <tommath_private.h>
+#ifdef BN_MP_N_ROOT_EX_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* find the n'th root of an integer
+ *
+ * Result found such that (c)**b <= a and (c+1)**b > a
+ *
+ * This algorithm uses Newton's approximation
+ * x[i+1] = x[i] - f(x[i])/f'(x[i])
+ * which will find the root in log(N) time where
+ * each step involves a fair bit.  This is not meant to
+ * find huge roots [square and cube, etc].
+ */
+int mp_n_root_ex (mp_int * a, mp_digit b, mp_int * c, int fast)
+{
+  mp_int  t1, t2, t3;
+  int     res, neg;
+
+  /* input must be positive if b is even */
+  if (((b & 1) == 0) && (a->sign == MP_NEG)) {
+    return MP_VAL;
+  }
+
+  if ((res = mp_init (&t1)) != MP_OKAY) {
+    return res;
+  }
+
+  if ((res = mp_init (&t2)) != MP_OKAY) {
+    goto LBL_T1;
+  }
+
+  if ((res = mp_init (&t3)) != MP_OKAY) {
+    goto LBL_T2;
+  }
+
+  /* if a is negative fudge the sign but keep track */
+  neg     = a->sign;
+  a->sign = MP_ZPOS;
+
+  /* t2 = 2 */
+  mp_set (&t2, 2);
+
+  do {
+    /* t1 = t2 */
+    if ((res = mp_copy (&t2, &t1)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+
+    /* t2 = t1 - ((t1**b - a) / (b * t1**(b-1))) */
+
+    /* t3 = t1**(b-1) */
+    if ((res = mp_expt_d_ex (&t1, b - 1, &t3, fast)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+
+    /* numerator */
+    /* t2 = t1**b */
+    if ((res = mp_mul (&t3, &t1, &t2)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+
+    /* t2 = t1**b - a */
+    if ((res = mp_sub (&t2, a, &t2)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+
+    /* denominator */
+    /* t3 = t1**(b-1) * b  */
+    if ((res = mp_mul_d (&t3, b, &t3)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+
+    /* t3 = (t1**b - a)/(b * t1**(b-1)) */
+    if ((res = mp_div (&t2, &t3, &t3, NULL)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+
+    if ((res = mp_sub (&t1, &t3, &t2)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+  }  while (mp_cmp (&t1, &t2) != MP_EQ);
+
+  /* result can be off by a few so check */
+  for (;;) {
+    if ((res = mp_expt_d_ex (&t1, b, &t2, fast)) != MP_OKAY) {
+      goto LBL_T3;
+    }
+
+    if (mp_cmp (&t2, a) == MP_GT) {
+      if ((res = mp_sub_d (&t1, 1, &t1)) != MP_OKAY) {
+         goto LBL_T3;
+      }
+    } else {
+      break;
+    }
+  }
+
+  /* reset the sign of a first */
+  a->sign = neg;
+
+  /* set the result */
+  mp_exch (&t1, c);
+
+  /* set the sign of the result */
+  c->sign = neg;
+
+  res = MP_OKAY;
+
+LBL_T3:mp_clear (&t3);
+LBL_T2:mp_clear (&t2);
+LBL_T1:mp_clear (&t1);
+  return res;
+}
+#endif
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */
diff --git a/libtommath/bn_mp_set_long.c b/libtommath/bn_mp_set_long.c
new file mode 100644
index 0000000..281fce7
--- /dev/null
+++ b/libtommath/bn_mp_set_long.c
@@ -0,0 +1,24 @@
+#include <tommath_private.h>
+#ifdef BN_MP_SET_LONG_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* set a platform dependent unsigned long int */
+MP_SET_XLONG(mp_set_long, unsigned long)
+#endif
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */
diff --git a/libtommath/bn_mp_set_long_long.c b/libtommath/bn_mp_set_long_long.c
new file mode 100644
index 0000000..3c4b01a
--- /dev/null
+++ b/libtommath/bn_mp_set_long_long.c
@@ -0,0 +1,24 @@
+#include <tommath_private.h>
+#ifdef BN_MP_SET_LONG_LONG_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://libtom.org
+ */
+
+/* set a platform dependent unsigned long long int */
+MP_SET_XLONG(mp_set_long_long, unsigned long long)
+#endif
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */
diff --git a/libtommath/bn_mp_sqrtmod_prime.c b/libtommath/bn_mp_sqrtmod_prime.c
new file mode 100644
index 0000000..968729e
--- /dev/null
+++ b/libtommath/bn_mp_sqrtmod_prime.c
@@ -0,0 +1,124 @@
+#include <tommath_private.h>
+#ifdef BN_MP_SQRTMOD_PRIME_C
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ */
+
+/* Tonelli-Shanks algorithm
+ * https://en.wikipedia.org/wiki/Tonelli%E2%80%93Shanks_algorithm
+ * https://gmplib.org/list-archives/gmp-discuss/2013-April/005300.html
+ *
+ */
+
+int mp_sqrtmod_prime(mp_int *n, mp_int *prime, mp_int *ret)
+{
+  int res, legendre;
+  mp_int t1, C, Q, S, Z, M, T, R, two;
+  mp_digit i;
+
+  /* first handle the simple cases */
+  if (mp_cmp_d(n, 0) == MP_EQ) {
+    mp_zero(ret);
+    return MP_OKAY;
+  }
+  if (mp_cmp_d(prime, 2) == MP_EQ)                              return MP_VAL; /* prime must be odd */
+  if ((res = mp_jacobi(n, prime, &legendre)) != MP_OKAY)        return res;
+  if (legendre == -1)                                           return MP_VAL; /* quadratic non-residue mod prime */
+
+  if ((res = mp_init_multi(&t1, &C, &Q, &S, &Z, &M, &T, &R, &two, NULL)) != MP_OKAY) {
+	return res;
+  }
+
+  /* SPECIAL CASE: if prime mod 4 == 3
+   * compute directly: res = n^(prime+1)/4 mod prime
+   * Handbook of Applied Cryptography algorithm 3.36
+   */
+  if ((res = mp_mod_d(prime, 4, &i)) != MP_OKAY)                goto cleanup;
+  if (i == 3) {
+    if ((res = mp_add_d(prime, 1, &t1)) != MP_OKAY)             goto cleanup;
+    if ((res = mp_div_2(&t1, &t1)) != MP_OKAY)                  goto cleanup;
+    if ((res = mp_div_2(&t1, &t1)) != MP_OKAY)                  goto cleanup;
+    if ((res = mp_exptmod(n, &t1, prime, ret)) != MP_OKAY)      goto cleanup;
+    res = MP_OKAY;
+    goto cleanup;
+  }
+
+  /* NOW: Tonelli-Shanks algorithm */
+
+  /* factor out powers of 2 from prime-1, defining Q and S as: prime-1 = Q*2^S */
+  if ((res = mp_copy(prime, &Q)) != MP_OKAY)                    goto cleanup;
+  if ((res = mp_sub_d(&Q, 1, &Q)) != MP_OKAY)                   goto cleanup;
+  /* Q = prime - 1 */
+  mp_zero(&S);
+  /* S = 0 */
+  while (mp_iseven(&Q) != MP_NO) {
+    if ((res = mp_div_2(&Q, &Q)) != MP_OKAY)                    goto cleanup;
+    /* Q = Q / 2 */
+    if ((res = mp_add_d(&S, 1, &S)) != MP_OKAY)                 goto cleanup;
+    /* S = S + 1 */
+  }
+
+  /* find a Z such that the Legendre symbol (Z|prime) == -1 */
+  if ((res = mp_set_int(&Z, 2)) != MP_OKAY)                     goto cleanup;
+  /* Z = 2 */
+  while(1) {
+    if ((res = mp_jacobi(&Z, prime, &legendre)) != MP_OKAY)     goto cleanup;
+    if (legendre == -1) break;
+    if ((res = mp_add_d(&Z, 1, &Z)) != MP_OKAY)                 goto cleanup;
+    /* Z = Z + 1 */
+  }
+
+  if ((res = mp_exptmod(&Z, &Q, prime, &C)) != MP_OKAY)         goto cleanup;
+  /* C = Z ^ Q mod prime */
+  if ((res = mp_add_d(&Q, 1, &t1)) != MP_OKAY)                  goto cleanup;
+  if ((res = mp_div_2(&t1, &t1)) != MP_OKAY)                    goto cleanup;
+  /* t1 = (Q + 1) / 2 */
+  if ((res = mp_exptmod(n, &t1, prime, &R)) != MP_OKAY)         goto cleanup;
+  /* R = n ^ ((Q + 1) / 2) mod prime */
+  if ((res = mp_exptmod(n, &Q, prime, &T)) != MP_OKAY)          goto cleanup;
+  /* T = n ^ Q mod prime */
+  if ((res = mp_copy(&S, &M)) != MP_OKAY)                       goto cleanup;
+  /* M = S */
+  if ((res = mp_set_int(&two, 2)) != MP_OKAY)                   goto cleanup;
+
+  res = MP_VAL;
+  while (1) {
+    if ((res = mp_copy(&T, &t1)) != MP_OKAY)                    goto cleanup;
+    i = 0;
+    while (1) {
+      if (mp_cmp_d(&t1, 1) == MP_EQ) break;
+      if ((res = mp_exptmod(&t1, &two, prime, &t1)) != MP_OKAY) goto cleanup;
+      i++;
+    }
+    if (i == 0) {
+      if ((res = mp_copy(&R, ret)) != MP_OKAY)                  goto cleanup;
+      res = MP_OKAY;
+      goto cleanup;
+    }
+    if ((res = mp_sub_d(&M, i, &t1)) != MP_OKAY)                goto cleanup;
+    if ((res = mp_sub_d(&t1, 1, &t1)) != MP_OKAY)               goto cleanup;
+    if ((res = mp_exptmod(&two, &t1, prime, &t1)) != MP_OKAY)   goto cleanup;
+    /* t1 = 2 ^ (M - i - 1) */
+    if ((res = mp_exptmod(&C, &t1, prime, &t1)) != MP_OKAY)     goto cleanup;
+    /* t1 = C ^ (2 ^ (M - i - 1)) mod prime */
+    if ((res = mp_sqrmod(&t1, prime, &C)) != MP_OKAY)           goto cleanup;
+    /* C = (t1 * t1) mod prime */
+    if ((res = mp_mulmod(&R, &t1, prime, &R)) != MP_OKAY)       goto cleanup;
+    /* R = (R * t1) mod prime */
+    if ((res = mp_mulmod(&T, &C, prime, &T)) != MP_OKAY)        goto cleanup;
+    /* T = (T * C) mod prime */
+    mp_set(&M, i);
+    /* M = i */
+  }
+
+cleanup:
+  mp_clear_multi(&t1, &C, &Q, &S, &Z, &M, &T, &R, &two, NULL);
+  return res;
+}
+
+#endif
diff --git a/libtommath/tommath_private.h b/libtommath/tommath_private.h
new file mode 100644
index 0000000..d23c333
--- /dev/null
+++ b/libtommath/tommath_private.h
@@ -0,0 +1,125 @@
+/* LibTomMath, multiple-precision integer library -- Tom St Denis
+ *
+ * LibTomMath is a library that provides multiple-precision
+ * integer arithmetic as well as number theoretic functionality.
+ *
+ * The library was designed directly after the MPI library by
+ * Michael Fromberger but has been written from scratch with
+ * additional optimizations in place.
+ *
+ * The library is free for all purposes without any express
+ * guarantee it works.
+ *
+ * Tom St Denis, tstdenis82@gmail.com, http://math.libtomcrypt.com
+ */
+#ifndef TOMMATH_PRIV_H_
+#define TOMMATH_PRIV_H_
+
+#include <tommath.h>
+#include <ctype.h>
+
+#ifndef MIN
+#define MIN(x,y) (((x) < (y)) ? (x) : (y))
+#endif
+
+#ifndef MAX
+#define MAX(x,y) (((x) > (y)) ? (x) : (y))
+#endif
+
+#ifdef __cplusplus
+extern "C" {
+
+/* C++ compilers don't like assigning void * to mp_digit * */
+#define  OPT_CAST(x)  (x *)
+
+#else
+
+/* C on the other hand doesn't care */
+#define  OPT_CAST(x)
+
+#endif
+
+/* define heap macros */
+#if 0
+#ifndef XMALLOC
+   /* default to libc stuff */
+   #define XMALLOC  malloc
+   #define XFREE    free
+   #define XREALLOC realloc
+   #define XCALLOC  calloc
+#else
+   /* prototypes for our heap functions */
+   extern void *XMALLOC(size_t n);
+   extern void *XREALLOC(void *p, size_t n);
+   extern void *XCALLOC(size_t n, size_t s);
+   extern void XFREE(void *p);
+#endif
+#endif
+
+/* lowlevel functions, do not call! */
+int s_mp_add(mp_int *a, mp_int *b, mp_int *c);
+int s_mp_sub(mp_int *a, mp_int *b, mp_int *c);
+#define s_mp_mul(a, b, c) s_mp_mul_digs(a, b, c, (a)->used + (b)->used + 1)
+int fast_s_mp_mul_digs(mp_int *a, mp_int *b, mp_int *c, int digs);
+int s_mp_mul_digs(mp_int *a, mp_int *b, mp_int *c, int digs);
+int fast_s_mp_mul_high_digs(mp_int *a, mp_int *b, mp_int *c, int digs);
+int s_mp_mul_high_digs(mp_int *a, mp_int *b, mp_int *c, int digs);
+int fast_s_mp_sqr(mp_int *a, mp_int *b);
+int s_mp_sqr(mp_int *a, mp_int *b);
+int mp_karatsuba_mul(mp_int *a, mp_int *b, mp_int *c);
+int mp_toom_mul(mp_int *a, mp_int *b, mp_int *c);
+int mp_karatsuba_sqr(mp_int *a, mp_int *b);
+int mp_toom_sqr(mp_int *a, mp_int *b);
+int fast_mp_invmod(mp_int *a, mp_int *b, mp_int *c);
+int mp_invmod_slow (mp_int * a, mp_int * b, mp_int * c);
+int fast_mp_montgomery_reduce(mp_int *x, mp_int *n, mp_digit rho);
+int mp_exptmod_fast(mp_int *G, mp_int *X, mp_int *P, mp_int *Y, int redmode);
+int s_mp_exptmod (mp_int * G, mp_int * X, mp_int * P, mp_int * Y, int redmode);
+void bn_reverse(unsigned char *s, int len);
+
+extern const char *mp_s_rmap;
+
+/* Fancy macro to set an MPI from another type.
+ * There are several things assumed:
+ *  x is the counter and unsigned
+ *  a is the pointer to the MPI
+ *  b is the original value that should be set in the MPI.
+ */
+#define MP_SET_XLONG(func_name, type)                    \
+int func_name (mp_int * a, type b)                       \
+{                                                        \
+  unsigned int  x;                                       \
+  int           res;                                     \
+                                                         \
+  mp_zero (a);                                           \
+                                                         \
+  /* set four bits at a time */                          \
+  for (x = 0; x < (sizeof(type) * 2u); x++) {            \
+    /* shift the number up four bits */                  \
+    if ((res = mp_mul_2d (a, 4, a)) != MP_OKAY) {        \
+      return res;                                        \
+    }                                                    \
+                                                         \
+    /* OR in the top four bits of the source */          \
+    a->dp[0] |= (b >> ((sizeof(type) * 8u) - 4u)) & 15u; \
+                                                         \
+    /* shift the source up to the next four bits */      \
+    b <<= 4;                                             \
+                                                         \
+    /* ensure that digits are not clamped off */         \
+    a->used += 1;                                        \
+  }                                                      \
+  mp_clamp (a);                                          \
+  return MP_OKAY;                                        \
+}
+
+#ifdef __cplusplus
+   }
+#endif
+
+#endif
+
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */