🚸 Ancillary and garbage handling

include/checker/dd/DDAlternatingChecker.hpp

@@ -39,6 +39,12 @@ class DDAlternatingChecker final
     j["checker"] = "decision_diagram_alternating";
   }
 
+  /// a function to determine whether the alternating checker can handle
+  /// checking both circuits. In particular, it checks whether both circuits
+  /// contain non-idle ancillaries.
+  static bool canHandle(const qc::QuantumComputation& qc1,
+                        const qc::QuantumComputation& qc2);
+
 private:
   qc::MatrixDD functionality{};
 

src/EquivalenceCheckingManager.cpp

@@ -196,6 +196,16 @@ EquivalenceCheckingManager::EquivalenceCheckingManager(
     fixOutputPermutationMismatch();
   }
 
+  // check whether the alternating checker is configured and can handle the
+  // circuits
+  if (configuration.execution.runAlternatingChecker &&
+      !DDAlternatingChecker::canHandle(this->qc1, this->qc2)) {
+    std::clog << "[QCEC] Warning: alternating checker cannot handle the "
+                 "circuits. Falling back to construction checker.\n";
+    this->configuration.execution.runAlternatingChecker  = false;
+    this->configuration.execution.runConstructionChecker = true;
+  }
+
   // initialize the stimuli generator
   stateGenerator = StateGenerator(configuration.simulation.seed);
 

src/checker/dd/DDAlternatingChecker.cpp

@@ -12,37 +12,26 @@ void DDAlternatingChecker::initialize() {
   functionality = dd->makeIdent(nqubits);
   dd->incRef(functionality);
 
-  // only count ancillaries that are present in but not acted upon in both of
-  // the circuits at the moment this is just to be on the safe side. It might be
-  // fine to also start with the reduced matrix for every ancillary without any
-  // restriction
-  // TODO: check whether the way ancillaries are handled here is theoretically
-  // sound
+  // Only count ancillaries that are present in but not acted upon in both of
+  // the circuits. Otherwise, the alternating checker must not be used.
+  // Counter-example: H |0><0| H^-1 = [[0.5, 0.5], [0.5, 0.5]] != |0><0|
   std::vector<bool> ancillary(nqubits);
   for (auto q = static_cast<dd::Qubit>(nqubits - 1U); q >= 0; --q) {
     if (qc1.logicalQubitIsAncillary(q) && qc2.logicalQubitIsAncillary(q)) {
-      bool found1  = false;
-      bool isIdle1 = false;
-      if (const auto it = std::find_if(
-              qc1.initialLayout.cbegin(), qc1.initialLayout.cend(),
-              [q](const auto& mapping) { return mapping.second == q; });
-          it != qc1.initialLayout.cend()) {
-        found1  = true;
-        isIdle1 = qc1.isIdleQubit(it->first);
-      }
-      bool found2  = false;
-      bool isIdle2 = false;
-      if (const auto it = std::find_if(
-              qc2.initialLayout.cbegin(), qc2.initialLayout.cend(),
-              [q](const auto& mapping) { return mapping.second == q; });
-          it != qc2.initialLayout.cend()) {
-        found2  = true;
-        isIdle2 = qc2.isIdleQubit(it->first);
-      }
+      const auto [found1, physical1] = qc1.containsLogicalQubit(q);
+      const auto isIdle1             = found1 && qc1.isIdleQubit(*physical1);
+
+      const auto [found2, physical2] = qc2.containsLogicalQubit(q);
+      const auto isIdle2             = found2 && qc2.isIdleQubit(*physical2);
 
       // qubit only really exists or is acted on in one of the circuits
       if ((found1 != found2) || (isIdle1 != isIdle2)) {
         ancillary[static_cast<std::size_t>(q)] = true;
+      } else {
+        throw std::invalid_argument(
+            "Alternating checker must not be used for "
+            "circuits that both have non-idle ancillary "
+            "qubits. Use the construction checker instead.");
       }
     }
   }
@@ -51,8 +40,8 @@ void DDAlternatingChecker::initialize() {
   // [1 0] if the qubit is no ancillary, or it is acted upon by both circuits
   // [0 1]
   //
-  // [1 0] for an ancillary that is present in one circuit and not acted upon in
-  // the other [0 0]
+  // [1 0] (= |0><0|) for an ancillary only acted on in one circuit
+  // [0 0]
   functionality = dd->reduceAncillae(functionality, ancillary);
 }
 
@@ -122,16 +111,6 @@ void DDAlternatingChecker::postprocess() {
   if (isDone()) {
     return;
   }
-
-  // TODO: check whether reducing ancillaries here is theoretically sound
-  taskManager1.reduceAncillae(functionality);
-  if (isDone()) {
-    return;
-  }
-  taskManager2.reduceAncillae(functionality);
-  if (isDone()) {
-    return;
-  }
 }
 
 EquivalenceCriterion DDAlternatingChecker::checkEquivalence() {
@@ -143,15 +122,14 @@ EquivalenceCriterion DDAlternatingChecker::checkEquivalence() {
   taskManager1.reduceGarbage(goalMatrix);
   taskManager2.reduceGarbage(goalMatrix);
 
-  // TODO: check whether reducing ancillaries here is theoretically sound
   taskManager1.reduceAncillae(goalMatrix);
   taskManager2.reduceAncillae(goalMatrix);
 
   // the resulting goal matrix is
   // [1 0] if the qubit is no ancillary
   // [0 1]
   //
-  // [1 0] for an ancillary that is present in either circuit
+  // [1 0] (= |0><0>|) for an ancillary that is present in either circuit
   // [0 0]
 
   // compare the obtained functionality to the goal matrix
@@ -169,4 +147,32 @@ EquivalenceCriterion DDAlternatingChecker::checkEquivalence() {
   return op1.equals(op2);
 }
 
+bool DDAlternatingChecker::canHandle(const qc::QuantumComputation& qc1,
+                                     const qc::QuantumComputation& qc2) {
+  assert(qc1.getNqubits() == qc2.getNqubits());
+  const auto nqubits = qc1.getNqubits();
+
+  for (auto q = static_cast<dd::Qubit>(nqubits - 1U); q >= 0; --q) {
+    if (qc1.logicalQubitIsAncillary(q) && qc2.logicalQubitIsAncillary(q)) {
+      const auto [found1, physical1] = qc1.containsLogicalQubit(q);
+      const auto [found2, physical2] = qc2.containsLogicalQubit(q);
+
+      // just continue, if a qubit is not found in both circuits
+      if (found1 != found2) {
+        continue;
+      }
+
+      const auto isIdle1 = found1 && qc1.isIdleQubit(*physical1);
+      const auto isIdle2 = found2 && qc2.isIdleQubit(*physical2);
+
+      // if an ancillary qubit is acted on in both circuits,
+      // the alternating checker cannot be used.
+      if (!isIdle1 && !isIdle2) {
+        return false;
+      }
+    }
+  }
+  return true;
+}
+
 } // namespace ec

src/checker/dd/DDEquivalenceChecker.cpp

@@ -52,14 +52,19 @@ DDEquivalenceChecker<DDType, DDPackage>::equals(const DDType& e,
     // product trace(U V^-1) and comparing it to some threshold. in a similar
     // fashion, we can simply compare U V^-1 with the identity, which results in
     // a much simpler check that is not prone to overflow.
-    bool isClose{};
-    if (e.p->isIdentity()) {
-      isClose =
-          dd->isCloseToIdentity(f, configuration.functionality.traceThreshold);
-    } else if (f.p->isIdentity()) {
-      isClose =
-          dd->isCloseToIdentity(e, configuration.functionality.traceThreshold);
+    bool       isClose{};
+    const bool eIsClose =
+        dd->isCloseToIdentity(e, configuration.functionality.traceThreshold);
+    const bool fIsClose =
+        dd->isCloseToIdentity(f, configuration.functionality.traceThreshold);
+    if (eIsClose || fIsClose) {
+      // if any of the DDs is close to the identity (structure), the result can
+      // be decided by whether both DDs are close enough to the identity.
+      isClose = eIsClose && fIsClose;
     } else {
+      // otherwise, one DD needs to be inverted before multiplying both of them
+      // together and checking whether the resulting DD is close enough to the
+      // identity.
       auto g = dd->multiply(e, dd->conjugateTranspose(f));
       isClose =
           dd->isCloseToIdentity(g, configuration.functionality.traceThreshold);

test/test_equality.cpp

@@ -125,3 +125,43 @@ TEST_F(EqualityTest, SimulationMoreThan64Qubits) {
   std::cout << ecm << std::endl;
   EXPECT_EQ(ecm.equivalence(), ec::EquivalenceCriterion::ProbablyEquivalent);
 }
+
+TEST_F(EqualityTest, AutomaticSwitchToConstructionChecker) {
+  // add ancillary qubits to both circuits
+  qc1.addAncillaryQubit(1, -1);
+  qc2.addAncillaryQubit(1, -1);
+
+  // perform the same action on both circuits' primary qubit
+  qc1.x(0);
+  qc2.x(0);
+
+  // perform a different action on the ancillary qubit
+  qc1.x(1);
+  qc2.z(1);
+
+  // setup default configuration
+  config = ec::Configuration{};
+  ec::EquivalenceCheckingManager ecm(qc1, qc2, config);
+
+  // this should notice that the alternating checker is not capable of running
+  // the circuit and should switch to the construction checker
+  const auto runConfig = ecm.getConfiguration();
+  EXPECT_TRUE(runConfig.execution.runConstructionChecker);
+  EXPECT_FALSE(runConfig.execution.runAlternatingChecker);
+
+  // run the equivalence checker
+  ecm.run();
+
+  // both circuits should be equivalent since their action only differs on an
+  // ancillary and garbage qubit
+  const auto result = ecm.equivalence();
+  EXPECT_EQ(result, ec::EquivalenceCriterion::Equivalent);
+
+  // Check an exception is raised for a checker configured after initialization.
+  // Note: this exception can only be caught in sequential mode since it is
+  // raised in a different thread otherwise.
+  ecm.reset();
+  ecm.setAlternatingChecker(true);
+  ecm.setParallel(false);
+  EXPECT_THROW(ecm.run(), std::invalid_argument);
+}