Allow measure sampling with initialize (#1219)

If initialize is the first instruction in a circuit, or is defined across the full width of a circuit, it is deterministic and hence is compatible with the measurement sampling optimization

releasenotes/notes/fix-initialize-sampling-b6948e39112f6a46.yaml

@@ -0,0 +1,7 @@
+---
+fixes:
+  - |
+    Fixes bug where the initialize instruction would disable measurement
+    sampling optimization for the statevector and matrix product state
+    simulation methods even when it was the first circuit instruction or
+    applied to all qubits and hence deterministic.

src/controllers/aer_controller.hpp

@@ -2034,21 +2034,26 @@ bool Controller::check_measure_sampling_opt(const Circuit &circ,
     return false;
   }
 
-  // If density matrix or superop method all noise instructions allow
-  // Sampling
-  if (method == Method::density_matrix || method == Method::unitary) {
+  // If density matrix, unitary, superop method all supported instructions
+  // allow sampling
+  if (method == Method::density_matrix ||
+      method == Method::superop ||
+      method == Method::unitary) {
     return true;
   }
+
+  // If circuit contains a non-initial initialize that is not a full width
+  // instruction we can't sample
+  if (circ.can_sample_initialize == false) {
+    return false;
+  }
 
   // Check if non-density matrix simulation and circuit contains
   // a stochastic instruction before measurement
-  // ie. initialize, reset, kraus, superop, conditional
+  // ie. reset, kraus, superop
   // TODO:
-  // * If initialize should be allowed if applied to product states (ie start of
-  // circuit)
   // * Resets should be allowed if applied to |0> state (no gates before).
   if (circ.opset().contains(Operations::OpType::reset) ||
-      circ.opset().contains(Operations::OpType::initialize) ||
       circ.opset().contains(Operations::OpType::kraus) ||
       circ.opset().contains(Operations::OpType::superop)) {
     return false;

src/controllers/qasm_controller.hpp

@@ -1216,20 +1216,27 @@ bool QasmController::check_measure_sampling_opt(const Circuit& circ,
     return false;
   }
 
+  // If density matrix method all noise instructions allowSampling
+  if (method == Method::density_matrix ||
+      method == Method::density_matrix_thrust_gpu ||
+      method == Method::density_matrix_thrust_cpu) {
+    return true;
+  }
+
+  // If circuit contains a non-initial initialize that is not a full width
+  // instruction we can't sample
+  if (circ.can_sample_initialize == false) {
+    return false;
+  }
+
   // Check if non-density matrix simulation and circuit contains
   // a stochastic instruction before measurement
-  // ie. initialize, reset, kraus, superop, conditional
+  // ie. reset, kraus, superop
   // TODO:
-  // * If initialize should be allowed if applied to product states (ie start of
-  // circuit)
   // * Resets should be allowed if applied to |0> state (no gates before).
-  bool density_mat = (method == Method::density_matrix ||
-                      method == Method::density_matrix_thrust_gpu ||
-                      method == Method::density_matrix_thrust_cpu);
-  if (!density_mat && (circ.opset().contains(Operations::OpType::reset) ||
-                       circ.opset().contains(Operations::OpType::initialize) ||
-                       circ.opset().contains(Operations::OpType::kraus) ||
-                       circ.opset().contains(Operations::OpType::superop))) {
+  if (circ.opset().contains(Operations::OpType::reset) ||
+      circ.opset().contains(Operations::OpType::kraus) ||
+      circ.opset().contains(Operations::OpType::superop)) {
     return false;
   }
   // Otherwise true

src/framework/circuit.hpp

@@ -42,9 +42,11 @@ class Circuit {
   uint_t num_registers = 0;     // maximum number of registers clbits needed for ops
 
   // Measurement params
-  bool has_conditional = false; // True if any ops are conditional
-  bool can_sample = true;       // True if circuit tail contains measure, roerror, barrier.
-  size_t first_measure_pos = 0; // Position of first measure instruction
+  bool has_conditional = false;      // True if any ops are conditional
+  bool can_sample = true;            // True if circuit tail contains measure, roerror, barrier.
+  size_t first_measure_pos = 0;      // Position of first measure instruction
+  bool can_sample_initialize = true; // True if circuit contains at most 1 initialize
+                                     // and it is the first instruction in the circuit
 
   // Circuit metadata constructed from json QobjExperiment
   uint_t shots = 1;
@@ -123,6 +125,7 @@ void Circuit::set_params() {
   // Memory size is loaded from qobj config
   // Also check if measure sampling is in principle possible
   bool first_measure = true;
+  size_t initialize_qubits = 0;
   for (size_t i = 0; i < ops.size(); ++i) {
     const auto &op = ops[i];
     has_conditional |= op.conditional;
@@ -140,6 +143,13 @@ void Circuit::set_params() {
       }
     }
 
+    // Keep track of minimum width of any non-initial initialize instructions
+    if (i > 0 && op.type == OpType::initialize) {
+      initialize_qubits = (initialize_qubits == 0)
+        ? op.qubits.size()
+        : std::min(op.qubits.size(), initialize_qubits);
+    }
+
     // Compute measure sampling check
     if (can_sample) {
       if (first_measure) {
@@ -162,6 +172,11 @@ void Circuit::set_params() {
   num_qubits = (qubitset_.empty()) ? 0 : 1 + *qubitset_.rbegin();
   num_memory = (memoryset_.empty()) ? 0 : 1 + *memoryset_.rbegin();
   num_registers = (registerset_.empty()) ? 0 : 1 + *registerset_.rbegin();
+
+  // Check that any non-initial initialize is defined on full width of qubits
+  if (initialize_qubits > 0 && initialize_qubits < num_qubits) {
+    can_sample_initialize = false;
+  }
 }
 
 Circuit::Circuit(const std::vector<Op> &_ops) : Circuit() {

test/terra/backends/qasm_simulator/qasm_initialize.py

@@ -13,6 +13,7 @@
 QasmSimulator Integration Tests
 """
 
+from qiskit import QuantumCircuit
 from test.terra.reference import ref_initialize
 from qiskit.compiler import assemble
 from qiskit.providers.aer import QasmSimulator
@@ -115,3 +116,32 @@ def test_initialize_entangled_qubits(self):
         result = self.SIMULATOR.run(qobj, **opts).result()
         self.assertSuccess(result)
         self.compare_counts(result, circuits, targets, delta=0.05 * shots)
+
+    def test_initialize_sampling_opt_enabled(self):
+        """Test sampling optimization"""
+        shots = 1000
+        circuit = QuantumCircuit(2)
+        circuit.initialize([0, 1], [1])
+        circuit.h([0, 1])
+        circuit.initialize([0, 0, 1, 0], [0, 1])
+        circuit.measure_all()
+        qobj = assemble(circuit, self.SIMULATOR, shots=shots)
+        opts = self.BACKEND_OPTS.copy()
+        result = self.SIMULATOR.run(qobj, **opts).result()
+        self.assertSuccess(result)
+        sampling = result.results[0].metadata.get('measure_sampling', None)
+        self.assertTrue(sampling)
+
+    def test_initialize_sampling_opt_disabled(self):
+        """Test sampling optimization"""
+        shots = 1000
+        circuit = QuantumCircuit(2)
+        circuit.h([0, 1])
+        circuit.initialize([0, 1], [1])
+        circuit.measure_all()
+        qobj = assemble(circuit, self.SIMULATOR, shots=shots)
+        opts = self.BACKEND_OPTS.copy()
+        result = self.SIMULATOR.run(qobj, **opts).result()
+        self.assertSuccess(result)
+        sampling = result.results[0].metadata.get('measure_sampling', None)
+        self.assertFalse(sampling)