MPS: Give the user control over the direction of the internal swaps (#1358)

qiskit/providers/aer/backends/aer_simulator.py

@@ -134,7 +134,7 @@ class AerSimulator(AerBackend):
     +--------------------------+---------------+
     | ``stabilizer``           | No            |
     +--------------------------+---------------+
-    | `"matrix_product_state`` | No            |
+    | ``matrix_product_state`` | No            |
     +--------------------------+---------------+
     | ``extended_stabilizer``  | No            |
     +--------------------------+---------------+
@@ -262,7 +262,7 @@ class AerSimulator(AerBackend):
       alongside setting extended_stabilizer_disable_measurement_opt
       to True (Default: 5000).
 
-    * ``"extended_stabilizer_approximation_error"`` (double): Set the error
+    * ``extended_stabilizer_approximation_error`` (double): Set the error
       in the approximation for the extended_stabilizer method. A
       smaller error needs more memory and computational time
       (Default: 0.05).
@@ -287,7 +287,7 @@ class AerSimulator(AerBackend):
       samples used to estimate probabilities in a probabilities snapshot
       (Default: 3000).
 
-    These backend options only apply when using the ``"matrix_product_state"``
+    These backend options only apply when using the ``matrix_product_state``
     simulation method:
 
     * ``matrix_product_state_max_bond_dimension`` (int): Sets a limit
@@ -303,13 +303,13 @@ class AerSimulator(AerBackend):
     * ``mps_sample_measure_algorithm`` (str): Choose which algorithm to use for
       ``"sample_measure"`` (Default: "mps_apply_measure").
 
-      - ``"mps_probabilities"``: This method first constructs the probability
+      - ``mps_probabilities``: This method first constructs the probability
         vector and then generates a sample per shot. It is more efficient for
         a large number of shots and a small number of qubits, with complexity
         O(2^n * n * D^2) to create the vector and O(1) per shot, where n is
         the number of qubits and D is the bond dimension.
 
-      - ``"mps_apply_measure"``: This method creates a copy of the mps structure
+      - ``mps_apply_measure``: This method creates a copy of the mps structure
         and measures directly on it. It is more efficient for a small number of
         shots, and a large number of qubits, with complexity around
         O(n * D^2) per shot.
@@ -318,6 +318,11 @@ class AerSimulator(AerBackend):
       structure: bond dimensions and values discarded during approximation.
       (Default: False)
 
+    * ``mps_swap_direction`` (str): Determine the direction of swapping the
+      qubits when internal swaps are inserted for a 2-qubit gate.
+      Possible values are "mps_swap_right" and "mps_swap_left".
+      (Default: "mps_swap_left")
+
     These backend options apply in circuit optimization passes:
 
     * ``fusion_enable`` (bool): Enable fusion optimization in circuit
@@ -488,6 +493,7 @@ class AerSimulator(AerBackend):
             matrix_product_state_max_bond_dimension=None,
             mps_sample_measure_algorithm='mps_heuristic',
             mps_log_data=False,
+            mps_swap_direction='mps_swap_left',
             chop_threshold=1e-8,
             mps_parallel_threshold=14,
             mps_omp_threads=1)

releasenotes/notes/mps_improve_swaps-e7a142ca3b8644cd.yaml

@@ -0,0 +1,6 @@
+other:
+  - |
+    Added an option that allows the user to determine the direction of 
+    swapping the qubits when internal swaps are inserted for a 2-qubit gate.
+    Possible values are "mps_swap_right" and "mps_swap_left".
+    The direction of the swaps may affect performance, depending on the circuit.

src/simulators/matrix_product_state/matrix_product_state.hpp

@@ -490,6 +490,15 @@ void State::set_config(const json_t &config) {
   bool mps_log_data;
   if (JSON::get_value(mps_log_data, "mps_log_data", config))
     MPS::set_mps_log_data(mps_log_data);
+
+// Set the direction for the internal swaps
+  std::string direction;
+  if (JSON::get_value(direction, "mps_swap_direction", config)) {
+    if (direction.compare("mps_swap_right") == 0)
+      MPS::set_mps_swap_direction(MPS_swap_direction::SWAP_RIGHT);
+    else
+      MPS::set_mps_swap_direction(MPS_swap_direction::SWAP_LEFT);
+  }
 }
 
 void State::add_metadata(ExperimentResult &result) const {

src/simulators/matrix_product_state/matrix_product_state_internal.cpp

@@ -41,6 +41,7 @@ static const cmatrix_t one_measure =
   uint_t MPS::omp_threads_ = 1;     
   uint_t MPS::omp_threshold_ = 14;  
   enum Sample_measure_alg MPS::sample_measure_alg_ = Sample_measure_alg::HEURISTIC; 
+  enum MPS_swap_direction MPS::mps_swap_direction_ = MPS_swap_direction::SWAP_LEFT;
   double MPS::json_chop_threshold_ = 1E-8;
   std::stringstream MPS::logging_str_;
   bool MPS::mps_log_data_ = 0;
@@ -589,26 +590,36 @@ void MPS::print_bond_dimensions() const {
 //    V is split by columns to yield two MPS_Tensors representing qubit B (in reshape_V_after_SVD),
 //    the diagonal of S becomes the Lambda-vector in between A and B.
 //-------------------------------------------------------------------------
-void MPS::apply_2_qubit_gate(uint_t index_A, uint_t index_B, Gates gate_type, const cmatrix_t &mat, bool is_diagonal)
+void MPS::apply_2_qubit_gate(uint_t index_A, uint_t index_B, 
+			     Gates gate_type, const cmatrix_t &mat, 
+			     bool is_diagonal)
 {
   // We first move the two qubits to be in consecutive positions
-  // If index_B > index_A, we move the qubit at index_B to index_A+1
-  // If index_B < index_A, we move the qubit at index_B to index_A-1, and then
-  // swap between the qubits
-  uint_t A = index_A;
+  // By default, the right qubit is moved to the position after the left qubit.
+  // However, the user can choose to move the left qubit to be in the position 
+  // before the right qubit by changing the MPS_swap_direction to SWAP_RIGHT.
+  // The direction of the swaps may affect performance, depending on the circuit.
 
   bool swapped = false;
+  uint_t low_qubit=0, high_qubit=0;
 
-  if (index_B > index_A+1) {
-    change_position(index_B, index_A+1);  // Move B to be right after A
-  } else if (index_A > 0 && index_B < index_A-1) {
-    change_position(index_B, index_A-1);  // Move B to be right before A
-  }
-  if (index_B < index_A) {
-    A = index_A - 1;
+  if (index_B > index_A) {
+    low_qubit = index_A;
+    high_qubit = index_B;
+  } else {
+    low_qubit = index_B;
+    high_qubit = index_A;
     swapped = true;
   }
-  common_apply_2_qubit_gate(A, gate_type, mat, swapped, is_diagonal);
+  if (mps_swap_direction_ == MPS_swap_direction::SWAP_LEFT) {
+    // Move high_qubit to be right after low_qubit
+      change_position(high_qubit, low_qubit+1);  
+  } else {  //mps_swap_right
+    // Move low_qubit to be right before high_qubit
+      change_position(low_qubit, high_qubit-1);  
+      low_qubit = high_qubit-1;
+    }
+  common_apply_2_qubit_gate(low_qubit, gate_type, mat, swapped, is_diagonal);
 }
 
 void MPS::common_apply_2_qubit_gate(uint_t A,  // the gate is applied to A and A+1

src/simulators/matrix_product_state/matrix_product_state_internal.hpp

@@ -37,6 +37,7 @@ enum Gates {
   //enum class Direction {RIGHT, LEFT};
 
   enum class Sample_measure_alg {APPLY_MEASURE, PROB, MEASURE_ALL, HEURISTIC};
+  enum class MPS_swap_direction {SWAP_LEFT, SWAP_RIGHT};
 
 //=========================================================================
 // MPS class
@@ -286,6 +287,10 @@ public:
     mps_log_data_ = mps_log_data;
   }
 
+  static void set_mps_swap_direction(MPS_swap_direction direction) {
+    mps_swap_direction_ = direction;
+  }
+
   static uint_t get_omp_threads() {
     return omp_threads_;
   }
@@ -307,6 +312,10 @@ public:
     return mps_log_data_;
   }
 
+  static MPS_swap_direction get_swap_direction() {
+    return mps_swap_direction_;
+  }
+
   //----------------------------------------------------------------
   // Function name: norm
   // Description: the norm is defined as <psi|A^dagger . A|psi>.
@@ -526,6 +535,7 @@ private:
   static bool enable_gate_opt_;      // allow optimizations on gates
   static std::stringstream logging_str_;
   static bool mps_log_data_;
+  static MPS_swap_direction mps_swap_direction_;
 };
 
 inline std::ostream &operator<<(std::ostream &out, const rvector_t &vec) {

test/terra/backends/aer_simulator/test_options.py

@@ -147,11 +147,14 @@ class TestOptions(SimulatorTestCase):
         value = sum(result.get_counts().values())
         self.assertEqual(value, shots)
 
-    def test_mps_approximation(self):
-        """Test MPS approximation"""
+    def test_mps_options(self):
+        """Test MPS options"""
         shots = 4000
         method="matrix_product_state"
-        backend_exact = self.backend(method=method)
+        backend_swap_left = self.backend(method=method,
+                                          mps_swap_direction='mps_swap_left')
+        backend_swap_right = self.backend(method=method,
+                                          mps_swap_direction='mps_swap_right')
         backend_approx = self.backend(method=method,
                                      matrix_product_state_max_bond_dimension=8)
         # The test must be large enough and entangled enough so that
@@ -165,13 +168,21 @@ class TestOptions(SimulatorTestCase):
                 circuit.cx(0, i)
         circuit.save_statevector('sv')
 
-        result_exact = backend_exact.run(circuit, shots=shots).result()
-        sv_exact = result_exact.data(0)['sv']
+        result_swap_left = backend_swap_left.run(circuit, shots=shots).result()
+        sv_left = result_swap_left.data(0)['sv']
+        
+        result_swap_right = backend_swap_right.run(circuit, shots=shots).result()
+        sv_right = result_swap_right.data(0)['sv']
+        
         result_approx = backend_approx.run(circuit, shots=shots).result()
         sv_approx = result_approx.data(0)['sv']
-        # Check that the fidelity is reasonable
-        self.assertGreaterEqual(state_fidelity(sv_exact, sv_approx), 0.80)
+
+        # swap_left and swap_right should give the same state vector
+        self.assertAlmostEqual(state_fidelity(sv_left, sv_right), 1.0)
+        
+        # Check that the fidelity of approximation is reasonable
+        self.assertGreaterEqual(state_fidelity(sv_left, sv_approx), 0.80)
 
         # Check that the approximated result is not identical to the exact
         # result, because that could mean there was actually no approximation
-        self.assertLessEqual(state_fidelity(sv_exact, sv_approx), 0.999)
+        self.assertLessEqual(state_fidelity(sv_left, sv_approx), 0.999)