QuEST.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 
18 # include "QuEST.h"
19 # include "QuEST_precision.h"
20 # include "QuEST_internal.h"
21 # include "QuEST_validation.h"
22 # include "QuEST_qasm.h"
23 
24 # include <stdlib.h>
25 # include <string.h>
26 
27 #ifdef __cplusplus
28 extern "C" {
29 #endif
30 
31 
32 /*
33  * state-vector management
34  */
35 
36 Qureg createQureg(int numQubits, QuESTEnv env) {
37  validateNumQubitsInQureg(numQubits, env.numRanks, __func__);
38 
39  Qureg qureg;
40  statevec_createQureg(&qureg, numQubits, env);
41  qureg.isDensityMatrix = 0;
42  qureg.numQubitsRepresented = numQubits;
43  qureg.numQubitsInStateVec = numQubits;
44 
45  qasm_setup(&qureg);
46  initZeroState(qureg); // safe call to public function
47  return qureg;
48 }
49 
50 Qureg createDensityQureg(int numQubits, QuESTEnv env) {
51  validateNumQubitsInQureg(2*numQubits, env.numRanks, __func__);
52 
53  Qureg qureg;
54  statevec_createQureg(&qureg, 2*numQubits, env);
55  qureg.isDensityMatrix = 1;
56  qureg.numQubitsRepresented = numQubits;
57  qureg.numQubitsInStateVec = 2*numQubits;
58 
59  qasm_setup(&qureg);
60  initZeroState(qureg); // safe call to public function
61  return qureg;
62 }
63 
65 
66  Qureg newQureg;
67  statevec_createQureg(&newQureg, qureg.numQubitsInStateVec, env);
68  newQureg.isDensityMatrix = qureg.isDensityMatrix;
71 
72  qasm_setup(&newQureg);
73  statevec_cloneQureg(newQureg, qureg);
74  return newQureg;
75 }
76 
77 void destroyQureg(Qureg qureg, QuESTEnv env) {
78  statevec_destroyQureg(qureg, env);
79  qasm_free(qureg);
80 }
81 
82 
83 /*
84  * QASM
85  */
86 
88  qasm_startRecording(qureg);
89 }
90 
91 void stopRecordingQASM(Qureg qureg) {
92  qasm_stopRecording(qureg);
93 }
94 
95 void clearRecordedQASM(Qureg qureg) {
96  qasm_clearRecorded(qureg);
97 }
98 
99 void printRecordedQASM(Qureg qureg) {
100  qasm_printRecorded(qureg);
101 }
102 
103 void writeRecordedQASMToFile(Qureg qureg, char* filename) {
104  int success = qasm_writeRecordedToFile(qureg, filename);
105  validateFileOpened(success, filename, __func__);
106 }
107 
108 
109 /*
110  * state initialisation
111  */
112 
113 void initZeroState(Qureg qureg) {
114  statevec_initZeroState(qureg); // valid for both statevec and density matrices
115 
116  qasm_recordInitZero(qureg);
117 }
118 
119 void initBlankState(Qureg qureg) {
121 
122  qasm_recordComment(qureg, "Here, the register was initialised to an unphysical all-zero-amplitudes 'state'.");
123 }
124 
125 void initPlusState(Qureg qureg) {
126  if (qureg.isDensityMatrix)
127  densmatr_initPlusState(qureg);
128  else
129  statevec_initPlusState(qureg);
130 
131  qasm_recordInitPlus(qureg);
132 }
133 
134 void initClassicalState(Qureg qureg, long long int stateInd) {
135  validateStateIndex(qureg, stateInd, __func__);
136 
137  if (qureg.isDensityMatrix)
138  densmatr_initClassicalState(qureg, stateInd);
139  else
140  statevec_initClassicalState(qureg, stateInd);
141 
142  qasm_recordInitClassical(qureg, stateInd);
143 }
144 
145 void initPureState(Qureg qureg, Qureg pure) {
146  validateSecondQuregStateVec(pure, __func__);
147  validateMatchingQuregDims(qureg, pure, __func__);
148 
149  if (qureg.isDensityMatrix)
150  densmatr_initPureState(qureg, pure);
151  else
152  statevec_cloneQureg(qureg, pure);
153 
154  qasm_recordComment(qureg, "Here, the register was initialised to an undisclosed given pure state.");
155 }
156 
157 void initStateFromAmps(Qureg qureg, qreal* reals, qreal* imags) {
158  validateStateVecQureg(qureg, __func__);
159 
160  statevec_setAmps(qureg, 0, reals, imags, qureg.numAmpsTotal);
161 
162  qasm_recordComment(qureg, "Here, the register was initialised to an undisclosed given pure state.");
163 }
164 
165 void cloneQureg(Qureg targetQureg, Qureg copyQureg) {
166  validateMatchingQuregTypes(targetQureg, copyQureg, __func__);
167  validateMatchingQuregDims(targetQureg, copyQureg, __func__);
168 
169  statevec_cloneQureg(targetQureg, copyQureg);
170 }
171 
172 
173 /*
174  * unitary gates
175  */
176 
177 void hadamard(Qureg qureg, int targetQubit) {
178  validateTarget(qureg, targetQubit, __func__);
179 
180  statevec_hadamard(qureg, targetQubit);
181  if (qureg.isDensityMatrix) {
182  statevec_hadamard(qureg, targetQubit+qureg.numQubitsRepresented);
183  }
184 
185  qasm_recordGate(qureg, GATE_HADAMARD, targetQubit);
186 }
187 
188 void rotateX(Qureg qureg, int targetQubit, qreal angle) {
189  validateTarget(qureg, targetQubit, __func__);
190 
191  statevec_rotateX(qureg, targetQubit, angle);
192  if (qureg.isDensityMatrix) {
193  statevec_rotateX(qureg, targetQubit+qureg.numQubitsRepresented, -angle);
194  }
195 
196  qasm_recordParamGate(qureg, GATE_ROTATE_X, targetQubit, angle);
197 }
198 
199 void rotateY(Qureg qureg, int targetQubit, qreal angle) {
200  validateTarget(qureg, targetQubit, __func__);
201 
202  statevec_rotateY(qureg, targetQubit, angle);
203  if (qureg.isDensityMatrix) {
204  statevec_rotateY(qureg, targetQubit+qureg.numQubitsRepresented, angle);
205  }
206 
207  qasm_recordParamGate(qureg, GATE_ROTATE_Y, targetQubit, angle);
208 }
209 
210 void rotateZ(Qureg qureg, int targetQubit, qreal angle) {
211  validateTarget(qureg, targetQubit, __func__);
212 
213  statevec_rotateZ(qureg, targetQubit, angle);
214  if (qureg.isDensityMatrix) {
215  statevec_rotateZ(qureg, targetQubit+qureg.numQubitsRepresented, -angle);
216  }
217 
218  qasm_recordParamGate(qureg, GATE_ROTATE_Z, targetQubit, angle);
219 }
220 
221 void controlledRotateX(Qureg qureg, int controlQubit, int targetQubit, qreal angle) {
222  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
223 
224  statevec_controlledRotateX(qureg, controlQubit, targetQubit, angle);
225  if (qureg.isDensityMatrix) {
226  int shift = qureg.numQubitsRepresented;
227  statevec_controlledRotateX(qureg, controlQubit+shift, targetQubit+shift, -angle);
228  }
229 
230  qasm_recordControlledParamGate(qureg, GATE_ROTATE_X, controlQubit, targetQubit, angle);
231 }
232 
233 void controlledRotateY(Qureg qureg, int controlQubit, int targetQubit, qreal angle) {
234  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
235 
236  statevec_controlledRotateY(qureg, controlQubit, targetQubit, angle);
237  if (qureg.isDensityMatrix) {
238  int shift = qureg.numQubitsRepresented;
239  statevec_controlledRotateY(qureg, controlQubit+shift, targetQubit+shift, angle); // rotateY is real
240  }
241 
242  qasm_recordControlledParamGate(qureg, GATE_ROTATE_Y, controlQubit, targetQubit, angle);
243 }
244 
245 void controlledRotateZ(Qureg qureg, int controlQubit, int targetQubit, qreal angle) {
246  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
247 
248  statevec_controlledRotateZ(qureg, controlQubit, targetQubit, angle);
249  if (qureg.isDensityMatrix) {
250  int shift = qureg.numQubitsRepresented;
251  statevec_controlledRotateZ(qureg, controlQubit+shift, targetQubit+shift, -angle);
252  }
253 
254  qasm_recordControlledParamGate(qureg, GATE_ROTATE_Z, controlQubit, targetQubit, angle);
255 }
256 
257 void twoQubitUnitary(Qureg qureg, int targetQubit1, int targetQubit2, ComplexMatrix4 u) {
258  validateMultiTargets(qureg, (int []) {targetQubit1, targetQubit2}, 2, __func__);
259  validateTwoQubitUnitaryMatrix(qureg, u, __func__);
260 
261  statevec_twoQubitUnitary(qureg, targetQubit1, targetQubit2, u);
262  if (qureg.isDensityMatrix) {
263  int shift = qureg.numQubitsRepresented;
264  statevec_twoQubitUnitary(qureg, targetQubit1+shift, targetQubit2+shift, getConjugateMatrix4(u));
265  }
266 
267  qasm_recordComment(qureg, "Here, an undisclosed 2-qubit unitary was applied.");
268 }
269 
270 void controlledTwoQubitUnitary(Qureg qureg, int controlQubit, int targetQubit1, int targetQubit2, ComplexMatrix4 u) {
271  validateMultiControlsMultiTargets(qureg, (int[]) {controlQubit}, 1, (int[]) {targetQubit1, targetQubit2}, 2, __func__);
272  validateTwoQubitUnitaryMatrix(qureg, u, __func__);
273 
274  statevec_controlledTwoQubitUnitary(qureg, controlQubit, targetQubit1, targetQubit2, u);
275  if (qureg.isDensityMatrix) {
276  int shift = qureg.numQubitsRepresented;
277  statevec_controlledTwoQubitUnitary(qureg, controlQubit+shift, targetQubit1+shift, targetQubit2+shift, getConjugateMatrix4(u));
278  }
279 
280  qasm_recordComment(qureg, "Here, an undisclosed controlled 2-qubit unitary was applied.");
281 }
282 
283 void multiControlledTwoQubitUnitary(Qureg qureg, int* controlQubits, int numControlQubits, int targetQubit1, int targetQubit2, ComplexMatrix4 u) {
284  validateMultiControlsMultiTargets(qureg, controlQubits, numControlQubits, (int[]) {targetQubit1, targetQubit2}, 2, __func__);
285  validateTwoQubitUnitaryMatrix(qureg, u, __func__);
286 
287  long long int ctrlQubitsMask = getQubitBitMask(controlQubits, numControlQubits);
288  statevec_multiControlledTwoQubitUnitary(qureg, ctrlQubitsMask, targetQubit1, targetQubit2, u);
289  if (qureg.isDensityMatrix) {
290  int shift = qureg.numQubitsRepresented;
291  statevec_multiControlledTwoQubitUnitary(qureg, ctrlQubitsMask<<shift, targetQubit1+shift, targetQubit2+shift, getConjugateMatrix4(u));
292  }
293 
294  qasm_recordComment(qureg, "Here, an undisclosed multi-controlled 2-qubit unitary was applied.");
295 }
296 
297 void multiQubitUnitary(Qureg qureg, int* targs, int numTargs, ComplexMatrixN u) {
298  validateMultiTargets(qureg, targs, numTargs, __func__);
299  validateMultiQubitUnitaryMatrix(qureg, u, numTargs, __func__);
300 
301  statevec_multiQubitUnitary(qureg, targs, numTargs, u);
302  if (qureg.isDensityMatrix) {
303  int shift = qureg.numQubitsRepresented;
304  shiftIndices(targs, numTargs, shift);
306  statevec_multiQubitUnitary(qureg, targs, numTargs, u);
307  shiftIndices(targs, numTargs, -shift);
309  }
310 
311  qasm_recordComment(qureg, "Here, an undisclosed multi-qubit unitary was applied.");
312 }
313 
314 void controlledMultiQubitUnitary(Qureg qureg, int ctrl, int* targs, int numTargs, ComplexMatrixN u) {
315  validateMultiControlsMultiTargets(qureg, (int[]) {ctrl}, 1, targs, numTargs, __func__);
316  validateMultiQubitUnitaryMatrix(qureg, u, numTargs, __func__);
317 
318  statevec_controlledMultiQubitUnitary(qureg, ctrl, targs, numTargs, u);
319  if (qureg.isDensityMatrix) {
320  int shift = qureg.numQubitsRepresented;
321  shiftIndices(targs, numTargs, shift);
323  statevec_controlledMultiQubitUnitary(qureg, ctrl+shift, targs, numTargs, u);
324  shiftIndices(targs, numTargs, -shift);
326  }
327 
328  qasm_recordComment(qureg, "Here, an undisclosed controlled multi-qubit unitary was applied.");
329 }
330 
331 void multiControlledMultiQubitUnitary(Qureg qureg, int* ctrls, int numCtrls, int* targs, int numTargs, ComplexMatrixN u) {
332  validateMultiControlsMultiTargets(qureg, ctrls, numCtrls, targs, numTargs, __func__);
333  validateMultiQubitUnitaryMatrix(qureg, u, numTargs, __func__);
334 
335  long long int ctrlMask = getQubitBitMask(ctrls, numCtrls);
336  statevec_multiControlledMultiQubitUnitary(qureg, ctrlMask, targs, numTargs, u);
337  if (qureg.isDensityMatrix) {
338  int shift = qureg.numQubitsRepresented;
339  shiftIndices(targs, numTargs, shift);
341  statevec_multiControlledMultiQubitUnitary(qureg, ctrlMask<<shift, targs, numTargs, u);
342  shiftIndices(targs, numTargs, -shift);
344  }
345 
346  qasm_recordComment(qureg, "Here, an undisclosed multi-controlled multi-qubit unitary was applied.");
347 }
348 
349 void unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u) {
350  validateTarget(qureg, targetQubit, __func__);
351  validateOneQubitUnitaryMatrix(u, __func__);
352 
353  statevec_unitary(qureg, targetQubit, u);
354  if (qureg.isDensityMatrix) {
355  statevec_unitary(qureg, targetQubit+qureg.numQubitsRepresented, getConjugateMatrix2(u));
356  }
357 
358  qasm_recordUnitary(qureg, u, targetQubit);
359 }
360 
361 void controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u) {
362  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
363  validateOneQubitUnitaryMatrix(u, __func__);
364 
365  statevec_controlledUnitary(qureg, controlQubit, targetQubit, u);
366  if (qureg.isDensityMatrix) {
367  int shift = qureg.numQubitsRepresented;
368  statevec_controlledUnitary(qureg, controlQubit+shift, targetQubit+shift, getConjugateMatrix2(u));
369  }
370 
371  qasm_recordControlledUnitary(qureg, u, controlQubit, targetQubit);
372 }
373 
374 void multiControlledUnitary(Qureg qureg, int* controlQubits, int numControlQubits, int targetQubit, ComplexMatrix2 u) {
375  validateMultiControlsTarget(qureg, controlQubits, numControlQubits, targetQubit, __func__);
376  validateOneQubitUnitaryMatrix(u, __func__);
377 
378  long long int ctrlQubitsMask = getQubitBitMask(controlQubits, numControlQubits);
379  long long int ctrlFlipMask = 0;
380  statevec_multiControlledUnitary(qureg, ctrlQubitsMask, ctrlFlipMask, targetQubit, u);
381  if (qureg.isDensityMatrix) {
382  int shift = qureg.numQubitsRepresented;
383  statevec_multiControlledUnitary(qureg, ctrlQubitsMask<<shift, ctrlFlipMask<<shift, targetQubit+shift, getConjugateMatrix2(u));
384  }
385 
386  qasm_recordMultiControlledUnitary(qureg, u, controlQubits, numControlQubits, targetQubit);
387 }
388 
389 void multiStateControlledUnitary(Qureg qureg, int* controlQubits, int* controlState, int numControlQubits, int targetQubit, ComplexMatrix2 u) {
390  validateMultiControlsTarget(qureg, controlQubits, numControlQubits, targetQubit, __func__);
391  validateOneQubitUnitaryMatrix(u, __func__);
392  validateControlState(controlState, numControlQubits, __func__);
393 
394  long long int ctrlQubitsMask = getQubitBitMask(controlQubits, numControlQubits);
395  long long int ctrlFlipMask = getControlFlipMask(controlQubits, controlState, numControlQubits);
396  statevec_multiControlledUnitary(qureg, ctrlQubitsMask, ctrlFlipMask, targetQubit, u);
397  if (qureg.isDensityMatrix) {
398  int shift = qureg.numQubitsRepresented;
399  statevec_multiControlledUnitary(qureg, ctrlQubitsMask<<shift, ctrlFlipMask<<shift, targetQubit+shift, getConjugateMatrix2(u));
400  }
401 
402  qasm_recordMultiStateControlledUnitary(qureg, u, controlQubits, controlState, numControlQubits, targetQubit);
403 }
404 
405 void compactUnitary(Qureg qureg, int targetQubit, Complex alpha, Complex beta) {
406  validateTarget(qureg, targetQubit, __func__);
407  validateUnitaryComplexPair(alpha, beta, __func__);
408 
409  statevec_compactUnitary(qureg, targetQubit, alpha, beta);
410  if (qureg.isDensityMatrix) {
411  int shift = qureg.numQubitsRepresented;
412  statevec_compactUnitary(qureg, targetQubit+shift, getConjugateScalar(alpha), getConjugateScalar(beta));
413  }
414 
415  qasm_recordCompactUnitary(qureg, alpha, beta, targetQubit);
416 }
417 
418 void controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta) {
419  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
420  validateUnitaryComplexPair(alpha, beta, __func__);
421 
422  statevec_controlledCompactUnitary(qureg, controlQubit, targetQubit, alpha, beta);
423  if (qureg.isDensityMatrix) {
424  int shift = qureg.numQubitsRepresented;
426  controlQubit+shift, targetQubit+shift,
427  getConjugateScalar(alpha), getConjugateScalar(beta));
428  }
429 
430  qasm_recordControlledCompactUnitary(qureg, alpha, beta, controlQubit, targetQubit);
431 }
432 
433 void pauliX(Qureg qureg, int targetQubit) {
434  validateTarget(qureg, targetQubit, __func__);
435 
436  statevec_pauliX(qureg, targetQubit);
437  if (qureg.isDensityMatrix) {
438  statevec_pauliX(qureg, targetQubit+qureg.numQubitsRepresented);
439  }
440 
441  qasm_recordGate(qureg, GATE_SIGMA_X, targetQubit);
442 }
443 
444 void pauliY(Qureg qureg, int targetQubit) {
445  validateTarget(qureg, targetQubit, __func__);
446 
447  statevec_pauliY(qureg, targetQubit);
448  if (qureg.isDensityMatrix) {
449  statevec_pauliYConj(qureg, targetQubit + qureg.numQubitsRepresented);
450  }
451 
452  qasm_recordGate(qureg, GATE_SIGMA_Y, targetQubit);
453 }
454 
455 void pauliZ(Qureg qureg, int targetQubit) {
456  validateTarget(qureg, targetQubit, __func__);
457 
458  statevec_pauliZ(qureg, targetQubit);
459  if (qureg.isDensityMatrix) {
460  statevec_pauliZ(qureg, targetQubit+qureg.numQubitsRepresented);
461  }
462 
463  qasm_recordGate(qureg, GATE_SIGMA_Z, targetQubit);
464 }
465 
466 void sGate(Qureg qureg, int targetQubit) {
467  validateTarget(qureg, targetQubit, __func__);
468 
469  statevec_sGate(qureg, targetQubit);
470  if (qureg.isDensityMatrix) {
471  statevec_sGateConj(qureg, targetQubit+qureg.numQubitsRepresented);
472  }
473 
474  qasm_recordGate(qureg, GATE_S, targetQubit);
475 }
476 
477 void tGate(Qureg qureg, int targetQubit) {
478  validateTarget(qureg, targetQubit, __func__);
479 
480  statevec_tGate(qureg, targetQubit);
481  if (qureg.isDensityMatrix) {
482  statevec_tGateConj(qureg, targetQubit+qureg.numQubitsRepresented);
483  }
484 
485  qasm_recordGate(qureg, GATE_T, targetQubit);
486 }
487 
488 void phaseShift(Qureg qureg, int targetQubit, qreal angle) {
489  validateTarget(qureg, targetQubit, __func__);
490 
491  statevec_phaseShift(qureg, targetQubit, angle);
492  if (qureg.isDensityMatrix) {
493  statevec_phaseShift(qureg, targetQubit+qureg.numQubitsRepresented, -angle);
494  }
495 
496  qasm_recordParamGate(qureg, GATE_PHASE_SHIFT, targetQubit, angle);
497 }
498 
499 void controlledPhaseShift(Qureg qureg, int idQubit1, int idQubit2, qreal angle) {
500  validateControlTarget(qureg, idQubit1, idQubit2, __func__);
501 
502  statevec_controlledPhaseShift(qureg, idQubit1, idQubit2, angle);
503  if (qureg.isDensityMatrix) {
504  int shift = qureg.numQubitsRepresented;
505  statevec_controlledPhaseShift(qureg, idQubit1+shift, idQubit2+shift, -angle);
506  }
507 
508  qasm_recordControlledParamGate(qureg, GATE_PHASE_SHIFT, idQubit1, idQubit2, angle);
509 }
510 
511 void multiControlledPhaseShift(Qureg qureg, int *controlQubits, int numControlQubits, qreal angle) {
512  validateMultiQubits(qureg, controlQubits, numControlQubits, __func__);
513 
514  statevec_multiControlledPhaseShift(qureg, controlQubits, numControlQubits, angle);
515  if (qureg.isDensityMatrix) {
516  int shift = qureg.numQubitsRepresented;
517  shiftIndices(controlQubits, numControlQubits, shift);
518  statevec_multiControlledPhaseShift(qureg, controlQubits, numControlQubits, -angle);
519  shiftIndices(controlQubits, numControlQubits, -shift);
520  }
521 
522  qasm_recordMultiControlledParamGate(qureg, GATE_PHASE_SHIFT, controlQubits, numControlQubits-1, controlQubits[numControlQubits-1], angle);
523 }
524 
525 void controlledNot(Qureg qureg, int controlQubit, int targetQubit) {
526  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
527 
528  statevec_controlledNot(qureg, controlQubit, targetQubit);
529  if (qureg.isDensityMatrix) {
530  int shift = qureg.numQubitsRepresented;
531  statevec_controlledNot(qureg, controlQubit+shift, targetQubit+shift);
532  }
533 
534  qasm_recordControlledGate(qureg, GATE_SIGMA_X, controlQubit, targetQubit);
535 }
536 
537 void controlledPauliY(Qureg qureg, int controlQubit, int targetQubit) {
538  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
539 
540  statevec_controlledPauliY(qureg, controlQubit, targetQubit);
541  if (qureg.isDensityMatrix) {
542  int shift = qureg.numQubitsRepresented;
543  statevec_controlledPauliYConj(qureg, controlQubit+shift, targetQubit+shift);
544  }
545 
546  qasm_recordControlledGate(qureg, GATE_SIGMA_Y, controlQubit, targetQubit);
547 }
548 
549 void controlledPhaseFlip(Qureg qureg, int idQubit1, int idQubit2) {
550  validateControlTarget(qureg, idQubit1, idQubit2, __func__);
551 
552  statevec_controlledPhaseFlip(qureg, idQubit1, idQubit2);
553  if (qureg.isDensityMatrix) {
554  int shift = qureg.numQubitsRepresented;
555  statevec_controlledPhaseFlip(qureg, idQubit1+shift, idQubit2+shift);
556  }
557 
558  qasm_recordControlledGate(qureg, GATE_SIGMA_Z, idQubit1, idQubit2);
559 }
560 
561 void multiControlledPhaseFlip(Qureg qureg, int *controlQubits, int numControlQubits) {
562  validateMultiQubits(qureg, controlQubits, numControlQubits, __func__);
563 
564  statevec_multiControlledPhaseFlip(qureg, controlQubits, numControlQubits);
565  if (qureg.isDensityMatrix) {
566  int shift = qureg.numQubitsRepresented;
567  shiftIndices(controlQubits, numControlQubits, shift);
568  statevec_multiControlledPhaseFlip(qureg, controlQubits, numControlQubits);
569  shiftIndices(controlQubits, numControlQubits, -shift);
570  }
571 
572  qasm_recordMultiControlledGate(qureg, GATE_SIGMA_Z, controlQubits, numControlQubits-1, controlQubits[numControlQubits-1]);
573 }
574 
575 void rotateAroundAxis(Qureg qureg, int rotQubit, qreal angle, Vector axis) {
576  validateTarget(qureg, rotQubit, __func__);
577  validateVector(axis, __func__);
578 
579  statevec_rotateAroundAxis(qureg, rotQubit, angle, axis);
580  if (qureg.isDensityMatrix) {
581  int shift = qureg.numQubitsRepresented;
582  statevec_rotateAroundAxisConj(qureg, rotQubit+shift, angle, axis);
583  }
584 
585  qasm_recordAxisRotation(qureg, angle, axis, rotQubit);
586 }
587 
588 void controlledRotateAroundAxis(Qureg qureg, int controlQubit, int targetQubit, qreal angle, Vector axis) {
589  validateControlTarget(qureg, controlQubit, targetQubit, __func__);
590  validateVector(axis, __func__);
591 
592  statevec_controlledRotateAroundAxis(qureg, controlQubit, targetQubit, angle, axis);
593  if (qureg.isDensityMatrix) {
594  int shift = qureg.numQubitsRepresented;
595  statevec_controlledRotateAroundAxisConj(qureg, controlQubit+shift, targetQubit+shift, angle, axis);
596  }
597 
598  qasm_recordControlledAxisRotation(qureg, angle, axis, controlQubit, targetQubit);
599 }
600 
601 void swapGate(Qureg qureg, int qb1, int qb2) {
602  validateUniqueTargets(qureg, qb1, qb2, __func__);
603 
604  statevec_swapQubitAmps(qureg, qb1, qb2);
605  if (qureg.isDensityMatrix) {
606  int shift = qureg.numQubitsRepresented;
607  statevec_swapQubitAmps(qureg, qb1+shift, qb2+shift);
608  }
609 
610  qasm_recordControlledGate(qureg, GATE_SWAP, qb1, qb2);
611 }
612 
613 void sqrtSwapGate(Qureg qureg, int qb1, int qb2) {
614  validateUniqueTargets(qureg, qb1, qb2, __func__);
615  validateMultiQubitMatrixFitsInNode(qureg, 2, __func__); // uses 2qb unitary in QuEST_common
616 
617  statevec_sqrtSwapGate(qureg, qb1, qb2);
618  if (qureg.isDensityMatrix) {
619  int shift = qureg.numQubitsRepresented;
620  statevec_sqrtSwapGateConj(qureg, qb1+shift, qb2+shift);
621  }
622 
623  qasm_recordControlledGate(qureg, GATE_SQRT_SWAP, qb1, qb2);
624 }
625 
626 void multiRotateZ(Qureg qureg, int* qubits, int numQubits, qreal angle) {
627  validateMultiTargets(qureg, qubits, numQubits, __func__);
628 
629  long long int mask = getQubitBitMask(qubits, numQubits);
630  statevec_multiRotateZ(qureg, mask, angle);
631  if (qureg.isDensityMatrix) {
632  int shift = qureg.numQubitsRepresented;
633  statevec_multiRotateZ(qureg, mask << shift, -angle);
634  }
635 
636  // @TODO: create actual QASM
637  qasm_recordComment(qureg,
638  "Here a %d-qubit multiRotateZ of angle %g was performed (QASM not yet implemented)",
639  numQubits, angle);
640 }
641 
642 void multiRotatePauli(Qureg qureg, int* targetQubits, enum pauliOpType* targetPaulis, int numTargets, qreal angle) {
643  validateMultiTargets(qureg, targetQubits, numTargets, __func__);
644  validatePauliCodes(targetPaulis, numTargets, __func__);
645 
646  int conj=0;
647  statevec_multiRotatePauli(qureg, targetQubits, targetPaulis, numTargets, angle, conj);
648  if (qureg.isDensityMatrix) {
649  conj = 1;
650  int shift = qureg.numQubitsRepresented;
651  shiftIndices(targetQubits, numTargets, shift);
652  statevec_multiRotatePauli(qureg, targetQubits, targetPaulis, numTargets, angle, conj);
653  shiftIndices(targetQubits, numTargets, -shift);
654  }
655 
656  // @TODO: create actual QASM
657  qasm_recordComment(qureg,
658  "Here a %d-qubit multiRotatePauli of angle %g was performed (QASM not yet implemented)",
659  numTargets, angle);
660 }
661 
662 
663 
664 /*
665  * register attributes
666  */
667 
668 int getNumQubits(Qureg qureg) {
669  return qureg.numQubitsRepresented;
670 }
671 
672 long long int getNumAmps(Qureg qureg) {
673  validateStateVecQureg(qureg, __func__);
674 
675  return qureg.numAmpsTotal;
676 }
677 
678 qreal getRealAmp(Qureg qureg, long long int index) {
679  validateStateVecQureg(qureg, __func__);
680  validateAmpIndex(qureg, index, __func__);
681 
682  return statevec_getRealAmp(qureg, index);
683 }
684 
685 qreal getImagAmp(Qureg qureg, long long int index) {
686  validateStateVecQureg(qureg, __func__);
687  validateAmpIndex(qureg, index, __func__);
688 
689  return statevec_getImagAmp(qureg, index);
690 }
691 
692 qreal getProbAmp(Qureg qureg, long long int index) {
693  validateStateVecQureg(qureg, __func__);
694  validateAmpIndex(qureg, index, __func__);
695 
696  return statevec_getProbAmp(qureg, index);
697 }
698 
699 Complex getAmp(Qureg qureg, long long int index) {
700  validateStateVecQureg(qureg, __func__);
701  validateAmpIndex(qureg, index, __func__);
702 
703  Complex amp;
704  amp.real = statevec_getRealAmp(qureg, index);
705  amp.imag = statevec_getImagAmp(qureg, index);
706  return amp;
707 }
708 
709 Complex getDensityAmp(Qureg qureg, long long int row, long long int col) {
710  validateDensityMatrQureg(qureg, __func__);
711  validateAmpIndex(qureg, row, __func__);
712  validateAmpIndex(qureg, col, __func__);
713 
714  long long ind = row + col*(1LL << qureg.numQubitsRepresented);
715  Complex amp;
716  amp.real = statevec_getRealAmp(qureg, ind);
717  amp.imag = statevec_getImagAmp(qureg, ind);
718  return amp;
719 }
720 
721 
722 /*
723  * non-unitary actions
724  */
725 
726 qreal collapseToOutcome(Qureg qureg, int measureQubit, int outcome) {
727  validateTarget(qureg, measureQubit, __func__);
728  validateOutcome(outcome, __func__);
729 
730  qreal outcomeProb;
731  if (qureg.isDensityMatrix) {
732  outcomeProb = densmatr_calcProbOfOutcome(qureg, measureQubit, outcome);
733  validateMeasurementProb(outcomeProb, __func__);
734  densmatr_collapseToKnownProbOutcome(qureg, measureQubit, outcome, outcomeProb);
735  } else {
736  outcomeProb = statevec_calcProbOfOutcome(qureg, measureQubit, outcome);
737  validateMeasurementProb(outcomeProb, __func__);
738  statevec_collapseToKnownProbOutcome(qureg, measureQubit, outcome, outcomeProb);
739  }
740 
741  qasm_recordMeasurement(qureg, measureQubit);
742  return outcomeProb;
743 }
744 
745 int measureWithStats(Qureg qureg, int measureQubit, qreal *outcomeProb) {
746  validateTarget(qureg, measureQubit, __func__);
747 
748  int outcome;
749  if (qureg.isDensityMatrix)
750  outcome = densmatr_measureWithStats(qureg, measureQubit, outcomeProb);
751  else
752  outcome = statevec_measureWithStats(qureg, measureQubit, outcomeProb);
753 
754  qasm_recordMeasurement(qureg, measureQubit);
755  return outcome;
756 }
757 
758 int measure(Qureg qureg, int measureQubit) {
759  validateTarget(qureg, measureQubit, __func__);
760 
761  int outcome;
762  qreal discardedProb;
763  if (qureg.isDensityMatrix)
764  outcome = densmatr_measureWithStats(qureg, measureQubit, &discardedProb);
765  else
766  outcome = statevec_measureWithStats(qureg, measureQubit, &discardedProb);
767 
768  qasm_recordMeasurement(qureg, measureQubit);
769  return outcome;
770 }
771 
772 void mixDensityMatrix(Qureg combineQureg, qreal otherProb, Qureg otherQureg) {
773  validateDensityMatrQureg(combineQureg, __func__);
774  validateDensityMatrQureg(otherQureg, __func__);
775  validateMatchingQuregDims(combineQureg, otherQureg, __func__);
776  validateProb(otherProb, __func__);
777 
778  densmatr_mixDensityMatrix(combineQureg, otherProb, otherQureg);
779 }
780 
781 void setAmps(Qureg qureg, long long int startInd, qreal* reals, qreal* imags, long long int numAmps) {
782  validateStateVecQureg(qureg, __func__);
783  validateNumAmps(qureg, startInd, numAmps, __func__);
784 
785  statevec_setAmps(qureg, startInd, reals, imags, numAmps);
786 
787  qasm_recordComment(qureg, "Here, some amplitudes in the statevector were manually edited.");
788 }
789 
790 void setDensityAmps(Qureg qureg, qreal* reals, qreal* imags) {
791  long long int numAmps = qureg.numAmpsTotal;
792  statevec_setAmps(qureg, 0, reals, imags, numAmps);
793 
794  qasm_recordComment(qureg, "Here, some amplitudes in the density matrix were manually edited.");
795 }
796 
797 void setWeightedQureg(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out) {
798  validateMatchingQuregTypes(qureg1, qureg2, __func__);
799  validateMatchingQuregTypes(qureg1, out, __func__);
800  validateMatchingQuregDims(qureg1, qureg2, __func__);
801  validateMatchingQuregDims(qureg1, out, __func__);
802 
803  statevec_setWeightedQureg(fac1, qureg1, fac2, qureg2, facOut, out);
804 
805  qasm_recordComment(out, "Here, the register was modified to an undisclosed and possibly unphysical state (setWeightedQureg).");
806 }
807 
808 void applyPauliSum(Qureg inQureg, enum pauliOpType* allPauliCodes, qreal* termCoeffs, int numSumTerms, Qureg outQureg) {
809  validateMatchingQuregTypes(inQureg, outQureg, __func__);
810  validateMatchingQuregDims(inQureg, outQureg, __func__);
811  validateNumPauliSumTerms(numSumTerms, __func__);
812  validatePauliCodes(allPauliCodes, numSumTerms*inQureg.numQubitsRepresented, __func__);
813 
814  statevec_applyPauliSum(inQureg, allPauliCodes, termCoeffs, numSumTerms, outQureg);
815 
816  qasm_recordComment(outQureg, "Here, the register was modified to an undisclosed and possibly unphysical state (applyPauliSum).");
817 }
818 
819 void applyPauliHamil(Qureg inQureg, PauliHamil hamil, Qureg outQureg) {
820  validateMatchingQuregTypes(inQureg, outQureg, __func__);
821  validateMatchingQuregDims(inQureg, outQureg, __func__);
822  validatePauliHamil(hamil, __func__);
823  validateMatchingQuregPauliHamilDims(inQureg, hamil, __func__);
824 
825  statevec_applyPauliSum(inQureg, hamil.pauliCodes, hamil.termCoeffs, hamil.numSumTerms, outQureg);
826 
827  qasm_recordComment(outQureg, "Here, the register was modified to an undisclosed and possibly unphysical state (applyPauliHamil).");
828 }
829 
830 void applyTrotterCircuit(Qureg qureg, PauliHamil hamil, qreal time, int order, int reps) {
831  validateTrotterParams(order, reps, __func__);
832  validatePauliHamil(hamil, __func__);
833  validateMatchingQuregPauliHamilDims(qureg, hamil, __func__);
834 
835  qasm_recordComment(qureg,
836  "Beginning of Trotter circuit (time %g, order %d, %d repetitions).",
837  time, order, reps);
838 
839  agnostic_applyTrotterCircuit(qureg, hamil, time, order, reps);
840 
841  qasm_recordComment(qureg, "End of Trotter circuit");
842 }
843 
844 void applyMatrix2(Qureg qureg, int targetQubit, ComplexMatrix2 u) {
845  validateTarget(qureg, targetQubit, __func__);
846 
847  // actually just left-multiplies any complex matrix
848  statevec_unitary(qureg, targetQubit, u);
849 
850  qasm_recordComment(qureg, "Here, an undisclosed 2-by-2 matrix (possibly non-unitary) was multiplied onto qubit %d", targetQubit);
851 }
852 
853 void applyMatrix4(Qureg qureg, int targetQubit1, int targetQubit2, ComplexMatrix4 u) {
854  validateMultiTargets(qureg, (int []) {targetQubit1, targetQubit2}, 2, __func__);
855  validateMultiQubitMatrixFitsInNode(qureg, 2, __func__);
856 
857  // actually just left-multiplies any complex matrix
858  statevec_twoQubitUnitary(qureg, targetQubit1, targetQubit2, u);
859 
860  qasm_recordComment(qureg, "Here, an undisclosed 4-by-4 matrix (possibly non-unitary) was multiplied onto qubits %d and %d", targetQubit1, targetQubit2);
861 }
862 
863 void applyMatrixN(Qureg qureg, int* targs, int numTargs, ComplexMatrixN u) {
864  validateMultiTargets(qureg, targs, numTargs, __func__);
865  validateMultiQubitMatrix(qureg, u, numTargs, __func__);
866 
867  // actually just left-multiplies any complex matrix
868  statevec_multiQubitUnitary(qureg, targs, numTargs, u);
869 
870  int dim = (1 << numTargs);
871  qasm_recordComment(qureg, "Here, an undisclosed %d-by-%d matrix (possibly non-unitary) was multiplied onto %d undisclosed qubits", dim, dim, numTargs);
872 }
873 
874 void applyMultiControlledMatrixN(Qureg qureg, int* ctrls, int numCtrls, int* targs, int numTargs, ComplexMatrixN u) {
875  validateMultiControlsMultiTargets(qureg, ctrls, numCtrls, targs, numTargs, __func__);
876  validateMultiQubitMatrix(qureg, u, numTargs, __func__);
877 
878  // actually just left-multiplies any complex matrix
879  long long int ctrlMask = getQubitBitMask(ctrls, numCtrls);
880  statevec_multiControlledMultiQubitUnitary(qureg, ctrlMask, targs, numTargs, u);
881 
882  int numTot = numTargs + numCtrls;
883  int dim = (1 << numTot );
884  qasm_recordComment(qureg, "Here, an undisclosed %d-by-%d matrix (possibly non-unitary, and including %d controlled qubits) was multiplied onto %d undisclosed qubits", dim, dim, numCtrls, numTot);
885 }
886 
888  validateDiagonalOp(qureg, op, __func__);
889 
890  if (qureg.isDensityMatrix)
891  densmatr_applyDiagonalOp(qureg, op);
892  else
893  statevec_applyDiagonalOp(qureg, op);
894 
895  qasm_recordComment(qureg, "Here, the register was modified to an undisclosed and possibly unphysical state (via applyDiagonalOp).");
896 }
897 
898 
899 /*
900  * calculations
901  */
902 
904  if (qureg.isDensityMatrix)
905  return densmatr_calcTotalProb(qureg);
906  else
907  return statevec_calcTotalProb(qureg);
908 }
909 
911  validateStateVecQureg(bra, __func__);
912  validateStateVecQureg(ket, __func__);
913  validateMatchingQuregDims(bra, ket, __func__);
914 
915  return statevec_calcInnerProduct(bra, ket);
916 }
917 
919  validateDensityMatrQureg(rho1, __func__);
920  validateDensityMatrQureg(rho2, __func__);
921  validateMatchingQuregDims(rho1, rho2, __func__);
922 
923  return densmatr_calcInnerProduct(rho1, rho2);
924 }
925 
926 qreal calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome) {
927  validateTarget(qureg, measureQubit, __func__);
928  validateOutcome(outcome, __func__);
929 
930  if (qureg.isDensityMatrix)
931  return densmatr_calcProbOfOutcome(qureg, measureQubit, outcome);
932  else
933  return statevec_calcProbOfOutcome(qureg, measureQubit, outcome);
934 }
935 
937  validateDensityMatrQureg(qureg, __func__);
938 
939  return densmatr_calcPurity(qureg);
940 }
941 
942 qreal calcFidelity(Qureg qureg, Qureg pureState) {
943  validateSecondQuregStateVec(pureState, __func__);
944  validateMatchingQuregDims(qureg, pureState, __func__);
945 
946  if (qureg.isDensityMatrix)
947  return densmatr_calcFidelity(qureg, pureState);
948  else
949  return statevec_calcFidelity(qureg, pureState);
950 }
951 
952 qreal calcExpecPauliProd(Qureg qureg, int* targetQubits, enum pauliOpType* pauliCodes, int numTargets, Qureg workspace) {
953  validateMultiTargets(qureg, targetQubits, numTargets, __func__);
954  validatePauliCodes(pauliCodes, numTargets, __func__);
955  validateMatchingQuregTypes(qureg, workspace, __func__);
956  validateMatchingQuregDims(qureg, workspace, __func__);
957 
958  return statevec_calcExpecPauliProd(qureg, targetQubits, pauliCodes, numTargets, workspace);
959 }
960 
961 qreal calcExpecPauliSum(Qureg qureg, enum pauliOpType* allPauliCodes, qreal* termCoeffs, int numSumTerms, Qureg workspace) {
962  validateNumPauliSumTerms(numSumTerms, __func__);
963  validatePauliCodes(allPauliCodes, numSumTerms*qureg.numQubitsRepresented, __func__);
964  validateMatchingQuregTypes(qureg, workspace, __func__);
965  validateMatchingQuregDims(qureg, workspace, __func__);
966 
967  return statevec_calcExpecPauliSum(qureg, allPauliCodes, termCoeffs, numSumTerms, workspace);
968 }
969 
970 qreal calcExpecPauliHamil(Qureg qureg, PauliHamil hamil, Qureg workspace) {
971  validateMatchingQuregTypes(qureg, workspace, __func__);
972  validateMatchingQuregDims(qureg, workspace, __func__);
973  validatePauliHamil(hamil, __func__);
974  validateMatchingQuregPauliHamilDims(qureg, hamil, __func__);
975 
976  return statevec_calcExpecPauliSum(qureg, hamil.pauliCodes, hamil.termCoeffs, hamil.numSumTerms, workspace);
977 }
978 
980  validateDiagonalOp(qureg, op, __func__);
981 
982  if (qureg.isDensityMatrix)
983  return densmatr_calcExpecDiagonalOp(qureg, op);
984  else
985  return statevec_calcExpecDiagonalOp(qureg, op);
986 }
987 
989  validateDensityMatrQureg(a, __func__);
990  validateDensityMatrQureg(b, __func__);
991  validateMatchingQuregDims(a, b, __func__);
992 
994 }
995 
996 
997 /*
998  * decoherence
999  */
1000 
1001 void mixDephasing(Qureg qureg, int targetQubit, qreal prob) {
1002  validateDensityMatrQureg(qureg, __func__);
1003  validateTarget(qureg, targetQubit, __func__);
1004  validateOneQubitDephaseProb(prob, __func__);
1005 
1006  densmatr_mixDephasing(qureg, targetQubit, 2*prob);
1007  qasm_recordComment(qureg,
1008  "Here, a phase (Z) error occured on qubit %d with probability %g", targetQubit, prob);
1009 }
1010 
1011 void mixTwoQubitDephasing(Qureg qureg, int qubit1, int qubit2, qreal prob) {
1012  validateDensityMatrQureg(qureg, __func__);
1013  validateUniqueTargets(qureg, qubit1, qubit2, __func__);
1014  validateTwoQubitDephaseProb(prob, __func__);
1015 
1016  ensureIndsIncrease(&qubit1, &qubit2);
1017  densmatr_mixTwoQubitDephasing(qureg, qubit1, qubit2, (4*prob)/3.0);
1018  qasm_recordComment(qureg,
1019  "Here, a phase (Z) error occured on either or both of qubits "
1020  "%d and %d with total probability %g", qubit1, qubit2, prob);
1021 }
1022 
1023 void mixDepolarising(Qureg qureg, int targetQubit, qreal prob) {
1024  validateDensityMatrQureg(qureg, __func__);
1025  validateTarget(qureg, targetQubit, __func__);
1026  validateOneQubitDepolProb(prob, __func__);
1027 
1028  densmatr_mixDepolarising(qureg, targetQubit, (4*prob)/3.0);
1029  qasm_recordComment(qureg,
1030  "Here, a homogeneous depolarising error (X, Y, or Z) occured on "
1031  "qubit %d with total probability %g", targetQubit, prob);
1032 }
1033 
1034 void mixDamping(Qureg qureg, int targetQubit, qreal prob) {
1035  validateDensityMatrQureg(qureg, __func__);
1036  validateTarget(qureg, targetQubit, __func__);
1037  validateOneQubitDampingProb(prob, __func__);
1038 
1039  densmatr_mixDamping(qureg, targetQubit, prob);
1040 }
1041 
1042 void mixTwoQubitDepolarising(Qureg qureg, int qubit1, int qubit2, qreal prob) {
1043  validateDensityMatrQureg(qureg, __func__);
1044  validateUniqueTargets(qureg, qubit1, qubit2, __func__);
1045  validateTwoQubitDepolProb(prob, __func__);
1046 
1047  ensureIndsIncrease(&qubit1, &qubit2);
1048  densmatr_mixTwoQubitDepolarising(qureg, qubit1, qubit2, (16*prob)/15.0);
1049  qasm_recordComment(qureg,
1050  "Here, a homogeneous depolarising error occured on qubits %d and %d "
1051  "with total probability %g", qubit1, qubit2, prob);
1052 }
1053 
1054 void mixPauli(Qureg qureg, int qubit, qreal probX, qreal probY, qreal probZ) {
1055  validateDensityMatrQureg(qureg, __func__);
1056  validateTarget(qureg, qubit, __func__);
1057  validateOneQubitPauliProbs(probX, probY, probZ, __func__);
1058 
1059  densmatr_mixPauli(qureg, qubit, probX, probY, probZ);
1060  qasm_recordComment(qureg,
1061  "Here, X, Y and Z errors occured on qubit %d with probabilities "
1062  "%g, %g and %g respectively", qubit, probX, probY, probZ);
1063 }
1064 
1065 void mixKrausMap(Qureg qureg, int target, ComplexMatrix2 *ops, int numOps) {
1066  validateDensityMatrQureg(qureg, __func__);
1067  validateTarget(qureg, target, __func__);
1068  validateOneQubitKrausMap(qureg, ops, numOps, __func__);
1069 
1070  densmatr_mixKrausMap(qureg, target, ops, numOps);
1071  qasm_recordComment(qureg,
1072  "Here, an undisclosed Kraus map was effected on qubit %d", target);
1073 }
1074 
1075 void mixTwoQubitKrausMap(Qureg qureg, int target1, int target2, ComplexMatrix4 *ops, int numOps) {
1076  validateDensityMatrQureg(qureg, __func__);
1077  validateMultiTargets(qureg, (int[]) {target1,target2}, 2, __func__);
1078  validateTwoQubitKrausMap(qureg, ops, numOps, __func__);
1079 
1080  densmatr_mixTwoQubitKrausMap(qureg, target1, target2, ops, numOps);
1081  qasm_recordComment(qureg,
1082  "Here, an undisclosed two-qubit Kraus map was effected on qubits %d and %d", target1, target2);
1083 }
1084 
1085 void mixMultiQubitKrausMap(Qureg qureg, int* targets, int numTargets, ComplexMatrixN* ops, int numOps) {
1086  validateDensityMatrQureg(qureg, __func__);
1087  validateMultiTargets(qureg, targets, numTargets, __func__);
1088  validateMultiQubitKrausMap(qureg, numTargets, ops, numOps, __func__);
1089 
1090  densmatr_mixMultiQubitKrausMap(qureg, targets, numTargets, ops, numOps);
1091  qasm_recordComment(qureg,
1092  "Here, an undisclosed %d-qubit Kraus map was applied to undisclosed qubits", numTargets);
1093 }
1094 
1095 /*
1096  * other data structures
1097  */
1098 
1100  validateNumQubitsInMatrix(numQubits, __func__);
1101 
1102  int numRows = 1 << numQubits;
1103 
1104  ComplexMatrixN m = {
1105  .numQubits = numQubits,
1106  .real = malloc(numRows * sizeof *m.real),
1107  .imag = malloc(numRows * sizeof *m.imag)};
1108 
1109  for (int n=0; n < 1<<numQubits; n++) {
1110  m.real[n] = calloc(numRows, sizeof **m.real);
1111  m.imag[n] = calloc(numRows, sizeof **m.imag);
1112  }
1113 
1114  // error if the ComplexMatrixN was not successfully malloc'ds
1115  validateMatrixInit(m, __func__);
1116 
1117  return m;
1118  }
1119 
1121  /* this checks m.real/imag != NULL, which is only ever set when the mallocs
1122  * in createComplexMatrixN fail, which already prompts an error. Hence
1123  * this check if useless
1124  */
1125  validateMatrixInit(m, __func__);
1126 
1127  int numRows = 1 << m.numQubits;
1128  for (int r=0; r < numRows; r++) {
1129  free(m.real[r]);
1130  free(m.imag[r]);
1131  }
1132  free(m.real);
1133  free(m.imag);
1134 }
1135 
1137  validateMatrixInit(m, __func__);
1138 
1139  int dim = 1 << m.numQubits;
1140  for (int i=0; i<dim; i++)
1141  for (int j=0; j<dim; j++) {
1142  m.real[i][j] = re[i][j];
1143  m.imag[i][j] = im[i][j];
1144  }
1145 }
1146 
1147 PauliHamil createPauliHamil(int numQubits, int numSumTerms) {
1148  validateHamilParams(numQubits, numSumTerms, __func__);
1149 
1150  PauliHamil h;
1151  h.numQubits = numQubits;
1152  h.numSumTerms = numSumTerms;
1153  h.termCoeffs = malloc(numSumTerms * sizeof *h.termCoeffs);
1154  h.pauliCodes = malloc(numQubits*numSumTerms * sizeof *h.pauliCodes);
1155 
1156  // initialise pauli codes to identity
1157  for (int i=0; i<numQubits*numSumTerms; i++)
1158  h.pauliCodes[i] = PAULI_I;
1159 
1160  return h;
1161 }
1162 
1164 
1165  free(h.termCoeffs);
1166  free(h.pauliCodes);
1167 }
1168 
1170 
1171  /* The validation in this function must close the file handle and free
1172  * allocated memory before raising an error (whether that's a C exit, or
1173  * an overriden C++ exception).
1174  */
1175 
1176  FILE* file = fopen(fn, "r");
1177  int success = (file != NULL);
1178  validateFileOpened(success, fn, __func__);
1179 
1180  /* file format: coeff {term} \n where {term} is #numQubits values of
1181  * 0 1 2 3 signifying I X Y Z acting on that qubit index
1182  */
1183 
1184  // count the number of qubits (ignore trailing whitespace)
1185  int numQubits = -1; // to exclude coeff at start
1186  char ch = getc(file);
1187  char prev = '0'; // anything not space
1188  while (ch != '\n' && ch != EOF) {
1189  if (ch == ' ' && prev != ' ') // skip multiple spaces
1190  numQubits++;
1191  prev = ch;
1192  ch = getc(file);
1193  }
1194  // edge-case: if we hit EOF/newline without a space
1195  if (prev != ' ')
1196  numQubits++;
1197 
1198  /* TODO:
1199  * The below code may break on Windows where newlines are multiple characters
1200  */
1201 
1202  // count the number of terms (being cautious of trailing newlines)
1203  int numTerms = 0;
1204  prev = '\n';
1205  rewind(file);
1206  while ((ch=getc(file)) != EOF) {
1207  if (ch == '\n' && prev != '\n')
1208  numTerms++;
1209  prev = ch;
1210  }
1211  // edge-case: if we hit EOF without a newline, count that line
1212  if (prev != '\n')
1213  numTerms++;
1214 
1215  // validate the inferred number of terms and qubits (closes file if error)
1216  validateHamilFileParams(numQubits, numTerms, file, fn, __func__);
1217 
1218  // allocate space for PauliHamil data
1219  PauliHamil h = createPauliHamil(numQubits, numTerms);
1220 
1221  // specifier for a qreal number then a space
1222  char strSpec[50];
1223  strcpy(strSpec, REAL_SPECIFIER);
1224  strcat(strSpec, " ");
1225 
1226  // collect coefficients and terms
1227  rewind(file);
1228  for (int t=0; t<numTerms; t++) {
1229 
1230  // record coefficient, and validate (closes file and frees h if error)
1231  success = fscanf(file, strSpec, &(h.termCoeffs[t])) == 1;
1232  validateHamilFileCoeffParsed(success, h, file, fn, __func__);
1233 
1234  // record Pauli operations, and validate (closes file and frees h if error)
1235  for (int q=0; q<numQubits; q++) {
1236  int i = t*numQubits + q;
1237 
1238  // verbose, to avoid type warnings
1239  int code;
1240  success = fscanf(file, "%d ", &code) == 1;
1241  h.pauliCodes[i] = (enum pauliOpType) code;
1242  validateHamilFilePauliParsed(success, h, file, fn, __func__);
1243  validateHamilFilePauliCode(h.pauliCodes[i], h, file, fn, __func__);
1244  }
1245 
1246  // the trailing newline is magically eaten
1247  }
1248 
1249  fclose(file);
1250  return h;
1251 }
1252 
1253 void initPauliHamil(PauliHamil hamil, qreal* coeffs, enum pauliOpType* codes) {
1254  validateHamilParams(hamil.numQubits, hamil.numSumTerms, __func__);
1255  validatePauliCodes(codes, hamil.numSumTerms*hamil.numQubits, __func__);
1256 
1257  int i=0;
1258  for (int t=0; t<hamil.numSumTerms; t++) {
1259  hamil.termCoeffs[t] = coeffs[t];
1260  for (int q=0; q<hamil.numQubits; q++) {
1261  hamil.pauliCodes[i] = codes[i];
1262  i++;
1263  }
1264  }
1265 }
1266 
1267 DiagonalOp createDiagonalOp(int numQubits, QuESTEnv env) {
1268  validateNumQubitsInDiagOp(numQubits, env.numRanks, __func__);
1269 
1270  return agnostic_createDiagonalOp(numQubits, env);
1271 }
1272 
1274  // env accepted for API consistency
1275  validateDiagOpInit(op, __func__);
1276 
1278 }
1279 
1281  validateDiagOpInit(op, __func__);
1282 
1284 }
1285 
1286 void initDiagonalOp(DiagonalOp op, qreal* real, qreal* imag) {
1287  validateDiagOpInit(op, __func__);
1288 
1289  agnostic_setDiagonalOpElems(op, 0, real, imag, 1LL << op.numQubits);
1290 }
1291 
1292 void setDiagonalOpElems(DiagonalOp op, long long int startInd, qreal* real, qreal* imag, long long int numElems) {
1293  validateDiagOpInit(op, __func__);
1294  validateNumElems(op, startInd, numElems, __func__);
1295 
1296  agnostic_setDiagonalOpElems(op, startInd, real, imag, numElems);
1297 }
1298 
1299 /*
1300  * debug
1301  */
1302 
1303 int compareStates(Qureg qureg1, Qureg qureg2, qreal precision) {
1304  validateMatchingQuregDims(qureg1, qureg2, __func__);
1305  return statevec_compareStates(qureg1, qureg2, precision);
1306 }
1307 
1308 void initDebugState(Qureg qureg) {
1309  statevec_initDebugState(qureg);
1310 }
1311 
1312 void initStateFromSingleFile(Qureg *qureg, char filename[200], QuESTEnv env) {
1313  int success = statevec_initStateFromSingleFile(qureg, filename, env);
1314  validateFileOpened(success, filename, __func__);
1315 }
1316 
1317 void initStateOfSingleQubit(Qureg *qureg, int qubitId, int outcome) {
1318  validateStateVecQureg(*qureg, __func__);
1319  validateTarget(*qureg, qubitId, __func__);
1320  validateOutcome(outcome, __func__);
1321  statevec_initStateOfSingleQubit(qureg, qubitId, outcome);
1322 }
1323 
1324 void reportStateToScreen(Qureg qureg, QuESTEnv env, int reportRank) {
1325  statevec_reportStateToScreen(qureg, env, reportRank);
1326 }
1327 
1329  validatePauliHamil(hamil, __func__);
1330 
1331  for (int t=0; t<hamil.numSumTerms; t++) {
1332  printf("%g\t", hamil.termCoeffs[t]);
1333  for (int q=0; q<hamil.numQubits; q++)
1334  printf("%d ", (int) hamil.pauliCodes[q+t*hamil.numQubits]);
1335  printf("\n");
1336  }
1337 }
1338 
1339 int getQuEST_PREC(void) {
1340  return sizeof(qreal)/4;
1341 }
1342 
1343 
1344 #ifdef __cplusplus
1345 }
1346 #endif
void agnostic_destroyDiagonalOp(DiagonalOp op)
Definition: QuEST_cpu.c:1357
void densmatr_mixDepolarising(Qureg qureg, int targetQubit, qreal depolLevel)
void controlledRotateZ(Qureg qureg, int controlQubit, int targetQubit, qreal angle)
Applies a controlled rotation by a given angle around the Z-axis of the Bloch-sphere.
Definition: QuEST.c:245
void qasm_printRecorded(Qureg qureg)
Definition: QuEST_qasm.c:484
void statevec_phaseShift(Qureg qureg, int targetQubit, qreal angle)
Definition: QuEST_common.c:251
void mixDephasing(Qureg qureg, int targetQubit, qreal prob)
Mixes a density matrix qureg to induce single-qubit dephasing noise.
Definition: QuEST.c:1001
int compareStates(Qureg qureg1, Qureg qureg2, qreal precision)
Return whether two given wavefunctions are equivalent within a given precision Global phase included ...
Definition: QuEST.c:1303
qreal getProbAmp(Qureg qureg, long long int index)
Get the probability of a state-vector at an index in the full state vector.
Definition: QuEST.c:692
void validateDensityMatrQureg(Qureg qureg, const char *caller)
void initBlankState(Qureg qureg)
Initialises a qureg to have all-zero-amplitudes.
Definition: QuEST.c:119
Represents a 3-vector of real numbers.
Definition: QuEST.h:148
void validateMultiControlsTarget(Qureg qureg, int *controlQubits, int numControlQubits, int targetQubit, const char *caller)
void statevec_sGate(Qureg qureg, int targetQubit)
Definition: QuEST_common.c:265
void statevec_pauliZ(Qureg qureg, int targetQubit)
Definition: QuEST_common.c:258
pauliOpType
Codes for specifying Pauli operators.
Definition: QuEST.h:96
void statevec_rotateX(Qureg qureg, int rotQubit, qreal angle)
Definition: QuEST_common.c:293
void validateMeasurementProb(qreal prob, const char *caller)
void validateTarget(Qureg qureg, int targetQubit, const char *caller)
void applyMultiControlledMatrixN(Qureg qureg, int *ctrls, int numCtrls, int *targs, int numTargs, ComplexMatrixN u)
Apply a general N-by-N matrix, which may be non-unitary, with additional controlled qubits.
Definition: QuEST.c:874
void twoQubitUnitary(Qureg qureg, int targetQubit1, int targetQubit2, ComplexMatrix4 u)
Apply a general two-qubit unitary (including a global phase factor).
Definition: QuEST.c:257
void controlledRotateX(Qureg qureg, int controlQubit, int targetQubit, qreal angle)
Applies a controlled rotation by a given angle around the X-axis of the Bloch-sphere.
Definition: QuEST.c:221
void initPauliHamil(PauliHamil hamil, qreal *coeffs, enum pauliOpType *codes)
Initialise a PauliHamil instance with the given term coefficients and Pauli codes (one for every qubi...
Definition: QuEST.c:1253
void validateOutcome(int outcome, const char *caller)
void densmatr_mixKrausMap(Qureg qureg, int target, ComplexMatrix2 *ops, int numOps)
Definition: QuEST_common.c:600
void initPureState(Qureg qureg, Qureg pure)
Initialise a set of qubits, which can be a state vector or density matrix, to a given pure state.
Definition: QuEST.c:145
void rotateAroundAxis(Qureg qureg, int rotQubit, qreal angle, Vector axis)
Rotate a single qubit by a given angle around a given Vector on the Bloch-sphere.
Definition: QuEST.c:575
void mixDepolarising(Qureg qureg, int targetQubit, qreal prob)
Mixes a density matrix qureg to induce single-qubit homogeneous depolarising noise.
Definition: QuEST.c:1023
void validateHamilFileCoeffParsed(int parsed, PauliHamil h, FILE *file, char *fn, const char *caller)
void qasm_recordParamGate(Qureg qureg, TargetGate gate, int targetQubit, qreal param)
Definition: QuEST_qasm.c:186
void validateStateIndex(Qureg qureg, long long int stateInd, const char *caller)
void densmatr_initPlusState(Qureg targetQureg)
Definition: QuEST_cpu.c:1154
qreal calcTotalProb(Qureg qureg)
A debugging function which calculates the probability of the qubits in qureg being in any state,...
Definition: QuEST.c:903
void mixDamping(Qureg qureg, int targetQubit, qreal prob)
Mixes a density matrix qureg to induce single-qubit amplitude damping (decay to 0 state).
Definition: QuEST.c:1034
void destroyComplexMatrixN(ComplexMatrixN m)
Destroy a ComplexMatrixN instance created with createComplexMatrixN()
Definition: QuEST.c:1120
void shiftIndices(int *indices, int numIndices, int shift)
Definition: QuEST_common.c:150
qreal densmatr_calcInnerProduct(Qureg a, Qureg b)
qreal statevec_calcExpecPauliProd(Qureg qureg, int *targetQubits, enum pauliOpType *pauliCodes, int numTargets, Qureg workspace)
Definition: QuEST_common.c:465
void reportPauliHamil(PauliHamil hamil)
Print the PauliHamil to screen.
Definition: QuEST.c:1328
void validateHamilParams(int numQubits, int numTerms, const char *caller)
void qasm_free(Qureg qureg)
Definition: QuEST_qasm.c:500
void statevec_controlledTwoQubitUnitary(Qureg qureg, int controlQubit, int targetQubit1, int targetQubit2, ComplexMatrix4 u)
Definition: QuEST_common.c:523
@ PAULI_I
Definition: QuEST.h:96
Complex getDensityAmp(Qureg qureg, long long int row, long long int col)
Get an amplitude from a density matrix at a given row and column.
Definition: QuEST.c:709
@ GATE_T
Definition: QuEST_qasm.h:24
@ GATE_PHASE_SHIFT
Definition: QuEST_qasm.h:32
void qasm_recordUnitary(Qureg qureg, ComplexMatrix2 u, int targetQubit)
Definition: QuEST_qasm.c:207
void validateStateVecQureg(Qureg qureg, const char *caller)
void qasm_recordInitZero(Qureg qureg)
Definition: QuEST_qasm.c:415
ComplexMatrixN createComplexMatrixN(int numQubits)
Create (dynamically) a square complex matrix which can be passed to the multi-qubit general unitary f...
Definition: QuEST.c:1099
qreal statevec_calcExpecPauliSum(Qureg qureg, enum pauliOpType *allCodes, qreal *termCoeffs, int numSumTerms, Qureg workspace)
Definition: QuEST_common.c:480
void qasm_recordMultiControlledGate(Qureg qureg, TargetGate gate, int *controlQubits, int numControlQubits, int targetQubit)
Definition: QuEST_qasm.c:316
qreal getImagAmp(Qureg qureg, long long int index)
Get the imaginary component of the complex probability amplitude at an index in the state vector.
Definition: QuEST.c:685
void validateMultiQubitMatrixFitsInNode(Qureg qureg, int numTargets, const char *caller)
void controlledRotateAroundAxis(Qureg qureg, int controlQubit, int targetQubit, qreal angle, Vector axis)
Applies a controlled rotation by a given angle around a given vector on the Bloch-sphere.
Definition: QuEST.c:588
void mixKrausMap(Qureg qureg, int target, ComplexMatrix2 *ops, int numOps)
Apply a general single-qubit Kraus map to a density matrix, as specified by at most four Kraus operat...
Definition: QuEST.c:1065
void statevec_twoQubitUnitary(Qureg qureg, int targetQubit1, int targetQubit2, ComplexMatrix4 u)
Definition: QuEST_common.c:517
void statevec_tGateConj(Qureg qureg, int targetQubit)
Definition: QuEST_common.c:286
qreal densmatr_calcPurity(Qureg qureg)
Computes the trace of the density matrix squared.
void syncDiagonalOp(DiagonalOp op)
Copy the elements in DiagonalOp op.real and op.imag to the persisent GPU memory.
Definition: QuEST.c:1280
void destroyDiagonalOp(DiagonalOp op, QuESTEnv env)
Destroys a DiagonalOp created with createDiagonalOp(), freeing its memory.
Definition: QuEST.c:1273
void multiControlledTwoQubitUnitary(Qureg qureg, int *controlQubits, int numControlQubits, int targetQubit1, int targetQubit2, ComplexMatrix4 u)
Apply a general multi-controlled two-qubit unitary (including a global phase factor).
Definition: QuEST.c:283
void validateNumAmps(Qureg qureg, long long int startInd, long long int numAmps, const char *caller)
void qasm_clearRecorded(Qureg qureg)
Definition: QuEST_qasm.c:477
void statevec_unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u)
void qasm_recordInitPlus(Qureg qureg)
Definition: QuEST_qasm.c:430
qreal collapseToOutcome(Qureg qureg, int measureQubit, int outcome)
Updates qureg to be consistent with measuring measureQubit in the given outcome (0 or 1),...
Definition: QuEST.c:726
ComplexMatrix4 getConjugateMatrix4(ComplexMatrix4 src)
Definition: QuEST_common.c:104
@ GATE_ROTATE_X
Definition: QuEST_qasm.h:27
void validateNumQubitsInQureg(int numQubits, int numRanks, const char *caller)
void validateMultiQubitUnitaryMatrix(Qureg qureg, ComplexMatrixN u, int numTargs, const char *caller)
qreal calcPurity(Qureg qureg)
Calculates the purity of a density matrix, by the trace of the density matrix squared.
Definition: QuEST.c:936
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 ...
void validateHamilFilePauliParsed(int parsed, PauliHamil h, FILE *file, char *fn, const char *caller)
void statevec_multiRotateZ(Qureg qureg, long long int mask, qreal angle)
Definition: QuEST_cpu.c:3109
qreal calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
Gives the probability of a specified qubit being measured in the given outcome (0 or 1).
Definition: QuEST.c:926
void qasm_recordControlledCompactUnitary(Qureg qureg, Complex alpha, Complex beta, int controlQubit, int targetQubit)
Definition: QuEST_qasm.c:264
void agnostic_syncDiagonalOp(DiagonalOp op)
Definition: QuEST_cpu.c:1362
void statevec_destroyQureg(Qureg qureg, QuESTEnv env)
Definition: QuEST_cpu.c:1317
qreal statevec_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
void clearRecordedQASM(Qureg qureg)
Clear all QASM so far recorded.
Definition: QuEST.c:95
int measure(Qureg qureg, int measureQubit)
Measures a single qubit, collapsing it randomly to 0 or 1.
Definition: QuEST.c:758
@ GATE_ROTATE_Z
Definition: QuEST_qasm.h:29
qreal calcFidelity(Qureg qureg, Qureg pureState)
Calculates the fidelity of qureg (a statevector or density matrix) against a reference pure state (ne...
Definition: QuEST.c:942
void validateHamilFilePauliCode(enum pauliOpType code, PauliHamil h, FILE *file, char *fn, const char *caller)
void unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Apply a general single-qubit unitary (including a global phase factor).
Definition: QuEST.c:349
void validateProb(qreal prob, const char *caller)
@ GATE_SIGMA_Z
Definition: QuEST_qasm.h:23
void statevec_controlledPauliY(Qureg qureg, int controlQubit, int targetQubit)
void mixTwoQubitDepolarising(Qureg qureg, int qubit1, int qubit2, qreal prob)
Mixes a density matrix qureg to induce two-qubit homogeneous depolarising noise.
Definition: QuEST.c:1042
void densmatr_mixMultiQubitKrausMap(Qureg qureg, int *targets, int numTargets, ComplexMatrixN *ops, int numOps)
Definition: QuEST_common.c:643
Complex getConjugateScalar(Complex scalar)
Definition: QuEST_common.c:85
void statevec_controlledPhaseShift(Qureg qureg, int idQubit1, int idQubit2, qreal angle)
Definition: QuEST_cpu.c:3019
Complex calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
Computes the expected value of the diagonal operator op for state qureg.
Definition: QuEST.c:979
@ GATE_HADAMARD
Definition: QuEST_qasm.h:26
ComplexMatrix2 getConjugateMatrix2(ComplexMatrix2 src)
Definition: QuEST_common.c:99
void statevec_initStateOfSingleQubit(Qureg *qureg, int qubitId, int outcome)
Initialise the state vector of probability amplitudes such that one qubit is set to 'outcome' and all...
Definition: QuEST_cpu.c:1545
void sGate(Qureg qureg, int targetQubit)
Apply the single-qubit S gate.
Definition: QuEST.c:466
void qasm_startRecording(Qureg qureg)
Definition: QuEST_qasm.c:84
void cloneQureg(Qureg targetQureg, Qureg copyQureg)
Set targetQureg to be a clone of copyQureg.
Definition: QuEST.c:165
Represents a 4x4 matrix of complex numbers.
Definition: QuEST.h:125
void rotateY(Qureg qureg, int targetQubit, qreal angle)
Rotate a single qubit by a given angle around the Y-axis of the Bloch-sphere.
Definition: QuEST.c:199
void densmatr_initClassicalState(Qureg qureg, long long int stateInd)
Definition: QuEST_cpu.c:1115
void applyPauliHamil(Qureg inQureg, PauliHamil hamil, Qureg outQureg)
Modifies outQureg to be the result of applying PauliHamil (a Hermitian but not necessarily unitary op...
Definition: QuEST.c:819
Information about the environment the program is running in.
Definition: QuEST.h:242
void setDensityAmps(Qureg qureg, qreal *reals, qreal *imags)
Set elements in the underlying state vector represenation of a density matrix.
Definition: QuEST.c:790
PauliHamil createPauliHamilFromFile(char *fn)
Create a PauliHamil instance, a real-weighted sum of products of Pauli operators, populated with the ...
Definition: QuEST.c:1169
void multiControlledPhaseShift(Qureg qureg, int *controlQubits, int numControlQubits, qreal angle)
Introduce a phase factor on state of the passed qubits.
Definition: QuEST.c:511
void setDiagonalOpElems(DiagonalOp op, long long int startInd, qreal *real, qreal *imag, long long int numElems)
Modifies a subset (starting at index startInd) of the elements in DiagonalOp op with the given elemen...
Definition: QuEST.c:1292
void statevec_multiControlledPhaseShift(Qureg qureg, int *controlQubits, int numControlQubits, qreal angle)
Definition: QuEST_cpu.c:3059
void statevec_pauliY(Qureg qureg, int targetQubit)
Represents a general 2^N by 2^N matrix of complex numbers.
Definition: QuEST.h:136
void validateTwoQubitDepolProb(qreal prob, const char *caller)
void statevec_controlledMultiQubitUnitary(Qureg qureg, int ctrl, int *targets, int numTargets, ComplexMatrixN u)
Definition: QuEST_common.c:535
void validateNumPauliSumTerms(int numTerms, const char *caller)
qreal statevec_calcFidelity(Qureg qureg, Qureg pureState)
Definition: QuEST_common.c:377
#define qreal
void validateMatrixInit(ComplexMatrixN matr, const char *caller)
void multiRotatePauli(Qureg qureg, int *targetQubits, enum pauliOpType *targetPaulis, int numTargets, qreal angle)
Apply a multi-qubit multi-Pauli rotation on a selected number of qubits.
Definition: QuEST.c:642
void validateFileOpened(int opened, char *fn, const char *caller)
void qasm_recordMultiStateControlledUnitary(Qureg qureg, ComplexMatrix2 u, int *controlQubits, int *controlState, int numControlQubits, int targetQubit)
Definition: QuEST_qasm.c:362
void stopRecordingQASM(Qureg qureg)
Disable QASM recording.
Definition: QuEST.c:91
void setAmps(Qureg qureg, long long int startInd, qreal *reals, qreal *imags, long long int numAmps)
Overwrites a subset of the amplitudes in qureg, with those passed in reals and imags.
Definition: QuEST.c:781
void multiControlledPhaseFlip(Qureg qureg, int *controlQubits, int numControlQubits)
Apply the multiple-qubit controlled phase flip gate, also known as the multiple-qubit controlled paul...
Definition: QuEST.c:561
void multiQubitUnitary(Qureg qureg, int *targs, int numTargs, ComplexMatrixN u)
Apply a general multi-qubit unitary (including a global phase factor) with any number of target qubit...
Definition: QuEST.c:297
void statevec_sqrtSwapGate(Qureg qureg, int qb1, int qb2)
Definition: QuEST_common.c:384
Complex getAmp(Qureg qureg, long long int index)
Get the complex amplitude at a given index in the state vector.
Definition: QuEST.c:699
qreal statevec_getProbAmp(Qureg qureg, long long int index)
Definition: QuEST_common.c:245
void controlledPauliY(Qureg qureg, int controlQubit, int targetQubit)
Apply the controlled pauliY (single control, single target) gate, also known as the c-Y and c-sigma-Y...
Definition: QuEST.c:537
DiagonalOp agnostic_createDiagonalOp(int numQubits, QuESTEnv env)
Definition: QuEST_cpu.c:1335
void initComplexMatrixN(ComplexMatrixN m, qreal re[][1<< m.numQubits], qreal im[][1<< m.numQubits])
Initialises a ComplexMatrixN instance to have the passed real and imag values.
Definition: QuEST.c:1136
void statevec_swapQubitAmps(Qureg qureg, int qb1, int qb2)
void statevec_controlledPauliYConj(Qureg qureg, int controlQubit, int targetQubit)
int measureWithStats(Qureg qureg, int measureQubit, qreal *outcomeProb)
Measures a single qubit, collapsing it randomly to 0 or 1, and additionally gives the probability of ...
Definition: QuEST.c:745
int numQubitsInStateVec
Number of qubits in the state-vector - this is double the number represented for mixed states.
Definition: QuEST.h:210
qreal calcExpecPauliHamil(Qureg qureg, PauliHamil hamil, Qureg workspace)
Computes the expected value of qureg under Hermitian operator hamil.
Definition: QuEST.c:970
void validateVector(Vector vec, const char *caller)
void multiControlledMultiQubitUnitary(Qureg qureg, int *ctrls, int numCtrls, int *targs, int numTargs, ComplexMatrixN u)
Apply a general multi-controlled multi-qubit unitary (including a global phase factor).
Definition: QuEST.c:331
qreal densmatr_calcTotalProb(Qureg qureg)
void rotateX(Qureg qureg, int targetQubit, qreal angle)
Rotate a single qubit by a given angle around the X-axis of the Bloch-sphere.
Definition: QuEST.c:188
void validateMultiQubits(Qureg qureg, int *qubits, int numQubits, const char *caller)
void densmatr_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal outcomeProb)
Renorms (/prob) every | * outcome * >< * outcome * | state, setting all others to zero.
Definition: QuEST_cpu.c:785
void statevec_controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
Complex statevec_calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
void setWeightedQureg(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out)
Modifies qureg out to the result of (facOut out + fac1 qureg1 + fac2 qureg2), imposing no constraints...
Definition: QuEST.c:797
void qasm_recordControlledGate(Qureg qureg, TargetGate gate, int controlQubit, int targetQubit)
Definition: QuEST_qasm.c:238
qreal calcHilbertSchmidtDistance(Qureg a, Qureg b)
Computes the Hilbert Schmidt distance between two density matrices a and b, defined as the Frobenius ...
Definition: QuEST.c:988
void validateControlTarget(Qureg qureg, int controlQubit, int targetQubit, const char *caller)
void qasm_recordAxisRotation(Qureg qureg, qreal angle, Vector axis, int targetQubit)
Definition: QuEST_qasm.c:223
void validateDiagonalOp(Qureg qureg, DiagonalOp op, const char *caller)
void densmatr_mixDensityMatrix(Qureg combineQureg, qreal otherProb, Qureg otherQureg)
Definition: QuEST_cpu.c:890
@ GATE_SQRT_SWAP
Definition: QuEST_qasm.h:34
void pauliZ(Qureg qureg, int targetQubit)
Apply the single-qubit Pauli-Z (also known as the Z, sigma-Z or phase-flip) gate.
Definition: QuEST.c:455
@ GATE_SIGMA_X
Definition: QuEST_qasm.h:21
void qasm_recordMeasurement(Qureg qureg, int measureQubit)
Definition: QuEST_qasm.c:398
void qasm_stopRecording(Qureg qureg)
Definition: QuEST_qasm.c:88
void statevec_initZeroState(Qureg qureg)
Definition: QuEST_cpu.c:1428
void applyDiagonalOp(Qureg qureg, DiagonalOp op)
Apply a diagonal complex operator, which is possibly non-unitary and non-Hermitian,...
Definition: QuEST.c:887
qreal * termCoeffs
The coefficient of each Pauli product. This is a length numSumTerms array.
Definition: QuEST.h:164
void applyMatrix4(Qureg qureg, int targetQubit1, int targetQubit2, ComplexMatrix4 u)
Apply a general 4-by-4 matrix, which may be non-unitary.
Definition: QuEST.c:853
void initStateOfSingleQubit(Qureg *qureg, int qubitId, int outcome)
Initialise the state vector of probability amplitudes such that one qubit is set to 'outcome' and all...
Definition: QuEST.c:1317
void agnostic_applyTrotterCircuit(Qureg qureg, PauliHamil hamil, qreal time, int order, int reps)
Definition: QuEST_common.c:773
void densmatr_mixDephasing(Qureg qureg, int targetQubit, qreal dephase)
Definition: QuEST_cpu.c:79
void validateOneQubitDephaseProb(qreal prob, const char *caller)
void swapGate(Qureg qureg, int qb1, int qb2)
Performs a SWAP gate between qubit1 and qubit2.
Definition: QuEST.c:601
enum pauliOpType * pauliCodes
The Pauli operators acting on each qubit, flattened over every operator.
Definition: QuEST.h:162
void statevec_tGate(Qureg qureg, int targetQubit)
Definition: QuEST_common.c:272
void statevec_initBlankState(Qureg qureg)
Definition: QuEST_cpu.c:1398
void initStateFromAmps(Qureg qureg, qreal *reals, qreal *imags)
Initialise qureg by specifying the complete statevector.
Definition: QuEST.c:157
void statevec_reportStateToScreen(Qureg qureg, QuESTEnv env, int reportRank)
Print the current state vector of probability amplitudes for a set of qubits to standard out.
Definition: QuEST_cpu.c:1366
void validateControlState(int *controlState, int numControlQubits, const char *caller)
void statevec_applyDiagonalOp(Qureg qureg, DiagonalOp op)
Definition: QuEST_cpu.c:3661
void validateMultiTargets(Qureg qureg, int *targetQubits, int numTargetQubits, const char *caller)
void qasm_recordInitClassical(Qureg qureg, long long int stateInd)
Definition: QuEST_qasm.c:458
void statevec_setAmps(Qureg qureg, long long int startInd, qreal *reals, qreal *imags, long long int numAmps)
Definition: QuEST_cpu.c:1237
int qasm_writeRecordedToFile(Qureg qureg, char *filename)
returns success of file write
Definition: QuEST_qasm.c:489
void qasm_recordControlledParamGate(Qureg qureg, TargetGate gate, int controlQubit, int targetQubit, qreal param)
Definition: QuEST_qasm.c:247
void statevec_cloneQureg(Qureg targetQureg, Qureg copyQureg)
works for both statevectors and density matrices
Definition: QuEST_cpu.c:1506
void statevec_rotateAroundAxis(Qureg qureg, int rotQubit, qreal angle, Vector axis)
Definition: QuEST_common.c:311
void writeRecordedQASMToFile(Qureg qureg, char *filename)
Writes recorded QASM to a file, throwing an error if inaccessible.
Definition: QuEST.c:103
int numRanks
Definition: QuEST.h:245
void densmatr_mixTwoQubitDephasing(Qureg qureg, int qubit1, int qubit2, qreal dephase)
Definition: QuEST_cpu.c:84
void setConjugateMatrixN(ComplexMatrixN m)
Definition: QuEST_common.c:109
void validateMatchingQuregDims(Qureg qureg1, Qureg qureg2, const char *caller)
qreal getRealAmp(Qureg qureg, long long int index)
Get the real component of the complex probability amplitude at an index in the state vector.
Definition: QuEST.c:678
void controlledPhaseFlip(Qureg qureg, int idQubit1, int idQubit2)
Apply the (two-qubit) controlled phase flip gate, also known as the controlled pauliZ gate.
Definition: QuEST.c:549
void densmatr_mixTwoQubitKrausMap(Qureg qureg, int target1, int target2, ComplexMatrix4 *ops, int numOps)
Definition: QuEST_common.c:635
void statevec_multiControlledTwoQubitUnitary(Qureg qureg, long long int ctrlMask, int targetQubit1, int targetQubit2, ComplexMatrix4 u)
This calls swapQubitAmps only when it would involve a distributed communication; if the qubit chunks ...
void statevec_initPlusState(Qureg qureg)
Definition: QuEST_cpu.c:1438
void qasm_recordControlledAxisRotation(Qureg qureg, qreal angle, Vector axis, int controlQubit, int targetQubit)
Definition: QuEST_qasm.c:300
int numQubits
The number of qubits this operator can act on (informing its size)
Definition: QuEST.h:181
int numSumTerms
The number of terms in the weighted sum, or the number of Pauli products.
Definition: QuEST.h:166
long long int getQubitBitMask(int *qubits, int numQubits)
Definition: QuEST_common.c:44
void statevec_compactUnitary(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
Represents a diagonal complex operator on the full Hilbert state of a Qureg.
Definition: QuEST.h:178
Represents a weighted sum of pauli products.
Definition: QuEST.h:158
int getNumQubits(Qureg qureg)
Get the number of qubits in a qureg object.
Definition: QuEST.c:668
void validateNumElems(DiagonalOp op, long long int startInd, long long int numElems, const char *caller)
void densmatr_mixDamping(Qureg qureg, int targetQubit, qreal damping)
Complex calcInnerProduct(Qureg bra, Qureg ket)
Computes the inner product of two equal-size state vectors, given by.
Definition: QuEST.c:910
void statevec_multiRotatePauli(Qureg qureg, int *targetQubits, enum pauliOpType *targetPaulis, int numTargets, qreal angle, int applyConj)
applyConj=1 will apply conjugate operation, else applyConj=0
Definition: QuEST_common.c:411
void statevec_multiControlledPhaseFlip(Qureg qureg, int *controlQubits, int numControlQubits)
Definition: QuEST_cpu.c:3331
void destroyQureg(Qureg qureg, QuESTEnv env)
Deallocate a Qureg object representing a set of qubits.
Definition: QuEST.c:77
void validateOneQubitKrausMap(Qureg qureg, ComplexMatrix2 *ops, int numOps, const char *caller)
void qasm_recordComment(Qureg qureg, char *comment,...)
Definition: QuEST_qasm.c:120
Complex densmatr_calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
void controlledMultiQubitUnitary(Qureg qureg, int ctrl, int *targs, int numTargs, ComplexMatrixN u)
Apply a general controlled multi-qubit unitary (including a global phase factor).
Definition: QuEST.c:314
void densmatr_mixTwoQubitDepolarising(Qureg qureg, int qubit1, int qubit2, qreal depolLevel)
void statevec_multiQubitUnitary(Qureg qureg, int *targets, int numTargets, ComplexMatrixN u)
Definition: QuEST_common.c:529
void mixTwoQubitKrausMap(Qureg qureg, int target1, int target2, ComplexMatrix4 *ops, int numOps)
Apply a general two-qubit Kraus map to a density matrix, as specified by at most sixteen Kraus operat...
Definition: QuEST.c:1075
qreal ** real
Definition: QuEST.h:139
void initDebugState(Qureg qureg)
Initialises qureg to be in the un-normalised, non-physical state with with n-th complex amplitude (2n...
Definition: QuEST.c:1308
void initStateFromSingleFile(Qureg *qureg, char filename[200], QuESTEnv env)
Initialises the wavefunction amplitudes according to those specified in a file.
Definition: QuEST.c:1312
void statevec_controlledRotateAroundAxis(Qureg qureg, int controlQubit, int targetQubit, qreal angle, Vector axis)
Definition: QuEST_common.c:327
void startRecordingQASM(Qureg qureg)
Enable QASM recording.
Definition: QuEST.c:87
void statevec_applyPauliSum(Qureg inQureg, enum pauliOpType *allCodes, qreal *termCoeffs, int numSumTerms, Qureg outQureg)
Definition: QuEST_common.c:494
void applyPauliSum(Qureg inQureg, enum pauliOpType *allPauliCodes, qreal *termCoeffs, int numSumTerms, Qureg outQureg)
Modifies outQureg to be the result of applying the weighted sum of Pauli products (a Hermitian but no...
Definition: QuEST.c:808
void applyMatrixN(Qureg qureg, int *targs, int numTargs, ComplexMatrixN u)
Apply a general N-by-N matrix, which may be non-unitary, on any number of target qubits.
Definition: QuEST.c:863
Complex statevec_calcInnerProduct(Qureg bra, Qureg ket)
Terrible code which unnecessarily individually computes and sums the real and imaginary components of...
void statevec_rotateAroundAxisConj(Qureg qureg, int rotQubit, qreal angle, Vector axis)
Definition: QuEST_common.c:318
void statevec_sGateConj(Qureg qureg, int targetQubit)
Definition: QuEST_common.c:279
void compactUnitary(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
Apply a single-qubit unitary parameterised by two given complex scalars.
Definition: QuEST.c:405
void initClassicalState(Qureg qureg, long long int stateInd)
Initialise a set of qubits to the classical state (also known as a "computational basis state") with...
Definition: QuEST.c:134
void validateOneQubitDepolProb(qreal prob, const char *caller)
int getQuEST_PREC(void)
Definition: QuEST.c:1339
void statevec_setWeightedQureg(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out)
Definition: QuEST_cpu.c:3619
void statevec_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal outcomeProb)
void validateTwoQubitKrausMap(Qureg qureg, ComplexMatrix4 *ops, int numOps, const char *caller)
void pauliY(Qureg qureg, int targetQubit)
Apply the single-qubit Pauli-Y (also known as the Y or sigma-Y) gate.
Definition: QuEST.c:444
long long int getControlFlipMask(int *controlQubits, int *controlState, int numControlQubits)
Definition: QuEST_common.c:54
void qasm_recordMultiControlledUnitary(Qureg qureg, ComplexMatrix2 u, int *controlQubits, int numControlQubits, int targetQubit)
additionally performs Rz on target to restore the global phase lost from u in QASM U(a,...
Definition: QuEST_qasm.c:341
qreal statevec_getImagAmp(Qureg qureg, long long int index)
Represents a system of qubits.
Definition: QuEST.h:203
void qasm_recordMultiControlledParamGate(Qureg qureg, TargetGate gate, int *controlQubits, int numControlQubits, int targetQubit, qreal param)
Definition: QuEST_qasm.c:324
void controlledPhaseShift(Qureg qureg, int idQubit1, int idQubit2, qreal angle)
Introduce a phase factor on state of qubits idQubit1 and idQubit2.
Definition: QuEST.c:499
void mixTwoQubitDephasing(Qureg qureg, int qubit1, int qubit2, qreal prob)
Mixes a density matrix qureg to induce two-qubit dephasing noise.
Definition: QuEST.c:1011
void controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
Apply a general controlled unitary (single control, single target), which can include a global phase ...
Definition: QuEST.c:361
int statevec_measureWithStats(Qureg qureg, int measureQubit, qreal *outcomeProb)
Definition: QuEST_common.c:361
void validateNumQubitsInMatrix(int numQubits, const char *caller)
void statevec_initDebugState(Qureg qureg)
Initialise the state vector of probability amplitudes to an (unphysical) state with each component of...
Definition: QuEST_cpu.c:1591
qreal ** imag
Definition: QuEST.h:140
void initDiagonalOp(DiagonalOp op, qreal *real, qreal *imag)
Updates the entire DiagonalOp op with the given elements, of which there must be 2^op....
Definition: QuEST.c:1286
void statevec_controlledNot(Qureg qureg, int controlQubit, int targetQubit)
void phaseShift(Qureg qureg, int targetQubit, qreal angle)
Shift the phase between and of a single qubit by a given angle.
Definition: QuEST.c:488
void validatePauliHamil(PauliHamil hamil, const char *caller)
void validateTwoQubitUnitaryMatrix(Qureg qureg, ComplexMatrix4 u, const char *caller)
void statevec_controlledRotateX(Qureg qureg, int controlQubit, int targetQubit, qreal angle)
Definition: QuEST_common.c:343
void validateDiagOpInit(DiagonalOp op, const char *caller)
void controlledNot(Qureg qureg, int controlQubit, int targetQubit)
Apply the controlled not (single control, single target) gate, also known as the c-X,...
Definition: QuEST.c:525
long long int getNumAmps(Qureg qureg)
Get the number of probability amplitudes in a qureg object, given by 2^numQubits.
Definition: QuEST.c:672
qreal densmatr_calcFidelity(Qureg qureg, Qureg pureState)
void statevec_initClassicalState(Qureg qureg, long long int stateInd)
Definition: QuEST_cpu.c:1470
void validateMatchingQuregPauliHamilDims(Qureg qureg, PauliHamil hamil, const char *caller)
int isDensityMatrix
Whether this instance is a density-state representation.
Definition: QuEST.h:206
void statevec_hadamard(Qureg qureg, int targetQubit)
void statevec_rotateY(Qureg qureg, int rotQubit, qreal angle)
Definition: QuEST_common.c:299
void validateMatchingQuregTypes(Qureg qureg1, Qureg qureg2, const char *caller)
int numQubits
Definition: QuEST.h:138
void validateAmpIndex(Qureg qureg, long long int ampInd, const char *caller)
void pauliX(Qureg qureg, int targetQubit)
Apply the single-qubit Pauli-X (also known as the X, sigma-X, NOT or bit-flip) gate.
Definition: QuEST.c:433
Qureg createCloneQureg(Qureg qureg, QuESTEnv env)
Create a new Qureg which is an exact clone of the passed qureg, which can be either a statevector or ...
Definition: QuEST.c:64
void statevec_multiControlledUnitary(Qureg qureg, long long int ctrlQubitsMask, long long int ctrlFlipMask, int targetQubit, ComplexMatrix2 u)
qreal calcExpecPauliSum(Qureg qureg, enum pauliOpType *allPauliCodes, qreal *termCoeffs, int numSumTerms, Qureg workspace)
Computes the expected value of a sum of products of Pauli operators.
Definition: QuEST.c:961
int statevec_initStateFromSingleFile(Qureg *qureg, char filename[200], QuESTEnv env)
Definition: QuEST_cpu.c:1625
void validateHamilFileParams(int numQubits, int numTerms, FILE *file, char *fn, const char *caller)
void validateOneQubitUnitaryMatrix(ComplexMatrix2 u, const char *caller)
void controlledTwoQubitUnitary(Qureg qureg, int controlQubit, int targetQubit1, int targetQubit2, ComplexMatrix4 u)
Apply a general controlled two-qubit unitary (including a global phase factor).
Definition: QuEST.c:270
void validateSecondQuregStateVec(Qureg qureg2, const char *caller)
int numQubits
The number of qubits for which this Hamiltonian is defined.
Definition: QuEST.h:168
void statevec_sqrtSwapGateConj(Qureg qureg, int qb1, int qb2)
Definition: QuEST_common.c:397
void statevec_controlledRotateZ(Qureg qureg, int controlQubit, int targetQubit, qreal angle)
Definition: QuEST_common.c:355
int numQubitsRepresented
The number of qubits represented in either the state-vector or density matrix.
Definition: QuEST.h:208
void statevec_pauliX(Qureg qureg, int targetQubit)
long long int numAmpsTotal
Total number of amplitudes, which are possibly distributed among machines.
Definition: QuEST.h:215
void validateMultiControlsMultiTargets(Qureg qureg, int *controlQubits, int numControlQubits, int *targetQubits, int numTargetQubits, const char *caller)
@ GATE_S
Definition: QuEST_qasm.h:25
@ GATE_SWAP
Definition: QuEST_qasm.h:33
qreal real
Definition: QuEST.h:105
void mixDensityMatrix(Qureg combineQureg, qreal otherProb, Qureg otherQureg)
Modifies combineQureg to become (1-prob)combineProb + prob otherQureg.
Definition: QuEST.c:772
Qureg createQureg(int numQubits, QuESTEnv env)
Create a Qureg object representing a set of qubits which will remain in a pure state.
Definition: QuEST.c:36
void destroyPauliHamil(PauliHamil h)
Destroy a PauliHamil instance, created with either createPauliHamil() or createPauliHamilFromFile().
Definition: QuEST.c:1163
void densmatr_mixPauli(Qureg qureg, int qubit, qreal probX, qreal probY, qreal probZ)
Definition: QuEST_common.c:676
void tGate(Qureg qureg, int targetQubit)
Apply the single-qubit T gate.
Definition: QuEST.c:477
void multiRotateZ(Qureg qureg, int *qubits, int numQubits, qreal angle)
Apply a multi-qubit Z rotation on a selected number of qubits.
Definition: QuEST.c:626
qreal imag
Definition: QuEST.h:106
void agnostic_setDiagonalOpElems(DiagonalOp op, long long int startInd, qreal *real, qreal *imag, long long int numElems)
Definition: QuEST_cpu.c:3842
void statevec_controlledPhaseFlip(Qureg qureg, int idQubit1, int idQubit2)
Definition: QuEST_cpu.c:3300
void validateMultiQubitKrausMap(Qureg qureg, int numTargs, ComplexMatrixN *ops, int numOps, const char *caller)
void validateMultiQubitMatrix(Qureg qureg, ComplexMatrixN u, int numTargs, const char *caller)
void validateNumQubitsInDiagOp(int numQubits, int numRanks, const char *caller)
@ GATE_SIGMA_Y
Definition: QuEST_qasm.h:22
void qasm_recordCompactUnitary(Qureg qureg, Complex alpha, Complex beta, int targetQubit)
Definition: QuEST_qasm.c:195
void mixPauli(Qureg qureg, int qubit, qreal probX, qreal probY, qreal probZ)
Mixes a density matrix qureg to induce general single-qubit Pauli noise.
Definition: QuEST.c:1054
void ensureIndsIncrease(int *ind1, int *ind2)
Definition: QuEST_common.c:64
qreal densmatr_calcHilbertSchmidtDistance(Qureg a, Qureg b)
int statevec_compareStates(Qureg mq1, Qureg mq2, qreal precision)
Definition: QuEST_cpu.c:1675
void qasm_recordControlledUnitary(Qureg qureg, ComplexMatrix2 u, int controlQubit, int targetQubit)
additionally performs Rz on target to restore the global phase lost from u in QASM U(a,...
Definition: QuEST_qasm.c:278
void hadamard(Qureg qureg, int targetQubit)
Apply the single-qubit Hadamard gate.
Definition: QuEST.c:177
Represents one complex number.
Definition: QuEST.h:103
void statevec_controlledRotateY(Qureg qureg, int controlQubit, int targetQubit, qreal angle)
Definition: QuEST_common.c:349
void applyMatrix2(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Apply a general 2-by-2 matrix, which may be non-unitary.
Definition: QuEST.c:844
void statevec_controlledRotateAroundAxisConj(Qureg qureg, int controlQubit, int targetQubit, qreal angle, Vector axis)
Definition: QuEST_common.c:334
void rotateZ(Qureg qureg, int targetQubit, qreal angle)
Rotate a single qubit by a given angle around the Z-axis of the Bloch-sphere (also known as a phase s...
Definition: QuEST.c:210
void multiControlledUnitary(Qureg qureg, int *controlQubits, int numControlQubits, int targetQubit, ComplexMatrix2 u)
Apply a general multiple-control single-target unitary, which can include a global phase factor.
Definition: QuEST.c:374
void initZeroState(Qureg qureg)
Initialise a set of qubits to the classical zero state .
Definition: QuEST.c:113
void validateTrotterParams(int order, int reps, const char *caller)
void qasm_setup(Qureg *qureg)
Definition: QuEST_qasm.c:60
void qasm_recordGate(Qureg qureg, TargetGate gate, int targetQubit)
Definition: QuEST_qasm.c:178
void validatePauliCodes(enum pauliOpType *pauliCodes, int numPauliCodes, const char *caller)
qreal calcDensityInnerProduct(Qureg rho1, Qureg rho2)
Computes the Hilbert-Schmidt scalar product (which is equivalent to the Frobenius inner product of ma...
Definition: QuEST.c:918
void statevec_createQureg(Qureg *qureg, int numQubits, QuESTEnv env)
Definition: QuEST_cpu.c:1279
qreal densmatr_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
Qureg createDensityQureg(int numQubits, QuESTEnv env)
Create a Qureg for qubits which are represented by a density matrix, and can be in mixed states.
Definition: QuEST.c:50
void validateUnitaryComplexPair(Complex alpha, Complex beta, const char *caller)
@ GATE_ROTATE_Y
Definition: QuEST_qasm.h:28
void statevec_rotateZ(Qureg qureg, int rotQubit, qreal angle)
Definition: QuEST_common.c:305
void applyTrotterCircuit(Qureg qureg, PauliHamil hamil, qreal time, int order, int reps)
Applies a trotterisation of unitary evolution to qureg.
Definition: QuEST.c:830
PauliHamil createPauliHamil(int numQubits, int numSumTerms)
Create a PauliHamil instance, which is a Hamiltonian expressed as a real-weighted sum of products of ...
Definition: QuEST.c:1147
void initPlusState(Qureg qureg)
Initialise a set of qubits to the plus state (and similarly for density matrices).
Definition: QuEST.c:125
qreal calcExpecPauliProd(Qureg qureg, int *targetQubits, enum pauliOpType *pauliCodes, int numTargets, Qureg workspace)
Computes the expected value of a product of Pauli operators.
Definition: QuEST.c:952
void validateOneQubitDampingProb(qreal prob, const char *caller)
void controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
Apply a controlled unitary (single control, single target) parameterised by two given complex scalars...
Definition: QuEST.c:418
void multiStateControlledUnitary(Qureg qureg, int *controlQubits, int *controlState, int numControlQubits, int targetQubit, ComplexMatrix2 u)
Apply a general multiple-control, conditioned on a specific bit sequence, single-target unitary,...
Definition: QuEST.c:389
void statevec_controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
qreal statevec_getRealAmp(Qureg qureg, long long int index)
void sqrtSwapGate(Qureg qureg, int qb1, int qb2)
Performs a sqrt SWAP gate between qubit1 and qubit2.
Definition: QuEST.c:613
void validateUniqueTargets(Qureg qureg, int qubit1, int qubit2, const char *caller)
void mixMultiQubitKrausMap(Qureg qureg, int *targets, int numTargets, ComplexMatrixN *ops, int numOps)
Apply a general N-qubit Kraus map to a density matrix, as specified by at most (2N)^2 Kraus operators...
Definition: QuEST.c:1085
void printRecordedQASM(Qureg qureg)
Print recorded QASM to stdout.
Definition: QuEST.c:99
void statevec_pauliYConj(Qureg qureg, int targetQubit)
void validateTwoQubitDephaseProb(qreal prob, const char *caller)
void reportStateToScreen(Qureg qureg, QuESTEnv env, int reportRank)
Print the current state vector of probability amplitudes for a set of qubits to standard out.
Definition: QuEST.c:1324
DiagonalOp createDiagonalOp(int numQubits, QuESTEnv env)
Creates a DiagonalOp representing a diagonal operator on the full Hilbert space of a Qureg.
Definition: QuEST.c:1267
int densmatr_measureWithStats(Qureg qureg, int measureQubit, qreal *outcomeProb)
Definition: QuEST_common.c:369
qreal statevec_calcTotalProb(Qureg qureg)
Represents a 2x2 matrix of complex numbers.
Definition: QuEST.h:114
void densmatr_initPureState(Qureg targetQureg, Qureg copyQureg)
void validateOneQubitPauliProbs(qreal probX, qreal probY, qreal probZ, const char *caller)
void controlledRotateY(Qureg qureg, int controlQubit, int targetQubit, qreal angle)
Applies a controlled rotation by a given angle around the Y-axis of the Bloch-sphere.
Definition: QuEST.c:233
void densmatr_applyDiagonalOp(Qureg qureg, DiagonalOp op)