Qrack  9.0
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 isSparse;
51  bool useTGadget;
60  int64_t devID;
64  std::vector<int64_t> deviceIDs;
65  std::vector<QInterfaceEngine> engineTypes;
66  std::vector<QInterfaceEngine> cloneEngineTypes;
67  std::vector<MpsShardPtr> shards;
68  std::map<bitCapInt, complex> stateMapCache;
69 
73 
74  void InvertBuffer(bitLenInt qubit);
75  void FlushH(bitLenInt qubit);
76  void FlushIfBlocked(bitLenInt control, bitLenInt target, bool isPhase = false);
78  bool TrimControls(const std::vector<bitLenInt>& lControls, std::vector<bitLenInt>& output, bool anti = false);
79  void CacheEigenstate(bitLenInt target);
80  void FlushBuffers();
81  void DumpBuffers()
82  {
83  for (size_t i = 0; i < shards.size(); ++i) {
84  shards[i] = NULL;
85  }
86  }
87  bool EitherIsBuffered(bool logical)
88  {
89  const size_t maxLcv = logical ? (size_t)qubitCount : shards.size();
90  for (size_t i = 0U; i < maxLcv; ++i) {
91  if (shards[i]) {
92  // We have a cached non-Clifford operation.
93  return true;
94  }
95  }
96 
97  return false;
98  }
99  bool IsBuffered() { return EitherIsBuffered(false); }
100  bool IsLogicalBuffered() { return EitherIsBuffered(true); }
101  bool EitherIsProbBuffered(bool logical)
102  {
103  const size_t maxLcv = logical ? (size_t)qubitCount : shards.size();
104  for (size_t i = 0U; i < maxLcv; ++i) {
105  MpsShardPtr shard = shards[i];
106  if (!shard) {
107  continue;
108  }
109  if (shard->IsHPhase() || shard->IsHInvert()) {
110  FlushH(i);
111  }
112  if (shard->IsInvert()) {
113  InvertBuffer(i);
114  }
115  if (!shard->IsPhase()) {
116  // We have a cached non-Clifford operation.
117  return true;
118  }
119  }
120 
121  return false;
122  }
123  bool IsProbBuffered() { return EitherIsProbBuffered(false); }
125 
126  std::unique_ptr<complex[]> GetQubitReducedDensityMatrix(bitLenInt qubit)
127  {
128  // Form the reduced density matrix of the single qubit.
129  const real1 z = (real1)(ONE_R1_F - 2 * stabilizer->Prob(qubit));
130  stabilizer->H(qubit);
131  const real1 x = (real1)(ONE_R1_F - 2 * stabilizer->Prob(qubit));
132  stabilizer->S(qubit);
133  const real1 y = (real1)(ONE_R1_F - 2 * stabilizer->Prob(qubit));
134  stabilizer->IS(qubit);
135  stabilizer->H(qubit);
136 
137  std::unique_ptr<complex[]> dMtrx(new complex[4]);
138  dMtrx[0] = (ONE_CMPLX + z) / complex((real1)2, ZERO_R1);
139  dMtrx[1] = x / complex((real1)2, ZERO_R1) - I_CMPLX * (y / complex((real1)2, ZERO_R1));
140  dMtrx[2] = x / complex((real1)2, ZERO_R1) + I_CMPLX * (y / complex((real1)2, ZERO_R1));
141  dMtrx[3] = (ONE_CMPLX + z) / complex((real1)2, ZERO_R1);
142  if (shards[qubit]) {
143  complex adj[4]{ std::conj(shards[qubit]->gate[0]), std::conj(shards[qubit]->gate[2]),
144  std::conj(shards[qubit]->gate[1]), std::conj(shards[qubit]->gate[3]) };
145  complex out[4];
146  mul2x2(dMtrx.get(), adj, out);
147  mul2x2(shards[qubit]->gate, out, dMtrx.get());
148  }
149 
150  return dMtrx;
151  }
152 
153  template <typename F>
154  void CheckShots(unsigned shots, bitCapInt m, real1_f partProb, const std::vector<bitCapInt>& qPowers,
155  std::vector<real1_f>& rng, F fn)
156  {
157  for (int64_t shot = rng.size() - 1U; shot >= 0; --shot) {
158  if (rng[shot] >= partProb) {
159  break;
160  }
161 
162  bitCapInt sample = 0U;
163  for (size_t i = 0U; i < qPowers.size(); ++i) {
164  if (m & qPowers[i]) {
165  sample |= pow2(i);
166  }
167  }
168  fn(sample, shot);
169 
170  rng.erase(rng.begin() + shot);
171  if (!rng.size()) {
172  break;
173  }
174  }
175  }
176 
177  std::vector<real1_f> GenerateShotProbs(unsigned shots)
178  {
179  std::vector<real1_f> rng;
180  rng.reserve(shots);
181  for (unsigned shot = 0U; shot < shots; ++shot) {
182  rng.push_back(Rand());
183  }
184  std::sort(rng.begin(), rng.end());
185  std::reverse(rng.begin(), rng.end());
186 
187  return rng;
188  }
189 
190  real1_f FractionalRzAngleWithFlush(bitLenInt i, real1_f angle, bool isGateSuppressed = false)
191  {
192  const real1_f sectorAngle = PI_R1 / 2;
193  const real1_f Period = 2 * PI_R1;
194  while (angle >= Period) {
195  angle -= Period;
196  }
197  while (angle < 0) {
198  angle += Period;
199  }
200 
201  int sector = std::round((real1_s)(angle / sectorAngle));
202  if (!isGateSuppressed) {
203  switch (sector) {
204  case 1U:
205  stabilizer->S(i);
206  break;
207  case 2U:
208  stabilizer->Z(i);
209  break;
210  case 3U:
211  stabilizer->IS(i);
212  break;
213  case 0U:
214  default:
215  break;
216  }
217  }
218 
219  real1_f correctionAngle = angle - (sector * sectorAngle);
220  if (correctionAngle > PI_R1) {
221  correctionAngle -= Period;
222  }
223  if (correctionAngle <= -PI_R1) {
224  correctionAngle += Period;
225  }
226 
227  return correctionAngle;
228  }
229 
231  {
232  for (size_t i = 0U; i < qubitCount; ++i) {
233  // Flush all buffers as close as possible to Clifforrd.
234  const MpsShardPtr& shard = shards[i];
235  if (!shard) {
236  continue;
237  }
238  if (shard->IsHPhase() || shard->IsHInvert()) {
239  FlushH(i);
240  }
241  if (shard->IsInvert()) {
242  InvertBuffer(i);
243  }
244  if (!shard->IsPhase()) {
245  // We have a cached non-phase operation.
246  continue;
247  }
248  const real1 angle = (real1)(FractionalRzAngleWithFlush(i, std::arg(shard->gate[3U] / shard->gate[0U])) / 2);
249  if ((2 * abs(angle) / PI_R1) <= FP_NORM_EPSILON) {
250  shards[i] = NULL;
251  continue;
252  }
253  const real1 angleCos = cos(angle);
254  const real1 angleSin = sin(angle);
255  shard->gate[0U] = complex(angleCos, -angleSin);
256  shard->gate[3U] = complex(angleCos, angleSin);
257  }
258 
259  RdmCloneFlush();
260  }
261 
262  void CombineAncillae();
263 
265  {
266  CombineAncillae();
267  QStabilizerHybridPtr clone = std::dynamic_pointer_cast<QStabilizerHybrid>(Clone());
268  clone->RdmCloneFlush(ONE_R1 / 2);
269 
270  return clone;
271  }
272  void RdmCloneFlush(real1_f threshold = FP_NORM_EPSILON);
273 
274  real1_f ExpectationFactorized(bool isFloat, const std::vector<bitLenInt>& bits, const std::vector<bitCapInt>& perms,
275  const std::vector<real1_f>& weights, bitCapInt offset, bool roundRz)
276  {
277  if (engine) {
278  return isFloat ? engine->ExpectationFloatsFactorizedRdm(roundRz, bits, weights)
279  : engine->ExpectationBitsFactorizedRdm(roundRz, bits, perms, offset);
280  ;
281  }
282 
283  CombineAncillae();
284 
285  if (!roundRz) {
286  return isFloat ? stabilizer->ExpectationFloatsFactorizedRdm(roundRz, bits, weights)
287  : stabilizer->ExpectationBitsFactorizedRdm(roundRz, bits, perms, offset);
288  }
289 
290  return isFloat ? RdmCloneHelper()->stabilizer->ExpectationFloatsFactorizedRdm(roundRz, bits, weights)
291  : RdmCloneHelper()->stabilizer->ExpectationBitsFactorizedRdm(roundRz, bits, perms, offset);
292  }
293 
295  {
296  if (stabilizer->TrySeparate(i)) {
297  stabilizer->Dispose(i, 1U);
298  shards.erase(shards.begin() + i);
299  } else {
300  const bitLenInt deadIndex = qubitCount + ancillaCount - 1U;
301  stabilizer->SetBit(i, false);
302  if (i != deadIndex) {
303  stabilizer->Swap(i, deadIndex);
304  shards[i].swap(shards[deadIndex]);
305  }
306  shards.erase(shards.begin() + deadIndex);
308  }
309  --ancillaCount;
310  }
311 
313  QStabilizerHybridPtr toCompare, bool isDiscreteBool, real1_f error_tol = TRYDECOMPOSE_EPSILON);
314 
315  void ISwapHelper(bitLenInt qubit1, bitLenInt qubit2, bool inverse);
316 
317  complex GetAmplitudeOrProb(bitCapInt perm, bool isProb = false);
318 
319 public:
320  QStabilizerHybrid(std::vector<QInterfaceEngine> eng, bitLenInt qBitCount, bitCapInt initState = 0U,
321  qrack_rand_gen_ptr rgp = nullptr, complex phaseFac = CMPLX_DEFAULT_ARG, bool doNorm = false,
322  bool randomGlobalPhase = true, bool useHostMem = false, int64_t deviceId = -1, bool useHardwareRNG = true,
323  bool useSparseStateVec = false, real1_f norm_thresh = REAL1_EPSILON, std::vector<int64_t> devList = {},
324  bitLenInt qubitThreshold = 0U, real1_f separation_thresh = FP_NORM_EPSILON_F);
325 
326  QStabilizerHybrid(bitLenInt qBitCount, bitCapInt initState = 0U, qrack_rand_gen_ptr rgp = nullptr,
327  complex phaseFac = CMPLX_DEFAULT_ARG, bool doNorm = false, bool randomGlobalPhase = true,
328  bool useHostMem = false, int64_t deviceId = -1, bool useHardwareRNG = true, bool useSparseStateVec = false,
329  real1_f norm_thresh = REAL1_EPSILON, std::vector<int64_t> devList = {}, bitLenInt qubitThreshold = 0U,
330  real1_f separation_thresh = FP_NORM_EPSILON_F)
331  : QStabilizerHybrid({ QINTERFACE_OPTIMAL_BASE }, qBitCount, initState, rgp, phaseFac, doNorm, randomGlobalPhase,
332  useHostMem, deviceId, useHardwareRNG, useSparseStateVec, norm_thresh, devList, qubitThreshold,
333  separation_thresh)
334  {
335  }
336 
337  void SetTInjection(bool useGadget) { useTGadget = useGadget; }
338  bool GetTInjection() { return useTGadget; }
339 
340  void Finish()
341  {
342  if (stabilizer) {
343  stabilizer->Finish();
344  } else {
345  engine->Finish();
346  }
347  };
348 
349  bool isFinished() { return (!stabilizer || stabilizer->isFinished()) && (!engine || engine->isFinished()); }
350 
351  void Dump()
352  {
353  if (stabilizer) {
354  stabilizer->Dump();
355  } else {
356  engine->Dump();
357  }
358  }
359 
360  void SetConcurrency(uint32_t threadCount)
361  {
362  QInterface::SetConcurrency(threadCount);
363  if (engine) {
365  }
366  }
367 
369  {
370  if (!ancillaCount || stabilizer->IsSeparable(qubit)) {
371  return Prob(qubit);
372  }
373 
374  std::unique_ptr<complex[]> dMtrx = GetQubitReducedDensityMatrix(qubit);
375  QRACK_CONST complex ONE_CMPLX_NEG = complex(-ONE_R1, ZERO_R1);
376  QRACK_CONST complex pauliZ[4]{ ONE_CMPLX, ZERO_CMPLX, ZERO_CMPLX, ONE_CMPLX_NEG };
377  complex pMtrx[4];
378  mul2x2(dMtrx.get(), pauliZ, pMtrx);
379  return (ONE_R1 - std::real(pMtrx[0] + pMtrx[1])) / 2;
380  }
381 
383  {
384  AntiCNOT(control, target);
385  const real1_f prob = ProbRdm(target);
386  AntiCNOT(control, target);
387 
388  return prob;
389  }
390 
392  {
393  CNOT(control, target);
394  const real1_f prob = ProbRdm(target);
395  CNOT(control, target);
396 
397  return prob;
398  }
399 
405  void SwitchToEngine();
406 
407  bool isClifford() { return !engine; }
408 
409  bool isClifford(bitLenInt qubit) { return !engine && !shards[qubit]; };
410 
411  bool isBinaryDecisionTree() { return engine && engine->isBinaryDecisionTree(); };
412 
413  using QInterface::Compose;
414  bitLenInt Compose(QStabilizerHybridPtr toCopy) { return ComposeEither(toCopy, false); };
415  bitLenInt Compose(QInterfacePtr toCopy) { return Compose(std::dynamic_pointer_cast<QStabilizerHybrid>(toCopy)); }
418  {
419  return Compose(std::dynamic_pointer_cast<QStabilizerHybrid>(toCopy), start);
420  }
421  bitLenInt ComposeNoClone(QStabilizerHybridPtr toCopy) { return ComposeEither(toCopy, true); };
423  {
424  return ComposeNoClone(std::dynamic_pointer_cast<QStabilizerHybrid>(toCopy));
425  }
426  bitLenInt ComposeEither(QStabilizerHybridPtr toCopy, bool willDestroy);
428  {
429  Decompose(start, std::dynamic_pointer_cast<QStabilizerHybrid>(dest));
430  }
431  void Decompose(bitLenInt start, QStabilizerHybridPtr dest);
433  void Dispose(bitLenInt start, bitLenInt length);
434  void Dispose(bitLenInt start, bitLenInt length, bitCapInt disposedPerm);
435  using QInterface::Allocate;
436  bitLenInt Allocate(bitLenInt start, bitLenInt length);
437 
438  void GetQuantumState(complex* outputState);
439  void GetProbs(real1* outputProbs);
440  complex GetAmplitude(bitCapInt perm) { return GetAmplitudeOrProb(perm, false); }
441  real1_f ProbAll(bitCapInt perm) { return (real1_f)norm(GetAmplitudeOrProb(perm, true)); }
442  void SetQuantumState(const complex* inputState);
444  {
445  SwitchToEngine();
446  engine->SetAmplitude(perm, amp);
447  }
448  void SetPermutation(bitCapInt perm, complex phaseFac = CMPLX_DEFAULT_ARG);
449 
450  void Swap(bitLenInt qubit1, bitLenInt qubit2);
451  void ISwap(bitLenInt qubit1, bitLenInt qubit2) { ISwapHelper(qubit1, qubit2, false); }
452  void IISwap(bitLenInt qubit1, bitLenInt qubit2) { ISwapHelper(qubit1, qubit2, true); }
453  void CSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
454  void CSqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
455  void AntiCSqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
456  void CISqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
457  void AntiCISqrtSwap(const std::vector<bitLenInt>& lControls, bitLenInt qubit1, bitLenInt qubit2);
458 
459  void XMask(bitCapInt mask);
460  void YMask(bitCapInt mask);
461  void ZMask(bitCapInt mask);
462 
463  real1_f Prob(bitLenInt qubit);
464 
465  bool ForceM(bitLenInt qubit, bool result, bool doForce = true, bool doApply = true);
466 
467  bitCapInt MAll();
468 
469  void Mtrx(const complex* mtrx, bitLenInt target);
470  void MCMtrx(const std::vector<bitLenInt>& controls, const complex* mtrx, bitLenInt target);
471  void MCPhase(const std::vector<bitLenInt>& controls, complex topLeft, complex bottomRight, bitLenInt target);
472  void MCInvert(const std::vector<bitLenInt>& controls, complex topRight, complex bottomLeft, bitLenInt target);
473  void MACMtrx(const std::vector<bitLenInt>& controls, const complex* mtrx, bitLenInt target);
474  void MACPhase(const std::vector<bitLenInt>& controls, complex topLeft, complex bottomRight, bitLenInt target);
475  void MACInvert(const std::vector<bitLenInt>& controls, complex topRight, complex bottomLeft, bitLenInt target);
476 
479  const std::vector<bitLenInt>& controls, bitLenInt qubitIndex, const complex* mtrxs);
480 
481  std::map<bitCapInt, int> MultiShotMeasureMask(const std::vector<bitCapInt>& qPowers, unsigned shots);
482  void MultiShotMeasureMask(const std::vector<bitCapInt>& qPowers, unsigned shots, unsigned long long* shotsArray);
483 
485  bool ForceMParity(bitCapInt mask, bool result, bool doForce = true);
486  void CUniformParityRZ(const std::vector<bitLenInt>& controls, bitCapInt mask, real1_f angle)
487  {
488  SwitchToEngine();
489  QINTERFACE_TO_QPARITY(engine)->CUniformParityRZ(controls, mask, angle);
490  }
491 
492 #if ENABLE_ALU
493  using QInterface::M;
494  bool M(bitLenInt q) { return QInterface::M(q); }
495  using QInterface::X;
496  void X(bitLenInt q) { QInterface::X(q); }
497  void CPhaseFlipIfLess(bitCapInt greaterPerm, bitLenInt start, bitLenInt length, bitLenInt flagIndex)
498  {
499  SwitchToEngine();
500  QINTERFACE_TO_QALU(engine)->CPhaseFlipIfLess(greaterPerm, start, length, flagIndex);
501  }
502  void PhaseFlipIfLess(bitCapInt greaterPerm, bitLenInt start, bitLenInt length)
503  {
504  SwitchToEngine();
505  QINTERFACE_TO_QALU(engine)->PhaseFlipIfLess(greaterPerm, start, length);
506  }
507 
508  void INC(bitCapInt toAdd, bitLenInt start, bitLenInt length)
509  {
510  if (stabilizer) {
511  QInterface::INC(toAdd, start, length);
512  return;
513  }
514 
515  engine->INC(toAdd, start, length);
516  }
517  void DEC(bitCapInt toSub, bitLenInt start, bitLenInt length)
518  {
519  if (stabilizer) {
520  QInterface::DEC(toSub, start, length);
521  return;
522  }
523 
524  engine->DEC(toSub, start, length);
525  }
526  void DECS(bitCapInt toSub, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
527  {
528  if (stabilizer) {
529  QInterface::DECS(toSub, start, length, overflowIndex);
530  return;
531  }
532 
533  engine->DECS(toSub, start, length, overflowIndex);
534  }
535  void CINC(bitCapInt toAdd, bitLenInt inOutStart, bitLenInt length, const std::vector<bitLenInt>& controls)
536  {
537  if (stabilizer) {
538  QInterface::CINC(toAdd, inOutStart, length, controls);
539  return;
540  }
541 
542  engine->CINC(toAdd, inOutStart, length, controls);
543  }
544  void INCS(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
545  {
546  if (stabilizer) {
547  QInterface::INCS(toAdd, start, length, overflowIndex);
548  return;
549  }
550 
551  engine->INCS(toAdd, start, length, overflowIndex);
552  }
553  void INCDECC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
554  {
555  if (stabilizer) {
556  QInterface::INCDECC(toAdd, start, length, carryIndex);
557  return;
558  }
559 
560  engine->INCDECC(toAdd, start, length, carryIndex);
561  }
562  void INCDECSC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex, bitLenInt carryIndex)
563  {
564  SwitchToEngine();
565  QINTERFACE_TO_QALU(engine)->INCDECSC(toAdd, start, length, overflowIndex, carryIndex);
566  }
567  void INCDECSC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
568  {
569  SwitchToEngine();
570  QINTERFACE_TO_QALU(engine)->INCDECSC(toAdd, start, length, carryIndex);
571  }
572 #if ENABLE_BCD
573  void INCBCD(bitCapInt toAdd, bitLenInt start, bitLenInt length)
574  {
575  SwitchToEngine();
576  QINTERFACE_TO_QALU(engine)->INCBCD(toAdd, start, length);
577  }
578  void INCDECBCDC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
579  {
580  SwitchToEngine();
581  QINTERFACE_TO_QALU(engine)->INCDECBCDC(toAdd, start, length, carryIndex);
582  }
583 #endif
584  void MUL(bitCapInt toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
585  {
586  SwitchToEngine();
587  QINTERFACE_TO_QALU(engine)->MUL(toMul, inOutStart, carryStart, length);
588  }
589  void DIV(bitCapInt toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
590  {
591  SwitchToEngine();
592  QINTERFACE_TO_QALU(engine)->DIV(toDiv, inOutStart, carryStart, length);
593  }
594  void MULModNOut(bitCapInt toMul, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
595  {
596  SwitchToEngine();
597  QINTERFACE_TO_QALU(engine)->MULModNOut(toMul, modN, inStart, outStart, length);
598  }
599  void IMULModNOut(bitCapInt toMul, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
600  {
601  SwitchToEngine();
602  QINTERFACE_TO_QALU(engine)->IMULModNOut(toMul, modN, inStart, outStart, length);
603  }
604  void POWModNOut(bitCapInt base, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
605  {
606  SwitchToEngine();
607  QINTERFACE_TO_QALU(engine)->POWModNOut(base, modN, inStart, outStart, length);
608  }
609  void CMUL(bitCapInt toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length,
610  const std::vector<bitLenInt>& controls)
611  {
612  SwitchToEngine();
613  QINTERFACE_TO_QALU(engine)->CMUL(toMul, inOutStart, carryStart, length, controls);
614  }
615  void CDIV(bitCapInt toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length,
616  const std::vector<bitLenInt>& controls)
617  {
618  SwitchToEngine();
619  QINTERFACE_TO_QALU(engine)->CDIV(toDiv, inOutStart, carryStart, length, controls);
620  }
621  void CMULModNOut(bitCapInt toMul, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length,
622  const std::vector<bitLenInt>& controls)
623  {
624  SwitchToEngine();
625  QINTERFACE_TO_QALU(engine)->CMULModNOut(toMul, modN, inStart, outStart, length, controls);
626  }
627  void CIMULModNOut(bitCapInt toMul, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length,
628  const std::vector<bitLenInt>& controls)
629  {
630  SwitchToEngine();
631  QINTERFACE_TO_QALU(engine)->CIMULModNOut(toMul, modN, inStart, outStart, length, controls);
632  }
633  void CPOWModNOut(bitCapInt base, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length,
634  const std::vector<bitLenInt>& controls)
635  {
636  SwitchToEngine();
637  QINTERFACE_TO_QALU(engine)->CPOWModNOut(base, modN, inStart, outStart, length, controls);
638  }
639 
640  bitCapInt IndexedLDA(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength,
641  const unsigned char* values, bool resetValue = true)
642  {
643  SwitchToEngine();
644  return QINTERFACE_TO_QALU(engine)->IndexedLDA(
645  indexStart, indexLength, valueStart, valueLength, values, resetValue);
646  }
647  bitCapInt IndexedADC(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength,
648  bitLenInt carryIndex, const unsigned char* values)
649  {
650  SwitchToEngine();
651  return QINTERFACE_TO_QALU(engine)->IndexedADC(
652  indexStart, indexLength, valueStart, valueLength, carryIndex, values);
653  }
654  bitCapInt IndexedSBC(bitLenInt indexStart, bitLenInt indexLength, bitLenInt valueStart, bitLenInt valueLength,
655  bitLenInt carryIndex, const unsigned char* values)
656  {
657  SwitchToEngine();
658  return QINTERFACE_TO_QALU(engine)->IndexedSBC(
659  indexStart, indexLength, valueStart, valueLength, carryIndex, values);
660  }
661  void Hash(bitLenInt start, bitLenInt length, const unsigned char* values)
662  {
663  SwitchToEngine();
664  QINTERFACE_TO_QALU(engine)->Hash(start, length, values);
665  }
666 #endif
667 
668  void PhaseFlip()
669  {
670  if (stabilizer) {
671  stabilizer->PhaseFlip();
672  } else {
673  engine->PhaseFlip();
674  }
675  }
676  void ZeroPhaseFlip(bitLenInt start, bitLenInt length)
677  {
678  SwitchToEngine();
679  engine->ZeroPhaseFlip(start, length);
680  }
681 
682  void SqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
683  {
684  if (stabilizer) {
685  QInterface::SqrtSwap(qubitIndex1, qubitIndex2);
686  return;
687  }
688 
689  SwitchToEngine();
690  engine->SqrtSwap(qubitIndex1, qubitIndex2);
691  }
692  void ISqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
693  {
694  if (stabilizer) {
695  QInterface::ISqrtSwap(qubitIndex1, qubitIndex2);
696  return;
697  }
698 
699  SwitchToEngine();
700  engine->ISqrtSwap(qubitIndex1, qubitIndex2);
701  }
702  void FSim(real1_f theta, real1_f phi, bitLenInt qubit1, bitLenInt qubit2)
703  {
704  const std::vector<bitLenInt> controls{ qubit1 };
705  real1 sinTheta = (real1)sin(theta);
706 
707  if ((sinTheta * sinTheta) <= FP_NORM_EPSILON) {
708  MCPhase(controls, ONE_CMPLX, exp(complex(ZERO_R1, (real1)phi)), qubit2);
709  return;
710  }
711 
712  const real1 sinThetaDiffNeg = ONE_R1 + sinTheta;
713  if ((sinThetaDiffNeg * sinThetaDiffNeg) <= FP_NORM_EPSILON) {
714  ISwap(qubit1, qubit2);
715  MCPhase(controls, ONE_CMPLX, exp(complex(ZERO_R1, (real1)phi)), qubit2);
716  return;
717  }
718 
719  const real1 sinThetaDiffPos = ONE_R1 - sinTheta;
720  if ((sinThetaDiffPos * sinThetaDiffPos) <= FP_NORM_EPSILON) {
721  IISwap(qubit1, qubit2);
722  MCPhase(controls, ONE_CMPLX, exp(complex(ZERO_R1, (real1)phi)), qubit2);
723  return;
724  }
725 
726  SwitchToEngine();
727  engine->FSim(theta, phi, qubit1, qubit2);
728  }
729 
730  real1_f ProbMask(bitCapInt mask, bitCapInt permutation)
731  {
732  SwitchToEngine();
733  return engine->ProbMask(mask, permutation);
734  }
735 
737  {
738  return ApproxCompareHelper(std::dynamic_pointer_cast<QStabilizerHybrid>(toCompare), false);
739  }
741  {
742  return error_tol >=
743  ApproxCompareHelper(std::dynamic_pointer_cast<QStabilizerHybrid>(toCompare), true, error_tol);
744  }
745 
747  {
748  if (engine) {
749  engine->UpdateRunningNorm(norm_thresh);
750  }
751  }
752 
753  void NormalizeState(
754  real1_f nrm = REAL1_DEFAULT_ARG, real1_f norm_thresh = REAL1_DEFAULT_ARG, real1_f phaseArg = ZERO_R1_F);
755 
756  real1_f ProbAllRdm(bool roundRz, bitCapInt fullRegister);
757  real1_f ProbMaskRdm(bool roundRz, bitCapInt mask, bitCapInt permutation);
758  real1_f ExpectationBitsAll(const std::vector<bitLenInt>& bits, bitCapInt offset = 0)
759  {
760  if (stabilizer) {
761  return QInterface::ExpectationBitsAll(bits, offset);
762  }
763 
764  return engine->ExpectationBitsAll(bits, offset);
765  }
766  real1_f ExpectationBitsAllRdm(bool roundRz, const std::vector<bitLenInt>& bits, bitCapInt offset = 0U)
767  {
768  if (engine) {
769  return engine->ExpectationBitsAllRdm(roundRz, bits, offset);
770  }
771 
772  CombineAncillae();
773 
774  if (!roundRz) {
775  return stabilizer->ExpectationBitsAll(bits, offset);
776  }
777 
778  return RdmCloneHelper()->stabilizer->ExpectationBitsAll(bits, offset);
779  }
781  const std::vector<bitLenInt>& bits, const std::vector<bitCapInt>& perms, bitCapInt offset = 0U)
782  {
783  if (stabilizer) {
784  return QInterface::ExpectationBitsFactorized(bits, perms, offset);
785  }
786 
787  return engine->ExpectationBitsFactorized(bits, perms, offset);
788  }
790  bool roundRz, const std::vector<bitLenInt>& bits, const std::vector<bitCapInt>& perms, bitCapInt offset = 0U)
791  {
792  return ExpectationFactorized(false, bits, perms, std::vector<real1_f>(), offset, roundRz);
793  }
794  real1_f ExpectationFloatsFactorized(const std::vector<bitLenInt>& bits, const std::vector<real1_f>& weights)
795  {
796  if (stabilizer) {
797  return QInterface::ExpectationFloatsFactorized(bits, weights);
798  }
799 
800  return engine->ExpectationFloatsFactorized(bits, weights);
801  }
803  bool roundRz, const std::vector<bitLenInt>& bits, const std::vector<real1_f>& weights)
804  {
805  return ExpectationFactorized(true, bits, std::vector<bitCapInt>(), weights, 0U, roundRz);
806  }
807 
808  bool TrySeparate(bitLenInt qubit);
809  bool TrySeparate(bitLenInt qubit1, bitLenInt qubit2);
810  bool TrySeparate(const std::vector<bitLenInt>& qubits, real1_f error_tol);
811 
813 
814  void SetDevice(int64_t dID)
815  {
816  devID = dID;
817  if (engine) {
818  engine->SetDevice(dID);
819  }
820  }
821 
822  int64_t GetDeviceID() { return devID; }
823 
825  {
826  if (stabilizer) {
827  return QInterface::GetMaxSize();
828  }
829 
830  return engine->GetMaxSize();
831  }
832 
833  friend std::ostream& operator<<(std::ostream& os, const QStabilizerHybridPtr s);
834  friend std::istream& operator>>(std::istream& is, const QStabilizerHybridPtr s);
835 };
836 } // namespace Qrack
unsigned GetConcurrencyLevel()
Definition: parallel_for.hpp:38
Definition: qalu.hpp:22
A "Qrack::QInterface" is an abstract interface exposing qubit permutation state vector with methods t...
Definition: qinterface.hpp:146
virtual void SetConcurrency(uint32_t threadsPerEngine)
Set the number of threads in parallel for loops, per component QEngine.
Definition: qinterface.hpp:249
virtual bitLenInt Allocate(bitLenInt length)
Allocate new "length" count of |0> state qubits at end of qubit index position.
Definition: qinterface.hpp:434
virtual bitLenInt Compose(QInterfacePtr toCopy)
Combine another QInterface with this one, after the last bit index of this one.
Definition: qinterface.hpp:338
bitLenInt qubitCount
Definition: qinterface.hpp:151
real1_f Rand()
Generate a random real number between 0 and 1.
Definition: qinterface.hpp:260
Definition: qparity.hpp:22
A "Qrack::QStabilizerHybrid" internally switched between Qrack::QStabilizer and Qrack::QEngine to max...
Definition: qstabilizerhybrid.hpp:43
bool TrimControls(const std::vector< bitLenInt > &lControls, std::vector< bitLenInt > &output, bool anti=false)
Definition: qstabilizerhybrid.cpp:270
complex GetAmplitude(bitCapInt perm)
Get the representational amplitude of a full permutation.
Definition: qstabilizerhybrid.hpp:440
void MACPhase(const std::vector< bitLenInt > &controls, complex topLeft, complex bottomRight, bitLenInt target)
Apply a single bit transformation that only effects phase, with arbitrary (anti-)control bits.
Definition: qstabilizerhybrid.cpp:1266
void MUL(bitCapInt toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
Multiply by integer.
Definition: qstabilizerhybrid.hpp:584
bitLenInt ComposeNoClone(QInterfacePtr toCopy)
Definition: qstabilizerhybrid.hpp:422
bitLenInt maxEngineQubitCount
Definition: qstabilizerhybrid.hpp:56
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:1240
void SetAmplitude(bitCapInt perm, complex amp)
Sets the representational amplitude of a full permutation.
Definition: qstabilizerhybrid.hpp:443
real1_f ProbAllRdm(bool roundRz, bitCapInt fullRegister)
Direct measure of full permutation probability, treating all ancillary qubits as post-selected T gate...
Definition: qstabilizerhybrid.cpp:378
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:959
bitLenInt ComposeEither(QStabilizerHybridPtr toCopy, bool willDestroy)
Definition: qstabilizerhybrid.cpp:483
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:452
void SetTInjection(bool useGadget)
Set the option to use T-injection gadgets (off by default)
Definition: qstabilizerhybrid.hpp:337
real1_f FractionalRzAngleWithFlush(bitLenInt i, real1_f angle, bool isGateSuppressed=false)
Definition: qstabilizerhybrid.hpp:190
int64_t devID
Definition: qstabilizerhybrid.hpp:60
real1_f ProbMaskRdm(bool roundRz, bitCapInt mask, bitCapInt permutation)
Direct measure of masked permutation probability, treating all ancillary qubits as post-selected T ga...
Definition: qstabilizerhybrid.cpp:391
bitLenInt maxStateMapCacheQubitCount
Definition: qstabilizerhybrid.hpp:58
void INC(bitCapInt toAdd, bitLenInt start, bitLenInt length)
Add integer (without sign)
Definition: qstabilizerhybrid.hpp:508
bool isFinished()
Returns "false" if asynchronous work is still running, and "true" if all previously dispatched asynch...
Definition: qstabilizerhybrid.hpp:349
QStabilizerHybridPtr RdmCloneHelper()
Definition: qstabilizerhybrid.hpp:264
void FlushBuffers()
Definition: qstabilizerhybrid.cpp:251
void XMask(bitCapInt mask)
Masked X gate.
Definition: qstabilizerhybrid.cpp:992
std::vector< int64_t > deviceIDs
Definition: qstabilizerhybrid.hpp:64
bool doNormalize
Definition: qstabilizerhybrid.hpp:49
friend std::ostream & operator<<(std::ostream &os, const QStabilizerHybridPtr s)
Definition: qstabilizerhybrid.cpp:2182
void SqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
Square root of Swap gate.
Definition: qstabilizerhybrid.hpp:682
void MULModNOut(bitCapInt toMul, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
Multiplication modulo N by integer, (out of place)
Definition: qstabilizerhybrid.hpp:594
virtual bitLenInt Allocate(bitLenInt length)
Allocate new "length" count of |0> state qubits at end of qubit index position.
Definition: qinterface.hpp:434
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:943
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:368
void Swap(bitLenInt qubit1, bitLenInt qubit2)
Swap values of two bits in register.
Definition: qstabilizerhybrid.cpp:897
void INCDECSC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCSC and DECSC (without overflow flag)
Definition: qstabilizerhybrid.hpp:567
void GetProbs(real1 *outputProbs)
Get the pure quantum state representation.
Definition: qstabilizerhybrid.cpp:670
void RdmCloneFlush(real1_f threshold=FP_NORM_EPSILON)
Flush non-Clifford phase gate gadgets with angle below a threshold.
Definition: qstabilizerhybrid.cpp:1936
void CSqrtSwap(const std::vector< bitLenInt > &lControls, bitLenInt qubit1, bitLenInt qubit2)
Apply a square root of swap with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:927
void PhaseFlip()
Phase flip always - equivalent to Z X Z X on any bit in the QInterface.
Definition: qstabilizerhybrid.hpp:668
void FlushIfBlocked(bitLenInt control, bitLenInt target, bool isPhase=false)
Definition: qstabilizerhybrid.cpp:139
bitCapInt MAll()
Measure permutation state of all coherent bits.
Definition: qstabilizerhybrid.cpp:1528
virtual bitLenInt Compose(QInterfacePtr toCopy)
Combine another QInterface with this one, after the last bit index of this one.
Definition: qinterface.hpp:338
real1_f ExpectationBitsFactorized(const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, bitCapInt offset=0U)
Get expectation value of bits, given an array of qubit weights.
Definition: qstabilizerhybrid.hpp:780
real1_f ProbAll(bitCapInt perm)
Direct measure of full permutation probability.
Definition: qstabilizerhybrid.hpp:441
bool isClifford()
Returns "true" if current state is identifiably within the Clifford set, or "false" if it is not or c...
Definition: qstabilizerhybrid.hpp:407
void FlushCliffordFromBuffers()
Definition: qstabilizerhybrid.hpp:230
void ClearAncilla(bitLenInt i)
Definition: qstabilizerhybrid.hpp:294
bool useTGadget
Definition: qstabilizerhybrid.hpp:51
QUnitCliffordPtr MakeStabilizer(bitCapInt perm=0U)
Definition: qstabilizerhybrid.cpp:98
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:647
void CDIV(bitCapInt toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Controlled division by power of integer.
Definition: qstabilizerhybrid.hpp:615
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:802
void DECS(bitCapInt toSub, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Add a classical integer to the register, with sign and without carry.
Definition: qstabilizerhybrid.hpp:526
void MACInvert(const std::vector< bitLenInt > &controls, complex topRight, complex bottomLeft, bitLenInt target)
Apply a single bit transformation that reverses bit probability and might effect phase,...
Definition: qstabilizerhybrid.cpp:1315
void CacheEigenstate(bitLenInt target)
Definition: qstabilizerhybrid.cpp:310
void InvertBuffer(bitLenInt qubit)
Definition: qstabilizerhybrid.cpp:120
bitLenInt ancillaCount
Definition: qstabilizerhybrid.hpp:54
void DumpBuffers()
Definition: qstabilizerhybrid.hpp:81
real1_f SumSqrDiff(QInterfacePtr toCompare)
Definition: qstabilizerhybrid.hpp:736
std::vector< MpsShardPtr > shards
Definition: qstabilizerhybrid.hpp:67
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:654
real1_f ACProbRdm(bitLenInt control, bitLenInt target)
Definition: qstabilizerhybrid.hpp:391
bool IsBuffered()
Definition: qstabilizerhybrid.hpp:99
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:1113
QInterfacePtr MakeEngine(bitCapInt perm=0U)
Definition: qstabilizerhybrid.cpp:103
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:746
bitLenInt maxAncillaCount
Definition: qstabilizerhybrid.hpp:57
virtual void UniformlyControlledSingleBit(const std::vector< bitLenInt > &controls, bitLenInt qubitIndex, const complex *mtrxs)
Apply a "uniformly controlled" arbitrary single bit unitary transformation.
Definition: qinterface.hpp:590
bool isSparse
Definition: qstabilizerhybrid.hpp:50
real1_f ProbParity(bitCapInt mask)
Overall probability of any odd permutation of the masked set of bits.
Definition: qstabilizerhybrid.cpp:1802
bitLenInt Compose(QInterfacePtr toCopy, bitLenInt start)
Definition: qstabilizerhybrid.hpp:417
bitLenInt ComposeNoClone(QStabilizerHybridPtr toCopy)
Definition: qstabilizerhybrid.hpp:421
bool TrySeparate(bitLenInt qubit)
Single-qubit TrySeparate()
Definition: qstabilizerhybrid.cpp:2141
QInterfacePtr engine
Definition: qstabilizerhybrid.hpp:62
real1_f ExpectationBitsAllRdm(bool roundRz, const std::vector< bitLenInt > &bits, bitCapInt offset=0U)
Get permutation expectation value of bits, treating all ancillary qubits as post-selected T gate gadg...
Definition: qstabilizerhybrid.hpp:766
void INCDECC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCC and DECC (without sign, with carry)
Definition: qstabilizerhybrid.hpp:553
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:451
void X(bitLenInt q)
Definition: qstabilizerhybrid.hpp:496
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:640
std::vector< QInterfaceEngine > cloneEngineTypes
Definition: qstabilizerhybrid.hpp:66
void ISwapHelper(bitLenInt qubit1, bitLenInt qubit2, bool inverse)
Definition: qstabilizerhybrid.cpp:2085
std::map< bitCapInt, int > MultiShotMeasureMask(const std::vector< bitCapInt > &qPowers, unsigned shots)
Statistical measure of masked permutation probability.
Definition: qstabilizerhybrid.cpp:1636
bitLenInt Compose(QStabilizerHybridPtr toCopy)
Definition: qstabilizerhybrid.hpp:414
void Mtrx(const complex *mtrx, bitLenInt target)
Apply an arbitrary single bit unitary transformation.
Definition: qstabilizerhybrid.cpp:1035
void MCInvert(const std::vector< bitLenInt > &controls, complex topRight, complex bottomLeft, bitLenInt target)
Apply a single bit transformation that reverses bit probability and might effect phase,...
Definition: qstabilizerhybrid.cpp:1192
bool EitherIsBuffered(bool logical)
Definition: qstabilizerhybrid.hpp:87
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:794
bool M(bitLenInt q)
Definition: qstabilizerhybrid.hpp:494
real1_f ExpectationBitsFactorizedRdm(bool roundRz, const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, bitCapInt offset=0U)
Get (reduced density matrix) expectation value of bits, given an array of qubit weights.
Definition: qstabilizerhybrid.hpp:789
bitCapIntOcl GetMaxSize()
Definition: qstabilizerhybrid.hpp:824
bool useHostRam
Definition: qstabilizerhybrid.hpp:48
void CUniformParityRZ(const std::vector< bitLenInt > &controls, 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:486
bool GetTInjection()
Get the option to use T-injection gadgets.
Definition: qstabilizerhybrid.hpp:338
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:975
void ZeroPhaseFlip(bitLenInt start, bitLenInt length)
Reverse the phase of the state where the register equals zero.
Definition: qstabilizerhybrid.hpp:676
QInterfacePtr Clone()
Clone this QInterface.
Definition: qstabilizerhybrid.cpp:350
void MCPhase(const std::vector< bitLenInt > &controls, complex topLeft, complex bottomRight, bitLenInt target)
Apply a single bit transformation that only effects phase, with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:1139
void INCS(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Add a classical integer to the register, with sign and without carry.
Definition: qstabilizerhybrid.hpp:544
bool IsLogicalBuffered()
Definition: qstabilizerhybrid.hpp:100
void CSwap(const std::vector< bitLenInt > &lControls, bitLenInt qubit1, bitLenInt qubit2)
Apply a swap with arbitrary control bits.
Definition: qstabilizerhybrid.cpp:911
void DIV(bitCapInt toDiv, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length)
Divide by integer.
Definition: qstabilizerhybrid.hpp:589
void SetDevice(int64_t dID)
Set the device index, if more than one device is available.
Definition: qstabilizerhybrid.hpp:814
void CPOWModNOut(bitCapInt base, 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:633
real1_f ApproxCompareHelper(QStabilizerHybridPtr toCompare, bool isDiscreteBool, real1_f error_tol=TRYDECOMPOSE_EPSILON)
Definition: qstabilizerhybrid.cpp:2002
void ISqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
Inverse square root of Swap gate.
Definition: qstabilizerhybrid.hpp:692
void CMULModNOut(bitCapInt toMul, 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:621
bool isBinaryDecisionTree()
Returns "true" if current state representation is definitely a binary decision tree,...
Definition: qstabilizerhybrid.hpp:411
bool EitherIsProbBuffered(bool logical)
Definition: qstabilizerhybrid.hpp:101
void CPhaseFlipIfLess(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:497
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:613
bool isRoundingFlushed
Definition: qstabilizerhybrid.hpp:52
real1_f ProbMask(bitCapInt mask, bitCapInt permutation)
Direct measure of masked permutation probability.
Definition: qstabilizerhybrid.hpp:730
real1_f separabilityThreshold
Definition: qstabilizerhybrid.hpp:59
complex GetAmplitudeOrProb(bitCapInt perm, bool isProb=false)
Definition: qstabilizerhybrid.cpp:686
void Decompose(bitLenInt start, QInterfacePtr dest)
Minimally decompose a set of contiguous bits from the separably composed unit, into "destination".
Definition: qstabilizerhybrid.hpp:427
void YMask(bitCapInt mask)
Masked Y gate.
Definition: qstabilizerhybrid.cpp:1006
real1_f ExpectationBitsAll(const std::vector< bitLenInt > &bits, bitCapInt offset=0)
Get permutation expectation value of bits.
Definition: qstabilizerhybrid.hpp:758
QUnitCliffordPtr stabilizer
Definition: qstabilizerhybrid.hpp:63
void SetQuantumState(const complex *inputState)
Set an arbitrary pure quantum state representation.
Definition: qstabilizerhybrid.cpp:841
void IMULModNOut(bitCapInt toMul, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
Inverse of multiplication modulo N by integer, (out of place)
Definition: qstabilizerhybrid.hpp:599
bitLenInt deadAncillaCount
Definition: qstabilizerhybrid.hpp:55
void CMUL(bitCapInt toMul, bitLenInt inOutStart, bitLenInt carryStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Controlled multiplication by integer.
Definition: qstabilizerhybrid.hpp:609
bool ForceMParity(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:1815
bitLenInt thresholdQubits
Definition: qstabilizerhybrid.hpp:53
QStabilizerHybrid(std::vector< QInterfaceEngine > eng, bitLenInt qBitCount, bitCapInt initState=0U, qrack_rand_gen_ptr rgp=nullptr, complex phaseFac=CMPLX_DEFAULT_ARG, bool doNorm=false, bool randomGlobalPhase=true, bool useHostMem=false, int64_t deviceId=-1, bool useHardwareRNG=true, bool useSparseStateVec=false, real1_f norm_thresh=REAL1_EPSILON, std::vector< int64_t > devList={}, bitLenInt qubitThreshold=0U, real1_f separation_thresh=FP_NORM_EPSILON_F)
Definition: qstabilizerhybrid.cpp:35
std::map< bitCapInt, complex > stateMapCache
Definition: qstabilizerhybrid.hpp:68
void SwitchToEngine()
Switches between CPU and GPU used modes.
Definition: qstabilizerhybrid.cpp:408
std::vector< QInterfaceEngine > engineTypes
Definition: qstabilizerhybrid.hpp:65
void CheckShots(unsigned shots, bitCapInt m, real1_f partProb, const std::vector< bitCapInt > &qPowers, std::vector< real1_f > &rng, F fn)
Definition: qstabilizerhybrid.hpp:154
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:2128
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:409
void Finish()
If asynchronous work is still running, block until it finishes.
Definition: qstabilizerhybrid.hpp:340
bool CollapseSeparableShard(bitLenInt qubit)
Definition: qstabilizerhybrid.cpp:227
bool IsLogicalProbBuffered()
Definition: qstabilizerhybrid.hpp:124
friend std::istream & operator>>(std::istream &is, const QStabilizerHybridPtr s)
Definition: qstabilizerhybrid.cpp:2210
void INCDECBCDC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCSC and DECSC (without overflow flag)
Definition: qstabilizerhybrid.hpp:578
void SetPermutation(bitCapInt perm, complex phaseFac=CMPLX_DEFAULT_ARG)
Set to a specific permutation of all qubits.
Definition: qstabilizerhybrid.cpp:880
real1_f Prob(bitLenInt qubit)
Direct measure of bit probability to be in |1> state.
Definition: qstabilizerhybrid.cpp:1363
std::unique_ptr< complex[]> GetQubitReducedDensityMatrix(bitLenInt qubit)
Definition: qstabilizerhybrid.hpp:126
void DEC(bitCapInt toSub, bitLenInt start, bitLenInt length)
Add integer (without sign)
Definition: qstabilizerhybrid.hpp:517
QStabilizerHybrid(bitLenInt qBitCount, bitCapInt initState=0U, qrack_rand_gen_ptr rgp=nullptr, complex phaseFac=CMPLX_DEFAULT_ARG, bool doNorm=false, bool randomGlobalPhase=true, bool useHostMem=false, int64_t deviceId=-1, bool useHardwareRNG=true, bool useSparseStateVec=false, real1_f norm_thresh=REAL1_EPSILON, std::vector< int64_t > devList={}, bitLenInt qubitThreshold=0U, real1_f separation_thresh=FP_NORM_EPSILON_F)
Definition: qstabilizerhybrid.hpp:326
real1_f ExpectationFactorized(bool isFloat, const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, const std::vector< real1_f > &weights, bitCapInt offset, bool roundRz)
Definition: qstabilizerhybrid.hpp:274
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:1447
void GetQuantumState(complex *outputState)
Get the pure quantum state representation.
Definition: qstabilizerhybrid.cpp:653
void INCBCD(bitCapInt toAdd, bitLenInt start, bitLenInt length)
Add classical BCD integer (without sign)
Definition: qstabilizerhybrid.hpp:573
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:661
void FlushH(bitLenInt qubit)
Definition: qstabilizerhybrid.cpp:129
void CombineAncillae()
Definition: qstabilizerhybrid.cpp:1831
void ZMask(bitCapInt mask)
Masked Z gate.
Definition: qstabilizerhybrid.cpp:1020
bool ApproxCompare(QInterfacePtr toCompare, real1_f error_tol=TRYDECOMPOSE_EPSILON)
Compare state vectors approximately, component by component, to determine whether this state vector i...
Definition: qstabilizerhybrid.hpp:740
void POWModNOut(bitCapInt base, bitCapInt modN, bitLenInt inStart, bitLenInt outStart, bitLenInt length)
Raise a classical base to a quantum power, modulo N, (out of place)
Definition: qstabilizerhybrid.hpp:604
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:702
void SetConcurrency(uint32_t threadCount)
Set the number of threads in parallel for loops, per component QEngine.
Definition: qstabilizerhybrid.hpp:360
void CIMULModNOut(bitCapInt toMul, 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:627
bitLenInt Compose(QInterfacePtr toCopy)
Combine another QInterface with this one, after the last bit index of this one.
Definition: qstabilizerhybrid.hpp:415
real1_f CProbRdm(bitLenInt control, bitLenInt target)
Definition: qstabilizerhybrid.hpp:382
void CINC(bitCapInt toAdd, bitLenInt inOutStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Add integer (without sign, with controls)
Definition: qstabilizerhybrid.hpp:535
void PhaseFlipIfLess(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:502
std::vector< real1_f > GenerateShotProbs(unsigned shots)
Definition: qstabilizerhybrid.hpp:177
void Dump()
If asynchronous work is still running, let the simulator know that it can be aborted.
Definition: qstabilizerhybrid.hpp:351
void INCDECSC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex, bitLenInt carryIndex)
Common driver method behind INCSC and DECSC (with overflow flag)
Definition: qstabilizerhybrid.hpp:562
int64_t GetDeviceID()
Definition: qstabilizerhybrid.hpp:822
bool IsProbBuffered()
Definition: qstabilizerhybrid.hpp:123
Half-precision floating-point type.
Definition: half.hpp:2222
virtual void INCDECC(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt carryIndex)
Common driver method behind INCC and DECC.
Definition: arithmetic.cpp:63
virtual void CINC(bitCapInt toAdd, bitLenInt inOutStart, bitLenInt length, const std::vector< bitLenInt > &controls)
Add integer (without sign, with controls)
Definition: arithmetic.cpp:114
virtual void INCS(bitCapInt toAdd, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Add a classical integer to the register, with sign and without carry.
Definition: arithmetic.cpp:167
virtual void DEC(bitCapInt toSub, bitLenInt start, bitLenInt length)
Subtract classical integer (without sign)
Definition: arithmetic.cpp:57
virtual void INC(bitCapInt toAdd, bitLenInt start, bitLenInt length)
Add integer (without sign)
Definition: arithmetic.cpp:24
virtual void DECS(bitCapInt toSub, bitLenInt start, bitLenInt length, bitLenInt overflowIndex)
Subtract a classical integer from the register, with sign and without carry.
Definition: arithmetic.cpp:182
virtual void CNOT(bitLenInt control, bitLenInt target)
Controlled NOT gate.
Definition: qinterface.hpp:672
virtual bool M(bitLenInt qubitIndex)
Measurement gate.
Definition: qinterface.hpp:976
virtual void UniformlyControlledSingleBit(const std::vector< bitLenInt > &controls, bitLenInt qubitIndex, const complex *mtrxs)
Apply a "uniformly controlled" arbitrary single bit unitary transformation.
Definition: qinterface.hpp:590
virtual void X(bitLenInt qubit)
X gate.
Definition: qinterface.hpp:1054
virtual void AntiCNOT(bitLenInt control, bitLenInt target)
Anti controlled NOT gate.
Definition: qinterface.hpp:683
virtual void U(bitLenInt target, real1_f theta, real1_f phi, real1_f lambda)
General unitary gate.
Definition: rotational.cpp:18
virtual void ISqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
Inverse square root of Swap gate.
Definition: gates.cpp:211
virtual void SqrtSwap(bitLenInt qubitIndex1, bitLenInt qubitIndex2)
Square root of Swap gate.
Definition: gates.cpp:188
virtual real1_f ExpectationBitsAll(const std::vector< bitLenInt > &bits, bitCapInt offset=0U)
Get permutation expectation value of bits.
Definition: qinterface.hpp:2394
bitCapIntOcl GetMaxSize()
Get maximum number of amplitudes that can be allocated on current device.
Definition: qinterface.hpp:2704
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:510
virtual real1_f ExpectationBitsFactorized(const std::vector< bitLenInt > &bits, const std::vector< bitCapInt > &perms, bitCapInt offset=0U)
Get expectation value of bits, given an array of qubit weights.
Definition: qinterface.cpp:469
Definition: complex16x2simd.hpp:25
std::complex< half_float::half > complex
Definition: qrack_types.hpp:62
@ QINTERFACE_OPTIMAL_BASE
Definition: qinterface.hpp:129
std::shared_ptr< QInterface > QInterfacePtr
Definition: qinterface.hpp:28
QRACK_CONST real1_f TRYDECOMPOSE_EPSILON
Definition: qrack_types.hpp:244
std::shared_ptr< QStabilizerHybrid > QStabilizerHybridPtr
Definition: qstabilizerhybrid.hpp:35
constexpr real1_f ZERO_R1_F
Definition: qrack_types.hpp:152
constexpr real1_f FP_NORM_EPSILON_F
Definition: qrack_types.hpp:245
const real1 ONE_R1
Definition: qrack_types.hpp:153
void mul2x2(complex const *left, complex const *right, complex *out)
Definition: functions.cpp:97
constexpr real1_f ONE_R1_F
Definition: qrack_types.hpp:154
half_float::half real1
Definition: qrack_types.hpp:63
QRACK_CONST real1 FP_NORM_EPSILON
Definition: qrack_types.hpp:243
bitCapInt pow2(const bitLenInt &p)
Definition: qrack_functions.hpp:22
void reverse(BidirectionalIterator first, BidirectionalIterator last, bitCapInt stride)
std::shared_ptr< QUnitClifford > QUnitCliffordPtr
Definition: qunitclifford.hpp:19
double norm(const complex2 &c)
Definition: complex16x2simd.hpp:101
const real1 REAL1_DEFAULT_ARG
Definition: qrack_types.hpp:155
const real1 PI_R1
Definition: qrack_types.hpp:158
QRACK_CONST complex ONE_CMPLX
Definition: qrack_types.hpp:239
const real1 ZERO_R1
Definition: qrack_types.hpp:151
float real1_f
Definition: qrack_types.hpp:64
float real1_s
Definition: qrack_types.hpp:65
QRACK_CONST complex CMPLX_DEFAULT_ARG
Definition: qrack_types.hpp:242
std::shared_ptr< MpsShard > MpsShardPtr
Definition: mpsshard.hpp:18
QRACK_CONST complex I_CMPLX
Definition: qrack_types.hpp:241
const real1 REAL1_EPSILON
Definition: qrack_types.hpp:157
QRACK_CONST complex ZERO_CMPLX
Definition: qrack_types.hpp:240
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
half round(half arg)
Nearest integer.
Definition: half.hpp:4495
half exp(half arg)
Exponential function.
Definition: half.hpp:3201
#define QRACK_CONST
Definition: qrack_types.hpp:150
#define bitLenInt
Definition: qrack_types.hpp:44
#define qrack_rand_gen_ptr
Definition: qrack_types.hpp:146
#define bitCapInt
Definition: qrack_types.hpp:105
#define bitCapIntOcl
Definition: qrack_types.hpp:91
#define QINTERFACE_TO_QALU(qReg)
Definition: qstabilizerhybrid.hpp:18
#define QINTERFACE_TO_QPARITY(qReg)
Definition: qstabilizerhybrid.hpp:19
Definition: qstabilizerhybrid.hpp:23
QUnitCliffordAmp(complex a, QUnitCliffordPtr s)
Definition: qstabilizerhybrid.hpp:27
QUnitCliffordPtr stabilizer
Definition: qstabilizerhybrid.hpp:25
complex amp
Definition: qstabilizerhybrid.hpp:24