Add support for ControlledQubitUnitary

pennylane_qulacs/qulacs_device.py

@@ -21,7 +21,7 @@
 import numpy as np
 
 from pennylane import QubitDevice, DeviceError
-from pennylane.ops import QubitStateVector, BasisState, QubitUnitary, CRZ, PhaseShift, Adjoint
+from pennylane.ops import QubitStateVector, BasisState, QubitUnitary, ControlledQubitUnitary, CRZ, PhaseShift, Adjoint
 
 import qulacs.gate as gate
 from qulacs import QuantumCircuit, QuantumState, Observable
@@ -88,6 +88,7 @@ class QulacsDevice(QubitDevice):
         "QubitStateVector": None,
         "BasisState": None,
         "QubitUnitary": None,
+        "ControlledQubitUnitary": None,
         "Toffoli": gate.TOFFOLI,
         "CSWAP": gate.FREDKIN,
         "CRZ": crz,
@@ -175,6 +176,11 @@ def apply_operations(self, operations):
                 self._apply_basis_state(op)
             elif isinstance(op, QubitUnitary):
                 self._apply_qubit_unitary(op, inverse)
+            elif isinstance(op, ControlledQubitUnitary):
+                if len(op.control_wires) == 0:
+                    self._apply_qubit_unitary(op, inverse)
+                else:
+                    self._apply_controlled_qubit_unitary(op, inverse)
             elif isinstance(op, (CRZ, PhaseShift)):
                 self._apply_matrix(op, inverse)
             else:
@@ -254,6 +260,34 @@ def _apply_qubit_unitary(self, op, inverse=False):
         self._circuit.add_gate(unitary_gate)
         unitary_gate.update_quantum_state(self._state)
 
+    def _apply_controlled_qubit_unitary(self, op, inverse=False):
+        """Apply controlled unitary to state"""
+        # translate op wire labels to consecutive wire labels used by the device
+        device_wires = self.map_wires(op.wires)
+        target_wires = self.map_wires(op.target_wires)
+        control_wires = self.map_wires(op.control_wires)
+        control_values = op.control_values 
+        par = op.parameters
+
+        if len(control_wires) + len(target_wires) != len(device_wires):
+            raise ValueError("len(device_wire) should be equal to len(control_wires) + len(target_wires)")
+
+        if len(par[0]) != 2 ** len(target_wires):
+            raise ValueError("Unitary matrix must be of shape (2**target_wires, 2**target_wires).")
+
+        if inverse:
+            par[0] = par[0].conj().T
+
+        # reverse wires (could also change par[0])
+        reverse_target_wire_labels = target_wires.tolist()[::-1]
+        reverse_control_wire_labels = control_wires.tolist()[::-1]
+        reverse_control_value = control_values[::-1]
+        unitary_gate = gate.DenseMatrix(reverse_target_wire_labels, par[0])
+        for i,j in enumerate(reverse_control_wire_labels):
+            unitary_gate.add_control_qubit(j, reverse_control_value[i])
+        self._circuit.add_gate(unitary_gate)
+        unitary_gate.update_quantum_state(self._state)
+
     def _apply_matrix(self, op, inverse=False):
         """Apply predefined gate-matrix to state (must follow qulacs convention)"""
         # translate op wire labels to consecutive wire labels used by the device

tests/test_apply.py

@@ -249,6 +249,33 @@ def test_qubit_unitary(self, init_state, mat, tol):
         expected = mat @ state
         assert np.allclose(res, expected, tol)
 
+    @pytest.mark.parametrize("control_wires", [[0], [0, 1]])
+    @pytest.mark.parametrize("mat", [U, U2])
+    def test_controlled_qubit_unitary(self, init_state, mat, control_wires, tol):
+        """Test ControlledQubitUnitary application"""
+
+        N = int(np.log2(len(mat)))
+        dev = QulacsDevice(N + len(control_wires))
+        state = init_state(N + len(control_wires))
+
+        op = qml.ControlledQubitUnitary(mat, wires=list(range(len(control_wires), N+len(control_wires))), 
+                                        control_wires=control_wires)
+        dev.apply([qml.QubitStateVector(state, wires=list(range(N + len(control_wires)))), op])
+        dev._obs_queue = []
+
+        res = dev.state
+        zero_op = np.array([[1, 0], [0, 0]])
+        one_op = np.array([[0, 0], [0, 1]])
+        if control_wires == [0]:
+            expected = (np.kron(zero_op, np.eye(2 ** N)) + np.kron(one_op, mat)) @ state
+        elif control_wires == [0, 1]:
+            expected = (np.kron(np.kron(zero_op, zero_op), np.eye(2 ** N)) + 
+                        np.kron(np.kron(zero_op, one_op), np.eye(2 ** N)) + 
+                        np.kron(np.kron(one_op, zero_op), np.eye(2 ** N)) + 
+                        np.kron(np.kron(one_op, one_op), mat)) @ state 
+
+        assert np.allclose(res, expected, tol)
+
     def test_invalid_qubit_state_unitary(self):
         """Test that an exception is raised if the
         unitary matrix is the wrong size"""