mirror of
https://github.com/aicodix/code.git
synced 2026-04-27 14:30:36 +00:00
added Galois field arithmetic
This commit is contained in:
Ahmet Inan
parent
266f299f86
commit
2a7fa2b7b7
3 changed files with 352 additions and 0 deletions
|
|
@ -31,3 +31,7 @@ Here a [Pseudorandom number generator](https://en.wikipedia.org/wiki/Pseudorando
|
|||
|
||||
When dealing with unaligned and arbitrary-bit-sized elements in a data stream, the bitwise stream container might help avoiding some headaches.
|
||||
|
||||
### [galois_field.hh](galois_field.hh)
|
||||
|
||||
We have to thank Évariste Galois' contribution of the [Finite field](https://en.wikipedia.org/wiki/Finite_field) to mathematics, which laid the cornerstone for a variety of applications that we take for granted today.
|
||||
|
||||
|
|
|
|||
302
galois_field.hh
Normal file
302
galois_field.hh
Normal file
|
|
@ -0,0 +1,302 @@
|
|||
/*
|
||||
Galois field arithmetic
|
||||
|
||||
Copyright 2018 Ahmet Inan <inan@aicodix.de>
|
||||
*/
|
||||
|
||||
#ifndef GALOIS_FIELD_HH
|
||||
#define GALOIS_FIELD_HH
|
||||
|
||||
#include <cassert>
|
||||
|
||||
namespace CODE {
|
||||
namespace GF {
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
struct Index;
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
struct Value
|
||||
{
|
||||
static const int Q = 1 << M, N = Q - 1;
|
||||
static_assert(M <= 8 * sizeof(TYPE), "TYPE not wide enough");
|
||||
static_assert(Q == (POLY & ~N), "POLY not of degree Q");
|
||||
TYPE v;
|
||||
Value() {}
|
||||
explicit Value(TYPE v) : v(v)
|
||||
{
|
||||
assert(v <= N);
|
||||
}
|
||||
explicit operator bool () const { return v; }
|
||||
explicit operator int () const { return v; }
|
||||
Value<M, POLY, TYPE> operator *= (Index<M, POLY, TYPE> a)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
return *this = *this * a;
|
||||
}
|
||||
Value<M, POLY, TYPE> operator *= (Value<M, POLY, TYPE> a)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
return *this = *this * a;
|
||||
}
|
||||
Value<M, POLY, TYPE> operator += (Value<M, POLY, TYPE> a)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
return *this = *this + a;
|
||||
}
|
||||
static const Value<M, POLY, TYPE> zero()
|
||||
{
|
||||
return Value<M, POLY, TYPE>(0);
|
||||
}
|
||||
};
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
struct Index
|
||||
{
|
||||
static const int Q = 1 << M, N = Q - 1;
|
||||
static_assert(M <= 8 * sizeof(TYPE), "TYPE not wide enough");
|
||||
static_assert(Q == (POLY & ~N), "POLY not of degree Q");
|
||||
TYPE i;
|
||||
Index() {}
|
||||
explicit Index(TYPE i) : i(i)
|
||||
{
|
||||
assert(i < modulus());
|
||||
}
|
||||
explicit operator int () const { return i; }
|
||||
Index<M, POLY, TYPE> operator *= (Index<M, POLY, TYPE> a)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(i < modulus());
|
||||
return *this = *this * a;
|
||||
}
|
||||
Index<M, POLY, TYPE> operator /= (Index<M, POLY, TYPE> a)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(i < modulus());
|
||||
return *this = *this / a;
|
||||
}
|
||||
static const TYPE modulus()
|
||||
{
|
||||
return N;
|
||||
}
|
||||
};
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
struct Tables
|
||||
{
|
||||
static const int Q = 1 << M, N = Q - 1;
|
||||
static_assert(M <= 8 * sizeof(TYPE), "TYPE not wide enough");
|
||||
static_assert(Q == (POLY & ~N), "POLY not of degree Q");
|
||||
static TYPE *LOG, *EXP;
|
||||
TYPE log_[Q], exp_[Q];
|
||||
static TYPE next(TYPE a)
|
||||
{
|
||||
return a & (TYPE)(Q >> 1) ? (a << 1) ^ (TYPE)POLY : a << 1;
|
||||
}
|
||||
static TYPE log(TYPE a)
|
||||
{
|
||||
assert(LOG != nullptr);
|
||||
assert(a <= N);
|
||||
return LOG[a];
|
||||
}
|
||||
static TYPE exp(TYPE a)
|
||||
{
|
||||
assert(EXP != nullptr);
|
||||
assert(a <= N);
|
||||
return EXP[a];
|
||||
}
|
||||
Tables()
|
||||
{
|
||||
assert(LOG == nullptr);
|
||||
LOG = log_;
|
||||
assert(EXP == nullptr);
|
||||
EXP = exp_;
|
||||
log_[exp_[N] = 0] = N;
|
||||
TYPE a = 1;
|
||||
for (int i = 0; i < N; ++i, a = next(a))
|
||||
log_[exp_[i] = a] = i;
|
||||
assert(1 == a);
|
||||
}
|
||||
~Tables()
|
||||
{
|
||||
assert(LOG != nullptr);
|
||||
LOG = nullptr;
|
||||
assert(EXP != nullptr);
|
||||
EXP = nullptr;
|
||||
}
|
||||
};
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
TYPE *Tables<M, POLY, TYPE>::LOG;
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
TYPE *Tables<M, POLY, TYPE>::EXP;
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Index<M, POLY, TYPE> index(Value<M, POLY, TYPE> a)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(a.v);
|
||||
return Index<M, POLY, TYPE>(Tables<M, POLY, TYPE>::log(a.v));
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> value(Index<M, POLY, TYPE> a) {
|
||||
assert(a.i < a.modulus());
|
||||
return Value<M, POLY, TYPE>(Tables<M, POLY, TYPE>::exp(a.i));
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
bool operator == (Value<M, POLY, TYPE> a, Value<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.v <= b.N);
|
||||
return a.v == b.v;
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
bool operator != (Value<M, POLY, TYPE> a, Value<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.v <= b.N);
|
||||
return a.v != b.v;
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> operator + (Value<M, POLY, TYPE> a, Value<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.v <= b.N);
|
||||
return Value<M, POLY, TYPE>(a.v ^ b.v);
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Index<M, POLY, TYPE> operator * (Index<M, POLY, TYPE> a, Index<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(b.i < b.modulus());
|
||||
TYPE tmp = a.i + b.i;
|
||||
return Index<M, POLY, TYPE>(a.modulus() - a.i <= b.i ? tmp - a.modulus() : tmp);
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> operator * (Value<M, POLY, TYPE> a, Value<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.v <= b.N);
|
||||
return (!a.v || !b.v) ? a.zero() : value(index(a) * index(b));
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> rcp(Value<M, POLY, TYPE> a)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(a.v);
|
||||
return value(Index<M, POLY, TYPE>(index(a).modulus() - index(a).i));
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Index<M, POLY, TYPE> operator / (Index<M, POLY, TYPE> a, Index<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(b.i < b.modulus());
|
||||
TYPE tmp = a.i - b.i;
|
||||
return Index<M, POLY, TYPE>(a.i < b.i ? tmp + a.modulus() : tmp);
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> operator / (Value<M, POLY, TYPE> a, Value<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.v <= b.N);
|
||||
assert(b.v);
|
||||
return !a.v ? a.zero() : value(index(a) / index(b));
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> operator / (Index<M, POLY, TYPE> a, Value<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(b.v <= b.N);
|
||||
assert(b.v);
|
||||
return value(a / index(b));
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> operator / (Value<M, POLY, TYPE> a, Index<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.i < b.modulus());
|
||||
return !a.v ? a.zero() : value(index(a) / b);
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> operator * (Index<M, POLY, TYPE> a, Value<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(b.v <= b.N);
|
||||
return !b.v ? b.zero() : value(a * index(b));
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> operator * (Value<M, POLY, TYPE> a, Index<M, POLY, TYPE> b)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.i < b.modulus());
|
||||
return !a.v ? a.zero() : value(index(a) * b);
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> fma(Index<M, POLY, TYPE> a, Index<M, POLY, TYPE> b, Value<M, POLY, TYPE> c)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(b.i < b.modulus());
|
||||
assert(c.v <= c.N);
|
||||
return value(a * b) + c;
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> fma(Index<M, POLY, TYPE> a, Value<M, POLY, TYPE> b, Value<M, POLY, TYPE> c)
|
||||
{
|
||||
assert(a.i < a.modulus());
|
||||
assert(b.v <= b.N);
|
||||
assert(c.v <= c.N);
|
||||
return !b.v ? c : (value(a * index(b)) + c);
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> fma(Value<M, POLY, TYPE> a, Index<M, POLY, TYPE> b, Value<M, POLY, TYPE> c)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.i < b.modulus());
|
||||
assert(c.v <= c.N);
|
||||
return !a.v ? c : (value(index(a) * b) + c);
|
||||
}
|
||||
|
||||
template <int M, int POLY, typename TYPE>
|
||||
Value<M, POLY, TYPE> fma(Value<M, POLY, TYPE> a, Value<M, POLY, TYPE> b, Value<M, POLY, TYPE> c)
|
||||
{
|
||||
assert(a.v <= a.N);
|
||||
assert(b.v <= b.N);
|
||||
assert(c.v <= c.N);
|
||||
return (!a.v || !b.v) ? c : (value(index(a) * index(b)) + c);
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
template <int WIDTH, int POLY, typename TYPE>
|
||||
class GaloisField
|
||||
{
|
||||
public:
|
||||
static const int M = WIDTH, Q = 1 << M, N = Q - 1;
|
||||
static_assert(M <= 8 * sizeof(TYPE), "TYPE not wide enough");
|
||||
static_assert(Q == (POLY & ~N), "POLY not of degree Q");
|
||||
typedef TYPE value_type;
|
||||
typedef GF::Value<M, POLY, TYPE> ValueType;
|
||||
typedef GF::Index<M, POLY, TYPE> IndexType;
|
||||
private:
|
||||
GF::Tables<M, POLY, TYPE> Tables;
|
||||
};
|
||||
|
||||
}
|
||||
#endif
|
||||
46
tests/gf_test.cc
Normal file
46
tests/gf_test.cc
Normal file
|
|
@ -0,0 +1,46 @@
|
|||
/*
|
||||
Test for the Galois field arithmetic
|
||||
|
||||
Copyright 2018 Ahmet Inan <inan@aicodix.de>
|
||||
*/
|
||||
|
||||
#include <cassert>
|
||||
#include <iostream>
|
||||
#include "galois_field.hh"
|
||||
|
||||
int main()
|
||||
{
|
||||
if (1) {
|
||||
// BBC WHP031 RS(15, 11) T=2
|
||||
typedef CODE::GaloisField<4, 0b10011, uint8_t> GF;
|
||||
GF instance;
|
||||
GF::ValueType a(3), b(7), c(15), d(6);
|
||||
assert(a * b + c == d);
|
||||
}
|
||||
if (1) {
|
||||
// DVB-T RS(255, 239) T=8
|
||||
typedef CODE::GaloisField<8, 0b100011101, uint8_t> GF;
|
||||
GF instance;
|
||||
GF::ValueType a(154), b(83), c(144), d(63);
|
||||
assert(fma(a, b, c) == d);
|
||||
}
|
||||
if (1) {
|
||||
// FUN RS(65535, 65471) T=32
|
||||
typedef CODE::GaloisField<16, 0b10001000000001011, uint16_t> GF;
|
||||
GF *instance = new GF();
|
||||
GF::ValueType a(42145), b(13346), c(40958), d(35941);
|
||||
assert(a / b + c == d);
|
||||
delete instance;
|
||||
}
|
||||
if (1) {
|
||||
// DVB-S2 FULL BCH(65535, 65343) T=12
|
||||
typedef CODE::GaloisField<16, 0b10000000000101101, uint16_t> GF;
|
||||
GF *instance = new GF();
|
||||
GF::ValueType a(38532), b(7932), c(34283), d(22281);
|
||||
assert(a / (rcp(b) + c) == d);
|
||||
delete instance;
|
||||
}
|
||||
std::cerr << "Galois field arithmetic test passed!" << std::endl;
|
||||
return 0;
|
||||
}
|
||||
|
||||
Loading…
Add table
Add a link
Reference in a new issue