Add example using PermuteQubits()

In the example 'communication_reduction_via_qubit_reordering.cpp',
we simulate a simple circuit (with gates on the last 10 qubits only)
using the trivial qubit ordering and the inverse one.
With increasing number of MPI processes, the relative cost of the
simulation with inverse qubit order (w.r.t. the simulation with
trivial qubit order) decreases.
This is expected since no MPI communication is needed for the inverse order.

examples/CMakeLists.txt

@@ -18,6 +18,9 @@ target_link_libraries(noisy_circuit_test.exe iqs)
 add_executable(test_of_custom_gates.exe test_of_custom_gates.cpp)
 target_link_libraries(test_of_custom_gates.exe iqs)
 
+add_executable(communication_reduction_via_qubit_reordering.exe communication_reduction_via_qubit_reordering.cpp)
+target_link_libraries(communication_reduction_via_qubit_reordering.exe iqs)
+
 ################################################################################
 
 set_target_properties( benchgates.exe
@@ -26,6 +29,7 @@ set_target_properties( benchgates.exe
                        heisenberg_dynamics.exe
                        noisy_circuit_test.exe
                        test_of_custom_gates.exe
+                       communication_reduction_via_qubit_reordering.exe
     PROPERTIES
     RUNTIME_OUTPUT_DIRECTORY "${CMAKE_CURRENT_SOURCE_DIR}/bin"
 )

examples/communication_reduction_via_qubit_reordering.cpp

@@ -0,0 +1,133 @@
+//------------------------------------------------------------------------------
+// Copyright (C) 2020 Intel Corporation 
+//
+// Licensed under the Apache License, Version 2.0 (the "License");
+// you may not use this file except in compliance with the License.
+// You may obtain a copy of the License at
+//
+// http://www.apache.org/licenses/LICENSE-2.0
+//
+// Unless required by applicable law or agreed to in writing, software
+// distributed under the License is distributed on an "AS IS" BASIS,
+// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+// See the License for the specific language governing permissions and
+// limitations under the License.
+//------------------------------------------------------------------------------
+
+#include <sys/time.h>
+
+#include "../include/qureg.hpp"
+
+#ifdef USE_MKL
+#include <mkl.h>
+#endif
+
+/// @file communication_reduction_via_qubit_reordering.cpp
+/// Simple program to illustrate how the communication overhead between MPI
+/// processes can be reduced by reordering the qubits.
+
+
+/////////////////////////////////////////////////////////////////////////////////////////
+/////////////////////////////////////////////////////////////////////////////////////////
+
+int main(int argc, char **argv)
+{
+  qhipster::mpi::Environment env(argc, argv);
+  if (env.IsUsefulRank() == false) return 0;
+  int myrank = env.GetStateRank();
+  if (myrank==0)
+      std::cout << " ------------------------------------------------------------ \n"
+                << " ---- Reducing communication overhead via qubit reordering -- \n"
+                << " ------------------------------------------------------------ \n";
+
+  int num_qubits = 22;
+  if (argc != 2)
+  {
+    if (myrank==0)
+        fprintf(stderr, "usage: %s <num_qubits>\n", argv[0]);
+    exit(1);
+  }
+  else
+  {
+    num_qubits = atoi(argv[1]);
+  }
+
+  if (num_qubits<22)
+  {
+    if (myrank==0)
+        fprintf(stderr, "simulate circuits of at least 22 qubits\n");
+    exit(1);
+  }
+
+
+#pragma omp parallel
+#pragma omp master
+  {
+    int nthreads = omp_get_num_threads();
+    if (myrank==0)
+      fprintf(stdout, "OMP number of threads = %d \n", nthreads);
+  }
+
+/////////////////////////////////////////////////////////////////////////////////////////
+
+  // Variables for the timing.
+  struct timeval time;
+  double start, end;
+  double duration_trivial_order;
+  double duration_inverse_order;
+
+  // We consider a simple circuit, with Hadamard, X, Y, Z acting on the last 10 qubits.
+  QubitRegister<ComplexDP> psi_trivial_order(num_qubits, "base", 0);
+  psi_trivial_order.ApplyPauliZ(num_qubits-1); // dummy operation to avoid the first MPI communication during timing
+
+  std::vector<std::size_t> original_map = psi_trivial_order.qubit_permutation->map;
+  std::vector<std::size_t> new_map(num_qubits);
+  for (std::size_t q = 0; q < num_qubits; q++)
+  {
+      // Verify that qubit order is trivial.
+      assert(original_map[q] == q);
+      // create inverse map.
+      new_map[num_qubits-q-1] = q;
+  }
+ 
+  // Simulate circuit.
+  qhipster::mpi::StateBarrier();
+  gettimeofday(&time, (struct timezone*)0);
+  start =  time.tv_sec + time.tv_usec * 1.0e-6;
+  for (int q=num_qubits-10; q<num_qubits; ++q)
+  {
+    psi_trivial_order.ApplyHadamard(q);
+    psi_trivial_order.ApplyPauliX(q);
+    psi_trivial_order.ApplyPauliY(q);
+    psi_trivial_order.ApplyPauliZ(q);
+  } 
+  qhipster::mpi::StateBarrier();
+  gettimeofday(&time, (struct timezone*)0);
+  end =  time.tv_sec + time.tv_usec * 1.0e-6;
+  if (myrank==0)
+      std::cout << "trivial qubit order --> Simulation time = " << end-start << "\n";
+  
+/////////////////////////////////////////////////////////////////////////////////////////
+
+  QubitRegister<ComplexDP> psi_inverse_order(num_qubits, "base", 0);
+  psi_inverse_order.PermuteQubits(new_map, "direct");
+  psi_trivial_order.ApplyPauliZ(0); // dummy operation to avoid the first MPI communication during timing
+  // Simulate circuit.
+  qhipster::mpi::StateBarrier();
+  gettimeofday(&time, (struct timezone*)0);
+  start =  time.tv_sec + time.tv_usec * 1.0e-6;
+  for (int q=num_qubits-10; q<num_qubits; ++q)
+  {
+    psi_inverse_order.ApplyHadamard(q);
+    psi_inverse_order.ApplyPauliX(q);
+    psi_inverse_order.ApplyPauliY(q);
+    psi_inverse_order.ApplyPauliZ(q);
+  } 
+  qhipster::mpi::StateBarrier();
+  gettimeofday(&time, (struct timezone*)0);
+  end =  time.tv_sec + time.tv_usec * 1.0e-6;
+  if (myrank==0)
+      std::cout << "inverse qubit order --> Simulation time = " << end-start << "\n";
+             
+  return 0;
+}