algebraic decoding

include/codecrypt.h

@@ -26,6 +26,8 @@ typedef int sint;
  * vector over GF(2). We rely on STL's vector<bool> == bit_vector
  * specialization for efficiency.
  */
+class polynomial;
+class gf2m;
 class bvector : public std::vector<bool>
 {
 protected:
@@ -34,6 +36,8 @@ public:
 	uint hamming_weight();
 	void add (const bvector&);
 	bool operator* (const bvector&); //dot product
+	bool zero() const;
+	void to_poly (polynomial&, gf2m&);
 };
 
 /*
@@ -115,6 +119,7 @@ public:
 	uint mult (uint, uint);
 	uint exp (uint, sint);
 	uint inv (uint);
+	uint sq_root (uint);
 };
 
 /*
@@ -129,13 +134,26 @@ public:
 	void strip();
 	int degree() const;
 	bool zero() const;
+	void shift (uint);
+
 	uint eval (uint, gf2m&) const;
 	void add (const polynomial&, gf2m&);
-	void mod (const polynomial&, gf2m&);
 	void mult (const polynomial&, gf2m&);
+	void add_mult (const polynomial&, uint mult, gf2m&);
+	void mod (const polynomial&, gf2m&);
+	void div (polynomial&, polynomial&, gf2m&);
+	void divmod (polynomial&, polynomial&, polynomial&, gf2m&);
+	void square (gf2m&);
+	void inv (polynomial&, gf2m&);
+	void make_monic (gf2m&);
+
+	void sqrt (vector<polynomial>&, gf2m&);
 	polynomial gcd (polynomial, gf2m&);
+	void mod_to_fracton (polynomial&, polynomial&, polynomial&, gf2m&);
+
 	bool is_irreducible (gf2m&) const;
 	void generate_random_irreducible (uint s, gf2m&, prng&);
+
 	bool compute_square_root_matrix (std::vector<polynomial>&, gf2m&);
 	void compute_goppa_check_matrix (matrix&, gf2m&);
 };

lib/bvector.cpp

@@ -25,3 +25,16 @@ bool bvector::operator* (const bvector&a)
 	return r;
 }
 
+bool bvector::zero() const
+{
+	for (uint i = 0; i < size(); ++i) if (item (i) ) return false;
+	return true;
+}
+
+void bvector::to_poly (polynomial&r, gf2m&fld)
+{
+	r.clear();
+	r.resize ( (size() % fld.m) ? 1 : 0 + (size() / fld.m), 0);
+	for (uint i = 0; i < size(); ++i)
+		if (item (i) ) r[i/fld.m] |= 1 << (i % fld.m);
+}

lib/decoding.cpp

@@ -0,0 +1,30 @@
+
+#include "decoding.h"
+
+void syndrome_decode (bvector&syndrome, gf2m&fld, polynomial& goppa,
+                      std::vector<polynomial>& sqInv, bvector&ev)
+
+{
+	ev.resize (fld.n, 0);
+	if (syndrome.zero() ) return;
+
+	polynomial p;
+	syndrome.to_poly (p, fld);
+	p.inv (goppa, fld); // p=S^-1 mod goppa
+
+	p[1] = fld.add (1, p[1]); //p is now tau
+	p.sqrt (sqInv, fld); //tau = sqrt(T+x) mod goppa
+
+	polynomial a, b;
+	p.mod_to_fracton (a, b, goppa, fld);
+	a.square (fld);
+	b.square (fld);
+	b.shift (1);
+	a.add (b, fld); //new a = a^2 + x b^2
+
+	a.make_monic (fld); //now it is the error locator.
+
+	for (uint i = 0; i < fld.n; ++i) {
+		if (0 == a.eval (i, fld) ) ev[i] = 1;
+	}
+}

lib/decoding.h

@@ -0,0 +1,14 @@
+
+#ifndef _decoding_h_
+#define _decoding_h_
+
+#include "codecrypt.h"
+
+using namespace ccr;
+void syndrome_decode (bvector&syndrome,
+                      gf2m&fld,
+                      polynomial& gp,
+                      std::vector<polynomial>& sqInv,
+                      bvector&ev);
+
+#endif

lib/gf2m.cpp

@@ -117,3 +117,9 @@ uint gf2m::inv (uint a)
 	return exp (a, n - 2);
 }
 
+uint gf2m::sq_root (uint a)
+{
+	for (uint i = 0; i < m; ++i)
+		a = mult (a, a);
+	return a;
+}

lib/matrix.cpp

@@ -62,8 +62,8 @@ bool matrix::compute_inversion (matrix&res)
 		if (m[i][i] != 1) {
 			for (j = i + 1; j < s; ++j) if (m[j][i] == 1) break;
 			if (j == s) return false; //noninvertible
-			m[i].add (m[j]);
-			r[i].add (r[j]);
+			m[i].swap (m[j]);
+			r[i].swap (r[j]);
 		}
 		//remove 1's below
 		for (j = i + 1; j < s; ++j) if (m[j][i]) {

lib/mce.cpp

@@ -69,7 +69,6 @@ int pubkey::encrypt (const bvector& in, bvector&out, prng&rng)
 
 int privkey::decrypt (const bvector&in, bvector&out)
 {
-
 	return -1; //TODO
 }
 

lib/polynomial.cpp

@@ -28,6 +28,14 @@ void polynomial::add (const polynomial&f, gf2m&fld)
 	for (int i = 0; i <= df; ++i) item (i) = fld.add (item (i), f[i]);
 }
 
+void polynomial::add_mult (const polynomial&f, uint mult, gf2m&fld)
+{
+	int df = f.degree();
+	if (df > degree() ) resize (df + 1);
+	for (int i = 0; i <= df; ++i)
+		item (i) = fld.add (item (i), fld.mult (mult, f[i]) );
+}
+
 void polynomial::mod (const polynomial&f, gf2m&fld)
 {
 	int df = f.degree();
@@ -226,3 +234,103 @@ void polynomial::compute_goppa_check_matrix (matrix&r, gf2m&fld)
 			r[i][j] = (h[j/fld.m][i] >> (j % fld.m) ) & 1;
 	}
 }
+
+void polynomial::make_monic (gf2m&fld)
+{
+	int d = degree();
+	if (d < 0) return;
+	uint m = fld.inv (item (d) );
+	for (uint i = 0; i <= d; ++i) item (i) = fld.mult (item (i), m);
+}
+
+void polynomial::shift (uint n)
+{
+	if (degree() < 0) return;
+	insert (begin(), n, 0);
+}
+
+void polynomial::square (gf2m&fld)
+{
+	polynomial a = *this;
+	this->mult (a, fld);
+}
+
+void polynomial::sqrt (vector<polynomial>& sqInv, gf2m&fld)
+{
+	polynomial a = *this;
+	clear();
+	for (uint i = 0; i < a.size(); ++i) add_mult (sqInv[i], a[i], fld);
+	for (uint i = 0; i < size(); ++i) item (i) = fld.sq_root (item (i) );
+}
+
+void polynomial::div (polynomial&p, polynomial&m, gf2m&fld)
+{
+	int degp = p.degree();
+	if (degp < 0) return;
+
+	uint headInv = fld.inv (p[degp]);
+	polynomial A = *this;
+	A.mod (m, fld);
+	clear();
+	int da;
+	while ( (da = A.degree() ) >= degp) {
+		int rp = da - degp;
+		if (size() < rp + 1) resize (rp + 1, 0);
+		item (rp) = fld.mult (headInv, A[da]);
+		for (uint i = 0; i <= degp; ++i)
+			A[i+rp] = fld.add (A[i+rp], fld.mult (item (rp), p[i]) );
+	}
+}
+
+void polynomial::divmod (polynomial&d, polynomial&res, polynomial&rem, gf2m&fld)
+{
+	int degd = d.degree();
+	if (degd < 0) return;
+
+	uint headInv = fld.inv (d[degd]);
+	rem = *this;
+	res.clear();
+	int t;
+	while ( (t = rem.degree() ) >= degd) {
+		int rp = t - degd;
+		if (res.size() < rp + 1) res.resize (rp + 1, 0);
+		res[rp] = fld.mult (headInv, rem[t]);
+		for (uint i = 0; i <= degd; ++i)
+			rem[i+rp] = fld.add (rem[i+rp], fld.mult (res[rp], d[i]) );
+	}
+}
+
+void polynomial::inv (polynomial&m, gf2m&fld)
+{
+	polynomial a = *this;
+	this->resize (2);
+	item (0) = 0;
+	item (1) = 1;
+	div (a, m, fld);
+}
+
+void polynomial::mod_to_fracton (polynomial&a, polynomial&b, polynomial&m, gf2m&fld)
+{
+	int deg = m.degree() / 2;
+	polynomial a0, a1, b0, b1, t1, t2;
+	a0 = m;
+	a0.make_monic (fld);
+	a1 = *this;
+	a1.mod (m, fld);
+	b0.resize (1, 0);
+	b1.resize (1, 1);
+	while (a1.degree() > deg) {
+
+		a0.divmod (a1, t1, t2, fld);
+		a0.swap (a1);
+		a1.swap (t2);
+
+		t1.mult (b1);
+		t1.mod (m);
+		t1.add (b0);
+		b0.swap (b1);
+		b1.swap (t1);
+	}
+	a = a1;
+	b = b1;
+}