Move GF(p) math code from pk/ecdsa to math/gfpmath

1 files changed, 1157 insertions, 0 deletions

diff --git a/src/math/gfpmath/point_gfp.cpp b/src/math/gfpmath/point_gfp.cpp
new file mode 100644
index 000000000..23b6d4518
--- /dev/null
+++ b/src/math/gfpmath/point_gfp.cpp
@@ -0,0 +1,1157 @@
+/******************************************************
+ * Arithmetic for point groups of elliptic curves     *
+ * over GF(p) (source file)                           *
+ *                                                    *
+ * (C) 2007 Martin Döring                             *
+ *          [email protected]    *
+ *          Christoph Ludwig                          *
+ *          [email protected]                        *
+ *          Falko Strenzke                            *
+ *          [email protected]                    *
+ ******************************************************/
+
+#include <botan/point_gfp.h>
+#include <botan/numthry.h>
+
+namespace Botan {
+
+// construct the point at infinity or a random point
+PointGFp::PointGFp(CurveGFp const& curve)
+   :	mC(curve),
+        mX(curve.get_p(), 0),
+        mY(curve.get_p(), 1),
+        mZ(curve.get_p(), 0),
+        mZpow2(curve.get_p(),0),
+        mZpow3(curve.get_p(),0),
+        mAZpow4(curve.get_p(),0),
+        mZpow2_set(false),
+        mZpow3_set(false),
+        mAZpow4_set(false)
+   {
+   // first set the point wide pointer
+
+   set_shrd_mod(mC.get_ptr_mod());
+
+   }
+
+
+// construct a point given its jacobian projective coordinates
+PointGFp::PointGFp(CurveGFp const& curve, GFpElement const& x,
+                   GFpElement const& y, GFpElement const& z)
+   :	mC(curve),
+        mX(x),
+        mY(y),
+        mZ(z),
+        mZpow2(curve.get_p(),0),
+        mZpow3(curve.get_p(),0),
+        mAZpow4(curve.get_p(),0),
+        mZpow2_set(false),
+        mZpow3_set(false),
+        mAZpow4_set(false)
+   {
+   set_shrd_mod(mC.get_ptr_mod());
+   }
+PointGFp::PointGFp ( CurveGFp const& curve, GFpElement const& x,
+                     GFpElement const& y )
+   :mC(curve),
+    mX(x),
+    mY(y),
+    mZ(curve.get_p(),1),
+    mZpow2(curve.get_p(),0),
+    mZpow3(curve.get_p(),0),
+    mAZpow4(curve.get_p(),0),
+    mZpow2_set(false),
+    mZpow3_set(false),
+    mAZpow4_set(false)
+   {
+   set_shrd_mod(mC.get_ptr_mod());
+   }
+
+// copy constructor
+PointGFp::PointGFp(PointGFp const& other)
+   :	mC(other.mC),
+        mX(other.mX),
+        mY(other.mY),
+        mZ(other.mZ),
+        mZpow2(other.mZpow2),
+        mZpow3(other.mZpow3),
+        mAZpow4(other.mAZpow4),
+        mZpow2_set(other.mZpow2_set),
+        mZpow3_set(other.mZpow3_set),
+        mAZpow4_set(other.mAZpow4_set)
+   {
+   set_shrd_mod(mC.get_ptr_mod());
+   }
+
+// assignment operator
+PointGFp const& PointGFp::operator=(PointGFp const& other)
+   {
+   mC = other.get_curve();
+   mX = other.get_jac_proj_x();
+   mY = other.get_jac_proj_y();
+   mZ = other.get_jac_proj_z();
+   mZpow2 = GFpElement(other.mZpow2);
+   mZpow3 = GFpElement(other.mZpow3);
+   mAZpow4 = GFpElement(other.mAZpow4);
+   mZpow2_set = other.mZpow2_set;
+   mZpow3_set = other.mZpow3_set;
+   mAZpow4_set = other.mAZpow4_set;
+   set_shrd_mod(mC.get_ptr_mod());
+   return *this;
+   }
+
+PointGFp const& PointGFp::assign_within_same_curve(PointGFp const& other)
+   {
+   mX = other.get_jac_proj_x();
+   mY = other.get_jac_proj_y();
+   mZ = other.get_jac_proj_z();
+   mZpow2_set = false;
+   mZpow3_set = false;
+   mAZpow4_set = false;
+   // the rest stays!
+   return *this;
+   }
+
+void PointGFp::set_shrd_mod(std::tr1::shared_ptr<GFpModulus> p_mod)
+   {
+   mX.set_shrd_mod(p_mod);
+   mY.set_shrd_mod(p_mod);
+   mZ.set_shrd_mod(p_mod);
+   mZpow2.set_shrd_mod(p_mod);
+   mZpow3.set_shrd_mod(p_mod);
+   mAZpow4.set_shrd_mod(p_mod);
+   }
+
+void PointGFp::ensure_worksp() const
+   {
+   if (mp_worksp_gfp_el.get() != 0)
+      {
+      if ((*mp_worksp_gfp_el).size() == GFPEL_WKSP_SIZE)
+         {
+         return;
+         }
+      else
+         {
+         throw Invalid_State("encountered incorrect size for PointGFp´s GFpElement workspace");
+         }
+      }
+
+   mp_worksp_gfp_el = std::tr1::shared_ptr<std::vector<GFpElement> >(new std::vector<GFpElement>);
+   mp_worksp_gfp_el->reserve(9);
+   for (u32bit i=0; i<GFPEL_WKSP_SIZE; i++)
+      {
+      mp_worksp_gfp_el->push_back(GFpElement(1,0));
+
+      }
+   }
+
+// arithmetic operators
+PointGFp& PointGFp::operator+=(PointGFp const& rhs)
+   {
+   if (is_zero())
+      {
+      *this = rhs;
+      return *this;
+      }
+   if (rhs.is_zero())
+      {
+      return *this;
+      }
+   ensure_worksp();
+
+   if (rhs.mZ == *(mC.get_mres_one()))
+      {
+      //U1 = mX;
+      (*mp_worksp_gfp_el)[0].share_assign(mX);
+
+      //S1 = mY;
+      (*mp_worksp_gfp_el)[2].share_assign(mY);
+      }
+   else
+      {
+      if ((!rhs.mZpow2_set) || (!rhs.mZpow3_set))
+         {
+         rhs.mZpow2 = rhs.mZ;
+         rhs.mZpow2 *= rhs.mZ;
+         rhs.mZpow3 = rhs.mZpow2;
+         rhs.mZpow3 *= rhs.mZ;
+
+         rhs.mZpow2_set = true;
+         rhs.mZpow3_set = true;
+         }
+      //U1 = mX * rhs.mZpow2;
+      (*mp_worksp_gfp_el)[0].share_assign(mX);
+      (*mp_worksp_gfp_el)[0] *= rhs.mZpow2;
+
+      //S1 = mY * rhs.mZpow3;
+      (*mp_worksp_gfp_el)[2].share_assign(mY);
+      (*mp_worksp_gfp_el)[2] *= rhs.mZpow3;
+
+      }
+   if (mZ == *(mC.get_mres_one()))
+      {
+      //U2 = rhs.mX;
+      (*mp_worksp_gfp_el)[1].share_assign(rhs.mX);
+
+      //S2 = rhs.mY;
+      (*mp_worksp_gfp_el)[3].share_assign(rhs.mY);
+      }
+   else
+      {
+      if ((!mZpow2_set) || (!mZpow3_set))
+         {
+         // precomputation can´t be used, because *this changes anyway
+         mZpow2 = mZ;
+         mZpow2 *= mZ;
+
+         mZpow3 = mZpow2;
+         mZpow3 *= mZ;
+         }
+      //U2 = rhs.mX * mZpow2;
+      (*mp_worksp_gfp_el)[1].share_assign(rhs.mX);
+      (*mp_worksp_gfp_el)[1] *= mZpow2;
+
+      //S2 = rhs.mY * mZpow3;
+      (*mp_worksp_gfp_el)[3].share_assign(rhs.mY);
+      (*mp_worksp_gfp_el)[3] *= mZpow3;
+
+      }
+   //GFpElement H(U2 - U1);
+
+   (*mp_worksp_gfp_el)[4].share_assign((*mp_worksp_gfp_el)[1]);
+   (*mp_worksp_gfp_el)[4] -= (*mp_worksp_gfp_el)[0];
+
+   //GFpElement r(S2 - S1);
+   (*mp_worksp_gfp_el)[5].share_assign((*mp_worksp_gfp_el)[3]);
+   (*mp_worksp_gfp_el)[5] -= (*mp_worksp_gfp_el)[2];
+
+   //if(H.is_zero())
+   if ((*mp_worksp_gfp_el)[4].is_zero())
+
+      {
+      if ((*mp_worksp_gfp_el)[5].is_zero())
+
+         {
+         mult2_in_place();
+         return *this;
+         }
+      *this = PointGFp(mC); // setting myself to zero
+      return *this;
+      }
+
+   //U2 = H * H;
+   (*mp_worksp_gfp_el)[1].share_assign((*mp_worksp_gfp_el)[4]);
+   (*mp_worksp_gfp_el)[1] *= (*mp_worksp_gfp_el)[4];
+
+   //S2 = U2 * H;
+   (*mp_worksp_gfp_el)[3].share_assign((*mp_worksp_gfp_el)[1]);
+   (*mp_worksp_gfp_el)[3] *= (*mp_worksp_gfp_el)[4];
+
+   //U2 *= U1;
+   (*mp_worksp_gfp_el)[1] *= (*mp_worksp_gfp_el)[0];
+
+   //GFpElement x(r*r - S2 - (U2+U2));
+   (*mp_worksp_gfp_el)[6].share_assign((*mp_worksp_gfp_el)[5]);
+   (*mp_worksp_gfp_el)[6] *= (*mp_worksp_gfp_el)[5];
+   (*mp_worksp_gfp_el)[6] -= (*mp_worksp_gfp_el)[3];
+   (*mp_worksp_gfp_el)[6] -= (*mp_worksp_gfp_el)[1];
+   (*mp_worksp_gfp_el)[6] -= (*mp_worksp_gfp_el)[1];
+
+   //GFpElement z(S1 * S2);
+   (*mp_worksp_gfp_el)[8].share_assign((*mp_worksp_gfp_el)[2]);
+   (*mp_worksp_gfp_el)[8] *= (*mp_worksp_gfp_el)[3];
+
+   //GFpElement y(r * (U2-x) - z);
+   (*mp_worksp_gfp_el)[7].share_assign((*mp_worksp_gfp_el)[1]);
+   (*mp_worksp_gfp_el)[7] -= (*mp_worksp_gfp_el)[6];
+   (*mp_worksp_gfp_el)[7] *= (*mp_worksp_gfp_el)[5];
+   (*mp_worksp_gfp_el)[7] -= (*mp_worksp_gfp_el)[8];
+
+   if (mZ == *(mC.get_mres_one()))
+      {
+      if (rhs.mZ != *(mC.get_mres_one()))
+         {
+         //z = rhs.mZ * H;
+         (*mp_worksp_gfp_el)[8].share_assign(rhs.mZ);
+         (*mp_worksp_gfp_el)[8] *= (*mp_worksp_gfp_el)[4];
+         }
+      else
+         {
+         //z = H;
+         (*mp_worksp_gfp_el)[8].share_assign((*mp_worksp_gfp_el)[4]);
+         }
+      }
+   else if (rhs.mZ != *(mC.get_mres_one()))
+      {
+      //U1 = mZ * rhs.mZ;
+      (*mp_worksp_gfp_el)[0].share_assign(mZ);
+      (*mp_worksp_gfp_el)[0] *= rhs.mZ;
+
+      //z = U1 * H;
+      (*mp_worksp_gfp_el)[8].share_assign((*mp_worksp_gfp_el)[0]);
+      (*mp_worksp_gfp_el)[8] *= (*mp_worksp_gfp_el)[4];
+
+      }
+   else
+      {
+      //z = mZ * H;
+      (*mp_worksp_gfp_el)[8].share_assign(mZ);
+      (*mp_worksp_gfp_el)[8] *= (*mp_worksp_gfp_el)[4];
+
+      }
+   mZpow2_set = false;
+   mZpow3_set = false;
+   mAZpow4_set = false;
+
+   mX = (*mp_worksp_gfp_el)[6];
+   mY = (*mp_worksp_gfp_el)[7];
+   mZ = (*mp_worksp_gfp_el)[8];
+
+   return *this;
+
+   }
+PointGFp& PointGFp::operator-=(PointGFp const& rhs)
+   {
+   PointGFp minus_rhs = PointGFp(rhs).negate();
+
+   if (is_zero())
+      {
+      *this = minus_rhs;
+      }
+   else
+      {
+      *this += minus_rhs;
+      }
+   return *this;
+   }
+
+PointGFp& PointGFp::mult_this_secure(const BigInt& scalar,
+                                     const BigInt& /*point_order*/,
+                                     const BigInt& /*max_secr*/)
+   {
+   // NOTE: FS: so far this is code duplication of op*=.
+   // we have to see how we deal with this.
+   // fact is that we will probably modify this function
+   // while evaluating the countermeasures
+   // whereas we probably will not start modifying the
+   // function operator*=.
+   // however, in the end both should be merged.
+
+   // use montgomery mult. in this operation
+   this->turn_on_sp_red_mul();
+
+   std::tr1::shared_ptr<PointGFp> H(new PointGFp(this->mC));
+   std::tr1::shared_ptr<PointGFp> tmp; // used for AADA
+
+   PointGFp P(*this);
+   BigInt m(scalar);
+
+   if (m < BigInt(0))
+      {
+      m = -m;
+      P.negate();
+      }
+   if (P.is_zero() || (m == BigInt(0)))
+      {
+      *this = *H;
+      return *this;
+      }
+   if (m == BigInt(1))
+      {
+      return *this;
+      }
+   //
+#ifdef CM_AADA
+#ifndef CM_RAND_EXP
+   int max_secr_bits = max_secr.bits();
+#endif
+#endif
+
+   int mul_bits = m.bits(); // this is used for a determined number of loop runs in
+   // the mult_loop where leading zero´s are padded if necessary.
+   // Here we assign the value that will be used when no countermeasures are specified
+#ifdef CM_RAND_EXP
+   u32bit rand_r_bit_len = 20; // Coron(99) proposes 20 bit for r
+
+#ifdef CM_AADA
+
+   BigInt r_max(1);
+
+#endif // CM_AADA
+
+   // use randomized exponent
+#ifdef TA_COLL_T
+   static BigInt r_randexp;
+   if (new_rand)
+      {
+      r_randexp = random_integer(rand_r_bit_len);
+      }
+   //assert(!r_randexp.is_zero());
+#else
+   BigInt r_randexp(random_integer(rand_r_bit_len));
+#endif
+
+   m += r_randexp * point_order;
+   // determine mul_bits...
+#ifdef CM_AADA
+   // AADA with rand. Exp.
+   //assert(rand_r_bit_len > 0);
+   r_max <<= rand_r_bit_len;
+   r_max -= 1;
+   //assert(r_max.bits() == rand_r_bit_len);
+   mul_bits = (max_secr + point_order * r_max).bits();
+#else
+   // rand. Exp. without AADA
+   mul_bits = m.bits();
+#endif // CM_AADA
+
+
+#endif // CM_RAND_EXP
+
+   // determine mul_bits...
+#if (CM_AADA == 1 && CM_RAND_EXP != 1)
+
+   mul_bits = max_secr_bits;
+#endif // CM_AADA without CM_RAND_EXP
+
+   //assert(mul_bits != 0);
+
+
+   H = mult_loop(mul_bits-1, m, H, tmp, P);
+
+   if (!H->is_zero()) // cannot convert if H == O
+      {
+      *this = H->get_z_to_one();
+      }else
+      {
+      *this = *H;
+      }
+   mX.turn_off_sp_red_mul();
+   mY.turn_off_sp_red_mul();
+   mZ.turn_off_sp_red_mul();
+   return *this;
+   }
+PointGFp& PointGFp::operator*=(const BigInt& scalar)
+   {
+   // use montgomery mult. in this operation
+
+   this->turn_on_sp_red_mul();
+
+   PointGFp H(this->mC); // create as zero
+   H.turn_on_sp_red_mul();
+   PointGFp P(*this);
+   P.turn_on_sp_red_mul();
+   BigInt m(scalar);
+   if (m < BigInt(0))
+      {
+      m = -m;
+      P.negate();
+      }
+   if (P.is_zero() || (m == BigInt(0)))
+      {
+      *this = H;
+      return *this;
+      }
+   if (m == BigInt(1))
+      {
+      //*this == P already
+      return *this;
+      }
+
+   const int l = m.bits() - 1;
+   for (int i=l; i >=0; i--)
+      {
+
+      H.mult2_in_place();
+      if (m.get_bit(i))
+         {
+         H += P;
+         }
+      }
+
+   if (!H.is_zero()) // cannot convert if H == O
+      {
+      *this = H.get_z_to_one();
+      }else
+      {
+      *this = H;
+      }
+   return *this;
+   }
+
+inline std::tr1::shared_ptr<PointGFp> PointGFp::mult_loop(
+   int l,
+   const BigInt& m,
+   std::tr1::shared_ptr<PointGFp> H,
+   std::tr1::shared_ptr<PointGFp> tmp,
+   PointGFp const& P)
+   {
+
+   //assert(l >= (int)m.bits()- 1);
+   tmp = H;
+   std::tr1::shared_ptr<PointGFp> to_add(new PointGFp(P)); // we just need some point
+   // so that we can use op=
+   // inside the loop
+   for (int i=l; i >=0; i--)
+      {
+      H->mult2_in_place();
+
+#ifndef CM_AADA
+
+      if (m.get_bit(i))
+         {
+         *H += P;
+         }
+#else // (CM_AADA is in)
+
+      if (H.get() == to_add.get())
+         {
+         to_add = tmp; // otherwise all pointers might point to the same object
+         // and we always need two objects to be able to switch around
+         }
+      to_add->assign_within_same_curve(*H);
+      tmp = H;
+      *tmp += P; // tmp already points to H
+
+      if (m.get_bit(i))
+         {
+         H = tmp; // NOTE: assign the pointer, not the value!
+         // (so that the operation is fast and thus as difficult
+         // to detect as possible)
+         }
+      else
+         {
+         H = to_add; // NOTE: this is necessary, because the assignment
+         // "*tmp = ..." already changed what H pointed to
+
+
+         }
+#endif // CM_AADA
+
+      }
+   return H;
+   }
+PointGFp& PointGFp::negate()
+   {
+   if (!is_zero())
+      {
+      mY.negate();
+      }
+   return *this;
+   }
+// *this *= 2
+PointGFp& PointGFp::mult2_in_place()
+   {
+   if (is_zero())
+      {
+      return *this;
+      }
+   if (mY.is_zero())
+      {
+
+      *this = PointGFp(mC); // setting myself to zero
+      return *this;
+      }
+   ensure_worksp();
+
+   (*mp_worksp_gfp_el)[0].share_assign(mY);
+   (*mp_worksp_gfp_el)[0] *= mY;
+
+   //GFpElement S(mX * z);
+   (*mp_worksp_gfp_el)[1].share_assign(mX);
+   (*mp_worksp_gfp_el)[1] *= (*mp_worksp_gfp_el)[0];
+
+   //GFpElement x(S + S);
+   (*mp_worksp_gfp_el)[2].share_assign((*mp_worksp_gfp_el)[1]);
+   (*mp_worksp_gfp_el)[2] += (*mp_worksp_gfp_el)[1];
+
+   //S = x + x;
+   (*mp_worksp_gfp_el)[1].share_assign((*mp_worksp_gfp_el)[2]);
+   (*mp_worksp_gfp_el)[1] += (*mp_worksp_gfp_el)[2];
+
+   if (!mAZpow4_set)
+      {
+      if (mZ == *(mC.get_mres_one()))
+         {
+         mAZpow4 = mC.get_mres_a();
+         mAZpow4_set = true;
+         }
+      else
+         {
+         if (!mZpow2_set)
+            {
+            mZpow2 = mZ;
+            mZpow2 *= mZ;
+
+            mZpow2_set = true;
+            }
+         //x = mZpow2 * mZpow2;
+         (*mp_worksp_gfp_el)[2].share_assign(mZpow2);
+         (*mp_worksp_gfp_el)[2] *= mZpow2;
+
+         //mAZpow4 = mC.get_mres_a() * x;
+         mAZpow4 = mC.get_mres_a();
+         mAZpow4 *= (*mp_worksp_gfp_el)[2];
+
+         }
+
+      }
+
+   //GFpElement y(mX * mX);
+   (*mp_worksp_gfp_el)[3].share_assign(mX);
+   (*mp_worksp_gfp_el)[3] *= mX;
+
+   //GFpElement M(y + y + y + mAZpow4);
+   (*mp_worksp_gfp_el)[4].share_assign((*mp_worksp_gfp_el)[3]);
+   (*mp_worksp_gfp_el)[4] += (*mp_worksp_gfp_el)[3];
+   (*mp_worksp_gfp_el)[4] += (*mp_worksp_gfp_el)[3];
+   (*mp_worksp_gfp_el)[4] += mAZpow4;
+
+   //x = M * M - (S+S);
+   (*mp_worksp_gfp_el)[2].share_assign((*mp_worksp_gfp_el)[4]);
+   (*mp_worksp_gfp_el)[2] *= (*mp_worksp_gfp_el)[4];
+   (*mp_worksp_gfp_el)[2] -= (*mp_worksp_gfp_el)[1];
+   (*mp_worksp_gfp_el)[2] -= (*mp_worksp_gfp_el)[1];
+
+   //y = z * z;
+   (*mp_worksp_gfp_el)[3].share_assign((*mp_worksp_gfp_el)[0]);
+   (*mp_worksp_gfp_el)[3] *= (*mp_worksp_gfp_el)[0];
+
+   //GFpElement U(y + y);
+   (*mp_worksp_gfp_el)[5].share_assign((*mp_worksp_gfp_el)[3]);
+   (*mp_worksp_gfp_el)[5] += (*mp_worksp_gfp_el)[3];
+
+   //z = U + U;
+   (*mp_worksp_gfp_el)[0].share_assign((*mp_worksp_gfp_el)[5]);
+   (*mp_worksp_gfp_el)[0] += (*mp_worksp_gfp_el)[5];
+
+   //U = z + z;
+   (*mp_worksp_gfp_el)[5].share_assign((*mp_worksp_gfp_el)[0]);
+   (*mp_worksp_gfp_el)[5] += (*mp_worksp_gfp_el)[0];
+
+   //y = M * (S - x) - U;
+   (*mp_worksp_gfp_el)[3].share_assign((*mp_worksp_gfp_el)[1]);
+   (*mp_worksp_gfp_el)[3] -= (*mp_worksp_gfp_el)[2];
+   (*mp_worksp_gfp_el)[3] *= (*mp_worksp_gfp_el)[4];
+   (*mp_worksp_gfp_el)[3] -= (*mp_worksp_gfp_el)[5];
+
+   if (mZ != *(mC.get_mres_one()))
+      {
+      //z = mY * mZ;
+      (*mp_worksp_gfp_el)[0].share_assign(mY);
+      (*mp_worksp_gfp_el)[0] *= mZ;
+
+      }
+   else
+      {
+      //z = mY;
+      (*mp_worksp_gfp_el)[0].share_assign(mY);
+
+      }
+   //z = z + z;
+   (*mp_worksp_gfp_el)[6].share_assign((*mp_worksp_gfp_el)[0]);
+   (*mp_worksp_gfp_el)[0] += (*mp_worksp_gfp_el)[6];
+
+   //mX = x;
+   //mY = y;
+   //mZ = z;
+   mX = (*mp_worksp_gfp_el)[2];
+   mY = (*mp_worksp_gfp_el)[3];
+   mZ = (*mp_worksp_gfp_el)[0];
+
+   mZpow2_set = false;
+   mZpow3_set = false;
+   mAZpow4_set = false;
+   return *this;
+
+   }
+void PointGFp::turn_on_sp_red_mul() const
+   {
+   mX.turn_on_sp_red_mul();
+   mY.turn_on_sp_red_mul();
+   mZ.turn_on_sp_red_mul();
+
+   // also pretransform, otherwise
+   // we might have bad results with respect to
+   // performance because
+   // additions/subtractions in mult2_in_place()
+   // and op+= spread untransformed GFpElements
+   mX.get_mres();
+   mY.get_mres();
+   mZ.get_mres();
+
+   mZpow2.turn_on_sp_red_mul();
+   mZpow3.turn_on_sp_red_mul();
+   mAZpow4.turn_on_sp_red_mul();
+   }
+// getters
+
+/**
+* returns a point equivalent to *this but were
+* Z has value one, i.e. x and y correspond to
+* their values in affine coordinates
+*/
+PointGFp const PointGFp::get_z_to_one() const
+   {
+   return PointGFp(*this).set_z_to_one();
+   }
+/**
+* changes the representation of *this so that
+* Z has value one, i.e. x and y correspond to
+* their values in affine coordinates.
+* returns *this.
+*/
+PointGFp const& PointGFp::set_z_to_one() const
+   {
+   if (!(mZ.get_value() == BigInt(1)) && !(mZ.get_value() == BigInt(0)))
+      {
+      GFpElement z = inverse(mZ);
+      GFpElement z2 = z * z;
+      z *= z2;
+      GFpElement x = mX * z2;
+      GFpElement y = mY * z;
+      mZ = GFpElement(mC.get_p(), BigInt(1));
+      mX = x;
+      mY = y;
+      }
+   else
+      {
+      if (mZ.get_value() == BigInt(0))
+         {
+         throw Illegal_Transformation("cannot convert Z to one");
+         }
+      }
+   return *this; // mZ = 1 already
+   }
+CurveGFp const PointGFp::get_curve() const
+   {
+   return mC;
+   }
+GFpElement const PointGFp::get_affine_x() const
+   {
+
+   if (is_zero())
+      {
+      throw Illegal_Transformation("cannot convert to affine");
+
+      }
+   /*if(!mZpow2_set)
+   {*/
+   mZpow2 = mZ * mZ;
+   mZpow2_set = true;
+   //}
+   //assert(mZpow2 == mZ*mZ);
+   GFpElement z2 = mZpow2;
+   return mX * z2.inverse_in_place();
+   }
+
+GFpElement const PointGFp::get_affine_y() const
+   {
+
+   if (is_zero())
+      {
+      throw Illegal_Transformation("cannot convert to affine");
+
+      }
+   /*if(!mZpow3_set )
+   {*/
+   mZpow3 = mZ * mZ * mZ;
+   mZpow3_set = true;
+   //}
+   //assert(mZpow3 == mZ * mZ *mZ);
+   GFpElement z3 = mZpow3;
+   return mY * z3.inverse_in_place();
+   }
+GFpElement const PointGFp::get_jac_proj_x() const
+   {
+   return GFpElement(mX);
+   }
+GFpElement const PointGFp::get_jac_proj_y() const
+   {
+   return GFpElement(mY);
+   }
+GFpElement const PointGFp::get_jac_proj_z() const
+   {
+   return GFpElement(mZ);
+   }
+
+// Is this the point at infinity?
+bool PointGFp::is_zero() const
+   {
+   return(mX.is_zero() && mZ.is_zero());
+   //NOTE: the calls to GFpElement::is_zero() instead of getting the value and
+   // and comparing it are import because they do not provoke backtransformations
+   // to the ordinary residue.
+   }
+
+// Is the point still on the curve??
+// (If everything is correct, the point is always on its curve; then the
+// function will return silently. If Oskar managed to corrupt this object's state,
+// then it will throw an exception.)
+void PointGFp::check_invariants() const
+   {
+   if (is_zero())
+      {
+      return;
+      }
+   const GFpElement y2 = mY * mY;
+   const GFpElement x3 = mX * mX * mX;
+
+   if (mZ.get_value() == BigInt(1))
+      {
+      GFpElement ax = mC.get_a() * mX;
+      if (y2 != (x3 + ax + mC.get_b()))
+         {
+         throw Illegal_Point();
+         }
+
+      }
+   /*if (!mZpow2_set)
+   {*/
+   mZpow2 = mZ * mZ;
+   mZpow2_set = true;
+   /*}
+   if (!mZpow3_set)
+   {*/
+   mZpow3 = mZpow2 * mZ;
+   mZpow3_set = true;
+   /*}
+   if(!mAZpow4_set)
+   {*/
+   mAZpow4 = mZpow3 * mZ * mC.get_a();
+   mAZpow4_set = true;
+   //}
+   const GFpElement aXZ4 = mAZpow4 * mX;
+   const GFpElement bZ6 = mC.get_b() * mZpow3 * mZpow3;
+
+   if (y2 != (x3 + aXZ4 + bZ6))
+      {
+      throw Illegal_Point();
+      }
+
+   }
+// swaps the states of *this and other, does not throw!
+void PointGFp::swap(PointGFp& other)
+   {
+   mC.swap(other.mC);
+   mX.swap(other.mX);
+   mY.swap(other.mY);
+   mZ.swap(other.mZ);
+   mZpow2.swap(other.mZpow2);
+   mZpow3.swap(other.mZpow3);
+   mAZpow4.swap(other.mAZpow4);
+   std::swap<bool>(mZpow2_set, other.mZpow2_set);
+   std::swap<bool>(mZpow3_set, other.mZpow3_set);
+   std::swap<bool>(mAZpow4_set, other.mAZpow4_set);
+   }
+
+PointGFp const mult2(PointGFp const& point)
+   {
+   return (PointGFp(point)).mult2_in_place();
+   }
+
+
+bool operator==(PointGFp const& lhs, PointGFp const& rhs)
+   {
+   if (lhs.is_zero() && rhs.is_zero())
+      {
+      return true;
+      }
+   if ((lhs.is_zero() && !rhs.is_zero()) || (!lhs.is_zero() && rhs.is_zero()))
+      {
+      return false;
+      }
+   // neither operand is zero, so we can call get_z_to_one()
+   //assert(!lhs.is_zero());
+   //assert(!rhs.is_zero());
+   PointGFp aff_lhs = lhs.get_z_to_one();
+   PointGFp aff_rhs = rhs.get_z_to_one();
+   return (aff_lhs.get_curve() == aff_rhs.get_curve() &&
+           aff_lhs.get_jac_proj_x() == aff_rhs.get_jac_proj_x() &&
+           aff_lhs.get_jac_proj_y() == aff_rhs.get_jac_proj_y());
+   }
+
+// arithmetic operators
+PointGFp operator+(PointGFp const& lhs, PointGFp const& rhs)
+   {
+   PointGFp tmp(lhs);
+   return tmp += rhs;
+   }
+
+PointGFp operator-(PointGFp const& lhs, PointGFp const& rhs)
+   {
+   PointGFp tmp(lhs);
+   return tmp -= rhs;
+   }
+
+PointGFp operator-(PointGFp const& lhs)
+   {
+   return PointGFp(lhs).negate();
+   }
+
+PointGFp operator*(const BigInt& scalar, PointGFp const& point)
+   {
+   PointGFp result(point);
+   return result *= scalar;
+   }
+
+PointGFp operator*(PointGFp const& point, const BigInt& scalar)
+   {
+   PointGFp result(point);
+   return result *= scalar;
+   }
+
+PointGFp mult_point_secure(PointGFp const& point, const BigInt& scalar, const BigInt& point_order, const BigInt& max_secret)
+   {
+   PointGFp result(point);
+   result.mult_this_secure(scalar, point_order, max_secret);
+   return result;
+   }
+
+// encoding and decoding
+SecureVector<byte> EC2OSP(PointGFp const& point, byte format)
+   {
+   SecureVector<byte> result;
+   if (format == PointGFp::UNCOMPRESSED)
+      {
+      result = encode_uncompressed(point);
+      }
+   else if (format == PointGFp::COMPRESSED)
+      {
+      result = encode_compressed(point);
+
+      }
+   else if (format == PointGFp::HYBRID)
+      {
+      result = encode_hybrid(point);
+      }
+   else
+      {
+      throw Format_Error("illegal point encoding format specification");
+      }
+   return result;
+   }
+SecureVector<byte> encode_compressed(PointGFp const& point)
+   {
+
+
+   if (point.is_zero())
+      {
+      SecureVector<byte> result (1);
+      result[0] = 0;
+      return result;
+
+      }
+   u32bit l = point.get_curve().get_p().bits();
+   int dummy = l & 7;
+   if (dummy != 0)
+      {
+      l += 8 - dummy;
+      }
+   l /= 8;
+   SecureVector<byte> result (l+1);
+   result[0] = 2;
+   BigInt x = point.get_affine_x().get_value();
+   SecureVector<byte> bX = BigInt::encode_1363(x, l);
+   result.copy(1, bX.begin(), bX.size());
+   BigInt y = point.get_affine_y().get_value();
+   if (y.get_bit(0))
+      {
+      result[0] |= 1;
+      }
+   return result;
+   }
+
+
+SecureVector<byte> encode_uncompressed(PointGFp const& point)
+   {
+   if (point.is_zero())
+      {
+      SecureVector<byte> result (1);
+      result[0] = 0;
+      return result;
+      }
+   u32bit l = point.get_curve().get_p().bits();
+   int dummy = l & 7;
+   if (dummy != 0)
+      {
+      l += 8 - dummy;
+      }
+   l /= 8;
+   SecureVector<byte> result (2*l+1);
+   result[0] = 4;
+   BigInt x = point.get_affine_x().get_value();
+   BigInt y = point.get_affine_y().get_value();
+   SecureVector<byte> bX = BigInt::encode_1363(x, l);
+   SecureVector<byte> bY = BigInt::encode_1363(y, l);
+   result.copy(1, bX.begin(), l);
+   result.copy(l+1, bY.begin(), l);
+   return result;
+
+   }
+
+SecureVector<byte> encode_hybrid(PointGFp const& point)
+   {
+   if (point.is_zero())
+      {
+      SecureVector<byte> result (1);
+      result[0] = 0;
+      return result;
+      }
+   u32bit l = point.get_curve().get_p().bits();
+   int dummy = l & 7;
+   if (dummy != 0)
+      {
+      l += 8 - dummy;
+      }
+   l /= 8;
+   SecureVector<byte> result (2*l+1);
+   result[0] = 6;
+   BigInt x = point.get_affine_x().get_value();
+   BigInt y = point.get_affine_y().get_value();
+   SecureVector<byte> bX = BigInt::encode_1363(x, l);
+   SecureVector<byte> bY = BigInt::encode_1363(y, l);
+   result.copy(1, bX.begin(), bX.size());
+   result.copy(l+1, bY.begin(), bY.size());
+   if (y.get_bit(0))
+      {
+      result[0] |= 1;
+      }
+   return result;
+   }
+
+PointGFp OS2ECP(MemoryRegion<byte> const& os, CurveGFp const& curve)
+   {
+   if (os.size() == 1 && os[0] == 0)
+      {
+      return PointGFp(curve); // return zero
+      }
+   SecureVector<byte> bX;
+   SecureVector<byte> bY;
+
+   GFpElement x(1,0);
+   GFpElement y(1,0);
+   GFpElement z(1,0);
+
+   const byte pc = os[0];
+   BigInt bi_dec_x;
+   BigInt bi_dec_y;
+   switch (pc)
+      {
+      case 2:
+      case 3:
+         //compressed form
+         bX = SecureVector<byte>(os.size() - 1);
+         bX.copy(os.begin()+1, os.size()-1);
+
+         /* Problem wäre, wenn decode() das erste bit als Vorzeichen interpretiert.
+         *---------------------
+         * AW(FS): decode() interpretiert das erste Bit nicht als Vorzeichen
+         */
+         bi_dec_x = BigInt::decode(bX, bX.size());
+         x = GFpElement(curve.get_p(), bi_dec_x);
+         bool yMod2;
+         yMod2 = (pc & 1) == 1;
+         y = PointGFp::decompress(yMod2, x, curve);
+         break;
+      case 4:
+         // uncompressed form
+         int l;
+         l = (os.size() -1)/2;
+         bX = SecureVector<byte>(l);
+         bY = SecureVector<byte>(l);
+         bX.copy(os.begin()+1, l);
+         bY.copy(os.begin()+1+l, l);
+         bi_dec_x = BigInt::decode(bX.begin(), bX.size());
+
+         bi_dec_y = BigInt::decode(bY.begin(),bY.size());
+         x = GFpElement(curve.get_p(), bi_dec_x);
+         y = GFpElement(curve.get_p(), bi_dec_y);
+         break;
+
+      case 6:
+      case 7:
+         //hybrid form
+         l = (os.size() - 1)/2;
+         bX = SecureVector<byte>(l);
+         bY = SecureVector<byte>(l);
+         bX.copy(os.begin() + 1, l);
+         bY.copy(os.begin()+1+l, l);
+         yMod2 = (pc & 0x01) == 1;
+         if (!(PointGFp::decompress(yMod2, x, curve) == y))
+            {
+            throw Illegal_Point("error during decoding hybrid format");
+            }
+         break;
+      default:
+         throw Format_Error("encountered illegal format specification while decoding point");
+      }
+   z = GFpElement(curve.get_p(), BigInt(1));
+   //assert((x.is_trf_to_mres() && x.is_use_montgm()) || !x.is_trf_to_mres());
+   //assert((y.is_trf_to_mres() && y.is_use_montgm()) || !y.is_trf_to_mres());
+   //assert((z.is_trf_to_mres() && z.is_use_montgm()) || !z.is_trf_to_mres());
+   PointGFp result(curve, x, y, z);
+   result.check_invariants();
+   //assert((result.get_jac_proj_x().is_trf_to_mres() && result.get_jac_proj_x().is_use_montgm()) || !result.get_jac_proj_x().is_trf_to_mres());
+   //assert((result.get_jac_proj_y().is_trf_to_mres() && result.get_jac_proj_y().is_use_montgm()) || !result.get_jac_proj_y().is_trf_to_mres());
+   //assert((result.get_jac_proj_z().is_trf_to_mres() && result.get_jac_proj_z().is_use_montgm()) || !result.get_jac_proj_z().is_trf_to_mres());
+   return result;
+   }
+
+GFpElement PointGFp::decompress(bool yMod2, GFpElement const& x,
+                                CurveGFp const& curve)
+   {
+   BigInt xVal = x.get_value();
+   BigInt xpow3 = xVal * xVal * xVal;
+   BigInt g = curve.get_a().get_value() * xVal;
+   g += xpow3;
+   g += curve.get_b().get_value();
+   g = g%curve.get_p();
+   BigInt z = ressol(g, curve.get_p());
+
+   if(z < 0)
+      throw Illegal_Point("error during decompression");
+
+   bool zMod2 = z.get_bit(0);
+   if ((zMod2 && ! yMod2) || (!zMod2 && yMod2))
+      {
+      z = curve.get_p() - z;
+      }
+   return GFpElement(curve.get_p(),z);
+   }
+
+PointGFp const create_random_point(RandomNumberGenerator& rng,
+                                   CurveGFp const& curve)
+   {
+
+   // create a random point
+   GFpElement mX(1,1);
+   GFpElement mY(1,1);
+   GFpElement mZ(1,1);
+   GFpElement minusOne(curve.get_p(), BigInt(BigInt::Negative,1));
+   mY = minusOne;
+   GFpElement y2(1,1);
+   GFpElement x(1,1);
+
+   while (mY == minusOne)
+      {
+      BigInt value(rng, curve.get_p().bits());
+      mX = GFpElement(curve.get_p(),value);
+      y2 = curve.get_a() * mX;
+      x = mX * mX;
+      x *= mX;
+      y2 += (x + curve.get_b());
+
+      value = ressol(y2.get_value(), curve.get_p());
+
+      if(value < 0)
+         mY = minusOne;
+      else
+         mY = GFpElement(curve.get_p(), value);
+      }
+   mZ = GFpElement(curve.get_p(), BigInt(1));
+
+   return PointGFp(curve, mX, mY, mZ);
+   }
+
+} // namespace Botan