Add rotation gates to stabilizer (#1938)

* Add rotation gates to stabilizer

* format tests

* set method=statevctor to sample tests

* set method=statevector to Estimator test, fix ry gate

* format test

Author: doichanj

qiskit_aer/backends/backend_utils.py

@@ -211,6 +211,9 @@
             "delay",
             "pauli",
             "ecr",
+            "rx",
+            "ry",
+            "rz",
         ]
     ),
     "extended_stabilizer": sorted(

releasenotes/notes/stabilizer_rotation-8ce2effd9578ee0a.yaml

@@ -0,0 +1,9 @@
+---
+upgrade:
+  - |
+    Adding support of rotation gates (rx, ry and rz gates) to stabilizer method
+    when input theta is multiple of pi/2.
+    If ``method=automatic`` is specified (this is default), if all the input
+    theta of rotation gates are multiple of pi/2 ``method=stabilizer``
+    is selected. Of when user sets ``method=stabilizer`` and any of theta
+    is not multiple of pi/2, Aer raises an exception.

src/simulators/circuit_executor.hpp

@@ -933,13 +933,20 @@ bool Executor<state_t>::validate_state(const Config &config,
 
   JSON::get_value(circ_name, "name", circ.header);
 
+  state.set_config(config);
   // Check if a circuit is valid for state ops
   bool circ_valid = state.opset().contains(circ.opset());
   if (throw_except && !circ_valid) {
     error_msg << "Circuit " << circ_name << " contains invalid instructions ";
     error_msg << state.opset().difference(circ.opset());
     error_msg << " for \"" << state.name() << "\" method.";
   }
+  // check parameters set inf ops
+  circ_valid &= state.validate_parameters(circ.ops);
+  if (throw_except && !circ_valid) {
+    error_msg << "Circuit " << circ_name << " contains invalid parameters ";
+    error_msg << " for \"" << state.name() << "\" method.";
+  }
 
   // Check if a noise model valid for state ops
   bool noise_valid = noise.is_ideal() || state.opset().contains(noise.opset());
@@ -952,8 +959,8 @@ bool Executor<state_t>::validate_state(const Config &config,
   // Validate memory requirements
   bool memory_valid = true;
   if (max_memory_mb_ > 0) {
-    size_t required_mb =
-        required_memory_mb(config, circ, noise) / num_process_per_experiment_;
+    size_t required_mb = state.required_memory_mb(circ.num_qubits, circ.ops) /
+                         num_process_per_experiment_;
     size_t mem_size = (sim_device_ == Device::GPU)
                           ? max_memory_mb_ + max_gpu_memory_mb_
                           : max_memory_mb_;

src/simulators/stabilizer/stabilizer_state.hpp

@@ -38,8 +38,8 @@ const Operations::OpSet StateOpSet(
      OpType::save_amps_sq, OpType::save_stabilizer, OpType::save_clifford,
      OpType::save_state, OpType::set_stabilizer, OpType::jump, OpType::mark},
     // Gates
-    {"CX", "cx", "cy", "cz", "swap", "id", "x", "y", "z", "h", "s", "sdg", "sx",
-     "sxdg", "delay", "pauli", "ecr"});
+    {"CX", "cx",  "cy", "cz",   "swap",  "id",    "x",   "y",  "z",  "h",
+     "s",  "sdg", "sx", "sxdg", "delay", "pauli", "ecr", "rx", "ry", "rz"});
 
 enum class Gates {
   id,
@@ -56,7 +56,10 @@ enum class Gates {
   cz,
   swap,
   pauli,
-  ecr
+  ecr,
+  rx,
+  ry,
+  rz
 };
 
 //============================================================================
@@ -101,6 +104,9 @@ class State : public QuantumState::State<Clifford::Clifford> {
   virtual std::vector<reg_t> sample_measure(const reg_t &qubits, uint_t shots,
                                             RngEngine &rng) override;
 
+  bool
+  validate_parameters(const std::vector<Operations::Op> &ops) const override;
+
 protected:
   //-----------------------------------------------------------------------
   // Apply instructions
@@ -203,7 +209,10 @@ const stringmap_t<Gates> State::gateset_({
     {"cz", Gates::cz},       // Controlled-Z gate
     {"swap", Gates::swap},   // SWAP gate
     {"pauli", Gates::pauli}, // Pauli gate
-    {"ecr", Gates::ecr}      // ECR gate
+    {"ecr", Gates::ecr},     // ECR gate
+    {"rx", Gates::rx},       // RX gate (only support k * pi/2 cases)
+    {"ry", Gates::ry},       // RY gate (only support k * pi/2 cases)
+    {"rz", Gates::rz}        // RZ gate (only support k * pi/2 cases)
 });
 
 //============================================================================
@@ -245,6 +254,23 @@ void State::set_config(const Config &config) {
   max_qubits_snapshot_probs_ = std::max<uint_t>(max_qubits_snapshot_probs_, 64);
 }
 
+bool State::validate_parameters(const std::vector<Operations::Op> &ops) const {
+  for (int_t i = 0; i < ops.size(); i++) {
+    if (ops[i].type == OpType::gate) {
+      // check parameter of R gates
+      if (ops[i].name == "rx" || ops[i].name == "ry" || ops[i].name == "rz") {
+        double pi2 = std::real(ops[i].params[0]) * 2.0 / M_PI;
+        double pi2_int = (double)std::round(pi2);
+
+        if (!AER::Linalg::almost_equal(pi2, pi2_int)) {
+          return false;
+        }
+      }
+    }
+  }
+  return true;
+}
+
 //=========================================================================
 // Implementation: apply operations
 //=========================================================================
@@ -298,6 +324,7 @@ void State::apply_op(const Operations::Op &op, ExperimentResult &result,
 }
 
 void State::apply_gate(const Operations::Op &op) {
+  int_t pi2;
   // Check Op is supported by State
   auto it = gateset_.find(op.name);
   if (it == gateset_.end())
@@ -369,6 +396,53 @@ void State::apply_gate(const Operations::Op &op) {
     BaseState::qreg_.append_x(op.qubits[0]);
     BaseState::qreg_.append_x(op.qubits[1]);
     break;
+  case Gates::rx:
+    pi2 = (int_t)std::round(std::real(op.params[0]) * 2.0 / M_PI) & 3;
+    if (pi2 == 1) {
+      // HSH
+      BaseState::qreg_.append_h(op.qubits[0]);
+      BaseState::qreg_.append_s(op.qubits[0]);
+      BaseState::qreg_.append_h(op.qubits[0]);
+    } else if (pi2 == 2) {
+      // X
+      BaseState::qreg_.append_x(op.qubits[0]);
+    } else if (pi2 == 3) {
+      // HSdgH
+      BaseState::qreg_.append_h(op.qubits[0]);
+      BaseState::qreg_.append_z(op.qubits[0]);
+      BaseState::qreg_.append_s(op.qubits[0]);
+      BaseState::qreg_.append_h(op.qubits[0]);
+    }
+    break;
+  case Gates::ry:
+    pi2 = (int_t)std::round(std::real(op.params[0]) * 2.0 / M_PI) & 3;
+    if (pi2 == 1) {
+      // HX
+      BaseState::qreg_.append_x(op.qubits[0]);
+      BaseState::qreg_.append_h(op.qubits[0]);
+    } else if (pi2 == 2) {
+      // Y
+      BaseState::qreg_.append_y(op.qubits[0]);
+    } else if (pi2 == 3) {
+      // Hdg
+      BaseState::qreg_.append_h(op.qubits[0]);
+      BaseState::qreg_.append_x(op.qubits[0]);
+    }
+    break;
+  case Gates::rz:
+    pi2 = (int_t)std::round(std::real(op.params[0]) * 2.0 / M_PI) & 3;
+    if (pi2 == 1) {
+      // S
+      BaseState::qreg_.append_s(op.qubits[0]);
+    } else if (pi2 == 2) {
+      // Z
+      BaseState::qreg_.append_z(op.qubits[0]);
+    } else if (pi2 == 3) {
+      // Sdg
+      BaseState::qreg_.append_z(op.qubits[0]);
+      BaseState::qreg_.append_s(op.qubits[0]);
+    }
+    break;
   default:
     // We shouldn't reach here unless there is a bug in gateset
     throw std::invalid_argument(

src/simulators/state.hpp

@@ -115,6 +115,12 @@ class Base {
   // Typically this is the n-qubit all |0> state
   virtual void initialize_qreg(uint_t num_qubits) = 0;
 
+  // validate parameters in input operations
+  virtual bool
+  validate_parameters(const std::vector<Operations::Op> &ops) const {
+    return true;
+  }
+
   //-----------------------------------------------------------------------
   // ClassicalRegister methods
   //-----------------------------------------------------------------------

test/terra/backends/aer_simulator/test_rotation.py

@@ -0,0 +1,76 @@
+# This code is part of Qiskit.
+#
+# (C) Copyright IBM 2018, 2019.
+#
+# This code is licensed under the Apache License, Version 2.0. You may
+# obtain a copy of this license in the LICENSE.txt file in the root directory
+# of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.
+#
+# Any modifications or derivative works of this code must retain this
+# copyright notice, and modified files need to carry a notice indicating
+# that they have been altered from the originals.
+"""
+AerSimulator Integration Tests
+"""
+from ddt import ddt
+from test.terra.reference import ref_rotation
+from qiskit import transpile
+from test.terra.backends.simulator_test_case import SimulatorTestCase, supported_methods
+
+SUPPORTED_METHODS = [
+    "automatic",
+    "stabilizer",
+    "statevector",
+    "density_matrix",
+    "matrix_product_state",
+    "tensor_network",
+]
+
+
+@ddt
+class TestRotation(SimulatorTestCase):
+    """AerSimulator Rotation gate tests"""
+
+    SEED = 12345
+
+    # ---------------------------------------------------------------------
+    # Test rx-gate
+    # ---------------------------------------------------------------------
+    @supported_methods(SUPPORTED_METHODS)
+    def test_rx_gate_deterministic(self, method, device):
+        """Test rx-gate circuits"""
+        backend = self.backend(method=method, device=device, seed_simulator=self.SEED)
+        shots = 1000
+        circuits = ref_rotation.rx_gate_circuits_deterministic(final_measure=True)
+        targets = ref_rotation.rx_gate_counts_deterministic(shots)
+        result = backend.run(circuits, shots=shots).result()
+        self.assertSuccess(result)
+        self.compare_counts(result, circuits, targets, delta=0.05 * shots)
+
+    # ---------------------------------------------------------------------
+    # Test rz-gate
+    # ---------------------------------------------------------------------
+    @supported_methods(SUPPORTED_METHODS)
+    def test_rz_gate_deterministic(self, method, device):
+        """Test rz-gate circuits"""
+        backend = self.backend(method=method, device=device, seed_simulator=self.SEED)
+        shots = 1000
+        circuits = ref_rotation.rz_gate_circuits_deterministic(final_measure=True)
+        targets = ref_rotation.rz_gate_counts_deterministic(shots)
+        result = backend.run(circuits, shots=shots).result()
+        self.assertSuccess(result)
+        self.compare_counts(result, circuits, targets, delta=0.05 * shots)
+
+    # ---------------------------------------------------------------------
+    # Test ry-gate
+    # ---------------------------------------------------------------------
+    @supported_methods(SUPPORTED_METHODS)
+    def test_ry_gate_deterministic(self, method, device):
+        """Test ry-gate circuits"""
+        backend = self.backend(method=method, device=device, seed_simulator=self.SEED)
+        shots = 1000
+        circuits = ref_rotation.ry_gate_circuits_deterministic(final_measure=True)
+        targets = ref_rotation.ry_gate_counts_deterministic(shots)
+        result = backend.run(circuits, shots=shots).result()
+        self.assertSuccess(result)
+        self.compare_counts(result, circuits, targets, delta=0.05 * shots)

test/terra/primitives/test_estimator.py

@@ -57,7 +57,9 @@ def test_estimator(self, abelian_grouping):
             with self.subTest("PauliSumOp"):
                 observable = PauliSumOp.from_list(lst)
                 ansatz = RealAmplitudes(num_qubits=2, reps=2)
-                est = Estimator(abelian_grouping=abelian_grouping)
+                est = Estimator(
+                    backend_options={"method": "statevector"}, abelian_grouping=abelian_grouping
+                )
                 result = est.run(
                     ansatz, observable, parameter_values=[[0, 1, 1, 2, 3, 5]], seed=15
                 ).result()
@@ -67,7 +69,9 @@ def test_estimator(self, abelian_grouping):
         with self.subTest("SparsePauliOp"):
             observable = SparsePauliOp.from_list(lst)
             ansatz = RealAmplitudes(num_qubits=2, reps=2)
-            est = Estimator(abelian_grouping=abelian_grouping)
+            est = Estimator(
+                backend_options={"method": "statevector"}, abelian_grouping=abelian_grouping
+            )
             result = est.run(
                 ansatz, observable, parameter_values=[[0, 1, 1, 2, 3, 5]], seed=15
             ).result()
@@ -84,7 +88,9 @@ def test_estimator(self, abelian_grouping):
                 ]
             )
             ansatz = RealAmplitudes(num_qubits=2, reps=2)
-            est = Estimator(abelian_grouping=abelian_grouping)
+            est = Estimator(
+                backend_options={"method": "statevector"}, abelian_grouping=abelian_grouping
+            )
             result = est.run(ansatz, observable, parameter_values=[[0] * 6], seed=15).result()
             self.assertIsInstance(result, EstimatorResult)
             np.testing.assert_allclose(result.values, [-0.4], rtol=0.02)

test/terra/primitives/test_sampler.py

@@ -110,7 +110,7 @@ def test_sampler_param_order(self):
         qc.measure(1, 1)
         qc.measure(2, 2)
 
-        sampler = Sampler(backend_options={"seed_simulator": 15})
+        sampler = Sampler(backend_options={"method": "statevector", "seed_simulator": 15})
         result = sampler.run([qc] * 4, [[0, 0], [0, 0], [np.pi / 2, 0], [0, np.pi / 2]]).result()
         self.assertIsInstance(result, SamplerResult)
         self.assertEqual(len(result.quasi_dists), 4)
@@ -140,7 +140,7 @@ def test_sampler_reverse_meas_order(self):
         qc.measure(1, 1)
         qc.measure(2, 0)
 
-        sampler = Sampler()
+        sampler = Sampler(backend_options={"method": "statevector"})
         result = sampler.run(
             [qc, qc, qc, qc], [[0, 0], [0, 0], [np.pi / 2, 0], [0, np.pi / 2]], seed=15
         ).result()