QuEST_cpu_local.c
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 
12 # include "QuEST.h"
13 # include "QuEST_internal.h"
14 # include "QuEST_precision.h"
15 # include "mt19937ar.h"
16 
17 # include "QuEST_cpu_internal.h"
18 
19 # include <stdlib.h>
20 # include <stdio.h>
21 # include <math.h>
22 # include <time.h>
23 # include <sys/types.h>
24 
25 # ifdef _OPENMP
26 # include <omp.h>
27 # endif
28 
29 
30 void densmatr_mixDepolarising(Qureg qureg, int targetQubit, qreal depolLevel) {
31  if (depolLevel == 0)
32  return;
33 
34  densmatr_mixDepolarisingLocal(qureg, targetQubit, depolLevel);
35 }
36 
37 void densmatr_mixDamping(Qureg qureg, int targetQubit, qreal damping) {
38  if (damping == 0)
39  return;
40 
41  densmatr_mixDampingLocal(qureg, targetQubit, damping);
42 }
43 
44 void densmatr_mixTwoQubitDepolarising(Qureg qureg, int qubit1, int qubit2, qreal depolLevel){
45  if (depolLevel == 0)
46  return;
47  qreal eta = 2/depolLevel;
48  qreal delta = eta - 1 - sqrt( (eta-1)*(eta-1) - 1 );
49  qreal gamma = 1+delta;
50  // TODO -- test delta too small
51 
52  gamma = 1/(gamma*gamma*gamma);
53  densmatr_mixTwoQubitDephasing(qureg, qubit1, qubit2, depolLevel);
54  densmatr_mixTwoQubitDepolarisingLocal(qureg, qubit1, qubit2, delta, gamma);
55 }
56 
58 
59  return densmatr_calcPurityLocal(qureg);
60 }
61 
63 
65  qreal dist = sqrt(distSquared);
66  return dist;
67 }
68 
70 
72  return scalar;
73 }
74 
76 
77  // save pointers to qureg's pair state
78  qreal* quregPairRePtr = qureg.pairStateVec.real;
79  qreal* quregPairImPtr = qureg.pairStateVec.imag;
80 
81  // populate qureg pair state with pure state (by repointing)
82  qureg.pairStateVec.real = pureState.stateVec.real;
83  qureg.pairStateVec.imag = pureState.stateVec.imag;
84 
85  // calculate fidelity using pairState
86  qreal fid = densmatr_calcFidelityLocal(qureg, pureState);
87 
88  // restore pointers
89  qureg.pairStateVec.real = quregPairRePtr;
90  qureg.pairStateVec.imag = quregPairImPtr;
91 
92  return fid;
93 }
94 
95 void densmatr_initPureState(Qureg qureg, Qureg pureState) {
96 
97  // save pointers to qureg's pair state
98  qreal* quregPairRePtr = qureg.pairStateVec.real;
99  qreal* quregPairImPtr = qureg.pairStateVec.imag;
100 
101  // populate qureg pair state with pure state (by repointing)
102  qureg.pairStateVec.real = pureState.stateVec.real;
103  qureg.pairStateVec.imag = pureState.stateVec.imag;
104 
105  // populate density matrix via it's pairState
106  densmatr_initPureStateLocal(qureg, pureState);
107 
108  // restore pointers
109  qureg.pairStateVec.real = quregPairRePtr;
110  qureg.pairStateVec.imag = quregPairImPtr;
111 }
112 
114 
115  return statevec_calcInnerProductLocal(bra, ket);
116 }
117 
119 
120  // computes the trace using Kahan summation
121  qreal pTotal=0;
122  qreal y, t, c;
123  c = 0;
124 
125  long long int numCols = 1LL << qureg.numQubitsRepresented;
126  long long diagIndex;
127 
128  for (int col=0; col< numCols; col++) {
129  diagIndex = col*(numCols + 1);
130  y = qureg.stateVec.real[diagIndex] - c;
131  t = pTotal + y;
132  c = ( t - pTotal ) - y; // brackets are important
133  pTotal = t;
134  }
135 
136  // does not check imaginary component, by design
137 
138  return pTotal;
139 }
140 
142  // implemented using Kahan summation for greater accuracy at a slight floating
143  // point operation overhead. For more details see https://en.wikipedia.org/wiki/Kahan_summation_algorithm
144  qreal pTotal=0;
145  qreal y, t, c;
146  long long int index;
147  long long int numAmpsPerRank = qureg.numAmpsPerChunk;
148  c = 0.0;
149  for (index=0; index<numAmpsPerRank; index++){
150  // Perform pTotal+=qureg.stateVec.real[index]*qureg.stateVec.real[index]; by Kahan
151 
152  y = qureg.stateVec.real[index]*qureg.stateVec.real[index] - c;
153  t = pTotal + y;
154  // Don't change the bracketing on the following line
155  c = ( t - pTotal ) - y;
156  pTotal = t;
157 
158  // Perform pTotal+=qureg.stateVec.imag[index]*qureg.stateVec.imag[index]; by Kahan
159 
160  y = qureg.stateVec.imag[index]*qureg.stateVec.imag[index] - c;
161  t = pTotal + y;
162  // Don't change the bracketing on the following line
163  c = ( t - pTotal ) - y;
164  pTotal = t;
165  }
166  return pTotal;
167 }
168 
169 
171  // init MPI environment
172 
173  QuESTEnv env;
174  env.rank=0;
175  env.numRanks=1;
176 
178 
179  return env;
180 }
181 
183  // MPI Barrier goes here in MPI version.
184 }
185 
186 int syncQuESTSuccess(int successCode){
187  return successCode;
188 }
189 
191  // MPI finalize goes here in MPI version. Call this function anyway for consistency
192 }
193 
195  printf("EXECUTION ENVIRONMENT:\n");
196  printf("Running locally on one node\n");
197  printf("Number of ranks is %d\n", env.numRanks);
198 # ifdef _OPENMP
199  printf("OpenMP enabled\n");
200  printf("Number of threads available is %d\n", omp_get_max_threads());
201 # else
202  printf("OpenMP disabled\n");
203 # endif
204  printf("Precision: size of qreal is %ld bytes\n", sizeof(qreal));
205 }
206 
207 qreal statevec_getRealAmp(Qureg qureg, long long int index){
208  return qureg.stateVec.real[index];
209 }
210 
211 qreal statevec_getImagAmp(Qureg qureg, long long int index){
212  return qureg.stateVec.imag[index];
213 }
214 
215 void statevec_compactUnitary(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
216 {
217  statevec_compactUnitaryLocal(qureg, targetQubit, alpha, beta);
218 }
219 
220 void statevec_unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u)
221 {
222  statevec_unitaryLocal(qureg, targetQubit, u);
223 }
224 
225 void statevec_controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
226 {
227  statevec_controlledCompactUnitaryLocal(qureg, controlQubit, targetQubit, alpha, beta);
228 }
229 
230 void statevec_controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
231 {
232  statevec_controlledUnitaryLocal(qureg, controlQubit, targetQubit, u);
233 }
234 
235 void statevec_multiControlledUnitary(Qureg qureg, long long int ctrlQubitsMask, long long int ctrlFlipMask, int targetQubit, ComplexMatrix2 u)
236 {
237  statevec_multiControlledUnitaryLocal(qureg, targetQubit, ctrlQubitsMask, ctrlFlipMask, u);
238 }
239 
240 void statevec_pauliX(Qureg qureg, int targetQubit)
241 {
242  statevec_pauliXLocal(qureg, targetQubit);
243 }
244 
245 void statevec_pauliY(Qureg qureg, int targetQubit)
246 {
247  int conjFac = 1;
248  statevec_pauliYLocal(qureg, targetQubit, conjFac);
249 }
250 
251 void statevec_pauliYConj(Qureg qureg, int targetQubit)
252 {
253  int conjFac = -1;
254  statevec_pauliYLocal(qureg, targetQubit, conjFac);
255 }
256 
257 void statevec_controlledPauliY(Qureg qureg, int controlQubit, int targetQubit)
258 {
259  int conjFac = 1;
260  statevec_controlledPauliYLocal(qureg, controlQubit, targetQubit, conjFac);
261 }
262 
263 void statevec_controlledPauliYConj(Qureg qureg, int controlQubit, int targetQubit)
264 {
265  int conjFac = -1;
266  statevec_controlledPauliYLocal(qureg, controlQubit, targetQubit, conjFac);
267 }
268 
269 void statevec_hadamard(Qureg qureg, int targetQubit)
270 {
271  statevec_hadamardLocal(qureg, targetQubit);
272 }
273 
274 void statevec_controlledNot(Qureg qureg, int controlQubit, int targetQubit)
275 {
276  statevec_controlledNotLocal(qureg, controlQubit, targetQubit);
277 }
278 
279 qreal statevec_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
280 {
281  qreal stateProb=0;
282  stateProb = statevec_findProbabilityOfZeroLocal(qureg, measureQubit);
283  if (outcome==1) stateProb = 1.0 - stateProb;
284  return stateProb;
285 }
286 
287 qreal densmatr_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome) {
288 
289  qreal outcomeProb = densmatr_findProbabilityOfZeroLocal(qureg, measureQubit);
290  if (outcome == 1)
291  outcomeProb = 1.0 - outcomeProb;
292  return outcomeProb;
293 }
294 
295 void statevec_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal stateProb)
296 {
297  statevec_collapseToKnownProbOutcomeLocal(qureg, measureQubit, outcome, stateProb);
298 }
299 
300 void seedQuESTDefault(void){
301  // init MT random number generator with three keys -- time and pid
302  // for the MPI version, it is ok that all procs will get the same seed as random numbers will only be
303  // used by the master process
304 
305  unsigned long int key[2];
307  init_by_array(key, 2);
308 }
309 
310 void statevec_multiControlledTwoQubitUnitary(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u)
311 {
312  statevec_multiControlledTwoQubitUnitaryLocal(qureg, ctrlMask, q1, q2, u);
313 }
314 
315 void statevec_multiControlledMultiQubitUnitary(Qureg qureg, long long int ctrlMask, int* targs, int numTargs, ComplexMatrixN u)
316 {
317  statevec_multiControlledMultiQubitUnitaryLocal(qureg, ctrlMask, targs, numTargs, u);
318 }
319 
320 void statevec_swapQubitAmps(Qureg qureg, int qb1, int qb2)
321 {
322  statevec_swapQubitAmpsLocal(qureg, qb1, qb2);
323 }
324 
326 
327  // we must preload qureg.pairStateVec with the elements of op.
328  // instead of needless cloning, we'll just temporarily swap the pointers
329  qreal* rePtr = qureg.pairStateVec.real;
330  qreal* imPtr = qureg.pairStateVec.imag;
331  qureg.pairStateVec.real = op.real;
332  qureg.pairStateVec.imag = op.imag;
333 
335 
336  qureg.pairStateVec.real = rePtr;
337  qureg.pairStateVec.imag = imPtr;
338 }
339 
341 
342  return statevec_calcExpecDiagonalOpLocal(qureg, op);
343 }
344 
346 
347  return densmatr_calcExpecDiagonalOpLocal(qureg, op);
348 }
qreal densmatr_calcTotalProb(Qureg qureg)
void statevec_pauliYConj(Qureg qureg, int targetQubit)
void destroyQuESTEnv(QuESTEnv env)
Destroy the QuEST environment.
qreal statevec_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
void init_by_array(unsigned long init_key[], int key_length)
Definition: mt19937ar.c:80
void syncQuESTEnv(QuESTEnv env)
Guarantees that all code up to the given point has been executed on all nodes (if running in distribu...
void statevec_controlledNotLocal(Qureg qureg, int controlQubit, int targetQubit)
Definition: QuEST_cpu.c:2584
void statevec_pauliYLocal(Qureg qureg, int targetQubit, int conjFac)
Definition: QuEST_cpu.c:2682
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_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal stateProb)
int rank
Definition: QuEST.h:244
void statevec_compactUnitary(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
ComplexArray pairStateVec
Temporary storage for a chunk of the state vector received from another process in the MPI version.
Definition: QuEST.h:224
qreal densmatr_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
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_multiControlledUnitary(Qureg qureg, long long int ctrlQubitsMask, long long int ctrlFlipMask, int targetQubit, ComplexMatrix2 u)
void getQuESTDefaultSeedKey(unsigned long int *key)
Definition: QuEST_common.c:182
void statevec_controlledPauliY(Qureg qureg, int controlQubit, int targetQubit)
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 statevec_hadamard(Qureg qureg, int targetQubit)
void densmatr_mixDepolarising(Qureg qureg, int targetQubit, qreal depolLevel)
void statevec_multiControlledTwoQubitUnitary(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u)
This calls swapQubitAmps only when it would involve a distributed communication; if the qubit chunks ...
void statevec_compactUnitaryLocal(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_cpu.c:1688
qreal densmatr_calcPurityLocal(Qureg qureg)
Definition: QuEST_cpu.c:861
Complex statevec_calcExpecDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:3738
void statevec_multiControlledMultiQubitUnitary(Qureg qureg, long long int ctrlMask, int *targs, int numTargs, ComplexMatrixN u)
This calls swapQubitAmps only when it would involve a distributed communication; if the qubit chunks ...
Represents a 4x4 matrix of complex numbers.
Definition: QuEST.h:125
Information about the environment the program is running in.
Definition: QuEST.h:242
void statevec_multiControlledMultiQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int *targs, int numTargs, ComplexMatrixN u)
Definition: QuEST_cpu.c:1846
void statevec_multiControlledTwoQubitUnitaryLocal(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u)
Definition: QuEST_cpu.c:1747
Represents a general 2^N by 2^N matrix of complex numbers.
Definition: QuEST.h:136
void statevec_controlledNot(Qureg qureg, int controlQubit, int targetQubit)
#define qreal
void statevec_controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
void statevec_controlledPauliYConj(Qureg qureg, int controlQubit, int targetQubit)
qreal * imag
The imaginary values of the 2^numQubits complex elements.
Definition: QuEST.h:191
qreal densmatr_calcFidelity(Qureg qureg, Qureg pureState)
void densmatr_initPureState(Qureg qureg, Qureg pureState)
long long int numAmpsPerChunk
Number of probability amplitudes held in stateVec by this process In the non-MPI version,...
Definition: QuEST.h:213
void densmatr_mixTwoQubitDepolarisingLocal(Qureg qureg, int qubit1, int qubit2, qreal delta, qreal gamma)
Definition: QuEST_cpu.c:387
void statevec_controlledPauliYLocal(Qureg qureg, int controlQubit, int targetQubit, int conjFac)
Definition: QuEST_cpu.c:2776
Complex statevec_calcInnerProduct(Qureg bra, Qureg ket)
void statevec_controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
qreal densmatr_calcPurity(Qureg qureg)
int numRanks
Definition: QuEST.h:245
void densmatr_mixTwoQubitDephasing(Qureg qureg, int qubit1, int qubit2, qreal dephase)
Definition: QuEST_cpu.c:84
Complex statevec_calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
void statevec_multiControlledUnitaryLocal(Qureg qureg, int targetQubit, long long int ctrlQubitsMask, long long int ctrlFlipMask, ComplexMatrix2 u)
Definition: QuEST_cpu.c:2173
Represents a diagonal complex operator on the full Hilbert state of a Qureg.
Definition: QuEST.h:178
Complex densmatr_calcExpecDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:3781
void statevec_pauliXLocal(Qureg qureg, int targetQubit)
Definition: QuEST_cpu.c:2498
void densmatr_initPureStateLocal(Qureg targetQureg, Qureg copyQureg)
Definition: QuEST_cpu.c:1184
void statevec_pauliY(Qureg qureg, int targetQubit)
void statevec_pauliX(Qureg qureg, int targetQubit)
void densmatr_mixDampingLocal(Qureg qureg, int targetQubit, qreal damping)
Definition: QuEST_cpu.c:174
Represents a system of qubits.
Definition: QuEST.h:203
qreal densmatr_calcInnerProduct(Qureg a, Qureg b)
qreal densmatr_calcInnerProductLocal(Qureg a, Qureg b)
computes Tr(conjTrans(a) b) = sum of (a_ij^* b_ij)
Definition: QuEST_cpu.c:958
void densmatr_applyDiagonalOp(Qureg qureg, DiagonalOp op)
qreal statevec_calcTotalProb(Qureg qureg)
void densmatr_mixTwoQubitDepolarising(Qureg qureg, int qubit1, int qubit2, qreal depolLevel)
ComplexArray stateVec
Computational state amplitudes - a subset thereof in the MPI version.
Definition: QuEST.h:222
qreal statevec_getImagAmp(Qureg qureg, long long int index)
void densmatr_mixDepolarisingLocal(Qureg qureg, int targetQubit, qreal depolLevel)
Definition: QuEST_cpu.c:125
void seedQuESTDefault(void)
Seed the Mersenne Twister used for random number generation in the QuEST environment with an example ...
void statevec_hadamardLocal(Qureg qureg, int targetQubit)
Definition: QuEST_cpu.c:2872
qreal densmatr_calcHilbertSchmidtDistance(Qureg a, Qureg b)
void reportQuESTEnv(QuESTEnv env)
Report information about the QuEST environment.
void densmatr_applyDiagonalOpLocal(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:3696
void statevec_controlledUnitaryLocal(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_cpu.c:2241
int numQubitsRepresented
The number of qubits represented in either the state-vector or density matrix.
Definition: QuEST.h:208
int syncQuESTSuccess(int successCode)
Performs a logical AND on all successCodes held by all processes.
qreal * real
The real values of the 2^numQubits complex elements.
Definition: QuEST.h:189
void densmatr_mixDamping(Qureg qureg, int targetQubit, qreal damping)
void statevec_unitaryLocal(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_cpu.c:1932
qreal densmatr_findProbabilityOfZeroLocal(Qureg qureg, int measureQubit)
Definition: QuEST_cpu.c:3151
qreal statevec_getRealAmp(Qureg qureg, long long int index)
Complex densmatr_calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
Represents one complex number.
Definition: QuEST.h:103
QuESTEnv createQuESTEnv(void)
Create the QuEST execution environment.
void statevec_controlledCompactUnitaryLocal(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_cpu.c:2101
Complex statevec_calcInnerProductLocal(Qureg bra, Qureg ket)
Definition: QuEST_cpu.c:1071
void statevec_unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Represents a 2x2 matrix of complex numbers.
Definition: QuEST.h:114
void statevec_swapQubitAmps(Qureg qureg, int qb1, int qb2)
qreal densmatr_calcFidelityLocal(Qureg qureg, Qureg pureState)
computes a few dens-columns-worth of (vec^*T) dens * vec
Definition: QuEST_cpu.c:990
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