graph.h

Go to the documentation of this file.

//
// Project NuriKit - Copyright 2023 SNU Compbio Lab.
// SPDX-License-Identifier: Apache-2.0
//
 
#ifndef NURI_CORE_GRAPH_GRAPH_H_
#define NURI_CORE_GRAPH_GRAPH_H_

#include <algorithm>
#include <initializer_list>
#include <iterator>
#include <queue>
#include <type_traits>
#include <utility>
#include <vector>
 
#include <absl/algorithm/container.h>
#include <absl/base/attributes.h>
#include <absl/base/optimization.h>
#include <absl/container/flat_hash_set.h>
#include <absl/log/absl_check.h>
#include <absl/log/absl_log.h>
#include <boost/iterator/iterator_facade.hpp>
 
#include "nuri/core/container/container_ext.h"
#include "nuri/core/graph/traits.h"
#include "nuri/iterator.h"
#include "nuri/meta.h"
#include "nuri/utils.h"
 
namespace nuri {
namespace internal {
  template <class Derived, class GT, class DT, bool is_const>
  class DataIteratorBase
      : public ProxyIterator<Derived, DT, std::random_access_iterator_tag, int> {
    using Traits =
        std::iterator_traits<typename DataIteratorBase::iterator_facade_>;
 
  public:
    using parent_type = const_if_t<is_const, GT>;
 
    using iterator_category = typename Traits::iterator_category;
    using value_type = typename Traits::value_type;
    using difference_type = typename Traits::difference_type;
    using pointer = typename Traits::pointer;
    using reference = typename Traits::reference;
 
    constexpr DataIteratorBase() noexcept = default;
 
    constexpr DataIteratorBase(parent_type &graph,
                               difference_type index) noexcept
        : graph_(&graph), index_(index) { }
 
  protected:
    using Parent = DataIteratorBase;
 
    template <class Other,
              std::enable_if_t<!std::is_same_v<Derived, Other>, int> = 0>
    constexpr DataIteratorBase(const Other &other) noexcept
        : graph_(other.graph_), index_(other.index_) { }
 
    template <class Other,
              std::enable_if_t<!std::is_same_v<Derived, Other>, int> = 0>
    constexpr DataIteratorBase &operator=(const Other &other) noexcept {
      graph_ = other.graph_;
      index_ = other.index_;
      return *this;
    }
 
    constexpr parent_type &graph() const noexcept { return *graph_; }
 
    constexpr difference_type index() const noexcept { return index_; }
 
    template <class Other,
              std::enable_if_t<std::is_convertible_v<Other, Derived>, int> = 0>
    constexpr bool equal(const Other &other) const noexcept {
      return index_ == other.index_;
    }
 
    template <class Other,
              std::enable_if_t<std::is_convertible_v<Other, Derived>, int> = 0>
    constexpr difference_type distance_to(const Other &other) const noexcept {
      return other.index_ - index_;
    }
 
    void increment() noexcept { ++index_; }
    void decrement() noexcept { --index_; }
    void advance(difference_type n) noexcept { index_ += n; }
 
  private:
    template <class, class, class, bool>
    friend class DataIteratorBase;
 
    friend class boost::iterator_core_access;
 
    parent_type *graph_;
    difference_type index_;
  };
 
  template <class GT, bool is_const>
  class NodeWrapper;
 
  template <class GT, bool is_const>
  class AdjWrapper {
  public:
    using edge_value_type = const_if_t<is_const, typename GT::edge_data_type>;
    using parent_type = const_if_t<is_const, GT>;
 
    template <bool other_const>
    using Other = AdjWrapper<GT, other_const>;
 
    constexpr AdjWrapper(parent_type &graph, int src, int dst, int eid) noexcept
        : src_(src), dst_(dst), eid_(eid), graph_(&graph) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr AdjWrapper(const Other<other_const> &other) noexcept
        : src_(other.src_), dst_(other.dst_), eid_(other.eid_),
          graph_(other.graph_) { }
 
    constexpr auto src() const noexcept { return graph_->node(src_); }
    constexpr auto dst() const noexcept { return graph_->node(dst_); }
 
    constexpr int eid() const noexcept { return eid_; }
    constexpr edge_value_type &edge_data() const noexcept {
      return graph_->edge_data(eid_);
    }
 
    constexpr Other<true> as_const() const noexcept { return *this; }
 
  private:
    template <class, bool>
    friend class AdjWrapper;
 
    friend GT;
 
    int src_;
    int dst_;
    int eid_;
    parent_type *graph_;
  };
 
  template <class GT, bool is_const>
  class AdjIterator
      : public DataIteratorBase<AdjIterator<GT, is_const>, GT,
                                AdjWrapper<GT, is_const>, is_const> {
    using Base = typename AdjIterator::Parent;
 
  public:
    using typename Base::difference_type;
    using typename Base::iterator_category;
    using typename Base::pointer;
    using typename Base::reference;
    using typename Base::value_type;
 
    using typename Base::parent_type;
 
    template <bool other_const>
    using Other = AdjIterator<GT, other_const>;
 
    constexpr AdjIterator() noexcept = default;
 
    constexpr AdjIterator(parent_type &graph, difference_type idx,
                          difference_type nid) noexcept
        : Base(graph, idx), nid_(nid) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr AdjIterator(const Other<other_const> &other) noexcept
        : Base(other), nid_(other.nid_) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr AdjIterator &operator=(const Other<other_const> &other) noexcept {
      Base::operator=(other);
      nid_ = other.nid_;
      return *this;
    }
 
    constexpr bool begin() const noexcept {
      return Base::equal(this->graph().adj_begin(nid_));
    }
 
    constexpr bool end() const noexcept {
      return Base::equal(this->graph().adj_end(nid_));
    }
 
  private:
    friend Base;
 
    friend class boost::iterator_core_access;
 
    template <class, bool>
    friend class AdjIterator;
 
    template <
        class Other,
        std::enable_if_t<std::is_convertible_v<Other, AdjIterator>, int> = 0>
    constexpr bool equal(const Other &other) const noexcept {
      return nid_ == other.nid_ && Base::equal(other);
    }
 
    constexpr reference dereference() const noexcept {
      return this->graph().adjacent(nid_, this->index());
    }
 
    int nid_;
  };
 
  template <class GT, bool is_const>
  class NodeWrapper {
  public:
    using DT = typename GT::node_data_type;
 
    using parent_type = const_if_t<is_const, GT>;
    using value_type = const_if_t<is_const, DT>;
 
    using adjacency_iterator =
        std::conditional_t<is_const, typename GT::const_adjacency_iterator,
                           typename GT::adjacency_iterator>;
 
    template <bool other_const>
    using Other = NodeWrapper<GT, other_const>;
 
    constexpr NodeWrapper(int nid, parent_type &graph) noexcept
        : nid_(nid), graph_(&graph) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr NodeWrapper(const Other<other_const> &other) noexcept
        : nid_(other.nid_), graph_(other.graph_) { }
 
    constexpr int id() const noexcept { return nid_; }
 
    constexpr value_type &data() const noexcept {
      return graph_->node_data(nid_);
    }
 
    constexpr int degree() const noexcept { return graph_->degree(nid_); }
 
    adjacency_iterator begin() const noexcept {
      return graph_->adj_begin(nid_);
    }
 
    adjacency_iterator end() const noexcept { return graph_->adj_end(nid_); }
 
    AdjWrapper<GT, is_const> neighbor(int idx) const noexcept {
      return graph_->adjacent(nid_, idx);
    }
 
    AdjWrapper<GT, is_const> operator[](int idx) const noexcept {
      return neighbor(idx);
    }
 
    adjacency_iterator find_adjacent(int aid) const noexcept {
      return graph_->find_adjacent(nid_, aid);
    }
 
    adjacency_iterator
    find_adjacent(NodeWrapper<GT, true> node) const noexcept {
      return find_adjacent(node.id());
    }
 
    constexpr Other<true> as_const() const noexcept { return *this; }
 
  private:
    template <class, bool>
    friend class NodeWrapper;
 
    int nid_;
    parent_type *graph_;
  };
 
  template <class GT, bool is_const>
  class NodeIterator
      : public DataIteratorBase<NodeIterator<GT, is_const>, GT,
                                NodeWrapper<GT, is_const>, is_const> {
    using Base = typename NodeIterator::Parent;
 
  public:
    using typename Base::difference_type;
    using typename Base::iterator_category;
    using typename Base::pointer;
    using typename Base::reference;
    using typename Base::value_type;
 
    using Base::Base;
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr NodeIterator(const NodeIterator<GT, other_const> &other) noexcept
        : Base(other) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr NodeIterator &
    operator=(const NodeIterator<GT, other_const> &other) noexcept {
      Base::operator=(other);
      return *this;
    }
 
  private:
    friend Base;
 
    friend class boost::iterator_core_access;
 
    template <class, bool other_const>
    friend class NodeIterator;
 
    constexpr reference dereference() const noexcept {
      return this->graph().node(this->index());
    }
  };
 
  template <class GT, bool is_const>
  class EdgeWrapper {
  public:
    using DT = typename GT::edge_data_type;
 
    using parent_type = const_if_t<is_const, GT>;
    using value_type = const_if_t<is_const, DT>;
 
    template <bool other_const>
    using Other = EdgeWrapper<GT, other_const>;
 
    constexpr EdgeWrapper(int eid, parent_type &graph) noexcept
        : eid_(eid), graph_(&graph) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr EdgeWrapper(const Other<other_const> &other) noexcept
        : eid_(other.eid_), graph_(other.graph_) { }
 
    constexpr int id() const noexcept { return eid_; }
 
    constexpr NodeWrapper<GT, is_const> src() const noexcept {
      return graph_->edge_src(eid_);
    }
 
    constexpr NodeWrapper<GT, is_const> dst() const noexcept {
      return graph_->edge_dst(eid_);
    }
 
    constexpr value_type &data() const noexcept {
      return graph_->edge_data(eid_);
    }
 
    constexpr Other<true> as_const() const noexcept { return *this; }
 
  private:
    template <class, bool>
    friend class EdgeWrapper;
 
    int eid_;
    parent_type *graph_;
  };
 
  template <class GT, bool is_const>
  class EdgeIterator
      : public DataIteratorBase<EdgeIterator<GT, is_const>, GT,
                                EdgeWrapper<GT, is_const>, is_const> {
    using Base = typename EdgeIterator::Parent;
 
  public:
    using typename Base::difference_type;
    using typename Base::iterator_category;
    using typename Base::pointer;
    using typename Base::reference;
    using typename Base::value_type;
 
    using Base::Base;
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr EdgeIterator(const EdgeIterator<GT, other_const> &other) noexcept
        : Base(other) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    constexpr EdgeIterator &
    operator=(const EdgeIterator<GT, other_const> &other) noexcept {
      Base::operator=(other);
      return *this;
    }
 
  private:
    friend Base;
 
    friend class boost::iterator_core_access;
 
    template <class, bool other_const>
    friend class EdgeIterator;
 
    constexpr reference dereference() const noexcept {
      return this->graph().edge(this->index());
    }
  };
 
  template <class GT, bool is_const>
  class EdgesWrapper {
  public:
    template <class GU = GT,
              std::enable_if_t<std::is_same_v<GU, GT> && !is_const, int> = 0>
    EdgesWrapper(GT &graph): graph_(&graph) { }
 
    EdgesWrapper(const GT &graph): graph_(&graph) { }
 
    template <bool other_const>
    using Other = EdgesWrapper<GT, other_const>;
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    EdgesWrapper(const Other<other_const> &other) noexcept
        : graph_(other.graph_) { }
 
    template <bool other_const,
              std::enable_if_t<is_const && !other_const, int> = 0>
    EdgesWrapper &operator=(const Other<other_const> &other) noexcept {
      graph_ = other.graph_;
      return *this;
    }
 
    auto operator[](int id) const { return graph_->edge(id); }
 
    auto begin() const { return graph_->edge_begin(); }
    auto end() const { return graph_->edge_end(); }
 
    auto cbegin() const { return graph_->edge_cbegin(); }
    auto cend() const { return graph_->edge_cend(); }
 
    EdgesWrapper<GT, true> as_const() const { return *this; }
 
    int size() const { return graph_->num_edges(); }
 
  private:
    const_if_t<is_const, GT> *graph_;
  };
 
  template <class, bool>
  class SubEdgeWrapper;
 
  template <class T>
  struct OffsetForwarder {
    template <class U>
    T operator()(const U &edge) const {
      return T { edge.src().id() + offset, edge.dst().id() + offset,
                 edge.data() };
    }
 
    int offset;
  };
}  // namespace internal
 
using NodesErased = std::pair<std::pair<int, std::vector<int>>,
                              std::pair<int, std::vector<int>>>;
 
template <class NT, class ET, bool is_const>
class Subgraph;

template <class NT, class ET>
class Graph {
private:
  struct AdjEntry {
    int dst;
    int eid;
  };
 
public:
  struct StoredEdge {
    int src;
    int dst;
    ET data;
  };
 
  using node_data_type = NT;
  using edge_data_type = ET;
 
  using iterator = internal::NodeIterator<Graph, false>;
  using node_iterator = iterator;
  using const_iterator = internal::NodeIterator<Graph, true>;
  using const_node_iterator = const_iterator;
  using NodeRef = internal::NodeWrapper<Graph, false>;
  using ConstNodeRef = internal::NodeWrapper<Graph, true>;
 
  using edge_iterator = internal::EdgeIterator<Graph, false>;
  using const_edge_iterator = internal::EdgeIterator<Graph, true>;
  using EdgeRef = internal::EdgeWrapper<Graph, false>;
  using ConstEdgeRef = internal::EdgeWrapper<Graph, true>;
 
  using adjacency_iterator = internal::AdjIterator<Graph, false>;
  using const_adjacency_iterator = internal::AdjIterator<Graph, true>;
  using AdjRef = internal::AdjWrapper<Graph, false>;
  using ConstAdjRef = internal::AdjWrapper<Graph, true>;
 
  static_assert(std::is_same_v<typename iterator::reference, NodeRef>);
  static_assert(
      std::is_same_v<typename const_iterator::reference, ConstNodeRef>);
  static_assert(std::is_same_v<typename edge_iterator::reference, EdgeRef>);
  static_assert(
      std::is_same_v<typename const_edge_iterator::reference, ConstEdgeRef>);
  static_assert(std::is_same_v<typename adjacency_iterator::reference, AdjRef>);
  static_assert(
      std::is_same_v<typename const_adjacency_iterator::reference, ConstAdjRef>);
 
  Graph() = default;
  Graph(const Graph &) = default;
  Graph(Graph &&) noexcept = default;
  Graph &operator=(const Graph &) = default;
  Graph &operator=(Graph &&) noexcept = default;
  ~Graph() noexcept = default;

  Graph(int num_nodes): offsets_(num_nodes + 1), nodes_(num_nodes) { }

  Graph(int num_nodes, const NT &data)
      : offsets_(num_nodes + 1), nodes_(num_nodes, data) { }
 
  bool empty() const {
    // GCOV_EXCL_START
    ABSL_DCHECK(num_nodes() > 0 || num_edges() == 0)
        << "The graph is empty (num_nodes() == 0) but num_edges() == "
        << num_edges();
    // GCOV_EXCL_STOP
    return size() == 0;
  }
  int size() const { return num_nodes(); }
  int num_nodes() const { return static_cast<int>(nodes_.size()); }
 
  bool edge_empty() const {
    // GCOV_EXCL_START
    ABSL_DCHECK(num_nodes() > 0 || num_edges() == 0)
        << "The graph is empty (num_nodes() == 0) but num_edges() == "
        << num_edges();
    // GCOV_EXCL_STOP
    return num_edges() == 0;
  }
  int num_edges() const { return static_cast<int>(edges_.size()); }
 
  int degree(int id) const { return offsets_[id + 1] - offsets_[id]; }
 
  void reserve(int num_nodes) {
    nodes_.reserve(num_nodes);
    offsets_.reserve(num_nodes + 1);
  }
 
  void reserve_edges(int num_edges) {
    edges_.reserve(num_edges);
    adj_list_.reserve(num_edges * 2);
  }

  int add_node(const NT &data) {
    int id = num_nodes();
    nodes_.push_back(data);
    offsets_.push_back(offsets_[id]);
    return id;
  }

  int add_node(NT &&data) noexcept {
    int id = num_nodes();
    nodes_.push_back(std::move(data));
    offsets_.push_back(offsets_[id]);
    return id;
  }

  template <class Iterator,
            internal::enable_if_compatible_iter_t<Iterator, NT> = 0>
  int add_nodes(Iterator begin, Iterator end) {
    const int first = num_nodes();
    nodes_.insert(nodes_.end(), begin, end);
    offsets_.resize(num_nodes() + 1, offsets_.back());
    return first;
  }

  ABSL_DEPRECATED("Slow, due to CSR recompilation.")
  int add_edge(int src, int dst, const ET &data) {
    ABSL_DCHECK_NE(src, dst) << "self-loop is not allowed";
 
    int eid = num_edges();
    edges_.push_back({ src, dst, data });
    publish_edges_from(eid);
    return eid;
  }

  ABSL_DEPRECATED("Slow, due to CSR recompilation.")
  int add_edge(int src, int dst, ET &&data) noexcept {
    ABSL_DCHECK_NE(src, dst) << "self-loop is not allowed";
 
    int eid = num_edges();
    edges_.push_back({ src, dst, std::move(data) });
    publish_edges_from(eid);
    return eid;
  }

  template <class Iterator,
            internal::enable_if_compatible_iter_t<Iterator, StoredEdge> = 0>
  int add_edges(Iterator begin, Iterator end) {
    const int first = num_edges();
    edges_.insert(edges_.end(), begin, end);
    publish_edges_from(first);
    return first;
  }
 
  NodeRef operator[](int id) { return node(id); }
  ConstNodeRef operator[](int id) const { return node(id); }
 
  NodeRef node(int id) { return { id, *this }; }
  ConstNodeRef node(int id) const { return { id, *this }; }
 
  NT &node_data(int id) { return nodes_[id]; }
  const NT &node_data(int id) const { return nodes_[id]; }
 
  void clear() {
    nodes_.clear();
    edges_.clear();
    adj_list_.clear();
    offsets_ = { 0 };
  }

  NodesErased erase_nodes(const_iterator begin, const_iterator end) {
    return erase_nodes(begin, end, [](auto /* ref */) { return true; });
  }

  template <class UnaryPred>
  NodesErased erase_nodes(const_iterator begin, const_iterator end,
                          UnaryPred pred);

  template <class Iterator,
            class = internal::enable_if_compatible_iter_t<Iterator, int>>
  NodesErased erase_nodes(Iterator begin, Iterator end);
 
  iterator begin() { return { *this, 0 }; }
  iterator end() { return { *this, num_nodes() }; }
  const_iterator begin() const { return cbegin(); }
  const_iterator end() const { return cend(); }
  const_iterator cbegin() const { return { *this, 0 }; }
  const_iterator cend() const { return { *this, num_nodes() }; }
 
  EdgeRef edge(int id) { return { id, *this }; }
  ConstEdgeRef edge(int id) const { return { id, *this }; }
 
  NodeRef edge_src(int id) { return this->node(edges_[id].src); }
  ConstNodeRef edge_src(int id) const { return this->node(edges_[id].src); }
 
  NodeRef edge_dst(int id) { return this->node(edges_[id].dst); }
  ConstNodeRef edge_dst(int id) const { return this->node(edges_[id].dst); }
 
  ET &edge_data(int id) { return edges_[id].data; }
  const ET &edge_data(int id) const { return edges_[id].data; }

  edge_iterator find_edge(int src, int dst) {
    return find_edge_helper(*this, src, dst);
  }

  const_edge_iterator find_edge(int src, int dst) const {
    return find_edge_helper(*this, src, dst);
  }

  edge_iterator find_edge(ConstNodeRef src, ConstNodeRef dst) {
    return find_edge(src.id(), dst.id());
  }

  const_edge_iterator find_edge(ConstNodeRef src, ConstNodeRef dst) const {
    return find_edge(src.id(), dst.id());
  }

  void clear_edges() {
    edges_.clear();
    adj_list_.clear();
    absl::c_fill(offsets_, 0);
  }

  void erase_edge(int id) {
    const StoredEdge &edge = edges_[id];
    auto srcit = find_adjacency_entry(edge.src, edge.dst),
         dstit = find_adjacency_entry(edge.dst, edge.src);
    erase_adjacency_entry(edge.src, srcit, edge.dst, dstit);
 
    erase_edge_common(id);
  }

  bool erase_edge_between(int src, int dst);

  bool erase_edge_between(ConstNodeRef src, ConstNodeRef dst) {
    return erase_edge_between(src.id(), dst.id());
  }

  std::pair<int, std::vector<int>> erase_edges(const_edge_iterator begin,
                                               const_edge_iterator end) {
    return erase_edges(begin, end, [](auto /* ref */) { return true; });
  }

  template <class UnaryPred>
  std::pair<int, std::vector<int>> erase_edges(const_edge_iterator begin,
                                               const_edge_iterator end,
                                               UnaryPred pred);

  template <class Iterator,
            class = internal::enable_if_compatible_iter_t<Iterator, int>>
  std::pair<int, std::vector<int>> erase_edges(Iterator begin, Iterator end);
 
  internal::EdgesWrapper<Graph, false> edges() { return { *this }; }
  internal::EdgesWrapper<Graph, true> edges() const { return { *this }; }
  internal::EdgesWrapper<Graph, true> cedges() const { return { *this }; }
 
  edge_iterator edge_begin() { return { *this, 0 }; }
  edge_iterator edge_end() { return { *this, num_edges() }; }
  const_edge_iterator edge_begin() const { return edge_cbegin(); }
  const_edge_iterator edge_end() const { return edge_cend(); }
  const_edge_iterator edge_cbegin() const { return { *this, 0 }; }
  const_edge_iterator edge_cend() const { return { *this, num_edges() }; }

  adjacency_iterator find_adjacent(int src, int dst) {
    return find_adj_helper(*this, src, dst);
  }

  const_adjacency_iterator find_adjacent(int src, int dst) const {
    return find_adj_helper(*this, src, dst);
  }

  adjacency_iterator find_adjacent(ConstNodeRef src, ConstNodeRef dst) {
    return find_adjacent(src.id(), dst.id());
  }

  const_adjacency_iterator find_adjacent(ConstNodeRef src,
                                         ConstNodeRef dst) const {
    return find_adjacent(src.id(), dst.id());
  }
 
  adjacency_iterator adj_begin(int nid) { return { *this, 0, nid }; }
  adjacency_iterator adj_end(int nid) { return { *this, degree(nid), nid }; }
  const_adjacency_iterator adj_begin(int nid) const { return adj_cbegin(nid); }
  const_adjacency_iterator adj_end(int nid) const { return adj_cend(nid); }
  const_adjacency_iterator adj_cbegin(int nid) const {
    return { *this, 0, nid };
  }
  const_adjacency_iterator adj_cend(int nid) const {
    return { *this, degree(nid), nid };
  }

  template <class GraphLike>
  void merge(const GraphLike &other) {
    const int offset = size();
    reserve(offset + other.size());
 
    for (auto node: other)
      add_node(node.data());
 
    auto edges = other.edges();
    reserve_edges(num_edges() + edges.size());
 
    internal::OffsetForwarder<StoredEdge> as_stored_edge { offset };
    add_edges(internal::make_transform_iterator(edges.begin(), as_stored_edge),
              internal::make_transform_iterator(edges.end(), as_stored_edge));
  }
 
private:
  template <class, bool>
  friend class internal::NodeWrapper;
 
  template <class, bool>
  friend class internal::EdgeIterator;
 
  template <class, bool>
  friend class internal::AdjIterator;
 
  friend class Subgraph<NT, ET, false>;
  friend class Subgraph<NT, ET, true>;
 
  template <class GT>
  static internal::AdjIterator<Graph, std::is_const_v<GT>>
  find_adj_helper(GT &graph, int src, int dst) {
    auto ret = graph.adj_begin(src);
 
    for (; !ret.end(); ++ret)
      if (ret->dst().id() == dst)
        break;
 
    return ret;
  }
 
  template <class GT>
  static internal::EdgeIterator<Graph, std::is_const_v<GT>>
  find_edge_helper(GT &graph, int src, int dst) {
    if (graph.degree(src) > graph.degree(dst))
      return find_edge_helper(graph, dst, src);
 
    auto ait = find_adj_helper(graph, src, dst);
    if (ait.end())
      return graph.edge_end();
 
    return { graph, ait->eid_ };
  }
 
  std::pair<int, std::vector<int>>
  erase_nodes_common(std::vector<int> &node_keep, int first_erased_id,
                     bool erase_trailing);
 
  void erase_edge_common(int id);
 
  void erase_edges_common(std::vector<int> &edge_keep, int first_erased_id,
                          bool erase_trailing);

  void publish_edges_from(int first_eid);
 
  auto find_adjacency_entry(int src, int dst) {
    return std::find_if(adj_list_.begin() + offsets_[src],
                        adj_list_.begin() + offsets_[src + 1],
                        [dst](const AdjEntry &adj) { return adj.dst == dst; });
  }
 
  void erase_adjacency_entry(int src,
                             typename std::vector<AdjEntry>::iterator srcit,
                             int dst,
                             typename std::vector<AdjEntry>::iterator dstit) {
    if (src > dst) {
      erase_adjacency_entry(dst, dstit, src, srcit);
      return;
    }
 
    // Close the two gaps left by the removed entries via a pair of back-to-back
    // left shifts on the flat buffer.
    std::move(srcit + 1, dstit, srcit);
    std::move(dstit + 1, adj_list_.end(), dstit - 1);
    adj_list_.resize(adj_list_.size() - 2);
 
    // Slices of nodes in (src, dst] shrink by 1 (one entry gone before them);
    // slices of nodes > dst shrink by 2.
    for (int v = src + 1; v <= dst; ++v)
      --offsets_[v];
    for (int v = dst + 1, n = static_cast<int>(offsets_.size()); v < n; ++v)
      offsets_[v] -= 2;
  }
 
  AdjRef adjacent(int nid, int idx) {
    const AdjEntry &adj = adj_list_[offsets_[nid] + idx];
    return { *this, nid, adj.dst, adj.eid };
  }
 
  ConstAdjRef adjacent(int nid, int idx) const {
    const AdjEntry &adj = adj_list_[offsets_[nid] + idx];
    return { *this, nid, adj.dst, adj.eid };
  }
 
  std::vector<int> offsets_ = { 0 };
  std::vector<AdjEntry> adj_list_;
 
  std::vector<NT> nodes_;
  std::vector<StoredEdge> edges_;
};
 
/* Out-of-line definitions */
 
template <class NT, class ET>
template <class Iterator, class>
NodesErased Graph<NT, ET>::erase_nodes(Iterator begin, Iterator end) {
  // Note: the time complexity notations are only for very sparse graphs, i.e.,
  // E = O(V).
  if (begin == end) {
    return {
      { num_nodes(), {} },
      { num_edges(), {} },
    };
  }
 
  // Phase I: mark nodes for removal, O(V)
  std::vector<int> node_keep(num_nodes(), 1);
  int first_erased_id = num_nodes();
  for (auto it = begin; it != end; ++it) {
    int nid = *it;
    node_keep[nid] = 0;
    first_erased_id = std::min(first_erased_id, nid);
  }
 
  bool erase_trailing = true;
  for (int i = num_nodes() - 1; i >= first_erased_id; --i) {
    if (node_keep[i] == 1) {
      erase_trailing = false;
      break;
    }
  }
 
  auto edge_info =
      erase_nodes_common(node_keep, first_erased_id, erase_trailing);
  return {
    { erase_trailing ? first_erased_id : -1, std::move(node_keep) },
    std::move(edge_info),
  };
}
 
template <class NT, class ET>
template <class UnaryPred>
NodesErased Graph<NT, ET>::erase_nodes(const const_iterator begin,
                                       const const_iterator end,
                                       UnaryPred pred) {
  // Note: the time complexity notations are only for very sparse graphs, i.e.,
  // E = O(V).
 
  // This will also handle size() == 0 case correctly.
  if (begin >= end) {
    return {
      { num_nodes(), {} },
      { num_edges(), {} },
    };
  }
 
  // Phase I: mark nodes for removal, O(V)
  std::vector<int> node_keep(num_nodes(), 1);
  int first_erased_id = -1;
  bool erase_trailing = end == this->end();
  for (auto it = begin; it != end; ++it) {
    if (pred(*it)) {
      int nid = it->id();
 
      // GCOV_EXCL_START
      ABSL_ASSUME(nid >= 0);
      // GCOV_EXCL_STOP
 
      if (first_erased_id < 0)
        first_erased_id = nid;
 
      node_keep[nid] = 0;
    } else if (first_erased_id >= 0) {
      erase_trailing = false;
    }
  }
 
  auto edge_info =
      erase_nodes_common(node_keep, first_erased_id, erase_trailing);
  return {
    { erase_trailing ? first_erased_id : -1, std::move(node_keep) },
    std::move(edge_info),
  };
}
 
template <class NT, class ET>
std::pair<int, std::vector<int>>
Graph<NT, ET>::erase_nodes_common(std::vector<int> &node_keep,
                                  const int first_erased_id,
                                  const bool erase_trailing) {
  // Fast path 1: no node is erased
  if (first_erased_id < 0 || first_erased_id >= num_nodes())
    return { num_edges(), {} };
 
  // Phase II: erase the edges, O(V+E)
  auto ret = erase_edges(edge_begin(), edge_end(), [&](ConstEdgeRef edge) {
    return node_keep[edge.src().id()] == 0 || node_keep[edge.dst().id()] == 0;
  });
 
  // Phase III: erase the nodes & adjacencies
  if (erase_trailing) {
    // Fast path 2: only trailing nodes with incident edges; the edge erase
    // above already compacted adj_list_, so just trim nodes_/offsets_.
    ABSL_DLOG(INFO) << "resizing offset table & node list";
    // O(1) operations
    nodes_.resize(first_erased_id);
    offsets_.resize(first_erased_id + 1);
    return ret;
  }
 
  // Erase unused nodes, O(V)
  int i = 0;
  erase_if(nodes_,
           [&](const NT & /* unused */) { return node_keep[i++] == 0; });
 
  // Phase IV: update the node numbers in adjacencies and edges, O(V+E)
  mask_to_map(node_keep);
 
  for (int v = 0; v < static_cast<int>(node_keep.size()); ++v) {
    int w = node_keep[v];
    if (w < 0)
      continue;
    offsets_[w + 1] = offsets_[v + 1];
  }
  offsets_.resize(nodes_.size() + 1);
 
  for (AdjEntry &adj: adj_list_)
    adj.dst = node_keep[adj.dst];
 
  for (StoredEdge &edge: edges_) {
    edge.src = node_keep[edge.src];
    edge.dst = node_keep[edge.dst];
  }
 
  return ret;
}
 
template <class NT, class ET>
void Graph<NT, ET>::publish_edges_from(const int first_eid) {
  const int total_edges = num_edges();
  if (first_eid == total_edges)
    return;
 
  const int n = num_nodes();
 
  // Phase 1: count per-node delta from just the new edges; track v_min.
  std::vector<int> delta(n, 0);
  for (int eid = first_eid; eid < total_edges; ++eid) {
    const StoredEdge &e = edges_[eid];
    ++delta[e.src];
    ++delta[e.dst];
  }
  const int v_min = static_cast<int>(
      absl::c_find_if(delta, [](int d) { return d > 0; }) - delta.begin());
 
  // Phase 2: grow the flat buffer by the total number of new endpoints.
  const int added = 2 * (total_edges - first_eid);
  adj_list_.resize(adj_list_.size() + added);
 
  // Phase 3: back-to-front shift of existing slices; update offsets_ in place.
  //
  // Invariant on entry to iteration `v`: `shift == sum(delta[v..n-1])`, which
  // equals `shift_{v+1}` (= the offset-shift that was applied to node v+1's
  // slice). Subtracting `delta[v]` yields `shift_v` for the current node.
  int shift = added;
  for (int i = n - 1; i > v_min; --i) {
    const int old_end = offsets_[i + 1], old_start = offsets_[i],
              old_deg = old_end - old_start;
 
    offsets_[i + 1] = old_end + shift;
    shift -= delta[i];
    const int new_start = old_start + shift;
    std::move_backward(adj_list_.begin() + old_start,
                       adj_list_.begin() + old_end,
                       adj_list_.begin() + new_start + old_deg);
  }
  offsets_[v_min + 1] += shift;
 
  // GCOV_EXCL_START
  ABSL_DCHECK_EQ(shift, delta[v_min])
      << "shift mismatch: " << shift << " vs " << delta[v_min];
  // GCOV_EXCL_STOP
 
  // Phase 4: place new adjacency entries in eid order.
  for (int eid = first_eid; eid < total_edges; ++eid) {
    const StoredEdge &e = edges_[eid];
    int sd = delta[e.src]--, dd = delta[e.dst]--;
    adj_list_[offsets_[e.src + 1] - sd] = { e.dst, eid };
    adj_list_[offsets_[e.dst + 1] - dd] = { e.src, eid };
  }
}
 
template <class NT, class ET>
bool Graph<NT, ET>::erase_edge_between(int src, int dst) {
  if (degree(src) > degree(dst))
    return erase_edge_between(dst, src);
 
  auto srcit = find_adjacency_entry(src, dst);
  if (srcit == adj_list_.begin() + offsets_[src + 1])
    return false;
 
  int eid = srcit->eid;
  // NOLINTNEXTLINE(readability-suspicious-call-argument)
  auto dstit = find_adjacency_entry(dst, src);
  erase_adjacency_entry(src, srcit, dst, dstit);
 
  erase_edge_common(eid);
 
  return true;
}
 
template <class NT, class ET>
void Graph<NT, ET>::erase_edge_common(int id) {
  const int orig_edges = num_edges();
  edges_.erase(edges_.begin() + id);
 
  if (id == orig_edges - 1)
    return;
 
  for (AdjEntry &adj: adj_list_) {
    if (adj.eid > id)
      --adj.eid;
  }
}
 
template <class NT, class ET>
template <class UnaryPred>
std::pair<int, std::vector<int>>
Graph<NT, ET>::erase_edges(const_edge_iterator begin, const_edge_iterator end,
                           UnaryPred pred) {
  // This will also handle num_edges() == 0 case correctly.
  if (begin >= end)
    return { num_edges(), {} };
 
  // Phase I: mark edges for removal, O(E)
  std::vector<int> edge_keep(num_edges(), 1);
  int first_erased_id = -1;
  bool erase_trailing = end == this->edge_end();
  for (auto it = begin; it != end; ++it) {
    if (pred(*it)) {
      int eid = it->id();
      // GCOV_EXCL_START
      ABSL_ASSUME(eid >= 0);
      // GCOV_EXCL_STOP
      if (first_erased_id < 0)
        first_erased_id = eid;
 
      edge_keep[eid] = 0;
    } else if (first_erased_id >= 0) {
      erase_trailing = false;
    }
  }
 
  erase_edges_common(edge_keep, first_erased_id, erase_trailing);
  return { erase_trailing ? first_erased_id : -1, std::move(edge_keep) };
}
 
template <class NT, class ET>
template <class Iterator, class>
std::pair<int, std::vector<int>> Graph<NT, ET>::erase_edges(Iterator begin,
                                                            Iterator end) {
  if (begin == end)
    return { num_edges(), {} };
 
  // Phase I: mark edges for removal, O(E)
  std::vector<int> edge_keep(num_edges(), 1);
  int first_erased_id = num_edges();
  for (auto it = begin; it != end; ++it) {
    int eid = *it;
    edge_keep[eid] = 0;
    first_erased_id = std::min(first_erased_id, eid);
  }
 
  bool erase_trailing = true;
  for (int i = num_edges() - 1; i >= first_erased_id; --i) {
    if (edge_keep[i] == 1) {
      erase_trailing = false;
      break;
    }
  }
 
  erase_edges_common(edge_keep, first_erased_id, erase_trailing);
  return { erase_trailing ? first_erased_id : -1, std::move(edge_keep) };
}
 
template <class NT, class ET>
void Graph<NT, ET>::erase_edges_common(std::vector<int> &edge_keep,
                                       int first_erased_id,
                                       bool erase_trailing) {
  // Fast path 1: no edge is erased
  if (first_erased_id < 0 || first_erased_id >= num_edges())
    return;
 
  // Phase II: compact adj_list_ and rewrite offsets_ in place, O(V+E).
  const int n = num_nodes();
  int p = 0, old_end = 0;
  for (int i = 0; i < n; ++i) {
    const int old_start = old_end;
    old_end = offsets_[i + 1];
    for (int q = old_start; q < old_end; ++q) {
      if (edge_keep[adj_list_[q].eid] != 0)
        adj_list_[p++] = adj_list_[q];
    }
    offsets_[i + 1] = p;
  }
  adj_list_.resize(p);
 
  // Fast path 2: if only trailing edges are erased, no edge number needs to
  // be updated.
  if (erase_trailing) {
    ABSL_DLOG(INFO) << "resizing edge list";
    // O(1) operation
    edges_.resize(first_erased_id);
    return;
  }
 
  // Phase III: erase the edges, O(E)
  int i = 0;
  erase_if(edges_, [&](const StoredEdge & /* unused */) {
    return edge_keep[i++] == 0;
  });
 
  // Phase IV: update the edge numbers in adjacencies, O(E)
  mask_to_map(edge_keep);
 
  for (AdjEntry &adj: adj_list_)
    adj.eid = edge_keep[adj.eid];
}
 
namespace internal {
  // NOLINTBEGIN(readability-identifier-naming)
  template <class NT, class ET>
  struct GraphTraits<Graph<NT, ET>> {
    constexpr static bool is_degree_constant_time = true;
  };
  // NOLINTEND(readability-identifier-naming)
 
  template <class GT, bool is_const>
  class EigenNeighborIndexer {
  public:
    EigenNeighborIndexer(NodeWrapper<GT, is_const> node): node_(node) { }
 
    int size() const { return node_.degree(); }
 
    int operator[](int i) const { return node_[i].dst().id(); }
 
  private:
    NodeWrapper<GT, is_const> node_;
  };
}  // namespace internal
 
template <class GT, bool is_const>
internal::EigenNeighborIndexer<GT, is_const>
as_index(internal::NodeWrapper<GT, is_const> node) {
  return { node };
}
 
namespace internal {
  template <class NT, class ET>
  void connected_components_impl(const Graph<NT, ET> &g,
                                 absl::flat_hash_set<int> &visited,
                                 std::queue<int> &q) {
    while (!q.empty()) {
      int atom = q.front();
      q.pop();
 
      auto [_, inserted] = visited.insert(atom);
      if (!inserted)
        continue;
 
      for (auto ait = g.adj_begin(atom); !ait.end(); ++ait) {
        int dst = ait->dst().id();
        q.push(dst);
      }
    }
  }
}  // namespace internal

template <class NT, class ET>
absl::flat_hash_set<int> connected_components(const Graph<NT, ET> &g,
                                              int begin) {
  absl::flat_hash_set<int> visited;
  std::queue<int> q { { begin } };
 
  internal::connected_components_impl(g, visited, q);
 
  return visited;
}

template <class NT, class ET>
absl::flat_hash_set<int> connected_components(const Graph<NT, ET> &g, int begin,
                                              int exclude) {
  absl::flat_hash_set<int> visited { begin };
  std::queue<int> q;
 
  for (auto ait = g.adj_begin(begin); !ait.end(); ++ait) {
    int dst = ait->dst().id();
    if (dst != exclude)
      q.push(dst);
  }
 
  internal::connected_components_impl(g, visited, q);
 
  if (visited.contains(exclude))
    return {};
 
  return visited;
}
}  // namespace nuri
 
#endif /* NURI_CORE_GRAPH_GRAPH_H_ */