1 files changed, 73 insertions, 4 deletions

diff --git a/libraries/integer-gmp/cbits/wrappers.c b/libraries/integer-gmp/cbits/wrappers.c
index c99c0176a4..11e5179323 100644
--- a/libraries/integer-gmp/cbits/wrappers.c
+++ b/libraries/integer-gmp/cbits/wrappers.c
@@ -11,6 +11,7 @@
 
 #include "HsFFI.h"
 #include "MachDeps.h"
+#include "HsIntegerGmp.h"
 
 #include <assert.h>
 #include <stdbool.h>
@@ -285,9 +286,9 @@ integer_gmp_mpn_gcd(mp_limb_t r[],
  * reconstructed).
  *
  * g must have space for exactly gn=min(xn,yn) limbs.
- * s must have space for at least xn limbs.
+ * s must have space for at least yn limbs.
  *
- * return value: signed 'sn' of {sp,sn}
+ * return value: signed 'sn' of {sp,sn} where |sn| >= 1
  */
 mp_size_t
 integer_gmp_gcdext(mp_limb_t s0[], mp_limb_t g0[],
@@ -304,15 +305,25 @@ integer_gmp_gcdext(mp_limb_t s0[], mp_limb_t g0[],
 
   mpz_gcdext (g, s, NULL, x, y);
 
+  // g must be positive (0 <= gn).
+  // According to the docs for mpz_gcdext(), we have:
+  //       g < min(|y|/2|s|, |x|/2|t|)
+  // -->   g < min(|y|, |x|)
+  // -->   gn <= min(yn, xn)
+  // <->   gn <= gn0
   const mp_size_t gn = g[0]._mp_size;
   assert(0 <= gn && gn <= gn0);
   memset(g0, 0, gn0*sizeof(mp_limb_t));
   memcpy(g0, g[0]._mp_d, gn*sizeof(mp_limb_t));
   mpz_clear (g);
 
+  // According to the docs for mpz_gcdext(), we have:
+  //       |s| < |y| / 2g
+  // -->   |s| < |y|         (note g > 0)
+  // -->   sn <= yn
   const mp_size_t ssn = s[0]._mp_size;
   const mp_size_t sn  = mp_size_abs(ssn);
-  assert(sn <= xn);
+  assert(sn <= mp_size_abs(yn));
   memcpy(s0, s[0]._mp_d, sn*sizeof(mp_limb_t));
   mpz_clear (s);
 
@@ -626,7 +637,7 @@ integer_gmp_powm(mp_limb_t rp[], // result
   }
 
   const mpz_t b = CONST_MPZ_INIT(bp, mp_limb_zero_p(bp,bn) ? 0 : bn);
-  const mpz_t e = CONST_MPZ_INIT(ep, mp_limb_zero_p(ep,en) ? 0 : en);
+  const mpz_t e = CONST_MPZ_INIT(ep, en);
   const mpz_t m = CONST_MPZ_INIT(mp, mn);
 
   mpz_t r;
@@ -687,6 +698,64 @@ integer_gmp_powm_word(const mp_limb_t b0, // base
   return integer_gmp_powm1(&b0, !!b0, &e0, !!e0, m0);
 }
 
+/* version of integer_gmp_powm() based on mpz_powm_sec
+ *
+ * With GMP 5.0 or later execution time depends on size of arguments
+ * and size of result.
+ *
+ * 'M' must be odd and 'E' non-negative.
+ */
+mp_size_t
+integer_gmp_powm_sec(mp_limb_t rp[], // result
+                     const mp_limb_t bp[], const mp_size_t bn, // base
+                     const mp_limb_t ep[], const mp_size_t en, // exponent
+                     const mp_limb_t mp[], const mp_size_t mn) // mod
+{
+  assert(!mp_limb_zero_p(mp,mn));
+  assert(mp[0] & 1);
+
+  if ((mn == 1 || mn == -1) && mp[0] == 1) {
+    rp[0] = 0;
+    return 1;
+  }
+
+  if (mp_limb_zero_p(ep,en)) {
+    rp[0] = 1;
+    return 1;
+  }
+
+  assert(en > 0);
+
+  const mpz_t b = CONST_MPZ_INIT(bp, mp_limb_zero_p(bp,bn) ? 0 : bn);
+  const mpz_t e = CONST_MPZ_INIT(ep, en);
+  const mpz_t m = CONST_MPZ_INIT(mp, mn);
+
+  mpz_t r;
+  mpz_init (r);
+
+#if HAVE_SECURE_POWM == 0
+  mpz_powm(r, b, e, m);
+#else
+  mpz_powm_sec(r, b, e, m);
+#endif
+
+  const mp_size_t rn = r[0]._mp_size;
+
+  if (rn) {
+    assert(0 < rn && rn <= mn);
+    memcpy(rp, r[0]._mp_d, rn*sizeof(mp_limb_t));
+  }
+
+  mpz_clear (r);
+
+  if (!rn) {
+    rp[0] = 0;
+    return 1;
+  }
+
+  return rn;
+}
+
 
 /* wrapper around mpz_invert()
  *