diff --git a/autogen.sh b/autogen.sh index 7d4207c..eb02c46 100755 --- a/autogen.sh +++ b/autogen.sh @@ -24,7 +24,7 @@ echo "noinst_HEADERS += `find lib/ -type f -name \*.h |tr \"\n\" \" \" `" >>$OUT echo "libcodecrypt_la_CPPFLAGS = -I\$(srcdir)/lib/ ${COMMON_CPPFLAGS}" >>$OUT echo "libcodecrypt_la_CFLAGS = ${COMMON_CFLAGS}" >>$OUT echo "libcodecrypt_la_LDFLAGS = ${COMMON_LDFLAGS}" >>$OUT -#echo "libcodecrypt_la_LDADD = ${COMMON_LDADD} " >>$OUT +echo "libcodecrypt_la_LIBADD = -lgmp ${COMMON_LDADD} " >>$OUT [ -f "lib/Makefile.am.extra" ] && cat "lib/Makefile.am.extra" >>$OUT echo "bin_PROGRAMS = ccr" >>$OUT diff --git a/include/codecrypt.h b/include/codecrypt.h index 2bccedd..d79626e 100644 --- a/include/codecrypt.h +++ b/include/codecrypt.h @@ -23,7 +23,7 @@ typedef unsigned int uint; /* * vector over GF(2). We rely on STL's vector == bit_vector - * specialization for efficiency. + * specialization for space efficiency. */ class polynomial; class gf2m; @@ -36,8 +36,12 @@ public: void add (const bvector&); bool operator* (const bvector&); //dot product bool zero() const; + void to_poly (polynomial&, gf2m&); void from_poly (const polynomial&, gf2m&); + + void colex_rank (bvector&); + void colex_unrank (bvector&, uint n, uint k); }; /* diff --git a/lib/bvector.cpp b/lib/bvector.cpp index 9630af7..7d66bc0 100644 --- a/lib/bvector.cpp +++ b/lib/bvector.cpp @@ -47,3 +47,115 @@ void bvector::from_poly (const polynomial&r, gf2m&fld) for (uint i = 0; i < size(); ++i) item (i) = (r[i/fld.m] >> (i % fld.m) ) & 1; } + +/* + * utility colex (un)ranking for niederreiter and workalikes. + * see Ruskey's Combinatorial Generation, algorithm 4.10 + * + * TODO use (external) cache for combination numbers to speed this up. + */ + +#include + +static void combination_number (uint n, uint k, mpz_t& r) +{ + mpz_t t; + if (k > n) { + mpz_set_ui (r, 0); + return; + } + + if (k * 2 > n) k = n - k; + + mpz_set_ui (r, 1); + mpz_init (t); + + //upper part n*(n-1)*(n-2)*...*(n-k+1) + for (uint i = n; i > n - k; --i) { + mpz_swap (t, r); + mpz_mul_ui (r, t, i); + } + //lower part (div k!) + for (uint i = k; i > 1; --i) { + mpz_swap (t, r); + mpz_tdiv_q_ui (r, t, i); + } + + mpz_clear (t); +} + +static void bvector_to_mpz (bvector&v, mpz_t&r) +{ + mpz_set_ui (r, 0); + mpz_realloc2 (r, v.size() ); + for (uint i = 0; i < v.size(); ++i) + if (v[i]) + mpz_setbit (r, i); + else mpz_clrbit (r, i); +} + +static void mpz_to_bvector (mpz_t&x, bvector&r) +{ + r.resize (mpz_sizeinbase (x, 2) ); + for (uint i = 0; i < r.size(); ++i) + r[i] = mpz_tstbit (x, i); +} + +void bvector::colex_rank (bvector&r) +{ + mpz_t res, t, t2; + mpz_init_set_ui (res, 0); + mpz_init (t); + mpz_init (t2); + + uint i, j; + j = 1; + for (i = 0; i < size(); ++i) + if (item (i) ) { + combination_number (i, j, t); + mpz_swap (t2, res); + mpz_add (res, t, t2); + ++j; + } + + mpz_to_bvector (res, r); + + mpz_clear (t); + mpz_clear (t2); + mpz_clear (res); +} + +#include +void bvector::colex_unrank (bvector&res, uint n, uint k) +{ + mpz_t r, t, t2; + mpz_init (r); + mpz_init (t); + mpz_init (t2); + + bvector_to_mpz (*this, r); + + res.clear(); + res.resize (n, 0); + + uint p; + for (uint i = k; i > 0; --i) { + p = i; + for (;;) { + combination_number (p, i, t); + + if (mpz_cmp (t, r) > 0) break; + ++p; + } + + combination_number (p - 1, i, t); + mpz_swap (t2, r); + mpz_sub (r, t2, t); + if (p > n) continue; //overflow protection + res[p-1] = 1; + } + + mpz_clear (r); + mpz_clear (t); + mpz_clear (t2); +}