#pragma once
#include <cassert>
#include <concepts>
#include <cstdint>
#include <bit>
#include <functional>
#include <stdexcept>
#include <vector>
namespace librscpp {
template <std::unsigned_integral T, std::integral HP>
class GFBase
{
static_assert(sizeof(T) <= sizeof(HP) && sizeof(HP) >= 2,
"HP must be able to hold values of type T and the result of their addition");
protected:
const int _size;
const int _fcr;
std::vector<T> _expTable, _logTable;
inline static T fast_mod(HP input, HP mod) { return input < mod ? input : input - mod; }
GFBase(int size, int fcr, std::function<T(T)> next) : _size(size), _fcr(fcr)
{
if (size <= 0 || size > (1 << (sizeof(T) * 8)))
throw std::invalid_argument("Invalid field size " + std::to_string(size));
if (fcr < 0 || fcr >= size)
throw std::invalid_argument("Invalid first consecutive root " + std::to_string(fcr) + " for field of size "
+ std::to_string(size));
#ifdef LIBRSCPP_SAVE_MEMORY
_expTable.resize(size, 0);
#else
_expTable.resize(size * 2, 0);
#endif
_logTable.resize(size, 0);
int x = 1;
for (int i = 0; i < size; ++i) {
_expTable[i] = x;
x = next(x);
}
#ifndef LIBRSCPP_SAVE_MEMORY
for (int i = size - 1; i < size * 2; ++i)
_expTable[i] = _expTable[i - (size - 1)];
#endif
for (int i = 0; i < size - 1; ++i)
_logTable[_expTable[i]] = i;
}
public:
using value_type = T;
T exp(T a) const { return _expTable.at((a)); }
T log(T a) const
{
assert(a != 0);
return _logTable.at((a));
}
T inv(T a) const { return _expTable[_size - 1 - log(a)]; }
T mul(T a, T b) const noexcept
{
if (a == 0 || b == 0)
return 0;
#ifdef LIBRSCPP_SAVE_MEMORY
return _expTable[fast_mod(HP(_logTable[a]) + HP(_logTable[b]), _size - 1)];
#else
return _expTable[HP(_logTable[a]) + HP(_logTable[b])];
#endif
}
int size() const noexcept { return _size; }
int fcr() const noexcept { return _fcr; }
};
template <std::unsigned_integral T = uint16_t, std::integral HP = uint32_t>
class GF2n : public GFBase<T, HP>
{
public:
GF2n(uint32_t poly, int fcr)
: GFBase<T, HP>(1 << (std::bit_width(poly) - 1), fcr, [poly, size = HP(1 << (std::bit_width(poly) - 1))](T x) {
HP ix = x;
ix *= 2; if (ix >= size) ix ^= poly;
return ix;
})
{}
T add(T a, T b) const { return a ^ b; }
T sub(T a, T b) const { return a ^ b; } T neg(T a) const { return a; } };
template <std::unsigned_integral T = uint16_t, std::integral HP = uint32_t>
class GFp : public GFBase<T, HP>
{
using Base = GFBase<T, HP>;
public:
GFp(uint32_t prime, uint32_t generator, int fcr)
: GFBase<T, HP>(prime, fcr, [prime, generator](T x) { return (x * generator) % prime; })
{}
T add(T a, T b) const { return Base::fast_mod(HP(a) + HP(b), Base::size()); }
T sub(T a, T b) const { return Base::fast_mod(Base::size() + HP(a) - HP(b), Base::size()); }
T neg(T a) const { return sub(0, a); }
};
}