The Quantum Exact Simulation Toolkit v4.2.0
Loading...
Searching...
No Matches
operations.cpp
1/** @file
2 * Unit tests of both multiplication and operations modules, since
3 * they use inextricable testing logic. Beware that because the
4 * operation functions have so much interface and test-semantic
5 * overlap (e.g. the logic of control qubits, of control-states,
6 * of density matrix variants), this file has opted to make
7 * extensive use of generics and metaprogramming, to avoid the
8 * code duplication of defining each function an independent test.
9 * This may have been a mistake, and the file now full of spaghetti
10 * comprised of advanced, misused C++ facilities. View at your own
11 * peril!
12 *
13 * @author Tyson Jones
14 *
15 * @defgroup unitops Operations
16 * @ingroup unittests
17 *
18 * @defgroup unitmult Multiplication
19 * @ingroup unittests
20 */
21
22#include "quest.h"
23
24#include <catch2/catch_test_macros.hpp>
25#include <catch2/catch_template_test_macros.hpp>
26#include <catch2/matchers/catch_matchers_string.hpp>
27#include <catch2/generators/catch_generators_range.hpp>
28
29#include "tests/utils/config.hpp"
30#include "tests/utils/cache.hpp"
31#include "tests/utils/qvector.hpp"
32#include "tests/utils/qmatrix.hpp"
33#include "tests/utils/compare.hpp"
34#include "tests/utils/convert.hpp"
35#include "tests/utils/evolve.hpp"
36#include "tests/utils/linalg.hpp"
37#include "tests/utils/lists.hpp"
38#include "tests/utils/measure.hpp"
39#include "tests/utils/macros.hpp"
40#include "tests/utils/random.hpp"
41
42#include <tuple>
43
44using std::tuple;
45using Catch::Matchers::ContainsSubstring;
46
47
48/*
49 * INTERNAL QUEST FUNCITONS
50 */
51
52extern int paulis_getPauliAt(PauliStr str, int ind);
53
54
55/*
56 * UTILITIES
57 */
58
59#define TEST_CATEGORY_OPS \
60 LABEL_UNIT_TAG "[operations]"
61
62#define TEST_CATEGORY_MULT \
63 LABEL_UNIT_TAG "[multiplication]"
64
65
66/*
67 * reference operator matrices used by testing
68 */
69
70namespace FixedMatrices {
71
72 qmatrix H = {
73 {1/std::sqrt(2), 1/std::sqrt(2)},
74 {1/std::sqrt(2), -1/std::sqrt(2)}};
75
76 qmatrix X = getPauliMatrix(1);
77 qmatrix Y = getPauliMatrix(2);
78 qmatrix Z = getPauliMatrix(3);
79
80 qreal PI = 3.14159265358979323846;
81 qmatrix T = {
82 {1, 0},
83 {0, std::exp(1_i * PI/4)}};
84
85 qmatrix S = {
86 {1, 0},
87 {0, 1_i}};
88
89 qmatrix SWAP = {
90 {1, 0, 0, 0},
91 {0, 0, 1, 0},
92 {0, 1, 0, 0},
93 {0, 0, 0, 1}};
94
95 qmatrix sqrtSWAP = {
96 {1, 0, 0, 0},
97 {0, (1+1_i)/2, (1-1_i)/2, 0},
98 {0, (1-1_i)/2, (1+1_i)/2, 0},
99 {0, 0, 0, 1}};
100}
101
102namespace ParameterisedMatrices {
103
104 auto Rx = [](qreal p) { return getExponentialOfPauliMatrix(p, FixedMatrices::X); };
105 auto Ry = [](qreal p) { return getExponentialOfPauliMatrix(p, FixedMatrices::Y); };
106 auto Rz = [](qreal p) { return getExponentialOfPauliMatrix(p, FixedMatrices::Z); };
107
108 auto PS = [](qreal p) { return qmatrix{{1, 0}, {0, std::exp(p*1_i)}}; };
109 auto PS2 = [](qreal p) { return getControlledMatrix(PS(p), 1); };
110}
111
112namespace VariableSizeMatrices {
113
114 auto X = [](int n) { return getKroneckerProduct(FixedMatrices::X, n); };
115 auto PF = [](int n) { return getControlledMatrix(FixedMatrices::Z, n - 1); };
116}
117
118namespace VariableSizeParameterisedMatrices {
119
120 auto Z = [](qreal p, int n) {
121 qmatrix m = getKroneckerProduct(FixedMatrices::Z, n);
122 return getExponentialOfPauliMatrix(p, m);
123 };
124
125 auto PS = [](qreal p, int n) {
126 qmatrix m = ParameterisedMatrices::PS(p);
127 return getControlledMatrix(m, n - 1);
128 };
129}
130
131
132/*
133 * execute 'function' upon each kind of cached qureg
134 * (e.g. distributed, GPU-accelerated, etc) and a
135 * reference state (T1 = qvector or qmatrix), when
136 * both are initialised in the debug state, and
137 * thereafter confirm they approximately agree. Each
138 * qureg deployment is featured in a separate test
139 * section, so are accounted distinctly.
140 */
141
142void TEST_ON_CACHED_QUREGS(quregCache quregs, auto& reference, auto& function) {
143
144 for (auto& [label, qureg]: quregs) {
145
146 DYNAMIC_SECTION( label ) {
147
148 // no need to validate whether qureg successfully
149 // enters the debug state here, because the below
150 // serial setToDebugState() is guaranteed to succeed
151 initDebugState(qureg);
152 setToDebugState(reference);
153
154 function(qureg, reference);
155 REQUIRE_AGREE( qureg, reference );
156 }
157 }
158}
159
160void TEST_ON_CACHED_QUREG_AND_MATRIX(quregCache quregs, matrixCache matrices, auto apiFunc, auto refState, auto refMatr, auto refFunc) {
161
162 for (auto& [labelA, qureg]: quregs) {
163 for (auto& [labelB, matrix]: matrices) {
164
165 // skip illegal (distributed matrix, local qureg) combo
166 if (matrix.isDistributed && ! qureg.isDistributed)
167 continue;
168
169 DYNAMIC_SECTION( labelA + LABEL_DELIMITER + labelB ) {
170
171 // set qureg and reference to debug
172 initDebugState(qureg);
173 setToDebugState(refState);
174
175 // set API matrix to pre-initialised ref matrix
176 setFullStateDiagMatr(matrix, 0, getDiagonals(refMatr));
177
178 // API and reference functions should produce agreeing states
179 apiFunc(qureg, matrix);
180 refFunc(refState, refMatr);
181 REQUIRE_AGREE( qureg, refState );
182 }
183 }
184 }
185}
186
187
188/*
189 * simply avoids boilerplate
190 */
191
192#define PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef ) \
193 int numQubits = getNumCachedQubits(); \
194 auto statevecQuregs = getCachedStatevecs(); \
195 auto densmatrQuregs = getCachedDensmatrs(); \
196 qvector statevecRef = getZeroVector(getPow2(numQubits)); \
197 qmatrix densmatrRef = getZeroMatrix(getPow2(numQubits));
198
199
200/*
201 * Template flags for specifying what kind of additional
202 * arguments (in addition to ctrls/states/targs below) are
203 * accepted by an API operation, when passing said operation
204 * to automated testing facilities. This is NOT consulted
205 * when generically invoking the API operation (we instead
206 * used variadic templates above for that), but IS used
207 * by the testing code to decide how to prepare inputs.
208 *
209 * For example:
210 * - applyHadamard: none
211 * - applyRotateX: scalar
212 * - applyRotateAroundAxis: axisrots
213 * - applyDiagMatr1: diagmatr
214 * - applyDiagMatrPower: diagpower
215 * - applyCompMatr: compmatr
216 * - applyPauliStr: paulistr
217 * - applyPauliGadget: pauligad
218*/
219
220enum ArgsFlag { none, scalar, axisrots, diagmatr, diagpower, compmatr, paulistr, pauligad };
221
222
223/*
224 * Template flags for specifying how many control and
225 * target qubits are accepted by an API operation.
226 * Value 'anystates' is reserved for control qubits,
227 * indicating that ctrls must accompany ctrl-states in
228 * the API signature.
229 *
230 * For example:
231 * - applyHadamard: Ctrls=zero, Targs=one
232 * - applyControlledSwap: Ctrls=one, Targs=two
233 * - applyMultiControlledCompMatr: Ctrls=any, Targs=any
234 * - applyMultiStateControlledT: Ctrls=anystates, Targs=one
235 */
236
237enum NumQubitsFlag { zero, one, two, any, anystates };
238
239void assertNumQubitsFlagsAreValid(NumQubitsFlag ctrlsFlag, NumQubitsFlag targsFlag) {
240
241 DEMAND(
242 ctrlsFlag == zero ||
243 ctrlsFlag == one ||
244 ctrlsFlag == any ||
245 ctrlsFlag == anystates );
246
247 DEMAND(
248 targsFlag == zero ||
249 targsFlag == one ||
250 targsFlag == two ||
251 targsFlag == any);
252}
253
254void assertNumQubitsFlagsValid(
255 NumQubitsFlag ctrlsFlag, NumQubitsFlag targsFlag,
256 vector<int> ctrls, vector<int> states, vector<int> targs
257) {
258 assertNumQubitsFlagsAreValid(ctrlsFlag, targsFlag);
259
260 // we explicitly permit targsFlag=zero while
261 // targs.size() != 0, which occurs when targs
262 // are supplied to the API through an alternate
263 // argument (e.g. a PauliStr)
264
265 if (targsFlag == one)
266 DEMAND( targs.size() == 1 );
267
268 if (targsFlag == two)
269 DEMAND( targs.size() == 2 );
270
271 if (ctrlsFlag == zero)
272 DEMAND( ctrls.size() == 0 );
273
274 if (ctrlsFlag == one)
275 DEMAND( ctrls.size() == 1 );
276
277 if (ctrlsFlag == anystates)
278 DEMAND( states.size() == ctrls.size() );
279 else
280 DEMAND( states.size() == 0 );
281}
282
283
284/*
285 * extract the runtime values of the number of control
286 * and target qubtits from their compile-time templates.
287 * When their number is permitted to be 'any', we use
288 * a generator to successively test all possible numbers.
289 * As such, this function is named similar to a Catch2
290 * macro so the caller recognises it is a generator.
291 */
292
293template <NumQubitsFlag Ctrls, NumQubitsFlag Targs, ArgsFlag Args>
294int GENERATE_NUM_TARGS(int numQuregQubits) {
295
296 assertNumQubitsFlagsAreValid(zero, Targs);
297 DEMAND( Targs != one || numQuregQubits >= 1 );
298 DEMAND( Targs != two || numQuregQubits >= 2 );
299 DEMAND( numQuregQubits > 0 );
300
301 // single choice if #targs is compile-time set
302 if constexpr (Targs == one)
303 return 1;
304 if constexpr (Targs == two)
305 return 2;
306
307 if constexpr (Targs == any) {
308
309 // we can target all qubits...
310 int maxNumTargs = numQuregQubits;
311
312 // unless we're applying CompMatr in distributed
313 // mode. Technically we can support all targets
314 // if the Qureg is a density matrix or not
315 // distributed, but we safely choose the min here
316 if (Args == compmatr)
317 maxNumTargs = numQuregQubits - getLog2(getQuESTEnv().numNodes);
318
319 // we must also ensure there is space left for a forced ctrl
320 if (Ctrls == one && maxNumTargs == numQuregQubits)
321 maxNumTargs -= 1;
322
323 return GENERATE_COPY( range(1, maxNumTargs+1) );
324 }
325}
326
327template <NumQubitsFlag Ctrls>
328int GENERATE_NUM_CTRLS(int numFreeQubits) {
329
330 assertNumQubitsFlagsAreValid(Ctrls, one);
331 DEMAND( Ctrls != one || numFreeQubits >= 1 );
332 DEMAND( numFreeQubits >= 0 );
333
334 if constexpr (Ctrls == zero)
335 return 0;
336
337 if constexpr (Ctrls == one)
338 return 1;
339
340 if constexpr (Ctrls == any || Ctrls == anystates)
341 return GENERATE_COPY( range(0, numFreeQubits+1) );
342}
343
344
345/*
346 * invoke an API operation (e.g. applyHadamard), passing
347 * any elements of (ctrls,states,targs) that it accepts
348 * (as informed by the template values) along with any
349 * additional arguments. We explicitly accept numCtrls
350 * and numTargs, rather than inferring them from ctrls
351 * and targs, such that invalid numbers (like -1) can be
352 * passed by input validation testing.
353 *
354 * This big, ugly bespoke function is necessary, rather
355 * than a simple variadic template, because the QuEST
356 * API accepts fixed numbers of qubits as individual
357 * arguments, rather than as lists/vectors/pointers. Note
358 * that our use of variadic templates (args) means we do
359 * not need to include ArgsFlag as a template parameter.
360 */
361
362template <NumQubitsFlag Ctrls, NumQubitsFlag Targs>
363void invokeApiOperation(
364 auto operation, Qureg qureg,
365 vector<int> ctrls, vector<int> states, int numCtrls,
366 vector<int> targs, int numTargs,
367 auto&... args
368) {
369 assertNumQubitsFlagsValid(Ctrls, Targs, ctrls, states, targs);
370
371 if constexpr (Ctrls == zero) {
372 if constexpr (Targs == zero) operation(qureg, args...);
373 if constexpr (Targs == one) operation(qureg, targs[0], args...);
374 if constexpr (Targs == two) operation(qureg, targs[0], targs[1], args...);
375 if constexpr (Targs == any) operation(qureg, targs.data(), numTargs, args...);
376 }
377 if constexpr (Ctrls == one) {
378 if constexpr (Targs == zero) operation(qureg, ctrls[0], args...);
379 if constexpr (Targs == one) operation(qureg, ctrls[0], targs[0], args...);
380 if constexpr (Targs == two) operation(qureg, ctrls[0], targs[0], targs[1], args...);
381 if constexpr (Targs == any) operation(qureg, ctrls[0], targs.data(), numTargs, args...);
382 }
383 if constexpr (Ctrls == any) {
384 if constexpr (Targs == zero) operation(qureg, ctrls.data(), numCtrls, args...);
385 if constexpr (Targs == one) operation(qureg, ctrls.data(), numCtrls, targs[0], args...);
386 if constexpr (Targs == two) operation(qureg, ctrls.data(), numCtrls, targs[0], targs[1], args...);
387 if constexpr (Targs == any) operation(qureg, ctrls.data(), numCtrls, targs.data(), numTargs, args...);
388 }
389 if constexpr (Ctrls == anystates) {
390 if constexpr (Targs == zero) operation(qureg, ctrls.data(), states.data(), numCtrls, args...);
391 if constexpr (Targs == one) operation(qureg, ctrls.data(), states.data(), numCtrls, targs[0], args...);
392 if constexpr (Targs == two) operation(qureg, ctrls.data(), states.data(), numCtrls, targs[0], targs[1], args...);
393 if constexpr (Targs == any) operation(qureg, ctrls.data(), states.data(), numCtrls, targs.data(), numTargs, args...);
394 }
395}
396
397// overload to avoid passing numCtrls and numTargs
398
399template <NumQubitsFlag Ctrls, NumQubitsFlag Targs>
400void invokeApiOperation(auto operation, Qureg qureg, vector<int> ctrls, vector<int> states, vector<int> targs, auto&... args) {
401 invokeApiOperation<Ctrls,Targs>(operation, qureg, ctrls, states, ctrls.size(), targs, targs.size(), args...);
402}
403
404
405/*
406 * prepare an API matrix (e.g. CompMatr1), as per
407 * the given template parameters. Depending on the
408 * elemsFlag, the elements can be initialised to
409 * zero (0), the identity matrix (1), or randomly (2).
410 * This is used for testing API functions which accept
411 * matrices, in a function-agnostic way.
412 */
413
414template <NumQubitsFlag Targs, ArgsFlag Args>
415auto getRandomOrIdentityApiMatrix(int numTargs, int elemsFlag) {
416
417 DEMAND(
418 Args == diagmatr ||
419 Args == diagpower ||
420 Args == compmatr );
421 DEMAND(
422 elemsFlag == 0 ||
423 elemsFlag == 1 ||
424 elemsFlag == 2 );
425
426 qmatrix qm;
427 if (elemsFlag == 0)
428 qm = getZeroMatrix(getPow2(numTargs));
429 if (elemsFlag == 1)
430 qm = getIdentityMatrix(getPow2(numTargs));
431 if (elemsFlag == 2)
432 qm = (Args == compmatr)?
433 getRandomUnitary(numTargs) :
434 getRandomDiagonalUnitary(numTargs);
435
436 if constexpr (Args == compmatr && Targs == one)
437 return getCompMatr1(qm);
438
439 if constexpr (Args == compmatr && Targs == two)
440 return getCompMatr2(qm);
441
442 if constexpr (Args == compmatr && Targs == any) {
443 CompMatr cm = createCompMatr(numTargs); // must be freed
444 setCompMatr(cm, qm);
445 return cm;
446 }
447
448 qvector dv = getDiagonals(qm);
449 constexpr bool diag = (Args == diagmatr || Args == diagpower);
450
451 if constexpr (diag && Targs == one)
452 return getDiagMatr1(dv);
453
454 if constexpr (diag && Targs == two)
455 return getDiagMatr2(dv);
456
457 if constexpr (diag && Targs == any) {
458 DiagMatr dm = createDiagMatr(numTargs); // must be freed
459 setDiagMatr(dm, dv);
460 return dm;
461 }
462}
463
464template <NumQubitsFlag Targs, ArgsFlag Args> auto getZeroApiMatrix (int numTargs) { return getRandomOrIdentityApiMatrix<Targs,Args>(numTargs, 0); }
465template <NumQubitsFlag Targs, ArgsFlag Args> auto getIdentityApiMatrix(int numTargs) { return getRandomOrIdentityApiMatrix<Targs,Args>(numTargs, 1); }
466template <NumQubitsFlag Targs, ArgsFlag Args> auto getRandomApiMatrix (int numTargs) { return getRandomOrIdentityApiMatrix<Targs,Args>(numTargs, 2); }
467
468
469/*
470 * chooses random values for the remaining arguments
471 * (after controls/states/targets) to API operations,
472 * with types informed by the template parameters.
473 * For example:
474 * - applyRotateX() accepts a scalar
475 * - applyCompMatr() accepts a CompMatr
476 * - applyDiagMatrPower() accepts a DiagMatr and qcomp
477 */
478
479template <NumQubitsFlag Targs, ArgsFlag Args>
480auto getRandomRemainingArgs(vector<int> targs) {
481
482 if constexpr (Args == none)
483 return tuple{ };
484
485 if constexpr (Args == scalar) {
486 qreal angle = getRandomPhase();
487 return tuple{ angle };
488 }
489
490 if constexpr (Args == axisrots) {
491 qreal angle = getRandomPhase();
492 qreal x = getRandomReal(-1, 1);
493 qreal y = getRandomReal(-1, 1);
494 qreal z = getRandomReal(-1, 1);
495 return tuple{ angle, x, y, z };
496 }
497
498 if constexpr (Args == compmatr || Args == diagmatr) {
499 auto matrix = getRandomApiMatrix<Targs,Args>(targs.size()); // allocates heap mem
500 return tuple{ matrix };
501 }
502
503 if constexpr (Args == diagpower) {
504 DiagMatr matrix = getRandomApiMatrix<Targs,Args>(targs.size()); // allocates heap mem
505 qcomp exponent = qcomp(getRandomReal(-3, 3), 0); // real for unitarity
506 return tuple{ matrix, exponent };
507 }
508
509 if constexpr (Args == paulistr) {
510 PauliStr str = getRandomPauliStr(targs);
511 return tuple{ str };
512 }
513
514 if constexpr (Args == pauligad) {
515 PauliStr str = getRandomPauliStr(targs);
516 qreal angle = getRandomPhase();
517 return tuple{ str, angle };
518 }
519}
520
521
522template <NumQubitsFlag Targs, ArgsFlag Args>
523void freeRemainingArgs(auto args) {
524
525 if constexpr (Targs == any && Args == compmatr)
526 destroyCompMatr(std::get<0>(args));
527
528 if constexpr (Targs == any && Args == diagmatr)
529 destroyDiagMatr(std::get<0>(args));
530
531 if constexpr (Targs == any && Args == diagpower)
532 destroyDiagMatr(std::get<0>(args));
533}
534
535
536/*
537 * unpack the given reference operator matrix (a qmatrix)
538 * for an API operation, which is passed to testOperation(),
539 * and which will be effected upon the reference state (a
540 * qvector or qmatrix). The type/form of matrixRefGen depends
541 * on the type of API operation, indicated by template parameter.
542 */
543
544template <NumQubitsFlag Targs, ArgsFlag Args>
545qmatrix getReferenceMatrix(auto matrixRefGen, vector<int> targs, auto additionalArgs) {
546
547 if constexpr (Args == none && Targs != any)
548 return matrixRefGen;
549
550 if constexpr (Args == none && Targs == any)
551 return matrixRefGen(targs.size());
552
553 if constexpr (Args == scalar && Targs != any) {
554 qreal angle = std::get<0>(additionalArgs);
555 return matrixRefGen(angle);
556 }
557
558 if constexpr (Args == scalar && Targs == any) {
559 qreal angle = std::get<0>(additionalArgs);
560 return matrixRefGen(angle, targs.size());
561 }
562
563 if constexpr (Args == axisrots) {
564 qreal angle = std::get<0>(additionalArgs);
565 qreal x = std::get<1>(additionalArgs);
566 qreal y = std::get<2>(additionalArgs);
567 qreal z = std::get<3>(additionalArgs);
568 return getExponentialOfNormalisedPauliVector(angle, x, y, z);
569 }
570
571 if constexpr (Args == compmatr || Args == diagmatr) {
572 auto apiMatrix = std::get<0>(additionalArgs);
573 return getMatrix(apiMatrix);
574 }
575
576 if constexpr (Args == diagpower) {
577 auto apiMatrix = std::get<0>(additionalArgs);
578 qmatrix diag = getMatrix(apiMatrix);
579 qcomp power = std::get<1>(additionalArgs);
580 return getPowerOfDiagonalMatrix(diag, power);
581 }
582
583 if constexpr (Args == paulistr) {
584 PauliStr str = std::get<0>(additionalArgs);
585 return getMatrix(str, targs);
586 }
587
588 if constexpr (Args == pauligad) {
589 PauliStr str = std::get<0>(additionalArgs);
590 qreal angle = std::get<1>(additionalArgs);
591 qmatrix matr = getMatrix(str, targs);
592 return getExponentialOfPauliMatrix(angle, matr);
593 }
594}
595
596
597/*
598 * Template parameters which specify how the reference
599 * operatorshould be applied upon the reference state.
600 * Let |psi> be a statevector, rho be a density matrix,
601 * and matr be an operator matrix. The options perform:
602 *
603 * apply: |psi> -> matr |psi>, rho -> matr rho adj(matr)
604 * leftapply: |psi> -> matr |psi>, rho -> matr rho
605 * rightapply: rho -> rho matr
606 *
607 * Note this is necessarily a template parameter (rather
608 * than just a runtime parameter) only because the
609 * rightapplyReferenceOperator() function is defined
610 * only upon qmatrix (for density matrices)
611 */
612
613enum ApplyFlag { apply, leftapply, rightapply };
614
615
616/*
617 * display all/only relevant inputs given to an
618 * API operation when its subsequent test fails.
619 * This is like a customisation of CAPTURE, although
620 * we must use UNSCOPED_INFO (and ergo re-implement
621 * some printers) because our branching makes scopes
622 * which end CAPTURE's lifetime.
623 */
624
625
626// @todo surely this should live somewhere else,
627// and/or re-use printer_ functions as much as possible
628
629std::string toString(vector<int> list) {
630
631 std::string out = "{ ";
632 for (int& e : list)
633 out += std::to_string(e) + " ";
634 out += "}";
635 return out;
636}
637
638std::string toString(PauliStr str, vector<int> targs) {
639
640 std::string labels = "IXYZ";
641
642 // ugly but adequate - like me (call me)
643 std::string out = "";
644 for (int i=targs.size()-1; i>=0; i--)
645 out += labels[paulis_getPauliAt(str, targs[i])];
646
647 return out;
648}
649
650template <NumQubitsFlag Ctrls, NumQubitsFlag Targs, ArgsFlag Args>
651void CAPTURE_RELEVANT( vector<int> ctrls, vector<int> states, vector<int> targs, auto& args ) {
652
653 // display ctrls
654 if (Ctrls == one)
655 UNSCOPED_INFO( "control := " << ctrls[0] );
656 if (Ctrls == any || Ctrls == anystates )
657 UNSCOPED_INFO( "controls := " << toString(ctrls) );
658
659 // display states
660 if (Ctrls == anystates)
661 UNSCOPED_INFO( "states := " << toString(states) );
662
663 // display targs
664 if (Targs == one)
665 UNSCOPED_INFO( "target := " << targs[0] );
666
667 if (Targs == two) {
668 UNSCOPED_INFO( "target A := " << targs[0] );
669 UNSCOPED_INFO( "target B := " << targs[1] );
670 }
671 if (Targs == any)
672 UNSCOPED_INFO( "targets := " << toString(targs) );
673
674 // display rotation angle
675 if constexpr (Args == scalar)
676 UNSCOPED_INFO( "angle := " << std::get<0>(args) );
677
678 // display rotation angle and axis
679 if constexpr (Args == axisrots) {
680 UNSCOPED_INFO( "angle := " << std::get<0>(args) );
681 UNSCOPED_INFO( "x := " << std::get<1>(args) );
682 UNSCOPED_INFO( "y := " << std::get<2>(args) );
683 UNSCOPED_INFO( "z := " << std::get<3>(args) );
684 }
685
686 // note but don't display API matrices
687 if constexpr (Args == compmatr || Args == diagmatr || Args == diagpower)
688 UNSCOPED_INFO( "matrix := (omitted)" );
689
690 // display exponent of diagonal matrix
691 if constexpr (Args == diagpower) {
692 qcomp p = std::get<1>(args);
693 UNSCOPED_INFO( "exponent := " << std::real(p) << " + (" << std::imag(p) << ")i" );
694 }
695
696 // display PauliStr
697 if constexpr (Args == paulistr || Args == pauligad)
698 UNSCOPED_INFO( "paulis := " << toString(std::get<0>(args), targs) );
699
700 // display PauliStr angle
701 if constexpr (Args == pauligad)
702 UNSCOPED_INFO( "angle := " << std::get<1>(args) );
703}
704
705
706/*
707 * test the correctness of an API operation. The
708 * template parameters are compile-time clues
709 * about what inputs to prepare and pass to the
710 * operation, and how its reference matrix (arg
711 * matrixRefGen) is formatted.
712 */
713
714template <NumQubitsFlag Ctrls, NumQubitsFlag Targs, ArgsFlag Args, ApplyFlag Apply>
715void testOperationCorrectness(auto operation, auto matrixRefGen) {
716
717 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
718
719 // try all possible number of ctrls and targs
720 int numTargs = GENERATE_NUM_TARGS<Ctrls,Targs,Args>(numQubits);
721 int numCtrls = GENERATE_NUM_CTRLS<Ctrls>(numQubits - numTargs);
722
723 // either try all possible ctrls and targs, or randomise them
724 auto ctrlsAndTargs = GENERATE_CTRLS_AND_TARGS( numQubits, numCtrls, numTargs );
725 vector<int> ctrls = std::get<0>(ctrlsAndTargs);
726 vector<int> targs = std::get<1>(ctrlsAndTargs);
727
728 // randomise control states (if operation accepts them)
729 vector<int> states = getRandomOutcomes(numCtrls * (Ctrls == anystates));
730
731 // randomise remaining operation parameters
732 auto primaryArgs = tuple{ ctrls, states, targs };
733 auto furtherArgs = getRandomRemainingArgs<Targs,Args>(targs); // may allocate heap memory
734
735 // obtain the reference matrix for this operation
736 qmatrix matrixRef = getReferenceMatrix<Targs,Args>(matrixRefGen, targs, furtherArgs);
737
738 // PauliStr arg replaces target qubit list in API operations
739 constexpr NumQubitsFlag RevTargs = (Args==paulistr||Args==pauligad)? zero : Targs;
740
741 // disabling unitary-check validation for compmatr, since it's hard to
742 // generate numerically-precise random unitaries upon many qubits, or
743 // upon few qubits are single-precision. So we disable completely until
744 // we re-implement 'input validation' checks which force us to fix thresholds
745 (Args == compmatr)?
748
749 // prepare test function which will receive both statevectors and density matrices
750 auto testFunc = [&](Qureg qureg, auto& stateRef) -> void {
751
752 // invoke API operation, passing all args (unpacking variadic)
753 auto apiFunc = [](auto&&... args) { return invokeApiOperation<Ctrls,RevTargs>(args...); };
754 auto allArgs = std::tuple_cat(tuple{operation, qureg}, primaryArgs, furtherArgs);
755 std::apply(apiFunc, allArgs);
756
757 // update reference state (ctrls & states happen to only ever be used by apply)
758 if constexpr (Apply == apply) applyReferenceOperator( stateRef, ctrls, states, targs, matrixRef);
759 if constexpr (Apply == leftapply) leftapplyReferenceOperator( stateRef, ctrls, states, targs, matrixRef);
760 if constexpr (Apply == rightapply) rightapplyReferenceOperator(stateRef, ctrls, states, targs, matrixRef);
761 };
762
763 // report operation's input parameters if any subsequent test fails
764 CAPTURE_RELEVANT<Ctrls,Targs,Args>( ctrls, states, targs, furtherArgs );
765
766 // test API operation on all available deployment combinations (e.g. OMP, MPI, MPI+GPU, etc),
767 // though the postMultiply*() functions do not accept statevectors
768 if constexpr (Apply != rightapply) {
769 SECTION( LABEL_STATEVEC ) {
770 TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc);
771 }
772 }
773 SECTION( LABEL_DENSMATR ) {
774 TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc);
775 }
776
777 // free any heap-alloated API matrices and restore epsilon
778 freeRemainingArgs<Targs,Args>(furtherArgs);
780}
781
782
783/*
784 * test the input valdiation of an API operation.
785 * The template parameters are compile-time clues
786 * about what inputs are accepted by the operation
787 */
788
789template <NumQubitsFlag Ctrls, NumQubitsFlag Targs, ArgsFlag Args>
790auto getFixedCtrlsStatesTargs(int numQubits) {
791
792 // default each to empty
793 vector<int> targs, ctrls, states;
794
795 // assign when non-empty
796 if constexpr (Targs == one) targs = {0};
797 if constexpr (Targs == two) targs = {0,1};
798 if constexpr (Targs == any) targs = {0,1,2};
799 if constexpr (Ctrls == one) ctrls = {3};
800 if constexpr (Ctrls == any) ctrls = {3,4};
801 if constexpr (Ctrls == anystates) ctrls = {3,4};
802 if constexpr (Ctrls == anystates) states = {0,0};
803
804 DEMAND( numQubits >= targs.size() + ctrls.size() );
805
806 return tuple{ ctrls, states, targs };
807}
808
809template <NumQubitsFlag Targs, ArgsFlag Args>
810auto getFixedRemainingArgs(vector<int> targs) {
811
812 // getPauliStr uses gives length-3 hardcoded string
813 if constexpr (Args == paulistr || Args == pauligad)
814 DEMAND( targs.size() == 3 );
815
816 if constexpr (Args == none) return tuple{ };
817 if constexpr (Args == scalar) return tuple{ 0 }; // angle
818 if constexpr (Args == axisrots) return tuple{ 0, 1,1,1 }; // (angle,x,y,z)
819 if constexpr (Args == compmatr) return tuple{ getIdentityApiMatrix<Targs,Args>(targs.size()) }; // id
820 if constexpr (Args == diagmatr) return tuple{ getIdentityApiMatrix<Targs,Args>(targs.size()) }; // id
821 if constexpr (Args == diagpower) return tuple{ getIdentityApiMatrix<Targs,Args>(targs.size()), qcomp(1,0) }; // (id, exponent)
822 if constexpr (Args == paulistr) return tuple{ getPauliStr("XXX", targs) }; // XXX
823 if constexpr (Args == pauligad) return tuple{ getPauliStr("XXX", targs), 0 }; // (XXX, angle)
824}
825
826template <NumQubitsFlag Ctrls, NumQubitsFlag Targs, ArgsFlag Args, ApplyFlag Apply>
827void testOperationValidation(auto operation) {
828
829 // use any cached Qureg (though postMultiply*() functions accept only density matrices)
830 Qureg qureg = getCachedDensmatrs().begin()->second;
831
832 // in lieu of preparing random inputs like testOperationCorrectness()
833 // above, we instead obtain simple, fixed, compatible inputs
834 auto [ctrls,states,targs] = getFixedCtrlsStatesTargs<Ctrls,Targs,Args>(qureg.numQubits);
835 auto furtherArgs = getFixedRemainingArgs<Targs,Args>(targs);
836
837 // calling apiFunc() will pass the above args with their call-time values
838 auto apiFunc = [&]() {
839 constexpr NumQubitsFlag RevTargs = (Args==paulistr||Args==pauligad)? zero : Targs;
840 auto func = [](auto&&... allArgs) { return invokeApiOperation<Ctrls,RevTargs>(allArgs...); };
841 std::apply(func, std::tuple_cat(tuple{operation, qureg, ctrls, states, targs}, furtherArgs));
842 };
843
844 // convenience vars
845 int numQubits = qureg.numQubits;
846 int numTargs = (int) targs.size();
847 int numCtrls = (int) ctrls.size();
848
849 /// @todo
850 /// below, we return from intendedly skipped SECTIONS which
851 /// appears to work (does not corrupt test statistics, and
852 /// does not attribute skipped tests to having passed the
853 /// section) but is an undocumented Catch2 trick. Check safe!
854
855 SECTION( "qureg uninitialised" ) {
856
857 // spoof uninitialised value
858 qureg.numQubits = -123;
859 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("invalid Qureg") );
860 }
861
862 SECTION( "invalid target" ) {
863
864 // not applicable (PauliStr already made from targs)
865 if (Args == paulistr || Args == pauligad)
866 return;
867
868 // sabotage a target
869 int ind = GENERATE_COPY( range(0,numTargs) );
870 int val = GENERATE_COPY( -1, numQubits );
871 targs[ind] = val;
872
873 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("Invalid target qubit") );
874 }
875
876 SECTION( "invalid non-Identity Pauli index" ) {
877
878 if (Args != paulistr && Args != pauligad)
879 return;
880
881 PauliStr badStr = getPauliStr("X", {numQubits + 1});
882 if constexpr (Args == paulistr) furtherArgs = tuple{ badStr };
883 if constexpr (Args == pauligad) furtherArgs = tuple{ badStr, 1 };
884
885 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("highest-index non-identity Pauli operator") && ContainsSubstring("exceeds the maximum target") );
886 }
887
888 SECTION( "invalid control" ) {
889
890 if (numCtrls == 0)
891 return;
892
893 // sabotage a ctrl
894 int ind = GENERATE_COPY( range(0,numCtrls) );
895 int val = GENERATE_COPY( -1, numQubits );
896 ctrls[ind] = val;
897
898 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("Invalid control qubit") );
899 }
900
901 SECTION( "control and target collision" ) {
902
903 if (numCtrls == 0)
904 return;
905
906 // sabotage a ctrl
907 int targInd = GENERATE_COPY( range(0,numTargs) );
908 int ctrlInd = GENERATE_COPY( range(0,numCtrls) );
909 ctrls[ctrlInd] = targs[targInd];
910
911 if (Args==paulistr||Args==pauligad)
912 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("control qubit overlaps a non-identity Pauli operator") );
913 else
914 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("qubit appeared among both the control and target qubits") );
915 }
916
917 SECTION( "control states" ) {
918
919 if (states.empty())
920 return;
921
922 int ind = GENERATE_COPY( range(0,numCtrls) );
923 int val = GENERATE( -1, 2 );
924 states[ind] = val;
925
926 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("invalid control-state") );
927 }
928
929 SECTION( "repetition in controls" ) {
930
931 if (numCtrls < 2)
932 return;
933
934 int ind = GENERATE_COPY( range(1,numCtrls) );
935 ctrls[ind] = ctrls[0];
936
937 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("control qubits contained duplicates") );
938 }
939
940 SECTION( "repetition in targets" ) {
941
942 // not applicable to Pauli functions (getPauliStr would throw)
943 if (Args==paulistr||Args==pauligad)
944 return;
945
946 if (numTargs < 2)
947 return;
948
949 int ind = GENERATE_COPY( range(1,numTargs) );
950 targs[ind] = targs[0];
951
952 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("target qubits contained duplicates") );
953 }
954
955 SECTION( "number of targets" ) {
956
957 // not applicable to Pauli functions (getPauliStr would run fine,
958 // and the runtime error would be about the non-identity index)
959 if (Args == paulistr || Args == pauligad)
960 return;
961
962 if (Targs != any)
963 return;
964
965 // too few (cannot test less than 0)
966 targs = {};
967 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("Must specify one or more targets") );
968
969 // too many; exceeds Qureg
970 targs = getRange(qureg.numQubits+1);
971 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("number of target qubits") && ContainsSubstring("exceeds the number of qubits in the Qureg") );
972
973 // note numTargs + numCtrls > numQubits is caught by
974 // invalid index or overlapping (ctrls,targs) validation
975 }
976
977 SECTION( "mismatching matrix size" ) {
978
979 // only relevant to variable-sized matrices
980 if (Targs != any)
981 return;
982 if (!(Args == compmatr || Args == diagmatr || Args == diagpower))
983 return;
984
985 DEMAND( numTargs > 1 );
986 DEMAND( numCtrls + numTargs < numQubits );
987
988 targs.push_back(numQubits - 1);
989 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("matrix has an inconsistent size") );
990
991 targs.pop_back();
992 targs.pop_back();
993 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("matrix has an inconsistent size") );
994 }
995
996 SECTION( "matrix unitarity" ) {
997
998 // only relevant to matrix functions...
999 if (Args != compmatr && Args != diagmatr && Args != diagpower)
1000 return;
1001
1002 // which enforce unitarity
1003 if (Apply != apply)
1004 return;
1005
1006 if constexpr (Args == compmatr || Args == diagmatr)
1007 furtherArgs = tuple{ getZeroApiMatrix<Targs,Args>(targs.size()) };
1008 if constexpr (Args == diagpower)
1009 furtherArgs = tuple{ getZeroApiMatrix<Targs,Args>(targs.size()), 1 };
1010
1011 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("unitary") );
1012 }
1013
1014 SECTION( "matrix uninitialised" ) {
1015
1016 // sabotage matrix struct field
1017 if constexpr (Args == compmatr || Args == diagmatr || Args == diagpower)
1018 std::get<0>(furtherArgs).numQubits = -1;
1019 else
1020 return;
1021
1022 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("Invalid") );
1023
1024 // must correct field so that subsequent destructor doesn't whine!
1025 if constexpr (Args == compmatr || Args == diagmatr || Args == diagpower)
1026 std::get<0>(furtherArgs).numQubits = numTargs;
1027 }
1028
1029 SECTION( "matrix unsynced" ) {
1030
1031 if (!getQuESTEnv().isGpuAccelerated)
1032 return;
1033
1034 // only relevant to variable-size matrix functions
1035 if constexpr (Targs == any && (Args == compmatr || Args == diagmatr || Args == diagpower))
1036 *(std::get<0>(furtherArgs).wasGpuSynced) = 0;
1037 else
1038 return; // avoid empty test
1039
1040 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("sync") );
1041 }
1042
1043 SECTION( "targeted amps fit in node" ) {
1044
1045 // simplest to trigger validation using a statevector
1046 qureg = getCachedStatevecs().begin()->second;
1047
1048 // can only be validated when environment AND qureg
1049 // are distributed (over more than 1 node, of course)
1050 if (qureg.numNodes < 2)
1051 return;
1052
1053 // can only be validated if forced ctrl qubits permit
1054 // enough remaining targets
1055 int minNumCtrls = (Ctrls == one)? 1 : 0;
1056 int minNumTargs = numQubits - qureg.logNumNodes + 1;
1057 int maxNumTargs = numQubits - minNumCtrls;
1058 if (minNumTargs > maxNumTargs)
1059 return;
1060
1061 // only relevant to >=2-targ dense matrices, and further
1062 // only testable with any-targ variable-size matrices, since
1063 // 4x4 matrices might always be permissable by Qureg distribution
1064 if constexpr (Args == compmatr && Targs == any) {
1065
1066 // free existing matrix to avoid leak
1067 destroyCompMatr(std::get<0>(furtherArgs));
1068
1069 // try all illegally sized matrices
1070 int numNewTargs = GENERATE_COPY( range(minNumTargs, maxNumTargs+1) );
1071 targs = getRange(numNewTargs);
1072
1073 // ensure no overlap with ctrls; just get rid of them, EXCEPT when the API
1074 // function expects explicitly one ctrl which we must always supply
1075 ctrls = vector<int>(minNumCtrls, numQubits - 1); // {} or {last}
1076 states = {};
1077
1078 // create the new illegaly-sized matrix, which will be destroyed at test-case end
1079 CompMatr matr = getIdentityApiMatrix<Targs,Args>(numNewTargs);
1080 furtherArgs = tuple{ matr };
1081
1082 CAPTURE( minNumCtrls, numNewTargs, numQubits - minNumCtrls, ctrls, targs );
1083 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("cannot simultaneously store") && ContainsSubstring("remote amplitudes") );
1084
1085 } else {
1086 return; // avoid empty test
1087 }
1088 }
1089
1090 SECTION( "non-unitary exponent" ) {
1091
1092 // not relevant for functions which do not assert unitarity
1093 if (Apply != apply)
1094 return;
1095
1096 if constexpr (Args == diagpower)
1097 furtherArgs = tuple{ std::get<0>(furtherArgs), qcomp(1,1) };
1098 else
1099 return; // avoid empty test
1100
1101 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("exponent was not approximately real") );
1102 }
1103
1104 SECTION( "diverging exponent" ) {
1105
1106 // when being applied as a unitary, abs(elem)=1 so there's no
1107 // possibility of divergence (we'd merely trigger isUnitary)
1108 if (Apply == apply)
1109 return;
1110
1111 if constexpr (Args == diagpower)
1112 furtherArgs = tuple{ getZeroApiMatrix<Targs,Args>(numTargs), qcomp(-1,0) };
1113 else
1114 return; // avoid empty test
1115
1116 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("divergences") );
1117 }
1118
1119 SECTION( "zero axis rotation" ) {
1120
1121 if constexpr (Args == axisrots) {
1122 furtherArgs = tuple{ 0, 0, 0, 0 }; // (angle,x,y,z)
1123 } else
1124 return; // avoid empty test
1125
1126 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("zero vector") );
1127 }
1128
1129 SECTION( "qureg type" ) {
1130
1131 // only postMultiply*() functions discriminate Qureg
1132 if (Apply != rightapply)
1133 return;
1134
1135 // use any statevector
1136 qureg = getCachedStatevecs().begin()->second;
1137
1138 REQUIRE_THROWS_WITH( apiFunc(), ContainsSubstring("Expected a density matrix") );
1139 }
1140
1141 freeRemainingArgs<Targs,Args>(furtherArgs);
1142}
1143
1144
1145/*
1146 * fully test an API operation, on compatible
1147 * inputs as indicated by the template flags
1148 */
1149
1150template <NumQubitsFlag Ctrls, NumQubitsFlag Targs, ArgsFlag Args, ApplyFlag Apply>
1151void testOperation(auto operation, auto matrixRefGen) {
1152
1153 assertNumQubitsFlagsAreValid(Ctrls, Targs);
1154
1155 SECTION( LABEL_CORRECTNESS ) {
1156 testOperationCorrectness<Ctrls,Targs,Args,Apply>(operation, matrixRefGen);
1157 }
1158
1159 SECTION( LABEL_VALIDATION ) {
1160 testOperationValidation<Ctrls,Targs,Args,Apply>(operation);
1161 }
1162}
1163
1164
1165/*
1166 * perform unit tests for the four distinctly-controlled
1167 * variants of the given API operation (with the specified
1168 * function name suffix). 'numtargs' indicates the number
1169 * of target qubits accepted by the operation, and 'argtype'
1170 * indicates the types of remaining arguments (if any exist).
1171 * 'matrixgen' is the matrix representation (of varying
1172 * formats) of the operation, against which it will be compared.
1173 */
1174
1175// when every API operation had NO overloads, it was sufficient to
1176// call the below simple macro. Alas, since C++ overloads have been
1177// added, the below passing of a function (e.g. applyControlledHadamard)
1178// is ambigious, not specifying whether it is the int* or vector<int>
1179// version. We now have to explicit cast the function to its specific
1180// C-compatible version. Alas, to your imminent horror, this is.. erm...
1181
1182// #define TEST_ALL_CTRL_OPERATIONS( namesuffix, numtargs, argtype, matrixgen ) \
1183// TEST_CASE( "apply" #namesuffix, TEST_CATEGORY_OPS ) { testOperation<zero, numtargs,argtype>( apply ## namesuffix, matrixgen); } \
1184// TEST_CASE( "applyControlled" #namesuffix, TEST_CATEGORY_OPS ) { testOperation<one, numtargs,argtype>( applyControlled ## namesuffix, matrixgen); } \
1185// TEST_CASE( "applyMultiControlled" #namesuffix, TEST_CATEGORY_OPS ) { testOperation<any, numtargs,argtype>( applyMultiControlled ## namesuffix, matrixgen); } \
1186// TEST_CASE( "applyMultiStateControlled" #namesuffix, TEST_CATEGORY_OPS ) { testOperation<anystates,numtargs,argtype>( applyMultiStateControlled ## namesuffix, matrixgen); }
1187
1188
1189/*
1190 * define macros to pre-processor-time generate the API function
1191 * signatures to give to static_cast to disambiguate the function
1192 * from its C++ overloads, as per the above comment. The absolute
1193 * macro nightmare below results from not being able to propagate
1194 * templates to static_cast<...> which is not actually templated!
1195 * May God forgive me for the misdeeds commited here below.
1196 */
1197
1198// produces the control qubit arguments of a function signature,
1199// with a trailing comma (since ctrls always preceed more arguments)
1200#define GET_FUNC_CTRL_SUB_SIG( numctrls ) GET_FUNC_CTRL_SUB_SIG_##numctrls
1201#define GET_FUNC_CTRL_SUB_SIG_zero
1202#define GET_FUNC_CTRL_SUB_SIG_one int,
1203#define GET_FUNC_CTRL_SUB_SIG_any int*,int,
1204#define GET_FUNC_CTRL_SUB_SIG_anystates int*,int*,int,
1205
1206// produces the target qubit arguments of a function signature, complicated
1207// by paulistr and pauligad functions not ever passign explicit targ lists.
1208// trailing comma attached only when targs exist, and final args exist.
1209// beware, macros must never have spaces between the name and open-paranthesis!
1210#define GET_FUNC_TARG_SUB_SIG( numtargs, argtype ) GET_FUNC_TARG_SUB_SIG_##argtype( numtargs )
1211#define GET_FUNC_TARG_SUB_SIG_none( numtargs ) GET_FUNC_TARG_SUB_SIG_##numtargs
1212#define GET_FUNC_TARG_SUB_SIG_scalar( numtargs ) GET_FUNC_TARG_SUB_SIG_##numtargs ,
1213#define GET_FUNC_TARG_SUB_SIG_axisrots( numtargs ) GET_FUNC_TARG_SUB_SIG_##numtargs ,
1214#define GET_FUNC_TARG_SUB_SIG_compmatr( numtargs ) GET_FUNC_TARG_SUB_SIG_##numtargs ,
1215#define GET_FUNC_TARG_SUB_SIG_diagmatr( numtargs ) GET_FUNC_TARG_SUB_SIG_##numtargs ,
1216#define GET_FUNC_TARG_SUB_SIG_diagpower( numtargs ) GET_FUNC_TARG_SUB_SIG_##numtargs ,
1217#define GET_FUNC_TARG_SUB_SIG_paulistr( numtargs )
1218#define GET_FUNC_TARG_SUB_SIG_pauligad( numtargs )
1219#define GET_FUNC_TARG_SUB_SIG_one int
1220#define GET_FUNC_TARG_SUB_SIG_two int,int
1221#define GET_FUNC_TARG_SUB_SIG_any int*,int
1222
1223// produces the final arguments of a function signature (no trailing comma).
1224#define GET_FUNC_ARGS_SUB_SIG( numtargs, argtype ) GET_FUNC_ARGS_SUB_SIG_##numtargs##_##argtype
1225#define GET_FUNC_ARGS_SUB_SIG_one_none
1226#define GET_FUNC_ARGS_SUB_SIG_two_none
1227#define GET_FUNC_ARGS_SUB_SIG_any_none
1228#define GET_FUNC_ARGS_SUB_SIG_one_scalar qreal
1229#define GET_FUNC_ARGS_SUB_SIG_two_scalar qreal
1230#define GET_FUNC_ARGS_SUB_SIG_any_scalar qreal
1231#define GET_FUNC_ARGS_SUB_SIG_one_compmatr CompMatr1
1232#define GET_FUNC_ARGS_SUB_SIG_two_compmatr CompMatr2
1233#define GET_FUNC_ARGS_SUB_SIG_any_compmatr CompMatr
1234#define GET_FUNC_ARGS_SUB_SIG_one_diagmatr DiagMatr1
1235#define GET_FUNC_ARGS_SUB_SIG_two_diagmatr DiagMatr2
1236#define GET_FUNC_ARGS_SUB_SIG_any_diagmatr DiagMatr
1237#define GET_FUNC_ARGS_SUB_SIG_one_diagpower DiagMatr1,qcomp
1238#define GET_FUNC_ARGS_SUB_SIG_two_diagpower DiagMatr2,qcomp
1239#define GET_FUNC_ARGS_SUB_SIG_any_diagpower DiagMatr, qcomp
1240#define GET_FUNC_ARGS_SUB_SIG_one_axisrots qreal,qreal,qreal,qreal
1241#define GET_FUNC_ARGS_SUB_SIG_any_paulistr PauliStr
1242#define GET_FUNC_ARGS_SUB_SIG_any_pauligad PauliStr,qreal
1243
1244// produces the control-qubit-related prefix of a function name
1245#define GET_FUNC_NAME_PREFIX( numctrls ) GET_FUNC_NAME_PREFIX_##numctrls
1246#define GET_FUNC_NAME_PREFIX_zero apply
1247#define GET_FUNC_NAME_PREFIX_one applyControlled
1248#define GET_FUNC_NAME_PREFIX_any applyMultiControlled
1249#define GET_FUNC_NAME_PREFIX_anystates applyMultiStateControlled
1250
1251// produces a function name from the control qubits and the suffix, e.g. (any,T) -> applyMultiControlledT
1252#define GET_FUNC_NAME(numctrls, suffix) GET_FUNC_NAME_INNER(GET_FUNC_NAME_PREFIX(numctrls), suffix)
1253#define GET_FUNC_NAME_INNER(A, B) GET_FUNC_NAME_INNER_INNER(A, B)
1254#define GET_FUNC_NAME_INNER_INNER(A, B) A##B
1255
1256// converts the output of GET_FUNC_NAME() to a string, e.g. (any,T) -> "applyMultiControlledT"
1257#define GET_FUNC_NAME_STR(numctrls, suffix) GET_FUNC_NAME_STR_INNER( GET_FUNC_NAME(numctrls,suffix) )
1258#define GET_FUNC_NAME_STR_INNER(expr) GET_FUNC_NAME_STR_INNER_INNER(expr)
1259#define GET_FUNC_NAME_STR_INNER_INNER(symbol) #symbol
1260
1261// produces the signature of a function, e.g. (any,one,diagmatr) -> void(*)(Qureg, int*,int, int,DiagMatr1),
1262// which is the signature of applyMultiControlled(Qureg qureg, int* ctrls, int numCtrls, int targ, DiagMatr1 m);
1263// NOTE:
1264// THIS CURRENT EXCLUDES PAULISTR AND PAULIGAD argtype
1265#define GET_FUNC_SIG( numctrls, numtargs, argtype ) \
1266 void(*) ( \
1267 Qureg, \
1268 GET_FUNC_CTRL_SUB_SIG( numctrls ) \
1269 GET_FUNC_TARG_SUB_SIG( numtargs, argtype ) \
1270 GET_FUNC_ARGS_SUB_SIG( numtargs, argtype ) \
1271 )
1272
1273// produces a function name, casted to its explicit C-argument form, disambiguated from its C++ overload.
1274// e.g. (DiagMatrPower, zero, any, diagpower) -> static_cast<void(*)(Qureg,int*,int,DiagMatr,qcomp)>(applyDiagMatrPower)>
1275#define GET_CASTED_FUNC( namesuffix, numctrls, numtargs, argtype ) \
1276 static_cast< GET_FUNC_SIG(numctrls, numtargs, argtype) > ( \
1277 GET_FUNC_NAME(numctrls, namesuffix) )
1278
1279// defines a Catch2 test-case for the implied function
1280#define TEST_CASE_OPERATION( namesuffix, numctrls, numtargs, argtype, matrixgen ) \
1281 TEST_CASE( GET_FUNC_NAME_STR(numctrls, namesuffix), TEST_CATEGORY_OPS ) { \
1282 testOperation<numctrls, numtargs, argtype, apply>( \
1283 GET_CASTED_FUNC(namesuffix, numctrls, numtargs, argtype), \
1284 matrixgen); \
1285 }
1286
1287// automate the testing of a apply*() function for all its controlled variants
1288#define TEST_ALL_CTRL_OPERATIONS( namesuffix, numtargs, argtype, matrixgen ) \
1289 TEST_CASE_OPERATION( namesuffix, zero, numtargs, argtype, matrixgen ) \
1290 TEST_CASE_OPERATION( namesuffix, one, numtargs, argtype, matrixgen ) \
1291 TEST_CASE_OPERATION( namesuffix, any, numtargs, argtype, matrixgen ) \
1292 TEST_CASE_OPERATION( namesuffix, anystates, numtargs, argtype, matrixgen )
1293
1294
1295
1296/**
1297 * OPERATOR TESTS
1298 *
1299 * @ingroup unitops
1300 * @{
1301 */
1302
1303
1304/*
1305 * controlled operations
1306 */
1307
1308TEST_ALL_CTRL_OPERATIONS( PauliStr, any, paulistr, nullptr );
1309TEST_ALL_CTRL_OPERATIONS( PauliGadget, any, pauligad, nullptr );
1310TEST_ALL_CTRL_OPERATIONS( CompMatr1, one, compmatr, nullptr );
1311TEST_ALL_CTRL_OPERATIONS( CompMatr2, two, compmatr, nullptr );
1312TEST_ALL_CTRL_OPERATIONS( CompMatr, any, compmatr, nullptr );
1313TEST_ALL_CTRL_OPERATIONS( DiagMatr1, one, diagmatr, nullptr );
1314TEST_ALL_CTRL_OPERATIONS( DiagMatr2, two, diagmatr, nullptr );
1315TEST_ALL_CTRL_OPERATIONS( DiagMatr, any, diagmatr, nullptr );
1316TEST_ALL_CTRL_OPERATIONS( DiagMatrPower, any, diagpower, nullptr );
1317TEST_ALL_CTRL_OPERATIONS( Hadamard, one, none, FixedMatrices::H );
1318TEST_ALL_CTRL_OPERATIONS( PauliX, one, none, FixedMatrices::X );
1319TEST_ALL_CTRL_OPERATIONS( PauliY, one, none, FixedMatrices::Y );
1320TEST_ALL_CTRL_OPERATIONS( PauliZ, one, none, FixedMatrices::Z );
1321TEST_ALL_CTRL_OPERATIONS( T, one, none, FixedMatrices::T );
1322TEST_ALL_CTRL_OPERATIONS( S, one, none, FixedMatrices::S );
1323TEST_ALL_CTRL_OPERATIONS( Swap, two, none, FixedMatrices::SWAP );
1324TEST_ALL_CTRL_OPERATIONS( SqrtSwap, two, none, FixedMatrices::sqrtSWAP );
1325TEST_ALL_CTRL_OPERATIONS( RotateX, one, scalar, ParameterisedMatrices::Rx );
1326TEST_ALL_CTRL_OPERATIONS( RotateY, one, scalar, ParameterisedMatrices::Ry );
1327TEST_ALL_CTRL_OPERATIONS( RotateZ, one, scalar, ParameterisedMatrices::Rz );
1328TEST_ALL_CTRL_OPERATIONS( RotateAroundAxis, one, axisrots, nullptr );
1329TEST_ALL_CTRL_OPERATIONS( MultiQubitNot, any, none, VariableSizeMatrices::X );
1330TEST_ALL_CTRL_OPERATIONS( PhaseGadget, any, scalar, VariableSizeParameterisedMatrices::Z );
1331
1332
1333/*
1334 * non-controlled operations with no C++ overloads
1335 */
1336
1337TEST_CASE( "applyPhaseFlip", TEST_CATEGORY_OPS ) { testOperation<zero,one,none,apply> (applyPhaseFlip, VariableSizeMatrices::PF(1)); }
1338TEST_CASE( "applyTwoQubitPhaseFlip", TEST_CATEGORY_OPS ) { testOperation<zero,two,none,apply> (applyTwoQubitPhaseFlip, VariableSizeMatrices::PF(2)); }
1339TEST_CASE( "applyPhaseShift", TEST_CATEGORY_OPS ) { testOperation<zero,one,scalar,apply>(applyPhaseShift, ParameterisedMatrices::PS ); }
1340TEST_CASE( "applyTwoQubitPhaseShift", TEST_CATEGORY_OPS ) { testOperation<zero,two,scalar,apply>(applyTwoQubitPhaseShift, ParameterisedMatrices::PS2 ); }
1341
1342
1343/*
1344 * non-controlled operations which have a C++ overload
1345 * (because they accept qubit lists which become vector),
1346 * and so which require explicit casting to resolve the
1347 * compiler ambiguity (spaghetti 4 lyf)
1348 */
1349
1350TEST_CASE( "applyMultiQubitPhaseFlip", TEST_CATEGORY_OPS ) {
1351 auto func = static_cast<void(*)(Qureg, int*, int)>(applyMultiQubitPhaseFlip);
1352 testOperation<zero,any,none,apply>(func, VariableSizeMatrices::PF);
1353}
1354
1355TEST_CASE( "applyMultiQubitPhaseShift", TEST_CATEGORY_OPS ) {
1356 auto func = static_cast<void(*)(Qureg, int*, int, qreal)>(applyMultiQubitPhaseShift);
1357 testOperation<zero,any,scalar,apply>(func, VariableSizeParameterisedMatrices::PS);
1358}
1359
1360
1361/*
1362 * operations which need custom logic
1363 */
1364
1365TEST_CASE( "applyQuantumFourierTransform", TEST_CATEGORY_OPS ) {
1366
1367 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1368
1369 SECTION( LABEL_CORRECTNESS ) {
1370
1371 int numTargs = GENERATE_COPY( range(1,numQubits+1) );
1372 auto targs = GENERATE_TARGS( numQubits, numTargs );
1373
1374 CAPTURE( targs );
1375
1376 SECTION( LABEL_STATEVEC ) {
1377
1378 auto testFunc = [&](Qureg qureg, qvector& ref) {
1379 applyQuantumFourierTransform(qureg, targs.data(), targs.size());
1380 ref = getDisceteFourierTransform(ref, targs);
1381 };
1382
1383 TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc);
1384 }
1385
1386 SECTION( LABEL_DENSMATR ) {
1387
1388 // prepare a random mixture
1389 auto states = getRandomOrthonormalStateVectors(numQubits, getRandomInt(1,10));
1390 auto probs = getRandomProbabilities(states.size());
1391
1392 auto testFunc = [&](Qureg qureg, qmatrix& ref) {
1393
1394 // overwrite the Qureg debug state set by caller to above mixture
1395 setQuregToReference(qureg, getMixture(states, probs));
1396 applyQuantumFourierTransform(qureg, targs.data(), targs.size());
1397
1398 ref = getZeroMatrix(ref.size());
1399 for (size_t i=0; i<states.size(); i++) {
1400 qvector vec = getDisceteFourierTransform(states[i], targs);
1401 ref += probs[i] * getOuterProduct(vec, vec);
1402 }
1403 };
1404
1405 TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc);
1406 }
1407 }
1408
1409 /// @todo input validation
1410}
1411
1412
1413TEST_CASE( "applyFullQuantumFourierTransform", TEST_CATEGORY_OPS ) {
1414
1415 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1416
1417 SECTION( LABEL_CORRECTNESS ) {
1418
1419 GENERATE( range(0,10) );
1420
1421 SECTION( LABEL_STATEVEC ) {
1422
1423 auto testFunc = [&](Qureg qureg, qvector& ref) {
1425 ref = getDisceteFourierTransform(ref);
1426 };
1427
1428 TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc);
1429 }
1430
1431 SECTION( LABEL_DENSMATR ) {
1432
1433 // prepare a random mixture
1434 auto states = getRandomOrthonormalStateVectors(numQubits, getRandomInt(1,10));
1435 auto probs = getRandomProbabilities(states.size());
1436
1437 auto testFunc = [&](Qureg qureg, qmatrix& ref) {
1438
1439 // overwrite the Qureg debug state set by caller to above mixture
1440 setQuregToReference(qureg, getMixture(states, probs));
1442
1443 ref = getZeroMatrix(ref.size());
1444 for (size_t i=0; i<states.size(); i++) {
1445 qvector vec = getDisceteFourierTransform(states[i]);
1446 ref += probs[i] * getOuterProduct(vec, vec);
1447 }
1448 };
1449
1450 TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc);
1451 }
1452 }
1453
1454 /// @todo input validation
1455}
1456
1457
1458TEST_CASE( "applyQubitProjector", TEST_CATEGORY_OPS ) {
1459
1460 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1461
1462 SECTION( LABEL_CORRECTNESS ) {
1463
1464 GENERATE( range(0,10) );
1465 int target = GENERATE_COPY( range(0,numQubits) );
1466 int outcome = GENERATE( 0, 1 );
1467
1468 qmatrix projector = getProjector(outcome);
1469
1470 auto testFunc = [&](Qureg qureg, auto& ref) {
1471 applyQubitProjector(qureg, target, outcome);
1472 applyReferenceOperator(ref, {target}, projector);
1473 };
1474
1475 CAPTURE( target, outcome );
1476 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1477 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1478 }
1479
1480 /// @todo input validation
1481}
1482
1483
1484TEST_CASE( "applyMultiQubitProjector", TEST_CATEGORY_OPS ) {
1485
1486 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1487
1488 SECTION( LABEL_CORRECTNESS ) {
1489
1490 int numTargs = GENERATE_COPY( range(1,numQubits+1) );
1491 auto targets = GENERATE_TARGS( numQubits, numTargs );
1492 auto outcomes = getRandomOutcomes(numTargs);
1493
1494 qmatrix projector = getProjector(targets, outcomes, numQubits);
1495
1496 auto testFunc = [&](Qureg qureg, auto& ref) {
1497 applyMultiQubitProjector(qureg, targets.data(), outcomes.data(), numTargs);
1498 applyReferenceOperator(ref, projector);
1499 };
1500
1501 CAPTURE( targets, outcomes );
1502 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1503 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1504 }
1505
1506 /// @todo input validation
1507}
1508
1509
1510TEST_CASE( "applyForcedQubitMeasurement", TEST_CATEGORY_OPS ) {
1511
1512 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1513
1514 SECTION( LABEL_CORRECTNESS ) {
1515
1516 GENERATE( range(0,10) );
1517 int target = GENERATE_COPY( range(0,numQubits) );
1518 int outcome = GENERATE( 0, 1 );
1519
1520 qmatrix projector = getProjector(outcome);
1521
1522 auto testFunc = [&](Qureg qureg, auto& ref) {
1523
1524 // overwrite caller's setting of initDebugState, since
1525 // that precludes outcomes=|0><0| due to zero-probability
1526 setToRandomState(ref);
1527 setQuregToReference(qureg, ref);
1528
1529 // compare the probabilities...
1530 qreal apiProb = applyForcedQubitMeasurement(qureg, target, outcome);
1531 qreal refProb = getReferenceProbability(ref, {target}, {outcome});
1532 REQUIRE_AGREE( apiProb, refProb );
1533
1534 // and the post-projection states (caller calls subsequent REQUIRE_AGREE)
1535 applyReferenceOperator(ref, {target}, projector);
1536 ref /= (qureg.isDensityMatrix)?
1537 refProb : std::sqrt(refProb);
1538 };
1539
1540 CAPTURE( target, outcome );
1541 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1542 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1543 }
1544
1545 /// @todo input validation
1546}
1547
1548
1549TEST_CASE( "applyForcedMultiQubitMeasurement", TEST_CATEGORY_OPS ) {
1550
1551 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1552
1553 SECTION( LABEL_CORRECTNESS ) {
1554
1555 int numTargs = GENERATE_COPY( range(1,numQubits+1) );
1556 auto targets = GENERATE_TARGS( numQubits, numTargs );
1557 auto outcomes = getRandomOutcomes(numTargs);
1558
1559 qmatrix projector = getProjector(targets, outcomes, numQubits);
1560
1561 // this test may randomly request a measurement outcome which
1562 // is illegally unlikely, triggering validation; we merely
1563 // disable such validation and hope divergences don't break the test!
1565
1566 auto testFunc = [&](Qureg qureg, auto& ref) {
1567
1568 // overwrite caller's setting of initDebugState, since
1569 // that precludes outcomes=|0><0| due to zero-probability
1570 setToRandomState(ref);
1571 setQuregToReference(qureg, ref);
1572
1573 // compare the probabilities...
1574 qreal apiProb = applyForcedMultiQubitMeasurement(qureg, targets.data(), outcomes.data(), numTargs);
1575 qreal refProb = getReferenceProbability(ref, targets, outcomes);
1576 REQUIRE_AGREE( apiProb, refProb );
1577
1578 // and the post-measurement states (caller calls subsequent REQUIRE_AGREE)
1579 applyReferenceOperator(ref, projector);
1580 ref /= (qureg.isDensityMatrix)?
1581 refProb : std::sqrt(refProb);
1582 };
1583
1584 CAPTURE( targets, outcomes );
1585 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1586 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1587
1589 }
1590
1591 /// @todo input validation
1592}
1593
1594
1595TEST_CASE( "applyMultiQubitMeasurement", TEST_CATEGORY_OPS ) {
1596
1597 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1598
1599 SECTION( LABEL_CORRECTNESS ) {
1600
1601 int numTargs = GENERATE_COPY( range(1,numQubits+1) );
1602 auto targets = GENERATE_TARGS( numQubits, numTargs );
1603
1604 auto testFunc = [&](Qureg qureg, auto& ref) {
1605
1606 // overwrite caller's setting of initDebugState, since
1607 // sampling requires the outcome probs are normalised
1608 setToRandomState(ref);
1609 setQuregToReference(qureg, ref);
1610
1611 // the output API state...
1612 qindex apiOut = applyMultiQubitMeasurement(qureg, targets.data(), numTargs);
1613
1614 // informs the projector which determines the post-measurement reference
1615 auto apiOutBits = getBits(apiOut, numTargs);
1616 qmatrix projector = getProjector(targets, apiOutBits, numQubits);
1617 applyReferenceOperator(ref, projector);
1618 qreal refProb = getReferenceProbability(ref, targets, apiOutBits);
1619 ref /= (qureg.isDensityMatrix)?
1620 refProb : std::sqrt(refProb);
1621 };
1622
1623 CAPTURE( targets );
1624 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1625 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1626 }
1627
1628 /// @todo input validation
1629}
1630
1631
1632TEST_CASE( "applyMultiQubitMeasurementAndGetProb", TEST_CATEGORY_OPS ) {
1633
1634 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1635
1636 SECTION( LABEL_CORRECTNESS ) {
1637
1638 int numTargs = GENERATE_COPY( range(1,numQubits+1) );
1639 auto targets = GENERATE_TARGS( numQubits, numTargs );
1640
1641 auto testFunc = [&](Qureg qureg, auto& ref) {
1642
1643 // overwrite caller's setting of initDebugState, since
1644 // sampling requires the outcome probs are normalised
1645 setToRandomState(ref);
1646 setQuregToReference(qureg, ref);
1647
1648 // compare the measurement probability...
1649 qreal apiProb = -1;
1650 qindex apiOut = applyMultiQubitMeasurementAndGetProb(qureg, targets.data(), numTargs, &apiProb);
1651 auto apiOutBits = getBits(apiOut, numTargs);
1652 qreal refProb = getReferenceProbability(ref, targets, apiOutBits);
1653 REQUIRE_AGREE( apiProb, refProb );
1654
1655 // and the post-measurement states (caller calls subsequent REQUIRE_AGREE)
1656 qmatrix projector = getProjector(targets, apiOutBits, numQubits);
1657 applyReferenceOperator(ref, projector);
1658 ref /= (qureg.isDensityMatrix)?
1659 refProb : std::sqrt(refProb);
1660 };
1661
1662 CAPTURE( targets );
1663 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1664 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1665 }
1666
1667 /// @todo input validation
1668}
1669
1670
1671TEST_CASE( "applyQubitMeasurement", TEST_CATEGORY_OPS ) {
1672
1673 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1674
1675 SECTION( LABEL_CORRECTNESS ) {
1676
1677 GENERATE( range(0,10) );
1678 int target = GENERATE_COPY( range(0,numQubits) );
1679
1680 auto testFunc = [&](Qureg qureg, auto& ref) {
1681
1682 // overwrite caller's setting of initDebugState, since
1683 // sampling requires the outcome probs are normalised
1684 setToRandomState(ref);
1685 setQuregToReference(qureg, ref);
1686
1687 // the output API state...
1688 int apiOut = applyQubitMeasurement(qureg, target);
1689
1690 // informs the projector which determines the post-measurement reference
1691 qmatrix projector = getProjector(apiOut);
1692 applyReferenceOperator(ref, {target}, projector);
1693 qreal refProb = getReferenceProbability(ref, {target}, {apiOut});
1694 ref /= (qureg.isDensityMatrix)?
1695 refProb : std::sqrt(refProb);
1696 };
1697
1698 CAPTURE( target );
1699 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1700 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1701 }
1702
1703 /// @todo input validation
1704}
1705
1706
1707TEST_CASE( "applyQubitMeasurementAndGetProb", TEST_CATEGORY_OPS ) {
1708
1709 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1710
1711 SECTION( LABEL_CORRECTNESS ) {
1712
1713 GENERATE( range(0,10) );
1714 int target = GENERATE_COPY( range(0,numQubits) );
1715
1716 auto testFunc = [&](Qureg qureg, auto& ref) {
1717
1718 // overwrite caller's setting of initDebugState, since
1719 // sampling requires the outcome probs are normalised
1720 setToRandomState(ref);
1721 setQuregToReference(qureg, ref);
1722
1723 // compare the measurement probability...
1724 qreal apiProb = -1;
1725 int apiOut = applyQubitMeasurementAndGetProb(qureg, target, &apiProb);
1726 qreal refProb = getReferenceProbability(ref, {target}, {apiOut});
1727 REQUIRE_AGREE( apiProb, refProb );
1728
1729 // and the post-measurement states (caller calls subsequent REQUIRE_AGREE)
1730 qmatrix projector = getProjector(apiOut);
1731 applyReferenceOperator(ref, {target}, projector);
1732 ref /= (qureg.isDensityMatrix)?
1733 refProb : std::sqrt(refProb);
1734 };
1735
1736 CAPTURE( target );
1737 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1738 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1739 }
1740
1741 /// @todo input validation
1742}
1743
1744
1745TEST_CASE( "applyFullStateDiagMatr", TEST_CATEGORY_OPS LABEL_MIXED_DEPLOY_TAG ) {
1746
1747 PREPARE_TEST( numQubits, cachedSV, cachedDM, refSV, refDM );
1748
1749 auto cachedMatrs = getCachedFullStateDiagMatrs();
1750
1751 SECTION( LABEL_CORRECTNESS ) {
1752
1753 qmatrix refMatr = getRandomDiagonalUnitary(numQubits);
1754 auto apiFunc = applyFullStateDiagMatr;
1755
1756 GENERATE( range(0, getNumTestedMixedDeploymentRepetitions()) );
1757
1758 SECTION( LABEL_STATEVEC ) {
1759
1760 auto refFunc = [&] (qvector& state, qmatrix matr) { applyReferenceOperator(state, matr); };
1761
1762 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedSV, cachedMatrs, apiFunc, refSV, refMatr, refFunc);
1763 }
1764
1765 SECTION( LABEL_DENSMATR ) {
1766
1767 auto refFunc = [&] (qmatrix& state, qmatrix matr) { applyReferenceOperator(state, matr); };
1768
1769 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedDM, cachedMatrs, apiFunc, refDM, refMatr, refFunc);
1770 }
1771 }
1772
1773 /// @todo input validation
1774}
1775
1776
1777TEST_CASE( "applyFullStateDiagMatrPower", TEST_CATEGORY_OPS LABEL_MIXED_DEPLOY_TAG ) {
1778
1779 PREPARE_TEST( numQubits, cachedSV, cachedDM, refSV, refDM );
1780
1781 auto cachedMatrs = getCachedFullStateDiagMatrs();
1782
1783 SECTION( LABEL_CORRECTNESS ) {
1784
1785 qmatrix refMatr = getRandomDiagonalUnitary(numQubits);
1786
1787 // supplying a complex exponent requires disabling
1788 // numerical validation to relax unitarity
1789 bool testRealExp = GENERATE( true, false );
1790 qcomp exponent = (testRealExp)?
1791 qcomp(getRandomReal(-2, 2), 0):
1793
1794 auto apiFunc = [&](Qureg qureg, FullStateDiagMatr matr) {
1795 return applyFullStateDiagMatrPower(qureg, matr, exponent);
1796 };
1797
1798 CAPTURE( exponent );
1799
1800 GENERATE( range(0, getNumTestedMixedDeploymentRepetitions()) );
1801
1802 if (!testRealExp)
1804
1805 SECTION( LABEL_STATEVEC ) {
1806
1807 auto refFunc = [&] (qvector& state, qmatrix matr) {
1808 matr = getPowerOfDiagonalMatrix(matr, exponent);
1809 applyReferenceOperator(state, matr);
1810 };
1811
1812 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedSV, cachedMatrs, apiFunc, refSV, refMatr, refFunc);
1813 }
1814
1815 SECTION( LABEL_DENSMATR ) {
1816
1817 auto refFunc = [&] (qmatrix& state, qmatrix matr) {
1818 matr = getPowerOfDiagonalMatrix(matr, exponent);
1819 applyReferenceOperator(state, matr);
1820 };
1821
1822 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedDM, cachedMatrs, apiFunc, refDM, refMatr, refFunc);
1823 }
1824
1826 }
1827
1828 /// @todo input validation
1829}
1830
1831
1832TEST_CASE( "applyNonUnitaryPauliGadget", TEST_CATEGORY_OPS ) {
1833
1834 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
1835
1836 SECTION( LABEL_CORRECTNESS ) {
1837
1838 // prepare a random Pauli string and angle
1839 int numTargs = GENERATE_COPY( range(1, numQubits+1) );
1840 auto targs = GENERATE_TARGS( numQubits, numTargs );
1841 PauliStr str = getRandomPauliStr(targs);
1842 qcomp angle = getRandomComplex();
1843
1844 // prepare the corresponding reference matrix exp(-i angle pauli)
1845 auto matrRef = getExponentialOfPauliMatrix(angle, getMatrix(str, numQubits));
1846
1847 auto testFunc = [&](Qureg qureg, auto& stateRef) {
1848 applyNonUnitaryPauliGadget(qureg, str, angle);
1849 applyReferenceOperator(stateRef, matrRef);
1850 };
1851
1852 CAPTURE( targs, angle );
1853 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
1854 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
1855 }
1856
1857 /// @todo input validation
1858}
1859
1860
1861/** @} (end defgroup) */
1862
1863
1864
1865/**
1866 * OPERATOR TESTS
1867 *
1868 * @ingroup unitmult
1869 * @{
1870 */
1871
1872
1873TEST_CASE( "leftapplySwap", TEST_CATEGORY_MULT ) { testOperation<zero,two,none,leftapply>(leftapplySwap, FixedMatrices::SWAP); }
1874TEST_CASE( "leftapplyPauliX", TEST_CATEGORY_MULT ) { testOperation<zero,one,none,leftapply>(leftapplyPauliX, FixedMatrices::X); }
1875TEST_CASE( "leftapplyPauliY", TEST_CATEGORY_MULT ) { testOperation<zero,one,none,leftapply>(leftapplyPauliY, FixedMatrices::Y); }
1876TEST_CASE( "leftapplyPauliZ", TEST_CATEGORY_MULT ) { testOperation<zero,one,none,leftapply>(leftapplyPauliZ, FixedMatrices::Z); }
1877TEST_CASE( "leftapplyPauliStr", TEST_CATEGORY_MULT ) { testOperation<zero,any,paulistr,leftapply>(leftapplyPauliStr, nullptr); }
1878TEST_CASE( "leftapplyPauliGadget", TEST_CATEGORY_MULT ) { testOperation<zero,any,pauligad,leftapply>(leftapplyPauliGadget, nullptr); }
1879TEST_CASE( "leftapplyCompMatr1", TEST_CATEGORY_MULT ) { testOperation<zero,one,compmatr,leftapply>(leftapplyCompMatr1, nullptr); }
1880TEST_CASE( "leftapplyCompMatr2", TEST_CATEGORY_MULT ) { testOperation<zero,two,compmatr,leftapply>(leftapplyCompMatr2, nullptr); }
1881TEST_CASE( "leftapplyDiagMatr1", TEST_CATEGORY_MULT ) { testOperation<zero,one,diagmatr,leftapply>(leftapplyDiagMatr1, nullptr); }
1882TEST_CASE( "leftapplyDiagMatr2", TEST_CATEGORY_MULT ) { testOperation<zero,two,diagmatr,leftapply>(leftapplyDiagMatr2, nullptr); }
1883
1884TEST_CASE( "rightapplySwap", TEST_CATEGORY_MULT ) { testOperation<zero,two,none,rightapply>(rightapplySwap, FixedMatrices::SWAP); }
1885TEST_CASE( "rightapplyPauliX", TEST_CATEGORY_MULT ) { testOperation<zero,one,none,rightapply>(rightapplyPauliX, FixedMatrices::X); }
1886TEST_CASE( "rightapplyPauliY", TEST_CATEGORY_MULT ) { testOperation<zero,one,none,rightapply>(rightapplyPauliY, FixedMatrices::Y); }
1887TEST_CASE( "rightapplyPauliZ", TEST_CATEGORY_MULT ) { testOperation<zero,one,none,rightapply>(rightapplyPauliZ, FixedMatrices::Z); }
1888TEST_CASE( "rightapplyPauliStr", TEST_CATEGORY_MULT ) { testOperation<zero,any,paulistr,rightapply>(rightapplyPauliStr, nullptr); }
1889TEST_CASE( "rightapplyPauliGadget", TEST_CATEGORY_MULT ) { testOperation<zero,any,pauligad,rightapply>(rightapplyPauliGadget, nullptr); }
1890TEST_CASE( "rightapplyCompMatr1", TEST_CATEGORY_MULT ) { testOperation<zero,one,compmatr,rightapply>(rightapplyCompMatr1, nullptr); }
1891TEST_CASE( "rightapplyCompMatr2", TEST_CATEGORY_MULT ) { testOperation<zero,two,compmatr,rightapply>(rightapplyCompMatr2, nullptr); }
1892TEST_CASE( "rightapplyDiagMatr1", TEST_CATEGORY_MULT ) { testOperation<zero,one,diagmatr,rightapply>(rightapplyDiagMatr1, nullptr); }
1893TEST_CASE( "rightapplyDiagMatr2", TEST_CATEGORY_MULT ) { testOperation<zero,two,diagmatr,rightapply>(rightapplyDiagMatr2, nullptr); }
1894
1895
1896/*
1897 * C++ overloads which accept qubit lists as vectors
1898 * and so which require explicit casting to resolve the
1899 * compiler ambiguity (spaghetti 4 lyf)
1900 */
1901
1902
1903TEST_CASE( "leftapplyCompMatr", TEST_CATEGORY_MULT ) {
1904 auto func = static_cast<void(*)(Qureg, int*, int, CompMatr)>(leftapplyCompMatr);
1905 testOperation<zero,any,compmatr,leftapply>(func, nullptr);
1906}
1907
1908TEST_CASE( "leftapplyDiagMatr", TEST_CATEGORY_MULT ) {
1909 auto func = static_cast<void(*)(Qureg, int*, int, DiagMatr)>(leftapplyDiagMatr);
1910 testOperation<zero,any,diagmatr,leftapply>(func, nullptr);
1911}
1912
1913TEST_CASE( "leftapplyDiagMatrPower", TEST_CATEGORY_MULT ) {
1914 auto func = static_cast<void(*)(Qureg, int*, int, DiagMatr, qcomp)>(leftapplyDiagMatrPower);
1915 testOperation<zero,any,diagpower,leftapply>(func, nullptr);
1916}
1917
1918TEST_CASE( "leftapplyMultiQubitNot", TEST_CATEGORY_MULT ) {
1919 auto func = static_cast<void(*)(Qureg, int*, int)>(leftapplyMultiQubitNot);
1920 testOperation<zero,any,none,leftapply>(func, VariableSizeMatrices::X);
1921}
1922
1923TEST_CASE( "leftapplyPhaseGadget", TEST_CATEGORY_MULT ) {
1924 auto func = static_cast<void(*)(Qureg, int*, int, qreal)>(leftapplyPhaseGadget);
1925 testOperation<zero,any,scalar,leftapply>(func, VariableSizeParameterisedMatrices::Z);
1926}
1927
1928
1929TEST_CASE( "rightapplyCompMatr", TEST_CATEGORY_MULT ) {
1930 auto func = static_cast<void(*)(Qureg, int*, int, CompMatr)>(rightapplyCompMatr);
1931 testOperation<zero,any,compmatr,rightapply>(func, nullptr);
1932}
1933
1934TEST_CASE( "rightapplyDiagMatr", TEST_CATEGORY_MULT ) {
1935 auto func = static_cast<void(*)(Qureg, int*, int, DiagMatr)>(rightapplyDiagMatr);
1936 testOperation<zero,any,diagmatr,rightapply>(func, nullptr);
1937}
1938
1939TEST_CASE( "rightapplyDiagMatrPower", TEST_CATEGORY_MULT ) {
1940 auto func = static_cast<void(*)(Qureg, int*, int, DiagMatr, qcomp)>(rightapplyDiagMatrPower);
1941 testOperation<zero,any,diagpower,rightapply>(func, nullptr);
1942}
1943
1944TEST_CASE( "rightapplyMultiQubitNot", TEST_CATEGORY_MULT ) {
1945 auto func = static_cast<void(*)(Qureg, int*, int)>(rightapplyMultiQubitNot);
1946 testOperation<zero,any,none,rightapply>(func, VariableSizeMatrices::X);
1947}
1948
1949TEST_CASE( "rightapplyPhaseGadget", TEST_CATEGORY_MULT ) {
1950 auto func = static_cast<void(*)(Qureg, int*, int, qreal)>(rightapplyPhaseGadget);
1951 testOperation<zero,any,scalar,rightapply>(func, VariableSizeParameterisedMatrices::Z);
1952}
1953
1954
1955/*
1956 * operations which need custom logic
1957 */
1958
1959
1960TEST_CASE( "leftapplyFullStateDiagMatr", TEST_CATEGORY_MULT LABEL_MIXED_DEPLOY_TAG ) {
1961
1962 PREPARE_TEST( numQubits, cachedSV, cachedDM, refSV, refDM );
1963
1964 auto cachedMatrs = getCachedFullStateDiagMatrs();
1965
1966 SECTION( LABEL_CORRECTNESS ) {
1967
1968 qmatrix refMatr = getRandomDiagonalMatrix(getPow2(numQubits));
1969 auto apiFunc = leftapplyFullStateDiagMatr;
1970
1971 GENERATE( range(0, getNumTestedMixedDeploymentRepetitions()) );
1972
1973 SECTION( LABEL_STATEVEC ) {
1974
1975 auto refFunc = [&] (qvector& state, qmatrix matr) { leftapplyReferenceOperator(state, matr); };
1976
1977 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedSV, cachedMatrs, apiFunc, refSV, refMatr, refFunc);
1978 }
1979
1980 SECTION( LABEL_DENSMATR ) {
1981
1982 auto refFunc = [&] (qmatrix& state, qmatrix matr) { leftapplyReferenceOperator(state, matr); };
1983
1984 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedDM, cachedMatrs, apiFunc, refDM, refMatr, refFunc);
1985 }
1986 }
1987
1988 /// @todo input validation
1989}
1990
1991
1992TEST_CASE( "rightapplyFullStateDiagMatr", TEST_CATEGORY_MULT LABEL_MIXED_DEPLOY_TAG ) {
1993
1994 PREPARE_TEST( numQubits, cachedSV, cachedDM, refSV, refDM );
1995
1996 auto cachedMatrs = getCachedFullStateDiagMatrs();
1997
1998 SECTION( LABEL_CORRECTNESS ) {
1999
2000 qmatrix refMatr = getRandomDiagonalMatrix(getPow2(numQubits));
2001 auto apiFunc = rightapplyFullStateDiagMatr;
2002
2003 GENERATE( range(0, getNumTestedMixedDeploymentRepetitions()) );
2004
2005 SECTION( LABEL_DENSMATR ) {
2006
2007 auto refFunc = [&] (qmatrix& state, qmatrix matr) { rightapplyReferenceOperator(state, matr); };
2008
2009 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedDM, cachedMatrs, apiFunc, refDM, refMatr, refFunc);
2010 }
2011 }
2012
2013 /// @todo input validation
2014}
2015
2016
2017TEST_CASE( "leftapplyFullStateDiagMatrPower", TEST_CATEGORY_MULT LABEL_MIXED_DEPLOY_TAG ) {
2018
2019 PREPARE_TEST( numQubits, cachedSV, cachedDM, refSV, refDM );
2020
2021 auto cachedMatrs = getCachedFullStateDiagMatrs();
2022
2023 SECTION( LABEL_CORRECTNESS ) {
2024
2025 qmatrix refMatr = getRandomDiagonalMatrix(getPow2(numQubits));
2026 qcomp exponent = getRandomComplex();
2027
2028 auto apiFunc = [&](Qureg qureg, FullStateDiagMatr matr) {
2029 return leftapplyFullStateDiagMatrPower(qureg, matr, exponent);
2030 };
2031
2032 CAPTURE( exponent );
2033
2034 GENERATE( range(0, getNumTestedMixedDeploymentRepetitions()) );
2035
2036 SECTION( LABEL_STATEVEC ) {
2037
2038 auto refFunc = [&] (qvector& state, qmatrix matr) {
2039 matr = getPowerOfDiagonalMatrix(matr, exponent);
2040 leftapplyReferenceOperator(state, matr);
2041 };
2042
2043 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedSV, cachedMatrs, apiFunc, refSV, refMatr, refFunc);
2044 }
2045
2046 SECTION( LABEL_DENSMATR ) {
2047
2048 auto refFunc = [&] (qmatrix& state, qmatrix matr) {
2049 matr = getPowerOfDiagonalMatrix(matr, exponent);
2050 leftapplyReferenceOperator(state, matr);
2051 };
2052
2053 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedDM, cachedMatrs, apiFunc, refDM, refMatr, refFunc);
2054 }
2055 }
2056
2057 /// @todo input validation
2058}
2059
2060
2061TEST_CASE( "rightapplyFullStateDiagMatrPower", TEST_CATEGORY_MULT LABEL_MIXED_DEPLOY_TAG ) {
2062
2063 PREPARE_TEST( numQubits, cachedSV, cachedDM, refSV, refDM );
2064
2065 auto cachedMatrs = getCachedFullStateDiagMatrs();
2066
2067 SECTION( LABEL_CORRECTNESS ) {
2068
2069 qmatrix refMatr = getRandomDiagonalMatrix(getPow2(numQubits));
2070 qcomp exponent = getRandomComplex();
2071
2072 auto apiFunc = [&](Qureg qureg, FullStateDiagMatr matr) {
2073 return rightapplyFullStateDiagMatrPower(qureg, matr, exponent);
2074 };
2075
2076 CAPTURE( exponent );
2077
2078 GENERATE( range(0, getNumTestedMixedDeploymentRepetitions()) );
2079
2080 SECTION( LABEL_DENSMATR ) {
2081
2082 auto refFunc = [&] (qmatrix& state, qmatrix matr) {
2083 matr = getPowerOfDiagonalMatrix(matr, exponent);
2084 rightapplyReferenceOperator(state, matr);
2085 };
2086
2087 TEST_ON_CACHED_QUREG_AND_MATRIX( cachedDM, cachedMatrs, apiFunc, refDM, refMatr, refFunc);
2088 }
2089 }
2090
2091 /// @todo input validation
2092}
2093
2094
2095TEST_CASE( "leftapplyQubitProjector", TEST_CATEGORY_OPS ) {
2096
2097 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
2098
2099 SECTION( LABEL_CORRECTNESS ) {
2100
2101 GENERATE( range(0,10) );
2102 int target = GENERATE_COPY( range(0,numQubits) );
2103 int outcome = GENERATE( 0, 1 );
2104
2105 qmatrix projector = getProjector(outcome);
2106
2107 auto testFunc = [&](Qureg qureg, auto& ref) {
2108 leftapplyQubitProjector(qureg, target, outcome);
2109 leftapplyReferenceOperator(ref, {target}, projector);
2110 };
2111
2112 CAPTURE( target, outcome );
2113 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
2114 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
2115 }
2116
2117 /// @todo input validation
2118}
2119
2120
2121TEST_CASE( "rightapplyQubitProjector", TEST_CATEGORY_OPS ) {
2122
2123 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
2124
2125 SECTION( LABEL_CORRECTNESS ) {
2126
2127 GENERATE( range(0,10) );
2128 int target = GENERATE_COPY( range(0,numQubits) );
2129 int outcome = GENERATE( 0, 1 );
2130
2131 qmatrix projector = getProjector(outcome);
2132
2133 auto testFunc = [&](Qureg qureg, auto& ref) {
2134 rightapplyQubitProjector(qureg, target, outcome);
2135 rightapplyReferenceOperator(ref, {target}, projector);
2136 };
2137
2138 CAPTURE( target, outcome );
2139 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
2140 }
2141
2142 /// @todo input validation
2143}
2144
2145
2146TEST_CASE( "leftapplyMultiQubitProjector", TEST_CATEGORY_OPS ) {
2147
2148 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
2149
2150 SECTION( LABEL_CORRECTNESS ) {
2151
2152 int numTargs = GENERATE_COPY( range(1,numQubits+1) );
2153 auto targets = GENERATE_TARGS( numQubits, numTargs );
2154 auto outcomes = getRandomOutcomes(numTargs);
2155
2156 qmatrix projector = getProjector(targets, outcomes, numQubits);
2157
2158 auto testFunc = [&](Qureg qureg, auto& ref) {
2159 leftapplyMultiQubitProjector(qureg, targets.data(), outcomes.data(), numTargs);
2160 leftapplyReferenceOperator(ref, projector);
2161 };
2162
2163 CAPTURE( targets, outcomes );
2164 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
2165 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
2166 }
2167
2168 /// @todo input validation
2169}
2170
2171
2172TEST_CASE( "rightapplyMultiQubitProjector", TEST_CATEGORY_OPS ) {
2173
2174 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
2175
2176 SECTION( LABEL_CORRECTNESS ) {
2177
2178 int numTargs = GENERATE_COPY( range(1,numQubits+1) );
2179 auto targets = GENERATE_TARGS( numQubits, numTargs );
2180 auto outcomes = getRandomOutcomes(numTargs);
2181
2182 qmatrix projector = getProjector(targets, outcomes, numQubits);
2183
2184 auto testFunc = [&](Qureg qureg, auto& ref) {
2185 rightapplyMultiQubitProjector(qureg, targets.data(), outcomes.data(), numTargs);
2186 rightapplyReferenceOperator(ref, projector);
2187 };
2188
2189 CAPTURE( targets, outcomes );
2190 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
2191 }
2192
2193 /// @todo input validation
2194}
2195
2196
2197TEST_CASE( "leftapplyPauliStrSum", TEST_CATEGORY_MULT LABEL_MIXED_DEPLOY_TAG ) {
2198
2199 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
2200
2201 SECTION( LABEL_CORRECTNESS ) {
2202
2203 int numQubits = getNumCachedQubits();
2204 int numTerms = GENERATE_COPY( 1, 2, 10 );
2205
2206 PauliStrSum sum = createRandomPauliStrSum(numQubits, numTerms);
2207
2208 auto testFunc = [&](Qureg qureg, auto& ref) {
2209
2210 // must use (and ergo make) an identically-deployed workspace
2211 Qureg workspace = createCloneQureg(qureg);
2212 leftapplyPauliStrSum(qureg, sum, workspace);
2213 destroyQureg(workspace);
2214
2215 ref = getMatrix(sum, numQubits) * ref;
2216 };
2217
2218 CAPTURE( numTerms );
2219 SECTION( LABEL_STATEVEC ) { TEST_ON_CACHED_QUREGS(statevecQuregs, statevecRef, testFunc); }
2220 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
2221 }
2222
2223 /// @todo input validation
2224}
2225
2226
2227TEST_CASE( "rightapplyPauliStrSum", TEST_CATEGORY_MULT LABEL_MIXED_DEPLOY_TAG ) {
2228
2229 PREPARE_TEST( numQubits, statevecQuregs, densmatrQuregs, statevecRef, densmatrRef );
2230
2231 SECTION( LABEL_CORRECTNESS ) {
2232
2233 int numQubits = getNumCachedQubits();
2234 int numTerms = GENERATE_COPY( 1, 2, 10 );
2235
2236 PauliStrSum sum = createRandomPauliStrSum(numQubits, numTerms);
2237
2238 auto testFunc = [&](Qureg qureg, auto& ref) {
2239
2240 // must use (and ergo make) an identically-deployed workspace
2241 Qureg workspace = createCloneQureg(qureg);
2242 rightapplyPauliStrSum(qureg, sum, workspace);
2243 destroyQureg(workspace);
2244
2245 ref = ref * getMatrix(sum, numQubits);
2246 };
2247
2248 CAPTURE( numTerms );
2249 SECTION( LABEL_DENSMATR ) { TEST_ON_CACHED_QUREGS(densmatrQuregs, densmatrRef, testFunc); }
2250 }
2251
2252 /// @todo input validation
2253}
2254
2255
2256/** @} (end defgroup) */
void setValidationEpsilonToDefault()
Definition debug.cpp:108
void setValidationEpsilon(qreal eps)
Definition debug.cpp:100
QuESTEnv getQuESTEnv()
void initDebugState(Qureg qureg)
CompMatr createCompMatr(int numQubits)
Definition matrices.cpp:213
DiagMatr createDiagMatr(int numQubits)
Definition matrices.cpp:248
void destroyDiagMatr(DiagMatr matrix)
Definition matrices.cpp:399
void destroyCompMatr(CompMatr matrix)
Definition matrices.cpp:398
static CompMatr2 getCompMatr2(qcomp **in)
Definition matrices.h:351
static CompMatr1 getCompMatr1(qcomp **in)
Definition matrices.h:325
static DiagMatr2 getDiagMatr2(qcomp *in)
Definition matrices.h:403
static DiagMatr1 getDiagMatr1(qcomp *in)
Definition matrices.h:378
void setDiagMatr(DiagMatr out, qcomp *in)
Definition matrices.cpp:431
void setCompMatr(CompMatr matr, qcomp **vals)
Definition matrices.cpp:423
void setFullStateDiagMatr(FullStateDiagMatr out, qindex startInd, qcomp *in, qindex numElems)
Definition matrices.cpp:443
void rightapplyCompMatr1(Qureg qureg, int target, CompMatr1 matrix)
void leftapplyCompMatr1(Qureg qureg, int target, CompMatr1 matrix)
void rightapplyCompMatr2(Qureg qureg, int target1, int target2, CompMatr2 matrix)
void leftapplyCompMatr2(Qureg qureg, int target1, int target2, CompMatr2 matr)
void leftapplyCompMatr(Qureg qureg, int *targets, int numTargets, CompMatr matrix)
void rightapplyCompMatr(Qureg qureg, int *targets, int numTargets, CompMatr matrix)
void leftapplyDiagMatr1(Qureg qureg, int target, DiagMatr1 matr)
void rightapplyDiagMatr1(Qureg qureg, int target, DiagMatr1 matrix)
void leftapplyDiagMatr2(Qureg qureg, int target1, int target2, DiagMatr2 matr)
void rightapplyDiagMatr2(Qureg qureg, int target1, int target2, DiagMatr2 matrix)
void rightapplyDiagMatr(Qureg qureg, int *targets, int numTargets, DiagMatr matrix)
void leftapplyDiagMatrPower(Qureg qureg, int *targets, int numTargets, DiagMatr matrix, qcomp exponent)
void rightapplyDiagMatrPower(Qureg qureg, int *targets, int numTargets, DiagMatr matrix, qcomp exponent)
void leftapplyDiagMatr(Qureg qureg, int *targets, int numTargets, DiagMatr matrix)
void rightapplyFullStateDiagMatr(Qureg qureg, FullStateDiagMatr matrix)
void leftapplyFullStateDiagMatr(Qureg qureg, FullStateDiagMatr matrix)
void rightapplyFullStateDiagMatrPower(Qureg qureg, FullStateDiagMatr matrix, qcomp exponent)
void leftapplyFullStateDiagMatrPower(Qureg qureg, FullStateDiagMatr matrix, qcomp exponent)
void leftapplyMultiQubitNot(Qureg qureg, int *targets, int numTargets)
void rightapplyMultiQubitNot(Qureg qureg, int *targets, int numTargets)
void leftapplyPauliX(Qureg qureg, int target)
void rightapplyPauliY(Qureg qureg, int target)
void leftapplyPauliY(Qureg qureg, int target)
void leftapplyPauliZ(Qureg qureg, int target)
void rightapplyPauliX(Qureg qureg, int target)
void rightapplyPauliZ(Qureg qureg, int target)
void leftapplyPauliGadget(Qureg qureg, PauliStr str, qreal angle)
void rightapplyPauliGadget(Qureg qureg, PauliStr str, qreal angle)
void rightapplyPauliStr(Qureg qureg, PauliStr str)
void leftapplyPauliStr(Qureg qureg, PauliStr str)
void rightapplyPauliStrSum(Qureg qureg, PauliStrSum sum, Qureg workspace)
void leftapplyPauliStrSum(Qureg qureg, PauliStrSum sum, Qureg workspace)
void rightapplyPhaseGadget(Qureg qureg, int *targets, int numTargets, qreal angle)
void leftapplyPhaseGadget(Qureg qureg, int *targets, int numTargets, qreal angle)
void leftapplyMultiQubitProjector(Qureg qureg, int *qubits, int *outcomes, int numQubits)
void rightapplyQubitProjector(Qureg qureg, int qubit, int outcome)
void rightapplyMultiQubitProjector(Qureg qureg, int *qubits, int *outcomes, int numQubits)
void leftapplyQubitProjector(Qureg qureg, int qubit, int outcome)
void leftapplySwap(Qureg qureg, int qubit1, int qubit2)
void rightapplySwap(Qureg qureg, int qubit1, int qubit2)
void applyFullStateDiagMatr(Qureg qureg, FullStateDiagMatr matrix)
void applyFullStateDiagMatrPower(Qureg qureg, FullStateDiagMatr matrix, qcomp exponent)
qreal applyForcedQubitMeasurement(Qureg qureg, int target, int outcome)
qindex applyMultiQubitMeasurement(Qureg qureg, int *qubits, int numQubits)
int applyQubitMeasurement(Qureg qureg, int target)
qreal applyForcedMultiQubitMeasurement(Qureg qureg, int *qubits, int *outcomes, int numQubits)
qindex applyMultiQubitMeasurementAndGetProb(Qureg qureg, int *qubits, int numQubits, qreal *probability)
int applyQubitMeasurementAndGetProb(Qureg qureg, int target, qreal *probability)
void applyNonUnitaryPauliGadget(Qureg qureg, PauliStr str, qcomp angle)
void applyMultiQubitPhaseShift(Qureg qureg, int *targets, int numTargets, qreal angle)
void applyTwoQubitPhaseShift(Qureg qureg, int target1, int target2, qreal angle)
void applyTwoQubitPhaseFlip(Qureg qureg, int target1, int target2)
void applyPhaseShift(Qureg qureg, int target, qreal angle)
void applyPhaseFlip(Qureg qureg, int target)
void applyMultiQubitPhaseFlip(Qureg qureg, int *targets, int numTargets)
void applyMultiQubitProjector(Qureg qureg, int *qubits, int *outcomes, int numQubits)
void applyQubitProjector(Qureg qureg, int target, int outcome)
void applyFullQuantumFourierTransform(Qureg qureg)
void applyQuantumFourierTransform(Qureg qureg, int *targets, int numTargets)
PauliStr getPauliStr(const char *paulis, int *indices, int numPaulis)
Definition paulis.cpp:76
Qureg createCloneQureg(Qureg qureg)
Definition qureg.cpp:319
void destroyQureg(Qureg qureg)
Definition qureg.cpp:334
qmatrix getKroneckerProduct(qmatrix a, qmatrix b)
Definition linalg.cpp:523
qmatrix getIdentityMatrix(size_t dim)
Definition qmatrix.cpp:30
qmatrix getZeroMatrix(size_t dim)
Definition qmatrix.cpp:18
qcomp getRandomComplex()
Definition random.cpp:107
qmatrix getRandomUnitary(int numQb)
Definition random.cpp:348
qreal getRandomReal(qreal min, qreal maxExcl)
Definition random.cpp:63
vector< qreal > getRandomProbabilities(int numProbs)
Definition random.cpp:160
int getRandomInt(int min, int maxExcl)
Definition random.cpp:90
TEST_CASE("calcExpecPauliStr", TEST_CATEGORY)
TEST_ALL_CTRL_OPERATIONS(PauliStr, any, paulistr, nullptr)
Definition qureg.h:49