Check qubit release/measurement status on release

src/Qir/Runtime/lib/Simulators/FullstateSimulator.cpp

@@ -24,6 +24,7 @@
 
 #include "FloatUtils.hpp"
 #include "QirTypes.hpp" // TODO: Consider removing dependency on this file.
+#include "QirRuntime.hpp"
 #include "QirRuntimeApi_I.hpp"
 #include "QSharpSimApi_I.hpp"
 #include "SimFactory.hpp"
@@ -164,6 +165,19 @@ namespace Quantum
             return proc;
         }
 
+        void UnmarkAsMeasuredSingleQubit(Qubit q)
+        {
+            isMeasured[GetQubitId(q)] = false;
+        }
+
+        void UnmarkAsMeasuredQubitList(long num, Qubit* qubit)
+        {
+            for (const auto& id : GetQubitIds(num, qubit))
+            {
+                isMeasured[id] = false;
+            }
+        }
+
       public:
         CFullstateSimulator(uint32_t userProvidedSeed = 0) : handle(LoadQuantumSimulator())
         {
@@ -220,16 +234,28 @@ namespace Quantum
             Qubit q     = qubitManager->Allocate(); // Allocate qubit in qubit manager.
             unsigned id = GetQubitId(q);            // Get its id.
             allocateQubit(this->simulatorId, id);   // Allocate it in the simulator.
+            if (isMeasured.size() < id + 1)
+            {
+                isMeasured.resize(id + 1, false);
+            }
             return q;
         }
 
         void ReleaseQubit(Qubit q) override
         {
-            typedef void (*TReleaseQubit)(unsigned, unsigned);
+            typedef bool (*TReleaseQubit)(unsigned, unsigned);
             static TReleaseQubit releaseQubit = reinterpret_cast<TReleaseQubit>(this->GetProc("release"));
 
-            releaseQubit(this->simulatorId, GetQubitId(q)); // Release qubit in the simulator.
-            qubitManager->Release(q);                       // Release it in the qubit manager.
+            // Release qubit in the simulator, checking to make sure that release was valid.
+            auto id = GetQubitId(q);
+            if (!releaseQubit(this->simulatorId, id) && !isMeasured[id])
+            {
+                // We reject the release of a qubit that is not in the ground state (releaseQubit returns false),
+                // and was not recently measured (ie: the last operation was not measurement). This means the
+                // state is not well known, and therefore the safety of release is not guaranteed.
+                quantum__rt__fail_cstr("Released qubit neither measured nor in ground state.");
+            }
+            qubitManager->Release(q); // Release it in the qubit manager.
         }
 
         Result Measure(long numBases, PauliId bases[], long numTargets, Qubit targets[]) override
@@ -238,6 +264,11 @@ namespace Quantum
             typedef unsigned (*TMeasure)(unsigned, unsigned, unsigned*, unsigned*);
             static TMeasure m         = reinterpret_cast<TMeasure>(this->GetProc("Measure"));
             std::vector<unsigned> ids = GetQubitIds(numTargets, targets);
+            if (ids.size() == 1)
+            {
+                // If measuring exactly one qubit, mark it as measured for tracking.
+                isMeasured[ids[0]] = true;
+            }
             return reinterpret_cast<Result>(
                 m(this->simulatorId, (unsigned)numBases, reinterpret_cast<unsigned*>(bases), ids.data()));
         }
@@ -272,71 +303,87 @@ namespace Quantum
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("X"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledX(long numControls, Qubit controls[], Qubit target) override
         {
             static TSingleQubitControlledGate op = reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCX"));
             std::vector<unsigned> ids            = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void Y(Qubit q) override
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("Y"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledY(long numControls, Qubit controls[], Qubit target) override
         {
             static TSingleQubitControlledGate op = reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCY"));
             std::vector<unsigned> ids            = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void Z(Qubit q) override
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("Z"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledZ(long numControls, Qubit controls[], Qubit target) override
         {
             static TSingleQubitControlledGate op = reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCZ"));
             std::vector<unsigned> ids            = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void H(Qubit q) override
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("H"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledH(long numControls, Qubit controls[], Qubit target) override
         {
             static TSingleQubitControlledGate op = reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCH"));
             std::vector<unsigned> ids            = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void S(Qubit q) override
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("S"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledS(long numControls, Qubit controls[], Qubit target) override
         {
             static TSingleQubitControlledGate op = reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCS"));
             std::vector<unsigned> ids            = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void AdjointS(Qubit q) override
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("AdjS"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledAdjointS(long numControls, Qubit controls[], Qubit target) override
@@ -345,25 +392,31 @@ namespace Quantum
                 reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCAdjS"));
             std::vector<unsigned> ids = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void T(Qubit q) override
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("T"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledT(long numControls, Qubit controls[], Qubit target) override
         {
             static TSingleQubitControlledGate op = reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCT"));
             std::vector<unsigned> ids            = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void AdjointT(Qubit q) override
         {
             static TSingleQubitGate op = reinterpret_cast<TSingleQubitGate>(this->GetProc("AdjT"));
             op(this->simulatorId, GetQubitId(q));
+            UnmarkAsMeasuredSingleQubit(q);
         }
 
         void ControlledAdjointT(long numControls, Qubit controls[], Qubit target) override
@@ -372,6 +425,8 @@ namespace Quantum
                 reinterpret_cast<TSingleQubitControlledGate>(this->GetProc("MCAdjT"));
             std::vector<unsigned> ids = GetQubitIds(numControls, controls);
             op(this->simulatorId, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void R(PauliId axis, Qubit target, double theta) override
@@ -380,6 +435,7 @@ namespace Quantum
             static TR r = reinterpret_cast<TR>(this->GetProc("R"));
 
             r(this->simulatorId, GetBasis(axis), theta, GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
         }
 
         void ControlledR(long numControls, Qubit controls[], PauliId axis, Qubit target, double theta) override
@@ -389,6 +445,8 @@ namespace Quantum
 
             std::vector<unsigned> ids = GetQubitIds(numControls, controls);
             cr(this->simulatorId, GetBasis(axis), theta, (unsigned)numControls, ids.data(), GetQubitId(target));
+            UnmarkAsMeasuredSingleQubit(target);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         void Exp(long numTargets, PauliId paulis[], Qubit targets[], double theta) override
@@ -397,6 +455,7 @@ namespace Quantum
             static TExp exp           = reinterpret_cast<TExp>(this->GetProc("Exp"));
             std::vector<unsigned> ids = GetQubitIds(numTargets, targets);
             exp(this->simulatorId, (unsigned)numTargets, reinterpret_cast<unsigned*>(paulis), theta, ids.data());
+            UnmarkAsMeasuredQubitList(numTargets, targets);
         }
 
         void ControlledExp(long numControls, Qubit controls[], long numTargets, PauliId paulis[], Qubit targets[],
@@ -408,6 +467,8 @@ namespace Quantum
             std::vector<unsigned> idsControls = GetQubitIds(numControls, controls);
             cexp(this->simulatorId, (unsigned)numTargets, reinterpret_cast<unsigned*>(paulis), theta,
                  (unsigned)numControls, idsControls.data(), idsTargets.data());
+            UnmarkAsMeasuredQubitList(numTargets, targets);
+            UnmarkAsMeasuredQubitList(numControls, controls);
         }
 
         bool Assert(long numTargets, PauliId* bases, Qubit* targets, Result result, const char* failureMessage) override
@@ -436,6 +497,11 @@ namespace Quantum
         void GetStateTo(TDumpLocation location, TDumpToLocationCallback callback);
         bool GetRegisterTo(TDumpLocation location, TDumpToLocationCallback callback, const QirArray* qubits);
 
+        // This bit std::vector tracks whether the last operation on a given qubit was Measure.
+        // Note that `std::vector<bool>` is already specialized to use an underlying bitfied to save space.
+        // See: https://www.cplusplus.com/reference/vector/vector-bool/
+        std::vector<bool> isMeasured;
+
       private:
         TDumpToLocationCallback const dumpToLocationCallback = [](size_t idx, double re, double im,
                                                                   TDumpLocation location) -> bool {

src/Qir/Tests/FullstateSimulator/FullstateSimulatorTests.cpp

@@ -145,6 +145,10 @@ TEST_CASE("Fullstate simulator: ZZ measure", "[fullstate_simulator]")
     Result rOne = iqa->Measure(2, paulis, 2, q);
     REQUIRE(Result_One == sim->GetResultValue(rOne));
 
+    iqa->X(q[1]);
+    iqa->ControlledX(1, &q[0], q[1]);
+    iqa->H(q[0]);
+
     sim->ReleaseQubit(q[0]);
     sim->ReleaseQubit(q[1]);
 }
@@ -173,6 +177,8 @@ TEST_CASE("Fullstate simulator: assert probability", "[fullstate_simulator]")
     REQUIRE(!idig->Assert(2, xi, qs, sim->UseZero(), ""));
     REQUIRE(!idig->Assert(2, xi, qs, sim->UseOne(), ""));
 
+    iqa->X(qs[0]);
+
     sim->ReleaseQubit(qs[0]);
     sim->ReleaseQubit(qs[1]);
 }
@@ -196,6 +202,9 @@ TEST_CASE("Fullstate simulator: toffoli", "[fullstate_simulator]")
     iqa->ControlledX(2, qs, qs[2]);
     REQUIRE(Result_One == sim->GetResultValue(MZ(iqa, qs[2])));
 
+    iqa->X(qs[1]);
+    iqa->X(qs[0]);
+
     for (int i = 0; i < 3; i++)
     {
         sim->ReleaseQubit(qs[i]);
@@ -307,6 +316,8 @@ TEST_CASE("Fullstate simulator: exponents", "[fullstate_simulator]")
     PauliId paulis[3] = {PauliId_X, PauliId_Y, PauliId_Z};
     iqa->Exp(2, paulis, qs, 0.42);
     iqa->ControlledExp(2, qs, 3, paulis, &qs[2], 0.17);
+    iqa->ControlledExp(2, qs, 3, paulis, &qs[2], -0.17);
+    iqa->Exp(2, paulis, qs, -0.42);
 
     // not crashes? consider it passing
     REQUIRE(true);
@@ -377,12 +388,32 @@ TEST_CASE("Fullstate simulator: get qubit state of Bell state", "[fullstate_simu
     REQUIRE(1.0 == Approx(norm).epsilon(0.0001));
     norm = 0.0;
 
+    iqa->Y(qs[2]);
+    iqa->ControlledX(1, &qs[0], qs[1]);
+    iqa->H(qs[0]);
+
     for (int i = 0; i < n; i++)
     {
         sim->ReleaseQubit(qs[i]);
     }
 }
 
+extern "C" int Microsoft__Quantum__Testing__QIR__InvalidRelease__Interop(); // NOLINT
+TEST_CASE("QIR: Simulator rejects unmeasured, non-zero release", "[fullstate_simulator]")
+{
+    std::unique_ptr<IRuntimeDriver> sim = CreateFullstateSimulator();
+    QirExecutionContext::Scoped qirctx(sim.get(), true /*trackAllocatedObjects*/);
+    REQUIRE_THROWS(Microsoft__Quantum__Testing__QIR__InvalidRelease__Interop());
+}
+
+extern "C" int Microsoft__Quantum__Testing__QIR__MeasureRelease__Interop(); // NOLINT
+TEST_CASE("QIR: Simulator accepts measured release", "[fullstate_simulator]")
+{
+    std::unique_ptr<IRuntimeDriver> sim = CreateFullstateSimulator();
+    QirExecutionContext::Scoped qirctx(sim.get(), true /*trackAllocatedObjects*/);
+    REQUIRE_NOTHROW(Microsoft__Quantum__Testing__QIR__MeasureRelease__Interop());
+}
+
 extern "C" int Microsoft__Quantum__Testing__QIR__Test_Simulator_QIS__Interop(); // NOLINT
 TEST_CASE("QIR: invoke all standard Q# gates against the fullstate simulator", "[fullstate_simulator]")
 {

src/Qir/Tests/FullstateSimulator/qsharp/qir-test-simulator.qs

@@ -71,6 +71,7 @@ namespace Microsoft.Quantum.Testing.QIR
             H(ctls[0]);
             H(ctls[1]);
             if (M(targets[0]) != Zero) { set res = 2; }
+            ResetAll(targets + ctls);
         }
         if (res != 0) { return 70 + res; }
 
@@ -79,7 +80,23 @@ namespace Microsoft.Quantum.Testing.QIR
             H(qs[0]);
             H(qs[2]);
             if (Measure([PauliX, PauliZ, PauliX], qs) != Zero) { set res = 80; }
+            ResetAll(qs);
         }
         return res;
     }
+
+    @EntryPoint()
+    operation InvalidRelease() : Unit {
+        use q = Qubit();
+        let _ = M(q);
+        X(q);
+    }
+
+    @EntryPoint()
+    operation MeasureRelease() : Unit {
+        use qs = Qubit[2];
+        X(qs[0]);
+        let _ = Measure([PauliX], [qs[1]]);
+        let _ = M(qs[0]);
+    }
 }