Commit c5838c3b authored by Nguyen, Thien Minh's avatar Nguyen, Thien Minh
Browse files

Added Hamiltonian stretching/translation for IQPE 



Also, recompute the energy from phase estimation.

Signed-off-by: default avatarThien Nguyen <nguyentm@ornl.gov>
parent 6fed858d

lib/qsim/impls/qsim_impl.cpp

+42 −15
Original line number Diff line number Diff line
@@ -132,15 +132,12 @@ bool IterativeQpeWorkflow::initialize(const HeterogeneousMap &params) { 
  return (num_steps >= 1) && (num_iters >= 1);

}



std::shared_ptr<CompositeInstruction>

IterativeQpeWorkflow::constructQpeCircuit(const QuantumSimulationModel &model,

                                          int k, double omega,

                                          bool measure) const {

std::shared_ptr<CompositeInstruction> IterativeQpeWorkflow::constructQpeCircuit(

    std::shared_ptr<Observable> obs, int k, double omega, bool measure) const {

  auto provider = xacc::getIRProvider("quantum");



  auto kernel = provider->createComposite("__TEMP__QPE__LOOP__");

  

  const auto nbQubits = model.observable->nBits();

  const auto nbQubits = obs->nBits();

  const size_t ancBit = nbQubits;

  // Hadamard on ancilla qubit

  kernel->addInstruction(provider->createInstruction("H", ancBit));
@@ -150,7 +147,7 @@ IterativeQpeWorkflow::constructQpeCircuit(const QuantumSimulationModel &model, 
  TrotterEvolution method;

  const double trotterStepSize = -2 * M_PI / num_steps;

  auto trotterCir =

      method.create_ansatz(model.observable, {{"dt", trotterStepSize}}).circuit;

      method.create_ansatz(obs.get(), {{"dt", trotterStepSize}}).circuit;

  // std::cout << "Trotter circ:\n" << trotterCir->toString() << "\n";

  // Controlled-U

  auto ctrlKernel = std::dynamic_pointer_cast<CompositeInstruction>(
@@ -174,9 +171,8 @@ IterativeQpeWorkflow::constructQpeCircuit(const QuantumSimulationModel &model, 


  // Rz on ancilla qubit

  // Global phase due to identity pauli

  if (model.observable->getIdentitySubTerm()) {

    const double idCoeff =

        model.observable->getIdentitySubTerm()->coefficient().real();

  if (obs->getIdentitySubTerm()) {

    const double idCoeff = obs->getIdentitySubTerm()->coefficient().real();

    const double globalPhase = 2 * M_PI * idCoeff * power;

    // std::cout << "Global phase = " << globalPhase << "\n";

    kernel->addInstruction(
@@ -195,8 +191,38 @@ IterativeQpeWorkflow::constructQpeCircuit(const QuantumSimulationModel &model, 
  return kernel;

}



void IterativeQpeWorkflow::HamOpConverter::fromObservable(Observable *obs) {

  translation = 0.0;



  for (auto &term : obs->getSubTerms()) {

    translation += std::abs(term->coefficient());

  }

  

  stretch = 0.5 / translation;

}



std::shared_ptr<Observable>

IterativeQpeWorkflow::HamOpConverter::stretchObservable(Observable *obs) const {

  PauliOperator *pauliCast = static_cast<PauliOperator *>(obs);

  if (pauliCast) {

    auto result = std::make_shared<PauliOperator>(translation);

    result->operator+=(*pauliCast);

    result->operator*=(stretch);

    return result;

  } else {

    return nullptr;

  }

}



double

IterativeQpeWorkflow::HamOpConverter::computeEnergy(double phaseVal) const {

  return phaseVal / stretch - translation;

}



QuantumSimulationResult

IterativeQpeWorkflow::execute(const QuantumSimulationModel &model) {

  ham_converter.fromObservable(model.observable);

  auto stretchedObs = ham_converter.stretchObservable(model.observable);

  auto provider = xacc::getIRProvider("quantum");

  // Iterative Quantum Phase Estimation:

  // We're using XACC IR construction API here, since using QCOR kernels here
@@ -214,11 +240,12 @@ IterativeQpeWorkflow::execute(const QuantumSimulationModel &model) { 
    omega_coef = omega_coef/2.0;

    // Construct the QPE circuit and append to the kernel:

    auto k = num_iters - iterIdx;

    auto iterQpe = constructQpeCircuit(model, k, -2 * M_PI * omega_coef);

    

    auto iterQpe = constructQpeCircuit(stretchedObs, k, -2 * M_PI * omega_coef);

    kernel->addInstruction(iterQpe);

    // Executes the iterative QPE algorithm:

    auto qpu = xacc::internal_compiler::get_qpu();

    auto temp_buffer = xacc::qalloc(model.observable->nBits() + 1);

    auto temp_buffer = xacc::qalloc(stretchedObs->nBits() + 1);

    // std::cout << "Kernel: \n" << kernel->toString() << "\n";



    qpu->execute(temp_buffer, kernel);
@@ -249,8 +276,8 @@ IterativeQpeWorkflow::execute(const QuantumSimulationModel &model) { 
    // std::cout << "Iter " << iterIdx << ": Result = " << bitResult << "; omega_coef = " << omega_coef << "\n";



  }

  // std::cout << "Final phase = " << omega_coef << "\n";

  return { {"phase", omega_coef}};

  

  return { {"phase", omega_coef}, {"energy", ham_converter.computeEnergy(omega_coef)}};

}



std::shared_ptr<QuantumSimulationWorkflow>

lib/qsim/impls/qsim_impl.hpp

+11 −1
Original line number Diff line number Diff line
@@ -77,6 +77,15 @@ private: 
// This can be integrated as a CostFuncEvaluator if needed.

class IterativeQpeWorkflow : public QuantumSimulationWorkflow {

public:

  // Translate/stretch the Hamiltonian operator for QPE.

  struct HamOpConverter {

    double translation;

    double stretch;

    void fromObservable(Observable *obs);

    std::shared_ptr<Observable> stretchObservable(Observable *obs) const;

    double computeEnergy(double phaseVal) const;

  };



  virtual bool initialize(const HeterogeneousMap &params) override;

  virtual QuantumSimulationResult

  execute(const QuantumSimulationModel &model) override;
@@ -85,7 +94,7 @@ public: 


private:

  std::shared_ptr<CompositeInstruction>

  constructQpeCircuit(const QuantumSimulationModel &model, int k, double omega,

  constructQpeCircuit(std::shared_ptr<Observable> obs, int k, double omega,

                      bool measure = true) const;



private:
@@ -93,6 +102,7 @@ private: 
  int num_steps;

  // Number of iterations (>=1)

  int num_iters;

  HamOpConverter ham_converter;

};



class DefaultObjFuncEval : public CostFunctionEvaluator {

lib/qsim/impls/tests/IterativeQpeWorkflow.cpp

+7 −4
Original line number Diff line number Diff line
@@ -9,14 +9,13 @@ 
// Ansatz to bring the state into an eigenvector state of the Hamiltonian.

__qpu__ void eigen_state_prep(qreg q) {

  using qcor::openqasm;

  u3(pi / 4, 0, 0) q[0];

  u3(0.95531662,0,0) q[0];

  u1(pi/4) q[0];

}



int main(int argc, char **argv) {

  // Create Hamiltonian:

  // Important: this must be a weighted pauli operator list.

  // TODO: add function to convert arbitrary Hamiltonian to the weighted form.

  auto H = 0.5 + 0.25 * X(0) + 0.25 * Z(0);

  auto H = 1.0 + X(0) + Y(0) +  Z(0);

  auto problemModel =

      qsim::ModelBuilder::createModel(eigen_state_prep, H, 1, 0);
@@ -26,6 +25,10 @@ int main(int argc, char **argv) { 


  auto result = workflow->execute(problemModel);

  const double phaseValue = result.get<double>("phase");

  const double energy = result.get<double>("energy");



  std::cout << "Final phase = " << phaseValue << "\n";

  std::cout << "Energy = " << energy << "\n";



  return 0;

}
 
 No newline at end of file