✨ Support for Neutral Atom Mapper

.gitignore

@@ -53,6 +53,7 @@ coverage.xml
 .hypothesis/
 .pytest_cache/
 cover/
+*.profraw
 
 # Translations
 *.mo

include/mqt-core/CircuitOptimizer.hpp

@@ -10,9 +10,6 @@
 #include <unordered_set>
 
 namespace qc {
-static constexpr std::array<qc::OpType, 10> DIAGONAL_GATES = {
-    qc::Barrier, qc::I,   qc::Z, qc::S,  qc::Sdg,
-    qc::T,       qc::Tdg, qc::P, qc::RZ, qc::RZZ};
 
 class CircuitOptimizer {
 protected:

include/mqt-core/Definitions.hpp

@@ -120,6 +120,13 @@ constexpr void hashCombine(std::size_t& hash, const std::size_t with) noexcept {
   hash = combineHash(hash, with);
 }
 
+/// Pairs do not provide a hash function by default, this is the replacement
+template <class T, class U> struct PairHash {
+  size_t operator()(const std::pair<T, U>& x) const {
+    return combineHash(std::hash<T>{}(x.first), std::hash<U>{}(x.second));
+  }
+};
+
 /**
  * @brief Function used to mark unreachable code
  * @details Uses compiler specific extensions if possible. Even if no extension

include/mqt-core/datastructures/DirectedAcyclicGraph.hpp

@@ -0,0 +1,107 @@
+#pragma once
+
+#include "Definitions.hpp"
+#include "DirectedGraph.hpp"
+
+#include <algorithm>
+#include <cassert>
+#include <iostream>
+#include <sstream>
+#include <stack>
+#include <stdexcept>
+
+namespace qc {
+
+template <class V> class DirectedAcyclicGraph final : public DirectedGraph<V> {
+protected:
+  // transitive closure matrix to detect cycles
+  std::vector<std::vector<bool>> closureMatrix;
+
+public:
+  auto addVertex(const V& v) -> void override {
+    DirectedGraph<V>::addVertex(v);
+    for (auto& row : closureMatrix) {
+      row.emplace_back(false);
+    }
+    closureMatrix.emplace_back(this->nVertices, false);
+    const auto i = this->mapping.at(v);
+    closureMatrix[i][i] = true;
+  }
+  auto addEdge(const V& u, const V& v) -> void override {
+    if (this->mapping.find(u) == this->mapping.end()) {
+      addVertex(u);
+    }
+    if (this->mapping.find(v) == this->mapping.end()) {
+      addVertex(v);
+    }
+    std::size_t const i = this->mapping.at(u);
+    std::size_t const j = this->mapping.at(v);
+    if (closureMatrix[j][i]) {
+      std::ostringstream oss;
+      oss << "Adding edge (" << u << ", " << v << ") would create a cycle.";
+      throw std::logic_error(oss.str());
+    }
+    DirectedGraph<V>::addEdge(u, v);
+    for (std::size_t k = 0; k < this->nVertices; ++k) {
+      if (closureMatrix[k][i]) {
+        closureMatrix[k][j] = true;
+      }
+      if (closureMatrix[j][k]) {
+        closureMatrix[i][k] = true;
+      }
+    }
+  }
+  [[nodiscard]] auto isReachable(const V& u, const V& v) const -> bool {
+    if (this->mapping.find(u) == this->mapping.end()) {
+      throw std::invalid_argument("Vertex u not in graph.");
+    }
+    if (this->mapping.find(v) == this->mapping.end()) {
+      throw std::invalid_argument("Vertex v not in graph.");
+    }
+    return closureMatrix[this->mapping.at(u)][this->mapping.at(v)];
+  }
+  /// Perform a depth-first search on the graph and return the nodes in a
+  /// topological order
+  [[nodiscard]] auto orderTopologically() const -> std::vector<V> {
+    std::stack<std::size_t> stack{};
+    std::vector<std::size_t> result;
+    result.reserve(this->nVertices);
+    std::vector visited(this->nVertices, false);
+    // visitedInDegree is used to count the incoming edges that have been
+    // visited already such that the resulting order of the nodes is one that
+    // satisfies a topological ordering
+    std::vector<std::size_t> visitedInDegree(this->nVertices, 0);
+    // Push nodes with 0 indegree onto the stack
+    for (std::size_t k = 0; k < this->nVertices; ++k) {
+      if (this->inDegrees[k] == 0) {
+        stack.push(k);
+        visited[k] = true;
+      }
+    }
+    // Perform DFS
+    while (!stack.empty()) {
+      const auto u = stack.top();
+      stack.pop();
+      result.emplace_back(u);
+
+      for (std::size_t k = 0; k < this->nVertices; ++k) {
+        if (this->adjacencyMatrix[u][k]) {
+          if (!visited[k]) {
+            if (++visitedInDegree[k] == this->inDegrees[k]) {
+              stack.push(k);
+              visited[k] = true;
+            }
+          }
+        }
+      }
+    }
+    // Otherwise graph has a cycle
+    assert(result.size() == this->nVertices);
+    std::vector<V> vertices;
+    vertices.reserve(this->nVertices);
+    std::transform(result.cbegin(), result.cend(), std::back_inserter(vertices),
+                   [&](const auto i) { return this->invMapping.at(i); });
+    return vertices;
+  }
+};
+} // namespace qc

include/mqt-core/datastructures/DirectedGraph.hpp

@@ -0,0 +1,138 @@
+#pragma once
+
+#include "Definitions.hpp"
+
+#include <numeric>
+#include <sstream>
+#include <unordered_set>
+
+namespace qc {
+
+/**
+ * Class representing generic directed graphs.
+ *
+ * @tparam V the type of the vertices in the graph. Must implement `operator<<`.
+ */
+template <class V> class DirectedGraph {
+  static_assert(
+      std::is_same<decltype(std::declval<std::ostream&>() << std::declval<V>()),
+                   std::ostream&>::value,
+      "V must support `operator<<`.");
+
+protected:
+  // the adjecency matrix works with indices
+  std::vector<std::vector<bool>> adjacencyMatrix;
+  // the mapping of vertices to indices in the graph are stored in a map
+  std::unordered_map<V, std::size_t> mapping;
+  // the inverse mapping is used to get the vertex from the index
+  std::vector<V> invMapping;
+  // the number of vertices in the graph
+  std::size_t nVertices = 0;
+  // the number of edges in the graph
+  std::size_t nEdges = 0;
+  // the in-degrees of the vertices in the graph
+  std::vector<std::size_t> inDegrees;
+  // the out-degrees of the vertices in the graph
+  std::vector<std::size_t> outDegrees;
+
+public:
+  virtual ~DirectedGraph() = default;
+  virtual auto addVertex(const V& v) -> void {
+    // check whether the vertex is already in the graph, if so do nothing
+    if (mapping.find(v) != mapping.end()) {
+      std::stringstream ss;
+      ss << "The vertex " << v << " is already in the graph.";
+      throw std::invalid_argument(ss.str());
+    }
+    mapping[v] = nVertices;
+    invMapping.emplace_back(v);
+    ++nVertices;
+    for (auto& row : adjacencyMatrix) {
+      row.emplace_back(false);
+    }
+    adjacencyMatrix.emplace_back(nVertices, false);
+    inDegrees.emplace_back(0);
+    outDegrees.emplace_back(0);
+  }
+  virtual auto addEdge(const V& u, const V& v) -> void {
+    if (mapping.find(u) == mapping.end()) {
+      addVertex(u);
+    }
+    if (mapping.find(v) == mapping.end()) {
+      addVertex(v);
+    }
+    const auto i = mapping.at(u);
+    const auto j = mapping.at(v);
+    if (!adjacencyMatrix[i][j]) {
+      adjacencyMatrix[i][j] = true;
+      ++outDegrees[i];
+      ++inDegrees[j];
+      ++nEdges;
+    }
+  }
+  [[nodiscard]] auto getNVertices() const -> std::size_t { return nVertices; }
+  [[nodiscard]] auto getNEdges() const -> std::size_t { return nEdges; }
+  [[nodiscard]] auto getInDegree(const V& v) const -> std::size_t {
+    if (mapping.find(v) == mapping.end()) {
+      std::stringstream ss;
+      ss << "The vertex " << v << " is not in the graph.";
+      throw std::invalid_argument(ss.str());
+    }
+    const auto i = mapping.at(v);
+    return inDegrees[i];
+  }
+  [[nodiscard]] auto getOutDegree(const V& v) const -> std::size_t {
+    if (mapping.find(v) == mapping.end()) {
+      std::stringstream ss;
+      ss << "The vertex " << v << " is not in the graph.";
+      throw std::invalid_argument(ss.str());
+    }
+    const auto i = mapping.at(v);
+    return outDegrees[i];
+  }
+  [[nodiscard]] auto getVertices() const -> std::unordered_set<V> {
+    return std::accumulate(mapping.cbegin(), mapping.cend(),
+                           std::unordered_set<V>(),
+                           [](auto& acc, const auto& v) {
+                             acc.emplace(v.first);
+                             return acc;
+                           });
+  }
+  [[nodiscard]] auto isEdge(const V& u, const V& v) const -> bool {
+    if (mapping.find(u) == mapping.end()) {
+      std::stringstream ss;
+      ss << "The vertex " << u << " is not in the graph.";
+      throw std::invalid_argument(ss.str());
+    }
+    if (mapping.find(v) == mapping.end()) {
+      std::stringstream ss;
+      ss << "The vertex " << v << " is not in the graph.";
+      throw std::invalid_argument(ss.str());
+    }
+    const auto i = mapping.at(u);
+    const auto j = mapping.at(v);
+    return adjacencyMatrix[i][j];
+  }
+  /// Outputs a string representation of the graph in the DOT format
+  [[nodiscard]] auto toString() const -> std::string {
+    std::stringstream ss;
+    ss << "digraph {\n";
+    for (const auto& [v, i] : mapping) {
+      ss << "  " << i << " [label=\"" << v << "\"];\n";
+    }
+    for (std::size_t i = 0; i < nVertices; ++i) {
+      for (std::size_t j = 0; j < nVertices; ++j) {
+        if (adjacencyMatrix[i][j]) {
+          ss << "  " << i << " -> " << j << ";\n";
+        }
+      }
+    }
+    ss << "}\n";
+    return ss.str();
+  }
+  friend auto operator<<(std::ostream& os,
+                         const DirectedGraph& g) -> std::ostream& {
+    return os << g.toString(); // Using toString() method
+  }
+};
+} // namespace qc

include/mqt-core/datastructures/DisjointSet.hpp

@@ -0,0 +1,39 @@
+#pragma once
+
+#include <algorithm>
+#include <unordered_map>
+
+namespace qc {
+template <class T> struct DisjointSet {
+  std::unordered_map<T, T> parent;
+  std::unordered_map<T, std::size_t> rank;
+
+  template <class Iterator>
+  explicit DisjointSet(const Iterator& begin, const Iterator& end) {
+    std::for_each(begin, end, [&](const auto& element) {
+      parent[element] = element;
+      rank[element] = 0;
+    });
+  }
+
+  T findSet(const T& v) {
+    if (parent[v] != v) {
+      parent[v] = findSet(parent[v]);
+    }
+    return parent[v];
+  }
+
+  void unionSet(const T& x, const T& y) {
+    const auto& xSet = findSet(x);
+    const auto& ySet = findSet(y);
+    if (rank[xSet] > rank[ySet]) {
+      parent[ySet] = xSet;
+    } else {
+      parent[xSet] = ySet;
+      if (rank[xSet] == rank[ySet]) {
+        ++rank[ySet];
+      }
+    }
+  }
+};
+} // namespace qc