Move headers into include folder and add single header to include all reachability algorithms

include/algorithm/breadth_first_search.h

@@ -0,0 +1,80 @@
+#ifndef BREADTH_FIRST_SEARCH_H_
+#define BREADTH_FIRST_SEARCH_H_
+
+#include <map>
+#include <set>
+#include <queue>
+
+using namespace graph;
+
+namespace algo {
+
+template<typename T>
+class BreadthFirstSearch {
+public:
+  BreadthFirstSearch() = default;
+
+  BreadthFirstSearch(std::unordered_map<T, std::unordered_set<T>> adjList)
+    : adjList(adjList) {}
+
+  // Traverse whole graph using the BFS search, and save the tree graph
+  // which is created when visiting new vertices (Breadth First Tree)
+  std::unordered_map<T, std::unordered_set<T>> execute(const T& root);
+
+  // Search if target vertex exists in graph
+  bool query(const T& root, const T& target);
+private:
+  // Represents the graph on which the algorithm will be executed
+  std::unordered_map<T, std::unordered_set<T>> adjList;
+}; 
+
+template<typename T>
+std::unordered_map<T, std::unordered_set<T>> BreadthFirstSearch<T>::execute(const T& root) {
+  std::unordered_map<T, std::unordered_set<T>> tree;
+  std::unordered_map<T, bool> visited;
+  std::queue<T> Q;
+  Q.push(root);
+  visited[root] = true;
+
+  while (!Q.empty()) {
+    const auto v = Q.front();
+    Q.pop();
+
+    for (const auto& u : adjList[v]) {
+      if (!visited[u]) {
+        visited[u] = true;
+        tree[v].insert(u);
+        tree[u];
+        Q.push(u);
+      }
+    }
+  }
+  return tree;
+}
+
+template<typename T>
+bool BreadthFirstSearch<T>::query(const T& root, const T& target) {
+  std::unordered_map<T, bool> visited;
+  std::queue<T> Q;
+  Q.push(root);
+  visited[root] = true;
+
+  while (!Q.empty()) {
+    const auto v = Q.front();
+    Q.pop();
+
+    if (v == target) return true;
+
+    for (const auto& u : adjList[v]) {
+      if (!visited[u]) {
+        visited[u] = true;
+        Q.push(u);
+      }
+    }
+  }
+  return false;
+}
+
+} // namespace algo
+
+#endif

include/algorithm/decremental_reachability.h

@@ -0,0 +1,28 @@
+#ifndef DECREMENTAL_REACHABILITY_H_
+#define DECREMENTAL_REACHABILITY_H_
+
+#include "graph/digraph.h"
+
+namespace algo {
+
+template<typename T>
+class DecrementalReachability {
+public:
+  virtual ~DecrementalReachability() {};
+
+  // This method is implemented and executed by all roditty and zwick
+  // algorithms, it constructs the data structures used in other operations
+  virtual void init() =0;
+
+  // Answer if vertex v is reachable from vertex u
+  virtual bool query(const T& u, const T& v) =0;
+
+  // Remove edge e(u,v) and maintain the transitive closure matrix
+  virtual void remove(const T& u, const T& v) =0;
+protected:
+  Digraph<T> G;
+};
+
+}; // namespace algo
+
+#endif

include/algorithm/dynamic_reachability.h

@@ -0,0 +1,18 @@
+#ifndef DYNAMIC_REACHABILITY_H_
+#define DYNAMIC_REACHABILITY_H_
+
+#include "algorithm/decremental_reachability.h"
+
+namespace algo {
+
+template<typename T>
+class DynamicReachability : public DecrementalReachability<T> {
+
+  // Insert edge e(u,v) and maintain the transitive closure matrix
+  virtual void insert(const T& u, const T& v) =0;
+};
+
+
+} // namespace algo
+
+#endif

include/algorithm/frigioni.h

@@ -0,0 +1,202 @@
+#ifndef FRIGIONI_H_
+#define FRIGIONI_H_
+
+#include "algorithm/decremental_reachability.h"
+#include "algorithm/roditty_zwick.h"
+
+#include <utility>
+
+using namespace graph;
+
+namespace algo {
+
+template<typename T>
+class Frigioni : public DecrementalReachability<T> {
+public:
+  Frigioni() = default;
+
+  Frigioni(Digraph<T> G) { this->G = G; }
+
+  // Initialize the decremental maintenance data structure for general graphs
+  void init() override;
+
+  // Execute reachability query q(u, v) from vertex u to vertex v
+  // in O(1) using the transitive closure matrix
+  bool query(const T& u, const T& v) override;
+
+  // Delete edge e(u, v)
+  void remove(const T& u, const T& v) override;
+
+  // Delete set of edges and explicitely maintain the transitive closure
+  void remove(const std::vector<std::pair<T, T>>& edges);
+private:
+  // Transitive closure matrix, used to answer reachability queries in O(1)
+  std::unordered_map<T, std::unordered_map<T, bool>> TC;
+
+  // For each SCC, store a reachability tree created as a BFS tree
+  std::unordered_map<T, BreadthFirstTree<T>> RT;
+
+  // For each SCC, store collections of incoming, outgoing and internal edges
+  struct Edges {
+    std::set<std::pair<T, T>> in;
+    std::set<std::pair<T, T>> inc;
+    std::set<std::pair<T, T>> out;
+  };
+  std::unordered_map<T, Edges> E;
+  
+  // Decremental maintenance of strongly connected components
+  RodittyZwick<T> rodittyZwick;
+};
+
+template<typename T>
+void Frigioni<T>::init() {
+  auto SCCs = Tarjan<T>(this->G.adjList).execute();
+  rodittyZwick = RodittyZwick<T>(this->G);
+  rodittyZwick.init();
+
+  for (auto& scc : SCCs) {
+    RT[scc.id] = BreadthFirstTree<T>(this->G, scc.id);
+    for (const auto& u : this->G.vertices()) {
+      for (const auto& v : this->G.adjList[u]) {
+        if (scc.member(u)) {
+          if (scc.member(v))
+            E[scc.id].in.insert(std::make_pair(u, v));  
+          else
+            E[scc.id].out.insert(std::make_pair(u, v));
+        } else if (scc.member(v)) {
+          E[scc.id].inc.insert(std::make_pair(u, v));
+        }
+        TC[u][v] = false;
+      }
+    }
+  }
+  for (auto& scc : SCCs) {
+    for (const auto& u : scc.vertices()) {
+      for (const auto& v : RT[scc.id].vertices()) {
+        TC[u][v] = true;
+      }
+    }
+  }
+}
+
+template<typename T>
+bool Frigioni<T>::query(const T& u, const T& v) {
+  return TC[u][v];
+}
+
+template<typename T>
+void Frigioni<T>::remove(const T& u, const T& v) {
+  
+}
+
+template<typename T>
+void Frigioni<T>::remove(const std::vector<std::pair<T, T>>& edges) {
+  std::vector<std::pair<T, T>> Eint;
+  std::vector<std::pair<T, T>> Eext;
+
+  for (const auto& [u, v] : edges) {
+    if (rodittyZwick.query(u, v))
+      Eint.push_back({ u, v });
+    else
+      Eext.push_back({ u, v });
+  }
+
+  std::unordered_map<T, std::stack<T>> H;
+  for (const auto& [u, v] : Eint) {
+    this->G.remove(u, v);
+    rodittyZwick.remove(u, v);
+
+    if (!rodittyZwick.query(u, v)) {
+      TC[u][v] = false;
+      auto C = rodittyZwick.getSCCs();
+
+      E[C[u].id].inc.erase({ u, v });
+      E[C[u].id].out.erase({ u, v });
+      E[C[u].id].in.erase({ u, v });
+      auto D = E[C[u].id];
+      E.erase(C[u].id);
+
+      for (const auto& id : std::views::keys(C)) {
+        if (RT[id].contains(id, v)) {
+          if (E[C[v].id].inc.size() > 1)
+            H[id].push(C[v].id);
+          else {
+            for (const auto& w : C[id].vertices())
+              TC[w][v] = false;
+            for (const auto& c : E[v].out)
+              H[id].push(C[c.second].id);
+          }
+          RT[id].adjList[u].erase(v);
+        }
+      }
+
+      for (const auto& [w, z] : D.inc) {
+        E[C[z].id].inc.insert({ w, z });
+      }
+
+      for (const auto& [w, z] : D.out) {
+        E[C[w].id].out.insert({ w, z });
+      }
+
+      for (const auto& [w, z] : D.in) {
+        if (C[w] == C[z])
+          E[C[w].id].in.insert({ w, z });
+        else {
+          E[C[w].id].out.insert({ w, z });
+          E[C[z].id].in.insert({ w, z });
+        }
+      }
+    } else {
+      E[u].in.erase({ u, v });
+    }
+  }
+
+  for (const auto& [u, v] : Eext) {
+    this->G.remove(u, v);
+    auto C = rodittyZwick.getSCCs();
+
+    for (const auto& id : std::views::keys(C)) {
+      if (RT[id].contains(id, v)) {
+        if (E[C[v].id].inc.size() > 1)
+          H[id].push(C[v].id);
+        else {
+          for (const auto& w : C[id].vertices())
+            TC[w][v] = false;
+          for (const auto& c : E[v].out)
+            H[id].push(C[c.second].id);
+        }
+        RT[id].adjList[u].erase(v);
+      }
+    }
+  }
+
+  auto C = rodittyZwick.getSCCs();
+  for (const auto& id : std::views::keys(rodittyZwick.getSCCs())) {
+    while (H[id].size() > 0) {
+      const auto& h = H[id].top();
+      bool found = false;
+      
+      for (const auto& [u, v] : E[h].inc) {
+        if (RT[id].contains(id, u)) {
+          RT[id].adjList[u].insert(h);
+          found = true;
+          break;
+        }
+      }
+      H[id].pop();
+      if (!found) {
+        auto C = rodittyZwick.getSCCs();
+        for (const auto& w : C[id].vertices())
+          TC[w][h] = false;
+
+        for (const auto& [u, v] : E[h].out) {
+          H[id].push(v);
+        }
+      }
+    }
+  }
+}
+
+} // namespace algo
+
+#endif

include/algorithm/henzinger_king.h

@@ -0,0 +1,113 @@
+#ifndef HENZINGER_KING_H_
+#define HENZINGER_KING_H_
+
+#include "algorithm/dynamic_reachability.h"
+#include "algorithm/frigioni.h"
+
+using namespace graph;
+
+constexpr int threshold = 5;
+
+namespace algo {
+
+template<typename T>
+class HenzingerKing : public DynamicReachability<T> {
+public:
+  HenzingerKing() = default;
+
+  HenzingerKing(Digraph<T> G) { this->G = G; }
+
+  // Initialize decremental maintenance data structure and the empty set S
+  void init() override;
+
+  // Execute query q(u, v) from vertex u to vertex v by querying the
+  // decremental reachability data structure at the start of the phase
+  // and then if necessary check the set S
+  bool query(const T& u, const T& v) override;
+
+  // Delete edge e(u, v)
+  void remove(const T& u, const T& v) override;
+
+  // Remove collection of edges from the decremental maintenance data structure
+  // and for every vertex in set S, rebuilt reachability trees from scratch
+  void remove(const std::vector<std::pair<T,T>>& edges);
+
+  // Add edge e(u, v)
+  void insert(const T& u, const T& v) override;
+
+  // Insert collection of edges in set S, if threshold is reached re-initialize
+  // algorithm, otherwise construct reachability trees for the vertex that is
+  // the center-of-insertions
+  void insert(const T& c, const std::vector<T>& vertices);
+private:
+  // Decremental maintenance data structure
+  Frigioni<T> frigioni;
+
+  // Collection of vertices that have been centers of insertions in this phase
+  std::set<T> S;
+
+  // Maintain in-out bfs trees
+  std::unordered_map<T, BreadthFirstTree<T>> In;
+  std::unordered_map<T, BreadthFirstTree<T>> Out;
+};
+
+template<typename T>
+void HenzingerKing<T>::init() {
+  S.clear();
+  frigioni = Frigioni<T>(this->G);
+  frigioni.init();
+}
+
+template<typename T>
+bool HenzingerKing<T>::query(const T& u, const T& v) {
+  if (frigioni.query(u, v))
+    return true;
+
+  return std::any_of(S.begin(), S.end(),
+                     [&](const T& w) {
+                       return In[w].adjList.contains(u) &&
+                         Out[w].adjList.contains(v);
+                     });
+}
+
+template<typename T>
+void HenzingerKing<T>::remove(const T& u, const T& v) {
+  this->G.remove(u, v);
+}
+
+template<typename T>
+void HenzingerKing<T>::remove(const std::vector<std::pair<T, T>>& edges) {
+  for (const auto& [u, v]: edges) {
+    remove(u, v);
+  }
+  frigioni.remove(edges);
+
+  for (const auto& w : S) {
+    In[w] = BreadthFirstTree(this->G.reverse(), w);
+    Out[w] = BreadthFirstTree(this->G, w);
+  }
+}
+
+template<typename T>
+void HenzingerKing<T>::insert(const T& u, const T& v) {
+  this->G.insert(u, v);
+}
+
+template<typename T>
+void HenzingerKing<T>::insert(const T& c, const std::vector<T>& vertices) {
+  for (const auto& w : vertices)
+    insert(c, w);
+
+  S.insert(c);
+  if (S.size() > threshold) {
+    init();
+    return;
+  }
+
+  In[c] = BreadthFirstTree(this->G.reverse(), c);
+  Out[c] = BreadthFirstTree(this->G, c);
+}
+
+} // namespace algo
+
+#endif

include/algorithm/italiano.h

@@ -0,0 +1,118 @@
+#ifndef ITALIANO_H_
+#define ITALIANO_H_
+
+#include "algorithm/decremental_reachability.h"
+#include "graph/breadth_first_tree.h"
+
+#include <stack>
+
+using namespace graph;
+
+namespace algo {
+
+template<typename T>
+class Italiano : public DecrementalReachability<T> {
+public:
+  Italiano() = default;
+
+  Italiano(Digraph<T> G) { this->G = G; }
+
+  // Initialize the decremental maintenance data structure for DAGs
+  void init() override;
+
+  // Execute reachability query q(u, v) from vertex u to vertex v
+  // in O(1) using the transitive closure matrix
+  bool query(const T& u, const T& v) override;
+
+  // Delete edge e(u, v) and explicitely maintain the transitive closure
+  void remove(const T& u, const T& v) override;
+private:
+  // Transitive closure matrix
+  std::unordered_map<T, std::unordered_map<T, bool>> TC;
+  
+  // For each vertex, store a reachability tree created as a BFS tree
+  std::unordered_map<T, BreadthFirstTree<T>> RT;
+
+  // For each vertex, store collections of its incoming and outgoing edges
+  struct Edges {
+    std::set<T> inc;
+    std::set<T> out;
+  };
+  std::unordered_map<T, Edges> E;
+
+  // Repair reachability trees after edge deletions
+  void repairTrees(std::unordered_map<T, std::stack<T>>& H);
+};
+
+template<typename T>
+void Italiano<T>::init() {
+  for (const auto& u : this->G.vertices()) {
+    for (const auto& v : this->G.adjList[u]) {
+      E[v].inc.insert(u);
+      E[u].out.insert(v);
+      TC[u][v] = false;
+    }
+    RT[u] = BreadthFirstTree<T>(this->G, u);
+    TC[u][u] = true;
+    for (const auto& v : RT[u].vertices())
+      TC[u][v] = true;
+  }
+}
+
+template<typename T>
+bool Italiano<T>::query(const T& u, const T& v) {
+  return TC[u][v];
+}
+
+template<typename T>
+void Italiano<T>::remove(const T& u, const T& v) {
+  if (!this->G.adjList[u].contains(v)) return;
+
+  std::unordered_map<T, std::stack<T>> H;
+  for (const auto& w : this->G.vertices()) {
+    if (RT[w].contains(u, v)) {
+      if (E[v].inc.size() > 1)
+        H[w].push(v);
+      else {
+        TC[w][v] = false;
+        for (const auto& c : E[v].out)
+          H[w].push(c);
+      }
+      RT[w].adjList[u].erase(v);
+    }
+  }
+  E[u].out.erase(v);
+  E[v].inc.erase(u);
+  this->G.remove(u, v);
+
+  repairTrees(H);
+}
+
+template<typename T>
+void Italiano<T>::repairTrees(std::unordered_map<T, std::stack<T>>& H) {
+  for (const auto& z : this->G.vertices()) {
+    while (H[z].size() > 0) {
+      const auto& h = H[z].top();
+      bool found = false;
+
+      for (const auto& i : E[h].inc) {
+        if (RT[z].contains(z, i)) {
+          RT[z].adjList[i].insert(h);
+          found = true;
+          break;
+        }
+      }
+      H[z].pop();
+      if (!found) {
+        TC[z][h] = false;
+        for (const auto& o : E[h].out) {
+          H[z].push(o);
+        }
+      }
+    }
+  }
+}
+
+} // namespace algo
+
+#endif

include/algorithm/king.h

@@ -0,0 +1,73 @@
+#ifndef KING_H_
+#define KING_H_
+
+#include "algorithm/dynamic_reachability.h"
+#include "algorithm/italiano.h"
+
+using namespace graph;
+
+namespace algo {
+
+template<typename T>
+class King : public DynamicReachability<T> {
+public:
+  King() = default;
+
+  King(Digraph<T> G) { this->G = G; }
+
+  // Initialize decremental maintenance data structures for DAGs for each
+  // vertex's in and out reachability trees
+  void init() override;
+
+  // Execute reachability query q(u, v) from vertex u to vertex v
+  // in O(n) using the stored decremental maintenance data structure for DAGs
+  bool query(const T& u, const T& v) override;
+
+  // Delete edge e(u, v) from all reachability trees and update each one of
+  // them using the decremental reachability algorithm for DAGs
+  void remove(const T& u, const T& v) override;
+
+  // Insert edge e(u, v) by reconstructing all reachability trees
+  void insert(const T& u, const T& v) override;
+private:
+  // Connect each reachabiliy tree with decremental maintenance data structure
+  std::unordered_map<T, Italiano<T>> In;
+  std::unordered_map<T, Italiano<T>> Out;
+};
+
+template<typename T>
+void King<T>::init() {
+  for (const auto& u : this->G.vertices()) {
+    In[u] = Italiano<T>(BreadthFirstTree<T>(this->G.reverse(), u));
+    In[u].init();
+    Out[u] = Italiano<T>(BreadthFirstTree<T>(this->G, u));
+    Out[u].init();
+  }
+}
+
+template<typename T>
+bool King<T>::query(const T& u, const T& v) {
+  return std::any_of(this->G.vertices().begin(), this->G.vertices().end(),
+                     [&](const T& w) {
+                       return In[w].query(w, u) && Out[w].query(w, v);
+                     });
+}
+
+template<typename T>
+void King<T>::remove(const T& u, const T& v) {
+  this->G.remove(u, v);
+  for (const auto& w : this->G.vertices()) {
+    In[w].remove(v, u);
+    Out[w].remove(u, v);
+  }
+}
+
+template<typename T>
+void King<T>::insert(const T& u, const T& v) {
+  this->G.insert(u, v);
+  init();
+}
+
+} // namespace algo
+
+#endif

include/algorithm/roditty_zwick.h

@@ -0,0 +1,96 @@
+#ifndef RODITTY_ZWICK_SCC_H_
+#define RODITTY_ZWICK_SCC_H_
+
+#include "algorithm/decremental_reachability.h"
+#include "algorithm/tarjan.h"
+#include "graph/breadth_first_tree.h"
+
+using namespace graph;
+
+namespace algo {
+
+template<typename T>
+class RodittyZwick : public DecrementalReachability<T> {
+public:
+  RodittyZwick() = default;
+
+  RodittyZwick(Digraph<T> G) { this->G = G; }
+
+  //
+  void init() override;
+
+  //
+  void findSCC();
+
+  // Return true if u and v are in the same SCC
+  bool query(const T& u, const T& v) override;
+
+  // Remove edge (u,v) and update A accordingly for fast checking query
+  void remove(const T& u, const T& v) override;
+
+  std::unordered_map<T, SCC<T>> getSCCs() { return C; }
+private:
+  // Array used to answer strong connectivity queries in O(1) time
+  std::unordered_map<T, T> A;
+
+  // Connect each representative with its SCC
+  std::unordered_map<T, SCC<T>> C;
+
+  // Maintain in-out bfs trees
+  std::unordered_map<T, BreadthFirstTree<T>> In;
+  std::unordered_map<T, BreadthFirstTree<T>> Out;
+};
+
+template<typename T>
+void RodittyZwick<T>::init() {
+  findSCC();
+}
+
+template<typename T>
+void RodittyZwick<T>::findSCC() {
+  auto SCCs = Tarjan<T>(this->G.adjList).execute();
+  
+  for (auto& SCC : SCCs) {
+    const auto& w = SCC.id;
+
+    for (const auto& v : SCC.vertices())
+      A[v] = w;
+    
+    Out[w] = BreadthFirstTree<T>(SCC, w);
+    In[w] = BreadthFirstTree<T>(SCC.reverse(), w);
+
+    C[w] = SCC;
+  }
+}
+
+template<typename T>
+bool RodittyZwick<T>::query(const T& u, const T& v) {
+  return A[u] == A[v];
+}
+
+template<typename T>
+void RodittyZwick<T>::remove(const T& u, const T& v) {
+  const auto& w = A[u];
+  C[w].remove(u, v);
+  this->G.remove(u, v);
+
+  // If u and v are not in the same SCC, do nothing
+  if (A[u] != A[v]) return;
+
+  // If edge (u,v) is not contained in both inTree and outTree do nothing
+  if (!In[w].adjList[u].contains(v) &&
+      !Out[w].adjList[u].contains(v))
+    return;
+
+  // Update In(w) and Out(w)
+  Out[w] = BreadthFirstTree<T>(C[w], w);
+  In[w] = BreadthFirstTree<T>(C[w].reverse(), w);
+
+  // If a SCC is broken, compute all SCCs again
+  if (!In[w].adjList.count(u) || !Out[w].adjList.count(v))
+    findSCC();
+}
+
+}; // namespace algo
+
+#endif

include/algorithm/tarjan.h

@@ -0,0 +1,85 @@
+#ifndef TARJAN_H_
+#define TARJAN_H_
+
+#include "graph/scc.h"
+
+#include <stack>
+#include <vector>
+#include <ranges>
+
+using namespace graph;
+
+namespace algo {
+
+template<typename T>
+class Tarjan {
+public:
+  Tarjan() = default;
+
+  Tarjan(std::unordered_map<T, std::unordered_set<T>> adjList) : adjList(adjList) {}
+
+  //
+  auto execute();
+
+  //
+  void strongConnect(const T& u);
+private:
+  std::unordered_map<T, std::unordered_set<T>> adjList;
+  std::stack<T> S;
+  std::int16_t index = 0;
+  std::vector<SCC<T>> SCCs;
+  T cid;
+
+  struct Vertex {
+    int index = -1;
+    int lowlink = -1;
+    bool onStack = false;
+  };
+  std::unordered_map<T, Vertex> vmap;
+};
+
+template<typename T>
+void Tarjan<T>::strongConnect(const T& u) {
+  vmap[u].index = vmap[u].lowlink = index++;
+  S.push(u);
+  vmap[u].onStack = true;
+
+  for (const auto& w : adjList[u]) {
+    if (vmap[w].index == -1) {
+      strongConnect(w);
+      vmap[u].lowlink = std::min(vmap[u].lowlink, vmap[w].lowlink);
+    } else if (vmap[w].onStack) {
+      vmap[u].lowlink = std::min(vmap[u].lowlink, vmap[w].index);
+    }
+  }
+
+  // If u is a root node, pop the stack and generate an SCC
+  if (vmap[u].lowlink == vmap[u].index) {
+    std::unordered_map<T, std::unordered_set<T>> scc;
+    bool finished = false;
+    cid = S.top();
+
+    do {
+      const auto w = S.top();
+      S.pop();
+      vmap[w].onStack = false;
+      scc[w] = adjList[w];
+      finished = (w == u);
+    } while (!finished);
+    
+    SCCs.push_back({ scc, static_cast<T>(cid) });
+  }
+}
+
+template<typename T>
+auto Tarjan<T>::execute() {
+  for (const auto& u : std::views::keys(adjList)) {
+    if (vmap[u].index == -1)
+      strongConnect(u);
+  }
+  return SCCs;
+}
+
+} // namespace algo
+
+#endif