Unit test utilities

Functions used in the unit testing. These are mostly unoptimised, analytic implementations of the complex linear algebra that QuEST ultimately effects on quantum states. These are not part of the QuEST API, and require C++14. More...

Typedefs

typedef std::vector< std::vector< qcomp > > QMatrix
 A complex square matrix. More...
 
typedef std::vector< qcompQVector
 A complex vector, which can be zero-initialised with QVector(numAmps). More...
 

Functions

void applyReferenceMatrix (QMatrix &state, int *ctrls, int numCtrls, int *targs, int numTargs, QMatrix op)
 Modifies the density matrix state to be the result of left-multiplying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively). More...
 
void applyReferenceMatrix (QVector &state, int *ctrls, int numCtrls, int *targs, int numTargs, QMatrix op)
 Modifies the state-vector state to be the result of left-multiplying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively). More...
 
void applyReferenceOp (QMatrix &state, int *ctrls, int numCtrls, int *targs, int numTargs, QMatrix op)
 Modifies the density matrix state to be the result of applying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively). More...
 
void applyReferenceOp (QMatrix &state, int *ctrls, int numCtrls, int targ1, int targ2, QMatrix op)
 Modifies the density matrix state to be the result of applying the two-target operator matrix op, with the specified control qubits (in ctrls). More...
 
void applyReferenceOp (QMatrix &state, int *ctrls, int numCtrls, int target, QMatrix op)
 Modifies the density matrix state to be the result of applying the single-target operator matrix op, with the specified control qubits (in ctrls). More...
 
void applyReferenceOp (QMatrix &state, int *targs, int numTargs, QMatrix op)
 Modifies the density matrix state to be the result of applying the multi-target operator matrix op, with no control qubits. More...
 
void applyReferenceOp (QMatrix &state, int ctrl, int *targs, int numTargs, QMatrix op)
 Modifies the density matrix state to be the result of applying the multi-target operator matrix op, with a single control qubit ctrl. More...
 
void applyReferenceOp (QMatrix &state, int ctrl, int targ, QMatrix op)
 Modifies the density matrix state to be the result of applying the single-control single-target operator matrix op. More...
 
void applyReferenceOp (QMatrix &state, int ctrl, int targ1, int targ2, QMatrix op)
 Modifies the density matrix state to be the result of applying the two-target operator matrix op, with a single control qubit ctrl. More...
 
void applyReferenceOp (QMatrix &state, int targ, QMatrix op)
 Modifies the density matrix state to be the result of applying the single-target operator matrix op, with no control qubit. More...
 
void applyReferenceOp (QVector &state, int *ctrls, int numCtrls, int *targs, int numTargs, QMatrix op)
 Modifies the state-vector state to be the result of applying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively). More...
 
void applyReferenceOp (QVector &state, int *ctrls, int numCtrls, int targ1, int targ2, QMatrix op)
 Modifies the state-vector state to be the result of applying the two-target operator matrix op, with the specified control qubits (in ctrls). More...
 
void applyReferenceOp (QVector &state, int *ctrls, int numCtrls, int target, QMatrix op)
 Modifies the state-vector state to be the result of applying the single-target operator matrix op, with the specified control qubits (in ctrls). More...
 
void applyReferenceOp (QVector &state, int *targs, int numTargs, QMatrix op)
 Modifies the state-vector state to be the result of applying the multi-target operator matrix op, with no contorl qubits. More...
 
void applyReferenceOp (QVector &state, int ctrl, int *targs, int numTargs, QMatrix op)
 Modifies the state-vector state to be the result of applying the multi-target operator matrix op, with a single control qubit (ctrl) This updates state under. More...
 
void applyReferenceOp (QVector &state, int ctrl, int targ, QMatrix op)
 Modifies the state-vector state to be the result of applying the single-target operator matrix op, with a single control qubit (ctrl). More...
 
void applyReferenceOp (QVector &state, int ctrl, int targ1, int targ2, QMatrix op)
 Modifies the state-vector state to be the result of applying the two-target operator matrix op, with a single control qubit (ctrl). More...
 
void applyReferenceOp (QVector &state, int targ, QMatrix op)
 Modifies the state-vector state to be the result of applying the single-target operator matrix op, with no contorl qubits. More...
 
bool areEqual (QMatrix a, QMatrix b)
 Returns true if the absolute value of the difference between every amplitude in matrices a and b is less than REAL_EPS. More...
 
bool areEqual (Qureg qureg, QMatrix matr)
 Performs a hardware-agnostic comparison of density-matrix qureg to matr, checking whether the difference between the real and imaginary components of every amplitude is smaller than the QuEST_PREC-specific REAL_EPS (defined in QuEST_precision) precision. More...
 
bool areEqual (Qureg qureg, QMatrix matr, qreal precision)
 Performs a hardware-agnostic comparison of density-matrix qureg to matr, checking whether the difference between the real and imaginary components of every amplitude is smaller than precision. More...
 
bool areEqual (Qureg qureg, QVector vec)
 Performs a hardware-agnostic comparison of state-vector qureg to vec, checking whether the difference between the real and imaginary components of every amplitude is smaller than the QuEST_PREC-specific REAL_EPS (defined in QuEST_precision) precision. More...
 
bool areEqual (Qureg qureg, QVector vec, qreal precision)
 Performs a hardware-agnostic comparison of state-vector qureg to vec, checking whether the difference between the real and imaginary components of every amplitude is smaller than precision. More...
 
bool areEqual (Qureg qureg1, Qureg qureg2)
 Performs a hardware-agnostic comparison of the given quregs, checking whether the difference between the real and imaginary components of every amplitude is smaller than the QuEST_PREC-specific REAL_EPS (defined in QuEST_precision) precision. More...
 
bool areEqual (Qureg qureg1, Qureg qureg2, qreal precision)
 Performs a hardware-agnostic comparison of the given quregs, checking whether the difference between the real and imaginary components of every amplitude is smaller than precision. More...
 
bool areEqual (QVector a, QVector b)
 Returns true if the absolute value of the difference between every amplitude in vectors a and b is less than REAL_EPS. More...
 
bool areEqual (QVector vec, qreal *reals)
 Returns true if the absolute value of the difference between every element in vec (which must be strictly real) and those implied by reals, is less than REAL_EPS. More...
 
bool areEqual (QVector vec, qreal *reals, qreal *imags)
 Returns true if the absolute value of the difference between every element in vec and those implied by reals and imags, is less than REAL_EPS. More...
 
CatchGen< int * > bitsets (int numBits)
 Returns a Catch2 generator of every numBits-length bit-set, in increasing lexographic order, where left-most (zero index) bit is treated as LEAST significant (opposite typical convention). More...
 
unsigned int calcLog2 (long unsigned int res)
 Returns log2 of numbers which must be gauranteed to be 2^n. More...
 
void deleteFilesWithPrefixSynch (char *prefix)
 Deletes all files with filename starting with prefix. More...
 
QMatrix getConjugateTranspose (QMatrix a)
 Returns the conjugate transpose of the complex square matrix a. More...
 
QVector getDFT (QVector in)
 Returns the discrete fourier transform of vector in. More...
 
QVector getDFT (QVector in, int *targs, int numTargs)
 Returns the discrete fourier transform of a sub-partition of the vector in. More...
 
QMatrix getExponentialOfDiagonalMatrix (QMatrix a)
 Returns the matrix exponential of a diagonal, square, complex matrix. More...
 
QMatrix getExponentialOfPauliMatrix (qreal angle, QMatrix a)
 Returns the matrix exponential of a kronecker product of pauli matrices (or of any involutory matrices), with exponent factor (-i angle / 2). More...
 
QMatrix getFullOperatorMatrix (int *ctrls, int numCtrls, int *targs, int numTargs, QMatrix op, int numQubits)
 Takes a 2^numTargs-by-2^numTargs matrix op and a returns a 2^numQubits-by-2^numQubits matrix where op is controlled on the given ctrls qubits. More...
 
QMatrix getIdentityMatrix (size_t dim)
 Returns a dim-by-dim identity matrix. More...
 
QMatrix getKetBra (QVector ket, QVector bra)
 Returns the matrix |ket><bra|, with ith-jth element ket(i) conj(bra(j)), since |ket><bra| = sum_i a_i|i> sum_j b_j* <j| = sum_{ij} a_i b_j* |i><j|. More...
 
QMatrix getKroneckerProduct (QMatrix a, QMatrix b)
 Returns the kronecker product of a and b, where a and b are square but possibly differently-sized complex matrices. More...
 
QVector getKroneckerProduct (QVector b, QVector a)
 Returns b (otimes) a. More...
 
QVector getMatrixDiagonal (QMatrix matr)
 Returns the diagonal vector of the given matrix. More...
 
QMatrix getMixedDensityMatrix (std::vector< qreal > probs, std::vector< QVector > states)
 Returns a mixed density matrix formed from mixing the given pure states, which are assumed normalised, but not necessarily orthogonal. More...
 
QVector getNormalised (QVector vec)
 Returns an L2-normalised copy of vec, using Kahan summation for improved accuracy. More...
 
QMatrix getPureDensityMatrix (QVector state)
 Returns a density matrix initialised into the given pure state. More...
 
qcomp getRandomComplex ()
 Returns a random complex number within the square closing (-1-i) and (1+i), from a distribution uniformly randomising the individual real and imaginary components in their domains. More...
 
QMatrix getRandomDensityMatrix (int numQb)
 Returns a random numQb-by-numQb density matrix, from an undisclosed distribution, in a very mixed state. More...
 
int getRandomInt (int min, int max)
 Returns a random integer between min (inclusive) and max (exclusive), from the uniform distribution. More...
 
std::vector< QMatrixgetRandomKrausMap (int numQb, int numOps)
 Returns a random Kraus map of #numOps 2^numQb-by-2^numQb operators, from an undisclosed distribution. More...
 
std::vector< QVectorgetRandomOrthonormalVectors (int numQb, int numStates)
 Returns a list of random orthonormal complex vectors, from an undisclosed distribution. More...
 
std::vector< qrealgetRandomProbabilities (int numProbs)
 Returns a list of random real scalars, each in [0, 1], which sum to unity. More...
 
QMatrix getRandomPureDensityMatrix (int numQb)
 Returns a random numQb-by-numQb density matrix, from an undisclosed distribution, which is pure (corresponds to a random state-vector) More...
 
QMatrix getRandomQMatrix (int dim)
 Returns a dim-by-dim complex matrix, where the real and imaginary value of each element are independently random, under the standard normal distribution (mean 0, standard deviation 1). More...
 
QVector getRandomQVector (int dim)
 Returns a dim-length vector with random complex amplitudes in the square joining {-1-i, 1+i}, of an undisclosed distribution. More...
 
qreal getRandomReal (qreal min, qreal max)
 Returns a random real between min (inclusive) and max (exclusive), from the uniform distribution. More...
 
QVector getRandomStateVector (int numQb)
 Returns a random numQb-length L2-normalised state-vector from an undisclosed distribution. More...
 
QMatrix getRandomUnitary (int numQb)
 Returns a uniformly random (under Haar) 2^numQb-by-2^numQb unitary matrix. More...
 
QMatrix getSwapMatrix (int qb1, int qb2, int numQb)
 Returns the 2^numQb-by-2^numQb unitary matrix which swaps qubits qb1 and qb2; the SWAP gate of not-necessarily-adjacent qubits. More...
 
long long int getTwosComplement (long long int decimal, int numBits)
 Returns the two's complement signed encoding of the unsigned number decimal, which must be a number between 0 and 2^numBits (exclusive). More...
 
long long int getUnsigned (long long int twosComp, int numBits)
 Return the unsigned value of a number, made of #numBits bits, which under two's complement, encodes the signed number twosComp. More...
 
long long int getValueOfTargets (long long int ind, int *targs, int numTargs)
 Returns the integer value of the targeted sub-register for the given full state index ind. More...
 
QMatrix getZeroMatrix (size_t dim)
 Returns a dim-by-dim square complex matrix, initialised to all zeroes. More...
 
CatchGen< pauliOpType * > pauliseqs (int numPaulis)
 Returns a Catch2 generator of every numPaulis-length set of Pauli-matrix types (or base-4 integers). More...
 
CatchGen< int * > sequences (int base, int numDigits)
 Returns a Catch2 generator of every numDigits-length sequence in the given base, in increasing lexographic order, where left-most (zero index) bit is treated as LEAST significant (opposite typical convention). More...
 
void setDiagMatrixOverrides (QMatrix &matr, int *numQubitsPerReg, int numRegs, enum bitEncoding encoding, long long int *overrideInds, qreal *overridePhases, int numOverrides)
 Modifies the given diagonal matrix such that the diagonal elements which correspond to the coordinates in overrideInds are replaced with exp(i phase), as prescribed by overridePhases. More...
 
void setRandomDiagPauliHamil (PauliHamil hamil)
 Populates hamil with random coefficients and a random amount number of PAULI_I and PAULI_Z operators. More...
 
void setRandomPauliSum (PauliHamil hamil)
 Populates hamil with random coefficients and pauli codes. More...
 
void setRandomPauliSum (qreal *coeffs, pauliOpType *codes, int numQubits, int numTerms)
 Populates the coeffs array with random qreals in (-5, 5), and populates codes with random Pauli codes. More...
 
void setSubMatrix (QMatrix &dest, QMatrix sub, size_t r, size_t c)
 Modifies dest by overwriting its submatrix (from top-left corner (r, c) to bottom-right corner (r + dest.size(), c + dest.size()) with the complete elements of sub. More...
 
void setUniqueFilename (char *outFn, char *prefix)
 Modifies outFn to be a filename of format prefix_NUM.txt where NUM is a new unique integer so far. More...
 
CatchGen< int * > sublists (CatchGen< int > &&gen, int numSamps, const int *exclude, int numExclude)
 Returns a Catch2 generator of every length-sublen sublist of the elements generated by gen, which exclude all elements in exclude, in increasing lexographic order. More...
 
CatchGen< int * > sublists (CatchGen< int > &&gen, int numSamps, int excluded)
 Returns a Catch2 generator of every length-sublen sublist of the elements generated by gen which exclude element excluded, in increasing lexographic order. More...
 
CatchGen< int * > sublists (CatchGen< int > &&gen, int sublen)
 Returns a Catch2 generator of every length-sublen sublist of the elements generated by gen, in increasing lexographic order. More...
 
CatchGen< int * > sublists (int *list, int len, int sublen)
 Returns a Catch2 generator of every length-sublen sublist of length-len list, in increasing lexographic order. More...
 
ComplexMatrix2 toComplexMatrix2 (QMatrix qm)
 Returns a ComplexMatrix2 copy of QMatix qm. More...
 
ComplexMatrix4 toComplexMatrix4 (QMatrix qm)
 Returns a ComplexMatrix4 copy of QMatix qm. More...
 
void toComplexMatrixN (QMatrix qm, ComplexMatrixN cm)
 Initialises cm with the values of qm. More...
 
QMatrix toDiagonalQMatrix (QVector vec)
 Returns a diagonal complex matrix formed by the given vector. More...
 
QMatrix toQMatrix (Complex alpha, Complex beta)
 Returns the matrix (where a=alpha, b=beta) {{a, -conj(b)}, {b, conj(a)}} using the qcomp complex type. More...
 
QMatrix toQMatrix (ComplexMatrix2 src)
 Returns a copy of the given 2-by-2 matrix. More...
 
QMatrix toQMatrix (ComplexMatrix4 src)
 Returns a copy of the given 4-by-4 matrix. More...
 
QMatrix toQMatrix (ComplexMatrixN src)
 Returns a copy of the given 2^N-by-2^N matrix. More...
 
QMatrix toQMatrix (DiagonalOp op)
 Returns a 2^N-by-2^N complex diagonal matrix form of the DiagonalOp. More...
 
QMatrix toQMatrix (PauliHamil hamil)
 Returns a 2^N-by-2^N Hermitian matrix form of the PauliHamil. More...
 
QMatrix toQMatrix (qreal *coeffs, pauliOpType *paulis, int numQubits, int numTerms)
 Returns a 2^N-by-2^N Hermitian matrix form of the specified weighted sum of Pauli products. More...
 
QMatrix toQMatrix (Qureg qureg)
 Returns an equal-size copy of the given density matrix qureg. More...
 
void toQureg (Qureg qureg, QMatrix mat)
 Initialises the density matrix qureg to have the same amplitudes as mat. More...
 
void toQureg (Qureg qureg, QVector vec)
 Initialises the state-vector qureg to have the same amplitudes as vec. More...
 
QVector toQVector (DiagonalOp op)
 Returns a vector with the same of the full diagonal operator, populated with op's elements. More...
 
QVector toQVector (Qureg qureg)
 Returns an equal-size copy of the given state-vector qureg. More...
 
void writeToFileSynch (char *fn, const string &contents)
 Writes contents to the file with filename fn, which is created and/or overwritten. More...
 

Detailed Description

Functions used in the unit testing. These are mostly unoptimised, analytic implementations of the complex linear algebra that QuEST ultimately effects on quantum states. These are not part of the QuEST API, and require C++14.

Author
Tyson Jones

Typedef Documentation

◆ QMatrix

typedef std::vector<std::vector<qcomp> > QMatrix

A complex square matrix.

Should be initialised with getZeroMatrix(). These have all the natural linear-algebra operator overloads, including left-multiplication onto a vector.

This data-structure is not partitioned between nodes in distributed mode. That is, every node has a complete copy, allowing for safe comparisons.

Author
Tyson Jones

Definition at line 49 of file utilities.hpp.

◆ QVector

typedef std::vector<qcomp> QVector

A complex vector, which can be zero-initialised with QVector(numAmps).

These have all the natural linear-algebra operator overloads.

This data-structure is not partitioned between nodes in distributed mode. That is, every node has a complete copy, allowing for safe comparisons.

Author
Tyson Jones

Definition at line 60 of file utilities.hpp.

Function Documentation

◆ applyReferenceMatrix() [1/2]

void applyReferenceMatrix ( QMatrix state,
int *  ctrls,
int  numCtrls,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the density matrix state to be the result of left-multiplying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively).

Here, op is treated like a simple matrix and is hence left-multiplied onto the state once.

Author
Tyson Jones

Definition at line 847 of file utilities.cpp.

849  {
850  // for density matrices, op is left-multiplied only
851  int numQubits = calcLog2(state.size());
852  QMatrix leftOp = getFullOperatorMatrix(ctrls, numCtrls, targs, numTargs, op, numQubits);
853  state = leftOp * state;
854 }

References calcLog2(), and getFullOperatorMatrix().

◆ applyReferenceMatrix() [2/2]

void applyReferenceMatrix ( QVector state,
int *  ctrls,
int  numCtrls,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the state-vector state to be the result of left-multiplying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively).

This is an alias of applyReferenceOp(), since operators are always left-multiplied as matrices onto state-vectors.

Author
Tyson Jones

Definition at line 841 of file utilities.cpp.

843  {
844  // for state-vectors, the op is always just left-multiplied
845  applyReferenceOp(state, ctrls, numCtrls, targs, numTargs, op);
846 }

References applyReferenceOp().

Referenced by TEST_CASE().

◆ applyReferenceOp() [1/16]

void applyReferenceOp ( QMatrix state,
int *  ctrls,
int  numCtrls,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 2^numTargs-by-2^numTargs matrix. Furthermore, every element of targs must not appear in ctrls (and vice-versa), though this is not explicitly checked. Elements of targs and ctrls should be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 784 of file utilities.cpp.

786  {
787  int numQubits = calcLog2(state.size());
788  QMatrix leftOp = getFullOperatorMatrix(ctrls, numCtrls, targs, numTargs, op, numQubits);
789  QMatrix rightOp = getConjugateTranspose(leftOp);
790  state = leftOp * state * rightOp;
791 }

References calcLog2(), getConjugateTranspose(), and getFullOperatorMatrix().

◆ applyReferenceOp() [2/16]

void applyReferenceOp ( QMatrix state,
int *  ctrls,
int  numCtrls,
int  targ1,
int  targ2,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the two-target operator matrix op, with the specified control qubits (in ctrls).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 4-by-4 matrix. Both targ1 and targ2 must not appear in ctrls, though this is not explicitly checked. Elements of ctrls, and targ1 and targ2, should be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 792 of file utilities.cpp.

794  {
795  int targs[2] = {targ1, targ2};
796  applyReferenceOp(state, ctrls, numCtrls, targs, 2, op);
797 }

References applyReferenceOp().

◆ applyReferenceOp() [3/16]

void applyReferenceOp ( QMatrix state,
int *  ctrls,
int  numCtrls,
int  target,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the single-target operator matrix op, with the specified control qubits (in ctrls).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 2-by-2 matrix. target must not appear in ctrls, though this is not explicitly checked.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 798 of file utilities.cpp.

800  {
801  int targs[1] = {target};
802  applyReferenceOp(state, ctrls, numCtrls, targs, 1, op);
803 }

References applyReferenceOp().

◆ applyReferenceOp() [4/16]

void applyReferenceOp ( QMatrix state,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the multi-target operator matrix op, with no control qubits.

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 2^numTargs-by-2^numTargs matrix. Every element in targs should be unique, though this is not explicitly checked.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 804 of file utilities.cpp.

806  {
807  applyReferenceOp(state, NULL, 0, targs, numTargs, op);
808 }

References applyReferenceOp().

◆ applyReferenceOp() [5/16]

void applyReferenceOp ( QMatrix state,
int  ctrl,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the multi-target operator matrix op, with a single control qubit ctrl.

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 2^numTargs-by-2^numTargs matrix, and ctrl must not appear in targs (though this is not explicitly checked).

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 816 of file utilities.cpp.

818  {
819  int ctrls[1] = {ctrl};
820  applyReferenceOp(state, ctrls, 1, targs, numTargs, op);
821 }

References applyReferenceOp().

◆ applyReferenceOp() [6/16]

void applyReferenceOp ( QMatrix state,
int  ctrl,
int  targ,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the single-control single-target operator matrix op.

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 2-by-2 matrix, and ctrl and targ should be different.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 809 of file utilities.cpp.

811  {
812  int ctrls[1] = {ctrl};
813  int targs[1] = {targ};
814  applyReferenceOp(state, ctrls, 1, targs, 1, op);
815 }

References applyReferenceOp().

◆ applyReferenceOp() [7/16]

void applyReferenceOp ( QMatrix state,
int  ctrl,
int  targ1,
int  targ2,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the two-target operator matrix op, with a single control qubit ctrl.

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 4-by-4 matrix, and ctrl, targ1 and targ2 must be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 822 of file utilities.cpp.

824  {
825  int ctrls[1] = {ctrl};
826  int targs[2] = {targ1, targ2};
827  applyReferenceOp(state, ctrls, 1, targs, 2, op);
828 }

References applyReferenceOp().

◆ applyReferenceOp() [8/16]

void applyReferenceOp ( QMatrix state,
int  targ,
QMatrix  op 
)

Modifies the density matrix state to be the result of applying the single-target operator matrix op, with no control qubit.

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \, \text{op}^\dagger \]

even if op is not unitary (which is useful for applying Kraus operators).

op must be a 2-by-2 matrix.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multipling it to state, then right-multiplying its conjugate transpose onto the result.

Author
Tyson Jones

Definition at line 829 of file utilities.cpp.

831  {
832  int targs[1] = {targ};
833  applyReferenceOp(state, NULL, 0, targs, 1, op);
834 }

References applyReferenceOp().

◆ applyReferenceOp() [9/16]

void applyReferenceOp ( QVector state,
int *  ctrls,
int  numCtrls,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the multi-target operator matrix op, with the specified control and target qubits (in ctrls and targs respectively).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 2^numTargs-by-2^numTargs matrix. Furthermore, every element of targs must not appear in ctrls (and vice-versa), though this is not explicitly checked. Elements of targs and ctrls should be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 728 of file utilities.cpp.

730  {
731  int numQubits = calcLog2(state.size());
732  QMatrix fullOp = getFullOperatorMatrix(ctrls, numCtrls, targs, numTargs, op, numQubits);
733  state = fullOp * state;
734 }

References calcLog2(), and getFullOperatorMatrix().

Referenced by applyReferenceMatrix(), applyReferenceOp(), and TEST_CASE().

◆ applyReferenceOp() [10/16]

void applyReferenceOp ( QVector state,
int *  ctrls,
int  numCtrls,
int  targ1,
int  targ2,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the two-target operator matrix op, with the specified control qubits (in ctrls).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 4-by-4 matrix. Furthermore, ctrls, targ1 and targ2 should all be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 735 of file utilities.cpp.

737  {
738  int targs[2] = {targ1, targ2};
739  applyReferenceOp(state, ctrls, numCtrls, targs, 2, op);
740 }

References applyReferenceOp().

◆ applyReferenceOp() [11/16]

void applyReferenceOp ( QVector state,
int *  ctrls,
int  numCtrls,
int  target,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the single-target operator matrix op, with the specified control qubits (in ctrls).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 2-by-2 matrix. Furthermore, elements in ctrls and target should all be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 741 of file utilities.cpp.

743  {
744  int targs[1] = {target};
745  applyReferenceOp(state, ctrls, numCtrls, targs, 1, op);
746 }

References applyReferenceOp().

◆ applyReferenceOp() [12/16]

void applyReferenceOp ( QVector state,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the multi-target operator matrix op, with no contorl qubits.

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 2^numTargs-by-2^numTargs matrix. Furthermore, elements in targs should be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 747 of file utilities.cpp.

749  {
750  applyReferenceOp(state, NULL, 0, targs, numTargs, op);
751 }

References applyReferenceOp().

◆ applyReferenceOp() [13/16]

void applyReferenceOp ( QVector state,
int  ctrl,
int *  targs,
int  numTargs,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the multi-target operator matrix op, with a single control qubit (ctrl) This updates state under.

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 2^numTargs-by-2^numTargs matrix. Furthermore, elements in targs and ctrl should be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 759 of file utilities.cpp.

761  {
762  int ctrls[1] = {ctrl};
763  applyReferenceOp(state, ctrls, 1, targs, numTargs, op);
764 }

References applyReferenceOp().

◆ applyReferenceOp() [14/16]

void applyReferenceOp ( QVector state,
int  ctrl,
int  targ,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the single-target operator matrix op, with a single control qubit (ctrl).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 2-by-2 matrix. Furthermore, ctrl and targ must be different.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 752 of file utilities.cpp.

754  {
755  int ctrls[1] = {ctrl};
756  int targs[1] = {targ};
757  applyReferenceOp(state, ctrls, 1, targs, 1, op);
758 }

References applyReferenceOp().

◆ applyReferenceOp() [15/16]

void applyReferenceOp ( QVector state,
int  ctrl,
int  targ1,
int  targ2,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the two-target operator matrix op, with a single control qubit (ctrl).

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 4-by-4 matrix. Furthermore, ctrl, targ1 and targ2 should all be unique.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 765 of file utilities.cpp.

767  {
768  int ctrls[1] = {ctrl};
769  int targs[2] = {targ1, targ2};
770  applyReferenceOp(state, ctrls, 1, targs, 2, op);
771 }

References applyReferenceOp().

◆ applyReferenceOp() [16/16]

void applyReferenceOp ( QVector state,
int  targ,
QMatrix  op 
)

Modifies the state-vector state to be the result of applying the single-target operator matrix op, with no contorl qubits.

This updates state under

\[ \text{state} \to \text{op} \, \text{state} \]

even if op is not unitary.

op must be a 2-by-2 matrix.

This function works by computing getFullOperatorMatrix() from the given arguments, and left-multiplying it onto state.

Author
Tyson Jones

Definition at line 772 of file utilities.cpp.

774  {
775  int targs[1] = {targ};
776  applyReferenceOp(state, NULL, 0, targs, 1, op);
777 }

References applyReferenceOp().

◆ areEqual() [1/10]

bool areEqual ( QMatrix  a,
QMatrix  b 
)

Returns true if the absolute value of the difference between every amplitude in matrices a and b is less than REAL_EPS.

Author
Tyson Jones

Definition at line 407 of file utilities.cpp.

407  {
408  DEMAND( a.size() == b.size() );
409 
410  for (size_t i=0; i<a.size(); i++)
411  for (size_t j=0; j<b.size(); j++)
412  if (abs(a[i][j] - b[i][j]) > REAL_EPS)
413  return false;
414  return true;
415 }

References DEMAND.

◆ areEqual() [2/10]

bool areEqual ( Qureg  qureg,
QMatrix  matr 
)

Performs a hardware-agnostic comparison of density-matrix qureg to matr, checking whether the difference between the real and imaginary components of every amplitude is smaller than the QuEST_PREC-specific REAL_EPS (defined in QuEST_precision) precision.

This function demands qureg is a density matrix, and that qureg and matr have equal dimensions.

In GPU mode, this function involves a GPU to CPU memory copy overhead. In distributed mode, it involves a all-to-all single-int broadcast.

Author
Tyson Jones

Definition at line 983 of file utilities.cpp.

983  {
984  return areEqual(qureg, matr, REAL_EPS);
985 }

References areEqual().

◆ areEqual() [3/10]

bool areEqual ( Qureg  qureg,
QMatrix  matr,
qreal  precision 
)

Performs a hardware-agnostic comparison of density-matrix qureg to matr, checking whether the difference between the real and imaginary components of every amplitude is smaller than precision.

This function demands qureg is a density matrix, and that qureg and matr have equal dimensions.

In GPU mode, this function involves a GPU to CPU memory copy overhead. In distributed mode, it involves a all-to-all single-int broadcast.

Author
Tyson Jones

Definition at line 928 of file utilities.cpp.

928  {
929  DEMAND( qureg.isDensityMatrix );
930  DEMAND( (long long int) (matr.size()*matr.size()) == qureg.numAmpsTotal );
931 
932  // ensure local qureg.stateVec is up to date
933  copyStateFromGPU(qureg);
935 
936  // the starting index in vec of this node's qureg partition.
937  long long int startInd = qureg.chunkId * qureg.numAmpsPerChunk;
938  long long int globalInd, row, col, i;
939  int ampsAgree;
940 
941  // compare each of this node's amplitude to the corresponding matr sub-matrix
942  for (i=0; i<qureg.numAmpsPerChunk; i++) {
943  globalInd = startInd + i;
944  row = globalInd % matr.size();
945  col = globalInd / matr.size();
946  qreal realDif = absReal(qureg.stateVec.real[i] - real(matr[row][col]));
947  qreal imagDif = absReal(qureg.stateVec.imag[i] - imag(matr[row][col]));
948  ampsAgree = (realDif < precision && imagDif < precision);
949 
950  // DEBUG
951  if (!ampsAgree) {
952 
953  // debug
954  char buff[200];
955  sprintf(buff, "[msg from utilities.cpp] node %d has a disagreement at (global) index %lld of (%s) + i(%s)\n",
956  qureg.chunkId, globalInd, REAL_STRING_FORMAT, REAL_STRING_FORMAT);
957  printf(buff, realDif, imagDif);
958  }
959 
960  // break loop as soon as amplitudes disagree
961  if (!ampsAgree)
962  break;
963 
964  /* TODO:
965  * of the nodes which disagree, the lowest-rank should send its
966  * disagreeing (i, row, col, stateVec[i]) to rank 0 which should
967  * report it immediately (before the impending DEMAND failure)
968  * using FAIL_CHECK, so users can determine nature of disagreement
969  * (e.g. numerical precision).
970  * Note FAIL_CHECK accepts << like cout, e.g.
971  * FAIL_CHECK( "Amp at (" << row << ", " << col ") disagreed" );
972  */
973  }
974 
975  // if one node's partition wasn't equal, all-nodes must report not-equal
976  int allAmpsAgree = ampsAgree;
977 #ifdef DISTRIBUTED_MODE
978  MPI_Allreduce(&ampsAgree, &allAmpsAgree, 1, MPI_INT, MPI_LAND, MPI_COMM_WORLD);
979 #endif
980 
981  return allAmpsAgree;
982 }

References Qureg::chunkId, copyStateFromGPU(), DEMAND, Qureg::isDensityMatrix, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, qreal, QUEST_ENV, Qureg::stateVec, and syncQuESTEnv().

◆ areEqual() [4/10]

bool areEqual ( Qureg  qureg,
QVector  vec 
)

Performs a hardware-agnostic comparison of state-vector qureg to vec, checking whether the difference between the real and imaginary components of every amplitude is smaller than the QuEST_PREC-specific REAL_EPS (defined in QuEST_precision) precision.

This function demands qureg is a state-vector, and that qureg and vec have the same number of amplitudes.

In GPU mode, this function involves a GPU to CPU memory copy overhead. In distributed mode, it involves a all-to-all single-int broadcast.

Author
Tyson Jones

Definition at line 924 of file utilities.cpp.

924  {
925  return areEqual(qureg, vec, REAL_EPS);
926 }

References areEqual().

◆ areEqual() [5/10]

bool areEqual ( Qureg  qureg,
QVector  vec,
qreal  precision 
)

Performs a hardware-agnostic comparison of state-vector qureg to vec, checking whether the difference between the real and imaginary components of every amplitude is smaller than precision.

This function demands qureg is a state-vector, and that qureg and vec have the same number of amplitudes.

In GPU mode, this function involves a GPU to CPU memory copy overhead. In distributed mode, it involves a all-to-all single-int broadcast.

Author
Tyson Jones

Definition at line 883 of file utilities.cpp.

883  {
884  DEMAND( !qureg.isDensityMatrix );
885  DEMAND( (int) vec.size() == qureg.numAmpsTotal );
886 
887  copyStateFromGPU(qureg);
889 
890  // the starting index in vec of this node's qureg partition.
891  long long int startInd = qureg.chunkId * qureg.numAmpsPerChunk;
892 
893  int ampsAgree = 1;
894  for (long long int i=0; i<qureg.numAmpsPerChunk; i++) {
895  qreal realDif = absReal(qureg.stateVec.real[i] - real(vec[startInd+i]));
896  qreal imagDif = absReal(qureg.stateVec.imag[i] - imag(vec[startInd+i]));
897 
898  if (realDif > precision || imagDif > precision) {
899  ampsAgree = 0;
900 
901  // debug
902  char buff[200];
903  sprintf(buff, "Disagreement at %lld of (%s) + i(%s):\n\t%s + i(%s) VS %s + i(%s)\n",
904  startInd+i,
905  REAL_STRING_FORMAT, REAL_STRING_FORMAT, REAL_STRING_FORMAT,
906  REAL_STRING_FORMAT, REAL_STRING_FORMAT, REAL_STRING_FORMAT);
907  printf(buff,
908  realDif, imagDif,
909  qureg.stateVec.real[i], qureg.stateVec.imag[i],
910  real(vec[startInd+i]), imag(vec[startInd+i]));
911 
912  break;
913  }
914  }
915 
916  // if one node's partition wasn't equal, all-nodes must report not-equal
917  int allAmpsAgree = ampsAgree;
918 #ifdef DISTRIBUTED_MODE
919  MPI_Allreduce(&ampsAgree, &allAmpsAgree, 1, MPI_INT, MPI_LAND, MPI_COMM_WORLD);
920 #endif
921 
922  return allAmpsAgree;
923 }

References Qureg::chunkId, copyStateFromGPU(), DEMAND, Qureg::isDensityMatrix, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, qreal, QUEST_ENV, Qureg::stateVec, and syncQuESTEnv().

◆ areEqual() [6/10]

bool areEqual ( Qureg  qureg1,
Qureg  qureg2 
)

Performs a hardware-agnostic comparison of the given quregs, checking whether the difference between the real and imaginary components of every amplitude is smaller than the QuEST_PREC-specific REAL_EPS (defined in QuEST_precision) precision.

This function demands that qureg1 and qureg2 are of the same type (i.e. both state-vectors or both density matrices), and of an equal number of qubits.

In GPU mode, this function involves a GPU to CPU memory copy overhead. In distributed mode, it involves a all-to-all single-int broadcast.

Author
Tyson Jones

Definition at line 879 of file utilities.cpp.

879  {
880  return areEqual(qureg1, qureg2, REAL_EPS);
881 }

References areEqual().

◆ areEqual() [7/10]

bool areEqual ( Qureg  qureg1,
Qureg  qureg2,
qreal  precision 
)

Performs a hardware-agnostic comparison of the given quregs, checking whether the difference between the real and imaginary components of every amplitude is smaller than precision.

This function demands that qureg1 and qureg2 are of the same type (i.e. both state-vectors or both density matrices), and of an equal number of qubits.

In GPU mode, this function involves a GPU to CPU memory copy overhead. In distributed mode, it involves a all-to-all single-int broadcast.

Author
Tyson Jones

Definition at line 856 of file utilities.cpp.

856  {
857  DEMAND( qureg1.isDensityMatrix == qureg2.isDensityMatrix );
858  DEMAND( qureg1.numAmpsTotal == qureg2.numAmpsTotal );
859 
860  copyStateFromGPU(qureg1);
861  copyStateFromGPU(qureg2);
863 
864  // loop terminates when areEqual = 0
865  int ampsAgree = 1;
866  for (long long int i=0; ampsAgree && i<qureg1.numAmpsPerChunk; i++)
867  ampsAgree = (
868  absReal(qureg1.stateVec.real[i] - qureg2.stateVec.real[i]) < precision
869  && absReal(qureg1.stateVec.imag[i] - qureg2.stateVec.imag[i]) < precision);
870 
871  // if one node's partition wasn't equal, all-nodes must report not-equal
872  int allAmpsAgree = ampsAgree;
873 #ifdef DISTRIBUTED_MODE
874  MPI_Allreduce(&ampsAgree, &allAmpsAgree, 1, MPI_INT, MPI_LAND, MPI_COMM_WORLD);
875 #endif
876 
877  return allAmpsAgree;
878 }

References copyStateFromGPU(), DEMAND, Qureg::isDensityMatrix, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, QUEST_ENV, Qureg::stateVec, and syncQuESTEnv().

◆ areEqual() [8/10]

bool areEqual ( QVector  a,
QVector  b 
)

Returns true if the absolute value of the difference between every amplitude in vectors a and b is less than REAL_EPS.

Author
Tyson Jones

Definition at line 398 of file utilities.cpp.

398  {
399  DEMAND( a.size() == b.size() );
400 
401  for (size_t i=0; i<a.size(); i++)
402  if (abs(a[i] - b[i]) > REAL_EPS)
403  return false;
404  return true;
405 }

References DEMAND.

Referenced by areEqual(), getRandomKrausMap(), getRandomUnitary(), and TEST_CASE().

◆ areEqual() [9/10]

bool areEqual ( QVector  vec,
qreal reals 
)

Returns true if the absolute value of the difference between every element in vec (which must be strictly real) and those implied by reals, is less than REAL_EPS.

Author
Tyson Jones

Definition at line 1001 of file utilities.cpp.

1001  {
1002  for (size_t i=0; i<vec.size(); i++) {
1003  DEMAND( imag(vec[i]) == 0. );
1004 
1005  qreal dif = abs(real(vec[i]) - reals[i]);
1006  if (dif > REAL_EPS)
1007  return false;
1008  }
1009  return true;
1010 }

References DEMAND, and qreal.

◆ areEqual() [10/10]

bool areEqual ( QVector  vec,
qreal reals,
qreal imags 
)

Returns true if the absolute value of the difference between every element in vec and those implied by reals and imags, is less than REAL_EPS.

Author
Tyson Jones

Definition at line 987 of file utilities.cpp.

987  {
988 
989  qreal dif;
990  for (size_t i=0; i<vec.size(); i++) {
991  dif = absReal(real(vec[i]) - reals[i]);
992  if (dif > REAL_EPS)
993  return false;
994  dif = absReal(imag(vec[i]) - imags[i]);
995  if (dif > REAL_EPS)
996  return false;
997  }
998  return true;
999 }

References qreal.

◆ bitsets()

CatchGen<int*> bitsets ( int  numBits)

Returns a Catch2 generator of every numBits-length bit-set, in increasing lexographic order, where left-most (zero index) bit is treated as LEAST significant (opposite typical convention).

Note that the produced bitset must not be modified during generation.

This function can be used like

int* bits = GENERATE( bitsets(3) );

to produce {0,0,0}, {1,0,0}, {0,1,0}, {1,1,0}, {0,0,1}, {1,0,1}, {0,1,1}, {1,1,1}.

Author
Tyson Jones

Definition at line 1557 of file utilities.cpp.

1557  {
1558  return Catch::Generators::GeneratorWrapper<int*>(
1559  std::unique_ptr<Catch::Generators::IGenerator<int*>>(
1560  new SequenceGenerator<int>(1, numBits)));
1561 }

Referenced by TEST_CASE().

◆ calcLog2()

unsigned int calcLog2 ( long unsigned int  num)

Returns log2 of numbers which must be gauranteed to be 2^n.

Author
Tyson Jones

Returns log2 of numbers which must be gauranteed to be 2^n.

Definition at line 328 of file QuEST_validation.c.

328  {
329  unsigned int l = 0;
330  while (num >>= 1)
331  l++;
332  return l;
333 }

Referenced by applyReferenceMatrix(), applyReferenceOp(), TEST_CASE(), validateDiagPauliHamilFromFile(), validateNumQubitsInDiagOp(), and validateNumQubitsInQureg().

◆ deleteFilesWithPrefixSynch()

void deleteFilesWithPrefixSynch ( char *  prefix)

Deletes all files with filename starting with prefix.

In distributed mode, the master node deletes while the other nodes wait until complete.

Author
Tyson Jones

Definition at line 1375 of file utilities.cpp.

1375  {
1376 
1377  // master node deletes all files
1378  if (QUEST_ENV.rank == 0) {
1379  char cmd[200];
1380  sprintf(cmd, "exec rm %s*", prefix);
1381  system(cmd);
1382  }
1383 
1384  // other nodes wait
1386 }

References QUEST_ENV, QuESTEnv::rank, and syncQuESTEnv().

Referenced by TEST_CASE().

◆ getConjugateTranspose()

QMatrix getConjugateTranspose ( QMatrix  a)

Returns the conjugate transpose of the complex square matrix a.

Author
Tyson Jones

Definition at line 189 of file utilities.cpp.

189  {
190  QMatrix b = a;
191  for (size_t r=0; r<a.size(); r++)
192  for (size_t c=0; c<a.size(); c++)
193  b[r][c] = conj(a[c][r]);
194  return b;
195 }

Referenced by applyReferenceOp(), getRandomKrausMap(), and getRandomUnitary().

◆ getDFT() [1/2]

QVector getDFT ( QVector  in)

Returns the discrete fourier transform of vector in.

Author
Tyson Jones

Definition at line 652 of file utilities.cpp.

652  {
653  REQUIRE( in.size() > 0 );
654 
655  size_t dim = in.size();
656  qreal ampFac = 1 / sqrt( dim );
657  qreal phaseFac = 2 * M_PI / dim;
658 
659  QVector dftVec = QVector(dim);
660 
661  for (size_t x=0; x<dim; x++) {
662  dftVec[x] = 0;
663  for (long long int y=0; y<dim; y++)
664  dftVec[x] += expI(phaseFac * x * y) * in[y];
665  dftVec[x] *= ampFac;
666  }
667  return dftVec;
668 }

References expI(), M_PI, and qreal.

Referenced by TEST_CASE().

◆ getDFT() [2/2]

QVector getDFT ( QVector  in,
int *  targs,
int  numTargs 
)

Returns the discrete fourier transform of a sub-partition of the vector in.

Author
Tyson Jones

Definition at line 700 of file utilities.cpp.

700  {
701 
702  QVector out = QVector(in.size());
703  long long int targDim = (1LL << numTargs);
704 
705  for (size_t j=0; j<in.size(); j++) {
706 
707  // |j> = |x> (x) |...>, but mixed (not separated)
708  long long int x = getValueOfTargets(j, targs, numTargs);
709 
710  for (long long int y=0; y<targDim; y++) {
711 
712  // modifies sum_y |y> (x) ...
713  long long int outInd = getIndexOfTargetValues(j, targs, numTargs, y);
714 
715  qcomp elem = (in[j] / sqrt(pow(2,numTargs))) * expI(2*M_PI * x * y / pow(2,numTargs));
716  out[outInd] += elem;
717  }
718  }
719 
720  return out;
721 }

References expI(), getIndexOfTargetValues(), getValueOfTargets(), M_PI, and qcomp.

◆ getExponentialOfDiagonalMatrix()

QMatrix getExponentialOfDiagonalMatrix ( QMatrix  a)

Returns the matrix exponential of a diagonal, square, complex matrix.

This method explicitly checks that the passed matrix a is diagonal.

Author
Tyson Jones

Definition at line 197 of file utilities.cpp.

197  {
198 
199  // ensure diagonal
200  for (size_t r=0; r<a.size(); r++)
201  for (size_t c=0; c<a.size(); c++) {
202  if (r == c)
203  continue;
204  DEMAND( real(a[r][c]) == 0. );
205  DEMAND( imag(a[r][c]) == 0. );
206  }
207 
208  // exp(diagonal) = diagonal(exp)
209  QMatrix diag = a;
210  for (size_t i=0; i<a.size(); i++)
211  diag[i][i] = exp(diag[i][i]);
212 
213  return diag;
214 }

References DEMAND.

Referenced by TEST_CASE().

◆ getExponentialOfPauliMatrix()

QMatrix getExponentialOfPauliMatrix ( qreal  angle,
QMatrix  a 
)

Returns the matrix exponential of a kronecker product of pauli matrices (or of any involutory matrices), with exponent factor (-i angle / 2).

This method will not explicitly check that the passed matrix a is kronecker product of involutory matrices, but will otherwise return an incorrect exponential.

Author
Tyson Jones

Definition at line 216 of file utilities.cpp.

216  {
217  QMatrix iden = getIdentityMatrix(a.size());
218  QMatrix expo = (cos(angle/2) * iden) + ((qcomp) -1i * sin(angle/2) * a);
219  return expo;
220 }

References getIdentityMatrix(), and qcomp.

Referenced by TEST_CASE().

◆ getFullOperatorMatrix()

QMatrix getFullOperatorMatrix ( int *  ctrls,
int  numCtrls,
int *  targs,
int  numTargs,
QMatrix  op,
int  numQubits 
)

Takes a 2^numTargs-by-2^numTargs matrix op and a returns a 2^numQubits-by-2^numQubits matrix where op is controlled on the given ctrls qubits.

The union of {ctrls} and {targs} must be unique (though this is not explicitly checked), and every element must be >= 0 (not checked). The passed {ctrls} and {targs} arrays are unmodified.

This funciton works by first swapping {targs} and {ctrls} (via swap unitaries) to be strictly increasing {0,1,...}, building controlled(op), tensoring it to the full Hilbert space, and then 'unswapping'. The returned matrix has form: swap1 ... swapN . controlled(op) . swapN ... swap1

Author
Tyson Jones

Definition at line 304 of file utilities.cpp.

306  {
307  DEMAND( numCtrls >= 0 );
308  DEMAND( numTargs >= 0 );
309  DEMAND( numQubits >= (numCtrls+numTargs) );
310  DEMAND( op.size() == (1u << numTargs) );
311 
312  // copy {ctrls} and {targs}to restore at end
313  std::vector<int> ctrlsCopy(ctrls, ctrls+numCtrls);
314  std::vector<int> targsCopy(targs, targs+numTargs);
315 
316  // full-state matrix of qubit swaps
317  QMatrix swaps = getIdentityMatrix(1 << numQubits);
318  QMatrix unswaps = getIdentityMatrix(1 << numQubits);
319  QMatrix matr;
320 
321  // swap targs to {0, ..., numTargs-1}
322  for (int i=0; i<numTargs; i++) {
323  if (i != targs[i]) {
324  matr = getSwapMatrix(i, targs[i], numQubits);
325  swaps = matr * swaps;
326  unswaps = unswaps * matr;
327 
328  // even if this is the last targ, ctrls might still need updating
330  i, targs[i], (i < numTargs-1)? &targs[i+1] : NULL,
331  numTargs-i-1, ctrls, numCtrls);
332  }
333  }
334 
335  // swap ctrls to {numTargs, ..., numTargs+numCtrls-1}
336  for (int i=0; i<numCtrls; i++) {
337  int newInd = numTargs+i;
338  if (newInd != ctrls[i]) {
339  matr = getSwapMatrix(newInd, ctrls[i], numQubits);
340  swaps = matr * swaps;
341  unswaps = unswaps * matr;
342 
343  // update remaining ctrls (if any exist)
344  if (i < numCtrls-1)
345  updateIndices(newInd, ctrls[i], NULL, 0, &ctrls[i+1], numCtrls-i-1);
346  }
347  }
348 
349  // construct controlled-op matrix for qubits {0, ..., numCtrls+numTargs-1}
350  size_t dim = 1 << (numCtrls+numTargs);
351  QMatrix fullOp = getIdentityMatrix(dim);
352  setSubMatrix(fullOp, op, dim-op.size(), dim-op.size());
353 
354  // create full-state controlled-op matrix (left-pad identities)
355  if (numQubits > numCtrls+numTargs) {
356  size_t pad = 1 << (numQubits - numCtrls - numTargs);
357  fullOp = getKroneckerProduct(getIdentityMatrix(pad), fullOp);
358  }
359 
360  // apply swap to either side (to swap qubits back and forth)
361  fullOp = unswaps * fullOp * swaps;
362 
363  // restore {ctrls and targs}
364  for (int i=0; i<numCtrls; i++)
365  ctrls[i] = ctrlsCopy[i];
366  for (int i=0; i<numTargs; i++)
367  targs[i] = targsCopy[i];
368 
369  return fullOp;
370 }

References DEMAND, getIdentityMatrix(), getKroneckerProduct(), getSwapMatrix(), setSubMatrix(), and updateIndices().

Referenced by applyReferenceMatrix(), applyReferenceOp(), and TEST_CASE().

◆ getIdentityMatrix()

QMatrix getIdentityMatrix ( size_t  dim)

Returns a dim-by-dim identity matrix.

Author
Tyson Jones

Definition at line 161 of file utilities.cpp.

161  {
162  DEMAND( dim > 1 );
163  QMatrix matr = getZeroMatrix(dim);
164  for (size_t i=0; i<dim; i++)
165  matr[i][i] = 1;
166  return matr;
167 }

References DEMAND, and getZeroMatrix().

Referenced by getExponentialOfPauliMatrix(), getFullOperatorMatrix(), getRandomKrausMap(), getRandomUnitary(), and getSwapMatrix().

◆ getKetBra()

QMatrix getKetBra ( QVector  ket,
QVector  bra 
)

Returns the matrix |ket><bra|, with ith-jth element ket(i) conj(bra(j)), since |ket><bra| = sum_i a_i|i> sum_j b_j* <j| = sum_{ij} a_i b_j* |i><j|.

The dimensions of bra and ket must agree, and the returned square complex matrix has dimensions size(bra) x size(bra).

Author
Tyson Jones

Definition at line 169 of file utilities.cpp.

169  {
170  DEMAND( ket.size() == bra.size() );
171  QMatrix mat = getZeroMatrix(ket.size());
172 
173  for (size_t r=0; r<ket.size(); r++)
174  for (size_t c=0; c<ket.size(); c++)
175  mat[r][c] = ket[r] * conj(bra[c]);
176  return mat;
177 }

References DEMAND, and getZeroMatrix().

Referenced by getPureDensityMatrix(), getRandomDensityMatrix(), and TEST_CASE().

◆ getKroneckerProduct() [1/2]

QMatrix getKroneckerProduct ( QMatrix  a,
QMatrix  b 
)

Returns the kronecker product of a and b, where a and b are square but possibly differently-sized complex matrices.

Author
Tyson Jones

Definition at line 179 of file utilities.cpp.

179  {
180  QMatrix prod = getZeroMatrix(a.size() * b.size());
181  for (size_t r=0; r<b.size(); r++)
182  for (size_t c=0; c<b.size(); c++)
183  for (size_t i=0; i<a.size(); i++)
184  for (size_t j=0; j<a.size(); j++)
185  prod[r+b.size()*i][c+b.size()*j] = a[i][j] * b[r][c];
186  return prod;
187 }

References getZeroMatrix().

◆ getKroneckerProduct() [2/2]

QVector getKroneckerProduct ( QVector  b,
QVector  a 
)

Returns b (otimes) a.

If b and a are state-vectors, the resulting kronecker product is the seperable state formed by joining the qubits in the state-vectors, producing |b>|a> (a is least significant)

Author
Tyson Jones

Definition at line 143 of file utilities.cpp.

143  {
144 
145  QVector prod = QVector(a.size() * b.size());
146 
147  for (size_t i=0; i<prod.size(); i++)
148  prod[i] = b[i / a.size()] * a[i % a.size()];
149 
150  return prod;
151 }

Referenced by getFullOperatorMatrix(), getSwapMatrix(), TEST_CASE(), and toQMatrix().

◆ getMatrixDiagonal()

QVector getMatrixDiagonal ( QMatrix  matr)

Returns the diagonal vector of the given matrix.

Author
Tyson Jones

Definition at line 517 of file utilities.cpp.

517  {
518 
519  QVector vec = QVector(matr.size());
520  for (size_t i=0; i<vec.size(); i++)
521  vec[i] = matr[i][i];
522 
523  return vec;
524 }

◆ getMixedDensityMatrix()

QMatrix getMixedDensityMatrix ( std::vector< qreal probs,
std::vector< QVector states 
)

Returns a mixed density matrix formed from mixing the given pure states, which are assumed normalised, but not necessarily orthogonal.

Author
Tyson Jones

Definition at line 640 of file utilities.cpp.

640  {
641  DEMAND( probs.size() == states.size() );
642  DEMAND( probs.size() >= 1 );
643 
644  QMatrix matr = getZeroMatrix(states[0].size());
645 
646  for (size_t i=0; i<probs.size(); i++)
647  matr += probs[i] * getPureDensityMatrix(states[i]);
648 
649  return matr;
650 }

References DEMAND, getPureDensityMatrix(), and getZeroMatrix().

Referenced by TEST_CASE().

◆ getNormalised()

QVector getNormalised ( QVector  vec)

Returns an L2-normalised copy of vec, using Kahan summation for improved accuracy.

Author
Tyson Jones

Definition at line 446 of file utilities.cpp.

446  {
447  qreal norm = 0;
448  qreal y, t, c;
449  c = 0;
450 
451  for (size_t i=0; i<vec.size(); i++) {
452  y = real(vec[i])*real(vec[i]) - c;
453  t = norm + y;
454  c = ( t - norm ) - y;
455  norm = t;
456 
457  y = imag(vec[i])*imag(vec[i]) - c;
458  t = norm + y;
459  c = ( t - norm ) - y;
460  norm = t;
461  }
462 
463  for (size_t i=0; i<vec.size(); i++)
464  vec[i] /= sqrt(norm);
465  return vec;
466 }

References qreal.

Referenced by getRandomOrthonormalVectors(), and getRandomStateVector().

◆ getPureDensityMatrix()

QMatrix getPureDensityMatrix ( QVector  state)

Returns a density matrix initialised into the given pure state.

Author
Tyson Jones

Definition at line 507 of file utilities.cpp.

507  {
508  return getKetBra(state, state);
509 }

References getKetBra().

Referenced by getMixedDensityMatrix(), getRandomPureDensityMatrix(), and TEST_CASE().

◆ getRandomComplex()

qcomp getRandomComplex ( )

Returns a random complex number within the square closing (-1-i) and (1+i), from a distribution uniformly randomising the individual real and imaginary components in their domains.

Author
Tyson Jones

Definition at line 431 of file utilities.cpp.

431  {
432  return getRandomReal(-1,1) + getRandomReal(-1,1) * (qcomp) 1i;
433 }

References getRandomReal(), and qcomp.

Referenced by getRandomQVector(), and TEST_CASE().

◆ getRandomDensityMatrix()

QMatrix getRandomDensityMatrix ( int  numQb)

Returns a random numQb-by-numQb density matrix, from an undisclosed distribution, in a very mixed state.

This function works by generating 2^numQb random pure states, and mixing them with random probabilities.

Author
Tyson Jones

Definition at line 490 of file utilities.cpp.

490  {
491  DEMAND( numQb > 0 );
492 
493  // generate random probabilities to weight random pure states
494  int dim = 1<<numQb;
495  std::vector<qreal> probs = getRandomProbabilities(dim);
496 
497  // add random pure states
498  QMatrix dens = getZeroMatrix(dim);
499  for (int i=0; i<dim; i++) {
500  QVector pure = getRandomStateVector(numQb);
501  dens += probs[i] * getKetBra(pure, pure);
502  }
503 
504  return dens;
505 }

References DEMAND, getKetBra(), getRandomProbabilities(), getRandomStateVector(), and getZeroMatrix().

Referenced by TEST_CASE().

◆ getRandomInt()

int getRandomInt ( int  min,
int  max 
)

Returns a random integer between min (inclusive) and max (exclusive), from the uniform distribution.

Demands that max > min.

Author
Tyson Jones

Definition at line 526 of file utilities.cpp.

526  {
527  return round(getRandomReal(min, max-1));
528 }

References getRandomReal().

Referenced by setRandomPauliSum(), and TEST_CASE().

◆ getRandomKrausMap()

std::vector<QMatrix> getRandomKrausMap ( int  numQb,
int  numOps 
)

Returns a random Kraus map of #numOps 2^numQb-by-2^numQb operators, from an undisclosed distribution.

Note this method is very simple and cannot generate all possible Kraus maps. It works by generating numOps random unitary matrices, and randomly re-normalising them, such that the sum of ops[j]^dagger ops[j] = 1

Author
Tyson Jones

Definition at line 578 of file utilities.cpp.

578  {
579  DEMAND( numOps >= 1 );
580  DEMAND( numOps <= 4*numQb*numQb );
581 
582  // generate random unitaries
583  std::vector<QMatrix> ops;
584  for (int i=0; i<numOps; i++)
585  ops.push_back(getRandomUnitary(numQb));
586 
587  // generate random weights
588  qreal weights[numOps];
589  for (int i=0; i<numOps; i++)
590  weights[i] = getRandomReal(0, 1);
591 
592  // normalise random weights
593  qreal weightSum = 0;
594  for (int i=0; i<numOps; i++)
595  weightSum += weights[i];
596  for (int i=0; i<numOps; i++)
597  weights[i] = sqrt(weights[i]/weightSum);
598 
599  // normalise ops
600  for (int i=0; i<numOps; i++)
601  ops[i] *= weights[i];
602 
603  // check what we produced was a valid Kraus map
604  QMatrix iden = getIdentityMatrix(1 << numQb);
605  QMatrix prodSum = getZeroMatrix(1 << numQb);
606  for (int i=0; i<numOps; i++)
607  prodSum += getConjugateTranspose(ops[i]) * ops[i];
608  DEMAND( areEqual(prodSum, iden) );
609 
610  return ops;
611 }

References areEqual(), DEMAND, getConjugateTranspose(), getIdentityMatrix(), getRandomReal(), getRandomUnitary(), getZeroMatrix(), and qreal.

Referenced by TEST_CASE().

◆ getRandomOrthonormalVectors()

std::vector<QVector> getRandomOrthonormalVectors ( int  numQb,
int  numStates 
)

Returns a list of random orthonormal complex vectors, from an undisclosed distribution.

Author
Tyson Jones

Definition at line 613 of file utilities.cpp.

613  {
614  DEMAND( numQb >= 1 );
615  DEMAND( numStates >= 1);
616 
617  // set of orthonormal vectors
618  std::vector<QVector> vecs;
619 
620  for (size_t n=0; n<numStates; n++) {
621 
622  QVector vec = getRandomStateVector(numQb);
623 
624  // orthogonalise by substracting projections of existing vectors
625  for (int m=0; m<n; m++) {
626  qcomp prod = vec * vecs[m];
627  vec -= (prod * vecs[m]);
628  }
629 
630  // renormalise
631  vec = getNormalised(vec);
632 
633  // add to orthonormal set
634  vecs.push_back(vec);
635  }
636 
637  return vecs;
638 }

References DEMAND, getNormalised(), getRandomStateVector(), and qcomp.

Referenced by TEST_CASE().

◆ getRandomProbabilities()

std::vector<qreal> getRandomProbabilities ( int  numProbs)

Returns a list of random real scalars, each in [0, 1], which sum to unity.

Author
Tyson Jones

Definition at line 472 of file utilities.cpp.

472  {
473 
474  // generate random unnormalised scalars
475  std::vector<qreal> probs;
476  qreal total = 0;
477  for (int i=0; i<numProbs; i++) {
478  qreal prob = getRandomReal(0, 1);
479  probs.push_back(prob);
480  total += prob;
481  }
482 
483  // normalise
484  for (int i=0; i<numProbs; i++)
485  probs[i] /= total;
486 
487  return probs;
488 }

References getRandomReal(), and qreal.

Referenced by getRandomDensityMatrix(), and TEST_CASE().

◆ getRandomPureDensityMatrix()

QMatrix getRandomPureDensityMatrix ( int  numQb)

Returns a random numQb-by-numQb density matrix, from an undisclosed distribution, which is pure (corresponds to a random state-vector)

Author
Tyson Jones

Definition at line 511 of file utilities.cpp.

511  {
512  QVector vec = getRandomStateVector(numQb);
513  QMatrix mat = getPureDensityMatrix(vec);
514  return mat;
515 }

References getPureDensityMatrix(), and getRandomStateVector().

◆ getRandomQMatrix()

QMatrix getRandomQMatrix ( int  dim)

Returns a dim-by-dim complex matrix, where the real and imaginary value of each element are independently random, under the standard normal distribution (mean 0, standard deviation 1).

Author
Tyson Jones

Definition at line 379 of file utilities.cpp.

379  {
380  DEMAND( dim > 1 );
381 
382  QMatrix matr = getZeroMatrix(dim);
383  for (int i=0; i<dim; i++) {
384  for (int j=0; j<dim; j++) {
385 
386  // generate 2 normally-distributed random numbers via Box-Muller
387  qreal a = rand()/(qreal) RAND_MAX;
388  qreal b = rand()/(qreal) RAND_MAX;
389  qreal r1 = sqrt(-2 * log(a)) * cos(2 * 3.14159265 * b);
390  qreal r2 = sqrt(-2 * log(a)) * sin(2 * 3.14159265 * b);
391 
392  matr[i][j] = r1 + r2 * (qcomp) 1i;
393  }
394  }
395  return matr;
396 }

References DEMAND, getZeroMatrix(), qcomp, and qreal.

Referenced by getRandomUnitary(), and TEST_CASE().

◆ getRandomQVector()

QVector getRandomQVector ( int  dim)

Returns a dim-length vector with random complex amplitudes in the square joining {-1-i, 1+i}, of an undisclosed distribution.

The resulting vector is NOT L2-normalised.

Author
Tyson Jones

Definition at line 435 of file utilities.cpp.

435  {
436  QVector vec = QVector(dim);
437  for (int i=0; i<dim; i++)
438  vec[i] = getRandomComplex();
439 
440  // check we didn't get the impossibly-unlikely zero-amplitude outcome
441  DEMAND( real(vec[0]) != 0 );
442 
443  return vec;
444 }

References DEMAND, and getRandomComplex().

Referenced by getRandomStateVector(), and TEST_CASE().

◆ getRandomReal()

qreal getRandomReal ( qreal  min,
qreal  max 
)

Returns a random real between min (inclusive) and max (exclusive), from the uniform distribution.

Demands that max > min.

Author
Tyson Jones

Definition at line 421 of file utilities.cpp.

421  {
422  DEMAND( min <= max );
423  qreal r = min + (max - min) * (rand() / (qreal) RAND_MAX);
424 
425  // check bounds satisfied
426  DEMAND( r >= min );
427  DEMAND( r <= max );
428  return r;
429 }

References DEMAND, and qreal.

Referenced by getRandomComplex(), getRandomInt(), getRandomKrausMap(), getRandomProbabilities(), setRandomDiagPauliHamil(), setRandomPauliSum(), and TEST_CASE().

◆ getRandomStateVector()

QVector getRandomStateVector ( int  numQb)

Returns a random numQb-length L2-normalised state-vector from an undisclosed distribution.

This function works by randomly generating each complex amplitude, then L2-normalising.

Author
Tyson Jones

Definition at line 468 of file utilities.cpp.

468  {
469  return getNormalised(getRandomQVector(1<<numQb));
470 }

References getNormalised(), and getRandomQVector().

Referenced by getRandomDensityMatrix(), getRandomOrthonormalVectors(), getRandomPureDensityMatrix(), and TEST_CASE().

◆ getRandomUnitary()

QMatrix getRandomUnitary ( int  numQb)

Returns a uniformly random (under Haar) 2^numQb-by-2^numQb unitary matrix.

This function works by first generating a complex matrix where each element is independently random; the real and imaginary component thereof are independent standard normally-distributed (mean 0, standard-dev 1). Then, the matrix is orthonormalised via the Gram Schmidt algorithm. The resulting unitary matrix MAY be uniformly distributed under the Haar measure, but we make no assurance. This routine may return an identity matrix if it was unable to sufficiently precisely produce a unitary of the given size.

Author
Tyson Jones

Definition at line 530 of file utilities.cpp.

530  {
531  DEMAND( numQb >= 1 );
532 
533  QMatrix matr = getRandomQMatrix(1 << numQb);
534 
535  for (size_t i=0; i<matr.size(); i++) {
536  QVector row = matr[i];
537 
538  // compute new orthogonal row by subtracting proj row onto prevs
539  for (int k=i-1; k>=0; k--) {
540 
541  // compute row . prev = sum_n row_n conj(prev_n)
542  qcomp prod = 0;
543  for (size_t n=0; n<row.size(); n++)
544  prod += row[n] * conj(matr[k][n]);
545 
546  // subtract (proj row onto prev) = (prod * prev) from final row
547  for (size_t n=0; n<row.size(); n++)
548  matr[i][n] -= prod * matr[k][n];
549  }
550 
551  // compute row magnitude
552  qreal mag = 0;
553  for (size_t j=0; j<row.size(); j++)
554  mag += pow(abs(matr[i][j]), 2);
555  mag = sqrt(mag);
556 
557  // normalise row
558  for (size_t j=0; j<row.size(); j++)
559  matr[i][j] /= mag;
560  }
561 
562  // ensure matrix is indeed unitary
563  QMatrix conjprod = matr * getConjugateTranspose(matr);
564  QMatrix iden = getIdentityMatrix(1 << numQb);
565 
566  // generating big unitary matrices is hard; if we fail, default to identity
567  if ( numQb >= 3 && !areEqual(conjprod, iden) ) {
568 
569  matr = getIdentityMatrix(1 << numQb);
570  conjprod = matr;
571  }
572  DEMAND( areEqual(conjprod, iden) );
573 
574  // return the new orthonormal matrix
575  return matr;
576 }

References areEqual(), DEMAND, getConjugateTranspose(), getIdentityMatrix(), getRandomQMatrix(), qcomp, and qreal.

Referenced by getRandomKrausMap(), and TEST_CASE().

◆ getSwapMatrix()

QMatrix getSwapMatrix ( int  qb1,
int  qb2,
int  numQb 
)

Returns the 2^numQb-by-2^numQb unitary matrix which swaps qubits qb1 and qb2; the SWAP gate of not-necessarily-adjacent qubits.

If qb1 == qb2, returns the identity matrix.

Author
Tyson Jones

Definition at line 230 of file utilities.cpp.

230  {
231  DEMAND( numQb > 1 );
232  DEMAND( (qb1 >= 0 && qb1 < numQb) );
233  DEMAND( (qb2 >= 0 && qb2 < numQb) );
234 
235  if (qb1 > qb2)
236  std::swap(qb1, qb2);
237 
238  if (qb1 == qb2)
239  return getIdentityMatrix(1 << numQb);
240 
241  QMatrix swap;
242 
243  if (qb2 == qb1 + 1) {
244  // qubits are adjacent
245  swap = QMatrix{{1,0,0,0},{0,0,1,0},{0,1,0,0},{0,0,0,1}};
246 
247  } else {
248  // qubits are distant
249  int block = 1 << (qb2 - qb1);
250  swap = getZeroMatrix(block*2);
251  QMatrix iden = getIdentityMatrix(block/2);
252 
253  // Lemma 3.1 of arxiv.org/pdf/1711.09765.pdf
254  QMatrix p0{{1,0},{0,0}};
255  QMatrix l0{{0,1},{0,0}};
256  QMatrix l1{{0,0},{1,0}};
257  QMatrix p1{{0,0},{0,1}};
258 
259  /* notating a^(n+1) = identity(1<<n) (otimes) a, we construct the matrix
260  * [ p0^(N) l1^N ]
261  * [ l0^(N) p1^N ]
262  * where N = qb2 - qb1 */
263  setSubMatrix(swap, getKroneckerProduct(iden, p0), 0, 0);
264  setSubMatrix(swap, getKroneckerProduct(iden, l0), block, 0);
265  setSubMatrix(swap, getKroneckerProduct(iden, l1), 0, block);
266  setSubMatrix(swap, getKroneckerProduct(iden, p1), block, block);
267  }
268 
269  // pad swap with outer identities
270  if (qb1 > 0)
271  swap = getKroneckerProduct(swap, getIdentityMatrix(1<<qb1));
272  if (qb2 < numQb-1)
273  swap = getKroneckerProduct(getIdentityMatrix(1<<(numQb-qb2-1)), swap);
274 
275  return swap;
276 }

References DEMAND, getIdentityMatrix(), getKroneckerProduct(), getZeroMatrix(), and setSubMatrix().

Referenced by getFullOperatorMatrix().

◆ getTwosComplement()

long long int getTwosComplement ( long long int  decimal,
int  numBits 
)

Returns the two's complement signed encoding of the unsigned number decimal, which must be a number between 0 and 2^numBits (exclusive).

The returned number lies in [-2^(numBits-1), 2^(numBits-1)-1]

Author
Tyson Jones

Definition at line 1286 of file utilities.cpp.

1286  {
1287  DEMAND( decimal >= 0 );
1288  DEMAND( numBits >= 2 );
1289  DEMAND( decimal < (1LL << numBits) );
1290 
1291  long long int maxMag = 1LL << (numBits-1);
1292  if (decimal >= maxMag)
1293  return -maxMag + (decimal - maxMag);
1294  else
1295  return decimal;
1296 }

References DEMAND.

Referenced by TEST_CASE().

◆ getUnsigned()

long long int getUnsigned ( long long int  twosComp,
int  numBits 
)

Return the unsigned value of a number, made of #numBits bits, which under two's complement, encodes the signed number twosComp.

The returned number lies in [0, 2^(numBits)-1]

Author
Tyson Jones

Definition at line 1298 of file utilities.cpp.

1298  {
1299  DEMAND( numBits >= 2 );
1300  DEMAND( twosComp < (1LL << (numBits-1)) );
1301  DEMAND( twosComp >= - (1LL << (numBits-1)) );
1302 
1303  if (twosComp >= 0)
1304  return twosComp;
1305  else
1306  return (1<<numBits) + twosComp;
1307 }

References DEMAND.

Referenced by setDiagMatrixOverrides().

◆ getValueOfTargets()

long long int getValueOfTargets ( long long int  ind,
int *  targs,
int  numTargs 
)

Returns the integer value of the targeted sub-register for the given full state index ind.

Author
Tyson Jones

Definition at line 670 of file utilities.cpp.

670  {
671  DEMAND( ind >= 0 );
672 
673  long long int val = 0;
674 
675  for (int t=0; t<numTargs; t++)
676  val += ((ind >> targs[t]) & 1) * (1LL << t);
677 
678  return val;
679 }

References DEMAND.

Referenced by getDFT().

◆ getZeroMatrix()

QMatrix getZeroMatrix ( size_t  dim)

Returns a dim-by-dim square complex matrix, initialised to all zeroes.

Author
Tyson Jones

Definition at line 153 of file utilities.cpp.

153  {
154  DEMAND( dim > 1 );
155  QMatrix matr = QMatrix(dim);
156  for (size_t i=0; i<dim; i++)
157  matr[i].resize(dim);
158  return matr;
159 }

References DEMAND.

Referenced by getIdentityMatrix(), getKetBra(), getKroneckerProduct(), getMixedDensityMatrix(), getRandomDensityMatrix(), getRandomKrausMap(), getRandomQMatrix(), getSwapMatrix(), TEST_CASE(), toDiagonalQMatrix(), and toQMatrix().

◆ pauliseqs()

CatchGen<pauliOpType*> pauliseqs ( int  numPaulis)

Returns a Catch2 generator of every numPaulis-length set of Pauli-matrix types (or base-4 integers).

Generates in increasing lexographic order, where the left-most (zero index) pauli is treated as LEAST significant. Note that the sequence must not be modified during generation.

This function can be used like

pauliOpType* set = GENERATE( pauliseqs(2) );

to produce {I,I}, {X,I}, {Y,I}, {Z,I}, {I,X}, {X,X}, {Y,X}, {Z,X}, {I,Y}, {X,Y}, {Y,Y}, {Z,Y}, {I,Z}, {X,Z}, {Y,Z}, {Z,Z}/

Author
Tyson Jones

Definition at line 1567 of file utilities.cpp.

1567  {
1568  return Catch::Generators::GeneratorWrapper<pauliOpType*>(
1569  std::unique_ptr<Catch::Generators::IGenerator<pauliOpType*>>(
1570  new SequenceGenerator<pauliOpType>(PAULI_Z, numPaulis)));
1571 }

References PAULI_Z.

◆ sequences()

CatchGen<int*> sequences ( int  base,
int  numDigits 
)

Returns a Catch2 generator of every numDigits-length sequence in the given base, in increasing lexographic order, where left-most (zero index) bit is treated as LEAST significant (opposite typical convention).

Note that the sequence must not be modified during generation.

This function can be used like

int base = 3;
int numDigits = 2;
int* seq = GENERATE_COPY( sequences(base, numDigits) );

to produce {0,0}, {1,0}, {2,0}, {0,1}, {1,1}, {2,1}, {0,2}, {1,2}, {2,2}.

Author
Tyson Jones

Definition at line 1562 of file utilities.cpp.

1562  {
1563  return Catch::Generators::GeneratorWrapper<int*>(
1564  std::unique_ptr<Catch::Generators::IGenerator<int*>>(
1565  new SequenceGenerator<int>(base-1, numDigits)));
1566 }

◆ setDiagMatrixOverrides()

void setDiagMatrixOverrides ( QMatrix matr,
int *  numQubitsPerReg,
int  numRegs,
enum bitEncoding  encoding,
long long int *  overrideInds,
qreal overridePhases,
int  numOverrides 
)

Modifies the given diagonal matrix such that the diagonal elements which correspond to the coordinates in overrideInds are replaced with exp(i phase), as prescribed by overridePhases.

This function assumes that the given registers are contiguous, are in order of increasing significance, and that the matrix is proportionately sized and structured to act on the space of all registers combined. Overrides can be repeated, and only the first encountered for a given index will be effected (much like applyMultiVarPhaseFuncOverrides()).

Author
Tyson Jones

Definition at line 1316 of file utilities.cpp.

1316  {
1317  DEMAND( (encoding == UNSIGNED || encoding == TWOS_COMPLEMENT) );
1318  DEMAND( numRegs > 0 );
1319  DEMAND( numOverrides >= 0 );
1320 
1321  int totalQb = 0;
1322  for (int r=0; r<numRegs; r++) {
1323  DEMAND( numQubitsPerReg[r] > 0 );
1324  totalQb += numQubitsPerReg[r];
1325  }
1326  DEMAND( matr.size() == (1 << totalQb) );
1327 
1328  // record whether a diagonal index has been already overriden
1329  int hasBeenOverriden[1 << totalQb];
1330  for (int i=0; i<(1 << totalQb); i++)
1331  hasBeenOverriden[i] = 0;
1332 
1333  int flatInd = 0;
1334  for (int v=0; v<numOverrides; v++) {
1335  int matrInd = 0;
1336  int numQubitsLeft = 0;
1337 
1338  for (int r=0; r<numRegs; r++) {
1339 
1340  if (encoding == UNSIGNED)
1341  matrInd += overrideInds[flatInd] * (1 << numQubitsLeft);
1342  else if (encoding == TWOS_COMPLEMENT)
1343  matrInd += getUnsigned(overrideInds[flatInd], numQubitsPerReg[r]) * (1 << numQubitsLeft);
1344 
1345  numQubitsLeft += numQubitsPerReg[r];
1346  flatInd += 1;
1347  }
1348 
1349  if (!hasBeenOverriden[matrInd]) {
1350  matr[matrInd][matrInd] = expI(overridePhases[v]);
1351  hasBeenOverriden[matrInd] = 1;
1352  }
1353  }
1354 }

References DEMAND, expI(), getUnsigned(), TWOS_COMPLEMENT, and UNSIGNED.

Referenced by TEST_CASE().

◆ setRandomDiagPauliHamil()

void setRandomDiagPauliHamil ( PauliHamil  hamil)

Populates hamil with random coefficients and a random amount number of PAULI_I and PAULI_Z operators.

Author
Tyson Jones

Definition at line 1241 of file utilities.cpp.

1241  {
1242  int i=0;
1243  for (int n=0; n<hamil.numSumTerms; n++) {
1244  hamil.termCoeffs[n] = getRandomReal(-5, 5);
1245  for (int q=0; q<hamil.numQubits; q++)
1246  if (getRandomReal(-1,1) > 0)
1247  hamil.pauliCodes[i++] = PAULI_Z;
1248  else
1249  hamil.pauliCodes[i++] = PAULI_I;
1250  }
1251 }

References getRandomReal(), PauliHamil::numQubits, PauliHamil::numSumTerms, PAULI_I, PAULI_Z, PauliHamil::pauliCodes, and PauliHamil::termCoeffs.

Referenced by TEST_CASE().

◆ setRandomPauliSum() [1/2]

void setRandomPauliSum ( PauliHamil  hamil)

Populates hamil with random coefficients and pauli codes.

Author
Tyson Jones

Definition at line 1237 of file utilities.cpp.

1237  {
1238  setRandomPauliSum(hamil.termCoeffs, hamil.pauliCodes, hamil.numQubits, hamil.numSumTerms);
1239 }

References PauliHamil::numQubits, PauliHamil::numSumTerms, PauliHamil::pauliCodes, setRandomPauliSum(), and PauliHamil::termCoeffs.

◆ setRandomPauliSum() [2/2]

void setRandomPauliSum ( qreal coeffs,
pauliOpType codes,
int  numQubits,
int  numTerms 
)

Populates the coeffs array with random qreals in (-5, 5), and populates codes with random Pauli codes.

Author
Tyson Jones

Definition at line 1229 of file utilities.cpp.

1229  {
1230  int i=0;
1231  for (int n=0; n<numTerms; n++) {
1232  coeffs[n] = getRandomReal(-5, 5);
1233  for (int q=0; q<numQubits; q++)
1234  codes[i++] = (pauliOpType) getRandomInt(0,4);
1235  }
1236 }

References getRandomInt(), and getRandomReal().

Referenced by setRandomPauliSum(), and TEST_CASE().

◆ setSubMatrix()

void setSubMatrix ( QMatrix dest,
QMatrix  sub,
size_t  r,
size_t  c 
)

Modifies dest by overwriting its submatrix (from top-left corner (r, c) to bottom-right corner (r + dest.size(), c + dest.size()) with the complete elements of sub.

This demands that dest.size() >= sub.size() + max(r,c).

Author
Tyson Jones

Definition at line 222 of file utilities.cpp.

222  {
223  DEMAND( sub.size() + r <= dest.size() );
224  DEMAND( sub.size() + c <= dest.size() );
225  for (size_t i=0; i<sub.size(); i++)
226  for (size_t j=0; j<sub.size(); j++)
227  dest[r+i][c+j] = sub[i][j];
228 }

References DEMAND.

Referenced by getFullOperatorMatrix(), and getSwapMatrix().

◆ setUniqueFilename()

void setUniqueFilename ( char *  outFn,
char *  prefix 
)

Modifies outFn to be a filename of format prefix_NUM.txt where NUM is a new unique integer so far.

This is useful for getting unique filenames for independent test cases of functions requiring reading/writing to file, to avoid IO locks (especially common in distributed mode).

Author
Tyson Jones

Definition at line 1358 of file utilities.cpp.

1358  {
1359  sprintf(outFn, "%s_%d.txt", prefix, fn_unique_suffix_id++);
1360 }

References fn_unique_suffix_id.

Referenced by TEST_CASE().

◆ sublists() [1/4]

CatchGen<int*> sublists ( CatchGen< int > &&  gen,
int  numSamps,
const int *  exclude,
int  numExclude 
)

Returns a Catch2 generator of every length-sublen sublist of the elements generated by gen, which exclude all elements in exclude, in increasing lexographic order.

This generates every fixed-length combination of gen's elements the nominated exclusions, and every permutation of each.

There is on need for the elements of exclude to actually appear in those of gen. sublen must less than or equal to the number of elements in gen, after the nominated exclusions.

Note that the sublist must not be modified, else further generation may break (QuEST's internal functions will indeed modify but restore the qubit index lists given to them, which is ok). Assumes list contains no duplicates, otherwise the generated sublists may be duplicated.

This function can be used like

int sublen = 2;
int exclude[2] = {3,4};
int* sublist = GENERATE_COPY( sublists(range(1,6), sublen, exclude, 2) );

to generate {1,2}, {1,5}, {2,1}, {2,5}, {5,1}, {5,2}

Author
Tyson Jones

◆ sublists() [2/4]

CatchGen<int*> sublists ( CatchGen< int > &&  gen,
int  numSamps,
int  excluded 
)

Returns a Catch2 generator of every length-sublen sublist of the elements generated by gen which exclude element excluded, in increasing lexographic order.

This generates every fixed-length combination of gen's elements the nominated exclusions, and every permutation of each.

sublen must less than or equal to the number of elements in gen, after the nominated exclusion. There is no need for excluded to actually appear in the elements of gen.

Note that the sublist must not be modified, else further generation may break (QuEST's internal functions will indeed modify but restore the qubit index lists given to them, which is ok). Assumes list contains no duplicates, otherwise the generated sublists may be duplicated.

This function can be used like

int sublen = 2;
int excluded = 1;
int* sublist = GENERATE_COPY( sublists(range(1,4), sublen, excluded) );

to generate {2,3}, {3,2}.

Author
Tyson Jones

◆ sublists() [3/4]

CatchGen<int*> sublists ( CatchGen< int > &&  gen,
int  sublen 
)

Returns a Catch2 generator of every length-sublen sublist of the elements generated by gen, in increasing lexographic order.

This generates every fixed-length combination of gen's elements, and every permutation of each. Note that the produced sublist must not be modified, else further generation may break (QuEST's internal functions will indeed modify but restore the qubit index lists given to them, which is ok). Assumes list contains no duplicates, otherwise the generated sublists may be duplicated.

This function can be used like

int sublen = 2;
int* sublist = GENERATE_COPY( sublists(list, 4, sublen) );

to generate {1,2}, {1,3}, {1,4}, {2,1}, {2,3}, {2,4}, {3,1}, {3,2}, {3, 4}, {4,1}, {4,2}, {4, 3}.

Author
Tyson Jones

◆ sublists() [4/4]

CatchGen<int*> sublists ( int *  list,
int  len,
int  sublen 
)

Returns a Catch2 generator of every length-sublen sublist of length-len list, in increasing lexographic order.

This generates every fixed-length combination of the given list and every permutation of each. & If the sublist length is the full list length, this generator produces every permutation correctly. Note that the sublist must not be modified, else further & generation may break (QuEST's internal functions will indeed modify but restore the qubit index lists given to them, which is ok). Assumes list contains no duplicates, otherwise the generated sublists may be duplicated.

This function can be used like

int list[4] = {1,2,3,4};
int sublen = 2;
int* sublist = GENERATE_COPY( sublists(list, 4, sublen) );

to generate {1,2}, {1,3}, {1,4}, {2,1}, {2,3}, {2,4}, {3,1}, {3,2}, {3, 4}, {4,1}, {4,2}, {4, 3}.

Author
Tyson Jones

Definition at line 1488 of file utilities.cpp.

1490  {
1491  return Catch::Generators::GeneratorWrapper<int*>(
1492  std::unique_ptr<Catch::Generators::IGenerator<int*>>(
1493  new SubListGenerator(list, len, sublen)));
1494 }

Referenced by TEST_CASE().

◆ toComplexMatrix2()

ComplexMatrix2 toComplexMatrix2 ( QMatrix  qm)

Returns a ComplexMatrix2 copy of QMatix qm.

Demands that qm is a 2-by-2 matrix.

Author
Tyson Jones

Definition at line 1021 of file utilities.cpp.

1021  {
1022  DEMAND( qm.size() == 2 );
1023  ComplexMatrix2 cm;
1024  macro_copyQMatrix(cm, qm);
1025  return cm;
1026 }

References DEMAND, and macro_copyQMatrix.

Referenced by TEST_CASE().

◆ toComplexMatrix4()

ComplexMatrix4 toComplexMatrix4 ( QMatrix  qm)

Returns a ComplexMatrix4 copy of QMatix qm.

Demands that qm is a 4-by-4 matrix.

Author
Tyson Jones

Definition at line 1027 of file utilities.cpp.

1027  {
1028  DEMAND( qm.size() == 4 );
1029  ComplexMatrix4 cm;
1030  macro_copyQMatrix(cm, qm);
1031  return cm;
1032 }

References DEMAND, and macro_copyQMatrix.

Referenced by TEST_CASE().

◆ toComplexMatrixN()

void toComplexMatrixN ( QMatrix  qm,
ComplexMatrixN  cm 
)

Initialises cm with the values of qm.

Demands that cm is a previously created ComplexMatrixN instance, with the same dimensions as qm.

Author
Tyson Jones

Definition at line 1033 of file utilities.cpp.

1033  {
1034  DEMAND( qm.size() == (1u<<cm.numQubits) );
1035  macro_copyQMatrix(cm, qm);
1036 }

References DEMAND, macro_copyQMatrix, and ComplexMatrixN::numQubits.

Referenced by TEST_CASE().

◆ toDiagonalQMatrix()

QMatrix toDiagonalQMatrix ( QVector  vec)

Returns a diagonal complex matrix formed by the given vector.

Author
Tyson Jones

Definition at line 1309 of file utilities.cpp.

1309  {
1310  QMatrix mat = getZeroMatrix(vec.size());
1311  for (size_t i=0; i<vec.size(); i++)
1312  mat[i][i] = vec[i];
1313  return mat;
1314 }

References getZeroMatrix().

◆ toQMatrix() [1/8]

QMatrix toQMatrix ( Complex  alpha,
Complex  beta 
)

Returns the matrix (where a=alpha, b=beta) {{a, -conj(b)}, {b, conj(a)}} using the qcomp complex type.

Author
Tyson Jones

Definition at line 1062 of file utilities.cpp.

1062  {
1063  qcomp a = qcomp(alpha.real, alpha.imag);
1064  qcomp b = qcomp(beta.real, beta.imag);
1065  QMatrix matr{
1066  {a, -conj(b)},
1067  {b, conj(a)}};
1068  return matr;
1069 }

References Complex::imag, qcomp, and Complex::real.

◆ toQMatrix() [2/8]

QMatrix toQMatrix ( ComplexMatrix2  src)

Returns a copy of the given 2-by-2 matrix.

Author
Tyson Jones

Definition at line 1044 of file utilities.cpp.

1044  {
1045  QMatrix dest = getZeroMatrix(2);
1046  macro_copyComplexMatrix(dest, src);
1047  return dest;
1048 }

References getZeroMatrix(), and macro_copyComplexMatrix.

Referenced by TEST_CASE(), and toQMatrix().

◆ toQMatrix() [3/8]

QMatrix toQMatrix ( ComplexMatrix4  src)

Returns a copy of the given 4-by-4 matrix.

Author
Tyson Jones

Definition at line 1049 of file utilities.cpp.

1049  {
1050  QMatrix dest = getZeroMatrix(4);
1051  macro_copyComplexMatrix(dest, src);
1052  return dest;
1053 }

References getZeroMatrix(), and macro_copyComplexMatrix.

◆ toQMatrix() [4/8]

QMatrix toQMatrix ( ComplexMatrixN  src)

Returns a copy of the given 2^N-by-2^N matrix.

Author
Tyson Jones

Definition at line 1054 of file utilities.cpp.

1054  {
1055  DEMAND( src.real != NULL );
1056  DEMAND( src.imag != NULL );
1057  QMatrix dest = getZeroMatrix(1 << src.numQubits);
1058  macro_copyComplexMatrix(dest, src);
1059  return dest;
1060 }

References DEMAND, getZeroMatrix(), ComplexMatrixN::imag, macro_copyComplexMatrix, ComplexMatrixN::numQubits, and ComplexMatrixN::real.

◆ toQMatrix() [5/8]

QMatrix toQMatrix ( DiagonalOp  op)

Returns a 2^N-by-2^N complex diagonal matrix form of the DiagonalOp.

Author
Tyson Jones

Definition at line 1193 of file utilities.cpp.

1193  {
1194  QVector vec = toQVector(op);
1195  QMatrix mat = getZeroMatrix(1LL << op.numQubits);
1196  for (size_t i=0; i<mat.size(); i++)
1197  mat[i][i] = vec[i];
1198  return mat;
1199 }

References getZeroMatrix(), DiagonalOp::numQubits, and toQVector().

◆ toQMatrix() [6/8]

QMatrix toQMatrix ( PauliHamil  hamil)

Returns a 2^N-by-2^N Hermitian matrix form of the PauliHamil.

Author
Tyson Jones

Definition at line 1282 of file utilities.cpp.

1282  {
1283  return toQMatrix(hamil.termCoeffs, hamil.pauliCodes, hamil.numQubits, hamil.numSumTerms);
1284 }

References PauliHamil::numQubits, PauliHamil::numSumTerms, PauliHamil::pauliCodes, PauliHamil::termCoeffs, and toQMatrix().

◆ toQMatrix() [7/8]

QMatrix toQMatrix ( qreal coeffs,
pauliOpType paulis,
int  numQubits,
int  numTerms 
)

Returns a 2^N-by-2^N Hermitian matrix form of the specified weighted sum of Pauli products.

Author
Tyson Jones

Definition at line 1253 of file utilities.cpp.

1253  {
1254 
1255  // produce a numTargs-big matrix 'pauliSum' by pauli-matrix tensoring and summing
1256  QMatrix iMatr{{1,0},{0,1}};
1257  QMatrix xMatr{{0,1},{1,0}};
1258  QMatrix yMatr{{0,-qcomp(0,1)},{qcomp(0,1),0}};
1259  QMatrix zMatr{{1,0},{0,-1}};
1260  QMatrix pauliSum = getZeroMatrix(1<<numQubits);
1261 
1262  for (int t=0; t<numTerms; t++) {
1263  QMatrix pauliProd = QMatrix{{1}};
1264 
1265  for (int q=0; q<numQubits; q++) {
1266  int i = q + t*numQubits;
1267 
1268  QMatrix fac;
1269  pauliOpType code = paulis[i];
1270  if (code == PAULI_I) fac = iMatr;
1271  if (code == PAULI_X) fac = xMatr;
1272  if (code == PAULI_Y) fac = yMatr;
1273  if (code == PAULI_Z) fac = zMatr;
1274  pauliProd = getKroneckerProduct(fac, pauliProd);
1275  }
1276  pauliSum += coeffs[t] * pauliProd;
1277  }
1278 
1279  // a now 2^numQubits by 2^numQubits Hermitian matrix
1280  return pauliSum;
1281 }

References getKroneckerProduct(), getZeroMatrix(), PAULI_I, PAULI_X, PAULI_Y, PAULI_Z, and qcomp.

◆ toQMatrix() [8/8]

QMatrix toQMatrix ( Qureg  qureg)

Returns an equal-size copy of the given density matrix qureg.

In GPU mode, this function involves a copy of qureg from GPU memory to RAM. In distributed mode, this involves an all-to-all broadcast of qureg.

Author
Tyson Jones

Definition at line 1071 of file utilities.cpp.

1071  {
1072  DEMAND( qureg.isDensityMatrix );
1073 #ifdef DISTRIBUTED_MODE
1074  DEMAND( qureg.numAmpsTotal < MPI_MAX_AMPS_IN_MSG );
1075 #endif
1076 
1077  // ensure local qureg.stateVec is up to date
1078  copyStateFromGPU(qureg);
1080 
1081  qreal* fullRe;
1082  qreal* fullIm;
1083 
1084  // in distributed mode, give every node the full state vector
1085 #ifdef DISTRIBUTED_MODE
1086  fullRe = (qreal*) malloc(qureg.numAmpsTotal * sizeof *fullRe);
1087  fullIm = (qreal*) malloc(qureg.numAmpsTotal * sizeof *fullIm);
1088  MPI_Allgather(
1089  qureg.stateVec.real, qureg.numAmpsPerChunk, MPI_QuEST_REAL,
1090  fullRe, qureg.numAmpsPerChunk, MPI_QuEST_REAL, MPI_COMM_WORLD);
1091  MPI_Allgather(
1092  qureg.stateVec.imag, qureg.numAmpsPerChunk, MPI_QuEST_REAL,
1093  fullIm, qureg.numAmpsPerChunk, MPI_QuEST_REAL, MPI_COMM_WORLD);
1094 #else
1095  fullRe = qureg.stateVec.real;
1096  fullIm = qureg.stateVec.imag;
1097 #endif
1098 
1099  // copy full state vector into a QVector
1100  long long int dim = (1 << qureg.numQubitsRepresented);
1101  QMatrix matr = getZeroMatrix(dim);
1102  for (long long int n=0; n<qureg.numAmpsTotal; n++)
1103  matr[n%dim][n/dim] = qcomp(fullRe[n], fullIm[n]);
1104 
1105  // clean up if we malloc'd the distributed array
1106 #ifdef DISTRIBUTED_MODE
1107  free(fullRe);
1108  free(fullIm);
1109 #endif
1110  return matr;
1111 }

References copyStateFromGPU(), DEMAND, getZeroMatrix(), Qureg::isDensityMatrix, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, Qureg::numQubitsRepresented, qcomp, qreal, QUEST_ENV, Qureg::stateVec, and syncQuESTEnv().

◆ toQureg() [1/2]

void toQureg ( Qureg  qureg,
QMatrix  mat 
)

Initialises the density matrix qureg to have the same amplitudes as mat.

Demands qureg is a density matrix of equal dimensions to mat. In GPU mode, this function involves a copy from RAM to GPU memory. This function has no communication cost in distributed mode.

Author
Tyson Jones

Definition at line 1214 of file utilities.cpp.

1214  {
1215  DEMAND( qureg.isDensityMatrix );
1216  DEMAND( (1 << qureg.numQubitsRepresented) == (long long int) mat.size() );
1217 
1219 
1220  int len = (1 << qureg.numQubitsRepresented);
1221  for (int i=0; i<qureg.numAmpsPerChunk; i++) {
1222  int ind = qureg.chunkId*qureg.numAmpsPerChunk + i;
1223  qureg.stateVec.real[i] = real(mat[ind%len][ind/len]);
1224  qureg.stateVec.imag[i] = imag(mat[ind%len][ind/len]);
1225  }
1226  copyStateToGPU(qureg);
1227 }

References Qureg::chunkId, copyStateToGPU(), DEMAND, Qureg::isDensityMatrix, Qureg::numAmpsPerChunk, Qureg::numQubitsRepresented, QUEST_ENV, Qureg::stateVec, and syncQuESTEnv().

◆ toQureg() [2/2]

void toQureg ( Qureg  qureg,
QVector  vec 
)

Initialises the state-vector qureg to have the same amplitudes as vec.

Demands qureg is a state-vector of an equal size to vec. In GPU mode, this function involves a copy from RAM to GPU memory. This function has no communication cost in distributed mode.

Author
Tyson Jones

Definition at line 1201 of file utilities.cpp.

1201  {
1202  DEMAND( !qureg.isDensityMatrix );
1203  DEMAND( qureg.numAmpsTotal == (long long int) vec.size() );
1204 
1206 
1207  for (int i=0; i<qureg.numAmpsPerChunk; i++) {
1208  int ind = qureg.chunkId*qureg.numAmpsPerChunk + i;
1209  qureg.stateVec.real[i] = real(vec[ind]);
1210  qureg.stateVec.imag[i] = imag(vec[ind]);
1211  }
1212  copyStateToGPU(qureg);
1213 }

References Qureg::chunkId, copyStateToGPU(), DEMAND, Qureg::isDensityMatrix, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, QUEST_ENV, Qureg::stateVec, and syncQuESTEnv().

Referenced by TEST_CASE().

◆ toQVector() [1/2]

QVector toQVector ( DiagonalOp  op)

Returns a vector with the same of the full diagonal operator, populated with op's elements.

In distributed mode, this involves an all-to-all broadcast of op.

Author
Tyson Jones

Definition at line 1155 of file utilities.cpp.

1155  {
1156  long long int totalElems = (1LL << op.numQubits);
1157 #ifdef DISTRIBUTED_MODE
1158  DEMAND( totalElems < MPI_MAX_AMPS_IN_MSG );
1159 #endif
1160 
1161  qreal* fullRe;
1162  qreal* fullIm;
1163 
1164  // in distributed mode, give every node the full diagonal operator
1165 #ifdef DISTRIBUTED_MODE
1166  fullRe = (qreal*) malloc(totalElems * sizeof *fullRe);
1167  fullIm = (qreal*) malloc(totalElems * sizeof *fullIm);
1168 
1169  MPI_Allgather(
1170  op.real, op.numElemsPerChunk, MPI_QuEST_REAL,
1171  fullRe, op.numElemsPerChunk, MPI_QuEST_REAL, MPI_COMM_WORLD);
1172  MPI_Allgather(
1173  op.imag, op.numElemsPerChunk, MPI_QuEST_REAL,
1174  fullIm, op.numElemsPerChunk, MPI_QuEST_REAL, MPI_COMM_WORLD);
1175 #else
1176  fullRe = op.real;
1177  fullIm = op.imag;
1178 #endif
1179 
1180  // copy full state vector into a QVector
1181  QVector vec = QVector(totalElems);
1182  for (long long int i=0; i<totalElems; i++)
1183  vec[i] = qcomp(fullRe[i], fullIm[i]);
1184 
1185  // clean up if we malloc'd distrib array
1186 #ifdef DISTRIBUTED_MODE
1187  free(fullRe);
1188  free(fullIm);
1189 #endif
1190  return vec;
1191 }

References DEMAND, DiagonalOp::imag, DiagonalOp::numElemsPerChunk, DiagonalOp::numQubits, qcomp, qreal, and DiagonalOp::real.

◆ toQVector() [2/2]

QVector toQVector ( Qureg  qureg)

Returns an equal-size copy of the given state-vector qureg.

In GPU mode, this function involves a copy of qureg from GPU memory to RAM. In distributed mode, this involves an all-to-all broadcast of qureg.

Author
Tyson Jones

Definition at line 1113 of file utilities.cpp.

1113  {
1114  DEMAND( !qureg.isDensityMatrix );
1115 #ifdef DISTRIBUTED_MODE
1116  DEMAND( qureg.numAmpsTotal < MPI_MAX_AMPS_IN_MSG );
1117 #endif
1118 
1119  // ensure local qureg.stateVec is up to date
1120  copyStateFromGPU(qureg);
1122 
1123  qreal* fullRe;
1124  qreal* fullIm;
1125 
1126  // in distributed mode, give every node the full state vector
1127 #ifdef DISTRIBUTED_MODE
1128  fullRe = (qreal*) malloc(qureg.numAmpsTotal * sizeof *fullRe);
1129  fullIm = (qreal*) malloc(qureg.numAmpsTotal * sizeof *fullIm);
1130 
1131  MPI_Allgather(
1132  qureg.stateVec.real, qureg.numAmpsPerChunk, MPI_QuEST_REAL,
1133  fullRe, qureg.numAmpsPerChunk, MPI_QuEST_REAL, MPI_COMM_WORLD);
1134  MPI_Allgather(
1135  qureg.stateVec.imag, qureg.numAmpsPerChunk, MPI_QuEST_REAL,
1136  fullIm, qureg.numAmpsPerChunk, MPI_QuEST_REAL, MPI_COMM_WORLD);
1137 #else
1138  fullRe = qureg.stateVec.real;
1139  fullIm = qureg.stateVec.imag;
1140 #endif
1141 
1142  // copy full state vector into a QVector
1143  QVector vec = QVector(qureg.numAmpsTotal);
1144  for (long long int i=0; i<qureg.numAmpsTotal; i++)
1145  vec[i] = qcomp(fullRe[i], fullIm[i]);
1146 
1147  // clean up if we malloc'd distrib array
1148 #ifdef DISTRIBUTED_MODE
1149  free(fullRe);
1150  free(fullIm);
1151 #endif
1152  return vec;
1153 }

References copyStateFromGPU(), DEMAND, Qureg::isDensityMatrix, Qureg::numAmpsPerChunk, Qureg::numAmpsTotal, qcomp, qreal, QUEST_ENV, Qureg::stateVec, and syncQuESTEnv().

Referenced by TEST_CASE(), and toQMatrix().

◆ writeToFileSynch()

void writeToFileSynch ( char *  fn,
const string &  contents 
)

Writes contents to the file with filename fn, which is created and/or overwritten.

In distributed mode, the master node writes while the other nodes wait until complete.

Author
Tyson Jones

Definition at line 1362 of file utilities.cpp.

1362  {
1363 
1364  // master node writes
1365  if (QUEST_ENV.rank == 0) {
1366  FILE* file = fopen(fn, "w");
1367  fputs(contents.c_str(), file);
1368  fclose(file);
1369  }
1370 
1371  // other nodes wait
1373 }

References QUEST_ENV, QuESTEnv::rank, and syncQuESTEnv().

Referenced by TEST_CASE().

pauliOpType
Codes for specifying Pauli operators.
Definition: QuEST.h:96
QMatrix getFullOperatorMatrix(int *ctrls, int numCtrls, int *targs, int numTargs, QMatrix op, int numQubits)
Takes a 2^numTargs-by-2^numTargs matrix op and a returns a 2^numQubits-by-2^numQubits matrix where op...
Definition: utilities.cpp:304
long long int getUnsigned(long long int twosComp, int numBits)
Return the unsigned value of a number, made of #numBits bits, which under two's complement,...
Definition: utilities.cpp:1298
void copyStateFromGPU(Qureg qureg)
In GPU mode, this copies the state-vector (or density matrix) from GPU memory (qureg....
Definition: QuEST_cpu.c:45
QuESTEnv QUEST_ENV
The global QuESTEnv instance, to be created and destroyed once in this main(), so that the MPI enviro...
Definition: main.cpp:20
void syncQuESTEnv(QuESTEnv env)
Guarantees that all code up to the given point has been executed on all nodes (if running in distribu...
@ PAULI_Z
Definition: QuEST.h:96
QVector getKroneckerProduct(QVector b, QVector a)
Returns b (otimes) a.
Definition: utilities.cpp:143
int rank
Definition: QuEST.h:364
#define macro_copyComplexMatrix(dest, src)
Copies ComplexMatrix structures into a QMatrix.
Definition: utilities.cpp:1039
@ PAULI_I
Definition: QuEST.h:96
QMatrix getConjugateTranspose(QMatrix a)
Returns the conjugate transpose of the complex square matrix a.
Definition: utilities.cpp:189
int getRandomInt(int min, int max)
Returns a random integer between min (inclusive) and max (exclusive), from the uniform distribution.
Definition: utilities.cpp:526
QMatrix getSwapMatrix(int qb1, int qb2, int numQb)
Returns the 2^numQb-by-2^numQb unitary matrix which swaps qubits qb1 and qb2; the SWAP gate of not-ne...
Definition: utilities.cpp:230
@ TWOS_COMPLEMENT
Definition: QuEST.h:269
QMatrix getRandomUnitary(int numQb)
Returns a uniformly random (under Haar) 2^numQb-by-2^numQb unitary matrix.
Definition: utilities.cpp:530
qreal getRandomReal(qreal min, qreal max)
Returns a random real between min (inclusive) and max (exclusive), from the uniform distribution.
Definition: utilities.cpp:421
QMatrix getKetBra(QVector ket, QVector bra)
Returns the matrix |ket><bra|, with ith-jth element ket(i) conj(bra(j)), since |ket><bra| = sum_i a_i...
Definition: utilities.cpp:169
Represents a 4x4 matrix of complex numbers.
Definition: QuEST.h:175
@ UNSIGNED
Definition: QuEST.h:269
#define qreal
QMatrix toQMatrix(ComplexMatrix2 src)
Returns a copy of the given 2-by-2 matrix.
Definition: utilities.cpp:1044
unsigned int calcLog2(long unsigned int res)
Returns log2 of numbers which must be gauranteed to be 2^n.
Definition: utilities.cpp:372
@ PAULI_X
Definition: QuEST.h:96
static int fn_unique_suffix_id
Definition: utilities.cpp:1356
void setRandomPauliSum(qreal *coeffs, pauliOpType *codes, int numQubits, int numTerms)
Populates the coeffs array with random qreals in (-5, 5), and populates codes with random Pauli codes...
Definition: utilities.cpp:1229
int chunkId
The position of the chunk of the state vector held by this process in the full state vector.
Definition: QuEST.h:336
qcomp expI(qreal phase)
Returns the unit-norm complex number exp(i*phase).
Definition: utilities.cpp:417
std::vector< qcomp > QVector
A complex vector, which can be zero-initialised with QVector(numAmps).
Definition: utilities.hpp:60
qreal * imag
The imaginary values of the 2^numQubits complex elements.
Definition: QuEST.h:310
void setSubMatrix(QMatrix &dest, QMatrix sub, size_t r, size_t c)
Modifies dest by overwriting its submatrix (from top-left corner (r, c) to bottom-right corner (r + d...
Definition: utilities.cpp:222
QVector toQVector(Qureg qureg)
Returns an equal-size copy of the given state-vector qureg.
Definition: utilities.cpp:1113
long long int numAmpsPerChunk
Number of probability amplitudes held in stateVec by this process In the non-MPI version,...
Definition: QuEST.h:332
qreal * termCoeffs
The real coefficient of each Pauli product. This is an array of length PauliHamil....
Definition: QuEST.h:283
#define qcomp
enum pauliOpType * pauliCodes
The Pauli operators acting on each qubit, flattened over every operator.
Definition: QuEST.h:281
void copyStateToGPU(Qureg qureg)
In GPU mode, this copies the state-vector (or density matrix) from RAM (qureg.stateVec) to VRAM / GPU...
Definition: QuEST_cpu.c:42
int numQubits
The number of qubits this operator can act on (informing its size)
Definition: QuEST.h:300
int numSumTerms
The number of terms in the weighted sum, or the number of Pauli products.
Definition: QuEST.h:285
@ PAULI_Y
Definition: QuEST.h:96
QVector getRandomStateVector(int numQb)
Returns a random numQb-length L2-normalised state-vector from an undisclosed distribution.
Definition: utilities.cpp:468
qreal ** real
Definition: QuEST.h:189
QMatrix getRandomQMatrix(int dim)
Returns a dim-by-dim complex matrix, where the real and imaginary value of each element are independe...
Definition: utilities.cpp:379
long long int getIndexOfTargetValues(long long int ref, int *targs, int numTargs, int targVal)
Definition: utilities.cpp:689
qreal ** imag
Definition: QuEST.h:190
ComplexArray stateVec
Computational state amplitudes - a subset thereof in the MPI version.
Definition: QuEST.h:341
long long int numElemsPerChunk
The number of the 2^numQubits amplitudes stored on each distributed node.
Definition: QuEST.h:302
int isDensityMatrix
Whether this instance is a density-state representation.
Definition: QuEST.h:325
int numQubits
Definition: QuEST.h:188
void updateIndices(int oldEl, int newEl, int *list1, int len1, int *list2, int len2)
Definition: utilities.cpp:287
std::vector< std::vector< qcomp > > QMatrix
A complex square matrix.
Definition: utilities.hpp:49
int numQubits
The number of qubits informing the Hilbert dimension of the Hamiltonian.
Definition: QuEST.h:287
QVector getNormalised(QVector vec)
Returns an L2-normalised copy of vec, using Kahan summation for improved accuracy.
Definition: utilities.cpp:446
int numQubitsRepresented
The number of qubits represented in either the state-vector or density matrix.
Definition: QuEST.h:327
std::vector< qreal > getRandomProbabilities(int numProbs)
Returns a list of random real scalars, each in [0, 1], which sum to unity.
Definition: utilities.cpp:472
long long int numAmpsTotal
Total number of amplitudes, which are possibly distributed among machines.
Definition: QuEST.h:334
qreal * real
The real values of the 2^numQubits complex elements.
Definition: QuEST.h:308
qreal real
Definition: QuEST.h:105
QMatrix getIdentityMatrix(size_t dim)
Returns a dim-by-dim identity matrix.
Definition: utilities.cpp:161
QMatrix getPureDensityMatrix(QVector state)
Returns a density matrix initialised into the given pure state.
Definition: utilities.cpp:507
qreal imag
Definition: QuEST.h:106
QMatrix getZeroMatrix(size_t dim)
Returns a dim-by-dim square complex matrix, initialised to all zeroes.
Definition: utilities.cpp:153
QVector getRandomQVector(int dim)
Returns a dim-length vector with random complex amplitudes in the square joining {-1-i,...
Definition: utilities.cpp:435
bool areEqual(QVector a, QVector b)
Returns true if the absolute value of the difference between every amplitude in vectors a and b is le...
Definition: utilities.cpp:398
void applyReferenceOp(QVector &state, int *ctrls, int numCtrls, int *targs, int numTargs, QMatrix op)
Modifies the state-vector state to be the result of applying the multi-target operator matrix op,...
Definition: utilities.cpp:728
#define M_PI
Definition: QuEST_common.c:41
qcomp getRandomComplex()
Returns a random complex number within the square closing (-1-i) and (1+i), from a distribution unifo...
Definition: utilities.cpp:431
#define DEMAND(cond)
Definition: utilities.cpp:24
#define macro_copyQMatrix(dest, src)
Definition: utilities.cpp:1013
Represents a 2x2 matrix of complex numbers.
Definition: QuEST.h:137
long long int getValueOfTargets(long long int ind, int *targs, int numTargs)
Returns the integer value of the targeted sub-register for the given full state index ind.
Definition: utilities.cpp:670