test_calculations.cpp
Go to the documentation of this file.
1 
2 #include "catch.hpp"
3 #include "QuEST.h"
4 #include "utilities.hpp"
5 
6 /* allows concise use of Contains in catch's REQUIRE_THROWS_WITH */
7 using Catch::Matchers::Contains;
8 
9 
10 
15 TEST_CASE( "calcDensityInnerProduct", "[calculations]" ) {
16 
17  Qureg mat1 = createDensityQureg(NUM_QUBITS, QUEST_ENV);
18  Qureg mat2 = createDensityQureg(NUM_QUBITS, QUEST_ENV);
19 
20  SECTION( "correctness" ) {
21 
22  // repeat these random tests 10 times
23  GENERATE( range(0,10) );
24 
25  SECTION( "density-matrix" ) {
26 
27  SECTION( "pure" ) {
28 
29  // mat1 = |r1><r1|
30  QVector r1 = getRandomStateVector(NUM_QUBITS);
31  toQureg(mat1, getKetBra(r1,r1));
32 
33  // mat2 = |r2><r2|
34  QVector r2 = getRandomStateVector(NUM_QUBITS);
35  toQureg(mat2, getKetBra(r2,r2));
36 
37  // prod( |r1><r1|, |r2><r2| ) = |<r1|r2>|^2
38  qcomp prod = 0;
39  for (size_t i=0; i<r1.size(); i++)
40  prod += conj(r1[i]) * r2[i];
41  qreal densProd = pow(abs(prod),2);
42 
43  REQUIRE( calcDensityInnerProduct(mat1,mat2) == Approx(densProd) );
44  }
45  SECTION( "mixed" ) {
46 
47  QMatrix ref1 = getRandomDensityMatrix(NUM_QUBITS);
48  QMatrix ref2 = getRandomDensityMatrix(NUM_QUBITS);
49  toQureg(mat1, ref1);
50  toQureg(mat2, ref2);
51 
52  // prod(mat1, mat2) = sum_{ij} conj(mat1_{ij}) * mat2_{ij}
53  qcomp refProd = 0;
54  for (size_t i=0; i<ref1.size(); i++)
55  for (size_t j=0; j<ref1.size(); j++)
56  refProd += conj(ref1[i][j]) * ref2[i][j];
57  REQUIRE( imag(refProd) == Approx(0).margin(REAL_EPS) );
58 
59  REQUIRE( calcDensityInnerProduct(mat1,mat2) == Approx(real(refProd)) );
60 
61  // should be invariant under ordering
62  REQUIRE( calcDensityInnerProduct(mat1,mat2) == Approx(calcDensityInnerProduct(mat2,mat1)) );
63  }
64  SECTION( "unnormalised" ) {
65 
66  // set both to random (non-Hermitian) complex matrices
67  QMatrix ref1 = getRandomQMatrix(1<<NUM_QUBITS);
68  QMatrix ref2 = getRandomQMatrix(1<<NUM_QUBITS);
69  toQureg(mat1, ref1);
70  toQureg(mat2, ref2);
71 
72  // prod(mat1, mat2) = real(sum_{ij} conj(mat1_{ij}) * mat2_{ij})
73  qcomp refProd = 0;
74  for (size_t i=0; i<ref1.size(); i++)
75  for (size_t j=0; j<ref1.size(); j++)
76  refProd += conj(ref1[i][j]) * ref2[i][j];
77 
78  REQUIRE( calcDensityInnerProduct(mat1,mat2) == Approx(real(refProd)) );
79  }
80  }
81  }
82  SECTION( "input validation" ) {
83 
84  SECTION( "dimensions" ) {
85 
86  Qureg mat3 = createDensityQureg(NUM_QUBITS + 1, QUEST_ENV);
87  REQUIRE_THROWS_WITH( calcDensityInnerProduct(mat1,mat3), Contains("Dimensions") && Contains("don't match") );
88  destroyQureg(mat3, QUEST_ENV);
89  }
90  SECTION( "state-vectors" ) {
91 
92  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
93 
94  REQUIRE_THROWS_WITH( calcDensityInnerProduct(mat1,vec), Contains("valid only for density matrices") );
95  REQUIRE_THROWS_WITH( calcDensityInnerProduct(vec,mat1), Contains("valid only for density matrices") );
96  REQUIRE_THROWS_WITH( calcDensityInnerProduct(vec,vec), Contains("valid only for density matrices") );
97 
98  destroyQureg(vec, QUEST_ENV);
99  }
100  }
101  destroyQureg(mat1, QUEST_ENV);
102  destroyQureg(mat2, QUEST_ENV);
103 }
104 
105 
106 
111 TEST_CASE( "calcExpecPauliProd", "[calculations]" ) {
112 
113  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
114  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
115  initDebugState(vec);
116  initDebugState(mat);
117  QVector vecRef = toQVector(vec);
118  QMatrix matRef = toQMatrix(mat);
119 
120  Qureg vecWork = createQureg(NUM_QUBITS, QUEST_ENV);
121  Qureg matWork = createDensityQureg(NUM_QUBITS, QUEST_ENV);
122 
123  SECTION( "correctness" ) {
124 
125  int numTargs = GENERATE( range(1,NUM_QUBITS+1) );
126  int* targs = GENERATE_COPY( sublists(range(0,NUM_QUBITS), numTargs) );
127 
128  /* it's too expensive to try ALL Pauli sequences, via
129  * pauliOpType* paulis = GENERATE_COPY( pauliseqs(numTargs) );.
130  * Furthermore, take(10, pauliseqs(numTargs)) will try the same pauli codes.
131  * Hence, we instead opt to repeatedlyrandomly generate pauliseqs
132  */
133  GENERATE( range(0,10) ); // gen 10 random pauli-codes for every targs
134  pauliOpType paulis[numTargs];
135  for (int i=0; i<numTargs; i++)
136  paulis[i] = (pauliOpType) getRandomInt(0,4);
137 
138  // produce a numTargs-big matrix 'pauliProd' by pauli-matrix tensoring
139  QMatrix iMatr{{1,0},{0,1}};
140  QMatrix xMatr{{0,1},{1,0}};
141  QMatrix yMatr{{0,-1i},{1i,0}};
142  QMatrix zMatr{{1,0},{0,-1}};
143  QMatrix pauliProd{{1}};
144  for (int i=0; i<numTargs; i++) {
145  QMatrix fac;
146  if (paulis[i] == PAULI_I) fac = iMatr;
147  if (paulis[i] == PAULI_X) fac = xMatr;
148  if (paulis[i] == PAULI_Y) fac = yMatr;
149  if (paulis[i] == PAULI_Z) fac = zMatr;
150  pauliProd = getKroneckerProduct(fac, pauliProd);
151  }
152 
153  SECTION( "state-vector" ) {
154 
155  /* calcExpecPauliProd calculates <qureg|pauliProd|qureg> */
156 
157  QVector prodRef = vecRef;
158  applyReferenceOp(prodRef, targs, numTargs, pauliProd);
159  qcomp prod = 0;
160  for (size_t i=0; i<vecRef.size(); i++)
161  prod += conj(vecRef[i]) * prodRef[i];
162  REQUIRE( imag(prod) == Approx(0).margin(REAL_EPS) );
163 
164  qreal res = calcExpecPauliProd(vec, targs, paulis, numTargs, vecWork);
165  REQUIRE( res == Approx(real(prod)).margin(REAL_EPS) );
166  }
167  SECTION( "density matrix" ) {
168 
169  /* calcExpecPauliProd calculates Trace( pauliProd * qureg ) */
170 
171  // produce (pauliProd * mat)
172  QMatrix fullOp = getFullOperatorMatrix(NULL, 0, targs, numTargs, pauliProd, NUM_QUBITS);
173  matRef = fullOp * matRef;
174 
175  // compute real(trace(pauliProd * mat))
176  qreal tr = 0;
177  for (size_t i=0; i<matRef.size(); i++)
178  tr += real(matRef[i][i]);
179  // (get real, since we start in a non-Hermitian state, hence diagonal isn't real)
180 
181  qreal res = calcExpecPauliProd(mat, targs, paulis, numTargs, matWork);
182  REQUIRE( res == Approx(tr).margin(10*REAL_EPS) );
183  }
184  }
185  SECTION( "validation" ) {
186 
187  SECTION( "number of targets" ) {
188 
189  int numTargs = GENERATE( -1, 0, NUM_QUBITS+1 );
190  REQUIRE_THROWS_WITH( calcExpecPauliProd(vec, NULL, NULL, numTargs, vecWork), Contains("Invalid number of target") );
191  }
192  SECTION( "target indices" ) {
193 
194  int numTargs = 3;
195  int targs[3] = {0, 1, 2};
196 
197  // make one index wrong
198  targs[GENERATE( range(0,3) )] = GENERATE( -1, NUM_QUBITS );
199  REQUIRE_THROWS_WITH( calcExpecPauliProd(vec, targs, NULL, numTargs, vecWork), Contains("Invalid target qubit") );
200  }
201  SECTION( "repetition in targets" ) {
202 
203  int numTargs = 3;
204  int targs[3] = {0, 1, 1};
205  REQUIRE_THROWS_WITH( calcExpecPauliProd(vec, targs, NULL, numTargs, vecWork), Contains("target qubits must be unique") );
206  }
207  SECTION( "pauli codes" ) {
208 
209  int numTargs = 3;
210  int targs[3] = {0, 1, 2};
211  pauliOpType codes[3] = {PAULI_X, PAULI_Y, PAULI_Z};
212 
213  // make one pauli wrong
214  codes[GENERATE( range(0,3) )] = (pauliOpType) GENERATE( -1, 4 );
215  REQUIRE_THROWS_WITH( calcExpecPauliProd(vec, targs, codes, numTargs, vecWork), Contains("Invalid Pauli code") );
216  }
217  SECTION( "workspace type" ) {
218 
219  int numTargs = 1;
220  int targs[1] = {0};
221  pauliOpType codes[1] = {PAULI_I};
222 
223  REQUIRE_THROWS_WITH( calcExpecPauliProd(vec, targs, codes, numTargs, matWork), Contains("Registers must both be state-vectors or both be density matrices") );
224  REQUIRE_THROWS_WITH( calcExpecPauliProd(mat, targs, codes, numTargs, vecWork), Contains("Registers must both be state-vectors or both be density matrices") );
225  }
226  SECTION( "workspace dimensions" ) {
227 
228  int numTargs = 1;
229  int targs[1] = {0};
230  pauliOpType codes[1] = {PAULI_I};
231 
232  Qureg vec2 = createQureg(NUM_QUBITS + 1, QUEST_ENV);
233  REQUIRE_THROWS_WITH( calcExpecPauliProd(vec, targs, codes, numTargs, vec2), Contains("Dimensions") && Contains("don't match") );
234  destroyQureg(vec2, QUEST_ENV);
235 
236  Qureg mat2 = createDensityQureg(NUM_QUBITS + 1, QUEST_ENV);
237  REQUIRE_THROWS_WITH( calcExpecPauliProd(mat, targs, codes, numTargs, mat2), Contains("Dimensions") && Contains("don't match") );
238  destroyQureg(mat2, QUEST_ENV);
239  }
240  }
241  destroyQureg(vec, QUEST_ENV);
242  destroyQureg(mat, QUEST_ENV);
243  destroyQureg(vecWork, QUEST_ENV);
244  destroyQureg(matWork, QUEST_ENV);
245 }
246 
247 
248 
253 TEST_CASE( "calcExpecPauliSum", "[calculations]" ) {
254 
255  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
256  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
257  initDebugState(vec);
258  initDebugState(mat);
259  QVector vecRef = toQVector(vec);
260  QMatrix matRef = toQMatrix(mat);
261 
262  Qureg vecWork = createQureg(NUM_QUBITS, QUEST_ENV);
263  Qureg matWork = createDensityQureg(NUM_QUBITS, QUEST_ENV);
264 
265  SECTION( "correctness" ) {
266 
267  int numSumTerms = GENERATE( 1, 2, 10, 15 );
268 
269  /* it's too expensive to try every possible Pauli configuration, so
270  * we'll try 10 random codes, and for each, random coefficients
271  */
272  GENERATE( range(0,10) );
273  int totNumCodes = numSumTerms*NUM_QUBITS;
274  pauliOpType paulis[totNumCodes];
275  for (int i=0; i<totNumCodes; i++)
276  paulis[i] = (pauliOpType) getRandomInt(0,4);
277 
278  // for every above param configuration, try random coefficients
279  qreal coeffs[numSumTerms];
280  for (int i=0; i<numSumTerms; i++)
281  coeffs[i] = getRandomReal(-5, 5);
282 
283  // produce a numTargs-big matrix 'pauliSum' by pauli-matrix tensoring and summing
284  QMatrix iMatr{{1,0},{0,1}};
285  QMatrix xMatr{{0,1},{1,0}};
286  QMatrix yMatr{{0,-1i},{1i,0}};
287  QMatrix zMatr{{1,0},{0,-1}};
288  QMatrix pauliSum = getZeroMatrix(1<<NUM_QUBITS);
289 
290  for (int t=0; t<numSumTerms; t++) {
291  QMatrix pauliProd = QMatrix{{1}};
292 
293  for (int q=0; q<NUM_QUBITS; q++) {
294  int i = q + t*NUM_QUBITS;
295 
296  QMatrix fac;
297  if (paulis[i] == PAULI_I) fac = iMatr;
298  if (paulis[i] == PAULI_X) fac = xMatr;
299  if (paulis[i] == PAULI_Y) fac = yMatr;
300  if (paulis[i] == PAULI_Z) fac = zMatr;
301  pauliProd = getKroneckerProduct(fac, pauliProd);
302  }
303  pauliSum += coeffs[t] * pauliProd;
304  }
305  SECTION( "state-vector" ) {
306 
307  /* calcExpecPauliSum calculates <qureg|pauliSum|qureg> */
308 
309  QVector sumRef = pauliSum * vecRef;
310  qcomp prod = 0;
311  for (size_t i=0; i<vecRef.size(); i++)
312  prod += conj(vecRef[i]) * sumRef[i];
313  REQUIRE( imag(prod) == Approx(0).margin(10*REAL_EPS) );
314 
315  qreal res = calcExpecPauliSum(vec, paulis, coeffs, numSumTerms, vecWork);
316  REQUIRE( res == Approx(real(prod)).margin(10*REAL_EPS) );
317  }
318  SECTION( "density matrix" ) {
319 
320  /* calcExpecPauliSum calculates Trace( pauliSum * qureg ) */
321  matRef = pauliSum * matRef;
322  qreal tr = 0;
323  for (size_t i=0; i<matRef.size(); i++)
324  tr += real(matRef[i][i]);
325  // (get real, since we start in a non-Hermitian state, hence diagonal isn't real)
326 
327  qreal res = calcExpecPauliSum(mat, paulis, coeffs, numSumTerms, matWork);
328  REQUIRE( res == Approx(tr).margin(1E2*REAL_EPS) );
329  }
330  }
331  SECTION( "validation" ) {
332 
333  SECTION( "number of sum terms" ) {
334 
335  int numSumTerms = GENERATE( -1, 0 );
336  REQUIRE_THROWS_WITH( calcExpecPauliSum(vec, NULL, NULL, numSumTerms, vecWork), Contains("Invalid number of terms in the Pauli sum") );
337  }
338  SECTION( "pauli codes" ) {
339 
340  // make valid params
341  int numSumTerms = 3;
342  qreal coeffs[numSumTerms];
343  pauliOpType codes[numSumTerms*NUM_QUBITS];
344  for (int i=0; i<numSumTerms*NUM_QUBITS; i++)
345  codes[i] = PAULI_I;
346 
347  // make one pauli wrong
348  codes[GENERATE_COPY( range(0,numSumTerms*NUM_QUBITS) )] = (pauliOpType) GENERATE( -1, 4 );
349  REQUIRE_THROWS_WITH( calcExpecPauliSum(vec, codes, coeffs, numSumTerms, vecWork), Contains("Invalid Pauli code") );
350  }
351  SECTION( "workspace type" ) {
352 
353  // make valid params
354  int numSumTerms = 1;
355  qreal coeffs[1] = {0};
356  pauliOpType codes[NUM_QUBITS];
357  for (int i=0; i<NUM_QUBITS; i++)
358  codes[i] = PAULI_I;
359 
360  REQUIRE_THROWS_WITH( calcExpecPauliSum(vec, codes, coeffs, numSumTerms, mat), Contains("Registers must both be state-vectors or both be density matrices") );
361  REQUIRE_THROWS_WITH( calcExpecPauliSum(mat, codes, coeffs, numSumTerms, vec), Contains("Registers must both be state-vectors or both be density matrices") );
362  }
363  SECTION( "workspace dimensions" ) {
364 
365  // make valid params
366  int numSumTerms = 1;
367  qreal coeffs[1] = {0};
368  pauliOpType codes[NUM_QUBITS];
369  for (int i=0; i<NUM_QUBITS; i++)
370  codes[i] = PAULI_I;
371 
372  Qureg vec2 = createQureg(NUM_QUBITS + 1, QUEST_ENV);
373  REQUIRE_THROWS_WITH( calcExpecPauliSum(vec, codes, coeffs, numSumTerms, vec2), Contains("Dimensions") && Contains("don't match") );
374  destroyQureg(vec2, QUEST_ENV);
375 
376  Qureg mat2 = createDensityQureg(NUM_QUBITS + 1, QUEST_ENV);
377  REQUIRE_THROWS_WITH( calcExpecPauliSum(mat, codes, coeffs, numSumTerms, mat2), Contains("Dimensions") && Contains("don't match") );
378  destroyQureg(mat2, QUEST_ENV);
379  }
380  }
381  destroyQureg(vec, QUEST_ENV);
382  destroyQureg(mat, QUEST_ENV);
383  destroyQureg(vecWork, QUEST_ENV);
384  destroyQureg(matWork, QUEST_ENV);
385 }
386 
387 
388 
393 TEST_CASE( "calcFidelity", "[calculations]" ) {
394 
395  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
396  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
397  Qureg pure = createQureg(NUM_QUBITS, QUEST_ENV);
398 
399  SECTION( "correctness" ) {
400 
401  // repeat the below random tests 10 times
402  GENERATE( range(0,10) );
403 
404  SECTION( "state-vector" ) {
405 
406  /* calcFidelity computes |<vec|pure>|^2 */
407 
408  SECTION( "normalised" ) {
409 
410  // random L2 vectors
411  QVector vecRef = getRandomStateVector(NUM_QUBITS);
412  QVector pureRef = getRandomStateVector(NUM_QUBITS);
413  toQureg(vec, vecRef);
414  toQureg(pure, pureRef);
415 
416  // |<vec|vec>|^2 = |1|^2 = 1
417  REQUIRE( calcFidelity(vec,vec) == Approx(1) );
418 
419  // |<vec|pure>|^2 = |sum_j conj(vec_j) * pure_j|^2
420  qcomp dotProd = 0;
421  for (size_t i=0; i<vecRef.size(); i++)
422  dotProd += conj(vecRef[i]) * pureRef[i];
423  qreal refFid = pow(abs(dotProd), 2);
424 
425  REQUIRE( calcFidelity(vec,pure) == Approx(refFid) );
426  }
427  SECTION( "unnormalised" ) {
428 
429  // random unnormalised vectors
430  QVector vecRef = getRandomQVector(1<<NUM_QUBITS);
431  QVector pureRef = getRandomQVector(1<<NUM_QUBITS);
432  toQureg(vec, vecRef);
433  toQureg(pure, pureRef);
434 
435  // Let nv be magnitude of vec, hence |unit-vec> = 1/sqrt(nv)|vec>
436  qreal nv = 0;
437  for (size_t i=0; i<vecRef.size(); i++)
438  nv += pow(abs(vecRef[i]), 2);
439  // then <vec|vec> = sqrt(nv)*sqrt(nv) <unit-vec|unit-vec> = nv,
440  // hence |<vec|vec>|^2 = nv*nv
441  REQUIRE( calcFidelity(vec,vec) == Approx( nv*nv ) );
442 
443  qcomp dotProd = 0;
444  for (size_t i=0; i<vecRef.size(); i++)
445  dotProd += conj(vecRef[i]) * pureRef[i];
446  qreal refFid = pow(abs(dotProd), 2);
447 
448  REQUIRE( calcFidelity(vec,pure) == Approx(refFid) );
449  }
450  }
451  SECTION( "density-matrix" ) {
452 
453  /* calcFidelity computes <pure|mat|pure> */
454 
455  SECTION( "pure" ) {
456 
457  QVector pureRef = getRandomStateVector(NUM_QUBITS);
458  toQureg(pure, pureRef);
459 
460  // test when density matrix is the same pure state
461  QMatrix matRef = getKetBra(pureRef, pureRef);
462  toQureg(mat, matRef);
463  REQUIRE( calcFidelity(mat,pure) == Approx(1) );
464 
465  // test when density matrix is a random pure state
466  QVector r1 = getRandomStateVector(NUM_QUBITS);
467  matRef = getKetBra(r1, r1); // actually pure |r1><r1|
468  toQureg(mat, matRef);
469 
470  // <pure|r1><r1|pure> = |<r1|pure>|^2 = |sum_j conj(r1_j) * pure_j|^2
471  qcomp dotProd = 0;
472  for (size_t i=0; i<r1.size(); i++)
473  dotProd += conj(r1[i]) * pureRef[i];
474  qreal refFid = pow(abs(dotProd), 2);
475 
476  REQUIRE( calcFidelity(mat,pure) == Approx(refFid) );
477  }
478  SECTION( "mixed" ) {
479 
480  QVector pureRef = getRandomStateVector(NUM_QUBITS);
481  toQureg(pure, pureRef);
482 
483  // test when density matrix is mixed
484  QMatrix matRef = getRandomDensityMatrix(NUM_QUBITS);
485  toQureg(mat, matRef);
486 
487  // <pure|mat|pure> = <pure| (Mat|pure>)
488  QVector rhs = matRef * pureRef;
489  qcomp dotProd = 0;
490  for (size_t i=0; i<rhs.size(); i++)
491  dotProd += conj(pureRef[i]) * rhs[i];
492 
493  REQUIRE( imag(dotProd) == Approx(0).margin(REAL_EPS) );
494  REQUIRE( calcFidelity(mat,pure) == Approx(real(dotProd)) );
495  }
496  SECTION( "unnormalised" ) {
497 
498  // test when both density matrix and pure state are unnormalised
499  QVector pureRef = getRandomQVector(1<<NUM_QUBITS);
500  QMatrix matRef = getRandomQMatrix(1<<NUM_QUBITS);
501  toQureg(pure, pureRef);
502  toQureg(mat, matRef);
503 
504  // real[ <pure|mat|pure> ] = real[ <pure| (Mat|pure>) ]
505  QVector rhs = matRef * pureRef;
506  qcomp dotProd = 0;
507  for (size_t i=0; i<rhs.size(); i++)
508  dotProd += conj(pureRef[i]) * rhs[i];
509 
510  REQUIRE( calcFidelity(mat,pure) == Approx(real(dotProd)) );
511  }
512  }
513  }
514  SECTION( "input validation" ) {
515 
516  SECTION( "dimensions" ) {
517 
518  // two state-vectors
519  Qureg vec2 = createQureg(vec.numQubitsRepresented + 1, QUEST_ENV);
520  REQUIRE_THROWS_WITH( calcFidelity(vec2,vec), Contains("Dimensions") && Contains("don't match") );
521  destroyQureg(vec2, QUEST_ENV);
522 
523  // density-matrix and state-vector
524  Qureg mat2 = createDensityQureg(vec.numQubitsRepresented + 1, QUEST_ENV);
525  REQUIRE_THROWS_WITH( calcFidelity(mat2,vec), Contains("Dimensions") && Contains("don't match") );
526  destroyQureg(mat2, QUEST_ENV);
527  }
528  SECTION( "density-matrices" ) {
529 
530  // two mixed statess
531  REQUIRE_THROWS_WITH( calcFidelity(mat,mat), Contains("Second argument must be a state-vector") );
532  }
533  }
534  destroyQureg(vec, QUEST_ENV);
535  destroyQureg(mat, QUEST_ENV);
536  destroyQureg(pure, QUEST_ENV);
537 }
538 
539 
540 
545 TEST_CASE( "calcHilbertSchmidtDistance", "[calculations]" ) {
546 
547  Qureg mat1 = createDensityQureg(NUM_QUBITS, QUEST_ENV);
548  Qureg mat2 = createDensityQureg(NUM_QUBITS, QUEST_ENV);
549 
550  SECTION( "correctness" ) {
551 
552  // perform these random tests 10 times
553  GENERATE( range(0,10) );
554 
555  SECTION( "density-matrix" ) {
556 
557  SECTION( "pure" ) {
558 
559  // create random |r1><r1| and |r2><r2| states
560  QVector r1 = getRandomStateVector(NUM_QUBITS);
561  QMatrix m1 = getKetBra(r1,r1);
562  toQureg(mat1, m1);
563  QVector r2 = getRandomStateVector(NUM_QUBITS);
564  QMatrix m2 = getKetBra(r2,r2);
565  toQureg(mat2, m2);
566 
567  // Tr{ (a-b)(a-b)^dagger } = sum_{ij} |a_{ij} - b_{ij}|^2
568  qreal tr = 0;
569  for (size_t i=0; i<m1.size(); i++)
570  for (size_t j=0; j<m1.size(); j++)
571  tr += pow(abs(m1[i][j] - m2[i][j]), 2);
572 
573  qreal res = calcHilbertSchmidtDistance(mat1, mat2);
574  REQUIRE( res == Approx(sqrt(tr)) );
575 
576  }
577  SECTION( "normalised" ) {
578 
579  QMatrix ref1 = getRandomDensityMatrix(NUM_QUBITS);
580  QMatrix ref2 = getRandomDensityMatrix(NUM_QUBITS);
581  toQureg(mat1, ref1);
582  toQureg(mat2, ref2);
583 
584  // Tr{ (a-b)(a-b)^dagger } = sum_{ij} |a_{ij} - b_{ij}|^2
585  qreal tr = 0;
586  for (size_t i=0; i<ref1.size(); i++)
587  for (size_t j=0; j<ref1.size(); j++)
588  tr += pow(abs(ref1[i][j] - ref2[i][j]), 2);
589 
590  qreal res = calcHilbertSchmidtDistance(mat1, mat2);
591  REQUIRE( res == Approx(sqrt(tr)) );
592  }
593  SECTION( "unnormalised" ) {
594 
595  // mat1 and mat2 are both random matrices
596  QMatrix ref1 = getRandomQMatrix(1<<NUM_QUBITS);
597  QMatrix ref2 = getRandomQMatrix(1<<NUM_QUBITS);
598  toQureg(mat1, ref1);
599  toQureg(mat2, ref2);
600 
601  // Tr{ (a-b)(a-b)^dagger } = sum_{ij} |a_{ij} - b_{ij}|^2
602  qreal tr = 0;
603  for (size_t i=0; i<ref1.size(); i++)
604  for (size_t j=0; j<ref1.size(); j++)
605  tr += pow(abs(ref1[i][j] - ref2[i][j]), 2);
606 
607  qreal res = calcHilbertSchmidtDistance(mat1, mat2);
608  REQUIRE( res == Approx(sqrt(tr)) );
609  }
610  }
611  }
612  SECTION( "input validation") {
613 
614  SECTION( "dimensions" ) {
615 
616  Qureg mat3 = createDensityQureg(NUM_QUBITS + 1, QUEST_ENV);
617  REQUIRE_THROWS_WITH( calcHilbertSchmidtDistance(mat1,mat3), Contains("Dimensions") && Contains("don't match") );
618  destroyQureg(mat3, QUEST_ENV);
619  }
620  SECTION( "state-vector" ) {
621 
622  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
623 
624  REQUIRE_THROWS_WITH( calcHilbertSchmidtDistance(vec,mat1), Contains("valid only for density matrices") );
625  REQUIRE_THROWS_WITH( calcHilbertSchmidtDistance(mat1,vec), Contains("valid only for density matrices") );
626  REQUIRE_THROWS_WITH( calcHilbertSchmidtDistance(vec,vec), Contains("valid only for density matrices") );
627 
628  destroyQureg(vec, QUEST_ENV);
629  }
630  }
631  destroyQureg(mat1, QUEST_ENV);
632  destroyQureg(mat2, QUEST_ENV);
633 }
634 
635 
636 
641 TEST_CASE( "calcInnerProduct", "[calculations]" ) {
642 
643  Qureg vec1 = createQureg(NUM_QUBITS, QUEST_ENV);
644  Qureg vec2 = createQureg(NUM_QUBITS, QUEST_ENV);
645 
646  SECTION( "correctness" ) {
647 
648  // perform these random tests 10 times
649  GENERATE( range(0,10) );
650 
651  SECTION( "state-vector" ) {
652 
653  SECTION( "normalised" ) {
654 
655  // <r1|r2> = sum_j conj(r1_j) * r2_j
656  QVector r1 = getRandomStateVector(NUM_QUBITS);
657  QVector r2 = getRandomStateVector(NUM_QUBITS);
658  qcomp prod = 0;
659  for (size_t i=0; i<r1.size(); i++)
660  prod += conj(r1[i]) * r2[i];
661 
662  toQureg(vec1, r1);
663  toQureg(vec2, r2);
664  Complex res = calcInnerProduct(vec1,vec2);
665 
666  REQUIRE( res.real == Approx(real(prod)) );
667  REQUIRE( res.imag == Approx(imag(prod)) );
668  }
669  SECTION( "unnormalised" ) {
670 
671  // <r1|r2> = sum_j conj(r1_j) * r2_j
672  QVector r1 = getRandomQVector(1<<NUM_QUBITS);
673  QVector r2 = getRandomQVector(1<<NUM_QUBITS);
674  qcomp prod = 0;
675  for (size_t i=0; i<r1.size(); i++)
676  prod += conj(r1[i]) * r2[i];
677 
678  toQureg(vec1, r1);
679  toQureg(vec2, r2);
680  Complex res = calcInnerProduct(vec1,vec2);
681 
682  REQUIRE( res.real == Approx(real(prod)) );
683  REQUIRE( res.imag == Approx(imag(prod)) );
684  }
685  }
686  }
687  SECTION( "input validation" ) {
688 
689  SECTION( "dimensions" ) {
690 
691  Qureg vec3 = createQureg(NUM_QUBITS + 1, QUEST_ENV);
692  REQUIRE_THROWS_WITH( calcInnerProduct(vec1,vec3), Contains("Dimensions") && Contains("don't match") );
693  destroyQureg(vec3, QUEST_ENV);
694  }
695  SECTION( "density-matrix" ) {
696 
697  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
698 
699  REQUIRE_THROWS_WITH( calcInnerProduct(vec1,mat), Contains("valid only for state-vectors") );
700  REQUIRE_THROWS_WITH( calcInnerProduct(mat,vec1), Contains("valid only for state-vectors") );
701  REQUIRE_THROWS_WITH( calcInnerProduct(mat,mat), Contains("valid only for state-vectors") );
702 
703  destroyQureg(mat, QUEST_ENV);
704  }
705  }
706  destroyQureg(vec1, QUEST_ENV);
707  destroyQureg(vec2, QUEST_ENV);
708 }
709 
710 
711 
716 TEST_CASE( "calcProbOfOutcome", "[calculations]" ) {
717 
718  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
719  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
720 
721  SECTION( "correctness" ) {
722 
723  int target = GENERATE( range(0,NUM_QUBITS) );
724  int outcome = GENERATE( 0, 1 );
725 
726  SECTION( "state-vector" ) {
727 
728  SECTION( "normalised" ) {
729 
730  QVector ref = getRandomStateVector(NUM_QUBITS);
731  toQureg(vec, ref);
732 
733  // prob is sum of |amp|^2 of amplitudes where target bit is outcome
734  qreal prob = 0;
735  for (size_t ind=0; ind<ref.size(); ind++) {
736  int bit = (ind >> target) & 1; // target-th bit
737  if (bit == outcome)
738  prob += pow(abs(ref[ind]), 2);
739  }
740 
741  REQUIRE( calcProbOfOutcome(vec, target, outcome) == Approx(prob) );
742  }
743  SECTION( "unnormalised" ) {
744 
745  QVector ref = getRandomQVector(1<<NUM_QUBITS);
746  toQureg(vec, ref);
747 
748  // prob is (sum of |amp|^2 of amplitudes where target bit is zero)
749  // or 1 - (this) if outcome == 1
750  qreal prob = 0;
751  for (size_t ind=0; ind<ref.size(); ind++) {
752  int bit = (ind >> target) & 1; // target-th bit
753  if (bit == 0)
754  prob += pow(abs(ref[ind]), 2);
755  }
756  if (outcome == 1)
757  prob = 1 - prob;
758 
759  REQUIRE( calcProbOfOutcome(vec, target, outcome) == Approx(prob) );
760  }
761  }
762  SECTION( "density-matrix" ) {
763 
764  SECTION( "pure" ) {
765 
766  // set mat to a random |r><r|
767  QVector ref = getRandomStateVector(NUM_QUBITS);
768  toQureg(mat, getKetBra(ref, ref));
769 
770  // calc prob of the state-vector
771  qreal prob = 0;
772  for (size_t ind=0; ind<ref.size(); ind++) {
773  int bit = (ind >> target) & 1; // target-th bit
774  if (bit == outcome)
775  prob += pow(abs(ref[ind]), 2);
776  }
777 
778  REQUIRE( calcProbOfOutcome(mat, target, outcome) == Approx(prob) );
779  }
780  SECTION( "mixed" ) {
781 
782  QMatrix ref = getRandomDensityMatrix(NUM_QUBITS);
783  toQureg(mat, ref);
784 
785  // prob is sum of diagonal amps (should be real) where target bit is outcome
786  qcomp tr = 0;
787  for (size_t ind=0; ind<ref.size(); ind++) {
788  int bit = (ind >> target) & 1; // target-th bit
789  if (bit == outcome)
790  tr += ref[ind][ind];
791  }
792  REQUIRE( imag(tr) == Approx(0).margin(REAL_EPS) );
793 
794  REQUIRE( calcProbOfOutcome(mat, target, outcome) == Approx(real(tr)) );
795  }
796  SECTION( "unnormalised" ) {
797 
798  QMatrix ref = getRandomQMatrix(1<<NUM_QUBITS);
799  toQureg(mat, ref);
800 
801  // prob is (sum of real of diagonal amps where target bit is outcome)
802  // or 1 - (this) if outcome == 1
803  qreal tr = 0;
804  for (size_t ind=0; ind<ref.size(); ind++) {
805  int bit = (ind >> target) & 1; // target-th bit
806  if (bit == 0)
807  tr += real(ref[ind][ind]);
808  }
809  if (outcome == 1)
810  tr = 1 - tr;
811 
812  REQUIRE( calcProbOfOutcome(mat, target, outcome) == Approx(tr) );
813  }
814  }
815  }
816  SECTION( "validation" ) {
817 
818  SECTION( "qubit indices" ) {
819 
820  int target = GENERATE( -1, NUM_QUBITS );
821  REQUIRE_THROWS_WITH( calcProbOfOutcome(vec, target, 0), Contains("Invalid target qubit") );
822  }
823  SECTION( "outcome value" ) {
824 
825  int outcome = GENERATE( -1, 2 );
826  REQUIRE_THROWS_WITH( calcProbOfOutcome(vec, 0, outcome), Contains("Invalid measurement outcome") );
827  }
828  }
829  destroyQureg(vec, QUEST_ENV);
830  destroyQureg(mat, QUEST_ENV);
831 }
832 
833 
834 
839 TEST_CASE( "calcPurity", "[calculations]" ) {
840 
841  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
842 
843  SECTION( "correctness" ) {
844 
845  // perform the following random tests 10 times
846  GENERATE( range(1,10) );
847 
848  SECTION( "density-matrix" ) {
849 
850  SECTION( "pure" ) {
851 
852  // pure states have unity purity
853  initZeroState(mat);
854  REQUIRE( calcPurity(mat) == 1 );
855 
856  // (try also a pure random L2-vector)
857  QVector r1 = getRandomStateVector(NUM_QUBITS); // |r>
858  QMatrix m1 = getKetBra(r1, r1); // |r><r|
859  toQureg(mat, m1);
860  REQUIRE( calcPurity(mat) == Approx(1) );
861 
862  }
863  SECTION( "mixed" ) {
864 
865  // mixed states have 1/2^N < purity < 1
866  QMatrix ref = getRandomDensityMatrix(NUM_QUBITS);
867  toQureg(mat, ref);
868  qreal purity = calcPurity(mat);
869  REQUIRE( purity < 1 );
870  REQUIRE( purity >= 1/pow(2.,NUM_QUBITS) );
871 
872  // compare to Tr(rho^2)
873  QMatrix prod = ref*ref;
874  qreal tr = 0;
875  for (size_t i=0; i<prod.size(); i++)
876  tr += real(prod[i][i]);
877  REQUIRE( purity == Approx(tr) );
878  }
879  SECTION( "unnormalised" ) {
880 
881  // unphysical states give sum_{ij} |rho_ij|^2
882  QMatrix ref = getRandomQMatrix(1<<NUM_QUBITS);
883  qreal tot = 0;
884  for (size_t i=0; i<ref.size(); i++)
885  for (size_t j=0; j<ref.size(); j++)
886  tot += pow(abs(ref[i][j]), 2);
887 
888  toQureg(mat, ref);
889  REQUIRE( calcPurity(mat) == Approx(tot) );
890  }
891  }
892  }
893  SECTION( "input validation" ) {
894 
895  SECTION( "state-vector" ) {
896 
897  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
898  REQUIRE_THROWS_WITH( calcPurity(vec), Contains("valid only for density matrices") );
899  destroyQureg(vec, QUEST_ENV);
900  }
901  }
902  destroyQureg(mat, QUEST_ENV);
903 }
904 
905 
906 
911 TEST_CASE( "calcTotalProb", "[calculations]" ) {
912 
913  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
914  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
915 
916  SECTION( "correctness" ) {
917 
918  SECTION( "state-vector" ) {
919 
920  // normalised: prob(vec) = 1
921  initPlusState(vec);
922  REQUIRE( calcTotalProb(vec) == Approx(1) );
923 
924  // zero norm: prob(vec) = 0
925  initBlankState(vec);
926  REQUIRE( calcTotalProb(vec) == 0 );
927 
928  // random L2 state: prob(vec) = 1
929  toQureg(vec, getRandomStateVector(NUM_QUBITS));
930  REQUIRE( calcTotalProb(vec) == Approx(1) );
931 
932  // unnormalised: prob(vec) = sum_i |vec_i|^2
933  initDebugState(vec);
934  QVector ref = toQVector(vec);
935  qreal refProb = 0;
936  for (size_t i=0; i<ref.size(); i++)
937  refProb += pow(abs(ref[i]), 2);
938  REQUIRE( calcTotalProb(vec) == Approx(refProb) );
939  }
940  SECTION( "density-matrix" ) {
941 
942  // normalised: prob(mat) = 1
943  initPlusState(mat);
944  REQUIRE( calcTotalProb(mat) == Approx(1) );
945 
946  // zero norm: prob(mat) = 0
947  initBlankState(mat);
948  REQUIRE( calcTotalProb(mat) == 0 );
949 
950  // random density matrix: prob(mat) = 1
951  toQureg(mat, getRandomDensityMatrix(NUM_QUBITS));
952  REQUIRE( calcTotalProb(mat) == Approx(1) );
953 
954  // unnormalised: prob(mat) = sum_i real(mat_{ii})
955  initDebugState(mat);
956  QMatrix ref = toQMatrix(mat);
957  qreal refProb = 0;
958  for (size_t i=0; i<ref.size(); i++)
959  refProb += real(ref[i][i]);
960  REQUIRE( calcTotalProb(mat) == Approx(refProb) );
961  }
962  }
963  SECTION( "input validation" ) {
964 
965  // no validation
966  SUCCEED();
967  }
968  destroyQureg(vec, QUEST_ENV);
969  destroyQureg(mat, QUEST_ENV);
970 }
971 
972 
973 
978 TEST_CASE( "getAmp", "[calculations]" ) {
979 
980  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
981 
982  SECTION( "correctness" ) {
983 
984  SECTION( "state-vector" ) {
985 
986  initDebugState(vec);
987  QVector ref = toQVector(vec);
988 
989  int ind = GENERATE( range(0,1<<NUM_QUBITS) );
990  Complex amp = getAmp(vec,ind);
991  REQUIRE( fromComplex(amp) == ref[ind] );
992  }
993  }
994  SECTION( "input validation" ) {
995 
996  SECTION( "state index" ) {
997 
998  int ind = GENERATE( -1, 1<<NUM_QUBITS );
999  REQUIRE_THROWS_WITH( getAmp(vec,ind), Contains("Invalid amplitude index") );
1000  }
1001  SECTION( "density-matrix" ) {
1002 
1003  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
1004  REQUIRE_THROWS_WITH( getAmp(mat,0), Contains("valid only for state-vectors") );
1005  destroyQureg(mat, QUEST_ENV);
1006  }
1007  }
1008  destroyQureg(vec, QUEST_ENV);
1009 }
1010 
1011 
1012 
1017 TEST_CASE( "getDensityAmp", "[calculations]" ) {
1018 
1019  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
1020 
1021  SECTION( "correctness" ) {
1022 
1023  SECTION( "density-matrix" ) {
1024 
1025  initDebugState(mat);
1026  QMatrix ref = toQMatrix(mat);
1027 
1028  int row = GENERATE( range(0,1<<NUM_QUBITS) );
1029  int col = GENERATE( range(0,1<<NUM_QUBITS) );
1030 
1031  Complex amp = getDensityAmp(mat,row,col);
1032  REQUIRE( fromComplex(amp) == ref[row][col] );
1033  }
1034  }
1035  SECTION( "input validation" ) {
1036 
1037  SECTION( "state index" ) {
1038 
1039  int ind = GENERATE( -1, 1<<NUM_QUBITS );
1040  REQUIRE_THROWS_WITH( getDensityAmp(mat,ind,0), Contains("Invalid amplitude index") );
1041  REQUIRE_THROWS_WITH( getDensityAmp(mat,0,ind), Contains("Invalid amplitude index") );
1042 
1043  }
1044  SECTION( "state-vector" ) {
1045 
1046  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
1047  REQUIRE_THROWS_WITH( getDensityAmp(vec,0,0), Contains("valid only for density matrices") );
1048  destroyQureg(vec, QUEST_ENV);
1049  }
1050  }
1051  destroyQureg(mat, QUEST_ENV);
1052 }
1053 
1054 
1055 
1060 TEST_CASE( "getImagAmp", "[calculations]" ) {
1061 
1062  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
1063 
1064  SECTION( "correctness" ) {
1065 
1066  SECTION( "state-vector" ) {
1067 
1068  initDebugState(vec);
1069  QVector ref = toQVector(vec);
1070 
1071  int ind = GENERATE( range(0,1<<NUM_QUBITS) );
1072  REQUIRE( getImagAmp(vec,ind) == imag(ref[ind]) );
1073  }
1074  }
1075  SECTION( "input validation" ) {
1076 
1077  SECTION( "state index" ) {
1078 
1079  int ind = GENERATE( -1, 1<<NUM_QUBITS );
1080  REQUIRE_THROWS_WITH( getImagAmp(vec,ind), Contains("Invalid amplitude index") );
1081  }
1082  SECTION( "density-matrix" ) {
1083 
1084  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
1085  REQUIRE_THROWS_WITH( getImagAmp(mat,0), Contains("valid only for state-vectors") );
1086  destroyQureg(mat, QUEST_ENV);
1087  }
1088  }
1089  destroyQureg(vec, QUEST_ENV);
1090 }
1091 
1092 
1093 
1098 TEST_CASE( "getNumAmps", "[calculations]" ) {
1099 
1100  SECTION( "correctness" ) {
1101 
1102  // test >= NUM_QUBITS so as not to limit distribution size
1103  int numQb = GENERATE( range(NUM_QUBITS, NUM_QUBITS+10) );
1104 
1105  SECTION( "state-vector" ) {
1106 
1107  Qureg vec = createQureg(numQb, QUEST_ENV);
1108  REQUIRE( getNumAmps(vec) == (1<<numQb) );
1109  destroyQureg(vec, QUEST_ENV);
1110  }
1111  }
1112  SECTION( "input validation" ) {
1113 
1114  SECTION( "density-matrix" ) {
1115  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
1116  REQUIRE_THROWS_WITH( getNumAmps(mat), Contains("valid only for state-vectors") );
1117  destroyQureg(mat, QUEST_ENV);
1118  }
1119  }
1120 }
1121 
1122 
1123 
1128 TEST_CASE( "getNumQubits", "[calculations]" ) {
1129 
1130  SECTION( "correctness" ) {
1131 
1132  // test >= NUM_QUBITS so as not to limit distribution size
1133  int numQb = GENERATE( range(NUM_QUBITS, NUM_QUBITS+10) );
1134 
1135  SECTION( "state-vector" ) {
1136 
1137  Qureg vec = createQureg(numQb, QUEST_ENV);
1138  REQUIRE( getNumQubits(vec) == numQb );
1139  destroyQureg(vec, QUEST_ENV);
1140  }
1141  SECTION( "density-matrix" ) {
1142 
1143  Qureg mat = createDensityQureg(numQb, QUEST_ENV);
1144  REQUIRE( getNumQubits(mat) == numQb );
1145  destroyQureg(mat, QUEST_ENV);
1146  }
1147  }
1148  SECTION( "input validation" ) {
1149 
1150  // no validation
1151  SUCCEED();
1152  }
1153 }
1154 
1155 
1156 
1161 TEST_CASE( "getProbAmp", "[calculations]" ) {
1162 
1163  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
1164 
1165  SECTION( "correctness" ) {
1166 
1167  SECTION( "state-vector" ) {
1168 
1169  initDebugState(vec);
1170  QVector ref = toQVector(vec);
1171 
1172  int ind = GENERATE( range(0,1<<NUM_QUBITS) );
1173  qreal refCalc = pow(abs(ref[ind]), 2);
1174  REQUIRE( getProbAmp(vec,ind) == Approx(refCalc) );
1175  }
1176  }
1177  SECTION( "input validation" ) {
1178 
1179  SECTION( "state index" ) {
1180 
1181  int ind = GENERATE( -1, 1<<NUM_QUBITS );
1182  REQUIRE_THROWS_WITH( getProbAmp(vec,ind), Contains("Invalid amplitude index") );
1183  }
1184  SECTION( "density-matrix" ) {
1185 
1186  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
1187  REQUIRE_THROWS_WITH( getProbAmp(mat,0), Contains("valid only for state-vectors") );
1188  destroyQureg(mat, QUEST_ENV);
1189  }
1190  }
1191  destroyQureg(vec, QUEST_ENV);
1192 }
1193 
1194 
1195 
1200 TEST_CASE( "getRealAmp", "[calculations]" ) {
1201 
1202  Qureg vec = createQureg(NUM_QUBITS, QUEST_ENV);
1203 
1204  SECTION( "correctness" ) {
1205 
1206  SECTION( "state-vector" ) {
1207 
1208  initDebugState(vec);
1209  QVector ref = toQVector(vec);
1210 
1211  int ind = GENERATE( range(0,1<<NUM_QUBITS) );
1212  REQUIRE( getRealAmp(vec,ind) == real(ref[ind]) );
1213  }
1214  }
1215  SECTION( "input validation" ) {
1216 
1217  SECTION( "state index" ) {
1218 
1219  int ind = GENERATE( -1, 1<<NUM_QUBITS );
1220  REQUIRE_THROWS_WITH( getRealAmp(vec,ind), Contains("Invalid amplitude index") );
1221  }
1222  SECTION( "density-matrix" ) {
1223 
1224  Qureg mat = createDensityQureg(NUM_QUBITS, QUEST_ENV);
1225  REQUIRE_THROWS_WITH( getRealAmp(mat,0), Contains("valid only for state-vectors") );
1226  destroyQureg(mat, QUEST_ENV);
1227  }
1228  }
1229  destroyQureg(vec, QUEST_ENV);
1230 }
1231 
1232 
getProbAmp
qreal getProbAmp(Qureg qureg, long long int index)
Get the probability of a state-vector at an index in the full state vector.
Definition: QuEST.c:690
initBlankState
void initBlankState(Qureg qureg)
Initialises a qureg to have all-zero-amplitudes.
Definition: QuEST.c:117
pauliOpType
pauliOpType
Codes for specifying Pauli operators.
Definition: QuEST.h:96
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...
Definition: utilities.cpp:293
QUEST_ENV
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
fromComplex
#define fromComplex(comp)
PAULI_Z
@ PAULI_Z
Definition: QuEST.h:96
calcTotalProb
qreal calcTotalProb(Qureg qureg)
A debugging function which calculates the probability of the qubits in qureg being in any state,...
Definition: QuEST.c:822
PAULI_I
@ PAULI_I
Definition: QuEST.h:96
getDensityAmp
Complex getDensityAmp(Qureg qureg, long long int row, long long int col)
Get an amplitude from a density matrix at a given row and column.
Definition: QuEST.c:707
getRandomInt
int getRandomInt(int min, int max)
Returns a random integer between min (inclusive) and max (exclusive), from the uniform distribution.
Definition: utilities.cpp:481
getImagAmp
qreal getImagAmp(Qureg qureg, long long int index)
Get the imaginary component of the complex probability amplitude at an index in the state vector.
Definition: QuEST.c:683
NUM_QUBITS
#define NUM_QUBITS
The default number of qubits in the registers created for unit testing (both statevectors and density...
Definition: utilities.hpp:36
calcPurity
qreal calcPurity(Qureg qureg)
Calculates the purity of a density matrix, by the trace of the density matrix squared.
Definition: QuEST.c:855
calcFidelity
qreal calcFidelity(Qureg qureg, Qureg pureState)
Calculates the fidelity of qureg (a statevector or density matrix) against a reference pure state (ne...
Definition: QuEST.c:861
getRandomReal
qreal getRandomReal(qreal min, qreal max)
Returns a random real between min (inclusive) and max (exclusive), from the uniform distribution.
Definition: utilities.cpp:410
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...
Definition: utilities.cpp:159
qreal
#define qreal
toQMatrix
QMatrix toQMatrix(ComplexMatrix2 src)
Returns a copy of the given 2-by-2 matrix.
Definition: utilities.cpp:835
toQureg
void toQureg(Qureg qureg, QVector vec)
Initialises the state-vector qureg to have the same amplitudes as vec.
Definition: utilities.cpp:946
PAULI_X
@ PAULI_X
Definition: QuEST.h:96
getAmp
Complex getAmp(Qureg qureg, long long int index)
Get the complex amplitude at a given index in the state vector.
Definition: QuEST.c:697
calcHilbertSchmidtDistance
qreal calcHilbertSchmidtDistance(Qureg a, Qureg b)
Computes the Hilbert Schmidt distance between two density matrices a and b, defined as the Frobenius ...
Definition: QuEST.c:889
QVector
std::vector< qcomp > QVector
A complex vector, which can be zero-initialised with QVector(numAmps).
Definition: utilities.hpp:60
toQVector
QVector toQVector(Qureg qureg)
Returns an equal-size copy of the given state-vector qureg.
Definition: utilities.cpp:904
qcomp
#define qcomp
TEST_CASE
TEST_CASE("calcDensityInnerProduct", "[calculations]")
Definition: test_calculations.cpp:15
getRealAmp
qreal getRealAmp(Qureg qureg, long long int index)
Get the real component of the complex probability amplitude at an index in the state vector.
Definition: QuEST.c:676
PAULI_Y
@ PAULI_Y
Definition: QuEST.h:96
getNumQubits
int getNumQubits(Qureg qureg)
Get the number of qubits in a qureg object.
Definition: QuEST.c:666
calcInnerProduct
Complex calcInnerProduct(Qureg bra, Qureg ket)
Computes the inner product of two equal-size state vectors, given by.
Definition: QuEST.c:829
destroyQureg
void destroyQureg(Qureg qureg, QuESTEnv env)
Deallocate a Qureg object representing a set of qubits.
Definition: QuEST.c:75
getRandomStateVector
QVector getRandomStateVector(int numQb)
Returns a random numQb-length L2-normalised state-vector from an undisclosed distribution.
Definition: utilities.cpp:453
initDebugState
void initDebugState(Qureg qureg)
Initialises qureg to be in the un-normalised, non-physical state with with n-th complex amplitude (2n...
Definition: QuEST.c:1057
getRandomQMatrix
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:368
Qureg
Represents a system of qubits.
Definition: QuEST.h:160
getNumAmps
long long int getNumAmps(Qureg qureg)
Get the number of probability amplitudes in a qureg object, given by 2^numQubits.
Definition: QuEST.c:670
QMatrix
std::vector< std::vector< qcomp > > QMatrix
A complex square matrix.
Definition: utilities.hpp:49
calcExpecPauliSum
qreal calcExpecPauliSum(Qureg qureg, enum pauliOpType *allPauliCodes, qreal *termCoeffs, int numSumTerms, Qureg workspace)
Computes the expected value of a sum of products of Pauli operators.
Definition: QuEST.c:880
sublists
Catch::Generators::GeneratorWrapper< int * > sublists(int *list, int len, int sublen)
Returns a Catch2 generator of every length-sublen sublist of length-len list, in increasing lexograph...
Definition: utilities.cpp:1074
Qureg::numQubitsRepresented
int numQubitsRepresented
The number of qubits represented in either the state-vector or density matrix.
Definition: QuEST.h:165
Complex::real
qreal real
Definition: QuEST.h:105
createQureg
Qureg createQureg(int numQubits, QuESTEnv env)
Create a Qureg object representing a set of qubits which will remain in a pure state.
Definition: QuEST.c:34
calcProbOfOutcome
qreal calcProbOfOutcome(Qureg qureg, const int measureQubit, int outcome)
Gives the probability of a specified qubit being measured in the given outcome (0 or 1).
Definition: QuEST.c:845
Complex::imag
qreal imag
Definition: QuEST.h:106
getKroneckerProduct
QMatrix getKroneckerProduct(QMatrix a, QMatrix b)
Returns the kronecker product of a and b, where a and b are square but possibly differently-sized com...
Definition: utilities.cpp:169
getRandomDensityMatrix
QMatrix getRandomDensityMatrix(int numQb)
Returns a random numQb-by-numQb density matrix, from an undisclosed distribution, in a very mixed sta...
Definition: utilities.cpp:457
getZeroMatrix
QMatrix getZeroMatrix(size_t dim)
Returns a dim-by-dim square complex matrix, initialised to all zeroes.
Definition: utilities.cpp:143
Complex
Represents one complex number.
Definition: QuEST.h:103
getRandomQVector
QVector getRandomQVector(int dim)
Returns a dim-length vector with random complex amplitudes in the square joining {-1-i,...
Definition: utilities.cpp:420
initZeroState
void initZeroState(Qureg qureg)
Initialise a set of qubits to the classical zero state .
Definition: QuEST.c:111
calcDensityInnerProduct
qreal calcDensityInnerProduct(Qureg rho1, Qureg rho2)
Computes the Hilbert-Schmidt scalar product (which is equivalent to the Frobenius inner product of ma...
Definition: QuEST.c:837
createDensityQureg
Qureg createDensityQureg(int numQubits, QuESTEnv env)
Create a Qureg for qubits which are represented by a density matrix, and can be in mixed states.
Definition: QuEST.c:48
applyReferenceOp
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:573
initPlusState
void initPlusState(Qureg qureg)
Initialise a set of qubits to the plus state (and similarly for density matrices).
Definition: QuEST.c:123
calcExpecPauliProd
qreal calcExpecPauliProd(Qureg qureg, int *targetQubits, enum pauliOpType *pauliCodes, int numTargets, Qureg workspace)
Computes the expected value of a product of Pauli operators.
Definition: QuEST.c:871