Add power analysis countermeasures for ECC point multiplications.

The plain PointGFp operator* now uses Montgomery ladder exclusively.

Adds a blinded point multiply algorithm which uses exponent and point
randomization, as well as a Montgomery ladder technique that takes a
random walk of the possible addition chains for k.

10 files changed, 411 insertions, 242 deletions

diff --git a/src/lib/math/ec_gfp/point_gfp.cpp b/src/lib/math/ec_gfp/point_gfp.cpp
index 2505e4d54..a319d8657 100644
--- a/src/lib/math/ec_gfp/point_gfp.cpp
+++ b/src/lib/math/ec_gfp/point_gfp.cpp
@@ -2,36 +2,55 @@
 * Point arithmetic on elliptic curves over GF(p)
 *
 * (C) 2007 Martin Doering, Christoph Ludwig, Falko Strenzke
-*     2008-2011,2012,2014 Jack Lloyd
+*     2008-2011,2012,2014,2015 Jack Lloyd
 *
 * Botan is released under the Simplified BSD License (see license.txt)
 */
 
 #include <botan/point_gfp.h>
 #include <botan/numthry.h>
+#include <botan/loadstor.h>
 
 namespace Botan {
 
+const size_t BOTAN_POINTGFP_MONTGOMERY_BLINDING_BITS = 20;
+const size_t BOTAN_POINTGFP_RANDOMIZE_BLINDING_BITS = 48;
+
+#define BOTAN_POINTGFP_BLINDED_MULTIPLY_USE_MONTGOMERY_LADDER 1
+
 PointGFp::PointGFp(const CurveGFp& curve) :
-   curve(curve),
-   coord_x(0),
-   coord_y(1),
-   coord_z(0)
+   m_curve(curve),
+   m_coord_x(0),
+   m_coord_y(1),
+   m_coord_z(0)
    {
-   curve.to_rep(coord_x, ws);
-   curve.to_rep(coord_y, ws);
-   curve.to_rep(coord_z, ws);
+   m_curve.to_rep(m_coord_x, m_monty_ws);
+   m_curve.to_rep(m_coord_y, m_monty_ws);
+   m_curve.to_rep(m_coord_z, m_monty_ws);
    }
 
 PointGFp::PointGFp(const CurveGFp& curve, const BigInt& x, const BigInt& y) :
-   curve(curve),
-   coord_x(x),
-   coord_y(y),
-   coord_z(1)
+   m_curve(curve),
+   m_coord_x(x),
+   m_coord_y(y),
+   m_coord_z(1)
+   {
+   m_curve.to_rep(m_coord_x, m_monty_ws);
+   m_curve.to_rep(m_coord_y, m_monty_ws);
+   m_curve.to_rep(m_coord_z, m_monty_ws);
+   }
+
+void PointGFp::randomize_repr(RandomNumberGenerator& rng)
    {
-   curve.to_rep(coord_x, ws);
-   curve.to_rep(coord_y, ws);
-   curve.to_rep(coord_z, ws);
+   BigInt mask(rng, BOTAN_POINTGFP_RANDOMIZE_BLINDING_BITS, false);
+
+   m_curve.to_rep(mask, m_monty_ws);
+   const BigInt mask2 = curve_mult(mask, mask);
+   const BigInt mask3 = curve_mult(mask2, mask);
+
+   m_coord_x = curve_mult(m_coord_x, mask2);
+   m_coord_y = curve_mult(m_coord_y, mask3);
+   m_coord_z = curve_mult(m_coord_z, mask);
    }
 
 // Point addition
@@ -39,15 +58,15 @@ void PointGFp::add(const PointGFp& rhs, std::vector<BigInt>& ws_bn)
    {
    if(is_zero())
       {
-      coord_x = rhs.coord_x;
-      coord_y = rhs.coord_y;
-      coord_z = rhs.coord_z;
+      m_coord_x = rhs.m_coord_x;
+      m_coord_y = rhs.m_coord_y;
+      m_coord_z = rhs.m_coord_z;
       return;
       }
    else if(rhs.is_zero())
       return;
 
-   const BigInt& p = curve.get_p();
+   const BigInt& p = m_curve.get_p();
 
    BigInt& rhs_z2 = ws_bn[0];
    BigInt& U1 = ws_bn[1];
@@ -64,13 +83,13 @@ void PointGFp::add(const PointGFp& rhs, std::vector<BigInt>& ws_bn)
    http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#addition-add-1998-cmo-2
    */
 
-   curve_sqr(rhs_z2, rhs.coord_z);
-   curve_mult(U1, coord_x, rhs_z2);
-   curve_mult(S1, coord_y, curve_mult(rhs.coord_z, rhs_z2));
+   curve_sqr(rhs_z2, rhs.m_coord_z);
+   curve_mult(U1, m_coord_x, rhs_z2);
+   curve_mult(S1, m_coord_y, curve_mult(rhs.m_coord_z, rhs_z2));
 
-   curve_sqr(lhs_z2, coord_z);
-   curve_mult(U2, rhs.coord_x, lhs_z2);
-   curve_mult(S2, rhs.coord_y, curve_mult(coord_z, lhs_z2));
+   curve_sqr(lhs_z2, m_coord_z);
+   curve_mult(U2, rhs.m_coord_x, lhs_z2);
+   curve_mult(S2, rhs.m_coord_y, curve_mult(m_coord_z, lhs_z2));
 
    H = U2;
    H -= U1;
@@ -90,7 +109,10 @@ void PointGFp::add(const PointGFp& rhs, std::vector<BigInt>& ws_bn)
          return;
          }
 
-      *this = PointGFp(curve); // setting myself to zero
+      // setting to zero:
+      m_coord_x = 0;
+      m_coord_y = 1;
+      m_coord_z = 0;
       return;
       }
 
@@ -100,22 +122,22 @@ void PointGFp::add(const PointGFp& rhs, std::vector<BigInt>& ws_bn)
 
    U2 = curve_mult(U1, U2);
 
-   curve_sqr(coord_x, r);
-   coord_x -= S2;
-   coord_x -= (U2 << 1);
-   while(coord_x.is_negative())
-      coord_x += p;
+   curve_sqr(m_coord_x, r);
+   m_coord_x -= S2;
+   m_coord_x -= (U2 << 1);
+   while(m_coord_x.is_negative())
+      m_coord_x += p;
 
-   U2 -= coord_x;
+   U2 -= m_coord_x;
    if(U2.is_negative())
       U2 += p;
 
-   curve_mult(coord_y, r, U2);
-   coord_y -= curve_mult(S1, S2);
-   if(coord_y.is_negative())
-      coord_y += p;
+   curve_mult(m_coord_y, r, U2);
+   m_coord_y -= curve_mult(S1, S2);
+   if(m_coord_y.is_negative())
+      m_coord_y += p;
 
-   curve_mult(coord_z, curve_mult(coord_z, rhs.coord_z), H);
+   curve_mult(m_coord_z, curve_mult(m_coord_z, rhs.m_coord_z), H);
    }
 
 // *this *= 2
@@ -123,9 +145,9 @@ void PointGFp::mult2(std::vector<BigInt>& ws_bn)
    {
    if(is_zero())
       return;
-   else if(coord_y.is_zero())
+   else if(m_coord_y.is_zero())
       {
-      *this = PointGFp(curve); // setting myself to zero
+      *this = PointGFp(m_curve); // setting myself to zero
       return;
       }
 
@@ -133,7 +155,7 @@ void PointGFp::mult2(std::vector<BigInt>& ws_bn)
    http://hyperelliptic.org/EFD/g1p/auto-shortw-jacobian-3.html#doubling-dbl-1986-cc
    */
 
-   const BigInt& p = curve.get_p();
+   const BigInt& p = m_curve.get_p();
 
    BigInt& y_2 = ws_bn[0];
    BigInt& S = ws_bn[1];
@@ -145,17 +167,17 @@ void PointGFp::mult2(std::vector<BigInt>& ws_bn)
    BigInt& y = ws_bn[7];
    BigInt& z = ws_bn[8];
 
-   curve_sqr(y_2, coord_y);
+   curve_sqr(y_2, m_coord_y);
 
-   curve_mult(S, coord_x, y_2);
+   curve_mult(S, m_coord_x, y_2);
    S <<= 2; // * 4
    while(S >= p)
       S -= p;
 
-   curve_sqr(z4, curve_sqr(coord_z));
-   curve_mult(a_z4, curve.get_a_rep(), z4);
+   curve_sqr(z4, curve_sqr(m_coord_z));
+   curve_mult(a_z4, m_curve.get_a_rep(), z4);
 
-   M = curve_sqr(coord_x);
+   M = curve_sqr(m_coord_x);
    M *= 3;
    M += a_z4;
    while(M >= p)
@@ -180,14 +202,14 @@ void PointGFp::mult2(std::vector<BigInt>& ws_bn)
    if(y.is_negative())
       y += p;
 
-   curve_mult(z, coord_y, coord_z);
+   curve_mult(z, m_coord_y, m_coord_z);
    z <<= 1;
    if(z >= p)
       z -= p;
 
-   coord_x = x;
-   coord_y = y;
-   coord_z = z;
+   m_coord_x = x;
+   m_coord_y = y;
+   m_coord_z = z;
    }
 
 // arithmetic operators
@@ -221,7 +243,7 @@ PointGFp multi_exponentiate(const PointGFp& p1, const BigInt& z1,
    {
    const PointGFp p3 = p1 + p2;
 
-   PointGFp H(p1.curve); // create as zero
+   PointGFp H(p1.get_curve()); // create as zero
    size_t bits_left = std::max(z1.bits(), z2.bits());
 
    std::vector<BigInt> ws(9);
@@ -251,22 +273,24 @@ PointGFp multi_exponentiate(const PointGFp& p1, const BigInt& z1,
 
 PointGFp operator*(const BigInt& scalar, const PointGFp& point)
    {
-   //BOTAN_ASSERT(point.on_the_curve(), "Input is valid");
+   //BOTAN_ASSERT(point.on_the_curve(), "Input is on the curve");
 
    const CurveGFp& curve = point.get_curve();
 
-   if(scalar.is_zero())
-      return PointGFp(curve); // zero point
+   const size_t scalar_bits = scalar.bits();
 
    std::vector<BigInt> ws(9);
 
-   if(scalar.abs() <= 2) // special cases for small values
+   if(scalar_bits <= 2)
       {
-      byte value = scalar.abs().byte_at(0);
+      const byte abs_val = scalar.byte_at(0);
+
+      if(abs_val == 0)
+         return PointGFp::zero_of(curve);
 
       PointGFp result = point;
 
-      if(value == 2)
+      if(abs_val == 2)
          result.mult2(ws);
 
       if(scalar.is_negative())
@@ -275,94 +299,164 @@ PointGFp operator*(const BigInt& scalar, const PointGFp& point)
       return result;
       }
 
-   const size_t scalar_bits = scalar.bits();
+   PointGFp R[2] = { PointGFp(curve), point };
 
-   PointGFp x1(curve); // zero
+   for(size_t i = scalar_bits; i > 0; i--)
+      {
+      const size_t b = scalar.get_bit(i - 1);
+      R[b ^ 1].add(R[b], ws);
+      R[b].mult2(ws);
+      }
 
-   size_t bits_left = scalar_bits;
+   if(scalar.is_negative())
+      R[0].negate();
 
-#if BOTAN_CURVE_GFP_USE_MONTGOMERY_LADDER
+   //BOTAN_ASSERT(R[0].on_the_curve(), "Output is on the curve");
 
-   PointGFp x2 = point;
-   while(bits_left)
+   return R[0];
+   }
+
+Blinded_Point_Multiply::Blinded_Point_Multiply(const PointGFp& base, const BigInt& order, size_t h) :
+   m_order(order), m_h(h ? h : 4), m_ws(9)
+   {
+   // Upper bound is a sanity check rather than hard limit
+   if(m_h < 1 || m_h > 8)
+      throw std::invalid_argument("Blinded_Point_Multiply invalid h param");
+
+   const CurveGFp& curve = base.get_curve();
+
+#if BOTAN_POINTGFP_BLINDED_MULTIPLY_USE_MONTGOMERY_LADDER
+
+   const PointGFp inv = -base;
+
+   m_U.resize(6*m_h + 3);
+
+   m_U[3*m_h+0] = inv;
+   m_U[3*m_h+1] = PointGFp::zero_of(curve);
+   m_U[3*m_h+2] = base;
+
+   for(size_t i = 1; i <= 3 * m_h + 1; ++i)
       {
-      if(scalar.get_bit(bits_left - 1))
-         {
-         x1.add(x2, ws);
-         x2.mult2(ws);
-         }
-      else
-         {
-         x2.add(x1, ws);
-         x1.mult2(ws);
-         }
+      m_U[3*m_h+1+i] = m_U[3*m_h+i];
+      m_U[3*m_h+1+i].add(base, m_ws);
 
-      --bits_left;
+      m_U[3*m_h+1-i] = m_U[3*m_h+2-i];
+      m_U[3*m_h+1-i].add(inv, m_ws);
       }
-
 #else
-   const size_t window_bits = 4;
+   m_U.resize(1 << m_h);
+   m_U[0] = PointGFp::zero_of(curve);
+   m_U[1] = base;
+
+   for(size_t i = 2; i < m_U.size(); ++i)
+      {
+      m_U[i] = m_U[i-1];
+      m_U[i].add(base, m_ws);
+      }
+#endif
+   }
 
-   std::vector<PointGFp> Ps(1 << window_bits);
-   Ps[0] = x1;
-   Ps[1] = point;
+PointGFp Blinded_Point_Multiply::blinded_multiply(const BigInt& scalar_in,
+                                                  RandomNumberGenerator& rng)
+   {
+   if(scalar_in.is_negative())
+      throw std::invalid_argument("Blinded_Point_Multiply scalar must be positive");
+
+   // Choose a small mask m and use k' = k + m*order (Coron's 1st countermeasure)
+   const u64bit mask = rng.gen_mask(BOTAN_POINTGFP_MONTGOMERY_BLINDING_BITS);
+   const BigInt scalar = scalar_in + m_order * mask;
+   const size_t scalar_bits = scalar.bits();
+
+   // Randomize each point representation (Coron's 3rd countermeasure)
+   for(size_t i = 0; i != m_U.size(); ++i)
+      m_U[i].randomize_repr(rng);
+
+#if BOTAN_POINTGFP_BLINDED_MULTIPLY_USE_MONTGOMERY_LADDER
+   PointGFp R = m_U.at(3*m_h + 2); // base point
+   int32_t alpha = 0;
+
+   R.randomize_repr(rng);
+
+   /*
+   Algorithm 7 from "Randomizing the Montgomery Powering Ladder"
+   Duc-Phong Le, Chik How Tan and Michael Tunstall
+   http://eprint.iacr.org/2015/657
 
-   for(size_t i = 2; i < Ps.size(); ++i)
+   It takes a random walk through (a subset of) the set of addition
+   chains that end in k.
+   */
+   for(size_t i = scalar_bits; i > 0; i--)
       {
-      Ps[i] = Ps[i-1];
-      Ps[i].add(point, ws);
+      const int32_t ki = scalar.get_bit(i);
+
+      // choose gamma from -h,...,h
+      const int32_t gamma = static_cast<int32_t>((rng.next_byte() % (2*m_h))) - m_h;
+      const int32_t l = gamma - 2*alpha + ki - (ki ^ 1);
+
+      R.mult2(m_ws);
+      R.add(m_U.at(3*m_h + 1 + l), m_ws);
+      alpha = gamma;
       }
 
-   while(bits_left >= window_bits)
+   const int32_t k0 = scalar.get_bit(0);
+   R.add(m_U[3*m_h + 1 - alpha - (k0 ^ 1)], m_ws);
+
+#else
+
+   size_t bits_left = scalar_bits;
+
+   PointGFp R = m_U[0];
+
+   while(bits_left >= m_h)
       {
-      for(size_t i = 0; i != window_bits; ++i)
-         x1.mult2(ws);
+      for(size_t i = 0; i != m_h; ++i)
+         R.mult2(m_ws);
 
-      const u32bit nibble = scalar.get_substring(bits_left - window_bits, window_bits);
-      x1.add(Ps[nibble], ws);
-      bits_left -= window_bits;
+      const u32bit nibble = scalar.get_substring(bits_left - m_h, m_h);
+      R.add(m_U[nibble], m_ws);
+      bits_left -= m_h;
       }
 
    while(bits_left)
       {
-      x1.mult2(ws);
+      R.mult2(m_ws);
       if(scalar.get_bit(bits_left-1))
-         x1.add(point, ws);
+         R.add(m_U[1], m_ws);
       --bits_left;
       }
-
 #endif
 
-   if(scalar.is_negative())
-      x1.negate();
+   //BOTAN_ASSERT(R.on_the_curve(), "Output is on the curve");
 
-   //BOTAN_ASSERT(x1.on_the_curve(), "Output is on the curve");
-
-   return x1;
+   return R;
    }
 
 BigInt PointGFp::get_affine_x() const
    {
    if(is_zero())
+      abort();
+   if(is_zero())
       throw Illegal_Transformation("Cannot convert zero point to affine");
 
-   BigInt z2 = curve_sqr(coord_z);
-   curve.from_rep(z2, ws);
-   z2 = inverse_mod(z2, curve.get_p());
+   BigInt z2 = curve_sqr(m_coord_z);
+   m_curve.from_rep(z2, m_monty_ws);
+   z2 = inverse_mod(z2, m_curve.get_p());
 
-   return curve_mult(z2, coord_x);
+   return curve_mult(z2, m_coord_x);
    }
 
 BigInt PointGFp::get_affine_y() const
    {
    if(is_zero())
+      abort();
+   if(is_zero())
       throw Illegal_Transformation("Cannot convert zero point to affine");
 
-   BigInt z3 = curve_mult(coord_z, curve_sqr(coord_z));
-   z3 = inverse_mod(z3, curve.get_p());
-   curve.to_rep(z3, ws);
+   BigInt z3 = curve_mult(m_coord_z, curve_sqr(m_coord_z));
+   z3 = inverse_mod(z3, m_curve.get_p());
+   m_curve.to_rep(z3, m_monty_ws);
 
-   return curve_mult(z3, coord_y);
+   return curve_mult(z3, m_coord_y);
    }
 
 bool PointGFp::on_the_curve() const
@@ -376,22 +470,22 @@ bool PointGFp::on_the_curve() const
    if(is_zero())
       return true;
 
-   const BigInt y2 = curve.from_rep(curve_sqr(coord_y), ws);
-   const BigInt x3 = curve_mult(coord_x, curve_sqr(coord_x));
-   const BigInt ax = curve_mult(coord_x, curve.get_a_rep());
-   const BigInt z2 = curve_sqr(coord_z);
+   const BigInt y2 = m_curve.from_rep(curve_sqr(m_coord_y), m_monty_ws);
+   const BigInt x3 = curve_mult(m_coord_x, curve_sqr(m_coord_x));
+   const BigInt ax = curve_mult(m_coord_x, m_curve.get_a_rep());
+   const BigInt z2 = curve_sqr(m_coord_z);
 
-   if(coord_z == z2) // Is z equal to 1 (in Montgomery form)?
+   if(m_coord_z == z2) // Is z equal to 1 (in Montgomery form)?
       {
-      if(y2 != curve.from_rep(x3 + ax + curve.get_b_rep(), ws))
+      if(y2 != m_curve.from_rep(x3 + ax + m_curve.get_b_rep(), m_monty_ws))
          return false;
       }
 
-   const BigInt z3 = curve_mult(coord_z, z2);
+   const BigInt z3 = curve_mult(m_coord_z, z2);
    const BigInt ax_z4 = curve_mult(ax, curve_sqr(z2));
-   const BigInt b_z6 = curve_mult(curve.get_b_rep(), curve_sqr(z3));
+   const BigInt b_z6 = curve_mult(m_curve.get_b_rep(), curve_sqr(z3));
 
-   if(y2 != curve.from_rep(x3 + ax_z4 + b_z6, ws))
+   if(y2 != m_curve.from_rep(x3 + ax_z4 + b_z6, m_monty_ws))
       return false;
 
    return true;
@@ -400,11 +494,11 @@ bool PointGFp::on_the_curve() const
 // swaps the states of *this and other, does not throw!
 void PointGFp::swap(PointGFp& other)
    {
-   curve.swap(other.curve);
-   coord_x.swap(other.coord_x);
-   coord_y.swap(other.coord_y);
-   coord_z.swap(other.coord_z);
-   ws.swap(other.ws);
+   m_curve.swap(other.m_curve);
+   m_coord_x.swap(other.m_coord_x);
+   m_coord_y.swap(other.m_coord_y);
+   m_coord_z.swap(other.m_coord_z);
+   m_monty_ws.swap(other.m_monty_ws);
    }
 
 bool PointGFp::operator==(const PointGFp& other) const
diff --git a/src/lib/math/ec_gfp/point_gfp.h b/src/lib/math/ec_gfp/point_gfp.h
index 813ead81e..bb3438697 100644
--- a/src/lib/math/ec_gfp/point_gfp.h
+++ b/src/lib/math/ec_gfp/point_gfp.h
@@ -2,7 +2,7 @@
 * Point arithmetic on elliptic curves over GF(p)
 *
 * (C) 2007 Martin Doering, Christoph Ludwig, Falko Strenzke
-*     2008-2011,2014 Jack Lloyd
+*     2008-2011,2014,2015 Jack Lloyd
 *
 * Botan is released under the Simplified BSD License (see license.txt)
 */
@@ -58,6 +58,11 @@ class BOTAN_DLL PointGFp
       */
       PointGFp(const CurveGFp& curve);
 
+      static PointGFp zero_of(const CurveGFp& curve)
+         {
+         return PointGFp(curve);
+         }
+
       /**
       * Copy constructor
       */
@@ -113,6 +118,7 @@ class BOTAN_DLL PointGFp
       * @param scalar the PointGFp to multiply with *this
       * @result resulting PointGFp
       */
+
       PointGFp& operator*=(const BigInt& scalar);
 
       /**
@@ -142,7 +148,7 @@ class BOTAN_DLL PointGFp
       PointGFp& negate()
          {
          if(!is_zero())
-            coord_y = curve.get_p() - coord_y;
+            m_coord_y = m_curve.get_p() - m_coord_y;
          return *this;
          }
 
@@ -150,7 +156,7 @@ class BOTAN_DLL PointGFp
       * Return base curve of this point
       * @result the curve over GF(p) of this point
       */
-      const CurveGFp& get_curve() const { return curve; }
+      const CurveGFp& get_curve() const { return m_curve; }
 
       /**
       * get affine x coordinate
@@ -169,7 +175,7 @@ class BOTAN_DLL PointGFp
       * @result true, if this point is at infinity, false otherwise.
       */
       bool is_zero() const
-         { return (coord_x.is_zero() && coord_z.is_zero()); }
+         { return (m_coord_x.is_zero() && m_coord_z.is_zero()); }
 
       /**
       * Checks whether the point is to be found on the underlying
@@ -185,33 +191,40 @@ class BOTAN_DLL PointGFp
       void swap(PointGFp& other);
 
       /**
+      * Randomize the point representation
+      * The actual value (get_affine_x, get_affine_y) does not change
+      */
+      void randomize_repr(RandomNumberGenerator& rng);
+
+      /**
       * Equality operator
       */
       bool operator==(const PointGFp& other) const;
    private:
+      friend class Blinded_Point_Multiply;
 
       BigInt curve_mult(const BigInt& x, const BigInt& y) const
          {
          BigInt z;
-         curve.mul(z, x, y, ws);
+         m_curve.mul(z, x, y, m_monty_ws);
          return z;
          }
 
       void curve_mult(BigInt& z, const BigInt& x, const BigInt& y) const
          {
-         curve.mul(z, x, y, ws);
+         m_curve.mul(z, x, y, m_monty_ws);
          }
 
       BigInt curve_sqr(const BigInt& x) const
          {
          BigInt z;
-         curve.sqr(z, x, ws);
+         m_curve.sqr(z, x, m_monty_ws);
          return z;
          }
 
       void curve_sqr(BigInt& z, const BigInt& x) const
          {
-         curve.sqr(z, x, ws);
+         m_curve.sqr(z, x, m_monty_ws);
          }
 
       /**
@@ -226,9 +239,9 @@ class BOTAN_DLL PointGFp
       */
       void mult2(std::vector<BigInt>& workspace);
 
-      CurveGFp curve;
-      BigInt coord_x, coord_y, coord_z;
-      mutable secure_vector<word> ws; // workspace for Montgomery
+      CurveGFp m_curve;
+      BigInt m_coord_x, m_coord_y, m_coord_z;
+      mutable secure_vector<word> m_monty_ws; // workspace for Montgomery
    };
 
 // relational operators
@@ -270,6 +283,22 @@ template<typename Alloc>
 PointGFp OS2ECP(const std::vector<byte, Alloc>& data, const CurveGFp& curve)
    { return OS2ECP(data.data(), data.size(), curve); }
 
+/**
+
+*/
+class BOTAN_DLL Blinded_Point_Multiply
+   {
+   public:
+      Blinded_Point_Multiply(const PointGFp& base, const BigInt& order, size_t h = 0);
+
+      PointGFp blinded_multiply(const BigInt& scalar, RandomNumberGenerator& rng);
+   private:
+      const BigInt& m_order;
+      const size_t m_h;
+      std::vector<BigInt> m_ws;
+      std::vector<PointGFp> m_U;
+   };
+
 }
 
 namespace std {
diff --git a/src/lib/pubkey/ecdsa/ecdsa.cpp b/src/lib/pubkey/ecdsa/ecdsa.cpp
index 2518a14fe..c327a37c3 100644
--- a/src/lib/pubkey/ecdsa/ecdsa.cpp
+++ b/src/lib/pubkey/ecdsa/ecdsa.cpp
@@ -2,14 +2,13 @@
 * ECDSA implemenation
 * (C) 2007 Manuel Hartl, FlexSecure GmbH
 *     2007 Falko Strenzke, FlexSecure GmbH
-*     2008-2010 Jack Lloyd
+*     2008-2010,2015 Jack Lloyd
 *
 * Botan is released under the Simplified BSD License (see license.txt)
 */
 
 #include <botan/internal/pk_utils.h>
 #include <botan/ecdsa.h>
-#include <botan/reducer.h>
 #include <botan/keypair.h>
 #include <botan/rfc6979.h>
 
@@ -41,9 +40,10 @@ class ECDSA_Signature_Operation : public PK_Ops::Signature_with_EMSA
                                 const std::string& emsa) :
          PK_Ops::Signature_with_EMSA(emsa),
          base_point(ecdsa.domain().get_base_point()),
-         order(ecdsa.domain().get_order()),
-         x(ecdsa.private_value()),
-         mod_order(order),
+         m_order(ecdsa.domain().get_order()),
+         m_base_point(ecdsa.domain().get_base_point(), m_order),
+         m_x(ecdsa.private_value()),
+         m_mod_order(m_order),
          m_hash(hash_for_deterministic_signature(emsa))
          {
          }
@@ -52,34 +52,36 @@ class ECDSA_Signature_Operation : public PK_Ops::Signature_with_EMSA
                                    RandomNumberGenerator& rng) override;
 
       size_t message_parts() const override { return 2; }
-      size_t message_part_size() const override { return order.bytes(); }
-      size_t max_input_bits() const override { return order.bits(); }
+      size_t message_part_size() const override { return m_order.bytes(); }
+      size_t max_input_bits() const override { return m_order.bits(); }
 
    private:
       const PointGFp& base_point;
-      const BigInt& order;
-      const BigInt& x;
-      Modular_Reducer mod_order;
+      const BigInt& m_order;
+      Blinded_Point_Multiply m_base_point;
+      const BigInt& m_x;
+      Modular_Reducer m_mod_order;
       std::string m_hash;
    };
 
 secure_vector<byte>
 ECDSA_Signature_Operation::raw_sign(const byte msg[], size_t msg_len,
-                                    RandomNumberGenerator&)
+                                    RandomNumberGenerator& rng)
    {
    const BigInt m(msg, msg_len);
 
-   const BigInt k = generate_rfc6979_nonce(x, order, m, m_hash);
+   const BigInt k = generate_rfc6979_nonce(m_x, m_order, m, m_hash);
 
-   const PointGFp k_times_P = base_point * k;
-   const BigInt r = mod_order.reduce(k_times_P.get_affine_x());
-   const BigInt s = mod_order.multiply(inverse_mod(k, order), mul_add(x, r, m));
+   //const PointGFp k_times_P = base_point * k;
+   const PointGFp k_times_P = m_base_point.blinded_multiply(k, rng);
+   const BigInt r = m_mod_order.reduce(k_times_P.get_affine_x());
+   const BigInt s = m_mod_order.multiply(inverse_mod(k, m_order), mul_add(m_x, r, m));
 
    // With overwhelming probability, a bug rather than actual zero r/s
    BOTAN_ASSERT(s != 0, "invalid s");
    BOTAN_ASSERT(r != 0, "invalid r");
 
-   secure_vector<byte> output(2*order.bytes());
+   secure_vector<byte> output(2*m_order.bytes());
    r.binary_encode(&output[output.size() / 2 - r.bytes()]);
    s.binary_encode(&output[output.size() - s.bytes()]);
    return output;
diff --git a/src/lib/pubkey/gost_3410/gost_3410.cpp b/src/lib/pubkey/gost_3410/gost_3410.cpp
index 9c3a0ef3c..f04692d12 100644
--- a/src/lib/pubkey/gost_3410/gost_3410.cpp
+++ b/src/lib/pubkey/gost_3410/gost_3410.cpp
@@ -2,7 +2,7 @@
 * GOST 34.10-2001 implemenation
 * (C) 2007 Falko Strenzke, FlexSecure GmbH
 *          Manuel Hartl, FlexSecure GmbH
-* (C) 2008-2010 Jack Lloyd
+* (C) 2008-2010,2015 Jack Lloyd
 *
 * Botan is released under the Simplified BSD License (see license.txt)
 */
@@ -16,7 +16,6 @@ namespace Botan {
 
 std::vector<byte> GOST_3410_PublicKey::x509_subject_public_key() const
    {
-   // Trust CryptoPro to come up with something obnoxious
    const BigInt x = public_point().get_affine_x();
    const BigInt y = public_point().get_affine_y();
 
@@ -53,7 +52,7 @@ GOST_3410_PublicKey::GOST_3410_PublicKey(const AlgorithmIdentifier& alg_id,
    {
    OID ecc_param_id;
 
-   // Also includes hash and cipher OIDs... brilliant design guys
+   // The parameters also includes hash and cipher OIDs
    BER_Decoder(alg_id.parameters).start_cons(SEQUENCE).decode(ecc_param_id);
 
    domain_params = EC_Group(ecc_param_id);
@@ -101,21 +100,23 @@ class GOST_3410_Signature_Operation : public PK_Ops::Signature_with_EMSA
       GOST_3410_Signature_Operation(const GOST_3410_PrivateKey& gost_3410,
                                     const std::string& emsa) :
          PK_Ops::Signature_with_EMSA(emsa),
-         base_point(gost_3410.domain().get_base_point()),
-         order(gost_3410.domain().get_order()),
-         x(gost_3410.private_value()) {}
+         m_order(gost_3410.domain().get_order()),
+         m_mod_order(m_order),
+         m_base_point(gost_3410.domain().get_base_point(), m_order),
+         m_x(gost_3410.private_value()) {}
 
       size_t message_parts() const override { return 2; }
-      size_t message_part_size() const override { return order.bytes(); }
-      size_t max_input_bits() const override { return order.bits(); }
+      size_t message_part_size() const override { return m_order.bytes(); }
+      size_t max_input_bits() const override { return m_order.bits(); }
 
       secure_vector<byte> raw_sign(const byte msg[], size_t msg_len,
                                    RandomNumberGenerator& rng) override;
 
    private:
-      const PointGFp& base_point;
-      const BigInt& order;
-      const BigInt& x;
+      const BigInt& m_order;
+      Modular_Reducer m_mod_order;
+      Blinded_Point_Multiply m_base_point;
+      const BigInt& m_x;
    };
 
 secure_vector<byte>
@@ -124,26 +125,25 @@ GOST_3410_Signature_Operation::raw_sign(const byte msg[], size_t msg_len,
    {
    BigInt k;
    do
-      k.randomize(rng, order.bits()-1);
-   while(k >= order);
+      k.randomize(rng, m_order.bits()-1);
+   while(k >= m_order);
 
    BigInt e = decode_le(msg, msg_len);
 
-   e %= order;
+   e = m_mod_order.reduce(e);
    if(e == 0)
       e = 1;
 
-   PointGFp k_times_P = base_point * k;
+   const PointGFp k_times_P = m_base_point.blinded_multiply(k, rng);
    BOTAN_ASSERT(k_times_P.on_the_curve(), "GOST 34.10 k*g is on the curve");
 
-   BigInt r = k_times_P.get_affine_x() % order;
-
-   BigInt s = (r*x + k*e) % order;
+   const BigInt r = m_mod_order.reduce(k_times_P.get_affine_x());
+   const BigInt s = m_mod_order.reduce(r*m_x + k*e);
 
    if(r == 0 || s == 0)
       throw Invalid_State("GOST 34.10: r == 0 || s == 0");
 
-   secure_vector<byte> output(2*order.bytes());
+   secure_vector<byte> output(2*m_order.bytes());
    s.binary_encode(&output[output.size() / 2 - s.bytes()]);
    r.binary_encode(&output[output.size() - r.bytes()]);
    return output;
diff --git a/src/lib/rng/rng.h b/src/lib/rng/rng.h
index b1f78f75d..6ee67f66f 100644
--- a/src/lib/rng/rng.h
+++ b/src/lib/rng/rng.h
@@ -47,17 +47,35 @@ class BOTAN_DLL RandomNumberGenerator
          }
 
       /**
-      * Return a random byte
-      * @return random byte
+      * Only usable with POD types, only useful with integers
+      * get_random<u64bit>()
       */
-      byte next_byte()
+      template<typename T> T get_random()
          {
-         byte out;
-         this->randomize(&out, 1);
-         return out;
+         T r;
+         this->randomize(reinterpret_cast<byte*>(&r), sizeof(r));
+         return r;
          }
 
       /**
+      * Return a value in range [0,2^bits)
+      */
+      u64bit gen_mask(size_t bits)
+         {
+         if(bits == 0 || bits > 64)
+            throw std::invalid_argument("RandomNumberGenerator::gen_mask invalid argument");
+
+         const u64bit mask = ((1 << bits) - 1);
+         return this->get_random<u64bit>() & mask;
+         }
+
+      /**
+      * Return a random byte
+      * @return random byte
+      */
+      byte next_byte() { return get_random<byte>(); }
+
+      /**
       * Check whether this RNG is seeded.
       * @return true if this RNG was already seeded, false otherwise.
       */
diff --git a/src/tests/test_bigint.cpp b/src/tests/test_bigint.cpp
index 7bac56bc7..e6aa4a434 100644
--- a/src/tests/test_bigint.cpp
+++ b/src/tests/test_bigint.cpp
@@ -31,55 +31,6 @@ using namespace Botan;
 
 namespace {
 
-class Test_State
-   {
-   public:
-      void started(const std::string& /*msg*/) { m_tests_run++; }
-
-      void test_ran(const char* msg);
-
-      void failure(const char* test, const std::string& what_failed)
-         {
-         std::cout << "FAIL " << test << " " << what_failed << "\n";
-         m_tests_failed++;
-         }
-
-      size_t ran() const { return m_tests_run; }
-      size_t failed() const { return m_tests_failed; }
-   private:
-      size_t m_tests_run = 0, m_tests_failed = 0;
-   };
-
-#define BOTAN_CONFIRM_NOTHROW(block) do {                              \
-   try { block }                                                        \
-   catch(std::exception& e) {                                           \
-      _test.failure(BOTAN_CURRENT_FUNCTION, e.what());  \
-   } } while(0)                                                         \
-
-#define BOTAN_TEST(lhs, rhs, msg) do {                                     \
-   _test.started(msg);                                                  \
-   BOTAN_CONFIRM_NOTHROW({                                              \
-      const auto lhs_val = lhs;                                         \
-      const auto rhs_val = rhs;                                         \
-      const bool cmp = lhs_val == rhs_val;                              \
-      if(!cmp)                                                          \
-         {                                                              \
-         std::ostringstream fmt;                                        \
-         fmt << "expr '" << #lhs << " == " << #rhs << "' false, "       \
-             << "actually " << lhs_val << " " << rhs_val                \
-             << " (" << msg << ")";                                     \
-         _test.failure(BOTAN_CURRENT_FUNCTION, fmt.str()); \
-         }                                                              \
-      });                                                               \
-   } while(0)
-
-#define BOTAN_TEST_CASE(name, descr, block) size_t test_ ## name() {     \
-   Test_State _test;                                                    \
-   BOTAN_CONFIRM_NOTHROW(block);                                           \
-   test_report(descr, _test.ran(), _test.failed());                     \
-   return _test.failed();                                   \
-   }
-
 BOTAN_TEST_CASE(bigint_to_u32bit, "BigInt to_u32bit", {
    for(size_t i = 0; i != 32; ++i)
       {
diff --git a/src/tests/test_ecdsa.cpp b/src/tests/test_ecdsa.cpp
index a2ec8d115..8d385b4bf 100644
--- a/src/tests/test_ecdsa.cpp
+++ b/src/tests/test_ecdsa.cpp
@@ -25,7 +25,6 @@ size_t ecdsa_sig_kat(const std::string& group_id,
                      const std::string& x,
                      const std::string& hash,
                      const std::string& msg,
-                     const std::string& nonce,
                      const std::string& signature)
    {
    auto& rng = test_rng();
@@ -39,7 +38,7 @@ size_t ecdsa_sig_kat(const std::string& group_id,
    PK_Signer sign(ecdsa, padding);
 
    return validate_signature(verify, sign, "ECDSA/" + group_id + '/' + hash,
-                             msg, rng, nonce, signature);
+                             msg, rng, signature);
    }
 
 }
@@ -53,7 +52,7 @@ size_t test_ecdsa()
    fails += run_tests_bb(ecdsa_sig, "ECDSA Signature", "Signature", false,
              [](std::map<std::string, std::string> m) -> size_t
              {
-             return ecdsa_sig_kat(m["Group"], m["X"], m["Hash"], m["Msg"], m["Nonce"], m["Signature"]);
+             return ecdsa_sig_kat(m["Group"], m["X"], m["Hash"], m["Msg"], m["Signature"]);
              });
 
    return fails;
diff --git a/src/tests/tests.cpp b/src/tests/tests.cpp
index 63e6761ac..763417209 100644
--- a/src/tests/tests.cpp
+++ b/src/tests/tests.cpp
@@ -304,6 +304,7 @@ int main(int argc, char* argv[])
    DEF_TEST(mceliece);
 
    DEF_TEST(ecc_unit);
+   DEF_TEST(ecc_randomized);
    DEF_TEST(ecdsa_unit);
    DEF_TEST(ecdh_unit);
    DEF_TEST(pk_keygen);
diff --git a/src/tests/tests.h b/src/tests/tests.h
index 88102f289..14ec5a17b 100644
--- a/src/tests/tests.h
+++ b/src/tests/tests.h
@@ -16,6 +16,7 @@
 #include <string>
 #include <vector>
 #include <iostream>
+#include <sstream>
 
 Botan::RandomNumberGenerator& test_rng();
 
@@ -45,7 +46,69 @@ typedef std::function<size_t ()> test_fn;
 size_t run_tests(const std::vector<std::pair<std::string, test_fn>>& tests);
 void test_report(const std::string& name, size_t ran, size_t failed);
 
-#define TEST(expr, msg) do { if(!(expr)) { ++fails; std::cout << msg; } while(0)
+class Test_State
+   {
+   public:
+      void started(const std::string& /*msg*/) { m_tests_run++; }
+
+      void test_ran(const char* msg);
+
+      void failure(const char* test, const std::string& what_failed)
+         {
+         std::cout << "FAIL " << test << " " << what_failed << "\n";
+         m_tests_failed++;
+         }
+
+      size_t ran() const { return m_tests_run; }
+      size_t failed() const { return m_tests_failed; }
+   private:
+      size_t m_tests_run = 0, m_tests_failed = 0;
+   };
+
+#define BOTAN_CONFIRM_NOTHROW(block) do {                              \
+   try { block }                                                        \
+   catch(std::exception& e) {                                           \
+      _test.failure(BOTAN_CURRENT_FUNCTION, e.what());  \
+   } } while(0)                                                         \
+
+#define BOTAN_TEST(lhs, rhs, msg) do {                                     \
+   _test.started(msg);                                                  \
+   BOTAN_CONFIRM_NOTHROW({                                              \
+      const auto lhs_val = lhs;                                         \
+      const auto rhs_val = rhs;                                         \
+      const bool cmp = lhs_val == rhs_val;                              \
+      if(!cmp)                                                          \
+         {                                                              \
+         std::ostringstream fmt;                                        \
+         fmt << "expr '" << #lhs << " == " << #rhs << "' false, "       \
+             << "actually " << lhs_val << " " << rhs_val                \
+             << " (" << msg << ")";                                     \
+         _test.failure(BOTAN_CURRENT_FUNCTION, fmt.str()); \
+         }                                                              \
+      });                                                               \
+   } while(0)
+
+#define BOTAN_CONFIRM(expr, msg) do {                                     \
+   _test.started(msg);                                                  \
+   BOTAN_CONFIRM_NOTHROW({                                              \
+      const bool expr_val = expr;                                         \
+      if(!expr_val)                                                     \
+         {                                                              \
+         std::ostringstream fmt;                                        \
+         fmt << "expr '" << #expr << " false (" << msg << ")";           \
+         _test.failure(BOTAN_CURRENT_FUNCTION, fmt.str());              \
+         }                                                              \
+      });                                                               \
+   } while(0)
+
+#define BOTAN_TEST_CASE(name, descr, block) size_t test_ ## name() {     \
+   Test_State _test;                                                    \
+   BOTAN_CONFIRM_NOTHROW(block);                                           \
+   test_report(descr, _test.ran(), _test.failed());                     \
+   return _test.failed();                                   \
+   }
+
+//#define TEST(expr, msg) do { if(!(expr)) { ++fails; std::cout << msg; } while(0)
 
 #define TEST_DATA_DIR     "src/tests/data"
 #define TEST_DATA_DIR_PK  "src/tests/data/pubkey"
@@ -75,6 +138,7 @@ size_t test_dh();
 size_t test_dlies();
 size_t test_elgamal();
 size_t test_ecc_pointmul();
+size_t test_ecc_random();
 size_t test_ecdsa();
 size_t test_gost_3410();
 size_t test_curve25519();
@@ -94,6 +158,7 @@ size_t test_pk_keygen();
 size_t test_bigint();
 
 size_t test_ecc_unit();
+size_t test_ecc_randomized();
 size_t test_ecdsa_unit();
 size_t test_ecdh_unit();
 
diff --git a/src/tests/unit_ecc.cpp b/src/tests/unit_ecc.cpp
index 8498e3b43..3cb1436d3 100644
--- a/src/tests/unit_ecc.cpp
+++ b/src/tests/unit_ecc.cpp
@@ -816,9 +816,9 @@ size_t test_curve_cp_ctor()
    return 0;
    }
 
-size_t ecc_randomized_test()
-   {
-   const std::vector<std::string> groups = {
+namespace {
+
+const std::vector<std::string> ec_groups = {
       "brainpool160r1",
       "brainpool192r1",
       "brainpool224r1",
@@ -849,11 +849,16 @@ size_t ecc_randomized_test()
       "x962_p239v3"
    };
 
+}
+
+}
+
+BOTAN_TEST_CASE(ecc_randomized, "ECC Randomized", {
    auto& rng = test_rng();
    size_t fails = 0;
    size_t tests = 0;
 
-   for(auto&& group_name : groups)
+   for(auto&& group_name : ec_groups)
       {
       EC_Group group(group_name);
 
@@ -861,8 +866,8 @@ size_t ecc_randomized_test()
       const BigInt& group_order = group.get_order();
 
       const PointGFp inf = base_point * group_order;
-      CHECK(inf.is_zero());
-      CHECK(inf.on_the_curve());
+      BOTAN_CONFIRM(inf.is_zero(), "Group math ok");
+      BOTAN_CONFIRM(inf.on_the_curve(), "Infinity still on the curve");
 
       try
          {
@@ -870,6 +875,9 @@ size_t ecc_randomized_test()
             {
             ++tests;
 
+            const size_t h = 1 + (rng.next_byte() % 8);
+            Blinded_Point_Multiply blind(base_point, group_order, h);
+
             const BigInt a = BigInt::random_integer(rng, 2, group_order);
             const BigInt b = BigInt::random_integer(rng, 2, group_order);
             const BigInt c = a + b;
@@ -878,16 +886,24 @@ size_t ecc_randomized_test()
             const PointGFp Q = base_point * b;
             const PointGFp R = base_point * c;
 
+            const PointGFp P1 = blind.blinded_multiply(a, rng);
+            const PointGFp Q1 = blind.blinded_multiply(b, rng);
+            const PointGFp R1 = blind.blinded_multiply(c, rng);
+
             const PointGFp A1 = P + Q;
             const PointGFp A2 = Q + P;
 
-            CHECK(A1 == R);
-            CHECK(A2 == R);
-            CHECK(P.on_the_curve());
-            CHECK(Q.on_the_curve());
-            CHECK(R.on_the_curve());
-            CHECK(A1.on_the_curve());
-            CHECK(A2.on_the_curve());
+            BOTAN_TEST(A1, R, "Addition");
+            BOTAN_TEST(A2, R, "Addition");
+            BOTAN_CONFIRM(P.on_the_curve(), "On the curve");
+            BOTAN_CONFIRM(Q.on_the_curve(), "On the curve");
+            BOTAN_CONFIRM(R.on_the_curve(), "On the curve");
+            BOTAN_CONFIRM(A1.on_the_curve(), "On the curve");
+            BOTAN_CONFIRM(A2.on_the_curve(), "On the curve");
+
+            BOTAN_TEST(P, P1, "P1");
+            BOTAN_TEST(Q, Q1, "Q1");
+            BOTAN_TEST(R, R1, "R1");
             }
          }
       catch(std::exception& e)
@@ -896,12 +912,8 @@ size_t ecc_randomized_test()
          ++fails;
          }
       }
+   });
 
-   test_report("ECC Randomized", tests, fails);
-   return fails;
-   }
-
-}
 
 size_t test_ecc_unit()
    {
@@ -934,8 +946,6 @@ size_t test_ecc_unit()
 
    test_report("ECC", 0, fails);
 
-   ecc_randomized_test();
-
    return fails;
    }