PathIntegralSimulator.h

Go to the documentation of this file.

 
#pragma once
 
#ifndef _PATH_INTEGRAL_SIMULATOR_H_
#define _PATH_INTEGRAL_SIMULATOR_H_
 
#include <optional>
 
#include "PathIntegral.h"
#include "../Circuit/Circuit.h"
 
namespace Simulators {
 
    class PathIntegralSimulator {
        public:
          void SetSeed(uint64_t seed) { simulator.SetSeed(seed); }
          void SetTrimValue(double val) { simulator.SetTrimValue(val); }
 
          double GetTrimValue() const { return simulator.GetTrimValue(); }
 
          void SetMaxDoublingsForBackwardPaths(size_t doublings) {
            simulator.SetMaxDoublingsForBackwardPaths(doublings);
          }
 
          size_t GetMaxDoublingsForBackwardPaths() const { return simulator.GetMaxDoublingsForBackwardPaths(); }
 
          void Reset() { simulator.Reset(); }
 
          bool SetCircuit(const std::shared_ptr<Circuits::Circuit<>>& circuit)
          {
            const auto convertedCircuit = ConvertCircuit(circuit);
            if (convertedCircuit.empty()) return false;
            simulator.SetCircuit(convertedCircuit);
            return true;
          }
 
          std::complex<double> AmplitudeFromZero(const std::vector<bool>& endState)
          {
            return simulator.Propagate(endState);
          }
 
          std::complex<double> Amplitude(const std::vector<bool>& startState, const std::vector<bool>& endState)
          {
            return simulator.Propagate(startState, endState);
          }
 
          std::optional<QC::Gates::AppliedGate<>> ConvertGate(const std::shared_ptr<Circuits::IQuantumGate<>>& gate)
          {
            if (!gate) return std::nullopt;
 
            switch (gate->GetGateType())
            {
              // single qubit gates
            case Circuits::QuantumGateType::kPhaseGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::PhaseGate<>>(gate);
              pgate.SetPhaseShift(g->GetLambda());
              return QC::Gates::AppliedGate<>(pgate.getRawOperatorMatrix(), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kXGateType:
              return QC::Gates::AppliedGate<>(xgate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kYGateType:
              return QC::Gates::AppliedGate<>(ygate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kZGateType:
              return QC::Gates::AppliedGate<>(zgate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kHadamardGateType:
              return QC::Gates::AppliedGate<>(h.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kSGateType:
              return QC::Gates::AppliedGate<>(sgate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kSdgGateType:
              return QC::Gates::AppliedGate<>(sdggate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kTGateType:
              return QC::Gates::AppliedGate<>(tgate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kTdgGateType:
              return QC::Gates::AppliedGate<>(tdggate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kSxGateType:
              return QC::Gates::AppliedGate<>(sxgate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kSxDagGateType:
              return QC::Gates::AppliedGate<>(sxdaggate.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kKGateType:
              return QC::Gates::AppliedGate<>(k.getRawOperatorMatrix(), gate->GetQubit(0));
            case Circuits::QuantumGateType::kRxGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::RxGate<>>(gate);
              rxgate.SetTheta(g->GetTheta());
              return QC::Gates::AppliedGate<>(rxgate.getRawOperatorMatrix(), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kRyGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::RyGate<>>(gate);
              rygate.SetTheta(g->GetTheta());
              return QC::Gates::AppliedGate<>(rygate.getRawOperatorMatrix(), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kRzGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::RzGate<>>(gate);
              rzgate.SetTheta(g->GetTheta());
              return QC::Gates::AppliedGate<>(rzgate.getRawOperatorMatrix(), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kUGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::UGate<>>(gate);
              ugate.SetParams(g->GetTheta(), g->GetPhi(), g->GetLambda(), g->GetGamma());
              return QC::Gates::AppliedGate<>(ugate.getRawOperatorMatrix(), gate->GetQubit(0));
            }
 
            // two qubit gates
            case Circuits::QuantumGateType::kSwapGateType:
              return QC::Gates::AppliedGate<>(swapgate.getRawOperatorMatrix(), gate->GetQubit(0), gate->GetQubit(1));
            case Circuits::QuantumGateType::kCXGateType:
              return QC::Gates::AppliedGate<>(cxgate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            case Circuits::QuantumGateType::kCYGateType:
              return QC::Gates::AppliedGate<>(cygate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            case Circuits::QuantumGateType::kCZGateType:
              return QC::Gates::AppliedGate<>(czgate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            case Circuits::QuantumGateType::kCPGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::CPGate<>>(gate);
              cpgate.SetPhaseShift(g->GetLambda());
              return QC::Gates::AppliedGate<>(cpgate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kCRxGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::CRxGate<>>(gate);
              crxgate.SetTheta(g->GetTheta());
              return QC::Gates::AppliedGate<>(crxgate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kCRyGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::CRyGate<>>(gate);
              crygate.SetTheta(g->GetTheta());
              return QC::Gates::AppliedGate<>(crygate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kCRzGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::CRzGate<>>(gate);
              crzgate.SetTheta(g->GetTheta());
              return QC::Gates::AppliedGate<>(crzgate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            }
            case Circuits::QuantumGateType::kCHGateType:
              return QC::Gates::AppliedGate<>(ch.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            case Circuits::QuantumGateType::kCSxGateType:
              return QC::Gates::AppliedGate<>(csx.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            case Circuits::QuantumGateType::kCSxDagGateType:
              return QC::Gates::AppliedGate<>(csxdag.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            case Circuits::QuantumGateType::kCUGateType:
            {
              const auto g = std::static_pointer_cast<Circuits::CUGate<>>(gate);
              cugate.SetParams(g->GetTheta(), g->GetPhi(), g->GetLambda(), g->GetGamma());
              return QC::Gates::AppliedGate<>(cugate.getRawOperatorMatrix(), gate->GetQubit(1), gate->GetQubit(0));
            }
 
            // three qubit gates
            case Circuits::QuantumGateType::kCCXGateType:
              return QC::Gates::AppliedGate<>(ccxgate.getRawOperatorMatrix(), gate->GetQubit(2), gate->GetQubit(1), gate->GetQubit(0));
            case Circuits::QuantumGateType::kCSwapGateType:
              return QC::Gates::AppliedGate<>(cswapgate.getRawOperatorMatrix(), gate->GetQubit(2), gate->GetQubit(1), gate->GetQubit(0));
 
            default:
              return std::nullopt;
            }
          }
 
          std::vector<QC::Gates::AppliedGate<>> ConvertCircuit(const std::shared_ptr<Circuits::Circuit<>>& circuit)
          {
            std::vector<QC::Gates::AppliedGate<>> result;
            if (!circuit) return result;
 
            const auto& operations = circuit->GetOperations();
            result.reserve(operations.size());
 
            for (const auto& op : operations)
            {
              if (!op || op->GetType() != Circuits::OperationType::kGate)
                continue;
 
              const auto& gate = std::static_pointer_cast<Circuits::IQuantumGate<>>(op);
              if (auto applied = ConvertGate(gate))
                result.emplace_back(std::move(*applied));
            }
 
            return result;
          }
 
          void PropagateStep(
              const QC::Gates::AppliedGate<>& gate,
              std::unordered_map<QC::PathIntegral::FastVectorBool, std::complex<double>,
              QC::PathIntegral::FastVectorBoolHash>& currentAmplitudes)
          {
            simulator.PropagateStep(gate, currentAmplitudes);
          }
 
          static size_t GetBranchingForQcsimCircuit(const std::vector<QC::Gates::AppliedGate<>>& circuit) {
            size_t doublings = 0;
            for (const auto& gate : circuit) 
              if (gate.isBranching()) ++doublings;
 
            return doublings;
          }
 
          size_t GetBranchingForMaestroCircuit(const std::shared_ptr<Circuits::Circuit<>>& circuit) {
            /*
            const auto circ = ConvertCircuit(circuit);
            return GetBranchingForQcsimCircuit(circ);
            */
           
            return GetBranchingForMaestroCircuitFromPos(circuit, 0);
          }
 
          size_t GetBranchingForMaestroCircuitFromPos(
              const std::shared_ptr<Circuits::Circuit<>>& circuit, size_t pos) {
            if (!circuit) return 0;
            
            size_t doublings = 0;
            const auto& operations = circuit->GetOperations();
            for (size_t i = pos; i < operations.size(); ++i)
                if (operations[i]->IsBranching()) ++doublings;
 
            return doublings;
          }
 
          double QubitProbability(size_t qubit, bool value = true) {
            return simulator.QubitProbability(qubit, value);
          }
 
          std::complex<double> AmplitudeForOutcome(size_t outcome) {
            const auto& amplitudes = simulator.GetAmplitudes();
            if (amplitudes.empty()) return 0.0;
 
            QC::PathIntegral::FastVectorBool state(
                amplitudes.begin()->first.size());
            for (size_t q = 0; q < state.size(); ++q)
              state.set(q, (outcome >> q) & 1);
            if (auto it = amplitudes.find(state); it != amplitudes.end())
              return it->second;
            return 0.0;
          }
 
          double Probability(size_t outcome) {
            return std::norm(AmplitudeForOutcome(outcome));
          }
 
          QC::PathIntegral::FastVectorBool MeasureNoCollapse() {
            return simulator.MeasureNoCollapse();
          }
 
          bool MeasureQubit(size_t qubit) {
            return simulator.MeasureQubit(qubit);
          }
 
          std::unique_ptr<PathIntegralSimulator> Clone() const {
            auto clone = std::make_unique<PathIntegralSimulator>();
            clone->simulator = simulator.Clone();
            return clone;
          }
 
          std::unordered_map<QC::PathIntegral::FastVectorBool, std::complex<double>,
                             QC::PathIntegral::FastVectorBoolHash>&
              Amplitudes()
          {
            return simulator.GetAmplitudes();
          }
 
          void SetStartZeroState(size_t numQubits) {
            auto& amplitudes = simulator.GetAmplitudes();
 
            amplitudes.clear();
            QC::PathIntegral::FastVectorBool zeroState(numQubits);
            amplitudes[zeroState] = 1.0;
          }
 
          double ExpectationValue(const std::string& pauliStringOrig) {
            // save current amplitudes
            const auto& currentAmplitudes = simulator.GetAmplitudes();
            auto amplitudes = currentAmplitudes;
 
            // now apply the pauli string as gates
            for (size_t i = 0; i < pauliStringOrig.size(); ++i)
            {
              const char c = pauliStringOrig[i];
              switch (c)
              {
              case 'X':
                  [[fallthrough]];
              case 'x':
                PropagateStep(QC::Gates::AppliedGate<>(xgate.getRawOperatorMatrix(), i), amplitudes);
                break;
              case 'Y':
                  [[fallthrough]];
              case 'y':
                PropagateStep(QC::Gates::AppliedGate<>(ygate.getRawOperatorMatrix(), i), amplitudes);
                break;
              case 'Z':
                  [[fallthrough]];
              case 'z':
                PropagateStep(QC::Gates::AppliedGate<>(zgate.getRawOperatorMatrix(), i), amplitudes);
                break;
              default:
                break;
              }
            }
 
            // use the new amplitudes to calculate the expectation value
            double expectation = 0.0;
 
            for (const auto& [state, amplitude] : amplitudes)
            {
              if (currentAmplitudes.find(state) == currentAmplitudes.end())
                continue;
 
              const auto& oldAmplitude = currentAmplitudes.at(state);
              expectation += std::real(std::conj(oldAmplitude) * amplitude);
            }
 
            return expectation;
          }
 
          void ApplyGate(const QC::Gates::AppliedGate<>& gate) {
            simulator.PropagateStep(gate, simulator.GetAmplitudes());
          }
 
          void SaveState() { simulator.SaveAmplitudes(); }
 
          void RestoreState() { simulator.RestoreAmplitudes(); }
 
        private:
          QC::PathIntegral::PathIntegralSimulator simulator;
 
          QC::Gates::PhaseShiftGate<> pgate;
          QC::Gates::PauliXGate<> xgate;
          QC::Gates::PauliYGate<> ygate;
          QC::Gates::PauliZGate<> zgate;
          QC::Gates::HadamardGate<> h;
          // QC::Gates::UGate<> ugate;
          QC::Gates::SGate<> sgate;
          QC::Gates::SDGGate<> sdggate;
          QC::Gates::TGate<> tgate;
          QC::Gates::TDGGate<> tdggate;
          QC::Gates::SquareRootNOTGate<> sxgate;
          QC::Gates::SquareRootNOTDagGate<> sxdaggate;
          QC::Gates::HyGate<> k;
          QC::Gates::RxGate<> rxgate;
          QC::Gates::RyGate<> rygate;
          QC::Gates::RzGate<> rzgate;
          QC::Gates::UGate<> ugate;
          QC::Gates::CNOTGate<> cxgate;
          QC::Gates::ControlledYGate<> cygate;
          QC::Gates::ControlledZGate<> czgate;
          QC::Gates::ControlledPhaseShiftGate<> cpgate;
          QC::Gates::ControlledRxGate<> crxgate;
          QC::Gates::ControlledRyGate<> crygate;
          QC::Gates::ControlledRzGate<> crzgate;
          QC::Gates::ControlledHadamardGate<> ch;
          QC::Gates::ControlledSquareRootNOTGate<> csx;
          QC::Gates::ControlledSquareRootNOTDagGate<> csxdag;
          QC::Gates::SwapGate<> swapgate;
          QC::Gates::ToffoliGate<> ccxgate;
          QC::Gates::FredkinGate<> cswapgate;
          QC::Gates::ControlledUGate<> cugate;
    };
 
}
 
#endif  // _PATH_INTEGRAL_SIMULATOR_H_