Functions
QuEST_cpu_internal.h File Reference
#include "QuEST_precision.h"

Go to the source code of this file.

Functions

qreal densmatr_calcFidelityLocal (Qureg qureg, Qureg pureState)
 computes a few dens-columns-worth of (vec^*T) dens * vec More...
 
qreal densmatr_calcHilbertSchmidtDistanceSquaredLocal (Qureg a, Qureg b)
 computes Tr((a-b) conjTrans(a-b)) = sum of abs values of (a-b) More...
 
qreal densmatr_calcInnerProductLocal (Qureg a, Qureg b)
 computes Tr(conjTrans(a) b) = sum of (a_ij^* b_ij) More...
 
qreal densmatr_calcPurityLocal (Qureg qureg)
 
qreal densmatr_findProbabilityOfZeroLocal (Qureg qureg, const int measureQubit)
 
void densmatr_initPureStateLocal (Qureg targetQureg, Qureg copyQureg)
 
void densmatr_mixDampingDistributed (Qureg qureg, const int targetQubit, qreal damping)
 
void densmatr_mixDampingLocal (Qureg qureg, const int targetQubit, qreal damping)
 
void densmatr_mixDepolarisingDistributed (Qureg qureg, const int targetQubit, qreal depolLevel)
 
void densmatr_mixDepolarisingLocal (Qureg qureg, const int targetQubit, qreal depolLevel)
 
void densmatr_mixTwoQubitDepolarisingDistributed (Qureg qureg, const int targetQubit, const int qubit2, qreal delta, qreal gamma)
 
void densmatr_mixTwoQubitDepolarisingLocal (Qureg qureg, int qubit1, int qubit2, qreal delta, qreal gamma)
 
void densmatr_mixTwoQubitDepolarisingLocalPart1 (Qureg qureg, int qubit1, int qubit2, qreal delta)
 
void densmatr_mixTwoQubitDepolarisingQ1LocalQ2DistributedPart3 (Qureg qureg, const int targetQubit, const int qubit2, qreal delta, qreal gamma)
 
static int extractBit (const int locationOfBitFromRight, const long long int theEncodedNumber)
 
static long long int flipBit (long long int number, int bitInd)
 
static long long int insertTwoZeroBits (long long int number, int bit1, int bit2)
 
static long long int insertZeroBit (long long int number, int index)
 
static int isOddParity (long long int number, int qb1, int qb2)
 
static int maskContainsBit (long long int mask, int bitInd)
 
Complex statevec_calcInnerProductLocal (Qureg bra, Qureg ket)
 
void statevec_collapseToKnownProbOutcomeDistributedRenorm (Qureg qureg, const int measureQubit, const qreal totalProbability)
 Renormalise parts of the state vector where measureQubit=0 or 1, based on the total probability of that qubit being in state 0 or 1. More...
 
void statevec_collapseToKnownProbOutcomeLocal (Qureg qureg, int measureQubit, int outcome, qreal totalProbability)
 Update the state vector to be consistent with measuring measureQubit=0 if outcome=0 and measureQubit=1 if outcome=1. More...
 
void statevec_collapseToOutcomeDistributedSetZero (Qureg qureg)
 Set all amplitudes in one chunk to 0. More...
 
void statevec_compactUnitaryDistributed (Qureg qureg, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
 Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha and beta, and a subset of the state vector with upper and lower block values stored seperately. More...
 
void statevec_compactUnitaryLocal (Qureg qureg, const int targetQubit, Complex alpha, Complex beta)
 
void statevec_controlledCompactUnitaryDistributed (Qureg qureg, const int controlQubit, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
 Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha and beta and a subset of the state vector with upper and lower block values stored seperately. More...
 
void statevec_controlledCompactUnitaryLocal (Qureg qureg, const int controlQubit, const int targetQubit, Complex alpha, Complex beta)
 
void statevec_controlledNotDistributed (Qureg qureg, const int controlQubit, ComplexArray stateVecIn, ComplexArray stateVecOut)
 Rotate a single qubit by {{0,1},{1,0}. More...
 
void statevec_controlledNotLocal (Qureg qureg, const int controlQubit, const int targetQubit)
 
void statevec_controlledPauliYDistributed (Qureg qureg, const int controlQubit, ComplexArray stateVecIn, ComplexArray stateVecOut, const int conjFactor)
 
void statevec_controlledPauliYLocal (Qureg qureg, const int controlQubit, const int targetQubit, const int conjFactor)
 
void statevec_controlledUnitaryDistributed (Qureg qureg, const int controlQubit, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
 Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha and beta and a subset of the state vector with upper and lower block values stored seperately. More...
 
void statevec_controlledUnitaryLocal (Qureg qureg, const int controlQubit, const int targetQubit, ComplexMatrix2 u)
 
qreal statevec_findProbabilityOfZeroDistributed (Qureg qureg)
 Measure the probability of a specified qubit being in the zero state across all amplitudes held in this chunk. More...
 
qreal statevec_findProbabilityOfZeroLocal (Qureg qureg, const int measureQubit)
 Measure the total probability of a specified qubit being in the zero state across all amplitudes in this chunk. More...
 
void statevec_hadamardDistributed (Qureg qureg, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut, int updateUpper)
 Rotate a single qubit by {{1,1},{1,-1}}/sqrt2. More...
 
void statevec_hadamardLocal (Qureg qureg, const int targetQubit)
 
void statevec_multiControlledMultiQubitUnitaryLocal (Qureg qureg, long long int ctrlMask, int *targs, const int numTargs, ComplexMatrixN u)
 
void statevec_multiControlledTwoQubitUnitaryLocal (Qureg qureg, long long int ctrlMask, const int q1, const int q2, ComplexMatrix4 u)
 
void statevec_multiControlledUnitaryDistributed (Qureg qureg, const int targetQubit, long long int ctrlQubitsMask, long long int ctrlFlipMask, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
 Apply a unitary operation to a single qubit in the state vector of probability amplitudes, given a subset of the state vector with upper and lower block values stored seperately. More...
 
void statevec_multiControlledUnitaryLocal (Qureg qureg, const int targetQubit, long long int ctrlQubitsMask, long long int ctrlFlipMask, ComplexMatrix2 u)
 
void statevec_pauliXDistributed (Qureg qureg, ComplexArray stateVecIn, ComplexArray stateVecOut)
 Rotate a single qubit by {{0,1},{1,0}. More...
 
void statevec_pauliXLocal (Qureg qureg, const int targetQubit)
 
void statevec_pauliYDistributed (Qureg qureg, ComplexArray stateVecIn, ComplexArray stateVecOut, int updateUpper, const int conjFac)
 Rotate a single qubit by +-{{0,-i},{i,0}. More...
 
void statevec_pauliYLocal (Qureg qureg, const int targetQubit, const int conjFac)
 
void statevec_swapQubitAmpsDistributed (Qureg qureg, int pairRank, int qb1, int qb2)
 qureg.pairStateVec contains the entire set of amplitudes of the paired node which includes the set of all amplitudes which need to be swapped between |..0..1..> and |..1..0..> More...
 
void statevec_swapQubitAmpsLocal (Qureg qureg, int qb1, int qb2)
 It is ensured that all amplitudes needing to be swapped are on this node. More...
 
void statevec_unitaryDistributed (Qureg qureg, Complex rot1, Complex rot2, ComplexArray stateVecUp, ComplexArray stateVecLo, ComplexArray stateVecOut)
 Apply a unitary operation to a single qubit given a subset of the state vector with upper and lower block values stored seperately. More...
 
void statevec_unitaryLocal (Qureg qureg, const int targetQubit, ComplexMatrix2 u)
 

Detailed Description

Internal functions used to implement the pure backend in ../QuEST_ops_pure.h. Do not call these functions directly. In general, qubits_cpu_local.c and qubits_cpu_mpi.c will implement the pure backend by choosing the correct function or combination of functions to use from those included here, which are defined in QuEST_cpu.c

Author
Ania Brown
Tyson Jones
Balint Koczor

Definition in file QuEST_cpu_internal.h.

Function Documentation

◆ densmatr_calcFidelityLocal()

qreal densmatr_calcFidelityLocal ( Qureg  qureg,
Qureg  pureState 
)

computes a few dens-columns-worth of (vec^*T) dens * vec

Definition at line 989 of file QuEST_cpu.c.

989  {
990 
991  /* Here, elements of pureState are not accessed (instead grabbed from qureg.pair).
992  * We only consult the attributes.
993  *
994  * qureg is a density matrix, and pureState is a statevector.
995  * Every node contains as many columns of qureg as amps by pureState.
996  * (each node contains an integer, exponent-of-2 number of whole columns of qureg)
997  * Ergo, this node contains columns:
998  * qureg.chunkID * pureState.numAmpsPerChunk to
999  * (qureg.chunkID + 1) * pureState.numAmpsPerChunk
1000  *
1001  * The first pureState.numAmpsTotal elements of qureg.pairStateVec are the
1002  * entire pureState state-vector
1003  */
1004 
1005  // unpack everything for OPENMP
1006  qreal* vecRe = qureg.pairStateVec.real;
1007  qreal* vecIm = qureg.pairStateVec.imag;
1008  qreal* densRe = qureg.stateVec.real;
1009  qreal* densIm = qureg.stateVec.imag;
1010 
1011  int row, col;
1012  int dim = (int) pureState.numAmpsTotal;
1013  int colsPerNode = (int) pureState.numAmpsPerChunk;
1014  // using only int, because density matrix has squared as many amps so its
1015  // iteration would be impossible if the pureStates numAmpsTotal didn't fit into int
1016 
1017  // starting GLOBAL column index of the qureg columns on this node
1018  int startCol = (int) (qureg.chunkId * pureState.numAmpsPerChunk);
1019 
1020  qreal densElemRe, densElemIm;
1021  qreal prefacRe, prefacIm;
1022  qreal rowSumRe, rowSumIm;
1023  qreal vecElemRe, vecElemIm;
1024 
1025  // quantity computed by this node
1026  qreal globalSumRe = 0; // imag-component is assumed zero
1027 
1028 # ifdef _OPENMP
1029 # pragma omp parallel \
1030  shared (vecRe,vecIm,densRe,densIm, dim,colsPerNode,startCol) \
1031  private (row,col, prefacRe,prefacIm, rowSumRe,rowSumIm, densElemRe,densElemIm, vecElemRe,vecElemIm) \
1032  reduction ( +:globalSumRe )
1033 # endif
1034  {
1035 # ifdef _OPENMP
1036 # pragma omp for schedule (static)
1037 # endif
1038  // indices of my GLOBAL row
1039  for (row=0; row < dim; row++) {
1040 
1041  // single element of conj(pureState)
1042  prefacRe = vecRe[row];
1043  prefacIm = - vecIm[row];
1044 
1045  rowSumRe = 0;
1046  rowSumIm = 0;
1047 
1048  // indices of my LOCAL column
1049  for (col=0; col < colsPerNode; col++) {
1050 
1051  // my local density element
1052  densElemRe = densRe[row + dim*col];
1053  densElemIm = densIm[row + dim*col];
1054 
1055  // state-vector element
1056  vecElemRe = vecRe[startCol + col];
1057  vecElemIm = vecIm[startCol + col];
1058 
1059  rowSumRe += densElemRe*vecElemRe - densElemIm*vecElemIm;
1060  rowSumIm += densElemRe*vecElemIm + densElemIm*vecElemRe;
1061  }
1062 
1063  globalSumRe += rowSumRe*prefacRe - rowSumIm*prefacIm;
1064  }
1065  }
1066 
1067  return globalSumRe;
1068 }

References Qureg::chunkId, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, Qureg::pairStateVec, qreal, and Qureg::stateVec.

Referenced by densmatr_calcFidelity().

◆ densmatr_calcHilbertSchmidtDistanceSquaredLocal()

qreal densmatr_calcHilbertSchmidtDistanceSquaredLocal ( Qureg  a,
Qureg  b 
)

computes Tr((a-b) conjTrans(a-b)) = sum of abs values of (a-b)

Definition at line 922 of file QuEST_cpu.c.

922  {
923 
924  long long int index;
925  long long int numAmps = a.numAmpsPerChunk;
926 
927  qreal *aRe = a.stateVec.real;
928  qreal *aIm = a.stateVec.imag;
929  qreal *bRe = b.stateVec.real;
930  qreal *bIm = b.stateVec.imag;
931 
932  qreal trace = 0;
933  qreal difRe, difIm;
934 
935 # ifdef _OPENMP
936 # pragma omp parallel \
937  shared (aRe,aIm, bRe,bIm, numAmps) \
938  private (index,difRe,difIm) \
939  reduction ( +:trace )
940 # endif
941  {
942 # ifdef _OPENMP
943 # pragma omp for schedule (static)
944 # endif
945  for (index=0LL; index<numAmps; index++) {
946 
947  difRe = aRe[index] - bRe[index];
948  difIm = aIm[index] - bIm[index];
949  trace += difRe*difRe + difIm*difIm;
950  }
951  }
952 
953  return trace;
954 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by densmatr_calcHilbertSchmidtDistance().

◆ densmatr_calcInnerProductLocal()

qreal densmatr_calcInnerProductLocal ( Qureg  a,
Qureg  b 
)

computes Tr(conjTrans(a) b) = sum of (a_ij^* b_ij)

Definition at line 957 of file QuEST_cpu.c.

957  {
958 
959  long long int index;
960  long long int numAmps = a.numAmpsPerChunk;
961 
962  qreal *aRe = a.stateVec.real;
963  qreal *aIm = a.stateVec.imag;
964  qreal *bRe = b.stateVec.real;
965  qreal *bIm = b.stateVec.imag;
966 
967  qreal trace = 0;
968 
969 # ifdef _OPENMP
970 # pragma omp parallel \
971  shared (aRe,aIm, bRe,bIm, numAmps) \
972  private (index) \
973  reduction ( +:trace )
974 # endif
975  {
976 # ifdef _OPENMP
977 # pragma omp for schedule (static)
978 # endif
979  for (index=0LL; index<numAmps; index++) {
980  trace += aRe[index]*bRe[index] + aIm[index]*bIm[index];
981  }
982  }
983 
984  return trace;
985 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by densmatr_calcInnerProduct().

◆ densmatr_calcPurityLocal()

qreal densmatr_calcPurityLocal ( Qureg  qureg)

Definition at line 860 of file QuEST_cpu.c.

860  {
861 
862  /* sum of qureg^2, which is sum_i |qureg[i]|^2 */
863  long long int index;
864  long long int numAmps = qureg.numAmpsPerChunk;
865 
866  qreal trace = 0;
867  qreal *vecRe = qureg.stateVec.real;
868  qreal *vecIm = qureg.stateVec.imag;
869 
870 # ifdef _OPENMP
871 # pragma omp parallel \
872  shared (vecRe, vecIm, numAmps) \
873  private (index) \
874  reduction ( +:trace )
875 # endif
876  {
877 # ifdef _OPENMP
878 # pragma omp for schedule (static)
879 # endif
880  for (index=0LL; index<numAmps; index++) {
881 
882  trace += vecRe[index]*vecRe[index] + vecIm[index]*vecIm[index];
883  }
884  }
885 
886  return trace;
887 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by densmatr_calcPurity().

◆ densmatr_findProbabilityOfZeroLocal()

qreal densmatr_findProbabilityOfZeroLocal ( Qureg  qureg,
const int  measureQubit 
)

Definition at line 3111 of file QuEST_cpu.c.

3111  {
3112 
3113  // computes first local index containing a diagonal element
3114  long long int localNumAmps = qureg.numAmpsPerChunk;
3115  long long int densityDim = (1LL << qureg.numQubitsRepresented);
3116  long long int diagSpacing = 1LL + densityDim;
3117  long long int maxNumDiagsPerChunk = 1 + localNumAmps / diagSpacing;
3118  long long int numPrevDiags = (qureg.chunkId>0)? 1+(qureg.chunkId*localNumAmps)/diagSpacing : 0;
3119  long long int globalIndNextDiag = diagSpacing * numPrevDiags;
3120  long long int localIndNextDiag = globalIndNextDiag % localNumAmps;
3121 
3122  // computes how many diagonals are contained in this chunk
3123  long long int numDiagsInThisChunk = maxNumDiagsPerChunk;
3124  if (localIndNextDiag + (numDiagsInThisChunk-1)*diagSpacing >= localNumAmps)
3125  numDiagsInThisChunk -= 1;
3126 
3127  long long int visitedDiags; // number of visited diagonals in this chunk so far
3128  long long int basisStateInd; // current diagonal index being considered
3129  long long int index; // index in the local chunk
3130 
3131  qreal zeroProb = 0;
3132  qreal *stateVecReal = qureg.stateVec.real;
3133 
3134 # ifdef _OPENMP
3135 # pragma omp parallel \
3136  shared (localIndNextDiag, numPrevDiags, diagSpacing, stateVecReal, numDiagsInThisChunk) \
3137  private (visitedDiags, basisStateInd, index) \
3138  reduction ( +:zeroProb )
3139 # endif
3140  {
3141 # ifdef _OPENMP
3142 # pragma omp for schedule (static)
3143 # endif
3144  // sums the diagonal elems of the density matrix where measureQubit=0
3145  for (visitedDiags = 0; visitedDiags < numDiagsInThisChunk; visitedDiags++) {
3146 
3147  basisStateInd = numPrevDiags + visitedDiags;
3148  index = localIndNextDiag + diagSpacing * visitedDiags;
3149 
3150  if (extractBit(measureQubit, basisStateInd) == 0)
3151  zeroProb += stateVecReal[index]; // assume imag[diagonls] ~ 0
3152 
3153  }
3154  }
3155 
3156  return zeroProb;
3157 }

References Qureg::chunkId, extractBit(), Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, qreal, and Qureg::stateVec.

Referenced by densmatr_calcProbOfOutcome().

◆ densmatr_initPureStateLocal()

void densmatr_initPureStateLocal ( Qureg  targetQureg,
Qureg  copyQureg 
)

Definition at line 1183 of file QuEST_cpu.c.

1183  {
1184 
1185  /* copyQureg amps aren't explicitly used - they're accessed through targetQureg.pair,
1186  * which contains the full pure statevector.
1187  * targetQureg has as many columns on node as copyQureg has amps
1188  */
1189 
1190  long long int colOffset = targetQureg.chunkId * copyQureg.numAmpsPerChunk;
1191  long long int colsPerNode = copyQureg.numAmpsPerChunk;
1192  long long int rowsPerNode = copyQureg.numAmpsTotal;
1193 
1194  // unpack vars for OpenMP
1195  qreal* vecRe = targetQureg.pairStateVec.real;
1196  qreal* vecIm = targetQureg.pairStateVec.imag;
1197  qreal* densRe = targetQureg.stateVec.real;
1198  qreal* densIm = targetQureg.stateVec.imag;
1199 
1200  long long int col, row, index;
1201 
1202  // a_i conj(a_j) |i><j|
1203  qreal ketRe, ketIm, braRe, braIm;
1204 
1205 # ifdef _OPENMP
1206 # pragma omp parallel \
1207  default (none) \
1208  shared (colOffset, colsPerNode,rowsPerNode, vecRe,vecIm,densRe,densIm) \
1209  private (col,row, ketRe,ketIm,braRe,braIm, index)
1210 # endif
1211  {
1212 # ifdef _OPENMP
1213 # pragma omp for schedule (static)
1214 # endif
1215  // local column
1216  for (col=0; col < colsPerNode; col++) {
1217 
1218  // global row
1219  for (row=0; row < rowsPerNode; row++) {
1220 
1221  // get pure state amps
1222  ketRe = vecRe[row];
1223  ketIm = vecIm[row];
1224  braRe = vecRe[col + colOffset];
1225  braIm = - vecIm[col + colOffset]; // minus for conjugation
1226 
1227  // update density matrix
1228  index = row + col*rowsPerNode; // local ind
1229  densRe[index] = ketRe*braRe - ketIm*braIm;
1230  densIm[index] = ketRe*braIm + ketIm*braRe;
1231  }
1232  }
1233  }
1234 }

References Qureg::chunkId, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, Qureg::pairStateVec, qreal, and Qureg::stateVec.

Referenced by densmatr_initPureState().

◆ densmatr_mixDampingDistributed()

void densmatr_mixDampingDistributed ( Qureg  qureg,
const int  targetQubit,
qreal  damping 
)

Definition at line 299 of file QuEST_cpu.c.

299  {
300  qreal retain=1-damping;
301  qreal dephase=sqrt(1-damping);
302 
303  // multiply the off-diagonal (|0><1| and |1><0|) terms by sqrt(1-damping)
304  densmatr_oneQubitDegradeOffDiagonal(qureg, targetQubit, dephase);
305 
306  // below, we modify the diagonals terms which require |1><1| to |0><0| communication
307 
308  long long int sizeInnerBlock, sizeInnerHalfBlock;
309  long long int sizeOuterColumn, sizeOuterHalfColumn;
310  long long int thisInnerBlock, // current block
311  thisOuterColumn, // current column in density matrix
312  thisIndex, // current index in (density matrix representation) state vector
313  thisIndexInOuterColumn,
314  thisIndexInInnerBlock;
315  int outerBit;
316  int stateBit;
317 
318  long long int thisTask;
319  const long long int numTasks=qureg.numAmpsPerChunk>>1;
320 
321  // set dimensions
322  sizeInnerHalfBlock = 1LL << targetQubit;
323  sizeInnerBlock = 2LL * sizeInnerHalfBlock;
324  sizeOuterColumn = 1LL << qureg.numQubitsRepresented;
325  sizeOuterHalfColumn = sizeOuterColumn >> 1;
326 
327 # ifdef _OPENMP
328 # pragma omp parallel \
329  default (none) \
330  shared (sizeInnerBlock,sizeInnerHalfBlock,sizeOuterColumn,sizeOuterHalfColumn,qureg,damping, retain, dephase) \
331  private (thisTask,thisInnerBlock,thisOuterColumn,thisIndex,thisIndexInOuterColumn, \
332  thisIndexInInnerBlock,outerBit, stateBit)
333 # endif
334  {
335 # ifdef _OPENMP
336 # pragma omp for schedule (static)
337 # endif
338  // thisTask iterates over half the elements in this process' chunk of the density matrix
339  // treat this as iterating over all columns, then iterating over half the values
340  // within one column.
341  // If this function has been called, this process' chunk contains half an
342  // outer block or less
343  for (thisTask=0; thisTask<numTasks; thisTask++) {
344  // we want to process all columns in the density matrix,
345  // updating the values for half of each column (one half of each inner block)
346  thisOuterColumn = thisTask / sizeOuterHalfColumn;
347  thisIndexInOuterColumn = thisTask&(sizeOuterHalfColumn-1); // thisTask % sizeOuterHalfColumn
348  thisInnerBlock = thisIndexInOuterColumn/sizeInnerHalfBlock;
349  // get index in state vector corresponding to upper inner block
350  thisIndexInInnerBlock = thisTask&(sizeInnerHalfBlock-1); // thisTask % sizeInnerHalfBlock
351  thisIndex = thisOuterColumn*sizeOuterColumn + thisInnerBlock*sizeInnerBlock
352  + thisIndexInInnerBlock;
353  // check if we are in the upper or lower half of an outer block
354  outerBit = extractBit(targetQubit, (thisIndex+qureg.numAmpsPerChunk*qureg.chunkId)>>qureg.numQubitsRepresented);
355  // if we are in the lower half of an outer block, shift to be in the lower half
356  // of the inner block as well (we want to dephase |0><0| and |1><1| only)
357  thisIndex += outerBit*(sizeInnerHalfBlock);
358 
359  // NOTE: at this point thisIndex should be the index of the element we want to
360  // dephase in the chunk of the state vector on this process, in the
361  // density matrix representation.
362  // thisTask is the index of the pair element in pairStateVec
363 
364  // Extract state bit, is 0 if thisIndex corresponds to a state with 0 in the target qubit
365  // and is 1 if thisIndex corresponds to a state with 1 in the target qubit
366  stateBit = extractBit(targetQubit, (thisIndex+qureg.numAmpsPerChunk*qureg.chunkId));
367 
368  // state[thisIndex] = (1-depolLevel)*state[thisIndex] + depolLevel*(state[thisIndex]
369  // + pair[thisTask])/2
370  if(stateBit == 0){
371  qureg.stateVec.real[thisIndex] = qureg.stateVec.real[thisIndex] +
372  damping*( qureg.pairStateVec.real[thisTask]);
373 
374  qureg.stateVec.imag[thisIndex] = qureg.stateVec.imag[thisIndex] +
375  damping*( qureg.pairStateVec.imag[thisTask]);
376  } else{
377  qureg.stateVec.real[thisIndex] = retain*qureg.stateVec.real[thisIndex];
378 
379  qureg.stateVec.imag[thisIndex] = retain*qureg.stateVec.imag[thisIndex];
380  }
381  }
382  }
383 }

References Qureg::chunkId, densmatr_oneQubitDegradeOffDiagonal(), extractBit(), Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, Qureg::pairStateVec, qreal, and Qureg::stateVec.

Referenced by densmatr_mixDamping().

◆ densmatr_mixDampingLocal()

void densmatr_mixDampingLocal ( Qureg  qureg,
const int  targetQubit,
qreal  damping 
)

Definition at line 174 of file QuEST_cpu.c.

174  {
175  qreal retain=1-damping;
176  qreal dephase=sqrt(retain);
177 
178  const long long int numTasks = qureg.numAmpsPerChunk;
179  long long int innerMask = 1LL << targetQubit;
180  long long int outerMask = 1LL << (targetQubit + (qureg.numQubitsRepresented));
181  long long int totMask = innerMask|outerMask;
182 
183  long long int thisTask;
184  long long int partner;
185  long long int thisPattern;
186 
187  //qreal realAv, imagAv;
188 
189 # ifdef _OPENMP
190 # pragma omp parallel \
191  default (none) \
192  shared (innerMask,outerMask,totMask,qureg,retain,damping,dephase) \
193  private (thisTask,partner,thisPattern)
194 # endif
195  {
196 # ifdef _OPENMP
197 # pragma omp for schedule (static)
198 # endif
199  for (thisTask=0; thisTask<numTasks; thisTask++){
200  thisPattern = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMask;
201  if ((thisPattern==innerMask) || (thisPattern==outerMask)){
202  // do dephase
203  // the lines below will degrade the off-diagonal terms |..0..><..1..| and |..1..><..0..|
204  qureg.stateVec.real[thisTask] = dephase*qureg.stateVec.real[thisTask];
205  qureg.stateVec.imag[thisTask] = dephase*qureg.stateVec.imag[thisTask];
206  } else {
207  if ((thisTask&totMask)==0){ //this element relates to targetQubit in state 0
208  // do depolarise
209  partner = thisTask | totMask;
210  //realAv = (qureg.stateVec.real[thisTask] + qureg.stateVec.real[partner]) /2 ;
211  //imagAv = (qureg.stateVec.imag[thisTask] + qureg.stateVec.imag[partner]) /2 ;
212 
213  qureg.stateVec.real[thisTask] = qureg.stateVec.real[thisTask] + damping*qureg.stateVec.real[partner];
214  qureg.stateVec.imag[thisTask] = qureg.stateVec.imag[thisTask] + damping*qureg.stateVec.imag[partner];
215 
216  qureg.stateVec.real[partner] = retain*qureg.stateVec.real[partner];
217  qureg.stateVec.imag[partner] = retain*qureg.stateVec.imag[partner];
218  }
219  }
220  }
221  }
222 }

References Qureg::chunkId, Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, qreal, and Qureg::stateVec.

Referenced by densmatr_mixDamping().

◆ densmatr_mixDepolarisingDistributed()

void densmatr_mixDepolarisingDistributed ( Qureg  qureg,
const int  targetQubit,
qreal  depolLevel 
)

Definition at line 224 of file QuEST_cpu.c.

224  {
225 
226  // first do dephase part.
227  // TODO -- this might be more efficient to do at the same time as the depolarise if we move to
228  // iterating over all elements in the state vector for the purpose of vectorisation
229  // TODO -- if we keep this split, move this function to densmatr_mixDepolarising()
230  densmatr_mixDephasing(qureg, targetQubit, depolLevel);
231 
232  long long int sizeInnerBlock, sizeInnerHalfBlock;
233  long long int sizeOuterColumn, sizeOuterHalfColumn;
234  long long int thisInnerBlock, // current block
235  thisOuterColumn, // current column in density matrix
236  thisIndex, // current index in (density matrix representation) state vector
237  thisIndexInOuterColumn,
238  thisIndexInInnerBlock;
239  int outerBit;
240 
241  long long int thisTask;
242  const long long int numTasks=qureg.numAmpsPerChunk>>1;
243 
244  // set dimensions
245  sizeInnerHalfBlock = 1LL << targetQubit;
246  sizeInnerBlock = 2LL * sizeInnerHalfBlock;
247  sizeOuterColumn = 1LL << qureg.numQubitsRepresented;
248  sizeOuterHalfColumn = sizeOuterColumn >> 1;
249 
250 # ifdef _OPENMP
251 # pragma omp parallel \
252  default (none) \
253  shared (sizeInnerBlock,sizeInnerHalfBlock,sizeOuterColumn,sizeOuterHalfColumn,qureg,depolLevel) \
254  private (thisTask,thisInnerBlock,thisOuterColumn,thisIndex,thisIndexInOuterColumn, \
255  thisIndexInInnerBlock,outerBit)
256 # endif
257  {
258 # ifdef _OPENMP
259 # pragma omp for schedule (static)
260 # endif
261  // thisTask iterates over half the elements in this process' chunk of the density matrix
262  // treat this as iterating over all columns, then iterating over half the values
263  // within one column.
264  // If this function has been called, this process' chunk contains half an
265  // outer block or less
266  for (thisTask=0; thisTask<numTasks; thisTask++) {
267  // we want to process all columns in the density matrix,
268  // updating the values for half of each column (one half of each inner block)
269  thisOuterColumn = thisTask / sizeOuterHalfColumn;
270  thisIndexInOuterColumn = thisTask&(sizeOuterHalfColumn-1); // thisTask % sizeOuterHalfColumn
271  thisInnerBlock = thisIndexInOuterColumn/sizeInnerHalfBlock;
272  // get index in state vector corresponding to upper inner block
273  thisIndexInInnerBlock = thisTask&(sizeInnerHalfBlock-1); // thisTask % sizeInnerHalfBlock
274  thisIndex = thisOuterColumn*sizeOuterColumn + thisInnerBlock*sizeInnerBlock
275  + thisIndexInInnerBlock;
276  // check if we are in the upper or lower half of an outer block
277  outerBit = extractBit(targetQubit, (thisIndex+qureg.numAmpsPerChunk*qureg.chunkId)>>qureg.numQubitsRepresented);
278  // if we are in the lower half of an outer block, shift to be in the lower half
279  // of the inner block as well (we want to dephase |0><0| and |1><1| only)
280  thisIndex += outerBit*(sizeInnerHalfBlock);
281 
282  // NOTE: at this point thisIndex should be the index of the element we want to
283  // dephase in the chunk of the state vector on this process, in the
284  // density matrix representation.
285  // thisTask is the index of the pair element in pairStateVec
286 
287 
288  // state[thisIndex] = (1-depolLevel)*state[thisIndex] + depolLevel*(state[thisIndex]
289  // + pair[thisTask])/2
290  qureg.stateVec.real[thisIndex] = (1-depolLevel)*qureg.stateVec.real[thisIndex] +
291  depolLevel*(qureg.stateVec.real[thisIndex] + qureg.pairStateVec.real[thisTask])/2;
292 
293  qureg.stateVec.imag[thisIndex] = (1-depolLevel)*qureg.stateVec.imag[thisIndex] +
294  depolLevel*(qureg.stateVec.imag[thisIndex] + qureg.pairStateVec.imag[thisTask])/2;
295  }
296  }
297 }

References Qureg::chunkId, densmatr_mixDephasing(), extractBit(), Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, Qureg::pairStateVec, and Qureg::stateVec.

Referenced by densmatr_mixDepolarising().

◆ densmatr_mixDepolarisingLocal()

void densmatr_mixDepolarisingLocal ( Qureg  qureg,
const int  targetQubit,
qreal  depolLevel 
)

Definition at line 125 of file QuEST_cpu.c.

125  {
126  qreal retain=1-depolLevel;
127 
128  const long long int numTasks = qureg.numAmpsPerChunk;
129  long long int innerMask = 1LL << targetQubit;
130  long long int outerMask = 1LL << (targetQubit + (qureg.numQubitsRepresented));
131  long long int totMask = innerMask|outerMask;
132 
133  long long int thisTask;
134  long long int partner;
135  long long int thisPattern;
136 
137  qreal realAv, imagAv;
138 
139 # ifdef _OPENMP
140 # pragma omp parallel \
141  default (none) \
142  shared (innerMask,outerMask,totMask,qureg,retain,depolLevel) \
143  private (thisTask,partner,thisPattern,realAv,imagAv)
144 # endif
145  {
146 # ifdef _OPENMP
147 # pragma omp for schedule (static)
148 # endif
149  for (thisTask=0; thisTask<numTasks; thisTask++){
150  thisPattern = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMask;
151  if ((thisPattern==innerMask) || (thisPattern==outerMask)){
152  // do dephase
153  // the lines below will degrade the off-diagonal terms |..0..><..1..| and |..1..><..0..|
154  qureg.stateVec.real[thisTask] = retain*qureg.stateVec.real[thisTask];
155  qureg.stateVec.imag[thisTask] = retain*qureg.stateVec.imag[thisTask];
156  } else {
157  if ((thisTask&totMask)==0){ //this element relates to targetQubit in state 0
158  // do depolarise
159  partner = thisTask | totMask;
160  realAv = (qureg.stateVec.real[thisTask] + qureg.stateVec.real[partner]) /2 ;
161  imagAv = (qureg.stateVec.imag[thisTask] + qureg.stateVec.imag[partner]) /2 ;
162 
163  qureg.stateVec.real[thisTask] = retain*qureg.stateVec.real[thisTask] + depolLevel*realAv;
164  qureg.stateVec.imag[thisTask] = retain*qureg.stateVec.imag[thisTask] + depolLevel*imagAv;
165 
166  qureg.stateVec.real[partner] = retain*qureg.stateVec.real[partner] + depolLevel*realAv;
167  qureg.stateVec.imag[partner] = retain*qureg.stateVec.imag[partner] + depolLevel*imagAv;
168  }
169  }
170  }
171  }
172 }

References Qureg::chunkId, Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, qreal, and Qureg::stateVec.

Referenced by densmatr_mixDepolarising().

◆ densmatr_mixTwoQubitDepolarisingDistributed()

void densmatr_mixTwoQubitDepolarisingDistributed ( Qureg  qureg,
const int  targetQubit,
const int  qubit2,
qreal  delta,
qreal  gamma 
)

Definition at line 540 of file QuEST_cpu.c.

541  {
542 
543  long long int sizeInnerBlockQ1, sizeInnerHalfBlockQ1;
544  long long int sizeInnerBlockQ2, sizeInnerHalfBlockQ2, sizeInnerQuarterBlockQ2;
545  long long int sizeOuterColumn, sizeOuterQuarterColumn;
546  long long int thisInnerBlockQ2,
547  thisOuterColumn, // current column in density matrix
548  thisIndex, // current index in (density matrix representation) state vector
549  thisIndexInOuterColumn,
550  thisIndexInInnerBlockQ1,
551  thisIndexInInnerBlockQ2,
552  thisInnerBlockQ1InInnerBlockQ2;
553  int outerBitQ1, outerBitQ2;
554 
555  long long int thisTask;
556  const long long int numTasks=qureg.numAmpsPerChunk>>2;
557 
558  // set dimensions
559  sizeInnerHalfBlockQ1 = 1LL << targetQubit;
560  sizeInnerHalfBlockQ2 = 1LL << qubit2;
561  sizeInnerQuarterBlockQ2 = sizeInnerHalfBlockQ2 >> 1;
562  sizeInnerBlockQ2 = sizeInnerHalfBlockQ2 << 1;
563  sizeInnerBlockQ1 = 2LL * sizeInnerHalfBlockQ1;
564  sizeOuterColumn = 1LL << qureg.numQubitsRepresented;
565  sizeOuterQuarterColumn = sizeOuterColumn >> 2;
566 
567 # ifdef _OPENMP
568 # pragma omp parallel \
569  default (none) \
570  shared (sizeInnerBlockQ1,sizeInnerHalfBlockQ1,sizeInnerBlockQ2,sizeInnerHalfBlockQ2,sizeInnerQuarterBlockQ2,\
571  sizeOuterColumn,sizeOuterQuarterColumn,qureg,delta,gamma) \
572  private (thisTask,thisInnerBlockQ2,thisInnerBlockQ1InInnerBlockQ2, \
573  thisOuterColumn,thisIndex,thisIndexInOuterColumn, \
574  thisIndexInInnerBlockQ1,thisIndexInInnerBlockQ2,outerBitQ1,outerBitQ2)
575 # endif
576  {
577 # ifdef _OPENMP
578 # pragma omp for schedule (static)
579 # endif
580  // thisTask iterates over half the elements in this process' chunk of the density matrix
581  // treat this as iterating over all columns, then iterating over half the values
582  // within one column.
583  // If this function has been called, this process' chunk contains half an
584  // outer block or less
585  for (thisTask=0; thisTask<numTasks; thisTask++) {
586  // we want to process all columns in the density matrix,
587  // updating the values for half of each column (one half of each inner block)
588  thisOuterColumn = thisTask / sizeOuterQuarterColumn;
589  // thisTask % sizeOuterQuarterColumn
590  thisIndexInOuterColumn = thisTask&(sizeOuterQuarterColumn-1);
591  thisInnerBlockQ2 = thisIndexInOuterColumn / sizeInnerQuarterBlockQ2;
592  // thisTask % sizeInnerQuarterBlockQ2;
593  thisIndexInInnerBlockQ2 = thisTask&(sizeInnerQuarterBlockQ2-1);
594  thisInnerBlockQ1InInnerBlockQ2 = thisIndexInInnerBlockQ2 / sizeInnerHalfBlockQ1;
595  // thisTask % sizeInnerHalfBlockQ1;
596  thisIndexInInnerBlockQ1 = thisTask&(sizeInnerHalfBlockQ1-1);
597 
598  // get index in state vector corresponding to upper inner block
599  thisIndex = thisOuterColumn*sizeOuterColumn + thisInnerBlockQ2*sizeInnerBlockQ2
600  + thisInnerBlockQ1InInnerBlockQ2*sizeInnerBlockQ1 + thisIndexInInnerBlockQ1;
601 
602  // check if we are in the upper or lower half of an outer block for Q1
603  outerBitQ1 = extractBit(targetQubit, (thisIndex+qureg.numAmpsPerChunk*qureg.chunkId)>>qureg.numQubitsRepresented);
604  // if we are in the lower half of an outer block, shift to be in the lower half
605  // of the inner block as well (we want to dephase |0><0| and |1><1| only)
606  thisIndex += outerBitQ1*(sizeInnerHalfBlockQ1);
607 
608  // check if we are in the upper or lower half of an outer block for Q2
609  outerBitQ2 = extractBit(qubit2, (thisIndex+qureg.numAmpsPerChunk*qureg.chunkId)>>qureg.numQubitsRepresented);
610  // if we are in the lower half of an outer block, shift to be in the lower half
611  // of the inner block as well (we want to dephase |0><0| and |1><1| only)
612  thisIndex += outerBitQ2*(sizeInnerQuarterBlockQ2<<1);
613 
614  // NOTE: at this point thisIndex should be the index of the element we want to
615  // dephase in the chunk of the state vector on this process, in the
616  // density matrix representation.
617  // thisTask is the index of the pair element in pairStateVec
618 
619 
620  // state[thisIndex] = (1-depolLevel)*state[thisIndex] + depolLevel*(state[thisIndex]
621  // + pair[thisTask])/2
622  // NOTE: must set gamma=1 if using this function for steps 1 or 2
623  qureg.stateVec.real[thisIndex] = gamma*(qureg.stateVec.real[thisIndex] +
624  delta*qureg.pairStateVec.real[thisTask]);
625  qureg.stateVec.imag[thisIndex] = gamma*(qureg.stateVec.imag[thisIndex] +
626  delta*qureg.pairStateVec.imag[thisTask]);
627  }
628  }
629 }

References Qureg::chunkId, extractBit(), Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, Qureg::pairStateVec, and Qureg::stateVec.

Referenced by densmatr_mixTwoQubitDepolarising().

◆ densmatr_mixTwoQubitDepolarisingLocal()

void densmatr_mixTwoQubitDepolarisingLocal ( Qureg  qureg,
int  qubit1,
int  qubit2,
qreal  delta,
qreal  gamma 
)

Definition at line 385 of file QuEST_cpu.c.

385  {
386  const long long int numTasks = qureg.numAmpsPerChunk;
387  long long int innerMaskQubit1 = 1LL << qubit1;
388  long long int outerMaskQubit1= 1LL << (qubit1 + qureg.numQubitsRepresented);
389  long long int totMaskQubit1 = innerMaskQubit1 | outerMaskQubit1;
390  long long int innerMaskQubit2 = 1LL << qubit2;
391  long long int outerMaskQubit2 = 1LL << (qubit2 + qureg.numQubitsRepresented);
392  long long int totMaskQubit2 = innerMaskQubit2 | outerMaskQubit2;
393 
394  long long int thisTask;
395  long long int partner;
396  long long int thisPatternQubit1, thisPatternQubit2;
397 
398  qreal real00, imag00;
399 
400 # ifdef _OPENMP
401 # pragma omp parallel \
402  default (none) \
403  shared (totMaskQubit1,totMaskQubit2,qureg,delta,gamma) \
404  private (thisTask,partner,thisPatternQubit1,thisPatternQubit2,real00,imag00)
405 # endif
406  {
407 
408 # ifdef _OPENMP
409 # pragma omp for schedule (static)
410 # endif
411  //--------------------------------------- STEP ONE ---------------------
412  for (thisTask=0; thisTask<numTasks; thisTask++){
413  thisPatternQubit1 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit1;
414  thisPatternQubit2 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit2;
415  if ((thisPatternQubit1==0) && ((thisPatternQubit2==0)
416  || (thisPatternQubit2==totMaskQubit2))){
417  //this element of form |...X...0...><...X...0...| for X either 0 or 1.
418  partner = thisTask | totMaskQubit1;
419  real00 = qureg.stateVec.real[thisTask];
420  imag00 = qureg.stateVec.imag[thisTask];
421 
422  qureg.stateVec.real[thisTask] = qureg.stateVec.real[thisTask]
423  + delta*qureg.stateVec.real[partner];
424  qureg.stateVec.imag[thisTask] = qureg.stateVec.imag[thisTask]
425  + delta*qureg.stateVec.imag[partner];
426 
427  qureg.stateVec.real[partner] = qureg.stateVec.real[partner] + delta*real00;
428  qureg.stateVec.imag[partner] = qureg.stateVec.imag[partner] + delta*imag00;
429 
430  }
431  }
432 # ifdef _OPENMP
433 # pragma omp for schedule (static)
434 # endif
435  //--------------------------------------- STEP TWO ---------------------
436  for (thisTask=0; thisTask<numTasks; thisTask++){
437  thisPatternQubit1 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit1;
438  thisPatternQubit2 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit2;
439  if ((thisPatternQubit2==0) && ((thisPatternQubit1==0)
440  || (thisPatternQubit1==totMaskQubit1))){
441  //this element of form |...0...X...><...0...X...| for X either 0 or 1.
442  partner = thisTask | totMaskQubit2;
443  real00 = qureg.stateVec.real[thisTask];
444  imag00 = qureg.stateVec.imag[thisTask];
445 
446  qureg.stateVec.real[thisTask] = qureg.stateVec.real[thisTask]
447  + delta*qureg.stateVec.real[partner];
448  qureg.stateVec.imag[thisTask] = qureg.stateVec.imag[thisTask]
449  + delta*qureg.stateVec.imag[partner];
450 
451  qureg.stateVec.real[partner] = qureg.stateVec.real[partner] + delta*real00;
452  qureg.stateVec.imag[partner] = qureg.stateVec.imag[partner] + delta*imag00;
453 
454  }
455  }
456 
457 # ifdef _OPENMP
458 # pragma omp for schedule (static)
459 # endif
460  //--------------------------------------- STEP THREE ---------------------
461  for (thisTask=0; thisTask<numTasks; thisTask++){
462  thisPatternQubit1 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit1;
463  thisPatternQubit2 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit2;
464  if ((thisPatternQubit2==0) && ((thisPatternQubit1==0)
465  || (thisPatternQubit1==totMaskQubit1))){
466  //this element of form |...0...X...><...0...X...| for X either 0 or 1.
467  partner = thisTask | totMaskQubit2;
468  partner = partner ^ totMaskQubit1;
469  real00 = qureg.stateVec.real[thisTask];
470  imag00 = qureg.stateVec.imag[thisTask];
471 
472  qureg.stateVec.real[thisTask] = gamma * (qureg.stateVec.real[thisTask]
473  + delta*qureg.stateVec.real[partner]);
474  qureg.stateVec.imag[thisTask] = gamma * (qureg.stateVec.imag[thisTask]
475  + delta*qureg.stateVec.imag[partner]);
476 
477  qureg.stateVec.real[partner] = gamma * (qureg.stateVec.real[partner]
478  + delta*real00);
479  qureg.stateVec.imag[partner] = gamma * (qureg.stateVec.imag[partner]
480  + delta*imag00);
481 
482  }
483  }
484  }
485 }

References Qureg::chunkId, Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, qreal, and Qureg::stateVec.

Referenced by densmatr_mixTwoQubitDepolarising().

◆ densmatr_mixTwoQubitDepolarisingLocalPart1()

void densmatr_mixTwoQubitDepolarisingLocalPart1 ( Qureg  qureg,
int  qubit1,
int  qubit2,
qreal  delta 
)

Definition at line 487 of file QuEST_cpu.c.

487  {
488  const long long int numTasks = qureg.numAmpsPerChunk;
489  long long int innerMaskQubit1 = 1LL << qubit1;
490  long long int outerMaskQubit1= 1LL << (qubit1 + qureg.numQubitsRepresented);
491  long long int totMaskQubit1 = innerMaskQubit1 | outerMaskQubit1;
492  long long int innerMaskQubit2 = 1LL << qubit2;
493  long long int outerMaskQubit2 = 1LL << (qubit2 + qureg.numQubitsRepresented);
494  long long int totMaskQubit2 = innerMaskQubit2 | outerMaskQubit2;
495  // correct for being in a particular chunk
496  //totMaskQubit2 = totMaskQubit2&(qureg.numAmpsPerChunk-1); // totMaskQubit2 % numAmpsPerChunk
497 
498 
499  long long int thisTask;
500  long long int partner;
501  long long int thisPatternQubit1, thisPatternQubit2;
502 
503  qreal real00, imag00;
504 
505 # ifdef _OPENMP
506 # pragma omp parallel \
507  default (none) \
508  shared (totMaskQubit1,totMaskQubit2,qureg,delta) \
509  private (thisTask,partner,thisPatternQubit1,thisPatternQubit2,real00,imag00)
510 # endif
511  {
512 
513 # ifdef _OPENMP
514 # pragma omp for schedule (static)
515 # endif
516  //--------------------------------------- STEP ONE ---------------------
517  for (thisTask=0; thisTask<numTasks; thisTask ++){
518  thisPatternQubit1 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit1;
519  thisPatternQubit2 = (thisTask+qureg.numAmpsPerChunk*qureg.chunkId)&totMaskQubit2;
520  if ((thisPatternQubit1==0) && ((thisPatternQubit2==0)
521  || (thisPatternQubit2==totMaskQubit2))){
522  //this element of form |...X...0...><...X...0...| for X either 0 or 1.
523  partner = thisTask | totMaskQubit1;
524  real00 = qureg.stateVec.real[thisTask];
525  imag00 = qureg.stateVec.imag[thisTask];
526 
527  qureg.stateVec.real[thisTask] = qureg.stateVec.real[thisTask]
528  + delta*qureg.stateVec.real[partner];
529  qureg.stateVec.imag[thisTask] = qureg.stateVec.imag[thisTask]
530  + delta*qureg.stateVec.imag[partner];
531 
532  qureg.stateVec.real[partner] = qureg.stateVec.real[partner] + delta*real00;
533  qureg.stateVec.imag[partner] = qureg.stateVec.imag[partner] + delta*imag00;
534 
535  }
536  }
537  }
538 }

References Qureg::chunkId, Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, qreal, and Qureg::stateVec.

Referenced by densmatr_mixTwoQubitDepolarising().

◆ densmatr_mixTwoQubitDepolarisingQ1LocalQ2DistributedPart3()

void densmatr_mixTwoQubitDepolarisingQ1LocalQ2DistributedPart3 ( Qureg  qureg,
const int  targetQubit,
const int  qubit2,
qreal  delta,
qreal  gamma 
)

Definition at line 631 of file QuEST_cpu.c.

632  {
633 
634  long long int sizeInnerBlockQ1, sizeInnerHalfBlockQ1;
635  long long int sizeInnerBlockQ2, sizeInnerHalfBlockQ2, sizeInnerQuarterBlockQ2;
636  long long int sizeOuterColumn, sizeOuterQuarterColumn;
637  long long int thisInnerBlockQ2,
638  thisOuterColumn, // current column in density matrix
639  thisIndex, // current index in (density matrix representation) state vector
640  thisIndexInPairVector,
641  thisIndexInOuterColumn,
642  thisIndexInInnerBlockQ1,
643  thisIndexInInnerBlockQ2,
644  thisInnerBlockQ1InInnerBlockQ2;
645  int outerBitQ1, outerBitQ2;
646 
647  long long int thisTask;
648  const long long int numTasks=qureg.numAmpsPerChunk>>2;
649 
650  // set dimensions
651  sizeInnerHalfBlockQ1 = 1LL << targetQubit;
652  sizeInnerHalfBlockQ2 = 1LL << qubit2;
653  sizeInnerQuarterBlockQ2 = sizeInnerHalfBlockQ2 >> 1;
654  sizeInnerBlockQ2 = sizeInnerHalfBlockQ2 << 1;
655  sizeInnerBlockQ1 = 2LL * sizeInnerHalfBlockQ1;
656  sizeOuterColumn = 1LL << qureg.numQubitsRepresented;
657  sizeOuterQuarterColumn = sizeOuterColumn >> 2;
658 
659 //# if 0
660 # ifdef _OPENMP
661 # pragma omp parallel \
662  default (none) \
663  shared (sizeInnerBlockQ1,sizeInnerHalfBlockQ1,sizeInnerBlockQ2,sizeInnerHalfBlockQ2,sizeInnerQuarterBlockQ2,\
664  sizeOuterColumn,sizeOuterQuarterColumn,qureg,delta,gamma) \
665  private (thisTask,thisInnerBlockQ2,thisInnerBlockQ1InInnerBlockQ2, \
666  thisOuterColumn,thisIndex,thisIndexInPairVector,thisIndexInOuterColumn, \
667  thisIndexInInnerBlockQ1,thisIndexInInnerBlockQ2,outerBitQ1,outerBitQ2)
668 # endif
669  {
670 # ifdef _OPENMP
671 # pragma omp for schedule (static)
672 # endif
673 //# endif
674  // thisTask iterates over half the elements in this process' chunk of the density matrix
675  // treat this as iterating over all columns, then iterating over half the values
676  // within one column.
677  // If this function has been called, this process' chunk contains half an
678  // outer block or less
679  for (thisTask=0; thisTask<numTasks; thisTask++) {
680  // we want to process all columns in the density matrix,
681  // updating the values for half of each column (one half of each inner block)
682  thisOuterColumn = thisTask / sizeOuterQuarterColumn;
683  // thisTask % sizeOuterQuarterColumn
684  thisIndexInOuterColumn = thisTask&(sizeOuterQuarterColumn-1);
685  thisInnerBlockQ2 = thisIndexInOuterColumn / sizeInnerQuarterBlockQ2;
686  // thisTask % sizeInnerQuarterBlockQ2;
687  thisIndexInInnerBlockQ2 = thisTask&(sizeInnerQuarterBlockQ2-1);
688  thisInnerBlockQ1InInnerBlockQ2 = thisIndexInInnerBlockQ2 / sizeInnerHalfBlockQ1;
689  // thisTask % sizeInnerHalfBlockQ1;
690  thisIndexInInnerBlockQ1 = thisTask&(sizeInnerHalfBlockQ1-1);
691 
692  // get index in state vector corresponding to upper inner block
693  thisIndex = thisOuterColumn*sizeOuterColumn + thisInnerBlockQ2*sizeInnerBlockQ2
694  + thisInnerBlockQ1InInnerBlockQ2*sizeInnerBlockQ1 + thisIndexInInnerBlockQ1;
695 
696  // check if we are in the upper or lower half of an outer block for Q1
697  outerBitQ1 = extractBit(targetQubit, (thisIndex+qureg.numAmpsPerChunk*qureg.chunkId)>>qureg.numQubitsRepresented);
698  // if we are in the lower half of an outer block, shift to be in the lower half
699  // of the inner block as well (we want to dephase |0><0| and |1><1| only)
700  thisIndex += outerBitQ1*(sizeInnerHalfBlockQ1);
701 
702  // For part 3 we need to match elements such that (my Q1 != pair Q1) AND (my Q2 != pair Q2)
703  // Find correct index in pairStateVector
704  thisIndexInPairVector = thisTask + (1-outerBitQ1)*sizeInnerHalfBlockQ1*sizeOuterQuarterColumn -
705  outerBitQ1*sizeInnerHalfBlockQ1*sizeOuterQuarterColumn;
706 
707  // check if we are in the upper or lower half of an outer block for Q2
708  outerBitQ2 = extractBit(qubit2, (thisIndex+qureg.numAmpsPerChunk*qureg.chunkId)>>qureg.numQubitsRepresented);
709  // if we are in the lower half of an outer block, shift to be in the lower half
710  // of the inner block as well (we want to dephase |0><0| and |1><1| only)
711  thisIndex += outerBitQ2*(sizeInnerQuarterBlockQ2<<1);
712 
713 
714  // NOTE: at this point thisIndex should be the index of the element we want to
715  // dephase in the chunk of the state vector on this process, in the
716  // density matrix representation.
717 
718 
719  // state[thisIndex] = (1-depolLevel)*state[thisIndex] + depolLevel*(state[thisIndex]
720  // + pair[thisIndexInPairVector])/2
721  qureg.stateVec.real[thisIndex] = gamma*(qureg.stateVec.real[thisIndex] +
722  delta*qureg.pairStateVec.real[thisIndexInPairVector]);
723 
724  qureg.stateVec.imag[thisIndex] = gamma*(qureg.stateVec.imag[thisIndex] +
725  delta*qureg.pairStateVec.imag[thisIndexInPairVector]);
726  }
727  }
728 
729 }

References Qureg::chunkId, extractBit(), Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, Qureg::pairStateVec, and Qureg::stateVec.

Referenced by densmatr_mixTwoQubitDepolarising().

◆ extractBit()

static int extractBit ( const int  locationOfBitFromRight,
const long long int  theEncodedNumber 
)
inlinestatic

Definition at line 26 of file QuEST_cpu_internal.h.

26  {
27  return (theEncodedNumber & ( 1LL << locationOfBitFromRight )) >> locationOfBitFromRight;
28 }

Referenced by isOddParity().

◆ flipBit()

static long long int flipBit ( long long int  number,
int  bitInd 
)
inlinestatic

Definition at line 30 of file QuEST_cpu_internal.h.

30  {
31  return (number ^ (1LL << bitInd));
32 }

◆ insertTwoZeroBits()

static long long int insertTwoZeroBits ( long long int  number,
int  bit1,
int  bit2 
)
inlinestatic

Definition at line 49 of file QuEST_cpu_internal.h.

49  {
50  int small = (bit1 < bit2)? bit1 : bit2;
51  int big = (bit1 < bit2)? bit2 : bit1;
52  return insertZeroBit(insertZeroBit(number, small), big);
53 }

References insertZeroBit().

◆ insertZeroBit()

static long long int insertZeroBit ( long long int  number,
int  index 
)
inlinestatic

Definition at line 42 of file QuEST_cpu_internal.h.

42  {
43  long long int left, right;
44  left = (number >> index) << index;
45  right = number - left;
46  return (left << 1) ^ right;
47 }

Referenced by insertTwoZeroBits().

◆ isOddParity()

static int isOddParity ( long long int  number,
int  qb1,
int  qb2 
)
inlinestatic

Definition at line 38 of file QuEST_cpu_internal.h.

38  {
39  return extractBit(qb1, number) != extractBit(qb2, number);
40 }

References extractBit().

Referenced by statevec_swapQubitAmpsDistributed().

◆ maskContainsBit()

static int maskContainsBit ( long long int  mask,
int  bitInd 
)
inlinestatic

Definition at line 34 of file QuEST_cpu_internal.h.

34  {
35  return mask & (1LL << bitInd);
36 }

Referenced by statevec_multiControlledMultiQubitUnitary(), and statevec_multiControlledTwoQubitUnitary().

◆ statevec_calcInnerProductLocal()

Complex statevec_calcInnerProductLocal ( Qureg  bra,
Qureg  ket 
)

Definition at line 1070 of file QuEST_cpu.c.

1070  {
1071 
1072  qreal innerProdReal = 0;
1073  qreal innerProdImag = 0;
1074 
1075  long long int index;
1076  long long int numAmps = bra.numAmpsPerChunk;
1077  qreal *braVecReal = bra.stateVec.real;
1078  qreal *braVecImag = bra.stateVec.imag;
1079  qreal *ketVecReal = ket.stateVec.real;
1080  qreal *ketVecImag = ket.stateVec.imag;
1081 
1082  qreal braRe, braIm, ketRe, ketIm;
1083 
1084 # ifdef _OPENMP
1085 # pragma omp parallel \
1086  shared (braVecReal, braVecImag, ketVecReal, ketVecImag, numAmps) \
1087  private (index, braRe, braIm, ketRe, ketIm) \
1088  reduction ( +:innerProdReal, innerProdImag )
1089 # endif
1090  {
1091 # ifdef _OPENMP
1092 # pragma omp for schedule (static)
1093 # endif
1094  for (index=0; index < numAmps; index++) {
1095  braRe = braVecReal[index];
1096  braIm = braVecImag[index];
1097  ketRe = ketVecReal[index];
1098  ketIm = ketVecImag[index];
1099 
1100  // conj(bra_i) * ket_i
1101  innerProdReal += braRe*ketRe + braIm*ketIm;
1102  innerProdImag += braRe*ketIm - braIm*ketRe;
1103  }
1104  }
1105 
1106  Complex innerProd;
1107  innerProd.real = innerProdReal;
1108  innerProd.imag = innerProdImag;
1109  return innerProd;
1110 }

References Complex::imag, Qureg::numAmpsPerChunk, qreal, Complex::real, and Qureg::stateVec.

Referenced by statevec_calcInnerProduct().

◆ statevec_collapseToKnownProbOutcomeDistributedRenorm()

void statevec_collapseToKnownProbOutcomeDistributedRenorm ( Qureg  qureg,
const int  measureQubit,
const qreal  totalProbability 
)

Renormalise parts of the state vector where measureQubit=0 or 1, based on the total probability of that qubit being in state 0 or 1.

Measure in Zero performs an irreversible change to the state vector: it updates the vector according to the event that the value 'outcome' has been measured on the qubit indicated by measureQubit (where this label starts from 0, of course). It achieves this by setting all inconsistent amplitudes to 0 and then renormalising based on the total probability of measuring measureQubit=0 if outcome=0 and measureQubit=1 if outcome=1. In the distributed version, one block (with measureQubit=0 in the first half of the block and measureQubit=1 in the second half of the block) is spread over multiple chunks, meaning that each chunks performs only renormalisation or only setting amplitudes to 0. This function handles the renormalisation.

Parameters
[in,out]quregobject representing the set of qubits
[in]measureQubitqubit to measure
[in]totalProbabilityprobability of qubit measureQubit being zero

Definition at line 3422 of file QuEST_cpu.c.

3423 {
3424  // ----- temp variables
3425  long long int thisTask;
3426  long long int numTasks=qureg.numAmpsPerChunk;
3427 
3428  qreal renorm=1/sqrt(totalProbability);
3429 
3430  qreal *stateVecReal = qureg.stateVec.real;
3431  qreal *stateVecImag = qureg.stateVec.imag;
3432 
3433 # ifdef _OPENMP
3434 # pragma omp parallel \
3435  shared (numTasks,stateVecReal,stateVecImag) \
3436  private (thisTask)
3437 # endif
3438  {
3439 # ifdef _OPENMP
3440 # pragma omp for schedule (static)
3441 # endif
3442  for (thisTask=0; thisTask<numTasks; thisTask++) {
3443  stateVecReal[thisTask] = stateVecReal[thisTask]*renorm;
3444  stateVecImag[thisTask] = stateVecImag[thisTask]*renorm;
3445  }
3446  }
3447 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_collapseToKnownProbOutcome().

◆ statevec_collapseToKnownProbOutcomeLocal()

void statevec_collapseToKnownProbOutcomeLocal ( Qureg  qureg,
int  measureQubit,
int  outcome,
qreal  totalProbability 
)

Update the state vector to be consistent with measuring measureQubit=0 if outcome=0 and measureQubit=1 if outcome=1.

Performs an irreversible change to the state vector: it updates the vector according to the event that an outcome have been measured on the qubit indicated by measureQubit (where this label starts from 0, of course). It achieves this by setting all inconsistent amplitudes to 0 and then renormalising based on the total probability of measuring measureQubit=0 or 1 according to the value of outcome. In the local version, one or more blocks (with measureQubit=0 in the first half of the block and measureQubit=1 in the second half of the block) fit entirely into one chunk.

Parameters
[in,out]quregobject representing the set of qubits
[in]measureQubitqubit to measure
[in]totalProbabilityprobability of qubit measureQubit being either zero or one
[in]outcometo measure the probability of and set the state to – either zero or one

Definition at line 3340 of file QuEST_cpu.c.

3341 {
3342  // ----- sizes
3343  long long int sizeBlock, // size of blocks
3344  sizeHalfBlock; // size of blocks halved
3345  // ----- indices
3346  long long int thisBlock, // current block
3347  index; // current index for first half block
3348  // ----- measured probability
3349  qreal renorm; // probability (returned) value
3350  // ----- temp variables
3351  long long int thisTask; // task based approach for expose loop with small granularity
3352  // (good for shared memory parallelism)
3353  long long int numTasks=qureg.numAmpsPerChunk>>1;
3354 
3355  // ---------------------------------------------------------------- //
3356  // dimensions //
3357  // ---------------------------------------------------------------- //
3358  sizeHalfBlock = 1LL << (measureQubit); // number of state vector elements to sum,
3359  // and then the number to skip
3360  sizeBlock = 2LL * sizeHalfBlock; // size of blocks (pairs of measure and skip entries)
3361 
3362  renorm=1/sqrt(totalProbability);
3363  qreal *stateVecReal = qureg.stateVec.real;
3364  qreal *stateVecImag = qureg.stateVec.imag;
3365 
3366 
3367 # ifdef _OPENMP
3368 # pragma omp parallel \
3369  default (none) \
3370  shared (numTasks,sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag,renorm,outcome) \
3371  private (thisTask,thisBlock,index)
3372 # endif
3373  {
3374  if (outcome==0){
3375  // measure qubit is 0
3376 # ifdef _OPENMP
3377 # pragma omp for schedule (static)
3378 # endif
3379  for (thisTask=0; thisTask<numTasks; thisTask++) {
3380  thisBlock = thisTask / sizeHalfBlock;
3381  index = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
3382  stateVecReal[index]=stateVecReal[index]*renorm;
3383  stateVecImag[index]=stateVecImag[index]*renorm;
3384 
3385  stateVecReal[index+sizeHalfBlock]=0;
3386  stateVecImag[index+sizeHalfBlock]=0;
3387  }
3388  } else {
3389  // measure qubit is 1
3390 # ifdef _OPENMP
3391 # pragma omp for schedule (static)
3392 # endif
3393  for (thisTask=0; thisTask<numTasks; thisTask++) {
3394  thisBlock = thisTask / sizeHalfBlock;
3395  index = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
3396  stateVecReal[index]=0;
3397  stateVecImag[index]=0;
3398 
3399  stateVecReal[index+sizeHalfBlock]=stateVecReal[index+sizeHalfBlock]*renorm;
3400  stateVecImag[index+sizeHalfBlock]=stateVecImag[index+sizeHalfBlock]*renorm;
3401  }
3402  }
3403  }
3404 
3405 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_collapseToKnownProbOutcome().

◆ statevec_collapseToOutcomeDistributedSetZero()

void statevec_collapseToOutcomeDistributedSetZero ( Qureg  qureg)

Set all amplitudes in one chunk to 0.

Measure in Zero performs an irreversible change to the state vector: it updates the vector according to the event that a zero have been measured on the qubit indicated by measureQubit (where this label starts from 0, of course). It achieves this by setting all inconsistent amplitudes to 0 and then renormalising based on the total probability of measuring measureQubit=0 or 1. In the distributed version, one block (with measureQubit=0 in the first half of the block and measureQubit=1 in the second half of the block) is spread over multiple chunks, meaning that each chunks performs only renormalisation or only setting amplitudes to 0. This function handles setting amplitudes to 0.

Parameters
[in,out]quregobject representing the set of qubits
[in]measureQubitqubit to measure

Definition at line 3461 of file QuEST_cpu.c.

3462 {
3463  // ----- temp variables
3464  long long int thisTask;
3465  long long int numTasks=qureg.numAmpsPerChunk;
3466 
3467  // ---------------------------------------------------------------- //
3468  // find probability //
3469  // ---------------------------------------------------------------- //
3470 
3471  qreal *stateVecReal = qureg.stateVec.real;
3472  qreal *stateVecImag = qureg.stateVec.imag;
3473 
3474 # ifdef _OPENMP
3475 # pragma omp parallel \
3476  shared (numTasks,stateVecReal,stateVecImag) \
3477  private (thisTask)
3478 # endif
3479  {
3480 # ifdef _OPENMP
3481 # pragma omp for schedule (static)
3482 # endif
3483  for (thisTask=0; thisTask<numTasks; thisTask++) {
3484  stateVecReal[thisTask] = 0;
3485  stateVecImag[thisTask] = 0;
3486  }
3487  }
3488 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_collapseToKnownProbOutcome().

◆ statevec_compactUnitaryDistributed()

void statevec_compactUnitaryDistributed ( Qureg  qureg,
Complex  rot1,
Complex  rot2,
ComplexArray  stateVecUp,
ComplexArray  stateVecLo,
ComplexArray  stateVecOut 
)

Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha and beta, and a subset of the state vector with upper and lower block values stored seperately.

Parameters
[in,out]quregobject representing the set of qubits
[in]rot1rotation angle
[in]rot2rotation angle
[in]stateVecUpprobability amplitudes in upper half of a block
[in]stateVecLoprobability amplitudes in lower half of a block
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 1969 of file QuEST_cpu.c.

1974 {
1975 
1976  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
1977  long long int thisTask;
1978  const long long int numTasks=qureg.numAmpsPerChunk;
1979 
1980  qreal rot1Real=rot1.real, rot1Imag=rot1.imag;
1981  qreal rot2Real=rot2.real, rot2Imag=rot2.imag;
1982  qreal *stateVecRealUp=stateVecUp.real, *stateVecImagUp=stateVecUp.imag;
1983  qreal *stateVecRealLo=stateVecLo.real, *stateVecImagLo=stateVecLo.imag;
1984  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
1985 
1986 # ifdef _OPENMP
1987 # pragma omp parallel \
1988  default (none) \
1989  shared (stateVecRealUp,stateVecImagUp,stateVecRealLo,stateVecImagLo,stateVecRealOut,stateVecImagOut, \
1990  rot1Real,rot1Imag, rot2Real,rot2Imag) \
1991  private (thisTask,stateRealUp,stateImagUp,stateRealLo,stateImagLo)
1992 # endif
1993  {
1994 # ifdef _OPENMP
1995 # pragma omp for schedule (static)
1996 # endif
1997  for (thisTask=0; thisTask<numTasks; thisTask++) {
1998  // store current state vector values in temp variables
1999  stateRealUp = stateVecRealUp[thisTask];
2000  stateImagUp = stateVecImagUp[thisTask];
2001 
2002  stateRealLo = stateVecRealLo[thisTask];
2003  stateImagLo = stateVecImagLo[thisTask];
2004 
2005  // state[indexUp] = alpha * state[indexUp] - conj(beta) * state[indexLo]
2006  stateVecRealOut[thisTask] = rot1Real*stateRealUp - rot1Imag*stateImagUp + rot2Real*stateRealLo + rot2Imag*stateImagLo;
2007  stateVecImagOut[thisTask] = rot1Real*stateImagUp + rot1Imag*stateRealUp + rot2Real*stateImagLo - rot2Imag*stateRealLo;
2008  }
2009  }
2010 }

References Complex::imag, Qureg::numAmpsPerChunk, qreal, and Complex::real.

Referenced by statevec_compactUnitary().

◆ statevec_compactUnitaryLocal()

void statevec_compactUnitaryLocal ( Qureg  qureg,
const int  targetQubit,
Complex  alpha,
Complex  beta 
)

Definition at line 1656 of file QuEST_cpu.c.

1657 {
1658  long long int sizeBlock, sizeHalfBlock;
1659  long long int thisBlock, // current block
1660  indexUp,indexLo; // current index and corresponding index in lower half block
1661 
1662  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
1663  long long int thisTask;
1664  const long long int numTasks=qureg.numAmpsPerChunk>>1;
1665 
1666  // set dimensions
1667  sizeHalfBlock = 1LL << targetQubit;
1668  sizeBlock = 2LL * sizeHalfBlock;
1669 
1670  // Can't use qureg.stateVec as a private OMP var
1671  qreal *stateVecReal = qureg.stateVec.real;
1672  qreal *stateVecImag = qureg.stateVec.imag;
1673  qreal alphaImag=alpha.imag, alphaReal=alpha.real;
1674  qreal betaImag=beta.imag, betaReal=beta.real;
1675 
1676 # ifdef _OPENMP
1677 # pragma omp parallel \
1678  default (none) \
1679  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag, alphaReal,alphaImag, betaReal,betaImag) \
1680  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp,stateRealLo,stateImagLo)
1681 # endif
1682  {
1683 # ifdef _OPENMP
1684 # pragma omp for schedule (static)
1685 # endif
1686  for (thisTask=0; thisTask<numTasks; thisTask++) {
1687 
1688  thisBlock = thisTask / sizeHalfBlock;
1689  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
1690  indexLo = indexUp + sizeHalfBlock;
1691 
1692  // store current state vector values in temp variables
1693  stateRealUp = stateVecReal[indexUp];
1694  stateImagUp = stateVecImag[indexUp];
1695 
1696  stateRealLo = stateVecReal[indexLo];
1697  stateImagLo = stateVecImag[indexLo];
1698 
1699  // state[indexUp] = alpha * state[indexUp] - conj(beta) * state[indexLo]
1700  stateVecReal[indexUp] = alphaReal*stateRealUp - alphaImag*stateImagUp
1701  - betaReal*stateRealLo - betaImag*stateImagLo;
1702  stateVecImag[indexUp] = alphaReal*stateImagUp + alphaImag*stateRealUp
1703  - betaReal*stateImagLo + betaImag*stateRealLo;
1704 
1705  // state[indexLo] = beta * state[indexUp] + conj(alpha) * state[indexLo]
1706  stateVecReal[indexLo] = betaReal*stateRealUp - betaImag*stateImagUp
1707  + alphaReal*stateRealLo + alphaImag*stateImagLo;
1708  stateVecImag[indexLo] = betaReal*stateImagUp + betaImag*stateRealUp
1709  + alphaReal*stateImagLo - alphaImag*stateRealLo;
1710  }
1711  }
1712 
1713 }

References Complex::imag, Qureg::numAmpsPerChunk, qreal, Complex::real, and Qureg::stateVec.

Referenced by statevec_compactUnitary().

◆ statevec_controlledCompactUnitaryDistributed()

void statevec_controlledCompactUnitaryDistributed ( Qureg  qureg,
const int  controlQubit,
Complex  rot1,
Complex  rot2,
ComplexArray  stateVecUp,
ComplexArray  stateVecLo,
ComplexArray  stateVecOut 
)

Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha and beta and a subset of the state vector with upper and lower block values stored seperately.

Only perform the rotation where the control qubit is one.

Parameters
[in,out]quregobject representing the set of qubits
[in]controlQubitqubit to determine whether or not to perform a rotation
[in]rot1rotation angle
[in]rot2rotation angle
[in]stateVecUpprobability amplitudes in upper half of a block
[in]stateVecLoprobability amplitudes in lower half of a block
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2285 of file QuEST_cpu.c.

2290 {
2291 
2292  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2293  long long int thisTask;
2294  const long long int numTasks=qureg.numAmpsPerChunk;
2295  const long long int chunkSize=qureg.numAmpsPerChunk;
2296  const long long int chunkId=qureg.chunkId;
2297 
2298  int controlBit;
2299 
2300  qreal rot1Real=rot1.real, rot1Imag=rot1.imag;
2301  qreal rot2Real=rot2.real, rot2Imag=rot2.imag;
2302  qreal *stateVecRealUp=stateVecUp.real, *stateVecImagUp=stateVecUp.imag;
2303  qreal *stateVecRealLo=stateVecLo.real, *stateVecImagLo=stateVecLo.imag;
2304  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2305 
2306 # ifdef _OPENMP
2307 # pragma omp parallel \
2308  default (none) \
2309  shared (stateVecRealUp,stateVecImagUp,stateVecRealLo,stateVecImagLo,stateVecRealOut,stateVecImagOut, \
2310  rot1Real,rot1Imag, rot2Real,rot2Imag) \
2311  private (thisTask,stateRealUp,stateImagUp,stateRealLo,stateImagLo,controlBit)
2312 # endif
2313  {
2314 # ifdef _OPENMP
2315 # pragma omp for schedule (static)
2316 # endif
2317  for (thisTask=0; thisTask<numTasks; thisTask++) {
2318  controlBit = extractBit (controlQubit, thisTask+chunkId*chunkSize);
2319  if (controlBit){
2320  // store current state vector values in temp variables
2321  stateRealUp = stateVecRealUp[thisTask];
2322  stateImagUp = stateVecImagUp[thisTask];
2323 
2324  stateRealLo = stateVecRealLo[thisTask];
2325  stateImagLo = stateVecImagLo[thisTask];
2326 
2327  // state[indexUp] = alpha * state[indexUp] - conj(beta) * state[indexLo]
2328  stateVecRealOut[thisTask] = rot1Real*stateRealUp - rot1Imag*stateImagUp + rot2Real*stateRealLo + rot2Imag*stateImagLo;
2329  stateVecImagOut[thisTask] = rot1Real*stateImagUp + rot1Imag*stateRealUp + rot2Real*stateImagLo - rot2Imag*stateRealLo;
2330  }
2331  }
2332  }
2333 }

References Qureg::chunkId, extractBit(), Complex::imag, Qureg::numAmpsPerChunk, qreal, and Complex::real.

Referenced by statevec_controlledCompactUnitary().

◆ statevec_controlledCompactUnitaryLocal()

void statevec_controlledCompactUnitaryLocal ( Qureg  qureg,
const int  controlQubit,
const int  targetQubit,
Complex  alpha,
Complex  beta 
)

Definition at line 2069 of file QuEST_cpu.c.

2071 {
2072  long long int sizeBlock, sizeHalfBlock;
2073  long long int thisBlock, // current block
2074  indexUp,indexLo; // current index and corresponding index in lower half block
2075 
2076  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2077  long long int thisTask;
2078  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2079  const long long int chunkSize=qureg.numAmpsPerChunk;
2080  const long long int chunkId=qureg.chunkId;
2081 
2082  int controlBit;
2083 
2084  // set dimensions
2085  sizeHalfBlock = 1LL << targetQubit;
2086  sizeBlock = 2LL * sizeHalfBlock;
2087 
2088  // Can't use qureg.stateVec as a private OMP var
2089  qreal *stateVecReal = qureg.stateVec.real;
2090  qreal *stateVecImag = qureg.stateVec.imag;
2091  qreal alphaImag=alpha.imag, alphaReal=alpha.real;
2092  qreal betaImag=beta.imag, betaReal=beta.real;
2093 
2094 # ifdef _OPENMP
2095 # pragma omp parallel \
2096  default (none) \
2097  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag, alphaReal,alphaImag, betaReal,betaImag) \
2098  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp,stateRealLo,stateImagLo,controlBit)
2099 # endif
2100  {
2101 # ifdef _OPENMP
2102 # pragma omp for schedule (static)
2103 # endif
2104  for (thisTask=0; thisTask<numTasks; thisTask++) {
2105 
2106  thisBlock = thisTask / sizeHalfBlock;
2107  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2108  indexLo = indexUp + sizeHalfBlock;
2109 
2110  controlBit = extractBit (controlQubit, indexUp+chunkId*chunkSize);
2111  if (controlBit){
2112  // store current state vector values in temp variables
2113  stateRealUp = stateVecReal[indexUp];
2114  stateImagUp = stateVecImag[indexUp];
2115 
2116  stateRealLo = stateVecReal[indexLo];
2117  stateImagLo = stateVecImag[indexLo];
2118 
2119  // state[indexUp] = alpha * state[indexUp] - conj(beta) * state[indexLo]
2120  stateVecReal[indexUp] = alphaReal*stateRealUp - alphaImag*stateImagUp
2121  - betaReal*stateRealLo - betaImag*stateImagLo;
2122  stateVecImag[indexUp] = alphaReal*stateImagUp + alphaImag*stateRealUp
2123  - betaReal*stateImagLo + betaImag*stateRealLo;
2124 
2125  // state[indexLo] = beta * state[indexUp] + conj(alpha) * state[indexLo]
2126  stateVecReal[indexLo] = betaReal*stateRealUp - betaImag*stateImagUp
2127  + alphaReal*stateRealLo + alphaImag*stateImagLo;
2128  stateVecImag[indexLo] = betaReal*stateImagUp + betaImag*stateRealUp
2129  + alphaReal*stateImagLo - alphaImag*stateRealLo;
2130  }
2131  }
2132  }
2133 
2134 }

References Qureg::chunkId, extractBit(), Complex::imag, Qureg::numAmpsPerChunk, qreal, Complex::real, and Qureg::stateVec.

Referenced by statevec_controlledCompactUnitary().

◆ statevec_controlledNotDistributed()

void statevec_controlledNotDistributed ( Qureg  qureg,
const int  controlQubit,
ComplexArray  stateVecIn,
ComplexArray  stateVecOut 
)

Rotate a single qubit by {{0,1},{1,0}.

Operate on a subset of the state vector with upper and lower block values stored seperately. This rotation is just swapping upper and lower values, and stateVecIn must already be the correct section for this chunk. Only perform the rotation for elements where controlQubit is one.

Parameters
[in,out]quregobject representing the set of qubits
[in]stateVecInprobability amplitudes in lower or upper half of a block depending on chunkId
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2612 of file QuEST_cpu.c.

2615 {
2616 
2617  long long int thisTask;
2618  const long long int numTasks=qureg.numAmpsPerChunk;
2619  const long long int chunkSize=qureg.numAmpsPerChunk;
2620  const long long int chunkId=qureg.chunkId;
2621 
2622  int controlBit;
2623 
2624  qreal *stateVecRealIn=stateVecIn.real, *stateVecImagIn=stateVecIn.imag;
2625  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2626 
2627 # ifdef _OPENMP
2628 # pragma omp parallel \
2629  default (none) \
2630  shared (stateVecRealIn,stateVecImagIn,stateVecRealOut,stateVecImagOut) \
2631  private (thisTask,controlBit)
2632 # endif
2633  {
2634 # ifdef _OPENMP
2635 # pragma omp for schedule (static)
2636 # endif
2637  for (thisTask=0; thisTask<numTasks; thisTask++) {
2638  controlBit = extractBit (controlQubit, thisTask+chunkId*chunkSize);
2639  if (controlBit){
2640  stateVecRealOut[thisTask] = stateVecRealIn[thisTask];
2641  stateVecImagOut[thisTask] = stateVecImagIn[thisTask];
2642  }
2643  }
2644  }
2645 }

References Qureg::chunkId, extractBit(), Qureg::numAmpsPerChunk, and qreal.

Referenced by statevec_controlledNot().

◆ statevec_controlledNotLocal()

void statevec_controlledNotLocal ( Qureg  qureg,
const int  controlQubit,
const int  targetQubit 
)

Definition at line 2550 of file QuEST_cpu.c.

2551 {
2552  long long int sizeBlock, sizeHalfBlock;
2553  long long int thisBlock, // current block
2554  indexUp,indexLo; // current index and corresponding index in lower half block
2555 
2556  qreal stateRealUp,stateImagUp;
2557  long long int thisTask;
2558  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2559  const long long int chunkSize=qureg.numAmpsPerChunk;
2560  const long long int chunkId=qureg.chunkId;
2561 
2562  int controlBit;
2563 
2564  // set dimensions
2565  sizeHalfBlock = 1LL << targetQubit;
2566  sizeBlock = 2LL * sizeHalfBlock;
2567 
2568  // Can't use qureg.stateVec as a private OMP var
2569  qreal *stateVecReal = qureg.stateVec.real;
2570  qreal *stateVecImag = qureg.stateVec.imag;
2571 
2572 # ifdef _OPENMP
2573 # pragma omp parallel \
2574  default (none) \
2575  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag) \
2576  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp,controlBit)
2577 # endif
2578  {
2579 # ifdef _OPENMP
2580 # pragma omp for schedule (static)
2581 # endif
2582  for (thisTask=0; thisTask<numTasks; thisTask++) {
2583  thisBlock = thisTask / sizeHalfBlock;
2584  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2585  indexLo = indexUp + sizeHalfBlock;
2586 
2587  controlBit = extractBit(controlQubit, indexUp+chunkId*chunkSize);
2588  if (controlBit){
2589  stateRealUp = stateVecReal[indexUp];
2590  stateImagUp = stateVecImag[indexUp];
2591 
2592  stateVecReal[indexUp] = stateVecReal[indexLo];
2593  stateVecImag[indexUp] = stateVecImag[indexLo];
2594 
2595  stateVecReal[indexLo] = stateRealUp;
2596  stateVecImag[indexLo] = stateImagUp;
2597  }
2598  }
2599  }
2600 }

References Qureg::chunkId, extractBit(), Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_controlledNot().

◆ statevec_controlledPauliYDistributed()

void statevec_controlledPauliYDistributed ( Qureg  qureg,
const int  controlQubit,
ComplexArray  stateVecIn,
ComplexArray  stateVecOut,
const int  conjFactor 
)

Definition at line 2793 of file QuEST_cpu.c.

2796 {
2797 
2798  long long int thisTask;
2799  const long long int numTasks=qureg.numAmpsPerChunk;
2800  const long long int chunkSize=qureg.numAmpsPerChunk;
2801  const long long int chunkId=qureg.chunkId;
2802 
2803  int controlBit;
2804 
2805  qreal *stateVecRealIn=stateVecIn.real, *stateVecImagIn=stateVecIn.imag;
2806  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2807 
2808 # ifdef _OPENMP
2809 # pragma omp parallel \
2810  default (none) \
2811  shared (stateVecRealIn,stateVecImagIn,stateVecRealOut,stateVecImagOut) \
2812  private (thisTask,controlBit)
2813 # endif
2814  {
2815 # ifdef _OPENMP
2816 # pragma omp for schedule (static)
2817 # endif
2818  for (thisTask=0; thisTask<numTasks; thisTask++) {
2819  controlBit = extractBit (controlQubit, thisTask+chunkId*chunkSize);
2820  if (controlBit){
2821  stateVecRealOut[thisTask] = conjFac * stateVecImagIn[thisTask];
2822  stateVecImagOut[thisTask] = conjFac * -stateVecRealIn[thisTask];
2823  }
2824  }
2825  }
2826 }

References Qureg::chunkId, extractBit(), Qureg::numAmpsPerChunk, and qreal.

Referenced by statevec_controlledPauliY(), and statevec_controlledPauliYConj().

◆ statevec_controlledPauliYLocal()

void statevec_controlledPauliYLocal ( Qureg  qureg,
const int  controlQubit,
const int  targetQubit,
const int  conjFactor 
)

Definition at line 2740 of file QuEST_cpu.c.

2741 {
2742  long long int sizeBlock, sizeHalfBlock;
2743  long long int thisBlock, // current block
2744  indexUp,indexLo; // current index and corresponding index in lower half block
2745 
2746  qreal stateRealUp,stateImagUp;
2747  long long int thisTask;
2748  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2749  const long long int chunkSize=qureg.numAmpsPerChunk;
2750  const long long int chunkId=qureg.chunkId;
2751 
2752  int controlBit;
2753 
2754  // set dimensions
2755  sizeHalfBlock = 1LL << targetQubit;
2756  sizeBlock = 2LL * sizeHalfBlock;
2757 
2758  // Can't use qureg.stateVec as a private OMP var
2759  qreal *stateVecReal = qureg.stateVec.real;
2760  qreal *stateVecImag = qureg.stateVec.imag;
2761 
2762 # ifdef _OPENMP
2763 # pragma omp parallel \
2764  default (none) \
2765  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag) \
2766  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp,controlBit)
2767 # endif
2768  {
2769 # ifdef _OPENMP
2770 # pragma omp for schedule (static)
2771 # endif
2772  for (thisTask=0; thisTask<numTasks; thisTask++) {
2773  thisBlock = thisTask / sizeHalfBlock;
2774  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2775  indexLo = indexUp + sizeHalfBlock;
2776 
2777  controlBit = extractBit(controlQubit, indexUp+chunkId*chunkSize);
2778  if (controlBit){
2779  stateRealUp = stateVecReal[indexUp];
2780  stateImagUp = stateVecImag[indexUp];
2781 
2782  // update under +-{{0, -i}, {i, 0}}
2783  stateVecReal[indexUp] = conjFac * stateVecImag[indexLo];
2784  stateVecImag[indexUp] = conjFac * -stateVecReal[indexLo];
2785  stateVecReal[indexLo] = conjFac * -stateImagUp;
2786  stateVecImag[indexLo] = conjFac * stateRealUp;
2787  }
2788  }
2789  }
2790 }

References Qureg::chunkId, extractBit(), Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_controlledPauliY(), and statevec_controlledPauliYConj().

◆ statevec_controlledUnitaryDistributed()

void statevec_controlledUnitaryDistributed ( Qureg  qureg,
const int  controlQubit,
Complex  rot1,
Complex  rot2,
ComplexArray  stateVecUp,
ComplexArray  stateVecLo,
ComplexArray  stateVecOut 
)

Rotate a single qubit in the state vector of probability amplitudes, given two complex numbers alpha and beta and a subset of the state vector with upper and lower block values stored seperately.

Only perform the rotation where the control qubit is one.

Parameters
[in,out]quregobject representing the set of qubits
[in]controlQubitqubit to determine whether or not to perform a rotation
[in]rot1rotation angle
[in]rot2rotation angle
[in]stateVecUpprobability amplitudes in upper half of a block
[in]stateVecLoprobability amplitudes in lower half of a block
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2347 of file QuEST_cpu.c.

2352 {
2353 
2354  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2355  long long int thisTask;
2356  const long long int numTasks=qureg.numAmpsPerChunk;
2357  const long long int chunkSize=qureg.numAmpsPerChunk;
2358  const long long int chunkId=qureg.chunkId;
2359 
2360  int controlBit;
2361 
2362  qreal rot1Real=rot1.real, rot1Imag=rot1.imag;
2363  qreal rot2Real=rot2.real, rot2Imag=rot2.imag;
2364  qreal *stateVecRealUp=stateVecUp.real, *stateVecImagUp=stateVecUp.imag;
2365  qreal *stateVecRealLo=stateVecLo.real, *stateVecImagLo=stateVecLo.imag;
2366  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2367 
2368 # ifdef _OPENMP
2369 # pragma omp parallel \
2370  default (none) \
2371  shared (stateVecRealUp,stateVecImagUp,stateVecRealLo,stateVecImagLo,stateVecRealOut,stateVecImagOut, \
2372  rot1Real,rot1Imag, rot2Real,rot2Imag) \
2373  private (thisTask,stateRealUp,stateImagUp,stateRealLo,stateImagLo,controlBit)
2374 # endif
2375  {
2376 # ifdef _OPENMP
2377 # pragma omp for schedule (static)
2378 # endif
2379  for (thisTask=0; thisTask<numTasks; thisTask++) {
2380  controlBit = extractBit (controlQubit, thisTask+chunkId*chunkSize);
2381  if (controlBit){
2382  // store current state vector values in temp variables
2383  stateRealUp = stateVecRealUp[thisTask];
2384  stateImagUp = stateVecImagUp[thisTask];
2385 
2386  stateRealLo = stateVecRealLo[thisTask];
2387  stateImagLo = stateVecImagLo[thisTask];
2388 
2389  stateVecRealOut[thisTask] = rot1Real*stateRealUp - rot1Imag*stateImagUp
2390  + rot2Real*stateRealLo - rot2Imag*stateImagLo;
2391  stateVecImagOut[thisTask] = rot1Real*stateImagUp + rot1Imag*stateRealUp
2392  + rot2Real*stateImagLo + rot2Imag*stateRealLo;
2393  }
2394  }
2395  }
2396 }

References Qureg::chunkId, extractBit(), Complex::imag, Qureg::numAmpsPerChunk, qreal, and Complex::real.

Referenced by statevec_controlledUnitary().

◆ statevec_controlledUnitaryLocal()

void statevec_controlledUnitaryLocal ( Qureg  qureg,
const int  controlQubit,
const int  targetQubit,
ComplexMatrix2  u 
)

Definition at line 2207 of file QuEST_cpu.c.

2209 {
2210  long long int sizeBlock, sizeHalfBlock;
2211  long long int thisBlock, // current block
2212  indexUp,indexLo; // current index and corresponding index in lower half block
2213 
2214  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2215  long long int thisTask;
2216  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2217  const long long int chunkSize=qureg.numAmpsPerChunk;
2218  const long long int chunkId=qureg.chunkId;
2219 
2220  int controlBit;
2221 
2222  // set dimensions
2223  sizeHalfBlock = 1LL << targetQubit;
2224  sizeBlock = 2LL * sizeHalfBlock;
2225 
2226  // Can't use qureg.stateVec as a private OMP var
2227  qreal *stateVecReal = qureg.stateVec.real;
2228  qreal *stateVecImag = qureg.stateVec.imag;
2229 
2230 # ifdef _OPENMP
2231 # pragma omp parallel \
2232  default (none) \
2233  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag, u) \
2234  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp,stateRealLo,stateImagLo,controlBit)
2235 # endif
2236  {
2237 # ifdef _OPENMP
2238 # pragma omp for schedule (static)
2239 # endif
2240  for (thisTask=0; thisTask<numTasks; thisTask++) {
2241 
2242  thisBlock = thisTask / sizeHalfBlock;
2243  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2244  indexLo = indexUp + sizeHalfBlock;
2245 
2246  controlBit = extractBit (controlQubit, indexUp+chunkId*chunkSize);
2247  if (controlBit){
2248  // store current state vector values in temp variables
2249  stateRealUp = stateVecReal[indexUp];
2250  stateImagUp = stateVecImag[indexUp];
2251 
2252  stateRealLo = stateVecReal[indexLo];
2253  stateImagLo = stateVecImag[indexLo];
2254 
2255 
2256  // state[indexUp] = u00 * state[indexUp] + u01 * state[indexLo]
2257  stateVecReal[indexUp] = u.real[0][0]*stateRealUp - u.imag[0][0]*stateImagUp
2258  + u.real[0][1]*stateRealLo - u.imag[0][1]*stateImagLo;
2259  stateVecImag[indexUp] = u.real[0][0]*stateImagUp + u.imag[0][0]*stateRealUp
2260  + u.real[0][1]*stateImagLo + u.imag[0][1]*stateRealLo;
2261 
2262  // state[indexLo] = u10 * state[indexUp] + u11 * state[indexLo]
2263  stateVecReal[indexLo] = u.real[1][0]*stateRealUp - u.imag[1][0]*stateImagUp
2264  + u.real[1][1]*stateRealLo - u.imag[1][1]*stateImagLo;
2265  stateVecImag[indexLo] = u.real[1][0]*stateImagUp + u.imag[1][0]*stateRealUp
2266  + u.real[1][1]*stateImagLo + u.imag[1][1]*stateRealLo;
2267  }
2268  }
2269  }
2270 
2271 }

References Qureg::chunkId, extractBit(), ComplexMatrix2::imag, Qureg::numAmpsPerChunk, qreal, ComplexMatrix2::real, and Qureg::stateVec.

Referenced by statevec_controlledUnitary().

◆ statevec_findProbabilityOfZeroDistributed()

qreal statevec_findProbabilityOfZeroDistributed ( Qureg  qureg)

Measure the probability of a specified qubit being in the zero state across all amplitudes held in this chunk.

Size of regions to skip is a multiple of chunkSize. The results are communicated and aggregated by the caller

Parameters
[in]quregobject representing the set of qubits
Returns
probability of qubit measureQubit being zero

Definition at line 3222 of file QuEST_cpu.c.

3222  {
3223  // ----- measured probability
3224  qreal totalProbability; // probability (returned) value
3225  // ----- temp variables
3226  long long int thisTask; // task based approach for expose loop with small granularity
3227  long long int numTasks=qureg.numAmpsPerChunk;
3228 
3229  // ---------------------------------------------------------------- //
3230  // find probability //
3231  // ---------------------------------------------------------------- //
3232 
3233  // initialise returned value
3234  totalProbability = 0.0;
3235 
3236  qreal *stateVecReal = qureg.stateVec.real;
3237  qreal *stateVecImag = qureg.stateVec.imag;
3238 
3239 # ifdef _OPENMP
3240 # pragma omp parallel \
3241  shared (numTasks,stateVecReal,stateVecImag) \
3242  private (thisTask) \
3243  reduction ( +:totalProbability )
3244 # endif
3245  {
3246 # ifdef _OPENMP
3247 # pragma omp for schedule (static)
3248 # endif
3249  for (thisTask=0; thisTask<numTasks; thisTask++) {
3250  totalProbability += stateVecReal[thisTask]*stateVecReal[thisTask]
3251  + stateVecImag[thisTask]*stateVecImag[thisTask];
3252  }
3253  }
3254 
3255  return totalProbability;
3256 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_calcProbOfOutcome().

◆ statevec_findProbabilityOfZeroLocal()

qreal statevec_findProbabilityOfZeroLocal ( Qureg  qureg,
const int  measureQubit 
)

Measure the total probability of a specified qubit being in the zero state across all amplitudes in this chunk.

Size of regions to skip is less than the size of one chunk.

Parameters
[in]quregobject representing the set of qubits
[in]measureQubitqubit to measure
Returns
probability of qubit measureQubit being zero

Definition at line 3166 of file QuEST_cpu.c.

3168 {
3169  // ----- sizes
3170  long long int sizeBlock, // size of blocks
3171  sizeHalfBlock; // size of blocks halved
3172  // ----- indices
3173  long long int thisBlock, // current block
3174  index; // current index for first half block
3175  // ----- measured probability
3176  qreal totalProbability; // probability (returned) value
3177  // ----- temp variables
3178  long long int thisTask;
3179  long long int numTasks=qureg.numAmpsPerChunk>>1;
3180 
3181  // ---------------------------------------------------------------- //
3182  // dimensions //
3183  // ---------------------------------------------------------------- //
3184  sizeHalfBlock = 1LL << (measureQubit); // number of state vector elements to sum,
3185  // and then the number to skip
3186  sizeBlock = 2LL * sizeHalfBlock; // size of blocks (pairs of measure and skip entries)
3187 
3188  // initialise returned value
3189  totalProbability = 0.0;
3190 
3191  qreal *stateVecReal = qureg.stateVec.real;
3192  qreal *stateVecImag = qureg.stateVec.imag;
3193 
3194 # ifdef _OPENMP
3195 # pragma omp parallel \
3196  shared (numTasks,sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag) \
3197  private (thisTask,thisBlock,index) \
3198  reduction ( +:totalProbability )
3199 # endif
3200  {
3201 # ifdef _OPENMP
3202 # pragma omp for schedule (static)
3203 # endif
3204  for (thisTask=0; thisTask<numTasks; thisTask++) {
3205  thisBlock = thisTask / sizeHalfBlock;
3206  index = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
3207 
3208  totalProbability += stateVecReal[index]*stateVecReal[index]
3209  + stateVecImag[index]*stateVecImag[index];
3210  }
3211  }
3212  return totalProbability;
3213 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_calcProbOfOutcome().

◆ statevec_hadamardDistributed()

void statevec_hadamardDistributed ( Qureg  qureg,
ComplexArray  stateVecUp,
ComplexArray  stateVecLo,
ComplexArray  stateVecOut,
int  updateUpper 
)

Rotate a single qubit by {{1,1},{1,-1}}/sqrt2.

Operate on a subset of the state vector with upper and lower block values stored seperately. This rotation is just swapping upper and lower values, and stateVecIn must already be the correct section for this chunk

Parameters
[in,out]quregobject representing the set of qubits
[in]stateVecInprobability amplitudes in lower or upper half of a block depending on chunkId
[in]updateUpperflag, 1: updating upper values, 0: updating lower values in block
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2894 of file QuEST_cpu.c.

2899 {
2900 
2901  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2902  long long int thisTask;
2903  const long long int numTasks=qureg.numAmpsPerChunk;
2904 
2905  int sign;
2906  if (updateUpper) sign=1;
2907  else sign=-1;
2908 
2909  qreal recRoot2 = 1.0/sqrt(2);
2910 
2911  qreal *stateVecRealUp=stateVecUp.real, *stateVecImagUp=stateVecUp.imag;
2912  qreal *stateVecRealLo=stateVecLo.real, *stateVecImagLo=stateVecLo.imag;
2913  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2914 
2915 # ifdef _OPENMP
2916 # pragma omp parallel \
2917  default (none) \
2918  shared (stateVecRealUp,stateVecImagUp,stateVecRealLo,stateVecImagLo,stateVecRealOut,stateVecImagOut, \
2919  recRoot2, sign) \
2920  private (thisTask,stateRealUp,stateImagUp,stateRealLo,stateImagLo)
2921 # endif
2922  {
2923 # ifdef _OPENMP
2924 # pragma omp for schedule (static)
2925 # endif
2926  for (thisTask=0; thisTask<numTasks; thisTask++) {
2927  // store current state vector values in temp variables
2928  stateRealUp = stateVecRealUp[thisTask];
2929  stateImagUp = stateVecImagUp[thisTask];
2930 
2931  stateRealLo = stateVecRealLo[thisTask];
2932  stateImagLo = stateVecImagLo[thisTask];
2933 
2934  stateVecRealOut[thisTask] = recRoot2*(stateRealUp + sign*stateRealLo);
2935  stateVecImagOut[thisTask] = recRoot2*(stateImagUp + sign*stateImagLo);
2936  }
2937  }
2938 }

References Qureg::numAmpsPerChunk, and qreal.

Referenced by statevec_hadamard().

◆ statevec_hadamardLocal()

void statevec_hadamardLocal ( Qureg  qureg,
const int  targetQubit 
)

Definition at line 2834 of file QuEST_cpu.c.

2835 {
2836  long long int sizeBlock, sizeHalfBlock;
2837  long long int thisBlock, // current block
2838  indexUp,indexLo; // current index and corresponding index in lower half block
2839 
2840  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2841  long long int thisTask;
2842  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2843 
2844  // set dimensions
2845  sizeHalfBlock = 1LL << targetQubit;
2846  sizeBlock = 2LL * sizeHalfBlock;
2847 
2848  // Can't use qureg.stateVec as a private OMP var
2849  qreal *stateVecReal = qureg.stateVec.real;
2850  qreal *stateVecImag = qureg.stateVec.imag;
2851 
2852  qreal recRoot2 = 1.0/sqrt(2);
2853 
2854 # ifdef _OPENMP
2855 # pragma omp parallel \
2856  default (none) \
2857  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag, recRoot2) \
2858  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp,stateRealLo,stateImagLo)
2859 # endif
2860  {
2861 # ifdef _OPENMP
2862 # pragma omp for schedule (static)
2863 # endif
2864  for (thisTask=0; thisTask<numTasks; thisTask++) {
2865  thisBlock = thisTask / sizeHalfBlock;
2866  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2867  indexLo = indexUp + sizeHalfBlock;
2868 
2869  stateRealUp = stateVecReal[indexUp];
2870  stateImagUp = stateVecImag[indexUp];
2871 
2872  stateRealLo = stateVecReal[indexLo];
2873  stateImagLo = stateVecImag[indexLo];
2874 
2875  stateVecReal[indexUp] = recRoot2*(stateRealUp + stateRealLo);
2876  stateVecImag[indexUp] = recRoot2*(stateImagUp + stateImagLo);
2877 
2878  stateVecReal[indexLo] = recRoot2*(stateRealUp - stateRealLo);
2879  stateVecImag[indexLo] = recRoot2*(stateImagUp - stateImagLo);
2880  }
2881  }
2882 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_hadamard().

◆ statevec_multiControlledMultiQubitUnitaryLocal()

void statevec_multiControlledMultiQubitUnitaryLocal ( Qureg  qureg,
long long int  ctrlMask,
int *  targs,
const int  numTargs,
ComplexMatrixN  u 
)

Definition at line 1814 of file QuEST_cpu.c.

1815 {
1816  // can't use qureg.stateVec as a private OMP var
1817  qreal *reVec = qureg.stateVec.real;
1818  qreal *imVec = qureg.stateVec.imag;
1819 
1820  long long int numTasks = qureg.numAmpsPerChunk >> numTargs; // kernel called on every 1 in 2^numTargs amplitudes
1821  long long int numTargAmps = 1 << u.numQubits; // num amps to be modified by each task
1822 
1823  // the global (between all nodes) index of this node's start index
1824  long long int globalIndStart = qureg.chunkId*qureg.numAmpsPerChunk;
1825 
1826  long long int thisTask;
1827  long long int thisInd00; // this thread's index of |..0..0..> (target qubits = 0)
1828  long long int thisGlobalInd00; // the global (between all nodes) index of this thread's |..0..0..> state
1829  long long int ind; // each thread's iteration of amplitudes to modify
1830  int i, t, r, c; // each thread's iteration of amps and targets
1831  qreal reElem, imElem; // each thread's iteration of u elements
1832 
1833  // each thread/task will record and modify numTargAmps amplitudes, privately
1834  // (of course, tasks eliminated by the ctrlMask won't edit their allocation)
1835  long long int ampInds[numTargAmps];
1836  qreal reAmps[numTargAmps];
1837  qreal imAmps[numTargAmps];
1838 
1839  // we need a sorted targets list to find thisInd00 for each task.
1840  // we can't modify targets, because the user-ordering of targets matters in u
1841  int sortedTargs[numTargs];
1842  for (int t=0; t < numTargs; t++)
1843  sortedTargs[t] = targs[t];
1844  qsort(sortedTargs, numTargs, sizeof(int), qsortComp);
1845 
1846 # ifdef _OPENMP
1847 # pragma omp parallel \
1848  default (none) \
1849  shared (reVec,imVec, numTasks,numTargAmps,globalIndStart, ctrlMask,targs,sortedTargs,u) \
1850  private (thisTask,thisInd00,thisGlobalInd00,ind,i,t,r,c,reElem,imElem, ampInds,reAmps,imAmps)
1851 # endif
1852  {
1853 # ifdef _OPENMP
1854 # pragma omp for schedule (static)
1855 # endif
1856  for (thisTask=0; thisTask<numTasks; thisTask++) {
1857 
1858  // find this task's start index (where all targs are 0)
1859  thisInd00 = thisTask;
1860  for (t=0; t < numTargs; t++)
1861  thisInd00 = insertZeroBit(thisInd00, sortedTargs[t]);
1862 
1863  // this task only modifies amplitudes if control qubits are 1 for this state
1864  thisGlobalInd00 = thisInd00 + globalIndStart;
1865  if (ctrlMask && ((ctrlMask & thisGlobalInd00) != ctrlMask))
1866  continue;
1867 
1868  // determine the indices and record values of this tasks's target amps
1869  for (i=0; i < numTargAmps; i++) {
1870 
1871  // get statevec index of current target qubit assignment
1872  ind = thisInd00;
1873  for (t=0; t < numTargs; t++)
1874  if (extractBit(t, i))
1875  ind = flipBit(ind, targs[t]);
1876 
1877  // update this tasks's private arrays
1878  ampInds[i] = ind;
1879  reAmps [i] = reVec[ind];
1880  imAmps [i] = imVec[ind];
1881  }
1882 
1883  // modify this tasks's target amplitudes
1884  for (r=0; r < numTargAmps; r++) {
1885  ind = ampInds[r];
1886  reVec[ind] = 0;
1887  imVec[ind] = 0;
1888 
1889  for (c=0; c < numTargAmps; c++) {
1890  reElem = u.real[r][c];
1891  imElem = u.imag[r][c];
1892  reVec[ind] += reAmps[c]*reElem - imAmps[c]*imElem;
1893  imVec[ind] += reAmps[c]*imElem + imAmps[c]*reElem;
1894  }
1895  }
1896  }
1897  }
1898 }

References Qureg::chunkId, extractBit(), flipBit(), ComplexMatrixN::imag, insertZeroBit(), Qureg::numAmpsPerChunk, ComplexMatrixN::numQubits, qreal, qsortComp(), ComplexMatrixN::real, and Qureg::stateVec.

Referenced by statevec_multiControlledMultiQubitUnitary().

◆ statevec_multiControlledTwoQubitUnitaryLocal()

void statevec_multiControlledTwoQubitUnitaryLocal ( Qureg  qureg,
long long int  ctrlMask,
const int  q1,
const int  q2,
ComplexMatrix4  u 
)

Definition at line 1715 of file QuEST_cpu.c.

1715  {
1716 
1717  // can't use qureg.stateVec as a private OMP var
1718  qreal *reVec = qureg.stateVec.real;
1719  qreal *imVec = qureg.stateVec.imag;
1720 
1721  // the global (between all nodes) index of this node's start index
1722  long long int globalIndStart = qureg.chunkId*qureg.numAmpsPerChunk;
1723 
1724  long long int numTasks = qureg.numAmpsPerChunk >> 2; // each iteration updates 4 amplitudes
1725  long long int thisTask;
1726  long long int thisGlobalInd00;
1727  long long int ind00, ind01, ind10, ind11;
1728  qreal re00, re01, re10, re11;
1729  qreal im00, im01, im10, im11;
1730 
1731 # ifdef _OPENMP
1732 # pragma omp parallel \
1733  default (none) \
1734  shared (reVec,imVec,globalIndStart,numTasks,ctrlMask,u) \
1735  private (thisTask, thisGlobalInd00, ind00,ind01,ind10,ind11, re00,re01,re10,re11, im00,im01,im10,im11)
1736 # endif
1737  {
1738 # ifdef _OPENMP
1739 # pragma omp for schedule (static)
1740 # endif
1741  for (thisTask=0; thisTask<numTasks; thisTask++) {
1742 
1743  // determine ind00 of |..0..0..>
1744  ind00 = insertTwoZeroBits(thisTask, q1, q2);
1745 
1746  // skip amplitude if controls aren't in 1 state (overloaded for speed)
1747  thisGlobalInd00 = ind00 + globalIndStart;
1748  if (ctrlMask && ((ctrlMask & thisGlobalInd00) != ctrlMask))
1749  continue;
1750 
1751  // inds of |..0..1..>, |..1..0..> and |..1..1..>
1752  ind01 = flipBit(ind00, q1);
1753  ind10 = flipBit(ind00, q2);
1754  ind11 = flipBit(ind01, q2);
1755 
1756  // extract statevec amplitudes
1757  re00 = reVec[ind00]; im00 = imVec[ind00];
1758  re01 = reVec[ind01]; im01 = imVec[ind01];
1759  re10 = reVec[ind10]; im10 = imVec[ind10];
1760  re11 = reVec[ind11]; im11 = imVec[ind11];
1761 
1762  // apply u * {amp00, amp01, amp10, amp11}
1763  reVec[ind00] =
1764  u.real[0][0]*re00 - u.imag[0][0]*im00 +
1765  u.real[0][1]*re01 - u.imag[0][1]*im01 +
1766  u.real[0][2]*re10 - u.imag[0][2]*im10 +
1767  u.real[0][3]*re11 - u.imag[0][3]*im11;
1768  imVec[ind00] =
1769  u.imag[0][0]*re00 + u.real[0][0]*im00 +
1770  u.imag[0][1]*re01 + u.real[0][1]*im01 +
1771  u.imag[0][2]*re10 + u.real[0][2]*im10 +
1772  u.imag[0][3]*re11 + u.real[0][3]*im11;
1773 
1774  reVec[ind01] =
1775  u.real[1][0]*re00 - u.imag[1][0]*im00 +
1776  u.real[1][1]*re01 - u.imag[1][1]*im01 +
1777  u.real[1][2]*re10 - u.imag[1][2]*im10 +
1778  u.real[1][3]*re11 - u.imag[1][3]*im11;
1779  imVec[ind01] =
1780  u.imag[1][0]*re00 + u.real[1][0]*im00 +
1781  u.imag[1][1]*re01 + u.real[1][1]*im01 +
1782  u.imag[1][2]*re10 + u.real[1][2]*im10 +
1783  u.imag[1][3]*re11 + u.real[1][3]*im11;
1784 
1785  reVec[ind10] =
1786  u.real[2][0]*re00 - u.imag[2][0]*im00 +
1787  u.real[2][1]*re01 - u.imag[2][1]*im01 +
1788  u.real[2][2]*re10 - u.imag[2][2]*im10 +
1789  u.real[2][3]*re11 - u.imag[2][3]*im11;
1790  imVec[ind10] =
1791  u.imag[2][0]*re00 + u.real[2][0]*im00 +
1792  u.imag[2][1]*re01 + u.real[2][1]*im01 +
1793  u.imag[2][2]*re10 + u.real[2][2]*im10 +
1794  u.imag[2][3]*re11 + u.real[2][3]*im11;
1795 
1796  reVec[ind11] =
1797  u.real[3][0]*re00 - u.imag[3][0]*im00 +
1798  u.real[3][1]*re01 - u.imag[3][1]*im01 +
1799  u.real[3][2]*re10 - u.imag[3][2]*im10 +
1800  u.real[3][3]*re11 - u.imag[3][3]*im11;
1801  imVec[ind11] =
1802  u.imag[3][0]*re00 + u.real[3][0]*im00 +
1803  u.imag[3][1]*re01 + u.real[3][1]*im01 +
1804  u.imag[3][2]*re10 + u.real[3][2]*im10 +
1805  u.imag[3][3]*re11 + u.real[3][3]*im11;
1806  }
1807  }
1808 }

References Qureg::chunkId, flipBit(), ComplexMatrix4::imag, insertTwoZeroBits(), Qureg::numAmpsPerChunk, qreal, ComplexMatrix4::real, and Qureg::stateVec.

Referenced by statevec_multiControlledTwoQubitUnitary().

◆ statevec_multiControlledUnitaryDistributed()

void statevec_multiControlledUnitaryDistributed ( Qureg  qureg,
const int  targetQubit,
long long int  ctrlQubitsMask,
long long int  ctrlFlipMask,
Complex  rot1,
Complex  rot2,
ComplexArray  stateVecUp,
ComplexArray  stateVecLo,
ComplexArray  stateVecOut 
)

Apply a unitary operation to a single qubit in the state vector of probability amplitudes, given a subset of the state vector with upper and lower block values stored seperately.

Only perform the rotation where all the control qubits are 1.

Parameters
[in,out]quregobject representing the set of qubits
[in]targetQubitqubit to rotate
[in]ctrlQubitsMaska bit mask indicating whether each qubit is a control (1) or not (0)
[in]ctrlFlipMaska bit mask indicating whether each qubit (only controls are relevant) should be flipped when checking the control condition
[in]rot1rotation angle
[in]rot2rotation angle
[in]stateVecUpprobability amplitudes in upper half of a block
[in]stateVecLoprobability amplitudes in lower half of a block
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2413 of file QuEST_cpu.c.

2421 {
2422 
2423  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2424  long long int thisTask;
2425  const long long int numTasks=qureg.numAmpsPerChunk;
2426  const long long int chunkSize=qureg.numAmpsPerChunk;
2427  const long long int chunkId=qureg.chunkId;
2428 
2429  qreal rot1Real=rot1.real, rot1Imag=rot1.imag;
2430  qreal rot2Real=rot2.real, rot2Imag=rot2.imag;
2431  qreal *stateVecRealUp=stateVecUp.real, *stateVecImagUp=stateVecUp.imag;
2432  qreal *stateVecRealLo=stateVecLo.real, *stateVecImagLo=stateVecLo.imag;
2433  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2434 
2435 # ifdef _OPENMP
2436 # pragma omp parallel \
2437  default (none) \
2438  shared (stateVecRealUp,stateVecImagUp,stateVecRealLo,stateVecImagLo,stateVecRealOut,stateVecImagOut, \
2439  rot1Real,rot1Imag, rot2Real,rot2Imag, ctrlQubitsMask,ctrlFlipMask) \
2440  private (thisTask,stateRealUp,stateImagUp,stateRealLo,stateImagLo)
2441 # endif
2442  {
2443 # ifdef _OPENMP
2444 # pragma omp for schedule (static)
2445 # endif
2446  for (thisTask=0; thisTask<numTasks; thisTask++) {
2447  if (ctrlQubitsMask == (ctrlQubitsMask & ((thisTask+chunkId*chunkSize) ^ ctrlFlipMask))) {
2448  // store current state vector values in temp variables
2449  stateRealUp = stateVecRealUp[thisTask];
2450  stateImagUp = stateVecImagUp[thisTask];
2451 
2452  stateRealLo = stateVecRealLo[thisTask];
2453  stateImagLo = stateVecImagLo[thisTask];
2454 
2455  stateVecRealOut[thisTask] = rot1Real*stateRealUp - rot1Imag*stateImagUp
2456  + rot2Real*stateRealLo - rot2Imag*stateImagLo;
2457  stateVecImagOut[thisTask] = rot1Real*stateImagUp + rot1Imag*stateRealUp
2458  + rot2Real*stateImagLo + rot2Imag*stateRealLo;
2459  }
2460  }
2461  }
2462 }

References Qureg::chunkId, Complex::imag, Qureg::numAmpsPerChunk, qreal, and Complex::real.

Referenced by statevec_multiControlledUnitary().

◆ statevec_multiControlledUnitaryLocal()

void statevec_multiControlledUnitaryLocal ( Qureg  qureg,
const int  targetQubit,
long long int  ctrlQubitsMask,
long long int  ctrlFlipMask,
ComplexMatrix2  u 
)

Definition at line 2140 of file QuEST_cpu.c.

2144 {
2145  long long int sizeBlock, sizeHalfBlock;
2146  long long int thisBlock, // current block
2147  indexUp,indexLo; // current index and corresponding index in lower half block
2148 
2149  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2150  long long int thisTask;
2151  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2152  const long long int chunkSize=qureg.numAmpsPerChunk;
2153  const long long int chunkId=qureg.chunkId;
2154 
2155  // set dimensions
2156  sizeHalfBlock = 1LL << targetQubit;
2157  sizeBlock = 2LL * sizeHalfBlock;
2158 
2159  // Can't use qureg.stateVec as a private OMP var
2160  qreal *stateVecReal = qureg.stateVec.real;
2161  qreal *stateVecImag = qureg.stateVec.imag;
2162 
2163 # ifdef _OPENMP
2164 # pragma omp parallel \
2165  default (none) \
2166  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag, u, ctrlQubitsMask,ctrlFlipMask) \
2167  private (thisTask,thisBlock, indexUp,indexLo, stateRealUp,stateImagUp,stateRealLo,stateImagLo)
2168 # endif
2169  {
2170 # ifdef _OPENMP
2171 # pragma omp for schedule (static)
2172 # endif
2173  for (thisTask=0; thisTask<numTasks; thisTask++) {
2174 
2175  thisBlock = thisTask / sizeHalfBlock;
2176  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2177  indexLo = indexUp + sizeHalfBlock;
2178 
2179 
2180  // take the basis index, flip the designated (XOR) 'control' bits, AND with the controls.
2181  // if this equals the control mask, the control qubits have the desired values in the basis index
2182  if (ctrlQubitsMask == (ctrlQubitsMask & ((indexUp+chunkId*chunkSize) ^ ctrlFlipMask))) {
2183  // store current state vector values in temp variables
2184  stateRealUp = stateVecReal[indexUp];
2185  stateImagUp = stateVecImag[indexUp];
2186 
2187  stateRealLo = stateVecReal[indexLo];
2188  stateImagLo = stateVecImag[indexLo];
2189 
2190  // state[indexUp] = u00 * state[indexUp] + u01 * state[indexLo]
2191  stateVecReal[indexUp] = u.real[0][0]*stateRealUp - u.imag[0][0]*stateImagUp
2192  + u.real[0][1]*stateRealLo - u.imag[0][1]*stateImagLo;
2193  stateVecImag[indexUp] = u.real[0][0]*stateImagUp + u.imag[0][0]*stateRealUp
2194  + u.real[0][1]*stateImagLo + u.imag[0][1]*stateRealLo;
2195 
2196  // state[indexLo] = u10 * state[indexUp] + u11 * state[indexLo]
2197  stateVecReal[indexLo] = u.real[1][0]*stateRealUp - u.imag[1][0]*stateImagUp
2198  + u.real[1][1]*stateRealLo - u.imag[1][1]*stateImagLo;
2199  stateVecImag[indexLo] = u.real[1][0]*stateImagUp + u.imag[1][0]*stateRealUp
2200  + u.real[1][1]*stateImagLo + u.imag[1][1]*stateRealLo;
2201  }
2202  }
2203  }
2204 
2205 }

References Qureg::chunkId, ComplexMatrix2::imag, Qureg::numAmpsPerChunk, qreal, ComplexMatrix2::real, and Qureg::stateVec.

Referenced by statevec_multiControlledUnitary().

◆ statevec_pauliXDistributed()

void statevec_pauliXDistributed ( Qureg  qureg,
ComplexArray  stateVecIn,
ComplexArray  stateVecOut 
)

Rotate a single qubit by {{0,1},{1,0}.

Operate on a subset of the state vector with upper and lower block values stored seperately. This rotation is just swapping upper and lower values, and stateVecIn must already be the correct section for this chunk

Remarks
Qubits are zero-based and the
the first qubit is the rightmost
Parameters
[in,out]quregobject representing the set of qubits
[in]stateVecInprobability amplitudes in lower or upper half of a block depending on chunkId
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2522 of file QuEST_cpu.c.

2525 {
2526 
2527  long long int thisTask;
2528  const long long int numTasks=qureg.numAmpsPerChunk;
2529 
2530  qreal *stateVecRealIn=stateVecIn.real, *stateVecImagIn=stateVecIn.imag;
2531  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2532 
2533 # ifdef _OPENMP
2534 # pragma omp parallel \
2535  default (none) \
2536  shared (stateVecRealIn,stateVecImagIn,stateVecRealOut,stateVecImagOut) \
2537  private (thisTask)
2538 # endif
2539  {
2540 # ifdef _OPENMP
2541 # pragma omp for schedule (static)
2542 # endif
2543  for (thisTask=0; thisTask<numTasks; thisTask++) {
2544  stateVecRealOut[thisTask] = stateVecRealIn[thisTask];
2545  stateVecImagOut[thisTask] = stateVecImagIn[thisTask];
2546  }
2547  }
2548 }

References Qureg::numAmpsPerChunk, and qreal.

Referenced by statevec_pauliX().

◆ statevec_pauliXLocal()

void statevec_pauliXLocal ( Qureg  qureg,
const int  targetQubit 
)

Definition at line 2464 of file QuEST_cpu.c.

2465 {
2466  long long int sizeBlock, sizeHalfBlock;
2467  long long int thisBlock, // current block
2468  indexUp,indexLo; // current index and corresponding index in lower half block
2469 
2470  qreal stateRealUp,stateImagUp;
2471  long long int thisTask;
2472  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2473 
2474  // set dimensions
2475  sizeHalfBlock = 1LL << targetQubit;
2476  sizeBlock = 2LL * sizeHalfBlock;
2477 
2478  // Can't use qureg.stateVec as a private OMP var
2479  qreal *stateVecReal = qureg.stateVec.real;
2480  qreal *stateVecImag = qureg.stateVec.imag;
2481 
2482 # ifdef _OPENMP
2483 # pragma omp parallel \
2484  default (none) \
2485  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag) \
2486  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp)
2487 # endif
2488  {
2489 # ifdef _OPENMP
2490 # pragma omp for schedule (static)
2491 # endif
2492  for (thisTask=0; thisTask<numTasks; thisTask++) {
2493  thisBlock = thisTask / sizeHalfBlock;
2494  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2495  indexLo = indexUp + sizeHalfBlock;
2496 
2497  stateRealUp = stateVecReal[indexUp];
2498  stateImagUp = stateVecImag[indexUp];
2499 
2500  stateVecReal[indexUp] = stateVecReal[indexLo];
2501  stateVecImag[indexUp] = stateVecImag[indexLo];
2502 
2503  stateVecReal[indexLo] = stateRealUp;
2504  stateVecImag[indexLo] = stateImagUp;
2505  }
2506  }
2507 
2508 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_pauliX().

◆ statevec_pauliYDistributed()

void statevec_pauliYDistributed ( Qureg  qureg,
ComplexArray  stateVecIn,
ComplexArray  stateVecOut,
int  updateUpper,
const int  conjFac 
)

Rotate a single qubit by +-{{0,-i},{i,0}.

Operate on a subset of the state vector with upper and lower block values stored seperately. This rotation is just swapping upper and lower values, and stateVecIn must already be the correct section for this chunk

Remarks
Qubits are zero-based and the
the first qubit is the rightmost
Parameters
[in,out]quregobject representing the set of qubits
[in]stateVecInprobability amplitudes in lower or upper half of a block depending on chunkId
[in]updateUpperflag, 1: updating upper values, 0: updating lower values in block
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2704 of file QuEST_cpu.c.

2708 {
2709 
2710  long long int thisTask;
2711  const long long int numTasks=qureg.numAmpsPerChunk;
2712 
2713  qreal *stateVecRealIn=stateVecIn.real, *stateVecImagIn=stateVecIn.imag;
2714  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2715 
2716  int realSign=1, imagSign=1;
2717  if (updateUpper) imagSign=-1;
2718  else realSign = -1;
2719 
2720 # ifdef _OPENMP
2721 # pragma omp parallel \
2722  default (none) \
2723  shared (stateVecRealIn,stateVecImagIn,stateVecRealOut,stateVecImagOut,realSign,imagSign) \
2724  private (thisTask)
2725 # endif
2726  {
2727 # ifdef _OPENMP
2728 # pragma omp for schedule (static)
2729 # endif
2730  for (thisTask=0; thisTask<numTasks; thisTask++) {
2731  stateVecRealOut[thisTask] = conjFac * realSign * stateVecImagIn[thisTask];
2732  stateVecImagOut[thisTask] = conjFac * imagSign * stateVecRealIn[thisTask];
2733  }
2734  }
2735 }

References Qureg::numAmpsPerChunk, and qreal.

Referenced by statevec_pauliY(), and statevec_pauliYConj().

◆ statevec_pauliYLocal()

void statevec_pauliYLocal ( Qureg  qureg,
const int  targetQubit,
const int  conjFac 
)

Definition at line 2647 of file QuEST_cpu.c.

2648 {
2649  long long int sizeBlock, sizeHalfBlock;
2650  long long int thisBlock, // current block
2651  indexUp,indexLo; // current index and corresponding index in lower half block
2652 
2653  qreal stateRealUp,stateImagUp;
2654  long long int thisTask;
2655  const long long int numTasks=qureg.numAmpsPerChunk>>1;
2656 
2657  // set dimensions
2658  sizeHalfBlock = 1LL << targetQubit;
2659  sizeBlock = 2LL * sizeHalfBlock;
2660 
2661  // Can't use qureg.stateVec as a private OMP var
2662  qreal *stateVecReal = qureg.stateVec.real;
2663  qreal *stateVecImag = qureg.stateVec.imag;
2664 
2665 # ifdef _OPENMP
2666 # pragma omp parallel \
2667  default (none) \
2668  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag) \
2669  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp)
2670 # endif
2671  {
2672 # ifdef _OPENMP
2673 # pragma omp for schedule (static)
2674 # endif
2675  for (thisTask=0; thisTask<numTasks; thisTask++) {
2676  thisBlock = thisTask / sizeHalfBlock;
2677  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
2678  indexLo = indexUp + sizeHalfBlock;
2679 
2680  stateRealUp = stateVecReal[indexUp];
2681  stateImagUp = stateVecImag[indexUp];
2682 
2683  stateVecReal[indexUp] = conjFac * stateVecImag[indexLo];
2684  stateVecImag[indexUp] = conjFac * -stateVecReal[indexLo];
2685  stateVecReal[indexLo] = conjFac * -stateImagUp;
2686  stateVecImag[indexLo] = conjFac * stateRealUp;
2687  }
2688  }
2689 }

References Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_pauliY(), and statevec_pauliYConj().

◆ statevec_swapQubitAmpsDistributed()

void statevec_swapQubitAmpsDistributed ( Qureg  qureg,
int  pairRank,
int  qb1,
int  qb2 
)

qureg.pairStateVec contains the entire set of amplitudes of the paired node which includes the set of all amplitudes which need to be swapped between |..0..1..> and |..1..0..>

Definition at line 3539 of file QuEST_cpu.c.

3539  {
3540 
3541  // can't use qureg.stateVec as a private OMP var
3542  qreal *reVec = qureg.stateVec.real;
3543  qreal *imVec = qureg.stateVec.imag;
3544  qreal *rePairVec = qureg.pairStateVec.real;
3545  qreal *imPairVec = qureg.pairStateVec.imag;
3546 
3547  long long int numLocalAmps = qureg.numAmpsPerChunk;
3548  long long int globalStartInd = qureg.chunkId * numLocalAmps;
3549  long long int pairGlobalStartInd = pairRank * numLocalAmps;
3550 
3551  long long int localInd, globalInd;
3552  long long int pairLocalInd, pairGlobalInd;
3553 
3554 # ifdef _OPENMP
3555 # pragma omp parallel \
3556  default (none) \
3557  shared (reVec,imVec,rePairVec,imPairVec,numLocalAmps,globalStartInd,pairGlobalStartInd,qb1,qb2) \
3558  private (localInd,globalInd, pairLocalInd,pairGlobalInd)
3559 # endif
3560  {
3561 # ifdef _OPENMP
3562 # pragma omp for schedule (static)
3563 # endif
3564  for (localInd=0; localInd < numLocalAmps; localInd++) {
3565 
3566  globalInd = globalStartInd + localInd;
3567  if (isOddParity(globalInd, qb1, qb2)) {
3568 
3569  pairGlobalInd = flipBit(flipBit(globalInd, qb1), qb2);
3570  pairLocalInd = pairGlobalInd - pairGlobalStartInd;
3571 
3572  reVec[localInd] = rePairVec[pairLocalInd];
3573  imVec[localInd] = imPairVec[pairLocalInd];
3574  }
3575  }
3576  }
3577 }

References Qureg::chunkId, flipBit(), isOddParity(), Qureg::numAmpsPerChunk, Qureg::pairStateVec, qreal, and Qureg::stateVec.

Referenced by statevec_swapQubitAmps().

◆ statevec_swapQubitAmpsLocal()

void statevec_swapQubitAmpsLocal ( Qureg  qureg,
int  qb1,
int  qb2 
)

It is ensured that all amplitudes needing to be swapped are on this node.

This means that amplitudes for |a 0..0..> to |a 1..1..> all exist on this node and each node has a different bit-string prefix "a". The prefix 'a' (and ergo, the chunkID) don't enter the calculations for the offset of |a 0..1..> and |a 1..0..> from |a 0..0..> and ergo are not featured below.

Definition at line 3496 of file QuEST_cpu.c.

3496  {
3497 
3498  // can't use qureg.stateVec as a private OMP var
3499  qreal *reVec = qureg.stateVec.real;
3500  qreal *imVec = qureg.stateVec.imag;
3501 
3502  long long int numTasks = qureg.numAmpsPerChunk >> 2; // each iteration updates 2 amps and skips 2 amps
3503  long long int thisTask;
3504  long long int ind00, ind01, ind10;
3505  qreal re01, re10;
3506  qreal im01, im10;
3507 
3508 # ifdef _OPENMP
3509 # pragma omp parallel \
3510  default (none) \
3511  shared (reVec,imVec,numTasks,qb1,qb2) \
3512  private (thisTask, ind00,ind01,ind10, re01,re10, im01,im10)
3513 # endif
3514  {
3515 # ifdef _OPENMP
3516 # pragma omp for schedule (static)
3517 # endif
3518  for (thisTask=0; thisTask<numTasks; thisTask++) {
3519  // determine ind00 of |..0..0..>, |..0..1..> and |..1..0..>
3520  ind00 = insertTwoZeroBits(thisTask, qb1, qb2);
3521  ind01 = flipBit(ind00, qb1);
3522  ind10 = flipBit(ind00, qb2);
3523 
3524  // extract statevec amplitudes
3525  re01 = reVec[ind01]; im01 = imVec[ind01];
3526  re10 = reVec[ind10]; im10 = imVec[ind10];
3527 
3528  // swap 01 and 10 amps
3529  reVec[ind01] = re10; reVec[ind10] = re01;
3530  imVec[ind01] = im10; imVec[ind10] = im01;
3531  }
3532  }
3533 }

References flipBit(), insertTwoZeroBits(), Qureg::numAmpsPerChunk, qreal, and Qureg::stateVec.

Referenced by statevec_swapQubitAmps().

◆ statevec_unitaryDistributed()

void statevec_unitaryDistributed ( Qureg  qureg,
Complex  rot1,
Complex  rot2,
ComplexArray  stateVecUp,
ComplexArray  stateVecLo,
ComplexArray  stateVecOut 
)

Apply a unitary operation to a single qubit given a subset of the state vector with upper and lower block values stored seperately.

Remarks
Qubits are zero-based and the first qubit is the rightmost
Parameters
[in,out]quregobject representing the set of qubits
[in]uunitary matrix to apply
[in]stateVecUpprobability amplitudes in upper half of a block
[in]stateVecLoprobability amplitudes in lower half of a block
[out]stateVecOutarray section to update (will correspond to either the lower or upper half of a block)

Definition at line 2024 of file QuEST_cpu.c.

2029 {
2030 
2031  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
2032  long long int thisTask;
2033  const long long int numTasks=qureg.numAmpsPerChunk;
2034 
2035  qreal rot1Real=rot1.real, rot1Imag=rot1.imag;
2036  qreal rot2Real=rot2.real, rot2Imag=rot2.imag;
2037  qreal *stateVecRealUp=stateVecUp.real, *stateVecImagUp=stateVecUp.imag;
2038  qreal *stateVecRealLo=stateVecLo.real, *stateVecImagLo=stateVecLo.imag;
2039  qreal *stateVecRealOut=stateVecOut.real, *stateVecImagOut=stateVecOut.imag;
2040 
2041 
2042 # ifdef _OPENMP
2043 # pragma omp parallel \
2044  default (none) \
2045  shared (stateVecRealUp,stateVecImagUp,stateVecRealLo,stateVecImagLo,stateVecRealOut,stateVecImagOut, \
2046  rot1Real, rot1Imag, rot2Real, rot2Imag) \
2047  private (thisTask,stateRealUp,stateImagUp,stateRealLo,stateImagLo)
2048 # endif
2049  {
2050 # ifdef _OPENMP
2051 # pragma omp for schedule (static)
2052 # endif
2053  for (thisTask=0; thisTask<numTasks; thisTask++) {
2054  // store current state vector values in temp variables
2055  stateRealUp = stateVecRealUp[thisTask];
2056  stateImagUp = stateVecImagUp[thisTask];
2057 
2058  stateRealLo = stateVecRealLo[thisTask];
2059  stateImagLo = stateVecImagLo[thisTask];
2060 
2061  stateVecRealOut[thisTask] = rot1Real*stateRealUp - rot1Imag*stateImagUp
2062  + rot2Real*stateRealLo - rot2Imag*stateImagLo;
2063  stateVecImagOut[thisTask] = rot1Real*stateImagUp + rot1Imag*stateRealUp
2064  + rot2Real*stateImagLo + rot2Imag*stateRealLo;
2065  }
2066  }
2067 }

References Complex::imag, Qureg::numAmpsPerChunk, qreal, and Complex::real.

Referenced by statevec_unitary().

◆ statevec_unitaryLocal()

void statevec_unitaryLocal ( Qureg  qureg,
const int  targetQubit,
ComplexMatrix2  u 
)

Definition at line 1900 of file QuEST_cpu.c.

1901 {
1902  long long int sizeBlock, sizeHalfBlock;
1903  long long int thisBlock, // current block
1904  indexUp,indexLo; // current index and corresponding index in lower half block
1905 
1906  qreal stateRealUp,stateRealLo,stateImagUp,stateImagLo;
1907  long long int thisTask;
1908  const long long int numTasks=qureg.numAmpsPerChunk>>1;
1909 
1910  // set dimensions
1911  sizeHalfBlock = 1LL << targetQubit;
1912  sizeBlock = 2LL * sizeHalfBlock;
1913 
1914  // Can't use qureg.stateVec as a private OMP var
1915  qreal *stateVecReal = qureg.stateVec.real;
1916  qreal *stateVecImag = qureg.stateVec.imag;
1917 
1918 # ifdef _OPENMP
1919 # pragma omp parallel \
1920  default (none) \
1921  shared (sizeBlock,sizeHalfBlock, stateVecReal,stateVecImag, u) \
1922  private (thisTask,thisBlock ,indexUp,indexLo, stateRealUp,stateImagUp,stateRealLo,stateImagLo)
1923 # endif
1924  {
1925 # ifdef _OPENMP
1926 # pragma omp for schedule (static)
1927 # endif
1928  for (thisTask=0; thisTask<numTasks; thisTask++) {
1929 
1930  thisBlock = thisTask / sizeHalfBlock;
1931  indexUp = thisBlock*sizeBlock + thisTask%sizeHalfBlock;
1932  indexLo = indexUp + sizeHalfBlock;
1933 
1934  // store current state vector values in temp variables
1935  stateRealUp = stateVecReal[indexUp];
1936  stateImagUp = stateVecImag[indexUp];
1937 
1938  stateRealLo = stateVecReal[indexLo];
1939  stateImagLo = stateVecImag[indexLo];
1940 
1941 
1942  // state[indexUp] = u00 * state[indexUp] + u01 * state[indexLo]
1943  stateVecReal[indexUp] = u.real[0][0]*stateRealUp - u.imag[0][0]*stateImagUp
1944  + u.real[0][1]*stateRealLo - u.imag[0][1]*stateImagLo;
1945  stateVecImag[indexUp] = u.real[0][0]*stateImagUp + u.imag[0][0]*stateRealUp
1946  + u.real[0][1]*stateImagLo + u.imag[0][1]*stateRealLo;
1947 
1948  // state[indexLo] = u10 * state[indexUp] + u11 * state[indexLo]
1949  stateVecReal[indexLo] = u.real[1][0]*stateRealUp - u.imag[1][0]*stateImagUp
1950  + u.real[1][1]*stateRealLo - u.imag[1][1]*stateImagLo;
1951  stateVecImag[indexLo] = u.real[1][0]*stateImagUp + u.imag[1][0]*stateRealUp
1952  + u.real[1][1]*stateImagLo + u.imag[1][1]*stateRealLo;
1953 
1954  }
1955  }
1956 }

References ComplexMatrix2::imag, Qureg::numAmpsPerChunk, qreal, ComplexMatrix2::real, and Qureg::stateVec.

Referenced by statevec_unitary().

qsortComp
int qsortComp(const void *a, const void *b)
Definition: QuEST_cpu.c:1810
ComplexMatrix4::real
qreal real[4][4]
Definition: QuEST.h:127
Qureg::pairStateVec
ComplexArray pairStateVec
Temporary storage for a chunk of the state vector received from another process in the MPI version.
Definition: QuEST.h:181
extractBit
__forceinline__ __device__ int extractBit(int locationOfBitFromRight, long long int theEncodedNumber)
Definition: QuEST_gpu.cu:82
insertZeroBit
static long long int insertZeroBit(long long int number, int index)
Definition: QuEST_cpu_internal.h:42
qreal
#define qreal
insertTwoZeroBits
__forceinline__ __device__ long long int insertTwoZeroBits(long long int number, int bit1, int bit2)
Definition: QuEST_gpu.cu:106
extractBit
static int extractBit(const int locationOfBitFromRight, const long long int theEncodedNumber)
Definition: QuEST_cpu_internal.h:26
Qureg::chunkId
int chunkId
The position of the chunk of the state vector held by this process in the full state vector.
Definition: QuEST.h:174
ComplexMatrix2::imag
qreal imag[2][2]
Definition: QuEST.h:117
Qureg::numAmpsPerChunk
long long int numAmpsPerChunk
Number of probability amplitudes held in stateVec by this process In the non-MPI version,...
Definition: QuEST.h:170
ComplexMatrix4::imag
qreal imag[4][4]
Definition: QuEST.h:128
ComplexMatrixN::real
qreal ** real
Definition: QuEST.h:139
flipBit
__forceinline__ __device__ long long int flipBit(long long int number, int bitInd)
Definition: QuEST_gpu.cu:95
ComplexMatrixN::imag
qreal ** imag
Definition: QuEST.h:140
Qureg::stateVec
ComplexArray stateVec
Computational state amplitudes - a subset thereof in the MPI version.
Definition: QuEST.h:179
ComplexMatrix2::real
qreal real[2][2]
Definition: QuEST.h:116
densmatr_mixDephasing
void densmatr_mixDephasing(Qureg qureg, const int targetQubit, qreal dephase)
Definition: QuEST_cpu.c:79
ComplexMatrixN::numQubits
int numQubits
Definition: QuEST.h:138
Qureg::numQubitsRepresented
int numQubitsRepresented
The number of qubits represented in either the state-vector or density matrix.
Definition: QuEST.h:165
Qureg::numAmpsTotal
long long int numAmpsTotal
Total number of amplitudes, which are possibly distributed among machines.
Definition: QuEST.h:172
Complex::real
qreal real
Definition: QuEST.h:105
isOddParity
static int isOddParity(long long int number, int qb1, int qb2)
Definition: QuEST_cpu_internal.h:38
Complex::imag
qreal imag
Definition: QuEST.h:106
insertZeroBit
__forceinline__ __device__ long long int insertZeroBit(long long int number, int index)
Definition: QuEST_gpu.cu:99
Complex
Represents one complex number.
Definition: QuEST.h:103
densmatr_oneQubitDegradeOffDiagonal
void densmatr_oneQubitDegradeOffDiagonal(Qureg qureg, const int targetQubit, qreal retain)
Definition: QuEST_cpu.c:48