Debugging

Utilities for seeding and debugging, such as state-logging. More...

Functions

void copyStateFromGPU (Qureg qureg)
 In GPU mode, this copies the state-vector (or density matrix) from GPU memory (qureg.deviceStateVec) to RAM (qureg.stateVec), where it can be accessed/modified by the user. More...
 
void copyStateToGPU (Qureg qureg)
 In GPU mode, this copies the state-vector (or density matrix) from RAM (qureg.stateVec) to VRAM / GPU-memory (qureg.deviceStateVec), which is the version operated upon by other calls to the API. More...
 
void getEnvironmentString (QuESTEnv env, Qureg qureg, char str[200])
 Sets str to a string containing the number of qubits in qureg, and the hardware facilities used (e.g. More...
 
void initDebugState (Qureg qureg)
 Initialises qureg to be in the un-normalised, non-physical state with with n-th complex amplitude (2n/10 + i(2n+1)/10). More...
 
void invalidQuESTInputError (const char *errMsg, const char *errFunc)
 An internal function called when invalid arguments are passed to a QuEST API call, which the user can optionally override by redefining. More...
 
void reportPauliHamil (PauliHamil hamil)
 Print the PauliHamil to screen. More...
 
void reportQuESTEnv (QuESTEnv env)
 Report information about the QuEST environment. More...
 
void reportQuregParams (Qureg qureg)
 Report metainformation about a set of qubits: number of qubits, number of probability amplitudes. More...
 
void reportState (Qureg qureg)
 Print the current state vector of probability amplitudes for a set of qubits to file. More...
 
void reportStateToScreen (Qureg qureg, QuESTEnv env, int reportRank)
 Print the current state vector of probability amplitudes for a set of qubits to standard out. More...
 
void seedQuEST (unsigned long int *seedArray, int numSeeds)
 Seed the Mersenne Twister used for random number generation in the QuEST environment with a user defined seed. More...
 
void seedQuESTDefault (void)
 Seed the Mersenne Twister used for random number generation in the QuEST environment with an example defualt seed. More...
 
void syncQuESTEnv (QuESTEnv env)
 Guarantees that all code up to the given point has been executed on all nodes (if running in distributed mode) More...
 
int syncQuESTSuccess (int successCode)
 Performs a logical AND on all successCodes held by all processes. More...
 

Detailed Description

Utilities for seeding and debugging, such as state-logging.

Author
Ania Brown
Tyson Jones
Balint Koczor
Nicolas Vogt of HQS (one-qubit damping)

Function Documentation

◆ copyStateFromGPU()

void copyStateFromGPU ( Qureg  qureg)

In GPU mode, this copies the state-vector (or density matrix) from GPU memory (qureg.deviceStateVec) to RAM (qureg.stateVec), where it can be accessed/modified by the user.

In CPU mode, this function has no effect. In conjunction with copyStateToGPU(), this allows a user to directly modify the state-vector in a harware agnostic way. Note though that users should instead use setAmps() if possible.

For example, to set the first real element to 1, one could do:

copyStateFromGPU(qureg);
qureg.stateVec.real[0] = 1;
copyStateToGPU(qureg);

Note users should never access qureg.deviceStateVec directly.

Parameters
[in,out]quregthe qureg of which to copy .deviceStateVec to .stateVec in GPU mode
Author
Ania Brown
Tyson Jones (doc)

Definition at line 39 of file QuEST_cpu.c.

39  {
40 }

References DEBUG, Qureg::deviceStateVec, Qureg::numAmpsPerChunk, and Qureg::stateVec.

Referenced by areEqual(), densmatr_calcTotalProb(), statevec_calcTotalProb(), statevec_compareStates(), statevec_reportStateToScreen(), TEST_CASE(), toQMatrix(), and toQVector().

◆ copyStateToGPU()

void copyStateToGPU ( Qureg  qureg)

In GPU mode, this copies the state-vector (or density matrix) from RAM (qureg.stateVec) to VRAM / GPU-memory (qureg.deviceStateVec), which is the version operated upon by other calls to the API.

In CPU mode, this function has no effect. In conjunction with copyStateFromGPU() (which should be called first), this allows a user to directly modify the state-vector in a harware agnostic way. Note though that users should instead use setAmps() if possible.

For example, to set the first real element to 1, one could do:

copyStateFromGPU(qureg);
qureg.stateVec.real[0] = 1;
copyStateToGPU(qureg);

Note users should never access qureg.deviceStateVec directly.

Parameters
[in,out]quregthe qureg of which to copy .stateVec to .deviceStateVec in GPU mode
Author
Ania Brown
Tyson Jones (doc)

Definition at line 36 of file QuEST_cpu.c.

36  {
37 }

References DEBUG, Qureg::deviceStateVec, Qureg::numAmpsPerChunk, and Qureg::stateVec.

Referenced by statevec_initStateFromSingleFile(), TEST_CASE(), and toQureg().

◆ getEnvironmentString()

void getEnvironmentString ( QuESTEnv  env,
Qureg  qureg,
char  str[200] 
)

Sets str to a string containing the number of qubits in qureg, and the hardware facilities used (e.g.

GPU, MPI and/or OMP).

Parameters
[in]envobject representing the execution environment. A single instance is used for each program
[in]quregthe qureg of which to query the simulating hardware
[out]strto be populated with the output string
Author
Ania Brown

Definition at line 447 of file QuEST_gpu.cu.

447  {
448  sprintf(str, "%dqubits_GPU_noMpi_noOMP", qureg.numQubitsInStateVec);
449 }

References Qureg::numQubitsInStateVec.

◆ initDebugState()

void initDebugState ( Qureg  qureg)

Initialises qureg to be in the un-normalised, non-physical state with with n-th complex amplitude (2n/10 + i(2n+1)/10).

This is used internally for debugging and testing.

Parameters
[in,out]quregthe register to have its amplitudes overwritten
Author
Ania Brown
Tyson Jones (doc)

Definition at line 1308 of file QuEST.c.

1308  {
1309  statevec_initDebugState(qureg);
1310 }

References statevec_initDebugState().

Referenced by TEST_CASE().

◆ invalidQuESTInputError()

void invalidQuESTInputError ( const char *  errMsg,
const char *  errFunc 
)

An internal function called when invalid arguments are passed to a QuEST API call, which the user can optionally override by redefining.

This function is a weak symbol, so that users can choose how input errors are handled, by redefining it in their own code. Users must ensure that the triggered API call does not continue (e.g. the user exits or throws an exception), else QuEST will continue with the valid input and likely trigger a seg-fault. This function is triggered before any internal state-change, hence it is safe to interrupt with exceptions.

E.g. in C

void invalidQuESTInputError(const char* errMsg, const char* errFunc) {
// log to file
printf("ERROR! Writing to file...\n");
FILE *fp = fopen("errorlog.txt", "w");
fprintf(fp, "incorrect usage of function '%s': %s", errFunc, errMsg);
fclose(fp);
// exit
exit(1);
}

This function is compatible with C++ exceptions, though note the user of extern "C":

extern "C" void invalidQuESTInputError(const char* errMsg, const char* errFunc) {
string err = "in function " + string(errFunc) + ": " + string(errMsg);
throw std::invalid_argument(err);
}
Parameters
[in]errMsga string describing the nature of the argument error
[in]errFuncthe name of the invalidly-called API function
Exceptions
invalidQuESTInputErrorunless overriden by the user
Author
Tyson Jones

An internal function called when invalid arguments are passed to a QuEST API call, which the user can optionally override by redefining.

an negative qubit index). This is redefined here to, in lieu of printing and exiting, throw a C++ exception which can be caught (and hence unit tested for) by Catch2

Definition at line 176 of file QuEST_validation.c.

176  {
177  exitWithError(errMsg, errFunc);
178 }

References exitWithError().

Referenced by QuESTAssert(), validateFileOpened(), validateHamilFileCoeffParsed(), validateHamilFileParams(), validateHamilFilePauliCode(), and validateHamilFilePauliParsed().

◆ reportPauliHamil()

void reportPauliHamil ( PauliHamil  hamil)

Print the PauliHamil to screen.

The output features a new line for each term, each with format

c p1 p2 p3 ... pN

where c is the real coefficient of the term, and p1 ... pN are numbers 0, 1, 2, 3 to indicate identity, pauliX, pauliY and pauliZ operators respectively, acting on qubits 0 through N-1 (all qubits). A tab character separates c and p1, but single spaces separate the Pauli operators.

Parameters
[in]hamilan instantiated PauliHamil
Exceptions
invalidQuESTInputErrorif the parameters of hamil are invalid, i.e. if numQubits <= 0, or if numSumTerms <= 0, or if pauliCodes contains an invalid Pauli code.
Author
Tyson Jones

Definition at line 1328 of file QuEST.c.

1328  {
1329  validatePauliHamil(hamil, __func__);
1330 
1331  for (int t=0; t<hamil.numSumTerms; t++) {
1332  printf("%g\t", hamil.termCoeffs[t]);
1333  for (int q=0; q<hamil.numQubits; q++)
1334  printf("%d ", (int) hamil.pauliCodes[q+t*hamil.numQubits]);
1335  printf("\n");
1336  }
1337 }

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

◆ reportQuESTEnv()

void reportQuESTEnv ( QuESTEnv  env)

Report information about the QuEST environment.

Parameters
[in]envobject representing the execution environment. A single instance is used for each program
Author
Ania Brown

Definition at line 179 of file QuEST_cpu_distributed.c.

179  {
180  if (env.rank==0){
181  printf("EXECUTION ENVIRONMENT:\n");
182  printf("Running distributed (MPI) version\n");
183  printf("Number of ranks is %d\n", env.numRanks);
184 # ifdef _OPENMP
185  printf("OpenMP enabled\n");
186  printf("Number of threads available is %d\n", omp_get_max_threads());
187 # else
188  printf("OpenMP disabled\n");
189 # endif
190  printf("Precision: size of qreal is %ld bytes\n", sizeof(qreal) );
191  }
192 }

References QuESTEnv::numRanks, qreal, and QuESTEnv::rank.

◆ reportQuregParams()

void reportQuregParams ( Qureg  qureg)

Report metainformation about a set of qubits: number of qubits, number of probability amplitudes.

Parameters
[in]quregobject representing the set of qubits
Author
Ania Brown

Definition at line 234 of file QuEST_common.c.

234  {
235  long long int numAmps = 1LL << qureg.numQubitsInStateVec;
236  long long int numAmpsPerRank = numAmps/qureg.numChunks;
237  if (qureg.chunkId==0){
238  printf("QUBITS:\n");
239  printf("Number of qubits is %d.\n", qureg.numQubitsInStateVec);
240  printf("Number of amps is %lld.\n", numAmps);
241  printf("Number of amps per rank is %lld.\n", numAmpsPerRank);
242  }
243 }

References Qureg::chunkId, Qureg::numChunks, and Qureg::numQubitsInStateVec.

◆ reportState()

void reportState ( Qureg  qureg)

Print the current state vector of probability amplitudes for a set of qubits to file.

File format:

real, imag
realComponent1, imagComponent1
realComponent2, imagComponent2
...
realComponentN, imagComponentN

File naming convention:

For each node that the program runs on, a file 'state_rank_[node_rank].csv' is generated. If there is more than one node, ranks after the first do not include the header

real, imag

so that files are easier to combine.

Parameters
[in,out]quregobject representing the set of qubits
Author
Ania Brown

Definition at line 216 of file QuEST_common.c.

216  {
217  FILE *state;
218  char filename[100];
219  long long int index;
220  sprintf(filename, "state_rank_%d.csv", qureg.chunkId);
221  state = fopen(filename, "w");
222  if (qureg.chunkId==0) fprintf(state, "real, imag\n");
223 
224  for(index=0; index<qureg.numAmpsPerChunk; index++){
225  # if QuEST_PREC==1 || QuEST_PREC==2
226  fprintf(state, "%.12f, %.12f\n", qureg.stateVec.real[index], qureg.stateVec.imag[index]);
227  # elif QuEST_PREC == 4
228  fprintf(state, "%.12Lf, %.12Lf\n", qureg.stateVec.real[index], qureg.stateVec.imag[index]);
229  #endif
230  }
231  fclose(state);
232 }

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

◆ reportStateToScreen()

void reportStateToScreen ( Qureg  qureg,
QuESTEnv  env,
int  reportRank 
)

Print the current state vector of probability amplitudes for a set of qubits to standard out.

For debugging purposes. Each rank should print output serially. Only print output for systems <= 5 qubits

Author
Ania Brown

Definition at line 1324 of file QuEST.c.

1324  {
1325  statevec_reportStateToScreen(qureg, env, reportRank);
1326 }

References statevec_reportStateToScreen().

◆ seedQuEST()

void seedQuEST ( unsigned long int *  seedArray,
int  numSeeds 
)

Seed the Mersenne Twister used for random number generation in the QuEST environment with a user defined seed.

This function uses the mt19937 init_by_array function with numSeeds keys supplied by the user. Subsequent calls to mt19937 genrand functions will use this seeding. For a multi process code, the same seed is given to all process, therefore this seeding is only appropriate to use for functions such as measure where all processes require the same random value.

For more information about the MT, see http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/emt19937ar.html

Parameters
[in]seedArrayArray of integers to use as seed. This allows the MT to be initialised with more than a 32-bit integer if required
[in]numSeedsLength of seedArray
Author
Ania Brown

Seed the Mersenne Twister used for random number generation in the QuEST environment with a user defined seed.

Definition at line 209 of file QuEST_common.c.

209  {
210  // init MT random number generator with user defined list of seeds
211  // for the MPI version, it is ok that all procs will get the same seed as random numbers will only be
212  // used by the master process
213  init_by_array(seedArray, numSeeds);
214 }

References init_by_array().

◆ seedQuESTDefault()

void seedQuESTDefault ( void  )

Seed the Mersenne Twister used for random number generation in the QuEST environment with an example defualt seed.

This default seeding function uses the mt19937 init_by_array function with two keys – time and pid. Subsequent calls to mt19937 genrand functions will use this seeding. For a multi process code, the same seed is given to all process, therefore this seeding is only appropriate to use for functions such as measure where all processes require the same random value.

For more information about the MT, see http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MT2002/emt19937ar.html

Author
Ania Brown
Balint Koczor (Windows compatibility)

Definition at line 1318 of file QuEST_cpu_distributed.c.

1318  {
1319  // init MT random number generator with three keys -- time and pid
1320  // for the MPI version, it is ok that all procs will get the same seed as random numbers will only be
1321  // used by the master process
1322 
1323  unsigned long int key[2];
1325  // this seed will be used to generate the same random number on all procs,
1326  // therefore we want to make sure all procs receive the same key
1327  MPI_Bcast(key, 2, MPI_UNSIGNED_LONG, 0, MPI_COMM_WORLD);
1328  init_by_array(key, 2);
1329 }

References getQuESTDefaultSeedKey(), and init_by_array().

Referenced by createQuESTEnv().

◆ syncQuESTEnv()

void syncQuESTEnv ( QuESTEnv  env)

Guarantees that all code up to the given point has been executed on all nodes (if running in distributed mode)

Parameters
[in]envobject representing the execution environment. A single instance is used for each program
Author
Ania Brown

Definition at line 162 of file QuEST_cpu_distributed.c.

162  {
163  MPI_Barrier(MPI_COMM_WORLD);
164 }

Referenced by areEqual(), statevec_initStateFromSingleFile(), statevec_reportStateToScreen(), toQMatrix(), toQureg(), and toQVector().

◆ syncQuESTSuccess()

int syncQuESTSuccess ( int  successCode)

Performs a logical AND on all successCodes held by all processes.

If any one process has a zero successCode all processes will return a zero success code.

Parameters
[in]successCode1 if process task succeeded, 0 if process task failed
Returns
1 if all processes succeeded, 0 if any one process failed
Author
Ania Brown

Definition at line 166 of file QuEST_cpu_distributed.c.

166  {
167  int totalSuccess;
168  MPI_Allreduce(&successCode, &totalSuccess, 1, MPI_INT, MPI_LAND, MPI_COMM_WORLD);
169  return totalSuccess;
170 }
void init_by_array(unsigned long init_key[], int key_length)
Definition: mt19937ar.c:80
int rank
Definition: QuEST.h:244
int numChunks
Number of chunks the state vector is broken up into – the number of MPI processes used.
Definition: QuEST.h:219
void getQuESTDefaultSeedKey(unsigned long int *key)
Definition: QuEST_common.c:182
#define qreal
int numQubitsInStateVec
Number of qubits in the state-vector - this is double the number represented for mixed states.
Definition: QuEST.h:210
void invalidQuESTInputError(const char *errMsg, const char *errFunc)
An internal function called when invalid arguments are passed to a QuEST API call,...
int chunkId
The position of the chunk of the state vector held by this process in the full state vector.
Definition: QuEST.h:217
void exitWithError(const char *msg, const char *func)
long long int numAmpsPerChunk
Number of probability amplitudes held in stateVec by this process In the non-MPI version,...
Definition: QuEST.h:213
qreal * termCoeffs
The coefficient of each Pauli product. This is a length numSumTerms array.
Definition: QuEST.h:164
enum pauliOpType * pauliCodes
The Pauli operators acting on each qubit, flattened over every operator.
Definition: QuEST.h:162
void statevec_reportStateToScreen(Qureg qureg, QuESTEnv env, int reportRank)
Print the current state vector of probability amplitudes for a set of qubits to standard out.
Definition: QuEST_cpu.c:1366
int numRanks
Definition: QuEST.h:245
int numSumTerms
The number of terms in the weighted sum, or the number of Pauli products.
Definition: QuEST.h:166
void statevec_initDebugState(Qureg qureg)
Initialise the state vector of probability amplitudes to an (unphysical) state with each component of...
Definition: QuEST_cpu.c:1591
void validatePauliHamil(PauliHamil hamil, const char *caller)
ComplexArray stateVec
Computational state amplitudes - a subset thereof in the MPI version.
Definition: QuEST.h:222
int numQubits
The number of qubits for which this Hamiltonian is defined.
Definition: QuEST.h:168