Add an implementation of McEliece encryption based on HyMES

(https://www.rocq.inria.fr/secret/CBCrypto/index.php?pg=hymes).
The original version is LGPL but cryptsource GmbH has secured
permission to release it under a BSD license. Also includes the
Overbeck CCA2 message encoding scheme.

23 files changed, 3366 insertions, 13 deletions

diff --git a/src/lib/pubkey/mce/binary_matrix.cpp b/src/lib/pubkey/mce/binary_matrix.cpp
new file mode 100644
index 000000000..4670a75c4
--- /dev/null
+++ b/src/lib/pubkey/mce/binary_matrix.cpp
@@ -0,0 +1,97 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/internal/binary_matrix.h>
+#include <botan/internal/xor_buf.h>
+
+namespace Botan {
+
+binary_matrix::binary_matrix (u32bit rown, u32bit coln)
+   {
+   m_coln = coln;
+   m_rown = rown;
+   m_rwdcnt = (1 + (m_coln - 1) / BITS_PER_U32);
+   m_alloc_size = m_rown * (*this).m_rwdcnt * sizeof (u32bit);
+   m_elem = std::vector<u32bit>((*this).m_alloc_size/4);
+   }
+
+void binary_matrix::row_xor(u32bit a, u32bit b)
+   {
+   u32bit i;
+   for(i=0;i<m_rwdcnt;i++)
+      {
+      m_elem[a*m_rwdcnt+i]^=m_elem[b*m_rwdcnt+i];
+      }
+   }
+
+//the matrix is reduced from LSB...(from right)
+secure_vector<int> binary_matrix::row_reduced_echelon_form()
+   {
+   u32bit i, failcnt, findrow, max=m_coln - 1;
+
+   secure_vector<int> perm(m_coln);
+   for(i=0;i<m_coln;i++)
+      {
+      perm[i]=i;//initialize permutation.
+      }
+   failcnt = 0;
+
+   for(i=0;i<m_rown;i++,max--)
+      {
+      findrow=0;
+      for(u32bit j=i;j<m_rown;j++)
+         {
+         if(coef(j,max))
+            {
+            if (i!=j)//not needed as ith row is 0 and jth row is 1.
+               row_xor(i,j);//xor to the row.(swap)?
+            findrow=1;
+            break;
+            }//largest value found (end if)
+         }
+
+      if(!findrow)//if no row with a 1 found then swap last column and the column with no 1 down.
+         {
+         perm[m_coln - m_rown - 1 - failcnt] = max;
+         failcnt++;
+         if (!max)
+            {
+            //CSEC_FREE_MEM_CHK_SET_NULL(*p_perm);
+            //CSEC_THR_RETURN();
+            perm.resize(0);
+            }
+         i--;
+         }
+      else
+         {
+         perm[i+m_coln - m_rown] = max;
+         for(u32bit j=i+1;j<m_rown;j++)//fill the column downwards with 0's
+            {
+            if(coef(j,(max)))
+               {
+               row_xor(j,i);//check the arg. order.
+               }
+            }
+
+         for(int j=i-1;j>=0;j--)//fill the column with 0's upwards too.
+            {
+            if(coef(j,(max)))
+               {
+               row_xor(j,i);
+               }
+            }
+         }
+      }//end for(i)
+   return perm;
+   }
+
+
+} // end namespace Botan
diff --git a/src/lib/pubkey/mce/binary_matrix.h b/src/lib/pubkey/mce/binary_matrix.h
new file mode 100644
index 000000000..a5438666d
--- /dev/null
+++ b/src/lib/pubkey/mce/binary_matrix.h
@@ -0,0 +1,61 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef BOTAN_BINARY_MATRIX_H__
+#define BOTAN_BINARY_MATRIX_H__
+
+#include <botan/secmem.h>
+
+namespace Botan {
+
+#define BITS_PER_U32 (8 * sizeof (u32bit))
+
+struct binary_matrix
+   {
+   public:
+      binary_matrix(u32bit m_rown, u32bit m_coln);
+
+      void row_xor(u32bit a, u32bit b);
+      secure_vector<int> row_reduced_echelon_form();
+
+      /**
+      * return the coefficient out of F_2
+      */
+      u32bit coef(u32bit i, u32bit j)
+         {
+         return (m_elem[(i) * m_rwdcnt + (j) / BITS_PER_U32] >> (j % BITS_PER_U32)) & 1;
+         };
+
+      void set_coef_to_one(u32bit i, u32bit j)
+         {
+         m_elem[(i) * m_rwdcnt + (j) / BITS_PER_U32] |= (1UL << ((j) % BITS_PER_U32)) ;
+         };
+
+      void toggle_coeff(u32bit i, u32bit j)
+         {
+         m_elem[(i) * m_rwdcnt + (j) / BITS_PER_U32] ^= (1UL << ((j) % BITS_PER_U32)) ;
+         }
+
+      void set_to_zero()
+         {
+         std::memset(&m_elem[0], 0, m_alloc_size);
+         }
+
+      u32bit m_rown;  // number of rows.
+      u32bit m_coln; // number of columns.
+      u32bit m_rwdcnt; // number of words in a row
+      u32bit m_alloc_size; // number of allocated bytes
+      std::vector<u32bit> m_elem;
+   };
+
+}
+
+#endif
diff --git a/src/lib/pubkey/mce/code_based_key_gen.cpp b/src/lib/pubkey/mce/code_based_key_gen.cpp
new file mode 100644
index 000000000..997e51685
--- /dev/null
+++ b/src/lib/pubkey/mce/code_based_key_gen.cpp
@@ -0,0 +1,185 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/internal/code_based_key_gen.h>
+#include <botan/code_based_util.h>
+#include <botan/gf2m_rootfind_dcmp.h>
+#include <botan/internal/binary_matrix.h>
+#include <botan/loadstor.h>
+#include <botan/polyn_gf2m.h>
+
+namespace Botan {
+
+namespace {
+
+void randomize_support(u32bit n, std::vector<gf2m> & L, RandomNumberGenerator & rng)
+   {
+   unsigned int i, j;
+   gf2m_small_m::gf2m tmp;
+
+   for (i = 0; i < n; ++i)
+      {
+
+      gf2m_small_m::gf2m rnd;
+      rng.randomize(reinterpret_cast<byte*>(&rnd), sizeof(rnd));
+      j = rnd % n; // no rejection sampling, but for useful code-based parameters with n <= 13 this seem tolerable
+
+      tmp = L[j];
+      L[j] = L[i];
+      L[i] = tmp;
+      }
+   }
+
+std::unique_ptr<binary_matrix> generate_R(std::vector<gf2m> &L, polyn_gf2m* g, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field, u32bit code_length, u32bit t )
+   {
+   //L- Support
+   //t- Number of errors
+   //n- Length of the Goppa code
+   //m- The extension degree of the GF
+   //g- The generator polynomial.
+   gf2m_small_m::gf2m x,y;
+   u32bit i,j,k,r,n;
+   std::vector<int> Laux(code_length);
+   n=code_length;
+   r=t*sp_field->get_extension_degree();
+
+   binary_matrix H(r, n) ;
+
+   for(i=0;i< n;i++)
+      {
+      x = g->eval(lex_to_gray(L[i]));//evaluate the polynomial at the point L[i].
+      x = sp_field->gf_inv(x);
+      y = x;
+      for(j=0;j<t;j++)
+         {
+         for(k=0;k<sp_field->get_extension_degree();k++)
+            {
+            if(y & (1<<k))
+               {
+               //the co-eff. are set in 2^0,...,2^11 ; 2^0,...,2^11 format along the rows/cols?
+               H.set_coef_to_one(j*sp_field->get_extension_degree()+ k,i);
+               }
+            }
+         y = sp_field->gf_mul(y,lex_to_gray(L[i]));
+         }
+      }//The H matrix is fed.
+
+   secure_vector<int> perm = H.row_reduced_echelon_form();
+   if (perm.size() == 0)
+      {
+      // result still is NULL
+      throw Invalid_State("could not bring matrix in row reduced echelon form");
+      }
+
+   std::unique_ptr<binary_matrix> result(new binary_matrix(n-r,r)) ;
+   for (i = 0; i < (*result).m_rown; ++i)
+      {
+      for (j = 0; j < (*result).m_coln; ++j)
+         {
+         if (H.coef(j,perm[i]))
+            {
+            result->toggle_coeff(i,j);
+            }
+         }
+      }
+   for (i = 0; i < code_length; ++i)
+      {
+      Laux[i] = L[perm[i]];
+      }
+   for (i = 0; i < code_length; ++i)
+      {
+      L[i] = Laux[i];
+      }
+   return result;
+   }
+}
+
+McEliece_PrivateKey generate_mceliece_key( RandomNumberGenerator & rng, u32bit ext_deg, u32bit code_length, u32bit t)
+   {
+   u32bit i, j, k, l;
+   std::unique_ptr<binary_matrix> R;
+
+   u32bit codimension = t * ext_deg;
+   if(code_length <= codimension)
+      {
+      throw Invalid_Argument("invalid McEliece parameters");
+      }
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field ( new Gf2m_Field(ext_deg ));
+
+   //pick the support.........
+   std::vector<gf2m> L(code_length);
+
+   for(i=0;i<code_length;i++)
+      {
+      L[i]=i;
+      }
+   randomize_support(code_length,L,rng);
+   polyn_gf2m g(sp_field); // create as zero
+   bool success = false;
+   do
+      {
+      // create a random irreducible polynomial
+      g = polyn_gf2m (t, rng, sp_field);
+
+      try{
+      R = generate_R(L,&g, sp_field, code_length, t);
+      success = true;
+      }
+      catch(const Invalid_State &)
+         {
+         }
+      } while (!success);
+
+   std::vector<polyn_gf2m> sqrtmod = polyn_gf2m::sqrt_mod_init( g);
+   std::vector<polyn_gf2m> F = syndrome_init(g, L, code_length);
+
+   // Each F[i] is the (precomputed) syndrome of the error vector with
+   // a single '1' in i-th position.
+   // We do not store the F[i] as polynomials of degree t , but
+   // as binary vectors of length ext_deg * t (this will
+   // speed up the syndrome computation)
+   //
+   //
+   std::vector<u32bit> H(bit_size_to_32bit_size(codimension) * code_length );
+   u32bit* sk = &H[0];
+   for (i = 0; i < code_length; ++i)
+      {
+      for (l = 0; l < t; ++l)
+         {
+         k = (l * ext_deg) / 32;
+         j = (l * ext_deg) % 32;
+         sk[k] ^= F[i].get_coef( l) << j;
+         if (j + ext_deg > 32)
+            {
+            sk[k + 1] ^= F[i].get_coef( l) >> (32 - j);
+            }
+         }
+      sk += bit_size_to_32bit_size(codimension);
+      }
+
+   // We need the support L for decoding (decryption). In fact the
+   // inverse is needed
+
+   std::vector<gf2m> Linv(code_length) ;
+   for (i = 0; i < code_length; ++i)
+      {
+      Linv[L[i]] = i;
+      }
+   std::vector<byte> pubmat (R->m_alloc_size);
+   for(i = 0; i < R->m_alloc_size/4; i++)
+      {
+      store_le(R->m_elem[i], &pubmat[i*4] );
+      }
+
+   return McEliece_PrivateKey(g, H, sqrtmod, Linv, pubmat);
+   }
+
+}
diff --git a/src/lib/pubkey/mce/code_based_key_gen.h b/src/lib/pubkey/mce/code_based_key_gen.h
new file mode 100644
index 000000000..80201ebf5
--- /dev/null
+++ b/src/lib/pubkey/mce/code_based_key_gen.h
@@ -0,0 +1,26 @@
+
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+#ifndef __code_based_key_gen__H_
+#define __code_based_key_gen__H_
+
+#include <botan/mceliece_key.h>
+
+namespace Botan {
+
+McEliece_PrivateKey generate_mceliece_key(RandomNumberGenerator &rng,
+                                          u32bit ext_deg,
+                                          u32bit code_length,
+                                          u32bit t);
+
+}
+
+#endif
diff --git a/src/lib/pubkey/mce/code_based_util.h b/src/lib/pubkey/mce/code_based_util.h
new file mode 100644
index 000000000..fd8bcaa49
--- /dev/null
+++ b/src/lib/pubkey/mce/code_based_util.h
@@ -0,0 +1,57 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef BOTAN_CODE_BASED_UTIL_H__
+#define BOTAN_CODE_BASED_UTIL_H__
+
+#include <botan/gf2m_small_m.h>
+
+namespace Botan {
+
+/**
+* Expand an input to a bit mask depending on it being being zero or non-zero
+* @ param tst the input
+* @return the mask 0xFFFF if tst is non-zero and 0 otherwise
+*/
+template<typename T>
+u16bit expand_mask_16bit(T tst)
+   {
+   const u16bit result = (tst != 0);
+   return ~(result - 1);
+   }
+
+inline gf2m_small_m::gf2m gray_to_lex(gf2m_small_m::gf2m gray)
+   {
+   gf2m_small_m::gf2m result = gray ^ (gray>>8);
+   result ^= (result >> 4);
+   result ^= (result >> 2);
+   result ^= (result >> 1);
+   return result;
+   }
+
+inline gf2m_small_m::gf2m lex_to_gray(gf2m_small_m::gf2m lex)
+   {
+   return (lex>>1) ^ lex;
+   }
+
+inline u32bit bit_size_to_byte_size(u32bit bit_size)
+   {
+   return (bit_size - 1) / 8 + 1;
+   }
+
+inline u32bit bit_size_to_32bit_size(u32bit bit_size)
+   {
+   return (bit_size - 1) / 32 + 1;
+   }
+
+}
+
+#endif
diff --git a/src/lib/pubkey/mce/gf2m_rootfind_dcmp.cpp b/src/lib/pubkey/mce/gf2m_rootfind_dcmp.cpp
new file mode 100644
index 000000000..f6ff8e498
--- /dev/null
+++ b/src/lib/pubkey/mce/gf2m_rootfind_dcmp.cpp
@@ -0,0 +1,317 @@
+/**
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/gf2m_rootfind_dcmp.h>
+#include <botan/gf2m_small_m.h>
+#include <botan/internal/bit_ops.h>
+#include <botan/code_based_util.h>
+
+namespace Botan
+{
+
+namespace {
+
+  u32bit patch_root_array(
+    gf2m* res_root_arr,
+    u32bit res_root_arr_len,
+    u32bit root_pos)
+  {
+    volatile u32bit i;
+    volatile gf2m patch_elem = 0x01;
+    volatile gf2m cond_mask = (root_pos == res_root_arr_len);
+    cond_mask = expand_mask_16bit(cond_mask);
+    cond_mask = ~cond_mask; /* now cond = 1 if not enough roots */
+    patch_elem &= cond_mask;
+    for(i = 0; i < res_root_arr_len; i++)
+    {
+
+      gf2m masked_patch_elem = (patch_elem++) & cond_mask;
+      res_root_arr[i] ^= masked_patch_elem++;
+    }
+    return res_root_arr_len;
+  }
+
+struct gf2m_decomp_rootfind_state
+{
+    gf2m_decomp_rootfind_state(const polyn_gf2m & p_polyn, u32bit code_length);
+
+  void calc_LiK(const polyn_gf2m & sigma);
+  gf2m calc_Fxj_j_neq_0( const polyn_gf2m & sigma, gf2m j_gray);
+  void calc_next_Aij();
+  void calc_Ai_zero(const polyn_gf2m & sigma);
+  secure_vector<gf2m> find_roots(const polyn_gf2m & sigma);
+   u32bit get_code_length() const { return code_length; };
+    u32bit code_length;
+  secure_vector<gf2m> m_Lik; // size is outer_summands * m
+  secure_vector<gf2m> m_Aij; // ...
+  u32bit m_outer_summands;
+  gf2m m_j;
+  gf2m m_j_gray;
+  gf2m m_sigma_3_l;
+  gf2m m_sigma_3_neq_0_mask;
+} ;
+
+
+/*
+ * !! Attention: assumes gf2m is 16bit !!
+ */
+gf2m brootf_decomp__gray_to_lex(gf2m gray)
+{
+  static_assert(sizeof(gf2m) == 2, "Expected size");
+  gf2m result = gray ^ (gray>>8);
+  result ^= (result >> 4);
+  result ^= (result >> 2);
+  result ^= (result >> 1);
+  return result;
+}
+
+/**
+ * calculates ceil((t-4)/5) = outer_summands - 1
+ */
+u32bit brootf_decomp__calc_sum_limit(u32bit t)
+{
+  u32bit result;
+  if(t < 4)
+  {
+    return 0;
+  }
+  result = t - 4;
+  result += 4;
+  result /= 5;
+  return result;
+}
+}
+
+secure_vector<gf2m> find_roots_gf2m_decomp(const polyn_gf2m & polyn, u32bit code_length)
+{
+  gf2m_decomp_rootfind_state state(polyn, code_length);
+  return state.find_roots(polyn);
+
+}
+
+gf2m_decomp_rootfind_state::gf2m_decomp_rootfind_state(const polyn_gf2m & polyn, u32bit the_code_length)
+  :code_length(the_code_length)
+{
+
+  gf2m coeff_3;
+  gf2m coeff_head;
+  std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field = polyn.get_sp_field();
+  int deg_sigma = polyn.get_degree();
+  if(deg_sigma <= 3)
+  {
+    throw std::exception();
+  }
+  this->m_j = 0;
+  coeff_3 = polyn.get_coef( 3);
+  coeff_head = polyn.get_coef( deg_sigma); /* dummy value for SCA CM */
+  if(coeff_3 != 0)
+  {
+    this->m_sigma_3_l = sp_field->gf_l_from_n(coeff_3);
+    this->m_sigma_3_neq_0_mask = 0xFFFF;
+  }
+  else
+  {
+    // dummy value needed for timing countermeasure
+    this->m_sigma_3_l = sp_field->gf_l_from_n(coeff_head);
+    this->m_sigma_3_neq_0_mask = 0 ;
+  }
+
+  this->m_outer_summands =  1 + brootf_decomp__calc_sum_limit(deg_sigma);
+  this->m_Lik.resize(this->m_outer_summands * sp_field->get_extension_degree());
+  this->m_Aij.resize(this->m_outer_summands);
+}
+
+void gf2m_decomp_rootfind_state::calc_Ai_zero(const polyn_gf2m & sigma)
+{
+  u32bit i;
+  /*
+   * this function assumes this the first gray code element is zero
+   */
+  for(i = 0; i < this->m_outer_summands; i++)
+  {
+    this->m_Aij[i] = sigma.get_coef(5*i);
+  }
+  this->m_j = 0;
+  this->m_j_gray = 0;
+}
+
+void gf2m_decomp_rootfind_state::calc_next_Aij()
+{
+  /*
+   * upon function entry, we have in the state j, Aij.
+   * first thing, we declare Aij Aij_minusone and increase j.
+   * Case j=0 upon function entry also included, then Aij contains A_{i,j=0}.
+   */
+  u32bit i;
+  gf2m diff, new_j_gray;
+  u32bit Lik_pos_base;
+
+  this->m_j++;
+
+  new_j_gray =  lex_to_gray(this->m_j);
+
+  if(this->m_j & 1)  /* half of the times */
+  {
+    Lik_pos_base = 0;
+  }
+  else if(this->m_j & 2) /* one quarter of the times */
+  {
+    Lik_pos_base = this->m_outer_summands;
+  }
+  else if( this->m_j & 4) /* one eighth of the times */
+  {
+    Lik_pos_base = this->m_outer_summands * 2;
+  }
+  else if( this->m_j & 8) /* one sixteenth of the times */
+  {
+    Lik_pos_base = this->m_outer_summands * 3;
+  }
+  else if( this->m_j & 16) /* ... */
+  {
+    Lik_pos_base = this->m_outer_summands * 4;
+  }
+  else
+  {
+    gf2m delta_offs = 5;
+    diff = this->m_j_gray ^ new_j_gray;
+    while(((((gf2m)1) << delta_offs) & diff) == 0)
+    {
+      delta_offs++;
+
+    }
+    Lik_pos_base = delta_offs * this->m_outer_summands;
+  }
+  this->m_j_gray = new_j_gray;
+
+  i = 0;
+  for(; i < this->m_outer_summands; i++)
+  {
+    this->m_Aij[i] ^= this->m_Lik[Lik_pos_base + i];
+  }
+
+}
+void gf2m_decomp_rootfind_state::calc_LiK(const polyn_gf2m & sigma)
+{
+  std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field = sigma.get_sp_field();
+  u32bit i, k, d;
+  d = sigma.get_degree();
+  for(k = 0; k < sp_field->get_extension_degree(); k++)
+  {
+    u32bit Lik_pos_base = k * this->m_outer_summands;
+    gf2m alpha_l_k_tt2_ttj[4];
+    alpha_l_k_tt2_ttj[0] = sp_field->gf_l_from_n( ((gf2m)1) << k);
+    alpha_l_k_tt2_ttj[1] = sp_field->gf_mul_rrr(alpha_l_k_tt2_ttj[0], alpha_l_k_tt2_ttj[0]);
+    alpha_l_k_tt2_ttj[2] = sp_field->gf_mul_rrr(alpha_l_k_tt2_ttj[1],alpha_l_k_tt2_ttj[1] );
+
+    alpha_l_k_tt2_ttj[3] = sp_field->gf_mul_rrr(alpha_l_k_tt2_ttj[2], alpha_l_k_tt2_ttj[2]);
+    for(i = 0; i < this->m_outer_summands; i++)
+    {
+      u32bit j;
+      u32bit five_i = 5*i;
+      u32bit Lik_pos = Lik_pos_base + i;
+      this->m_Lik[Lik_pos] = 0;
+      for(j = 0; j <= 3; j++)
+      {
+        gf2m f, x;
+        u32bit f_ind = five_i + ((u32bit)1<<j);
+        if(f_ind > d)
+        {
+          break;
+        }
+        f = sigma.get_coef( f_ind);
+
+        x = sp_field->gf_mul_zrz(alpha_l_k_tt2_ttj[j], f);
+        this->m_Lik[Lik_pos] ^= x;
+      }
+    }
+  }
+}
+
+gf2m gf2m_decomp_rootfind_state::calc_Fxj_j_neq_0( const polyn_gf2m & sigma, gf2m j_gray)
+{
+  //needs the A_{ij} to compute F(x)_j
+  gf2m sum = 0;
+  u32bit i;
+  std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field = sigma.get_sp_field();
+  gf2m xl_j_tt_5i, xl_j_tt_5,  xl_gray_tt_3;
+  gf2m jl_gray;
+  jl_gray = sp_field->gf_l_from_n(j_gray);
+  xl_j_tt_5 = sp_field->gf_square_rr(jl_gray);
+  xl_gray_tt_3 = sp_field->gf_mul_rrr(xl_j_tt_5, jl_gray);
+  xl_j_tt_5 = sp_field->gf_mul_rrr(xl_j_tt_5, xl_gray_tt_3);
+
+
+  sum = sp_field->gf_mul_nrr(xl_gray_tt_3, this->m_sigma_3_l);
+  sum &= this->m_sigma_3_neq_0_mask;
+  /* here, we rely on compiler to be unable to optimize
+   * for the state->sigma_3_neq_0_mask value
+   */
+  /* treat i = 0 special: */
+  sum ^= this->m_Aij[0];
+  /* treat i = 1 special also */
+  if(this->m_outer_summands > 1)
+  {
+    gf2m x;
+    xl_j_tt_5i = xl_j_tt_5;
+    x = sp_field->gf_mul_zrz(xl_j_tt_5, this->m_Aij[1]); /* x_j^{5i} A_i^j */
+    sum ^= x;
+  }
+  for(i = 2; i < this->m_outer_summands; i++)
+  {
+    gf2m x;
+    xl_j_tt_5i = sp_field->gf_mul_rrr(xl_j_tt_5i, xl_j_tt_5);
+    // now x_j_tt_5i lives up to its name
+    x = sp_field->gf_mul_zrz(xl_j_tt_5i, this->m_Aij[i]); /* x_j^{5i} A_i^(j) */
+    sum ^= x;
+  }
+  return sum;
+}
+
+
+
+
+secure_vector<gf2m> gf2m_decomp_rootfind_state::find_roots(const polyn_gf2m & sigma)
+{
+
+  secure_vector<gf2m> result(sigma.get_degree());
+  u32bit root_pos = 0;
+
+  this->calc_Ai_zero(sigma);
+  this->calc_LiK(sigma);
+  do
+  {
+    gf2m eval_result;
+    if(this->m_j_gray == 0)
+    {
+      eval_result = sigma.get_coef( 0);
+    }
+    else
+    {
+      eval_result = this->calc_Fxj_j_neq_0(sigma, this->m_j_gray);
+    }
+
+    if(eval_result == 0)
+    {
+
+      result[root_pos] = this->m_j_gray;
+      root_pos++;
+
+    }
+    if(this->m_j + static_cast<u32bit>(1) == this->get_code_length())
+    {
+      break;
+    }
+    this->calc_next_Aij();
+  }while(1);
+
+  // side channel / fault attack countermeasure:
+  root_pos = patch_root_array(&result[0], result.size(), root_pos);
+  result.resize(root_pos);
+ return result;
+}
+} // end namespace Botan
diff --git a/src/lib/pubkey/mce/gf2m_rootfind_dcmp.h b/src/lib/pubkey/mce/gf2m_rootfind_dcmp.h
new file mode 100644
index 000000000..8e1e1cd75
--- /dev/null
+++ b/src/lib/pubkey/mce/gf2m_rootfind_dcmp.h
@@ -0,0 +1,24 @@
+/**
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef __gf2m_rootfind_dcmp__H_
+#define __gf2m_rootfind_dcmp__H_
+
+#include <botan/polyn_gf2m.h>
+namespace Botan
+{
+  /**
+   * Find the roots of a polynomial over GF(2^m) using the method by Federenko
+   * et al.
+   */
+  secure_vector<gf2m> find_roots_gf2m_decomp(const polyn_gf2m & polyn, u32bit code_length);
+
+} // end namespace Botan
+
+#endif /* h-guard */
diff --git a/src/lib/pubkey/mce/gf2m_small_m.cpp b/src/lib/pubkey/mce/gf2m_small_m.cpp
new file mode 100644
index 000000000..76e5b0efd
--- /dev/null
+++ b/src/lib/pubkey/mce/gf2m_small_m.cpp
@@ -0,0 +1,133 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ */
+
+
+#include <botan/gf2m_small_m.h>
+#include <botan/code_based_util.h>
+#include <string>
+
+namespace Botan {
+
+namespace gf2m_small_m {
+
+#define MAX_EXT_DEG 16
+
+namespace {
+
+unsigned int prim_poly[MAX_EXT_DEG + 1] = {
+   01,		/* extension degree 0 (!) never used */
+   03,		/* extension degree 1 (!) never used */
+   07, 		/* extension degree 2 */
+   013, 		/* extension degree 3 */
+   023, 		/* extension degree 4 */
+   045, 		/* extension degree 5 */
+   0103, 		/* extension degree 6 */
+   0203, 		/* extension degree 7 */
+   0435, 		/* extension degree 8 */
+   01041, 		/* extension degree 9 */
+   02011,		/* extension degree 10 */
+   04005,		/* extension degree 11 */
+   010123,		/* extension degree 12 */
+   020033,		/* extension degree 13 */
+   042103,		/* extension degree 14 */
+   0100003,		/* extension degree 15 */
+   0210013		/* extension degree 16 */
+};
+
+}
+
+u32bit encode_gf2m(gf2m to_enc, byte* mem)
+   {
+   mem[0] = to_enc >> 8;
+   mem[1] = to_enc & 0xFF;
+   return sizeof(to_enc);
+   }
+
+gf2m decode_gf2m(const byte* mem)
+   {
+   gf2m result;
+   result = mem[0] << 8;
+   result |= mem[1];
+   return result;
+   }
+
+// construct the table gf_exp[i]=alpha^i
+void Gf2m_Field::init_exp()
+   {
+   m_gf_exp_table.resize(1 << get_extension_degree());
+
+   m_gf_exp_table[0] = 1;
+   for(size_t i = 1; i < gf_ord(); ++i)
+      {
+      m_gf_exp_table[i] = m_gf_exp_table[i - 1] << 1;
+      if (m_gf_exp_table[i - 1] & (1 << (get_extension_degree()-1)))
+         {
+         m_gf_exp_table[i] ^= prim_poly[get_extension_degree()];
+         }
+      }
+
+   // hack for the multiplication
+   m_gf_exp_table[gf_ord()] = 1;
+   }
+
+// construct the table gf_log[alpha^i]=i
+void Gf2m_Field::init_log()
+   {
+   m_gf_log_table.resize(1 << get_extension_degree());
+
+   m_gf_log_table[0] = gf_ord(); // log of 0 par convention
+   for (size_t i = 0; i < gf_ord() ; ++i)
+      {
+      m_gf_log_table[m_gf_exp_table[i]] = i;
+      }
+   }
+
+
+Gf2m_Field::Gf2m_Field(size_t extdeg)
+   {
+   if(extdeg < 2 || extdeg > MAX_EXT_DEG)
+      throw std::runtime_error("Gf2m_Field does not support degree " + std::to_string(extdeg));
+
+   m_gf_extension_degree = extdeg;
+   m_gf_cardinality = 1 << extdeg;
+   m_gf_multiplicative_order = m_gf_cardinality - 1;
+
+   init_exp();
+   init_log();
+   }
+
+gf2m Gf2m_Field::gf_div(gf2m x, gf2m y)
+   {
+   s32bit sub_res = ((s32bit)m_gf_log_table[x]) - ((s32bit) m_gf_log_table[y]);
+   s32bit modq_res = ((s32bit)_gf_modq_1(sub_res));
+   s32bit div_res = (((s32bit)x) ? ((s32bit) m_gf_exp_table[modq_res]) : 0  );
+   return (gf2m) div_res;
+   }
+
+// we suppose i >= 0. Par convention 0^0 = 1
+gf2m Gf2m_Field::gf_pow(gf2m x, int i)
+   {
+   if (i == 0)
+      return 1;
+   else if (x == 0)
+      return 0;
+   else
+      {
+      // i mod (q-1)
+      while (i >> get_extension_degree())
+         i = (i & (gf_ord())) + (i >> get_extension_degree());
+      i *= m_gf_log_table[x];
+      while (i >> get_extension_degree())
+         i = (i & (gf_ord())) + (i >> get_extension_degree());
+      return m_gf_exp_table[i];
+      }
+   }
+
+}
+
+}
diff --git a/src/lib/pubkey/mce/gf2m_small_m.h b/src/lib/pubkey/mce/gf2m_small_m.h
new file mode 100644
index 000000000..9fc42f1fc
--- /dev/null
+++ b/src/lib/pubkey/mce/gf2m_small_m.h
@@ -0,0 +1,236 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef BOTAN_GF2M_SMALL_M_H__
+#define BOTAN_GF2M_SMALL_M_H__
+
+#include <vector>
+#include <botan/types.h>
+
+namespace Botan {
+
+namespace gf2m_small_m {
+
+typedef u16bit gf2m;
+
+class Gf2m_Field
+   {
+   public:
+      Gf2m_Field(size_t extdeg);
+
+      gf2m gf_mul(gf2m x, gf2m y)
+         {
+         return ((x) ? gf_mul_fast(x, y) : 0);
+         }
+
+      gf2m gf_square(gf2m x)
+         {
+         return ((x) ? m_gf_exp_table[_gf_modq_1(m_gf_log_table[x] << 1)] : 0);
+         }
+
+      gf2m square_rr(gf2m x)
+         {
+         return _gf_modq_1(x << 1);
+         }
+
+      // naming convention of GF(2^m) field operations:
+      // l logarithmic, unreduced
+      // r logarithmic, reduced
+      // n normal, non-zero
+      // z normal, might be zero
+      //
+      inline gf2m gf_mul_lll(gf2m a, gf2m b);
+      inline gf2m gf_mul_rrr(gf2m a, gf2m b);
+      inline gf2m gf_mul_nrr(gf2m a, gf2m b);
+      inline gf2m gf_mul_rrn(gf2m a, gf2m y);
+      inline gf2m gf_mul_lnn(gf2m x, gf2m y);
+      inline gf2m gf_mul_rnn(gf2m x, gf2m y);
+      inline gf2m gf_mul_nrn(gf2m a, gf2m y);
+      inline gf2m gf_mul_rnr(gf2m y, gf2m a);
+      inline gf2m gf_mul_zrz(gf2m a, gf2m y);
+      inline gf2m gf_mul_zzr(gf2m a, gf2m y);
+      inline gf2m gf_mul_nnr(gf2m y, gf2m a);
+      inline gf2m gf_sqrt(gf2m x) ;
+      gf2m gf_div(gf2m x, gf2m y);
+      inline gf2m gf_div_rnn(gf2m x, gf2m y);
+      inline gf2m gf_div_rnr(gf2m x, gf2m b);
+      inline gf2m gf_div_nrr(gf2m a, gf2m b);
+      inline gf2m gf_div_zzr(gf2m x, gf2m b);
+      inline gf2m gf_inv(gf2m x);
+      inline gf2m gf_inv_rn(gf2m x);
+      inline gf2m gf_square_ln(gf2m x);
+      inline gf2m gf_square_rr(gf2m a) ;
+      inline gf2m gf_l_from_n(gf2m x);
+
+      inline gf2m gf_mul_fast(gf2m a, gf2m b);
+
+      gf2m gf_exp(gf2m i)
+         {
+         return m_gf_exp_table[i]; /* alpha^i */
+         }
+
+      gf2m gf_log(gf2m i)
+         {
+         return m_gf_log_table[i]; /* return i when x=alpha^i */
+         }
+
+      inline gf2m gf_ord() const
+         {
+         return m_gf_multiplicative_order;
+         }
+
+      inline gf2m get_extension_degree() const
+         {
+         return m_gf_extension_degree;
+         }
+
+      inline gf2m get_cardinality() const
+         {
+         return m_gf_cardinality;
+         }
+
+      gf2m gf_pow(gf2m x, int i) ;
+
+   private:
+      gf2m m_gf_extension_degree, m_gf_cardinality, m_gf_multiplicative_order;
+      std::vector<gf2m> m_gf_log_table;
+      std::vector<gf2m> m_gf_exp_table;
+
+      inline gf2m _gf_modq_1(s32bit d);
+      void init_log();
+      void init_exp();
+   };
+
+gf2m Gf2m_Field::_gf_modq_1(s32bit d)
+   {
+   return  (((d) & gf_ord()) + ((d) >> m_gf_extension_degree));
+   }
+
+gf2m Gf2m_Field::gf_mul_fast(gf2m x, gf2m y)
+   {
+   return ((y) ? m_gf_exp_table[_gf_modq_1(m_gf_log_table[x] + m_gf_log_table[y])] : 0);
+   }
+
+gf2m Gf2m_Field::gf_mul_lll(gf2m a, gf2m b)
+   {
+   return  (a + b);
+   }
+
+gf2m Gf2m_Field::gf_mul_rrr(gf2m a, gf2m b)
+   {
+   return (_gf_modq_1(gf_mul_lll(a, b)));
+   }
+
+gf2m Gf2m_Field::gf_mul_nrr(gf2m a, gf2m b)
+   {
+   return (gf_exp(gf_mul_rrr(a, b)));
+   }
+
+gf2m Gf2m_Field::gf_mul_rrn(gf2m a, gf2m y)
+   {
+   return _gf_modq_1(gf_mul_lll(a, gf_log(y)));
+   }
+
+gf2m Gf2m_Field::gf_mul_rnr(gf2m y, gf2m a)
+   {
+   return gf_mul_rrn(a, y);
+   }
+
+gf2m Gf2m_Field::gf_mul_lnn(gf2m x, gf2m y)
+   {
+   return (m_gf_log_table[x] + m_gf_log_table[y]);
+   }
+gf2m Gf2m_Field::gf_mul_rnn(gf2m x, gf2m y)
+   {
+   return _gf_modq_1(gf_mul_lnn(x, y));
+   }
+
+gf2m Gf2m_Field::gf_mul_nrn(gf2m a, gf2m y)
+   {
+   return m_gf_exp_table[_gf_modq_1((a) + m_gf_log_table[y])];
+   }
+
+/**
+* zero operand allowed
+*/
+gf2m Gf2m_Field::gf_mul_zrz(gf2m a, gf2m y)
+   {
+   return ( (y == 0) ? 0 : gf_mul_nrn(a, y) );
+   }
+
+gf2m Gf2m_Field::gf_mul_zzr(gf2m a, gf2m y)
+   {
+   return gf_mul_zrz(y, a);
+   }
+/**
+* non-zero operand
+*/
+gf2m Gf2m_Field::gf_mul_nnr(gf2m y, gf2m a)
+   {
+   return gf_mul_nrn( a, y);
+   }
+
+gf2m Gf2m_Field::gf_sqrt(gf2m x)
+   {
+   return ((x) ? m_gf_exp_table[_gf_modq_1(m_gf_log_table[x] << (m_gf_extension_degree-1))] : 0);
+   }
+
+gf2m Gf2m_Field::gf_div_rnn(gf2m x, gf2m y)
+   {
+   return _gf_modq_1(m_gf_log_table[x] - m_gf_log_table[y]);
+   }
+gf2m Gf2m_Field::gf_div_rnr(gf2m x, gf2m b)
+   {
+   return _gf_modq_1(m_gf_log_table[x] - b);
+   }
+gf2m Gf2m_Field::gf_div_nrr(gf2m a, gf2m b)
+   {
+   return m_gf_exp_table[_gf_modq_1(a - b)];
+   }
+
+gf2m Gf2m_Field::gf_div_zzr(gf2m x, gf2m b)
+   {
+   return ((x) ? m_gf_exp_table[_gf_modq_1(m_gf_log_table[x] - b)] : 0);
+   }
+
+gf2m Gf2m_Field::gf_inv(gf2m x)
+   {
+   return m_gf_exp_table[gf_ord() - m_gf_log_table[x]];
+   }
+gf2m Gf2m_Field::gf_inv_rn(gf2m x)
+   {
+   return (gf_ord() - m_gf_log_table[x]);
+   }
+
+gf2m Gf2m_Field::gf_square_ln(gf2m x)
+   {
+   return m_gf_log_table[x] << 1;
+   }
+
+gf2m Gf2m_Field::gf_square_rr(gf2m a)
+   {
+   return a << 1;
+   }
+
+gf2m Gf2m_Field::gf_l_from_n(gf2m x)
+   {
+   return m_gf_log_table[x];
+   }
+
+u32bit encode_gf2m(gf2m to_enc, byte* mem);
+
+gf2m decode_gf2m(const byte* mem);
+
+}
+
+}
+
+#endif
diff --git a/src/lib/pubkey/mce/goppa_code.cpp b/src/lib/pubkey/mce/goppa_code.cpp
new file mode 100644
index 000000000..bc8fb7b32
--- /dev/null
+++ b/src/lib/pubkey/mce/goppa_code.cpp
@@ -0,0 +1,203 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/goppa_code.h>
+
+#include <botan/gf2m_rootfind_dcmp.h>
+#include <botan/code_based_util.h>
+
+namespace Botan {
+
+namespace {
+
+void matrix_arr_mul(std::vector<u32bit> matrix,
+                    u32bit numo_rows,
+                    u32bit words_per_row,
+                    const byte* input_vec,
+                    u32bit* output_vec, u32bit output_vec_len)
+   {
+   for(size_t j = 0; j < numo_rows; j++)
+      {
+      if((input_vec[j / 8] >> (j % 8)) & 1)
+         {
+         for(size_t i = 0; i < output_vec_len; i ++)
+            {
+            output_vec[i] ^= matrix[ j * (words_per_row) + i];
+            }
+         }
+      }
+   }
+
+/**
+* returns the error vector to the syndrome
+*/
+secure_vector<gf2m> goppa_decode(const polyn_gf2m & syndrom_polyn,
+                                 const polyn_gf2m & g,
+                                 const std::vector<polyn_gf2m> & sqrtmod,
+                                 const std::vector<gf2m> & Linv)
+   {
+   gf2m a;
+   u32bit code_length = Linv.size();
+   u32bit t = g.get_degree();
+
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field = g.get_sp_field();
+
+   std::pair<polyn_gf2m, polyn_gf2m> h__aux = polyn_gf2m::eea_with_coefficients( syndrom_polyn, g, 1);
+   polyn_gf2m & h = h__aux.first;
+   polyn_gf2m & aux = h__aux.second;
+   a = sp_field->gf_inv(aux.get_coef(0));
+   gf2m log_a = sp_field->gf_log(a);
+   for(int i = 0; i <= h.get_degree(); ++i)
+      {
+      h.set_coef(i,sp_field->gf_mul_zrz(log_a,h.get_coef(i)));
+      }
+
+   //  compute h(z) += z
+   h.add_to_coef( 1, 1);
+   // compute S square root of h (using sqrtmod)
+   polyn_gf2m S(t - 1, g.get_sp_field());
+
+   for(u32bit i=0;i<t;i++)
+      {
+      a = sp_field->gf_sqrt(h.get_coef(i));
+
+      if(i & 1)
+         {
+         for(u32bit j=0;j<t;j++)
+            {
+            S.add_to_coef( j, sp_field->gf_mul(a, sqrtmod[i/2].get_coef(j)));
+            }
+         }
+      else
+         {
+         S.add_to_coef( i/2, a);
+         }
+      } /* end for loop (i) */
+
+
+   S.get_degree();
+
+   std::pair<polyn_gf2m, polyn_gf2m> v__u = polyn_gf2m::eea_with_coefficients(S, g, t/2+1);
+   polyn_gf2m & u = v__u.second;
+   polyn_gf2m & v = v__u.first;
+
+   // sigma = u^2+z*v^2
+   polyn_gf2m sigma ( t , g.get_sp_field());
+
+   const size_t u_deg = u.get_degree();
+   for(size_t i = 0; i <= u_deg; ++i)
+      {
+      sigma.set_coef(2*i, sp_field->gf_square(u.get_coef(i)));
+      }
+
+   const size_t v_deg = v.get_degree();
+   for(size_t i = 0; i <= v_deg; ++i)
+      {
+      sigma.set_coef(2*i+1, sp_field->gf_square(v.get_coef(i)));
+      }
+
+   secure_vector<gf2m> res = find_roots_gf2m_decomp(sigma, code_length);
+   size_t d = res.size();
+
+   secure_vector<gf2m> result(d);
+   for(u32bit i = 0; i < d; ++i)
+      {
+      gf2m current = res[i];
+
+      gf2m tmp;
+      tmp = gray_to_lex(current);
+      if(tmp >= code_length) /* invalid root */
+         {
+         result[i] = i;
+         }
+      result[i] = Linv[tmp];
+      }
+
+   return result;
+   }
+}
+
+
+/**
+* @param p_err_pos_len must point to the available length of err_pos on input, the
+* function will set it to the actual number of errors returned in the err_pos
+* array */
+secure_vector<byte> mceliece_decrypt(
+   secure_vector<gf2m> & error_pos,
+   const byte *ciphertext, u32bit ciphertext_len,
+   const McEliece_PrivateKey & key)
+   {
+
+   u32bit dimension = key.get_dimension();
+   u32bit codimension = key.get_codimension();
+   u32bit t = key.get_goppa_polyn().get_degree();
+   polyn_gf2m syndrome_polyn(key.get_goppa_polyn().get_sp_field()); // init as zero polyn
+   const unsigned unused_pt_bits = dimension % 8;
+   const unsigned char unused_pt_bits_mask = (1 << unused_pt_bits) - 1;
+
+   if(ciphertext_len != (key.get_code_length()+7)/8)
+      {
+      throw Invalid_Argument("wrong size of McEliece ciphertext");
+      }
+   u32bit cleartext_len = (key.get_message_word_bit_length()+7)/8;
+
+   if(cleartext_len != bit_size_to_byte_size(dimension))
+      {
+      throw Invalid_Argument("mce-decryption: wrong length of cleartext buffer");
+      }
+
+
+   secure_vector<u32bit> syndrome_vec(bit_size_to_32bit_size(codimension));
+   matrix_arr_mul(
+      key.get_H_coeffs(),
+      key.get_code_length(),
+      bit_size_to_32bit_size(codimension),
+      ciphertext,
+      &syndrome_vec[0], syndrome_vec.size() );
+   secure_vector<byte> syndrome_byte_vec(bit_size_to_byte_size(codimension));
+   u32bit syndrome_byte_vec_size = syndrome_byte_vec.size();
+   for(u32bit i = 0; i < syndrome_byte_vec_size; i++)
+      {
+      syndrome_byte_vec[i] = syndrome_vec[i/4] >> (8* (i % 4));
+      }
+   syndrome_polyn = polyn_gf2m( t-1, &syndrome_byte_vec[0],bit_size_to_byte_size(codimension), key.get_goppa_polyn().get_sp_field());
+
+
+
+   syndrome_polyn.get_degree();
+   error_pos = goppa_decode(syndrome_polyn, key.get_goppa_polyn(), key.get_sqrtmod(), key.get_Linv());
+
+   u32bit nb_err = error_pos.size();
+
+
+   secure_vector<byte> cleartext(cleartext_len);
+   std::memcpy(&cleartext[0], ciphertext, cleartext_len);
+
+   for(u32bit i = 0; i < nb_err; i++)
+      {
+      gf2m current = error_pos[i];
+
+      if(current >= cleartext_len * 8)
+         {
+         // an invalid position, this shouldn't happen
+         continue;
+         }
+      cleartext[current / 8] ^= (1 << (current % 8));
+      }
+   if(unused_pt_bits)
+      {
+      cleartext[cleartext_len - 1] &= unused_pt_bits_mask;
+      }
+
+   return cleartext;
+   }
+
+}
diff --git a/src/lib/pubkey/mce/goppa_code.h b/src/lib/pubkey/mce/goppa_code.h
new file mode 100644
index 000000000..1c8a163b7
--- /dev/null
+++ b/src/lib/pubkey/mce/goppa_code.h
@@ -0,0 +1,32 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef __goppa_code__H_
+#define __goppa_code__H_
+
+#include <botan/polyn_gf2m.h>
+#include <botan/mceliece_key.h>
+
+
+
+namespace Botan
+{
+
+  std::vector<byte> mceliece_encrypt( const secure_vector<byte> & cleartext, std::vector<byte> const& public_matrix, const secure_vector<gf2m> & err_pos, u32bit code_length);
+
+
+secure_vector<byte> mceliece_decrypt(
+    secure_vector<gf2m> & error_pos,
+    const byte *ciphertext, u32bit ciphertext_len,
+    const McEliece_PrivateKey & key);
+} //end namepace Botan
+
+#endif /* h-guard */
diff --git a/src/lib/pubkey/mce/info.txt b/src/lib/pubkey/mce/info.txt
new file mode 100644
index 000000000..6c0da2199
--- /dev/null
+++ b/src/lib/pubkey/mce/info.txt
@@ -0,0 +1,17 @@
+define MCELIECE 20141124
+
+<header:public>
+code_based_util.h
+gf2m_rootfind_dcmp.h
+gf2m_small_m.h
+goppa_code.h
+mce_overbeck_cca2.h
+mceliece.h
+mceliece_key.h
+polyn_gf2m.h
+</header:public>
+
+<header:internal>
+binary_matrix.h
+code_based_key_gen.h
+</header:internal>
diff --git a/src/lib/pubkey/mce/mce_overbeck_cca2.cpp b/src/lib/pubkey/mce/mce_overbeck_cca2.cpp
new file mode 100644
index 000000000..7edd2e2a3
--- /dev/null
+++ b/src/lib/pubkey/mce/mce_overbeck_cca2.cpp
@@ -0,0 +1,182 @@
+/**
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/mce_overbeck_cca2.h>
+#include <botan/sha2_64.h>
+#include <botan/mceliece.h>
+#include <botan/internal/xor_buf.h>
+
+namespace Botan
+{
+
+  McEliece_Overbeck_CCA2_Public_Operation::McEliece_Overbeck_CCA2_Public_Operation(const McEliece_PublicKey& public_key)
+    :m_raw_pub_op(public_key, public_key.get_code_length())
+  {
+    if(public_key.get_message_word_bit_length() < 1024)
+    {
+      // k is smaller than the minimum required length for Overbeck conversion
+      // using SHA-512
+      throw Invalid_Argument("McEliece parameters are too small to support the Overbeck conversion with SHA-512, the dimension of the code must be at least 1024");
+    }
+  }
+
+
+  secure_vector<byte> McEliece_Overbeck_CCA2_Public_Operation::encrypt(const byte msg[], size_t msg_len, RandomNumberGenerator& rng)
+  {
+    const u32bit k = m_raw_pub_op.max_input_bits();
+    const u32bit l = 512; // output length of SHA-512
+    const u32bit l_bytes = l/8;
+    const u32bit u1_bit_length = k - l;
+    const u32bit u1_length_ceil = (u1_bit_length + 7)/8; // valid lengths ensured already during construction
+    const McEliece_PublicKey & key = m_raw_pub_op.get_key();
+    const u32bit n = key.get_code_length();
+    const u32bit n_bytes_ceil = (n+7)/8;
+    const u32bit k_bytes_ceil = (k+7)/8;
+
+    if(msg_len != l_bytes)
+    {
+      throw Invalid_Argument("McEliece/Overbeck message lengtth must be 64 bytes");
+    }
+     secure_vector<byte> u1(u1_length_ceil);
+    rng.randomize(&u1[0], u1.size());
+   // unused bits of final byte of u1 must be set to zero
+    u32bit used = u1_bit_length % 8;
+   if(used)
+   {
+      byte mask = (1 << used) - 1;
+
+      u1[u1.size() - 1] &= mask;
+   }
+
+    secure_vector<byte> u2(l_bytes);
+    rng.randomize(&u2[0], u2.size());
+
+    // compute the hash of m||u1:
+    SHA_512 hash;
+
+    hash.update(msg, msg_len);
+    hash.update(u2);
+    secure_vector<byte> hash_m_u2 = hash.final();
+
+    //std::cout << "enc hash_m_u2 " << hex_encode(hash_m_u2) << "\n";
+
+    secure_vector<byte> mce_msg(k_bytes_ceil);
+    std::memcpy(&mce_msg[0], &hash_m_u2[0], hash_m_u2.size());
+    std::memcpy(&mce_msg[hash_m_u2.size()], &u1[0], u1.size());
+
+// create the error vector
+    secure_vector<gf2m> err_pos = create_random_error_positions(n, key.get_t(), rng);
+
+    secure_vector<byte> err_vec = mceliece_message_parts::error_vector_from_error_positions(&err_pos[0], err_pos.size(), n);
+
+    mceliece_message_parts parts(err_pos, mce_msg, n);
+
+    secure_vector<byte> message_and_error_input = parts.get_concat();
+
+    //std::cout << "enc msg_and_error " << hex_encode(message_and_error_input) << "\n";
+    //std::cout << "enc h(msg_and_error) " << hex_encode(hash.process(message_and_error_input)) << "\n";
+
+    secure_vector<byte> mce_ct = m_raw_pub_op.encrypt(&message_and_error_input[0], message_and_error_input.size(), rng);
+
+    secure_vector<byte> result(n_bytes_ceil + 2*l_bytes);
+
+    BOTAN_ASSERT(mce_ct.size() == (key.get_code_length()+7)/8, "Expected size");
+
+    std::memcpy(&result[0], &mce_ct[0], mce_ct.size());
+
+
+    // z2 part of the ciphertext
+    SHA_512 hash2;
+    secure_vector<byte> hash_u1 = hash2.process(u1);
+
+    //std::cout << "enc hash_u1 " << hex_encode(hash_u1) << "\n";
+
+    xor_buf(&result[mce_ct.size()], &hash_u1[0], &msg[0], l_bytes);
+
+    // 3rd part of the overbeck ct
+    SHA_512 hash3;
+    secure_vector<byte> err_hash = hash3.process(err_vec);
+
+    //std::cout << "enc err_hash " << hex_encode(err_hash) << "\n";
+
+    const u32bit z3_offs = n_bytes_ceil + l_bytes;
+    xor_buf(&result[z3_offs], &u2[0], &err_hash[0], l_bytes);
+
+    return result;
+  }
+
+  McEliece_Overbeck_CCA2_Private_Operation::McEliece_Overbeck_CCA2_Private_Operation(const McEliece_PrivateKey& mce_key)
+    :m_raw_priv_op(mce_key)
+  {
+    if(mce_key.get_dimension() < 1024)
+    {
+      // k is smaller than the minimum required length for Overbeck conversion
+      // using SHA-512
+      throw Invalid_Argument("McEliece parameters are too small to support the Overbeck conversion with SHA-512, the dimension of the code must be at least 1024");
+    }
+  }
+
+  secure_vector<byte> McEliece_Overbeck_CCA2_Private_Operation::decrypt(const byte msg[], size_t msg_len)
+  {
+
+    const McEliece_PrivateKey& key = m_raw_priv_op.get_key();
+    const u32bit k = key.get_dimension();
+    const u32bit l = 512; // output length of SHA-512
+    const u32bit l_bytes = l/8;
+    const u32bit r_length_ceil = (k - l + 7)/8; // valid lengths ensured already during construction
+    const u32bit n = key.get_code_length();
+    const u32bit n_bytes_ceil = (n+7)/8;
+
+    const u32bit z2_offs = n_bytes_ceil;
+    const u32bit z3_offs = n_bytes_ceil + l_bytes;
+
+    if(msg_len != (max_input_bits()+7)/8)
+    {
+      throw Invalid_Argument("wrong length of McEliece/Overbeck ciphertext");
+    }
+    secure_vector<byte> mce_pt_and_err  = m_raw_priv_op.decrypt(msg, n_bytes_ceil);
+
+    SHA_512 hash;
+    //std::cout << "dec msg_and_error " << hex_encode(mce_pt_and_err) << "\n";
+    //std::cout << "dec h(msg_and_error) " << hex_encode(hash.process(mce_pt_and_err)) << "\n";
+
+    mceliece_message_parts parts(&mce_pt_and_err[0], mce_pt_and_err.size(), n);
+
+    secure_vector<byte> mce_pt = parts.get_message_word();
+    secure_vector<byte> err_vec = parts.get_error_vector();
+
+    secure_vector<byte> h(l_bytes);
+    std::memcpy(&h[0], &mce_pt[0], l_bytes);
+    secure_vector<byte> r(r_length_ceil);
+    std::memcpy(&r[0], &mce_pt[l_bytes], r.size());
+
+    secure_vector<byte> hash_r = hash.process(r);
+    //std::cout << "dec hash_r " << hex_encode(hash_r) << "\n";
+
+    secure_vector<byte> m(l_bytes);
+    xor_buf(&m[0], &msg[z2_offs], &hash_r[0], l_bytes);
+
+    SHA_512 hash2;
+    secure_vector<byte> hash_e = hash2.process(err_vec);
+    //std::cout << "dec hash_e " << hex_encode(hash_e) << "\n";
+    xor_buf(&hash_e[0], &msg[z3_offs], l_bytes);
+    // hash_e now is H(e) ^ z3 = u2
+
+    SHA_512 hash3;
+    hash3.update(m);
+    hash3.update(hash_e);
+    secure_vector<byte> h_cmp = hash3.final();
+
+    //std::cout << "dec hash_cmp " << hex_encode(h_cmp) << "\n";
+    if(h_cmp != h)
+       throw Integrity_Failure("McEliece/Overbeck CCA2 check failed");
+    return m;
+
+  }
+
+}
diff --git a/src/lib/pubkey/mce/mce_overbeck_cca2.h b/src/lib/pubkey/mce/mce_overbeck_cca2.h
new file mode 100644
index 000000000..1b9439753
--- /dev/null
+++ b/src/lib/pubkey/mce/mce_overbeck_cca2.h
@@ -0,0 +1,46 @@
+/**
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef __mce_overbeck_cca2__H_
+#define __mce_overbeck_cca2__H_
+
+#include <botan/pk_ops.h>
+#include <botan/mceliece_key.h>
+#include <botan/mceliece.h>
+#include <botan/types.h>
+#include <botan/secmem.h>
+
+namespace Botan
+{
+class BOTAN_DLL McEliece_Overbeck_CCA2_Public_Operation : public PK_Ops::Encryption
+{
+  public:
+  McEliece_Overbeck_CCA2_Public_Operation(const McEliece_PublicKey& public_key);
+
+  size_t max_input_bits() const { return 512; };
+  secure_vector<byte> encrypt(const byte msg[], size_t msg_len, RandomNumberGenerator&);
+
+  private:
+  McEliece_Public_Operation m_raw_pub_op;
+};
+
+class BOTAN_DLL McEliece_Overbeck_CCA2_Private_Operation : public PK_Ops::Decryption
+  {
+    public:
+      McEliece_Overbeck_CCA2_Private_Operation(const McEliece_PrivateKey& mce_key);
+
+      size_t max_input_bits() const { return (m_raw_priv_op.max_input_bits()+7)/8*8 + 2*512; };
+
+
+secure_vector<byte> decrypt(const byte msg[], size_t msg_len);
+    private:
+      McEliece_Private_Operation m_raw_priv_op;
+  };
+}
+
+#endif /* h-guard */
diff --git a/src/lib/pubkey/mce/mceliece.cpp b/src/lib/pubkey/mce/mceliece.cpp
new file mode 100644
index 000000000..69d5eb008
--- /dev/null
+++ b/src/lib/pubkey/mce/mceliece.cpp
@@ -0,0 +1,172 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/mceliece.h>
+#include <botan/mceliece_key.h>
+#include <botan/internal/code_based_key_gen.h>
+#include <botan/polyn_gf2m.h>
+#include <botan/code_based_util.h>
+#include <botan/goppa_code.h>
+#include <botan/internal/bit_ops.h>
+#include <botan/internal/xor_buf.h>
+#include <iostream>
+
+namespace Botan {
+
+namespace {
+void concat_vectors(unsigned char* x, const unsigned char* a, const unsigned char* b, u32bit dimension, u32bit codimension)
+   {
+   if(dimension % 8 == 0)
+      {
+      std::memcpy(x, a, bit_size_to_byte_size(dimension));
+      std::memcpy(((unsigned char *) x) + bit_size_to_byte_size(dimension), b, bit_size_to_byte_size(codimension));
+      }
+   else
+      {
+      u32bit i, j, k, l;
+      i = dimension - 8 * (dimension/ 8);
+      j = 8 - i;
+      l = dimension / 8;
+      std::memcpy(x, a, 1 * (dimension / 8));
+      x[l] = ((byte) (a[l] & ((1 << i) - 1)));
+
+      for(k = 0; k < codimension / 8; ++k)
+         {
+         x[l] ^= ((byte) ( b[k] << i));
+         ++l;
+         x[l] = ((byte) (b[k] >> j));
+         }
+      x[l] ^= ((byte) ( b[k] << i));
+      }
+   }
+
+std::vector<byte> mult_by_pubkey(const byte *cleartext,
+                                 std::vector<byte> const& public_matrix,
+                                 u32bit code_length, u32bit t)
+   {
+   std::vector<byte> ciphertext(code_length);
+   u32bit i, j;
+   u32bit ext_deg = ceil_log2(code_length);
+   u32bit codimension = ext_deg * t;
+   u32bit dimension = code_length - codimension;
+   std::vector<byte> cR(bit_size_to_32bit_size(codimension)* sizeof(u32bit));
+
+   const byte* pt = &public_matrix[0];
+
+   for(i = 0; i < dimension / 8; ++i)
+      {
+      for(j = 0; j < 8; ++j)
+         {
+         if(cleartext[i] & (1 << j))
+            {
+            xor_buf(&cR[0], pt, cR.size());
+            }
+         pt += bit_size_to_32bit_size(codimension) * sizeof(u32bit);
+         }
+      }
+
+   for(j = 0; j < dimension % 8 ; ++j)
+      {
+      if(cleartext[i] & (1 << j))
+         {
+         xor_buf(&cR[0], pt, bit_size_to_byte_size(codimension));
+         }
+      pt += bit_size_to_32bit_size(codimension) * sizeof(u32bit);
+      }
+
+   concat_vectors( &ciphertext[0],  cleartext, &cR[0], dimension, codimension);
+   return ciphertext;
+   }
+
+}
+
+secure_vector<gf2m> create_random_error_positions(unsigned code_length,
+                                                  unsigned error_weight,
+                                                  RandomNumberGenerator& rng)
+   {
+   secure_vector<gf2m> result(error_weight);
+   gf2m i;
+   for(i = 0; i < result.size(); i++)
+      {
+      unsigned j;
+      char try_again = 0;
+      do
+         {
+         try_again = 0;
+         gf2m new_pos = random_code_element(code_length, rng);
+         for(j = 0; j < i; j++)
+            {
+            if(new_pos == result[j])
+               {
+               try_again = 1;
+               break;
+               }
+            }
+         result[i] = new_pos;
+         } while(try_again);
+      }
+   return result;
+   }
+
+McEliece_Private_Operation::McEliece_Private_Operation(const McEliece_PrivateKey& private_key)
+   :m_priv_key(private_key)
+   {
+   }
+
+secure_vector<byte> McEliece_Private_Operation::decrypt(const byte msg[], size_t msg_len)
+   {
+   secure_vector<gf2m> err_pos;
+
+   secure_vector<byte> plaintext = mceliece_decrypt(
+      err_pos,
+      msg, msg_len,
+      m_priv_key
+      );
+
+   return mceliece_message_parts(err_pos, plaintext, m_priv_key.get_code_length()).get_concat();
+   }
+
+McEliece_Public_Operation::McEliece_Public_Operation(const McEliece_PublicKey& public_key, u32bit the_code_length)
+   :m_pub_key(public_key),
+    m_code_length(the_code_length)
+   {}
+
+secure_vector<byte> McEliece_Public_Operation::encrypt(const byte msg[], size_t msg_len, RandomNumberGenerator&)
+   {
+   mceliece_message_parts parts(msg, msg_len, m_pub_key.get_code_length());
+   secure_vector<gf2m> err_pos = parts.get_error_positions();
+   secure_vector<byte> message_word = parts.get_message_word();
+   secure_vector<byte> ciphertext((m_pub_key.get_code_length()+7)/8);
+
+
+   std::vector<byte> ciphertext_tmp = mceliece_encrypt( message_word,  m_pub_key.get_public_matrix(), err_pos, m_code_length);
+
+   std::memcpy(&ciphertext[0], &ciphertext_tmp[0], ciphertext.size());
+   return ciphertext;
+   }
+
+std::vector<byte> mceliece_encrypt(const secure_vector<byte> & cleartext,
+                                   std::vector<byte> const& public_matrix,
+                                   const secure_vector<gf2m> & err_pos,
+                                   u32bit code_length)
+   {
+   std::vector<byte> ciphertext = mult_by_pubkey(&cleartext[0], public_matrix, code_length, err_pos.size());
+
+   // flip t error positions
+   for(size_t i = 0; i < err_pos.size(); ++i)
+      {
+      ciphertext[err_pos[i] / 8] ^= (1 << (err_pos[i] % 8));
+      }
+
+   return ciphertext;
+   }
+
+}
diff --git a/src/lib/pubkey/mce/mceliece.h b/src/lib/pubkey/mce/mceliece.h
new file mode 100644
index 000000000..8705d8132
--- /dev/null
+++ b/src/lib/pubkey/mce/mceliece.h
@@ -0,0 +1,149 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+
+
+#include <botan/secmem.h>
+#include <botan/types.h>
+#include <botan/pk_ops.h>
+#include <botan/mceliece_key.h>
+
+#ifndef BOTAN_MCELIECE_H__
+#define BOTAN_MCELIECE_H__
+
+#define MASK_LOG2_BYTE ((1 << 3) - 1)
+#define _BITP_TO_BYTEP(__bit_pos) (__bit_pos >> 3)
+#define _BITP_TO_BYTEOFFS(__bit_pos) (__bit_pos & MASK_LOG2_BYTE)
+
+namespace Botan {
+
+
+
+  secure_vector<gf2m> create_random_error_positions(unsigned code_length, unsigned error_weight, RandomNumberGenerator& rng);
+
+
+class mceliece_message_parts
+{
+  public:
+
+    mceliece_message_parts(const secure_vector<gf2m>& err_pos, const byte* message, u32bit message_length, u32bit code_length)
+    :m_error_vector(error_vector_from_error_positions(&err_pos[0], err_pos.size(), code_length)),
+    m_code_length(code_length)
+  {
+    m_message_word.resize(message_length);
+    std::memcpy(&m_message_word[0], message, message_length);
+  };
+
+    mceliece_message_parts(const secure_vector<gf2m>& err_pos, const secure_vector<byte>& message, unsigned code_length)
+    :m_error_vector(error_vector_from_error_positions(&err_pos[0], err_pos.size(), code_length)),
+    m_message_word(message),
+    m_code_length(code_length)
+  {};
+    static secure_vector<byte> error_vector_from_error_positions(const gf2m* err_pos, size_t err_pos_len, size_t code_length)
+    {
+      secure_vector<byte> result((code_length+7)/8);
+      for(unsigned i = 0; i < err_pos_len; i++)
+      {
+        u16bit pos = err_pos[i];
+        u32bit byte_pos = _BITP_TO_BYTEP(pos);
+        if(byte_pos > result.size())
+        {
+          throw Invalid_Argument("error position larger than code size");
+        }
+       result[byte_pos] |= (1 << _BITP_TO_BYTEOFFS(pos));
+      }
+      return result;
+    };
+    mceliece_message_parts(const byte* message_concat_errors, size_t message_concat_errors_len, unsigned code_length)
+      :m_code_length(code_length)
+    {
+      size_t err_vec_len = (code_length+7)/8;
+      if(message_concat_errors_len < err_vec_len )
+      {
+        throw Invalid_Argument("cannot split McEliece message parts");
+      }
+      size_t err_vec_start_pos = message_concat_errors_len - err_vec_len;
+      m_message_word = secure_vector<byte>(err_vec_start_pos );
+      std::memcpy(&m_message_word[0], &message_concat_errors[0], err_vec_start_pos);
+      m_error_vector = secure_vector<byte>(err_vec_len );
+      std::memcpy(&m_error_vector[0],  &message_concat_errors[err_vec_start_pos], err_vec_len);
+
+    };
+    secure_vector<byte> get_concat() const
+    {
+      secure_vector<byte> result(m_error_vector.size() + m_message_word.size());
+      std::memcpy(&result[0], &m_message_word[0], m_message_word.size());
+      std::memcpy(&result[m_message_word.size()], &m_error_vector[0], m_error_vector.size());
+      return result;
+    };
+    secure_vector<gf2m> get_error_positions() const
+    {
+      secure_vector<gf2m> result;
+      for(unsigned i = 0; i < m_code_length; i++)
+      {
+
+        if ( i >= m_code_length)
+        {
+          throw Invalid_Argument("index out of range in get_error_positions()");
+        }
+        if((m_error_vector[_BITP_TO_BYTEP(i)] >> _BITP_TO_BYTEOFFS(i)) & 1)
+        {
+          result.push_back(i);
+        }
+      }
+      return result;
+    };
+    secure_vector<byte> get_error_vector() const { return m_error_vector; };
+    secure_vector<byte> get_message_word() const { return m_message_word; };
+  private:
+    secure_vector<byte> m_error_vector;
+    secure_vector<byte> m_message_word;
+    unsigned m_code_length;
+};
+
+class BOTAN_DLL McEliece_Private_Operation : public PK_Ops::Decryption
+  {
+    public:
+      McEliece_Private_Operation(const McEliece_PrivateKey& mce_key);
+
+      size_t max_input_bits() const {
+        return m_priv_key.max_input_bits();
+
+      };
+
+
+secure_vector<byte> decrypt(const byte msg[], size_t msg_len);
+
+  McEliece_PrivateKey const& get_key() const { return m_priv_key; };
+
+    private:
+      const McEliece_PrivateKey m_priv_key;
+  };
+
+class BOTAN_DLL McEliece_Public_Operation : public PK_Ops::Encryption
+{
+  public:
+  McEliece_Public_Operation(const McEliece_PublicKey& public_key, u32bit code_length);
+
+  size_t max_input_bits() const { return m_pub_key.max_input_bits(); };
+  secure_vector<byte> encrypt(const byte msg[], size_t msg_len, RandomNumberGenerator&);
+
+  McEliece_PublicKey const& get_key() const { return m_pub_key; };
+
+  private:
+  McEliece_PublicKey m_pub_key;
+  u32bit m_code_length;
+};
+
+}
+
+
+#endif
diff --git a/src/lib/pubkey/mce/mceliece_key.cpp b/src/lib/pubkey/mce/mceliece_key.cpp
new file mode 100644
index 000000000..71e0ec4e8
--- /dev/null
+++ b/src/lib/pubkey/mce/mceliece_key.cpp
@@ -0,0 +1,281 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/mceliece_key.h>
+#include <botan/internal/bit_ops.h>
+#include <botan/gf2m_small_m.h>
+#include <botan/mceliece.h>
+#include <botan/internal/code_based_key_gen.h>
+#include <botan/alg_id.h>
+#include <botan/der_enc.h>
+#include <botan/ber_dec.h>
+#include <botan/oids.h>
+
+namespace Botan {
+
+McEliece_PrivateKey::McEliece_PrivateKey(polyn_gf2m const& goppa_polyn,
+                                         std::vector<u32bit> const& parity_check_matrix_coeffs,
+                                         std::vector<polyn_gf2m> const& square_root_matrix,
+                                         std::vector<gf2m> const& inverse_support,
+                                         std::vector<byte> const& public_matrix) :
+   McEliece_PublicKey(public_matrix, goppa_polyn.get_degree(), inverse_support.size()),
+   m_g(goppa_polyn),
+   m_sqrtmod(square_root_matrix),
+   m_Linv(inverse_support),
+   m_coeffs(parity_check_matrix_coeffs),
+   m_codimension(ceil_log2(inverse_support.size()) * goppa_polyn.get_degree()),
+   m_dimension(inverse_support.size() - m_codimension)
+   {
+   };
+
+McEliece_PrivateKey::McEliece_PrivateKey(RandomNumberGenerator& rng, size_t code_length, size_t t)
+   {
+   u32bit ext_deg = ceil_log2(code_length);
+   *this = generate_mceliece_key(rng, ext_deg, code_length, t);
+   }
+
+unsigned McEliece_PublicKey::get_message_word_bit_length() const
+   {
+   u32bit codimension = ceil_log2(m_code_length) * m_t;
+   return m_code_length - codimension;
+   }
+
+AlgorithmIdentifier McEliece_PublicKey::algorithm_identifier() const
+   {
+   return AlgorithmIdentifier(get_oid(), std::vector<byte>());
+   }
+
+std::vector<byte> McEliece_PublicKey::x509_subject_public_key() const
+   {
+   // encode the public key
+   return unlock(DER_Encoder()
+                 .start_cons(SEQUENCE)
+                 .start_cons(SEQUENCE)
+                 .encode(static_cast<size_t>(get_code_length()))
+                 .encode(static_cast<size_t>(get_t()))
+                 .end_cons()
+                 .encode(m_public_matrix, OCTET_STRING)
+                 .end_cons()
+                 .get_contents());
+   }
+
+McEliece_PublicKey::McEliece_PublicKey(const McEliece_PublicKey & other) :
+   m_public_matrix(other.m_public_matrix),
+   m_t(other.m_t),
+   m_code_length(other.m_code_length)
+   {
+   }
+
+McEliece_PublicKey::McEliece_PublicKey(const std::vector<byte>& key_bits)
+   {
+   BER_Decoder dec(key_bits);
+   size_t n;
+   size_t t;
+   dec.start_cons(SEQUENCE)
+      .start_cons(SEQUENCE)
+      .decode(n)
+      .decode(t)
+      .end_cons()
+      .decode(m_public_matrix, OCTET_STRING)
+      .end_cons();
+   m_t = t;
+   m_code_length = n;
+   }
+
+secure_vector<byte> McEliece_PrivateKey::pkcs8_private_key() const
+   {
+   DER_Encoder enc;
+   enc.start_cons(SEQUENCE)
+      .start_cons(SEQUENCE)
+      .encode(static_cast<size_t>(get_code_length()))
+      .encode(static_cast<size_t>(get_t()))
+      .end_cons()
+      .encode(m_public_matrix, OCTET_STRING)
+      .encode(m_g.encode(), OCTET_STRING); // g as octet string
+   enc.start_cons(SEQUENCE);
+   for(u32bit i = 0; i < m_sqrtmod.size(); i++)
+      {
+      enc.encode(m_sqrtmod[i].encode(), OCTET_STRING);
+      }
+   enc.end_cons();
+   secure_vector<byte> enc_support;
+   for(u32bit i = 0; i < m_Linv.size(); i++)
+      {
+      enc_support.push_back(m_Linv[i] >> 8);
+      enc_support.push_back(m_Linv[i]);
+      }
+   enc.encode(enc_support, OCTET_STRING);
+   secure_vector<byte> enc_H;
+   for(u32bit i = 0; i < m_coeffs.size(); i++)
+      {
+      enc_H.push_back(m_coeffs[i] >> 24);
+      enc_H.push_back(m_coeffs[i] >> 16);
+      enc_H.push_back(m_coeffs[i] >> 8);
+      enc_H.push_back(m_coeffs[i]);
+      }
+   enc.encode(enc_H, OCTET_STRING);
+   enc.end_cons();
+   return enc.get_contents();
+   }
+
+bool McEliece_PrivateKey::check_key(RandomNumberGenerator& rng, bool) const
+   {
+   const McEliece_PublicKey* p_pk = dynamic_cast<const McEliece_PublicKey*>(this);
+
+   McEliece_Private_Operation priv_op(*this);
+   McEliece_Public_Operation pub_op(*p_pk, p_pk->get_code_length() );
+
+   secure_vector<byte> plaintext((p_pk->get_message_word_bit_length()+7)/8);
+   rng.randomize(&plaintext[0], plaintext.size() - 1);
+   secure_vector<gf2m> err_pos = create_random_error_positions(p_pk->get_code_length(), p_pk->get_t(), rng);
+
+   mceliece_message_parts parts(err_pos, plaintext, p_pk->get_code_length());
+   secure_vector<byte> message_and_error_input = parts.get_concat();
+   secure_vector<byte> ciphertext = pub_op.encrypt(&message_and_error_input[0], message_and_error_input.size(), rng);
+   secure_vector<byte> message_and_error_output = priv_op.decrypt(&ciphertext[0], ciphertext.size() );
+
+   return (message_and_error_input == message_and_error_output);
+   }
+
+McEliece_PrivateKey::McEliece_PrivateKey(const secure_vector<byte>& key_bits)
+   {
+   size_t n, t;
+   secure_vector<byte> g_enc;
+   BER_Decoder dec_base(key_bits);
+   BER_Decoder dec = dec_base.start_cons(SEQUENCE)
+      .start_cons(SEQUENCE)
+      .decode(n)
+      .decode(t)
+      .end_cons()
+      .decode(m_public_matrix, OCTET_STRING)
+      .decode(g_enc, OCTET_STRING);
+
+   if(t == 0 || n == 0)
+      throw Decoding_Error("invalid McEliece parameters");
+
+   u32bit ext_deg = ceil_log2(n);
+   m_code_length = n;
+   m_t = t;
+   m_codimension = (ext_deg * t);
+   m_dimension = (n - m_codimension);
+
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field(new gf2m_small_m::Gf2m_Field(ext_deg));
+   m_g = polyn_gf2m(g_enc, sp_field);
+   if(m_g.get_degree() != static_cast<int>(t))
+      {
+      throw Decoding_Error("degree of decoded Goppa polynomial is incorrect");
+      }
+   BER_Decoder dec2 = dec.start_cons(SEQUENCE);
+   for(u32bit i = 0; i < t/2; i++)
+      {
+      secure_vector<byte> sqrt_enc;
+      dec2.decode(sqrt_enc, OCTET_STRING);
+      while(sqrt_enc.size() < (t*2))
+         {
+         // ensure that the length is always t
+         sqrt_enc.push_back(0);
+         sqrt_enc.push_back(0);
+         }
+      if(sqrt_enc.size() != t*2)
+         {
+         throw Decoding_Error("length of square root polynomial entry is too large");
+         }
+      m_sqrtmod.push_back(polyn_gf2m(sqrt_enc, sp_field));
+      }
+   secure_vector<byte> enc_support;
+   BER_Decoder dec3 = dec2.end_cons()
+      .decode(enc_support, OCTET_STRING);
+   if(enc_support.size() % 2)
+      {
+      throw Decoding_Error("encoded support has odd length");
+      }
+   if(enc_support.size() / 2 != n)
+      {
+      throw Decoding_Error("encoded support has length different from code length");
+      }
+   for(u32bit i = 0; i < n*2; i+=2)
+      {
+      gf2m el = (enc_support[i] << 8) |  enc_support[i+1];
+      m_Linv.push_back(el);
+      }
+   secure_vector<byte> enc_H;
+   dec3.decode(enc_H, OCTET_STRING)
+      .end_cons();
+   if(enc_H.size() % 4)
+      {
+      throw Decoding_Error("encoded parity check matrix has length which is not a multiple of four");
+      }
+   if(enc_H.size()/4 != bit_size_to_32bit_size(m_codimension) * m_code_length )
+      {
+      throw Decoding_Error("encoded parity check matrix has wrong length");
+      }
+
+   for(u32bit i = 0; i < enc_H.size(); i+=4)
+      {
+      u32bit coeff = (enc_H[i] << 24) | (enc_H[i+1] << 16) | (enc_H[i+2] << 8) | enc_H[i+3];
+      m_coeffs.push_back(coeff);
+      }
+
+   }
+
+
+bool McEliece_PrivateKey::operator==(const McEliece_PrivateKey & other) const
+   {
+   if(*static_cast<const McEliece_PublicKey*>(this) != *static_cast<const McEliece_PublicKey*>(&other))
+      {
+      return false;
+      }
+   if(m_g != other.m_g)
+      {
+      return false;
+      }
+
+   if( m_sqrtmod != other.m_sqrtmod)
+      {
+      return false;
+      }
+   if( m_Linv != other.m_Linv)
+      {
+      return false;
+      }
+   if( m_coeffs != other.m_coeffs)
+      {
+      return false;
+      }
+
+   if(m_codimension != other.m_codimension || m_dimension != other.m_dimension)
+      {
+      return false;
+      }
+
+   return true;
+   }
+
+bool McEliece_PublicKey::operator==(const McEliece_PublicKey& other) const
+   {
+   if(m_public_matrix != other.m_public_matrix)
+      {
+      return false;
+      }
+   if(m_t != other.m_t )
+      {
+      return false;
+      }
+   if( m_code_length != other.m_code_length)
+      {
+      return false;
+      }
+   return true;
+   }
+
+}
+
+
diff --git a/src/lib/pubkey/mce/mceliece_key.h b/src/lib/pubkey/mce/mceliece_key.h
new file mode 100644
index 000000000..c51745bba
--- /dev/null
+++ b/src/lib/pubkey/mce/mceliece_key.h
@@ -0,0 +1,127 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef BOTAN_MCELIECE_KEY_H_
+#define BOTAN_MCELIECE_KEY_H_
+
+#include <botan/exceptn.h>
+#include <botan/pk_keys.h>
+#include <botan/polyn_gf2m.h>
+#include <botan/code_based_util.h>
+
+namespace Botan {
+
+class BOTAN_DLL McEliece_PublicKey : public virtual Public_Key
+   {
+   public:
+
+      McEliece_PublicKey(const std::vector<byte>& key_bits);
+
+      McEliece_PublicKey(std::vector<byte> const& pub_matrix, u32bit the_t, u32bit the_code_length) :
+         m_public_matrix(pub_matrix),
+         m_t(the_t),
+         m_code_length(the_code_length)
+            {}
+
+      McEliece_PublicKey(const McEliece_PublicKey & other);
+
+      std::string algo_name() const { return "McEliece/BIGGF2M"; }
+
+      /**
+      * Get the maximum number of bits allowed to be fed to this key.
+      * This is the bitlength of the order of the base point.
+      * @result the maximum number of input bits
+      */
+      size_t max_input_bits() const
+         {
+         return get_message_word_bit_length();
+         };
+
+      AlgorithmIdentifier algorithm_identifier() const;
+
+      size_t estimated_strength() const { return 0; }
+
+      std::vector<byte> x509_subject_public_key() const;
+
+      bool check_key(RandomNumberGenerator&, bool) const
+         { return true; }
+
+      u32bit get_t() const { return m_t; }
+      u32bit get_code_length() const { return m_code_length; }
+      u32bit get_message_word_bit_length() const;
+      std::vector<byte> const& get_public_matrix() const { return m_public_matrix; }
+
+      bool operator==(const McEliece_PublicKey& other) const;
+      bool operator!=(const McEliece_PublicKey& other) const { return !(*this == other); }
+
+   protected:
+      McEliece_PublicKey() {}
+
+      std::vector<byte> m_public_matrix;
+      u32bit m_t;
+      u32bit m_code_length;
+   };
+
+class BOTAN_DLL McEliece_PrivateKey : public virtual McEliece_PublicKey,
+                                      public virtual Private_Key
+   {
+   public:
+      /**
+      * Get the maximum number of bits allowed to be fed to this key.
+      * This is the bitlength of the order of the base point.
+      * @result the maximum number of input bits
+      */
+      size_t max_input_bits() const {
+      return m_Linv.size();
+      };
+
+      McEliece_PrivateKey(const secure_vector<byte>& key_bits);
+
+      McEliece_PrivateKey(polyn_gf2m const& goppa_polyn,
+                          std::vector<u32bit> const& parity_check_matrix_coeffs,
+                          std::vector<polyn_gf2m> const& square_root_matrix,
+                          std::vector<gf2m> const& inverse_support,
+                          std::vector<byte> const& public_matrix );
+
+      McEliece_PrivateKey(RandomNumberGenerator& rng, size_t code_length, size_t t);
+      bool check_key(RandomNumberGenerator& rng, bool strong) const;
+
+      polyn_gf2m const& get_goppa_polyn() const { return m_g; };
+      std::vector<u32bit> const& get_H_coeffs() const { return m_coeffs; };
+      std::vector<gf2m> const& get_Linv() const { return m_Linv; };
+      std::vector<polyn_gf2m> const& get_sqrtmod() const { return m_sqrtmod; };
+
+      inline u32bit get_dimension() const
+         { return m_dimension; };
+
+      inline u32bit get_codimension() const
+         { return m_codimension; };
+
+
+      secure_vector<byte> pkcs8_private_key() const;
+
+      bool operator==(const McEliece_PrivateKey & other) const;
+
+      bool operator!=(const McEliece_PrivateKey& other) const { return !(*this == other); };
+
+   private:
+      polyn_gf2m m_g;
+      std::vector<polyn_gf2m> m_sqrtmod;
+      std::vector<gf2m> m_Linv;
+      std::vector<u32bit> m_coeffs;
+
+      u32bit m_codimension;
+      u32bit m_dimension;
+   };
+
+}
+
+#endif
diff --git a/src/lib/pubkey/mce/polyn_gf2m.cpp b/src/lib/pubkey/mce/polyn_gf2m.cpp
new file mode 100644
index 000000000..a0da9b479
--- /dev/null
+++ b/src/lib/pubkey/mce/polyn_gf2m.cpp
@@ -0,0 +1,804 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#include <botan/polyn_gf2m.h>
+#include <cstring>
+#include <iostream>
+
+#include <botan/gf2m_rootfind_dcmp.h>
+#include <botan/code_based_util.h>
+#include <botan/gf2m_small_m.h>
+#include <botan/internal/bit_ops.h>
+
+namespace Botan {
+
+using namespace Botan::gf2m_small_m;
+
+namespace {
+
+gf2m generate_gf2m_mask(gf2m a)
+   {
+   gf2m result =  (a != 0);
+   return ~(result - 1);
+   }
+
+unsigned nlz_16bit(u16bit x)
+   {
+   unsigned n;
+   if(x == 0) return 16;
+   n = 0;
+   if(x <= 0x00FF) {n = n + 8; x = x << 8;}
+   if(x <= 0x0FFF) {n = n + 4; x = x << 4;}
+   if(x <= 0x3FFF) {n = n + 2; x = x << 2;}
+   if(x <= 0x7FFF) {n = n + 1;}
+   return n;
+   }
+}
+
+using namespace Botan::gf2m_small_m;
+
+int polyn_gf2m::calc_degree_secure() const
+   {
+   int i = this->coeff.size() - 1;
+   int result = 0;
+   u32bit found_mask = 0;
+   u32bit tracker_mask = 0xffff;
+   for( ; i >= 0; i--)
+      {
+      found_mask = expand_mask_16bit(this->coeff[i]);
+      result |= i & found_mask & tracker_mask;
+      // tracker mask shall become zero once found mask is set
+      // it shall remain zero from then on
+      tracker_mask = tracker_mask & ~found_mask;
+      }
+   const_cast<polyn_gf2m*>(this)->m_deg = result;
+   return result;
+   }
+/**
+* number of leading zeros
+*/
+
+gf2m random_code_element(unsigned code_length, Botan::RandomNumberGenerator& rng)
+   {
+   if(code_length == 0)
+      {
+      throw Invalid_Argument("random_code_element() was supplied a code length of zero");
+      }
+   unsigned nlz = nlz_16bit(code_length-1);
+   gf2m mask = (1 << (16-nlz)) -1;
+   gf2m result;
+   do
+      {
+      rng.randomize((byte*)&result, sizeof(result));
+
+      result &= mask;
+      } while(result >= code_length); // rejection sampling
+   return result;
+   }
+
+polyn_gf2m::polyn_gf2m(polyn_gf2m const& other)
+   :m_deg(other.m_deg),
+    coeff(other.coeff),
+    msp_field(other.msp_field)
+   { }
+
+polyn_gf2m::polyn_gf2m(   int d, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field)
+   :m_deg(-1),
+    coeff(d+1),
+    msp_field(sp_field)
+   {
+   }
+
+std::string polyn_gf2m::to_string() const
+   {
+   int d = get_degree();
+   std::string result;
+   for(int i = 0; i < d + 1; i ++)
+      {
+      result += std::to_string(this->coeff[i]);
+      if(i != d)
+         {
+         result += ", ";
+         }
+      }
+   return result;
+   }
+/**
+* doesn't save coefficients:
+*/
+void polyn_gf2m::realloc(u32bit new_size)
+   {
+   this->coeff = secure_vector<gf2m>(new_size);
+   }
+
+polyn_gf2m::polyn_gf2m(const byte* mem, u32bit mem_len, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field)
+   :msp_field(sp_field)
+   {
+   if(mem_len % sizeof(gf2m))
+      {
+      throw new Botan::Decoding_Error("illegal length of memory to decode ");
+      }
+
+   u32bit size = (mem_len / sizeof(this->coeff[0])) ;
+   this->coeff = secure_vector<gf2m>(size);
+   this->m_deg = -1;
+   for(u32bit i = 0; i < size; i++)
+      {
+      this->coeff[i] = gf2m_small_m::decode_gf2m(mem);
+      mem += sizeof(this->coeff[0]);
+      }
+   for(u32bit i = 0; i < size; i++)
+      {
+      if(this->coeff[i] >= (1 << sp_field->get_extension_degree()))
+         {
+         throw Botan::Decoding_Error("error decoding polynomial");
+         }
+      }
+   this->get_degree();
+   }
+
+
+polyn_gf2m::polyn_gf2m( std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field )
+   : m_deg(-1),
+     coeff(1),
+     msp_field(sp_field)
+   {};
+
+polyn_gf2m::polyn_gf2m(int degree, const unsigned  char* mem, u32bit mem_byte_len, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field)
+   :msp_field(sp_field)
+   {
+   u32bit j, k, l;
+   gf2m a;
+   u32bit polyn_size;
+   polyn_size = degree + 1;
+   if(polyn_size * sp_field->get_extension_degree() > 8 * mem_byte_len)
+      {
+      throw Botan::Decoding_Error("memory vector for polynomial has wrong size");
+      }
+   this->coeff = secure_vector<gf2m>(degree+1);
+   gf2m ext_deg = this->msp_field->get_extension_degree();
+   for (l = 0; l < polyn_size; l++)
+      {
+      k = (l * ext_deg) / 8;
+
+      j = (l * ext_deg) % 8;
+      a = mem[k] >> j;
+      if (j + ext_deg > 8)
+         {
+         a ^= mem[k + 1] << (8- j);
+         }
+      if(j + ext_deg > 16)
+         {
+         a ^= mem[k + 2] << (16- j);
+         }
+      a &= ((1 << ext_deg) - 1);
+      (*this).set_coef( l, a);
+      }
+
+   this->get_degree();
+   }
+
+#if 0
+void polyn_gf2m::encode(u32bit min_numo_coeffs, byte* mem, u32bit mem_len) const
+   {
+   u32bit i;
+   u32bit numo_coeffs, needed_size;
+   this->get_degree();
+   numo_coeffs = (min_numo_coeffs > static_cast<u32bit>(this->m_deg+1)) ? min_numo_coeffs : this->m_deg+1;
+   needed_size = sizeof(this->coeff[0]) * numo_coeffs;
+   if(mem_len < needed_size)
+      {
+      Invalid_Argument("provided memory too small to encode polynomial");
+      }
+
+   for(i = 0; i < numo_coeffs; i++)
+      {
+      gf2m to_enc;
+      if(i >= static_cast<u32bit>(this->m_deg+1))
+         {
+         /* encode a zero */
+         to_enc = 0;
+         }
+      else
+         {
+         to_enc = this->coeff[i];
+         }
+      mem += encode_gf2m(to_enc, mem);
+      }
+   }
+#endif
+
+
+void polyn_gf2m::set_to_zero()
+   {
+   memset(&this->coeff[0], 0, this->coeff.size() * sizeof (gf2m));
+   this->m_deg = -1;
+   }
+
+int polyn_gf2m::get_degree() const
+   {
+   int d = this->coeff.size() - 1;
+   while ((d >= 0) && (this->coeff[d] == 0))
+      --d;
+   const_cast<polyn_gf2m*>(this)->m_deg = d;
+   return d;
+   }
+
+
+static gf2m eval_aux(const gf2m * /*restrict*/ coeff, gf2m a, int d, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field)
+   {
+   gf2m b;
+   b = coeff[d--];
+   for (; d >= 0; --d)
+      if (b != 0)
+         {
+         b = sp_field->gf_mul(b, a) ^ coeff[d];
+         }
+      else
+         {
+         b = coeff[d];
+         }
+   return b;
+   }
+
+gf2m polyn_gf2m::eval(gf2m a)
+   {
+   return eval_aux(&this->coeff[0], a, this->m_deg, this->msp_field);
+   }
+
+
+// p will contain it's remainder modulo g
+void polyn_gf2m::remainder(polyn_gf2m &p, const polyn_gf2m & g)
+   {
+   int i, j, d;
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> msp_field = g.msp_field;
+   d = p.get_degree() - g.get_degree();
+   if (d >= 0) {
+   gf2m la = msp_field->gf_inv_rn(g.get_lead_coef());
+
+   for (i = p.get_degree(); d >= 0; --i, --d) {
+   if (p[i] != 0) {
+   gf2m lb = msp_field->gf_mul_rrn(la, p[i]);
+   for (j = 0; j < g.get_degree(); ++j)
+      {
+      p[j+d] ^= msp_field->gf_mul_zrz(lb, g[j]);
+      }
+   (*&p).set_coef( i, 0);
+   }
+   }
+   p.set_degree( g.get_degree() - 1);
+   while ((p.get_degree() >= 0) && (p[p.get_degree()] == 0))
+      p.set_degree( p.get_degree() - 1);
+   }
+   }
+
+std::vector<polyn_gf2m> polyn_gf2m::sqmod_init(const polyn_gf2m & g)
+   {
+   std::vector<polyn_gf2m> sq;
+   int signed_deg = g.get_degree();
+   if(signed_deg <= 0)
+      {
+      throw Invalid_Argument("cannot compute sqmod for such low degree");
+      }
+   u32bit d = (u32bit) signed_deg;
+   u32bit t = g.m_deg;
+   // create t zero polynomials
+   u32bit i;
+   for (i = 0; i < t; ++i)
+      {
+      sq.push_back(polyn_gf2m(t+1, g.get_sp_field()));
+      }
+   for (i = 0; i < d / 2; ++i)
+      {
+      sq[i].set_degree( 2 * i);
+      (*&sq[i]).set_coef( 2 * i, 1);
+      }
+
+   for (; i < d; ++i)
+      {
+      memset(&sq[i].coeff[0], 0, 2 * sizeof (gf2m));
+      memcpy(&sq[i].coeff[0] + 2, &sq[i - 1].coeff[0], d * sizeof (gf2m));
+      sq[i].set_degree( sq[i - 1].get_degree() + 2);
+      polyn_gf2m::remainder(sq[i], g);
+      }
+   return sq;
+   }
+
+/*Modulo p square of a certain polynomial g, sq[] contains the square
+Modulo g of the base canonical polynomials of degree < d, where d is
+the degree of G. The table sq[] will be calculated by polyn_gf2m__sqmod_init*/
+polyn_gf2m polyn_gf2m::sqmod( const std::vector<polyn_gf2m> & sq, int d)
+   {
+   int i, j;
+   gf2m la;
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field = this->msp_field;
+
+   polyn_gf2m result(d - 1, sp_field);
+   // terms of low degree
+   for (i = 0; i < d / 2; ++i)
+      {
+      (*&result).set_coef( i * 2, sp_field->gf_square((*this)[i]));
+      }
+
+   // terms of high degree
+   for (; i < d; ++i)
+      {
+      gf2m lpi = (*this)[i];
+      if (lpi != 0)
+         {
+         lpi = sp_field->gf_log(lpi);
+         la = sp_field->gf_mul_rrr(lpi, lpi);
+         for (j = 0; j < d; ++j)
+            {
+            result[j] ^= sp_field->gf_mul_zrz(la, sq[i][j]);
+            }
+         }
+      }
+
+   // Update degre
+   result.set_degree( d - 1);
+   while ((result.get_degree() >= 0) && (result[result.get_degree()] == 0))
+      result.set_degree( result.get_degree() - 1);
+   return result;
+   }
+
+
+// destructive
+polyn_gf2m polyn_gf2m::gcd_aux(polyn_gf2m& p1, polyn_gf2m& p2)
+   {
+   if (p2.get_degree() == -1)
+      return p1;
+   else {
+   polyn_gf2m::remainder(p1, p2);
+   return polyn_gf2m::gcd_aux(p2, p1);
+   }
+   }
+
+
+polyn_gf2m polyn_gf2m::gcd(polyn_gf2m const& p1, polyn_gf2m const& p2)
+   {
+   polyn_gf2m a(p1);
+   polyn_gf2m b(p2);
+   if (a.get_degree() < b.get_degree())
+      {
+      return polyn_gf2m(polyn_gf2m::gcd_aux(b, a));
+      }
+   else
+      {
+      return polyn_gf2m(polyn_gf2m::gcd_aux(a, b));
+      }
+   }
+
+
+
+
+
+// Returns the degree of the smallest factor
+void polyn_gf2m::degppf(const polyn_gf2m & g, int* p_result)
+   {
+   int i, d;
+   polyn_gf2m s(g.get_sp_field());
+
+   d = g.get_degree();
+   std::vector<polyn_gf2m> u = polyn_gf2m::sqmod_init(g);
+
+   polyn_gf2m p( d - 1, g.msp_field);
+
+   p.set_degree( 1);
+   (*&p).set_coef( 1, 1);
+   (*p_result) = d;
+   for (i = 1; i <= (d / 2) * g.msp_field->get_extension_degree(); ++i)
+      {
+      polyn_gf2m r = p.sqmod(u, d);
+      if ((i % g.msp_field->get_extension_degree()) == 0)
+         {
+         r[1] ^= 1;
+         r.get_degree(); // The degree may change
+         s = polyn_gf2m::gcd( g, r);
+
+         if (s.get_degree() > 0)
+            {
+            (*p_result) = i / g.msp_field->get_extension_degree();
+            break;
+            }
+         r[1] ^= 1;
+         r.get_degree(); // The degree may change
+         }
+      // No need for the exchange s
+      s = p;
+      p = r;
+      r = s;
+      }
+
+
+   }
+
+void polyn_gf2m::patchup_deg_secure( u32bit trgt_deg, volatile gf2m patch_elem)
+   {
+   u32bit i;
+   if(this->coeff.size() < trgt_deg)
+      {
+      return;
+      }
+   for(i = 0; i < this->coeff.size(); i++)
+      {
+      u32bit equal, equal_mask;
+      this->coeff[i] |= patch_elem;
+      equal = (i == trgt_deg);
+      equal_mask = expand_mask_16bit(equal);
+      patch_elem &= ~equal_mask;
+      }
+   this->calc_degree_secure();
+   }
+// We suppose m_deg(g) >= m_deg(p)
+// v is the problem
+std::pair<polyn_gf2m, polyn_gf2m> polyn_gf2m::eea_with_coefficients( const polyn_gf2m & p, const polyn_gf2m & g, int break_deg)
+   {
+
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> msp_field = g.msp_field;
+   int i, j, dr, du, delta;
+   gf2m a;
+   polyn_gf2m aux;
+
+   // initialisation of the local variables
+   // r0 <- g, r1 <- p, u0 <- 0, u1 <- 1
+   dr = g.get_degree();
+
+   polyn_gf2m r0(dr, g.msp_field);
+   polyn_gf2m r1(dr - 1, g.msp_field);
+   polyn_gf2m u0(dr - 1, g.msp_field);
+   polyn_gf2m u1(dr - 1, g.msp_field);
+
+   r0 = g;
+   r1 = p;
+   u0.set_to_zero();
+   u1.set_to_zero();
+   (*&u1).set_coef( 0, 1);
+   u1.set_degree( 0);
+
+
+   // invariants:
+   // r1 = u1 * p + v1 * g
+   // r0 = u0 * p + v0 * g
+   // and m_deg(u1) = m_deg(g) - m_deg(r0)
+   // It stops when m_deg (r1) <t (m_deg (r0)> = t)
+   // And therefore m_deg (u1) = m_deg (g) - m_deg (r0) <m_deg (g) - break_deg
+   du = 0;
+   dr = r1.get_degree();
+   delta = r0.get_degree() - dr;
+
+
+   while (dr >= break_deg)
+      {
+
+      for (j = delta; j >= 0; --j)
+         {
+         a = msp_field->gf_div(r0[dr + j], r1[dr]);
+         if (a != 0)
+            {
+            gf2m la = msp_field->gf_log(a);
+            // u0(z) <- u0(z) + a * u1(z) * z^j
+            for (i = 0; i <= du; ++i)
+               {
+               u0[i + j] ^= msp_field->gf_mul_zrz(la, u1[i]);
+               }
+            // r0(z) <- r0(z) + a * r1(z) * z^j
+            for (i = 0; i <= dr; ++i)
+               {
+               r0[i + j] ^= msp_field->gf_mul_zrz(la, r1[i]);
+               }
+            }
+         } // end loop over j
+
+      if(break_deg != 1) /* key eq. solving */
+         {
+         /* [ssms_icisc09] Countermeasure
+         * d_break from paper equals break_deg - 1
+         * */
+
+         volatile gf2m fake_elem = 0x01;
+         volatile gf2m cond1, cond2;
+         int trgt_deg = r1.get_degree() - 1;
+         r0.calc_degree_secure();
+         u0.calc_degree_secure();
+         if(!(g.get_degree() % 2))
+            {
+            /* t even */
+            cond1 = r0.get_degree() < break_deg - 1;
+            }
+         else
+            {
+            /* t odd */
+            cond1 =  r0.get_degree() <= break_deg - 1;
+            cond2 =  u0.get_degree() < break_deg - 1;
+            cond1 &= cond2;
+            }
+         /* expand cond1 to a full mask */
+         //CSEC_MASK__GEN_MASK_16B(cond1, mask);
+         gf2m mask = generate_gf2m_mask(cond1);
+         fake_elem &= mask;
+         r0.patchup_deg_secure(trgt_deg, fake_elem);
+         }
+      if(break_deg == 1) /* syndrome inversion */
+         {
+         volatile gf2m fake_elem = 0x00;
+         volatile u32bit trgt_deg = 0;
+         r0.calc_degree_secure();
+         u0.calc_degree_secure();
+         /**
+         * countermeasure against the low weight attacks for w=4, w=6 and w=8.
+         * Higher values are not covered since for w=8 we already have a
+         * probability for a positive of 1/n^3 from random ciphertexts with the
+         * given weight. For w = 10 it would be 1/n^4 and so on. Thus attacks
+         * based on such high values of w are considered impractical.
+         *
+         * The outer test for the degree of u ( Omega in the paper ) needs not to
+         * be disguised. Each of the three is performed at most once per EEA
+         * (syndrome inversion) execution, the attacker knows this already when
+         * preparing the ciphertext with the given weight. Inside these three
+         * cases however, we must use timing neutral (branch free) operations to
+         * implement the condition detection and the counteractions.
+         *
+         */
+         if(u0.get_degree() == 4)
+            {
+            u32bit mask = 0;
+            /**
+            * Condition that the EEA would break now
+            */
+            int cond_r = r0.get_degree() == 0;
+            /**
+            * Now come the conditions for all odd coefficients of this sigma
+            * candiate. If they are all fulfilled, then we know that we have a low
+            * weight error vector, since the key-equation solving EEA is skipped if
+            * the degree of tau^2 is low (=m_deg(u0)) and all its odd cofficients are
+            * zero (they would cause "full-lenght" contributions from the square
+            * root computation).
+            */
+            // Condition for the coefficient to Y to be cancelled out by the
+            // addition of Y before the square root computation:
+            int cond_u1 = msp_field->gf_mul(u0.coeff[1], msp_field->gf_inv(r0.coeff[0])) == 1;
+
+            // Condition sigma_3 = 0:
+            int cond_u3 = u0.coeff[3] == 0;
+            // combine the conditions:
+            cond_r &= (cond_u1 & cond_u3);
+            // mask generation:
+            mask = expand_mask_16bit(cond_r);
+            trgt_deg = 2 & mask;
+            fake_elem = 1 & mask;
+            }
+         else if(u0.get_degree() == 6)
+            {
+            u32bit mask = 0;
+            int cond_r= r0.get_degree() == 0;
+            int cond_u1 = msp_field->gf_mul(u0.coeff[1], msp_field->gf_inv(r0.coeff[0])) == 1;
+            int cond_u3 = u0.coeff[3] == 0;
+
+            int cond_u5 = u0.coeff[5] == 0;
+
+            cond_r &= (cond_u1 & cond_u3 & cond_u5);
+            mask = expand_mask_16bit(cond_r);
+            trgt_deg = 4 & mask;
+            fake_elem = 1 & mask;
+            }
+         else if(u0.get_degree() == 8)
+            {
+            u32bit mask = 0;
+            int cond_r= r0.get_degree() == 0;
+            int cond_u1 = msp_field->gf_mul(u0[1], msp_field->gf_inv(r0[0])) == 1;
+            int cond_u3 = u0.coeff[3] == 0;
+
+            int cond_u5 = u0.coeff[5] == 0;
+
+            int cond_u7 = u0.coeff[7] == 0;
+
+            cond_r &= (cond_u1 & cond_u3 & cond_u5 & cond_u7);
+            mask = expand_mask_16bit(cond_r);
+            trgt_deg = 6 & mask;
+            fake_elem = 1 & mask;
+            }
+         r0.patchup_deg_secure(trgt_deg, fake_elem);
+         }
+      // exchange
+      aux = r0; r0 = r1; r1 = aux;
+      aux = u0; u0 = u1; u1 = aux;
+
+      du = du + delta;
+      delta = 1;
+      while (r1[dr - delta] == 0)
+         {
+         delta++;
+         }
+
+
+      dr -= delta;
+      } /* end  while loop (dr >= break_deg) */
+
+
+   u1.set_degree( du);
+   r1.set_degree( dr);
+   //return u1 and r1;
+   return std::make_pair(u1,r1); // coefficients u,v
+   }
+
+polyn_gf2m::polyn_gf2m(int t, Botan::RandomNumberGenerator& rng, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field)
+   :m_deg(t),
+    coeff(t+1),
+    msp_field(sp_field)
+   {
+   int i;
+   (*this).set_coef( t, 1);
+   i = 0;
+   int m_deg;
+   do
+      {
+      for (i = 0; i < t; ++i)
+         {
+         (*this).set_coef( i, random_code_element(sp_field->get_cardinality(), rng));
+         }
+      polyn_gf2m::degppf(*this, &m_deg);
+      }
+   while (m_deg < t);
+   }
+
+
+void polyn_gf2m::poly_shiftmod( const polyn_gf2m & g)
+   {
+   int i, t;
+   gf2m a;
+
+   if(g.get_degree() <= 0)
+      {
+      throw Invalid_Argument("shiftmod cannot be called on polynomials of degree 0 or less");
+      }
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> msp_field = g.msp_field;
+
+   t = g.get_degree();
+   a = msp_field->gf_div(this->coeff[t-1], g.coeff[t]);
+   for (i = t - 1; i > 0; --i)
+      {
+      this->coeff[i] = this->coeff[i - 1] ^ this->msp_field->gf_mul(a, g.coeff[i]);
+      }
+   this->coeff[0] = msp_field->gf_mul(a, g.coeff[0]);
+   }
+
+std::vector<polyn_gf2m> polyn_gf2m::sqrt_mod_init(const polyn_gf2m & g)
+   {
+   u32bit i, t;
+   u32bit nb_polyn_sqrt_mat;
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> msp_field = g.msp_field;
+   std::vector<polyn_gf2m> result;
+   t = g.get_degree();
+   nb_polyn_sqrt_mat = t/2;
+
+   std::vector<polyn_gf2m> sq_aux = polyn_gf2m::sqmod_init(g);
+
+
+   polyn_gf2m p( t - 1, g.get_sp_field());
+   p.set_degree( 1);
+
+   (*&p).set_coef( 1, 1);
+   // q(z) = 0, p(z) = z
+   for (i = 0; i < t * msp_field->get_extension_degree() - 1; ++i)
+      {
+      // q(z) <- p(z)^2 mod g(z)
+      polyn_gf2m q = p.sqmod(sq_aux, t);
+      // q(z) <-> p(z)
+      polyn_gf2m aux = q;
+      q = p;
+      p = aux;
+      }
+   // p(z) = z^(2^(tm-1)) mod g(z) = sqrt(z) mod g(z)
+
+   for (i = 0; i < nb_polyn_sqrt_mat; ++i)
+      {
+      result.push_back(polyn_gf2m(t - 1, g.get_sp_field()));
+      }
+
+   result[0] = p;
+   result[0].get_degree();
+   for(i = 1; i < nb_polyn_sqrt_mat; i++)
+      {
+      result[i] = result[i - 1];
+      result[i].poly_shiftmod(g),
+         result[i].get_degree();
+      }
+
+   return result;
+   }
+
+std::vector<polyn_gf2m> syndrome_init(polyn_gf2m const& generator, std::vector<gf2m> const& support, int n)
+   {
+   int i,j,t;
+   gf2m a;
+
+
+   std::shared_ptr<gf2m_small_m::Gf2m_Field> msp_field = generator.msp_field;
+
+   std::vector<polyn_gf2m> result;
+   t = generator.get_degree();
+
+   //g(z)=g_t+g_(t-1).z^(t-1)+......+g_1.z+g_0
+   //f(z)=f_(t-1).z^(t-1)+......+f_1.z+f_0
+
+   for(j=0;j<n;j++)
+      {
+      result.push_back(polyn_gf2m( t-1, msp_field));
+
+      (*&result[j]).set_coef(t-1,1);
+      for(i=t-2;i>=0;i--)
+         {
+         (*&result[j]).set_coef(i, (generator)[i+1]  ^
+                                msp_field->gf_mul(lex_to_gray(support[j]),result[j][i+1]));
+         }
+      a = ((generator)[0] ^ msp_field->gf_mul(lex_to_gray(support[j]),result[j][0]));
+      for(i=0;i<t;i++)
+         {
+         (*&result[j]).set_coef(i, msp_field->gf_div(result[j][i],a));
+         }
+      }
+   return result;
+   }
+
+polyn_gf2m::polyn_gf2m(const secure_vector<byte>& encoded, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field )
+   :msp_field(sp_field)
+   {
+   if(encoded.size() % 2)
+      {
+      throw Decoding_Error("encoded polynomial has odd length");
+      }
+   for(u32bit i = 0; i < encoded.size(); i += 2)
+      {
+      gf2m el = (encoded[i] << 8) | encoded[i + 1];
+      coeff.push_back(el);
+      }
+   get_degree();
+
+   }
+secure_vector<byte> polyn_gf2m::encode() const
+   {
+   secure_vector<byte> result;
+
+   if(m_deg < 1)
+      {
+      result.push_back(0);
+      result.push_back(0);
+      return result;
+      }
+
+   u32bit len = m_deg+1;
+   for(unsigned i = 0; i < len; i++)
+      {
+      // "big endian" encoding of the GF(2^m) elements
+      result.push_back(coeff[i] >> 8);
+      result.push_back(coeff[i]);
+      }
+   return result;
+   }
+
+void polyn_gf2m::swap(polyn_gf2m& other)
+   {
+   std::swap(this->m_deg, other.m_deg);
+   std::swap(this->msp_field, other.msp_field);
+   std::swap(this->coeff, other.coeff);
+   }
+
+bool polyn_gf2m::operator==(const polyn_gf2m & other) const
+   {
+   if(m_deg != other.m_deg || coeff != other.coeff)
+      {
+      return false;
+      }
+   return true;
+   }
+
+}
diff --git a/src/lib/pubkey/mce/polyn_gf2m.h b/src/lib/pubkey/mce/polyn_gf2m.h
new file mode 100644
index 000000000..9dff549e0
--- /dev/null
+++ b/src/lib/pubkey/mce/polyn_gf2m.h
@@ -0,0 +1,161 @@
+/**
+ * (C) Copyright Projet SECRET, INRIA, Rocquencourt
+ * (C) Bhaskar Biswas and  Nicolas Sendrier
+ *
+ * (C) 2014 cryptosource GmbH
+ * (C) 2014 Falko Strenzke [email protected]
+ *
+ * Distributed under the terms of the Botan license
+ *
+ */
+
+#ifndef BOTAN_POLYN_GF2M_H__
+#define BOTAN_POLYN_GF2M_H__
+
+#include <botan/gf2m_small_m.h>
+#include <botan/rng.h>
+#include <memory>
+#include <utility>
+
+namespace Botan {
+
+using namespace gf2m_small_m;
+
+struct polyn_gf2m
+   {
+   public:
+      /**
+      * create a zero polynomial:
+      */
+      polyn_gf2m( std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field );
+
+      polyn_gf2m()
+         :m_deg(-1)
+         {};
+
+      polyn_gf2m(const secure_vector<byte>& encoded, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field );
+
+      polyn_gf2m& operator=(const polyn_gf2m&) = default;
+
+      bool operator==(const polyn_gf2m & other) const ;
+
+      bool operator!=(const polyn_gf2m & other) const { return !(*this == other); };
+
+      polyn_gf2m(polyn_gf2m&& other)
+         {
+         this->swap(other);
+         };
+
+      polyn_gf2m & operator=(polyn_gf2m&& other)
+         {
+         if(this != &other)
+            {
+            this->swap(other);
+            }
+         return *this;
+         }
+
+      void swap(polyn_gf2m& other);
+
+      secure_vector<byte> encode() const;
+      /**
+      * create zero polynomial with reservation of space for a degree d polynomial
+      */
+      polyn_gf2m(int d, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field);
+
+      polyn_gf2m(polyn_gf2m const& other);
+      /**
+      * create zero polynomial with allocated size determined by specified degree d:
+      */
+
+      /**
+      * random irreducible polynomial of degree t
+      */
+      polyn_gf2m(int t, RandomNumberGenerator& rng, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field);
+
+      std::shared_ptr<gf2m_small_m::Gf2m_Field> get_sp_field() const
+         { return msp_field; };
+
+      gf2m& operator[](size_t i) { return coeff[i]; };
+
+      gf2m operator[](size_t i) const { return coeff[i]; }
+
+      gf2m get_lead_coef() const { return coeff[m_deg]; }
+
+      gf2m get_coef(u32bit i) const { return coeff[i]; }
+
+      inline void set_coef(u32bit i, gf2m v)
+         {
+         coeff[i] = v;
+         };
+
+      inline void add_to_coef(u32bit i, gf2m v)
+         {
+         coeff[i] = coeff[i] ^ v;
+         }
+
+      std::string to_string() const;
+
+      /** decode a polynomial from memory: **/
+      polyn_gf2m(const byte* mem, u32bit mem_len, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field);
+      // remove one! ^v!
+      /**
+      *  create a polynomial from memory area (encoded)
+      */
+      polyn_gf2m(int degree, const unsigned  char* mem, u32bit mem_byte_len, std::shared_ptr<gf2m_small_m::Gf2m_Field> sp_field);
+
+      void encode(u32bit min_numo_coeffs, byte* mem, u32bit mem_len) const;
+
+      int get_degree() const;
+
+      /**
+      * determine the degree in a timing secure manner. the timing of this function
+      * only depends on the number of allocated coefficients, not on the actual
+      * degree
+      */
+      int calc_degree_secure() const;
+
+      void degppf(const polyn_gf2m & g, int* p_result);
+
+      static std::vector<polyn_gf2m> sqmod_init(const polyn_gf2m & g);
+
+      static std::vector<polyn_gf2m> sqrt_mod_init(const polyn_gf2m & g);
+
+
+      polyn_gf2m sqmod(const std::vector<polyn_gf2m> & sq, int d);
+      void set_to_zero();
+      gf2m eval(gf2m a);
+
+      static std::pair<polyn_gf2m, polyn_gf2m> eea_with_coefficients(const polyn_gf2m & p,
+                                                                     const polyn_gf2m & g,
+                                                                     int break_deg);
+
+      void patchup_deg_secure( u32bit trgt_deg, volatile gf2m patch_elem);
+
+   private:
+
+      void set_degree(int d) { m_deg = d; }
+
+      void poly_shiftmod( const polyn_gf2m & g);
+      void realloc(u32bit new_size);
+      static polyn_gf2m gcd(polyn_gf2m const& p1, polyn_gf2m const& p2);
+
+      /**
+      * destructive:
+      */
+      static void remainder(polyn_gf2m & p, const polyn_gf2m & g);
+
+      static polyn_gf2m gcd_aux(polyn_gf2m& p1, polyn_gf2m& p2);
+   public:
+      int m_deg;
+      secure_vector<gf2m> coeff;
+      std::shared_ptr<gf2m_small_m::Gf2m_Field> msp_field;
+   };
+
+gf2m random_code_element(unsigned code_length, RandomNumberGenerator& rng);
+
+std::vector<polyn_gf2m> syndrome_init(polyn_gf2m const& generator, std::vector<gf2m> const& support, int n);
+
+}
+
+#endif
diff --git a/src/lib/utils/bit_ops.h b/src/lib/utils/bit_ops.h
index 0072fde71..75a7584ad 100644
--- a/src/lib/utils/bit_ops.h
+++ b/src/lib/utils/bit_ops.h
@@ -1,6 +1,10 @@
 /*
 * Bit/Word Operations
 * (C) 1999-2008 Jack Lloyd
+* (C) Copyright Projet SECRET, INRIA, Rocquencourt
+* (C) Bhaskar Biswas and  Nicolas Sendrier
+* (C) 2014 cryptosource GmbH
+* (C) 2014 Falko Strenzke [email protected]
 *
 * Distributed under the terms of the Botan license
 */
@@ -98,6 +102,24 @@ inline size_t ctz(T n)
    return 8*sizeof(T);
    }
 
+template<typename T>
+size_t ceil_log2(T x)
+   {
+   if(x >> (sizeof(T)*8-1))
+      return sizeof(T)*8;
+
+   size_t result = 0;
+   T compare = 1;
+
+   while(compare < x)
+      {
+      compare <<= 1;
+      result++;
+      }
+
+   return result;
+   }
+
 }
 
 #endif
diff --git a/src/lib/utils/ta_utils.cpp b/src/lib/utils/ta_utils.cpp
index 86cf25969..9a2c0df49 100644
--- a/src/lib/utils/ta_utils.cpp
+++ b/src/lib/utils/ta_utils.cpp
@@ -22,30 +22,42 @@ namespace TA_CM {
 * anywhere.
 */
 
-u32bit gen_mask_u32bit(u32bit in)
+namespace {
+
+template<typename T>
+T expand_mask(T x)
+   {
+   volatile T r = x;
+   for(size_t i = 1; i != sizeof(T) * 8; i *= 2)
+      r |= r >> i;
+   r &= 1;
+   r = ~(r - 1);
+   return r;
+   }
+
+}
+
+u32bit expand_mask_u32bit(u32bit in)
+   {
+   return expand_mask<u32bit>(in);
+   }
+
+u16bit expand_mask_u16bit(u16bit in)
    {
-   volatile u32bit result = in;
-   result |= result >> 1;
-   result |= result >> 2;
-   result |= result >> 4;
-   result |= result >> 8;
-   result |= result >> 16;
-   result &= 1;
-   result = ~(result - 1);
-   return result;
+   return expand_mask<u16bit>(in);
    }
 
 u32bit max_32(u32bit a, u32bit b)
    {
    const u32bit a_larger = b - a; /* negative if a larger */
-   const u32bit mask = gen_mask_u32bit(a_larger >> 31);
+   const u32bit mask = expand_mask<u32bit>(a_larger >> 31);
    return (a & mask) | (b & ~mask);
    }
 
 u32bit min_32(u32bit a, u32bit b)
    {
    const u32bit a_larger = b - a; /* negative if a larger */
-   const u32bit mask = gen_mask_u32bit(a_larger >> 31);
+   const u32bit mask = expand_mask<u32bit>(a_larger >> 31);
    return (a & ~mask) | (b & mask);
    }
 
diff --git a/src/lib/utils/ta_utils.h b/src/lib/utils/ta_utils.h
index 36ee551cc..866a4cc37 100644
--- a/src/lib/utils/ta_utils.h
+++ b/src/lib/utils/ta_utils.h
@@ -21,7 +21,16 @@ namespace TA_CM {
 * @param in an integer
 * @return 0 if in == 0 else 0xFFFFFFFF
 */
-u32bit gen_mask_u32bit(u32bit in);
+u32bit expand_mask_u32bit(u32bit in);
+
+
+/**
+ * Expand an input to a bit mask depending on it being being zero or
+ * non-zero
+ * @ param in the input
+ * @return the mask 0xFFFF if tst is non-zero and 0 otherwise
+ */
+u16bit expand_mask_u16bit(u16bit in);
 
 /**
 * Branch-free maximum