Qrack  9.13
General classical-emulating-quantum development framework
qstabilizerhybrid.hpp
Go to the documentation of this file.
1 //
3 // (C) Daniel Strano and the Qrack contributors 2017-2023. All rights reserved.
4 //
5 // This is a multithreaded, universal quantum register simulation, allowing
6 // (nonphysical) register cloning and direct measurement of probability and
7 // phase, to leverage what advantages classical emulation of qubits can have.
8 //
9 // Licensed under the GNU Lesser General Public License V3.
10 // See LICENSE.md in the project root or https://www.gnu.org/licenses/lgpl-3.0.en.html
11 // for details.
12 #pragma once
13 
14 #include "mpsshard.hpp"
15 #include "qengine.hpp"
16 #include "qunitclifford.hpp"
17 
18 #define QINTERFACE_TO_QALU(qReg) std::dynamic_pointer_cast<QAlu>(qReg)
19 #define QINTERFACE_TO_QPARITY(qReg) std::dynamic_pointer_cast<QParity>(qReg)
20 
21 namespace Qrack {
22 
26 
28  : amp(a)
29  , stabilizer(s)
30  {
31  // Intentionally left blank.
32  }
33 };
34 
35 class QStabilizerHybrid;
36 typedef std::shared_ptr<QStabilizerHybrid> QStabilizerHybridPtr;
37 
42 #if ENABLE_ALU
43 class QStabilizerHybrid : public QAlu, public QParity, public QInterface {
44 #else
45 class QStabilizerHybrid : public QParity, public QInterface {
46 #endif
47 protected:
48  bool useHostRam;
50  bool useTGadget;
60  int64_t devID;
62  double logFidelity;
65  std::vector<int64_t> deviceIDs;
66  std::vector<QInterfaceEngine> engineTypes;
67  std::vector<QInterfaceEngine> cloneEngineTypes;
68  std::vector<MpsShardPtr> shards;
70 
73  QInterfacePtr MakeEngine(const bitCapInt& perm, bitLenInt qbCount);
74 
75  void InvertBuffer(bitLenInt qubit);
76  void FlushH(bitLenInt qubit);
77  void FlushIfBlocked(bitLenInt control, bitLenInt target, bool isPhase = false);
79  bool TrimControls(const std::vector<bitLenInt>& lControls, std::vector<bitLenInt>& output, bool anti = false);
80  void CacheEigenstate(bitLenInt target);
81  void FlushBuffers();
82  void DumpBuffers()
83  {
84  for (size_t i = 0U; i < shards.size(); ++i) {
85  shards[i] = NULL;
86  }
87  }
88  bool EitherIsBuffered(bool logical)
89  {
90  const size_t maxLcv = logical ? (size_t)qubitCount : shards.size();
91  for (size_t i = 0U; i < maxLcv; ++i) {
92  if (shards[i]) {
93  // We have a cached non-Clifford operation.
94  return true;
95  }
96  }
97 
98  return false;
99  }
100  bool IsBuffered() { return EitherIsBuffered(false); }
101  bool IsLogicalBuffered() { return EitherIsBuffered(true); }
102  bool EitherIsProbBuffered(bool logical)
103  {
104  const size_t maxLcv = logical ? (size_t)qubitCount : shards.size();
105  for (size_t i = 0U; i < maxLcv; ++i) {
106  MpsShardPtr shard = shards[i];
107  if (!shard) {
108  continue;
109  }
110  if (shard->IsHPhase() || shard->IsHInvert()) {
111  FlushH(i);
112  }
113  if (shard->IsInvert()) {
114  InvertBuffer(i);
115  }
116  if (!shard->IsPhase()) {
117  // We have a cached non-Clifford operation.
118  return true;
119  }
120  }
121 
122  return false;
123  }
124  bool IsProbBuffered() { return EitherIsProbBuffered(false); }
126 
127  std::unique_ptr<complex[]> GetQubitReducedDensityMatrix(bitLenInt qubit)
128  {
129  // Form the reduced density matrix of the single qubit.
130  const real1 z = (real1)(ONE_R1_F - 2 * stabilizer->Prob(qubit));
131  stabilizer->H(qubit);
132  const real1 x = (real1)(ONE_R1_F - 2 * stabilizer->Prob(qubit));
133  stabilizer->S(qubit);
134  const real1 y = (real1)(ONE_R1_F - 2 * stabilizer->Prob(qubit));
135  stabilizer->IS(qubit);
136  stabilizer->H(qubit);
137 
138  std::unique_ptr<complex[]> dMtrx(new complex[4]);
139  dMtrx[0] = (ONE_CMPLX + z) / complex((real1)2, ZERO_R1);
140  dMtrx[1] = x / complex((real1)2, ZERO_R1) - I_CMPLX * (y / complex((real1)2, ZERO_R1));
141  dMtrx[2] = x / complex((real1)2, ZERO_R1) + I_CMPLX * (y / complex((real1)2, ZERO_R1));
142  dMtrx[3] = (ONE_CMPLX + z) / complex((real1)2, ZERO_R1);
143  if (shards[qubit]) {
144  const complex adj[4]{ std::conj(shards[qubit]->gate[0]), std::conj(shards[qubit]->gate[2]),
145  std::conj(shards[qubit]->gate[1]), std::conj(shards[qubit]->gate[3]) };
146  complex out[4];
147  mul2x2(dMtrx.get(), adj, out);
148  mul2x2(shards[qubit]->gate, out, dMtrx.get());
149  }
150 
151  return dMtrx;
152  }
153 
154  template <typename F>
155  void CheckShots(unsigned shots, const bitCapInt& m, real1_f partProb, const std::vector<bitCapInt>& qPowers,
156  std::vector<real1_f>& rng, F fn)
157  {
158  for (int64_t shot = rng.size() - 1U; shot >= 0; --shot) {
159  if (rng[shot] >= partProb) {
160  break;
161  }
162 
163  bitCapInt sample = ZERO_BCI;
164  for (size_t i = 0U; i < qPowers.size(); ++i) {
165  if (bi_compare_0(m & qPowers[i]) != 0) {
166  bi_or_ip(&sample, pow2(i));
167  }
168  }
169  fn(sample, (unsigned int)shot);
170 
171  rng.erase(rng.begin() + shot);
172  if (rng.empty()) {
173  break;
174  }
175  }
176  }
177 
178  std::vector<real1_f> GenerateShotProbs(unsigned shots)
179  {
180  std::vector<real1_f> rng;
181  rng.reserve(shots);
182  for (unsigned shot = 0U; shot < shots; ++shot) {
183  rng.push_back(Rand());
184  }
185  std::sort(rng.begin(), rng.end());
186  std::reverse(rng.begin(), rng.end());
187 
188  return rng;
189  }
190 
191  real1_f FractionalRzAngleWithFlush(bitLenInt i, real1_f angle, bool isGateSuppressed = false)
192  {
193  const real1_f sectorAngle = PI_R1 / 2;
194  const real1_f Period = 2 * PI_R1;
195  while (angle >= Period) {
196  angle -= Period;
197  }
198  while (angle < 0U) {
199  angle += Period;
200  }
201 
202  const long sector = std::lround((real1_s)(angle / sectorAngle));
203  if (!isGateSuppressed) {
204  switch (sector) {
205  case 1:
206  stabilizer->S(i);
207  break;
208  case 2:
209  stabilizer->Z(i);
210  break;
211  case 3:
212  stabilizer->IS(i);
213  break;
214  case 0:
215  default:
216  break;
217  }
218  }
219 
220  angle -= (sector * sectorAngle);
221  if (angle > PI_R1) {
222  angle -= Period;
223  }
224  if (angle <= -PI_R1) {
225  angle += Period;
226  }
227 
228  return angle;
229  }
230 
232  {
233  for (size_t i = 0U; i < qubitCount; ++i) {
234  // Flush all buffers as close as possible to Clifford.
235  const MpsShardPtr& shard = shards[i];
236  if (!shard) {
237  continue;
238  }
239  if (shard->IsHPhase() || shard->IsHInvert()) {
240  FlushH(i);
241  }
242  if (shard->IsInvert()) {
243  InvertBuffer(i);
244  }
245  if (!shard->IsPhase()) {
246  // We have a cached non-phase operation.
247  continue;
248  }
249  const real1 angle = (real1)(FractionalRzAngleWithFlush(i, std::arg(shard->gate[3U] / shard->gate[0U])) / 2);
250  if ((2 * abs(angle) / PI_R1) <= FP_NORM_EPSILON) {
251  shards[i] = NULL;
252  continue;
253  }
254  const real1 angleCos = cos(angle);
255  const real1 angleSin = sin(angle);
256  shard->gate[0U] = complex(angleCos, -angleSin);
257  shard->gate[3U] = complex(angleCos, angleSin);
258  }
259 
260  RdmCloneFlush();
261  }
262 
264  {
265  QStabilizerHybridPtr clone = std::dynamic_pointer_cast<QStabilizerHybrid>(Clone());
266  clone->RdmCloneFlush(HALF_R1);
267 
268  return clone;
269  }
270  void RdmCloneFlush(real1_f threshold = FP_NORM_EPSILON);
271 
272  real1_f ExpVarFactorized(bool isExp, bool isFloat, const std::vector<bitLenInt>& bits,
273  const std::vector<bitCapInt>& perms, const std::vector<real1_f>& weights, const bitCapInt& offset, bool roundRz)
274  {
275  if (engine) {
276  return isExp ? isFloat ? engine->ExpectationFloatsFactorizedRdm(roundRz, bits, weights)
277  : engine->ExpectationBitsFactorizedRdm(roundRz, bits, perms, offset)
278  : isFloat ? engine->VarianceFloatsFactorizedRdm(roundRz, bits, weights)
279  : engine->VarianceBitsFactorizedRdm(roundRz, bits, perms, offset);
280  }
281 
282  if (!roundRz) {
283  return isExp ? isFloat ? stabilizer->ExpectationFloatsFactorizedRdm(roundRz, bits, weights)
284  : stabilizer->ExpectationBitsFactorizedRdm(roundRz, bits, perms, offset)
285  : isFloat ? stabilizer->VarianceFloatsFactorizedRdm(roundRz, bits, weights)
286  : stabilizer->VarianceBitsFactorizedRdm(roundRz, bits, perms, offset);
287  }
288 
289  return isExp ? isFloat
290  ? RdmCloneHelper()->stabilizer->ExpectationFloatsFactorizedRdm(roundRz, bits, weights)
291  : RdmCloneHelper()->stabilizer->ExpectationBitsFactorizedRdm(roundRz, bits, perms, offset)
292  : isFloat ? RdmCloneHelper()->stabilizer->VarianceFloatsFactorizedRdm(roundRz, bits, weights)
293  : RdmCloneHelper()->stabilizer->VarianceBitsFactorizedRdm(roundRz, bits, perms, offset);
294  }
295 
297  {
298  if (stabilizer->TrySeparate(i)) {
299  stabilizer->Dispose(i, 1U);
300  shards.erase(shards.begin() + i);
301  } else {
302  const bitLenInt deadIndex = qubitCount + ancillaCount - 1U;
303  stabilizer->SetBit(i, false);
304  if (i != deadIndex) {
305  stabilizer->Swap(i, deadIndex);
306  shards[i].swap(shards[deadIndex]);
307  }
308  shards.erase(shards.begin() + deadIndex);
310  }
311  --ancillaCount;
312  }
313 
315  QStabilizerHybridPtr toCompare, bool isDiscreteBool, real1_f error_tol = TRYDECOMPOSE_EPSILON);
316 
317  void ISwapHelper(bitLenInt qubit1, bitLenInt qubit2, bool inverse);
318 
319  complex GetAmplitudeOrProb(const bitCapInt& perm, bool isProb = false);
320 
321  QInterfacePtr CloneBody(bool isCopy);
322  using QInterface::Copy;
323  void Copy(QInterfacePtr orig) { Copy(std::dynamic_pointer_cast<QStabilizerHybrid>(orig)); }
325  {
326  QInterface::Copy(std::dynamic_pointer_cast<QInterface>(orig));
327  useHostRam = orig->useHostRam;
328  doNormalize = orig->doNormalize;
329  useTGadget = orig->useTGadget;
330  isRoundingFlushed = orig->isRoundingFlushed;
331  thresholdQubits = orig->thresholdQubits;
332  ancillaCount = orig->ancillaCount;
333  deadAncillaCount = orig->deadAncillaCount;
334  maxEngineQubitCount = orig->maxEngineQubitCount;
335  maxAncillaCount = orig->maxAncillaCount;
336  maxStateMapCacheQubitCount = orig->maxStateMapCacheQubitCount;
337  separabilityThreshold = orig->separabilityThreshold;
338  roundingThreshold = orig->roundingThreshold;
339  devID = orig->devID;
340  phaseFactor = orig->phaseFactor;
341  logFidelity = orig->logFidelity;
342  engine = orig->engine;
343  stabilizer = orig->stabilizer;
344  deviceIDs = orig->deviceIDs;
345  engineTypes = orig->engineTypes;
346  cloneEngineTypes = orig->cloneEngineTypes;
347  shards = orig->shards;
348  stateMapCache = orig->stateMapCache;
349  }
350 
351 public:
352  QStabilizerHybrid(std::vector<QInterfaceEngine> eng, bitLenInt qBitCount, const bitCapInt& initState = ZERO_BCI,
353  qrack_rand_gen_ptr rgp = nullptr, const complex& phaseFac = CMPLX_DEFAULT_ARG, bool doNorm = false,
354  bool randomGlobalPhase = true, bool useHostMem = false, int64_t deviceId = -1, bool useHardwareRNG = true,
355  bool ignored = false, real1_f norm_thresh = REAL1_EPSILON, std::vector<int64_t> devList = {},
356  bitLenInt qubitThreshold = 0U, real1_f separation_thresh = _qrack_qunit_sep_thresh);
357 
358  QStabilizerHybrid(bitLenInt qBitCount, const bitCapInt& initState = ZERO_BCI, qrack_rand_gen_ptr rgp = nullptr,
359  const complex& phaseFac = CMPLX_DEFAULT_ARG, bool doNorm = false, bool randomGlobalPhase = true,
360  bool useHostMem = false, int64_t deviceId = -1, bool useHardwareRNG = true, bool ignored = false,
361  real1_f norm_thresh = REAL1_EPSILON, std::vector<int64_t> devList = {}, bitLenInt qubitThreshold = 0U,
362  real1_f separation_thresh = _qrack_qunit_sep_thresh)
363  : QStabilizerHybrid({ QINTERFACE_HYBRID }, qBitCount, initState, rgp, phaseFac, doNorm, randomGlobalPhase,
364  useHostMem, deviceId, useHardwareRNG, ignored, norm_thresh, devList, qubitThreshold, separation_thresh)
365  {
366  }
367 
368  void SetNcrp(real1_f ncrp) { roundingThreshold = ncrp; };
369  void SetTInjection(bool useGadget) { useTGadget = useGadget; }
370  bool GetTInjection() { return useTGadget; }
371  double GetUnitaryFidelity() { return exp(logFidelity); }
373 
374  void Finish()
375  {
376  if (stabilizer) {
377  stabilizer->Finish();
378  } else {
379  engine->Finish();
380  }
381  };
382 
383  bool isFinished() { return (!stabilizer || stabilizer->isFinished()) && (!engine || engine->isFinished()); }
384 
385  void Dump()
386  {
387  if (stabilizer) {
388  stabilizer->Dump();
389  } else {
390  engine->Dump();
391  }
392  }
393 
394  void SetConcurrency(uint32_t threadCount)
395  {
396  QInterface::SetConcurrency(threadCount);
397  if (engine) {
399  }
400  }
401 
403  {
404  if (!ancillaCount || stabilizer->IsSeparable(qubit)) {
405  return Prob(qubit);
406  }
407 
408  std::unique_ptr<complex[]> dMtrx = GetQubitReducedDensityMatrix(qubit);
409  QRACK_CONST complex ONE_CMPLX_NEG = complex(-ONE_R1, ZERO_R1);
410  QRACK_CONST complex pauliZ[4]{ ONE_CMPLX, ZERO_CMPLX, ZERO_CMPLX, ONE_CMPLX_NEG };
411  complex pMtrx[4];
412  mul2x2(dMtrx.get(), pauliZ, pMtrx);
413  return (ONE_R1 - std::real(pMtrx[0]) + std::real(pMtrx[1])) / 2;
414  }
415 
417  {
418  AntiCNOT(control, target);
419  const real1_f prob = ProbRdm(target);
420  AntiCNOT(control, target);
421 
422  return prob;
423  }
424 
426  {
427  CNOT(control, target);
428  const real1_f prob = ProbRdm(target);
429  CNOT(control, target);
430 
431  return prob;
432  }
433 
439  void SwitchToEngine();
440 
441  bool isClifford() { return !engine; }
442 
443  bool isClifford(bitLenInt qubit) { return !engine && !shards[qubit]; };
444 
445  bool isBinaryDecisionTree() { return engine && engine->isBinaryDecisionTree(); };
446 
447  using QInterface::Compose;
448  bitLenInt Compose(QStabilizerHybridPtr toCopy) { return ComposeEither(toCopy, false); };
449  bitLenInt Compose(QInterfacePtr toCopy) { return Compose(std::dynamic_pointer_cast<QStabilizerHybrid>(toCopy)); }
452  {
453  return Compose(std::dynamic_pointer_cast<QStabilizerHybrid>(toCopy), start);
454  }
455  bitLenInt ComposeNoClone(QStabilizerHybridPtr toCopy) { return ComposeEither(toCopy, true); };
457  {
458  return ComposeNoClone(std::dynamic_pointer_cast<QStabilizerHybrid>(toCopy));
459  }
460  bitLenInt ComposeEither(QStabilizerHybridPtr toCopy, bool willDestroy);
462  {
463  Decompose(start, std::dynamic_pointer_cast<QStabilizerHybrid>(dest));
464  }
465  void Decompose(bitLenInt start, QStabilizerHybridPtr dest);
467  void Dispose(bitLenInt start, bitLenInt length);
468  void Dispose(bitLenInt start, bitLenInt length, const bitCapInt& disposedPerm);
469  using QInterface::Allocate;
470  bitLenInt Allocate(bitLenInt start, bitLenInt length);
471 
472  void GetQuantumState(complex* outputState);
473  void GetProbs(real1* outputProbs);
474  complex GetAmplitude(const bitCapInt& perm) { return GetAmplitudeOrProb(perm, false); }
475  real1_f ProbAll(const bitCapInt& perm) { return (real1_f)norm(GetAmplitudeOrProb(perm, true)); }
476  void SetQuantumState(const complex* inputState);
477  void SetAmplitude(const bitCapInt& perm, const complex& amp)
478  {
479  SwitchToEngine();
480  engine->SetAmplitude(perm, amp);
481  }
482  void SetPermutation(const bitCapInt& perm, const complex& phaseFac = CMPLX_DEFAULT_ARG);
483 
484  void Swap(bitLenInt qubit1, bitLenInt qubit2);
485  void ISwap(bitLenInt qubit1, bitLenInt qubit2) { ISwapHelper(qubit1, qubit2, false); }
486  void IISwap(bitLenInt qubit1, bitLenInt qubit2) { ISwapHelper(qubit1, qubit2, true); }
487  void CSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
488  void CSqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
489  void AntiCSqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
490  void CISqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
491  void AntiCISqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
492 
493  void XMask(const bitCapInt& mask);
494  void YMask(const bitCapInt& mask);
495  void ZMask(const bitCapInt& mask);
496 
497  real1_f Prob(bitLenInt qubit);
498 
499  bool ForceM(bitLenInt qubit, bool result, bool doForce = true, bool doApply = true);
500 
501  bitCapInt MAll();
502 
503  void Mtrx(const complex* mtrx, bitLenInt target);
504  void MCMtrx(const std::vector<bitLenInt>& controls, const complex* mtrx, bitLenInt target);
505  void MCPhase(
506  const std::vector<bitLenInt>& controls, const complex& topLeft, const complex& bottomRight, bitLenInt target);
507  void MCInvert(
508  const std::vector<bitLenInt>& controls, const complex& topRight, const complex& bottomLeft, bitLenInt target);
509  void MACMtrx(const std::vector<bitLenInt>& controls, const complex* mtrx, bitLenInt target);
510  void MACPhase(
511  const std::vector<bitLenInt>& controls, const complex& topLeft, const complex& bottomRight, bitLenInt target);
512  void MACInvert(
513  const std::vector<bitLenInt>& controls, const complex& topRight, const complex& bottomLeft, bitLenInt target);
514 
517  const std::vector<bitLenInt>& controls, bitLenInt qubitIndex, const complex* mtrxs);
518 
519  std::map<bitCapInt, int> MultiShotMeasureMask(const std::vector<bitCapInt>& qPowers, unsigned shots);
520  void MultiShotMeasureMask(const std::vector<bitCapInt>& qPowers, unsigned shots, unsigned long long* shotsArray);
521 
522  real1_f ProbParity(const bitCapInt& mask);
523  bool ForceMParity(const bitCapInt& mask, bool result, bool doForce = true);
524  void CUniformParityRZ(const std::vector<bitLenInt>& controls, const bitCapInt& mask, real1_f angle)
525  {
526  SwitchToEngine();
527  QINTERFACE_TO_QPARITY(engine)->CUniformParityRZ(controls, mask, angle);
528  }
529 
530 #if ENABLE_ALU
531  using QInterface::M;
532  bool M(bitLenInt q) { return QInterface::M(q); }
533  using QInterface::X;
534  void X(bitLenInt q) { QInterface::X(q); }
535  void CPhaseFlipIfLess(const bitCapInt& greaterPerm, bitLenInt start, bitLenInt length, bitLenInt flagIndex)
536  {
537  SwitchToEngine();
538  QINTERFACE_TO_QALU(engine)->CPhaseFlipIfLess(greaterPerm, start, length, flagIndex);
539  }
540  void PhaseFlipIfLess(const bitCapInt& greaterPerm, bitLenInt start, bitLenInt length)
541  {
542  SwitchToEngine();
543  QINTERFACE_TO_QALU(engine)->PhaseFlipIfLess(greaterPerm, start, length);
544  }
545 
546  void INC(const bitCapInt& toAdd, bitLenInt start, bitLenInt length)
547  {
548  if (stabilizer) {
549  QInterface::INC(toAdd, start, length);
550  return;
551  }
552 
553  engine->INC(toAdd, start, length);
554  }
555  void DEC(const bitCapInt& toSub, bitLenInt start, bitLenInt length)
556  {
557  if (stabilizer) {
558  QInterface::DEC(toSub, start, length);
559  return;
560  }
561 
562  engine->DEC(toSub, start, length);
563  }
564  void DECS(const bitCapInt& toSub, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
565  {
566  if (stabilizer) {
567  QInterface::DECS(toSub, start, length, overflowIndex);
568  return;
569  }
570 
571  engine->DECS(toSub, start, length, overflowIndex);
572  }
573  void CINC(const bitCapInt& toAdd, bitLenInt inOutStart, bitLenInt length, const std::vector<bitLenInt>& controls)
574  {
575  if (stabilizer) {
576  QInterface::CINC(toAdd, inOutStart, length, controls);
577  return;
578  }
579 
580  engine->CINC(toAdd, inOutStart, length, controls);
581  }
582  void INCS(const bitCapInt& toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
583  {
584  if (stabilizer) {
585  QInterface::INCS(toAdd, start, length, overflowIndex);
586  return;
587  }
588 
589  engine->INCS(toAdd, start, length, overflowIndex);
590  }
591  void INCDECC(const bitCapInt& toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
592  {
593  if (stabilizer) {
594  QInterface::INCDECC(toAdd, start, length, carryIndex);
595  return;
596  }
597 
598  engine->INCDECC(toAdd, start, length, carryIndex);
599  }
600  void INCDECSC(
601  const bitCapInt& toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex, bitLenInt carryIndex)
602  {
603  SwitchToEngine();
604  QINTERFACE_TO_QALU(engine)->INCDECSC(toAdd, start, length, overflowIndex, carryIndex);
605  }
606  void INCDECSC(const bitCapInt& toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
607  {
608  SwitchToEngine();
609  QINTERFACE_TO_QALU(engine)->INCDECSC(toAdd, start, length, carryIndex);
610  }
611 #if ENABLE_BCD
612  void INCBCD(const bitCapInt& toAdd, bitLenInt start, bitLenInt length)
613  {
614  SwitchToEngine();
615  QINTERFACE_TO_QALU(engine)->INCBCD(toAdd, start, length);
616  }
617  void INCDECBCDC(const bitCapInt& toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
618  {
619  SwitchToEngine();
620  QINTERFACE_TO_QALU(engine)->INCDECBCDC(toAdd, start, length, carryIndex);
621  }
622 #endif
623  void MUL(const bitCapInt& toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
624  {
625  SwitchToEngine();
626  QINTERFACE_TO_QALU(engine)->MUL(toMul, inOutStart, carryStart, length);
627  }
628  void DIV(const bitCapInt& toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
629  {
630  SwitchToEngine();
631  QINTERFACE_TO_QALU(engine)->DIV(toDiv, inOutStart, carryStart, length);
632  }
634  const bitCapInt& toMul, const bitCapInt& modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
635  {
636  SwitchToEngine();
637  QINTERFACE_TO_QALU(engine)->MULModNOut(toMul, modN, inStart, outStart, length);
638  }
640  const bitCapInt& toMul, const bitCapInt& modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
641  {
642  SwitchToEngine();
643  QINTERFACE_TO_QALU(engine)->IMULModNOut(toMul, modN, inStart, outStart, length);
644  }
646  const bitCapInt& base, const bitCapInt& modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
647  {
648  SwitchToEngine();
649  QINTERFACE_TO_QALU(engine)->POWModNOut(base, modN, inStart, outStart, length);
650  }
651  void CMUL(const bitCapInt& toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length,
652  const std::vector<bitLenInt>& controls)
653  {
654  SwitchToEngine();
655  QINTERFACE_TO_QALU(engine)->CMUL(toMul, inOutStart, carryStart, length, controls);
656  }
657  void CDIV(const bitCapInt& toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length,
658  const std::vector<bitLenInt>& controls)
659  {
660  SwitchToEngine();
661  QINTERFACE_TO_QALU(engine)->CDIV(toDiv, inOutStart, carryStart, length, controls);
662  }
663  void CMULModNOut(const bitCapInt& toMul, const bitCapInt& modN, bitLenInt inStart, bitLenInt outStart,
664  bitLenInt length, const std::vector<bitLenInt>& controls)
665  {
666  SwitchToEngine();
667  QINTERFACE_TO_QALU(engine)->CMULModNOut(toMul, modN, inStart, outStart, length, controls);
668  }
669  void CIMULModNOut(const bitCapInt& toMul, const bitCapInt& modN, bitLenInt inStart, bitLenInt outStart,
670  bitLenInt length, const std::vector<bitLenInt>& controls)
671  {
672  SwitchToEngine();
673  QINTERFACE_TO_QALU(engine)->CIMULModNOut(toMul, modN, inStart, outStart, length, controls);
674  }
675  void CPOWModNOut(const bitCapInt& base, const bitCapInt& modN, bitLenInt inStart, bitLenInt outStart,
676  bitLenInt length, const std::vector<bitLenInt>& controls)
677  {
678  SwitchToEngine();
679  QINTERFACE_TO_QALU(engine)->CPOWModNOut(base, modN, inStart, outStart, length, controls);
680  }
681 
682  bitCapInt IndexedLDA(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength,
683  const unsigned char* values, bool resetValue = true)
684  {
685  SwitchToEngine();
686  return QINTERFACE_TO_QALU(engine)->IndexedLDA(
687  indexStart, indexLength, valueStart, valueLength, values, resetValue);
688  }
689  bitCapInt IndexedADC(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength,
690  bitLenInt carryIndex, const unsigned char* values)
691  {
692  SwitchToEngine();
693  return QINTERFACE_TO_QALU(engine)->IndexedADC(
694  indexStart, indexLength, valueStart, valueLength, carryIndex, values);
695  }
696  bitCapInt IndexedSBC(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength,
697  bitLenInt carryIndex, const unsigned char* values)
698  {
699  SwitchToEngine();
700  return QINTERFACE_TO_QALU(engine)->IndexedSBC(
701  indexStart, indexLength, valueStart, valueLength, carryIndex, values);
702  }
703  void Hash(bitLenInt start, bitLenInt length, const unsigned char* values)
704  {
705  SwitchToEngine();
706  QINTERFACE_TO_QALU(engine)->Hash(start, length, values);
707  }
708 #endif
709 
710  void PhaseFlip()
711  {
712  if (stabilizer) {
713  stabilizer->PhaseFlip();
714  } else {
715  engine->PhaseFlip();
716  }
717  }
718  void ZeroPhaseFlip(bitLenInt start, bitLenInt length)
719  {
720  SwitchToEngine();
721  engine->ZeroPhaseFlip(start, length);
722  }
723 
724  void SqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
725  {
726  if (stabilizer) {
727  QInterface::SqrtSwap(qubitIndex1, qubitIndex2);
728  return;
729  }
730 
731  SwitchToEngine();
732  engine->SqrtSwap(qubitIndex1, qubitIndex2);
733  }
734  void ISqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
735  {
736  if (stabilizer) {
737  QInterface::ISqrtSwap(qubitIndex1, qubitIndex2);
738  return;
739  }
740 
741  SwitchToEngine();
742  engine->ISqrtSwap(qubitIndex1, qubitIndex2);
743  }
744  void FSim(real1_f theta, real1_f phi, bitLenInt qubit1, bitLenInt qubit2)
745  {
746  const std::vector<bitLenInt> controls{ qubit1 };
747  const real1 sinTheta = (real1)sin(theta);
748 
749  if ((sinTheta * sinTheta) <= FP_NORM_EPSILON) {
750  MCPhase(controls, ONE_CMPLX, exp(complex(ZERO_R1, (real1)phi)), qubit2);
751  return;
752  }
753 
754  const real1 sinThetaDiffNeg = ONE_R1 + sinTheta;
755  if ((sinThetaDiffNeg * sinThetaDiffNeg) <= FP_NORM_EPSILON) {
756  ISwap(qubit1, qubit2);
757  MCPhase(controls, ONE_CMPLX, exp(complex(ZERO_R1, (real1)phi)), qubit2);
758  return;
759  }
760 
761  const real1 sinThetaDiffPos = ONE_R1 - sinTheta;
762  if ((sinThetaDiffPos * sinThetaDiffPos) <= FP_NORM_EPSILON) {
763  IISwap(qubit1, qubit2);
764  MCPhase(controls, ONE_CMPLX, exp(complex(ZERO_R1, (real1)phi)), qubit2);
765  return;
766  }
767 
768  SwitchToEngine();
769  engine->FSim(theta, phi, qubit1, qubit2);
770  }
771 
772  real1_f ProbMask(const bitCapInt& mask, const bitCapInt& permutation)
773  {
774  SwitchToEngine();
775  return engine->ProbMask(mask, permutation);
776  }
777 
779  {
780  return ApproxCompareHelper(std::dynamic_pointer_cast<QStabilizerHybrid>(toCompare), false);
781  }
783  {
784  return error_tol >=
785  ApproxCompareHelper(std::dynamic_pointer_cast<QStabilizerHybrid>(toCompare), true, error_tol);
786  }
787 
789  {
790  if (engine) {
791  engine->UpdateRunningNorm(norm_thresh);
792  }
793  }
794 
795  void NormalizeState(
796  real1_f nrm = REAL1_DEFAULT_ARG, real1_f norm_thresh = REAL1_DEFAULT_ARG, real1_f phaseArg = ZERO_R1_F);
797 
798  real1_f ProbAllRdm(bool roundRz, const bitCapInt& fullRegister);
799  real1_f ProbMaskRdm(bool roundRz, const bitCapInt& mask, const bitCapInt& permutation);
800  real1_f ExpectationBitsAll(const std::vector<bitLenInt>& bits, const bitCapInt& offset = ZERO_BCI)
801  {
802  if (stabilizer) {
803  return QInterface::ExpectationBitsAll(bits, offset);
804  }
805 
806  return engine->ExpectationBitsAll(bits, offset);
807  }
808  real1_f ExpectationBitsAllRdm(bool roundRz, const std::vector<bitLenInt>& bits, const bitCapInt& offset = ZERO_BCI)
809  {
810  if (engine) {
811  return engine->ExpectationBitsAllRdm(roundRz, bits, offset);
812  }
813 
814  if (!roundRz) {
815  return stabilizer->ExpectationBitsAll(bits, offset);
816  }
817 
818  return RdmCloneHelper()->stabilizer->ExpectationBitsAll(bits, offset);
819  }
821  const std::vector<bitLenInt>& bits, const std::vector<bitCapInt>& perms, const bitCapInt& offset = ZERO_BCI)
822  {
823  if (stabilizer) {
824  return QInterface::ExpectationBitsFactorized(bits, perms, offset);
825  }
826 
827  return engine->ExpectationBitsFactorized(bits, perms, offset);
828  }
829  real1_f ExpectationBitsFactorizedRdm(bool roundRz, const std::vector<bitLenInt>& bits,
830  const std::vector<bitCapInt>& perms, const bitCapInt& offset = ZERO_BCI)
831  {
832  return ExpVarFactorized(true, false, bits, perms, std::vector<real1_f>(), offset, roundRz);
833  }
834  real1_f ExpectationFloatsFactorized(const std::vector<bitLenInt>& bits, const std::vector<real1_f>& weights)
835  {
836  if (stabilizer) {
837  return QInterface::ExpectationFloatsFactorized(bits, weights);
838  }
839 
840  return engine->ExpectationFloatsFactorized(bits, weights);
841  }
843  bool roundRz, const std::vector<bitLenInt>& bits, const std::vector<real1_f>& weights)
844  {
845  return ExpVarFactorized(true, true, bits, std::vector<bitCapInt>(), weights, ZERO_BCI, roundRz);
846  }
847  real1_f VarianceBitsAll(const std::vector<bitLenInt>& bits, const bitCapInt& offset = ZERO_BCI)
848  {
849  if (stabilizer) {
850  return QInterface::VarianceBitsAll(bits, offset);
851  }
852 
853  return engine->VarianceBitsAll(bits, offset);
854  }
855  real1_f VarianceBitsAllRdm(bool roundRz, const std::vector<bitLenInt>& bits, const bitCapInt& offset = ZERO_BCI)
856  {
857  if (engine) {
858  return engine->VarianceBitsAllRdm(roundRz, bits, offset);
859  }
860 
861  if (!roundRz) {
862  return stabilizer->VarianceBitsAll(bits, offset);
863  }
864 
865  return RdmCloneHelper()->stabilizer->VarianceBitsAll(bits, offset);
866  }
868  const std::vector<bitLenInt>& bits, const std::vector<bitCapInt>& perms, const bitCapInt& offset = ZERO_BCI)
869  {
870  if (stabilizer) {
871  return QInterface::VarianceBitsFactorized(bits, perms, offset);
872  }
873 
874  return engine->VarianceBitsFactorized(bits, perms, offset);
875  }
876  real1_f VarianceBitsFactorizedRdm(bool roundRz, const std::vector<bitLenInt>& bits,
877  const std::vector<bitCapInt>& perms, const bitCapInt& offset = ZERO_BCI)
878  {
879  return ExpVarFactorized(true, false, bits, perms, std::vector<real1_f>(), offset, roundRz);
880  }
881  real1_f VarianceFloatsFactorized(const std::vector<bitLenInt>& bits, const std::vector<real1_f>& weights)
882  {
883  if (stabilizer) {
884  return QInterface::VarianceFloatsFactorized(bits, weights);
885  }
886 
887  return engine->VarianceFloatsFactorized(bits, weights);
888  }
890  bool roundRz, const std::vector<bitLenInt>& bits, const std::vector<real1_f>& weights)
891  {
892  return ExpVarFactorized(true, true, bits, std::vector<bitCapInt>(), weights, ZERO_BCI, roundRz);
893  }
894 
895  bool TrySeparate(bitLenInt qubit);
896  bool TrySeparate(bitLenInt qubit1, bitLenInt qubit2);
897  bool TrySeparate(const std::vector<bitLenInt>& qubits, real1_f error_tol);
898 
899  QInterfacePtr Clone() { return CloneBody(false); }
900  QInterfacePtr Copy() { return CloneBody(true); }
901 
902  void SetDevice(int64_t dID)
903  {
904  devID = dID;
905  if (engine) {
906  engine->SetDevice(dID);
907  }
908  }
909 
910  int64_t GetDeviceID() { return devID; }
911 
913  {
914  if (stabilizer) {
915  return QInterface::GetMaxSize();
916  }
917 
918  return engine->GetMaxSize();
919  }
920 
921  friend std::ostream& operator<<(std::ostream& os, const QStabilizerHybridPtr s);
922  friend std::istream& operator>>(std::istream& is, const QStabilizerHybridPtr s);
923 };
924 } // namespace Qrack
void bi_or_ip(BigInteger *left, const BigInteger &right)
Definition: big_integer.hpp:429
int bi_compare_0(const BigInteger &left)
Definition: big_integer.hpp:134
unsigned GetConcurrencyLevel()
Definition: parallel_for.hpp:41
Definition: qalu.hpp:22
A "Qrack::QInterface" is an abstract interface exposing qubit permutation state vector with methods t...
Definition: qinterface.hpp:141
virtual void SetConcurrency(uint32_t threadsPerEngine)
Set the number of threads in parallel for loops, per component QEngine.
Definition: qinterface.hpp:275
virtual bitLenInt Allocate(bitLenInt length)
Allocate new "length" count of |0> state qubits at end of qubit index position.
Definition: qinterface.hpp:470
virtual bitLenInt Compose(QInterfacePtr toCopy)
Combine another QInterface with this one, after the last bit index of this one.
Definition: qinterface.hpp:364
bitLenInt qubitCount
Definition: qinterface.hpp:146
real1_f Rand()
Generate a random real number between 0 and 1.
Definition: qinterface.hpp:286
Definition: qparity.hpp:22
A "Qrack::QStabilizerHybrid" internally switched between Qrack::QStabilizer and Qrack::QEngine to max...
Definition: qstabilizerhybrid.hpp:43
real1_f VarianceBitsFactorizedRdm(bool roundRz, const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const bitCapInt &offset=ZERO_BCI)
Get (reduced density matrix) expectation value of bits, given an array of qubit weights.
Definition: qstabilizerhybrid.hpp:876
bool TrimControls(const std::vector< bitLenInt > &lControls, std::vector< bitLenInt > &output, bool anti=false)
Definition: qstabilizerhybrid.cpp:284
real1_f ExpectationBitsAll(const std::vector< bitLenInt > &bits, const bitCapInt &offset=ZERO_BCI)
Get permutation expectation value of bits.
Definition: qstabilizerhybrid.hpp:800
void CPhaseFlipIfLess(const bitCapInt &greaterPerm, bitLenInt start, bitLenInt length, bitLenInt flagIndex)
The 6502 uses its carry flag also as a greater-than/less-than flag, for the CMP operation.
Definition: qstabilizerhybrid.hpp:535
void CMUL(const bitCapInt &toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Controlled multiplication by integer.
Definition: qstabilizerhybrid.hpp:651
void ZMask(const bitCapInt &mask)
Masked Z gate.
Definition: qstabilizerhybrid.cpp:1069
void SetPermutation(const bitCapInt &perm, const complex &phaseFac=CMPLX_DEFAULT_ARG)
Set to a specific permutation of all qubits.
Definition: qstabilizerhybrid.cpp:927
real1_f ProbAll(const bitCapInt &perm)
Direct measure of full permutation probability.
Definition: qstabilizerhybrid.hpp:475
bitLenInt ComposeNoClone(QInterfacePtr toCopy)
This is a variant of Compose() for a toCopy argument that will definitely not be reused once "Compose...
Definition: qstabilizerhybrid.hpp:456
bitLenInt maxEngineQubitCount
Definition: qstabilizerhybrid.hpp:55
void MACMtrx(const std::vector< bitLenInt > &controls, const complex *mtrx, bitLenInt target)
Apply an arbitrary single bit unitary transformation, with arbitrary (anti-)control bits.
Definition: qstabilizerhybrid.cpp:1288
QStabilizerHybrid(std::vector< QInterfaceEngine > eng, bitLenInt qBitCount, const bitCapInt &initState=ZERO_BCI, qrack_rand_gen_ptr rgp=nullptr, const complex &phaseFac=CMPLX_DEFAULT_ARG, bool doNorm=false, bool randomGlobalPhase=true, bool useHostMem=false, int64_t deviceId=-1, bool useHardwareRNG=true, bool ignored=false, real1_f norm_thresh=REAL1_EPSILON, std::vector< int64_t > devList={}, bitLenInt qubitThreshold=0U, real1_f separation_thresh=_qrack_qunit_sep_thresh)
Definition: qstabilizerhybrid.cpp:43
void CISqrtSwap(const std::vector< bitLenInt > &lControls, bitLenInt qubit1, bitLenInt qubit2)
Apply an inverse square root of swap with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:1006
bitLenInt ComposeEither(QStabilizerHybridPtr toCopy, bool willDestroy)
Definition: qstabilizerhybrid.cpp:497
void IISwap(bitLenInt qubit1, bitLenInt qubit2)
Inverse ISwap - Swap values of two bits in register, and apply phase factor of -i if bits are differe...
Definition: qstabilizerhybrid.hpp:486
void SetTInjection(bool useGadget)
Set the option to use T-injection gadgets (off by default)
Definition: qstabilizerhybrid.hpp:369
real1_f FractionalRzAngleWithFlush(bitLenInt i, real1_f angle, bool isGateSuppressed=false)
Definition: qstabilizerhybrid.hpp:191
int64_t devID
Definition: qstabilizerhybrid.hpp:60
void MACInvert(const std::vector< bitLenInt > &controls, const complex &topRight, const complex &bottomLeft, bitLenInt target)
Apply a single bit transformation that reverses bit probability and might effect phase,...
Definition: qstabilizerhybrid.cpp:1363
bitLenInt maxStateMapCacheQubitCount
Definition: qstabilizerhybrid.hpp:57
bool isFinished()
Returns "false" if asynchronous work is still running, and "true" if all previously dispatched asynch...
Definition: qstabilizerhybrid.hpp:383
QStabilizerHybridPtr RdmCloneHelper()
Definition: qstabilizerhybrid.hpp:263
void FlushBuffers()
Definition: qstabilizerhybrid.cpp:265
real1_f VarianceBitsAll(const std::vector< bitLenInt > &bits, const bitCapInt &offset=ZERO_BCI)
Direct measure of variance of listed permutation probability.
Definition: qstabilizerhybrid.hpp:847
std::vector< int64_t > deviceIDs
Definition: qstabilizerhybrid.hpp:65
bool doNormalize
Definition: qstabilizerhybrid.hpp:49
friend std::ostream & operator<<(std::ostream &os, const QStabilizerHybridPtr s)
Definition: qstabilizerhybrid.cpp:2146
void SqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
Square root of Swap gate.
Definition: qstabilizerhybrid.hpp:724
virtual bitLenInt Allocate(bitLenInt length)
Allocate new "length" count of |0> state qubits at end of qubit index position.
Definition: qinterface.hpp:470
void AntiCSqrtSwap(const std::vector< bitLenInt > &lControls, bitLenInt qubit1, bitLenInt qubit2)
Apply a square root of swap with arbitrary (anti) control bits.
Definition: qstabilizerhybrid.cpp:990
real1_f ProbRdm(bitLenInt qubit)
Direct measure of bit probability to be in |1> state, treating all ancillary qubits as post-selected ...
Definition: qstabilizerhybrid.hpp:402
void SetNcrp(real1_f ncrp)
Set the "Near-clifford rounding parameter" value, (between 0 and 1)
Definition: qstabilizerhybrid.hpp:368
void Swap(bitLenInt qubit1, bitLenInt qubit2)
Swap values of two bits in register.
Definition: qstabilizerhybrid.cpp:944
QInterfacePtr CloneBody(bool isCopy)
Definition: qstabilizerhybrid.cpp:364
void INCDECBCDC(const bitCapInt &toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCSC and DECSC (without overflow flag)
Definition: qstabilizerhybrid.hpp:617
void GetProbs(real1 *outputProbs)
Get the pure quantum state representation.
Definition: qstabilizerhybrid.cpp:684
void INCBCD(const bitCapInt &toAdd, bitLenInt start, bitLenInt length)
Add classical BCD integer (without sign)
Definition: qstabilizerhybrid.hpp:612
void RdmCloneFlush(real1_f threshold=FP_NORM_EPSILON)
Flush non-Clifford phase gate gadgets with angle below a threshold.
Definition: qstabilizerhybrid.cpp:1897
void CSqrtSwap(const std::vector< bitLenInt > &lControls, bitLenInt qubit1, bitLenInt qubit2)
Apply a square root of swap with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:974
virtual void UniformlyControlledSingleBit(const std::vector< bitLenInt > &controls, bitLenInt qubit, const complex *mtrxs)
Apply a "uniformly controlled" arbitrary single bit unitary transformation.
Definition: qinterface.hpp:627
void PhaseFlip()
Phase flip always - equivalent to Z X Z X on any bit in the QInterface.
Definition: qstabilizerhybrid.hpp:710
void FlushIfBlocked(bitLenInt control, bitLenInt target, bool isPhase=false)
Definition: qstabilizerhybrid.cpp:155
void Copy(QStabilizerHybridPtr orig)
Definition: qstabilizerhybrid.hpp:324
bitCapInt MAll()
Measure permutation state of all coherent bits.
Definition: qstabilizerhybrid.cpp:1577
virtual bitLenInt Compose(QInterfacePtr toCopy)
Combine another QInterface with this one, after the last bit index of this one.
Definition: qinterface.hpp:364
void INCS(const bitCapInt &toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Add a classical integer to the register, with sign and without carry.
Definition: qstabilizerhybrid.hpp:582
bool isClifford()
Returns "true" if current state is identifiably within the Clifford set, or "false" if it is not or c...
Definition: qstabilizerhybrid.hpp:441
double GetUnitaryFidelity()
When "Schmidt-decomposition rounding parameter" ("SDRP") is being used, starting from initial 1....
Definition: qstabilizerhybrid.hpp:371
void DEC(const bitCapInt &toSub, bitLenInt start, bitLenInt length)
Add integer (without sign)
Definition: qstabilizerhybrid.hpp:555
void FlushCliffordFromBuffers()
Definition: qstabilizerhybrid.hpp:231
void ClearAncilla(bitLenInt i)
Definition: qstabilizerhybrid.hpp:296
bool useTGadget
Definition: qstabilizerhybrid.hpp:50
bitCapInt IndexedADC(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength, bitLenInt carryIndex, const unsigned char *values)
Add to entangled 8 bit register state with a superposed index-offset-based read from classical memory...
Definition: qstabilizerhybrid.hpp:689
real1_f ProbMask(const bitCapInt &mask, const bitCapInt &permutation)
Direct measure of masked permutation probability.
Definition: qstabilizerhybrid.hpp:772
real1_f ExpectationFloatsFactorizedRdm(bool roundRz, const std::vector< bitLenInt > &bits, const std::vector< real1_f > &weights)
Get (reduced density matrix) expectation value of bits, given a (floating-point) array of qubit weigh...
Definition: qstabilizerhybrid.hpp:842
real1_f VarianceFloatsFactorized(const std::vector< bitLenInt > &bits, const std::vector< real1_f > &weights)
Direct measure of variance of listed bit string probability.
Definition: qstabilizerhybrid.hpp:881
complex GetAmplitude(const bitCapInt &perm)
Get the representational amplitude of a full permutation.
Definition: qstabilizerhybrid.hpp:474
void CINC(const bitCapInt &toAdd, bitLenInt inOutStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Add integer (without sign, with controls)
Definition: qstabilizerhybrid.hpp:573
void INCDECSC(const bitCapInt &toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCSC and DECSC (without overflow flag)
Definition: qstabilizerhybrid.hpp:606
void CacheEigenstate(bitLenInt target)
Definition: qstabilizerhybrid.cpp:324
void InvertBuffer(bitLenInt qubit)
Definition: qstabilizerhybrid.cpp:137
void INCDECC(const bitCapInt &toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCC and DECC (without sign, with carry)
Definition: qstabilizerhybrid.hpp:591
bitLenInt ancillaCount
Definition: qstabilizerhybrid.hpp:53
void DumpBuffers()
Definition: qstabilizerhybrid.hpp:82
real1_f SumSqrDiff(QInterfacePtr toCompare)
Calculates (1 - <\psi_e|\psi_c>) between states |\psi_c> and |\psi_e>.
Definition: qstabilizerhybrid.hpp:778
std::vector< MpsShardPtr > shards
Definition: qstabilizerhybrid.hpp:68
bitCapInt IndexedSBC(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength, bitLenInt carryIndex, const unsigned char *values)
Subtract from an entangled 8 bit register state with a superposed index-offset-based read from classi...
Definition: qstabilizerhybrid.hpp:696
void YMask(const bitCapInt &mask)
Masked Y gate.
Definition: qstabilizerhybrid.cpp:1054
real1_f ACProbRdm(bitLenInt control, bitLenInt target)
Definition: qstabilizerhybrid.hpp:425
bool IsBuffered()
Definition: qstabilizerhybrid.hpp:100
void CMULModNOut(const bitCapInt &toMul, const bitCapInt &modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Controlled multiplication modulo N by integer, (out of place)
Definition: qstabilizerhybrid.hpp:663
complex phaseFactor
Definition: qstabilizerhybrid.hpp:61
void MCMtrx(const std::vector< bitLenInt > &controls, const complex *mtrx, bitLenInt target)
Apply an arbitrary single bit unitary transformation, with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:1161
void INCDECSC(const bitCapInt &toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex, bitLenInt carryIndex)
Common driver method behind INCSC and DECSC (with overflow flag)
Definition: qstabilizerhybrid.hpp:600
void UpdateRunningNorm(real1_f norm_thresh=REAL1_DEFAULT_ARG)
Force a calculation of the norm of the state vector, in order to make it unit length before the next ...
Definition: qstabilizerhybrid.hpp:788
bitLenInt maxAncillaCount
Definition: qstabilizerhybrid.hpp:56
void DECS(const bitCapInt &toSub, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Add a classical integer to the register, with sign and without carry.
Definition: qstabilizerhybrid.hpp:564
bitLenInt Compose(QInterfacePtr toCopy, bitLenInt start)
Compose() a QInterface peer, inserting its qubit into index order at start index.
Definition: qstabilizerhybrid.hpp:451
real1_f roundingThreshold
Definition: qstabilizerhybrid.hpp:59
bitLenInt ComposeNoClone(QStabilizerHybridPtr toCopy)
Definition: qstabilizerhybrid.hpp:455
void PhaseFlipIfLess(const bitCapInt &greaterPerm, bitLenInt start, bitLenInt length)
This is an expedient for an adaptive Grover's search for a function's global minimum.
Definition: qstabilizerhybrid.hpp:540
bool TrySeparate(bitLenInt qubit)
Single-qubit TrySeparate()
Definition: qstabilizerhybrid.cpp:2105
QInterfacePtr engine
Definition: qstabilizerhybrid.hpp:63
real1_f VarianceBitsAllRdm(bool roundRz, const std::vector< bitLenInt > &bits, const bitCapInt &offset=ZERO_BCI)
Direct measure of (reduced density matrix) variance of listed permutation probability.
Definition: qstabilizerhybrid.hpp:855
double logFidelity
Definition: qstabilizerhybrid.hpp:62
void ISwap(bitLenInt qubit1, bitLenInt qubit2)
Swap values of two bits in register, and apply phase factor of i if bits are different.
Definition: qstabilizerhybrid.hpp:485
real1_f ProbMaskRdm(bool roundRz, const bitCapInt &mask, const bitCapInt &permutation)
Direct measure of masked permutation probability, treating all ancillary qubits as post-selected T ga...
Definition: qstabilizerhybrid.cpp:405
void X(bitLenInt q)
Definition: qstabilizerhybrid.hpp:534
bitCapInt IndexedLDA(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength, const unsigned char *values, bool resetValue=true)
Set 8 bit register bits by a superposed index-offset-based read from classical memory.
Definition: qstabilizerhybrid.hpp:682
std::vector< QInterfaceEngine > cloneEngineTypes
Definition: qstabilizerhybrid.hpp:67
void ISwapHelper(bitLenInt qubit1, bitLenInt qubit2, bool inverse)
Definition: qstabilizerhybrid.cpp:2049
std::map< bitCapInt, int > MultiShotMeasureMask(const std::vector< bitCapInt > &qPowers, unsigned shots)
Statistical measure of masked permutation probability.
Definition: qstabilizerhybrid.cpp:1684
void CUniformParityRZ(const std::vector< bitLenInt > &controls, const bitCapInt &mask, real1_f angle)
If the controls are set and the target qubit set parity is odd, this applies a phase factor of .
Definition: qstabilizerhybrid.hpp:524
bitLenInt Compose(QStabilizerHybridPtr toCopy)
Definition: qstabilizerhybrid.hpp:448
real1_f ExpectationBitsFactorized(const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const bitCapInt &offset=ZERO_BCI)
Get expectation value of bits, given an array of qubit weights.
Definition: qstabilizerhybrid.hpp:820
void Mtrx(const complex *mtrx, bitLenInt target)
Apply an arbitrary single bit unitary transformation.
Definition: qstabilizerhybrid.cpp:1085
real1_f VarianceBitsFactorized(const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const bitCapInt &offset=ZERO_BCI)
Get expectation value of bits, given an array of qubit weights.
Definition: qstabilizerhybrid.hpp:867
bool EitherIsBuffered(bool logical)
Definition: qstabilizerhybrid.hpp:88
real1_f ExpectationFloatsFactorized(const std::vector< bitLenInt > &bits, const std::vector< real1_f > &weights)
Get expectation value of bits, given a (floating-point) array of qubit weights.
Definition: qstabilizerhybrid.hpp:834
bool M(bitLenInt q)
Definition: qstabilizerhybrid.hpp:532
void IMULModNOut(const bitCapInt &toMul, const bitCapInt &modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
Inverse of multiplication modulo N by integer, (out of place)
Definition: qstabilizerhybrid.hpp:639
bitCapIntOcl GetMaxSize()
Definition: qstabilizerhybrid.hpp:912
bool useHostRam
Definition: qstabilizerhybrid.hpp:48
void MULModNOut(const bitCapInt &toMul, const bitCapInt &modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
Multiplication modulo N by integer, (out of place)
Definition: qstabilizerhybrid.hpp:633
bool GetTInjection()
Get the option to use T-injection gadgets.
Definition: qstabilizerhybrid.hpp:370
void AntiCISqrtSwap(const std::vector< bitLenInt > &lControls, bitLenInt qubit1, bitLenInt qubit2)
Apply an inverse square root of swap with arbitrary (anti) control bits.
Definition: qstabilizerhybrid.cpp:1022
real1_f ProbAllRdm(bool roundRz, const bitCapInt &fullRegister)
Direct measure of full permutation probability, treating all ancillary qubits as post-selected T gate...
Definition: qstabilizerhybrid.cpp:392
void ZeroPhaseFlip(bitLenInt start, bitLenInt length)
Reverse the phase of the state where the register equals zero.
Definition: qstabilizerhybrid.hpp:718
QInterfacePtr Clone()
Clone this QInterface.
Definition: qstabilizerhybrid.hpp:899
bool IsLogicalBuffered()
Definition: qstabilizerhybrid.hpp:101
real1_f ExpectationBitsFactorizedRdm(bool roundRz, const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const bitCapInt &offset=ZERO_BCI)
Get (reduced density matrix) expectation value of bits, given an array of qubit weights.
Definition: qstabilizerhybrid.hpp:829
void CSwap(const std::vector< bitLenInt > &lControls, bitLenInt qubit1, bitLenInt qubit2)
Apply a swap with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:958
void SetDevice(int64_t dID)
Set the device index, if more than one device is available.
Definition: qstabilizerhybrid.hpp:902
void MCPhase(const std::vector< bitLenInt > &controls, const complex &topLeft, const complex &bottomRight, bitLenInt target)
Apply a single bit transformation that only effects phase, with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:1187
real1_f ApproxCompareHelper(QStabilizerHybridPtr toCompare, bool isDiscreteBool, real1_f error_tol=TRYDECOMPOSE_EPSILON)
Definition: qstabilizerhybrid.cpp:1966
void ISqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
Inverse square root of Swap gate.
Definition: qstabilizerhybrid.hpp:734
bool isBinaryDecisionTree()
Returns "true" if current state representation is definitely a binary decision tree,...
Definition: qstabilizerhybrid.hpp:445
bool EitherIsProbBuffered(bool logical)
Definition: qstabilizerhybrid.hpp:102
void SetAmplitude(const bitCapInt &perm, const complex &amp)
Sets the representational amplitude of a full permutation.
Definition: qstabilizerhybrid.hpp:477
void Dispose(bitLenInt start, bitLenInt length)
Minimally decompose a set of contiguous bits from the separably composed unit, and discard the separa...
Definition: qstabilizerhybrid.cpp:627
bool isRoundingFlushed
Definition: qstabilizerhybrid.hpp:51
real1_f separabilityThreshold
Definition: qstabilizerhybrid.hpp:58
complex GetAmplitudeOrProb(const bitCapInt &perm, bool isProb=false)
Definition: qstabilizerhybrid.cpp:701
QStabilizerHybrid(bitLenInt qBitCount, const bitCapInt &initState=ZERO_BCI, qrack_rand_gen_ptr rgp=nullptr, const complex &phaseFac=CMPLX_DEFAULT_ARG, bool doNorm=false, bool randomGlobalPhase=true, bool useHostMem=false, int64_t deviceId=-1, bool useHardwareRNG=true, bool ignored=false, real1_f norm_thresh=REAL1_EPSILON, std::vector< int64_t > devList={}, bitLenInt qubitThreshold=0U, real1_f separation_thresh=_qrack_qunit_sep_thresh)
Definition: qstabilizerhybrid.hpp:358
void Decompose(bitLenInt start, QInterfacePtr dest)
Minimally decompose a set of contiguous bits from the separably composed unit, into "destination".
Definition: qstabilizerhybrid.hpp:461
QInterfacePtr Copy()
Copy this QInterface.
Definition: qstabilizerhybrid.hpp:900
QUnitCliffordPtr stabilizer
Definition: qstabilizerhybrid.hpp:64
void SetQuantumState(const complex *inputState)
Set an arbitrary pure quantum state representation.
Definition: qstabilizerhybrid.cpp:888
bitLenInt deadAncillaCount
Definition: qstabilizerhybrid.hpp:54
void CIMULModNOut(const bitCapInt &toMul, const bitCapInt &modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Inverse of controlled multiplication modulo N by integer, (out of place)
Definition: qstabilizerhybrid.hpp:669
void INC(const bitCapInt &toAdd, bitLenInt start, bitLenInt length)
Add integer (without sign)
Definition: qstabilizerhybrid.hpp:546
void Copy(QInterfacePtr orig)
Definition: qstabilizerhybrid.hpp:323
bitLenInt thresholdQubits
Definition: qstabilizerhybrid.hpp:52
real1_f ExpectationBitsAllRdm(bool roundRz, const std::vector< bitLenInt > &bits, const bitCapInt &offset=ZERO_BCI)
Get permutation expectation value of bits, treating all ancillary qubits as post-selected T gate gadg...
Definition: qstabilizerhybrid.hpp:808
void SwitchToEngine()
Switches between CPU and GPU used modes.
Definition: qstabilizerhybrid.cpp:422
void DIV(const bitCapInt &toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
Divide by integer.
Definition: qstabilizerhybrid.hpp:628
std::vector< QInterfaceEngine > engineTypes
Definition: qstabilizerhybrid.hpp:66
void NormalizeState(real1_f nrm=REAL1_DEFAULT_ARG, real1_f norm_thresh=REAL1_DEFAULT_ARG, real1_f phaseArg=ZERO_R1_F)
Apply the normalization factor found by UpdateRunningNorm() or on the fly by a single bit gate.
Definition: qstabilizerhybrid.cpp:2092
bool isClifford(bitLenInt qubit)
Returns "true" if current qubit state is identifiably within the Clifford set, or "false" if it is no...
Definition: qstabilizerhybrid.hpp:443
void Finish()
If asynchronous work is still running, block until it finishes.
Definition: qstabilizerhybrid.hpp:374
bool CollapseSeparableShard(bitLenInt qubit)
Definition: qstabilizerhybrid.cpp:241
bool IsLogicalProbBuffered()
Definition: qstabilizerhybrid.hpp:125
bool ForceMParity(const bitCapInt &mask, bool result, bool doForce=true)
Act as if is a measurement of parity of the masked set of qubits was applied, except force the (usual...
Definition: qstabilizerhybrid.cpp:1880
friend std::istream & operator>>(std::istream &is, const QStabilizerHybridPtr s)
Definition: qstabilizerhybrid.cpp:2169
void ResetUnitaryFidelity()
Reset the internal fidelity calculation tracker to 1.0.
Definition: qstabilizerhybrid.hpp:372
void MACPhase(const std::vector< bitLenInt > &controls, const complex &topLeft, const complex &bottomRight, bitLenInt target)
Apply a single bit transformation that only effects phase, with arbitrary (anti-)control bits.
Definition: qstabilizerhybrid.cpp:1314
void CheckShots(unsigned shots, const bitCapInt &m, real1_f partProb, const std::vector< bitCapInt > &qPowers, std::vector< real1_f > &rng, F fn)
Definition: qstabilizerhybrid.hpp:155
real1_f Prob(bitLenInt qubit)
Direct measure of bit probability to be in |1> state.
Definition: qstabilizerhybrid.cpp:1411
std::unique_ptr< complex[]> GetQubitReducedDensityMatrix(bitLenInt qubit)
Definition: qstabilizerhybrid.hpp:127
void CPOWModNOut(const bitCapInt &base, const bitCapInt &modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Controlled, raise a classical base to a quantum power, modulo N, (out of place)
Definition: qstabilizerhybrid.hpp:675
void CDIV(const bitCapInt &toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Controlled division by power of integer.
Definition: qstabilizerhybrid.hpp:657
bool ForceM(bitLenInt qubit, bool result, bool doForce=true, bool doApply=true)
Act as if is a measurement was applied, except force the (usually random) result.
Definition: qstabilizerhybrid.cpp:1496
real1_f ProbParity(const bitCapInt &mask)
Overall probability of any odd permutation of the masked set of bits.
Definition: qstabilizerhybrid.cpp:1867
void GetQuantumState(complex *outputState)
Get the pure quantum state representation.
Definition: qstabilizerhybrid.cpp:667
QInterfacePtr MakeEngine(const bitCapInt &perm=ZERO_BCI)
Definition: qstabilizerhybrid.cpp:120
void XMask(const bitCapInt &mask)
Masked X gate.
Definition: qstabilizerhybrid.cpp:1039
void MUL(const bitCapInt &toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
Multiply by integer.
Definition: qstabilizerhybrid.hpp:623
void Hash(bitLenInt start, bitLenInt length, const unsigned char *values)
Transform a length of qubit register via lookup through a hash table.
Definition: qstabilizerhybrid.hpp:703
void FlushH(bitLenInt qubit)
Definition: qstabilizerhybrid.cpp:146
bool ApproxCompare(QInterfacePtr toCompare, real1_f error_tol=TRYDECOMPOSE_EPSILON)
Compare state vectors approximately, to determine whether this state vector is the same as the target...
Definition: qstabilizerhybrid.hpp:782
void POWModNOut(const bitCapInt &base, const bitCapInt &modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
Raise a classical base to a quantum power, modulo N, (out of place)
Definition: qstabilizerhybrid.hpp:645
QUnitStateVectorPtr stateMapCache
Definition: qstabilizerhybrid.hpp:69
void FSim(real1_f theta, real1_f phi, bitLenInt qubit1, bitLenInt qubit2)
The 2-qubit "fSim" gate, (useful in the simulation of particles with fermionic statistics)
Definition: qstabilizerhybrid.hpp:744
void SetConcurrency(uint32_t threadCount)
Set the number of threads in parallel for loops, per component QEngine.
Definition: qstabilizerhybrid.hpp:394
real1_f VarianceFloatsFactorizedRdm(bool roundRz, const std::vector< bitLenInt > &bits, const std::vector< real1_f > &weights)
Direct measure of (reduced density matrix) variance of bits, given an array of qubit weights.
Definition: qstabilizerhybrid.hpp:889
QUnitCliffordPtr MakeStabilizer(const bitCapInt &perm=ZERO_BCI)
Definition: qstabilizerhybrid.cpp:115
real1_f ExpVarFactorized(bool isExp, bool isFloat, const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const std::vector< real1_f > &weights, const bitCapInt &offset, bool roundRz)
Definition: qstabilizerhybrid.hpp:272
bitLenInt Compose(QInterfacePtr toCopy)
Combine another QInterface with this one, after the last bit index of this one.
Definition: qstabilizerhybrid.hpp:449
real1_f CProbRdm(bitLenInt control, bitLenInt target)
Definition: qstabilizerhybrid.hpp:416
std::vector< real1_f > GenerateShotProbs(unsigned shots)
Definition: qstabilizerhybrid.hpp:178
void Dump()
If asynchronous work is still running, let the simulator know that it can be aborted.
Definition: qstabilizerhybrid.hpp:385
int64_t GetDeviceID()
Definition: qstabilizerhybrid.hpp:910
void MCInvert(const std::vector< bitLenInt > &controls, const complex &topRight, const complex &bottomLeft, bitLenInt target)
Apply a single bit transformation that reverses bit probability and might effect phase,...
Definition: qstabilizerhybrid.cpp:1240
bool IsProbBuffered()
Definition: qstabilizerhybrid.hpp:124
Half-precision floating-point type.
Definition: half.hpp:2222
virtual void DECS(const bitCapInt &toSub, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Subtract a classical integer from the register, with sign and without carry.
Definition: qinterface.hpp:2193
virtual void INCDECC(const bitCapInt &toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCC and DECC.
Definition: arithmetic.cpp:52
virtual void CINC(const bitCapInt &toAdd, bitLenInt inOutStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Add integer (without sign, with controls)
Definition: arithmetic.cpp:78
virtual void INCS(const bitCapInt &toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Add a classical integer to the register, with sign and without carry.
Definition: qinterface.hpp:2182
virtual void DEC(const bitCapInt &toSub, bitLenInt start, bitLenInt length)
Subtract classical integer (without sign)
Definition: qinterface.hpp:2134
virtual void INC(const bitCapInt &toAdd, bitLenInt start, bitLenInt length)
Add integer (without sign)
Definition: arithmetic.cpp:20
virtual void CNOT(bitLenInt control, bitLenInt target)
Controlled NOT gate.
Definition: qinterface.hpp:709
virtual void UniformlyControlledSingleBit(const std::vector< bitLenInt > &controls, bitLenInt qubit, const complex *mtrxs)
Apply a "uniformly controlled" arbitrary single bit unitary transformation.
Definition: qinterface.hpp:627
virtual void X(bitLenInt qubit)
X gate.
Definition: qinterface.hpp:1084
virtual void AntiCNOT(bitLenInt control, bitLenInt target)
Anti controlled NOT gate.
Definition: qinterface.hpp:720
virtual void U(bitLenInt target, real1_f theta, real1_f phi, real1_f lambda)
General unitary gate.
Definition: rotational.cpp:18
virtual bool M(bitLenInt qubit)
Measurement gate.
Definition: qinterface.hpp:1013
virtual void ISqrtSwap(bitLenInt qubit1, bitLenInt qubit2)
Inverse square root of Swap gate.
Definition: gates.cpp:225
virtual void SqrtSwap(bitLenInt qubit1, bitLenInt qubit2)
Square root of Swap gate.
Definition: gates.cpp:202
virtual real1_f VarianceBitsFactorized(const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const bitCapInt &offset=ZERO_BCI)
Get expectation value of bits, given an array of qubit weights.
Definition: qinterface.cpp:578
virtual real1_f VarianceFloatsFactorized(const std::vector< bitLenInt > &bits, const std::vector< real1_f > &weights)
Direct measure of variance of listed bit string probability.
Definition: qinterface.cpp:618
virtual QInterfacePtr Copy()
Copy this QInterface.
Definition: qinterface.hpp:2965
virtual real1_f ExpectationBitsAll(const std::vector< bitLenInt > &bits, const bitCapInt &offset=ZERO_BCI)
Get permutation expectation value of bits.
Definition: qinterface.hpp:2607
bitCapIntOcl GetMaxSize()
Get maximum number of amplitudes that can be allocated on current device.
Definition: qinterface.hpp:2980
virtual real1_f VarianceBitsAll(const std::vector< bitLenInt > &bits, const bitCapInt &offset=ZERO_BCI)
Direct measure of variance of listed permutation probability.
Definition: qinterface.hpp:2491
virtual real1_f ExpectationBitsFactorized(const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const bitCapInt &offset=ZERO_BCI)
Get expectation value of bits, given an array of qubit weights.
Definition: qinterface.cpp:541
virtual real1_f ExpectationFloatsFactorized(const std::vector< bitLenInt > &bits, const std::vector< real1_f > &weights)
Get expectation value of bits, given a (floating-point) array of qubit weights.
Definition: qinterface.cpp:768
GLOSSARY: bitLenInt - "bit-length integer" - unsigned integer ID of qubit position in register bitCap...
Definition: complex16x2simd.hpp:25
@ QINTERFACE_HYBRID
Create a QHybrid, switching between QEngineCPU and QEngineOCL as efficient.
Definition: qinterface.hpp:57
std::shared_ptr< QInterface > QInterfacePtr
Definition: qinterface.hpp:29
const real1_f _qrack_qunit_sep_thresh
Definition: qrack_functions.hpp:235
QRACK_CONST real1_f TRYDECOMPOSE_EPSILON
Definition: qrack_types.hpp:260
std::shared_ptr< QStabilizerHybrid > QStabilizerHybridPtr
Definition: qstabilizerhybrid.hpp:35
void mul2x2(const complex *left, const complex *right, complex *out)
Definition: functions.cpp:111
QRACK_CONST real1 HALF_R1
Definition: qrack_types.hpp:184
half_float::half real1
Definition: qrack_types.hpp:94
std::complex< real1 > complex
Definition: qrack_types.hpp:128
std::shared_ptr< QUnitStateVector > QUnitStateVectorPtr
Definition: qunitstatevector.hpp:17
QRACK_CONST real1 FP_NORM_EPSILON
Definition: qrack_types.hpp:258
bitCapInt pow2(const bitLenInt &p)
Definition: qrack_functions.hpp:136
std::shared_ptr< QUnitClifford > QUnitCliffordPtr
Definition: qunitclifford.hpp:20
double norm(const complex2 &c)
Definition: complex16x2simd.hpp:101
QRACK_CONST real1 REAL1_EPSILON
Definition: qrack_types.hpp:200
QRACK_CONST complex ONE_CMPLX
Definition: qrack_types.hpp:252
QRACK_CONST real1 ONE_R1
Definition: qrack_types.hpp:185
QRACK_CONST real1 ZERO_R1
Definition: qrack_types.hpp:183
float real1_f
Definition: qrack_types.hpp:95
float real1_s
Definition: qrack_types.hpp:96
QRACK_CONST complex CMPLX_DEFAULT_ARG
Definition: qrack_types.hpp:257
std::shared_ptr< MpsShard > MpsShardPtr
Definition: mpsshard.hpp:18
QRACK_CONST complex I_CMPLX
Definition: qrack_types.hpp:254
QRACK_CONST complex ZERO_CMPLX
Definition: qrack_types.hpp:253
QRACK_CONST real1 PI_R1
Definition: qrack_types.hpp:178
void reverse(BidirectionalIterator first, BidirectionalIterator last, const bitCapInt &stride)
const bitCapInt ZERO_BCI
Definition: qrack_types.hpp:130
HALF_CONSTEXPR half abs(half arg)
Absolute value.
Definition: half.hpp:2975
half sin(half arg)
Sine function.
Definition: half.hpp:3885
half cos(half arg)
Cosine function.
Definition: half.hpp:3922
long lround(half arg)
Nearest integer.
Definition: half.hpp:4505
half exp(half arg)
Exponential function.
Definition: half.hpp:3201
#define REAL1_DEFAULT_ARG
Definition: qrack_types.hpp:177
#define QRACK_CONST
Definition: qrack_types.hpp:174
#define bitLenInt
Definition: qrack_types.hpp:38
#define ZERO_R1_F
Definition: qrack_types.hpp:160
#define qrack_rand_gen_ptr
Definition: qrack_types.hpp:156
#define bitCapInt
Definition: qrack_types.hpp:62
#define bitCapIntOcl
Definition: qrack_types.hpp:50
#define ONE_R1_F
Definition: qrack_types.hpp:163
#define QINTERFACE_TO_QALU(qReg)
Definition: qstabilizerhybrid.hpp:18
#define QINTERFACE_TO_QPARITY(qReg)
Definition: qstabilizerhybrid.hpp:19
Definition: qstabilizerhybrid.hpp:23
QUnitCliffordAmp(const complex &a, QUnitCliffordPtr s)
Definition: qstabilizerhybrid.hpp:27
QUnitCliffordPtr stabilizer
Definition: qstabilizerhybrid.hpp:25
complex amp
Definition: qstabilizerhybrid.hpp:24