--- /dev/null
+/* TomsFastMath, a fast ISO C bignum library.
+ *
+ * This project is meant to fill in where LibTomMath
+ * falls short. That is speed ;-)
+ *
+ * This project is public domain and free for all purposes.
+ *
+ * Tom St Denis, tomstdenis@gmail.com
+ */
+#include <tfm_private.h>
+
+/* setups the montgomery reduction */
+int fp_montgomery_setup(fp_int *a, fp_digit *rho)
+{
+ fp_digit x, b;
+
+/* fast inversion mod 2**k
+ *
+ * Based on the fact that
+ *
+ * XA = 1 (mod 2**n) => (X(2-XA)) A = 1 (mod 2**2n)
+ * => 2*X*A - X*X*A*A = 1
+ * => 2*(1) - (1) = 1
+ */
+ b = a->dp[0];
+
+ if ((b & 1) == 0) {
+ return FP_VAL;
+ }
+
+ x = (((b + 2) & 4) << 1) + b; /* here x*a==1 mod 2**4 */
+ x *= 2 - b * x; /* here x*a==1 mod 2**8 */
+ x *= 2 - b * x; /* here x*a==1 mod 2**16 */
+ x *= 2 - b * x; /* here x*a==1 mod 2**32 */
+#ifdef FP_64BIT
+ x *= 2 - b * x; /* here x*a==1 mod 2**64 */
+#endif
+
+ /* rho = -1/m mod b */
+ *rho = (((fp_word) 1 << ((fp_word) DIGIT_BIT)) - ((fp_word)x));
+
+ return FP_OKAY;
+}
+
+
+/* $Source$ */
+/* $Revision$ */
+/* $Date$ */