State.h

Go to the documentation of this file.

 
#pragma once
 
#ifndef _SIMULATOR_STATE_H_
#define _SIMULATOR_STATE_H_
 
#include <Eigen/Eigen>
#include <cmath>
#include <complex>
#include <stdexcept>
#include <unordered_map>
#include <unordered_set>
#include <utility>
#include <vector>
 
#ifndef NO_QISKIT_AER
#include "framework/linalg/vector.hpp"
#endif
 
#include "QuantumChannel.h"
#include "SimulatorObserver.h"
 
namespace Circuits {
template <typename Time> class IOperation;
}
 
namespace Simulators {

/*
template<class T> class avoid_init_allocator : public std::allocator<T>
{
public:
        using std::allocator<T>::allocator;
 
        template <class U, class... Args> void construct(U*, Args&&...) {}
};
*/

enum class SimulatorType : int {
#ifndef NO_QISKIT_AER
  kQiskitAer, 
#endif
  kQCSim, 
#ifndef NO_QISKIT_AER
  kCompositeQiskitAer, 
#endif
  kCompositeQCSim, 
  kGpuSim,         
  kQuestSim,       
  kDistGpuSim,     
  kDistMpiGpuSim   
};
 
inline bool IsDistributedGpuSimulator(SimulatorType type) {
  return type == SimulatorType::kDistGpuSim || type == SimulatorType::kDistMpiGpuSim;
}
inline bool IsGpuSimulator(SimulatorType type) {
  return type == SimulatorType::kGpuSim || IsDistributedGpuSimulator(type);
}

enum class SimulationType : int {
  kStatevector,        
  kMatrixProductState, 
  kStabilizer,         
  kTensorNetwork,      
  kPauliPropagator,    
  kExtendedStabilizer, 
  kPathIntegral,       
  kDensityMatrix,      
  kMatrixProductOperator, 
  kOther 
};

class IState {
 public:
  virtual ~IState() = default;

  virtual void SetSeed(uint64_t seed) {
    const std::string value = std::to_string(seed);
    Configure("seed", value.c_str());
  }
 
  static uint64_t DeriveSeed(uint64_t seed, uint64_t stream) {
    uint64_t value = seed + 0x9e3779b97f4a7c15ULL * (stream + 1);
    value = (value ^ (value >> 30)) * 0xbf58476d1ce4e5b9ULL;
    value = (value ^ (value >> 27)) * 0x94d049bb133111ebULL;
    return value ^ (value >> 31);
  }

  virtual void Initialize() = 0;

  virtual void InitializeState(
      size_t num_qubits, std::vector<std::complex<double>> &amplitudes) = 0;

#ifndef NO_QISKIT_AER
  virtual void InitializeState(
      size_t num_qubits, AER::Vector<std::complex<double>> &amplitudes) = 0;
#endif

  virtual void InitializeState(size_t num_qubits,
                               Eigen::VectorXcd &amplitudes) = 0;

  virtual void InitializeToBasisState(size_t num_qubits,
                                      Types::qubit_t basisState) {
    (void)num_qubits;
    (void)basisState;
    throw std::runtime_error(
        "This simulator does not support initialization to a computational "
        "basis state");
  }

  virtual void InitializeToBasisState(size_t num_qubits,
                                      const std::vector<bool> &basisState) {
    (void)num_qubits;
    (void)basisState;
    throw std::runtime_error(
        "This simulator does not support initialization to a computational "
        "basis state");
  }

  virtual void InitializeToMixtureOfBasisStates(
      size_t num_qubits,
      const std::vector<std::pair<Types::qubit_t, double>> &mixture) {
    (void)num_qubits;
    (void)mixture;
    throw std::runtime_error(
        "This simulator does not support initialization to a mixture of "
        "computational basis states");
  }

  virtual void InitializeToMixtureOfBasisStates(
      size_t num_qubits,
      const std::vector<std::pair<std::vector<bool>, double>> &mixture) {
    (void)num_qubits;
    (void)mixture;
    throw std::runtime_error(
        "This simulator does not support initialization to a mixture of "
        "computational basis states");
  }

  virtual void Reset() = 0;

  virtual bool SupportsMPSSwapOptimization() const { return false; }

  virtual void SetInitialQubitsMap(
      const std::vector<long long int> &initialMap) {}

  virtual void SetUseOptimalMeetingPosition(bool /*enable*/) {}

  virtual void SetLookaheadDepth(int /*depth*/) {}

  virtual void SetLookaheadDepthWithHeuristic(int /*depth*/) {}

  virtual void SetUpcomingGates(
      const std::vector<std::shared_ptr<Circuits::IOperation<double>>> & /*gates*/) {}

  virtual long long int GetGatesCounter() const { return 0; }

  virtual void SetGatesCounter(long long int /*counter*/) {}

  virtual void IncrementGatesCounter() {}
 
  // the following four functions are also for MPS swaps optimizations, might be removed in the future
 
  virtual double getGrowthFactorSwap() const { return 0.; }
  virtual double getGrowthFactorGate() const { return 0.; }
 
  virtual void setGrowthFactorSwap(double factor) { }
  virtual void setGrowthFactorGate(double factor) { }

  virtual void Configure(const char *key, const char *value) = 0;

  virtual std::string GetConfiguration(const char *key) const = 0;

  virtual size_t AllocateQubits(size_t num_qubits) = 0;

  virtual size_t GetNumberOfQubits() const = 0;

  virtual void Clear() = 0;

  virtual size_t Measure(const Types::qubits_vector &qubits) = 0;

  virtual std::vector<bool> MeasureMany(const Types::qubits_vector &qubits) = 0;

  virtual void ApplyReset(const Types::qubits_vector &qubits) = 0;

  virtual bool SupportsQuantumChannels() const { return false; }

  virtual void ApplyQuantumChannel(const Types::qubits_vector &targets,
                                   const QuantumChannel &channel) {
    (void)targets;
    (void)channel;
    throw std::runtime_error(
        "This simulator does not support exact quantum-channel evolution");
  }

  virtual std::complex<double> DensityMatrixTrace() const {
    throw std::runtime_error("This simulator does not expose a density-matrix trace");
  }
  virtual double DensityMatrixPurity() const {
    throw std::runtime_error("This simulator does not expose density-matrix purity");
  }
  virtual std::complex<double> DensityMatrixTraceOfSquare() const {
    throw std::runtime_error("This simulator does not expose Tr(rho^2)");
  }
  virtual std::complex<double> DensityMatrixOverlap(const IState &) const {
    throw std::runtime_error("This simulator does not support density-matrix overlap");
  }
  virtual double DensityMatrixHermiticityResidual() const {
    throw std::runtime_error("This simulator does not expose a Hermiticity residual");
  }
  virtual bool IsDensityMatrixHermitian(double = 1e-10) const {
    throw std::runtime_error("This simulator does not expose a Hermiticity test");
  }
  virtual Eigen::MatrixXcd PartialTrace(const Types::qubits_vector &) const {
    throw std::runtime_error("This simulator does not support partial trace");
  }
  virtual double FidelityWithStatevector(const Eigen::VectorXcd &) const {
    throw std::runtime_error("This simulator does not support mixed-state fidelity");
  }
  virtual void RestoreDensityMatrixTrace() {
    throw std::runtime_error("This simulator cannot restore density-matrix trace");
  }
  virtual void HermitizeDensityMatrix() {
    throw std::runtime_error("This simulator cannot hermitize its density matrix");
  }
  virtual void TrimMatrixProductOperator() {
    throw std::runtime_error("This simulator is not a matrix-product operator");
  }
  virtual void ReCanonicalizeMatrixProductOperator() {
    throw std::runtime_error("This simulator is not a matrix-product operator");
  }

  void ApplyKrausChannel(
      const Types::qubits_vector &targets,
      const QuantumChannel::KrausOperators &krausOperators) {
    ApplyQuantumChannel(targets, QuantumChannel(krausOperators));
  }

  void ApplyChannel(
      const Types::qubits_vector &targets,
      const QuantumChannel::KrausOperators &krausOperators) {
    ApplyKrausChannel(targets, krausOperators);
  }

  void ApplyPauliChannel(const Types::qubits_vector &targets,
                         const std::vector<double> &probabilities) {
    ApplyQuantumChannel(targets, QuantumChannel::Pauli(probabilities));
  }

  void ApplyPauliChannel(Types::qubit_t qubit, double px, double py,
                         double pz) {
    ApplyQuantumChannel({qubit}, QuantumChannel::Pauli(px, py, pz));
  }
 
  void ApplyBitFlipNoise(Types::qubit_t qubit, double probability) {
    ApplyQuantumChannel({qubit}, QuantumChannel::BitFlip(probability));
  }
 
  void ApplyBitPhaseFlipNoise(Types::qubit_t qubit, double probability) {
    ApplyQuantumChannel({qubit},
                        QuantumChannel::BitPhaseFlip(probability));
  }
 
  void ApplyPhaseFlipNoise(Types::qubit_t qubit, double probability) {
    ApplyQuantumChannel({qubit}, QuantumChannel::PhaseFlip(probability));
  }

  void ApplyDephasingNoise(Types::qubit_t qubit, double probability) {
    ApplyPhaseFlipNoise(qubit, probability);
  }

  void ApplyDepolarizingNoise(Types::qubit_t qubit,
                              double errorProbability) {
    ApplyQuantumChannel({qubit},
                        QuantumChannel::Depolarizing(errorProbability));
  }

  void ApplyDepolarizingMixingNoise(Types::qubit_t qubit,
                                    double mixingProbability) {
    ApplyQuantumChannel(
        {qubit}, QuantumChannel::DepolarizingMixing(mixingProbability));
  }
 
  void ApplyAmplitudeDamping(Types::qubit_t qubit, double gamma) {
    ApplyQuantumChannel({qubit}, QuantumChannel::AmplitudeDamping(gamma));
  }

  void ApplyT1Relaxation(Types::qubit_t qubit, double gamma) {
    ApplyAmplitudeDamping(qubit, gamma);
  }

  void ApplyT1RelaxationFromTime(Types::qubit_t qubit, double duration,
                                 double t1) {
    if (!std::isfinite(duration) || duration < 0.0)
      throw std::invalid_argument(
          "T1-relaxation duration must be finite and nonnegative");
    if (std::isnan(t1) || t1 <= 0.0)
      throw std::invalid_argument(
          "T1 must be positive (infinity is allowed)");
    const double gamma =
        std::isinf(t1) ? 0.0 : -std::expm1(-duration / t1);
    ApplyAmplitudeDamping(qubit, gamma);
  }

  void ApplyPhaseDamping(Types::qubit_t qubit, double gamma) {
    ApplyQuantumChannel({qubit}, QuantumChannel::PhaseDamping(gamma));
  }

  void ApplyPhaseDampingFromTime(Types::qubit_t qubit, double duration,
                                 double tPhi) {
    if (!std::isfinite(duration) || duration < 0.0)
      throw std::invalid_argument(
          "Phase-damping duration must be finite and nonnegative");
    if (std::isnan(tPhi) || tPhi <= 0.0)
      throw std::invalid_argument(
          "T_phi must be positive (infinity is allowed)");
    const double gamma =
        std::isinf(tPhi) ? 0.0 : -std::expm1(-2.0 * duration / tPhi);
    ApplyPhaseDamping(qubit, gamma);
  }
 
  void ApplyGeneralizedAmplitudeDamping(Types::qubit_t qubit, double gamma,
                                        double excitedStatePopulation) {
    ApplyQuantumChannel(
        {qubit}, QuantumChannel::GeneralizedAmplitudeDamping(
                      gamma, excitedStatePopulation));
  }
 
  void ApplyThermalRelaxation(Types::qubit_t qubit, double duration,
                              double t1, double t2,
                              double excitedStatePopulation = 0.0) {
    ApplyQuantumChannel(
        {qubit}, QuantumChannel::ThermalRelaxation(
                      duration, t1, t2, excitedStatePopulation));
  }
 
  void ApplyCorrelatedPhaseFlipNoise(Types::qubit_t qubit0,
                                     Types::qubit_t qubit1,
                                     double probability) {
    ApplyQuantumChannel(
        {qubit0, qubit1},
        QuantumChannel::CorrelatedPhaseFlip(probability));
  }
 
  void ApplyCorrelatedPhaseFlipNoise(Types::qubit_t qubit0,
                                     Types::qubit_t qubit1,
                                     double probability,
                                     double correlation) {
    ApplyQuantumChannel(
        {qubit0, qubit1},
        QuantumChannel::CorrelatedPhaseFlip(probability, correlation));
  }
 
  void ApplyTwoQubitDepolarizingNoise(Types::qubit_t qubit0,
                                      Types::qubit_t qubit1,
                                      double errorProbability) {
    ApplyQuantumChannel(
        {qubit0, qubit1},
        QuantumChannel::TwoQubitDepolarizing(errorProbability));
  }
 
  void ApplyTwoQubitDepolarizingMixingNoise(
      Types::qubit_t qubit0, Types::qubit_t qubit1,
      double mixingProbability) {
    ApplyQuantumChannel(
        {qubit0, qubit1},
        QuantumChannel::TwoQubitDepolarizingMixing(mixingProbability));
  }

  virtual double Probability(Types::qubit_t outcome) = 0;

  virtual std::complex<double> Amplitude(Types::qubit_t outcome) = 0;

  virtual std::complex<double> ProjectOnZero() = 0;

  virtual std::vector<double> AllProbabilities() = 0;

  virtual std::vector<double> Probabilities(
      const Types::qubits_vector &qubits) = 0;

  virtual std::unordered_map<Types::qubit_t, Types::qubit_t> SampleCounts(
      const Types::qubits_vector &qubits, size_t shots = 1000) = 0;

  virtual std::unordered_map<std::vector<bool>, Types::qubit_t>
  SampleCountsMany(const Types::qubits_vector &qubits, size_t shots = 1000) = 0;

  virtual double ExpectationValue(const std::string &pauliString) = 0;

  void RegisterObserver(const std::shared_ptr<ISimulatorObserver> &observer) {
    observers.insert(observer);
  }

  void UnregisterObserver(const std::shared_ptr<ISimulatorObserver> &observer) {
    observers.erase(observer);
  }

  void ClearObservers() { observers.clear(); }

  virtual SimulatorType GetType() const = 0;

  virtual SimulationType GetSimulationType() const = 0;

  virtual void Flush() = 0;

  virtual void SaveStateToInternalDestructive() = 0;

  virtual void RestoreInternalDestructiveSavedState() = 0;

  virtual void SaveState() = 0;

  virtual void RestoreState() = 0;

  virtual std::complex<double> AmplitudeRaw(Types::qubit_t outcome) = 0;

  virtual void SetMultithreading(bool multithreading = true) = 0;

  virtual bool GetMultithreading() const = 0;

  virtual bool IsQcsim() const = 0;

  virtual Types::qubit_t MeasureNoCollapse() = 0;

  virtual std::vector<bool> MeasureNoCollapseMany() = 0;

  virtual size_t GetCurrentMaxBondDimension() const { return 0; }
 
  virtual const std::unordered_map<std::string, std::string>& GetConfigMap() const = 0;
 
 protected:
  void DontNotify() { notifyObservers = false; }

  void Notify() { notifyObservers = true; }

  void NotifyObservers(const Types::qubits_vector &affectedQubits) {
    if (!notifyObservers) return;
 
    for (auto &observer : observers) {
      observer->Update(affectedQubits);
    }
  }
 
 
 private:
  std::unordered_set<std::shared_ptr<ISimulatorObserver>>
      observers; 
  bool notifyObservers =
      true; 
};
 
}  // namespace Simulators
 
#endif  // !_SIMULATOR_STATE_H_