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 
95 /*
96  * state vector operations
97  */
98 
100 
101 void statevec_compactUnitaryLocal (Qureg qureg, int targetQubit, Complex alpha, Complex beta);
102 
104  Complex rot1, Complex rot2,
105  ComplexArray stateVecUp,
106  ComplexArray stateVecLo,
107  ComplexArray stateVecOut);
108 
109 void statevec_unitaryLocal(Qureg qureg, int targetQubit, ComplexMatrix2 u);
110 
112  Complex rot1, Complex rot2,
113  ComplexArray stateVecUp,
114  ComplexArray stateVecLo,
115  ComplexArray stateVecOut);
116 
117 void statevec_controlledCompactUnitaryLocal (Qureg qureg, int controlQubit, int targetQubit,
118  Complex alpha, Complex beta);
119 
120 void statevec_controlledCompactUnitaryDistributed (Qureg qureg, int controlQubit,
121  Complex rot1, Complex rot2,
122  ComplexArray stateVecUp,
123  ComplexArray stateVecLo,
124  ComplexArray stateVecOut);
125 
126 void statevec_controlledUnitaryLocal(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u);
127 
128 void statevec_controlledUnitaryDistributed (Qureg qureg, int controlQubit,
129  Complex rot1, Complex rot2,
130  ComplexArray stateVecUp,
131  ComplexArray stateVecLo,
132  ComplexArray stateVecOut);
133 
134 void statevec_multiControlledUnitaryLocal(Qureg qureg, int targetQubit,
135  long long int ctrlQubitsMask, long long int ctrlFlipMask, ComplexMatrix2 u);
136 
138  int targetQubit,
139  long long int ctrlQubitsMask, long long int ctrlFlipMask,
140  Complex rot1, Complex rot2,
141  ComplexArray stateVecUp,
142  ComplexArray stateVecLo,
143  ComplexArray stateVecOut);
144 
145 void statevec_pauliXLocal(Qureg qureg, int targetQubit);
146 
148  ComplexArray stateVecIn,
149  ComplexArray stateVecOut);
150 
151 void statevec_pauliYLocal(Qureg qureg, int targetQubit, int conjFac);
152 
154  ComplexArray stateVecIn,
155  ComplexArray stateVecOut,
156  int updateUpper, int conjFac);
157 
158 void statevec_controlledPauliYLocal(Qureg qureg, int controlQubit, int targetQubit, int conjFactor);
159 
160 void statevec_controlledPauliYDistributed(Qureg qureg, int controlQubit,
161  ComplexArray stateVecIn,
162  ComplexArray stateVecOut, int conjFactor);
163 
164 void statevec_hadamardLocal (Qureg qureg, int targetQubit);
165 
167  ComplexArray stateVecUp,
168  ComplexArray stateVecLo,
169  ComplexArray stateVecOut, int updateUpper);
170 
171 void statevec_controlledNotLocal(Qureg qureg, int controlQubit, int targetQubit);
172 
173 void statevec_controlledNotDistributed (Qureg qureg, int controlQubit,
174  ComplexArray stateVecIn,
175  ComplexArray stateVecOut);
176 
177 qreal statevec_findProbabilityOfZeroLocal (Qureg qureg, int measureQubit);
178 
180 
181 void statevec_collapseToKnownProbOutcomeLocal(Qureg qureg, int measureQubit, int outcome, qreal totalProbability);
182 
183 void statevec_collapseToKnownProbOutcomeDistributedRenorm (Qureg qureg, int measureQubit, qreal totalProbability);
184 
186 
187 void statevec_swapQubitAmpsLocal(Qureg qureg, int qb1, int qb2);
188 
189 void statevec_swapQubitAmpsDistributed(Qureg qureg, int pairRank, int qb1, int qb2);
190 
191 void statevec_multiControlledTwoQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u);
192 
193 void statevec_multiControlledMultiQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int* targs, int numTargs, ComplexMatrixN u);
194 
196 
197 
198 # endif // QUEST_CPU_INTERNAL_H
qreal densmatr_calcPurityLocal(Qureg qureg)
Definition: QuEST_cpu.c:861
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:923
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:3380
void densmatr_mixTwoQubitDepolarisingDistributed(Qureg qureg, int targetQubit, int qubit2, qreal delta, qreal gamma)
Definition: QuEST_cpu.c:541
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:3579
void statevec_controlledUnitaryLocal(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_cpu.c:2241
void densmatr_mixDepolarisingLocal(Qureg qureg, int targetQubit, qreal depolLevel)
Definition: QuEST_cpu.c:125
void statevec_controlledPauliYDistributed(Qureg qureg, int controlQubit, ComplexArray stateVecIn, ComplexArray stateVecOut, int conjFactor)
Definition: QuEST_cpu.c:2830
Represents a 4x4 matrix of complex numbers.
Definition: QuEST.h:125
void densmatr_initPureStateLocal(Qureg targetQureg, Qureg copyQureg)
Definition: QuEST_cpu.c:1184
void densmatr_mixDampingDistributed(Qureg qureg, int targetQubit, qreal damping)
Definition: QuEST_cpu.c:300
Represents a general 2^N by 2^N matrix of complex numbers.
Definition: QuEST.h:136
#define qreal
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:2056
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:2646
void statevec_controlledPauliYLocal(Qureg qureg, int controlQubit, int targetQubit, int conjFactor)
Definition: QuEST_cpu.c:2776
void statevec_multiControlledMultiQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int *targs, int numTargs, ComplexMatrixN u)
Definition: QuEST_cpu.c:1846
static int extractBit(const int locationOfBitFromRight, const long long int theEncodedNumber)
Complex statevec_calcInnerProductLocal(Qureg bra, Qureg ket)
Definition: QuEST_cpu.c:1071
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:2932
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:3262
void densmatr_applyDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:3696
void statevec_unitaryLocal(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_cpu.c:1932
void statevec_collapseToOutcomeDistributedSetZero(Qureg qureg)
Set all amplitudes in one chunk to 0.
Definition: QuEST_cpu.c:3501
qreal densmatr_calcFidelityLocal(Qureg qureg, Qureg pureState)
computes a few dens-columns-worth of (vec^*T) dens * vec
Definition: QuEST_cpu.c:990
void densmatr_mixDampingLocal(Qureg qureg, int targetQubit, qreal damping)
Definition: QuEST_cpu.c:174
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:3462
void densmatr_mixDepolarisingDistributed(Qureg qureg, int targetQubit, qreal depolLevel)
Definition: QuEST_cpu.c:224
qreal densmatr_findProbabilityOfZeroLocal(Qureg qureg, int measureQubit)
Definition: QuEST_cpu.c:3151
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:2381
void statevec_hadamardLocal(Qureg qureg, int targetQubit)
Definition: QuEST_cpu.c:2872
Complex densmatr_calcExpecDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:3781
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:2447
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:2739
Represents a diagonal complex operator on the full Hilbert state of a Qureg.
Definition: QuEST.h:178
void densmatr_mixTwoQubitDepolarisingLocal(Qureg qureg, int qubit1, int qubit2, qreal delta, qreal gamma)
Definition: QuEST_cpu.c:387
Complex statevec_calcExpecDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:3738
void statevec_multiControlledUnitaryLocal(Qureg qureg, int targetQubit, long long int ctrlQubitsMask, long long int ctrlFlipMask, ComplexMatrix2 u)
Definition: QuEST_cpu.c:2173
void densmatr_mixTwoQubitDepolarisingQ1LocalQ2DistributedPart3(Qureg qureg, int targetQubit, int qubit2, qreal delta, qreal gamma)
Definition: QuEST_cpu.c:632
void statevec_pauliXLocal(Qureg qureg, int targetQubit)
Definition: QuEST_cpu.c:2498
qreal densmatr_calcInnerProductLocal(Qureg a, Qureg b)
computes Tr(conjTrans(a) b) = sum of (a_ij^* b_ij)
Definition: QuEST_cpu.c:958
Represents a system of qubits.
Definition: QuEST.h:203
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:1747
void statevec_pauliYLocal(Qureg qureg, int targetQubit, int conjFac)
Definition: QuEST_cpu.c:2682
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:3536
void statevec_compactUnitaryLocal(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_cpu.c:1688
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:2319
void densmatr_mixTwoQubitDepolarisingLocalPart1(Qureg qureg, int qubit1, int qubit2, qreal delta)
Definition: QuEST_cpu.c:488
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:2101
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:2556
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:2001
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:3206
void statevec_controlledNotLocal(Qureg qureg, int controlQubit, int targetQubit)
Definition: QuEST_cpu.c:2584
Represents a 2x2 matrix of complex numbers.
Definition: QuEST.h:114