Add native N-control gate support to lightning.gpu

.github/CHANGELOG.md

@@ -2,6 +2,9 @@
 
 ### New features since last release
 
+* Add native N-controlled gates support to `lightning.gpu`'s single-GPU backend.
+  [(#938)](https://github.com/PennyLaneAI/pennylane-lightning/pull/938)
+
 * Add `mid-circuit measurements` support to `lightning.gpu`'s single-GPU backend.
   [(#931)](https://github.com/PennyLaneAI/pennylane-lightning/pull/931)
 

pennylane_lightning/core/_version.py

@@ -16,4 +16,4 @@
    Version number (major.minor.patch[-label])
 """
 
-__version__ = "0.39.0-dev51"
+__version__ = "0.39.0-dev52"

pennylane_lightning/core/src/simulators/lightning_gpu/bindings/LGPUBindings.hpp

@@ -51,9 +51,55 @@ using StateVectorBackends =
                               StateVectorCudaManaged<double>, void>;
 
 /**
- * @brief Get a gate kernel map for a statevector.
+ * @brief Register controlled matrix kernel.
  */
+template <class StateVectorT>
+void applyControlledMatrix(
+    StateVectorT &st,
+    const py::array_t<std::complex<typename StateVectorT::PrecisionT>,
+                      py::array::c_style | py::array::forcecast> &matrix,
+    const std::vector<std::size_t> &controlled_wires,
+    const std::vector<bool> &controlled_values,
+    const std::vector<std::size_t> &wires, bool inverse = false) {
+    using ComplexT = typename StateVectorT::ComplexT;
+    st.applyControlledMatrix(
+        static_cast<const ComplexT *>(matrix.request().ptr),
+        static_cast<const std::size_t>(matrix.request().size), controlled_wires,
+        controlled_values, wires, inverse);
+}
+
+template <class StateVectorT, class PyClass>
+void registerControlledGate(PyClass &pyclass) {
+    using PrecisionT =
+        typename StateVectorT::PrecisionT; // Statevector's precision
+    using ParamT = PrecisionT;             // Parameter's data precision
+
+    using Pennylane::Gates::ControlledGateOperation;
+    using Pennylane::Util::for_each_enum;
+    namespace Constant = Pennylane::Gates::Constant;
 
+    for_each_enum<ControlledGateOperation>(
+        [&pyclass](ControlledGateOperation gate_op) {
+            using Pennylane::Util::lookup;
+            const auto gate_name =
+                std::string(lookup(Constant::controlled_gate_names, gate_op));
+            const std::string doc = "Apply the " + gate_name + " gate.";
+            auto func = [gate_name](
+                            StateVectorT &sv,
+                            const std::vector<std::size_t> &controlled_wires,
+                            const std::vector<bool> &controlled_values,
+                            const std::vector<std::size_t> &wires, bool inverse,
+                            const std::vector<ParamT> &params) {
+                sv.applyOperation(gate_name, controlled_wires,
+                                  controlled_values, wires, inverse, params);
+            };
+            pyclass.def(gate_name.c_str(), func, doc.c_str());
+        });
+}
+
+/**
+ * @brief Get a gate kernel map for a statevector.
+ */
 template <class StateVectorT, class PyClass>
 void registerBackendClassSpecificBindings(PyClass &pyclass) {
     using PrecisionT =
@@ -65,6 +111,7 @@ void registerBackendClassSpecificBindings(PyClass &pyclass) {
                                  py::array::c_style | py::array::forcecast>;
 
     registerGatesForStateVector<StateVectorT>(pyclass);
+    registerControlledGate<StateVectorT>(pyclass);
 
     pyclass
         .def(py::init<std::size_t>())              // qubits, device
@@ -96,6 +143,8 @@ void registerBackendClassSpecificBindings(PyClass &pyclass) {
             },
             "Set State Vector on GPU with values for the state vector and "
             "wires on the host memory.")
+        .def("applyControlledMatrix", &applyControlledMatrix<StateVectorT>,
+             "Apply controlled operation")
         .def(
             "DeviceToDevice",
             [](StateVectorT &sv, const StateVectorT &other, bool async) {

pennylane_lightning/core/src/simulators/lightning_gpu/gates/tests/Test_StateVectorCudaManaged_NonParam.cpp

@@ -1105,4 +1105,257 @@ TEMPLATE_TEST_CASE("StateVectorCudaManaged::SetIthStates",
 
         CHECK(expected_state == Pennylane::Util::approx(sv.getDataVector()));
     }
-}
+}
+
+TEMPLATE_TEST_CASE("StateVectorCudaManaged::applyOperation non-param "
+                   "one-qubit with controls",
+                   "[StateVectorCudaManaged]", float, double) {
+    using PrecisionT = TestType;
+    std::mt19937 re{1337};
+    const int num_qubits = 4;
+    const auto margin = PrecisionT{1e-5};
+    const std::size_t control = GENERATE(0, 1, 2, 3);
+    const std::size_t wire = GENERATE(0, 1, 2, 3);
+    auto st0 = createRandomStateVectorData<PrecisionT>(re, num_qubits);
+    StateVectorCudaManaged<PrecisionT> sv0(num_qubits);
+    StateVectorCudaManaged<PrecisionT> sv1(num_qubits);
+
+    sv0.CopyHostDataToGpu(st0.data(), st0.size());
+    sv1.CopyHostDataToGpu(st0.data(), st0.size());
+
+    DYNAMIC_SECTION("N-controlled PauliX - "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire) {
+            sv0.applyOperation("CNOT", {control, wire});
+            sv1.applyOperation("PauliX", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        } else {
+            PL_REQUIRE_THROWS_MATCHES(
+                sv0.applyOperation("PauliX", std::vector<std::size_t>{control},
+                                   std::vector<bool>{true},
+                                   std::vector<std::size_t>{wire}),
+                LightningException,
+                "`controlled_wires` and target wires must be disjoint.");
+        }
+
+        if (control != 0 && wire != 0 && control != wire) {
+            sv0.applyOperation("Toffoli", {0, control, wire});
+            sv1.applyOperation("PauliX", std::vector<std::size_t>{0, control},
+                               std::vector<bool>{true, true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+
+            sv0.applyOperation("Toffoli", {control, 0, wire});
+            sv1.applyOperation("PauliX", std::vector<std::size_t>{control, 0},
+                               std::vector<bool>{true, true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+
+    DYNAMIC_SECTION("N-controlled PauliY - "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire) {
+            sv0.applyOperation("CY", {control, wire});
+            sv1.applyOperation("PauliY", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+
+    DYNAMIC_SECTION("N-controlled PauliZ - "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire) {
+            sv0.applyOperation("CZ", {control, wire});
+            sv1.applyOperation("PauliZ", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+
+    DYNAMIC_SECTION("N-controlled Hadamard - "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire) {
+            const auto matrix = getHadamard<std::complex, PrecisionT>();
+
+            sv0.applyControlledMatrix(
+                matrix.data(), matrix.size(), std::vector<std::size_t>{control},
+                std::vector<bool>{true}, std::vector<std::size_t>{wire});
+            sv1.applyOperation("Hadamard", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+    DYNAMIC_SECTION("N-controlled S - "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire) {
+            const auto matrix = getS<std::complex, PrecisionT>();
+
+            sv0.applyControlledMatrix(
+                matrix.data(), matrix.size(), std::vector<std::size_t>{control},
+                std::vector<bool>{true}, std::vector<std::size_t>{wire});
+            sv1.applyOperation("S", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+
+    DYNAMIC_SECTION("N-controlled T - "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire) {
+            const std::vector<std::complex<PrecisionT>> matrix =
+                getT<std::complex, PrecisionT>();
+
+            sv0.applyControlledMatrix(
+                matrix.data(), matrix.size(), std::vector<std::size_t>{control},
+                std::vector<bool>{true}, std::vector<std::size_t>{wire});
+            sv1.applyOperation("T", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+}
+
+TEMPLATE_TEST_CASE("StateVectorCudaManaged::applyOperation non-param "
+                   "two-qubit with controls",
+                   "[StateVectorCudaManaged]", float, double) {
+    using PrecisionT = TestType;
+    std::mt19937 re{1337};
+    const int num_qubits = 4;
+    const auto margin = PrecisionT{1e-5};
+    const std::size_t control = GENERATE(0, 1, 2, 3);
+    const std::size_t wire0 = GENERATE(0, 1, 2, 3);
+    const std::size_t wire1 = GENERATE(0, 1, 2, 3);
+    auto st0 = createRandomStateVectorData<PrecisionT>(re, num_qubits);
+    StateVectorCudaManaged<PrecisionT> sv0(num_qubits);
+    StateVectorCudaManaged<PrecisionT> sv1(num_qubits);
+
+    sv0.CopyHostDataToGpu(st0.data(), st0.size());
+    sv1.CopyHostDataToGpu(st0.data(), st0.size());
+
+    DYNAMIC_SECTION("N-controlled SWAP - "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire0 << ", " << wire1 << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire0 && control != wire1 && wire0 != wire1) {
+            sv0.applyOperation("CSWAP", {control, wire0, wire1});
+            sv1.applyOperation("SWAP", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire0, wire1});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+
+    DYNAMIC_SECTION("N-controlled SWAP with matrix- "
+                    << "controls = {" << control << "} "
+                    << ", wires = {" << wire0 << ", " << wire1 << "} - "
+                    << PrecisionToName<PrecisionT>::value) {
+        if (control != wire0 && control != wire1 && wire0 != wire1) {
+            const std::vector<std::complex<PrecisionT>> matrix =
+                getSWAP<std::complex, PrecisionT>();
+            sv0.applyControlledMatrix(matrix.data(), matrix.size(),
+                                      std::vector<std::size_t>{control},
+                                      std::vector<bool>{true},
+                                      std::vector<std::size_t>{wire0, wire1});
+            sv1.applyOperation("SWAP", std::vector<std::size_t>{control},
+                               std::vector<bool>{true},
+                               std::vector<std::size_t>{wire0, wire1});
+            REQUIRE(sv0.getDataVector() ==
+                    approx(sv1.getDataVector()).margin(margin));
+        }
+    }
+}
+
+TEMPLATE_TEST_CASE("StateVectorCudaManaged::controlled Toffoli",
+                   "[StateVectorCudaManaged]", float, double) {
+    using PrecisionT = TestType;
+    std::mt19937 re{1337};
+    const int num_qubits = 6;
+    const auto margin = PrecisionT{1e-5};
+    const std::size_t control = GENERATE(0, 1, 2);
+    auto st0 = createRandomStateVectorData<PrecisionT>(re, num_qubits);
+
+    StateVectorCudaManaged<PrecisionT> sv0(num_qubits);
+    StateVectorCudaManaged<PrecisionT> sv1(num_qubits);
+
+    sv0.CopyHostDataToGpu(st0.data(), st0.size());
+    sv1.CopyHostDataToGpu(st0.data(), st0.size());
+
+    const std::vector<std::complex<PrecisionT>> matrix =
+        getToffoli<std::complex, PrecisionT>();
+    sv0.applyControlledMatrix(
+        matrix.data(), matrix.size(), std::vector<std::size_t>{control},
+        std::vector<bool>{true}, std::vector<std::size_t>{3, 4, 5});
+    sv1.applyOperation("PauliX", std::vector<std::size_t>{control, 3, 4},
+                       std::vector<bool>{true, true, true},
+                       std::vector<std::size_t>{5});
+    REQUIRE(sv0.getDataVector() == approx(sv1.getDataVector()).margin(margin));
+}
+
+TEMPLATE_TEST_CASE(
+    "StateVectorCudaManaged::controlled_gate_direct_matrix_offload",
+    "[StateVectorCudaManaged]", float, double) {
+    using PrecisionT = TestType;
+    std::mt19937 re{1337};
+    const int num_qubits = 6;
+    const auto margin = PrecisionT{1e-5};
+    const std::size_t control = GENERATE(0, 1, 2);
+    auto st0 = createRandomStateVectorData<PrecisionT>(re, num_qubits);
+
+    StateVectorCudaManaged<PrecisionT> sv0(num_qubits);
+    StateVectorCudaManaged<PrecisionT> sv1(num_qubits);
+
+    sv0.CopyHostDataToGpu(st0.data(), st0.size());
+    sv1.CopyHostDataToGpu(st0.data(), st0.size());
+
+    SECTION("Catch failures caused by unsupported named gates") {
+        PL_CHECK_THROWS_MATCHES(
+            sv0.applyOperation("paulix",
+                               std::vector<std::size_t>{control, 3, 4},
+                               std::vector<bool>{true, true, true},
+                               std::vector<std::size_t>{5}),
+            LightningException, "Currently unsupported gate: paulix");
+    }
+
+    SECTION("direct base matrix offload") {
+        const std::vector<std::complex<PrecisionT>> matrix = {
+            {0.0, 0.0}, {1.0, 0.0}, {1.0, 0.0}, {0.0, 0.0}};
+        sv0.applyControlledMatrix(matrix.data(), matrix.size(),
+                                  std::vector<std::size_t>{control, 3, 4},
+                                  std::vector<bool>{true, true, true},
+                                  std::vector<std::size_t>{5});
+        sv1.applyOperation("Paulix", std::vector<std::size_t>{control, 3, 4},
+                           std::vector<bool>{true, true, true},
+                           std::vector<std::size_t>{5}, false, {0.0}, matrix);
+        REQUIRE(sv0.getDataVector() ==
+                approx(sv1.getDataVector()).margin(margin));
+    }
+}