diff --git a/projectq/backends/_qracksim/_cpp/qracksimulator.hpp b/projectq/backends/_qracksim/_cpp/qracksimulator.hpp index 8608271e5..0e9acf96d 100755 --- a/projectq/backends/_qracksim/_cpp/qracksimulator.hpp +++ b/projectq/backends/_qracksim/_cpp/qracksimulator.hpp @@ -40,6 +40,10 @@ class QrackSimulator{ using StateVector = std::vector, aligned_allocator,64>>; using Map = std::map; using RndEngine = qrack_rand_gen; + using Term = std::vector>; + using TermsDict = std::vector>; + using ComplexTermsDict = std::vector>>; + using Matrix = std::vector, aligned_allocator, 64>>>; enum Qrack::QInterfaceEngine QrackEngine = Qrack::QINTERFACE_QUNIT; enum Qrack::QInterfaceEngine QrackSubengine1 = Qrack::QINTERFACE_QFUSION; #if ENABLE_OPENCL @@ -422,6 +426,31 @@ class QrackSimulator{ delete[] substateVec; } + void apply_qubit_operator(ComplexTermsDict const& td, std::vector const& ids){ + for (auto const& term : td){ + apply_term(term.first, term.second, ids, {}); + } + } + + calc_type get_expectation_value(TermsDict const& td, std::vector const& ids){ + calc_type expectation = 0; + + std::size_t mask = 0; + for (unsigned i = 0; i < ids.size(); i++){ + mask |= 1UL << map_[ids[i]]; + } + + run(); + + Qrack::QInterfacePtr qRegOrig = qReg->Clone(); + for (auto const& term : td){ + expectation += diagonalize(term.first, term.second, ids); + qReg = qRegOrig->Clone(); + } + + return expectation; + } + std::tuple cheat(){ if (qReg == NULL) { StateVector vec(1, 0.0); @@ -556,6 +585,51 @@ class QrackSimulator{ } } + void apply_term(Term const& term, std::complex coeff, std::vector const& ids, + std::vector const& ctrl){ + + std::complex I(0., 1.); + Matrix X = {{0., 1.}, {1., 0.}}; + Matrix Y = {{0., -I}, {I, 0.}}; + Matrix Z = {{1., 0.}, {0., -1.}}; + std::vector gates = {X, Y, Z}; + + for (auto const& local_op : term){ + unsigned id = ids[local_op.first]; + auto temp = gates[local_op.second - 'X']; + for (unsigned i = 0; i < 4; i++) { + temp[i / 2][i % 2] *= -coeff; + } + apply_controlled_gate(temp, {id}, ctrl); + } + } + + calc_type diagonalize(Term const& term, std::complex coeff, std::vector const& ids){ + calc_type expectation = 1; + calc_type angle = arg(coeff); + calc_type len = abs(coeff); + std::complex phaseFac(cos(angle), sin(angle)); + bitCapInt idPower; + + std::complex I(0., 1.); + Matrix X = {{M_SQRT1_2, M_SQRT1_2}, {M_SQRT1_2, -M_SQRT1_2}}; + Matrix Y = {{M_SQRT1_2, -M_SQRT1_2 * I}, {M_SQRT1_2, M_SQRT1_2 * I}}; + Matrix Z = {{1., 0.}, {0., 1.}}; + std::vector gates = {X, Y, Z}; + + for (auto const& local_op : term){ + unsigned id = ids[local_op.first]; + auto temp = gates[local_op.second - 'X']; + for (unsigned i = 0; i < 4; i++) { + temp[i / 2][i % 2] *= phaseFac; + } + apply_controlled_gate(temp, {id}, {}); + idPower = 1U << map_[id]; + expectation *= (qReg->ProbMask(idPower, 0) - qReg->ProbMask(idPower, idPower)); + } + return len * expectation; + } + Map map_; std::shared_ptr rnd_eng_; Qrack::QInterfacePtr qReg; diff --git a/projectq/backends/_qracksim/_qracksim.cpp b/projectq/backends/_qracksim/_qracksim.cpp index 5f16d09e6..89426da65 100755 --- a/projectq/backends/_qracksim/_qracksim.cpp +++ b/projectq/backends/_qracksim/_qracksim.cpp @@ -51,6 +51,8 @@ PYBIND11_PLUGIN(_qracksim) { .def("apply_controlled_dec", &QrackSimulator::apply_controlled_dec) .def("apply_controlled_mul", &QrackSimulator::apply_controlled_mul) .def("apply_controlled_div", &QrackSimulator::apply_controlled_div) + .def("get_expectation_value", &QrackSimulator::get_expectation_value) + .def("apply_qubit_operator", &QrackSimulator::apply_qubit_operator) .def("get_probability", &QrackSimulator::get_probability) .def("get_amplitude", &QrackSimulator::get_amplitude) .def("set_wavefunction", &QrackSimulator::set_wavefunction) diff --git a/projectq/backends/_qracksim/_simulator.py b/projectq/backends/_qracksim/_simulator.py index 0324f8bbc..dd821e0d1 100755 --- a/projectq/backends/_qracksim/_simulator.py +++ b/projectq/backends/_qracksim/_simulator.py @@ -1,4 +1,3 @@ -# Copyright 2017 ProjectQ-Framework (www.projectq.ch) # # Licensed under the Apache License, Version 2.0 (the "License"); # you may not use this file except in compliance with the License. @@ -31,7 +30,8 @@ Deallocate, UniformlyControlledRy, UniformlyControlledRz, - StatePreparation) + StatePreparation, + QubitOperator) from projectq.libs.math import (AddConstant, AddConstantModN, MultiplyByConstantModN) @@ -119,6 +119,8 @@ def is_available(self, cmd): if (isinstance(cmd.gate, StatePreparation) and not cmd.control_qubits): # Qrack has inexpensive ways of preparing a partial state, without controls. return True + elif (isinstance(cmd.gate, QubitOperator) and not np.isclose(1, np.absolute(cmd.gate.coefficient))): + return True except: pass @@ -154,12 +156,81 @@ def _convert_logical_to_mapped_qureg(self, qureg): return qureg def get_expectation_value(self, qubit_operator, qureg): - # To maintain compatibility with default Simulator, for the moment. - pass + """ + Get the expectation value of qubit_operator w.r.t. the current wave + function represented by the supplied quantum register. + + Args: + qubit_operator (projectq.ops.QubitOperator): Operator to measure. + qureg (list[Qubit],Qureg): Quantum bits to measure. + + Returns: + Expectation value + + Note: + Make sure all previous commands (especially allocations) have + passed through the compilation chain (call main_engine.flush() to + make sure). + + Note: + If there is a mapper present in the compiler, this function + automatically converts from logical qubits to mapped qubits for + the qureg argument. + + Raises: + Exception: If `qubit_operator` acts on more qubits than present in + the `qureg` argument. + """ + qureg = self._convert_logical_to_mapped_qureg(qureg) + num_qubits = len(qureg) + for term, _ in qubit_operator.terms.items(): + if not term == () and term[-1][0] >= num_qubits: + raise Exception("qubit_operator acts on more qubits than " + "contained in the qureg.") + operator = [(list(term), coeff) for (term, coeff) + in qubit_operator.terms.items()] + return self._simulator.get_expectation_value(operator, + [qb.id for qb in qureg]) def apply_qubit_operator(self, qubit_operator, qureg): - # To maintain compatibility with default Simulator, for the moment. - pass + """ + Apply a (possibly non-unitary) qubit_operator to the current wave + function represented by the supplied quantum register. + + Args: + qubit_operator (projectq.ops.QubitOperator): Operator to apply. + qureg (list[Qubit],Qureg): Quantum bits to which to apply the + operator. + + Raises: + Exception: If `qubit_operator` acts on more qubits than present in + the `qureg` argument. + + Warning: + This function allows applying non-unitary gates and it will not + re-normalize the wave function! It is for numerical experiments + only and should not be used for other purposes. + + Note: + Make sure all previous commands (especially allocations) have + passed through the compilation chain (call main_engine.flush() to + make sure). + + Note: + If there is a mapper present in the compiler, this function + automatically converts from logical qubits to mapped qubits for + the qureg argument. + """ + qureg = self._convert_logical_to_mapped_qureg(qureg) + num_qubits = len(qureg) + for term, _ in qubit_operator.terms.items(): + if not term == () and term[-1][0] >= num_qubits: + raise Exception("qubit_operator acts on more qubits than " + "contained in the qureg.") + operator = [(list(term), coeff) for (term, coeff) + in qubit_operator.terms.items()] + return self._simulator.apply_qubit_operator(operator, + [qb.id for qb in qureg]) def get_probability(self, bit_string, qureg): """ @@ -357,6 +428,9 @@ def _handle(self, cmd): self._simulator.apply_controlled_sqrtswap(ids1, ids2, [qb.id for qb in cmd.control_qubits]) + elif isinstance(cmd.gate, QubitOperator): + ids = [qb.id for qb in cmd.qubits[0]] + self.apply_qubit_operator(cmd.gate, ids) elif isinstance(cmd.gate, AddConstant) or isinstance(cmd.gate, AddConstantModN): #Unless there's a carry, the only unitary addition is mod (2^len(ids)) ids = [qb.id for qr in cmd.qubits for qb in qr] diff --git a/projectq/backends/_qracksim/_simulator_test.py b/projectq/backends/_qracksim/_simulator_test.py index 0cba59575..b907b6be8 100755 --- a/projectq/backends/_qracksim/_simulator_test.py +++ b/projectq/backends/_qracksim/_simulator_test.py @@ -45,6 +45,9 @@ from projectq.backends import Simulator +tolerance = 1e-6 + + def test_is_qrack_simulator_present(): _qracksim = pytest.importorskip("projectq.backends._qracksim._qracksim") import projectq.backends._qracksim._qracksim as _ @@ -640,3 +643,64 @@ def test_uniformly_controlled_r(sim, gate_classes): All(Measure) | correct_qb + correct_ctrl_qureg test_eng.flush(deallocate_qubits=True) correct_eng.flush(deallocate_qubits=True) + +def test_get_expectation_value(sim): + num_qubits = 2 + test_eng = MainEngine(sim) + test_qureg = test_eng.allocate_qureg(num_qubits) + test_eng.flush() + + qubit_op = QubitOperator("X0", 1) + test_eng.backend.set_wavefunction([1 / math.sqrt(2), 1 / math.sqrt(2), 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(1, rel=tolerance, abs=tolerance)) + test_eng.backend.set_wavefunction([1 / math.sqrt(2), -1 / math.sqrt(2), 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(-1, rel=tolerance, abs=tolerance)) + + qubit_op = QubitOperator("Y0", 1) + test_eng.backend.set_wavefunction([1 / math.sqrt(2), 1j / math.sqrt(2), 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(1, rel=tolerance, abs=tolerance)) + test_eng.backend.set_wavefunction([1 / math.sqrt(2), -1j / math.sqrt(2), 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(-1, rel=tolerance, abs=tolerance)) + + qubit_op = QubitOperator("Z0", 1) + test_eng.backend.set_wavefunction([1, 0, 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(1, rel=tolerance, abs=tolerance)) + test_eng.backend.set_wavefunction([0, 1, 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(-1, rel=tolerance, abs=tolerance)) + + qubit_op = QubitOperator("Z0", 1 / 4) + test_eng.backend.set_wavefunction([1, 0, 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(1 / 4, rel=tolerance, abs=tolerance)) + test_eng.backend.set_wavefunction([0, 1, 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(-1 / 4, rel=tolerance, abs=tolerance)) + + qubit_op = QubitOperator("Z0 Z1", 1) + test_eng.backend.set_wavefunction([1, 0, 0, 0], + test_qureg) + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(1, rel=tolerance, abs=tolerance)) + X | test_qureg[0] + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(-1, rel=tolerance, abs=tolerance)) + X | test_qureg[1] + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(1, rel=tolerance, abs=tolerance)) + X | test_qureg[0] + test_eng.flush() + assert(sim.get_expectation_value(qubit_op, test_qureg) == pytest.approx(-1, rel=tolerance, abs=tolerance)) diff --git a/projectq/setups/decompositions/qubitop2onequbit_test.py b/projectq/setups/decompositions/qubitop2onequbit_test.py index eddf4087b..1e8bb85f4 100644 --- a/projectq/setups/decompositions/qubitop2onequbit_test.py +++ b/projectq/setups/decompositions/qubitop2onequbit_test.py @@ -18,7 +18,7 @@ from projectq import MainEngine # Qrack simulator does not yet support QubitOperator, so always use the default simulator: -from projectq.backends._sim import Simulator +from projectq.backends import Simulator from projectq.cengines import (AutoReplacer, DecompositionRuleSet, DummyEngine, InstructionFilter) from projectq.meta import Control