QuEST.h
Go to the documentation of this file.
1 // Distributed under MIT licence. See https://github.com/QuEST-Kit/QuEST/blob/master/LICENCE.txt for details
2 
32 # ifndef QUEST_H
33 # define QUEST_H
34 
35 # include "QuEST_precision.h"
36 
37 // prevent C++ name mangling
38 #ifdef __cplusplus
39 extern "C" {
40 #endif
41 
42 
43 /*
44  * private structures
45  */
46 
47 // hide these from doxygen
49 
55 enum phaseGateType {SIGMA_Z=0, S_GATE=1, T_GATE=2};
56 
62 typedef struct {
63 
64  char* buffer; // generated QASM string
65  int bufferSize; // maximum number of chars before overflow
66  int bufferFill; // number of chars currently in buffer
67  int isLogging; // whether gates are being added to buffer
68 
69 } QASMLogger;
70 
77 typedef struct ComplexArray
78 {
79  qreal *real;
80  qreal *imag;
81 } ComplexArray;
82 
84 
85 
86 
87 /*
88  * public structures
89  */
90 
97 
103 typedef struct Complex
104 {
107 } Complex;
108 
114 typedef struct ComplexMatrix2
115 {
116  qreal real[2][2];
117  qreal imag[2][2];
119 
125 typedef struct ComplexMatrix4
126 {
127  qreal real[4][4];
128  qreal imag[4][4];
130 
136 typedef struct ComplexMatrixN
137 {
142 
148 typedef struct Vector
149 {
150  qreal x, y, z;
151 } Vector;
152 
158 typedef struct PauliHamil
159 {
169 } PauliHamil;
170 
178 typedef struct DiagonalOp
179 {
183  long long int numElemsPerChunk;
187  int chunkId;
193  ComplexArray deviceOperator;
194 } DiagonalOp;
195 
203 typedef struct Qureg
204 {
213  long long int numAmpsPerChunk;
215  long long int numAmpsTotal;
217  int chunkId;
220 
222  ComplexArray stateVec;
224  ComplexArray pairStateVec;
225 
227  ComplexArray deviceStateVec;
230 
232  QASMLogger* qasmLog;
233 
234 } Qureg;
235 
242 typedef struct QuESTEnv
243 {
244  int rank;
245  int numRanks;
246 } QuESTEnv;
247 
248 
249 
250 /*
251  * public functions
252  */
253 
268 Qureg createQureg(int numQubits, QuESTEnv env);
269 
284 Qureg createDensityQureg(int numQubits, QuESTEnv env);
285 
298 
308 void destroyQureg(Qureg qureg, QuESTEnv env);
309 
326 ComplexMatrixN createComplexMatrixN(int numQubits);
327 
339 
340 #ifndef __cplusplus
341 
362 void initComplexMatrixN(ComplexMatrixN m, qreal real[][1<<m.numQubits], qreal imag[][1<<m.numQubits]);
363 #endif
364 
383 PauliHamil createPauliHamil(int numQubits, int numSumTerms);
384 
391 void destroyPauliHamil(PauliHamil hamil);
392 
422 
439 void initPauliHamil(PauliHamil hamil, qreal* coeffs, enum pauliOpType* codes);
440 
483 DiagonalOp createDiagonalOp(int numQubits, QuESTEnv env);
484 
494 
509 void syncDiagonalOp(DiagonalOp op);
510 
527 void initDiagonalOp(DiagonalOp op, qreal* real, qreal* imag);
528 
565 void setDiagonalOpElems(DiagonalOp op, long long int startInd, qreal* real, qreal* imag, long long int numElems);
566 
577 void applyDiagonalOp(Qureg qureg, DiagonalOp op);
578 
602 
626 void reportState(Qureg qureg);
627 
635 void reportStateToScreen(Qureg qureg, QuESTEnv env, int reportRank);
636 
643 void reportQuregParams(Qureg qureg);
644 
662 void reportPauliHamil(PauliHamil hamil);
663 
669 int getNumQubits(Qureg qureg);
670 
676 long long int getNumAmps(Qureg qureg);
677 
686 void initBlankState(Qureg qureg);
687 
696 void initZeroState(Qureg qureg);
697 
712 void initPlusState(Qureg qureg);
713 
738 void initClassicalState(Qureg qureg, long long int stateInd);
739 
751 void initPureState(Qureg qureg, Qureg pure);
752 
762 void initDebugState(Qureg qureg);
763 
782 void initStateFromAmps(Qureg qureg, qreal* reals, qreal* imags);
783 
802 void setAmps(Qureg qureg, long long int startInd, qreal* reals, qreal* imags, long long int numAmps);
803 
814 void cloneQureg(Qureg targetQureg, Qureg copyQureg);
815 
848 void phaseShift(Qureg qureg, int targetQubit, qreal angle);
849 
891 void controlledPhaseShift(Qureg qureg, int idQubit1, int idQubit2, qreal angle);
892 
932 void multiControlledPhaseShift(Qureg qureg, int *controlQubits, int numControlQubits, qreal angle);
933 
969 void controlledPhaseFlip (Qureg qureg, int idQubit1, int idQubit2);
970 
1017 void multiControlledPhaseFlip(Qureg qureg, int *controlQubits, int numControlQubits);
1018 
1049 void sGate(Qureg qureg, int targetQubit);
1050 
1081 void tGate(Qureg qureg, int targetQubit);
1082 
1093 QuESTEnv createQuESTEnv(void);
1094 
1103 void destroyQuESTEnv(QuESTEnv env);
1104 
1111 void syncQuESTEnv(QuESTEnv env);
1112 
1121 int syncQuESTSuccess(int successCode);
1122 
1129 void reportQuESTEnv(QuESTEnv env);
1130 
1140 void getEnvironmentString(QuESTEnv env, Qureg qureg, char str[200]);
1141 
1163 void copyStateToGPU(Qureg qureg);
1164 
1186 void copyStateFromGPU(Qureg qureg);
1187 
1199 Complex getAmp(Qureg qureg, long long int index);
1200 
1212 qreal getRealAmp(Qureg qureg, long long int index);
1213 
1225 qreal getImagAmp(Qureg qureg, long long int index);
1226 
1238 qreal getProbAmp(Qureg qureg, long long int index);
1239 
1252 Complex getDensityAmp(Qureg qureg, long long int row, long long int col);
1253 
1279 qreal calcTotalProb(Qureg qureg);
1280 
1318 void compactUnitary(Qureg qureg, int targetQubit, Complex alpha, Complex beta);
1319 
1350 void unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u);
1351 
1383 void rotateX(Qureg qureg, int rotQubit, qreal angle);
1384 
1417 void rotateY(Qureg qureg, int rotQubit, qreal angle);
1418 
1451 void rotateZ(Qureg qureg, int rotQubit, qreal angle);
1452 
1471 void rotateAroundAxis(Qureg qureg, int rotQubit, qreal angle, Vector axis);
1472 
1473 
1505 void controlledRotateX(Qureg qureg, int controlQubit, int targetQubit, qreal angle);
1506 
1538 void controlledRotateY(Qureg qureg, int controlQubit, int targetQubit, qreal angle);
1539 
1571 void controlledRotateZ(Qureg qureg, int controlQubit, int targetQubit, qreal angle);
1572 
1609 void controlledRotateAroundAxis(Qureg qureg, int controlQubit, int targetQubit, qreal angle, Vector axis);
1610 
1656 void controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta);
1657 
1701 void controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u);
1702 
1761 void multiControlledUnitary(Qureg qureg, int* controlQubits, int numControlQubits, int targetQubit, ComplexMatrix2 u);
1762 
1793 void pauliX(Qureg qureg, int targetQubit);
1794 
1826 void pauliY(Qureg qureg, int targetQubit);
1827 
1859 void pauliZ(Qureg qureg, int targetQubit);
1860 
1894 void hadamard(Qureg qureg, int targetQubit);
1895 
1936 void controlledNot(Qureg qureg, int controlQubit, int targetQubit);
1937 
1980 void controlledPauliY(Qureg qureg, int controlQubit, int targetQubit);
1981 
2005 qreal calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome);
2006 
2032 qreal collapseToOutcome(Qureg qureg, int measureQubit, int outcome);
2033 
2047 int measure(Qureg qureg, int measureQubit);
2048 
2064 int measureWithStats(Qureg qureg, int measureQubit, qreal *outcomeProb);
2065 
2089 
2138 
2152 void seedQuESTDefault(void);
2153 
2169 void seedQuEST(unsigned long int *seedArray, int numSeeds);
2170 
2179 void startRecordingQASM(Qureg qureg);
2180 
2190 void stopRecordingQASM(Qureg qureg);
2191 
2198 void clearRecordedQASM(Qureg qureg);
2199 
2207 void printRecordedQASM(Qureg qureg);
2208 
2217 void writeRecordedQASMToFile(Qureg qureg, char* filename);
2218 
2240 void mixDephasing(Qureg qureg, int targetQubit, qreal prob);
2241 
2269 void mixTwoQubitDephasing(Qureg qureg, int qubit1, int qubit2, qreal prob);
2270 
2304 void mixDepolarising(Qureg qureg, int targetQubit, qreal prob);
2305 
2335 void mixDamping(Qureg qureg, int targetQubit, qreal prob);
2336 
2399 void mixTwoQubitDepolarising(Qureg qureg, int qubit1, int qubit2, qreal prob);
2400 
2434 void mixPauli(Qureg qureg, int targetQubit, qreal probX, qreal probY, qreal probZ);
2435 
2449 void mixDensityMatrix(Qureg combineQureg, qreal prob, Qureg otherQureg);
2450 
2472 qreal calcPurity(Qureg qureg);
2473 
2500 qreal calcFidelity(Qureg qureg, Qureg pureState);
2501 
2543 void swapGate(Qureg qureg, int qubit1, int qubit2);
2544 
2545 
2590 void sqrtSwapGate(Qureg qureg, int qb1, int qb2);
2591 
2644  Qureg qureg, int* controlQubits, int* controlState, int numControlQubits,
2645  int targetQubit, ComplexMatrix2 u
2646 );
2647 
2671 void multiRotateZ(Qureg qureg, int* qubits, int numQubits, qreal angle);
2672 
2719 void multiRotatePauli(Qureg qureg, int* targetQubits, enum pauliOpType* targetPaulis, int numTargets, qreal angle);
2720 
2765 qreal calcExpecPauliProd(Qureg qureg, int* targetQubits, enum pauliOpType* pauliCodes, int numTargets, Qureg workspace);
2766 
2815 qreal calcExpecPauliSum(Qureg qureg, enum pauliOpType* allPauliCodes, qreal* termCoeffs, int numSumTerms, Qureg workspace);
2816 
2850 qreal calcExpecPauliHamil(Qureg qureg, PauliHamil hamil, Qureg workspace);
2851 
2913 void twoQubitUnitary(Qureg qureg, int targetQubit1, int targetQubit2, ComplexMatrix4 u);
2914 
2977 void controlledTwoQubitUnitary(Qureg qureg, int controlQubit, int targetQubit1, int targetQubit2, ComplexMatrix4 u);
2978 
3052 void multiControlledTwoQubitUnitary(Qureg qureg, int* controlQubits, int numControlQubits, int targetQubit1, int targetQubit2, ComplexMatrix4 u);
3053 
3129 void multiQubitUnitary(Qureg qureg, int* targs, int numTargs, ComplexMatrixN u);
3130 
3200 void controlledMultiQubitUnitary(Qureg qureg, int ctrl, int* targs, int numTargs, ComplexMatrixN u);
3201 
3277 void multiControlledMultiQubitUnitary(Qureg qureg, int* ctrls, int numCtrls, int* targs, int numTargs, ComplexMatrixN u);
3278 
3312 void mixKrausMap(Qureg qureg, int target, ComplexMatrix2 *ops, int numOps);
3313 
3346 void mixTwoQubitKrausMap(Qureg qureg, int target1, int target2, ComplexMatrix4 *ops, int numOps);
3347 
3390 void mixMultiQubitKrausMap(Qureg qureg, int* targets, int numTargets, ComplexMatrixN* ops, int numOps);
3391 
3419 
3443 void setWeightedQureg(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out);
3444 
3498 void applyPauliSum(Qureg inQureg, enum pauliOpType* allPauliCodes, qreal* termCoeffs, int numSumTerms, Qureg outQureg);
3499 
3536 void applyPauliHamil(Qureg inQureg, PauliHamil hamil, Qureg outQureg);
3537 
3595 void applyTrotterCircuit(Qureg qureg, PauliHamil hamil, qreal time, int order, int reps);
3596 
3612 void applyMatrix2(Qureg qureg, int targetQubit, ComplexMatrix2 u);
3613 
3660 void applyMatrix4(Qureg qureg, int targetQubit1, int targetQubit2, ComplexMatrix4 u);
3661 
3723 void applyMatrixN(Qureg qureg, int* targs, int numTargs, ComplexMatrixN u);
3724 
3776 void applyMultiControlledMatrixN(Qureg qureg, int* ctrls, int numCtrls, int* targs, int numTargs, ComplexMatrixN u);
3777 
3816 void invalidQuESTInputError(const char* errMsg, const char* errFunc);
3817 
3818 #ifndef __cplusplus
3819  // hide this function from doxygen
3821 
3860  int numQubits, qreal re[][1<<numQubits], qreal im[][1<<numQubits],
3861  qreal** reStorage, qreal** imStorage);
3862 #endif
3863 
3865 // hide this function from doxygen
3867 #define UNPACK_ARR(...) __VA_ARGS__
3868 
3870 #ifndef __cplusplus
3871 
3910 #define getStaticComplexMatrixN(numQubits, re, im) \
3911  bindArraysToStackComplexMatrixN( \
3912  numQubits, \
3913  (qreal[1<<numQubits][1<<numQubits]) UNPACK_ARR re, \
3914  (qreal[1<<numQubits][1<<numQubits]) UNPACK_ARR im, \
3915  (double*[1<<numQubits]) {NULL}, (double*[1<<numQubits]) {NULL} \
3916  )
3917 #endif
3918 
3919 
3920 // end prevention of C++ name mangling
3921 #ifdef __cplusplus
3922 }
3923 #endif
3924 
3925 #endif // QUEST_H
3926 
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 mixDephasing(Qureg qureg, int targetQubit, qreal prob)
Mixes a density matrix qureg to induce single-qubit dephasing noise.
Definition: QuEST.c:1001
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 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 destroyQuESTEnv(QuESTEnv env)
Destroy the QuEST environment.
pauliOpType
Codes for specifying Pauli operators.
Definition: QuEST.h:96
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 reportQuregParams(Qureg qureg)
Report metainformation about a set of qubits: number of qubits, number of probability amplitudes.
Definition: QuEST_common.c:234
void copyStateFromGPU(Qureg qureg)
In GPU mode, this copies the state-vector (or density matrix) from GPU memory (qureg....
Definition: QuEST_cpu.c:39
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
qreal real[4][4]
Definition: QuEST.h:127
void syncQuESTEnv(QuESTEnv env)
Guarantees that all code up to the given point has been executed on all nodes (if running in distribu...
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
@ PAULI_Z
Definition: QuEST.h:96
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
int rank
Definition: QuEST.h:244
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 matr)
Destroy a ComplexMatrixN instance created with createComplexMatrixN()
Definition: QuEST.c:1120
int numChunks
The number of nodes between which the elements of this operator are split.
Definition: QuEST.h:185
void reportPauliHamil(PauliHamil hamil)
Print the PauliHamil to screen.
Definition: QuEST.c:1328
ComplexArray pairStateVec
Temporary storage for a chunk of the state vector received from another process in the MPI version.
Definition: QuEST.h:224
@ 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
void getEnvironmentString(QuESTEnv env, Qureg qureg, char str[200])
Sets str to a string containing the number of qubits in qureg, and the hardware facilities used (e....
Definition: QuEST_gpu.cu:447
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 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 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 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
qreal z
Definition: QuEST.h:150
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
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
int numChunks
Number of chunks the state vector is broken up into – the number of MPI processes used.
Definition: QuEST.h:219
qreal calcPurity(Qureg qureg)
Calculates the purity of a density matrix, by the trace of the density matrix squared.
Definition: QuEST.c:936
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
ComplexArray deviceOperator
A copy of the elements stored persistently on the GPU.
Definition: QuEST.h:193
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
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 unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Apply a general single-qubit unitary (including a global phase factor).
Definition: QuEST.c:349
int chunkId
The position of the chunk of the operator held by this process in the full operator.
Definition: QuEST.h:187
ComplexArray deviceStateVec
Storage for wavefunction amplitudes in the GPU version.
Definition: QuEST.h:227
void seedQuEST(unsigned long int *seedArray, int numSeeds)
Seed the Mersenne Twister used for random number generation in the QuEST environment with a user defi...
Definition: QuEST_common.c:209
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
Complex calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
Computes the expected value of the diagonal operator op for state qureg.
Definition: QuEST.c:979
void sGate(Qureg qureg, int targetQubit)
Apply the single-qubit S gate.
Definition: QuEST.c:466
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 rotQubit, qreal angle)
Rotate a single qubit by a given angle around the Y-axis of the Bloch-sphere.
Definition: QuEST.c:199
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
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
ComplexMatrixN bindArraysToStackComplexMatrixN(int numQubits, qreal re[][1<< numQubits], qreal im[][1<< numQubits], qreal **reStorage, qreal **imStorage)
Definition: QuEST_common.c:607
Represents a general 2^N by 2^N matrix of complex numbers.
Definition: QuEST.h:136
#define qreal
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 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
@ PAULI_X
Definition: QuEST.h:96
Complex getAmp(Qureg qureg, long long int index)
Get the complex amplitude at a given index in the state vector.
Definition: QuEST.c:699
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
void reportState(Qureg qureg)
Print the current state vector of probability amplitudes for a set of qubits to file.
Definition: QuEST_common.c:216
void initComplexMatrixN(ComplexMatrixN m, qreal real[][1<< m.numQubits], qreal imag[][1<< m.numQubits])
Initialises a ComplexMatrixN instance to have the passed real and imag values.
Definition: QuEST.c:1136
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 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
void rotateX(Qureg qureg, int rotQubit, qreal angle)
Rotate a single qubit by a given angle around the X-axis of the Bloch-sphere.
Definition: QuEST.c:188
void invalidQuESTInputError(const char *errMsg, const char *errFunc)
An internal function called when invalid arguments are passed to a QuEST API call,...
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
int chunkId
The position of the chunk of the state vector held by this process in the full state vector.
Definition: QuEST.h:217
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
qreal y
Definition: QuEST.h:150
qreal imag[2][2]
Definition: QuEST.h:117
qreal * imag
The imaginary values of the 2^numQubits complex elements.
Definition: QuEST.h:191
qreal x
Definition: QuEST.h:150
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
long long int numAmpsPerChunk
Number of probability amplitudes held in stateVec by this process In the non-MPI version,...
Definition: QuEST.h:213
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 swapGate(Qureg qureg, int qubit1, int qubit2)
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 initStateFromAmps(Qureg qureg, qreal *reals, qreal *imags)
Initialise qureg by specifying the complete statevector.
Definition: QuEST.c:157
void copyStateToGPU(Qureg qureg)
In GPU mode, this copies the state-vector (or density matrix) from RAM (qureg.stateVec) to VRAM / GPU...
Definition: QuEST_cpu.c:36
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
qreal imag[4][4]
Definition: QuEST.h:128
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
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
Represents a diagonal complex operator on the full Hilbert state of a Qureg.
Definition: QuEST.h:178
@ PAULI_Y
Definition: QuEST.h:96
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
Complex calcInnerProduct(Qureg bra, Qureg ket)
Computes the inner product of two equal-size state vectors, given by.
Definition: QuEST.c:910
void destroyQureg(Qureg qureg, QuESTEnv env)
Deallocate a Qureg object representing a set of qubits.
Definition: QuEST.c:77
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 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
QASMLogger * qasmLog
Storage for generated QASM output.
Definition: QuEST.h:232
qreal * secondLevelReduction
Definition: QuEST.h:229
void startRecordingQASM(Qureg qureg)
Enable QASM recording.
Definition: QuEST.c:87
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
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 pauliY(Qureg qureg, int targetQubit)
Apply the single-qubit Pauli-Y (also known as the Y or sigma-Y) gate.
Definition: QuEST.c:444
Represents a system of qubits.
Definition: QuEST.h:203
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
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 phaseShift(Qureg qureg, int targetQubit, qreal angle)
Shift the phase between and of a single qubit by a given angle.
Definition: QuEST.c:488
ComplexArray stateVec
Computational state amplitudes - a subset thereof in the MPI version.
Definition: QuEST.h:222
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
qreal real[2][2]
Definition: QuEST.h:116
long long int getNumAmps(Qureg qureg)
Get the number of probability amplitudes in a qureg object, given by 2^numQubits.
Definition: QuEST.c:672
void seedQuESTDefault(void)
Seed the Mersenne Twister used for random number generation in the QuEST environment with an example ...
long long int numElemsPerChunk
The number of the 2^numQubits amplitudes stored on each distributed node.
Definition: QuEST.h:183
int isDensityMatrix
Whether this instance is a density-state representation.
Definition: QuEST.h:206
void reportQuESTEnv(QuESTEnv env)
Report information about the QuEST environment.
int numQubits
Definition: QuEST.h:138
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
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
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
int numQubits
The number of qubits for which this Hamiltonian is defined.
Definition: QuEST.h:168
int numQubitsRepresented
The number of qubits represented in either the state-vector or density matrix.
Definition: QuEST.h:208
long long int numAmpsTotal
Total number of amplitudes, which are possibly distributed among machines.
Definition: QuEST.h:215
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
qreal real
Definition: QuEST.h:105
void mixDensityMatrix(Qureg combineQureg, qreal prob, 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 hamil)
Destroy a PauliHamil instance, created with either createPauliHamil() or createPauliHamilFromFile().
Definition: QuEST.c:1163
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 mixPauli(Qureg qureg, int targetQubit, qreal probX, qreal probY, qreal probZ)
Mixes a density matrix qureg to induce general single-qubit Pauli noise.
Definition: QuEST.c:1054
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 applyMatrix2(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Apply a general 2-by-2 matrix, which may be non-unitary.
Definition: QuEST.c:844
QuESTEnv createQuESTEnv(void)
Create the QuEST execution environment.
void rotateZ(Qureg qureg, int rotQubit, 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
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
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
qreal * firstLevelReduction
Storage for reduction of probabilities on GPU.
Definition: QuEST.h:229
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 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 sqrtSwapGate(Qureg qureg, int qb1, int qb2)
Performs a sqrt SWAP gate between qubit1 and qubit2.
Definition: QuEST.c:613
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 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
Represents a 2x2 matrix of complex numbers.
Definition: QuEST.h:114
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