add mp_get_double, mp_set_double

bn_mp_get_double.c

@@ -0,0 +1,32 @@
+#include "tommath_private.h"
+#ifdef BN_MP_GET_DOUBLE_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.
+ */
+
+double mp_get_double(const mp_int *a)
+{
+   int i;
+   double d = 0, fac = 1;
+   for (i = 0; i < DIGIT_BIT; ++i) {
+      fac *= 2;
+   }
+   for (i = USED(a); i --> 0;) {
+      d = d * fac + (double)DIGIT(a, i);
+   }
+   return mp_isneg(a) ? -d : d;
+}
+#endif
+
+/* ref:         $Format:%D$ */
+/* git commit:  $Format:%H$ */
+/* commit time: $Format:%ai$ */

bn_mp_set_double.c

@@ -0,0 +1,54 @@
+#include "tommath_private.h"
+#ifdef BN_MP_SET_DOUBLE_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.
+ */
+
+#if defined(__STDC_IEC_559__) || defined(__GCC_IEC_559)
+int mp_set_double(mp_int *a, double d)
+{
+   uint64_t frac;
+   int exp, res;
+   union {
+      double   dbl;
+      uint64_t bits;
+   } cast;
+   cast.dbl = d;
+
+   exp = (int)(cast.bits >> 52) & 0x7FF;
+   frac = (cast.bits & ((1ULL << 52) - 1)) | (1ULL << 52);
+
+   if (exp == 0x7FF) { /* +-inf, NaN */
+      return MP_VAL;
+   }
+   exp -= 1023 + 52;
+
+   res = mp_set_long_long(a, frac);
+   if (res != MP_OKAY) {
+      return res;
+   }
+
+   res = exp < 0 ? mp_div_2d(a, -exp, a, 0) : mp_mul_2d(a, exp, a);
+   if ((cast.bits >> 63) && !mp_iszero(a)) {
+      SIGN(a) = MP_NEG;
+   }
+
+   return MP_OKAY;
+}
+#else
+#  warning "mp_set_double implementation is only available on platforms with IEEE754 floating point format"
+#endif
+#endif
+
+/* ref:         $Format:%D$ */
+/* git commit:  $Format:%H$ */
+/* commit time: $Format:%ai$ */

demo/demo.c

@@ -413,6 +413,48 @@ int main(void)
       }
    }
 
+   // test mp_get_double/mp_set_double
+#if defined(__STDC_IEC_559__) || defined(__GCC_IEC_559)
+   printf("\n\nTesting: mp_get_double");
+   if (mp_set_double(&a, +1.0/0.0) != MP_VAL) {
+      printf("\nmp_set_double should return MP_VAL for +inf");
+      return EXIT_FAILURE;
+   }
+   if (mp_set_double(&a, -1.0/0.0) != MP_VAL) {
+      printf("\nmp_set_double should return MP_VAL for -inf");
+      return EXIT_FAILURE;
+   }
+   if (mp_set_double(&a, +0.0/0.0) != MP_VAL) {
+      printf("\nmp_set_double should return MP_VAL for NaN");
+      return EXIT_FAILURE;
+   }
+   if (mp_set_double(&a, -0.0/0.0) != MP_VAL) {
+      printf("\nmp_set_double should return MP_VAL for NaN");
+      return EXIT_FAILURE;
+   }
+
+   for (i = 0; i < 1000; ++i) {
+      int tmp = rand();
+      double dbl = (double)tmp * rand() + 1;
+      if (mp_set_double(&a, dbl) != MP_OKAY) {
+         printf("\nmp_set_double() failed");
+         return EXIT_FAILURE;
+      }
+      if (dbl != mp_get_double(&a)) {
+         printf("\nmp_get_double() bad result!");
+         return EXIT_FAILURE;
+      }
+      if (mp_set_double(&a, -dbl) != MP_OKAY) {
+         printf("\nmp_set_double() failed");
+         return EXIT_FAILURE;
+      }
+      if (-dbl != mp_get_double(&a)) {
+         printf("\nmp_get_double() bad result!");
+         return EXIT_FAILURE;
+      }
+   }
+#endif
+
    // test mp_get_int
    printf("\n\nTesting: mp_get_int");
    for (i = 0; i < 1000; ++i) {

tommath.h

@@ -201,6 +201,9 @@ void mp_zero(mp_int *a);
 /* set to a digit */
 void mp_set(mp_int *a, mp_digit b);
 
+/* set a double */
+int mp_set_double(mp_int *a, double b);
+
 /* set a 32-bit const */
 int mp_set_int(mp_int *a, unsigned long b);
 
@@ -210,6 +213,9 @@ int mp_set_long(mp_int *a, unsigned long b);
 /* set a platform dependent unsigned long long value */
 int mp_set_long_long(mp_int *a, unsigned long long b);
 
+/* get a double */
+double mp_get_double(const mp_int *a);
+
 /* get a 32-bit value */
 unsigned long mp_get_int(const mp_int *a);