diff --git a/CMakeLists.txt b/CMakeLists.txt index 29ed6c252..b89ca98b2 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -302,7 +302,6 @@ install (FILES include/qbdt.hpp include/qbdt_node.hpp include/qbdt_node_interface.hpp - include/qbdt_qengine_node.hpp DESTINATION ${CMAKE_INSTALL_INCLUDEDIR}/qrack ) diff --git a/cmake/Qbdt.cmake b/cmake/Qbdt.cmake index 9c23aeb08..fa38c8222 100644 --- a/cmake/Qbdt.cmake +++ b/cmake/Qbdt.cmake @@ -4,7 +4,6 @@ if (ENABLE_QBDT) target_sources (qrack PRIVATE src/qbdt/node_interface.cpp src/qbdt/node.cpp - src/qbdt/qengine_node.cpp src/qbdt/tree.cpp ) endif (ENABLE_QBDT) diff --git a/include/qbdt.hpp b/include/qbdt.hpp index 0c024080b..532bab628 100644 --- a/include/qbdt.hpp +++ b/include/qbdt.hpp @@ -16,10 +16,9 @@ #pragma once -#include "qbdt_qengine_node.hpp" +#include "qbdt_node.hpp" #include "qengine.hpp" -#define NODE_TO_QENGINE(leaf) (std::dynamic_pointer_cast(leaf)->qReg) #define QINTERFACE_TO_QALU(qReg) std::dynamic_pointer_cast(qReg) #define QINTERFACE_TO_QPARITY(qReg) std::dynamic_pointer_cast(qReg) @@ -34,9 +33,6 @@ class QBdt : public QAlu, public QParity, public QInterface { class QBdt : public QParity, public QInterface { #endif protected: - bitLenInt attachedQubitCount; - bitLenInt bdtQubitCount; - bitLenInt maxPageQubits; bitLenInt bdtStride; int64_t devID; QBdtNodeInterfacePtr root; @@ -44,54 +40,61 @@ class QBdt : public QParity, public QInterface { std::vector deviceIDs; std::vector engines; - void SetQubitCount(bitLenInt qb, bitLenInt aqb) - { - attachedQubitCount = aqb; - SetQubitCount(qb); - } + QEnginePtr MakeQEngine(bitLenInt qbCount, bitCapInt perm = 0U); - void SetQubitCount(bitLenInt qb) + template void GetTraversal(Fn getLambda) { - QInterface::SetQubitCount(qb); - bdtQubitCount = qubitCount - attachedQubitCount; - bdtMaxQPower = pow2(bdtQubitCount); + Finish(); + + for (bitCapInt i = 0U; i < maxQPower; ++i) { + QBdtNodeInterfacePtr leaf = root; + complex scale = leaf->scale; + for (bitLenInt j = 0U; j < qubitCount; ++j) { + if (norm(leaf->scale) <= _qrack_qbdt_sep_thresh) { + break; + } + leaf = leaf->branches[SelectBit(i, j)]; + scale *= leaf->scale; + } + + getLambda((bitCapIntOcl)i, scale); + } } + template void SetTraversal(Fn setLambda) + { + Dump(); - QBdtQEngineNodePtr MakeQEngineNode(complex scale, bitLenInt qbCount, bitCapInt perm = 0U); + root = std::make_shared(); + root->Branch(qubitCount); - QInterfacePtr MakeTempStateVector() - { - QInterfacePtr copyPtr = NODE_TO_QENGINE(MakeQEngineNode(ONE_R1, qubitCount)); - Finish(); - GetQuantumState(copyPtr); + _par_for(maxQPower, [&](const bitCapInt& i, const unsigned& cpu) { + QBdtNodeInterfacePtr prevLeaf = root; + QBdtNodeInterfacePtr leaf = root; + for (bitLenInt j = 0U; j < qubitCount; ++j) { + prevLeaf = leaf; + leaf = leaf->branches[SelectBit(i, j)]; + } - // If the calling function fully deferences our return, it's automatically freed. - return copyPtr; - } + setLambda((bitCapIntOcl)i, leaf); + }); - template void GetTraversal(Fn getLambda); - template void SetTraversal(Fn setLambda); + root->PopStateVector(qubitCount); + root->Prune(qubitCount); + } template void ExecuteAsStateVector(Fn operation) { - if (!bdtQubitCount) { - operation(NODE_TO_QENGINE(root)); - return; - } - - SetStateVector(); - operation(NODE_TO_QENGINE(root)); - ResetStateVector(); + QInterfacePtr qReg = MakeQEngine(qubitCount); + GetQuantumState(qReg); + operation(qReg); + SetQuantumState(qReg); } template bitCapInt BitCapIntAsStateVector(Fn operation) { - if (!bdtQubitCount) { - return operation(NODE_TO_QENGINE(root)); - } - - SetStateVector(); - bitCapInt toRet = operation(NODE_TO_QENGINE(root)); - ResetStateVector(); + QInterfacePtr qReg = MakeQEngine(qubitCount); + GetQuantumState(qReg); + const bitCapInt toRet = operation(qReg); + SetQuantumState(qReg); return toRet; } @@ -133,29 +136,9 @@ class QBdt : public QParity, public QInterface { { } - QBdt(QEnginePtr enginePtr, std::vector eng, bitLenInt qBitCount, bitCapInt ignored = 0U, - qrack_rand_gen_ptr rgp = nullptr, complex phaseFac = CMPLX_DEFAULT_ARG, bool doNorm = false, - bool randomGlobalPhase = true, bool useHostMem = false, int64_t deviceId = -1, bool useHardwareRNG = true, - bool useSparseStateVec = false, real1_f norm_thresh = REAL1_EPSILON, std::vector devList = {}, - bitLenInt qubitThreshold = 0U, real1_f separation_thresh = FP_NORM_EPSILON_F); - - QEnginePtr ReleaseEngine() - { - if (bdtQubitCount) { - throw std::domain_error("Cannot release QEngine from QBdt with BDT qubits!"); - } - - return NODE_TO_QENGINE(root); - } - - void LockEngine(QEnginePtr eng) { root = std::make_shared(ONE_CMPLX, eng); } - bool isBinaryDecisionTree() { return true; }; - void SetStateVector(); - void ResetStateVector(bitLenInt aqb = 0U); - - void SetDevice(int64_t dID); + void SetDevice(int64_t dID) { devID = dID; } void UpdateRunningNorm(real1_f norm_thresh = REAL1_DEFAULT_ARG) { @@ -165,7 +148,7 @@ class QBdt : public QParity, public QInterface { void NormalizeState( real1_f nrm = REAL1_DEFAULT_ARG, real1_f norm_thresh = REAL1_DEFAULT_ARG, real1_f phaseArg = ZERO_R1_F) { - root->Normalize(bdtQubitCount); + root->Normalize(qubitCount); } real1_f SumSqrDiff(QInterfacePtr toCompare) { return SumSqrDiff(std::dynamic_pointer_cast(toCompare)); } @@ -175,11 +158,26 @@ class QBdt : public QParity, public QInterface { QInterfacePtr Clone(); - void GetQuantumState(complex* state); - void GetQuantumState(QInterfacePtr eng); - void SetQuantumState(const complex* state); - void SetQuantumState(QInterfacePtr eng); - void GetProbs(real1* outputProbs); + void GetQuantumState(complex* state) + { + GetTraversal([state](bitCapIntOcl i, complex scale) { state[i] = scale; }); + } + void GetQuantumState(QInterfacePtr eng) + { + GetTraversal([eng](bitCapIntOcl i, complex scale) { eng->SetAmplitude(i, scale); }); + } + void SetQuantumState(const complex* state) + { + SetTraversal([state](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { leaf->scale = state[i]; }); + } + void SetQuantumState(QInterfacePtr eng) + { + SetTraversal([eng](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { leaf->scale = eng->GetAmplitude(i); }); + } + void GetProbs(real1* outputProbs) + { + GetTraversal([outputProbs](bitCapIntOcl i, complex scale) { outputProbs[i] = norm(scale); }); + } complex GetAmplitude(bitCapInt perm); void SetAmplitude(bitCapInt perm, complex amp) @@ -196,41 +194,14 @@ class QBdt : public QParity, public QInterface { void Decompose(bitLenInt start, QInterfacePtr dest) { QBdtPtr d = std::dynamic_pointer_cast(dest); - if (!bdtQubitCount) { - d->root = d->MakeQEngineNode(ONE_CMPLX, d->qubitCount, 0U); - NODE_TO_QENGINE(root)->Decompose(start, NODE_TO_QENGINE(d->root)); - d->SetQubitCount(d->qubitCount, d->qubitCount); - SetQubitCount(qubitCount - d->qubitCount, qubitCount - d->qubitCount); - - return; - } - DecomposeDispose(start, dest->GetQubitCount(), d); } QInterfacePtr Decompose(bitLenInt start, bitLenInt length); - void Dispose(bitLenInt start, bitLenInt length) - { - if (!bdtQubitCount) { - NODE_TO_QENGINE(root)->Dispose(start, length); - SetQubitCount(qubitCount - length, qubitCount - length); - - return; - } - - DecomposeDispose(start, length, NULL); - } + void Dispose(bitLenInt start, bitLenInt length) { DecomposeDispose(start, length, NULL); } void Dispose(bitLenInt start, bitLenInt length, bitCapInt disposedPerm) { - if (!bdtQubitCount) { - NODE_TO_QENGINE(root)->Dispose(start, length, disposedPerm); - SetQubitCount(qubitCount - length, qubitCount - length); - - return; - } - ForceMReg(start, length, disposedPerm); - DecomposeDispose(start, length, NULL); } @@ -262,8 +233,8 @@ class QBdt : public QParity, public QInterface { } real1_f toRet; - ExecuteAsStateVector( - [&](QInterfacePtr eng) { toRet = QINTERFACE_TO_QPARITY(NODE_TO_QENGINE(root))->ProbParity(mask); }); + ExecuteAsStateVector([&](QInterfacePtr eng) { toRet = QINTERFACE_TO_QPARITY(eng)->ProbParity(mask); }); + return toRet; } void CUniformParityRZ(const std::vector& controls, bitCapInt mask, real1_f angle) @@ -283,8 +254,11 @@ class QBdt : public QParity, public QInterface { return ForceM(log2(mask), result, doForce); } - SetStateVector(); - return QINTERFACE_TO_QPARITY(NODE_TO_QENGINE(root))->ForceMParity(mask, result, doForce); + bool toRet; + ExecuteAsStateVector( + [&](QInterfacePtr eng) { toRet = QINTERFACE_TO_QPARITY(eng)->ForceMParity(mask, result, doForce); }); + + return toRet; } #if ENABLE_ALU diff --git a/include/qbdt_qengine_node.hpp b/include/qbdt_qengine_node.hpp deleted file mode 100644 index 1044eb70d..000000000 --- a/include/qbdt_qengine_node.hpp +++ /dev/null @@ -1,105 +0,0 @@ -////////////////////////////////////////////////////////////////////////////////////// -// -// (C) Daniel Strano and the Qrack contributors 2017-2021. All rights reserved. -// -// QBinaryDecision tree is an alternative approach to quantum state representation, as -// opposed to state vector representation. This is a compressed form that can be -// operated directly on while compressed. Inspiration for the Qrack implementation was -// taken from JKQ DDSIM, maintained by the Institute for Integrated Circuits at the -// Johannes Kepler University Linz: -// -// https://github.com/iic-jku/ddsim -// -// Licensed under the GNU Lesser General Public License V3. -// See LICENSE.md in the project root or https://www.gnu.org/licenses/lgpl-3.0.en.html -// for details. - -#pragma once - -#include "qbdt_node_interface.hpp" -#include "qengine.hpp" - -namespace Qrack { - -class QBdtQEngineNode; -typedef std::shared_ptr QBdtQEngineNodePtr; - -class QBdtQEngineNode : public QBdtNodeInterface { -protected: -#if ENABLE_COMPLEX_X2 - virtual void PushStateVector(const complex2& mtrxCol1, const complex2& mtrxCol2, const complex2& mtrxColShuff1, - const complex2& mtrxColShuff2, QBdtNodeInterfacePtr& b0, QBdtNodeInterfacePtr& b1, bitLenInt depth, - bitLenInt parDepth = 1U) -#else - virtual void PushStateVector(complex const* mtrx, QBdtNodeInterfacePtr& b0, QBdtNodeInterfacePtr& b1, - bitLenInt depth, bitLenInt parDepth = 1U) -#endif - { - throw std::out_of_range("QBdtQEngineNode::PushStateVector() not implemented!"); - } - -public: - QEnginePtr qReg; - - QBdtQEngineNode() - : QBdtNodeInterface(ZERO_CMPLX) - , qReg(NULL) - { - // Intentionally left blank. - } - - QBdtQEngineNode(complex scl, QEnginePtr q) - : QBdtNodeInterface(scl) - , qReg(q) - { - // Intentionally left blank - } - - virtual ~QBdtQEngineNode() - { - // Virtual destructor for inheritance - } - - virtual void SetZero() - { - QBdtNodeInterface::SetZero(); - qReg = NULL; - } - - virtual QBdtNodeInterfacePtr ShallowClone() { return std::make_shared(scale, qReg); } - - virtual bool isEqual(QBdtNodeInterfacePtr r); - - virtual bool isEqualUnder(QBdtNodeInterfacePtr r); - - virtual void Normalize(bitLenInt depth = 1U); - - virtual void Branch(bitLenInt depth = 1U, bitLenInt parDepth = 1U); - - virtual void Prune(bitLenInt depth = 1U, bitLenInt parDepth = 1U); - - virtual void InsertAtDepth(QBdtNodeInterfacePtr b, bitLenInt depth, const bitLenInt& size, bitLenInt parDepth = 1U); - - virtual QBdtNodeInterfacePtr RemoveSeparableAtDepth( - bitLenInt depth, const bitLenInt& size, bitLenInt parDepth = 1U); - - virtual void PopStateVector(bitLenInt depth = 1U, bitLenInt parDepth = 1U); - -#if ENABLE_COMPLEX_X2 - virtual void Apply2x2(const complex2& mtrxCol1, const complex2& mtrxCol2, const complex2& mtrxColShuff1, - const complex2& mtrxColShuff2, bitLenInt depth) -#else - virtual void Apply2x2(complex const* mtrx, bitLenInt depth) -#endif - { - throw std::out_of_range("QBdtQEngineNode::Apply2x2() not implemented!"); - } -#if ENABLE_COMPLEX_X2 - virtual void PushSpecial(const complex2& mtrxCol1, const complex2& mtrxCol2, const complex2& mtrxColShuff1, - const complex2& mtrxColShuff2, QBdtNodeInterfacePtr& b1); -#else - virtual void PushSpecial(complex const* mtrx, QBdtNodeInterfacePtr& b1); -#endif -}; - -} // namespace Qrack diff --git a/src/qbdt/qengine_node.cpp b/src/qbdt/qengine_node.cpp deleted file mode 100644 index fb188fe80..000000000 --- a/src/qbdt/qengine_node.cpp +++ /dev/null @@ -1,224 +0,0 @@ -////////////////////////////////////////////////////////////////////////////////////// -// -// (C) Daniel Strano and the Qrack contributors 2017-2021. All rights reserved. -// -// QBinaryDecision tree is an alternative approach to quantum state representation, as -// opposed to state vector representation. This is a compressed form that can be -// operated directly on while compressed. Inspiration for the Qrack implementation was -// taken from JKQ DDSIM, maintained by the Institute for Integrated Circuits at the -// Johannes Kepler University Linz: -// -// https://github.com/iic-jku/ddsim -// -// Licensed under the GNU Lesser General Public License V3. -// See LICENSE.md in the project root or https://www.gnu.org/licenses/lgpl-3.0.en.html -// for details. - -#include "qbdt_qengine_node.hpp" - -#define IS_NODE_0(c) (norm(c) <= _qrack_qbdt_sep_thresh) -#define IS_SAME_AMP(a, b) (abs((a) - (b)) <= REAL1_EPSILON) - -namespace Qrack { -bool QBdtQEngineNode::isEqual(QBdtNodeInterfacePtr r) -{ - if (!r) { - return false; - } - - if (this == r.get()) { - return true; - } - - if (!IS_SAME_AMP(scale, r->scale)) { - return false; - } - - if (IS_NODE_0(scale)) { - return true; - } - - QEnginePtr rReg = std::dynamic_pointer_cast(r)->qReg; - - if (qReg.get() == rReg.get()) { - return true; - } - - if (qReg->ApproxCompare(rReg)) { - qReg = rReg; - return true; - } - - return false; -} - -bool QBdtQEngineNode::isEqualUnder(QBdtNodeInterfacePtr r) -{ - if (!r) { - return false; - } - - if (this == r.get()) { - return true; - } - - if (IS_NODE_0(scale)) { - return IS_NODE_0(r->scale); - } - - QEnginePtr rReg = std::dynamic_pointer_cast(r)->qReg; - - if (qReg.get() == rReg.get()) { - return true; - } - - if (qReg->ApproxCompare(rReg)) { - qReg = rReg; - return true; - } - - return false; -} - -void QBdtQEngineNode::Normalize(bitLenInt depth) -{ - if (!depth) { - return; - } - - if (IS_NODE_0(scale)) { - return; - } - - if (qReg) { - qReg->UpdateRunningNorm(); - qReg->NormalizeState(); - } -} - -void QBdtQEngineNode::Branch(bitLenInt depth, bitLenInt parDepth) -{ - if (!depth) { - return; - } - - if (IS_NODE_0(scale)) { - SetZero(); - return; - } - - if (qReg) { - qReg = std::dynamic_pointer_cast(qReg->Clone()); - } -} - -void QBdtQEngineNode::Prune(bitLenInt depth, bitLenInt unused) -{ - if (IS_NODE_0(scale)) { - SetZero(); - return; - } - - const real1_f phaseArg = qReg->FirstNonzeroPhase(); - qReg->NormalizeState(REAL1_DEFAULT_ARG, REAL1_DEFAULT_ARG, -phaseArg); - scale *= std::polar(ONE_R1, (real1)phaseArg); -} - -void QBdtQEngineNode::InsertAtDepth(QBdtNodeInterfacePtr b, bitLenInt depth, const bitLenInt& size, bitLenInt parDepth) -{ - if (IS_NODE_0(scale)) { - SetZero(); - return; - } - - QBdtQEngineNodePtr bEng = std::dynamic_pointer_cast(b); - qReg->Compose(bEng->qReg, depth); -} - -QBdtNodeInterfacePtr QBdtQEngineNode::RemoveSeparableAtDepth(bitLenInt depth, const bitLenInt& size, bitLenInt parDepth) -{ - if (!size) { - return NULL; - } - - if (IS_NODE_0(scale)) { - SetZero(); - return NULL; - } - - QBdtQEngineNodePtr toRet = std::dynamic_pointer_cast(ShallowClone()); - toRet->scale /= abs(toRet->scale); - - if (!qReg) { - return toRet; - } - - toRet->qReg = std::dynamic_pointer_cast(qReg->Decompose(depth, size)); - - return toRet; -} - -#if ENABLE_COMPLEX_X2 -void QBdtQEngineNode::PushSpecial(const complex2& mtrxCol1, const complex2& mtrxCol2, const complex2& mtrxColShuff1, - const complex2& mtrxColShuff2, QBdtNodeInterfacePtr& b1) -{ - const complex mtrx[4U]{ mtrxCol1.c(0U), mtrxCol2.c(0U), mtrxCol1.c(1U), mtrxCol2.c(1U) }; -#else -void QBdtQEngineNode::PushSpecial(complex const* mtrx, QBdtNodeInterfacePtr& b1) -{ -#endif - const bool is0Zero = IS_NODE_0(scale); - const bool is1Zero = IS_NODE_0(b1->scale); - - if (is0Zero && is1Zero) { - SetZero(); - b1->SetZero(); - - return; - } - - QEnginePtr qReg0 = qReg; - QEnginePtr qReg1 = std::dynamic_pointer_cast(b1)->qReg; - - if (is0Zero) { - qReg0 = std::dynamic_pointer_cast(b1)->qReg->CloneEmpty(); - } else if (is1Zero) { - qReg1 = qReg->CloneEmpty(); - } - - qReg0->NormalizeState((real1_f)(ONE_R1 / norm(scale)), REAL1_DEFAULT_ARG, (real1_f)std::arg(scale)); - qReg1->NormalizeState((real1_f)(ONE_R1 / norm(b1->scale)), REAL1_DEFAULT_ARG, (real1_f)std::arg(b1->scale)); - - scale = SQRT1_2_R1; - b1->scale = SQRT1_2_R1; - - qReg0->ShuffleBuffers(qReg1); - - qReg0->Mtrx(mtrx, qReg0->GetQubitCount() - 1U); - qReg1->Mtrx(mtrx, qReg1->GetQubitCount() - 1U); - - qReg0->ShuffleBuffers(qReg1); - - qReg = qReg0; - std::dynamic_pointer_cast(b1)->qReg = qReg1; -} - -void QBdtQEngineNode::PopStateVector(bitLenInt depth, bitLenInt unused) -{ - if (IS_NODE_0(scale)) { - SetZero(); - return; - } - - qReg->UpdateRunningNorm(); - const real1_f nrm = qReg->GetRunningNorm(); - - if (nrm <= _qrack_qbdt_sep_thresh) { - SetZero(); - return; - } - - qReg->NormalizeState(); - scale = std::polar((real1)sqrt(nrm), ZERO_R1); -} -} // namespace Qrack diff --git a/src/qbdt/tree.cpp b/src/qbdt/tree.cpp index d0d0f0a21..13db7701a 100644 --- a/src/qbdt/tree.cpp +++ b/src/qbdt/tree.cpp @@ -14,7 +14,7 @@ // See LICENSE.md in the project root or https://www.gnu.org/licenses/lgpl-3.0.en.html // for details. -#include "qbdt_node.hpp" +#include "qbdt.hpp" #include "qfactory.hpp" #define IS_NODE_0(c) (norm(c) <= _qrack_qbdt_sep_thresh) @@ -26,7 +26,6 @@ QBdt::QBdt(std::vector eng, bitLenInt qBitCount, bitCapInt ini bool useSparseStateVec, real1_f norm_thresh, std::vector devIds, bitLenInt qubitThreshold, real1_f sep_thresh) : QInterface(qBitCount, rgp, doNorm, useHardwareRNG, randomGlobalPhase, doNorm ? norm_thresh : ZERO_R1_F) - , maxPageQubits(0U) , devID(deviceId) , root(NULL) , deviceIDs(devIds) @@ -34,29 +33,9 @@ QBdt::QBdt(std::vector eng, bitLenInt qBitCount, bitCapInt ini { Init(); - SetQubitCount(qBitCount, (maxPageQubits < qBitCount) ? maxPageQubits : qBitCount); - SetPermutation(initState); } -QBdt::QBdt(QEnginePtr enginePtr, std::vector eng, bitLenInt qBitCount, bitCapInt initState, - qrack_rand_gen_ptr rgp, complex phaseFac, bool doNorm, bool randomGlobalPhase, bool useHostMem, int64_t deviceId, - bool useHardwareRNG, bool useSparseStateVec, real1_f norm_thresh, std::vector devIds, - bitLenInt qubitThreshold, real1_f sep_thresh) - : QInterface(qBitCount, rgp, doNorm, useHardwareRNG, randomGlobalPhase, doNorm ? norm_thresh : ZERO_R1_F) - , maxPageQubits(0U) - , devID(deviceId) - , root(NULL) - , deviceIDs(devIds) - , engines(eng) -{ - Init(); - - SetQubitCount(qBitCount, qBitCount); - - LockEngine(enginePtr); -} - void QBdt::Init() { #if ENABLE_PTHREAD @@ -78,18 +57,12 @@ void QBdt::Init() ++engineLevel; rootEngine = engines[engineLevel]; } - - if (getenv("QRACK_QBDT_THRESHOLD_QB")) { - maxPageQubits = (bitLenInt)std::stoi(std::string(getenv("QRACK_QBDT_THRESHOLD_QB"))); - } } -QBdtQEngineNodePtr QBdt::MakeQEngineNode(complex scale, bitLenInt qbCount, bitCapInt perm) +QEnginePtr QBdt::MakeQEngine(bitLenInt qbCount, bitCapInt perm) { - return std::make_shared(scale, - std::dynamic_pointer_cast( - CreateQuantumInterface(engines, qbCount, perm, rand_generator, ONE_CMPLX, doNormalize, false, false, devID, - hardware_rand_generator != NULL, false, (real1_f)amplitudeFloor, deviceIDs))); + return std::dynamic_pointer_cast(CreateQuantumInterface(engines, qbCount, perm, rand_generator, ONE_CMPLX, + doNormalize, false, false, devID, hardware_rand_generator != NULL, false, (real1_f)amplitudeFloor, deviceIDs)); } void QBdt::par_for_qbdt(const bitCapInt& end, bitLenInt maxQubit, BdtFunc fn) @@ -226,27 +199,14 @@ void QBdt::SetPermutation(bitCapInt initState, complex phaseFac) } } - if (!bdtQubitCount) { - root = MakeQEngineNode(phaseFac, attachedQubitCount, initState); - - return; - } - - const bitLenInt maxQubit = attachedQubitCount ? (bdtQubitCount - 1U) : bdtQubitCount; root = std::make_shared(phaseFac); QBdtNodeInterfacePtr leaf = root; - for (bitLenInt qubit = 0U; qubit < maxQubit; ++qubit) { + for (bitLenInt qubit = 0U; qubit < qubitCount; ++qubit) { const size_t bit = SelectBit(initState, qubit); leaf->branches[bit] = std::make_shared(ONE_CMPLX); leaf->branches[bit ^ 1U] = std::make_shared(ZERO_CMPLX); leaf = leaf->branches[bit]; } - - if (attachedQubitCount) { - const size_t bit = SelectBit(initState, maxQubit); - leaf->branches[bit] = MakeQEngineNode(ONE_CMPLX, attachedQubitCount, initState >> bdtQubitCount); - leaf->branches[bit ^ 1U] = std::make_shared(); - } } QInterfacePtr QBdt::Clone() @@ -257,192 +217,11 @@ QInterfacePtr QBdt::Clone() false, -1, (hardware_rand_generator == NULL) ? false : true, false, (real1_f)amplitudeFloor); copyPtr->root = root ? root->ShallowClone() : NULL; - copyPtr->SetQubitCount(qubitCount, attachedQubitCount); - - if (!attachedQubitCount) { - return copyPtr; - } - - if (!bdtQubitCount) { - QBdtQEngineNodePtr eLeaf = std::dynamic_pointer_cast(copyPtr->root); - if (eLeaf->qReg) { - eLeaf->qReg = std::dynamic_pointer_cast(eLeaf->qReg->Clone()); - } - return copyPtr; - } - - std::map qis; - - copyPtr->SetTraversal([&qis](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { - QBdtQEngineNodePtr qenp = std::dynamic_pointer_cast(leaf); - QEnginePtr qi = NODE_TO_QENGINE(qenp); - if (qis.find(qi) == qis.end()) { - qis[qi] = std::dynamic_pointer_cast(qi->Clone()); - } - NODE_TO_QENGINE(qenp) = qis[qi]; - }); - copyPtr->root->Prune(bdtQubitCount); + copyPtr->SetQubitCount(qubitCount); return copyPtr; } -template void QBdt::GetTraversal(Fn getLambda) -{ - Finish(); - - for (bitCapInt i = 0U; i < maxQPower; ++i) { - QBdtNodeInterfacePtr leaf = root; - complex scale = leaf->scale; - for (bitLenInt j = 0U; j < bdtQubitCount; ++j) { - if (IS_NODE_0(leaf->scale)) { - break; - } - leaf = leaf->branches[SelectBit(i, j)]; - scale *= leaf->scale; - } - - if (!IS_NODE_0(leaf->scale) && attachedQubitCount) { - scale *= NODE_TO_QENGINE(leaf)->GetAmplitude(i >> bdtQubitCount); - } - - getLambda((bitCapIntOcl)i, scale); - } -} -template void QBdt::SetTraversal(Fn setLambda) -{ - Dump(); - - root = std::make_shared(); - root->Branch(bdtQubitCount); - - _par_for(maxQPower, [&](const bitCapInt& i, const unsigned& cpu) { - QBdtNodeInterfacePtr prevLeaf = root; - QBdtNodeInterfacePtr leaf = root; - for (bitLenInt j = 0U; j < bdtQubitCount; ++j) { - prevLeaf = leaf; - leaf = leaf->branches[SelectBit(i, j)]; - } - - if (attachedQubitCount) { - leaf = MakeQEngineNode(ONE_CMPLX, attachedQubitCount, 0U); - prevLeaf->branches[SelectBit(i, bdtQubitCount - 1U)] = leaf; - } - - setLambda((bitCapIntOcl)i, leaf); - }); - - root->PopStateVector(bdtQubitCount); - root->Prune(bdtQubitCount); -} -void QBdt::GetQuantumState(complex* state) -{ - GetTraversal([state](bitCapIntOcl i, complex scale) { state[i] = scale; }); -} -void QBdt::GetQuantumState(QInterfacePtr eng) -{ - GetTraversal([eng](bitCapIntOcl i, complex scale) { eng->SetAmplitude(i, scale); }); -} -void QBdt::SetQuantumState(const complex* state) -{ - if (!bdtQubitCount) { - NODE_TO_QENGINE(root)->SetQuantumState(state); - return; - } - - if (attachedQubitCount) { - const bitLenInt qbCount = bdtQubitCount; - SetTraversal([qbCount, state](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { - NODE_TO_QENGINE(leaf)->SetAmplitude(i >> qbCount, state[i]); - }); - } else { - SetTraversal([state](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { leaf->scale = state[i]; }); - } -} -void QBdt::SetQuantumState(QInterfacePtr eng) -{ - eng->Finish(); - - if (!bdtQubitCount) { - NODE_TO_QENGINE(root) = std::dynamic_pointer_cast(eng->Clone()); - return; - } - - if (attachedQubitCount) { - const bitLenInt qbCount = bdtQubitCount; - SetTraversal([qbCount, eng](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { - NODE_TO_QENGINE(leaf)->SetAmplitude(i >> qbCount, eng->GetAmplitude(i)); - }); - } else { - SetTraversal([eng](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { leaf->scale = eng->GetAmplitude(i); }); - } -} -void QBdt::GetProbs(real1* outputProbs) -{ - GetTraversal([outputProbs](bitCapIntOcl i, complex scale) { outputProbs[i] = norm(scale); }); -} - -void QBdt::SetStateVector() -{ - Finish(); - - if (!bdtQubitCount) { - return; - } - - QBdtQEngineNodePtr nRoot = MakeQEngineNode(ONE_R1, qubitCount); - GetQuantumState(NODE_TO_QENGINE(nRoot)); - root = nRoot; - SetQubitCount(qubitCount, qubitCount); -} -void QBdt::ResetStateVector(bitLenInt aqb) -{ - if (attachedQubitCount <= aqb) { - return; - } - - Finish(); - - if (!bdtQubitCount) { - QBdtQEngineNodePtr oRoot = std::dynamic_pointer_cast(root); - SetQubitCount(qubitCount, aqb); - SetQuantumState(NODE_TO_QENGINE(oRoot)); - } - - const bitLenInt length = attachedQubitCount - aqb; - const bitLenInt oBdtQubitCount = bdtQubitCount; - QBdtPtr nQubits = std::make_shared(engines, length, 0U, rand_generator, ONE_CMPLX, doNormalize, - randGlobalPhase, false, -1, (hardware_rand_generator == NULL) ? false : true, false, (real1_f)amplitudeFloor); - nQubits->SetQubitCount(length, 0U); - nQubits->SetPermutation(0U); - root->InsertAtDepth(nQubits->root, oBdtQubitCount, length); - SetQubitCount(qubitCount + length, attachedQubitCount); - for (bitLenInt i = 0U; i < length; ++i) { - Swap(oBdtQubitCount + i, oBdtQubitCount + length + i); - } - root->RemoveSeparableAtDepth(qubitCount - length, length); - SetQubitCount(qubitCount - length, 0U); -} - -void QBdt::SetDevice(int64_t dID) -{ - if (devID == dID) { - return; - } - - devID = dID; - - if (!attachedQubitCount) { - return; - } - - if (!bdtQubitCount) { - NODE_TO_QENGINE(root)->SetDevice(dID); - return; - } - - SetTraversal([dID](bitCapIntOcl i, QBdtNodeInterfacePtr leaf) { NODE_TO_QENGINE(leaf)->SetDevice(dID); }); -} - real1_f QBdt::SumSqrDiff(QBdtPtr toCompare) { if (this == toCompare.get()) { @@ -489,7 +268,7 @@ complex QBdt::GetAmplitude(bitCapInt perm) QBdtNodeInterfacePtr leaf = root; complex scale = leaf->scale; - for (bitLenInt j = 0U; j < bdtQubitCount; ++j) { + for (bitLenInt j = 0U; j < qubitCount; ++j) { if (IS_NODE_0(leaf->scale)) { break; } @@ -497,10 +276,6 @@ complex QBdt::GetAmplitude(bitCapInt perm) scale *= leaf->scale; } - if (!IS_NODE_0(leaf->scale) && attachedQubitCount) { - scale *= NODE_TO_QENGINE(leaf)->GetAmplitude(perm >> bdtQubitCount); - } - return scale; } @@ -514,58 +289,18 @@ bitLenInt QBdt::Compose(QBdtPtr toCopy, bitLenInt start) return start; } - if (maxPageQubits < (attachedQubitCount + toCopy->attachedQubitCount)) { - const bitLenInt diff = (attachedQubitCount + toCopy->attachedQubitCount) - maxPageQubits; - ResetStateVector((diff < qubitCount) ? (qubitCount - diff) : 0U); - - if (maxPageQubits < (attachedQubitCount + toCopy->attachedQubitCount)) { - const bitLenInt diff = (attachedQubitCount + toCopy->attachedQubitCount) - maxPageQubits; - if (toCopy->qubitCount < diff) { - throw std::runtime_error("Too many attached qubits to compose in QBdt::Compose()!"); - } - toCopy->ResetStateVector(toCopy->qubitCount - diff); - } - } - - if (bdtQubitCount && (attachedQubitCount || toCopy->attachedQubitCount)) { - if (start < bdtQubitCount) { - const bitLenInt offset = bdtQubitCount - start; - ROR(qubitCount - offset, 0U, qubitCount); - Compose(toCopy, offset); - ROL(qubitCount - offset, 0U, qubitCount); - - return start; - } - - if (start > bdtQubitCount) { - const bitLenInt offset = start - bdtQubitCount; - ROR(offset, 0U, qubitCount); - Compose(toCopy, qubitCount - offset); - ROL(offset, 0U, qubitCount); - - return start; - } - } - Finish(); - if (!bdtQubitCount && !toCopy->bdtQubitCount) { - NODE_TO_QENGINE(root)->Compose(NODE_TO_QENGINE(toCopy->root), start); - SetQubitCount(qubitCount + toCopy->qubitCount, qubitCount + toCopy->qubitCount); - - return start; - } - root->InsertAtDepth(toCopy->root->ShallowClone(), start, toCopy->qubitCount); - SetQubitCount(qubitCount + toCopy->qubitCount, attachedQubitCount + toCopy->attachedQubitCount); + SetQubitCount(qubitCount + toCopy->qubitCount); return start; } QInterfacePtr QBdt::Decompose(bitLenInt start, bitLenInt length) { - QBdtPtr dest = std::make_shared(engines, bdtQubitCount, length, rand_generator, ONE_CMPLX, doNormalize, - randGlobalPhase, false, -1, (hardware_rand_generator == NULL) ? false : true, false, (real1_f)amplitudeFloor); + QBdtPtr dest = std::make_shared(engines, length, 0U, rand_generator, ONE_CMPLX, doNormalize, randGlobalPhase, + false, -1, (hardware_rand_generator == NULL) ? false : true, false, (real1_f)amplitudeFloor); Decompose(start, dest); @@ -582,32 +317,6 @@ void QBdt::DecomposeDispose(bitLenInt start, bitLenInt length, QBdtPtr dest) return; } - if (length > bdtQubitCount) { - if (dest) { - ExecuteAsStateVector([&](QInterfacePtr eng) { - dest->SetStateVector(); - eng->Decompose(start, NODE_TO_QENGINE(dest->root)); - SetQubitCount(qubitCount - length); - dest->ResetStateVector(); - }); - } else { - ExecuteAsStateVector([&](QInterfacePtr eng) { - eng->Dispose(start, length); - SetQubitCount(qubitCount - length); - }); - } - - return; - } - - if (start && bdtQubitCount && attachedQubitCount) { - ROR(start, 0U, qubitCount); - DecomposeDispose(0U, length, dest); - ROL(start, 0U, qubitCount); - - return; - } - Finish(); if (dest) { @@ -616,8 +325,8 @@ void QBdt::DecomposeDispose(bitLenInt start, bitLenInt length, QBdtPtr dest) } else { root->RemoveSeparableAtDepth(start, length); } - SetQubitCount(qubitCount - length, attachedQubitCount); - root->Prune(bdtQubitCount); + SetQubitCount(qubitCount - length); + root->Prune(qubitCount); } bitLenInt QBdt::Allocate(bitLenInt start, bitLenInt length) @@ -630,11 +339,10 @@ bitLenInt QBdt::Allocate(bitLenInt start, bitLenInt length) QBdtPtr nQubits = std::make_shared(engines, length, 0U, rand_generator, ONE_CMPLX, doNormalize, randGlobalPhase, false, -1, (hardware_rand_generator == NULL) ? false : true, false, (real1_f)amplitudeFloor); - nQubits->SetQubitCount(length, 0U); nQubits->SetPermutation(0U); nQubits->root->InsertAtDepth(root, length, qubitCount); root = nQubits->root; - SetQubitCount(qubitCount + length, attachedQubitCount); + SetQubitCount(qubitCount + length); ROR(length, 0U, start + length); return start; @@ -645,14 +353,9 @@ real1_f QBdt::Prob(bitLenInt qubit) if (qubit >= qubitCount) { throw std::invalid_argument("QBdt::Prob qubit index parameter must be within allocated qubit bounds!"); } - - const bool isKet = (qubit >= bdtQubitCount); - const bitLenInt maxQubit = isKet ? bdtQubitCount : qubit; - const bitCapInt qPower = pow2(maxQubit); - - std::map qiProbs; - + const bitCapInt qPower = pow2(qubit); const unsigned numCores = GetConcurrencyLevel(); + std::map qiProbs; std::unique_ptr oneChanceBuff(new real1[numCores]()); Finish(); @@ -660,7 +363,7 @@ real1_f QBdt::Prob(bitLenInt qubit) _par_for(qPower, [&](const bitCapInt& i, const unsigned& cpu) { QBdtNodeInterfacePtr leaf = root; complex scale = leaf->scale; - for (bitLenInt j = 0U; j < maxQubit; ++j) { + for (bitLenInt j = 0U; j < qubit; ++j) { if (IS_NODE_0(leaf->scale)) { break; } @@ -672,17 +375,6 @@ real1_f QBdt::Prob(bitLenInt qubit) return; } - if (isKet) { - // Phase effects don't matter, for probability expectation. - QEnginePtr qi = NODE_TO_QENGINE(leaf); - if (qiProbs.find(qi) == qiProbs.end()) { - qiProbs[qi] = sqrt(NODE_TO_QENGINE(leaf)->Prob(qubit - bdtQubitCount)); - } - oneChanceBuff[cpu] += norm(scale * qiProbs[qi]); - - return; - } - oneChanceBuff[cpu] += norm(scale * leaf->branches[1U]->scale); }); @@ -701,7 +393,7 @@ real1_f QBdt::ProbAll(bitCapInt perm) QBdtNodeInterfacePtr leaf = root; complex scale = leaf->scale; - for (bitLenInt j = 0U; j < bdtQubitCount; ++j) { + for (bitLenInt j = 0U; j < qubitCount; ++j) { if (IS_NODE_0(leaf->scale)) { break; } @@ -709,10 +401,6 @@ real1_f QBdt::ProbAll(bitCapInt perm) scale *= leaf->scale; } - if (!IS_NODE_0(leaf->scale) && attachedQubitCount) { - scale *= NODE_TO_QENGINE(leaf)->GetAmplitude(perm >> bdtQubitCount); - } - return clampProb((real1_f)norm(scale)); } @@ -735,14 +423,12 @@ bool QBdt::ForceM(bitLenInt qubit, bool result, bool doForce, bool doApply) return result; } - const bool isKet = (qubit >= bdtQubitCount); - const bitLenInt maxQubit = isKet ? bdtQubitCount : qubit; - const bitCapInt qPower = pow2(maxQubit); + const bitCapInt qPower = pow2(qubit); root->scale = GetNonunitaryPhase(); _par_for(qPower, [&](const bitCapInt& i, const unsigned& cpu) { QBdtNodeInterfacePtr leaf = root; - for (bitLenInt j = 0U; j < maxQubit; ++j) { + for (bitLenInt j = 0U; j < qubit; ++j) { if (IS_NODE_0(leaf->scale)) { break; } @@ -758,11 +444,6 @@ bool QBdt::ForceM(bitLenInt qubit, bool result, bool doForce, bool doApply) leaf->Branch(); - if (isKet) { - NODE_TO_QENGINE(leaf)->ForceM(qubit - bdtQubitCount, result, true, true); - return; - } - QBdtNodeInterfacePtr& b0 = leaf->branches[0U]; QBdtNodeInterfacePtr& b1 = leaf->branches[1U]; @@ -783,7 +464,7 @@ bool QBdt::ForceM(bitLenInt qubit, bool result, bool doForce, bool doApply) } }); - root->Prune(maxQubit); + root->Prune(qubit); return result; } @@ -795,7 +476,7 @@ bitCapInt QBdt::MAll() Finish(); - for (bitLenInt i = 0U; i < bdtQubitCount; ++i) { + for (bitLenInt i = 0U; i < qubitCount; ++i) { leaf->Branch(); real1_f oneChance = clampProb((real1_f)norm(leaf->branches[1U]->scale)); bool bitResult; @@ -819,11 +500,6 @@ bitCapInt QBdt::MAll() } } - if (bdtQubitCount < qubitCount) { - // Theoretically, there's only 1 copy of this leaf left, so no need to branch. - result |= NODE_TO_QENGINE(leaf)->MAll() << bdtQubitCount; - } - return result; } @@ -833,19 +509,12 @@ void QBdt::ApplySingle(const complex* mtrx, bitLenInt target) throw std::invalid_argument("QBdt::ApplySingle target parameter must be within allocated qubit bounds!"); } - if (!bdtQubitCount) { - NODE_TO_QENGINE(root)->Mtrx(mtrx, target); - return; - } - if (IS_NORM_0(mtrx[1U]) && IS_NORM_0(mtrx[2U]) && IS_NORM_0(mtrx[0U] - mtrx[3U]) && (randGlobalPhase || IS_NORM_0(ONE_CMPLX - mtrx[0U]))) { return; } - const bool isKet = (target >= bdtQubitCount); - const bitLenInt maxQubit = isKet ? bdtQubitCount : target; - const bitCapInt qPower = pow2(maxQubit); + const bitCapInt qPower = pow2(target); #if ENABLE_COMPLEX_X2 const complex2 mtrxCol1(mtrx[0U], mtrx[2U]); @@ -855,21 +524,20 @@ void QBdt::ApplySingle(const complex* mtrx, bitLenInt target) const complex2 mtrxCol2Shuff = mtrxColShuff(mtrxCol2); #endif - par_for_qbdt(qPower, maxQubit, + par_for_qbdt(qPower, target, #if ENABLE_COMPLEX_X2 - [this, maxQubit, target, mtrx, &mtrxCol1, &mtrxCol2, &mtrxCol1Shuff, &mtrxCol2Shuff, isKet]( - const bitCapInt& i) { + [this, target, mtrx, &mtrxCol1, &mtrxCol2, &mtrxCol1Shuff, &mtrxCol2Shuff](const bitCapInt& i) { #else - [this, maxQubit, target, mtrx, isKet](const bitCapInt& i) { + [this, target, mtrx](const bitCapInt& i) { #endif QBdtNodeInterfacePtr leaf = root; // Iterate to qubit depth. - for (bitLenInt j = 0U; j < maxQubit; ++j) { + for (bitLenInt j = 0U; j < target; ++j) { if (IS_NODE_0(leaf->scale)) { // WARNING: Mutates loop control variable! - return (bitCapInt)(pow2(maxQubit - j) - ONE_BCI); + return (bitCapInt)(pow2(target - j) - ONE_BCI); } - leaf = leaf->branches[SelectBit(i, maxQubit - (j + 1U))]; + leaf = leaf->branches[SelectBit(i, target - (j + 1U))]; } std::lock_guard lock(leaf->mtx); @@ -878,16 +546,11 @@ void QBdt::ApplySingle(const complex* mtrx, bitLenInt target) return (bitCapInt)0U; } - if (isKet) { - leaf->Branch(); - NODE_TO_QENGINE(leaf)->Mtrx(mtrx, target - bdtQubitCount); - } else { #if ENABLE_COMPLEX_X2 - leaf->Apply2x2(mtrxCol1, mtrxCol2, mtrxCol1Shuff, mtrxCol2Shuff, bdtQubitCount - target); + leaf->Apply2x2(mtrxCol1, mtrxCol2, mtrxCol1Shuff, mtrxCol2Shuff, qubitCount - target); #else - leaf->Apply2x2(mtrx, bdtQubitCount - target); + leaf->Apply2x2(mtrx, qubitCount - target); #endif - } return (bitCapInt)0U; }); @@ -904,15 +567,6 @@ void QBdt::ApplyControlledSingle( ThrowIfQbIdArrayIsBad(controls, qubitCount, "QBdt::ApplyControlledSingle parameter controls array values must be within allocated qubit bounds!"); - if (!bdtQubitCount) { - if (isAnti) { - NODE_TO_QENGINE(root)->MACMtrx(controls, mtrx, target); - } else { - NODE_TO_QENGINE(root)->MCMtrx(controls, mtrx, target); - } - return; - } - if (IS_NORM_0(mtrx[1U]) && IS_NORM_0(mtrx[2U]) && IS_NORM_0(ONE_CMPLX - mtrx[0U]) && IS_NORM_0(ONE_CMPLX - mtrx[3U])) { return; @@ -920,26 +574,19 @@ void QBdt::ApplyControlledSingle( std::vector controlVec(controls.begin(), controls.end()); std::sort(controlVec.begin(), controlVec.end()); - const bool isSwapped = (target < controlVec.back()) && (target < bdtQubitCount); + const bool isSwapped = target < controlVec.back(); if (isSwapped) { Swap(target, controlVec.back()); std::swap(target, controlVec.back()); } - const bool isKet = (target >= bdtQubitCount); - const bitLenInt maxQubit = isKet ? bdtQubitCount : target; - const bitCapInt qPower = pow2(maxQubit); - std::vector ketControlsVec; - bitCapInt lowControlMask = 0U; + const bitCapInt qPower = pow2(target); + bitCapInt controlMask = 0U; for (size_t c = 0U; c < controls.size(); ++c) { const bitLenInt control = controlVec[c]; - if (control < bdtQubitCount) { - lowControlMask |= pow2(maxQubit - (control + 1U)); - } else { - ketControlsVec.push_back(control - bdtQubitCount); - } + controlMask |= pow2(target - (control + 1U)); } - bitCapInt lowControlPerm = isAnti ? 0U : lowControlMask; + const bitCapInt controlPerm = isAnti ? 0U : controlMask; #if ENABLE_COMPLEX_X2 const complex2 mtrxCol1(mtrx[0U], mtrx[2U]); @@ -949,25 +596,25 @@ void QBdt::ApplyControlledSingle( const complex2 mtrxCol2Shuff = mtrxColShuff(mtrxCol2); #endif - par_for_qbdt(qPower, maxQubit, + par_for_qbdt(qPower, target, #if ENABLE_COMPLEX_X2 - [this, lowControlMask, lowControlPerm, maxQubit, target, mtrx, &mtrxCol1, &mtrxCol2, &mtrxCol1Shuff, - &mtrxCol2Shuff, isKet, isAnti, ketControlsVec](const bitCapInt& i) { + [this, controlMask, controlPerm, target, mtrx, &mtrxCol1, &mtrxCol2, &mtrxCol1Shuff, &mtrxCol2Shuff, isAnti]( + const bitCapInt& i) { #else - [this, lowControlMask, lowControlPerm, maxQubit, target, mtrx, isKet, isAnti, ketControlsVec](const bitCapInt& i) { + [this, controlMask, controlPerm, target, mtrx, isAnti](const bitCapInt& i) { #endif - if ((i & lowControlMask) != lowControlPerm) { - return (bitCapInt)(lowControlMask - ONE_BCI); + if ((i & controlMask) != controlPerm) { + return (bitCapInt)(controlMask - ONE_BCI); } QBdtNodeInterfacePtr leaf = root; // Iterate to qubit depth. - for (bitLenInt j = 0U; j < maxQubit; ++j) { + for (bitLenInt j = 0U; j < target; ++j) { if (IS_NODE_0(leaf->scale)) { // WARNING: Mutates loop control variable! - return (bitCapInt)(pow2(maxQubit - j) - ONE_BCI); + return (bitCapInt)(pow2(target - j) - ONE_BCI); } - leaf = leaf->branches[SelectBit(i, maxQubit - (j + 1U))]; + leaf = leaf->branches[SelectBit(i, target - (j + 1U))]; } std::lock_guard lock(leaf->mtx); @@ -976,21 +623,11 @@ void QBdt::ApplyControlledSingle( return (bitCapInt)0U; } - if (isKet) { - leaf->Branch(); - QEnginePtr qi = NODE_TO_QENGINE(leaf); - if (isAnti) { - qi->MACMtrx(ketControlsVec, mtrx, target - bdtQubitCount); - } else { - qi->MCMtrx(ketControlsVec, mtrx, target - bdtQubitCount); - } - } else { #if ENABLE_COMPLEX_X2 - leaf->Apply2x2(mtrxCol1, mtrxCol2, mtrxCol1Shuff, mtrxCol2Shuff, bdtQubitCount - target); + leaf->Apply2x2(mtrxCol1, mtrxCol2, mtrxCol1Shuff, mtrxCol2Shuff, qubitCount - target); #else - leaf->Apply2x2(mtrx, bdtQubitCount - target); + leaf->Apply2x2(mtrx, qubitCount - target); #endif - } return (bitCapInt)0U; }); @@ -1072,13 +709,15 @@ void QBdt::MCInvert(const std::vector& controls, complex topRight, co } std::vector lControls(controls); - lControls.push_back(target); std::sort(lControls.begin(), lControls.end()); - const bitLenInt lTarget = lControls.back(); - lControls.pop_back(); + + if (lControls[controls.size() - 1U] < target) { + ApplyControlledSingle(mtrx, lControls, target, false); + return; + } H(target); - MCPhase(lControls, ONE_CMPLX, -ONE_CMPLX, lTarget); + MCPhase(lControls, ONE_CMPLX, -ONE_CMPLX, target); H(target); } diff --git a/src/qstabilizerhybrid.cpp b/src/qstabilizerhybrid.cpp index 9506a0100..f498b603a 100644 --- a/src/qstabilizerhybrid.cpp +++ b/src/qstabilizerhybrid.cpp @@ -402,14 +402,8 @@ void QStabilizerHybrid::SwitchToEngine() return; } - const bool isBdt = engineTypes.size() && (engineTypes[0] == QINTERFACE_BDT); - if ((qubitCount + ancillaCount) > maxEngineQubitCount) { QInterfacePtr e = MakeEngine(0); - if (isBdt) { - std::dynamic_pointer_cast(e)->SetStateVector(); - } - #if ENABLE_QUNIT_CPU_PARALLEL && ENABLE_PTHREAD const unsigned numCores = GetConcurrencyLevel(); std::vector clones; @@ -446,11 +440,6 @@ void QStabilizerHybrid::SwitchToEngine() if (!doNormalize) { engine->NormalizeState(); } - - if (isBdt) { - std::dynamic_pointer_cast(engine)->ResetStateVector(); - } - // We have extra "gate fusion" shards leftover. shards.erase(shards.begin() + qubitCount, shards.end()); // We're done with ancillae. @@ -460,25 +449,16 @@ void QStabilizerHybrid::SwitchToEngine() } engine = MakeEngine(0, stabilizer->GetQubitCount()); - if (isBdt) { - std::dynamic_pointer_cast(engine)->SetStateVector(); - } stabilizer->GetQuantumState(engine); stabilizer = NULL; FlushBuffers(); if (!ancillaCount) { - if (isBdt) { - std::dynamic_pointer_cast(engine)->ResetStateVector(); - } return; } // When we measure, we act postselection on reverse T-gadgets. engine->ForceMReg(qubitCount, ancillaCount, 0, true, true); - if (isBdt) { - std::dynamic_pointer_cast(engine)->ResetStateVector(); - } // Ancillae are separable after measurement. engine->Dispose(qubitCount, ancillaCount); // We have extra "gate fusion" shards leftover.