7 #ifndef _MINISKETCH_SKETCH_IMPL_H_
8 #define _MINISKETCH_SKETCH_IMPL_H_
18 void PolyMod(
const std::vector<typename F::Elem>& mod, std::vector<typename F::Elem>& val,
const F& field) {
19 size_t modsize = mod.size();
21 if (val.size() < modsize)
return;
23 while (val.size() >= modsize) {
24 auto term = val.back();
27 typename F::Multiplier mul(field, term);
28 for (
size_t x = 0; x < mod.size() - 1; ++x) {
29 val[val.size() - modsize + 1 + x] ^= mul(mod[x]);
33 while (val.size() > 0 && val.back() == 0) val.pop_back();
38 void DivMod(
const std::vector<typename F::Elem>& mod, std::vector<typename F::Elem>& val, std::vector<typename F::Elem>& div,
const F& field) {
39 size_t modsize = mod.size();
41 if (val.size() < mod.size()) {
46 div.resize(val.size() - mod.size() + 1);
47 while (val.size() >= modsize) {
48 auto term = val.back();
49 div[val.size() - modsize] = term;
52 typename F::Multiplier mul(field, term);
53 for (
size_t x = 0; x < mod.size() - 1; ++x) {
54 val[val.size() - modsize + 1 + x] ^= mul(mod[x]);
62 typename F::Elem
MakeMonic(std::vector<typename F::Elem>& a,
const F& field) {
64 if (a.back() == 1)
return 0;
65 auto inv = field.Inv(a.back());
66 typename F::Multiplier mul(field, inv);
68 for (
size_t i = 0; i < a.size() - 1; ++i) {
76 void GCD(std::vector<typename F::Elem>& a, std::vector<typename F::Elem>& b,
const F& field) {
77 if (a.size() < b.size()) std::swap(a, b);
78 while (b.size() > 0) {
92 void Sqr(std::vector<typename F::Elem>& poly,
const F& field) {
93 if (poly.size() == 0)
return;
94 poly.resize(poly.size() * 2 - 1);
95 for (
int x = poly.size() - 1; x >= 0; --x) {
96 poly[x] = (x & 1) ? 0 : field.Sqr(poly[x / 2]);
102 void TraceMod(
const std::vector<typename F::Elem>& mod, std::vector<typename F::Elem>&
out,
const typename F::Elem& param,
const F& field) {
103 out.reserve(mod.size() * 2);
108 for (
int i = 0; i < field.Bits() - 1; ++i) {
110 if (
out.size() < 2)
out.resize(2);
128 bool RecFindRoots(std::vector<std::vector<typename F::Elem>>& stack,
size_t pos, std::vector<typename F::Elem>& roots,
bool fully_factorizable,
int depth,
typename F::Elem randv,
const F& field) {
129 auto& ppoly = stack[pos];
134 CHECK_SAFE(ppoly.size() > 1 && ppoly.back() == 1);
136 if (ppoly.size() == 2) {
137 roots.push_back(ppoly[0]);
141 if (ppoly.size() == 3) {
143 auto input = field.Mul(ppoly[0], field.Sqr(field.Inv(ppoly[1])));
144 auto root = field.Qrt(input);
145 if ((field.Sqr(root) ^ root) != input) {
149 auto sol = field.Mul(root, ppoly[1]);
150 roots.push_back(sol);
151 roots.push_back(sol ^ ppoly[1]);
155 if (pos + 3 > stack.size()) {
157 stack.resize((pos + 3) * 2);
159 auto& poly = stack[pos];
160 auto& tmp = stack[pos + 1];
161 auto& trace = stack[pos + 2];
164 for (
int iter = 0;; ++iter) {
167 TraceMod(poly, trace, randv, field);
169 if (iter >= 1 && !fully_factorizable) {
205 for (
size_t i = 0; i < trace.size(); ++i) {
208 while (tmp.size() && tmp.back() == 0) tmp.pop_back();
215 if (tmp.size() != 0)
return false;
216 fully_factorizable =
true;
219 if (fully_factorizable) {
224 CHECK_RETURN(field.Bits() - depth >= std::numeric_limits<decltype(poly.size())>::digits ||
225 (poly.size() - 2) >> (field.Bits() - depth) == 0,
false);
231 randv = field.Mul2(randv);
233 GCD(trace, tmp, field);
234 if (trace.size() != poly.size() && trace.size() > 1)
break;
237 DivMod(trace, poly, tmp, field);
241 std::swap(poly, trace);
244 if (!
RecFindRoots(stack, pos + 1, roots, fully_factorizable, depth, randv, field))
return false;
263 std::vector<typename F::Elem>
FindRoots(
const std::vector<typename F::Elem>& poly,
typename F::Elem basis,
const F& field) {
264 std::vector<typename F::Elem> roots;
267 if (poly.size() == 1)
return roots;
268 roots.reserve(poly.size() - 1);
269 std::vector<std::vector<typename F::Elem>> stack = {poly};
272 if (!
RecFindRoots(stack, 0, roots,
false, 0, basis, field)) {
281 std::vector<typename F::Elem>
BerlekampMassey(
const std::vector<typename F::Elem>& syndromes,
size_t max_degree,
const F& field) {
282 std::vector<typename F::Multiplier> table;
283 std::vector<typename F::Elem> current, prev, tmp;
284 current.reserve(syndromes.size() / 2 + 1);
285 prev.reserve(syndromes.size() / 2 + 1);
286 tmp.reserve(syndromes.size() / 2 + 1);
291 typename F::Elem b = 1, b_inv = 1;
292 bool b_have_inv =
true;
293 table.reserve(syndromes.size());
295 for (
size_t n = 0; n != syndromes.size(); ++n) {
296 table.emplace_back(field, syndromes[n]);
297 auto discrepancy = syndromes[n];
298 for (
size_t i = 1; i < current.size(); ++i) discrepancy ^= table[n - i](current[i]);
299 if (discrepancy != 0) {
300 int x = n + 1 - (current.size() - 1) - (prev.size() - 1);
302 b_inv = field.Inv(b);
305 bool swap = 2 * (current.size() - 1) <= n;
307 if (prev.size() + x - 1 > max_degree)
return {};
309 current.resize(prev.size() + x);
311 typename F::Multiplier mul(field, field.Mul(discrepancy, b_inv));
312 for (
size_t i = 0; i < prev.size(); ++i) current[i + x] ^= mul(prev[i]);
314 std::swap(prev, tmp);
320 CHECK_RETURN(current.size() && current.back() != 0, {});
326 std::vector<typename F::Elem> all_syndromes;
327 all_syndromes.resize(odd_syndromes.size() * 2);
328 for (
size_t i = 0; i < odd_syndromes.size(); ++i) {
329 all_syndromes[i * 2] = odd_syndromes[i];
330 all_syndromes[i * 2 + 1] = field.Sqr(all_syndromes[i]);
332 return all_syndromes;
336 void AddToOddSyndromes(std::vector<typename F::Elem>& osyndromes,
typename F::Elem data,
const F& field) {
337 auto sqr = field.Sqr(data);
338 typename F::Multiplier mul(field, sqr);
339 for (
auto& osyndrome : osyndromes) {
346 std::vector<typename F::Elem>
FullDecode(
const std::vector<typename F::Elem>& osyndromes,
const F& field) {
347 auto asyndromes = ReconstructAllSyndromes<typename F::Elem>(osyndromes, field);
349 std::reverse(poly.begin(), poly.end());
361 template<
typename... Args>
363 std::random_device rng;
364 std::uniform_int_distribution<uint64_t> dist;
371 void Add(uint64_t val)
override
373 auto elem =
m_field.FromUint64(val);
381 m_field.Serialize(writer, val);
390 val =
m_field.Deserialize(reader);
398 if (poly.size() == 0)
return -1;
399 if (poly.size() == 1)
return 0;
400 if ((
int)poly.size() > 1 + max_count)
return -1;
401 std::reverse(poly.begin(), poly.end());
403 if (roots.size() == 0)
return -1;
405 for (
const auto& root : roots) {
425 if (seed == (uint64_t)-1) {
Abstract class for internal representation of a minisketch object.
SketchImpl(int implementation, int bits, const Args &... args)
void Serialize(unsigned char *ptr) const override
size_t Merge(const Sketch *other_sketch) override
int Decode(int max_count, uint64_t *out) const override
void Add(uint64_t val) override
void Deserialize(const unsigned char *ptr) override
size_t Syndromes() const override
void Init(int count) override
std::vector< typename F::Elem > m_syndromes
void SetSeed(uint64_t seed) override
#define CHECK_RETURN(cond, rvar)
Check a condition and return on failure in non-verify builds, crash in verify builds.
#define CHECK_SAFE(cond)
Check macro that does nothing in normal non-verify builds but crashes in verify builds.
void AddToOddSyndromes(std::vector< typename F::Elem > &osyndromes, typename F::Elem data, const F &field)
void TraceMod(const std::vector< typename F::Elem > &mod, std::vector< typename F::Elem > &out, const typename F::Elem ¶m, const F &field)
Compute the trace map of (param*x) modulo mod, putting the result in out.
bool RecFindRoots(std::vector< std::vector< typename F::Elem >> &stack, size_t pos, std::vector< typename F::Elem > &roots, bool fully_factorizable, int depth, typename F::Elem randv, const F &field)
One step of the root finding algorithm; finds roots of stack[pos] and adds them to roots.
void PolyMod(const std::vector< typename F::Elem > &mod, std::vector< typename F::Elem > &val, const F &field)
Compute the remainder of a polynomial division of val by mod, putting the result in mod.
void GCD(std::vector< typename F::Elem > &a, std::vector< typename F::Elem > &b, const F &field)
Compute the GCD of two polynomials, putting the result in a.
F::Elem MakeMonic(std::vector< typename F::Elem > &a, const F &field)
Make a polynomial monic.
void Sqr(std::vector< typename F::Elem > &poly, const F &field)
Square a polynomial.
std::vector< typename F::Elem > ReconstructAllSyndromes(const std::vector< typename F::Elem > &odd_syndromes, const F &field)
std::vector< typename F::Elem > BerlekampMassey(const std::vector< typename F::Elem > &syndromes, size_t max_degree, const F &field)
void DivMod(const std::vector< typename F::Elem > &mod, std::vector< typename F::Elem > &val, std::vector< typename F::Elem > &div, const F &field)
Compute the quotient of a polynomial division of val by mod, putting the quotient in div and the rema...
std::vector< typename F::Elem > FullDecode(const std::vector< typename F::Elem > &osyndromes, const F &field)
std::vector< typename F::Elem > FindRoots(const std::vector< typename F::Elem > &poly, typename F::Elem basis, const F &field)
Returns the roots of a fully factorizable polynomial.