QuEST_gpu.cu
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
82 __forceinline__ __device__ int extractBit (const int locationOfBitFromRight, const long long int theEncodedNumber) {
99 __forceinline__ __device__ long long int insertZeroBit(const long long int number, const int index) {
106 __forceinline__ __device__ long long int insertTwoZeroBits(const long long int number, const int bit1, const int bit2) {
112 __forceinline__ __device__ long long int insertZeroBits(long long int number, int* inds, const int numInds) {
153 void statevec_setAmps(Qureg qureg, long long int startInd, qreal* reals, qreal* imags, long long int numAmps) {
216 __global__ void densmatr_initPlusStateKernel(long long int stateVecSize, qreal probFactor, qreal *stateVecReal, qreal *stateVecImag){
307 cudaMalloc(&(qureg->deviceStateVec.real), qureg->numAmpsPerChunk*sizeof(*(qureg->deviceStateVec.real)));
308 cudaMalloc(&(qureg->deviceStateVec.imag), qureg->numAmpsPerChunk*sizeof(*(qureg->deviceStateVec.imag)));
309 cudaMalloc(&(qureg->firstLevelReduction), ceil(qureg->numAmpsPerChunk/(qreal)REDUCE_SHARED_SIZE)*sizeof(qreal));
310 cudaMalloc(&(qureg->secondLevelReduction), ceil(qureg->numAmpsPerChunk/(qreal)(REDUCE_SHARED_SIZE*REDUCE_SHARED_SIZE))*
492 printf(REAL_STRING_FORMAT ", " REAL_STRING_FORMAT "\n", qureg.stateVec.real[index], qureg.stateVec.imag[index]);
515 __global__ void statevec_initBlankStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag){
536 __global__ void statevec_initZeroStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag){
563 __global__ void statevec_initPlusStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag){
585 __global__ void statevec_initClassicalStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag, long long int stateInd){
611 __global__ void statevec_initDebugStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag){
632 __global__ void statevec_initStateOfSingleQubitKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag, int qubitId, int outcome){
655 statevec_initStateOfSingleQubitKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg->numAmpsPerChunk, qureg->deviceStateVec.real, qureg->deviceStateVec.imag, qubitId, outcome);
721 __global__ void statevec_compactUnitaryKernel (Qureg qureg, int rotQubit, Complex alpha, Complex beta){
781 statevec_compactUnitaryKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, targetQubit, alpha, beta);
784 __global__ void statevec_controlledCompactUnitaryKernel (Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta){
843 void statevec_controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
848 statevec_controlledCompactUnitaryKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, controlQubit, targetQubit, alpha, beta);
909 statevec_unitaryKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, targetQubit, argifyMatrix2(u));
914 qreal* uRe, qreal* uIm, long long int* ampInds, qreal* reAmps, qreal* imAmps, long long int numTargAmps)
919 long long int numTasks = qureg.numAmpsPerChunk >> numTargs; // kernel called on every 1 in 2^numTargs amplitudes
971 void statevec_multiControlledMultiQubitUnitary(Qureg qureg, long long int ctrlMask, int* targs, int numTargs, ComplexMatrixN u)
994 // allocate device space for global u.real and u.imag (flatten by concatenating rows) and populate
1028 __global__ void statevec_multiControlledTwoQubitUnitaryKernel(Qureg qureg, long long int ctrlMask, int q1, int q2, ArgMatrix4 u){
1032 long long int numTasks = qureg.numAmpsPerChunk >> 2; // kernel called on every 1 in 4 amplitudes
1104 void statevec_multiControlledTwoQubitUnitary(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u)
1107 int CUDABlocks = ceil((qreal)(qureg.numAmpsPerChunk>>2)/threadsPerCUDABlock); // one kernel eval for every 4 amplitudes
1108 statevec_multiControlledTwoQubitUnitaryKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, ctrlMask, q1, q2, argifyMatrix4(u));
1111 __global__ void statevec_controlledUnitaryKernel(Qureg qureg, int controlQubit, int targetQubit, ArgMatrix2 u){
1169 void statevec_controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
1174 statevec_controlledUnitaryKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, controlQubit, targetQubit, argifyMatrix2(u));
1341 __global__ void statevec_controlledPauliYKernel(Qureg qureg, int controlQubit, int targetQubit, int conjFac)
1383 statevec_controlledPauliYKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, controlQubit, targetQubit, conjFactor);
1392 statevec_controlledPauliYKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, controlQubit, targetQubit, conjFactor);
1395 __global__ void statevec_phaseShiftByTermKernel(Qureg qureg, int targetQubit, qreal cosAngle, qreal sinAngle) {
1431 statevec_phaseShiftByTermKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, targetQubit, cosAngle, sinAngle);
1434 __global__ void statevec_controlledPhaseShiftKernel(Qureg qureg, int idQubit1, int idQubit2, qreal cosAngle, qreal sinAngle)
1467 statevec_controlledPhaseShiftKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, idQubit1, idQubit2, cosAngle, sinAngle);
1470 __global__ void statevec_multiControlledPhaseShiftKernel(Qureg qureg, long long int mask, qreal cosAngle, qreal sinAngle) {
1490 void statevec_multiControlledPhaseShift(Qureg qureg, int *controlQubits, int numControlQubits, qreal angle)
1500 statevec_multiControlledPhaseShiftKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, mask, cosAngle, sinAngle);
1503 __global__ void statevec_multiRotateZKernel(Qureg qureg, long long int mask, qreal cosAngle, qreal sinAngle) {
1528 statevec_multiRotateZKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, mask, cosAngle, sinAngle);
1556 point operation overhead. For more details see https://en.wikipedia.org/wiki/Kahan_summation_algorithm */
1612 statevec_controlledPhaseFlipKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, idQubit1, idQubit2);
1647 long long int numTasks = qureg.numAmpsPerChunk >> 2; // each iteration updates 2 amps and skips 2 amps
1777 statevec_controlledNotKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, controlQubit, targetQubit);
2211 __global__ void densmatr_calcFidelityKernel(Qureg dens, Qureg vec, long long int dim, qreal* reducedArray) {
2347 // spawn threads to sum the probs in each block (store reduction temp values in a's reduction array)
2349 densmatr_calcHilbertSchmidtDistanceSquaredKernel<<<numCUDABlocks, valuesPerCUDABlock, sharedMemSize>>>(
2375 __global__ void densmatr_calcPurityKernel(qreal* vecReal, qreal* vecImag, long long int numAmpsToSum, qreal *reducedArray) {
2439 cudaMemcpy(&traceDensSquared, qureg.firstLevelReduction, sizeof(qreal), cudaMemcpyDeviceToHost);
2443 __global__ void statevec_collapseToKnownProbOutcomeKernel(Qureg qureg, int measureQubit, int outcome, qreal totalProbability)
2500 void statevec_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal outcomeProb)
2505 statevec_collapseToKnownProbOutcomeKernel<<<CUDABlocks, threadsPerCUDABlock>>>(qureg, measureQubit, outcome, outcomeProb);
2535 void densmatr_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal outcomeProb) {
2564 __global__ void densmatr_mixDensityMatrixKernel(Qureg combineQureg, qreal otherProb, Qureg otherQureg, long long int numAmpsToVisit) {
2645 qreal fac, qreal* vecReal, qreal *vecImag, long long int numBackgroundStates, long long int numAmpsToVisit,
2646 long long int part1, long long int part2, long long int part3, long long int part4, long long int part5,
2657 long long int shift = rowBit2*((meta>>3)%2) + rowBit1*((meta>>2)%2) + colBit2*((meta>>1)%2) + colBit1*(meta%2);
2661 (scanInd&part1) + ((scanInd&part2)<<1) + ((scanInd&part3)<<2) + ((scanInd&part4)<<3) + ((scanInd&part5)<<4));
2667 // @TODO is separating these 12 amplitudes really faster than letting every 16th base modify 12 elems?
2700 dephFac, qureg.deviceStateVec.real, qureg.deviceStateVec.imag, numBackgroundStates, numAmpsToVisit,
2816 long long int ind00 = (scanInd&part1) + ((scanInd&part2)<<1) + ((scanInd&part3)<<2) + ((scanInd&part4)<<3) + ((scanInd&part5)<<4);
2873 __global__ void statevec_setWeightedQuregKernel(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out) {
2901 vecReOut[index] = (facReOut*reOut - facImOut*imOut) + (facRe1*re1 - facIm1*im1) + (facRe2*re2 - facIm2*im2);
2902 vecImOut[index] = (facReOut*imOut + facImOut*reOut) + (facRe1*im1 + facIm1*re1) + (facRe2*im2 + facIm2*re2);
2905 void statevec_setWeightedQureg(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out) {
3233 void agnostic_setDiagonalOpElems(DiagonalOp op, long long int startInd, qreal* real, qreal* imag, long long int numElems) {
void densmatr_mixDamping(Qureg qureg, int targetQubit, qreal damping)
Definition: QuEST_gpu.cu:2778
__global__ void statevec_multiControlledUnitaryKernel(Qureg qureg, long long int ctrlQubitsMask, long long int ctrlFlipMask, int targetQubit, ArgMatrix2 u)
Definition: QuEST_gpu.cu:1177
__global__ void statevec_initPlusStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag)
Definition: QuEST_gpu.cu:563
__global__ void statevec_setWeightedQuregKernel(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out)
Definition: QuEST_gpu.cu:2873
__global__ void densmatr_collapseToKnownProbOutcomeKernel(qreal outcomeProb, qreal *vecReal, qreal *vecImag, long long int numBasesToVisit, long long int part1, long long int part2, long long int part3, long long int rowBit, long long int colBit, long long int desired, long long int undesired)
Maps thread ID to a |..0..><..0..| state and then locates |0><1|, |1><0| and |1><1|.
Definition: QuEST_gpu.cu:2509
void statevec_swapQubitAmps(Qureg qureg, int qb1, int qb2)
Definition: QuEST_gpu.cu:1668
void copyStateFromGPU(Qureg qureg)
In GPU mode, this copies the state-vector (or density matrix) from GPU memory (qureg....
Definition: QuEST_gpu.cu:461
__global__ void statevec_controlledUnitaryKernel(Qureg qureg, int controlQubit, int targetQubit, ArgMatrix2 u)
Definition: QuEST_gpu.cu:1111
void statevec_multiRotateZ(Qureg qureg, long long int mask, qreal angle)
Definition: QuEST_gpu.cu:1520
void syncQuESTEnv(QuESTEnv env)
Guarantees that all code up to the given point has been executed on all nodes (if running in distribu...
Definition: QuEST_gpu.cu:423
__global__ void densmatr_initPureStateKernel(long long int numPureAmps, qreal *targetVecReal, qreal *targetVecImag, qreal *copyVecReal, qreal *copyVecImag)
Definition: QuEST_gpu.cu:186
void statevec_setAmps(Qureg qureg, long long int startInd, qreal *reals, qreal *imags, long long int numAmps)
Definition: QuEST_gpu.cu:153
__global__ void densmatr_mixTwoQubitDepolarisingKernel(qreal depolLevel, qreal *vecReal, qreal *vecImag, long long int numAmpsToVisit, long long int part1, long long int part2, long long int part3, long long int part4, long long int part5, long long int rowCol1, long long int rowCol2)
Called once for every 16 amplitudes.
Definition: QuEST_gpu.cu:2806
void statevec_multiControlledPhaseShift(Qureg qureg, int *controlQubits, int numControlQubits, qreal angle)
Definition: QuEST_gpu.cu:1490
__global__ void copySharedReduceBlock(qreal *arrayIn, qreal *reducedArray, int length)
Definition: QuEST_gpu.cu:1806
__global__ void statevec_multiControlledMultiQubitUnitaryKernel(Qureg qureg, long long int ctrlMask, int *targs, int numTargs, qreal *uRe, qreal *uIm, long long int *ampInds, qreal *reAmps, qreal *imAmps, long long int numTargAmps)
Definition: QuEST_gpu.cu:912
__global__ void densmatr_initClassicalStateKernel(long long int densityNumElems, qreal *densityReal, qreal *densityImag, long long int densityInd)
Definition: QuEST_gpu.cu:239
int statevec_initStateFromSingleFile(Qureg *qureg, char filename[200], QuESTEnv env)
Definition: QuEST_gpu.cu:659
int numChunks
The number of nodes between which the elements of this operator are split.
Definition: QuEST.h:185
ComplexArray pairStateVec
Temporary storage for a chunk of the state vector received from another process in the MPI version.
Definition: QuEST.h:224
qreal statevec_findProbabilityOfZero(Qureg qureg, int measureQubit)
Definition: QuEST_gpu.cu:1967
void statevec_controlledCompactUnitary(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_gpu.cu:843
__global__ void densmatr_calcExpecDiagonalOpKernel(int getRealComp, qreal *matReal, qreal *matImag, qreal *opReal, qreal *opImag, int numQubits, long long int numTermsToSum, qreal *reducedArray)
Definition: QuEST_gpu.cu:3097
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
__global__ void statevec_controlledPhaseShiftKernel(Qureg qureg, int idQubit1, int idQubit2, qreal cosAngle, qreal sinAngle)
Definition: QuEST_gpu.cu:1434
__global__ void statevec_multiRotateZKernel(Qureg qureg, long long int mask, qreal cosAngle, qreal sinAngle)
Definition: QuEST_gpu.cu:1503
int numChunks
Number of chunks the state vector is broken up into – the number of MPI processes used.
Definition: QuEST.h:219
__global__ void densmatr_applyDiagonalOpKernel(Qureg qureg, DiagonalOp op)
Definition: QuEST_gpu.cu:2947
void getQuESTDefaultSeedKey(unsigned long int *key)
Definition: QuEST_common.c:182
ComplexArray deviceOperator
A copy of the elements stored persistently on the GPU.
Definition: QuEST.h:193
void statevec_multiControlledTwoQubitUnitary(Qureg qureg, long long int ctrlMask, int q1, int q2, ComplexMatrix4 u)
This calls swapQubitAmps only when it would involve a distributed communication; if the qubit chunks ...
Definition: QuEST_gpu.cu:1104
__global__ void densmatr_calcHilbertSchmidtDistanceSquaredKernel(qreal *aRe, qreal *aIm, qreal *bRe, qreal *bIm, long long int numAmpsToSum, qreal *reducedArray)
Definition: QuEST_gpu.cu:2299
void statevec_controlledPhaseFlip(Qureg qureg, int idQubit1, int idQubit2)
Definition: QuEST_gpu.cu:1607
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
qreal statevec_getRealAmp(Qureg qureg, long long int index)
Definition: QuEST_gpu.cu:501
void statevec_createQureg(Qureg *qureg, int numQubits, QuESTEnv env)
Definition: QuEST_gpu.cu:275
qreal densmatr_calcInnerProduct(Qureg a, Qureg b)
Definition: QuEST_gpu.cu:2043
__global__ void densmatr_mixDampingKernel(qreal damping, qreal *vecReal, qreal *vecImag, long long int numAmpsToVisit, long long int part1, long long int part2, long long int part3, long long int bothBits)
Works like mixDephasing but modifies every other element, and elements are averaged in pairs.
Definition: QuEST_gpu.cu:2731
__global__ void statevec_calcExpecDiagonalOpKernel(int getRealComp, qreal *vecReal, qreal *vecImag, qreal *opReal, qreal *opImag, long long int numTermsToSum, qreal *reducedArray)
computes either a real or imag term of |vec_i|^2 op_i
Definition: QuEST_gpu.cu:2979
void densmatr_oneQubitDegradeOffDiagonal(Qureg qureg, int targetQubit, qreal dephFac)
Definition: QuEST_gpu.cu:2609
__global__ void densmatr_mixTwoQubitDephasingKernel(qreal fac, qreal *vecReal, qreal *vecImag, long long int numBackgroundStates, long long int numAmpsToVisit, long long int part1, long long int part2, long long int part3, long long int part4, long long int part5, long long int colBit1, long long int rowBit1, long long int colBit2, long long int rowBit2)
Called 12 times for every 16 amplitudes in density matrix Each sums from the |..0....
Definition: QuEST_gpu.cu:2644
void statevec_initClassicalState(Qureg qureg, long long int stateInd)
Definition: QuEST_gpu.cu:600
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 ...
Definition: QuEST_gpu.cu:971
void statevec_multiControlledPhaseFlip(Qureg qureg, int *controlQubits, int numControlQubits)
Definition: QuEST_gpu.cu:1633
__global__ void statevec_controlledCompactUnitaryKernel(Qureg qureg, int controlQubit, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_gpu.cu:784
#define qreal
void statevec_unitary(Qureg qureg, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_gpu.cu:904
void statevec_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal outcomeProb)
Definition: QuEST_gpu.cu:2500
void densmatr_initClassicalState(Qureg qureg, long long int stateInd)
Definition: QuEST_gpu.cu:258
__forceinline__ __device__ long long int flipBit(const long long int number, const int bitInd)
Definition: QuEST_gpu.cu:95
__global__ void statevec_unitaryKernel(Qureg qureg, int targetQubit, ArgMatrix2 u)
Definition: QuEST_gpu.cu:851
__global__ void statevec_multiControlledPhaseFlipKernel(Qureg qureg, long long int mask)
Definition: QuEST_gpu.cu:1615
__global__ void statevec_findProbabilityOfZeroKernel(Qureg qureg, int measureQubit, qreal *reducedArray)
Definition: QuEST_gpu.cu:1853
__global__ void statevec_applyDiagonalOpKernel(Qureg qureg, DiagonalOp op)
Definition: QuEST_gpu.cu:2917
int numQubitsInStateVec
Number of qubits in the state-vector - this is double the number represented for mixed states.
Definition: QuEST.h:210
__global__ void statevec_initBlankStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag)
Definition: QuEST_gpu.cu:515
qreal statevec_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
Definition: QuEST_gpu.cu:2005
__global__ void densmatr_calcInnerProductKernel(Qureg a, Qureg b, long long int numTermsToSum, qreal *reducedArray)
computes Tr(conjTrans(a) b) = sum of (a_ij^* b_ij), which is a real number
Definition: QuEST_gpu.cu:2022
int chunkId
The position of the chunk of the state vector held by this process in the full state vector.
Definition: QuEST.h:217
void statevec_phaseShiftByTerm(Qureg qureg, int targetQubit, Complex term)
Definition: QuEST_gpu.cu:1423
int statevec_compareStates(Qureg mq1, Qureg mq2, qreal precision)
Definition: QuEST_gpu.cu:703
__forceinline__ __device__ long long int insertZeroBit(const long long int number, const int index)
Definition: QuEST_gpu.cu:99
void densmatr_mixDensityMatrix(Qureg combineQureg, qreal otherProb, Qureg otherQureg)
Definition: QuEST_gpu.cu:2576
qreal densmatr_findProbabilityOfZero(Qureg qureg, int measureQubit)
Definition: QuEST_gpu.cu:1919
__global__ void statevec_pauliXKernel(Qureg qureg, int targetQubit)
Definition: QuEST_gpu.cu:1248
__forceinline__ __device__ int getBitMaskParity(long long int mask)
Definition: QuEST_gpu.cu:86
long long int numAmpsPerChunk
Number of probability amplitudes held in stateVec by this process In the non-MPI version,...
Definition: QuEST.h:213
void statevec_applyDiagonalOp(Qureg qureg, DiagonalOp op)
Definition: QuEST_gpu.cu:2939
void densmatr_initPureState(Qureg targetQureg, Qureg copyQureg)
Definition: QuEST_gpu.cu:205
__global__ void statevec_phaseShiftByTermKernel(Qureg qureg, int targetQubit, qreal cosAngle, qreal sinAngle)
Definition: QuEST_gpu.cu:1395
__global__ void statevec_initZeroStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag)
Definition: QuEST_gpu.cu:536
void copyStateToGPU(Qureg qureg)
In GPU mode, this copies the state-vector (or density matrix) from RAM (qureg.stateVec) to VRAM / GPU...
Definition: QuEST_gpu.cu:451
__global__ void densmatr_calcPurityKernel(qreal *vecReal, qreal *vecImag, long long int numAmpsToSum, qreal *reducedArray)
Definition: QuEST_gpu.cu:2375
DiagonalOp agnostic_createDiagonalOp(int numQubits, QuESTEnv env)
Definition: QuEST_gpu.cu:338
qreal statevec_getImagAmp(Qureg qureg, long long int index)
Definition: QuEST_gpu.cu:508
qreal densmatr_calcProbOfOutcome(Qureg qureg, int measureQubit, int outcome)
Definition: QuEST_gpu.cu:2013
void statevec_compactUnitary(Qureg qureg, int targetQubit, Complex alpha, Complex beta)
Definition: QuEST_gpu.cu:776
int numQubits
The number of qubits this operator can act on (informing its size)
Definition: QuEST.h:181
long long int getQubitBitMask(int *qubits, int numQubits)
Definition: QuEST_common.c:44
Represents a diagonal complex operator on the full Hilbert state of a Qureg.
Definition: QuEST.h:178
__forceinline__ __device__ long long int insertZeroBits(long long int number, int *inds, const int numInds)
Definition: QuEST_gpu.cu:112
Complex densmatr_calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
Definition: QuEST_gpu.cu:3142
__global__ void densmatr_mixDephasingKernel(qreal fac, qreal *vecReal, qreal *vecImag, long long int numAmpsToVisit, long long int part1, long long int part2, long long int part3, long long int colBit, long long int rowBit)
Called once for every 4 amplitudes in density matrix Works by establishing the |.....
Definition: QuEST_gpu.cu:2593
int getNumReductionLevels(long long int numValuesToReduce, int numReducedPerLevel)
Definition: QuEST_gpu.cu:1903
__global__ void statevec_hadamardKernel(Qureg qureg, int targetQubit)
Definition: QuEST_gpu.cu:1676
__global__ void statevec_controlledPhaseFlipKernel(Qureg qureg, int idQubit1, int idQubit2)
Definition: QuEST_gpu.cu:1586
__global__ void densmatr_initPlusStateKernel(long long int stateVecSize, qreal probFactor, qreal *stateVecReal, qreal *stateVecImag)
Definition: QuEST_gpu.cu:216
__forceinline__ __device__ int extractBit(const int locationOfBitFromRight, const long long int theEncodedNumber)
Definition: QuEST_gpu.cu:82
__global__ void densmatr_mixDepolarisingKernel(qreal depolLevel, qreal *vecReal, qreal *vecImag, long long int numAmpsToVisit, long long int part1, long long int part2, long long int part3, long long int bothBits)
Works like mixDephasing but modifies every other element, and elements are averaged in pairs.
Definition: QuEST_gpu.cu:2705
void densmatr_mixTwoQubitDepolarising(Qureg qureg, int qubit1, int qubit2, qreal depolLevel)
Definition: QuEST_gpu.cu:2838
void statevec_controlledPhaseShift(Qureg qureg, int idQubit1, int idQubit2, qreal angle)
Definition: QuEST_gpu.cu:1459
qreal densmatr_calcHilbertSchmidtDistance(Qureg a, Qureg b)
Definition: QuEST_gpu.cu:2323
void densmatr_mixTwoQubitDephasing(Qureg qureg, int qubit1, int qubit2, qreal dephase)
Definition: QuEST_gpu.cu:2668
void statevec_controlledPauliY(Qureg qureg, int controlQubit, int targetQubit)
Definition: QuEST_gpu.cu:1377
__global__ void statevec_initDebugStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag)
Definition: QuEST_gpu.cu:611
__forceinline__ __device__ long long int insertTwoZeroBits(const long long int number, const int bit1, const int bit2)
Definition: QuEST_gpu.cu:106
Complex statevec_calcInnerProduct(Qureg bra, Qureg ket)
Terrible code which unnecessarily individually computes and sums the real and imaginary components of...
Definition: QuEST_gpu.cu:2123
void densmatr_collapseToKnownProbOutcome(Qureg qureg, int measureQubit, int outcome, qreal outcomeProb)
This involves finding |...i...><...j...| states and killing those where i!=j.
Definition: QuEST_gpu.cu:2535
__global__ void statevec_pauliYKernel(Qureg qureg, int targetQubit, int conjFac)
Definition: QuEST_gpu.cu:1300
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_gpu.cu:475
__global__ void statevec_initClassicalStateKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag, long long int stateInd)
Definition: QuEST_gpu.cu:585
ComplexArray stateVec
Computational state amplitudes - a subset thereof in the MPI version.
Definition: QuEST.h:222
__global__ void statevec_calcInnerProductKernel(int getRealComp, qreal *vecReal1, qreal *vecImag1, qreal *vecReal2, qreal *vecImag2, long long int numTermsToSum, qreal *reducedArray)
computes either a real or imag term in the inner product
Definition: QuEST_gpu.cu:2093
__global__ void densmatr_findProbabilityOfZeroKernel(Qureg qureg, int measureQubit, qreal *reducedArray)
Definition: QuEST_gpu.cu:1815
void seedQuESTDefault()
Seed the Mersenne Twister used for random number generation in the QuEST environment with an example ...
Definition: QuEST_gpu.cu:3252
long long int numElemsPerChunk
The number of the 2^numQubits amplitudes stored on each distributed node.
Definition: QuEST.h:183
__global__ void statevec_controlledPauliYKernel(Qureg qureg, int controlQubit, int targetQubit, int conjFac)
Definition: QuEST_gpu.cu:1341
void reportQuESTEnv(QuESTEnv env)
Report information about the QuEST environment.
Definition: QuEST_gpu.cu:435
void densmatr_mixDephasing(Qureg qureg, int targetQubit, qreal dephase)
Definition: QuEST_gpu.cu:2629
void statevec_controlledUnitary(Qureg qureg, int controlQubit, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_gpu.cu:1169
void statevec_pauliYConj(Qureg qureg, int targetQubit)
Definition: QuEST_gpu.cu:1333
void statevec_destroyQureg(Qureg qureg, QuESTEnv env)
Definition: QuEST_gpu.cu:321
void statevec_controlledNot(Qureg qureg, int controlQubit, int targetQubit)
Definition: QuEST_gpu.cu:1772
int numQubitsRepresented
The number of qubits represented in either the state-vector or density matrix.
Definition: QuEST.h:208
void statevec_cloneQureg(Qureg targetQureg, Qureg copyQureg)
works for both statevectors and density matrices
Definition: QuEST_gpu.cu:170
void agnostic_destroyDiagonalOp(DiagonalOp op)
Definition: QuEST_gpu.cu:375
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.
Definition: QuEST_gpu.cu:427
__device__ void reduceBlock(qreal *arrayIn, qreal *reducedArray, int length)
Definition: QuEST_gpu.cu:1787
__global__ void statevec_initStateOfSingleQubitKernel(long long int stateVecSize, qreal *stateVecReal, qreal *stateVecImag, int qubitId, int outcome)
Definition: QuEST_gpu.cu:632
__global__ void statevec_collapseToKnownProbOutcomeKernel(Qureg qureg, int measureQubit, int outcome, qreal totalProbability)
Definition: QuEST_gpu.cu:2443
void densmatr_applyDiagonalOp(Qureg qureg, DiagonalOp op)
Definition: QuEST_gpu.cu:2970
void statevec_initDebugState(Qureg qureg)
Initialise the state vector of probability amplitudes to an (unphysical) state with each component of...
Definition: QuEST_gpu.cu:621
void statevec_multiControlledUnitary(Qureg qureg, long long int ctrlQubitsMask, long long int ctrlFlipMask, int targetQubit, ComplexMatrix2 u)
Definition: QuEST_gpu.cu:1237
qreal densmatr_calcPurity(Qureg qureg)
Computes the trace of the density matrix squared.
Definition: QuEST_gpu.cu:2394
__global__ void statevec_multiControlledTwoQubitUnitaryKernel(Qureg qureg, long long int ctrlMask, int q1, int q2, ArgMatrix4 u)
Definition: QuEST_gpu.cu:1028
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_gpu.cu:650
void statevec_setWeightedQureg(Complex fac1, Qureg qureg1, Complex fac2, Qureg qureg2, Complex facOut, Qureg out)
Definition: QuEST_gpu.cu:2905
qreal densmatr_calcFidelity(Qureg qureg, Qureg pureState)
Definition: QuEST_gpu.cu:2249
__global__ void statevec_compactUnitaryKernel(Qureg qureg, int rotQubit, Complex alpha, Complex beta)
Definition: QuEST_gpu.cu:721
__global__ void statevec_controlledNotKernel(Qureg qureg, int controlQubit, int targetQubit)
Definition: QuEST_gpu.cu:1733
Complex statevec_calcExpecDiagonalOp(Qureg qureg, DiagonalOp op)
Definition: QuEST_gpu.cu:3006
void agnostic_setDiagonalOpElems(DiagonalOp op, long long int startInd, qreal *real, qreal *imag, long long int numElems)
Definition: QuEST_gpu.cu:3233
__global__ void densmatr_mixDensityMatrixKernel(Qureg combineQureg, qreal otherProb, Qureg otherQureg, long long int numAmpsToVisit)
Definition: QuEST_gpu.cu:2564
void densmatr_mixDepolarising(Qureg qureg, int targetQubit, qreal depolLevel)
Definition: QuEST_gpu.cu:2752
__global__ void statevec_swapQubitAmpsKernel(Qureg qureg, int qb1, int qb2)
Definition: QuEST_gpu.cu:1642
void statevec_controlledPauliYConj(Qureg qureg, int controlQubit, int targetQubit)
Definition: QuEST_gpu.cu:1386
__global__ void densmatr_calcFidelityKernel(Qureg dens, Qureg vec, long long int dim, qreal *reducedArray)
computes one term of (vec^*T) dens * vec
Definition: QuEST_gpu.cu:2211
__global__ void statevec_multiControlledPhaseShiftKernel(Qureg qureg, long long int mask, qreal cosAngle, qreal sinAngle)
Definition: QuEST_gpu.cu:1470