QuEST_cpu_internal.h
Go to the documentation of this file.
1 // Distributed under MIT licence. See https://github.com/QuEST-Kit/QuEST/blob/master/LICENCE.txt for details
2 
13 # ifndef QUEST_CPU_INTERNAL_H
14 # define QUEST_CPU_INTERNAL_H
15 
16 # include "QuEST_precision.h"
17 
18 
19 /*
20 * Bit twiddling functions are defined seperately here in the CPU backend,
21 * since the GPU backend needs a device-specific redefinition to be callable
22 * from GPU kernels. These are called in both QuEST_cpu and QuEST_cpu_distributed
23 * and defined in here since public inline methods in C must go in the header
24 */
25 
26 static inline int extractBit (const int locationOfBitFromRight, const long long int theEncodedNumber) {
27  return (theEncodedNumber & ( 1LL << locationOfBitFromRight )) >> locationOfBitFromRight;
28 }
29 
30 static inline long long int flipBit(const long long int number, const int bitInd) {
31  return (number ^ (1LL << bitInd));
32 }
33 
34 static inline int maskContainsBit(const long long int mask, const int bitInd) {
35  return mask & (1LL << bitInd);
36 }
37 
38 static inline int isOddParity(const long long int number, const int qb1, const int qb2) {
39  return extractBit(qb1, number) != extractBit(qb2, number);
40 }
41 
42 static inline long long int insertZeroBit(const long long int number, const int index) {
43  long long int left, right;
44  left = (number >> index) << index;
45  right = number - left;
46  return (left << 1) ^ right;
47 }
48 
49 static inline long long int insertTwoZeroBits(const long long int number, const int bit1, const int bit2) {
50  int small = (bit1 < bit2)? bit1 : bit2;
51  int big = (bit1 < bit2)? bit2 : bit1;
52  return insertZeroBit(insertZeroBit(number, small), big);
53 }
54 
55 
56 /*
57  * density matrix operations
58  */
59 
61 
62 void densmatr_initPureStateLocal(Qureg targetQureg, Qureg copyQureg);
63 
65 
67 
69 
70 qreal densmatr_findProbabilityOfZeroLocal(Qureg qureg, int measureQubit);
71 
72 void densmatr_mixDepolarisingLocal(Qureg qureg, int targetQubit, qreal depolLevel);
73 
74 void densmatr_mixDepolarisingDistributed(Qureg qureg, int targetQubit, qreal depolLevel);
75 
76 void densmatr_mixDampingLocal(Qureg qureg, int targetQubit, qreal damping);
77 
78 void densmatr_mixDampingDistributed(Qureg qureg, int targetQubit, qreal damping);
79 
80 void densmatr_mixTwoQubitDepolarisingLocal(Qureg qureg, int qubit1, int qubit2, qreal delta, qreal gamma);
81 
82 void densmatr_mixTwoQubitDepolarisingLocalPart1(Qureg qureg, int qubit1, int qubit2, qreal delta);
83 
84 void densmatr_mixTwoQubitDepolarisingDistributed(Qureg qureg, int targetQubit,
85  int qubit2, qreal delta, qreal gamma);
86 
88  int qubit2, qreal delta, qreal gamma);
89 
91 
93 
94 void densmatr_calcProbOfAllOutcomesLocal(qreal* retProbs, Qureg qureg, int* qubits, int numQubits);
95 
96 
97 /*
98  * state vector operations
99  */
100 
102 
103 void statevec_compactUnitaryLocal (Qureg qureg, int targetQubit, Complex alpha, Complex beta);
104 
106  Complex rot1, Complex rot2,
107  ComplexArray stateVecUp,
108  ComplexArray stateVecLo,
109  ComplexArray stateVecOut);
110 
111 void statevec_unitaryLocal(Qureg qureg, int targetQubit, ComplexMatrix2 u);
112 
114  Complex rot1, Complex rot2,
115  ComplexArray stateVecUp,
116  ComplexArray stateVecLo,
117  ComplexArray stateVecOut);
118 
119 void statevec_controlledCompactUnitaryLocal (Qureg qureg, int controlQubit, int targetQubit,
120  Complex alpha, Complex beta);
121 
122 void statevec_controlledCompactUnitaryDistributed (Qureg qureg, int controlQubit,
123  Complex rot1, Complex rot2,
124  ComplexArray stateVecUp,
125  ComplexArray stateVecLo,
126  ComplexArray stateVecOut);
127 
128 void statevec_controlledUnitaryLocal(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u);
129 
130 void statevec_controlledUnitaryDistributed (Qureg qureg, int controlQubit,
131  Complex rot1, Complex rot2,
132  ComplexArray stateVecUp,
133  ComplexArray stateVecLo,
134  ComplexArray stateVecOut);
135 
136 void statevec_multiControlledUnitaryLocal(Qureg qureg, int targetQubit,
137  long long int ctrlQubitsMask, long long int ctrlFlipMask, ComplexMatrix2 u);
138 
140  int targetQubit,
141  long long int ctrlQubitsMask, long long int ctrlFlipMask,
142  Complex rot1, Complex rot2,
143  ComplexArray stateVecUp,
144  ComplexArray stateVecLo,
145  ComplexArray stateVecOut);
146 
147 void statevec_pauliXLocal(Qureg qureg, int targetQubit);
148 
150  ComplexArray stateVecIn,
151  ComplexArray stateVecOut);
152 
153 void statevec_pauliYLocal(Qureg qureg, int targetQubit, int conjFac);
154 
156  ComplexArray stateVecIn,
157  ComplexArray stateVecOut,
158  int updateUpper, int conjFac);
159 
160 void statevec_controlledPauliYLocal(Qureg qureg, int controlQubit, int targetQubit, int conjFactor);
161 
162 void statevec_controlledPauliYDistributed(Qureg qureg, int controlQubit,
163  ComplexArray stateVecIn,
164  ComplexArray stateVecOut, int conjFactor);
165 
166 void statevec_hadamardLocal (Qureg qureg, int targetQubit);
167 
169  ComplexArray stateVecUp,
170  ComplexArray stateVecLo,
171  ComplexArray stateVecOut, int updateUpper);
172 
173 void statevec_controlledNotLocal(Qureg qureg, int controlQubit, int targetQubit);
174 
175 void statevec_controlledNotDistributed (Qureg qureg, int controlQubit,
176  ComplexArray stateVecIn,
177  ComplexArray stateVecOut);
178 
179 void statevec_multiControlledMultiQubitNotLocal(Qureg qureg, int ctrlMask, int targMask);
180 
181 void statevec_multiControlledMultiQubitNotDistributed(Qureg qureg, int ctrlMask, int targMask,
182  ComplexArray stateVecIn,
183  ComplexArray stateVecOut);
184 
185 qreal statevec_findProbabilityOfZeroLocal (Qureg qureg, int measureQubit);
186 
188 
189 void statevec_collapseToKnownProbOutcomeLocal(Qureg qureg, int measureQubit, int outcome, qreal totalProbability);
190 
191 void statevec_collapseToKnownProbOutcomeDistributedRenorm (Qureg qureg, int measureQubit, qreal totalProbability);
192 
194 
195 void statevec_swapQubitAmpsLocal(Qureg qureg, int qb1, int qb2);
196 
197 void statevec_swapQubitAmpsDistributed(Qureg qureg, int pairRank, int qb1, int qb2);
198 
199 void statevec_multiControlledTwoQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u);
200 
201 void statevec_multiControlledMultiQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int* targs, int numTargs, ComplexMatrixN u);
202 
204 
205 void statevec_calcProbOfAllOutcomesLocal(qreal* retProbs, Qureg qureg, int* qubits, int numQubits);
206 
207 
208 # endif // QUEST_CPU_INTERNAL_H
qreal densmatr_calcPurityLocal(Qureg qureg)
Definition: QuEST_cpu.c:872
qreal densmatr_calcHilbertSchmidtDistanceSquaredLocal(Qureg a, Qureg b)
computes Tr((a-b) conjTrans(a-b)) = sum of abs values of (a-b)
Definition: QuEST_cpu.c:934
void statevec_collapseToKnownProbOutcomeLocal(Qureg qureg, int measureQubit, int outcome, qreal totalProbability)
Update the state vector to be consistent with measuring measureQubit=0 if outcome=0 and measureQubit=...
Definition: QuEST_cpu.c:3767
void densmatr_mixTwoQubitDepolarisingDistributed(Qureg qureg, int targetQubit, int qubit2, qreal delta, qreal gamma)
Definition: QuEST_cpu.c:547
void statevec_swapQubitAmpsDistributed(Qureg qureg, int pairRank, int qb1, int qb2)
qureg.pairStateVec contains the entire set of amplitudes of the paired node which includes the set of...
Definition: QuEST_cpu.c:3965
void statevec_controlledUnitaryLocal(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_cpu.c:2336
void densmatr_mixDepolarisingLocal(Qureg qureg, int targetQubit, qreal depolLevel)
Definition: QuEST_cpu.c:131
void statevec_controlledPauliYDistributed(Qureg qureg, int controlQubit, ComplexArray stateVecIn, ComplexArray stateVecOut, int conjFactor)
Definition: QuEST_cpu.c:3036
Represents a 4x4 matrix of complex numbers.
Definition: QuEST.h:175
void densmatr_initPureStateLocal(Qureg targetQureg, Qureg copyQureg)
Definition: QuEST_cpu.c:1195
void densmatr_mixDampingDistributed(Qureg qureg, int targetQubit, qreal damping)
Definition: QuEST_cpu.c:306
Represents a general 2^N by 2^N matrix of complex numbers.
Definition: QuEST.h:186
#define qreal
void statevec_multiControlledMultiQubitNotLocal(Qureg qureg, int ctrlMask, int targMask)
Definition: QuEST_cpu.c:2778
void statevec_unitaryDistributed(Qureg qureg, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
Apply a unitary operation to a single qubit given a subset of the state vector with upper and lower b...
Definition: QuEST_cpu.c:2151
void densmatr_calcProbOfAllOutcomesLocal(qreal *retProbs, Qureg qureg, int *qubits, int numQubits)
Definition: QuEST_cpu.c:3616
static long long int flipBit(const long long int number, const int bitInd)
void statevec_controlledNotDistributed(Qureg qureg, int controlQubit, ComplexArray stateVecIn, ComplexArray stateVecOut)
Rotate a single qubit by {{0,1},{1,0}.
Definition: QuEST_cpu.c:2742
void statevec_controlledPauliYLocal(Qureg qureg, int controlQubit, int targetQubit, int conjFactor)
Definition: QuEST_cpu.c:2982
void statevec_multiControlledMultiQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int *targs, int numTargs, ComplexMatrixN u)
Definition: QuEST_cpu.c:1912
static int extractBit(const int locationOfBitFromRight, const long long int theEncodedNumber)
Complex statevec_calcInnerProductLocal(Qureg bra, Qureg ket)
Definition: QuEST_cpu.c:1082
void statevec_hadamardDistributed(Qureg qureg, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut, int updateUpper)
Rotate a single qubit by {{1,1},{1,-1}}/sqrt2.
Definition: QuEST_cpu.c:3139
qreal statevec_findProbabilityOfZeroDistributed(Qureg qureg)
Measure the probability of a specified qubit being in the zero state across all amplitudes held in th...
Definition: QuEST_cpu.c:3513
void densmatr_applyDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:4082
void statevec_unitaryLocal(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_cpu.c:2026
void statevec_calcProbOfAllOutcomesLocal(qreal *retProbs, Qureg qureg, int *qubits, int numQubits)
Definition: QuEST_cpu.c:3549
void statevec_collapseToOutcomeDistributedSetZero(Qureg qureg)
Definition: QuEST_cpu.c:3887
qreal densmatr_calcFidelityLocal(Qureg qureg, Qureg pureState)
computes a few dens-columns-worth of (vec^*T) dens * vec
Definition: QuEST_cpu.c:1001
void densmatr_mixDampingLocal(Qureg qureg, int targetQubit, qreal damping)
Definition: QuEST_cpu.c:180
void statevec_collapseToKnownProbOutcomeDistributedRenorm(Qureg qureg, int measureQubit, qreal totalProbability)
Renormalise parts of the state vector where measureQubit=0 or 1, based on the total probability of th...
Definition: QuEST_cpu.c:3849
void densmatr_mixDepolarisingDistributed(Qureg qureg, int targetQubit, qreal depolLevel)
Definition: QuEST_cpu.c:230
qreal densmatr_findProbabilityOfZeroLocal(Qureg qureg, int measureQubit)
Definition: QuEST_cpu.c:3402
void statevec_controlledUnitaryDistributed(Qureg qureg, int controlQubit, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha ...
Definition: QuEST_cpu.c:2476
void statevec_hadamardLocal(Qureg qureg, int targetQubit)
Definition: QuEST_cpu.c:3078
Complex densmatr_calcExpecDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:4167
void statevec_multiControlledUnitaryDistributed(Qureg qureg, int targetQubit, long long int ctrlQubitsMask, long long int ctrlFlipMask, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
Apply a unitary operation to a single qubit in the state vector of probability amplitudes,...
Definition: QuEST_cpu.c:2542
void statevec_pauliYDistributed(Qureg qureg, ComplexArray stateVecIn, ComplexArray stateVecOut, int updateUpper, int conjFac)
Rotate a single qubit by +-{{0,-i},{i,0}.
Definition: QuEST_cpu.c:2945
Represents a diagonal complex operator on the full Hilbert state of a Qureg.
Definition: QuEST.h:297
void densmatr_mixTwoQubitDepolarisingLocal(Qureg qureg, int qubit1, int qubit2, qreal delta, qreal gamma)
Definition: QuEST_cpu.c:393
Complex statevec_calcExpecDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:4124
void statevec_multiControlledUnitaryLocal(Qureg qureg, int targetQubit, long long int ctrlQubitsMask, long long int ctrlFlipMask, ComplexMatrix2 u)
Definition: QuEST_cpu.c:2268
void densmatr_mixTwoQubitDepolarisingQ1LocalQ2DistributedPart3(Qureg qureg, int targetQubit, int qubit2, qreal delta, qreal gamma)
Definition: QuEST_cpu.c:638
void statevec_pauliXLocal(Qureg qureg, int targetQubit)
Definition: QuEST_cpu.c:2593
qreal densmatr_calcInnerProductLocal(Qureg a, Qureg b)
computes Tr(conjTrans(a) b) = sum of (a_ij^* b_ij)
Definition: QuEST_cpu.c:969
Represents a system of qubits.
Definition: QuEST.h:322
void statevec_multiControlledMultiQubitNotDistributed(Qureg qureg, int ctrlMask, int targMask, ComplexArray stateVecIn, ComplexArray stateVecOut)
Definition: QuEST_cpu.c:2834
static long long int insertZeroBit(const long long int number, const int index)
void statevec_multiControlledTwoQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u)
Definition: QuEST_cpu.c:1813
void statevec_pauliYLocal(Qureg qureg, int targetQubit, int conjFac)
Definition: QuEST_cpu.c:2887
static int isOddParity(const long long int number, const int qb1, const int qb2)
void statevec_swapQubitAmpsLocal(Qureg qureg, int qb1, int qb2)
It is ensured that all amplitudes needing to be swapped are on this node.
Definition: QuEST_cpu.c:3922
void statevec_compactUnitaryLocal(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_cpu.c:1754
void statevec_controlledCompactUnitaryDistributed(Qureg qureg, int controlQubit, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha ...
Definition: QuEST_cpu.c:2414
void densmatr_mixTwoQubitDepolarisingLocalPart1(Qureg qureg, int qubit1, int qubit2, qreal delta)
Definition: QuEST_cpu.c:494
Represents one complex number.
Definition: QuEST.h:103
void statevec_controlledCompactUnitaryLocal(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_cpu.c:2196
static int maskContainsBit(const long long int mask, const int bitInd)
void statevec_pauliXDistributed(Qureg qureg, ComplexArray stateVecIn, ComplexArray stateVecOut)
Rotate a single qubit by {{0,1},{1,0}.
Definition: QuEST_cpu.c:2651
void statevec_compactUnitaryDistributed(Qureg qureg, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha ...
Definition: QuEST_cpu.c:2095
static long long int insertTwoZeroBits(const long long int number, const int bit1, const int bit2)
qreal statevec_findProbabilityOfZeroLocal(Qureg qureg, int measureQubit)
Measure the total probability of a specified qubit being in the zero state across all amplitudes in t...
Definition: QuEST_cpu.c:3457
void statevec_controlledNotLocal(Qureg qureg, int controlQubit, int targetQubit)
Definition: QuEST_cpu.c:2679
Represents a 2x2 matrix of complex numbers.
Definition: QuEST.h:137