Merge GH #1750 Improve BigInt const time behavior

8 files changed, 349 insertions, 157 deletions

diff --git a/src/lib/math/bigint/big_ops2.cpp b/src/lib/math/bigint/big_ops2.cpp
index df672d035..2fd775ccf 100644
--- a/src/lib/math/bigint/big_ops2.cpp
+++ b/src/lib/math/bigint/big_ops2.cpp
@@ -114,9 +114,41 @@ BigInt& BigInt::mod_add(const BigInt& s, const BigInt& mod, secure_vector<word>&
    if(this->is_negative() || s.is_negative() || mod.is_negative())
       throw Invalid_Argument("BigInt::mod_add expects all arguments are positive");
 
-   // TODO add optimized version of this
-   *this += s;
-   this->reduce_below(mod, ws);
+   BOTAN_DEBUG_ASSERT(*this < mod);
+   BOTAN_DEBUG_ASSERT(s < mod);
+
+   /*
+   t + s or t + s - p == t - (p - s)
+
+   So first compute ws = p - s
+
+   Then compute t + s and t - ws
+
+   If t - ws does not borrow, then that is the correct valued
+   */
+
+   const size_t mod_sw = mod.sig_words();
+   BOTAN_ARG_CHECK(mod_sw > 0, "BigInt::mod_add modulus must be positive");
+
+   this->grow_to(mod_sw);
+   s.grow_to(mod_sw);
+
+   // First mod_sw for p - s, 2*mod_sw for bigint_addsub workspace
+   if(ws.size() < 3*mod_sw)
+      ws.resize(3*mod_sw);
+
+   word borrow = bigint_sub3(&ws[0], mod.data(), mod_sw, s.data(), mod_sw);
+   CT::unpoison(borrow);
+   BOTAN_ASSERT_NOMSG(borrow == 0);
+
+   // Compute t - ws
+   borrow = bigint_sub3(&ws[mod_sw], this->data(), mod_sw, &ws[0], mod_sw);
+
+   // Compute t + s
+   bigint_add3_nc(&ws[mod_sw*2], this->data(), mod_sw, s.data(), mod_sw);
+
+   CT::conditional_copy_mem(borrow, &ws[0], &ws[mod_sw*2], &ws[mod_sw], mod_sw);
+   set_words(&ws[0], mod_sw);
 
    return (*this);
    }
@@ -126,50 +158,30 @@ BigInt& BigInt::mod_sub(const BigInt& s, const BigInt& mod, secure_vector<word>&
    if(this->is_negative() || s.is_negative() || mod.is_negative())
       throw Invalid_Argument("BigInt::mod_sub expects all arguments are positive");
 
-   const size_t mod_sw = mod.sig_words();
-
+   // We are assuming in this function that *this and s are no more than mod_sw words long
    BOTAN_DEBUG_ASSERT(*this < mod);
    BOTAN_DEBUG_ASSERT(s < mod);
 
-   // We are assuming here that *this and s are no more than mod_sw words long
-   const size_t t_w = std::min(mod_sw, size());
-   const size_t s_w = std::min(mod_sw, s.size());
-
-   /*
-   TODO make this const time
-   */
-
-   int32_t relative_size = bigint_cmp(data(), t_w, s.data(), s_w);
+   const size_t mod_sw = mod.sig_words();
 
-   if(relative_size >= 0)
-      {
-      /*
-      this >= s in which case just subtract
+   this->grow_to(mod_sw);
+   s.grow_to(mod_sw);
 
-      Here s_w might be > t_w because these values are just based on
-      the size of the buffer. But we know that because *this < s, then
-      this->sig_words() must be <= s.sig_words() so set the size of s
-      to the minimum of t and s words.
-      */
-      BOTAN_DEBUG_ASSERT(sig_words() <= s.sig_words());
-      bigint_sub2(mutable_data(), t_w, s.data(), std::min(t_w, s_w));
-      }
-   else
-      {
-       // Otherwise we must sub s and then add p (or add (p - s) as here)
+   if(ws.size() < mod_sw)
+      ws.resize(mod_sw);
 
-      if(ws.size() < mod_sw)
-         ws.resize(mod_sw);
+   // is t < s or not?
+   const word is_lt = bigint_ct_is_lt(data(), mod_sw, s.data(), mod_sw);
 
-      word borrow = bigint_sub3(ws.data(), mod.data(), mod_sw, s.data(), s_w);
-      BOTAN_ASSERT_NOMSG(borrow == 0);
+   // ws = p - s
+   word borrow = bigint_sub3(ws.data(), mod.data(), mod_sw, s.data(), mod_sw);
+   CT::unpoison(borrow);
+   BOTAN_ASSERT_NOMSG(borrow == 0);
 
-      if(size() < mod_sw)
-         grow_to(mod_sw);
-
-      word carry = bigint_add2_nc(mutable_data(), size(), ws.data(), mod_sw);
-      BOTAN_ASSERT_NOMSG(carry == 0);
-      }
+   // Compute either (t - s) or (t + (p - s)) depending on mask
+   word carry = bigint_cnd_addsub(is_lt, mutable_data(), ws.data(), s.data(), mod_sw);
+   CT::unpoison(carry);
+   BOTAN_ASSERT_NOMSG(carry == 0);
 
    return (*this);
    }
@@ -185,24 +197,18 @@ BigInt& BigInt::rev_sub(const word y[], size_t y_sw, secure_vector<word>& ws)
 
    const size_t x_sw = this->sig_words();
 
-   const int32_t relative_size = bigint_cmp(y, y_sw, this->data(), x_sw);
+   // TODO use bigint_sub_abs or a new variant of it
 
    ws.resize(std::max(y_sw, x_sw) + 1);
    clear_mem(ws.data(), ws.size());
 
-   if(relative_size < 0)
+   word borrow = bigint_sub3(ws.data(), y, y_sw, this->data(), x_sw);
+
+   if(borrow)
       {
       bigint_sub3(ws.data(), this->data(), x_sw, y, y_sw);
       this->flip_sign();
       }
-   else if(relative_size == 0)
-      {
-      ws.clear();
-      }
-   else if(relative_size > 0)
-      {
-      bigint_sub3(ws.data(), y, y_sw, this->data(), x_sw);
-      }
 
    this->swap_reg(ws);
 
diff --git a/src/lib/math/bigint/bigint.cpp b/src/lib/math/bigint/bigint.cpp
index bd89be7fb..2cb9394ce 100644
--- a/src/lib/math/bigint/bigint.cpp
+++ b/src/lib/math/bigint/bigint.cpp
@@ -140,6 +140,33 @@ int32_t BigInt::cmp(const BigInt& other, bool check_signs) const
                      other.data(), other.sig_words());
    }
 
+bool BigInt::is_equal(const BigInt& other) const
+   {
+   if(this->sign() != other.sign())
+      return false;
+
+   return bigint_ct_is_eq(this->data(), this->sig_words(),
+                          other.data(), other.sig_words());
+   }
+
+bool BigInt::is_less_than(const BigInt& other) const
+   {
+   if(this->is_negative() && other.is_positive())
+      return true;
+
+   if(this->is_positive() && other.is_negative())
+      return false;
+
+   if(other.is_negative() && this->is_negative())
+      {
+      return !bigint_ct_is_lt(other.data(), other.sig_words(),
+                              this->data(), this->sig_words(), true);
+      }
+
+   return bigint_ct_is_lt(this->data(), this->sig_words(),
+                          other.data(), other.sig_words());
+   }
+
 void BigInt::encode_words(word out[], size_t size) const
    {
    const size_t words = sig_words();
@@ -155,25 +182,21 @@ size_t BigInt::Data::calc_sig_words() const
    {
    size_t sig = m_reg.size();
 
-#if 0
-   // Const time, but slower ...
-
-   word seen_only_zeros = MP_WORD_MASK;
    word sub = 1;
 
    for(size_t i = 0; i != m_reg.size(); ++i)
       {
       const word w = m_reg[m_reg.size() - i - 1];
-      seen_only_zeros &= CT::is_zero(w);
-      sub &= seen_only_zeros;
-
+      sub &= CT::is_zero(w);
       sig -= sub;
       }
 
-#else
-   while(sig && (m_reg[sig-1] == 0))
-      sig--;
-#endif
+   /*
+   * This depends on the data so is poisoned, but unpoison it here as
+   * later conditionals are made on the size.
+   */
+   CT::unpoison(sig);
+
    return sig;
    }
 
@@ -221,10 +244,14 @@ uint32_t BigInt::to_u32bit() const
 void BigInt::set_bit(size_t n)
    {
    const size_t which = n / BOTAN_MP_WORD_BITS;
-   const word mask = static_cast<word>(1) << (n % BOTAN_MP_WORD_BITS);
-   if(which >= size()) grow_to(which + 1);
 
-   m_data.set_word_at(which, m_data.get_word_at(which) | mask);
+   if(which >= size())
+      {
+      grow_to(which + 1);
+      }
+
+   const word mask = static_cast<word>(1) << (n % BOTAN_MP_WORD_BITS);
+   m_data.set_word_at(which, word_at(which) | mask);
    }
 
 /*
@@ -233,9 +260,12 @@ void BigInt::set_bit(size_t n)
 void BigInt::clear_bit(size_t n)
    {
    const size_t which = n / BOTAN_MP_WORD_BITS;
-   const word mask = ~(static_cast<word>(1) << (n % BOTAN_MP_WORD_BITS));
+
    if(which < size())
-      m_data.set_word_at(which, m_data.get_word_at(which) & mask);
+      {
+      const word mask = ~(static_cast<word>(1) << (n % BOTAN_MP_WORD_BITS));
+      m_data.set_word_at(which, word_at(which) & mask);
+      }
    }
 
 size_t BigInt::bytes() const
@@ -254,7 +284,10 @@ size_t BigInt::bits() const
       return 0;
 
    const size_t full_words = words - 1;
-   return (full_words * BOTAN_MP_WORD_BITS + high_bit(word_at(full_words)));
+   const size_t bits = (full_words * BOTAN_MP_WORD_BITS + high_bit(word_at(full_words)));
+   // Need to unpoison due to high_bit not being const time
+   CT::unpoison(bits);
+   return bits;
    }
 
 /*
@@ -303,6 +336,7 @@ void BigInt::reduce_below(const BigInt& p, secure_vector<word>& ws)
       {
       word borrow = bigint_sub3(ws.data(), data(), p_words + 1, p.data(), p_words);
 
+      //CT::unpoison(borrow); // fixme
       if(borrow)
          break;
 
diff --git a/src/lib/math/bigint/bigint.h b/src/lib/math/bigint/bigint.h
index 64e408798..1de3f7bc5 100644
--- a/src/lib/math/bigint/bigint.h
+++ b/src/lib/math/bigint/bigint.h
@@ -314,7 +314,7 @@ class BOTAN_PUBLIC_API(2,0) BigInt final
      *
      * Assumes that *this is (if anything) only slightly larger than
      * mod and performs repeated subtractions. It should not be used if
-     * *this is much larger than mod, instead of modulo operator.
+     * *this is much larger than mod, instead use modulo operator.
      */
      void reduce_below(const BigInt& mod, secure_vector<word> &ws);
 
@@ -334,6 +334,20 @@ class BOTAN_PUBLIC_API(2,0) BigInt final
      int32_t cmp(const BigInt& n, bool check_signs = true) const;
 
      /**
+     * Compare this to another BigInt
+     * @param n the BigInt value to compare with
+     * @result true if this == n or false otherwise
+     */
+     bool is_equal(const BigInt& n) const;
+
+     /**
+     * Compare this to another BigInt
+     * @param n the BigInt value to compare with
+     * @result true if this < n or false otherwise
+     */
+     bool is_less_than(const BigInt& n) const;
+
+     /**
      * Compare this to an integer
      * @param n the value to compare with
      * @result if (this<n) return -1, if (this>n) return 1, if both
@@ -562,7 +576,7 @@ class BOTAN_PUBLIC_API(2,0) BigInt final
      * Increase internal register buffer to at least n words
      * @param n new size of register
      */
-     void grow_to(size_t n) { m_data.grow_to(n); }
+     void grow_to(size_t n) const { m_data.grow_to(n); }
 
      /**
      * Resize the vector to the minimum word size to hold the integer, or
@@ -896,7 +910,7 @@ class BOTAN_PUBLIC_API(2,0) BigInt final
                  }
               }
 
-           void grow_to(size_t n)
+           void grow_to(size_t n) const
               {
               if(n > size())
                  {
@@ -954,7 +968,7 @@ class BOTAN_PUBLIC_API(2,0) BigInt final
 
            size_t calc_sig_words() const;
 
-           secure_vector<word> m_reg;
+           mutable secure_vector<word> m_reg;
            mutable size_t m_sig_words = sig_words_npos;
         };
 
@@ -986,17 +1000,17 @@ BigInt BOTAN_PUBLIC_API(2,0) operator>>(const BigInt& x, size_t n);
 * Comparison Operators
 */
 inline bool operator==(const BigInt& a, const BigInt& b)
-   { return (a.cmp(b) == 0); }
+   { return a.is_equal(b); }
 inline bool operator!=(const BigInt& a, const BigInt& b)
-   { return (a.cmp(b) != 0); }
+   { return !a.is_equal(b); }
 inline bool operator<=(const BigInt& a, const BigInt& b)
    { return (a.cmp(b) <= 0); }
 inline bool operator>=(const BigInt& a, const BigInt& b)
    { return (a.cmp(b) >= 0); }
 inline bool operator<(const BigInt& a, const BigInt& b)
-   { return (a.cmp(b) < 0); }
+   { return a.is_less_than(b); }
 inline bool operator>(const BigInt& a, const BigInt& b)
-   { return (a.cmp(b) > 0); }
+   { return b.is_less_than(a); }
 
 inline bool operator==(const BigInt& a, word b)
    { return (a.cmp_word(b) == 0); }
diff --git a/src/lib/math/mp/mp_core.h b/src/lib/math/mp/mp_core.h
index 40ff5fadd..9a19a46be 100644
--- a/src/lib/math/mp/mp_core.h
+++ b/src/lib/math/mp/mp_core.h
@@ -15,6 +15,7 @@
 #include <botan/mem_ops.h>
 #include <botan/internal/mp_asmi.h>
 #include <botan/internal/ct_utils.h>
+#include <algorithm>
 
 namespace Botan {
 
@@ -170,6 +171,46 @@ inline void bigint_cnd_addsub(word mask, word x[], const word y[], size_t size)
    }
 
 /*
+* Equivalent to
+*   bigint_cnd_add( mask, x, size, y, size);
+*   bigint_cnd_sub(~mask, x, size, z, size);
+*
+* Mask must be either 0 or all 1 bits
+*
+* Returns the carry or borrow resp
+*/
+inline word bigint_cnd_addsub(word mask, word x[],
+                              const word y[], const word z[],
+                              size_t size)
+   {
+   const size_t blocks = size - (size % 8);
+
+   word carry = 0;
+   word borrow = 0;
+
+   word t0[8] = { 0 };
+   word t1[8] = { 0 };
+
+   for(size_t i = 0; i != blocks; i += 8)
+      {
+      carry = word8_add3(t0, x + i, y + i, carry);
+      borrow = word8_sub3(t1, x + i, z + i, borrow);
+
+      for(size_t j = 0; j != 8; ++j)
+         x[i+j] = CT::select(mask, t0[j], t1[j]);
+      }
+
+   for(size_t i = blocks; i != size; ++i)
+      {
+      t0[0] = word_add(x[i], y[i], &carry);
+      t1[0] = word_sub(x[i], z[i], &borrow);
+      x[i] = CT::select(mask, t0[0], t1[0]);
+      }
+
+   return CT::select(mask, carry, borrow);
+   }
+
+/*
 * 2s complement absolute value
 * If cond > 0 sets x to ~x + 1
 * Runs in constant time
@@ -331,7 +372,7 @@ inline word bigint_sub3(word z[],
 * Otherwise compute z = y - x
 * No borrow is possible since the result is always >= 0
 *
-* Returns 1 if x >= y or 0 if x < y
+* Returns ~0 if x >= y or 0 if x < y
 * @param z output array of at least N words
 * @param x input array of N words
 * @param y input array of N words
@@ -364,9 +405,7 @@ inline word bigint_sub_abs(word z[],
       ws1[i] = word_sub(y[i], x[i], &borrow1);
       }
 
-   word mask = CT::conditional_copy_mem(borrow1, z, ws0, ws1, N);
-
-   return CT::select<word>(mask, 0, 1);
+   return CT::conditional_copy_mem(borrow0, z, ws1, ws0, N);
    }
 
 /*
@@ -495,23 +534,6 @@ inline void bigint_linmul3(word z[], const word x[], size_t x_size, word y)
    }
 
 /**
-* Montgomery Reduction
-* @param z integer to reduce, of size exactly 2*(p_size+1).
-           Output is in the first p_size+1 words, higher
-           words are set to zero.
-* @param p modulus
-* @param p_size size of p
-* @param p_dash Montgomery value
-* @param workspace array of at least 2*(p_size+1) words
-* @param ws_size size of workspace in words
-*/
-void bigint_monty_redc(word z[],
-                       const word p[], size_t p_size,
-                       word p_dash,
-                       word workspace[],
-                       size_t ws_size);
-
-/**
 * Compare x and y
 * Return -1 if x < y
 * Return 0 if x == y
@@ -520,24 +542,116 @@ void bigint_monty_redc(word z[],
 inline int32_t bigint_cmp(const word x[], size_t x_size,
                           const word y[], size_t y_size)
    {
-   if(x_size < y_size) { return (-bigint_cmp(y, y_size, x, x_size)); }
+   static_assert(sizeof(word) >= sizeof(uint32_t), "Size assumption");
+
+   const uint32_t LT = static_cast<uint32_t>(-1);
+   const uint32_t EQ = 0;
+   const uint32_t GT = 1;
+
+   const size_t common_elems = std::min(x_size, y_size);
+
+   uint32_t result = EQ; // until found otherwise
+
+   for(size_t i = 0; i != common_elems; i++)
+      {
+      const word is_eq = CT::is_equal(x[i], y[i]);
+      const word is_lt = CT::is_less(x[i], y[i]);
+      result = CT::select<uint32_t>(is_eq, result, CT::select<uint32_t>(is_lt, LT, GT));
+      }
+
+   if(x_size < y_size)
+      {
+      word mask = 0;
+      for(size_t i = x_size; i != y_size; i++)
+         mask |= y[i];
+
+      // If any bits were set in high part of y, then x < y
+      result = CT::select<uint32_t>(CT::is_zero(mask), result, LT);
+      }
+   else if(y_size < x_size)
+      {
+      word mask = 0;
+      for(size_t i = y_size; i != x_size; i++)
+         mask |= x[i];
+
+      // If any bits were set in high part of x, then x > y
+      result = CT::select<uint32_t>(CT::is_zero(mask), result, GT);
+      }
+
+   CT::unpoison(result);
+   BOTAN_DEBUG_ASSERT(result == LT || result == GT || result == EQ);
+   return static_cast<int32_t>(result);
+   }
+
+/**
+* Compare x and y
+* Return ~0 if x[0:x_size] < y[0:y_size] or 0 otherwise
+* If lt_or_equal is true, returns ~0 also for x == y
+*/
+inline word bigint_ct_is_lt(const word x[], size_t x_size,
+                            const word y[], size_t y_size,
+                            bool lt_or_equal = false)
+   {
+   const size_t common_elems = std::min(x_size, y_size);
+
+   word is_lt = CT::expand_mask<word>(lt_or_equal);
 
-   while(x_size > y_size)
+   for(size_t i = 0; i != common_elems; i++)
       {
-      if(x[x_size-1])
-         return 1;
-      x_size--;
+      const word eq = CT::is_equal(x[i], y[i]);
+      const word lt = CT::is_less(x[i], y[i]);
+      is_lt = CT::select(eq, is_lt, lt);
       }
 
-   for(size_t i = x_size; i > 0; --i)
+   if(x_size < y_size)
       {
-      if(x[i-1] > y[i-1])
-         return 1;
-      if(x[i-1] < y[i-1])
-         return -1;
+      word mask = 0;
+      for(size_t i = x_size; i != y_size; i++)
+         mask |= y[i];
+      // If any bits were set in high part of y, then is_lt should be forced true
+      is_lt |= ~CT::is_zero(mask);
       }
+   else if(y_size < x_size)
+      {
+      word mask = 0;
+      for(size_t i = y_size; i != x_size; i++)
+         mask |= x[i];
 
-   return 0;
+      // If any bits were set in high part of x, then is_lt should be false
+      is_lt &= CT::is_zero(mask);
+      }
+
+   CT::unpoison(is_lt);
+   return is_lt;
+   }
+
+inline word bigint_ct_is_eq(const word x[], size_t x_size,
+                            const word y[], size_t y_size)
+   {
+   const size_t common_elems = std::min(x_size, y_size);
+
+   word diff = 0;
+
+   for(size_t i = 0; i != common_elems; i++)
+      {
+      diff |= (x[i] ^ y[i]);
+      }
+
+   // If any bits were set in high part of x/y, then they are not equal
+   if(x_size < y_size)
+      {
+      for(size_t i = x_size; i != y_size; i++)
+         diff |= y[i];
+      }
+   else if(y_size < x_size)
+      {
+      for(size_t i = y_size; i != x_size; i++)
+         diff |= x[i];
+      }
+
+   const word is_equal = CT::is_zero(diff);
+   CT::unpoison(is_equal);
+   return is_equal;
    }
 
 /**
@@ -578,6 +692,9 @@ inline word bigint_divop(word n1, word n0, word d)
 */
 inline word bigint_modop(word n1, word n0, word d)
    {
+   if(d == 0)
+      throw Invalid_Argument("bigint_modop divide by zero");
+
 #if defined(BOTAN_HAS_MP_DWORD)
    return ((static_cast<dword>(n1) << BOTAN_MP_WORD_BITS) | n0) % d;
 #else
@@ -605,6 +722,23 @@ void bigint_comba_sqr9(word out[18], const word in[9]);
 void bigint_comba_sqr16(word out[32], const word in[16]);
 void bigint_comba_sqr24(word out[48], const word in[24]);
 
+/**
+* Montgomery Reduction
+* @param z integer to reduce, of size exactly 2*(p_size+1).
+           Output is in the first p_size+1 words, higher
+           words are set to zero.
+* @param p modulus
+* @param p_size size of p
+* @param p_dash Montgomery value
+* @param workspace array of at least 2*(p_size+1) words
+* @param ws_size size of workspace in words
+*/
+void bigint_monty_redc(word z[],
+                       const word p[], size_t p_size,
+                       word p_dash,
+                       word workspace[],
+                       size_t ws_size);
+
 /*
 * High Level Multiplication/Squaring Interfaces
 */
diff --git a/src/lib/math/mp/mp_karat.cpp b/src/lib/math/mp/mp_karat.cpp
index 69e2ae665..c0fc5304b 100644
--- a/src/lib/math/mp/mp_karat.cpp
+++ b/src/lib/math/mp/mp_karat.cpp
@@ -83,18 +83,21 @@ void karatsuba_mul(word z[], const word x[], const word y[], size_t N,
    {
    if(N < KARATSUBA_MULTIPLY_THRESHOLD || N % 2)
       {
-      if(N == 6)
-         return bigint_comba_mul6(z, x, y);
-      else if(N == 8)
-         return bigint_comba_mul8(z, x, y);
-      else if(N == 9)
-         return bigint_comba_mul9(z, x, y);
-      else if(N == 16)
-         return bigint_comba_mul16(z, x, y);
-      else if(N == 24)
-         return bigint_comba_mul24(z, x, y);
-      else
-         return basecase_mul(z, 2*N, x, N, y, N);
+      switch(N)
+         {
+         case 6:
+            return bigint_comba_mul6(z, x, y);
+         case 8:
+            return bigint_comba_mul8(z, x, y);
+         case 9:
+            return bigint_comba_mul9(z, x, y);
+         case 16:
+            return bigint_comba_mul16(z, x, y);
+         case 24:
+            return bigint_comba_mul24(z, x, y);
+         default:
+            return basecase_mul(z, 2*N, x, N, y, N);
+         }
       }
 
    const size_t N2 = N / 2;
@@ -123,6 +126,7 @@ void karatsuba_mul(word z[], const word x[], const word y[], size_t N,
    // First compute (X_lo - X_hi)*(Y_hi - Y_lo)
    const word cmp0 = bigint_sub_abs(z0, x0, x1, N2, workspace);
    const word cmp1 = bigint_sub_abs(z1, y1, y0, N2, workspace);
+   const word neg_mask = ~(cmp0 ^ cmp1);
 
    karatsuba_mul(ws0, z0, z1, N2, ws1);
 
@@ -140,8 +144,6 @@ void karatsuba_mul(word z[], const word x[], const word y[], size_t N,
 
    clear_mem(workspace + N, N2);
 
-   const word neg_mask = CT::is_equal<word>(cmp0, cmp1);
-
    bigint_cnd_addsub(neg_mask, z + N2, workspace, 2*N-N2);
    }
 
@@ -152,18 +154,21 @@ void karatsuba_sqr(word z[], const word x[], size_t N, word workspace[])
    {
    if(N < KARATSUBA_SQUARE_THRESHOLD || N % 2)
       {
-      if(N == 6)
-         return bigint_comba_sqr6(z, x);
-      else if(N == 8)
-         return bigint_comba_sqr8(z, x);
-      else if(N == 9)
-         return bigint_comba_sqr9(z, x);
-      else if(N == 16)
-         return bigint_comba_sqr16(z, x);
-      else if(N == 24)
-         return bigint_comba_sqr24(z, x);
-      else
-         return basecase_sqr(z, 2*N, x, N);
+      switch(N)
+         {
+         case 6:
+            return bigint_comba_sqr6(z, x);
+         case 8:
+            return bigint_comba_sqr8(z, x);
+         case 9:
+            return bigint_comba_sqr9(z, x);
+         case 16:
+            return bigint_comba_sqr16(z, x);
+         case 24:
+            return bigint_comba_sqr24(z, x);
+         default:
+            return basecase_sqr(z, 2*N, x, N);
+         }
       }
 
    const size_t N2 = N / 2;
diff --git a/src/lib/math/numbertheory/nistp_redc.cpp b/src/lib/math/numbertheory/nistp_redc.cpp
index a005e1abb..a74c4de9f 100644
--- a/src/lib/math/numbertheory/nistp_redc.cpp
+++ b/src/lib/math/numbertheory/nistp_redc.cpp
@@ -25,18 +25,14 @@ void redc_p521(BigInt& x, secure_vector<word>& ws)
    const size_t p_top_bits = 521 % BOTAN_MP_WORD_BITS;
    const size_t p_words = p_full_words + 1;
 
-   const size_t x_sw = x.sig_words();
-
-   if(x_sw < p_words)
-      return; // already smaller
-
    if(ws.size() < p_words + 1)
       ws.resize(p_words + 1);
 
    clear_mem(ws.data(), ws.size());
-   bigint_shr2(ws.data(), x.data(), x_sw, p_full_words, p_top_bits);
+   bigint_shr2(ws.data(), x.data(), std::min(x.size(), 2*p_words), p_full_words, p_top_bits);
 
    x.mask_bits(521);
+   x.grow_to(p_words);
 
    // Word-level carry will be zero
    word carry = bigint_add3_nc(x.mutable_data(), x.data(), p_words, ws.data(), p_words);
@@ -45,13 +41,6 @@ void redc_p521(BigInt& x, secure_vector<word>& ws)
    // Now find the actual carry in bit 522
    const word bit_522_set = x.word_at(p_full_words) >> p_top_bits;
 
-#if (BOTAN_MP_WORD_BITS == 64)
-   static const word p521_words[9] = {
-      0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF,
-      0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF,
-      0x1FF };
-#endif
-
    /*
    * If bit 522 is set then we overflowed and must reduce. Otherwise, if the
    * top bit is set, it is possible we have x == 2**521 - 1 so check for that.
@@ -59,7 +48,12 @@ void redc_p521(BigInt& x, secure_vector<word>& ws)
    if(bit_522_set)
       {
 #if (BOTAN_MP_WORD_BITS == 64)
-      bigint_sub2(x.mutable_data(), x.size(), p521_words, 9);
+      static const word p521_words[9] = {
+         0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF,
+         0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF, 0xFFFFFFFFFFFFFFFF,
+         0x1FF };
+
+      bigint_sub2(x.mutable_data(), p_words, p521_words, 9);
 #else
       x -= prime_p521();
 #endif
@@ -141,6 +135,7 @@ void redc_p192(BigInt& x, secure_vector<word>& ws)
    const uint64_t S5 = X05 + X09 + X11;
 
    x.mask_bits(192);
+   x.resize(p192_limbs + 1);
 
    uint64_t S = 0;
    uint32_t R0 = 0, R1 = 0;
@@ -194,9 +189,10 @@ void redc_p192(BigInt& x, secure_vector<word>& ws)
 #endif
    };
 
-   word borrow = bigint_sub2(x.mutable_data(), x.size(), p192_mults[S], p192_limbs);
+   BOTAN_ASSERT_NOMSG(x.size() == p192_limbs + 1);
+   word borrow = bigint_sub2(x.mutable_data(), p192_limbs + 1, p192_mults[S], p192_limbs);
    BOTAN_DEBUG_ASSERT(borrow == 0 || borrow == 1);
-   bigint_cnd_add(borrow, x.mutable_data(), x.size(), p192_mults[0], p192_limbs);
+   bigint_cnd_add(borrow, x.mutable_data(), p192_limbs + 1, p192_mults[0], p192_limbs);
    }
 
 const BigInt& prime_p224()
@@ -293,9 +289,10 @@ void redc_p224(BigInt& x, secure_vector<word>& ws)
 
    };
 
-   word borrow = bigint_sub2(x.mutable_data(), x.size(), p224_mults[S], p224_limbs);
+   BOTAN_ASSERT_NOMSG(x.size() == p224_limbs + 1);
+   word borrow = bigint_sub2(x.mutable_data(), p224_limbs + 1, p224_mults[S], p224_limbs);
    BOTAN_DEBUG_ASSERT(borrow == 0 || borrow == 1);
-   bigint_cnd_add(borrow, x.mutable_data(), x.size(), p224_mults[0], p224_limbs);
+   bigint_cnd_add(borrow, x.mutable_data(), p224_limbs + 1, p224_mults[0], p224_limbs);
    }
 
 const BigInt& prime_p256()
@@ -420,9 +417,10 @@ void redc_p256(BigInt& x, secure_vector<word>& ws)
 
    CT::unpoison(S);
 
-   word borrow = bigint_sub2(x.mutable_data(), x.size(), p256_mults[S], p256_limbs);
+   BOTAN_ASSERT_NOMSG(x.size() == p256_limbs + 1);
+   word borrow = bigint_sub2(x.mutable_data(), p256_limbs + 1, p256_mults[S], p256_limbs);
    BOTAN_DEBUG_ASSERT(borrow == 0 || borrow == 1);
-   bigint_cnd_add(borrow, x.mutable_data(), x.size(), p256_mults[0], p256_limbs);
+   bigint_cnd_add(borrow, x.mutable_data(), p256_limbs + 1, p256_mults[0], p256_limbs);
    }
 
 const BigInt& prime_p384()
@@ -570,9 +568,10 @@ void redc_p384(BigInt& x, secure_vector<word>& ws)
 #endif
    };
 
-   word borrow = bigint_sub2(x.mutable_data(), x.size(), p384_mults[S], p384_limbs);
+   BOTAN_ASSERT_NOMSG(x.size() == p384_limbs + 1);
+   word borrow = bigint_sub2(x.mutable_data(), p384_limbs + 1, p384_mults[S], p384_limbs);
    BOTAN_DEBUG_ASSERT(borrow == 0 || borrow == 1);
-   bigint_cnd_add(borrow, x.mutable_data(), x.size(), p384_mults[0], p384_limbs);
+   bigint_cnd_add(borrow, x.mutable_data(), p384_limbs + 1, p384_mults[0], p384_limbs);
    }
 
 #endif
diff --git a/src/lib/math/numbertheory/numthry.cpp b/src/lib/math/numbertheory/numthry.cpp
index a8dfe8eaf..399a49cea 100644
--- a/src/lib/math/numbertheory/numthry.cpp
+++ b/src/lib/math/numbertheory/numthry.cpp
@@ -292,9 +292,8 @@ BigInt inverse_mod(const BigInt& n, const BigInt& mod)
       throw BigInt::DivideByZero();
    if(mod.is_negative() || n.is_negative())
       throw Invalid_Argument("inverse_mod: arguments must be non-negative");
-
-   if(n.is_zero() || (n.is_even() && mod.is_even()))
-      return 0; // fast fail checks
+   if(n.is_zero())
+      return 0;
 
    if(mod.is_odd() && n < mod)
       return ct_inverse_mod_odd_modulus(n, mod);
diff --git a/src/lib/pubkey/ec_group/point_gfp.cpp b/src/lib/pubkey/ec_group/point_gfp.cpp
index 77803de78..7bc6c4975 100644
--- a/src/lib/pubkey/ec_group/point_gfp.cpp
+++ b/src/lib/pubkey/ec_group/point_gfp.cpp
@@ -341,6 +341,7 @@ void PointGFp::mult2(std::vector<BigInt>& ws_bn)
 
       m_curve.sqr(T4, m_coord_x, ws); // x^2
       T4 *= 3; // 3*x^2
+      T4.reduce_below(p, sub_ws);
       T4.mod_add(T3, p, sub_ws); // 3*x^2 + a*z^4
       }