From 2ab527b7e068df4b0e6a82b674c6fcf8d24935fd Mon Sep 17 00:00:00 2001
From: Daniel Mendler <mail@daniel-mendler.de>
Date: Mon, 10 Sep 2018 18:51:26 +0200
Subject: [PATCH] add mp_get_double, mp_set_double

---
 bn_mp_get_double.c | 32 +++++++++++++++++++++++++++
 bn_mp_set_double.c | 54 ++++++++++++++++++++++++++++++++++++++++++++++
 demo/demo.c        | 42 ++++++++++++++++++++++++++++++++++++
 tommath.h          |  6 ++++++
 4 files changed, 134 insertions(+)
 create mode 100644 bn_mp_get_double.c
 create mode 100644 bn_mp_set_double.c

diff --git a/bn_mp_get_double.c b/bn_mp_get_double.c
new file mode 100644
index 0000000..542993d
--- /dev/null
+++ b/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$ */
diff --git a/bn_mp_set_double.c b/bn_mp_set_double.c
new file mode 100644
index 0000000..0a5f771
--- /dev/null
+++ b/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$ */
diff --git a/demo/demo.c b/demo/demo.c
index b5af727..b62954f 100644
--- a/demo/demo.c
+++ b/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) {
diff --git a/tommath.h b/tommath.h
index 9cec473..04274eb 100644
--- a/tommath.h
+++ b/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);