cleanup on qpt integration with benchmarking suite, fixed bug in setting physical qubits for qpt

Signed-off-by: Alex McCaskey <mccaskeyaj@ornl.gov>

python/benchmark/qpt/benchmarks/ibm_poughkeepsie_hadamard_qpt_1.ini

+[XACC]
+accelerator = ibm:ibmq_poughkeepsie
+verbose = True
+
+[Benchmark]
+name = qpt
+analysis = ['fidelity', 'heat-maps']
+chi-theoretical-real = [0., 0., 0., 0., 0., 1., 0., 1., 0., 0., 0., 0., 0., 1., 0., 1.]
+
+[Circuit]
+# Logical circuit source code
+source = .compiler xasm
+    .circuit hadamard
+    .qbit q
+    H(q[0]);
+
+# Can specify physical qubit to run on
+qubit-map = [1]
+

python/benchmark/qpt/qpt_benchmark.py

@@ -18,6 +18,7 @@ class QPT(Benchmark):
         self.qpt = None
         self.nq = None
         self.qpu = None
+        self.qubit_map = []
 
     def execute(self, inputParams):
         xacc_opts = inputParams['XACC']
@@ -54,7 +55,11 @@ class QPT(Benchmark):
 
         opts = {'circuit':ansatz, 'accelerator':self.qpu}
         if 'qubit-map' in inputParams['Circuit']:
-            opts['qubit-map'] = ast.literal_eval(inputParams['Circuit']['qubit-map'])
+            raw_qbit_map = inputParams['Circuit']['qubit-map']
+            if not isinstance(raw_qbit_map, list):
+                raw_qbit_map = ast.literal_eval(raw_qbit_map)
+            self.qubit_map = raw_qbit_map
+            opts['qubit-map'] = self.qubit_map
 
         self.qpt = xacc.getAlgorithm('qpt', opts)
 
@@ -69,6 +74,9 @@ class QPT(Benchmark):
     def analyze(self, buffer, inputParams):
         import numpy as np, ast
 
+        if self.qubit_map:
+            buffer.addExtraInfo('qubit-map', self.qubit_map)
+
         # Can get the chi process matrix and print it
         chi_real_part = np.array(buffer['chi-real'])
         chi_imag_part = np.array(buffer['chi-imag'])
@@ -76,13 +84,17 @@ class QPT(Benchmark):
         print('\nComputed Chi Process:\n', chi)
 
         if 'fidelity' in inputParams['Benchmark']['analysis']:
-            opts = {'chi-theoretical-real':ast.literal_eval(inputParams['Benchmark']['chi-theoretical-real'])}
+            raw_chi_real = inputParams['Benchmark']['chi-theoretical-real']
+            if not isinstance(raw_chi_real, list):
+                raw_chi_real = ast.literal_eval(inputParams['Benchmark']['chi-theoretical-real'])
+            opts = {'chi-theoretical-real':raw_chi_real}
             if 'chi-theoretical-imag' in inputParams['Benchmark']:
                 opts['chi-theoretical-imag'] = ast.literal_eval(inputParams['Benchmark']['chi-theoretical-imag'])
 
             # Compute the fidelity with respect
             F = self.qpt.calculate('fidelity', buffer, opts)
             print('\nComputed Process Fidelity: ', F)
+            buffer.addExtraInfo('fidelity', F)
 
         if 'heat-maps' in inputParams['Benchmark']['analysis']:
             import matplotlib.pyplot as plt
@@ -95,9 +107,12 @@ class QPT(Benchmark):
             ax1, ax2 = axes
             rc('text', usetex=True)
 
+            raw_chi_real = inputParams['Benchmark']['chi-theoretical-real']
+            if not isinstance(raw_chi_real, list):
+                raw_chi_real = ast.literal_eval(inputParams['Benchmark']['chi-theoretical-real'])
             # Heat maps.
             im1 = ax1.matshow(np.real(chi), cmap='viridis',interpolation='nearest')
-            im2 = ax2.matshow(np.abs(np.real(chi) - np.reshape(np.array(ast.literal_eval(inputParams['Benchmark']['chi-theoretical-real'])), (2**(2*self.nq), 2**(2*self.nq)))), cmap='viridis')
+            im2 = ax2.matshow(np.abs(np.real(chi) - np.reshape(np.array(raw_chi_real), (2**(2*self.nq), 2**(2*self.nq)))), cmap='viridis')
 
             # Formatting for heat map 1.
             ax1.set_xticks(range(2**(2*self.nq)))
@@ -124,6 +139,6 @@ class QPT(Benchmark):
             from datetime import datetime
             now = datetime.now() # current date and time
             date_time = now.strftime("%m_%d_%Y_%H_%M_%S")
-            plt.savefig('qpt_{}_{}.png'.format(self.circuit_name, date_time))
-            print('\nHeat-maps output to qpt_{}_{}.png\n'.format(self.circuit_name, date_time))
+            plt.savefig('qpt_{}_{}_{}.png'.format(self.circuit_name, date_time, '_'.join([str(q) for q in self.qubit_map])))
+            print('\nHeat-maps output to qpt_{}_{}_{}.png\n'.format(self.circuit_name, date_time, '_'.join([str(q) for q in self.qubit_map])))
 

python/examples/cx_qpt.py

+import xacc
+import numpy as np
+xacc.set_verbose(True)
+# Choose the QPU on which to
+# characterize the process matrix for a Hadamard
+qpu = xacc.getAccelerator('ibm:ibmq_poughkeepsie')
+
+# Create the CompositeInstruction containing a
+# single Hadamard instruction
+provider = xacc.getIRProvider('quantum')
+circuit = provider.createComposite('U')
+cnot = provider.createInstruction('CX', [0,1])
+circuit.addInstruction(cnot)
+
+# Create the Algorithm, give it the circuit
+# to characterize and the backend to target
+qpt = xacc.getAlgorithm('qpt', {'circuit':circuit, 'accelerator':qpu, 'qubit-map':[1,2]})
+
+# Allocate a qubit, this will
+# store our tomography results
+buffer = xacc.qalloc(2)
+
+# Execute
+qpt.execute(buffer)
+
+theory_cx = [1., 0., 0., 1., 1., 0., 0., -1., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0.,
+        1., 1., 0., 0., -1., 0., 0., 0., 0., 0., 0., 0., 0., 1., 0., 0., 1., 1.,
+        0., 0., -1., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., -1., 0., 0., -1., -1., 0., 0., 1., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+        0., 0., 0., 0., 0., 0., 0., 0.]
+
+# Compute the fidelity with respect to
+# the ideal hadamard process
+F = qpt.calculate('fidelity', buffer, {'chi-theoretical-real':theory_cx})
+print('\nFidelity: ', F)
+
+# Can get the chi process matrix and print it
+chi_real_part = np.array(buffer['chi-real'])
+chi_imag_part = np.array(buffer['chi-imag'])
+chi = np.reshape(chi_real_part + 1j*chi_imag_part, (16,16))
+print('\nChi:\n', chi)
+

python/examples/hadamard_qpt.py

 import xacc
 import numpy as np
-
+xacc.set_verbose(True)
 # Choose the QPU on which to
 # characterize the process matrix for a Hadamard
-qpu = xacc.getAccelerator('aer')
+qpu = xacc.getAccelerator('ibm:ibmq_poughkeepsie')
 
 # Create the CompositeInstruction containing a
 # single Hadamard instruction
@@ -14,7 +14,7 @@ circuit.addInstruction(hadamard)
 
 # Create the Algorithm, give it the circuit
 # to characterize and the backend to target
-qpt = xacc.getAlgorithm('qpt', {'circuit':circuit, 'accelerator':qpu})
+qpt = xacc.getAlgorithm('qpt', {'circuit':circuit, 'accelerator':qpu}) # map logical 0 to physical 1
 
 # Allocate a qubit, this will
 # store our tomography results

python/xacc.py

@@ -362,18 +362,7 @@ if not pelix.framework.FrameworkFactory.is_framework_running(None):
     serviceRegistry = PyServiceRegistry()
     serviceRegistry.initialize()
 
-
-def benchmark(opts):
-    if opts.benchmark is not None:
-        inputfile = opts.benchmark
-        config = configparser.ConfigParser()
-        config.read(inputfile)
-        xacc_settings = {section: dict(config.items(section)) for section in config.sections()}
-    else:
-        error('Must provide input file for benchmark.')
-        return
-
-    Initialize()
+def benchmark(xacc_settings):
 
     if 'benchmark' not in serviceRegistry.registry:
         error('No benchmarks available')
@@ -389,6 +378,8 @@ def benchmark(opts):
               xacc_settings['benchmark'] + " is not installed.")
         exit(1)
 
+    inputfile = xacc_settings['input-file-name'] if 'input-file-name' in xacc_settings else None
+    
     starttime = time.time()
     buffer = _benchmark.execute(xacc_settings)
     elapsedtime = time.time() - starttime
@@ -396,16 +387,33 @@ def benchmark(opts):
     # for k, v in xacc_settings.items():
     #     buffer.addExtraInfo(k, v)
     # Analyze the buffer
-    head, tail = os.path.split(inputfile)
-    buffer.addExtraInfo('file-name', tail)
-    _benchmark.analyze(buffer, xacc_settings)
     timestr = time.strftime("%Y%m%d-%H%M%S")
+    tail = '{}_{}'.format(xacc_settings['Benchmark']['name'], xacc_settings['XACC']['accelerator'].replace(':','_'))
+    if inputfile is not None:
+        head, tail = os.path.split(inputfile)
+        buffer.addExtraInfo('file-name', tail)
+
+    _benchmark.analyze(buffer, xacc_settings)
+
     results_name = "%s_%s_out" % (os.path.splitext(
         tail)[0], timestr)
     f = open(results_name+".ab", 'w')
     f.write(str(buffer))
     f.close()
 
+def benchmark_from_cmd_line(opts):
+    if opts.benchmark is not None:
+        inputfile = opts.benchmark
+        config = configparser.ConfigParser()
+        config.read(inputfile)
+        xacc_settings = {section: dict(config.items(section)) for section in config.sections()}
+        xacc_settings['input-file-name'] = inputfile
+    else:
+        error('Must provide input file for benchmark.')
+        return
+
+    Initialize()
+    benchmark(xacc_settings)
     Finalize()
 
 
@@ -429,7 +437,7 @@ def main(argv=None):
         setCredentials(opts)
 
     if not opts.benchmark == None:
-        benchmark(opts)
+        benchmark_from_cmd_line(opts)
 
 initialize()
 

quantum/gate/ir/Circuit.hpp

@@ -93,11 +93,19 @@ public:
 
   const std::string getTag() override {return "";}
   void setTag(const std::string& tag) override {return;}
-  
+
   void mapBits(std::vector<std::size_t> bitMap) override {
-    for (auto &inst : instructions) {
-      inst->mapBits(bitMap);
+    InstructionIterator iter(shared_from_this());
+    while(iter.hasNext()) {
+        auto next = iter.next();
+        if (!next->isComposite()) {
+            next->mapBits(bitMap);
+        }
     }
+
+    // for (auto &inst : instructions) {
+    // #  inst->mapBits(bitMap);
+    // }
   }
   void setBits(const std::vector<std::size_t> bits) override {}
   const std::vector<std::size_t> bits() override { return {}; }
@@ -363,7 +371,12 @@ public:
 
   DEFINE_VISITABLE()
   std::shared_ptr<Instruction> clone() override {
-    return std::make_shared<Circuit>(*this);
+    auto cloned = std::make_shared<Circuit>(name(), variables);
+
+    for (auto i : instructions) {
+        cloned->addInstruction(i->clone());
+    }
+    return cloned;//std::make_shared<Circuit>(*this);
   }
 
   virtual ~Circuit() {}

quantum/gate/ir/Gate.hpp

@@ -75,7 +75,7 @@ public:
 
 #define DEFINE_CLONE(CLASS)                                                    \
   std::shared_ptr<Instruction> clone() override {                              \
-    return std::make_shared<CLASS>();                                          \
+    return std::make_shared<CLASS>(*this);                                     \
   }
 };
 

quantum/observable/pauli/PauliOperator.cpp

@@ -110,7 +110,7 @@ PauliOperator::observe(std::shared_ptr<CompositeInstruction> function) {
     gateFunction->setCoefficient(spinInst.coeff());
 
     if (function->hasChildren()) {
-      gateFunction->addInstruction(function);
+      gateFunction->addInstruction(function->clone());
     }
 
     // Loop over all terms in the Spin Instruction

quantum/plugins/algorithms/qpt/qpt.cpp

@@ -209,10 +209,15 @@ void QPT::execute(const std::shared_ptr<AcceleratorBuffer> buffer) const {
       if (i > 0) {
         prep->setBits({i});
       }
+
+
       prep_on_all_qbits->addInstructions(prep->getInstructions());
     }
 
-    prep_on_all_qbits->addInstruction(circuit_as_shared);
+    // Add circuit to be characterized
+    for (auto& inst : circuit_as_shared->getInstructions()) {
+       prep_on_all_qbits->addInstruction(inst->clone());
+    }
 
     // Create Prep Operator
     Eigen::MatrixXcd prep;
@@ -230,19 +235,6 @@ void QPT::execute(const std::shared_ptr<AcceleratorBuffer> buffer) const {
     for (auto &basis_element : basis_observables) {
       auto prep_on_all_and_meas = basis_element.observe(prep_on_all_qbits)[0];
 
-      // Observable::observe will change basis to {X,Y,Z}, we want to
-      // measure Y.
-      for (int i = 0; i < prep_on_all_and_meas->nInstructions(); i++) {
-        auto inst = prep_on_all_and_meas->getInstruction(i);
-        if (inst->name() == "Rx") {
-          auto sdg = provider->createInstruction("Sdg", inst->bits());
-          auto h = provider->createInstruction("H", inst->bits());
-          auto c = provider->createComposite("__tmp_sdg_h");
-          c->addInstructions({sdg, h});
-          prep_on_all_and_meas->replaceInstruction(i, c);
-        }
-      }
-
       std::stringstream ss;
       std::vector<std::string> measure_label;
       std::vector<int> idxs;
@@ -256,6 +248,7 @@ void QPT::execute(const std::shared_ptr<AcceleratorBuffer> buffer) const {
       ss << measure_label;
 
       names.push_back(comp_name + "_" + ss.str());
+      prep_on_all_and_meas->setName(comp_name + "_"+ss.str());
 
       all_circuits.push_back(prep_on_all_and_meas);
 
@@ -307,6 +300,12 @@ void QPT::execute(const std::shared_ptr<AcceleratorBuffer> buffer) const {
 
   Eigen::MatrixXcd basis_matrix = VStack(blocks);
 
+  // Perform default circuit optimization
+//   auto optimizer = xacc::getIRTransformation("circuit-optimizer");
+//   for (auto circuit : all_circuits) {
+//     optimizer->apply(circuit, nullptr);
+//   }
+
   // Map to physical qubits if specified
   if (!qubit_map.empty()) {
       auto placement = xacc::getIRTransformation("default-placement");
@@ -320,13 +319,18 @@ void QPT::execute(const std::shared_ptr<AcceleratorBuffer> buffer) const {
   auto children = buffer->getChildren();
 
   std::vector<std::complex<double>> data_vec;
+  auto xasm = xacc::getCompiler("xasm");
   int i = 0;
   for (auto &child : children) {
+    auto circuit_src = xasm->translate(all_circuits[i]);
+    child->addExtraInfo("xasm-src", circuit_src);
+
     for (auto bit_string : all_bit_strings) {
       // FIXME NEEDED TO REVERSE HERE TO REPRODUCE RESULTS
       if (qpu->getBitOrder() == Accelerator::BitOrder::MSB) {
           std::reverse(bit_string.begin(), bit_string.end());
       }
+    //   std::cout << names[i] << ", " << child->getMeasurementCounts() << ", " << bit_string << "\n";
       auto prob = child->computeMeasurementProbability(bit_string);
       data_vec.push_back(prob);
     }

quantum/plugins/algorithms/qpt/tests/QPTTester.cpp

@@ -89,7 +89,7 @@ TEST(QPTTester, checkHadamard) {
 
 TEST(QPTTester, checkCX) {
   if (xacc::hasAccelerator("aer")) {
-    auto acc = xacc::getAccelerator("qpp", {std::make_pair("shots", 1024)});
+    auto acc = xacc::getAccelerator("aer", {std::make_pair("shots", 1024)});
     auto buffer = xacc::qalloc(2);
 
     auto compiler = xacc::getCompiler("xasm");

quantum/plugins/ibm/accelerator/IBMAccelerator.cpp

@@ -227,6 +227,8 @@ void IBMAccelerator::execute(
   auto job_id = reserve_response_json["id"].get<std::string>();
   currentJobId = job_id;
 
+  buffer->addExtraInfo("ibm-job-id", job_id);
+  
   // Now we ask IBM for an upload URL for the QObj
   auto job_upload_url_response =
       get(IBM_API_URL, IBM_CREDENTIALS_PATH + "/Jobs/" + job_id +