Bravyi-Kitaev

quantum/observable/fermion/FermionOperator.cpp

@@ -29,16 +29,15 @@ FermionTerm &FermionTerm::operator*=(const FermionTerm &v) noexcept {
     auto site = kv.first;
     auto c_or_a = kv.second;
     Operators o = ops();
-    for (auto oo : o ) {
+    for (auto oo : o) {
       if (oo.first == site) {
         if (oo.second == c_or_a) {
           ops().clear();
         } else {
           ops().push_back({site, c_or_a});
         }
         break;
-      }
-      else {
+      } else {
         ops().push_back({site, c_or_a});
         break;
       }
@@ -87,10 +86,18 @@ FermionOperator::FermionOperator(Operators operators, double coeff,
 
 void FermionOperator::clear() { terms.clear(); }
 
-std::vector<std::shared_ptr<CompositeInstruction>> FermionOperator::observe(
-    std::shared_ptr<CompositeInstruction> function) {
-  auto transform = xacc::getService<ObservableTransform>("jw");
-  return transform->transform(xacc::as_shared_ptr(this))->observe(function);
+std::vector<std::shared_ptr<CompositeInstruction>>
+FermionOperator::observe(std::shared_ptr<CompositeInstruction> function,
+                         std::string transformName) {
+
+  auto mapping = xacc::getService<ObservableTransform>(transformName);
+  return mapping->transform(xacc::as_shared_ptr(this))->observe(function);
+}
+
+std::vector<std::shared_ptr<CompositeInstruction>>
+FermionOperator::observe(std::shared_ptr<CompositeInstruction> function) {
+  auto mapping = xacc::getService<ObservableTransform>("jw");
+  return mapping->transform(xacc::as_shared_ptr(this))->observe(function);
 }
 
 const std::string FermionOperator::toString() {
@@ -135,7 +142,8 @@ void FermionOperator::fromString(const std::string str) {
 }
 const int FermionOperator::nBits() {
   auto maxInt = 0;
-  if (terms.empty()) return 0;
+  if (terms.empty())
+    return 0;
 
   for (auto &kv : terms) {
     auto ops = kv.second.ops();
@@ -148,8 +156,8 @@ const int FermionOperator::nBits() {
   return maxInt + 1;
 }
 
-FermionOperator &FermionOperator::operator+=(
-    const FermionOperator &v) noexcept {
+FermionOperator &
+FermionOperator::operator+=(const FermionOperator &v) noexcept {
   FermionOperator vv = v;
   for (auto &kv : v.terms) {
     auto termId = kv.first;
@@ -171,13 +179,13 @@ FermionOperator &FermionOperator::operator+=(
   return *this;
 }
 
-FermionOperator &FermionOperator::operator-=(
-    const FermionOperator &v) noexcept {
+FermionOperator &
+FermionOperator::operator-=(const FermionOperator &v) noexcept {
   return operator+=(-1.0 * v);
 }
 
-FermionOperator &FermionOperator::operator*=(
-    const FermionOperator &v) noexcept {
+FermionOperator &
+FermionOperator::operator*=(const FermionOperator &v) noexcept {
   std::unordered_map<std::string, FermionTerm> newTerms;
   for (auto &kv : terms) {
     for (auto &vkv : v.terms) {
@@ -236,16 +244,16 @@ FermionOperator &FermionOperator::operator*=(const double v) noexcept {
   return operator*=(std::complex<double>(v, 0));
 }
 
-FermionOperator &FermionOperator::operator*=(
-    const std::complex<double> v) noexcept {
+FermionOperator &
+FermionOperator::operator*=(const std::complex<double> v) noexcept {
   for (auto &kv : terms) {
     std::get<0>(kv.second) *= v;
   }
   return *this;
 }
 
-std::shared_ptr<Observable> FermionOperator::commutator(
-    std::shared_ptr<Observable> op) {
+std::shared_ptr<Observable>
+FermionOperator::commutator(std::shared_ptr<Observable> op) {
   FermionOperator &A = *std::dynamic_pointer_cast<FermionOperator>(op);
   std::shared_ptr<FermionOperator> commutatorHA =
       std::make_shared<FermionOperator>((*this) * A - A * (*this));
@@ -279,7 +287,7 @@ void FermionOperator::normalize() {
   double norm = 0.0;
 
   for (auto &kv : terms) {
-    norm += std::pow(std::norm(std::get<0>(kv.second)), 2);
+    norm += std::norm(std::get<0>(kv.second));
   }
   norm = std::sqrt(norm);
 
@@ -293,12 +301,19 @@ void FermionOperator::normalize() {
 double FermionOperator::postProcess(std::shared_ptr<AcceleratorBuffer> buffer,
                                     const std::string &postProcessTask,
                                     const HeterogeneousMap &extra_data) {
-  auto transform = xacc::getService<ObservableTransform>("jw");
-  return transform->transform(xacc::as_shared_ptr(this))->postProcess(buffer, postProcessTask, extra_data);
+
+  std::string mapping = "jw";
+  if (extra_data.stringExists("transform")) {
+    mapping = extra_data.getString("transform");
+  }
+  auto transform = xacc::getService<ObservableTransform>(mapping);
+
+  return transform->transform(xacc::as_shared_ptr(this))
+      ->postProcess(buffer, postProcessTask, extra_data);
 }
 
 } // namespace quantum
-}  // namespace xacc
+} // namespace xacc
 
 bool operator==(const xacc::quantum::FermionOperator &lhs,
                 const xacc::quantum::FermionOperator &rhs) {

quantum/observable/fermion/FermionOperator.hpp

@@ -34,10 +34,12 @@ using Operator = std::pair<int, bool>;
 using Operators = std::vector<Operator>;
 // using SiteMap = std::map<Operators, std::complex<double>>;
 
-using FermionTermTuple = std::tuple<std::complex<double>, Operators, std::string>;
-class FermionTerm : public FermionTermTuple,
-             public tao::operators::commutative_multipliable<FermionTerm>,
-             public tao::operators::equality_comparable<FermionTerm> {
+using FermionTermTuple =
+    std::tuple<std::complex<double>, Operators, std::string>;
+class FermionTerm
+    : public FermionTermTuple,
+      public tao::operators::commutative_multipliable<FermionTerm>,
+      public tao::operators::equality_comparable<FermionTerm> {
 
 public:
   FermionTerm() {
@@ -58,11 +60,10 @@ class FermionTerm : public FermionTermTuple,
     std::get<2>(*this) = "";
   }
 
- FermionTerm(double c) {
-    std::get<0>(*this) = std::complex<double>(c,0.0);
+  FermionTerm(double c) {
+    std::get<0>(*this) = std::complex<double>(c, 0.0);
     std::get<1>(*this) = {};
     std::get<2>(*this) = "";
-
   }
 
   FermionTerm(std::complex<double> c, Operators ops) {
@@ -71,7 +72,7 @@ class FermionTerm : public FermionTermTuple,
     std::get<2>(*this) = "";
   }
 
- FermionTerm(std::complex<double> c, Operators ops, std::string var) {
+  FermionTerm(std::complex<double> c, Operators ops, std::string var) {
     std::get<0>(*this) = c;
     std::get<1>(*this) = ops;
     std::get<2>(*this) = var;
@@ -121,7 +122,7 @@ class FermionTerm : public FermionTermTuple,
   Operators &ops() { return std::get<1>(*this); }
 
   std::complex<double> &coeff() { return std::get<0>(*this); }
-  std::string var() {return std::get<2>(*this);}
+  std::string var() { return std::get<2>(*this); }
 
   FermionTerm &operator*=(const FermionTerm &v) noexcept;
 
@@ -131,7 +132,8 @@ class FermionTerm : public FermionTermTuple,
 };
 
 class FermionOperator
-    : public Observable, public Cloneable<Observable>,
+    : public Observable,
+      public Cloneable<Observable>,
       public std::enable_shared_from_this<FermionOperator>,
       public tao::operators::commutative_ring<FermionOperator>,
       public tao::operators::equality_comparable<FermionOperator>,
@@ -144,7 +146,7 @@ class FermionOperator
 
 public:
   std::shared_ptr<Observable> clone() override {
-      return std::make_shared<FermionOperator>();
+    return std::make_shared<FermionOperator>();
   }
   FermionOperator();
   FermionOperator(std::complex<double> c);
@@ -158,12 +160,19 @@ class FermionOperator
 
   std::vector<std::shared_ptr<CompositeInstruction>>
   observe(std::shared_ptr<CompositeInstruction> function) override;
+
+  std::vector<std::shared_ptr<CompositeInstruction>>
+  observe(std::shared_ptr<CompositeInstruction> function,
+          std::string transformName) override;
+
   const std::string toString() override;
   void fromString(const std::string str) override;
   const int nBits() override;
 
   void clear();
-  std::unordered_map<std::string, FermionTerm> getTerms() const { return terms; }
+  std::unordered_map<std::string, FermionTerm> getTerms() const {
+    return terms;
+  }
 
   FermionOperator &operator+=(const FermionOperator &v) noexcept;
   FermionOperator &operator-=(const FermionOperator &v) noexcept;
@@ -172,16 +181,10 @@ class FermionOperator
   FermionOperator &operator*=(const double v) noexcept;
   FermionOperator &operator*=(const std::complex<double> v) noexcept;
 
-  const std::string name() const override {
-      return "fermion";
-  }
+  const std::string name() const override { return "fermion"; }
 
-  const std::string description() const override {
-      return "";
-  }
-  void fromOptions(const HeterogeneousMap& options) override {
-      return;
-  }
+  const std::string description() const override { return ""; }
+  void fromOptions(const HeterogeneousMap &options) override { return; }
 
   std::shared_ptr<Observable>
   commutator(std::shared_ptr<Observable> obs) override;

quantum/observable/pauli/PauliOperator.cpp

@@ -800,7 +800,7 @@ void PauliOperator::normalize() {
 
   double norm = 0.0;
   for (auto &kv : terms) {
-    norm += std::pow(std::norm(std::get<0>(kv.second)), 2);
+    norm += std::norm(std::get<0>(kv.second));
   }
   norm = std::sqrt(norm);
 

quantum/observable/transforms/BravyiKitaev/BK.cpp

@@ -0,0 +1,106 @@
+/*******************************************************************************
+ * Copyright (c) 2021 UT-Battelle, LLC.
+ * All rights reserved. This program and the accompanying materials
+ * are made available under the terms of the Eclipse Public License v1.0
+ * and Eclipse Distribution License v1.0 which accompanies this
+ * distribution. The Eclipse Public License is available at
+ * http://www.eclipse.org/legal/epl-v10.html and the Eclipse Distribution
+ *License is available at https://eclipse.org/org/documents/edl-v10.php
+ *
+ * Contributors:
+ *    Alexander J. McCaskey - initial API and implementation
+ *    Daniel Claudino - porting
+ ******************************************************************************/
+#include "BK.hpp"
+#include "FermionOperator.hpp"
+#include "PauliOperator.hpp"
+#include "Fenwick.hpp"
+#include "xacc_observable.hpp"
+
+namespace xacc {
+namespace quantum {
+
+template <typename T>
+bool ptr_is_a(std::shared_ptr<Observable> ptr) {
+  return std::dynamic_pointer_cast<T>(ptr) != nullptr;
+}
+
+std::shared_ptr<Observable>
+BK::transform(std::shared_ptr<Observable> Hptr_input) {
+
+  // First we pre-process the observable to a FermionOperator
+  std::shared_ptr<Observable> observable;
+  auto obs_str = Hptr_input->toString();
+  if (ptr_is_a<FermionOperator>(Hptr_input)) {
+    observable = Hptr_input;
+  } else if (obs_str.find("^") != std::string::npos) {
+    observable = xacc::quantum::getObservable("fermion", obs_str);
+  } else {
+    XACCLogger::instance()->error("[Jordan Wigner] Error, cannot cast incoming Observable ptr to something we can process.");
+  }
+  auto fermionObservable =
+      std::dynamic_pointer_cast<FermionOperator>(observable);
+
+  if (!fermionObservable) {
+    XACCLogger::instance()->info(
+        "Cannot execute Bravyi Kitaev on a non-fermion observable.");
+    return observable;
+  }
+
+	int nQubits = observable->nBits();
+
+	FenwickTree tree(nQubits);
+
+  auto result = std::make_shared<PauliOperator>();
+
+  auto terms = fermionObservable->getTerms();
+
+  // Loop over all Fermionic terms...
+  for (auto &kv : terms) {
+
+    auto var = kv.second.var();
+    auto coeff = kv.second.coeff();
+    auto operators = kv.second.ops();
+
+    PauliOperator current(coeff, var);
+    for (auto &kv2 : operators) {
+
+      auto isCreation = kv2.second;
+      int index = kv2.first;
+
+			auto paritySet = tree.getParitySet(index);
+			auto ancestors = tree.getUpdateSet(index);
+			auto ancestorChildren = tree.getRemainderSet(index);
+
+			std::complex<double> dcoeff, ccoeff(.5,0);
+			if (isCreation) {
+				dcoeff = std::complex<double>(0,-.5);
+			} else {
+				dcoeff = std::complex<double>(0,.5);
+			}
+
+			std::map<int, std::string> dTerms{{index, "Y"}}, cTerms{{index, "X"}};
+			for (auto ac : ancestorChildren) {
+				dTerms[ac->index] = "Z";
+			}
+			for (auto p : paritySet) {
+				cTerms[p->index] = "Z";
+			}
+			for (auto a : ancestors) {
+				dTerms[a->index] = "X";
+				cTerms[a->index] = "X";
+			}
+
+			PauliOperator d(dTerms,dcoeff), c(cTerms,ccoeff);
+
+			current *= (c+d);
+    }
+
+    result->operator+=(current);
+  }
+
+  return result;
+}
+
+} // namespace quantum
+} // namespace xacc

quantum/observable/transforms/BravyiKitaev/BK.hpp

@@ -0,0 +1,31 @@
+/*******************************************************************************
+ * Copyright (c) 2021 UT-Battelle, LLC.
+ * All rights reserved. This program and the accompanying materials
+ * are made available under the terms of the Eclipse Public License v1.0
+ * and Eclipse Distribution License v1.0 which accompanies this
+ * distribution. The Eclipse Public License is available at
+ * http://www.eclipse.org/legal/epl-v10.html and the Eclipse Distribution
+ *License is available at https://eclipse.org/org/documents/edl-v10.php
+ *
+ * Contributors:
+ *    Alexander J. McCaskey - initial API and implementation
+ *    Daniel Claudino - porting
+ ******************************************************************************/
+#ifndef XACC_IR_OBSERVABLETRANSFORM_BK_HPP_
+#define XACC_IR_OBSERVABLETRANSFORM_BK_HPP_
+#include "ObservableTransform.hpp"
+
+namespace xacc {
+namespace quantum {
+
+class BK : public ObservableTransform {
+public:
+  std::shared_ptr<Observable>
+  transform(std::shared_ptr<Observable> obs) override;
+  const std::string name() const override { return "bk"; }
+
+  const std::string description() const override { return ""; }
+};
+} // namespace quantum
+} // namespace xacc
+#endif