books and problem

src/main-bob-test-01.cpp

@@ -0,0 +1,400 @@
+// Minesweeper next-move solver from problems/book.pdf: compute the probability
+// that the safest optimal hidden-cell move is not a bomb.
+
+#include <algorithm>
+#include <cstdio>
+#include <functional>
+#include <iostream>
+#include <numeric>
+#include <string>
+#include <unordered_map>
+#include <utility>
+#include <vector>
+
+namespace {
+
+struct Constraint {
+  std::vector<int> hidden;
+  int need = 0;
+};
+
+struct LocalConstraint {
+  std::vector<int> cells;
+  int need = 0;
+};
+
+struct ComponentResult {
+  std::vector<int> globalHiddenIds;
+  std::vector<long double> ways;
+  std::vector<std::vector<long double>> bombCount;
+};
+
+struct DisjointSet {
+  std::vector<int> parent;
+  std::vector<int> rank;
+
+  explicit DisjointSet(int n) : parent(n), rank(n, 0) {
+    std::iota(parent.begin(), parent.end(), 0);
+  }
+
+  int find(int x) {
+    if (parent[x] != x) {
+      parent[x] = find(parent[x]);
+    }
+    return parent[x];
+  }
+
+  void unite(int a, int b) {
+    a = find(a);
+    b = find(b);
+    if (a == b) {
+      return;
+    }
+    if (rank[a] < rank[b]) {
+      std::swap(a, b);
+    }
+    parent[b] = a;
+    if (rank[a] == rank[b]) {
+      ++rank[a];
+    }
+  }
+};
+
+long double choose(int n, int k) {
+  if (k < 0 || k > n) {
+    return 0.0L;
+  }
+  if (k == 0 || k == n) {
+    return 1.0L;
+  }
+  k = std::min(k, n - k);
+  long double result = 1.0L;
+  for (int i = 1; i <= k; ++i) {
+    result *= static_cast<long double>(n - k + i);
+    result /= static_cast<long double>(i);
+  }
+  return result;
+}
+
+std::vector<long double> convolve(const std::vector<long double>& a,
+                                  const std::vector<long double>& b) {
+  std::vector<long double> out(a.size() + b.size() - 1, 0.0L);
+  for (int i = 0; i < static_cast<int>(a.size()); ++i) {
+    if (a[i] == 0.0L) {
+      continue;
+    }
+    for (int j = 0; j < static_cast<int>(b.size()); ++j) {
+      if (b[j] == 0.0L) {
+        continue;
+      }
+      out[i + j] += a[i] * b[j];
+    }
+  }
+  return out;
+}
+
+ComponentResult enumerateComponent(const std::vector<int>& globalHiddenIds,
+                                   const std::vector<LocalConstraint>& constraints) {
+  const int n = static_cast<int>(globalHiddenIds.size());
+  ComponentResult result;
+  result.globalHiddenIds = globalHiddenIds;
+  result.ways.assign(n + 1, 0.0L);
+  result.bombCount.assign(n, std::vector<long double>(n + 1, 0.0L));
+
+  std::vector<std::vector<int>> incident(n);
+  for (int ci = 0; ci < static_cast<int>(constraints.size()); ++ci) {
+    for (int cell : constraints[ci].cells) {
+      incident[cell].push_back(ci);
+    }
+  }
+
+  std::vector<int> order(n);
+  std::iota(order.begin(), order.end(), 0);
+  std::sort(order.begin(), order.end(), [&](int lhs, int rhs) {
+    return incident[lhs].size() > incident[rhs].size();
+  });
+
+  std::vector<int> assigned(constraints.size(), 0);
+  std::vector<int> unassigned(constraints.size(), 0);
+  for (int ci = 0; ci < static_cast<int>(constraints.size()); ++ci) {
+    unassigned[ci] = static_cast<int>(constraints[ci].cells.size());
+  }
+
+  std::vector<int> current(n, 0);
+
+  std::function<void(int, int)> dfs = [&](int pos, int bombsUsed) {
+    if (pos == n) {
+      for (int ci = 0; ci < static_cast<int>(constraints.size()); ++ci) {
+        if (assigned[ci] != constraints[ci].need) {
+          return;
+        }
+      }
+      result.ways[bombsUsed] += 1.0L;
+      for (int cell = 0; cell < n; ++cell) {
+        if (current[cell] == 1) {
+          result.bombCount[cell][bombsUsed] += 1.0L;
+        }
+      }
+      return;
+    }
+
+    const int cell = order[pos];
+    for (int value = 0; value <= 1; ++value) {
+      bool ok = true;
+      current[cell] = value;
+      for (int ci : incident[cell]) {
+        --unassigned[ci];
+        assigned[ci] += value;
+        if (assigned[ci] > constraints[ci].need ||
+            assigned[ci] + unassigned[ci] < constraints[ci].need) {
+          ok = false;
+        }
+      }
+      if (ok) {
+        dfs(pos + 1, bombsUsed + value);
+      }
+      for (int ci : incident[cell]) {
+        assigned[ci] -= value;
+        ++unassigned[ci];
+      }
+      current[cell] = 0;
+    }
+  };
+
+  dfs(0, 0);
+  return result;
+}
+
+}  // namespace
+
+int main() {
+  std::ios::sync_with_stdio(false);
+  std::cin.tie(nullptr);
+
+  int r = 0;
+  int c = 0;
+  int b = 0;
+  if (!(std::cin >> r >> c >> b)) {
+    return 0;
+  }
+
+  std::vector<std::string> grid(r);
+  for (int i = 0; i < r; ++i) {
+    if (!(std::cin >> grid[i]) || static_cast<int>(grid[i].size()) != c) {
+      std::printf("0.000\n");
+      return 0;
+    }
+  }
+
+  std::vector<std::vector<int>> hiddenId(r, std::vector<int>(c, -1));
+  int hiddenCount = 0;
+  bool invalidInput = false;
+  for (int i = 0; i < r; ++i) {
+    for (int j = 0; j < c; ++j) {
+      const char ch = grid[i][j];
+      if (ch == '#') {
+        hiddenId[i][j] = hiddenCount++;
+      } else if (ch < '0' || ch > '8') {
+        invalidInput = true;
+      }
+    }
+  }
+
+  if (invalidInput || hiddenCount == 0) {
+    std::printf("0.000\n");
+    return 0;
+  }
+
+  constexpr int dr[8] = {-1, -1, -1, 0, 0, 1, 1, 1};
+  constexpr int dc[8] = {-1, 0, 1, -1, 1, -1, 0, 1};
+
+  std::vector<Constraint> constraints;
+  std::vector<char> isFrontier(hiddenCount, 0);
+  bool impossible = false;
+
+  for (int i = 0; i < r; ++i) {
+    for (int j = 0; j < c; ++j) {
+      if (grid[i][j] == '#') {
+        continue;
+      }
+      Constraint constraint;
+      constraint.need = grid[i][j] - '0';
+      for (int dir = 0; dir < 8; ++dir) {
+        const int ni = i + dr[dir];
+        const int nj = j + dc[dir];
+        if (ni < 0 || ni >= r || nj < 0 || nj >= c) {
+          continue;
+        }
+        if (grid[ni][nj] == '#') {
+          const int id = hiddenId[ni][nj];
+          constraint.hidden.push_back(id);
+          isFrontier[id] = 1;
+        }
+      }
+      if (constraint.hidden.empty()) {
+        if (constraint.need != 0) {
+          impossible = true;
+        }
+      } else {
+        if (constraint.need < 0 || constraint.need > static_cast<int>(constraint.hidden.size())) {
+          impossible = true;
+        }
+        constraints.push_back(std::move(constraint));
+      }
+    }
+  }
+
+  if (impossible || b < 0 || b > hiddenCount) {
+    std::printf("0.000\n");
+    return 0;
+  }
+
+  std::vector<int> frontierHiddenIds;
+  frontierHiddenIds.reserve(hiddenCount);
+  for (int id = 0; id < hiddenCount; ++id) {
+    if (isFrontier[id]) {
+      frontierHiddenIds.push_back(id);
+    }
+  }
+
+  const int frontierCount = static_cast<int>(frontierHiddenIds.size());
+  const int interiorCount = hiddenCount - frontierCount;
+
+  if (frontierCount == 0) {
+    const long double answer = 1.0L - static_cast<long double>(b) / static_cast<long double>(hiddenCount);
+    std::printf("%.3Lf\n", std::clamp(answer, 0.0L, 1.0L));
+    return 0;
+  }
+
+  std::vector<int> frontierPos(hiddenCount, -1);
+  for (int i = 0; i < frontierCount; ++i) {
+    frontierPos[frontierHiddenIds[i]] = i;
+  }
+
+  DisjointSet dsu(frontierCount);
+  for (const auto& constraint : constraints) {
+    if (constraint.hidden.empty()) {
+      continue;
+    }
+    const int first = frontierPos[constraint.hidden.front()];
+    for (int idx = 1; idx < static_cast<int>(constraint.hidden.size()); ++idx) {
+      dsu.unite(first, frontierPos[constraint.hidden[idx]]);
+    }
+  }
+
+  std::unordered_map<int, int> componentIndex;
+  std::vector<std::vector<int>> componentCells;
+  for (int fp = 0; fp < frontierCount; ++fp) {
+    const int root = dsu.find(fp);
+    auto [it, inserted] = componentIndex.emplace(root, static_cast<int>(componentCells.size()));
+    if (inserted) {
+      componentCells.push_back({});
+    }
+    componentCells[it->second].push_back(frontierHiddenIds[fp]);
+  }
+
+  std::vector<std::vector<LocalConstraint>> componentConstraints(componentCells.size());
+  for (const auto& constraint : constraints) {
+    const int root = dsu.find(frontierPos[constraint.hidden.front()]);
+    const int cid = componentIndex[root];
+    std::unordered_map<int, int> localIndex;
+    const auto& cells = componentCells[cid];
+    for (int i = 0; i < static_cast<int>(cells.size()); ++i) {
+      localIndex.emplace(cells[i], i);
+    }
+    LocalConstraint local;
+    local.need = constraint.need;
+    local.cells.reserve(constraint.hidden.size());
+    for (int globalId : constraint.hidden) {
+      local.cells.push_back(localIndex[globalId]);
+    }
+    componentConstraints[cid].push_back(std::move(local));
+  }
+
+  std::vector<ComponentResult> components;
+  components.reserve(componentCells.size());
+  for (int cid = 0; cid < static_cast<int>(componentCells.size()); ++cid) {
+    components.push_back(enumerateComponent(componentCells[cid], componentConstraints[cid]));
+  }
+
+  const int componentCount = static_cast<int>(components.size());
+  std::vector<std::vector<long double>> prefix(componentCount + 1);
+  std::vector<std::vector<long double>> suffix(componentCount + 1);
+  prefix[0] = {1.0L};
+  for (int i = 0; i < componentCount; ++i) {
+    prefix[i + 1] = convolve(prefix[i], components[i].ways);
+  }
+  suffix[componentCount] = {1.0L};
+  for (int i = componentCount - 1; i >= 0; --i) {
+    suffix[i] = convolve(components[i].ways, suffix[i + 1]);
+  }
+
+  const std::vector<long double>& waysFrontier = prefix[componentCount];
+  std::vector<long double> interiorWays(waysFrontier.size(), 0.0L);
+  long double totalWays = 0.0L;
+  for (int k = 0; k < static_cast<int>(waysFrontier.size()); ++k) {
+    const int remaining = b - k;
+    if (remaining < 0 || remaining > interiorCount) {
+      continue;
+    }
+    interiorWays[k] = choose(interiorCount, remaining);
+    totalWays += waysFrontier[k] * interiorWays[k];
+  }
+
+  if (totalWays == 0.0L) {
+    std::printf("0.000\n");
+    return 0;
+  }
+
+  std::vector<long double> bombWays(hiddenCount, 0.0L);
+  for (int cid = 0; cid < componentCount; ++cid) {
+    const std::vector<long double> waysOther = convolve(prefix[cid], suffix[cid + 1]);
+    const auto& component = components[cid];
+    for (int localCell = 0; localCell < static_cast<int>(component.globalHiddenIds.size()); ++localCell) {
+      long double totalBombWays = 0.0L;
+      for (int kc = 0; kc < static_cast<int>(component.ways.size()); ++kc) {
+        const long double localBombWays = component.bombCount[localCell][kc];
+        if (localBombWays == 0.0L) {
+          continue;
+        }
+        for (int ko = 0; ko < static_cast<int>(waysOther.size()); ++ko) {
+          if (waysOther[ko] == 0.0L) {
+            continue;
+          }
+          const int totalFrontierBombs = kc + ko;
+          if (totalFrontierBombs >= static_cast<int>(interiorWays.size()) ||
+              interiorWays[totalFrontierBombs] == 0.0L) {
+            continue;
+          }
+          totalBombWays += localBombWays * waysOther[ko] * interiorWays[totalFrontierBombs];
+        }
+      }
+      bombWays[component.globalHiddenIds[localCell]] = totalBombWays;
+    }
+  }
+
+  if (interiorCount > 0) {
+    long double interiorBombWays = 0.0L;
+    for (int k = 0; k < static_cast<int>(waysFrontier.size()); ++k) {
+      const int remaining = b - k;
+      if (remaining <= 0 || remaining > interiorCount) {
+        continue;
+      }
+      interiorBombWays += waysFrontier[k] * choose(interiorCount - 1, remaining - 1);
+    }
+    for (int id = 0; id < hiddenCount; ++id) {
+      if (!isFrontier[id]) {
+        bombWays[id] = interiorBombWays;
+      }
+    }
+  }
+
+  long double minBombProbability = 1.0L;
+  for (int id = 0; id < hiddenCount; ++id) {
+    minBombProbability = std::min(minBombProbability, bombWays[id] / totalWays);
+  }
+
+  const long double answer = std::clamp(1.0L - minBombProbability, 0.0L, 1.0L);
+  std::printf("%.3Lf\n", answer);
+  return 0;
+}