GCC Code Coverage Report

Directory: ./
File: Circuit/AssertionSynthesis.cpp
Date: 2022-10-15 05:10:18
Warnings: 1 unchecked decisions!
Exec Total Coverage
Lines: 155 162 95.7%
Functions: 10 10 100.0%
Branches: 216 389 55.5%
Decisions: 57 62 91.9%
Line Branch Decision Exec Source
1 // Copyright 2019-2022 Cambridge Quantum Computing
2 //
3 // Licensed under the Apache License, Version 2.0 (the "License");
4 // you may not use this file except in compliance with the License.
5 // You may obtain a copy of the License at
6 //
7 // http://www.apache.org/licenses/LICENSE-2.0
8 //
9 // Unless required by applicable law or agreed to in writing, software
10 // distributed under the License is distributed on an "AS IS" BASIS,
11 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12 // See the License for the specific language governing permissions and
13 // limitations under the License.
14
15 #include "AssertionSynthesis.hpp"
16
17 #include <array>
18 #include <cmath>
19 #include <complex>
20 #include <optional>
21 #include <stdexcept>
22 #include <tkassert/Assert.hpp>
23
24 #include "CircUtils.hpp"
25 #include "Circuit.hpp"
26 #include "OpType/OpType.hpp"
27 #include "Utils/Constants.hpp"
28 #include "Utils/CosSinDecomposition.hpp"
29 #include "Utils/EigenConfig.hpp"
30 #include "Utils/MatrixAnalysis.hpp"
31 #include "Utils/UnitID.hpp"
32
33 namespace tket {
34
35 /**
36 * Diagonalise a projector matrix P such that P = U*D*U.dag, where
37 * U is a unitary and D is a diagonal binary matrix.
38 * The resulting diagonal matrix always have 1s on the left.
39 * The benefit of permuting 1s to the left is that the
40 * diagonal matrix can always be factorised into |0><0|, |1><1| or I if the
41 * projector has a valid rank. However, this is not optimised.
42 *
43 * @param P projector matrix
44 *
45 * @return D as a boolean vector, U, rank
46 */
47 9 static std::tuple<VectorXb, Eigen::MatrixXcd, int> projector_diagonalisation(
48 const Eigen::MatrixXcd &P) {
49 // Solve the eigen values
50
1/2
✓ Branch 1 taken 9 times.
✗ Branch 2 not taken.
 9 Eigen::SelfAdjointEigenSolver<Eigen::MatrixXcd> eigen_solver(P);
51 // Make a permutation matrix to reorder D and U_dag
52
1/2
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
 9 Eigen::PermutationMatrix<Eigen::Dynamic> perm(P.rows());
53 9 Eigen::MatrixXcd P1 = eigen_solver.eigenvectors() *
54
3/6
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 9 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 9 times.
✗ Branch 9 not taken.
 18 eigen_solver.eigenvalues().asDiagonal() *
55
2/4
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 9 times.
✗ Branch 6 not taken.
 18 eigen_solver.eigenvectors().adjoint();
56 TKET_ASSERT(P1.isApprox(P));
57 // Cast eigenvalues to booleans
58
1/2
✓ Branch 3 taken 9 times.
✗ Branch 4 not taken.
 9 VectorXb eigenvalues(eigen_solver.eigenvalues().rows());
59
2/2
✓ Branch 1 taken 54 times.
✓ Branch 2 taken 9 times.
2/2
✓ Decision 'true' taken 54 times.
✓ Decision 'false' taken 9 times.
 63 for (unsigned i = 0; i < eigenvalues.rows(); i++) {
60
3/4
✓ Branch 2 taken 54 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 30 times.
✓ Branch 6 taken 24 times.
2/2
✓ Decision 'true' taken 30 times.
✓ Decision 'false' taken 24 times.
 54 if (std::abs(eigen_solver.eigenvalues()[i]) < EPS) {
61
1/2
✓ Branch 1 taken 30 times.
✗ Branch 2 not taken.
 30 eigenvalues[i] = 0;
62 } else {
63
1/2
✓ Branch 1 taken 24 times.
✗ Branch 2 not taken.
 24 eigenvalues[i] = 1;
64 }
65 }
66 // Move 1 to the front
67 9 unsigned j = 0;
68
2/2
✓ Branch 1 taken 54 times.
✓ Branch 2 taken 9 times.
2/2
✓ Decision 'true' taken 54 times.
✓ Decision 'false' taken 9 times.
 63 for (unsigned i = 0; i < P.rows(); i++) {
69
3/4
✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 24 times.
✓ Branch 4 taken 30 times.
2/2
✓ Decision 'true' taken 24 times.
✓ Decision 'false' taken 30 times.
 54 if (eigenvalues[i] == 1) {
70
1/2
✓ Branch 2 taken 24 times.
✗ Branch 3 not taken.
 24 perm.indices()[j] = i;
71 24 j++;
72 }
73 }
74 // Assign rank
75 9 unsigned rank = j;
76 // Pad 0 to the end
77
2/2
✓ Branch 1 taken 54 times.
✓ Branch 2 taken 9 times.
2/2
✓ Decision 'true' taken 54 times.
✓ Decision 'false' taken 9 times.
 63 for (unsigned i = 0; i < P.rows(); i++) {
78
3/4
✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 30 times.
✓ Branch 4 taken 24 times.
2/2
✓ Decision 'true' taken 30 times.
✓ Decision 'false' taken 24 times.
 54 if (eigenvalues[i] == 0) {
79
1/2
✓ Branch 2 taken 30 times.
✗ Branch 3 not taken.
 30 perm.indices()[j] = i;
80 30 j++;
81 }
82 }
83
3/6
✓ Branch 1 taken 9 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 9 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 9 times.
✗ Branch 8 not taken.
 9 VectorXb D = perm.transpose() * eigenvalues;
84
2/4
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 9 times.
✗ Branch 6 not taken.
 9 Eigen::MatrixXcd U = eigen_solver.eigenvectors() * perm;
85 TKET_ASSERT((U * D.cast<double>().asDiagonal() * U.adjoint()).isApprox(P));
86
1/2
✓ Branch 1 taken 9 times.
✗ Branch 2 not taken.
 18 return {D, U, rank};
87 9 }
88
89 // Get the n_th bit of an integer
90 89 static bool get_bit(unsigned value, unsigned index) {
91 89 return (value & (1 << index)) >> index;
92 }
93
94 /**
95 * Tensor-factorize a binary diagonal matrix with dimension 2^n.
96 * Assume D has a rank <= 2^(n-1), and the rank is power of 2
97 *
98 * @param D a vector representation of a diagonal binary matrix
99 *
100 * @return A list of 2x2 diagonal binary matrices D_0,...,D_m (as vectors)
101 * such that tensor(D_0,...,D_m) == D
102 */
103 10 static std::vector<Vector2b> tensor_factorization(const VectorXb &D) {
104 10 unsigned dim = D.rows();
105 10 unsigned log_dim = (unsigned)log2(dim);
106 10 std::vector<Vector2b> factorisations;
107
2/2
✓ Branch 0 taken 26 times.
✓ Branch 1 taken 10 times.
2/2
✓ Decision 'true' taken 26 times.
✓ Decision 'false' taken 10 times.
 36 for (unsigned i = 0; i < log_dim; i++) {
108
3/6
✓ Branch 1 taken 26 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 26 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 26 times.
✗ Branch 8 not taken.
 26 factorisations.push_back(Vector2b::Zero());
109 }
110 // The diagonal matrix D can be factorised into a tensorproduct of
111 // log_dim 2x2 diagonal boolean matrices in the set {\0><0|, |1><1|, I}
112 // The factorisation can be determined by the diagonals in D
113 // For example, if D[0,0] == 1, then all the 2x2 matices should have 1 in
114 // their top left entry.
115
2/2
✓ Branch 0 taken 66 times.
✓ Branch 1 taken 10 times.
2/2
✓ Decision 'true' taken 66 times.
✓ Decision 'false' taken 10 times.
 76 for (unsigned i = 0; i < dim; i++) {
116
3/4
✓ Branch 1 taken 66 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 31 times.
✓ Branch 4 taken 35 times.
0/1
? Decision couldn't be analyzed.
 66 if (D(i) == 1) {
117
2/2
✓ Branch 0 taken 89 times.
✓ Branch 1 taken 31 times.
2/2
✓ Decision 'true' taken 89 times.
✓ Decision 'false' taken 31 times.
 120 for (unsigned j = 0; j < log_dim; j++) {
118
2/2
✓ Branch 1 taken 28 times.
✓ Branch 2 taken 61 times.
2/2
✓ Decision 'true' taken 28 times.
✓ Decision 'false' taken 61 times.
 89 if (get_bit(i, j)) {
119
1/2
✓ Branch 2 taken 28 times.
✗ Branch 3 not taken.
 28 factorisations[j](1) = 1;
120 } else {
121
1/2
✓ Branch 2 taken 61 times.
✗ Branch 3 not taken.
 61 factorisations[j](0) = 1;
122 }
123 }
124 }
125 }
126 10 return factorisations;
127 }
128
129 20 static void apply_unitary(Circuit *circ, const Eigen::MatrixXcd &U) {
130
2/2
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 18 times.
2/2
✓ Decision 'true' taken 2 times.
✓ Decision 'false' taken 18 times.
 20 if (U.rows() == 2) {
131 // Add a Unitary1qBox
132
2/4
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
 2 Unitary1qBox ubox(U);
133
3/6
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
 2 circ->add_op<unsigned>(std::make_shared<Unitary1qBox>(ubox), {0});
134
2/2
✓ Branch 2 taken 4 times.
✓ Branch 3 taken 14 times.
2/2
✓ Decision 'true' taken 4 times.
✓ Decision 'false' taken 16 times.
 20 } else if (U.rows() == 4) {
135 // Add a Unitary2qBox
136
2/4
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4 times.
✗ Branch 5 not taken.
 4 Unitary2qBox ubox(U);
137
3/6
✓ Branch 3 taken 4 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 4 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 4 times.
✗ Branch 11 not taken.
 4 circ->add_op<unsigned>(std::make_shared<Unitary2qBox>(ubox), {0, 1});
138
1/2
✓ Branch 2 taken 14 times.
✗ Branch 3 not taken.
2/2
✓ Decision 'true' taken 14 times.
✓ Decision 'false' taken 4 times.
 18 } else if (U.rows() == 8) {
139 // Add a Unitary3qBox
140
2/4
✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
 14 Unitary3qBox ubox(U);
141
3/6
✓ Branch 3 taken 14 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 14 times.
✗ Branch 7 not taken.
✓ Branch 10 taken 14 times.
✗ Branch 11 not taken.
 14 circ->add_op<unsigned>(std::make_shared<Unitary3qBox>(ubox), {0, 1, 2});
142
143 14 } else {
144 ✗ throw CircuitInvalidity("Only 2x2, 4x4, and 8x8 projectors are supported");
145 }
146 20 }
147
148 /**
149 * Apply measurements and update the expectations
150 *
151 * @param circ the circuit to apply measurements to
152 * @param z_projectors a list of 2x2 diagonal binary matrices
153 * @param debug_bit_index the next available index of the default classical
154 * register
155 * @param expected_readouts a list of expectations. e.g.
156 * expected_readouts[i]=True indicates that the ith bit in the default register
157 * expects the readout to be 1.
158 *
159 * @return the next available index of the default classical register
160 */
161 10 static unsigned apply_z_measurements(
162 Circuit &circ, const std::vector<Vector2b> z_projectors,
163 unsigned debug_bit_index, std::vector<bool> &expected_readouts) {
164
2/2
✓ Branch 1 taken 26 times.
✓ Branch 2 taken 10 times.
2/2
✓ Decision 'true' taken 26 times.
✓ Decision 'false' taken 10 times.
 36 for (unsigned i = 0; i < z_projectors.size(); i++) {
165
7/10
✓ Branch 2 taken 26 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 26 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 26 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 14 times.
✓ Branch 11 taken 12 times.
✓ Branch 12 taken 14 times.
✓ Branch 13 taken 12 times.
2/2
✓ Decision 'true' taken 14 times.
✓ Decision 'false' taken 12 times.
 26 if (z_projectors[i](0) == 1 && z_projectors[i](1) == 1) {
166 14 continue;
167 }
168
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 Qubit q(i);
169
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 Bit b(debug_bit_index++);
170
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 circ.add_bit(b);
171
4/8
✓ Branch 5 taken 12 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 12 times.
✗ Branch 9 not taken.
✓ Branch 12 taken 24 times.
✓ Branch 13 taken 12 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
 36 circ.add_op<UnitID>(OpType::Measure, {q, b});
172 const bool is_1_projector =
173
2/8
✓ Branch 2 taken 12 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 12 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 12 z_projectors[i](0) == 0 && z_projectors[i](1) == 1;
174
1/2
✓ Branch 1 taken 12 times.
✗ Branch 2 not taken.
 12 expected_readouts.push_back(is_1_projector);
175 12 }
176 10 return debug_bit_index;
177 }
178
179 /**
180 * Split a projector into two 2^(n-1)-ranked projectors that project the same
181 * subspace when combined. Assume D has 1s on the left. Suppose D.rows() = 16,
182 * rank = 5, then the diagonal of D is 11111000|00000000. D can't be factorised
183 * into Z basis projectors, so we create two rank=2^(n-1) projectors that
184 * project the same subspace when combined (intersection of two subspaces).
185 * e.g. D1: 11111111|00000000, D2: 11111000|11100000.
186 * After permuting D2, both diagonal matrices should have diagonal entires:
187 * 11111111|00000000.
188 * Similar to the projector_diagonalisation function, the permutation of the
189 * eigen decomposition is not optimised.
190 *
191 * @param D a vector representation of a diagonal binary matrix
192 * @param U a unitary matrix
193 * @param rank the rank of D
194 *
195 * @return D1, U1, D2, U2 such that the intersection of the subspaces projected
196 * by P1=U1*D1*U1.dag and P2=U2*D2*U2.dag is the subspace projected by
197 * P=U*D*U.dag.
198 */
199 static std::tuple<VectorXb, Eigen::MatrixXcd, VectorXb, Eigen::MatrixXcd>
200 3 projector_split(
201 const VectorXb &D, const Eigen::MatrixXcd &U, const unsigned &rank) {
202
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 VectorXb D_new(D);
203
2/2
✓ Branch 1 taken 3 times.
✓ Branch 2 taken 3 times.
2/2
✓ Decision 'true' taken 3 times.
✓ Decision 'false' taken 3 times.
 6 for (unsigned i = rank; i < D_new.rows() / 2; i++) {
204
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 D_new(i) = 1;
205 }
206 // U1
207
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 Eigen::MatrixXcd U1(U);
208 // U2
209
1/2
✓ Branch 2 taken 3 times.
✗ Branch 3 not taken.
 3 Eigen::PermutationMatrix<Eigen::Dynamic> perm(D_new.rows());
210
2/2
✓ Branch 1 taken 24 times.
✓ Branch 2 taken 3 times.
2/2
✓ Decision 'true' taken 24 times.
✓ Decision 'false' taken 3 times.
 27 for (unsigned i = 0; i < D_new.rows(); i++) {
211
6/6
✓ Branch 0 taken 15 times.
✓ Branch 1 taken 9 times.
✓ Branch 3 taken 3 times.
✓ Branch 4 taken 12 times.
✓ Branch 5 taken 3 times.
✓ Branch 6 taken 21 times.
2/2
✓ Decision 'true' taken 3 times.
✓ Decision 'false' taken 21 times.
 24 if (i >= rank && i < D_new.rows() / 2) {
212
1/2
✓ Branch 3 taken 3 times.
✗ Branch 4 not taken.
 3 perm.indices()[i] = i + D_new.rows() / 2 - rank;
213
6/6
✓ Branch 1 taken 12 times.
✓ Branch 2 taken 9 times.
✓ Branch 4 taken 3 times.
✓ Branch 5 taken 9 times.
✓ Branch 6 taken 3 times.
✓ Branch 7 taken 18 times.
2/2
✓ Decision 'true' taken 3 times.
✓ Decision 'false' taken 18 times.
 21 } else if (i >= D_new.rows() / 2 && i < D_new.rows() - rank) {
214
1/2
✓ Branch 3 taken 3 times.
✗ Branch 4 not taken.
 3 perm.indices()[i] = i - (D_new.rows() / 2 - rank);
215 } else {
216
1/2
✓ Branch 2 taken 18 times.
✗ Branch 3 not taken.
 18 perm.indices()[i] = i;
217 }
218 }
219
2/4
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
 3 Eigen::MatrixXcd U2 = U * perm;
220
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 6 return {D_new, U1, D_new, U2};
221 3 }
222
223 8 std::tuple<Circuit, std::vector<bool>> projector_assertion_synthesis(
224 const Eigen::MatrixXcd &P) {
225 // Diagonoalise the projector P
226
1/2
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
 8 Eigen::MatrixXcd projector(P);
227 8 unsigned n_qubits = (unsigned)log2(projector.rows());
228
1/2
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
 8 VectorXb D;
229
1/2
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
 8 Eigen::MatrixXcd U;
230 unsigned rank;
231
1/2
✓ Branch 1 taken 8 times.
✗ Branch 2 not taken.
 8 std::tie(D, U, rank) = projector_diagonalisation(projector);
232
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 8 times.
1/2
✗ Decision 'true' not taken.
✓ Decision 'false' taken 8 times.
 8 if (rank == 0) {
233 ✗ throw CircuitInvalidity("The projector must have none-zero rank");
234 }
235 8 std::vector<bool> expected_readouts;
236
1/2
✓ Branch 2 taken 8 times.
✗ Branch 3 not taken.
 8 Circuit circ(n_qubits);
237
2/2
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 6 times.
2/2
✓ Decision 'true' taken 2 times.
✓ Decision 'false' taken 6 times.
 8 if (rank > projector.rows() / 2) {
238
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
2/2
✓ Decision 'true' taken 1 times.
✓ Decision 'false' taken 1 times.
 2 if (n_qubits >= 3) {
239
2/4
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 1 times.
✗ Branch 6 not taken.
 1 throw CircuitInvalidity(
240 2 "8x8 projector that requires an ancilla is not supported");
241 }
242 // Add an auxilary qubit to the end
243
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 1 Qubit ancilla(n_qubits);
244
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 1 circ.add_qubit(ancilla, false);
245
4/8
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1 times.
✗ Branch 8 not taken.
✓ Branch 11 taken 1 times.
✓ Branch 12 taken 1 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 2 circ.add_op<Qubit>(OpType::Reset, {ancilla});
246
2/4
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
 1 Eigen::MatrixXcd zero_projector = Eigen::MatrixXcd::Zero(2, 2);
247
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 1 zero_projector(0, 0) = 1;
248 Eigen::MatrixXcd new_projector =
249
2/4
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
 1 Eigen::KroneckerProduct(projector, zero_projector);
250
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 1 std::tie(D, U, rank) = projector_diagonalisation(new_projector);
251 1 }
252
253 // Initialise debug bit index
254 7 unsigned debug_bit_index = 0;
255 // Check if the rank is power of 2
256 // If a number is power of 2 then the count of binary 1s will be 1.
257
3/4
✓ Branch 0 taken 7 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 4 times.
✓ Branch 3 taken 3 times.
2/2
✓ Decision 'true' taken 4 times.
✓ Decision 'false' taken 3 times.
 7 if (rank != 0 && (rank & (rank - 1)) == 0) {
258 // Implement the projection
259 // Apply the transformation before the measurements
260
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 std::vector<Vector2b> tensors = tensor_factorization(D);
261
3/6
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 4 times.
✗ Branch 8 not taken.
 4 apply_unitary(&circ, U.adjoint());
262
2/4
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4 times.
✗ Branch 5 not taken.
 4 apply_z_measurements(circ, tensors, debug_bit_index, expected_readouts);
263
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
 4 apply_unitary(&circ, U);
264 4 }
265 // Transform this into two projectors
266 else {
267
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 auto [D1, U1, D2, U2] = projector_split(D, U, rank);
268
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 std::vector<Vector2b> tensors1 = tensor_factorization(D1);
269
3/6
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 3 times.
✗ Branch 8 not taken.
 3 apply_unitary(&circ, U1.adjoint());
270
2/4
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
 3 debug_bit_index = apply_z_measurements(
271 circ, tensors1, debug_bit_index, expected_readouts);
272
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 apply_unitary(&circ, U1);
273
274
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 std::vector<Vector2b> tensors2 = tensor_factorization(D2);
275
3/6
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 3 times.
✗ Branch 8 not taken.
 3 apply_unitary(&circ, U2.adjoint());
276
2/4
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
 3 apply_z_measurements(circ, tensors2, debug_bit_index, expected_readouts);
277
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 apply_unitary(&circ, U2);
278 3 }
279
1/2
✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
 14 return {circ, expected_readouts};
280 12 }
281
282 5 static unsigned get_n_qubits_from_stabilisers(
283 const PauliStabiliserList &paulis) {
284
2/2
✓ Branch 1 taken 1 times.
✓ Branch 2 taken 4 times.
2/2
✓ Decision 'true' taken 1 times.
✓ Decision 'false' taken 4 times.
 5 if (paulis.size() == 0) {
285
2/4
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
 1 throw CircuitInvalidity("Stabilisers cannot be empty");
286 }
287 4 unsigned stabiliser_len = paulis[0].string.size();
288
2/2
✓ Branch 1 taken 5 times.
✓ Branch 2 taken 3 times.
2/2
✓ Decision 'true' taken 5 times.
✓ Decision 'false' taken 3 times.
 8 for (unsigned i = 1; i < paulis.size(); i++) {
289
2/2
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
2/2
✓ Decision 'true' taken 1 times.
✓ Decision 'false' taken 4 times.
 5 if (paulis[i].string.size() != stabiliser_len) {
290
2/4
✓ Branch 3 taken 1 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 1 times.
✗ Branch 7 not taken.
 1 throw CircuitInvalidity("Stabilisers have unequal lengths");
291 }
292 }
293 3 return stabiliser_len;
294 }
295
296 /**
297 * Apply assertion operators and update the expectations
298 *
299 * @param circ the circuit to apply assertion operators to
300 * @param pauli a Pauli stabiliser
301 * @param debug_bit_index the next available index of the default classical
302 * register
303 * @param expected_readouts a list of expectations. e.g.
304 * expected_readouts[i]=True indicates that the ith bit in the default register
305 * expects the readout to be 1.
306 *
307 * @return the next available index of the default classical register
308 */
309 7 static unsigned add_assertion_operator(
310 Circuit &circ, const PauliStabiliser &pauli, unsigned debug_bit_index,
311 const Qubit &ancilla, std::vector<bool> &expected_readouts) {
312
4/8
✓ Branch 4 taken 7 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 7 times.
✗ Branch 8 not taken.
✓ Branch 11 taken 7 times.
✓ Branch 12 taken 7 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 14 circ.add_op<Qubit>(OpType::Reset, {ancilla});
313
4/8
✓ Branch 4 taken 7 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 7 times.
✗ Branch 8 not taken.
✓ Branch 11 taken 7 times.
✓ Branch 12 taken 7 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 14 circ.add_op<Qubit>(OpType::H, {ancilla});
314
2/2
✓ Branch 1 taken 14 times.
✓ Branch 2 taken 7 times.
2/2
✓ Decision 'true' taken 14 times.
✓ Decision 'false' taken 7 times.
 21 for (unsigned i = 0; i < pauli.string.size(); i++) {
315
2/5
✗ Branch 1 not taken.
✓ Branch 2 taken 6 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 8 times.
✗ Branch 5 not taken.
 14 switch (pauli.string[i]) {
316
0/1
✗ Decision 'true' not taken.
 ✗ case Pauli::I:
317 ✗ break;
318
1/1
✓ Decision 'true' taken 18 times.
 6 case Pauli::X:
319
5/10
✓ Branch 3 taken 6 times.
✗ Branch 4 not taken.
✓ Branch 7 taken 6 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 6 times.
✗ Branch 11 not taken.
✓ Branch 14 taken 12 times.
✓ Branch 15 taken 6 times.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
 18 circ.add_op<Qubit>(OpType::CX, {ancilla, Qubit(i)});
320 6 break;
321
0/1
✗ Decision 'true' not taken.
 ✗ case Pauli::Y:
322 ✗ circ.add_op<Qubit>(OpType::CY, {ancilla, Qubit(i)});
323 ✗ break;
324
1/1
✓ Decision 'true' taken 24 times.
 8 case Pauli::Z:
325
5/10
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✓ Branch 7 taken 8 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 8 times.
✗ Branch 11 not taken.
✓ Branch 14 taken 16 times.
✓ Branch 15 taken 8 times.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
 24 circ.add_op<Qubit>(OpType::CZ, {ancilla, Qubit(i)});
326 8 break;
327 }
328 }
329
4/8
✓ Branch 4 taken 7 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 7 times.
✗ Branch 8 not taken.
✓ Branch 11 taken 7 times.
✓ Branch 12 taken 7 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 14 circ.add_op<Qubit>(OpType::H, {ancilla});
330
1/2
✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
 7 Bit b(debug_bit_index++);
331
1/2
✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
 7 circ.add_bit(b);
332
4/8
✓ Branch 5 taken 7 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 7 times.
✗ Branch 9 not taken.
✓ Branch 12 taken 14 times.
✓ Branch 13 taken 7 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
 21 circ.add_op<UnitID>(OpType::Measure, {ancilla, b});
333
1/2
✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
 7 expected_readouts.push_back(!pauli.coeff);
334 14 return debug_bit_index;
335 }
336
337 // Assume all Pauli stabilisers have equal lengths and no identity
338 5 std::tuple<Circuit, std::vector<bool>> stabiliser_assertion_synthesis(
339 const PauliStabiliserList &paulis) {
340
2/2
✓ Branch 1 taken 3 times.
✓ Branch 2 taken 2 times.
 5 unsigned n_qubits = get_n_qubits_from_stabilisers(paulis);
341 3 std::vector<bool> expected_readouts;
342
1/2
✓ Branch 2 taken 3 times.
✗ Branch 3 not taken.
 3 Circuit circ(n_qubits);
343
344 // Initialise debug bit index
345 3 unsigned debug_bit_index = 0;
346 // Add an ancilla
347
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 Qubit ancilla(n_qubits);
348
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 circ.add_qubit(ancilla, false);
349
350
2/2
✓ Branch 4 taken 7 times.
✓ Branch 5 taken 3 times.
2/2
✓ Decision 'true' taken 7 times.
✓ Decision 'false' taken 3 times.
 10 for (auto &pauli : paulis) {
351
1/2
✓ Branch 1 taken 7 times.
✗ Branch 2 not taken.
 7 debug_bit_index = add_assertion_operator(
352 circ, pauli, debug_bit_index, ancilla, expected_readouts);
353 }
354
355
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 6 return {circ, expected_readouts};
356 3 }
357
358 } // namespace tket
359