implement fmap randomization while the conversion, bmap -> fmap

conv/bmap_fmap.cc

@@ -23,6 +23,9 @@ HelpFlag help {
 ValueFlag<uint32_t> fnum {
   parser, "64", "number of feature kinds", {"fnum"}, 64,
 };
+ValueFlag<uint32_t> seed {
+  parser, "0", "random seed to randomize selection of combination (0=no randomize)", {"seed"}, 0,
+};
 
 }  // namespace param
 
@@ -45,6 +48,13 @@ static bool NextCombination(auto begin, auto end) noexcept {
   return true;
 }
 
+static uint32_t XorShift(uint32_t* v) noexcept {
+  *v ^= *v << 13;
+  *v ^= *v >> 17;
+  *v ^= *v << 5;
+  return *v;
+}
+
 static void Exec() {
   const auto bmap = ReadMatrix<uint32_t>(std::cin);
   const auto fnum  = args::get(param::fnum);
@@ -54,15 +64,47 @@ static void Exec() {
     const auto  bnum   = static_cast<uint32_t>(blocks.size());
     Enforce(fnum < pow(bnum, 2), "number of blocks is too less");
 
+    // calc skip vector
+    std::vector<uint32_t> skip(fnum);
+    uint32_t nCr  = 1;
+    for (uint32_t f = 0, r = 1; f < fnum; ++r) {
+      nCr *= bnum - (r-1);
+      nCr /= r;
+
+      f += nCr;
+      if (f <= fnum || 0 == args::get(param::seed)) {
+        const auto n = fnum - (f-nCr);
+        auto begin = skip.begin() + f-nCr;
+        std::fill(begin, begin+n, 1);
+      } else {
+        uint32_t seed = args::get(param::seed);
+        uint64_t sum  = 0;
+        for (uint32_t i = f-nCr; i < fnum; ++i) {
+          skip[i] = XorShift(&seed);
+          sum += skip[i];
+        }
+        sum += XorShift(&seed);
+
+        const auto coe = nCr*1. / sum;
+        for (uint32_t i = f-nCr; i < fnum; ++i) {
+          skip[i] = std::max(uint32_t {1}, static_cast<uint32_t>(skip[i] * coe));
+        }
+      }
+    }
+
     std::vector<bool> C(bnum);
     uint32_t r = 0;
     for (uint32_t f = 0; f < fnum; ++f) {
-      if (!NextCombination(C.begin(), C.end())) {
-        ++r;
-        assert(r <= bnum);
-        std::fill(C.begin(), C.begin()+r, 1);
-        std::fill(C.begin()+r, C.end(), 0);
+      for (uint32_t s = 0; s < skip[f]; ++s) {
+        if (!NextCombination(C.begin(), C.end())) {
+          ++r;
+          assert(r <= bnum);
+          std::fill(C.begin(), C.begin()+r, 1);
+          std::fill(C.begin()+r, C.end(), 0);
+          s = 0;  // s will be 1 on next iteration
+        }
       }
+
       auto itr = C.begin();
       for (uint32_t i = 0; i < r; ++i, ++itr) {
         itr = std::find(itr, C.end(), true);